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