PostGIS 3.7.0dev-r@@SVN_REVISION@@
Loading...
Searching...
No Matches

◆ geography_substring()

LWGEOM * geography_substring ( const LWLINE lwline,
const SPHEROID s,
double  from,
double  to,
double  tolerance 
)

Return the part of a line between two fractional locations.

Definition at line 104 of file lwgeodetic_measures.c.

109{
110 POINTARRAY *dpa;
111 POINTARRAY *ipa = lwline->points;
112 LWGEOM *lwresult;
113 POINT4D pt;
114 POINT4D p1, p2;
115 GEOGRAPHIC_POINT g1, g2;
116 int nsegs, i;
117 double length, slength, tlength;
118 int state = 0; /* 0 = before, 1 = inside */
119 uint32_t srid = lwline->srid;
120
121 /*
122 * Create a dynamic pointarray with an initial capacity
123 * equal to full copy of input points
124 */
126 (char) FLAGS_GET_Z(ipa->flags),
127 (char) FLAGS_GET_M(ipa->flags),
128 ipa->npoints);
129
130 /* Compute total line length */
131 length = ptarray_length_spheroid(ipa, s);
132
133 /* Get 'from' and 'to' lengths */
134 from = length * from;
135 to = length * to;
136 tlength = 0;
137 getPoint4d_p(ipa, 0, &p1);
138 geographic_point_init(p1.x, p1.y, &g1);
139 nsegs = ipa->npoints - 1;
140 for (i = 0; i < nsegs; i++)
141 {
142 double dseg;
143 getPoint4d_p(ipa, (uint32_t) i+1, &p2);
144 geographic_point_init(p2.x, p2.y, &g2);
145
146 /* Find the length of this segment */
147 /* Special sphere case */
148 if ( s->a == s->b )
149 slength = s->radius * sphere_distance(&g1, &g2);
150 /* Spheroid case */
151 else
152 slength = spheroid_distance(&g1, &g2, s);
153
154 /* We are before requested start. */
155 if (state == 0) /* before */
156 {
157 if (fabs ( from - ( tlength + slength ) ) <= tolerance)
158 {
159 /* Second point is our start */
161 state = 1; /* we're inside now */
162 goto END;
163 }
164 else if (fabs(from - tlength) <= tolerance)
165 {
166 /* First point is our start */
168 /*
169 * We're inside now, but will check
170 * 'to' point as well
171 */
172 state = 1;
173 }
174 /*
175 * Didn't reach the 'from' point,
176 * nothing to do
177 */
178 else if (from > tlength + slength)
179 {
180 goto END;
181 }
182 else /* tlength < from < tlength+slength */
183 {
184 /* Our start is between first and second point */
185 dseg = (from - tlength) / slength;
186 interpolate_point4d_spheroid(&p1, &p2, &pt, s, dseg);
188 /* We're inside now, but will check 'to' point as well */
189 state = 1;
190 }
191 }
192
193 if (state == 1) /* inside */
194 {
195 /* 'to' point is our second point. */
196 if (fabs(to - ( tlength + slength ) ) <= tolerance )
197 {
199 break; /* substring complete */
200 }
201 /*
202 * 'to' point is our first point.
203 * (should only happen if 'to' is 0)
204 */
205 else if (fabs(to - tlength) <= tolerance)
206 {
208 break; /* substring complete */
209 }
210 /*
211 * Didn't reach the 'end' point,
212 * just copy second point
213 */
214 else if (to > tlength + slength)
215 {
217 goto END;
218 }
219 /*
220 * 'to' point falls on this segment
221 * Interpolate and break.
222 */
223 else if (to < tlength + slength )
224 {
225 dseg = (to - tlength) / slength;
226 interpolate_point4d_spheroid(&p1, &p2, &pt, s, dseg);
228 break;
229 }
230 }
231
232 END:
233 tlength += slength;
234 memcpy(&p1, &p2, sizeof(POINT4D));
235 geographic_point_init(p1.x, p1.y, &g1);
236 }
237
238 if (dpa->npoints <= 1) {
239 lwresult = lwpoint_as_lwgeom(lwpoint_construct(srid, NULL, dpa));
240 }
241 else {
242 lwresult = lwline_as_lwgeom(lwline_construct(srid, NULL, dpa));
243 }
244
245 return lwresult;
246}
char * s
Definition cu_in_wkt.c:23
LWGEOM * lwpoint_as_lwgeom(const LWPOINT *obj)
Definition lwgeom.c:372
#define LW_FALSE
Definition liblwgeom.h:94
LWPOINT * lwpoint_construct(int32_t srid, GBOX *bbox, POINTARRAY *point)
Definition lwpoint.c:129
POINTARRAY * ptarray_construct_empty(char hasz, char hasm, uint32_t maxpoints)
Create a new POINTARRAY with no points.
Definition ptarray.c:59
#define FLAGS_GET_Z(flags)
Definition liblwgeom.h:165
LWLINE * lwline_construct(int32_t srid, GBOX *bbox, POINTARRAY *points)
Definition lwline.c:42
LWGEOM * lwline_as_lwgeom(const LWLINE *obj)
Definition lwgeom.c:367
#define FLAGS_GET_M(flags)
Definition liblwgeom.h:166
int getPoint4d_p(const POINTARRAY *pa, uint32_t n, POINT4D *point)
Definition lwgeom_api.c:125
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:147
void geographic_point_init(double lon, double lat, GEOGRAPHIC_POINT *g)
Initialize a geographic point.
Definition lwgeodetic.c:180
double ptarray_length_spheroid(const POINTARRAY *pa, const SPHEROID *s)
double sphere_distance(const GEOGRAPHIC_POINT *s, const GEOGRAPHIC_POINT *e)
Given two points on a unit sphere, calculate their distance apart in radians.
Definition lwgeodetic.c:896
double spheroid_distance(const GEOGRAPHIC_POINT *a, const GEOGRAPHIC_POINT *b, const SPHEROID *spheroid)
Computes the shortest distance along the surface of the spheroid between two points,...
Definition lwspheroid.c:79
static void interpolate_point4d_spheroid(const POINT4D *p1, const POINT4D *p2, POINT4D *p, const SPHEROID *s, double f)
Find interpolation point p between geography points p1 and p2 so that the len(p1,p) == len(p1,...
Point in spherical coordinates on the world.
Definition lwgeodetic.h:54
POINTARRAY * points
Definition liblwgeom.h:483
int32_t srid
Definition liblwgeom.h:484
double x
Definition liblwgeom.h:414
double y
Definition liblwgeom.h:414
lwflags_t flags
Definition liblwgeom.h:431
uint32_t npoints
Definition liblwgeom.h:427

References POINTARRAY::flags, FLAGS_GET_M, FLAGS_GET_Z, geographic_point_init(), getPoint4d_p(), interpolate_point4d_spheroid(), LW_FALSE, lwline_as_lwgeom(), lwline_construct(), lwpoint_as_lwgeom(), lwpoint_construct(), POINTARRAY::npoints, LWLINE::points, ptarray_append_point(), ptarray_construct_empty(), ptarray_length_spheroid(), s, sphere_distance(), spheroid_distance(), LWLINE::srid, POINT4D::x, and POINT4D::y.

Referenced by geography_line_substring(), and test_geography_substring().

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