[GRASS-SVN] r56542 - in sandbox/martinl: . v.ten

svn_grass at osgeo.org svn_grass at osgeo.org
Sat Jun 1 15:28:53 PDT 2013


Author: martinl
Date: 2013-06-01 15:28:53 -0700 (Sat, 01 Jun 2013)
New Revision: 56542

Added:
   sandbox/martinl/v.ten/
   sandbox/martinl/v.ten/Makefile
   sandbox/martinl/v.ten/local_proto.h
   sandbox/martinl/v.ten/main.cpp
   sandbox/martinl/v.ten/read.cpp
   sandbox/martinl/v.ten/v.ten.html
   sandbox/martinl/v.ten/write.cpp
Log:
sandbox: experimental module for 3D triangulation based on CGAL

Added: sandbox/martinl/v.ten/Makefile
===================================================================
--- sandbox/martinl/v.ten/Makefile	                        (rev 0)
+++ sandbox/martinl/v.ten/Makefile	2013-06-01 22:28:53 UTC (rev 56542)
@@ -0,0 +1,18 @@
+MODULE_TOPDIR = ../..
+
+PGM = v.ten
+
+LIBES = $(VECTORLIB) $(GISLIB)
+DEPENDENCIES = $(VECTORDEP) $(GISDEP)
+EXTRA_INC = $(VECT_INC)
+EXTRA_CFLAGS = $(VECT_CFLAGS)
+
+include $(MODULE_TOPDIR)/include/Make/Module.make
+
+EXTRA_CFLAGS = -frounding-math
+EXTRA_LDFLAGS = -lCGAL -lgmp
+LINK = $(CXX)
+
+ifneq ($(strip $(CXX)),)
+default: cmd
+endif

Added: sandbox/martinl/v.ten/local_proto.h
===================================================================
--- sandbox/martinl/v.ten/local_proto.h	                        (rev 0)
+++ sandbox/martinl/v.ten/local_proto.h	2013-06-01 22:28:53 UTC (rev 56542)
@@ -0,0 +1,16 @@
+#ifndef __LOCAL_PROTO_H__
+
+#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
+#include <CGAL/Triangulation_3.h>
+
+typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
+
+typedef CGAL::Triangulation_3<K>      Triangulation;
+typedef Triangulation::Point          Point;
+
+/* read.cpp */
+int read_points(struct Map_info *, int, std::vector<Point>&);
+
+/* write.cpp */
+void write_lines(struct Map_info *, int, const Triangulation &);
+#endif

Added: sandbox/martinl/v.ten/main.cpp
===================================================================
--- sandbox/martinl/v.ten/main.cpp	                        (rev 0)
+++ sandbox/martinl/v.ten/main.cpp	2013-06-01 22:28:53 UTC (rev 56542)
@@ -0,0 +1,98 @@
+#include <cstdlib>
+#include <vector>
+
+extern "C" {
+#include <grass/vector.h>
+#include <grass/glocale.h>
+}
+
+#include "local_proto.h"
+
+int main(int argc, char *argv[])
+{
+    int type;  /* line or face */
+    int field;
+    unsigned int npoints, nvertices;
+    
+    struct GModule *module;
+    
+    struct {
+        struct Option *input, *field, *output;
+    } opt;
+    struct {
+        struct Flag *line;
+    } flag;
+    
+    struct Map_info In, Out;
+    
+    std::vector<Point> points;
+
+    G_gisinit(argv[0]);
+    
+    module = G_define_module();
+    G_add_keyword(_("vector"));
+    G_add_keyword(_("geometry"));
+    G_add_keyword(_("3D triangulation"));
+    module->description = _("Creates a 3D triangulation from an input vector map containing points or centroids.");
+
+    opt.input = G_define_standard_option(G_OPT_V_INPUT);
+
+    opt.field = G_define_standard_option(G_OPT_V_FIELD_ALL);
+
+    opt.output = G_define_standard_option(G_OPT_V_OUTPUT);
+
+    flag.line = G_define_flag();
+    flag.line->key = 'l';
+    flag.line->description =
+	_("Output triangulation as a graph (lines), not faces");
+
+    if (G_parser(argc, argv)) {
+        exit(EXIT_FAILURE);
+    }
+
+    if (flag.line->answer)
+	type = GV_LINE;
+    else
+	type = GV_BOUNDARY;
+
+    /* open input map */
+    Vect_open_old2(&In, opt.input->answer, "", opt.field->answer);
+    Vect_set_error_handler_io(&In, &Out);
+    /* check if the map is 3D */
+    if (!Vect_is_3d(&In))
+        G_fatal_error(_("Vector map <%s> is not 3D"), Vect_get_full_name(&In));
+
+    field = Vect_get_field_number(&In, opt.field->answer);
+    
+    /* create output */
+    Vect_open_new(&Out, opt.output->answer, WITH_Z); /* output is always 3D */
+    Vect_hist_copy(&In, &Out);
+    Vect_hist_command(&Out);
+
+    /* read points */
+    npoints = read_points(&In, field, points);
+    Vect_close(&In);
+    
+    /* do 3D triangulation */
+    G_message(_("Creating TEN..."));
+    Triangulation T(points.begin(), points.end());
+
+    nvertices = T.number_of_vertices();
+    if (nvertices != npoints)
+        G_fatal_error(_("Invalid number of vertices %d (%d)"), 
+                      nvertices, npoints);
+    
+    G_message(_("Number of vertices: %d"), nvertices);
+    G_message(_("Number of edges: %d"), T.number_of_edges());
+    G_message(_("Number of faces: %d"), T.number_of_facets());
+    G_message(_("Number of tetrahedrons: %d"), T.number_of_cells());
+    
+    G_message(_("Writing output features..."));
+
+    write_lines(&Out, type, T);
+
+    Vect_build(&Out);
+    Vect_close(&Out);
+
+    exit(EXIT_SUCCESS);
+}

