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