1 : /******************************************************************************
2 : * $Id: ogrgeometryfactory.cpp 19992 2010-07-08 18:24:09Z rouault $
3 : *
4 : * Project: OpenGIS Simple Features Reference Implementation
5 : * Purpose: Factory for converting geometry to and from well known binary
6 : * format.
7 : * Author: Frank Warmerdam, warmerdam@pobox.com
8 : *
9 : ******************************************************************************
10 : * Copyright (c) 1999, Frank Warmerdam
11 : *
12 : * Permission is hereby granted, free of charge, to any person obtaining a
13 : * copy of this software and associated documentation files (the "Software"),
14 : * to deal in the Software without restriction, including without limitation
15 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
16 : * and/or sell copies of the Software, and to permit persons to whom the
17 : * Software is furnished to do so, subject to the following conditions:
18 : *
19 : * The above copyright notice and this permission notice shall be included
20 : * in all copies or substantial portions of the Software.
21 : *
22 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
23 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
24 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
25 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
26 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
27 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
28 : * DEALINGS IN THE SOFTWARE.
29 : ****************************************************************************/
30 :
31 : #include "ogr_geometry.h"
32 : #include "ogr_api.h"
33 : #include "ogr_p.h"
34 : #include "ogr_geos.h"
35 :
36 : CPL_CVSID("$Id: ogrgeometryfactory.cpp 19992 2010-07-08 18:24:09Z rouault $");
37 :
38 : #ifndef PI
39 : #define PI 3.14159265358979323846
40 : #endif
41 :
42 : /************************************************************************/
43 : /* createFromWkb() */
44 : /************************************************************************/
45 :
46 : /**
47 : * \brief Create a geometry object of the appropriate type from it's well known binary representation.
48 : *
49 : * Note that if nBytes is passed as zero, no checking can be done on whether
50 : * the pabyData is sufficient. This can result in a crash if the input
51 : * data is corrupt. This function returns no indication of the number of
52 : * bytes from the data source actually used to represent the returned
53 : * geometry object. Use OGRGeometry::WkbSize() on the returned geometry to
54 : * establish the number of bytes it required in WKB format.
55 : *
56 : * Also note that this is a static method, and that there
57 : * is no need to instantiate an OGRGeometryFactory object.
58 : *
59 : * The C function OGR_G_CreateFromWkb() is the same as this method.
60 : *
61 : * @param pabyData pointer to the input BLOB data.
62 : * @param poSR pointer to the spatial reference to be assigned to the
63 : * created geometry object. This may be NULL.
64 : * @param ppoReturn the newly created geometry object will be assigned to the
65 : * indicated pointer on return. This will be NULL in case
66 : * of failure.
67 : * @param nBytes the number of bytes available in pabyData, or -1 if it isn't
68 : * known.
69 : *
70 : * @return OGRERR_NONE if all goes well, otherwise any of
71 : * OGRERR_NOT_ENOUGH_DATA, OGRERR_UNSUPPORTED_GEOMETRY_TYPE, or
72 : * OGRERR_CORRUPT_DATA may be returned.
73 : */
74 :
75 : OGRErr OGRGeometryFactory::createFromWkb(unsigned char *pabyData,
76 : OGRSpatialReference * poSR,
77 : OGRGeometry **ppoReturn,
78 1768 : int nBytes )
79 :
80 : {
81 : OGRwkbGeometryType eGeometryType;
82 : OGRwkbByteOrder eByteOrder;
83 : OGRErr eErr;
84 : OGRGeometry *poGeom;
85 :
86 1768 : *ppoReturn = NULL;
87 :
88 1768 : if( nBytes < 9 && nBytes != -1 )
89 1020 : return OGRERR_NOT_ENOUGH_DATA;
90 :
91 : /* -------------------------------------------------------------------- */
92 : /* Get the byte order byte. The extra tests are to work around */
93 : /* bug sin the WKB of DB2 v7.2 as identified by Safe Software. */
94 : /* -------------------------------------------------------------------- */
95 748 : eByteOrder = DB2_V72_FIX_BYTE_ORDER((OGRwkbByteOrder) *pabyData);
96 :
97 :
98 748 : if( eByteOrder != wkbXDR && eByteOrder != wkbNDR )
99 : {
100 : CPLDebug( "OGR",
101 : "OGRGeometryFactory::createFromWkb() - got corrupt data.\n"
102 : "%02X%02X%02X%02X%02X%02X%02X%02X%02X\n",
103 : pabyData[0],
104 : pabyData[1],
105 : pabyData[2],
106 : pabyData[3],
107 : pabyData[4],
108 : pabyData[5],
109 : pabyData[6],
110 : pabyData[7],
111 0 : pabyData[8]);
112 0 : return OGRERR_CORRUPT_DATA;
113 : }
114 :
115 : /* -------------------------------------------------------------------- */
116 : /* Get the geometry feature type. For now we assume that */
117 : /* geometry type is between 0 and 255 so we only have to fetch */
118 : /* one byte. */
119 : /* -------------------------------------------------------------------- */
120 748 : if( eByteOrder == wkbNDR )
121 668 : eGeometryType = (OGRwkbGeometryType) pabyData[1];
122 : else
123 80 : eGeometryType = (OGRwkbGeometryType) pabyData[4];
124 :
125 : /* -------------------------------------------------------------------- */
126 : /* Instantiate a geometry of the appropriate type, and */
127 : /* initialize from the input stream. */
128 : /* -------------------------------------------------------------------- */
129 748 : poGeom = createGeometry( eGeometryType );
130 :
131 748 : if( poGeom == NULL )
132 19 : return OGRERR_UNSUPPORTED_GEOMETRY_TYPE;
133 :
134 : /* -------------------------------------------------------------------- */
135 : /* Import from binary. */
136 : /* -------------------------------------------------------------------- */
137 729 : eErr = poGeom->importFromWkb( pabyData, nBytes );
138 :
139 : /* -------------------------------------------------------------------- */
140 : /* Assign spatial reference system. */
141 : /* -------------------------------------------------------------------- */
142 729 : if( eErr == OGRERR_NONE )
143 : {
144 729 : poGeom->assignSpatialReference( poSR );
145 729 : *ppoReturn = poGeom;
146 : }
147 : else
148 : {
149 0 : delete poGeom;
150 : }
151 :
152 729 : return eErr;
153 : }
154 :
155 : /************************************************************************/
156 : /* OGR_G_CreateFromWkb() */
157 : /************************************************************************/
158 : /**
159 : * \brief Create a geometry object of the appropriate type from it's well known binary representation.
160 : *
161 : * Note that if nBytes is passed as zero, no checking can be done on whether
162 : * the pabyData is sufficient. This can result in a crash if the input
163 : * data is corrupt. This function returns no indication of the number of
164 : * bytes from the data source actually used to represent the returned
165 : * geometry object. Use OGR_G_WkbSize() on the returned geometry to
166 : * establish the number of bytes it required in WKB format.
167 : *
168 : * The OGRGeometryFactory::createFromWkb() CPP method is the same as this
169 : * function.
170 : *
171 : * @param pabyData pointer to the input BLOB data.
172 : * @param hSRS handle to the spatial reference to be assigned to the
173 : * created geometry object. This may be NULL.
174 : * @param phGeometry the newly created geometry object will
175 : * be assigned to the indicated handle on return. This will be NULL in case
176 : * of failure.
177 : * @param nBytes the number of bytes of data available in pabyData, or -1
178 : * if it is not known, but assumed to be sufficient.
179 : *
180 : * @return OGRERR_NONE if all goes well, otherwise any of
181 : * OGRERR_NOT_ENOUGH_DATA, OGRERR_UNSUPPORTED_GEOMETRY_TYPE, or
182 : * OGRERR_CORRUPT_DATA may be returned.
183 : */
184 :
185 : OGRErr CPL_DLL OGR_G_CreateFromWkb( unsigned char *pabyData,
186 : OGRSpatialReferenceH hSRS,
187 : OGRGeometryH *phGeometry,
188 62 : int nBytes )
189 :
190 : {
191 : return OGRGeometryFactory::createFromWkb( pabyData,
192 : (OGRSpatialReference *) hSRS,
193 : (OGRGeometry **) phGeometry,
194 62 : nBytes );
195 : }
196 :
197 : /************************************************************************/
198 : /* createFromWkt() */
199 : /************************************************************************/
200 :
201 : /**
202 : * \brief Create a geometry object of the appropriate type from it's well known text representation.
203 : *
204 : * The C function OGR_G_CreateFromWkt() is the same as this method.
205 : *
206 : * @param ppszData input zero terminated string containing well known text
207 : * representation of the geometry to be created. The pointer
208 : * is updated to point just beyond that last character consumed.
209 : * @param poSR pointer to the spatial reference to be assigned to the
210 : * created geometry object. This may be NULL.
211 : * @param ppoReturn the newly created geometry object will be assigned to the
212 : * indicated pointer on return. This will be NULL if the
213 : * method fails.
214 : *
215 : * <b>Example:</b>
216 : *
217 : * <pre>
218 : * const char* wkt= "POINT(0 0)";
219 : *
220 : * // cast because OGR_G_CreateFromWkt will move the pointer
221 : * char* pszWkt = (char*) wkt.c_str();
222 : * OGRSpatialReferenceH ref = OSRNewSpatialReference(NULL);
223 : * OGRGeometryH new_geom;
224 : * OGRErr err = OGR_G_CreateFromWkt(&pszWkt, ref, &new_geom);
225 : * </pre>
226 : *
227 : *
228 : *
229 : * @return OGRERR_NONE if all goes well, otherwise any of
230 : * OGRERR_NOT_ENOUGH_DATA, OGRERR_UNSUPPORTED_GEOMETRY_TYPE, or
231 : * OGRERR_CORRUPT_DATA may be returned.
232 : */
233 :
234 : OGRErr OGRGeometryFactory::createFromWkt(char **ppszData,
235 : OGRSpatialReference * poSR,
236 16913 : OGRGeometry **ppoReturn )
237 :
238 : {
239 : OGRErr eErr;
240 : char szToken[OGR_WKT_TOKEN_MAX];
241 16913 : char *pszInput = *ppszData;
242 : OGRGeometry *poGeom;
243 :
244 16913 : *ppoReturn = NULL;
245 :
246 : /* -------------------------------------------------------------------- */
247 : /* Get the first token, which should be the geometry type. */
248 : /* -------------------------------------------------------------------- */
249 16913 : if( OGRWktReadToken( pszInput, szToken ) == NULL )
250 0 : return OGRERR_CORRUPT_DATA;
251 :
252 : /* -------------------------------------------------------------------- */
253 : /* Instantiate a geometry of the appropriate type. */
254 : /* -------------------------------------------------------------------- */
255 16913 : if( EQUAL(szToken,"POINT") )
256 : {
257 15058 : poGeom = new OGRPoint();
258 : }
259 :
260 1855 : else if( EQUAL(szToken,"LINESTRING") )
261 : {
262 182 : poGeom = new OGRLineString();
263 : }
264 :
265 1673 : else if( EQUAL(szToken,"POLYGON") )
266 : {
267 272 : poGeom = new OGRPolygon();
268 : }
269 :
270 1401 : else if( EQUAL(szToken,"GEOMETRYCOLLECTION") )
271 : {
272 68 : poGeom = new OGRGeometryCollection();
273 : }
274 :
275 1333 : else if( EQUAL(szToken,"MULTIPOLYGON") )
276 : {
277 126 : poGeom = new OGRMultiPolygon();
278 : }
279 :
280 1207 : else if( EQUAL(szToken,"MULTIPOINT") )
281 : {
282 85 : poGeom = new OGRMultiPoint();
283 : }
284 :
285 1122 : else if( EQUAL(szToken,"MULTILINESTRING") )
286 : {
287 95 : poGeom = new OGRMultiLineString();
288 : }
289 :
290 : else
291 : {
292 1027 : return OGRERR_UNSUPPORTED_GEOMETRY_TYPE;
293 : }
294 :
295 : /* -------------------------------------------------------------------- */
296 : /* Do the import. */
297 : /* -------------------------------------------------------------------- */
298 15886 : eErr = poGeom->importFromWkt( &pszInput );
299 :
300 : /* -------------------------------------------------------------------- */
301 : /* Assign spatial reference system. */
302 : /* -------------------------------------------------------------------- */
303 15886 : if( eErr == OGRERR_NONE )
304 : {
305 15712 : poGeom->assignSpatialReference( poSR );
306 15712 : *ppoReturn = poGeom;
307 15712 : *ppszData = pszInput;
308 : }
309 : else
310 : {
311 174 : delete poGeom;
312 : }
313 :
314 15886 : return eErr;
315 : }
316 :
317 : /************************************************************************/
318 : /* OGR_G_CreateFromWkt() */
319 : /************************************************************************/
320 : /**
321 : * \brief Create a geometry object of the appropriate type from it's well known text representation.
322 : *
323 : * The OGRGeometryFactory::createFromWkt CPP method is the same as this
324 : * function.
325 : *
326 : * @param ppszData input zero terminated string containing well known text
327 : * representation of the geometry to be created. The pointer
328 : * is updated to point just beyond that last character consumed.
329 : * @param hSRS handle to the spatial reference to be assigned to the
330 : * created geometry object. This may be NULL.
331 : * @param phGeometry the newly created geometry object will be assigned to the
332 : * indicated handle on return. This will be NULL if the
333 : * method fails.
334 : *
335 : * @return OGRERR_NONE if all goes well, otherwise any of
336 : * OGRERR_NOT_ENOUGH_DATA, OGRERR_UNSUPPORTED_GEOMETRY_TYPE, or
337 : * OGRERR_CORRUPT_DATA may be returned.
338 : */
339 :
340 : OGRErr CPL_DLL OGR_G_CreateFromWkt( char **ppszData,
341 : OGRSpatialReferenceH hSRS,
342 15628 : OGRGeometryH *phGeometry )
343 :
344 : {
345 : return OGRGeometryFactory::createFromWkt( ppszData,
346 : (OGRSpatialReference *) hSRS,
347 15628 : (OGRGeometry **) phGeometry );
348 : }
349 :
350 : /************************************************************************/
351 : /* createGeometry() */
352 : /************************************************************************/
353 :
354 : /**
355 : * \brief Create an empty geometry of desired type.
356 : *
357 : * This is equivelent to allocating the desired geometry with new, but
358 : * the allocation is guaranteed to take place in the context of the
359 : * GDAL/OGR heap.
360 : *
361 : * This method is the same as the C function OGR_G_CreateGeometry().
362 : *
363 : * @param eGeometryType the type code of the geometry class to be instantiated.
364 : *
365 : * @return the newly create geometry or NULL on failure.
366 : */
367 :
368 : OGRGeometry *
369 2081 : OGRGeometryFactory::createGeometry( OGRwkbGeometryType eGeometryType )
370 :
371 : {
372 2081 : switch( wkbFlatten(eGeometryType) )
373 : {
374 : case wkbPoint:
375 166 : return new OGRPoint();
376 :
377 : case wkbLineString:
378 904 : return new OGRLineString();
379 :
380 : case wkbPolygon:
381 554 : return new OGRPolygon();
382 :
383 : case wkbGeometryCollection:
384 70 : return new OGRGeometryCollection();
385 :
386 : case wkbMultiPolygon:
387 46 : return new OGRMultiPolygon();
388 :
389 : case wkbMultiPoint:
390 23 : return new OGRMultiPoint();
391 :
392 : case wkbMultiLineString:
393 30 : return new OGRMultiLineString();
394 :
395 : case wkbLinearRing:
396 269 : return new OGRLinearRing();
397 :
398 : default:
399 19 : return NULL;
400 : }
401 : }
402 :
403 : /************************************************************************/
404 : /* OGR_G_CreateGeometry() */
405 : /************************************************************************/
406 : /**
407 : * \brief Create an empty geometry of desired type.
408 : *
409 : * This is equivelent to allocating the desired geometry with new, but
410 : * the allocation is guaranteed to take place in the context of the
411 : * GDAL/OGR heap.
412 : *
413 : * This function is the same as the CPP method
414 : * OGRGeometryFactory::createGeometry.
415 : *
416 : * @param eGeometryType the type code of the geometry to be created.
417 : *
418 : * @return handle to the newly create geometry or NULL on failure.
419 : */
420 :
421 1327 : OGRGeometryH OGR_G_CreateGeometry( OGRwkbGeometryType eGeometryType )
422 :
423 : {
424 1327 : return (OGRGeometryH) OGRGeometryFactory::createGeometry( eGeometryType );
425 : }
426 :
427 :
428 : /************************************************************************/
429 : /* destroyGeometry() */
430 : /************************************************************************/
431 :
432 : /**
433 : * \brief Destroy geometry object.
434 : *
435 : * Equivalent to invoking delete on a geometry, but it guaranteed to take
436 : * place within the context of the GDAL/OGR heap.
437 : *
438 : * This method is the same as the C function OGR_G_DestroyGeometry().
439 : *
440 : * @param poGeom the geometry to deallocate.
441 : */
442 :
443 17801 : void OGRGeometryFactory::destroyGeometry( OGRGeometry *poGeom )
444 :
445 : {
446 17801 : delete poGeom;
447 17801 : }
448 :
449 :
450 : /************************************************************************/
451 : /* OGR_G_DestroyGeometry() */
452 : /************************************************************************/
453 : /**
454 : * \brief Destroy geometry object.
455 : *
456 : * Equivalent to invoking delete on a geometry, but it guaranteed to take
457 : * place within the context of the GDAL/OGR heap.
458 : *
459 : * This function is the same as the CPP method
460 : * OGRGeometryFactory::destroyGeometry.
461 : *
462 : * @param hGeom handle to the geometry to delete.
463 : */
464 :
465 17672 : void OGR_G_DestroyGeometry( OGRGeometryH hGeom )
466 :
467 : {
468 17672 : OGRGeometryFactory::destroyGeometry( (OGRGeometry *) hGeom );
469 17672 : }
470 :
471 : /************************************************************************/
472 : /* forceToPolygon() */
473 : /************************************************************************/
474 :
475 : /**
476 : * \brief Convert to polygon.
477 : *
478 : * Tries to force the provided geometry to be a polygon. Currently
479 : * this just effects a change on multipolygons. The passed in geometry is
480 : * consumed and a new one returned (or potentially the same one).
481 : *
482 : * @param poGeom the input geometry - ownership is passed to the method.
483 : * @return new geometry.
484 : */
485 :
486 45 : OGRGeometry *OGRGeometryFactory::forceToPolygon( OGRGeometry *poGeom )
487 :
488 : {
489 45 : if( poGeom == NULL )
490 0 : return NULL;
491 :
492 45 : OGRwkbGeometryType eGeomType = wkbFlatten(poGeom->getGeometryType());
493 :
494 45 : if( eGeomType != wkbGeometryCollection
495 : && eGeomType != wkbMultiPolygon )
496 31 : return poGeom;
497 :
498 : // build an aggregated polygon from all the polygon rings in the container.
499 14 : OGRPolygon *poPolygon = new OGRPolygon();
500 14 : OGRGeometryCollection *poGC = (OGRGeometryCollection *) poGeom;
501 : int iGeom;
502 :
503 30 : for( iGeom = 0; iGeom < poGC->getNumGeometries(); iGeom++ )
504 : {
505 16 : if( wkbFlatten(poGC->getGeometryRef(iGeom)->getGeometryType())
506 : != wkbPolygon )
507 9 : continue;
508 :
509 7 : OGRPolygon *poOldPoly = (OGRPolygon *) poGC->getGeometryRef(iGeom);
510 : int iRing;
511 :
512 7 : if( poOldPoly->getExteriorRing() == NULL )
513 2 : continue;
514 :
515 5 : poPolygon->addRing( poOldPoly->getExteriorRing() );
516 :
517 9 : for( iRing = 0; iRing < poOldPoly->getNumInteriorRings(); iRing++ )
518 4 : poPolygon->addRing( poOldPoly->getInteriorRing( iRing ) );
519 : }
520 :
521 14 : delete poGC;
522 :
523 14 : return poPolygon;
524 : }
525 :
526 : /************************************************************************/
527 : /* OGR_G_ForceToPolygon() */
528 : /************************************************************************/
529 :
530 : /**
531 : * \brief Convert to polygon.
532 : *
533 : * This function is the same as the C++ method
534 : * OGRGeometryFactory::forceToPolygon().
535 : *
536 : * @param hGeom handle to the geometry to convert (ownership surrendered).
537 : * @return the converted geometry (ownership to caller).
538 : *
539 : * @since GDAL/OGR 1.8.0
540 : */
541 :
542 35 : OGRGeometryH OGR_G_ForceToPolygon( OGRGeometryH hGeom )
543 :
544 : {
545 : return (OGRGeometryH)
546 35 : OGRGeometryFactory::forceToPolygon( (OGRGeometry *) hGeom );
547 : }
548 :
549 : /************************************************************************/
550 : /* forceToMultiPolygon() */
551 : /************************************************************************/
552 :
553 : /**
554 : * \brief Convert to multipolygon.
555 : *
556 : * Tries to force the provided geometry to be a multipolygon. Currently
557 : * this just effects a change on polygons. The passed in geometry is
558 : * consumed and a new one returned (or potentially the same one).
559 : *
560 : * @return new geometry.
561 : */
562 :
563 34 : OGRGeometry *OGRGeometryFactory::forceToMultiPolygon( OGRGeometry *poGeom )
564 :
565 : {
566 34 : if( poGeom == NULL )
567 0 : return NULL;
568 :
569 34 : OGRwkbGeometryType eGeomType = wkbFlatten(poGeom->getGeometryType());
570 :
571 : /* -------------------------------------------------------------------- */
572 : /* If this is already a MultiPolygon, nothing to do */
573 : /* -------------------------------------------------------------------- */
574 34 : if( eGeomType == wkbMultiPolygon )
575 : {
576 2 : return poGeom;
577 : }
578 :
579 : /* -------------------------------------------------------------------- */
580 : /* Check for the case of a geometrycollection that can be */
581 : /* promoted to MultiPolygon. */
582 : /* -------------------------------------------------------------------- */
583 32 : if( eGeomType == wkbGeometryCollection )
584 : {
585 : int iGeom;
586 12 : int bAllPoly = TRUE;
587 12 : OGRGeometryCollection *poGC = (OGRGeometryCollection *) poGeom;
588 :
589 25 : for( iGeom = 0; iGeom < poGC->getNumGeometries(); iGeom++ )
590 : {
591 13 : if( wkbFlatten(poGC->getGeometryRef(iGeom)->getGeometryType())
592 : != wkbPolygon )
593 9 : bAllPoly = FALSE;
594 : }
595 :
596 12 : if( !bAllPoly )
597 8 : return poGeom;
598 :
599 4 : OGRMultiPolygon *poMP = new OGRMultiPolygon();
600 :
601 12 : while( poGC->getNumGeometries() > 0 )
602 : {
603 4 : poMP->addGeometryDirectly( poGC->getGeometryRef(0) );
604 4 : poGC->removeGeometry( 0, FALSE );
605 : }
606 :
607 4 : delete poGC;
608 :
609 4 : return poMP;
610 : }
611 :
612 : /* -------------------------------------------------------------------- */
613 : /* Eventually we should try to split the polygon into component */
614 : /* island polygons. But thats alot of work and can be put off. */
615 : /* -------------------------------------------------------------------- */
616 20 : if( eGeomType != wkbPolygon )
617 8 : return poGeom;
618 :
619 12 : OGRMultiPolygon *poMP = new OGRMultiPolygon();
620 12 : poMP->addGeometryDirectly( poGeom );
621 :
622 12 : return poMP;
623 : }
624 :
625 : /************************************************************************/
626 : /* OGR_G_ForceToMultiPolygon() */
627 : /************************************************************************/
628 :
629 : /**
630 : * \brief Convert to multipolygon.
631 : *
632 : * This function is the same as the C++ method
633 : * OGRGeometryFactory::forceToMultiPolygon().
634 : *
635 : * @param hGeom handle to the geometry to convert (ownership surrendered).
636 : * @return the converted geometry (ownership to caller).
637 : *
638 : * @since GDAL/OGR 1.8.0
639 : */
640 :
641 30 : OGRGeometryH OGR_G_ForceToMultiPolygon( OGRGeometryH hGeom )
642 :
643 : {
644 : return (OGRGeometryH)
645 30 : OGRGeometryFactory::forceToMultiPolygon( (OGRGeometry *) hGeom );
646 : }
647 :
648 : /************************************************************************/
649 : /* forceToMultiPoint() */
650 : /************************************************************************/
651 :
652 : /**
653 : * \brief Convert to multipoint.
654 : *
655 : * Tries to force the provided geometry to be a multipoint. Currently
656 : * this just effects a change on points. The passed in geometry is
657 : * consumed and a new one returned (or potentially the same one).
658 : *
659 : * @return new geometry.
660 : */
661 :
662 26 : OGRGeometry *OGRGeometryFactory::forceToMultiPoint( OGRGeometry *poGeom )
663 :
664 : {
665 26 : if( poGeom == NULL )
666 0 : return NULL;
667 :
668 26 : OGRwkbGeometryType eGeomType = wkbFlatten(poGeom->getGeometryType());
669 :
670 : /* -------------------------------------------------------------------- */
671 : /* If this is already a MultiPoint, nothing to do */
672 : /* -------------------------------------------------------------------- */
673 26 : if( eGeomType == wkbMultiPoint )
674 : {
675 2 : return poGeom;
676 : }
677 :
678 : /* -------------------------------------------------------------------- */
679 : /* Check for the case of a geometrycollection that can be */
680 : /* promoted to MultiPoint. */
681 : /* -------------------------------------------------------------------- */
682 24 : if( eGeomType == wkbGeometryCollection )
683 : {
684 : int iGeom;
685 12 : int bAllPoint = TRUE;
686 12 : OGRGeometryCollection *poGC = (OGRGeometryCollection *) poGeom;
687 :
688 26 : for( iGeom = 0; iGeom < poGC->getNumGeometries(); iGeom++ )
689 : {
690 14 : if( wkbFlatten(poGC->getGeometryRef(iGeom)->getGeometryType())
691 : != wkbPoint )
692 10 : bAllPoint = FALSE;
693 : }
694 :
695 12 : if( !bAllPoint )
696 8 : return poGeom;
697 :
698 4 : OGRMultiPoint *poMP = new OGRMultiPoint();
699 :
700 12 : while( poGC->getNumGeometries() > 0 )
701 : {
702 4 : poMP->addGeometryDirectly( poGC->getGeometryRef(0) );
703 4 : poGC->removeGeometry( 0, FALSE );
704 : }
705 :
706 4 : delete poGC;
707 :
708 4 : return poMP;
709 : }
710 :
711 12 : if( eGeomType != wkbPoint )
712 9 : return poGeom;
713 :
714 3 : OGRMultiPoint *poMP = new OGRMultiPoint();
715 3 : poMP->addGeometryDirectly( poGeom );
716 :
717 3 : return poMP;
718 : }
719 :
720 : /************************************************************************/
721 : /* OGR_G_ForceToMultiPoint() */
722 : /************************************************************************/
723 :
724 : /**
725 : * \brief Convert to multipoint.
726 : *
727 : * This function is the same as the C++ method
728 : * OGRGeometryFactory::forceToMultiPoint().
729 : *
730 : * @param hGeom handle to the geometry to convert (ownership surrendered).
731 : * @return the converted geometry (ownership to caller).
732 : *
733 : * @since GDAL/OGR 1.8.0
734 : */
735 :
736 26 : OGRGeometryH OGR_G_ForceToMultiPoint( OGRGeometryH hGeom )
737 :
738 : {
739 : return (OGRGeometryH)
740 26 : OGRGeometryFactory::forceToMultiPoint( (OGRGeometry *) hGeom );
741 : }
742 :
743 : /************************************************************************/
744 : /* forceToMultiLinestring() */
745 : /************************************************************************/
746 :
747 : /**
748 : * \brief Convert to multilinestring.
749 : *
750 : * Tries to force the provided geometry to be a multilinestring.
751 : *
752 : * - linestrings are placed in a multilinestring.
753 : * - geometry collections will be converted to multilinestring if they only
754 : * contain linestrings.
755 : * - polygons will be changed to a collection of linestrings (one per ring).
756 : *
757 : * The passed in geometry is
758 : * consumed and a new one returned (or potentially the same one).
759 : *
760 : * @return new geometry.
761 : */
762 :
763 47 : OGRGeometry *OGRGeometryFactory::forceToMultiLineString( OGRGeometry *poGeom )
764 :
765 : {
766 47 : if( poGeom == NULL )
767 0 : return NULL;
768 :
769 47 : OGRwkbGeometryType eGeomType = wkbFlatten(poGeom->getGeometryType());
770 :
771 : /* -------------------------------------------------------------------- */
772 : /* If this is already a MultiLineString, nothing to do */
773 : /* -------------------------------------------------------------------- */
774 47 : if( eGeomType == wkbMultiLineString )
775 : {
776 8 : return poGeom;
777 : }
778 :
779 : /* -------------------------------------------------------------------- */
780 : /* Check for the case of a geometrycollection that can be */
781 : /* promoted to MultiLineString. */
782 : /* -------------------------------------------------------------------- */
783 39 : if( eGeomType == wkbGeometryCollection )
784 : {
785 : int iGeom;
786 12 : int bAllLines = TRUE;
787 12 : OGRGeometryCollection *poGC = (OGRGeometryCollection *) poGeom;
788 :
789 26 : for( iGeom = 0; iGeom < poGC->getNumGeometries(); iGeom++ )
790 : {
791 14 : if( wkbFlatten(poGC->getGeometryRef(iGeom)->getGeometryType())
792 : != wkbLineString )
793 9 : bAllLines = FALSE;
794 : }
795 :
796 12 : if( !bAllLines )
797 8 : return poGeom;
798 :
799 4 : OGRMultiLineString *poMP = new OGRMultiLineString();
800 :
801 13 : while( poGC->getNumGeometries() > 0 )
802 : {
803 5 : poMP->addGeometryDirectly( poGC->getGeometryRef(0) );
804 5 : poGC->removeGeometry( 0, FALSE );
805 : }
806 :
807 4 : delete poGC;
808 :
809 4 : return poMP;
810 : }
811 :
812 : /* -------------------------------------------------------------------- */
813 : /* Turn a linestring into a multilinestring. */
814 : /* -------------------------------------------------------------------- */
815 27 : if( eGeomType == wkbLineString )
816 : {
817 3 : OGRMultiLineString *poMP = new OGRMultiLineString();
818 3 : poMP->addGeometryDirectly( poGeom );
819 3 : return poMP;
820 : }
821 :
822 : /* -------------------------------------------------------------------- */
823 : /* Convert polygons into a multilinestring. */
824 : /* -------------------------------------------------------------------- */
825 24 : if( eGeomType == wkbPolygon )
826 : {
827 5 : OGRMultiLineString *poMP = new OGRMultiLineString();
828 5 : OGRPolygon *poPoly = (OGRPolygon *) poGeom;
829 : int iRing;
830 :
831 13 : for( iRing = 0; iRing < poPoly->getNumInteriorRings()+1; iRing++ )
832 : {
833 : OGRLineString *poNewLS, *poLR;
834 :
835 9 : if( iRing == 0 )
836 : {
837 5 : poLR = poPoly->getExteriorRing();
838 5 : if( poLR == NULL )
839 1 : break;
840 : }
841 : else
842 4 : poLR = poPoly->getInteriorRing(iRing-1);
843 :
844 8 : if (poLR == NULL || poLR->getNumPoints() == 0)
845 2 : continue;
846 :
847 6 : poNewLS = new OGRLineString();
848 6 : poNewLS->addSubLineString( poLR );
849 6 : poMP->addGeometryDirectly( poNewLS );
850 : }
851 :
852 5 : delete poPoly;
853 :
854 5 : return poMP;
855 : }
856 :
857 : /* -------------------------------------------------------------------- */
858 : /* Convert multi-polygons into a multilinestring. */
859 : /* -------------------------------------------------------------------- */
860 19 : if( eGeomType == wkbMultiPolygon )
861 : {
862 3 : OGRMultiLineString *poMP = new OGRMultiLineString();
863 3 : OGRMultiPolygon *poMPoly = (OGRMultiPolygon *) poGeom;
864 : int iPoly;
865 :
866 8 : for( iPoly = 0; iPoly < poMPoly->getNumGeometries(); iPoly++ )
867 : {
868 5 : OGRPolygon *poPoly = (OGRPolygon*) poMPoly->getGeometryRef(iPoly);
869 : int iRing;
870 :
871 12 : for( iRing = 0; iRing < poPoly->getNumInteriorRings()+1; iRing++ )
872 : {
873 : OGRLineString *poNewLS, *poLR;
874 :
875 8 : if( iRing == 0 )
876 : {
877 5 : poLR = poPoly->getExteriorRing();
878 5 : if( poLR == NULL )
879 1 : break;
880 : }
881 : else
882 3 : poLR = poPoly->getInteriorRing(iRing-1);
883 :
884 7 : if (poLR == NULL || poLR->getNumPoints() == 0)
885 1 : continue;
886 :
887 6 : poNewLS = new OGRLineString();
888 6 : poNewLS->addSubLineString( poLR );
889 6 : poMP->addGeometryDirectly( poNewLS );
890 : }
891 : }
892 3 : delete poMPoly;
893 :
894 3 : return poMP;
895 : }
896 :
897 16 : return poGeom;
898 : }
899 :
900 : /************************************************************************/
901 : /* OGR_G_ForceToMultiLineString() */
902 : /************************************************************************/
903 :
904 : /**
905 : * \brief Convert to multilinestring.
906 : *
907 : * This function is the same as the C++ method
908 : * OGRGeometryFactory::forceToMultiLineString().
909 : *
910 : * @param hGeom handle to the geometry to convert (ownership surrendered).
911 : * @return the converted geometry (ownership to caller).
912 : *
913 : * @since GDAL/OGR 1.8.0
914 : */
915 :
916 28 : OGRGeometryH OGR_G_ForceToMultiLineString( OGRGeometryH hGeom )
917 :
918 : {
919 : return (OGRGeometryH)
920 28 : OGRGeometryFactory::forceToMultiLineString( (OGRGeometry *) hGeom );
921 : }
922 :
923 : /************************************************************************/
924 : /* organizePolygons() */
925 : /************************************************************************/
926 :
927 : typedef struct _sPolyExtended sPolyExtended;
928 :
929 : struct _sPolyExtended
930 594 : {
931 : OGRPolygon* poPolygon;
932 : OGREnvelope sEnvelope;
933 : OGRLinearRing* poExteriorRing;
934 : OGRPoint poAPoint;
935 : int nInitialIndex;
936 : int bIsTopLevel;
937 : OGRPolygon* poEnclosingPolygon;
938 : double dfArea;
939 : int bIsCW;
940 : };
941 :
942 341 : static int OGRGeometryFactoryCompareArea(const void* p1, const void* p2)
943 : {
944 341 : const sPolyExtended* psPoly1 = (const sPolyExtended*) p1;
945 341 : const sPolyExtended* psPoly2 = (const sPolyExtended*) p2;
946 341 : if (psPoly2->dfArea < psPoly1->dfArea)
947 189 : return -1;
948 152 : else if (psPoly2->dfArea > psPoly1->dfArea)
949 136 : return 1;
950 : else
951 16 : return 0;
952 : }
953 :
954 352 : static int OGRGeometryFactoryCompareByIndex(const void* p1, const void* p2)
955 : {
956 352 : const sPolyExtended* psPoly1 = (const sPolyExtended*) p1;
957 352 : const sPolyExtended* psPoly2 = (const sPolyExtended*) p2;
958 352 : if (psPoly1->nInitialIndex < psPoly2->nInitialIndex)
959 210 : return -1;
960 142 : else if (psPoly1->nInitialIndex > psPoly2->nInitialIndex)
961 142 : return 1;
962 : else
963 0 : return 0;
964 : }
965 :
966 : #define N_CRITICAL_PART_NUMBER 100
967 :
968 : typedef enum
969 : {
970 : METHOD_NORMAL,
971 : METHOD_SKIP,
972 : METHOD_ONLY_CCW
973 : } OrganizePolygonMethod;
974 :
975 : /**
976 : * \brief Organize polygons based on geometries.
977 : *
978 : * Analyse a set of rings (passed as simple polygons), and based on a
979 : * geometric analysis convert them into a polygon with inner rings,
980 : * or a MultiPolygon if dealing with more than one polygon.
981 : *
982 : * All the input geometries must be OGRPolygons with only a valid exterior
983 : * ring (at least 4 points) and no interior rings.
984 : *
985 : * The passed in geometries become the responsibility of the method, but the
986 : * papoPolygons "pointer array" remains owned by the caller.
987 : *
988 : * For faster computation, a polygon is considered to be inside
989 : * another one if a single point of its external ring is included into the other one.
990 : * (unless 'OGR_DEBUG_ORGANIZE_POLYGONS' configuration option is set to TRUE.
991 : * In that case, a slower algorithm that tests exact topological relationships
992 : * is used if GEOS is available.)
993 : *
994 : * In cases where a big number of polygons is passed to this function, the default processing
995 : * may be really slow. You can skip the processing by adding METHOD=SKIP
996 : * to the option list (the result of the function will be a multi-polygon with all polygons
997 : * as toplevel polygons) or only make it analyze counterclockwise polygons by adding
998 : * METHOD=ONLY_CCW to the option list if you can assume that the outline
999 : * of holes is counterclockwise defined (this is the convention for shapefiles e.g.)
1000 : *
1001 : * If the OGR_ORGANIZE_POLYGONS configuration option is defined, its value will override
1002 : * the value of the METHOD option of papszOptions (usefull to modify the behaviour of the
1003 : * shapefile driver)
1004 : *
1005 : * @param papoPolygons array of geometry pointers - should all be OGRPolygons.
1006 : * Ownership of the geometries is passed, but not of the array itself.
1007 : * @param nPolygonCount number of items in papoPolygons
1008 : * @param pbIsValidGeometry value will be set TRUE if result is valid or
1009 : * FALSE otherwise.
1010 : * @param papszOptions a list of strings for passing options
1011 : *
1012 : * @return a single resulting geometry (either OGRPolygon or OGRMultiPolygon).
1013 : */
1014 :
1015 : OGRGeometry* OGRGeometryFactory::organizePolygons( OGRGeometry **papoPolygons,
1016 : int nPolygonCount,
1017 : int *pbIsValidGeometry,
1018 77 : const char** papszOptions )
1019 : {
1020 : int bUseFastVersion;
1021 : int i, j;
1022 77 : OGRGeometry* geom = NULL;
1023 77 : OrganizePolygonMethod method = METHOD_NORMAL;
1024 :
1025 : /* -------------------------------------------------------------------- */
1026 : /* Trivial case of a single polygon. */
1027 : /* -------------------------------------------------------------------- */
1028 77 : if (nPolygonCount == 1)
1029 : {
1030 0 : geom = papoPolygons[0];
1031 0 : papoPolygons[0] = NULL;
1032 :
1033 0 : if( pbIsValidGeometry )
1034 0 : *pbIsValidGeometry = TRUE;
1035 :
1036 0 : return geom;
1037 : }
1038 :
1039 77 : if (CSLTestBoolean(CPLGetConfigOption("OGR_DEBUG_ORGANIZE_POLYGONS", "NO")))
1040 : {
1041 : /* -------------------------------------------------------------------- */
1042 : /* A wee bit of a warning. */
1043 : /* -------------------------------------------------------------------- */
1044 : static int firstTime = 1;
1045 0 : if (!haveGEOS() && firstTime)
1046 : {
1047 : CPLDebug("OGR",
1048 : "In OGR_DEBUG_ORGANIZE_POLYGONS mode, GDAL should be built with GEOS support enabled in order "
1049 0 : "OGRGeometryFactory::organizePolygons to provide reliable results on complex polygons.");
1050 0 : firstTime = 0;
1051 : }
1052 0 : bUseFastVersion = !haveGEOS();
1053 : }
1054 : else
1055 : {
1056 77 : bUseFastVersion = TRUE;
1057 : }
1058 :
1059 : /* -------------------------------------------------------------------- */
1060 : /* Setup per polygon envelope and area information. */
1061 : /* -------------------------------------------------------------------- */
1062 77 : sPolyExtended* asPolyEx = new sPolyExtended[nPolygonCount];
1063 :
1064 77 : int go_on = TRUE;
1065 77 : int bMixedUpGeometries = FALSE;
1066 77 : int bNonPolygon = FALSE;
1067 77 : int bFoundCCW = FALSE;
1068 :
1069 77 : const char* pszMethodValue = CSLFetchNameValue( (char**)papszOptions, "METHOD" );
1070 77 : const char* pszMethodValueOption = CPLGetConfigOption("OGR_ORGANIZE_POLYGONS", NULL);
1071 140 : if (pszMethodValueOption != NULL && pszMethodValueOption[0] != '\0')
1072 1 : pszMethodValue = pszMethodValueOption;
1073 :
1074 77 : if (pszMethodValue != NULL)
1075 : {
1076 24 : if (EQUAL(pszMethodValue, "SKIP"))
1077 : {
1078 0 : method = METHOD_SKIP;
1079 0 : bMixedUpGeometries = TRUE;
1080 : }
1081 24 : else if (EQUAL(pszMethodValue, "ONLY_CCW"))
1082 : {
1083 23 : method = METHOD_ONLY_CCW;
1084 : }
1085 1 : else if (!EQUAL(pszMethodValue, "DEFAULT"))
1086 : {
1087 : CPLError(CE_Warning, CPLE_AppDefined,
1088 0 : "Unrecognized value for METHOD option : %s", pszMethodValue);
1089 : }
1090 : }
1091 :
1092 77 : int nCountCWPolygon = 0;
1093 77 : int indexOfCWPolygon = -1;
1094 :
1095 374 : for(i=0;i<nPolygonCount;i++)
1096 : {
1097 297 : asPolyEx[i].nInitialIndex = i;
1098 297 : asPolyEx[i].poPolygon = (OGRPolygon*)papoPolygons[i];
1099 297 : papoPolygons[i]->getEnvelope(&asPolyEx[i].sEnvelope);
1100 :
1101 297 : if( wkbFlatten(papoPolygons[i]->getGeometryType()) == wkbPolygon
1102 : && ((OGRPolygon *) papoPolygons[i])->getNumInteriorRings() == 0
1103 : && ((OGRPolygon *)papoPolygons[i])->getExteriorRing()->getNumPoints() >= 4)
1104 : {
1105 297 : asPolyEx[i].dfArea = asPolyEx[i].poPolygon->get_Area();
1106 297 : asPolyEx[i].poExteriorRing = asPolyEx[i].poPolygon->getExteriorRing();
1107 297 : asPolyEx[i].poExteriorRing->getPoint(0, &asPolyEx[i].poAPoint);
1108 297 : asPolyEx[i].bIsCW = asPolyEx[i].poExteriorRing->isClockwise();
1109 297 : if (asPolyEx[i].bIsCW)
1110 : {
1111 224 : indexOfCWPolygon = i;
1112 224 : nCountCWPolygon ++;
1113 : }
1114 297 : if (!bFoundCCW)
1115 177 : bFoundCCW = ! (asPolyEx[i].bIsCW);
1116 : }
1117 : else
1118 : {
1119 0 : if( !bMixedUpGeometries )
1120 : {
1121 : CPLError(
1122 : CE_Warning, CPLE_AppDefined,
1123 : "organizePolygons() received an unexpected geometry.\n"
1124 : "Either a polygon with interior rings, or a polygon with less than 4 points,\n"
1125 0 : "or a non-Polygon geometry. Return arguments as a collection." );
1126 0 : bMixedUpGeometries = TRUE;
1127 : }
1128 0 : if( wkbFlatten(papoPolygons[i]->getGeometryType()) != wkbPolygon )
1129 0 : bNonPolygon = TRUE;
1130 : }
1131 : }
1132 :
1133 : /* If we are in ONLY_CCW mode and that we have found that there is only one outer ring, */
1134 : /* then it is pretty easy : we can assume that all other rings are inside */
1135 77 : if (method == METHOD_ONLY_CCW && nCountCWPolygon == 1 && bUseFastVersion)
1136 : {
1137 6 : geom = asPolyEx[indexOfCWPolygon].poPolygon;
1138 19 : for(i=0; i<nPolygonCount; i++)
1139 : {
1140 13 : if (i != indexOfCWPolygon)
1141 : {
1142 7 : ((OGRPolygon*)geom)->addRing(asPolyEx[i].poPolygon->getExteriorRing());
1143 7 : delete asPolyEx[i].poPolygon;
1144 : }
1145 : }
1146 6 : delete [] asPolyEx;
1147 6 : if (pbIsValidGeometry)
1148 6 : *pbIsValidGeometry = TRUE;
1149 6 : return geom;
1150 : }
1151 :
1152 : /* Emits a warning if the number of parts is sufficiently big to anticipate for */
1153 : /* very long computation time, and the user didn't specify an explicit method */
1154 71 : if (nPolygonCount > N_CRITICAL_PART_NUMBER && method == METHOD_NORMAL && pszMethodValue == NULL)
1155 : {
1156 : static int firstTime = 1;
1157 0 : if (firstTime)
1158 : {
1159 0 : if (bFoundCCW)
1160 : {
1161 : CPLError( CE_Warning, CPLE_AppDefined,
1162 : "organizePolygons() received a polygon with more than %d parts. "
1163 : "The processing may be really slow.\n"
1164 : "You can skip the processing by setting METHOD=SKIP, "
1165 : "or only make it analyze counter-clock wise parts by setting "
1166 : "METHOD=ONLY_CCW if you can assume that the "
1167 0 : "outline of holes is counter-clock wise defined", N_CRITICAL_PART_NUMBER);
1168 : }
1169 : else
1170 : {
1171 : CPLError( CE_Warning, CPLE_AppDefined,
1172 : "organizePolygons() received a polygon with more than %d parts. "
1173 : "The processing may be really slow.\n"
1174 : "You can skip the processing by setting METHOD=SKIP.",
1175 0 : N_CRITICAL_PART_NUMBER);
1176 : }
1177 0 : firstTime = 0;
1178 : }
1179 : }
1180 :
1181 :
1182 : /* This a several steps algorithm :
1183 : 1) Sort polygons by descending areas
1184 : 2) For each polygon of rank i, find its smallest enclosing polygon
1185 : among the polygons of rank [i-1 ... 0]. If there are no such polygon,
1186 : this is a toplevel polygon. Otherwise, depending on if the enclosing
1187 : polygon is toplevel or not, we can decide if we are toplevel or not
1188 : 3) Re-sort the polygons to retrieve their inital order (nicer for some applications)
1189 : 4) For each non toplevel polygon (= inner ring), add it to its outer ring
1190 : 5) Add the toplevel polygons to the multipolygon
1191 :
1192 : Complexity : O(nPolygonCount^2)
1193 : */
1194 :
1195 : /* Compute how each polygon relate to the other ones
1196 : To save a bit of computation we always begin the computation by a test
1197 : on the enveloppe. We also take into account the areas to avoid some
1198 : useless tests. (A contains B implies envelop(A) contains envelop(B)
1199 : and area(A) > area(B)) In practise, we can hope that few full geometry
1200 : intersection of inclusion test is done:
1201 : * if the polygons are well separated geographically (a set of islands
1202 : for example), no full geometry intersection or inclusion test is done.
1203 : (the envelopes don't intersect each other)
1204 :
1205 : * if the polygons are 'lake inside an island inside a lake inside an
1206 : area' and that each polygon is much smaller than its enclosing one,
1207 : their bounding boxes are stricly contained into each oter, and thus,
1208 : no full geometry intersection or inclusion test is done
1209 : */
1210 :
1211 71 : if (!bMixedUpGeometries)
1212 : {
1213 : /* STEP 1 : Sort polygons by descending area */
1214 71 : qsort(asPolyEx, nPolygonCount, sizeof(sPolyExtended), OGRGeometryFactoryCompareArea);
1215 : }
1216 71 : papoPolygons = NULL; /* just to use to avoid it afterwards */
1217 :
1218 : /* -------------------------------------------------------------------- */
1219 : /* Compute relationships, if things seem well structured. */
1220 : /* -------------------------------------------------------------------- */
1221 :
1222 : /* The first (largest) polygon is necessarily top-level */
1223 71 : asPolyEx[0].bIsTopLevel = TRUE;
1224 71 : asPolyEx[0].poEnclosingPolygon = NULL;
1225 :
1226 71 : int nCountTopLevel = 1;
1227 :
1228 : /* STEP 2 */
1229 284 : for(i=1; !bMixedUpGeometries && go_on && i<nPolygonCount; i++)
1230 : {
1231 213 : if (method == METHOD_ONLY_CCW && asPolyEx[i].bIsCW)
1232 : {
1233 24 : nCountTopLevel ++;
1234 24 : asPolyEx[i].bIsTopLevel = TRUE;
1235 24 : asPolyEx[i].poEnclosingPolygon = NULL;
1236 24 : continue;
1237 : }
1238 :
1239 409 : for(j=i-1; go_on && j>=0;j--)
1240 : {
1241 373 : int b_i_inside_j = FALSE;
1242 :
1243 373 : if (method == METHOD_ONLY_CCW && asPolyEx[j].bIsCW == FALSE)
1244 : {
1245 : /* In that mode, i which is CCW if we reach here can only be */
1246 : /* included in a CW polygon */
1247 5 : continue;
1248 : }
1249 :
1250 368 : if (asPolyEx[j].sEnvelope.Contains(asPolyEx[i].sEnvelope))
1251 : {
1252 154 : if (bUseFastVersion)
1253 : {
1254 : /* Note that isPointInRing only test strict inclusion in the ring */
1255 154 : if (asPolyEx[j].poExteriorRing->isPointInRing(&asPolyEx[i].poAPoint, FALSE))
1256 : {
1257 152 : b_i_inside_j = TRUE;
1258 : }
1259 2 : else if (asPolyEx[j].poExteriorRing->isPointOnRingBoundary(&asPolyEx[i].poAPoint, FALSE))
1260 : {
1261 : /* If the point of i is on the boundary of j, we will iterate over the other points of i */
1262 2 : int k, nPoints = asPolyEx[i].poExteriorRing->getNumPoints();
1263 7 : for(k=1;k<nPoints;k++)
1264 : {
1265 6 : OGRPoint point;
1266 6 : asPolyEx[i].poExteriorRing->getPoint(k, &point);
1267 6 : if (asPolyEx[j].poExteriorRing->isPointInRing(&point, FALSE))
1268 : {
1269 : /* If then point is strictly included in j, then i is considered inside j */
1270 1 : b_i_inside_j = TRUE;
1271 1 : break;
1272 : }
1273 5 : else if (asPolyEx[j].poExteriorRing->isPointOnRingBoundary(&point, FALSE))
1274 : {
1275 : /* If it is on the boundary of j, iterate again */
1276 : }
1277 : else
1278 : {
1279 : /* If it is outside, then i cannot be inside j */
1280 : break;
1281 : }
1282 : }
1283 : }
1284 : }
1285 0 : else if (asPolyEx[j].poPolygon->Contains(asPolyEx[i].poPolygon))
1286 : {
1287 0 : b_i_inside_j = TRUE;
1288 : }
1289 : }
1290 :
1291 :
1292 368 : if (b_i_inside_j)
1293 : {
1294 153 : if (asPolyEx[j].bIsTopLevel)
1295 : {
1296 : /* We are a lake */
1297 119 : asPolyEx[i].bIsTopLevel = FALSE;
1298 119 : asPolyEx[i].poEnclosingPolygon = asPolyEx[j].poPolygon;
1299 : }
1300 : else
1301 : {
1302 : /* We are included in a something not toplevel (a lake), */
1303 : /* so in OGCSF we are considered as toplevel too */
1304 34 : nCountTopLevel ++;
1305 34 : asPolyEx[i].bIsTopLevel = TRUE;
1306 34 : asPolyEx[i].poEnclosingPolygon = NULL;
1307 : }
1308 153 : break;
1309 : }
1310 : /* We use Overlaps instead of Intersects to be more
1311 : tolerant about touching polygons */
1312 215 : else if ( bUseFastVersion || !asPolyEx[i].sEnvelope.Intersects(asPolyEx[j].sEnvelope)
1313 : || !asPolyEx[i].poPolygon->Overlaps(asPolyEx[j].poPolygon) )
1314 : {
1315 :
1316 : }
1317 : else
1318 : {
1319 : /* Bad... The polygons are intersecting but no one is
1320 : contained inside the other one. This is a really broken
1321 : case. We just make a multipolygon with the whole set of
1322 : polygons */
1323 0 : go_on = FALSE;
1324 : #ifdef DEBUG
1325 : char* wkt1;
1326 : char* wkt2;
1327 0 : asPolyEx[i].poPolygon->exportToWkt(&wkt1);
1328 0 : asPolyEx[j].poPolygon->exportToWkt(&wkt2);
1329 : CPLDebug( "OGR",
1330 : "Bad intersection for polygons %d and %d\n"
1331 : "geom %d: %s\n"
1332 : "geom %d: %s",
1333 0 : i, j, i, wkt1, j, wkt2 );
1334 0 : CPLFree(wkt1);
1335 0 : CPLFree(wkt2);
1336 : #endif
1337 : }
1338 : }
1339 :
1340 189 : if (j < 0)
1341 : {
1342 : /* We come here because we are not included in anything */
1343 : /* We are toplevel */
1344 36 : nCountTopLevel ++;
1345 36 : asPolyEx[i].bIsTopLevel = TRUE;
1346 36 : asPolyEx[i].poEnclosingPolygon = NULL;
1347 : }
1348 : }
1349 :
1350 71 : if (pbIsValidGeometry)
1351 71 : *pbIsValidGeometry = go_on && !bMixedUpGeometries;
1352 :
1353 : /* -------------------------------------------------------------------- */
1354 : /* Things broke down - just turn everything into a multipolygon. */
1355 : /* -------------------------------------------------------------------- */
1356 71 : if ( !go_on || bMixedUpGeometries )
1357 : {
1358 0 : if( bNonPolygon )
1359 0 : geom = new OGRGeometryCollection();
1360 : else
1361 0 : geom = new OGRMultiPolygon();
1362 :
1363 0 : for( i=0; i < nPolygonCount; i++ )
1364 : {
1365 : ((OGRGeometryCollection*)geom)->
1366 0 : addGeometryDirectly( asPolyEx[i].poPolygon );
1367 : }
1368 : }
1369 :
1370 : /* -------------------------------------------------------------------- */
1371 : /* Try to turn into one or more polygons based on the ring */
1372 : /* relationships. */
1373 : /* -------------------------------------------------------------------- */
1374 : else
1375 : {
1376 : /* STEP 3: Resort in initial order */
1377 71 : qsort(asPolyEx, nPolygonCount, sizeof(sPolyExtended), OGRGeometryFactoryCompareByIndex);
1378 :
1379 : /* STEP 4: Add holes as rings of their enclosing polygon */
1380 355 : for(i=0;i<nPolygonCount;i++)
1381 : {
1382 284 : if (asPolyEx[i].bIsTopLevel == FALSE)
1383 : {
1384 : asPolyEx[i].poEnclosingPolygon->addRing(
1385 119 : asPolyEx[i].poPolygon->getExteriorRing());
1386 119 : delete asPolyEx[i].poPolygon;
1387 : }
1388 165 : else if (nCountTopLevel == 1)
1389 : {
1390 19 : geom = asPolyEx[i].poPolygon;
1391 : }
1392 : }
1393 :
1394 : /* STEP 5: Add toplevel polygons */
1395 71 : if (nCountTopLevel > 1)
1396 : {
1397 296 : for(i=0;i<nPolygonCount;i++)
1398 : {
1399 244 : if (asPolyEx[i].bIsTopLevel)
1400 : {
1401 146 : if (geom == NULL)
1402 : {
1403 52 : geom = new OGRMultiPolygon();
1404 52 : ((OGRMultiPolygon*)geom)->addGeometryDirectly(asPolyEx[i].poPolygon);
1405 : }
1406 : else
1407 : {
1408 94 : ((OGRMultiPolygon*)geom)->addGeometryDirectly(asPolyEx[i].poPolygon);
1409 : }
1410 : }
1411 : }
1412 : }
1413 : }
1414 :
1415 71 : delete[] asPolyEx;
1416 :
1417 71 : return geom;
1418 : }
1419 :
1420 : /************************************************************************/
1421 : /* createFromGML() */
1422 : /************************************************************************/
1423 :
1424 : /**
1425 : * \brief Create geometry from GML.
1426 : *
1427 : * This method translates a fragment of GML containing only the geometry
1428 : * portion into a corresponding OGRGeometry. There are many limitations
1429 : * on the forms of GML geometries supported by this parser, but they are
1430 : * too numerous to list here.
1431 : *
1432 : * The C function OGR_G_CreateFromGML() is the same as this method.
1433 : *
1434 : * @param pszData The GML fragment for the geometry.
1435 : *
1436 : * @return a geometry on succes, or NULL on error.
1437 : */
1438 :
1439 57 : OGRGeometry *OGRGeometryFactory::createFromGML( const char *pszData )
1440 :
1441 : {
1442 : OGRGeometryH hGeom;
1443 :
1444 57 : hGeom = OGR_G_CreateFromGML( pszData );
1445 :
1446 57 : return (OGRGeometry *) hGeom;
1447 : }
1448 :
1449 : /************************************************************************/
1450 : /* createFromGEOS() */
1451 : /************************************************************************/
1452 :
1453 : OGRGeometry *
1454 115 : OGRGeometryFactory::createFromGEOS( GEOSGeom geosGeom )
1455 :
1456 : {
1457 : #ifndef HAVE_GEOS
1458 :
1459 : CPLError( CE_Failure, CPLE_NotSupported,
1460 : "GEOS support not enabled." );
1461 : return NULL;
1462 :
1463 : #else
1464 :
1465 115 : size_t nSize = 0;
1466 115 : unsigned char *pabyBuf = NULL;
1467 115 : OGRGeometry *poGeometry = NULL;
1468 :
1469 : /* Special case as POINT EMPTY cannot be translated to WKB */
1470 115 : if (GEOSGeomTypeId(geosGeom) == GEOS_POINT &&
1471 : GEOSisEmpty(geosGeom))
1472 0 : return new OGRPoint();
1473 :
1474 : #if GEOS_VERSION_MAJOR > 3 || (GEOS_VERSION_MAJOR == 3 && GEOS_VERSION_MINOR >= 3)
1475 : int nCoordDim = GEOSGeom_getCoordinateDimension(geosGeom);
1476 : GEOSWKBWriter* wkbwriter = GEOSWKBWriter_create();
1477 : GEOSWKBWriter_setOutputDimension(wkbwriter, nCoordDim);
1478 : pabyBuf = GEOSWKBWriter_write( wkbwriter, geosGeom, &nSize );
1479 : GEOSWKBWriter_destroy(wkbwriter);
1480 : #else
1481 115 : pabyBuf = GEOSGeomToWKB_buf( geosGeom, &nSize );
1482 : #endif
1483 115 : if( pabyBuf == NULL || nSize == 0 )
1484 : {
1485 0 : return NULL;
1486 : }
1487 :
1488 115 : if( OGRGeometryFactory::createFromWkb( (unsigned char *) pabyBuf,
1489 : NULL, &poGeometry, (int) nSize )
1490 : != OGRERR_NONE )
1491 : {
1492 0 : poGeometry = NULL;
1493 : }
1494 :
1495 115 : if( pabyBuf != NULL )
1496 : {
1497 : #if GEOS_CAPI_VERSION_MAJOR >= 2 || (GEOS_CAPI_VERSION_MAJOR == 1 && GEOS_CAPI_VERSION_MINOR >= 6)
1498 : GEOSFree( pabyBuf );
1499 : #else
1500 115 : free( pabyBuf );
1501 : #endif
1502 : }
1503 :
1504 115 : return poGeometry;
1505 :
1506 : #endif /* HAVE_GEOS */
1507 : }
1508 :
1509 : /************************************************************************/
1510 : /* getGEOSGeometryFactory() */
1511 : /************************************************************************/
1512 :
1513 0 : void *OGRGeometryFactory::getGEOSGeometryFactory()
1514 :
1515 : {
1516 : // XXX - mloskot - What to do with this call
1517 : // after GEOS C++ API has been stripped?
1518 0 : return NULL;
1519 : }
1520 :
1521 : /************************************************************************/
1522 : /* haveGEOS() */
1523 : /************************************************************************/
1524 :
1525 : /**
1526 : * \brief Test if GEOS enabled.
1527 : *
1528 : * This static method returns TRUE if GEOS support is built into OGR,
1529 : * otherwise it returns FALSE.
1530 : *
1531 : * @return TRUE if available, otherwise FALSE.
1532 : */
1533 :
1534 51 : int OGRGeometryFactory::haveGEOS()
1535 :
1536 : {
1537 : #ifndef HAVE_GEOS
1538 : return FALSE;
1539 : #else
1540 51 : return TRUE;
1541 : #endif
1542 : }
1543 :
1544 : /************************************************************************/
1545 : /* createFromFgf() */
1546 : /************************************************************************/
1547 :
1548 : /**
1549 : * \brief Create a geometry object of the appropriate type from it's FGF (FDO Geometry Format) binary representation.
1550 : *
1551 : * Also note that this is a static method, and that there
1552 : * is no need to instantiate an OGRGeometryFactory object.
1553 : *
1554 : * The C function OGR_G_CreateFromFgf() is the same as this method.
1555 : *
1556 : * @param pabyData pointer to the input BLOB data.
1557 : * @param poSR pointer to the spatial reference to be assigned to the
1558 : * created geometry object. This may be NULL.
1559 : * @param ppoReturn the newly created geometry object will be assigned to the
1560 : * indicated pointer on return. This will be NULL in case
1561 : * of failure.
1562 : * @param nBytes the number of bytes available in pabyData.
1563 : * @param pnBytesConsumed if not NULL, it will be set to the number of bytes
1564 : * consumed (at most nBytes).
1565 : *
1566 : * @return OGRERR_NONE if all goes well, otherwise any of
1567 : * OGRERR_NOT_ENOUGH_DATA, OGRERR_UNSUPPORTED_GEOMETRY_TYPE, or
1568 : * OGRERR_CORRUPT_DATA may be returned.
1569 : */
1570 :
1571 : OGRErr OGRGeometryFactory::createFromFgf( unsigned char *pabyData,
1572 : OGRSpatialReference * poSR,
1573 : OGRGeometry **ppoReturn,
1574 : int nBytes,
1575 10 : int *pnBytesConsumed )
1576 :
1577 : {
1578 10 : OGRErr eErr = OGRERR_NONE;
1579 10 : OGRGeometry *poGeom = NULL;
1580 : GInt32 nGType, nGDim;
1581 10 : int nTupleSize = 0;
1582 10 : int iOrdinal = 0;
1583 :
1584 : (void) iOrdinal;
1585 :
1586 10 : *ppoReturn = NULL;
1587 :
1588 10 : if( nBytes < 4 )
1589 0 : return OGRERR_NOT_ENOUGH_DATA;
1590 :
1591 : /* -------------------------------------------------------------------- */
1592 : /* Decode the geometry type. */
1593 : /* -------------------------------------------------------------------- */
1594 10 : memcpy( &nGType, pabyData + 0, 4 );
1595 : CPL_LSBPTR32( &nGType );
1596 :
1597 10 : if( nGType < 0 || nGType > 13 )
1598 0 : return OGRERR_UNSUPPORTED_GEOMETRY_TYPE;
1599 :
1600 : /* -------------------------------------------------------------------- */
1601 : /* Decode the dimentionality if appropriate. */
1602 : /* -------------------------------------------------------------------- */
1603 10 : switch( nGType )
1604 : {
1605 : case 1: // Point
1606 : case 2: // LineString
1607 : case 3: // Polygon
1608 :
1609 8 : if( nBytes < 8 )
1610 0 : return OGRERR_NOT_ENOUGH_DATA;
1611 :
1612 8 : memcpy( &nGDim, pabyData + 4, 4 );
1613 : CPL_LSBPTR32( &nGDim );
1614 :
1615 8 : if( nGDim < 0 || nGDim > 3 )
1616 0 : return OGRERR_CORRUPT_DATA;
1617 :
1618 8 : nTupleSize = 2;
1619 8 : if( nGDim & 0x01 ) // Z
1620 1 : nTupleSize++;
1621 8 : if( nGDim & 0x02 ) // M
1622 0 : nTupleSize++;
1623 :
1624 : break;
1625 :
1626 : default:
1627 : break;
1628 : }
1629 :
1630 : /* -------------------------------------------------------------------- */
1631 : /* None */
1632 : /* -------------------------------------------------------------------- */
1633 10 : if( nGType == 0 )
1634 : {
1635 0 : if( pnBytesConsumed )
1636 0 : *pnBytesConsumed = 4;
1637 : }
1638 :
1639 : /* -------------------------------------------------------------------- */
1640 : /* Point */
1641 : /* -------------------------------------------------------------------- */
1642 10 : else if( nGType == 1 )
1643 : {
1644 : double adfTuple[4];
1645 :
1646 2 : if( nBytes < nTupleSize * 8 + 8 )
1647 0 : return OGRERR_NOT_ENOUGH_DATA;
1648 :
1649 2 : memcpy( adfTuple, pabyData + 8, nTupleSize*8 );
1650 : #ifdef CPL_MSB
1651 : for( iOrdinal = 0; iOrdinal < nTupleSize; iOrdinal++ )
1652 : CPL_SWAP64PTR( adfTuple + iOrdinal );
1653 : #endif
1654 2 : if( nTupleSize > 2 )
1655 1 : poGeom = new OGRPoint( adfTuple[0], adfTuple[1], adfTuple[2] );
1656 : else
1657 1 : poGeom = new OGRPoint( adfTuple[0], adfTuple[1] );
1658 :
1659 2 : if( pnBytesConsumed )
1660 0 : *pnBytesConsumed = 8 + nTupleSize * 8;
1661 : }
1662 :
1663 : /* -------------------------------------------------------------------- */
1664 : /* LineString */
1665 : /* -------------------------------------------------------------------- */
1666 8 : else if( nGType == 2 )
1667 : {
1668 : double adfTuple[4];
1669 : GInt32 nPointCount;
1670 : int iPoint;
1671 : OGRLineString *poLS;
1672 :
1673 2 : if( nBytes < 12 )
1674 0 : return OGRERR_NOT_ENOUGH_DATA;
1675 :
1676 2 : memcpy( &nPointCount, pabyData + 8, 4 );
1677 : CPL_LSBPTR32( &nPointCount );
1678 :
1679 2 : if (nPointCount < 0 || nPointCount > INT_MAX / (nTupleSize * 8))
1680 0 : return OGRERR_CORRUPT_DATA;
1681 :
1682 2 : if( nBytes - 12 < nTupleSize * 8 * nPointCount )
1683 0 : return OGRERR_NOT_ENOUGH_DATA;
1684 :
1685 2 : poGeom = poLS = new OGRLineString();
1686 2 : poLS->setNumPoints( nPointCount );
1687 :
1688 4 : for( iPoint = 0; iPoint < nPointCount; iPoint++ )
1689 : {
1690 : memcpy( adfTuple, pabyData + 12 + 8*nTupleSize*iPoint,
1691 2 : nTupleSize*8 );
1692 : #ifdef CPL_MSB
1693 : for( iOrdinal = 0; iOrdinal < nTupleSize; iOrdinal++ )
1694 : CPL_SWAP64PTR( adfTuple + iOrdinal );
1695 : #endif
1696 2 : if( nTupleSize > 2 )
1697 0 : poLS->setPoint( iPoint, adfTuple[0], adfTuple[1], adfTuple[2] );
1698 : else
1699 2 : poLS->setPoint( iPoint, adfTuple[0], adfTuple[1] );
1700 : }
1701 :
1702 2 : if( pnBytesConsumed )
1703 0 : *pnBytesConsumed = 12 + nTupleSize * 8 * nPointCount;
1704 : }
1705 :
1706 : /* -------------------------------------------------------------------- */
1707 : /* Polygon */
1708 : /* -------------------------------------------------------------------- */
1709 6 : else if( nGType == 3 )
1710 : {
1711 : double adfTuple[4];
1712 : GInt32 nPointCount;
1713 : GInt32 nRingCount;
1714 : int iPoint, iRing;
1715 : OGRLinearRing *poLR;
1716 : OGRPolygon *poPoly;
1717 : int nNextByte;
1718 :
1719 4 : if( nBytes < 12 )
1720 0 : return OGRERR_NOT_ENOUGH_DATA;
1721 :
1722 4 : memcpy( &nRingCount, pabyData + 8, 4 );
1723 : CPL_LSBPTR32( &nRingCount );
1724 :
1725 4 : if (nRingCount < 0 || nRingCount > INT_MAX / 4)
1726 0 : return OGRERR_CORRUPT_DATA;
1727 :
1728 : /* Each ring takes at least 4 bytes */
1729 4 : if (nBytes - 12 < nRingCount * 4)
1730 0 : return OGRERR_NOT_ENOUGH_DATA;
1731 :
1732 4 : nNextByte = 12;
1733 :
1734 4 : poGeom = poPoly = new OGRPolygon();
1735 :
1736 10 : for( iRing = 0; iRing < nRingCount; iRing++ )
1737 : {
1738 6 : if( nBytes - nNextByte < 4 )
1739 0 : return OGRERR_NOT_ENOUGH_DATA;
1740 :
1741 6 : memcpy( &nPointCount, pabyData + nNextByte, 4 );
1742 : CPL_LSBPTR32( &nPointCount );
1743 :
1744 6 : if (nPointCount < 0 || nPointCount > INT_MAX / (nTupleSize * 8))
1745 0 : return OGRERR_CORRUPT_DATA;
1746 :
1747 6 : nNextByte += 4;
1748 :
1749 6 : if( nBytes - nNextByte < nTupleSize * 8 * nPointCount )
1750 0 : return OGRERR_NOT_ENOUGH_DATA;
1751 :
1752 6 : poLR = new OGRLinearRing();
1753 6 : poLR->setNumPoints( nPointCount );
1754 :
1755 12 : for( iPoint = 0; iPoint < nPointCount; iPoint++ )
1756 : {
1757 6 : memcpy( adfTuple, pabyData + nNextByte, nTupleSize*8 );
1758 6 : nNextByte += nTupleSize * 8;
1759 :
1760 : #ifdef CPL_MSB
1761 : for( iOrdinal = 0; iOrdinal < nTupleSize; iOrdinal++ )
1762 : CPL_SWAP64PTR( adfTuple + iOrdinal );
1763 : #endif
1764 6 : if( nTupleSize > 2 )
1765 0 : poLR->setPoint( iPoint, adfTuple[0], adfTuple[1], adfTuple[2] );
1766 : else
1767 6 : poLR->setPoint( iPoint, adfTuple[0], adfTuple[1] );
1768 : }
1769 :
1770 6 : poPoly->addRingDirectly( poLR );
1771 : }
1772 :
1773 4 : if( pnBytesConsumed )
1774 2 : *pnBytesConsumed = nNextByte;
1775 : }
1776 :
1777 : /* -------------------------------------------------------------------- */
1778 : /* GeometryCollections of various kinds. */
1779 : /* -------------------------------------------------------------------- */
1780 4 : else if( nGType == 4 // MultiPoint
1781 : || nGType == 5 // MultiLineString
1782 : || nGType == 6 // MultiPolygon
1783 : || nGType == 7 ) // MultiGeometry
1784 : {
1785 2 : OGRGeometryCollection *poGC = NULL;
1786 : GInt32 nGeomCount;
1787 : int iGeom, nBytesUsed;
1788 :
1789 2 : if( nGType == 4 )
1790 0 : poGC = new OGRMultiPoint();
1791 2 : else if( nGType == 5 )
1792 0 : poGC = new OGRMultiLineString();
1793 2 : else if( nGType == 6 )
1794 0 : poGC = new OGRMultiPolygon();
1795 2 : else if( nGType == 7 )
1796 2 : poGC = new OGRGeometryCollection();
1797 :
1798 2 : if( nBytes < 8 )
1799 0 : return OGRERR_NOT_ENOUGH_DATA;
1800 :
1801 2 : memcpy( &nGeomCount, pabyData + 4, 4 );
1802 : CPL_LSBPTR32( &nGeomCount );
1803 :
1804 2 : if (nGeomCount < 0 || nGeomCount > INT_MAX / 4)
1805 0 : return OGRERR_CORRUPT_DATA;
1806 :
1807 : /* Each geometry takes at least 4 bytes */
1808 2 : if (nBytes - 8 < 4 * nGeomCount)
1809 0 : return OGRERR_NOT_ENOUGH_DATA;
1810 :
1811 2 : nBytesUsed = 8;
1812 :
1813 4 : for( iGeom = 0; iGeom < nGeomCount; iGeom++ )
1814 : {
1815 : int nThisGeomSize;
1816 2 : OGRGeometry *poThisGeom = NULL;
1817 :
1818 : eErr = createFromFgf( pabyData + nBytesUsed, poSR, &poThisGeom,
1819 2 : nBytes - nBytesUsed, &nThisGeomSize);
1820 2 : if( eErr != OGRERR_NONE )
1821 : {
1822 0 : delete poGC;
1823 0 : return eErr;
1824 : }
1825 :
1826 2 : nBytesUsed += nThisGeomSize;
1827 2 : eErr = poGC->addGeometryDirectly( poThisGeom );
1828 2 : if( eErr != OGRERR_NONE )
1829 : {
1830 0 : delete poGC;
1831 0 : return eErr;
1832 : }
1833 : }
1834 :
1835 2 : poGeom = poGC;
1836 2 : if( pnBytesConsumed )
1837 0 : *pnBytesConsumed = nBytesUsed;
1838 : }
1839 :
1840 : /* -------------------------------------------------------------------- */
1841 : /* Currently unsupported geometry. */
1842 : /* */
1843 : /* We need to add 10/11/12/13 curve types in some fashion. */
1844 : /* -------------------------------------------------------------------- */
1845 : else
1846 : {
1847 0 : return OGRERR_UNSUPPORTED_GEOMETRY_TYPE;
1848 : }
1849 :
1850 : /* -------------------------------------------------------------------- */
1851 : /* Assign spatial reference system. */
1852 : /* -------------------------------------------------------------------- */
1853 10 : if( eErr == OGRERR_NONE )
1854 : {
1855 10 : if( poGeom != NULL && poSR )
1856 0 : poGeom->assignSpatialReference( poSR );
1857 10 : *ppoReturn = poGeom;
1858 : }
1859 : else
1860 : {
1861 0 : delete poGeom;
1862 : }
1863 :
1864 10 : return eErr;
1865 : }
1866 :
1867 : /************************************************************************/
1868 : /* OGR_G_CreateFromFgf() */
1869 : /************************************************************************/
1870 :
1871 : OGRErr CPL_DLL OGR_G_CreateFromFgf( unsigned char *pabyData,
1872 : OGRSpatialReferenceH hSRS,
1873 : OGRGeometryH *phGeometry,
1874 0 : int nBytes, int *pnBytesConsumed )
1875 :
1876 : {
1877 : return OGRGeometryFactory::createFromFgf( pabyData,
1878 : (OGRSpatialReference *) hSRS,
1879 : (OGRGeometry **) phGeometry,
1880 0 : nBytes, pnBytesConsumed );
1881 : }
1882 :
1883 : /************************************************************************/
1884 : /* SplitLineStringAtDateline() */
1885 : /************************************************************************/
1886 :
1887 : #define SWAP_DBL(a,b) do { double tmp = a; a = b; b = tmp; } while(0)
1888 :
1889 : static void SplitLineStringAtDateline(OGRGeometryCollection* poMulti,
1890 1 : const OGRLineString* poLS)
1891 : {
1892 : int i;
1893 1 : int bIs3D = poLS->getCoordinateDimension() == 3;
1894 1 : OGRLineString* poNewLS = new OGRLineString();
1895 1 : poMulti->addGeometryDirectly(poNewLS);
1896 12 : for(i=0;i<poLS->getNumPoints();i++)
1897 : {
1898 11 : double dfX = poLS->getX(i);
1899 11 : if (i > 0 && fabs(dfX - poLS->getX(i-1)) > 350)
1900 : {
1901 2 : double dfX1 = poLS->getX(i-1);
1902 2 : double dfY1 = poLS->getY(i-1);
1903 2 : double dfZ1 = poLS->getY(i-1);
1904 2 : double dfX2 = poLS->getX(i);
1905 2 : double dfY2 = poLS->getY(i);
1906 2 : double dfZ2 = poLS->getY(i);
1907 :
1908 2 : if (dfX1 > -180 && dfX1 < -170 && dfX2 == 180 &&
1909 : i+1 < poLS->getNumPoints() &&
1910 : poLS->getX(i+1) > -180 && poLS->getX(i+1) < -170)
1911 : {
1912 0 : if( bIs3D )
1913 0 : poNewLS->addPoint(-180, poLS->getY(i), poLS->getZ(i));
1914 : else
1915 0 : poNewLS->addPoint(-180, poLS->getY(i));
1916 :
1917 0 : i++;
1918 :
1919 0 : if( bIs3D )
1920 0 : poNewLS->addPoint(poLS->getX(i), poLS->getY(i), poLS->getZ(i));
1921 : else
1922 0 : poNewLS->addPoint(poLS->getX(i), poLS->getY(i));
1923 0 : continue;
1924 : }
1925 2 : else if (dfX1 > 170 && dfX1 < 180 && dfX2 == -180 &&
1926 : i+1 < poLS->getNumPoints() &&
1927 : poLS->getX(i+1) > 170 && poLS->getX(i+1) < 180)
1928 : {
1929 0 : if( bIs3D )
1930 0 : poNewLS->addPoint(180, poLS->getY(i), poLS->getZ(i));
1931 : else
1932 0 : poNewLS->addPoint(180, poLS->getY(i));
1933 :
1934 0 : i++;
1935 :
1936 0 : if( bIs3D )
1937 0 : poNewLS->addPoint(poLS->getX(i), poLS->getY(i), poLS->getZ(i));
1938 : else
1939 0 : poNewLS->addPoint(poLS->getX(i), poLS->getY(i));
1940 0 : continue;
1941 : }
1942 :
1943 2 : if (dfX1 < -170 && dfX2 > 170)
1944 : {
1945 1 : SWAP_DBL(dfX1, dfX2);
1946 1 : SWAP_DBL(dfY1, dfY2);
1947 1 : SWAP_DBL(dfZ1, dfZ2);
1948 : }
1949 2 : if (dfX1 > 170 && dfX2 < -170)
1950 2 : dfX2 += 360;
1951 :
1952 4 : if (dfX1 <= 180 && dfX2 >= 180 && dfX1 < dfX2)
1953 : {
1954 2 : double dfRatio = (180 - dfX1) / (dfX2 - dfX1);
1955 2 : double dfY = dfRatio * dfY2 + (1 - dfRatio) * dfY1;
1956 2 : double dfZ = dfRatio * dfZ2 + (1 - dfRatio) * dfZ1;
1957 2 : if( bIs3D )
1958 0 : poNewLS->addPoint(poLS->getX(i-1) > 170 ? 180 : -180, dfY, dfZ);
1959 : else
1960 2 : poNewLS->addPoint(poLS->getX(i-1) > 170 ? 180 : -180, dfY);
1961 2 : poNewLS = new OGRLineString();
1962 2 : if( bIs3D )
1963 0 : poNewLS->addPoint(poLS->getX(i-1) > 170 ? -180 : 180, dfY, dfZ);
1964 : else
1965 2 : poNewLS->addPoint(poLS->getX(i-1) > 170 ? -180 : 180, dfY);
1966 2 : poMulti->addGeometryDirectly(poNewLS);
1967 : }
1968 : else
1969 : {
1970 0 : poNewLS = new OGRLineString();
1971 0 : poMulti->addGeometryDirectly(poNewLS);
1972 : }
1973 : }
1974 11 : if( bIs3D )
1975 0 : poNewLS->addPoint(dfX, poLS->getY(i), poLS->getZ(i));
1976 : else
1977 11 : poNewLS->addPoint(dfX, poLS->getY(i));
1978 : }
1979 1 : }
1980 :
1981 : /************************************************************************/
1982 : /* FixPolygonCoordinatesAtDateLine() */
1983 : /************************************************************************/
1984 :
1985 3 : static void FixPolygonCoordinatesAtDateLine(OGRPolygon* poPoly)
1986 : {
1987 : int i, iPart;
1988 6 : for(iPart = 0; iPart < 1 + poPoly->getNumInteriorRings(); iPart++)
1989 : {
1990 : OGRLineString* poLS = (iPart == 0) ? poPoly->getExteriorRing() :
1991 3 : poPoly->getInteriorRing(iPart-1);
1992 3 : int bGoEast = FALSE;
1993 3 : int bIs3D = poLS->getCoordinateDimension() == 3;
1994 31 : for(i=1;i<poLS->getNumPoints();i++)
1995 : {
1996 28 : double dfX = poLS->getX(i);
1997 28 : double dfPrevX = poLS->getX(i-1);
1998 28 : double dfDiffLong = fabs(dfX - dfPrevX);
1999 28 : if (dfDiffLong > 350)
2000 : {
2001 29 : if ((dfPrevX > 170 && dfX < -170) || (dfX < 0 && bGoEast))
2002 : {
2003 14 : dfX += 360;
2004 14 : bGoEast = TRUE;
2005 14 : if( bIs3D )
2006 0 : poLS->setPoint(i, dfX, poLS->getY(i), poLS->getZ(i));
2007 : else
2008 14 : poLS->setPoint(i, dfX, poLS->getY(i));
2009 : }
2010 2 : else if (dfPrevX < -170 && dfX > 170)
2011 : {
2012 : int j;
2013 6 : for(j=i-1;j>=0;j--)
2014 : {
2015 5 : dfX = poLS->getX(j);
2016 5 : if (dfX < 0)
2017 : {
2018 5 : if( bIs3D )
2019 0 : poLS->setPoint(j, dfX + 360, poLS->getY(j), poLS->getZ(j));
2020 : else
2021 5 : poLS->setPoint(j, dfX + 360, poLS->getY(j));
2022 : }
2023 : }
2024 1 : bGoEast = FALSE;
2025 : }
2026 : else
2027 : {
2028 0 : bGoEast = FALSE;
2029 : }
2030 : }
2031 : }
2032 : }
2033 3 : }
2034 :
2035 : /************************************************************************/
2036 : /* Sub360ToLon() */
2037 : /************************************************************************/
2038 :
2039 6 : static void Sub360ToLon( OGRGeometry* poGeom )
2040 : {
2041 6 : switch (wkbFlatten(poGeom->getGeometryType()))
2042 : {
2043 : case wkbPolygon:
2044 : case wkbMultiLineString:
2045 : case wkbMultiPolygon:
2046 : case wkbGeometryCollection:
2047 : {
2048 3 : int nSubGeomCount = OGR_G_GetGeometryCount((OGRGeometryH)poGeom);
2049 6 : for( int iGeom = 0; iGeom < nSubGeomCount; iGeom++ )
2050 : {
2051 3 : Sub360ToLon((OGRGeometry*)OGR_G_GetGeometryRef((OGRGeometryH)poGeom, iGeom));
2052 : }
2053 :
2054 3 : break;
2055 : }
2056 :
2057 : case wkbLineString:
2058 : case wkbLinearRing:
2059 : {
2060 3 : OGRLineString* poLineString = (OGRLineString* )poGeom;
2061 3 : int nPointCount = poLineString->getNumPoints();
2062 3 : int nCoordDim = poLineString->getCoordinateDimension();
2063 30 : for( int iPoint = 0; iPoint < nPointCount; iPoint++)
2064 : {
2065 27 : if (nCoordDim == 2)
2066 : poLineString->setPoint(iPoint,
2067 : poLineString->getX(iPoint) - 360,
2068 27 : poLineString->getY(iPoint));
2069 : else
2070 : poLineString->setPoint(iPoint,
2071 : poLineString->getX(iPoint) - 360,
2072 : poLineString->getY(iPoint),
2073 0 : poLineString->getZ(iPoint));
2074 : }
2075 : break;
2076 : }
2077 :
2078 : default:
2079 : break;
2080 : }
2081 6 : }
2082 :
2083 : /************************************************************************/
2084 : /* AddSimpleGeomToMulti() */
2085 : /************************************************************************/
2086 :
2087 : static void AddSimpleGeomToMulti(OGRGeometryCollection* poMulti,
2088 6 : const OGRGeometry* poGeom)
2089 : {
2090 6 : switch (wkbFlatten(poGeom->getGeometryType()))
2091 : {
2092 : case wkbPolygon:
2093 : case wkbLineString:
2094 6 : poMulti->addGeometry(poGeom);
2095 6 : break;
2096 :
2097 : case wkbMultiLineString:
2098 : case wkbMultiPolygon:
2099 : case wkbGeometryCollection:
2100 : {
2101 0 : int nSubGeomCount = OGR_G_GetGeometryCount((OGRGeometryH)poGeom);
2102 0 : for( int iGeom = 0; iGeom < nSubGeomCount; iGeom++ )
2103 : {
2104 : OGRGeometry* poSubGeom =
2105 0 : (OGRGeometry*)OGR_G_GetGeometryRef((OGRGeometryH)poGeom, iGeom);
2106 0 : AddSimpleGeomToMulti(poMulti, poSubGeom);
2107 : }
2108 : break;
2109 : }
2110 :
2111 : default:
2112 : break;
2113 : }
2114 6 : }
2115 :
2116 : /************************************************************************/
2117 : /* CutGeometryOnDateLineAndAddToMulti() */
2118 : /************************************************************************/
2119 :
2120 : static void CutGeometryOnDateLineAndAddToMulti(OGRGeometryCollection* poMulti,
2121 4 : const OGRGeometry* poGeom)
2122 : {
2123 4 : OGRwkbGeometryType eGeomType = wkbFlatten(poGeom->getGeometryType());
2124 4 : switch (eGeomType)
2125 : {
2126 : case wkbPolygon:
2127 : case wkbLineString:
2128 : {
2129 4 : int bWrapDateline = FALSE;
2130 4 : int bSplitLineStringAtDateline = FALSE;
2131 4 : OGREnvelope oEnvelope;
2132 :
2133 4 : poGeom->getEnvelope(&oEnvelope);
2134 :
2135 : /* Naive heuristics... Place to improvement... */
2136 4 : OGRGeometry* poDupGeom = NULL;
2137 :
2138 4 : if (oEnvelope.MinX > 170 && oEnvelope.MaxX > 180)
2139 : {
2140 : #ifndef HAVE_GEOS
2141 : CPLError( CE_Failure, CPLE_NotSupported,
2142 : "GEOS support not enabled." );
2143 : #else
2144 0 : bWrapDateline = TRUE;
2145 : #endif
2146 : }
2147 : else
2148 : {
2149 : OGRLineString* poLS;
2150 4 : if (eGeomType == wkbPolygon)
2151 3 : poLS = ((OGRPolygon*)poGeom)->getExteriorRing();
2152 : else
2153 1 : poLS = (OGRLineString*)poGeom;
2154 4 : if (poLS)
2155 : {
2156 : int i;
2157 4 : double dfMaxSmallDiffLong = 0;
2158 4 : int bHasBigDiff = FALSE;
2159 : /* Detect big gaps in longitude */
2160 42 : for(i=1;i<poLS->getNumPoints();i++)
2161 : {
2162 38 : double dfPrevX = poLS->getX(i-1);
2163 38 : double dfX = poLS->getX(i);
2164 38 : double dfDiffLong = fabs(dfX - dfPrevX);
2165 46 : if (dfDiffLong > 350 &&
2166 : ((dfX > 170 && dfPrevX < -170) || (dfPrevX > 170 && dfX < -170)))
2167 8 : bHasBigDiff = TRUE;
2168 30 : else if (dfDiffLong > dfMaxSmallDiffLong)
2169 7 : dfMaxSmallDiffLong = dfDiffLong;
2170 : }
2171 4 : if (bHasBigDiff && dfMaxSmallDiffLong < 10)
2172 : {
2173 4 : if (eGeomType == wkbLineString)
2174 1 : bSplitLineStringAtDateline = TRUE;
2175 : else
2176 : {
2177 : #ifndef HAVE_GEOS
2178 : CPLError( CE_Failure, CPLE_NotSupported,
2179 : "GEOS support not enabled." );
2180 : #else
2181 3 : bWrapDateline = TRUE;
2182 3 : poDupGeom = poGeom->clone();
2183 3 : FixPolygonCoordinatesAtDateLine((OGRPolygon*)poDupGeom);
2184 : #endif
2185 : }
2186 : }
2187 : }
2188 : }
2189 :
2190 4 : if (bSplitLineStringAtDateline)
2191 : {
2192 1 : SplitLineStringAtDateline(poMulti, (OGRLineString*)poGeom);
2193 : }
2194 3 : else if (bWrapDateline)
2195 : {
2196 3 : const OGRGeometry* poWorkGeom = (poDupGeom) ? poDupGeom : poGeom;
2197 3 : OGRGeometry* poRectangle1 = NULL;
2198 3 : OGRGeometry* poRectangle2 = NULL;
2199 3 : const char* pszWKT1 = "POLYGON((0 90,180 90,180 -90,0 -90,0 90))";
2200 3 : const char* pszWKT2 = "POLYGON((180 90,360 90,360 -90,180 -90,180 90))";
2201 3 : OGRGeometryFactory::createFromWkt((char**)&pszWKT1, NULL, &poRectangle1);
2202 3 : OGRGeometryFactory::createFromWkt((char**)&pszWKT2, NULL, &poRectangle2);
2203 3 : OGRGeometry* poGeom1 = poWorkGeom->Intersection(poRectangle1);
2204 3 : OGRGeometry* poGeom2 = poWorkGeom->Intersection(poRectangle2);
2205 3 : delete poRectangle1;
2206 3 : delete poRectangle2;
2207 :
2208 6 : if (poGeom1 != NULL && poGeom2 != NULL)
2209 : {
2210 3 : AddSimpleGeomToMulti(poMulti, poGeom1);
2211 3 : Sub360ToLon(poGeom2);
2212 3 : AddSimpleGeomToMulti(poMulti, poGeom2);
2213 : }
2214 : else
2215 : {
2216 0 : AddSimpleGeomToMulti(poMulti, poGeom);
2217 : }
2218 :
2219 3 : delete poGeom1;
2220 3 : delete poGeom2;
2221 3 : delete poDupGeom;
2222 : }
2223 : else
2224 : {
2225 0 : poMulti->addGeometry(poGeom);
2226 : }
2227 4 : break;
2228 : }
2229 :
2230 : case wkbMultiLineString:
2231 : case wkbMultiPolygon:
2232 : case wkbGeometryCollection:
2233 : {
2234 0 : int nSubGeomCount = OGR_G_GetGeometryCount((OGRGeometryH)poGeom);
2235 0 : for( int iGeom = 0; iGeom < nSubGeomCount; iGeom++ )
2236 : {
2237 : OGRGeometry* poSubGeom =
2238 0 : (OGRGeometry*)OGR_G_GetGeometryRef((OGRGeometryH)poGeom, iGeom);
2239 0 : CutGeometryOnDateLineAndAddToMulti(poMulti, poSubGeom);
2240 : }
2241 : break;
2242 : }
2243 :
2244 : default:
2245 : break;
2246 : }
2247 4 : }
2248 :
2249 : /************************************************************************/
2250 : /* transformWithOptions() */
2251 : /************************************************************************/
2252 :
2253 : OGRGeometry* OGRGeometryFactory::transformWithOptions( const OGRGeometry* poSrcGeom,
2254 : OGRCoordinateTransformation *poCT,
2255 14 : char** papszOptions )
2256 : {
2257 14 : OGRGeometry* poDstGeom = poSrcGeom->clone();
2258 14 : if (poCT != NULL)
2259 : {
2260 11 : OGRErr eErr = poDstGeom->transform(poCT);
2261 11 : if (eErr != OGRERR_NONE)
2262 : {
2263 0 : delete poDstGeom;
2264 0 : return NULL;
2265 : }
2266 : }
2267 :
2268 14 : if (CSLTestBoolean(CSLFetchNameValueDef(papszOptions, "WRAPDATELINE", "NO")))
2269 : {
2270 4 : OGRwkbGeometryType eType = wkbFlatten(poSrcGeom->getGeometryType());
2271 : OGRwkbGeometryType eNewType;
2272 7 : if (eType == wkbPolygon || eType == wkbMultiPolygon)
2273 3 : eNewType = wkbMultiPolygon;
2274 2 : else if (eType == wkbLineString || eType == wkbMultiLineString)
2275 1 : eNewType = wkbMultiLineString;
2276 : else
2277 0 : eNewType = wkbGeometryCollection;
2278 :
2279 : OGRGeometryCollection* poMulti =
2280 4 : (OGRGeometryCollection* )createGeometry(eNewType);
2281 :
2282 4 : CutGeometryOnDateLineAndAddToMulti(poMulti, poDstGeom);
2283 :
2284 4 : if (poMulti->getNumGeometries() == 0)
2285 : {
2286 0 : delete poMulti;
2287 : }
2288 4 : else if (poMulti->getNumGeometries() == 1)
2289 : {
2290 0 : delete poDstGeom;
2291 0 : poDstGeom = poMulti->getGeometryRef(0)->clone();
2292 0 : delete poMulti;
2293 : }
2294 : else
2295 : {
2296 4 : delete poDstGeom;
2297 4 : poDstGeom = poMulti;
2298 : }
2299 : }
2300 :
2301 14 : return poDstGeom;
2302 : }
2303 :
2304 : /************************************************************************/
2305 : /* approximateArcAngles() */
2306 : /************************************************************************/
2307 :
2308 : /**
2309 : * Stroke arc to linestring.
2310 : *
2311 : * Stroke an arc of a circle to a linestring based on a center
2312 : * point, radius, start angle and end angle, all angles in degrees.
2313 : *
2314 : * If the dfMaxAngleStepSizeDegrees is zero, then a default value will be
2315 : * used. This is currently 4 degrees unless the user has overridden the
2316 : * value with the OGR_ARC_STEPSIZE configuration variable.
2317 : *
2318 : * @see CPLSetConfigOption()
2319 : *
2320 : * @param dfCenterX center X
2321 : * @param dfCenterY center Y
2322 : * @param dfZ center Z
2323 : * @param dfPrimaryRadius X radius of ellipse.
2324 : * @param dfSecondaryRadius Y radius of ellipse.
2325 : * @param dfRotation rotation of the ellipse clockwise.
2326 : * @param dfStartAngle angle to first point on arc (clockwise of X-positive)
2327 : * @param dfEndAngle angle to last point on arc (clockwise of X-positive)
2328 : * @param dfMaxAngleStepSizeDegrees the largest step in degrees along the
2329 : * arc, zero to use the default setting.
2330 : *
2331 : * @return OGRLineString geometry representing an approximation of the arc.
2332 : */
2333 :
2334 : OGRGeometry* OGRGeometryFactory::approximateArcAngles(
2335 : double dfCenterX, double dfCenterY, double dfZ,
2336 : double dfPrimaryRadius, double dfSecondaryRadius, double dfRotation,
2337 : double dfStartAngle, double dfEndAngle,
2338 15 : double dfMaxAngleStepSizeDegrees )
2339 :
2340 : {
2341 : double dfSlice;
2342 : int iPoint, nVertexCount;
2343 15 : OGRLineString *poLine = new OGRLineString();
2344 15 : double dfRotationRadians = dfRotation * PI / 180.0;
2345 :
2346 : // support default arc step setting.
2347 15 : if( dfMaxAngleStepSizeDegrees == 0 )
2348 : {
2349 : dfMaxAngleStepSizeDegrees =
2350 14 : atof(CPLGetConfigOption("OGR_ARC_STEPSIZE","4"));
2351 : }
2352 :
2353 : // switch direction
2354 15 : dfStartAngle *= -1;
2355 15 : dfEndAngle *= -1;
2356 :
2357 : // Figure out the number of slices to make this into.
2358 : nVertexCount = (int)
2359 15 : ceil(fabs(dfEndAngle - dfStartAngle)/dfMaxAngleStepSizeDegrees) + 1;
2360 15 : nVertexCount = MAX(2,nVertexCount);
2361 15 : dfSlice = (dfEndAngle-dfStartAngle)/(nVertexCount-1);
2362 :
2363 : /* -------------------------------------------------------------------- */
2364 : /* Compute the interpolated points. */
2365 : /* -------------------------------------------------------------------- */
2366 723 : for( iPoint=0; iPoint < nVertexCount; iPoint++ )
2367 : {
2368 : double dfAngleOnEllipse;
2369 : double dfArcX, dfArcY;
2370 : double dfEllipseX, dfEllipseY;
2371 :
2372 708 : dfAngleOnEllipse = (dfStartAngle + iPoint * dfSlice) * PI / 180.0;
2373 :
2374 : // Compute position on the unrotated ellipse.
2375 708 : dfEllipseX = cos(dfAngleOnEllipse) * dfPrimaryRadius;
2376 708 : dfEllipseY = sin(dfAngleOnEllipse) * dfSecondaryRadius;
2377 :
2378 : // Rotate this position around the center of the ellipse.
2379 : dfArcX = dfCenterX
2380 : + dfEllipseX * cos(dfRotationRadians)
2381 708 : + dfEllipseY * sin(dfRotationRadians);
2382 : dfArcY = dfCenterY
2383 : - dfEllipseX * sin(dfRotationRadians)
2384 708 : + dfEllipseY * cos(dfRotationRadians);
2385 :
2386 708 : poLine->setPoint( iPoint, dfArcX, dfArcY, dfZ );
2387 : }
2388 :
2389 15 : return poLine;
2390 : }
2391 :
2392 : /************************************************************************/
2393 : /* OGR_G_ApproximateArcAngles() */
2394 : /************************************************************************/
2395 :
2396 : OGRGeometryH CPL_DLL
2397 : OGR_G_ApproximateArcAngles(
2398 : double dfCenterX, double dfCenterY, double dfZ,
2399 : double dfPrimaryRadius, double dfSecondaryAxis, double dfRotation,
2400 : double dfStartAngle, double dfEndAngle,
2401 1 : double dfMaxAngleStepSizeDegrees )
2402 :
2403 : {
2404 : return (OGRGeometryH) OGRGeometryFactory::approximateArcAngles(
2405 : dfCenterX, dfCenterY, dfZ,
2406 : dfPrimaryRadius, dfSecondaryAxis, dfRotation,
2407 1 : dfStartAngle, dfEndAngle, dfMaxAngleStepSizeDegrees );
2408 : }
|