Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: GDAL
4 : * Purpose: gdal "raster compare" subcommand
5 : * Author: Even Rouault <even dot rouault at spatialys.com>
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2025, Even Rouault <even dot rouault at spatialys.com>
9 : *
10 : * SPDX-License-Identifier: MIT
11 : ****************************************************************************/
12 :
13 : #include "gdalalg_raster_compare.h"
14 :
15 : #include "cpl_conv.h"
16 : #include "cpl_vsi_virtual.h"
17 : #include "gdal_alg.h"
18 : #include "gdal_priv.h"
19 :
20 : #include <algorithm>
21 : #include <array>
22 : #include <cmath>
23 :
24 : //! @cond Doxygen_Suppress
25 :
26 : #ifndef _
27 : #define _(x) (x)
28 : #endif
29 :
30 : /************************************************************************/
31 : /* GDALRasterCompareAlgorithm::GDALRasterCompareAlgorithm() */
32 : /************************************************************************/
33 :
34 407 : GDALRasterCompareAlgorithm::GDALRasterCompareAlgorithm(bool standaloneStep)
35 : : GDALRasterPipelineStepAlgorithm(
36 : NAME, DESCRIPTION, HELP_URL,
37 0 : ConstructorOptions()
38 407 : .SetStandaloneStep(standaloneStep)
39 407 : .SetInputDatasetMaxCount(1)
40 814 : .SetInputDatasetHelpMsg(_("Input raster dataset"))
41 1221 : .SetAddDefaultArguments(false))
42 : {
43 407 : AddProgressArg();
44 :
45 : auto &referenceDatasetArg = AddArg("reference", 0, _("Reference dataset"),
46 814 : &m_referenceDataset, GDAL_OF_RASTER)
47 407 : .SetPositional()
48 407 : .SetRequired();
49 :
50 407 : SetAutoCompleteFunctionForFilename(referenceDatasetArg, GDAL_OF_RASTER);
51 :
52 407 : if (standaloneStep)
53 : {
54 187 : AddRasterInputArgs(/* openForMixedRasterVector = */ false,
55 : /* hiddenForCLI = */ false);
56 : }
57 : else
58 : {
59 220 : AddRasterHiddenInputDatasetArg();
60 : }
61 :
62 : AddArg("skip-all-optional", 0, _("Skip all optional comparisons"),
63 407 : &m_skipAllOptional);
64 407 : AddArg("skip-binary", 0, _("Skip binary file comparison"), &m_skipBinary);
65 407 : AddArg("skip-crs", 0, _("Skip CRS comparison"), &m_skipCRS);
66 : AddArg("skip-geotransform", 0, _("Skip geotransform comparison"),
67 407 : &m_skipGeotransform);
68 407 : AddArg("skip-overview", 0, _("Skip overview comparison"), &m_skipOverview);
69 407 : AddArg("skip-metadata", 0, _("Skip metadata comparison"), &m_skipMetadata);
70 407 : AddArg("skip-rpc", 0, _("Skip RPC metadata comparison"), &m_skipRPC);
71 : AddArg("skip-geolocation", 0, _("Skip Geolocation metadata comparison"),
72 407 : &m_skipGeolocation);
73 : AddArg("skip-subdataset", 0, _("Skip subdataset comparison"),
74 407 : &m_skipSubdataset);
75 :
76 407 : AddOutputStringArg(&m_output);
77 :
78 814 : AddArg("return-code", 0, _("Return code"), &m_retCode)
79 407 : .SetHiddenForCLI()
80 407 : .SetIsInput(false)
81 407 : .SetIsOutput(true);
82 407 : }
83 :
84 : /************************************************************************/
85 : /* GDALRasterCompareAlgorithm::BinaryComparison() */
86 : /************************************************************************/
87 :
88 12 : bool GDALRasterCompareAlgorithm::BinaryComparison(
89 : std::vector<std::string> &aosReport, GDALDataset *poRefDS,
90 : GDALDataset *poInputDS)
91 : {
92 12 : if (poRefDS->GetDescription()[0] == 0)
93 : {
94 1 : ReportError(
95 : CE_Warning, CPLE_AppDefined,
96 : "Reference dataset has no name. Skipping binary file comparison");
97 1 : return false;
98 : }
99 :
100 11 : auto poRefDrv = poRefDS->GetDriver();
101 11 : if (poRefDrv && EQUAL(poRefDrv->GetDescription(), "MEM"))
102 : {
103 1 : ReportError(CE_Warning, CPLE_AppDefined,
104 : "Reference dataset is a in-memory dataset. Skipping binary "
105 : "file comparison");
106 1 : return false;
107 : }
108 :
109 10 : if (poInputDS->GetDescription()[0] == 0)
110 : {
111 1 : ReportError(
112 : CE_Warning, CPLE_AppDefined,
113 : "Input dataset has no name. Skipping binary file comparison");
114 1 : return false;
115 : }
116 :
117 9 : auto poInputDrv = poInputDS->GetDriver();
118 9 : if (poInputDrv && EQUAL(poInputDrv->GetDescription(), "MEM"))
119 : {
120 1 : ReportError(CE_Warning, CPLE_AppDefined,
121 : "Input dataset is a in-memory dataset. Skipping binary "
122 : "file comparison");
123 1 : return false;
124 : }
125 :
126 16 : VSIVirtualHandleUniquePtr fpRef(VSIFOpenL(poRefDS->GetDescription(), "rb"));
127 : VSIVirtualHandleUniquePtr fpInput(
128 16 : VSIFOpenL(poInputDS->GetDescription(), "rb"));
129 8 : if (!fpRef)
130 : {
131 1 : ReportError(CE_Warning, CPLE_AppDefined,
132 : "Reference dataset '%s' is not a file. Skipping binary "
133 : "file comparison",
134 1 : poRefDS->GetDescription());
135 1 : return false;
136 : }
137 :
138 7 : if (!fpInput)
139 : {
140 1 : ReportError(
141 : CE_Warning, CPLE_AppDefined,
142 : "Input dataset '%s' is not a file. Skipping binary file comparison",
143 1 : poInputDS->GetDescription());
144 1 : return false;
145 : }
146 :
147 6 : fpRef->Seek(0, SEEK_END);
148 6 : fpInput->Seek(0, SEEK_END);
149 6 : const auto nRefSize = fpRef->Tell();
150 6 : const auto nInputSize = fpInput->Tell();
151 6 : if (nRefSize != nInputSize)
152 : {
153 4 : aosReport.push_back("Reference file has size " +
154 6 : std::to_string(nRefSize) +
155 4 : " bytes, whereas input file has size " +
156 8 : std::to_string(nInputSize) + " bytes.");
157 :
158 2 : return false;
159 : }
160 :
161 4 : constexpr size_t BUF_SIZE = 1024 * 1024;
162 8 : std::vector<GByte> abyRef(BUF_SIZE);
163 8 : std::vector<GByte> abyInput(BUF_SIZE);
164 :
165 4 : fpRef->Seek(0, SEEK_SET);
166 4 : fpInput->Seek(0, SEEK_SET);
167 :
168 0 : do
169 : {
170 4 : const size_t nRefRead = fpRef->Read(abyRef.data(), 1, BUF_SIZE);
171 4 : const size_t nInputRead = fpInput->Read(abyInput.data(), 1, BUF_SIZE);
172 :
173 4 : if (nRefRead != BUF_SIZE && fpRef->Tell() != nRefSize)
174 : {
175 0 : aosReport.push_back("Failed to fully read reference file");
176 0 : return false;
177 : }
178 :
179 4 : if (nInputRead != BUF_SIZE && fpInput->Tell() != nRefSize)
180 : {
181 0 : aosReport.push_back("Failed to fully read input file");
182 0 : return false;
183 : }
184 :
185 4 : if (abyRef != abyInput)
186 : {
187 1 : aosReport.push_back(
188 : "Reference file and input file differ at the binary level.");
189 1 : return false;
190 : }
191 3 : } while (fpRef->Tell() < nRefSize);
192 :
193 3 : return true;
194 : }
195 :
196 : /************************************************************************/
197 : /* GDALRasterCompareAlgorithm::CRSComparison() */
198 : /************************************************************************/
199 :
200 175 : void GDALRasterCompareAlgorithm::CRSComparison(
201 : std::vector<std::string> &aosReport, GDALDataset *poRefDS,
202 : GDALDataset *poInputDS)
203 : {
204 175 : const auto poRefCRS = poRefDS->GetSpatialRef();
205 175 : const auto poInputCRS = poInputDS->GetSpatialRef();
206 :
207 175 : if (poRefCRS == nullptr)
208 : {
209 157 : if (poInputCRS)
210 : {
211 1 : aosReport.push_back(
212 : "Reference dataset has no CRS, but input dataset has one.");
213 : }
214 174 : return;
215 : }
216 :
217 18 : if (poInputCRS == nullptr)
218 : {
219 1 : aosReport.push_back(
220 : "Reference dataset has a CRS, but input dataset has none.");
221 1 : return;
222 : }
223 :
224 17 : if (poRefCRS->IsSame(poInputCRS))
225 16 : return;
226 :
227 1 : const char *apszOptions[] = {"FORMAT=WKT2_2019", nullptr};
228 2 : const auto poRefWKT = poRefCRS->exportToWkt(apszOptions);
229 1 : const auto poInputWKT = poInputCRS->exportToWkt(apszOptions);
230 1 : aosReport.push_back(
231 2 : "Reference and input CRS are not equivalent. Reference one is '" +
232 2 : poRefWKT + "'. Input one is '" + poInputWKT + "'");
233 : }
234 :
235 : /************************************************************************/
236 : /* GDALRasterCompareAlgorithm::GeotransformComparison() */
237 : /************************************************************************/
238 :
239 175 : void GDALRasterCompareAlgorithm::GeoTransformComparison(
240 : std::vector<std::string> &aosReport, GDALDataset *poRefDS,
241 : GDALDataset *poInputDS)
242 : {
243 175 : GDALGeoTransform refGT;
244 175 : CPLErr eErr1 = poRefDS->GetGeoTransform(refGT);
245 175 : GDALGeoTransform inputGT;
246 175 : CPLErr eErr2 = poInputDS->GetGeoTransform(inputGT);
247 175 : if (eErr1 == CE_Failure && eErr2 == CE_Failure)
248 154 : return;
249 :
250 24 : if (eErr1 == CE_Failure && eErr2 == CE_None)
251 : {
252 1 : aosReport.push_back(
253 : "Reference dataset has no geotransform, but input one has one.");
254 1 : return;
255 : }
256 :
257 23 : if (eErr1 == CE_None && eErr2 == CE_Failure)
258 : {
259 1 : aosReport.push_back(
260 : "Reference dataset has a geotransform, but input one has none.");
261 1 : return;
262 : }
263 :
264 153 : for (int i = 0; i < 6; ++i)
265 : {
266 132 : if ((refGT[i] != 0 &&
267 263 : std::fabs(refGT[i] - inputGT[i]) > 1e-10 * std::fabs(refGT[i])) ||
268 131 : (refGT[i] == 0 && std::fabs(refGT[i] - inputGT[i]) > 1e-10))
269 : {
270 : std::string s = "Geotransform of reference and input dataset are "
271 1 : "not equivalent. Reference geotransform is (";
272 7 : for (int j = 0; j < 6; ++j)
273 : {
274 6 : if (j > 0)
275 5 : s += ',';
276 6 : s += std::to_string(refGT[j]);
277 : }
278 1 : s += "). Input geotransform is (";
279 7 : for (int j = 0; j < 6; ++j)
280 : {
281 6 : if (j > 0)
282 5 : s += ',';
283 6 : s += std::to_string(inputGT[j]);
284 : }
285 1 : s += ')';
286 1 : aosReport.push_back(std::move(s));
287 1 : return;
288 : }
289 : }
290 : }
291 :
292 : /************************************************************************/
293 : /* Diff() */
294 : /************************************************************************/
295 :
296 16880 : template <class T> CPL_NOSANITIZE_UNSIGNED_INT_OVERFLOW static T Diff(T a, T b)
297 : {
298 16880 : return a - b;
299 : }
300 :
301 : /************************************************************************/
302 : /* CompareVectors() */
303 : /************************************************************************/
304 :
305 : template <class T, class Tdiff, bool bIsComplex>
306 1254 : static void CompareVectors(size_t nValCount, const T *refValues,
307 : const T *inputValues, uint64_t &countDiffPixels,
308 : Tdiff &maxDiffValue)
309 : {
310 1254 : constexpr bool bIsFloatingPoint =
311 : std::is_same_v<T, float> || std::is_same_v<T, double>;
312 : if constexpr (bIsComplex)
313 : {
314 374 : for (size_t i = 0; i < nValCount; ++i)
315 : {
316 : if constexpr (bIsFloatingPoint)
317 : {
318 : static_assert(std::is_same_v<T, Tdiff>);
319 121 : if (std::isnan(refValues[2 * i]) &&
320 4 : std::isnan(inputValues[2 * i]) &&
321 123 : std::isnan(refValues[2 * i + 1]) &&
322 2 : std::isnan(inputValues[2 * i + 1]))
323 : {
324 2 : continue;
325 : }
326 : }
327 :
328 185 : if (refValues[2 * i] != inputValues[2 * i] ||
329 175 : refValues[2 * i + 1] != inputValues[2 * i + 1])
330 : {
331 : const Tdiff diff =
332 10 : std::hypot(static_cast<Tdiff>(refValues[2 * i]) -
333 10 : static_cast<Tdiff>(inputValues[2 * i]),
334 10 : static_cast<Tdiff>(refValues[2 * i + 1]) -
335 10 : static_cast<Tdiff>(inputValues[2 * i + 1]));
336 10 : ++countDiffPixels;
337 10 : if (diff > maxDiffValue)
338 10 : maxDiffValue = diff;
339 : }
340 : }
341 : }
342 : else
343 : {
344 : static_assert(sizeof(Tdiff) == sizeof(T));
345 18277 : for (size_t i = 0; i < nValCount; ++i)
346 : {
347 : if constexpr (bIsFloatingPoint)
348 : {
349 : static_assert(std::is_same_v<T, Tdiff>);
350 342 : if (std::isnan(refValues[i]))
351 : {
352 8 : if (!std::isnan(inputValues[i]))
353 : {
354 2 : ++countDiffPixels;
355 : }
356 8 : continue;
357 : }
358 334 : else if (std::isnan(inputValues[i]))
359 : {
360 2 : ++countDiffPixels;
361 2 : continue;
362 : }
363 332 : else if (refValues[i] == inputValues[i])
364 : {
365 320 : continue;
366 : }
367 : }
368 :
369 33760 : const Tdiff diff = refValues[i] >= inputValues[i]
370 16880 : ? Diff(static_cast<Tdiff>(refValues[i]),
371 16858 : static_cast<Tdiff>(inputValues[i]))
372 22 : : Diff(static_cast<Tdiff>(inputValues[i]),
373 22 : static_cast<Tdiff>(refValues[i]));
374 16880 : if (diff > 0)
375 : {
376 45 : ++countDiffPixels;
377 45 : if (diff > maxDiffValue)
378 45 : maxDiffValue = diff;
379 : }
380 : }
381 : }
382 1254 : }
383 :
384 : /************************************************************************/
385 : /* DatasetPixelComparison() */
386 : /************************************************************************/
387 :
388 : template <class T, class Tdiff, bool bIsComplex>
389 61 : static void DatasetPixelComparison(std::vector<std::string> &aosReport,
390 : GDALDataset *poRefDS, GDALDataset *poInputDS,
391 : GDALDataType eReqDT,
392 : GDALProgressFunc pfnProgress,
393 : void *pProgressData)
394 : {
395 122 : std::vector<T> refValues;
396 122 : std::vector<T> inputValues;
397 :
398 61 : CPLAssert(GDALDataTypeIsComplex(eReqDT) == bIsComplex);
399 :
400 61 : const uint64_t nTotalPixels =
401 61 : static_cast<uint64_t>(poRefDS->GetRasterXSize()) *
402 61 : poRefDS->GetRasterYSize();
403 61 : uint64_t nIterPixels = 0;
404 :
405 61 : constexpr int nValPerPixel = bIsComplex ? 2 : 1;
406 61 : const int nBands = poRefDS->GetRasterCount();
407 :
408 122 : std::vector<Tdiff> maxDiffValue(nBands, 0);
409 122 : std::vector<uint64_t> countDiffPixels(nBands, 0);
410 :
411 181 : for (const auto &window : poRefDS->GetRasterBand(1)->IterateWindows())
412 : {
413 61 : const size_t nValCount =
414 61 : static_cast<size_t>(window.nXSize) * window.nYSize;
415 61 : refValues.resize(nValCount * nValPerPixel * nBands);
416 61 : inputValues.resize(nValCount * nValPerPixel * nBands);
417 :
418 122 : if (poRefDS->RasterIO(GF_Read, window.nXOff, window.nYOff,
419 61 : window.nXSize, window.nYSize, refValues.data(),
420 61 : window.nXSize, window.nYSize, eReqDT, nBands,
421 122 : nullptr, 0, 0, 0, nullptr) == CE_None &&
422 : poInputDS->RasterIO(
423 61 : GF_Read, window.nXOff, window.nYOff, window.nXSize,
424 122 : window.nYSize, inputValues.data(), window.nXSize, window.nYSize,
425 61 : eReqDT, nBands, nullptr, 0, 0, 0, nullptr) == CE_None)
426 : {
427 1037 : for (int i = 0; i < nBands; ++i)
428 : {
429 976 : CompareVectors<T, Tdiff, bIsComplex>(
430 976 : nValCount, refValues.data() + i * nValCount * nValPerPixel,
431 976 : inputValues.data() + i * nValCount * nValPerPixel,
432 976 : countDiffPixels[i], maxDiffValue[i]);
433 : }
434 : }
435 : else
436 : {
437 0 : aosReport.push_back("I/O error when comparing pixel values");
438 : }
439 :
440 61 : if (pfnProgress)
441 : {
442 2 : nIterPixels += nValCount;
443 2 : if (!pfnProgress(static_cast<double>(nIterPixels) /
444 : static_cast<double>(nTotalPixels),
445 : "", pProgressData))
446 : {
447 1 : CPLError(CE_Failure, CPLE_UserInterrupt, "Interrupted by user");
448 1 : break;
449 : }
450 : }
451 : }
452 1037 : for (int i = 0; i < nBands; ++i)
453 : {
454 976 : if (countDiffPixels[i])
455 : {
456 29 : aosReport.push_back(
457 : "Band " + std::to_string(i + 1) +
458 29 : ": pixels differing: " + std::to_string(countDiffPixels[i]));
459 29 : aosReport.push_back("Band " + std::to_string(i + 1) +
460 : ": maximum pixel value difference: " +
461 29 : std::to_string(maxDiffValue[i]));
462 : }
463 : }
464 61 : }
465 :
466 : /************************************************************************/
467 : /* GDALRasterCompareAlgorithm::DatasetComparison() */
468 : /************************************************************************/
469 :
470 184 : void GDALRasterCompareAlgorithm::DatasetComparison(
471 : std::vector<std::string> &aosReport, GDALDataset *poRefDS,
472 : GDALDataset *poInputDS, GDALProgressFunc pfnProgress, void *pProgressData)
473 : {
474 184 : if (!m_skipCRS)
475 : {
476 175 : CRSComparison(aosReport, poRefDS, poInputDS);
477 : }
478 :
479 184 : if (!m_skipGeotransform)
480 : {
481 175 : GeoTransformComparison(aosReport, poRefDS, poInputDS);
482 : }
483 :
484 184 : bool ret = true;
485 184 : if (poRefDS->GetRasterCount() != poInputDS->GetRasterCount())
486 : {
487 2 : aosReport.push_back("Reference dataset has " +
488 3 : std::to_string(poRefDS->GetRasterCount()) +
489 2 : " band(s), but input dataset has " +
490 2 : std::to_string(poInputDS->GetRasterCount()));
491 1 : ret = false;
492 : }
493 :
494 184 : if (poRefDS->GetRasterXSize() != poInputDS->GetRasterXSize())
495 : {
496 2 : aosReport.push_back("Reference dataset width is " +
497 3 : std::to_string(poRefDS->GetRasterXSize()) +
498 2 : ", but input dataset width is " +
499 2 : std::to_string(poInputDS->GetRasterXSize()));
500 1 : ret = false;
501 : }
502 :
503 184 : if (poRefDS->GetRasterYSize() != poInputDS->GetRasterYSize())
504 : {
505 2 : aosReport.push_back("Reference dataset height is " +
506 3 : std::to_string(poRefDS->GetRasterYSize()) +
507 2 : ", but input dataset height is " +
508 2 : std::to_string(poInputDS->GetRasterYSize()));
509 1 : ret = false;
510 : }
511 :
512 184 : if (!m_skipMetadata)
513 : {
514 175 : MetadataComparison(aosReport, "(dataset default metadata domain)",
515 175 : poRefDS->GetMetadata(), poInputDS->GetMetadata());
516 : }
517 :
518 184 : if (!m_skipRPC)
519 : {
520 175 : MetadataComparison(aosReport, "RPC", poRefDS->GetMetadata("RPC"),
521 175 : poInputDS->GetMetadata("RPC"));
522 : }
523 :
524 184 : if (!m_skipGeolocation)
525 : {
526 175 : MetadataComparison(aosReport, "GEOLOCATION",
527 175 : poRefDS->GetMetadata("GEOLOCATION"),
528 175 : poInputDS->GetMetadata("GEOLOCATION"));
529 : }
530 :
531 184 : if (!ret)
532 3 : return;
533 :
534 181 : const int nBands = poRefDS->GetRasterCount();
535 :
536 181 : bool doBandBasedPixelComparison = true;
537 : // Do not do band-by-band pixel difference if there are too many interleaved
538 : // bands as this could be extremely slow
539 181 : if (nBands > 10)
540 : {
541 : const char *pszRefInterleave =
542 61 : poRefDS->GetMetadataItem("INTERLEAVE", "IMAGE_STRUCTURE");
543 : const char *pszInputInterleave =
544 61 : poInputDS->GetMetadataItem("INTERLEAVE", "IMAGE_STRUCTURE");
545 61 : if ((pszRefInterleave && EQUAL(pszRefInterleave, "PIXEL")) ||
546 0 : (pszInputInterleave && EQUAL(pszInputInterleave, "PIXEL")))
547 : {
548 61 : doBandBasedPixelComparison = false;
549 : }
550 : }
551 :
552 1273 : for (int i = 0; i < nBands; ++i)
553 : {
554 2184 : void *pScaledProgress = GDALCreateScaledProgress(
555 1092 : static_cast<double>(i) / nBands,
556 1092 : static_cast<double>(i + 1) / nBands, pfnProgress, pProgressData);
557 2184 : BandComparison(
558 2184 : aosReport, std::to_string(i + 1), doBandBasedPixelComparison,
559 : poRefDS->GetRasterBand(i + 1), poInputDS->GetRasterBand(i + 1),
560 : pScaledProgress ? GDALScaledProgress : nullptr, pScaledProgress);
561 1092 : GDALDestroyScaledProgress(pScaledProgress);
562 : }
563 :
564 181 : if (!doBandBasedPixelComparison)
565 : {
566 : const auto eReqDT =
567 61 : GDALDataTypeUnion(poRefDS->GetRasterBand(1)->GetRasterDataType(),
568 : poInputDS->GetRasterBand(1)->GetRasterDataType());
569 61 : switch (eReqDT)
570 : {
571 5 : case GDT_UInt8:
572 5 : DatasetPixelComparison<uint8_t, uint8_t, false>(
573 : aosReport, poRefDS, poInputDS, eReqDT, pfnProgress,
574 : pProgressData);
575 5 : break;
576 4 : case GDT_Int8:
577 4 : DatasetPixelComparison<int8_t, uint8_t, false>(
578 : aosReport, poRefDS, poInputDS, eReqDT, pfnProgress,
579 : pProgressData);
580 4 : break;
581 3 : case GDT_UInt16:
582 3 : DatasetPixelComparison<uint16_t, uint16_t, false>(
583 : aosReport, poRefDS, poInputDS, eReqDT, pfnProgress,
584 : pProgressData);
585 3 : break;
586 4 : case GDT_Int16:
587 4 : DatasetPixelComparison<int16_t, uint16_t, false>(
588 : aosReport, poRefDS, poInputDS, eReqDT, pfnProgress,
589 : pProgressData);
590 4 : break;
591 3 : case GDT_UInt32:
592 3 : DatasetPixelComparison<uint32_t, uint32_t, false>(
593 : aosReport, poRefDS, poInputDS, eReqDT, pfnProgress,
594 : pProgressData);
595 3 : break;
596 4 : case GDT_Int32:
597 4 : DatasetPixelComparison<int32_t, uint32_t, false>(
598 : aosReport, poRefDS, poInputDS, eReqDT, pfnProgress,
599 : pProgressData);
600 4 : break;
601 3 : case GDT_UInt64:
602 3 : DatasetPixelComparison<uint64_t, uint64_t, false>(
603 : aosReport, poRefDS, poInputDS, eReqDT, pfnProgress,
604 : pProgressData);
605 3 : break;
606 4 : case GDT_Int64:
607 4 : DatasetPixelComparison<int64_t, uint64_t, false>(
608 : aosReport, poRefDS, poInputDS, eReqDT, pfnProgress,
609 : pProgressData);
610 4 : break;
611 14 : case GDT_Float16:
612 : case GDT_Float32:
613 14 : DatasetPixelComparison<float, float, false>(
614 : aosReport, poRefDS, poInputDS, GDT_Float32, pfnProgress,
615 : pProgressData);
616 14 : break;
617 6 : case GDT_Float64:
618 6 : DatasetPixelComparison<double, double, false>(
619 : aosReport, poRefDS, poInputDS, eReqDT, pfnProgress,
620 : pProgressData);
621 6 : break;
622 2 : case GDT_CInt16:
623 2 : DatasetPixelComparison<int16_t, float, true>(
624 : aosReport, poRefDS, poInputDS, eReqDT, pfnProgress,
625 : pProgressData);
626 2 : break;
627 2 : case GDT_CInt32:
628 2 : DatasetPixelComparison<int32_t, double, true>(
629 : aosReport, poRefDS, poInputDS, eReqDT, pfnProgress,
630 : pProgressData);
631 2 : break;
632 5 : case GDT_CFloat16:
633 : case GDT_CFloat32:
634 5 : DatasetPixelComparison<float, float, true>(
635 : aosReport, poRefDS, poInputDS, GDT_CFloat32, pfnProgress,
636 : pProgressData);
637 5 : break;
638 2 : case GDT_CFloat64:
639 2 : DatasetPixelComparison<double, double, true>(
640 : aosReport, poRefDS, poInputDS, eReqDT, pfnProgress,
641 : pProgressData);
642 2 : break;
643 0 : case GDT_Unknown:
644 : case GDT_TypeCount:
645 0 : break;
646 : }
647 : }
648 : }
649 :
650 : /************************************************************************/
651 : /* ComparePixels() */
652 : /************************************************************************/
653 :
654 : template <class T, class Tdiff, bool bIsComplex>
655 119 : static void ComparePixels(std::vector<std::string> &aosReport,
656 : const std::string &bandId, GDALRasterBand *poRefBand,
657 : GDALRasterBand *poInputBand, GDALDataType eReqDT,
658 : GDALProgressFunc pfnProgress, void *pProgressData)
659 : {
660 238 : std::vector<T> refValues;
661 238 : std::vector<T> inputValues;
662 119 : Tdiff maxDiffValue = 0;
663 119 : uint64_t countDiffPixels = 0;
664 :
665 119 : CPLAssert(GDALDataTypeIsComplex(eReqDT) == bIsComplex);
666 119 : const uint64_t nTotalPixels =
667 119 : static_cast<uint64_t>(poRefBand->GetXSize()) * poRefBand->GetYSize();
668 119 : uint64_t nIterPixels = 0;
669 :
670 119 : constexpr int nValPerPixel = bIsComplex ? 2 : 1;
671 :
672 673 : for (const auto &window : poRefBand->IterateWindows())
673 : {
674 278 : const size_t nValCount =
675 278 : static_cast<size_t>(window.nXSize) * window.nYSize;
676 278 : refValues.resize(nValCount * nValPerPixel);
677 278 : inputValues.resize(nValCount * nValPerPixel);
678 :
679 556 : if (poRefBand->RasterIO(GF_Read, window.nXOff, window.nYOff,
680 278 : window.nXSize, window.nYSize, refValues.data(),
681 278 : window.nXSize, window.nYSize, eReqDT, 0, 0,
682 556 : nullptr) == CE_None &&
683 : poInputBand->RasterIO(
684 278 : GF_Read, window.nXOff, window.nYOff, window.nXSize,
685 556 : window.nYSize, inputValues.data(), window.nXSize, window.nYSize,
686 278 : eReqDT, 0, 0, nullptr) == CE_None)
687 : {
688 278 : CompareVectors<T, Tdiff, bIsComplex>(nValCount, refValues.data(),
689 : inputValues.data(),
690 : countDiffPixels, maxDiffValue);
691 : }
692 : else
693 : {
694 0 : aosReport.push_back("I/O error when comparing pixel values");
695 : }
696 :
697 278 : if (pfnProgress)
698 : {
699 6 : nIterPixels += nValCount;
700 6 : if (!pfnProgress(static_cast<double>(nIterPixels) /
701 : static_cast<double>(nTotalPixels),
702 : "", pProgressData))
703 : {
704 1 : CPLError(CE_Failure, CPLE_UserInterrupt, "Interrupted by user");
705 1 : break;
706 : }
707 : }
708 : }
709 119 : if (countDiffPixels)
710 : {
711 30 : aosReport.push_back(
712 : bandId + ": pixels differing: " + std::to_string(countDiffPixels));
713 30 : aosReport.push_back(bandId + ": maximum pixel value difference: " +
714 : std::to_string(maxDiffValue));
715 : }
716 119 : }
717 :
718 : /************************************************************************/
719 : /* ComparePixels() */
720 : /************************************************************************/
721 :
722 119 : static void ComparePixels(std::vector<std::string> &aosReport,
723 : const std::string &bandId, GDALRasterBand *poRefBand,
724 : GDALRasterBand *poInputBand,
725 : GDALProgressFunc pfnProgress, void *pProgressData)
726 : {
727 119 : const auto eReqDT = GDALDataTypeUnion(poRefBand->GetRasterDataType(),
728 : poInputBand->GetRasterDataType());
729 119 : switch (eReqDT)
730 : {
731 59 : case GDT_UInt8:
732 59 : ComparePixels<uint8_t, uint8_t, false>(aosReport, bandId, poRefBand,
733 : poInputBand, eReqDT,
734 : pfnProgress, pProgressData);
735 59 : break;
736 4 : case GDT_Int8:
737 4 : ComparePixels<int8_t, uint8_t, false>(aosReport, bandId, poRefBand,
738 : poInputBand, eReqDT,
739 : pfnProgress, pProgressData);
740 4 : break;
741 4 : case GDT_UInt16:
742 4 : ComparePixels<uint16_t, uint16_t, false>(
743 : aosReport, bandId, poRefBand, poInputBand, eReqDT, pfnProgress,
744 : pProgressData);
745 4 : break;
746 5 : case GDT_Int16:
747 5 : ComparePixels<int16_t, uint16_t, false>(
748 : aosReport, bandId, poRefBand, poInputBand, eReqDT, pfnProgress,
749 : pProgressData);
750 5 : break;
751 3 : case GDT_UInt32:
752 3 : ComparePixels<uint32_t, uint32_t, false>(
753 : aosReport, bandId, poRefBand, poInputBand, eReqDT, pfnProgress,
754 : pProgressData);
755 3 : break;
756 4 : case GDT_Int32:
757 4 : ComparePixels<int32_t, uint32_t, false>(
758 : aosReport, bandId, poRefBand, poInputBand, eReqDT, pfnProgress,
759 : pProgressData);
760 4 : break;
761 3 : case GDT_UInt64:
762 3 : ComparePixels<uint64_t, uint64_t, false>(
763 : aosReport, bandId, poRefBand, poInputBand, eReqDT, pfnProgress,
764 : pProgressData);
765 3 : break;
766 4 : case GDT_Int64:
767 4 : ComparePixels<int64_t, uint64_t, false>(
768 : aosReport, bandId, poRefBand, poInputBand, eReqDT, pfnProgress,
769 : pProgressData);
770 4 : break;
771 16 : case GDT_Float16:
772 : case GDT_Float32:
773 16 : ComparePixels<float, float, false>(aosReport, bandId, poRefBand,
774 : poInputBand, GDT_Float32,
775 : pfnProgress, pProgressData);
776 16 : break;
777 6 : case GDT_Float64:
778 6 : ComparePixels<double, double, false>(aosReport, bandId, poRefBand,
779 : poInputBand, eReqDT,
780 : pfnProgress, pProgressData);
781 6 : break;
782 2 : case GDT_CInt16:
783 2 : ComparePixels<int16_t, float, true>(aosReport, bandId, poRefBand,
784 : poInputBand, eReqDT,
785 : pfnProgress, pProgressData);
786 2 : break;
787 2 : case GDT_CInt32:
788 2 : ComparePixels<int32_t, double, true>(aosReport, bandId, poRefBand,
789 : poInputBand, eReqDT,
790 : pfnProgress, pProgressData);
791 2 : break;
792 5 : case GDT_CFloat16:
793 : case GDT_CFloat32:
794 5 : ComparePixels<float, float, true>(aosReport, bandId, poRefBand,
795 : poInputBand, GDT_CFloat32,
796 : pfnProgress, pProgressData);
797 5 : break;
798 2 : case GDT_CFloat64:
799 2 : ComparePixels<double, double, true>(aosReport, bandId, poRefBand,
800 : poInputBand, eReqDT,
801 : pfnProgress, pProgressData);
802 2 : break;
803 0 : case GDT_Unknown:
804 : case GDT_TypeCount:
805 0 : break;
806 : }
807 119 : }
808 :
809 : /************************************************************************/
810 : /* GDALRasterCompareAlgorithm::BandComparison() */
811 : /************************************************************************/
812 :
813 1096 : void GDALRasterCompareAlgorithm::BandComparison(
814 : std::vector<std::string> &aosReport, const std::string &bandId,
815 : bool doBandBasedPixelComparison, GDALRasterBand *poRefBand,
816 : GDALRasterBand *poInputBand, GDALProgressFunc pfnProgress,
817 : void *pProgressData)
818 : {
819 1096 : bool ret = true;
820 :
821 1096 : if (poRefBand->GetXSize() != poInputBand->GetXSize())
822 : {
823 2 : aosReport.push_back("Reference band width is " +
824 3 : std::to_string(poRefBand->GetXSize()) +
825 2 : ", but input band width is " +
826 2 : std::to_string(poInputBand->GetXSize()));
827 1 : ret = false;
828 : }
829 :
830 1096 : if (poRefBand->GetYSize() != poInputBand->GetYSize())
831 : {
832 2 : aosReport.push_back("Reference band height is " +
833 3 : std::to_string(poRefBand->GetYSize()) +
834 2 : ", but input band height is " +
835 2 : std::to_string(poInputBand->GetYSize()));
836 1 : ret = false;
837 : }
838 :
839 1096 : if (strcmp(poRefBand->GetDescription(), poInputBand->GetDescription()) != 0)
840 : {
841 3 : aosReport.push_back("Reference band " + bandId + " has description " +
842 4 : std::string(poRefBand->GetDescription()) +
843 2 : ", but input band has description " +
844 2 : std::string(poInputBand->GetDescription()));
845 : }
846 :
847 1096 : if (poRefBand->GetRasterDataType() != poInputBand->GetRasterDataType())
848 : {
849 2 : aosReport.push_back(
850 4 : "Reference band " + bandId + " has data type " +
851 8 : std::string(GDALGetDataTypeName(poRefBand->GetRasterDataType())) +
852 4 : ", but input band has data type " +
853 4 : std::string(GDALGetDataTypeName(poInputBand->GetRasterDataType())));
854 : }
855 :
856 1096 : int bRefHasNoData = false;
857 1096 : const double dfRefNoData = poRefBand->GetNoDataValue(&bRefHasNoData);
858 1096 : int bInputHasNoData = false;
859 1096 : const double dfInputNoData = poInputBand->GetNoDataValue(&bInputHasNoData);
860 1096 : if (!bRefHasNoData && !bInputHasNoData)
861 : {
862 : // ok
863 : }
864 6 : else if (bRefHasNoData && !bInputHasNoData)
865 : {
866 3 : aosReport.push_back("Reference band " + bandId + " has nodata value " +
867 4 : std::to_string(dfRefNoData) +
868 : ", but input band has none.");
869 : }
870 5 : else if (!bRefHasNoData && bInputHasNoData)
871 : {
872 3 : aosReport.push_back("Reference band " + bandId +
873 2 : " has no nodata value, " +
874 2 : "but input band has no data value " +
875 4 : std::to_string(dfInputNoData) + ".");
876 : }
877 4 : else if ((std::isnan(dfRefNoData) && std::isnan(dfInputNoData)) ||
878 : dfRefNoData == dfInputNoData)
879 : {
880 : // ok
881 : }
882 : else
883 : {
884 6 : aosReport.push_back("Reference band " + bandId + " has nodata value " +
885 8 : std::to_string(dfRefNoData) +
886 4 : ", but input band has no data value " +
887 8 : std::to_string(dfInputNoData) + ".");
888 : }
889 :
890 1096 : if (poRefBand->GetColorInterpretation() !=
891 1096 : poInputBand->GetColorInterpretation())
892 : {
893 3 : aosReport.push_back("Reference band " + bandId +
894 2 : " has color interpretation " +
895 3 : std::string(GDALGetColorInterpretationName(
896 3 : poRefBand->GetColorInterpretation())) +
897 2 : ", but input band has color interpretation " +
898 3 : std::string(GDALGetColorInterpretationName(
899 1 : poInputBand->GetColorInterpretation())));
900 : }
901 :
902 1096 : if (!ret)
903 1 : return;
904 :
905 : const uint64_t nBasePixels =
906 1095 : static_cast<uint64_t>(poRefBand->GetXSize()) * poRefBand->GetYSize();
907 1095 : uint64_t nTotalPixels = nBasePixels;
908 1095 : const int nOvrCount = poRefBand->GetOverviewCount();
909 1097 : for (int i = 0; i < nOvrCount; ++i)
910 : {
911 2 : auto poOvrBand = poRefBand->GetOverview(i);
912 : const uint64_t nOvrPixels =
913 2 : static_cast<uint64_t>(poOvrBand->GetXSize()) *
914 2 : poOvrBand->GetYSize();
915 2 : nTotalPixels += nOvrPixels;
916 : }
917 :
918 1095 : if (doBandBasedPixelComparison)
919 : {
920 : void *pScaledProgress =
921 238 : GDALCreateScaledProgress(0.0,
922 119 : static_cast<double>(nBasePixels) /
923 119 : static_cast<double>(nTotalPixels),
924 : pfnProgress, pProgressData);
925 119 : ComparePixels(aosReport, bandId, poRefBand, poInputBand,
926 : pScaledProgress ? GDALScaledProgress : nullptr,
927 : pScaledProgress);
928 119 : GDALDestroyScaledProgress(pScaledProgress);
929 : }
930 :
931 1095 : if (!m_skipOverview)
932 : {
933 1087 : if (nOvrCount != poInputBand->GetOverviewCount())
934 : {
935 1 : aosReport.push_back(
936 2 : "Reference band " + bandId + " has " +
937 4 : std::to_string(nOvrCount) +
938 2 : " overview band(s), but input band has " +
939 2 : std::to_string(poInputBand->GetOverviewCount()));
940 : }
941 : else
942 : {
943 1086 : uint64_t nIterPixels = nBasePixels;
944 :
945 1088 : for (int i = 0; i < nOvrCount; ++i)
946 : {
947 2 : GDALRasterBand *poOvrBand = poRefBand->GetOverview(i);
948 : const uint64_t nOvrPixels =
949 2 : static_cast<uint64_t>(poOvrBand->GetXSize()) *
950 2 : poOvrBand->GetYSize();
951 4 : void *pScaledProgress = GDALCreateScaledProgress(
952 2 : static_cast<double>(nIterPixels) /
953 2 : static_cast<double>(nTotalPixels),
954 2 : static_cast<double>(nIterPixels + nOvrPixels) /
955 2 : static_cast<double>(nTotalPixels),
956 : pfnProgress, pProgressData);
957 4 : BandComparison(aosReport, "overview of band " + bandId,
958 : doBandBasedPixelComparison, poOvrBand,
959 2 : poInputBand->GetOverview(i),
960 : pScaledProgress ? GDALScaledProgress : nullptr,
961 : pScaledProgress);
962 2 : GDALDestroyScaledProgress(pScaledProgress);
963 2 : nIterPixels += nOvrPixels;
964 : }
965 : }
966 : }
967 :
968 1095 : if (poRefBand->GetMaskFlags() != poInputBand->GetMaskFlags())
969 : {
970 6 : aosReport.push_back("Reference band " + bandId + " has mask flags = " +
971 8 : std::to_string(poRefBand->GetMaskFlags()) +
972 4 : " , but input band has mask flags = " +
973 4 : std::to_string(poInputBand->GetMaskFlags()));
974 : }
975 1093 : else if (poRefBand->GetMaskFlags() == GMF_PER_DATASET)
976 : {
977 2 : BandComparison(aosReport, "mask of band " + bandId, true,
978 2 : poRefBand->GetMaskBand(), poInputBand->GetMaskBand(),
979 : nullptr, nullptr);
980 : }
981 :
982 1095 : if (!m_skipMetadata)
983 : {
984 1087 : MetadataComparison(aosReport, "(band default metadata domain)",
985 1087 : poRefBand->GetMetadata(),
986 1087 : poInputBand->GetMetadata());
987 : }
988 : }
989 :
990 : /************************************************************************/
991 : /* GDALRasterCompareAlgorithm::MetadataComparison() */
992 : /************************************************************************/
993 :
994 1612 : void GDALRasterCompareAlgorithm::MetadataComparison(
995 : std::vector<std::string> &aosReport, const std::string &metadataDomain,
996 : CSLConstList aosRef, CSLConstList aosInput)
997 : {
998 3224 : std::map<std::string, std::string> oMapRef;
999 3224 : std::map<std::string, std::string> oMapInput;
1000 :
1001 1612 : std::array<const char *, 3> ignoredKeys = {
1002 : "backend", // from gdalcompare.py. Not sure why
1003 : "ERR_BIAS", // RPC optional key
1004 : "ERR_RAND", // RPC optional key
1005 : };
1006 :
1007 1641 : for (const auto &[key, value] : cpl::IterateNameValue(aosRef))
1008 : {
1009 29 : const char *pszKey = key;
1010 87 : const auto eq = [pszKey](const char *s)
1011 87 : { return strcmp(pszKey, s) == 0; };
1012 29 : auto it = std::find_if(ignoredKeys.begin(), ignoredKeys.end(), eq);
1013 29 : if (it == ignoredKeys.end())
1014 : {
1015 28 : oMapRef[key] = value;
1016 : }
1017 : }
1018 :
1019 1641 : for (const auto &[key, value] : cpl::IterateNameValue(aosInput))
1020 : {
1021 29 : const char *pszKey = key;
1022 86 : const auto eq = [pszKey](const char *s)
1023 86 : { return strcmp(pszKey, s) == 0; };
1024 29 : auto it = std::find_if(ignoredKeys.begin(), ignoredKeys.end(), eq);
1025 29 : if (it == ignoredKeys.end())
1026 : {
1027 28 : oMapInput[key] = value;
1028 : }
1029 : }
1030 :
1031 4 : const auto strip = [](const std::string &s)
1032 : {
1033 4 : const auto posBegin = s.find_first_not_of(' ');
1034 4 : if (posBegin == std::string::npos)
1035 0 : return std::string();
1036 4 : const auto posEnd = s.find_last_not_of(' ');
1037 4 : return s.substr(posBegin, posEnd - posBegin + 1);
1038 : };
1039 :
1040 1640 : for (const auto &sKeyValuePair : oMapRef)
1041 : {
1042 28 : const auto oIter = oMapInput.find(sKeyValuePair.first);
1043 28 : if (oIter == oMapInput.end())
1044 : {
1045 6 : aosReport.push_back("Reference metadata " + metadataDomain +
1046 6 : " contains key '" + sKeyValuePair.first +
1047 : "' but input metadata does not.");
1048 : }
1049 : else
1050 : {
1051 : // this will always have the current date set
1052 26 : if (sKeyValuePair.first == "NITF_FDT")
1053 1 : continue;
1054 :
1055 50 : std::string ref = oIter->second;
1056 50 : std::string input = sKeyValuePair.second;
1057 25 : if (metadataDomain == "RPC")
1058 : {
1059 : // _RPC.TXT files and in-file have a difference
1060 : // in white space that is not otherwise meaningful.
1061 2 : ref = strip(ref);
1062 2 : input = strip(input);
1063 : }
1064 25 : if (ref != input)
1065 : {
1066 2 : aosReport.push_back(
1067 4 : "Reference metadata " + metadataDomain + " has value '" +
1068 6 : ref + "' for key '" + sKeyValuePair.first +
1069 4 : "' but input metadata has value '" + input + "'.");
1070 : }
1071 : }
1072 : }
1073 :
1074 1640 : for (const auto &sKeyValuePair : oMapInput)
1075 : {
1076 28 : if (!cpl::contains(oMapRef, sKeyValuePair.first))
1077 : {
1078 6 : aosReport.push_back("Input metadata " + metadataDomain +
1079 6 : " contains key '" + sKeyValuePair.first +
1080 : "' but reference metadata does not.");
1081 : }
1082 : }
1083 1612 : }
1084 :
1085 : /************************************************************************/
1086 : /* GDALRasterCompareAlgorithm::RunStep() */
1087 : /************************************************************************/
1088 :
1089 184 : bool GDALRasterCompareAlgorithm::RunStep(GDALPipelineStepRunContext &ctxt)
1090 : {
1091 184 : auto poRefDS = m_referenceDataset.GetDatasetRef();
1092 184 : CPLAssert(poRefDS);
1093 :
1094 184 : CPLAssert(m_inputDataset.size() == 1);
1095 184 : auto poInputDS = m_inputDataset[0].GetDatasetRef();
1096 184 : CPLAssert(poInputDS);
1097 :
1098 184 : if (m_skipAllOptional)
1099 : {
1100 8 : m_skipBinary = true;
1101 8 : m_skipCRS = true;
1102 8 : m_skipGeotransform = true;
1103 8 : m_skipOverview = true;
1104 8 : m_skipMetadata = true;
1105 8 : m_skipRPC = true;
1106 8 : m_skipGeolocation = true;
1107 8 : m_skipSubdataset = true;
1108 : }
1109 :
1110 368 : std::vector<std::string> aosReport;
1111 :
1112 184 : if (!m_skipBinary)
1113 : {
1114 12 : if (BinaryComparison(aosReport, poRefDS, poInputDS))
1115 : {
1116 3 : return true;
1117 : }
1118 : }
1119 :
1120 : CSLConstList papszSubDSRef =
1121 181 : m_skipSubdataset ? nullptr : poRefDS->GetMetadata("SUBDATASETS");
1122 181 : const int nCountRef = CSLCount(papszSubDSRef) / 2;
1123 : CSLConstList papszSubDSInput =
1124 181 : m_skipSubdataset ? nullptr : poInputDS->GetMetadata("SUBDATASETS");
1125 181 : const int nCountInput = CSLCount(papszSubDSInput) / 2;
1126 :
1127 181 : if (!m_skipSubdataset)
1128 : {
1129 172 : if (nCountRef != nCountInput)
1130 : {
1131 2 : aosReport.push_back("Reference dataset has " +
1132 3 : std::to_string(nCountRef) +
1133 2 : " subdataset(s) whereas input dataset has " +
1134 4 : std::to_string(nCountInput) + " one(s).");
1135 1 : m_skipSubdataset = true;
1136 : }
1137 : }
1138 :
1139 : // Compute total number of pixels, including in subdatasets
1140 : const uint64_t nBasePixels =
1141 181 : static_cast<uint64_t>(poRefDS->GetRasterXSize()) *
1142 181 : poRefDS->GetRasterYSize() * poRefDS->GetRasterCount();
1143 181 : uint64_t nTotalPixels = nBasePixels;
1144 181 : if (ctxt.m_pfnProgress && !m_skipSubdataset)
1145 : {
1146 7 : for (int i = 0; i < nCountRef; ++i)
1147 : {
1148 1 : const char *pszRef = CSLFetchNameValue(
1149 : papszSubDSRef, CPLSPrintf("SUBDATASET_%d_NAME", i + 1));
1150 1 : const char *pszInput = CSLFetchNameValue(
1151 : papszSubDSInput, CPLSPrintf("SUBDATASET_%d_NAME", i + 1));
1152 1 : if (pszRef && pszInput)
1153 : {
1154 : auto poSubRef = std::unique_ptr<GDALDataset>(
1155 2 : GDALDataset::Open(pszRef, GDAL_OF_RASTER));
1156 : auto poSubInput = std::unique_ptr<GDALDataset>(
1157 2 : GDALDataset::Open(pszInput, GDAL_OF_RASTER));
1158 1 : if (poSubRef && poSubInput)
1159 : {
1160 : const uint64_t nSubDSPixels =
1161 1 : static_cast<uint64_t>(poSubRef->GetRasterXSize()) *
1162 1 : poSubRef->GetRasterYSize() * poSubRef->GetRasterCount();
1163 1 : nTotalPixels += nSubDSPixels;
1164 : }
1165 : }
1166 : }
1167 : }
1168 :
1169 : {
1170 : void *pScaledProgress =
1171 362 : GDALCreateScaledProgress(0.0,
1172 181 : static_cast<double>(nBasePixels) /
1173 181 : static_cast<double>(nTotalPixels),
1174 : ctxt.m_pfnProgress, ctxt.m_pProgressData);
1175 181 : DatasetComparison(aosReport, poRefDS, poInputDS,
1176 : pScaledProgress ? GDALScaledProgress : nullptr,
1177 : pScaledProgress);
1178 181 : GDALDestroyScaledProgress(pScaledProgress);
1179 : }
1180 :
1181 181 : if (!m_skipSubdataset)
1182 : {
1183 171 : uint64_t nIterPixels = nBasePixels;
1184 174 : for (int i = 0; i < nCountRef; ++i)
1185 : {
1186 3 : const char *pszRef = CSLFetchNameValue(
1187 : papszSubDSRef, CPLSPrintf("SUBDATASET_%d_NAME", i + 1));
1188 3 : const char *pszInput = CSLFetchNameValue(
1189 : papszSubDSInput, CPLSPrintf("SUBDATASET_%d_NAME", i + 1));
1190 3 : if (pszRef && pszInput)
1191 : {
1192 : auto poSubRef = std::unique_ptr<GDALDataset>(GDALDataset::Open(
1193 6 : pszRef, GDAL_OF_RASTER | GDAL_OF_VERBOSE_ERROR));
1194 : auto poSubInput =
1195 : std::unique_ptr<GDALDataset>(GDALDataset::Open(
1196 6 : pszInput, GDAL_OF_RASTER | GDAL_OF_VERBOSE_ERROR));
1197 3 : if (poSubRef && poSubInput)
1198 : {
1199 : const uint64_t nSubDSPixels =
1200 3 : static_cast<uint64_t>(poSubRef->GetRasterXSize()) *
1201 3 : poSubRef->GetRasterYSize() * poSubRef->GetRasterCount();
1202 6 : void *pScaledProgress = GDALCreateScaledProgress(
1203 3 : static_cast<double>(nIterPixels) /
1204 3 : static_cast<double>(nTotalPixels),
1205 3 : static_cast<double>(nIterPixels + nSubDSPixels) /
1206 3 : static_cast<double>(nTotalPixels),
1207 : ctxt.m_pfnProgress, ctxt.m_pProgressData);
1208 3 : DatasetComparison(
1209 : aosReport, poSubRef.get(), poSubInput.get(),
1210 : pScaledProgress ? GDALScaledProgress : nullptr,
1211 : pScaledProgress);
1212 3 : GDALDestroyScaledProgress(pScaledProgress);
1213 3 : nIterPixels += nSubDSPixels;
1214 : }
1215 : }
1216 : }
1217 : }
1218 :
1219 331 : for (const auto &s : aosReport)
1220 : {
1221 150 : m_output += s;
1222 150 : m_output += '\n';
1223 : }
1224 :
1225 181 : m_retCode = static_cast<int>(aosReport.size());
1226 :
1227 181 : return true;
1228 : }
1229 :
1230 : /************************************************************************/
1231 : /* ~GDALRasterCompareAlgorithmStandalone() */
1232 : /************************************************************************/
1233 :
1234 : GDALRasterCompareAlgorithmStandalone::~GDALRasterCompareAlgorithmStandalone() =
1235 : default;
1236 :
1237 : //! @endcond
|
