Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: Virtual GDAL Datasets
4 : * Purpose: Implementation of VRTProcessedDataset.
5 : * Author: Even Rouault <even.rouault at spatialys.com>
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2024, Even Rouault <even.rouault at spatialys.com>
9 : *
10 : * SPDX-License-Identifier: MIT
11 : ****************************************************************************/
12 :
13 : #include "cpl_minixml.h"
14 : #include "cpl_string.h"
15 : #include "gdal_utils.h"
16 : #include "vrtdataset.h"
17 :
18 : #include <algorithm>
19 : #include <limits>
20 : #include <map>
21 : #include <vector>
22 :
23 : /************************************************************************/
24 : /* VRTProcessedDatasetFunc */
25 : /************************************************************************/
26 :
27 : //! Structure holding information for a VRTProcessedDataset function.
28 : struct VRTProcessedDatasetFunc
29 : {
30 : //! Processing function name
31 : std::string osFuncName{};
32 :
33 : //! User data to provide to pfnInit, pfnFree, pfnProcess callbacks.
34 : void *pUserData = nullptr;
35 :
36 : //! Whether XML metadata has been specified
37 : bool bMetadataSpecified = false;
38 :
39 : //! Map of (constant argument name, constant value)
40 : std::map<std::string, std::string> oMapConstantArguments{};
41 :
42 : //! Set of builtin argument names (e.g "offset", "scale", "nodata")
43 : std::set<std::string> oSetBuiltinArguments{};
44 :
45 : //! Arguments defined in the VRT
46 : struct OtherArgument
47 : {
48 : std::string osType{};
49 : bool bRequired = false;
50 : };
51 :
52 : std::map<std::string, OtherArgument> oOtherArguments{};
53 :
54 : //! Requested input data type.
55 : GDALDataType eRequestedInputDT = GDT_Unknown;
56 :
57 : //! List of supported input datatypes. Empty if no restriction.
58 : std::vector<GDALDataType> aeSupportedInputDT{};
59 :
60 : //! List of supported input band counts. Empty if no restriction.
61 : std::vector<int> anSupportedInputBandCount{};
62 :
63 : //! Optional initialization function
64 : GDALVRTProcessedDatasetFuncInit pfnInit = nullptr;
65 :
66 : //! Optional free function
67 : GDALVRTProcessedDatasetFuncFree pfnFree = nullptr;
68 :
69 : //! Required processing function
70 : GDALVRTProcessedDatasetFuncProcess pfnProcess = nullptr;
71 : };
72 :
73 : /************************************************************************/
74 : /* GetGlobalMapProcessedDatasetFunc() */
75 : /************************************************************************/
76 :
77 : /** Return the registry of VRTProcessedDatasetFunc functions */
78 : static std::map<std::string, VRTProcessedDatasetFunc> &
79 7082 : GetGlobalMapProcessedDatasetFunc()
80 : {
81 7082 : static std::map<std::string, VRTProcessedDatasetFunc> goMap;
82 7082 : return goMap;
83 : }
84 :
85 : /************************************************************************/
86 : /* Step::~Step() */
87 : /************************************************************************/
88 :
89 : /*! @cond Doxygen_Suppress */
90 :
91 : /** Step destructor */
92 128 : VRTProcessedDataset::Step::~Step()
93 : {
94 128 : deinit();
95 128 : }
96 :
97 : /************************************************************************/
98 : /* Step::deinit() */
99 : /************************************************************************/
100 :
101 : /** Free pWorkingData */
102 128 : void VRTProcessedDataset::Step::deinit()
103 : {
104 128 : if (pWorkingData)
105 : {
106 48 : const auto &oMapFunctions = GetGlobalMapProcessedDatasetFunc();
107 48 : const auto oIterFunc = oMapFunctions.find(osAlgorithm);
108 48 : if (oIterFunc != oMapFunctions.end())
109 : {
110 48 : if (oIterFunc->second.pfnFree)
111 : {
112 96 : oIterFunc->second.pfnFree(osAlgorithm.c_str(),
113 48 : oIterFunc->second.pUserData,
114 : pWorkingData);
115 : }
116 : }
117 : else
118 : {
119 0 : CPLAssert(false);
120 : }
121 48 : pWorkingData = nullptr;
122 : }
123 128 : }
124 :
125 : /************************************************************************/
126 : /* Step::Step(Step&& other) */
127 : /************************************************************************/
128 :
129 : /** Move constructor */
130 51 : VRTProcessedDataset::Step::Step(Step &&other)
131 51 : : osAlgorithm(std::move(other.osAlgorithm)),
132 102 : aosArguments(std::move(other.aosArguments)), eInDT(other.eInDT),
133 51 : eOutDT(other.eOutDT), nInBands(other.nInBands),
134 51 : nOutBands(other.nOutBands), adfInNoData(other.adfInNoData),
135 51 : adfOutNoData(other.adfOutNoData), pWorkingData(other.pWorkingData)
136 : {
137 51 : other.pWorkingData = nullptr;
138 51 : }
139 :
140 : /************************************************************************/
141 : /* Step operator=(Step&& other) */
142 : /************************************************************************/
143 :
144 : /** Move assignment operator */
145 0 : VRTProcessedDataset::Step &VRTProcessedDataset::Step::operator=(Step &&other)
146 : {
147 0 : if (&other != this)
148 : {
149 0 : deinit();
150 0 : osAlgorithm = std::move(other.osAlgorithm);
151 0 : aosArguments = std::move(other.aosArguments);
152 0 : eInDT = other.eInDT;
153 0 : eOutDT = other.eOutDT;
154 0 : nInBands = other.nInBands;
155 0 : nOutBands = other.nOutBands;
156 0 : adfInNoData = std::move(other.adfInNoData);
157 0 : adfOutNoData = std::move(other.adfOutNoData);
158 0 : std::swap(pWorkingData, other.pWorkingData);
159 : }
160 0 : return *this;
161 : }
162 :
163 : /************************************************************************/
164 : /* VRTProcessedDataset() */
165 : /************************************************************************/
166 :
167 : /** Constructor */
168 87 : VRTProcessedDataset::VRTProcessedDataset(int nXSize, int nYSize)
169 87 : : VRTDataset(nXSize, nYSize)
170 : {
171 87 : }
172 :
173 : /************************************************************************/
174 : /* ~VRTProcessedDataset() */
175 : /************************************************************************/
176 :
177 174 : VRTProcessedDataset::~VRTProcessedDataset()
178 :
179 : {
180 87 : VRTProcessedDataset::FlushCache(true);
181 87 : VRTProcessedDataset::CloseDependentDatasets();
182 174 : }
183 :
184 : /************************************************************************/
185 : /* XMLInit() */
186 : /************************************************************************/
187 :
188 : /** Instantiate object from XML tree */
189 84 : CPLErr VRTProcessedDataset::XMLInit(const CPLXMLNode *psTree,
190 : const char *pszVRTPathIn)
191 :
192 : {
193 84 : if (Init(psTree, pszVRTPathIn, nullptr, nullptr, -1) != CE_None)
194 41 : return CE_Failure;
195 :
196 43 : const auto poSrcFirstBand = m_poSrcDS->GetRasterBand(1);
197 43 : const int nOvrCount = poSrcFirstBand->GetOverviewCount();
198 46 : for (int i = 0; i < nOvrCount; ++i)
199 : {
200 3 : auto poOvrDS = std::make_unique<VRTProcessedDataset>(0, 0);
201 3 : if (poOvrDS->Init(psTree, pszVRTPathIn, this, m_poSrcDS.get(), i) !=
202 : CE_None)
203 0 : break;
204 3 : m_apoOverviewDatasets.emplace_back(std::move(poOvrDS));
205 : }
206 :
207 43 : return CE_None;
208 : }
209 :
210 72 : static bool HasScaleOffset(GDALDataset &oSrcDS)
211 : {
212 281 : for (int i = 1; i <= oSrcDS.GetRasterCount(); i++)
213 : {
214 : int pbSuccess;
215 211 : GDALRasterBand &oBand = *oSrcDS.GetRasterBand(i);
216 211 : double scale = oBand.GetScale(&pbSuccess);
217 211 : if (pbSuccess && scale != 1)
218 : {
219 2 : return true;
220 : }
221 209 : double offset = oBand.GetOffset(&pbSuccess);
222 209 : if (pbSuccess && offset != 0)
223 : {
224 0 : return true;
225 : }
226 : }
227 :
228 70 : return false;
229 : }
230 :
231 : /** Instantiate object from XML tree */
232 87 : CPLErr VRTProcessedDataset::Init(const CPLXMLNode *psTree,
233 : const char *pszVRTPathIn,
234 : const VRTProcessedDataset *poParentDS,
235 : GDALDataset *poParentSrcDS, int iOvrLevel)
236 :
237 : {
238 87 : const CPLXMLNode *psInput = CPLGetXMLNode(psTree, "Input");
239 87 : if (!psInput)
240 : {
241 1 : CPLError(CE_Failure, CPLE_AppDefined, "Input element missing");
242 1 : return CE_Failure;
243 : }
244 :
245 86 : if (pszVRTPathIn)
246 12 : m_osVRTPath = pszVRTPathIn;
247 :
248 86 : if (poParentSrcDS)
249 : {
250 3 : m_poSrcDS.reset(
251 : GDALCreateOverviewDataset(poParentSrcDS, iOvrLevel, true));
252 : }
253 83 : else if (const CPLXMLNode *psSourceFileNameNode =
254 83 : CPLGetXMLNode(psInput, "SourceFilename"))
255 : {
256 81 : const bool bRelativeToVRT = CPL_TO_BOOL(
257 : atoi(CPLGetXMLValue(psSourceFileNameNode, "relativetoVRT", "0")));
258 : const std::string osFilename = VRTDataset::BuildSourceFilename(
259 : CPLGetXMLValue(psInput, "SourceFilename", ""), pszVRTPathIn,
260 162 : bRelativeToVRT);
261 81 : m_poSrcDS.reset(GDALDataset::Open(
262 : osFilename.c_str(), GDAL_OF_RASTER | GDAL_OF_VERBOSE_ERROR, nullptr,
263 : nullptr, nullptr));
264 : }
265 2 : else if (const CPLXMLNode *psVRTDataset =
266 2 : CPLGetXMLNode(psInput, "VRTDataset"))
267 : {
268 1 : CPLXMLNode sVRTDatasetTmp = *psVRTDataset;
269 1 : sVRTDatasetTmp.psNext = nullptr;
270 1 : char *pszXML = CPLSerializeXMLTree(&sVRTDatasetTmp);
271 1 : m_poSrcDS.reset(VRTDataset::OpenXML(pszXML, pszVRTPathIn, GA_ReadOnly));
272 1 : CPLFree(pszXML);
273 : }
274 : else
275 : {
276 1 : CPLError(
277 : CE_Failure, CPLE_AppDefined,
278 : "Input element should have a SourceFilename or VRTDataset element");
279 1 : return CE_Failure;
280 : }
281 :
282 85 : if (!m_poSrcDS)
283 2 : return CE_Failure;
284 :
285 83 : const char *pszUnscale = CPLGetXMLValue(psInput, "unscale", "AUTO");
286 83 : bool bUnscale = false;
287 83 : if (EQUAL(pszUnscale, "AUTO"))
288 : {
289 72 : if (HasScaleOffset(*m_poSrcDS))
290 : {
291 2 : bUnscale = true;
292 : }
293 : }
294 11 : else if (EQUAL(pszUnscale, "YES") || EQUAL(pszUnscale, "ON") ||
295 11 : EQUAL(pszUnscale, "TRUE") || EQUAL(pszUnscale, "1"))
296 : {
297 6 : bUnscale = true;
298 : }
299 5 : else if (!(EQUAL(pszUnscale, "NO") || EQUAL(pszUnscale, "OFF") ||
300 5 : EQUAL(pszUnscale, "FALSE") || EQUAL(pszUnscale, "0")))
301 : {
302 1 : CPLError(CE_Failure, CPLE_AppDefined, "Invalid value of 'unscale'");
303 1 : return CE_Failure;
304 : }
305 :
306 82 : if (bUnscale)
307 : {
308 8 : CPLStringList oArgs;
309 8 : oArgs.AddString("-unscale");
310 8 : oArgs.AddString("-ot");
311 8 : oArgs.AddString("Float64");
312 8 : oArgs.AddString("-of");
313 8 : oArgs.AddString("VRT");
314 8 : oArgs.AddString("-a_nodata");
315 8 : oArgs.AddString("nan");
316 8 : auto *poArgs = GDALTranslateOptionsNew(oArgs.List(), nullptr);
317 : int pbUsageError;
318 8 : CPLAssert(poArgs);
319 8 : m_poVRTSrcDS.reset(m_poSrcDS.release());
320 : // https://trac.cppcheck.net/ticket/11325
321 : // cppcheck-suppress accessMoved
322 8 : m_poSrcDS.reset(GDALDataset::FromHandle(
323 8 : GDALTranslate("", m_poVRTSrcDS.get(), poArgs, &pbUsageError)));
324 8 : GDALTranslateOptionsFree(poArgs);
325 :
326 8 : if (pbUsageError || !m_poSrcDS)
327 : {
328 0 : return CE_Failure;
329 : }
330 : }
331 :
332 82 : if (nRasterXSize == 0 && nRasterYSize == 0)
333 : {
334 80 : nRasterXSize = m_poSrcDS->GetRasterXSize();
335 80 : nRasterYSize = m_poSrcDS->GetRasterYSize();
336 : }
337 2 : else if (nRasterXSize != m_poSrcDS->GetRasterXSize() ||
338 0 : nRasterYSize != m_poSrcDS->GetRasterYSize())
339 : {
340 2 : CPLError(CE_Failure, CPLE_AppDefined,
341 : "Inconsistent declared VRT dimensions with input dataset");
342 2 : return CE_Failure;
343 : }
344 :
345 80 : if (m_poSrcDS->GetRasterCount() == 0)
346 0 : return CE_Failure;
347 :
348 : // Inherit SRS from source if not explicitly defined in VRT
349 80 : if (!CPLGetXMLNode(psTree, "SRS"))
350 : {
351 80 : const OGRSpatialReference *poSRS = m_poSrcDS->GetSpatialRef();
352 80 : if (poSRS)
353 : {
354 10 : m_poSRS.reset(poSRS->Clone());
355 : }
356 : }
357 :
358 : // Inherit GeoTransform from source if not explicitly defined in VRT
359 80 : if (iOvrLevel < 0 && !CPLGetXMLNode(psTree, "GeoTransform"))
360 : {
361 77 : if (m_poSrcDS->GetGeoTransform(m_adfGeoTransform) == CE_None)
362 40 : m_bGeoTransformSet = true;
363 : }
364 :
365 : /* -------------------------------------------------------------------- */
366 : /* Initialize blocksize before calling sub-init so that the */
367 : /* band initializers can get it from the dataset object when */
368 : /* they are created. */
369 : /* -------------------------------------------------------------------- */
370 :
371 80 : const auto poSrcFirstBand = m_poSrcDS->GetRasterBand(1);
372 80 : poSrcFirstBand->GetBlockSize(&m_nBlockXSize, &m_nBlockYSize);
373 80 : bool bUserBlockSize = false;
374 80 : if (const char *pszBlockXSize =
375 80 : CPLGetXMLValue(psTree, "BlockXSize", nullptr))
376 : {
377 0 : bUserBlockSize = true;
378 0 : m_nBlockXSize = atoi(pszBlockXSize);
379 0 : if (m_nBlockXSize <= 1)
380 : {
381 0 : CPLError(CE_Failure, CPLE_AppDefined, "Invalid BlockXSize");
382 0 : return CE_Failure;
383 : }
384 : }
385 80 : if (const char *pszBlockYSize =
386 80 : CPLGetXMLValue(psTree, "BlockYSize", nullptr))
387 : {
388 0 : bUserBlockSize = true;
389 0 : m_nBlockYSize = atoi(pszBlockYSize);
390 0 : if (m_nBlockYSize <= 1)
391 : {
392 0 : CPLError(CE_Failure, CPLE_AppDefined, "Invalid BlockYSize");
393 0 : return CE_Failure;
394 : }
395 : }
396 :
397 : // Initialize all the general VRT stuff.
398 80 : if (VRTDataset::XMLInit(psTree, pszVRTPathIn) != CE_None)
399 : {
400 0 : return CE_Failure;
401 : }
402 :
403 : // Use geotransform from parent for overviews
404 80 : if (iOvrLevel >= 0 && poParentDS->m_bGeoTransformSet)
405 : {
406 1 : m_bGeoTransformSet = true;
407 1 : m_adfGeoTransform[0] = poParentDS->m_adfGeoTransform[0];
408 1 : m_adfGeoTransform[1] = poParentDS->m_adfGeoTransform[1];
409 1 : m_adfGeoTransform[2] = poParentDS->m_adfGeoTransform[2];
410 1 : m_adfGeoTransform[3] = poParentDS->m_adfGeoTransform[3];
411 1 : m_adfGeoTransform[4] = poParentDS->m_adfGeoTransform[4];
412 1 : m_adfGeoTransform[5] = poParentDS->m_adfGeoTransform[5];
413 :
414 1 : m_adfGeoTransform[1] *=
415 1 : static_cast<double>(poParentDS->GetRasterXSize()) / nRasterXSize;
416 1 : m_adfGeoTransform[2] *=
417 1 : static_cast<double>(poParentDS->GetRasterYSize()) / nRasterYSize;
418 1 : m_adfGeoTransform[4] *=
419 1 : static_cast<double>(poParentDS->GetRasterXSize()) / nRasterXSize;
420 1 : m_adfGeoTransform[5] *=
421 1 : static_cast<double>(poParentDS->GetRasterYSize()) / nRasterYSize;
422 : }
423 :
424 80 : const CPLXMLNode *psOutputBands = CPLGetXMLNode(psTree, "OutputBands");
425 80 : if (psOutputBands)
426 : {
427 14 : if (const char *pszCount =
428 14 : CPLGetXMLValue(psOutputBands, "count", nullptr))
429 : {
430 14 : if (EQUAL(pszCount, "FROM_LAST_STEP"))
431 : {
432 8 : m_outputBandCountProvenance = ValueProvenance::FROM_LAST_STEP;
433 : }
434 6 : else if (!EQUAL(pszCount, "FROM_SOURCE"))
435 : {
436 4 : if (CPLGetValueType(pszCount) == CPL_VALUE_INTEGER)
437 : {
438 3 : m_outputBandCountProvenance =
439 : ValueProvenance::USER_PROVIDED;
440 3 : m_outputBandCountValue = atoi(pszCount);
441 3 : if (!GDALCheckBandCount(m_outputBandCountValue,
442 : /* bIsZeroAllowed = */ false))
443 : {
444 : // Error emitted by above function
445 1 : return CE_Failure;
446 : }
447 : }
448 : else
449 : {
450 1 : CPLError(CE_Failure, CPLE_AppDefined,
451 : "Invalid value for OutputBands.count");
452 1 : return CE_Failure;
453 : }
454 : }
455 : }
456 :
457 12 : if (const char *pszType =
458 12 : CPLGetXMLValue(psOutputBands, "dataType", nullptr))
459 : {
460 12 : if (EQUAL(pszType, "FROM_LAST_STEP"))
461 : {
462 4 : m_outputBandDataTypeProvenance =
463 : ValueProvenance::FROM_LAST_STEP;
464 : }
465 8 : else if (!EQUAL(pszType, "FROM_SOURCE"))
466 : {
467 6 : m_outputBandDataTypeProvenance = ValueProvenance::USER_PROVIDED;
468 6 : m_outputBandDataTypeValue = GDALGetDataTypeByName(pszType);
469 6 : if (m_outputBandDataTypeValue == GDT_Unknown)
470 : {
471 1 : CPLError(CE_Failure, CPLE_AppDefined,
472 : "Invalid value for OutputBands.dataType");
473 1 : return CE_Failure;
474 : }
475 : }
476 : }
477 : }
478 66 : else if (CPLGetXMLNode(psTree, "VRTRasterBand"))
479 : {
480 5 : m_outputBandCountProvenance = ValueProvenance::FROM_VRTRASTERBAND;
481 5 : m_outputBandDataTypeProvenance = ValueProvenance::FROM_VRTRASTERBAND;
482 : }
483 :
484 77 : int nOutputBandCount = 0;
485 77 : switch (m_outputBandCountProvenance)
486 : {
487 1 : case ValueProvenance::USER_PROVIDED:
488 1 : nOutputBandCount = m_outputBandCountValue;
489 1 : break;
490 63 : case ValueProvenance::FROM_SOURCE:
491 63 : nOutputBandCount = m_poSrcDS->GetRasterCount();
492 63 : break;
493 5 : case ValueProvenance::FROM_VRTRASTERBAND:
494 5 : nOutputBandCount = nBands;
495 5 : break;
496 8 : case ValueProvenance::FROM_LAST_STEP:
497 8 : break;
498 : }
499 :
500 : const CPLXMLNode *psProcessingSteps =
501 77 : CPLGetXMLNode(psTree, "ProcessingSteps");
502 77 : if (!psProcessingSteps)
503 : {
504 1 : CPLError(CE_Failure, CPLE_AppDefined,
505 : "ProcessingSteps element missing");
506 1 : return CE_Failure;
507 : }
508 :
509 76 : const auto eInDT = poSrcFirstBand->GetRasterDataType();
510 219 : for (int i = 1; i < m_poSrcDS->GetRasterCount(); ++i)
511 : {
512 143 : const auto eDT = m_poSrcDS->GetRasterBand(i + 1)->GetRasterDataType();
513 143 : if (eDT != eInDT)
514 : {
515 0 : CPLError(CE_Warning, CPLE_AppDefined,
516 : "Not all bands of the input dataset have the same data "
517 : "type. The data type of the first band will be used as "
518 : "the reference one.");
519 0 : break;
520 : }
521 : }
522 :
523 76 : GDALDataType eCurrentDT = eInDT;
524 76 : int nCurrentBandCount = m_poSrcDS->GetRasterCount();
525 :
526 152 : std::vector<double> adfNoData;
527 295 : for (int i = 1; i <= nCurrentBandCount; ++i)
528 : {
529 219 : int bHasVal = FALSE;
530 : const double dfVal =
531 219 : m_poSrcDS->GetRasterBand(i)->GetNoDataValue(&bHasVal);
532 : adfNoData.emplace_back(
533 219 : bHasVal ? dfVal : std::numeric_limits<double>::quiet_NaN());
534 : }
535 :
536 76 : int nStepCount = 0;
537 153 : for (const CPLXMLNode *psStep = psProcessingSteps->psChild; psStep;
538 77 : psStep = psStep->psNext)
539 : {
540 77 : if (psStep->eType == CXT_Element &&
541 77 : strcmp(psStep->pszValue, "Step") == 0)
542 : {
543 77 : ++nStepCount;
544 : }
545 : }
546 :
547 76 : int iStep = 0;
548 124 : for (const CPLXMLNode *psStep = psProcessingSteps->psChild; psStep;
549 48 : psStep = psStep->psNext)
550 : {
551 77 : if (psStep->eType == CXT_Element &&
552 77 : strcmp(psStep->pszValue, "Step") == 0)
553 : {
554 77 : ++iStep;
555 77 : const bool bIsFinalStep = (iStep == nStepCount);
556 77 : std::vector<double> adfOutNoData;
557 77 : if (bIsFinalStep)
558 : {
559 : // Initialize adfOutNoData with nodata value of *output* bands
560 : // for final step
561 75 : if (m_outputBandCountProvenance ==
562 : ValueProvenance::FROM_VRTRASTERBAND)
563 : {
564 13 : for (int i = 1; i <= nBands; ++i)
565 : {
566 8 : int bHasVal = FALSE;
567 : const double dfVal =
568 8 : GetRasterBand(i)->GetNoDataValue(&bHasVal);
569 : adfOutNoData.emplace_back(
570 8 : bHasVal ? dfVal
571 8 : : std::numeric_limits<double>::quiet_NaN());
572 : }
573 : }
574 70 : else if (m_outputBandCountProvenance ==
575 : ValueProvenance::FROM_SOURCE)
576 : {
577 213 : for (int i = 1; i <= m_poSrcDS->GetRasterCount(); ++i)
578 : {
579 152 : int bHasVal = FALSE;
580 : const double dfVal =
581 152 : m_poSrcDS->GetRasterBand(i)->GetNoDataValue(
582 152 : &bHasVal);
583 : adfOutNoData.emplace_back(
584 182 : bHasVal ? dfVal
585 182 : : std::numeric_limits<double>::quiet_NaN());
586 : }
587 : }
588 9 : else if (m_outputBandCountProvenance ==
589 : ValueProvenance::USER_PROVIDED)
590 : {
591 1 : adfOutNoData.resize(
592 1 : m_outputBandCountValue,
593 1 : std::numeric_limits<double>::quiet_NaN());
594 : }
595 : }
596 77 : if (!ParseStep(psStep, bIsFinalStep, eCurrentDT, nCurrentBandCount,
597 : adfNoData, adfOutNoData))
598 29 : return CE_Failure;
599 48 : adfNoData = std::move(adfOutNoData);
600 : }
601 : }
602 :
603 47 : if (m_aoSteps.empty())
604 : {
605 1 : CPLError(CE_Failure, CPLE_AppDefined,
606 : "At least one step should be defined");
607 1 : return CE_Failure;
608 : }
609 :
610 46 : int nLargestInDTSizeTimesBand = 1;
611 46 : int nLargestOutDTSizeTimesBand = 1;
612 94 : for (const auto &oStep : m_aoSteps)
613 : {
614 : const int nInDTSizeTimesBand =
615 48 : GDALGetDataTypeSizeBytes(oStep.eInDT) * oStep.nInBands;
616 48 : nLargestInDTSizeTimesBand =
617 48 : std::max(nLargestInDTSizeTimesBand, nInDTSizeTimesBand);
618 : const int nOutDTSizeTimesBand =
619 48 : GDALGetDataTypeSizeBytes(oStep.eOutDT) * oStep.nOutBands;
620 48 : nLargestOutDTSizeTimesBand =
621 48 : std::max(nLargestOutDTSizeTimesBand, nOutDTSizeTimesBand);
622 : }
623 46 : m_nWorkingBytesPerPixel =
624 46 : nLargestInDTSizeTimesBand + nLargestOutDTSizeTimesBand;
625 :
626 : // Use only up to 40% of RAM to acquire source bands and generate the output
627 : // buffer.
628 46 : m_nAllowedRAMUsage = CPLGetUsablePhysicalRAM() / 10 * 4;
629 : // Only for tests now
630 46 : const char *pszMAX_RAM = "VRT_PROCESSED_DATASET_ALLOWED_RAM_USAGE";
631 46 : if (const char *pszVal = CPLGetConfigOption(pszMAX_RAM, nullptr))
632 : {
633 3 : CPL_IGNORE_RET_VAL(
634 3 : CPLParseMemorySize(pszVal, &m_nAllowedRAMUsage, nullptr));
635 : }
636 :
637 46 : if (m_nAllowedRAMUsage > 0)
638 : {
639 46 : bool bBlockSizeModified = false;
640 78 : while ((m_nBlockXSize >= 2 || m_nBlockYSize >= 2) &&
641 71 : static_cast<GIntBig>(m_nBlockXSize) * m_nBlockYSize >
642 71 : m_nAllowedRAMUsage / m_nWorkingBytesPerPixel)
643 : {
644 32 : if ((m_nBlockXSize == nRasterXSize ||
645 30 : m_nBlockYSize >= m_nBlockXSize) &&
646 17 : m_nBlockYSize >= 2)
647 : {
648 16 : m_nBlockYSize /= 2;
649 : }
650 : else
651 : {
652 16 : m_nBlockXSize /= 2;
653 : }
654 32 : bBlockSizeModified = true;
655 : }
656 46 : if (bBlockSizeModified)
657 : {
658 3 : if (bUserBlockSize)
659 : {
660 0 : CPLError(
661 : CE_Warning, CPLE_AppDefined,
662 : "Reducing block size to %d x %d to avoid consuming too "
663 : "much RAM",
664 : m_nBlockXSize, m_nBlockYSize);
665 : }
666 : else
667 : {
668 3 : CPLDebug(
669 : "VRT",
670 : "Reducing block size to %d x %d to avoid consuming too "
671 : "much RAM",
672 : m_nBlockXSize, m_nBlockYSize);
673 : }
674 : }
675 : }
676 :
677 46 : if (m_outputBandCountProvenance == ValueProvenance::FROM_LAST_STEP)
678 : {
679 7 : nOutputBandCount = nCurrentBandCount;
680 : }
681 39 : else if (nOutputBandCount != nCurrentBandCount)
682 : {
683 : // Should not happen frequently as pixel init functions are expected
684 : // to validate that they can accept the number of output bands provided
685 : // to them
686 0 : CPLError(
687 : CE_Failure, CPLE_AppDefined,
688 : "Number of output bands of last step (%d) is not consistent with "
689 : "number of VRTProcessedRasterBand's (%d)",
690 : nCurrentBandCount, nBands);
691 0 : return CE_Failure;
692 : }
693 :
694 46 : if (m_outputBandDataTypeProvenance == ValueProvenance::FROM_LAST_STEP)
695 : {
696 3 : m_outputBandDataTypeValue = eCurrentDT;
697 : }
698 :
699 53 : if (nBands != 0 &&
700 7 : (nBands != nOutputBandCount ||
701 6 : (m_outputBandDataTypeProvenance == ValueProvenance::FROM_LAST_STEP &&
702 1 : m_outputBandDataTypeValue != papoBands[0]->GetRasterDataType())))
703 : {
704 2 : for (int i = 0; i < nBands; ++i)
705 1 : delete papoBands[i];
706 1 : CPLFree(papoBands);
707 1 : papoBands = nullptr;
708 1 : nBands = 0;
709 : }
710 :
711 176 : const auto GetOutputBandType = [this, eCurrentDT](GDALDataType eSourceDT)
712 : {
713 87 : if (m_outputBandDataTypeProvenance == ValueProvenance::FROM_LAST_STEP)
714 3 : return eCurrentDT;
715 84 : else if (m_outputBandDataTypeProvenance ==
716 : ValueProvenance::USER_PROVIDED)
717 5 : return m_outputBandDataTypeValue;
718 : else
719 79 : return eSourceDT;
720 46 : };
721 :
722 46 : if (m_outputBandCountProvenance == ValueProvenance::FROM_SOURCE)
723 : {
724 112 : for (int i = 0; i < m_poSrcDS->GetRasterCount(); ++i)
725 : {
726 79 : const auto poSrcBand = m_poSrcDS->GetRasterBand(i + 1);
727 : const GDALDataType eOutputBandType =
728 79 : GetOutputBandType(poSrcBand->GetRasterDataType());
729 : auto poBand =
730 79 : new VRTProcessedRasterBand(this, i + 1, eOutputBandType);
731 79 : poBand->CopyCommonInfoFrom(poSrcBand);
732 79 : SetBand(i + 1, poBand);
733 : }
734 : }
735 13 : else if (m_outputBandCountProvenance != ValueProvenance::FROM_VRTRASTERBAND)
736 : {
737 8 : const GDALDataType eOutputBandType = GetOutputBandType(eInDT);
738 35 : for (int i = 0; i < nOutputBandCount; ++i)
739 : {
740 : auto poBand =
741 27 : new VRTProcessedRasterBand(this, i + 1, eOutputBandType);
742 27 : SetBand(i + 1, poBand);
743 : }
744 : }
745 :
746 46 : if (nBands > 1)
747 37 : SetMetadataItem("INTERLEAVE", "PIXEL", "IMAGE_STRUCTURE");
748 :
749 46 : m_oXMLTree.reset(CPLCloneXMLTree(psTree));
750 :
751 46 : return CE_None;
752 : }
753 :
754 : /************************************************************************/
755 : /* ParseStep() */
756 : /************************************************************************/
757 :
758 : /** Parse the current Step node and create a corresponding entry in m_aoSteps.
759 : *
760 : * @param psStep Step node
761 : * @param bIsFinalStep Whether this is the final step.
762 : * @param[in,out] eCurrentDT Input data type for this step.
763 : * Updated to output data type at end of method.
764 : * @param[in,out] nCurrentBandCount Input band count for this step.
765 : * Updated to output band cout at end of
766 : * method.
767 : * @param adfInNoData Input nodata values
768 : * @param[in,out] adfOutNoData Output nodata values, to be filled by this
769 : * method. When bIsFinalStep, this is also an
770 : * input parameter.
771 : * @return true on success.
772 : */
773 77 : bool VRTProcessedDataset::ParseStep(const CPLXMLNode *psStep, bool bIsFinalStep,
774 : GDALDataType &eCurrentDT,
775 : int &nCurrentBandCount,
776 : std::vector<double> &adfInNoData,
777 : std::vector<double> &adfOutNoData)
778 : {
779 77 : const char *pszStepName = CPLGetXMLValue(
780 77 : psStep, "name", CPLSPrintf("nr %d", 1 + int(m_aoSteps.size())));
781 77 : const char *pszAlgorithm = CPLGetXMLValue(psStep, "Algorithm", nullptr);
782 77 : if (!pszAlgorithm)
783 : {
784 0 : CPLError(CE_Failure, CPLE_AppDefined,
785 : "Step '%s' lacks a Algorithm element", pszStepName);
786 0 : return false;
787 : }
788 :
789 77 : const auto &oMapFunctions = GetGlobalMapProcessedDatasetFunc();
790 77 : const auto oIterFunc = oMapFunctions.find(pszAlgorithm);
791 77 : if (oIterFunc == oMapFunctions.end())
792 : {
793 0 : CPLError(CE_Failure, CPLE_AppDefined,
794 : "Step '%s' uses unregistered algorithm '%s'", pszStepName,
795 : pszAlgorithm);
796 0 : return false;
797 : }
798 :
799 77 : const auto &oFunc = oIterFunc->second;
800 :
801 77 : if (!oFunc.aeSupportedInputDT.empty())
802 : {
803 0 : if (std::find(oFunc.aeSupportedInputDT.begin(),
804 : oFunc.aeSupportedInputDT.end(),
805 0 : eCurrentDT) == oFunc.aeSupportedInputDT.end())
806 : {
807 0 : CPLError(CE_Failure, CPLE_AppDefined,
808 : "Step '%s' (using algorithm '%s') does not "
809 : "support input data type = '%s'",
810 : pszStepName, pszAlgorithm,
811 : GDALGetDataTypeName(eCurrentDT));
812 0 : return false;
813 : }
814 : }
815 :
816 77 : if (!oFunc.anSupportedInputBandCount.empty())
817 : {
818 0 : if (std::find(oFunc.anSupportedInputBandCount.begin(),
819 : oFunc.anSupportedInputBandCount.end(),
820 0 : nCurrentBandCount) ==
821 0 : oFunc.anSupportedInputBandCount.end())
822 : {
823 0 : CPLError(CE_Failure, CPLE_AppDefined,
824 : "Step '%s' (using algorithm '%s') does not "
825 : "support input band count = %d",
826 : pszStepName, pszAlgorithm, nCurrentBandCount);
827 0 : return false;
828 : }
829 : }
830 :
831 154 : Step oStep;
832 77 : oStep.osAlgorithm = pszAlgorithm;
833 154 : oStep.eInDT = oFunc.eRequestedInputDT != GDT_Unknown
834 77 : ? oFunc.eRequestedInputDT
835 : : eCurrentDT;
836 77 : oStep.nInBands = nCurrentBandCount;
837 :
838 : // Unless modified by pfnInit...
839 77 : oStep.eOutDT = oStep.eInDT;
840 :
841 77 : oStep.adfInNoData = adfInNoData;
842 77 : oStep.adfOutNoData = bIsFinalStep ? adfOutNoData : adfInNoData;
843 :
844 : // Deal with constant arguments
845 77 : for (const auto &nameValuePair : oFunc.oMapConstantArguments)
846 : {
847 : oStep.aosArguments.AddNameValue(nameValuePair.first.c_str(),
848 0 : nameValuePair.second.c_str());
849 : }
850 :
851 : // Deal with built-in arguments
852 77 : if (oFunc.oSetBuiltinArguments.find("nodata") !=
853 154 : oFunc.oSetBuiltinArguments.end())
854 : {
855 0 : int bHasVal = false;
856 0 : const auto poSrcFirstBand = m_poSrcDS->GetRasterBand(1);
857 0 : const double dfVal = poSrcFirstBand->GetNoDataValue(&bHasVal);
858 0 : if (bHasVal)
859 : {
860 : oStep.aosArguments.AddNameValue("nodata",
861 0 : CPLSPrintf("%.17g", dfVal));
862 : }
863 : }
864 :
865 77 : if (oFunc.oSetBuiltinArguments.find("offset_{band}") !=
866 154 : oFunc.oSetBuiltinArguments.end())
867 : {
868 0 : for (int i = 1; i <= m_poSrcDS->GetRasterCount(); ++i)
869 : {
870 0 : int bHasVal = false;
871 0 : const double dfVal = GetRasterBand(i)->GetOffset(&bHasVal);
872 : oStep.aosArguments.AddNameValue(
873 : CPLSPrintf("offset_%d", i),
874 0 : CPLSPrintf("%.17g", bHasVal ? dfVal : 0.0));
875 : }
876 : }
877 :
878 77 : if (oFunc.oSetBuiltinArguments.find("scale_{band}") !=
879 154 : oFunc.oSetBuiltinArguments.end())
880 : {
881 0 : for (int i = 1; i <= m_poSrcDS->GetRasterCount(); ++i)
882 : {
883 0 : int bHasVal = false;
884 0 : const double dfVal = GetRasterBand(i)->GetScale(&bHasVal);
885 : oStep.aosArguments.AddNameValue(
886 : CPLSPrintf("scale_%d", i),
887 0 : CPLSPrintf("%.17g", bHasVal ? dfVal : 1.0));
888 : }
889 : }
890 :
891 : // Parse arguments specified in VRT
892 154 : std::set<std::string> oFoundArguments;
893 :
894 421 : for (const CPLXMLNode *psStepChild = psStep->psChild; psStepChild;
895 344 : psStepChild = psStepChild->psNext)
896 : {
897 344 : if (psStepChild->eType == CXT_Element &&
898 320 : strcmp(psStepChild->pszValue, "Argument") == 0)
899 : {
900 : const char *pszParamName =
901 243 : CPLGetXMLValue(psStepChild, "name", nullptr);
902 243 : if (!pszParamName)
903 : {
904 0 : CPLError(CE_Failure, CPLE_AppDefined,
905 : "Step '%s' has a Argument without a name attribute",
906 : pszStepName);
907 0 : return false;
908 : }
909 243 : const char *pszValue = CPLGetXMLValue(psStepChild, nullptr, "");
910 : auto oOtherArgIter =
911 243 : oFunc.oOtherArguments.find(CPLString(pszParamName).tolower());
912 486 : if (!oFunc.oOtherArguments.empty() &&
913 486 : oOtherArgIter == oFunc.oOtherArguments.end())
914 : {
915 : // If we got a parameter name like 'coefficients_1',
916 : // try to fetch the generic 'coefficients_{band}'
917 296 : std::string osParamName(pszParamName);
918 148 : const auto nPos = osParamName.rfind('_');
919 148 : if (nPos != std::string::npos)
920 : {
921 148 : osParamName.resize(nPos + 1);
922 148 : osParamName += "{band}";
923 : oOtherArgIter = oFunc.oOtherArguments.find(
924 148 : CPLString(osParamName).tolower());
925 : }
926 : }
927 243 : if (oOtherArgIter != oFunc.oOtherArguments.end())
928 : {
929 243 : oFoundArguments.insert(oOtherArgIter->first);
930 :
931 243 : const std::string &osType = oOtherArgIter->second.osType;
932 243 : if (osType == "boolean")
933 : {
934 0 : if (!EQUAL(pszValue, "true") && !EQUAL(pszValue, "false"))
935 : {
936 0 : CPLError(CE_Failure, CPLE_NotSupported,
937 : "Step '%s' has a Argument '%s' whose "
938 : "value '%s' is not a boolean",
939 : pszStepName, pszParamName, pszValue);
940 0 : return false;
941 : }
942 : }
943 243 : else if (osType == "integer")
944 : {
945 44 : if (CPLGetValueType(pszValue) != CPL_VALUE_INTEGER)
946 : {
947 0 : CPLError(CE_Failure, CPLE_NotSupported,
948 : "Step '%s' has a Argument '%s' whose "
949 : "value '%s' is not a integer",
950 : pszStepName, pszParamName, pszValue);
951 0 : return false;
952 : }
953 : }
954 199 : else if (osType == "double")
955 : {
956 49 : const auto eType = CPLGetValueType(pszValue);
957 49 : if (eType != CPL_VALUE_INTEGER && eType != CPL_VALUE_REAL)
958 : {
959 0 : CPLError(CE_Failure, CPLE_NotSupported,
960 : "Step '%s' has a Argument '%s' whose "
961 : "value '%s' is not a double",
962 : pszStepName, pszParamName, pszValue);
963 0 : return false;
964 : }
965 : }
966 150 : else if (osType == "double_list")
967 : {
968 : const CPLStringList aosTokens(
969 88 : CSLTokenizeString2(pszValue, ",", 0));
970 405 : for (int i = 0; i < aosTokens.size(); ++i)
971 : {
972 317 : const auto eType = CPLGetValueType(aosTokens[i]);
973 317 : if (eType != CPL_VALUE_INTEGER &&
974 : eType != CPL_VALUE_REAL)
975 : {
976 0 : CPLError(CE_Failure, CPLE_NotSupported,
977 : "Step '%s' has a Argument '%s' "
978 : "whose value '%s' is not a "
979 : "comma-separated list of doubles",
980 : pszStepName, pszParamName, pszValue);
981 0 : return false;
982 : }
983 : }
984 : }
985 62 : else if (osType != "string")
986 : {
987 0 : CPLDebug("VRT", "Unhandled argument type '%s'",
988 : osType.c_str());
989 0 : CPLAssert(0);
990 : }
991 : }
992 0 : else if (oFunc.bMetadataSpecified &&
993 0 : oFunc.oSetBuiltinArguments.find(
994 0 : CPLString(pszParamName).tolower()) ==
995 0 : oFunc.oSetBuiltinArguments.end() &&
996 0 : oFunc.oMapConstantArguments.find(
997 0 : CPLString(pszParamName).tolower()) ==
998 0 : oFunc.oMapConstantArguments.end())
999 : {
1000 0 : CPLError(CE_Warning, CPLE_NotSupported,
1001 : "Step '%s' has a Argument '%s' which is not "
1002 : "supported",
1003 : pszStepName, pszParamName);
1004 : }
1005 :
1006 243 : oStep.aosArguments.AddNameValue(pszParamName, pszValue);
1007 : }
1008 : }
1009 :
1010 : // Check that required arguments have been specified
1011 619 : for (const auto &oIter : oFunc.oOtherArguments)
1012 : {
1013 693 : if (oIter.second.bRequired &&
1014 693 : oFoundArguments.find(oIter.first) == oFoundArguments.end())
1015 : {
1016 3 : CPLError(CE_Failure, CPLE_AppDefined,
1017 : "Step '%s' lacks required Argument '%s'", pszStepName,
1018 : oIter.first.c_str());
1019 3 : return false;
1020 : }
1021 : }
1022 :
1023 74 : if (oFunc.pfnInit)
1024 : {
1025 74 : double *padfOutNoData = nullptr;
1026 74 : if (bIsFinalStep && !adfOutNoData.empty())
1027 : {
1028 64 : oStep.nOutBands = static_cast<int>(adfOutNoData.size());
1029 64 : padfOutNoData = static_cast<double *>(
1030 64 : CPLMalloc(adfOutNoData.size() * sizeof(double)));
1031 64 : memcpy(padfOutNoData, adfOutNoData.data(),
1032 64 : adfOutNoData.size() * sizeof(double));
1033 : }
1034 : else
1035 : {
1036 10 : oStep.nOutBands = 0;
1037 : }
1038 :
1039 148 : if (oFunc.pfnInit(pszAlgorithm, oFunc.pUserData,
1040 74 : oStep.aosArguments.List(), oStep.nInBands,
1041 : oStep.eInDT, adfInNoData.data(), &(oStep.nOutBands),
1042 : &(oStep.eOutDT), &padfOutNoData, m_osVRTPath.c_str(),
1043 74 : &(oStep.pWorkingData)) != CE_None)
1044 : {
1045 26 : CPLError(CE_Failure, CPLE_AppDefined,
1046 : "Step '%s' (using algorithm '%s') init() function "
1047 : "failed",
1048 : pszStepName, pszAlgorithm);
1049 26 : CPLFree(padfOutNoData);
1050 26 : return false;
1051 : }
1052 :
1053 : // Input nodata values may have been modified by pfnInit()
1054 48 : oStep.adfInNoData = adfInNoData;
1055 :
1056 48 : if (padfOutNoData)
1057 : {
1058 88 : adfOutNoData = std::vector<double>(padfOutNoData,
1059 44 : padfOutNoData + oStep.nOutBands);
1060 : }
1061 : else
1062 : {
1063 12 : adfOutNoData = std::vector<double>(
1064 8 : oStep.nOutBands, std::numeric_limits<double>::quiet_NaN());
1065 : }
1066 48 : CPLFree(padfOutNoData);
1067 :
1068 48 : oStep.adfOutNoData = adfOutNoData;
1069 : }
1070 : else
1071 : {
1072 0 : oStep.nOutBands = oStep.nInBands;
1073 0 : adfOutNoData = oStep.adfOutNoData;
1074 : }
1075 :
1076 48 : eCurrentDT = oStep.eOutDT;
1077 48 : nCurrentBandCount = oStep.nOutBands;
1078 :
1079 48 : m_aoSteps.emplace_back(std::move(oStep));
1080 :
1081 48 : return true;
1082 : }
1083 :
1084 : /************************************************************************/
1085 : /* SerializeToXML() */
1086 : /************************************************************************/
1087 :
1088 1 : CPLXMLNode *VRTProcessedDataset::SerializeToXML(const char *pszVRTPathIn)
1089 :
1090 : {
1091 1 : CPLXMLNode *psTree = CPLCloneXMLTree(m_oXMLTree.get());
1092 1 : if (psTree == nullptr)
1093 0 : return psTree;
1094 :
1095 : /* -------------------------------------------------------------------- */
1096 : /* Remove VRTRasterBand nodes from the original tree and find the */
1097 : /* last child. */
1098 : /* -------------------------------------------------------------------- */
1099 1 : CPLXMLNode *psLastChild = psTree->psChild;
1100 1 : CPLXMLNode *psPrevChild = nullptr;
1101 5 : while (psLastChild)
1102 : {
1103 5 : CPLXMLNode *psNextChild = psLastChild->psNext;
1104 5 : if (psLastChild->eType == CXT_Element &&
1105 3 : strcmp(psLastChild->pszValue, "VRTRasterBand") == 0)
1106 : {
1107 1 : if (psPrevChild)
1108 1 : psPrevChild->psNext = psNextChild;
1109 : else
1110 0 : psTree->psChild = psNextChild;
1111 1 : psLastChild->psNext = nullptr;
1112 1 : CPLDestroyXMLNode(psLastChild);
1113 1 : psLastChild = psPrevChild ? psPrevChild : psTree->psChild;
1114 : }
1115 4 : else if (!psNextChild)
1116 : {
1117 1 : break;
1118 : }
1119 : else
1120 : {
1121 3 : psPrevChild = psLastChild;
1122 3 : psLastChild = psNextChild;
1123 : }
1124 : }
1125 1 : CPLAssert(psLastChild); // we have at least Input
1126 :
1127 : /* -------------------------------------------------------------------- */
1128 : /* Serialize bands. */
1129 : /* -------------------------------------------------------------------- */
1130 1 : bool bHasWarnedAboutRAMUsage = false;
1131 1 : size_t nAccRAMUsage = 0;
1132 2 : for (int iBand = 0; iBand < nBands; iBand++)
1133 : {
1134 : CPLXMLNode *psBandTree =
1135 1 : static_cast<VRTRasterBand *>(papoBands[iBand])
1136 2 : ->SerializeToXML(pszVRTPathIn, bHasWarnedAboutRAMUsage,
1137 1 : nAccRAMUsage);
1138 :
1139 1 : if (psBandTree != nullptr)
1140 : {
1141 1 : psLastChild->psNext = psBandTree;
1142 1 : psLastChild = psBandTree;
1143 : }
1144 : }
1145 :
1146 1 : return psTree;
1147 : }
1148 :
1149 : /************************************************************************/
1150 : /* SerializeToXML() */
1151 : /************************************************************************/
1152 :
1153 : CPLXMLNode *
1154 1 : VRTProcessedRasterBand::SerializeToXML(const char *pszVRTPathIn,
1155 : bool &bHasWarnedAboutRAMUsage,
1156 : size_t &nAccRAMUsage)
1157 :
1158 : {
1159 1 : CPLXMLNode *psTree = VRTRasterBand::SerializeToXML(
1160 : pszVRTPathIn, bHasWarnedAboutRAMUsage, nAccRAMUsage);
1161 :
1162 : /* -------------------------------------------------------------------- */
1163 : /* Set subclass. */
1164 : /* -------------------------------------------------------------------- */
1165 1 : CPLCreateXMLNode(CPLCreateXMLNode(psTree, CXT_Attribute, "subClass"),
1166 : CXT_Text, "VRTProcessedRasterBand");
1167 :
1168 1 : return psTree;
1169 : }
1170 :
1171 : /************************************************************************/
1172 : /* GetBlockSize() */
1173 : /************************************************************************/
1174 :
1175 : /** Return block size */
1176 119 : void VRTProcessedDataset::GetBlockSize(int *pnBlockXSize,
1177 : int *pnBlockYSize) const
1178 :
1179 : {
1180 119 : *pnBlockXSize = m_nBlockXSize;
1181 119 : *pnBlockYSize = m_nBlockYSize;
1182 119 : }
1183 :
1184 : /************************************************************************/
1185 : /* ProcessRegion() */
1186 : /************************************************************************/
1187 :
1188 : /** Compute pixel values for the specified region.
1189 : *
1190 : * The output is stored in m_abyInput in a pixel-interleaved way.
1191 : */
1192 60 : bool VRTProcessedDataset::ProcessRegion(int nXOff, int nYOff, int nBufXSize,
1193 : int nBufYSize,
1194 : GDALProgressFunc pfnProgress,
1195 : void *pProgressData)
1196 : {
1197 :
1198 60 : CPLAssert(!m_aoSteps.empty());
1199 :
1200 60 : const size_t nPixelCount = static_cast<size_t>(nBufXSize) * nBufYSize;
1201 :
1202 60 : const int nFirstBandCount = m_aoSteps.front().nInBands;
1203 60 : CPLAssert(nFirstBandCount == m_poSrcDS->GetRasterCount());
1204 60 : const GDALDataType eFirstDT = m_aoSteps.front().eInDT;
1205 60 : const int nFirstDTSize = GDALGetDataTypeSizeBytes(eFirstDT);
1206 60 : auto &abyInput = m_abyInput;
1207 60 : auto &abyOutput = m_abyOutput;
1208 :
1209 : const char *pszInterleave =
1210 60 : m_poSrcDS->GetMetadataItem("INTERLEAVE", "IMAGE_STRUCTURE");
1211 60 : if (nFirstBandCount > 1 && (!pszInterleave || EQUAL(pszInterleave, "BAND")))
1212 : {
1213 : // If there are several bands and the source dataset organization
1214 : // is apparently band interleaved, then first acquire data in
1215 : // a BSQ organization in the abyInput array use in the native
1216 : // data type.
1217 : // And then transpose it and convert it to the expected data type
1218 : // of the first step.
1219 1 : const auto eSrcDT = m_poSrcDS->GetRasterBand(1)->GetRasterDataType();
1220 : try
1221 : {
1222 2 : abyInput.resize(nPixelCount * nFirstBandCount *
1223 1 : GDALGetDataTypeSizeBytes(eSrcDT));
1224 : }
1225 0 : catch (const std::bad_alloc &)
1226 : {
1227 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
1228 : "Out of memory allocating working buffer");
1229 0 : return false;
1230 : }
1231 :
1232 : try
1233 : {
1234 1 : abyOutput.resize(nPixelCount * nFirstBandCount * nFirstDTSize);
1235 : }
1236 0 : catch (const std::bad_alloc &)
1237 : {
1238 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
1239 : "Out of memory allocating working buffer");
1240 0 : return false;
1241 : }
1242 :
1243 : GDALRasterIOExtraArg sArg;
1244 1 : INIT_RASTERIO_EXTRA_ARG(sArg);
1245 1 : sArg.pfnProgress = GDALScaledProgress;
1246 1 : sArg.pProgressData =
1247 1 : GDALCreateScaledProgress(0, 0.5, pfnProgress, pProgressData);
1248 1 : if (sArg.pProgressData == nullptr)
1249 1 : sArg.pfnProgress = nullptr;
1250 :
1251 1 : CPLDebugOnly("VRT", "ProcessRegion(): start RasterIO()");
1252 : const bool bOK =
1253 2 : m_poSrcDS->RasterIO(GF_Read, nXOff, nYOff, nBufXSize, nBufYSize,
1254 1 : abyInput.data(), nBufXSize, nBufYSize, eSrcDT,
1255 : nFirstBandCount, nullptr, 0, 0, 0,
1256 1 : &sArg) == CE_None;
1257 1 : CPLDebugOnly("VRT", "ProcessRegion(): end RasterIO()");
1258 1 : GDALDestroyScaledProgress(sArg.pProgressData);
1259 1 : if (!bOK)
1260 0 : return false;
1261 :
1262 1 : CPLDebugOnly("VRT", "ProcessRegion(): start GDALTranspose2D()");
1263 1 : GDALTranspose2D(abyInput.data(), eSrcDT, abyOutput.data(), eFirstDT,
1264 1 : static_cast<size_t>(nBufXSize) * nBufYSize,
1265 : nFirstBandCount);
1266 1 : CPLDebugOnly("VRT", "ProcessRegion(): end GDALTranspose2D()");
1267 :
1268 : // Swap arrays
1269 1 : std::swap(abyInput, abyOutput);
1270 : }
1271 : else
1272 : {
1273 : try
1274 : {
1275 59 : abyInput.resize(nPixelCount * nFirstBandCount * nFirstDTSize);
1276 : }
1277 0 : catch (const std::bad_alloc &)
1278 : {
1279 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
1280 : "Out of memory allocating working buffer");
1281 0 : return false;
1282 : }
1283 :
1284 : GDALRasterIOExtraArg sArg;
1285 59 : INIT_RASTERIO_EXTRA_ARG(sArg);
1286 59 : sArg.pfnProgress = GDALScaledProgress;
1287 59 : sArg.pProgressData =
1288 59 : GDALCreateScaledProgress(0, 0.5, pfnProgress, pProgressData);
1289 59 : if (sArg.pProgressData == nullptr)
1290 59 : sArg.pfnProgress = nullptr;
1291 :
1292 : const bool bOK =
1293 118 : m_poSrcDS->RasterIO(
1294 59 : GF_Read, nXOff, nYOff, nBufXSize, nBufYSize, abyInput.data(),
1295 : nBufXSize, nBufYSize, eFirstDT, nFirstBandCount, nullptr,
1296 59 : static_cast<GSpacing>(nFirstDTSize) * nFirstBandCount,
1297 59 : static_cast<GSpacing>(nFirstDTSize) * nFirstBandCount *
1298 59 : nBufXSize,
1299 59 : nFirstDTSize, &sArg) == CE_None;
1300 :
1301 59 : GDALDestroyScaledProgress(sArg.pProgressData);
1302 59 : if (!bOK)
1303 3 : return false;
1304 : }
1305 :
1306 57 : const double dfSrcXOff = nXOff;
1307 57 : const double dfSrcYOff = nYOff;
1308 57 : const double dfSrcXSize = nBufXSize;
1309 57 : const double dfSrcYSize = nBufYSize;
1310 :
1311 : double adfSrcGT[6];
1312 57 : if (m_poSrcDS->GetGeoTransform(adfSrcGT) != CE_None)
1313 : {
1314 38 : adfSrcGT[0] = 0;
1315 38 : adfSrcGT[1] = 1;
1316 38 : adfSrcGT[2] = 0;
1317 38 : adfSrcGT[3] = 0;
1318 38 : adfSrcGT[4] = 0;
1319 38 : adfSrcGT[5] = 1;
1320 : }
1321 :
1322 57 : GDALDataType eLastDT = eFirstDT;
1323 57 : const auto &oMapFunctions = GetGlobalMapProcessedDatasetFunc();
1324 :
1325 57 : int iStep = 0;
1326 116 : for (const auto &oStep : m_aoSteps)
1327 : {
1328 59 : const auto oIterFunc = oMapFunctions.find(oStep.osAlgorithm);
1329 59 : CPLAssert(oIterFunc != oMapFunctions.end());
1330 :
1331 : // Data type adaptation
1332 59 : if (eLastDT != oStep.eInDT)
1333 : {
1334 : try
1335 : {
1336 0 : abyOutput.resize(nPixelCount * oStep.nInBands *
1337 0 : GDALGetDataTypeSizeBytes(oStep.eInDT));
1338 : }
1339 0 : catch (const std::bad_alloc &)
1340 : {
1341 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
1342 : "Out of memory allocating working buffer");
1343 0 : return false;
1344 : }
1345 :
1346 0 : GDALCopyWords64(abyInput.data(), eLastDT,
1347 0 : GDALGetDataTypeSizeBytes(eLastDT), abyOutput.data(),
1348 0 : oStep.eInDT, GDALGetDataTypeSizeBytes(oStep.eInDT),
1349 0 : nPixelCount * oStep.nInBands);
1350 :
1351 0 : std::swap(abyInput, abyOutput);
1352 : }
1353 :
1354 : try
1355 : {
1356 118 : abyOutput.resize(nPixelCount * oStep.nOutBands *
1357 59 : GDALGetDataTypeSizeBytes(oStep.eOutDT));
1358 : }
1359 0 : catch (const std::bad_alloc &)
1360 : {
1361 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
1362 : "Out of memory allocating working buffer");
1363 0 : return false;
1364 : }
1365 :
1366 59 : const auto &oFunc = oIterFunc->second;
1367 236 : if (oFunc.pfnProcess(
1368 59 : oStep.osAlgorithm.c_str(), oFunc.pUserData, oStep.pWorkingData,
1369 : oStep.aosArguments.List(), nBufXSize, nBufYSize,
1370 59 : abyInput.data(), abyInput.size(), oStep.eInDT, oStep.nInBands,
1371 59 : oStep.adfInNoData.data(), abyOutput.data(), abyOutput.size(),
1372 59 : oStep.eOutDT, oStep.nOutBands, oStep.adfOutNoData.data(),
1373 : dfSrcXOff, dfSrcYOff, dfSrcXSize, dfSrcYSize, adfSrcGT,
1374 : m_osVRTPath.c_str(),
1375 59 : /*papszExtra=*/nullptr) != CE_None)
1376 : {
1377 0 : return false;
1378 : }
1379 :
1380 59 : std::swap(abyInput, abyOutput);
1381 59 : eLastDT = oStep.eOutDT;
1382 :
1383 59 : ++iStep;
1384 59 : if (pfnProgress &&
1385 0 : !pfnProgress(0.5 + 0.5 * iStep / static_cast<int>(m_aoSteps.size()),
1386 : "", pProgressData))
1387 0 : return false;
1388 : }
1389 :
1390 57 : return true;
1391 : }
1392 :
1393 : /************************************************************************/
1394 : /* VRTProcessedRasterBand() */
1395 : /************************************************************************/
1396 :
1397 : /** Constructor */
1398 119 : VRTProcessedRasterBand::VRTProcessedRasterBand(VRTProcessedDataset *poDSIn,
1399 : int nBandIn,
1400 119 : GDALDataType eDataTypeIn)
1401 : {
1402 119 : Initialize(poDSIn->GetRasterXSize(), poDSIn->GetRasterYSize());
1403 :
1404 119 : poDS = poDSIn;
1405 119 : nBand = nBandIn;
1406 119 : eAccess = GA_Update;
1407 119 : eDataType = eDataTypeIn;
1408 :
1409 119 : poDSIn->GetBlockSize(&nBlockXSize, &nBlockYSize);
1410 119 : }
1411 :
1412 : /************************************************************************/
1413 : /* GetOverviewCount() */
1414 : /************************************************************************/
1415 :
1416 6 : int VRTProcessedRasterBand::GetOverviewCount()
1417 : {
1418 6 : auto poVRTDS = cpl::down_cast<VRTProcessedDataset *>(poDS);
1419 6 : return static_cast<int>(poVRTDS->m_apoOverviewDatasets.size());
1420 : }
1421 :
1422 : /************************************************************************/
1423 : /* GetOverview() */
1424 : /************************************************************************/
1425 :
1426 11 : GDALRasterBand *VRTProcessedRasterBand::GetOverview(int iOvr)
1427 : {
1428 11 : auto poVRTDS = cpl::down_cast<VRTProcessedDataset *>(poDS);
1429 22 : if (iOvr < 0 ||
1430 11 : iOvr >= static_cast<int>(poVRTDS->m_apoOverviewDatasets.size()))
1431 0 : return nullptr;
1432 11 : return poVRTDS->m_apoOverviewDatasets[iOvr]->GetRasterBand(nBand);
1433 : }
1434 :
1435 : /************************************************************************/
1436 : /* IReadBlock() */
1437 : /************************************************************************/
1438 :
1439 39 : CPLErr VRTProcessedRasterBand::IReadBlock(int nBlockXOff, int nBlockYOff,
1440 : void *pImage)
1441 :
1442 : {
1443 39 : auto poVRTDS = cpl::down_cast<VRTProcessedDataset *>(poDS);
1444 :
1445 39 : int nBufXSize = 0;
1446 39 : int nBufYSize = 0;
1447 39 : GetActualBlockSize(nBlockXOff, nBlockYOff, &nBufXSize, &nBufYSize);
1448 :
1449 39 : const int nXPixelOff = nBlockXOff * nBlockXSize;
1450 39 : const int nYPixelOff = nBlockYOff * nBlockYSize;
1451 39 : if (!poVRTDS->ProcessRegion(nXPixelOff, nYPixelOff, nBufXSize, nBufYSize,
1452 : nullptr, nullptr))
1453 : {
1454 0 : return CE_Failure;
1455 : }
1456 :
1457 39 : const int nOutBands = poVRTDS->m_aoSteps.back().nOutBands;
1458 39 : CPLAssert(nOutBands == poVRTDS->GetRasterCount());
1459 39 : const auto eLastDT = poVRTDS->m_aoSteps.back().eOutDT;
1460 39 : const int nLastDTSize = GDALGetDataTypeSizeBytes(eLastDT);
1461 39 : const int nDTSize = GDALGetDataTypeSizeBytes(eDataType);
1462 :
1463 : // Dispatch final output buffer to cached blocks of output bands
1464 114 : for (int iDstBand = 0; iDstBand < nOutBands; ++iDstBand)
1465 : {
1466 75 : GDALRasterBlock *poBlock = nullptr;
1467 : GByte *pDst;
1468 75 : if (iDstBand + 1 == nBand)
1469 : {
1470 39 : pDst = static_cast<GByte *>(pImage);
1471 : }
1472 : else
1473 : {
1474 36 : auto poOtherBand = poVRTDS->papoBands[iDstBand];
1475 36 : poBlock = poOtherBand->TryGetLockedBlockRef(nBlockXOff, nBlockYOff);
1476 36 : if (poBlock)
1477 : {
1478 0 : poBlock->DropLock();
1479 0 : continue;
1480 : }
1481 72 : poBlock = poOtherBand->GetLockedBlockRef(
1482 36 : nBlockXOff, nBlockYOff, /* bJustInitialized = */ true);
1483 36 : if (!poBlock)
1484 0 : continue;
1485 36 : pDst = static_cast<GByte *>(poBlock->GetDataRef());
1486 : }
1487 847 : for (int iY = 0; iY < nBufYSize; ++iY)
1488 : {
1489 1544 : GDALCopyWords64(poVRTDS->m_abyInput.data() +
1490 772 : (iDstBand + static_cast<size_t>(iY) *
1491 772 : nBufXSize * nOutBands) *
1492 772 : nLastDTSize,
1493 : eLastDT, nLastDTSize * nOutBands,
1494 772 : pDst +
1495 772 : static_cast<size_t>(iY) * nBlockXSize * nDTSize,
1496 : eDataType, nDTSize, nBufXSize);
1497 : }
1498 75 : if (poBlock)
1499 36 : poBlock->DropLock();
1500 : }
1501 :
1502 39 : return CE_None;
1503 : }
1504 :
1505 : /************************************************************************/
1506 : /* VRTProcessedDataset::IRasterIO() */
1507 : /************************************************************************/
1508 :
1509 52 : CPLErr VRTProcessedDataset::IRasterIO(
1510 : GDALRWFlag eRWFlag, int nXOff, int nYOff, int nXSize, int nYSize,
1511 : void *pData, int nBufXSize, int nBufYSize, GDALDataType eBufType,
1512 : int nBandCount, BANDMAP_TYPE panBandMap, GSpacing nPixelSpace,
1513 : GSpacing nLineSpace, GSpacing nBandSpace, GDALRasterIOExtraArg *psExtraArg)
1514 : {
1515 : // Try to pass the request to the most appropriate overview dataset.
1516 52 : if (nBufXSize < nXSize && nBufYSize < nYSize)
1517 : {
1518 2 : int bTried = FALSE;
1519 2 : const CPLErr eErr = TryOverviewRasterIO(
1520 : eRWFlag, nXOff, nYOff, nXSize, nYSize, pData, nBufXSize, nBufYSize,
1521 : eBufType, nBandCount, panBandMap, nPixelSpace, nLineSpace,
1522 : nBandSpace, psExtraArg, &bTried);
1523 2 : if (bTried)
1524 2 : return eErr;
1525 : }
1526 :
1527 : // Optimize reading of all bands at nominal resolution for BIP-like or
1528 : // BSQ-like buffer spacing.
1529 50 : if (eRWFlag == GF_Read && nXSize == nBufXSize && nYSize == nBufYSize &&
1530 49 : nBandCount == nBands)
1531 : {
1532 474 : const auto IsSequentialBandMap = [panBandMap, nBandCount]()
1533 : {
1534 236 : for (int i = 0; i < nBandCount; ++i)
1535 : {
1536 190 : if (panBandMap[i] != i + 1)
1537 : {
1538 2 : return false;
1539 : }
1540 : }
1541 46 : return true;
1542 48 : };
1543 :
1544 48 : const int nBufTypeSize = GDALGetDataTypeSizeBytes(eBufType);
1545 : const bool bIsBIPLike =
1546 9 : nBandSpace == nBufTypeSize && nPixelSpace == nBandSpace * nBands &&
1547 57 : nLineSpace >= nPixelSpace * nBufXSize && IsSequentialBandMap();
1548 91 : const bool bIsBSQLike = nPixelSpace == nBufTypeSize &&
1549 43 : nLineSpace >= nPixelSpace * nBufXSize &&
1550 134 : nBandSpace >= nLineSpace * nBufYSize &&
1551 43 : IsSequentialBandMap();
1552 48 : if (bIsBIPLike || bIsBSQLike)
1553 : {
1554 46 : GByte *pabyData = static_cast<GByte *>(pData);
1555 : // If acquiring the region of interest in a single time is going
1556 : // to consume too much RAM, split in halves.
1557 46 : if (m_nAllowedRAMUsage > 0 &&
1558 46 : static_cast<GIntBig>(nBufXSize) * nBufYSize >
1559 46 : m_nAllowedRAMUsage / m_nWorkingBytesPerPixel)
1560 : {
1561 25 : if ((nBufXSize == nRasterXSize || nBufYSize >= nBufXSize) &&
1562 : nBufYSize >= 2)
1563 : {
1564 : GDALRasterIOExtraArg sArg;
1565 12 : INIT_RASTERIO_EXTRA_ARG(sArg);
1566 12 : const int nHalfHeight = nBufYSize / 2;
1567 :
1568 12 : sArg.pfnProgress = GDALScaledProgress;
1569 12 : sArg.pProgressData = GDALCreateScaledProgress(
1570 : 0, 0.5, psExtraArg->pfnProgress,
1571 : psExtraArg->pProgressData);
1572 12 : if (sArg.pProgressData == nullptr)
1573 12 : sArg.pfnProgress = nullptr;
1574 : bool bOK =
1575 12 : IRasterIO(eRWFlag, nXOff, nYOff, nBufXSize, nHalfHeight,
1576 : pabyData, nBufXSize, nHalfHeight, eBufType,
1577 : nBandCount, panBandMap, nPixelSpace,
1578 12 : nLineSpace, nBandSpace, &sArg) == CE_None;
1579 12 : GDALDestroyScaledProgress(sArg.pProgressData);
1580 :
1581 12 : if (bOK)
1582 : {
1583 2 : sArg.pfnProgress = GDALScaledProgress;
1584 2 : sArg.pProgressData = GDALCreateScaledProgress(
1585 : 0.5, 1, psExtraArg->pfnProgress,
1586 : psExtraArg->pProgressData);
1587 2 : if (sArg.pProgressData == nullptr)
1588 2 : sArg.pfnProgress = nullptr;
1589 4 : bOK = IRasterIO(eRWFlag, nXOff, nYOff + nHalfHeight,
1590 : nBufXSize, nBufYSize - nHalfHeight,
1591 2 : pabyData + nHalfHeight * nLineSpace,
1592 : nBufXSize, nBufYSize - nHalfHeight,
1593 : eBufType, nBandCount, panBandMap,
1594 : nPixelSpace, nLineSpace, nBandSpace,
1595 : &sArg) == CE_None;
1596 2 : GDALDestroyScaledProgress(sArg.pProgressData);
1597 : }
1598 12 : return bOK ? CE_None : CE_Failure;
1599 : }
1600 13 : else if (nBufXSize >= 2)
1601 : {
1602 : GDALRasterIOExtraArg sArg;
1603 13 : INIT_RASTERIO_EXTRA_ARG(sArg);
1604 13 : const int nHalfWidth = nBufXSize / 2;
1605 :
1606 13 : sArg.pfnProgress = GDALScaledProgress;
1607 13 : sArg.pProgressData = GDALCreateScaledProgress(
1608 : 0, 0.5, psExtraArg->pfnProgress,
1609 : psExtraArg->pProgressData);
1610 13 : if (sArg.pProgressData == nullptr)
1611 13 : sArg.pfnProgress = nullptr;
1612 : bool bOK =
1613 13 : IRasterIO(eRWFlag, nXOff, nYOff, nHalfWidth, nBufYSize,
1614 : pabyData, nHalfWidth, nBufYSize, eBufType,
1615 : nBandCount, panBandMap, nPixelSpace,
1616 13 : nLineSpace, nBandSpace, &sArg) == CE_None;
1617 13 : GDALDestroyScaledProgress(sArg.pProgressData);
1618 :
1619 13 : if (bOK)
1620 : {
1621 3 : sArg.pfnProgress = GDALScaledProgress;
1622 3 : sArg.pProgressData = GDALCreateScaledProgress(
1623 : 0.5, 1, psExtraArg->pfnProgress,
1624 : psExtraArg->pProgressData);
1625 3 : if (sArg.pProgressData == nullptr)
1626 3 : sArg.pfnProgress = nullptr;
1627 6 : bOK = IRasterIO(eRWFlag, nXOff + nHalfWidth, nYOff,
1628 : nBufXSize - nHalfWidth, nBufYSize,
1629 3 : pabyData + nHalfWidth * nPixelSpace,
1630 : nBufXSize - nHalfWidth, nBufYSize,
1631 : eBufType, nBandCount, panBandMap,
1632 : nPixelSpace, nLineSpace, nBandSpace,
1633 : &sArg) == CE_None;
1634 3 : GDALDestroyScaledProgress(sArg.pProgressData);
1635 : }
1636 13 : return bOK ? CE_None : CE_Failure;
1637 : }
1638 : }
1639 :
1640 21 : if (!ProcessRegion(nXOff, nYOff, nBufXSize, nBufYSize,
1641 : psExtraArg->pfnProgress,
1642 : psExtraArg->pProgressData))
1643 : {
1644 3 : return CE_Failure;
1645 : }
1646 18 : const auto eLastDT = m_aoSteps.back().eOutDT;
1647 18 : const int nLastDTSize = GDALGetDataTypeSizeBytes(eLastDT);
1648 18 : if (bIsBIPLike)
1649 : {
1650 10 : for (int iY = 0; iY < nBufYSize; ++iY)
1651 : {
1652 21 : GDALCopyWords64(
1653 14 : m_abyInput.data() + static_cast<size_t>(iY) * nBands *
1654 7 : nBufXSize * nLastDTSize,
1655 7 : eLastDT, nLastDTSize, pabyData + iY * nLineSpace,
1656 : eBufType, GDALGetDataTypeSizeBytes(eBufType),
1657 7 : static_cast<size_t>(nBufXSize) * nBands);
1658 : }
1659 : }
1660 : else
1661 : {
1662 15 : CPLAssert(bIsBSQLike);
1663 79 : for (int iBand = 0; iBand < nBands; ++iBand)
1664 : {
1665 148 : for (int iY = 0; iY < nBufYSize; ++iY)
1666 : {
1667 168 : GDALCopyWords64(
1668 168 : m_abyInput.data() +
1669 84 : (static_cast<size_t>(iY) * nBands * nBufXSize +
1670 84 : iBand) *
1671 84 : nLastDTSize,
1672 84 : eLastDT, nLastDTSize * nBands,
1673 84 : pabyData + iBand * nBandSpace + iY * nLineSpace,
1674 : eBufType, nBufTypeSize, nBufXSize);
1675 : }
1676 : }
1677 : }
1678 18 : return CE_None;
1679 : }
1680 : }
1681 :
1682 4 : return VRTDataset::IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize, pData,
1683 : nBufXSize, nBufYSize, eBufType, nBandCount,
1684 : panBandMap, nPixelSpace, nLineSpace,
1685 4 : nBandSpace, psExtraArg);
1686 : }
1687 :
1688 : /*! @endcond */
1689 :
1690 : /************************************************************************/
1691 : /* GDALVRTRegisterProcessedDatasetFunc() */
1692 : /************************************************************************/
1693 :
1694 : /** Register a function to be used by VRTProcessedDataset.
1695 :
1696 : An example of content for pszXMLMetadata is:
1697 : \verbatim
1698 : <ProcessedDatasetFunctionArgumentsList>
1699 : <Argument name='src_nodata' type='double' description='Override input nodata value'/>
1700 : <Argument name='dst_nodata' type='double' description='Override output nodata value'/>
1701 : <Argument name='replacement_nodata' description='value to substitute to a valid computed value that would be nodata' type='double'/>
1702 : <Argument name='dst_intended_datatype' type='string' description='Intented datatype of output (which might be different than the working data type)'/>
1703 : <Argument name='coefficients_{band}' description='Comma-separated coefficients for combining bands. First one is constant term' type='double_list' required='true'/>
1704 : </ProcessedDatasetFunctionArgumentsList>
1705 : \endverbatim
1706 :
1707 : @param pszFuncName Function name. Must be unique and not null.
1708 : @param pUserData User data. May be nullptr. Must remain valid during the
1709 : lifetime of GDAL.
1710 : @param pszXMLMetadata XML metadata describing the function arguments. May be
1711 : nullptr if there are no arguments.
1712 : @param eRequestedInputDT If the pfnProcess callback only supports a single
1713 : data type, it should be specified in this parameter.
1714 : Otherwise set it to GDT_Unknown.
1715 : @param paeSupportedInputDT List of supported input data types. May be nullptr
1716 : if all are supported or if eRequestedInputDT is
1717 : set to a non GDT_Unknown value.
1718 : @param nSupportedInputDTSize Size of paeSupportedInputDT
1719 : @param panSupportedInputBandCount List of supported band count. May be nullptr
1720 : if any source band count is supported.
1721 : @param nSupportedInputBandCountSize Size of panSupportedInputBandCount
1722 : @param pfnInit Initialization function called when a VRTProcessedDataset
1723 : step uses the register function. This initialization function
1724 : will return the output data type, output band count and
1725 : potentially initialize a working structure, typically parsing
1726 : arguments. May be nullptr.
1727 : If not specified, it will be assumed that the input and output
1728 : data types are the same, and that the input number of bands
1729 : and output number of bands are the same.
1730 : @param pfnFree Free function that will free the working structure allocated
1731 : by pfnInit. May be nullptr.
1732 : @param pfnProcess Processing function called to compute pixel values. Must
1733 : not be nullptr.
1734 : @param papszOptions Unused currently. Must be nullptr.
1735 : @return CE_None in case of success, error otherwise.
1736 : @since 3.9
1737 : */
1738 6900 : CPLErr GDALVRTRegisterProcessedDatasetFunc(
1739 : const char *pszFuncName, void *pUserData, const char *pszXMLMetadata,
1740 : GDALDataType eRequestedInputDT, const GDALDataType *paeSupportedInputDT,
1741 : size_t nSupportedInputDTSize, const int *panSupportedInputBandCount,
1742 : size_t nSupportedInputBandCountSize,
1743 : GDALVRTProcessedDatasetFuncInit pfnInit,
1744 : GDALVRTProcessedDatasetFuncFree pfnFree,
1745 : GDALVRTProcessedDatasetFuncProcess pfnProcess,
1746 : CPL_UNUSED CSLConstList papszOptions)
1747 : {
1748 6900 : if (pszFuncName == nullptr || pszFuncName[0] == '\0')
1749 : {
1750 0 : CPLError(CE_Failure, CPLE_AppDefined,
1751 : "pszFuncName should be non-empty");
1752 0 : return CE_Failure;
1753 : }
1754 :
1755 6900 : auto &oMap = GetGlobalMapProcessedDatasetFunc();
1756 6900 : if (oMap.find(pszFuncName) != oMap.end())
1757 : {
1758 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s already registered",
1759 : pszFuncName);
1760 0 : return CE_Failure;
1761 : }
1762 :
1763 6900 : if (!pfnProcess)
1764 : {
1765 0 : CPLError(CE_Failure, CPLE_AppDefined, "pfnProcess should not be null");
1766 0 : return CE_Failure;
1767 : }
1768 :
1769 13800 : VRTProcessedDatasetFunc oFunc;
1770 6900 : oFunc.osFuncName = pszFuncName;
1771 6900 : oFunc.pUserData = pUserData;
1772 6900 : if (pszXMLMetadata)
1773 : {
1774 6900 : oFunc.bMetadataSpecified = true;
1775 6900 : auto psTree = CPLXMLTreeCloser(CPLParseXMLString(pszXMLMetadata));
1776 6900 : if (!psTree)
1777 : {
1778 0 : CPLError(CE_Failure, CPLE_AppDefined,
1779 : "Cannot parse pszXMLMetadata=%s for %s", pszXMLMetadata,
1780 : pszFuncName);
1781 0 : return CE_Failure;
1782 : }
1783 6900 : const CPLXMLNode *psRoot = CPLGetXMLNode(
1784 : psTree.get(), "=ProcessedDatasetFunctionArgumentsList");
1785 6900 : if (!psRoot)
1786 : {
1787 0 : CPLError(CE_Failure, CPLE_AppDefined,
1788 : "No root ProcessedDatasetFunctionArgumentsList element in "
1789 : "pszXMLMetadata=%s for %s",
1790 : pszXMLMetadata, pszFuncName);
1791 0 : return CE_Failure;
1792 : }
1793 52440 : for (const CPLXMLNode *psIter = psRoot->psChild; psIter;
1794 45540 : psIter = psIter->psNext)
1795 : {
1796 45540 : if (psIter->eType == CXT_Element &&
1797 45540 : strcmp(psIter->pszValue, "Argument") == 0)
1798 : {
1799 45540 : const char *pszName = CPLGetXMLValue(psIter, "name", nullptr);
1800 45540 : if (!pszName)
1801 : {
1802 0 : CPLError(CE_Failure, CPLE_AppDefined,
1803 : "Missing Argument.name attribute in "
1804 : "pszXMLMetadata=%s for %s",
1805 : pszXMLMetadata, pszFuncName);
1806 0 : return CE_Failure;
1807 : }
1808 45540 : const char *pszType = CPLGetXMLValue(psIter, "type", nullptr);
1809 45540 : if (!pszType)
1810 : {
1811 0 : CPLError(CE_Failure, CPLE_AppDefined,
1812 : "Missing Argument.type attribute in "
1813 : "pszXMLMetadata=%s for %s",
1814 : pszXMLMetadata, pszFuncName);
1815 0 : return CE_Failure;
1816 : }
1817 45540 : if (strcmp(pszType, "constant") == 0)
1818 : {
1819 : const char *pszValue =
1820 0 : CPLGetXMLValue(psIter, "value", nullptr);
1821 0 : if (!pszValue)
1822 : {
1823 0 : CPLError(CE_Failure, CPLE_AppDefined,
1824 : "Missing Argument.value attribute in "
1825 : "pszXMLMetadata=%s for %s",
1826 : pszXMLMetadata, pszFuncName);
1827 0 : return CE_Failure;
1828 : }
1829 0 : oFunc.oMapConstantArguments[CPLString(pszName).tolower()] =
1830 0 : pszValue;
1831 : }
1832 45540 : else if (strcmp(pszType, "builtin") == 0)
1833 : {
1834 0 : if (EQUAL(pszName, "nodata") ||
1835 0 : EQUAL(pszName, "offset_{band}") ||
1836 0 : EQUAL(pszName, "scale_{band}"))
1837 : {
1838 0 : oFunc.oSetBuiltinArguments.insert(
1839 0 : CPLString(pszName).tolower());
1840 : }
1841 : else
1842 : {
1843 0 : CPLError(CE_Failure, CPLE_NotSupported,
1844 : "Unsupported builtin parameter name %s in "
1845 : "pszXMLMetadata=%s for %s. Only nodata, "
1846 : "offset_{band} and scale_{band} are supported",
1847 : pszName, pszXMLMetadata, pszFuncName);
1848 0 : return CE_Failure;
1849 : }
1850 : }
1851 45540 : else if (strcmp(pszType, "boolean") == 0 ||
1852 42780 : strcmp(pszType, "string") == 0 ||
1853 33120 : strcmp(pszType, "integer") == 0 ||
1854 24840 : strcmp(pszType, "double") == 0 ||
1855 1380 : strcmp(pszType, "double_list") == 0)
1856 : {
1857 45540 : VRTProcessedDatasetFunc::OtherArgument otherArgument;
1858 45540 : otherArgument.bRequired = CPLTestBool(
1859 : CPLGetXMLValue(psIter, "required", "false"));
1860 45540 : otherArgument.osType = pszType;
1861 91080 : oFunc.oOtherArguments[CPLString(pszName).tolower()] =
1862 136620 : std::move(otherArgument);
1863 : }
1864 : else
1865 : {
1866 0 : CPLError(CE_Failure, CPLE_NotSupported,
1867 : "Unsupported type for parameter %s in "
1868 : "pszXMLMetadata=%s for %s. Only boolean, string, "
1869 : "integer, double and double_list are supported",
1870 : pszName, pszXMLMetadata, pszFuncName);
1871 0 : return CE_Failure;
1872 : }
1873 : }
1874 : }
1875 : }
1876 6900 : oFunc.eRequestedInputDT = eRequestedInputDT;
1877 6900 : if (nSupportedInputDTSize)
1878 : {
1879 : oFunc.aeSupportedInputDT.insert(
1880 0 : oFunc.aeSupportedInputDT.end(), paeSupportedInputDT,
1881 0 : paeSupportedInputDT + nSupportedInputDTSize);
1882 : }
1883 6900 : if (nSupportedInputBandCountSize)
1884 : {
1885 : oFunc.anSupportedInputBandCount.insert(
1886 0 : oFunc.anSupportedInputBandCount.end(), panSupportedInputBandCount,
1887 0 : panSupportedInputBandCount + nSupportedInputBandCountSize);
1888 : }
1889 6900 : oFunc.pfnInit = pfnInit;
1890 6900 : oFunc.pfnFree = pfnFree;
1891 6900 : oFunc.pfnProcess = pfnProcess;
1892 :
1893 6900 : oMap[pszFuncName] = std::move(oFunc);
1894 :
1895 6900 : return CE_None;
1896 : }
|