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