Added: sandbox/martinl/v.ten/read.cpp
===================================================================
--- sandbox/martinl/v.ten/read.cpp	                        (rev 0)
+++ sandbox/martinl/v.ten/read.cpp	2013-06-01 22:28:53 UTC (rev 56542)
@@ -0,0 +1,51 @@
+extern "C" {
+#include <grass/vector.h>
+#include <grass/glocale.h>
+}
+
+#include "local_proto.h"
+
+int read_points(struct Map_info *Map, int field, std::vector<Point>& OutPoints)
+{
+    int npoints;
+    double x, y, z;
+    
+    struct line_pnts *Points;
+    
+    Points = Vect_new_line_struct();
+
+    /* set constraints */
+    Vect_set_constraint_type(Map, GV_POINTS);
+    if (field > 0)
+        Vect_set_constraint_field(Map, field);
+    
+    /* read points */
+    npoints = 0;
+    G_message(_("Reading points..."));
+    while(TRUE) {
+        if (Vect_read_next_line(Map, Points, NULL) < 0)
+            break;
+
+        G_progress(npoints, 1e3);
+        
+        if (Points->n_points != 1) {
+            G_warning(_("Invalid point skipped"));
+            continue;
+        }
+        
+        x = Points->x[0];
+        y = Points->y[0];
+        z = Points->z[0];
+        
+        G_debug(3, "new point added: %f, %f, %f", x, y, z);
+        OutPoints.push_back(Point(x, y, z));
+        npoints++;
+    }
+    G_progress(1, 1);
+
+    Vect_destroy_line_struct(Points);
+
+    G_debug(1, "read_points(): %d", npoints);
+    
+    return npoints;
+}

Added: sandbox/martinl/v.ten/v.ten.html
===================================================================
--- sandbox/martinl/v.ten/v.ten.html	                        (rev 0)
+++ sandbox/martinl/v.ten/v.ten.html	2013-06-01 22:28:53 UTC (rev 56542)
@@ -0,0 +1,3 @@
+<h2>DESCRIPTION</h2>
+
+TODO

Added: sandbox/martinl/v.ten/write.cpp
===================================================================
--- sandbox/martinl/v.ten/write.cpp	                        (rev 0)
+++ sandbox/martinl/v.ten/write.cpp	2013-06-01 22:28:53 UTC (rev 56542)
@@ -0,0 +1,37 @@
+extern "C" {
+#include <grass/vector.h>
+#include <grass/glocale.h>
+}
+
+#include "local_proto.h"
+
+void write_lines(struct Map_info *Out, int do_lines, const Triangulation &T)
+{
+     Triangulation::Finite_edges_iterator eit;
+     Point pt1, pt2;
+     
+     struct line_pnts *Points;
+     
+     Points = Vect_new_line_struct();
+     
+     if (do_lines) {
+         /* write edges as lines */
+         for (eit = T.finite_edges_begin(); eit != T.finite_edges_end(); ++eit) {
+             pt1 = eit->first->vertex(eit->second)->point();
+             pt2 = eit->first->vertex(eit->third)->point();
+             G_debug(3, "edge: %f %f %f | %f %f %f", pt1.x(), pt1.y(), pt1.z(),
+                     pt2.x(), pt2.y(), pt2.z());
+             
+             Vect_reset_line(Points);
+             Vect_append_point(Points, pt1.x(), pt1.y(), pt1.z());
+             Vect_append_point(Points, pt2.x(), pt2.y(), pt2.z());
+             
+             Vect_write_line(Out, GV_LINE, Points, NULL);
+         }
+     }
+     else {
+         /* TODO */
+     }
+
+     Vect_destroy_line_struct(Points);
+}



More information about the grass-commit mailing list