Ticket #390: addmeasure.patch

liblwgeom/lwline.c

@@ -524 +524 @@
 {
         setPoint4d(line->points, index, newpoint);
 }
+
+/**
+* Re-write the measure ordinate (or add one, if it isn't already there) interpolating
+* the measure between the supplied start and end values. 
+*/
+LWLINE* 
+lwline_measured_from_lwline(const LWLINE *lwline, double m_start, double m_end)
+{
+        int i = 0;
+        int hasm = 0, hasz = 0;
+        int npoints = 0;
+        double length = 0.0;
+        double length_so_far = 0.0;
+        double m_range = m_end - m_start;
+        double m;
+        POINTARRAY *pa = NULL;
+        POINT3DZ p1, p2;
+        
+        if( TYPE_GETTYPE(lwline->type) != LINETYPE )
+        {
+                lwerror("lwmline_construct_from_lwline: only line types supported");
+                return NULL;
+        }
+        
+        hasz = TYPE_HASZ(lwline->type);
+        hasm = 1;
+        
+        /* Null points or npoints == 0 will result in empty return geometry */
+        if( lwline->points )
+        {
+                npoints = lwline->points->npoints;
+                length = lwgeom_pointarray_length2d(lwline->points);
+                getPoint3dz_p(lwline->points, 0, &p1);
+        }
+        
+        pa = ptarray_construct(hasz, hasm, npoints);
+        
+        for( i = 0; i < npoints; i++ )
+        {               
+                POINT4D q;
+                getPoint3dz_p(lwline->points, i, &p2);
+                length_so_far += distance2d_pt_pt((POINT2D*)&p1, (POINT2D*)&p2);
+                if ( length > 0.0 )
+                        m = m_start + m_range * length_so_far / length;
+                else
+                        m = 0.0;
+                q.x = p2.x;
+                q.y = p2.y;
+                q.z = p2.z;
+                q.m = m;
+                setPoint4d(pa, i, &q);
+                p1 = p2;
+        }
+        
+        return lwline_construct(lwline->SRID, NULL, pa);
+}

liblwgeom/liblwgeom.h

@@ -1379 +1379 @@
 extern POINTARRAY *ptarray_substring(POINTARRAY *, double, double);
 extern double ptarray_locate_point(POINTARRAY *, POINT2D *);
 extern void closest_point_on_segment(POINT2D *p, POINT2D *A, POINT2D *B, POINT2D *ret);
+extern LWLINE *lwline_measured_from_lwline(const LWLINE *lwline, double m_start, double m_end);
 
 /*
  * Ensure every segment is at most 'dist' long.

postgis/lwgeom_functions_lrs.c

@@ -18 +18 @@
 
 
 Datum LWGEOM_locate_between_m(PG_FUNCTION_ARGS);
+Datum ST_AddMeasure(PG_FUNCTION_ARGS);
 
 typedef struct
 {
@@ -535 +536 @@
         PG_RETURN_POINTER(gout);
 }
 
+
+/*
+* CREATE OR REPLACE FUNCTION ST_AddMeasure(geometry, float8, float8) 
+* RETURNS geometry 
+* AS '$libdir/postgis-1.5', 'ST_AddMeasure' 
+* LANGUAGE 'C' IMMUTABLE STRICT;
+*/
+PG_FUNCTION_INFO_V1(ST_AddMeasure);
+Datum ST_AddMeasure(PG_FUNCTION_ARGS)
+{
+        PG_LWGEOM *gin = (PG_LWGEOM *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+        PG_LWGEOM *gout;
+        double start_measure = PG_GETARG_FLOAT8(1);
+        double end_measure = PG_GETARG_FLOAT8(2);
+        LWGEOM *lwin, *lwout;
+
+        /* Raise an error if input is not a linestring */
+        if ( TYPE_GETTYPE(gin->type) != LINETYPE )
+        {
+                lwerror("Only LINESTRING is supported");
+                PG_RETURN_NULL();
+        }
+
+        lwin = pglwgeom_deserialize(gin);
+        lwout = (LWGEOM*)lwline_measured_from_lwline((LWLINE*)lwin, start_measure, end_measure);
+        lwgeom_release(lwin);
+
+        if ( lwout == NULL )
+                PG_RETURN_NULL();
+
+        gout = pglwgeom_serialize(lwout);
+        lwgeom_release(lwout);
+
+        PG_RETURN_POINTER(gout);
+}
+