1 : /******************************************************************************
2 : * $Id: ograssemblepolygon.cpp 17556 2009-08-21 17:12:30Z rouault $
3 : *
4 : * Project: S-57 Reader
5 : * Purpose: Implements polygon assembly from a bunch of arcs.
6 : * Author: Frank Warmerdam, warmerdam@pobox.com
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 1999, Frank Warmerdam
10 : *
11 : * Permission is hereby granted, free of charge, to any person obtaining a
12 : * copy of this software and associated documentation files (the "Software"),
13 : * to deal in the Software without restriction, including without limitation
14 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
15 : * and/or sell copies of the Software, and to permit persons to whom the
16 : * Software is furnished to do so, subject to the following conditions:
17 : *
18 : * The above copyright notice and this permission notice shall be included
19 : * in all copies or substantial portions of the Software.
20 : *
21 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
22 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
23 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
24 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
25 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
26 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
27 : * DEALINGS IN THE SOFTWARE.
28 : ****************************************************************************/
29 :
30 : #include "ogr_geometry.h"
31 : #include "ogr_api.h"
32 : #include "cpl_conv.h"
33 :
34 : CPL_CVSID("$Id: ograssemblepolygon.cpp 17556 2009-08-21 17:12:30Z rouault $");
35 :
36 : /************************************************************************/
37 : /* CheckPoints() */
38 : /* */
39 : /* Check if two points are closer than the current best */
40 : /* distance. Update the current best distance if they are. */
41 : /************************************************************************/
42 :
43 266 : static int CheckPoints( OGRLineString *poLine1, int iPoint1,
44 : OGRLineString *poLine2, int iPoint2,
45 : double *pdfDistance )
46 :
47 : {
48 : double dfDeltaX, dfDeltaY, dfDistance;
49 :
50 266 : if( pdfDistance == NULL || *pdfDistance == 0 )
51 : return poLine1->getX(iPoint1) == poLine2->getX(iPoint2)
52 266 : && poLine1->getY(iPoint1) == poLine2->getY(iPoint2);
53 :
54 0 : dfDeltaX = poLine1->getX(iPoint1) - poLine2->getX(iPoint2);
55 0 : dfDeltaY = poLine1->getY(iPoint1) - poLine2->getY(iPoint2);
56 :
57 0 : dfDeltaX = ABS(dfDeltaX);
58 0 : dfDeltaY = ABS(dfDeltaY);
59 :
60 0 : if( dfDeltaX > *pdfDistance || dfDeltaY > *pdfDistance )
61 0 : return FALSE;
62 :
63 0 : dfDistance = sqrt(dfDeltaX*dfDeltaX + dfDeltaY*dfDeltaY);
64 :
65 0 : if( dfDistance < *pdfDistance )
66 : {
67 0 : *pdfDistance = dfDistance;
68 0 : return TRUE;
69 : }
70 : else
71 0 : return FALSE;
72 : }
73 :
74 : /************************************************************************/
75 : /* AddEdgeToRing() */
76 : /************************************************************************/
77 :
78 25 : static void AddEdgeToRing( OGRLinearRing * poRing, OGRLineString * poLine,
79 : int bReverse )
80 :
81 : {
82 : /* -------------------------------------------------------------------- */
83 : /* Establish order and range of traverse. */
84 : /* -------------------------------------------------------------------- */
85 25 : int iStart=0, iEnd=0, iStep=0;
86 25 : int nVertToAdd = poLine->getNumPoints();
87 :
88 25 : if( !bReverse )
89 : {
90 19 : iStart = 0;
91 19 : iEnd = nVertToAdd - 1;
92 19 : iStep = 1;
93 : }
94 : else
95 : {
96 6 : iStart = nVertToAdd - 1;
97 6 : iEnd = 0;
98 6 : iStep = -1;
99 : }
100 :
101 : /* -------------------------------------------------------------------- */
102 : /* Should we skip a repeating vertex? */
103 : /* -------------------------------------------------------------------- */
104 25 : if( poRing->getNumPoints() > 0
105 : && CheckPoints( poRing, poRing->getNumPoints()-1,
106 : poLine, iStart, NULL ) )
107 : {
108 22 : iStart += iStep;
109 : }
110 :
111 25 : poRing->addSubLineString( poLine, iStart, iEnd );
112 25 : }
113 :
114 : /************************************************************************/
115 : /* OGRBuildPolygonFromEdges() */
116 : /************************************************************************/
117 :
118 : /**
119 : * Build a ring from a bunch of arcs.
120 : *
121 : * @param hLines handle to an OGRGeometryCollection (or OGRMultiLineString) containing the line string geometries to be built into rings.
122 : * @param bBestEffort not yet implemented???.
123 : * @param bAutoClose indicates if the ring should be close when first and
124 : * last points of the ring are the same.
125 : * @param dfTolerance tolerance into which two arcs are considered
126 : * close enough to be joined.
127 : * @param peErr OGRERR_NONE on success, or OGRERR_FAILURE on failure.
128 : * @return an handle to the new geometry, a polygon.
129 : *
130 : */
131 :
132 5 : OGRGeometryH OGRBuildPolygonFromEdges( OGRGeometryH hLines,
133 : int bBestEffort,
134 : int bAutoClose,
135 : double dfTolerance,
136 : OGRErr * peErr )
137 :
138 : {
139 5 : if( hLines == NULL )
140 : {
141 0 : if (peErr != NULL)
142 0 : *peErr = OGRERR_NONE;
143 0 : return NULL;
144 : }
145 :
146 : /* -------------------------------------------------------------------- */
147 : /* Check for the case of a geometrycollection that can be */
148 : /* promoted to MultiLineString. */
149 : /* -------------------------------------------------------------------- */
150 5 : OGRGeometry* poGeom = (OGRGeometry*) hLines;
151 5 : if( wkbFlatten(poGeom->getGeometryType()) == wkbGeometryCollection )
152 : {
153 : int iGeom;
154 3 : OGRGeometryCollection *poGC = (OGRGeometryCollection *) poGeom;
155 :
156 23 : for( iGeom = 0; iGeom < poGC->getNumGeometries(); iGeom++ )
157 : {
158 21 : if( wkbFlatten(poGC->getGeometryRef(iGeom)->getGeometryType())
159 : != wkbLineString )
160 : {
161 1 : if (peErr != NULL)
162 1 : *peErr = OGRERR_FAILURE;
163 : CPLError(CE_Failure, CPLE_NotSupported,
164 1 : "The geometry collection contains non line string geometries");
165 1 : return NULL;
166 : }
167 : }
168 : }
169 2 : else if (wkbFlatten(poGeom->getGeometryType()) != wkbMultiLineString )
170 : {
171 1 : if (peErr != NULL)
172 1 : *peErr = OGRERR_FAILURE;
173 : CPLError(CE_Failure, CPLE_NotSupported,
174 : "The passed geometry is not an OGRGeometryCollection (or OGRMultiLineString) "
175 1 : "containing line string geometries");
176 1 : return NULL;
177 : }
178 :
179 3 : int bSuccess = TRUE;
180 3 : OGRGeometryCollection *poLines = (OGRGeometryCollection *) hLines;
181 3 : OGRPolygon *poPolygon = new OGRPolygon();
182 :
183 : (void) bBestEffort;
184 :
185 : /* -------------------------------------------------------------------- */
186 : /* Setup array of line markers indicating if they have been */
187 : /* added to a ring yet. */
188 : /* -------------------------------------------------------------------- */
189 3 : int nEdges = poLines->getNumGeometries();
190 3 : int *panEdgeConsumed, nRemainingEdges = nEdges;
191 :
192 3 : panEdgeConsumed = (int *) CPLCalloc(sizeof(int),nEdges);
193 :
194 : /* ==================================================================== */
195 : /* Loop generating rings. */
196 : /* ==================================================================== */
197 9 : while( nRemainingEdges > 0 )
198 : {
199 : int iEdge;
200 : OGRLineString *poLine;
201 :
202 : /* -------------------------------------------------------------------- */
203 : /* Find the first unconsumed edge. */
204 : /* -------------------------------------------------------------------- */
205 3 : for( iEdge = 0; panEdgeConsumed[iEdge]; iEdge++ ) {}
206 :
207 3 : poLine = (OGRLineString *) poLines->getGeometryRef(iEdge);
208 :
209 : /* -------------------------------------------------------------------- */
210 : /* Start a new ring, copying in the current line directly */
211 : /* -------------------------------------------------------------------- */
212 3 : OGRLinearRing *poRing = new OGRLinearRing();
213 :
214 3 : AddEdgeToRing( poRing, poLine, FALSE );
215 :
216 3 : panEdgeConsumed[iEdge] = TRUE;
217 3 : nRemainingEdges--;
218 :
219 : /* ==================================================================== */
220 : /* Loop adding edges to this ring until we make a whole pass */
221 : /* within finding anything to add. */
222 : /* ==================================================================== */
223 3 : int bWorkDone = TRUE;
224 3 : double dfBestDist = dfTolerance;
225 :
226 28 : while( !CheckPoints(poRing,0,poRing,poRing->getNumPoints()-1,NULL)
227 : && nRemainingEdges > 0
228 : && bWorkDone )
229 : {
230 22 : int iBestEdge = -1, bReverse = FALSE;
231 :
232 22 : bWorkDone = FALSE;
233 22 : dfBestDist = dfTolerance;
234 :
235 : // We consider linking the end to the beginning. If this is
236 : // closer than any other option we will just close the loop.
237 :
238 : //CheckPoints(poRing,0,poRing,poRing->getNumPoints()-1,&dfBestDist);
239 :
240 : // Find unused edge with end point closest to our loose end.
241 238 : for( iEdge = 0; iEdge < nEdges; iEdge++ )
242 : {
243 216 : if( panEdgeConsumed[iEdge] )
244 108 : continue;
245 :
246 108 : poLine = (OGRLineString *) poLines->getGeometryRef(iEdge);
247 :
248 108 : if( CheckPoints(poLine,0,poRing,poRing->getNumPoints()-1,
249 : &dfBestDist) )
250 : {
251 16 : iBestEdge = iEdge;
252 16 : bReverse = FALSE;
253 : }
254 108 : if( CheckPoints(poLine,poLine->getNumPoints()-1,
255 : poRing,poRing->getNumPoints()-1,
256 : &dfBestDist) )
257 : {
258 6 : iBestEdge = iEdge;
259 6 : bReverse = TRUE;
260 : }
261 : }
262 :
263 : // We found one within tolerance - add it.
264 22 : if( iBestEdge != -1 )
265 : {
266 : poLine = (OGRLineString *)
267 22 : poLines->getGeometryRef(iBestEdge);
268 :
269 22 : AddEdgeToRing( poRing, poLine, bReverse );
270 :
271 22 : panEdgeConsumed[iBestEdge] = TRUE;
272 22 : nRemainingEdges--;
273 22 : bWorkDone = TRUE;
274 : }
275 : }
276 :
277 : /* -------------------------------------------------------------------- */
278 : /* Did we fail to complete the ring? */
279 : /* -------------------------------------------------------------------- */
280 3 : dfBestDist = dfTolerance;
281 :
282 3 : if( !CheckPoints(poRing,0,poRing,poRing->getNumPoints()-1,
283 : &dfBestDist) )
284 : {
285 : CPLDebug( "OGR",
286 : "Failed to close ring %d.\n"
287 : "End Points are: (%.8f,%.7f) and (%.7f,%.7f)\n",
288 : poPolygon->getNumInteriorRings()+1,
289 : poRing->getX(0), poRing->getY(0),
290 : poRing->getX(poRing->getNumPoints()-1),
291 0 : poRing->getY(poRing->getNumPoints()-1) );
292 :
293 0 : bSuccess = FALSE;
294 : }
295 :
296 : /* -------------------------------------------------------------------- */
297 : /* Do we need to auto-close this ring? */
298 : /* -------------------------------------------------------------------- */
299 3 : if( bAutoClose
300 : && !CheckPoints(poRing,0,poRing,poRing->getNumPoints()-1,NULL) )
301 : {
302 : poRing->addPoint( poRing->getX(0),
303 : poRing->getY(0),
304 0 : poRing->getZ(0));
305 : }
306 :
307 3 : poPolygon->addRingDirectly( poRing );
308 : } /* next ring */
309 :
310 : /* -------------------------------------------------------------------- */
311 : /* Cleanup. */
312 : /* -------------------------------------------------------------------- */
313 3 : CPLFree( panEdgeConsumed );
314 :
315 : // Eventually we should at least identify the external ring properly,
316 : // perhaps even ordering the direction of rings, though this isn't
317 : // required by the OGC geometry model.
318 :
319 3 : if( peErr != NULL )
320 : {
321 3 : if( bSuccess )
322 3 : *peErr = OGRERR_NONE;
323 : else
324 0 : *peErr = OGRERR_FAILURE;
325 : }
326 :
327 3 : return (OGRGeometryH) poPolygon;
328 : }
329 :
330 :
331 :
|