[postgis-tickets] r16486 - Add function ST_FilterByM #4056

nicklas.aven at jordogskog.no nicklas.aven at jordogskog.no
Fri Mar 23 07:11:53 PDT 2018


Author: nicklas
Date: 2018-03-23 07:11:53 -0700 (Fri, 23 Mar 2018)
New Revision: 16486

Added:
   trunk/liblwgeom/cunit/cu_filterm.c
   trunk/liblwgeom/lwmval.c
   trunk/regress/filterm.sql
   trunk/regress/filterm_expected
Modified:
   trunk/NEWS
   trunk/doc/reference_processing.xml
   trunk/liblwgeom/Makefile.in
   trunk/liblwgeom/cunit/Makefile.in
   trunk/liblwgeom/cunit/cu_tester.c
   trunk/liblwgeom/effectivearea.c
   trunk/liblwgeom/liblwgeom.h.in
   trunk/postgis/lwgeom_functions_basic.c
   trunk/postgis/postgis.sql.in
   trunk/regress/Makefile.in
Log:
Add function ST_FilterByM #4056

Modified: trunk/NEWS
===================================================================
--- trunk/NEWS	2018-03-23 09:24:00 UTC (rev 16485)
+++ trunk/NEWS	2018-03-23 14:11:53 UTC (rev 16486)
@@ -2,6 +2,7 @@
 2018/xx/xx
 
 * New Features *
+  - #4056, ST_FilterByM (Nicklas Avén)
   - #4050, ST_ChaikinSmoothing (Nicklas Avén)
   - #3989, ST_Buffer single sided option (Stephen Knox)
   - #3876, ST_Angle function (Rémi Cura)
@@ -13,6 +14,7 @@
   - #4029, Add ST_QuantizeCoordinates (Dan Baston)
 
 * Breaking Changes *
+  - #4054, ST_SimplifyVW changed from > tolerance to >= tolerance
   - #3885, version number removed from address_standardize lib file
   - #3893, raster support functions can only be loaded in the same schema
            with core PostGIS functions.

Modified: trunk/doc/reference_processing.xml
===================================================================
--- trunk/doc/reference_processing.xml	2018-03-23 09:24:00 UTC (rev 16485)
+++ trunk/doc/reference_processing.xml	2018-03-23 14:11:53 UTC (rev 16486)
@@ -3119,7 +3119,7 @@
 FROM (SELECT  'LINESTRING(5 2, 3 8, 6 20, 7 25, 10 10)'::geometry geom) As foo;
 -result
  simplified
------------+-------------------+
+------------------------------
 LINESTRING(5 2,7 25,10 10)
 
 				</programlisting>
@@ -3130,7 +3130,7 @@
 		  </refsection>
 	</refentry>
 
-    <refentry id="ST_ChaikinSmoothing">
+	<refentry id="ST_ChaikinSmoothing">
 	  <refnamediv>
 		<refname>ST_ChaikinSmoothing</refname>
 		<refpurpose>Returns a "smoothed" version of the given geometry using the Chaikin algorithm</refpurpose>
@@ -3174,7 +3174,7 @@
 FROM (SELECT  'LINESTRING(0 0, 8 8, 0 16)'::geometry geom) As foo;
 -result
  smoothed
------------+-------------------+
+------------------------------
 LINESTRING(0 0,6 6,6 10,0 16)
 
 				</programlisting>
@@ -3184,7 +3184,64 @@
 			<para><xref linkend="ST_Simplify" />, <xref linkend="ST_SimplifyVW" /></para>
 		  </refsection>
 	</refentry>
