Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: GDAL TileDB Driver
4 : * Purpose: Implement GDAL TileDB Support based on https://www.tiledb.io
5 : * Author: TileDB, Inc
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2019, TileDB, Inc
9 : *
10 : * SPDX-License-Identifier: MIT
11 : ****************************************************************************/
12 :
13 : #include <cassert>
14 : #include <cmath>
15 : #include <cinttypes>
16 : #include <limits>
17 :
18 : #include "gdal_priv_templates.hpp"
19 : #include "tiledbheaders.h"
20 :
21 : // XML element inside _gdal XML metadata to store the number of overview levels
22 : constexpr const char *OVERVIEW_COUNT_KEY = "tiledb:OverviewCount";
23 :
24 : /************************************************************************/
25 : /* ==================================================================== */
26 : /* TileDBRasterBand */
27 : /* ==================================================================== */
28 : /************************************************************************/
29 :
30 : class TileDBRasterBand final : public GDALPamRasterBand
31 : {
32 : friend class TileDBRasterDataset;
33 :
34 : protected:
35 : TileDBRasterDataset *poGDS;
36 : bool bStats;
37 : CPLString osAttrName;
38 : double m_dfNoData = 0;
39 : bool m_bNoDataSet = false;
40 :
41 : CPL_DISALLOW_COPY_ASSIGN(TileDBRasterBand)
42 :
43 : public:
44 : TileDBRasterBand(TileDBRasterDataset *, int,
45 : const std::string &osAttr = TILEDB_VALUES);
46 : virtual CPLErr IReadBlock(int, int, void *) override;
47 : virtual CPLErr IWriteBlock(int, int, void *) override;
48 : virtual CPLErr IRasterIO(GDALRWFlag, int, int, int, int, void *, int, int,
49 : GDALDataType, GSpacing, GSpacing,
50 : GDALRasterIOExtraArg *psExtraArg) override;
51 : virtual GDALColorInterp GetColorInterpretation() override;
52 :
53 : double GetNoDataValue(int *pbHasNoData) override;
54 : CPLErr SetNoDataValue(double dfNoData) override;
55 :
56 : int GetOverviewCount() override;
57 : GDALRasterBand *GetOverview(int nIdx) override;
58 : };
59 :
60 226 : static CPLErr option_to_index_type(const char *pszIndexingType,
61 : TILEDB_INTERLEAVE_MODE &eMode)
62 : {
63 226 : CPLErr eErr = CE_None;
64 :
65 226 : if (pszIndexingType)
66 : {
67 139 : if EQUAL (pszIndexingType, "BAND")
68 108 : eMode = BAND;
69 31 : else if EQUAL (pszIndexingType, "ATTRIBUTES")
70 8 : eMode = ATTRIBUTES;
71 23 : else if EQUAL (pszIndexingType, "PIXEL")
72 23 : eMode = PIXEL;
73 : else
74 : {
75 0 : eErr = CE_Failure;
76 0 : CPLError(eErr, CPLE_AppDefined,
77 : "Unable to identify TileDB index mode %s.",
78 : pszIndexingType);
79 : }
80 : }
81 : else
82 : {
83 87 : eMode = BAND;
84 : }
85 :
86 226 : return eErr;
87 : }
88 :
89 106 : static const char *index_type_name(TILEDB_INTERLEAVE_MODE eMode)
90 : {
91 106 : switch (eMode)
92 : {
93 8 : case PIXEL:
94 8 : return "PIXEL";
95 5 : case ATTRIBUTES:
96 5 : return "ATTRIBUTES";
97 93 : case BAND:
98 93 : return "BAND";
99 0 : default:
100 0 : return nullptr;
101 : }
102 : }
103 :
104 : /************************************************************************/
105 : /* SetBuffer() */
106 : /************************************************************************/
107 :
108 402 : static CPLErr SetBuffer(tiledb::Query *poQuery, GDALDataType eType,
109 : const CPLString &osAttrName, void *pImage, size_t nSize)
110 : {
111 402 : switch (eType)
112 : {
113 200 : case GDT_Byte:
114 : poQuery->set_data_buffer(
115 200 : osAttrName, reinterpret_cast<unsigned char *>(pImage), nSize);
116 200 : break;
117 2 : case GDT_Int8:
118 : poQuery->set_data_buffer(osAttrName,
119 2 : reinterpret_cast<int8_t *>(pImage), nSize);
120 2 : break;
121 3 : case GDT_UInt16:
122 : poQuery->set_data_buffer(
123 3 : osAttrName, reinterpret_cast<unsigned short *>(pImage), nSize);
124 3 : break;
125 3 : case GDT_UInt32:
126 : poQuery->set_data_buffer(
127 3 : osAttrName, reinterpret_cast<unsigned int *>(pImage), nSize);
128 3 : break;
129 2 : case GDT_UInt64:
130 : poQuery->set_data_buffer(
131 2 : osAttrName, reinterpret_cast<uint64_t *>(pImage), nSize);
132 2 : break;
133 3 : case GDT_Int16:
134 : poQuery->set_data_buffer(osAttrName,
135 3 : reinterpret_cast<short *>(pImage), nSize);
136 3 : break;
137 9 : case GDT_Int32:
138 : poQuery->set_data_buffer(osAttrName,
139 9 : reinterpret_cast<int *>(pImage), nSize);
140 9 : break;
141 2 : case GDT_Int64:
142 : poQuery->set_data_buffer(
143 2 : osAttrName, reinterpret_cast<int64_t *>(pImage), nSize);
144 2 : break;
145 159 : case GDT_Float32:
146 : poQuery->set_data_buffer(osAttrName,
147 159 : reinterpret_cast<float *>(pImage), nSize);
148 159 : break;
149 3 : case GDT_Float64:
150 : poQuery->set_data_buffer(osAttrName,
151 3 : reinterpret_cast<double *>(pImage), nSize);
152 3 : break;
153 3 : case GDT_CInt16:
154 : poQuery->set_data_buffer(
155 3 : osAttrName, reinterpret_cast<short *>(pImage), nSize * 2);
156 3 : break;
157 3 : case GDT_CInt32:
158 : poQuery->set_data_buffer(
159 3 : osAttrName, reinterpret_cast<int *>(pImage), nSize * 2);
160 3 : break;
161 3 : case GDT_CFloat32:
162 : poQuery->set_data_buffer(
163 3 : osAttrName, reinterpret_cast<float *>(pImage), nSize * 2);
164 3 : break;
165 7 : case GDT_CFloat64:
166 : poQuery->set_data_buffer(
167 7 : osAttrName, reinterpret_cast<double *>(pImage), nSize * 2);
168 7 : break;
169 0 : default:
170 0 : return CE_Failure;
171 : }
172 402 : return CE_None;
173 : }
174 :
175 : /************************************************************************/
176 : /* TileDBRasterBand() */
177 : /************************************************************************/
178 :
179 928 : TileDBRasterBand::TileDBRasterBand(TileDBRasterDataset *poDSIn, int nBandIn,
180 928 : const std::string &osAttr)
181 928 : : poGDS(poDSIn), bStats(poDSIn->bStats), osAttrName(osAttr)
182 : {
183 928 : poDS = poDSIn;
184 928 : nBand = nBandIn;
185 928 : eDataType = poGDS->eDataType;
186 928 : if (eDataType == GDT_Unknown)
187 : {
188 : try
189 : {
190 14 : auto attr = (poGDS->m_roArray ? poGDS->m_roArray : poGDS->m_array)
191 14 : ->schema()
192 28 : .attribute(osAttr);
193 14 : switch (attr.type())
194 : {
195 1 : case TILEDB_INT8:
196 1 : eDataType = GDT_Int8;
197 1 : break;
198 1 : case TILEDB_UINT8:
199 1 : eDataType = GDT_Byte;
200 1 : break;
201 2 : case TILEDB_INT16:
202 2 : eDataType =
203 2 : attr.cell_val_num() == 2 ? GDT_CInt16 : GDT_Int16;
204 2 : break;
205 1 : case TILEDB_UINT16:
206 1 : eDataType = GDT_UInt16;
207 1 : break;
208 2 : case TILEDB_INT32:
209 2 : eDataType =
210 2 : attr.cell_val_num() == 2 ? GDT_CInt32 : GDT_Int32;
211 2 : break;
212 1 : case TILEDB_UINT32:
213 1 : eDataType = GDT_UInt32;
214 1 : break;
215 1 : case TILEDB_INT64:
216 1 : eDataType = GDT_Int64;
217 1 : break;
218 1 : case TILEDB_UINT64:
219 1 : eDataType = GDT_UInt64;
220 1 : break;
221 2 : case TILEDB_FLOAT32:
222 2 : eDataType =
223 2 : attr.cell_val_num() == 2 ? GDT_CFloat32 : GDT_Float32;
224 2 : break;
225 2 : case TILEDB_FLOAT64:
226 2 : eDataType =
227 2 : attr.cell_val_num() == 2 ? GDT_CFloat64 : GDT_Float64;
228 2 : break;
229 0 : default:
230 : {
231 0 : const char *pszTypeName = "";
232 0 : tiledb_datatype_to_str(attr.type(), &pszTypeName);
233 0 : CPLError(CE_Failure, CPLE_NotSupported,
234 : "Unhandled TileDB data type: %s", pszTypeName);
235 0 : break;
236 : }
237 : }
238 : }
239 0 : catch (const tiledb::TileDBError &e)
240 : {
241 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
242 : }
243 : }
244 928 : eAccess = poGDS->eAccess;
245 928 : nRasterXSize = poGDS->nRasterXSize;
246 928 : nRasterYSize = poGDS->nRasterYSize;
247 928 : nBlockXSize = poGDS->nBlockXSize;
248 928 : nBlockYSize = poGDS->nBlockYSize;
249 928 : }
250 :
251 : /************************************************************************/
252 : /* IRasterIO() */
253 : /************************************************************************/
254 :
255 344 : CPLErr TileDBRasterBand::IRasterIO(GDALRWFlag eRWFlag, int nXOff, int nYOff,
256 : int nXSize, int nYSize, void *pData,
257 : int nBufXSize, int nBufYSize,
258 : GDALDataType eBufType, GSpacing nPixelSpace,
259 : GSpacing nLineSpace,
260 : GDALRasterIOExtraArg *psExtraArg)
261 : {
262 344 : if (!poGDS->m_bDeferredCreateHasRun)
263 5 : poGDS->DeferredCreate(/* bCreateArray = */ true);
264 344 : if (!poGDS->m_bDeferredCreateHasBeenSuccessful)
265 0 : return CE_Failure;
266 344 : if (!poGDS->m_array)
267 : {
268 1 : CPLError(CE_Failure, CPLE_AppDefined, "Dataset has been closed");
269 1 : return CE_Failure;
270 : }
271 :
272 343 : if (poGDS->eIndexMode == ATTRIBUTES && eRWFlag == GF_Write)
273 : {
274 0 : CPLError(CE_Failure, CPLE_NoWriteAccess,
275 : "Unable to write using band ordered IRasterIO when using "
276 : "interleave 'ATTRIBUTES'.\n");
277 0 : return CE_Failure;
278 : }
279 :
280 343 : const int nBufferDTSize(GDALGetDataTypeSizeBytes(eBufType));
281 :
282 343 : if (eBufType == eDataType && nXSize == nBufXSize && nYSize == nBufYSize &&
283 280 : nBufferDTSize > 0 && nPixelSpace == nBufferDTSize &&
284 280 : nLineSpace == nBufXSize * nPixelSpace)
285 : {
286 280 : const uint64_t nBandIdx = poGDS->nBandStart + nBand - 1;
287 : std::vector<uint64_t> oaSubarray = {
288 : nBandIdx, nBandIdx,
289 280 : (uint64_t)nYOff, (uint64_t)nYOff + nYSize - 1,
290 560 : (uint64_t)nXOff, (uint64_t)nXOff + nXSize - 1};
291 280 : if (poGDS->eIndexMode == PIXEL)
292 84 : std::rotate(oaSubarray.begin(), oaSubarray.begin() + 2,
293 84 : oaSubarray.end());
294 :
295 : try
296 : {
297 280 : tiledb::Context *ctx = poGDS->m_ctx.get();
298 280 : const auto &oArray = poGDS->GetArray(eRWFlag == GF_Write, ctx);
299 :
300 560 : auto poQuery = std::make_unique<tiledb::Query>(*ctx, oArray);
301 560 : tiledb::Subarray subarray(*ctx, oArray);
302 280 : if (poGDS->m_array->schema().domain().ndim() == 3)
303 : {
304 213 : subarray.set_subarray(oaSubarray);
305 : }
306 : else
307 : {
308 268 : subarray.set_subarray(std::vector<uint64_t>(
309 201 : oaSubarray.cbegin() + 2, oaSubarray.cend()));
310 : }
311 280 : poQuery->set_subarray(subarray);
312 :
313 280 : const size_t nValues = static_cast<size_t>(nBufXSize) * nBufYSize;
314 280 : SetBuffer(poQuery.get(), eDataType, osAttrName, pData, nValues);
315 :
316 : // write additional co-registered values
317 560 : std::vector<std::unique_ptr<void, decltype(&VSIFree)>> aBlocks;
318 :
319 280 : if (poGDS->m_lpoAttributeDS.size() > 0)
320 : {
321 12 : const int nXSizeToRead = nXOff + nXSize > nRasterXSize
322 6 : ? nRasterXSize - nXOff
323 : : nXSize;
324 12 : const int nYSizeToRead = nYOff + nYSize > nRasterYSize
325 6 : ? nRasterYSize - nYOff
326 : : nYSize;
327 18 : for (auto const &poAttrDS : poGDS->m_lpoAttributeDS)
328 : {
329 12 : GDALRasterBand *poAttrBand = poAttrDS->GetRasterBand(nBand);
330 12 : GDALDataType eAttrType = poAttrBand->GetRasterDataType();
331 12 : int nBytes = GDALGetDataTypeSizeBytes(eAttrType);
332 12 : void *pAttrBlock = VSI_MALLOC_VERBOSE(nBytes * nValues);
333 :
334 12 : if (pAttrBlock == nullptr)
335 : {
336 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
337 : "Cannot allocate attribute buffer");
338 0 : return CE_Failure;
339 : }
340 12 : aBlocks.emplace_back(pAttrBlock, &VSIFree);
341 :
342 12 : poAttrBand->AdviseRead(nXOff, nYOff, nXSizeToRead,
343 : nYSizeToRead, nXSizeToRead,
344 12 : nYSizeToRead, eAttrType, nullptr);
345 :
346 12 : CPLErr eErr = poAttrBand->RasterIO(
347 : GF_Read, nXOff, nYOff, nXSizeToRead, nYSizeToRead,
348 : pAttrBlock, nXSizeToRead, nYSizeToRead, eAttrType,
349 : nPixelSpace, nLineSpace, psExtraArg);
350 :
351 12 : if (eErr == CE_None)
352 : {
353 : CPLString osName =
354 24 : CPLGetBasenameSafe(poAttrDS->GetDescription());
355 :
356 12 : SetBuffer(poQuery.get(), eAttrType, osName, pAttrBlock,
357 : nValues);
358 : }
359 : else
360 : {
361 0 : return eErr;
362 : }
363 : }
364 : }
365 :
366 280 : if (bStats)
367 0 : tiledb::Stats::enable();
368 :
369 280 : auto status = poQuery->submit();
370 :
371 280 : if (bStats)
372 : {
373 0 : tiledb::Stats::dump(stdout);
374 0 : tiledb::Stats::disable();
375 : }
376 :
377 280 : if (status == tiledb::Query::Status::FAILED)
378 0 : return CE_Failure;
379 : else
380 280 : return CE_None;
381 : }
382 0 : catch (const tiledb::TileDBError &e)
383 : {
384 0 : CPLError(CE_Failure, CPLE_AppDefined,
385 : "TileDB: TileDBRasterBand::IRasterIO() failed: %s",
386 0 : e.what());
387 0 : return CE_Failure;
388 : }
389 : }
390 :
391 63 : return GDALPamRasterBand::IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize,
392 : pData, nBufXSize, nBufYSize, eBufType,
393 63 : nPixelSpace, nLineSpace, psExtraArg);
394 : }
395 :
396 140 : CPLErr TileDBRasterBand::IReadBlock(int nBlockXOff, int nBlockYOff,
397 : void *pImage)
398 : {
399 140 : const int nXOff = nBlockXOff * nBlockXSize;
400 140 : const int nYOff = nBlockYOff * nBlockYSize;
401 140 : const int nDTSize = GDALGetDataTypeSizeBytes(eDataType);
402 280 : return IRasterIO(GF_Read, nXOff, nYOff, nBlockXSize, nBlockYSize, pImage,
403 : nBlockXSize, nBlockYSize, eDataType, nDTSize,
404 140 : nDTSize * nBlockXSize, nullptr);
405 : }
406 :
407 : /************************************************************************/
408 : /* IWriteBlock() */
409 : /************************************************************************/
410 :
411 15 : CPLErr TileDBRasterBand::IWriteBlock(int nBlockXOff, int nBlockYOff,
412 : void *pImage)
413 :
414 : {
415 15 : if (eAccess == GA_ReadOnly)
416 : {
417 0 : CPLError(CE_Failure, CPLE_NoWriteAccess,
418 : "Unable to write block, dataset is opened read only.\n");
419 0 : return CE_Failure;
420 : }
421 :
422 15 : CPLAssert(poGDS != nullptr && nBlockXOff >= 0 && nBlockYOff >= 0 &&
423 : pImage != nullptr);
424 :
425 15 : int nStartX = nBlockXSize * nBlockXOff;
426 15 : int nStartY = nBlockYSize * nBlockYOff;
427 :
428 15 : const int nDTSize = GDALGetDataTypeSizeBytes(eDataType);
429 30 : return IRasterIO(GF_Write, nStartX, nStartY, nBlockXSize, nBlockYSize,
430 : pImage, nBlockXSize, nBlockYSize, eDataType, nDTSize,
431 15 : nDTSize * nBlockXSize, nullptr);
432 : }
433 :
434 : /************************************************************************/
435 : /* GetColorInterpretation() */
436 : /************************************************************************/
437 :
438 70 : GDALColorInterp TileDBRasterBand::GetColorInterpretation()
439 :
440 : {
441 70 : if (poGDS->nBands == 1)
442 20 : return GCI_GrayIndex;
443 :
444 50 : if (nBand == 1)
445 20 : return GCI_RedBand;
446 :
447 30 : else if (nBand == 2)
448 15 : return GCI_GreenBand;
449 :
450 15 : else if (nBand == 3)
451 15 : return GCI_BlueBand;
452 :
453 0 : return GCI_AlphaBand;
454 : }
455 :
456 : /************************************************************************/
457 : /* GetNoDataValue() */
458 : /************************************************************************/
459 :
460 79 : double TileDBRasterBand::GetNoDataValue(int *pbHasNoData)
461 : {
462 79 : if (pbHasNoData)
463 76 : *pbHasNoData = false;
464 79 : if (m_bNoDataSet)
465 : {
466 11 : if (pbHasNoData)
467 11 : *pbHasNoData = true;
468 11 : return m_dfNoData;
469 : }
470 68 : if (!poGDS->m_bDeferredCreateHasBeenSuccessful)
471 0 : return 0.0;
472 68 : if (!poGDS->m_array)
473 0 : return 0.0;
474 68 : double dfNoData = 0.0;
475 : try
476 : {
477 68 : const void *value = nullptr;
478 68 : uint64_t size = 0;
479 : // Caution: 2 below statements must not be combined in a single one,
480 : // as the lifetime of value is linked to the return value of
481 : // attribute()
482 68 : auto attr = (poGDS->m_roArray ? poGDS->m_roArray : poGDS->m_array)
483 68 : ->schema()
484 136 : .attribute(osAttrName);
485 68 : attr.get_fill_value(&value, &size);
486 68 : if (value &&
487 68 : size == static_cast<uint64_t>(GDALGetDataTypeSizeBytes(eDataType)))
488 : {
489 68 : switch (eDataType)
490 : {
491 51 : case GDT_Byte:
492 51 : dfNoData = *static_cast<const uint8_t *>(value);
493 51 : break;
494 1 : case GDT_Int8:
495 1 : dfNoData = *static_cast<const int8_t *>(value);
496 1 : break;
497 1 : case GDT_UInt16:
498 1 : dfNoData = *static_cast<const uint16_t *>(value);
499 1 : break;
500 2 : case GDT_Int16:
501 : case GDT_CInt16:
502 2 : dfNoData = *static_cast<const int16_t *>(value);
503 2 : break;
504 1 : case GDT_UInt32:
505 1 : dfNoData = *static_cast<const uint32_t *>(value);
506 1 : break;
507 2 : case GDT_Int32:
508 : case GDT_CInt32:
509 2 : dfNoData = *static_cast<const int32_t *>(value);
510 2 : break;
511 0 : case GDT_UInt64:
512 0 : dfNoData = static_cast<double>(
513 0 : *static_cast<const uint64_t *>(value));
514 0 : break;
515 0 : case GDT_Int64:
516 0 : dfNoData = static_cast<double>(
517 0 : *static_cast<const int64_t *>(value));
518 0 : break;
519 0 : case GDT_Float16:
520 : case GDT_CFloat16:
521 : // tileDB does not support float16
522 0 : CPLAssert(false);
523 : break;
524 5 : case GDT_Float32:
525 : case GDT_CFloat32:
526 5 : dfNoData = *static_cast<const float *>(value);
527 5 : break;
528 5 : case GDT_Float64:
529 : case GDT_CFloat64:
530 5 : dfNoData = *static_cast<const double *>(value);
531 5 : break;
532 0 : case GDT_Unknown:
533 : case GDT_TypeCount:
534 0 : break;
535 : }
536 68 : if (pbHasNoData)
537 65 : *pbHasNoData = true;
538 : }
539 : }
540 0 : catch (const tiledb::TileDBError &e)
541 : {
542 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
543 : }
544 68 : return dfNoData;
545 : }
546 :
547 : /************************************************************************/
548 : /* IsValidNoData() */
549 : /************************************************************************/
550 :
551 37 : template <class T> static bool IsValidNoData(double dfNoData)
552 : {
553 72 : return GDALIsValueInRange<T>(dfNoData) &&
554 72 : dfNoData == static_cast<double>(static_cast<T>(dfNoData));
555 : }
556 :
557 : /************************************************************************/
558 : /* SetNoDataValue() */
559 : /************************************************************************/
560 :
561 45 : CPLErr TileDBRasterBand::SetNoDataValue(double dfNoData)
562 : {
563 45 : if (poGDS->m_bDeferredCreateHasBeenSuccessful)
564 : {
565 1 : CPLError(CE_Failure, CPLE_NotSupported,
566 : "TileDBRasterBand::SetNoDataValue(): cannot be called after "
567 : "pixel values have been set");
568 1 : return CE_Failure;
569 : }
570 :
571 56 : if (nBand != 1 &&
572 : dfNoData !=
573 12 : cpl::down_cast<TileDBRasterBand *>(poGDS->papoBands[0])->m_dfNoData)
574 : {
575 1 : CPLError(CE_Failure, CPLE_NotSupported,
576 : "TileDBRasterBand::SetNoDataValue(): all bands should have "
577 : "the same nodata value");
578 1 : return CE_Failure;
579 : }
580 :
581 43 : bool bIsValid = false;
582 43 : switch (eDataType)
583 : {
584 26 : case GDT_Byte:
585 26 : bIsValid = IsValidNoData<uint8_t>(dfNoData);
586 26 : break;
587 1 : case GDT_Int8:
588 1 : bIsValid = IsValidNoData<int8_t>(dfNoData);
589 1 : break;
590 1 : case GDT_UInt16:
591 1 : bIsValid = IsValidNoData<uint16_t>(dfNoData);
592 1 : break;
593 2 : case GDT_Int16:
594 : case GDT_CInt16:
595 2 : bIsValid = IsValidNoData<int16_t>(dfNoData);
596 2 : break;
597 1 : case GDT_UInt32:
598 1 : bIsValid = IsValidNoData<uint32_t>(dfNoData);
599 1 : break;
600 2 : case GDT_Int32:
601 : case GDT_CInt32:
602 2 : bIsValid = IsValidNoData<int32_t>(dfNoData);
603 2 : break;
604 0 : case GDT_UInt64:
605 0 : bIsValid = IsValidNoData<uint64_t>(dfNoData);
606 0 : break;
607 0 : case GDT_Int64:
608 0 : bIsValid = IsValidNoData<int64_t>(dfNoData);
609 0 : break;
610 0 : case GDT_Float16:
611 : case GDT_CFloat16:
612 : // tileDB does not support float16
613 0 : bIsValid = CPLIsNan(dfNoData) || IsValidNoData<float>(dfNoData);
614 0 : break;
615 5 : case GDT_Float32:
616 : case GDT_CFloat32:
617 5 : bIsValid = CPLIsNan(dfNoData) || IsValidNoData<float>(dfNoData);
618 5 : break;
619 5 : case GDT_Float64:
620 : case GDT_CFloat64:
621 5 : bIsValid = true;
622 5 : break;
623 0 : case GDT_Unknown:
624 : case GDT_TypeCount:
625 0 : break;
626 : }
627 43 : if (!bIsValid)
628 : {
629 3 : CPLError(CE_Failure, CPLE_NotSupported,
630 : "TileDBRasterBand::SetNoDataValue(): nodata value cannot be "
631 : "stored in band data type");
632 3 : return CE_Failure;
633 : }
634 :
635 40 : m_dfNoData = dfNoData;
636 40 : m_bNoDataSet = true;
637 40 : return CE_None;
638 : }
639 :
640 : /************************************************************************/
641 : /* GetOverviewCount() */
642 : /************************************************************************/
643 :
644 15 : int TileDBRasterBand::GetOverviewCount()
645 : {
646 15 : if (poGDS->m_nOverviewCountFromMetadata > 0)
647 : {
648 7 : poGDS->LoadOverviews();
649 7 : return static_cast<int>(poGDS->m_apoOverviewDS.size());
650 : }
651 8 : return GDALPamRasterBand::GetOverviewCount();
652 : }
653 :
654 : /************************************************************************/
655 : /* GetOverview() */
656 : /************************************************************************/
657 :
658 22 : GDALRasterBand *TileDBRasterBand::GetOverview(int nIdx)
659 : {
660 22 : if (poGDS->m_nOverviewCountFromMetadata > 0)
661 : {
662 21 : poGDS->LoadOverviews();
663 21 : if (nIdx >= 0 && nIdx < static_cast<int>(poGDS->m_apoOverviewDS.size()))
664 : {
665 17 : return poGDS->m_apoOverviewDS[nIdx]->GetRasterBand(nBand);
666 : }
667 4 : return nullptr;
668 : }
669 1 : return GDALPamRasterBand::GetOverview(nIdx);
670 : }
671 :
672 : /************************************************************************/
673 : /* ==================================================================== */
674 : /* TileDBRasterDataset */
675 : /* ==================================================================== */
676 : /************************************************************************/
677 :
678 : /************************************************************************/
679 : /* ~TileDBRasterDataset() */
680 : /************************************************************************/
681 :
682 498 : TileDBRasterDataset::~TileDBRasterDataset()
683 :
684 : {
685 249 : TileDBRasterDataset::CloseDependentDatasets();
686 249 : TileDBRasterDataset::Close();
687 498 : }
688 :
689 : /************************************************************************/
690 : /* GetArray() */
691 : /************************************************************************/
692 :
693 307 : tiledb::Array &TileDBRasterDataset::GetArray(bool bForWrite,
694 : tiledb::Context *&ctx)
695 : {
696 307 : if (eAccess == GA_Update && !bForWrite)
697 : {
698 105 : if (!m_roArray)
699 : {
700 34 : m_roCtx = std::make_unique<tiledb::Context>(m_ctx->config());
701 34 : if (nTimestamp)
702 : {
703 16 : m_roArray = std::make_unique<tiledb::Array>(
704 8 : *m_roCtx, m_osArrayURI, TILEDB_READ,
705 24 : tiledb::TemporalPolicy(tiledb::TimeTravel, nTimestamp));
706 : }
707 : else
708 : {
709 52 : m_roArray = std::make_unique<tiledb::Array>(
710 78 : *m_roCtx, m_osArrayURI, TILEDB_READ);
711 : }
712 : }
713 :
714 105 : ctx = m_roCtx.get();
715 105 : return *(m_roArray.get());
716 : }
717 :
718 202 : ctx = m_ctx.get();
719 202 : return *(m_array.get());
720 : }
721 :
722 : /************************************************************************/
723 : /* IRasterio() */
724 : /************************************************************************/
725 :
726 85 : CPLErr TileDBRasterDataset::IRasterIO(
727 : GDALRWFlag eRWFlag, int nXOff, int nYOff, int nXSize, int nYSize,
728 : void *pData, int nBufXSize, int nBufYSize, GDALDataType eBufType,
729 : int nBandCount, BANDMAP_TYPE panBandMap, GSpacing nPixelSpace,
730 : GSpacing nLineSpace, GSpacing nBandSpace, GDALRasterIOExtraArg *psExtraArg)
731 :
732 : {
733 85 : if (!m_bDeferredCreateHasRun)
734 48 : DeferredCreate(/* bCreateArray = */ true);
735 85 : if (!m_bDeferredCreateHasBeenSuccessful)
736 0 : return CE_Failure;
737 85 : if (!m_array)
738 : {
739 1 : CPLError(CE_Failure, CPLE_AppDefined, "Dataset has been closed");
740 1 : return CE_Failure;
741 : }
742 :
743 : // support special case of writing attributes for bands, all attributes have to be set at once
744 84 : const int nBufferDTSize(GDALGetDataTypeSizeBytes(eBufType));
745 :
746 84 : if (eIndexMode == ATTRIBUTES && nBandCount == nBands &&
747 4 : eBufType == eDataType && nXSize == nBufXSize && nYSize == nBufYSize &&
748 4 : nBufferDTSize > 0 && nPixelSpace == nBufferDTSize &&
749 4 : nLineSpace == nBufXSize * nPixelSpace &&
750 1 : ((nBandSpace == 0 && nBandCount == 1) ||
751 3 : nBandSpace == nBufYSize * nLineSpace))
752 : {
753 : std::vector<uint64_t> oaSubarray = {
754 4 : (uint64_t)nYOff, (uint64_t)nYOff + nYSize - 1, (uint64_t)nXOff,
755 8 : (uint64_t)nXOff + nXSize - 1};
756 :
757 : try
758 : {
759 4 : tiledb::Context *ctx = m_ctx.get();
760 4 : const auto &oArray = GetArray(eRWFlag == GF_Write, ctx);
761 :
762 8 : auto poQuery = std::make_unique<tiledb::Query>(*ctx, oArray);
763 8 : tiledb::Subarray subarray(*ctx, oArray);
764 4 : subarray.set_subarray(oaSubarray);
765 4 : poQuery->set_subarray(subarray);
766 :
767 14 : for (int b = 0; b < nBandCount; b++)
768 : {
769 10 : TileDBRasterBand *poBand = cpl::down_cast<TileDBRasterBand *>(
770 10 : GetRasterBand(panBandMap[b]));
771 10 : const size_t nRegionSize =
772 10 : static_cast<size_t>(nBufXSize) * nBufYSize * nBufferDTSize;
773 10 : SetBuffer(poQuery.get(), eDataType, poBand->osAttrName,
774 10 : static_cast<GByte *>(pData) + b * nRegionSize,
775 : nRegionSize);
776 : }
777 :
778 4 : if (bStats)
779 0 : tiledb::Stats::enable();
780 :
781 4 : auto status = poQuery->submit();
782 :
783 4 : if (bStats)
784 : {
785 0 : tiledb::Stats::dump(stdout);
786 0 : tiledb::Stats::disable();
787 : }
788 :
789 4 : if (status == tiledb::Query::Status::FAILED)
790 0 : return CE_Failure;
791 : else
792 4 : return CE_None;
793 : }
794 0 : catch (const tiledb::TileDBError &e)
795 : {
796 0 : CPLError(CE_Failure, CPLE_AppDefined,
797 : "TileDB: TileDBRasterDataset::IRasterIO() failed: %s",
798 0 : e.what());
799 0 : return CE_Failure;
800 : }
801 : }
802 :
803 80 : return GDALPamDataset::IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize,
804 : pData, nBufXSize, nBufYSize, eBufType,
805 : nBandCount, panBandMap, nPixelSpace,
806 80 : nLineSpace, nBandSpace, psExtraArg);
807 : }
808 :
809 : /************************************************************************/
810 : /* AddDimensions() */
811 : /************************************************************************/
812 :
813 117 : CPLErr TileDBRasterDataset::AddDimensions(tiledb::Domain &domain,
814 : const char *pszAttrName,
815 : tiledb::Dimension &y,
816 : tiledb::Dimension &x,
817 : tiledb::Dimension *poBands)
818 :
819 : {
820 117 : CPLErr eErr = CE_None;
821 :
822 : const bool bHasFillValue =
823 117 : nBands ? cpl::down_cast<TileDBRasterBand *>(papoBands[0])->m_bNoDataSet
824 117 : : false;
825 : const double dfFillValue =
826 117 : nBands ? cpl::down_cast<TileDBRasterBand *>(papoBands[0])->m_dfNoData
827 117 : : 0.0;
828 :
829 117 : switch (eIndexMode)
830 : {
831 6 : case ATTRIBUTES:
832 6 : domain.add_dimensions(y, x);
833 6 : eErr = CreateAttribute(eDataType, pszAttrName, nBands,
834 : bHasFillValue, dfFillValue);
835 6 : break;
836 8 : case PIXEL:
837 8 : assert(poBands);
838 8 : domain.add_dimensions(y, x, *poBands);
839 8 : eErr = CreateAttribute(eDataType, pszAttrName, 1, bHasFillValue,
840 : dfFillValue);
841 8 : break;
842 103 : default: // BAND
843 103 : assert(poBands);
844 103 : domain.add_dimensions(*poBands, y, x);
845 103 : eErr = CreateAttribute(eDataType, pszAttrName, 1, bHasFillValue,
846 : dfFillValue);
847 103 : break;
848 : }
849 :
850 117 : return eErr;
851 : }
852 :
853 : /************************************************************************/
854 : /* Close() */
855 : /************************************************************************/
856 :
857 487 : CPLErr TileDBRasterDataset::Close()
858 : {
859 487 : CPLErr eErr = CE_None;
860 487 : if (nOpenFlags != OPEN_FLAGS_CLOSED)
861 : {
862 249 : if (TileDBRasterDataset::FlushCache(true) != CE_None)
863 0 : eErr = CE_Failure;
864 :
865 249 : if (!m_bDeferredCreateHasRun)
866 62 : DeferredCreate(/* bCreateArray = */ true);
867 :
868 249 : if (nPamFlags & GPF_DIRTY)
869 135 : TrySaveXML();
870 :
871 : try
872 : {
873 249 : if (m_array)
874 : {
875 246 : m_array->close();
876 246 : m_array.reset();
877 : }
878 : }
879 0 : catch (const tiledb::TileDBError &e)
880 : {
881 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
882 : }
883 : try
884 : {
885 249 : if (m_roArray)
886 : {
887 34 : m_roArray->close();
888 34 : m_roArray.reset();
889 : }
890 : }
891 0 : catch (const tiledb::TileDBError &e)
892 : {
893 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
894 : }
895 :
896 249 : if (GDALPamDataset::Close() != CE_None)
897 0 : eErr = CE_Failure;
898 : }
899 487 : return eErr;
900 : }
901 :
902 : /************************************************************************/
903 : /* CloseDependentDatasets() */
904 : /************************************************************************/
905 :
906 249 : int TileDBRasterDataset::CloseDependentDatasets()
907 : {
908 249 : int bDroppedRef = GDALPamDataset::CloseDependentDatasets();
909 249 : if (!m_apoOverviewDS.empty())
910 6 : bDroppedRef = TRUE;
911 249 : m_apoOverviewDS.clear();
912 249 : return bDroppedRef;
913 : }
914 :
915 : /************************************************************************/
916 : /* FlushCache() */
917 : /************************************************************************/
918 :
919 315 : CPLErr TileDBRasterDataset::FlushCache(bool bAtClosing)
920 :
921 : {
922 315 : return BlockBasedFlushCache(bAtClosing);
923 : }
924 :
925 : /************************************************************************/
926 : /* TrySaveXML() */
927 : /************************************************************************/
928 :
929 135 : CPLErr TileDBRasterDataset::TrySaveXML()
930 :
931 : {
932 135 : if (m_array == nullptr)
933 0 : return CE_None;
934 :
935 135 : CPLXMLNode *psTree = nullptr;
936 : try
937 : {
938 135 : nPamFlags &= ~GPF_DIRTY;
939 :
940 135 : if (psPam == nullptr || (nPamFlags & GPF_NOSAVE))
941 0 : return CE_None;
942 :
943 : /* --------------------------------------------------------------------
944 : */
945 : /* Make sure we know the filename we want to store in. */
946 : /* --------------------------------------------------------------------
947 : */
948 135 : if (!BuildPamFilename())
949 0 : return CE_None;
950 :
951 : /* --------------------------------------------------------------------
952 : */
953 : /* Build the XML representation of the auxiliary metadata. */
954 : /* --------------------------------------------------------------------
955 : */
956 135 : psTree = SerializeToXML(nullptr);
957 :
958 135 : if (psTree == nullptr)
959 : {
960 : /* If we have unset all metadata, we have to delete the PAM file */
961 0 : m_array->delete_metadata(GDAL_ATTRIBUTE_NAME);
962 :
963 0 : if (m_bDatasetInGroup)
964 : {
965 0 : tiledb::Group group(*m_ctx, GetDescription(), TILEDB_WRITE);
966 0 : group.delete_metadata(GDAL_ATTRIBUTE_NAME);
967 : }
968 :
969 0 : return CE_None;
970 : }
971 :
972 135 : if (m_poSubDatasetsTree)
973 : {
974 3 : CPLAddXMLChild(psTree,
975 3 : CPLCloneXMLTree(m_poSubDatasetsTree->psChild));
976 : }
977 :
978 : /* --------------------------------------------------------------------
979 : */
980 : /* If we are working with a subdataset, we need to integrate */
981 : /* the subdataset tree within the whole existing pam tree, */
982 : /* after removing any old version of the same subdataset. */
983 : /* --------------------------------------------------------------------
984 : */
985 135 : if (!psPam->osSubdatasetName.empty())
986 : {
987 : CPLXMLNode *psOldTree, *psSubTree;
988 :
989 0 : CPLErrorReset();
990 : {
991 0 : CPLErrorHandlerPusher oQuietError(CPLQuietErrorHandler);
992 0 : psOldTree = CPLParseXMLFile(psPam->pszPamFilename);
993 : }
994 :
995 0 : if (psOldTree == nullptr)
996 : psOldTree =
997 0 : CPLCreateXMLNode(nullptr, CXT_Element, "PAMDataset");
998 :
999 0 : for (psSubTree = psOldTree->psChild; psSubTree != nullptr;
1000 0 : psSubTree = psSubTree->psNext)
1001 : {
1002 0 : if (psSubTree->eType != CXT_Element ||
1003 0 : !EQUAL(psSubTree->pszValue, "Subdataset"))
1004 0 : continue;
1005 :
1006 0 : if (!EQUAL(CPLGetXMLValue(psSubTree, "name", ""),
1007 : psPam->osSubdatasetName))
1008 0 : continue;
1009 :
1010 0 : break;
1011 : }
1012 :
1013 0 : if (psSubTree == nullptr)
1014 : {
1015 : psSubTree =
1016 0 : CPLCreateXMLNode(psOldTree, CXT_Element, "Subdataset");
1017 0 : CPLCreateXMLNode(
1018 : CPLCreateXMLNode(psSubTree, CXT_Attribute, "name"),
1019 0 : CXT_Text, psPam->osSubdatasetName);
1020 : }
1021 :
1022 : CPLXMLNode *psOldPamDataset =
1023 0 : CPLGetXMLNode(psSubTree, "PAMDataset");
1024 0 : if (psOldPamDataset != nullptr)
1025 : {
1026 0 : CPLRemoveXMLChild(psSubTree, psOldPamDataset);
1027 0 : CPLDestroyXMLNode(psOldPamDataset);
1028 : }
1029 :
1030 0 : CPLAddXMLChild(psSubTree, psTree);
1031 0 : psTree = psOldTree;
1032 : }
1033 :
1034 135 : if (m_nOverviewCountFromMetadata)
1035 : {
1036 4 : CPLCreateXMLElementAndValue(
1037 : psTree, OVERVIEW_COUNT_KEY,
1038 : CPLSPrintf("%d", m_nOverviewCountFromMetadata));
1039 : }
1040 :
1041 : /* --------------------------------------------------------------------
1042 : */
1043 : /* Try saving the auxiliary metadata. */
1044 : /* --------------------------------------------------------------------
1045 : */
1046 270 : CPLErrorHandlerPusher oQuietError(CPLQuietErrorHandler);
1047 135 : char *pszTree = CPLSerializeXMLTree(psTree);
1048 :
1049 0 : std::unique_ptr<tiledb::Array> poArrayToClose;
1050 135 : tiledb::Array *poArray = m_array.get();
1051 135 : if (eAccess == GA_ReadOnly)
1052 : {
1053 1 : if (nTimestamp)
1054 : {
1055 0 : poArrayToClose = std::make_unique<tiledb::Array>(
1056 0 : *m_ctx, m_array->uri(), TILEDB_WRITE,
1057 0 : tiledb::TemporalPolicy(tiledb::TimeTravel, nTimestamp));
1058 : }
1059 : else
1060 : {
1061 2 : poArrayToClose = std::make_unique<tiledb::Array>(
1062 3 : *m_ctx, m_array->uri(), TILEDB_WRITE);
1063 : }
1064 1 : poArray = poArrayToClose.get();
1065 : }
1066 :
1067 135 : poArray->put_metadata(DATASET_TYPE_ATTRIBUTE_NAME, TILEDB_STRING_UTF8,
1068 : static_cast<int>(strlen(RASTER_DATASET_TYPE)),
1069 : RASTER_DATASET_TYPE);
1070 135 : poArray->put_metadata(GDAL_ATTRIBUTE_NAME, TILEDB_UINT8,
1071 135 : static_cast<int>(strlen(pszTree)), pszTree);
1072 :
1073 135 : if (poArrayToClose)
1074 1 : poArrayToClose->close();
1075 :
1076 135 : if (m_bDatasetInGroup)
1077 : {
1078 369 : tiledb::Group group(*m_ctx, GetDescription(), TILEDB_WRITE);
1079 123 : group.put_metadata(GDAL_ATTRIBUTE_NAME, TILEDB_STRING_UTF8,
1080 123 : static_cast<int>(strlen(pszTree)), pszTree);
1081 123 : group.close();
1082 : }
1083 :
1084 135 : CPLFree(pszTree);
1085 :
1086 : /* --------------------------------------------------------------------
1087 : */
1088 : /* Cleanup */
1089 : /* --------------------------------------------------------------------
1090 : */
1091 135 : if (psTree)
1092 135 : CPLDestroyXMLNode(psTree);
1093 :
1094 135 : return CE_None;
1095 : }
1096 0 : catch (const std::exception &e)
1097 : {
1098 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
1099 0 : if (psTree)
1100 0 : CPLDestroyXMLNode(psTree);
1101 0 : return CE_Failure;
1102 : }
1103 : }
1104 :
1105 : /************************************************************************/
1106 : /* TryLoadXML() */
1107 : /************************************************************************/
1108 :
1109 129 : CPLErr TileDBRasterDataset::TryLoadXML(CSLConstList papszSiblingFiles)
1110 :
1111 : {
1112 129 : return TryLoadCachedXML(papszSiblingFiles, true);
1113 : }
1114 :
1115 : /************************************************************************/
1116 : /* TryLoadCachedXML() */
1117 : /************************************************************************/
1118 :
1119 257 : CPLErr TileDBRasterDataset::TryLoadCachedXML(CSLConstList /*papszSiblingFiles*/,
1120 : bool bReload)
1121 :
1122 : {
1123 257 : CPLXMLNode *psTree = nullptr;
1124 : try
1125 : {
1126 257 : PamInitialize();
1127 514 : tiledb::VFS vfs(*m_ctx, m_ctx->config());
1128 :
1129 : /* --------------------------------------------------------------------
1130 : */
1131 : /* Clear dirty flag. Generally when we get to this point is */
1132 : /* from a call at the end of the Open() method, and some calls */
1133 : /* may have already marked the PAM info as dirty (for instance */
1134 : /* setting metadata), but really everything to this point is */
1135 : /* reproducible, and so the PAM info should not really be */
1136 : /* thought of as dirty. */
1137 : /* --------------------------------------------------------------------
1138 : */
1139 257 : nPamFlags &= ~GPF_DIRTY;
1140 :
1141 : /* --------------------------------------------------------------------
1142 : */
1143 : /* Try reading the file. */
1144 : /* --------------------------------------------------------------------
1145 : */
1146 257 : if (!BuildPamFilename())
1147 0 : return CE_None;
1148 :
1149 : /* --------------------------------------------------------------------
1150 : */
1151 : /* In case the PAM filename is a .aux.xml file next to the */
1152 : /* physical file and we have a siblings list, then we can skip */
1153 : /* stat'ing the filesystem. */
1154 : /* --------------------------------------------------------------------
1155 : */
1156 257 : CPLErr eLastErr = CPLGetLastErrorType();
1157 257 : int nLastErrNo = CPLGetLastErrorNo();
1158 514 : CPLString osLastErrorMsg = CPLGetLastErrorMsg();
1159 :
1160 257 : CPLErrorReset();
1161 : {
1162 514 : CPLErrorHandlerPusher oQuietError(CPLQuietErrorHandler);
1163 :
1164 257 : if (bReload)
1165 : {
1166 129 : tiledb_datatype_t v_type =
1167 : TILEDB_UINT8; // CPLSerializeXMLTree returns char*
1168 129 : const void *v_r = nullptr;
1169 129 : uint32_t v_num = 0;
1170 23 : auto &oArray = ((eAccess == GA_Update) && (m_roArray))
1171 : ? m_roArray
1172 152 : : m_array;
1173 129 : oArray->get_metadata(GDAL_ATTRIBUTE_NAME, &v_type, &v_num,
1174 : &v_r);
1175 128 : if (v_r && (v_type == TILEDB_UINT8 || v_type == TILEDB_CHAR ||
1176 0 : v_type == TILEDB_STRING_ASCII ||
1177 129 : v_type == TILEDB_STRING_UTF8))
1178 : {
1179 : osMetaDoc =
1180 128 : CPLString(static_cast<const char *>(v_r), v_num);
1181 : }
1182 129 : psTree = CPLParseXMLString(osMetaDoc);
1183 : }
1184 :
1185 258 : if (bReload && psTree == nullptr &&
1186 258 : vfs.is_file(psPam->pszPamFilename))
1187 : {
1188 1 : auto nBytes = vfs.file_size(psPam->pszPamFilename);
1189 2 : tiledb::VFS::filebuf fbuf(vfs);
1190 1 : fbuf.open(psPam->pszPamFilename, std::ios::in);
1191 1 : std::istream is(&fbuf);
1192 1 : osMetaDoc.resize(nBytes);
1193 1 : is.read((char *)osMetaDoc.data(), nBytes);
1194 1 : fbuf.close();
1195 1 : psTree = CPLParseXMLString(osMetaDoc);
1196 : }
1197 :
1198 257 : if (!bReload)
1199 : {
1200 128 : psTree = CPLParseXMLString(osMetaDoc);
1201 : }
1202 : }
1203 257 : CPLErrorReset();
1204 :
1205 257 : if (eLastErr != CE_None)
1206 0 : CPLErrorSetState(eLastErr, nLastErrNo, osLastErrorMsg.c_str());
1207 :
1208 : /* --------------------------------------------------------------------
1209 : */
1210 : /* If we are looking for a subdataset, search for its subtree not.
1211 : */
1212 : /* --------------------------------------------------------------------
1213 : */
1214 257 : if (psTree && !psPam->osSubdatasetName.empty())
1215 : {
1216 11 : CPLXMLNode *psSubTree = psTree->psChild;
1217 :
1218 65 : for (; psSubTree != nullptr; psSubTree = psSubTree->psNext)
1219 : {
1220 64 : if (psSubTree->eType != CXT_Element ||
1221 64 : !EQUAL(psSubTree->pszValue, "Subdataset"))
1222 38 : continue;
1223 :
1224 26 : if (!EQUAL(CPLGetXMLValue(psSubTree, "name", ""),
1225 : psPam->osSubdatasetName))
1226 16 : continue;
1227 :
1228 10 : psSubTree = CPLGetXMLNode(psSubTree, "PAMDataset");
1229 10 : break;
1230 : }
1231 :
1232 11 : if (psSubTree != nullptr)
1233 10 : psSubTree = CPLCloneXMLTree(psSubTree);
1234 :
1235 11 : CPLDestroyXMLNode(psTree);
1236 11 : psTree = psSubTree;
1237 : }
1238 257 : if (psTree == nullptr)
1239 1 : return CE_Failure;
1240 :
1241 256 : m_nOverviewCountFromMetadata =
1242 256 : atoi(CPLGetXMLValue(psTree, OVERVIEW_COUNT_KEY, "0"));
1243 256 : if (m_nOverviewCountFromMetadata)
1244 : {
1245 10 : CPLDebugOnly("TileDB", "OverviewCount = %d",
1246 : m_nOverviewCountFromMetadata);
1247 : }
1248 :
1249 : /* --------------------------------------------------------------------
1250 : */
1251 : /* Initialize ourselves from this XML tree. */
1252 : /* --------------------------------------------------------------------
1253 : */
1254 :
1255 256 : CPLString osPath(CPLGetPathSafe(psPam->pszPamFilename));
1256 256 : const CPLErr eErr = XMLInit(psTree, osPath);
1257 :
1258 256 : CPLDestroyXMLNode(psTree);
1259 :
1260 256 : if (eErr != CE_None)
1261 0 : PamClear();
1262 :
1263 256 : return eErr;
1264 : }
1265 0 : catch (const tiledb::TileDBError &e)
1266 : {
1267 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
1268 0 : if (psTree)
1269 0 : CPLDestroyXMLNode(psTree);
1270 0 : return CE_Failure;
1271 : }
1272 : }
1273 :
1274 : /************************************************************************/
1275 : /* GetMetadata() */
1276 : /************************************************************************/
1277 :
1278 216 : char **TileDBRasterDataset::GetMetadata(const char *pszDomain)
1279 :
1280 : {
1281 216 : if (pszDomain != nullptr && EQUAL(pszDomain, "SUBDATASETS"))
1282 : {
1283 6 : if (m_aosSubdatasetMD.empty())
1284 : {
1285 5 : CSLConstList papszMD = GDALPamDataset::GetMetadata(pszDomain);
1286 8 : for (const auto &[pszKey, pszValue] :
1287 13 : cpl::IterateNameValue(papszMD))
1288 : {
1289 4 : if (STARTS_WITH(pszKey, "SUBDATASET_") &&
1290 4 : EQUAL(pszKey + strlen(pszKey) - strlen("_NAME"), "_NAME") &&
1291 2 : !STARTS_WITH(pszValue, "TILEDB:"))
1292 : {
1293 2 : m_aosSubdatasetMD.AddNameValue(
1294 2 : pszKey, CPLString().Printf("TILEDB:\"%s\":%s",
1295 2 : GetDescription(), pszValue));
1296 : }
1297 : else
1298 : {
1299 2 : m_aosSubdatasetMD.AddNameValue(pszKey, pszValue);
1300 : }
1301 : }
1302 : }
1303 6 : return m_aosSubdatasetMD.List();
1304 : }
1305 : else
1306 : {
1307 210 : return GDALPamDataset::GetMetadata(pszDomain);
1308 : }
1309 : }
1310 :
1311 : /************************************************************************/
1312 : /* Open() */
1313 : /************************************************************************/
1314 :
1315 129 : GDALDataset *TileDBRasterDataset::Open(GDALOpenInfo *poOpenInfo,
1316 : tiledb::Object::Type objectType)
1317 :
1318 : {
1319 : try
1320 : {
1321 129 : return OpenInternal(poOpenInfo, objectType);
1322 : }
1323 0 : catch (const tiledb::TileDBError &e)
1324 : {
1325 0 : CPLError(CE_Failure, CPLE_AppDefined, "TileDB: Open(%s) failed: %s",
1326 0 : poOpenInfo->pszFilename, e.what());
1327 0 : return nullptr;
1328 : }
1329 : }
1330 :
1331 : /************************************************************************/
1332 : /* OpenInternal() */
1333 : /************************************************************************/
1334 :
1335 129 : GDALDataset *TileDBRasterDataset::OpenInternal(GDALOpenInfo *poOpenInfo,
1336 : tiledb::Object::Type objectType)
1337 :
1338 : {
1339 258 : auto poDS = std::make_unique<TileDBRasterDataset>();
1340 :
1341 129 : poDS->m_bDeferredCreateHasRun = true;
1342 129 : poDS->m_bDeferredCreateHasBeenSuccessful = true;
1343 :
1344 : const char *pszConfig =
1345 129 : CSLFetchNameValue(poOpenInfo->papszOpenOptions, "TILEDB_CONFIG");
1346 :
1347 : const char *pszTimestamp =
1348 129 : CSLFetchNameValue(poOpenInfo->papszOpenOptions, "TILEDB_TIMESTAMP");
1349 :
1350 258 : poDS->bStats =
1351 129 : CSLFetchBoolean(poOpenInfo->papszOpenOptions, "STATS", FALSE);
1352 :
1353 129 : if (pszConfig != nullptr)
1354 : {
1355 0 : poDS->m_osConfigFilename = pszConfig;
1356 0 : tiledb::Config cfg(pszConfig);
1357 0 : poDS->m_ctx.reset(new tiledb::Context(cfg));
1358 : }
1359 : else
1360 : {
1361 129 : poDS->m_ctx.reset(new tiledb::Context());
1362 : }
1363 129 : if (pszTimestamp)
1364 18 : poDS->nTimestamp = std::strtoull(pszTimestamp, nullptr, 10);
1365 :
1366 258 : CPLString osURI;
1367 258 : CPLString osSubdataset;
1368 :
1369 258 : std::string osAttrNameTmp;
1370 : const char *pszAttr =
1371 129 : CSLFetchNameValue(poOpenInfo->papszOpenOptions, "TILEDB_ATTRIBUTE");
1372 :
1373 129 : if (STARTS_WITH_CI(poOpenInfo->pszFilename, "TILEDB:") &&
1374 2 : !STARTS_WITH_CI(poOpenInfo->pszFilename, "TILEDB://"))
1375 : {
1376 : // form required read attributes and open file
1377 : // Create a corresponding GDALDataset.
1378 : CPLStringList apszName(
1379 2 : CSLTokenizeString2(poOpenInfo->pszFilename, ":",
1380 2 : CSLT_HONOURSTRINGS | CSLT_PRESERVEESCAPES));
1381 :
1382 2 : if (apszName.size() != 3)
1383 : {
1384 0 : return nullptr;
1385 : }
1386 :
1387 2 : osURI = TileDBDataset::VSI_to_tiledb_uri(apszName[1]);
1388 2 : osSubdataset = apszName[2];
1389 4 : poDS->SetSubdatasetName(osSubdataset.c_str());
1390 : }
1391 : else
1392 : {
1393 127 : if (pszAttr != nullptr)
1394 : {
1395 4 : poDS->SetSubdatasetName(pszAttr);
1396 : }
1397 :
1398 127 : osURI = TileDBDataset::VSI_to_tiledb_uri(poOpenInfo->pszFilename);
1399 : }
1400 :
1401 258 : CPLString osAux = CPLGetBasenameSafe(osURI);
1402 129 : osAux += ".tdb";
1403 :
1404 : // aux file is in array folder
1405 258 : poDS->SetPhysicalFilename(
1406 258 : CPLFormFilenameSafe(osURI, osAux, nullptr).c_str());
1407 : // Initialize any PAM information.
1408 129 : poDS->SetDescription(osURI);
1409 :
1410 129 : poDS->m_osArrayURI = osURI;
1411 129 : if (objectType == tiledb::Object::Type::Group)
1412 : {
1413 116 : poDS->m_bDatasetInGroup = true;
1414 :
1415 : // Full resolution raster array
1416 116 : poDS->m_osArrayURI = CPLFormFilenameSafe(osURI.c_str(), "l_0", nullptr);
1417 : }
1418 :
1419 129 : const tiledb_query_type_t eMode =
1420 129 : poOpenInfo->eAccess == GA_Update ? TILEDB_WRITE : TILEDB_READ;
1421 :
1422 129 : if (poDS->nTimestamp)
1423 : {
1424 36 : poDS->m_array = std::make_unique<tiledb::Array>(
1425 18 : *poDS->m_ctx, poDS->m_osArrayURI, eMode,
1426 54 : tiledb::TemporalPolicy(tiledb::TimeTravel, poDS->nTimestamp));
1427 : }
1428 : else
1429 : {
1430 222 : poDS->m_array = std::make_unique<tiledb::Array>(
1431 222 : *poDS->m_ctx, poDS->m_osArrayURI, eMode);
1432 : }
1433 :
1434 129 : poDS->eAccess = poOpenInfo->eAccess;
1435 :
1436 : // Force opening read-only dataset
1437 129 : if (eMode == TILEDB_WRITE)
1438 : {
1439 23 : tiledb::Context *ctx = nullptr;
1440 23 : poDS->GetArray(false, ctx);
1441 : }
1442 :
1443 : // dependent on PAM metadata for information about array
1444 129 : poDS->TryLoadXML();
1445 :
1446 258 : tiledb::ArraySchema schema = poDS->m_array->schema();
1447 :
1448 129 : char **papszStructMeta = poDS->GetMetadata("IMAGE_STRUCTURE");
1449 129 : const char *pszXSize = CSLFetchNameValue(papszStructMeta, "X_SIZE");
1450 129 : if (pszXSize)
1451 : {
1452 113 : poDS->nRasterXSize = atoi(pszXSize);
1453 113 : if (poDS->nRasterXSize <= 0)
1454 : {
1455 0 : CPLError(CE_Failure, CPLE_AppDefined,
1456 : "Width %i should be greater than zero.",
1457 0 : poDS->nRasterXSize);
1458 0 : return nullptr;
1459 : }
1460 : }
1461 :
1462 129 : const char *pszYSize = CSLFetchNameValue(papszStructMeta, "Y_SIZE");
1463 129 : if (pszYSize)
1464 : {
1465 113 : poDS->nRasterYSize = atoi(pszYSize);
1466 113 : if (poDS->nRasterYSize <= 0)
1467 : {
1468 0 : CPLError(CE_Failure, CPLE_AppDefined,
1469 : "Height %i should be greater than zero.",
1470 0 : poDS->nRasterYSize);
1471 0 : return nullptr;
1472 : }
1473 : }
1474 :
1475 129 : const char *pszNBits = CSLFetchNameValue(papszStructMeta, "NBITS");
1476 129 : if (pszNBits)
1477 : {
1478 113 : poDS->nBitsPerSample = atoi(pszNBits);
1479 : }
1480 :
1481 129 : const char *pszDataType = CSLFetchNameValue(papszStructMeta, "DATA_TYPE");
1482 129 : if (pszDataType)
1483 : {
1484 : // handle the case where arrays have been written with int type
1485 : // (2.5.0)
1486 113 : GDALDataType eDT = GDALGetDataTypeByName(pszDataType);
1487 113 : if (eDT == GDT_Unknown)
1488 : {
1489 1 : int t = atoi(pszDataType);
1490 1 : if ((t > 0) && (t < GDT_TypeCount))
1491 1 : poDS->eDataType = static_cast<GDALDataType>(atoi(pszDataType));
1492 : else
1493 : {
1494 0 : CPLError(CE_Failure, CPLE_AppDefined, "Unknown data type %s.",
1495 : pszDataType);
1496 0 : return nullptr;
1497 : }
1498 : }
1499 : else
1500 : {
1501 112 : poDS->eDataType = eDT;
1502 : }
1503 : }
1504 : else
1505 : {
1506 16 : if (!pszAttr)
1507 : {
1508 16 : if (schema.attribute_num() == 1)
1509 : {
1510 14 : osAttrNameTmp = schema.attribute(0).name();
1511 14 : pszAttr = osAttrNameTmp.c_str();
1512 : }
1513 : }
1514 : }
1515 :
1516 129 : const char *pszIndexMode = CSLFetchNameValue(papszStructMeta, "INTERLEAVE");
1517 :
1518 129 : if (pszIndexMode)
1519 107 : option_to_index_type(pszIndexMode, poDS->eIndexMode);
1520 :
1521 258 : std::vector<tiledb::Dimension> dims = schema.domain().dimensions();
1522 :
1523 129 : int iYDim = 0;
1524 129 : int iXDim = 1;
1525 129 : if ((dims.size() == 2) || (dims.size() == 3))
1526 : {
1527 129 : if (dims.size() == 3)
1528 : {
1529 141 : if ((pszAttr != nullptr) &&
1530 141 : (schema.attributes().count(pszAttr) == 0))
1531 : {
1532 0 : CPLError(CE_Failure, CPLE_NotSupported,
1533 : "%s attribute is not found in TileDB schema.",
1534 : pszAttr);
1535 0 : return nullptr;
1536 : }
1537 :
1538 123 : if (poDS->eIndexMode == PIXEL)
1539 15 : std::rotate(dims.begin(), dims.begin() + 2, dims.end());
1540 :
1541 123 : if (dims[0].type() != TILEDB_UINT64)
1542 : {
1543 0 : const char *pszTypeName = "";
1544 0 : tiledb_datatype_to_str(dims[0].type(), &pszTypeName);
1545 0 : CPLError(CE_Failure, CPLE_NotSupported,
1546 : "Unsupported BAND dimension type: %s", pszTypeName);
1547 0 : return nullptr;
1548 : }
1549 123 : poDS->nBandStart = dims[0].domain<uint64_t>().first;
1550 123 : const uint64_t nBandEnd = dims[0].domain<uint64_t>().second;
1551 246 : if (nBandEnd < poDS->nBandStart ||
1552 123 : nBandEnd - poDS->nBandStart >
1553 123 : static_cast<uint64_t>(std::numeric_limits<int>::max() - 1))
1554 : {
1555 0 : CPLError(CE_Failure, CPLE_NotSupported,
1556 : "Invalid bounds for BAND dimension.");
1557 0 : return nullptr;
1558 : }
1559 123 : poDS->nBands = static_cast<int>(nBandEnd - poDS->nBandStart + 1);
1560 123 : iYDim = 1;
1561 123 : iXDim = 2;
1562 : }
1563 : else
1564 : {
1565 : const char *pszBands =
1566 6 : poDS->GetMetadataItem("NUM_BANDS", "IMAGE_STRUCTURE");
1567 6 : if (pszBands)
1568 : {
1569 1 : poDS->nBands = atoi(pszBands);
1570 : }
1571 :
1572 6 : poDS->eIndexMode = ATTRIBUTES;
1573 : }
1574 : }
1575 : else
1576 : {
1577 0 : CPLError(CE_Failure, CPLE_AppDefined,
1578 : "Wrong number of dimensions %d: expected 2 or 3.",
1579 0 : static_cast<int>(dims.size()));
1580 0 : return nullptr;
1581 : }
1582 :
1583 258 : if (!GDALCheckDatasetDimensions(poDS->nRasterXSize, poDS->nRasterYSize) ||
1584 129 : !GDALCheckBandCount(poDS->nBands, /*bIsZeroAllowed=*/true))
1585 : {
1586 0 : return nullptr;
1587 : }
1588 :
1589 129 : if (dims[iYDim].type() != TILEDB_UINT64)
1590 : {
1591 0 : const char *pszTypeName = "";
1592 0 : tiledb_datatype_to_str(dims[0].type(), &pszTypeName);
1593 0 : CPLError(CE_Failure, CPLE_NotSupported,
1594 : "Unsupported Y dimension type: %s", pszTypeName);
1595 0 : return nullptr;
1596 : }
1597 129 : if (!pszYSize)
1598 : {
1599 16 : const uint64_t nStart = dims[iYDim].domain<uint64_t>().first;
1600 16 : const uint64_t nEnd = dims[iYDim].domain<uint64_t>().second;
1601 32 : if (nStart != 0 ||
1602 16 : nEnd > static_cast<uint64_t>(std::numeric_limits<int>::max() - 1))
1603 : {
1604 0 : CPLError(CE_Failure, CPLE_NotSupported,
1605 : "Invalid bounds for Y dimension.");
1606 0 : return nullptr;
1607 : }
1608 16 : poDS->nRasterYSize = static_cast<int>(nEnd - nStart + 1);
1609 : }
1610 129 : const uint64_t nBlockYSize = dims[iYDim].tile_extent<uint64_t>();
1611 129 : if (nBlockYSize > static_cast<uint64_t>(std::numeric_limits<int>::max()))
1612 : {
1613 0 : CPLError(CE_Failure, CPLE_NotSupported, "Too large block Y size.");
1614 0 : return nullptr;
1615 : }
1616 129 : poDS->nBlockYSize = static_cast<int>(nBlockYSize);
1617 :
1618 129 : if (dims[iXDim].type() != TILEDB_UINT64)
1619 : {
1620 0 : const char *pszTypeName = "";
1621 0 : tiledb_datatype_to_str(dims[0].type(), &pszTypeName);
1622 0 : CPLError(CE_Failure, CPLE_NotSupported,
1623 : "Unsupported Y dimension type: %s", pszTypeName);
1624 0 : return nullptr;
1625 : }
1626 129 : if (!pszXSize)
1627 : {
1628 16 : const uint64_t nStart = dims[iXDim].domain<uint64_t>().first;
1629 16 : const uint64_t nEnd = dims[iXDim].domain<uint64_t>().second;
1630 32 : if (nStart != 0 ||
1631 16 : nEnd > static_cast<uint64_t>(std::numeric_limits<int>::max() - 1))
1632 : {
1633 0 : CPLError(CE_Failure, CPLE_NotSupported,
1634 : "Invalid bounds for X dimension.");
1635 0 : return nullptr;
1636 : }
1637 16 : poDS->nRasterXSize = static_cast<int>(nEnd - nStart + 1);
1638 : }
1639 129 : const uint64_t nBlockXSize = dims[iXDim].tile_extent<uint64_t>();
1640 129 : if (nBlockXSize > static_cast<uint64_t>(std::numeric_limits<int>::max()))
1641 : {
1642 0 : CPLError(CE_Failure, CPLE_NotSupported, "Too large block X size.");
1643 0 : return nullptr;
1644 : }
1645 129 : poDS->nBlockXSize = static_cast<int>(nBlockXSize);
1646 :
1647 129 : poDS->nBlocksX = DIV_ROUND_UP(poDS->nRasterXSize, poDS->nBlockXSize);
1648 129 : poDS->nBlocksY = DIV_ROUND_UP(poDS->nRasterYSize, poDS->nBlockYSize);
1649 :
1650 129 : if (dims.size() == 3)
1651 : {
1652 : // Create band information objects.
1653 844 : for (int i = 1; i <= poDS->nBands; ++i)
1654 : {
1655 721 : if (pszAttr == nullptr)
1656 189 : poDS->SetBand(i, new TileDBRasterBand(poDS.get(), i));
1657 : else
1658 1064 : poDS->SetBand(
1659 1064 : i, new TileDBRasterBand(poDS.get(), i, CPLString(pszAttr)));
1660 : }
1661 : }
1662 : else // subdatasets or only attributes
1663 : {
1664 7 : if ((poOpenInfo->eAccess == GA_Update) &&
1665 1 : (poDS->GetMetadata("SUBDATASETS") != nullptr))
1666 : {
1667 0 : CPLError(CE_Failure, CPLE_NotSupported,
1668 : "The TileDB driver does not support update access "
1669 : "to subdatasets.");
1670 0 : return nullptr;
1671 : }
1672 :
1673 6 : if (!osSubdataset.empty())
1674 : {
1675 : // do we have the attribute in the schema
1676 2 : if (schema.attributes().count(osSubdataset))
1677 : {
1678 0 : poDS->SetBand(
1679 0 : 1, new TileDBRasterBand(poDS.get(), 1, osSubdataset));
1680 : }
1681 : else
1682 : {
1683 2 : if (schema.attributes().count(osSubdataset + "_1"))
1684 : {
1685 : // Create band information objects.
1686 2 : for (int i = 1; i <= poDS->nBands; ++i)
1687 : {
1688 1 : CPLString osAttr = CPLString().Printf(
1689 2 : "%s_%d", osSubdataset.c_str(), i);
1690 2 : poDS->SetBand(
1691 1 : i, new TileDBRasterBand(poDS.get(), i, osAttr));
1692 : }
1693 : }
1694 : else
1695 : {
1696 1 : CPLError(CE_Failure, CPLE_NotSupported,
1697 : "%s attribute is not found in TileDB schema.",
1698 : osSubdataset.c_str());
1699 1 : return nullptr;
1700 : }
1701 : }
1702 : }
1703 : else
1704 : {
1705 4 : char **papszMeta = poDS->GetMetadata("SUBDATASETS");
1706 4 : if (papszMeta != nullptr)
1707 : {
1708 1 : if ((CSLCount(papszMeta) / 2) == 1)
1709 : {
1710 : const char *pszSubDSName =
1711 0 : poDS->m_aosSubdatasetMD.FetchNameValueDef(
1712 : "SUBDATASET_1_NAME", "");
1713 0 : return GDALDataset::FromHandle(
1714 0 : GDALOpen(pszSubDSName, poOpenInfo->eAccess));
1715 : }
1716 : }
1717 3 : else if (poDS->eIndexMode == ATTRIBUTES)
1718 : {
1719 3 : poDS->nBands = schema.attribute_num();
1720 :
1721 : // Create band information objects.
1722 11 : for (int i = 1; i <= poDS->nBands; ++i)
1723 : {
1724 : CPLString osAttr =
1725 24 : TILEDB_VALUES + CPLString().Printf("_%i", i);
1726 16 : poDS->SetBand(i,
1727 8 : new TileDBRasterBand(poDS.get(), i, osAttr));
1728 : }
1729 : }
1730 : else
1731 : {
1732 0 : CPLError(CE_Failure, CPLE_NotSupported,
1733 : "%s is missing required TileDB subdataset metadata.",
1734 : osURI.c_str());
1735 0 : return nullptr;
1736 : }
1737 : }
1738 : }
1739 :
1740 : // reload metadata now that bands are created to populate band metadata
1741 128 : poDS->TryLoadCachedXML(nullptr, false);
1742 :
1743 256 : tiledb::VFS vfs(*poDS->m_ctx, poDS->m_ctx->config());
1744 :
1745 128 : if (!STARTS_WITH_CI(osURI, "TILEDB:") && vfs.is_dir(osURI))
1746 128 : poDS->oOvManager.Initialize(poDS.get(), ":::VIRTUAL:::");
1747 : else
1748 0 : CPLError(CE_Warning, CPLE_AppDefined,
1749 : "Overviews not supported for network writes.");
1750 :
1751 128 : return poDS.release();
1752 : }
1753 :
1754 : /************************************************************************/
1755 : /* CreateAttribute() */
1756 : /************************************************************************/
1757 :
1758 : template <class T, class NoDataT = T>
1759 129 : static tiledb::Attribute CreateAttribute(tiledb::Context &ctx,
1760 : const std::string &osAttrName,
1761 : tiledb::FilterList &filterList,
1762 : bool bHasFillValue, double dfFillValue)
1763 : {
1764 129 : auto attr = tiledb::Attribute::create<T>(ctx, osAttrName, filterList);
1765 129 : if (bHasFillValue && GDALIsValueInRange<NoDataT>(dfFillValue))
1766 : {
1767 30 : const auto nVal = static_cast<NoDataT>(dfFillValue);
1768 30 : if (dfFillValue == static_cast<double>(nVal))
1769 : {
1770 : if constexpr (sizeof(T) == sizeof(NoDataT))
1771 : {
1772 26 : attr.set_fill_value(&nVal, sizeof(nVal));
1773 : }
1774 : else
1775 : {
1776 4 : T aVal = {nVal, nVal};
1777 4 : attr.set_fill_value(&aVal, sizeof(aVal));
1778 : }
1779 : }
1780 : }
1781 129 : return attr;
1782 : }
1783 :
1784 : /************************************************************************/
1785 : /* CreateAttribute() */
1786 : /************************************************************************/
1787 :
1788 123 : CPLErr TileDBRasterDataset::CreateAttribute(GDALDataType eType,
1789 : const CPLString &osAttrName,
1790 : const int nSubRasterCount,
1791 : bool bHasFillValue,
1792 : double dfFillValue)
1793 : {
1794 : try
1795 : {
1796 252 : for (int i = 0; i < nSubRasterCount; ++i)
1797 : {
1798 129 : CPLString osName(osAttrName);
1799 : // a few special cases
1800 : // remove any leading slashes or
1801 : // additional slashes as in the case of hdf5
1802 129 : if STARTS_WITH (osName, "//")
1803 : {
1804 2 : osName = osName.substr(2);
1805 : }
1806 :
1807 129 : osName.replaceAll("/", "_");
1808 129 : CPLString osPrettyName = osName;
1809 :
1810 129 : if ((eIndexMode == ATTRIBUTES) ||
1811 115 : ((m_bHasSubDatasets) && (nSubRasterCount > 1)))
1812 : {
1813 14 : osName = CPLString().Printf("%s_%d", osName.c_str(), i + 1);
1814 : }
1815 :
1816 129 : switch (eType)
1817 : {
1818 53 : case GDT_Byte:
1819 : {
1820 159 : m_schema->add_attribute(::CreateAttribute<unsigned char>(
1821 53 : *m_ctx, osName, *m_filterList, bHasFillValue,
1822 106 : dfFillValue));
1823 53 : nBitsPerSample = 8;
1824 53 : break;
1825 : }
1826 4 : case GDT_Int8:
1827 : {
1828 4 : m_schema->add_attribute(
1829 8 : ::CreateAttribute<int8_t>(*m_ctx, osName, *m_filterList,
1830 8 : bHasFillValue, dfFillValue));
1831 4 : nBitsPerSample = 8;
1832 4 : break;
1833 : }
1834 5 : case GDT_UInt16:
1835 : {
1836 15 : m_schema->add_attribute(::CreateAttribute<uint16_t>(
1837 5 : *m_ctx, osName, *m_filterList, bHasFillValue,
1838 10 : dfFillValue));
1839 5 : nBitsPerSample = 16;
1840 5 : break;
1841 : }
1842 5 : case GDT_UInt32:
1843 : {
1844 15 : m_schema->add_attribute(::CreateAttribute<uint32_t>(
1845 5 : *m_ctx, osName, *m_filterList, bHasFillValue,
1846 10 : dfFillValue));
1847 5 : nBitsPerSample = 32;
1848 5 : break;
1849 : }
1850 4 : case GDT_UInt64:
1851 : {
1852 12 : m_schema->add_attribute(::CreateAttribute<uint64_t>(
1853 4 : *m_ctx, osName, *m_filterList, bHasFillValue,
1854 8 : dfFillValue));
1855 4 : nBitsPerSample = 64;
1856 4 : break;
1857 : }
1858 5 : case GDT_Int16:
1859 : {
1860 15 : m_schema->add_attribute(::CreateAttribute<int16_t>(
1861 5 : *m_ctx, osName, *m_filterList, bHasFillValue,
1862 10 : dfFillValue));
1863 5 : nBitsPerSample = 16;
1864 5 : break;
1865 : }
1866 7 : case GDT_Int32:
1867 : {
1868 21 : m_schema->add_attribute(::CreateAttribute<int32_t>(
1869 7 : *m_ctx, osName, *m_filterList, bHasFillValue,
1870 14 : dfFillValue));
1871 7 : nBitsPerSample = 32;
1872 7 : break;
1873 : }
1874 4 : case GDT_Int64:
1875 : {
1876 12 : m_schema->add_attribute(::CreateAttribute<int64_t>(
1877 4 : *m_ctx, osName, *m_filterList, bHasFillValue,
1878 8 : dfFillValue));
1879 4 : nBitsPerSample = 64;
1880 4 : break;
1881 : }
1882 12 : case GDT_Float32:
1883 : {
1884 12 : m_schema->add_attribute(
1885 24 : ::CreateAttribute<float>(*m_ctx, osName, *m_filterList,
1886 24 : bHasFillValue, dfFillValue));
1887 12 : nBitsPerSample = 32;
1888 12 : break;
1889 : }
1890 8 : case GDT_Float64:
1891 : {
1892 8 : m_schema->add_attribute(
1893 16 : ::CreateAttribute<double>(*m_ctx, osName, *m_filterList,
1894 16 : bHasFillValue, dfFillValue));
1895 8 : nBitsPerSample = 64;
1896 8 : break;
1897 : }
1898 5 : case GDT_CInt16:
1899 : {
1900 5 : m_schema->add_attribute(
1901 10 : ::CreateAttribute<int16_t[2], int16_t>(
1902 5 : *m_ctx, osName, *m_filterList, bHasFillValue,
1903 10 : dfFillValue));
1904 5 : nBitsPerSample = 16;
1905 5 : break;
1906 : }
1907 5 : case GDT_CInt32:
1908 : {
1909 5 : m_schema->add_attribute(
1910 10 : ::CreateAttribute<int32_t[2], int32_t>(
1911 5 : *m_ctx, osName, *m_filterList, bHasFillValue,
1912 10 : dfFillValue));
1913 5 : nBitsPerSample = 32;
1914 5 : break;
1915 : }
1916 5 : case GDT_CFloat32:
1917 : {
1918 15 : m_schema->add_attribute(::CreateAttribute<float[2], float>(
1919 5 : *m_ctx, osName, *m_filterList, bHasFillValue,
1920 10 : dfFillValue));
1921 5 : nBitsPerSample = 32;
1922 5 : break;
1923 : }
1924 7 : case GDT_CFloat64:
1925 : {
1926 7 : m_schema->add_attribute(
1927 14 : ::CreateAttribute<double[2], double>(
1928 7 : *m_ctx, osName, *m_filterList, bHasFillValue,
1929 14 : dfFillValue));
1930 7 : nBitsPerSample = 64;
1931 7 : break;
1932 : }
1933 0 : default:
1934 0 : return CE_Failure;
1935 : }
1936 :
1937 129 : if ((m_bHasSubDatasets) && (i == 0))
1938 : {
1939 8 : CPLString osDim;
1940 8 : switch (nSubRasterCount)
1941 : {
1942 0 : case 2:
1943 0 : osDim.Printf("%dx%d", nRasterXSize, nRasterYSize);
1944 0 : break;
1945 8 : default:
1946 : osDim.Printf("%dx%dx%d", nSubRasterCount, nRasterXSize,
1947 8 : nRasterYSize);
1948 8 : break;
1949 : }
1950 :
1951 8 : const int nSubDataCount = 1 + m_aosSubdatasetMD.size() / 2;
1952 : m_aosSubdatasetMD.SetNameValue(
1953 16 : CPLString().Printf("SUBDATASET_%d_NAME", nSubDataCount),
1954 24 : CPLString().Printf("%s", osPrettyName.c_str()));
1955 :
1956 : m_aosSubdatasetMD.SetNameValue(
1957 16 : CPLString().Printf("SUBDATASET_%d_DESC", nSubDataCount),
1958 8 : CPLString().Printf("[%s] %s (%s)", osDim.c_str(),
1959 : osPrettyName.c_str(),
1960 16 : GDALGetDataTypeName(eType)));
1961 :
1962 : // add to PAM metadata
1963 8 : if (!m_poSubDatasetsTree)
1964 : {
1965 3 : m_poSubDatasetsTree.reset(
1966 : CPLCreateXMLNode(nullptr, CXT_Element, "PAMDataset"));
1967 : }
1968 :
1969 8 : CPLXMLNode *psSubNode = CPLCreateXMLNode(
1970 : m_poSubDatasetsTree.get(), CXT_Element, "Subdataset");
1971 8 : CPLAddXMLAttributeAndValue(psSubNode, "name",
1972 : osPrettyName.c_str());
1973 :
1974 8 : CPLXMLNode *psMetaNode = CPLCreateXMLNode(
1975 : CPLCreateXMLNode(psSubNode, CXT_Element, "PAMDataset"),
1976 : CXT_Element, "Metadata");
1977 8 : CPLAddXMLAttributeAndValue(psMetaNode, "domain",
1978 : "IMAGE_STRUCTURE");
1979 :
1980 8 : CPLAddXMLAttributeAndValue(
1981 : CPLCreateXMLElementAndValue(
1982 : psMetaNode, "MDI",
1983 16 : CPLString().Printf("%d", nRasterXSize)),
1984 : "KEY", "X_SIZE");
1985 :
1986 8 : CPLAddXMLAttributeAndValue(
1987 : CPLCreateXMLElementAndValue(
1988 : psMetaNode, "MDI",
1989 16 : CPLString().Printf("%d", nRasterYSize)),
1990 : "KEY", "Y_SIZE");
1991 :
1992 8 : CPLAddXMLAttributeAndValue(
1993 : CPLCreateXMLElementAndValue(
1994 : psMetaNode, "MDI",
1995 16 : CPLString().Printf("%s", GDALGetDataTypeName(eType))),
1996 : "KEY", "DATA_TYPE");
1997 :
1998 8 : if (m_lpoAttributeDS.size() > 0)
1999 : {
2000 6 : CPLAddXMLAttributeAndValue(
2001 : CPLCreateXMLElementAndValue(
2002 : psMetaNode, "MDI",
2003 12 : CPLString().Printf("%d", nBands)),
2004 : "KEY", "NUM_BANDS");
2005 : }
2006 : else
2007 : {
2008 2 : CPLAddXMLAttributeAndValue(
2009 : CPLCreateXMLElementAndValue(
2010 : psMetaNode, "MDI",
2011 4 : CPLString().Printf("%d", nSubRasterCount)),
2012 : "KEY", "NUM_BANDS");
2013 : }
2014 :
2015 8 : CPLAddXMLAttributeAndValue(
2016 : CPLCreateXMLElementAndValue(
2017 : psMetaNode, "MDI",
2018 16 : CPLString().Printf("%d", nBitsPerSample)),
2019 : "KEY", "NBITS");
2020 : }
2021 : }
2022 123 : return CE_None;
2023 : }
2024 0 : catch (const tiledb::TileDBError &e)
2025 : {
2026 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
2027 0 : return CE_Failure;
2028 : }
2029 : }
2030 :
2031 : /************************************************************************/
2032 : /* SetBlockSize() */
2033 : /************************************************************************/
2034 :
2035 1 : void TileDBRasterDataset::SetBlockSize(GDALRasterBand *poBand,
2036 : CPLStringList &aosOptions)
2037 :
2038 : {
2039 1 : int nX = 0;
2040 1 : int nY = 0;
2041 1 : poBand->GetBlockSize(&nX, &nY);
2042 :
2043 1 : if (!aosOptions.FetchNameValue("BLOCKXSIZE"))
2044 : {
2045 1 : aosOptions.SetNameValue("BLOCKXSIZE", CPLString().Printf("%d", nX));
2046 : }
2047 :
2048 1 : if (!aosOptions.FetchNameValue("BLOCKYSIZE"))
2049 : {
2050 1 : aosOptions.SetNameValue("BLOCKYSIZE", CPLString().Printf("%d", nY));
2051 : }
2052 1 : }
2053 :
2054 : /************************************************************************/
2055 : /* CreateLL() */
2056 : /* */
2057 : /* Shared functionality between TileDBDataset::Create() and */
2058 : /* TileDBDataset::CreateCopy() for creating TileDB array based on */
2059 : /* a set of options and a configuration. */
2060 : /************************************************************************/
2061 :
2062 120 : TileDBRasterDataset *TileDBRasterDataset::CreateLL(const char *pszFilename,
2063 : int nXSize, int nYSize,
2064 : int nBandsIn,
2065 : GDALDataType eType,
2066 : CSLConstList papszOptions)
2067 : {
2068 : try
2069 : {
2070 120 : if ((nXSize <= 0 && nYSize <= 0))
2071 : {
2072 0 : return nullptr;
2073 : }
2074 :
2075 240 : auto poDS = std::make_unique<TileDBRasterDataset>();
2076 120 : poDS->nRasterXSize = nXSize;
2077 120 : poDS->nRasterYSize = nYSize;
2078 120 : poDS->eDataType = eType;
2079 120 : poDS->nBands = nBandsIn;
2080 120 : poDS->eAccess = GA_Update;
2081 :
2082 120 : if (poDS->nBands == 0) // subdatasets
2083 : {
2084 1 : poDS->eIndexMode = ATTRIBUTES;
2085 : }
2086 : else
2087 : {
2088 : const char *pszIndexMode =
2089 119 : CSLFetchNameValue(papszOptions, "INTERLEAVE");
2090 :
2091 119 : if (option_to_index_type(pszIndexMode, poDS->eIndexMode))
2092 0 : return nullptr;
2093 : }
2094 :
2095 : const char *pszConfig =
2096 120 : CSLFetchNameValue(papszOptions, "TILEDB_CONFIG");
2097 120 : if (pszConfig != nullptr)
2098 : {
2099 0 : poDS->m_osConfigFilename = pszConfig;
2100 0 : tiledb::Config cfg(pszConfig);
2101 0 : poDS->m_ctx.reset(new tiledb::Context(cfg));
2102 : }
2103 : else
2104 : {
2105 120 : poDS->m_ctx.reset(new tiledb::Context());
2106 : }
2107 :
2108 120 : if (CPLTestBool(
2109 : CSLFetchNameValueDef(papszOptions, "CREATE_GROUP", "YES")))
2110 : {
2111 110 : poDS->m_bDatasetInGroup = true;
2112 116 : tiledb::create_group(*(poDS->m_ctx.get()), pszFilename);
2113 :
2114 : {
2115 107 : tiledb::Group group(*(poDS->m_ctx.get()), pszFilename,
2116 428 : TILEDB_WRITE);
2117 107 : group.put_metadata(
2118 : DATASET_TYPE_ATTRIBUTE_NAME, TILEDB_STRING_UTF8,
2119 : static_cast<int>(strlen(RASTER_DATASET_TYPE)),
2120 : RASTER_DATASET_TYPE);
2121 :
2122 107 : group.close();
2123 : }
2124 :
2125 : // Full resolution raster array
2126 107 : poDS->m_osArrayURI =
2127 214 : CPLFormFilenameSafe(pszFilename, "l_0", nullptr);
2128 : }
2129 : else
2130 : {
2131 10 : poDS->m_osArrayURI = pszFilename;
2132 : }
2133 :
2134 : const char *pszCompression =
2135 117 : CSLFetchNameValue(papszOptions, "COMPRESSION");
2136 : const char *pszCompressionLevel =
2137 117 : CSLFetchNameValue(papszOptions, "COMPRESSION_LEVEL");
2138 :
2139 : const char *pszBlockXSize =
2140 117 : CSLFetchNameValue(papszOptions, "BLOCKXSIZE");
2141 117 : poDS->nBlockXSize = (pszBlockXSize) ? atoi(pszBlockXSize) : 256;
2142 : const char *pszBlockYSize =
2143 117 : CSLFetchNameValue(papszOptions, "BLOCKYSIZE");
2144 117 : poDS->nBlockYSize = (pszBlockYSize) ? atoi(pszBlockYSize) : 256;
2145 117 : poDS->bStats = CSLFetchBoolean(papszOptions, "STATS", FALSE);
2146 :
2147 : const char *pszTimestamp =
2148 117 : CSLFetchNameValue(papszOptions, "TILEDB_TIMESTAMP");
2149 117 : if (pszTimestamp != nullptr)
2150 2 : poDS->nTimestamp = std::strtoull(pszTimestamp, nullptr, 10);
2151 :
2152 : // set dimensions and attribute type for schema
2153 234 : poDS->m_schema.reset(
2154 117 : new tiledb::ArraySchema(*poDS->m_ctx, TILEDB_DENSE));
2155 117 : poDS->m_schema->set_tile_order(TILEDB_ROW_MAJOR);
2156 117 : poDS->m_schema->set_cell_order(TILEDB_ROW_MAJOR);
2157 :
2158 117 : poDS->m_filterList.reset(new tiledb::FilterList(*poDS->m_ctx));
2159 :
2160 117 : if (pszCompression != nullptr)
2161 : {
2162 0 : int nLevel = (pszCompressionLevel) ? atoi(pszCompressionLevel) : -1;
2163 0 : if (TileDBDataset::AddFilter(*(poDS->m_ctx.get()),
2164 0 : *(poDS->m_filterList.get()),
2165 0 : pszCompression, nLevel) == CE_None)
2166 : {
2167 0 : poDS->SetMetadataItem("COMPRESSION", pszCompression,
2168 0 : "IMAGE_STRUCTURE");
2169 0 : poDS->m_schema->set_coords_filter_list(*poDS->m_filterList);
2170 : }
2171 : }
2172 :
2173 234 : CPLString osAux = CPLGetBasenameSafe(pszFilename);
2174 117 : osAux += ".tdb";
2175 :
2176 234 : poDS->SetPhysicalFilename(
2177 234 : CPLFormFilenameSafe(pszFilename, osAux.c_str(), nullptr).c_str());
2178 :
2179 : // Initialize PAM information.
2180 117 : poDS->SetDescription(pszFilename);
2181 :
2182 : // Note the dimension bounds are inclusive and are expanded to the match
2183 : // the block size
2184 117 : poDS->nBlocksX = DIV_ROUND_UP(nXSize, poDS->nBlockXSize);
2185 117 : poDS->nBlocksY = DIV_ROUND_UP(nYSize, poDS->nBlockYSize);
2186 :
2187 : // register additional attributes to the pixel value, these will be
2188 : // be reported as subdatasets on future reads
2189 : CSLConstList papszAttributes =
2190 117 : CSLFetchNameValueMultiple(papszOptions, "TILEDB_ATTRIBUTE");
2191 121 : for (const char *pszAttribute : cpl::Iterate(papszAttributes))
2192 : {
2193 : // modeling additional attributes as subdatasets
2194 4 : poDS->m_bHasSubDatasets = true;
2195 : // check each attribute is a GDAL source
2196 : std::unique_ptr<GDALDataset> poAttrDS(
2197 8 : GDALDataset::Open(pszAttribute, GA_ReadOnly));
2198 :
2199 4 : if (poAttrDS != nullptr)
2200 : {
2201 : // check each is co-registered
2202 : // candidate band
2203 4 : int nAttrBands = poAttrDS->GetRasterCount();
2204 4 : if (nAttrBands > 0)
2205 : {
2206 4 : GDALRasterBand *poAttrBand = poAttrDS->GetRasterBand(1);
2207 :
2208 4 : if ((poAttrBand->GetXSize() == poDS->nRasterXSize) &&
2209 8 : (poAttrBand->GetYSize() == poDS->nRasterYSize) &&
2210 4 : (poDS->nBands == nAttrBands))
2211 : {
2212 : // could check geotransform, but it is sufficient
2213 : // that cartesian dimensions are equal
2214 4 : poDS->m_lpoAttributeDS.push_back(std::move(poAttrDS));
2215 : }
2216 : else
2217 : {
2218 0 : CPLError(
2219 : CE_Warning, CPLE_AppDefined,
2220 : "Skipping %s as it has a different dimension\n",
2221 : pszAttribute);
2222 : }
2223 : }
2224 : else
2225 : {
2226 0 : CPLError(CE_Warning, CPLE_AppDefined,
2227 : "Skipping %s as it doesn't have any bands\n",
2228 : pszAttribute);
2229 : }
2230 : }
2231 : else
2232 : {
2233 0 : CPLError(CE_Warning, CPLE_AppDefined,
2234 : "Skipping %s, not recognized as a GDAL dataset\n",
2235 : pszAttribute);
2236 : }
2237 : }
2238 :
2239 117 : return poDS.release();
2240 : }
2241 3 : catch (const tiledb::TileDBError &e)
2242 : {
2243 3 : CPLError(CE_Failure, CPLE_AppDefined, "TileDB: %s", e.what());
2244 3 : return nullptr;
2245 : }
2246 : }
2247 :
2248 : /************************************************************************/
2249 : /* DeferredCreate() */
2250 : /* */
2251 : /* Create dimension, domains and attributes. and optionally the array */
2252 : /************************************************************************/
2253 :
2254 120 : bool TileDBRasterDataset::DeferredCreate(bool bCreateArray)
2255 : {
2256 120 : CPLAssert(!m_bDeferredCreateHasRun);
2257 120 : m_bDeferredCreateHasRun = true;
2258 120 : m_bDeferredCreateHasBeenSuccessful = false;
2259 :
2260 : try
2261 : {
2262 : // this driver enforces that all subdatasets are the same size
2263 240 : tiledb::Domain domain(*m_ctx);
2264 :
2265 120 : uint64_t w = (uint64_t)nBlocksX * nBlockXSize - 1;
2266 120 : uint64_t h = (uint64_t)nBlocksY * nBlockYSize - 1;
2267 :
2268 120 : auto d1 = tiledb::Dimension::create<uint64_t>(*m_ctx, "X", {0, w},
2269 483 : uint64_t(nBlockXSize));
2270 117 : auto d2 = tiledb::Dimension::create<uint64_t>(*m_ctx, "Y", {0, h},
2271 468 : uint64_t(nBlockYSize));
2272 :
2273 : {
2274 : CPLErr eErr;
2275 : // Only used for unit test purposes (to check ability of GDAL to read
2276 : // an arbitrary array)
2277 : const char *pszAttrName =
2278 117 : CPLGetConfigOption("TILEDB_ATTRIBUTE", TILEDB_VALUES);
2279 117 : if ((nBands == 0) || (eIndexMode == ATTRIBUTES))
2280 : {
2281 6 : eErr = AddDimensions(domain, pszAttrName, d2, d1, nullptr);
2282 : }
2283 : else
2284 : {
2285 : auto d3 = tiledb::Dimension::create<uint64_t>(
2286 222 : *m_ctx, "BANDS", {1, uint64_t(nBands)}, 1);
2287 111 : eErr = AddDimensions(domain, pszAttrName, d2, d1, &d3);
2288 : }
2289 117 : if (eErr != CE_None)
2290 0 : return false;
2291 : }
2292 :
2293 117 : m_schema->set_domain(domain).set_order(
2294 117 : {{TILEDB_ROW_MAJOR, TILEDB_ROW_MAJOR}});
2295 :
2296 : // register additional attributes to the pixel value, these will be
2297 : // be reported as subdatasets on future reads
2298 121 : for (const auto &poAttrDS : m_lpoAttributeDS)
2299 : {
2300 : const std::string osAttrName =
2301 4 : CPLGetBasenameSafe(poAttrDS->GetDescription());
2302 4 : GDALRasterBand *poAttrBand = poAttrDS->GetRasterBand(1);
2303 4 : int bHasNoData = false;
2304 4 : const double dfNoData = poAttrBand->GetNoDataValue(&bHasNoData);
2305 4 : CreateAttribute(poAttrBand->GetRasterDataType(), osAttrName.c_str(),
2306 4 : 1, CPL_TO_BOOL(bHasNoData), dfNoData);
2307 : }
2308 :
2309 117 : if (bCreateArray)
2310 : {
2311 116 : CreateArray();
2312 : }
2313 :
2314 117 : m_bDeferredCreateHasBeenSuccessful = true;
2315 117 : return true;
2316 : }
2317 3 : catch (const tiledb::TileDBError &e)
2318 : {
2319 3 : CPLError(CE_Failure, CPLE_AppDefined, "TileDB: %s", e.what());
2320 3 : return false;
2321 : }
2322 : }
2323 :
2324 : /************************************************************************/
2325 : /* CreateArray() */
2326 : /************************************************************************/
2327 :
2328 117 : void TileDBRasterDataset::CreateArray()
2329 : {
2330 117 : tiledb::Array::create(m_osArrayURI, *m_schema);
2331 :
2332 117 : if (m_bDatasetInGroup)
2333 : {
2334 321 : tiledb::Group group(*m_ctx, GetDescription(), TILEDB_WRITE);
2335 107 : group.add_member(m_osArrayURI, false);
2336 107 : group.close();
2337 : }
2338 :
2339 117 : if (nTimestamp)
2340 4 : m_array.reset(new tiledb::Array(
2341 2 : *m_ctx, m_osArrayURI, TILEDB_WRITE,
2342 6 : tiledb::TemporalPolicy(tiledb::TimeTravel, nTimestamp)));
2343 : else
2344 115 : m_array.reset(new tiledb::Array(*m_ctx, m_osArrayURI, TILEDB_WRITE));
2345 117 : }
2346 :
2347 : /************************************************************************/
2348 : /* CopySubDatasets() */
2349 : /* */
2350 : /* Copy SubDatasets from src to a TileDBRasterDataset */
2351 : /* */
2352 : /************************************************************************/
2353 :
2354 1 : CPLErr TileDBRasterDataset::CopySubDatasets(GDALDataset *poSrcDS,
2355 : TileDBRasterDataset *poDstDS,
2356 : GDALProgressFunc pfnProgress,
2357 : void *pProgressData)
2358 :
2359 : {
2360 : try
2361 : {
2362 2 : std::vector<std::unique_ptr<GDALDataset>> apoDatasets;
2363 1 : poDstDS->m_bHasSubDatasets = true;
2364 1 : CSLConstList papszSrcSubDatasets = poSrcDS->GetMetadata("SUBDATASETS");
2365 1 : if (!papszSrcSubDatasets)
2366 0 : return CE_Failure;
2367 : const char *pszSubDSName =
2368 1 : CSLFetchNameValue(papszSrcSubDatasets, "SUBDATASET_1_NAME");
2369 1 : if (!pszSubDSName)
2370 0 : return CE_Failure;
2371 :
2372 : CPLStringList apszTokens(CSLTokenizeString2(
2373 2 : pszSubDSName, ":", CSLT_HONOURSTRINGS | CSLT_PRESERVEESCAPES));
2374 : // FIXME? this is tailored for HDF5-like subdataset names
2375 : // HDF5:foo.hdf5:attrname.
2376 1 : if (apszTokens.size() != 3)
2377 : {
2378 0 : CPLError(CE_Failure, CPLE_AppDefined,
2379 : "Cannot guess attribute name in %s", pszSubDSName);
2380 0 : return CE_Failure;
2381 : }
2382 :
2383 : std::unique_ptr<GDALDataset> poSubDataset(
2384 2 : GDALDataset::Open(pszSubDSName));
2385 2 : if (poSubDataset.get() == nullptr ||
2386 1 : poSubDataset->GetRasterCount() == 0)
2387 : {
2388 0 : return CE_Failure;
2389 : }
2390 :
2391 1 : uint64_t nSubXSize = poSubDataset->GetRasterXSize();
2392 1 : uint64_t nSubYSize = poSubDataset->GetRasterYSize();
2393 :
2394 1 : const char *pszAttrName = apszTokens[2];
2395 :
2396 1 : auto poFirstSubDSBand = poSubDataset->GetRasterBand(1);
2397 1 : int bFirstSubDSBandHasNoData = FALSE;
2398 : const double dfFirstSubDSBandNoData =
2399 1 : poFirstSubDSBand->GetNoDataValue(&bFirstSubDSBandHasNoData);
2400 1 : poDstDS->CreateAttribute(poFirstSubDSBand->GetRasterDataType(),
2401 : pszAttrName, poSubDataset->GetRasterCount(),
2402 1 : CPL_TO_BOOL(bFirstSubDSBandHasNoData),
2403 : dfFirstSubDSBandNoData);
2404 1 : apoDatasets.push_back(std::move(poSubDataset));
2405 :
2406 8 : for (const auto &[pszKey, pszValue] :
2407 9 : cpl::IterateNameValue(papszSrcSubDatasets))
2408 : {
2409 4 : if (EQUAL(pszKey, "SUBDATASET_1_NAME") || !strstr(pszKey, "_NAME"))
2410 : {
2411 3 : continue;
2412 : }
2413 1 : pszSubDSName = pszValue;
2414 : apszTokens = CSLTokenizeString2(
2415 1 : pszSubDSName, ":", CSLT_HONOURSTRINGS | CSLT_PRESERVEESCAPES);
2416 1 : if (apszTokens.size() != 3)
2417 : {
2418 0 : CPLError(CE_Failure, CPLE_AppDefined,
2419 : "Cannot guess attribute name in %s", pszSubDSName);
2420 0 : continue;
2421 : }
2422 :
2423 : std::unique_ptr<GDALDataset> poSubDS(
2424 2 : GDALDataset::Open(pszSubDSName));
2425 1 : if ((poSubDS != nullptr) && poSubDS->GetRasterCount() > 0)
2426 : {
2427 1 : GDALRasterBand *poBand = poSubDS->GetRasterBand(1);
2428 : int nBlockXSize, nBlockYSize;
2429 1 : poBand->GetBlockSize(&nBlockXSize, &nBlockYSize);
2430 :
2431 1 : int bHasNoData = FALSE;
2432 1 : const double dfNoData = poBand->GetNoDataValue(&bHasNoData);
2433 :
2434 1 : if ((poSubDS->GetRasterXSize() != (int)nSubXSize) ||
2435 1 : (poSubDS->GetRasterYSize() != (int)nSubYSize) ||
2436 3 : (nBlockXSize != poDstDS->nBlockXSize) ||
2437 1 : (nBlockYSize != poDstDS->nBlockYSize))
2438 : {
2439 0 : CPLError(CE_Warning, CPLE_AppDefined,
2440 : "Sub-datasets must have the same dimension,"
2441 : " and block sizes, skipping %s",
2442 : pszSubDSName);
2443 : }
2444 : else
2445 : {
2446 1 : pszAttrName = apszTokens[2];
2447 1 : poDstDS->CreateAttribute(
2448 : poSubDS->GetRasterBand(1)->GetRasterDataType(),
2449 : pszAttrName, poSubDS->GetRasterCount(),
2450 1 : CPL_TO_BOOL(bHasNoData), dfNoData);
2451 1 : apoDatasets.push_back(std::move(poSubDS));
2452 : }
2453 : }
2454 : else
2455 : {
2456 0 : CPLError(
2457 : CE_Warning, CPLE_AppDefined,
2458 : "Sub-datasets must be not null and contain data in bands,"
2459 : "skipping %s\n",
2460 : pszSubDSName);
2461 : }
2462 : }
2463 :
2464 1 : poDstDS->SetMetadata(poDstDS->m_aosSubdatasetMD.List(), "SUBDATASETS");
2465 :
2466 1 : poDstDS->CreateArray();
2467 :
2468 : /* -------------------------------------------------------- */
2469 : /* Report preliminary (0) progress. */
2470 : /* --------------------------------------------------------- */
2471 1 : if (!pfnProgress(0.0, nullptr, pProgressData))
2472 : {
2473 0 : CPLError(CE_Failure, CPLE_UserInterrupt,
2474 : "User terminated CreateCopy()");
2475 0 : return CE_Failure;
2476 : }
2477 :
2478 : // copy over subdatasets by block
2479 2 : tiledb::Query query(*poDstDS->m_ctx, *poDstDS->m_array);
2480 1 : query.set_layout(TILEDB_GLOBAL_ORDER);
2481 1 : const uint64_t nTotalBlocks =
2482 1 : static_cast<uint64_t>(poDstDS->nBlocksX) * poDstDS->nBlocksY;
2483 1 : uint64_t nBlockCounter = 0;
2484 :
2485 : // row-major
2486 6 : for (int j = 0; j < poDstDS->nBlocksY; ++j)
2487 : {
2488 55 : for (int i = 0; i < poDstDS->nBlocksX; ++i)
2489 : {
2490 50 : std::vector<std::unique_ptr<void, decltype(&VSIFree)>> aBlocks;
2491 : // have to write set all tiledb attributes on write
2492 50 : int iAttr = 0;
2493 150 : for (auto &poSubDS : apoDatasets)
2494 : {
2495 : const GDALDataType eDT =
2496 100 : poSubDS->GetRasterBand(1)->GetRasterDataType();
2497 :
2498 200 : for (int b = 1; b <= poSubDS->GetRasterCount(); ++b)
2499 : {
2500 100 : const int nDTSize = GDALGetDataTypeSizeBytes(eDT);
2501 100 : const size_t nValues =
2502 100 : static_cast<size_t>(poDstDS->nBlockXSize) *
2503 100 : poDstDS->nBlockYSize;
2504 100 : void *pBlock = VSI_MALLOC_VERBOSE(nDTSize * nValues);
2505 100 : if (!pBlock)
2506 0 : return CE_Failure;
2507 100 : aBlocks.emplace_back(pBlock, &VSIFree);
2508 100 : GDALRasterBand *poBand = poSubDS->GetRasterBand(b);
2509 100 : if (poBand->ReadBlock(i, j, pBlock) == CE_None)
2510 : {
2511 100 : SetBuffer(
2512 : &query, eDT,
2513 200 : poDstDS->m_schema->attribute(iAttr).name(),
2514 : pBlock, nValues);
2515 : }
2516 100 : ++iAttr;
2517 : }
2518 : }
2519 :
2520 50 : if (poDstDS->bStats)
2521 0 : tiledb::Stats::enable();
2522 :
2523 50 : auto status = query.submit();
2524 :
2525 50 : if (poDstDS->bStats)
2526 : {
2527 0 : tiledb::Stats::dump(stdout);
2528 0 : tiledb::Stats::disable();
2529 : }
2530 :
2531 50 : if (status == tiledb::Query::Status::FAILED)
2532 : {
2533 0 : return CE_Failure;
2534 : }
2535 :
2536 50 : ++nBlockCounter;
2537 50 : if (!pfnProgress(static_cast<double>(nBlockCounter) /
2538 50 : static_cast<double>(nTotalBlocks),
2539 : nullptr, pProgressData))
2540 : {
2541 0 : CPLError(CE_Failure, CPLE_UserInterrupt,
2542 : "User terminated CreateCopy()");
2543 0 : return CE_Failure;
2544 : }
2545 : }
2546 : }
2547 :
2548 1 : query.finalize();
2549 :
2550 1 : return CE_None;
2551 : }
2552 0 : catch (const tiledb::TileDBError &e)
2553 : {
2554 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
2555 0 : return CE_Failure;
2556 : }
2557 : }
2558 :
2559 : /************************************************************************/
2560 : /* Create() */
2561 : /************************************************************************/
2562 :
2563 119 : TileDBRasterDataset *TileDBRasterDataset::Create(const char *pszFilename,
2564 : int nXSize, int nYSize,
2565 : int nBandsIn,
2566 : GDALDataType eType,
2567 : char **papszOptions)
2568 :
2569 : {
2570 238 : CPLString osArrayPath = TileDBDataset::VSI_to_tiledb_uri(pszFilename);
2571 :
2572 : std::unique_ptr<TileDBRasterDataset> poDS(TileDBRasterDataset::CreateLL(
2573 238 : osArrayPath, nXSize, nYSize, nBandsIn, eType, papszOptions));
2574 :
2575 119 : if (!poDS)
2576 3 : return nullptr;
2577 :
2578 : const char *pszAttrName =
2579 116 : CPLGetConfigOption("TILEDB_ATTRIBUTE", TILEDB_VALUES);
2580 314 : for (int i = 0; i < poDS->nBands; i++)
2581 : {
2582 198 : if (poDS->eIndexMode == ATTRIBUTES)
2583 24 : poDS->SetBand(
2584 : i + 1, new TileDBRasterBand(
2585 12 : poDS.get(), i + 1,
2586 24 : TILEDB_VALUES + CPLString().Printf("_%i", i + 1)));
2587 : else
2588 558 : poDS->SetBand(i + 1,
2589 372 : new TileDBRasterBand(poDS.get(), i + 1, pszAttrName));
2590 : }
2591 :
2592 : // TILEDB_WRITE_IMAGE_STRUCTURE=NO only used for unit test purposes (to
2593 : // check ability of GDAL to read an arbitrary array)
2594 116 : if (CPLTestBool(CPLGetConfigOption("TILEDB_WRITE_IMAGE_STRUCTURE", "YES")))
2595 : {
2596 204 : CPLStringList aosImageStruct;
2597 : aosImageStruct.SetNameValue(
2598 102 : "NBITS", CPLString().Printf("%d", poDS->nBitsPerSample));
2599 : aosImageStruct.SetNameValue(
2600 : "DATA_TYPE",
2601 102 : CPLString().Printf("%s", GDALGetDataTypeName(poDS->eDataType)));
2602 : aosImageStruct.SetNameValue(
2603 102 : "X_SIZE", CPLString().Printf("%d", poDS->nRasterXSize));
2604 : aosImageStruct.SetNameValue(
2605 102 : "Y_SIZE", CPLString().Printf("%d", poDS->nRasterYSize));
2606 : aosImageStruct.SetNameValue("INTERLEAVE",
2607 102 : index_type_name(poDS->eIndexMode));
2608 102 : aosImageStruct.SetNameValue("DATASET_TYPE", RASTER_DATASET_TYPE);
2609 :
2610 102 : if (poDS->m_lpoAttributeDS.size() > 0)
2611 : {
2612 2 : int i = 0;
2613 6 : for (auto const &poAttrDS : poDS->m_lpoAttributeDS)
2614 : {
2615 4 : aosImageStruct.SetNameValue(
2616 8 : CPLString().Printf("TILEDB_ATTRIBUTE_%i", ++i),
2617 12 : CPLGetBasenameSafe(poAttrDS->GetDescription()).c_str());
2618 : }
2619 : }
2620 102 : poDS->SetMetadata(aosImageStruct.List(), "IMAGE_STRUCTURE");
2621 : }
2622 :
2623 116 : return poDS.release();
2624 : }
2625 :
2626 : /************************************************************************/
2627 : /* CreateCopy() */
2628 : /************************************************************************/
2629 :
2630 52 : GDALDataset *TileDBRasterDataset::CreateCopy(const char *pszFilename,
2631 : GDALDataset *poSrcDS, int bStrict,
2632 : char **papszOptions,
2633 : GDALProgressFunc pfnProgress,
2634 : void *pProgressData)
2635 :
2636 : {
2637 104 : CPLStringList aosOptions(CSLDuplicate(papszOptions));
2638 104 : CPLString osArrayPath = TileDBDataset::VSI_to_tiledb_uri(pszFilename);
2639 :
2640 52 : std::unique_ptr<TileDBRasterDataset> poDstDS;
2641 :
2642 52 : if (CSLFetchNameValue(papszOptions, "APPEND_SUBDATASET"))
2643 : {
2644 : // TileDB schemas are fixed
2645 0 : CPLError(CE_Failure, CPLE_NotSupported,
2646 : "TileDB driver does not support "
2647 : "appending to an existing schema.");
2648 0 : return nullptr;
2649 : }
2650 :
2651 52 : char **papszSrcSubDatasets = poSrcDS->GetMetadata("SUBDATASETS");
2652 :
2653 52 : if (papszSrcSubDatasets == nullptr)
2654 : {
2655 51 : const int nBands = poSrcDS->GetRasterCount();
2656 :
2657 51 : if (nBands > 0)
2658 : {
2659 51 : GDALRasterBand *poBand = poSrcDS->GetRasterBand(1);
2660 51 : GDALDataType eType = poBand->GetRasterDataType();
2661 :
2662 85 : for (int i = 2; i <= nBands; ++i)
2663 : {
2664 34 : if (eType != poSrcDS->GetRasterBand(i)->GetRasterDataType())
2665 : {
2666 0 : CPLError(CE_Failure, CPLE_NotSupported,
2667 : "TileDB driver does not support "
2668 : "source dataset with different band data types.");
2669 0 : return nullptr;
2670 : }
2671 : }
2672 :
2673 51 : poDstDS.reset(TileDBRasterDataset::Create(
2674 : osArrayPath, poSrcDS->GetRasterXSize(),
2675 : poSrcDS->GetRasterYSize(), nBands, eType, papszOptions));
2676 :
2677 51 : if (!poDstDS)
2678 : {
2679 3 : return nullptr;
2680 : }
2681 :
2682 130 : for (int i = 1; i <= nBands; ++i)
2683 : {
2684 82 : int bHasNoData = FALSE;
2685 : const double dfNoData =
2686 82 : poSrcDS->GetRasterBand(i)->GetNoDataValue(&bHasNoData);
2687 82 : if (bHasNoData)
2688 3 : poDstDS->GetRasterBand(i)->SetNoDataValue(dfNoData);
2689 : }
2690 :
2691 : CPLErr eErr =
2692 48 : GDALDatasetCopyWholeRaster(poSrcDS, poDstDS.get(), papszOptions,
2693 : pfnProgress, pProgressData);
2694 :
2695 48 : if (eErr != CE_None)
2696 : {
2697 0 : CPLError(eErr, CPLE_AppDefined,
2698 : "Error copying raster to TileDB.");
2699 0 : return nullptr;
2700 : }
2701 : }
2702 : else
2703 : {
2704 0 : CPLError(CE_Failure, CPLE_NotSupported,
2705 : "TileDB driver does not support "
2706 : "source dataset with zero bands.");
2707 0 : return nullptr;
2708 : }
2709 : }
2710 : else
2711 : {
2712 1 : if (bStrict)
2713 : {
2714 0 : CPLError(CE_Failure, CPLE_NotSupported,
2715 : "TileDB driver does not support copying "
2716 : "subdatasets in strict mode.");
2717 0 : return nullptr;
2718 : }
2719 :
2720 2 : if (CSLFetchNameValue(papszOptions, "BLOCKXSIZE") ||
2721 1 : CSLFetchNameValue(papszOptions, "BLOCKYSIZE"))
2722 : {
2723 0 : CPLError(CE_Failure, CPLE_NotSupported,
2724 : "Changing block size is not supported when copying "
2725 : "subdatasets.");
2726 0 : return nullptr;
2727 : }
2728 :
2729 1 : const int nSubDatasetCount = CSLCount(papszSrcSubDatasets) / 2;
2730 : const int nMaxFiles =
2731 1 : atoi(CPLGetConfigOption("GDAL_READDIR_LIMIT_ON_OPEN", "1000"));
2732 :
2733 1 : aosOptions.SetNameValue("CREATE_GROUP", "NO");
2734 :
2735 1 : if (nSubDatasetCount <= nMaxFiles)
2736 : {
2737 : const char *pszSource =
2738 1 : CSLFetchNameValue(papszSrcSubDatasets, "SUBDATASET_1_NAME");
2739 1 : if (pszSource)
2740 : {
2741 : std::unique_ptr<GDALDataset> poSubDataset(
2742 1 : GDALDataset::Open(pszSource));
2743 1 : if (poSubDataset && poSubDataset->GetRasterCount() > 0)
2744 : {
2745 1 : GDALRasterBand *poBand = poSubDataset->GetRasterBand(1);
2746 :
2747 1 : TileDBRasterDataset::SetBlockSize(poBand, aosOptions);
2748 1 : poDstDS.reset(TileDBRasterDataset::CreateLL(
2749 : osArrayPath, poBand->GetXSize(), poBand->GetYSize(), 0,
2750 1 : poBand->GetRasterDataType(), aosOptions.List()));
2751 :
2752 1 : if (poDstDS)
2753 : {
2754 1 : if (!poDstDS->DeferredCreate(
2755 : /* bCreateArray = */ false))
2756 0 : return nullptr;
2757 :
2758 1 : if (TileDBRasterDataset::CopySubDatasets(
2759 : poSrcDS, poDstDS.get(), pfnProgress,
2760 1 : pProgressData) != CE_None)
2761 : {
2762 0 : poDstDS.reset();
2763 0 : CPLError(CE_Failure, CPLE_AppDefined,
2764 : "Unable to copy subdatasets.");
2765 : }
2766 : }
2767 : }
2768 : }
2769 : }
2770 : else
2771 : {
2772 0 : CPLError(CE_Failure, CPLE_AppDefined,
2773 : "Please increase GDAL_READDIR_LIMIT_ON_OPEN variable.");
2774 : }
2775 : }
2776 :
2777 : // TODO Supporting mask bands is a possible future task
2778 49 : if (poDstDS != nullptr)
2779 : {
2780 49 : int nCloneFlags = GCIF_PAM_DEFAULT & ~GCIF_MASK;
2781 49 : poDstDS->CloneInfo(poSrcDS, nCloneFlags);
2782 :
2783 49 : if (poDstDS->FlushCache(false) != CE_None)
2784 : {
2785 0 : CPLError(CE_Failure, CPLE_AppDefined, "FlushCache() failed");
2786 0 : return nullptr;
2787 : }
2788 :
2789 49 : return poDstDS.release();
2790 : }
2791 0 : return nullptr;
2792 : }
2793 :
2794 : /************************************************************************/
2795 : /* LoadOverviews() */
2796 : /************************************************************************/
2797 :
2798 28 : void TileDBRasterDataset::LoadOverviews()
2799 : {
2800 28 : if (m_bLoadOverviewsDone)
2801 23 : return;
2802 5 : m_bLoadOverviewsDone = true;
2803 :
2804 : // NOTE: read overview_model.rst for a high level explanation of overviews
2805 : // are stored.
2806 :
2807 5 : if (!m_bDatasetInGroup)
2808 0 : return;
2809 :
2810 5 : CPLStringList aosOpenOptions;
2811 5 : if (nTimestamp)
2812 : {
2813 : aosOpenOptions.SetNameValue("TILEDB_TIMESTAMP",
2814 0 : CPLSPrintf("%" PRIu64, nTimestamp));
2815 : }
2816 5 : if (!m_osConfigFilename.empty())
2817 : {
2818 : aosOpenOptions.SetNameValue("TILEDB_CONFIG",
2819 0 : m_osConfigFilename.c_str());
2820 : }
2821 10 : for (int i = 0; i < m_nOverviewCountFromMetadata; ++i)
2822 : {
2823 5 : const std::string osArrayName = CPLSPrintf("l_%d", 1 + i);
2824 : const std::string osOvrDatasetName =
2825 5 : CPLFormFilenameSafe(GetDescription(), osArrayName.c_str(), nullptr);
2826 :
2827 5 : GDALOpenInfo oOpenInfo(osOvrDatasetName.c_str(), eAccess);
2828 5 : oOpenInfo.papszOpenOptions = aosOpenOptions.List();
2829 : auto poOvrDS = std::unique_ptr<GDALDataset>(
2830 5 : Open(&oOpenInfo, tiledb::Object::Type::Array));
2831 5 : if (!poOvrDS)
2832 0 : return;
2833 5 : if (poOvrDS->GetRasterCount() != nBands)
2834 : {
2835 0 : CPLError(CE_Failure, CPLE_AppDefined,
2836 : "Overview %s has not the same number of bands as full "
2837 : "resolution dataset",
2838 : osOvrDatasetName.c_str());
2839 0 : return;
2840 : }
2841 5 : m_apoOverviewDS.emplace_back(std::move(poOvrDS));
2842 : }
2843 : }
2844 :
2845 : /************************************************************************/
2846 : /* IBuildOverviews() */
2847 : /************************************************************************/
2848 :
2849 13 : CPLErr TileDBRasterDataset::IBuildOverviews(
2850 : const char *pszResampling, int nOverviews, const int *panOverviewList,
2851 : int nListBands, const int *panBandList, GDALProgressFunc pfnProgress,
2852 : void *pProgressData, CSLConstList papszOptions)
2853 : {
2854 : // NOTE: read overview_model.rst for a high level explanation of overviews
2855 : // are stored.
2856 :
2857 13 : if (eAccess == GA_ReadOnly)
2858 : {
2859 3 : if (!CPLTestBool(
2860 : CPLGetConfigOption("TILEDB_GEOTIFF_OVERVIEWS", "FALSE")))
2861 : {
2862 2 : ReportError(
2863 : CE_Failure, CPLE_NotSupported,
2864 : "Cannot %s overviews in TileDB format in read-only mode",
2865 : nOverviews ? "create" : "delete");
2866 2 : return CE_Failure;
2867 : }
2868 :
2869 : // GeoTIFF overviews. This used to be supported before GDAL 3.10
2870 : // although likely not desirable.
2871 1 : return GDALPamDataset::IBuildOverviews(
2872 : pszResampling, nOverviews, panOverviewList, nListBands, panBandList,
2873 1 : pfnProgress, pProgressData, papszOptions);
2874 : }
2875 :
2876 10 : if (nBands == 0)
2877 : {
2878 0 : return CE_Failure;
2879 : }
2880 :
2881 : // If we already have PAM overview (i.e. GeoTIFF based), go through PAM
2882 10 : if (cpl::down_cast<GDALPamRasterBand *>(GetRasterBand(1))
2883 10 : ->GDALPamRasterBand::GetOverviewCount() > 0)
2884 : {
2885 1 : return GDALPamDataset::IBuildOverviews(
2886 : pszResampling, nOverviews, panOverviewList, nListBands, panBandList,
2887 1 : pfnProgress, pProgressData, papszOptions);
2888 : }
2889 :
2890 9 : if (!m_bDatasetInGroup)
2891 : {
2892 2 : ReportError(CE_Failure, CPLE_NotSupported,
2893 : "IBuildOverviews() only supported for datasets created "
2894 : "with CREATE_GROUP=YES");
2895 2 : return CE_Failure;
2896 : }
2897 :
2898 : /* -------------------------------------------------------------------- */
2899 : /* Our overview support currently only works safely if all */
2900 : /* bands are handled at the same time. */
2901 : /* -------------------------------------------------------------------- */
2902 7 : if (nListBands != nBands)
2903 : {
2904 0 : ReportError(CE_Failure, CPLE_NotSupported,
2905 : "Generation of TileDB overviews currently only "
2906 : "supported when operating on all bands. "
2907 : "Operation failed.");
2908 0 : return CE_Failure;
2909 : }
2910 :
2911 : // Force loading existing overviews
2912 7 : if (m_nOverviewCountFromMetadata)
2913 4 : GetRasterBand(1)->GetOverviewCount();
2914 7 : m_bLoadOverviewsDone = true;
2915 :
2916 : /* -------------------------------------------------------------------- */
2917 : /* Deletes existing overviews if requested. */
2918 : /* -------------------------------------------------------------------- */
2919 7 : if (nOverviews == 0)
2920 : {
2921 2 : CPLErr eErr = CE_None;
2922 :
2923 : // Unlink arrays from he group
2924 : try
2925 : {
2926 6 : tiledb::Group group(*m_ctx, GetDescription(), TILEDB_WRITE);
2927 4 : for (auto &&poODS : m_apoOverviewDS)
2928 : {
2929 2 : group.remove_member(poODS->GetDescription());
2930 : }
2931 2 : group.close();
2932 : }
2933 0 : catch (const tiledb::TileDBError &e)
2934 : {
2935 0 : eErr = CE_Failure;
2936 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
2937 : }
2938 :
2939 2 : tiledb::VFS vfs(*m_ctx, m_ctx->config());
2940 :
2941 : // Delete arrays
2942 4 : for (auto &&poODS : m_apoOverviewDS)
2943 : {
2944 : try
2945 : {
2946 2 : CPL_IGNORE_RET_VAL(poODS->Close());
2947 2 : tiledb::Array::delete_array(*m_ctx, poODS->GetDescription());
2948 2 : if (vfs.is_dir(poODS->GetDescription()))
2949 : {
2950 2 : vfs.remove_dir(poODS->GetDescription());
2951 : }
2952 : }
2953 0 : catch (const tiledb::TileDBError &e)
2954 : {
2955 0 : eErr = CE_Failure;
2956 0 : CPLError(CE_Failure, CPLE_AppDefined,
2957 : "Array::delete_array(%s) failed: %s",
2958 0 : poODS->GetDescription(), e.what());
2959 : }
2960 2 : m_apoOverviewDSRemoved.emplace_back(std::move(poODS));
2961 : }
2962 :
2963 2 : m_apoOverviewDS.clear();
2964 2 : m_nOverviewCountFromMetadata = 0;
2965 2 : MarkPamDirty();
2966 2 : return eErr;
2967 : }
2968 :
2969 : /* -------------------------------------------------------------------- */
2970 : /* Establish which of the overview levels we already have, and */
2971 : /* which are new. */
2972 : /* -------------------------------------------------------------------- */
2973 10 : std::vector<bool> abRequireNewOverview(nOverviews, true);
2974 10 : for (int i = 0; i < nOverviews; ++i)
2975 : {
2976 5 : const int nOXSize = DIV_ROUND_UP(GetRasterXSize(), panOverviewList[i]);
2977 5 : const int nOYSize = DIV_ROUND_UP(GetRasterYSize(), panOverviewList[i]);
2978 :
2979 6 : for (const auto &poODS : m_apoOverviewDS)
2980 : {
2981 : const int nOvFactor =
2982 2 : GDALComputeOvFactor(poODS->GetRasterXSize(), GetRasterXSize(),
2983 : poODS->GetRasterYSize(), GetRasterYSize());
2984 :
2985 : // If we already have a 1x1 overview and this new one would result
2986 : // in it too, then don't create it.
2987 2 : if (poODS->GetRasterXSize() == 1 && poODS->GetRasterYSize() == 1 &&
2988 2 : nOXSize == 1 && nOYSize == 1)
2989 : {
2990 0 : abRequireNewOverview[i] = false;
2991 0 : break;
2992 : }
2993 :
2994 3 : if (nOvFactor == panOverviewList[i] ||
2995 1 : nOvFactor == GDALOvLevelAdjust2(panOverviewList[i],
2996 : GetRasterXSize(),
2997 : GetRasterYSize()))
2998 : {
2999 1 : abRequireNewOverview[i] = false;
3000 1 : break;
3001 : }
3002 : }
3003 :
3004 5 : if (abRequireNewOverview[i])
3005 : {
3006 4 : CPLStringList aosCreationOptions;
3007 4 : aosCreationOptions.SetNameValue("CREATE_GROUP", "NO");
3008 : aosCreationOptions.SetNameValue(
3009 4 : "NBITS", CPLString().Printf("%d", nBitsPerSample));
3010 : aosCreationOptions.SetNameValue("INTERLEAVE",
3011 4 : index_type_name(eIndexMode));
3012 4 : if (nTimestamp)
3013 : {
3014 : aosCreationOptions.SetNameValue(
3015 0 : "TILEDB_TIMESTAMP", CPLSPrintf("%" PRIu64, nTimestamp));
3016 : }
3017 4 : if (!m_osConfigFilename.empty())
3018 : {
3019 : aosCreationOptions.SetNameValue("TILEDB_CONFIG",
3020 0 : m_osConfigFilename.c_str());
3021 : }
3022 :
3023 : const std::string osArrayName =
3024 4 : CPLSPrintf("l_%d", 1 + int(m_apoOverviewDS.size()));
3025 : const std::string osOvrDatasetName = CPLFormFilenameSafe(
3026 4 : GetDescription(), osArrayName.c_str(), nullptr);
3027 :
3028 : auto poOvrDS = std::unique_ptr<TileDBRasterDataset>(
3029 : Create(osOvrDatasetName.c_str(), nOXSize, nOYSize, nBands,
3030 : GetRasterBand(1)->GetRasterDataType(),
3031 4 : aosCreationOptions.List()));
3032 4 : if (!poOvrDS)
3033 0 : return CE_Failure;
3034 :
3035 : // Apply nodata from main dataset
3036 16 : for (int j = 0; j < nBands; ++j)
3037 : {
3038 12 : int bHasNoData = FALSE;
3039 : const double dfNoData =
3040 12 : GetRasterBand(j + 1)->GetNoDataValue(&bHasNoData);
3041 12 : if (bHasNoData)
3042 12 : poOvrDS->GetRasterBand(j + 1)->SetNoDataValue(dfNoData);
3043 : }
3044 :
3045 : // Apply georeferencing from main dataset
3046 4 : poOvrDS->SetSpatialRef(GetSpatialRef());
3047 : double adfGeoTransform[6];
3048 4 : if (GetGeoTransform(adfGeoTransform) == CE_None)
3049 : {
3050 4 : adfGeoTransform[1] *=
3051 4 : static_cast<double>(nRasterXSize) / nOXSize;
3052 4 : adfGeoTransform[2] *=
3053 4 : static_cast<double>(nRasterXSize) / nOXSize;
3054 4 : adfGeoTransform[4] *=
3055 4 : static_cast<double>(nRasterYSize) / nOYSize;
3056 4 : adfGeoTransform[5] *=
3057 4 : static_cast<double>(nRasterYSize) / nOYSize;
3058 4 : poOvrDS->SetGeoTransform(adfGeoTransform);
3059 : }
3060 :
3061 4 : poOvrDS->DeferredCreate(/* bCreateArray = */ true);
3062 :
3063 : try
3064 : {
3065 12 : tiledb::Group group(*m_ctx, GetDescription(), TILEDB_WRITE);
3066 4 : group.add_member(osOvrDatasetName, false);
3067 4 : group.close();
3068 : }
3069 0 : catch (const tiledb::TileDBError &e)
3070 : {
3071 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
3072 : }
3073 :
3074 4 : m_apoOverviewDS.emplace_back(std::move(poOvrDS));
3075 : }
3076 : }
3077 :
3078 5 : m_nOverviewCountFromMetadata = static_cast<int>(m_apoOverviewDS.size());
3079 5 : MarkPamDirty();
3080 :
3081 : /* -------------------------------------------------------------------- */
3082 : /* Refresh/generate overviews that are listed. */
3083 : /* -------------------------------------------------------------------- */
3084 10 : std::vector<GDALRasterBand *> apoSrcBands;
3085 10 : std::vector<std::vector<GDALRasterBand *>> aapoOverviewBands;
3086 5 : CPLErr eErr = CE_None;
3087 : const auto osNormalizedResampling =
3088 5 : GDALGetNormalizedOvrResampling(pszResampling);
3089 20 : for (int iBand = 0; eErr == CE_None && iBand < nBands; iBand++)
3090 : {
3091 15 : apoSrcBands.push_back(GetRasterBand(iBand + 1));
3092 30 : std::vector<GDALRasterBand *> apoOverviewBands;
3093 :
3094 : std::vector<bool> abAlreadyUsedOverviewBand(m_apoOverviewDS.size(),
3095 30 : false);
3096 :
3097 30 : for (int i = 0; i < nOverviews; i++)
3098 : {
3099 15 : bool bFound = false;
3100 18 : for (size_t j = 0; j < m_apoOverviewDS.size(); ++j)
3101 : {
3102 18 : if (!abAlreadyUsedOverviewBand[j])
3103 : {
3104 18 : auto &poODS = m_apoOverviewDS[j];
3105 18 : int nOvFactor = GDALComputeOvFactor(
3106 : poODS->GetRasterXSize(), nRasterXSize,
3107 : poODS->GetRasterYSize(), nRasterYSize);
3108 :
3109 21 : if (nOvFactor == panOverviewList[i] ||
3110 3 : nOvFactor == GDALOvLevelAdjust2(panOverviewList[i],
3111 : nRasterXSize,
3112 : nRasterYSize))
3113 : {
3114 15 : abAlreadyUsedOverviewBand[j] = true;
3115 15 : auto poOvrBand = poODS->GetRasterBand(iBand + 1);
3116 15 : if (!osNormalizedResampling.empty())
3117 : {
3118 : // Store resampling method in band metadata, as it
3119 : // can be used by the gdaladdo utilities to refresh
3120 : // existing overviews with the method previously
3121 : // used
3122 15 : poOvrBand->SetMetadataItem(
3123 15 : "RESAMPLING", osNormalizedResampling.c_str());
3124 : }
3125 15 : apoOverviewBands.push_back(poOvrBand);
3126 15 : bFound = true;
3127 15 : break;
3128 : }
3129 : }
3130 : }
3131 15 : if (!bFound)
3132 : {
3133 0 : CPLError(CE_Failure, CPLE_AppDefined,
3134 : "Could not find dataset corresponding to ov factor %d",
3135 0 : panOverviewList[i]);
3136 0 : eErr = CE_Failure;
3137 : }
3138 : }
3139 15 : if (iBand > 0)
3140 : {
3141 10 : CPLAssert(apoOverviewBands.size() == aapoOverviewBands[0].size());
3142 : }
3143 15 : aapoOverviewBands.emplace_back(std::move(apoOverviewBands));
3144 : }
3145 :
3146 5 : if (eErr == CE_None)
3147 : {
3148 5 : eErr = GDALRegenerateOverviewsMultiBand(apoSrcBands, aapoOverviewBands,
3149 : pszResampling, pfnProgress,
3150 : pProgressData, papszOptions);
3151 : }
3152 :
3153 5 : return eErr;
3154 : }
|