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