PostGIS  3.7.0dev-r@@SVN_REVISION@@

◆ lwline_interpolate_points()

POINTARRAY* lwline_interpolate_points ( const LWLINE line,
double  length_fraction,
char  repeat 
)

Interpolate one or more points along a line.

Definition at line 538 of file lwline.c.

538  {
539  POINT4D pt;
540  uint32_t i;
541  uint32_t points_to_interpolate;
542  uint32_t points_found = 0;
543  double length;
544  double length_fraction_increment = length_fraction;
545  double length_fraction_consumed = 0;
546  char has_z = (char) lwgeom_has_z(lwline_as_lwgeom(line));
547  char has_m = (char) lwgeom_has_m(lwline_as_lwgeom(line));
548  const POINTARRAY* ipa = line->points;
549  POINTARRAY* opa;
550 
551  /* Empty.InterpolatePoint == Point Empty */
552  if ( lwline_is_empty(line) )
553  {
554  return ptarray_construct_empty(has_z, has_m, 0);
555  }
556 
557  /* If distance is one of the two extremes, return the point on that
558  * end rather than doing any computations
559  */
560  if ( length_fraction == 0.0 || length_fraction == 1.0 )
561  {
562  if ( length_fraction == 0.0 )
563  getPoint4d_p(ipa, 0, &pt);
564  else
565  getPoint4d_p(ipa, ipa->npoints-1, &pt);
566 
567  opa = ptarray_construct(has_z, has_m, 1);
568  ptarray_set_point4d(opa, 0, &pt);
569 
570  return opa;
571  }
572 
573  /* Interpolate points along the line */
574  length = ptarray_length_2d(ipa);
575  points_to_interpolate = repeat ? (uint32_t) floor(1 / length_fraction) : 1;
576  opa = ptarray_construct(has_z, has_m, points_to_interpolate);
577 
578  const POINT2D* p1 = getPoint2d_cp(ipa, 0);
579  for ( i = 0; i < ipa->npoints - 1 && points_found < points_to_interpolate; i++ )
580  {
581  const POINT2D* p2 = getPoint2d_cp(ipa, i+1);
582  double segment_length_frac = distance2d_pt_pt(p1, p2) / length;
583 
584  /* If our target distance is before the total length we've seen
585  * so far. create a new point some distance down the current
586  * segment.
587  */
588  while ( length_fraction < length_fraction_consumed + segment_length_frac && points_found < points_to_interpolate )
589  {
590  POINT4D p1_4d = getPoint4d(ipa, i);
591  POINT4D p2_4d = getPoint4d(ipa, i+1);
592 
593  double segment_fraction = (length_fraction - length_fraction_consumed) / segment_length_frac;
594  interpolate_point4d(&p1_4d, &p2_4d, &pt, segment_fraction);
595  ptarray_set_point4d(opa, points_found++, &pt);
596  length_fraction += length_fraction_increment;
597  }
598 
599  length_fraction_consumed += segment_length_frac;
600 
601  p1 = p2;
602  }
603 
604  /* Return the last point on the line. This shouldn't happen, but
605  * could if there's some floating point rounding errors. */
606  if (points_found < points_to_interpolate) {
607  getPoint4d_p(ipa, ipa->npoints - 1, &pt);
608  ptarray_set_point4d(opa, points_found, &pt);
609  }
610 
611  return opa;
612 }
POINT4D getPoint4d(const POINTARRAY *pa, uint32_t n)
Definition: lwgeom_api.c:107
LWGEOM * lwline_as_lwgeom(const LWLINE *obj)
Definition: lwgeom.c:367
double distance2d_pt_pt(const POINT2D *p1, const POINT2D *p2)
Definition: measures.c:2344
void interpolate_point4d(const POINT4D *A, const POINT4D *B, POINT4D *I, double F)
Find interpolation point I between point A and point B so that the len(AI) == len(AB)*F and I falls o...
Definition: lwgeom_api.c:646
POINTARRAY * ptarray_construct(char hasz, char hasm, uint32_t npoints)
Construct an empty pointarray, allocating storage and setting the npoints, but not filling in any inf...
Definition: ptarray.c:51
double ptarray_length_2d(const POINTARRAY *pts)
Find the 2d length of the given POINTARRAY (even if it's 3d)
Definition: ptarray.c:1975
int lwgeom_has_z(const LWGEOM *geom)
Return LW_TRUE if geometry has Z ordinates.
Definition: lwgeom.c:962
POINTARRAY * ptarray_construct_empty(char hasz, char hasm, uint32_t maxpoints)
Create a new POINTARRAY with no points.
Definition: ptarray.c:59
int getPoint4d_p(const POINTARRAY *pa, uint32_t n, POINT4D *point)
Definition: lwgeom_api.c:125
int lwgeom_has_m(const LWGEOM *geom)
Return LW_TRUE if geometry has M ordinates.
Definition: lwgeom.c:969
void ptarray_set_point4d(POINTARRAY *pa, uint32_t n, const POINT4D *p4d)
Definition: lwgeom_api.c:369
int lwline_is_empty(const LWLINE *line)
static const POINT2D * getPoint2d_cp(const POINTARRAY *pa, uint32_t n)
Returns a POINT2D pointer into the POINTARRAY serialized_ptlist, suitable for reading from.
Definition: lwinline.h:97
POINTARRAY * points
Definition: liblwgeom.h:483
uint32_t npoints
Definition: liblwgeom.h:427

References distance2d_pt_pt(), getPoint2d_cp(), getPoint4d(), getPoint4d_p(), interpolate_point4d(), lwgeom_has_m(), lwgeom_has_z(), lwline_as_lwgeom(), lwline_is_empty(), POINTARRAY::npoints, LWLINE::points, ptarray_construct(), ptarray_construct_empty(), ptarray_length_2d(), and ptarray_set_point4d().

Referenced by LWGEOM_line_interpolate_point(), and test_lwline_interpolate_points().

Here is the call graph for this function:
Here is the caller graph for this function: