PostGIS  3.4.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 528 of file lwline.c.

528  {
529  POINT4D pt;
530  uint32_t i;
531  uint32_t points_to_interpolate;
532  uint32_t points_found = 0;
533  double length;
534  double length_fraction_increment = length_fraction;
535  double length_fraction_consumed = 0;
536  char has_z = (char) lwgeom_has_z(lwline_as_lwgeom(line));
537  char has_m = (char) lwgeom_has_m(lwline_as_lwgeom(line));
538  const POINTARRAY* ipa = line->points;
539  POINTARRAY* opa;
540 
541  /* Empty.InterpolatePoint == Point Empty */
542  if ( lwline_is_empty(line) )
543  {
544  return ptarray_construct_empty(has_z, has_m, 0);
545  }
546 
547  /* If distance is one of the two extremes, return the point on that
548  * end rather than doing any computations
549  */
550  if ( length_fraction == 0.0 || length_fraction == 1.0 )
551  {
552  if ( length_fraction == 0.0 )
553  getPoint4d_p(ipa, 0, &pt);
554  else
555  getPoint4d_p(ipa, ipa->npoints-1, &pt);
556 
557  opa = ptarray_construct(has_z, has_m, 1);
558  ptarray_set_point4d(opa, 0, &pt);
559 
560  return opa;
561  }
562 
563  /* Interpolate points along the line */
564  length = ptarray_length_2d(ipa);
565  points_to_interpolate = repeat ? (uint32_t) floor(1 / length_fraction) : 1;
566  opa = ptarray_construct(has_z, has_m, points_to_interpolate);
567 
568  const POINT2D* p1 = getPoint2d_cp(ipa, 0);
569  for ( i = 0; i < ipa->npoints - 1 && points_found < points_to_interpolate; i++ )
570  {
571  const POINT2D* p2 = getPoint2d_cp(ipa, i+1);
572  double segment_length_frac = distance2d_pt_pt(p1, p2) / length;
573 
574  /* If our target distance is before the total length we've seen
575  * so far. create a new point some distance down the current
576  * segment.
577  */
578  while ( length_fraction < length_fraction_consumed + segment_length_frac && points_found < points_to_interpolate )
579  {
580  POINT4D p1_4d = getPoint4d(ipa, i);
581  POINT4D p2_4d = getPoint4d(ipa, i+1);
582 
583  double segment_fraction = (length_fraction - length_fraction_consumed) / segment_length_frac;
584  interpolate_point4d(&p1_4d, &p2_4d, &pt, segment_fraction);
585  ptarray_set_point4d(opa, points_found++, &pt);
586  length_fraction += length_fraction_increment;
587  }
588 
589  length_fraction_consumed += segment_length_frac;
590 
591  p1 = p2;
592  }
593 
594  /* Return the last point on the line. This shouldn't happen, but
595  * could if there's some floating point rounding errors. */
596  if (points_found < points_to_interpolate) {
597  getPoint4d_p(ipa, ipa->npoints - 1, &pt);
598  ptarray_set_point4d(opa, points_found, &pt);
599  }
600 
601  return opa;
602 }
POINT4D getPoint4d(const POINTARRAY *pa, uint32_t n)
Definition: lwgeom_api.c:108
LWGEOM * lwline_as_lwgeom(const LWLINE *obj)
Definition: lwgeom.c:339
double distance2d_pt_pt(const POINT2D *p1, const POINT2D *p2)
Definition: measures.c:2398
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:649
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:1832
int lwgeom_has_z(const LWGEOM *geom)
Return LW_TRUE if geometry has Z ordinates.
Definition: lwgeom.c:934
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:941
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:101
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: