]> granicus.if.org Git - postgis/commitdiff
Add ST_LineClipZ(geometry, from, to) SQL and C functions.
authorPaul Ramsey <pramsey@cleverelephant.ca>
Sun, 21 Dec 2008 06:37:16 +0000 (06:37 +0000)
committerPaul Ramsey <pramsey@cleverelephant.ca>
Sun, 21 Dec 2008 06:37:16 +0000 (06:37 +0000)
git-svn-id: http://svn.osgeo.org/postgis/trunk@3456 b70326c6-7e19-0410-871a-916f4a2858ee

liblwgeom/cunit/cu_algorithm.c
liblwgeom/cunit/cu_algorithm.h
liblwgeom/lwalgorithm.c
lwgeom/lwgeom_functions_analytic.c
lwgeom/lwpostgis.sql.in.c

index e6da538a7c1161a2dbe7cf3b0583c6a8a3b9c611..a9a8b6fadc0827996d2cc06975ab2752e1f86173 100644 (file)
@@ -33,7 +33,8 @@ CU_pSuite register_cg_suite(void)
            (NULL == CU_add_test(pSuite, "test_lwpoint_set_ordinate()", test_lwpoint_set_ordinate)) || 
            (NULL == CU_add_test(pSuite, "test_lwpoint_get_ordinate()", test_lwpoint_get_ordinate)) ||
            (NULL == CU_add_test(pSuite, "test_lwpoint_interpolate()", test_lwpoint_interpolate)) ||
-           (NULL == CU_add_test(pSuite, "test_lwline_clip()", test_lwline_clip)) 
+           (NULL == CU_add_test(pSuite, "test_lwline_clip()", test_lwline_clip)) ||
+           (NULL == CU_add_test(pSuite, "test_lwline_clip_big()", test_lwline_clip_big)) 
        )
        {
                CU_cleanup_registry();
@@ -674,3 +675,36 @@ void test_lwline_clip(void)
 
 }
 
+void test_lwline_clip_big(void)
+{
+       POINTARRAY *pa = ptarray_construct(1, 0, 3);
+       LWLINE *line = lwline_construct(-1, NULL, pa);
+       LWCOLLECTION *c;
+       char *ewkt;
+       int rv = 0;
+
+       p->x = 0.0;
+       p->y = 0.0;
+       p->z = 0.0;
+       setPoint4d(pa, 0, p);
+       
+       p->x = 1.0;
+       p->y = 1.0;
+       p->z = 1.0;
+       setPoint4d(pa, 1, p);
+       
+       p->x = 2.0;
+       p->y = 2.0;
+       p->z = 2.0;
+       setPoint4d(pa, 2, p);
+       
+       c = lwline_clip_to_ordinate_range(line, 2, 0.5, 1.5);
+       ewkt = lwgeom_to_ewkt((LWGEOM*)c,0);
+       //printf("c = %s\n", ewkt);
+       CU_ASSERT_STRING_EQUAL(ewkt, "MULTILINESTRING((0.5 0.5 0.5,1 1 1,1.5 1.5 1.5))" );
+
+       lwfree(ewkt);
+       lwfree_collection(c);
+       lwfree_line(line);
+}
+
index 0ee5979019a796ea7e56bad792e9250d40d2112f..057351086b5039bbea4ad74fb0ec08e90553d935 100644 (file)
@@ -35,3 +35,5 @@ void test_lwpoint_set_ordinate(void);
 void test_lwpoint_get_ordinate(void);
 void test_lwpoint_interpolate(void);
 void test_lwline_clip(void);