-    
+
+	<refentry id="ST_FilterByM">
+	  <refnamediv>
+		<refname>ST_FilterByM</refname>
+		<refpurpose>Filters vertex points based on their m-value</refpurpose>
+	  </refnamediv>
+
+	  <refsynopsisdiv>
+		<funcsynopsis>
+		  <funcprototype>
+			<funcdef>geometry <function>ST_FilterByM</function></funcdef>
+			<paramdef><type>geometry</type> <parameter>geom</parameter></paramdef>
+			<paramdef><type>double precision</type> <parameter>min</parameter></paramdef>
+			<paramdef><type>double precision</type> <parameter>max = null</parameter></paramdef>
+			<paramdef><type>boolean</type> <parameter>returnM = false</parameter></paramdef>
+		  </funcprototype>
+		</funcsynopsis>
+	  </refsynopsisdiv>
+
+	  <refsection>
+		<title>Description</title>
+		<para>Filters away vertex points based on their m-value. Returns a geometry with only 
+			vertex points that have a m-value larger or equal to the min value and smaller or equal to
+			the max value. If max-value argument is left out only min value is considered. If fourth argument is left out the m-value 
+			will not be in the resulting geoemtry. If resulting geometry have too few vertex points left for its geometry type an empty 
+			geoemtry will be returned. In a geometry collection 
+			geometries without enough points will just be left out silently. If </para>
+		<para>This function is mainly intended to be used in conjunction with ST_SetEffectiveArea. ST_EffectiveArea sets the effective area 
+			of a vertex in it's m-value. With ST_FilterByM it then is possible to get a simplified version of the geoemtry without any calculations, just by filtering</para>
+		
+		<note><para>There is a difference in what ST_SimplifyVW returns when not enough points meets the creterias compared to ST_FilterByM.
+				ST_SimplifyVW returns the geometry with enough points while ST_FilterByM returns an empty geometry</para></note>
+		<note><para>Note that the retuned geometry might be invalid</para></note>
+		<note><para>This function returns all dimensions, also the z and m-value</para></note>
+		<para>Availability: 2.5.0</para>
+	  </refsection>
+
+		  <refsection>
+			<title>Examples</title>
+			<para>A linestring is filtered</para>
+				<programlisting>
+
+
+SELECT ST_AsText(ST_FilterByM(geom,30)) simplified
+FROM (SELECT  ST_SetEffectiveArea('LINESTRING(5 2, 3 8, 6 20, 7 25, 10 10)'::geometry) geom) As foo;
+-result
+         simplified
+----------------------------
+ LINESTRING(5 2,7 25,10 10)
+
+				</programlisting>
+		  </refsection>
+		  <refsection>
+			<title>See Also</title>
+			<para><xref linkend="ST_SetEffectiveArea" />, <xref linkend="ST_SimplifyVW" /></para>
+		  </refsection>
+	</refentry>
+
     <refentry id="ST_SetEffectiveArea">
 	  <refnamediv>
 		<refname>ST_SetEffectiveArea</refname>

Modified: trunk/liblwgeom/Makefile.in
===================================================================
--- trunk/liblwgeom/Makefile.in	2018-03-23 09:24:00 UTC (rev 16485)
+++ trunk/liblwgeom/Makefile.in	2018-03-23 14:11:53 UTC (rev 16486)
@@ -115,6 +115,7 @@
 	lwunionfind.o \
 	effectivearea.o \
 	lwchaikins.o \
+	lwmval.o \
 	lwkmeans.o \
 	varint.o
 

Modified: trunk/liblwgeom/cunit/Makefile.in
===================================================================
--- trunk/liblwgeom/cunit/Makefile.in	2018-03-23 09:24:00 UTC (rev 16485)
+++ trunk/liblwgeom/cunit/Makefile.in	2018-03-23 14:11:53 UTC (rev 16486)
@@ -37,6 +37,7 @@
 	cu_measures.o \
 	cu_effectivearea.o \
 	cu_chaikin.o \
+	cu_filterm.o \
 	cu_node.o \
 	cu_clip_by_rect.o \
 	cu_libgeom.o \

