Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: GDAL TileDB Driver
4 : * Purpose: Implement GDAL TileDB multidimensional support based on https://www.tiledb.io
5 : * Author: TileDB, Inc
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2023, TileDB, Inc
9 : *
10 : * SPDX-License-Identifier: MIT
11 : ****************************************************************************/
12 :
13 : #include "tiledbmultidim.h"
14 :
15 : #include <algorithm>
16 : #include <limits>
17 :
18 : /************************************************************************/
19 : /* TileDBArray::TileDBArray() */
20 : /************************************************************************/
21 :
22 75 : TileDBArray::TileDBArray(
23 : const std::shared_ptr<TileDBSharedResource> &poSharedResource,
24 : const std::string &osParentName, const std::string &osName,
25 : const std::vector<std::shared_ptr<GDALDimension>> &aoDims,
26 75 : const GDALExtendedDataType &oType, const std::string &osPath)
27 : : GDALAbstractMDArray(osParentName, osName),
28 : GDALMDArray(osParentName, osName), m_poSharedResource(poSharedResource),
29 : m_aoDims(aoDims), m_oType(oType), m_osPath(osPath),
30 75 : m_bStats(poSharedResource->GetDumpStats())
31 : {
32 75 : }
33 :
34 : /************************************************************************/
35 : /* TileDBArray::Create() */
36 : /************************************************************************/
37 :
38 75 : /*static*/ std::shared_ptr<TileDBArray> TileDBArray::Create(
39 : const std::shared_ptr<TileDBSharedResource> &poSharedResource,
40 : const std::string &osParentName, const std::string &osName,
41 : const std::vector<std::shared_ptr<GDALDimension>> &aoDims,
42 : const GDALExtendedDataType &oType, const std::string &osPath)
43 : {
44 : auto poArray = std::shared_ptr<TileDBArray>(new TileDBArray(
45 75 : poSharedResource, osParentName, osName, aoDims, oType, osPath));
46 75 : poArray->SetSelf(poArray);
47 75 : return poArray;
48 : }
49 :
50 : /************************************************************************/
51 : /* TileDBArray::~TileDBArray() */
52 : /************************************************************************/
53 :
54 150 : TileDBArray::~TileDBArray()
55 : {
56 75 : if (!m_bFinalized)
57 5 : Finalize();
58 150 : }
59 :
60 : /************************************************************************/
61 : /* BuildDimensionLabelName() */
62 : /************************************************************************/
63 :
64 : static std::string
65 77 : BuildDimensionLabelName(const std::shared_ptr<GDALDimension> &poDim)
66 : {
67 77 : return poDim->GetName() + "_label";
68 : }
69 :
70 : /************************************************************************/
71 : /* TileDBDataTypeToGDALDataType() */
72 : /************************************************************************/
73 :
74 : /*static*/ GDALDataType
75 71 : TileDBArray::TileDBDataTypeToGDALDataType(tiledb_datatype_t tiledb_dt)
76 : {
77 71 : GDALDataType eDT = GDT_Unknown;
78 71 : switch (tiledb_dt)
79 : {
80 7 : case TILEDB_UINT8:
81 7 : eDT = GDT_Byte;
82 7 : break;
83 :
84 1 : case TILEDB_INT8:
85 1 : eDT = GDT_Int8;
86 1 : break;
87 :
88 1 : case TILEDB_UINT16:
89 1 : eDT = GDT_UInt16;
90 1 : break;
91 :
92 2 : case TILEDB_INT16:
93 2 : eDT = GDT_Int16;
94 2 : break;
95 :
96 1 : case TILEDB_UINT32:
97 1 : eDT = GDT_UInt32;
98 1 : break;
99 :
100 7 : case TILEDB_INT32:
101 7 : eDT = GDT_Int32;
102 7 : break;
103 :
104 1 : case TILEDB_UINT64:
105 1 : eDT = GDT_UInt64;
106 1 : break;
107 :
108 1 : case TILEDB_INT64:
109 1 : eDT = GDT_Int64;
110 1 : break;
111 :
112 12 : case TILEDB_FLOAT32:
113 12 : eDT = GDT_Float32;
114 12 : break;
115 :
116 37 : case TILEDB_FLOAT64:
117 37 : eDT = GDT_Float64;
118 37 : break;
119 :
120 1 : case TILEDB_CHAR:
121 : case TILEDB_STRING_ASCII:
122 : case TILEDB_STRING_UTF8:
123 : case TILEDB_STRING_UTF16:
124 : case TILEDB_STRING_UTF32:
125 : case TILEDB_STRING_UCS2:
126 : case TILEDB_STRING_UCS4:
127 : case TILEDB_ANY:
128 : case TILEDB_DATETIME_YEAR:
129 : case TILEDB_DATETIME_MONTH:
130 : case TILEDB_DATETIME_WEEK:
131 : case TILEDB_DATETIME_DAY:
132 : case TILEDB_DATETIME_HR:
133 : case TILEDB_DATETIME_MIN:
134 : case TILEDB_DATETIME_SEC:
135 : case TILEDB_DATETIME_MS:
136 : case TILEDB_DATETIME_US:
137 : case TILEDB_DATETIME_NS:
138 : case TILEDB_DATETIME_PS:
139 : case TILEDB_DATETIME_FS:
140 : case TILEDB_DATETIME_AS:
141 : case TILEDB_TIME_HR:
142 : case TILEDB_TIME_MIN:
143 : case TILEDB_TIME_SEC:
144 : case TILEDB_TIME_MS:
145 : case TILEDB_TIME_US:
146 : case TILEDB_TIME_NS:
147 : case TILEDB_TIME_PS:
148 : case TILEDB_TIME_FS:
149 : case TILEDB_TIME_AS:
150 : case TILEDB_BLOB:
151 : case TILEDB_BOOL:
152 : #ifdef HAS_TILEDB_GEOM_WKB_WKT
153 : case TILEDB_GEOM_WKB:
154 : case TILEDB_GEOM_WKT:
155 : #endif
156 : {
157 1 : break;
158 : }
159 : }
160 71 : return eDT;
161 : }
162 :
163 : /************************************************************************/
164 : /* TileDBArray::Finalize() */
165 : /************************************************************************/
166 :
167 30 : bool TileDBArray::Finalize() const
168 : {
169 30 : if (m_bFinalized)
170 0 : return m_poTileDBArray != nullptr;
171 :
172 30 : m_bFinalized = true;
173 :
174 30 : CPLAssert(m_poSchema);
175 30 : CPLAssert(m_poAttr);
176 :
177 : try
178 : {
179 : // TODO: set nodata as fill_value
180 :
181 30 : m_poSchema->add_attribute(*(m_poAttr.get()));
182 :
183 30 : tiledb::Array::create(m_osPath, *m_poSchema);
184 :
185 30 : bool bAdded = false;
186 60 : auto poGroup = m_poParent.lock();
187 30 : if (!poGroup)
188 : {
189 : // Temporarily instantiate a TileDBGroup to call AddMember() on it
190 15 : poGroup = TileDBGroup::OpenFromDisk(
191 5 : m_poSharedResource,
192 10 : /* osParentName = */ std::string(),
193 10 : CPLGetFilename(m_osParentPath.c_str()), m_osParentPath);
194 : }
195 30 : if (poGroup)
196 : {
197 30 : bAdded = poGroup->AddMember(m_osPath, m_osName);
198 : }
199 30 : if (!bAdded)
200 : {
201 0 : CPLError(CE_Failure, CPLE_AppDefined,
202 : "Could not add array %s as a member of group %s",
203 0 : m_osName.c_str(), m_osParentPath.c_str());
204 : }
205 :
206 30 : auto &ctx = m_poSharedResource->GetCtx();
207 : m_poTileDBArray =
208 30 : std::make_unique<tiledb::Array>(ctx, m_osPath, TILEDB_READ);
209 30 : if (m_nTimestamp > 0)
210 0 : m_poTileDBArray->set_open_timestamp_end(m_nTimestamp);
211 : m_poSchema =
212 30 : std::make_unique<tiledb::ArraySchema>(m_poTileDBArray->schema());
213 30 : m_poAttr.reset();
214 :
215 : // Write dimension label values
216 71 : for (const auto &poDim : m_aoDims)
217 : {
218 82 : auto poVar = poDim->GetIndexingVariable();
219 41 : if (poVar)
220 : {
221 16 : const std::string osLabelName(BuildDimensionLabelName(poDim));
222 8 : if (tiledb::ArraySchemaExperimental::has_dimension_label(
223 8 : ctx, *(m_poSchema.get()), osLabelName))
224 : {
225 : auto label =
226 : tiledb::ArraySchemaExperimental::dimension_label(
227 16 : ctx, *(m_poSchema.get()), osLabelName);
228 16 : tiledb::Array labelArray(ctx, label.uri(), TILEDB_WRITE);
229 16 : auto label_attr = labelArray.schema().attribute(0);
230 : const auto eDT =
231 8 : TileDBDataTypeToGDALDataType(label_attr.type());
232 8 : if (eDT != GDT_Unknown)
233 : {
234 16 : std::vector<GByte> abyVals;
235 8 : abyVals.resize(static_cast<size_t>(
236 8 : poVar->GetDimensions()[0]->GetSize() *
237 8 : GDALGetDataTypeSizeBytes(eDT)));
238 8 : GUInt64 anStart[1] = {0};
239 : size_t anCount[1] = {static_cast<size_t>(
240 8 : poVar->GetDimensions()[0]->GetSize())};
241 16 : if (poVar->Read(anStart, anCount, nullptr, nullptr,
242 16 : GDALExtendedDataType::Create(eDT),
243 8 : abyVals.data()))
244 : {
245 16 : tiledb::Query query(ctx, labelArray);
246 : query.set_data_buffer(
247 16 : label_attr.name(),
248 8 : static_cast<void *>(abyVals.data()),
249 16 : anCount[0]);
250 8 : if (query.submit() !=
251 : tiledb::Query::Status::COMPLETE)
252 : {
253 0 : CPLError(CE_Failure, CPLE_AppDefined,
254 : "Could not write values for dimension "
255 : "label %s",
256 : osLabelName.c_str());
257 : }
258 :
259 8 : if (!poDim->GetType().empty())
260 : {
261 6 : labelArray.put_metadata(
262 : DIM_TYPE_ATTRIBUTE_NAME, TILEDB_STRING_UTF8,
263 : static_cast<uint32_t>(
264 6 : poDim->GetType().size()),
265 6 : poDim->GetType().c_str());
266 : }
267 :
268 8 : if (!poDim->GetDirection().empty())
269 : {
270 4 : labelArray.put_metadata(
271 : DIM_DIRECTION_ATTRIBUTE_NAME,
272 : TILEDB_STRING_UTF8,
273 : static_cast<uint32_t>(
274 4 : poDim->GetDirection().size()),
275 4 : poDim->GetDirection().c_str());
276 : }
277 : }
278 : }
279 : }
280 : }
281 : }
282 : }
283 0 : catch (const std::exception &e)
284 : {
285 0 : CPLError(CE_Failure, CPLE_AppDefined,
286 0 : "Array %s creation failed with: %s", m_osName.c_str(),
287 0 : e.what());
288 0 : return false;
289 : }
290 :
291 30 : return true;
292 : }
293 :
294 : /************************************************************************/
295 : /* TileDBArray::OpenFromDisk() */
296 : /************************************************************************/
297 :
298 : /* static */
299 45 : std::shared_ptr<TileDBArray> TileDBArray::OpenFromDisk(
300 : const std::shared_ptr<TileDBSharedResource> &poSharedResource,
301 : const std::shared_ptr<GDALGroup> &poParent, const std::string &osParentName,
302 : const std::string &osName, const std::string &osAttributeName,
303 : const std::string &osPath, CSLConstList papszOptions)
304 : {
305 : try
306 : {
307 45 : auto &ctx = poSharedResource->GetCtx();
308 45 : uint64_t nTimestamp = poSharedResource->GetTimestamp();
309 : const char *pszTimestamp =
310 45 : CSLFetchNameValue(papszOptions, "TILEDB_TIMESTAMP");
311 45 : if (pszTimestamp)
312 0 : nTimestamp = std::strtoull(pszTimestamp, nullptr, 10);
313 :
314 : auto poTileDBArray =
315 90 : std::make_unique<tiledb::Array>(ctx, osPath, TILEDB_READ);
316 45 : if (nTimestamp > 0)
317 0 : poTileDBArray->set_open_timestamp_end(nTimestamp);
318 :
319 90 : auto schema = poTileDBArray->schema();
320 :
321 45 : if (schema.attribute_num() != 1 && osAttributeName.empty())
322 : {
323 0 : CPLError(CE_Failure, CPLE_NotSupported,
324 : "Array %s has %u attributes. "
325 : "osAttributeName must be specified",
326 : osName.c_str(), schema.attribute_num());
327 0 : return nullptr;
328 : }
329 :
330 45 : const auto &attr = osAttributeName.empty()
331 : ? schema.attribute(0)
332 90 : : schema.attribute(osAttributeName);
333 45 : GDALDataType eDT = TileDBDataTypeToGDALDataType(attr.type());
334 45 : if (attr.type() == TILEDB_CHAR)
335 1 : eDT = GDT_Byte;
336 45 : if (eDT == GDT_Unknown)
337 : {
338 0 : const char *pszTypeName = "";
339 0 : tiledb_datatype_to_str(attr.type(), &pszTypeName);
340 0 : CPLError(CE_Failure, CPLE_NotSupported,
341 : "Array %s has type %s, which is unsupported",
342 : osName.c_str(), pszTypeName);
343 0 : return nullptr;
344 : }
345 :
346 45 : if (attr.variable_sized())
347 : {
348 0 : CPLError(CE_Failure, CPLE_NotSupported,
349 : "Variable sized attribute not supported");
350 0 : return nullptr;
351 : }
352 45 : if (attr.cell_val_num() == 2)
353 : {
354 4 : if (attr.type() == TILEDB_INT16)
355 1 : eDT = GDT_CInt16;
356 3 : else if (attr.type() == TILEDB_INT32)
357 1 : eDT = GDT_CInt32;
358 2 : else if (attr.type() == TILEDB_FLOAT32)
359 1 : eDT = GDT_CFloat32;
360 1 : else if (attr.type() == TILEDB_FLOAT64)
361 1 : eDT = GDT_CFloat64;
362 : else
363 : {
364 0 : const char *pszTypeName = "";
365 0 : tiledb_datatype_to_str(attr.type(), &pszTypeName);
366 0 : CPLError(CE_Failure, CPLE_NotSupported,
367 : "Attribute with number of values per cell = %u not "
368 : "supported for type %s",
369 : attr.cell_val_num(), pszTypeName);
370 0 : return nullptr;
371 : }
372 : }
373 41 : else if (attr.cell_val_num() != 1)
374 : {
375 0 : CPLError(
376 : CE_Failure, CPLE_NotSupported,
377 : "Attribute with number of values per cell = %u not supported",
378 : attr.cell_val_num());
379 0 : return nullptr;
380 : }
381 :
382 : // Compatibility with the 2D raster side: extract X_SIZE, Y_SIZE, SRS
383 : // and geotransform
384 45 : int nXSize = 0;
385 45 : int nYSize = 0;
386 45 : std::shared_ptr<OGRSpatialReference> poSRS;
387 45 : double adfGeoTransform[6] = {0};
388 45 : bool bHasGeoTransform = false;
389 : {
390 45 : tiledb_datatype_t value_type = TILEDB_ANY;
391 45 : uint32_t value_num = 0;
392 45 : const void *value = nullptr;
393 45 : poTileDBArray->get_metadata(GDAL_ATTRIBUTE_NAME, &value_type,
394 : &value_num, &value);
395 49 : if (value && value_num && value_type == TILEDB_UINT8 &&
396 4 : CPLIsUTF8(static_cast<const char *>(value), value_num))
397 : {
398 8 : std::string osXML;
399 4 : osXML.assign(static_cast<const char *>(value), value_num);
400 4 : CPLXMLNode *psRoot = CPLParseXMLString(osXML.c_str());
401 4 : if (psRoot)
402 : {
403 : const CPLXMLNode *psDataset =
404 4 : CPLGetXMLNode(psRoot, "=PAMDataset");
405 4 : if (psDataset)
406 : {
407 4 : for (const CPLXMLNode *psIter = psDataset->psChild;
408 20 : psIter; psIter = psIter->psNext)
409 : {
410 48 : if (psIter->eType == CXT_Element &&
411 24 : strcmp(psIter->pszValue, "Metadata") == 0 &&
412 8 : strcmp(CPLGetXMLValue(psIter, "domain", ""),
413 : "IMAGE_STRUCTURE") == 0)
414 : {
415 4 : for (const CPLXMLNode *psIter2 =
416 : psIter->psChild;
417 32 : psIter2; psIter2 = psIter2->psNext)
418 : {
419 80 : if (psIter2->eType == CXT_Element &&
420 52 : strcmp(psIter2->pszValue, "MDI") == 0 &&
421 24 : strcmp(
422 : CPLGetXMLValue(psIter2, "key", ""),
423 : "X_SIZE") == 0)
424 : {
425 4 : nXSize = atoi(CPLGetXMLValue(
426 : psIter2, nullptr, "0"));
427 : }
428 68 : else if (psIter2->eType == CXT_Element &&
429 20 : strcmp(psIter2->pszValue, "MDI") ==
430 44 : 0 &&
431 20 : strcmp(CPLGetXMLValue(psIter2,
432 : "key", ""),
433 : "Y_SIZE") == 0)
434 : {
435 4 : nYSize = atoi(CPLGetXMLValue(
436 : psIter2, nullptr, "0"));
437 : }
438 : }
439 : }
440 : }
441 :
442 : const char *pszSRS =
443 4 : CPLGetXMLValue(psDataset, "SRS", nullptr);
444 4 : if (pszSRS)
445 : {
446 4 : poSRS = std::make_shared<OGRSpatialReference>();
447 4 : poSRS->SetAxisMappingStrategy(
448 : OAMS_TRADITIONAL_GIS_ORDER);
449 4 : if (poSRS->importFromWkt(pszSRS) != OGRERR_NONE)
450 : {
451 0 : poSRS.reset();
452 : }
453 : }
454 :
455 : const char *pszGeoTransform =
456 4 : CPLGetXMLValue(psDataset, "GeoTransform", nullptr);
457 4 : if (pszGeoTransform)
458 : {
459 : const CPLStringList aosTokens(
460 8 : CSLTokenizeString2(pszGeoTransform, ", ", 0));
461 4 : if (aosTokens.size() == 6)
462 : {
463 4 : bHasGeoTransform = true;
464 28 : for (int i = 0; i < 6; ++i)
465 24 : adfGeoTransform[i] = CPLAtof(aosTokens[i]);
466 : }
467 : }
468 : }
469 :
470 4 : CPLDestroyXMLNode(psRoot);
471 : }
472 : }
473 : }
474 :
475 : // Read CRS from _CRS attribute otherwise
476 45 : if (!poSRS)
477 : {
478 41 : tiledb_datatype_t value_type = TILEDB_ANY;
479 41 : uint32_t value_num = 0;
480 41 : const void *value = nullptr;
481 41 : poTileDBArray->get_metadata(CRS_ATTRIBUTE_NAME, &value_type,
482 : &value_num, &value);
483 41 : if (value && value_num &&
484 2 : (value_type == TILEDB_STRING_ASCII ||
485 2 : value_type == TILEDB_STRING_UTF8))
486 : {
487 4 : std::string osStr;
488 2 : osStr.assign(static_cast<const char *>(value), value_num);
489 2 : poSRS = std::make_shared<OGRSpatialReference>();
490 2 : poSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
491 2 : if (poSRS->SetFromUserInput(
492 : osStr.c_str(),
493 : OGRSpatialReference::
494 2 : SET_FROM_USER_INPUT_LIMITATIONS_get()) !=
495 : OGRERR_NONE)
496 : {
497 0 : poSRS.reset();
498 : }
499 : }
500 : }
501 :
502 : // Read unit
503 90 : std::string osUnit;
504 : {
505 45 : tiledb_datatype_t value_type = TILEDB_ANY;
506 45 : uint32_t value_num = 0;
507 45 : const void *value = nullptr;
508 45 : poTileDBArray->get_metadata(UNIT_ATTRIBUTE_NAME, &value_type,
509 : &value_num, &value);
510 45 : if (value && value_num &&
511 5 : (value_type == TILEDB_STRING_ASCII ||
512 5 : value_type == TILEDB_STRING_UTF8))
513 : {
514 5 : osUnit.assign(static_cast<const char *>(value), value_num);
515 : }
516 : }
517 :
518 : // Read dimensions
519 90 : std::vector<std::shared_ptr<GDALDimension>> aoDims;
520 90 : const auto dims = schema.domain().dimensions();
521 90 : std::vector<GUInt64> anBlockSize;
522 90 : std::vector<uint64_t> anStartDimOffset;
523 : const std::string osArrayFullName(
524 90 : (osParentName == "/" ? std::string() : osParentName) + "/" +
525 90 : osName);
526 106 : for (size_t i = 0; i < dims.size(); ++i)
527 : {
528 61 : const auto &dim = dims[i];
529 61 : if (dim.type() != TILEDB_UINT64)
530 : {
531 0 : const char *pszTypeName = "";
532 0 : tiledb_datatype_to_str(dim.type(), &pszTypeName);
533 0 : CPLError(CE_Failure, CPLE_NotSupported,
534 : "Dimension %s of array %s has type %s, which is "
535 : "unsupported. Only UInt64 is supported",
536 0 : dim.name().c_str(), osName.c_str(), pszTypeName);
537 0 : return nullptr;
538 : }
539 61 : const auto domain = dim.domain<uint64_t>();
540 61 : anStartDimOffset.push_back(domain.first);
541 118 : const uint64_t nSize = (i + 2 == dims.size() && nYSize > 0) ? nYSize
542 102 : : (i + 1 == dims.size() && nXSize > 0)
543 61 : ? nXSize
544 53 : : domain.second - domain.first + 1;
545 122 : std::string osType;
546 122 : std::string osDirection;
547 : auto poDim = std::make_shared<TileDBDimension>(
548 122 : osArrayFullName, dim.name(), osType, osDirection, nSize);
549 :
550 61 : const std::string osLabelName(BuildDimensionLabelName(poDim));
551 61 : if (tiledb::ArraySchemaExperimental::has_dimension_label(
552 : ctx, schema, osLabelName))
553 : {
554 : auto label = tiledb::ArraySchemaExperimental::dimension_label(
555 16 : ctx, schema, osLabelName);
556 : auto poIndexingVar = OpenFromDisk(
557 : poSharedResource, nullptr, osArrayFullName,
558 8 : poDim->GetName(),
559 16 : /* osAttributeName = */ std::string(), label.uri(),
560 24 : /* papszOptions= */ nullptr);
561 8 : if (poIndexingVar)
562 : {
563 : auto poAttr =
564 16 : poIndexingVar->GetAttribute(DIM_TYPE_ATTRIBUTE_NAME);
565 14 : if (poAttr &&
566 14 : poAttr->GetDataType().GetClass() == GEDTC_STRING)
567 : {
568 6 : const char *pszVal = poAttr->ReadAsString();
569 6 : if (pszVal)
570 6 : osType = pszVal;
571 : }
572 :
573 16 : poAttr = poIndexingVar->GetAttribute(
574 8 : DIM_DIRECTION_ATTRIBUTE_NAME);
575 12 : if (poAttr &&
576 12 : poAttr->GetDataType().GetClass() == GEDTC_STRING)
577 : {
578 4 : const char *pszVal = poAttr->ReadAsString();
579 4 : if (pszVal)
580 4 : osDirection = pszVal;
581 : }
582 :
583 8 : if (!osType.empty() || !osDirection.empty())
584 : {
585 : // Recreate dimension with type and/or direction info
586 12 : poDim = std::make_shared<TileDBDimension>(
587 12 : osArrayFullName, dim.name(), osType, osDirection,
588 6 : nSize);
589 : }
590 :
591 8 : poDim->SetIndexingVariableOneTime(poIndexingVar);
592 : }
593 : }
594 79 : if (bHasGeoTransform && !poDim->GetIndexingVariable() &&
595 87 : i + 2 >= dims.size() && adfGeoTransform[2] == 0 &&
596 8 : adfGeoTransform[4] == 0)
597 : {
598 : // Recreate dimension with type and/or direction info
599 8 : if (i + 2 == dims.size())
600 : {
601 4 : osType = GDAL_DIM_TYPE_HORIZONTAL_Y;
602 4 : osDirection = "NORTH";
603 : }
604 : else /* if( i + 1 == dims.size()) */
605 : {
606 4 : osType = GDAL_DIM_TYPE_HORIZONTAL_X;
607 4 : osDirection = "EAST";
608 : }
609 16 : poDim = std::make_shared<TileDBDimension>(
610 24 : osArrayFullName, dim.name(), osType, osDirection, nSize);
611 : // Do not create indexing variable with poDim, otherwise
612 : // both dimension and indexing variable will have a shared_ptr
613 : // to each other, causing memory leak
614 : auto poDimTmp = std::make_shared<GDALDimension>(
615 16 : std::string(), dim.name(), /* osType = */ std::string(),
616 16 : /* osDirection = */ std::string(), nSize);
617 : const double dfStart =
618 8 : (i + 2 == dims.size())
619 8 : ? adfGeoTransform[3] + adfGeoTransform[5] / 2
620 4 : : adfGeoTransform[0] + adfGeoTransform[1] / 2;
621 8 : const double dfStep = (i + 2 == dims.size())
622 8 : ? adfGeoTransform[5]
623 8 : : adfGeoTransform[1];
624 16 : poDim->SetIndexingVariableOneTime(
625 16 : GDALMDArrayRegularlySpaced::Create(
626 8 : osArrayFullName, poDim->GetName(), poDimTmp, dfStart,
627 : dfStep, 0));
628 : }
629 :
630 61 : if (poParent && dims.size() >= 2)
631 : {
632 69 : for (const auto &osOtherArray : poParent->GetMDArrayNames())
633 : {
634 46 : if (osOtherArray != osName)
635 : {
636 23 : auto poOtherArray = poParent->OpenMDArray(osOtherArray);
637 44 : if (poOtherArray &&
638 40 : poOtherArray->GetDimensionCount() == 1 &&
639 19 : poOtherArray->GetDataType().GetClass() ==
640 44 : GEDTC_NUMERIC &&
641 42 : poOtherArray->GetAttribute(
642 42 : std::string("__tiledb_attr.")
643 19 : .append(poDim->GetName())
644 38 : .append(".data.standard_name")))
645 : {
646 2 : if (dim.name() == "x")
647 : {
648 1 : osType = GDAL_DIM_TYPE_HORIZONTAL_X;
649 1 : osDirection = "EAST";
650 : }
651 1 : else if (dim.name() == "y")
652 : {
653 1 : osType = GDAL_DIM_TYPE_HORIZONTAL_Y;
654 1 : osDirection = "NORTH";
655 : }
656 2 : if (!osType.empty())
657 : {
658 4 : poDim = std::make_shared<TileDBDimension>(
659 4 : osArrayFullName, dim.name(), osType,
660 2 : osDirection, nSize);
661 : }
662 2 : poDim->SetIndexingVariableOneTime(poOtherArray);
663 2 : break;
664 : }
665 : }
666 : }
667 : }
668 :
669 61 : aoDims.emplace_back(std::move(poDim));
670 61 : anBlockSize.push_back(dim.tile_extent<uint64_t>());
671 : }
672 :
673 90 : GDALExtendedDataType oType = GDALExtendedDataType::Create(eDT);
674 : auto poArray = Create(poSharedResource, osParentName, osName, aoDims,
675 90 : oType, osPath);
676 45 : poArray->m_poTileDBArray = std::move(poTileDBArray);
677 90 : poArray->m_poSchema = std::make_unique<tiledb::ArraySchema>(
678 135 : poArray->m_poTileDBArray->schema());
679 45 : poArray->m_anBlockSize = std::move(anBlockSize);
680 45 : poArray->m_anStartDimOffset = std::move(anStartDimOffset);
681 45 : poArray->m_osAttrName = attr.name();
682 45 : poArray->m_osUnit = std::move(osUnit);
683 45 : poArray->m_nTimestamp = nTimestamp;
684 :
685 : // Try to get SRS from CF-1 conventions, if dataset has been generated
686 : // with https://github.com/TileDB-Inc/TileDB-CF-Py
687 45 : if (poParent && !poSRS)
688 : {
689 : const auto ENDS_WITH_CI =
690 9 : [](const char *pszStr, const char *pszNeedle)
691 : {
692 9 : const size_t nLenStr = strlen(pszStr);
693 9 : const size_t nLenNeedle = strlen(pszNeedle);
694 18 : return nLenStr >= nLenNeedle &&
695 9 : memcmp(pszStr + (nLenStr - nLenNeedle), pszNeedle,
696 9 : nLenNeedle) == 0;
697 : };
698 :
699 : const auto GetSRSFromGridMappingArray =
700 1 : [](const std::shared_ptr<GDALMDArray> &poOtherArray,
701 : const std::string &osGMPrefix)
702 : {
703 2 : CPLStringList aosGridMappingKeyValues;
704 11 : for (const auto &poGMAttr : poOtherArray->GetAttributes())
705 : {
706 10 : if (STARTS_WITH(poGMAttr->GetName().c_str(),
707 : osGMPrefix.c_str()))
708 : {
709 : const std::string osKey =
710 20 : poGMAttr->GetName().c_str() + osGMPrefix.size();
711 10 : if (poGMAttr->GetDataType().GetClass() == GEDTC_STRING)
712 : {
713 2 : const char *pszValue = poGMAttr->ReadAsString();
714 2 : if (pszValue)
715 : aosGridMappingKeyValues.AddNameValue(
716 2 : osKey.c_str(), pszValue);
717 : }
718 8 : else if (poGMAttr->GetDataType().GetClass() ==
719 : GEDTC_NUMERIC)
720 : {
721 : const auto aosValues =
722 16 : poGMAttr->ReadAsDoubleArray();
723 16 : std::string osVal;
724 17 : for (double dfVal : aosValues)
725 : {
726 9 : if (!osVal.empty())
727 1 : osVal += ',';
728 9 : osVal += CPLSPrintf("%.17g", dfVal);
729 : }
730 : aosGridMappingKeyValues.AddNameValue(osKey.c_str(),
731 8 : osVal.c_str());
732 : }
733 : }
734 : }
735 1 : auto l_poSRS = std::make_shared<OGRSpatialReference>();
736 1 : l_poSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
737 1 : if (l_poSRS->importFromCF1(aosGridMappingKeyValues.List(),
738 1 : nullptr) != OGRERR_NONE)
739 : {
740 0 : l_poSRS.reset();
741 : }
742 2 : return l_poSRS;
743 : };
744 :
745 62 : const auto poAttributes = poArray->GetAttributes();
746 57 : for (const auto &poAttr : poAttributes)
747 : {
748 27 : if (poAttr->GetDataType().GetClass() == GEDTC_STRING &&
749 19 : STARTS_WITH_CI(poAttr->GetName().c_str(),
750 46 : "__tiledb_attr.") &&
751 9 : ENDS_WITH_CI(poAttr->GetName().c_str(), ".grid_mapping"))
752 : {
753 1 : const char *pszGridMapping = poAttr->ReadAsString();
754 1 : if (pszGridMapping)
755 : {
756 3 : for (const auto &osOtherArray :
757 8 : poParent->GetMDArrayNames())
758 : {
759 4 : if (osOtherArray != osName)
760 : {
761 : auto poOtherArray =
762 3 : poParent->OpenMDArray(osOtherArray);
763 3 : if (poOtherArray)
764 : {
765 : const std::string osGMPrefix =
766 6 : std::string("__tiledb_attr.")
767 3 : .append(pszGridMapping)
768 3 : .append(".");
769 : auto poGridMappingNameAttr =
770 3 : poOtherArray->GetAttribute(
771 6 : std::string(osGMPrefix)
772 3 : .append("grid_mapping_name")
773 9 : .c_str());
774 3 : if (poGridMappingNameAttr)
775 : {
776 2 : poSRS = GetSRSFromGridMappingArray(
777 1 : poOtherArray, osGMPrefix);
778 1 : break;
779 : }
780 : }
781 : }
782 : }
783 : }
784 1 : break;
785 : }
786 : }
787 : }
788 :
789 : // Set SRS DataAxisToSRSAxisMapping
790 45 : if (poSRS)
791 : {
792 7 : int iDimX = 0;
793 7 : int iDimY = 0;
794 7 : int iCount = 1;
795 22 : for (const auto &poDim : aoDims)
796 : {
797 15 : if (poDim->GetType() == GDAL_DIM_TYPE_HORIZONTAL_X)
798 5 : iDimX = iCount;
799 10 : else if (poDim->GetType() == GDAL_DIM_TYPE_HORIZONTAL_Y)
800 5 : iDimY = iCount;
801 15 : iCount++;
802 : }
803 7 : if ((iDimX == 0 || iDimY == 0) && aoDims.size() >= 2)
804 : {
805 2 : iDimX = static_cast<int>(aoDims.size());
806 2 : iDimY = iDimX - 1;
807 : }
808 7 : if (iDimX > 0 && iDimY > 0)
809 : {
810 7 : if (poSRS->GetDataAxisToSRSAxisMapping() ==
811 14 : std::vector<int>{2, 1})
812 5 : poSRS->SetDataAxisToSRSAxisMapping({iDimY, iDimX});
813 2 : else if (poSRS->GetDataAxisToSRSAxisMapping() ==
814 4 : std::vector<int>{1, 2})
815 2 : poSRS->SetDataAxisToSRSAxisMapping({iDimX, iDimY});
816 : }
817 : }
818 :
819 45 : poArray->m_poSRS = std::move(poSRS);
820 :
821 90 : const auto filters = attr.filter_list();
822 90 : std::string osFilters;
823 46 : for (uint32_t j = 0; j < filters.nfilters(); ++j)
824 : {
825 1 : const auto filter = filters.filter(j);
826 1 : if (j > 0)
827 0 : osFilters += ',';
828 1 : osFilters += tiledb::Filter::to_str(filter.filter_type());
829 : }
830 45 : if (!osFilters.empty())
831 : {
832 1 : poArray->m_aosStructuralInfo.SetNameValue("FILTER_LIST",
833 1 : osFilters.c_str());
834 : }
835 :
836 45 : return poArray;
837 : }
838 0 : catch (const std::exception &e)
839 : {
840 0 : CPLError(CE_Failure, CPLE_AppDefined, "OpenFromDisk() failed with: %s",
841 0 : e.what());
842 0 : return nullptr;
843 : }
844 : }
845 :
846 : /************************************************************************/
847 : /* TileDBArray::EnsureOpenAs() */
848 : /************************************************************************/
849 :
850 225 : bool TileDBArray::EnsureOpenAs(tiledb_query_type_t mode) const
851 : {
852 225 : if (!m_bFinalized && !Finalize())
853 0 : return false;
854 225 : if (!m_poTileDBArray)
855 0 : return false;
856 225 : if (m_poTileDBArray->query_type() == mode && m_poTileDBArray->is_open())
857 182 : return true;
858 : try
859 : {
860 43 : m_poTileDBArray->close();
861 43 : m_poTileDBArray->open(mode);
862 : }
863 0 : catch (const std::exception &e)
864 : {
865 0 : CPLError(CE_Failure, CPLE_AppDefined, "%s", e.what());
866 0 : m_poTileDBArray.reset();
867 0 : return false;
868 : }
869 43 : return true;
870 : }
871 :
872 : /************************************************************************/
873 : /* TileDBArray::IRead() */
874 : /************************************************************************/
875 :
876 55 : bool TileDBArray::IRead(const GUInt64 *arrayStartIdx, const size_t *count,
877 : const GInt64 *arrayStep, const GPtrDiff_t *bufferStride,
878 : const GDALExtendedDataType &bufferDataType,
879 : void *pDstBuffer) const
880 : {
881 55 : if (!EnsureOpenAs(TILEDB_READ))
882 0 : return false;
883 :
884 55 : if (!IsStepOneContiguousRowMajorOrderedSameDataType(
885 : count, arrayStep, bufferStride, bufferDataType))
886 : {
887 7 : return ReadUsingContiguousIRead(arrayStartIdx, count, arrayStep,
888 : bufferStride, bufferDataType,
889 7 : pDstBuffer);
890 : }
891 : else
892 : {
893 48 : std::vector<uint64_t> anSubArray;
894 48 : const auto nDims = m_aoDims.size();
895 48 : anSubArray.reserve(2 * nDims);
896 : size_t nBufferSize =
897 48 : GDALDataTypeIsComplex(m_oType.GetNumericDataType()) ? 2 : 1;
898 113 : for (size_t i = 0; i < nDims; ++i)
899 : {
900 65 : anSubArray.push_back(m_anStartDimOffset[i] + arrayStartIdx[i]);
901 65 : anSubArray.push_back(m_anStartDimOffset[i] + arrayStartIdx[i] +
902 65 : count[i] - 1);
903 65 : nBufferSize *= count[i];
904 : }
905 : try
906 : {
907 48 : tiledb::Query query(m_poSharedResource->GetCtx(),
908 144 : *(m_poTileDBArray.get()));
909 48 : tiledb::Subarray subarray(m_poSharedResource->GetCtx(),
910 96 : *(m_poTileDBArray.get()));
911 48 : subarray.set_subarray(anSubArray);
912 48 : query.set_subarray(subarray);
913 48 : query.set_data_buffer(m_osAttrName, pDstBuffer, nBufferSize);
914 :
915 48 : if (m_bStats)
916 0 : tiledb::Stats::enable();
917 :
918 48 : const auto ret = query.submit();
919 :
920 48 : if (m_bStats)
921 : {
922 0 : tiledb::Stats::dump(stdout);
923 0 : tiledb::Stats::disable();
924 : }
925 :
926 48 : return ret == tiledb::Query::Status::COMPLETE;
927 : }
928 0 : catch (const std::exception &e)
929 : {
930 0 : CPLError(CE_Failure, CPLE_AppDefined, "Read() failed with %s",
931 0 : e.what());
932 : }
933 : }
934 :
935 0 : return false;
936 : }
937 :
938 : /************************************************************************/
939 : /* TileDBArray::IWrite() */
940 : /************************************************************************/
941 :
942 23 : bool TileDBArray::IWrite(const GUInt64 *arrayStartIdx, const size_t *count,
943 : const GInt64 *arrayStep,
944 : const GPtrDiff_t *bufferStride,
945 : const GDALExtendedDataType &bufferDataType,
946 : const void *pSrcBuffer)
947 : {
948 23 : if (!IsWritable())
949 : {
950 0 : CPLError(CE_Failure, CPLE_NotSupported,
951 : "Dataset not open in update mode");
952 0 : return false;
953 : }
954 :
955 23 : if (!EnsureOpenAs(TILEDB_WRITE))
956 0 : return false;
957 :
958 23 : if (!IsStepOneContiguousRowMajorOrderedSameDataType(
959 : count, arrayStep, bufferStride, bufferDataType))
960 : {
961 1 : CPLError(CE_Failure, CPLE_NotSupported,
962 : "Write parameters not supported");
963 1 : return false;
964 : }
965 : else
966 : {
967 22 : std::vector<uint64_t> anSubArray;
968 22 : const auto nDims = m_aoDims.size();
969 22 : anSubArray.reserve(2 * nDims);
970 : size_t nBufferSize =
971 22 : GDALDataTypeIsComplex(m_oType.GetNumericDataType()) ? 2 : 1;
972 50 : for (size_t i = 0; i < nDims; ++i)
973 : {
974 28 : anSubArray.push_back(m_anStartDimOffset[i] + arrayStartIdx[i]);
975 28 : anSubArray.push_back(m_anStartDimOffset[i] + arrayStartIdx[i] +
976 28 : count[i] - 1);
977 28 : nBufferSize *= count[i];
978 : }
979 : try
980 : {
981 22 : tiledb::Query query(m_poSharedResource->GetCtx(),
982 66 : *(m_poTileDBArray.get()));
983 22 : tiledb::Subarray subarray(m_poSharedResource->GetCtx(),
984 44 : *(m_poTileDBArray.get()));
985 22 : subarray.set_subarray(anSubArray);
986 22 : query.set_subarray(subarray);
987 22 : query.set_data_buffer(m_osAttrName, const_cast<void *>(pSrcBuffer),
988 22 : nBufferSize);
989 :
990 22 : return query.submit() == tiledb::Query::Status::COMPLETE;
991 : }
992 0 : catch (const std::exception &e)
993 : {
994 0 : CPLError(CE_Failure, CPLE_AppDefined, "Write() failed with %s",
995 0 : e.what());
996 : }
997 : }
998 :
999 0 : return false;
1000 : }
1001 :
1002 : /************************************************************************/
1003 : /* TileDBArray::GetRawNoDataValue() */
1004 : /************************************************************************/
1005 :
1006 9 : const void *TileDBArray::GetRawNoDataValue() const
1007 : {
1008 9 : if (!m_bFinalized)
1009 0 : return nullptr;
1010 :
1011 9 : if (m_abyNoData.empty())
1012 : {
1013 9 : const void *value = nullptr;
1014 9 : uint64_t size = 0;
1015 : // Caution: 2 below statements must not be combined in a single one,
1016 : // as the lifetime of value is linked to the return value of
1017 : // attribute()
1018 18 : auto attr = m_poSchema->attribute(m_osAttrName);
1019 9 : attr.get_fill_value(&value, &size);
1020 9 : if (size == m_oType.GetSize())
1021 : {
1022 9 : m_abyNoData.resize(size);
1023 9 : memcpy(m_abyNoData.data(), value, size);
1024 : }
1025 : }
1026 :
1027 9 : return m_abyNoData.empty() ? nullptr : m_abyNoData.data();
1028 : }
1029 :
1030 : /************************************************************************/
1031 : /* TileDBArray::SetRawNoDataValue() */
1032 : /************************************************************************/
1033 :
1034 3 : bool TileDBArray::SetRawNoDataValue(const void *pRawNoData)
1035 : {
1036 3 : if (m_bFinalized)
1037 : {
1038 1 : CPLError(CE_Failure, CPLE_NotSupported,
1039 : "SetRawNoDataValue() not supported after array has been "
1040 : "finalized.");
1041 1 : return false;
1042 : }
1043 :
1044 2 : if (pRawNoData)
1045 : {
1046 2 : CPLAssert(m_poAttr);
1047 2 : m_poAttr->set_fill_value(pRawNoData, m_oType.GetSize());
1048 2 : m_abyNoData.resize(m_oType.GetSize());
1049 2 : memcpy(m_abyNoData.data(), pRawNoData, m_oType.GetSize());
1050 : }
1051 :
1052 2 : Finalize();
1053 :
1054 2 : return true;
1055 : }
1056 :
1057 : /************************************************************************/
1058 : /* TileDBArray::CreateAttribute() */
1059 : /************************************************************************/
1060 :
1061 12 : std::shared_ptr<GDALAttribute> TileDBArray::CreateAttribute(
1062 : const std::string &osName, const std::vector<GUInt64> &anDimensions,
1063 : const GDALExtendedDataType &oDataType, CSLConstList papszOptions)
1064 : {
1065 12 : return CreateAttributeImpl(osName, anDimensions, oDataType, papszOptions);
1066 : }
1067 :
1068 : /************************************************************************/
1069 : /* TileDBArray::GetAttribute() */
1070 : /************************************************************************/
1071 :
1072 : std::shared_ptr<GDALAttribute>
1073 39 : TileDBArray::GetAttribute(const std::string &osName) const
1074 : {
1075 39 : return GetAttributeImpl(osName);
1076 : }
1077 :
1078 : /************************************************************************/
1079 : /* TileDBArray::GetAttributes() */
1080 : /************************************************************************/
1081 :
1082 : std::vector<std::shared_ptr<GDALAttribute>>
1083 40 : TileDBArray::GetAttributes(CSLConstList papszOptions) const
1084 : {
1085 40 : return GetAttributesImpl(papszOptions);
1086 : }
1087 :
1088 : /************************************************************************/
1089 : /* TileDBArray::DeleteAttribute() */
1090 : /************************************************************************/
1091 :
1092 0 : bool TileDBArray::DeleteAttribute(const std::string &osName,
1093 : CSLConstList papszOptions)
1094 : {
1095 0 : return DeleteAttributeImpl(osName, papszOptions);
1096 : }
1097 :
1098 : /************************************************************************/
1099 : /* TileDBArray::SetSpatialRef() */
1100 : /************************************************************************/
1101 :
1102 2 : bool TileDBArray::SetSpatialRef(const OGRSpatialReference *poSRS)
1103 : {
1104 2 : if (!IsWritable())
1105 : {
1106 0 : CPLError(CE_Failure, CPLE_NotSupported,
1107 : "Dataset not open in update mode");
1108 0 : return false;
1109 : }
1110 :
1111 2 : if (!EnsureOpenAs(TILEDB_WRITE))
1112 0 : return false;
1113 :
1114 : try
1115 : {
1116 2 : if (m_poSRS && !poSRS)
1117 0 : m_poTileDBArray->delete_metadata(CRS_ATTRIBUTE_NAME);
1118 :
1119 2 : m_poSRS.reset();
1120 2 : if (poSRS)
1121 : {
1122 2 : m_poSRS.reset(poSRS->Clone());
1123 :
1124 2 : char *pszPROJJSON = nullptr;
1125 2 : if (m_poSRS->exportToPROJJSON(&pszPROJJSON, nullptr) ==
1126 4 : OGRERR_NONE &&
1127 2 : pszPROJJSON != nullptr)
1128 : {
1129 4 : m_poTileDBArray->put_metadata(
1130 : CRS_ATTRIBUTE_NAME, TILEDB_STRING_UTF8,
1131 2 : static_cast<uint32_t>(strlen(pszPROJJSON)), pszPROJJSON);
1132 2 : CPLFree(pszPROJJSON);
1133 : }
1134 : else
1135 : {
1136 0 : CPLFree(pszPROJJSON);
1137 0 : return false;
1138 : }
1139 : }
1140 2 : return true;
1141 : }
1142 0 : catch (const std::exception &e)
1143 : {
1144 0 : CPLError(CE_Failure, CPLE_AppDefined, "SetSpatialRef() failed with: %s",
1145 0 : e.what());
1146 0 : return false;
1147 : }
1148 : }
1149 :
1150 : /************************************************************************/
1151 : /* TileDBArray::SetUnit() */
1152 : /************************************************************************/
1153 :
1154 5 : bool TileDBArray::SetUnit(const std::string &osUnit)
1155 : {
1156 5 : if (!IsWritable())
1157 : {
1158 0 : CPLError(CE_Failure, CPLE_NotSupported,
1159 : "Dataset not open in update mode");
1160 0 : return false;
1161 : }
1162 :
1163 5 : if (!EnsureOpenAs(TILEDB_WRITE))
1164 0 : return false;
1165 :
1166 : try
1167 : {
1168 5 : if (!m_osUnit.empty() && osUnit.empty())
1169 0 : m_poTileDBArray->delete_metadata(UNIT_ATTRIBUTE_NAME);
1170 :
1171 5 : m_osUnit = osUnit;
1172 5 : if (!osUnit.empty())
1173 : {
1174 10 : m_poTileDBArray->put_metadata(
1175 : UNIT_ATTRIBUTE_NAME, TILEDB_STRING_UTF8,
1176 5 : static_cast<uint32_t>(osUnit.size()), osUnit.data());
1177 : }
1178 5 : return true;
1179 : }
1180 0 : catch (const std::exception &e)
1181 : {
1182 0 : CPLError(CE_Failure, CPLE_AppDefined, "SetUnit() failed with: %s",
1183 0 : e.what());
1184 0 : return false;
1185 : }
1186 : }
1187 :
1188 : /************************************************************************/
1189 : /* FillBlockSize() */
1190 : /************************************************************************/
1191 :
1192 : static bool
1193 30 : FillBlockSize(const std::vector<std::shared_ptr<GDALDimension>> &aoDimensions,
1194 : const GDALExtendedDataType &oDataType,
1195 : std::vector<GUInt64> &anBlockSize, CSLConstList papszOptions)
1196 : {
1197 30 : const auto nDims = aoDimensions.size();
1198 30 : anBlockSize.resize(nDims);
1199 71 : for (size_t i = 0; i < nDims; ++i)
1200 41 : anBlockSize[i] = 1;
1201 30 : if (nDims >= 2)
1202 : {
1203 18 : anBlockSize[nDims - 2] =
1204 18 : std::min(std::max<GUInt64>(1, aoDimensions[nDims - 2]->GetSize()),
1205 18 : static_cast<GUInt64>(256));
1206 18 : anBlockSize[nDims - 1] =
1207 18 : std::min(std::max<GUInt64>(1, aoDimensions[nDims - 1]->GetSize()),
1208 27 : static_cast<GUInt64>(256));
1209 : }
1210 21 : else if (nDims == 1)
1211 : {
1212 21 : anBlockSize[0] = std::max<GUInt64>(1, aoDimensions[0]->GetSize());
1213 : }
1214 :
1215 30 : const char *pszBlockSize = CSLFetchNameValue(papszOptions, "BLOCKSIZE");
1216 30 : if (pszBlockSize)
1217 : {
1218 : const auto aszTokens(
1219 1 : CPLStringList(CSLTokenizeString2(pszBlockSize, ",", 0)));
1220 1 : if (static_cast<size_t>(aszTokens.size()) != nDims)
1221 : {
1222 0 : CPLError(CE_Failure, CPLE_AppDefined,
1223 : "Invalid number of values in BLOCKSIZE");
1224 0 : return false;
1225 : }
1226 1 : size_t nBlockSize = oDataType.GetSize();
1227 3 : for (size_t i = 0; i < nDims; ++i)
1228 : {
1229 2 : anBlockSize[i] = static_cast<GUInt64>(CPLAtoGIntBig(aszTokens[i]));
1230 2 : if (anBlockSize[i] == 0)
1231 : {
1232 0 : CPLError(CE_Failure, CPLE_AppDefined,
1233 : "Values in BLOCKSIZE should be > 0");
1234 0 : return false;
1235 : }
1236 2 : if (anBlockSize[i] >
1237 2 : std::numeric_limits<size_t>::max() / nBlockSize)
1238 : {
1239 0 : CPLError(CE_Failure, CPLE_AppDefined,
1240 : "Too large values in BLOCKSIZE");
1241 0 : return false;
1242 : }
1243 2 : nBlockSize *= static_cast<size_t>(anBlockSize[i]);
1244 : }
1245 : }
1246 30 : return true;
1247 : }
1248 :
1249 : /************************************************************************/
1250 : /* TileDBArray::GDALDataTypeToTileDB() */
1251 : /************************************************************************/
1252 :
1253 52 : /*static*/ bool TileDBArray::GDALDataTypeToTileDB(GDALDataType dt,
1254 : tiledb_datatype_t &tiledb_dt)
1255 : {
1256 52 : switch (dt)
1257 : {
1258 4 : case GDT_Byte:
1259 4 : tiledb_dt = TILEDB_UINT8;
1260 4 : break;
1261 1 : case GDT_Int8:
1262 1 : tiledb_dt = TILEDB_INT8;
1263 1 : break;
1264 1 : case GDT_UInt16:
1265 1 : tiledb_dt = TILEDB_UINT16;
1266 1 : break;
1267 2 : case GDT_CInt16:
1268 : case GDT_Int16:
1269 2 : tiledb_dt = TILEDB_INT16;
1270 2 : break;
1271 1 : case GDT_UInt32:
1272 1 : tiledb_dt = TILEDB_UINT32;
1273 1 : break;
1274 7 : case GDT_CInt32:
1275 : case GDT_Int32:
1276 7 : tiledb_dt = TILEDB_INT32;
1277 7 : break;
1278 1 : case GDT_UInt64:
1279 1 : tiledb_dt = TILEDB_UINT64;
1280 1 : break;
1281 1 : case GDT_Int64:
1282 1 : tiledb_dt = TILEDB_INT64;
1283 1 : break;
1284 8 : case GDT_CFloat32:
1285 : case GDT_Float32:
1286 8 : tiledb_dt = TILEDB_FLOAT32;
1287 8 : break;
1288 26 : case GDT_CFloat64:
1289 : case GDT_Float64:
1290 26 : tiledb_dt = TILEDB_FLOAT64;
1291 26 : break;
1292 :
1293 0 : case GDT_Unknown:
1294 : case GDT_TypeCount:
1295 : {
1296 0 : CPLError(CE_Failure, CPLE_NotSupported, "Unsupported data type: %s",
1297 : GDALGetDataTypeName(dt));
1298 0 : return false;
1299 : }
1300 : }
1301 52 : return true;
1302 : }
1303 :
1304 : /************************************************************************/
1305 : /* IsIncreasingOrDecreasing1DVar() */
1306 : /************************************************************************/
1307 :
1308 : static void
1309 8 : IsIncreasingOrDecreasing1DVar(const std::shared_ptr<GDALMDArray> &poVar,
1310 : bool &bIncreasing, bool &bDecreasing)
1311 : {
1312 8 : bIncreasing = false;
1313 8 : bDecreasing = false;
1314 16 : std::vector<double> adfVals;
1315 : try
1316 : {
1317 8 : adfVals.resize(
1318 8 : static_cast<size_t>(poVar->GetDimensions()[0]->GetSize()));
1319 : }
1320 0 : catch (const std::exception &)
1321 : {
1322 : }
1323 8 : if (adfVals.size() > 1)
1324 : {
1325 8 : GUInt64 anStart[1] = {0};
1326 8 : size_t anCount[1] = {adfVals.size()};
1327 16 : if (poVar->Read(anStart, anCount, nullptr, nullptr,
1328 16 : GDALExtendedDataType::Create(GDT_Float64),
1329 8 : adfVals.data()))
1330 : {
1331 8 : if (adfVals[1] > adfVals[0])
1332 5 : bIncreasing = true;
1333 3 : else if (adfVals[1] < adfVals[0])
1334 3 : bDecreasing = true;
1335 8 : if (bIncreasing || bDecreasing)
1336 : {
1337 34 : for (size_t i = 2; i < adfVals.size(); ++i)
1338 : {
1339 26 : if (bIncreasing)
1340 : {
1341 16 : if (!(adfVals[i] > adfVals[i - 1]))
1342 : {
1343 0 : bIncreasing = false;
1344 0 : break;
1345 : }
1346 : }
1347 : else
1348 : {
1349 10 : if (!(adfVals[i] < adfVals[i - 1]))
1350 : {
1351 0 : bDecreasing = false;
1352 0 : break;
1353 : }
1354 : }
1355 : }
1356 : }
1357 : }
1358 : }
1359 8 : }
1360 :
1361 : /************************************************************************/
1362 : /* TileDBArray::CreateOnDisk() */
1363 : /************************************************************************/
1364 :
1365 : /* static */
1366 34 : std::shared_ptr<TileDBArray> TileDBArray::CreateOnDisk(
1367 : const std::shared_ptr<TileDBSharedResource> &poSharedResource,
1368 : const std::shared_ptr<TileDBGroup> &poParent, const std::string &osName,
1369 : const std::vector<std::shared_ptr<GDALDimension>> &aoDimensions,
1370 : const GDALExtendedDataType &oDataType, CSLConstList papszOptions)
1371 : {
1372 34 : if (aoDimensions.empty())
1373 : {
1374 1 : CPLError(CE_Failure, CPLE_NotSupported,
1375 : "Zero-dimensions arrays are not supported by TileDB");
1376 1 : return nullptr;
1377 : }
1378 :
1379 33 : tiledb_datatype_t tiledb_dt = TILEDB_ANY;
1380 33 : if (oDataType.GetClass() != GEDTC_NUMERIC)
1381 : {
1382 1 : CPLError(CE_Failure, CPLE_NotSupported,
1383 : "Only numeric data types are supported");
1384 1 : return nullptr;
1385 : }
1386 32 : if (!GDALDataTypeToTileDB(oDataType.GetNumericDataType(), tiledb_dt))
1387 0 : return nullptr;
1388 :
1389 : try
1390 : {
1391 : const auto osSanitizedName =
1392 64 : TileDBSharedResource::SanitizeNameForPath(osName);
1393 63 : if (osSanitizedName.empty() || STARTS_WITH(osName.c_str(), "./") ||
1394 31 : STARTS_WITH(osName.c_str(), "../") ||
1395 94 : STARTS_WITH(osName.c_str(), ".\\") ||
1396 31 : STARTS_WITH(osName.c_str(), "..\\"))
1397 : {
1398 1 : CPLError(CE_Failure, CPLE_AppDefined, "Invalid array name");
1399 1 : return nullptr;
1400 : }
1401 62 : std::string osArrayPath = poParent->GetPath() + "/" + osSanitizedName;
1402 31 : const char *pszURI = CSLFetchNameValue(papszOptions, "URI");
1403 31 : if (pszURI)
1404 0 : osArrayPath = pszURI;
1405 :
1406 31 : auto &ctx = poSharedResource->GetCtx();
1407 62 : tiledb::VFS vfs(ctx);
1408 31 : if (vfs.is_dir(osArrayPath))
1409 : {
1410 1 : CPLError(CE_Failure, CPLE_AppDefined, "Path %s already exists",
1411 : osArrayPath.c_str());
1412 1 : return nullptr;
1413 : }
1414 :
1415 60 : std::vector<GUInt64> anBlockSize;
1416 : #if defined(__GNUC__)
1417 : #pragma GCC diagnostic push
1418 : #pragma GCC diagnostic ignored "-Wnull-dereference"
1419 : #endif
1420 30 : if (!FillBlockSize(aoDimensions, oDataType, anBlockSize, papszOptions))
1421 0 : return nullptr;
1422 : #if defined(__GNUC__)
1423 : #pragma GCC diagnostic pop
1424 : #endif
1425 : auto poSchema =
1426 60 : std::make_unique<tiledb::ArraySchema>(ctx, TILEDB_DENSE);
1427 30 : poSchema->set_tile_order(TILEDB_ROW_MAJOR);
1428 30 : poSchema->set_cell_order(TILEDB_ROW_MAJOR);
1429 :
1430 60 : tiledb::FilterList filterList(ctx);
1431 : const char *pszCompression =
1432 30 : CSLFetchNameValue(papszOptions, "COMPRESSION");
1433 : const char *pszCompressionLevel =
1434 30 : CSLFetchNameValue(papszOptions, "COMPRESSION_LEVEL");
1435 :
1436 30 : if (pszCompression != nullptr)
1437 : {
1438 1 : int nLevel = (pszCompressionLevel) ? atoi(pszCompressionLevel) : -1;
1439 1 : if (TileDBDataset::AddFilter(ctx, filterList, pszCompression,
1440 1 : nLevel) != CE_None)
1441 : {
1442 0 : return nullptr;
1443 : }
1444 : }
1445 30 : poSchema->set_coords_filter_list(filterList);
1446 :
1447 60 : tiledb::Domain domain(ctx);
1448 71 : for (size_t i = 0; i < aoDimensions.size(); ++i)
1449 : {
1450 41 : const auto &poDim = aoDimensions[i];
1451 41 : if (poDim->GetSize() == 0)
1452 : {
1453 0 : CPLError(CE_Failure, CPLE_AppDefined, "Invalid dim size: 0");
1454 0 : return nullptr;
1455 : }
1456 82 : std::string osDimName(poDim->GetName());
1457 41 : if (poDim->GetName() == osName)
1458 5 : osDimName += "_dim";
1459 : auto dim = tiledb::Dimension::create<uint64_t>(
1460 82 : ctx, osDimName, {0, poDim->GetSize() - 1}, anBlockSize[i]);
1461 41 : domain.add_dimension(std::move(dim));
1462 : }
1463 :
1464 30 : poSchema->set_domain(domain);
1465 :
1466 60 : std::vector<std::shared_ptr<GDALMDArray>> apoIndexingVariables;
1467 71 : for (size_t i = 0; i < aoDimensions.size(); ++i)
1468 : {
1469 41 : const auto &poDim = aoDimensions[i];
1470 82 : auto poIndexingVar = poDim->GetIndexingVariable();
1471 41 : bool bDimLabelCreated = false;
1472 41 : tiledb_datatype_t dim_label_tiledb_dt = TILEDB_ANY;
1473 57 : if (poIndexingVar && poIndexingVar->GetDimensionCount() == 1 &&
1474 16 : poIndexingVar->GetDataType().GetClass() == GEDTC_NUMERIC &&
1475 8 : poIndexingVar->GetDimensions()[0]->GetName() ==
1476 16 : poDim->GetName() &&
1477 8 : poIndexingVar->GetDimensions()[0]->GetSize() <
1478 8 : 10 * 1024 * 1024 &&
1479 8 : !GDALDataTypeIsComplex(
1480 57 : poIndexingVar->GetDataType().GetNumericDataType()) &&
1481 8 : GDALDataTypeToTileDB(
1482 8 : poIndexingVar->GetDataType().GetNumericDataType(),
1483 : dim_label_tiledb_dt))
1484 : {
1485 8 : bool bIncreasing = false;
1486 8 : bool bDecreasing = false;
1487 8 : IsIncreasingOrDecreasing1DVar(poIndexingVar, bIncreasing,
1488 : bDecreasing);
1489 8 : if (bIncreasing || bDecreasing)
1490 : {
1491 8 : bDimLabelCreated = true;
1492 8 : apoIndexingVariables.push_back(poIndexingVar);
1493 16 : tiledb::ArraySchemaExperimental::add_dimension_label(
1494 8 : ctx, *(poSchema.get()), static_cast<uint32_t>(i),
1495 16 : BuildDimensionLabelName(poDim),
1496 : bIncreasing ? TILEDB_INCREASING_DATA
1497 : : TILEDB_DECREASING_DATA,
1498 : dim_label_tiledb_dt,
1499 16 : std::optional<tiledb::FilterList>(filterList));
1500 : }
1501 : }
1502 41 : if (poIndexingVar && !bDimLabelCreated)
1503 : {
1504 0 : CPLDebug("TILEDB",
1505 : "Dimension %s has indexing variable %s, "
1506 : "but not compatible of a dimension label",
1507 0 : poDim->GetName().c_str(),
1508 0 : poIndexingVar->GetName().c_str());
1509 : }
1510 : }
1511 :
1512 : auto attr = std::make_unique<tiledb::Attribute>(
1513 60 : tiledb::Attribute::create(ctx, osName, tiledb_dt));
1514 30 : if (GDALDataTypeIsComplex(oDataType.GetNumericDataType()))
1515 4 : attr->set_cell_val_num(2);
1516 30 : attr->set_filter_list(filterList);
1517 :
1518 : // Implement a deferred TileDB array creation given that we might
1519 : // need to set the fill value of the attribute from the nodata value
1520 30 : auto poArray = Create(poSharedResource, poParent->GetFullName(), osName,
1521 60 : aoDimensions, oDataType, osArrayPath);
1522 30 : poArray->m_bFinalized = false;
1523 30 : poArray->m_poParent = poParent;
1524 30 : poArray->m_osParentPath = poParent->GetPath();
1525 30 : poArray->m_poSchema = std::move(poSchema);
1526 30 : poArray->m_osAttrName = attr->name();
1527 30 : poArray->m_poAttr = std::move(attr);
1528 30 : poArray->m_anBlockSize = std::move(anBlockSize);
1529 30 : poArray->m_anStartDimOffset.resize(aoDimensions.size());
1530 : // To keep a reference on the indexing variables, so they are still
1531 : // alive at Finalize() time
1532 30 : poArray->m_apoIndexingVariables = std::move(apoIndexingVariables);
1533 30 : if (CPLTestBool(CSLFetchNameValueDef(papszOptions, "STATS", "FALSE")))
1534 0 : poArray->m_bStats = true;
1535 :
1536 30 : uint64_t nTimestamp = poSharedResource->GetTimestamp();
1537 : const char *pszTimestamp =
1538 30 : CSLFetchNameValue(papszOptions, "TILEDB_TIMESTAMP");
1539 30 : if (pszTimestamp)
1540 0 : nTimestamp = std::strtoull(pszTimestamp, nullptr, 10);
1541 30 : poArray->m_nTimestamp = nTimestamp;
1542 :
1543 30 : return poArray;
1544 : }
1545 0 : catch (const std::exception &e)
1546 : {
1547 0 : CPLError(CE_Failure, CPLE_AppDefined, "CreateMDArray() failed with: %s",
1548 0 : e.what());
1549 0 : return nullptr;
1550 : }
1551 : }
1552 :
1553 : /************************************************************************/
1554 : /* TileDBArray::GetStructuralInfo() */
1555 : /************************************************************************/
1556 :
1557 7 : CSLConstList TileDBArray::GetStructuralInfo() const
1558 : {
1559 7 : return m_aosStructuralInfo.List();
1560 : }
|
