Coverage Report

Created: 2026-08-11 08:26

next uncovered line (L), next uncovered region (R), next uncovered branch (B)
/src/gdal/frmts/zarr/zarrdriver.cpp
Line
Count
Source
1
/******************************************************************************
2
 *
3
 * Project:  GDAL
4
 * Purpose:  Zarr driver
5
 * Author:   Even Rouault <even dot rouault at spatialys.com>
6
 *
7
 ******************************************************************************
8
 * Copyright (c) 2021, Even Rouault <even dot rouault at spatialys.com>
9
 *
10
 * SPDX-License-Identifier: MIT
11
 ****************************************************************************/
12
13
#include "zarr.h"
14
#include "zarrdrivercore.h"
15
#include "vsikerchunk.h"
16
17
#include "cpl_minixml.h"
18
19
#include "gdalalgorithm.h"
20
#include "gdal_frmts.h"
21
22
#include <algorithm>
23
#include <cassert>
24
#include <cinttypes>
25
#include <cmath>
26
#include <limits>
27
#include <future>
28
#include <mutex>
29
30
#ifdef HAVE_BLOSC
31
#include <blosc.h>
32
#endif
33
34
/************************************************************************/
35
/*                            ZarrDataset()                             */
36
/************************************************************************/
37
38
ZarrDataset::ZarrDataset(const std::shared_ptr<ZarrGroupBase> &poRootGroup)
39
76
    : m_poRootGroup(poRootGroup)
40
76
{
41
76
}
42
43
/************************************************************************/
44
/*                            OpenMultidim()                            */
45
/************************************************************************/
46
47
GDALDataset *ZarrDataset::OpenMultidim(const char *pszFilename,
48
                                       bool bUpdateMode,
49
                                       CSLConstList papszOpenOptionsIn)
50
273
{
51
273
    CPLString osFilename(pszFilename);
52
273
    if (osFilename.back() == '/')
53
12
        osFilename.pop_back();
54
55
273
    auto poSharedResource = ZarrSharedResource::Create(osFilename, bUpdateMode);
56
273
    poSharedResource->SetOpenOptions(papszOpenOptionsIn);
57
58
273
    auto poRG = poSharedResource->GetRootGroup();
59
273
    if (!poRG)
60
273
    {
61
        // Kerchunk Parquet auto-detection: OpenRootGroup found a
62
        // .zmetadata with record_size, signaling a redirect.
63
273
        const auto &osKerchunkPath = poSharedResource->GetKerchunkParquetPath();
64
273
        if (!osKerchunkPath.empty())
65
0
            return OpenMultidim(osKerchunkPath.c_str(), bUpdateMode,
66
0
                                papszOpenOptionsIn);
67
273
        return nullptr;
68
273
    }
69
0
    return new ZarrDataset(poRG);
70
273
}
71
72
/************************************************************************/
73
/*                            ExploreGroup()                            */
74
/************************************************************************/
75
76
static bool ExploreGroup(const std::shared_ptr<GDALGroup> &poGroup,
77
                         std::vector<std::string> &aosArrays, int nRecCount)
78
0
{
79
0
    if (nRecCount == 32)
80
0
    {
81
0
        CPLError(CE_Failure, CPLE_NotSupported,
82
0
                 "Too deep recursion level in ExploreGroup()");
83
0
        return false;
84
0
    }
85
0
    const auto aosGroupArrayNames = poGroup->GetMDArrayNames();
86
0
    for (const auto &osArrayName : aosGroupArrayNames)
87
0
    {
88
0
        std::string osArrayFullname = poGroup->GetFullName();
89
0
        if (osArrayName != "/")
90
0
        {
91
0
            if (osArrayFullname != "/")
92
0
                osArrayFullname += '/';
93
0
            osArrayFullname += osArrayName;
94
0
        }
95
0
        aosArrays.emplace_back(std::move(osArrayFullname));
96
0
        if (aosArrays.size() == 10000)
97
0
        {
98
0
            CPLError(CE_Failure, CPLE_NotSupported,
99
0
                     "Too many arrays found by ExploreGroup()");
100
0
            return false;
101
0
        }
102
0
    }
103
104
0
    const auto aosSubGroups = poGroup->GetGroupNames();
105
0
    for (const auto &osSubGroup : aosSubGroups)
106
0
    {
107
0
        const auto poSubGroup = poGroup->OpenGroup(osSubGroup);
108
0
        if (poSubGroup)
109
0
        {
110
0
            if (!ExploreGroup(poSubGroup, aosArrays, nRecCount + 1))
111
0
                return false;
112
0
        }
113
0
    }
114
0
    return true;
115
0
}
116
117
/************************************************************************/
118
/*                          GetMetadataItem()                           */
119
/************************************************************************/
120
121
const char *ZarrDataset::GetMetadataItem(const char *pszName,
122
                                         const char *pszDomain)
123
76
{
124
76
    if (pszDomain != nullptr && EQUAL(pszDomain, GDAL_MDD_SUBDATASETS))
125
0
        return m_aosSubdatasets.FetchNameValue(pszName);
126
76
    if (pszDomain != nullptr && EQUAL(pszDomain, GDAL_MDD_IMAGE_STRUCTURE))
127
76
        return GDALDataset::GetMetadataItem(pszName, pszDomain);
128
0
    return nullptr;
129
76
}
130
131
/************************************************************************/
132
/*                            GetMetadata()                             */
133
/************************************************************************/
134
135
CSLConstList ZarrDataset::GetMetadata(const char *pszDomain)
136
0
{
137
0
    if (pszDomain != nullptr && EQUAL(pszDomain, GDAL_MDD_SUBDATASETS))
138
0
        return m_aosSubdatasets.List();
139
0
    if (pszDomain != nullptr && EQUAL(pszDomain, GDAL_MDD_IMAGE_STRUCTURE))
140
0
        return GDALDataset::GetMetadata(pszDomain);
141
0
    return nullptr;
142
0
}
143
144
/************************************************************************/
145
/*                          IBuildOverviews()                           */
146
/************************************************************************/
147
148
CPLErr ZarrDataset::IBuildOverviews(const char *pszResampling, int nOverviews,
149
                                    const int *panOverviewList,
150
                                    int /* nListBands */,
151
                                    const int * /*panBandList*/,
152
                                    GDALProgressFunc pfnProgress,
153
                                    void *pProgressData,
154
                                    CSLConstList papszOptions)
155
0
{
156
0
    for (int i = 0; i < nBands; ++i)
157
0
    {
158
0
        auto poBand = cpl::down_cast<ZarrRasterBand *>(papoBands[i]);
159
0
        for (auto &[_, value] : poBand->m_oMapOverview)
160
0
        {
161
0
            poBand->m_aoOverviewOld.push_back(std::move(value));
162
0
        }
163
0
        poBand->m_oMapOverview.clear();
164
0
    }
165
166
0
    if (m_poSingleArray)
167
0
    {
168
0
        return m_poSingleArray->BuildOverviews(pszResampling, nOverviews,
169
0
                                               panOverviewList, pfnProgress,
170
0
                                               pProgressData, papszOptions);
171
0
    }
172
0
    else
173
0
    {
174
0
        CPLErr eErr = CE_None;
175
0
        for (int i = 0; i < nBands && eErr == CE_None; ++i)
176
0
        {
177
0
            std::unique_ptr<void, decltype(&GDALDestroyScaledProgress)>
178
0
                pScaledProgress(
179
0
                    GDALCreateScaledProgress(static_cast<double>(i) /
180
0
                                                 static_cast<double>(nBands),
181
0
                                             static_cast<double>(i + 1) /
182
0
                                                 static_cast<double>(nBands),
183
0
                                             pfnProgress, pProgressData),
184
0
                    GDALDestroyScaledProgress);
185
0
            eErr = cpl::down_cast<ZarrRasterBand *>(papoBands[i])
186
0
                       ->m_poArray->BuildOverviews(
187
0
                           pszResampling, nOverviews, panOverviewList,
188
0
                           pScaledProgress ? GDALScaledProgress : nullptr,
189
0
                           pScaledProgress.get(), papszOptions);
190
0
        }
191
0
        return eErr;
192
0
    }
193
0
}
194
195
/************************************************************************/
196
/*                       GetXYDimensionIndices()                        */
197
/************************************************************************/
198
199
static void GetXYDimensionIndices(const std::shared_ptr<GDALMDArray> &poArray,
200
                                  const GDALOpenInfo *poOpenInfo, size_t &iXDim,
201
                                  size_t &iYDim)
202
0
{
203
0
    const size_t nDims = poArray->GetDimensionCount();
204
0
    iYDim = nDims >= 2 ? nDims - 2 : 0;
205
0
    iXDim = nDims >= 2 ? nDims - 1 : 0;
206
207
0
    if (nDims >= 2)
208
0
    {
209
0
        const char *pszDimX =
210
0
            CSLFetchNameValue(poOpenInfo->papszOpenOptions, "DIM_X");
211
0
        const char *pszDimY =
212
0
            CSLFetchNameValue(poOpenInfo->papszOpenOptions, "DIM_Y");
213
0
        bool bFoundX = false;
214
0
        bool bFoundY = false;
215
0
        const auto &apoDims = poArray->GetDimensions();
216
0
        for (size_t i = 0; i < nDims; ++i)
217
0
        {
218
0
            if (pszDimX && apoDims[i]->GetName() == pszDimX)
219
0
            {
220
0
                bFoundX = true;
221
0
                iXDim = i;
222
0
            }
223
0
            else if (pszDimY && apoDims[i]->GetName() == pszDimY)
224
0
            {
225
0
                bFoundY = true;
226
0
                iYDim = i;
227
0
            }
228
0
            else if (!pszDimX &&
229
0
                     (apoDims[i]->GetType() == GDAL_DIM_TYPE_HORIZONTAL_X ||
230
0
                      apoDims[i]->GetName() == "X"))
231
0
                iXDim = i;
232
0
            else if (!pszDimY &&
233
0
                     (apoDims[i]->GetType() == GDAL_DIM_TYPE_HORIZONTAL_Y ||
234
0
                      apoDims[i]->GetName() == "Y"))
235
0
                iYDim = i;
236
0
        }
237
0
        if (pszDimX)
238
0
        {
239
0
            if (!bFoundX && CPLGetValueType(pszDimX) == CPL_VALUE_INTEGER)
240
0
            {
241
0
                const int nTmp = atoi(pszDimX);
242
0
                if (nTmp >= 0 && nTmp <= static_cast<int>(nDims))
243
0
                {
244
0
                    iXDim = nTmp;
245
0
                    bFoundX = true;
246
0
                }
247
0
            }
248
0
            if (!bFoundX)
249
0
            {
250
0
                CPLError(CE_Warning, CPLE_AppDefined,
251
0
                         "Cannot find dimension DIM_X=%s", pszDimX);
252
0
            }
253
0
        }
254
0
        if (pszDimY)
255
0
        {
256
0
            if (!bFoundY && CPLGetValueType(pszDimY) == CPL_VALUE_INTEGER)
257
0
            {
258
0
                const int nTmp = atoi(pszDimY);
259
0
                if (nTmp >= 0 && nTmp <= static_cast<int>(nDims))
260
0
                {
261
0
                    iYDim = nTmp;
262
0
                    bFoundY = true;
263
0
                }
264
0
            }
265
0
            if (!bFoundY)
266
0
            {
267
0
                CPLError(CE_Warning, CPLE_AppDefined,
268
0
                         "Cannot find dimension DIM_Y=%s", pszDimY);
269
0
            }
270
0
        }
271
0
    }
272
0
}
273
274
/************************************************************************/
275
/*                       GetExtraDimSampleCount()                       */
276
/************************************************************************/
277
278
static uint64_t
279
GetExtraDimSampleCount(const std::shared_ptr<GDALMDArray> &poArray,
280
                       size_t iXDim, size_t iYDim)
281
0
{
282
0
    uint64_t nExtraDimSamples = 1;
283
0
    const auto &apoDims = poArray->GetDimensions();
284
0
    for (size_t i = 0; i < apoDims.size(); ++i)
285
0
    {
286
0
        if (i != iXDim && i != iYDim)
287
0
            nExtraDimSamples *= apoDims[i]->GetSize();
288
0
    }
289
0
    return nExtraDimSamples;
290
0
}
291
292
/************************************************************************/
293
/*                        PrefetchCoordArrays()                         */
294
/************************************************************************/
295
296
// Warm g_oCoordCache by reading X and Y coordinate arrays in parallel.
297
// For remote datasets this avoids sequential HTTP round-trips in
298
// GuessGeoTransform() (can save ~800 ms).  Each ZarrArray has its own
299
// mutex and VSI opens independent handles, so sibling reads are safe.
300
static void PrefetchCoordArrays(const std::shared_ptr<GDALMDArray> &poArray,
301
                                size_t iXDim, size_t iYDim)
302
0
{
303
0
    const auto nDimCount = poArray->GetDimensionCount();
304
0
    if (nDimCount < 2 || iXDim >= nDimCount || iYDim >= nDimCount)
305
0
        return;
306
0
    const auto &dims = poArray->GetDimensions();
307
0
    auto poVarX = dims[iXDim]->GetIndexingVariable();
308
0
    auto poVarY = dims[iYDim]->GetIndexingVariable();
309
0
    if (!poVarX || poVarX->GetDimensionCount() != 1 || !poVarY ||
310
0
        poVarY->GetDimensionCount() != 1)
311
0
        return;
312
0
    if (VSIIsLocal(poVarX->GetFilename().c_str()))
313
0
        return;
314
315
0
    double dfXStart = 0, dfXSpacing = 0, dfYStart = 0, dfYSpacing = 0;
316
0
    auto futureX =
317
0
        std::async(std::launch::async, [&poVarX, &dfXStart, &dfXSpacing]()
318
0
                   { return poVarX->IsRegularlySpaced(dfXStart, dfXSpacing); });
319
0
    CPL_IGNORE_RET_VAL(poVarY->IsRegularlySpaced(dfYStart, dfYSpacing));
320
0
    CPL_IGNORE_RET_VAL(futureX.get());
321
0
}
322
323
/************************************************************************/
324
/*                                Open()                                */
325
/************************************************************************/
326
327
GDALDataset *ZarrDataset::Open(GDALOpenInfo *poOpenInfo)
328
274
{
329
274
    if (!ZARRDriverIdentify(poOpenInfo))
330
0
    {
331
0
        return nullptr;
332
0
    }
333
334
    // Used by gdal_translate kerchunk_ref.json kerchunk_parq.parq -of ZARR -co CONVERT_TO_KERCHUNK_PARQUET_REFERENCE=YES
335
274
    if (STARTS_WITH(poOpenInfo->pszFilename, "ZARR_DUMMY:"))
336
0
    {
337
0
        class ZarrDummyDataset final : public GDALDataset
338
0
        {
339
0
          public:
340
0
            ZarrDummyDataset()
341
0
            {
342
0
                nRasterXSize = 0;
343
0
                nRasterYSize = 0;
344
0
            }
345
0
        };
346
347
0
        auto poDS = std::make_unique<ZarrDummyDataset>();
348
0
        poDS->SetDescription(poOpenInfo->pszFilename + strlen("ZARR_DUMMY:"));
349
0
        return poDS.release();
350
0
    }
351
352
274
    const bool bKerchunkCached = CPLFetchBool(poOpenInfo->papszOpenOptions,
353
274
                                              "CACHE_KERCHUNK_JSON", false);
354
355
274
    if (ZARRIsLikelyKerchunkJSONRef(poOpenInfo))
356
0
    {
357
0
        GDALOpenInfo oOpenInfo(std::string("ZARR:\"")
358
0
                                   .append(bKerchunkCached
359
0
                                               ? JSON_REF_CACHED_FS_PREFIX
360
0
                                               : JSON_REF_FS_PREFIX)
361
0
                                   .append("{")
362
0
                                   .append(poOpenInfo->pszFilename)
363
0
                                   .append("}\"")
364
0
                                   .c_str(),
365
0
                               GA_ReadOnly);
366
0
        oOpenInfo.nOpenFlags = poOpenInfo->nOpenFlags;
367
0
        oOpenInfo.papszOpenOptions = poOpenInfo->papszOpenOptions;
368
0
        return Open(&oOpenInfo);
369
0
    }
370
274
    else if (STARTS_WITH(poOpenInfo->pszFilename, JSON_REF_FS_PREFIX) ||
371
273
             STARTS_WITH(poOpenInfo->pszFilename, JSON_REF_CACHED_FS_PREFIX))
372
1
    {
373
1
        GDALOpenInfo oOpenInfo(
374
1
            std::string("ZARR:").append(poOpenInfo->pszFilename).c_str(),
375
1
            GA_ReadOnly);
376
1
        oOpenInfo.nOpenFlags = poOpenInfo->nOpenFlags;
377
1
        oOpenInfo.papszOpenOptions = poOpenInfo->papszOpenOptions;
378
1
        return Open(&oOpenInfo);
379
1
    }
380
381
273
    CPLString osFilename(poOpenInfo->pszFilename);
382
273
    if (!poOpenInfo->bIsDirectory)
383
273
    {
384
273
        osFilename = CPLGetPathSafe(osFilename);
385
273
    }
386
273
    CPLString osArrayOfInterest;
387
273
    std::vector<uint64_t> anExtraDimIndices;
388
273
    if (STARTS_WITH(poOpenInfo->pszFilename, "ZARR:"))
389
74
    {
390
74
        const CPLStringList aosTokens(CSLTokenizeString2(
391
74
            poOpenInfo->pszFilename, ":", CSLT_HONOURSTRINGS));
392
74
        if (aosTokens.size() < 2)
393
0
            return nullptr;
394
74
        osFilename = aosTokens[1];
395
396
74
        if (!cpl::starts_with(osFilename, JSON_REF_FS_PREFIX) &&
397
73
            !cpl::starts_with(osFilename, JSON_REF_CACHED_FS_PREFIX) &&
398
73
            CPLGetExtensionSafe(osFilename) == "json")
399
0
        {
400
0
            VSIStatBufL sStat;
401
0
            if (VSIStatL(osFilename.c_str(), &sStat) == 0 &&
402
0
                !VSI_ISDIR(sStat.st_mode))
403
0
            {
404
0
                osFilename =
405
0
                    std::string(bKerchunkCached ? JSON_REF_CACHED_FS_PREFIX
406
0
                                                : JSON_REF_FS_PREFIX)
407
0
                        .append("{")
408
0
                        .append(osFilename)
409
0
                        .append("}");
410
0
            }
411
0
        }
412
413
74
        std::string osErrorMsg;
414
74
        if (osFilename == "http" || osFilename == "https")
415
0
        {
416
0
            osErrorMsg = "There is likely a quoting error of the whole "
417
0
                         "connection string, and the filename should "
418
0
                         "likely be prefixed with /vsicurl/";
419
0
        }
420
74
        else if (osFilename == "/vsicurl/http" ||
421
74
                 osFilename == "/vsicurl/https")
422
0
        {
423
0
            osErrorMsg = "There is likely a quoting error of the whole "
424
0
                         "connection string.";
425
0
        }
426
74
        else if (STARTS_WITH(osFilename.c_str(), "http://") ||
427
74
                 STARTS_WITH(osFilename.c_str(), "https://"))
428
0
        {
429
0
            osErrorMsg =
430
0
                "The filename should likely be prefixed with /vsicurl/";
431
0
        }
432
74
        if (!osErrorMsg.empty())
433
0
        {
434
0
            CPLError(CE_Failure, CPLE_AppDefined, "%s", osErrorMsg.c_str());
435
0
            return nullptr;
436
0
        }
437
74
        if (aosTokens.size() >= 3)
438
31
        {
439
31
            osArrayOfInterest = aosTokens[2];
440
136
            for (int i = 3; i < aosTokens.size(); ++i)
441
105
            {
442
105
                anExtraDimIndices.push_back(
443
105
                    static_cast<uint64_t>(CPLAtoGIntBig(aosTokens[i])));
444
105
            }
445
31
        }
446
74
    }
447
448
273
    auto poDSMultiDim = std::unique_ptr<GDALDataset>(
449
273
        OpenMultidim(osFilename.c_str(), poOpenInfo->eAccess == GA_Update,
450
273
                     poOpenInfo->papszOpenOptions));
451
273
    if (poDSMultiDim == nullptr ||
452
0
        (poOpenInfo->nOpenFlags & GDAL_OF_MULTIDIM_RASTER) != 0)
453
273
    {
454
273
        return poDSMultiDim.release();
455
273
    }
456
457
0
    auto poRG = poDSMultiDim->GetRootGroup();
458
459
0
    auto poDS = std::make_unique<ZarrDataset>(nullptr);
460
0
    std::shared_ptr<GDALMDArray> poMainArray;
461
0
    std::vector<std::string> aosArrays;
462
0
    std::string osMainArray;
463
0
    const bool bMultiband = CPLTestBool(
464
0
        CSLFetchNameValueDef(poOpenInfo->papszOpenOptions, "MULTIBAND", "YES"));
465
0
    size_t iXDim = 0;
466
0
    size_t iYDim = 0;
467
468
0
    if (!osArrayOfInterest.empty())
469
0
    {
470
0
        poMainArray = osArrayOfInterest == "/"
471
0
                          ? poRG->OpenMDArray("/")
472
0
                          : poRG->OpenMDArrayFromFullname(osArrayOfInterest);
473
0
        if (poMainArray == nullptr)
474
0
            return nullptr;
475
0
        GetXYDimensionIndices(poMainArray, poOpenInfo, iXDim, iYDim);
476
477
0
        if (poMainArray->GetDimensionCount() > 2)
478
0
        {
479
0
            if (anExtraDimIndices.empty())
480
0
            {
481
0
                const uint64_t nExtraDimSamples =
482
0
                    GetExtraDimSampleCount(poMainArray, iXDim, iYDim);
483
0
                if (bMultiband)
484
0
                {
485
0
                    if (nExtraDimSamples > 65536)  // arbitrary limit
486
0
                    {
487
0
                        if (poMainArray->GetDimensionCount() == 3)
488
0
                        {
489
0
                            CPLError(CE_Warning, CPLE_AppDefined,
490
0
                                     "Too many samples along the > 2D "
491
0
                                     "dimensions of %s. "
492
0
                                     "Use ZARR:\"%s\":%s:{i} syntax",
493
0
                                     osArrayOfInterest.c_str(),
494
0
                                     osFilename.c_str(),
495
0
                                     osArrayOfInterest.c_str());
496
0
                        }
497
0
                        else
498
0
                        {
499
0
                            CPLError(CE_Warning, CPLE_AppDefined,
500
0
                                     "Too many samples along the > 2D "
501
0
                                     "dimensions of %s. "
502
0
                                     "Use ZARR:\"%s\":%s:{i}:{j} syntax",
503
0
                                     osArrayOfInterest.c_str(),
504
0
                                     osFilename.c_str(),
505
0
                                     osArrayOfInterest.c_str());
506
0
                        }
507
0
                        return nullptr;
508
0
                    }
509
0
                }
510
0
                else if (nExtraDimSamples != 1)
511
0
                {
512
0
                    CPLError(CE_Failure, CPLE_AppDefined,
513
0
                             "Indices of extra dimensions must be specified");
514
0
                    return nullptr;
515
0
                }
516
0
            }
517
0
            else if (anExtraDimIndices.size() !=
518
0
                     poMainArray->GetDimensionCount() - 2)
519
0
            {
520
0
                CPLError(CE_Failure, CPLE_AppDefined,
521
0
                         "Wrong number of indices of extra dimensions");
522
0
                return nullptr;
523
0
            }
524
0
            else
525
0
            {
526
0
                for (const auto idx : anExtraDimIndices)
527
0
                {
528
0
                    poMainArray = poMainArray->at(idx);
529
0
                    if (poMainArray == nullptr)
530
0
                        return nullptr;
531
0
                }
532
0
                GetXYDimensionIndices(poMainArray, poOpenInfo, iXDim, iYDim);
533
0
            }
534
0
        }
535
0
        else if (!anExtraDimIndices.empty())
536
0
        {
537
0
            CPLError(CE_Failure, CPLE_AppDefined, "Unexpected extra indices");
538
0
            return nullptr;
539
0
        }
540
0
    }
541
0
    else
542
0
    {
543
0
        ExploreGroup(poRG, aosArrays, 0);
544
0
        if (aosArrays.empty())
545
0
            return nullptr;
546
547
0
        const bool bListAllArrays = CPLTestBool(CSLFetchNameValueDef(
548
0
            poOpenInfo->papszOpenOptions, "LIST_ALL_ARRAYS", "NO"));
549
550
0
        if (!bListAllArrays)
551
0
        {
552
0
            if (aosArrays.size() == 1)
553
0
            {
554
0
                poMainArray = poRG->OpenMDArrayFromFullname(aosArrays[0]);
555
0
                if (poMainArray)
556
0
                    osMainArray = poMainArray->GetFullName();
557
0
            }
558
0
            else  // at least 2 arrays
559
0
            {
560
0
                for (const auto &osArrayName : aosArrays)
561
0
                {
562
0
                    auto poArray = poRG->OpenMDArrayFromFullname(osArrayName);
563
0
                    if (poArray && poArray->GetDimensionCount() >= 2 &&
564
0
                        osArrayName.find("/ovr_") == std::string::npos)
565
0
                    {
566
0
                        if (osMainArray.empty())
567
0
                        {
568
0
                            poMainArray = std::move(poArray);
569
0
                            osMainArray = osArrayName;
570
0
                        }
571
0
                        else
572
0
                        {
573
0
                            poMainArray.reset();
574
0
                            osMainArray.clear();
575
0
                            break;
576
0
                        }
577
0
                    }
578
0
                }
579
0
            }
580
581
0
            if (poMainArray)
582
0
                GetXYDimensionIndices(poMainArray, poOpenInfo, iXDim, iYDim);
583
0
        }
584
585
0
        int iCountSubDS = 1;
586
587
0
        if (poMainArray && poMainArray->GetDimensionCount() > 2)
588
0
        {
589
0
            const auto &apoDims = poMainArray->GetDimensions();
590
0
            const uint64_t nExtraDimSamples =
591
0
                GetExtraDimSampleCount(poMainArray, iXDim, iYDim);
592
0
            if (nExtraDimSamples > 65536)  // arbitrary limit
593
0
            {
594
0
                if (apoDims.size() == 3)
595
0
                {
596
0
                    CPLError(
597
0
                        CE_Warning, CPLE_AppDefined,
598
0
                        "Too many samples along the > 2D dimensions of %s. "
599
0
                        "Use ZARR:\"%s\":%s:{i} syntax",
600
0
                        osMainArray.c_str(), osFilename.c_str(),
601
0
                        osMainArray.c_str());
602
0
                }
603
0
                else
604
0
                {
605
0
                    CPLError(
606
0
                        CE_Warning, CPLE_AppDefined,
607
0
                        "Too many samples along the > 2D dimensions of %s. "
608
0
                        "Use ZARR:\"%s\":%s:{i}:{j} syntax",
609
0
                        osMainArray.c_str(), osFilename.c_str(),
610
0
                        osMainArray.c_str());
611
0
                }
612
0
            }
613
0
            else if (nExtraDimSamples > 1 && bMultiband)
614
0
            {
615
                // nothing to do
616
0
            }
617
0
            else if (nExtraDimSamples > 1 && apoDims.size() == 3)
618
0
            {
619
0
                for (int i = 0; i < static_cast<int>(nExtraDimSamples); ++i)
620
0
                {
621
0
                    poDS->m_aosSubdatasets.AddString(CPLSPrintf(
622
0
                        "SUBDATASET_%d_NAME=ZARR:\"%s\":%s:%d", iCountSubDS,
623
0
                        osFilename.c_str(), osMainArray.c_str(), i));
624
0
                    poDS->m_aosSubdatasets.AddString(CPLSPrintf(
625
0
                        "SUBDATASET_%d_DESC=Array %s at index %d of %s",
626
0
                        iCountSubDS, osMainArray.c_str(), i,
627
0
                        apoDims[0]->GetName().c_str()));
628
0
                    ++iCountSubDS;
629
0
                }
630
0
            }
631
0
            else if (nExtraDimSamples > 1)
632
0
            {
633
0
                int nDimIdxI = 0;
634
0
                int nDimIdxJ = 0;
635
0
                for (int i = 0; i < static_cast<int>(nExtraDimSamples); ++i)
636
0
                {
637
0
                    poDS->m_aosSubdatasets.AddString(
638
0
                        CPLSPrintf("SUBDATASET_%d_NAME=ZARR:\"%s\":%s:%d:%d",
639
0
                                   iCountSubDS, osFilename.c_str(),
640
0
                                   osMainArray.c_str(), nDimIdxI, nDimIdxJ));
641
0
                    poDS->m_aosSubdatasets.AddString(
642
0
                        CPLSPrintf("SUBDATASET_%d_DESC=Array %s at "
643
0
                                   "index %d of %s and %d of %s",
644
0
                                   iCountSubDS, osMainArray.c_str(), nDimIdxI,
645
0
                                   apoDims[0]->GetName().c_str(), nDimIdxJ,
646
0
                                   apoDims[1]->GetName().c_str()));
647
0
                    ++iCountSubDS;
648
0
                    ++nDimIdxJ;
649
0
                    if (nDimIdxJ == static_cast<int>(apoDims[1]->GetSize()))
650
0
                    {
651
0
                        nDimIdxJ = 0;
652
0
                        ++nDimIdxI;
653
0
                    }
654
0
                }
655
0
            }
656
0
        }
657
658
0
        if (bListAllArrays || aosArrays.size() >= 2)
659
0
        {
660
0
            for (size_t i = 0; i < aosArrays.size(); ++i)
661
0
            {
662
0
                auto poArray = poRG->OpenMDArrayFromFullname(aosArrays[i]);
663
0
                if (poArray && (bListAllArrays || aosArrays[i].find("/ovr_") ==
664
0
                                                      std::string::npos))
665
0
                {
666
0
                    bool bAddSubDS = false;
667
0
                    if (bListAllArrays)
668
0
                    {
669
0
                        bAddSubDS = true;
670
0
                    }
671
0
                    else if (poArray->GetDimensionCount() >= 2)
672
0
                    {
673
0
                        bAddSubDS = true;
674
0
                    }
675
0
                    if (bAddSubDS)
676
0
                    {
677
0
                        std::string osDim;
678
0
                        const auto &apoDims = poArray->GetDimensions();
679
0
                        for (const auto &poDim : apoDims)
680
0
                        {
681
0
                            if (!osDim.empty())
682
0
                                osDim += "x";
683
0
                            osDim += CPLSPrintf(
684
0
                                "%" PRIu64,
685
0
                                static_cast<uint64_t>(poDim->GetSize()));
686
0
                        }
687
688
0
                        std::string osDataType;
689
0
                        if (poArray->GetDataType().GetClass() == GEDTC_STRING)
690
0
                        {
691
0
                            osDataType = "string type";
692
0
                        }
693
0
                        else if (poArray->GetDataType().GetClass() ==
694
0
                                 GEDTC_NUMERIC)
695
0
                        {
696
0
                            osDataType = GDALGetDataTypeName(
697
0
                                poArray->GetDataType().GetNumericDataType());
698
0
                        }
699
0
                        else
700
0
                        {
701
0
                            osDataType = "compound type";
702
0
                        }
703
704
0
                        poDS->m_aosSubdatasets.AddString(CPLSPrintf(
705
0
                            "SUBDATASET_%d_NAME=ZARR:\"%s\":%s", iCountSubDS,
706
0
                            osFilename.c_str(), aosArrays[i].c_str()));
707
0
                        poDS->m_aosSubdatasets.AddString(CPLSPrintf(
708
0
                            "SUBDATASET_%d_DESC=[%s] %s (%s)", iCountSubDS,
709
0
                            osDim.c_str(), aosArrays[i].c_str(),
710
0
                            osDataType.c_str()));
711
0
                        ++iCountSubDS;
712
0
                    }
713
0
                }
714
0
            }
715
0
        }
716
0
    }
717
718
0
    if (poMainArray && (bMultiband || poMainArray->GetDimensionCount() <= 2))
719
0
    {
720
0
        PrefetchCoordArrays(poMainArray, iXDim, iYDim);
721
722
        // Pass papszOpenOptions for LOAD_EXTRA_DIM_METADATA_DELAY
723
0
        auto poNewDS =
724
0
            std::unique_ptr<GDALDataset>(poMainArray->AsClassicDataset(
725
0
                iXDim, iYDim, poRG, poOpenInfo->papszOpenOptions));
726
0
        if (!poNewDS)
727
0
            return nullptr;
728
729
0
        if (poMainArray->GetDimensionCount() >= 2)
730
0
        {
731
            // If we have 3 arrays, check that the 2 ones that are not the main
732
            // 2D array are indexing variables of its dimensions. If so, don't
733
            // expose them as subdatasets
734
0
            if (aosArrays.size() == 3)
735
0
            {
736
0
                std::vector<std::string> aosOtherArrays;
737
0
                for (size_t i = 0; i < aosArrays.size(); ++i)
738
0
                {
739
0
                    if (aosArrays[i] != osMainArray)
740
0
                    {
741
0
                        aosOtherArrays.emplace_back(aosArrays[i]);
742
0
                    }
743
0
                }
744
0
                bool bMatchFound[] = {false, false};
745
0
                for (int i = 0; i < 2; i++)
746
0
                {
747
0
                    auto poIndexingVar =
748
0
                        poMainArray->GetDimensions()[i == 0 ? iXDim : iYDim]
749
0
                            ->GetIndexingVariable();
750
0
                    if (poIndexingVar)
751
0
                    {
752
0
                        for (int j = 0; j < 2; j++)
753
0
                        {
754
0
                            if (aosOtherArrays[j] ==
755
0
                                poIndexingVar->GetFullName())
756
0
                            {
757
0
                                bMatchFound[i] = true;
758
0
                                break;
759
0
                            }
760
0
                        }
761
0
                    }
762
0
                }
763
0
                if (bMatchFound[0] && bMatchFound[1])
764
0
                {
765
0
                    poDS->m_aosSubdatasets.Clear();
766
0
                }
767
0
            }
768
0
        }
769
0
        if (!poDS->m_aosSubdatasets.empty())
770
0
        {
771
0
            poNewDS->SetMetadata(poDS->m_aosSubdatasets.List(),
772
0
                                 GDAL_MDD_SUBDATASETS);
773
0
        }
774
0
        return poNewDS.release();
775
0
    }
776
777
0
    return poDS.release();
778
0
}
779
780
/************************************************************************/
781
/*                         ZarrDatasetDelete()                          */
782
/************************************************************************/
783
784
static CPLErr ZarrDatasetDelete(const char *pszFilename)
785
27
{
786
27
    if (STARTS_WITH(pszFilename, "ZARR:"))
787
0
    {
788
0
        CPLError(CE_Failure, CPLE_AppDefined,
789
0
                 "Delete() only supported on ZARR connection names "
790
0
                 "not starting with the ZARR: prefix");
791
0
        return CE_Failure;
792
0
    }
793
27
    return VSIRmdirRecursive(pszFilename) == 0 ? CE_None : CE_Failure;
794
27
}
795
796
/************************************************************************/
797
/*                         ZarrDatasetRename()                          */
798
/************************************************************************/
799
800
static CPLErr ZarrDatasetRename(const char *pszNewName, const char *pszOldName)
801
0
{
802
0
    if (STARTS_WITH(pszNewName, "ZARR:") || STARTS_WITH(pszOldName, "ZARR:"))
803
0
    {
804
0
        CPLError(CE_Failure, CPLE_AppDefined,
805
0
                 "Rename() only supported on ZARR connection names "
806
0
                 "not starting with the ZARR: prefix");
807
0
        return CE_Failure;
808
0
    }
809
0
    return VSIRename(pszOldName, pszNewName) == 0 ? CE_None : CE_Failure;
810
0
}
811
812
/************************************************************************/
813
/*                        ZarrDatasetCopyFiles()                        */
814
/************************************************************************/
815
816
static CPLErr ZarrDatasetCopyFiles(const char *pszNewName,
817
                                   const char *pszOldName)
