Line data    Source code

       1             : /******************************************************************************
       2             :  * Name:     gdalmultidim_subsetdimension.cpp
       3             :  * Project:  GDAL Core
       4             :  * Purpose:  GDALGroup::SubsetDimensionFromSelection() implementation
       5             :  * Author:   Even Rouault <even.rouault at spatialys.com>
       6             :  *
       7             :  ******************************************************************************
       8             :  * Copyright (c) 2023, Even Rouault <even.rouault at spatialys.com>
       9             :  *
      10             :  * Permission is hereby granted, free of charge, to any person obtaining a
      11             :  * copy of this software and associated documentation files (the "Software"),
      12             :  * to deal in the Software without restriction, including without limitation
      13             :  * the rights to use, copy, modify, merge, publish, distribute, sublicense,
      14             :  * and/or sell copies of the Software, and to permit persons to whom the
      15             :  * Software is furnished to do so, subject to the following conditions:
      16             :  *
      17             :  * The above copyright notice and this permission notice shall be included
      18             :  * in all copies or substantial portions of the Software.
      19             :  *
      20             :  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
      21             :  * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
      22             :  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
      23             :  * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
      24             :  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
      25             :  * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
      26             :  * DEALINGS IN THE SOFTWARE.
      27             :  ****************************************************************************/
      28             : 
      29             : #include "gdal_priv.h"
      30             : #include "gdal_pam.h"
      31             : 
      32             : #include <algorithm>
      33             : 
      34             : /************************************************************************/
      35             : /*                           GetParentName()                            */
      36             : /************************************************************************/
      37             : 
      38          42 : static std::string GetParentName(const std::string &osPath)
      39             : {
      40          42 :     if (osPath == "/" || osPath.rfind('/') == 0)
      41          36 :         return "/";
      42           6 :     return osPath.substr(0, osPath.rfind('/'));
      43             : }
      44             : 
      45             : /************************************************************************/
      46             : /*                   GDALSubsetGroupSharedResources                     */
      47             : /************************************************************************/
      48             : 
      49             : struct GDALSubsetGroupSharedResources
      50             : {
      51             :     std::shared_ptr<GDALGroup> m_poRootGroup{};  // may be nullptr
      52             :     std::string m_osDimFullName{};
      53             :     std::vector<int> m_anMapNewDimToOldDim{};
      54             :     std::string m_osSelection{};
      55             :     std::shared_ptr<GDALDimension> m_poNewDim{};
      56             :     std::shared_ptr<GDALMDArray> m_poNewIndexingVar{};
      57             : };
      58             : 
      59             : /************************************************************************/
      60             : /*                          CreateContext()                             */
      61             : /************************************************************************/
      62             : 
      63             : static std::string
      64          14 : CreateContext(const std::string &osParentContext,
      65             :               const std::shared_ptr<GDALSubsetGroupSharedResources> &poShared)
      66             : {
      67          14 :     std::string osRet(osParentContext);
      68          14 :     if (!osRet.empty())
      69           0 :         osRet += ". ";
      70          14 :     osRet += "Selection ";
      71          14 :     osRet += poShared->m_osSelection;
      72          14 :     return osRet;
      73             : }
      74             : 
      75             : /************************************************************************/
      76             : /*                           GDALSubsetGroup                            */
      77             : /************************************************************************/
      78             : 
      79             : class GDALSubsetGroup final : public GDALGroup
      80             : {
      81             :     std::shared_ptr<GDALGroup> m_poParent{};
      82             :     std::shared_ptr<GDALSubsetGroupSharedResources> m_poShared{};
      83             : 
      84           9 :     GDALSubsetGroup(
      85             :         const std::shared_ptr<GDALGroup> &poParent,
      86             :         const std::shared_ptr<GDALSubsetGroupSharedResources> &poShared)
      87          18 :         : GDALGroup(GetParentName(poParent->GetFullName()), poParent->GetName(),
      88          18 :                     CreateContext(poParent->GetContext(), poShared)),
      89          36 :           m_poParent(std::move(poParent)), m_poShared(std::move(poShared))
      90             :     {
      91           9 :     }
      92             : 
      93             :   public:
      94             :     static std::shared_ptr<GDALGroup>
      95           9 :     Create(const std::shared_ptr<GDALGroup> &poParent,
      96             :            const std::shared_ptr<GDALSubsetGroupSharedResources> &poShared)
      97             :     {
      98             :         auto poGroup = std::shared_ptr<GDALSubsetGroup>(
      99          18 :             new GDALSubsetGroup(poParent, poShared));
     100           9 :         poGroup->SetSelf(poGroup);
     101          18 :         return poGroup;
     102             :     }
     103             : 
     104             :     std::vector<std::string>
     105           2 :     GetMDArrayNames(CSLConstList papszOptions = nullptr) const override
     106             :     {
     107           2 :         return m_poParent->GetMDArrayNames(papszOptions);
     108             :     }
     109             : 
     110             :     std::shared_ptr<GDALMDArray>
     111             :     OpenMDArray(const std::string &osName,
     112             :                 CSLConstList papszOptions = nullptr) const override;
     113             : 
     114             :     std::vector<std::string>
     115           1 :     GetGroupNames(CSLConstList papszOptions = nullptr) const override
     116             :     {
     117           1 :         return m_poParent->GetGroupNames(papszOptions);
     118             :     }
     119             : 
     120             :     std::shared_ptr<GDALGroup>
     121             :     OpenGroup(const std::string &osName,
     122             :               CSLConstList papszOptions = nullptr) const override;
     123             : 
     124             :     std::vector<std::shared_ptr<GDALDimension>>
     125             :     GetDimensions(CSLConstList papszOptions = nullptr) const override;
     126             : 
     127             :     std::shared_ptr<GDALAttribute>
     128           2 :     GetAttribute(const std::string &osName) const override
     129             :     {
     130           2 :         return m_poParent->GetAttribute(osName);
     131             :     }
     132             : 
     133             :     std::vector<std::shared_ptr<GDALAttribute>>
     134           1 :     GetAttributes(CSLConstList papszOptions = nullptr) const override
     135             :     {
     136           1 :         return m_poParent->GetAttributes(papszOptions);
     137             :     }
     138             : };
     139             : 
     140             : /************************************************************************/
     141             : /*                           GDALSubsetArray                            */
     142             : /************************************************************************/
     143             : 
     144             : class GDALSubsetArray final : public GDALPamMDArray
     145             : {
     146             :   private:
     147             :     std::shared_ptr<GDALMDArray> m_poParent{};
     148             :     std::shared_ptr<GDALSubsetGroupSharedResources> m_poShared{};
     149             :     std::vector<std::shared_ptr<GDALDimension>> m_apoDims{};
     150             :     std::vector<bool> m_abPatchedDim{};
     151             :     bool m_bPatchedDimIsFirst = false;
     152             : 
     153             :   protected:
     154          14 :     GDALSubsetArray(
     155             :         const std::shared_ptr<GDALMDArray> &poParent,
     156             :         const std::shared_ptr<GDALSubsetGroupSharedResources> &poShared,
     157             :         const std::string &osContext)
     158          28 :         : GDALAbstractMDArray(GetParentName(poParent->GetFullName()),
     159          14 :                               poParent->GetName()),
     160          28 :           GDALPamMDArray(GetParentName(poParent->GetFullName()),
     161          28 :                          poParent->GetName(), GDALPamMultiDim::GetPAM(poParent),
     162             :                          osContext),
     163          84 :           m_poParent(std::move(poParent)), m_poShared(std::move(poShared))
     164             :     {
     165          14 :         m_apoDims = m_poParent->GetDimensions();
     166          39 :         for (size_t i = 0; i < m_apoDims.size(); ++i)
     167             :         {
     168          25 :             auto &poDim = m_apoDims[i];
     169          25 :             if (poDim->GetFullName() == m_poShared->m_osDimFullName)
     170             :             {
     171          15 :                 m_bPatchedDimIsFirst = (i == 0);
     172          15 :                 poDim = m_poShared->m_poNewDim;
     173          15 :                 m_abPatchedDim.push_back(true);
     174             :             }
     175             :             else
     176             :             {
     177          10 :                 m_abPatchedDim.push_back(false);
     178             :             }
     179             :         }
     180          14 :     }
     181             : 
     182             :     bool IRead(const GUInt64 *arrayStartIdx, const size_t *count,
     183             :                const GInt64 *arrayStep, const GPtrDiff_t *bufferStride,
     184             :                const GDALExtendedDataType &bufferDataType,
     185             :                void *pDstBuffer) const override;
     186             : 
     187             :   public:
     188             :     static std::shared_ptr<GDALSubsetArray>
     189          14 :     Create(const std::shared_ptr<GDALMDArray> &poParent,
     190             :            const std::shared_ptr<GDALSubsetGroupSharedResources> &poShared,
     191             :            const std::string &osContext)
     192             :     {
     193             :         auto newAr(std::shared_ptr<GDALSubsetArray>(
     194          14 :             new GDALSubsetArray(poParent, poShared, osContext)));
     195          14 :         newAr->SetSelf(newAr);
     196          14 :         return newAr;
     197             :     }
     198             : 
     199           2 :     bool IsWritable() const override
     200             :     {
     201           2 :         return false;
     202             :     }
     203             : 
     204          19 :     const std::string &GetFilename() const override
     205             :     {
     206          19 :         return m_poParent->GetFilename();
     207             :     }
     208             : 
     209             :     const std::vector<std::shared_ptr<GDALDimension>> &
     210         108 :     GetDimensions() const override
     211             :     {
     212         108 :         return m_apoDims;
     213             :     }
     214             : 
     215          62 :     const GDALExtendedDataType &GetDataType() const override
     216             :     {
     217          62 :         return m_poParent->GetDataType();
     218             :     }
     219             : 
     220           3 :     const std::string &GetUnit() const override
     221             :     {
     222           3 :         return m_poParent->GetUnit();
     223             :     }
     224             : 
     225           1 :     std::shared_ptr<OGRSpatialReference> GetSpatialRef() const override
     226             :     {
     227           1 :         return m_poParent->GetSpatialRef();
     228             :     }
     229             : 
     230           8 :     const void *GetRawNoDataValue() const override
     231             :     {
     232           8 :         return m_poParent->GetRawNoDataValue();
     233             :     }
     234             : 
     235           4 :     std::vector<GUInt64> GetBlockSize() const override
     236             :     {
     237           4 :         std::vector<GUInt64> ret(m_poParent->GetBlockSize());
     238          13 :         for (size_t i = 0; i < m_apoDims.size(); ++i)
     239             :         {
     240           9 :             if (m_abPatchedDim[i])
     241           4 :                 ret[1] = 1;
     242             :         }
     243           4 :         return ret;
     244             :     }
     245             : 
     246             :     std::shared_ptr<GDALAttribute>
     247           8 :     GetAttribute(const std::string &osName) const override
     248             :     {
     249           8 :         return m_poParent->GetAttribute(osName);
     250             :     }
     251             : 
     252             :     std::vector<std::shared_ptr<GDALAttribute>>
     253           3 :     GetAttributes(CSLConstList papszOptions = nullptr) const override
     254             :     {
     255           3 :         return m_poParent->GetAttributes(papszOptions);
     256             :     }
     257             : 
     258           1 :     std::shared_ptr<GDALGroup> GetRootGroup() const override
     259             :     {
     260           1 :         if (m_poShared->m_poRootGroup)
     261             :         {
     262           2 :             return GDALSubsetGroup::Create(m_poShared->m_poRootGroup,
     263           1 :                                            m_poShared);
     264             :         }
     265           0 :         return nullptr;
     266             :     }
     267             : };
     268             : 
     269             : /************************************************************************/
     270             : /*                            OpenMDArray()                             */
     271             : /************************************************************************/
     272             : 
     273             : std::shared_ptr<GDALMDArray>
     274          12 : GDALSubsetGroup::OpenMDArray(const std::string &osName,
     275             :                              CSLConstList papszOptions) const
     276             : {
     277          24 :     auto poArray = m_poParent->OpenMDArray(osName, papszOptions);
     278          12 :     if (poArray)
     279             :     {
     280          22 :         for (const auto &poDim : poArray->GetDimensions())
     281             :         {
     282          20 :             if (poDim->GetFullName() == m_poShared->m_osDimFullName)
     283             :             {
     284          18 :                 return GDALSubsetArray::Create(poArray, m_poShared,
     285           9 :                                                GetContext());
     286             :             }
     287             :         }
     288             :     }
     289           3 :     return poArray;
     290             : }
     291             : 
     292             : /************************************************************************/
     293             : /*                             OpenGroup()                              */
     294             : /************************************************************************/
     295             : 
     296             : std::shared_ptr<GDALGroup>
     297           4 : GDALSubsetGroup::OpenGroup(const std::string &osName,
     298             :                            CSLConstList papszOptions) const
     299             : {
     300           4 :     auto poSubGroup = m_poParent->OpenGroup(osName, papszOptions);
     301           4 :     if (poSubGroup)
     302             :     {
     303           3 :         poSubGroup = GDALSubsetGroup::Create(poSubGroup, m_poShared);
     304             :     }
     305           4 :     return poSubGroup;
     306             : }
     307             : 
     308             : /************************************************************************/
     309             : /*                             GetDimensions()                          */
     310             : /************************************************************************/
     311             : 
     312             : std::vector<std::shared_ptr<GDALDimension>>
     313           1 : GDALSubsetGroup::GetDimensions(CSLConstList papszOptions) const
     314             : {
     315           1 :     auto apoDims = m_poParent->GetDimensions(papszOptions);
     316           5 :     for (auto &poDim : apoDims)
     317             :     {
     318           4 :         if (poDim->GetFullName() == m_poShared->m_osDimFullName)
     319             :         {
     320           1 :             poDim = m_poShared->m_poNewDim;
     321             :         }
     322             :     }
     323           1 :     return apoDims;
     324             : }
     325             : 
     326             : /************************************************************************/
     327             : /*                             IRead()                                  */
     328             : /************************************************************************/
     329             : 
     330          28 : bool GDALSubsetArray::IRead(const GUInt64 *arrayStartIdx, const size_t *count,
     331             :                             const GInt64 *arrayStep,
     332             :                             const GPtrDiff_t *bufferStride,
     333             :                             const GDALExtendedDataType &bufferDataType,
     334             :                             void *pDstBuffer) const
     335             : {
     336          28 :     const auto nDims = m_apoDims.size();
     337          56 :     std::vector<GUInt64> newArrayStartIdx(nDims);
     338             :     // the +1 in nDims + 1 is to make happy -Werror=null-dereference when
     339             :     // doing newCount[0] = 1 and newArrayStep[0] = 1
     340          56 :     std::vector<size_t> newCount(nDims + 1, 1);
     341          56 :     std::vector<GInt64> newArrayStep(nDims + 1, 1);
     342          28 :     const size_t nBufferDTSize = bufferDataType.GetSize();
     343             : 
     344          28 :     if (m_bPatchedDimIsFirst)
     345             :     {
     346             :         // Optimized case when the only patched dimension is the first one.
     347           3 :         std::copy_n(arrayStartIdx, nDims, newArrayStartIdx.data());
     348           3 :         std::copy_n(count, nDims, newCount.data());
     349           3 :         std::copy_n(arrayStep, nDims, newArrayStep.data());
     350           3 :         GUInt64 arrayIdx = arrayStartIdx[0];
     351           3 :         GByte *pabyDstBuffer = static_cast<GByte *>(pDstBuffer);
     352             : #if defined(__GNUC__)
     353             : #pragma GCC diagnostic push
     354             : #pragma GCC diagnostic ignored "-Wnull-dereference"
     355             : #endif
     356           3 :         newCount[0] = 1;
     357           3 :         newArrayStep[0] = 1;
     358             : #if defined(__GNUC__)
     359             : #pragma GCC diagnostic pop
     360             : #endif
     361           9 :         for (size_t i = 0; i < count[0]; ++i)
     362             :         {
     363           6 :             if (i > 0)
     364             :             {
     365           3 :                 if (arrayStep[0] > 0)
     366           2 :                     arrayIdx += arrayStep[0];
     367             :                 else
     368           1 :                     arrayIdx -= static_cast<GUInt64>(-arrayStep[0]);
     369           3 :                 pabyDstBuffer += bufferStride[0] * nBufferDTSize;
     370             :             }
     371           6 :             newArrayStartIdx[0] =
     372           6 :                 m_poShared->m_anMapNewDimToOldDim[static_cast<int>(arrayIdx)];
     373          12 :             if (!m_poParent->Read(newArrayStartIdx.data(), newCount.data(),
     374           6 :                                   newArrayStep.data(), bufferStride,
     375             :                                   bufferDataType, pabyDstBuffer))
     376             :             {
     377           0 :                 return false;
     378             :             }
     379             :         }
     380           3 :         return true;
     381             :     }
     382             : 
     383             :     // Slow/unoptimized case
     384          50 :     std::vector<size_t> anStackIter(nDims);
     385          50 :     std::vector<GUInt64> anStackArrayIdx(nDims);
     386          50 :     std::vector<GByte *> pabyDstBufferStack(nDims + 1);
     387             : #if defined(__GNUC__)
     388             : #pragma GCC diagnostic push
     389             : #pragma GCC diagnostic ignored "-Wnull-dereference"
     390             : #endif
     391          25 :     pabyDstBufferStack[0] = static_cast<GByte *>(pDstBuffer);
     392             : #if defined(__GNUC__)
     393             : #pragma GCC diagnostic pop
     394             : #endif
     395          25 :     size_t iDim = 0;
     396         183 : lbl_next_depth:
     397         183 :     if (iDim == nDims)
     398             :     {
     399         160 :         if (!m_poParent->Read(newArrayStartIdx.data(), newCount.data(),
     400          80 :                               newArrayStep.data(), bufferStride, bufferDataType,
     401          80 :                               pabyDstBufferStack[iDim]))
     402             :         {
     403           0 :             return false;
     404             :         }
     405             :     }
     406             :     else
     407             :     {
     408         103 :         anStackIter[iDim] = 0;
     409         103 :         anStackArrayIdx[iDim] = arrayStartIdx[iDim];
     410             :         while (true)
     411             :         {
     412         158 :             if (m_abPatchedDim[iDim])
     413             :             {
     414          84 :                 newArrayStartIdx[iDim] =
     415          84 :                     m_poShared->m_anMapNewDimToOldDim[static_cast<int>(
     416          84 :                         anStackArrayIdx[iDim])];
     417             :             }
     418             :             else
     419             :             {
     420          74 :                 newArrayStartIdx[iDim] = anStackArrayIdx[iDim];
     421             :             }
     422         158 :             ++iDim;
     423         158 :             pabyDstBufferStack[iDim] = pabyDstBufferStack[iDim - 1];
     424         158 :             goto lbl_next_depth;
     425         158 :         lbl_return_to_caller_in_loop:
     426         158 :             --iDim;
     427         158 :             ++anStackIter[iDim];
     428         158 :             if (anStackIter[iDim] == count[iDim])
     429         103 :                 break;
     430          55 :             if (arrayStep[iDim] > 0)
     431          54 :                 anStackArrayIdx[iDim] += arrayStep[iDim];
     432             :             else
     433           1 :                 anStackArrayIdx[iDim] -= -arrayStep[iDim];
     434          55 :             pabyDstBufferStack[iDim] += bufferStride[iDim] * nBufferDTSize;
     435             :         }
     436             :     }
     437         183 :     if (iDim > 0)
     438         158 :         goto lbl_return_to_caller_in_loop;
     439             : 
     440          25 :     return true;
     441             : }
     442             : 
     443             : /************************************************************************/
     444             : /*                   SubsetDimensionFromSelection()                     */
     445             : /************************************************************************/
     446             : 
     447             : /** Return a virtual group whose one dimension has been subset according to a
     448             :  * selection.
     449             :  *
     450             :  * The selection criterion is currently restricted to the form
     451             :  * "/path/to/array=numeric_value" (no spaces around equal)
     452             :  *
     453             :  * This is similar to XArray indexing by name and label on a XArray Dataset
     454             :  * using the sel() method.
     455             :  * Cf https://docs.xarray.dev/en/latest/user-guide/indexing.html#quick-overview
     456             :  *
     457             :  * For example on a EMIT L2A product
     458             :  * (https://github.com/nasa/EMIT-Data-Resources/blob/main/python/tutorials/Exploring_EMIT_L2A_Reflectance.ipynb),
     459             :  * this can be used to keep only valid bands with
     460             :  * SubsetDimensionFromSelection("/sensor_band_parameters/good_wavelengths=1")
     461             :  *
     462             :  * This is the same as the C function GDALGroupSubsetDimensionFromSelection().
     463             :  *
     464             :  * @param osSelection Selection criterion.
     465             :  * @return a virtual group, or nullptr in case of error
     466             :  * @since 3.8
     467             :  */
     468             : std::shared_ptr<GDALGroup>
     469          14 : GDALGroup::SubsetDimensionFromSelection(const std::string &osSelection) const
     470             : {
     471          28 :     auto self = std::dynamic_pointer_cast<GDALGroup>(m_pSelf.lock());
     472          14 :     if (!self)
     473             :     {
     474           0 :         CPLError(CE_Failure, CPLE_AppDefined,
     475             :                  "Driver implementation issue: m_pSelf not set !");
     476           0 :         return nullptr;
     477             :     }
     478             : 
     479          14 :     const auto nEqualPos = osSelection.find('=');
     480          14 :     if (nEqualPos == std::string::npos)
     481             :     {
     482           2 :         CPLError(CE_Failure, CPLE_AppDefined, "Invalid value for selection");
     483           2 :         return nullptr;
     484             :     }
     485          24 :     const auto osArrayName = osSelection.substr(0, nEqualPos);
     486          24 :     const auto osValue = osSelection.substr(nEqualPos + 1);
     487          13 :     if (CPLGetValueType(osValue.c_str()) != CPL_VALUE_INTEGER &&
     488           1 :         CPLGetValueType(osValue.c_str()) != CPL_VALUE_REAL)
     489             :     {
     490           1 :         CPLError(CE_Failure, CPLE_AppDefined,
     491             :                  "Non-numeric value in selection criterion");
     492           1 :         return nullptr;
     493             :     }
     494          22 :     auto poArray = OpenMDArrayFromFullname(osArrayName);
     495          11 :     if (!poArray)
     496             :     {
     497           1 :         CPLError(CE_Failure, CPLE_AppDefined, "Cannot find array %s",
     498             :                  osArrayName.c_str());
     499           1 :         return nullptr;
     500             :     }
     501          10 :     if (poArray->GetDimensionCount() != 1)
     502             :     {
     503           1 :         CPLError(CE_Failure, CPLE_AppDefined,
     504             :                  "Array %s is not single dimensional", osArrayName.c_str());
     505           1 :         return nullptr;
     506             :     }
     507           9 :     if (poArray->GetDataType().GetClass() != GEDTC_NUMERIC)
     508             :     {
     509           1 :         CPLError(CE_Failure, CPLE_AppDefined, "Array %s is not of numeric type",
     510             :                  osArrayName.c_str());
     511           1 :         return nullptr;
     512             :     }
     513             : 
     514           8 :     const auto nElts = poArray->GetTotalElementsCount();
     515           8 :     if (nElts > 10 * 1024 * 1024)
     516             :     {
     517           1 :         CPLError(CE_Failure, CPLE_AppDefined, "Too many values in %s",
     518             :                  osArrayName.c_str());
     519           1 :         return nullptr;
     520             :     }
     521          14 :     std::vector<double> values;
     522             :     try
     523             :     {
     524           7 :         values.resize(static_cast<size_t>(nElts));
     525             :     }
     526           0 :     catch (const std::bad_alloc &e)
     527             :     {
     528           0 :         CPLError(CE_Failure, CPLE_AppDefined, "Out of memory: %s", e.what());
     529           0 :         return nullptr;
     530             :     }
     531           7 :     const GUInt64 startIdx[1] = {0};
     532           7 :     const size_t count[1] = {values.size()};
     533          14 :     if (!poArray->Read(startIdx, count, nullptr, nullptr,
     534          14 :                        GDALExtendedDataType::Create(GDT_Float64), &values[0],
     535           7 :                        values.data(), values.size() * sizeof(values[0])))
     536             :     {
     537           0 :         return nullptr;
     538             :     }
     539           7 :     const double dfSelectionValue = CPLAtof(osValue.c_str());
     540          14 :     std::vector<int> anMapNewDimToOldDim;
     541          77 :     for (int i = 0; i < static_cast<int>(nElts); ++i)
     542             :     {
     543          70 :         if (values[i] == dfSelectionValue)
     544           8 :             anMapNewDimToOldDim.push_back(i);
     545             :     }
     546           7 :     if (anMapNewDimToOldDim.empty())
     547             :     {
     548           1 :         CPLError(CE_Failure, CPLE_AppDefined, "No value in %s matching %f",
     549             :                  osArrayName.c_str(), dfSelectionValue);
     550           1 :         return nullptr;
     551             :     }
     552           6 :     if (anMapNewDimToOldDim.size() == nElts)
     553             :     {
     554           1 :         return self;
     555             :     }
     556             : 
     557          10 :     auto poDim = poArray->GetDimensions()[0];
     558          10 :     auto poShared = std::make_shared<GDALSubsetGroupSharedResources>();
     559           5 :     if (GetFullName() == "/")
     560           5 :         poShared->m_poRootGroup = self;
     561           5 :     poShared->m_osSelection = osSelection;
     562           5 :     poShared->m_osDimFullName = poArray->GetDimensions()[0]->GetFullName();
     563           5 :     poShared->m_anMapNewDimToOldDim = std::move(anMapNewDimToOldDim);
     564             : 
     565             :     // Create a modified dimension of reduced size
     566             :     auto poNewDim = std::make_shared<GDALDimensionWeakIndexingVar>(
     567           5 :         GetParentName(poDim->GetFullName()), poDim->GetName(), poDim->GetType(),
     568          15 :         poDim->GetDirection(), poShared->m_anMapNewDimToOldDim.size());
     569           5 :     poShared->m_poNewDim = poNewDim;
     570             : 
     571          10 :     auto poIndexingVar = poDim->GetIndexingVariable();
     572           5 :     if (poIndexingVar)
     573             :     {
     574             :         // poNewIndexingVar must be created with a different GDALSubsetGroupSharedResources
     575             :         // instance than poShared, to avoid cross reference, that would result in
     576             :         // objects not being freed !
     577             :         auto poSpecificShared =
     578          10 :             std::make_shared<GDALSubsetGroupSharedResources>();
     579           5 :         poSpecificShared->m_poRootGroup = poShared->m_poRootGroup;
     580           5 :         poSpecificShared->m_osSelection = osSelection;
     581           5 :         poSpecificShared->m_osDimFullName =
     582          10 :             poArray->GetDimensions()[0]->GetFullName();
     583           5 :         poSpecificShared->m_anMapNewDimToOldDim =
     584           5 :             poShared->m_anMapNewDimToOldDim;
     585           5 :         poSpecificShared->m_poNewDim = poNewDim;
     586             :         auto poNewIndexingVar =
     587             :             GDALSubsetArray::Create(poIndexingVar, poSpecificShared,
     588          10 :                                     CreateContext(GetContext(), poShared));
     589           5 :         poNewDim->SetIndexingVariable(poNewIndexingVar);
     590           5 :         poShared->m_poNewIndexingVar = poNewIndexingVar;
     591             :     }
     592             : 
     593           5 :     return GDALSubsetGroup::Create(self, poShared);
     594             : }