Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: GDAL
4 : * Purpose: Vector polygon rasterization code.
5 : * Author: Frank Warmerdam, warmerdam@pobox.com
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2000, Frank Warmerdam <warmerdam@pobox.com>
9 : * Copyright (c) 2011, Even Rouault <even dot rouault at spatialys.com>
10 : *
11 : * SPDX-License-Identifier: MIT
12 : ****************************************************************************/
13 :
14 : #include "cpl_port.h"
15 : #include "gdal_alg_priv.h"
16 :
17 : #include <cmath>
18 : #include <cstdlib>
19 : #include <cstring>
20 :
21 : #include <algorithm>
22 : #include <set>
23 : #include <utility>
24 : #include <vector>
25 :
26 : #include "gdal_alg.h"
27 :
28 : /************************************************************************/
29 : /* dllImageFilledPolygon() */
30 : /* */
31 : /* Perform scanline conversion of the passed multi-ring */
32 : /* polygon. Note the polygon does not need to be explicitly */
33 : /* closed. The scanline function will be called with */
34 : /* horizontal scanline chunks which may not be entirely */
35 : /* contained within the valid raster area (in the X */
36 : /* direction). */
37 : /* */
38 : /* NEW: Nodes' coordinate are kept as double in order */
39 : /* to compute accurately the intersections with the lines */
40 : /* */
41 : /* A pixel is considered inside a polygon if its center */
42 : /* falls inside the polygon. This is robust unless */
43 : /* the nodes are placed in the center of the pixels in which */
44 : /* case, due to numerical inaccuracies, it's hard to predict */
45 : /* if the pixel will be considered inside or outside the shape. */
46 : /************************************************************************/
47 :
48 : /*
49 : * NOTE: This code was originally adapted from the gdImageFilledPolygon()
50 : * function in libgd.
51 : *
52 : * http://www.boutell.com/gd/
53 : *
54 : * It was later adapted for direct inclusion in GDAL and relicensed under
55 : * the GDAL MIT license (pulled from the OpenEV distribution).
56 : */
57 :
58 150 : void GDALdllImageFilledPolygon(int nRasterXSize, int nRasterYSize,
59 : int nPartCount, const int *panPartSize,
60 : const double *padfX, const double *padfY,
61 : const double *dfVariant,
62 : llScanlineFunc pfnScanlineFunc,
63 : GDALRasterizeInfo *pCBData,
64 : bool bAvoidBurningSamePoints)
65 : {
66 150 : if (!nPartCount)
67 : {
68 0 : return;
69 : }
70 :
71 150 : int n = 0;
72 311 : for (int part = 0; part < nPartCount; part++)
73 161 : n += panPartSize[part];
74 :
75 300 : std::vector<int> polyInts(n);
76 300 : std::vector<int> polyInts2;
77 150 : if (bAvoidBurningSamePoints)
78 15 : polyInts2.resize(n);
79 :
80 150 : double dminy = padfY[0];
81 150 : double dmaxy = padfY[0];
82 6886 : for (int i = 1; i < n; i++)
83 : {
84 6736 : if (padfY[i] < dminy)
85 : {
86 952 : dminy = padfY[i];
87 : }
88 6736 : if (padfY[i] > dmaxy)
89 : {
90 1322 : dmaxy = padfY[i];
91 : }
92 : }
93 150 : int miny = static_cast<int>(dminy);
94 150 : int maxy = static_cast<int>(dmaxy);
95 :
96 150 : if (miny < 0)
97 10 : miny = 0;
98 150 : if (maxy >= nRasterYSize)
99 16 : maxy = nRasterYSize - 1;
100 :
101 150 : int minx = 0;
102 150 : const int maxx = nRasterXSize - 1;
103 :
104 : // Fix in 1.3: count a vertex only once.
105 9465 : for (int y = miny; y <= maxy; y++)
106 : {
107 9315 : int partoffset = 0;
108 :
109 9315 : const double dy = y + 0.5; // Center height of line.
110 :
111 9315 : int part = 0;
112 9315 : int ints = 0;
113 9315 : int ints2 = 0;
114 :
115 6965280 : for (int i = 0; i < n; i++)
116 : {
117 6955970 : if (i == partoffset + panPartSize[part])
118 : {
119 150 : partoffset += panPartSize[part];
120 150 : part++;
121 : }
122 :
123 6955970 : int ind1 = 0;
124 6955970 : int ind2 = 0;
125 6955970 : if (i == partoffset)
126 : {
127 9465 : ind1 = partoffset + panPartSize[part] - 1;
128 9465 : ind2 = partoffset;
129 : }
130 : else
131 : {
132 6946500 : ind1 = i - 1;
133 6946500 : ind2 = i;
134 : }
135 :
136 6955970 : double dy1 = padfY[ind1];
137 6955970 : double dy2 = padfY[ind2];
138 :
139 6955970 : if ((dy1 < dy && dy2 < dy) || (dy1 > dy && dy2 > dy))
140 6937340 : continue;
141 :
142 18750 : double dx1 = 0.0;
143 18750 : double dx2 = 0.0;
144 18750 : if (dy1 < dy2)
145 : {
146 9317 : dx1 = padfX[ind1];
147 9317 : dx2 = padfX[ind2];
148 : }
149 9433 : else if (dy1 > dy2)
150 : {
151 9317 : std::swap(dy1, dy2);
152 9317 : dx2 = padfX[ind1];
153 9317 : dx1 = padfX[ind2];
154 : }
155 : else // if( fabs(dy1-dy2) < 1.e-6 )
156 : {
157 :
158 : // AE: DO NOT skip bottom horizontal segments
159 : // -Fill them separately-
160 116 : if (padfX[ind1] > padfX[ind2])
161 : {
162 36 : const int horizontal_x1 =
163 36 : static_cast<int>(floor(padfX[ind2] + 0.5));
164 36 : const int horizontal_x2 =
165 36 : static_cast<int>(floor(padfX[ind1] + 0.5));
166 :
167 36 : if ((horizontal_x1 > maxx) || (horizontal_x2 <= minx))
168 0 : continue;
169 :
170 : // Fill the horizontal segment (separately from the rest).
171 36 : if (bAvoidBurningSamePoints)
172 : {
173 10 : polyInts2[ints2++] = horizontal_x1;
174 10 : polyInts2[ints2++] = horizontal_x2;
175 : }
176 : else
177 : {
178 26 : pfnScanlineFunc(
179 : pCBData, y, horizontal_x1, horizontal_x2 - 1,
180 : dfVariant == nullptr ? 0 : dfVariant[0]);
181 : }
182 : }
183 : // else: Skip top horizontal segments.
184 : // They are already filled in the regular loop.
185 116 : continue;
186 : }
187 :
188 18634 : if (dy < dy2 && dy >= dy1)
189 : {
190 18554 : const double intersect =
191 18554 : (dy - dy1) * (dx2 - dx1) / (dy2 - dy1) + dx1;
192 :
193 18554 : polyInts[ints++] = static_cast<int>(floor(intersect + 0.5));
194 : }
195 : }
196 :
197 9315 : std::sort(polyInts.begin(), polyInts.begin() + ints);
198 9315 : std::sort(polyInts2.begin(), polyInts2.begin() + ints2);
199 :
200 18592 : for (int i = 0; i + 1 < ints; i += 2)
201 : {
202 9277 : if (polyInts[i] <= maxx && polyInts[i + 1] > minx)
203 : {
204 8609 : pfnScanlineFunc(pCBData, y, polyInts[i], polyInts[i + 1] - 1,
205 : dfVariant == nullptr ? 0 : dfVariant[0]);
206 : }
207 : }
208 :
209 9325 : for (int i2 = 0, i = 0; i2 + 1 < ints2; i2 += 2)
210 : {
211 10 : if (polyInts2[i2] <= maxx && polyInts2[i2 + 1] > minx)
212 : {
213 : // "synchronize" polyInts[i] with polyInts2[i2]
214 10 : while (i + 1 < ints && polyInts[i] < polyInts2[i2])
215 0 : i += 2;
216 : // Only burn if we don't have a common segment between
217 : // polyInts[] and polyInts2[]
218 10 : if (i + 1 >= ints || polyInts[i] != polyInts2[i2])
219 : {
220 4 : pfnScanlineFunc(pCBData, y, polyInts2[i2],
221 2 : polyInts2[i2 + 1] - 1,
222 : dfVariant == nullptr ? 0 : dfVariant[0]);
223 : }
224 : }
225 : }
226 : }
227 : }
228 :
229 : /************************************************************************/
230 : /* GDALdllImagePoint() */
231 : /************************************************************************/
232 :
233 5 : void GDALdllImagePoint(int nRasterXSize, int nRasterYSize, int nPartCount,
234 : const int * /*panPartSize*/, const double *padfX,
235 : const double *padfY, const double *padfVariant,
236 : llPointFunc pfnPointFunc, GDALRasterizeInfo *pCBData)
237 : {
238 10 : for (int i = 0; i < nPartCount; i++)
239 : {
240 5 : const int nX = static_cast<int>(floor(padfX[i]));
241 5 : const int nY = static_cast<int>(floor(padfY[i]));
242 5 : double dfVariant = 0.0;
243 5 : if (padfVariant != nullptr)
244 0 : dfVariant = padfVariant[i];
245 :
246 5 : if (0 <= nX && nX < nRasterXSize && 0 <= nY && nY < nRasterYSize)
247 5 : pfnPointFunc(pCBData, nY, nX, dfVariant);
248 : }
249 5 : }
250 :
251 : /************************************************************************/
252 : /* GDALdllImageLine() */
253 : /************************************************************************/
254 :
255 1874 : void GDALdllImageLine(int nRasterXSize, int nRasterYSize, int nPartCount,
256 : const int *panPartSize, const double *padfX,
257 : const double *padfY, const double *padfVariant,
258 : llPointFunc pfnPointFunc, GDALRasterizeInfo *pCBData)
259 : {
260 1874 : if (!nPartCount)
261 0 : return;
262 :
263 3748 : for (int i = 0, n = 0; i < nPartCount; n += panPartSize[i++])
264 : {
265 48667 : for (int j = 1; j < panPartSize[i]; j++)
266 : {
267 46793 : int iX = static_cast<int>(floor(padfX[n + j - 1]));
268 46793 : int iY = static_cast<int>(floor(padfY[n + j - 1]));
269 :
270 46793 : const int iX1 = static_cast<int>(floor(padfX[n + j]));
271 46793 : const int iY1 = static_cast<int>(floor(padfY[n + j]));
272 :
273 46793 : double dfVariant = 0.0;
274 46793 : double dfVariant1 = 0.0;
275 46793 : if (padfVariant != nullptr &&
276 46758 : pCBData->eBurnValueSource != GBV_UserBurnValue)
277 : {
278 46758 : dfVariant = padfVariant[n + j - 1];
279 46758 : dfVariant1 = padfVariant[n + j];
280 : }
281 :
282 46793 : int nDeltaX = std::abs(iX1 - iX);
283 46793 : int nDeltaY = std::abs(iY1 - iY);
284 :
285 : // Step direction depends on line direction.
286 46793 : const int nXStep = (iX > iX1) ? -1 : 1;
287 46793 : const int nYStep = (iY > iY1) ? -1 : 1;
288 :
289 : // Determine the line slope.
290 46793 : if (nDeltaX >= nDeltaY)
291 : {
292 38876 : const int nXError = nDeltaY << 1;
293 38876 : const int nYError = nXError - (nDeltaX << 1);
294 38876 : int nError = nXError - nDeltaX;
295 : // == 0 makes clang -fcatch-undefined-behavior -ftrapv happy,
296 : // but if it is == 0, dfDeltaVariant is not really used, so any
297 : // value is okay.
298 38876 : const double dfDeltaVariant =
299 38876 : nDeltaX == 0 ? 0.0 : (dfVariant1 - dfVariant) / nDeltaX;
300 :
301 : // Do not burn the end point, unless we are in the last
302 : // segment. This is to avoid burning twice intermediate points,
303 : // which causes artifacts in Add mode
304 38876 : if (j != panPartSize[i] - 1)
305 : {
306 37185 : nDeltaX--;
307 : }
308 :
309 67162 : while (nDeltaX-- >= 0)
310 : {
311 28286 : if (0 <= iX && iX < nRasterXSize && 0 <= iY &&
312 : iY < nRasterYSize)
313 28089 : pfnPointFunc(pCBData, iY, iX, dfVariant);
314 :
315 28286 : dfVariant += dfDeltaVariant;
316 28286 : iX += nXStep;
317 28286 : if (nError > 0)
318 : {
319 11024 : iY += nYStep;
320 11024 : nError += nYError;
321 : }
322 : else
323 : {
324 17262 : nError += nXError;
325 : }
326 : }
327 : }
328 : else
329 : {
330 7917 : const int nXError = nDeltaX << 1;
331 7917 : const int nYError = nXError - (nDeltaY << 1);
332 7917 : int nError = nXError - nDeltaY;
333 : // == 0 makes clang -fcatch-undefined-behavior -ftrapv happy,
334 : // but if it is == 0, dfDeltaVariant is not really used, so any
335 : // value is okay.
336 7917 : double dfDeltaVariant =
337 7917 : nDeltaY == 0 ? 0.0 : (dfVariant1 - dfVariant) / nDeltaY;
338 :
339 : // Do not burn the end point, unless we are in the last
340 : // segment. This is to avoid burning twice intermediate points,
341 : // which causes artifacts in Add mode
342 7917 : if (j != panPartSize[i] - 1)
343 : {
344 7734 : nDeltaY--;
345 : }
346 :
347 16207 : while (nDeltaY-- >= 0)
348 : {
349 8290 : if (0 <= iX && iX < nRasterXSize && 0 <= iY &&
350 : iY < nRasterYSize)
351 8201 : pfnPointFunc(pCBData, iY, iX, dfVariant);
352 :
353 8290 : dfVariant += dfDeltaVariant;
354 8290 : iY += nYStep;
355 8290 : if (nError > 0)
356 : {
357 35 : iX += nXStep;
358 35 : nError += nYError;
359 : }
360 : else
361 : {
362 8255 : nError += nXError;
363 : }
364 : }
365 : }
366 : }
367 : }
368 : }
369 :
370 : /************************************************************************/
371 : /* GDALdllImageLineAllTouched() */
372 : /* */
373 : /* This alternate line drawing algorithm attempts to ensure */
374 : /* that every pixel touched at all by the line will get set. */
375 : /* @param padfVariant should contain the values that are to be */
376 : /* added to the burn value. The values along the line between the */
377 : /* points will be linearly interpolated. These values are used */
378 : /* only if pCBData->eBurnValueSource is set to something other */
379 : /* than GBV_UserBurnValue. If NULL is passed, a monotonous line */
380 : /* will be drawn with the burn value. */
381 : /************************************************************************/
382 :
383 84 : void GDALdllImageLineAllTouched(
384 : int nRasterXSize, int nRasterYSize, int nPartCount, const int *panPartSize,
385 : const double *padfX, const double *padfY, const double *padfVariant,
386 : llPointFunc pfnPointFunc, GDALRasterizeInfo *pCBData,
387 : bool bAvoidBurningSamePoints, bool bIntersectOnly)
388 :
389 : {
390 : // This is an epsilon to detect geometries that are aligned with pixel
391 : // coordinates. Hard to find the right value. We put it to that value
392 : // to satisfy the scenarios of https://github.com/OSGeo/gdal/issues/7523
393 : // and https://github.com/OSGeo/gdal/issues/6414
394 84 : constexpr double EPSILON_INTERSECT_ONLY = 1e-4;
395 :
396 84 : if (!nPartCount)
397 0 : return;
398 :
399 171 : for (int i = 0, n = 0; i < nPartCount; n += panPartSize[i++])
400 : {
401 174 : std::set<std::pair<int, int>> lastBurntPoints;
402 174 : std::set<std::pair<int, int>> newBurntPoints;
403 :
404 412 : for (int j = 1; j < panPartSize[i]; j++)
405 : {
406 325 : lastBurntPoints = std::move(newBurntPoints);
407 325 : newBurntPoints.clear();
408 :
409 325 : double dfX = padfX[n + j - 1];
410 325 : double dfY = padfY[n + j - 1];
411 :
412 325 : double dfXEnd = padfX[n + j];
413 325 : double dfYEnd = padfY[n + j];
414 :
415 325 : double dfVariant = 0.0;
416 325 : double dfVariantEnd = 0.0;
417 325 : if (padfVariant != nullptr &&
418 11 : static_cast<GDALRasterizeInfo *>(pCBData)->eBurnValueSource !=
419 : GBV_UserBurnValue)
420 : {
421 11 : dfVariant = padfVariant[n + j - 1];
422 11 : dfVariantEnd = padfVariant[n + j];
423 : }
424 :
425 : // Skip segments that are off the target region.
426 325 : if ((dfY < 0.0 && dfYEnd < 0.0) ||
427 323 : (dfY > nRasterYSize && dfYEnd > nRasterYSize) ||
428 321 : (dfX < 0.0 && dfXEnd < 0.0) ||
429 319 : (dfX > nRasterXSize && dfXEnd > nRasterXSize))
430 280 : continue;
431 :
432 : // Swap if needed so we can proceed from left2right (X increasing)
433 316 : if (dfX > dfXEnd)
434 : {
435 104 : std::swap(dfX, dfXEnd);
436 104 : std::swap(dfY, dfYEnd);
437 104 : std::swap(dfVariant, dfVariantEnd);
438 : }
439 :
440 : // Special case for vertical lines.
441 316 : if (floor(dfX) == floor(dfXEnd) || fabs(dfX - dfXEnd) < .01)
442 : {
443 137 : if (bIntersectOnly)
444 : {
445 128 : if (std::abs(dfX - std::round(dfX)) <
446 139 : EPSILON_INTERSECT_ONLY &&
447 11 : std::abs(dfXEnd - std::round(dfXEnd)) <
448 : EPSILON_INTERSECT_ONLY)
449 11 : continue;
450 : }
451 :
452 126 : if (dfYEnd < dfY)
453 : {
454 69 : std::swap(dfY, dfYEnd);
455 69 : std::swap(dfVariant, dfVariantEnd);
456 : }
457 :
458 126 : const int iX = static_cast<int>(floor(dfXEnd));
459 126 : int iY = static_cast<int>(floor(dfY));
460 126 : int iYEnd =
461 126 : static_cast<int>(floor(dfYEnd - EPSILON_INTERSECT_ONLY));
462 :
463 126 : if (iX < 0 || iX >= nRasterXSize)
464 0 : continue;
465 :
466 126 : double dfDeltaVariant = 0.0;
467 126 : if (dfYEnd - dfY > 0.0)
468 126 : dfDeltaVariant = (dfVariantEnd - dfVariant) /
469 126 : (dfYEnd - dfY); // Per unit change in iY.
470 :
471 : // Clip to the borders of the target region.
472 126 : if (iY < 0)
473 0 : iY = 0;
474 126 : if (iYEnd >= nRasterYSize)
475 0 : iYEnd = nRasterYSize - 1;
476 126 : dfVariant += dfDeltaVariant * (iY - dfY);
477 :
478 126 : if (padfVariant == nullptr)
479 : {
480 1235 : for (; iY <= iYEnd; iY++)
481 : {
482 1114 : if (bAvoidBurningSamePoints)
483 : {
484 149 : auto yx = std::pair<int, int>(iY, iX);
485 149 : if (lastBurntPoints.find(yx) !=
486 298 : lastBurntPoints.end())
487 : {
488 15 : continue;
489 : }
490 134 : newBurntPoints.insert(yx);
491 : }
492 1099 : pfnPointFunc(pCBData, iY, iX, 0.0);
493 : }
494 : }
495 : else
496 : {
497 18 : for (; iY <= iYEnd; iY++, dfVariant += dfDeltaVariant)
498 : {
499 13 : if (bAvoidBurningSamePoints)
500 : {
501 0 : auto yx = std::pair<int, int>(iY, iX);
502 0 : if (lastBurntPoints.find(yx) !=
503 0 : lastBurntPoints.end())
504 : {
505 0 : continue;
506 : }
507 0 : newBurntPoints.insert(yx);
508 : }
509 13 : pfnPointFunc(pCBData, iY, iX, dfVariant);
510 : }
511 : }
512 :
513 126 : continue; // Next segment.
514 : }
515 :
516 179 : const double dfDeltaVariant =
517 179 : (dfVariantEnd - dfVariant) /
518 179 : (dfXEnd - dfX); // Per unit change in iX.
519 :
520 : // Special case for horizontal lines.
521 179 : if (floor(dfY) == floor(dfYEnd) || fabs(dfY - dfYEnd) < .01)
522 : {
523 134 : if (bIntersectOnly)
524 : {
525 125 : if (std::abs(dfY - std::round(dfY)) <
526 133 : EPSILON_INTERSECT_ONLY &&
527 8 : std::abs(dfYEnd - std::round(dfYEnd)) <
528 : EPSILON_INTERSECT_ONLY)
529 8 : continue;
530 : }
531 :
532 126 : if (dfXEnd < dfX)
533 : {
534 0 : std::swap(dfX, dfXEnd);
535 0 : std::swap(dfVariant, dfVariantEnd);
536 : }
537 :
538 126 : int iX = static_cast<int>(floor(dfX));
539 126 : const int iY = static_cast<int>(floor(dfY));
540 126 : int iXEnd =
541 126 : static_cast<int>(floor(dfXEnd - EPSILON_INTERSECT_ONLY));
542 :
543 126 : if (iY < 0 || iY >= nRasterYSize)
544 0 : continue;
545 :
546 : // Clip to the borders of the target region.
547 126 : if (iX < 0)
548 0 : iX = 0;
549 126 : if (iXEnd >= nRasterXSize)
550 0 : iXEnd = nRasterXSize - 1;
551 126 : dfVariant += dfDeltaVariant * (iX - dfX);
552 :
553 126 : if (padfVariant == nullptr)
554 : {
555 1000 : for (; iX <= iXEnd; iX++)
556 : {
557 879 : if (bAvoidBurningSamePoints)
558 : {
559 147 : auto yx = std::pair<int, int>(iY, iX);
560 147 : if (lastBurntPoints.find(yx) !=
561 294 : lastBurntPoints.end())
562 : {
563 13 : continue;
564 : }
565 134 : newBurntPoints.insert(yx);
566 : }
567 866 : pfnPointFunc(pCBData, iY, iX, 0.0);
568 : }
569 : }
570 : else
571 : {
572 18 : for (; iX <= iXEnd; iX++, dfVariant += dfDeltaVariant)
573 : {
574 13 : if (bAvoidBurningSamePoints)
575 : {
576 0 : auto yx = std::pair<int, int>(iY, iX);
577 0 : if (lastBurntPoints.find(yx) !=
578 0 : lastBurntPoints.end())
579 : {
580 0 : continue;
581 : }
582 0 : newBurntPoints.insert(yx);
583 : }
584 13 : pfnPointFunc(pCBData, iY, iX, dfVariant);
585 : }
586 : }
587 :
588 126 : continue; // Next segment.
589 : }
590 :
591 : /* --------------------------------------------------------------------
592 : */
593 : /* General case - left to right sloped. */
594 : /* --------------------------------------------------------------------
595 : */
596 45 : const double dfSlope = (dfYEnd - dfY) / (dfXEnd - dfX);
597 :
598 : // Clip segment in X.
599 45 : if (dfXEnd > nRasterXSize)
600 : {
601 0 : dfYEnd -= (dfXEnd - nRasterXSize) * dfSlope;
602 0 : dfXEnd = nRasterXSize;
603 : }
604 45 : if (dfX < 0.0)
605 : {
606 1 : dfY += (0.0 - dfX) * dfSlope;
607 1 : dfVariant += dfDeltaVariant * (0.0 - dfX);
608 1 : dfX = 0.0;
609 : }
610 :
611 : // Clip segment in Y.
612 45 : if (dfYEnd > dfY)
613 : {
614 6 : if (dfY < 0.0)
615 : {
616 0 : const double dfDiffX = (0.0 - dfY) / dfSlope;
617 0 : dfX += dfDiffX;
618 0 : dfVariant += dfDeltaVariant * dfDiffX;
619 0 : dfY = 0.0;
620 : }
621 6 : if (dfYEnd >= nRasterYSize)
622 : {
623 2 : dfXEnd += (dfYEnd - nRasterYSize) / dfSlope;
624 2 : if (dfXEnd > nRasterXSize)
625 2 : dfXEnd = nRasterXSize;
626 : // dfYEnd is no longer used afterwards, but for
627 : // consistency it should be:
628 : // dfYEnd = nRasterXSize;
629 : }
630 : }
631 : else
632 : {
633 39 : if (dfY >= nRasterYSize)
634 : {
635 0 : const double dfDiffX = (nRasterYSize - dfY) / dfSlope;
636 0 : dfX += dfDiffX;
637 0 : dfVariant += dfDeltaVariant * dfDiffX;
638 0 : dfY = nRasterYSize;
639 : }
640 39 : if (dfYEnd < 0.0)
641 : {
642 0 : dfXEnd -= (dfYEnd - 0) / dfSlope;
643 : // dfYEnd is no longer used afterwards, but for
644 : // consistency it should be:
645 : // dfYEnd = 0.0;
646 : }
647 : }
648 :
649 : // Step from pixel to pixel.
650 4690 : while (dfX >= 0.0 && dfX < dfXEnd)
651 : {
652 4645 : const int iX = static_cast<int>(floor(dfX));
653 4645 : const int iY = static_cast<int>(floor(dfY));
654 :
655 : // Burn in the current point.
656 : // We should be able to drop the Y check because we clipped
657 : // in Y, but there may be some error with all the small steps.
658 4645 : if (iY >= 0 && iY < nRasterYSize)
659 : {
660 4644 : if (bAvoidBurningSamePoints)
661 : {
662 34 : auto yx = std::pair<int, int>(iY, iX);
663 67 : if (lastBurntPoints.find(yx) == lastBurntPoints.end() &&
664 67 : newBurntPoints.find(yx) == newBurntPoints.end())
665 : {
666 24 : newBurntPoints.insert(yx);
667 24 : pfnPointFunc(pCBData, iY, iX, dfVariant);
668 : }
669 : }
670 : else
671 : {
672 4610 : pfnPointFunc(pCBData, iY, iX, dfVariant);
673 : }
674 : }
675 :
676 4645 : double dfStepX = floor(dfX + 1.0) - dfX;
677 4645 : double dfStepY = dfStepX * dfSlope;
678 :
679 : // Step to right pixel without changing scanline?
680 4645 : if (static_cast<int>(floor(dfY + dfStepY)) == iY)
681 : {
682 2992 : dfX += dfStepX;
683 2992 : dfY += dfStepY;
684 2992 : dfVariant += dfDeltaVariant * dfStepX;
685 : }
686 1653 : else if (dfSlope < 0)
687 : {
688 1122 : dfStepY = iY - dfY;
689 1122 : if (dfStepY > -0.000000001)
690 562 : dfStepY = -0.000000001;
691 :
692 1122 : dfStepX = dfStepY / dfSlope;
693 1122 : dfX += dfStepX;
694 1122 : dfY += dfStepY;
695 1122 : dfVariant += dfDeltaVariant * dfStepX;
696 : }
697 : else
698 : {
699 531 : dfStepY = (iY + 1) - dfY;
700 531 : if (dfStepY < 0.000000001)
701 0 : dfStepY = 0.000000001;
702 :
703 531 : dfStepX = dfStepY / dfSlope;
704 531 : dfX += dfStepX;
705 531 : dfY += dfStepY;
706 531 : dfVariant += dfDeltaVariant * dfStepX;
707 : }
708 : } // Next step along segment.
709 : } // Next segment.
710 : } // Next part.
711 : }
|
