Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: OpenGIS Simple Features Reference Implementation
4 : * Purpose: Implements writing of FileGDB indices
5 : * Author: Even Rouault, <even dot rouault at spatialys.com>
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2022, Even Rouault <even dot rouault at spatialys.com>
9 : *
10 : * SPDX-License-Identifier: MIT
11 : ****************************************************************************/
12 :
13 : #include "cpl_port.h"
14 :
15 : #include "filegdbtable.h"
16 : #include "filegdbtable_priv.h"
17 :
18 : #include <cctype>
19 : #include <cstdint>
20 : #include <algorithm>
21 : #include <limits>
22 :
23 : #include "cpl_string.h"
24 :
25 : namespace OpenFileGDB
26 : {
27 :
28 : /************************************************************************/
29 : /* RemoveIndices() */
30 : /************************************************************************/
31 :
32 3097 : void FileGDBTable::RemoveIndices()
33 : {
34 3097 : if (!m_bUpdate)
35 0 : return;
36 :
37 3097 : CPLString osUCGeomFieldName;
38 3097 : if (m_iGeomField >= 0)
39 : {
40 881 : osUCGeomFieldName = m_apoFields[m_iGeomField]->GetName();
41 881 : osUCGeomFieldName.toupper();
42 : }
43 :
44 3097 : GetIndexCount();
45 3681 : for (const auto &poIndex : m_apoIndexes)
46 : {
47 1168 : if (m_iObjectIdField >= 0 &&
48 584 : m_apoFields[m_iObjectIdField]->m_poIndex == poIndex.get())
49 : {
50 321 : continue;
51 : }
52 :
53 526 : CPLString osUCIndexFieldName(poIndex->GetExpression());
54 263 : osUCIndexFieldName.toupper();
55 263 : if (osUCIndexFieldName == osUCGeomFieldName)
56 : {
57 243 : VSIUnlink(
58 486 : CPLResetExtensionSafe(m_osFilename.c_str(), "spx").c_str());
59 : }
60 : else
61 : {
62 20 : VSIUnlink(CPLResetExtensionSafe(
63 : m_osFilename.c_str(),
64 40 : (poIndex->GetIndexName() + ".atx").c_str())
65 : .c_str());
66 : }
67 : }
68 :
69 3097 : m_nHasSpatialIndex = false;
70 : }
71 :
72 : /************************************************************************/
73 : /* RefreshIndices() */
74 : /************************************************************************/
75 :
76 3052 : void FileGDBTable::RefreshIndices()
77 : {
78 3052 : if (!m_bUpdate)
79 0 : return;
80 :
81 3052 : RemoveIndices();
82 :
83 3540 : for (const auto &poIndex : m_apoIndexes)
84 : {
85 976 : if (m_iObjectIdField >= 0 &&
86 488 : m_apoFields[m_iObjectIdField]->m_poIndex == poIndex.get())
87 : {
88 276 : continue;
89 : }
90 :
91 210 : if (m_iGeomField >= 0 &&
92 422 : m_apoFields[m_iGeomField]->m_poIndex == poIndex.get() &&
93 199 : m_eTableGeomType != FGTGT_MULTIPATCH)
94 : {
95 193 : CreateSpatialIndex();
96 : }
97 : else
98 : {
99 38 : const std::string osFieldName = poIndex->GetFieldName();
100 19 : const int iField = GetFieldIdx(osFieldName);
101 19 : if (iField >= 0)
102 : {
103 19 : const auto eFieldType = m_apoFields[iField]->GetType();
104 19 : if (eFieldType == FGFT_INT16 || eFieldType == FGFT_INT32 ||
105 16 : eFieldType == FGFT_FLOAT32 || eFieldType == FGFT_FLOAT64 ||
106 13 : eFieldType == FGFT_STRING || eFieldType == FGFT_DATETIME)
107 : {
108 7 : CreateAttributeIndex(poIndex.get());
109 : }
110 : }
111 : }
112 : }
113 : }
114 :
115 : /************************************************************************/
116 : /* CreateIndex() */
117 : /************************************************************************/
118 :
119 581 : bool FileGDBTable::CreateIndex(const std::string &osIndexName,
120 : const std::string &osExpression)
121 : {
122 581 : if (!m_bUpdate)
123 0 : return false;
124 :
125 1162 : if (osIndexName.empty() ||
126 581 : !((osIndexName[0] >= 'a' && osIndexName[0] <= 'z') ||
127 564 : (osIndexName[0] >= 'A' && osIndexName[0] <= 'Z')))
128 : {
129 1 : CPLError(CE_Failure, CPLE_AppDefined,
130 : "Invalid index name: must start with a letter");
131 1 : return false;
132 : }
133 :
134 6814 : for (const char ch : osIndexName)
135 : {
136 6235 : if (!isalnum(static_cast<unsigned char>(ch)) && ch != '_')
137 : {
138 1 : CPLError(CE_Failure, CPLE_AppDefined,
139 : "Invalid index name: must contain only alpha numeric "
140 : "character or _");
141 1 : return false;
142 : }
143 : }
144 :
145 579 : if (osIndexName.size() > 16)
146 : {
147 1 : CPLError(CE_Failure, CPLE_AppDefined,
148 : "Invalid index name: cannot be greater than 16 characters");
149 1 : return false;
150 : }
151 :
152 578 : GetIndexCount();
153 883 : for (const auto &poIndex : m_apoIndexes)
154 : {
155 306 : if (EQUAL(poIndex->GetIndexName().c_str(), osIndexName.c_str()))
156 : {
157 1 : CPLError(CE_Failure, CPLE_AppDefined,
158 : "An index with same name already exists");
159 1 : return false;
160 : }
161 : }
162 :
163 : const std::string osFieldName =
164 1154 : FileGDBIndex::GetFieldNameFromExpression(osExpression);
165 577 : const int iField = GetFieldIdx(osFieldName);
166 577 : if (iField < 0)
167 : {
168 1 : CPLError(CE_Failure, CPLE_AppDefined, "Cannot find field %s",
169 : osFieldName.c_str());
170 1 : return false;
171 : }
172 :
173 576 : if (m_apoFields[iField]->m_poIndex != nullptr)
174 : {
175 2 : CPLError(CE_Failure, CPLE_AppDefined,
176 : "Field %s has already a registered index",
177 : osFieldName.c_str());
178 2 : return false;
179 : }
180 :
181 574 : const auto eFieldType = m_apoFields[iField]->GetType();
182 574 : if (eFieldType != FGFT_OBJECTID && eFieldType != FGFT_GEOMETRY &&
183 10 : eFieldType != FGFT_INT16 && eFieldType != FGFT_INT32 &&
184 8 : eFieldType != FGFT_FLOAT32 && eFieldType != FGFT_FLOAT64 &&
185 5 : eFieldType != FGFT_STRING && eFieldType != FGFT_DATETIME &&
186 3 : eFieldType != FGFT_INT64 && eFieldType != FGFT_DATE &&
187 1 : eFieldType != FGFT_TIME && eFieldType != FGFT_DATETIME_WITH_OFFSET)
188 : {
189 : // FGFT_GUID could potentially be added (cf a00000007.gdbindexes /
190 : // GDBItemRelationshipTypes ) Not sure about FGFT_GLOBALID, FGFT_XML or
191 : // FGFT_RASTER
192 0 : CPLError(CE_Failure, CPLE_AppDefined,
193 : "Unsupported field type for index creation");
194 0 : return false;
195 : }
196 :
197 574 : m_bDirtyGdbIndexesFile = true;
198 :
199 1148 : auto poIndex = std::make_unique<FileGDBIndex>();
200 574 : poIndex->m_osIndexName = osIndexName;
201 574 : poIndex->m_osExpression = osExpression;
202 :
203 574 : if (iField != m_iObjectIdField && iField != m_iGeomField)
204 : {
205 11 : if (!CreateAttributeIndex(poIndex.get()))
206 0 : return false;
207 : }
208 :
209 574 : m_apoFields[iField]->m_poIndex = poIndex.get();
210 :
211 574 : m_apoIndexes.push_back(std::move(poIndex));
212 :
213 574 : return true;
214 : }
215 :
216 : /************************************************************************/
217 : /* CreateGdbIndexesFile() */
218 : /************************************************************************/
219 :
220 329 : void FileGDBTable::CreateGdbIndexesFile()
221 : {
222 329 : std::vector<GByte> abyBuffer;
223 :
224 329 : WriteUInt32(abyBuffer, static_cast<uint32_t>(m_apoIndexes.size()));
225 923 : for (const auto &poIndex : m_apoIndexes)
226 : {
227 594 : const FileGDBField *poField = nullptr;
228 947 : for (size_t i = 0; i < m_apoFields.size(); i++)
229 : {
230 947 : if (CPLString(poIndex->GetFieldName()).toupper() ==
231 1894 : CPLString(m_apoFields[i]->GetName()).toupper())
232 : {
233 594 : poField = m_apoFields[i].get();
234 594 : break;
235 : }
236 : }
237 594 : if (poField == nullptr)
238 : {
239 0 : CPLError(CE_Failure, CPLE_AppDefined,
240 : "Cannot find field corresponding to index field name %s",
241 0 : poIndex->GetFieldName().c_str());
242 0 : return;
243 : }
244 :
245 594 : WriteUTF16String(abyBuffer, poIndex->GetIndexName().c_str(),
246 : NUMBER_OF_CHARS_ON_UINT32);
247 594 : WriteUInt16(abyBuffer, 0); // unknown semantics
248 594 : if (poField->GetType() == FGFT_OBJECTID)
249 : {
250 329 : WriteUInt32(abyBuffer, 16); // unknown semantics
251 329 : WriteUInt16(abyBuffer, 0xFFFF); // unknown semantics
252 : }
253 265 : else if (poField->GetType() == FGFT_GEOMETRY)
254 : {
255 241 : WriteUInt32(abyBuffer, 4); // unknown semantics
256 241 : WriteUInt16(abyBuffer, 0); // unknown semantics
257 : }
258 : else
259 : {
260 24 : WriteUInt32(abyBuffer, 2); // unknown semantics
261 24 : WriteUInt16(abyBuffer, 0); // unknown semantics
262 : }
263 594 : WriteUInt32(abyBuffer, 1); // unknown semantics
264 594 : WriteUTF16String(abyBuffer, poIndex->GetExpression().c_str(),
265 : NUMBER_OF_CHARS_ON_UINT32);
266 594 : WriteUInt16(abyBuffer, 0); // unknown semantics
267 : }
268 :
269 329 : VSILFILE *fp = VSIFOpenL(
270 658 : CPLResetExtensionSafe(m_osFilename.c_str(), "gdbindexes").c_str(),
271 : "wb");
272 329 : if (fp == nullptr)
273 0 : return;
274 329 : CPL_IGNORE_RET_VAL(VSIFWriteL(abyBuffer.data(), abyBuffer.size(), 1, fp));
275 329 : VSIFCloseL(fp);
276 : }
277 :
278 : /************************************************************************/
279 : /* ComputeOptimalSpatialIndexGridResolution() */
280 : /************************************************************************/
281 :
282 193 : void FileGDBTable::ComputeOptimalSpatialIndexGridResolution()
283 : {
284 386 : if (m_nValidRecordCount == 0 || m_iGeomField < 0 ||
285 193 : m_adfSpatialIndexGridResolution.size() != 1)
286 : {
287 0 : return;
288 : }
289 :
290 : auto poGeomField =
291 193 : cpl::down_cast<FileGDBGeomField *>(m_apoFields[m_iGeomField].get());
292 193 : if (m_eTableGeomType == FGTGT_POINT)
293 : {
294 : // For point, use the density as the grid resolution
295 75 : int nValid = 0;
296 3927 : for (int64_t iCurFeat = 0; iCurFeat < m_nTotalRecordCount; ++iCurFeat)
297 : {
298 3852 : iCurFeat = GetAndSelectNextNonEmptyRow(iCurFeat);
299 3852 : if (iCurFeat < 0)
300 0 : break;
301 3852 : const OGRField *psField = GetFieldValue(m_iGeomField);
302 3852 : if (psField != nullptr)
303 3806 : nValid++;
304 : }
305 75 : if (nValid > 0)
306 : {
307 : const double dfArea =
308 58 : (poGeomField->GetXMax() - poGeomField->GetXMin()) *
309 58 : (poGeomField->GetYMax() - poGeomField->GetYMin());
310 58 : if (dfArea != 0)
311 : {
312 33 : m_adfSpatialIndexGridResolution[0] = sqrt(dfArea / nValid);
313 : }
314 25 : else if (poGeomField->GetXMax() > poGeomField->GetXMin())
315 : {
316 0 : m_adfSpatialIndexGridResolution[0] =
317 0 : (poGeomField->GetXMax() - poGeomField->GetXMin()) / nValid;
318 : }
319 25 : else if (poGeomField->GetYMax() > poGeomField->GetYMin())
320 : {
321 0 : m_adfSpatialIndexGridResolution[0] =
322 0 : (poGeomField->GetYMax() - poGeomField->GetYMin()) / nValid;
323 : }
324 : else
325 : {
326 25 : return;
327 : }
328 33 : m_bDirtyGeomFieldSpatialIndexGridRes = true;
329 : poGeomField->m_adfSpatialIndexGridResolution =
330 33 : m_adfSpatialIndexGridResolution;
331 : }
332 : }
333 :
334 118 : else if (m_eTableGeomType == FGTGT_MULTIPOINT)
335 : {
336 : // For multipoint, use the density as the grid resolution
337 13 : int64_t nValid = 0;
338 : auto poGeomConverter = std::unique_ptr<FileGDBOGRGeometryConverter>(
339 13 : FileGDBOGRGeometryConverter::BuildConverter(poGeomField));
340 46 : for (int64_t iCurFeat = 0; iCurFeat < m_nTotalRecordCount; ++iCurFeat)
341 : {
342 33 : iCurFeat = GetAndSelectNextNonEmptyRow(iCurFeat);
343 33 : if (iCurFeat < 0)
344 0 : break;
345 33 : const OGRField *psField = GetFieldValue(m_iGeomField);
346 33 : if (psField != nullptr)
347 : {
348 : auto poGeom = std::unique_ptr<OGRGeometry>(
349 58 : poGeomConverter->GetAsGeometry(psField));
350 58 : if (poGeom != nullptr &&
351 29 : wkbFlatten(poGeom->getGeometryType()) == wkbMultiPoint)
352 : {
353 29 : nValid += poGeom->toMultiPoint()->getNumGeometries();
354 : }
355 : }
356 : }
357 13 : if (nValid > 0)
358 : {
359 : const double dfArea =
360 13 : (poGeomField->GetXMax() - poGeomField->GetXMin()) *
361 13 : (poGeomField->GetYMax() - poGeomField->GetYMin());
362 13 : if (dfArea != 0)
363 : {
364 8 : m_adfSpatialIndexGridResolution[0] = sqrt(dfArea / nValid);
365 : }
366 5 : else if (poGeomField->GetXMax() > poGeomField->GetXMin())
367 : {
368 0 : m_adfSpatialIndexGridResolution[0] =
369 0 : (poGeomField->GetXMax() - poGeomField->GetXMin()) / nValid;
370 : }
371 5 : else if (poGeomField->GetYMax() > poGeomField->GetYMin())
372 : {
373 0 : m_adfSpatialIndexGridResolution[0] =
374 0 : (poGeomField->GetYMax() - poGeomField->GetYMin()) / nValid;
375 : }
376 : else
377 : {
378 5 : return;
379 : }
380 8 : m_bDirtyGeomFieldSpatialIndexGridRes = true;
381 : poGeomField->m_adfSpatialIndexGridResolution =
382 8 : m_adfSpatialIndexGridResolution;
383 : }
384 : }
385 :
386 : else
387 : {
388 105 : CPLDebug("OpenFileGDB", "Computing optimal grid size...");
389 :
390 : // For other types of geometries, just take the maximum extent along x/y
391 : // of all geometries
392 105 : double dfMaxSize = 0;
393 105 : OGREnvelope sEnvelope;
394 391 : for (int64_t iCurFeat = 0; iCurFeat < m_nTotalRecordCount; ++iCurFeat)
395 : {
396 286 : iCurFeat = GetAndSelectNextNonEmptyRow(iCurFeat);
397 286 : if (iCurFeat < 0)
398 0 : break;
399 286 : const OGRField *psField = GetFieldValue(m_iGeomField);
400 286 : if (psField != nullptr)
401 : {
402 243 : if (GetFeatureExtent(psField, &sEnvelope))
403 : {
404 243 : dfMaxSize =
405 243 : std::max(dfMaxSize, sEnvelope.MaxX - sEnvelope.MinX);
406 243 : dfMaxSize =
407 243 : std::max(dfMaxSize, sEnvelope.MaxY - sEnvelope.MinY);
408 : }
409 : }
410 : }
411 105 : CPLDebug("OpenFileGDB", "Optimal grid size = %f", dfMaxSize);
412 :
413 105 : if (dfMaxSize > 0)
414 : {
415 86 : m_bDirtyGeomFieldSpatialIndexGridRes = true;
416 86 : m_adfSpatialIndexGridResolution[0] = dfMaxSize;
417 : poGeomField->m_adfSpatialIndexGridResolution =
418 86 : m_adfSpatialIndexGridResolution;
419 : }
420 : }
421 : }
422 :
423 : /************************************************************************/
424 : /* SortByAscendingValuesAndOID() */
425 : /************************************************************************/
426 :
427 : // recent libc++ std::sort() involve unsigned integer overflow in some
428 : // situation
429 : template <class ValueOIDPair>
430 : CPL_NOSANITIZE_UNSIGNED_INT_OVERFLOW static void
431 145 : SortByAscendingValuesAndOID(std::vector<ValueOIDPair> &asValues)
432 : {
433 145 : std::sort(asValues.begin(), asValues.end(),
434 47153 : [](const ValueOIDPair &a, const ValueOIDPair &b) {
435 65332 : return a.first < b.first ||
436 65332 : (a.first == b.first && a.second < b.second);
437 : });
438 145 : }
439 :
440 : /************************************************************************/
441 : /* WriteIndex() */
442 : /************************************************************************/
443 :
444 : template <class ValueOIDPair>
445 146 : static bool WriteIndex(
446 : VSILFILE *fp, std::vector<ValueOIDPair> &asValues,
447 : void (*writeValueFunc)(std::vector<GByte> &abyPage,
448 : const typename ValueOIDPair::first_type &value,
449 : int maxStrSize),
450 : int &nDepth, int maxStrSize = 0)
451 : {
452 146 : constexpr int IDX_PAGE_SIZE = 4096;
453 146 : constexpr int HEADER_SIZE_PAGE_REFERENCING_FEATURES = 12; // 3 * int32
454 146 : constexpr int SIZEOF_FEATURE_ID = 4; // sizeof(int)
455 146 : const int SIZEOF_INDEXED_VALUE =
456 4 : maxStrSize ? sizeof(uint16_t) * maxStrSize
457 : : sizeof(typename ValueOIDPair::first_type);
458 146 : const int NUM_MAX_FEATURES_PER_PAGE =
459 : (IDX_PAGE_SIZE - HEADER_SIZE_PAGE_REFERENCING_FEATURES) /
460 146 : (SIZEOF_FEATURE_ID + SIZEOF_INDEXED_VALUE);
461 : // static_assert(NUM_MAX_FEATURES_PER_PAGE == 340,
462 : // "NUM_MAX_FEATURES_PER_PAGE == 340");
463 146 : const int OFFSET_FIRST_VAL_IN_PAGE =
464 146 : HEADER_SIZE_PAGE_REFERENCING_FEATURES +
465 : NUM_MAX_FEATURES_PER_PAGE * SIZEOF_FEATURE_ID;
466 :
467 : // Configurable only for debugging & autotest purposes
468 292 : const int numMaxFeaturesPerPage = [NUM_MAX_FEATURES_PER_PAGE]()
469 : {
470 146 : const int nVal = atoi(
471 : CPLGetConfigOption("OPENFILEGDB_MAX_FEATURES_PER_SPX_PAGE",
472 : CPLSPrintf("%d", NUM_MAX_FEATURES_PER_PAGE)));
473 146 : if (nVal < 2)
474 0 : return 2;
475 146 : if (nVal > NUM_MAX_FEATURES_PER_PAGE)
476 0 : return NUM_MAX_FEATURES_PER_PAGE;
477 146 : return nVal;
478 146 : }();
479 :
480 292 : if (asValues.size() > static_cast<size_t>(INT_MAX) ||
481 : // Maximum number of values for depth == 4: this evaluates to ~ 13
482 : // billion values (~ features)
483 146 : asValues.size() > (((static_cast<uint64_t>(numMaxFeaturesPerPage) + 1) *
484 146 : numMaxFeaturesPerPage +
485 146 : 1) *
486 146 : numMaxFeaturesPerPage +
487 146 : 1) *
488 146 : numMaxFeaturesPerPage)
489 : {
490 1 : CPLError(CE_Failure, CPLE_NotSupported,
491 : "More values in spatial index than can be handled");
492 1 : return false;
493 : }
494 :
495 : // Sort by ascending values, and for same value by ascending OID
496 145 : SortByAscendingValuesAndOID(asValues);
497 :
498 145 : bool bRet = true;
499 290 : std::vector<GByte> abyPage;
500 145 : abyPage.reserve(IDX_PAGE_SIZE);
501 :
502 145 : const auto WriteRootPageNonLeaf =
503 240 : [=, &bRet, &asValues, &abyPage](int nNumDirectChildren,
504 : int nSubPageIdxToFeatIdxMultiplier)
505 : {
506 : // Write root page (level 1)
507 18 : WriteUInt32(abyPage, 0); // id of next page at same level
508 36 : WriteUInt32(abyPage,
509 18 : nNumDirectChildren == 1 ? 1 : nNumDirectChildren - 1);
510 :
511 68 : for (int i = 0; i < nNumDirectChildren; i++)
512 : {
513 50 : WriteUInt32(abyPage, 2 + i); // id of subpage
514 : }
515 :
516 : // Add padding
517 18 : abyPage.resize(OFFSET_FIRST_VAL_IN_PAGE);
518 :
519 18 : if (nNumDirectChildren == 1)
520 : {
521 : // Should only happen if OPENFILEGDB_FORCE_SPX_DEPTH is forced
522 0 : writeValueFunc(abyPage, asValues.back().first, maxStrSize);
523 : }
524 : else
525 : {
526 50 : for (int i = 0; i < nNumDirectChildren - 1; i++)
527 : {
528 32 : const int nFeatIdx =
529 32 : (i + 1) * nSubPageIdxToFeatIdxMultiplier - 1;
530 32 : writeValueFunc(abyPage, asValues[nFeatIdx].first, maxStrSize);
531 : }
532 : }
533 18 : abyPage.resize(IDX_PAGE_SIZE);
534 18 : bRet &= VSIFWriteL(abyPage.data(), abyPage.size(), 1, fp) == 1;
535 : };
536 :
537 13339 : const auto WriteLeafPages = [=, &bRet, &asValues, &abyPage](
538 : int pageBaseOffset, int nNumFeaturePages)
539 : {
540 : // Write leaf pages
541 128 : for (int i = 0; i < nNumFeaturePages; ++i)
542 : {
543 110 : abyPage.clear();
544 110 : int nNumFeaturesInPage = numMaxFeaturesPerPage;
545 110 : if (i + 1 < nNumFeaturePages)
546 : {
547 92 : WriteUInt32(abyPage, pageBaseOffset + i +
548 : 1); // id of next page at same level
549 : }
550 : else
551 : {
552 18 : WriteUInt32(abyPage, 0);
553 18 : nNumFeaturesInPage = static_cast<int>(asValues.size()) -
554 18 : i * numMaxFeaturesPerPage;
555 : }
556 110 : CPLAssert(nNumFeaturesInPage > 0 &&
557 : nNumFeaturesInPage <= NUM_MAX_FEATURES_PER_PAGE);
558 110 : WriteUInt32(abyPage, nNumFeaturesInPage);
559 110 : WriteUInt32(abyPage, 0); // unknown semantics
560 :
561 : // Write features' ID
562 3143 : for (int j = 0; j < nNumFeaturesInPage; j++)
563 : {
564 3033 : WriteUInt32(
565 : abyPage,
566 : static_cast<uint32_t>(
567 3033 : asValues[i * numMaxFeaturesPerPage + j].second));
568 : }
569 :
570 : // Add padding
571 110 : abyPage.resize(OFFSET_FIRST_VAL_IN_PAGE);
572 :
573 : // Write features' spatial index value
574 3143 : for (int j = 0; j < nNumFeaturesInPage; j++)
575 : {
576 3033 : writeValueFunc(abyPage,
577 3033 : asValues[i * numMaxFeaturesPerPage + j].first,
578 : maxStrSize);
579 : }
580 :
581 110 : abyPage.resize(IDX_PAGE_SIZE);
582 110 : bRet &= VSIFWriteL(abyPage.data(), abyPage.size(), 1, fp) == 1;
583 : }
584 : };
585 :
586 145 : const auto WriteIntermediatePages =
587 804 : [=, &bRet, &asValues,
588 : &abyPage](int pageBaseOffset, int nNumPagesThisLevel,
589 : int nNumPagesNextLevel, int nSubPageIdxToFeatIdxMultiplier)
590 : {
591 68 : for (int i = 0; i < nNumPagesThisLevel; ++i)
592 : {
593 52 : abyPage.clear();
594 52 : int nNumItemsInPage = numMaxFeaturesPerPage;
595 52 : if (i + 1 < nNumPagesThisLevel)
596 : {
597 36 : WriteUInt32(abyPage, pageBaseOffset + i +
598 : 1); // id of next page at same level
599 : }
600 : else
601 : {
602 16 : WriteUInt32(abyPage, 0);
603 16 : nNumItemsInPage =
604 16 : nNumPagesNextLevel - i * numMaxFeaturesPerPage;
605 16 : CPLAssert(nNumItemsInPage > 1 &&
606 : nNumItemsInPage <= NUM_MAX_FEATURES_PER_PAGE + 1);
607 16 : nNumItemsInPage--;
608 : }
609 52 : CPLAssert(nNumItemsInPage > 0 &&
610 : nNumItemsInPage <= NUM_MAX_FEATURES_PER_PAGE);
611 52 : WriteUInt32(abyPage, nNumItemsInPage);
612 :
613 : // Write subpages' ID
614 200 : for (int j = 0; j < 1 + nNumItemsInPage; j++)
615 : {
616 148 : WriteUInt32(abyPage, pageBaseOffset + nNumPagesThisLevel +
617 148 : i * numMaxFeaturesPerPage + j);
618 : }
619 :
620 : // Add padding
621 52 : abyPage.resize(OFFSET_FIRST_VAL_IN_PAGE);
622 :
623 : // Write features' spatial index value
624 148 : for (int j = 0; j < nNumItemsInPage; j++)
625 : {
626 96 : const int nFeatIdx = (i * numMaxFeaturesPerPage + j + 1) *
627 : nSubPageIdxToFeatIdxMultiplier -
628 : 1;
629 96 : writeValueFunc(abyPage, asValues[nFeatIdx].first, maxStrSize);
630 : }
631 :
632 52 : abyPage.resize(IDX_PAGE_SIZE);
633 52 : bRet &= VSIFWriteL(abyPage.data(), abyPage.size(), 1, fp) == 1;
634 : }
635 : };
636 :
637 145 : const auto WriteLastTwoLevelPages =
638 24 : [numMaxFeaturesPerPage, WriteIntermediatePages,
639 : WriteLeafPages](int pageBaseOffset, int nNumPagesBeforeLastLevel,
640 : int nNumFeaturePages)
641 : {
642 : // Write pages at level depth-1 (referencing pages of level depth)
643 12 : WriteIntermediatePages(pageBaseOffset, nNumPagesBeforeLastLevel,
644 : nNumFeaturePages, numMaxFeaturesPerPage);
645 :
646 : // Write leaf pages
647 12 : WriteLeafPages(pageBaseOffset + nNumPagesBeforeLastLevel,
648 : nNumFeaturePages);
649 : };
650 :
651 285 : if (asValues.empty() || nDepth == 1 ||
652 140 : (nDepth == 0 &&
653 140 : static_cast<int>(asValues.size()) <= numMaxFeaturesPerPage))
654 : {
655 127 : nDepth = 1;
656 :
657 127 : WriteUInt32(abyPage, 0); // id of next page
658 127 : WriteUInt32(abyPage, static_cast<uint32_t>(asValues.size()));
659 127 : WriteUInt32(abyPage, 0); // unknown semantics
660 :
661 : // Write features' ID
662 1468 : for (const auto &pair : asValues)
663 1341 : WriteUInt32(abyPage, static_cast<uint32_t>(pair.second));
664 :
665 : // Add padding
666 127 : abyPage.resize(OFFSET_FIRST_VAL_IN_PAGE);
667 :
668 : // Write features' spatial index value
669 1468 : for (const auto &pair : asValues)
670 1341 : writeValueFunc(abyPage, pair.first, maxStrSize);
671 :
672 127 : abyPage.resize(IDX_PAGE_SIZE);
673 127 : bRet &= VSIFWriteL(abyPage.data(), abyPage.size(), 1, fp) == 1;
674 : }
675 36 : else if (nDepth == 2 || (nDepth == 0 && static_cast<int>(asValues.size()) <=
676 18 : (numMaxFeaturesPerPage + 1) *
677 : numMaxFeaturesPerPage))
678 : {
679 6 : nDepth = 2;
680 :
681 18 : const int nNumFeaturePages = static_cast<int>(
682 6 : DIV_ROUND_UP(asValues.size(), numMaxFeaturesPerPage));
683 6 : CPLAssert(nNumFeaturePages - 1 <= NUM_MAX_FEATURES_PER_PAGE);
684 :
685 : // Write root page (level 1)
686 6 : WriteRootPageNonLeaf(nNumFeaturePages, numMaxFeaturesPerPage);
687 :
688 : // Write leaf pages (level 2)
689 6 : WriteLeafPages(2, nNumFeaturePages);
690 : }
691 24 : else if (nDepth == 3 ||
692 12 : (nDepth == 0 &&
693 12 : static_cast<int>(asValues.size()) <=
694 12 : ((numMaxFeaturesPerPage + 1) * numMaxFeaturesPerPage + 1) *
695 : numMaxFeaturesPerPage))
696 : {
697 8 : nDepth = 3;
698 :
699 : // imagine simpler case: NUM_MAX_FEATURES_PER_PAGE = 2 and 9 values
700 : // ==> nNumFeaturePages = ceil(9 / 2) = 5
701 : // ==> nNumPagesLevel2 = ceil((5-1) / 2) = 2
702 : // level 1:
703 : // page 1: point to page 2(, 3)
704 : // level 2:
705 : // page 2: point to page 4, 5(, 6)
706 : // page 3: point to page 6, 7(, 8)
707 : // level 3:
708 : // page 4: point to feature 1, 2
709 : // page 5: point to feature 3, 4
710 : // page 6: point to feature 5, 6
711 : // page 7: point to feature 7, 8
712 : // page 8: point to feature 9
713 :
714 : // or NUM_MAX_FEATURES_PER_PAGE = 2 and 11 values
715 : // ==> nNumFeaturePages = ceil(11 / 2) = 6
716 : // ==> nNumPagesLevel2 = ceil((6-1) / 2) = 3
717 : // level 1:
718 : // page 1: point to page 2, 3(, 4)
719 : // level 2:
720 : // page 2: point to page 5, 6(, 7)
721 : // page 3: point to page 7, 8(, 9)
722 : // page 4: point to page 9(, 10)
723 : // level 3:
724 : // page 5: point to feature 1, 2
725 : // page 6: point to feature 3, 4
726 : // page 7: point to feature 5, 6
727 : // page 8: point to feature 7, 8
728 : // page 9: point to feature 9, 10
729 : // page 10: point to feature 11
730 :
731 : // or NUM_MAX_FEATURES_PER_PAGE = 2 and 14 values
732 : // ==> nNumFeaturePages = ceil(14 / 2) = 7
733 : // ==> nNumPagesLevel2 = ceil((7-1) / 2) = 3
734 : // level 1:
735 : // page 1: point to page 2, 3(, 4)
736 : // level 2:
737 : // page 2: point to page 5, 6(, 7)
738 : // page 3: point to page 7, 8(, 9)
739 : // page 4: point to page 9, 10(, 11)
740 : // level 3:
741 : // page 5: point to feature 1, 2
742 : // page 6: point to feature 3, 4
743 : // page 7: point to feature 5, 6
744 : // page 8: point to feature 7, 8
745 : // page 9: point to feature 9, 10
746 : // page 10: point to feature 11, 12
747 : // page 11: point to feature 13, 14
748 :
749 24 : const int nNumFeaturePages = static_cast<int>(
750 8 : DIV_ROUND_UP(asValues.size(), numMaxFeaturesPerPage));
751 16 : const int nNumPagesLevel2 =
752 : nNumFeaturePages == 1
753 : ? 1
754 8 : : DIV_ROUND_UP(nNumFeaturePages - 1, numMaxFeaturesPerPage);
755 8 : CPLAssert(nNumPagesLevel2 - 1 <= NUM_MAX_FEATURES_PER_PAGE);
756 :
757 : // Write root page (level 1)
758 8 : WriteRootPageNonLeaf(nNumPagesLevel2,
759 : numMaxFeaturesPerPage * numMaxFeaturesPerPage);
760 :
761 : // Write level 2 and level 3 pages
762 8 : WriteLastTwoLevelPages(2, nNumPagesLevel2, nNumFeaturePages);
763 : }
764 : else
765 : {
766 4 : nDepth = 4;
767 :
768 12 : const int nNumFeaturePages = static_cast<int>(
769 4 : DIV_ROUND_UP(asValues.size(), numMaxFeaturesPerPage));
770 8 : const int nNumPagesLevel3 =
771 : nNumFeaturePages == 1
772 : ? 1
773 4 : : DIV_ROUND_UP(nNumFeaturePages - 1, numMaxFeaturesPerPage);
774 8 : const int nNumPagesLevel2 =
775 : nNumPagesLevel3 == 1
776 : ? 1
777 4 : : DIV_ROUND_UP(nNumPagesLevel3 - 1, numMaxFeaturesPerPage);
778 4 : CPLAssert(nNumPagesLevel2 - 1 <= NUM_MAX_FEATURES_PER_PAGE);
779 :
780 : // Write root page (level 1)
781 4 : WriteRootPageNonLeaf(nNumPagesLevel2, numMaxFeaturesPerPage *
782 : numMaxFeaturesPerPage *
783 : numMaxFeaturesPerPage);
784 :
785 : // Write pages at level 2 (referencing pages of level 3)
786 4 : WriteIntermediatePages(2, nNumPagesLevel2, nNumPagesLevel3,
787 : numMaxFeaturesPerPage * numMaxFeaturesPerPage);
788 :
789 : // Write pages at level 3 and 4
790 4 : WriteLastTwoLevelPages(2 + nNumPagesLevel2, nNumPagesLevel3,
791 : nNumFeaturePages);
792 : }
793 :
794 : // Write trailer
795 145 : std::vector<GByte> abyTrailer;
796 145 : CPLAssert(SIZEOF_INDEXED_VALUE <= 255);
797 145 : WriteUInt8(abyTrailer, static_cast<uint8_t>(SIZEOF_INDEXED_VALUE));
798 145 : WriteUInt8(abyTrailer, maxStrSize ? 0x20 : 0x40); // unknown semantics
799 145 : WriteUInt32(abyTrailer, 1); // unknown semantics
800 145 : WriteUInt32(abyTrailer, nDepth); // index depth
801 145 : WriteUInt32(abyTrailer, static_cast<uint32_t>(asValues.size()));
802 145 : WriteUInt32(abyTrailer, 0); // unknown semantics
803 145 : WriteUInt32(abyTrailer, 1); // unknown semantics
804 145 : bRet &= VSIFWriteL(abyTrailer.data(), abyTrailer.size(), 1, fp) == 1;
805 :
806 145 : return bRet;
807 : }
808 :
809 : /************************************************************************/
810 : /* CreateSpatialIndex() */
811 : /************************************************************************/
812 :
813 193 : bool FileGDBTable::CreateSpatialIndex()
814 : {
815 386 : if (m_iGeomField < 0 || m_adfSpatialIndexGridResolution.empty() ||
816 193 : m_adfSpatialIndexGridResolution.size() > 3)
817 : {
818 0 : return false;
819 : }
820 :
821 193 : if (m_eTableGeomType == FGTGT_MULTIPATCH)
822 : {
823 0 : CPLError(CE_Failure, CPLE_NotSupported,
824 : "Multipatch not supported for spatial index generation");
825 0 : return false;
826 : }
827 :
828 : auto poGeomField =
829 193 : cpl::down_cast<FileGDBGeomField *>(m_apoFields[m_iGeomField].get());
830 193 : if (m_adfSpatialIndexGridResolution.size() == 1)
831 : {
832 : // Debug only
833 : const char *pszGridSize =
834 193 : CPLGetConfigOption("OPENFILEGDB_GRID_SIZE", nullptr);
835 193 : if (pszGridSize)
836 : {
837 0 : m_bDirtyGeomFieldSpatialIndexGridRes = true;
838 0 : m_adfSpatialIndexGridResolution[0] = CPLAtof(pszGridSize);
839 : poGeomField->m_adfSpatialIndexGridResolution =
840 0 : m_adfSpatialIndexGridResolution;
841 : }
842 : else
843 : {
844 193 : ComputeOptimalSpatialIndexGridResolution();
845 193 : if (m_adfSpatialIndexGridResolution[0] == 0)
846 65 : return false;
847 : }
848 : }
849 : auto poGeomConverter = std::unique_ptr<FileGDBOGRGeometryConverter>(
850 256 : FileGDBOGRGeometryConverter::BuildConverter(poGeomField));
851 : typedef std::pair<int64_t, int64_t> ValueOIDPair;
852 256 : std::vector<ValueOIDPair> asValues;
853 :
854 128 : const double dfGridStep = m_adfSpatialIndexGridResolution.back();
855 : const double dfShift =
856 128 : (1 << 29) / (dfGridStep / m_adfSpatialIndexGridResolution[0]);
857 :
858 : double dfYMinClamped;
859 : double dfYMaxClamped;
860 128 : GetMinMaxProjYForSpatialIndex(dfYMinClamped, dfYMaxClamped);
861 :
862 : const auto AddPointToIndex =
863 3813 : [this, dfGridStep, dfShift, dfYMinClamped, dfYMaxClamped](
864 7626 : double dfX, double dfY, std::vector<int64_t> &aSetValues)
865 : {
866 3813 : dfY = std::min(std::max(dfY, dfYMinClamped), dfYMaxClamped);
867 :
868 3813 : const double dfXShifted = dfX / dfGridStep + dfShift;
869 3813 : const double dfYShifted = dfY / dfGridStep + dfShift;
870 : // Each value must fit on 31 bit (sign included)
871 7626 : if (std::abs(dfXShifted) < (1 << 30) &&
872 3813 : std::abs(dfYShifted) < (1 << 30))
873 : {
874 3813 : const unsigned nX =
875 3813 : static_cast<unsigned>(static_cast<int>(std::floor(dfXShifted)));
876 3813 : const unsigned nY =
877 3813 : static_cast<unsigned>(static_cast<int>(std::floor(dfYShifted)));
878 : const uint64_t nVal =
879 3813 : ((static_cast<uint64_t>(m_adfSpatialIndexGridResolution.size() -
880 : 1))
881 3813 : << 62) |
882 3813 : ((static_cast<uint64_t>(nX)) << 31) | nY;
883 3813 : aSetValues.push_back(static_cast<int64_t>(nVal));
884 3813 : return true;
885 : }
886 0 : return false;
887 128 : };
888 :
889 : // Adapted from GDALdllImageLineAllTouched() of alg/llrasterize.cpp
890 : const auto AddLineStringToIndex =
891 269 : [this, dfGridStep, dfShift, dfYMinClamped, dfYMaxClamped](
892 16542 : const OGRLineString *poLS, std::vector<int64_t> &aSetValues)
893 : {
894 269 : const int nNumPoints = poLS->getNumPoints();
895 269 : if (nNumPoints < 2)
896 0 : return;
897 269 : OGREnvelope sEnvelope;
898 269 : poLS->getEnvelope(&sEnvelope);
899 269 : double dfYShift = 0;
900 269 : if (sEnvelope.MaxY > dfYMaxClamped)
901 30 : dfYShift = dfYMaxClamped - sEnvelope.MaxY;
902 239 : else if (sEnvelope.MinY < dfYMinClamped)
903 0 : dfYShift = dfYMinClamped - sEnvelope.MinY;
904 3382 : for (int i = 0; i < nNumPoints - 1; i++)
905 : {
906 3113 : double dfX = poLS->getX(i) / dfGridStep + dfShift;
907 3113 : double dfY = (poLS->getY(i) + dfYShift) / dfGridStep + dfShift;
908 3113 : double dfXEnd = poLS->getX(i + 1) / dfGridStep + dfShift;
909 : double dfYEnd =
910 3113 : (poLS->getY(i + 1) + dfYShift) / dfGridStep + dfShift;
911 6226 : if (!(std::abs(dfX) < (1 << 30) && std::abs(dfY) < (1 << 30) &&
912 3113 : std::abs(dfXEnd) < (1 << 30) && std::abs(dfYEnd) < (1 << 30)))
913 : {
914 0 : return;
915 : }
916 :
917 : // Swap if needed so we can proceed from left2right (X increasing)
918 3113 : if (dfX > dfXEnd)
919 : {
920 1312 : std::swap(dfX, dfXEnd);
921 1312 : std::swap(dfY, dfYEnd);
922 : }
923 :
924 : // Special case for vertical lines.
925 3113 : if (floor(dfX) == floor(dfXEnd) || fabs(dfX - dfXEnd) < .01)
926 : {
927 2784 : if (dfYEnd < dfY)
928 : {
929 1093 : std::swap(dfY, dfYEnd);
930 : }
931 :
932 2784 : const int iX = static_cast<int>(floor(dfXEnd));
933 2784 : int iY = static_cast<int>(floor(dfY));
934 2784 : int iYEnd = static_cast<int>(floor(dfYEnd));
935 :
936 5694 : for (; iY <= iYEnd; iY++)
937 : {
938 2910 : const unsigned nX = static_cast<unsigned>(iX);
939 2910 : const unsigned nY = static_cast<unsigned>(iY);
940 : const uint64_t nVal =
941 : ((static_cast<uint64_t>(
942 2910 : m_adfSpatialIndexGridResolution.size() - 1))
943 2910 : << 62) |
944 2910 : ((static_cast<uint64_t>(nX)) << 31) | nY;
945 2910 : aSetValues.push_back(static_cast<int64_t>(nVal));
946 : }
947 :
948 3049 : continue; // Next segment.
949 : }
950 :
951 : // Special case for horizontal lines.
952 329 : if (floor(dfY) == floor(dfYEnd) || fabs(dfY - dfYEnd) < .01)
953 : {
954 265 : if (dfXEnd < dfX)
955 : {
956 0 : std::swap(dfX, dfXEnd);
957 : }
958 :
959 265 : int iX = static_cast<int>(floor(dfX));
960 265 : const int iY = static_cast<int>(floor(dfY));
961 265 : int iXEnd = static_cast<int>(floor(dfXEnd));
962 :
963 795 : for (; iX <= iXEnd; iX++)
964 : {
965 530 : const unsigned nX = static_cast<unsigned>(iX);
966 530 : const unsigned nY = static_cast<unsigned>(iY);
967 : const uint64_t nVal =
968 : ((static_cast<uint64_t>(
969 530 : m_adfSpatialIndexGridResolution.size() - 1))
970 530 : << 62) |
971 530 : ((static_cast<uint64_t>(nX)) << 31) | nY;
972 530 : aSetValues.push_back(static_cast<int64_t>(nVal));
973 : }
974 :
975 265 : continue; // Next segment.
976 : }
977 :
978 : /* --------------------------------------------------------------------
979 : */
980 : /* General case - left to right sloped. */
981 : /* --------------------------------------------------------------------
982 : */
983 :
984 : // Recenter coordinates to avoid numeric precision issues
985 : // particularly the tests against a small epsilon below that could
986 : // lead to infinite looping otherwise.
987 64 : const int nXShift = static_cast<int>(floor(dfX));
988 64 : const int nYShift = static_cast<int>(floor(dfY));
989 64 : dfX -= nXShift;
990 64 : dfY -= nYShift;
991 64 : dfXEnd -= nXShift;
992 64 : dfYEnd -= nYShift;
993 :
994 64 : const double dfSlope = (dfYEnd - dfY) / (dfXEnd - dfX);
995 :
996 : // Step from pixel to pixel.
997 176 : while (dfX < dfXEnd)
998 : {
999 112 : const int iX = static_cast<int>(floor(dfX));
1000 112 : const int iY = static_cast<int>(floor(dfY));
1001 :
1002 : // Burn in the current point.
1003 112 : const unsigned nX = static_cast<unsigned>(iX + nXShift);
1004 112 : const unsigned nY = static_cast<unsigned>(iY + nYShift);
1005 : const uint64_t nVal =
1006 : ((static_cast<uint64_t>(
1007 112 : m_adfSpatialIndexGridResolution.size() - 1))
1008 112 : << 62) |
1009 112 : ((static_cast<uint64_t>(nX)) << 31) | nY;
1010 112 : aSetValues.push_back(static_cast<int64_t>(nVal));
1011 :
1012 112 : double dfStepX = floor(dfX + 1.0) - dfX;
1013 112 : double dfStepY = dfStepX * dfSlope;
1014 :
1015 : // Step to right pixel without changing scanline?
1016 112 : if (static_cast<int>(floor(dfY + dfStepY)) == iY)
1017 : {
1018 40 : dfX += dfStepX;
1019 40 : dfY += dfStepY;
1020 : }
1021 72 : else if (dfSlope < 0)
1022 : {
1023 2 : dfStepY = iY - dfY;
1024 2 : if (dfStepY > -0.000000001)
1025 2 : dfStepY = -0.000000001;
1026 :
1027 2 : dfStepX = dfStepY / dfSlope;
1028 2 : dfX += dfStepX;
1029 2 : dfY += dfStepY;
1030 : }
1031 : else
1032 : {
1033 70 : dfStepY = (iY + 1) - dfY;
1034 70 : if (dfStepY < 0.000000001)
1035 0 : dfStepY = 0.000000001;
1036 :
1037 70 : dfStepX = dfStepY / dfSlope;
1038 70 : dfX += dfStepX;
1039 70 : dfY += dfStepY;
1040 : }
1041 : } // Next step along segment.
1042 : }
1043 128 : };
1044 :
1045 : // Adapted from GDALdllImageFilledPolygon() of alg/llrasterize.cpp
1046 : const auto AddPolygonToIndex =
1047 184 : [this, dfGridStep, dfShift, dfYMinClamped, dfYMaxClamped,
1048 : AddLineStringToIndex](const OGRPolygon *poPoly,
1049 9707 : std::vector<int64_t> &aSetValues)
1050 : {
1051 184 : if (poPoly->IsEmpty())
1052 0 : return;
1053 :
1054 : // Burn contour of exterior ring, because burning the interior
1055 : // can often result in nothing
1056 184 : AddLineStringToIndex(poPoly->getExteriorRing(), aSetValues);
1057 :
1058 184 : OGREnvelope sEnvelope;
1059 184 : poPoly->getEnvelope(&sEnvelope);
1060 :
1061 184 : double dfYShift = 0;
1062 184 : if (sEnvelope.MaxY > dfYMaxClamped)
1063 30 : dfYShift = dfYMaxClamped - sEnvelope.MaxY;
1064 154 : else if (sEnvelope.MinY < dfYMinClamped)
1065 0 : dfYShift = dfYMinClamped - sEnvelope.MinY;
1066 :
1067 184 : const int miny = static_cast<int>(
1068 184 : floor((sEnvelope.MinY + dfYShift) / dfGridStep + dfShift));
1069 184 : const int maxy = static_cast<int>(
1070 184 : floor((sEnvelope.MaxY + dfYShift) / dfGridStep + dfShift));
1071 368 : std::vector<double> intersections;
1072 :
1073 : // Burn interior of polygon
1074 438 : for (int iY = miny; iY <= maxy; iY++)
1075 : {
1076 254 : const double dy = iY + 0.5;
1077 254 : intersections.clear();
1078 539 : for (const auto *poRing : *poPoly)
1079 : {
1080 285 : const int nNumPoints = poRing->getNumPoints();
1081 4513 : for (int i = 0; i < nNumPoints - 1; ++i)
1082 : {
1083 : double dy1 =
1084 4228 : (poRing->getY(i) + dfYShift) / dfGridStep + dfShift;
1085 : double dy2 =
1086 4228 : (poRing->getY(i + 1) + dfYShift) / dfGridStep + dfShift;
1087 4228 : if ((dy1 < dy && dy2 < dy) || (dy1 > dy && dy2 > dy))
1088 3974 : continue;
1089 :
1090 256 : double dx1 = 0.0;
1091 256 : double dx2 = 0.0;
1092 256 : if (dy1 < dy2)
1093 : {
1094 126 : dx1 = poRing->getX(i) / dfGridStep + dfShift;
1095 126 : dx2 = poRing->getX(i + 1) / dfGridStep + dfShift;
1096 : }
1097 130 : else if (dy1 > dy2)
1098 : {
1099 128 : std::swap(dy1, dy2);
1100 128 : dx2 = poRing->getX(i) / dfGridStep + dfShift;
1101 128 : dx1 = poRing->getX(i + 1) / dfGridStep + dfShift;
1102 : }
1103 : else // if( fabs(dy1-dy2) < 1.e-6 )
1104 : {
1105 2 : const int iX1 = static_cast<int>(floor(
1106 2 : std::min(poRing->getX(i), poRing->getX(i + 1)) /
1107 : dfGridStep +
1108 2 : dfShift));
1109 2 : const int iX2 = static_cast<int>(floor(
1110 2 : std::max(poRing->getX(i), poRing->getX(i + 1)) /
1111 : dfGridStep +
1112 2 : dfShift));
1113 :
1114 : // Fill the horizontal segment (separately from the
1115 : // rest).
1116 6 : for (int iX = iX1; iX <= iX2; ++iX)
1117 : {
1118 4 : const unsigned nX = static_cast<unsigned>(iX);
1119 4 : const unsigned nY = static_cast<unsigned>(iY);
1120 : const uint64_t nVal =
1121 : ((static_cast<uint64_t>(
1122 4 : m_adfSpatialIndexGridResolution.size() -
1123 : 1))
1124 4 : << 62) |
1125 4 : ((static_cast<uint64_t>(nX)) << 31) | nY;
1126 4 : aSetValues.push_back(static_cast<int64_t>(nVal));
1127 : }
1128 2 : continue;
1129 : }
1130 :
1131 254 : if (dy < dy2 && dy >= dy1)
1132 : {
1133 250 : const double intersect =
1134 250 : (dy - dy1) * (dx2 - dx1) / (dy2 - dy1) + dx1;
1135 250 : intersections.emplace_back(intersect);
1136 : }
1137 : }
1138 : }
1139 :
1140 254 : std::sort(intersections.begin(), intersections.end());
1141 :
1142 379 : for (size_t i = 0; i + 1 < intersections.size(); i += 2)
1143 : {
1144 125 : const int iX1 = static_cast<int>(floor(intersections[i]));
1145 125 : const int iX2 = static_cast<int>(floor(intersections[i + 1]));
1146 :
1147 308 : for (int iX = iX1; iX <= iX2; ++iX)
1148 : {
1149 183 : const unsigned nX = static_cast<unsigned>(iX);
1150 183 : const unsigned nY = static_cast<unsigned>(iY);
1151 : const uint64_t nVal =
1152 : ((static_cast<uint64_t>(
1153 183 : m_adfSpatialIndexGridResolution.size() - 1))
1154 183 : << 62) |
1155 183 : ((static_cast<uint64_t>(nX)) << 31) | nY;
1156 183 : aSetValues.push_back(static_cast<int64_t>(nVal));
1157 : }
1158 : }
1159 : }
1160 128 : };
1161 :
1162 256 : std::vector<int64_t> aSetValues;
1163 128 : int64_t iLastReported = 0;
1164 128 : const auto nReportIncrement = m_nTotalRecordCount / 20;
1165 : try
1166 : {
1167 4174 : for (int64_t iCurFeat = 0; iCurFeat < m_nTotalRecordCount; ++iCurFeat)
1168 : {
1169 4046 : if (m_nTotalRecordCount > 10000 &&
1170 0 : (iCurFeat + 1 == m_nTotalRecordCount ||
1171 0 : iCurFeat - iLastReported >= nReportIncrement))
1172 : {
1173 0 : CPLDebug("OpenFileGDB", "Spatial index building: %02.2f %%",
1174 0 : 100 * double(iCurFeat + 1) /
1175 0 : double(m_nTotalRecordCount));
1176 0 : iLastReported = iCurFeat;
1177 : }
1178 4046 : iCurFeat = GetAndSelectNextNonEmptyRow(iCurFeat);
1179 4046 : if (iCurFeat < 0)
1180 0 : break;
1181 4046 : const OGRField *psField = GetFieldValue(m_iGeomField);
1182 4046 : if (psField != nullptr)
1183 : {
1184 : auto poGeom = std::unique_ptr<OGRGeometry>(
1185 8042 : poGeomConverter->GetAsGeometry(psField));
1186 4021 : if (poGeom != nullptr && !poGeom->IsEmpty())
1187 : {
1188 4017 : aSetValues.clear();
1189 : const auto eGeomType =
1190 4017 : wkbFlatten(poGeom->getGeometryType());
1191 4017 : if (eGeomType == wkbPoint)
1192 : {
1193 3761 : const auto poPoint = poGeom->toPoint();
1194 3761 : AddPointToIndex(poPoint->getX(), poPoint->getY(),
1195 : aSetValues);
1196 : }
1197 256 : else if (eGeomType == wkbMultiPoint)
1198 : {
1199 76 : for (const auto poPoint : *(poGeom->toMultiPoint()))
1200 : {
1201 52 : AddPointToIndex(poPoint->getX(), poPoint->getY(),
1202 : aSetValues);
1203 : }
1204 : }
1205 232 : else if (eGeomType == wkbLineString)
1206 : {
1207 46 : AddLineStringToIndex(poGeom->toLineString(),
1208 : aSetValues);
1209 : }
1210 186 : else if (eGeomType == wkbMultiLineString)
1211 : {
1212 48 : for (const auto poLS : *(poGeom->toMultiLineString()))
1213 : {
1214 32 : AddLineStringToIndex(poLS, aSetValues);
1215 : }
1216 : }
1217 170 : else if (eGeomType == wkbCircularString ||
1218 : eGeomType == wkbCompoundCurve)
1219 : {
1220 4 : poGeom.reset(poGeom->getLinearGeometry());
1221 4 : if (poGeom)
1222 4 : AddLineStringToIndex(poGeom->toLineString(),
1223 : aSetValues);
1224 : }
1225 166 : else if (eGeomType == wkbMultiCurve)
1226 : {
1227 1 : poGeom.reset(poGeom->getLinearGeometry());
1228 1 : if (poGeom)
1229 : {
1230 3 : for (const auto poLS :
1231 4 : *(poGeom->toMultiLineString()))
1232 : {
1233 3 : AddLineStringToIndex(poLS, aSetValues);
1234 : }
1235 : }
1236 : }
1237 165 : else if (eGeomType == wkbPolygon)
1238 : {
1239 140 : AddPolygonToIndex(poGeom->toPolygon(), aSetValues);
1240 : }
1241 25 : else if (eGeomType == wkbCurvePolygon)
1242 : {
1243 6 : poGeom.reset(poGeom->getLinearGeometry());
1244 6 : if (poGeom)
1245 6 : AddPolygonToIndex(poGeom->toPolygon(), aSetValues);
1246 : }
1247 19 : else if (eGeomType == wkbMultiPolygon)
1248 : {
1249 54 : for (const auto poPoly : *(poGeom->toMultiPolygon()))
1250 : {
1251 36 : AddPolygonToIndex(poPoly, aSetValues);
1252 : }
1253 : }
1254 1 : else if (eGeomType == wkbMultiSurface)
1255 : {
1256 1 : poGeom.reset(poGeom->getLinearGeometry());
1257 1 : if (poGeom)
1258 : {
1259 2 : for (const auto poPoly :
1260 3 : *(poGeom->toMultiPolygon()))
1261 : {
1262 2 : AddPolygonToIndex(poPoly, aSetValues);
1263 : }
1264 : }
1265 : }
1266 :
1267 4017 : std::sort(aSetValues.begin(), aSetValues.end());
1268 :
1269 4017 : int64_t nLastVal = std::numeric_limits<int64_t>::min();
1270 11569 : for (auto nVal : aSetValues)
1271 : {
1272 7552 : if (nVal != nLastVal)
1273 : {
1274 4376 : asValues.push_back(
1275 4376 : ValueOIDPair(nVal, iCurFeat + 1));
1276 4376 : nLastVal = nVal;
1277 : }
1278 : }
1279 : }
1280 : }
1281 : }
1282 : }
1283 0 : catch (const std::exception &e)
1284 : {
1285 0 : CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
1286 0 : return false;
1287 : }
1288 :
1289 : const std::string osSPXFilename(
1290 256 : CPLResetExtensionSafe(m_osFilename.c_str(), "spx"));
1291 128 : VSILFILE *fp = VSIFOpenL(osSPXFilename.c_str(), "wb");
1292 128 : if (fp == nullptr)
1293 0 : return false;
1294 :
1295 : // Configurable only for debugging purposes
1296 128 : int nDepth = atoi(CPLGetConfigOption("OPENFILEGDB_FORCE_SPX_DEPTH", "0"));
1297 :
1298 128 : const auto writeValueFunc =
1299 4473 : +[](std::vector<GByte> &abyPage,
1300 : const typename ValueOIDPair::first_type &nval, int /* maxStrSize */)
1301 4473 : { WriteUInt64(abyPage, static_cast<uint64_t>(nval)); };
1302 :
1303 128 : bool bRet = WriteIndex(fp, asValues, writeValueFunc, nDepth);
1304 :
1305 128 : CPLDebug("OpenFileGDB", "Spatial index of depth %d", nDepth);
1306 :
1307 128 : VSIFCloseL(fp);
1308 :
1309 128 : if (!bRet)
1310 : {
1311 1 : CPLError(CE_Failure, CPLE_FileIO, "Write error during .spx generation");
1312 1 : VSIUnlink(osSPXFilename.c_str());
1313 : }
1314 :
1315 128 : return bRet;
1316 : }
1317 :
1318 : /************************************************************************/
1319 : /* CreateAttributeIndex() */
1320 : /************************************************************************/
1321 :
1322 18 : bool FileGDBTable::CreateAttributeIndex(const FileGDBIndex *poIndex)
1323 : {
1324 36 : const std::string osFieldName = poIndex->GetFieldName();
1325 18 : const int iField = GetFieldIdx(osFieldName);
1326 18 : if (iField < 0)
1327 : {
1328 0 : CPLError(CE_Failure, CPLE_AppDefined, "Cannot find field %s",
1329 : osFieldName.c_str());
1330 0 : return false;
1331 : }
1332 :
1333 : const std::string osIdxFilename(CPLResetExtensionSafe(
1334 36 : m_osFilename.c_str(), (poIndex->GetIndexName() + ".atx").c_str()));
1335 18 : VSILFILE *fp = VSIFOpenL(osIdxFilename.c_str(), "wb");
1336 18 : if (fp == nullptr)
1337 0 : return false;
1338 :
1339 : bool bRet;
1340 18 : int nDepth = 0;
1341 :
1342 : try
1343 : {
1344 18 : const auto eFieldType = m_apoFields[iField]->GetType();
1345 18 : if (eFieldType == FGFT_INT16)
1346 : {
1347 : typedef std::pair<int16_t, int64_t> ValueOIDPair;
1348 2 : std::vector<ValueOIDPair> asValues;
1349 8 : for (int64_t iCurFeat = 0; iCurFeat < m_nTotalRecordCount;
1350 : ++iCurFeat)
1351 : {
1352 6 : iCurFeat = GetAndSelectNextNonEmptyRow(iCurFeat);
1353 6 : if (iCurFeat < 0)
1354 0 : break;
1355 6 : const OGRField *psField = GetFieldValue(iField);
1356 6 : if (psField != nullptr)
1357 : {
1358 2 : asValues.push_back(ValueOIDPair(
1359 4 : static_cast<int16_t>(psField->Integer), iCurFeat + 1));
1360 : }
1361 : }
1362 :
1363 2 : const auto writeValueFunc =
1364 2 : +[](std::vector<GByte> &abyPage,
1365 : const typename ValueOIDPair::first_type &val,
1366 : int /* maxStrSize */)
1367 2 : { WriteUInt16(abyPage, static_cast<uint16_t>(val)); };
1368 :
1369 2 : bRet = WriteIndex(fp, asValues, writeValueFunc, nDepth);
1370 : }
1371 16 : else if (eFieldType == FGFT_INT32)
1372 : {
1373 : typedef std::pair<int32_t, int64_t> ValueOIDPair;
1374 2 : std::vector<ValueOIDPair> asValues;
1375 8 : for (int64_t iCurFeat = 0; iCurFeat < m_nTotalRecordCount;
1376 : ++iCurFeat)
1377 : {
1378 6 : iCurFeat = GetAndSelectNextNonEmptyRow(iCurFeat);
1379 6 : if (iCurFeat < 0)
1380 0 : break;
1381 6 : const OGRField *psField = GetFieldValue(iField);
1382 6 : if (psField != nullptr)
1383 : {
1384 2 : asValues.push_back(
1385 4 : ValueOIDPair(psField->Integer, iCurFeat + 1));
1386 : }
1387 : }
1388 :
1389 2 : const auto writeValueFunc =
1390 2 : +[](std::vector<GByte> &abyPage,
1391 : const typename ValueOIDPair::first_type &val,
1392 : int /* maxStrSize */)
1393 2 : { WriteUInt32(abyPage, static_cast<uint32_t>(val)); };
1394 :
1395 2 : bRet = WriteIndex(fp, asValues, writeValueFunc, nDepth);
1396 : }
1397 14 : else if (eFieldType == FGFT_INT64)
1398 : {
1399 : typedef std::pair<int64_t, int64_t> ValueOIDPair;
1400 1 : std::vector<ValueOIDPair> asValues;
1401 3 : for (int64_t iCurFeat = 0; iCurFeat < m_nTotalRecordCount;
1402 : ++iCurFeat)
1403 : {
1404 2 : iCurFeat = GetAndSelectNextNonEmptyRow(iCurFeat);
1405 2 : if (iCurFeat < 0)
1406 0 : break;
1407 2 : const OGRField *psField = GetFieldValue(iField);
1408 2 : if (psField != nullptr)
1409 : {
1410 2 : asValues.push_back(
1411 4 : ValueOIDPair(psField->Integer64, iCurFeat + 1));
1412 : }
1413 : }
1414 :
1415 1 : const auto writeValueFunc =
1416 2 : +[](std::vector<GByte> &abyPage,
1417 : const typename ValueOIDPair::first_type &val,
1418 : int /* maxStrSize */)
1419 2 : { WriteUInt64(abyPage, static_cast<uint64_t>(val)); };
1420 :
1421 1 : bRet = WriteIndex(fp, asValues, writeValueFunc, nDepth);
1422 : }
1423 13 : else if (eFieldType == FGFT_FLOAT32)
1424 : {
1425 : typedef std::pair<float, int64_t> ValueOIDPair;
1426 2 : std::vector<ValueOIDPair> asValues;
1427 8 : for (int64_t iCurFeat = 0; iCurFeat < m_nTotalRecordCount;
1428 : ++iCurFeat)
1429 : {
1430 6 : iCurFeat = GetAndSelectNextNonEmptyRow(iCurFeat);
1431 6 : if (iCurFeat < 0)
1432 0 : break;
1433 6 : const OGRField *psField = GetFieldValue(iField);
1434 6 : if (psField != nullptr)
1435 : {
1436 2 : asValues.push_back(ValueOIDPair(
1437 4 : static_cast<float>(psField->Real), iCurFeat + 1));
1438 : }
1439 : }
1440 :
1441 2 : const auto writeValueFunc =
1442 2 : +[](std::vector<GByte> &abyPage,
1443 : const typename ValueOIDPair::first_type &val,
1444 2 : int /* maxStrSize */) { WriteFloat32(abyPage, val); };
1445 :
1446 2 : bRet = WriteIndex(fp, asValues, writeValueFunc, nDepth);
1447 : }
1448 11 : else if (eFieldType == FGFT_FLOAT64 || eFieldType == FGFT_DATETIME ||
1449 6 : eFieldType == FGFT_DATE || eFieldType == FGFT_TIME ||
1450 : eFieldType == FGFT_DATETIME_WITH_OFFSET)
1451 : {
1452 : typedef std::pair<double, int64_t> ValueOIDPair;
1453 7 : std::vector<ValueOIDPair> asValues;
1454 : // Hack to force reading DateTime as double
1455 7 : m_apoFields[iField]->m_bReadAsDouble = true;
1456 28 : for (int64_t iCurFeat = 0; iCurFeat < m_nTotalRecordCount;
1457 : ++iCurFeat)
1458 : {
1459 21 : iCurFeat = GetAndSelectNextNonEmptyRow(iCurFeat);
1460 21 : if (iCurFeat < 0)
1461 0 : break;
1462 21 : const OGRField *psField = GetFieldValue(iField);
1463 21 : if (psField != nullptr)
1464 : {
1465 13 : asValues.push_back(
1466 26 : ValueOIDPair(psField->Real, iCurFeat + 1));
1467 : }
1468 : }
1469 7 : m_apoFields[iField]->m_bReadAsDouble = false;
1470 :
1471 7 : const auto writeValueFunc =
1472 13 : +[](std::vector<GByte> &abyPage,
1473 : const typename ValueOIDPair::first_type &val,
1474 13 : int /* maxStrSize */) { WriteFloat64(abyPage, val); };
1475 :
1476 7 : bRet = WriteIndex(fp, asValues, writeValueFunc, nDepth);
1477 : }
1478 4 : else if (eFieldType == FGFT_STRING)
1479 : {
1480 : typedef std::pair<std::vector<std::uint16_t>, int64_t> ValueOIDPair;
1481 8 : std::vector<ValueOIDPair> asValues;
1482 4 : bRet = true;
1483 : const bool bIsLower =
1484 4 : STARTS_WITH_CI(poIndex->GetExpression().c_str(), "LOWER(");
1485 4 : int maxStrSize = 0;
1486 16 : for (int64_t iCurFeat = 0; iCurFeat < m_nTotalRecordCount;
1487 : ++iCurFeat)
1488 : {
1489 12 : iCurFeat = GetAndSelectNextNonEmptyRow(iCurFeat);
1490 12 : if (iCurFeat < 0)
1491 0 : break;
1492 12 : const OGRField *psField = GetFieldValue(iField);
1493 12 : if (psField != nullptr)
1494 : {
1495 16 : wchar_t *pWide = CPLRecodeToWChar(
1496 8 : psField->String, CPL_ENC_UTF8, CPL_ENC_UCS2);
1497 8 : if (pWide == nullptr)
1498 : {
1499 0 : bRet = false;
1500 0 : break;
1501 : }
1502 8 : int nCount = 0;
1503 8 : std::vector<std::uint16_t> asUTF16Str;
1504 352 : while (nCount < MAX_CAR_COUNT_INDEXED_STR &&
1505 348 : pWide[nCount] != 0)
1506 : {
1507 344 : nCount++;
1508 : }
1509 8 : if (nCount > maxStrSize)
1510 6 : maxStrSize = nCount;
1511 8 : asUTF16Str.reserve(nCount);
1512 352 : for (int i = 0; i < nCount; ++i)
1513 : {
1514 344 : if (bIsLower && pWide[i] >= 'A' && pWide[i] <= 'Z')
1515 10 : asUTF16Str.push_back(
1516 10 : static_cast<uint16_t>(pWide[i] + ('a' - 'A')));
1517 : else
1518 334 : asUTF16Str.push_back(
1519 334 : static_cast<uint16_t>(pWide[i]));
1520 : }
1521 8 : CPLFree(pWide);
1522 :
1523 8 : asValues.push_back(ValueOIDPair(asUTF16Str, iCurFeat + 1));
1524 : }
1525 : }
1526 4 : if (maxStrSize < MAX_CAR_COUNT_INDEXED_STR)
1527 2 : maxStrSize++;
1528 :
1529 4 : const auto writeValueFunc =
1530 8 : +[](std::vector<GByte> &abyPage,
1531 : const typename ValueOIDPair::first_type &val,
1532 : int l_maxStrSize)
1533 : {
1534 352 : for (size_t i = 0; i < val.size(); ++i)
1535 344 : WriteUInt16(abyPage, val[i]);
1536 158 : for (size_t i = val.size();
1537 158 : i < static_cast<size_t>(l_maxStrSize); ++i)
1538 150 : WriteUInt16(abyPage, 32); // space
1539 8 : };
1540 :
1541 4 : if (bRet)
1542 : {
1543 4 : bRet = WriteIndex(fp, asValues, writeValueFunc, nDepth,
1544 : maxStrSize);
1545 : }
1546 : }
1547 : else
1548 : {
1549 0 : CPLError(CE_Failure, CPLE_AppDefined,
1550 : "CreateAttributeIndex(%s): "
1551 : "Unsupported field type for index creation",
1552 : osFieldName.c_str());
1553 0 : bRet = false;
1554 : }
1555 : }
1556 0 : catch (const std::exception &e)
1557 : {
1558 0 : CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
1559 0 : bRet = false;
1560 : }
1561 :
1562 18 : VSIFCloseL(fp);
1563 :
1564 18 : if (!bRet)
1565 : {
1566 0 : CPLError(CE_Failure, CPLE_FileIO, "Write error during %s generation",
1567 : osIdxFilename.c_str());
1568 0 : VSIUnlink(osIdxFilename.c_str());
1569 : }
1570 :
1571 18 : return bRet;
1572 : }
1573 :
1574 : } /* namespace OpenFileGDB */
|
