PostGIS  2.1.10dev-r@@SVN_REVISION@@
lwgeom_functions_lrs.c
Go to the documentation of this file.
1 /**********************************************************************
2  *
3  * PostGIS - Spatial Types for PostgreSQL
4  * http://postgis.net
5  *
6  * Copyright (C) 2001-2005 Refractions Research Inc.
7  *
8  * This is free software; you can redistribute and/or modify it under
9  * the terms of the GNU General Public Licence. See the COPYING file.
10  *
11  **********************************************************************/
12 
13 #include <math.h>
14 
15 #include "postgres.h"
16 #include "fmgr.h"
17 
18 #include "../postgis_config.h"
19 #include "liblwgeom.h"
20 #include "lwgeom_pg.h"
21 
22 /*
23 * Add a measure dimension to a line, interpolating linearly from the
24 * start value to the end value.
25 * ST_AddMeasure(Geometry, StartMeasure, EndMeasure) returns Geometry
26 */
27 Datum ST_AddMeasure(PG_FUNCTION_ARGS);
29 Datum ST_AddMeasure(PG_FUNCTION_ARGS)
30 {
31  GSERIALIZED *gin = (GSERIALIZED *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
32  GSERIALIZED *gout;
33  double start_measure = PG_GETARG_FLOAT8(1);
34  double end_measure = PG_GETARG_FLOAT8(2);
35  LWGEOM *lwin, *lwout;
36  int type = gserialized_get_type(gin);
37 
38  /* Raise an error if input is not a linestring or multilinestring */
39  if ( type != LINETYPE && type != MULTILINETYPE )
40  {
41  lwerror("Only LINESTRING and MULTILINESTRING are supported");
42  PG_RETURN_NULL();
43  }
44 
45  lwin = lwgeom_from_gserialized(gin);
46  if ( type == LINETYPE )
47  lwout = (LWGEOM*)lwline_measured_from_lwline((LWLINE*)lwin, start_measure, end_measure);
48  else
49  lwout = (LWGEOM*)lwmline_measured_from_lwmline((LWMLINE*)lwin, start_measure, end_measure);
50 
51  lwgeom_free(lwin);
52 
53  if ( lwout == NULL )
54  PG_RETURN_NULL();
55 
56  gout = geometry_serialize(lwout);
57  lwgeom_free(lwout);
58 
59  PG_RETURN_POINTER(gout);
60 }
61 
62 
63 /*
64 * Locate a point along a feature based on a measure value.
65 * ST_LocateAlong(Geometry, Measure, [Offset]) returns Geometry
66 */
67 Datum ST_LocateAlong(PG_FUNCTION_ARGS);
69 Datum ST_LocateAlong(PG_FUNCTION_ARGS)
70 {
71  GSERIALIZED *gin = (GSERIALIZED *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
72  GSERIALIZED *gout;
73  LWGEOM *lwin = NULL, *lwout = NULL;
74  double measure = PG_GETARG_FLOAT8(1);
75  double offset = PG_GETARG_FLOAT8(2);;
76 
77  lwin = lwgeom_from_gserialized(gin);
78  lwout = lwgeom_locate_along(lwin, measure, offset);
79  lwgeom_free(lwin);
80  PG_FREE_IF_COPY(gin, 0);
81 
82  if ( ! lwout )
83  PG_RETURN_NULL();
84 
85  gout = geometry_serialize(lwout);
86  lwgeom_free(lwout);
87 
88  PG_RETURN_POINTER(gout);
89 }
90 
91 
92 /*
93 * Locate the portion of a line between the specified measures
94 */
95 Datum ST_LocateBetween(PG_FUNCTION_ARGS);
97 Datum ST_LocateBetween(PG_FUNCTION_ARGS)
98 {
99  GSERIALIZED *geom_in = (GSERIALIZED *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
100  double from = PG_GETARG_FLOAT8(1);
101  double to = PG_GETARG_FLOAT8(2);
102  double offset = PG_GETARG_FLOAT8(3);
103  LWCOLLECTION *geom_out = NULL;
104  LWGEOM *line_in = NULL;
105  static char ordinate = 'M'; /* M */
106 
107  if ( ! gserialized_has_m(geom_in) )
108  {
109  elog(ERROR,"This function only accepts geometries that have an M dimension.");
110  PG_RETURN_NULL();
111  }
112 
113  /* This should be a call to ST_LocateAlong! */
114  if ( to == from )
115  {
116  PG_RETURN_DATUM(DirectFunctionCall3(ST_LocateAlong, PG_GETARG_DATUM(0), PG_GETARG_DATUM(1), PG_GETARG_DATUM(3)));
117  }
118 
119  line_in = lwgeom_from_gserialized(geom_in);
120  geom_out = lwgeom_clip_to_ordinate_range(line_in, ordinate, from, to, offset);
121  lwgeom_free(line_in);
122  PG_FREE_IF_COPY(geom_in, 0);
123 
124  if ( ! geom_out )
125  {
126  elog(ERROR,"lwline_clip_to_ordinate_range returned null");
127  PG_RETURN_NULL();
128  }
129 
130  PG_RETURN_POINTER(geometry_serialize((LWGEOM*)geom_out));
131 }
132 
133 /*
134 * Locate the portion of a line between the specified elevations
135 */
136 Datum ST_LocateBetweenElevations(PG_FUNCTION_ARGS);
138 Datum ST_LocateBetweenElevations(PG_FUNCTION_ARGS)
139 {
140  GSERIALIZED *geom_in = (GSERIALIZED *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
141  double from = PG_GETARG_FLOAT8(1);
142  double to = PG_GETARG_FLOAT8(2);
143  LWCOLLECTION *geom_out = NULL;
144  LWGEOM *line_in = NULL;
145  static char ordinate = 'Z'; /* Z */
146  static double offset = 0.0;
147 
148  if ( ! gserialized_has_z(geom_in) )
149  {
150  elog(ERROR,"This function only accepts LINESTRING or MULTILINESTRING with Z dimensions.");
151  PG_RETURN_NULL();
152  }
153 
154  line_in = lwgeom_from_gserialized(geom_in);
155  geom_out = lwgeom_clip_to_ordinate_range(line_in, ordinate, from, to, offset);
156  lwgeom_free(line_in);
157  PG_FREE_IF_COPY(geom_in, 0);
158 
159  if ( ! geom_out )
160  {
161  elog(ERROR,"lwline_clip_to_ordinate_range returned null");
162  PG_RETURN_NULL();
163  }
164 
165  PG_RETURN_POINTER(geometry_serialize((LWGEOM*)geom_out));
166 }
167 
168 
169 Datum ST_InterpolatePoint(PG_FUNCTION_ARGS);
171 Datum ST_InterpolatePoint(PG_FUNCTION_ARGS)
172 {
173  GSERIALIZED *gser_line = (GSERIALIZED *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
174  GSERIALIZED *gser_point = (GSERIALIZED *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
175  LWGEOM *lwline;
176  LWPOINT *lwpoint;
177 
178  if ( gserialized_get_type(gser_line) != LINETYPE )
179  {
180  elog(ERROR,"ST_InterpolatePoint: 1st argument isn't a line");
181  PG_RETURN_NULL();
182  }
183  if ( gserialized_get_type(gser_point) != POINTTYPE )
184  {
185  elog(ERROR,"ST_InterpolatePoint: 2st argument isn't a point");
186  PG_RETURN_NULL();
187  }
188  if ( gserialized_get_srid(gser_line) != gserialized_get_srid(gser_point) )
189  {
190  elog(ERROR, "Operation on two geometries with different SRIDs");
191  PG_RETURN_NULL();
192  }
193  if ( ! gserialized_has_m(gser_line) )
194  {
195  elog(ERROR,"ST_InterpolatePoint only accepts geometries that have an M dimension");
196  PG_RETURN_NULL();
197  }
198 
199  lwpoint = lwgeom_as_lwpoint(lwgeom_from_gserialized(gser_point));
200  lwline = lwgeom_from_gserialized(gser_line);
201 
202  PG_RETURN_FLOAT8(lwgeom_interpolate_point(lwline, lwpoint));
203 }
204 
205 
206 Datum LWGEOM_line_locate_point(PG_FUNCTION_ARGS);
208 Datum LWGEOM_line_locate_point(PG_FUNCTION_ARGS)
209 {
210  GSERIALIZED *geom1 = (GSERIALIZED *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
211  GSERIALIZED *geom2 = (GSERIALIZED *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
212  LWLINE *lwline;
213  LWPOINT *lwpoint;
214  POINTARRAY *pa;
215  POINT4D p, p_proj;
216  double ret;
217 
218  if ( gserialized_get_type(geom1) != LINETYPE )
219  {
220  elog(ERROR,"line_locate_point: 1st arg isnt a line");
221  PG_RETURN_NULL();
222  }
223  if ( gserialized_get_type(geom2) != POINTTYPE )
224  {
225  elog(ERROR,"line_locate_point: 2st arg isnt a point");
226  PG_RETURN_NULL();
227  }
228  if ( gserialized_get_srid(geom1) != gserialized_get_srid(geom2) )
229  {
230  elog(ERROR, "Operation on two geometries with different SRIDs");
231  PG_RETURN_NULL();
232  }
233 
234  lwline = lwgeom_as_lwline(lwgeom_from_gserialized(geom1));
235  lwpoint = lwgeom_as_lwpoint(lwgeom_from_gserialized(geom2));
236 
237  pa = lwline->points;
238  lwpoint_getPoint4d_p(lwpoint, &p);
239 
240  ret = ptarray_locate_point(pa, &p, NULL, &p_proj);
241 
242  PG_RETURN_FLOAT8(ret);
243 }
244 
245 
246 /***********************************************************************
247 * LEGACY SUPPORT FOR locate_between_measures and locate_along_measure
248 * Deprecated at PostGIS 2.0. To be removed.
249 */
250 
251 
252 typedef struct
253 {
256 }
258 
260  POINTARRAY *ipa, double m0, double m1);
261 
263  LWCOLLECTION *lwcoll, double m0, double m1);
264 
266  LWGEOM *lwin, double m0, double m1);
267 
269  LWLINE *lwline_in, double m0, double m1);
270 
272  LWPOINT *lwpoint, double m0, double m1);
273 
274 static int clip_seg_by_m_range(
275  POINT4D *p1, POINT4D *p2, double m0, double m1);
276 
277 
278 /*
279  * Clip a segment by a range of measures.
280  * Z and M values are interpolated in case of clipping.
281  *
282  * Returns a bitfield, flags being:
283  * 0x0001 : segment intersects the range
284  * 0x0010 : first point is modified
285  * 0x0100 : second point is modified
286  *
287  * Values:
288  * - 0 segment fully outside the range, no modifications
289  * - 1 segment fully inside the range, no modifications
290  * - 7 segment crosses the range, both points modified.
291  * - 3 first point out, second in, first point modified
292  * - 5 first point in, second out, second point modified
293  */
294 static int
295 clip_seg_by_m_range(POINT4D *p1, POINT4D *p2, double m0, double m1)
296 {
297  double dM0, dM1, dX, dY, dZ;
298  POINT4D *tmp;
299  int swapped=0;
300  int ret=0;
301 
302  POSTGIS_DEBUGF(3, "m0: %g m1: %g", m0, m1);
303 
304  /* Handle corner case of m values being the same */
305  if ( p1->m == p2->m )
306  {
307  /* out of range, no clipping */
308  if ( p1->m < m0 || p1->m > m1 )
309  return 0;
310 
311  /* inside range, no clipping */
312  return 1;
313  }
314 
315  /*
316  * Order points so that p1 has the smaller M
317  */
318  if ( p1->m > p2->m )
319  {
320  tmp=p2;
321  p2=p1;
322  p1=tmp;
323  swapped=1;
324  }
325 
326  /*
327  * The M range is not intersected, segment
328  * fully out of range, no clipping.
329  */
330  if ( p2->m < m0 || p1->m > m1 )
331  return 0;
332 
333  /*
334  * The segment is fully inside the range,
335  * no clipping.
336  */
337  if ( p1->m >= m0 && p2->m <= m1 )
338  return 1;
339 
340  /*
341  * Segment intersects range, lets compute
342  * the proportional location of the two
343  * measures wrt p1/p2 m range.
344  *
345  * if p1 and p2 have the same measure
346  * this should never be reached (either
347  * both inside or both outside)
348  *
349  */
350  dM0=(m0-p1->m)/(p2->m-p1->m); /* delta-M0 */
351  dM1=(m1-p2->m)/(p2->m-p1->m); /* delta-M1 */
352  dX=p2->x-p1->x;
353  dY=p2->y-p1->y;
354  dZ=p2->z-p1->z;
355 
356  POSTGIS_DEBUGF(3, "dM0:%g dM1:%g", dM0, dM1);
357  POSTGIS_DEBUGF(3, "dX:%g dY:%g dZ:%g", dX, dY, dZ);
358  POSTGIS_DEBUGF(3, "swapped: %d", swapped);
359 
360  /*
361  * First point out of range, project
362  * it on the range
363  */
364  if ( p1->m < m0 )
365  {
366  /*
367  * To prevent rounding errors, then if m0==m1 and p2 lies within the range, copy
368  * p1 as a direct copy of p2
369  */
370  if (m0 == m1 && p2->m <= m1)
371  {
372  memcpy(p1, p2, sizeof(POINT4D));
373 
374  POSTGIS_DEBUG(3, "Projected p1 on range (as copy of p2)");
375  }
376  else
377  {
378  /* Otherwise interpolate coordinates */
379  p1->x += (dX*dM0);
380  p1->y += (dY*dM0);
381  p1->z += (dZ*dM0);
382  p1->m = m0;
383 
384  POSTGIS_DEBUG(3, "Projected p1 on range");
385  }
386 
387  if ( swapped ) ret |= 0x0100;
388  else ret |= 0x0010;
389  }
390 
391  /*
392  * Second point out of range, project
393  * it on the range
394  */
395  if ( p2->m > m1 )
396  {
397  /*
398  * To prevent rounding errors, then if m0==m1 and p1 lies within the range, copy
399  * p2 as a direct copy of p1
400  */
401  if (m0 == m1 && p1->m >= m0)
402  {
403  memcpy(p2, p1, sizeof(POINT4D));
404 
405  POSTGIS_DEBUG(3, "Projected p2 on range (as copy of p1)");
406  }
407  else
408  {
409  /* Otherwise interpolate coordinates */
410  p2->x += (dX*dM1);
411  p2->y += (dY*dM1);
412  p2->z += (dZ*dM1);
413  p2->m = m1;
414 
415  POSTGIS_DEBUG(3, "Projected p2 on range");
416  }
417 
418  if ( swapped ) ret |= 0x0010;
419  else ret |= 0x0100;
420  }
421 
422  /* Clipping occurred */
423  return ret;
424 
425 }
426 
427 static POINTARRAYSET
428 ptarray_locate_between_m(POINTARRAY *ipa, double m0, double m1)
429 {
430  POINTARRAYSET ret;
431  POINTARRAY *dpa=NULL;
432  int i;
433 
434  ret.nptarrays=0;
435 
436  /*
437  * We allocate space for as many pointarray as
438  * segments in the input POINTARRAY, as worst
439  * case is for each segment to cross the M range
440  * window.
441  * TODO: rework this to reduce used memory
442  */
443  ret.ptarrays=lwalloc(sizeof(POINTARRAY *)*ipa->npoints-1);
444 
445  POSTGIS_DEBUGF(2, "ptarray_locate...: called for pointarray %x, m0:%g, m1:%g",
446  ipa, m0, m1);
447 
448 
449  for (i=1; i<ipa->npoints; i++)
450  {
451  POINT4D p1, p2;
452  int clipval;
453 
454  getPoint4d_p(ipa, i-1, &p1);
455  getPoint4d_p(ipa, i, &p2);
456 
457  POSTGIS_DEBUGF(3, " segment %d-%d [ %g %g %g %g - %g %g %g %g ]",
458  i-1, i,
459  p1.x, p1.y, p1.z, p1.m,
460  p2.x, p2.y, p2.z, p2.m);
461 
462  clipval = clip_seg_by_m_range(&p1, &p2, m0, m1);
463 
464  /* segment completely outside, nothing to do */
465  if (! clipval ) continue;
466 
467  POSTGIS_DEBUGF(3, " clipped to: [ %g %g %g %g - %g %g %g %g ] clipval: %x", p1.x, p1.y, p1.z, p1.m,
468  p2.x, p2.y, p2.z, p2.m, clipval);
469 
470  /* If no points have been accumulated so far, then if clipval != 0 the first point must be the
471  start of the intersection */
472  if (dpa == NULL)
473  {
474  POSTGIS_DEBUGF(3, " 1 creating new POINTARRAY with first point %g,%g,%g,%g", p1.x, p1.y, p1.z, p1.m);
475 
477  ptarray_append_point(dpa, &p1, LW_TRUE);
478  }
479 
480  /* Otherwise always add the next point, avoiding duplicates */
481  if (dpa)
482  ptarray_append_point(dpa, &p2, LW_FALSE);
483 
484  /*
485  * second point has been clipped
486  */
487  if ( clipval & 0x0100 || i == ipa->npoints-1 )
488  {
489  POSTGIS_DEBUGF(3, " closing pointarray %x with %d points", dpa, dpa->npoints);
490 
491  ret.ptarrays[ret.nptarrays++] = dpa;
492  dpa = NULL;
493  }
494  }
495 
496  /*
497  * if dpa!=NULL it means we didn't close it yet.
498  * this should never happen.
499  */
500  if ( dpa != NULL ) lwerror("Something wrong with algorithm");
501 
502  return ret;
503 }
504 
505 /*
506  * Point is assumed to have an M value.
507  * Return NULL if point is not in the given range (inclusive)
508  * Return an LWPOINT *copy* otherwise.
509  */
510 static LWGEOM *
511 lwpoint_locate_between_m(LWPOINT *lwpoint, double m0, double m1)
512 {
513  POINT3DM p3dm;
514 
515  POSTGIS_DEBUGF(2, "lwpoint_locate_between called for lwpoint %x", lwpoint);
516 
517  lwpoint_getPoint3dm_p(lwpoint, &p3dm);
518  if ( p3dm.m >= m0 && p3dm.m <= m1)
519  {
520  POSTGIS_DEBUG(3, " lwpoint... returning a clone of input");
521 
522  return (LWGEOM *)lwpoint_clone(lwpoint);
523  }
524  else
525  {
526  POSTGIS_DEBUG(3, " lwpoint... returning a clone of input");
527 
528  return NULL;
529  }
530 }
531 
532 /*
533  * Line is assumed to have an M value.
534  *
535  * Return NULL if no parts of the line are in the given range (inclusive)
536  *
537  * Return an LWCOLLECTION with LWLINES and LWPOINT being consecutive
538  * and isolated points on the line falling in the range.
539  *
540  * X,Y and Z (if present) ordinates are interpolated.
541  *
542  */
543 static LWGEOM *
544 lwline_locate_between_m(LWLINE *lwline_in, double m0, double m1)
545 {
546  POINTARRAY *ipa=lwline_in->points;
547  int i;
548  LWGEOM **geoms;
549  int ngeoms;
550  int outtype;
551  int typeflag=0; /* see flags below */
552  const int pointflag=0x01;
553  const int lineflag=0x10;
554  POINTARRAYSET paset=ptarray_locate_between_m(ipa, m0, m1);
555 
556  POSTGIS_DEBUGF(2, "lwline_locate_between called for lwline %x", lwline_in);
557 
558  POSTGIS_DEBUGF(3, " ptarray_locate... returned %d pointarrays",
559  paset.nptarrays);
560 
561  if ( paset.nptarrays == 0 )
562  {
563  return NULL;
564  }
565 
566  ngeoms=paset.nptarrays;
567  /* TODO: rework this to reduce used memory */
568  geoms=lwalloc(sizeof(LWGEOM *)*ngeoms);
569  for (i=0; i<ngeoms; i++)
570  {
571  LWPOINT *lwpoint;
572  LWLINE *lwline;
573 
574  POINTARRAY *pa=paset.ptarrays[i];
575 
576  /* This is a point */
577  if ( pa->npoints == 1 )
578  {
579  lwpoint = lwpoint_construct(lwline_in->srid, NULL, pa);
580  geoms[i]=(LWGEOM *)lwpoint;
581  typeflag|=pointflag;
582  }
583 
584  /* This is a line */
585  else if ( pa->npoints > 1 )
586  {
587  lwline = lwline_construct(lwline_in->srid, NULL, pa);
588  geoms[i]=(LWGEOM *)lwline;
589  typeflag|=lineflag;
590  }
591 
592  /* This is a bug */
593  else
594  {
595  lwerror("ptarray_locate_between_m returned a POINARRAY set containing POINTARRAY with 0 points");
596  }
597 
598  }
599 
600  if ( ngeoms == 1 )
601  {
602  return geoms[0];
603  }
604  else
605  {
606  /* Choose best type */
607  if ( typeflag == 1 ) outtype=MULTIPOINTTYPE;
608  else if ( typeflag == 2 ) outtype=MULTILINETYPE;
609  else outtype = COLLECTIONTYPE;
610 
611  return (LWGEOM *)lwcollection_construct(outtype,
612  lwline_in->srid, NULL, ngeoms, geoms);
613  }
614 }
615 
616 /*
617  * Return a fully new allocated LWCOLLECTION
618  * always tagged as COLLECTIONTYPE.
619  */
620 static LWGEOM *
621 lwcollection_locate_between_m(LWCOLLECTION *lwcoll, double m0, double m1)
622 {
623  int i;
624  int ngeoms=0;
625  LWGEOM **geoms;
626 
627  POSTGIS_DEBUGF(2, "lwcollection_locate_between_m called for lwcoll %x", lwcoll);
628 
629  geoms=lwalloc(sizeof(LWGEOM *)*lwcoll->ngeoms);
630  for (i=0; i<lwcoll->ngeoms; i++)
631  {
632  LWGEOM *sub=lwgeom_locate_between_m(lwcoll->geoms[i],
633  m0, m1);
634  if ( sub != NULL )
635  geoms[ngeoms++] = sub;
636  }
637 
638  if ( ngeoms == 0 ) return NULL;
639 
641  lwcoll->srid, NULL, ngeoms, geoms);
642 }
643 
644 /*
645  * Return a fully allocated LWGEOM containing elements
646  * intersected/interpolated with the given M range.
647  * Return NULL if none of the elements fall in the range.
648  *
649  * m0 is assumed to be less-or-equal to m1.
650  * input LWGEOM is assumed to contain an M value.
651  *
652  */
653 static LWGEOM *
654 lwgeom_locate_between_m(LWGEOM *lwin, double m0, double m1)
655 {
656  POSTGIS_DEBUGF(2, "lwgeom_locate_between called for lwgeom %x", lwin);
657 
658  switch (lwin->type)
659  {
660  case POINTTYPE:
662  (LWPOINT *)lwin, m0, m1);
663  case LINETYPE:
665  (LWLINE *)lwin, m0, m1);
666 
667  case MULTIPOINTTYPE:
668  case MULTILINETYPE:
669  case COLLECTIONTYPE:
671  (LWCOLLECTION *)lwin, m0, m1);
672 
673  /* Polygon types are not supported */
674  case POLYGONTYPE:
675  case MULTIPOLYGONTYPE:
676  lwerror("Areal geometries are not supported by locate_between_measures");
677  return NULL;
678  }
679 
680  lwerror("Unkonwn geometry type (%s:%d)", __FILE__, __LINE__);
681  return NULL;
682 }
683 
684 /*
685  * Return a derived geometry collection value with elements that match
686  * the specified range of measures inclusively.
687  *
688  * Implements SQL/MM ST_LocateBetween(measure, measure) method.
689  * See ISO/IEC CD 13249-3:200x(E)
690  *
691  */
692 Datum LWGEOM_locate_between_m(PG_FUNCTION_ARGS);
694 Datum LWGEOM_locate_between_m(PG_FUNCTION_ARGS)
695 {
696  GSERIALIZED *gin = (GSERIALIZED *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
697  GSERIALIZED *gout;
698  double start_measure = PG_GETARG_FLOAT8(1);
699  double end_measure = PG_GETARG_FLOAT8(2);
700  LWGEOM *lwin, *lwout;
701  int hasz = gserialized_has_z(gin);
702  int hasm = gserialized_has_m(gin);
703  int type;
704 
705  elog(NOTICE,"ST_Locate_Between_Measures and ST_Locate_Along_Measure are deprecated. Use ST_LocateAlong and ST_LocateBetween.");
706 
707  if ( end_measure < start_measure )
708  {
709  lwerror("locate_between_m: 2nd arg must be bigger then 1st arg");
710  PG_RETURN_NULL();
711  }
712 
713  /*
714  * Return error if input doesn't have a measure
715  */
716  if ( ! hasm )
717  {
718  lwerror("Geometry argument does not have an 'M' ordinate");
719  PG_RETURN_NULL();
720  }
721 
722  /*
723  * Raise an error if input is a polygon, a multipolygon
724  * or a collection
725  */
726  type = gserialized_get_type(gin);
727 
728  if ( type == POLYGONTYPE || type == MULTIPOLYGONTYPE || type == COLLECTIONTYPE )
729  {
730  lwerror("Areal or Collection types are not supported");
731  PG_RETURN_NULL();
732  }
733 
734  lwin = lwgeom_from_gserialized(gin);
735 
736  lwout = lwgeom_locate_between_m(lwin,
737  start_measure, end_measure);
738 
739  lwgeom_free(lwin);
740 
741  if ( lwout == NULL )
742  {
744  gserialized_get_srid(gin), hasz, hasm);
745  }
746 
747  gout = geometry_serialize(lwout);
748  lwgeom_free(lwout);
749 
750  PG_RETURN_POINTER(gout);
751 }
752 
double x
Definition: liblwgeom.h:308
#define LINETYPE
Definition: liblwgeom.h:61
uint32_t gserialized_get_type(const GSERIALIZED *s)
Extract the geometry type from the serialized form (it hides in the anonymous data area...
Definition: g_serialized.c:56
LWLINE * lwline_measured_from_lwline(const LWLINE *lwline, double m_start, double m_end)
Add a measure dimension to a line, interpolating linearly from the start to the end value...
Definition: lwline.c:366
static POINTARRAYSET ptarray_locate_between_m(POINTARRAY *ipa, double m0, double m1)
static LWGEOM * lwcollection_locate_between_m(LWCOLLECTION *lwcoll, double m0, double m1)
double m
Definition: liblwgeom.h:308
PG_FUNCTION_INFO_V1(ST_AddMeasure)
LWCOLLECTION * lwcollection_construct(uint8_t type, int srid, GBOX *bbox, uint32_t ngeoms, LWGEOM **geoms)
Definition: lwcollection.c:30
LWGEOM * lwgeom_from_gserialized(const GSERIALIZED *g)
Allocate a new LWGEOM from a GSERIALIZED.
int npoints
Definition: liblwgeom.h:327
int gserialized_has_m(const GSERIALIZED *gser)
Check if a GSERIALIZED has an M ordinate.
Definition: g_serialized.c:30
Datum ST_LocateBetweenElevations(PG_FUNCTION_ARGS)
#define POLYGONTYPE
Definition: liblwgeom.h:62
POINTARRAY * ptarray_construct_empty(char hasz, char hasm, uint32_t maxpoints)
Create a new POINTARRAY with no points.
Definition: ptarray.c:57
void lwgeom_free(LWGEOM *geom)
Definition: lwgeom.c:1006
#define MULTIPOINTTYPE
Definition: liblwgeom.h:63
POINTARRAY ** ptarrays
Datum ST_LocateBetween(PG_FUNCTION_ARGS)
Datum LWGEOM_locate_between_m(PG_FUNCTION_ARGS)
Datum ST_InterpolatePoint(PG_FUNCTION_ARGS)
LWPOINT * lwgeom_as_lwpoint(const LWGEOM *lwgeom)
Definition: lwgeom.c:80
int32_t srid
Definition: liblwgeom.h:377
int gserialized_has_z(const GSERIALIZED *gser)
Check if a GSERIALIZED has a Z ordinate.
Definition: g_serialized.c:25
static LWGEOM * lwline_locate_between_m(LWLINE *lwline_in, double m0, double m1)
double lwgeom_interpolate_point(const LWGEOM *lwin, const LWPOINT *lwpt)
Find the measure value at the location on the line closest to the point.
static LWGEOM * lwgeom_locate_between_m(LWGEOM *lwin, double m0, double m1)
void lwerror(const char *fmt,...)
Write a notice out to the error handler.
Definition: lwutil.c:67
double ptarray_locate_point(const POINTARRAY *pa, const POINT4D *pt, double *dist, POINT4D *p_located)
Definition: ptarray.c:1267
unsigned int uint32
Definition: shpopen.c:274
int ptarray_append_point(POINTARRAY *pa, const POINT4D *pt, int allow_duplicates)
Append a point to the end of an existing POINTARRAY If allow_duplicate is LW_TRUE, then a duplicate point will not be added.
Definition: ptarray.c:141
#define LW_FALSE
Definition: liblwgeom.h:52
double m
Definition: liblwgeom.h:302
uint8_t flags
Definition: liblwgeom.h:325
Datum ST_LocateAlong(PG_FUNCTION_ARGS)
LWPOINT * lwpoint_clone(const LWPOINT *lwgeom)
Definition: lwpoint.c:206
#define LW_TRUE
Return types for functions with status returns.
Definition: liblwgeom.h:51
LWLINE * lwline_construct(int srid, GBOX *bbox, POINTARRAY *points)
Definition: lwline.c:29
LWGEOM ** geoms
Definition: liblwgeom.h:465
int lwpoint_getPoint4d_p(const LWPOINT *point, POINT4D *out)
Definition: lwpoint.c:42
static LWGEOM * lwpoint_locate_between_m(LWPOINT *lwpoint, double m0, double m1)
int32_t srid
Definition: liblwgeom.h:462
#define FLAGS_GET_Z(flags)
Macros for manipulating the 'flags' byte.
Definition: liblwgeom.h:106
LWCOLLECTION * lwgeom_clip_to_ordinate_range(const LWGEOM *lwin, char ordinate, double from, double to, double offset)
Given a geometry clip based on the from/to range of one of its ordinates (x, y, z, m).
double z
Definition: liblwgeom.h:308
LWMLINE * lwmline_measured_from_lwmline(const LWMLINE *lwmline, double m_start, double m_end)
Re-write the measure ordinate (or add one, if it isn't already there) interpolating the measure betwe...
Definition: lwmline.c:43
LWLINE * lwgeom_as_lwline(const LWGEOM *lwgeom)
Definition: lwgeom.c:89
LWGEOM * lwgeom_locate_along(const LWGEOM *lwin, double m, double offset)
Determine the location(s) along a measured line where m occurs and return as a multipoint.
Datum LWGEOM_line_locate_point(PG_FUNCTION_ARGS)
#define MULTIPOLYGONTYPE
Definition: liblwgeom.h:65
int lwpoint_getPoint3dm_p(const LWPOINT *point, POINT3DM *out)
Definition: lwpoint.c:37
Datum ST_AddMeasure(PG_FUNCTION_ARGS)
GSERIALIZED * geometry_serialize(LWGEOM *lwgeom)
#define POINTTYPE
LWTYPE numbers, used internally by PostGIS.
Definition: liblwgeom.h:60
#define FLAGS_GET_M(flags)
Definition: liblwgeom.h:107
uint8_t type
Definition: liblwgeom.h:352
LWPOINT * lwpoint_construct(int srid, GBOX *bbox, POINTARRAY *point)
Definition: lwpoint.c:96
void * lwalloc(size_t size)
Definition: lwutil.c:175
double y
Definition: liblwgeom.h:308
#define MULTILINETYPE
Definition: liblwgeom.h:64
LWCOLLECTION * lwcollection_construct_empty(uint8_t type, int srid, char hasz, char hasm)
Definition: lwcollection.c:81
int32_t gserialized_get_srid(const GSERIALIZED *s)
Extract the SRID from the serialized form (it is packed into three bytes so this is a handy function)...
Definition: g_serialized.c:70
int getPoint4d_p(const POINTARRAY *pa, int n, POINT4D *point)
Definition: lwgeom_api.c:217
#define COLLECTIONTYPE
Definition: liblwgeom.h:66
static int clip_seg_by_m_range(POINT4D *p1, POINT4D *p2, double m0, double m1)
This library is the generic geometry handling section of PostGIS.
POINTARRAY * points
Definition: liblwgeom.h:378