Added: trunk/liblwgeom/cunit/cu_filterm.c
===================================================================
--- trunk/liblwgeom/cunit/cu_filterm.c	                        (rev 0)
+++ trunk/liblwgeom/cunit/cu_filterm.c	2018-03-23 14:11:53 UTC (rev 16486)
@@ -0,0 +1,64 @@
+/**********************************************************************
+ *
+ * PostGIS - Spatial Types for PostgreSQL
+ * http://postgis.net
+ *
+ * Copyright 2018 Nicklas Avén <nicklas.aven at jordogskog.no>
+ *
+ * This is free software; you can redistribute and/or modify it under
+ * the terms of the GNU General Public Licence. See the COPYING file.
+ *
+ **********************************************************************/
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include "CUnit/Basic.h"
+
+#include "liblwgeom_internal.h"
+#include "cu_tester.h"
+
+static void do_test_filterm(char *geom_txt,char *expected, double min, double max)
+{
+	LWGEOM *geom_in, *geom_out;
+	char *out_txt;
+	geom_in = lwgeom_from_wkt(geom_txt, LW_PARSER_CHECK_NONE);
+	geom_out = lwgeom_filter_m(geom_in,min, max, 1);
+	out_txt = lwgeom_to_wkt(geom_out, WKT_EXTENDED, 3, NULL);
+	if(strcmp(expected, out_txt))
+		printf("%s is not equal to %s\n", expected, out_txt);
+	CU_ASSERT_STRING_EQUAL(expected, out_txt)
+	lwfree(out_txt);
+	lwgeom_free(geom_in);
+	lwgeom_free(geom_out);
+	return;
+}
+
+static void do_test_filterm_single_geometries(void)
+{
+	/*Point*/
+	do_test_filterm("POINT(0 0 0 0)","POINT(0 0 0 0)", 0, 11);
+	do_test_filterm("POINT(0 0 0 0)","POINT EMPTY", 1, 11);
+	/*Line*/
+	do_test_filterm("LINESTRING(0 0 0 0,1 1 0 2,2 2 0 5,3 1 0 3)","LINESTRING(1 1 0 2,3 1 0 3)",2,3);
+	do_test_filterm("LINESTRING(0 0 0 0,1 1 0 2,2 2 0 5,3 1 0 3)","LINESTRING EMPTY",10, 11);
+	/*Polygon*/
+	do_test_filterm("POLYGON((0 0 0 0,1 1 0 2,5 5 0 5,3 1 0 3,0 0 0 0))","POLYGON((0 0 0 0,1 1 0 2,3 1 0 3,0 0 0 0))",0,4);
+	do_test_filterm("POLYGON((0 0 0 0,1 1 0 2,5 5 0 5,3 1 0 3,0 0 0 0))","POLYGON EMPTY",10, 11);
+	return;
+}
+
+static void do_test_filterm_collections(void)
+{
+	do_test_filterm("GEOMETRYCOLLECTION(POINT(1 1 1 1), LINESTRING(1 1 1 1, 1 2 1 4, 2 2 1 3), POLYGON((0 0 0 4,1 1 0 2,5 5 0 5,3 1 0 3,0 0 0 4)))","GEOMETRYCOLLECTION(POINT(1 1 1 1),LINESTRING(1 1 1 1,1 2 1 4,2 2 1 3),POLYGON((0 0 0 4,1 1 0 2,3 1 0 3,0 0 0 4)))",0,4);
+	do_test_filterm("GEOMETRYCOLLECTION(POINT(1 1 1 1), LINESTRING(1 1 1 1, 1 2 1 4, 2 2 1 3), POLYGON((0 0 0 4,1 1 0 2,5 5 0 5,3 1 0 3,0 0 0 4)))","GEOMETRYCOLLECTION(LINESTRING(1 2 1 4,2 2 1 3),POLYGON((0 0 0 4,1 1 0 2,3 1 0 3,0 0 0 4)))",2,4);
+	return;
+}
+
+void filterm_suite_setup(void);
+void filterm_suite_setup(void)
+{
+	CU_pSuite suite = CU_add_suite("filterm",NULL,NULL);
+	PG_ADD_TEST(suite, do_test_filterm_single_geometries);
+	PG_ADD_TEST(suite, do_test_filterm_collections);
+}

Modified: trunk/liblwgeom/cunit/cu_tester.c
===================================================================
--- trunk/liblwgeom/cunit/cu_tester.c	2018-03-23 09:24:00 UTC (rev 16485)
+++ trunk/liblwgeom/cunit/cu_tester.c	2018-03-23 14:11:53 UTC (rev 16486)
@@ -48,6 +48,7 @@
 extern void measures_suite_setup(void);
 extern void effectivearea_suite_setup(void);
 extern void chaikin_suite_setup(void);
+extern void filterm_suite_setup(void);
 extern void minimum_bounding_circle_suite_setup(void);
 extern void misc_suite_setup(void);
 extern void node_suite_setup(void);
@@ -99,6 +100,7 @@
 	measures_suite_setup,
 	effectivearea_suite_setup,
 	chaikin_suite_setup,
+	filterm_suite_setup,
 	minimum_bounding_circle_suite_setup,
 	misc_suite_setup,
 	node_suite_setup,

Modified: trunk/liblwgeom/effectivearea.c
===================================================================
--- trunk/liblwgeom/effectivearea.c	2018-03-23 09:24:00 UTC (rev 16485)
+++ trunk/liblwgeom/effectivearea.c	2018-03-23 14:11:53 UTC (rev 16486)
@@ -332,7 +332,7 @@
 		ea->initial_arealist[after_current].prev = ea->initial_arealist[current].prev;
 
 		/*Check if we are finished*/