+void test_lwline_clip_big(void);
+
index 352faadd032c4e8c573e71c322130dae8ee894e6..d44ac9da45827f7c84bb34d50f79e0c9f29a427c 100644 (file)
@@ -403,9 +403,9 @@ LWCOLLECTION *lwline_clip_to_ordinate_range(LWLINE *line, int ordinate, double f
        double ordinate_value_p = 0.0, ordinate_value_q = 0.0;
        char hasz = TYPE_HASZ(line->type);
        char hasm = TYPE_HASM(line->type);
+       char dims = TYPE_NDIMS(line->type);
        char hassrid = TYPE_HASSRID(line->type);
 
-
        LWDEBUGF(5, "hassrid = %d", hassrid);
 
        /* Null input, nothing we can do. */
@@ -423,12 +423,10 @@ LWCOLLECTION *lwline_clip_to_ordinate_range(LWLINE *line, int ordinate, double f
                to = t;
        }
        
-       int ndims = TYPE_NDIMS(line->type);
-       
        /* Asking for an ordinate we don't have. Error. */
-       if( ordinate >= ndims ) 
+       if( ordinate >= dims ) 
        {
-               lwerror("Cannot clip on ordinate %d in a %d-d geometry.", ordinate, ndims);
+               lwerror("Cannot clip on ordinate %d in a %d-d geometry.", ordinate, dims);
                return NULL;
        }
        
@@ -471,7 +469,7 @@ LWCOLLECTION *lwline_clip_to_ordinate_range(LWLINE *line, int ordinate, double f
                        if( ! added_last_point ) 
                        {       
                                if( dp ) lwfree(dp);
-                               dp = dynptarray_create(64, ndims);
+                               dp = dynptarray_create(64, line->type);
                        }
                        rv = dynptarray_addPoint4d(dp, p, 0);
             added_last_point = 2; /* Added on a boundary. */
@@ -488,7 +486,7 @@ LWCOLLECTION *lwline_clip_to_ordinate_range(LWLINE *line, int ordinate, double f
                                /* We didn't add the previous point, so this is a new segment.
                                *  Make a new point array. */
                                if( dp ) lwfree(dp);
-                               dp = dynptarray_create(64, ndims);
+                               dp = dynptarray_create(64, line->type);
 
                        /* We're transiting into the range so add an interpolated 
                                *  point at the range boundary. */
@@ -496,13 +494,13 @@ LWCOLLECTION *lwline_clip_to_ordinate_range(LWLINE *line, int ordinate, double f
                    {
                     double interpolation_value;
                     (ordinate_value_q > to) ? (interpolation_value = to) : (interpolation_value = from);
-                    rv = lwpoint_interpolate(q, p, r, ndims, ordinate, interpolation_value);
+                    rv = lwpoint_interpolate(q, p, r, dims, ordinate, interpolation_value);
                     rv = dynptarray_addPoint4d(dp, r, 0);
                    }
                        }
                    /* Add the current vertex to the point array. */
                rv = dynptarray_addPoint4d(dp, p, 0);
-            added_last_point = 1;
+            added_last_point = 1; /* Added inside range. */
                } 
                /* Is this point inside the ordinate range? No. */
                else 
@@ -513,7 +511,7 @@ LWCOLLECTION *lwline_clip_to_ordinate_range(LWLINE *line, int ordinate, double f
                                *  to the point array at the range boundary. */
                 double interpolation_value;
                 (ordinate_value_p > to) ? (interpolation_value = to) : (interpolation_value = from);
-                rv = lwpoint_interpolate(q, p, r, ndims, ordinate, interpolation_value);
+                rv = lwpoint_interpolate(q, p, r, dims, ordinate, interpolation_value);
                 rv = dynptarray_addPoint4d(dp, r, 0);
                }
                        else if ( added_last_point == 2 ) 
@@ -526,10 +524,10 @@ LWCOLLECTION *lwline_clip_to_ordinate_range(LWLINE *line, int ordinate, double f
                                *  so we need to add *two* interpolated points! */
                                pa_out = ptarray_construct(hasz, hasm, 2);
                                /* Interpolate lower point. */
-                               rv = lwpoint_interpolate(p, q, r, ndims, ordinate, from);
+                               rv = lwpoint_interpolate(p, q, r, dims, ordinate, from);
                                setPoint4d(pa_out, 0, r);
                                /* Interpolate upper point. */
-                               rv = lwpoint_interpolate(p, q, r, ndims, ordinate, to);
+                               rv = lwpoint_interpolate(p, q, r, dims, ordinate, to);
                                setPoint4d(pa_out, 1, r);
                        }
                        else if ( ordinate_value_q > to && ordinate_value_p < from ) {
@@ -537,10 +535,10 @@ LWCOLLECTION *lwline_clip_to_ordinate_range(LWLINE *line, int ordinate, double f
                                *  so we need to add *two* interpolated points! */
                                pa_out = ptarray_construct(hasz, hasm, 2);
                                /* Interpolate upper point. */
-                               rv = lwpoint_interpolate(p, q, r, ndims, ordinate, to);
+                               rv = lwpoint_interpolate(p, q, r, dims, ordinate, to);
                                setPoint4d(pa_out, 0, r);
                                /* Interpolate lower point. */
-                               rv = lwpoint_interpolate(p, q, r, ndims, ordinate, from);
+                               rv = lwpoint_interpolate(p, q, r, dims, ordinate, from);
                                setPoint4d(pa_out, 1, r);
                        }
                        /* We have an extant point-array, save it out to a multi-line. */
@@ -553,11 +551,13 @@ LWCOLLECTION *lwline_clip_to_ordinate_range(LWLINE *line, int ordinate, double f
                                        if( dp->pa->npoints == 1 ) 
                                        {
                                                oline = (LWGEOM*)lwpoint_construct(line->SRID, NULL, dp->pa);
+                                               oline->type = lwgeom_makeType(hasz, hasm, hassrid, POINTTYPE);
                                                lwgeom_out->type = lwgeom_makeType(hasz, hasm, hassrid, COLLECTIONTYPE);
                                        }
                                        else
                                        {
                                                oline = (LWGEOM*)lwline_construct(line->SRID, NULL, dp->pa);
+                                               oline->type = lwgeom_makeType(hasz, hasm, hassrid, LINETYPE);
                                        }
                                }
                                else 
@@ -565,7 +565,14 @@ LWCOLLECTION *lwline_clip_to_ordinate_range(LWLINE *line, int ordinate, double f
                                        oline = (LWGEOM*)lwline_construct(line->SRID, NULL, pa_out);
                                }
                                lwgeom_out->ngeoms++;
-                               lwgeom_out->geoms = lwrealloc(lwgeom_out->geoms, sizeof(LWGEOM*) * lwgeom_out->ngeoms);
+                               if( lwgeom_out->geoms ) /* We can't just realloc, since repalloc chokes on a starting null ptr. */
+                               {
+                                       lwgeom_out->geoms = lwrealloc(lwgeom_out->geoms, sizeof(LWGEOM*) * lwgeom_out->ngeoms);
+                               }
+                               else
+                               {
+                                       lwgeom_out->geoms = lwalloc(sizeof(LWGEOM*) * lwgeom_out->ngeoms);
+                               }
                                lwgeom_out->geoms[lwgeom_out->ngeoms - 1] = oline;                              
                                lwgeom_dropBBOX((LWGEOM*)lwgeom_out);
                                lwgeom_addBBOX((LWGEOM*)lwgeom_out);
@@ -582,8 +589,16 @@ LWCOLLECTION *lwline_clip_to_ordinate_range(LWLINE *line, int ordinate, double f
        if( dp && dp->pa->npoints > 0 ) {
                LWGEOM *oline;
                oline = (LWGEOM*)lwline_construct(line->SRID, NULL, dp->pa);
+               oline->type = lwgeom_makeType(hasz, hasm, hassrid, LINETYPE);
                lwgeom_out->ngeoms++;
-               lwgeom_out->geoms = lwrealloc(lwgeom_out->geoms, sizeof(LWGEOM*) * lwgeom_out->ngeoms);
+               if( lwgeom_out->geoms ) /* We can't just realloc, since repalloc chokes on a starting null ptr. */
+               {
+                       lwgeom_out->geoms = lwrealloc(lwgeom_out->geoms, sizeof(LWGEOM*) * lwgeom_out->ngeoms);
+               }
+               else
+               {
+                       lwgeom_out->geoms = lwalloc(sizeof(LWGEOM*) * lwgeom_out->ngeoms);
+               }
                lwgeom_out->geoms[lwgeom_out->ngeoms - 1] = oline;                              
                lwgeom_dropBBOX((LWGEOM*)lwgeom_out);
                lwgeom_addBBOX((LWGEOM*)lwgeom_out);
index 415ecd39beea4547337ac00ef88d92f9282b749c..6700782c0f497769ac42ac92e14172bb7aadb25c 100644 (file)
@@ -37,6 +37,7 @@ LWCOLLECTION *simplify2d_collection(const LWCOLLECTION *igeom, double dist);
 LWGEOM *simplify2d_lwgeom(const LWGEOM *igeom, double dist);
 Datum LWGEOM_simplify2d(PG_FUNCTION_ARGS);
 Datum crossingDirection(PG_FUNCTION_ARGS);
+Datum ST_LineClipZ(PG_FUNCTION_ARGS);
 
 double determineSide(POINT2D *seg1, POINT2D *seg2, POINT2D *point);
 int isOnSegment(POINT2D *seg1, POINT2D *seg2, POINT2D *point);
@@ -993,6 +994,45 @@ Datum crossingDirection(PG_FUNCTION_ARGS)
 
 }
 
+PG_FUNCTION_INFO_V1(ST_LineClipZ);
+Datum ST_LineClipZ(PG_FUNCTION_ARGS)
+{
+       PG_LWGEOM *geom_in = (PG_LWGEOM *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+       double from = PG_GETARG_FLOAT8(1);
+       double to = PG_GETARG_FLOAT8(2);
+       LWCOLLECTION *geom_out = NULL;
+       LWLINE *line_in = NULL;
+       uchar type = (uchar)SERIALIZED_FORM(geom_in)[0];
+       char geomtype = TYPE_GETTYPE(type);
+       char hasz = TYPE_HASZ(type);
+       static int ordinate = 2; /* Z */
+       
+       if ( geomtype != LINETYPE )
+       {
+               elog(ERROR,"This function only accepts LINESTRING as arguments.");
+               PG_RETURN_NULL();
+       }
+
+       if( ! hasz ) 
+       {
+               elog(ERROR,"This function only accepts LINESTRING with Z values as arguments.");
+               PG_RETURN_NULL();
+       }
+
+       line_in = lwline_deserialize(SERIALIZED_FORM(geom_in));
+       geom_out = lwline_clip_to_ordinate_range(line_in, ordinate, from, to);
+       lwfree_line(line_in);
+
+       if( ! geom_out ) {
+               elog(ERROR,"The lwline_clip_to_ordinate_range returned null.");
+               PG_RETURN_NULL();
+       }
+
+       PG_FREE_IF_COPY(geom_in, 0);
+       PG_RETURN_POINTER(pglwgeom_serialize((LWGEOM*)geom_out));
+       
+}
+
 /***********************************************************************
  * --strk@keybit.net
  ***********************************************************************/
index 9f3a7ecaf2c76279f607236812a74a4474ff7e67..4668876f3b3f52bd23037b0024f105f8074d846b 100644 (file)
@@ -3972,6 +3972,12 @@ CREATEFUNCTION ST_CrossingDirection(geometry, geometry)
     AS 'MODULE_PATHNAME', 'crossingDirection'
     LANGUAGE 'C' _IMMUTABLE_STRICT; -- WITH (isstrict,iscachable);
 
+-- Only accepts LINESTRING as parameters.
+-- Availability: 1.4.0
+CREATEFUNCTION ST_LineClipZ(geometry, float8, float8)
+    RETURNS geometry
+    AS 'MODULE_PATHNAME', 'ST_LineClipZ'
+    LANGUAGE 'C' _IMMUTABLE_STRICT; -- WITH (isstrict,iscachable);
 
 #if POSTGIS_GEOS_VERSION >= 30
 -- Requires GEOS >= 3.0.0