Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: GDAL Core
4 : * Purpose: Implementation of GDALNoDataMaskBand, a class implementing all
5 : * a default band mask based on nodata values.
6 : * Author: Frank Warmerdam, warmerdam@pobox.com
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 2007, Frank Warmerdam
10 : * Copyright (c) 2008-2012, Even Rouault <even dot rouault at spatialys.com>
11 : *
12 : * SPDX-License-Identifier: MIT
13 : ****************************************************************************/
14 :
15 : #include "cpl_port.h"
16 : #include "gdal_priv.h"
17 :
18 : #include <algorithm>
19 : #include <cmath>
20 : #include <cstring>
21 : #include <utility>
22 :
23 : #include "cpl_conv.h"
24 : #include "cpl_error.h"
25 : #include "cpl_vsi.h"
26 : #include "gdal.h"
27 : #include "gdal_priv_templates.hpp"
28 :
29 : //! @cond Doxygen_Suppress
30 : /************************************************************************/
31 : /* GDALNoDataMaskBand() */
32 : /************************************************************************/
33 :
34 797 : GDALNoDataMaskBand::GDALNoDataMaskBand(GDALRasterBand *poParentIn)
35 797 : : m_poParent(poParentIn)
36 : {
37 797 : poDS = nullptr;
38 797 : nBand = 0;
39 :
40 797 : nRasterXSize = m_poParent->GetXSize();
41 797 : nRasterYSize = m_poParent->GetYSize();
42 :
43 797 : eDataType = GDT_Byte;
44 797 : m_poParent->GetBlockSize(&nBlockXSize, &nBlockYSize);
45 :
46 797 : const auto eParentDT = m_poParent->GetRasterDataType();
47 797 : if (eParentDT == GDT_Int64)
48 13 : m_nNoDataValueInt64 = m_poParent->GetNoDataValueAsInt64();
49 784 : else if (eParentDT == GDT_UInt64)
50 12 : m_nNoDataValueUInt64 = m_poParent->GetNoDataValueAsUInt64();
51 : else
52 772 : m_dfNoDataValue = m_poParent->GetNoDataValue();
53 797 : }
54 :
55 : /************************************************************************/
56 : /* GDALNoDataMaskBand() */
57 : /************************************************************************/
58 :
59 1 : GDALNoDataMaskBand::GDALNoDataMaskBand(GDALRasterBand *poParentIn,
60 1 : double dfNoDataValue)
61 1 : : m_poParent(poParentIn)
62 : {
63 1 : poDS = nullptr;
64 1 : nBand = 0;
65 :
66 1 : nRasterXSize = m_poParent->GetXSize();
67 1 : nRasterYSize = m_poParent->GetYSize();
68 :
69 1 : eDataType = GDT_Byte;
70 1 : m_poParent->GetBlockSize(&nBlockXSize, &nBlockYSize);
71 :
72 1 : const auto eParentDT = m_poParent->GetRasterDataType();
73 1 : if (eParentDT == GDT_Int64)
74 0 : m_nNoDataValueInt64 = static_cast<int64_t>(dfNoDataValue);
75 1 : else if (eParentDT == GDT_UInt64)
76 0 : m_nNoDataValueUInt64 = static_cast<uint64_t>(dfNoDataValue);
77 : else
78 1 : m_dfNoDataValue = dfNoDataValue;
79 1 : }
80 :
81 : /************************************************************************/
82 : /* ~GDALNoDataMaskBand() */
83 : /************************************************************************/
84 :
85 : GDALNoDataMaskBand::~GDALNoDataMaskBand() = default;
86 :
87 : /************************************************************************/
88 : /* GetWorkDataType() */
89 : /************************************************************************/
90 :
91 10304 : static GDALDataType GetWorkDataType(GDALDataType eDataType)
92 : {
93 10304 : GDALDataType eWrkDT = GDT_Unknown;
94 10304 : switch (eDataType)
95 : {
96 793 : case GDT_Byte:
97 793 : eWrkDT = GDT_Byte;
98 793 : break;
99 :
100 1741 : case GDT_Int16:
101 1741 : eWrkDT = GDT_Int16;
102 1741 : break;
103 :
104 48 : case GDT_UInt16:
105 48 : eWrkDT = GDT_UInt16;
106 48 : break;
107 :
108 29 : case GDT_UInt32:
109 29 : eWrkDT = GDT_UInt32;
110 29 : break;
111 :
112 101 : case GDT_Int8:
113 : case GDT_Int32:
114 : case GDT_CInt16:
115 : case GDT_CInt32:
116 101 : eWrkDT = GDT_Int32;
117 101 : break;
118 :
119 7493 : case GDT_Float32:
120 : case GDT_CFloat32:
121 7493 : eWrkDT = GDT_Float32;
122 7493 : break;
123 :
124 71 : case GDT_Float64:
125 : case GDT_CFloat64:
126 71 : eWrkDT = GDT_Float64;
127 71 : break;
128 :
129 28 : case GDT_Int64:
130 : case GDT_UInt64:
131 28 : eWrkDT = eDataType;
132 28 : break;
133 :
134 0 : case GDT_Unknown:
135 : case GDT_TypeCount:
136 0 : CPLAssert(false);
137 : eWrkDT = GDT_Float64;
138 : break;
139 : }
140 10304 : return eWrkDT;
141 : }
142 :
143 : /************************************************************************/
144 : /* IsNoDataInRange() */
145 : /************************************************************************/
146 :
147 804 : bool GDALNoDataMaskBand::IsNoDataInRange(double dfNoDataValue,
148 : GDALDataType eDataTypeIn)
149 : {
150 804 : GDALDataType eWrkDT = GetWorkDataType(eDataTypeIn);
151 804 : switch (eWrkDT)
152 : {
153 345 : case GDT_Byte:
154 : {
155 345 : return GDALIsValueInRange<GByte>(dfNoDataValue);
156 : }
157 :
158 69 : case GDT_Int16:
159 : {
160 69 : return GDALIsValueInRange<GInt16>(dfNoDataValue);
161 : }
162 :
163 38 : case GDT_UInt16:
164 : {
165 38 : return GDALIsValueInRange<GUInt16>(dfNoDataValue);
166 : }
167 :
168 20 : case GDT_UInt32:
169 : {
170 20 : return GDALIsValueInRange<GUInt32>(dfNoDataValue);
171 : }
172 24 : case GDT_Int32:
173 : {
174 24 : return GDALIsValueInRange<GInt32>(dfNoDataValue);
175 : }
176 :
177 0 : case GDT_UInt64:
178 : {
179 0 : return GDALIsValueInRange<uint64_t>(dfNoDataValue);
180 : }
181 :
182 0 : case GDT_Int64:
183 : {
184 0 : return GDALIsValueInRange<int64_t>(dfNoDataValue);
185 : }
186 :
187 262 : case GDT_Float32:
188 : {
189 510 : return std::isnan(dfNoDataValue) || std::isinf(dfNoDataValue) ||
190 510 : GDALIsValueInRange<float>(dfNoDataValue);
191 : }
192 :
193 46 : case GDT_Float64:
194 : {
195 46 : return true;
196 : }
197 :
198 0 : default:
199 0 : CPLAssert(false);
200 : return false;
201 : }
202 : }
203 :
204 : /************************************************************************/
205 : /* IReadBlock() */
206 : /************************************************************************/
207 :
208 26 : CPLErr GDALNoDataMaskBand::IReadBlock(int nXBlockOff, int nYBlockOff,
209 : void *pImage)
210 :
211 : {
212 26 : const int nXOff = nXBlockOff * nBlockXSize;
213 26 : const int nXSizeRequest = std::min(nBlockXSize, nRasterXSize - nXOff);
214 26 : const int nYOff = nYBlockOff * nBlockYSize;
215 26 : const int nYSizeRequest = std::min(nBlockYSize, nRasterYSize - nYOff);
216 :
217 26 : if (nBlockXSize != nXSizeRequest || nBlockYSize != nYSizeRequest)
218 : {
219 0 : memset(pImage, 0, static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize);
220 : }
221 :
222 : GDALRasterIOExtraArg sExtraArg;
223 26 : INIT_RASTERIO_EXTRA_ARG(sExtraArg);
224 26 : return IRasterIO(GF_Read, nXOff, nYOff, nXSizeRequest, nYSizeRequest,
225 : pImage, nXSizeRequest, nYSizeRequest, GDT_Byte, 1,
226 52 : nBlockXSize, &sExtraArg);
227 : }
228 :
229 : /************************************************************************/
230 : /* SetZeroOr255() */
231 : /************************************************************************/
232 :
233 : #if (defined(__GNUC__) && !defined(__clang__))
234 : __attribute__((optimize("tree-vectorize")))
235 : #endif
236 : static void
237 387 : SetZeroOr255(GByte *pabyDestAndSrc, size_t nBufSize, GByte byNoData)
238 : {
239 4074150 : for (size_t i = 0; i < nBufSize; ++i)
240 : {
241 4073760 : pabyDestAndSrc[i] = (pabyDestAndSrc[i] == byNoData) ? 0 : 255;
242 : }
243 387 : }
244 :
245 : template <class T>
246 : #if (defined(__GNUC__) && !defined(__clang__))
247 : __attribute__((optimize("tree-vectorize")))
248 : #endif
249 : static void
250 1776 : SetZeroOr255(GByte *pabyDest, const T *panSrc, size_t nBufSize, T nNoData)
251 : {
252 821389 : for (size_t i = 0; i < nBufSize; ++i)
253 : {
254 819613 : pabyDest[i] = (panSrc[i] == nNoData) ? 0 : 255;
255 : }
256 1776 : }
257 :
258 : template <class T>
259 1776 : static void SetZeroOr255(GByte *pabyDest, const T *panSrc, int nBufXSize,
260 : int nBufYSize, GSpacing nPixelSpace,
261 : GSpacing nLineSpace, T nNoData)
262 : {
263 1776 : if (nPixelSpace == 1 && nLineSpace == nBufXSize)
264 : {
265 1776 : const size_t nBufSize = static_cast<size_t>(nBufXSize) * nBufYSize;
266 1776 : SetZeroOr255(pabyDest, panSrc, nBufSize, nNoData);
267 : }
268 0 : else if (nPixelSpace == 1)
269 : {
270 0 : for (int iY = 0; iY < nBufYSize; iY++)
271 : {
272 0 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nNoData);
273 0 : pabyDest += nLineSpace;
274 0 : panSrc += nBufXSize;
275 : }
276 : }
277 : else
278 : {
279 0 : size_t i = 0;
280 0 : for (int iY = 0; iY < nBufYSize; iY++)
281 : {
282 0 : GByte *pabyLineDest = pabyDest + iY * nLineSpace;
283 0 : for (int iX = 0; iX < nBufXSize; iX++)
284 : {
285 0 : *pabyLineDest = (panSrc[i] == nNoData) ? 0 : 255;
286 0 : ++i;
287 0 : pabyLineDest += nPixelSpace;
288 : }
289 : }
290 : }
291 1776 : }
292 :
293 : /************************************************************************/
294 : /* IRasterIO() */
295 : /************************************************************************/
296 :
297 9500 : CPLErr GDALNoDataMaskBand::IRasterIO(GDALRWFlag eRWFlag, int nXOff, int nYOff,
298 : int nXSize, int nYSize, void *pData,
299 : int nBufXSize, int nBufYSize,
300 : GDALDataType eBufType,
301 : GSpacing nPixelSpace, GSpacing nLineSpace,
302 : GDALRasterIOExtraArg *psExtraArg)
303 : {
304 9500 : if (eRWFlag != GF_Read)
305 : {
306 0 : return CE_Failure;
307 : }
308 9500 : const auto eParentDT = m_poParent->GetRasterDataType();
309 9500 : const GDALDataType eWrkDT = GetWorkDataType(eParentDT);
310 :
311 : // Optimization in common use case (#4488).
312 : // This avoids triggering the block cache on this band, which helps
313 : // reducing the global block cache consumption.
314 9500 : if (eBufType == GDT_Byte && eWrkDT == GDT_Byte)
315 : {
316 387 : const CPLErr eErr = m_poParent->RasterIO(
317 : GF_Read, nXOff, nYOff, nXSize, nYSize, pData, nBufXSize, nBufYSize,
318 : eBufType, nPixelSpace, nLineSpace, psExtraArg);
319 387 : if (eErr != CE_None)
320 0 : return eErr;
321 :
322 387 : GByte *pabyData = static_cast<GByte *>(pData);
323 387 : const GByte byNoData = static_cast<GByte>(m_dfNoDataValue);
324 :
325 387 : if (nPixelSpace == 1 && nLineSpace == nBufXSize)
326 : {
327 387 : const size_t nBufSize = static_cast<size_t>(nBufXSize) * nBufYSize;
328 387 : SetZeroOr255(pabyData, nBufSize, byNoData);
329 : }
330 : else
331 : {
332 0 : SetZeroOr255(pabyData, pabyData, nBufXSize, nBufYSize, nPixelSpace,
333 : nLineSpace, byNoData);
334 : }
335 387 : return CE_None;
336 : }
337 :
338 : const auto AllocTempBufferOrFallback =
339 9113 : [this, eRWFlag, nXOff, nYOff, nXSize, nYSize, pData, nBufXSize,
340 : nBufYSize, eBufType, nPixelSpace, nLineSpace,
341 18282 : psExtraArg](int nWrkDTSize) -> std::pair<CPLErr, void *>
342 : {
343 9113 : auto poParentDS = m_poParent->GetDataset();
344 : // Check if we must simulate a memory allocation failure
345 : // Before checking the env variable, which is slightly expensive,
346 : // check first for a special dataset name, which is a cheap test.
347 : const char *pszOptVal =
348 9113 : poParentDS && strcmp(poParentDS->GetDescription(), "__debug__") == 0
349 18226 : ? CPLGetConfigOption(
350 : "GDAL_SIMUL_MEM_ALLOC_FAILURE_NODATA_MASK_BAND", "NO")
351 9113 : : "NO";
352 : const bool bSimulMemAllocFailure =
353 18222 : EQUAL(pszOptVal, "ALWAYS") ||
354 9109 : (CPLTestBool(pszOptVal) &&
355 14 : GDALMajorObject::GetMetadataItem(__func__, "__INTERNAL__") ==
356 9113 : nullptr);
357 9113 : void *pTemp = nullptr;
358 9113 : if (!bSimulMemAllocFailure)
359 : {
360 9101 : CPLErrorStateBackuper oErrorStateBackuper(CPLQuietErrorHandler);
361 9101 : pTemp = VSI_MALLOC3_VERBOSE(nWrkDTSize, nBufXSize, nBufYSize);
362 : }
363 9113 : if (!pTemp)
364 : {
365 : const bool bAllocHasAlreadyFailed =
366 12 : GDALMajorObject::GetMetadataItem(__func__, "__INTERNAL__") !=
367 12 : nullptr;
368 12 : CPLError(bAllocHasAlreadyFailed ? CE_Failure : CE_Warning,
369 : CPLE_OutOfMemory,
370 : "GDALNoDataMaskBand::IRasterIO(): cannot allocate %d x %d "
371 : "x %d bytes%s",
372 : nBufXSize, nBufYSize, nWrkDTSize,
373 : bAllocHasAlreadyFailed
374 : ? ""
375 : : ". Falling back to block-based approach");
376 12 : if (bAllocHasAlreadyFailed)
377 2 : return std::pair(CE_Failure, nullptr);
378 : // Sets a metadata item to prevent potential infinite recursion
379 10 : GDALMajorObject::SetMetadataItem(__func__, "IN", "__INTERNAL__");
380 10 : const CPLErr eErr = GDALRasterBand::IRasterIO(
381 : eRWFlag, nXOff, nYOff, nXSize, nYSize, pData, nBufXSize,
382 10 : nBufYSize, eBufType, nPixelSpace, nLineSpace, psExtraArg);
383 10 : GDALMajorObject::SetMetadataItem(__func__, nullptr, "__INTERNAL__");
384 10 : return std::pair(eErr, nullptr);
385 : }
386 9101 : return std::pair(CE_None, pTemp);
387 9113 : };
388 :
389 9113 : if (eBufType == GDT_Byte)
390 : {
391 9034 : const int nWrkDTSize = GDALGetDataTypeSizeBytes(eWrkDT);
392 9034 : auto [eErr, pTemp] = AllocTempBufferOrFallback(nWrkDTSize);
393 9034 : if (!pTemp)
394 12 : return eErr;
395 :
396 18044 : eErr = m_poParent->RasterIO(
397 : GF_Read, nXOff, nYOff, nXSize, nYSize, pTemp, nBufXSize, nBufYSize,
398 9022 : eWrkDT, nWrkDTSize, static_cast<GSpacing>(nBufXSize) * nWrkDTSize,
399 : psExtraArg);
400 9022 : if (eErr != CE_None)
401 : {
402 0 : VSIFree(pTemp);
403 0 : return eErr;
404 : }
405 :
406 9022 : const bool bIsNoDataNan = std::isnan(m_dfNoDataValue) != 0;
407 9022 : GByte *pabyDest = static_cast<GByte *>(pData);
408 :
409 : /* --------------------------------------------------------------------
410 : */
411 : /* Process different cases. */
412 : /* --------------------------------------------------------------------
413 : */
414 9022 : switch (eWrkDT)
415 : {
416 1670 : case GDT_Int16:
417 : {
418 1670 : const auto nNoData = static_cast<int16_t>(m_dfNoDataValue);
419 1670 : const auto *panSrc = static_cast<const int16_t *>(pTemp);
420 1670 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nBufYSize,
421 : nPixelSpace, nLineSpace, nNoData);
422 : }
423 1670 : break;
424 :
425 8 : case GDT_UInt16:
426 : {
427 8 : const auto nNoData = static_cast<uint16_t>(m_dfNoDataValue);
428 8 : const auto *panSrc = static_cast<const uint16_t *>(pTemp);
429 8 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nBufYSize,
430 : nPixelSpace, nLineSpace, nNoData);
431 : }
432 8 : break;
433 :
434 7 : case GDT_UInt32:
435 : {
436 7 : const auto nNoData = static_cast<GUInt32>(m_dfNoDataValue);
437 7 : const auto *panSrc = static_cast<const GUInt32 *>(pTemp);
438 7 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nBufYSize,
439 : nPixelSpace, nLineSpace, nNoData);
440 : }
441 7 : break;
442 :
443 75 : case GDT_Int32:
444 : {
445 75 : const auto nNoData = static_cast<GInt32>(m_dfNoDataValue);
446 75 : const auto *panSrc = static_cast<const GInt32 *>(pTemp);
447 75 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nBufYSize,
448 : nPixelSpace, nLineSpace, nNoData);
449 : }
450 75 : break;
451 :
452 7226 : case GDT_Float32:
453 : {
454 7226 : const float fNoData = static_cast<float>(m_dfNoDataValue);
455 7226 : const float *pafSrc = static_cast<const float *>(pTemp);
456 :
457 7226 : size_t i = 0;
458 14516 : for (int iY = 0; iY < nBufYSize; iY++)
459 : {
460 7290 : GByte *pabyLineDest = pabyDest + iY * nLineSpace;
461 1962120 : for (int iX = 0; iX < nBufXSize; iX++)
462 : {
463 1954830 : const float fVal = pafSrc[i];
464 1954830 : if (bIsNoDataNan && std::isnan(fVal))
465 20 : *pabyLineDest = 0;
466 1954810 : else if (ARE_REAL_EQUAL(fVal, fNoData))
467 876112 : *pabyLineDest = 0;
468 : else
469 1078690 : *pabyLineDest = 255;
470 1954830 : ++i;
471 1954830 : pabyLineDest += nPixelSpace;
472 : }
473 : }
474 : }
475 7226 : break;
476 :
477 20 : case GDT_Float64:
478 : {
479 20 : const double *padfSrc = static_cast<const double *>(pTemp);
480 :
481 20 : size_t i = 0;
482 64 : for (int iY = 0; iY < nBufYSize; iY++)
483 : {
484 44 : GByte *pabyLineDest = pabyDest + iY * nLineSpace;
485 209 : for (int iX = 0; iX < nBufXSize; iX++)
486 : {
487 165 : const double dfVal = padfSrc[i];
488 165 : if (bIsNoDataNan && std::isnan(dfVal))
489 20 : *pabyLineDest = 0;
490 145 : else if (ARE_REAL_EQUAL(dfVal, m_dfNoDataValue))
491 82 : *pabyLineDest = 0;
492 : else
493 63 : *pabyLineDest = 255;
494 165 : ++i;
495 165 : pabyLineDest += nPixelSpace;
496 : }
497 : }
498 : }
499 20 : break;
500 :
501 8 : case GDT_Int64:
502 : {
503 8 : const auto *panSrc = static_cast<const int64_t *>(pTemp);
504 8 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nBufYSize,
505 : nPixelSpace, nLineSpace, m_nNoDataValueInt64);
506 : }
507 8 : break;
508 :
509 8 : case GDT_UInt64:
510 : {
511 8 : const auto *panSrc = static_cast<const uint64_t *>(pTemp);
512 8 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nBufYSize,
513 : nPixelSpace, nLineSpace, m_nNoDataValueUInt64);
514 : }
515 8 : break;
516 :
517 0 : default:
518 0 : CPLAssert(false);
519 : break;
520 : }
521 :
522 9022 : VSIFree(pTemp);
523 9022 : return CE_None;
524 : }
525 :
526 : // Output buffer is non-Byte. Ask for Byte and expand to user requested
527 : // type
528 79 : auto [eErr, pTemp] = AllocTempBufferOrFallback(sizeof(GByte));
529 79 : if (!pTemp)
530 0 : return eErr;
531 :
532 158 : eErr = IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize, pTemp, nBufXSize,
533 79 : nBufYSize, GDT_Byte, 1, nBufXSize, psExtraArg);
534 79 : if (eErr != CE_None)
535 : {
536 0 : VSIFree(pTemp);
537 0 : return eErr;
538 : }
539 :
540 186 : for (int iY = 0; iY < nBufYSize; iY++)
541 : {
542 107 : GDALCopyWords(
543 107 : static_cast<GByte *>(pTemp) + static_cast<size_t>(iY) * nBufXSize,
544 107 : GDT_Byte, 1, static_cast<GByte *>(pData) + iY * nLineSpace,
545 : eBufType, static_cast<int>(nPixelSpace), nBufXSize);
546 : }
547 79 : VSIFree(pTemp);
548 79 : return CE_None;
549 : }
550 :
551 : //! @endcond
|