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