[GRASS-SVN] r71769 - in grass/trunk/vector/v.profile: . testsuite
svn_grass at osgeo.org
svn_grass at osgeo.org
Sat Nov 18 15:37:09 PST 2017
Author: marisn
Date: 2017-11-18 15:37:09 -0800 (Sat, 18 Nov 2017)
New Revision: 71769
Modified:
grass/trunk/vector/v.profile/Makefile
grass/trunk/vector/v.profile/main.c
grass/trunk/vector/v.profile/testsuite/test_v_profile.py
grass/trunk/vector/v.profile/v.profile.html
Log:
v.profile: Generate buffers with GEOS.
Current state of both GRASS native buffer generation methods
is so bad, that they should be avoided at any cost.
This commit can be reverted as soon as native buffers are fixed.
Modified: grass/trunk/vector/v.profile/Makefile
===================================================================
--- grass/trunk/vector/v.profile/Makefile 2017-11-18 20:38:49 UTC (rev 71768)
+++ grass/trunk/vector/v.profile/Makefile 2017-11-18 23:37:09 UTC (rev 71769)
@@ -3,10 +3,10 @@
PGM = v.profile
-LIBES = $(VECTORLIB) $(DBMILIB) $(GISLIB)
+LIBES = $(VECTORLIB) $(DBMILIB) $(GISLIB) $(GEOSLIBS)
DEPENDENCIES = $(VECTORDEP) $(DBMIDEP) $(GISDEP)
EXTRA_INC = $(VECT_INC)
-EXTRA_CFLAGS = $(VECT_CFLAGS)
+EXTRA_CFLAGS = $(VECT_CFLAGS) $(GEOSCFLAGS)
include $(MODULE_TOPDIR)/include/Make/Module.make
Modified: grass/trunk/vector/v.profile/main.c
===================================================================
--- grass/trunk/vector/v.profile/main.c 2017-11-18 20:38:49 UTC (rev 71768)
+++ grass/trunk/vector/v.profile/main.c 2017-11-18 23:37:09 UTC (rev 71769)
@@ -37,6 +37,78 @@
#include "local_proto.h"
+#if HAVE_GEOS
+#include <float.h>
+
+/* A copy of ring2pts from v.buffer geos.c
+ * It is a terrible approach to copy/pasta functions around,
+ * but all this GEOS stuff is just a temporary solution before native
+ * buffers are fixed. (Will happen at some point after next ten years
+ * or so as there's nothing more permanent than a temporary solution.)
+ * 2017-11-19
+ */
+static int ring2pts(const GEOSGeometry *geom, struct line_pnts *Points)
+{
+ int i, ncoords;
+ double x, y, z;
+ const GEOSCoordSequence *seq = NULL;
+
+ G_debug(3, "ring2pts()");
+
+ Vect_reset_line(Points);
+ if (!geom) {
+ G_warning(_("Invalid GEOS geometry!"));
+ return 0;
+ }
+ z = 0.0;
+ ncoords = GEOSGetNumCoordinates(geom);
+ if (!ncoords) {
+ G_warning(_("No coordinates in GEOS geometry (can be ok for negative distance)!"));
+ return 0;
+ }
+ seq = GEOSGeom_getCoordSeq(geom);
+ for (i = 0; i < ncoords; i++) {
+ GEOSCoordSeq_getX(seq, i, &x);
+ GEOSCoordSeq_getY(seq, i, &y);
+ if (x != x || x > DBL_MAX || x < -DBL_MAX)
+ G_fatal_error(_("Invalid x coordinate %f"), x);
+ if (y != y || y > DBL_MAX || y < -DBL_MAX)
+ G_fatal_error(_("Invalid y coordinate %f"), y);
+ Vect_append_point(Points, x, y, z);
+ }
+
+ return 1;
+}
+
+/* Helper for converting multipoligons to GRASS poligons */
+static int add_poly(const GEOSGeometry *OGeom, struct line_pnts *Buffer) {
+ const GEOSGeometry *geom2;
+ static struct line_pnts *gPoints;
+ int i, nrings;
+
+ gPoints = Vect_new_line_struct();
+
+ geom2 = GEOSGetExteriorRing(OGeom);
+ if (!ring2pts(geom2, gPoints)) {
+ G_fatal_error(_("Corrupt GEOS geometry"));
+ }
+
+ Vect_append_points(Buffer, gPoints, GV_FORWARD);
+ Vect_reset_line(gPoints);
+
+ nrings = GEOSGetNumInteriorRings(OGeom);
+
+ for (i = 0; i < nrings; i++) {
+ geom2 = GEOSGetInteriorRingN(OGeom, i);
+ if (!ring2pts(geom2, gPoints)) {
+ G_fatal_error(_("Corrupt GEOS geometry"));
+ }
+ Vect_append_points(Buffer, gPoints, GV_FORWARD);
+ Vect_reset_line(gPoints);
+ }
+}
+#endif
+
static int compdist(const void *, const void *);
Result *resultset;
@@ -197,6 +269,16 @@
if (G_parser(argc, argv))
exit(EXIT_FAILURE);
+
+#if HAVE_GEOS
+#if (GEOS_VERSION_MAJOR < 3 || (GEOS_VERSION_MAJOR >= 3 && GEOS_VERSION_MINOR < 3))
+ G_fatal_error("This module requires GEOS >= 3.3");
+#endif
+ initGEOS(G_message, G_fatal_error);
+#else
+ G_fatal_error("GRASS native buffering functions are known to return incorrect results.\n"
+ "Till those errors are fixed, this module requires GRASS to be compiled with GEOS support.");
+#endif
/* Start with simple input validation and then move to more complex ones */
otype = Vect_option_to_types(type_opt);
@@ -404,8 +486,47 @@
/* Create a buffer around profile line for point sampling
Tolerance is calculated in such way that buffer will have flat end and no cap. */
+ /* Native buffering is known to fail.
Vect_line_buffer(Profil, bufsize, 1 - (bufsize * cos((2 * M_PI) / 2)),
Buffer);
+ */
+#ifdef HAVE_GEOS
+ /* Code lifted from v.buffer geos.c (with modifications) */
+ GEOSGeometry *IGeom;
+ GEOSGeometry *OGeom = NULL;
+ const GEOSGeometry *geom2 = NULL;
+
+ IGeom = Vect_line_to_geos(Profil, GV_LINE, 0);
+ if (!IGeom) {
+ G_fatal_error(_("Failed to convert GRASS line to GEOS line"));
+ }
+
+ GEOSBufferParams* geos_params = GEOSBufferParams_create();
+ GEOSBufferParams_setEndCapStyle(geos_params, GEOSBUF_CAP_FLAT);
+ OGeom = GEOSBufferWithParams(IGeom, geos_params, bufsize);
+ GEOSBufferParams_destroy(geos_params);
+ if (!OGeom) {
+ G_fatal_error(_("Buffering failed"));
+ }
+
+ if (GEOSGeomTypeId(OGeom) == GEOS_MULTIPOLYGON) {
+ int ngeoms = GEOSGetNumGeometries(OGeom);
+ for (i = 0; i < ngeoms; i++) {
+ geom2 = GEOSGetGeometryN(OGeom, i);
+ add_poly(geom2, Buffer);
+ }
+ }
+ else {
+ add_poly(OGeom, Buffer);
+ }
+
+ if (IGeom)
+ GEOSGeom_destroy(IGeom);
+ if (OGeom)
+ GEOSGeom_destroy(OGeom);
+ finishGEOS();
+#endif
+
Vect_cat_set(Cats, 1, 1);
/* Should we store used buffer for later examination? */
Modified: grass/trunk/vector/v.profile/testsuite/test_v_profile.py
===================================================================
--- grass/trunk/vector/v.profile/testsuite/test_v_profile.py 2017-11-18 20:38:49 UTC (rev 71768)
+++ grass/trunk/vector/v.profile/testsuite/test_v_profile.py 2017-11-18 23:37:09 UTC (rev 71769)
@@ -36,13 +36,40 @@
2|24.34|1029|999647|"Greshams Lake Dam"|"Dam"|"NC"|37|"Wake"|183|"078:34:31W"|"35:52:43N"|35.878484|-78.57528|""|""|""|""|77|"Wake Forest"
"""
+buf_points = """\
+626382.68026139|228917.44816672|1
+626643.91393428|228738.2879083|2
+626907.14939778|228529.10079092|3
+"""
+buf_output = """\
+Number,Distance,cat,dbl_1,dbl_2,int_1
+1,2102.114,3,626907.14939778,228529.10079092,3
+2,2854.300,2,626643.91393428,228738.2879083,2
+3,2960.822,1,626382.68026139,228917.44816672,1
+"""
class TestProfiling(TestCase):
+ to_remove = []
+ points = 'test_v_profile_points'
in_points = 'points_of_interest'
in_map = 'roadsmajor'
where = "cat='354'"
prof_ponts = (647952, 236176, 647950, 236217)
outfile = 'test_out.csv'
+
+ @classmethod
+ def setUpClass(cls):
+ """Create vector map with points for sampling"""
+ cls.runModule('v.in.ascii', input='-', output=cls.points,
+ stdin=buf_points)
+ cls.to_remove.append(cls.points)
+
+ @classmethod
+ def tearDownClass(cls):
+ """Remove vector maps"""
+ if cls.to_remove:
+ cls.runModule('g.remove', flags='f', type='vector',
+ name=','.join(cls.to_remove), verbose=True)
def testParsing(self):
"""Test input parameter parsing"""
@@ -99,6 +126,13 @@
vpro = SimpleModule('v.profile', input=self.in_points, coordinates=self.prof_ponts, buffer=200)
vpro.run()
self.assertLooksLike(reference=output_coords, actual=vpro.outputs.stdout)
+
+ def testBuffering(self):
+ """Test against errors in buffering implementation"""
+ vpro = SimpleModule('v.profile', input=self.points, separator='comma', dp=3,
+ buffer=500, profile_map='roadsmajor at PERMANENT', profile_where='cat=193')
+ vpro.run()
+
if __name__ == '__main__':
Modified: grass/trunk/vector/v.profile/v.profile.html
===================================================================
--- grass/trunk/vector/v.profile/v.profile.html 2017-11-18 20:38:49 UTC (rev 71768)
+++ grass/trunk/vector/v.profile/v.profile.html 2017-11-18 23:37:09 UTC (rev 71769)
@@ -18,10 +18,18 @@
<h2>NOTES</h2>
+<p>
Currently the module can profile only points and lines (including 3D ones).
Areas and other complex features are not supported. If in future users can
provide reasonable examples how area sampling should work and why it is
important, area (or any other feature type) sampling can be added.
+</p>
+<p>
+Due to bugs in GRASS native buffering algorithms, this module for now
+depends on GEOS and will not function if GRASS is compiled without GEOS.
+This restriction will be removed as soon as GRASS native buffer generation
+is fixed.
+</p>
<h2>EXAMPLES</h2>
More information about the grass-commit
mailing list