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