818
0
{
819
0
    if (STARTS_WITH(pszNewName, "ZARR:") || STARTS_WITH(pszOldName, "ZARR:"))
820
0
    {
821
0
        CPLError(CE_Failure, CPLE_AppDefined,
822
0
                 "CopyFiles() only supported on ZARR connection names "
823
0
                 "not starting with the ZARR: prefix");
824
0
        return CE_Failure;
825
0
    }
826
    // VSISync() returns true in case of success
827
0
    return VSISync((std::string(pszOldName) + '/').c_str(), pszNewName, nullptr,
828
0
                   nullptr, nullptr, nullptr)
829
0
               ? CE_None
830
0
               : CE_Failure;
831
0
}
832
833
/************************************************************************/
834
/*                       ZarrDriverClearCaches()                        */
835
/************************************************************************/
836
837
static void ZarrDriverClearCaches(GDALDriver *)
838
0
{
839
0
    ZarrClearCoordinateCache();
840
0
    ZarrClearShardIndexCache();
841
0
}
842
843
/************************************************************************/
844
/*                             ZarrDriver()                             */
845
/************************************************************************/
846
847
class ZarrDriver final : public GDALDriver
848
{
849
    std::recursive_mutex m_oMutex{};
850
    bool m_bMetadataInitialized = false;
851
    void InitMetadata();
852
853
  public:
854
    const char *GetMetadataItem(const char *pszName,
855
                                const char *pszDomain) override;
856
857
    CSLConstList GetMetadata(const char *pszDomain) override
858
79
    {
859
79
        std::lock_guard oLock(m_oMutex);
860
79
        InitMetadata();
861
79
        return GDALDriver::GetMetadata(pszDomain);
862
79
    }
863
};
864
865
const char *ZarrDriver::GetMetadataItem(const char *pszName,
866
                                        const char *pszDomain)
867
1.11M
{
868
1.11M
    std::lock_guard oLock(m_oMutex);
869
1.11M
    if (EQUAL(pszName, "COMPRESSORS") || EQUAL(pszName, "BLOSC_COMPRESSORS") ||
870
1.11M
        EQUAL(pszName, GDAL_DMD_CREATIONOPTIONLIST) ||
871
1.11M
        EQUAL(pszName, GDAL_DMD_MULTIDIM_ARRAY_CREATIONOPTIONLIST))
872
184
    {
873
184
        InitMetadata();
874
184
    }
875
1.11M
    return GDALDriver::GetMetadataItem(pszName, pszDomain);
876
1.11M
}
877
878
void ZarrDriver::InitMetadata()
879
263
{
880
263
    if (m_bMetadataInitialized)
881
262
        return;
882
1
    m_bMetadataInitialized = true;
883
884
1
    {
885
        // A bit of a hack. Normally GetMetadata() should also return it,
886
        // but as this is only used for tests, just make GetMetadataItem()
887
        // handle it
888
1
        std::string osCompressors;
889
1
        std::string osFilters;
890
1
        char **decompressors = CPLGetDecompressors();
891
6
        for (auto iter = decompressors; iter && *iter; ++iter)
892
5
        {
893
5
            const auto psCompressor = CPLGetCompressor(*iter);
894
5
            if (psCompressor)
895
5
            {
896
5
                if (psCompressor->eType == CCT_COMPRESSOR)
897
4
                {
898
4
                    if (!osCompressors.empty())
899
3
                        osCompressors += ',';
900
4
                    osCompressors += *iter;
901
4
                }
902
1
                else if (psCompressor->eType == CCT_FILTER)
903
1
                {
904
1
                    if (!osFilters.empty())
905
0
                        osFilters += ',';
906
1
                    osFilters += *iter;
907
1
                }
908
5
            }
909
5
        }
910
1
        CSLDestroy(decompressors);
911
1
        GDALDriver::SetMetadataItem("COMPRESSORS", osCompressors.c_str());
912
1
        GDALDriver::SetMetadataItem("FILTERS", osFilters.c_str());
913
1
    }
914
#ifdef HAVE_BLOSC
915
    {
916
        GDALDriver::SetMetadataItem("BLOSC_COMPRESSORS",
917
                                    blosc_list_compressors());
918
    }
919
#endif
920
921
1
    {
922
1
        CPLXMLTreeCloser oTree(
923
1
            CPLCreateXMLNode(nullptr, CXT_Element, "CreationOptionList"));
924
1
        char **compressors = CPLGetCompressors();
925
926
1
        auto psCompressNode =
927
1
            CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
928
1
        CPLAddXMLAttributeAndValue(psCompressNode, "name", "COMPRESS");
929
1
        CPLAddXMLAttributeAndValue(psCompressNode, "type", "string-select");
930
1
        CPLAddXMLAttributeAndValue(psCompressNode, "description",
931
1
                                   "Compression method");
932
1
        CPLAddXMLAttributeAndValue(psCompressNode, "default", "NONE");
933
1
        {
934
1
            auto poValueNode =
935
1
                CPLCreateXMLNode(psCompressNode, CXT_Element, "Value");
936
1
            CPLCreateXMLNode(poValueNode, CXT_Text, "NONE");
937
1
        }
938
939
1
        auto psFilterNode =
940
1
            CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
941
1
        CPLAddXMLAttributeAndValue(psFilterNode, "name", "FILTER");
942
1
        CPLAddXMLAttributeAndValue(psFilterNode, "type", "string-select");
943
1
        CPLAddXMLAttributeAndValue(psFilterNode, "description",
944
1
                                   "Filter method (only for ZARR_V2)");
945
1
        CPLAddXMLAttributeAndValue(psFilterNode, "default", "NONE");
946
1
        {
947
1
            auto poValueNode =
948
1
                CPLCreateXMLNode(psFilterNode, CXT_Element, "Value");
949
1
            CPLCreateXMLNode(poValueNode, CXT_Text, "NONE");
950
1
        }
951
952
1
        auto psBlockSizeNode =
953
1
            CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
954
1
        CPLAddXMLAttributeAndValue(psBlockSizeNode, "name", "BLOCKSIZE");
955
1
        CPLAddXMLAttributeAndValue(psBlockSizeNode, "type", "string");
956
1
        CPLAddXMLAttributeAndValue(
957
1
            psBlockSizeNode, "description",
958
1
            "Comma separated list of chunk size along each dimension");
959
960
1
        auto psChunkMemoryLayout =
961
1
            CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
962
1
        CPLAddXMLAttributeAndValue(psChunkMemoryLayout, "name",
963
1
                                   "CHUNK_MEMORY_LAYOUT");
964
1
        CPLAddXMLAttributeAndValue(psChunkMemoryLayout, "type",
965
1
                                   "string-select");
966
1
        CPLAddXMLAttributeAndValue(psChunkMemoryLayout, "description",
967
1
                                   "Whether to use C (row-major) order or F "
968
1
                                   "(column-major) order in chunks");
969
1
        CPLAddXMLAttributeAndValue(psChunkMemoryLayout, "default", "C");
970
1
        {
971
1
            auto poValueNode =
972
1
                CPLCreateXMLNode(psChunkMemoryLayout, CXT_Element, "Value");
973
1
            CPLCreateXMLNode(poValueNode, CXT_Text, "C");
974
1
        }
975
1
        {
976
1
            auto poValueNode =
977
1
                CPLCreateXMLNode(psChunkMemoryLayout, CXT_Element, "Value");
978
1
            CPLCreateXMLNode(poValueNode, CXT_Text, "F");
979
1
        }
980
981
1
        auto psStringFormat =
982
1
            CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
983
1
        CPLAddXMLAttributeAndValue(psStringFormat, "name", "STRING_FORMAT");
984
1
        CPLAddXMLAttributeAndValue(psStringFormat, "type", "string-select");
985
1
        CPLAddXMLAttributeAndValue(psStringFormat, "default", "STRING");
986
1
        {
987
1
            auto poValueNode =
988
1
                CPLCreateXMLNode(psStringFormat, CXT_Element, "Value");
989
1
            CPLCreateXMLNode(poValueNode, CXT_Text, "STRING");
990
1
        }
991
1
        {
992
1
            auto poValueNode =
993
1
                CPLCreateXMLNode(psStringFormat, CXT_Element, "Value");
994
1
            CPLCreateXMLNode(poValueNode, CXT_Text, "UNICODE");
995
1
        }
996
997
1
        auto psDimSeparatorNode =
998
1
            CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
999
1
        CPLAddXMLAttributeAndValue(psDimSeparatorNode, "name", "DIM_SEPARATOR");
1000
1
        CPLAddXMLAttributeAndValue(psDimSeparatorNode, "type", "string");
1001
1
        CPLAddXMLAttributeAndValue(
1002
1
            psDimSeparatorNode, "description",
1003
1
            "Dimension separator in chunk filenames. Default to decimal point "
1004
1
            "for ZarrV2 and slash for ZarrV3");
1005
1006
6
        for (auto iter = compressors; iter && *iter; ++iter)
1007
5
        {
1008
5
            const auto psCompressor = CPLGetCompressor(*iter);
1009
5
            if (psCompressor)
1010
5
            {
1011
5
                auto poValueNode = CPLCreateXMLNode(
1012
5
                    (psCompressor->eType == CCT_COMPRESSOR) ? psCompressNode
1013
5
                                                            : psFilterNode,
1014
5
                    CXT_Element, "Value");
1015
5
                CPLCreateXMLNode(poValueNode, CXT_Text,
1016
5
                                 CPLString(*iter).toupper().c_str());
1017
1018
5
                const char *pszOptions =
1019
5
                    CSLFetchNameValue(psCompressor->papszMetadata, "OPTIONS");
1020
5
                if (pszOptions)
1021
5
                {
1022
5
                    CPLXMLTreeCloser oTreeCompressor(
1023
5
                        CPLParseXMLString(pszOptions));
1024
5
                    const auto psRoot =
1025
5
                        oTreeCompressor.get()
1026
5
                            ? CPLGetXMLNode(oTreeCompressor.get(), "=Options")
1027
5
                            : nullptr;
1028
5
                    if (psRoot)
1029
5
                    {
1030
5
                        for (CPLXMLNode *psNode = psRoot->psChild;
1031
12
                             psNode != nullptr; psNode = psNode->psNext)
1032
7
                        {
1033
7
                            if (psNode->eType == CXT_Element)
1034
7
                            {
1035
7
                                const char *pszName =
1036
7
                                    CPLGetXMLValue(psNode, "name", nullptr);
1037
7
                                if (pszName &&
1038
7
                                    !EQUAL(pszName, "TYPESIZE")   // Blosc
1039
7
                                    && !EQUAL(pszName, "HEADER")  // LZ4
1040
7
                                )
1041
7
                                {
1042
7
                                    CPLXMLNode *psNext = psNode->psNext;
1043
7
                                    psNode->psNext = nullptr;
1044
7
                                    CPLXMLNode *psOption =
1045
7
                                        CPLCloneXMLTree(psNode);
1046
1047
7
                                    CPLXMLNode *psName =
1048
7
                                        CPLGetXMLNode(psOption, "name");
1049
7
                                    if (psName &&
1050
7
                                        psName->eType == CXT_Attribute &&
1051
7
                                        psName->psChild &&
1052
7
                                        psName->psChild->pszValue)
1053
7
                                    {
1054
7
                                        CPLString osNewValue(*iter);
1055
7
                                        osNewValue = osNewValue.toupper();
1056
7
                                        osNewValue += '_';
1057
7
                                        osNewValue += psName->psChild->pszValue;
1058
7
                                        CPLFree(psName->psChild->pszValue);
1059
7
                                        psName->psChild->pszValue =
1060
7
                                            CPLStrdup(osNewValue.c_str());
1061
7
                                    }
1062
1063
7
                                    CPLXMLNode *psDescription =
1064
7
                                        CPLGetXMLNode(psOption, "description");
1065
7
                                    if (psDescription &&
1066
7
                                        psDescription->eType == CXT_Attribute &&
1067
7
                                        psDescription->psChild &&
1068
7
                                        psDescription->psChild->pszValue)
1069
7
                                    {
1070
7
                                        std::string osNewValue(
1071
7
                                            psDescription->psChild->pszValue);
1072
7
                                        if (psCompressor->eType ==
1073
7
                                            CCT_COMPRESSOR)
1074
6
                                        {
1075
6
                                            osNewValue +=
1076
6
                                                ". Only used when COMPRESS=";
1077
6
                                        }
1078
1
                                        else
1079
1
                                        {
1080
1
                                            osNewValue +=
1081
1
                                                ". Only used when FILTER=";
1082
1
                                        }
1083
7
                                        osNewValue +=
1084
7
                                            CPLString(*iter).toupper();
1085
7
                                        CPLFree(
1086
7
                                            psDescription->psChild->pszValue);
1087
7
                                        psDescription->psChild->pszValue =
1088
7
                                            CPLStrdup(osNewValue.c_str());
1089
7
                                    }
1090
1091
7
                                    CPLAddXMLChild(oTree.get(), psOption);
1092
7
                                    psNode->psNext = psNext;
1093
7
                                }
1094
7
                            }
1095
7
                        }
1096
5
                    }
1097
5
                }
1098
5
            }
1099
5
        }
1100
1
        CSLDestroy(compressors);
1101
1102
1
        auto psGeoreferencingConvention =
1103
1
            CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
1104
1
        CPLAddXMLAttributeAndValue(psGeoreferencingConvention, "name",
1105
1
                                   "GEOREFERENCING_CONVENTION");
1106
1
        CPLAddXMLAttributeAndValue(psGeoreferencingConvention, "type",
1107
1
                                   "string-select");
1108
1
        CPLAddXMLAttributeAndValue(psGeoreferencingConvention, "default",
1109
1
                                   "GDAL");
1110
1
        CPLAddXMLAttributeAndValue(psGeoreferencingConvention, "description",
1111
1
                                   "Georeferencing convention to use");
1112
1113
1
        {
1114
1
            auto poValueNode = CPLCreateXMLNode(psGeoreferencingConvention,
1115
1
                                                CXT_Element, "Value");
1116
1
            CPLCreateXMLNode(poValueNode, CXT_Text, "GDAL");
1117
1
        }
1118
1
        {
1119
1
            auto poValueNode = CPLCreateXMLNode(psGeoreferencingConvention,
1120
1
                                                CXT_Element, "Value");
1121
1
            CPLCreateXMLNode(poValueNode, CXT_Text, "SPATIAL_PROJ");
1122
1
        }
1123
1124
1
        {
1125
1
            char *pszXML = CPLSerializeXMLTree(oTree.get());
1126
1
            GDALDriver::SetMetadataItem(
1127
1
                GDAL_DMD_MULTIDIM_ARRAY_CREATIONOPTIONLIST,
1128
1
                CPLString(pszXML)
1129
1
                    .replaceAll("CreationOptionList",
1130
1
                                "MultiDimArrayCreationOptionList")
1131
1
                    .c_str());
1132
1
            CPLFree(pszXML);
1133
1
        }
1134
1135
1
        {
1136
1
            auto psArrayNameOption =
1137
1
                CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
1138
1
            CPLAddXMLAttributeAndValue(psArrayNameOption, "name", "ARRAY_NAME");
1139
1
            CPLAddXMLAttributeAndValue(psArrayNameOption, "type", "string");
1140
1
            CPLAddXMLAttributeAndValue(
1141
1
                psArrayNameOption, "description",
1142
1
                "Array name. If not specified, deduced from the filename");
1143
1144
1
            auto psAppendSubDSOption =
1145
1
                CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
1146
1
            CPLAddXMLAttributeAndValue(psAppendSubDSOption, "name",
1147
1
                                       "APPEND_SUBDATASET");
1148
1
            CPLAddXMLAttributeAndValue(psAppendSubDSOption, "type", "boolean");
1149
1
            CPLAddXMLAttributeAndValue(psAppendSubDSOption, "description",
1150
1
                                       "Whether to append the new dataset to "
1151
1
                                       "an existing Zarr hierarchy");
1152
1
            CPLAddXMLAttributeAndValue(psAppendSubDSOption, "default", "NO");
1153
1154
1
            auto psFormat =
1155
1
                CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
1156
1
            CPLAddXMLAttributeAndValue(psFormat, "name", "FORMAT");
1157
1
            CPLAddXMLAttributeAndValue(psFormat, "type", "string-select");
1158
1
            CPLAddXMLAttributeAndValue(psFormat, "default", "ZARR_V3");
1159
1
            {
1160
1
                auto poValueNode =
1161
1
                    CPLCreateXMLNode(psFormat, CXT_Element, "Value");
1162
1
                CPLCreateXMLNode(poValueNode, CXT_Text, "ZARR_V2");
1163
1
            }
1164
1
            {
1165
1
                auto poValueNode =
1166
1
                    CPLCreateXMLNode(psFormat, CXT_Element, "Value");
1167
1
                CPLCreateXMLNode(poValueNode, CXT_Text, "ZARR_V3");
1168
1
            }
1169
1170
1
            auto psCreateZMetadata =
1171
1
                CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
1172
1
            CPLAddXMLAttributeAndValue(psCreateZMetadata, "name",
1173
1
                                       "CREATE_CONSOLIDATED_METADATA");
1174
1
            CPLAddXMLAttributeAndValue(psCreateZMetadata, "alias",
1175
1
                                       "CREATE_ZMETADATA");
1176
1
            CPLAddXMLAttributeAndValue(psCreateZMetadata, "type", "boolean");
1177
1
            CPLAddXMLAttributeAndValue(
1178
1
                psCreateZMetadata, "description",
1179
1
                "Whether to create consolidated metadata");
1180
1
            CPLAddXMLAttributeAndValue(psCreateZMetadata, "default", "YES");
1181
1182
1
            auto psSingleArrayNode =
1183
1
                CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
1184
1
            CPLAddXMLAttributeAndValue(psSingleArrayNode, "name",
1185
1
                                       "SINGLE_ARRAY");
1186
1
            CPLAddXMLAttributeAndValue(psSingleArrayNode, "type", "boolean");
1187
1
            CPLAddXMLAttributeAndValue(
1188
1
                psSingleArrayNode, "description",
1189
1
                "Whether to write a multi-band dataset as a single array, or "
1190
1
                "one array per band");
1191
1
            CPLAddXMLAttributeAndValue(psSingleArrayNode, "default", "YES");
1192
1193
1
            auto psInterleaveNode =
1194
1
                CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
1195
1
            CPLAddXMLAttributeAndValue(psInterleaveNode, "name",
1196
1
                                       GDALMD_INTERLEAVE);
1197
1
            CPLAddXMLAttributeAndValue(psInterleaveNode, "type",
1198
1
                                       "string-select");
1199
1
            CPLAddXMLAttributeAndValue(psInterleaveNode, "default", "BAND");
1200
1
            {
1201
1
                auto poValueNode =
1202
1
                    CPLCreateXMLNode(psInterleaveNode, CXT_Element, "Value");
1203
1
                CPLCreateXMLNode(poValueNode, CXT_Text, "BAND");
1204
1
            }
1205
1
            {
1206
1
                auto poValueNode =
1207
1
                    CPLCreateXMLNode(psInterleaveNode, CXT_Element, "Value");
1208
1
                CPLCreateXMLNode(poValueNode, CXT_Text, "PIXEL");
1209
1
            }
1210
1211
1
            auto psConvertToParquet =
1212
1
                CPLCreateXMLNode(oTree.get(), CXT_Element, "Option");
1213
1
            CPLAddXMLAttributeAndValue(psConvertToParquet, "name",
1214
1
                                       "CONVERT_TO_KERCHUNK_PARQUET_REFERENCE");
1215
1
            CPLAddXMLAttributeAndValue(psConvertToParquet, "type", "boolean");
1216
1
            CPLAddXMLAttributeAndValue(
1217
1
                psConvertToParquet, "description",
1218
1
                "Whether to convert a Kerchunk JSON reference store to a "
1219
1
                "Kerchunk Parquet reference store. (CreateCopy() only)");
1220
1221
1
            char *pszXML = CPLSerializeXMLTree(oTree.get());
1222
1
            GDALDriver::SetMetadataItem(GDAL_DMD_CREATIONOPTIONLIST, pszXML);
1223
1
            CPLFree(pszXML);
1224
1
        }
1225
1
    }
1226
1
}
1227
1228
/************************************************************************/
1229
/*                       CreateMultiDimensional()                       */
1230
/************************************************************************/
1231
1232
GDALDataset *
1233
ZarrDataset::CreateMultiDimensional(const char *pszFilename,
1234
                                    CSLConstList /*papszRootGroupOptions*/,
1235
                                    CSLConstList papszOptions)
1236
0
{
1237
0
    const char *pszFormat =
1238
0
        CSLFetchNameValueDef(papszOptions, "FORMAT", "ZARR_V3");
1239
0
    std::shared_ptr<ZarrGroupBase> poRG;
1240
0
    auto poSharedResource =
1241
0
        ZarrSharedResource::Create(pszFilename, /*bUpdatable=*/true);
1242
0
    const bool bCreateZMetadata = CPLTestBool(CSLFetchNameValueDef(
1243
0
        papszOptions, "CREATE_CONSOLIDATED_METADATA",
1244
0
        CSLFetchNameValueDef(papszOptions, "CREATE_ZMETADATA", "YES")));
1245
0
    if (bCreateZMetadata)
1246
0
    {
1247
0
        poSharedResource->EnableConsolidatedMetadata(
1248
0
            EQUAL(pszFormat, "ZARR_V3")
1249
0
                ? ZarrSharedResource::ConsolidatedMetadataKind::INTERNAL
1250
0
                : ZarrSharedResource::ConsolidatedMetadataKind::EXTERNAL);
1251
0
    }
1252
0
    if (EQUAL(pszFormat, "ZARR_V3"))
1253
0
    {
1254
0
        poRG = ZarrV3Group::CreateOnDisk(poSharedResource, std::string(), "/",
1255
0
                                         pszFilename);
1256
0
    }
1257
0
    else
1258
0
    {
1259
0
        poRG = ZarrV2Group::CreateOnDisk(poSharedResource, std::string(), "/",
1260
0
                                         pszFilename);
1261
0
    }
1262
0
    if (!poRG)
1263
0
        return nullptr;
1264
1265
0
    auto poDS = new ZarrDataset(poRG);
1266
0
    poDS->SetDescription(pszFilename);
1267
0
    return poDS;
1268
0
}
1269
1270
/************************************************************************/
1271
/*                               Create()                               */
1272
/************************************************************************/
1273
1274
GDALDataset *ZarrDataset::Create(const char *pszName, int nXSize, int nYSize,
1275
                                 int nBandsIn, GDALDataType eType,
1276
                                 CSLConstList papszOptions)
1277
77
{
1278
    // To avoid any issue with short-lived string that would be passed to us
1279
77
    const std::string osName = pszName;
1280
77
    pszName = osName.c_str();
1281
1282
77
    if (nBandsIn <= 0 || nXSize <= 0 || nYSize <= 0)
1283
0
    {
1284
0
        CPLError(CE_Failure, CPLE_NotSupported,
1285
0
                 "nBands, nXSize, nYSize should be > 0");
1286
0
        return nullptr;
1287
0
    }
1288
1289
77
    const bool bAppendSubDS = CPLTestBool(
1290
77
        CSLFetchNameValueDef(papszOptions, "APPEND_SUBDATASET", "NO"));
1291
77
    const char *pszArrayName = CSLFetchNameValue(papszOptions, "ARRAY_NAME");
1292
1293
77
    std::shared_ptr<ZarrGroupBase> poRG;
1294
77
    if (bAppendSubDS)
1295
1
    {
1296
1
        if (pszArrayName == nullptr)
1297
1
        {
1298
1
            CPLError(CE_Failure, CPLE_AppDefined,
1299
1
                     "ARRAY_NAME should be provided when "
1300
1
                     "APPEND_SUBDATASET is set to YES");
1301
1
            return nullptr;
1302
1
        }
1303
0
        auto poDS =
1304
0
            std::unique_ptr<GDALDataset>(OpenMultidim(pszName, true, nullptr));
1305
0
        if (poDS == nullptr)
1306
0
        {
1307
0
            CPLError(CE_Failure, CPLE_AppDefined, "Cannot open %s", pszName);
1308
0
            return nullptr;
1309
0
        }
1310
0
        poRG = std::dynamic_pointer_cast<ZarrGroupBase>(poDS->GetRootGroup());
1311
0
    }
1312
76
    else
1313
76
    {
1314
76
        VSIStatBufL sStat;
1315
76
        const bool bExists = VSIStatL(pszName, &sStat) == 0;
1316
76
        const char *pszObjType = nullptr;
1317
76
        if (bExists && !VSI_ISDIR(sStat.st_mode))
1318
0
            pszObjType = "File";
1319
76
        else if ((bExists /* && VSI_ISDIR(sStat.st_mode)*/) ||
1320
76
                 !CPLStringList(VSIReadDirEx(pszName, 1)).empty())
1321
0
            pszObjType = "Directory";
1322
76
        if (pszObjType)
1323
0
        {
1324
0
            CPLError(CE_Failure, CPLE_FileIO, "%s %s already exists.",
1325
0
                     pszObjType, pszName);
1326
0
            return nullptr;
1327
0
        }
1328
1329
76
        const char *pszFormat =
1330
76
            CSLFetchNameValueDef(papszOptions, "FORMAT", "ZARR_V3");
1331
76
        auto poSharedResource =
1332
76
            ZarrSharedResource::Create(pszName, /*bUpdatable=*/true);
1333
76
        const bool bCreateZMetadata = CPLTestBool(CSLFetchNameValueDef(
1334
76
            papszOptions, "CREATE_CONSOLIDATED_METADATA",
1335
76
            CSLFetchNameValueDef(papszOptions, "CREATE_ZMETADATA", "YES")));
1336
76
        if (bCreateZMetadata)
1337
76
        {
1338
76
            poSharedResource->EnableConsolidatedMetadata(
1339
76
                EQUAL(pszFormat, "ZARR_V3")
1340
76
                    ? ZarrSharedResource::ConsolidatedMetadataKind::INTERNAL
1341
76
                    : ZarrSharedResource::ConsolidatedMetadataKind::EXTERNAL);
1342
76
        }
1343
76
        if (EQUAL(pszFormat, "ZARR_V3"))
1344
76
        {
1345
76
            poRG = ZarrV3Group::CreateOnDisk(poSharedResource, std::string(),
1346
76
                                             "/", pszName);
1347
76
        }
1348
0
        else
1349
0
        {
1350
0
            poRG = ZarrV2Group::CreateOnDisk(poSharedResource, std::string(),
1351
0
                                             "/", pszName);
1352
0
        }
1353
76
        poSharedResource->SetRootGroup(poRG);
1354
76
    }
1355
76
    if (!poRG)
1356
0
        return nullptr;
1357
1358
76
    auto poDS = std::make_unique<ZarrDataset>(poRG);
1359
76
    poDS->SetDescription(pszName);
1360
76
    poDS->nRasterYSize = nYSize;
1361
76
    poDS->nRasterXSize = nXSize;
1362
76
    poDS->eAccess = GA_Update;
1363
1364
76
    const auto CleanupCreatedFiles =
1365
76
        [bAppendSubDS, pszName, pszArrayName, &poRG, &poDS]()
1366
76
    {
1367
        // Make sure all objects are released so that ZarrSharedResource
1368
        // is finalized and all files are serialized.
1369
0
        poRG.reset();
1370
0
        poDS.reset();
1371
1372
0
        if (bAppendSubDS)
1373
0
        {
1374
0
            VSIRmdir(
1375
0
                CPLFormFilenameSafe(pszName, pszArrayName, nullptr).c_str());
1376
0
        }
1377
0
        else
1378
0
        {
1379
            // Be a bit careful before wiping too much stuff...
1380
            // At most 5 files expected for ZARR_V2: .zgroup, .zmetadata,
1381
            // one (empty) subdir, . and ..
1382
            // and for ZARR_V3: zarr.json, one (empty) subdir, . and ..
1383
0
            const CPLStringList aosFiles(VSIReadDirEx(pszName, 6));
1384
0
            if (aosFiles.size() < 6)
1385
0
            {
1386
0
                for (const char *pszFile : aosFiles)
1387
0
                {
1388
0
                    if (pszArrayName && strcmp(pszFile, pszArrayName) == 0)
1389
0
                    {
1390
0
                        VSIRmdir(CPLFormFilenameSafe(pszName, pszFile, nullptr)
1391
0
                                     .c_str());
1392
0
                    }
1393
0
                    else if (!pszArrayName &&
1394
0
                             strcmp(pszFile,
1395
0
                                    CPLGetBasenameSafe(pszName).c_str()) == 0)
1396
0
                    {
1397
0
                        VSIRmdir(CPLFormFilenameSafe(pszName, pszFile, nullptr)
1398
0
                                     .c_str());
1399
0
                    }
1400
0
                    else if (strcmp(pszFile, ".zgroup") == 0 ||
1401
0
                             strcmp(pszFile, ".zmetadata") == 0 ||
1402
0
                             strcmp(pszFile, "zarr.json") == 0)
1403
0
                    {
1404
0
                        VSIUnlink(CPLFormFilenameSafe(pszName, pszFile, nullptr)
1405
0
                                      .c_str());
1406
0
                    }
1407
0
                }
1408
0
                VSIRmdir(pszName);
1409
0
            }
1410
0
        }
1411
0
    };
1412
1413
76
    std::string osDimXType, osDimYType;
1414
76
    if (CPLTestBool(
1415
76
            CSLFetchNameValueDef(papszOptions, "@HAS_GEOTRANSFORM", "NO")))
1416
37
    {
1417
37
        osDimXType = GDAL_DIM_TYPE_HORIZONTAL_X;
1418
37
        osDimYType = GDAL_DIM_TYPE_HORIZONTAL_Y;
1419
37
    }
1420
76
    poDS->m_bSpatialProjConvention = EQUAL(
1421
76
        CSLFetchNameValueDef(papszOptions, "GEOREFERENCING_CONVENTION", "GDAL"),
1422
76
        "SPATIAL_PROJ");
1423
1424
76
    if (bAppendSubDS)
1425
0
    {
1426
0
        auto aoDims = poRG->GetDimensions();
1427
0
        for (const auto &poDim : aoDims)
1428
0
        {
1429
0
            if (poDim->GetName() == "Y" &&
1430
0
                poDim->GetSize() == static_cast<GUInt64>(nYSize))
1431
0
            {
1432
0
                poDS->m_poDimY = poDim;
1433
0
            }
1434
0
            else if (poDim->GetName() == "X" &&
1435
0
                     poDim->GetSize() == static_cast<GUInt64>(nXSize))
1436
0
            {
1437
0
                poDS->m_poDimX = poDim;
1438
0
            }
1439
0
        }
1440
0
        if (poDS->m_poDimY == nullptr)
1441
0
        {
1442
0
            poDS->m_poDimY =
1443
0
                poRG->CreateDimension(std::string(pszArrayName) + "_Y",
1444
0
                                      osDimYType, std::string(), nYSize);
1445
0
        }
1446
0
        if (poDS->m_poDimX == nullptr)
1447
0
        {
1448
0
            poDS->m_poDimX =
1449
0
                poRG->CreateDimension(std::string(pszArrayName) + "_X",
1450
0
                                      osDimXType, std::string(), nXSize);
1451
0
        }
1452
0
    }
1453
76
    else
1454
76
    {
1455
76
        poDS->m_poDimY =
1456
76
            poRG->CreateDimension("Y", osDimYType, std::string(), nYSize);
1457
76
        poDS->m_poDimX =
1458
76
            poRG->CreateDimension("X", osDimXType, std::string(), nXSize);
1459
76
    }
1460
76
    if (poDS->m_poDimY == nullptr || poDS->m_poDimX == nullptr)
1461
0
    {
1462
0
        CleanupCreatedFiles();
1463
0
        return nullptr;
1464
0
    }
1465
1466
76
    const bool bSingleArray =
1467
76
        CPLTestBool(CSLFetchNameValueDef(papszOptions, "SINGLE_ARRAY", "YES"));
1468
76
    const bool bBandInterleave = EQUAL(
1469
76
        CSLFetchNameValueDef(papszOptions, GDALMD_INTERLEAVE, "BAND"), "BAND");
1470
76
    std::shared_ptr<GDALDimension> poBandDim(
1471
76
        (bSingleArray && nBandsIn > 1)
1472
76
            ? poRG->CreateDimension("Band", std::string(), std::string(),
1473
75
                                    nBandsIn)
1474
76
            : nullptr);
1475
1476
76
    const std::string osNonNullArrayName =
1477
76
        pszArrayName ? std::string(pszArrayName) : CPLGetBasenameSafe(pszName);
1478
76
    if (poBandDim)
1479
75
    {
1480
75
        std::vector<std::shared_ptr<GDALDimension>> apoDims;
1481
75
        if (bBandInterleave)
1482
55
        {
1483
55
            apoDims = std::vector<std::shared_ptr<GDALDimension>>{
1484
55
                poBandDim, poDS->m_poDimY, poDS->m_poDimX};
1485
55
        }
1486
20
        else
1487
20
        {
1488
20
            apoDims = std::vector<std::shared_ptr<GDALDimension>>{
1489
20
                poDS->m_poDimY, poDS->m_poDimX, poBandDim};
1490
20
        }
1491
75
        CPL_IGNORE_RET_VAL(poBandDim);
1492
75
        poDS->m_poSingleArray =
1493
75
            std::dynamic_pointer_cast<ZarrArray>(poRG->CreateMDArray(
1494
75
                osNonNullArrayName.c_str(), apoDims,
1495
75
                GDALExtendedDataType::Create(eType), papszOptions));
1496
75
        if (!poDS->m_poSingleArray)
1497
0
        {
1498
0
            CleanupCreatedFiles();
1499
0
            return nullptr;
1500
0
        }
1501
75
        poDS->SetMetadataItem(GDALMD_INTERLEAVE,
1502
75
                              bBandInterleave ? "BAND" : "PIXEL",
1503
75
                              GDAL_MDD_IMAGE_STRUCTURE);
1504
75
        if (bBandInterleave)
1505
55
        {
1506
55
            const char *pszBlockSize =
1507
55
                CSLFetchNameValue(papszOptions, "BLOCKSIZE");
1508
55
            if (pszBlockSize)
1509
0
            {
1510
0
                const CPLStringList aosTokens(
1511
0
                    CSLTokenizeString2(pszBlockSize, ",", 0));
1512
0
                if (aosTokens.size() == 3 && atoi(aosTokens[0]) == nBandsIn)
1513
0
                {
1514
                    // Actually expose as pixel interleaved
1515
0
                    poDS->SetMetadataItem(GDALMD_INTERLEAVE, "PIXEL",
1516
0
                                          GDAL_MDD_IMAGE_STRUCTURE);
1517
0
                }
1518
0
            }
1519
55
        }
1520
12.2k
        for (int i = 0; i < nBandsIn; i++)
1521
12.1k
        {
1522
12.1k
            const std::string viewDef =
1523
12.1k
                CPLSPrintf(bBandInterleave ? "[%d,::,::]" : "[::,::,%d]", i);
1524
12.1k
            auto poSlicedArray = poDS->m_poSingleArray->GetView(viewDef);
1525
12.1k
            poDS->SetBand(i + 1, std::make_unique<ZarrRasterBand>(poSlicedArray,
1526
12.1k
                                                                  viewDef));
1527
12.1k
        }
1528
75
    }
1529
1
    else
1530
1
    {
1531
1
        const auto apoDims = std::vector<std::shared_ptr<GDALDimension>>{
1532
1
            poDS->m_poDimY, poDS->m_poDimX};
1533
2
        for (int i = 0; i < nBandsIn; i++)
1534
1
        {
1535
1
            auto poArray = poRG->CreateMDArray(
1536
1
                nBandsIn == 1  ? osNonNullArrayName.c_str()
1537
1
                : pszArrayName ? CPLSPrintf("%s_band%d", pszArrayName, i + 1)
1538
0
                               : CPLSPrintf("Band%d", i + 1),
1539
1
                apoDims, GDALExtendedDataType::Create(eType), papszOptions);
1540
1
            if (poArray == nullptr)
1541
0
            {
1542
0
                CleanupCreatedFiles();
1543
0
                return nullptr;
1544
0
            }
1545
1
            poDS->SetBand(i + 1, std::make_unique<ZarrRasterBand>(poArray));
1546
1
        }
1547
1
    }
1548
1549
76
    return poDS.release();
1550
76
}
1551
1552
/************************************************************************/
1553
/*                            ~ZarrDataset()                            */
1554
/************************************************************************/
1555
1556
ZarrDataset::~ZarrDataset()
1557
76
{
1558
76
    ZarrDataset::FlushCache(true);
1559
76
}
1560
1561
/************************************************************************/
1562
/*                             FlushCache()                             */
1563
/************************************************************************/
1564
1565
CPLErr ZarrDataset::FlushCache(bool bAtClosing)
1566
176
{
1567
176
    CPLErr eErr = GDALDataset::FlushCache(bAtClosing);
1568
1569
176
    if (m_poSingleArray && !m_poSingleArray->Flush())
1570
0
    {
1571
0
        eErr = CE_Failure;
1572
0
    }
1573
1574
176
    if (bAtClosing && m_poSingleArray)
1575
75
    {
1576
75
        auto poFirstBand = cpl::down_cast<ZarrRasterBand *>(papoBands[0]);
1577
75
        if (poFirstBand->m_dfNoData.has_value())
1578
27
        {
1579
27
            bool bSameValue = true;
1580
3.57k
            for (int i = 1; bSameValue && i < nBands; ++i)
1581
3.54k
            {
1582
3.54k
                auto poBand = cpl::down_cast<ZarrRasterBand *>(papoBands[i]);
1583
3.54k
                bSameValue = poBand->m_dfNoData.has_value() &&
1584
3.53k
                             ((std::isnan(poFirstBand->m_dfNoData.value()) &&
1585
0
                               std::isnan(poBand->m_dfNoData.value())) ||
1586
3.53k
                              (poFirstBand->m_dfNoData.value() ==
1587
3.53k
                               poBand->m_dfNoData.value()));
1588
3.54k
            }
1589
27
            if (!bSameValue)
1590
7
            {
1591
7
                CPLError(CE_Failure, CPLE_NotSupported,
1592
7
                         "Not all bands have the same nodata value. It will be "
1593
7
                         "ignored as the array can only have a single nodata "
1594
7
                         "value for all bands.");
1595
7
                eErr = CE_Failure;
1596
7
            }
1597
20
            else
1598
20
            {
1599
20
                m_poSingleArray->SetNoDataValue(
1600
20
                    poFirstBand->m_dfNoData.value());
1601
20
            }
1602
27
        }
1603
48
        else if (poFirstBand->m_nNoDataInt64.has_value())
1604
0
        {
1605
0
            bool bSameValue = true;
1606
0
            for (int i = 1; bSameValue && i < nBands; ++i)
1607
0
            {
1608
0
                auto poBand = cpl::down_cast<ZarrRasterBand *>(papoBands[i]);
1609
0
                bSameValue = poBand->m_nNoDataInt64.has_value() &&
1610
0
                             poFirstBand->m_nNoDataInt64.value() ==
1611
0
                                 poBand->m_nNoDataInt64.value();
1612
0
            }
1613
0
            if (!bSameValue)
1614
0
            {
1615
0
                CPLError(CE_Failure, CPLE_NotSupported,
1616
0
                         "Not all bands have the same nodata value. It will be "
1617
0
                         "ignored as the array can only have a single nodata "
1618
0
                         "value for all bands.");
1619
0
                eErr = CE_Failure;
1620
0
            }
1621
0
            else
1622
0
            {
1623
0
                m_poSingleArray->SetNoDataValue(
1624
0
                    poFirstBand->m_nNoDataInt64.value());
1625
0
            }
1626
0
        }
1627
48
        else if (poFirstBand->m_nNoDataUInt64.has_value())
1628
0
        {
1629
0
            bool bSameValue = true;
1630
0
            for (int i = 1; bSameValue && i < nBands; ++i)
1631
0
            {
1632
0
                auto poBand = cpl::down_cast<ZarrRasterBand *>(papoBands[i]);
1633
0
                bSameValue = poBand->m_nNoDataUInt64.has_value() &&
1634
0
                             poFirstBand->m_nNoDataUInt64.value() ==
1635
0
                                 poBand->m_nNoDataUInt64.value();
1636
0
            }
1637
0
            if (!bSameValue)
1638
0
            {
1639
0
                CPLError(CE_Failure, CPLE_NotSupported,
1640
0
                         "Not all bands have the same nodata value. It will be "
1641
0
                         "ignored as the array can only have a single nodata "
1642
0
                         "value for all bands.");
1643
0
                eErr = CE_Failure;
1644
0
            }
1645
0
            else
1646
0
            {
1647
0
                m_poSingleArray->SetNoDataValue(
1648
0
                    poFirstBand->m_nNoDataUInt64.value());
1649
0
            }
1650
0
        }
1651
1652
75
        if (poFirstBand->m_dfOffset.has_value())
1653
16
        {
1654
16
            bool bSameValue = true;
1655
495
            for (int i = 1; bSameValue && i < nBands; ++i)
1656
479
            {
1657
479
                auto poBand = cpl::down_cast<ZarrRasterBand *>(papoBands[i]);
1658
479
                bSameValue = poBand->m_dfOffset.has_value() &&
1659
464
                             poFirstBand->m_dfOffset.value() ==
1660
464
                                 poBand->m_dfOffset.value();
1661
479
            }
1662
16
            if (!bSameValue)
1663
16
            {
1664
16
                CPLError(CE_Failure, CPLE_NotSupported,
1665
16
                         "Not all bands have the same offset value. It will be "
1666
16
                         "ignored as the array can only have a single offset "
1667
16
                         "value for all bands.");
1668
16
                eErr = CE_Failure;
1669
16
            }
1670
0
            else
1671
0
            {
1672
0
                m_poSingleArray->SetOffset(poFirstBand->m_dfOffset.value());
1673
0
            }
1674
16
        }
1675
1676
75
        if (poFirstBand->m_dfScale.has_value())
1677
53
        {
1678
53
            bool bSameValue = true;
1679
2.63k
            for (int i = 1; bSameValue && i < nBands; ++i)
1680
2.58k
            {
1681
2.58k
                auto poBand = cpl::down_cast<ZarrRasterBand *>(papoBands[i]);
1682
2.58k
                bSameValue =
1683
2.58k
                    poBand->m_dfScale.has_value() &&
1684
2.53k
                    poFirstBand->m_dfScale.value() == poBand->m_dfScale.value();
1685
2.58k
            }
1686
53
            if (!bSameValue)
1687
53
            {
1688
53
                CPLError(CE_Failure, CPLE_NotSupported,
1689
53
                         "Not all bands have the same scale value. It will be "
1690
53
                         "ignored as the array can only have a single scale "
1691
53
                         "value for all bands.");
1692
53
                eErr = CE_Failure;
1693
53
            }
1694
0
            else
1695
0
            {
1696
0
                m_poSingleArray->SetScale(poFirstBand->m_dfScale.value());
1697
0
            }
1698
53
        }
1699
1700
75
        bool bFoundColorInterp = false;
1701
12.2k
        for (int i = 0; i < nBands; ++i)
1702
12.1k
        {
1703
12.1k
            if (papoBands[i]->GetColorInterpretation() != GCI_Undefined)
1704
9.02k
                bFoundColorInterp = true;
1705
12.1k
        }
1706
75
        if (bFoundColorInterp)
1707
56
        {
1708
56
            const auto oStringDT = GDALExtendedDataType::CreateString();
1709
56
            auto poAttr = m_poSingleArray->GetAttribute("COLOR_INTERPRETATION");
1710
56
            if (!poAttr)
1711
56
                poAttr = m_poSingleArray->CreateAttribute(
1712
56
                    "COLOR_INTERPRETATION", {static_cast<GUInt64>(nBands)},
1713
56
                    oStringDT);
1714
56
            if (poAttr)
1715
56
            {
1716
56
                const GUInt64 nStartIndex = 0;
1717
56
                const size_t nCount = nBands;
1718
56
                const GInt64 arrayStep = 1;
1719
56
                const GPtrDiff_t bufferStride = 1;
1720
56
                std::vector<const char *> apszValues;
1721
10.1k
                for (int i = 0; i < nBands; ++i)
1722
10.1k
                {
1723
10.1k
                    const auto eColorInterp =
1724
10.1k
                        papoBands[i]->GetColorInterpretation();
1725
10.1k
                    apszValues.push_back(
1726
10.1k
                        GDALGetColorInterpretationName(eColorInterp));
1727
10.1k
                }
1728
56
                poAttr->Write(&nStartIndex, &nCount, &arrayStep, &bufferStride,
1729
56
                              oStringDT, apszValues.data());
1730
56
            }
1731
56
        }
1732
75
    }
1733
1734
176
    if (m_poRootGroup)
1735
176
    {
1736
176
        if (bAtClosing)
1737
76
        {
1738
76
            if (!m_poRootGroup->Close())
1739
0
                eErr = CE_Failure;
1740
76
        }
1741
100
        else
1742
100
        {
1743
100
            if (!m_poRootGroup->Flush())
1744
0
                eErr = CE_Failure;
1745
100
        }
1746
176
    }
1747
1748
176
    return eErr;
1749
176
}
1750
1751
/************************************************************************/
1752
/*                            GetRootGroup()                            */
1753
/************************************************************************/
1754
1755
std::shared_ptr<GDALGroup> ZarrDataset::GetRootGroup() const
1756
0
{
1757
0
    return m_poRootGroup;
1758
0
}
1759
1760
/************************************************************************/
1761
/*                           GetSpatialRef()                            */
1762
/************************************************************************/
1763
1764
const OGRSpatialReference *ZarrDataset::GetSpatialRef() const
1765
0
{
1766
0
    if (m_poSingleArray)
1767
0
    {
1768
0
        return m_poSingleArray->GetSpatialRef().get();
1769
0
    }
1770
0
    else if (nBands >= 1)
1771
0
    {
1772
0
        return cpl::down_cast<ZarrRasterBand *>(papoBands[0])
1773
0
            ->m_poArray->GetSpatialRef()
1774
0
            .get();
1775
0
    }
1776
0
    else
1777
0
    {
1778
0
        return nullptr;
1779
0
    }
1780
0
}
1781
1782
/************************************************************************/
1783
/*                           SetSpatialRef()                            */
1784
/************************************************************************/
1785
1786
CPLErr ZarrDataset::SetSpatialRef(const OGRSpatialReference *poSRS)
1787
46
{
1788
46
    if (m_poSingleArray)
1789
45
    {
1790
45
        m_poSingleArray->SetSpatialRef(poSRS);
1791
45
    }
1792
1
    else
1793
1
    {
1794
2
        for (int i = 0; i < nBands; ++i)
1795
1
        {
1796
1
            cpl::down_cast<ZarrRasterBand *>(papoBands[i])
1797
1
                ->m_poArray->SetSpatialRef(poSRS);
1798
1
        }
1799
1
    }
1800
46
    return CE_None;
1801
46
}
1802
1803
/************************************************************************/
1804
/*                          GetGeoTransform()                           */
1805
/************************************************************************/
1806
1807
CPLErr ZarrDataset::GetGeoTransform(GDALGeoTransform &gt) const
1808
0
{
1809
0
    gt = m_gt;
1810
0
    return m_bHasGT ? CE_None : CE_Failure;
1811
0
}
1812
1813
/************************************************************************/
1814
/*                          SetGeoTransform()                           */
1815
/************************************************************************/
1816
1817
CPLErr ZarrDataset::SetGeoTransform(const GDALGeoTransform &gt)
1818
37
{
1819
37
    const bool bHasRotatedTerms = (gt.xrot != 0 || gt.yrot != 0);
1820
1821
37
    if (bHasRotatedTerms)
1822
0
    {
1823
0
        if (!m_bSpatialProjConvention)
1824
0
        {
1825
0
            CPLError(CE_Failure, CPLE_NotSupported,
1826
0
                     "Geotransform with rotated terms not supported with "
1827
0
                     "GEOREFERENCING_CONVENTION=GDAL, but would be with "
1828
0
                     "SPATIAL_PROJ");
1829
0
            return CE_Failure;
1830
0
        }
1831
0
    }
1832
37
    else if (m_poDimX == nullptr || m_poDimY == nullptr)
1833
0
    {
1834
0
        CPLError(CE_Failure, CPLE_AppDefined,
1835
0
                 "SetGeoTransform() failed because of missing X/Y dimension");
1836
0
        return CE_Failure;
1837
0
    }
1838
1839
37
    m_gt = gt;
1840
37
    m_bHasGT = true;
1841
1842
37
    if (m_bSpatialProjConvention)
1843
0
    {
1844
0
        const auto bSingleArray = m_poSingleArray != nullptr;
1845
0
        const int nIters = bSingleArray ? 1 : nBands;
1846
0
        for (int i = 0; i < nIters; ++i)
1847
0
        {
1848
0
            auto *poArray = bSingleArray
1849
0
                                ? m_poSingleArray.get()
1850
0
                                : cpl::down_cast<ZarrRasterBand *>(papoBands[i])
1851
0
                                      ->m_poArray.get();
1852
0
            auto oAttrDT = GDALExtendedDataType::Create(GDT_Float64);
1853
0
            auto poAttr =
1854
0
                poArray->CreateAttribute("gdal:geotransform", {6}, oAttrDT);
1855
0
            if (poAttr)
1856
0
            {
1857
0
                const GUInt64 nStartIndex = 0;
1858
0
                const size_t nCount = 6;
1859
0
                const GInt64 arrayStep = 1;
1860
0
                const GPtrDiff_t bufferStride = 1;
1861
0
                poAttr->Write(&nStartIndex, &nCount, &arrayStep, &bufferStride,
1862
0
                              oAttrDT, m_gt.data());
1863
0
            }
1864
0
        }
1865
0
    }
1866
1867
37
    if (!bHasRotatedTerms)
1868
37
    {
1869
37
        CPLAssert(m_poDimX);
1870
37
        CPLAssert(m_poDimY);
1871
1872
37
        const auto oDTFloat64 = GDALExtendedDataType::Create(GDT_Float64);
1873
37
        {
1874
37
            auto poX = m_poRootGroup->OpenMDArray(m_poDimX->GetName());
1875
37
            if (!poX)
1876
37
                poX = m_poRootGroup->CreateMDArray(
1877
37
                    m_poDimX->GetName(), {m_poDimX}, oDTFloat64, nullptr);
1878
37
            if (!poX)
1879
0
                return CE_Failure;
1880
37
            m_poDimX->SetIndexingVariable(poX);
1881
37
            std::vector<double> adfX;
1882
37
            try
1883
37
            {
1884
37
                adfX.reserve(nRasterXSize);
1885
3.37k
                for (int i = 0; i < nRasterXSize; ++i)
1886
3.33k
                    adfX.emplace_back(m_gt.xorig + m_gt.xscale * (i + 0.5));
1887
37
            }
1888
37
            catch (const std::exception &)
1889
37
            {
1890
0
                CPLError(CE_Failure, CPLE_OutOfMemory,
1891
0
                         "Out of memory when allocating X array");
1892
0
                return CE_Failure;
1893
0
            }
1894
37
            const GUInt64 nStartIndex = 0;
1895
37
            const size_t nCount = adfX.size();
1896
37
            const GInt64 arrayStep = 1;
1897
37
            const GPtrDiff_t bufferStride = 1;
1898
37
            if (!poX->Write(&nStartIndex, &nCount, &arrayStep, &bufferStride,
1899
37
                            poX->GetDataType(), adfX.data()))
1900
0
            {
1901
0
                return CE_Failure;
1902
0
            }
1903
37
        }
1904
1905
37
        auto poY = m_poRootGroup->OpenMDArray(m_poDimY->GetName());
1906
37
        if (!poY)
1907
37
            poY = m_poRootGroup->CreateMDArray(m_poDimY->GetName(), {m_poDimY},
1908
37
                                               oDTFloat64, nullptr);
1909
37
        if (!poY)
1910
0
            return CE_Failure;
1911
37
        m_poDimY->SetIndexingVariable(poY);
1912
37
        std::vector<double> adfY;
1913
37
        try
1914
37
        {
1915
37
            adfY.reserve(nRasterYSize);
1916
4.29k
            for (int i = 0; i < nRasterYSize; ++i)
1917
4.25k
                adfY.emplace_back(m_gt.yorig + m_gt.yscale * (i + 0.5));
1918
37
        }
1919
37
        catch (const std::exception &)
1920
37
        {
1921
0
            CPLError(CE_Failure, CPLE_OutOfMemory,
1922
0
                     "Out of memory when allocating Y array");
1923
0
            return CE_Failure;
1924
0
        }
1925
37
        const GUInt64 nStartIndex = 0;
1926
37
        const size_t nCount = adfY.size();
1927
37
        const GInt64 arrayStep = 1;
1928
37
        const GPtrDiff_t bufferStride = 1;
1929
37
        if (!poY->Write(&nStartIndex, &nCount, &arrayStep, &bufferStride,
1930
37
                        poY->GetDataType(), adfY.data()))
1931
0
        {
1932
0
            return CE_Failure;
1933
0
        }
1934
37
    }
1935
1936
37
    return CE_None;
1937
37
}
1938
1939
/************************************************************************/
1940
/*                            SetMetadata()                             */
1941
/************************************************************************/
1942
1943
CPLErr ZarrDataset::SetMetadata(CSLConstList papszMetadata,
1944
                                const char *pszDomain)
1945
64
{
1946
64
    if (nBands >= 1 && (pszDomain == nullptr || pszDomain[0] == '\0'))
1947
64
    {
1948
64
        const auto oStringDT = GDALExtendedDataType::CreateString();
1949
64
        const auto bSingleArray = m_poSingleArray != nullptr;
1950
64
        const int nIters = bSingleArray ? 1 : nBands;
1951
128
        for (int i = 0; i < nIters; ++i)
1952
64
        {
1953
64
            auto *poArray = bSingleArray
1954
64
                                ? m_poSingleArray.get()
1955
64
                                : cpl::down_cast<ZarrRasterBand *>(papoBands[i])
1956
1
                                      ->m_poArray.get();
1957
1.53k
            for (auto iter = papszMetadata; iter && *iter; ++iter)
1958
1.46k
            {
1959
1.46k
                char *pszKey = nullptr;
1960
1.46k
                const char *pszValue = CPLParseNameValue(*iter, &pszKey);
1961
1.46k
                if (pszKey && pszValue)
1962
1.46k
                {
1963
1.46k
                    auto poAttr =
1964
1.46k
                        poArray->CreateAttribute(pszKey, {}, oStringDT);
1965
1.46k
                    if (poAttr)
1966
1.46k
                    {
1967
1.46k
                        const GUInt64 nStartIndex = 0;
1968
1.46k
                        const size_t nCount = 1;
1969
1.46k
                        const GInt64 arrayStep = 1;
1970
1.46k
                        const GPtrDiff_t bufferStride = 1;
1971
1.46k
                        poAttr->Write(&nStartIndex, &nCount, &arrayStep,
1972
1.46k
                                      &bufferStride, oStringDT, &pszValue);
1973
1.46k
                    }
1974
1.46k
                }
1975
1.46k
                CPLFree(pszKey);
1976
1.46k
            }
1977
64
        }
1978
64
    }
1979
64
    return GDALDataset::SetMetadata(papszMetadata, pszDomain);
1980
64
}
1981
1982
/************************************************************************/
1983
/*                   ZarrRasterBand::ZarrRasterBand()                   */
1984
/************************************************************************/
1985
1986
ZarrRasterBand::ZarrRasterBand(const std::shared_ptr<GDALMDArray> &poArray,
1987
                               const std::string &viewDef)
1988
12.1k
    : m_poArray(poArray), m_osViewDef(viewDef)
1989
12.1k
{
1990
12.1k
    assert(poArray->GetDimensionCount() == 2);
1991
12.1k
    eDataType = poArray->GetDataType().GetNumericDataType();
1992
12.1k
    nRasterXSize = static_cast<int>(poArray->GetDimensions()[1]->GetSize());
1993
12.1k
    nRasterYSize = static_cast<int>(poArray->GetDimensions()[0]->GetSize());
1994
12.1k
    nBlockXSize = static_cast<int>(poArray->GetBlockSize()[1]);
1995
12.1k
    nBlockYSize = static_cast<int>(poArray->GetBlockSize()[0]);
1996
12.1k
}
1997
1998
/************************************************************************/
1999
/*                           GetNoDataValue()                           */
2000
/************************************************************************/
2001
2002
double ZarrRasterBand::GetNoDataValue(int *pbHasNoData)
2003
0
{
2004
0
    if (m_dfNoData.has_value())
2005
0
    {
2006
0
        if (pbHasNoData)
2007
0
            *pbHasNoData = true;
2008
0
        return m_dfNoData.value();
2009
0
    }
2010
0
    bool bHasNodata = false;
2011
0
    const auto res = m_poArray->GetNoDataValueAsDouble(&bHasNodata);
2012
0
    if (pbHasNoData)
2013
0
        *pbHasNoData = bHasNodata;
2014
0
    return res;
2015
0
}
2016
2017
/************************************************************************/
2018
/*                       GetNoDataValueAsInt64()                        */
2019
/************************************************************************/
2020
2021
int64_t ZarrRasterBand::GetNoDataValueAsInt64(int *pbHasNoData)
2022
0
{
2023
0
    if (m_nNoDataInt64.has_value())
2024
0
    {
2025
0
        if (pbHasNoData)
2026
0
            *pbHasNoData = true;
2027
0
        return m_nNoDataInt64.value();
2028
0
    }
2029
0
    bool bHasNodata = false;
2030
0
    const auto res = m_poArray->GetNoDataValueAsInt64(&bHasNodata);
2031
0
    if (pbHasNoData)
2032
0
        *pbHasNoData = bHasNodata;
2033
0
    return res;
2034
0
}
2035
2036
/************************************************************************/
2037
/*                       GetNoDataValueAsUInt64()                       */
2038
/************************************************************************/
2039
2040
uint64_t ZarrRasterBand::GetNoDataValueAsUInt64(int *pbHasNoData)
2041
0
{
2042
0
    if (m_nNoDataUInt64.has_value())
2043
0
    {
2044
0
        if (pbHasNoData)
2045
0
            *pbHasNoData = true;
2046
0
        return m_nNoDataUInt64.value();
2047
0
    }
2048
0
    bool bHasNodata = false;
2049
0
    const auto res = m_poArray->GetNoDataValueAsUInt64(&bHasNodata);
2050
0
    if (pbHasNoData)
2051
0
        *pbHasNoData = bHasNodata;
2052
0
    return res;
2053
0
}
2054
2055
/************************************************************************/
2056
/*                           SetNoDataValue()                           */
2057
/************************************************************************/
2058
2059
CPLErr ZarrRasterBand::SetNoDataValue(double dfNoData)
2060
4.48k
{
2061
4.48k
    auto poGDS = cpl::down_cast<ZarrDataset *>(poDS);
2062
4.48k
    if (!poGDS->m_poSingleArray)
2063
0
    {
2064
0
        return m_poArray->SetNoDataValue(dfNoData) ? CE_None : CE_Failure;
2065
0
    }
2066
4.48k
    m_dfNoData = dfNoData;
2067
4.48k
    return CE_None;
2068
4.48k
}
2069
2070
/************************************************************************/
2071
/*                       SetNoDataValueAsInt64()                        */
2072
/************************************************************************/
2073
2074
CPLErr ZarrRasterBand::SetNoDataValueAsInt64(int64_t nNoData)
2075
0
{
2076
0
    auto poGDS = cpl::down_cast<ZarrDataset *>(poDS);
2077
0
    if (!poGDS->m_poSingleArray)
2078
0
    {
2079
0
        return m_poArray->SetNoDataValue(nNoData) ? CE_None : CE_Failure;
2080
0
    }
2081
0
    m_nNoDataInt64 = nNoData;
2082
0
    return CE_None;
2083
0
}
2084
2085
/************************************************************************/
2086
/*                       SetNoDataValueAsUInt64()                       */
2087
/************************************************************************/
2088
2089
CPLErr ZarrRasterBand::SetNoDataValueAsUInt64(uint64_t nNoData)
2090
0
{
2091
0
    auto poGDS = cpl::down_cast<ZarrDataset *>(poDS);
2092
0
    if (!poGDS->m_poSingleArray)
2093
0
    {
2094
0
        return m_poArray->SetNoDataValue(nNoData) ? CE_None : CE_Failure;
2095
0
    }
2096
0
    m_nNoDataUInt64 = nNoData;
2097
0
    return CE_None;
2098
0
}
2099
2100
/************************************************************************/
2101
/*                             GetOffset()                              */
2102
/************************************************************************/
2103
2104
double ZarrRasterBand::GetOffset(int *pbSuccess)
2105
0
{
2106
0
    if (m_dfOffset.has_value())
2107
0
    {
2108
0
        if (pbSuccess)
2109
0
            *pbSuccess = true;
2110
0
        return m_dfOffset.value();
2111
0
    }
2112
0
    bool bHasValue = false;
2113
0
    double dfRet = m_poArray->GetOffset(&bHasValue);
2114
0
    if (pbSuccess)
2115
0
        *pbSuccess = bHasValue ? TRUE : FALSE;
2116
0
    return dfRet;
2117
0
}
2118
2119
/************************************************************************/
2120
/*                             SetOffset()                              */
2121
/************************************************************************/
2122
2123
CPLErr ZarrRasterBand::SetOffset(double dfNewOffset)
2124
787
{
2125
787
    auto poGDS = cpl::down_cast<ZarrDataset *>(poDS);
2126
787
    if (!poGDS->m_poSingleArray)
2127
0
    {
2128
0
        return m_poArray->SetOffset(dfNewOffset) ? CE_None : CE_Failure;
2129
0
    }
2130
787
    m_dfOffset = dfNewOffset;
2131
787
    return CE_None;
2132
787
}
2133
2134
/************************************************************************/
2135
/*                              GetScale()                              */
2136
/************************************************************************/
2137
2138
double ZarrRasterBand::GetScale(int *pbSuccess)
2139
0
{
2140
0
    if (m_dfScale.has_value())
2141
0
    {
2142
0
        if (pbSuccess)
2143
0
            *pbSuccess = true;
2144
0
        return m_dfScale.value();
2145
0
    }
2146
0
    bool bHasValue = false;
2147
0
    double dfRet = m_poArray->GetScale(&bHasValue);
2148
0
    if (pbSuccess)
2149
0
        *pbSuccess = bHasValue ? TRUE : FALSE;
2150
0
    return dfRet;
2151
0
}
2152
2153
/************************************************************************/
2154
/*                              SetScale()                              */
2155
/************************************************************************/
2156
2157
CPLErr ZarrRasterBand::SetScale(double dfNewScale)
2158
3.29k
{
2159
3.29k
    auto poGDS = cpl::down_cast<ZarrDataset *>(poDS);
2160
3.29k
    if (!poGDS->m_poSingleArray)
2161
0
    {
2162
0
        return m_poArray->SetScale(dfNewScale) ? CE_None : CE_Failure;
2163
0
    }
2164
3.29k
    m_dfScale = dfNewScale;
2165
3.29k
    return CE_None;
2166
3.29k
}
2167
2168
/************************************************************************/
2169
/*                            GetUnitType()                             */
2170
/************************************************************************/
2171
2172
const char *ZarrRasterBand::GetUnitType()
2173
0
{
2174
0
    return m_poArray->GetUnit().c_str();
2175
0
}
2176
2177
/************************************************************************/
2178
/*                            SetUnitType()                             */
2179
/************************************************************************/
2180
2181
CPLErr ZarrRasterBand::SetUnitType(const char *pszNewValue)
2182
0
{
2183
0
    return m_poArray->SetUnit(pszNewValue ? pszNewValue : "") ? CE_None
2184
0
                                                              : CE_Failure;
2185
0
}
2186
2187
/************************************************************************/
2188
/*                       GetColorInterpretation()                       */
2189
/************************************************************************/
2190
2191
GDALColorInterp ZarrRasterBand::GetColorInterpretation()
2192
31.3k
{
2193
31.3k
    return m_eColorInterp;
2194
31.3k
}
2195
2196
/************************************************************************/
2197
/*                       SetColorInterpretation()                       */
2198
/************************************************************************/
2199
2200
CPLErr ZarrRasterBand::SetColorInterpretation(GDALColorInterp eColorInterp)
2201
9.02k
{
2202
9.02k
    auto poGDS = cpl::down_cast<ZarrDataset *>(poDS);
2203
9.02k
    m_eColorInterp = eColorInterp;
2204
9.02k
    if (!poGDS->m_poSingleArray)
2205
0
    {
2206
0
        const auto oStringDT = GDALExtendedDataType::CreateString();
2207
0
        auto poAttr = m_poArray->GetAttribute("COLOR_INTERPRETATION");
2208
0
        if (poAttr && (poAttr->GetDimensionCount() != 0 ||
2209
0
                       poAttr->GetDataType().GetClass() != GEDTC_STRING))
2210
0
            return CE_None;
2211
0
        if (!poAttr)
2212
0
            poAttr = m_poArray->CreateAttribute("COLOR_INTERPRETATION", {},
2213
0
                                                oStringDT);
2214
0
        if (poAttr)
2215
0
        {
2216
0
            const GUInt64 nStartIndex = 0;
2217
0
            const size_t nCount = 1;
2218
0
            const GInt64 arrayStep = 1;
2219
0
            const GPtrDiff_t bufferStride = 1;
2220
0
            const char *pszValue = GDALGetColorInterpretationName(eColorInterp);
2221
0
            poAttr->Write(&nStartIndex, &nCount, &arrayStep, &bufferStride,
2222
0
                          oStringDT, &pszValue);
2223
0
        }
2224
0
    }
2225
9.02k
    return CE_None;
2226
9.02k
}
2227
2228
/************************************************************************/
2229
/*                          GetOverviewCount()                          */
2230
/************************************************************************/
2231
2232
int ZarrRasterBand::GetOverviewCount()
2233
0
{
2234
0
    auto poGDS = cpl::down_cast<ZarrDataset *>(poDS);
2235
0
    if (poGDS->m_poSingleArray)
2236
0
        return poGDS->m_poSingleArray->GetOverviewCount();
2237
0
    return m_poArray->GetOverviewCount();
2238
0
}
2239
2240
/************************************************************************/
2241
/*                            GetOverview()                             */
2242
/************************************************************************/
2243
2244
GDALRasterBand *ZarrRasterBand::GetOverview(int idx)
2245
0
{
2246
0
    auto oIter = m_oMapOverview.find(idx);
2247
0
    if (oIter != m_oMapOverview.end())
2248
0
        return oIter->second.get();
2249
0
    auto poGDS = cpl::down_cast<ZarrDataset *>(poDS);
2250
0
    if (poGDS->m_poSingleArray)
2251
0
    {
2252
0
        auto ovrArray = poGDS->m_poSingleArray->GetOverview(idx);
2253
0
        if (!ovrArray)
2254
0
            return nullptr;
2255
0
        auto ovrArrayView = ovrArray->GetView(m_osViewDef);
2256
0
        if (!ovrArrayView)
2257
0
            return nullptr;  // not supposed to happen
2258
0
        return m_oMapOverview
2259
0
            .insert({idx, std::make_unique<ZarrRasterBand>(ovrArrayView,
2260
0
                                                           m_osViewDef)})
2261
0
            .first->second.get();
2262
0
    }
2263
0
    else
2264
0
    {
2265
0
        auto ovrArray = m_poArray->GetOverview(idx);
2266
0
        if (!ovrArray)
2267
0
            return nullptr;  // not supposed to happen
2268
0
        return m_oMapOverview
2269
0
            .insert({idx, std::make_unique<ZarrRasterBand>(ovrArray)})
2270
0
            .first->second.get();
2271
0
    }
2272
0
}
2273
2274
/************************************************************************/
2275
/*                     ZarrRasterBand::IReadBlock()                     */
2276
/************************************************************************/
2277
2278
CPLErr ZarrRasterBand::IReadBlock(int nBlockXOff, int nBlockYOff, void *pData)
2279
0
{
2280
2281
0
    const int nXOff = nBlockXOff * nBlockXSize;
2282
0
    const int nYOff = nBlockYOff * nBlockYSize;
2283
0
    const int nReqXSize = std::min(nRasterXSize - nXOff, nBlockXSize);
2284
0
    const int nReqYSize = std::min(nRasterYSize - nYOff, nBlockYSize);
2285
0
    GUInt64 arrayStartIdx[] = {static_cast<GUInt64>(nYOff),
2286
0
                               static_cast<GUInt64>(nXOff)};
2287
0
    size_t count[] = {static_cast<size_t>(nReqYSize),
2288
0
                      static_cast<size_t>(nReqXSize)};
2289
0
    constexpr GInt64 arrayStep[] = {1, 1};
2290
0
    GPtrDiff_t bufferStride[] = {nBlockXSize, 1};
2291
0
    return m_poArray->Read(arrayStartIdx, count, arrayStep, bufferStride,
2292
0
                           m_poArray->GetDataType(), pData)
2293
0
               ? CE_None
2294
0
               : CE_Failure;
2295
0
}
2296
2297
/************************************************************************/
2298
/*                    ZarrRasterBand::IWriteBlock()                     */
2299
/************************************************************************/
2300
2301
CPLErr ZarrRasterBand::IWriteBlock(int nBlockXOff, int nBlockYOff, void *pData)
2302
0
{
2303
0
    const int nXOff = nBlockXOff * nBlockXSize;
2304
0
    const int nYOff = nBlockYOff * nBlockYSize;
2305
0
    const int nReqXSize = std::min(nRasterXSize - nXOff, nBlockXSize);
2306
0
    const int nReqYSize = std::min(nRasterYSize - nYOff, nBlockYSize);
2307
0
    GUInt64 arrayStartIdx[] = {static_cast<GUInt64>(nYOff),
2308
0
                               static_cast<GUInt64>(nXOff)};
2309
0
    size_t count[] = {static_cast<size_t>(nReqYSize),
2310
0
                      static_cast<size_t>(nReqXSize)};
2311
0
    constexpr GInt64 arrayStep[] = {1, 1};
2312
0
    GPtrDiff_t bufferStride[] = {nBlockXSize, 1};
2313
0
    return m_poArray->Write(arrayStartIdx, count, arrayStep, bufferStride,
2314
0
                            m_poArray->GetDataType(), pData)
2315
0
               ? CE_None
2316
0
               : CE_Failure;
2317
0
}
2318
2319
/************************************************************************/
2320
/*                             IRasterIO()                              */
2321
/************************************************************************/
2322
2323
CPLErr ZarrRasterBand::IRasterIO(GDALRWFlag eRWFlag, int nXOff, int nYOff,
2324
                                 int nXSize, int nYSize, void *pData,
2325
                                 int nBufXSize, int nBufYSize,
2326
                                 GDALDataType eBufType, GSpacing nPixelSpaceBuf,
2327
                                 GSpacing nLineSpaceBuf,
2328
                                 GDALRasterIOExtraArg *psExtraArg)
2329
9.34k
{
2330
9.34k
    const int nBufferDTSize(GDALGetDataTypeSizeBytes(eBufType));
2331
    // If reading/writing at full resolution and with proper stride, go
2332
    // directly to the array, but, for performance reasons,
2333
    // only if exactly on chunk boundaries, otherwise go through the block cache.
2334
9.34k
    if (nXSize == nBufXSize && nYSize == nBufYSize && nBufferDTSize > 0 &&
2335
9.34k
        (nPixelSpaceBuf % nBufferDTSize) == 0 &&
2336
9.34k
        (nLineSpaceBuf % nBufferDTSize) == 0 && (nXOff % nBlockXSize) == 0 &&
2337
9.34k
        (nYOff % nBlockYSize) == 0 &&
2338
9.34k
        ((nXSize % nBlockXSize) == 0 || nXOff + nXSize == nRasterXSize) &&
2339
9.34k
        ((nYSize % nBlockYSize) == 0 || nYOff + nYSize == nRasterYSize))
2340
9.34k
    {
2341
9.34k
        GUInt64 arrayStartIdx[] = {static_cast<GUInt64>(nYOff),
2342
9.34k
                                   static_cast<GUInt64>(nXOff)};
2343
9.34k
        size_t count[] = {static_cast<size_t>(nYSize),
2344
9.34k
                          static_cast<size_t>(nXSize)};
2345
9.34k
        constexpr GInt64 arrayStep[] = {1, 1};
2346
9.34k
        GPtrDiff_t bufferStride[] = {
2347
9.34k
            static_cast<GPtrDiff_t>(nLineSpaceBuf / nBufferDTSize),
2348
9.34k
            static_cast<GPtrDiff_t>(nPixelSpaceBuf / nBufferDTSize)};
2349
2350
9.34k
        if (eRWFlag == GF_Read)
2351
0
        {
2352
0
            return m_poArray->Read(
2353
0
                       arrayStartIdx, count, arrayStep, bufferStride,
2354
0
                       GDALExtendedDataType::Create(eBufType), pData)
2355
0
                       ? CE_None
2356
0
                       : CE_Failure;
2357
0
        }
2358
9.34k
        else
2359
9.34k
        {
2360
9.34k
            return m_poArray->Write(
2361
9.34k
                       arrayStartIdx, count, arrayStep, bufferStride,
2362
9.34k
                       GDALExtendedDataType::Create(eBufType), pData)
2363
9.34k
                       ? CE_None
2364
9.34k
                       : CE_Failure;
2365
9.34k
        }
2366
9.34k
    }
2367
2368
0
    return GDALRasterBand::IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize,
2369
0
                                     pData, nBufXSize, nBufYSize, eBufType,
2370
0
                                     nPixelSpaceBuf, nLineSpaceBuf, psExtraArg);
2371
9.34k
}
2372
2373
/************************************************************************/
2374
/*                       ZarrDataset::IRasterIO()                       */
2375
/************************************************************************/
2376
2377
CPLErr ZarrDataset::IRasterIO(GDALRWFlag eRWFlag, int nXOff, int nYOff,
2378
                              int nXSize, int nYSize, void *pData,
2379
                              int nBufXSize, int nBufYSize,
2380
                              GDALDataType eBufType, int nBandCount,
2381
                              BANDMAP_TYPE panBandMap, GSpacing nPixelSpace,
2382
                              GSpacing nLineSpace, GSpacing nBandSpace,
2383
                              GDALRasterIOExtraArg *psExtraArg)
2384
9.35k
{
2385
9.35k
    const int nBufferDTSize(GDALGetDataTypeSizeBytes(eBufType));
2386
    // If reading/writing at full resolution and with proper stride, go
2387
    // directly to the array, but, for performance reasons,
2388
    // only if exactly on chunk boundaries, otherwise go through the block cache.
2389
9.35k
    int nBlockXSize = 0, nBlockYSize = 0;
2390
9.35k
    papoBands[0]->GetBlockSize(&nBlockXSize, &nBlockYSize);
2391
9.35k
    if (m_poSingleArray && nXSize == nBufXSize && nYSize == nBufYSize &&
2392
9.35k
        nBufferDTSize > 0 && (nPixelSpace % nBufferDTSize) == 0 &&
2393
9.35k
        (nLineSpace % nBufferDTSize) == 0 &&
2394
9.35k
        (nBandSpace % nBufferDTSize) == 0 && (nXOff % nBlockXSize) == 0 &&
2395
9.35k
        (nYOff % nBlockYSize) == 0 &&
2396
9.35k
        ((nXSize % nBlockXSize) == 0 || nXOff + nXSize == nRasterXSize) &&
2397
9.35k
        ((nYSize % nBlockYSize) == 0 || nYOff + nYSize == nRasterYSize) &&
2398
9.35k
        IsAllBands(nBandCount, panBandMap))
2399
3
    {
2400
3
        CPLAssert(m_poSingleArray->GetDimensionCount() == 3);
2401
3
        if (m_poSingleArray->GetDimensions().back().get() == m_poDimX.get())
2402
0
        {
2403
0
            GUInt64 arrayStartIdx[] = {0, static_cast<GUInt64>(nYOff),
2404
0
                                       static_cast<GUInt64>(nXOff)};
2405
0
            size_t count[] = {static_cast<size_t>(nBands),
2406
0
                              static_cast<size_t>(nYSize),
2407
0
                              static_cast<size_t>(nXSize)};
2408
0
            constexpr GInt64 arrayStep[] = {1, 1, 1};
2409
0
            GPtrDiff_t bufferStride[] = {
2410
0
                static_cast<GPtrDiff_t>(nBandSpace / nBufferDTSize),
2411
0
                static_cast<GPtrDiff_t>(nLineSpace / nBufferDTSize),
2412
0
                static_cast<GPtrDiff_t>(nPixelSpace / nBufferDTSize)};
2413
2414
0
            if (eRWFlag == GF_Read)
2415
0
            {
2416
0
                return m_poSingleArray->Read(
2417
0
                           arrayStartIdx, count, arrayStep, bufferStride,
2418
0
                           GDALExtendedDataType::Create(eBufType), pData)
2419
0
                           ? CE_None
2420
0
                           : CE_Failure;
2421
0
            }
2422
0
            else
2423
0
            {
2424
0
                return m_poSingleArray->Write(
2425
0
                           arrayStartIdx, count, arrayStep, bufferStride,
2426
0
                           GDALExtendedDataType::Create(eBufType), pData)
2427
0
                           ? CE_None
2428
0
                           : CE_Failure;
2429
0
            }
2430
0
        }
2431
3
        else
2432
3
        {
2433
3
            GUInt64 arrayStartIdx[] = {static_cast<GUInt64>(nYOff),
2434
3
                                       static_cast<GUInt64>(nXOff), 0};
2435
3
            size_t count[] = {static_cast<size_t>(nYSize),
2436
3
                              static_cast<size_t>(nXSize),
2437
3
                              static_cast<size_t>(nBands)};
2438
3
            constexpr GInt64 arrayStep[] = {1, 1, 1};
2439
3
            GPtrDiff_t bufferStride[] = {
2440
3
                static_cast<GPtrDiff_t>(nLineSpace / nBufferDTSize),
2441
3
                static_cast<GPtrDiff_t>(nPixelSpace / nBufferDTSize),
2442
3
                static_cast<GPtrDiff_t>(nBandSpace / nBufferDTSize),
2443
3
            };
2444
2445
3
            if (eRWFlag == GF_Read)
2446
0
            {
2447
0
                return m_poSingleArray->Read(
2448
0
                           arrayStartIdx, count, arrayStep, bufferStride,
2449
0
                           GDALExtendedDataType::Create(eBufType), pData)
2450
0
                           ? CE_None
2451
0
                           : CE_Failure;
2452
0
            }
2453
3
            else
2454
3
            {
2455
3
                return m_poSingleArray->Write(
2456
3
                           arrayStartIdx, count, arrayStep, bufferStride,
2457
3
                           GDALExtendedDataType::Create(eBufType), pData)
2458
3
                           ? CE_None
2459
3
                           : CE_Failure;
2460
3
            }
2461
3
        }
2462
3
    }
2463
2464
9.34k
    return GDALDataset::IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize, pData,
2465
9.34k
                                  nBufXSize, nBufYSize, eBufType, nBandCount,
2466
9.34k
                                  panBandMap, nPixelSpace, nLineSpace,
2467
9.34k
                                  nBandSpace, psExtraArg);
2468
9.35k
}
2469
2470
/************************************************************************/
2471
/*                      ZarrDataset::CreateCopy()                       */
2472
/************************************************************************/
2473
2474
/* static */
2475
GDALDataset *ZarrDataset::CreateCopy(const char *pszFilename,
2476
                                     GDALDataset *poSrcDS, int bStrict,
2477
                                     CSLConstList papszOptions,
2478
                                     GDALProgressFunc pfnProgress,
2479
                                     void *pProgressData)
2480
77
{
2481
77
    if (CPLFetchBool(papszOptions, "CONVERT_TO_KERCHUNK_PARQUET_REFERENCE",
2482
77
                     false))
2483
0
    {
2484
0
        if (VSIKerchunkConvertJSONToParquet(poSrcDS->GetDescription(),
2485
0
                                            pszFilename, pfnProgress,
2486
0
                                            pProgressData))
2487
0
        {
2488
0
            GDALOpenInfo oOpenInfo(
2489
0
                std::string("ZARR:\"").append(pszFilename).append("\"").c_str(),
2490
0
                GA_ReadOnly);
2491
0
            return Open(&oOpenInfo);
2492
0
        }
2493
0
    }
2494
77
    else
2495
77
    {
2496
77
        auto poDriver = GetGDALDriverManager()->GetDriverByName(DRIVER_NAME);
2497
77
        CPLStringList aosCreationOptions(
2498
77
            const_cast<CSLConstList>(papszOptions));
2499
77
        GDALGeoTransform gt;
2500
77
        if (poSrcDS->GetGeoTransform(gt) == CE_None)
2501
38
        {
2502
38
            aosCreationOptions.SetNameValue("@HAS_GEOTRANSFORM", "YES");
2503
38
        }
2504
77
        auto poDS = std::unique_ptr<GDALDataset>(poDriver->DefaultCreateCopy(
2505
77
            pszFilename, poSrcDS, bStrict, aosCreationOptions.List(),
2506
77
            pfnProgress, pProgressData));
2507
77
        if (poDS)
2508
50
        {
2509
50
            if (poDS->FlushCache() != CE_None)
2510
0
                poDS.reset();
2511
50
        }
2512
77
        return poDS.release();
2513
77
    }
2514
0
    return nullptr;
2515
77
}
2516
2517
/************************************************************************/
2518
/*               ZARRAddGeoreferencingConventionAlgorithm               */
2519
/************************************************************************/
2520
2521
#ifndef _
2522
0
#define _(x) (x)
2523
#endif
2524
2525
namespace
2526
{
2527
class ZARRAddGeoreferencingConventionAlgorithm final : public GDALAlgorithm
2528
{
2529
  public:
2530
    ZARRAddGeoreferencingConventionAlgorithm()
2531
0
        : GDALAlgorithm(
2532
0
              "add-georeferencing-convention",
2533
0
              std::string("Add a georeferencing convention to an existing ZARR "
2534
0
                          "dataset"),
2535
0
              "/programs/gdal_driver_zarr_add_georeferencing_convention.html")
2536
0
    {
2537
0
        AddProgressArg(/* hidden = */ true);
2538
0
        AddInputDatasetArg(&m_dataset,
2539
0
                           GDAL_OF_MULTIDIM_RASTER | GDAL_OF_UPDATE);
2540
0
        AddArg("convention", 0, _("Georeferencing convention"),
2541
0
               &m_georeferencingConvention)
2542
0
            .SetRequired()
2543
0
            .SetPositional()
2544
0
            .SetChoices("GDAL", "spatial_proj");
2545
0
    }
2546
2547
  protected:
2548
    bool RunImpl(GDALProgressFunc, void *) override;
2549
2550
  private:
2551
    GDALArgDatasetValue m_dataset{};
2552
    std::string m_georeferencingConvention{};
2553
};
2554
2555
bool ZARRAddGeoreferencingConventionAlgorithm::RunImpl(GDALProgressFunc, void *)
2556
0
{
2557
0
    auto poDS = dynamic_cast<ZarrDataset *>(m_dataset.GetDatasetRef());
2558
0
    if (!poDS)
2559
0
    {
2560
0
        ReportError(CE_Failure, CPLE_AppDefined, "%s is not a ZARR dataset",
2561
0
                    m_dataset.GetName().c_str());
2562
0
        return false;
2563
0
    }
2564
2565
0
    auto poRG = poDS->GetRootGroup();
2566
0
    CPLAssert(poRG);
2567
2568
0
    poRG->RecursivelyVisitArrays(
2569
0
        [this](const std::shared_ptr<GDALMDArray> &poArray)
2570
0
        {
2571
0
            ZarrArray *poZarrArray = dynamic_cast<ZarrArray *>(poArray.get());
2572
0
            if (poZarrArray && poZarrArray->GetSpatialRef())
2573
0
            {
2574
0
                CPLStringList aosOptions;
2575
0
                aosOptions.SetNameValue("GEOREFERENCING_CONVENTION",
2576
0
                                        m_georeferencingConvention.c_str());
2577
0
                poZarrArray->SetCreationOptions(aosOptions.List());
2578
0
                poZarrArray->InvalidateGeoreferencing();
2579
0
            }
2580
0
        });
2581
2582
0
    return true;
2583
0
}
2584
}  // namespace
2585
2586
/************************************************************************/
2587
/*                   ZarrDriverInstantiateAlgorithm()                   */
2588
/************************************************************************/
2589
2590
static GDALAlgorithm *
2591
ZarrDriverInstantiateAlgorithm(const std::vector<std::string> &aosPath)
2592
0
{
2593
0
    if (aosPath.size() == 1 && aosPath[0] == "add-georeferencing-convention")
2594
0
    {
2595
0
        return std::make_unique<ZARRAddGeoreferencingConventionAlgorithm>()
2596
0
            .release();
2597
0
    }
2598
0
    else
2599
0
    {
2600
0
        return nullptr;
2601
0
    }
2602
0
}
2603
2604
/************************************************************************/
2605
/*                         GDALRegister_Zarr()                          */
2606
/************************************************************************/
2607
2608
void GDALRegister_Zarr()
2609
2610
22
{
2611
22
    if (GDALGetDriverByName(DRIVER_NAME) != nullptr)
2612
0
        return;
2613
2614
22
    VSIInstallKerchunkFileSystems();
2615
2616
22
    GDALDriver *poDriver = new ZarrDriver();
2617
22
    ZARRDriverSetCommonMetadata(poDriver);
2618
2619
#ifdef HAVE_PCODEC
2620
    poDriver->SetMetadataItem("HAVE_PCODEC", "YES");
2621
#endif
2622
2623
22
    poDriver->pfnOpen = ZarrDataset::Open;
2624
22
    poDriver->pfnCreateMultiDimensional = ZarrDataset::CreateMultiDimensional;
2625
22
    poDriver->pfnCreate = ZarrDataset::Create;
2626
22
    poDriver->pfnCreateCopy = ZarrDataset::CreateCopy;
2627
22
    poDriver->pfnDelete = ZarrDatasetDelete;
2628
22
    poDriver->pfnRename = ZarrDatasetRename;
2629
22
    poDriver->pfnCopyFiles = ZarrDatasetCopyFiles;
2630
22
    poDriver->pfnClearCaches = ZarrDriverClearCaches;
2631
22
    poDriver->pfnInstantiateAlgorithm = ZarrDriverInstantiateAlgorithm;
2632
2633
22
    GetGDALDriverManager()->RegisterDriver(poDriver);
2634
22
}