Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: SDTS Translator
4 : * Purpose: Implementation of SDTSPolygonReader and SDTSRawPolygon classes.
5 : * Author: Frank Warmerdam, warmerdam@pobox.com
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 1999, Frank Warmerdam
9 : *
10 : * Permission is hereby granted, free of charge, to any person obtaining a
11 : * copy of this software and associated documentation files (the "Software"),
12 : * to deal in the Software without restriction, including without limitation
13 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
14 : * and/or sell copies of the Software, and to permit persons to whom the
15 : * Software is furnished to do so, subject to the following conditions:
16 : *
17 : * The above copyright notice and this permission notice shall be included
18 : * in all copies or substantial portions of the Software.
19 : *
20 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
21 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
22 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
23 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
24 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
25 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
26 : * DEALINGS IN THE SOFTWARE.
27 : ****************************************************************************/
28 :
29 : #include "sdts_al.h"
30 :
31 : #include <cmath>
32 :
33 : /************************************************************************/
34 : /* ==================================================================== */
35 : /* SDTSRawPolygon */
36 : /* */
37 : /* This is a simple class for holding the data related with a */
38 : /* polygon feature. */
39 : /* ==================================================================== */
40 : /************************************************************************/
41 :
42 : /************************************************************************/
43 : /* SDTSRawPolygon() */
44 : /************************************************************************/
45 :
46 35 : SDTSRawPolygon::SDTSRawPolygon()
47 : : nEdges(0), papoEdges(nullptr), nRings(0), nVertices(0),
48 35 : panRingStart(nullptr), padfX(nullptr), padfY(nullptr), padfZ(nullptr)
49 : {
50 35 : nAttributes = 0;
51 35 : }
52 :
53 : /************************************************************************/
54 : /* ~SDTSRawPolygon() */
55 : /************************************************************************/
56 :
57 70 : SDTSRawPolygon::~SDTSRawPolygon()
58 :
59 : {
60 35 : CPLFree(papoEdges);
61 35 : CPLFree(panRingStart);
62 35 : CPLFree(padfX);
63 35 : CPLFree(padfY);
64 35 : CPLFree(padfZ);
65 70 : }
66 :
67 : /************************************************************************/
68 : /* Read() */
69 : /* */
70 : /* Read a record from the passed SDTSPolygonReader, and assign the */
71 : /* values from that record to this object. This is the bulk of */
72 : /* the work in this whole file. */
73 : /************************************************************************/
74 :
75 35 : int SDTSRawPolygon::Read(DDFRecord *poRecord)
76 :
77 : {
78 : /* ==================================================================== */
79 : /* Loop over fields in this record, looking for those we */
80 : /* recognise, and need. */
81 : /* ==================================================================== */
82 105 : for (int iField = 0; iField < poRecord->GetFieldCount(); iField++)
83 : {
84 70 : DDFField *poField = poRecord->GetField(iField);
85 70 : if (poField == nullptr)
86 0 : return FALSE;
87 70 : DDFFieldDefn *poFieldDefn = poField->GetFieldDefn();
88 70 : if (poFieldDefn == nullptr)
89 0 : return FALSE;
90 :
91 70 : const char *pszFieldName = poFieldDefn->GetName();
92 :
93 70 : if (EQUAL(pszFieldName, "POLY"))
94 : {
95 35 : oModId.Set(poField);
96 : }
97 :
98 35 : else if (EQUAL(pszFieldName, "ATID"))
99 : {
100 0 : ApplyATID(poField);
101 : }
102 : }
103 :
104 35 : return TRUE;
105 : }
106 :
107 : /************************************************************************/
108 : /* AddEdge() */
109 : /************************************************************************/
110 :
111 50 : void SDTSRawPolygon::AddEdge(SDTSRawLine *poNewLine)
112 :
113 : {
114 50 : nEdges++;
115 :
116 50 : papoEdges = reinterpret_cast<SDTSRawLine **>(
117 50 : CPLRealloc(papoEdges, sizeof(void *) * nEdges));
118 50 : papoEdges[nEdges - 1] = poNewLine;
119 50 : }
120 :
121 : /************************************************************************/
122 : /* AddEdgeToRing() */
123 : /************************************************************************/
124 :
125 52 : void SDTSRawPolygon::AddEdgeToRing(int nVertToAdd, double *padfXToAdd,
126 : double *padfYToAdd, double *padfZToAdd,
127 : int bReverse, int bDropVertex)
128 :
129 : {
130 52 : int iStart = 0;
131 52 : int iEnd = nVertToAdd - 1;
132 52 : int iStep = 1;
133 :
134 52 : if (bDropVertex && bReverse)
135 : {
136 4 : iStart = nVertToAdd - 2;
137 4 : iEnd = 0;
138 4 : iStep = -1;
139 : }
140 48 : else if (bDropVertex && !bReverse)
141 : {
142 29 : iStart = 1;
143 29 : iEnd = nVertToAdd - 1;
144 29 : iStep = 1;
145 : }
146 19 : else if (!bDropVertex && !bReverse)
147 : {
148 19 : iStart = 0;
149 19 : iEnd = nVertToAdd - 1;
150 19 : iStep = 1;
151 : }
152 0 : else if (!bDropVertex && bReverse)
153 : {
154 0 : iStart = nVertToAdd - 1;
155 0 : iEnd = 0;
156 0 : iStep = -1;
157 : }
158 :
159 1211 : for (int i = iStart; i != (iEnd + iStep); i += iStep)
160 : {
161 1159 : padfX[nVertices] = padfXToAdd[i];
162 1159 : padfY[nVertices] = padfYToAdd[i];
163 1159 : padfZ[nVertices] = padfZToAdd[i];
164 :
165 1159 : nVertices++;
166 : }
167 52 : }
168 :
169 : /************************************************************************/
170 : /* AssembleRings() */
171 : /************************************************************************/
172 :
173 : /**
174 : * Form border lines (arcs) into outer and inner rings.
175 : *
176 : * See SDTSPolygonReader::AssemblePolygons() for a simple one step process
177 : * to assembling geometry for all polygons in a transfer.
178 : *
179 : * This method will assemble the lines attached to a polygon into
180 : * an outer ring, and zero or more inner rings. Before calling it is
181 : * necessary that all the lines associated with this polygon have already
182 : * been attached. Normally this is accomplished by calling
183 : * SDTSLineReader::AttachToPolygons() on all line layers that might
184 : * contain edges related to this layer.
185 : *
186 : * This method then forms the lines into rings. Rings are formed by:
187 : * <ol>
188 : * <li> Take a previously unconsumed line, and start a ring with it. Mark
189 : * it as consumed, and keep track of its start and end node ids as
190 : * being the start and end node ids of the ring.
191 : * <li> If the rings start id is the same as the end node id then this ring
192 : * is completely formed, return to step 1.
193 : * <li> Search all unconsumed lines for a line with the same start or end
194 : * node id as the rings current node id. If none are found then the
195 : * assembly has failed. Return to step 1 but report failure on
196 : * completion.
197 : * <li> Once found, add the line to the current ring, dropping the duplicated
198 : * vertex and reverse order if necessary. Mark the line as consumed,
199 : * and update the rings end node id accordingly.
200 : * <li> go to step 2.
201 : * </ol>
202 : *
203 : * Once ring assembly from lines is complete, another pass is made to
204 : * order the rings such that the exterior ring is first, the first ring
205 : * has counter-clockwise vertex ordering and the inner rings have clockwise
206 : * vertex ordering. This is accomplished based on the assumption that the
207 : * outer ring has the largest area, and using the +/- sign of area to establish
208 : * direction of rings.
209 : *
210 : * @return TRUE if all rings assembled without problems or FALSE if a problem
211 : * occurred. If a problem occurs rings are still formed from all lines, but
212 : * some of the rings will not be closed, and rings will have no particular
213 : * order or direction.
214 : */
215 :
216 35 : int SDTSRawPolygon::AssembleRings()
217 :
218 : {
219 35 : if (nRings > 0)
220 0 : return TRUE;
221 :
222 35 : if (nEdges == 0)
223 22 : return FALSE;
224 :
225 : /* -------------------------------------------------------------------- */
226 : /* Setup array of line markers indicating if they have been */
227 : /* added to a ring yet. */
228 : /* -------------------------------------------------------------------- */
229 13 : int nRemainingEdges = nEdges;
230 :
231 : int *panEdgeConsumed =
232 13 : reinterpret_cast<int *>(CPLCalloc(sizeof(int), nEdges));
233 :
234 : /* -------------------------------------------------------------------- */
235 : /* Allocate ring arrays. */
236 : /* -------------------------------------------------------------------- */
237 13 : panRingStart = reinterpret_cast<int *>(CPLMalloc(sizeof(int) * nEdges));
238 :
239 13 : nVertices = 0;
240 63 : for (int iEdge = 0; iEdge < nEdges; iEdge++)
241 : {
242 50 : if (papoEdges[iEdge]->nVertices < 2)
243 : {
244 0 : panEdgeConsumed[iEdge] = TRUE;
245 0 : nRemainingEdges--;
246 : }
247 : else
248 : {
249 50 : nVertices += papoEdges[iEdge]->nVertices;
250 : }
251 : }
252 :
253 13 : padfX = reinterpret_cast<double *>(CPLMalloc(sizeof(double) * nVertices));
254 13 : padfY = reinterpret_cast<double *>(CPLMalloc(sizeof(double) * nVertices));
255 13 : padfZ = reinterpret_cast<double *>(CPLMalloc(sizeof(double) * nVertices));
256 :
257 13 : nVertices = 0;
258 :
259 : /* ==================================================================== */
260 : /* Loop generating rings. */
261 : /* ==================================================================== */
262 13 : bool bSuccess = true;
263 :
264 30 : while (nRemainingEdges > 0)
265 : {
266 : /* --------------------------------------------------------------------
267 : */
268 : /* Find the first unconsumed edge. */
269 : /* --------------------------------------------------------------------
270 : */
271 17 : int iEdge = 0;
272 25 : for (; panEdgeConsumed[iEdge]; iEdge++)
273 : {
274 : }
275 :
276 17 : SDTSRawLine *poEdge = papoEdges[iEdge];
277 :
278 : /* --------------------------------------------------------------------
279 : */
280 : /* Start a new ring, copying in the current line directly */
281 : /* --------------------------------------------------------------------
282 : */
283 17 : panRingStart[nRings++] = nVertices;
284 :
285 17 : AddEdgeToRing(poEdge->nVertices, poEdge->padfX, poEdge->padfY,
286 : poEdge->padfZ, FALSE, FALSE);
287 :
288 17 : panEdgeConsumed[iEdge] = TRUE;
289 17 : nRemainingEdges--;
290 :
291 17 : const int nStartNode = poEdge->oStartNode.nRecord;
292 17 : int nLinkNode = poEdge->oEndNode.nRecord;
293 :
294 : /* ====================================================================
295 : */
296 : /* Loop adding edges to this ring until we make a whole pass */
297 : /* within finding anything to add. */
298 : /* ====================================================================
299 : */
300 17 : bool bWorkDone = true;
301 :
302 38 : while (nLinkNode != nStartNode && nRemainingEdges > 0 && bWorkDone)
303 : {
304 21 : bWorkDone = false;
305 :
306 319 : for (iEdge = 0; iEdge < nEdges; iEdge++)
307 : {
308 298 : if (panEdgeConsumed[iEdge])
309 182 : continue;
310 :
311 116 : poEdge = papoEdges[iEdge];
312 116 : if (poEdge->oStartNode.nRecord == nLinkNode)
313 : {
314 29 : AddEdgeToRing(poEdge->nVertices, poEdge->padfX,
315 : poEdge->padfY, poEdge->padfZ, FALSE, TRUE);
316 29 : nLinkNode = poEdge->oEndNode.nRecord;
317 : }
318 87 : else if (poEdge->oEndNode.nRecord == nLinkNode)
319 : {
320 4 : AddEdgeToRing(poEdge->nVertices, poEdge->padfX,
321 : poEdge->padfY, poEdge->padfZ, TRUE, TRUE);
322 4 : nLinkNode = poEdge->oStartNode.nRecord;
323 : }
324 : else
325 : {
326 83 : continue;
327 : }
328 :
329 33 : panEdgeConsumed[iEdge] = TRUE;
330 33 : nRemainingEdges--;
331 33 : bWorkDone = true;
332 : }
333 : }
334 :
335 : /* --------------------------------------------------------------------
336 : */
337 : /* Did we fail to complete the ring? */
338 : /* --------------------------------------------------------------------
339 : */
340 17 : if (nLinkNode != nStartNode)
341 15 : bSuccess = false;
342 : } /* next ring */
343 :
344 13 : CPLFree(panEdgeConsumed);
345 :
346 13 : if (!bSuccess)
347 11 : return bSuccess;
348 :
349 : /* ==================================================================== */
350 : /* Compute the area of each ring. The sign will be positive */
351 : /* for counter clockwise rings, otherwise negative. */
352 : /* */
353 : /* The algorithm used in this function was taken from _Graphics */
354 : /* Gems II_, James Arvo, 1991, Academic Press, Inc., section 1.1, */
355 : /* "The Area of a Simple Polygon", Jon Rokne, pp. 5-6. */
356 : /* ==================================================================== */
357 2 : double dfMaxArea = 0.0;
358 2 : int iBiggestRing = -1;
359 :
360 : double *padfRingArea =
361 2 : reinterpret_cast<double *>(CPLCalloc(sizeof(double), nRings));
362 :
363 4 : for (int iRing = 0; iRing < nRings; iRing++)
364 : {
365 : int nRingVertices;
366 2 : if (iRing == nRings - 1)
367 2 : nRingVertices = nVertices - panRingStart[iRing];
368 : else
369 0 : nRingVertices = panRingStart[iRing + 1] - panRingStart[iRing];
370 :
371 2 : double dfSum1 = 0.0;
372 2 : double dfSum2 = 0.0;
373 2 : for (int i = panRingStart[iRing];
374 382 : i < panRingStart[iRing] + nRingVertices - 1; i++)
375 : {
376 380 : dfSum1 += padfX[i] * padfY[i + 1];
377 380 : dfSum2 += padfY[i] * padfX[i + 1];
378 : }
379 :
380 2 : padfRingArea[iRing] = (dfSum1 - dfSum2) / 2;
381 :
382 2 : if (std::abs(padfRingArea[iRing]) > dfMaxArea)
383 : {
384 2 : dfMaxArea = std::abs(padfRingArea[iRing]);
385 2 : iBiggestRing = iRing;
386 : }
387 : }
388 :
389 2 : if (iBiggestRing < 0)
390 : {
391 0 : CPLFree(padfRingArea);
392 0 : return FALSE;
393 : }
394 :
395 : /* ==================================================================== */
396 : /* Make a new set of vertices, and copy the largest ring into */
397 : /* it, adjusting the direction if necessary to ensure that this */
398 : /* outer ring is counter clockwise. */
399 : /* ==================================================================== */
400 2 : double *padfXRaw = padfX;
401 2 : double *padfYRaw = padfY;
402 2 : double *padfZRaw = padfZ;
403 2 : int *panRawRingStart = panRingStart;
404 2 : int nRawVertices = nVertices;
405 2 : int nRawRings = nRings;
406 :
407 2 : padfX = reinterpret_cast<double *>(CPLMalloc(sizeof(double) * nVertices));
408 2 : padfY = reinterpret_cast<double *>(CPLMalloc(sizeof(double) * nVertices));
409 2 : padfZ = reinterpret_cast<double *>(CPLMalloc(sizeof(double) * nVertices));
410 2 : panRingStart = reinterpret_cast<int *>(CPLMalloc(sizeof(int) * nRawRings));
411 2 : nVertices = 0;
412 2 : nRings = 0;
413 :
414 : int nRingVertices;
415 2 : if (iBiggestRing == nRawRings - 1)
416 2 : nRingVertices = nRawVertices - panRawRingStart[iBiggestRing];
417 : else
418 0 : nRingVertices =
419 0 : panRawRingStart[iBiggestRing + 1] - panRawRingStart[iBiggestRing];
420 :
421 2 : panRingStart[nRings++] = 0;
422 2 : AddEdgeToRing(nRingVertices, padfXRaw + panRawRingStart[iBiggestRing],
423 2 : padfYRaw + panRawRingStart[iBiggestRing],
424 2 : padfZRaw + panRawRingStart[iBiggestRing],
425 2 : padfRingArea[iBiggestRing] < 0.0, FALSE);
426 :
427 : /* ==================================================================== */
428 : /* Add the rest of the rings, which must be holes, in clockwise */
429 : /* order. */
430 : /* ==================================================================== */
431 4 : for (int iRing = 0; iRing < nRawRings; iRing++)
432 : {
433 2 : if (iRing == iBiggestRing)
434 2 : continue;
435 :
436 0 : if (iRing == nRawRings - 1)
437 0 : nRingVertices = nRawVertices - panRawRingStart[iRing];
438 : else
439 0 : nRingVertices = panRawRingStart[iRing + 1] - panRawRingStart[iRing];
440 :
441 0 : panRingStart[nRings++] = nVertices;
442 0 : AddEdgeToRing(nRingVertices, padfXRaw + panRawRingStart[iRing],
443 0 : padfYRaw + panRawRingStart[iRing],
444 0 : padfZRaw + panRawRingStart[iRing],
445 0 : padfRingArea[iRing] > 0.0, FALSE);
446 : }
447 :
448 : /* -------------------------------------------------------------------- */
449 : /* Cleanup */
450 : /* -------------------------------------------------------------------- */
451 2 : CPLFree(padfXRaw);
452 2 : CPLFree(padfYRaw);
453 2 : CPLFree(padfZRaw);
454 2 : CPLFree(padfRingArea);
455 2 : CPLFree(panRawRingStart);
456 :
457 2 : CPLFree(papoEdges);
458 2 : papoEdges = nullptr;
459 2 : nEdges = 0;
460 :
461 2 : return TRUE;
462 : }
463 :
464 : /************************************************************************/
465 : /* Dump() */
466 : /************************************************************************/
467 :
468 0 : void SDTSRawPolygon::Dump(FILE *fp)
469 :
470 : {
471 0 : fprintf(fp, "SDTSRawPolygon %s: ", oModId.GetName());
472 :
473 0 : for (int i = 0; i < nAttributes; i++)
474 0 : fprintf(fp, " ATID[%d]=%s", i, paoATID[i].GetName());
475 :
476 0 : fprintf(fp, "\n");
477 0 : }
478 :
479 : /************************************************************************/
480 : /* ==================================================================== */
481 : /* SDTSPolygonReader */
482 : /* */
483 : /* This is the class used to read a Polygon module. */
484 : /* ==================================================================== */
485 : /************************************************************************/
486 :
487 : /************************************************************************/
488 : /* SDTSPolygonReader() */
489 : /************************************************************************/
490 :
491 1 : SDTSPolygonReader::SDTSPolygonReader() : bRingsAssembled(FALSE)
492 : {
493 1 : }
494 :
495 : /************************************************************************/
496 : /* ~SDTSPolygonReader() */
497 : /************************************************************************/
498 :
499 2 : SDTSPolygonReader::~SDTSPolygonReader()
500 : {
501 2 : }
502 :
503 : /************************************************************************/
504 : /* Close() */
505 : /************************************************************************/
506 :
507 0 : void SDTSPolygonReader::Close()
508 :
509 : {
510 0 : oDDFModule.Close();
511 0 : }
512 :
513 : /************************************************************************/
514 : /* Open() */
515 : /* */
516 : /* Open the requested line file, and prepare to start reading */
517 : /* data records. */
518 : /************************************************************************/
519 :
520 1 : int SDTSPolygonReader::Open(const char *pszFilename)
521 :
522 : {
523 1 : return oDDFModule.Open(pszFilename);
524 : }
525 :
526 : /************************************************************************/
527 : /* GetNextPolygon() */
528 : /* */
529 : /* Fetch the next feature as an STDSRawPolygon. */
530 : /************************************************************************/
531 :
532 36 : SDTSRawPolygon *SDTSPolygonReader::GetNextPolygon()
533 :
534 : {
535 : /* -------------------------------------------------------------------- */
536 : /* Read a record. */
537 : /* -------------------------------------------------------------------- */
538 36 : if (oDDFModule.GetFP() == nullptr)
539 0 : return nullptr;
540 :
541 36 : DDFRecord *poRecord = oDDFModule.ReadRecord();
542 :
543 36 : if (poRecord == nullptr)
544 1 : return nullptr;
545 :
546 : /* -------------------------------------------------------------------- */
547 : /* Transform into a Polygon feature. */
548 : /* -------------------------------------------------------------------- */
549 35 : SDTSRawPolygon *poRawPolygon = new SDTSRawPolygon();
550 :
551 35 : if (poRawPolygon->Read(poRecord))
552 : {
553 35 : return poRawPolygon;
554 : }
555 :
556 0 : delete poRawPolygon;
557 0 : return nullptr;
558 : }
559 :
560 : /************************************************************************/
561 : /* AssembleRings() */
562 : /************************************************************************/
563 :
564 : /**
565 : * Assemble geometry for a polygon transfer.
566 : *
567 : * This method takes care of attaching lines from all the line layers in
568 : * this transfer to this polygon layer, assembling the lines into rings on
569 : * the polygons, and then cleaning up unnecessary intermediate results.
570 : *
571 : * Currently this method will leave the line layers rewound to the beginning
572 : * but indexed, and the polygon layer rewound but indexed. In the future
573 : * it may restore reading positions, and possibly flush line indexes if they
574 : * were not previously indexed.
575 : *
576 : * This method does nothing if the rings have already been assembled on
577 : * this layer using this method.
578 : *
579 : * See SDTSRawPolygon::AssembleRings() for more information on how the lines
580 : * are assembled into rings.
581 : *
582 : * @param poTransfer the SDTSTransfer that this reader is a part of. Used
583 : * to get a list of line layers that might be needed.
584 : * @param iPolyLayer the polygon reader instance number, used to avoid
585 : * processing lines for other layers.
586 : */
587 :
588 38 : void SDTSPolygonReader::AssembleRings(SDTSTransfer *poTransfer, int iPolyLayer)
589 :
590 : {
591 38 : if (bRingsAssembled)
592 37 : return;
593 :
594 1 : bRingsAssembled = TRUE;
595 :
596 : /* -------------------------------------------------------------------- */
597 : /* To write polygons we need to build them from their related */
598 : /* arcs. We don't know off hand which arc (line) layers */
599 : /* contribute so we process all line layers, attaching them to */
600 : /* polygons as appropriate. */
601 : /* -------------------------------------------------------------------- */
602 9 : for (int iLineLayer = 0; iLineLayer < poTransfer->GetLayerCount();
603 : iLineLayer++)
604 : {
605 8 : if (poTransfer->GetLayerType(iLineLayer) != SLTLine)
606 7 : continue;
607 :
608 : SDTSLineReader *poLineReader = reinterpret_cast<SDTSLineReader *>(
609 1 : poTransfer->GetLayerIndexedReader(iLineLayer));
610 1 : if (poLineReader == nullptr)
611 0 : continue;
612 :
613 1 : poLineReader->AttachToPolygons(poTransfer, iPolyLayer);
614 1 : poLineReader->Rewind();
615 : }
616 :
617 1 : if (!IsIndexed())
618 0 : return;
619 :
620 : /* -------------------------------------------------------------------- */
621 : /* Scan all polygons indexed on this reader, and assemble their */
622 : /* rings. */
623 : /* -------------------------------------------------------------------- */
624 1 : Rewind();
625 :
626 1 : SDTSFeature *poFeature = nullptr;
627 36 : while ((poFeature = GetNextFeature()) != nullptr)
628 : {
629 35 : SDTSRawPolygon *poPoly = reinterpret_cast<SDTSRawPolygon *>(poFeature);
630 :
631 35 : poPoly->AssembleRings();
632 : }
633 :
634 1 : Rewind();
635 : }
|
