From e239d4e975ef4ae1161ae8a41a27874d780824a5 Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Thu, 5 Oct 2023 13:02:37 +0200 Subject: [PATCH] Add a GDALPamDataset::SetDerivedDatasetName() method, and use it to be able to save statistics of datasets returned by GDALMDArray::AsClassicDataset() --- autotest/gdrivers/netcdf_multidim.py | 70 ++++++++++++++++++++ gcore/gdal_pam.h | 2 + gcore/gdalmultidim.cpp | 16 ++++- gcore/gdalpamdataset.cpp | 97 ++++++++++++++++++++++------ 4 files changed, 161 insertions(+), 24 deletions(-) diff --git a/autotest/gdrivers/netcdf_multidim.py b/autotest/gdrivers/netcdf_multidim.py index 06ac484ca6a0..49e390ea947d 100755 --- a/autotest/gdrivers/netcdf_multidim.py +++ b/autotest/gdrivers/netcdf_multidim.py @@ -3740,3 +3740,73 @@ def test_netcdf_multidim_getresampled_with_geoloc_EMIT(): )[0] == -10 ) + + +############################################################################### + + +@gdaltest.enable_exceptions() +def test_netcdf_multidim_serialize_statistics_asclassicdataset(tmp_path): + + filename = str( + tmp_path / "test_netcdf_multidim_serialize_statistics_asclassicdataset.nc" + ) + shutil.copy("data/netcdf/byte_no_cf.nc", filename) + + def test(): + ds = gdal.OpenEx(filename, gdal.OF_MULTIDIM_RASTER | gdal.OF_UPDATE) + rg = ds.GetRootGroup() + ar = rg.OpenMDArray("Band1") + + view = ar.GetView("[0:10,...]") + classic_ds = view.AsClassicDataset(1, 0) + assert classic_ds.GetRasterBand(1).GetStatistics(False, False) == [ + 0.0, + 0.0, + 0.0, + -1.0, + ] + classic_ds.GetRasterBand(1).ComputeStatistics(False) + + view = ar.GetView("[10:20,...]") + classic_ds = view.AsClassicDataset(1, 0) + classic_ds.GetRasterBand(1).ComputeStatistics(False) + + def reopen(): + + aux_xml = open(filename + ".aux.xml", "rb").read().decode("UTF-8") + assert ( + '' + in aux_xml + ) + assert ( + '' + in aux_xml + ) + + ds = gdal.OpenEx(filename, gdal.OF_MULTIDIM_RASTER) + rg = ds.GetRootGroup() + ar = rg.OpenMDArray("Band1") + + view = ar.GetView("[0:10,...]") + classic_ds = view.AsClassicDataset(1, 0) + assert classic_ds.GetRasterBand(1).GetStatistics(False, False) == pytest.approx( + [74.0, 255.0, 126.82, 26.729713803182] + ) + + view = ar.GetView("[10:20,...]") + classic_ds = view.AsClassicDataset(1, 0) + assert classic_ds.GetRasterBand(1).GetStatistics(False, False) == pytest.approx( + [99.0, 206.0, 126.71, 18.356086184152] + ) + + classic_ds = ar.AsClassicDataset(1, 0) + assert classic_ds.GetRasterBand(1).GetStatistics(False, False) == [ + 0.0, + 0.0, + 0.0, + -1.0, + ] + + test() + reopen() diff --git a/gcore/gdal_pam.h b/gcore/gdal_pam.h index 193a434b14cd..5b23e3d0e2f1 100644 --- a/gcore/gdal_pam.h +++ b/gcore/gdal_pam.h @@ -103,6 +103,7 @@ class GDALDatasetPamInfo CPLString osPhysicalFilename{}; CPLString osSubdatasetName{}; + CPLString osDerivedDatasetName{}; CPLString osAuxFilename{}; int bHasMetadata = false; @@ -145,6 +146,7 @@ class CPL_DLL GDALPamDataset : public GDALDataset const char *GetPhysicalFilename(); void SetSubdatasetName(const char *); const char *GetSubdatasetName(); + void SetDerivedDatasetName(const char *); //! @endcond public: diff --git a/gcore/gdalmultidim.cpp b/gcore/gdalmultidim.cpp index d2fe0e48cc7e..abb5583d9b23 100644 --- a/gcore/gdalmultidim.cpp +++ b/gcore/gdalmultidim.cpp @@ -8249,7 +8249,7 @@ struct MetadataItem }; } // namespace -class GDALRasterBandFromArray final : public GDALRasterBand +class GDALRasterBandFromArray final : public GDALPamRasterBand { std::vector m_anOffset{}; std::vector m_anCount{}; @@ -8281,7 +8281,7 @@ class GDALRasterBandFromArray final : public GDALRasterBand GDALColorInterp GetColorInterpretation() override; }; -class GDALDatasetFromArray final : public GDALDataset +class GDALDatasetFromArray final : public GDALPamDataset { friend class GDALRasterBandFromArray; @@ -8315,7 +8315,8 @@ class GDALDatasetFromArray final : public GDALDataset CPLErr eErr = CE_None; if (nOpenFlags != OPEN_FLAGS_CLOSED) { - if (GDALDatasetFromArray::FlushCache() != CE_None) + if (GDALDatasetFromArray::FlushCache(/*bAtClosing=*/true) != + CE_None) eErr = CE_Failure; m_poArray.reset(); } @@ -9138,6 +9139,15 @@ GDALDatasetFromArray *GDALDatasetFromArray::Create( if (iDim > 0) goto lbl_return_to_caller; + if (!array->GetFilename().empty()) + { + poDS->SetPhysicalFilename(array->GetFilename().c_str()); + poDS->SetDerivedDatasetName( + CPLSPrintf("AsClassicDataset(%d,%d) view of %s", int(iXDim), + int(iYDim), array->GetFullName().c_str())); + poDS->TryLoadXML(); + } + return poDS.release(); } diff --git a/gcore/gdalpamdataset.cpp b/gcore/gdalpamdataset.cpp index 0651f074d34d..49fc78caa575 100644 --- a/gcore/gdalpamdataset.cpp +++ b/gcore/gdalpamdataset.cpp @@ -137,6 +137,18 @@ * * poDS->TryLoadXML(); * \endcode + * + * In some situations where a derived dataset (e.g. used by + * GDALMDArray::AsClassicDataset()) is linked to a physical file, the name of + * the derived dataset is set with the SetDerivedSubdatasetName() method. + * + * \code + * poDS->SetDescription( poOpenInfo->pszFilename ); + * poDS->SetPhysicalFilename( poDS->pszFilename ); + * poDS->SetDerivedDatasetName( osDerivedDatasetName ); + * + * poDS->TryLoadXML(); + * \endcode */ class GDALPamDataset; @@ -761,6 +773,7 @@ const char *GDALPamDataset::GetPhysicalFilename() /* SetSubdatasetName() */ /************************************************************************/ +/* Mutually exclusive with SetDerivedDatasetName() */ void GDALPamDataset::SetSubdatasetName(const char *pszSubdataset) { @@ -770,6 +783,20 @@ void GDALPamDataset::SetSubdatasetName(const char *pszSubdataset) psPam->osSubdatasetName = pszSubdataset; } +/************************************************************************/ +/* SetDerivedDatasetName() */ +/************************************************************************/ + +/* Mutually exclusive with SetSubdatasetName() */ +void GDALPamDataset::SetDerivedDatasetName(const char *pszDerivedDataset) + +{ + PamInitialize(); + + if (psPam) + psPam->osDerivedDatasetName = pszDerivedDataset; +} + /************************************************************************/ /* GetSubdatasetName() */ /************************************************************************/ @@ -916,29 +943,44 @@ CPLErr GDALPamDataset::TryLoadXML(char **papszSiblingFiles) /* -------------------------------------------------------------------- */ /* If we are looking for a subdataset, search for its subtree now. */ /* -------------------------------------------------------------------- */ - if (psTree && !psPam->osSubdatasetName.empty()) + if (psTree) { - CPLXMLNode *psSubTree = psTree->psChild; - - for (; psSubTree != nullptr; psSubTree = psSubTree->psNext) + std::string osSubNode; + std::string osSubNodeValue; + if (!psPam->osSubdatasetName.empty()) { - if (psSubTree->eType != CXT_Element || - !EQUAL(psSubTree->pszValue, "Subdataset")) - continue; + osSubNode = "Subdataset"; + osSubNodeValue = psPam->osSubdatasetName; + } + else if (!psPam->osDerivedDatasetName.empty()) + { + osSubNode = "DerivedDataset"; + osSubNodeValue = psPam->osDerivedDatasetName; + } + if (!osSubNode.empty()) + { + CPLXMLNode *psSubTree = psTree->psChild; - if (!EQUAL(CPLGetXMLValue(psSubTree, "name", ""), - psPam->osSubdatasetName)) - continue; + for (; psSubTree != nullptr; psSubTree = psSubTree->psNext) + { + if (psSubTree->eType != CXT_Element || + !EQUAL(psSubTree->pszValue, osSubNode.c_str())) + continue; - psSubTree = CPLGetXMLNode(psSubTree, "PAMDataset"); - break; - } + if (!EQUAL(CPLGetXMLValue(psSubTree, "name", ""), + osSubNodeValue.c_str())) + continue; - if (psSubTree != nullptr) - psSubTree = CPLCloneXMLTree(psSubTree); + psSubTree = CPLGetXMLNode(psSubTree, "PAMDataset"); + break; + } - CPLDestroyXMLNode(psTree); - psTree = psSubTree; + if (psSubTree != nullptr) + psSubTree = CPLCloneXMLTree(psSubTree); + + CPLDestroyXMLNode(psTree); + psTree = psSubTree; + } } /* -------------------------------------------------------------------- */ @@ -1000,7 +1042,19 @@ CPLErr GDALPamDataset::TrySaveXML() /* the subdataset tree within the whole existing pam tree, */ /* after removing any old version of the same subdataset. */ /* -------------------------------------------------------------------- */ + std::string osSubNode; + std::string osSubNodeValue; if (!psPam->osSubdatasetName.empty()) + { + osSubNode = "Subdataset"; + osSubNodeValue = psPam->osSubdatasetName; + } + else if (!psPam->osDerivedDatasetName.empty()) + { + osSubNode = "DerivedDataset"; + osSubNodeValue = psPam->osDerivedDatasetName; + } + if (!osSubNode.empty()) { CPLXMLNode *psOldTree = nullptr; @@ -1022,11 +1076,11 @@ CPLErr GDALPamDataset::TrySaveXML() psSubTree = psSubTree->psNext) { if (psSubTree->eType != CXT_Element || - !EQUAL(psSubTree->pszValue, "Subdataset")) + !EQUAL(psSubTree->pszValue, osSubNode.c_str())) continue; if (!EQUAL(CPLGetXMLValue(psSubTree, "name", ""), - psPam->osSubdatasetName)) + osSubNodeValue.c_str())) continue; break; @@ -1034,9 +1088,10 @@ CPLErr GDALPamDataset::TrySaveXML() if (psSubTree == nullptr) { - psSubTree = CPLCreateXMLNode(psOldTree, CXT_Element, "Subdataset"); + psSubTree = + CPLCreateXMLNode(psOldTree, CXT_Element, osSubNode.c_str()); CPLCreateXMLNode(CPLCreateXMLNode(psSubTree, CXT_Attribute, "name"), - CXT_Text, psPam->osSubdatasetName); + CXT_Text, osSubNodeValue.c_str()); } CPLXMLNode *psOldPamDataset = CPLGetXMLNode(psSubTree, "PAMDataset");