/********************************************************************** * $Id: postgis_fn.c,v 1.32 2004/02/12 10:34:49 cspoerri Exp $ * * PostGIS - Spatial Types for PostgreSQL * http://postgis.refractions.net * Copyright 2001-2003 Refractions Research Inc. * * This is free software; you can redistribute and/or modify it under * the terms of hte GNU General Public Licence. See the COPYING file. * * ********************************************************************** * $Log: postgis_lrs.c,v $ * Revision 1.0.1 2005/07/22 20:00:00 cspoerri * Revision 1.0.2 2005/08/26 20:00:00 cspoerri * Revision 1.0.3 2005/10/30 20:00:00 cspoerri * **********************************************************************/ /* Installation: 1. copy this file into the lwgeom directory of the postgis distribution 2. add 'lwgeom_functions_lrs.o' to the Makefile (in the lwgeom directory) at the line where it says 'OBJS=$(SA_OBJS) lib.....' 3. Note: If you don't want the debug info comment out the line #define DEBUG in this file 4. recompile the postgis and install 5. run the following SQL to install/register the new functions: CREATE OR REPLACE FUNCTION lrs_create(geometry, int4) RETURNS geometry AS '$libdir/liblwgeom.so.1','lrs_create' LANGUAGE 'C' IMMUTABLE STRICT; -- WITH (isstrict,iscachable); CREATE OR REPLACE FUNCTION lrs_create_cust(geometry, int4, int4) RETURNS geometry AS '$libdir/liblwgeom.so.1','lrs_create_cust' LANGUAGE 'C' IMMUTABLE STRICT; -- WITH (isstrict,iscachable); CREATE OR REPLACE FUNCTION lrs_getmeasure(geometry, geometry) RETURNS float8 AS '$libdir/liblwgeom.so.1','lrs_getmeasure' LANGUAGE 'C' IMMUTABLE STRICT; -- WITH (isstrict,iscachable); CREATE OR REPLACE FUNCTION lrs_getgeom(geometry, float8) RETURNS geometry AS '$libdir/liblwgeom.so.1','lrs_getgeom_point' LANGUAGE 'C' IMMUTABLE STRICT; -- WITH (isstrict,iscachable); CREATE OR REPLACE FUNCTION lrs_getgeom(geometry, float8, float8) RETURNS geometry AS '$libdir/liblwgeom.so.1','lrs_getgeom_line' LANGUAGE 'C' IMMUTABLE STRICT; -- WITH (isstrict,iscachable); CREATE OR REPLACE FUNCTION lrs_getgeom_offset(geometry, float8, float8, float8) RETURNS geometry AS '$libdir/liblwgeom.so.1','lrs_getgeom_pointoffset' LANGUAGE 'C' IMMUTABLE STRICT; -- WITH (isstrict,iscachable); 5. Done Functions: lrs_create(geom, int): returns a new LRS geometry int = 0: LRS uses line length for line measures, the measure for a given point is the same as the distance of this point from the start point of the line (Note: distance is measured along the line) int = 1: LRS uses % as line measure, the measure for a given point is the same as the ratio = (distance point -> start point)/total line length lrs_create_cust(geom, int1, int2): returns a new LRS geometry geom must be a LINESTRING or MULTILINESTRING (2D, 3DM/D, or 4D), if 3DM/4D the existing M values are overwritten int1: start measure for the M values int2: end measures for the M values lrs_getmeasure(geom1, geom2): returns the measure of the closet 'point' on geom1 from geom2 geom1 must have a M value geom2 must be a point geom lrs_getgeom(geom, double): returns a point along the LRS geom must be a LINESTRING or MULTILINESTRING with M values double is the measure along the geom for which you would like a point returned, should be within the measure range of geom lrs_getgeom(geom, double1, double2): return a new linestring along the LRS geom must be a LINESTRING or MULTILINESTRING with M values double1 is the start measure/point of the new linestring, should be within the measure range of geom double2 is the end measure/point of the new linestring, should be within the measure range of geom Note: if double1 or double2 are out of range, the LRS end points are assumed and used instead lrs_getgeom_offset(geom, double1, double2, double3): returns a new point offset from LRS, geom must be a LINESTRING or MULTILINESTRING with M values double1 is the measure along the geom for which you would like a point returned, should be within the measure range of geom double2 is the radian (0 -> 2* Pi) measure to calculate offset. 0 => new point is offset to the right, Pi = > new point is offset to the left, etc. double3 is the distance by which the new point should be offset TODO: - if a Point is calculated for a MULTILINESTRING, and the measure is one that falls on two points (happens for points that are not connected, but are in sequence), we get a point that doesn't lay on the line, but between the lines. - make sure that the SRID in the returned geometry is the same as in the input geometry - original Z values are set to 0, should be set to the same value as Z value of input geometry - for 3D (2.5D) features, calculations should be bassed on x, y and z. Comments: - 3D features are not supported yet, all calculations are done in 2D space - FLT_EPSILON was used to compare small numbers, since DBL_EPSILON didn't provide the correct results (not sure why). Yet, FLT_EPSILON should be fine for most applications. History: Oct 30, 2005: - started to code offset function for points (not available yet through database) - bug fix in lrs_create function Oct 12, 2005: - code clean up - fixed a bug with getmeasure function (returned measures were fixed at range 0 - 100). Now the function returns any measures, which means that it can interpolate beyond the line - added argument hasZ to addlastPoint() Oct 4, 2005: - renamed functions to make make consistent, the all start with 'lrs_' now - combined the getpoint and getline function to lrs_getgeom() for simplicity, they differe now by signature */ #define DEBUG #include #include #include #include #include #include "postgres.h" #include "fmgr.h" #include "utils/elog.h" #include "utils/array.h" #include "utils/geo_decls.h" #include "liblwgeom.h" #include "lwgeom_pg.h" #include "profile.h" #define MAX_MEAS 100 #define MIN_MEAS 0 Datum lrs_create(PG_FUNCTION_ARGS); Datum lrs_create_cust(PG_FUNCTION_ARGS); Datum lrs_getmeasure(PG_FUNCTION_ARGS); Datum lrs_getgeom_point(PG_FUNCTION_ARGS); Datum lrs_getgeom_line(PG_FUNCTION_ARGS); Datum lrs_getgeom_pointoffset(PG_FUNCTION_ARGS); LWGEOM *lrs_create_common(char *geom, int32 measureType, double from, double to); LWGEOM *calcMeasures3DM(char *new_geom, int32 measureType, double from, double mRatio); bool calcMeasures4D(char *new_geom, int32 measureType, double from, double mRatio); LWGEOM *getLRSPoint4D (char * lwgeom1, double m); LWGEOM *getLRSPoint3DM (char * lwgeom1, double m); LWGEOM *getLRSLine4D (const char * lwgeom1, double frm_m, double to_m, char isReversed); LWGEOM *getLRSLine3DM (const char * lwgeom1, double frm_m, double to_m, char isReversed); LWGEOM *addLastPoint(char *opts , int npts, char* pt1, int ptsize, int hasZ, LWGEOM *resLine); PG_FUNCTION_INFO_V1(lrs_create_cust); Datum lrs_create_cust(PG_FUNCTION_ARGS) { PG_LWGEOM *geom = (PG_LWGEOM *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0)); double from = PG_GETARG_FLOAT8(1); double to = PG_GETARG_FLOAT8(2); LWGEOM *res; if (!((TYPE_GETTYPE(geom->type) == LINETYPE) || (TYPE_GETTYPE(geom->type) == MULTILINETYPE) )) { elog(ERROR,"create_lrs: 1st arg isnt of type LINE or MULTILINE"); PG_RETURN_NULL(); } LWGEOM *in = lwgeom_deserialize(SERIALIZED_FORM(geom)); // res = lrs_create_common(in, 1, from, to); PG_RETURN_POINTER(pglwgeom_serialize(res)); } PG_FUNCTION_INFO_V1(lrs_create); Datum lrs_create(PG_FUNCTION_ARGS) { PG_LWGEOM *geom = (PG_LWGEOM *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0)); int32 measureType = PG_GETARG_INT32(1); LWGEOM *res; if (!((TYPE_GETTYPE(geom->type) == LINETYPE) || (TYPE_GETTYPE(geom->type) == MULTILINETYPE) )) { elog(ERROR,"create_lrs: 1st arg isnt of type LINE or MULTILINE"); PG_RETURN_NULL();; } if (measureType) { res = lrs_create_common(SERIALIZED_FORM(geom), 1, 0, 100); } else { res = lrs_create_common(SERIALIZED_FORM(geom), 0, 0, 0); } PG_RETURN_POINTER(pglwgeom_serialize(res)); } LWGEOM *lrs_create_common(char *geom, int32 measureType, double from, double to) { /* * Purpose: creates a new linear referencing system (calculates the M values). * The M values can be a percentage of the line length (0 - 100) at which the given * point occurs. Or it can be the length from the start point to the current point. * In case we have a Multilinestring, the entire length of the line is considered * for the calculations, yet the gaps between the line segements is NOT considered. * Input: * 1. Geometry of type Linestring or Multilinestring (2,3,4D) * 2. Int if 1 => uses percentage for measure, otherwise uses length * Output: Geometry with M values, in case of 3DM and 4D the existing M values are * overwritten * Calls: * calcMeasures3DM * calcMeasures4D */ //LWGEOM *new_geom; // the new geometry, used during the processing size_t size = 0; // used in calc of size of new geometry int32 i; //used for looping double totalDist = 0.0; //the total distance of the line int olddims; // the dimension of the original geometry char *srlOut; // the output geometry returned from the function char *serialized_form; double mRatio = 0; // the ratio used to calculate the measure in case lenght is not used olddims = lwgeom_ndims(geom[0]); #ifdef DEBUG elog(INFO,"checking dimension and ZM values"); #endif if ( olddims < 4) { // allocate double as memory a larger for safety //size = sizeof(POINT3D * geom->points; size = lwgeom_size(geom) * 2; } else { //size = geom->size; size = lwgeom_size(geom); } #ifdef DEBUG elog(INFO,"allocate new geometry"); #endif srlOut = lwalloc(size); if (olddims < 3 || TYPE_HASM(geom[0])) { lwgeom_force3dm_recursive(geom,srlOut, &size); } else { lwgeom_force4d_recursive(geom,srlOut, &size); } // result->size = size + 4; //new_geom = lwgeom_deserialize(SERIALIZED_FORM(result)); // serialized_form = geom; //SERIALIZED_FORM(geom); // we need the entire length of the line in case we should create percent measures #ifdef DEBUG elog(INFO,"calculate new length"); #endif int totPoints = 0; if (measureType) { #ifdef DEBUG elog(INFO,"getting distance for perc/custom measures"); // printType(geom->type); elog(INFO,"num of subgeoma is %i",lwgeom_getnumgeometries(srlOut)); #endif for (i=0; ipoints); } #ifdef DEBUG elog(INFO,"totalDist = %f",totalDist); #endif mRatio = abs((from - to) )/ totalDist; } if (olddims < 3 || TYPE_HASM(geom[0])) { // if ((result=calcMeasures3DM(serialized_form, measureType, from, mRatio))==NULL) if (!calcMeasures3DM(srlOut, measureType, from, mRatio)) { elog(ERROR,"create_lrs: the shape is not of LINE or MULTILINE type."); return NULL; } } else { if (!calcMeasures4D(srlOut, measureType, from, mRatio)) { elog(ERROR,"create_lrs: the shape is not of LINE or MULTILINE type."); return NULL; } } #ifdef DEBUG elog(INFO,"test A"); #endif return lwgeom_deserialize(srlOut); } LWGEOM *calcMeasures3DM(char *new_geom, int32 measureType, double startM, double mRatio) { /* * Purpose: calculates the measure for a 3D geometry. * Input: * 1. Geometry of type Linestring or Multilinestring (2,3,4D) * 2. Int if 1 => uses percentage for measure, otherwise uses length * 3. Double which is the total length of the geometry * Output: Geometry with M values, existing M values are overwritten * Called by: * create_lrs * Author: Christoph Spoerri */ int j, i, npts = 0; double dist1 = 0.0; POINTARRAY *ipa, *opa; POINT3DM *frm3M, *to3M; LWGEOM *g1; char *outPts; int ptsize = sizeof(POINT3DM); LWLINE *resLine; LWGEOM **geomColl; int numGeom=0; #ifdef DEBUG elog(INFO,"calculate new measures for 3D line"); #endif //geomColl = lwalloc(sizeof(LWGEOM *)*lwgeom_getnumgeometries(new_geom)); // now go through the entire line and assign M values to all the points for (i=0; itype) == LINETYPE) { ipa = ((LWLINE *)g1)->points; } else { return FALSE; } #ifdef DEBUG elog(INFO,"test3"); #endif npts = 0; frm3M = (POINT3DM *) getPoint_internal(ipa,npts); npts++; #ifdef DEBUG elog(INFO,"test4"); #endif if (measureType) { frm3M->m = startM + dist1 * mRatio; } else { frm3M->m = dist1; } #ifdef DEBUG elog(INFO,"test5"); #endif // for (j=1;jnpoints;j++) do { to3M = (POINT3DM *) getPoint_internal(ipa,npts); // memcpy(outPts + ptsize * npts, getPoint_internal(ipa,npts), ptsize); // to3M = (POINT3DM *) outPts[ptsize*npts]; npts++; dist1 += sqrt( ( (frm3M->x - to3M->x)*(frm3M->x - to3M->x) ) + ((frm3M->y - to3M->y)*(frm3M->y - to3M->y) )); #ifdef DEBUG elog(INFO,"test6"); #endif frm3M = to3M; #ifdef DEBUG elog(INFO,"test7"); #endif if (measureType) { frm3M->m = startM + dist1 * mRatio; } else { frm3M->m = dist1; } } while (npts < ipa->npoints); } #ifdef DEBUG elog(INFO,"test10"); #endif return TRUE; } bool calcMeasures4D(char *new_geom, int32 measureType, double startM, double mRatio) { /* * Purpose: calculates the measure for a 4D geometry. * Input: * 1. Geometry of type Linestring or Multilinestring (2,3,4D) * 2. Int if 1 => uses percentage for measure, otherwise uses length * 3. Double which is the total length of the geometry * Output: Geometry with M values, existing M values are overwritten * Called by: * create_lrs * Author: Christoph Spoerri */ int j, i; double dist1 = 0.0; POINTARRAY *ipa; POINT4D *frm4D, *to4D; LWGEOM *g1; #ifdef DEBUG elog(INFO,"calculate new measures"); #endif // now go through the entire line and assign M values to all the points for (i=0; itype) == LINETYPE) { ipa = ((LWLINE *)g1)->points; } else { return FALSE; } frm4D = (POINT4D *) getPoint_internal(ipa,0); if (measureType) { frm4D->m = startM + dist1 * mRatio; } else { frm4D->m = dist1; } for (j=1;jnpoints;j++) { to4D = (POINT4D *) getPoint_internal(ipa,j); dist1 += sqrt( ( (frm4D->x - to4D->x)*(frm4D->x - to4D->x) ) + ((frm4D->y - to4D->y)*(frm4D->y - to4D->y) )); frm4D = to4D; if (measureType) { frm4D->m = startM + dist1 * mRatio; } else { frm4D->m = dist1; } } } return TRUE; } PG_FUNCTION_INFO_V1(lrs_getmeasure); Datum lrs_getmeasure(PG_FUNCTION_ARGS) { /* * Purpose: calculates the M value for any point along the line. In case the * point does not lay on the, any point on the line closest to the given point * is used to determine the M value * Input: * 1. Geometry that contains the LRS system (M values) * 2. Geometry of type POINT which is used as reference to calculate M * Output: * Double which is the calculated M value * Author: Christoph Spoerri */ PG_LWGEOM *geom1 = (PG_LWGEOM *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0)); PG_LWGEOM *geom2 = (PG_LWGEOM *) PG_DETOAST_DATUM(PG_GETARG_DATUM(1)); int32 j,i; double dist = 0.0, meas = 0.0, mindist = 99999999999.99; char *serform_geom1; LWPOINT *lwpt; POINT2D pt; POINTARRAY *ipa, *minLine; POINT2D frm, to; int minPoint; LWGEOM *g1; POINT3DM lwMinFrm, lwMinTo; LWPOINT *minFrmPt, *minToPt; #ifdef DEBUG elog(INFO,"start getmeasure stuff"); #endif if (!((TYPE_GETTYPE(geom1->type) == LINETYPE) || (TYPE_GETTYPE(geom1->type) == MULTILINETYPE) )) { elog(ERROR,"getmeasure: 1st arg isnt of type LINE or MULTILINE"); PG_RETURN_FLOAT8(-1); } if (!TYPE_HASM(geom1->type)) { elog(ERROR,"getmeasure: line doesn't contain any M values"); PG_RETURN_FLOAT8(-1); } serform_geom1 = SERIALIZED_FORM(geom1); if (TYPE_GETTYPE(geom2->type) != POINTTYPE) { elog(ERROR,"getmeasure: 2nd arg isnt of type POINT"); PG_RETURN_FLOAT8(-1); } lwpt = lwgeom_getpoint(SERIALIZED_FORM(geom2),0); pt = getPoint2d(lwpt->point,0); #ifdef DEBUG elog(INFO,"start searching"); #endif for (i=0; i< lwgeom_getnumgeometries(serform_geom1); i++) //for each object in geom1 { g1 = (LWGEOM *)lwgeom_deserialize(lwgeom_getsubgeometry(serform_geom1,i)); ipa = ((LWLINE *)g1)->points; frm = getPoint2d(ipa,0); for (j=1; j < ipa->npoints; j++) { to = getPoint2d(ipa,j); dist = distance2d_pt_seg(&pt, &frm, &to); #ifdef DEBUG elog(INFO,"dist = %f",dist); #endif if (dist < mindist) { mindist = dist; minLine = ipa; minPoint = j; } frm = to; } } #ifdef DEBUG elog(INFO,"mindist = %0.20f",mindist); #endif if (TYPE_NDIMS(geom1->type) == 4) { POINT4D frm4, to4; frm4 = getPoint4d(minLine,minPoint-1); lwMinFrm.x = frm4.x; lwMinFrm.y = frm4.y; lwMinFrm.m = frm4.m; to4 = getPoint4d(minLine,minPoint); lwMinTo.x = to4.x; lwMinTo.y = to4.y; lwMinTo.m = to4.m; } else { lwMinFrm = getPoint3dm(minLine,minPoint-1); lwMinTo = getPoint3dm(minLine,minPoint); } minFrmPt = make_lwpoint2d(-1,lwMinFrm.x, lwMinFrm.y); minToPt = make_lwpoint2d(-1, lwMinTo.x, lwMinTo.y); #ifdef DEBUG elog(INFO,"meas DBL_EPSILON = %0.20f",DBL_EPSILON); elog(INFO,"meas FLT_EPSILON = %0.20f",FLT_EPSILON); #endif // SRID needs to be corrected !!! if (fabs(mindist - distance2d_point_point(lwpt,minFrmPt)) < FLT_EPSILON) { meas = lwMinFrm.m; #ifdef DEBUG elog(INFO,"meas (start) = %f", meas); #endif } // SRID needs to be corrected !!! else if (fabs(mindist - distance2d_point_point(lwpt, minToPt)) < FLT_EPSILON) { meas = lwMinTo.m; #ifdef DEBUG elog(INFO,"meas (end) = %f", meas); #endif } else { double c, d1, d2; //get distance between point and from_point c = (lwMinFrm.x-pt.x)*(lwMinFrm.x-pt.x) + (lwMinFrm.y-pt.y)*(lwMinFrm.y-pt.y); #ifdef DEBUG elog(INFO,"c (1) = %f", c); #endif //get distnace between 'intersection' point and from point d1 = sqrt(c - mindist * mindist); d2 = sqrt((lwMinFrm.x-lwMinTo.x)*(lwMinFrm.x-lwMinTo.x) + (lwMinFrm.y-lwMinTo.y)*(lwMinFrm.y-lwMinTo.y)); #ifdef DEBUG elog(INFO,"d1 (1a) = %f , mindist = %f", d1,mindist); #endif meas = d1/d2 * (lwMinTo.m - lwMinFrm.m); #ifdef DEBUG elog(INFO,"meas (2) = %f, pt.x = %f , min_frm.m = %f , min_to.m = %f", meas, pt.x, lwMinFrm.m, lwMinTo.m); #endif // now calc the new measure value for the 'intersection' meas = lwMinFrm.m + meas; } #ifdef DEBUG elog(INFO,"meas (3) = %f", meas); #endif /* // make sure that we do not have measures larger or smaller than possible if (meas > MAX_MEAS) { meas = MAX_MEAS; } else if (meas < MIN_MEAS) { meas = MIN_MEAS; } */ PG_RETURN_FLOAT8(meas); } PG_FUNCTION_INFO_V1(lrs_getgeom_pointoffset); Datum lrs_getgeom_pointoffset(PG_FUNCTION_ARGS) { PG_LWGEOM *geom1 = (PG_LWGEOM *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0)); double m = PG_GETARG_FLOAT8(1); double angle = PG_GETARG_FLOAT8(2); // angle of offset in degrees double d = PG_GETARG_FLOAT8(3); // distance for offset LWGEOM *res = NULL; LWPOINT *lwpnt = NULL; POINT2D *pt = NULL; #ifdef DEBUG elog(INFO,"start getting point from measures with offset"); #endif /* xNew = x + cos(angle * 3.14 / 180) * d; yNew = y + sin(angle * 3.14 / 180) * d; */ if (!((TYPE_GETTYPE(geom1->type) == LINETYPE) || (TYPE_GETTYPE(geom1->type) == MULTILINETYPE) )) { elog(ERROR,"getmeasure: 1st arg isnt of type LINE or MULTILINE"); PG_RETURN_FLOAT8(-1); } if (!TYPE_HASM(geom1->type)) { elog(ERROR,"getmeasure: line doesn't contain any M values"); PG_RETURN_FLOAT8(-1); } if (TYPE_NDIMS(geom1->type) == 4) { res = getLRSPoint4D(SERIALIZED_FORM(geom1), m); } else { res = getLRSPoint3DM(SERIALIZED_FORM(geom1), m); } #ifdef DEBUG elog(INFO,"calculate offset point"); #endif lwpnt = lwgeom_as_lwpoint(res); lwpnt = (LWPOINT*) res; if (lwpnt == NULL) { elog(ERROR,"can't cast geom to lwpoint"); PG_RETURN_FLOAT8(-1); } #ifdef DEBUG elog(INFO,"test1"); #endif pt = (POINT2D*)getPoint_internal(lwpnt->point,0); #ifdef DEBUG elog(INFO,"test2"); #endif pt->x = pt->x + cos(angle * 3.14 / 180) * d; #ifdef DEBUG elog(INFO,"test3"); #endif pt->y = pt->y + sin(angle * 3.14 / 180) * d; #ifdef DEBUG elog(INFO,"test4"); #endif PG_RETURN_POINTER(pglwgeom_serialize(res)); } PG_FUNCTION_INFO_V1(lrs_getgeom_point); Datum lrs_getgeom_point(PG_FUNCTION_ARGS) { /* * Purpose: calculates a point along the line given a M value. In case the M * value is smaller/larger than the first/last point in the LRS, the first/last point * is returned. * Input: * 1. LRS system (type Geometry) * 2. M value (type Double) * Output: * Point at M (type Geometry) * Author: Christoph Spoerri * Calls: * getLRSPoint4D * getLRSPoint3DM */ PG_LWGEOM *geom1 = (PG_LWGEOM *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0)); double m = PG_GETARG_FLOAT8(1); LWGEOM *res = NULL; if (!((TYPE_GETTYPE(geom1->type) == LINETYPE) || (TYPE_GETTYPE(geom1->type) == MULTILINETYPE) )) { elog(ERROR,"getmeasure: 1st arg isnt of type LINE or MULTILINE"); PG_RETURN_FLOAT8(-1); } if (!TYPE_HASM(geom1->type)) { elog(ERROR,"getmeasure: line doesn't contain any M values"); PG_RETURN_FLOAT8(-1); } if (TYPE_NDIMS(geom1->type) == 4) { res = pglwgeom_serialize(getLRSPoint4D(SERIALIZED_FORM(geom1), m)); } else { res = pglwgeom_serialize(getLRSPoint3DM(SERIALIZED_FORM(geom1), m)); } PG_RETURN_POINTER(res); } LWGEOM *getLRSPoint4D (char *lwgeom1, double m) { /* * Purpose: returns a point along a 4D LRS based on a M value * Input: * 1. LRS geometry (type Geometry) * 2. M value (type double) * Output: * Point at location M (type Geometry/POINT) * Author: Christoph Spoerri * Called by: * lrs_getgeom_point */ int32 geomCount, endPtNum, ptCount, endIdx, incVar, invFlag; LWPOINT *pt1; LWGEOM *g1; double r; POINTARRAY *ipa; POINT4D frm, to; #ifdef DEBUG elog(INFO,"start getting point from measures"); #endif // determine the order of the measures, so that the points can be traversed accordently g1 = (LWGEOM *)lwgeom_deserialize(lwgeom_getsubgeometry(lwgeom1,0)); ipa = ((LWLINE *)g1)->points; frm = getPoint4d(ipa,0); to = getPoint4d(ipa,1); if (to.m > frm.m) { geomCount = 0; endIdx = lwgeom_getnumgeometries(lwgeom1); incVar = 1; invFlag = 0; if (frm.m > m) { elog(WARNING,"Measure value %f is out of range. Returning first point in line geometry.",m); return (LWGEOM *) make_lwpoint4d(-1, frm.x, frm.y, frm.z, frm.m); } } else { geomCount = lwgeom_getnumgeometries(lwgeom1)-1; endIdx = -1; incVar = -1; invFlag = 1; g1 = (LWGEOM *)lwgeom_deserialize(lwgeom_getsubgeometry(lwgeom1,0)); ipa = ((LWLINE *)g1)->points; to = getPoint4d(ipa,ipa->npoints-1); if (to.m > m) { elog(WARNING,"Measure value %f is out of range. Returning first point in line geometry.",m); return (LWGEOM *)make_lwpoint4d(-1, to.x, to.y, to.z, to.m); } } while (geomCount != endIdx) //for each object in geom1 { #ifdef DEBUG elog(INFO,"getting line %i",geomCount); #endif g1 = (LWGEOM *)lwgeom_deserialize(lwgeom_getsubgeometry(lwgeom1,geomCount)); ipa = ((LWLINE *)g1)->points; if (!invFlag) { ptCount = 0; endPtNum = ipa->npoints; } else { ptCount = ipa->npoints - 1; endPtNum = -1; } frm = getPoint4d(ipa,ptCount); ptCount = ptCount + incVar; while (ptCount != endPtNum) { to = getPoint4d(ipa,ptCount); if (to.m >= m) { #ifdef DEBUG elog(INFO,"return point: m = %f",m); elog(INFO,"return point: frm.m = %f : to.m = %f",frm.m, to.m); #endif r = (m - frm.m) / (to.m - frm.m); // SRID needs to be corrected !!! #ifdef DEBUG elog(INFO,"return point: r = %f",r); elog(INFO,"return point: x = %f : y = %f",frm.x + (to.x - frm.x)*r, frm.y + (to.y - frm.y)*r); #endif pt1 = make_lwpoint4d(-1, frm.x + (to.x - frm.x)*r, frm.y + (to.y - frm.y)*r, 0, m); return (LWGEOM *)pt1; } frm = to; ptCount = ptCount + incVar; } geomCount = geomCount + incVar; } #ifdef DEBUG elog(INFO,"didn't find a point"); #endif // looks like measure of more than existing measures in line // return the last point in the line elog(WARNING,"Measure value %f is out of range. Returning last point in line geometry.",m); return (LWGEOM *)make_lwpoint4d(-1, to.x, to.y, to.z, to.m); } LWGEOM *getLRSPoint3DM (char *lwgeom1, double m) { /* * Purpose: returns a point along a 3DM LRS based on a M value * Input: * 1. LRS geometry (type Geometry) * 2. M value (type double) * Output: * Point at location M (type Geometry/POINT) * Author: Christoph Spoerri * Called by: * lrs_getgeom_point */ int32 geomCount, endPtNum, ptCount, endIdx, incVar, invFlag; LWPOINT *pt1; LWGEOM *g1; double r; POINTARRAY *ipa; POINT3DM frm, to; #ifdef DEBUG elog(INFO,"start getting point from measures (3d)"); #endif // determine the order of the measures, so that the points can be traversed accordently g1 = (LWGEOM *)lwgeom_deserialize(lwgeom_getsubgeometry(lwgeom1,0)); ipa = ((LWLINE *)g1)->points; frm = getPoint3dm(ipa,0); to = getPoint3dm(ipa,1); if (to.m > frm.m) { geomCount = 0; endIdx = lwgeom_getnumgeometries(lwgeom1); incVar = 1; invFlag = 0; if (frm.m > m) { elog(WARNING,"Measure value %f is out of range. Returning first point in line geometry.",m); return lwpoint_as_lwgeom(make_lwpoint3dm(-1, frm.x, frm.y, frm.m)); } } else { geomCount = lwgeom_getnumgeometries(lwgeom1)-1; endIdx = -1; incVar = -1; invFlag = 1; g1 = (LWGEOM *)lwgeom_deserialize(lwgeom_getsubgeometry(lwgeom1,0)); ipa = ((LWLINE *)g1)->points; to = getPoint3dm(ipa,ipa->npoints-1); if (to.m > m) { elog(WARNING,"Measure value %f is out of range. Returning first point in line geometry.",m); return lwpoint_as_lwgeom(make_lwpoint3dm(-1, to.x, to.y, to.m)); } } while (geomCount != endIdx) //for each object in geom1 { #ifdef DEBUG elog(INFO,"getting line %i",geomCount); #endif g1 = (LWGEOM *)lwgeom_deserialize(lwgeom_getsubgeometry(lwgeom1,geomCount)); ipa = ((LWLINE *)g1)->points; if (!invFlag) { ptCount = 0; endPtNum = ipa->npoints; } else { ptCount = ipa->npoints - 1; endPtNum = -1; } frm = getPoint3dm(ipa,ptCount); ptCount = ptCount + incVar; while (ptCount != endPtNum) { // frm = getPoint3dm(ipa,0); // for (j=1; j < ipa->npoints; j++) { to = getPoint3dm(ipa,ptCount); #ifdef DEBUG elog(INFO,"to.m = %f",to.m); #endif if (to.m >= m) { #ifdef DEBUG elog(NOTICE,"found a point"); #endif #ifdef DEBUG elog(NOTICE,"return point: m = %f",m); elog(NOTICE,"return point: frm.m = %f : to.m = %f",frm.m, to.m); #endif r = (m - frm.m) / (to.m - frm.m); // SRID needs to be corrected !!! pt1 = make_lwpoint3dm(-1, frm.x + (to.x - frm.x)*r, frm.y + (to.y - frm.y)*r, m); #ifdef DEBUG printLWPOINT(pt1); #endif return lwpoint_as_lwgeom(pt1); } frm = to; ptCount = ptCount + incVar; } geomCount = geomCount + incVar; } #ifdef DEBUG elog(INFO,"didn't find a point"); #endif // looks like measure of more than existing measures in line // return the last point in the line elog(WARNING,"Measure value %f is out of range. Returning last point in line geometry.",m); return lwpoint_as_lwgeom(make_lwpoint3dm(-1, to.x, to.y, to.m)); } PG_FUNCTION_INFO_V1(lrs_getgeom_line); Datum lrs_getgeom_line(PG_FUNCTION_ARGS) { /* * Purpose: calculates a line LRS system based on the two input M values. The first * M value is the start of the new line, the second M value is the end point. In * case of a MULTILINESTRING the lines are appended as necessary. If either M value * is smaller/large than the first/last point of the LRS, the first/last point will * be used * Input: * 1. LRS geometry (type GEOMETRY) * 2. first M value (type DOUBLE) * 3. second M value (type DOUBLE) * Output: * new LRS system (type GEOMETRY) * Author: Christoph Spoerri * Calls: * getLRSLine4D * getLRSLine3DM */ PG_LWGEOM *geom1 = (PG_LWGEOM *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0)); double frm_m = PG_GETARG_FLOAT8(1); double to_m = PG_GETARG_FLOAT8(2); int has4Dim; LWGEOM *res = NULL; char *lwgeom1; char isReversed; // used to indicate with LRS system increases (0) or decreases (1 - reversed) LWGEOM_INSPECTED *inspected; LWLINE *startline = NULL, *endLine = NULL; #ifdef DEBUG elog(INFO,"lsr_getgeom_line: start"); #endif if (!((TYPE_GETTYPE(geom1->type) == LINETYPE) || (TYPE_GETTYPE(geom1->type) == MULTILINETYPE) )) { elog(ERROR,"getmeasure: 1st arg isnt of type LINE or MULTILINE"); PG_RETURN_FLOAT8(-1); } if (!TYPE_HASM(geom1->type)) { elog(ERROR,"getmeasure: line doesn't contain any M values"); PG_RETURN_FLOAT8(-1); } has4Dim = TYPE_NDIMS(geom1->type)-3; if (frm_m == to_m) { if (has4Dim) { res = lwgeom_serialize(getLRSPoint4D(SERIALIZED_FORM(geom1), frm_m)); } else { res = lwgeom_serialize(getLRSPoint3DM(SERIALIZED_FORM(geom1), frm_m)); } elog(WARNING,"From and To measure the same. Returning a single point."); PG_RETURN_POINTER(res); } lwgeom1 = SERIALIZED_FORM(geom1); inspected = lwgeom_inspect(lwgeom1); startline = lwgeom_getline(lwgeom1, 0); if ( startline == NULL ) { PG_FREE_IF_COPY(geom1, 0); PG_RETURN_NULL(); } endLine = lwgeom_getline_inspected(inspected,inspected->ngeometries-1); if ( endLine == NULL ) { PG_FREE_IF_COPY(geom1, 0); PG_RETURN_NULL(); } pfree_inspected(inspected); if (has4Dim) { #ifdef DEBUG elog(INFO,"4D geom"); #endif POINT4D *startPt = NULL, *endPt = NULL; startPt = (POINT4D *)getPoint_internal(startline->points, 0); endPt = (POINT4D *)getPoint_internal(endLine->points, endLine->points->npoints-1); /* in case the from and to measure don't 'point' the same way as the measure of the LRS we want to switch the input measures in order to get a correct result */ #ifdef DEBUG elog(INFO,"to_m = %f, frm_m = %f",to_m, frm_m); #endif isReversed = 0; if (endPt->m < startPt->m) isReversed = 1; if (((to_m > frm_m) && (endPt->m > startPt->m)) || ((to_m < frm_m) && (endPt->m < startPt->m))) { res = getLRSLine4D(lwgeom1, frm_m, to_m, isReversed); } else { #ifdef DEBUG elog(INFO,"reversed"); #endif res = getLRSLine4D(lwgeom1, to_m, frm_m, isReversed); } #ifdef DEBUG elog(INFO,"get new line"); #endif } else { #ifdef DEBUG elog(INFO,"3D geom"); #endif POINT3DM *startPt = NULL; POINT3DM *endPt = NULL; startPt = (POINT3DM *)getPoint_internal(startline->points, 0); endPt = (POINT3DM *)getPoint_internal(endLine->points, endLine->points->npoints-1); /* in case the from and to measure don't 'point' the same way as the measure of the LRS we want to switch the given measures in order to get a correct result */ isReversed = 0; if (endPt->m < startPt->m) isReversed = 1; if (((to_m > frm_m) && (endPt->m > startPt->m)) || ((to_m < frm_m) && (endPt->m < startPt->m))) { res = getLRSLine3DM(lwgeom1, frm_m, to_m, isReversed); } else { res = getLRSLine3DM(lwgeom1, to_m, frm_m, isReversed); } } #ifdef DEBUG elog(INFO,"checking for null"); #endif if (res != NULL) { #ifdef DEBUG elog(INFO,"return the geom"); #endif PG_RETURN_POINTER(pglwgeom_serialize(res)); } elog(WARNING,"getlrsline: no shape could be calculated based on input data (from = %f, to = %f.)",frm_m,to_m); PG_RETURN_NULL(); } LWGEOM *getLRSLine4D(const char *lwgeom1, double frm_m, double to_m, char isReversed) { /* * Purpose: returns a new LRS along a 4D LRS based on two M values * Input: * 1. LRS geometry (type Geometry) * 2. first M value (type DOUBLE) * 3. second M value (type DOUBLE) * Output: * new LRS system from first M to second M (type Geometry/MULTILINESTRING) * Author: Christoph Spoerri * Called by: * lrs_getgeom_line */ char found; int32 type1,j,i; int ptsize, npts; LWGEOM *g1; double r; POINTARRAY *ipa, *opa; POINT4D frm, to, pt1; LWMLINE *res = NULL; LWLINE *resLine = NULL; POINT4D pt4d; char *opts = NULL; double **t = NULL, **q = NULL; //pointer to the pointer of the provided M values double *user_m = NULL, *point_m = NULL; //pointer to the actuall M values LWGEOM *whatLine; // used for creating MultiLine LWMLINE *toLine; // used for creating MultiLine LWGEOM **geomColl; int numGeom=0; // set the pointers to the initial value address point_m = &frm_m; user_m = &frm_m; if (isReversed) { q = &point_m; t = &user_m; } else { t = &point_m; q = &user_m; } found = false; geomColl = lwalloc(sizeof(LWGEOM *)*lwgeom_getnumgeometries(lwgeom1)); for (i=0; i < lwgeom_getnumgeometries(lwgeom1); i++) //for each object in geom1 { #ifdef DEBUG elog(INFO,"getting line %i",i); #endif g1 = (LWGEOM *)lwgeom_deserialize(lwgeom_getsubgeometry(lwgeom1,i)); ipa = ((LWLINE *)g1)->points; //create pointarray to hold new points ptsize = pointArray_ptsize(ipa); // in case we have a multipart line, additional parts will be // appended if (i != 0) { // append new points/line if necessary to existing line as multiline // if (resLine == NULL) if (numGeom) { opa = pointArray_construct(opts, 1, 1, npts); // whatLine = (LWGEOM*) lwline_construct(-1, NULL, opa); // toLine = (LWMLINE *) resLine; resLine = lwline_construct(-1, NULL, opa); //destroy to and what, since they are not used anymore // lwgeom_release((LWGEOM*) toLine); // lwgeom_release((LWGEOM*) whatLine); geomColl[numGeom++]=(LWGEOM*)resLine; } else { if (npts > 0) { opa = pointArray_construct(opts, 1, 1, npts); resLine = lwline_construct(-1, NULL, opa); geomColl[numGeom++]=(LWGEOM*)resLine; } } } // create new point array to continue //opa = palloc(sizeof(POINTARRAY)); opts = (char *)palloc(ptsize*ipa->npoints); npts = 0; frm = getPoint4d(ipa,0); // in case the first point is already found, we need to automatically add // the first point of the new linesegment if (found) { // adding the points for the current line memcpy(opts + ptsize * npts, &frm, ptsize); npts++; } for (j=1; j < ipa->npoints; j++) { to = getPoint4d(ipa,j); if (found) { point_m = &to.m; user_m = &to_m; if (**t >= **q) { r = (to_m - frm.m) / (to.m - frm.m); pt1.x = frm.x + (to.x - frm.x)*r; pt1.y = frm.y + (to.y - frm.y)*r; pt1.z = 0; pt1.m = to_m; #ifdef DEBUG elog(INFO,"lrs_getgeom_line: last point %i (m = %f)",npts, to.m); #endif // adding the new point memcpy(opts + ptsize * npts, &pt1, ptsize); npts++; // adjusting the size of the 'point' memory to the exact # of points // opts = (char *)repalloc(opts, ptsize*npts); opa=pointArray_construct(opts,1,1,npts); resLine = lwline_construct(-1, NULL, opa); // if (!numGeom) // { // return (LWGEOM*)resLine; // } geomColl[numGeom++]=resLine; return lwcollection_as_lwgeom(lwcollection_construct(MULTILINETYPE, -1,NULL, numGeom, geomColl)); } else { // adding the points for the current line memcpy(opts + ptsize * npts, &to, ptsize); npts++; #ifdef DEBUG elog(INFO,"lrs_getgeom_line: next point %i (m = %f)",npts, to.m); #endif } } else { point_m = &to.m; user_m = &frm_m; if (**t >= **q) { r = (frm_m - frm.m) / (to.m - frm.m); pt4d.x = frm.x + (to.x - frm.x)*r; pt4d.y = frm.y + (to.y - frm.y)*r; pt4d.m = frm_m; pt4d.z = 0; // adding the new point memcpy(opts + ptsize * npts, &pt4d, ptsize); npts++; // adding the 'to' point memcpy(opts + ptsize * npts, &to, ptsize); npts++; found = true; #ifdef DEBUG elog(INFO,"lrs_getgeom_line: found first point"); #endif point_m = &to.m; user_m = &to_m; if (**t >= **q) { r = (to_m - frm.m) / (to.m - frm.m); pt1.x = frm.x + (to.x - frm.x)*r; pt1.y = frm.y + (to.y - frm.y)*r; pt1.z = 0; pt1.m = to_m; #ifdef DEBUG elog(INFO,"lrs_getgeom_line: last point %i",npts); #endif // adding the new point memcpy(opts + ptsize * npts, &to, ptsize); npts++; // adjusting the size of the 'point' memory to the exact # of points // opts = (char *)repalloc(opts, ptsize*npts); opa=pointArray_construct(opts,1,1,npts); resLine = lwline_construct(-1, NULL, opa); // if (!numGeom) // { return (LWGEOM*)resLine; // } // geomColl[numGeom++]=resLine; // return lwcollection_as_lwgeom(lwcollection_construct(MULTILINETYPE, -1,NULL, numGeom, geomColl)); } } } frm = to; } #ifdef DEBUG elog(INFO,"lrs_getgeom_line: checking next line seg."); #endif } // in case we found the first point but not the second/end point, we just assume that the last point of the // LRS line will end where the 'parent' line ends if (found) { elog(WARNING,"getlrsline: could not find last point (%f), assume end point of LRS line",to_m); pt1.x = frm.x; pt1.y = frm.y; pt1.z = 0; pt1.m = frm.m; // adding the new point memcpy(opts + ptsize * npts, &to, ptsize); npts++; // adjusting the size of the 'point' memory to the exact # of points // opts = (char *)repalloc(opts, ptsize*npts); opa=pointArray_construct(opts,1,1,npts); resLine = lwline_construct(-1, NULL, opa); // if (!numGeom) // { // return (LWGEOM*)resLine; // } geomColl[numGeom++]=resLine; return lwcollection_as_lwgeom(lwcollection_construct(MULTILINETYPE, -1,NULL, numGeom, geomColl)); } return NULL; } LWGEOM *getLRSLine3DM(const char *lwgeom1, double frm_m, double to_m, char isReversed) { /* * Purpose: returns a new LRS along a 3DM LRS based on two M values * Input: * 1. LRS geometry (type Geometry) * 2. first M value (type DOUBLE) * 3. second M value (type DOUBLE) * Output: * new LRS system from first M to second M (type Geometry/MULTILINESTRING) * Author: Christoph Spoerri * Called by: * lrs_getgeom_line */ char found; int32 type1,j,i; int ptsize, npts; LWGEOM *g1; double r; POINTARRAY *ipa, *opa; POINT3DM frm, to, pt1; LWMLINE *res = NULL; LWGEOM *resLine = NULL; char *opts = NULL; double **t = NULL, **q = NULL; //pointer to the pointer of the values double *user_m = NULL, *point_m = NULL; //pointer to the actuall values to compare LWGEOM *whatLine; // used for creating MultiLine LWMLINE *toLine; // used for creating MultiLine LWGEOM **geomColl; int numGeom=0; point_m = &frm_m; user_m = &frm_m; if (isReversed) { q = &point_m; t = &user_m; } else { t = &point_m; q = &user_m; } found = false; geomColl = lwalloc(sizeof(LWGEOM *)*lwgeom_getnumgeometries(lwgeom1)); for (i=0; i< lwgeom_getnumgeometries(lwgeom1); i++) //for each object in geom1 { #ifdef DEBUG elog(INFO,"getting line %i",i); #endif g1 = (LWGEOM *)lwgeom_deserialize(lwgeom_getsubgeometry(lwgeom1,i)); ipa = ((LWLINE *)g1)->points; //create pointarray to hold new points ptsize = pointArray_ptsize(ipa); // in case we have a multipart line, additional parts will be // appended if (i != 0) { // append new points/line if necessary to existing line as multiline // if (resLine == NULL) if (numGeom) { opa = pointArray_construct(opts, 1, 1, npts); // whatLine = (LWGEOM*) lwline_construct(-1, NULL, opa); // toLine = (LWMLINE *) resLine; resLine = lwline_construct(-1, NULL, opa); //destroy to and what, since they are not used anymore // lwgeom_release((LWGEOM*) toLine); // lwgeom_release((LWGEOM*) whatLine); geomColl[numGeom++]=(LWGEOM*)resLine; } else { if (npts > 0) { opa = pointArray_construct(opts, 1, 1, npts); resLine = lwline_construct(-1, NULL, opa); geomColl[numGeom++]=(LWGEOM*)resLine; } } } // create new point array to continue //opa = palloc(sizeof(POINTARRAY)); opts = (char *)palloc(ptsize*ipa->npoints); npts = 0; frm = getPoint3dm(ipa,0); // in case the first point is already found, we need to automatically add // the first point of the new linesegment if (found) { // adding the points for the current line memcpy(opts + ptsize * npts, &frm, ptsize); npts++; } for (j=1; j < ipa->npoints; j++) { to = getPoint3dm(ipa,j); if (found) { point_m = &to.m; user_m = &to_m; if (**t >= **q) { r = (to_m - frm.m) / (to.m - frm.m); pt1.x = frm.x + (to.x - frm.x)*r; pt1.y = frm.y + (to.y - frm.y)*r; pt1.m = to_m; #ifdef DEBUG elog(INFO,"lrs_getgeom_line: last point %i",npts); #endif // adding the new point memcpy(opts + ptsize * npts, &pt1, ptsize); npts++; // adjusting the size of the 'point' memory to the exact # of points // opts = (char *)repalloc(opts, ptsize*npts); opa=pointArray_construct(opts,0,1,npts); resLine = lwline_construct(-1, NULL, opa); // if (!numGeom) // { // return (LWGEOM*)resLine; // } geomColl[numGeom++]=resLine; return lwcollection_as_lwgeom(lwcollection_construct(MULTILINETYPE, -1,NULL, numGeom, geomColl)); } else { // adding the points for the current line memcpy(opts + ptsize * npts, &to, ptsize); npts++; #ifdef DEBUG elog(INFO,"lrs_getgeom_line: next point %i (m = %f)",npts, to.m); #endif } } else { point_m = &to.m; user_m = &frm_m; if (**t >= **q) { r = (frm_m - frm.m) / (to.m - frm.m); pt1.x = frm.x + (to.x - frm.x)*r; pt1.y = frm.y + (to.y - frm.y)*r; pt1.m = frm_m; // adding the new point memcpy(opts + ptsize * npts, &pt1, ptsize); npts++; // adding the 'to' point memcpy(opts + ptsize * npts, &to, ptsize); npts++; found = true; #ifdef DEBUG elog(INFO,"lrs_getgeom_line: found first point"); #endif point_m = &to.m; user_m = &frm_m; if (**t >= **q) { r = (to_m - frm.m) / (to.m - frm.m); pt1.x = frm.x + (to.x - frm.x)*r; pt1.y = frm.y + (to.y - frm.y)*r; pt1.m = to_m; #ifdef DEBUG elog(INFO,"lrs_getgeom_line: last point %i",npts); #endif // adding the new point memcpy(opts + ptsize * npts, &to, ptsize); npts++; // adjusting the size of the 'point' memory to the exact # of points // opts = (char *)repalloc(opts, ptsize*npts); opa=pointArray_construct(opts,0,1,npts); resLine = lwline_construct(-1, NULL, opa); // if (!numGeom) // { return (LWGEOM*)resLine; // } // geomColl[numGeom++]=resLine; // return lwcollection_as_lwgeom(lwcollection_construct(MULTILINETYPE, -1,NULL, numGeom, geomColl)); } } } frm = to; } } // in case we found the first point but not the second/end point, we just assume that the last point of the // LRS line will end where the 'parent' line ends if (found) { elog(WARNING,"getlrsline: could not find last point (%f), assume end point of LRS line",to_m); pt1.x = frm.x; pt1.y = frm.y; pt1.m = frm.m; // adding the new point memcpy(opts + ptsize * npts, &to, ptsize); npts++; // adjusting the size of the 'point' memory to the exact # of points // opts = (char *)repalloc(opts, ptsize*npts); opa=pointArray_construct(opts,0,1,npts); resLine = lwline_construct(-1, NULL, opa); // if (!numGeom) // { // return (LWGEOM*)resLine; // } geomColl[numGeom++]=resLine; return lwcollection_as_lwgeom(lwcollection_construct(MULTILINETYPE, -1,NULL, numGeom, geomColl)); } return NULL; } LWGEOM *addLastPoint(char *opts , int npts, char *pt1, int ptsize, int hasZ, LWGEOM *resLine) /* Purpose: adds the list of collected points and the last point (pt1) to the already collected line. It's basically used to finish up the creation of the new geometry Input: - *opts: the current list of points - npts: the current number of points - *pt1: the last point to be added - hasZ: specifies if the collection has a Z value or not (Note: it always has a M value) - *resLine: the current line Output: - the final line (or NULL on error) Called by: - getLRSLine3DM - getLRSLine4D */ { POINTARRAY *opa; LWMLINE *to; LWGEOM *what; // adding the last point to the 'point' memory memcpy(opts + ptsize * npts, pt1, ptsize); npts++; // adjusting the size of the 'point' memory to the exact # of points opts = (char *)repalloc(opts, ptsize*npts); if ( opts == NULL ) { elog(ERROR, "Out of virtual memory"); return NULL; } opa = pointArray_construct(opts, hasZ, 1, npts); what = (LWGEOM*) lwline_construct(-1, NULL, opa); if (resLine == NULL) { //resLine = (LWGEOM*) what; //resLine = lwcollection_as_lwgeom( lwcollection_construct(MULTILINETYPE,-1,NULL,1, &what)); //return (LWGEOM*) lwcollection_construct(MULTILINETYPE,-1,NULL,1, &what); } else { to = (LWMLINE *) resLine; resLine = lwmline_add(to,-1,what); //destroy to and what, since they are not used anymore lwgeom_release((LWGEOM*)to); lwgeom_release((LWGEOM*)what); } return resLine; }