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 : * Permission is hereby granted, free of charge, to any person obtaining a
11 : * copy of this software and associated documentation files (the "Software"),
12 : * to deal in the Software without restriction, including without limitation
13 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
14 : * and/or sell copies of the Software, and to permit persons to whom the
15 : * Software is furnished to do so, subject to the following conditions:
16 : *
17 : * The above copyright notice and this permission notice shall be included
18 : * in all copies or substantial portions of the Software.
19 : *
20 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
21 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
22 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
23 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
24 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
25 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
26 : * DEALINGS IN THE SOFTWARE.
27 : ****************************************************************************/
28 :
29 : #include <cassert>
30 : #include <limits>
31 :
32 : #include "tiledbheaders.h"
33 :
34 : constexpr const char *RASTER_DATASET_TYPE = "raster";
35 :
36 : /************************************************************************/
37 : /* ==================================================================== */
38 : /* TileDBRasterBand */
39 : /* ==================================================================== */
40 : /************************************************************************/
41 :
42 : class TileDBRasterBand final : public GDALPamRasterBand
43 : {
44 : friend class TileDBRasterDataset;
45 :
46 : protected:
47 : TileDBRasterDataset *poGDS;
48 : bool bStats;
49 : CPLString osAttrName;
50 :
51 : public:
52 : TileDBRasterBand(TileDBRasterDataset *, int,
53 : const std::string &osAttr = TILEDB_VALUES);
54 : virtual CPLErr IReadBlock(int, int, void *) override;
55 : virtual CPLErr IWriteBlock(int, int, void *) override;
56 : virtual CPLErr IRasterIO(GDALRWFlag, int, int, int, int, void *, int, int,
57 : GDALDataType, GSpacing, GSpacing,
58 : GDALRasterIOExtraArg *psExtraArg) override;
59 : virtual GDALColorInterp GetColorInterpretation() override;
60 : };
61 :
62 153 : static CPLErr option_to_index_type(const char *pszIndexingType,
63 : TILEDB_INTERLEAVE_MODE &eMode)
64 : {
65 153 : CPLErr eErr = CE_None;
66 :
67 153 : if (pszIndexingType)
68 : {
69 80 : if EQUAL (pszIndexingType, "BAND")
70 53 : eMode = BAND;
71 27 : else if EQUAL (pszIndexingType, "ATTRIBUTES")
72 6 : eMode = ATTRIBUTES;
73 21 : else if EQUAL (pszIndexingType, "PIXEL")
74 21 : eMode = PIXEL;
75 : else
76 : {
77 0 : eErr = CE_Failure;
78 0 : CPLError(eErr, CPLE_AppDefined,
79 : "Unable to identify TileDB index mode %s.",
80 : pszIndexingType);
81 : }
82 : }
83 : else
84 : {
85 73 : eMode = BAND;
86 : }
87 :
88 153 : return eErr;
89 : }
90 :
91 65 : static const char *index_type_name(TILEDB_INTERLEAVE_MODE eMode)
92 : {
93 65 : switch (eMode)
94 : {
95 7 : case PIXEL:
96 7 : return "PIXEL";
97 4 : case ATTRIBUTES:
98 4 : return "ATTRIBUTES";
99 54 : case BAND:
100 54 : return "BAND";
101 0 : default:
102 0 : return nullptr;
103 : }
104 : }
105 :
106 : /************************************************************************/
107 : /* SetBuffer() */
108 : /************************************************************************/
109 :
110 332 : static CPLErr SetBuffer(tiledb::Query *poQuery, GDALDataType eType,
111 : const CPLString &osAttrName, void *pImage, int nSize)
112 : {
113 332 : switch (eType)
114 : {
115 136 : case GDT_Byte:
116 : poQuery->set_data_buffer(
117 136 : osAttrName, reinterpret_cast<unsigned char *>(pImage), nSize);
118 136 : break;
119 2 : case GDT_Int8:
120 : poQuery->set_data_buffer(osAttrName,
121 2 : reinterpret_cast<int8_t *>(pImage), nSize);
122 2 : break;
123 3 : case GDT_UInt16:
124 : poQuery->set_data_buffer(
125 3 : osAttrName, reinterpret_cast<unsigned short *>(pImage), nSize);
126 3 : break;
127 3 : case GDT_UInt32:
128 : poQuery->set_data_buffer(
129 3 : osAttrName, reinterpret_cast<unsigned int *>(pImage), nSize);
130 3 : break;
131 2 : case GDT_UInt64:
132 : poQuery->set_data_buffer(
133 2 : osAttrName, reinterpret_cast<uint64_t *>(pImage), nSize);
134 2 : break;
135 3 : case GDT_Int16:
136 : poQuery->set_data_buffer(osAttrName,
137 3 : reinterpret_cast<short *>(pImage), nSize);
138 3 : break;
139 6 : case GDT_Int32:
140 : poQuery->set_data_buffer(osAttrName,
141 6 : reinterpret_cast<int *>(pImage), nSize);
142 6 : break;
143 2 : case GDT_Int64:
144 : poQuery->set_data_buffer(
145 2 : osAttrName, reinterpret_cast<int64_t *>(pImage), nSize);
146 2 : break;
147 156 : case GDT_Float32:
148 : poQuery->set_data_buffer(osAttrName,
149 156 : reinterpret_cast<float *>(pImage), nSize);
150 156 : break;
151 3 : case GDT_Float64:
152 : poQuery->set_data_buffer(osAttrName,
153 3 : reinterpret_cast<double *>(pImage), nSize);
154 3 : break;
155 3 : case GDT_CInt16:
156 : poQuery->set_data_buffer(
157 3 : osAttrName, reinterpret_cast<short *>(pImage), nSize * 2);
158 3 : break;
159 3 : case GDT_CInt32:
160 : poQuery->set_data_buffer(
161 3 : osAttrName, reinterpret_cast<int *>(pImage), nSize * 2);
162 3 : break;
163 3 : case GDT_CFloat32:
164 : poQuery->set_data_buffer(
165 3 : osAttrName, reinterpret_cast<float *>(pImage), nSize * 2);
166 3 : break;
167 7 : case GDT_CFloat64:
168 : poQuery->set_data_buffer(
169 7 : osAttrName, reinterpret_cast<double *>(pImage), nSize * 2);
170 7 : break;
171 0 : default:
172 0 : return CE_Failure;
173 : }
174 332 : return CE_None;
175 : }
176 :
177 : /************************************************************************/
178 : /* TileDBRasterBand() */
179 : /************************************************************************/
180 :
181 792 : TileDBRasterBand::TileDBRasterBand(TileDBRasterDataset *poDSIn, int nBandIn,
182 792 : const std::string &osAttr)
183 792 : : poGDS(poDSIn), bStats(poDSIn->bStats), osAttrName(osAttr)
184 : {
185 792 : poDS = poDSIn;
186 792 : nBand = nBandIn;
187 792 : eDataType = poGDS->eDataType;
188 792 : if (eDataType == GDT_Unknown)
189 : {
190 : try
191 : {
192 14 : auto attr = (poGDS->m_roArray ? poGDS->m_roArray : poGDS->m_array)
193 14 : ->schema()
194 28 : .attribute(osAttr);
195 14 : switch (attr.type())
196 : {
197 1 : case TILEDB_INT8:
198 1 : eDataType = GDT_Int8;
199 1 : break;
200 1 : case TILEDB_UINT8:
201 1 : eDataType = GDT_Byte;
202 1 : break;
203 2 : case TILEDB_INT16:
204 2 : eDataType =
205 2 : attr.cell_val_num() == 2 ? GDT_CInt16 : GDT_Int16;
206 2 : break;
207 1 : case TILEDB_UINT16:
208 1 : eDataType = GDT_UInt16;
209 1 : break;
210 2 : case TILEDB_INT32:
211 2 : eDataType =
212 2 : attr.cell_val_num() == 2 ? GDT_CInt32 : GDT_Int32;
213 2 : break;
214 1 : case TILEDB_UINT32:
215 1 : eDataType = GDT_UInt32;
216 1 : break;
217 1 : case TILEDB_INT64:
218 1 : eDataType = GDT_Int64;
219 1 : break;
220 1 : case TILEDB_UINT64:
221 1 : eDataType = GDT_UInt64;
222 1 : break;
223 2 : case TILEDB_FLOAT32:
224 2 : eDataType =
225 2 : attr.cell_val_num() == 2 ? GDT_CFloat32 : GDT_Float32;
226 2 : break;
227 2 : case TILEDB_FLOAT64:
228 2 : eDataType =
229 2 : attr.cell_val_num() == 2 ? GDT_CFloat64 : GDT_Float64;
230 2 : break;
231 0 : default:
232 : {
233 0 : const char *pszTypeName = "";
234 0 : tiledb_datatype_to_str(attr.type(), &pszTypeName);
235 0 : CPLError(CE_Failure, CPLE_NotSupported,
236 : "Unhandled TileDB data type: %s", pszTypeName);
237 0 : break;
238 : }
239 : }
240 : }
241 0 : catch (const tiledb::TileDBError &e)
242 : {
243 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
244 : }
245 : }
246 792 : eAccess = poGDS->eAccess;
247 792 : nRasterXSize = poGDS->nRasterXSize;
248 792 : nRasterYSize = poGDS->nRasterYSize;
249 792 : nBlockXSize = poGDS->nBlockXSize;
250 792 : nBlockYSize = poGDS->nBlockYSize;
251 792 : }
252 :
253 : /************************************************************************/
254 : /* IRasterIO() */
255 : /************************************************************************/
256 :
257 264 : CPLErr TileDBRasterBand::IRasterIO(GDALRWFlag eRWFlag, int nXOff, int nYOff,
258 : int nXSize, int nYSize, void *pData,
259 : int nBufXSize, int nBufYSize,
260 : GDALDataType eBufType, GSpacing nPixelSpace,
261 : GSpacing nLineSpace,
262 : GDALRasterIOExtraArg *psExtraArg)
263 : {
264 264 : if (poGDS->eIndexMode == ATTRIBUTES && eRWFlag == GF_Write)
265 : {
266 0 : CPLError(CE_Failure, CPLE_NoWriteAccess,
267 : "Unable to write using band ordered IRasterIO when using "
268 : "interleave 'ATTRIBUTES'.\n");
269 0 : return CE_Failure;
270 : }
271 :
272 264 : const int nBufferDTSize(GDALGetDataTypeSizeBytes(eBufType));
273 :
274 264 : if (eBufType == eDataType && nXSize == nBufXSize && nYSize == nBufYSize &&
275 216 : nBufferDTSize > 0 && (nPixelSpace % nBufferDTSize) == 0 &&
276 216 : (nLineSpace % nBufferDTSize) == 0)
277 : {
278 216 : const uint64_t nBandIdx = poGDS->nBandStart + nBand - 1;
279 : std::vector<uint64_t> oaSubarray = {
280 : nBandIdx, nBandIdx,
281 216 : (uint64_t)nYOff, (uint64_t)nYOff + nYSize - 1,
282 432 : (uint64_t)nXOff, (uint64_t)nXOff + nXSize - 1};
283 216 : if (poGDS->eIndexMode == PIXEL)
284 78 : std::rotate(oaSubarray.begin(), oaSubarray.begin() + 2,
285 78 : oaSubarray.end());
286 :
287 : const bool bUseReadOnlyObjs =
288 270 : ((eRWFlag == GF_Read) && (eAccess == GA_Update) &&
289 54 : (poGDS->m_roArray));
290 216 : const auto &oCtxt = bUseReadOnlyObjs ? *poGDS->m_roCtx : *poGDS->m_ctx;
291 : const auto &oArray =
292 216 : bUseReadOnlyObjs ? *poGDS->m_roArray : *poGDS->m_array;
293 :
294 432 : auto poQuery = std::make_unique<tiledb::Query>(oCtxt, oArray);
295 432 : tiledb::Subarray subarray(oCtxt, oArray);
296 216 : if (poGDS->m_array->schema().domain().ndim() == 3)
297 : {
298 149 : subarray.set_subarray(oaSubarray);
299 : }
300 : else
301 : {
302 134 : subarray.set_subarray(std::vector<uint64_t>(oaSubarray.cbegin() + 2,
303 67 : oaSubarray.cend()));
304 : }
305 216 : poQuery->set_subarray(subarray);
306 :
307 216 : SetBuffer(poQuery.get(), eDataType, osAttrName, pData, nXSize * nYSize);
308 :
309 : // write additional co-registered values
310 432 : std::vector<std::unique_ptr<void, decltype(&VSIFree)>> aBlocks;
311 :
312 216 : if (poGDS->lpoAttributeDS.size() > 0)
313 : {
314 10 : for (auto const &poAttrDS : poGDS->lpoAttributeDS)
315 : {
316 7 : GDALRasterBand *poAttrBand = poAttrDS->GetRasterBand(nBand);
317 7 : GDALDataType eAttrType = poAttrBand->GetRasterDataType();
318 7 : int nBytes = GDALGetDataTypeSizeBytes(eAttrType);
319 7 : size_t nValues = static_cast<size_t>(nBufXSize) * nBufYSize;
320 7 : void *pAttrBlock = VSIMalloc(nBytes * nValues);
321 7 : aBlocks.emplace_back(pAttrBlock, &VSIFree);
322 :
323 7 : if (pAttrBlock == nullptr)
324 : {
325 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
326 : "Cannot allocate attribute buffer");
327 1 : return CE_Failure;
328 : }
329 :
330 7 : poAttrBand->AdviseRead(nXOff, nYOff, nXSize, nYSize, nBufXSize,
331 7 : nBufYSize, eAttrType, nullptr);
332 :
333 7 : CPLErr eErr = poAttrBand->RasterIO(
334 : GF_Read, nXOff, nYOff, nXSize, nYSize, pAttrBlock,
335 : nBufXSize, nBufYSize, eAttrType, nPixelSpace, nLineSpace,
336 : psExtraArg);
337 :
338 7 : if (eErr == CE_None)
339 : {
340 6 : CPLString osName = CPLString().Printf(
341 12 : "%s", CPLGetBasename(poAttrDS->GetDescription()));
342 :
343 6 : SetBuffer(poQuery.get(), eAttrType, osName, pAttrBlock,
344 : nBufXSize * nBufYSize);
345 : }
346 : else
347 : {
348 1 : return eErr;
349 : }
350 : }
351 : }
352 :
353 215 : if (bStats)
354 0 : tiledb::Stats::enable();
355 :
356 215 : auto status = poQuery->submit();
357 :
358 215 : if (bStats)
359 : {
360 0 : tiledb::Stats::dump(stdout);
361 0 : tiledb::Stats::disable();
362 : }
363 :
364 215 : if (status == tiledb::Query::Status::FAILED)
365 0 : return CE_Failure;
366 : else
367 215 : return CE_None;
368 : }
369 :
370 48 : return GDALPamRasterBand::IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize,
371 : pData, nBufXSize, nBufYSize, eBufType,
372 48 : nPixelSpace, nLineSpace, psExtraArg);
373 : }
374 :
375 124 : CPLErr TileDBRasterBand::IReadBlock(int nBlockXOff, int nBlockYOff,
376 : void *pImage)
377 : {
378 124 : const int nXOff = nBlockXOff * nBlockXSize;
379 124 : const int nYOff = nBlockYOff * nBlockYSize;
380 124 : const int nDTSize = GDALGetDataTypeSizeBytes(eDataType);
381 248 : return IRasterIO(GF_Read, nXOff, nYOff, nBlockXSize, nBlockYSize, pImage,
382 : nBlockXSize, nBlockYSize, eDataType, nDTSize,
383 124 : nDTSize * nBlockXSize, nullptr);
384 : }
385 :
386 : /************************************************************************/
387 : /* IWriteBlock() */
388 : /************************************************************************/
389 :
390 21 : CPLErr TileDBRasterBand::IWriteBlock(int nBlockXOff, int nBlockYOff,
391 : void *pImage)
392 :
393 : {
394 21 : if (eAccess == GA_ReadOnly)
395 : {
396 0 : CPLError(CE_Failure, CPLE_NoWriteAccess,
397 : "Unable to write block, dataset is opened read only.\n");
398 0 : return CE_Failure;
399 : }
400 :
401 21 : CPLAssert(poGDS != nullptr && nBlockXOff >= 0 && nBlockYOff >= 0 &&
402 : pImage != nullptr);
403 :
404 21 : int nStartX = nBlockXSize * nBlockXOff;
405 21 : int nStartY = nBlockYSize * nBlockYOff;
406 :
407 21 : const int nDTSize = GDALGetDataTypeSizeBytes(eDataType);
408 42 : return IRasterIO(GF_Write, nStartX, nStartY, nBlockXSize, nBlockYSize,
409 : pImage, nBlockXSize, nBlockYSize, eDataType, nDTSize,
410 21 : nDTSize * nBlockXSize, nullptr);
411 : }
412 :
413 : /************************************************************************/
414 : /* GetColorInterpretation() */
415 : /************************************************************************/
416 :
417 49 : GDALColorInterp TileDBRasterBand::GetColorInterpretation()
418 :
419 : {
420 49 : if (poGDS->nBands == 1)
421 19 : return GCI_GrayIndex;
422 :
423 30 : if (nBand == 1)
424 10 : return GCI_RedBand;
425 :
426 20 : else if (nBand == 2)
427 10 : return GCI_GreenBand;
428 :
429 10 : else if (nBand == 3)
430 10 : return GCI_BlueBand;
431 :
432 0 : return GCI_AlphaBand;
433 : }
434 :
435 : /************************************************************************/
436 : /* ==================================================================== */
437 : /* TileDBRasterDataset */
438 : /* ==================================================================== */
439 : /************************************************************************/
440 :
441 : /************************************************************************/
442 : /* ~TileDBRasterDataset() */
443 : /************************************************************************/
444 :
445 350 : TileDBRasterDataset::~TileDBRasterDataset()
446 :
447 : {
448 175 : TileDBRasterDataset::FlushCache(true);
449 :
450 : try
451 : {
452 175 : if (m_array)
453 : {
454 162 : m_array->close();
455 : }
456 : }
457 0 : catch (const tiledb::TileDBError &e)
458 : {
459 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
460 : }
461 175 : CPLDestroyXMLNode(psSubDatasetsTree);
462 175 : CSLDestroy(papszSubDatasets);
463 175 : CSLDestroy(papszAttributes);
464 350 : }
465 :
466 : /************************************************************************/
467 : /* IRasterio() */
468 : /************************************************************************/
469 :
470 76 : CPLErr TileDBRasterDataset::IRasterIO(
471 : GDALRWFlag eRWFlag, int nXOff, int nYOff, int nXSize, int nYSize,
472 : void *pData, int nBufXSize, int nBufYSize, GDALDataType eBufType,
473 : int nBandCount, int *panBandMap, GSpacing nPixelSpace, GSpacing nLineSpace,
474 : GSpacing nBandSpace, CPL_UNUSED GDALRasterIOExtraArg *psExtraArg)
475 :
476 : {
477 : // support special case of writing attributes for bands, all attributes have to be set at once
478 76 : const int nBufferDTSize(GDALGetDataTypeSizeBytes(eBufType));
479 :
480 76 : if (eIndexMode == ATTRIBUTES && nBandCount == nBands &&
481 4 : eBufType == eDataType && nXSize == nBufXSize && nYSize == nBufYSize &&
482 4 : nBufferDTSize > 0 && (nPixelSpace % nBufferDTSize) == 0 &&
483 4 : (nLineSpace % nBufferDTSize) == 0)
484 : {
485 : std::vector<uint64_t> oaSubarray = {
486 4 : (uint64_t)nYOff, (uint64_t)nYOff + nYSize - 1, (uint64_t)nXOff,
487 8 : (uint64_t)nXOff + nXSize - 1};
488 :
489 : const bool bUseReadOnlyObjs =
490 4 : ((eRWFlag == GF_Read) && (eAccess == GA_Update) && (m_roArray));
491 4 : const auto &oCtxt = bUseReadOnlyObjs ? *m_roCtx : *m_ctx;
492 4 : const auto &oArray = bUseReadOnlyObjs ? *m_roArray : *m_array;
493 :
494 8 : auto poQuery = std::make_unique<tiledb::Query>(oCtxt, oArray);
495 8 : tiledb::Subarray subarray(oCtxt, oArray);
496 4 : subarray.set_subarray(oaSubarray);
497 4 : poQuery->set_subarray(subarray);
498 :
499 14 : for (int b = 0; b < nBandCount; b++)
500 : {
501 : TileDBRasterBand *poBand =
502 10 : (TileDBRasterBand *)GetRasterBand(panBandMap[b]);
503 10 : int nRegionSize = nBufXSize * nBufYSize * nBufferDTSize;
504 10 : SetBuffer(poQuery.get(), eDataType, poBand->osAttrName,
505 10 : ((GByte *)pData) + b * nRegionSize, nRegionSize);
506 : }
507 :
508 4 : if (bStats)
509 0 : tiledb::Stats::enable();
510 :
511 4 : auto status = poQuery->submit();
512 :
513 4 : if (bStats)
514 : {
515 0 : tiledb::Stats::dump(stdout);
516 0 : tiledb::Stats::disable();
517 : }
518 :
519 4 : if (status == tiledb::Query::Status::FAILED)
520 0 : return CE_Failure;
521 : else
522 4 : return CE_None;
523 : }
524 :
525 72 : return GDALPamDataset::IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize,
526 : pData, nBufXSize, nBufYSize, eBufType,
527 : nBandCount, panBandMap, nPixelSpace,
528 72 : nLineSpace, nBandSpace, psExtraArg);
529 : }
530 :
531 : /************************************************************************/
532 : /* AddDimensions() */
533 : /************************************************************************/
534 :
535 93 : CPLErr TileDBRasterDataset::AddDimensions(tiledb::Domain &domain,
536 : const char *pszAttrName,
537 : tiledb::Dimension &y,
538 : tiledb::Dimension &x,
539 : tiledb::Dimension *poBands)
540 :
541 : {
542 93 : CPLErr eErr = CE_None;
543 93 : switch (eIndexMode)
544 : {
545 5 : case ATTRIBUTES:
546 5 : domain.add_dimensions(y, x);
547 5 : CreateAttribute(eDataType, pszAttrName, nBands);
548 5 : break;
549 7 : case PIXEL:
550 7 : assert(poBands);
551 7 : domain.add_dimensions(y, x, *poBands);
552 7 : CreateAttribute(eDataType, pszAttrName, 1);
553 7 : break;
554 81 : default: // BAND
555 81 : assert(poBands);
556 81 : domain.add_dimensions(*poBands, y, x);
557 81 : CreateAttribute(eDataType, pszAttrName, 1);
558 81 : break;
559 : }
560 :
561 93 : return eErr;
562 : }
563 :
564 : /************************************************************************/
565 : /* FlushCache() */
566 : /************************************************************************/
567 :
568 196 : CPLErr TileDBRasterDataset::FlushCache(bool bAtClosing)
569 :
570 : {
571 196 : const CPLErr eErr = BlockBasedFlushCache(bAtClosing);
572 :
573 196 : if (nPamFlags & GPF_DIRTY)
574 91 : TrySaveXML();
575 196 : return eErr;
576 : }
577 :
578 : /************************************************************************/
579 : /* TrySaveXML() */
580 : /************************************************************************/
581 :
582 91 : CPLErr TileDBRasterDataset::TrySaveXML()
583 :
584 : {
585 91 : if (m_array == nullptr)
586 0 : return CE_None;
587 :
588 91 : CPLXMLNode *psTree = nullptr;
589 : try
590 : {
591 182 : tiledb::VFS vfs(*m_ctx, m_ctx->config());
592 :
593 91 : nPamFlags &= ~GPF_DIRTY;
594 :
595 91 : if (psPam == nullptr || (nPamFlags & GPF_NOSAVE))
596 0 : return CE_None;
597 :
598 : /* --------------------------------------------------------------------
599 : */
600 : /* Make sure we know the filename we want to store in. */
601 : /* --------------------------------------------------------------------
602 : */
603 91 : if (!BuildPamFilename())
604 0 : return CE_None;
605 :
606 : /* --------------------------------------------------------------------
607 : */
608 : /* Build the XML representation of the auxiliary metadata. */
609 : /* --------------------------------------------------------------------
610 : */
611 91 : psTree = SerializeToXML(nullptr);
612 :
613 91 : if (psTree == nullptr)
614 : {
615 : /* If we have unset all metadata, we have to delete the PAM file */
616 0 : m_array->delete_metadata(GDAL_ATTRIBUTE_NAME);
617 0 : return CE_None;
618 : }
619 :
620 91 : if (psSubDatasetsTree != nullptr)
621 : {
622 3 : CPLAddXMLChild(psTree, CPLCloneXMLTree(psSubDatasetsTree->psChild));
623 : }
624 :
625 : /* --------------------------------------------------------------------
626 : */
627 : /* If we are working with a subdataset, we need to integrate */
628 : /* the subdataset tree within the whole existing pam tree, */
629 : /* after removing any old version of the same subdataset. */
630 : /* --------------------------------------------------------------------
631 : */
632 91 : if (!psPam->osSubdatasetName.empty())
633 : {
634 : CPLXMLNode *psOldTree, *psSubTree;
635 :
636 0 : CPLErrorReset();
637 : {
638 0 : CPLErrorHandlerPusher oQuietError(CPLQuietErrorHandler);
639 0 : psOldTree = CPLParseXMLFile(psPam->pszPamFilename);
640 : }
641 :
642 0 : if (psOldTree == nullptr)
643 : psOldTree =
644 0 : CPLCreateXMLNode(nullptr, CXT_Element, "PAMDataset");
645 :
646 0 : for (psSubTree = psOldTree->psChild; psSubTree != nullptr;
647 0 : psSubTree = psSubTree->psNext)
648 : {
649 0 : if (psSubTree->eType != CXT_Element ||
650 0 : !EQUAL(psSubTree->pszValue, "Subdataset"))
651 0 : continue;
652 :
653 0 : if (!EQUAL(CPLGetXMLValue(psSubTree, "name", ""),
654 : psPam->osSubdatasetName))
655 0 : continue;
656 :
657 0 : break;
658 : }
659 :
660 0 : if (psSubTree == nullptr)
661 : {
662 : psSubTree =
663 0 : CPLCreateXMLNode(psOldTree, CXT_Element, "Subdataset");
664 0 : CPLCreateXMLNode(
665 : CPLCreateXMLNode(psSubTree, CXT_Attribute, "name"),
666 0 : CXT_Text, psPam->osSubdatasetName);
667 : }
668 :
669 : CPLXMLNode *psOldPamDataset =
670 0 : CPLGetXMLNode(psSubTree, "PAMDataset");
671 0 : if (psOldPamDataset != nullptr)
672 : {
673 0 : CPLRemoveXMLChild(psSubTree, psOldPamDataset);
674 0 : CPLDestroyXMLNode(psOldPamDataset);
675 : }
676 :
677 0 : CPLAddXMLChild(psSubTree, psTree);
678 0 : psTree = psOldTree;
679 : }
680 :
681 : /* --------------------------------------------------------------------
682 : */
683 : /* Try saving the auxiliary metadata. */
684 : /* --------------------------------------------------------------------
685 : */
686 91 : CPLErrorHandlerPusher oQuietError(CPLQuietErrorHandler);
687 91 : char *pszTree = CPLSerializeXMLTree(psTree);
688 :
689 91 : if (eAccess == GA_ReadOnly)
690 : {
691 40 : if (nTimestamp)
692 : {
693 : auto oMeta = std::unique_ptr<tiledb::Array>(new tiledb::Array(
694 0 : *m_ctx, m_array->uri(), TILEDB_WRITE,
695 0 : tiledb::TemporalPolicy(tiledb::TimeTravel, nTimestamp)));
696 0 : oMeta->put_metadata(GDAL_ATTRIBUTE_NAME, TILEDB_UINT8,
697 0 : static_cast<int>(strlen(pszTree)), pszTree);
698 0 : oMeta->close();
699 : }
700 : else
701 : {
702 : auto oMeta = std::unique_ptr<tiledb::Array>(
703 80 : new tiledb::Array(*m_ctx, m_array->uri(), TILEDB_WRITE));
704 80 : oMeta->put_metadata(GDAL_ATTRIBUTE_NAME, TILEDB_UINT8,
705 40 : static_cast<int>(strlen(pszTree)), pszTree);
706 40 : oMeta->close();
707 : }
708 : }
709 : else
710 : {
711 51 : m_array->put_metadata("dataset_type", TILEDB_STRING_UTF8,
712 : static_cast<int>(strlen(RASTER_DATASET_TYPE)),
713 : RASTER_DATASET_TYPE);
714 102 : m_array->put_metadata(GDAL_ATTRIBUTE_NAME, TILEDB_UINT8,
715 51 : static_cast<int>(strlen(pszTree)), pszTree);
716 : }
717 :
718 91 : CPLFree(pszTree);
719 :
720 : /* --------------------------------------------------------------------
721 : */
722 : /* Cleanup */
723 : /* --------------------------------------------------------------------
724 : */
725 91 : if (psTree)
726 91 : CPLDestroyXMLNode(psTree);
727 :
728 91 : return CE_None;
729 : }
730 0 : catch (const std::exception &e)
731 : {
732 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
733 0 : if (psTree)
734 0 : CPLDestroyXMLNode(psTree);
735 0 : return CE_Failure;
736 : }
737 : }
738 :
739 : /************************************************************************/
740 : /* TryLoadXML() */
741 : /************************************************************************/
742 :
743 82 : CPLErr TileDBRasterDataset::TryLoadXML(char **papszSiblingFiles)
744 :
745 : {
746 82 : return TryLoadCachedXML(papszSiblingFiles, true);
747 : }
748 :
749 : /************************************************************************/
750 : /* TryLoadCachedXML() */
751 : /************************************************************************/
752 :
753 163 : CPLErr TileDBRasterDataset::TryLoadCachedXML(char ** /*papszSiblingFiles*/,
754 : bool bReload)
755 :
756 : {
757 163 : CPLXMLNode *psTree = nullptr;
758 : try
759 : {
760 163 : PamInitialize();
761 326 : tiledb::VFS vfs(*m_ctx, m_ctx->config());
762 :
763 : /* --------------------------------------------------------------------
764 : */
765 : /* Clear dirty flag. Generally when we get to this point is */
766 : /* from a call at the end of the Open() method, and some calls */
767 : /* may have already marked the PAM info as dirty (for instance */
768 : /* setting metadata), but really everything to this point is */
769 : /* reproducible, and so the PAM info should not really be */
770 : /* thought of as dirty. */
771 : /* --------------------------------------------------------------------
772 : */
773 163 : nPamFlags &= ~GPF_DIRTY;
774 :
775 : /* --------------------------------------------------------------------
776 : */
777 : /* Try reading the file. */
778 : /* --------------------------------------------------------------------
779 : */
780 163 : if (!BuildPamFilename())
781 0 : return CE_None;
782 :
783 : /* --------------------------------------------------------------------
784 : */
785 : /* In case the PAM filename is a .aux.xml file next to the */
786 : /* physical file and we have a siblings list, then we can skip */
787 : /* stat'ing the filesystem. */
788 : /* --------------------------------------------------------------------
789 : */
790 163 : CPLErr eLastErr = CPLGetLastErrorType();
791 163 : int nLastErrNo = CPLGetLastErrorNo();
792 326 : CPLString osLastErrorMsg = CPLGetLastErrorMsg();
793 :
794 163 : CPLErrorReset();
795 : {
796 326 : CPLErrorHandlerPusher oQuietError(CPLQuietErrorHandler);
797 :
798 163 : if (bReload)
799 : {
800 82 : tiledb_datatype_t v_type =
801 : TILEDB_UINT8; // CPLSerializeXMLTree returns char*
802 82 : const void *v_r = nullptr;
803 82 : uint32_t v_num = 0;
804 82 : if ((eAccess == GA_Update) && (m_roArray))
805 : {
806 13 : m_roArray->get_metadata(GDAL_ATTRIBUTE_NAME, &v_type,
807 : &v_num, &v_r);
808 13 : if (v_r)
809 : {
810 : osMetaDoc =
811 13 : CPLString(static_cast<const char *>(v_r), v_num);
812 : }
813 : }
814 : else
815 : {
816 69 : m_array->get_metadata(GDAL_ATTRIBUTE_NAME, &v_type, &v_num,
817 : &v_r);
818 :
819 69 : if (v_r)
820 : {
821 : osMetaDoc =
822 68 : CPLString(static_cast<const char *>(v_r), v_num);
823 : }
824 : }
825 82 : psTree = CPLParseXMLString(osMetaDoc);
826 : }
827 :
828 164 : if (bReload && psTree == nullptr &&
829 164 : vfs.is_file(psPam->pszPamFilename))
830 : {
831 1 : auto nBytes = vfs.file_size(psPam->pszPamFilename);
832 2 : tiledb::VFS::filebuf fbuf(vfs);
833 1 : fbuf.open(psPam->pszPamFilename, std::ios::in);
834 1 : std::istream is(&fbuf);
835 1 : osMetaDoc.resize(nBytes);
836 1 : is.read((char *)osMetaDoc.data(), nBytes);
837 1 : fbuf.close();
838 1 : psTree = CPLParseXMLString(osMetaDoc);
839 : }
840 :
841 163 : if (!bReload)
842 : {
843 81 : psTree = CPLParseXMLString(osMetaDoc);
844 : }
845 : }
846 163 : CPLErrorReset();
847 :
848 163 : if (eLastErr != CE_None)
849 0 : CPLErrorSetState(eLastErr, nLastErrNo, osLastErrorMsg.c_str());
850 :
851 : /* --------------------------------------------------------------------
852 : */
853 : /* If we are looking for a subdataset, search for its subtree not.
854 : */
855 : /* --------------------------------------------------------------------
856 : */
857 163 : if (psTree && !psPam->osSubdatasetName.empty())
858 : {
859 11 : CPLXMLNode *psSubTree = psTree->psChild;
860 :
861 57 : for (; psSubTree != nullptr; psSubTree = psSubTree->psNext)
862 : {
863 56 : if (psSubTree->eType != CXT_Element ||
864 56 : !EQUAL(psSubTree->pszValue, "Subdataset"))
865 38 : continue;
866 :
867 18 : if (!EQUAL(CPLGetXMLValue(psSubTree, "name", ""),
868 : psPam->osSubdatasetName))
869 8 : continue;
870 :
871 10 : psSubTree = CPLGetXMLNode(psSubTree, "PAMDataset");
872 10 : break;
873 : }
874 :
875 11 : if (psSubTree != nullptr)
876 10 : psSubTree = CPLCloneXMLTree(psSubTree);
877 :
878 11 : CPLDestroyXMLNode(psTree);
879 11 : psTree = psSubTree;
880 : }
881 163 : if (psTree == nullptr)
882 1 : return CE_Failure;
883 :
884 : /* --------------------------------------------------------------------
885 : */
886 : /* Initialize ourselves from this XML tree. */
887 : /* --------------------------------------------------------------------
888 : */
889 :
890 162 : CPLString osPath(CPLGetPath(psPam->pszPamFilename));
891 162 : const CPLErr eErr = XMLInit(psTree, osPath);
892 :
893 162 : CPLDestroyXMLNode(psTree);
894 :
895 162 : if (eErr != CE_None)
896 0 : PamClear();
897 :
898 162 : return eErr;
899 : }
900 0 : catch (const tiledb::TileDBError &e)
901 : {
902 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
903 0 : if (psTree)
904 0 : CPLDestroyXMLNode(psTree);
905 0 : return CE_Failure;
906 : }
907 : }
908 :
909 : /************************************************************************/
910 : /* GetMetadata() */
911 : /************************************************************************/
912 :
913 164 : char **TileDBRasterDataset::GetMetadata(const char *pszDomain)
914 :
915 : {
916 164 : if (pszDomain != nullptr && EQUAL(pszDomain, "SUBDATASETS"))
917 : {
918 3 : char **papszMeta = CSLDuplicate(GDALPamDataset::GetMetadata(pszDomain));
919 3 : for (int i = 0; papszMeta && papszMeta[i]; i++)
920 : {
921 0 : if (STARTS_WITH(papszMeta[i], "SUBDATASET_") &&
922 0 : strstr(papszMeta[i], "_NAME="))
923 : {
924 0 : char *pszKey = nullptr;
925 0 : const char *pszAttr = CPLParseNameValue(papszMeta[i], &pszKey);
926 0 : if (pszAttr && !STARTS_WITH(pszAttr, "TILEDB:"))
927 : {
928 0 : CPLString osAttr(pszAttr);
929 0 : CPLFree(papszMeta[i]);
930 0 : papszMeta[i] = CPLStrdup(
931 0 : CPLString().Printf("%s=TILEDB:\"%s\":%s", pszKey,
932 0 : GetDescription(), osAttr.c_str()));
933 : }
934 0 : CPLFree(pszKey);
935 : }
936 : }
937 3 : m_osSubdatasetMD.Assign(papszMeta);
938 3 : return m_osSubdatasetMD.List();
939 : }
940 : else
941 : {
942 161 : return GDALPamDataset::GetMetadata(pszDomain);
943 : }
944 : }
945 :
946 : /************************************************************************/
947 : /* Open() */
948 : /************************************************************************/
949 :
950 82 : GDALDataset *TileDBRasterDataset::Open(GDALOpenInfo *poOpenInfo)
951 :
952 : {
953 164 : auto poDS = std::make_unique<TileDBRasterDataset>();
954 :
955 : const char *pszConfig =
956 82 : CSLFetchNameValue(poOpenInfo->papszOpenOptions, "TILEDB_CONFIG");
957 :
958 : const char *pszTimestamp =
959 82 : CSLFetchNameValue(poOpenInfo->papszOpenOptions, "TILEDB_TIMESTAMP");
960 :
961 164 : poDS->bStats =
962 82 : CSLFetchBoolean(poOpenInfo->papszOpenOptions, "STATS", FALSE);
963 :
964 82 : if (pszConfig != nullptr)
965 : {
966 0 : tiledb::Config cfg(pszConfig);
967 0 : poDS->m_ctx.reset(new tiledb::Context(cfg));
968 : }
969 : else
970 : {
971 82 : poDS->m_ctx.reset(new tiledb::Context());
972 : }
973 82 : if (pszTimestamp)
974 18 : poDS->nTimestamp = std::strtoull(pszTimestamp, nullptr, 10);
975 :
976 164 : CPLString osArrayPath;
977 164 : CPLString osAux;
978 164 : CPLString osSubdataset;
979 :
980 164 : std::string osAttrNameTmp;
981 : const char *pszAttr =
982 82 : CSLFetchNameValue(poOpenInfo->papszOpenOptions, "TILEDB_ATTRIBUTE");
983 :
984 82 : if (STARTS_WITH_CI(poOpenInfo->pszFilename, "TILEDB:") &&
985 2 : !STARTS_WITH_CI(poOpenInfo->pszFilename, "TILEDB://"))
986 : {
987 : // form required read attributes and open file
988 : // Create a corresponding GDALDataset.
989 : CPLStringList apszName(
990 2 : CSLTokenizeString2(poOpenInfo->pszFilename, ":",
991 2 : CSLT_HONOURSTRINGS | CSLT_PRESERVEESCAPES));
992 :
993 2 : if (apszName.size() != 3)
994 : {
995 0 : return nullptr;
996 : }
997 :
998 2 : osArrayPath = TileDBDataset::VSI_to_tiledb_uri(apszName[1]);
999 2 : osSubdataset = apszName[2];
1000 4 : poDS->SetSubdatasetName(osSubdataset.c_str());
1001 : }
1002 : else
1003 : {
1004 80 : if (pszAttr != nullptr)
1005 : {
1006 4 : poDS->SetSubdatasetName(pszAttr);
1007 : }
1008 :
1009 80 : osArrayPath = TileDBDataset::VSI_to_tiledb_uri(poOpenInfo->pszFilename);
1010 : }
1011 :
1012 82 : const char *pszArrayName = CPLGetBasename(osArrayPath);
1013 82 : osAux.Printf("%s.tdb", pszArrayName);
1014 :
1015 : // aux file is in array folder
1016 82 : poDS->SetPhysicalFilename(CPLFormFilename(osArrayPath, osAux, nullptr));
1017 : // Initialize any PAM information.
1018 82 : poDS->SetDescription(osArrayPath);
1019 :
1020 82 : tiledb_query_type_t eMode = TILEDB_READ;
1021 82 : if (poOpenInfo->eAccess == GA_Update)
1022 : {
1023 13 : eMode = TILEDB_WRITE;
1024 13 : poDS->m_roCtx.reset(new tiledb::Context(poDS->m_ctx->config()));
1025 39 : poDS->m_roArray.reset(
1026 26 : new tiledb::Array(*poDS->m_roCtx, osArrayPath, TILEDB_READ));
1027 : }
1028 :
1029 82 : if (poDS->nTimestamp)
1030 : {
1031 18 : if (eMode == TILEDB_READ)
1032 : {
1033 30 : poDS->m_array.reset(new tiledb::Array(
1034 10 : *poDS->m_ctx, osArrayPath, TILEDB_READ,
1035 30 : tiledb::TemporalPolicy(tiledb::TimeTravel, poDS->nTimestamp)));
1036 : }
1037 : else
1038 : {
1039 24 : poDS->m_array.reset(new tiledb::Array(
1040 8 : *poDS->m_ctx, osArrayPath, TILEDB_WRITE,
1041 24 : tiledb::TemporalPolicy(tiledb::TimeTravel, poDS->nTimestamp)));
1042 : }
1043 : }
1044 : else
1045 192 : poDS->m_array.reset(
1046 128 : new tiledb::Array(*poDS->m_ctx, osArrayPath, eMode));
1047 :
1048 82 : poDS->eAccess = poOpenInfo->eAccess;
1049 :
1050 : // dependent on PAM metadata for information about array
1051 82 : poDS->TryLoadXML();
1052 :
1053 164 : tiledb::ArraySchema schema = poDS->m_array->schema();
1054 :
1055 82 : char **papszStructMeta = poDS->GetMetadata("IMAGE_STRUCTURE");
1056 82 : const char *pszXSize = CSLFetchNameValue(papszStructMeta, "X_SIZE");
1057 82 : if (pszXSize)
1058 : {
1059 67 : poDS->nRasterXSize = atoi(pszXSize);
1060 67 : if (poDS->nRasterXSize <= 0)
1061 : {
1062 0 : CPLError(CE_Failure, CPLE_AppDefined,
1063 : "Width %i should be greater than zero.",
1064 0 : poDS->nRasterXSize);
1065 0 : return nullptr;
1066 : }
1067 : }
1068 :
1069 82 : const char *pszYSize = CSLFetchNameValue(papszStructMeta, "Y_SIZE");
1070 82 : if (pszYSize)
1071 : {
1072 67 : poDS->nRasterYSize = atoi(pszYSize);
1073 67 : if (poDS->nRasterYSize <= 0)
1074 : {
1075 0 : CPLError(CE_Failure, CPLE_AppDefined,
1076 : "Height %i should be greater than zero.",
1077 0 : poDS->nRasterYSize);
1078 0 : return nullptr;
1079 : }
1080 : }
1081 :
1082 82 : const char *pszNBits = CSLFetchNameValue(papszStructMeta, "NBITS");
1083 82 : if (pszNBits)
1084 : {
1085 67 : poDS->nBitsPerSample = atoi(pszNBits);
1086 : }
1087 :
1088 82 : const char *pszDataType = CSLFetchNameValue(papszStructMeta, "DATA_TYPE");
1089 82 : if (pszDataType)
1090 : {
1091 : // handle the case where arrays have been written with int type
1092 : // (2.5.0)
1093 67 : GDALDataType eDT = GDALGetDataTypeByName(pszDataType);
1094 67 : if (eDT == GDT_Unknown)
1095 : {
1096 1 : int t = atoi(pszDataType);
1097 1 : if ((t > 0) && (t < GDT_TypeCount))
1098 1 : poDS->eDataType = static_cast<GDALDataType>(atoi(pszDataType));
1099 : else
1100 : {
1101 0 : CPLError(CE_Failure, CPLE_AppDefined, "Unknown data type %s.",
1102 : pszDataType);
1103 0 : return nullptr;
1104 : }
1105 : }
1106 : else
1107 : {
1108 66 : poDS->eDataType = eDT;
1109 : }
1110 : }
1111 : else
1112 : {
1113 15 : if (!pszAttr)
1114 : {
1115 15 : if (schema.attribute_num() == 1)
1116 : {
1117 14 : osAttrNameTmp = schema.attribute(0).name();
1118 14 : pszAttr = osAttrNameTmp.c_str();
1119 : }
1120 : }
1121 : }
1122 :
1123 82 : const char *pszIndexMode = CSLFetchNameValue(papszStructMeta, "INTERLEAVE");
1124 :
1125 82 : if (pszIndexMode)
1126 61 : option_to_index_type(pszIndexMode, poDS->eIndexMode);
1127 :
1128 164 : std::vector<tiledb::Dimension> dims = schema.domain().dimensions();
1129 :
1130 82 : int iYDim = 0;
1131 82 : int iXDim = 1;
1132 82 : if ((dims.size() == 2) || (dims.size() == 3))
1133 : {
1134 82 : if (dims.size() == 3)
1135 : {
1136 96 : if ((pszAttr != nullptr) &&
1137 96 : (schema.attributes().count(pszAttr) == 0))
1138 : {
1139 0 : CPLError(CE_Failure, CPLE_NotSupported,
1140 : "%s attribute is not found in TileDB schema.",
1141 : pszAttr);
1142 0 : return nullptr;
1143 : }
1144 :
1145 78 : if (poDS->eIndexMode == PIXEL)
1146 14 : std::rotate(dims.begin(), dims.begin() + 2, dims.end());
1147 :
1148 78 : if (dims[0].type() != TILEDB_UINT64)
1149 : {
1150 0 : const char *pszTypeName = "";
1151 0 : tiledb_datatype_to_str(dims[0].type(), &pszTypeName);
1152 0 : CPLError(CE_Failure, CPLE_NotSupported,
1153 : "Unsupported BAND dimension type: %s", pszTypeName);
1154 0 : return nullptr;
1155 : }
1156 78 : poDS->nBandStart = dims[0].domain<uint64_t>().first;
1157 78 : const uint64_t nBandEnd = dims[0].domain<uint64_t>().second;
1158 156 : if (nBandEnd < poDS->nBandStart ||
1159 78 : nBandEnd - poDS->nBandStart >
1160 78 : static_cast<uint64_t>(std::numeric_limits<int>::max() - 1))
1161 : {
1162 0 : CPLError(CE_Failure, CPLE_NotSupported,
1163 : "Invalid bounds for BAND dimension.");
1164 0 : return nullptr;
1165 : }
1166 78 : poDS->nBands = static_cast<int>(nBandEnd - poDS->nBandStart + 1);
1167 78 : iYDim = 1;
1168 78 : iXDim = 2;
1169 : }
1170 : else
1171 : {
1172 : const char *pszBands =
1173 4 : poDS->GetMetadataItem("NUM_BANDS", "IMAGE_STRUCTURE");
1174 4 : if (pszBands)
1175 : {
1176 1 : poDS->nBands = atoi(pszBands);
1177 : }
1178 :
1179 4 : poDS->eIndexMode = ATTRIBUTES;
1180 : }
1181 : }
1182 : else
1183 : {
1184 0 : CPLError(CE_Failure, CPLE_AppDefined,
1185 : "Wrong number of dimensions %d: expected 2 or 3.",
1186 0 : static_cast<int>(dims.size()));
1187 0 : return nullptr;
1188 : }
1189 :
1190 164 : if (!GDALCheckDatasetDimensions(poDS->nRasterXSize, poDS->nRasterYSize) ||
1191 82 : !GDALCheckBandCount(poDS->nBands, /*bIsZeroAllowed=*/true))
1192 : {
1193 0 : return nullptr;
1194 : }
1195 :
1196 82 : if (dims[iYDim].type() != TILEDB_UINT64)
1197 : {
1198 0 : const char *pszTypeName = "";
1199 0 : tiledb_datatype_to_str(dims[0].type(), &pszTypeName);
1200 0 : CPLError(CE_Failure, CPLE_NotSupported,
1201 : "Unsupported Y dimension type: %s", pszTypeName);
1202 0 : return nullptr;
1203 : }
1204 82 : if (!pszYSize)
1205 : {
1206 15 : const uint64_t nStart = dims[iYDim].domain<uint64_t>().first;
1207 15 : const uint64_t nEnd = dims[iYDim].domain<uint64_t>().second;
1208 30 : if (nStart != 0 ||
1209 15 : nEnd > static_cast<uint64_t>(std::numeric_limits<int>::max() - 1))
1210 : {
1211 0 : CPLError(CE_Failure, CPLE_NotSupported,
1212 : "Invalid bounds for Y dimension.");
1213 0 : return nullptr;
1214 : }
1215 15 : poDS->nRasterYSize = static_cast<int>(nEnd - nStart + 1);
1216 : }
1217 82 : const uint64_t nBlockYSize = dims[iYDim].tile_extent<uint64_t>();
1218 82 : if (nBlockYSize > static_cast<uint64_t>(std::numeric_limits<int>::max()))
1219 : {
1220 0 : CPLError(CE_Failure, CPLE_NotSupported, "Too large block Y size.");
1221 0 : return nullptr;
1222 : }
1223 82 : poDS->nBlockYSize = static_cast<int>(nBlockYSize);
1224 :
1225 82 : if (dims[iXDim].type() != TILEDB_UINT64)
1226 : {
1227 0 : const char *pszTypeName = "";
1228 0 : tiledb_datatype_to_str(dims[0].type(), &pszTypeName);
1229 0 : CPLError(CE_Failure, CPLE_NotSupported,
1230 : "Unsupported Y dimension type: %s", pszTypeName);
1231 0 : return nullptr;
1232 : }
1233 82 : if (!pszXSize)
1234 : {
1235 15 : const uint64_t nStart = dims[iXDim].domain<uint64_t>().first;
1236 15 : const uint64_t nEnd = dims[iXDim].domain<uint64_t>().second;
1237 30 : if (nStart != 0 ||
1238 15 : nEnd > static_cast<uint64_t>(std::numeric_limits<int>::max() - 1))
1239 : {
1240 0 : CPLError(CE_Failure, CPLE_NotSupported,
1241 : "Invalid bounds for X dimension.");
1242 0 : return nullptr;
1243 : }
1244 15 : poDS->nRasterXSize = static_cast<int>(nEnd - nStart + 1);
1245 : }
1246 82 : const uint64_t nBlockXSize = dims[iXDim].tile_extent<uint64_t>();
1247 82 : if (nBlockXSize > static_cast<uint64_t>(std::numeric_limits<int>::max()))
1248 : {
1249 0 : CPLError(CE_Failure, CPLE_NotSupported, "Too large block X size.");
1250 0 : return nullptr;
1251 : }
1252 82 : poDS->nBlockXSize = static_cast<int>(nBlockXSize);
1253 :
1254 82 : poDS->nBlocksX = DIV_ROUND_UP(poDS->nRasterXSize, poDS->nBlockXSize);
1255 82 : poDS->nBlocksY = DIV_ROUND_UP(poDS->nRasterYSize, poDS->nBlockYSize);
1256 :
1257 82 : if (dims.size() == 3)
1258 : {
1259 : // Create band information objects.
1260 718 : for (int i = 1; i <= poDS->nBands; ++i)
1261 : {
1262 640 : if (pszAttr == nullptr)
1263 108 : poDS->SetBand(i, new TileDBRasterBand(poDS.get(), i));
1264 : else
1265 1064 : poDS->SetBand(
1266 1064 : i, new TileDBRasterBand(poDS.get(), i, CPLString(pszAttr)));
1267 : }
1268 : }
1269 : else // subdatasets or only attributes
1270 : {
1271 5 : if ((poOpenInfo->eAccess == GA_Update) &&
1272 1 : (poDS->GetMetadata("SUBDATASETS") != nullptr))
1273 : {
1274 0 : CPLError(CE_Failure, CPLE_NotSupported,
1275 : "The TileDB driver does not support update access "
1276 : "to subdatasets.");
1277 0 : return nullptr;
1278 : }
1279 :
1280 4 : if (!osSubdataset.empty())
1281 : {
1282 : // do we have the attribute in the schema
1283 2 : if (schema.attributes().count(osSubdataset))
1284 : {
1285 0 : poDS->SetBand(
1286 0 : 1, new TileDBRasterBand(poDS.get(), 1, osSubdataset));
1287 : }
1288 : else
1289 : {
1290 2 : if (schema.attributes().count(osSubdataset + "_1"))
1291 : {
1292 : // Create band information objects.
1293 2 : for (int i = 1; i <= poDS->nBands; ++i)
1294 : {
1295 1 : CPLString osAttr = CPLString().Printf(
1296 2 : "%s_%d", osSubdataset.c_str(), i);
1297 2 : poDS->SetBand(
1298 1 : i, new TileDBRasterBand(poDS.get(), i, osAttr));
1299 : }
1300 : }
1301 : else
1302 : {
1303 1 : CPLError(CE_Failure, CPLE_NotSupported,
1304 : "%s attribute is not found in TileDB schema.",
1305 : osSubdataset.c_str());
1306 1 : return nullptr;
1307 : }
1308 : }
1309 : }
1310 : else
1311 : {
1312 2 : char **papszMeta = poDS->GetMetadata("SUBDATASETS");
1313 2 : if (papszMeta != nullptr)
1314 : {
1315 0 : if ((CSLCount(papszMeta) / 2) == 1)
1316 : {
1317 : CPLString osDSName = CSLFetchNameValueDef(
1318 0 : poDS->papszSubDatasets, "SUBDATASET_1_NAME", "");
1319 0 : return (GDALDataset *)GDALOpen(osDSName,
1320 0 : poOpenInfo->eAccess);
1321 : }
1322 : }
1323 2 : else if (poDS->eIndexMode == ATTRIBUTES)
1324 : {
1325 2 : poDS->nBands = schema.attribute_num();
1326 :
1327 : // Create band information objects.
1328 8 : for (int i = 1; i <= poDS->nBands; ++i)
1329 : {
1330 : CPLString osAttr =
1331 18 : TILEDB_VALUES + CPLString().Printf("_%i", i);
1332 12 : poDS->SetBand(i,
1333 6 : new TileDBRasterBand(poDS.get(), i, osAttr));
1334 : }
1335 : }
1336 : else
1337 : {
1338 0 : CPLError(CE_Failure, CPLE_NotSupported,
1339 : "%s is missing required TileDB subdataset metadata.",
1340 : osArrayPath.c_str());
1341 0 : return nullptr;
1342 : }
1343 : }
1344 : }
1345 :
1346 : // reload metadata now that bands are created to populate band metadata
1347 81 : poDS->TryLoadCachedXML(nullptr, false);
1348 :
1349 162 : tiledb::VFS vfs(*poDS->m_ctx, poDS->m_ctx->config());
1350 :
1351 81 : if (!STARTS_WITH_CI(osArrayPath, "TILEDB:") && vfs.is_dir(osArrayPath))
1352 81 : poDS->oOvManager.Initialize(poDS.get(), ":::VIRTUAL:::");
1353 : else
1354 0 : CPLError(CE_Warning, CPLE_AppDefined,
1355 : "Overviews not supported for network writes.");
1356 :
1357 81 : return poDS.release();
1358 : }
1359 :
1360 : /************************************************************************/
1361 : /* CreateAttribute() */
1362 : /************************************************************************/
1363 :
1364 99 : CPLErr TileDBRasterDataset::CreateAttribute(GDALDataType eType,
1365 : const CPLString &osAttrName,
1366 : const int nSubRasterCount)
1367 : {
1368 : try
1369 : {
1370 203 : for (int i = 0; i < nSubRasterCount; ++i)
1371 : {
1372 104 : CPLString osName(osAttrName);
1373 : // a few special cases
1374 : // remove any leading slashes or
1375 : // additional slashes as in the case of hdf5
1376 104 : if STARTS_WITH (osName, "//")
1377 : {
1378 2 : osName = osName.substr(2);
1379 : }
1380 :
1381 104 : osName.replaceAll("/", "_");
1382 104 : CPLString osPrettyName = osName;
1383 :
1384 104 : if ((eIndexMode == ATTRIBUTES) ||
1385 92 : ((bHasSubDatasets) && (nSubRasterCount > 1)))
1386 : {
1387 12 : osName = CPLString().Printf("%s_%d", osName.c_str(), i + 1);
1388 : }
1389 :
1390 104 : switch (eType)
1391 : {
1392 47 : case GDT_Byte:
1393 : {
1394 47 : m_schema->add_attribute(
1395 94 : tiledb::Attribute::create<unsigned char>(
1396 94 : *m_ctx, osName, *m_filterList));
1397 47 : nBitsPerSample = 8;
1398 47 : break;
1399 : }
1400 3 : case GDT_Int8:
1401 : {
1402 9 : m_schema->add_attribute(tiledb::Attribute::create<int8_t>(
1403 6 : *m_ctx, osName, *m_filterList));
1404 3 : nBitsPerSample = 8;
1405 3 : break;
1406 : }
1407 4 : case GDT_UInt16:
1408 : {
1409 4 : m_schema->add_attribute(
1410 8 : tiledb::Attribute::create<unsigned short>(
1411 8 : *m_ctx, osName, *m_filterList));
1412 4 : nBitsPerSample = 16;
1413 4 : break;
1414 : }
1415 4 : case GDT_UInt32:
1416 : {
1417 4 : m_schema->add_attribute(
1418 8 : tiledb::Attribute::create<unsigned int>(*m_ctx, osName,
1419 8 : *m_filterList));
1420 4 : nBitsPerSample = 32;
1421 4 : break;
1422 : }
1423 3 : case GDT_UInt64:
1424 : {
1425 9 : m_schema->add_attribute(tiledb::Attribute::create<uint64_t>(
1426 6 : *m_ctx, osName, *m_filterList));
1427 3 : nBitsPerSample = 64;
1428 3 : break;
1429 : }
1430 4 : case GDT_Int16:
1431 : {
1432 12 : m_schema->add_attribute(tiledb::Attribute::create<short>(
1433 8 : *m_ctx, osName, *m_filterList));
1434 4 : nBitsPerSample = 16;
1435 4 : break;
1436 : }
1437 6 : case GDT_Int32:
1438 : {
1439 18 : m_schema->add_attribute(tiledb::Attribute::create<int>(
1440 12 : *m_ctx, osName, *m_filterList));
1441 6 : nBitsPerSample = 32;
1442 6 : break;
1443 : }
1444 3 : case GDT_Int64:
1445 : {
1446 9 : m_schema->add_attribute(tiledb::Attribute::create<int64_t>(
1447 6 : *m_ctx, osName, *m_filterList));
1448 3 : nBitsPerSample = 64;
1449 3 : break;
1450 : }
1451 8 : case GDT_Float32:
1452 : {
1453 24 : m_schema->add_attribute(tiledb::Attribute::create<float>(
1454 16 : *m_ctx, osName, *m_filterList));
1455 8 : nBitsPerSample = 32;
1456 8 : break;
1457 : }
1458 4 : case GDT_Float64:
1459 : {
1460 12 : m_schema->add_attribute(tiledb::Attribute::create<double>(
1461 8 : *m_ctx, osName, *m_filterList));
1462 4 : nBitsPerSample = 64;
1463 4 : break;
1464 : }
1465 4 : case GDT_CInt16:
1466 : {
1467 12 : m_schema->add_attribute(tiledb::Attribute::create<short[2]>(
1468 8 : *m_ctx, osName, *m_filterList));
1469 4 : nBitsPerSample = 16;
1470 4 : break;
1471 : }
1472 4 : case GDT_CInt32:
1473 : {
1474 12 : m_schema->add_attribute(tiledb::Attribute::create<int[2]>(
1475 8 : *m_ctx, osName, *m_filterList));
1476 4 : nBitsPerSample = 32;
1477 4 : break;
1478 : }
1479 4 : case GDT_CFloat32:
1480 : {
1481 12 : m_schema->add_attribute(tiledb::Attribute::create<float[2]>(
1482 8 : *m_ctx, osName, *m_filterList));
1483 4 : nBitsPerSample = 32;
1484 4 : break;
1485 : }
1486 6 : case GDT_CFloat64:
1487 : {
1488 6 : m_schema->add_attribute(
1489 12 : tiledb::Attribute::create<double[2]>(*m_ctx, osName,
1490 12 : *m_filterList));
1491 6 : nBitsPerSample = 64;
1492 6 : break;
1493 : }
1494 0 : default:
1495 0 : return CE_Failure;
1496 : }
1497 :
1498 104 : if ((bHasSubDatasets) && (i == 0))
1499 : {
1500 6 : ++nSubDataCount;
1501 :
1502 6 : CPLString osDim;
1503 6 : switch (nSubRasterCount)
1504 : {
1505 0 : case 2:
1506 0 : osDim.Printf("%dx%d", nRasterXSize, nRasterYSize);
1507 0 : break;
1508 6 : default:
1509 : osDim.Printf("%dx%dx%d", nSubRasterCount, nRasterXSize,
1510 6 : nRasterYSize);
1511 6 : break;
1512 : }
1513 :
1514 6 : papszSubDatasets = CSLSetNameValue(
1515 : papszSubDatasets,
1516 12 : CPLString().Printf("SUBDATASET_%d_NAME", nSubDataCount),
1517 12 : CPLString().Printf("%s", osPrettyName.c_str()));
1518 :
1519 6 : papszSubDatasets = CSLSetNameValue(
1520 : papszSubDatasets,
1521 12 : CPLString().Printf("SUBDATASET_%d_DESC", nSubDataCount),
1522 12 : CPLString().Printf("[%s] %s (%s)", osDim.c_str(),
1523 : osPrettyName.c_str(),
1524 6 : GDALGetDataTypeName(eType)));
1525 :
1526 : // add to PAM metadata
1527 6 : if (psSubDatasetsTree == nullptr)
1528 : {
1529 3 : psSubDatasetsTree =
1530 3 : CPLCreateXMLNode(nullptr, CXT_Element, "PAMDataset");
1531 : }
1532 :
1533 6 : CPLXMLNode *psSubNode = CPLCreateXMLNode(
1534 : psSubDatasetsTree, CXT_Element, "Subdataset");
1535 6 : CPLAddXMLAttributeAndValue(psSubNode, "name",
1536 : osPrettyName.c_str());
1537 :
1538 6 : CPLXMLNode *psMetaNode = CPLCreateXMLNode(
1539 : CPLCreateXMLNode(psSubNode, CXT_Element, "PAMDataset"),
1540 : CXT_Element, "Metadata");
1541 6 : CPLAddXMLAttributeAndValue(psMetaNode, "domain",
1542 : "IMAGE_STRUCTURE");
1543 :
1544 6 : CPLAddXMLAttributeAndValue(
1545 : CPLCreateXMLElementAndValue(
1546 : psMetaNode, "MDI",
1547 12 : CPLString().Printf("%d", nRasterXSize)),
1548 : "KEY", "X_SIZE");
1549 :
1550 6 : CPLAddXMLAttributeAndValue(
1551 : CPLCreateXMLElementAndValue(
1552 : psMetaNode, "MDI",
1553 12 : CPLString().Printf("%d", nRasterYSize)),
1554 : "KEY", "Y_SIZE");
1555 :
1556 6 : CPLAddXMLAttributeAndValue(
1557 : CPLCreateXMLElementAndValue(
1558 : psMetaNode, "MDI",
1559 12 : CPLString().Printf("%s", GDALGetDataTypeName(eType))),
1560 : "KEY", "DATA_TYPE");
1561 :
1562 6 : if (lpoAttributeDS.size() > 0)
1563 : {
1564 4 : CPLAddXMLAttributeAndValue(
1565 : CPLCreateXMLElementAndValue(
1566 : psMetaNode, "MDI",
1567 8 : CPLString().Printf("%d", nBands)),
1568 : "KEY", "NUM_BANDS");
1569 : }
1570 : else
1571 : {
1572 2 : CPLAddXMLAttributeAndValue(
1573 : CPLCreateXMLElementAndValue(
1574 : psMetaNode, "MDI",
1575 4 : CPLString().Printf("%d", nSubRasterCount)),
1576 : "KEY", "NUM_BANDS");
1577 : }
1578 :
1579 6 : CPLAddXMLAttributeAndValue(
1580 : CPLCreateXMLElementAndValue(
1581 : psMetaNode, "MDI",
1582 12 : CPLString().Printf("%d", nBitsPerSample)),
1583 : "KEY", "NBITS");
1584 : }
1585 : }
1586 99 : return CE_None;
1587 : }
1588 0 : catch (const tiledb::TileDBError &e)
1589 : {
1590 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
1591 0 : return CE_Failure;
1592 : }
1593 : }
1594 :
1595 : /************************************************************************/
1596 : /* SetBlockSize() */
1597 : /************************************************************************/
1598 :
1599 1 : void TileDBRasterDataset::SetBlockSize(GDALRasterBand *poBand,
1600 : char **&papszOptions)
1601 :
1602 : {
1603 1 : int nX = 0;
1604 1 : int nY = 0;
1605 1 : poBand->GetBlockSize(&nX, &nY);
1606 :
1607 1 : if (CSLFetchNameValue(papszOptions, "BLOCKXSIZE") == nullptr)
1608 : {
1609 1 : papszOptions = CSLSetNameValue(papszOptions, "BLOCKXSIZE",
1610 2 : CPLString().Printf("%d", nX));
1611 : }
1612 :
1613 1 : if (CSLFetchNameValue(papszOptions, "BLOCKYSIZE") == nullptr)
1614 : {
1615 1 : papszOptions = CSLSetNameValue(papszOptions, "BLOCKYSIZE",
1616 2 : CPLString().Printf("%d", nY));
1617 : }
1618 1 : }
1619 :
1620 : /************************************************************************/
1621 : /* CreateLL() */
1622 : /* */
1623 : /* Shared functionality between TileDBDataset::Create() and */
1624 : /* TileDBDataset::CreateCopy() for creating TileDB array based on */
1625 : /* a set of options and a configuration. */
1626 : /************************************************************************/
1627 :
1628 93 : TileDBRasterDataset *TileDBRasterDataset::CreateLL(const char *pszFilename,
1629 : int nXSize, int nYSize,
1630 : int nBandsIn,
1631 : GDALDataType eType,
1632 : char **papszOptions)
1633 : {
1634 : try
1635 : {
1636 93 : if ((nXSize <= 0 && nYSize <= 0))
1637 : {
1638 0 : return nullptr;
1639 : }
1640 :
1641 186 : auto poDS = std::make_unique<TileDBRasterDataset>();
1642 93 : poDS->nRasterXSize = nXSize;
1643 93 : poDS->nRasterYSize = nYSize;
1644 93 : poDS->eDataType = eType;
1645 93 : poDS->nBands = nBandsIn;
1646 93 : poDS->eAccess = GA_Update;
1647 :
1648 93 : if (poDS->nBands == 0) // subdatasets
1649 : {
1650 1 : poDS->eIndexMode = ATTRIBUTES;
1651 : }
1652 : else
1653 : {
1654 : const char *pszIndexMode =
1655 92 : CSLFetchNameValue(papszOptions, "INTERLEAVE");
1656 :
1657 92 : if (option_to_index_type(pszIndexMode, poDS->eIndexMode))
1658 0 : return nullptr;
1659 : }
1660 :
1661 : const char *pszConfig =
1662 93 : CSLFetchNameValue(papszOptions, "TILEDB_CONFIG");
1663 93 : if (pszConfig != nullptr)
1664 : {
1665 0 : tiledb::Config cfg(pszConfig);
1666 0 : poDS->m_ctx.reset(new tiledb::Context(cfg));
1667 : }
1668 : else
1669 : {
1670 93 : poDS->m_ctx.reset(new tiledb::Context());
1671 : }
1672 :
1673 : const char *pszCompression =
1674 93 : CSLFetchNameValue(papszOptions, "COMPRESSION");
1675 : const char *pszCompressionLevel =
1676 93 : CSLFetchNameValue(papszOptions, "COMPRESSION_LEVEL");
1677 :
1678 : const char *pszBlockXSize =
1679 93 : CSLFetchNameValue(papszOptions, "BLOCKXSIZE");
1680 93 : poDS->nBlockXSize = (pszBlockXSize) ? atoi(pszBlockXSize) : 256;
1681 : const char *pszBlockYSize =
1682 93 : CSLFetchNameValue(papszOptions, "BLOCKYSIZE");
1683 93 : poDS->nBlockYSize = (pszBlockYSize) ? atoi(pszBlockYSize) : 256;
1684 93 : poDS->bStats = CSLFetchBoolean(papszOptions, "STATS", FALSE);
1685 :
1686 : const char *pszTimestamp =
1687 93 : CSLFetchNameValue(papszOptions, "TILEDB_TIMESTAMP");
1688 93 : if (pszTimestamp != nullptr)
1689 2 : poDS->nTimestamp = std::strtoull(pszTimestamp, nullptr, 10);
1690 :
1691 : // set dimensions and attribute type for schema
1692 186 : poDS->m_schema.reset(
1693 93 : new tiledb::ArraySchema(*poDS->m_ctx, TILEDB_DENSE));
1694 93 : poDS->m_schema->set_tile_order(TILEDB_ROW_MAJOR);
1695 93 : poDS->m_schema->set_cell_order(TILEDB_ROW_MAJOR);
1696 :
1697 93 : poDS->m_filterList.reset(new tiledb::FilterList(*poDS->m_ctx));
1698 :
1699 93 : if (pszCompression != nullptr)
1700 : {
1701 0 : int nLevel = (pszCompressionLevel) ? atoi(pszCompressionLevel) : -1;
1702 0 : if (TileDBDataset::AddFilter(*(poDS->m_ctx.get()),
1703 0 : *(poDS->m_filterList.get()),
1704 0 : pszCompression, nLevel) == CE_None)
1705 : {
1706 0 : poDS->SetMetadataItem("COMPRESSION", pszCompression,
1707 0 : "IMAGE_STRUCTURE");
1708 0 : poDS->m_schema->set_coords_filter_list(*poDS->m_filterList);
1709 : }
1710 : }
1711 :
1712 186 : CPLString osAux;
1713 93 : const char *pszArrayName = CPLGetBasename(pszFilename);
1714 93 : osAux.Printf("%s.tdb", pszArrayName);
1715 :
1716 93 : poDS->SetPhysicalFilename(
1717 : CPLFormFilename(pszFilename, osAux.c_str(), nullptr));
1718 :
1719 : // Initialize PAM information.
1720 93 : poDS->SetDescription(pszFilename);
1721 :
1722 : // this driver enforces that all subdatasets are the same size
1723 186 : tiledb::Domain domain(*poDS->m_ctx);
1724 :
1725 : // Note the dimension bounds are inclusive and are expanded to the match
1726 : // the block size
1727 93 : poDS->nBlocksX = DIV_ROUND_UP(nXSize, poDS->nBlockXSize);
1728 93 : poDS->nBlocksY = DIV_ROUND_UP(nYSize, poDS->nBlockYSize);
1729 :
1730 93 : uint64_t w = (uint64_t)poDS->nBlocksX * poDS->nBlockXSize - 1;
1731 93 : uint64_t h = (uint64_t)poDS->nBlocksY * poDS->nBlockYSize - 1;
1732 :
1733 : auto d1 = tiledb::Dimension::create<uint64_t>(
1734 279 : *poDS->m_ctx, "X", {0, w}, uint64_t(poDS->nBlockXSize));
1735 : auto d2 = tiledb::Dimension::create<uint64_t>(
1736 279 : *poDS->m_ctx, "Y", {0, h}, uint64_t(poDS->nBlockYSize));
1737 :
1738 : {
1739 : // Only used for unit test purposes (to check ability of GDAL to read
1740 : // an arbitrary array)
1741 : const char *pszAttrName =
1742 93 : CPLGetConfigOption("TILEDB_ATTRIBUTE", TILEDB_VALUES);
1743 93 : if ((poDS->nBands == 0) || (poDS->eIndexMode == ATTRIBUTES))
1744 : {
1745 5 : poDS->AddDimensions(domain, pszAttrName, d2, d1, nullptr);
1746 : }
1747 : else
1748 : {
1749 : auto d3 = tiledb::Dimension::create<uint64_t>(
1750 264 : *poDS->m_ctx, "BANDS", {1, uint64_t(poDS->nBands)}, 1);
1751 88 : poDS->AddDimensions(domain, pszAttrName, d2, d1, &d3);
1752 : }
1753 : }
1754 :
1755 93 : poDS->m_schema->set_domain(domain).set_order(
1756 93 : {{TILEDB_ROW_MAJOR, TILEDB_ROW_MAJOR}});
1757 : ;
1758 :
1759 : // register additional attributes to the pixel value, these will be
1760 : // be reported as subdatasets on future reads
1761 186 : poDS->papszAttributes =
1762 93 : CSLFetchNameValueMultiple(papszOptions, "TILEDB_ATTRIBUTE");
1763 :
1764 103 : for (int i = 0; poDS->papszAttributes != nullptr &&
1765 6 : poDS->papszAttributes[i] != nullptr;
1766 : i++)
1767 : {
1768 : // modeling additional attributes as subdatasets
1769 4 : poDS->bHasSubDatasets = true;
1770 : // check each attribute is a GDAL source
1771 : std::unique_ptr<GDALDataset> poAttrDS(
1772 8 : GDALDataset::Open(poDS->papszAttributes[i], GA_ReadOnly));
1773 :
1774 4 : if (poAttrDS != nullptr)
1775 : {
1776 : const char *pszAttrName =
1777 4 : CPLGetBasename(poAttrDS->GetDescription());
1778 : // check each is co-registered
1779 : // candidate band
1780 4 : int nAttrBands = poAttrDS->GetRasterCount();
1781 4 : if (nAttrBands > 0)
1782 : {
1783 4 : GDALRasterBand *poAttrBand = poAttrDS->GetRasterBand(1);
1784 :
1785 4 : if ((poAttrBand->GetXSize() == poDS->nRasterXSize) &&
1786 8 : (poAttrBand->GetYSize() == poDS->nRasterYSize) &&
1787 4 : (poDS->nBands == nAttrBands))
1788 : {
1789 : // could check geotransform, but it is sufficient
1790 : // that cartesian dimensions are equal
1791 4 : poDS->lpoAttributeDS.push_back(std::move(poAttrDS));
1792 4 : poDS->CreateAttribute(poAttrBand->GetRasterDataType(),
1793 : pszAttrName, 1);
1794 : }
1795 : else
1796 : {
1797 0 : CPLError(
1798 : CE_Warning, CPLE_AppDefined,
1799 : "Skipping %s as it has a different dimension\n",
1800 0 : poDS->papszAttributes[i]);
1801 : }
1802 : }
1803 : else
1804 : {
1805 0 : CPLError(CE_Warning, CPLE_AppDefined,
1806 : "Skipping %s as it doesn't have any bands\n",
1807 0 : poDS->papszAttributes[i]);
1808 : }
1809 : }
1810 : else
1811 : {
1812 0 : CPLError(CE_Warning, CPLE_AppDefined,
1813 : "Skipping %s, not recognized as a GDAL dataset\n",
1814 0 : poDS->papszAttributes[i]);
1815 : }
1816 : }
1817 :
1818 93 : return poDS.release();
1819 : }
1820 0 : catch (const tiledb::TileDBError &e)
1821 : {
1822 0 : CPLError(CE_Failure, CPLE_AppDefined, "TileDB: %s", e.what());
1823 0 : return nullptr;
1824 : }
1825 : }
1826 :
1827 : /************************************************************************/
1828 : /* CopySubDatasets() */
1829 : /* */
1830 : /* Copy SubDatasets from src to a TileDBRasterDataset */
1831 : /* */
1832 : /************************************************************************/
1833 :
1834 1 : CPLErr TileDBRasterDataset::CopySubDatasets(GDALDataset *poSrcDS,
1835 : TileDBRasterDataset *poDstDS,
1836 : GDALProgressFunc pfnProgress,
1837 : void *pProgressData)
1838 :
1839 : {
1840 : try
1841 : {
1842 2 : std::vector<std::unique_ptr<GDALDataset>> apoDatasets;
1843 1 : poDstDS->bHasSubDatasets = true;
1844 1 : char **papszSrcSubDatasets = poSrcDS->GetMetadata("SUBDATASETS");
1845 1 : if (!papszSrcSubDatasets)
1846 0 : return CE_Failure;
1847 : const char *pszSubDSName =
1848 1 : CSLFetchNameValue(papszSrcSubDatasets, "SUBDATASET_1_NAME");
1849 1 : if (!pszSubDSName)
1850 0 : return CE_Failure;
1851 :
1852 : CPLStringList apszTokens(CSLTokenizeString2(
1853 2 : pszSubDSName, ":", CSLT_HONOURSTRINGS | CSLT_PRESERVEESCAPES));
1854 : // FIXME? this is tailored for HDF5-like subdataset names
1855 : // HDF5:foo.hdf5:attrname.
1856 1 : if (apszTokens.size() != 3)
1857 : {
1858 0 : CPLError(CE_Failure, CPLE_AppDefined,
1859 : "Cannot guess attribute name in %s", pszSubDSName);
1860 0 : return CE_Failure;
1861 : }
1862 :
1863 : std::unique_ptr<GDALDataset> poSubDataset(
1864 2 : GDALDataset::Open(pszSubDSName));
1865 2 : if (poSubDataset.get() == nullptr ||
1866 1 : poSubDataset->GetRasterCount() == 0)
1867 : {
1868 0 : return CE_Failure;
1869 : }
1870 :
1871 1 : uint64_t nSubXSize = poSubDataset->GetRasterXSize();
1872 1 : uint64_t nSubYSize = poSubDataset->GetRasterYSize();
1873 :
1874 1 : const char *pszAttrName = apszTokens[2];
1875 :
1876 1 : poDstDS->CreateAttribute(
1877 : poSubDataset->GetRasterBand(1)->GetRasterDataType(), pszAttrName,
1878 : poSubDataset->GetRasterCount());
1879 1 : apoDatasets.push_back(std::move(poSubDataset));
1880 :
1881 5 : for (int i = 0; papszSrcSubDatasets[i] != nullptr; i++)
1882 : {
1883 4 : if (STARTS_WITH_CI(papszSrcSubDatasets[i], "SUBDATASET_1_NAME=") ||
1884 3 : strstr(papszSrcSubDatasets[i], "_DESC="))
1885 : {
1886 3 : continue;
1887 : }
1888 1 : pszSubDSName = CPLParseNameValue(papszSrcSubDatasets[i], nullptr);
1889 : apszTokens = CSLTokenizeString2(
1890 1 : pszSubDSName, ":", CSLT_HONOURSTRINGS | CSLT_PRESERVEESCAPES);
1891 1 : if (apszTokens.size() != 3)
1892 : {
1893 0 : CPLError(CE_Failure, CPLE_AppDefined,
1894 : "Cannot guess attribute name in %s", pszSubDSName);
1895 0 : continue;
1896 : }
1897 :
1898 : std::unique_ptr<GDALDataset> poSubDS(
1899 2 : GDALDataset::Open(pszSubDSName));
1900 1 : if ((poSubDS != nullptr) && poSubDS->GetRasterCount() > 0)
1901 : {
1902 1 : GDALRasterBand *poBand = poSubDS->GetRasterBand(1);
1903 : int nBlockXSize, nBlockYSize;
1904 1 : poBand->GetBlockSize(&nBlockXSize, &nBlockYSize);
1905 :
1906 1 : if ((poSubDS->GetRasterXSize() != (int)nSubXSize) ||
1907 1 : (poSubDS->GetRasterYSize() != (int)nSubYSize) ||
1908 3 : (nBlockXSize != poDstDS->nBlockXSize) ||
1909 1 : (nBlockYSize != poDstDS->nBlockYSize))
1910 : {
1911 0 : CPLError(CE_Warning, CPLE_AppDefined,
1912 : "Sub-datasets must have the same dimension,"
1913 : " and block sizes, skipping %s\n",
1914 : pszSubDSName);
1915 : }
1916 : else
1917 : {
1918 1 : pszAttrName = apszTokens[2];
1919 1 : poDstDS->CreateAttribute(
1920 : poSubDS->GetRasterBand(1)->GetRasterDataType(),
1921 : pszAttrName, poSubDS->GetRasterCount());
1922 1 : apoDatasets.push_back(std::move(poSubDS));
1923 : }
1924 : }
1925 : else
1926 : {
1927 0 : CPLError(
1928 : CE_Warning, CPLE_AppDefined,
1929 : "Sub-datasets must be not null and contain data in bands,"
1930 : "skipping %s\n",
1931 : pszSubDSName);
1932 : }
1933 : }
1934 :
1935 1 : poDstDS->SetMetadata(poDstDS->papszSubDatasets, "SUBDATASETS");
1936 1 : tiledb::Array::create(poDstDS->GetDescription(), *poDstDS->m_schema);
1937 :
1938 1 : if (poDstDS->nTimestamp)
1939 : {
1940 0 : poDstDS->m_array.reset(new tiledb::Array(
1941 0 : *poDstDS->m_ctx, poDstDS->GetDescription(), TILEDB_WRITE,
1942 0 : tiledb::TemporalPolicy(tiledb::TimeTravel,
1943 0 : poDstDS->nTimestamp)));
1944 : }
1945 : else
1946 3 : poDstDS->m_array.reset(new tiledb::Array(
1947 3 : *poDstDS->m_ctx, poDstDS->GetDescription(), TILEDB_WRITE));
1948 :
1949 : /* -------------------------------------------------------- */
1950 : /* Report preliminary (0) progress. */
1951 : /* --------------------------------------------------------- */
1952 1 : if (!pfnProgress(0.0, nullptr, pProgressData))
1953 : {
1954 0 : CPLError(CE_Failure, CPLE_UserInterrupt,
1955 : "User terminated CreateCopy()");
1956 0 : return CE_Failure;
1957 : }
1958 :
1959 : // copy over subdatasets by block
1960 2 : tiledb::Query query(*poDstDS->m_ctx, *poDstDS->m_array);
1961 1 : query.set_layout(TILEDB_GLOBAL_ORDER);
1962 1 : int nTotalBlocks = poDstDS->nBlocksX * poDstDS->nBlocksY;
1963 :
1964 : // row-major
1965 6 : for (int j = 0; j < poDstDS->nBlocksY; ++j)
1966 : {
1967 55 : for (int i = 0; i < poDstDS->nBlocksX; ++i)
1968 : {
1969 50 : std::vector<std::unique_ptr<void, decltype(&VSIFree)>> aBlocks;
1970 : // have to write set all tiledb attributes on write
1971 50 : int iAttr = 0;
1972 150 : for (auto &poSubDS : apoDatasets)
1973 : {
1974 : GDALDataType eDT =
1975 100 : poSubDS->GetRasterBand(1)->GetRasterDataType();
1976 :
1977 200 : for (int b = 1; b <= poSubDS->GetRasterCount(); ++b)
1978 : {
1979 100 : int nBytes = GDALGetDataTypeSizeBytes(eDT);
1980 100 : int nValues = nBytes * poDstDS->nBlockXSize *
1981 100 : poDstDS->nBlockYSize;
1982 100 : void *pBlock = VSIMalloc(nBytes * nValues);
1983 100 : aBlocks.emplace_back(pBlock, &VSIFree);
1984 100 : GDALRasterBand *poBand = poSubDS->GetRasterBand(b);
1985 100 : if (poBand->ReadBlock(i, j, pBlock) == CE_None)
1986 : {
1987 100 : SetBuffer(
1988 : &query, eDT,
1989 200 : poDstDS->m_schema->attribute(iAttr++).name(),
1990 : pBlock,
1991 100 : poDstDS->nBlockXSize * poDstDS->nBlockYSize);
1992 : }
1993 : }
1994 : }
1995 :
1996 50 : if (poDstDS->bStats)
1997 0 : tiledb::Stats::enable();
1998 :
1999 50 : auto status = query.submit();
2000 :
2001 50 : if (poDstDS->bStats)
2002 : {
2003 0 : tiledb::Stats::dump(stdout);
2004 0 : tiledb::Stats::disable();
2005 : }
2006 :
2007 50 : if (status == tiledb::Query::Status::FAILED)
2008 : {
2009 0 : return CE_Failure;
2010 : }
2011 :
2012 50 : int nBlocks = ((j + 1) * poDstDS->nBlocksX);
2013 :
2014 50 : if (!pfnProgress(nBlocks / static_cast<double>(nTotalBlocks),
2015 : nullptr, pProgressData))
2016 : {
2017 0 : CPLError(CE_Failure, CPLE_UserInterrupt,
2018 : "User terminated CreateCopy()");
2019 0 : return CE_Failure;
2020 : }
2021 : }
2022 : }
2023 :
2024 1 : query.finalize();
2025 :
2026 1 : return CE_None;
2027 : }
2028 0 : catch (const tiledb::TileDBError &e)
2029 : {
2030 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
2031 0 : return CE_Failure;
2032 : }
2033 : }
2034 :
2035 : /************************************************************************/
2036 : /* Create() */
2037 : /************************************************************************/
2038 :
2039 92 : GDALDataset *TileDBRasterDataset::Create(const char *pszFilename, int nXSize,
2040 : int nYSize, int nBandsIn,
2041 : GDALDataType eType,
2042 : char **papszOptions)
2043 :
2044 : {
2045 184 : CPLString osArrayPath = TileDBDataset::VSI_to_tiledb_uri(pszFilename);
2046 :
2047 : std::unique_ptr<TileDBRasterDataset> poDS(TileDBRasterDataset::CreateLL(
2048 184 : osArrayPath, nXSize, nYSize, nBandsIn, eType, papszOptions));
2049 :
2050 92 : if (!poDS)
2051 0 : return nullptr;
2052 :
2053 92 : tiledb::Array::create(osArrayPath, *poDS->m_schema);
2054 :
2055 79 : if (poDS->nTimestamp)
2056 6 : poDS->m_array.reset(new tiledb::Array(
2057 2 : *poDS->m_ctx, osArrayPath, TILEDB_WRITE,
2058 6 : tiledb::TemporalPolicy(tiledb::TimeTravel, poDS->nTimestamp)));
2059 : else
2060 231 : poDS->m_array.reset(
2061 154 : new tiledb::Array(*poDS->m_ctx, osArrayPath, TILEDB_WRITE));
2062 :
2063 : const char *pszAttrName =
2064 79 : CPLGetConfigOption("TILEDB_ATTRIBUTE", TILEDB_VALUES);
2065 224 : for (int i = 0; i < poDS->nBands; i++)
2066 : {
2067 145 : if (poDS->eIndexMode == ATTRIBUTES)
2068 20 : poDS->SetBand(
2069 : i + 1, new TileDBRasterBand(
2070 10 : poDS.get(), i + 1,
2071 20 : TILEDB_VALUES + CPLString().Printf("_%i", i + 1)));
2072 : else
2073 405 : poDS->SetBand(i + 1,
2074 270 : new TileDBRasterBand(poDS.get(), i + 1, pszAttrName));
2075 : }
2076 :
2077 : // Only used for unit test purposes (to check ability of GDAL to read
2078 : // an arbitrary array)
2079 79 : if (CPLTestBool(CPLGetConfigOption("TILEDB_WRITE_IMAGE_STRUCTURE", "YES")))
2080 : {
2081 65 : char **papszImageStruct = nullptr;
2082 : papszImageStruct =
2083 65 : CSLAddNameValue(papszImageStruct, "NBITS",
2084 130 : CPLString().Printf("%d", poDS->nBitsPerSample));
2085 65 : papszImageStruct = CSLAddNameValue(
2086 : papszImageStruct, "DATA_TYPE",
2087 130 : CPLString().Printf("%s", GDALGetDataTypeName(poDS->eDataType)));
2088 : papszImageStruct =
2089 65 : CSLAddNameValue(papszImageStruct, "X_SIZE",
2090 130 : CPLString().Printf("%d", poDS->nRasterXSize));
2091 : papszImageStruct =
2092 65 : CSLAddNameValue(papszImageStruct, "Y_SIZE",
2093 130 : CPLString().Printf("%d", poDS->nRasterYSize));
2094 65 : papszImageStruct = CSLAddNameValue(papszImageStruct, "INTERLEAVE",
2095 65 : index_type_name(poDS->eIndexMode));
2096 65 : papszImageStruct = CSLAddNameValue(papszImageStruct, "DATASET_TYPE",
2097 : RASTER_DATASET_TYPE);
2098 :
2099 65 : if (poDS->lpoAttributeDS.size() > 0)
2100 : {
2101 2 : int i = 0;
2102 6 : for (auto const &poAttrDS : poDS->lpoAttributeDS)
2103 : {
2104 8 : papszImageStruct = CSLAddNameValue(
2105 : papszImageStruct,
2106 8 : CPLString().Printf("TILEDB_ATTRIBUTE_%i", ++i),
2107 4 : CPLGetBasename(poAttrDS->GetDescription()));
2108 : }
2109 : }
2110 65 : poDS->SetMetadata(papszImageStruct, "IMAGE_STRUCTURE");
2111 :
2112 65 : CSLDestroy(papszImageStruct);
2113 : }
2114 :
2115 79 : return poDS.release();
2116 : }
2117 :
2118 : /************************************************************************/
2119 : /* CreateCopy() */
2120 : /************************************************************************/
2121 :
2122 58 : GDALDataset *TileDBRasterDataset::CreateCopy(const char *pszFilename,
2123 : GDALDataset *poSrcDS, int bStrict,
2124 : char **papszOptions,
2125 : GDALProgressFunc pfnProgress,
2126 : void *pProgressData)
2127 :
2128 : {
2129 58 : char **papszCopyOptions = CSLDuplicate(papszOptions);
2130 116 : CPLString osArrayPath = TileDBDataset::VSI_to_tiledb_uri(pszFilename);
2131 :
2132 58 : std::unique_ptr<TileDBRasterDataset> poDstDS;
2133 :
2134 58 : if (CSLFetchNameValue(papszOptions, "APPEND_SUBDATASET"))
2135 : {
2136 : // TileDB schemas are fixed
2137 0 : CPLError(CE_Failure, CPLE_NotSupported,
2138 : "TileDB driver does not support "
2139 : "appending to an existing schema.");
2140 0 : CSLDestroy(papszCopyOptions);
2141 0 : return nullptr;
2142 : }
2143 :
2144 58 : char **papszSrcSubDatasets = poSrcDS->GetMetadata("SUBDATASETS");
2145 :
2146 58 : if (papszSrcSubDatasets == nullptr)
2147 : {
2148 57 : const int nBands = poSrcDS->GetRasterCount();
2149 :
2150 57 : if (nBands > 0)
2151 : {
2152 57 : GDALRasterBand *poBand = poSrcDS->GetRasterBand(1);
2153 57 : GDALDataType eType = poBand->GetRasterDataType();
2154 :
2155 87 : for (int i = 2; i <= nBands; ++i)
2156 : {
2157 30 : if (eType != poSrcDS->GetRasterBand(i)->GetRasterDataType())
2158 : {
2159 0 : CPLError(CE_Failure, CPLE_NotSupported,
2160 : "TileDB driver does not support "
2161 : "source dataset with different band data types.");
2162 0 : CSLDestroy(papszCopyOptions);
2163 0 : return nullptr;
2164 : }
2165 : }
2166 :
2167 44 : poDstDS.reset(
2168 57 : static_cast<TileDBRasterDataset *>(TileDBRasterDataset::Create(
2169 : osArrayPath, poSrcDS->GetRasterXSize(),
2170 : poSrcDS->GetRasterYSize(), nBands, eType, papszOptions)));
2171 :
2172 44 : if (!poDstDS)
2173 : {
2174 0 : CSLDestroy(papszCopyOptions);
2175 0 : return nullptr;
2176 : }
2177 :
2178 : CPLErr eErr =
2179 44 : GDALDatasetCopyWholeRaster(poSrcDS, poDstDS.get(), papszOptions,
2180 : pfnProgress, pProgressData);
2181 :
2182 44 : if (eErr != CE_None)
2183 : {
2184 0 : CPLError(eErr, CPLE_AppDefined,
2185 : "Error copying raster to TileDB.");
2186 : }
2187 : }
2188 : else
2189 : {
2190 0 : CPLError(CE_Failure, CPLE_NotSupported,
2191 : "TileDB driver does not support "
2192 : "source dataset with zero bands.");
2193 : }
2194 : }
2195 : else
2196 : {
2197 1 : if (bStrict)
2198 : {
2199 0 : CPLError(CE_Failure, CPLE_NotSupported,
2200 : "TileDB driver does not support copying "
2201 : "subdatasets in strict mode.");
2202 0 : CSLDestroy(papszCopyOptions);
2203 0 : return nullptr;
2204 : }
2205 :
2206 2 : if (CSLFetchNameValue(papszOptions, "BLOCKXSIZE") ||
2207 1 : CSLFetchNameValue(papszOptions, "BLOCKYSIZE"))
2208 : {
2209 0 : CPLError(CE_Failure, CPLE_NotSupported,
2210 : "Changing block size is not supported when copying "
2211 : "subdatasets.");
2212 0 : CSLDestroy(papszCopyOptions);
2213 0 : return nullptr;
2214 : }
2215 :
2216 1 : const int nSubDatasetCount = CSLCount(papszSrcSubDatasets) / 2;
2217 : const int nMaxFiles =
2218 1 : atoi(CPLGetConfigOption("GDAL_READDIR_LIMIT_ON_OPEN", "1000"));
2219 :
2220 1 : if (nSubDatasetCount <= nMaxFiles)
2221 : {
2222 : const char *pszSource =
2223 1 : CSLFetchNameValue(papszSrcSubDatasets, "SUBDATASET_1_NAME");
2224 1 : if (pszSource)
2225 : {
2226 : std::unique_ptr<GDALDataset> poSubDataset(
2227 2 : GDALDataset::Open(pszSource));
2228 1 : if (poSubDataset && poSubDataset->GetRasterCount() > 0)
2229 : {
2230 1 : GDALRasterBand *poBand = poSubDataset->GetRasterBand(1);
2231 :
2232 1 : TileDBRasterDataset::SetBlockSize(poBand, papszCopyOptions);
2233 1 : poDstDS.reset(TileDBRasterDataset::CreateLL(
2234 : osArrayPath, poBand->GetXSize(), poBand->GetYSize(), 0,
2235 : poBand->GetRasterDataType(), papszCopyOptions));
2236 :
2237 1 : if (poDstDS && TileDBRasterDataset::CopySubDatasets(
2238 : poSrcDS, poDstDS.get(), pfnProgress,
2239 1 : pProgressData) != CE_None)
2240 : {
2241 0 : poDstDS.reset();
2242 0 : CPLError(CE_Failure, CPLE_AppDefined,
2243 : "Unable to copy subdatasets.");
2244 : }
2245 : }
2246 : }
2247 : }
2248 : else
2249 : {
2250 0 : CPLError(CE_Failure, CPLE_AppDefined,
2251 : "Please increase GDAL_READDIR_LIMIT_ON_OPEN variable.");
2252 : }
2253 : }
2254 :
2255 45 : CSLDestroy(papszCopyOptions);
2256 :
2257 : // TODO Supporting mask bands is a possible future task
2258 45 : if (poDstDS != nullptr)
2259 : {
2260 45 : int nCloneFlags = GCIF_PAM_DEFAULT & ~GCIF_MASK;
2261 45 : poDstDS->CloneInfo(poSrcDS, nCloneFlags);
2262 :
2263 45 : if (poDstDS->eIndexMode == ATTRIBUTES)
2264 : {
2265 5 : poDstDS->FlushCache(false);
2266 : }
2267 :
2268 45 : poDstDS->m_array->close();
2269 45 : poDstDS->eAccess = GA_ReadOnly;
2270 45 : poDstDS->m_array->open(TILEDB_READ);
2271 :
2272 45 : return poDstDS.release();
2273 : }
2274 0 : return nullptr;
2275 : }
|
