PostGIS  3.0.6dev-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 525 of file lwline.c.

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

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: