Line data Source code
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 "ucs4_utf8.hpp"
15 :
16 : #include "cpl_float.h"
17 : #include "cpl_multiproc.h"
18 : #include "cpl_vsi_virtual.h"
19 :
20 : #include "netcdf_cf_constants.h" // for CF_UNITS, etc
21 :
22 : #include <algorithm>
23 : #include <cassert>
24 : #include <cmath>
25 : #include <cstdlib>
26 : #include <limits>
27 : #include <map>
28 : #include <set>
29 :
30 : #if defined(__clang__) || defined(_MSC_VER)
31 : #define COMPILER_WARNS_ABOUT_ABSTRACT_VBASE_INIT
32 : #endif
33 :
34 : namespace
35 : {
36 :
37 4 : inline std::vector<GByte> UTF8ToUCS4(const char *pszStr, bool needByteSwap)
38 : {
39 4 : const size_t nLen = strlen(pszStr);
40 : // Worst case if that we need 4 more bytes than the UTF-8 one
41 : // (when the content is pure ASCII)
42 4 : if (nLen > std::numeric_limits<size_t>::max() / sizeof(uint32_t))
43 0 : throw std::bad_alloc();
44 4 : std::vector<GByte> ret(nLen * sizeof(uint32_t));
45 4 : size_t outPos = 0;
46 27 : for (size_t i = 0; i < nLen; outPos += sizeof(uint32_t))
47 : {
48 23 : uint32_t ucs4 = 0;
49 23 : int consumed = FcUtf8ToUcs4(
50 23 : reinterpret_cast<const uint8_t *>(pszStr + i), &ucs4, nLen - i);
51 23 : if (consumed <= 0)
52 : {
53 0 : ret.resize(outPos);
54 : }
55 23 : if (needByteSwap)
56 : {
57 1 : CPL_SWAP32PTR(&ucs4);
58 : }
59 23 : memcpy(&ret[outPos], &ucs4, sizeof(uint32_t));
60 23 : i += consumed;
61 : }
62 4 : ret.resize(outPos);
63 4 : return ret;
64 : }
65 :
66 6 : inline char *UCS4ToUTF8(const uint8_t *ucs4Ptr, size_t nSize, bool needByteSwap)
67 : {
68 : // A UCS4 char can require up to 6 bytes in UTF8.
69 6 : if (nSize > (std::numeric_limits<size_t>::max() - 1) / 6 * 4)
70 0 : return nullptr;
71 6 : const size_t nOutSize = nSize / 4 * 6 + 1;
72 6 : char *ret = static_cast<char *>(VSI_MALLOC_VERBOSE(nOutSize));
73 6 : if (ret == nullptr)
74 0 : return nullptr;
75 6 : size_t outPos = 0;
76 30 : for (size_t i = 0; i + sizeof(uint32_t) - 1 < nSize; i += sizeof(uint32_t))
77 : {
78 : uint32_t ucs4;
79 24 : memcpy(&ucs4, ucs4Ptr + i, sizeof(uint32_t));
80 24 : if (needByteSwap)
81 : {
82 2 : CPL_SWAP32PTR(&ucs4);
83 : }
84 : int written =
85 24 : FcUcs4ToUtf8(ucs4, reinterpret_cast<uint8_t *>(ret + outPos));
86 24 : outPos += written;
87 : }
88 6 : ret[outPos] = 0;
89 6 : return ret;
90 : }
91 :
92 : } // namespace
93 :
94 : /************************************************************************/
95 : /* ZarrArray::ParseChunkSize() */
96 : /************************************************************************/
97 :
98 872 : /* static */ bool ZarrArray::ParseChunkSize(const CPLJSONArray &oChunks,
99 : const GDALExtendedDataType &oType,
100 : std::vector<GUInt64> &anBlockSize)
101 : {
102 872 : size_t nBlockSize = oType.GetSize();
103 2217 : for (const auto &item : oChunks)
104 : {
105 1348 : const auto nSize = static_cast<GUInt64>(item.ToLong());
106 1348 : if (nSize == 0)
107 : {
108 1 : CPLError(CE_Failure, CPLE_AppDefined, "Invalid content for chunks");
109 3 : return false;
110 : }
111 1347 : if (nBlockSize > std::numeric_limits<size_t>::max() / nSize)
112 : {
113 2 : CPLError(CE_Failure, CPLE_AppDefined, "Too large chunks");
114 2 : return false;
115 : }
116 1345 : nBlockSize *= static_cast<size_t>(nSize);
117 1345 : anBlockSize.emplace_back(nSize);
118 : }
119 :
120 869 : return true;
121 : }
122 :
123 : /************************************************************************/
124 : /* ZarrArray::ComputeTileCount() */
125 : /************************************************************************/
126 :
127 1278 : /* static */ uint64_t ZarrArray::ComputeTileCount(
128 : const std::string &osName,
129 : const std::vector<std::shared_ptr<GDALDimension>> &aoDims,
130 : const std::vector<GUInt64> &anBlockSize)
131 : {
132 1278 : uint64_t nTotalTileCount = 1;
133 3226 : for (size_t i = 0; i < aoDims.size(); ++i)
134 : {
135 : uint64_t nTileThisDim =
136 1950 : (aoDims[i]->GetSize() / anBlockSize[i]) +
137 1950 : (((aoDims[i]->GetSize() % anBlockSize[i]) != 0) ? 1 : 0);
138 3900 : if (nTileThisDim != 0 &&
139 : nTotalTileCount >
140 1950 : std::numeric_limits<uint64_t>::max() / nTileThisDim)
141 : {
142 2 : CPLError(
143 : CE_Failure, CPLE_NotSupported,
144 : "Array %s has more than 2^64 tiles. This is not supported.",
145 : osName.c_str());
146 2 : return 0;
147 : }
148 1948 : nTotalTileCount *= nTileThisDim;
149 : }
150 1276 : return nTotalTileCount;
151 : }
152 :
153 : /************************************************************************/
154 : /* ZarrArray::ZarrArray() */
155 : /************************************************************************/
156 :
157 1278 : ZarrArray::ZarrArray(
158 : const std::shared_ptr<ZarrSharedResource> &poSharedResource,
159 : const std::string &osParentName, const std::string &osName,
160 : const std::vector<std::shared_ptr<GDALDimension>> &aoDims,
161 : const GDALExtendedDataType &oType, const std::vector<DtypeElt> &aoDtypeElts,
162 0 : const std::vector<GUInt64> &anBlockSize)
163 : :
164 : #if !defined(COMPILER_WARNS_ABOUT_ABSTRACT_VBASE_INIT)
165 : GDALAbstractMDArray(osParentName, osName),
166 : #endif
167 : GDALPamMDArray(osParentName, osName, poSharedResource->GetPAM()),
168 : m_poSharedResource(poSharedResource), m_aoDims(aoDims), m_oType(oType),
169 : m_aoDtypeElts(aoDtypeElts), m_anBlockSize(anBlockSize),
170 1278 : m_oAttrGroup(m_osFullName, /*bContainerIsGroup=*/false)
171 : {
172 1278 : m_nTotalTileCount = ComputeTileCount(osName, aoDims, anBlockSize);
173 1278 : if (m_nTotalTileCount == 0)
174 2 : return;
175 :
176 : // Compute individual tile size
177 : const size_t nSourceSize =
178 1276 : m_aoDtypeElts.back().nativeOffset + m_aoDtypeElts.back().nativeSize;
179 1276 : m_nTileSize = nSourceSize;
180 3220 : for (const auto &nBlockSize : m_anBlockSize)
181 : {
182 1944 : m_nTileSize *= static_cast<size_t>(nBlockSize);
183 : }
184 :
185 1276 : m_bUseOptimizedCodePaths = CPLTestBool(
186 : CPLGetConfigOption("GDAL_ZARR_USE_OPTIMIZED_CODE_PATHS", "YES"));
187 : }
188 :
189 : /************************************************************************/
190 : /* ~ZarrArray() */
191 : /************************************************************************/
192 :
193 1278 : ZarrArray::~ZarrArray()
194 : {
195 1278 : if (m_pabyNoData)
196 : {
197 206 : m_oType.FreeDynamicMemory(&m_pabyNoData[0]);
198 206 : CPLFree(m_pabyNoData);
199 : }
200 :
201 1278 : DeallocateDecodedTileData();
202 1278 : }
203 :
204 : /************************************************************************/
205 : /* ZarrArray::SerializeSpecialAttributes() */
206 : /************************************************************************/
207 :
208 390 : CPLJSONObject ZarrArray::SerializeSpecialAttributes()
209 : {
210 390 : m_bSRSModified = false;
211 390 : m_oAttrGroup.UnsetModified();
212 :
213 390 : auto oAttrs = m_oAttrGroup.Serialize();
214 :
215 390 : if (m_poSRS)
216 : {
217 40 : CPLJSONObject oCRS;
218 40 : const char *const apszOptions[] = {"FORMAT=WKT2_2019", nullptr};
219 40 : char *pszWKT = nullptr;
220 40 : if (m_poSRS->exportToWkt(&pszWKT, apszOptions) == OGRERR_NONE)
221 : {
222 40 : oCRS.Add("wkt", pszWKT);
223 : }
224 40 : CPLFree(pszWKT);
225 :
226 : {
227 80 : CPLErrorStateBackuper oErrorStateBackuper(CPLQuietErrorHandler);
228 40 : char *projjson = nullptr;
229 80 : if (m_poSRS->exportToPROJJSON(&projjson, nullptr) == OGRERR_NONE &&
230 40 : projjson != nullptr)
231 : {
232 80 : CPLJSONDocument oDocProjJSON;
233 40 : if (oDocProjJSON.LoadMemory(std::string(projjson)))
234 : {
235 40 : oCRS.Add("projjson", oDocProjJSON.GetRoot());
236 : }
237 : }
238 40 : CPLFree(projjson);
239 : }
240 :
241 40 : const char *pszAuthorityCode = m_poSRS->GetAuthorityCode(nullptr);
242 40 : const char *pszAuthorityName = m_poSRS->GetAuthorityName(nullptr);
243 40 : if (pszAuthorityCode && pszAuthorityName &&
244 7 : EQUAL(pszAuthorityName, "EPSG"))
245 : {
246 7 : oCRS.Add("url",
247 14 : std::string("http://www.opengis.net/def/crs/EPSG/0/") +
248 : pszAuthorityCode);
249 : }
250 :
251 40 : oAttrs.Add(CRS_ATTRIBUTE_NAME, oCRS);
252 : }
253 :
254 390 : if (m_osUnit.empty())
255 : {
256 381 : if (m_bUnitModified)
257 0 : oAttrs.Delete(CF_UNITS);
258 : }
259 : else
260 : {
261 9 : oAttrs.Set(CF_UNITS, m_osUnit);
262 : }
263 390 : m_bUnitModified = false;
264 :
265 390 : if (!m_bHasOffset)
266 : {
267 387 : oAttrs.Delete(CF_ADD_OFFSET);
268 : }
269 : else
270 : {
271 3 : oAttrs.Set(CF_ADD_OFFSET, m_dfOffset);
272 : }
273 390 : m_bOffsetModified = false;
274 :
275 390 : if (!m_bHasScale)
276 : {
277 387 : oAttrs.Delete(CF_SCALE_FACTOR);
278 : }
279 : else
280 : {
281 3 : oAttrs.Set(CF_SCALE_FACTOR, m_dfScale);
282 : }
283 390 : m_bScaleModified = false;
284 :
285 390 : return oAttrs;
286 : }
287 :
288 : /************************************************************************/
289 : /* FillBlockSize() */
290 : /************************************************************************/
291 :
292 : /* static */
293 446 : bool ZarrArray::FillBlockSize(
294 : const std::vector<std::shared_ptr<GDALDimension>> &aoDimensions,
295 : const GDALExtendedDataType &oDataType, std::vector<GUInt64> &anBlockSize,
296 : CSLConstList papszOptions)
297 : {
298 446 : const auto nDims = aoDimensions.size();
299 446 : anBlockSize.resize(nDims);
300 1108 : for (size_t i = 0; i < nDims; ++i)
301 662 : anBlockSize[i] = 1;
302 446 : if (nDims >= 2)
303 : {
304 436 : anBlockSize[nDims - 2] =
305 436 : std::min(std::max<GUInt64>(1, aoDimensions[nDims - 2]->GetSize()),
306 436 : static_cast<GUInt64>(256));
307 436 : anBlockSize[nDims - 1] =
308 436 : std::min(std::max<GUInt64>(1, aoDimensions[nDims - 1]->GetSize()),
309 654 : static_cast<GUInt64>(256));
310 : }
311 228 : else if (nDims == 1)
312 : {
313 181 : anBlockSize[0] = std::max<GUInt64>(1, aoDimensions[0]->GetSize());
314 : }
315 :
316 446 : const char *pszBlockSize = CSLFetchNameValue(papszOptions, "BLOCKSIZE");
317 446 : if (pszBlockSize)
318 : {
319 : const auto aszTokens(
320 18 : CPLStringList(CSLTokenizeString2(pszBlockSize, ",", 0)));
321 18 : if (static_cast<size_t>(aszTokens.size()) != nDims)
322 : {
323 5 : CPLError(CE_Failure, CPLE_AppDefined,
324 : "Invalid number of values in BLOCKSIZE");
325 5 : return false;
326 : }
327 13 : size_t nBlockSize = oDataType.GetSize();
328 40 : for (size_t i = 0; i < nDims; ++i)
329 : {
330 28 : const auto v = static_cast<GUInt64>(CPLAtoGIntBig(aszTokens[i]));
331 28 : if (v > 0)
332 : {
333 28 : anBlockSize[i] = v;
334 : }
335 28 : if (anBlockSize[i] >
336 28 : std::numeric_limits<size_t>::max() / nBlockSize)
337 : {
338 1 : CPLError(CE_Failure, CPLE_AppDefined,
339 : "Too large values in BLOCKSIZE");
340 1 : return false;
341 : }
342 27 : nBlockSize *= static_cast<size_t>(anBlockSize[i]);
343 : }
344 : }
345 440 : return true;
346 : }
347 :
348 : /************************************************************************/
349 : /* DeallocateDecodedTileData() */
350 : /************************************************************************/
351 :
352 2781 : void ZarrArray::DeallocateDecodedTileData()
353 : {
354 2781 : if (!m_abyDecodedTileData.empty())
355 : {
356 223 : const size_t nDTSize = m_oType.GetSize();
357 223 : GByte *pDst = &m_abyDecodedTileData[0];
358 223 : const size_t nValues = m_abyDecodedTileData.size() / nDTSize;
359 454 : for (const auto &elt : m_aoDtypeElts)
360 : {
361 231 : if (elt.nativeType == DtypeElt::NativeType::STRING_ASCII ||
362 230 : elt.nativeType == DtypeElt::NativeType::STRING_UNICODE)
363 : {
364 2 : for (size_t i = 0; i < nValues; i++, pDst += nDTSize)
365 : {
366 : char *ptr;
367 1 : char **pptr =
368 1 : reinterpret_cast<char **>(pDst + elt.gdalOffset);
369 1 : memcpy(&ptr, pptr, sizeof(ptr));
370 1 : VSIFree(ptr);
371 : }
372 : }
373 : }
374 : }
375 2781 : }
376 :
377 : /************************************************************************/
378 : /* EncodeElt() */
379 : /************************************************************************/
380 :
381 : /* Encode from GDAL raw type to Zarr native type */
382 : /*static*/
383 121 : void ZarrArray::EncodeElt(const std::vector<DtypeElt> &elts, const GByte *pSrc,
384 : GByte *pDst)
385 : {
386 243 : for (const auto &elt : elts)
387 : {
388 122 : if (elt.nativeType == DtypeElt::NativeType::STRING_UNICODE)
389 : {
390 0 : const char *pStr =
391 0 : *reinterpret_cast<const char *const *>(pSrc + elt.gdalOffset);
392 0 : if (pStr)
393 : {
394 : try
395 : {
396 0 : const auto ucs4 = UTF8ToUCS4(pStr, elt.needByteSwapping);
397 0 : const auto ucs4Len = ucs4.size();
398 0 : memcpy(pDst + elt.nativeOffset, ucs4.data(),
399 0 : std::min(ucs4Len, elt.nativeSize));
400 0 : if (ucs4Len > elt.nativeSize)
401 : {
402 0 : CPLError(CE_Warning, CPLE_AppDefined,
403 : "Too long string truncated");
404 : }
405 0 : else if (ucs4Len < elt.nativeSize)
406 : {
407 0 : memset(pDst + elt.nativeOffset + ucs4Len, 0,
408 0 : elt.nativeSize - ucs4Len);
409 : }
410 : }
411 0 : catch (const std::exception &)
412 : {
413 0 : memset(pDst + elt.nativeOffset, 0, elt.nativeSize);
414 : }
415 : }
416 : else
417 : {
418 0 : memset(pDst + elt.nativeOffset, 0, elt.nativeSize);
419 : }
420 : }
421 122 : else if (elt.needByteSwapping)
422 : {
423 120 : if (elt.nativeSize == 2)
424 : {
425 24 : if (elt.gdalTypeIsApproxOfNative)
426 : {
427 0 : CPLAssert(elt.nativeType == DtypeElt::NativeType::IEEEFP);
428 0 : CPLAssert(elt.gdalType.GetNumericDataType() == GDT_Float32);
429 0 : const uint32_t uint32Val =
430 0 : *reinterpret_cast<const uint32_t *>(pSrc +
431 0 : elt.gdalOffset);
432 0 : bool bHasWarned = false;
433 : uint16_t uint16Val =
434 0 : CPL_SWAP16(CPLFloatToHalf(uint32Val, bHasWarned));
435 0 : memcpy(pDst + elt.nativeOffset, &uint16Val,
436 : sizeof(uint16Val));
437 : }
438 : else
439 : {
440 24 : const uint16_t val =
441 24 : CPL_SWAP16(*reinterpret_cast<const uint16_t *>(
442 : pSrc + elt.gdalOffset));
443 24 : memcpy(pDst + elt.nativeOffset, &val, sizeof(val));
444 : }
445 : }
446 96 : else if (elt.nativeSize == 4)
447 : {
448 36 : const uint32_t val = CPL_SWAP32(
449 : *reinterpret_cast<const uint32_t *>(pSrc + elt.gdalOffset));
450 36 : memcpy(pDst + elt.nativeOffset, &val, sizeof(val));
451 : }
452 60 : else if (elt.nativeSize == 8)
453 : {
454 48 : if (elt.nativeType == DtypeElt::NativeType::COMPLEX_IEEEFP)
455 : {
456 12 : uint32_t val =
457 12 : CPL_SWAP32(*reinterpret_cast<const uint32_t *>(
458 : pSrc + elt.gdalOffset));
459 12 : memcpy(pDst + elt.nativeOffset, &val, sizeof(val));
460 12 : val = CPL_SWAP32(*reinterpret_cast<const uint32_t *>(
461 : pSrc + elt.gdalOffset + 4));
462 12 : memcpy(pDst + elt.nativeOffset + 4, &val, sizeof(val));
463 : }
464 : else
465 : {
466 36 : const uint64_t val =
467 36 : CPL_SWAP64(*reinterpret_cast<const uint64_t *>(
468 : pSrc + elt.gdalOffset));
469 36 : memcpy(pDst + elt.nativeOffset, &val, sizeof(val));
470 : }
471 : }
472 12 : else if (elt.nativeSize == 16)
473 : {
474 12 : uint64_t val = CPL_SWAP64(
475 : *reinterpret_cast<const uint64_t *>(pSrc + elt.gdalOffset));
476 12 : memcpy(pDst + elt.nativeOffset, &val, sizeof(val));
477 12 : val = CPL_SWAP64(*reinterpret_cast<const uint64_t *>(
478 : pSrc + elt.gdalOffset + 8));
479 12 : memcpy(pDst + elt.nativeOffset + 8, &val, sizeof(val));
480 : }
481 : else
482 : {
483 0 : CPLAssert(false);
484 : }
485 : }
486 2 : else if (elt.gdalTypeIsApproxOfNative)
487 : {
488 0 : if (elt.nativeType == DtypeElt::NativeType::IEEEFP &&
489 0 : elt.nativeSize == 2)
490 : {
491 0 : CPLAssert(elt.gdalType.GetNumericDataType() == GDT_Float32);
492 0 : const uint32_t uint32Val =
493 0 : *reinterpret_cast<const uint32_t *>(pSrc + elt.gdalOffset);
494 0 : bool bHasWarned = false;
495 : const uint16_t uint16Val =
496 0 : CPLFloatToHalf(uint32Val, bHasWarned);
497 0 : memcpy(pDst + elt.nativeOffset, &uint16Val, sizeof(uint16Val));
498 : }
499 : else
500 : {
501 0 : CPLAssert(false);
502 : }
503 : }
504 2 : else if (elt.nativeType == DtypeElt::NativeType::STRING_ASCII)
505 : {
506 0 : const char *pStr =
507 0 : *reinterpret_cast<const char *const *>(pSrc + elt.gdalOffset);
508 0 : if (pStr)
509 : {
510 0 : const size_t nLen = strlen(pStr);
511 0 : memcpy(pDst + elt.nativeOffset, pStr,
512 0 : std::min(nLen, elt.nativeSize));
513 0 : if (nLen < elt.nativeSize)
514 0 : memset(pDst + elt.nativeOffset + nLen, 0,
515 0 : elt.nativeSize - nLen);
516 : }
517 : else
518 : {
519 0 : memset(pDst + elt.nativeOffset, 0, elt.nativeSize);
520 : }
521 : }
522 : else
523 : {
524 2 : CPLAssert(elt.nativeSize == elt.gdalSize);
525 2 : memcpy(pDst + elt.nativeOffset, pSrc + elt.gdalOffset,
526 2 : elt.nativeSize);
527 : }
528 : }
529 121 : }
530 :
531 : /************************************************************************/
532 : /* ZarrArray::SerializeNumericNoData() */
533 : /************************************************************************/
534 :
535 16 : void ZarrArray::SerializeNumericNoData(CPLJSONObject &oRoot) const
536 : {
537 16 : if (m_oType.GetNumericDataType() == GDT_Int64)
538 : {
539 0 : const auto nVal = GetNoDataValueAsInt64();
540 0 : oRoot.Add("fill_value", static_cast<GInt64>(nVal));
541 : }
542 16 : else if (m_oType.GetNumericDataType() == GDT_UInt64)
543 : {
544 0 : const auto nVal = GetNoDataValueAsUInt64();
545 0 : oRoot.Add("fill_value", static_cast<uint64_t>(nVal));
546 : }
547 : else
548 : {
549 16 : const double dfVal = GetNoDataValueAsDouble();
550 16 : if (std::isnan(dfVal))
551 2 : oRoot.Add("fill_value", "NaN");
552 14 : else if (dfVal == std::numeric_limits<double>::infinity())
553 2 : oRoot.Add("fill_value", "Infinity");
554 12 : else if (dfVal == -std::numeric_limits<double>::infinity())
555 2 : oRoot.Add("fill_value", "-Infinity");
556 10 : else if (GDALDataTypeIsInteger(m_oType.GetNumericDataType()))
557 8 : oRoot.Add("fill_value", static_cast<GInt64>(dfVal));
558 : else
559 2 : oRoot.Add("fill_value", dfVal);
560 : }
561 16 : }
562 :
563 : /************************************************************************/
564 : /* ZarrArray::GetSpatialRef() */
565 : /************************************************************************/
566 :
567 22 : std::shared_ptr<OGRSpatialReference> ZarrArray::GetSpatialRef() const
568 : {
569 22 : if (!CheckValidAndErrorOutIfNot())
570 0 : return nullptr;
571 :
572 22 : if (m_poSRS)
573 12 : return m_poSRS;
574 10 : return GDALPamMDArray::GetSpatialRef();
575 : }
576 :
577 : /************************************************************************/
578 : /* SetRawNoDataValue() */
579 : /************************************************************************/
580 :
581 18 : bool ZarrArray::SetRawNoDataValue(const void *pRawNoData)
582 : {
583 18 : if (!CheckValidAndErrorOutIfNot())
584 0 : return false;
585 :
586 18 : if (!m_bUpdatable)
587 : {
588 0 : CPLError(CE_Failure, CPLE_AppDefined, "Array opened in read-only mode");
589 0 : return false;
590 : }
591 18 : m_bDefinitionModified = true;
592 18 : RegisterNoDataValue(pRawNoData);
593 18 : return true;
594 : }
595 :
596 : /************************************************************************/
597 : /* RegisterNoDataValue() */
598 : /************************************************************************/
599 :
600 206 : void ZarrArray::RegisterNoDataValue(const void *pNoData)
601 : {
602 206 : if (m_pabyNoData)
603 : {
604 0 : m_oType.FreeDynamicMemory(&m_pabyNoData[0]);
605 : }
606 :
607 206 : if (pNoData == nullptr)
608 : {
609 0 : CPLFree(m_pabyNoData);
610 0 : m_pabyNoData = nullptr;
611 : }
612 : else
613 : {
614 206 : const auto nSize = m_oType.GetSize();
615 206 : if (m_pabyNoData == nullptr)
616 : {
617 206 : m_pabyNoData = static_cast<GByte *>(CPLMalloc(nSize));
618 : }
619 206 : memset(m_pabyNoData, 0, nSize);
620 206 : GDALExtendedDataType::CopyValue(pNoData, m_oType, m_pabyNoData,
621 206 : m_oType);
622 : }
623 206 : }
624 :
625 : /************************************************************************/
626 : /* ZarrArray::BlockTranspose() */
627 : /************************************************************************/
628 :
629 66 : void ZarrArray::BlockTranspose(const ZarrByteVectorQuickResize &abySrc,
630 : ZarrByteVectorQuickResize &abyDst,
631 : bool bDecode) const
632 : {
633 : // Perform transposition
634 66 : const size_t nDims = m_anBlockSize.size();
635 : const size_t nSourceSize =
636 66 : m_aoDtypeElts.back().nativeOffset + m_aoDtypeElts.back().nativeSize;
637 :
638 : struct Stack
639 : {
640 : size_t nIters = 0;
641 : const GByte *src_ptr = nullptr;
642 : GByte *dst_ptr = nullptr;
643 : size_t src_inc_offset = 0;
644 : size_t dst_inc_offset = 0;
645 : };
646 :
647 132 : std::vector<Stack> stack(nDims);
648 : stack.emplace_back(
649 66 : Stack()); // to make gcc 9.3 -O2 -Wnull-dereference happy
650 :
651 66 : if (bDecode)
652 : {
653 46 : stack[0].src_inc_offset = nSourceSize;
654 118 : for (size_t i = 1; i < nDims; ++i)
655 : {
656 144 : stack[i].src_inc_offset = stack[i - 1].src_inc_offset *
657 72 : static_cast<size_t>(m_anBlockSize[i - 1]);
658 : }
659 :
660 46 : stack[nDims - 1].dst_inc_offset = nSourceSize;
661 118 : for (size_t i = nDims - 1; i > 0;)
662 : {
663 72 : --i;
664 144 : stack[i].dst_inc_offset = stack[i + 1].dst_inc_offset *
665 72 : static_cast<size_t>(m_anBlockSize[i + 1]);
666 : }
667 : }
668 : else
669 : {
670 20 : stack[0].dst_inc_offset = nSourceSize;
671 60 : for (size_t i = 1; i < nDims; ++i)
672 : {
673 80 : stack[i].dst_inc_offset = stack[i - 1].dst_inc_offset *
674 40 : static_cast<size_t>(m_anBlockSize[i - 1]);
675 : }
676 :
677 20 : stack[nDims - 1].src_inc_offset = nSourceSize;
678 60 : for (size_t i = nDims - 1; i > 0;)
679 : {
680 40 : --i;
681 80 : stack[i].src_inc_offset = stack[i + 1].src_inc_offset *
682 40 : static_cast<size_t>(m_anBlockSize[i + 1]);
683 : }
684 : }
685 :
686 66 : stack[0].src_ptr = abySrc.data();
687 66 : stack[0].dst_ptr = &abyDst[0];
688 :
689 66 : size_t dimIdx = 0;
690 1058 : lbl_next_depth:
691 1058 : if (dimIdx == nDims)
692 : {
693 744 : void *dst_ptr = stack[nDims].dst_ptr;
694 744 : const void *src_ptr = stack[nDims].src_ptr;
695 744 : if (nSourceSize == 1)
696 264 : *stack[nDims].dst_ptr = *stack[nDims].src_ptr;
697 480 : else if (nSourceSize == 2)
698 120 : *static_cast<uint16_t *>(dst_ptr) =
699 120 : *static_cast<const uint16_t *>(src_ptr);
700 360 : else if (nSourceSize == 4)
701 168 : *static_cast<uint32_t *>(dst_ptr) =
702 168 : *static_cast<const uint32_t *>(src_ptr);
703 192 : else if (nSourceSize == 8)
704 168 : *static_cast<uint64_t *>(dst_ptr) =
705 168 : *static_cast<const uint64_t *>(src_ptr);
706 : else
707 24 : memcpy(dst_ptr, src_ptr, nSourceSize);
708 : }
709 : else
710 : {
711 314 : stack[dimIdx].nIters = static_cast<size_t>(m_anBlockSize[dimIdx]);
712 : while (true)
713 : {
714 992 : dimIdx++;
715 992 : stack[dimIdx].src_ptr = stack[dimIdx - 1].src_ptr;
716 992 : stack[dimIdx].dst_ptr = stack[dimIdx - 1].dst_ptr;
717 992 : goto lbl_next_depth;
718 992 : lbl_return_to_caller:
719 992 : dimIdx--;
720 992 : if ((--stack[dimIdx].nIters) == 0)
721 314 : break;
722 678 : stack[dimIdx].src_ptr += stack[dimIdx].src_inc_offset;
723 678 : stack[dimIdx].dst_ptr += stack[dimIdx].dst_inc_offset;
724 : }
725 : }
726 1058 : if (dimIdx > 0)
727 992 : goto lbl_return_to_caller;
728 66 : }
729 :
730 : /************************************************************************/
731 : /* DecodeSourceElt() */
732 : /************************************************************************/
733 :
734 : /* static */
735 1936 : void ZarrArray::DecodeSourceElt(const std::vector<DtypeElt> &elts,
736 : const GByte *pSrc, GByte *pDst)
737 : {
738 3904 : for (const auto &elt : elts)
739 : {
740 1968 : if (elt.nativeType == DtypeElt::NativeType::STRING_UNICODE)
741 : {
742 : char *ptr;
743 0 : char **pDstPtr = reinterpret_cast<char **>(pDst + elt.gdalOffset);
744 0 : memcpy(&ptr, pDstPtr, sizeof(ptr));
745 0 : VSIFree(ptr);
746 :
747 0 : char *pDstStr = UCS4ToUTF8(pSrc + elt.nativeOffset, elt.nativeSize,
748 0 : elt.needByteSwapping);
749 0 : memcpy(pDstPtr, &pDstStr, sizeof(pDstStr));
750 : }
751 1968 : else if (elt.needByteSwapping)
752 : {
753 1922 : if (elt.nativeSize == 2)
754 : {
755 : uint16_t val;
756 386 : memcpy(&val, pSrc + elt.nativeOffset, sizeof(val));
757 386 : if (elt.gdalTypeIsApproxOfNative)
758 : {
759 0 : CPLAssert(elt.nativeType == DtypeElt::NativeType::IEEEFP);
760 0 : CPLAssert(elt.gdalType.GetNumericDataType() == GDT_Float32);
761 0 : uint32_t uint32Val = CPLHalfToFloat(CPL_SWAP16(val));
762 0 : memcpy(pDst + elt.gdalOffset, &uint32Val,
763 : sizeof(uint32Val));
764 : }
765 : else
766 : {
767 386 : *reinterpret_cast<uint16_t *>(pDst + elt.gdalOffset) =
768 386 : CPL_SWAP16(val);
769 : }
770 : }
771 1536 : else if (elt.nativeSize == 4)
772 : {
773 : uint32_t val;
774 576 : memcpy(&val, pSrc + elt.nativeOffset, sizeof(val));
775 576 : *reinterpret_cast<uint32_t *>(pDst + elt.gdalOffset) =
776 576 : CPL_SWAP32(val);
777 : }
778 960 : else if (elt.nativeSize == 8)
779 : {
780 768 : if (elt.nativeType == DtypeElt::NativeType::COMPLEX_IEEEFP)
781 : {
782 : uint32_t val;
783 192 : memcpy(&val, pSrc + elt.nativeOffset, sizeof(val));
784 192 : *reinterpret_cast<uint32_t *>(pDst + elt.gdalOffset) =
785 192 : CPL_SWAP32(val);
786 192 : memcpy(&val, pSrc + elt.nativeOffset + 4, sizeof(val));
787 192 : *reinterpret_cast<uint32_t *>(pDst + elt.gdalOffset + 4) =
788 192 : CPL_SWAP32(val);
789 : }
790 : else
791 : {
792 : uint64_t val;
793 576 : memcpy(&val, pSrc + elt.nativeOffset, sizeof(val));
794 576 : *reinterpret_cast<uint64_t *>(pDst + elt.gdalOffset) =
795 576 : CPL_SWAP64(val);
796 : }
797 : }
798 192 : else if (elt.nativeSize == 16)
799 : {
800 : uint64_t val;
801 192 : memcpy(&val, pSrc + elt.nativeOffset, sizeof(val));
802 192 : *reinterpret_cast<uint64_t *>(pDst + elt.gdalOffset) =
803 192 : CPL_SWAP64(val);
804 192 : memcpy(&val, pSrc + elt.nativeOffset + 8, sizeof(val));
805 192 : *reinterpret_cast<uint64_t *>(pDst + elt.gdalOffset + 8) =
806 192 : CPL_SWAP64(val);
807 : }
808 : else
809 : {
810 0 : CPLAssert(false);
811 : }
812 : }
813 46 : else if (elt.gdalTypeIsApproxOfNative)
814 : {
815 0 : if (elt.nativeType == DtypeElt::NativeType::IEEEFP &&
816 0 : elt.nativeSize == 2)
817 : {
818 0 : CPLAssert(elt.gdalType.GetNumericDataType() == GDT_Float32);
819 : uint16_t uint16Val;
820 0 : memcpy(&uint16Val, pSrc + elt.nativeOffset, sizeof(uint16Val));
821 0 : uint32_t uint32Val = CPLHalfToFloat(uint16Val);
822 0 : memcpy(pDst + elt.gdalOffset, &uint32Val, sizeof(uint32Val));
823 : }
824 : else
825 : {
826 0 : CPLAssert(false);
827 : }
828 : }
829 46 : else if (elt.nativeType == DtypeElt::NativeType::STRING_ASCII)
830 : {
831 : char *ptr;
832 3 : char **pDstPtr = reinterpret_cast<char **>(pDst + elt.gdalOffset);
833 3 : memcpy(&ptr, pDstPtr, sizeof(ptr));
834 3 : VSIFree(ptr);
835 :
836 3 : char *pDstStr = static_cast<char *>(CPLMalloc(elt.nativeSize + 1));
837 3 : memcpy(pDstStr, pSrc + elt.nativeOffset, elt.nativeSize);
838 3 : pDstStr[elt.nativeSize] = 0;
839 3 : memcpy(pDstPtr, &pDstStr, sizeof(pDstStr));
840 : }
841 : else
842 : {
843 43 : CPLAssert(elt.nativeSize == elt.gdalSize);
844 43 : memcpy(pDst + elt.gdalOffset, pSrc + elt.nativeOffset,
845 43 : elt.nativeSize);
846 : }
847 : }
848 1936 : }
849 :
850 : /************************************************************************/
851 : /* ZarrArray::IAdviseReadCommon() */
852 : /************************************************************************/
853 :
854 12 : bool ZarrArray::IAdviseReadCommon(const GUInt64 *arrayStartIdx,
855 : const size_t *count,
856 : CSLConstList papszOptions,
857 : std::vector<uint64_t> &anIndicesCur,
858 : int &nThreadsMax,
859 : std::vector<uint64_t> &anReqTilesIndices,
860 : size_t &nReqTiles) const
861 : {
862 12 : if (!CheckValidAndErrorOutIfNot())
863 0 : return false;
864 :
865 12 : const size_t nDims = m_aoDims.size();
866 12 : anIndicesCur.resize(nDims);
867 24 : std::vector<uint64_t> anIndicesMin(nDims);
868 24 : std::vector<uint64_t> anIndicesMax(nDims);
869 :
870 : // Compute min and max tile indices in each dimension, and the total
871 : // number of tiles this represents.
872 12 : nReqTiles = 1;
873 36 : for (size_t i = 0; i < nDims; ++i)
874 : {
875 24 : anIndicesMin[i] = arrayStartIdx[i] / m_anBlockSize[i];
876 24 : anIndicesMax[i] = (arrayStartIdx[i] + count[i] - 1) / m_anBlockSize[i];
877 : // Overflow on number of tiles already checked in Create()
878 24 : nReqTiles *= static_cast<size_t>(anIndicesMax[i] - anIndicesMin[i] + 1);
879 : }
880 :
881 : // Find available cache size
882 12 : const size_t nCacheSize = [papszOptions]()
883 : {
884 : size_t nCacheSizeTmp;
885 : const char *pszCacheSize =
886 12 : CSLFetchNameValue(papszOptions, "CACHE_SIZE");
887 12 : if (pszCacheSize)
888 : {
889 4 : const auto nCacheSizeBig = CPLAtoGIntBig(pszCacheSize);
890 8 : if (nCacheSizeBig < 0 || static_cast<uint64_t>(nCacheSizeBig) >
891 4 : std::numeric_limits<size_t>::max() / 2)
892 : {
893 0 : CPLError(CE_Failure, CPLE_OutOfMemory, "Too big CACHE_SIZE");
894 0 : return std::numeric_limits<size_t>::max();
895 : }
896 4 : nCacheSizeTmp = static_cast<size_t>(nCacheSizeBig);
897 : }
898 : else
899 : {
900 : // Arbitrarily take half of remaining cache size
901 8 : nCacheSizeTmp = static_cast<size_t>(std::min(
902 16 : static_cast<uint64_t>(
903 8 : (GDALGetCacheMax64() - GDALGetCacheUsed64()) / 2),
904 16 : static_cast<uint64_t>(std::numeric_limits<size_t>::max() / 2)));
905 8 : CPLDebug(ZARR_DEBUG_KEY, "Using implicit CACHE_SIZE=" CPL_FRMT_GUIB,
906 : static_cast<GUIntBig>(nCacheSizeTmp));
907 : }
908 12 : return nCacheSizeTmp;
909 12 : }();
910 12 : if (nCacheSize == std::numeric_limits<size_t>::max())
911 0 : return false;
912 :
913 : // Check that cache size is sufficient to hold all needed tiles.
914 : // Also check that anReqTilesIndices size computation won't overflow.
915 12 : if (nReqTiles > nCacheSize / std::max(m_nTileSize, nDims))
916 : {
917 4 : CPLError(
918 : CE_Failure, CPLE_OutOfMemory,
919 : "CACHE_SIZE=" CPL_FRMT_GUIB " is not big enough to cache "
920 : "all needed tiles. "
921 : "At least " CPL_FRMT_GUIB " bytes would be needed",
922 : static_cast<GUIntBig>(nCacheSize),
923 4 : static_cast<GUIntBig>(nReqTiles * std::max(m_nTileSize, nDims)));
924 4 : return false;
925 : }
926 :
927 8 : const char *pszNumThreads = CSLFetchNameValueDef(
928 : papszOptions, "NUM_THREADS",
929 : CPLGetConfigOption("GDAL_NUM_THREADS", "ALL_CPUS"));
930 8 : if (EQUAL(pszNumThreads, "ALL_CPUS"))
931 8 : nThreadsMax = CPLGetNumCPUs();
932 : else
933 0 : nThreadsMax = std::max(1, atoi(pszNumThreads));
934 8 : if (nThreadsMax > 1024)
935 0 : nThreadsMax = 1024;
936 8 : if (nThreadsMax <= 1)
937 0 : return true;
938 8 : CPLDebug(ZARR_DEBUG_KEY, "IAdviseRead(): Using up to %d threads",
939 : nThreadsMax);
940 :
941 8 : m_oMapTileIndexToCachedTile.clear();
942 :
943 : // Overflow checked above
944 : try
945 : {
946 8 : anReqTilesIndices.resize(nDims * nReqTiles);
947 : }
948 0 : catch (const std::bad_alloc &e)
949 : {
950 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
951 0 : "Cannot allocate anReqTilesIndices: %s", e.what());
952 0 : return false;
953 : }
954 :
955 8 : size_t dimIdx = 0;
956 8 : size_t nTileIter = 0;
957 21608 : lbl_next_depth:
958 21608 : if (dimIdx == nDims)
959 : {
960 21344 : if (nDims == 2)
961 : {
962 : // optimize in common case
963 21344 : memcpy(&anReqTilesIndices[nTileIter * nDims], anIndicesCur.data(),
964 : sizeof(uint64_t) * 2);
965 : }
966 0 : else if (nDims == 3)
967 : {
968 : // optimize in common case
969 0 : memcpy(&anReqTilesIndices[nTileIter * nDims], anIndicesCur.data(),
970 : sizeof(uint64_t) * 3);
971 : }
972 : else
973 : {
974 0 : memcpy(&anReqTilesIndices[nTileIter * nDims], anIndicesCur.data(),
975 0 : sizeof(uint64_t) * nDims);
976 : }
977 21344 : nTileIter++;
978 : }
979 : else
980 : {
981 : // This level of loop loops over blocks
982 264 : anIndicesCur[dimIdx] = anIndicesMin[dimIdx];
983 : while (true)
984 : {
985 21600 : dimIdx++;
986 21600 : goto lbl_next_depth;
987 21600 : lbl_return_to_caller:
988 21600 : dimIdx--;
989 21600 : if (anIndicesCur[dimIdx] == anIndicesMax[dimIdx])
990 264 : break;
991 21336 : ++anIndicesCur[dimIdx];
992 : }
993 : }
994 21608 : if (dimIdx > 0)
995 21600 : goto lbl_return_to_caller;
996 8 : assert(nTileIter == nReqTiles);
997 :
998 8 : return true;
999 : }
1000 :
1001 : /************************************************************************/
1002 : /* ZarrArray::IRead() */
1003 : /************************************************************************/
1004 :
1005 1634 : bool ZarrArray::IRead(const GUInt64 *arrayStartIdx, const size_t *count,
1006 : const GInt64 *arrayStep, const GPtrDiff_t *bufferStride,
1007 : const GDALExtendedDataType &bufferDataType,
1008 : void *pDstBuffer) const
1009 : {
1010 1634 : if (!CheckValidAndErrorOutIfNot())
1011 0 : return false;
1012 :
1013 1634 : if (!AllocateWorkingBuffers())
1014 3 : return false;
1015 :
1016 : // Need to be kept in top-level scope
1017 3262 : std::vector<GUInt64> arrayStartIdxMod;
1018 3262 : std::vector<GInt64> arrayStepMod;
1019 3262 : std::vector<GPtrDiff_t> bufferStrideMod;
1020 :
1021 1631 : const size_t nDims = m_aoDims.size();
1022 1631 : bool negativeStep = false;
1023 4505 : for (size_t i = 0; i < nDims; ++i)
1024 : {
1025 3076 : if (arrayStep[i] < 0)
1026 : {
1027 202 : negativeStep = true;
1028 202 : break;
1029 : }
1030 : }
1031 :
1032 : // const auto eBufferDT = bufferDataType.GetNumericDataType();
1033 1631 : const auto nBufferDTSize = static_cast<int>(bufferDataType.GetSize());
1034 :
1035 : // Make sure that arrayStep[i] are positive for sake of simplicity
1036 1631 : if (negativeStep)
1037 : {
1038 : #if defined(__GNUC__)
1039 : #pragma GCC diagnostic push
1040 : #pragma GCC diagnostic ignored "-Wnull-dereference"
1041 : #endif
1042 202 : arrayStartIdxMod.resize(nDims);
1043 202 : arrayStepMod.resize(nDims);
1044 202 : bufferStrideMod.resize(nDims);
1045 : #if defined(__GNUC__)
1046 : #pragma GCC diagnostic pop
1047 : #endif
1048 606 : for (size_t i = 0; i < nDims; ++i)
1049 : {
1050 404 : if (arrayStep[i] < 0)
1051 : {
1052 808 : arrayStartIdxMod[i] =
1053 404 : arrayStartIdx[i] - (count[i] - 1) * (-arrayStep[i]);
1054 404 : arrayStepMod[i] = -arrayStep[i];
1055 404 : bufferStrideMod[i] = -bufferStride[i];
1056 404 : pDstBuffer =
1057 : static_cast<GByte *>(pDstBuffer) +
1058 404 : bufferStride[i] *
1059 404 : static_cast<GPtrDiff_t>(nBufferDTSize * (count[i] - 1));
1060 : }
1061 : else
1062 : {
1063 0 : arrayStartIdxMod[i] = arrayStartIdx[i];
1064 0 : arrayStepMod[i] = arrayStep[i];
1065 0 : bufferStrideMod[i] = bufferStride[i];
1066 : }
1067 : }
1068 202 : arrayStartIdx = arrayStartIdxMod.data();
1069 202 : arrayStep = arrayStepMod.data();
1070 202 : bufferStride = bufferStrideMod.data();
1071 : }
1072 :
1073 3262 : std::vector<uint64_t> indicesOuterLoop(nDims + 1);
1074 3262 : std::vector<GByte *> dstPtrStackOuterLoop(nDims + 1);
1075 :
1076 3262 : std::vector<uint64_t> indicesInnerLoop(nDims + 1);
1077 3262 : std::vector<GByte *> dstPtrStackInnerLoop(nDims + 1);
1078 :
1079 3262 : std::vector<GPtrDiff_t> dstBufferStrideBytes;
1080 4909 : for (size_t i = 0; i < nDims; ++i)
1081 : {
1082 3278 : dstBufferStrideBytes.push_back(bufferStride[i] *
1083 3278 : static_cast<GPtrDiff_t>(nBufferDTSize));
1084 : }
1085 1631 : dstBufferStrideBytes.push_back(0);
1086 :
1087 1631 : const auto nDTSize = m_oType.GetSize();
1088 :
1089 3262 : std::vector<uint64_t> tileIndices(nDims);
1090 : const size_t nSourceSize =
1091 1631 : m_aoDtypeElts.back().nativeOffset + m_aoDtypeElts.back().nativeSize;
1092 :
1093 3262 : std::vector<size_t> countInnerLoopInit(nDims + 1, 1);
1094 3262 : std::vector<size_t> countInnerLoop(nDims);
1095 :
1096 3245 : const bool bBothAreNumericDT = m_oType.GetClass() == GEDTC_NUMERIC &&
1097 1614 : bufferDataType.GetClass() == GEDTC_NUMERIC;
1098 : const bool bSameNumericDT =
1099 3245 : bBothAreNumericDT &&
1100 1614 : m_oType.GetNumericDataType() == bufferDataType.GetNumericDataType();
1101 1631 : const auto nSameDTSize = bSameNumericDT ? m_oType.GetSize() : 0;
1102 : const bool bSameCompoundAndNoDynamicMem =
1103 1634 : m_oType.GetClass() == GEDTC_COMPOUND && m_oType == bufferDataType &&
1104 3 : !m_oType.NeedsFreeDynamicMemory();
1105 3262 : std::vector<GByte> abyTargetNoData;
1106 1631 : bool bNoDataIsZero = false;
1107 :
1108 1631 : size_t dimIdx = 0;
1109 1631 : dstPtrStackOuterLoop[0] = static_cast<GByte *>(pDstBuffer);
1110 29989 : lbl_next_depth:
1111 29989 : if (dimIdx == nDims)
1112 : {
1113 25360 : size_t dimIdxSubLoop = 0;
1114 25360 : dstPtrStackInnerLoop[0] = dstPtrStackOuterLoop[nDims];
1115 25360 : bool bEmptyTile = false;
1116 :
1117 25360 : const GByte *pabySrcTile = m_abyDecodedTileData.empty()
1118 24369 : ? m_abyRawTileData.data()
1119 25360 : : m_abyDecodedTileData.data();
1120 25360 : bool bMatchFoundInMapTileIndexToCachedTile = false;
1121 :
1122 : // Use cache built by IAdviseRead() if possible
1123 25360 : if (!m_oMapTileIndexToCachedTile.empty())
1124 : {
1125 21352 : uint64_t nTileIdx = 0;
1126 64056 : for (size_t j = 0; j < nDims; ++j)
1127 : {
1128 42704 : if (j > 0)
1129 21352 : nTileIdx *= m_aoDims[j - 1]->GetSize();
1130 42704 : nTileIdx += tileIndices[j];
1131 : }
1132 21352 : const auto oIter = m_oMapTileIndexToCachedTile.find(nTileIdx);
1133 21352 : if (oIter != m_oMapTileIndexToCachedTile.end())
1134 : {
1135 21344 : bMatchFoundInMapTileIndexToCachedTile = true;
1136 21344 : if (oIter->second.abyDecoded.empty())
1137 : {
1138 4 : bEmptyTile = true;
1139 : }
1140 : else
1141 : {
1142 21340 : pabySrcTile = oIter->second.abyDecoded.data();
1143 : }
1144 : }
1145 : else
1146 : {
1147 8 : CPLDebugOnly(ZARR_DEBUG_KEY,
1148 : "Cache miss for tile " CPL_FRMT_GUIB,
1149 : static_cast<GUIntBig>(nTileIdx));
1150 : }
1151 : }
1152 :
1153 25360 : if (!bMatchFoundInMapTileIndexToCachedTile)
1154 : {
1155 4016 : if (!tileIndices.empty() && tileIndices == m_anCachedTiledIndices)
1156 : {
1157 221 : if (!m_bCachedTiledValid)
1158 7 : return false;
1159 221 : bEmptyTile = m_bCachedTiledEmpty;
1160 : }
1161 : else
1162 : {
1163 3795 : if (!FlushDirtyTile())
1164 0 : return false;
1165 :
1166 3795 : m_anCachedTiledIndices = tileIndices;
1167 3795 : m_bCachedTiledValid =
1168 3795 : LoadTileData(tileIndices.data(), bEmptyTile);
1169 3795 : if (!m_bCachedTiledValid)
1170 : {
1171 7 : return false;
1172 : }
1173 3788 : m_bCachedTiledEmpty = bEmptyTile;
1174 : }
1175 :
1176 8018 : pabySrcTile = m_abyDecodedTileData.empty()
1177 3018 : ? m_abyRawTileData.data()
1178 991 : : m_abyDecodedTileData.data();
1179 : }
1180 : const size_t nSrcDTSize =
1181 25353 : m_abyDecodedTileData.empty() ? nSourceSize : nDTSize;
1182 :
1183 76182 : for (size_t i = 0; i < nDims; ++i)
1184 : {
1185 50829 : countInnerLoopInit[i] = 1;
1186 50829 : if (arrayStep[i] != 0)
1187 : {
1188 : const auto nextBlockIdx =
1189 50127 : std::min((1 + indicesOuterLoop[i] / m_anBlockSize[i]) *
1190 50127 : m_anBlockSize[i],
1191 100254 : arrayStartIdx[i] + count[i] * arrayStep[i]);
1192 50127 : countInnerLoopInit[i] = static_cast<size_t>(
1193 50127 : (nextBlockIdx - indicesOuterLoop[i] + arrayStep[i] - 1) /
1194 50127 : arrayStep[i]);
1195 : }
1196 : }
1197 :
1198 25353 : if (bEmptyTile && bBothAreNumericDT && abyTargetNoData.empty())
1199 : {
1200 603 : abyTargetNoData.resize(nBufferDTSize);
1201 603 : if (m_pabyNoData)
1202 : {
1203 156 : GDALExtendedDataType::CopyValue(
1204 156 : m_pabyNoData, m_oType, &abyTargetNoData[0], bufferDataType);
1205 156 : bNoDataIsZero = true;
1206 1308 : for (size_t i = 0; i < abyTargetNoData.size(); ++i)
1207 : {
1208 1152 : if (abyTargetNoData[i] != 0)
1209 334 : bNoDataIsZero = false;
1210 : }
1211 : }
1212 : else
1213 : {
1214 447 : bNoDataIsZero = true;
1215 447 : GByte zero = 0;
1216 447 : GDALCopyWords(&zero, GDT_Byte, 0, &abyTargetNoData[0],
1217 : bufferDataType.GetNumericDataType(), 0, 1);
1218 : }
1219 : }
1220 :
1221 24750 : lbl_next_depth_inner_loop:
1222 467811 : if (nDims == 0 || dimIdxSubLoop == nDims - 1)
1223 : {
1224 442331 : indicesInnerLoop[dimIdxSubLoop] = indicesOuterLoop[dimIdxSubLoop];
1225 442331 : void *dst_ptr = dstPtrStackInnerLoop[dimIdxSubLoop];
1226 :
1227 439844 : if (m_bUseOptimizedCodePaths && bEmptyTile && bBothAreNumericDT &&
1228 882175 : bNoDataIsZero &&
1229 1135 : nBufferDTSize == dstBufferStrideBytes[dimIdxSubLoop])
1230 : {
1231 1069 : memset(dst_ptr, 0,
1232 1069 : nBufferDTSize * countInnerLoopInit[dimIdxSubLoop]);
1233 1069 : goto end_inner_loop;
1234 : }
1235 438775 : else if (m_bUseOptimizedCodePaths && bEmptyTile &&
1236 880412 : !abyTargetNoData.empty() && bBothAreNumericDT &&
1237 375 : dstBufferStrideBytes[dimIdxSubLoop] <
1238 375 : std::numeric_limits<int>::max())
1239 : {
1240 750 : GDALCopyWords64(
1241 375 : abyTargetNoData.data(), bufferDataType.GetNumericDataType(),
1242 : 0, dst_ptr, bufferDataType.GetNumericDataType(),
1243 375 : static_cast<int>(dstBufferStrideBytes[dimIdxSubLoop]),
1244 375 : static_cast<GPtrDiff_t>(countInnerLoopInit[dimIdxSubLoop]));
1245 375 : goto end_inner_loop;
1246 : }
1247 440887 : else if (bEmptyTile)
1248 : {
1249 3875 : for (size_t i = 0; i < countInnerLoopInit[dimIdxSubLoop];
1250 2528 : ++i, dst_ptr = static_cast<uint8_t *>(dst_ptr) +
1251 2528 : dstBufferStrideBytes[dimIdxSubLoop])
1252 : {
1253 2528 : if (bNoDataIsZero)
1254 : {
1255 1930 : if (nBufferDTSize == 1)
1256 : {
1257 36 : *static_cast<uint8_t *>(dst_ptr) = 0;
1258 : }
1259 1894 : else if (nBufferDTSize == 2)
1260 : {
1261 48 : *static_cast<uint16_t *>(dst_ptr) = 0;
1262 : }
1263 1846 : else if (nBufferDTSize == 4)
1264 : {
1265 72 : *static_cast<uint32_t *>(dst_ptr) = 0;
1266 : }
1267 1774 : else if (nBufferDTSize == 8)
1268 : {
1269 1534 : *static_cast<uint64_t *>(dst_ptr) = 0;
1270 : }
1271 240 : else if (nBufferDTSize == 16)
1272 : {
1273 240 : static_cast<uint64_t *>(dst_ptr)[0] = 0;
1274 240 : static_cast<uint64_t *>(dst_ptr)[1] = 0;
1275 : }
1276 : else
1277 : {
1278 0 : CPLAssert(false);
1279 : }
1280 : }
1281 598 : else if (m_pabyNoData)
1282 : {
1283 597 : if (bBothAreNumericDT)
1284 : {
1285 588 : const void *src_ptr_v = abyTargetNoData.data();
1286 588 : if (nBufferDTSize == 1)
1287 24 : *static_cast<uint8_t *>(dst_ptr) =
1288 24 : *static_cast<const uint8_t *>(src_ptr_v);
1289 564 : else if (nBufferDTSize == 2)
1290 0 : *static_cast<uint16_t *>(dst_ptr) =
1291 0 : *static_cast<const uint16_t *>(src_ptr_v);
1292 564 : else if (nBufferDTSize == 4)
1293 60 : *static_cast<uint32_t *>(dst_ptr) =
1294 60 : *static_cast<const uint32_t *>(src_ptr_v);
1295 504 : else if (nBufferDTSize == 8)
1296 504 : *static_cast<uint64_t *>(dst_ptr) =
1297 504 : *static_cast<const uint64_t *>(src_ptr_v);
1298 0 : else if (nBufferDTSize == 16)
1299 : {
1300 0 : static_cast<uint64_t *>(dst_ptr)[0] =
1301 0 : static_cast<const uint64_t *>(src_ptr_v)[0];
1302 0 : static_cast<uint64_t *>(dst_ptr)[1] =
1303 0 : static_cast<const uint64_t *>(src_ptr_v)[1];
1304 : }
1305 : else
1306 : {
1307 0 : CPLAssert(false);
1308 : }
1309 : }
1310 : else
1311 : {
1312 9 : GDALExtendedDataType::CopyValue(
1313 9 : m_pabyNoData, m_oType, dst_ptr, bufferDataType);
1314 : }
1315 : }
1316 : else
1317 : {
1318 1 : memset(dst_ptr, 0, nBufferDTSize);
1319 : }
1320 : }
1321 :
1322 1347 : goto end_inner_loop;
1323 : }
1324 :
1325 439540 : size_t nOffset = 0;
1326 1331060 : for (size_t i = 0; i < nDims; i++)
1327 : {
1328 891516 : nOffset = static_cast<size_t>(
1329 891516 : nOffset * m_anBlockSize[i] +
1330 1783030 : (indicesInnerLoop[i] - tileIndices[i] * m_anBlockSize[i]));
1331 : }
1332 439540 : const GByte *src_ptr = pabySrcTile + nOffset * nSrcDTSize;
1333 439540 : const auto step = nDims == 0 ? 0 : arrayStep[dimIdxSubLoop];
1334 :
1335 438390 : if (m_bUseOptimizedCodePaths && bBothAreNumericDT &&
1336 438366 : step <= static_cast<GIntBig>(std::numeric_limits<int>::max() /
1337 877930 : nDTSize) &&
1338 438366 : dstBufferStrideBytes[dimIdxSubLoop] <=
1339 438366 : std::numeric_limits<int>::max())
1340 : {
1341 438366 : GDALCopyWords64(
1342 : src_ptr, m_oType.GetNumericDataType(),
1343 : static_cast<int>(step * nDTSize), dst_ptr,
1344 : bufferDataType.GetNumericDataType(),
1345 438366 : static_cast<int>(dstBufferStrideBytes[dimIdxSubLoop]),
1346 438366 : static_cast<GPtrDiff_t>(countInnerLoopInit[dimIdxSubLoop]));
1347 :
1348 438366 : goto end_inner_loop;
1349 : }
1350 :
1351 3391 : for (size_t i = 0; i < countInnerLoopInit[dimIdxSubLoop];
1352 2217 : ++i, src_ptr += step * nSrcDTSize,
1353 2217 : dst_ptr = static_cast<uint8_t *>(dst_ptr) +
1354 2217 : dstBufferStrideBytes[dimIdxSubLoop])
1355 : {
1356 2217 : if (bSameNumericDT)
1357 : {
1358 814 : const void *src_ptr_v = src_ptr;
1359 814 : if (nSameDTSize == 1)
1360 70 : *static_cast<uint8_t *>(dst_ptr) =
1361 70 : *static_cast<const uint8_t *>(src_ptr_v);
1362 744 : else if (nSameDTSize == 2)
1363 : {
1364 56 : *static_cast<uint16_t *>(dst_ptr) =
1365 56 : *static_cast<const uint16_t *>(src_ptr_v);
1366 : }
1367 688 : else if (nSameDTSize == 4)
1368 : {
1369 154 : *static_cast<uint32_t *>(dst_ptr) =
1370 154 : *static_cast<const uint32_t *>(src_ptr_v);
1371 : }
1372 534 : else if (nSameDTSize == 8)
1373 : {
1374 482 : *static_cast<uint64_t *>(dst_ptr) =
1375 482 : *static_cast<const uint64_t *>(src_ptr_v);
1376 : }
1377 52 : else if (nSameDTSize == 16)
1378 : {
1379 52 : static_cast<uint64_t *>(dst_ptr)[0] =
1380 52 : static_cast<const uint64_t *>(src_ptr_v)[0];
1381 52 : static_cast<uint64_t *>(dst_ptr)[1] =
1382 52 : static_cast<const uint64_t *>(src_ptr_v)[1];
1383 : }
1384 : else
1385 : {
1386 0 : CPLAssert(false);
1387 : }
1388 : }
1389 1403 : else if (bSameCompoundAndNoDynamicMem)
1390 : {
1391 4 : memcpy(dst_ptr, src_ptr, nDTSize);
1392 : }
1393 1399 : else if (m_oType.GetClass() == GEDTC_STRING)
1394 : {
1395 26 : if (m_aoDtypeElts.back().nativeType ==
1396 : DtypeElt::NativeType::STRING_UNICODE)
1397 : {
1398 : char *pDstStr =
1399 12 : UCS4ToUTF8(src_ptr, nSourceSize,
1400 6 : m_aoDtypeElts.back().needByteSwapping);
1401 6 : char **pDstPtr = static_cast<char **>(dst_ptr);
1402 6 : memcpy(pDstPtr, &pDstStr, sizeof(pDstStr));
1403 : }
1404 : else
1405 : {
1406 : char *pDstStr =
1407 20 : static_cast<char *>(CPLMalloc(nSourceSize + 1));
1408 20 : memcpy(pDstStr, src_ptr, nSourceSize);
1409 20 : pDstStr[nSourceSize] = 0;
1410 20 : char **pDstPtr = static_cast<char **>(dst_ptr);
1411 20 : memcpy(pDstPtr, &pDstStr, sizeof(char *));
1412 : }
1413 : }
1414 : else
1415 : {
1416 1373 : GDALExtendedDataType::CopyValue(src_ptr, m_oType, dst_ptr,
1417 : bufferDataType);
1418 : }
1419 1174 : }
1420 : }
1421 : else
1422 : {
1423 : // This level of loop loops over individual samples, within a
1424 : // block
1425 25480 : indicesInnerLoop[dimIdxSubLoop] = indicesOuterLoop[dimIdxSubLoop];
1426 25480 : countInnerLoop[dimIdxSubLoop] = countInnerLoopInit[dimIdxSubLoop];
1427 : while (true)
1428 : {
1429 442458 : dimIdxSubLoop++;
1430 442458 : dstPtrStackInnerLoop[dimIdxSubLoop] =
1431 442458 : dstPtrStackInnerLoop[dimIdxSubLoop - 1];
1432 442458 : goto lbl_next_depth_inner_loop;
1433 442458 : lbl_return_to_caller_inner_loop:
1434 442458 : dimIdxSubLoop--;
1435 442458 : --countInnerLoop[dimIdxSubLoop];
1436 442458 : if (countInnerLoop[dimIdxSubLoop] == 0)
1437 : {
1438 25480 : break;
1439 : }
1440 416978 : indicesInnerLoop[dimIdxSubLoop] += arrayStep[dimIdxSubLoop];
1441 416978 : dstPtrStackInnerLoop[dimIdxSubLoop] +=
1442 416978 : dstBufferStrideBytes[dimIdxSubLoop];
1443 : }
1444 : }
1445 467811 : end_inner_loop:
1446 467811 : if (dimIdxSubLoop > 0)
1447 442458 : goto lbl_return_to_caller_inner_loop;
1448 : }
1449 : else
1450 : {
1451 : // This level of loop loops over blocks
1452 4629 : indicesOuterLoop[dimIdx] = arrayStartIdx[dimIdx];
1453 4629 : tileIndices[dimIdx] = indicesOuterLoop[dimIdx] / m_anBlockSize[dimIdx];
1454 : while (true)
1455 : {
1456 28358 : dimIdx++;
1457 28358 : dstPtrStackOuterLoop[dimIdx] = dstPtrStackOuterLoop[dimIdx - 1];
1458 28358 : goto lbl_next_depth;
1459 28348 : lbl_return_to_caller:
1460 28348 : dimIdx--;
1461 28348 : if (count[dimIdx] == 1 || arrayStep[dimIdx] == 0)
1462 : break;
1463 :
1464 : size_t nIncr;
1465 27255 : if (static_cast<GUInt64>(arrayStep[dimIdx]) < m_anBlockSize[dimIdx])
1466 : {
1467 : // Compute index at next block boundary
1468 : auto newIdx =
1469 26847 : indicesOuterLoop[dimIdx] +
1470 26847 : (m_anBlockSize[dimIdx] -
1471 26847 : (indicesOuterLoop[dimIdx] % m_anBlockSize[dimIdx]));
1472 : // And round up compared to arrayStartIdx, arrayStep
1473 26847 : nIncr = static_cast<size_t>((newIdx - indicesOuterLoop[dimIdx] +
1474 26847 : arrayStep[dimIdx] - 1) /
1475 26847 : arrayStep[dimIdx]);
1476 : }
1477 : else
1478 : {
1479 408 : nIncr = 1;
1480 : }
1481 27255 : indicesOuterLoop[dimIdx] += nIncr * arrayStep[dimIdx];
1482 27255 : if (indicesOuterLoop[dimIdx] >
1483 27255 : arrayStartIdx[dimIdx] + (count[dimIdx] - 1) * arrayStep[dimIdx])
1484 3526 : break;
1485 23729 : dstPtrStackOuterLoop[dimIdx] +=
1486 23729 : bufferStride[dimIdx] *
1487 23729 : static_cast<GPtrDiff_t>(nIncr * nBufferDTSize);
1488 47458 : tileIndices[dimIdx] =
1489 23729 : indicesOuterLoop[dimIdx] / m_anBlockSize[dimIdx];
1490 23729 : }
1491 : }
1492 29972 : if (dimIdx > 0)
1493 28348 : goto lbl_return_to_caller;
1494 :
1495 1624 : return true;
1496 : }
1497 :
1498 : /************************************************************************/
1499 : /* ZarrArray::IWrite() */
1500 : /************************************************************************/
1501 :
1502 494 : bool ZarrArray::IWrite(const GUInt64 *arrayStartIdx, const size_t *count,
1503 : const GInt64 *arrayStep, const GPtrDiff_t *bufferStride,
1504 : const GDALExtendedDataType &bufferDataType,
1505 : const void *pSrcBuffer)
1506 : {
1507 494 : if (!CheckValidAndErrorOutIfNot())
1508 0 : return false;
1509 :
1510 494 : if (!AllocateWorkingBuffers())
1511 0 : return false;
1512 :
1513 494 : m_oMapTileIndexToCachedTile.clear();
1514 :
1515 : // Need to be kept in top-level scope
1516 988 : std::vector<GUInt64> arrayStartIdxMod;
1517 988 : std::vector<GInt64> arrayStepMod;
1518 988 : std::vector<GPtrDiff_t> bufferStrideMod;
1519 :
1520 494 : const size_t nDims = m_aoDims.size();
1521 494 : bool negativeStep = false;
1522 494 : bool bWriteWholeTileInit = true;
1523 1384 : for (size_t i = 0; i < nDims; ++i)
1524 : {
1525 890 : if (arrayStep[i] < 0)
1526 : {
1527 132 : negativeStep = true;
1528 132 : if (arrayStep[i] != -1 && count[i] > 1)
1529 0 : bWriteWholeTileInit = false;
1530 : }
1531 758 : else if (arrayStep[i] != 1 && count[i] > 1)
1532 0 : bWriteWholeTileInit = false;
1533 : }
1534 :
1535 494 : const auto nBufferDTSize = static_cast<int>(bufferDataType.GetSize());
1536 :
1537 : // Make sure that arrayStep[i] are positive for sake of simplicity
1538 494 : if (negativeStep)
1539 : {
1540 : #if defined(__GNUC__)
1541 : #pragma GCC diagnostic push
1542 : #pragma GCC diagnostic ignored "-Wnull-dereference"
1543 : #endif
1544 66 : arrayStartIdxMod.resize(nDims);
1545 66 : arrayStepMod.resize(nDims);
1546 66 : bufferStrideMod.resize(nDims);
1547 : #if defined(__GNUC__)
1548 : #pragma GCC diagnostic pop
1549 : #endif
1550 198 : for (size_t i = 0; i < nDims; ++i)
1551 : {
1552 132 : if (arrayStep[i] < 0)
1553 : {
1554 264 : arrayStartIdxMod[i] =
1555 132 : arrayStartIdx[i] - (count[i] - 1) * (-arrayStep[i]);
1556 132 : arrayStepMod[i] = -arrayStep[i];
1557 132 : bufferStrideMod[i] = -bufferStride[i];
1558 132 : pSrcBuffer =
1559 : static_cast<const GByte *>(pSrcBuffer) +
1560 132 : bufferStride[i] *
1561 132 : static_cast<GPtrDiff_t>(nBufferDTSize * (count[i] - 1));
1562 : }
1563 : else
1564 : {
1565 0 : arrayStartIdxMod[i] = arrayStartIdx[i];
1566 0 : arrayStepMod[i] = arrayStep[i];
1567 0 : bufferStrideMod[i] = bufferStride[i];
1568 : }
1569 : }
1570 66 : arrayStartIdx = arrayStartIdxMod.data();
1571 66 : arrayStep = arrayStepMod.data();
1572 66 : bufferStride = bufferStrideMod.data();
1573 : }
1574 :
1575 988 : std::vector<uint64_t> indicesOuterLoop(nDims + 1);
1576 988 : std::vector<const GByte *> srcPtrStackOuterLoop(nDims + 1);
1577 :
1578 988 : std::vector<size_t> offsetDstBuffer(nDims + 1);
1579 988 : std::vector<const GByte *> srcPtrStackInnerLoop(nDims + 1);
1580 :
1581 988 : std::vector<GPtrDiff_t> srcBufferStrideBytes;
1582 1384 : for (size_t i = 0; i < nDims; ++i)
1583 : {
1584 890 : srcBufferStrideBytes.push_back(bufferStride[i] *
1585 890 : static_cast<GPtrDiff_t>(nBufferDTSize));
1586 : }
1587 494 : srcBufferStrideBytes.push_back(0);
1588 :
1589 494 : const auto nDTSize = m_oType.GetSize();
1590 :
1591 988 : std::vector<uint64_t> tileIndices(nDims);
1592 : const size_t nNativeSize =
1593 494 : m_aoDtypeElts.back().nativeOffset + m_aoDtypeElts.back().nativeSize;
1594 :
1595 988 : std::vector<size_t> countInnerLoopInit(nDims + 1, 1);
1596 988 : std::vector<size_t> countInnerLoop(nDims);
1597 :
1598 984 : const bool bBothAreNumericDT = m_oType.GetClass() == GEDTC_NUMERIC &&
1599 490 : bufferDataType.GetClass() == GEDTC_NUMERIC;
1600 : const bool bSameNumericDT =
1601 984 : bBothAreNumericDT &&
1602 490 : m_oType.GetNumericDataType() == bufferDataType.GetNumericDataType();
1603 494 : const auto nSameDTSize = bSameNumericDT ? m_oType.GetSize() : 0;
1604 : const bool bSameCompoundAndNoDynamicMem =
1605 494 : m_oType.GetClass() == GEDTC_COMPOUND && m_oType == bufferDataType &&
1606 0 : !m_oType.NeedsFreeDynamicMemory();
1607 :
1608 494 : size_t dimIdx = 0;
1609 494 : size_t dimIdxForCopy = nDims == 0 ? 0 : nDims - 1;
1610 494 : if (nDims)
1611 : {
1612 562 : while (dimIdxForCopy > 0 && count[dimIdxForCopy] == 1)
1613 68 : --dimIdxForCopy;
1614 : }
1615 :
1616 494 : srcPtrStackOuterLoop[0] = static_cast<const GByte *>(pSrcBuffer);
1617 24644 : lbl_next_depth:
1618 24644 : if (dimIdx == nDims)
1619 : {
1620 22820 : bool bWriteWholeTile = bWriteWholeTileInit;
1621 22820 : bool bPartialTile = false;
1622 68623 : for (size_t i = 0; i < nDims; ++i)
1623 : {
1624 45803 : countInnerLoopInit[i] = 1;
1625 45803 : if (arrayStep[i] != 0)
1626 : {
1627 : const auto nextBlockIdx =
1628 45520 : std::min((1 + indicesOuterLoop[i] / m_anBlockSize[i]) *
1629 45520 : m_anBlockSize[i],
1630 91040 : arrayStartIdx[i] + count[i] * arrayStep[i]);
1631 45520 : countInnerLoopInit[i] = static_cast<size_t>(
1632 45520 : (nextBlockIdx - indicesOuterLoop[i] + arrayStep[i] - 1) /
1633 45520 : arrayStep[i]);
1634 : }
1635 45803 : if (bWriteWholeTile)
1636 : {
1637 : const bool bWholePartialTileThisDim =
1638 46689 : indicesOuterLoop[i] == 0 &&
1639 1834 : countInnerLoopInit[i] == m_aoDims[i]->GetSize();
1640 44855 : bWriteWholeTile = (countInnerLoopInit[i] == m_anBlockSize[i] ||
1641 : bWholePartialTileThisDim);
1642 44855 : if (bWholePartialTileThisDim)
1643 : {
1644 587 : bPartialTile = true;
1645 : }
1646 : }
1647 : }
1648 :
1649 22820 : size_t dimIdxSubLoop = 0;
1650 22820 : srcPtrStackInnerLoop[0] = srcPtrStackOuterLoop[nDims];
1651 : const size_t nCacheDTSize =
1652 22820 : m_abyDecodedTileData.empty() ? nNativeSize : nDTSize;
1653 22820 : auto &abyTile = m_abyDecodedTileData.empty() ? m_abyRawTileData
1654 22820 : : m_abyDecodedTileData;
1655 :
1656 22820 : if (!tileIndices.empty() && tileIndices == m_anCachedTiledIndices)
1657 : {
1658 4 : if (!m_bCachedTiledValid)
1659 0 : return false;
1660 : }
1661 : else
1662 : {
1663 22816 : if (!FlushDirtyTile())
1664 0 : return false;
1665 :
1666 22816 : m_anCachedTiledIndices = tileIndices;
1667 22816 : m_bCachedTiledValid = true;
1668 :
1669 22816 : if (bWriteWholeTile)
1670 : {
1671 21292 : if (bPartialTile)
1672 : {
1673 288 : DeallocateDecodedTileData();
1674 288 : memset(&abyTile[0], 0, abyTile.size());
1675 : }
1676 : }
1677 : else
1678 : {
1679 : // If we don't write the whole tile, we need to fetch a
1680 : // potentially existing one.
1681 1524 : bool bEmptyTile = false;
1682 1524 : m_bCachedTiledValid =
1683 1524 : LoadTileData(tileIndices.data(), bEmptyTile);
1684 1524 : if (!m_bCachedTiledValid)
1685 : {
1686 0 : return false;
1687 : }
1688 :
1689 1524 : if (bEmptyTile)
1690 : {
1691 1215 : DeallocateDecodedTileData();
1692 :
1693 1215 : if (m_pabyNoData == nullptr)
1694 : {
1695 495 : memset(&abyTile[0], 0, abyTile.size());
1696 : }
1697 : else
1698 : {
1699 720 : const size_t nElts = abyTile.size() / nCacheDTSize;
1700 720 : GByte *dstPtr = &abyTile[0];
1701 720 : if (m_oType.GetClass() == GEDTC_NUMERIC)
1702 : {
1703 720 : GDALCopyWords64(
1704 720 : m_pabyNoData, m_oType.GetNumericDataType(), 0,
1705 : dstPtr, m_oType.GetNumericDataType(),
1706 720 : static_cast<int>(m_oType.GetSize()),
1707 : static_cast<GPtrDiff_t>(nElts));
1708 : }
1709 : else
1710 : {
1711 0 : for (size_t i = 0; i < nElts; ++i)
1712 : {
1713 0 : GDALExtendedDataType::CopyValue(
1714 0 : m_pabyNoData, m_oType, dstPtr, m_oType);
1715 0 : dstPtr += nCacheDTSize;
1716 : }
1717 : }
1718 : }
1719 : }
1720 : }
1721 : }
1722 22820 : m_bDirtyTile = true;
1723 22820 : m_bCachedTiledEmpty = false;
1724 22820 : if (nDims)
1725 22820 : offsetDstBuffer[0] = static_cast<size_t>(
1726 22820 : indicesOuterLoop[0] - tileIndices[0] * m_anBlockSize[0]);
1727 :
1728 22820 : GByte *pabyTile = &abyTile[0];
1729 :
1730 453424 : lbl_next_depth_inner_loop:
1731 453424 : if (dimIdxSubLoop == dimIdxForCopy)
1732 : {
1733 430666 : size_t nOffset = offsetDstBuffer[dimIdxSubLoop];
1734 430666 : GInt64 step = nDims == 0 ? 0 : arrayStep[dimIdxSubLoop];
1735 430891 : for (size_t i = dimIdxSubLoop + 1; i < nDims; ++i)
1736 : {
1737 225 : nOffset = static_cast<size_t>(
1738 225 : nOffset * m_anBlockSize[i] +
1739 450 : (indicesOuterLoop[i] - tileIndices[i] * m_anBlockSize[i]));
1740 225 : step *= m_anBlockSize[i];
1741 : }
1742 430666 : const void *src_ptr = srcPtrStackInnerLoop[dimIdxSubLoop];
1743 430666 : GByte *dst_ptr = pabyTile + nOffset * nCacheDTSize;
1744 :
1745 430666 : if (m_bUseOptimizedCodePaths && bBothAreNumericDT)
1746 : {
1747 429936 : if (countInnerLoopInit[dimIdxSubLoop] == 1 && bSameNumericDT)
1748 : {
1749 86 : void *dst_ptr_v = dst_ptr;
1750 86 : if (nSameDTSize == 1)
1751 8 : *static_cast<uint8_t *>(dst_ptr_v) =
1752 8 : *static_cast<const uint8_t *>(src_ptr);
1753 78 : else if (nSameDTSize == 2)
1754 : {
1755 2 : *static_cast<uint16_t *>(dst_ptr_v) =
1756 2 : *static_cast<const uint16_t *>(src_ptr);
1757 : }
1758 76 : else if (nSameDTSize == 4)
1759 : {
1760 2 : *static_cast<uint32_t *>(dst_ptr_v) =
1761 2 : *static_cast<const uint32_t *>(src_ptr);
1762 : }
1763 74 : else if (nSameDTSize == 8)
1764 : {
1765 52 : *static_cast<uint64_t *>(dst_ptr_v) =
1766 52 : *static_cast<const uint64_t *>(src_ptr);
1767 : }
1768 22 : else if (nSameDTSize == 16)
1769 : {
1770 22 : static_cast<uint64_t *>(dst_ptr_v)[0] =
1771 22 : static_cast<const uint64_t *>(src_ptr)[0];
1772 22 : static_cast<uint64_t *>(dst_ptr_v)[1] =
1773 22 : static_cast<const uint64_t *>(src_ptr)[1];
1774 : }
1775 : else
1776 : {
1777 0 : CPLAssert(false);
1778 : }
1779 : }
1780 429850 : else if (step <=
1781 429850 : static_cast<GIntBig>(
1782 859700 : std::numeric_limits<int>::max() / nDTSize) &&
1783 429850 : srcBufferStrideBytes[dimIdxSubLoop] <=
1784 429850 : std::numeric_limits<int>::max())
1785 : {
1786 859700 : GDALCopyWords64(
1787 : src_ptr, bufferDataType.GetNumericDataType(),
1788 429850 : static_cast<int>(srcBufferStrideBytes[dimIdxSubLoop]),
1789 : dst_ptr, m_oType.GetNumericDataType(),
1790 : static_cast<int>(step * nDTSize),
1791 : static_cast<GPtrDiff_t>(
1792 429850 : countInnerLoopInit[dimIdxSubLoop]));
1793 : }
1794 429936 : goto end_inner_loop;
1795 : }
1796 :
1797 2188 : for (size_t i = 0; i < countInnerLoopInit[dimIdxSubLoop];
1798 1458 : ++i, dst_ptr += step * nCacheDTSize,
1799 1458 : src_ptr = static_cast<const uint8_t *>(src_ptr) +
1800 1458 : srcBufferStrideBytes[dimIdxSubLoop])
1801 : {
1802 1458 : if (bSameNumericDT)
1803 : {
1804 300 : void *dst_ptr_v = dst_ptr;
1805 300 : if (nSameDTSize == 1)
1806 0 : *static_cast<uint8_t *>(dst_ptr_v) =
1807 0 : *static_cast<const uint8_t *>(src_ptr);
1808 300 : else if (nSameDTSize == 2)
1809 : {
1810 0 : *static_cast<uint16_t *>(dst_ptr_v) =
1811 0 : *static_cast<const uint16_t *>(src_ptr);
1812 : }
1813 300 : else if (nSameDTSize == 4)
1814 : {
1815 0 : *static_cast<uint32_t *>(dst_ptr_v) =
1816 0 : *static_cast<const uint32_t *>(src_ptr);
1817 : }
1818 300 : else if (nSameDTSize == 8)
1819 : {
1820 220 : *static_cast<uint64_t *>(dst_ptr_v) =
1821 220 : *static_cast<const uint64_t *>(src_ptr);
1822 : }
1823 80 : else if (nSameDTSize == 16)
1824 : {
1825 80 : static_cast<uint64_t *>(dst_ptr_v)[0] =
1826 80 : static_cast<const uint64_t *>(src_ptr)[0];
1827 80 : static_cast<uint64_t *>(dst_ptr_v)[1] =
1828 80 : static_cast<const uint64_t *>(src_ptr)[1];
1829 : }
1830 : else
1831 : {
1832 0 : CPLAssert(false);
1833 : }
1834 : }
1835 1158 : else if (bSameCompoundAndNoDynamicMem)
1836 : {
1837 0 : memcpy(dst_ptr, src_ptr, nDTSize);
1838 : }
1839 1158 : else if (m_oType.GetClass() == GEDTC_STRING)
1840 : {
1841 6 : const char *pSrcStr =
1842 : *static_cast<const char *const *>(src_ptr);
1843 6 : if (pSrcStr)
1844 : {
1845 6 : const size_t nLen = strlen(pSrcStr);
1846 6 : if (m_aoDtypeElts.back().nativeType ==
1847 : DtypeElt::NativeType::STRING_UNICODE)
1848 : {
1849 : try
1850 : {
1851 : const auto ucs4 = UTF8ToUCS4(
1852 : pSrcStr,
1853 8 : m_aoDtypeElts.back().needByteSwapping);
1854 4 : const auto ucs4Len = ucs4.size();
1855 4 : memcpy(dst_ptr, ucs4.data(),
1856 4 : std::min(ucs4Len, nNativeSize));
1857 4 : if (ucs4Len > nNativeSize)
1858 : {
1859 1 : CPLError(CE_Warning, CPLE_AppDefined,
1860 : "Too long string truncated");
1861 : }
1862 3 : else if (ucs4Len < nNativeSize)
1863 : {
1864 1 : memset(dst_ptr + ucs4Len, 0,
1865 1 : nNativeSize - ucs4Len);
1866 : }
1867 : }
1868 0 : catch (const std::exception &)
1869 : {
1870 0 : memset(dst_ptr, 0, nNativeSize);
1871 : }
1872 : }
1873 : else
1874 : {
1875 2 : memcpy(dst_ptr, pSrcStr,
1876 2 : std::min(nLen, nNativeSize));
1877 2 : if (nLen < nNativeSize)
1878 1 : memset(dst_ptr + nLen, 0, nNativeSize - nLen);
1879 : }
1880 : }
1881 : else
1882 : {
1883 0 : memset(dst_ptr, 0, nNativeSize);
1884 : }
1885 : }
1886 : else
1887 : {
1888 1152 : if (m_oType.NeedsFreeDynamicMemory())
1889 0 : m_oType.FreeDynamicMemory(dst_ptr);
1890 1152 : GDALExtendedDataType::CopyValue(src_ptr, bufferDataType,
1891 1152 : dst_ptr, m_oType);
1892 : }
1893 : }
1894 : }
1895 : else
1896 : {
1897 : // This level of loop loops over individual samples, within a
1898 : // block
1899 22758 : countInnerLoop[dimIdxSubLoop] = countInnerLoopInit[dimIdxSubLoop];
1900 : while (true)
1901 : {
1902 430604 : dimIdxSubLoop++;
1903 430604 : srcPtrStackInnerLoop[dimIdxSubLoop] =
1904 430604 : srcPtrStackInnerLoop[dimIdxSubLoop - 1];
1905 861208 : offsetDstBuffer[dimIdxSubLoop] =
1906 1291810 : static_cast<size_t>(offsetDstBuffer[dimIdxSubLoop - 1] *
1907 430604 : m_anBlockSize[dimIdxSubLoop] +
1908 430604 : (indicesOuterLoop[dimIdxSubLoop] -
1909 861208 : tileIndices[dimIdxSubLoop] *
1910 430604 : m_anBlockSize[dimIdxSubLoop]));
1911 430604 : goto lbl_next_depth_inner_loop;
1912 430604 : lbl_return_to_caller_inner_loop:
1913 430604 : dimIdxSubLoop--;
1914 430604 : --countInnerLoop[dimIdxSubLoop];
1915 430604 : if (countInnerLoop[dimIdxSubLoop] == 0)
1916 : {
1917 22758 : break;
1918 : }
1919 407846 : srcPtrStackInnerLoop[dimIdxSubLoop] +=
1920 407846 : srcBufferStrideBytes[dimIdxSubLoop];
1921 407846 : offsetDstBuffer[dimIdxSubLoop] +=
1922 407846 : static_cast<size_t>(arrayStep[dimIdxSubLoop]);
1923 : }
1924 : }
1925 430666 : end_inner_loop:
1926 453424 : if (dimIdxSubLoop > 0)
1927 430604 : goto lbl_return_to_caller_inner_loop;
1928 : }
1929 : else
1930 : {
1931 : // This level of loop loops over blocks
1932 1824 : indicesOuterLoop[dimIdx] = arrayStartIdx[dimIdx];
1933 1824 : tileIndices[dimIdx] = indicesOuterLoop[dimIdx] / m_anBlockSize[dimIdx];
1934 : while (true)
1935 : {
1936 24150 : dimIdx++;
1937 24150 : srcPtrStackOuterLoop[dimIdx] = srcPtrStackOuterLoop[dimIdx - 1];
1938 24150 : goto lbl_next_depth;
1939 24150 : lbl_return_to_caller:
1940 24150 : dimIdx--;
1941 24150 : if (count[dimIdx] == 1 || arrayStep[dimIdx] == 0)
1942 : break;
1943 :
1944 : size_t nIncr;
1945 23897 : if (static_cast<GUInt64>(arrayStep[dimIdx]) < m_anBlockSize[dimIdx])
1946 : {
1947 : // Compute index at next block boundary
1948 : auto newIdx =
1949 23705 : indicesOuterLoop[dimIdx] +
1950 23705 : (m_anBlockSize[dimIdx] -
1951 23705 : (indicesOuterLoop[dimIdx] % m_anBlockSize[dimIdx]));
1952 : // And round up compared to arrayStartIdx, arrayStep
1953 23705 : nIncr = static_cast<size_t>((newIdx - indicesOuterLoop[dimIdx] +
1954 23705 : arrayStep[dimIdx] - 1) /
1955 23705 : arrayStep[dimIdx]);
1956 : }
1957 : else
1958 : {
1959 192 : nIncr = 1;
1960 : }
1961 23897 : indicesOuterLoop[dimIdx] += nIncr * arrayStep[dimIdx];
1962 23897 : if (indicesOuterLoop[dimIdx] >
1963 23897 : arrayStartIdx[dimIdx] + (count[dimIdx] - 1) * arrayStep[dimIdx])
1964 1571 : break;
1965 22326 : srcPtrStackOuterLoop[dimIdx] +=
1966 22326 : bufferStride[dimIdx] *
1967 22326 : static_cast<GPtrDiff_t>(nIncr * nBufferDTSize);
1968 44652 : tileIndices[dimIdx] =
1969 22326 : indicesOuterLoop[dimIdx] / m_anBlockSize[dimIdx];
1970 22326 : }
1971 : }
1972 24644 : if (dimIdx > 0)
1973 24150 : goto lbl_return_to_caller;
1974 :
1975 494 : return true;
1976 : }
1977 :
1978 : /************************************************************************/
1979 : /* ZarrArray::IsEmptyTile() */
1980 : /************************************************************************/
1981 :
1982 22820 : bool ZarrArray::IsEmptyTile(const ZarrByteVectorQuickResize &abyTile) const
1983 : {
1984 44459 : if (m_pabyNoData == nullptr || (m_oType.GetClass() == GEDTC_NUMERIC &&
1985 21639 : GetNoDataValueAsDouble() == 0.0))
1986 : {
1987 22512 : const size_t nBytes = abyTile.size();
1988 22512 : size_t i = 0;
1989 25762 : for (; i + (sizeof(size_t) - 1) < nBytes; i += sizeof(size_t))
1990 : {
1991 25115 : if (*reinterpret_cast<const size_t *>(abyTile.data() + i) != 0)
1992 : {
1993 21865 : return false;
1994 : }
1995 : }
1996 1576 : for (; i < nBytes; ++i)
1997 : {
1998 994 : if (abyTile[i] != 0)
1999 : {
2000 65 : return false;
2001 : }
2002 : }
2003 582 : return true;
2004 : }
2005 616 : else if (m_oType.GetClass() == GEDTC_NUMERIC &&
2006 308 : !GDALDataTypeIsComplex(m_oType.GetNumericDataType()))
2007 : {
2008 308 : const int nDTSize = static_cast<int>(m_oType.GetSize());
2009 308 : const size_t nElts = abyTile.size() / nDTSize;
2010 308 : const auto eDT = m_oType.GetNumericDataType();
2011 308 : return GDALBufferHasOnlyNoData(abyTile.data(), GetNoDataValueAsDouble(),
2012 : nElts, // nWidth
2013 : 1, // nHeight
2014 : nElts, // nLineStride
2015 : 1, // nComponents
2016 : nDTSize * 8, // nBitsPerSample
2017 308 : GDALDataTypeIsInteger(eDT)
2018 112 : ? (GDALDataTypeIsSigned(eDT)
2019 112 : ? GSF_SIGNED_INT
2020 : : GSF_UNSIGNED_INT)
2021 308 : : GSF_FLOATING_POINT);
2022 : }
2023 0 : return false;
2024 : }
2025 :
2026 : /************************************************************************/
2027 : /* ZarrArray::OpenTilePresenceCache() */
2028 : /************************************************************************/
2029 :
2030 : std::shared_ptr<GDALMDArray>
2031 26678 : ZarrArray::OpenTilePresenceCache(bool bCanCreate) const
2032 : {
2033 26678 : if (m_bHasTriedCacheTilePresenceArray)
2034 26228 : return m_poCacheTilePresenceArray;
2035 450 : m_bHasTriedCacheTilePresenceArray = true;
2036 :
2037 450 : if (m_nTotalTileCount == 1)
2038 241 : return nullptr;
2039 :
2040 418 : std::string osCacheFilename;
2041 418 : auto poRGCache = GetCacheRootGroup(bCanCreate, osCacheFilename);
2042 209 : if (!poRGCache)
2043 199 : return nullptr;
2044 :
2045 10 : const std::string osTilePresenceArrayName(MassageName(GetFullName()) +
2046 20 : "_tile_presence");
2047 20 : auto poTilePresenceArray = poRGCache->OpenMDArray(osTilePresenceArrayName);
2048 20 : const auto eByteDT = GDALExtendedDataType::Create(GDT_Byte);
2049 10 : if (poTilePresenceArray)
2050 : {
2051 8 : bool ok = true;
2052 8 : const auto &apoDimsCache = poTilePresenceArray->GetDimensions();
2053 16 : if (poTilePresenceArray->GetDataType() != eByteDT ||
2054 8 : apoDimsCache.size() != m_aoDims.size())
2055 : {
2056 0 : ok = false;
2057 : }
2058 : else
2059 : {
2060 24 : for (size_t i = 0; i < m_aoDims.size(); i++)
2061 : {
2062 : const auto nExpectedDimSize =
2063 16 : (m_aoDims[i]->GetSize() + m_anBlockSize[i] - 1) /
2064 16 : m_anBlockSize[i];
2065 16 : if (apoDimsCache[i]->GetSize() != nExpectedDimSize)
2066 : {
2067 0 : ok = false;
2068 0 : break;
2069 : }
2070 : }
2071 : }
2072 8 : if (!ok)
2073 : {
2074 0 : CPLError(CE_Failure, CPLE_NotSupported,
2075 : "Array %s in %s has not expected characteristics",
2076 : osTilePresenceArrayName.c_str(), osCacheFilename.c_str());
2077 0 : return nullptr;
2078 : }
2079 :
2080 8 : if (!poTilePresenceArray->GetAttribute("filling_status") && !bCanCreate)
2081 : {
2082 0 : CPLDebug(ZARR_DEBUG_KEY,
2083 : "Cache tile presence array for %s found, but filling not "
2084 : "finished",
2085 0 : GetFullName().c_str());
2086 0 : return nullptr;
2087 : }
2088 :
2089 8 : CPLDebug(ZARR_DEBUG_KEY, "Using cache tile presence for %s",
2090 8 : GetFullName().c_str());
2091 : }
2092 2 : else if (bCanCreate)
2093 : {
2094 2 : int idxDim = 0;
2095 2 : std::string osBlockSize;
2096 2 : std::vector<std::shared_ptr<GDALDimension>> apoNewDims;
2097 6 : for (const auto &poDim : m_aoDims)
2098 : {
2099 4 : auto poNewDim = poRGCache->CreateDimension(
2100 8 : osTilePresenceArrayName + '_' + std::to_string(idxDim),
2101 8 : std::string(), std::string(),
2102 4 : (poDim->GetSize() + m_anBlockSize[idxDim] - 1) /
2103 12 : m_anBlockSize[idxDim]);
2104 4 : if (!poNewDim)
2105 0 : return nullptr;
2106 4 : apoNewDims.emplace_back(poNewDim);
2107 :
2108 4 : if (!osBlockSize.empty())
2109 2 : osBlockSize += ',';
2110 4 : constexpr GUInt64 BLOCKSIZE = 256;
2111 : osBlockSize +=
2112 4 : std::to_string(std::min(poNewDim->GetSize(), BLOCKSIZE));
2113 :
2114 4 : idxDim++;
2115 : }
2116 :
2117 2 : CPLStringList aosOptionsTilePresence;
2118 2 : aosOptionsTilePresence.SetNameValue("BLOCKSIZE", osBlockSize.c_str());
2119 : poTilePresenceArray =
2120 6 : poRGCache->CreateMDArray(osTilePresenceArrayName, apoNewDims,
2121 4 : eByteDT, aosOptionsTilePresence.List());
2122 2 : if (!poTilePresenceArray)
2123 : {
2124 0 : CPLError(CE_Failure, CPLE_NotSupported, "Cannot create %s in %s",
2125 : osTilePresenceArrayName.c_str(), osCacheFilename.c_str());
2126 0 : return nullptr;
2127 : }
2128 2 : poTilePresenceArray->SetNoDataValue(0);
2129 : }
2130 : else
2131 : {
2132 0 : return nullptr;
2133 : }
2134 :
2135 10 : m_poCacheTilePresenceArray = poTilePresenceArray;
2136 :
2137 10 : return poTilePresenceArray;
2138 : }
2139 :
2140 : /************************************************************************/
2141 : /* ZarrArray::CacheTilePresence() */
2142 : /************************************************************************/
2143 :
2144 4 : bool ZarrArray::CacheTilePresence()
2145 : {
2146 4 : if (m_nTotalTileCount == 1)
2147 0 : return true;
2148 :
2149 8 : const std::string osDirectoryName = GetDataDirectory();
2150 :
2151 : auto psDir = std::unique_ptr<VSIDIR, decltype(&VSICloseDir)>(
2152 8 : VSIOpenDir(osDirectoryName.c_str(), -1, nullptr), VSICloseDir);
2153 4 : if (!psDir)
2154 0 : return false;
2155 :
2156 8 : auto poTilePresenceArray = OpenTilePresenceCache(true);
2157 4 : if (!poTilePresenceArray)
2158 : {
2159 0 : return false;
2160 : }
2161 :
2162 4 : if (poTilePresenceArray->GetAttribute("filling_status"))
2163 : {
2164 2 : CPLDebug(ZARR_DEBUG_KEY,
2165 : "CacheTilePresence(): %s already filled. Nothing to do",
2166 2 : poTilePresenceArray->GetName().c_str());
2167 2 : return true;
2168 : }
2169 :
2170 4 : std::vector<GUInt64> anTileIdx(m_aoDims.size());
2171 4 : const std::vector<size_t> anCount(m_aoDims.size(), 1);
2172 4 : const std::vector<GInt64> anArrayStep(m_aoDims.size(), 0);
2173 4 : const std::vector<GPtrDiff_t> anBufferStride(m_aoDims.size(), 0);
2174 2 : const auto &apoDimsCache = poTilePresenceArray->GetDimensions();
2175 4 : const auto eByteDT = GDALExtendedDataType::Create(GDT_Byte);
2176 :
2177 2 : CPLDebug(ZARR_DEBUG_KEY,
2178 : "CacheTilePresence(): Iterating over %s to find which tiles are "
2179 : "present...",
2180 : osDirectoryName.c_str());
2181 2 : uint64_t nCounter = 0;
2182 : const char chSrcFilenameDirSeparator =
2183 2 : VSIGetDirectorySeparator(osDirectoryName.c_str())[0];
2184 14 : while (const VSIDIREntry *psEntry = VSIGetNextDirEntry(psDir.get()))
2185 : {
2186 12 : if (!VSI_ISDIR(psEntry->nMode))
2187 : {
2188 : const CPLStringList aosTokens = GetTileIndicesFromFilename(
2189 0 : CPLString(psEntry->pszName)
2190 9 : .replaceAll(chSrcFilenameDirSeparator, '/')
2191 18 : .c_str());
2192 9 : if (aosTokens.size() == static_cast<int>(m_aoDims.size()))
2193 : {
2194 : // Get tile indices from filename
2195 5 : bool unexpectedIndex = false;
2196 15 : for (int i = 0; i < aosTokens.size(); ++i)
2197 : {
2198 10 : if (CPLGetValueType(aosTokens[i]) != CPL_VALUE_INTEGER)
2199 : {
2200 2 : unexpectedIndex = true;
2201 : }
2202 20 : anTileIdx[i] =
2203 10 : static_cast<GUInt64>(CPLAtoGIntBig(aosTokens[i]));
2204 10 : if (anTileIdx[i] >= apoDimsCache[i]->GetSize())
2205 : {
2206 0 : unexpectedIndex = true;
2207 : }
2208 : }
2209 5 : if (unexpectedIndex)
2210 : {
2211 1 : continue;
2212 : }
2213 :
2214 4 : nCounter++;
2215 4 : if ((nCounter % 1000) == 0)
2216 : {
2217 0 : CPLDebug(ZARR_DEBUG_KEY,
2218 : "CacheTilePresence(): Listing in progress "
2219 : "(last examined %s, at least %.02f %% completed)",
2220 0 : psEntry->pszName,
2221 0 : 100.0 * double(nCounter) /
2222 0 : double(m_nTotalTileCount));
2223 : }
2224 4 : constexpr GByte byOne = 1;
2225 : // CPLDebugOnly(ZARR_DEBUG_KEY, "Marking %s has present",
2226 : // psEntry->pszName);
2227 8 : if (!poTilePresenceArray->Write(
2228 4 : anTileIdx.data(), anCount.data(), anArrayStep.data(),
2229 : anBufferStride.data(), eByteDT, &byOne))
2230 : {
2231 0 : return false;
2232 : }
2233 : }
2234 : }
2235 12 : }
2236 2 : CPLDebug(ZARR_DEBUG_KEY, "CacheTilePresence(): finished");
2237 :
2238 : // Write filling_status attribute
2239 2 : auto poAttr = poTilePresenceArray->CreateAttribute(
2240 4 : "filling_status", {}, GDALExtendedDataType::CreateString(), nullptr);
2241 2 : if (poAttr)
2242 : {
2243 2 : if (nCounter == 0)
2244 0 : poAttr->Write("no_tile_present");
2245 2 : else if (nCounter == m_nTotalTileCount)
2246 0 : poAttr->Write("all_tiles_present");
2247 : else
2248 2 : poAttr->Write("some_tiles_missing");
2249 : }
2250 :
2251 : // Force closing
2252 2 : m_poCacheTilePresenceArray = nullptr;
2253 2 : m_bHasTriedCacheTilePresenceArray = false;
2254 :
2255 2 : return true;
2256 : }
2257 :
2258 : /************************************************************************/
2259 : /* ZarrArray::CreateAttribute() */
2260 : /************************************************************************/
2261 :
2262 115 : std::shared_ptr<GDALAttribute> ZarrArray::CreateAttribute(
2263 : const std::string &osName, const std::vector<GUInt64> &anDimensions,
2264 : const GDALExtendedDataType &oDataType, CSLConstList papszOptions)
2265 : {
2266 115 : if (!CheckValidAndErrorOutIfNot())
2267 0 : return nullptr;
2268 :
2269 115 : if (!m_bUpdatable)
2270 : {
2271 2 : CPLError(CE_Failure, CPLE_NotSupported,
2272 : "Dataset not open in update mode");
2273 2 : return nullptr;
2274 : }
2275 113 : if (anDimensions.size() >= 2)
2276 : {
2277 2 : CPLError(CE_Failure, CPLE_NotSupported,
2278 : "Cannot create attributes of dimension >= 2");
2279 2 : return nullptr;
2280 : }
2281 : return m_oAttrGroup.CreateAttribute(osName, anDimensions, oDataType,
2282 111 : papszOptions);
2283 : }
2284 :
2285 : /************************************************************************/
2286 : /* ZarrGroupBase::DeleteAttribute() */
2287 : /************************************************************************/
2288 :
2289 18 : bool ZarrArray::DeleteAttribute(const std::string &osName, CSLConstList)
2290 : {
2291 18 : if (!CheckValidAndErrorOutIfNot())
2292 0 : return false;
2293 :
2294 18 : if (!m_bUpdatable)
2295 : {
2296 6 : CPLError(CE_Failure, CPLE_NotSupported,
2297 : "Dataset not open in update mode");
2298 6 : return false;
2299 : }
2300 :
2301 12 : return m_oAttrGroup.DeleteAttribute(osName);
2302 : }
2303 :
2304 : /************************************************************************/
2305 : /* ZarrArray::SetSpatialRef() */
2306 : /************************************************************************/
2307 :
2308 42 : bool ZarrArray::SetSpatialRef(const OGRSpatialReference *poSRS)
2309 : {
2310 42 : if (!CheckValidAndErrorOutIfNot())
2311 0 : return false;
2312 :
2313 42 : if (!m_bUpdatable)
2314 : {
2315 2 : return GDALPamMDArray::SetSpatialRef(poSRS);
2316 : }
2317 40 : m_poSRS.reset();
2318 40 : if (poSRS)
2319 40 : m_poSRS.reset(poSRS->Clone());
2320 40 : m_bSRSModified = true;
2321 40 : return true;
2322 : }
2323 :
2324 : /************************************************************************/
2325 : /* ZarrArray::SetUnit() */
2326 : /************************************************************************/
2327 :
2328 9 : bool ZarrArray::SetUnit(const std::string &osUnit)
2329 : {
2330 9 : if (!CheckValidAndErrorOutIfNot())
2331 0 : return false;
2332 :
2333 9 : if (!m_bUpdatable)
2334 : {
2335 0 : CPLError(CE_Failure, CPLE_NotSupported,
2336 : "Dataset not open in update mode");
2337 0 : return false;
2338 : }
2339 9 : m_osUnit = osUnit;
2340 9 : m_bUnitModified = true;
2341 9 : return true;
2342 : }
2343 :
2344 : /************************************************************************/
2345 : /* ZarrArray::GetOffset() */
2346 : /************************************************************************/
2347 :
2348 79 : double ZarrArray::GetOffset(bool *pbHasOffset,
2349 : GDALDataType *peStorageType) const
2350 : {
2351 79 : if (pbHasOffset)
2352 79 : *pbHasOffset = m_bHasOffset;
2353 79 : if (peStorageType)
2354 0 : *peStorageType = GDT_Unknown;
2355 79 : return m_dfOffset;
2356 : }
2357 :
2358 : /************************************************************************/
2359 : /* ZarrArray::GetScale() */
2360 : /************************************************************************/
2361 :
2362 77 : double ZarrArray::GetScale(bool *pbHasScale, GDALDataType *peStorageType) const
2363 : {
2364 77 : if (pbHasScale)
2365 77 : *pbHasScale = m_bHasScale;
2366 77 : if (peStorageType)
2367 0 : *peStorageType = GDT_Unknown;
2368 77 : return m_dfScale;
2369 : }
2370 :
2371 : /************************************************************************/
2372 : /* ZarrArray::SetOffset() */
2373 : /************************************************************************/
2374 :
2375 3 : bool ZarrArray::SetOffset(double dfOffset, GDALDataType /* eStorageType */)
2376 : {
2377 3 : if (!CheckValidAndErrorOutIfNot())
2378 0 : return false;
2379 :
2380 3 : m_dfOffset = dfOffset;
2381 3 : m_bHasOffset = true;
2382 3 : m_bOffsetModified = true;
2383 3 : return true;
2384 : }
2385 :
2386 : /************************************************************************/
2387 : /* ZarrArray::SetScale() */
2388 : /************************************************************************/
2389 :
2390 3 : bool ZarrArray::SetScale(double dfScale, GDALDataType /* eStorageType */)
2391 : {
2392 3 : if (!CheckValidAndErrorOutIfNot())
2393 0 : return false;
2394 :
2395 3 : m_dfScale = dfScale;
2396 3 : m_bHasScale = true;
2397 3 : m_bScaleModified = true;
2398 3 : return true;
2399 : }
2400 :
2401 : /************************************************************************/
2402 : /* GetDimensionTypeDirection() */
2403 : /************************************************************************/
2404 :
2405 : /* static */
2406 154 : void ZarrArray::GetDimensionTypeDirection(CPLJSONObject &oAttributes,
2407 : std::string &osType,
2408 : std::string &osDirection)
2409 : {
2410 308 : std::string osUnit;
2411 462 : const auto unit = oAttributes[CF_UNITS];
2412 154 : if (unit.GetType() == CPLJSONObject::Type::String)
2413 : {
2414 52 : osUnit = unit.ToString();
2415 : }
2416 :
2417 462 : const auto oStdName = oAttributes[CF_STD_NAME];
2418 154 : if (oStdName.GetType() == CPLJSONObject::Type::String)
2419 : {
2420 156 : const auto osStdName = oStdName.ToString();
2421 52 : if (osStdName == CF_PROJ_X_COORD || osStdName == CF_LONGITUDE_STD_NAME)
2422 : {
2423 26 : osType = GDAL_DIM_TYPE_HORIZONTAL_X;
2424 26 : oAttributes.Delete(CF_STD_NAME);
2425 26 : if (osUnit == CF_DEGREES_EAST)
2426 : {
2427 24 : osDirection = "EAST";
2428 : }
2429 : }
2430 50 : else if (osStdName == CF_PROJ_Y_COORD ||
2431 24 : osStdName == CF_LATITUDE_STD_NAME)
2432 : {
2433 26 : osType = GDAL_DIM_TYPE_HORIZONTAL_Y;
2434 26 : oAttributes.Delete(CF_STD_NAME);
2435 26 : if (osUnit == CF_DEGREES_NORTH)
2436 : {
2437 24 : osDirection = "NORTH";
2438 : }
2439 : }
2440 0 : else if (osStdName == "time")
2441 : {
2442 0 : osType = GDAL_DIM_TYPE_TEMPORAL;
2443 0 : oAttributes.Delete(CF_STD_NAME);
2444 : }
2445 : }
2446 :
2447 462 : const auto osAxis = oAttributes[CF_AXIS].ToString();
2448 154 : if (osAxis == "Z")
2449 : {
2450 0 : osType = GDAL_DIM_TYPE_VERTICAL;
2451 0 : const auto osPositive = oAttributes["positive"].ToString();
2452 0 : if (osPositive == "up")
2453 : {
2454 0 : osDirection = "UP";
2455 0 : oAttributes.Delete("positive");
2456 : }
2457 0 : else if (osPositive == "down")
2458 : {
2459 0 : osDirection = "DOWN";
2460 0 : oAttributes.Delete("positive");
2461 : }
2462 0 : oAttributes.Delete(CF_AXIS);
2463 : }
2464 154 : }
2465 :
2466 : /************************************************************************/
2467 : /* GetCoordinateVariables() */
2468 : /************************************************************************/
2469 :
2470 : std::vector<std::shared_ptr<GDALMDArray>>
2471 2 : ZarrArray::GetCoordinateVariables() const
2472 : {
2473 2 : if (!CheckValidAndErrorOutIfNot())
2474 0 : return {};
2475 :
2476 4 : std::vector<std::shared_ptr<GDALMDArray>> ret;
2477 6 : const auto poCoordinates = GetAttribute("coordinates");
2478 1 : if (poCoordinates &&
2479 3 : poCoordinates->GetDataType().GetClass() == GEDTC_STRING &&
2480 1 : poCoordinates->GetDimensionCount() == 0)
2481 : {
2482 1 : const char *pszCoordinates = poCoordinates->ReadAsString();
2483 1 : if (pszCoordinates)
2484 : {
2485 2 : auto poGroup = m_poGroupWeak.lock();
2486 1 : if (!poGroup)
2487 : {
2488 0 : CPLError(CE_Failure, CPLE_AppDefined,
2489 : "Cannot access coordinate variables of %s has "
2490 : "belonging group has gone out of scope",
2491 0 : GetName().c_str());
2492 : }
2493 : else
2494 : {
2495 : const CPLStringList aosNames(
2496 2 : CSLTokenizeString2(pszCoordinates, " ", 0));
2497 3 : for (int i = 0; i < aosNames.size(); i++)
2498 : {
2499 6 : auto poCoordinateVar = poGroup->OpenMDArray(aosNames[i]);
2500 2 : if (poCoordinateVar)
2501 : {
2502 2 : ret.emplace_back(poCoordinateVar);
2503 : }
2504 : else
2505 : {
2506 0 : CPLError(CE_Warning, CPLE_AppDefined,
2507 : "Cannot find variable corresponding to "
2508 : "coordinate %s",
2509 : aosNames[i]);
2510 : }
2511 : }
2512 : }
2513 : }
2514 : }
2515 :
2516 2 : return ret;
2517 : }
2518 :
2519 : /************************************************************************/
2520 : /* Resize() */
2521 : /************************************************************************/
2522 :
2523 16 : bool ZarrArray::Resize(const std::vector<GUInt64> &anNewDimSizes,
2524 : CSLConstList /* papszOptions */)
2525 : {
2526 16 : if (!CheckValidAndErrorOutIfNot())
2527 0 : return false;
2528 :
2529 16 : if (!IsWritable())
2530 : {
2531 3 : CPLError(CE_Failure, CPLE_AppDefined,
2532 : "Resize() not supported on read-only file");
2533 3 : return false;
2534 : }
2535 :
2536 13 : const auto nDimCount = GetDimensionCount();
2537 13 : if (anNewDimSizes.size() != nDimCount)
2538 : {
2539 0 : CPLError(CE_Failure, CPLE_IllegalArg,
2540 : "Not expected number of values in anNewDimSizes.");
2541 0 : return false;
2542 : }
2543 :
2544 13 : auto &dims = GetDimensions();
2545 26 : std::vector<size_t> anGrownDimIdx;
2546 26 : std::map<GDALDimension *, GUInt64> oMapDimToSize;
2547 27 : for (size_t i = 0; i < nDimCount; ++i)
2548 : {
2549 21 : auto oIter = oMapDimToSize.find(dims[i].get());
2550 21 : if (oIter != oMapDimToSize.end() && oIter->second != anNewDimSizes[i])
2551 : {
2552 2 : CPLError(CE_Failure, CPLE_AppDefined,
2553 : "Cannot resize a dimension referenced several times "
2554 : "to different sizes");
2555 7 : return false;
2556 : }
2557 19 : if (anNewDimSizes[i] != dims[i]->GetSize())
2558 : {
2559 14 : if (anNewDimSizes[i] < dims[i]->GetSize())
2560 : {
2561 5 : CPLError(CE_Failure, CPLE_NotSupported,
2562 : "Resize() does not support shrinking the array.");
2563 5 : return false;
2564 : }
2565 :
2566 9 : oMapDimToSize[dims[i].get()] = anNewDimSizes[i];
2567 9 : anGrownDimIdx.push_back(i);
2568 : }
2569 : else
2570 : {
2571 5 : oMapDimToSize[dims[i].get()] = dims[i]->GetSize();
2572 : }
2573 : }
2574 6 : if (!anGrownDimIdx.empty())
2575 : {
2576 6 : m_bDefinitionModified = true;
2577 13 : for (size_t dimIdx : anGrownDimIdx)
2578 : {
2579 14 : auto dim = std::dynamic_pointer_cast<ZarrDimension>(dims[dimIdx]);
2580 7 : if (dim)
2581 : {
2582 7 : dim->SetSize(anNewDimSizes[dimIdx]);
2583 7 : if (dim->GetName() != dim->GetFullName())
2584 : {
2585 : // This is not a local dimension
2586 7 : m_poSharedResource->UpdateDimensionSize(dim);
2587 : }
2588 : }
2589 : else
2590 : {
2591 0 : CPLAssert(false);
2592 : }
2593 : }
2594 : }
2595 6 : return true;
2596 : }
2597 :
2598 : /************************************************************************/
2599 : /* NotifyChildrenOfRenaming() */
2600 : /************************************************************************/
2601 :
2602 15 : void ZarrArray::NotifyChildrenOfRenaming()
2603 : {
2604 15 : m_oAttrGroup.ParentRenamed(m_osFullName);
2605 15 : }
2606 :
2607 : /************************************************************************/
2608 : /* ParentRenamed() */
2609 : /************************************************************************/
2610 :
2611 9 : void ZarrArray::ParentRenamed(const std::string &osNewParentFullName)
2612 : {
2613 9 : GDALMDArray::ParentRenamed(osNewParentFullName);
2614 :
2615 9 : auto poParent = m_poGroupWeak.lock();
2616 : // The parent necessarily exist, since it notified us
2617 9 : CPLAssert(poParent);
2618 :
2619 27 : m_osFilename = CPLFormFilenameSafe(
2620 18 : CPLFormFilenameSafe(poParent->GetDirectoryName().c_str(),
2621 9 : m_osName.c_str(), nullptr)
2622 : .c_str(),
2623 9 : CPLGetFilename(m_osFilename.c_str()), nullptr);
2624 9 : }
2625 :
2626 : /************************************************************************/
2627 : /* Rename() */
2628 : /************************************************************************/
2629 :
2630 21 : bool ZarrArray::Rename(const std::string &osNewName)
2631 : {
2632 21 : if (!CheckValidAndErrorOutIfNot())
2633 0 : return false;
2634 :
2635 21 : if (!m_bUpdatable)
2636 : {
2637 6 : CPLError(CE_Failure, CPLE_NotSupported,
2638 : "Dataset not open in update mode");
2639 6 : return false;
2640 : }
2641 15 : if (!ZarrGroupBase::IsValidObjectName(osNewName))
2642 : {
2643 3 : CPLError(CE_Failure, CPLE_NotSupported, "Invalid array name");
2644 3 : return false;
2645 : }
2646 :
2647 24 : auto poParent = m_poGroupWeak.lock();
2648 12 : if (poParent)
2649 : {
2650 12 : if (!poParent->CheckArrayOrGroupWithSameNameDoesNotExist(osNewName))
2651 6 : return false;
2652 : }
2653 :
2654 : const std::string osRootDirectoryName(
2655 12 : CPLGetDirnameSafe(CPLGetDirnameSafe(m_osFilename.c_str()).c_str()));
2656 : const std::string osOldDirectoryName = CPLFormFilenameSafe(
2657 12 : osRootDirectoryName.c_str(), m_osName.c_str(), nullptr);
2658 : const std::string osNewDirectoryName = CPLFormFilenameSafe(
2659 12 : osRootDirectoryName.c_str(), osNewName.c_str(), nullptr);
2660 :
2661 6 : if (VSIRename(osOldDirectoryName.c_str(), osNewDirectoryName.c_str()) != 0)
2662 : {
2663 0 : CPLError(CE_Failure, CPLE_AppDefined, "Renaming of %s to %s failed",
2664 : osOldDirectoryName.c_str(), osNewDirectoryName.c_str());
2665 0 : return false;
2666 : }
2667 :
2668 6 : m_poSharedResource->RenameZMetadataRecursive(osOldDirectoryName,
2669 : osNewDirectoryName);
2670 :
2671 : m_osFilename =
2672 12 : CPLFormFilenameSafe(osNewDirectoryName.c_str(),
2673 6 : CPLGetFilename(m_osFilename.c_str()), nullptr);
2674 :
2675 6 : if (poParent)
2676 : {
2677 6 : poParent->NotifyArrayRenamed(m_osName, osNewName);
2678 : }
2679 :
2680 6 : BaseRename(osNewName);
2681 :
2682 6 : return true;
2683 : }
2684 :
2685 : /************************************************************************/
2686 : /* NotifyChildrenOfDeletion() */
2687 : /************************************************************************/
2688 :
2689 8 : void ZarrArray::NotifyChildrenOfDeletion()
2690 : {
2691 8 : m_oAttrGroup.ParentDeleted();
2692 8 : }
2693 :
2694 : /************************************************************************/
2695 : /* ParseSpecialAttributes() */
2696 : /************************************************************************/
2697 :
2698 839 : void ZarrArray::ParseSpecialAttributes(
2699 : const std::shared_ptr<GDALGroup> &poGroup, CPLJSONObject &oAttributes)
2700 : {
2701 2517 : const auto crs = oAttributes[CRS_ATTRIBUTE_NAME];
2702 839 : std::shared_ptr<OGRSpatialReference> poSRS;
2703 839 : if (crs.GetType() == CPLJSONObject::Type::Object)
2704 : {
2705 48 : for (const char *key : {"url", "wkt", "projjson"})
2706 : {
2707 96 : const auto item = crs[key];
2708 48 : if (item.IsValid())
2709 : {
2710 29 : poSRS = std::make_shared<OGRSpatialReference>();
2711 29 : poSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
2712 58 : if (poSRS->SetFromUserInput(
2713 58 : item.ToString().c_str(),
2714 : OGRSpatialReference::
2715 29 : SET_FROM_USER_INPUT_LIMITATIONS_get()) ==
2716 : OGRERR_NONE)
2717 : {
2718 29 : oAttributes.Delete(CRS_ATTRIBUTE_NAME);
2719 29 : break;
2720 : }
2721 0 : poSRS.reset();
2722 : }
2723 : }
2724 : }
2725 : else
2726 : {
2727 : // Check if SRS is using CF-1 conventions
2728 2430 : const auto gridMapping = oAttributes["grid_mapping"];
2729 810 : if (gridMapping.GetType() == CPLJSONObject::Type::String)
2730 : {
2731 : const auto gridMappingArray =
2732 6 : poGroup->OpenMDArray(gridMapping.ToString());
2733 2 : if (gridMappingArray)
2734 : {
2735 2 : poSRS = std::make_shared<OGRSpatialReference>();
2736 2 : poSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
2737 4 : CPLStringList aosKeyValues;
2738 22 : for (const auto &poAttr : gridMappingArray->GetAttributes())
2739 : {
2740 20 : if (poAttr->GetDataType().GetClass() == GEDTC_STRING)
2741 : {
2742 4 : aosKeyValues.SetNameValue(poAttr->GetName().c_str(),
2743 8 : poAttr->ReadAsString());
2744 : }
2745 16 : else if (poAttr->GetDataType().GetClass() == GEDTC_NUMERIC)
2746 : {
2747 32 : std::string osVal;
2748 32 : for (double val : poAttr->ReadAsDoubleArray())
2749 : {
2750 16 : if (!osVal.empty())
2751 0 : osVal += ',';
2752 16 : osVal += CPLSPrintf("%.17g", val);
2753 : }
2754 16 : aosKeyValues.SetNameValue(poAttr->GetName().c_str(),
2755 32 : osVal.c_str());
2756 : }
2757 : }
2758 2 : if (poSRS->importFromCF1(aosKeyValues.List(), nullptr) !=
2759 : OGRERR_NONE)
2760 : {
2761 0 : poSRS.reset();
2762 : }
2763 : }
2764 : }
2765 : }
2766 :
2767 : // For EOPF Sentinel Zarr Samples Service datasets, read attributes from
2768 : // the STAC Proj extension attributes to get the CRS.
2769 839 : if (!poSRS)
2770 : {
2771 2424 : const auto oProjEPSG = oAttributes["proj:epsg"];
2772 808 : if (oProjEPSG.GetType() == CPLJSONObject::Type::Integer)
2773 : {
2774 1 : poSRS = std::make_shared<OGRSpatialReference>();
2775 1 : poSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
2776 1 : if (poSRS->importFromEPSG(oProjEPSG.ToInteger()) != OGRERR_NONE)
2777 : {
2778 0 : poSRS.reset();
2779 : }
2780 : }
2781 : else
2782 : {
2783 2421 : const auto oProjWKT2 = oAttributes["proj:wkt2"];
2784 807 : if (oProjWKT2.GetType() == CPLJSONObject::Type::String)
2785 : {
2786 1 : poSRS = std::make_shared<OGRSpatialReference>();
2787 1 : poSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
2788 1 : if (poSRS->importFromWkt(oProjWKT2.ToString().c_str()) !=
2789 : OGRERR_NONE)
2790 : {
2791 0 : poSRS.reset();
2792 : }
2793 : }
2794 : }
2795 :
2796 : // There is also a "proj:transform" attribute, but we don't need to
2797 : // use it since the x and y dimensions are already associated with a
2798 : // 1-dimensional array with the values.
2799 : }
2800 :
2801 839 : if (poSRS)
2802 : {
2803 33 : int iDimX = 0;
2804 33 : int iDimY = 0;
2805 33 : int iCount = 1;
2806 99 : for (const auto &poDim : GetDimensions())
2807 : {
2808 66 : if (poDim->GetType() == GDAL_DIM_TYPE_HORIZONTAL_X)
2809 2 : iDimX = iCount;
2810 64 : else if (poDim->GetType() == GDAL_DIM_TYPE_HORIZONTAL_Y)
2811 2 : iDimY = iCount;
2812 66 : iCount++;
2813 : }
2814 33 : if ((iDimX == 0 || iDimY == 0) && GetDimensionCount() >= 2)
2815 : {
2816 31 : iDimX = static_cast<int>(GetDimensionCount());
2817 31 : iDimY = iDimX - 1;
2818 : }
2819 33 : if (iDimX > 0 && iDimY > 0)
2820 : {
2821 33 : const auto &oMapping = poSRS->GetDataAxisToSRSAxisMapping();
2822 92 : if (oMapping == std::vector<int>{2, 1} ||
2823 59 : oMapping == std::vector<int>{2, 1, 3})
2824 7 : poSRS->SetDataAxisToSRSAxisMapping({iDimY, iDimX});
2825 53 : else if (oMapping == std::vector<int>{1, 2} ||
2826 27 : oMapping == std::vector<int>{1, 2, 3})
2827 26 : poSRS->SetDataAxisToSRSAxisMapping({iDimX, iDimY});
2828 : }
2829 :
2830 33 : SetSRS(poSRS);
2831 : }
2832 :
2833 2517 : const auto unit = oAttributes[CF_UNITS];
2834 839 : if (unit.GetType() == CPLJSONObject::Type::String)
2835 : {
2836 165 : std::string osUnit = unit.ToString();
2837 55 : oAttributes.Delete(CF_UNITS);
2838 55 : RegisterUnit(osUnit);
2839 : }
2840 :
2841 2517 : const auto offset = oAttributes[CF_ADD_OFFSET];
2842 839 : const auto offsetType = offset.GetType();
2843 839 : if (offsetType == CPLJSONObject::Type::Integer ||
2844 839 : offsetType == CPLJSONObject::Type::Long ||
2845 : offsetType == CPLJSONObject::Type::Double)
2846 : {
2847 3 : double dfOffset = offset.ToDouble();
2848 3 : oAttributes.Delete(CF_ADD_OFFSET);
2849 3 : RegisterOffset(dfOffset);
2850 : }
2851 :
2852 2517 : const auto scale = oAttributes[CF_SCALE_FACTOR];
2853 839 : const auto scaleType = scale.GetType();
2854 839 : if (scaleType == CPLJSONObject::Type::Integer ||
2855 839 : scaleType == CPLJSONObject::Type::Long ||
2856 : scaleType == CPLJSONObject::Type::Double)
2857 : {
2858 3 : double dfScale = scale.ToDouble();
2859 3 : oAttributes.Delete(CF_SCALE_FACTOR);
2860 3 : RegisterScale(dfScale);
2861 : }
2862 839 : }
2863 :
2864 : /************************************************************************/
2865 : /* SetStatistics() */
2866 : /************************************************************************/
2867 :
2868 1 : bool ZarrArray::SetStatistics(bool bApproxStats, double dfMin, double dfMax,
2869 : double dfMean, double dfStdDev,
2870 : GUInt64 nValidCount, CSLConstList papszOptions)
2871 : {
2872 2 : if (!bApproxStats && m_bUpdatable &&
2873 1 : CPLTestBool(
2874 : CSLFetchNameValueDef(papszOptions, "UPDATE_METADATA", "NO")))
2875 : {
2876 3 : auto poAttr = GetAttribute("actual_range");
2877 1 : if (!poAttr)
2878 : {
2879 : poAttr =
2880 1 : CreateAttribute("actual_range", {2}, GetDataType(), nullptr);
2881 : }
2882 1 : if (poAttr)
2883 : {
2884 2 : std::vector<GUInt64> startIdx = {0};
2885 2 : std::vector<size_t> count = {2};
2886 1 : std::vector<double> values = {dfMin, dfMax};
2887 2 : poAttr->Write(startIdx.data(), count.data(), nullptr, nullptr,
2888 2 : GDALExtendedDataType::Create(GDT_Float64),
2889 1 : values.data(), nullptr, 0);
2890 : }
2891 : }
2892 1 : return GDALPamMDArray::SetStatistics(bApproxStats, dfMin, dfMax, dfMean,
2893 1 : dfStdDev, nValidCount, papszOptions);
2894 : }
2895 :
2896 : /************************************************************************/
2897 : /* ZarrArray::IsTileMissingFromCacheInfo() */
2898 : /************************************************************************/
2899 :
2900 26674 : bool ZarrArray::IsTileMissingFromCacheInfo(const std::string &osFilename,
2901 : const uint64_t *tileIndices) const
2902 : {
2903 26674 : CPL_IGNORE_RET_VAL(osFilename);
2904 53348 : auto poTilePresenceArray = OpenTilePresenceCache(false);
2905 26674 : if (poTilePresenceArray)
2906 : {
2907 36 : std::vector<GUInt64> anTileIdx(m_aoDims.size());
2908 36 : const std::vector<size_t> anCount(m_aoDims.size(), 1);
2909 36 : const std::vector<GInt64> anArrayStep(m_aoDims.size(), 0);
2910 36 : const std::vector<GPtrDiff_t> anBufferStride(m_aoDims.size(), 0);
2911 36 : const auto eByteDT = GDALExtendedDataType::Create(GDT_Byte);
2912 108 : for (size_t i = 0; i < m_aoDims.size(); ++i)
2913 : {
2914 72 : anTileIdx[i] = static_cast<GUInt64>(tileIndices[i]);
2915 : }
2916 36 : GByte byValue = 0;
2917 36 : if (poTilePresenceArray->Read(anTileIdx.data(), anCount.data(),
2918 : anArrayStep.data(), anBufferStride.data(),
2919 72 : eByteDT, &byValue) &&
2920 36 : byValue == 0)
2921 : {
2922 26 : CPLDebugOnly(ZARR_DEBUG_KEY, "Tile %s missing (=nodata)",
2923 : osFilename.c_str());
2924 26 : return true;
2925 : }
2926 : }
2927 26648 : return false;
2928 : }
2929 :
2930 : /************************************************************************/
2931 : /* ZarrArray::GetRawBlockInfo() */
2932 : /************************************************************************/
2933 :
2934 13 : bool ZarrArray::GetRawBlockInfo(const uint64_t *panBlockCoordinates,
2935 : GDALMDArrayRawBlockInfo &info) const
2936 : {
2937 13 : info.clear();
2938 26 : for (size_t i = 0; i < m_anBlockSize.size(); ++i)
2939 : {
2940 15 : const auto nBlockSize = m_anBlockSize[i];
2941 : const auto nBlockCount =
2942 15 : cpl::div_round_up(m_aoDims[i]->GetSize(), nBlockSize);
2943 15 : if (panBlockCoordinates[i] >= nBlockCount)
2944 : {
2945 2 : CPLError(CE_Failure, CPLE_AppDefined,
2946 : "GetRawBlockInfo() failed: array %s: "
2947 : "invalid block coordinate (%u) for dimension %u",
2948 2 : GetName().c_str(),
2949 2 : static_cast<unsigned>(panBlockCoordinates[i]),
2950 : static_cast<unsigned>(i));
2951 2 : return false;
2952 : }
2953 : }
2954 :
2955 22 : std::string osFilename = BuildTileFilename(panBlockCoordinates);
2956 :
2957 : // For network file systems, get the streaming version of the filename,
2958 : // as we don't need arbitrary seeking in the file
2959 11 : osFilename = VSIFileManager::GetHandler(osFilename.c_str())
2960 11 : ->GetStreamingFilename(osFilename);
2961 : {
2962 11 : std::lock_guard<std::mutex> oLock(m_oMutex);
2963 11 : if (IsTileMissingFromCacheInfo(osFilename, panBlockCoordinates))
2964 0 : return true;
2965 : }
2966 :
2967 11 : VSILFILE *fp = nullptr;
2968 : // This is the number of files returned in a S3 directory listing operation
2969 11 : const char *const apszOpenOptions[] = {"IGNORE_FILENAME_RESTRICTIONS=YES",
2970 : nullptr};
2971 11 : const auto nErrorBefore = CPLGetErrorCounter();
2972 : {
2973 : // Avoid issuing ReadDir() when a lot of files are expected
2974 : CPLConfigOptionSetter optionSetter("GDAL_DISABLE_READDIR_ON_OPEN",
2975 11 : "YES", true);
2976 11 : fp = VSIFOpenEx2L(osFilename.c_str(), "rb", 0, apszOpenOptions);
2977 : }
2978 11 : if (fp == nullptr)
2979 : {
2980 1 : if (nErrorBefore != CPLGetErrorCounter())
2981 : {
2982 0 : return false;
2983 : }
2984 : else
2985 : {
2986 : // Missing files are OK and indicate nodata_value
2987 1 : return true;
2988 : }
2989 : }
2990 10 : VSIFSeekL(fp, 0, SEEK_END);
2991 10 : const auto nFileSize = VSIFTellL(fp);
2992 10 : VSIFCloseL(fp);
2993 :
2994 : // For Kerchunk files, get information on the actual location
2995 : const CPLStringList aosMetadata(
2996 10 : VSIGetFileMetadata(osFilename.c_str(), "CHUNK_INFO", nullptr));
2997 10 : if (!aosMetadata.empty())
2998 : {
2999 2 : const char *pszFilename = aosMetadata.FetchNameValue("FILENAME");
3000 2 : if (pszFilename)
3001 1 : info.pszFilename = CPLStrdup(pszFilename);
3002 2 : info.nOffset = std::strtoull(
3003 : aosMetadata.FetchNameValueDef("OFFSET", "0"), nullptr, 10);
3004 2 : info.nSize = std::strtoull(aosMetadata.FetchNameValueDef("SIZE", "0"),
3005 : nullptr, 10);
3006 2 : const char *pszBase64 = aosMetadata.FetchNameValue("BASE64");
3007 2 : if (pszBase64)
3008 : {
3009 1 : const size_t nSizeBase64 = strlen(pszBase64) + 1;
3010 1 : info.pabyInlineData = static_cast<GByte *>(CPLMalloc(nSizeBase64));
3011 1 : memcpy(info.pabyInlineData, pszBase64, nSizeBase64);
3012 : const int nDecodedSize =
3013 1 : CPLBase64DecodeInPlace(info.pabyInlineData);
3014 1 : CPLAssert(static_cast<size_t>(nDecodedSize) ==
3015 : static_cast<size_t>(info.nSize));
3016 1 : CPL_IGNORE_RET_VAL(nDecodedSize);
3017 : }
3018 : }
3019 : else
3020 : {
3021 8 : info.pszFilename = CPLStrdup(osFilename.c_str());
3022 8 : info.nOffset = 0;
3023 8 : info.nSize = nFileSize;
3024 : }
3025 :
3026 10 : info.papszInfo = CSLDuplicate(GetRawBlockInfoInfo().List());
3027 :
3028 10 : return true;
3029 : }
|