1 /**********************************************************************
2 *
3 * rttopo - topology library
4 * http://git.osgeo.org/gitea/rttopo/librttopo
5 *
6 * rttopo 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 * rttopo 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 rttopo. If not, see <http://www.gnu.org/licenses/>.
18 *
19 **********************************************************************
20 *
21 * Copyright (C) 2001-2006 Refractions Research Inc.
22 *
23 **********************************************************************/
24
25
26
27 #include "rttopo_config.h"
28 #include <stdio.h>
29 #include <stdlib.h>
30 #include <stdarg.h>
31
32 #include "librttopo_geom_internal.h"
33 #include "rtgeom_log.h"
34
35
36 /** Force Right-hand-rule on RTGEOM polygons **/
37 void
rtgeom_force_clockwise(const RTCTX * ctx,RTGEOM * rtgeom)38 rtgeom_force_clockwise(const RTCTX *ctx, RTGEOM *rtgeom)
39 {
40 RTCOLLECTION *coll;
41 int i;
42
43 switch (rtgeom->type)
44 {
45 case RTPOLYGONTYPE:
46 rtpoly_force_clockwise(ctx, (RTPOLY *)rtgeom);
47 return;
48
49 case RTTRIANGLETYPE:
50 rttriangle_force_clockwise(ctx, (RTTRIANGLE *)rtgeom);
51 return;
52
53 /* Not handle POLYHEDRALSURFACE and TIN
54 as they are supposed to be well oriented */
55 case RTMULTIPOLYGONTYPE:
56 case RTCOLLECTIONTYPE:
57 coll = (RTCOLLECTION *)rtgeom;
58 for (i=0; i<coll->ngeoms; i++)
59 rtgeom_force_clockwise(ctx, coll->geoms[i]);
60 return;
61 }
62 }
63
64 /** Reverse vertex order of RTGEOM **/
65 void
rtgeom_reverse(const RTCTX * ctx,RTGEOM * rtgeom)66 rtgeom_reverse(const RTCTX *ctx, RTGEOM *rtgeom)
67 {
68 int i;
69 RTCOLLECTION *col;
70
71 switch (rtgeom->type)
72 {
73 case RTLINETYPE:
74 rtline_reverse(ctx, (RTLINE *)rtgeom);
75 return;
76 case RTPOLYGONTYPE:
77 rtpoly_reverse(ctx, (RTPOLY *)rtgeom);
78 return;
79 case RTTRIANGLETYPE:
80 rttriangle_reverse(ctx, (RTTRIANGLE *)rtgeom);
81 return;
82 case RTMULTILINETYPE:
83 case RTMULTIPOLYGONTYPE:
84 case RTPOLYHEDRALSURFACETYPE:
85 case RTTINTYPE:
86 case RTCOLLECTIONTYPE:
87 col = (RTCOLLECTION *)rtgeom;
88 for (i=0; i<col->ngeoms; i++)
89 rtgeom_reverse(ctx, col->geoms[i]);
90 return;
91 }
92 }
93
94 RTPOINT *
rtgeom_as_rtpoint(const RTCTX * ctx,const RTGEOM * rtgeom)95 rtgeom_as_rtpoint(const RTCTX *ctx, const RTGEOM *rtgeom)
96 {
97 if ( rtgeom == NULL ) return NULL;
98 if ( rtgeom->type == RTPOINTTYPE )
99 return (RTPOINT *)rtgeom;
100 else return NULL;
101 }
102
103 RTLINE *
rtgeom_as_rtline(const RTCTX * ctx,const RTGEOM * rtgeom)104 rtgeom_as_rtline(const RTCTX *ctx, const RTGEOM *rtgeom)
105 {
106 if ( rtgeom == NULL ) return NULL;
107 if ( rtgeom->type == RTLINETYPE )
108 return (RTLINE *)rtgeom;
109 else return NULL;
110 }
111
112 RTCIRCSTRING *
rtgeom_as_rtcircstring(const RTCTX * ctx,const RTGEOM * rtgeom)113 rtgeom_as_rtcircstring(const RTCTX *ctx, const RTGEOM *rtgeom)
114 {
115 if ( rtgeom == NULL ) return NULL;
116 if ( rtgeom->type == RTCIRCSTRINGTYPE )
117 return (RTCIRCSTRING *)rtgeom;
118 else return NULL;
119 }
120
121 RTCOMPOUND *
rtgeom_as_rtcompound(const RTCTX * ctx,const RTGEOM * rtgeom)122 rtgeom_as_rtcompound(const RTCTX *ctx, const RTGEOM *rtgeom)
123 {
124 if ( rtgeom == NULL ) return NULL;
125 if ( rtgeom->type == RTCOMPOUNDTYPE )
126 return (RTCOMPOUND *)rtgeom;
127 else return NULL;
128 }
129
130 RTCURVEPOLY *
rtgeom_as_rtcurvepoly(const RTCTX * ctx,const RTGEOM * rtgeom)131 rtgeom_as_rtcurvepoly(const RTCTX *ctx, const RTGEOM *rtgeom)
132 {
133 if ( rtgeom == NULL ) return NULL;
134 if ( rtgeom->type == RTCURVEPOLYTYPE )
135 return (RTCURVEPOLY *)rtgeom;
136 else return NULL;
137 }
138
139 RTPOLY *
rtgeom_as_rtpoly(const RTCTX * ctx,const RTGEOM * rtgeom)140 rtgeom_as_rtpoly(const RTCTX *ctx, const RTGEOM *rtgeom)
141 {
142 if ( rtgeom == NULL ) return NULL;
143 if ( rtgeom->type == RTPOLYGONTYPE )
144 return (RTPOLY *)rtgeom;
145 else return NULL;
146 }
147
148 RTTRIANGLE *
rtgeom_as_rttriangle(const RTCTX * ctx,const RTGEOM * rtgeom)149 rtgeom_as_rttriangle(const RTCTX *ctx, const RTGEOM *rtgeom)
150 {
151 if ( rtgeom == NULL ) return NULL;
152 if ( rtgeom->type == RTTRIANGLETYPE )
153 return (RTTRIANGLE *)rtgeom;
154 else return NULL;
155 }
156
157 RTCOLLECTION *
rtgeom_as_rtcollection(const RTCTX * ctx,const RTGEOM * rtgeom)158 rtgeom_as_rtcollection(const RTCTX *ctx, const RTGEOM *rtgeom)
159 {
160 if ( rtgeom == NULL ) return NULL;
161 if ( rtgeom_is_collection(ctx, rtgeom) )
162 return (RTCOLLECTION*)rtgeom;
163 else return NULL;
164 }
165
166 RTMPOINT *
rtgeom_as_rtmpoint(const RTCTX * ctx,const RTGEOM * rtgeom)167 rtgeom_as_rtmpoint(const RTCTX *ctx, const RTGEOM *rtgeom)
168 {
169 if ( rtgeom == NULL ) return NULL;
170 if ( rtgeom->type == RTMULTIPOINTTYPE )
171 return (RTMPOINT *)rtgeom;
172 else return NULL;
173 }
174
175 RTMLINE *
rtgeom_as_rtmline(const RTCTX * ctx,const RTGEOM * rtgeom)176 rtgeom_as_rtmline(const RTCTX *ctx, const RTGEOM *rtgeom)
177 {
178 if ( rtgeom == NULL ) return NULL;
179 if ( rtgeom->type == RTMULTILINETYPE )
180 return (RTMLINE *)rtgeom;
181 else return NULL;
182 }
183
184 RTMPOLY *
rtgeom_as_rtmpoly(const RTCTX * ctx,const RTGEOM * rtgeom)185 rtgeom_as_rtmpoly(const RTCTX *ctx, const RTGEOM *rtgeom)
186 {
187 if ( rtgeom == NULL ) return NULL;
188 if ( rtgeom->type == RTMULTIPOLYGONTYPE )
189 return (RTMPOLY *)rtgeom;
190 else return NULL;
191 }
192
193 RTPSURFACE *
rtgeom_as_rtpsurface(const RTCTX * ctx,const RTGEOM * rtgeom)194 rtgeom_as_rtpsurface(const RTCTX *ctx, const RTGEOM *rtgeom)
195 {
196 if ( rtgeom->type == RTPOLYHEDRALSURFACETYPE )
197 return (RTPSURFACE *)rtgeom;
198 else return NULL;
199 }
200
201 RTTIN *
rtgeom_as_rttin(const RTCTX * ctx,const RTGEOM * rtgeom)202 rtgeom_as_rttin(const RTCTX *ctx, const RTGEOM *rtgeom)
203 {
204 if ( rtgeom->type == RTTINTYPE )
205 return (RTTIN *)rtgeom;
206 else return NULL;
207 }
208
rttin_as_rtgeom(const RTCTX * ctx,const RTTIN * obj)209 RTGEOM * rttin_as_rtgeom(const RTCTX *ctx, const RTTIN *obj)
210 {
211 return (RTGEOM *)obj;
212 }
213
rtpsurface_as_rtgeom(const RTCTX * ctx,const RTPSURFACE * obj)214 RTGEOM * rtpsurface_as_rtgeom(const RTCTX *ctx, const RTPSURFACE *obj)
215 {
216 return (RTGEOM *)obj;
217 }
218
rtmpoly_as_rtgeom(const RTCTX * ctx,const RTMPOLY * obj)219 RTGEOM * rtmpoly_as_rtgeom(const RTCTX *ctx, const RTMPOLY *obj)
220 {
221 if ( obj == NULL ) return NULL;
222 return (RTGEOM *)obj;
223 }
rtmline_as_rtgeom(const RTCTX * ctx,const RTMLINE * obj)224 RTGEOM * rtmline_as_rtgeom(const RTCTX *ctx, const RTMLINE *obj)
225 {
226 if ( obj == NULL ) return NULL;
227 return (RTGEOM *)obj;
228 }
rtmpoint_as_rtgeom(const RTCTX * ctx,const RTMPOINT * obj)229 RTGEOM * rtmpoint_as_rtgeom(const RTCTX *ctx, const RTMPOINT *obj)
230 {
231 if ( obj == NULL ) return NULL;
232 return (RTGEOM *)obj;
233 }
rtcollection_as_rtgeom(const RTCTX * ctx,const RTCOLLECTION * obj)234 RTGEOM * rtcollection_as_rtgeom(const RTCTX *ctx, const RTCOLLECTION *obj)
235 {
236 if ( obj == NULL ) return NULL;
237 return (RTGEOM *)obj;
238 }
rtcircstring_as_rtgeom(const RTCTX * ctx,const RTCIRCSTRING * obj)239 RTGEOM * rtcircstring_as_rtgeom(const RTCTX *ctx, const RTCIRCSTRING *obj)
240 {
241 if ( obj == NULL ) return NULL;
242 return (RTGEOM *)obj;
243 }
rtcurvepoly_as_rtgeom(const RTCTX * ctx,const RTCURVEPOLY * obj)244 RTGEOM * rtcurvepoly_as_rtgeom(const RTCTX *ctx, const RTCURVEPOLY *obj)
245 {
246 if ( obj == NULL ) return NULL;
247 return (RTGEOM *)obj;
248 }
rtcompound_as_rtgeom(const RTCTX * ctx,const RTCOMPOUND * obj)249 RTGEOM * rtcompound_as_rtgeom(const RTCTX *ctx, const RTCOMPOUND *obj)
250 {
251 if ( obj == NULL ) return NULL;
252 return (RTGEOM *)obj;
253 }
rtpoly_as_rtgeom(const RTCTX * ctx,const RTPOLY * obj)254 RTGEOM * rtpoly_as_rtgeom(const RTCTX *ctx, const RTPOLY *obj)
255 {
256 if ( obj == NULL ) return NULL;
257 return (RTGEOM *)obj;
258 }
rttriangle_as_rtgeom(const RTCTX * ctx,const RTTRIANGLE * obj)259 RTGEOM * rttriangle_as_rtgeom(const RTCTX *ctx, const RTTRIANGLE *obj)
260 {
261 if ( obj == NULL ) return NULL;
262 return (RTGEOM *)obj;
263 }
rtline_as_rtgeom(const RTCTX * ctx,const RTLINE * obj)264 RTGEOM * rtline_as_rtgeom(const RTCTX *ctx, const RTLINE *obj)
265 {
266 if ( obj == NULL ) return NULL;
267 return (RTGEOM *)obj;
268 }
rtpoint_as_rtgeom(const RTCTX * ctx,const RTPOINT * obj)269 RTGEOM * rtpoint_as_rtgeom(const RTCTX *ctx, const RTPOINT *obj)
270 {
271 if ( obj == NULL ) return NULL;
272 return (RTGEOM *)obj;
273 }
274
275
276 /**
277 ** Look-up for the correct MULTI* type promotion for singleton types.
278 */
279 uint8_t RTMULTITYPE[RTNUMTYPES] =
280 {
281 0,
282 RTMULTIPOINTTYPE, /* 1 */
283 RTMULTILINETYPE, /* 2 */
284 RTMULTIPOLYGONTYPE, /* 3 */
285 0,0,0,0,
286 RTMULTICURVETYPE, /* 8 */
287 RTMULTICURVETYPE, /* 9 */
288 RTMULTISURFACETYPE, /* 10 */
289 RTPOLYHEDRALSURFACETYPE, /* 11 */
290 0, 0,
291 RTTINTYPE, /* 14 */
292 0
293 };
294
295 /**
296 * Create a new RTGEOM of the appropriate MULTI* type.
297 */
298 RTGEOM *
rtgeom_as_multi(const RTCTX * ctx,const RTGEOM * rtgeom)299 rtgeom_as_multi(const RTCTX *ctx, const RTGEOM *rtgeom)
300 {
301 RTGEOM **ogeoms;
302 RTGEOM *ogeom = NULL;
303 RTGBOX *box = NULL;
304 int type;
305
306 type = rtgeom->type;
307
308 if ( ! RTMULTITYPE[type] ) return rtgeom_clone(ctx, rtgeom);
309
310 if( rtgeom_is_empty(ctx, rtgeom) )
311 {
312 ogeom = (RTGEOM *)rtcollection_construct_empty(ctx,
313 RTMULTITYPE[type],
314 rtgeom->srid,
315 RTFLAGS_GET_Z(rtgeom->flags),
316 RTFLAGS_GET_M(rtgeom->flags)
317 );
318 }
319 else
320 {
321 ogeoms = rtalloc(ctx, sizeof(RTGEOM*));
322 ogeoms[0] = rtgeom_clone(ctx, rtgeom);
323
324 /* Sub-geometries are not allowed to have bboxes or SRIDs, move the bbox to the collection */
325 box = ogeoms[0]->bbox;
326 ogeoms[0]->bbox = NULL;
327 ogeoms[0]->srid = SRID_UNKNOWN;
328
329 ogeom = (RTGEOM *)rtcollection_construct(ctx, RTMULTITYPE[type], rtgeom->srid, box, 1, ogeoms);
330 }
331
332 return ogeom;
333 }
334
335 /**
336 * Create a new RTGEOM of the appropriate CURVE* type.
337 */
338 RTGEOM *
rtgeom_as_curve(const RTCTX * ctx,const RTGEOM * rtgeom)339 rtgeom_as_curve(const RTCTX *ctx, const RTGEOM *rtgeom)
340 {
341 RTGEOM *ogeom;
342 int type = rtgeom->type;
343 /*
344 int hasz = RTFLAGS_GET_Z(rtgeom->flags);
345 int hasm = RTFLAGS_GET_M(rtgeom->flags);
346 int srid = rtgeom->srid;
347 */
348
349 switch(type)
350 {
351 case RTLINETYPE:
352 /* turn to COMPOUNDCURVE */
353 ogeom = (RTGEOM*)rtcompound_construct_from_rtline(ctx, (RTLINE*)rtgeom);
354 break;
355 case RTPOLYGONTYPE:
356 ogeom = (RTGEOM*)rtcurvepoly_construct_from_rtpoly(ctx, rtgeom_as_rtpoly(ctx, rtgeom));
357 break;
358 case RTMULTILINETYPE:
359 /* turn to MULTICURVE */
360 ogeom = rtgeom_clone(ctx, rtgeom);
361 ogeom->type = RTMULTICURVETYPE;
362 break;
363 case RTMULTIPOLYGONTYPE:
364 /* turn to MULTISURFACE */
365 ogeom = rtgeom_clone(ctx, rtgeom);
366 ogeom->type = RTMULTISURFACETYPE;
367 break;
368 case RTCOLLECTIONTYPE:
369 default:
370 ogeom = rtgeom_clone(ctx, rtgeom);
371 break;
372 }
373
374 /* TODO: copy bbox from input geom ? */
375
376 return ogeom;
377 }
378
379
380 /**
381 * Free the containing RTGEOM and the associated BOX. Leave the underlying
382 * geoms/points/point objects intact. Useful for functions that are stripping
383 * out subcomponents of complex objects, or building up new temporary objects
384 * on top of subcomponents.
385 */
386 void
rtgeom_release(const RTCTX * ctx,RTGEOM * rtgeom)387 rtgeom_release(const RTCTX *ctx, RTGEOM *rtgeom)
388 {
389 if ( ! rtgeom )
390 rterror(ctx, "rtgeom_release: someone called on 0x0");
391
392 RTDEBUGF(ctx, 3, "releasing type %s", rttype_name(ctx, rtgeom->type));
393
394 /* Drop bounding box (artays a copy) */
395 if ( rtgeom->bbox )
396 {
397 RTDEBUGF(ctx, 3, "rtgeom_release: releasing bbox. %p", rtgeom->bbox);
398 rtfree(ctx, rtgeom->bbox);
399 }
400 rtfree(ctx, rtgeom);
401
402 }
403
404
405 /* @brief Clone RTGEOM object. Serialized point lists are not copied.
406 *
407 * @see ptarray_clone
408 */
409 RTGEOM *
rtgeom_clone(const RTCTX * ctx,const RTGEOM * rtgeom)410 rtgeom_clone(const RTCTX *ctx, const RTGEOM *rtgeom)
411 {
412 RTDEBUGF(ctx, 2, "rtgeom_clone called with %p, %s",
413 rtgeom, rttype_name(ctx, rtgeom->type));
414
415 switch (rtgeom->type)
416 {
417 case RTPOINTTYPE:
418 return (RTGEOM *)rtpoint_clone(ctx, (RTPOINT *)rtgeom);
419 case RTLINETYPE:
420 return (RTGEOM *)rtline_clone(ctx, (RTLINE *)rtgeom);
421 case RTCIRCSTRINGTYPE:
422 return (RTGEOM *)rtcircstring_clone(ctx, (RTCIRCSTRING *)rtgeom);
423 case RTPOLYGONTYPE:
424 return (RTGEOM *)rtpoly_clone(ctx, (RTPOLY *)rtgeom);
425 case RTTRIANGLETYPE:
426 return (RTGEOM *)rttriangle_clone(ctx, (RTTRIANGLE *)rtgeom);
427 case RTCOMPOUNDTYPE:
428 case RTCURVEPOLYTYPE:
429 case RTMULTICURVETYPE:
430 case RTMULTISURFACETYPE:
431 case RTMULTIPOINTTYPE:
432 case RTMULTILINETYPE:
433 case RTMULTIPOLYGONTYPE:
434 case RTPOLYHEDRALSURFACETYPE:
435 case RTTINTYPE:
436 case RTCOLLECTIONTYPE:
437 return (RTGEOM *)rtcollection_clone(ctx, (RTCOLLECTION *)rtgeom);
438 default:
439 rterror(ctx, "rtgeom_clone: Unknown geometry type: %s", rttype_name(ctx, rtgeom->type));
440 return NULL;
441 }
442 }
443
444 /**
445 * Deep-clone an #RTGEOM object. #RTPOINTARRAY <em>are</em> copied.
446 */
447 RTGEOM *
rtgeom_clone_deep(const RTCTX * ctx,const RTGEOM * rtgeom)448 rtgeom_clone_deep(const RTCTX *ctx, const RTGEOM *rtgeom)
449 {
450 RTDEBUGF(ctx, 2, "rtgeom_clone called with %p, %s",
451 rtgeom, rttype_name(ctx, rtgeom->type));
452
453 switch (rtgeom->type)
454 {
455 case RTPOINTTYPE:
456 case RTLINETYPE:
457 case RTCIRCSTRINGTYPE:
458 case RTTRIANGLETYPE:
459 return (RTGEOM *)rtline_clone_deep(ctx, (RTLINE *)rtgeom);
460 case RTPOLYGONTYPE:
461 return (RTGEOM *)rtpoly_clone_deep(ctx, (RTPOLY *)rtgeom);
462 case RTCOMPOUNDTYPE:
463 case RTCURVEPOLYTYPE:
464 case RTMULTICURVETYPE:
465 case RTMULTISURFACETYPE:
466 case RTMULTIPOINTTYPE:
467 case RTMULTILINETYPE:
468 case RTMULTIPOLYGONTYPE:
469 case RTPOLYHEDRALSURFACETYPE:
470 case RTTINTYPE:
471 case RTCOLLECTIONTYPE:
472 return (RTGEOM *)rtcollection_clone_deep(ctx, (RTCOLLECTION *)rtgeom);
473 default:
474 rterror(ctx, "rtgeom_clone_deep: Unknown geometry type: %s", rttype_name(ctx, rtgeom->type));
475 return NULL;
476 }
477 }
478
479
480 /**
481 * Return an alloced string
482 */
483 char*
rtgeom_to_ewkt(const RTCTX * ctx,const RTGEOM * rtgeom)484 rtgeom_to_ewkt(const RTCTX *ctx, const RTGEOM *rtgeom)
485 {
486 char* wkt = NULL;
487 size_t wkt_size = 0;
488
489 wkt = rtgeom_to_wkt(ctx, rtgeom, RTWKT_EXTENDED, 12, &wkt_size);
490
491 if ( ! wkt )
492 {
493 rterror(ctx, "Error writing geom %p to RTWKT", rtgeom);
494 }
495
496 return wkt;
497 }
498
499 /**
500 * @brief geom1 same as geom2
501 * iff
502 * + have same type
503 * + have same # objects
504 * + have same bvol
505 * + each object in geom1 has a corresponding object in geom2 (see above)
506 * @param rtgeom1
507 * @param rtgeom2
508 */
509 char
rtgeom_same(const RTCTX * ctx,const RTGEOM * rtgeom1,const RTGEOM * rtgeom2)510 rtgeom_same(const RTCTX *ctx, const RTGEOM *rtgeom1, const RTGEOM *rtgeom2)
511 {
512 RTDEBUGF(ctx, 2, "rtgeom_same(ctx, %s, %s) called",
513 rttype_name(ctx, rtgeom1->type),
514 rttype_name(ctx, rtgeom2->type));
515
516 if ( rtgeom1->type != rtgeom2->type )
517 {
518 RTDEBUG(ctx, 3, " type differ");
519
520 return RT_FALSE;
521 }
522
523 if ( RTFLAGS_GET_ZM(rtgeom1->flags) != RTFLAGS_GET_ZM(rtgeom2->flags) )
524 {
525 RTDEBUG(ctx, 3, " ZM flags differ");
526
527 return RT_FALSE;
528 }
529
530 /* Check boxes if both already computed */
531 if ( rtgeom1->bbox && rtgeom2->bbox )
532 {
533 /*rtnotice(ctx, "bbox1:%p, bbox2:%p", rtgeom1->bbox, rtgeom2->bbox);*/
534 if ( ! gbox_same(ctx, rtgeom1->bbox, rtgeom2->bbox) )
535 {
536 RTDEBUG(ctx, 3, " bounding boxes differ");
537
538 return RT_FALSE;
539 }
540 }
541
542 /* geoms have same type, invoke type-specific function */
543 switch (rtgeom1->type)
544 {
545 case RTPOINTTYPE:
546 return rtpoint_same(ctx, (RTPOINT *)rtgeom1,
547 (RTPOINT *)rtgeom2);
548 case RTLINETYPE:
549 return rtline_same(ctx, (RTLINE *)rtgeom1,
550 (RTLINE *)rtgeom2);
551 case RTPOLYGONTYPE:
552 return rtpoly_same(ctx, (RTPOLY *)rtgeom1,
553 (RTPOLY *)rtgeom2);
554 case RTTRIANGLETYPE:
555 return rttriangle_same(ctx, (RTTRIANGLE *)rtgeom1,
556 (RTTRIANGLE *)rtgeom2);
557 case RTCIRCSTRINGTYPE:
558 return rtcircstring_same(ctx, (RTCIRCSTRING *)rtgeom1,
559 (RTCIRCSTRING *)rtgeom2);
560 case RTMULTIPOINTTYPE:
561 case RTMULTILINETYPE:
562 case RTMULTIPOLYGONTYPE:
563 case RTMULTICURVETYPE:
564 case RTMULTISURFACETYPE:
565 case RTCOMPOUNDTYPE:
566 case RTCURVEPOLYTYPE:
567 case RTPOLYHEDRALSURFACETYPE:
568 case RTTINTYPE:
569 case RTCOLLECTIONTYPE:
570 return rtcollection_same(ctx, (RTCOLLECTION *)rtgeom1,
571 (RTCOLLECTION *)rtgeom2);
572 default:
573 rterror(ctx, "rtgeom_same: unsupported geometry type: %s",
574 rttype_name(ctx, rtgeom1->type));
575 return RT_FALSE;
576 }
577
578 }
579
580 int
rtpoint_inside_circle(const RTCTX * ctx,const RTPOINT * p,double cx,double cy,double rad)581 rtpoint_inside_circle(const RTCTX *ctx, const RTPOINT *p, double cx, double cy, double rad)
582 {
583 const RTPOINT2D *pt;
584 RTPOINT2D center;
585
586 if ( ! p || ! p->point )
587 return RT_FALSE;
588
589 pt = rt_getPoint2d_cp(ctx, p->point, 0);
590
591 center.x = cx;
592 center.y = cy;
593
594 if ( distance2d_pt_pt(ctx, pt, ¢er) < rad )
595 return RT_TRUE;
596
597 return RT_FALSE;
598 }
599
600 void
rtgeom_drop_bbox(const RTCTX * ctx,RTGEOM * rtgeom)601 rtgeom_drop_bbox(const RTCTX *ctx, RTGEOM *rtgeom)
602 {
603 if ( rtgeom->bbox ) rtfree(ctx, rtgeom->bbox);
604 rtgeom->bbox = NULL;
605 RTFLAGS_SET_BBOX(rtgeom->flags, 0);
606 }
607
608 /**
609 * Ensure there's a box in the RTGEOM.
610 * If the box is already there just return,
611 * else compute it.
612 */
613 void
rtgeom_add_bbox(const RTCTX * ctx,RTGEOM * rtgeom)614 rtgeom_add_bbox(const RTCTX *ctx, RTGEOM *rtgeom)
615 {
616 /* an empty RTGEOM has no bbox */
617 if ( rtgeom_is_empty(ctx, rtgeom) ) return;
618
619 if ( rtgeom->bbox ) return;
620 RTFLAGS_SET_BBOX(rtgeom->flags, 1);
621 rtgeom->bbox = gbox_new(ctx, rtgeom->flags);
622 rtgeom_calculate_gbox(ctx, rtgeom, rtgeom->bbox);
623 }
624
625 void
rtgeom_add_bbox_deep(const RTCTX * ctx,RTGEOM * rtgeom,RTGBOX * gbox)626 rtgeom_add_bbox_deep(const RTCTX *ctx, RTGEOM *rtgeom, RTGBOX *gbox)
627 {
628 if ( rtgeom_is_empty(ctx, rtgeom) ) return;
629
630 RTFLAGS_SET_BBOX(rtgeom->flags, 1);
631
632 if ( ! ( gbox || rtgeom->bbox ) )
633 {
634 rtgeom->bbox = gbox_new(ctx, rtgeom->flags);
635 rtgeom_calculate_gbox(ctx, rtgeom, rtgeom->bbox);
636 }
637 else if ( gbox && ! rtgeom->bbox )
638 {
639 rtgeom->bbox = gbox_clone(ctx, gbox);
640 }
641
642 if ( rtgeom_is_collection(ctx, rtgeom) )
643 {
644 int i;
645 RTCOLLECTION *rtcol = (RTCOLLECTION*)rtgeom;
646
647 for ( i = 0; i < rtcol->ngeoms; i++ )
648 {
649 rtgeom_add_bbox_deep(ctx, rtcol->geoms[i], rtgeom->bbox);
650 }
651 }
652 }
653
654 const RTGBOX *
rtgeom_get_bbox(const RTCTX * ctx,const RTGEOM * rtg)655 rtgeom_get_bbox(const RTCTX *ctx, const RTGEOM *rtg)
656 {
657 /* add it if not already there */
658 rtgeom_add_bbox(ctx, (RTGEOM *)rtg);
659 return rtg->bbox;
660 }
661
662
663 /**
664 * Calculate the gbox for this goemetry, a cartesian box or
665 * geodetic box, depending on how it is flagged.
666 */
rtgeom_calculate_gbox(const RTCTX * ctx,const RTGEOM * rtgeom,RTGBOX * gbox)667 int rtgeom_calculate_gbox(const RTCTX *ctx, const RTGEOM *rtgeom, RTGBOX *gbox)
668 {
669 gbox->flags = rtgeom->flags;
670 if( RTFLAGS_GET_GEODETIC(rtgeom->flags) )
671 return rtgeom_calculate_gbox_geodetic(ctx, rtgeom, gbox);
672 else
673 return rtgeom_calculate_gbox_cartesian(ctx, rtgeom, gbox);
674 }
675
676 void
rtgeom_drop_srid(const RTCTX * ctx,RTGEOM * rtgeom)677 rtgeom_drop_srid(const RTCTX *ctx, RTGEOM *rtgeom)
678 {
679 rtgeom->srid = SRID_UNKNOWN; /* TODO: To be changed to SRID_UNKNOWN */
680 }
681
682 RTGEOM *
rtgeom_segmentize2d(const RTCTX * ctx,RTGEOM * rtgeom,double dist)683 rtgeom_segmentize2d(const RTCTX *ctx, RTGEOM *rtgeom, double dist)
684 {
685 switch (rtgeom->type)
686 {
687 case RTLINETYPE:
688 return (RTGEOM *)rtline_segmentize2d(ctx, (RTLINE *)rtgeom,
689 dist);
690 case RTPOLYGONTYPE:
691 return (RTGEOM *)rtpoly_segmentize2d(ctx, (RTPOLY *)rtgeom,
692 dist);
693 case RTMULTILINETYPE:
694 case RTMULTIPOLYGONTYPE:
695 case RTCOLLECTIONTYPE:
696 return (RTGEOM *)rtcollection_segmentize2d(ctx,
697 (RTCOLLECTION *)rtgeom, dist);
698
699 default:
700 return rtgeom_clone(ctx, rtgeom);
701 }
702 }
703
704 RTGEOM*
rtgeom_force_2d(const RTCTX * ctx,const RTGEOM * geom)705 rtgeom_force_2d(const RTCTX *ctx, const RTGEOM *geom)
706 {
707 return rtgeom_force_dims(ctx, geom, 0, 0);
708 }
709
710 RTGEOM*
rtgeom_force_3dz(const RTCTX * ctx,const RTGEOM * geom)711 rtgeom_force_3dz(const RTCTX *ctx, const RTGEOM *geom)
712 {
713 return rtgeom_force_dims(ctx, geom, 1, 0);
714 }
715
716 RTGEOM*
rtgeom_force_3dm(const RTCTX * ctx,const RTGEOM * geom)717 rtgeom_force_3dm(const RTCTX *ctx, const RTGEOM *geom)
718 {
719 return rtgeom_force_dims(ctx, geom, 0, 1);
720 }
721
722 RTGEOM*
rtgeom_force_4d(const RTCTX * ctx,const RTGEOM * geom)723 rtgeom_force_4d(const RTCTX *ctx, const RTGEOM *geom)
724 {
725 return rtgeom_force_dims(ctx, geom, 1, 1);
726 }
727
728 RTGEOM*
rtgeom_force_dims(const RTCTX * ctx,const RTGEOM * geom,int hasz,int hasm)729 rtgeom_force_dims(const RTCTX *ctx, const RTGEOM *geom, int hasz, int hasm)
730 {
731 switch(geom->type)
732 {
733 case RTPOINTTYPE:
734 return rtpoint_as_rtgeom(ctx, rtpoint_force_dims(ctx, (RTPOINT*)geom, hasz, hasm));
735 case RTCIRCSTRINGTYPE:
736 case RTLINETYPE:
737 case RTTRIANGLETYPE:
738 return rtline_as_rtgeom(ctx, rtline_force_dims(ctx, (RTLINE*)geom, hasz, hasm));
739 case RTPOLYGONTYPE:
740 return rtpoly_as_rtgeom(ctx, rtpoly_force_dims(ctx, (RTPOLY*)geom, hasz, hasm));
741 case RTCOMPOUNDTYPE:
742 case RTCURVEPOLYTYPE:
743 case RTMULTICURVETYPE:
744 case RTMULTISURFACETYPE:
745 case RTMULTIPOINTTYPE:
746 case RTMULTILINETYPE:
747 case RTMULTIPOLYGONTYPE:
748 case RTPOLYHEDRALSURFACETYPE:
749 case RTTINTYPE:
750 case RTCOLLECTIONTYPE:
751 return rtcollection_as_rtgeom(ctx, rtcollection_force_dims(ctx, (RTCOLLECTION*)geom, hasz, hasm));
752 default:
753 rterror(ctx, "rtgeom_force_2d: unsupported geom type: %s", rttype_name(ctx, geom->type));
754 return NULL;
755 }
756 }
757
758 RTGEOM*
rtgeom_force_sfs(const RTCTX * ctx,RTGEOM * geom,int version)759 rtgeom_force_sfs(const RTCTX *ctx, RTGEOM *geom, int version)
760 {
761 RTCOLLECTION *col;
762 int i;
763 RTGEOM *g;
764
765 /* SFS 1.2 version */
766 if (version == 120)
767 {
768 switch(geom->type)
769 {
770 /* SQL/MM types */
771 case RTCIRCSTRINGTYPE:
772 case RTCOMPOUNDTYPE:
773 case RTCURVEPOLYTYPE:
774 case RTMULTICURVETYPE:
775 case RTMULTISURFACETYPE:
776 return rtgeom_stroke(ctx, geom, 32);
777
778 case RTCOLLECTIONTYPE:
779 col = (RTCOLLECTION*)geom;
780 for ( i = 0; i < col->ngeoms; i++ )
781 col->geoms[i] = rtgeom_force_sfs(ctx, (RTGEOM*)col->geoms[i], version);
782
783 return rtcollection_as_rtgeom(ctx, (RTCOLLECTION*)geom);
784
785 default:
786 return (RTGEOM *)geom;
787 }
788 }
789
790
791 /* SFS 1.1 version */
792 switch(geom->type)
793 {
794 /* SQL/MM types */
795 case RTCIRCSTRINGTYPE:
796 case RTCOMPOUNDTYPE:
797 case RTCURVEPOLYTYPE:
798 case RTMULTICURVETYPE:
799 case RTMULTISURFACETYPE:
800 return rtgeom_stroke(ctx, geom, 32);
801
802 /* SFS 1.2 types */
803 case RTTRIANGLETYPE:
804 g = rtpoly_as_rtgeom(ctx, rtpoly_from_rtlines(ctx, (RTLINE*)geom, 0, NULL));
805 rtgeom_free(ctx, geom);
806 return g;
807
808 case RTTINTYPE:
809 col = (RTCOLLECTION*) geom;
810 for ( i = 0; i < col->ngeoms; i++ )
811 {
812 g = rtpoly_as_rtgeom(ctx, rtpoly_from_rtlines(ctx, (RTLINE*)col->geoms[i], 0, NULL));
813 rtgeom_free(ctx, col->geoms[i]);
814 col->geoms[i] = g;
815 }
816 col->type = RTCOLLECTIONTYPE;
817 return rtmpoly_as_rtgeom(ctx, (RTMPOLY*)geom);
818
819 case RTPOLYHEDRALSURFACETYPE:
820 geom->type = RTCOLLECTIONTYPE;
821 return (RTGEOM *)geom;
822
823 /* Collection */
824 case RTCOLLECTIONTYPE:
825 col = (RTCOLLECTION*)geom;
826 for ( i = 0; i < col->ngeoms; i++ )
827 col->geoms[i] = rtgeom_force_sfs(ctx, (RTGEOM*)col->geoms[i], version);
828
829 return rtcollection_as_rtgeom(ctx, (RTCOLLECTION*)geom);
830
831 default:
832 return (RTGEOM *)geom;
833 }
834 }
835
836 int32_t
rtgeom_get_srid(const RTCTX * ctx,const RTGEOM * geom)837 rtgeom_get_srid(const RTCTX *ctx, const RTGEOM *geom)
838 {
839 if ( ! geom ) return SRID_UNKNOWN;
840 return geom->srid;
841 }
842
843 uint32_t
rtgeom_get_type(const RTCTX * ctx,const RTGEOM * geom)844 rtgeom_get_type(const RTCTX *ctx, const RTGEOM *geom)
845 {
846 if ( ! geom ) return 0;
847 return geom->type;
848 }
849
850 int
rtgeom_has_z(const RTCTX * ctx,const RTGEOM * geom)851 rtgeom_has_z(const RTCTX *ctx, const RTGEOM *geom)
852 {
853 if ( ! geom ) return RT_FALSE;
854 return RTFLAGS_GET_Z(geom->flags);
855 }
856
857 int
rtgeom_has_m(const RTCTX * ctx,const RTGEOM * geom)858 rtgeom_has_m(const RTCTX *ctx, const RTGEOM *geom)
859 {
860 if ( ! geom ) return RT_FALSE;
861 return RTFLAGS_GET_M(geom->flags);
862 }
863
864 int
rtgeom_ndims(const RTCTX * ctx,const RTGEOM * geom)865 rtgeom_ndims(const RTCTX *ctx, const RTGEOM *geom)
866 {
867 if ( ! geom ) return 0;
868 return RTFLAGS_NDIMS(geom->flags);
869 }
870
871
872 void
rtgeom_set_geodetic(const RTCTX * ctx,RTGEOM * geom,int value)873 rtgeom_set_geodetic(const RTCTX *ctx, RTGEOM *geom, int value)
874 {
875 RTPOINT *pt;
876 RTLINE *ln;
877 RTPOLY *ply;
878 RTCOLLECTION *col;
879 int i;
880
881 RTFLAGS_SET_GEODETIC(geom->flags, value);
882 if ( geom->bbox )
883 RTFLAGS_SET_GEODETIC(geom->bbox->flags, value);
884
885 switch(geom->type)
886 {
887 case RTPOINTTYPE:
888 pt = (RTPOINT*)geom;
889 if ( pt->point )
890 RTFLAGS_SET_GEODETIC(pt->point->flags, value);
891 break;
892 case RTLINETYPE:
893 ln = (RTLINE*)geom;
894 if ( ln->points )
895 RTFLAGS_SET_GEODETIC(ln->points->flags, value);
896 break;
897 case RTPOLYGONTYPE:
898 ply = (RTPOLY*)geom;
899 for ( i = 0; i < ply->nrings; i++ )
900 RTFLAGS_SET_GEODETIC(ply->rings[i]->flags, value);
901 break;
902 case RTMULTIPOINTTYPE:
903 case RTMULTILINETYPE:
904 case RTMULTIPOLYGONTYPE:
905 case RTCOLLECTIONTYPE:
906 col = (RTCOLLECTION*)geom;
907 for ( i = 0; i < col->ngeoms; i++ )
908 rtgeom_set_geodetic(ctx, col->geoms[i], value);
909 break;
910 default:
911 rterror(ctx, "rtgeom_set_geodetic: unsupported geom type: %s", rttype_name(ctx, geom->type));
912 return;
913 }
914 }
915
916 void
rtgeom_longitude_shift(const RTCTX * ctx,RTGEOM * rtgeom)917 rtgeom_longitude_shift(const RTCTX *ctx, RTGEOM *rtgeom)
918 {
919 int i;
920 switch (rtgeom->type)
921 {
922 RTPOINT *point;
923 RTLINE *line;
924 RTPOLY *poly;
925 RTTRIANGLE *triangle;
926 RTCOLLECTION *coll;
927
928 case RTPOINTTYPE:
929 point = (RTPOINT *)rtgeom;
930 ptarray_longitude_shift(ctx, point->point);
931 return;
932 case RTLINETYPE:
933 line = (RTLINE *)rtgeom;
934 ptarray_longitude_shift(ctx, line->points);
935 return;
936 case RTPOLYGONTYPE:
937 poly = (RTPOLY *)rtgeom;
938 for (i=0; i<poly->nrings; i++)
939 ptarray_longitude_shift(ctx, poly->rings[i]);
940 return;
941 case RTTRIANGLETYPE:
942 triangle = (RTTRIANGLE *)rtgeom;
943 ptarray_longitude_shift(ctx, triangle->points);
944 return;
945 case RTMULTIPOINTTYPE:
946 case RTMULTILINETYPE:
947 case RTMULTIPOLYGONTYPE:
948 case RTPOLYHEDRALSURFACETYPE:
949 case RTTINTYPE:
950 case RTCOLLECTIONTYPE:
951 coll = (RTCOLLECTION *)rtgeom;
952 for (i=0; i<coll->ngeoms; i++)
953 rtgeom_longitude_shift(ctx, coll->geoms[i]);
954 return;
955 default:
956 rterror(ctx, "rtgeom_longitude_shift: unsupported geom type: %s",
957 rttype_name(ctx, rtgeom->type));
958 }
959 }
960
961 int
rtgeom_is_closed(const RTCTX * ctx,const RTGEOM * geom)962 rtgeom_is_closed(const RTCTX *ctx, const RTGEOM *geom)
963 {
964 int type = geom->type;
965
966 if( rtgeom_is_empty(ctx, geom) )
967 return RT_FALSE;
968
969 /* Test linear types for closure */
970 switch (type)
971 {
972 case RTLINETYPE:
973 return rtline_is_closed(ctx, (RTLINE*)geom);
974 case RTPOLYGONTYPE:
975 return rtpoly_is_closed(ctx, (RTPOLY*)geom);
976 case RTCIRCSTRINGTYPE:
977 return rtcircstring_is_closed(ctx, (RTCIRCSTRING*)geom);
978 case RTCOMPOUNDTYPE:
979 return rtcompound_is_closed(ctx, (RTCOMPOUND*)geom);
980 case RTTINTYPE:
981 return rttin_is_closed(ctx, (RTTIN*)geom);
982 case RTPOLYHEDRALSURFACETYPE:
983 return rtpsurface_is_closed(ctx, (RTPSURFACE*)geom);
984 }
985
986 /* Recurse into collections and see if anything is not closed */
987 if ( rtgeom_is_collection(ctx, geom) )
988 {
989 RTCOLLECTION *col = rtgeom_as_rtcollection(ctx, geom);
990 int i;
991 int closed;
992 for ( i = 0; i < col->ngeoms; i++ )
993 {
994 closed = rtgeom_is_closed(ctx, col->geoms[i]);
995 if ( ! closed )
996 return RT_FALSE;
997 }
998 return RT_TRUE;
999 }
1000
1001 /* All non-linear non-collection types we will call closed */
1002 return RT_TRUE;
1003 }
1004
1005 int
rtgeom_is_collection(const RTCTX * ctx,const RTGEOM * geom)1006 rtgeom_is_collection(const RTCTX *ctx, const RTGEOM *geom)
1007 {
1008 if( ! geom ) return RT_FALSE;
1009 return rttype_is_collection(ctx, geom->type);
1010 }
1011
1012 /** Return TRUE if the geometry may contain sub-geometries, i.e. it is a MULTI* or COMPOUNDCURVE */
1013 int
rttype_is_collection(const RTCTX * ctx,uint8_t type)1014 rttype_is_collection(const RTCTX *ctx, uint8_t type)
1015 {
1016
1017 switch (type)
1018 {
1019 case RTMULTIPOINTTYPE:
1020 case RTMULTILINETYPE:
1021 case RTMULTIPOLYGONTYPE:
1022 case RTCOLLECTIONTYPE:
1023 case RTCURVEPOLYTYPE:
1024 case RTCOMPOUNDTYPE:
1025 case RTMULTICURVETYPE:
1026 case RTMULTISURFACETYPE:
1027 case RTPOLYHEDRALSURFACETYPE:
1028 case RTTINTYPE:
1029 return RT_TRUE;
1030 break;
1031
1032 default:
1033 return RT_FALSE;
1034 }
1035 }
1036
1037 /**
1038 * Given an rttype number, what homogeneous collection can hold it?
1039 */
1040 int
rttype_get_collectiontype(const RTCTX * ctx,uint8_t type)1041 rttype_get_collectiontype(const RTCTX *ctx, uint8_t type)
1042 {
1043 switch (type)
1044 {
1045 case RTPOINTTYPE:
1046 return RTMULTIPOINTTYPE;
1047 case RTLINETYPE:
1048 return RTMULTILINETYPE;
1049 case RTPOLYGONTYPE:
1050 return RTMULTIPOLYGONTYPE;
1051 case RTCIRCSTRINGTYPE:
1052 return RTMULTICURVETYPE;
1053 case RTCOMPOUNDTYPE:
1054 return RTMULTICURVETYPE;
1055 case RTCURVEPOLYTYPE:
1056 return RTMULTISURFACETYPE;
1057 case RTTRIANGLETYPE:
1058 return RTTINTYPE;
1059 default:
1060 return RTCOLLECTIONTYPE;
1061 }
1062 }
1063
1064
rtgeom_free(const RTCTX * ctx,RTGEOM * rtgeom)1065 void rtgeom_free(const RTCTX *ctx, RTGEOM *rtgeom)
1066 {
1067
1068 /* There's nothing here to free... */
1069 if( ! rtgeom ) return;
1070
1071 RTDEBUGF(ctx, 5,"freeing a %s",rttype_name(ctx, rtgeom->type));
1072
1073 switch (rtgeom->type)
1074 {
1075 case RTPOINTTYPE:
1076 rtpoint_free(ctx, (RTPOINT *)rtgeom);
1077 break;
1078 case RTLINETYPE:
1079 rtline_free(ctx, (RTLINE *)rtgeom);
1080 break;
1081 case RTPOLYGONTYPE:
1082 rtpoly_free(ctx, (RTPOLY *)rtgeom);
1083 break;
1084 case RTCIRCSTRINGTYPE:
1085 rtcircstring_free(ctx, (RTCIRCSTRING *)rtgeom);
1086 break;
1087 case RTTRIANGLETYPE:
1088 rttriangle_free(ctx, (RTTRIANGLE *)rtgeom);
1089 break;
1090 case RTMULTIPOINTTYPE:
1091 rtmpoint_free(ctx, (RTMPOINT *)rtgeom);
1092 break;
1093 case RTMULTILINETYPE:
1094 rtmline_free(ctx, (RTMLINE *)rtgeom);
1095 break;
1096 case RTMULTIPOLYGONTYPE:
1097 rtmpoly_free(ctx, (RTMPOLY *)rtgeom);
1098 break;
1099 case RTPOLYHEDRALSURFACETYPE:
1100 rtpsurface_free(ctx, (RTPSURFACE *)rtgeom);
1101 break;
1102 case RTTINTYPE:
1103 rttin_free(ctx, (RTTIN *)rtgeom);
1104 break;
1105 case RTCURVEPOLYTYPE:
1106 case RTCOMPOUNDTYPE:
1107 case RTMULTICURVETYPE:
1108 case RTMULTISURFACETYPE:
1109 case RTCOLLECTIONTYPE:
1110 rtcollection_free(ctx, (RTCOLLECTION *)rtgeom);
1111 break;
1112 default:
1113 rterror(ctx, "rtgeom_free called with unknown type (%d) %s", rtgeom->type, rttype_name(ctx, rtgeom->type));
1114 }
1115 return;
1116 }
1117
rtgeom_needs_bbox(const RTCTX * ctx,const RTGEOM * geom)1118 int rtgeom_needs_bbox(const RTCTX *ctx, const RTGEOM *geom)
1119 {
1120 assert(geom);
1121 if ( geom->type == RTPOINTTYPE )
1122 {
1123 return RT_FALSE;
1124 }
1125 else if ( geom->type == RTLINETYPE )
1126 {
1127 if ( rtgeom_count_vertices(ctx, geom) <= 2 )
1128 return RT_FALSE;
1129 else
1130 return RT_TRUE;
1131 }
1132 else if ( geom->type == RTMULTIPOINTTYPE )
1133 {
1134 if ( ((RTCOLLECTION*)geom)->ngeoms == 1 )
1135 return RT_FALSE;
1136 else
1137 return RT_TRUE;
1138 }
1139 else if ( geom->type == RTMULTILINETYPE )
1140 {
1141 if ( ((RTCOLLECTION*)geom)->ngeoms == 1 && rtgeom_count_vertices(ctx, geom) <= 2 )
1142 return RT_FALSE;
1143 else
1144 return RT_TRUE;
1145 }
1146 else
1147 {
1148 return RT_TRUE;
1149 }
1150 }
1151
1152 /**
1153 * Count points in an #RTGEOM.
1154 */
rtgeom_count_vertices(const RTCTX * ctx,const RTGEOM * geom)1155 int rtgeom_count_vertices(const RTCTX *ctx, const RTGEOM *geom)
1156 {
1157 int result = 0;
1158
1159 /* Null? Zero. */
1160 if( ! geom ) return 0;
1161
1162 RTDEBUGF(ctx, 4, "rtgeom_count_vertices got type %s",
1163 rttype_name(ctx, geom->type));
1164
1165 /* Empty? Zero. */
1166 if( rtgeom_is_empty(ctx, geom) ) return 0;
1167
1168 switch (geom->type)
1169 {
1170 case RTPOINTTYPE:
1171 result = 1;
1172 break;
1173 case RTTRIANGLETYPE:
1174 case RTCIRCSTRINGTYPE:
1175 case RTLINETYPE:
1176 result = rtline_count_vertices(ctx, (RTLINE *)geom);
1177 break;
1178 case RTPOLYGONTYPE:
1179 result = rtpoly_count_vertices(ctx, (RTPOLY *)geom);
1180 break;
1181 case RTCOMPOUNDTYPE:
1182 case RTCURVEPOLYTYPE:
1183 case RTMULTICURVETYPE:
1184 case RTMULTISURFACETYPE:
1185 case RTMULTIPOINTTYPE:
1186 case RTMULTILINETYPE:
1187 case RTMULTIPOLYGONTYPE:
1188 case RTPOLYHEDRALSURFACETYPE:
1189 case RTTINTYPE:
1190 case RTCOLLECTIONTYPE:
1191 result = rtcollection_count_vertices(ctx, (RTCOLLECTION *)geom);
1192 break;
1193 default:
1194 rterror(ctx, "%s: unsupported input geometry type: %s",
1195 __func__, rttype_name(ctx, geom->type));
1196 break;
1197 }
1198 RTDEBUGF(ctx, 3, "counted %d vertices", result);
1199 return result;
1200 }
1201
1202 /**
1203 * For an #RTGEOM, returns 0 for points, 1 for lines,
1204 * 2 for polygons, 3 for volume, and the max dimension
1205 * of a collection.
1206 */
rtgeom_dimension(const RTCTX * ctx,const RTGEOM * geom)1207 int rtgeom_dimension(const RTCTX *ctx, const RTGEOM *geom)
1208 {
1209
1210 /* Null? Zero. */
1211 if( ! geom ) return -1;
1212
1213 RTDEBUGF(ctx, 4, "rtgeom_dimension got type %s",
1214 rttype_name(ctx, geom->type));
1215
1216 /* Empty? Zero. */
1217 /* if( rtgeom_is_empty(ctx, geom) ) return 0; */
1218
1219 switch (geom->type)
1220 {
1221 case RTPOINTTYPE:
1222 case RTMULTIPOINTTYPE:
1223 return 0;
1224 case RTCIRCSTRINGTYPE:
1225 case RTLINETYPE:
1226 case RTCOMPOUNDTYPE:
1227 case RTMULTICURVETYPE:
1228 case RTMULTILINETYPE:
1229 return 1;
1230 case RTTRIANGLETYPE:
1231 case RTPOLYGONTYPE:
1232 case RTCURVEPOLYTYPE:
1233 case RTMULTISURFACETYPE:
1234 case RTMULTIPOLYGONTYPE:
1235 case RTTINTYPE:
1236 return 2;
1237 case RTPOLYHEDRALSURFACETYPE:
1238 {
1239 /* A closed polyhedral surface contains a volume. */
1240 int closed = rtpsurface_is_closed(ctx, (RTPSURFACE*)geom);
1241 return ( closed ? 3 : 2 );
1242 }
1243 case RTCOLLECTIONTYPE:
1244 {
1245 int maxdim = 0, i;
1246 RTCOLLECTION *col = (RTCOLLECTION*)geom;
1247 for( i = 0; i < col->ngeoms; i++ )
1248 {
1249 int dim = rtgeom_dimension(ctx, col->geoms[i]);
1250 maxdim = ( dim > maxdim ? dim : maxdim );
1251 }
1252 return maxdim;
1253 }
1254 default:
1255 rterror(ctx, "%s: unsupported input geometry type: %s",
1256 __func__, rttype_name(ctx, geom->type));
1257 }
1258 return -1;
1259 }
1260
1261 /**
1262 * Count rings in an #RTGEOM.
1263 */
rtgeom_count_rings(const RTCTX * ctx,const RTGEOM * geom)1264 int rtgeom_count_rings(const RTCTX *ctx, const RTGEOM *geom)
1265 {
1266 int result = 0;
1267
1268 /* Null? Empty? Zero. */
1269 if( ! geom || rtgeom_is_empty(ctx, geom) )
1270 return 0;
1271
1272 switch (geom->type)
1273 {
1274 case RTPOINTTYPE:
1275 case RTCIRCSTRINGTYPE:
1276 case RTCOMPOUNDTYPE:
1277 case RTMULTICURVETYPE:
1278 case RTMULTIPOINTTYPE:
1279 case RTMULTILINETYPE:
1280 case RTLINETYPE:
1281 result = 0;
1282 break;
1283 case RTTRIANGLETYPE:
1284 result = 1;
1285 break;
1286 case RTPOLYGONTYPE:
1287 result = ((RTPOLY *)geom)->nrings;
1288 break;
1289 case RTCURVEPOLYTYPE:
1290 result = ((RTCURVEPOLY *)geom)->nrings;
1291 break;
1292 case RTMULTISURFACETYPE:
1293 case RTMULTIPOLYGONTYPE:
1294 case RTPOLYHEDRALSURFACETYPE:
1295 case RTTINTYPE:
1296 case RTCOLLECTIONTYPE:
1297 {
1298 RTCOLLECTION *col = (RTCOLLECTION*)geom;
1299 int i = 0;
1300 for( i = 0; i < col->ngeoms; i++ )
1301 result += rtgeom_count_rings(ctx, col->geoms[i]);
1302 break;
1303 }
1304 default:
1305 rterror(ctx, "rtgeom_count_rings: unsupported input geometry type: %s", rttype_name(ctx, geom->type));
1306 break;
1307 }
1308 RTDEBUGF(ctx, 3, "counted %d rings", result);
1309 return result;
1310 }
1311
rtgeom_is_empty(const RTCTX * ctx,const RTGEOM * geom)1312 int rtgeom_is_empty(const RTCTX *ctx, const RTGEOM *geom)
1313 {
1314 int result = RT_FALSE;
1315 RTDEBUGF(ctx, 4, "rtgeom_is_empty: got type %s",
1316 rttype_name(ctx, geom->type));
1317
1318 switch (geom->type)
1319 {
1320 case RTPOINTTYPE:
1321 return rtpoint_is_empty(ctx, (RTPOINT*)geom);
1322 break;
1323 case RTLINETYPE:
1324 return rtline_is_empty(ctx, (RTLINE*)geom);
1325 break;
1326 case RTCIRCSTRINGTYPE:
1327 return rtcircstring_is_empty(ctx, (RTCIRCSTRING*)geom);
1328 break;
1329 case RTPOLYGONTYPE:
1330 return rtpoly_is_empty(ctx, (RTPOLY*)geom);
1331 break;
1332 case RTTRIANGLETYPE:
1333 return rttriangle_is_empty(ctx, (RTTRIANGLE*)geom);
1334 break;
1335 case RTMULTIPOINTTYPE:
1336 case RTMULTILINETYPE:
1337 case RTMULTIPOLYGONTYPE:
1338 case RTCOMPOUNDTYPE:
1339 case RTCURVEPOLYTYPE:
1340 case RTMULTICURVETYPE:
1341 case RTMULTISURFACETYPE:
1342 case RTPOLYHEDRALSURFACETYPE:
1343 case RTTINTYPE:
1344 case RTCOLLECTIONTYPE:
1345 return rtcollection_is_empty(ctx, (RTCOLLECTION *)geom);
1346 break;
1347 default:
1348 rterror(ctx, "rtgeom_is_empty: unsupported input geometry type: %s",
1349 rttype_name(ctx, geom->type));
1350 break;
1351 }
1352 return result;
1353 }
1354
rtgeom_has_srid(const RTCTX * ctx,const RTGEOM * geom)1355 int rtgeom_has_srid(const RTCTX *ctx, const RTGEOM *geom)
1356 {
1357 if ( geom->srid != SRID_UNKNOWN )
1358 return RT_TRUE;
1359
1360 return RT_FALSE;
1361 }
1362
1363
rtcollection_dimensionality(const RTCTX * ctx,RTCOLLECTION * col)1364 static int rtcollection_dimensionality(const RTCTX *ctx, RTCOLLECTION *col)
1365 {
1366 int i;
1367 int dimensionality = 0;
1368 for ( i = 0; i < col->ngeoms; i++ )
1369 {
1370 int d = rtgeom_dimensionality(ctx, col->geoms[i]);
1371 if ( d > dimensionality )
1372 dimensionality = d;
1373 }
1374 return dimensionality;
1375 }
1376
rtgeom_dimensionality(const RTCTX * ctx,RTGEOM * geom)1377 extern int rtgeom_dimensionality(const RTCTX *ctx, RTGEOM *geom)
1378 {
1379 int dim;
1380
1381 RTDEBUGF(ctx, 3, "rtgeom_dimensionality got type %s",
1382 rttype_name(ctx, geom->type));
1383
1384 switch (geom->type)
1385 {
1386 case RTPOINTTYPE:
1387 case RTMULTIPOINTTYPE:
1388 return 0;
1389 break;
1390 case RTLINETYPE:
1391 case RTCIRCSTRINGTYPE:
1392 case RTMULTILINETYPE:
1393 case RTCOMPOUNDTYPE:
1394 case RTMULTICURVETYPE:
1395 return 1;
1396 break;
1397 case RTPOLYGONTYPE:
1398 case RTTRIANGLETYPE:
1399 case RTCURVEPOLYTYPE:
1400 case RTMULTIPOLYGONTYPE:
1401 case RTMULTISURFACETYPE:
1402 return 2;
1403 break;
1404
1405 case RTPOLYHEDRALSURFACETYPE:
1406 case RTTINTYPE:
1407 dim = rtgeom_is_closed(ctx, geom)?3:2;
1408 return dim;
1409 break;
1410
1411 case RTCOLLECTIONTYPE:
1412 return rtcollection_dimensionality(ctx, (RTCOLLECTION *)geom);
1413 break;
1414 default:
1415 rterror(ctx, "rtgeom_dimensionality: unsupported input geometry type: %s",
1416 rttype_name(ctx, geom->type));
1417 break;
1418 }
1419 return 0;
1420 }
1421
rtgeom_remove_repeated_points(const RTCTX * ctx,const RTGEOM * in,double tolerance)1422 extern RTGEOM* rtgeom_remove_repeated_points(const RTCTX *ctx, const RTGEOM *in, double tolerance)
1423 {
1424 RTDEBUGF(ctx, 4, "rtgeom_remove_repeated_points got type %s",
1425 rttype_name(ctx, in->type));
1426
1427 if(rtgeom_is_empty(ctx, in))
1428 {
1429 return rtgeom_clone_deep(ctx, in);
1430 }
1431
1432 switch (in->type)
1433 {
1434 case RTMULTIPOINTTYPE:
1435 return rtmpoint_remove_repeated_points(ctx, (RTMPOINT*)in, tolerance);
1436 break;
1437 case RTLINETYPE:
1438 return rtline_remove_repeated_points(ctx, (RTLINE*)in, tolerance);
1439
1440 case RTMULTILINETYPE:
1441 case RTCOLLECTIONTYPE:
1442 case RTMULTIPOLYGONTYPE:
1443 case RTPOLYHEDRALSURFACETYPE:
1444 return rtcollection_remove_repeated_points(ctx, (RTCOLLECTION *)in, tolerance);
1445
1446 case RTPOLYGONTYPE:
1447 return rtpoly_remove_repeated_points(ctx, (RTPOLY *)in, tolerance);
1448 break;
1449
1450 case RTPOINTTYPE:
1451 case RTTRIANGLETYPE:
1452 case RTTINTYPE:
1453 /* No point is repeated for a single point, or for Triangle or TIN */
1454 return rtgeom_clone_deep(ctx, in);
1455
1456 case RTCIRCSTRINGTYPE:
1457 case RTCOMPOUNDTYPE:
1458 case RTMULTICURVETYPE:
1459 case RTCURVEPOLYTYPE:
1460 case RTMULTISURFACETYPE:
1461 /* Dunno how to handle these, will return untouched */
1462 return rtgeom_clone_deep(ctx, in);
1463
1464 default:
1465 rtnotice(ctx, "%s: unsupported geometry type: %s",
1466 __func__, rttype_name(ctx, in->type));
1467 return rtgeom_clone_deep(ctx, in);
1468 break;
1469 }
1470 return 0;
1471 }
1472
rtgeom_flip_coordinates(const RTCTX * ctx,RTGEOM * in)1473 RTGEOM* rtgeom_flip_coordinates(const RTCTX *ctx, RTGEOM *in)
1474 {
1475 rtgeom_swap_ordinates(ctx, in,RTORD_X,RTORD_Y);
1476 return in;
1477 }
1478
rtgeom_swap_ordinates(const RTCTX * ctx,RTGEOM * in,RTORD o1,RTORD o2)1479 void rtgeom_swap_ordinates(const RTCTX *ctx, RTGEOM *in, RTORD o1, RTORD o2)
1480 {
1481 RTCOLLECTION *col;
1482 RTPOLY *poly;
1483 int i;
1484
1485 #if PARANOIA_LEVEL > 0
1486 assert(o1 < 4);
1487 assert(o2 < 4);
1488 #endif
1489
1490 if ( (!in) || rtgeom_is_empty(ctx, in) ) return;
1491
1492 /* TODO: check for rtgeom NOT having the specified dimension ? */
1493
1494 RTDEBUGF(ctx, 4, "rtgeom_flip_coordinates, got type: %s",
1495 rttype_name(ctx, in->type));
1496
1497 switch (in->type)
1498 {
1499 case RTPOINTTYPE:
1500 ptarray_swap_ordinates(ctx, rtgeom_as_rtpoint(ctx, in)->point, o1, o2);
1501 break;
1502
1503 case RTLINETYPE:
1504 ptarray_swap_ordinates(ctx, rtgeom_as_rtline(ctx, in)->points, o1, o2);
1505 break;
1506
1507 case RTCIRCSTRINGTYPE:
1508 ptarray_swap_ordinates(ctx, rtgeom_as_rtcircstring(ctx, in)->points, o1, o2);
1509 break;
1510
1511 case RTPOLYGONTYPE:
1512 poly = (RTPOLY *) in;
1513 for (i=0; i<poly->nrings; i++)
1514 {
1515 ptarray_swap_ordinates(ctx, poly->rings[i], o1, o2);
1516 }
1517 break;
1518
1519 case RTTRIANGLETYPE:
1520 ptarray_swap_ordinates(ctx, rtgeom_as_rttriangle(ctx, in)->points, o1, o2);
1521 break;
1522
1523 case RTMULTIPOINTTYPE:
1524 case RTMULTILINETYPE:
1525 case RTMULTIPOLYGONTYPE:
1526 case RTCOLLECTIONTYPE:
1527 case RTCOMPOUNDTYPE:
1528 case RTCURVEPOLYTYPE:
1529 case RTMULTISURFACETYPE:
1530 case RTMULTICURVETYPE:
1531 case RTPOLYHEDRALSURFACETYPE:
1532 case RTTINTYPE:
1533 col = (RTCOLLECTION *) in;
1534 for (i=0; i<col->ngeoms; i++)
1535 {
1536 rtgeom_swap_ordinates(ctx, col->geoms[i], o1, o2);
1537 }
1538 break;
1539
1540 default:
1541 rterror(ctx, "rtgeom_swap_ordinates: unsupported geometry type: %s",
1542 rttype_name(ctx, in->type));
1543 return;
1544 }
1545
1546 /* only refresh bbox if X or Y changed */
1547 if ( in->bbox && (o1 < 2 || o2 < 2) )
1548 {
1549 rtgeom_drop_bbox(ctx, in);
1550 rtgeom_add_bbox(ctx, in);
1551 }
1552 }
1553
rtgeom_set_srid(const RTCTX * ctx,RTGEOM * geom,int32_t srid)1554 void rtgeom_set_srid(const RTCTX *ctx, RTGEOM *geom, int32_t srid)
1555 {
1556 int i;
1557
1558 RTDEBUGF(ctx, 4,"entered with srid=%d",srid);
1559
1560 geom->srid = srid;
1561
1562 if ( rtgeom_is_collection(ctx, geom) )
1563 {
1564 /* All the children are set to the same SRID value */
1565 RTCOLLECTION *col = rtgeom_as_rtcollection(ctx, geom);
1566 for ( i = 0; i < col->ngeoms; i++ )
1567 {
1568 rtgeom_set_srid(ctx, col->geoms[i], srid);
1569 }
1570 }
1571 }
1572
rtgeom_simplify(const RTCTX * ctx,const RTGEOM * igeom,double dist,int preserve_collapsed)1573 RTGEOM* rtgeom_simplify(const RTCTX *ctx, const RTGEOM *igeom, double dist, int preserve_collapsed)
1574 {
1575 switch (igeom->type)
1576 {
1577 case RTPOINTTYPE:
1578 case RTMULTIPOINTTYPE:
1579 return rtgeom_clone(ctx, igeom);
1580 case RTLINETYPE:
1581 return (RTGEOM*)rtline_simplify(ctx, (RTLINE*)igeom, dist, preserve_collapsed);
1582 case RTPOLYGONTYPE:
1583 return (RTGEOM*)rtpoly_simplify(ctx, (RTPOLY*)igeom, dist, preserve_collapsed);
1584 case RTMULTILINETYPE:
1585 case RTMULTIPOLYGONTYPE:
1586 case RTCOLLECTIONTYPE:
1587 return (RTGEOM*)rtcollection_simplify(ctx, (RTCOLLECTION *)igeom, dist, preserve_collapsed);
1588 default:
1589 rterror(ctx, "%s: unsupported geometry type: %s", __func__, rttype_name(ctx, igeom->type));
1590 }
1591 return NULL;
1592 }
1593
rtgeom_area(const RTCTX * ctx,const RTGEOM * geom)1594 double rtgeom_area(const RTCTX *ctx, const RTGEOM *geom)
1595 {
1596 int type = geom->type;
1597
1598 if ( type == RTPOLYGONTYPE )
1599 return rtpoly_area(ctx, (RTPOLY*)geom);
1600 else if ( type == RTCURVEPOLYTYPE )
1601 return rtcurvepoly_area(ctx, (RTCURVEPOLY*)geom);
1602 else if (type == RTTRIANGLETYPE )
1603 return rttriangle_area(ctx, (RTTRIANGLE*)geom);
1604 else if ( rtgeom_is_collection(ctx, geom) )
1605 {
1606 double area = 0.0;
1607 int i;
1608 RTCOLLECTION *col = (RTCOLLECTION*)geom;
1609 for ( i = 0; i < col->ngeoms; i++ )
1610 area += rtgeom_area(ctx, col->geoms[i]);
1611 return area;
1612 }
1613 else
1614 return 0.0;
1615 }
1616
rtgeom_perimeter(const RTCTX * ctx,const RTGEOM * geom)1617 double rtgeom_perimeter(const RTCTX *ctx, const RTGEOM *geom)
1618 {
1619 int type = geom->type;
1620 if ( type == RTPOLYGONTYPE )
1621 return rtpoly_perimeter(ctx, (RTPOLY*)geom);
1622 else if ( type == RTCURVEPOLYTYPE )
1623 return rtcurvepoly_perimeter(ctx, (RTCURVEPOLY*)geom);
1624 else if ( type == RTTRIANGLETYPE )
1625 return rttriangle_perimeter(ctx, (RTTRIANGLE*)geom);
1626 else if ( rtgeom_is_collection(ctx, geom) )
1627 {
1628 double perimeter = 0.0;
1629 int i;
1630 RTCOLLECTION *col = (RTCOLLECTION*)geom;
1631 for ( i = 0; i < col->ngeoms; i++ )
1632 perimeter += rtgeom_perimeter(ctx, col->geoms[i]);
1633 return perimeter;
1634 }
1635 else
1636 return 0.0;
1637 }
1638
rtgeom_perimeter_2d(const RTCTX * ctx,const RTGEOM * geom)1639 double rtgeom_perimeter_2d(const RTCTX *ctx, const RTGEOM *geom)
1640 {
1641 int type = geom->type;
1642 if ( type == RTPOLYGONTYPE )
1643 return rtpoly_perimeter_2d(ctx, (RTPOLY*)geom);
1644 else if ( type == RTCURVEPOLYTYPE )
1645 return rtcurvepoly_perimeter_2d(ctx, (RTCURVEPOLY*)geom);
1646 else if ( type == RTTRIANGLETYPE )
1647 return rttriangle_perimeter_2d(ctx, (RTTRIANGLE*)geom);
1648 else if ( rtgeom_is_collection(ctx, geom) )
1649 {
1650 double perimeter = 0.0;
1651 int i;
1652 RTCOLLECTION *col = (RTCOLLECTION*)geom;
1653 for ( i = 0; i < col->ngeoms; i++ )
1654 perimeter += rtgeom_perimeter_2d(ctx, col->geoms[i]);
1655 return perimeter;
1656 }
1657 else
1658 return 0.0;
1659 }
1660
rtgeom_length(const RTCTX * ctx,const RTGEOM * geom)1661 double rtgeom_length(const RTCTX *ctx, const RTGEOM *geom)
1662 {
1663 int type = geom->type;
1664 if ( type == RTLINETYPE )
1665 return rtline_length(ctx, (RTLINE*)geom);
1666 else if ( type == RTCIRCSTRINGTYPE )
1667 return rtcircstring_length(ctx, (RTCIRCSTRING*)geom);
1668 else if ( type == RTCOMPOUNDTYPE )
1669 return rtcompound_length(ctx, (RTCOMPOUND*)geom);
1670 else if ( rtgeom_is_collection(ctx, geom) )
1671 {
1672 double length = 0.0;
1673 int i;
1674 RTCOLLECTION *col = (RTCOLLECTION*)geom;
1675 for ( i = 0; i < col->ngeoms; i++ )
1676 length += rtgeom_length(ctx, col->geoms[i]);
1677 return length;
1678 }
1679 else
1680 return 0.0;
1681 }
1682
rtgeom_length_2d(const RTCTX * ctx,const RTGEOM * geom)1683 double rtgeom_length_2d(const RTCTX *ctx, const RTGEOM *geom)
1684 {
1685 int type = geom->type;
1686 if ( type == RTLINETYPE )
1687 return rtline_length_2d(ctx, (RTLINE*)geom);
1688 else if ( type == RTCIRCSTRINGTYPE )
1689 return rtcircstring_length_2d(ctx, (RTCIRCSTRING*)geom);
1690 else if ( type == RTCOMPOUNDTYPE )
1691 return rtcompound_length_2d(ctx, (RTCOMPOUND*)geom);
1692 else if ( rtgeom_is_collection(ctx, geom) )
1693 {
1694 double length = 0.0;
1695 int i;
1696 RTCOLLECTION *col = (RTCOLLECTION*)geom;
1697 for ( i = 0; i < col->ngeoms; i++ )
1698 length += rtgeom_length_2d(ctx, col->geoms[i]);
1699 return length;
1700 }
1701 else
1702 return 0.0;
1703 }
1704
1705 void
rtgeom_affine(const RTCTX * ctx,RTGEOM * geom,const RTAFFINE * affine)1706 rtgeom_affine(const RTCTX *ctx, RTGEOM *geom, const RTAFFINE *affine)
1707 {
1708 int type = geom->type;
1709 int i;
1710
1711 switch(type)
1712 {
1713 /* Take advantage of fact tht pt/ln/circ/tri have same memory structure */
1714 case RTPOINTTYPE:
1715 case RTLINETYPE:
1716 case RTCIRCSTRINGTYPE:
1717 case RTTRIANGLETYPE:
1718 {
1719 RTLINE *l = (RTLINE*)geom;
1720 ptarray_affine(ctx, l->points, affine);
1721 break;
1722 }
1723 case RTPOLYGONTYPE:
1724 {
1725 RTPOLY *p = (RTPOLY*)geom;
1726 for( i = 0; i < p->nrings; i++ )
1727 ptarray_affine(ctx, p->rings[i], affine);
1728 break;
1729 }
1730 case RTCURVEPOLYTYPE:
1731 {
1732 RTCURVEPOLY *c = (RTCURVEPOLY*)geom;
1733 for( i = 0; i < c->nrings; i++ )
1734 rtgeom_affine(ctx, c->rings[i], affine);
1735 break;
1736 }
1737 default:
1738 {
1739 if( rtgeom_is_collection(ctx, geom) )
1740 {
1741 RTCOLLECTION *c = (RTCOLLECTION*)geom;
1742 for( i = 0; i < c->ngeoms; i++ )
1743 {
1744 rtgeom_affine(ctx, c->geoms[i], affine);
1745 }
1746 }
1747 else
1748 {
1749 rterror(ctx, "rtgeom_affine: unable to handle type '%s'", rttype_name(ctx, type));
1750 }
1751 }
1752 }
1753
1754 }
1755
1756 void
rtgeom_scale(const RTCTX * ctx,RTGEOM * geom,const RTPOINT4D * factor)1757 rtgeom_scale(const RTCTX *ctx, RTGEOM *geom, const RTPOINT4D *factor)
1758 {
1759 int type = geom->type;
1760 int i;
1761
1762 switch(type)
1763 {
1764 /* Take advantage of fact tht pt/ln/circ/tri have same memory structure */
1765 case RTPOINTTYPE:
1766 case RTLINETYPE:
1767 case RTCIRCSTRINGTYPE:
1768 case RTTRIANGLETYPE:
1769 {
1770 RTLINE *l = (RTLINE*)geom;
1771 ptarray_scale(ctx, l->points, factor);
1772 break;
1773 }
1774 case RTPOLYGONTYPE:
1775 {
1776 RTPOLY *p = (RTPOLY*)geom;
1777 for( i = 0; i < p->nrings; i++ )
1778 ptarray_scale(ctx, p->rings[i], factor);
1779 break;
1780 }
1781 case RTCURVEPOLYTYPE:
1782 {
1783 RTCURVEPOLY *c = (RTCURVEPOLY*)geom;
1784 for( i = 0; i < c->nrings; i++ )
1785 rtgeom_scale(ctx, c->rings[i], factor);
1786 break;
1787 }
1788 default:
1789 {
1790 if( rtgeom_is_collection(ctx, geom) )
1791 {
1792 RTCOLLECTION *c = (RTCOLLECTION*)geom;
1793 for( i = 0; i < c->ngeoms; i++ )
1794 {
1795 rtgeom_scale(ctx, c->geoms[i], factor);
1796 }
1797 }
1798 else
1799 {
1800 rterror(ctx, "rtgeom_scale: unable to handle type '%s'", rttype_name(ctx, type));
1801 }
1802 }
1803 }
1804
1805 /* Recompute bbox if needed */
1806
1807 if ( geom->bbox )
1808 {
1809 /* TODO: expose a gbox_scale function */
1810 geom->bbox->xmin *= factor->x;
1811 geom->bbox->xmax *= factor->x;
1812 geom->bbox->ymin *= factor->y;
1813 geom->bbox->ymax *= factor->y;
1814 geom->bbox->zmin *= factor->z;
1815 geom->bbox->zmax *= factor->z;
1816 geom->bbox->mmin *= factor->m;
1817 geom->bbox->mmax *= factor->m;
1818 }
1819 }
1820
1821 RTGEOM*
rtgeom_construct_empty(const RTCTX * ctx,uint8_t type,int srid,char hasz,char hasm)1822 rtgeom_construct_empty(const RTCTX *ctx, uint8_t type, int srid, char hasz, char hasm)
1823 {
1824 switch(type)
1825 {
1826 case RTPOINTTYPE:
1827 return rtpoint_as_rtgeom(ctx, rtpoint_construct_empty(ctx, srid, hasz, hasm));
1828 case RTLINETYPE:
1829 return rtline_as_rtgeom(ctx, rtline_construct_empty(ctx, srid, hasz, hasm));
1830 case RTPOLYGONTYPE:
1831 return rtpoly_as_rtgeom(ctx, rtpoly_construct_empty(ctx, srid, hasz, hasm));
1832 case RTCURVEPOLYTYPE:
1833 return rtcurvepoly_as_rtgeom(ctx, rtcurvepoly_construct_empty(ctx, srid, hasz, hasm));
1834 case RTCIRCSTRINGTYPE:
1835 return rtcircstring_as_rtgeom(ctx, rtcircstring_construct_empty(ctx, srid, hasz, hasm));
1836 case RTTRIANGLETYPE:
1837 return rttriangle_as_rtgeom(ctx, rttriangle_construct_empty(ctx, srid, hasz, hasm));
1838 case RTCOMPOUNDTYPE:
1839 case RTMULTIPOINTTYPE:
1840 case RTMULTILINETYPE:
1841 case RTMULTIPOLYGONTYPE:
1842 case RTCOLLECTIONTYPE:
1843 return rtcollection_as_rtgeom(ctx, rtcollection_construct_empty(ctx, type, srid, hasz, hasm));
1844 default:
1845 rterror(ctx, "rtgeom_construct_empty: unsupported geometry type: %s",
1846 rttype_name(ctx, type));
1847 return NULL;
1848 }
1849 }
1850
1851 int
rtgeom_startpoint(const RTCTX * ctx,const RTGEOM * rtgeom,RTPOINT4D * pt)1852 rtgeom_startpoint(const RTCTX *ctx, const RTGEOM* rtgeom, RTPOINT4D* pt)
1853 {
1854 if ( ! rtgeom )
1855 return RT_FAILURE;
1856
1857 switch( rtgeom->type )
1858 {
1859 case RTPOINTTYPE:
1860 return ptarray_startpoint(ctx, ((RTPOINT*)rtgeom)->point, pt);
1861 case RTTRIANGLETYPE:
1862 case RTCIRCSTRINGTYPE:
1863 case RTLINETYPE:
1864 return ptarray_startpoint(ctx, ((RTLINE*)rtgeom)->points, pt);
1865 case RTPOLYGONTYPE:
1866 return rtpoly_startpoint(ctx, (RTPOLY*)rtgeom, pt);
1867 case RTCURVEPOLYTYPE:
1868 case RTCOMPOUNDTYPE:
1869 case RTMULTIPOINTTYPE:
1870 case RTMULTILINETYPE:
1871 case RTMULTIPOLYGONTYPE:
1872 case RTCOLLECTIONTYPE:
1873 return rtcollection_startpoint(ctx, (RTCOLLECTION*)rtgeom, pt);
1874 default:
1875 rterror(ctx, "int: unsupported geometry type: %s",
1876 rttype_name(ctx, rtgeom->type));
1877 return RT_FAILURE;
1878 }
1879 }
1880
1881
1882 RTGEOM *
rtgeom_grid(const RTCTX * ctx,const RTGEOM * rtgeom,const gridspec * grid)1883 rtgeom_grid(const RTCTX *ctx, const RTGEOM *rtgeom, const gridspec *grid)
1884 {
1885 switch ( rtgeom->type )
1886 {
1887 case RTPOINTTYPE:
1888 return (RTGEOM *)rtpoint_grid(ctx, (RTPOINT *)rtgeom, grid);
1889 case RTLINETYPE:
1890 return (RTGEOM *)rtline_grid(ctx, (RTLINE *)rtgeom, grid);
1891 case RTPOLYGONTYPE:
1892 return (RTGEOM *)rtpoly_grid(ctx, (RTPOLY *)rtgeom, grid);
1893 case RTMULTIPOINTTYPE:
1894 case RTMULTILINETYPE:
1895 case RTMULTIPOLYGONTYPE:
1896 case RTCOLLECTIONTYPE:
1897 case RTCOMPOUNDTYPE:
1898 return (RTGEOM *)rtcollection_grid(ctx, (RTCOLLECTION *)rtgeom, grid);
1899 case RTCIRCSTRINGTYPE:
1900 return (RTGEOM *)rtcircstring_grid(ctx, (RTCIRCSTRING *)rtgeom, grid);
1901 default:
1902 rterror(ctx, "rtgeom_grid: Unsupported geometry type: %s",
1903 rttype_name(ctx, rtgeom->type));
1904 return NULL;
1905 }
1906 }
1907
1908
1909 /* Prototype for recursion */
1910 static int
1911 rtgeom_subdivide_recursive(const RTCTX *ctx, const RTGEOM *geom, int maxvertices, int depth, RTCOLLECTION *col, const RTGBOX *clip);
1912
1913 static int
rtgeom_subdivide_recursive(const RTCTX * ctx,const RTGEOM * geom,int maxvertices,int depth,RTCOLLECTION * col,const RTGBOX * clip)1914 rtgeom_subdivide_recursive(const RTCTX *ctx, const RTGEOM *geom, int maxvertices, int depth, RTCOLLECTION *col, const RTGBOX *clip)
1915 {
1916 const int maxdepth = 50;
1917 int nvertices = 0;
1918 int i, n = 0;
1919 double width = clip->xmax - clip->xmin;
1920 double height = clip->ymax - clip->ymin;
1921 RTGBOX subbox1, subbox2;
1922 RTGEOM *clipped1, *clipped2;
1923
1924 if ( geom->type == RTPOLYHEDRALSURFACETYPE || geom->type == RTTINTYPE )
1925 {
1926 rterror(ctx, "%s: unsupported geometry type '%s'", __func__, rttype_name(ctx, geom->type));
1927 }
1928
1929 if ( width == 0.0 && height == 0.0 )
1930 return 0;
1931
1932 /* Artays just recurse into collections */
1933 if ( rtgeom_is_collection(ctx, geom) )
1934 {
1935 RTCOLLECTION *incol = (RTCOLLECTION*)geom;
1936 int n = 0;
1937 for ( i = 0; i < incol->ngeoms; i++ )
1938 {
1939 /* Don't increment depth yet, since we aren't actually subdividing geomtries yet */
1940 n += rtgeom_subdivide_recursive(ctx, incol->geoms[i], maxvertices, depth, col, clip);
1941 }
1942 return n;
1943 }
1944
1945 /* But don't go too far. 2^25 = 33M, that's enough subdivision */
1946 /* Signal callers above that we depth'ed out with a negative */
1947 /* return value */
1948 if ( depth > maxdepth )
1949 {
1950 return 0;
1951 }
1952
1953 nvertices = rtgeom_count_vertices(ctx, geom);
1954 /* Skip empties entirely */
1955 if ( nvertices == 0 )
1956 {
1957 return 0;
1958 }
1959
1960 /* If it is under the vertex tolerance, just add it, we're done */
1961 if ( nvertices < maxvertices )
1962 {
1963 rtcollection_add_rtgeom(ctx, col, rtgeom_clone_deep(ctx, geom));
1964 return 1;
1965 }
1966
1967 subbox1 = subbox2 = *clip;
1968 if ( width > height )
1969 {
1970 subbox1.xmax = subbox2.xmin = (clip->xmin + clip->xmax)/2;
1971 }
1972 else
1973 {
1974 subbox1.ymax = subbox2.ymin = (clip->ymin + clip->ymax)/2;
1975 }
1976
1977 if ( height == 0 )
1978 {
1979 subbox1.ymax += FP_TOLERANCE;
1980 subbox2.ymax += FP_TOLERANCE;
1981 subbox1.ymin -= FP_TOLERANCE;
1982 subbox2.ymin -= FP_TOLERANCE;
1983 }
1984
1985 if ( width == 0 )
1986 {
1987 subbox1.xmax += FP_TOLERANCE;
1988 subbox2.xmax += FP_TOLERANCE;
1989 subbox1.xmin -= FP_TOLERANCE;
1990 subbox2.xmin -= FP_TOLERANCE;
1991 }
1992
1993 clipped1 = rtgeom_clip_by_rect(ctx, geom, subbox1.xmin, subbox1.ymin, subbox1.xmax, subbox1.ymax);
1994 clipped2 = rtgeom_clip_by_rect(ctx, geom, subbox2.xmin, subbox2.ymin, subbox2.xmax, subbox2.ymax);
1995
1996 if ( clipped1 )
1997 {
1998 n += rtgeom_subdivide_recursive(ctx, clipped1, maxvertices, ++depth, col, &subbox1);
1999 rtgeom_free(ctx, clipped1);
2000 }
2001
2002 if ( clipped2 )
2003 {
2004 n += rtgeom_subdivide_recursive(ctx, clipped2, maxvertices, ++depth, col, &subbox2);
2005 rtgeom_free(ctx, clipped2);
2006 }
2007
2008 return n;
2009
2010 }
2011
2012 RTCOLLECTION *
rtgeom_subdivide(const RTCTX * ctx,const RTGEOM * geom,int maxvertices)2013 rtgeom_subdivide(const RTCTX *ctx, const RTGEOM *geom, int maxvertices)
2014 {
2015 static int startdepth = 0;
2016 static int minmaxvertices = 8;
2017 RTCOLLECTION *col;
2018 RTGBOX clip;
2019
2020 col = rtcollection_construct_empty(ctx, RTCOLLECTIONTYPE, geom->srid, rtgeom_has_z(ctx, geom), rtgeom_has_m(ctx, geom));
2021
2022 if ( rtgeom_is_empty(ctx, geom) )
2023 return col;
2024
2025 if ( maxvertices < minmaxvertices )
2026 {
2027 rtcollection_free(ctx, col);
2028 rterror(ctx, "%s: cannot subdivide to fewer than %d vertices per output", __func__, minmaxvertices);
2029 }
2030
2031 clip = *(rtgeom_get_bbox(ctx, geom));
2032 rtgeom_subdivide_recursive(ctx, geom, maxvertices, startdepth, col, &clip);
2033 rtgeom_set_srid(ctx, (RTGEOM*)col, geom->srid);
2034 return col;
2035 }
2036
2037
2038 int
rtgeom_is_trajectory(const RTCTX * ctx,const RTGEOM * geom)2039 rtgeom_is_trajectory(const RTCTX *ctx, const RTGEOM *geom)
2040 {
2041 int type = geom->type;
2042
2043 if( type != RTLINETYPE )
2044 {
2045 rtnotice(ctx, "Geometry is not a LINESTRING");
2046 return RT_FALSE;
2047 }
2048 return rtline_is_trajectory(ctx, (RTLINE*)geom);
2049 }
2050
2051