-		if((!set_area && ea->res_arealist[current]>trshld) || (ea->initial_arealist[0].next==(npoints-1)))
+		if((!set_area && ea->res_arealist[current]>=trshld) || (ea->initial_arealist[0].next==(npoints-1)))
 			go_on=0;
 
 		i++;
@@ -422,7 +422,7 @@
 		/*Only return points with an effective area above the threshold*/
 		for (p=0;p<ea->inpts->npoints;p++)
 		{
-			if(ea->res_arealist[p]>trshld)
+			if(ea->res_arealist[p]>=trshld)
 			{
 				pt=getPoint4d(ea->inpts, p);
 				pt.m=ea->res_arealist[p];
@@ -435,7 +435,7 @@
 		/*Only return points with an effective area above the threshold*/
 		for (p=0;p<ea->inpts->npoints;p++)
 		{
-			if(ea->res_arealist[p]>trshld)
+			if(ea->res_arealist[p]>=trshld)
 			{
 				pt=getPoint4d(ea->inpts, p);
 				ptarray_append_point(opts, &pt, LW_TRUE);
@@ -505,7 +505,7 @@
 
 	opoly->type = ipoly->type;
 
-	if( lwpoly_is_empty(opoly) )
+	if( lwpoly_is_empty(opoly))
 		return NULL;
 
 	return opoly;

Modified: trunk/liblwgeom/liblwgeom.h.in
===================================================================
--- trunk/liblwgeom/liblwgeom.h.in	2018-03-23 09:24:00 UTC (rev 16485)
+++ trunk/liblwgeom/liblwgeom.h.in	2018-03-23 14:11:53 UTC (rev 16486)
@@ -991,6 +991,7 @@
 
 extern LWGEOM* lwgeom_set_effective_area(const LWGEOM *igeom, int set_area, double area);
 extern LWGEOM* lwgeom_chaikin(const LWGEOM *igeom, int n_iterations, int preserve_endpoint);
+extern LWGEOM* lwgeom_filter_m(LWGEOM *geom, double min, double max, int returnm);
 
 /*
  * Force to use SFS 1.1 geometry type

Added: trunk/liblwgeom/lwmval.c
===================================================================
--- trunk/liblwgeom/lwmval.c	                        (rev 0)
+++ trunk/liblwgeom/lwmval.c	2018-03-23 14:11:53 UTC (rev 16486)
@@ -0,0 +1,269 @@
+
+/**********************************************************************
+ *
+ * PostGIS - Spatial Types for PostgreSQL
+ * http://postgis.net
+ *
+ * PostGIS is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 2 of the License, or
+ * (at your option) any later version.
+ *
+ * PostGIS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with PostGIS.  If not, see <http://www.gnu.org/licenses/>.
+ *
+ **********************************************************************
+ *
+ * Copyright 2018 Nicklas Avén
+ *
+ **********************************************************************/
+
+
+#include "liblwgeom_internal.h"
+
+static LWGEOM* lwgeom_filter_m_ignore_null(LWGEOM *geom,double min,double max, int returnm);
+
+static POINTARRAY* ptarray_filterm(POINTARRAY *pa,double min, double max, int returnm)
+{
+	LWDEBUGF(2, "Entered %s", __func__);
+
+	/*Check if M exists
+	* This should already be tested earlier so we don't handle it properly.
+	* If this happens it is because the function is used in another context than filterM
+	and we throw an error*/
+	if(!FLAGS_GET_M(pa->flags))
+		lwerror("missing m-value in function %s\n",__func__);
+
+	/*Dimensions in input geometry*/
+	int ndims = FLAGS_NDIMS(pa->flags);
+
+	/*Dimensions in result*/
+	int res_ndims;
+	if(returnm)
+		res_ndims = ndims;
+	else
+		res_ndims = ndims-1;
+
+	int pointsize = res_ndims * sizeof(double);
+
+	double m_val;
+
+	//M-value will always be the last dimension
+	int m_pos = ndims-1;
+
+	int i, counter=0;
+	for(i=0;i<pa->npoints;i++)
+	{
+		m_val = *((double*)pa->serialized_pointlist + i*ndims + m_pos);
+		if(m_val >= min && m_val <= max)
+			counter++;
+	}
+
+	POINTARRAY *pa_res = ptarray_construct(FLAGS_GET_Z(pa->flags),returnm * FLAGS_GET_M(pa->flags), counter);
+
+	double *res_cursor = (double*) pa_res->serialized_pointlist;
+	for(i=0;i<pa->npoints;i++)
+	{
+		m_val = *((double*)pa->serialized_pointlist + i*ndims + m_pos);
+		if(m_val >= min && m_val <= max)
+		{
+			memcpy(res_cursor, (double*) pa->serialized_pointlist + i*ndims, pointsize);
+			res_cursor+=res_ndims;
+		}
+	}
+
+	return pa_res;
+
+}
+static LWPOINT* lwpoint_filterm(LWPOINT  *pt,double min,double max, int returnm)
+{
+	LWDEBUGF(2, "Entered %s", __func__);
+
+	POINTARRAY *pa;
+
+	pa = ptarray_filterm(pt->point, min, max,returnm);
+	if(pa->npoints < 1)
+	{
+		ptarray_free(pa);
+		return NULL;
+	}
+	return lwpoint_construct(pt->srid, NULL, pa);
+
+}
+
+static LWLINE* lwline_filterm(LWLINE  *line,double min,double max, int returnm)
+{
+	LWDEBUGF(2, "Entered %s", __func__);
+
+	POINTARRAY *pa;
+	pa = ptarray_filterm(line->points, min, max,returnm);
+
+	if(pa->npoints < 2 )
+	{
+		ptarray_free(pa);
+		return NULL;
+	}
+
+	return lwline_construct(line->srid, NULL, pa);
+}
+
+static LWPOLY* lwpoly_filterm(LWPOLY  *poly,double min,double max, int returnm)
+{
+	LWDEBUGF(2, "Entered %s", __func__);
+	int i, nrings;
+	LWPOLY *poly_res = lwpoly_construct_empty(poly->srid, FLAGS_GET_Z(poly->flags),returnm * FLAGS_GET_M(poly->flags));
+
+	nrings = poly->nrings;
+	for( i = 0; i < nrings; i++ )
+	{
+		/* Ret number of points */
+		POINTARRAY *pa = ptarray_filterm(poly->rings[i], min, max,returnm);
+
+		/* Skip empty rings */
+		if( pa == NULL )
+			continue;
+
+		if(pa->npoints>=4)
+		{
+			if (lwpoly_add_ring(poly_res, pa) == LW_FAILURE )
+			{
+				LWDEBUG(2, "Unable to add ring to polygon");
+				lwerror("Unable to add ring to polygon");
+			}
+		}
+		else if (i==0)
+		{
+			ptarray_free(pa);
+			lwpoly_free(poly_res);
+			return NULL;
+		}
+		else
+			ptarray_free(pa);
+	}
+	return poly_res;
+}
+
+
+
+static LWCOLLECTION* lwcollection_filterm(const LWCOLLECTION *igeom,double min,double max, int returnm)
+{
+	LWDEBUGF(2, "Entered %s",__func__);
+
+	uint32_t i;
+
+	LWCOLLECTION *out = lwcollection_construct_empty(igeom->type, igeom->srid, FLAGS_GET_Z(igeom->flags),returnm * FLAGS_GET_M(igeom->flags));
+
+	if( lwcollection_is_empty(igeom) )
+		return out;
+
+	for( i = 0; i < igeom->ngeoms; i++ )
+	{
+		LWGEOM *ngeom = lwgeom_filter_m_ignore_null(igeom->geoms[i],min, max,returnm);
+		if ( ngeom ) out = lwcollection_add_lwgeom(out, ngeom);
+	}
+
+	return out;
+}
+
+static LWGEOM* lwgeom_filter_m_ignore_null(LWGEOM *geom,double min,double max, int returnm)
+{
+	LWDEBUGF(2, "Entered %s",__func__);
+
+	LWGEOM *geom_out = NULL;
+
+	if(!FLAGS_GET_M(geom->flags))
+		return geom;
+	switch ( geom->type )
+	{
+		case POINTTYPE:
+		{
+			LWDEBUGF(4,"Type found is Point, %d", geom->type);
+			geom_out = lwpoint_as_lwgeom(lwpoint_filterm((LWPOINT*) geom, min, max,returnm));
+			break;
+		}
+		case LINETYPE:
+		{
+			LWDEBUGF(4,"Type found is Linestring, %d", geom->type);
+			geom_out = lwline_as_lwgeom(lwline_filterm((LWLINE*) geom, min,max,returnm));
+			break;
+		}
+		/* Polygon has 'nrings' and 'rings' elements */
+		case POLYGONTYPE:
+		{
+			LWDEBUGF(4,"Type found is Polygon, %d", geom->type);
+			geom_out = lwpoly_as_lwgeom(lwpoly_filterm((LWPOLY*)geom, min, max,returnm));
+			break;
+		}
+
+		/* All these Collection types have 'ngeoms' and 'geoms' elements */
+		case MULTIPOINTTYPE:
+		case MULTILINETYPE:
+		case MULTIPOLYGONTYPE:
+		case COLLECTIONTYPE:
+		{
+			LWDEBUGF(4,"Type found is collection, %d", geom->type);
+			geom_out = (LWGEOM*) lwcollection_filterm((LWCOLLECTION*) geom, min, max,returnm);
+			break;
+		}
+		/* Unknown type! */
+		default:
+			lwerror("Unsupported geometry type: %s [%d] in function %s", lwtype_name((geom)->type), (geom)->type, __func__);
+	}
+	return geom_out;
+
+}
+
+LWGEOM* lwgeom_filter_m(LWGEOM *geom,double min,double max, int returnm)
+{
+	LWDEBUGF(2, "Entered %s",__func__);
+
+	LWGEOM *ngeom = lwgeom_filter_m_ignore_null(geom,min, max,returnm);
+
+	if(ngeom)
+		return ngeom;
+	else
+	{
+		switch ( geom->type )
+		{
+			case POINTTYPE:
+			{
+				return (LWGEOM*) lwpoint_construct_empty(geom->srid, FLAGS_GET_Z(geom->flags), returnm * FLAGS_GET_M(geom->flags));
+				break;
+			}
+			case LINETYPE:
+			{
+				return (LWGEOM*) lwline_construct_empty(geom->srid, FLAGS_GET_Z(geom->flags), returnm * FLAGS_GET_M(geom->flags));
+				break;
+			}
+			/* Polygon has 'nrings' and 'rings' elements */
+			case POLYGONTYPE:
+			{
+				return (LWGEOM*) lwpoly_construct_empty(geom->srid, FLAGS_GET_Z(geom->flags), returnm * FLAGS_GET_M(geom->flags));
+				break;
+			}
+
+			/* All these Collection types have 'ngeoms' and 'geoms' elements */
+			case MULTIPOINTTYPE:
+			case MULTILINETYPE:
+			case MULTIPOLYGONTYPE:
+			case COLLECTIONTYPE:
+			{
+				return (LWGEOM*) lwcollection_construct_empty(geom->type, geom->srid, FLAGS_GET_Z(geom->flags), returnm * FLAGS_GET_M(geom->flags));
+				break;
+			}
+			/* Unknown type! */
+			default:
+				lwerror("Unsupported geometry type: %s [%d] in function %s", lwtype_name((geom)->type), (geom)->type, __func__);
+		}
+	}
+	/*Shouldn't be possible*/
+	return NULL;
+}
+
+
+

Modified: trunk/postgis/lwgeom_functions_basic.c
===================================================================
--- trunk/postgis/lwgeom_functions_basic.c	2018-03-23 09:24:00 UTC (rev 16485)
+++ trunk/postgis/lwgeom_functions_basic.c	2018-03-23 14:11:53 UTC (rev 16486)
@@ -115,8 +115,8 @@
 Datum ST_IsCollection(PG_FUNCTION_ARGS);
 Datum ST_QuantizeCoordinates(PG_FUNCTION_ARGS);
 Datum ST_WrapX(PG_FUNCTION_ARGS);
+Datum LWGEOM_FilterByM(PG_FUNCTION_ARGS);
 
-
 /*------------------------------------------------------------------*/
 
 /** find the size of geometry */
@@ -3023,3 +3023,68 @@
 	PG_FREE_IF_COPY(input, 0);
 	PG_RETURN_POINTER(result);
 }
+
+
+/*
+ * ST_FilterByM(in geometry, val double precision)
+ */
+PG_FUNCTION_INFO_V1(LWGEOM_FilterByM);
+Datum LWGEOM_FilterByM(PG_FUNCTION_ARGS)
+{
+	GSERIALIZED *geom_in;
+	GSERIALIZED *geom_out;
+	LWGEOM *lwgeom_in;
+	LWGEOM *lwgeom_out;
+	double min, max;
+	int returnm;
+	if ( PG_NARGS() > 0 && ! PG_ARGISNULL(0))
+	{
+		geom_in = PG_GETARG_GSERIALIZED_P(0);
+	}
+	else
+	{
+		PG_RETURN_NULL();
+	}
+
+	if ( PG_NARGS() > 1 && ! PG_ARGISNULL(1))
+		min = PG_GETARG_FLOAT8(1);
+	else
+	{
+		min = DBL_MIN;
+	}
+	if ( PG_NARGS() > 2  && ! PG_ARGISNULL(2))
+		max = PG_GETARG_FLOAT8(2);
+	else
+	{
+		max = DBL_MAX;
+	}
+	if ( PG_NARGS() > 3  && ! PG_ARGISNULL(3) && PG_GETARG_BOOL(3))
+		returnm = 1;
+	else
+	{
+		returnm=0;
+	}
+
+	if(min>max)
+	{
+		elog(ERROR,"Min-value cannot be larger than Max value\n");
+		PG_RETURN_NULL();
+	}
+
+	lwgeom_in = lwgeom_from_gserialized(geom_in);
+
+	int hasm = FLAGS_GET_M(lwgeom_in->flags);
+
+	if(!hasm)
+	{
+		elog(NOTICE,"No M-value, No vertex removed\n");
+		PG_RETURN_POINTER(geom_in);
+	}
+
+	lwgeom_out = lwgeom_filter_m(lwgeom_in, min, max,returnm);
+
+	geom_out = geometry_serialize(lwgeom_out);
+	lwgeom_free(lwgeom_out);
+	PG_RETURN_POINTER(geom_out);
+
+}

Modified: trunk/postgis/postgis.sql.in
===================================================================
--- trunk/postgis/postgis.sql.in	2018-03-23 09:24:00 UTC (rev 16485)
+++ trunk/postgis/postgis.sql.in	2018-03-23 14:11:53 UTC (rev 16486)
@@ -3146,6 +3146,13 @@
 	COST 1; -- reset cost, see #3675
 	
 -- Availability: 2.5.0
+CREATE OR REPLACE FUNCTION ST_FilterByM(geometry, double precision, double precision default null, boolean default false)
+	RETURNS geometry
+	AS 'MODULE_PATHNAME', 'LWGEOM_FilterByM'
+	LANGUAGE 'c' IMMUTABLE _PARALLEL
+	COST 1; -- reset cost, see #3675
+	
+-- Availability: 2.5.0
 CREATE OR REPLACE FUNCTION ST_ChaikinSmoothing(geometry, integer default 1, boolean default true)
 	RETURNS geometry
 	AS 'MODULE_PATHNAME', 'LWGEOM_ChaikinSmoothing'

Modified: trunk/regress/Makefile.in
===================================================================
--- trunk/regress/Makefile.in	2018-03-23 09:24:00 UTC (rev 16485)
+++ trunk/regress/Makefile.in	2018-03-23 14:11:53 UTC (rev 16486)
@@ -81,6 +81,7 @@
 	binary \
 	boundary \
 	chaikin \
+	filterm \
 	cluster \
 	concave_hull\
 	ctors \

Added: trunk/regress/filterm.sql
===================================================================
--- trunk/regress/filterm.sql	                        (rev 0)
+++ trunk/regress/filterm.sql	2018-03-23 14:11:53 UTC (rev 16486)
@@ -0,0 +1,2 @@
+SELECT '1', ST_astext(ST_FilterByM('LINESTRING(0 0 0 0, 8 8 5 2, 0 16 0 5, 16 16 0 2)',2,null,'t'));
+SELECT '2', ST_astext(ST_FilterByM('LINESTRING(0 0 0 0, 8 8 5 2, 0 16 0 5, 16 16 0 2)',2,4,'t'));

Added: trunk/regress/filterm_expected
===================================================================
--- trunk/regress/filterm_expected	                        (rev 0)
+++ trunk/regress/filterm_expected	2018-03-23 14:11:53 UTC (rev 16486)
@@ -0,0 +1,2 @@
+1|LINESTRING ZM (8 8 5 2,0 16 0 5,16 16 0 2)
+2|LINESTRING ZM (8 8 5 2,16 16 0 2)



More information about the postgis-tickets mailing list