Line data Source code
1 : /******************************************************************************
2 : * Project: GDAL
3 : * Purpose: Implements The Two-Arm Chains EdgeTracing Algorithm
4 : * Author: kikitte.lee
5 : *
6 : ******************************************************************************
7 : * Copyright (c) 2023, kikitte.lee <kikitte.lee@gmail.com>
8 : *
9 : * SPDX-License-Identifier: MIT
10 : ****************************************************************************/
11 :
12 : /*! @cond Doxygen_Suppress */
13 :
14 : #include "polygonize_polygonizer.h"
15 :
16 : #include <algorithm>
17 :
18 : namespace gdal
19 : {
20 : namespace polygonizer
21 : {
22 41588 : IndexedArc RPolygon::newArc(bool bFollowRighthand)
23 : {
24 41588 : const std::size_t iArcIndex = oArcs.size();
25 : const auto &oNewArc =
26 41588 : oArcs.emplace_back(static_cast<unsigned>(iArcIndex), bFollowRighthand);
27 41588 : return IndexedArc{oNewArc.poArc.get(), iArcIndex};
28 : }
29 :
30 41588 : void RPolygon::setArcConnection(const IndexedArc &oArc,
31 : const IndexedArc &oNextArc)
32 : {
33 41588 : oArcs[oArc.iIndex].nConnection = static_cast<unsigned>(oNextArc.iIndex);
34 41588 : }
35 :
36 431423 : void RPolygon::updateBottomRightPos(IndexType iRow, IndexType iCol)
37 : {
38 431423 : iBottomRightRow = iRow;
39 431423 : iBottomRightCol = iCol;
40 431423 : }
41 :
42 : /**
43 : * Process different kinds of Arm connections.
44 : */
45 431423 : static void ProcessArmConnections(TwoArm *poCurrent, TwoArm *poAbove,
46 : TwoArm *poLeft)
47 : {
48 431423 : poCurrent->poPolyInside->updateBottomRightPos(poCurrent->iRow,
49 : poCurrent->iCol);
50 431423 : poCurrent->bSolidVertical = poCurrent->poPolyInside != poLeft->poPolyInside;
51 431423 : poCurrent->bSolidHorizontal =
52 431423 : poCurrent->poPolyInside != poAbove->poPolyInside;
53 431423 : poCurrent->poPolyAbove = poAbove->poPolyInside;
54 431423 : poCurrent->poPolyLeft = poLeft->poPolyInside;
55 :
56 431423 : constexpr int BIT_CUR_HORIZ = 0;
57 431423 : constexpr int BIT_CUR_VERT = 1;
58 431423 : constexpr int BIT_LEFT = 2;
59 431423 : constexpr int BIT_ABOVE = 3;
60 :
61 431423 : const int nArmConnectionType =
62 431423 : (static_cast<int>(poAbove->bSolidVertical) << BIT_ABOVE) |
63 431423 : (static_cast<int>(poLeft->bSolidHorizontal) << BIT_LEFT) |
64 431423 : (static_cast<int>(poCurrent->bSolidVertical) << BIT_CUR_VERT) |
65 431423 : (static_cast<int>(poCurrent->bSolidHorizontal) << BIT_CUR_HORIZ);
66 :
67 431423 : constexpr int VIRTUAL = 0;
68 431423 : constexpr int SOLID = 1;
69 :
70 431423 : constexpr int ABOVE_VIRTUAL = VIRTUAL << BIT_ABOVE;
71 431423 : constexpr int ABOVE_SOLID = SOLID << BIT_ABOVE;
72 :
73 431423 : constexpr int LEFT_VIRTUAL = VIRTUAL << BIT_LEFT;
74 431423 : constexpr int LEFT_SOLID = SOLID << BIT_LEFT;
75 :
76 431423 : constexpr int CUR_VERT_VIRTUAL = VIRTUAL << BIT_CUR_VERT;
77 431423 : constexpr int CUR_VERT_SOLID = SOLID << BIT_CUR_VERT;
78 :
79 431423 : constexpr int CUR_HORIZ_VIRTUAL = VIRTUAL << BIT_CUR_HORIZ;
80 431423 : constexpr int CUR_HORIZ_SOLID = SOLID << BIT_CUR_HORIZ;
81 :
82 : /**
83 : * There are 12 valid connection types depending on the arm types(virtual or solid)
84 : * The following diagram illustrates these kinds of connection types, ⇢⇣ means virtual arm, →↓ means solid arm.
85 : * ⇣ ⇣ ⇣ ⇣ ↓
86 : * ⇢ → → → → ⇢ → → ⇢ →
87 : * ↓ ⇣ ↓ ↓ ⇣
88 : * type=3 type=5 type=6 type=7 type=9
89 : *
90 : * ↓ ↓ ↓ ↓ ↓
91 : * ⇢ ⇢ ⇢ → → ⇢ → → → ⇢
92 : * ↓ ↓ ⇣ ⇣ ↓
93 : * type=10 type=11 type=12 type=13 type=14
94 : *
95 : * ↓ ⇣
96 : * → → ⇢ ⇢
97 : * ↓ ⇣
98 : * type=15 type=0
99 : *
100 : * For each connection type, we may create new arc, ,
101 : * Depending on the connection type, we may do the following things:
102 : * 1. Create new arc. If the arc is closed to the inner polygon, it is called "Inner Arc", otherwise "Outer Arc"
103 : * 2. Pass an arc to the next arm.
104 : * 3. "Close" two arcs. If two arcs meet at the bottom right corner of a cell, close them by recording the arc connection.
105 : * 4. Add grid position(row, col) to an arc.
106 : */
107 :
108 431423 : switch (nArmConnectionType)
109 : {
110 361885 : case ABOVE_VIRTUAL | LEFT_VIRTUAL | CUR_VERT_VIRTUAL |
111 : CUR_HORIZ_VIRTUAL: // 0
112 : // nothing to do
113 361885 : break;
114 :
115 7294 : case ABOVE_VIRTUAL | LEFT_VIRTUAL | CUR_VERT_SOLID |
116 : CUR_HORIZ_SOLID: // 3
117 : // add inner arcs
118 7294 : poCurrent->oArcVerInner = poCurrent->poPolyInside->newArc(true);
119 7294 : poCurrent->oArcHorInner = poCurrent->poPolyInside->newArc(false);
120 7294 : poCurrent->poPolyInside->setArcConnection(poCurrent->oArcHorInner,
121 7294 : poCurrent->oArcVerInner);
122 7294 : poCurrent->oArcVerInner.poArc->emplace_back(
123 7294 : Point{poCurrent->iRow, poCurrent->iCol});
124 :
125 : // add outer arcs
126 7294 : poCurrent->oArcHorOuter = poAbove->poPolyInside->newArc(true);
127 7294 : poCurrent->oArcVerOuter = poAbove->poPolyInside->newArc(false);
128 7294 : poAbove->poPolyInside->setArcConnection(poCurrent->oArcVerOuter,
129 7294 : poCurrent->oArcHorOuter);
130 7294 : poCurrent->oArcHorOuter.poArc->push_back(
131 7294 : Point{poCurrent->iRow, poCurrent->iCol});
132 :
133 7294 : break;
134 16372 : case ABOVE_VIRTUAL | LEFT_SOLID | CUR_VERT_VIRTUAL |
135 : CUR_HORIZ_SOLID: // 5
136 : // pass arcs
137 16372 : poCurrent->oArcHorInner = poLeft->oArcHorInner;
138 16372 : poCurrent->oArcHorOuter = poLeft->oArcHorOuter;
139 :
140 16372 : break;
141 8805 : case ABOVE_VIRTUAL | LEFT_SOLID | CUR_VERT_SOLID |
142 : CUR_HORIZ_VIRTUAL: // 6
143 : // pass arcs
144 8805 : poCurrent->oArcVerInner = poLeft->oArcHorOuter;
145 8805 : poCurrent->oArcVerOuter = poLeft->oArcHorInner;
146 8805 : poCurrent->oArcVerInner.poArc->push_back(
147 8805 : Point{poCurrent->iRow, poCurrent->iCol});
148 8805 : poCurrent->oArcVerOuter.poArc->push_back(
149 8805 : Point{poCurrent->iRow, poCurrent->iCol});
150 :
151 8805 : break;
152 1328 : case ABOVE_VIRTUAL | LEFT_SOLID | CUR_VERT_SOLID |
153 : CUR_HORIZ_SOLID: // 7
154 : // pass arcs
155 1328 : poCurrent->oArcHorOuter = poLeft->oArcHorOuter;
156 1328 : poCurrent->oArcVerOuter = poLeft->oArcHorInner;
157 1328 : poLeft->oArcHorInner.poArc->push_back(
158 1328 : Point{poCurrent->iRow, poCurrent->iCol});
159 :
160 : // add inner arcs
161 1328 : poCurrent->oArcVerInner = poCurrent->poPolyInside->newArc(true);
162 1328 : poCurrent->oArcHorInner = poCurrent->poPolyInside->newArc(false);
163 1328 : poCurrent->poPolyInside->setArcConnection(poCurrent->oArcHorInner,
164 1328 : poCurrent->oArcVerInner);
165 1328 : poCurrent->oArcVerInner.poArc->push_back(
166 1328 : Point{poCurrent->iRow, poCurrent->iCol});
167 :
168 1328 : break;
169 8769 : case ABOVE_SOLID | LEFT_VIRTUAL | CUR_VERT_VIRTUAL |
170 : CUR_HORIZ_SOLID: // 9
171 : // pass arcs
172 8769 : poCurrent->oArcHorOuter = poAbove->oArcVerInner;
173 8769 : poCurrent->oArcHorInner = poAbove->oArcVerOuter;
174 8769 : poCurrent->oArcHorOuter.poArc->push_back(
175 8769 : Point{poCurrent->iRow, poCurrent->iCol});
176 8769 : poCurrent->oArcHorInner.poArc->push_back(
177 8769 : Point{poCurrent->iRow, poCurrent->iCol});
178 :
179 8769 : break;
180 12056 : case ABOVE_SOLID | LEFT_VIRTUAL | CUR_VERT_SOLID |
181 : CUR_HORIZ_VIRTUAL: // 10
182 : // pass arcs
183 12056 : poCurrent->oArcVerInner = poAbove->oArcVerInner;
184 12056 : poCurrent->oArcVerOuter = poAbove->oArcVerOuter;
185 :
186 12056 : break;
187 1485 : case ABOVE_SOLID | LEFT_VIRTUAL | CUR_VERT_SOLID |
188 : CUR_HORIZ_SOLID: // 11
189 : // pass arcs
190 1485 : poCurrent->oArcHorOuter = poAbove->oArcVerInner;
191 1485 : poCurrent->oArcVerOuter = poAbove->oArcVerOuter;
192 1485 : poCurrent->oArcHorOuter.poArc->push_back(
193 1485 : Point{poCurrent->iRow, poCurrent->iCol});
194 : // add inner arcs
195 1485 : poCurrent->oArcVerInner = poCurrent->poPolyInside->newArc(true);
196 1485 : poCurrent->oArcHorInner = poCurrent->poPolyInside->newArc(false);
197 1485 : poCurrent->poPolyInside->setArcConnection(poCurrent->oArcHorInner,
198 1485 : poCurrent->oArcVerInner);
199 1485 : poCurrent->oArcVerInner.poArc->push_back(
200 1485 : Point{poCurrent->iRow, poCurrent->iCol});
201 :
202 1485 : break;
203 7323 : case ABOVE_SOLID | LEFT_SOLID | CUR_VERT_VIRTUAL |
204 : CUR_HORIZ_VIRTUAL: // 12
205 : // close arcs
206 7323 : poLeft->oArcHorOuter.poArc->push_back(
207 7323 : Point{poCurrent->iRow, poCurrent->iCol});
208 7323 : poLeft->poPolyAbove->setArcConnection(poLeft->oArcHorOuter,
209 7323 : poAbove->oArcVerOuter);
210 : // close arcs
211 7323 : poAbove->oArcVerInner.poArc->push_back(
212 7323 : Point{poCurrent->iRow, poCurrent->iCol});
213 7323 : poCurrent->poPolyInside->setArcConnection(poAbove->oArcVerInner,
214 7323 : poLeft->oArcHorInner);
215 :
216 7323 : break;
217 1335 : case ABOVE_SOLID | LEFT_SOLID | CUR_VERT_VIRTUAL |
218 : CUR_HORIZ_SOLID: // 13
219 : // close arcs
220 1335 : poLeft->oArcHorOuter.poArc->push_back(
221 1335 : Point{poCurrent->iRow, poCurrent->iCol});
222 1335 : poLeft->poPolyAbove->setArcConnection(poLeft->oArcHorOuter,
223 1335 : poAbove->oArcVerOuter);
224 : // pass arcs
225 1335 : poCurrent->oArcHorOuter = poAbove->oArcVerInner;
226 1335 : poCurrent->oArcHorInner = poLeft->oArcHorInner;
227 1335 : poCurrent->oArcHorOuter.poArc->push_back(
228 1335 : Point{poCurrent->iRow, poCurrent->iCol});
229 :
230 1335 : break;
231 1420 : case ABOVE_SOLID | LEFT_SOLID | CUR_VERT_SOLID |
232 : CUR_HORIZ_VIRTUAL: // 14
233 : // close arcs
234 1420 : poLeft->oArcHorOuter.poArc->push_back(
235 1420 : Point{poCurrent->iRow, poCurrent->iCol});
236 1420 : poLeft->poPolyAbove->setArcConnection(poLeft->oArcHorOuter,
237 1420 : poAbove->oArcVerOuter);
238 : // pass arcs
239 1420 : poCurrent->oArcVerInner = poAbove->oArcVerInner;
240 1420 : poCurrent->oArcVerOuter = poLeft->oArcHorInner;
241 1420 : poCurrent->oArcVerOuter.poArc->push_back(
242 1420 : Point{poCurrent->iRow, poCurrent->iCol});
243 :
244 1420 : break;
245 3351 : case ABOVE_SOLID | LEFT_SOLID | CUR_VERT_SOLID | CUR_HORIZ_SOLID: // 15
246 : // Tow pixels of the main diagonal belong to the same polygon
247 3351 : if (poAbove->poPolyLeft == poCurrent->poPolyInside)
248 : {
249 : // pass arcs
250 62 : poCurrent->oArcVerInner = poLeft->oArcHorOuter;
251 62 : poCurrent->oArcHorInner = poAbove->oArcVerOuter;
252 62 : poCurrent->oArcVerInner.poArc->push_back(
253 62 : Point{poCurrent->iRow, poCurrent->iCol});
254 62 : poCurrent->oArcHorInner.poArc->push_back(
255 62 : Point{poCurrent->iRow, poCurrent->iCol});
256 : }
257 : else
258 : {
259 : // close arcs
260 3289 : poLeft->oArcHorOuter.poArc->push_back(
261 3289 : Point{poCurrent->iRow, poCurrent->iCol});
262 3289 : poLeft->poPolyAbove->setArcConnection(poLeft->oArcHorOuter,
263 3289 : poAbove->oArcVerOuter);
264 : // add inner arcs
265 3289 : poCurrent->oArcVerInner = poCurrent->poPolyInside->newArc(true);
266 : poCurrent->oArcHorInner =
267 3289 : poCurrent->poPolyInside->newArc(false);
268 3289 : poCurrent->poPolyInside->setArcConnection(
269 3289 : poCurrent->oArcHorInner, poCurrent->oArcVerInner);
270 3289 : poCurrent->oArcVerInner.poArc->push_back(
271 3289 : Point{poCurrent->iRow, poCurrent->iCol});
272 : }
273 :
274 : // Tow pixels of the secondary diagonal belong to the same polygon
275 3351 : if (poAbove->poPolyInside == poLeft->poPolyInside)
276 : {
277 : // close arcs
278 104 : poAbove->poPolyInside->setArcConnection(poAbove->oArcVerInner,
279 104 : poLeft->oArcHorInner);
280 104 : poAbove->oArcVerInner.poArc->push_back(
281 104 : Point{poCurrent->iRow, poCurrent->iCol});
282 : // add outer arcs
283 104 : poCurrent->oArcHorOuter = poAbove->poPolyInside->newArc(true);
284 104 : poCurrent->oArcVerOuter = poAbove->poPolyInside->newArc(false);
285 104 : poCurrent->oArcHorOuter.poArc->push_back(
286 104 : Point{poCurrent->iRow, poCurrent->iCol});
287 104 : poAbove->poPolyInside->setArcConnection(
288 104 : poCurrent->oArcVerOuter, poCurrent->oArcHorOuter);
289 : }
290 : else
291 : {
292 : // pass arcs
293 3247 : poCurrent->oArcHorOuter = poAbove->oArcVerInner;
294 3247 : poCurrent->oArcVerOuter = poLeft->oArcHorInner;
295 3247 : poCurrent->oArcHorOuter.poArc->push_back(
296 3247 : Point{poCurrent->iRow, poCurrent->iCol});
297 3247 : poCurrent->oArcVerOuter.poArc->push_back(
298 3247 : Point{poCurrent->iRow, poCurrent->iCol});
299 : }
300 :
301 3351 : break;
302 :
303 0 : case ABOVE_VIRTUAL | LEFT_VIRTUAL | CUR_VERT_VIRTUAL |
304 : CUR_HORIZ_SOLID: // 1
305 : case ABOVE_VIRTUAL | LEFT_VIRTUAL | CUR_VERT_SOLID |
306 : CUR_HORIZ_VIRTUAL: // 2
307 : case ABOVE_VIRTUAL | LEFT_SOLID | CUR_VERT_VIRTUAL |
308 : CUR_HORIZ_VIRTUAL: // 4
309 : default:
310 : // Impossible case
311 0 : CPLAssert(false);
312 : break;
313 : }
314 431423 : }
315 :
316 : template <typename PolyIdType, typename DataType>
317 90 : Polygonizer<PolyIdType, DataType>::Polygonizer(
318 : PolyIdType nInvalidPolyId, PolygonReceiver<DataType> *poPolygonReceiver)
319 90 : : nInvalidPolyId_(nInvalidPolyId), poPolygonReceiver_(poPolygonReceiver)
320 : {
321 90 : poTheOuterPolygon_ = createPolygon(THE_OUTER_POLYGON_ID);
322 90 : }
323 :
324 : template <typename PolyIdType, typename DataType>
325 90 : Polygonizer<PolyIdType, DataType>::~Polygonizer()
326 : {
327 : // cppcheck-suppress constVariableReference
328 180 : for (auto &pair : oPolygonMap_)
329 : {
330 90 : delete pair.second;
331 : }
332 90 : }
333 :
334 : template <typename PolyIdType, typename DataType>
335 429503 : RPolygon *Polygonizer<PolyIdType, DataType>::getPolygon(PolyIdType nPolygonId)
336 : {
337 429503 : const auto oIter = oPolygonMap_.find(nPolygonId);
338 429503 : if (oIter == oPolygonMap_.end())
339 : {
340 6104 : return createPolygon(nPolygonId);
341 : }
342 : else
343 : {
344 423399 : return oIter->second;
345 : }
346 : }
347 :
348 : template <typename PolyIdType, typename DataType>
349 : RPolygon *
350 6194 : Polygonizer<PolyIdType, DataType>::createPolygon(PolyIdType nPolygonId)
351 : {
352 6194 : auto polygon = new RPolygon();
353 6194 : oPolygonMap_[nPolygonId] = polygon;
354 6194 : return polygon;
355 : }
356 :
357 : template <typename PolyIdType, typename DataType>
358 6104 : void Polygonizer<PolyIdType, DataType>::destroyPolygon(PolyIdType nPolygonId)
359 : {
360 6104 : const auto oIter = oPolygonMap_.find(nPolygonId);
361 6104 : CPLAssert(oIter != oPolygonMap_.end());
362 6104 : delete oIter->second;
363 6104 : oPolygonMap_.erase(oIter);
364 6104 : }
365 :
366 : template <typename PolyIdType, typename DataType>
367 1920 : bool Polygonizer<PolyIdType, DataType>::processLine(
368 : const PolyIdType *panThisLineId, const DataType *panLastLineVal,
369 : TwoArm *poThisLineArm, TwoArm *poLastLineArm, const IndexType nCurrentRow,
370 : const IndexType nCols)
371 : {
372 : TwoArm *poCurrent, *poAbove, *poLeft;
373 :
374 : try
375 : {
376 1920 : poCurrent = poThisLineArm + 1;
377 1920 : poCurrent->iRow = nCurrentRow;
378 1920 : poCurrent->iCol = 0;
379 1920 : poCurrent->poPolyInside = getPolygon(panThisLineId[0]);
380 1920 : poAbove = poLastLineArm + 1;
381 1920 : poLeft = poThisLineArm;
382 1920 : poLeft->poPolyInside = poTheOuterPolygon_;
383 1920 : ProcessArmConnections(poCurrent, poAbove, poLeft);
384 429503 : for (IndexType col = 1; col < nCols; ++col)
385 : {
386 427583 : IndexType iArmIndex = col + 1;
387 427583 : poCurrent = poThisLineArm + iArmIndex;
388 427583 : poCurrent->iRow = nCurrentRow;
389 427583 : poCurrent->iCol = col;
390 427583 : poCurrent->poPolyInside = getPolygon(panThisLineId[col]);
391 427583 : poAbove = poLastLineArm + iArmIndex;
392 427583 : poLeft = poThisLineArm + iArmIndex - 1;
393 427583 : ProcessArmConnections(poCurrent, poAbove, poLeft);
394 : }
395 1920 : poCurrent = poThisLineArm + nCols + 1;
396 1920 : poCurrent->iRow = nCurrentRow;
397 1920 : poCurrent->iCol = nCols;
398 1920 : poCurrent->poPolyInside = poTheOuterPolygon_;
399 1920 : poAbove = poLastLineArm + nCols + 1;
400 1920 : poAbove->poPolyInside = poTheOuterPolygon_;
401 1920 : poLeft = poThisLineArm + nCols;
402 1920 : ProcessArmConnections(poCurrent, poAbove, poLeft);
403 :
404 : /**
405 : *
406 : * Find those polygons haven't been processed on this line as we can be sure they are completed
407 : *
408 : */
409 1920 : std::vector<PolygonMapEntry> oCompletedPolygons;
410 40762 : for (auto &entry : oPolygonMap_)
411 : {
412 38842 : RPolygon *poPolygon = entry.second;
413 :
414 38842 : if (poPolygon->iBottomRightRow + 1 == nCurrentRow)
415 : {
416 6104 : oCompletedPolygons.push_back(entry);
417 : }
418 : }
419 : // cppcheck-suppress constVariableReference
420 8024 : for (auto &entry : oCompletedPolygons)
421 : {
422 6104 : PolyIdType nPolyId = entry.first;
423 6104 : RPolygon *poPolygon = entry.second;
424 :
425 : // emit valid polygon only
426 6104 : if (nPolyId != nInvalidPolyId_)
427 : {
428 6075 : poPolygonReceiver_->receive(
429 6075 : poPolygon, panLastLineVal[poPolygon->iBottomRightCol]);
430 : }
431 :
432 6104 : destroyPolygon(nPolyId);
433 : }
434 1920 : return true;
435 : }
436 0 : catch (const std::bad_alloc &)
437 : {
438 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
439 : "Out of memory in Polygonizer::processLine");
440 0 : return false;
441 : }
442 : }
443 :
444 : template <typename DataType>
445 90 : OGRPolygonWriter<DataType>::OGRPolygonWriter(OGRLayerH hOutLayer,
446 : int iPixValField,
447 : double *padfGeoTransform)
448 90 : : PolygonReceiver<DataType>(), poOutLayer_(OGRLayer::FromHandle(hOutLayer)),
449 90 : iPixValField_(iPixValField), padfGeoTransform_(padfGeoTransform)
450 : {
451 90 : poFeature_ = std::make_unique<OGRFeature>(poOutLayer_->GetLayerDefn());
452 90 : poPolygon_ = new OGRPolygon();
453 90 : poFeature_->SetGeometryDirectly(poPolygon_);
454 90 : }
455 :
456 : template <typename DataType>
457 6075 : void OGRPolygonWriter<DataType>::receive(RPolygon *poPolygon,
458 : DataType nPolygonCellValue)
459 : {
460 6075 : std::vector<bool> oAccessedArc(poPolygon->oArcs.size(), false);
461 6075 : double *padfGeoTransform = padfGeoTransform_;
462 :
463 6075 : OGRLinearRing *poFirstRing = poPolygon_->getExteriorRing();
464 6075 : if (poFirstRing && poPolygon_->getNumInteriorRings() == 0)
465 : {
466 5974 : poFirstRing->empty();
467 : }
468 : else
469 : {
470 101 : poFirstRing = nullptr;
471 101 : poPolygon_->empty();
472 : }
473 :
474 6075 : auto AddRingToPolygon =
475 90303 : [this, &poPolygon, &oAccessedArc,
476 : padfGeoTransform](std::size_t iFirstArcIndex, OGRLinearRing *poRing)
477 : {
478 6099 : std::unique_ptr<OGRLinearRing> poNewRing;
479 6099 : if (!poRing)
480 : {
481 125 : poNewRing = std::make_unique<OGRLinearRing>();
482 125 : poRing = poNewRing.get();
483 : }
484 :
485 6099 : auto AddArcToRing =
486 575030 : [&poPolygon, poRing, padfGeoTransform](std::size_t iArcIndex)
487 : {
488 38990 : const auto &oArc = poPolygon->oArcs[iArcIndex];
489 38990 : const bool bArcFollowRighthand = oArc.bFollowRighthand;
490 38990 : const int nArcPointCount = static_cast<int>(oArc.poArc->size());
491 38990 : int nDstPointIdx = poRing->getNumPoints();
492 38990 : poRing->setNumPoints(nDstPointIdx + nArcPointCount,
493 : /* bZeroizeNewContent = */ false);
494 38990 : if (poRing->getNumPoints() < nDstPointIdx + nArcPointCount)
495 : {
496 0 : return false;
497 : }
498 122804 : for (int i = 0; i < nArcPointCount; ++i)
499 : {
500 106099 : const Point &oPixel =
501 83814 : (*oArc.poArc)[bArcFollowRighthand
502 : ? i
503 22285 : : (nArcPointCount - i - 1)];
504 :
505 251442 : const double dfX = padfGeoTransform[0] +
506 83814 : oPixel[1] * padfGeoTransform[1] +
507 83814 : oPixel[0] * padfGeoTransform[2];
508 251442 : const double dfY = padfGeoTransform[3] +
509 83814 : oPixel[1] * padfGeoTransform[4] +
510 83814 : oPixel[0] * padfGeoTransform[5];
511 :
512 83814 : poRing->setPoint(nDstPointIdx, dfX, dfY);
513 83814 : ++nDstPointIdx;
514 : }
515 38990 : return true;
516 : };
517 :
518 6099 : if (!AddArcToRing(iFirstArcIndex))
519 : {
520 0 : return false;
521 : }
522 :
523 6099 : std::size_t iArcIndex = iFirstArcIndex;
524 6099 : std::size_t iNextArcIndex = poPolygon->oArcs[iArcIndex].nConnection;
525 6099 : oAccessedArc[iArcIndex] = true;
526 38990 : while (iNextArcIndex != iFirstArcIndex)
527 : {
528 32891 : if (!AddArcToRing(iNextArcIndex))
529 : {
530 0 : return false;
531 : }
532 32891 : iArcIndex = iNextArcIndex;
533 32891 : iNextArcIndex = poPolygon->oArcs[iArcIndex].nConnection;
534 32891 : oAccessedArc[iArcIndex] = true;
535 : }
536 :
537 : // close ring manually
538 6099 : poRing->closeRings();
539 :
540 6099 : if (poNewRing)
541 125 : poPolygon_->addRingDirectly(poNewRing.release());
542 6099 : return true;
543 : };
544 :
545 45065 : for (size_t i = 0; i < oAccessedArc.size(); ++i)
546 : {
547 38990 : if (!oAccessedArc[i])
548 : {
549 6099 : if (!AddRingToPolygon(i, poFirstRing))
550 : {
551 0 : eErr_ = CE_Failure;
552 0 : return;
553 : }
554 6099 : poFirstRing = nullptr;
555 : }
556 : }
557 :
558 : // Create the feature object
559 6075 : poFeature_->SetFID(OGRNullFID);
560 6075 : if (iPixValField_ >= 0)
561 5452 : poFeature_->SetField(iPixValField_,
562 : static_cast<double>(nPolygonCellValue));
563 :
564 : // Write the to the layer.
565 6075 : if (poOutLayer_->CreateFeature(poFeature_.get()) != OGRERR_NONE)
566 0 : eErr_ = CE_Failure;
567 :
568 : // Shouldn't happen for well behaved drivers, but better check...
569 6075 : else if (poFeature_->GetGeometryRef() != poPolygon_)
570 : {
571 0 : poPolygon_ = new OGRPolygon();
572 0 : poFeature_->SetGeometryDirectly(poPolygon_);
573 : }
574 : }
575 :
576 : } // namespace polygonizer
577 : } // namespace gdal
578 :
579 : /*! @endcond */
|
