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 870 : GDALNoDataMaskBand::GDALNoDataMaskBand(GDALRasterBand *poParentIn)
35 870 : : m_poParent(poParentIn)
36 : {
37 870 : poDS = nullptr;
38 870 : nBand = 0;
39 :
40 870 : nRasterXSize = m_poParent->GetXSize();
41 870 : nRasterYSize = m_poParent->GetYSize();
42 :
43 870 : eDataType = GDT_Byte;
44 870 : m_poParent->GetBlockSize(&nBlockXSize, &nBlockYSize);
45 :
46 870 : const auto eParentDT = m_poParent->GetRasterDataType();
47 870 : if (eParentDT == GDT_Int64)
48 13 : m_nNoDataValueInt64 = m_poParent->GetNoDataValueAsInt64();
49 857 : else if (eParentDT == GDT_UInt64)
50 12 : m_nNoDataValueUInt64 = m_poParent->GetNoDataValueAsUInt64();
51 : else
52 845 : m_dfNoDataValue = m_poParent->GetNoDataValue();
53 870 : }
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 10370 : static GDALDataType GetWorkDataType(GDALDataType eDataType)
92 : {
93 10370 : GDALDataType eWrkDT = GDT_Unknown;
94 10370 : switch (eDataType)
95 : {
96 795 : case GDT_Byte:
97 795 : eWrkDT = GDT_Byte;
98 795 : break;
99 :
100 1759 : case GDT_Int16:
101 1759 : eWrkDT = GDT_Int16;
102 1759 : break;
103 :
104 52 : case GDT_UInt16:
105 52 : eWrkDT = GDT_UInt16;
106 52 : 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 7511 : case GDT_Float16:
120 : case GDT_CFloat16:
121 : case GDT_Float32:
122 : case GDT_CFloat32:
123 7511 : eWrkDT = GDT_Float32;
124 7511 : break;
125 :
126 95 : case GDT_Float64:
127 : case GDT_CFloat64:
128 95 : eWrkDT = GDT_Float64;
129 95 : break;
130 :
131 28 : case GDT_Int64:
132 : case GDT_UInt64:
133 28 : eWrkDT = eDataType;
134 28 : break;
135 :
136 0 : case GDT_Unknown:
137 : case GDT_TypeCount:
138 0 : CPLAssert(false);
139 : eWrkDT = GDT_Float64;
140 : break;
141 : }
142 10370 : return eWrkDT;
143 : }
144 :
145 : /************************************************************************/
146 : /* IsNoDataInRange() */
147 : /************************************************************************/
148 :
149 845 : bool GDALNoDataMaskBand::IsNoDataInRange(double dfNoDataValue,
150 : GDALDataType eDataTypeIn)
151 : {
152 845 : GDALDataType eWrkDT = GetWorkDataType(eDataTypeIn);
153 845 : switch (eWrkDT)
154 : {
155 346 : case GDT_Byte:
156 : {
157 346 : return GDALIsValueInRange<GByte>(dfNoDataValue);
158 : }
159 :
160 81 : case GDT_Int16:
161 : {
162 81 : return GDALIsValueInRange<GInt16>(dfNoDataValue);
163 : }
164 :
165 42 : case GDT_UInt16:
166 : {
167 42 : return GDALIsValueInRange<GUInt16>(dfNoDataValue);
168 : }
169 :
170 20 : case GDT_UInt32:
171 : {
172 20 : return GDALIsValueInRange<GUInt32>(dfNoDataValue);
173 : }
174 24 : case GDT_Int32:
175 : {
176 24 : return GDALIsValueInRange<GInt32>(dfNoDataValue);
177 : }
178 :
179 0 : case GDT_UInt64:
180 : {
181 0 : return GDALIsValueInRange<uint64_t>(dfNoDataValue);
182 : }
183 :
184 0 : case GDT_Int64:
185 : {
186 0 : return GDALIsValueInRange<int64_t>(dfNoDataValue);
187 : }
188 :
189 0 : case GDT_Float16:
190 : {
191 0 : return isnan(dfNoDataValue) || isinf(dfNoDataValue) ||
192 0 : GDALIsValueInRange<GFloat16>(dfNoDataValue);
193 : }
194 :
195 274 : case GDT_Float32:
196 : {
197 522 : return std::isnan(dfNoDataValue) || std::isinf(dfNoDataValue) ||
198 522 : GDALIsValueInRange<float>(dfNoDataValue);
199 : }
200 :
201 58 : case GDT_Float64:
202 : {
203 58 : return true;
204 : }
205 :
206 0 : default:
207 0 : CPLAssert(false);
208 : return false;
209 : }
210 : }
211 :
212 : /************************************************************************/
213 : /* IReadBlock() */
214 : /************************************************************************/
215 :
216 26 : CPLErr GDALNoDataMaskBand::IReadBlock(int nXBlockOff, int nYBlockOff,
217 : void *pImage)
218 :
219 : {
220 26 : const int nXOff = nXBlockOff * nBlockXSize;
221 26 : const int nXSizeRequest = std::min(nBlockXSize, nRasterXSize - nXOff);
222 26 : const int nYOff = nYBlockOff * nBlockYSize;
223 26 : const int nYSizeRequest = std::min(nBlockYSize, nRasterYSize - nYOff);
224 :
225 26 : if (nBlockXSize != nXSizeRequest || nBlockYSize != nYSizeRequest)
226 : {
227 0 : memset(pImage, 0, static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize);
228 : }
229 :
230 : GDALRasterIOExtraArg sExtraArg;
231 26 : INIT_RASTERIO_EXTRA_ARG(sExtraArg);
232 26 : return IRasterIO(GF_Read, nXOff, nYOff, nXSizeRequest, nYSizeRequest,
233 : pImage, nXSizeRequest, nYSizeRequest, GDT_Byte, 1,
234 52 : nBlockXSize, &sExtraArg);
235 : }
236 :
237 : /************************************************************************/
238 : /* SetZeroOr255() */
239 : /************************************************************************/
240 :
241 : #if (defined(__GNUC__) && !defined(__clang__))
242 : __attribute__((optimize("tree-vectorize")))
243 : #endif
244 : static void
245 388 : SetZeroOr255(GByte *pabyDestAndSrc, size_t nBufSize, GByte byNoData)
246 : {
247 4074160 : for (size_t i = 0; i < nBufSize; ++i)
248 : {
249 4073770 : pabyDestAndSrc[i] = (pabyDestAndSrc[i] == byNoData) ? 0 : 255;
250 : }
251 388 : }
252 :
253 : template <class T>
254 : #if (defined(__GNUC__) && !defined(__clang__))
255 : __attribute__((optimize("tree-vectorize")))
256 : #endif
257 : static void
258 1782 : SetZeroOr255(GByte *pabyDest, const T *panSrc, size_t nBufSize, T nNoData)
259 : {
260 821431 : for (size_t i = 0; i < nBufSize; ++i)
261 : {
262 819649 : pabyDest[i] = (panSrc[i] == nNoData) ? 0 : 255;
263 : }
264 1782 : }
265 :
266 : template <class T>
267 1782 : static void SetZeroOr255(GByte *pabyDest, const T *panSrc, int nBufXSize,
268 : int nBufYSize, GSpacing nPixelSpace,
269 : GSpacing nLineSpace, T nNoData)
270 : {
271 1782 : if (nPixelSpace == 1 && nLineSpace == nBufXSize)
272 : {
273 1782 : const size_t nBufSize = static_cast<size_t>(nBufXSize) * nBufYSize;
274 1782 : SetZeroOr255(pabyDest, panSrc, nBufSize, nNoData);
275 : }
276 0 : else if (nPixelSpace == 1)
277 : {
278 0 : for (int iY = 0; iY < nBufYSize; iY++)
279 : {
280 0 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nNoData);
281 0 : pabyDest += nLineSpace;
282 0 : panSrc += nBufXSize;
283 : }
284 : }
285 : else
286 : {
287 0 : size_t i = 0;
288 0 : for (int iY = 0; iY < nBufYSize; iY++)
289 : {
290 0 : GByte *pabyLineDest = pabyDest + iY * nLineSpace;
291 0 : for (int iX = 0; iX < nBufXSize; iX++)
292 : {
293 0 : *pabyLineDest = (panSrc[i] == nNoData) ? 0 : 255;
294 0 : ++i;
295 0 : pabyLineDest += nPixelSpace;
296 : }
297 : }
298 : }
299 1782 : }
300 :
301 : /************************************************************************/
302 : /* IRasterIO() */
303 : /************************************************************************/
304 :
305 9525 : CPLErr GDALNoDataMaskBand::IRasterIO(GDALRWFlag eRWFlag, int nXOff, int nYOff,
306 : int nXSize, int nYSize, void *pData,
307 : int nBufXSize, int nBufYSize,
308 : GDALDataType eBufType,
309 : GSpacing nPixelSpace, GSpacing nLineSpace,
310 : GDALRasterIOExtraArg *psExtraArg)
311 : {
312 9525 : if (eRWFlag != GF_Read)
313 : {
314 0 : return CE_Failure;
315 : }
316 9525 : const auto eParentDT = m_poParent->GetRasterDataType();
317 9525 : const GDALDataType eWrkDT = GetWorkDataType(eParentDT);
318 :
319 : // Optimization in common use case (#4488).
320 : // This avoids triggering the block cache on this band, which helps
321 : // reducing the global block cache consumption.
322 9525 : if (eBufType == GDT_Byte && eWrkDT == GDT_Byte)
323 : {
324 388 : const CPLErr eErr = m_poParent->RasterIO(
325 : GF_Read, nXOff, nYOff, nXSize, nYSize, pData, nBufXSize, nBufYSize,
326 : eBufType, nPixelSpace, nLineSpace, psExtraArg);
327 388 : if (eErr != CE_None)
328 0 : return eErr;
329 :
330 388 : GByte *pabyData = static_cast<GByte *>(pData);
331 388 : const GByte byNoData = static_cast<GByte>(m_dfNoDataValue);
332 :
333 388 : if (nPixelSpace == 1 && nLineSpace == nBufXSize)
334 : {
335 388 : const size_t nBufSize = static_cast<size_t>(nBufXSize) * nBufYSize;
336 388 : SetZeroOr255(pabyData, nBufSize, byNoData);
337 : }
338 : else
339 : {
340 0 : SetZeroOr255(pabyData, pabyData, nBufXSize, nBufYSize, nPixelSpace,
341 : nLineSpace, byNoData);
342 : }
343 388 : return CE_None;
344 : }
345 :
346 : const auto AllocTempBufferOrFallback =
347 9137 : [this, eRWFlag, nXOff, nYOff, nXSize, nYSize, pData, nBufXSize,
348 : nBufYSize, eBufType, nPixelSpace, nLineSpace,
349 18330 : psExtraArg](int nWrkDTSize) -> std::pair<CPLErr, void *>
350 : {
351 9137 : auto poParentDS = m_poParent->GetDataset();
352 : // Check if we must simulate a memory allocation failure
353 : // Before checking the env variable, which is slightly expensive,
354 : // check first for a special dataset name, which is a cheap test.
355 : const char *pszOptVal =
356 9137 : poParentDS && strcmp(poParentDS->GetDescription(), "__debug__") == 0
357 18274 : ? CPLGetConfigOption(
358 : "GDAL_SIMUL_MEM_ALLOC_FAILURE_NODATA_MASK_BAND", "NO")
359 9137 : : "NO";
360 : const bool bSimulMemAllocFailure =
361 18270 : EQUAL(pszOptVal, "ALWAYS") ||
362 9133 : (CPLTestBool(pszOptVal) &&
363 14 : GDALMajorObject::GetMetadataItem(__func__, "__INTERNAL__") ==
364 9137 : nullptr);
365 9137 : void *pTemp = nullptr;
366 9137 : if (!bSimulMemAllocFailure)
367 : {
368 9125 : CPLErrorStateBackuper oErrorStateBackuper(CPLQuietErrorHandler);
369 9125 : pTemp = VSI_MALLOC3_VERBOSE(nWrkDTSize, nBufXSize, nBufYSize);
370 : }
371 9137 : if (!pTemp)
372 : {
373 : const bool bAllocHasAlreadyFailed =
374 12 : GDALMajorObject::GetMetadataItem(__func__, "__INTERNAL__") !=
375 12 : nullptr;
376 12 : CPLError(bAllocHasAlreadyFailed ? CE_Failure : CE_Warning,
377 : CPLE_OutOfMemory,
378 : "GDALNoDataMaskBand::IRasterIO(): cannot allocate %d x %d "
379 : "x %d bytes%s",
380 : nBufXSize, nBufYSize, nWrkDTSize,
381 : bAllocHasAlreadyFailed
382 : ? ""
383 : : ". Falling back to block-based approach");
384 12 : if (bAllocHasAlreadyFailed)
385 2 : return std::pair(CE_Failure, nullptr);
386 : // Sets a metadata item to prevent potential infinite recursion
387 10 : GDALMajorObject::SetMetadataItem(__func__, "IN", "__INTERNAL__");
388 10 : const CPLErr eErr = GDALRasterBand::IRasterIO(
389 : eRWFlag, nXOff, nYOff, nXSize, nYSize, pData, nBufXSize,
390 10 : nBufYSize, eBufType, nPixelSpace, nLineSpace, psExtraArg);
391 10 : GDALMajorObject::SetMetadataItem(__func__, nullptr, "__INTERNAL__");
392 10 : return std::pair(eErr, nullptr);
393 : }
394 9125 : return std::pair(CE_None, pTemp);
395 9137 : };
396 :
397 9137 : if (eBufType == GDT_Byte)
398 : {
399 9058 : const int nWrkDTSize = GDALGetDataTypeSizeBytes(eWrkDT);
400 9058 : auto [eErr, pTemp] = AllocTempBufferOrFallback(nWrkDTSize);
401 9058 : if (!pTemp)
402 12 : return eErr;
403 :
404 18092 : eErr = m_poParent->RasterIO(
405 : GF_Read, nXOff, nYOff, nXSize, nYSize, pTemp, nBufXSize, nBufYSize,
406 9046 : eWrkDT, nWrkDTSize, static_cast<GSpacing>(nBufXSize) * nWrkDTSize,
407 : psExtraArg);
408 9046 : if (eErr != CE_None)
409 : {
410 0 : VSIFree(pTemp);
411 0 : return eErr;
412 : }
413 :
414 9046 : const bool bIsNoDataNan = std::isnan(m_dfNoDataValue) != 0;
415 9046 : GByte *pabyDest = static_cast<GByte *>(pData);
416 :
417 : /* --------------------------------------------------------------------
418 : */
419 : /* Process different cases. */
420 : /* --------------------------------------------------------------------
421 : */
422 9046 : switch (eWrkDT)
423 : {
424 1676 : case GDT_Int16:
425 : {
426 1676 : const auto nNoData = static_cast<int16_t>(m_dfNoDataValue);
427 1676 : const auto *panSrc = static_cast<const int16_t *>(pTemp);
428 1676 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nBufYSize,
429 : nPixelSpace, nLineSpace, nNoData);
430 : }
431 1676 : break;
432 :
433 8 : case GDT_UInt16:
434 : {
435 8 : const auto nNoData = static_cast<uint16_t>(m_dfNoDataValue);
436 8 : const auto *panSrc = static_cast<const uint16_t *>(pTemp);
437 8 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nBufYSize,
438 : nPixelSpace, nLineSpace, nNoData);
439 : }
440 8 : break;
441 :
442 7 : case GDT_UInt32:
443 : {
444 7 : const auto nNoData = static_cast<GUInt32>(m_dfNoDataValue);
445 7 : const auto *panSrc = static_cast<const GUInt32 *>(pTemp);
446 7 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nBufYSize,
447 : nPixelSpace, nLineSpace, nNoData);
448 : }
449 7 : break;
450 :
451 75 : case GDT_Int32:
452 : {
453 75 : const auto nNoData = static_cast<GInt32>(m_dfNoDataValue);
454 75 : const auto *panSrc = static_cast<const GInt32 *>(pTemp);
455 75 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nBufYSize,
456 : nPixelSpace, nLineSpace, nNoData);
457 : }
458 75 : break;
459 :
460 7232 : case GDT_Float32:
461 : {
462 7232 : const float fNoData = static_cast<float>(m_dfNoDataValue);
463 7232 : const float *pafSrc = static_cast<const float *>(pTemp);
464 :
465 7232 : size_t i = 0;
466 14540 : for (int iY = 0; iY < nBufYSize; iY++)
467 : {
468 7308 : GByte *pabyLineDest = pabyDest + iY * nLineSpace;
469 1962170 : for (int iX = 0; iX < nBufXSize; iX++)
470 : {
471 1954860 : const float fVal = pafSrc[i];
472 1954860 : if (bIsNoDataNan && std::isnan(fVal))
473 26 : *pabyLineDest = 0;
474 1954840 : else if (ARE_REAL_EQUAL(fVal, fNoData))
475 876112 : *pabyLineDest = 0;
476 : else
477 1078720 : *pabyLineDest = 255;
478 1954860 : ++i;
479 1954860 : pabyLineDest += nPixelSpace;
480 : }
481 : }
482 : }
483 7232 : break;
484 :
485 32 : case GDT_Float64:
486 : {
487 32 : const double *padfSrc = static_cast<const double *>(pTemp);
488 :
489 32 : size_t i = 0;
490 112 : for (int iY = 0; iY < nBufYSize; iY++)
491 : {
492 80 : GByte *pabyLineDest = pabyDest + iY * nLineSpace;
493 317 : for (int iX = 0; iX < nBufXSize; iX++)
494 : {
495 237 : const double dfVal = padfSrc[i];
496 237 : if (bIsNoDataNan && std::isnan(dfVal))
497 32 : *pabyLineDest = 0;
498 205 : else if (ARE_REAL_EQUAL(dfVal, m_dfNoDataValue))
499 82 : *pabyLineDest = 0;
500 : else
501 123 : *pabyLineDest = 255;
502 237 : ++i;
503 237 : pabyLineDest += nPixelSpace;
504 : }
505 : }
506 : }
507 32 : break;
508 :
509 8 : case GDT_Int64:
510 : {
511 8 : const auto *panSrc = static_cast<const int64_t *>(pTemp);
512 8 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nBufYSize,
513 : nPixelSpace, nLineSpace, m_nNoDataValueInt64);
514 : }
515 8 : break;
516 :
517 8 : case GDT_UInt64:
518 : {
519 8 : const auto *panSrc = static_cast<const uint64_t *>(pTemp);
520 8 : SetZeroOr255(pabyDest, panSrc, nBufXSize, nBufYSize,
521 : nPixelSpace, nLineSpace, m_nNoDataValueUInt64);
522 : }
523 8 : break;
524 :
525 0 : default:
526 0 : CPLAssert(false);
527 : break;
528 : }
529 :
530 9046 : VSIFree(pTemp);
531 9046 : return CE_None;
532 : }
533 :
534 : // Output buffer is non-Byte. Ask for Byte and expand to user requested
535 : // type
536 79 : auto [eErr, pTemp] = AllocTempBufferOrFallback(sizeof(GByte));
537 79 : if (!pTemp)
538 0 : return eErr;
539 :
540 158 : eErr = IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize, pTemp, nBufXSize,
541 79 : nBufYSize, GDT_Byte, 1, nBufXSize, psExtraArg);
542 79 : if (eErr != CE_None)
543 : {
544 0 : VSIFree(pTemp);
545 0 : return eErr;
546 : }
547 :
548 186 : for (int iY = 0; iY < nBufYSize; iY++)
549 : {
550 107 : GDALCopyWords64(
551 107 : static_cast<GByte *>(pTemp) + static_cast<size_t>(iY) * nBufXSize,
552 107 : GDT_Byte, 1, static_cast<GByte *>(pData) + iY * nLineSpace,
553 : eBufType, static_cast<int>(nPixelSpace), nBufXSize);
554 : }
555 79 : VSIFree(pTemp);
556 79 : return CE_None;
557 : }
558 :
559 : //! @endcond
|
