Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: Virtual GDAL Datasets
4 : * Purpose: Implementation of VRTSimpleSource, VRTFuncSource and
5 : * VRTAveragedSource.
6 : * Author: Frank Warmerdam <warmerdam@pobox.com>
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 2001, Frank Warmerdam <warmerdam@pobox.com>
10 : * Copyright (c) 2008-2013, Even Rouault <even dot rouault at spatialys.com>
11 : *
12 : * SPDX-License-Identifier: MIT
13 : ****************************************************************************/
14 :
15 : #include "gdal_vrt.h"
16 : #include "vrtdataset.h"
17 :
18 : #include <cassert>
19 : #include <climits>
20 : #include <cmath>
21 : #include <cstddef>
22 : #include <cstdio>
23 : #include <cstdlib>
24 : #include <cstring>
25 : #include <algorithm>
26 : #include <limits>
27 : #include <string>
28 :
29 : #include "cpl_conv.h"
30 : #include "cpl_error.h"
31 : #include "cpl_hash_set.h"
32 : #include "cpl_minixml.h"
33 : #include "cpl_progress.h"
34 : #include "cpl_string.h"
35 : #include "cpl_vsi.h"
36 : #include "gdal.h"
37 : #include "gdal_priv.h"
38 : #include "gdal_proxy.h"
39 : #include "gdal_priv_templates.hpp"
40 :
41 : /*! @cond Doxygen_Suppress */
42 :
43 : // #define DEBUG_VERBOSE 1
44 :
45 : // See #5459
46 : #ifdef isnan
47 : #define HAS_ISNAN_MACRO
48 : #endif
49 : #include <algorithm>
50 : #if defined(HAS_ISNAN_MACRO) && !defined(isnan)
51 : #define isnan std::isnan
52 : #endif
53 :
54 : /************************************************************************/
55 : /* ==================================================================== */
56 : /* VRTSource */
57 : /* ==================================================================== */
58 : /************************************************************************/
59 :
60 13333 : VRTSource::~VRTSource()
61 : {
62 13333 : }
63 :
64 : /************************************************************************/
65 : /* GetFileList() */
66 : /************************************************************************/
67 :
68 0 : void VRTSource::GetFileList(char *** /* ppapszFileList */, int * /* pnSize */,
69 : int * /* pnMaxSize */, CPLHashSet * /* hSetFiles */)
70 : {
71 0 : }
72 :
73 : /************************************************************************/
74 : /* ==================================================================== */
75 : /* VRTSimpleSource */
76 : /* ==================================================================== */
77 : /************************************************************************/
78 :
79 : /************************************************************************/
80 : /* VRTSimpleSource() */
81 : /************************************************************************/
82 :
83 : VRTSimpleSource::VRTSimpleSource() = default;
84 :
85 : /************************************************************************/
86 : /* VRTSimpleSource() */
87 : /************************************************************************/
88 :
89 76 : VRTSimpleSource::VRTSimpleSource(const VRTSimpleSource *poSrcSource,
90 76 : double dfXDstRatio, double dfYDstRatio)
91 76 : : m_poMapSharedSources(poSrcSource->m_poMapSharedSources),
92 76 : m_poRasterBand(poSrcSource->m_poRasterBand),
93 76 : m_poMaskBandMainBand(poSrcSource->m_poMaskBandMainBand),
94 76 : m_aosOpenOptionsOri(poSrcSource->m_aosOpenOptionsOri),
95 76 : m_aosOpenOptions(poSrcSource->m_aosOpenOptions),
96 76 : m_bSrcDSNameFromVRT(poSrcSource->m_bSrcDSNameFromVRT),
97 76 : m_nBand(poSrcSource->m_nBand),
98 76 : m_bGetMaskBand(poSrcSource->m_bGetMaskBand),
99 76 : m_dfSrcXOff(poSrcSource->m_dfSrcXOff),
100 76 : m_dfSrcYOff(poSrcSource->m_dfSrcYOff),
101 76 : m_dfSrcXSize(poSrcSource->m_dfSrcXSize),
102 76 : m_dfSrcYSize(poSrcSource->m_dfSrcYSize),
103 76 : m_nMaxValue(poSrcSource->m_nMaxValue), m_bRelativeToVRTOri(-1),
104 76 : m_nExplicitSharedStatus(poSrcSource->m_nExplicitSharedStatus),
105 76 : m_osSrcDSName(poSrcSource->m_osSrcDSName),
106 76 : m_bDropRefOnSrcBand(poSrcSource->m_bDropRefOnSrcBand)
107 : {
108 76 : if (!poSrcSource->IsSrcWinSet() && !poSrcSource->IsDstWinSet() &&
109 0 : (dfXDstRatio != 1.0 || dfYDstRatio != 1.0))
110 : {
111 1 : auto l_band = GetRasterBand();
112 1 : if (l_band)
113 : {
114 1 : m_dfSrcXOff = 0;
115 1 : m_dfSrcYOff = 0;
116 1 : m_dfSrcXSize = l_band->GetXSize();
117 1 : m_dfSrcYSize = l_band->GetYSize();
118 1 : m_dfDstXOff = 0;
119 1 : m_dfDstYOff = 0;
120 1 : m_dfDstXSize = l_band->GetXSize() * dfXDstRatio;
121 1 : m_dfDstYSize = l_band->GetYSize() * dfYDstRatio;
122 : }
123 : }
124 75 : else if (poSrcSource->IsDstWinSet())
125 : {
126 75 : m_dfDstXOff = poSrcSource->m_dfDstXOff * dfXDstRatio;
127 75 : m_dfDstYOff = poSrcSource->m_dfDstYOff * dfYDstRatio;
128 75 : m_dfDstXSize = poSrcSource->m_dfDstXSize * dfXDstRatio;
129 75 : m_dfDstYSize = poSrcSource->m_dfDstYSize * dfYDstRatio;
130 : }
131 76 : }
132 :
133 : /************************************************************************/
134 : /* ~VRTSimpleSource() */
135 : /************************************************************************/
136 :
137 27942 : VRTSimpleSource::~VRTSimpleSource()
138 :
139 : {
140 13289 : if (!m_bDropRefOnSrcBand)
141 396 : return;
142 :
143 12893 : if (m_poMaskBandMainBand != nullptr)
144 : {
145 51 : if (m_poMaskBandMainBand->GetDataset() != nullptr)
146 : {
147 51 : m_poMaskBandMainBand->GetDataset()->ReleaseRef();
148 : }
149 : }
150 23144 : else if (m_poRasterBand != nullptr &&
151 10302 : m_poRasterBand->GetDataset() != nullptr)
152 : {
153 10302 : m_poRasterBand->GetDataset()->ReleaseRef();
154 : }
155 25962 : }
156 :
157 : /************************************************************************/
158 : /* GetTypeStatic() */
159 : /************************************************************************/
160 :
161 35433 : const char *VRTSimpleSource::GetTypeStatic()
162 : {
163 : static const char *TYPE = "SimpleSource";
164 35433 : return TYPE;
165 : }
166 :
167 : /************************************************************************/
168 : /* GetType() */
169 : /************************************************************************/
170 :
171 16123 : const char *VRTSimpleSource::GetType() const
172 : {
173 16123 : return GetTypeStatic();
174 : }
175 :
176 : /************************************************************************/
177 : /* FlushCache() */
178 : /************************************************************************/
179 :
180 13349 : CPLErr VRTSimpleSource::FlushCache(bool bAtClosing)
181 :
182 : {
183 13349 : if (m_poMaskBandMainBand != nullptr)
184 : {
185 5 : return m_poMaskBandMainBand->FlushCache(bAtClosing);
186 : }
187 13344 : else if (m_poRasterBand != nullptr)
188 : {
189 10793 : return m_poRasterBand->FlushCache(bAtClosing);
190 : }
191 2551 : return CE_None;
192 : }
193 :
194 : /************************************************************************/
195 : /* UnsetPreservedRelativeFilenames() */
196 : /************************************************************************/
197 :
198 5 : void VRTSimpleSource::UnsetPreservedRelativeFilenames()
199 : {
200 9 : if (!STARTS_WITH(m_osSourceFileNameOri.c_str(), "http://") &&
201 4 : !STARTS_WITH(m_osSourceFileNameOri.c_str(), "https://"))
202 : {
203 4 : m_bRelativeToVRTOri = -1;
204 4 : m_osSourceFileNameOri = "";
205 : }
206 5 : }
207 :
208 : /************************************************************************/
209 : /* SetSrcBand() */
210 : /************************************************************************/
211 :
212 446 : void VRTSimpleSource::SetSrcBand(const char *pszFilename, int nBand)
213 :
214 : {
215 446 : m_nBand = nBand;
216 446 : m_osSrcDSName = pszFilename;
217 446 : }
218 :
219 : /************************************************************************/
220 : /* SetSrcBand() */
221 : /************************************************************************/
222 :
223 9292 : void VRTSimpleSource::SetSrcBand(GDALRasterBand *poNewSrcBand)
224 :
225 : {
226 9292 : m_poRasterBand = poNewSrcBand;
227 9292 : m_nBand = m_poRasterBand->GetBand();
228 9292 : auto poDS = poNewSrcBand->GetDataset();
229 9292 : if (poDS != nullptr)
230 : {
231 9292 : m_osSrcDSName = poDS->GetDescription();
232 9292 : m_aosOpenOptions = CSLDuplicate(poDS->GetOpenOptions());
233 9292 : m_aosOpenOptionsOri = m_aosOpenOptions;
234 : }
235 9292 : }
236 :
237 : /************************************************************************/
238 : /* SetSrcMaskBand() */
239 : /************************************************************************/
240 :
241 : // poSrcBand is not the mask band, but the band from which the mask band is
242 : // taken.
243 32 : void VRTSimpleSource::SetSrcMaskBand(GDALRasterBand *poNewSrcBand)
244 :
245 : {
246 32 : m_poRasterBand = poNewSrcBand->GetMaskBand();
247 32 : m_poMaskBandMainBand = poNewSrcBand;
248 32 : m_nBand = poNewSrcBand->GetBand();
249 32 : auto poDS = poNewSrcBand->GetDataset();
250 32 : if (poDS != nullptr)
251 : {
252 32 : m_osSrcDSName = poDS->GetDescription();
253 32 : m_aosOpenOptions = CSLDuplicate(poDS->GetOpenOptions());
254 32 : m_aosOpenOptionsOri = m_aosOpenOptions;
255 : }
256 32 : m_bGetMaskBand = true;
257 32 : }
258 :
259 : /************************************************************************/
260 : /* RoundIfCloseToInt() */
261 : /************************************************************************/
262 :
263 298566 : static double RoundIfCloseToInt(double dfValue)
264 : {
265 298566 : double dfClosestInt = floor(dfValue + 0.5);
266 298566 : return (fabs(dfValue - dfClosestInt) < 1e-3) ? dfClosestInt : dfValue;
267 : }
268 :
269 : /************************************************************************/
270 : /* SetSrcWindow() */
271 : /************************************************************************/
272 :
273 12667 : void VRTSimpleSource::SetSrcWindow(double dfNewXOff, double dfNewYOff,
274 : double dfNewXSize, double dfNewYSize)
275 :
276 : {
277 12667 : m_dfSrcXOff = RoundIfCloseToInt(dfNewXOff);
278 12667 : m_dfSrcYOff = RoundIfCloseToInt(dfNewYOff);
279 12667 : m_dfSrcXSize = RoundIfCloseToInt(dfNewXSize);
280 12667 : m_dfSrcYSize = RoundIfCloseToInt(dfNewYSize);
281 12667 : }
282 :
283 : /************************************************************************/
284 : /* SetDstWindow() */
285 : /************************************************************************/
286 :
287 12666 : void VRTSimpleSource::SetDstWindow(double dfNewXOff, double dfNewYOff,
288 : double dfNewXSize, double dfNewYSize)
289 :
290 : {
291 12666 : m_dfDstXOff = RoundIfCloseToInt(dfNewXOff);
292 12666 : m_dfDstYOff = RoundIfCloseToInt(dfNewYOff);
293 12666 : m_dfDstXSize = RoundIfCloseToInt(dfNewXSize);
294 12666 : m_dfDstYSize = RoundIfCloseToInt(dfNewYSize);
295 12666 : }
296 :
297 : /************************************************************************/
298 : /* GetDstWindow() */
299 : /************************************************************************/
300 :
301 517 : void VRTSimpleSource::GetDstWindow(double &dfDstXOff, double &dfDstYOff,
302 : double &dfDstXSize, double &dfDstYSize) const
303 : {
304 517 : dfDstXOff = m_dfDstXOff;
305 517 : dfDstYOff = m_dfDstYOff;
306 517 : dfDstXSize = m_dfDstXSize;
307 517 : dfDstYSize = m_dfDstYSize;
308 517 : }
309 :
310 : /************************************************************************/
311 : /* DstWindowIntersects() */
312 : /************************************************************************/
313 :
314 1073 : bool VRTSimpleSource::DstWindowIntersects(double dfXOff, double dfYOff,
315 : double dfXSize, double dfYSize) const
316 : {
317 2137 : return IsDstWinSet() && m_dfDstXOff + m_dfDstXSize > dfXOff &&
318 1064 : m_dfDstYOff + m_dfDstYSize > dfYOff &&
319 2137 : m_dfDstXOff < dfXOff + dfXSize && m_dfDstYOff < dfYOff + dfYSize;
320 : }
321 :
322 : /************************************************************************/
323 : /* IsSlowSource() */
324 : /************************************************************************/
325 :
326 2505 : static bool IsSlowSource(const char *pszSrcName)
327 : {
328 5010 : return strstr(pszSrcName, "/vsicurl/http") != nullptr ||
329 5010 : strstr(pszSrcName, "/vsicurl/ftp") != nullptr ||
330 2505 : (strstr(pszSrcName, "/vsicurl?") != nullptr &&
331 2505 : strstr(pszSrcName, "&url=http") != nullptr);
332 : }
333 :
334 : /************************************************************************/
335 : /* AddSourceFilenameNode() */
336 : /************************************************************************/
337 :
338 2537 : void VRTSimpleSource::AddSourceFilenameNode(const char *pszVRTPath,
339 : CPLXMLNode *psSrc)
340 : {
341 :
342 : VSIStatBufL sStat;
343 2537 : int bRelativeToVRT = FALSE; // TODO(schwehr): Make this a bool?
344 5074 : std::string osSourceFilename;
345 :
346 2537 : if (m_bRelativeToVRTOri >= 0)
347 : {
348 32 : osSourceFilename = m_osSourceFileNameOri;
349 32 : bRelativeToVRT = m_bRelativeToVRTOri;
350 : }
351 2505 : else if (IsSlowSource(m_osSrcDSName))
352 : {
353 : // Testing the existence of remote resources can be excruciating
354 : // slow, so let's just suppose they exist.
355 0 : osSourceFilename = m_osSrcDSName;
356 0 : bRelativeToVRT = FALSE;
357 : }
358 : // If this isn't actually a file, don't even try to know if it is a
359 : // relative path. It can't be !, and unfortunately CPLIsFilenameRelative()
360 : // can only work with strings that are filenames To be clear
361 : // NITF_TOC_ENTRY:CADRG_JOG-A_250K_1_0:some_path isn't a relative file
362 : // path.
363 2505 : else if (VSIStatExL(m_osSrcDSName, &sStat, VSI_STAT_EXISTS_FLAG) != 0)
364 : {
365 58 : osSourceFilename = m_osSrcDSName;
366 58 : bRelativeToVRT = FALSE;
367 :
368 : // Try subdatasetinfo API first
369 : // Note: this will become the only branch when subdatasetinfo will become
370 : // available for NITF_IM, RASTERLITE and TILEDB
371 58 : const auto oSubDSInfo{GDALGetSubdatasetInfo(osSourceFilename.c_str())};
372 58 : if (oSubDSInfo && !oSubDSInfo->GetPathComponent().empty())
373 : {
374 10 : auto path{oSubDSInfo->GetPathComponent()};
375 : std::string relPath{CPLExtractRelativePath(pszVRTPath, path.c_str(),
376 10 : &bRelativeToVRT)};
377 5 : osSourceFilename = oSubDSInfo->ModifyPathComponent(relPath);
378 5 : GDALDestroySubdatasetInfo(oSubDSInfo);
379 : }
380 : else
381 : {
382 303 : for (const char *pszSyntax : VRTDataset::apszSpecialSyntax)
383 : {
384 253 : CPLString osPrefix(pszSyntax);
385 253 : osPrefix.resize(strchr(pszSyntax, ':') - pszSyntax + 1);
386 253 : if (pszSyntax[osPrefix.size()] == '"')
387 50 : osPrefix += '"';
388 253 : if (EQUALN(osSourceFilename.c_str(), osPrefix, osPrefix.size()))
389 : {
390 3 : if (STARTS_WITH_CI(pszSyntax + osPrefix.size(), "{ANY}"))
391 : {
392 : const char *pszLastPart =
393 3 : strrchr(osSourceFilename.c_str(), ':') + 1;
394 : // CSV:z:/foo.xyz
395 1 : if ((pszLastPart[0] == '/' || pszLastPart[0] == '\\') &&
396 6 : pszLastPart - osSourceFilename.c_str() >= 3 &&
397 2 : pszLastPart[-3] == ':')
398 0 : pszLastPart -= 2;
399 3 : CPLString osPrefixFilename(osSourceFilename);
400 3 : osPrefixFilename.resize(pszLastPart -
401 3 : osSourceFilename.c_str());
402 : osSourceFilename = CPLExtractRelativePath(
403 3 : pszVRTPath, pszLastPart, &bRelativeToVRT);
404 3 : osSourceFilename = osPrefixFilename + osSourceFilename;
405 : }
406 0 : else if (STARTS_WITH_CI(pszSyntax + osPrefix.size(),
407 : "{FILENAME}"))
408 : {
409 0 : CPLString osFilename(osSourceFilename.c_str() +
410 0 : osPrefix.size());
411 0 : size_t nPos = 0;
412 0 : if (osFilename.size() >= 3 && osFilename[1] == ':' &&
413 0 : (osFilename[2] == '\\' || osFilename[2] == '/'))
414 0 : nPos = 2;
415 0 : nPos = osFilename.find(
416 0 : pszSyntax[osPrefix.size() + strlen("{FILENAME}")],
417 : nPos);
418 0 : if (nPos != std::string::npos)
419 : {
420 0 : const CPLString osSuffix = osFilename.substr(nPos);
421 0 : osFilename.resize(nPos);
422 : osSourceFilename = CPLExtractRelativePath(
423 0 : pszVRTPath, osFilename, &bRelativeToVRT);
424 : osSourceFilename =
425 0 : osPrefix + osSourceFilename + osSuffix;
426 : }
427 : }
428 3 : break;
429 : }
430 : }
431 : }
432 : }
433 : else
434 : {
435 4894 : std::string osVRTFilename = pszVRTPath;
436 4894 : std::string osSourceDataset = m_osSrcDSName;
437 2447 : char *pszCurDir = CPLGetCurrentDir();
438 2447 : if (CPLIsFilenameRelative(osSourceDataset.c_str()) &&
439 2447 : !CPLIsFilenameRelative(osVRTFilename.c_str()) &&
440 : pszCurDir != nullptr)
441 : {
442 128 : osSourceDataset = CPLFormFilenameSafe(
443 64 : pszCurDir, osSourceDataset.c_str(), nullptr);
444 : }
445 2383 : else if (!CPLIsFilenameRelative(osSourceDataset.c_str()) &&
446 2383 : CPLIsFilenameRelative(osVRTFilename.c_str()) &&
447 : pszCurDir != nullptr)
448 : {
449 : osVRTFilename =
450 94 : CPLFormFilenameSafe(pszCurDir, osVRTFilename.c_str(), nullptr);
451 : }
452 2447 : CPLFree(pszCurDir);
453 : osSourceFilename = CPLExtractRelativePath(
454 2447 : osVRTFilename.c_str(), osSourceDataset.c_str(), &bRelativeToVRT);
455 : }
456 :
457 2537 : CPLSetXMLValue(psSrc, "SourceFilename", osSourceFilename.c_str());
458 :
459 2537 : CPLCreateXMLNode(CPLCreateXMLNode(CPLGetXMLNode(psSrc, "SourceFilename"),
460 : CXT_Attribute, "relativeToVRT"),
461 2537 : CXT_Text, bRelativeToVRT ? "1" : "0");
462 :
463 : // Determine if we must write the shared attribute. The config option
464 : // will override the m_nExplicitSharedStatus value
465 2537 : const char *pszShared = CPLGetConfigOption("VRT_SHARED_SOURCE", nullptr);
466 2537 : if ((pszShared == nullptr && m_nExplicitSharedStatus == 0) ||
467 0 : (pszShared != nullptr && !CPLTestBool(pszShared)))
468 : {
469 0 : CPLCreateXMLNode(
470 : CPLCreateXMLNode(CPLGetXMLNode(psSrc, "SourceFilename"),
471 : CXT_Attribute, "shared"),
472 : CXT_Text, "0");
473 : }
474 2537 : }
475 :
476 : /************************************************************************/
477 : /* SerializeToXML() */
478 : /************************************************************************/
479 :
480 2538 : CPLXMLNode *VRTSimpleSource::SerializeToXML(const char *pszVRTPath)
481 :
482 : {
483 : CPLXMLNode *const psSrc =
484 2538 : CPLCreateXMLNode(nullptr, CXT_Element, GetTypeStatic());
485 :
486 2538 : if (!m_osResampling.empty())
487 : {
488 24 : CPLCreateXMLNode(CPLCreateXMLNode(psSrc, CXT_Attribute, "resampling"),
489 : CXT_Text, m_osResampling.c_str());
490 : }
491 :
492 2538 : if (m_bSrcDSNameFromVRT)
493 : {
494 1 : CPLAddXMLChild(psSrc, CPLParseXMLString(m_osSrcDSName.c_str()));
495 : }
496 : else
497 : {
498 2537 : AddSourceFilenameNode(pszVRTPath, psSrc);
499 : }
500 :
501 2538 : GDALSerializeOpenOptionsToXML(psSrc, m_aosOpenOptionsOri.List());
502 :
503 2538 : if (m_bGetMaskBand)
504 18 : CPLSetXMLValue(psSrc, "SourceBand", CPLSPrintf("mask,%d", m_nBand));
505 : else
506 2520 : CPLSetXMLValue(psSrc, "SourceBand", CPLSPrintf("%d", m_nBand));
507 :
508 : // TODO: in a later version, no longer emit SourceProperties, which
509 : // is no longer used by GDAL 3.4
510 2538 : if (m_poRasterBand)
511 : {
512 : /* Write a few additional useful properties of the dataset */
513 : /* so that we can use a proxy dataset when re-opening. See XMLInit() */
514 : /* below */
515 2430 : CPLSetXMLValue(psSrc, "SourceProperties.#RasterXSize",
516 2430 : CPLSPrintf("%d", m_poRasterBand->GetXSize()));
517 2430 : CPLSetXMLValue(psSrc, "SourceProperties.#RasterYSize",
518 2430 : CPLSPrintf("%d", m_poRasterBand->GetYSize()));
519 2430 : CPLSetXMLValue(
520 : psSrc, "SourceProperties.#DataType",
521 2430 : GDALGetDataTypeName(m_poRasterBand->GetRasterDataType()));
522 :
523 2430 : int nBlockXSize = 0;
524 2430 : int nBlockYSize = 0;
525 2430 : m_poRasterBand->GetBlockSize(&nBlockXSize, &nBlockYSize);
526 :
527 2430 : CPLSetXMLValue(psSrc, "SourceProperties.#BlockXSize",
528 : CPLSPrintf("%d", nBlockXSize));
529 2430 : CPLSetXMLValue(psSrc, "SourceProperties.#BlockYSize",
530 : CPLSPrintf("%d", nBlockYSize));
531 : }
532 :
533 2538 : if (IsSrcWinSet())
534 : {
535 2525 : CPLSetXMLValue(psSrc, "SrcRect.#xOff",
536 : CPLSPrintf("%.15g", m_dfSrcXOff));
537 2525 : CPLSetXMLValue(psSrc, "SrcRect.#yOff",
538 : CPLSPrintf("%.15g", m_dfSrcYOff));
539 2525 : CPLSetXMLValue(psSrc, "SrcRect.#xSize",
540 : CPLSPrintf("%.15g", m_dfSrcXSize));
541 2525 : CPLSetXMLValue(psSrc, "SrcRect.#ySize",
542 : CPLSPrintf("%.15g", m_dfSrcYSize));
543 : }
544 :
545 2538 : if (IsDstWinSet())
546 : {
547 2525 : CPLSetXMLValue(psSrc, "DstRect.#xOff",
548 : CPLSPrintf("%.15g", m_dfDstXOff));
549 2525 : CPLSetXMLValue(psSrc, "DstRect.#yOff",
550 : CPLSPrintf("%.15g", m_dfDstYOff));
551 2525 : CPLSetXMLValue(psSrc, "DstRect.#xSize",
552 : CPLSPrintf("%.15g", m_dfDstXSize));
553 2525 : CPLSetXMLValue(psSrc, "DstRect.#ySize",
554 : CPLSPrintf("%.15g", m_dfDstYSize));
555 : }
556 :
557 2538 : return psSrc;
558 : }
559 :
560 : /************************************************************************/
561 : /* XMLInit() */
562 : /************************************************************************/
563 :
564 3097 : CPLErr VRTSimpleSource::XMLInit(const CPLXMLNode *psSrc, const char *pszVRTPath,
565 : VRTMapSharedResources &oMapSharedSources)
566 :
567 : {
568 3097 : m_poMapSharedSources = &oMapSharedSources;
569 :
570 3097 : m_osResampling = CPLGetXMLValue(psSrc, "resampling", "");
571 3097 : m_osName = CPLGetXMLValue(psSrc, "name", "");
572 :
573 : /* -------------------------------------------------------------------- */
574 : /* Prepare filename. */
575 : /* -------------------------------------------------------------------- */
576 : const CPLXMLNode *psSourceFileNameNode =
577 3097 : CPLGetXMLNode(psSrc, "SourceFilename");
578 3097 : const CPLXMLNode *psSourceVRTDataset = CPLGetXMLNode(psSrc, "VRTDataset");
579 : const char *pszFilename =
580 3097 : psSourceFileNameNode ? CPLGetXMLValue(psSourceFileNameNode, nullptr, "")
581 3097 : : "";
582 :
583 3097 : if (pszFilename[0] == '\0' && !psSourceVRTDataset)
584 : {
585 1 : CPLError(CE_Warning, CPLE_AppDefined,
586 : "Missing <SourceFilename> or <VRTDataset> element in <%s>.",
587 1 : psSrc->pszValue);
588 1 : return CE_Failure;
589 : }
590 :
591 : // Backup original filename and relativeToVRT so as to be able to
592 : // serialize them identically again (#5985)
593 3096 : m_osSourceFileNameOri = pszFilename;
594 3096 : if (pszFilename[0])
595 : {
596 3093 : m_bRelativeToVRTOri =
597 3093 : atoi(CPLGetXMLValue(psSourceFileNameNode, "relativetoVRT", "0"));
598 : const char *pszShared =
599 3093 : CPLGetXMLValue(psSourceFileNameNode, "shared", nullptr);
600 3093 : if (pszShared == nullptr)
601 : {
602 3091 : pszShared = CPLGetConfigOption("VRT_SHARED_SOURCE", nullptr);
603 : }
604 3093 : if (pszShared != nullptr)
605 : {
606 8 : m_nExplicitSharedStatus = CPLTestBool(pszShared);
607 : }
608 :
609 6186 : m_osSrcDSName = VRTDataset::BuildSourceFilename(
610 6186 : pszFilename, pszVRTPath, CPL_TO_BOOL(m_bRelativeToVRTOri));
611 : }
612 3 : else if (psSourceVRTDataset)
613 : {
614 : CPLXMLNode sNode;
615 3 : sNode.eType = psSourceVRTDataset->eType;
616 3 : sNode.pszValue = psSourceVRTDataset->pszValue;
617 3 : sNode.psNext = nullptr;
618 3 : sNode.psChild = psSourceVRTDataset->psChild;
619 3 : char *pszXML = CPLSerializeXMLTree(&sNode);
620 3 : if (pszXML)
621 : {
622 3 : m_bSrcDSNameFromVRT = true;
623 3 : m_osSrcDSName = pszXML;
624 3 : CPLFree(pszXML);
625 : }
626 : }
627 :
628 3096 : const char *pszSourceBand = CPLGetXMLValue(psSrc, "SourceBand", "1");
629 3096 : m_bGetMaskBand = false;
630 3096 : if (STARTS_WITH_CI(pszSourceBand, "mask"))
631 : {
632 22 : m_bGetMaskBand = true;
633 22 : if (pszSourceBand[4] == ',')
634 22 : m_nBand = atoi(pszSourceBand + 5);
635 : else
636 0 : m_nBand = 1;
637 : }
638 : else
639 : {
640 3074 : m_nBand = atoi(pszSourceBand);
641 : }
642 3096 : if (!GDALCheckBandCount(m_nBand, 0))
643 : {
644 0 : CPLError(CE_Warning, CPLE_AppDefined,
645 : "Invalid <SourceBand> element in VRTRasterBand.");
646 0 : return CE_Failure;
647 : }
648 :
649 3096 : m_aosOpenOptions = GDALDeserializeOpenOptionsFromXML(psSrc);
650 3096 : m_aosOpenOptionsOri = m_aosOpenOptions;
651 3096 : if (strstr(m_osSrcDSName.c_str(), "<VRTDataset") != nullptr)
652 5 : m_aosOpenOptions.SetNameValue("ROOT_PATH", pszVRTPath);
653 :
654 3096 : return ParseSrcRectAndDstRect(psSrc);
655 : }
656 :
657 : /************************************************************************/
658 : /* ParseSrcRectAndDstRect() */
659 : /************************************************************************/
660 :
661 3110 : CPLErr VRTSimpleSource::ParseSrcRectAndDstRect(const CPLXMLNode *psSrc)
662 : {
663 23188 : const auto GetAttrValue = [](const CPLXMLNode *psNode,
664 : const char *pszAttrName, double dfDefaultVal)
665 : {
666 23188 : if (const char *pszVal = CPLGetXMLValue(psNode, pszAttrName, nullptr))
667 23188 : return CPLAtof(pszVal);
668 : else
669 0 : return dfDefaultVal;
670 : };
671 :
672 : /* -------------------------------------------------------------------- */
673 : /* Set characteristics. */
674 : /* -------------------------------------------------------------------- */
675 3110 : const CPLXMLNode *const psSrcRect = CPLGetXMLNode(psSrc, "SrcRect");
676 3110 : if (psSrcRect)
677 : {
678 2899 : double xOff = GetAttrValue(psSrcRect, "xOff", UNINIT_WINDOW);
679 2899 : double yOff = GetAttrValue(psSrcRect, "yOff", UNINIT_WINDOW);
680 2899 : double xSize = GetAttrValue(psSrcRect, "xSize", UNINIT_WINDOW);
681 2899 : double ySize = GetAttrValue(psSrcRect, "ySize", UNINIT_WINDOW);
682 : // Test written that way to catch NaN values
683 2899 : if (!(xOff >= INT_MIN && xOff <= INT_MAX) ||
684 2899 : !(yOff >= INT_MIN && yOff <= INT_MAX) ||
685 2899 : !(xSize > 0 || xSize == UNINIT_WINDOW) || xSize > INT_MAX ||
686 2898 : !(ySize > 0 || ySize == UNINIT_WINDOW) || ySize > INT_MAX)
687 : {
688 1 : CPLError(CE_Failure, CPLE_AppDefined, "Wrong values in SrcRect");
689 1 : return CE_Failure;
690 : }
691 2898 : SetSrcWindow(xOff, yOff, xSize, ySize);
692 : }
693 : else
694 : {
695 211 : m_dfSrcXOff = UNINIT_WINDOW;
696 211 : m_dfSrcYOff = UNINIT_WINDOW;
697 211 : m_dfSrcXSize = UNINIT_WINDOW;
698 211 : m_dfSrcYSize = UNINIT_WINDOW;
699 : }
700 :
701 3109 : const CPLXMLNode *const psDstRect = CPLGetXMLNode(psSrc, "DstRect");
702 3109 : if (psDstRect)
703 : {
704 2898 : double xOff = GetAttrValue(psDstRect, "xOff", UNINIT_WINDOW);
705 : ;
706 2898 : double yOff = GetAttrValue(psDstRect, "yOff", UNINIT_WINDOW);
707 2898 : double xSize = GetAttrValue(psDstRect, "xSize", UNINIT_WINDOW);
708 : ;
709 2898 : double ySize = GetAttrValue(psDstRect, "ySize", UNINIT_WINDOW);
710 : // Test written that way to catch NaN values
711 2898 : if (!(xOff >= INT_MIN && xOff <= INT_MAX) ||
712 2898 : !(yOff >= INT_MIN && yOff <= INT_MAX) ||
713 2898 : !(xSize > 0 || xSize == UNINIT_WINDOW) || xSize > INT_MAX ||
714 2898 : !(ySize > 0 || ySize == UNINIT_WINDOW) || ySize > INT_MAX)
715 : {
716 1 : CPLError(CE_Failure, CPLE_AppDefined, "Wrong values in DstRect");
717 1 : return CE_Failure;
718 : }
719 2897 : SetDstWindow(xOff, yOff, xSize, ySize);
720 : }
721 : else
722 : {
723 211 : m_dfDstXOff = UNINIT_WINDOW;
724 211 : m_dfDstYOff = UNINIT_WINDOW;
725 211 : m_dfDstXSize = UNINIT_WINDOW;
726 211 : m_dfDstYSize = UNINIT_WINDOW;
727 : }
728 :
729 3108 : return CE_None;
730 : }
731 :
732 : /************************************************************************/
733 : /* GetFileList() */
734 : /************************************************************************/
735 :
736 112 : void VRTSimpleSource::GetFileList(char ***ppapszFileList, int *pnSize,
737 : int *pnMaxSize, CPLHashSet *hSetFiles)
738 : {
739 112 : if (!m_osSrcDSName.empty() && !m_bSrcDSNameFromVRT)
740 : {
741 112 : const char *pszFilename = m_osSrcDSName.c_str();
742 :
743 : /* --------------------------------------------------------------------
744 : */
745 : /* Is it already in the list ? */
746 : /* --------------------------------------------------------------------
747 : */
748 112 : if (CPLHashSetLookup(hSetFiles, pszFilename) != nullptr)
749 15 : return;
750 :
751 : /* --------------------------------------------------------------------
752 : */
753 : /* Grow array if necessary */
754 : /* --------------------------------------------------------------------
755 : */
756 97 : if (*pnSize + 1 >= *pnMaxSize)
757 : {
758 66 : *pnMaxSize = std::max(*pnSize + 2, 2 + 2 * (*pnMaxSize));
759 66 : *ppapszFileList = static_cast<char **>(
760 66 : CPLRealloc(*ppapszFileList, sizeof(char *) * (*pnMaxSize)));
761 : }
762 :
763 : /* --------------------------------------------------------------------
764 : */
765 : /* Add the string to the list */
766 : /* --------------------------------------------------------------------
767 : */
768 97 : (*ppapszFileList)[*pnSize] = CPLStrdup(pszFilename);
769 97 : (*ppapszFileList)[(*pnSize + 1)] = nullptr;
770 97 : CPLHashSetInsert(hSetFiles, (*ppapszFileList)[*pnSize]);
771 :
772 97 : (*pnSize)++;
773 : }
774 : }
775 :
776 : /************************************************************************/
777 : /* OpenSource() */
778 : /************************************************************************/
779 :
780 1229 : void VRTSimpleSource::OpenSource() const
781 : {
782 1229 : CPLAssert(m_poRasterBand == nullptr);
783 :
784 : /* ----------------------------------------------------------------- */
785 : /* Create a proxy dataset */
786 : /* ----------------------------------------------------------------- */
787 1229 : GDALProxyPoolDataset *proxyDS = nullptr;
788 1229 : std::string osKeyMapSharedSources;
789 1229 : if (m_poMapSharedSources)
790 : {
791 1182 : osKeyMapSharedSources = m_osSrcDSName;
792 1191 : for (int i = 0; i < m_aosOpenOptions.size(); ++i)
793 : {
794 9 : osKeyMapSharedSources += "||";
795 9 : osKeyMapSharedSources += m_aosOpenOptions[i];
796 : }
797 :
798 1182 : proxyDS = cpl::down_cast<GDALProxyPoolDataset *>(
799 1182 : m_poMapSharedSources->Get(osKeyMapSharedSources));
800 : }
801 :
802 1229 : if (proxyDS == nullptr)
803 : {
804 1058 : int bShared = true;
805 1058 : if (m_nExplicitSharedStatus != -1)
806 8 : bShared = m_nExplicitSharedStatus;
807 :
808 1058 : const CPLString osUniqueHandle(CPLSPrintf("%p", m_poMapSharedSources));
809 1058 : proxyDS = GDALProxyPoolDataset::Create(
810 : m_osSrcDSName, m_aosOpenOptions.List(), GA_ReadOnly, bShared,
811 : osUniqueHandle.c_str());
812 1058 : if (proxyDS == nullptr)
813 213 : return;
814 : }
815 : else
816 : {
817 171 : proxyDS->Reference();
818 : }
819 :
820 1016 : if (m_bGetMaskBand)
821 : {
822 : GDALProxyPoolRasterBand *poMaskBand =
823 15 : cpl::down_cast<GDALProxyPoolRasterBand *>(
824 15 : proxyDS->GetRasterBand(m_nBand));
825 15 : poMaskBand->AddSrcMaskBandDescriptionFromUnderlying();
826 : }
827 :
828 : /* -------------------------------------------------------------------- */
829 : /* Get the raster band. */
830 : /* -------------------------------------------------------------------- */
831 :
832 1016 : m_poRasterBand = proxyDS->GetRasterBand(m_nBand);
833 1016 : if (m_poRasterBand == nullptr || !ValidateOpenedBand(m_poRasterBand))
834 : {
835 2 : proxyDS->ReleaseRef();
836 2 : return;
837 : }
838 :
839 1014 : if (m_bGetMaskBand)
840 : {
841 15 : m_poRasterBand = m_poRasterBand->GetMaskBand();
842 15 : if (m_poRasterBand == nullptr)
843 : {
844 0 : proxyDS->ReleaseRef();
845 0 : return;
846 : }
847 15 : m_poMaskBandMainBand = m_poRasterBand;
848 : }
849 :
850 1014 : if (m_poMapSharedSources)
851 : {
852 967 : m_poMapSharedSources->Insert(osKeyMapSharedSources, proxyDS);
853 : }
854 : }
855 :
856 : /************************************************************************/
857 : /* GetRasterBand() */
858 : /************************************************************************/
859 :
860 104856 : GDALRasterBand *VRTSimpleSource::GetRasterBand() const
861 : {
862 104856 : if (m_poRasterBand == nullptr)
863 1227 : OpenSource();
864 104855 : return m_poRasterBand;
865 : }
866 :
867 : /************************************************************************/
868 : /* GetMaskBandMainBand() */
869 : /************************************************************************/
870 :
871 1245 : GDALRasterBand *VRTSimpleSource::GetMaskBandMainBand()
872 : {
873 1245 : if (m_poRasterBand == nullptr)
874 2 : OpenSource();
875 1245 : return m_poMaskBandMainBand;
876 : }
877 :
878 : /************************************************************************/
879 : /* IsSameExceptBandNumber() */
880 : /************************************************************************/
881 :
882 971 : bool VRTSimpleSource::IsSameExceptBandNumber(
883 : const VRTSimpleSource *poOtherSource) const
884 : {
885 1942 : return m_dfSrcXOff == poOtherSource->m_dfSrcXOff &&
886 971 : m_dfSrcYOff == poOtherSource->m_dfSrcYOff &&
887 971 : m_dfSrcXSize == poOtherSource->m_dfSrcXSize &&
888 971 : m_dfSrcYSize == poOtherSource->m_dfSrcYSize &&
889 971 : m_dfDstXOff == poOtherSource->m_dfDstXOff &&
890 971 : m_dfDstYOff == poOtherSource->m_dfDstYOff &&
891 971 : m_dfDstXSize == poOtherSource->m_dfDstXSize &&
892 971 : m_dfDstYSize == poOtherSource->m_dfDstYSize &&
893 2913 : !m_osSrcDSName.empty() &&
894 1942 : m_osSrcDSName == poOtherSource->m_osSrcDSName;
895 : }
896 :
897 : /************************************************************************/
898 : /* SrcToDst() */
899 : /* */
900 : /* Note: this is a no-op if the both src and dst windows are unset */
901 : /************************************************************************/
902 :
903 6510 : void VRTSimpleSource::SrcToDst(double dfX, double dfY, double &dfXOut,
904 : double &dfYOut) const
905 :
906 : {
907 6510 : dfXOut = ((dfX - m_dfSrcXOff) / m_dfSrcXSize) * m_dfDstXSize + m_dfDstXOff;
908 6510 : dfYOut = ((dfY - m_dfSrcYOff) / m_dfSrcYSize) * m_dfDstYSize + m_dfDstYOff;
909 6510 : }
910 :
911 : /************************************************************************/
912 : /* DstToSrc() */
913 : /* */
914 : /* Note: this is a no-op if the both src and dst windows are unset */
915 : /************************************************************************/
916 :
917 4022310 : void VRTSimpleSource::DstToSrc(double dfX, double dfY, double &dfXOut,
918 : double &dfYOut) const
919 :
920 : {
921 4022310 : dfXOut = ((dfX - m_dfDstXOff) / m_dfDstXSize) * m_dfSrcXSize + m_dfSrcXOff;
922 4022310 : dfYOut = ((dfY - m_dfDstYOff) / m_dfDstYSize) * m_dfSrcYSize + m_dfSrcYOff;
923 4022310 : }
924 :
925 : /************************************************************************/
926 : /* GetSrcDstWindow() */
927 : /************************************************************************/
928 :
929 69734 : int VRTSimpleSource::GetSrcDstWindow(
930 : double dfXOff, double dfYOff, double dfXSize, double dfYSize, int nBufXSize,
931 : int nBufYSize, double *pdfReqXOff, double *pdfReqYOff, double *pdfReqXSize,
932 : double *pdfReqYSize, int *pnReqXOff, int *pnReqYOff, int *pnReqXSize,
933 : int *pnReqYSize, int *pnOutXOff, int *pnOutYOff, int *pnOutXSize,
934 : int *pnOutYSize, bool &bErrorOut)
935 :
936 : {
937 69734 : bErrorOut = false;
938 :
939 69734 : if (m_dfSrcXSize == 0.0 || m_dfSrcYSize == 0.0 || m_dfDstXSize == 0.0 ||
940 69735 : m_dfDstYSize == 0.0)
941 : {
942 0 : return FALSE;
943 : }
944 :
945 69735 : const bool bDstWinSet = IsDstWinSet();
946 :
947 : #ifdef DEBUG
948 69734 : const bool bSrcWinSet = IsSrcWinSet();
949 :
950 69734 : if (bSrcWinSet != bDstWinSet)
951 : {
952 0 : return FALSE;
953 : }
954 : #endif
955 :
956 : /* -------------------------------------------------------------------- */
957 : /* If the input window completely misses the portion of the */
958 : /* virtual dataset provided by this source we have nothing to do. */
959 : /* -------------------------------------------------------------------- */
960 69734 : if (bDstWinSet)
961 : {
962 69464 : if (dfXOff >= m_dfDstXOff + m_dfDstXSize ||
963 59035 : dfYOff >= m_dfDstYOff + m_dfDstYSize ||
964 58499 : dfXOff + dfXSize <= m_dfDstXOff || dfYOff + dfYSize <= m_dfDstYOff)
965 20311 : return FALSE;
966 : }
967 :
968 : /* -------------------------------------------------------------------- */
969 : /* This request window corresponds to the whole output buffer. */
970 : /* -------------------------------------------------------------------- */
971 49423 : *pnOutXOff = 0;
972 49423 : *pnOutYOff = 0;
973 49423 : *pnOutXSize = nBufXSize;
974 49423 : *pnOutYSize = nBufYSize;
975 :
976 : /* -------------------------------------------------------------------- */
977 : /* If the input window extents outside the portion of the on */
978 : /* the virtual file that this source can set, then clip down */
979 : /* the requested window. */
980 : /* -------------------------------------------------------------------- */
981 49423 : bool bModifiedX = false;
982 49423 : bool bModifiedY = false;
983 49423 : double dfRXOff = dfXOff;
984 49423 : double dfRYOff = dfYOff;
985 49423 : double dfRXSize = dfXSize;
986 49423 : double dfRYSize = dfYSize;
987 :
988 49423 : if (bDstWinSet)
989 : {
990 49153 : if (dfRXOff < m_dfDstXOff)
991 : {
992 1860 : dfRXSize = dfRXSize + dfRXOff - m_dfDstXOff;
993 1860 : dfRXOff = m_dfDstXOff;
994 1860 : bModifiedX = true;
995 : }
996 :
997 49153 : if (dfRYOff < m_dfDstYOff)
998 : {
999 166 : dfRYSize = dfRYSize + dfRYOff - m_dfDstYOff;
1000 166 : dfRYOff = m_dfDstYOff;
1001 166 : bModifiedY = true;
1002 : }
1003 :
1004 49153 : if (dfRXOff + dfRXSize > m_dfDstXOff + m_dfDstXSize)
1005 : {
1006 1935 : dfRXSize = m_dfDstXOff + m_dfDstXSize - dfRXOff;
1007 1935 : bModifiedX = true;
1008 : }
1009 :
1010 49153 : if (dfRYOff + dfRYSize > m_dfDstYOff + m_dfDstYSize)
1011 : {
1012 265 : dfRYSize = m_dfDstYOff + m_dfDstYSize - dfRYOff;
1013 265 : bModifiedY = true;
1014 : }
1015 : }
1016 :
1017 : /* -------------------------------------------------------------------- */
1018 : /* Translate requested region in virtual file into the source */
1019 : /* band coordinates. */
1020 : /* -------------------------------------------------------------------- */
1021 49423 : const double dfScaleX = m_dfSrcXSize / m_dfDstXSize;
1022 49423 : const double dfScaleY = m_dfSrcYSize / m_dfDstYSize;
1023 :
1024 49423 : *pdfReqXOff = (dfRXOff - m_dfDstXOff) * dfScaleX + m_dfSrcXOff;
1025 49423 : *pdfReqYOff = (dfRYOff - m_dfDstYOff) * dfScaleY + m_dfSrcYOff;
1026 49423 : *pdfReqXSize = dfRXSize * dfScaleX;
1027 49423 : *pdfReqYSize = dfRYSize * dfScaleY;
1028 :
1029 98846 : if (!std::isfinite(*pdfReqXOff) || !std::isfinite(*pdfReqYOff) ||
1030 49424 : !std::isfinite(*pdfReqXSize) || !std::isfinite(*pdfReqYSize) ||
1031 148269 : *pdfReqXOff > INT_MAX || *pdfReqYOff > INT_MAX || *pdfReqXSize < 0 ||
1032 49423 : *pdfReqYSize < 0)
1033 : {
1034 0 : return FALSE;
1035 : }
1036 :
1037 : /* -------------------------------------------------------------------- */
1038 : /* Clamp within the bounds of the available source data. */
1039 : /* -------------------------------------------------------------------- */
1040 49423 : if (*pdfReqXOff < 0)
1041 : {
1042 6 : *pdfReqXSize += *pdfReqXOff;
1043 6 : *pdfReqXOff = 0;
1044 6 : bModifiedX = true;
1045 : }
1046 49423 : if (*pdfReqYOff < 0)
1047 : {
1048 6 : *pdfReqYSize += *pdfReqYOff;
1049 6 : *pdfReqYOff = 0;
1050 6 : bModifiedY = true;
1051 : }
1052 :
1053 49423 : *pnReqXOff = static_cast<int>(floor(*pdfReqXOff));
1054 49423 : *pnReqYOff = static_cast<int>(floor(*pdfReqYOff));
1055 :
1056 49423 : constexpr double EPS = 1e-3;
1057 49423 : constexpr double ONE_MINUS_EPS = 1.0 - EPS;
1058 49423 : if (*pdfReqXOff - *pnReqXOff > ONE_MINUS_EPS)
1059 : {
1060 2 : (*pnReqXOff)++;
1061 2 : *pdfReqXOff = *pnReqXOff;
1062 : }
1063 49423 : if (*pdfReqYOff - *pnReqYOff > ONE_MINUS_EPS)
1064 : {
1065 18 : (*pnReqYOff)++;
1066 18 : *pdfReqYOff = *pnReqYOff;
1067 : }
1068 :
1069 49423 : if (*pdfReqXSize > INT_MAX)
1070 0 : *pnReqXSize = INT_MAX;
1071 : else
1072 49423 : *pnReqXSize = static_cast<int>(floor(*pdfReqXSize + 0.5));
1073 :
1074 49423 : if (*pdfReqYSize > INT_MAX)
1075 0 : *pnReqYSize = INT_MAX;
1076 : else
1077 49423 : *pnReqYSize = static_cast<int>(floor(*pdfReqYSize + 0.5));
1078 :
1079 : /* -------------------------------------------------------------------- */
1080 : /* Clamp within the bounds of the available source data. */
1081 : /* -------------------------------------------------------------------- */
1082 :
1083 49423 : if (*pnReqXSize == 0)
1084 675 : *pnReqXSize = 1;
1085 49423 : if (*pnReqYSize == 0)
1086 22678 : *pnReqYSize = 1;
1087 :
1088 49423 : auto l_band = GetRasterBand();
1089 49423 : if (!l_band)
1090 : {
1091 106 : bErrorOut = true;
1092 106 : return FALSE;
1093 : }
1094 98634 : if (*pnReqXSize > INT_MAX - *pnReqXOff ||
1095 49317 : *pnReqXOff + *pnReqXSize > l_band->GetXSize())
1096 : {
1097 44 : *pnReqXSize = l_band->GetXSize() - *pnReqXOff;
1098 44 : bModifiedX = true;
1099 : }
1100 49317 : if (*pdfReqXOff + *pdfReqXSize > l_band->GetXSize())
1101 : {
1102 46 : *pdfReqXSize = l_band->GetXSize() - *pdfReqXOff;
1103 44 : bModifiedX = true;
1104 : }
1105 :
1106 98634 : if (*pnReqYSize > INT_MAX - *pnReqYOff ||
1107 49317 : *pnReqYOff + *pnReqYSize > l_band->GetYSize())
1108 : {
1109 44 : *pnReqYSize = l_band->GetYSize() - *pnReqYOff;
1110 44 : bModifiedY = true;
1111 : }
1112 49317 : if (*pdfReqYOff + *pdfReqYSize > l_band->GetYSize())
1113 : {
1114 56 : *pdfReqYSize = l_band->GetYSize() - *pdfReqYOff;
1115 56 : bModifiedY = true;
1116 : }
1117 :
1118 : /* -------------------------------------------------------------------- */
1119 : /* Don't do anything if the requesting region is completely off */
1120 : /* the source image. */
1121 : /* -------------------------------------------------------------------- */
1122 98630 : if (*pnReqXOff >= l_band->GetXSize() || *pnReqYOff >= l_band->GetYSize() ||
1123 98630 : *pnReqXSize <= 0 || *pnReqYSize <= 0)
1124 : {
1125 9 : return FALSE;
1126 : }
1127 :
1128 : /* -------------------------------------------------------------------- */
1129 : /* If we haven't had to modify the source rectangle, then the */
1130 : /* destination rectangle must be the whole region. */
1131 : /* -------------------------------------------------------------------- */
1132 49308 : if (bModifiedX || bModifiedY)
1133 : {
1134 : /* --------------------------------------------------------------------
1135 : */
1136 : /* Now transform this possibly reduced request back into the */
1137 : /* destination buffer coordinates in case the output region is */
1138 : /* less than the whole buffer. */
1139 : /* --------------------------------------------------------------------
1140 : */
1141 3255 : double dfDstULX = 0.0;
1142 3255 : double dfDstULY = 0.0;
1143 3255 : double dfDstLRX = 0.0;
1144 3255 : double dfDstLRY = 0.0;
1145 :
1146 3255 : SrcToDst(*pdfReqXOff, *pdfReqYOff, dfDstULX, dfDstULY);
1147 3255 : SrcToDst(*pdfReqXOff + *pdfReqXSize, *pdfReqYOff + *pdfReqYSize,
1148 : dfDstLRX, dfDstLRY);
1149 : #if DEBUG_VERBOSE
1150 : CPLDebug("VRT", "dfDstULX=%g dfDstULY=%g dfDstLRX=%g dfDstLRY=%g",
1151 : dfDstULX, dfDstULY, dfDstLRX, dfDstLRY);
1152 : #endif
1153 :
1154 3255 : if (bModifiedX)
1155 : {
1156 3211 : const double dfScaleWinToBufX = nBufXSize / dfXSize;
1157 :
1158 3211 : const double dfOutXOff = (dfDstULX - dfXOff) * dfScaleWinToBufX;
1159 3211 : if (dfOutXOff <= 0)
1160 1352 : *pnOutXOff = 0;
1161 1859 : else if (dfOutXOff > INT_MAX)
1162 0 : *pnOutXOff = INT_MAX;
1163 : else
1164 1859 : *pnOutXOff = static_cast<int>(dfOutXOff + EPS);
1165 :
1166 : // Apply correction on floating-point source window
1167 : {
1168 3211 : double dfDstDeltaX =
1169 3211 : (dfOutXOff - *pnOutXOff) / dfScaleWinToBufX;
1170 3211 : double dfSrcDeltaX = dfDstDeltaX / m_dfDstXSize * m_dfSrcXSize;
1171 3211 : *pdfReqXOff -= dfSrcDeltaX;
1172 6422 : *pdfReqXSize = std::min(*pdfReqXSize + dfSrcDeltaX,
1173 3211 : static_cast<double>(INT_MAX));
1174 : }
1175 :
1176 3211 : double dfOutRightXOff = (dfDstLRX - dfXOff) * dfScaleWinToBufX;
1177 3211 : if (dfOutRightXOff < dfOutXOff)
1178 0 : return FALSE;
1179 3211 : if (dfOutRightXOff > INT_MAX)
1180 0 : dfOutRightXOff = INT_MAX;
1181 3211 : const int nOutRightXOff =
1182 3211 : static_cast<int>(ceil(dfOutRightXOff - EPS));
1183 3211 : *pnOutXSize = nOutRightXOff - *pnOutXOff;
1184 :
1185 3211 : if (*pnOutXSize > INT_MAX - *pnOutXOff ||
1186 3211 : *pnOutXOff + *pnOutXSize > nBufXSize)
1187 0 : *pnOutXSize = nBufXSize - *pnOutXOff;
1188 :
1189 : // Apply correction on floating-point source window
1190 : {
1191 3211 : double dfDstDeltaX =
1192 3211 : (nOutRightXOff - dfOutRightXOff) / dfScaleWinToBufX;
1193 3211 : double dfSrcDeltaX = dfDstDeltaX / m_dfDstXSize * m_dfSrcXSize;
1194 6422 : *pdfReqXSize = std::min(*pdfReqXSize + dfSrcDeltaX,
1195 3211 : static_cast<double>(INT_MAX));
1196 : }
1197 : }
1198 :
1199 3255 : if (bModifiedY)
1200 : {
1201 369 : const double dfScaleWinToBufY = nBufYSize / dfYSize;
1202 :
1203 369 : const double dfOutYOff = (dfDstULY - dfYOff) * dfScaleWinToBufY;
1204 369 : if (dfOutYOff <= 0)
1205 192 : *pnOutYOff = 0;
1206 177 : else if (dfOutYOff > INT_MAX)
1207 0 : *pnOutYOff = INT_MAX;
1208 : else
1209 177 : *pnOutYOff = static_cast<int>(dfOutYOff + EPS);
1210 :
1211 : // Apply correction on floating-point source window
1212 : {
1213 369 : double dfDstDeltaY =
1214 369 : (dfOutYOff - *pnOutYOff) / dfScaleWinToBufY;
1215 369 : double dfSrcDeltaY = dfDstDeltaY / m_dfDstYSize * m_dfSrcYSize;
1216 369 : *pdfReqYOff -= dfSrcDeltaY;
1217 738 : *pdfReqYSize = std::min(*pdfReqYSize + dfSrcDeltaY,
1218 369 : static_cast<double>(INT_MAX));
1219 : }
1220 :
1221 369 : double dfOutTopYOff = (dfDstLRY - dfYOff) * dfScaleWinToBufY;
1222 369 : if (dfOutTopYOff < dfOutYOff)
1223 0 : return FALSE;
1224 369 : if (dfOutTopYOff > INT_MAX)
1225 0 : dfOutTopYOff = INT_MAX;
1226 369 : const int nOutTopYOff = static_cast<int>(ceil(dfOutTopYOff - EPS));
1227 369 : *pnOutYSize = nOutTopYOff - *pnOutYOff;
1228 :
1229 369 : if (*pnOutYSize > INT_MAX - *pnOutYOff ||
1230 369 : *pnOutYOff + *pnOutYSize > nBufYSize)
1231 0 : *pnOutYSize = nBufYSize - *pnOutYOff;
1232 :
1233 : // Apply correction on floating-point source window
1234 : {
1235 369 : double dfDstDeltaY =
1236 369 : (nOutTopYOff - dfOutTopYOff) / dfScaleWinToBufY;
1237 369 : double dfSrcDeltaY = dfDstDeltaY / m_dfDstYSize * m_dfSrcYSize;
1238 738 : *pdfReqYSize = std::min(*pdfReqYSize + dfSrcDeltaY,
1239 369 : static_cast<double>(INT_MAX));
1240 : }
1241 : }
1242 :
1243 3255 : if (*pnOutXSize < 1 || *pnOutYSize < 1)
1244 0 : return FALSE;
1245 : }
1246 :
1247 49308 : *pdfReqXOff = RoundIfCloseToInt(*pdfReqXOff);
1248 49308 : *pdfReqYOff = RoundIfCloseToInt(*pdfReqYOff);
1249 49309 : *pdfReqXSize = RoundIfCloseToInt(*pdfReqXSize);
1250 49308 : *pdfReqYSize = RoundIfCloseToInt(*pdfReqYSize);
1251 :
1252 49309 : return TRUE;
1253 : }
1254 :
1255 : /************************************************************************/
1256 : /* NeedMaxValAdjustment() */
1257 : /************************************************************************/
1258 :
1259 41986 : int VRTSimpleSource::NeedMaxValAdjustment() const
1260 : {
1261 41986 : if (!m_nMaxValue)
1262 41964 : return FALSE;
1263 :
1264 22 : auto l_band = GetRasterBand();
1265 22 : if (!l_band)
1266 0 : return FALSE;
1267 22 : const char *pszNBITS = l_band->GetMetadataItem("NBITS", "IMAGE_STRUCTURE");
1268 22 : const int nBits = (pszNBITS) ? atoi(pszNBITS) : 0;
1269 22 : if (nBits >= 1 && nBits <= 31)
1270 : {
1271 0 : const int nBandMaxValue = static_cast<int>((1U << nBits) - 1);
1272 0 : return nBandMaxValue > m_nMaxValue;
1273 : }
1274 22 : return TRUE;
1275 : }
1276 :
1277 : /************************************************************************/
1278 : /* RasterIO() */
1279 : /************************************************************************/
1280 :
1281 50213 : CPLErr VRTSimpleSource::RasterIO(GDALDataType eVRTBandDataType, int nXOff,
1282 : int nYOff, int nXSize, int nYSize, void *pData,
1283 : int nBufXSize, int nBufYSize,
1284 : GDALDataType eBufType, GSpacing nPixelSpace,
1285 : GSpacing nLineSpace,
1286 : GDALRasterIOExtraArg *psExtraArgIn,
1287 : WorkingState & /*oWorkingState*/)
1288 :
1289 : {
1290 : GDALRasterIOExtraArg sExtraArg;
1291 50213 : INIT_RASTERIO_EXTRA_ARG(sExtraArg);
1292 50213 : GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
1293 :
1294 50213 : double dfXOff = nXOff;
1295 50213 : double dfYOff = nYOff;
1296 50213 : double dfXSize = nXSize;
1297 50213 : double dfYSize = nYSize;
1298 50213 : if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
1299 : {
1300 62 : dfXOff = psExtraArgIn->dfXOff;
1301 62 : dfYOff = psExtraArgIn->dfYOff;
1302 62 : dfXSize = psExtraArgIn->dfXSize;
1303 62 : dfYSize = psExtraArgIn->dfYSize;
1304 : }
1305 :
1306 : // The window we will actually request from the source raster band.
1307 50213 : double dfReqXOff = 0.0;
1308 50213 : double dfReqYOff = 0.0;
1309 50213 : double dfReqXSize = 0.0;
1310 50213 : double dfReqYSize = 0.0;
1311 50213 : int nReqXOff = 0;
1312 50213 : int nReqYOff = 0;
1313 50213 : int nReqXSize = 0;
1314 50213 : int nReqYSize = 0;
1315 :
1316 : // The window we will actual set _within_ the pData buffer.
1317 50213 : int nOutXOff = 0;
1318 50213 : int nOutYOff = 0;
1319 50213 : int nOutXSize = 0;
1320 50213 : int nOutYSize = 0;
1321 :
1322 50213 : bool bError = false;
1323 50213 : if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
1324 : &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
1325 : &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
1326 : &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
1327 : {
1328 9313 : return bError ? CE_Failure : CE_None;
1329 : }
1330 : #if DEBUG_VERBOSE
1331 : CPLDebug("VRT",
1332 : "nXOff=%d, nYOff=%d, nXSize=%d, nYSize=%d, nBufXSize=%d, "
1333 : "nBufYSize=%d,\n"
1334 : "dfReqXOff=%g, dfReqYOff=%g, dfReqXSize=%g, dfReqYSize=%g,\n"
1335 : "nReqXOff=%d, nReqYOff=%d, nReqXSize=%d, nReqYSize=%d,\n"
1336 : "nOutXOff=%d, nOutYOff=%d, nOutXSize=%d, nOutYSize=%d",
1337 : nXOff, nYOff, nXSize, nYSize, nBufXSize, nBufYSize, dfReqXOff,
1338 : dfReqYOff, dfReqXSize, dfReqYSize, nReqXOff, nReqYOff, nReqXSize,
1339 : nReqYSize, nOutXOff, nOutYOff, nOutXSize, nOutYSize);
1340 : #endif
1341 :
1342 : /* -------------------------------------------------------------------- */
1343 : /* Actually perform the IO request. */
1344 : /* -------------------------------------------------------------------- */
1345 40900 : if (!m_osResampling.empty())
1346 : {
1347 2050 : psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
1348 : }
1349 38850 : else if (psExtraArgIn != nullptr)
1350 : {
1351 38850 : psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
1352 : }
1353 40900 : psExtraArg->bFloatingPointWindowValidity = TRUE;
1354 40900 : psExtraArg->dfXOff = dfReqXOff;
1355 40900 : psExtraArg->dfYOff = dfReqYOff;
1356 40900 : psExtraArg->dfXSize = dfReqXSize;
1357 40900 : psExtraArg->dfYSize = dfReqYSize;
1358 40900 : if (psExtraArgIn)
1359 : {
1360 40900 : psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
1361 40900 : psExtraArg->pProgressData = psExtraArgIn->pProgressData;
1362 : }
1363 :
1364 40900 : GByte *pabyOut = static_cast<unsigned char *>(pData) +
1365 40900 : nOutXOff * nPixelSpace +
1366 40900 : static_cast<GPtrDiff_t>(nOutYOff) * nLineSpace;
1367 :
1368 40900 : auto l_band = GetRasterBand();
1369 40900 : if (!l_band)
1370 0 : return CE_Failure;
1371 :
1372 40900 : CPLErr eErr = CE_Failure;
1373 40900 : if (GDALDataTypeIsConversionLossy(l_band->GetRasterDataType(),
1374 40899 : eVRTBandDataType))
1375 : {
1376 49 : const int nBandDTSize = GDALGetDataTypeSizeBytes(eVRTBandDataType);
1377 49 : void *pTemp = VSI_MALLOC3_VERBOSE(nOutXSize, nOutYSize, nBandDTSize);
1378 49 : if (pTemp)
1379 : {
1380 49 : eErr = l_band->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize,
1381 : nReqYSize, pTemp, nOutXSize, nOutYSize,
1382 : eVRTBandDataType, 0, 0, psExtraArg);
1383 49 : if (eErr == CE_None)
1384 : {
1385 49 : GByte *pabyTemp = static_cast<GByte *>(pTemp);
1386 779 : for (int iY = 0; iY < nOutYSize; iY++)
1387 : {
1388 730 : GDALCopyWords(
1389 730 : pabyTemp +
1390 730 : static_cast<size_t>(iY) * nBandDTSize * nOutXSize,
1391 : eVRTBandDataType, nBandDTSize,
1392 730 : pabyOut + static_cast<GPtrDiff_t>(iY * nLineSpace),
1393 : eBufType, static_cast<int>(nPixelSpace), nOutXSize);
1394 : }
1395 : }
1396 49 : VSIFree(pTemp);
1397 : }
1398 : }
1399 : else
1400 : {
1401 40850 : eErr = l_band->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize,
1402 : nReqYSize, pabyOut, nOutXSize, nOutYSize,
1403 : eBufType, nPixelSpace, nLineSpace, psExtraArg);
1404 : }
1405 :
1406 40900 : if (NeedMaxValAdjustment())
1407 : {
1408 234 : for (int j = 0; j < nOutYSize; j++)
1409 : {
1410 4620 : for (int i = 0; i < nOutXSize; i++)
1411 : {
1412 4400 : int nVal = 0;
1413 4400 : GDALCopyWords(pabyOut + j * nLineSpace + i * nPixelSpace,
1414 : eBufType, 0, &nVal, GDT_Int32, 0, 1);
1415 4400 : if (nVal > m_nMaxValue)
1416 800 : nVal = m_nMaxValue;
1417 4400 : GDALCopyWords(&nVal, GDT_Int32, 0,
1418 4400 : pabyOut + j * nLineSpace + i * nPixelSpace,
1419 : eBufType, 0, 1);
1420 : }
1421 : }
1422 : }
1423 :
1424 40900 : if (psExtraArg->pfnProgress)
1425 2928 : psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
1426 :
1427 40900 : return eErr;
1428 : }
1429 :
1430 : /************************************************************************/
1431 : /* GetMinimum() */
1432 : /************************************************************************/
1433 :
1434 52 : double VRTSimpleSource::GetMinimum(int nXSize, int nYSize, int *pbSuccess)
1435 : {
1436 : // The window we will actually request from the source raster band.
1437 52 : double dfReqXOff = 0.0;
1438 52 : double dfReqYOff = 0.0;
1439 52 : double dfReqXSize = 0.0;
1440 52 : double dfReqYSize = 0.0;
1441 52 : int nReqXOff = 0;
1442 52 : int nReqYOff = 0;
1443 52 : int nReqXSize = 0;
1444 52 : int nReqYSize = 0;
1445 :
1446 : // The window we will actual set _within_ the pData buffer.
1447 52 : int nOutXOff = 0;
1448 52 : int nOutYOff = 0;
1449 52 : int nOutXSize = 0;
1450 52 : int nOutYSize = 0;
1451 :
1452 52 : bool bError = false;
1453 52 : auto l_band = GetRasterBand();
1454 104 : if (!l_band ||
1455 52 : !GetSrcDstWindow(0, 0, nXSize, nYSize, nXSize, nYSize, &dfReqXOff,
1456 : &dfReqYOff, &dfReqXSize, &dfReqYSize, &nReqXOff,
1457 : &nReqYOff, &nReqXSize, &nReqYSize, &nOutXOff,
1458 52 : &nOutYOff, &nOutXSize, &nOutYSize, bError) ||
1459 156 : nReqXOff != 0 || nReqYOff != 0 || nReqXSize != l_band->GetXSize() ||
1460 52 : nReqYSize != l_band->GetYSize())
1461 : {
1462 0 : *pbSuccess = FALSE;
1463 0 : return 0;
1464 : }
1465 :
1466 52 : const double dfVal = l_band->GetMinimum(pbSuccess);
1467 52 : if (NeedMaxValAdjustment() && dfVal > m_nMaxValue)
1468 2 : return m_nMaxValue;
1469 50 : return dfVal;
1470 : }
1471 :
1472 : /************************************************************************/
1473 : /* GetMaximum() */
1474 : /************************************************************************/
1475 :
1476 51 : double VRTSimpleSource::GetMaximum(int nXSize, int nYSize, int *pbSuccess)
1477 : {
1478 : // The window we will actually request from the source raster band.
1479 51 : double dfReqXOff = 0.0;
1480 51 : double dfReqYOff = 0.0;
1481 51 : double dfReqXSize = 0.0;
1482 51 : double dfReqYSize = 0.0;
1483 51 : int nReqXOff = 0;
1484 51 : int nReqYOff = 0;
1485 51 : int nReqXSize = 0;
1486 51 : int nReqYSize = 0;
1487 :
1488 : // The window we will actual set _within_ the pData buffer.
1489 51 : int nOutXOff = 0;
1490 51 : int nOutYOff = 0;
1491 51 : int nOutXSize = 0;
1492 51 : int nOutYSize = 0;
1493 :
1494 51 : bool bError = false;
1495 51 : auto l_band = GetRasterBand();
1496 102 : if (!l_band ||
1497 51 : !GetSrcDstWindow(0, 0, nXSize, nYSize, nXSize, nYSize, &dfReqXOff,
1498 : &dfReqYOff, &dfReqXSize, &dfReqYSize, &nReqXOff,
1499 : &nReqYOff, &nReqXSize, &nReqYSize, &nOutXOff,
1500 51 : &nOutYOff, &nOutXSize, &nOutYSize, bError) ||
1501 153 : nReqXOff != 0 || nReqYOff != 0 || nReqXSize != l_band->GetXSize() ||
1502 51 : nReqYSize != l_band->GetYSize())
1503 : {
1504 0 : *pbSuccess = FALSE;
1505 0 : return 0;
1506 : }
1507 :
1508 51 : const double dfVal = l_band->GetMaximum(pbSuccess);
1509 51 : if (NeedMaxValAdjustment() && dfVal > m_nMaxValue)
1510 2 : return m_nMaxValue;
1511 49 : return dfVal;
1512 : }
1513 :
1514 : /************************************************************************/
1515 : /* GetHistogram() */
1516 : /************************************************************************/
1517 :
1518 4 : CPLErr VRTSimpleSource::GetHistogram(int nXSize, int nYSize, double dfMin,
1519 : double dfMax, int nBuckets,
1520 : GUIntBig *panHistogram,
1521 : int bIncludeOutOfRange, int bApproxOK,
1522 : GDALProgressFunc pfnProgress,
1523 : void *pProgressData)
1524 : {
1525 : // The window we will actually request from the source raster band.
1526 4 : double dfReqXOff = 0.0;
1527 4 : double dfReqYOff = 0.0;
1528 4 : double dfReqXSize = 0.0;
1529 4 : double dfReqYSize = 0.0;
1530 4 : int nReqXOff = 0;
1531 4 : int nReqYOff = 0;
1532 4 : int nReqXSize = 0;
1533 4 : int nReqYSize = 0;
1534 :
1535 : // The window we will actual set _within_ the pData buffer.
1536 4 : int nOutXOff = 0;
1537 4 : int nOutYOff = 0;
1538 4 : int nOutXSize = 0;
1539 4 : int nOutYSize = 0;
1540 :
1541 4 : bool bError = false;
1542 4 : auto l_band = GetRasterBand();
1543 4 : if (!l_band || NeedMaxValAdjustment() ||
1544 4 : !GetSrcDstWindow(0, 0, nXSize, nYSize, nXSize, nYSize, &dfReqXOff,
1545 : &dfReqYOff, &dfReqXSize, &dfReqYSize, &nReqXOff,
1546 : &nReqYOff, &nReqXSize, &nReqYSize, &nOutXOff,
1547 4 : &nOutYOff, &nOutXSize, &nOutYSize, bError) ||
1548 12 : nReqXOff != 0 || nReqYOff != 0 || nReqXSize != l_band->GetXSize() ||
1549 4 : nReqYSize != l_band->GetYSize())
1550 : {
1551 0 : return CE_Failure;
1552 : }
1553 :
1554 4 : return l_band->GetHistogram(dfMin, dfMax, nBuckets, panHistogram,
1555 : bIncludeOutOfRange, bApproxOK, pfnProgress,
1556 4 : pProgressData);
1557 : }
1558 :
1559 : /************************************************************************/
1560 : /* DatasetRasterIO() */
1561 : /************************************************************************/
1562 :
1563 992 : CPLErr VRTSimpleSource::DatasetRasterIO(
1564 : GDALDataType eVRTBandDataType, int nXOff, int nYOff, int nXSize, int nYSize,
1565 : void *pData, int nBufXSize, int nBufYSize, GDALDataType eBufType,
1566 : int nBandCount, const int *panBandMap, GSpacing nPixelSpace,
1567 : GSpacing nLineSpace, GSpacing nBandSpace,
1568 : GDALRasterIOExtraArg *psExtraArgIn)
1569 : {
1570 992 : if (GetType() != VRTSimpleSource::GetTypeStatic())
1571 : {
1572 0 : CPLError(CE_Failure, CPLE_NotSupported,
1573 0 : "DatasetRasterIO() not implemented for %s", GetType());
1574 0 : return CE_Failure;
1575 : }
1576 :
1577 : GDALRasterIOExtraArg sExtraArg;
1578 991 : INIT_RASTERIO_EXTRA_ARG(sExtraArg);
1579 991 : GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
1580 :
1581 991 : double dfXOff = nXOff;
1582 991 : double dfYOff = nYOff;
1583 991 : double dfXSize = nXSize;
1584 991 : double dfYSize = nYSize;
1585 991 : if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
1586 : {
1587 11 : dfXOff = psExtraArgIn->dfXOff;
1588 11 : dfYOff = psExtraArgIn->dfYOff;
1589 11 : dfXSize = psExtraArgIn->dfXSize;
1590 11 : dfYSize = psExtraArgIn->dfYSize;
1591 : }
1592 :
1593 : // The window we will actually request from the source raster band.
1594 991 : double dfReqXOff = 0.0;
1595 991 : double dfReqYOff = 0.0;
1596 991 : double dfReqXSize = 0.0;
1597 991 : double dfReqYSize = 0.0;
1598 991 : int nReqXOff = 0;
1599 991 : int nReqYOff = 0;
1600 991 : int nReqXSize = 0;
1601 991 : int nReqYSize = 0;
1602 :
1603 : // The window we will actual set _within_ the pData buffer.
1604 991 : int nOutXOff = 0;
1605 991 : int nOutYOff = 0;
1606 991 : int nOutXSize = 0;
1607 991 : int nOutYSize = 0;
1608 :
1609 991 : bool bError = false;
1610 991 : if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
1611 : &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
1612 : &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
1613 : &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
1614 : {
1615 71 : return bError ? CE_Failure : CE_None;
1616 : }
1617 :
1618 921 : auto l_band = GetRasterBand();
1619 920 : if (!l_band)
1620 0 : return CE_Failure;
1621 :
1622 920 : GDALDataset *poDS = l_band->GetDataset();
1623 920 : if (poDS == nullptr)
1624 0 : return CE_Failure;
1625 :
1626 920 : if (!m_osResampling.empty())
1627 : {
1628 6 : psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
1629 : }
1630 914 : else if (psExtraArgIn != nullptr)
1631 : {
1632 915 : psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
1633 : }
1634 920 : psExtraArg->bFloatingPointWindowValidity = TRUE;
1635 920 : psExtraArg->dfXOff = dfReqXOff;
1636 920 : psExtraArg->dfYOff = dfReqYOff;
1637 920 : psExtraArg->dfXSize = dfReqXSize;
1638 920 : psExtraArg->dfYSize = dfReqYSize;
1639 920 : if (psExtraArgIn)
1640 : {
1641 921 : psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
1642 921 : psExtraArg->pProgressData = psExtraArgIn->pProgressData;
1643 : }
1644 :
1645 920 : GByte *pabyOut = static_cast<unsigned char *>(pData) +
1646 920 : nOutXOff * nPixelSpace +
1647 920 : static_cast<GPtrDiff_t>(nOutYOff) * nLineSpace;
1648 :
1649 920 : CPLErr eErr = CE_Failure;
1650 :
1651 920 : if (GDALDataTypeIsConversionLossy(l_band->GetRasterDataType(),
1652 920 : eVRTBandDataType))
1653 : {
1654 1 : const int nBandDTSize = GDALGetDataTypeSizeBytes(eVRTBandDataType);
1655 1 : void *pTemp = VSI_MALLOC3_VERBOSE(
1656 : nOutXSize, nOutYSize, cpl::fits_on<int>(nBandDTSize * nBandCount));
1657 1 : if (pTemp)
1658 : {
1659 1 : eErr = poDS->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize,
1660 : nReqYSize, pTemp, nOutXSize, nOutYSize,
1661 : eVRTBandDataType, nBandCount, panBandMap, 0,
1662 : 0, 0, psExtraArg);
1663 1 : if (eErr == CE_None)
1664 : {
1665 1 : GByte *pabyTemp = static_cast<GByte *>(pTemp);
1666 1 : const size_t nSrcBandSpace =
1667 1 : static_cast<size_t>(nOutYSize) * nOutXSize * nBandDTSize;
1668 2 : for (int iBand = 0; iBand < nBandCount; iBand++)
1669 : {
1670 21 : for (int iY = 0; iY < nOutYSize; iY++)
1671 : {
1672 20 : GDALCopyWords(
1673 20 : pabyTemp + iBand * nSrcBandSpace +
1674 20 : static_cast<size_t>(iY) * nBandDTSize *
1675 20 : nOutXSize,
1676 : eVRTBandDataType, nBandDTSize,
1677 20 : pabyOut + static_cast<GPtrDiff_t>(
1678 20 : iY * nLineSpace + iBand * nBandSpace),
1679 : eBufType, static_cast<int>(nPixelSpace), nOutXSize);
1680 : }
1681 : }
1682 : }
1683 1 : VSIFree(pTemp);
1684 : }
1685 : }
1686 : else
1687 : {
1688 919 : eErr = poDS->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
1689 : pabyOut, nOutXSize, nOutYSize, eBufType,
1690 : nBandCount, panBandMap, nPixelSpace, nLineSpace,
1691 : nBandSpace, psExtraArg);
1692 : }
1693 :
1694 921 : if (NeedMaxValAdjustment())
1695 : {
1696 0 : for (int k = 0; k < nBandCount; k++)
1697 : {
1698 0 : for (int j = 0; j < nOutYSize; j++)
1699 : {
1700 0 : for (int i = 0; i < nOutXSize; i++)
1701 : {
1702 0 : int nVal = 0;
1703 0 : GDALCopyWords(pabyOut + k * nBandSpace + j * nLineSpace +
1704 0 : i * nPixelSpace,
1705 : eBufType, 0, &nVal, GDT_Int32, 0, 1);
1706 :
1707 0 : if (nVal > m_nMaxValue)
1708 0 : nVal = m_nMaxValue;
1709 :
1710 0 : GDALCopyWords(&nVal, GDT_Int32, 0,
1711 0 : pabyOut + k * nBandSpace + j * nLineSpace +
1712 0 : i * nPixelSpace,
1713 : eBufType, 0, 1);
1714 : }
1715 : }
1716 : }
1717 : }
1718 :
1719 921 : if (psExtraArg->pfnProgress)
1720 460 : psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
1721 :
1722 921 : return eErr;
1723 : }
1724 :
1725 : /************************************************************************/
1726 : /* SetResampling() */
1727 : /************************************************************************/
1728 :
1729 4069 : void VRTSimpleSource::SetResampling(const char *pszResampling)
1730 : {
1731 4069 : m_osResampling = (pszResampling) ? pszResampling : "";
1732 4069 : }
1733 :
1734 : /************************************************************************/
1735 : /* ==================================================================== */
1736 : /* VRTAveragedSource */
1737 : /* ==================================================================== */
1738 : /************************************************************************/
1739 :
1740 : /************************************************************************/
1741 : /* VRTAveragedSource() */
1742 : /************************************************************************/
1743 :
1744 16 : VRTAveragedSource::VRTAveragedSource()
1745 : {
1746 16 : }
1747 :
1748 : /************************************************************************/
1749 : /* GetTypeStatic() */
1750 : /************************************************************************/
1751 :
1752 3094 : const char *VRTAveragedSource::GetTypeStatic()
1753 : {
1754 : static const char *TYPE = "AveragedSource";
1755 3094 : return TYPE;
1756 : }
1757 :
1758 : /************************************************************************/
1759 : /* GetType() */
1760 : /************************************************************************/
1761 :
1762 11 : const char *VRTAveragedSource::GetType() const
1763 : {
1764 11 : return GetTypeStatic();
1765 : }
1766 :
1767 : /************************************************************************/
1768 : /* SerializeToXML() */
1769 : /************************************************************************/
1770 :
1771 0 : CPLXMLNode *VRTAveragedSource::SerializeToXML(const char *pszVRTPath)
1772 :
1773 : {
1774 0 : CPLXMLNode *const psSrc = VRTSimpleSource::SerializeToXML(pszVRTPath);
1775 :
1776 0 : if (psSrc == nullptr)
1777 0 : return nullptr;
1778 :
1779 0 : CPLFree(psSrc->pszValue);
1780 0 : psSrc->pszValue = CPLStrdup(GetTypeStatic());
1781 :
1782 0 : return psSrc;
1783 : }
1784 :
1785 : /************************************************************************/
1786 : /* SetNoDataValue() */
1787 : /************************************************************************/
1788 :
1789 0 : void VRTAveragedSource::SetNoDataValue(double dfNewNoDataValue)
1790 :
1791 : {
1792 0 : if (dfNewNoDataValue == VRT_NODATA_UNSET)
1793 : {
1794 0 : m_bNoDataSet = FALSE;
1795 0 : m_dfNoDataValue = VRT_NODATA_UNSET;
1796 0 : return;
1797 : }
1798 :
1799 0 : m_bNoDataSet = TRUE;
1800 0 : m_dfNoDataValue = dfNewNoDataValue;
1801 : }
1802 :
1803 : /************************************************************************/
1804 : /* RasterIO() */
1805 : /************************************************************************/
1806 :
1807 33 : CPLErr VRTAveragedSource::RasterIO(GDALDataType /*eVRTBandDataType*/, int nXOff,
1808 : int nYOff, int nXSize, int nYSize,
1809 : void *pData, int nBufXSize, int nBufYSize,
1810 : GDALDataType eBufType, GSpacing nPixelSpace,
1811 : GSpacing nLineSpace,
1812 : GDALRasterIOExtraArg *psExtraArgIn,
1813 : WorkingState & /*oWorkingState*/)
1814 :
1815 : {
1816 : GDALRasterIOExtraArg sExtraArg;
1817 33 : INIT_RASTERIO_EXTRA_ARG(sExtraArg);
1818 33 : GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
1819 :
1820 33 : double dfXOff = nXOff;
1821 33 : double dfYOff = nYOff;
1822 33 : double dfXSize = nXSize;
1823 33 : double dfYSize = nYSize;
1824 33 : if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
1825 : {
1826 0 : dfXOff = psExtraArgIn->dfXOff;
1827 0 : dfYOff = psExtraArgIn->dfYOff;
1828 0 : dfXSize = psExtraArgIn->dfXSize;
1829 0 : dfYSize = psExtraArgIn->dfYSize;
1830 : }
1831 :
1832 : // The window we will actually request from the source raster band.
1833 33 : double dfReqXOff = 0.0;
1834 33 : double dfReqYOff = 0.0;
1835 33 : double dfReqXSize = 0.0;
1836 33 : double dfReqYSize = 0.0;
1837 33 : int nReqXOff = 0;
1838 33 : int nReqYOff = 0;
1839 33 : int nReqXSize = 0;
1840 33 : int nReqYSize = 0;
1841 :
1842 : // The window we will actual set _within_ the pData buffer.
1843 33 : int nOutXOff = 0;
1844 33 : int nOutYOff = 0;
1845 33 : int nOutXSize = 0;
1846 33 : int nOutYSize = 0;
1847 :
1848 33 : bool bError = false;
1849 33 : if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
1850 : &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
1851 : &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
1852 : &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
1853 : {
1854 0 : return bError ? CE_Failure : CE_None;
1855 : }
1856 :
1857 33 : auto l_band = GetRasterBand();
1858 33 : if (!l_band)
1859 0 : return CE_Failure;
1860 :
1861 : /* -------------------------------------------------------------------- */
1862 : /* Allocate a temporary buffer to whole the full resolution */
1863 : /* data from the area of interest. */
1864 : /* -------------------------------------------------------------------- */
1865 : float *const pafSrc = static_cast<float *>(
1866 33 : VSI_MALLOC3_VERBOSE(sizeof(float), nReqXSize, nReqYSize));
1867 33 : if (pafSrc == nullptr)
1868 : {
1869 0 : return CE_Failure;
1870 : }
1871 :
1872 : /* -------------------------------------------------------------------- */
1873 : /* Load it. */
1874 : /* -------------------------------------------------------------------- */
1875 33 : if (!m_osResampling.empty())
1876 : {
1877 28 : psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
1878 : }
1879 5 : else if (psExtraArgIn != nullptr)
1880 : {
1881 5 : psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
1882 : }
1883 :
1884 33 : psExtraArg->bFloatingPointWindowValidity = TRUE;
1885 33 : psExtraArg->dfXOff = dfReqXOff;
1886 33 : psExtraArg->dfYOff = dfReqYOff;
1887 33 : psExtraArg->dfXSize = dfReqXSize;
1888 33 : psExtraArg->dfYSize = dfReqYSize;
1889 33 : if (psExtraArgIn)
1890 : {
1891 33 : psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
1892 33 : psExtraArg->pProgressData = psExtraArgIn->pProgressData;
1893 : }
1894 :
1895 33 : const CPLErr eErr = l_band->RasterIO(
1896 : GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize, pafSrc, nReqXSize,
1897 : nReqYSize, GDT_Float32, 0, 0, psExtraArg);
1898 :
1899 33 : if (eErr != CE_None)
1900 : {
1901 0 : VSIFree(pafSrc);
1902 0 : return eErr;
1903 : }
1904 :
1905 : /* -------------------------------------------------------------------- */
1906 : /* Do the averaging. */
1907 : /* -------------------------------------------------------------------- */
1908 5956 : for (int iBufLine = nOutYOff; iBufLine < nOutYOff + nOutYSize; iBufLine++)
1909 : {
1910 5923 : const double dfYDst =
1911 5923 : (iBufLine / static_cast<double>(nBufYSize)) * nYSize + nYOff;
1912 :
1913 2017080 : for (int iBufPixel = nOutXOff; iBufPixel < nOutXOff + nOutXSize;
1914 : iBufPixel++)
1915 : {
1916 : double dfXSrcStart, dfXSrcEnd, dfYSrcStart, dfYSrcEnd;
1917 : int iXSrcStart, iYSrcStart, iXSrcEnd, iYSrcEnd;
1918 :
1919 2011160 : const double dfXDst =
1920 2011160 : (iBufPixel / static_cast<double>(nBufXSize)) * nXSize + nXOff;
1921 :
1922 : // Compute the source image rectangle needed for this pixel.
1923 2011160 : DstToSrc(dfXDst, dfYDst, dfXSrcStart, dfYSrcStart);
1924 2011160 : DstToSrc(dfXDst + 1.0, dfYDst + 1.0, dfXSrcEnd, dfYSrcEnd);
1925 :
1926 : // Convert to integers, assuming that the center of the source
1927 : // pixel must be in our rect to get included.
1928 2011160 : if (dfXSrcEnd >= dfXSrcStart + 1)
1929 : {
1930 1049560 : iXSrcStart = static_cast<int>(floor(dfXSrcStart + 0.5));
1931 1049560 : iXSrcEnd = static_cast<int>(floor(dfXSrcEnd + 0.5));
1932 : }
1933 : else
1934 : {
1935 : /* If the resampling factor is less than 100%, the distance */
1936 : /* between the source pixel is < 1, so we stick to nearest */
1937 : /* neighbour */
1938 961600 : iXSrcStart = static_cast<int>(floor(dfXSrcStart));
1939 961600 : iXSrcEnd = iXSrcStart + 1;
1940 : }
1941 2011160 : if (dfYSrcEnd >= dfYSrcStart + 1)
1942 : {
1943 1049560 : iYSrcStart = static_cast<int>(floor(dfYSrcStart + 0.5));
1944 1049560 : iYSrcEnd = static_cast<int>(floor(dfYSrcEnd + 0.5));
1945 : }
1946 : else
1947 : {
1948 961600 : iYSrcStart = static_cast<int>(floor(dfYSrcStart));
1949 961600 : iYSrcEnd = iYSrcStart + 1;
1950 : }
1951 :
1952 : // Transform into the coordinate system of the source *buffer*
1953 2011160 : iXSrcStart -= nReqXOff;
1954 2011160 : iYSrcStart -= nReqYOff;
1955 2011160 : iXSrcEnd -= nReqXOff;
1956 2011160 : iYSrcEnd -= nReqYOff;
1957 :
1958 2011160 : double dfSum = 0.0;
1959 2011160 : int nPixelCount = 0;
1960 :
1961 4022510 : for (int iY = iYSrcStart; iY < iYSrcEnd; iY++)
1962 : {
1963 2011360 : if (iY < 0 || iY >= nReqYSize)
1964 0 : continue;
1965 :
1966 4023130 : for (int iX = iXSrcStart; iX < iXSrcEnd; iX++)
1967 : {
1968 2011780 : if (iX < 0 || iX >= nReqXSize)
1969 0 : continue;
1970 :
1971 2011780 : const float fSampledValue =
1972 2011780 : pafSrc[iX + static_cast<size_t>(iY) * nReqXSize];
1973 2011780 : if (std::isnan(fSampledValue))
1974 0 : continue;
1975 :
1976 0 : if (m_bNoDataSet &&
1977 2011780 : GDALIsValueInRange<float>(m_dfNoDataValue) &&
1978 0 : ARE_REAL_EQUAL(fSampledValue,
1979 0 : static_cast<float>(m_dfNoDataValue)))
1980 0 : continue;
1981 :
1982 2011780 : nPixelCount++;
1983 2011780 : dfSum += pafSrc[iX + static_cast<size_t>(iY) * nReqXSize];
1984 : }
1985 : }
1986 :
1987 2011160 : if (nPixelCount == 0)
1988 0 : continue;
1989 :
1990 : // Compute output value.
1991 2011160 : const float dfOutputValue = static_cast<float>(dfSum / nPixelCount);
1992 :
1993 : // Put it in the output buffer.
1994 2011160 : GByte *pDstLocation =
1995 2011160 : static_cast<GByte *>(pData) + nPixelSpace * iBufPixel +
1996 2011160 : static_cast<GPtrDiff_t>(nLineSpace) * iBufLine;
1997 :
1998 2011160 : if (eBufType == GDT_Byte)
1999 2008660 : *pDstLocation = static_cast<GByte>(
2000 2008660 : std::min(255.0, std::max(0.0, dfOutputValue + 0.5)));
2001 : else
2002 2500 : GDALCopyWords(&dfOutputValue, GDT_Float32, 4, pDstLocation,
2003 : eBufType, 8, 1);
2004 : }
2005 : }
2006 :
2007 33 : VSIFree(pafSrc);
2008 :
2009 33 : if (psExtraArg->pfnProgress)
2010 0 : psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
2011 :
2012 33 : return CE_None;
2013 : }
2014 :
2015 : /************************************************************************/
2016 : /* GetMinimum() */
2017 : /************************************************************************/
2018 :
2019 0 : double VRTAveragedSource::GetMinimum(int /* nXSize */, int /* nYSize */,
2020 : int *pbSuccess)
2021 : {
2022 0 : *pbSuccess = FALSE;
2023 0 : return 0.0;
2024 : }
2025 :
2026 : /************************************************************************/
2027 : /* GetMaximum() */
2028 : /************************************************************************/
2029 :
2030 0 : double VRTAveragedSource::GetMaximum(int /* nXSize */, int /* nYSize */,
2031 : int *pbSuccess)
2032 : {
2033 0 : *pbSuccess = FALSE;
2034 0 : return 0.0;
2035 : }
2036 :
2037 : /************************************************************************/
2038 : /* GetHistogram() */
2039 : /************************************************************************/
2040 :
2041 0 : CPLErr VRTAveragedSource::GetHistogram(
2042 : int /* nXSize */, int /* nYSize */, double /* dfMin */, double /* dfMax */,
2043 : int /* nBuckets */, GUIntBig * /* panHistogram */,
2044 : int /* bIncludeOutOfRange */, int /* bApproxOK */,
2045 : GDALProgressFunc /* pfnProgress */, void * /* pProgressData */)
2046 : {
2047 0 : return CE_Failure;
2048 : }
2049 :
2050 : /************************************************************************/
2051 : /* ==================================================================== */
2052 : /* VRTNoDataFromMaskSource */
2053 : /* ==================================================================== */
2054 : /************************************************************************/
2055 :
2056 : /************************************************************************/
2057 : /* VRTNoDataFromMaskSource() */
2058 : /************************************************************************/
2059 :
2060 23 : VRTNoDataFromMaskSource::VRTNoDataFromMaskSource()
2061 : {
2062 23 : }
2063 :
2064 : /************************************************************************/
2065 : /* XMLInit() */
2066 : /************************************************************************/
2067 :
2068 : CPLErr
2069 8 : VRTNoDataFromMaskSource::XMLInit(const CPLXMLNode *psSrc,
2070 : const char *pszVRTPath,
2071 : VRTMapSharedResources &oMapSharedSources)
2072 :
2073 : {
2074 : /* -------------------------------------------------------------------- */
2075 : /* Do base initialization. */
2076 : /* -------------------------------------------------------------------- */
2077 : {
2078 : const CPLErr eErr =
2079 8 : VRTSimpleSource::XMLInit(psSrc, pszVRTPath, oMapSharedSources);
2080 8 : if (eErr != CE_None)
2081 0 : return eErr;
2082 : }
2083 :
2084 8 : if (const char *pszNODATA = CPLGetXMLValue(psSrc, "NODATA", nullptr))
2085 : {
2086 8 : m_bNoDataSet = true;
2087 8 : m_dfNoDataValue = CPLAtofM(pszNODATA);
2088 : }
2089 :
2090 8 : m_dfMaskValueThreshold =
2091 8 : CPLAtofM(CPLGetXMLValue(psSrc, "MaskValueThreshold", "0"));
2092 :
2093 8 : if (const char *pszRemappedValue =
2094 8 : CPLGetXMLValue(psSrc, "RemappedValue", nullptr))
2095 : {
2096 0 : m_bHasRemappedValue = true;
2097 0 : m_dfRemappedValue = CPLAtofM(pszRemappedValue);
2098 : }
2099 :
2100 8 : return CE_None;
2101 : }
2102 :
2103 : /************************************************************************/
2104 : /* GetTypeStatic() */
2105 : /************************************************************************/
2106 :
2107 27 : const char *VRTNoDataFromMaskSource::GetTypeStatic()
2108 : {
2109 : static const char *TYPE = "NoDataFromMaskSource";
2110 27 : return TYPE;
2111 : }
2112 :
2113 : /************************************************************************/
2114 : /* GetType() */
2115 : /************************************************************************/
2116 :
2117 11 : const char *VRTNoDataFromMaskSource::GetType() const
2118 : {
2119 11 : return GetTypeStatic();
2120 : }
2121 :
2122 : /************************************************************************/
2123 : /* SerializeToXML() */
2124 : /************************************************************************/
2125 :
2126 8 : CPLXMLNode *VRTNoDataFromMaskSource::SerializeToXML(const char *pszVRTPath)
2127 :
2128 : {
2129 8 : CPLXMLNode *const psSrc = VRTSimpleSource::SerializeToXML(pszVRTPath);
2130 :
2131 8 : if (psSrc == nullptr)
2132 0 : return nullptr;
2133 :
2134 8 : CPLFree(psSrc->pszValue);
2135 8 : psSrc->pszValue = CPLStrdup(GetTypeStatic());
2136 :
2137 8 : if (m_bNoDataSet)
2138 : {
2139 8 : CPLSetXMLValue(psSrc, "MaskValueThreshold",
2140 : CPLSPrintf("%.17g", m_dfMaskValueThreshold));
2141 :
2142 8 : GDALDataType eBandDT = GDT_Unknown;
2143 8 : double dfNoDataValue = m_dfNoDataValue;
2144 8 : const auto kMaxFloat = std::numeric_limits<float>::max();
2145 8 : if (std::fabs(std::fabs(m_dfNoDataValue) - kMaxFloat) <
2146 : 1e-10 * kMaxFloat)
2147 : {
2148 0 : auto l_band = GetRasterBand();
2149 0 : if (l_band)
2150 : {
2151 0 : eBandDT = l_band->GetRasterDataType();
2152 0 : if (eBandDT == GDT_Float32)
2153 : {
2154 : dfNoDataValue =
2155 0 : GDALAdjustNoDataCloseToFloatMax(m_dfNoDataValue);
2156 : }
2157 : }
2158 : }
2159 8 : CPLSetXMLValue(psSrc, "NODATA",
2160 16 : VRTSerializeNoData(dfNoDataValue, eBandDT, 18).c_str());
2161 : }
2162 :
2163 8 : if (m_bHasRemappedValue)
2164 : {
2165 0 : CPLSetXMLValue(psSrc, "RemappedValue",
2166 : CPLSPrintf("%.17g", m_dfRemappedValue));
2167 : }
2168 :
2169 8 : return psSrc;
2170 : }
2171 :
2172 : /************************************************************************/
2173 : /* SetParameters() */
2174 : /************************************************************************/
2175 :
2176 15 : void VRTNoDataFromMaskSource::SetParameters(double dfNoDataValue,
2177 : double dfMaskValueThreshold)
2178 : {
2179 15 : m_bNoDataSet = true;
2180 15 : m_dfNoDataValue = dfNoDataValue;
2181 15 : m_dfMaskValueThreshold = dfMaskValueThreshold;
2182 15 : if (!m_bHasRemappedValue)
2183 15 : m_dfRemappedValue = m_dfNoDataValue;
2184 15 : }
2185 :
2186 : /************************************************************************/
2187 : /* SetParameters() */
2188 : /************************************************************************/
2189 :
2190 0 : void VRTNoDataFromMaskSource::SetParameters(double dfNoDataValue,
2191 : double dfMaskValueThreshold,
2192 : double dfRemappedValue)
2193 : {
2194 0 : SetParameters(dfNoDataValue, dfMaskValueThreshold);
2195 0 : m_bHasRemappedValue = true;
2196 0 : m_dfRemappedValue = dfRemappedValue;
2197 0 : }
2198 :
2199 : /************************************************************************/
2200 : /* RasterIO() */
2201 : /************************************************************************/
2202 :
2203 13 : CPLErr VRTNoDataFromMaskSource::RasterIO(
2204 : GDALDataType eVRTBandDataType, int nXOff, int nYOff, int nXSize, int nYSize,
2205 : void *pData, int nBufXSize, int nBufYSize, GDALDataType eBufType,
2206 : GSpacing nPixelSpace, GSpacing nLineSpace,
2207 : GDALRasterIOExtraArg *psExtraArgIn, WorkingState &oWorkingState)
2208 :
2209 : {
2210 13 : if (!m_bNoDataSet)
2211 : {
2212 0 : return VRTSimpleSource::RasterIO(eVRTBandDataType, nXOff, nYOff, nXSize,
2213 : nYSize, pData, nBufXSize, nBufYSize,
2214 : eBufType, nPixelSpace, nLineSpace,
2215 0 : psExtraArgIn, oWorkingState);
2216 : }
2217 :
2218 : GDALRasterIOExtraArg sExtraArg;
2219 13 : INIT_RASTERIO_EXTRA_ARG(sExtraArg);
2220 13 : GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
2221 :
2222 13 : double dfXOff = nXOff;
2223 13 : double dfYOff = nYOff;
2224 13 : double dfXSize = nXSize;
2225 13 : double dfYSize = nYSize;
2226 13 : if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
2227 : {
2228 0 : dfXOff = psExtraArgIn->dfXOff;
2229 0 : dfYOff = psExtraArgIn->dfYOff;
2230 0 : dfXSize = psExtraArgIn->dfXSize;
2231 0 : dfYSize = psExtraArgIn->dfYSize;
2232 : }
2233 :
2234 : // The window we will actually request from the source raster band.
2235 13 : double dfReqXOff = 0.0;
2236 13 : double dfReqYOff = 0.0;
2237 13 : double dfReqXSize = 0.0;
2238 13 : double dfReqYSize = 0.0;
2239 13 : int nReqXOff = 0;
2240 13 : int nReqYOff = 0;
2241 13 : int nReqXSize = 0;
2242 13 : int nReqYSize = 0;
2243 :
2244 : // The window we will actual set _within_ the pData buffer.
2245 13 : int nOutXOff = 0;
2246 13 : int nOutYOff = 0;
2247 13 : int nOutXSize = 0;
2248 13 : int nOutYSize = 0;
2249 :
2250 13 : bool bError = false;
2251 13 : if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
2252 : &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
2253 : &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
2254 : &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
2255 : {
2256 0 : return bError ? CE_Failure : CE_None;
2257 : }
2258 :
2259 13 : auto l_band = GetRasterBand();
2260 13 : if (!l_band)
2261 0 : return CE_Failure;
2262 :
2263 : /* -------------------------------------------------------------------- */
2264 : /* Allocate temporary buffer(s). */
2265 : /* -------------------------------------------------------------------- */
2266 13 : const auto eSrcBandDT = l_band->GetRasterDataType();
2267 13 : const int nSrcBandDTSize = GDALGetDataTypeSizeBytes(eSrcBandDT);
2268 13 : const auto eSrcMaskBandDT = l_band->GetMaskBand()->GetRasterDataType();
2269 13 : const int nSrcMaskBandDTSize = GDALGetDataTypeSizeBytes(eSrcMaskBandDT);
2270 13 : double dfRemappedValue = m_dfRemappedValue;
2271 13 : if (!m_bHasRemappedValue)
2272 : {
2273 19 : if (eSrcBandDT == GDT_Byte &&
2274 6 : m_dfNoDataValue >= std::numeric_limits<GByte>::min() &&
2275 25 : m_dfNoDataValue <= std::numeric_limits<GByte>::max() &&
2276 6 : static_cast<int>(m_dfNoDataValue) == m_dfNoDataValue)
2277 : {
2278 6 : if (m_dfNoDataValue == std::numeric_limits<GByte>::max())
2279 1 : dfRemappedValue = m_dfNoDataValue - 1;
2280 : else
2281 5 : dfRemappedValue = m_dfNoDataValue + 1;
2282 : }
2283 10 : else if (eSrcBandDT == GDT_UInt16 &&
2284 3 : m_dfNoDataValue >= std::numeric_limits<uint16_t>::min() &&
2285 13 : m_dfNoDataValue <= std::numeric_limits<uint16_t>::max() &&
2286 3 : static_cast<int>(m_dfNoDataValue) == m_dfNoDataValue)
2287 : {
2288 3 : if (m_dfNoDataValue == std::numeric_limits<uint16_t>::max())
2289 1 : dfRemappedValue = m_dfNoDataValue - 1;
2290 : else
2291 2 : dfRemappedValue = m_dfNoDataValue + 1;
2292 : }
2293 6 : else if (eSrcBandDT == GDT_Int16 &&
2294 2 : m_dfNoDataValue >= std::numeric_limits<int16_t>::min() &&
2295 8 : m_dfNoDataValue <= std::numeric_limits<int16_t>::max() &&
2296 2 : static_cast<int>(m_dfNoDataValue) == m_dfNoDataValue)
2297 : {
2298 2 : if (m_dfNoDataValue == std::numeric_limits<int16_t>::max())
2299 1 : dfRemappedValue = m_dfNoDataValue - 1;
2300 : else
2301 1 : dfRemappedValue = m_dfNoDataValue + 1;
2302 : }
2303 : else
2304 : {
2305 2 : constexpr double EPS = 1e-3;
2306 2 : if (m_dfNoDataValue == 0)
2307 1 : dfRemappedValue = EPS;
2308 : else
2309 1 : dfRemappedValue = m_dfNoDataValue * (1 + EPS);
2310 : }
2311 : }
2312 13 : const bool bByteOptim =
2313 6 : (eSrcBandDT == GDT_Byte && eBufType == GDT_Byte &&
2314 5 : eSrcMaskBandDT == GDT_Byte && m_dfMaskValueThreshold >= 0 &&
2315 5 : m_dfMaskValueThreshold <= 255 &&
2316 5 : static_cast<int>(m_dfMaskValueThreshold) == m_dfMaskValueThreshold &&
2317 4 : m_dfNoDataValue >= 0 && m_dfNoDataValue <= 255 &&
2318 4 : static_cast<int>(m_dfNoDataValue) == m_dfNoDataValue &&
2319 23 : dfRemappedValue >= 0 && dfRemappedValue <= 255 &&
2320 4 : static_cast<int>(dfRemappedValue) == dfRemappedValue);
2321 : GByte *pabyWrkBuffer;
2322 : try
2323 : {
2324 13 : if (bByteOptim && nOutXOff == 0 && nOutYOff == 0 &&
2325 4 : nOutXSize == nBufXSize && nOutYSize == nBufYSize &&
2326 4 : eSrcBandDT == eBufType && nPixelSpace == nSrcBandDTSize &&
2327 4 : nLineSpace == nPixelSpace * nBufXSize)
2328 : {
2329 4 : pabyWrkBuffer = static_cast<GByte *>(pData);
2330 : }
2331 : else
2332 : {
2333 9 : oWorkingState.m_abyWrkBuffer.resize(static_cast<size_t>(nOutXSize) *
2334 9 : nOutYSize * nSrcBandDTSize);
2335 : pabyWrkBuffer =
2336 9 : reinterpret_cast<GByte *>(oWorkingState.m_abyWrkBuffer.data());
2337 : }
2338 13 : oWorkingState.m_abyWrkBufferMask.resize(static_cast<size_t>(nOutXSize) *
2339 13 : nOutYSize * nSrcMaskBandDTSize);
2340 : }
2341 0 : catch (const std::exception &)
2342 : {
2343 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
2344 : "Out of memory when allocating buffers");
2345 0 : return CE_Failure;
2346 : }
2347 :
2348 : /* -------------------------------------------------------------------- */
2349 : /* Load data. */
2350 : /* -------------------------------------------------------------------- */
2351 13 : if (!m_osResampling.empty())
2352 : {
2353 0 : psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
2354 : }
2355 13 : else if (psExtraArgIn != nullptr)
2356 : {
2357 13 : psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
2358 : }
2359 :
2360 13 : psExtraArg->bFloatingPointWindowValidity = TRUE;
2361 13 : psExtraArg->dfXOff = dfReqXOff;
2362 13 : psExtraArg->dfYOff = dfReqYOff;
2363 13 : psExtraArg->dfXSize = dfReqXSize;
2364 13 : psExtraArg->dfYSize = dfReqYSize;
2365 13 : if (psExtraArgIn)
2366 : {
2367 13 : psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
2368 13 : psExtraArg->pProgressData = psExtraArgIn->pProgressData;
2369 : }
2370 :
2371 13 : if (l_band->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
2372 : pabyWrkBuffer, nOutXSize, nOutYSize, eSrcBandDT, 0, 0,
2373 13 : psExtraArg) != CE_None)
2374 : {
2375 0 : return CE_Failure;
2376 : }
2377 :
2378 26 : if (l_band->GetMaskBand()->RasterIO(
2379 : GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
2380 13 : oWorkingState.m_abyWrkBufferMask.data(), nOutXSize, nOutYSize,
2381 13 : eSrcMaskBandDT, 0, 0, psExtraArg) != CE_None)
2382 : {
2383 0 : return CE_Failure;
2384 : }
2385 :
2386 : /* -------------------------------------------------------------------- */
2387 : /* Do the processing. */
2388 : /* -------------------------------------------------------------------- */
2389 :
2390 13 : GByte *const pabyOut = static_cast<GByte *>(pData) +
2391 13 : nPixelSpace * nOutXOff +
2392 13 : static_cast<GPtrDiff_t>(nLineSpace) * nOutYOff;
2393 13 : if (bByteOptim)
2394 : {
2395 : // Special case when everything fits on Byte
2396 4 : const GByte nMaskValueThreshold =
2397 4 : static_cast<GByte>(m_dfMaskValueThreshold);
2398 4 : const GByte nNoDataValue = static_cast<GByte>(m_dfNoDataValue);
2399 4 : const GByte nRemappedValue = static_cast<GByte>(dfRemappedValue);
2400 4 : size_t nSrcIdx = 0;
2401 8 : for (int iY = 0; iY < nOutYSize; iY++)
2402 : {
2403 4 : GSpacing nDstOffset = iY * nLineSpace;
2404 12 : for (int iX = 0; iX < nOutXSize; iX++)
2405 : {
2406 : const GByte nMaskVal =
2407 8 : oWorkingState.m_abyWrkBufferMask[nSrcIdx];
2408 8 : if (nMaskVal <= nMaskValueThreshold)
2409 : {
2410 4 : pabyOut[static_cast<GPtrDiff_t>(nDstOffset)] = nNoDataValue;
2411 : }
2412 : else
2413 : {
2414 4 : if (pabyWrkBuffer[nSrcIdx] == nNoDataValue)
2415 : {
2416 2 : pabyOut[static_cast<GPtrDiff_t>(nDstOffset)] =
2417 : nRemappedValue;
2418 : }
2419 : else
2420 : {
2421 2 : pabyOut[static_cast<GPtrDiff_t>(nDstOffset)] =
2422 2 : pabyWrkBuffer[nSrcIdx];
2423 : }
2424 : }
2425 8 : nDstOffset += nPixelSpace;
2426 8 : nSrcIdx++;
2427 : }
2428 : }
2429 : }
2430 : else
2431 : {
2432 9 : size_t nSrcIdx = 0;
2433 9 : double dfMaskVal = 0;
2434 9 : const int nBufDTSize = GDALGetDataTypeSizeBytes(eBufType);
2435 18 : std::vector<GByte> abyDstNoData(nBufDTSize);
2436 9 : GDALCopyWords(&m_dfNoDataValue, GDT_Float64, 0, abyDstNoData.data(),
2437 : eBufType, 0, 1);
2438 18 : std::vector<GByte> abyRemappedValue(nBufDTSize);
2439 9 : GDALCopyWords(&dfRemappedValue, GDT_Float64, 0, abyRemappedValue.data(),
2440 : eBufType, 0, 1);
2441 18 : for (int iY = 0; iY < nOutYSize; iY++)
2442 : {
2443 9 : GSpacing nDstOffset = iY * nLineSpace;
2444 28 : for (int iX = 0; iX < nOutXSize; iX++)
2445 : {
2446 19 : if (eSrcMaskBandDT == GDT_Byte)
2447 : {
2448 19 : dfMaskVal = oWorkingState.m_abyWrkBufferMask[nSrcIdx];
2449 : }
2450 : else
2451 : {
2452 0 : GDALCopyWords(oWorkingState.m_abyWrkBufferMask.data() +
2453 0 : nSrcIdx * nSrcMaskBandDTSize,
2454 : eSrcMaskBandDT, 0, &dfMaskVal, GDT_Float64, 0,
2455 : 1);
2456 : }
2457 19 : void *const pDst =
2458 19 : pabyOut + static_cast<GPtrDiff_t>(nDstOffset);
2459 19 : if (!(dfMaskVal > m_dfMaskValueThreshold))
2460 : {
2461 9 : memcpy(pDst, abyDstNoData.data(), nBufDTSize);
2462 : }
2463 : else
2464 : {
2465 10 : const void *const pSrc =
2466 10 : pabyWrkBuffer + nSrcIdx * nSrcBandDTSize;
2467 10 : if (eSrcBandDT == eBufType)
2468 : {
2469 : // coverity[overrun-buffer-arg]
2470 8 : memcpy(pDst, pSrc, nBufDTSize);
2471 : }
2472 : else
2473 : {
2474 2 : GDALCopyWords(pSrc, eSrcBandDT, 0, pDst, eBufType, 0,
2475 : 1);
2476 : }
2477 10 : if (memcmp(pDst, abyDstNoData.data(), nBufDTSize) == 0)
2478 9 : memcpy(pDst, abyRemappedValue.data(), nBufDTSize);
2479 : }
2480 19 : nDstOffset += nPixelSpace;
2481 19 : nSrcIdx++;
2482 : }
2483 : }
2484 : }
2485 :
2486 13 : if (psExtraArg->pfnProgress)
2487 0 : psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
2488 :
2489 13 : return CE_None;
2490 : }
2491 :
2492 : /************************************************************************/
2493 : /* GetMinimum() */
2494 : /************************************************************************/
2495 :
2496 0 : double VRTNoDataFromMaskSource::GetMinimum(int /* nXSize */, int /* nYSize */,
2497 : int *pbSuccess)
2498 : {
2499 0 : *pbSuccess = FALSE;
2500 0 : return 0.0;
2501 : }
2502 :
2503 : /************************************************************************/
2504 : /* GetMaximum() */
2505 : /************************************************************************/
2506 :
2507 0 : double VRTNoDataFromMaskSource::GetMaximum(int /* nXSize */, int /* nYSize */,
2508 : int *pbSuccess)
2509 : {
2510 0 : *pbSuccess = FALSE;
2511 0 : return 0.0;
2512 : }
2513 :
2514 : /************************************************************************/
2515 : /* GetHistogram() */
2516 : /************************************************************************/
2517 :
2518 0 : CPLErr VRTNoDataFromMaskSource::GetHistogram(
2519 : int /* nXSize */, int /* nYSize */, double /* dfMin */, double /* dfMax */,
2520 : int /* nBuckets */, GUIntBig * /* panHistogram */,
2521 : int /* bIncludeOutOfRange */, int /* bApproxOK */,
2522 : GDALProgressFunc /* pfnProgress */, void * /* pProgressData */)
2523 : {
2524 0 : return CE_Failure;
2525 : }
2526 :
2527 : /************************************************************************/
2528 : /* ==================================================================== */
2529 : /* VRTComplexSource */
2530 : /* ==================================================================== */
2531 : /************************************************************************/
2532 :
2533 : /************************************************************************/
2534 : /* VRTComplexSource() */
2535 : /************************************************************************/
2536 :
2537 5 : VRTComplexSource::VRTComplexSource(const VRTComplexSource *poSrcSource,
2538 5 : double dfXDstRatio, double dfYDstRatio)
2539 : : VRTSimpleSource(poSrcSource, dfXDstRatio, dfYDstRatio),
2540 5 : m_nProcessingFlags(poSrcSource->m_nProcessingFlags),
2541 5 : m_dfNoDataValue(poSrcSource->m_dfNoDataValue),
2542 5 : m_osNoDataValueOri(poSrcSource->m_osNoDataValueOri),
2543 5 : m_dfScaleOff(poSrcSource->m_dfScaleOff),
2544 5 : m_dfScaleRatio(poSrcSource->m_dfScaleRatio),
2545 5 : m_bSrcMinMaxDefined(poSrcSource->m_bSrcMinMaxDefined),
2546 5 : m_dfSrcMin(poSrcSource->m_dfSrcMin), m_dfSrcMax(poSrcSource->m_dfSrcMax),
2547 5 : m_dfDstMin(poSrcSource->m_dfDstMin), m_dfDstMax(poSrcSource->m_dfDstMax),
2548 5 : m_dfExponent(poSrcSource->m_dfExponent),
2549 5 : m_nColorTableComponent(poSrcSource->m_nColorTableComponent),
2550 5 : m_adfLUTInputs(poSrcSource->m_adfLUTInputs),
2551 5 : m_adfLUTOutputs(poSrcSource->m_adfLUTOutputs)
2552 : {
2553 5 : }
2554 :
2555 : /************************************************************************/
2556 : /* GetTypeStatic() */
2557 : /************************************************************************/
2558 :
2559 13523 : const char *VRTComplexSource::GetTypeStatic()
2560 : {
2561 : static const char *TYPE = "ComplexSource";
2562 13523 : return TYPE;
2563 : }
2564 :
2565 : /************************************************************************/
2566 : /* GetType() */
2567 : /************************************************************************/
2568 :
2569 3843 : const char *VRTComplexSource::GetType() const
2570 : {
2571 3843 : return GetTypeStatic();
2572 : }
2573 :
2574 : /************************************************************************/
2575 : /* SetNoDataValue() */
2576 : /************************************************************************/
2577 :
2578 65 : void VRTComplexSource::SetNoDataValue(double dfNewNoDataValue)
2579 :
2580 : {
2581 65 : if (dfNewNoDataValue == VRT_NODATA_UNSET)
2582 : {
2583 0 : m_nProcessingFlags &= ~PROCESSING_FLAG_NODATA;
2584 0 : m_dfNoDataValue = VRT_NODATA_UNSET;
2585 0 : return;
2586 : }
2587 :
2588 65 : m_nProcessingFlags |= PROCESSING_FLAG_NODATA;
2589 65 : m_dfNoDataValue = dfNewNoDataValue;
2590 : }
2591 :
2592 : /************************************************************************/
2593 : /* GetAdjustedNoDataValue() */
2594 : /************************************************************************/
2595 :
2596 4728 : double VRTComplexSource::GetAdjustedNoDataValue() const
2597 : {
2598 4728 : if ((m_nProcessingFlags & PROCESSING_FLAG_NODATA) != 0)
2599 : {
2600 36 : auto l_band = GetRasterBand();
2601 36 : if (l_band && l_band->GetRasterDataType() == GDT_Float32)
2602 : {
2603 15 : return GDALAdjustNoDataCloseToFloatMax(m_dfNoDataValue);
2604 : }
2605 : }
2606 4713 : return m_dfNoDataValue;
2607 : }
2608 :
2609 : /************************************************************************/
2610 : /* SerializeToXML() */
2611 : /************************************************************************/
2612 :
2613 65 : CPLXMLNode *VRTComplexSource::SerializeToXML(const char *pszVRTPath)
2614 :
2615 : {
2616 65 : CPLXMLNode *psSrc = VRTSimpleSource::SerializeToXML(pszVRTPath);
2617 :
2618 65 : if (psSrc == nullptr)
2619 0 : return nullptr;
2620 :
2621 65 : CPLFree(psSrc->pszValue);
2622 65 : psSrc->pszValue = CPLStrdup(GetTypeStatic());
2623 :
2624 65 : if ((m_nProcessingFlags & PROCESSING_FLAG_USE_MASK_BAND) != 0)
2625 : {
2626 30 : CPLSetXMLValue(psSrc, "UseMaskBand", "true");
2627 : }
2628 :
2629 65 : if ((m_nProcessingFlags & PROCESSING_FLAG_NODATA) != 0)
2630 : {
2631 25 : if (!m_osNoDataValueOri.empty() && GetRasterBandNoOpen() == nullptr)
2632 : {
2633 2 : CPLSetXMLValue(psSrc, "NODATA", m_osNoDataValueOri.c_str());
2634 : }
2635 : else
2636 : {
2637 23 : GDALDataType eBandDT = GDT_Unknown;
2638 23 : double dfNoDataValue = m_dfNoDataValue;
2639 23 : const auto kMaxFloat = std::numeric_limits<float>::max();
2640 23 : if (std::fabs(std::fabs(m_dfNoDataValue) - kMaxFloat) <
2641 : 1e-10 * kMaxFloat)
2642 : {
2643 1 : auto l_band = GetRasterBand();
2644 1 : if (l_band)
2645 : {
2646 1 : dfNoDataValue = GetAdjustedNoDataValue();
2647 1 : eBandDT = l_band->GetRasterDataType();
2648 : }
2649 : }
2650 23 : CPLSetXMLValue(
2651 : psSrc, "NODATA",
2652 46 : VRTSerializeNoData(dfNoDataValue, eBandDT, 18).c_str());
2653 : }
2654 : }
2655 :
2656 65 : if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_LINEAR) != 0)
2657 : {
2658 1 : CPLSetXMLValue(psSrc, "ScaleOffset", CPLSPrintf("%g", m_dfScaleOff));
2659 1 : CPLSetXMLValue(psSrc, "ScaleRatio", CPLSPrintf("%g", m_dfScaleRatio));
2660 : }
2661 64 : else if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_EXPONENTIAL) != 0)
2662 : {
2663 0 : CPLSetXMLValue(psSrc, "Exponent", CPLSPrintf("%g", m_dfExponent));
2664 0 : if (m_bSrcMinMaxDefined)
2665 : {
2666 0 : CPLSetXMLValue(psSrc, "SrcMin", CPLSPrintf("%g", m_dfSrcMin));
2667 0 : CPLSetXMLValue(psSrc, "SrcMax", CPLSPrintf("%g", m_dfSrcMax));
2668 : }
2669 0 : CPLSetXMLValue(psSrc, "DstMin", CPLSPrintf("%g", m_dfDstMin));
2670 0 : CPLSetXMLValue(psSrc, "DstMax", CPLSPrintf("%g", m_dfDstMax));
2671 : }
2672 :
2673 65 : if (!m_adfLUTInputs.empty())
2674 : {
2675 : // Make sure we print with sufficient precision to address really close
2676 : // entries (#6422).
2677 0 : CPLString osLUT;
2678 0 : if (m_adfLUTInputs.size() >= 2 &&
2679 0 : CPLString().Printf("%g", m_adfLUTInputs[0]) ==
2680 0 : CPLString().Printf("%g", m_adfLUTInputs[1]))
2681 : {
2682 0 : osLUT = CPLString().Printf("%.17g:%g", m_adfLUTInputs[0],
2683 0 : m_adfLUTOutputs[0]);
2684 : }
2685 : else
2686 : {
2687 0 : osLUT = CPLString().Printf("%g:%g", m_adfLUTInputs[0],
2688 0 : m_adfLUTOutputs[0]);
2689 : }
2690 0 : for (size_t i = 1; i < m_adfLUTInputs.size(); i++)
2691 : {
2692 0 : if (CPLString().Printf("%g", m_adfLUTInputs[i]) ==
2693 0 : CPLString().Printf("%g", m_adfLUTInputs[i - 1]) ||
2694 0 : (i + 1 < m_adfLUTInputs.size() &&
2695 0 : CPLString().Printf("%g", m_adfLUTInputs[i]) ==
2696 0 : CPLString().Printf("%g", m_adfLUTInputs[i + 1])))
2697 : {
2698 : // TODO(schwehr): An explanation of the 18 would be helpful.
2699 : // Can someone distill the issue down to a quick comment?
2700 : // https://trac.osgeo.org/gdal/ticket/6422
2701 0 : osLUT += CPLString().Printf(",%.17g:%g", m_adfLUTInputs[i],
2702 0 : m_adfLUTOutputs[i]);
2703 : }
2704 : else
2705 : {
2706 0 : osLUT += CPLString().Printf(",%g:%g", m_adfLUTInputs[i],
2707 0 : m_adfLUTOutputs[i]);
2708 : }
2709 : }
2710 0 : CPLSetXMLValue(psSrc, "LUT", osLUT);
2711 : }
2712 :
2713 65 : if (m_nColorTableComponent)
2714 : {
2715 7 : CPLSetXMLValue(psSrc, "ColorTableComponent",
2716 : CPLSPrintf("%d", m_nColorTableComponent));
2717 : }
2718 :
2719 65 : return psSrc;
2720 : }
2721 :
2722 : /************************************************************************/
2723 : /* XMLInit() */
2724 : /************************************************************************/
2725 :
2726 217 : CPLErr VRTComplexSource::XMLInit(const CPLXMLNode *psSrc,
2727 : const char *pszVRTPath,
2728 : VRTMapSharedResources &oMapSharedSources)
2729 :
2730 : {
2731 : /* -------------------------------------------------------------------- */
2732 : /* Do base initialization. */
2733 : /* -------------------------------------------------------------------- */
2734 : {
2735 : const CPLErr eErr =
2736 217 : VRTSimpleSource::XMLInit(psSrc, pszVRTPath, oMapSharedSources);
2737 217 : if (eErr != CE_None)
2738 0 : return eErr;
2739 : }
2740 :
2741 : /* -------------------------------------------------------------------- */
2742 : /* Complex parameters. */
2743 : /* -------------------------------------------------------------------- */
2744 217 : const char *pszScaleOffset = CPLGetXMLValue(psSrc, "ScaleOffset", nullptr);
2745 217 : const char *pszScaleRatio = CPLGetXMLValue(psSrc, "ScaleRatio", nullptr);
2746 217 : if (pszScaleOffset || pszScaleRatio)
2747 : {
2748 18 : m_nProcessingFlags |= PROCESSING_FLAG_SCALING_LINEAR;
2749 18 : if (pszScaleOffset)
2750 18 : m_dfScaleOff = CPLAtof(pszScaleOffset);
2751 18 : if (pszScaleRatio)
2752 15 : m_dfScaleRatio = CPLAtof(pszScaleRatio);
2753 : }
2754 199 : else if (CPLGetXMLValue(psSrc, "Exponent", nullptr) != nullptr &&
2755 200 : CPLGetXMLValue(psSrc, "DstMin", nullptr) != nullptr &&
2756 1 : CPLGetXMLValue(psSrc, "DstMax", nullptr) != nullptr)
2757 : {
2758 1 : m_nProcessingFlags |= PROCESSING_FLAG_SCALING_EXPONENTIAL;
2759 1 : m_dfExponent = CPLAtof(CPLGetXMLValue(psSrc, "Exponent", "1.0"));
2760 :
2761 1 : const char *pszSrcMin = CPLGetXMLValue(psSrc, "SrcMin", nullptr);
2762 1 : const char *pszSrcMax = CPLGetXMLValue(psSrc, "SrcMax", nullptr);
2763 1 : if (pszSrcMin && pszSrcMax)
2764 : {
2765 0 : m_dfSrcMin = CPLAtof(pszSrcMin);
2766 0 : m_dfSrcMax = CPLAtof(pszSrcMax);
2767 0 : m_bSrcMinMaxDefined = true;
2768 : }
2769 :
2770 1 : m_dfDstMin = CPLAtof(CPLGetXMLValue(psSrc, "DstMin", "0.0"));
2771 1 : m_dfDstMax = CPLAtof(CPLGetXMLValue(psSrc, "DstMax", "0.0"));
2772 : }
2773 :
2774 217 : if (const char *pszNODATA = CPLGetXMLValue(psSrc, "NODATA", nullptr))
2775 : {
2776 50 : m_nProcessingFlags |= PROCESSING_FLAG_NODATA;
2777 50 : m_osNoDataValueOri = pszNODATA;
2778 50 : m_dfNoDataValue = CPLAtofM(m_osNoDataValueOri.c_str());
2779 : }
2780 :
2781 217 : const char *pszUseMaskBand = CPLGetXMLValue(psSrc, "UseMaskBand", nullptr);
2782 217 : if (pszUseMaskBand && CPLTestBool(pszUseMaskBand))
2783 : {
2784 38 : m_nProcessingFlags |= PROCESSING_FLAG_USE_MASK_BAND;
2785 : }
2786 :
2787 217 : const char *pszLUT = CPLGetXMLValue(psSrc, "LUT", nullptr);
2788 217 : if (pszLUT)
2789 : {
2790 : const CPLStringList aosValues(
2791 59 : CSLTokenizeString2(pszLUT, ",:", CSLT_ALLOWEMPTYTOKENS));
2792 :
2793 59 : const int nLUTItemCount = aosValues.size() / 2;
2794 : try
2795 : {
2796 59 : m_adfLUTInputs.resize(nLUTItemCount);
2797 59 : m_adfLUTOutputs.resize(nLUTItemCount);
2798 : }
2799 0 : catch (const std::bad_alloc &e)
2800 : {
2801 0 : CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
2802 0 : m_adfLUTInputs.clear();
2803 0 : m_adfLUTOutputs.clear();
2804 0 : return CE_Failure;
2805 : }
2806 :
2807 544 : for (int nIndex = 0; nIndex < nLUTItemCount; nIndex++)
2808 : {
2809 485 : m_adfLUTInputs[nIndex] = CPLAtof(aosValues[nIndex * 2]);
2810 485 : m_adfLUTOutputs[nIndex] = CPLAtof(aosValues[nIndex * 2 + 1]);
2811 :
2812 : // Enforce the requirement that the LUT input array is
2813 : // monotonically non-decreasing.
2814 485 : if (std::isnan(m_adfLUTInputs[nIndex]) && nIndex != 0)
2815 : {
2816 0 : CPLError(CE_Failure, CPLE_AppDefined,
2817 : "A Not-A-Number (NaN) source value should be the "
2818 : "first one of the LUT.");
2819 0 : m_adfLUTInputs.clear();
2820 0 : m_adfLUTOutputs.clear();
2821 0 : return CE_Failure;
2822 : }
2823 911 : else if (nIndex > 0 &&
2824 426 : m_adfLUTInputs[nIndex] < m_adfLUTInputs[nIndex - 1])
2825 : {
2826 0 : CPLError(CE_Failure, CPLE_AppDefined,
2827 : "Source values of the LUT are not listed in a "
2828 : "monotonically non-decreasing order");
2829 0 : m_adfLUTInputs.clear();
2830 0 : m_adfLUTOutputs.clear();
2831 0 : return CE_Failure;
2832 : }
2833 : }
2834 :
2835 59 : m_nProcessingFlags |= PROCESSING_FLAG_LUT;
2836 : }
2837 :
2838 : const char *pszColorTableComponent =
2839 217 : CPLGetXMLValue(psSrc, "ColorTableComponent", nullptr);
2840 217 : if (pszColorTableComponent)
2841 : {
2842 15 : m_nColorTableComponent = atoi(pszColorTableComponent);
2843 15 : m_nProcessingFlags |= PROCESSING_FLAG_COLOR_TABLE_EXPANSION;
2844 : }
2845 :
2846 217 : return CE_None;
2847 : }
2848 :
2849 : /************************************************************************/
2850 : /* LookupValue() */
2851 : /************************************************************************/
2852 :
2853 583336 : double VRTComplexSource::LookupValue(double dfInput)
2854 : {
2855 583336 : auto beginIter = m_adfLUTInputs.begin();
2856 583336 : auto endIter = m_adfLUTInputs.end();
2857 583336 : size_t offset = 0;
2858 583336 : if (std::isnan(m_adfLUTInputs[0]))
2859 : {
2860 6 : if (std::isnan(dfInput) || m_adfLUTInputs.size() == 1)
2861 1 : return m_adfLUTOutputs[0];
2862 5 : ++beginIter;
2863 5 : offset = 1;
2864 : }
2865 :
2866 : // Find the index of the first element in the LUT input array that
2867 : // is not smaller than the input value.
2868 : const size_t i =
2869 : offset +
2870 583335 : std::distance(beginIter, std::lower_bound(beginIter, endIter, dfInput));
2871 :
2872 583335 : if (i == offset)
2873 129039 : return m_adfLUTOutputs[offset];
2874 :
2875 : // If the index is beyond the end of the LUT input array, the input
2876 : // value is larger than all the values in the array.
2877 454296 : if (i == m_adfLUTInputs.size())
2878 8 : return m_adfLUTOutputs.back();
2879 :
2880 454288 : if (m_adfLUTInputs[i] == dfInput)
2881 179297 : return m_adfLUTOutputs[i];
2882 :
2883 : // Otherwise, interpolate.
2884 274991 : return m_adfLUTOutputs[i - 1] +
2885 274991 : (dfInput - m_adfLUTInputs[i - 1]) *
2886 274991 : ((m_adfLUTOutputs[i] - m_adfLUTOutputs[i - 1]) /
2887 274991 : (m_adfLUTInputs[i] - m_adfLUTInputs[i - 1]));
2888 : }
2889 :
2890 : /************************************************************************/
2891 : /* SetLinearScaling() */
2892 : /************************************************************************/
2893 :
2894 96 : void VRTComplexSource::SetLinearScaling(double dfOffset, double dfScale)
2895 : {
2896 96 : m_nProcessingFlags &= ~PROCESSING_FLAG_SCALING_EXPONENTIAL;
2897 96 : m_nProcessingFlags |= PROCESSING_FLAG_SCALING_LINEAR;
2898 96 : m_dfScaleOff = dfOffset;
2899 96 : m_dfScaleRatio = dfScale;
2900 96 : }
2901 :
2902 : /************************************************************************/
2903 : /* SetPowerScaling() */
2904 : /************************************************************************/
2905 :
2906 8 : void VRTComplexSource::SetPowerScaling(double dfExponentIn, double dfSrcMinIn,
2907 : double dfSrcMaxIn, double dfDstMinIn,
2908 : double dfDstMaxIn)
2909 : {
2910 8 : m_nProcessingFlags &= ~PROCESSING_FLAG_SCALING_LINEAR;
2911 8 : m_nProcessingFlags |= PROCESSING_FLAG_SCALING_EXPONENTIAL;
2912 8 : m_dfExponent = dfExponentIn;
2913 8 : m_dfSrcMin = dfSrcMinIn;
2914 8 : m_dfSrcMax = dfSrcMaxIn;
2915 8 : m_dfDstMin = dfDstMinIn;
2916 8 : m_dfDstMax = dfDstMaxIn;
2917 8 : m_bSrcMinMaxDefined = true;
2918 8 : }
2919 :
2920 : /************************************************************************/
2921 : /* SetColorTableComponent() */
2922 : /************************************************************************/
2923 :
2924 188 : void VRTComplexSource::SetColorTableComponent(int nComponent)
2925 : {
2926 188 : m_nProcessingFlags |= PROCESSING_FLAG_COLOR_TABLE_EXPANSION;
2927 188 : m_nColorTableComponent = nComponent;
2928 188 : }
2929 :
2930 : /************************************************************************/
2931 : /* RasterIO() */
2932 : /************************************************************************/
2933 :
2934 15705 : CPLErr VRTComplexSource::RasterIO(GDALDataType eVRTBandDataType, int nXOff,
2935 : int nYOff, int nXSize, int nYSize,
2936 : void *pData, int nBufXSize, int nBufYSize,
2937 : GDALDataType eBufType, GSpacing nPixelSpace,
2938 : GSpacing nLineSpace,
2939 : GDALRasterIOExtraArg *psExtraArgIn,
2940 : WorkingState &oWorkingState)
2941 :
2942 : {
2943 : GDALRasterIOExtraArg sExtraArg;
2944 15705 : INIT_RASTERIO_EXTRA_ARG(sExtraArg);
2945 15705 : GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
2946 :
2947 15705 : double dfXOff = nXOff;
2948 15705 : double dfYOff = nYOff;
2949 15705 : double dfXSize = nXSize;
2950 15705 : double dfYSize = nYSize;
2951 15705 : if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
2952 : {
2953 87 : dfXOff = psExtraArgIn->dfXOff;
2954 87 : dfYOff = psExtraArgIn->dfYOff;
2955 87 : dfXSize = psExtraArgIn->dfXSize;
2956 87 : dfYSize = psExtraArgIn->dfYSize;
2957 : }
2958 :
2959 : // The window we will actually request from the source raster band.
2960 15705 : double dfReqXOff = 0.0;
2961 15705 : double dfReqYOff = 0.0;
2962 15705 : double dfReqXSize = 0.0;
2963 15705 : double dfReqYSize = 0.0;
2964 15705 : int nReqXOff = 0;
2965 15705 : int nReqYOff = 0;
2966 15705 : int nReqXSize = 0;
2967 15705 : int nReqYSize = 0;
2968 :
2969 : // The window we will actual set _within_ the pData buffer.
2970 15705 : int nOutXOff = 0;
2971 15705 : int nOutYOff = 0;
2972 15705 : int nOutXSize = 0;
2973 15705 : int nOutYSize = 0;
2974 :
2975 15705 : bool bError = false;
2976 15705 : if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
2977 : &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
2978 : &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
2979 : &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
2980 : {
2981 10946 : return bError ? CE_Failure : CE_None;
2982 : }
2983 : #if DEBUG_VERBOSE
2984 : CPLDebug("VRT",
2985 : "nXOff=%d, nYOff=%d, nXSize=%d, nYSize=%d, nBufXSize=%d, "
2986 : "nBufYSize=%d,\n"
2987 : "dfReqXOff=%g, dfReqYOff=%g, dfReqXSize=%g, dfReqYSize=%g,\n"
2988 : "nReqXOff=%d, nReqYOff=%d, nReqXSize=%d, nReqYSize=%d,\n"
2989 : "nOutXOff=%d, nOutYOff=%d, nOutXSize=%d, nOutYSize=%d",
2990 : nXOff, nYOff, nXSize, nYSize, nBufXSize, nBufYSize, dfReqXOff,
2991 : dfReqYOff, dfReqXSize, dfReqYSize, nReqXOff, nReqYOff, nReqXSize,
2992 : nReqYSize, nOutXOff, nOutYOff, nOutXSize, nOutYSize);
2993 : #endif
2994 :
2995 4759 : auto poSourceBand = GetRasterBand();
2996 4759 : if (!poSourceBand)
2997 0 : return CE_Failure;
2998 :
2999 4759 : if (!m_osResampling.empty())
3000 : {
3001 28 : psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
3002 : }
3003 4731 : else if (psExtraArgIn != nullptr)
3004 : {
3005 4731 : psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
3006 : }
3007 4759 : psExtraArg->bFloatingPointWindowValidity = TRUE;
3008 4759 : psExtraArg->dfXOff = dfReqXOff;
3009 4759 : psExtraArg->dfYOff = dfReqYOff;
3010 4759 : psExtraArg->dfXSize = dfReqXSize;
3011 4759 : psExtraArg->dfYSize = dfReqYSize;
3012 4759 : if (psExtraArgIn)
3013 : {
3014 4759 : psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
3015 4759 : psExtraArg->pProgressData = psExtraArgIn->pProgressData;
3016 : }
3017 :
3018 4759 : GByte *const pabyOut = static_cast<GByte *>(pData) +
3019 4759 : nPixelSpace * nOutXOff +
3020 4759 : static_cast<GPtrDiff_t>(nLineSpace) * nOutYOff;
3021 4759 : const auto eSourceType = poSourceBand->GetRasterDataType();
3022 4759 : if (m_nProcessingFlags == PROCESSING_FLAG_NODATA)
3023 : {
3024 : // Optimization if doing only nodata processing
3025 65 : if (eSourceType == GDT_Byte)
3026 : {
3027 36 : if (!GDALIsValueInRange<GByte>(m_dfNoDataValue))
3028 : {
3029 1 : return VRTSimpleSource::RasterIO(
3030 : eVRTBandDataType, nXOff, nYOff, nXSize, nYSize, pData,
3031 : nBufXSize, nBufYSize, eBufType, nPixelSpace, nLineSpace,
3032 1 : psExtraArgIn, oWorkingState);
3033 : }
3034 :
3035 35 : return RasterIOProcessNoData<GByte, GDT_Byte>(
3036 : poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3037 : nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3038 35 : nLineSpace, psExtraArg, oWorkingState);
3039 : }
3040 29 : else if (eSourceType == GDT_Int16)
3041 : {
3042 2 : if (!GDALIsValueInRange<int16_t>(m_dfNoDataValue))
3043 : {
3044 1 : return VRTSimpleSource::RasterIO(
3045 : eVRTBandDataType, nXOff, nYOff, nXSize, nYSize, pData,
3046 : nBufXSize, nBufYSize, eBufType, nPixelSpace, nLineSpace,
3047 1 : psExtraArgIn, oWorkingState);
3048 : }
3049 :
3050 1 : return RasterIOProcessNoData<int16_t, GDT_Int16>(
3051 : poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3052 : nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3053 1 : nLineSpace, psExtraArg, oWorkingState);
3054 : }
3055 27 : else if (eSourceType == GDT_UInt16)
3056 : {
3057 8 : if (!GDALIsValueInRange<uint16_t>(m_dfNoDataValue))
3058 : {
3059 1 : return VRTSimpleSource::RasterIO(
3060 : eVRTBandDataType, nXOff, nYOff, nXSize, nYSize, pData,
3061 : nBufXSize, nBufYSize, eBufType, nPixelSpace, nLineSpace,
3062 1 : psExtraArgIn, oWorkingState);
3063 : }
3064 :
3065 7 : return RasterIOProcessNoData<uint16_t, GDT_UInt16>(
3066 : poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3067 : nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3068 7 : nLineSpace, psExtraArg, oWorkingState);
3069 : }
3070 : }
3071 :
3072 4713 : const bool bIsComplex = CPL_TO_BOOL(GDALDataTypeIsComplex(eBufType));
3073 : CPLErr eErr;
3074 : // For Int32, float32 isn't sufficiently precise as working data type
3075 4713 : if (eVRTBandDataType == GDT_CInt32 || eVRTBandDataType == GDT_CFloat64 ||
3076 4711 : eVRTBandDataType == GDT_Int32 || eVRTBandDataType == GDT_UInt32 ||
3077 4709 : eVRTBandDataType == GDT_Int64 || eVRTBandDataType == GDT_UInt64 ||
3078 4690 : eVRTBandDataType == GDT_Float64 || eSourceType == GDT_Int32 ||
3079 4689 : eSourceType == GDT_UInt32 || eSourceType == GDT_Int64 ||
3080 : eSourceType == GDT_UInt64)
3081 : {
3082 25 : eErr = RasterIOInternal<double>(
3083 : poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3084 : nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3085 : nLineSpace, psExtraArg, bIsComplex ? GDT_CFloat64 : GDT_Float64,
3086 : oWorkingState);
3087 : }
3088 : else
3089 : {
3090 4688 : eErr = RasterIOInternal<float>(
3091 : poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3092 : nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3093 : nLineSpace, psExtraArg, bIsComplex ? GDT_CFloat32 : GDT_Float32,
3094 : oWorkingState);
3095 : }
3096 :
3097 4713 : if (psExtraArg->pfnProgress)
3098 44 : psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
3099 :
3100 4713 : return eErr;
3101 : }
3102 :
3103 : /************************************************************************/
3104 : /* hasZeroByte() */
3105 : /************************************************************************/
3106 :
3107 : CPL_NOSANITIZE_UNSIGNED_INT_OVERFLOW
3108 3863300 : static inline bool hasZeroByte(uint32_t v)
3109 : {
3110 : // Cf https://graphics.stanford.edu/~seander/bithacks.html#ZeroInWord
3111 3863300 : return (((v)-0x01010101U) & ~(v)&0x80808080U) != 0;
3112 : }
3113 :
3114 : /************************************************************************/
3115 : /* CopyWord() */
3116 : /************************************************************************/
3117 :
3118 : template <class SrcType>
3119 6418944 : static void CopyWord(const SrcType *pSrcVal, GDALDataType eSrcType,
3120 : void *pDstVal, GDALDataType eDstType)
3121 : {
3122 6418944 : switch (eDstType)
3123 : {
3124 2145060 : case GDT_Byte:
3125 2145060 : GDALCopyWord(*pSrcVal, *static_cast<uint8_t *>(pDstVal));
3126 2145060 : break;
3127 80 : case GDT_Int8:
3128 80 : GDALCopyWord(*pSrcVal, *static_cast<int8_t *>(pDstVal));
3129 80 : break;
3130 1348720 : case GDT_UInt16:
3131 1348720 : GDALCopyWord(*pSrcVal, *static_cast<uint16_t *>(pDstVal));
3132 1348720 : break;
3133 2920180 : case GDT_Int16:
3134 2920180 : GDALCopyWord(*pSrcVal, *static_cast<int16_t *>(pDstVal));
3135 2920180 : break;
3136 20 : case GDT_UInt32:
3137 20 : GDALCopyWord(*pSrcVal, *static_cast<uint32_t *>(pDstVal));
3138 20 : break;
3139 638 : case GDT_Int32:
3140 638 : GDALCopyWord(*pSrcVal, *static_cast<int32_t *>(pDstVal));
3141 638 : break;
3142 20 : case GDT_UInt64:
3143 20 : GDALCopyWord(*pSrcVal, *static_cast<uint64_t *>(pDstVal));
3144 20 : break;
3145 40 : case GDT_Int64:
3146 40 : GDALCopyWord(*pSrcVal, *static_cast<int64_t *>(pDstVal));
3147 40 : break;
3148 1650 : case GDT_Float32:
3149 1650 : GDALCopyWord(*pSrcVal, *static_cast<float *>(pDstVal));
3150 1650 : break;
3151 1525 : case GDT_Float64:
3152 1525 : GDALCopyWord(*pSrcVal, *static_cast<double *>(pDstVal));
3153 1525 : break;
3154 1016 : default:
3155 1016 : GDALCopyWords(pSrcVal, eSrcType, 0, pDstVal, eDstType, 0, 1);
3156 1016 : break;
3157 : }
3158 6418944 : }
3159 :
3160 : /************************************************************************/
3161 : /* RasterIOProcessNoData() */
3162 : /************************************************************************/
3163 :
3164 : // This method is an optimization of the generic RasterIOInternal()
3165 : // that deals with a VRTComplexSource with only a NODATA value in it, and
3166 : // no other processing flags.
3167 :
3168 : // nReqXOff, nReqYOff, nReqXSize, nReqYSize are expressed in source band
3169 : // referential.
3170 : template <class SourceDT, GDALDataType eSourceType>
3171 43 : CPLErr VRTComplexSource::RasterIOProcessNoData(
3172 : GDALRasterBand *poSourceBand, GDALDataType eVRTBandDataType, int nReqXOff,
3173 : int nReqYOff, int nReqXSize, int nReqYSize, void *pData, int nOutXSize,
3174 : int nOutYSize, GDALDataType eBufType, GSpacing nPixelSpace,
3175 : GSpacing nLineSpace, GDALRasterIOExtraArg *psExtraArg,
3176 : WorkingState &oWorkingState)
3177 : {
3178 43 : CPLAssert(m_nProcessingFlags == PROCESSING_FLAG_NODATA);
3179 43 : CPLAssert(GDALIsValueInRange<SourceDT>(m_dfNoDataValue));
3180 :
3181 : /* -------------------------------------------------------------------- */
3182 : /* Read into a temporary buffer. */
3183 : /* -------------------------------------------------------------------- */
3184 : try
3185 : {
3186 : // Cannot overflow since pData should at least have that number of
3187 : // elements
3188 43 : const size_t nPixelCount = static_cast<size_t>(nOutXSize) * nOutYSize;
3189 43 : if (nPixelCount >
3190 43 : static_cast<size_t>(std::numeric_limits<ptrdiff_t>::max()) /
3191 : sizeof(SourceDT))
3192 : {
3193 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
3194 : "Too large temporary buffer");
3195 0 : return CE_Failure;
3196 : }
3197 43 : oWorkingState.m_abyWrkBuffer.resize(sizeof(SourceDT) * nPixelCount);
3198 : }
3199 0 : catch (const std::bad_alloc &e)
3200 : {
3201 0 : CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
3202 0 : return CE_Failure;
3203 : }
3204 : const auto paSrcData =
3205 43 : reinterpret_cast<const SourceDT *>(oWorkingState.m_abyWrkBuffer.data());
3206 :
3207 43 : const GDALRIOResampleAlg eResampleAlgBack = psExtraArg->eResampleAlg;
3208 43 : if (!m_osResampling.empty())
3209 : {
3210 0 : psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
3211 : }
3212 :
3213 43 : const CPLErr eErr = poSourceBand->RasterIO(
3214 : GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
3215 43 : oWorkingState.m_abyWrkBuffer.data(), nOutXSize, nOutYSize, eSourceType,
3216 8 : sizeof(SourceDT), sizeof(SourceDT) * static_cast<GSpacing>(nOutXSize),
3217 : psExtraArg);
3218 43 : if (!m_osResampling.empty())
3219 0 : psExtraArg->eResampleAlg = eResampleAlgBack;
3220 :
3221 43 : if (eErr != CE_None)
3222 : {
3223 0 : return eErr;
3224 : }
3225 :
3226 43 : const auto nNoDataValue = static_cast<SourceDT>(m_dfNoDataValue);
3227 43 : size_t idxBuffer = 0;
3228 74 : if (eSourceType == eBufType &&
3229 31 : !GDALDataTypeIsConversionLossy(eSourceType, eVRTBandDataType))
3230 : {
3231 : // Most optimized case: the output type is the same as the source type,
3232 : // and conversion from the source type to the VRT band data type is
3233 : // not lossy
3234 14 : for (int iY = 0; iY < nOutYSize; iY++)
3235 : {
3236 19591 : GByte *pDstLocation = static_cast<GByte *>(pData) +
3237 19591 : static_cast<GPtrDiff_t>(nLineSpace) * iY;
3238 :
3239 19591 : int iX = 0;
3240 19579 : if (sizeof(SourceDT) == 1 && nPixelSpace == 1)
3241 : {
3242 : // Optimization to detect more quickly if source pixels are
3243 : // at nodata.
3244 19583 : const GByte byNoDataValue = static_cast<GByte>(nNoDataValue);
3245 19583 : const uint32_t wordNoData =
3246 19583 : (static_cast<uint32_t>(byNoDataValue) << 24) |
3247 19583 : (byNoDataValue << 16) | (byNoDataValue << 8) |
3248 19583 : byNoDataValue;
3249 :
3250 : // Warning: hasZeroByte() assumes WORD_SIZE = 4
3251 19583 : constexpr int WORD_SIZE = 4;
3252 3918250 : for (; iX < nOutXSize - (WORD_SIZE - 1); iX += WORD_SIZE)
3253 : {
3254 : uint32_t v;
3255 : static_assert(sizeof(v) == WORD_SIZE,
3256 : "sizeof(v) == WORD_SIZE");
3257 3950310 : memcpy(&v, paSrcData + idxBuffer, sizeof(v));
3258 : // Cf https://graphics.stanford.edu/~seander/bithacks.html#ValueInWord
3259 3950310 : if (!hasZeroByte(v ^ wordNoData))
3260 : {
3261 : // No bytes are at nodata
3262 3789080 : memcpy(pDstLocation, &v, WORD_SIZE);
3263 3789080 : idxBuffer += WORD_SIZE;
3264 3789080 : pDstLocation += WORD_SIZE;
3265 : }
3266 109587 : else if (v == wordNoData)
3267 : {
3268 : // All bytes are at nodata
3269 131075 : idxBuffer += WORD_SIZE;
3270 131075 : pDstLocation += WORD_SIZE;
3271 : }
3272 : else
3273 : {
3274 : // There are both bytes at nodata and valid bytes
3275 0 : for (int k = 0; k < WORD_SIZE; ++k)
3276 : {
3277 48 : if (paSrcData[idxBuffer] != nNoDataValue)
3278 : {
3279 36 : memcpy(pDstLocation, &paSrcData[idxBuffer],
3280 : sizeof(SourceDT));
3281 : }
3282 48 : idxBuffer++;
3283 48 : pDstLocation += nPixelSpace;
3284 : }
3285 : }
3286 : }
3287 : }
3288 :
3289 72 : for (; iX < nOutXSize;
3290 3709 : iX++, pDstLocation += nPixelSpace, idxBuffer++)
3291 : {
3292 3709 : if (paSrcData[idxBuffer] != nNoDataValue)
3293 : {
3294 3679 : memcpy(pDstLocation, &paSrcData[idxBuffer],
3295 : sizeof(SourceDT));
3296 : }
3297 : }
3298 : }
3299 : }
3300 13 : else if (!GDALDataTypeIsConversionLossy(eSourceType, eVRTBandDataType))
3301 : {
3302 : // Conversion from the source type to the VRT band data type is
3303 : // not lossy, so we can directly convert from the source type to
3304 : // the the output type
3305 106 : for (int iY = 0; iY < nOutYSize; iY++)
3306 : {
3307 94 : GByte *pDstLocation = static_cast<GByte *>(pData) +
3308 94 : static_cast<GPtrDiff_t>(nLineSpace) * iY;
3309 :
3310 1556 : for (int iX = 0; iX < nOutXSize;
3311 1462 : iX++, pDstLocation += nPixelSpace, idxBuffer++)
3312 : {
3313 1462 : if (paSrcData[idxBuffer] != nNoDataValue)
3314 : {
3315 658 : CopyWord(&paSrcData[idxBuffer], eSourceType, pDstLocation,
3316 : eBufType);
3317 : }
3318 : }
3319 : }
3320 : }
3321 : else
3322 : {
3323 : GByte abyTemp[2 * sizeof(double)];
3324 7 : for (int iY = 0; iY < nOutYSize; iY++)
3325 : {
3326 6 : GByte *pDstLocation = static_cast<GByte *>(pData) +
3327 6 : static_cast<GPtrDiff_t>(nLineSpace) * iY;
3328 :
3329 36 : for (int iX = 0; iX < nOutXSize;
3330 30 : iX++, pDstLocation += nPixelSpace, idxBuffer++)
3331 : {
3332 30 : if (paSrcData[idxBuffer] != nNoDataValue)
3333 : {
3334 : // Convert first to the VRTRasterBand data type
3335 : // to get its clamping, before outputting to buffer data type
3336 20 : CopyWord(&paSrcData[idxBuffer], eSourceType, abyTemp,
3337 : eVRTBandDataType);
3338 20 : GDALCopyWords(abyTemp, eVRTBandDataType, 0, pDstLocation,
3339 : eBufType, 0, 1);
3340 : }
3341 : }
3342 : }
3343 : }
3344 :
3345 8 : return CE_None;
3346 : }
3347 :
3348 : /************************************************************************/
3349 : /* RasterIOInternal() */
3350 : /************************************************************************/
3351 :
3352 : // nReqXOff, nReqYOff, nReqXSize, nReqYSize are expressed in source band
3353 : // referential.
3354 : template <class WorkingDT>
3355 4727 : CPLErr VRTComplexSource::RasterIOInternal(
3356 : GDALRasterBand *poSourceBand, GDALDataType eVRTBandDataType, int nReqXOff,
3357 : int nReqYOff, int nReqXSize, int nReqYSize, void *pData, int nOutXSize,
3358 : int nOutYSize, GDALDataType eBufType, GSpacing nPixelSpace,
3359 : GSpacing nLineSpace, GDALRasterIOExtraArg *psExtraArg,
3360 : GDALDataType eWrkDataType, WorkingState &oWorkingState)
3361 : {
3362 4727 : const GDALColorTable *poColorTable = nullptr;
3363 4727 : const bool bIsComplex = CPL_TO_BOOL(GDALDataTypeIsComplex(eBufType));
3364 4727 : const int nWordSize = GDALGetDataTypeSizeBytes(eWrkDataType);
3365 4727 : assert(nWordSize != 0);
3366 :
3367 : // If no explicit <NODATA> is set, but UseMaskBand is set, and the band
3368 : // has a nodata value, then use it as if it was set as <NODATA>
3369 4727 : int bNoDataSet = (m_nProcessingFlags & PROCESSING_FLAG_NODATA) != 0;
3370 4727 : double dfNoDataValue = GetAdjustedNoDataValue();
3371 :
3372 4830 : if ((m_nProcessingFlags & PROCESSING_FLAG_USE_MASK_BAND) != 0 &&
3373 103 : poSourceBand->GetMaskFlags() == GMF_NODATA)
3374 : {
3375 0 : dfNoDataValue = poSourceBand->GetNoDataValue(&bNoDataSet);
3376 : }
3377 :
3378 4727 : const bool bNoDataSetIsNan = bNoDataSet && std::isnan(dfNoDataValue);
3379 4727 : const bool bNoDataSetAndNotNan =
3380 4754 : bNoDataSet && !std::isnan(dfNoDataValue) &&
3381 27 : GDALIsValueInRange<WorkingDT>(dfNoDataValue);
3382 4727 : const auto fWorkingDataTypeNoData = static_cast<WorkingDT>(dfNoDataValue);
3383 :
3384 4727 : const GByte *pabyMask = nullptr;
3385 4727 : const WorkingDT *pafData = nullptr;
3386 4727 : if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_LINEAR) != 0 &&
3387 4431 : m_dfScaleRatio == 0 && bNoDataSet == FALSE &&
3388 4041 : (m_nProcessingFlags & PROCESSING_FLAG_USE_MASK_BAND) == 0)
3389 : {
3390 : /* ------------------------------------------------------------------ */
3391 : /* Optimization when writing a constant value */
3392 : /* (used by the -addalpha option of gdalbuildvrt) */
3393 : /* ------------------------------------------------------------------ */
3394 : // Already set to NULL when defined.
3395 : // pafData = NULL;
3396 : }
3397 : else
3398 : {
3399 : /* ---------------------------------------------------------------- */
3400 : /* Read into a temporary buffer. */
3401 : /* ---------------------------------------------------------------- */
3402 686 : const size_t nPixelCount = static_cast<size_t>(nOutXSize) * nOutYSize;
3403 : try
3404 : {
3405 : // Cannot overflow since pData should at least have that number of
3406 : // elements
3407 686 : if (nPixelCount >
3408 686 : static_cast<size_t>(std::numeric_limits<ptrdiff_t>::max()) /
3409 686 : static_cast<size_t>(nWordSize))
3410 : {
3411 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
3412 : "Too large temporary buffer");
3413 0 : return CE_Failure;
3414 : }
3415 686 : oWorkingState.m_abyWrkBuffer.resize(nWordSize * nPixelCount);
3416 : }
3417 0 : catch (const std::bad_alloc &e)
3418 : {
3419 0 : CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
3420 0 : return CE_Failure;
3421 : }
3422 : pafData = reinterpret_cast<const WorkingDT *>(
3423 686 : oWorkingState.m_abyWrkBuffer.data());
3424 :
3425 686 : const GDALRIOResampleAlg eResampleAlgBack = psExtraArg->eResampleAlg;
3426 686 : if (!m_osResampling.empty())
3427 : {
3428 28 : psExtraArg->eResampleAlg =
3429 28 : GDALRasterIOGetResampleAlg(m_osResampling);
3430 : }
3431 :
3432 686 : const CPLErr eErr = poSourceBand->RasterIO(
3433 : GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
3434 686 : oWorkingState.m_abyWrkBuffer.data(), nOutXSize, nOutYSize,
3435 : eWrkDataType, nWordSize,
3436 686 : nWordSize * static_cast<GSpacing>(nOutXSize), psExtraArg);
3437 686 : if (!m_osResampling.empty())
3438 28 : psExtraArg->eResampleAlg = eResampleAlgBack;
3439 :
3440 686 : if (eErr != CE_None)
3441 : {
3442 0 : return eErr;
3443 : }
3444 :
3445 : // Allocate and read mask band if needed
3446 2023 : if (!bNoDataSet &&
3447 789 : (m_nProcessingFlags & PROCESSING_FLAG_USE_MASK_BAND) != 0 &&
3448 103 : (poSourceBand->GetMaskFlags() != GMF_ALL_VALID ||
3449 51 : poSourceBand->GetColorInterpretation() == GCI_AlphaBand ||
3450 31 : GetMaskBandMainBand() != nullptr))
3451 : {
3452 : try
3453 : {
3454 103 : oWorkingState.m_abyWrkBufferMask.resize(nPixelCount);
3455 : }
3456 0 : catch (const std::exception &)
3457 : {
3458 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
3459 : "Out of memory when allocating mask buffer");
3460 0 : return CE_Failure;
3461 : }
3462 : pabyMask = reinterpret_cast<const GByte *>(
3463 103 : oWorkingState.m_abyWrkBufferMask.data());
3464 51 : auto poMaskBand =
3465 103 : (poSourceBand->GetColorInterpretation() == GCI_AlphaBand ||
3466 83 : GetMaskBandMainBand() != nullptr)
3467 : ? poSourceBand
3468 52 : : poSourceBand->GetMaskBand();
3469 103 : if (poMaskBand->RasterIO(
3470 : GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
3471 103 : oWorkingState.m_abyWrkBufferMask.data(), nOutXSize,
3472 : nOutYSize, GDT_Byte, 1, static_cast<GSpacing>(nOutXSize),
3473 103 : psExtraArg) != CE_None)
3474 : {
3475 0 : return CE_Failure;
3476 : }
3477 : }
3478 :
3479 686 : if (m_nColorTableComponent != 0)
3480 : {
3481 64 : poColorTable = poSourceBand->GetColorTable();
3482 64 : if (poColorTable == nullptr)
3483 : {
3484 0 : CPLError(CE_Failure, CPLE_AppDefined,
3485 : "Source band has no color table.");
3486 0 : return CE_Failure;
3487 : }
3488 : }
3489 : }
3490 :
3491 : /* -------------------------------------------------------------------- */
3492 : /* Selectively copy into output buffer with nodata masking, */
3493 : /* and/or scaling. */
3494 : /* -------------------------------------------------------------------- */
3495 :
3496 4727 : const bool bTwoStepDataTypeConversion =
3497 4727 : CPL_TO_BOOL(
3498 9397 : GDALDataTypeIsConversionLossy(eWrkDataType, eVRTBandDataType)) &&
3499 4670 : !CPL_TO_BOOL(GDALDataTypeIsConversionLossy(eVRTBandDataType, eBufType));
3500 :
3501 4727 : size_t idxBuffer = 0;
3502 75184 : for (int iY = 0; iY < nOutYSize; iY++)
3503 : {
3504 70457 : GByte *pDstLocation = static_cast<GByte *>(pData) +
3505 70457 : static_cast<GPtrDiff_t>(nLineSpace) * iY;
3506 :
3507 11977242 : for (int iX = 0; iX < nOutXSize;
3508 11906848 : iX++, pDstLocation += nPixelSpace, idxBuffer++)
3509 : {
3510 : WorkingDT afResult[2];
3511 11906848 : if (pafData && !bIsComplex)
3512 : {
3513 6297146 : WorkingDT fResult = pafData[idxBuffer];
3514 6297146 : if (bNoDataSetIsNan && std::isnan(fResult))
3515 1188632 : continue;
3516 6298386 : if (bNoDataSetAndNotNan &&
3517 1266 : ARE_REAL_EQUAL(fResult, fWorkingDataTypeNoData))
3518 260 : continue;
3519 6296854 : if (pabyMask && pabyMask[idxBuffer] == 0)
3520 1188350 : continue;
3521 :
3522 5108514 : if (poColorTable)
3523 : {
3524 : const GDALColorEntry *poEntry =
3525 1925900 : poColorTable->GetColorEntry(static_cast<int>(fResult));
3526 1925900 : if (poEntry)
3527 : {
3528 1925900 : if (m_nColorTableComponent == 1)
3529 683585 : fResult = poEntry->c1;
3530 1242320 : else if (m_nColorTableComponent == 2)
3531 529008 : fResult = poEntry->c2;
3532 713309 : else if (m_nColorTableComponent == 3)
3533 529008 : fResult = poEntry->c3;
3534 184301 : else if (m_nColorTableComponent == 4)
3535 184301 : fResult = poEntry->c4;
3536 : }
3537 : else
3538 : {
3539 0 : CPLErrorOnce(CE_Failure, CPLE_AppDefined,
3540 : "No entry %d.", static_cast<int>(fResult));
3541 0 : continue;
3542 : }
3543 : }
3544 :
3545 5108514 : if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_LINEAR) != 0)
3546 : {
3547 1650010 : fResult = static_cast<WorkingDT>(fResult * m_dfScaleRatio +
3548 1650010 : m_dfScaleOff);
3549 : }
3550 3458504 : else if ((m_nProcessingFlags &
3551 : PROCESSING_FLAG_SCALING_EXPONENTIAL) != 0)
3552 : {
3553 3543 : if (!m_bSrcMinMaxDefined)
3554 : {
3555 1 : int bSuccessMin = FALSE;
3556 1 : int bSuccessMax = FALSE;
3557 2 : double adfMinMax[2] = {
3558 1 : poSourceBand->GetMinimum(&bSuccessMin),
3559 1 : poSourceBand->GetMaximum(&bSuccessMax)};
3560 2 : if ((bSuccessMin && bSuccessMax) ||
3561 1 : poSourceBand->ComputeRasterMinMax(
3562 : TRUE, adfMinMax) == CE_None)
3563 : {
3564 1 : m_dfSrcMin = adfMinMax[0];
3565 1 : m_dfSrcMax = adfMinMax[1];
3566 1 : m_bSrcMinMaxDefined = true;
3567 : }
3568 : else
3569 : {
3570 0 : CPLError(CE_Failure, CPLE_AppDefined,
3571 : "Cannot determine source min/max value");
3572 0 : return CE_Failure;
3573 : }
3574 : }
3575 :
3576 3543 : double dfPowVal =
3577 3543 : (fResult - m_dfSrcMin) / (m_dfSrcMax - m_dfSrcMin);
3578 3543 : if (dfPowVal < 0.0)
3579 0 : dfPowVal = 0.0;
3580 3543 : else if (dfPowVal > 1.0)
3581 699 : dfPowVal = 1.0;
3582 3543 : fResult =
3583 3543 : static_cast<WorkingDT>((m_dfDstMax - m_dfDstMin) *
3584 3543 : pow(dfPowVal, m_dfExponent) +
3585 3543 : m_dfDstMin);
3586 : }
3587 :
3588 5108514 : if (!m_adfLUTInputs.empty())
3589 583336 : fResult = static_cast<WorkingDT>(LookupValue(fResult));
3590 :
3591 5108514 : if (m_nMaxValue != 0 && fResult > m_nMaxValue)
3592 800 : fResult = static_cast<WorkingDT>(m_nMaxValue);
3593 :
3594 5108514 : afResult[0] = fResult;
3595 5108514 : afResult[1] = 0;
3596 : }
3597 5609652 : else if (pafData && bIsComplex)
3598 : {
3599 1015 : afResult[0] = pafData[2 * idxBuffer];
3600 1015 : afResult[1] = pafData[2 * idxBuffer + 1];
3601 :
3602 : // Do not use color table.
3603 1015 : if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_LINEAR) != 0)
3604 : {
3605 1 : afResult[0] = static_cast<WorkingDT>(
3606 1 : afResult[0] * m_dfScaleRatio + m_dfScaleOff);
3607 1 : afResult[1] = static_cast<WorkingDT>(
3608 1 : afResult[1] * m_dfScaleRatio + m_dfScaleOff);
3609 : }
3610 :
3611 : /* Do not use LUT */
3612 : }
3613 : else
3614 : {
3615 5608642 : afResult[0] = static_cast<WorkingDT>(m_dfScaleOff);
3616 5608642 : afResult[1] = 0;
3617 :
3618 5608642 : if (!m_adfLUTInputs.empty())
3619 0 : afResult[0] =
3620 0 : static_cast<WorkingDT>(LookupValue(afResult[0]));
3621 :
3622 5608642 : if (m_nMaxValue != 0 && afResult[0] > m_nMaxValue)
3623 0 : afResult[0] = static_cast<WorkingDT>(m_nMaxValue);
3624 : }
3625 :
3626 10718146 : if (eBufType == GDT_Byte && eVRTBandDataType == GDT_Byte)
3627 : {
3628 4299880 : *pDstLocation = static_cast<GByte>(std::min(
3629 8599760 : 255.0f,
3630 4299880 : std::max(0.0f, static_cast<float>(afResult[0]) + 0.5f)));
3631 : }
3632 6418286 : else if (!bTwoStepDataTypeConversion)
3633 : {
3634 4021 : CopyWord<WorkingDT>(afResult, eWrkDataType, pDstLocation,
3635 : eBufType);
3636 : }
3637 6414270 : else if (eVRTBandDataType == GDT_Float32 && eBufType == GDT_Float64)
3638 : {
3639 : // Particular case of the below 2-step conversion.
3640 : // Helps a bit for some geolocation based warping with Sentinel3
3641 : // data where the longitude/latitude arrays are Int32 bands,
3642 : // rescaled in VRT as Float32 and requested as Float64
3643 : float fVal;
3644 20 : GDALCopyWord(afResult[0], fVal);
3645 20 : *reinterpret_cast<double *>(pDstLocation) = fVal;
3646 : }
3647 : else
3648 : {
3649 : GByte abyTemp[2 * sizeof(double)];
3650 : // Convert first to the VRTRasterBand data type
3651 : // to get its clamping, before outputting to buffer data type
3652 6414250 : CopyWord<WorkingDT>(afResult, eWrkDataType, abyTemp,
3653 : eVRTBandDataType);
3654 6414250 : GDALCopyWords(abyTemp, eVRTBandDataType, 0, pDstLocation,
3655 : eBufType, 0, 1);
3656 : }
3657 : }
3658 : }
3659 :
3660 4727 : return CE_None;
3661 : }
3662 :
3663 : // Explicitly instantiate template method, as it is used in another file.
3664 : template CPLErr VRTComplexSource::RasterIOInternal<float>(
3665 : GDALRasterBand *poSourceBand, GDALDataType eVRTBandDataType, int nReqXOff,
3666 : int nReqYOff, int nReqXSize, int nReqYSize, void *pData, int nOutXSize,
3667 : int nOutYSize, GDALDataType eBufType, GSpacing nPixelSpace,
3668 : GSpacing nLineSpace, GDALRasterIOExtraArg *psExtraArg,
3669 : GDALDataType eWrkDataType, WorkingState &oWorkingState);
3670 :
3671 : /************************************************************************/
3672 : /* AreValuesUnchanged() */
3673 : /************************************************************************/
3674 :
3675 9 : bool VRTComplexSource::AreValuesUnchanged() const
3676 : {
3677 10 : return m_dfScaleOff == 0.0 && m_dfScaleRatio == 1.0 &&
3678 24 : m_adfLUTInputs.empty() && m_nColorTableComponent == 0 &&
3679 14 : (m_nProcessingFlags & PROCESSING_FLAG_SCALING_EXPONENTIAL) == 0;
3680 : }
3681 :
3682 : /************************************************************************/
3683 : /* GetMinimum() */
3684 : /************************************************************************/
3685 :
3686 2 : double VRTComplexSource::GetMinimum(int nXSize, int nYSize, int *pbSuccess)
3687 : {
3688 2 : if (AreValuesUnchanged())
3689 : {
3690 1 : return VRTSimpleSource::GetMinimum(nXSize, nYSize, pbSuccess);
3691 : }
3692 :
3693 1 : *pbSuccess = FALSE;
3694 1 : return 0;
3695 : }
3696 :
3697 : /************************************************************************/
3698 : /* GetMaximum() */
3699 : /************************************************************************/
3700 :
3701 2 : double VRTComplexSource::GetMaximum(int nXSize, int nYSize, int *pbSuccess)
3702 : {
3703 2 : if (AreValuesUnchanged())
3704 : {
3705 1 : return VRTSimpleSource::GetMaximum(nXSize, nYSize, pbSuccess);
3706 : }
3707 :
3708 1 : *pbSuccess = FALSE;
3709 1 : return 0;
3710 : }
3711 :
3712 : /************************************************************************/
3713 : /* GetHistogram() */
3714 : /************************************************************************/
3715 :
3716 0 : CPLErr VRTComplexSource::GetHistogram(int nXSize, int nYSize, double dfMin,
3717 : double dfMax, int nBuckets,
3718 : GUIntBig *panHistogram,
3719 : int bIncludeOutOfRange, int bApproxOK,
3720 : GDALProgressFunc pfnProgress,
3721 : void *pProgressData)
3722 : {
3723 0 : if (AreValuesUnchanged())
3724 : {
3725 0 : return VRTSimpleSource::GetHistogram(
3726 : nXSize, nYSize, dfMin, dfMax, nBuckets, panHistogram,
3727 0 : bIncludeOutOfRange, bApproxOK, pfnProgress, pProgressData);
3728 : }
3729 :
3730 0 : return CE_Failure;
3731 : }
3732 :
3733 : /************************************************************************/
3734 : /* ==================================================================== */
3735 : /* VRTFuncSource */
3736 : /* ==================================================================== */
3737 : /************************************************************************/
3738 :
3739 : /************************************************************************/
3740 : /* VRTFuncSource() */
3741 : /************************************************************************/
3742 :
3743 12 : VRTFuncSource::VRTFuncSource()
3744 : : pfnReadFunc(nullptr), pCBData(nullptr), eType(GDT_Byte),
3745 12 : fNoDataValue(static_cast<float>(VRT_NODATA_UNSET))
3746 : {
3747 12 : }
3748 :
3749 : /************************************************************************/
3750 : /* ~VRTFuncSource() */
3751 : /************************************************************************/
3752 :
3753 24 : VRTFuncSource::~VRTFuncSource()
3754 : {
3755 24 : }
3756 :
3757 : /************************************************************************/
3758 : /* GetType() */
3759 : /************************************************************************/
3760 :
3761 0 : const char *VRTFuncSource::GetType() const
3762 : {
3763 : static const char *TYPE = "FuncSource";
3764 0 : return TYPE;
3765 : }
3766 :
3767 : /************************************************************************/
3768 : /* SerializeToXML() */
3769 : /************************************************************************/
3770 :
3771 0 : CPLXMLNode *VRTFuncSource::SerializeToXML(CPL_UNUSED const char *pszVRTPath)
3772 : {
3773 0 : return nullptr;
3774 : }
3775 :
3776 : /************************************************************************/
3777 : /* RasterIO() */
3778 : /************************************************************************/
3779 :
3780 10 : CPLErr VRTFuncSource::RasterIO(GDALDataType /*eVRTBandDataType*/, int nXOff,
3781 : int nYOff, int nXSize, int nYSize, void *pData,
3782 : int nBufXSize, int nBufYSize,
3783 : GDALDataType eBufType, GSpacing nPixelSpace,
3784 : GSpacing nLineSpace,
3785 : GDALRasterIOExtraArg * /* psExtraArg */,
3786 : WorkingState & /* oWorkingState */)
3787 : {
3788 10 : if (nPixelSpace * 8 == GDALGetDataTypeSize(eBufType) &&
3789 10 : nLineSpace == nPixelSpace * nXSize && nBufXSize == nXSize &&
3790 20 : nBufYSize == nYSize && eBufType == eType)
3791 : {
3792 10 : return pfnReadFunc(pCBData, nXOff, nYOff, nXSize, nYSize, pData);
3793 : }
3794 : else
3795 : {
3796 0 : CPLError(CE_Failure, CPLE_AppDefined,
3797 : "VRTFuncSource::RasterIO() - Irregular request.");
3798 0 : CPLDebug("VRT", "Irregular request: %d,%d %d,%d, %d,%d %d,%d %d,%d",
3799 : static_cast<int>(nPixelSpace) * 8,
3800 : GDALGetDataTypeSize(eBufType), static_cast<int>(nLineSpace),
3801 : static_cast<int>(nPixelSpace) * nXSize, nBufXSize, nXSize,
3802 : nBufYSize, nYSize, static_cast<int>(eBufType),
3803 0 : static_cast<int>(eType));
3804 :
3805 0 : return CE_Failure;
3806 : }
3807 : }
3808 :
3809 : /************************************************************************/
3810 : /* GetMinimum() */
3811 : /************************************************************************/
3812 :
3813 0 : double VRTFuncSource::GetMinimum(int /* nXSize */, int /* nYSize */,
3814 : int *pbSuccess)
3815 : {
3816 0 : *pbSuccess = FALSE;
3817 0 : return 0;
3818 : }
3819 :
3820 : /************************************************************************/
3821 : /* GetMaximum() */
3822 : /************************************************************************/
3823 :
3824 0 : double VRTFuncSource::GetMaximum(int /* nXSize */, int /* nYSize */,
3825 : int *pbSuccess)
3826 : {
3827 0 : *pbSuccess = FALSE;
3828 0 : return 0;
3829 : }
3830 :
3831 : /************************************************************************/
3832 : /* GetHistogram() */
3833 : /************************************************************************/
3834 :
3835 0 : CPLErr VRTFuncSource::GetHistogram(
3836 : int /* nXSize */, int /* nYSize */, double /* dfMin */, double /* dfMax */,
3837 : int /* nBuckets */, GUIntBig * /* panHistogram */,
3838 : int /* bIncludeOutOfRange */, int /* bApproxOK */,
3839 : GDALProgressFunc /* pfnProgress */, void * /* pProgressData */)
3840 : {
3841 0 : return CE_Failure;
3842 : }
3843 :
3844 : /************************************************************************/
3845 : /* VRTParseCoreSources() */
3846 : /************************************************************************/
3847 :
3848 3083 : VRTSource *VRTParseCoreSources(const CPLXMLNode *psChild,
3849 : const char *pszVRTPath,
3850 : VRTMapSharedResources &oMapSharedSources)
3851 :
3852 : {
3853 3083 : VRTSource *poSource = nullptr;
3854 :
3855 6159 : if (EQUAL(psChild->pszValue, VRTAveragedSource::GetTypeStatic()) ||
3856 3076 : (EQUAL(psChild->pszValue, VRTSimpleSource::GetTypeStatic()) &&
3857 2865 : STARTS_WITH_CI(CPLGetXMLValue(psChild, "Resampling", "Nearest"),
3858 : "Aver")))
3859 : {
3860 16 : poSource = new VRTAveragedSource();
3861 : }
3862 3067 : else if (EQUAL(psChild->pszValue, VRTSimpleSource::GetTypeStatic()))
3863 : {
3864 2856 : poSource = new VRTSimpleSource();
3865 : }
3866 211 : else if (EQUAL(psChild->pszValue, VRTComplexSource::GetTypeStatic()))
3867 : {
3868 203 : poSource = new VRTComplexSource();
3869 : }
3870 8 : else if (EQUAL(psChild->pszValue, VRTNoDataFromMaskSource::GetTypeStatic()))
3871 : {
3872 8 : poSource = new VRTNoDataFromMaskSource();
3873 : }
3874 : else
3875 : {
3876 0 : CPLError(CE_Failure, CPLE_AppDefined,
3877 : "VRTParseCoreSources() - Unknown source : %s",
3878 0 : psChild->pszValue);
3879 0 : return nullptr;
3880 : }
3881 :
3882 3083 : if (poSource->XMLInit(psChild, pszVRTPath, oMapSharedSources) == CE_None)
3883 3080 : return poSource;
3884 :
3885 3 : delete poSource;
3886 3 : return nullptr;
3887 : }
3888 :
3889 : /*! @endcond */
|
