Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: GeoPackage Translator
4 : * Purpose: Implements GDALGeoPackageRasterBand class
5 : * Author: Even Rouault <even dot rouault at spatialys dot com>
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2014, Even Rouault <even dot rouault at spatialys dot com>
9 : *
10 : * SPDX-License-Identifier: MIT
11 : ****************************************************************************/
12 :
13 : #include "ogr_geopackage.h"
14 : #include "memdataset.h"
15 : #include "gdal_alg_priv.h"
16 : #include "ogrsqlitevfs.h"
17 : #include "cpl_error.h"
18 : #include "cpl_float.h"
19 :
20 : #include <algorithm>
21 : #include <cmath>
22 : #include <limits>
23 : #include <set>
24 : #include <utility>
25 :
26 : #if !defined(DEBUG_VERBOSE) && defined(DEBUG_VERBOSE_GPKG)
27 : #define DEBUG_VERBOSE
28 : #endif
29 :
30 : /************************************************************************/
31 : /* GDALGPKGMBTilesLikePseudoDataset() */
32 : /************************************************************************/
33 :
34 2649 : GDALGPKGMBTilesLikePseudoDataset::GDALGPKGMBTilesLikePseudoDataset()
35 : : m_bForceTempDBCompaction(
36 2649 : CPLTestBool(CPLGetConfigOption("GPKG_FORCE_TEMPDB_COMPACTION", "NO")))
37 : {
38 13245 : for (int i = 0; i < 4; i++)
39 : {
40 10596 : m_asCachedTilesDesc[i].nRow = -1;
41 10596 : m_asCachedTilesDesc[i].nCol = -1;
42 10596 : m_asCachedTilesDesc[i].nIdxWithinTileData = -1;
43 10596 : m_asCachedTilesDesc[i].abBandDirty[0] = FALSE;
44 10596 : m_asCachedTilesDesc[i].abBandDirty[1] = FALSE;
45 10596 : m_asCachedTilesDesc[i].abBandDirty[2] = FALSE;
46 10596 : m_asCachedTilesDesc[i].abBandDirty[3] = FALSE;
47 : }
48 2649 : }
49 :
50 : /************************************************************************/
51 : /* ~GDALGPKGMBTilesLikePseudoDataset() */
52 : /************************************************************************/
53 :
54 2649 : GDALGPKGMBTilesLikePseudoDataset::~GDALGPKGMBTilesLikePseudoDataset()
55 : {
56 2649 : if (m_poParentDS == nullptr && m_hTempDB != nullptr)
57 : {
58 48 : sqlite3_close(m_hTempDB);
59 48 : m_hTempDB = nullptr;
60 48 : VSIUnlink(m_osTempDBFilename);
61 48 : if (m_pMyVFS)
62 : {
63 29 : sqlite3_vfs_unregister(m_pMyVFS);
64 29 : CPLFree(m_pMyVFS->pAppData);
65 29 : CPLFree(m_pMyVFS);
66 : }
67 : }
68 2649 : CPLFree(m_pabyCachedTiles);
69 2649 : delete m_poCT;
70 2649 : CPLFree(m_pabyHugeColorArray);
71 2649 : }
72 :
73 : /************************************************************************/
74 : /* SetDataType() */
75 : /************************************************************************/
76 :
77 319 : void GDALGPKGMBTilesLikePseudoDataset::SetDataType(GDALDataType eDT)
78 : {
79 319 : CPLAssert(eDT == GDT_Byte || eDT == GDT_Int16 || eDT == GDT_UInt16 ||
80 : eDT == GDT_Float32);
81 319 : m_eDT = eDT;
82 319 : m_nDTSize = GDALGetDataTypeSizeBytes(m_eDT);
83 319 : }
84 :
85 : /************************************************************************/
86 : /* SetGlobalOffsetScale() */
87 : /************************************************************************/
88 :
89 75 : void GDALGPKGMBTilesLikePseudoDataset::SetGlobalOffsetScale(double dfOffset,
90 : double dfScale)
91 : {
92 75 : m_dfOffset = dfOffset;
93 75 : m_dfScale = dfScale;
94 75 : }
95 :
96 : /************************************************************************/
97 : /* GDALGPKGMBTilesLikeRasterBand() */
98 : /************************************************************************/
99 :
100 : // Recent GCC versions complain about null dereference of m_poTPD in
101 : // the constructor body
102 : #ifdef __GNUC__
103 : #pragma GCC diagnostic push
104 : #pragma GCC diagnostic ignored "-Wnull-dereference"
105 : #endif
106 :
107 2551 : GDALGPKGMBTilesLikeRasterBand::GDALGPKGMBTilesLikeRasterBand(
108 2551 : GDALGPKGMBTilesLikePseudoDataset *poTPD, int nTileWidth, int nTileHeight)
109 2551 : : m_poTPD(poTPD)
110 : {
111 2551 : eDataType = m_poTPD->m_eDT;
112 2551 : m_nDTSize = m_poTPD->m_nDTSize;
113 2551 : nBlockXSize = nTileWidth;
114 2551 : nBlockYSize = nTileHeight;
115 2551 : }
116 :
117 : #ifdef __GNUC__
118 : #pragma GCC diagnostic pop
119 : #endif
120 :
121 : /************************************************************************/
122 : /* FlushCache() */
123 : /************************************************************************/
124 :
125 5477 : CPLErr GDALGPKGMBTilesLikeRasterBand::FlushCache(bool bAtClosing)
126 : {
127 5477 : m_poTPD->m_nLastSpaceCheckTimestamp = -1; // disable partial flushes
128 5477 : CPLErr eErr = GDALPamRasterBand::FlushCache(bAtClosing);
129 5477 : if (eErr == CE_None)
130 5474 : eErr = m_poTPD->IFlushCacheWithErrCode(bAtClosing);
131 5477 : m_poTPD->m_nLastSpaceCheckTimestamp = 0;
132 5477 : return eErr;
133 : }
134 :
135 : /************************************************************************/
136 : /* FlushTiles() */
137 : /************************************************************************/
138 :
139 3383 : CPLErr GDALGPKGMBTilesLikePseudoDataset::FlushTiles()
140 : {
141 3383 : CPLErr eErr = CE_None;
142 3383 : GDALGPKGMBTilesLikePseudoDataset *poMainDS =
143 3383 : m_poParentDS ? m_poParentDS : this;
144 3383 : if (poMainDS->m_nTileInsertionCount < 0)
145 7 : return CE_Failure;
146 :
147 3376 : if (IGetUpdate())
148 : {
149 2017 : if (m_nShiftXPixelsMod || m_nShiftYPixelsMod)
150 : {
151 607 : eErr = FlushRemainingShiftedTiles(false /* total flush*/);
152 : }
153 : else
154 : {
155 1410 : eErr = WriteTile();
156 : }
157 : }
158 :
159 3376 : if (poMainDS->m_nTileInsertionCount > 0)
160 : {
161 188 : if (poMainDS->ICommitTransaction() != OGRERR_NONE)
162 : {
163 0 : poMainDS->m_nTileInsertionCount = -1;
164 0 : eErr = CE_Failure;
165 : }
166 : else
167 : {
168 188 : poMainDS->m_nTileInsertionCount = 0;
169 : }
170 : }
171 3376 : return eErr;
172 : }
173 :
174 : /************************************************************************/
175 : /* GetColorTable() */
176 : /************************************************************************/
177 :
178 345 : GDALColorTable *GDALGPKGMBTilesLikeRasterBand::GetColorTable()
179 : {
180 345 : if (poDS->GetRasterCount() != 1)
181 127 : return nullptr;
182 :
183 218 : if (!m_poTPD->m_bTriedEstablishingCT)
184 : {
185 64 : m_poTPD->m_bTriedEstablishingCT = true;
186 64 : if (m_poTPD->m_poParentDS != nullptr)
187 : {
188 16 : m_poTPD->m_poCT =
189 8 : m_poTPD->m_poParentDS->IGetRasterBand(1)->GetColorTable();
190 8 : if (m_poTPD->m_poCT)
191 4 : m_poTPD->m_poCT = m_poTPD->m_poCT->Clone();
192 8 : return m_poTPD->m_poCT;
193 : }
194 :
195 68 : for (int i = 0; i < 2; i++)
196 : {
197 62 : bool bRetry = false;
198 62 : char *pszSQL = nullptr;
199 62 : if (i == 0)
200 : {
201 56 : pszSQL = sqlite3_mprintf("SELECT tile_data FROM \"%w\" "
202 : "WHERE zoom_level = %d LIMIT 1",
203 56 : m_poTPD->m_osRasterTable.c_str(),
204 56 : m_poTPD->m_nZoomLevel);
205 : }
206 : else
207 : {
208 : // Try a tile in the middle of the raster
209 6 : pszSQL = sqlite3_mprintf(
210 : "SELECT tile_data FROM \"%w\" "
211 : "WHERE zoom_level = %d AND tile_column = %d AND tile_row = "
212 : "%d",
213 6 : m_poTPD->m_osRasterTable.c_str(), m_poTPD->m_nZoomLevel,
214 6 : m_poTPD->m_nShiftXTiles + nRasterXSize / 2 / nBlockXSize,
215 6 : m_poTPD->GetRowFromIntoTopConvention(
216 6 : m_poTPD->m_nShiftYTiles +
217 6 : nRasterYSize / 2 / nBlockYSize));
218 : }
219 62 : sqlite3_stmt *hStmt = nullptr;
220 62 : int rc = SQLPrepareWithError(m_poTPD->IGetDB(), pszSQL, -1, &hStmt,
221 : nullptr);
222 62 : if (rc == SQLITE_OK)
223 : {
224 62 : rc = sqlite3_step(hStmt);
225 79 : if (rc == SQLITE_ROW &&
226 17 : sqlite3_column_type(hStmt, 0) == SQLITE_BLOB)
227 : {
228 17 : const int nBytes = sqlite3_column_bytes(hStmt, 0);
229 : GByte *pabyRawData = reinterpret_cast<GByte *>(
230 17 : const_cast<void *>(sqlite3_column_blob(hStmt, 0)));
231 : const CPLString osMemFileName(
232 34 : VSIMemGenerateHiddenFilename("gpkg_read_tile"));
233 17 : VSILFILE *fp = VSIFileFromMemBuffer(
234 : osMemFileName.c_str(), pabyRawData, nBytes, FALSE);
235 17 : VSIFCloseL(fp);
236 :
237 : // Only PNG can have color table.
238 17 : const char *const apszDrivers[] = {"PNG", nullptr};
239 : auto poDSTile = std::unique_ptr<GDALDataset>(
240 : GDALDataset::Open(osMemFileName.c_str(),
241 : GDAL_OF_RASTER | GDAL_OF_INTERNAL,
242 34 : apszDrivers, nullptr, nullptr));
243 17 : if (poDSTile != nullptr)
244 : {
245 17 : if (poDSTile->GetRasterCount() == 1)
246 : {
247 10 : m_poTPD->m_poCT =
248 5 : poDSTile->GetRasterBand(1)->GetColorTable();
249 5 : if (m_poTPD->m_poCT != nullptr)
250 1 : m_poTPD->m_poCT = m_poTPD->m_poCT->Clone();
251 : }
252 : else
253 : {
254 12 : bRetry = true;
255 : }
256 : }
257 : else
258 0 : bRetry = true;
259 :
260 17 : VSIUnlink(osMemFileName);
261 : }
262 : }
263 62 : sqlite3_free(pszSQL);
264 62 : sqlite3_finalize(hStmt);
265 62 : if (!bRetry)
266 50 : break;
267 : }
268 : }
269 :
270 210 : return m_poTPD->m_poCT;
271 : }
272 :
273 : /************************************************************************/
274 : /* SetColorTable() */
275 : /************************************************************************/
276 :
277 15 : CPLErr GDALGPKGMBTilesLikeRasterBand::SetColorTable(GDALColorTable *poCT)
278 : {
279 15 : if (m_poTPD->m_eDT != GDT_Byte)
280 0 : return CE_Failure;
281 15 : if (poDS->GetRasterCount() != 1)
282 : {
283 1 : CPLError(CE_Failure, CPLE_NotSupported,
284 : "SetColorTable() only supported for a single band dataset");
285 1 : return CE_Failure;
286 : }
287 14 : if (!m_poTPD->m_bNew || m_poTPD->m_bTriedEstablishingCT)
288 : {
289 2 : CPLError(CE_Failure, CPLE_NotSupported,
290 : "SetColorTable() only supported on a newly created dataset");
291 2 : return CE_Failure;
292 : }
293 :
294 12 : AssignColorTable(poCT);
295 12 : return CE_None;
296 : }
297 :
298 : /************************************************************************/
299 : /* AssignColorTable() */
300 : /************************************************************************/
301 :
302 15 : void GDALGPKGMBTilesLikeRasterBand::AssignColorTable(const GDALColorTable *poCT)
303 : {
304 15 : m_poTPD->m_bTriedEstablishingCT = true;
305 15 : delete m_poTPD->m_poCT;
306 15 : if (poCT != nullptr)
307 14 : m_poTPD->m_poCT = poCT->Clone();
308 : else
309 1 : m_poTPD->m_poCT = nullptr;
310 15 : }
311 :
312 : /************************************************************************/
313 : /* GetColorInterpretation() */
314 : /************************************************************************/
315 :
316 209 : GDALColorInterp GDALGPKGMBTilesLikeRasterBand::GetColorInterpretation()
317 : {
318 209 : if (m_poTPD->m_eDT != GDT_Byte)
319 22 : return GCI_Undefined;
320 187 : if (poDS->GetRasterCount() == 1)
321 49 : return GetColorTable() ? GCI_PaletteIndex : GCI_GrayIndex;
322 138 : else if (poDS->GetRasterCount() == 2)
323 5 : return (nBand == 1) ? GCI_GrayIndex : GCI_AlphaBand;
324 : else
325 133 : return static_cast<GDALColorInterp>(GCI_RedBand + (nBand - 1));
326 : }
327 :
328 : /************************************************************************/
329 : /* SetColorInterpretation() */
330 : /************************************************************************/
331 :
332 : CPLErr
333 31 : GDALGPKGMBTilesLikeRasterBand::SetColorInterpretation(GDALColorInterp eInterp)
334 : {
335 31 : if (eInterp == GCI_Undefined)
336 1 : return CE_None;
337 32 : if (poDS->GetRasterCount() == 1 &&
338 2 : (eInterp == GCI_GrayIndex || eInterp == GCI_PaletteIndex))
339 23 : return CE_None;
340 11 : if (poDS->GetRasterCount() == 2 &&
341 4 : ((nBand == 1 && eInterp == GCI_GrayIndex) ||
342 3 : (nBand == 2 && eInterp == GCI_AlphaBand)))
343 2 : return CE_None;
344 5 : if (poDS->GetRasterCount() >= 3 && eInterp == GCI_RedBand + nBand - 1)
345 1 : return CE_None;
346 4 : CPLError(CE_Warning, CPLE_NotSupported,
347 : "%s color interpretation not supported. Will be ignored",
348 : GDALGetColorInterpretationName(eInterp));
349 4 : return CE_Warning;
350 : }
351 :
352 : /************************************************************************/
353 : /* GPKGFindBestEntry() */
354 : /************************************************************************/
355 :
356 1 : static int GPKGFindBestEntry(GDALColorTable *poCT, GByte c1, GByte c2, GByte c3,
357 : GByte c4, int nTileBandCount)
358 : {
359 1 : const int nEntries = std::min(256, poCT->GetColorEntryCount());
360 1 : int iBestIdx = 0;
361 1 : int nBestDistance = 4 * 256 * 256;
362 257 : for (int i = 0; i < nEntries; i++)
363 : {
364 256 : const GDALColorEntry *psEntry = poCT->GetColorEntry(i);
365 256 : int nDistance = (psEntry->c1 - c1) * (psEntry->c1 - c1) +
366 256 : (psEntry->c2 - c2) * (psEntry->c2 - c2) +
367 256 : (psEntry->c3 - c3) * (psEntry->c3 - c3);
368 256 : if (nTileBandCount == 4)
369 256 : nDistance += (psEntry->c4 - c4) * (psEntry->c4 - c4);
370 256 : if (nDistance < nBestDistance)
371 : {
372 10 : iBestIdx = i;
373 10 : nBestDistance = nDistance;
374 : }
375 : }
376 1 : return iBestIdx;
377 : }
378 :
379 : /************************************************************************/
380 : /* FillBuffer() */
381 : /************************************************************************/
382 :
383 14254 : void GDALGPKGMBTilesLikePseudoDataset::FillBuffer(GByte *pabyData,
384 : size_t nPixels)
385 : {
386 14254 : int bHasNoData = FALSE;
387 14254 : const double dfNoDataValue = IGetRasterBand(1)->GetNoDataValue(&bHasNoData);
388 14254 : if (!bHasNoData || dfNoDataValue == 0.0)
389 : {
390 : // cppcheck-suppress nullPointer
391 13952 : memset(pabyData, 0, nPixels * m_nDTSize);
392 : }
393 : else
394 : {
395 302 : GDALCopyWords64(&dfNoDataValue, GDT_Float64, 0, pabyData, m_eDT,
396 : m_nDTSize, nPixels);
397 : }
398 14254 : }
399 :
400 : /************************************************************************/
401 : /* FillEmptyTile() */
402 : /************************************************************************/
403 :
404 2955 : void GDALGPKGMBTilesLikePseudoDataset::FillEmptyTile(GByte *pabyData)
405 : {
406 : int nBlockXSize, nBlockYSize;
407 2955 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
408 2955 : const int nBands = IGetRasterCount();
409 2955 : const size_t nPixels =
410 2955 : static_cast<size_t>(nBands) * nBlockXSize * nBlockYSize;
411 2955 : FillBuffer(pabyData, nPixels);
412 2955 : }
413 :
414 : /************************************************************************/
415 : /* FillEmptyTileSingleBand() */
416 : /************************************************************************/
417 :
418 1481 : void GDALGPKGMBTilesLikePseudoDataset::FillEmptyTileSingleBand(GByte *pabyData)
419 : {
420 : int nBlockXSize, nBlockYSize;
421 1481 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
422 1481 : const size_t nPixels = static_cast<size_t>(nBlockXSize) * nBlockYSize;
423 1481 : FillBuffer(pabyData, nPixels);
424 1481 : }
425 :
426 : /************************************************************************/
427 : /* ReadTile() */
428 : /************************************************************************/
429 :
430 2582 : CPLErr GDALGPKGMBTilesLikePseudoDataset::ReadTile(
431 : const CPLString &osMemFileName, GByte *pabyTileData, double dfTileOffset,
432 : double dfTileScale, bool *pbIsLossyFormat)
433 : {
434 2582 : const char *const apszDriversByte[] = {"JPEG", "PNG", "WEBP", nullptr};
435 2582 : const char *const apszDriversInt[] = {"PNG", nullptr};
436 2582 : const char *const apszDriversFloat[] = {"GTiff", nullptr};
437 : int nBlockXSize, nBlockYSize;
438 2582 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
439 2582 : const int nBands = IGetRasterCount();
440 : auto poDSTile = std::unique_ptr<GDALDataset>(GDALDataset::Open(
441 : osMemFileName.c_str(), GDAL_OF_RASTER | GDAL_OF_INTERNAL,
442 2582 : (m_eDT == GDT_Byte) ? apszDriversByte
443 43 : : (m_eTF == GPKG_TF_TIFF_32BIT_FLOAT) ? apszDriversFloat
444 : : apszDriversInt,
445 5164 : nullptr, nullptr));
446 2582 : if (poDSTile == nullptr)
447 : {
448 2 : CPLError(CE_Failure, CPLE_AppDefined, "Cannot parse tile data");
449 2 : FillEmptyTile(pabyTileData);
450 2 : return CE_Failure;
451 : }
452 :
453 2580 : const int nTileBandCount = poDSTile->GetRasterCount();
454 :
455 5160 : if (!(poDSTile->GetRasterXSize() == nBlockXSize &&
456 2580 : poDSTile->GetRasterYSize() == nBlockYSize &&
457 7740 : (nTileBandCount >= 1 && nTileBandCount <= 4)) ||
458 2580 : (m_eDT != GDT_Byte && nTileBandCount != 1))
459 : {
460 0 : CPLError(CE_Failure, CPLE_AppDefined,
461 : "Inconsistent tiles characteristics");
462 0 : FillEmptyTile(pabyTileData);
463 0 : return CE_Failure;
464 : }
465 :
466 2580 : GDALDataType eRequestDT = GDT_Byte;
467 2580 : if (m_eTF == GPKG_TF_PNG_16BIT)
468 : {
469 40 : CPLAssert(m_eDT == GDT_Int16 || m_eDT == GDT_UInt16 ||
470 : m_eDT == GDT_Float32);
471 40 : eRequestDT = GDT_UInt16;
472 : }
473 2540 : else if (m_eTF == GPKG_TF_TIFF_32BIT_FLOAT)
474 : {
475 3 : CPLAssert(m_eDT == GDT_Float32);
476 3 : eRequestDT = GDT_Float32;
477 : }
478 :
479 2580 : if (poDSTile->RasterIO(GF_Read, 0, 0, nBlockXSize, nBlockYSize,
480 : pabyTileData, nBlockXSize, nBlockYSize, eRequestDT,
481 : poDSTile->GetRasterCount(), nullptr, 0, 0, 0,
482 2580 : nullptr) != CE_None)
483 : {
484 0 : FillEmptyTile(pabyTileData);
485 0 : return CE_Failure;
486 : }
487 :
488 2580 : if (m_eDT != GDT_Byte)
489 : {
490 43 : int bHasNoData = FALSE;
491 : const double dfNoDataValue =
492 43 : IGetRasterBand(1)->GetNoDataValue(&bHasNoData);
493 43 : if (m_eDT == GDT_Int16)
494 : {
495 29 : CPLAssert(eRequestDT == GDT_UInt16);
496 29 : for (size_t i = 0;
497 1900570 : i < static_cast<size_t>(nBlockXSize) * nBlockYSize; i++)
498 : {
499 : GUInt16 usVal;
500 1900540 : memcpy(&usVal, pabyTileData + i * sizeof(GUInt16),
501 : sizeof(usVal));
502 : double dfVal =
503 1900540 : floor((usVal * dfTileScale + dfTileOffset) * m_dfScale +
504 1900540 : m_dfOffset + 0.5);
505 1900540 : if (bHasNoData && usVal == m_usGPKGNull)
506 195804 : dfVal = dfNoDataValue;
507 1900540 : if (dfVal > 32767)
508 0 : dfVal = 32767;
509 1900540 : else if (dfVal < -32768)
510 0 : dfVal = -32768;
511 1900540 : GInt16 sVal = static_cast<GInt16>(dfVal);
512 1900540 : memcpy(pabyTileData + i * sizeof(GUInt16), &sVal, sizeof(sVal));
513 : }
514 : }
515 14 : else if (m_eDT == GDT_UInt16 &&
516 6 : (m_dfOffset != 0.0 || m_dfScale != 1.0 ||
517 6 : dfTileOffset != 0.0 || dfTileScale != 1.0))
518 : {
519 0 : CPLAssert(eRequestDT == GDT_UInt16);
520 0 : GUInt16 *psVal = reinterpret_cast<GUInt16 *>(pabyTileData);
521 0 : for (size_t i = 0;
522 0 : i < static_cast<size_t>(nBlockXSize) * nBlockYSize; i++)
523 : {
524 0 : const GUInt16 nVal = psVal[i];
525 : double dfVal =
526 0 : floor((nVal * dfTileScale + dfTileOffset) * m_dfScale +
527 0 : m_dfOffset + 0.5);
528 0 : if (bHasNoData && nVal == m_usGPKGNull)
529 0 : dfVal = dfNoDataValue;
530 0 : if (dfVal > 65535)
531 0 : dfVal = 65535;
532 0 : else if (dfVal < 0)
533 0 : dfVal = 0;
534 0 : psVal[i] = static_cast<GUInt16>(dfVal);
535 0 : }
536 : }
537 14 : else if (m_eDT == GDT_Float32 && eRequestDT == GDT_UInt16)
538 : {
539 : // Due to non identical data type size, we need to start from the
540 : // end of the buffer.
541 327685 : for (GPtrDiff_t i =
542 5 : static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize - 1;
543 327685 : i >= 0; i--)
544 : {
545 : // Use memcpy() and not reinterpret_cast<GUInt16*> and
546 : // reinterpret_cast<float*>, otherwise compilers such as ICC
547 : // may (ab)use rules about aliasing to generate wrong code!
548 : GUInt16 usVal;
549 327680 : memcpy(&usVal, pabyTileData + i * sizeof(GUInt16),
550 : sizeof(usVal));
551 327680 : double dfVal =
552 327680 : (usVal * dfTileScale + dfTileOffset) * m_dfScale +
553 327680 : m_dfOffset;
554 327680 : if (m_dfPrecision == 1.0)
555 327680 : dfVal = floor(dfVal + 0.5);
556 327680 : if (bHasNoData && usVal == m_usGPKGNull)
557 196206 : dfVal = dfNoDataValue;
558 327680 : const float fVal = static_cast<float>(dfVal);
559 327680 : memcpy(pabyTileData + i * sizeof(float), &fVal, sizeof(fVal));
560 : }
561 : }
562 :
563 43 : return CE_None;
564 : }
565 :
566 2537 : GDALColorTable *poCT = nullptr;
567 2537 : if (nBands == 1 || nTileBandCount == 1)
568 : {
569 153 : poCT = poDSTile->GetRasterBand(1)->GetColorTable();
570 153 : IGetRasterBand(1)->GetColorTable();
571 : }
572 :
573 2537 : if (pbIsLossyFormat)
574 9 : *pbIsLossyFormat =
575 9 : !EQUAL(poDSTile->GetDriver()->GetDescription(), "PNG") ||
576 0 : (poCT != nullptr && poCT->GetColorEntryCount() == 256) /* PNG8 */;
577 :
578 : /* Map RGB(A) tile to single-band color indexed */
579 2537 : const GPtrDiff_t nBlockPixels =
580 2537 : static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize;
581 2537 : if (nBands == 1 && m_poCT != nullptr && nTileBandCount != 1)
582 : {
583 1 : std::map<GUInt32, int> oMapEntryToIndex;
584 1 : const int nEntries = std::min(256, m_poCT->GetColorEntryCount());
585 257 : for (int i = 0; i < nEntries; i++)
586 : {
587 256 : const GDALColorEntry *psEntry = m_poCT->GetColorEntry(i);
588 256 : GByte c1 = static_cast<GByte>(psEntry->c1);
589 256 : GByte c2 = static_cast<GByte>(psEntry->c2);
590 256 : GByte c3 = static_cast<GByte>(psEntry->c3);
591 256 : GUInt32 nVal = c1 + (c2 << 8) + (c3 << 16);
592 256 : if (nTileBandCount == 4)
593 256 : nVal += (static_cast<GUInt32>(psEntry->c4) << 24);
594 256 : oMapEntryToIndex[nVal] = i;
595 : }
596 : int iBestEntryFor0 =
597 1 : GPKGFindBestEntry(m_poCT, 0, 0, 0, 0, nTileBandCount);
598 10001 : for (GPtrDiff_t i = 0; i < nBlockPixels; i++)
599 : {
600 10000 : const GByte c1 = pabyTileData[i];
601 10000 : const GByte c2 = pabyTileData[i + nBlockPixels];
602 10000 : const GByte c3 = pabyTileData[i + 2 * nBlockPixels];
603 10000 : const GByte c4 = pabyTileData[i + 3 * nBlockPixels];
604 10000 : GUInt32 nVal = c1 + (c2 << 8) + (c3 << 16);
605 10000 : if (nTileBandCount == 4)
606 10000 : nVal += (c4 << 24);
607 10000 : if (nVal == 0)
608 : // In most cases we will reach that point at partial tiles.
609 5000 : pabyTileData[i] = static_cast<GByte>(iBestEntryFor0);
610 : else
611 : {
612 : std::map<GUInt32, int>::iterator oMapEntryToIndexIter =
613 5000 : oMapEntryToIndex.find(nVal);
614 5000 : if (oMapEntryToIndexIter == oMapEntryToIndex.end())
615 : /* Could happen with JPEG tiles */
616 0 : pabyTileData[i] = static_cast<GByte>(GPKGFindBestEntry(
617 : m_poCT, c1, c2, c3, c4, nTileBandCount));
618 : else
619 5000 : pabyTileData[i] =
620 5000 : static_cast<GByte>(oMapEntryToIndexIter->second);
621 : }
622 : }
623 1 : return CE_None;
624 : }
625 :
626 31 : if (nBands == 1 && nTileBandCount == 1 && poCT != nullptr &&
627 2567 : m_poCT != nullptr && !poCT->IsSame(m_poCT))
628 : {
629 0 : CPLError(CE_Warning, CPLE_NotSupported,
630 : "Different color tables. Unhandled for now");
631 : }
632 2536 : else if ((nBands == 1 && nTileBandCount >= 3) ||
633 31 : (nBands == 1 && nTileBandCount == 1 && m_poCT != nullptr &&
634 2536 : poCT == nullptr) ||
635 2536 : ((nBands == 1 || nBands == 2) && nTileBandCount == 1 &&
636 35 : m_poCT == nullptr && poCT != nullptr))
637 : {
638 0 : CPLError(CE_Failure, CPLE_AppDefined,
639 : "Inconsistent dataset and tiles band characteristics");
640 : }
641 :
642 2536 : if (nBands == 2)
643 : {
644 : // assuming that the RGB is Grey,Grey,Grey
645 307 : if (nTileBandCount == 1 || nTileBandCount == 3)
646 : {
647 : /* Create fully opaque alpha */
648 97 : memset(pabyTileData + 1 * nBlockPixels, 255, nBlockPixels);
649 : }
650 210 : else if (nTileBandCount == 4)
651 : {
652 : /* Transfer alpha band */
653 67 : memcpy(pabyTileData + 1 * nBlockPixels,
654 67 : pabyTileData + 3 * nBlockPixels, nBlockPixels);
655 : }
656 : }
657 2229 : else if (nTileBandCount == 2)
658 : {
659 : /* Do Grey+Alpha -> RGBA */
660 362 : memcpy(pabyTileData + 3 * nBlockPixels, pabyTileData + 1 * nBlockPixels,
661 : nBlockPixels);
662 362 : memcpy(pabyTileData + 1 * nBlockPixels, pabyTileData, nBlockPixels);
663 362 : memcpy(pabyTileData + 2 * nBlockPixels, pabyTileData, nBlockPixels);
664 : }
665 1867 : else if (nTileBandCount == 1 && !(nBands == 1 && m_poCT != nullptr))
666 : {
667 : /* Expand color indexed to RGB(A) */
668 122 : if (poCT != nullptr)
669 : {
670 : GByte abyCT[4 * 256];
671 15 : const int nEntries = std::min(256, poCT->GetColorEntryCount());
672 3600 : for (int i = 0; i < nEntries; i++)
673 : {
674 3585 : const GDALColorEntry *psEntry = poCT->GetColorEntry(i);
675 3585 : abyCT[4 * i] = static_cast<GByte>(psEntry->c1);
676 3585 : abyCT[4 * i + 1] = static_cast<GByte>(psEntry->c2);
677 3585 : abyCT[4 * i + 2] = static_cast<GByte>(psEntry->c3);
678 3585 : abyCT[4 * i + 3] = static_cast<GByte>(psEntry->c4);
679 : }
680 270 : for (int i = nEntries; i < 256; i++)
681 : {
682 255 : abyCT[4 * i] = 0;
683 255 : abyCT[4 * i + 1] = 0;
684 255 : abyCT[4 * i + 2] = 0;
685 255 : abyCT[4 * i + 3] = 0;
686 : }
687 967175 : for (GPtrDiff_t i = 0; i < nBlockPixels; i++)
688 : {
689 967160 : const GByte byVal = pabyTileData[i];
690 967160 : pabyTileData[i] = abyCT[4 * byVal];
691 967160 : pabyTileData[i + 1 * nBlockPixels] = abyCT[4 * byVal + 1];
692 967160 : pabyTileData[i + 2 * nBlockPixels] = abyCT[4 * byVal + 2];
693 967160 : pabyTileData[i + 3 * nBlockPixels] = abyCT[4 * byVal + 3];
694 : }
695 : }
696 : else
697 : {
698 107 : memcpy(pabyTileData + 1 * nBlockPixels, pabyTileData, nBlockPixels);
699 107 : memcpy(pabyTileData + 2 * nBlockPixels, pabyTileData, nBlockPixels);
700 107 : if (nBands == 4)
701 : {
702 96 : memset(pabyTileData + 3 * nBlockPixels, 255, nBlockPixels);
703 : }
704 122 : }
705 : }
706 1745 : else if (nTileBandCount == 3 && nBands == 4)
707 : {
708 : /* Create fully opaque alpha */
709 224 : memset(pabyTileData + 3 * nBlockPixels, 255, nBlockPixels);
710 : }
711 :
712 2536 : return CE_None;
713 : }
714 :
715 : /************************************************************************/
716 : /* ReadTile() */
717 : /************************************************************************/
718 :
719 7542 : GByte *GDALGPKGMBTilesLikePseudoDataset::ReadTile(int nRow, int nCol)
720 : {
721 : int nBlockXSize, nBlockYSize;
722 7542 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
723 7542 : const int nBands = IGetRasterCount();
724 7542 : const size_t nBandBlockSize =
725 7542 : static_cast<size_t>(nBlockXSize) * nBlockYSize * m_nDTSize;
726 7542 : const int nTileBands = m_eDT == GDT_Byte ? 4 : 1;
727 7542 : if (m_nShiftXPixelsMod || m_nShiftYPixelsMod)
728 : {
729 4803 : GByte *pabyData = nullptr;
730 4803 : int i = 0;
731 11981 : for (; i < 4; i++)
732 : {
733 11981 : if (m_asCachedTilesDesc[i].nRow == nRow &&
734 7181 : m_asCachedTilesDesc[i].nCol == nCol)
735 : {
736 4803 : if (m_asCachedTilesDesc[i].nIdxWithinTileData >= 0)
737 : {
738 2040 : return m_pabyCachedTiles +
739 2040 : nBandBlockSize *
740 2040 : m_asCachedTilesDesc[i].nIdxWithinTileData *
741 2040 : nTileBands;
742 : }
743 : else
744 : {
745 2763 : if (i == 0)
746 194 : m_asCachedTilesDesc[i].nIdxWithinTileData =
747 194 : (m_asCachedTilesDesc[1].nIdxWithinTileData == 0)
748 194 : ? 1
749 : : 0;
750 2569 : else if (i == 1)
751 1189 : m_asCachedTilesDesc[i].nIdxWithinTileData =
752 1189 : (m_asCachedTilesDesc[0].nIdxWithinTileData == 0)
753 1189 : ? 1
754 : : 0;
755 1380 : else if (i == 2)
756 191 : m_asCachedTilesDesc[i].nIdxWithinTileData =
757 191 : (m_asCachedTilesDesc[3].nIdxWithinTileData == 2)
758 191 : ? 3
759 : : 2;
760 : else
761 1189 : m_asCachedTilesDesc[i].nIdxWithinTileData =
762 1189 : (m_asCachedTilesDesc[2].nIdxWithinTileData == 2)
763 1189 : ? 3
764 : : 2;
765 2763 : pabyData = m_pabyCachedTiles +
766 2763 : nBandBlockSize *
767 2763 : m_asCachedTilesDesc[i].nIdxWithinTileData *
768 2763 : nTileBands;
769 2763 : break;
770 : }
771 : }
772 : }
773 2763 : CPLAssert(i < 4);
774 2763 : return ReadTile(nRow, nCol, pabyData);
775 : }
776 : else
777 : {
778 2739 : GByte *pabyDest = m_pabyCachedTiles + 2 * nTileBands * nBandBlockSize;
779 2739 : bool bAllNonDirty = true;
780 11227 : for (int i = 0; i < nBands; i++)
781 : {
782 8488 : if (m_asCachedTilesDesc[0].abBandDirty[i])
783 2 : bAllNonDirty = false;
784 : }
785 2739 : if (bAllNonDirty)
786 : {
787 2737 : return ReadTile(nRow, nCol, pabyDest);
788 : }
789 :
790 : /* If some bands of the blocks are dirty/written we need to fetch */
791 : /* the tile in a temporary buffer in order not to override dirty bands*/
792 2 : GByte *pabyTemp = m_pabyCachedTiles + 3 * nTileBands * nBandBlockSize;
793 2 : if (ReadTile(nRow, nCol, pabyTemp) != nullptr)
794 : {
795 8 : for (int i = 0; i < nBands; i++)
796 : {
797 6 : if (!m_asCachedTilesDesc[0].abBandDirty[i])
798 : {
799 4 : memcpy(pabyDest + i * nBandBlockSize,
800 4 : pabyTemp + i * nBandBlockSize, nBandBlockSize);
801 : }
802 : }
803 : }
804 2 : return pabyDest;
805 : }
806 : }
807 :
808 : /************************************************************************/
809 : /* GetTileOffsetAndScale() */
810 : /************************************************************************/
811 :
812 2582 : void GDALGPKGMBTilesLikePseudoDataset::GetTileOffsetAndScale(
813 : GIntBig nTileId, double &dfTileOffset, double &dfTileScale)
814 : {
815 2582 : dfTileOffset = 0.0;
816 2582 : dfTileScale = 1.0;
817 :
818 2582 : if (m_eTF == GPKG_TF_PNG_16BIT)
819 : {
820 40 : char *pszSQL = sqlite3_mprintf(
821 : "SELECT offset, scale FROM gpkg_2d_gridded_tile_ancillary WHERE "
822 : "tpudt_name = '%q' AND tpudt_id = ?",
823 : m_osRasterTable.c_str());
824 40 : sqlite3_stmt *hStmt = nullptr;
825 40 : int rc = SQLPrepareWithError(IGetDB(), pszSQL, -1, &hStmt, nullptr);
826 40 : if (rc == SQLITE_OK)
827 : {
828 40 : sqlite3_bind_int64(hStmt, 1, nTileId);
829 40 : rc = sqlite3_step(hStmt);
830 40 : if (rc == SQLITE_ROW)
831 : {
832 40 : if (sqlite3_column_type(hStmt, 0) == SQLITE_FLOAT)
833 40 : dfTileOffset = sqlite3_column_double(hStmt, 0);
834 40 : if (sqlite3_column_type(hStmt, 1) == SQLITE_FLOAT)
835 40 : dfTileScale = sqlite3_column_double(hStmt, 1);
836 : }
837 40 : sqlite3_finalize(hStmt);
838 : }
839 40 : sqlite3_free(pszSQL);
840 : }
841 2582 : }
842 :
843 : /************************************************************************/
844 : /* ReadTile() */
845 : /************************************************************************/
846 :
847 5520 : GByte *GDALGPKGMBTilesLikePseudoDataset::ReadTile(int nRow, int nCol,
848 : GByte *pabyData,
849 : bool *pbIsLossyFormat)
850 : {
851 5520 : int nBlockXSize = 0;
852 5520 : int nBlockYSize = 0;
853 5520 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
854 5520 : const int nBands = IGetRasterCount();
855 :
856 5520 : if (pbIsLossyFormat)
857 18 : *pbIsLossyFormat = false;
858 :
859 5520 : const size_t nBandBlockSize =
860 5520 : static_cast<size_t>(nBlockXSize) * nBlockYSize * m_nDTSize;
861 5520 : if (nRow < 0 || nCol < 0 || nRow >= m_nTileMatrixHeight ||
862 5463 : nCol >= m_nTileMatrixWidth)
863 : {
864 73 : FillEmptyTile(pabyData);
865 73 : return pabyData;
866 : }
867 :
868 : #ifdef DEBUG_VERBOSE
869 : CPLDebug("GPKG", "ReadTile(row=%d, col=%d)", nRow, nCol);
870 : #endif
871 :
872 16341 : char *pszSQL = sqlite3_mprintf(
873 : "SELECT tile_data%s FROM \"%w\" "
874 : "WHERE zoom_level = %d AND tile_row = %d AND tile_column = %d%s",
875 5447 : m_eDT != GDT_Byte ? ", id" : "", // MBTiles do not have an id
876 : m_osRasterTable.c_str(), m_nZoomLevel,
877 5447 : GetRowFromIntoTopConvention(nRow), nCol,
878 5447 : !m_osWHERE.empty() ? CPLSPrintf(" AND (%s)", m_osWHERE.c_str()) : "");
879 :
880 : #ifdef DEBUG_VERBOSE
881 : CPLDebug("GPKG", "%s", pszSQL);
882 : #endif
883 :
884 5447 : sqlite3_stmt *hStmt = nullptr;
885 5447 : int rc = SQLPrepareWithError(IGetDB(), pszSQL, -1, &hStmt, nullptr);
886 5447 : sqlite3_free(pszSQL);
887 5447 : if (rc != SQLITE_OK)
888 : {
889 0 : return nullptr;
890 : }
891 5447 : rc = sqlite3_step(hStmt);
892 :
893 5447 : if (rc == SQLITE_ROW && sqlite3_column_type(hStmt, 0) == SQLITE_BLOB)
894 : {
895 2565 : const int nBytes = sqlite3_column_bytes(hStmt, 0);
896 : GIntBig nTileId =
897 2565 : (m_eDT == GDT_Byte) ? 0 : sqlite3_column_int64(hStmt, 1);
898 : GByte *pabyRawData = static_cast<GByte *>(
899 2565 : const_cast<void *>(sqlite3_column_blob(hStmt, 0)));
900 : const CPLString osMemFileName(
901 5130 : VSIMemGenerateHiddenFilename("gpkg_read_tile"));
902 2565 : VSILFILE *fp = VSIFileFromMemBuffer(osMemFileName.c_str(), pabyRawData,
903 : nBytes, FALSE);
904 2565 : VSIFCloseL(fp);
905 :
906 2565 : double dfTileOffset = 0.0;
907 2565 : double dfTileScale = 1.0;
908 2565 : GetTileOffsetAndScale(nTileId, dfTileOffset, dfTileScale);
909 2565 : ReadTile(osMemFileName, pabyData, dfTileOffset, dfTileScale,
910 : pbIsLossyFormat);
911 2565 : VSIUnlink(osMemFileName);
912 2565 : sqlite3_finalize(hStmt);
913 : }
914 2882 : else if (rc == SQLITE_BUSY)
915 : {
916 0 : FillEmptyTile(pabyData);
917 0 : CPLError(CE_Failure, CPLE_AppDefined,
918 : "sqlite3_step(%s) failed (SQLITE_BUSY): %s",
919 0 : sqlite3_sql(hStmt), sqlite3_errmsg(IGetDB()));
920 0 : sqlite3_finalize(hStmt);
921 0 : return pabyData;
922 : }
923 : else
924 : {
925 2882 : sqlite3_finalize(hStmt);
926 2882 : hStmt = nullptr;
927 :
928 2882 : if (m_hTempDB && (m_nShiftXPixelsMod || m_nShiftYPixelsMod))
929 : {
930 602 : const char *pszSQLNew = CPLSPrintf(
931 : "SELECT partial_flag, tile_data_band_1, tile_data_band_2, "
932 : "tile_data_band_3, tile_data_band_4 FROM partial_tiles WHERE "
933 : "zoom_level = %d AND tile_row = %d AND tile_column = %d",
934 : m_nZoomLevel, nRow, nCol);
935 :
936 : #ifdef DEBUG_VERBOSE
937 : CPLDebug("GPKG", "%s", pszSQLNew);
938 : #endif
939 :
940 602 : rc = SQLPrepareWithError(m_hTempDB, pszSQLNew, -1, &hStmt, nullptr);
941 602 : if (rc != SQLITE_OK)
942 : {
943 0 : FillEmptyTile(pabyData);
944 0 : return pabyData;
945 : }
946 :
947 602 : rc = sqlite3_step(hStmt);
948 602 : if (rc == SQLITE_ROW)
949 : {
950 2 : const int nPartialFlag = sqlite3_column_int(hStmt, 0);
951 10 : for (int iBand = 1; iBand <= nBands; iBand++)
952 : {
953 8 : GByte *pabyDestBand =
954 8 : pabyData + (iBand - 1) * nBandBlockSize;
955 8 : if (nPartialFlag & (((1 << 4) - 1) << (4 * (iBand - 1))))
956 : {
957 2 : CPLAssert(sqlite3_column_bytes(hStmt, iBand) ==
958 : static_cast<int>(nBandBlockSize));
959 2 : memcpy(pabyDestBand, sqlite3_column_blob(hStmt, iBand),
960 : nBandBlockSize);
961 : }
962 : else
963 : {
964 6 : FillEmptyTileSingleBand(pabyDestBand);
965 : }
966 : }
967 : }
968 : else
969 : {
970 600 : FillEmptyTile(pabyData);
971 : }
972 602 : sqlite3_finalize(hStmt);
973 : }
974 : else
975 : {
976 2280 : FillEmptyTile(pabyData);
977 : }
978 : }
979 :
980 5447 : return pabyData;
981 : }
982 :
983 : /************************************************************************/
984 : /* IReadBlock() */
985 : /************************************************************************/
986 :
987 3950 : CPLErr GDALGPKGMBTilesLikeRasterBand::IReadBlock(int nBlockXOff, int nBlockYOff,
988 : void *pData)
989 : {
990 : #ifdef DEBUG_VERBOSE
991 : CPLDebug("GPKG",
992 : "IReadBlock(nBand=%d,nBlockXOff=%d,nBlockYOff=%d,m_nZoomLevel=%d)",
993 : nBand, nBlockXOff, nBlockYOff, m_poTPD->m_nZoomLevel);
994 : #endif
995 :
996 3950 : if (m_poTPD->m_pabyCachedTiles == nullptr)
997 0 : return CE_Failure;
998 :
999 3950 : const int nRowMin = nBlockYOff + m_poTPD->m_nShiftYTiles;
1000 3950 : int nRowMax = nRowMin;
1001 3950 : if (m_poTPD->m_nShiftYPixelsMod)
1002 1211 : nRowMax++;
1003 :
1004 3950 : const int nColMin = nBlockXOff + m_poTPD->m_nShiftXTiles;
1005 3950 : int nColMax = nColMin;
1006 3950 : if (m_poTPD->m_nShiftXPixelsMod)
1007 1189 : nColMax++;
1008 :
1009 2761 : retry:
1010 : /* Optimize for left to right reading at constant row */
1011 3953 : if (m_poTPD->m_nShiftXPixelsMod || m_poTPD->m_nShiftYPixelsMod)
1012 : {
1013 1214 : if (nRowMin == m_poTPD->m_asCachedTilesDesc[0].nRow &&
1014 1036 : nColMin == m_poTPD->m_asCachedTilesDesc[0].nCol + 1 &&
1015 1030 : m_poTPD->m_asCachedTilesDesc[0].nIdxWithinTileData >= 0)
1016 : {
1017 1030 : CPLAssert(nRowMin == m_poTPD->m_asCachedTilesDesc[1].nRow);
1018 1030 : CPLAssert(nColMin == m_poTPD->m_asCachedTilesDesc[1].nCol);
1019 1030 : CPLAssert(m_poTPD->m_asCachedTilesDesc[0].nIdxWithinTileData == 0 ||
1020 : m_poTPD->m_asCachedTilesDesc[0].nIdxWithinTileData == 1);
1021 :
1022 : /* 0 1 --> 1 -1 */
1023 : /* 2 3 3 -1 */
1024 : /* or */
1025 : /* 1 0 --> 0 -1 */
1026 : /* 3 2 2 -1 */
1027 1030 : m_poTPD->m_asCachedTilesDesc[0].nIdxWithinTileData =
1028 1030 : m_poTPD->m_asCachedTilesDesc[1].nIdxWithinTileData;
1029 1030 : m_poTPD->m_asCachedTilesDesc[2].nIdxWithinTileData =
1030 1030 : m_poTPD->m_asCachedTilesDesc[3].nIdxWithinTileData;
1031 : }
1032 : else
1033 : {
1034 184 : m_poTPD->m_asCachedTilesDesc[0].nIdxWithinTileData = -1;
1035 184 : m_poTPD->m_asCachedTilesDesc[2].nIdxWithinTileData = -1;
1036 : }
1037 1214 : m_poTPD->m_asCachedTilesDesc[0].nRow = nRowMin;
1038 1214 : m_poTPD->m_asCachedTilesDesc[0].nCol = nColMin;
1039 1214 : m_poTPD->m_asCachedTilesDesc[1].nRow = nRowMin;
1040 1214 : m_poTPD->m_asCachedTilesDesc[1].nCol = nColMin + 1;
1041 1214 : m_poTPD->m_asCachedTilesDesc[2].nRow = nRowMin + 1;
1042 1214 : m_poTPD->m_asCachedTilesDesc[2].nCol = nColMin;
1043 1214 : m_poTPD->m_asCachedTilesDesc[3].nRow = nRowMin + 1;
1044 1214 : m_poTPD->m_asCachedTilesDesc[3].nCol = nColMin + 1;
1045 1214 : m_poTPD->m_asCachedTilesDesc[1].nIdxWithinTileData = -1;
1046 1214 : m_poTPD->m_asCachedTilesDesc[3].nIdxWithinTileData = -1;
1047 : }
1048 :
1049 9114 : for (int nRow = nRowMin; nRow <= nRowMax; nRow++)
1050 : {
1051 12703 : for (int nCol = nColMin; nCol <= nColMax; nCol++)
1052 : {
1053 7542 : if (m_poTPD->m_nShiftXPixelsMod == 0 &&
1054 2786 : m_poTPD->m_nShiftYPixelsMod == 0)
1055 : {
1056 2739 : if (!(nRow == m_poTPD->m_asCachedTilesDesc[0].nRow &&
1057 7 : nCol == m_poTPD->m_asCachedTilesDesc[0].nCol &&
1058 3 : m_poTPD->m_asCachedTilesDesc[0].nIdxWithinTileData == 0))
1059 : {
1060 2737 : if (m_poTPD->WriteTile() != CE_None)
1061 0 : return CE_Failure;
1062 : }
1063 : }
1064 :
1065 7542 : GByte *pabyTileData = m_poTPD->ReadTile(nRow, nCol);
1066 7542 : if (pabyTileData == nullptr)
1067 0 : return CE_Failure;
1068 :
1069 34600 : for (int iBand = 1; iBand <= poDS->GetRasterCount(); iBand++)
1070 : {
1071 27061 : GDALRasterBlock *poBlock = nullptr;
1072 27061 : GByte *pabyDest = nullptr;
1073 27061 : if (iBand == nBand)
1074 : {
1075 7542 : pabyDest = static_cast<GByte *>(pData);
1076 : }
1077 : else
1078 : {
1079 19519 : poBlock = poDS->GetRasterBand(iBand)->GetLockedBlockRef(
1080 19519 : nBlockXOff, nBlockYOff, TRUE);
1081 19519 : if (poBlock == nullptr)
1082 0 : continue;
1083 19519 : if (poBlock->GetDirty())
1084 : {
1085 2 : poBlock->DropLock();
1086 2 : continue;
1087 : }
1088 : /* if we are short of GDAL cache max and there are dirty
1089 : * blocks */
1090 : /* of our dataset, the above GetLockedBlockRef() might have
1091 : * reset */
1092 : /* (at least part of) the 4 tiles we want to cache and have
1093 : */
1094 : /* already read */
1095 : // FIXME this is way too fragile.
1096 19517 : if ((m_poTPD->m_nShiftXPixelsMod != 0 ||
1097 5857 : m_poTPD->m_nShiftYPixelsMod != 0) &&
1098 13770 : (m_poTPD->m_asCachedTilesDesc[0].nRow != nRowMin ||
1099 13767 : m_poTPD->m_asCachedTilesDesc[0].nCol != nColMin))
1100 : {
1101 3 : poBlock->DropLock();
1102 3 : goto retry;
1103 : }
1104 19514 : pabyDest = static_cast<GByte *>(poBlock->GetDataRef());
1105 : }
1106 :
1107 : // Composite tile data into block data
1108 27056 : if (m_poTPD->m_nShiftXPixelsMod == 0 &&
1109 8640 : m_poTPD->m_nShiftYPixelsMod == 0)
1110 : {
1111 8486 : const size_t nBandBlockSize =
1112 8486 : static_cast<size_t>(nBlockXSize) * nBlockYSize *
1113 8486 : m_nDTSize;
1114 8486 : memcpy(pabyDest,
1115 8486 : pabyTileData + (iBand - 1) * nBandBlockSize,
1116 8486 : nBandBlockSize);
1117 : #ifdef DEBUG_VERBOSE
1118 : if (eDataType == GDT_Byte &&
1119 : (nBlockXOff + 1) * nBlockXSize <= nRasterXSize &&
1120 : (nBlockYOff + 1) * nBlockYSize > nRasterYSize)
1121 : {
1122 : bool bFoundNonZero = false;
1123 : for (int y = nRasterYSize - nBlockYOff * nBlockYSize;
1124 : y < nBlockYSize; y++)
1125 : {
1126 : for (int x = 0; x < nBlockXSize; x++)
1127 : {
1128 : if (pabyDest[static_cast<GPtrDiff_t>(y) *
1129 : nBlockXSize +
1130 : x] != 0 &&
1131 : !bFoundNonZero)
1132 : {
1133 : CPLDebug("GPKG",
1134 : "IReadBlock(): Found non-zero "
1135 : "content in ghost part of "
1136 : "tile(nBand=%d,nBlockXOff=%d,"
1137 : "nBlockYOff=%d,m_nZoomLevel=%d)\n",
1138 : iBand, nBlockXOff, nBlockYOff,
1139 : m_poTPD->m_nZoomLevel);
1140 : bFoundNonZero = true;
1141 : }
1142 : }
1143 : }
1144 : }
1145 : #endif
1146 : }
1147 : else
1148 : {
1149 18570 : int nSrcXOffset = 0;
1150 18570 : int nSrcXSize = 0;
1151 18570 : int nDstXOffset = 0;
1152 18570 : if (nCol == nColMin)
1153 : {
1154 9362 : nSrcXOffset = m_poTPD->m_nShiftXPixelsMod;
1155 9362 : nSrcXSize = nBlockXSize - m_poTPD->m_nShiftXPixelsMod;
1156 9362 : nDstXOffset = 0;
1157 : }
1158 : else
1159 : {
1160 9208 : nSrcXOffset = 0;
1161 9208 : nSrcXSize = m_poTPD->m_nShiftXPixelsMod;
1162 9208 : nDstXOffset = nBlockXSize - m_poTPD->m_nShiftXPixelsMod;
1163 : }
1164 18570 : int nSrcYOffset = 0;
1165 18570 : int nSrcYSize = 0;
1166 18570 : int nDstYOffset = 0;
1167 18570 : if (nRow == nRowMin)
1168 : {
1169 9288 : nSrcYOffset = m_poTPD->m_nShiftYPixelsMod;
1170 9288 : nSrcYSize = nBlockYSize - m_poTPD->m_nShiftYPixelsMod;
1171 9288 : nDstYOffset = 0;
1172 : }
1173 : else
1174 : {
1175 9282 : nSrcYOffset = 0;
1176 9282 : nSrcYSize = m_poTPD->m_nShiftYPixelsMod;
1177 9282 : nDstYOffset = nBlockYSize - m_poTPD->m_nShiftYPixelsMod;
1178 : }
1179 :
1180 : #ifdef DEBUG_VERBOSE
1181 : CPLDebug("GPKG",
1182 : "Copy source tile x=%d,w=%d,y=%d,h=%d into "
1183 : "buffer at x=%d,y=%d",
1184 : nSrcXOffset, nSrcXSize, nSrcYOffset, nSrcYSize,
1185 : nDstXOffset, nDstYOffset);
1186 : #endif
1187 :
1188 2368610 : for (GPtrDiff_t y = 0; y < nSrcYSize; y++)
1189 : {
1190 2350040 : GByte *pSrc =
1191 : pabyTileData +
1192 2350040 : (static_cast<GPtrDiff_t>(iBand - 1) * nBlockXSize *
1193 2350040 : nBlockYSize +
1194 2350040 : (y + nSrcYOffset) * nBlockXSize + nSrcXOffset) *
1195 2350040 : m_nDTSize;
1196 2350040 : GByte *pDst =
1197 : pabyDest +
1198 2350040 : ((y + nDstYOffset) * nBlockXSize + nDstXOffset) *
1199 2350040 : m_nDTSize;
1200 2350040 : GDALCopyWords(pSrc, eDataType, m_nDTSize, pDst,
1201 : eDataType, m_nDTSize, nSrcXSize);
1202 : }
1203 : }
1204 :
1205 27056 : if (poBlock)
1206 19514 : poBlock->DropLock();
1207 : }
1208 : }
1209 : }
1210 :
1211 3950 : return CE_None;
1212 : }
1213 :
1214 : /************************************************************************/
1215 : /* IGetDataCoverageStatus() */
1216 : /************************************************************************/
1217 :
1218 8 : int GDALGPKGMBTilesLikeRasterBand::IGetDataCoverageStatus(int nXOff, int nYOff,
1219 : int nXSize,
1220 : int nYSize,
1221 : int nMaskFlagStop,
1222 : double *pdfDataPct)
1223 : {
1224 8 : if (eAccess == GA_Update)
1225 0 : FlushCache(false);
1226 :
1227 8 : const int iColMin = nXOff / nBlockXSize + m_poTPD->m_nShiftXTiles;
1228 16 : const int iColMax = (nXOff + nXSize - 1) / nBlockXSize +
1229 8 : m_poTPD->m_nShiftXTiles +
1230 8 : (m_poTPD->m_nShiftXPixelsMod ? 1 : 0);
1231 8 : const int iRowMin = nYOff / nBlockYSize + m_poTPD->m_nShiftYTiles;
1232 16 : const int iRowMax = (nYOff + nYSize - 1) / nBlockYSize +
1233 8 : m_poTPD->m_nShiftYTiles +
1234 8 : (m_poTPD->m_nShiftYPixelsMod ? 1 : 0);
1235 :
1236 8 : int iDBRowMin = m_poTPD->GetRowFromIntoTopConvention(iRowMin);
1237 8 : int iDBRowMax = m_poTPD->GetRowFromIntoTopConvention(iRowMax);
1238 8 : if (iDBRowMin > iDBRowMax)
1239 0 : std::swap(iDBRowMin, iDBRowMax);
1240 :
1241 16 : char *pszSQL = sqlite3_mprintf(
1242 : "SELECT tile_row, tile_column FROM \"%w\" "
1243 : "WHERE zoom_level = %d AND "
1244 : "(tile_row BETWEEN %d AND %d) AND "
1245 : "(tile_column BETWEEN %d AND %d)"
1246 : "%s",
1247 8 : m_poTPD->m_osRasterTable.c_str(), m_poTPD->m_nZoomLevel, iDBRowMin,
1248 : iDBRowMax, iColMin, iColMax,
1249 8 : !m_poTPD->m_osWHERE.empty()
1250 0 : ? CPLSPrintf(" AND (%s)", m_poTPD->m_osWHERE.c_str())
1251 : : "");
1252 :
1253 : #ifdef DEBUG_VERBOSE
1254 : CPLDebug("GPKG", "%s", pszSQL);
1255 : #endif
1256 :
1257 8 : sqlite3_stmt *hStmt = nullptr;
1258 : int rc =
1259 8 : SQLPrepareWithError(m_poTPD->IGetDB(), pszSQL, -1, &hStmt, nullptr);
1260 8 : if (rc != SQLITE_OK)
1261 : {
1262 0 : sqlite3_free(pszSQL);
1263 : return GDAL_DATA_COVERAGE_STATUS_UNIMPLEMENTED |
1264 0 : GDAL_DATA_COVERAGE_STATUS_DATA;
1265 : }
1266 8 : sqlite3_free(pszSQL);
1267 8 : rc = sqlite3_step(hStmt);
1268 16 : std::set<std::pair<int, int>> oSetTiles; // (row, col)
1269 12 : while (rc == SQLITE_ROW)
1270 : {
1271 0 : oSetTiles.insert(std::pair(
1272 4 : m_poTPD->GetRowFromIntoTopConvention(sqlite3_column_int(hStmt, 0)),
1273 8 : sqlite3_column_int(hStmt, 1)));
1274 4 : rc = sqlite3_step(hStmt);
1275 : }
1276 8 : sqlite3_finalize(hStmt);
1277 8 : if (rc != SQLITE_DONE)
1278 : {
1279 : return GDAL_DATA_COVERAGE_STATUS_UNIMPLEMENTED |
1280 0 : GDAL_DATA_COVERAGE_STATUS_DATA;
1281 : }
1282 8 : if (oSetTiles.empty())
1283 : {
1284 4 : if (pdfDataPct)
1285 4 : *pdfDataPct = 0.0;
1286 4 : return GDAL_DATA_COVERAGE_STATUS_EMPTY;
1287 : }
1288 :
1289 8 : if (m_poTPD->m_nShiftXPixelsMod == 0 && m_poTPD->m_nShiftYPixelsMod == 0 &&
1290 4 : oSetTiles.size() == static_cast<size_t>(iRowMax - iRowMin + 1) *
1291 4 : (iColMax - iColMin + 1))
1292 : {
1293 2 : if (pdfDataPct)
1294 2 : *pdfDataPct = 100.0;
1295 2 : return GDAL_DATA_COVERAGE_STATUS_DATA;
1296 : }
1297 :
1298 2 : if (m_poTPD->m_nShiftXPixelsMod == 0 && m_poTPD->m_nShiftYPixelsMod == 0)
1299 : {
1300 2 : int nStatus = 0;
1301 2 : GIntBig nPixelsData = 0;
1302 7 : for (int iY = iRowMin; iY <= iRowMax; ++iY)
1303 : {
1304 21 : for (int iX = iColMin; iX <= iColMax; ++iX)
1305 : {
1306 : // coverity[swapped_arguments]
1307 16 : if (oSetTiles.find(std::pair(iY, iX)) == oSetTiles.end())
1308 : {
1309 14 : nStatus |= GDAL_DATA_COVERAGE_STATUS_EMPTY;
1310 : }
1311 : else
1312 : {
1313 2 : const int iXGDAL = iX - m_poTPD->m_nShiftXTiles;
1314 2 : const int iYGDAL = iY - m_poTPD->m_nShiftYTiles;
1315 2 : const int nXBlockRight =
1316 2 : (iXGDAL * nBlockXSize > INT_MAX - nBlockXSize)
1317 2 : ? INT_MAX
1318 2 : : (iXGDAL + 1) * nBlockXSize;
1319 2 : const int nYBlockBottom =
1320 2 : (iYGDAL * nBlockYSize > INT_MAX - nBlockYSize)
1321 2 : ? INT_MAX
1322 2 : : (iYGDAL + 1) * nBlockYSize;
1323 :
1324 4 : nPixelsData += (static_cast<GIntBig>(std::min(
1325 2 : nXBlockRight, nXOff + nXSize)) -
1326 2 : std::max(iXGDAL * nBlockXSize, nXOff)) *
1327 2 : (std::min(nYBlockBottom, nYOff + nYSize) -
1328 2 : std::max(iYGDAL * nBlockYSize, nYOff));
1329 2 : nStatus |= GDAL_DATA_COVERAGE_STATUS_DATA;
1330 : }
1331 16 : if (nMaskFlagStop != 0 && (nMaskFlagStop & nStatus) != 0)
1332 : {
1333 0 : if (pdfDataPct)
1334 0 : *pdfDataPct = -1.0;
1335 0 : return nStatus;
1336 : }
1337 : }
1338 : }
1339 :
1340 2 : if (pdfDataPct)
1341 : {
1342 2 : *pdfDataPct =
1343 2 : 100.0 * nPixelsData / (static_cast<GIntBig>(nXSize) * nYSize);
1344 : }
1345 2 : return nStatus;
1346 : }
1347 :
1348 : return GDAL_DATA_COVERAGE_STATUS_UNIMPLEMENTED |
1349 0 : GDAL_DATA_COVERAGE_STATUS_DATA;
1350 : }
1351 :
1352 : /************************************************************************/
1353 : /* WEBPSupports4Bands() */
1354 : /************************************************************************/
1355 :
1356 83 : static bool WEBPSupports4Bands()
1357 : {
1358 : static int bRes = -1;
1359 83 : if (bRes < 0)
1360 : {
1361 1 : GDALDriver *poDrv = GDALDriver::FromHandle(GDALGetDriverByName("WEBP"));
1362 2 : if (poDrv == nullptr ||
1363 1 : CPLTestBool(CPLGetConfigOption("GPKG_SIMUL_WEBP_3BAND", "FALSE")))
1364 0 : bRes = false;
1365 : else
1366 : {
1367 : // LOSSLESS and RGBA support appeared in the same version
1368 1 : bRes = strstr(poDrv->GetMetadataItem(GDAL_DMD_CREATIONOPTIONLIST),
1369 1 : "LOSSLESS") != nullptr;
1370 : }
1371 1 : if (poDrv != nullptr && !bRes)
1372 : {
1373 0 : CPLError(
1374 : CE_Warning, CPLE_AppDefined,
1375 : "The version of WEBP available does not support 4-band RGBA");
1376 : }
1377 : }
1378 83 : return CPL_TO_BOOL(bRes);
1379 : }
1380 :
1381 : /************************************************************************/
1382 : /* GetTileId() */
1383 : /************************************************************************/
1384 :
1385 58 : GIntBig GDALGPKGMBTilesLikePseudoDataset::GetTileId(int nRow, int nCol)
1386 : {
1387 : char *pszSQL =
1388 58 : sqlite3_mprintf("SELECT id FROM \"%w\" WHERE zoom_level = %d AND "
1389 : "tile_row = %d AND tile_column = %d",
1390 : m_osRasterTable.c_str(), m_nZoomLevel,
1391 58 : GetRowFromIntoTopConvention(nRow), nCol);
1392 58 : GIntBig nRes = SQLGetInteger64(IGetDB(), pszSQL, nullptr);
1393 58 : sqlite3_free(pszSQL);
1394 58 : return nRes;
1395 : }
1396 :
1397 : /************************************************************************/
1398 : /* DeleteTile() */
1399 : /************************************************************************/
1400 :
1401 460 : bool GDALGPKGMBTilesLikePseudoDataset::DeleteTile(int nRow, int nCol)
1402 : {
1403 : char *pszSQL =
1404 460 : sqlite3_mprintf("DELETE FROM \"%w\" "
1405 : "WHERE zoom_level = %d AND tile_row = %d AND "
1406 : "tile_column = %d",
1407 : m_osRasterTable.c_str(), m_nZoomLevel,
1408 460 : GetRowFromIntoTopConvention(nRow), nCol);
1409 : #ifdef DEBUG_VERBOSE
1410 : CPLDebug("GPKG", "%s", pszSQL);
1411 : #endif
1412 460 : char *pszErrMsg = nullptr;
1413 460 : int rc = sqlite3_exec(IGetDB(), pszSQL, nullptr, nullptr, &pszErrMsg);
1414 460 : if (rc != SQLITE_OK)
1415 : {
1416 0 : CPLError(CE_Failure, CPLE_AppDefined,
1417 : "Failure when deleting tile (row=%d,col=%d) "
1418 : "at zoom_level=%d : %s",
1419 0 : GetRowFromIntoTopConvention(nRow), nCol, m_nZoomLevel,
1420 0 : pszErrMsg ? pszErrMsg : "");
1421 : }
1422 460 : sqlite3_free(pszSQL);
1423 460 : sqlite3_free(pszErrMsg);
1424 460 : return rc == SQLITE_OK;
1425 : }
1426 :
1427 : /************************************************************************/
1428 : /* DeleteFromGriddedTileAncillary() */
1429 : /************************************************************************/
1430 :
1431 58 : bool GDALGPKGMBTilesLikePseudoDataset::DeleteFromGriddedTileAncillary(
1432 : GIntBig nTileId)
1433 : {
1434 : char *pszSQL =
1435 58 : sqlite3_mprintf("DELETE FROM gpkg_2d_gridded_tile_ancillary WHERE "
1436 : "tpudt_name = '%q' AND tpudt_id = ?",
1437 : m_osRasterTable.c_str());
1438 58 : sqlite3_stmt *hStmt = nullptr;
1439 58 : int rc = SQLPrepareWithError(IGetDB(), pszSQL, -1, &hStmt, nullptr);
1440 58 : if (rc == SQLITE_OK)
1441 : {
1442 58 : sqlite3_bind_int64(hStmt, 1, nTileId);
1443 58 : rc = sqlite3_step(hStmt);
1444 58 : sqlite3_finalize(hStmt);
1445 : }
1446 58 : sqlite3_free(pszSQL);
1447 58 : return rc == SQLITE_OK;
1448 : }
1449 :
1450 : /************************************************************************/
1451 : /* ProcessInt16UInt16Tile() */
1452 : /************************************************************************/
1453 :
1454 : template <class T>
1455 48 : static void ProcessInt16UInt16Tile(
1456 : const void *pabyData, GPtrDiff_t nPixels, bool bIsInt16, bool bHasNoData,
1457 : double dfNoDataValue, GUInt16 usGPKGNull, double m_dfOffset,
1458 : double m_dfScale, GUInt16 *pTempTileBuffer, double &dfTileOffset,
1459 : double &dfTileScale, double &dfTileMin, double &dfTileMax,
1460 : double &dfTileMean, double &dfTileStdDev, GPtrDiff_t &nValidPixels)
1461 : {
1462 48 : const T *pSrc = reinterpret_cast<const T *>(pabyData);
1463 48 : T nMin = 0;
1464 48 : T nMax = 0;
1465 48 : double dfM2 = 0.0;
1466 3145777 : for (int i = 0; i < nPixels; i++)
1467 : {
1468 3145726 : const T nVal = pSrc[i];
1469 3145726 : if (bHasNoData && nVal == dfNoDataValue)
1470 718161 : continue;
1471 :
1472 2427565 : if (nValidPixels == 0)
1473 : {
1474 47 : nMin = nVal;
1475 47 : nMax = nVal;
1476 : }
1477 : else
1478 : {
1479 2427515 : nMin = std::min(nMin, nVal);
1480 2427515 : nMax = std::max(nMax, nVal);
1481 : }
1482 2427565 : nValidPixels++;
1483 2427565 : const double dfDelta = nVal - dfTileMean;
1484 2427565 : dfTileMean += dfDelta / nValidPixels;
1485 2427565 : dfM2 += dfDelta * (nVal - dfTileMean);
1486 : }
1487 48 : dfTileMin = nMin;
1488 48 : dfTileMax = nMax;
1489 48 : if (nValidPixels)
1490 47 : dfTileStdDev = sqrt(dfM2 / nValidPixels);
1491 :
1492 48 : double dfGlobalMin = (nMin - m_dfOffset) / m_dfScale;
1493 48 : double dfGlobalMax = (nMax - m_dfOffset) / m_dfScale;
1494 48 : double dfRange = 65535.0;
1495 48 : if (bHasNoData && usGPKGNull == 65535 &&
1496 8 : dfGlobalMax - dfGlobalMin >= dfRange)
1497 : {
1498 0 : dfRange = 65534.0;
1499 : }
1500 :
1501 48 : if (dfGlobalMax - dfGlobalMin > dfRange)
1502 : {
1503 0 : dfTileScale = (dfGlobalMax - dfGlobalMin) / dfRange;
1504 : }
1505 48 : if (dfGlobalMin < 0.0)
1506 : {
1507 0 : dfTileOffset = -dfGlobalMin;
1508 : }
1509 48 : else if (dfGlobalMax / dfTileScale > dfRange)
1510 : {
1511 0 : dfTileOffset = dfGlobalMax - dfRange * dfTileScale;
1512 : }
1513 :
1514 48 : if (bHasNoData && std::numeric_limits<T>::min() == 0 && m_dfOffset == 0.0 &&
1515 : m_dfScale == 1.0)
1516 : {
1517 7 : dfTileOffset = 0.0;
1518 7 : dfTileScale = 1.0;
1519 : }
1520 41 : else if (bHasNoData && bIsInt16 && dfNoDataValue == -32768.0 &&
1521 1 : usGPKGNull == 65535 && m_dfOffset == -32768.0 && m_dfScale == 1.0)
1522 : {
1523 1 : dfTileOffset = 1.0;
1524 1 : dfTileScale = 1.0;
1525 : }
1526 :
1527 3145777 : for (GPtrDiff_t i = 0; i < nPixels; i++)
1528 : {
1529 3145726 : const T nVal = pSrc[i];
1530 3145726 : if (bHasNoData && nVal == dfNoDataValue)
1531 718161 : pTempTileBuffer[i] = usGPKGNull;
1532 : else
1533 : {
1534 2427565 : double dfVal =
1535 2427565 : ((nVal - m_dfOffset) / m_dfScale - dfTileOffset) / dfTileScale;
1536 2427565 : CPLAssert(dfVal >= 0.0 && dfVal < 65535.5);
1537 2427565 : pTempTileBuffer[i] = static_cast<GUInt16>(dfVal + 0.5);
1538 2427565 : if (bHasNoData && pTempTileBuffer[i] == usGPKGNull)
1539 : {
1540 0 : if (usGPKGNull > 0)
1541 0 : pTempTileBuffer[i]--;
1542 : else
1543 0 : pTempTileBuffer[i]++;
1544 : ;
1545 : }
1546 : }
1547 : }
1548 48 : }
1549 :
1550 : /************************************************************************/
1551 : /* WriteTile() */
1552 : /************************************************************************/
1553 :
1554 13456 : CPLErr GDALGPKGMBTilesLikePseudoDataset::WriteTile()
1555 : {
1556 13456 : GDALGPKGMBTilesLikePseudoDataset *poMainDS =
1557 13456 : m_poParentDS ? m_poParentDS : this;
1558 13456 : if (poMainDS->m_nTileInsertionCount < 0)
1559 501 : return CE_Failure;
1560 :
1561 12955 : if (m_bInWriteTile)
1562 : {
1563 : // Shouldn't happen in practice, but #7022 shows that the unexpected
1564 : // can happen sometimes.
1565 0 : CPLError(
1566 : CE_Failure, CPLE_AppDefined,
1567 : "Recursive call to GDALGPKGMBTilesLikePseudoDataset::WriteTile()");
1568 0 : return CE_Failure;
1569 : }
1570 12955 : GDALRasterBlock::EnterDisableDirtyBlockFlush();
1571 12955 : m_bInWriteTile = true;
1572 12955 : CPLErr eErr = WriteTileInternal();
1573 12955 : m_bInWriteTile = false;
1574 12955 : GDALRasterBlock::LeaveDisableDirtyBlockFlush();
1575 12955 : return eErr;
1576 : }
1577 :
1578 : /* should only be called by WriteTile() */
1579 12955 : CPLErr GDALGPKGMBTilesLikePseudoDataset::WriteTileInternal()
1580 : {
1581 17579 : if (!(IGetUpdate() && m_asCachedTilesDesc[0].nRow >= 0 &&
1582 4624 : m_asCachedTilesDesc[0].nCol >= 0 &&
1583 4624 : m_asCachedTilesDesc[0].nIdxWithinTileData == 0))
1584 8331 : return CE_None;
1585 :
1586 4624 : int nRow = m_asCachedTilesDesc[0].nRow;
1587 4624 : int nCol = m_asCachedTilesDesc[0].nCol;
1588 :
1589 4624 : bool bAllDirty = true;
1590 4624 : bool bAllNonDirty = true;
1591 4624 : const int nBands = IGetRasterCount();
1592 18006 : for (int i = 0; i < nBands; i++)
1593 : {
1594 13382 : if (m_asCachedTilesDesc[0].abBandDirty[i])
1595 13349 : bAllNonDirty = false;
1596 : else
1597 33 : bAllDirty = false;
1598 : }
1599 4624 : if (bAllNonDirty)
1600 0 : return CE_None;
1601 :
1602 : int nBlockXSize, nBlockYSize;
1603 4624 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
1604 :
1605 : /* If all bands for that block are not dirty/written, we need to */
1606 : /* fetch the missing ones if the tile exists */
1607 4624 : bool bIsLossyFormat = false;
1608 4624 : const size_t nBandBlockSize =
1609 4624 : static_cast<size_t>(nBlockXSize) * nBlockYSize * m_nDTSize;
1610 4624 : if (!bAllDirty)
1611 : {
1612 72 : for (int i = 1; i <= 3; i++)
1613 : {
1614 54 : m_asCachedTilesDesc[i].nRow = -1;
1615 54 : m_asCachedTilesDesc[i].nCol = -1;
1616 54 : m_asCachedTilesDesc[i].nIdxWithinTileData = -1;
1617 : }
1618 18 : const int nTileBands = m_eDT == GDT_Byte ? 4 : 1;
1619 18 : GByte *pabyTemp = m_pabyCachedTiles + nTileBands * nBandBlockSize;
1620 18 : ReadTile(nRow, nCol, pabyTemp, &bIsLossyFormat);
1621 82 : for (int i = 0; i < nBands; i++)
1622 : {
1623 64 : if (!m_asCachedTilesDesc[0].abBandDirty[i])
1624 : {
1625 33 : memcpy(m_pabyCachedTiles + i * nBandBlockSize,
1626 33 : pabyTemp + i * nBandBlockSize, nBandBlockSize);
1627 : }
1628 : }
1629 : }
1630 :
1631 : /* Compute origin of tile in GDAL raster space */
1632 4624 : int nXOff = (nCol - m_nShiftXTiles) * nBlockXSize - m_nShiftXPixelsMod;
1633 4624 : int nYOff = (nRow - m_nShiftYTiles) * nBlockYSize - m_nShiftYPixelsMod;
1634 :
1635 : /* Assert that the tile at least intersects some of the GDAL raster space */
1636 4624 : CPLAssert(nXOff > -nBlockXSize);
1637 4624 : CPLAssert(nYOff > -nBlockYSize);
1638 : /* Can happen if the tile of the raster is less than the block size */
1639 4624 : const int nRasterXSize = IGetRasterBand(1)->GetXSize();
1640 4624 : const int nRasterYSize = IGetRasterBand(1)->GetYSize();
1641 4624 : if (nXOff >= nRasterXSize || nYOff >= nRasterYSize)
1642 0 : return CE_None;
1643 :
1644 : #ifdef DEBUG_VERBOSE
1645 : if (m_nShiftXPixelsMod == 0 && m_nShiftYPixelsMod == 0 && m_eDT == GDT_Byte)
1646 : {
1647 : int nBlockXOff = nCol;
1648 : int nBlockYOff = nRow;
1649 : if (nBlockXOff * nBlockXSize <= nRasterXSize - nBlockXSize &&
1650 : nBlockYOff * nBlockYSize > nRasterYSize - nBlockYSize)
1651 : {
1652 : for (int i = 0; i < nBands; i++)
1653 : {
1654 : bool bFoundNonZero = false;
1655 : for (GPtrDiff_t y = nRasterYSize - nBlockYOff * nBlockYSize;
1656 : y < nBlockYSize; y++)
1657 : {
1658 : for (int x = 0; x < nBlockXSize; x++)
1659 : {
1660 : if (m_pabyCachedTiles[y * nBlockXSize + x +
1661 : i * nBandBlockSize] != 0 &&
1662 : !bFoundNonZero)
1663 : {
1664 : CPLDebug("GPKG",
1665 : "WriteTileInternal(): Found non-zero "
1666 : "content in ghost part of "
1667 : "tile(band=%d,nBlockXOff=%d,nBlockYOff=%d,"
1668 : "m_nZoomLevel=%d)\n",
1669 : i + 1, nBlockXOff, nBlockYOff,
1670 : m_nZoomLevel);
1671 : bFoundNonZero = true;
1672 : }
1673 : }
1674 : }
1675 : }
1676 : }
1677 : }
1678 : #endif
1679 :
1680 : /* Validity area of tile data in intra-tile coordinate space */
1681 4624 : int iXOff = 0;
1682 4624 : GPtrDiff_t iYOff = 0;
1683 4624 : int iXCount = nBlockXSize;
1684 4624 : int iYCount = nBlockYSize;
1685 :
1686 4624 : bool bPartialTile = false;
1687 4624 : int nAlphaBand = (nBands == 2) ? 2 : (nBands == 4) ? 4 : 0;
1688 4624 : if (nAlphaBand == 0)
1689 : {
1690 3742 : if (nXOff < 0)
1691 : {
1692 4 : bPartialTile = true;
1693 4 : iXOff = -nXOff;
1694 4 : iXCount += nXOff;
1695 : }
1696 3742 : if (nXOff > nRasterXSize - nBlockXSize)
1697 : {
1698 93 : bPartialTile = true;
1699 93 : iXCount -= static_cast<int>(static_cast<GIntBig>(nXOff) +
1700 93 : nBlockXSize - nRasterXSize);
1701 : }
1702 3742 : if (nYOff < 0)
1703 : {
1704 8 : bPartialTile = true;
1705 8 : iYOff = -nYOff;
1706 8 : iYCount += nYOff;
1707 : }
1708 3742 : if (nYOff > nRasterYSize - nBlockYSize)
1709 : {
1710 107 : bPartialTile = true;
1711 107 : iYCount -= static_cast<int>(static_cast<GIntBig>(nYOff) +
1712 107 : nBlockYSize - nRasterYSize);
1713 : }
1714 3742 : CPLAssert(iXOff >= 0);
1715 3742 : CPLAssert(iYOff >= 0);
1716 3742 : CPLAssert(iXCount > 0);
1717 3742 : CPLAssert(iYCount > 0);
1718 3742 : CPLAssert(iXOff + iXCount <= nBlockXSize);
1719 3742 : CPLAssert(iYOff + iYCount <= nBlockYSize);
1720 : }
1721 :
1722 4624 : m_asCachedTilesDesc[0].nRow = -1;
1723 4624 : m_asCachedTilesDesc[0].nCol = -1;
1724 4624 : m_asCachedTilesDesc[0].nIdxWithinTileData = -1;
1725 4624 : m_asCachedTilesDesc[0].abBandDirty[0] = false;
1726 4624 : m_asCachedTilesDesc[0].abBandDirty[1] = false;
1727 4624 : m_asCachedTilesDesc[0].abBandDirty[2] = false;
1728 4624 : m_asCachedTilesDesc[0].abBandDirty[3] = false;
1729 :
1730 4624 : CPLErr eErr = CE_Failure;
1731 :
1732 4624 : int bHasNoData = FALSE;
1733 4624 : double dfNoDataValue = IGetRasterBand(1)->GetNoDataValue(&bHasNoData);
1734 4624 : const bool bHasNanNoData = bHasNoData && std::isnan(dfNoDataValue);
1735 :
1736 4624 : bool bAllOpaque = true;
1737 : // Detect fully transparent tiles, but only if all bands are dirty (that is
1738 : // the user wrote content into them) !
1739 4624 : if (bAllDirty && m_eDT == GDT_Byte && m_poCT == nullptr && nAlphaBand != 0)
1740 : {
1741 872 : GByte byFirstAlphaVal =
1742 872 : m_pabyCachedTiles[(nAlphaBand - 1) * nBlockXSize * nBlockYSize];
1743 872 : GPtrDiff_t i = 1;
1744 24105900 : for (; i < static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize; i++)
1745 : {
1746 24105400 : if (m_pabyCachedTiles[static_cast<GPtrDiff_t>(nAlphaBand - 1) *
1747 24105400 : nBlockXSize * nBlockYSize +
1748 24105400 : i] != byFirstAlphaVal)
1749 367 : break;
1750 : }
1751 872 : if (i == static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize)
1752 : {
1753 : // If tile is fully transparent, don't serialize it and remove it if
1754 : // it exists
1755 505 : if (byFirstAlphaVal == 0)
1756 : {
1757 440 : DeleteTile(nRow, nCol);
1758 :
1759 440 : return CE_None;
1760 : }
1761 65 : bAllOpaque = (byFirstAlphaVal == 255);
1762 : }
1763 : else
1764 432 : bAllOpaque = false;
1765 : }
1766 3752 : else if (bAllDirty && m_eDT == GDT_Byte && m_poCT == nullptr &&
1767 3654 : (!bHasNoData || dfNoDataValue == 0.0))
1768 : {
1769 3654 : bool bAllEmpty = true;
1770 3654 : const auto nPixels =
1771 3654 : static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize * nBands;
1772 1207860 : for (GPtrDiff_t i = 0; i < nPixels; i++)
1773 : {
1774 1207840 : if (m_pabyCachedTiles[i] != 0)
1775 : {
1776 3637 : bAllEmpty = false;
1777 3637 : break;
1778 : }
1779 : }
1780 3654 : if (bAllEmpty)
1781 : {
1782 : // If tile is fully transparent, don't serialize it and remove it if
1783 : // it exists
1784 17 : DeleteTile(nRow, nCol);
1785 :
1786 17 : return CE_None;
1787 3637 : }
1788 : }
1789 98 : else if (bAllDirty && m_eDT == GDT_Float32)
1790 : {
1791 12 : const float *pSrc = reinterpret_cast<float *>(m_pabyCachedTiles);
1792 : GPtrDiff_t i;
1793 12 : const float fNoDataValueOrZero =
1794 12 : bHasNoData ? static_cast<float>(dfNoDataValue) : 0.0f;
1795 131084 : for (i = 0; i < static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize; i++)
1796 : {
1797 131082 : const float fVal = pSrc[i];
1798 131082 : if (bHasNanNoData)
1799 : {
1800 0 : if (std::isnan(fVal))
1801 0 : continue;
1802 : }
1803 131082 : else if (fVal == fNoDataValueOrZero)
1804 : {
1805 131072 : continue;
1806 : }
1807 10 : break;
1808 : }
1809 12 : if (i == static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize)
1810 : {
1811 : // If tile is fully transparent, don't serialize it and remove it if
1812 : // it exists
1813 2 : DeleteTile(nRow, nCol);
1814 :
1815 2 : return CE_None;
1816 : }
1817 : }
1818 :
1819 4165 : if (bIsLossyFormat)
1820 : {
1821 1 : CPLDebug("GPKG",
1822 : "Had to read tile (row=%d,col=%d) at zoom_level=%d, "
1823 : "stored in a lossy format, before rewriting it, causing "
1824 : "potential extra quality loss",
1825 : nRow, nCol, m_nZoomLevel);
1826 : }
1827 :
1828 : const CPLString osMemFileName(
1829 8330 : VSIMemGenerateHiddenFilename("gpkg_write_tile"));
1830 4165 : const char *pszDriverName = "PNG";
1831 4165 : CPL_IGNORE_RET_VAL(pszDriverName); // Make CSA happy
1832 4165 : bool bTileDriverSupports1Band = false;
1833 4165 : bool bTileDriverSupports2Bands = false;
1834 4165 : bool bTileDriverSupports4Bands = false;
1835 4165 : bool bTileDriverSupportsCT = false;
1836 :
1837 4165 : if (nBands == 1 && m_eDT == GDT_Byte)
1838 76 : IGetRasterBand(1)->GetColorTable();
1839 :
1840 4165 : GDALDataType eTileDT = GDT_Byte;
1841 : // If not all bands are dirty, then (temporarily) use a lossless format
1842 4165 : if (m_eTF == GPKG_TF_PNG_JPEG || (!bAllDirty && m_eTF == GPKG_TF_JPEG))
1843 : {
1844 97 : bTileDriverSupports1Band = true;
1845 97 : if (bPartialTile || !bAllDirty || (nBands == 2 && !bAllOpaque) ||
1846 19 : (nBands == 4 && !bAllOpaque) || m_poCT != nullptr)
1847 : {
1848 89 : pszDriverName = "PNG";
1849 89 : bTileDriverSupports2Bands = m_bPNGSupports2Bands;
1850 89 : bTileDriverSupports4Bands = true;
1851 89 : bTileDriverSupportsCT = m_bPNGSupportsCT;
1852 : }
1853 : else
1854 8 : pszDriverName = "JPEG";
1855 : }
1856 4068 : else if (m_eTF == GPKG_TF_PNG || m_eTF == GPKG_TF_PNG8)
1857 : {
1858 3845 : pszDriverName = "PNG";
1859 3845 : bTileDriverSupports1Band = true;
1860 3845 : bTileDriverSupports2Bands = m_bPNGSupports2Bands;
1861 3845 : bTileDriverSupports4Bands = true;
1862 3845 : bTileDriverSupportsCT = m_bPNGSupportsCT;
1863 : }
1864 223 : else if (m_eTF == GPKG_TF_JPEG)
1865 : {
1866 82 : pszDriverName = "JPEG";
1867 82 : bTileDriverSupports1Band = true;
1868 : }
1869 141 : else if (m_eTF == GPKG_TF_WEBP)
1870 : {
1871 83 : pszDriverName = "WEBP";
1872 83 : bTileDriverSupports4Bands = WEBPSupports4Bands();
1873 : }
1874 58 : else if (m_eTF == GPKG_TF_PNG_16BIT)
1875 : {
1876 53 : pszDriverName = "PNG";
1877 53 : eTileDT = GDT_UInt16;
1878 53 : bTileDriverSupports1Band = true;
1879 : }
1880 5 : else if (m_eTF == GPKG_TF_TIFF_32BIT_FLOAT)
1881 : {
1882 5 : pszDriverName = "GTiff";
1883 5 : eTileDT = GDT_Float32;
1884 5 : bTileDriverSupports1Band = true;
1885 : }
1886 : else
1887 : {
1888 0 : CPLAssert(false);
1889 : }
1890 :
1891 : GDALDriver *l_poDriver =
1892 4165 : GDALDriver::FromHandle(GDALGetDriverByName(pszDriverName));
1893 4165 : if (l_poDriver != nullptr)
1894 : {
1895 4164 : auto poMEMDS = MEMDataset::Create("", nBlockXSize, nBlockYSize, 0,
1896 : eTileDT, nullptr);
1897 4164 : int nTileBands = nBands;
1898 4164 : if (bPartialTile && nBands == 1 && m_poCT == nullptr &&
1899 : bTileDriverSupports2Bands)
1900 31 : nTileBands = 2;
1901 4133 : else if (bPartialTile && bTileDriverSupports4Bands)
1902 40 : nTileBands = 4;
1903 : // only use (somewhat lossy) PNG8 if all bands are dirty
1904 4093 : else if (bAllDirty && m_eTF == GPKG_TF_PNG8 && nBands >= 3 &&
1905 10 : bAllOpaque && !bPartialTile)
1906 10 : nTileBands = 1;
1907 4083 : else if (nBands == 2)
1908 : {
1909 383 : if (bAllOpaque)
1910 : {
1911 51 : if (bTileDriverSupports2Bands)
1912 30 : nTileBands = 1;
1913 : else
1914 21 : nTileBands = 3;
1915 : }
1916 332 : else if (!bTileDriverSupports2Bands)
1917 : {
1918 135 : if (bTileDriverSupports4Bands)
1919 68 : nTileBands = 4;
1920 : else
1921 67 : nTileBands = 3;
1922 : }
1923 : }
1924 3700 : else if (nBands == 4 && (bAllOpaque || !bTileDriverSupports4Bands))
1925 21 : nTileBands = 3;
1926 3679 : else if (nBands == 1 && m_poCT != nullptr && !bTileDriverSupportsCT)
1927 : {
1928 1 : nTileBands = 3;
1929 1 : if (bTileDriverSupports4Bands)
1930 : {
1931 0 : for (int i = 0; i < m_poCT->GetColorEntryCount(); i++)
1932 : {
1933 0 : const GDALColorEntry *psEntry = m_poCT->GetColorEntry(i);
1934 0 : if (psEntry->c4 == 0)
1935 : {
1936 0 : nTileBands = 4;
1937 0 : break;
1938 : }
1939 : }
1940 1 : }
1941 : }
1942 3678 : else if (nBands == 1 && m_poCT == nullptr && !bTileDriverSupports1Band)
1943 1 : nTileBands = 3;
1944 :
1945 4164 : if (bPartialTile && (nTileBands == 2 || nTileBands == 4))
1946 : {
1947 71 : int nTargetAlphaBand = nTileBands;
1948 71 : memset(m_pabyCachedTiles + (nTargetAlphaBand - 1) * nBandBlockSize,
1949 : 0, nBandBlockSize);
1950 4356 : for (GPtrDiff_t iY = iYOff; iY < iYOff + iYCount; iY++)
1951 : {
1952 4285 : memset(m_pabyCachedTiles +
1953 4285 : (static_cast<size_t>(nTargetAlphaBand - 1) *
1954 4285 : nBlockYSize +
1955 4285 : iY) *
1956 4285 : nBlockXSize +
1957 4285 : iXOff,
1958 : 255, iXCount);
1959 : }
1960 : }
1961 :
1962 4164 : GUInt16 *pTempTileBuffer = nullptr;
1963 4164 : GPtrDiff_t nValidPixels = 0;
1964 4164 : double dfTileMin = 0.0;
1965 4164 : double dfTileMax = 0.0;
1966 4164 : double dfTileMean = 0.0;
1967 4164 : double dfTileStdDev = 0.0;
1968 4164 : double dfTileOffset = 0.0;
1969 4164 : double dfTileScale = 1.0;
1970 4164 : if (m_eTF == GPKG_TF_PNG_16BIT)
1971 : {
1972 : pTempTileBuffer = static_cast<GUInt16 *>(
1973 53 : VSI_MALLOC3_VERBOSE(2, nBlockXSize, nBlockYSize));
1974 :
1975 53 : if (m_eDT == GDT_Int16)
1976 : {
1977 74 : ProcessInt16UInt16Tile<GInt16>(
1978 37 : m_pabyCachedTiles,
1979 37 : static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize, true,
1980 37 : CPL_TO_BOOL(bHasNoData), dfNoDataValue, m_usGPKGNull,
1981 : m_dfOffset, m_dfScale, pTempTileBuffer, dfTileOffset,
1982 : dfTileScale, dfTileMin, dfTileMax, dfTileMean, dfTileStdDev,
1983 : nValidPixels);
1984 : }
1985 16 : else if (m_eDT == GDT_UInt16)
1986 : {
1987 22 : ProcessInt16UInt16Tile<GUInt16>(
1988 11 : m_pabyCachedTiles,
1989 11 : static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize, false,
1990 11 : CPL_TO_BOOL(bHasNoData), dfNoDataValue, m_usGPKGNull,
1991 : m_dfOffset, m_dfScale, pTempTileBuffer, dfTileOffset,
1992 : dfTileScale, dfTileMin, dfTileMax, dfTileMean, dfTileStdDev,
1993 : nValidPixels);
1994 : }
1995 5 : else if (m_eDT == GDT_Float32)
1996 : {
1997 5 : const float *pSrc =
1998 : reinterpret_cast<float *>(m_pabyCachedTiles);
1999 5 : float fMin = 0.0f;
2000 5 : float fMax = 0.0f;
2001 5 : double dfM2 = 0.0;
2002 327685 : for (GPtrDiff_t i = 0;
2003 327685 : i < static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize;
2004 : i++)
2005 : {
2006 327680 : const float fVal = pSrc[i];
2007 327680 : if (bHasNanNoData)
2008 : {
2009 0 : if (std::isnan(fVal))
2010 196206 : continue;
2011 : }
2012 720494 : else if (bHasNoData &&
2013 327680 : fVal == static_cast<float>(dfNoDataValue))
2014 : {
2015 196206 : continue;
2016 : }
2017 131474 : if (std::isinf(fVal))
2018 0 : continue;
2019 :
2020 131474 : if (nValidPixels == 0)
2021 : {
2022 5 : fMin = fVal;
2023 5 : fMax = fVal;
2024 : }
2025 : else
2026 : {
2027 131469 : fMin = std::min(fMin, fVal);
2028 131469 : fMax = std::max(fMax, fVal);
2029 : }
2030 131474 : nValidPixels++;
2031 131474 : const double dfDelta = fVal - dfTileMean;
2032 131474 : dfTileMean += dfDelta / nValidPixels;
2033 131474 : dfM2 += dfDelta * (fVal - dfTileMean);
2034 : }
2035 5 : dfTileMin = fMin;
2036 5 : dfTileMax = fMax;
2037 5 : if (nValidPixels)
2038 5 : dfTileStdDev = sqrt(dfM2 / nValidPixels);
2039 :
2040 5 : double dfGlobalMin = (fMin - m_dfOffset) / m_dfScale;
2041 5 : double dfGlobalMax = (fMax - m_dfOffset) / m_dfScale;
2042 5 : if (dfGlobalMax > dfGlobalMin)
2043 : {
2044 4 : if (bHasNoData && m_usGPKGNull == 65535 &&
2045 2 : dfGlobalMax - dfGlobalMin >= 65534.0)
2046 : {
2047 1 : dfTileOffset = dfGlobalMin;
2048 1 : dfTileScale = (dfGlobalMax - dfGlobalMin) / 65534.0;
2049 : }
2050 3 : else if (bHasNoData && m_usGPKGNull == 0 &&
2051 0 : (dfNoDataValue - m_dfOffset) / m_dfScale != 0)
2052 : {
2053 0 : dfTileOffset =
2054 0 : (65535.0 * dfGlobalMin - dfGlobalMax) / 65534.0;
2055 0 : dfTileScale = dfGlobalMin - dfTileOffset;
2056 : }
2057 : else
2058 : {
2059 3 : dfTileOffset = dfGlobalMin;
2060 3 : dfTileScale = (dfGlobalMax - dfGlobalMin) / 65535.0;
2061 : }
2062 : }
2063 : else
2064 : {
2065 1 : dfTileOffset = dfGlobalMin;
2066 : }
2067 :
2068 327685 : for (GPtrDiff_t i = 0;
2069 327685 : i < static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize;
2070 : i++)
2071 : {
2072 327680 : const float fVal = pSrc[i];
2073 327680 : if (bHasNanNoData)
2074 : {
2075 0 : if (std::isnan(fVal))
2076 : {
2077 0 : pTempTileBuffer[i] = m_usGPKGNull;
2078 0 : continue;
2079 : }
2080 : }
2081 327680 : else if (bHasNoData)
2082 : {
2083 196608 : if (fVal == static_cast<float>(dfNoDataValue))
2084 : {
2085 196206 : pTempTileBuffer[i] = m_usGPKGNull;
2086 196206 : continue;
2087 : }
2088 : }
2089 : double dfVal =
2090 131474 : std::isfinite(fVal)
2091 131474 : ? ((fVal - m_dfOffset) / m_dfScale - dfTileOffset) /
2092 : dfTileScale
2093 : : (fVal > 0) ? 65535
2094 131474 : : 0;
2095 131474 : CPLAssert(dfVal >= 0.0 && dfVal < 65535.5);
2096 131474 : pTempTileBuffer[i] = static_cast<GUInt16>(dfVal + 0.5);
2097 131474 : if (bHasNoData && pTempTileBuffer[i] == m_usGPKGNull)
2098 : {
2099 1 : if (m_usGPKGNull > 0)
2100 1 : pTempTileBuffer[i]--;
2101 : else
2102 0 : pTempTileBuffer[i]++;
2103 : }
2104 : }
2105 : }
2106 :
2107 53 : auto hBand = MEMCreateRasterBandEx(
2108 : poMEMDS, 1, reinterpret_cast<GByte *>(pTempTileBuffer),
2109 : GDT_UInt16, 0, 0, false);
2110 53 : poMEMDS->AddMEMBand(hBand);
2111 : }
2112 4111 : else if (m_eTF == GPKG_TF_TIFF_32BIT_FLOAT)
2113 : {
2114 5 : const float *pSrc = reinterpret_cast<float *>(m_pabyCachedTiles);
2115 5 : float fMin = 0.0f;
2116 5 : float fMax = 0.0f;
2117 5 : double dfM2 = 0.0;
2118 327685 : for (GPtrDiff_t i = 0;
2119 327685 : i < static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize; i++)
2120 : {
2121 327680 : const float fVal = pSrc[i];
2122 327680 : if (bHasNanNoData)
2123 : {
2124 0 : if (std::isnan(fVal))
2125 65136 : continue;
2126 : }
2127 458352 : else if (bHasNoData &&
2128 327680 : fVal == static_cast<float>(dfNoDataValue))
2129 : {
2130 65136 : continue;
2131 : }
2132 :
2133 262544 : if (nValidPixels == 0)
2134 : {
2135 5 : fMin = fVal;
2136 5 : fMax = fVal;
2137 : }
2138 : else
2139 : {
2140 262539 : fMin = std::min(fMin, fVal);
2141 262539 : fMax = std::max(fMax, fVal);
2142 : }
2143 262544 : nValidPixels++;
2144 262544 : const double dfDelta = fVal - dfTileMean;
2145 262544 : dfTileMean += dfDelta / nValidPixels;
2146 262544 : dfM2 += dfDelta * (fVal - dfTileMean);
2147 : }
2148 5 : dfTileMin = fMin;
2149 5 : dfTileMax = fMax;
2150 5 : if (nValidPixels)
2151 5 : dfTileStdDev = sqrt(dfM2 / nValidPixels);
2152 :
2153 5 : auto hBand = MEMCreateRasterBandEx(poMEMDS, 1, m_pabyCachedTiles,
2154 : GDT_Float32, 0, 0, false);
2155 5 : poMEMDS->AddMEMBand(hBand);
2156 : }
2157 : else
2158 : {
2159 4106 : CPLAssert(m_eDT == GDT_Byte);
2160 16190 : for (int i = 0; i < nTileBands; i++)
2161 : {
2162 12084 : int iSrc = i;
2163 12084 : if (nBands == 1 && m_poCT == nullptr && nTileBands == 3)
2164 3 : iSrc = 0;
2165 12081 : else if (nBands == 1 && m_poCT == nullptr && bPartialTile &&
2166 : nTileBands == 4)
2167 8 : iSrc = (i < 3) ? 0 : 3;
2168 12073 : else if (nBands == 2 && nTileBands >= 3)
2169 536 : iSrc = (i < 3) ? 0 : 1;
2170 :
2171 24168 : auto hBand = MEMCreateRasterBandEx(
2172 : poMEMDS, i + 1,
2173 12084 : m_pabyCachedTiles + iSrc * nBlockXSize * nBlockYSize,
2174 : GDT_Byte, 0, 0, false);
2175 12084 : poMEMDS->AddMEMBand(hBand);
2176 :
2177 12084 : if (i == 0 && nTileBands == 1 && m_poCT != nullptr)
2178 15 : poMEMDS->GetRasterBand(1)->SetColorTable(m_poCT);
2179 : }
2180 : }
2181 :
2182 4164 : if ((m_eTF == GPKG_TF_PNG_16BIT || m_eTF == GPKG_TF_TIFF_32BIT_FLOAT) &&
2183 58 : nValidPixels == 0)
2184 : {
2185 : // If tile is fully transparent, don't serialize it and remove
2186 : // it if it exists.
2187 1 : GIntBig nId = GetTileId(nRow, nCol);
2188 1 : if (nId > 0)
2189 : {
2190 1 : DeleteTile(nRow, nCol);
2191 :
2192 1 : DeleteFromGriddedTileAncillary(nId);
2193 : }
2194 :
2195 1 : CPLFree(pTempTileBuffer);
2196 1 : delete poMEMDS;
2197 2 : return CE_None;
2198 : }
2199 :
2200 4163 : if (m_eTF == GPKG_TF_PNG8 && nTileBands == 1 && nBands >= 3)
2201 : {
2202 10 : CPLAssert(bAllDirty);
2203 :
2204 10 : auto poMEM_RGB_DS = MEMDataset::Create("", nBlockXSize, nBlockYSize,
2205 : 0, GDT_Byte, nullptr);
2206 40 : for (int i = 0; i < 3; i++)
2207 : {
2208 60 : auto hBand = MEMCreateRasterBandEx(
2209 30 : poMEMDS, i + 1, m_pabyCachedTiles + i * nBandBlockSize,
2210 : GDT_Byte, 0, 0, false);
2211 30 : poMEM_RGB_DS->AddMEMBand(hBand);
2212 : }
2213 :
2214 10 : if (m_pabyHugeColorArray == nullptr)
2215 : {
2216 5 : if (nBlockXSize <= 65536 / nBlockYSize)
2217 4 : m_pabyHugeColorArray =
2218 4 : VSIMalloc(MEDIAN_CUT_AND_DITHER_BUFFER_SIZE_65536);
2219 : else
2220 1 : m_pabyHugeColorArray =
2221 1 : VSIMalloc2(256 * 256 * 256, sizeof(GUInt32));
2222 : }
2223 :
2224 10 : GDALColorTable *poCT = new GDALColorTable();
2225 20 : GDALComputeMedianCutPCTInternal(
2226 10 : poMEM_RGB_DS->GetRasterBand(1), poMEM_RGB_DS->GetRasterBand(2),
2227 10 : poMEM_RGB_DS->GetRasterBand(3),
2228 : /*NULL, NULL, NULL,*/
2229 10 : m_pabyCachedTiles, m_pabyCachedTiles + nBandBlockSize,
2230 10 : m_pabyCachedTiles + 2 * nBandBlockSize, nullptr,
2231 : 256, /* max colors */
2232 : 8, /* bit depth */
2233 : static_cast<GUInt32 *>(
2234 10 : m_pabyHugeColorArray), /* preallocated histogram */
2235 : poCT, nullptr, nullptr);
2236 :
2237 20 : GDALDitherRGB2PCTInternal(
2238 10 : poMEM_RGB_DS->GetRasterBand(1), poMEM_RGB_DS->GetRasterBand(2),
2239 10 : poMEM_RGB_DS->GetRasterBand(3), poMEMDS->GetRasterBand(1), poCT,
2240 : 8, /* bit depth */
2241 : static_cast<GInt16 *>(
2242 10 : m_pabyHugeColorArray), /* pasDynamicColorMap */
2243 10 : m_bDither, nullptr, nullptr);
2244 10 : poMEMDS->GetRasterBand(1)->SetColorTable(poCT);
2245 10 : delete poCT;
2246 10 : GDALClose(poMEM_RGB_DS);
2247 : }
2248 4153 : else if (nBands == 1 && m_poCT != nullptr && nTileBands > 1)
2249 : {
2250 : GByte abyCT[4 * 256];
2251 7 : const int nEntries = std::min(256, m_poCT->GetColorEntryCount());
2252 1799 : for (int i = 0; i < nEntries; i++)
2253 : {
2254 1792 : const GDALColorEntry *psEntry = m_poCT->GetColorEntry(i);
2255 1792 : abyCT[4 * i] = static_cast<GByte>(psEntry->c1);
2256 1792 : abyCT[4 * i + 1] = static_cast<GByte>(psEntry->c2);
2257 1792 : abyCT[4 * i + 2] = static_cast<GByte>(psEntry->c3);
2258 1792 : abyCT[4 * i + 3] = static_cast<GByte>(psEntry->c4);
2259 : }
2260 7 : for (int i = nEntries; i < 256; i++)
2261 : {
2262 0 : abyCT[4 * i] = 0;
2263 0 : abyCT[4 * i + 1] = 0;
2264 0 : abyCT[4 * i + 2] = 0;
2265 0 : abyCT[4 * i + 3] = 0;
2266 : }
2267 7 : if (iYOff > 0)
2268 : {
2269 0 : memset(m_pabyCachedTiles + 0 * nBandBlockSize, 0,
2270 0 : nBlockXSize * iYOff);
2271 0 : memset(m_pabyCachedTiles + 1 * nBandBlockSize, 0,
2272 0 : nBlockXSize * iYOff);
2273 0 : memset(m_pabyCachedTiles + 2 * nBandBlockSize, 0,
2274 0 : nBlockXSize * iYOff);
2275 0 : memset(m_pabyCachedTiles + 3 * nBandBlockSize, 0,
2276 0 : nBlockXSize * iYOff);
2277 : }
2278 1057 : for (GPtrDiff_t iY = iYOff; iY < iYOff + iYCount; iY++)
2279 : {
2280 1050 : if (iXOff > 0)
2281 : {
2282 0 : const GPtrDiff_t i = iY * nBlockXSize;
2283 0 : memset(m_pabyCachedTiles + 0 * nBandBlockSize + i, 0,
2284 : iXOff);
2285 0 : memset(m_pabyCachedTiles + 1 * nBandBlockSize + i, 0,
2286 : iXOff);
2287 0 : memset(m_pabyCachedTiles + 2 * nBandBlockSize + i, 0,
2288 : iXOff);
2289 0 : memset(m_pabyCachedTiles + 3 * nBandBlockSize + i, 0,
2290 : iXOff);
2291 : }
2292 183250 : for (int iX = iXOff; iX < iXOff + iXCount; iX++)
2293 : {
2294 182200 : const GPtrDiff_t i = iY * nBlockXSize + iX;
2295 182200 : GByte byVal = m_pabyCachedTiles[i];
2296 182200 : m_pabyCachedTiles[i] = abyCT[4 * byVal];
2297 182200 : m_pabyCachedTiles[i + 1 * nBandBlockSize] =
2298 182200 : abyCT[4 * byVal + 1];
2299 182200 : m_pabyCachedTiles[i + 2 * nBandBlockSize] =
2300 182200 : abyCT[4 * byVal + 2];
2301 182200 : m_pabyCachedTiles[i + 3 * nBandBlockSize] =
2302 182200 : abyCT[4 * byVal + 3];
2303 : }
2304 1050 : if (iXOff + iXCount < nBlockXSize)
2305 : {
2306 800 : const GPtrDiff_t i = iY * nBlockXSize + iXOff + iXCount;
2307 800 : memset(m_pabyCachedTiles + 0 * nBandBlockSize + i, 0,
2308 800 : nBlockXSize - (iXOff + iXCount));
2309 800 : memset(m_pabyCachedTiles + 1 * nBandBlockSize + i, 0,
2310 800 : nBlockXSize - (iXOff + iXCount));
2311 800 : memset(m_pabyCachedTiles + 2 * nBandBlockSize + i, 0,
2312 800 : nBlockXSize - (iXOff + iXCount));
2313 800 : memset(m_pabyCachedTiles + 3 * nBandBlockSize + i, 0,
2314 800 : nBlockXSize - (iXOff + iXCount));
2315 : }
2316 : }
2317 7 : if (iYOff + iYCount < nBlockYSize)
2318 : {
2319 7 : const GPtrDiff_t i = (iYOff + iYCount) * nBlockXSize;
2320 7 : memset(m_pabyCachedTiles + 0 * nBandBlockSize + i, 0,
2321 7 : nBlockXSize * (nBlockYSize - (iYOff + iYCount)));
2322 7 : memset(m_pabyCachedTiles + 1 * nBandBlockSize + i, 0,
2323 7 : nBlockXSize * (nBlockYSize - (iYOff + iYCount)));
2324 7 : memset(m_pabyCachedTiles + 2 * nBandBlockSize + i, 0,
2325 7 : nBlockXSize * (nBlockYSize - (iYOff + iYCount)));
2326 7 : memset(m_pabyCachedTiles + 3 * nBandBlockSize + i, 0,
2327 7 : nBlockXSize * (nBlockYSize - (iYOff + iYCount)));
2328 : }
2329 : }
2330 :
2331 : char **papszDriverOptions =
2332 4163 : CSLSetNameValue(nullptr, "_INTERNAL_DATASET", "YES");
2333 4163 : if (EQUAL(pszDriverName, "JPEG") || EQUAL(pszDriverName, "WEBP"))
2334 : {
2335 : // If not all bands are dirty, then use lossless WEBP
2336 172 : if (!bAllDirty && EQUAL(pszDriverName, "WEBP"))
2337 : {
2338 2 : papszDriverOptions =
2339 2 : CSLSetNameValue(papszDriverOptions, "LOSSLESS", "YES");
2340 : }
2341 : else
2342 : {
2343 : papszDriverOptions =
2344 170 : CSLSetNameValue(papszDriverOptions, "QUALITY",
2345 : CPLSPrintf("%d", m_nQuality));
2346 : }
2347 : }
2348 3991 : else if (EQUAL(pszDriverName, "PNG"))
2349 : {
2350 3986 : papszDriverOptions = CSLSetNameValue(papszDriverOptions, "ZLEVEL",
2351 : CPLSPrintf("%d", m_nZLevel));
2352 : }
2353 5 : else if (EQUAL(pszDriverName, "GTiff"))
2354 : {
2355 : papszDriverOptions =
2356 5 : CSLSetNameValue(papszDriverOptions, "COMPRESS", "LZW");
2357 5 : if (nBlockXSize * nBlockYSize <= 512 * 512)
2358 : {
2359 : // If tile is not too big, create it as single-strip TIFF
2360 : papszDriverOptions =
2361 5 : CSLSetNameValue(papszDriverOptions, "BLOCKYSIZE",
2362 : CPLSPrintf("%d", nBlockYSize));
2363 : }
2364 : }
2365 : #ifdef DEBUG
2366 : VSIStatBufL sStat;
2367 4163 : CPLAssert(VSIStatL(osMemFileName, &sStat) != 0);
2368 : #endif
2369 : GDALDataset *poOutDS =
2370 4163 : l_poDriver->CreateCopy(osMemFileName, poMEMDS, FALSE,
2371 : papszDriverOptions, nullptr, nullptr);
2372 4163 : CSLDestroy(papszDriverOptions);
2373 4163 : CPLFree(pTempTileBuffer);
2374 :
2375 4163 : if (poOutDS)
2376 : {
2377 4163 : GDALClose(poOutDS);
2378 4163 : vsi_l_offset nBlobSize = 0;
2379 : GByte *pabyBlob =
2380 4163 : VSIGetMemFileBuffer(osMemFileName, &nBlobSize, TRUE);
2381 :
2382 : /* Create or commit and recreate transaction */
2383 4163 : GDALGPKGMBTilesLikePseudoDataset *poMainDS =
2384 4163 : m_poParentDS ? m_poParentDS : this;
2385 4163 : if (poMainDS->m_nTileInsertionCount == 0)
2386 : {
2387 189 : poMainDS->IStartTransaction();
2388 : }
2389 3974 : else if (poMainDS->m_nTileInsertionCount == 1000)
2390 : {
2391 3 : if (poMainDS->ICommitTransaction() != OGRERR_NONE)
2392 : {
2393 1 : poMainDS->m_nTileInsertionCount = -1;
2394 1 : CPLFree(pabyBlob);
2395 1 : VSIUnlink(osMemFileName);
2396 1 : delete poMEMDS;
2397 1 : return CE_Failure;
2398 : }
2399 2 : poMainDS->IStartTransaction();
2400 2 : poMainDS->m_nTileInsertionCount = 0;
2401 : }
2402 4162 : poMainDS->m_nTileInsertionCount++;
2403 :
2404 : char *pszSQL =
2405 4162 : sqlite3_mprintf("INSERT OR REPLACE INTO \"%w\" "
2406 : "(zoom_level, tile_row, tile_column, "
2407 : "tile_data) VALUES (%d, %d, %d, ?)",
2408 : m_osRasterTable.c_str(), m_nZoomLevel,
2409 4162 : GetRowFromIntoTopConvention(nRow), nCol);
2410 : #ifdef DEBUG_VERBOSE
2411 : CPLDebug("GPKG", "%s", pszSQL);
2412 : #endif
2413 4162 : sqlite3_stmt *hStmt = nullptr;
2414 4162 : int rc = SQLPrepareWithError(IGetDB(), pszSQL, -1, &hStmt, nullptr);
2415 4162 : if (rc != SQLITE_OK)
2416 : {
2417 0 : CPLFree(pabyBlob);
2418 : }
2419 : else
2420 : {
2421 4162 : sqlite3_bind_blob(hStmt, 1, pabyBlob,
2422 : static_cast<int>(nBlobSize), CPLFree);
2423 4162 : rc = sqlite3_step(hStmt);
2424 4162 : if (rc == SQLITE_DONE)
2425 4162 : eErr = CE_None;
2426 : else
2427 : {
2428 0 : CPLError(CE_Failure, CPLE_AppDefined,
2429 : "Failure when inserting tile (row=%d,col=%d) at "
2430 : "zoom_level=%d : %s",
2431 0 : GetRowFromIntoTopConvention(nRow), nCol,
2432 0 : m_nZoomLevel, sqlite3_errmsg(IGetDB()));
2433 : }
2434 : }
2435 4162 : sqlite3_finalize(hStmt);
2436 4162 : sqlite3_free(pszSQL);
2437 :
2438 4162 : if (m_eTF == GPKG_TF_PNG_16BIT || m_eTF == GPKG_TF_TIFF_32BIT_FLOAT)
2439 : {
2440 57 : GIntBig nTileId = GetTileId(nRow, nCol);
2441 57 : if (nTileId == 0)
2442 0 : eErr = CE_Failure;
2443 : else
2444 : {
2445 57 : DeleteFromGriddedTileAncillary(nTileId);
2446 :
2447 57 : pszSQL = sqlite3_mprintf(
2448 : "INSERT INTO gpkg_2d_gridded_tile_ancillary "
2449 : "(tpudt_name, tpudt_id, scale, offset, min, max, "
2450 : "mean, std_dev) VALUES "
2451 : "('%q', ?, %.17g, %.17g, ?, ?, ?, ?)",
2452 : m_osRasterTable.c_str(), dfTileScale, dfTileOffset);
2453 : #ifdef DEBUG_VERBOSE
2454 : CPLDebug("GPKG", "%s", pszSQL);
2455 : #endif
2456 57 : hStmt = nullptr;
2457 57 : rc = SQLPrepareWithError(IGetDB(), pszSQL, -1, &hStmt,
2458 : nullptr);
2459 57 : if (rc != SQLITE_OK)
2460 : {
2461 0 : eErr = CE_Failure;
2462 : }
2463 : else
2464 : {
2465 57 : sqlite3_bind_int64(hStmt, 1, nTileId);
2466 57 : sqlite3_bind_double(hStmt, 2, dfTileMin);
2467 57 : sqlite3_bind_double(hStmt, 3, dfTileMax);
2468 57 : sqlite3_bind_double(hStmt, 4, dfTileMean);
2469 57 : sqlite3_bind_double(hStmt, 5, dfTileStdDev);
2470 57 : rc = sqlite3_step(hStmt);
2471 57 : if (rc == SQLITE_DONE)
2472 : {
2473 57 : eErr = CE_None;
2474 : }
2475 : else
2476 : {
2477 0 : CPLError(CE_Failure, CPLE_AppDefined,
2478 : "Cannot insert into "
2479 : "gpkg_2d_gridded_tile_ancillary");
2480 0 : eErr = CE_Failure;
2481 : }
2482 : }
2483 57 : sqlite3_finalize(hStmt);
2484 57 : sqlite3_free(pszSQL);
2485 : }
2486 : }
2487 : }
2488 :
2489 4162 : VSIUnlink(osMemFileName);
2490 4162 : delete poMEMDS;
2491 : }
2492 : else
2493 : {
2494 1 : CPLError(CE_Failure, CPLE_NotSupported, "Cannot find driver %s",
2495 : pszDriverName);
2496 : }
2497 :
2498 4163 : return eErr;
2499 : }
2500 :
2501 : /************************************************************************/
2502 : /* FlushRemainingShiftedTiles() */
2503 : /************************************************************************/
2504 :
2505 : CPLErr
2506 609 : GDALGPKGMBTilesLikePseudoDataset::FlushRemainingShiftedTiles(bool bPartialFlush)
2507 : {
2508 609 : if (m_hTempDB == nullptr)
2509 460 : return CE_None;
2510 :
2511 745 : for (int i = 0; i <= 3; i++)
2512 : {
2513 596 : m_asCachedTilesDesc[i].nRow = -1;
2514 596 : m_asCachedTilesDesc[i].nCol = -1;
2515 596 : m_asCachedTilesDesc[i].nIdxWithinTileData = -1;
2516 : }
2517 :
2518 : int nBlockXSize, nBlockYSize;
2519 149 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
2520 149 : const int nBands = IGetRasterCount();
2521 149 : const int nRasterXSize = IGetRasterBand(1)->GetXSize();
2522 149 : const int nRasterYSize = IGetRasterBand(1)->GetYSize();
2523 149 : const int nXBlocks = DIV_ROUND_UP(nRasterXSize, nBlockXSize);
2524 149 : const int nYBlocks = DIV_ROUND_UP(nRasterYSize, nBlockYSize);
2525 :
2526 149 : int nPartialActiveTiles = 0;
2527 149 : if (bPartialFlush)
2528 : {
2529 2 : sqlite3_stmt *hStmt = nullptr;
2530 4 : CPLString osSQL;
2531 : osSQL.Printf("SELECT COUNT(*) FROM partial_tiles WHERE zoom_level = %d "
2532 : "AND partial_flag != 0",
2533 2 : m_nZoomLevel);
2534 2 : if (SQLPrepareWithError(m_hTempDB, osSQL.c_str(), -1, &hStmt,
2535 2 : nullptr) == SQLITE_OK)
2536 : {
2537 2 : if (sqlite3_step(hStmt) == SQLITE_ROW)
2538 : {
2539 2 : nPartialActiveTiles = sqlite3_column_int(hStmt, 0);
2540 2 : CPLDebug("GPKG", "Active partial tiles before flush: %d",
2541 : nPartialActiveTiles);
2542 : }
2543 2 : sqlite3_finalize(hStmt);
2544 : }
2545 : }
2546 :
2547 298 : CPLString osSQL = "SELECT tile_row, tile_column, partial_flag";
2548 616 : for (int nBand = 1; nBand <= nBands; nBand++)
2549 : {
2550 467 : osSQL += CPLSPrintf(", tile_data_band_%d", nBand);
2551 : }
2552 : osSQL += CPLSPrintf(" FROM partial_tiles WHERE "
2553 : "zoom_level = %d AND partial_flag != 0",
2554 149 : m_nZoomLevel);
2555 149 : if (bPartialFlush)
2556 : {
2557 2 : osSQL += " ORDER BY age";
2558 : }
2559 149 : const char *pszSQL = osSQL.c_str();
2560 :
2561 : #ifdef DEBUG_VERBOSE
2562 : CPLDebug("GPKG", "%s", pszSQL);
2563 : #endif
2564 149 : sqlite3_stmt *hStmt = nullptr;
2565 149 : int rc = SQLPrepareWithError(m_hTempDB, pszSQL, -1, &hStmt, nullptr);
2566 149 : if (rc != SQLITE_OK)
2567 : {
2568 0 : return CE_Failure;
2569 : }
2570 :
2571 149 : CPLErr eErr = CE_None;
2572 149 : bool bGotPartialTiles = false;
2573 149 : int nCountFlushedTiles = 0;
2574 149 : const size_t nBandBlockSize =
2575 149 : static_cast<size_t>(nBlockXSize) * nBlockYSize * m_nDTSize;
2576 175 : do
2577 : {
2578 324 : rc = sqlite3_step(hStmt);
2579 324 : if (rc == SQLITE_ROW)
2580 : {
2581 177 : bGotPartialTiles = true;
2582 :
2583 177 : int nRow = sqlite3_column_int(hStmt, 0);
2584 177 : int nCol = sqlite3_column_int(hStmt, 1);
2585 177 : int nPartialFlags = sqlite3_column_int(hStmt, 2);
2586 :
2587 177 : if (bPartialFlush)
2588 : {
2589 : // This method assumes that there are no dirty blocks alive
2590 : // so check this assumption.
2591 : // When called with bPartialFlush = false, FlushCache() has
2592 : // already been called, so no need to check.
2593 16 : bool bFoundDirtyBlock = false;
2594 16 : int nBlockXOff = nCol - m_nShiftXTiles;
2595 16 : int nBlockYOff = nRow - m_nShiftYTiles;
2596 96 : for (int iX = 0; !bFoundDirtyBlock &&
2597 48 : iX < ((m_nShiftXPixelsMod != 0) ? 2 : 1);
2598 : iX++)
2599 : {
2600 32 : if (nBlockXOff + iX < 0 || nBlockXOff + iX >= nXBlocks)
2601 12 : continue;
2602 120 : for (int iY = 0; !bFoundDirtyBlock &&
2603 60 : iY < ((m_nShiftYPixelsMod != 0) ? 2 : 1);
2604 : iY++)
2605 : {
2606 40 : if (nBlockYOff + iY < 0 || nBlockYOff + iY >= nYBlocks)
2607 17 : continue;
2608 69 : for (int iBand = 1;
2609 69 : !bFoundDirtyBlock && iBand <= nBands; iBand++)
2610 : {
2611 : GDALRasterBlock *poBlock =
2612 : cpl::down_cast<GDALGPKGMBTilesLikeRasterBand *>(
2613 46 : IGetRasterBand(iBand))
2614 46 : ->AccessibleTryGetLockedBlockRef(
2615 : nBlockXOff + iX, nBlockYOff + iY);
2616 46 : if (poBlock)
2617 : {
2618 0 : if (poBlock->GetDirty())
2619 0 : bFoundDirtyBlock = true;
2620 0 : poBlock->DropLock();
2621 : }
2622 : }
2623 : }
2624 : }
2625 16 : if (bFoundDirtyBlock)
2626 : {
2627 : #ifdef DEBUG_VERBOSE
2628 : CPLDebug("GPKG",
2629 : "Skipped flushing tile row = %d, column = %d "
2630 : "because it has dirty block(s) in GDAL cache",
2631 : nRow, nCol);
2632 : #endif
2633 0 : continue;
2634 : }
2635 : }
2636 :
2637 177 : nCountFlushedTiles++;
2638 177 : if (bPartialFlush && nCountFlushedTiles >= nPartialActiveTiles / 2)
2639 : {
2640 2 : CPLDebug("GPKG", "Flushed %d tiles", nCountFlushedTiles);
2641 2 : break;
2642 : }
2643 :
2644 705 : for (int nBand = 1; nBand <= nBands; nBand++)
2645 : {
2646 530 : if (nPartialFlags & (((1 << 4) - 1) << (4 * (nBand - 1))))
2647 : {
2648 498 : CPLAssert(sqlite3_column_bytes(hStmt, 2 + nBand) ==
2649 : static_cast<int>(nBandBlockSize));
2650 498 : memcpy(m_pabyCachedTiles + (nBand - 1) * nBandBlockSize,
2651 : sqlite3_column_blob(hStmt, 2 + nBand),
2652 : nBandBlockSize);
2653 : }
2654 : else
2655 : {
2656 32 : FillEmptyTileSingleBand(m_pabyCachedTiles +
2657 32 : (nBand - 1) * nBandBlockSize);
2658 : }
2659 : }
2660 :
2661 175 : int nFullFlags = (1 << (4 * nBands)) - 1;
2662 :
2663 : // In case the partial flags indicate that there's some quadrant
2664 : // missing, check in the main database if there is already a tile
2665 : // In which case, use the parts of that tile that aren't in the
2666 : // temporary database
2667 175 : if (nPartialFlags != nFullFlags)
2668 : {
2669 525 : char *pszNewSQL = sqlite3_mprintf(
2670 : "SELECT tile_data%s FROM \"%w\" "
2671 : "WHERE zoom_level = %d AND tile_row = %d AND tile_column = "
2672 : "%d%s",
2673 175 : m_eDT != GDT_Byte ? ", id"
2674 : : "", // MBTiles do not have an id
2675 : m_osRasterTable.c_str(), m_nZoomLevel,
2676 175 : GetRowFromIntoTopConvention(nRow), nCol,
2677 175 : !m_osWHERE.empty()
2678 0 : ? CPLSPrintf(" AND (%s)", m_osWHERE.c_str())
2679 : : "");
2680 : #ifdef DEBUG_VERBOSE
2681 : CPLDebug("GPKG", "%s", pszNewSQL);
2682 : #endif
2683 175 : sqlite3_stmt *hNewStmt = nullptr;
2684 175 : rc = SQLPrepareWithError(IGetDB(), pszNewSQL, -1, &hNewStmt,
2685 : nullptr);
2686 175 : if (rc == SQLITE_OK)
2687 : {
2688 175 : rc = sqlite3_step(hNewStmt);
2689 192 : if (rc == SQLITE_ROW &&
2690 17 : sqlite3_column_type(hNewStmt, 0) == SQLITE_BLOB)
2691 : {
2692 17 : const int nBytes = sqlite3_column_bytes(hNewStmt, 0);
2693 : GIntBig nTileId =
2694 17 : (m_eDT == GDT_Byte)
2695 17 : ? 0
2696 0 : : sqlite3_column_int64(hNewStmt, 1);
2697 : GByte *pabyRawData =
2698 : const_cast<GByte *>(static_cast<const GByte *>(
2699 17 : sqlite3_column_blob(hNewStmt, 0)));
2700 : const CPLString osMemFileName(
2701 34 : VSIMemGenerateHiddenFilename("gpkg_read_tile"));
2702 17 : VSILFILE *fp = VSIFileFromMemBuffer(
2703 : osMemFileName.c_str(), pabyRawData, nBytes, FALSE);
2704 17 : VSIFCloseL(fp);
2705 :
2706 17 : double dfTileOffset = 0.0;
2707 17 : double dfTileScale = 1.0;
2708 17 : GetTileOffsetAndScale(nTileId, dfTileOffset,
2709 : dfTileScale);
2710 17 : const int nTileBands = m_eDT == GDT_Byte ? 4 : 1;
2711 17 : GByte *pabyTemp =
2712 17 : m_pabyCachedTiles + nTileBands * nBandBlockSize;
2713 17 : ReadTile(osMemFileName, pabyTemp, dfTileOffset,
2714 : dfTileScale);
2715 17 : VSIUnlink(osMemFileName);
2716 :
2717 17 : int iYQuadrantMax = (m_nShiftYPixelsMod) ? 1 : 0;
2718 17 : int iXQuadrantMax = (m_nShiftXPixelsMod) ? 1 : 0;
2719 51 : for (int iYQuadrant = 0; iYQuadrant <= iYQuadrantMax;
2720 : iYQuadrant++)
2721 : {
2722 100 : for (int iXQuadrant = 0;
2723 100 : iXQuadrant <= iXQuadrantMax; iXQuadrant++)
2724 : {
2725 226 : for (int nBand = 1; nBand <= nBands; nBand++)
2726 : {
2727 160 : int iQuadrantFlag = 0;
2728 160 : if (iXQuadrant == 0 && iYQuadrant == 0)
2729 42 : iQuadrantFlag |= (1 << 0);
2730 160 : if (iXQuadrant == iXQuadrantMax &&
2731 : iYQuadrant == 0)
2732 42 : iQuadrantFlag |= (1 << 1);
2733 160 : if (iXQuadrant == 0 &&
2734 : iYQuadrant == iYQuadrantMax)
2735 42 : iQuadrantFlag |= (1 << 2);
2736 160 : if (iXQuadrant == iXQuadrantMax &&
2737 : iYQuadrant == iYQuadrantMax)
2738 42 : iQuadrantFlag |= (1 << 3);
2739 160 : int nLocalFlag = iQuadrantFlag
2740 160 : << (4 * (nBand - 1));
2741 160 : if (!(nPartialFlags & nLocalFlag))
2742 : {
2743 : int nXOff, nYOff, nXSize, nYSize;
2744 118 : if (iXQuadrant == 0 &&
2745 60 : m_nShiftXPixelsMod != 0)
2746 : {
2747 56 : nXOff = 0;
2748 56 : nXSize = m_nShiftXPixelsMod;
2749 : }
2750 : else
2751 : {
2752 62 : nXOff = m_nShiftXPixelsMod;
2753 62 : nXSize = nBlockXSize -
2754 62 : m_nShiftXPixelsMod;
2755 : }
2756 118 : if (iYQuadrant == 0 &&
2757 63 : m_nShiftYPixelsMod != 0)
2758 : {
2759 63 : nYOff = 0;
2760 63 : nYSize = m_nShiftYPixelsMod;
2761 : }
2762 : else
2763 : {
2764 55 : nYOff = m_nShiftYPixelsMod;
2765 55 : nYSize = nBlockYSize -
2766 55 : m_nShiftYPixelsMod;
2767 : }
2768 118 : for (int iY = nYOff;
2769 13888 : iY < nYOff + nYSize; iY++)
2770 : {
2771 13770 : memcpy(m_pabyCachedTiles +
2772 13770 : ((static_cast<size_t>(
2773 13770 : nBand - 1) *
2774 13770 : nBlockYSize +
2775 13770 : iY) *
2776 13770 : nBlockXSize +
2777 13770 : nXOff) *
2778 13770 : m_nDTSize,
2779 13770 : pabyTemp +
2780 13770 : ((static_cast<size_t>(
2781 13770 : nBand - 1) *
2782 13770 : nBlockYSize +
2783 13770 : iY) *
2784 13770 : nBlockXSize +
2785 13770 : nXOff) *
2786 13770 : m_nDTSize,
2787 13770 : static_cast<size_t>(nXSize) *
2788 13770 : m_nDTSize);
2789 : }
2790 : }
2791 : }
2792 : }
2793 : }
2794 : }
2795 158 : else if (rc != SQLITE_DONE)
2796 : {
2797 0 : CPLError(CE_Failure, CPLE_AppDefined,
2798 : "sqlite3_step(%s) failed: %s", pszNewSQL,
2799 : sqlite3_errmsg(m_hTempDB));
2800 : }
2801 175 : sqlite3_finalize(hNewStmt);
2802 : }
2803 175 : sqlite3_free(pszNewSQL);
2804 : }
2805 :
2806 175 : m_asCachedTilesDesc[0].nRow = nRow;
2807 175 : m_asCachedTilesDesc[0].nCol = nCol;
2808 175 : m_asCachedTilesDesc[0].nIdxWithinTileData = 0;
2809 175 : m_asCachedTilesDesc[0].abBandDirty[0] = true;
2810 175 : m_asCachedTilesDesc[0].abBandDirty[1] = true;
2811 175 : m_asCachedTilesDesc[0].abBandDirty[2] = true;
2812 175 : m_asCachedTilesDesc[0].abBandDirty[3] = true;
2813 :
2814 175 : eErr = WriteTile();
2815 :
2816 175 : if (eErr == CE_None && bPartialFlush)
2817 : {
2818 : pszSQL =
2819 14 : CPLSPrintf("DELETE FROM partial_tiles WHERE zoom_level = "
2820 : "%d AND tile_row = %d AND tile_column = %d",
2821 : m_nZoomLevel, nRow, nCol);
2822 : #ifdef DEBUG_VERBOSE
2823 : CPLDebug("GPKG", "%s", pszSQL);
2824 : #endif
2825 14 : if (SQLCommand(m_hTempDB, pszSQL) != OGRERR_NONE)
2826 0 : eErr = CE_None;
2827 : }
2828 : }
2829 : else
2830 : {
2831 147 : if (rc != SQLITE_DONE)
2832 : {
2833 0 : CPLError(CE_Failure, CPLE_AppDefined,
2834 : "sqlite3_step(%s) failed: %s", pszSQL,
2835 : sqlite3_errmsg(m_hTempDB));
2836 : }
2837 147 : break;
2838 : }
2839 175 : } while (eErr == CE_None);
2840 :
2841 149 : sqlite3_finalize(hStmt);
2842 :
2843 149 : if (bPartialFlush && nCountFlushedTiles < nPartialActiveTiles / 2)
2844 : {
2845 0 : CPLDebug("GPKG", "Flushed %d tiles. Target was %d", nCountFlushedTiles,
2846 : nPartialActiveTiles / 2);
2847 : }
2848 :
2849 149 : if (bGotPartialTiles && !bPartialFlush)
2850 : {
2851 : #ifdef DEBUG_VERBOSE
2852 : pszSQL = CPLSPrintf("SELECT p1.id, p1.tile_row, p1.tile_column FROM "
2853 : "partial_tiles p1, partial_tiles p2 "
2854 : "WHERE p1.zoom_level = %d AND p2.zoom_level = %d "
2855 : "AND p1.tile_row = p2.tile_row AND p1.tile_column "
2856 : "= p2.tile_column AND p2.partial_flag != 0",
2857 : -1 - m_nZoomLevel, m_nZoomLevel);
2858 : rc = SQLPrepareWithError(m_hTempDB, pszSQL, -1, &hStmt, nullptr);
2859 : CPLAssert(rc == SQLITE_OK);
2860 : while ((rc = sqlite3_step(hStmt)) == SQLITE_ROW)
2861 : {
2862 : CPLDebug("GPKG",
2863 : "Conflict: existing id = %d, tile_row = %d, tile_column = "
2864 : "%d, zoom_level = %d",
2865 : sqlite3_column_int(hStmt, 0), sqlite3_column_int(hStmt, 1),
2866 : sqlite3_column_int(hStmt, 2), m_nZoomLevel);
2867 : }
2868 : sqlite3_finalize(hStmt);
2869 : #endif
2870 :
2871 110 : pszSQL = CPLSPrintf("UPDATE partial_tiles SET zoom_level = %d, "
2872 : "partial_flag = 0, age = -1 WHERE zoom_level = %d "
2873 : "AND partial_flag != 0",
2874 55 : -1 - m_nZoomLevel, m_nZoomLevel);
2875 : #ifdef DEBUG_VERBOSE
2876 : CPLDebug("GPKG", "%s", pszSQL);
2877 : #endif
2878 55 : SQLCommand(m_hTempDB, pszSQL);
2879 : }
2880 :
2881 149 : return eErr;
2882 : }
2883 :
2884 : /************************************************************************/
2885 : /* DoPartialFlushOfPartialTilesIfNecessary() */
2886 : /************************************************************************/
2887 :
2888 : CPLErr
2889 4556 : GDALGPKGMBTilesLikePseudoDataset::DoPartialFlushOfPartialTilesIfNecessary()
2890 : {
2891 4556 : time_t nCurTimeStamp = time(nullptr);
2892 4556 : if (m_nLastSpaceCheckTimestamp == 0)
2893 3 : m_nLastSpaceCheckTimestamp = nCurTimeStamp;
2894 4556 : if (m_nLastSpaceCheckTimestamp > 0 &&
2895 81 : (m_bForceTempDBCompaction ||
2896 13 : nCurTimeStamp - m_nLastSpaceCheckTimestamp > 10))
2897 : {
2898 68 : m_nLastSpaceCheckTimestamp = nCurTimeStamp;
2899 : GIntBig nFreeSpace =
2900 68 : VSIGetDiskFreeSpace(CPLGetDirnameSafe(m_osTempDBFilename).c_str());
2901 68 : bool bTryFreeing = false;
2902 68 : if (nFreeSpace >= 0 && nFreeSpace < 1024 * 1024 * 1024)
2903 : {
2904 0 : CPLDebug("GPKG",
2905 : "Free space below 1GB. Flushing part of partial tiles");
2906 0 : bTryFreeing = true;
2907 : }
2908 : else
2909 : {
2910 : VSIStatBufL sStat;
2911 68 : if (VSIStatL(m_osTempDBFilename, &sStat) == 0)
2912 : {
2913 68 : GIntBig nTempSpace = sStat.st_size;
2914 136 : if (VSIStatL((m_osTempDBFilename + "-journal").c_str(),
2915 68 : &sStat) == 0)
2916 0 : nTempSpace += sStat.st_size;
2917 136 : else if (VSIStatL((m_osTempDBFilename + "-wal").c_str(),
2918 68 : &sStat) == 0)
2919 0 : nTempSpace += sStat.st_size;
2920 :
2921 : int nBlockXSize, nBlockYSize;
2922 68 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
2923 68 : const int nBands = IGetRasterCount();
2924 :
2925 68 : if (nTempSpace >
2926 68 : 4 * static_cast<GIntBig>(IGetRasterBand(1)->GetXSize()) *
2927 68 : nBlockYSize * nBands * m_nDTSize)
2928 : {
2929 2 : CPLDebug("GPKG",
2930 : "Partial tiles DB is " CPL_FRMT_GIB
2931 : " bytes. Flushing part of partial tiles",
2932 : nTempSpace);
2933 2 : bTryFreeing = true;
2934 : }
2935 : }
2936 : }
2937 68 : if (bTryFreeing)
2938 : {
2939 2 : if (FlushRemainingShiftedTiles(true /* partial flush*/) != CE_None)
2940 : {
2941 0 : return CE_Failure;
2942 : }
2943 2 : SQLCommand(m_hTempDB,
2944 : "DELETE FROM partial_tiles WHERE zoom_level < 0");
2945 2 : SQLCommand(m_hTempDB, "VACUUM");
2946 : }
2947 : }
2948 4556 : return CE_None;
2949 : }
2950 :
2951 : /************************************************************************/
2952 : /* WriteShiftedTile() */
2953 : /************************************************************************/
2954 :
2955 4556 : CPLErr GDALGPKGMBTilesLikePseudoDataset::WriteShiftedTile(
2956 : int nRow, int nCol, int nBand, int nDstXOffset, int nDstYOffset,
2957 : int nDstXSize, int nDstYSize)
2958 : {
2959 4556 : CPLAssert(m_nShiftXPixelsMod || m_nShiftYPixelsMod);
2960 4556 : CPLAssert(nRow >= 0);
2961 4556 : CPLAssert(nCol >= 0);
2962 4556 : CPLAssert(nRow < m_nTileMatrixHeight);
2963 4556 : CPLAssert(nCol < m_nTileMatrixWidth);
2964 :
2965 4556 : if (m_hTempDB == nullptr &&
2966 55 : (m_poParentDS == nullptr || m_poParentDS->m_hTempDB == nullptr))
2967 : {
2968 : const char *pszBaseFilename =
2969 48 : m_poParentDS ? m_poParentDS->IGetFilename() : IGetFilename();
2970 : m_osTempDBFilename =
2971 48 : CPLResetExtensionSafe(pszBaseFilename, "partial_tiles.db");
2972 48 : CPLPushErrorHandler(CPLQuietErrorHandler);
2973 48 : VSIUnlink(m_osTempDBFilename);
2974 48 : CPLPopErrorHandler();
2975 48 : m_hTempDB = nullptr;
2976 48 : int rc = 0;
2977 48 : if (STARTS_WITH(m_osTempDBFilename, "/vsi"))
2978 : {
2979 33 : m_pMyVFS = OGRSQLiteCreateVFS(nullptr, nullptr);
2980 33 : sqlite3_vfs_register(m_pMyVFS, 0);
2981 33 : rc = sqlite3_open_v2(m_osTempDBFilename, &m_hTempDB,
2982 : SQLITE_OPEN_READWRITE | SQLITE_OPEN_CREATE |
2983 : SQLITE_OPEN_NOMUTEX,
2984 33 : m_pMyVFS->zName);
2985 : }
2986 : else
2987 : {
2988 15 : rc = sqlite3_open(m_osTempDBFilename, &m_hTempDB);
2989 : }
2990 48 : if (rc != SQLITE_OK || m_hTempDB == nullptr)
2991 : {
2992 0 : CPLError(CE_Failure, CPLE_AppDefined,
2993 : "Cannot create temporary database %s",
2994 : m_osTempDBFilename.c_str());
2995 0 : return CE_Failure;
2996 : }
2997 48 : SQLCommand(m_hTempDB, "PRAGMA synchronous = OFF");
2998 48 : SQLCommand(m_hTempDB, "PRAGMA journal_mode = OFF");
2999 48 : SQLCommand(m_hTempDB, "CREATE TABLE partial_tiles("
3000 : "id INTEGER PRIMARY KEY AUTOINCREMENT,"
3001 : "zoom_level INTEGER NOT NULL,"
3002 : "tile_column INTEGER NOT NULL,"
3003 : "tile_row INTEGER NOT NULL,"
3004 : "tile_data_band_1 BLOB,"
3005 : "tile_data_band_2 BLOB,"
3006 : "tile_data_band_3 BLOB,"
3007 : "tile_data_band_4 BLOB,"
3008 : "partial_flag INTEGER NOT NULL,"
3009 : "age INTEGER NOT NULL,"
3010 : "UNIQUE (zoom_level, tile_column, tile_row))");
3011 48 : SQLCommand(m_hTempDB, "CREATE INDEX partial_tiles_partial_flag_idx "
3012 : "ON partial_tiles(partial_flag)");
3013 48 : SQLCommand(m_hTempDB, "CREATE INDEX partial_tiles_age_idx "
3014 : "ON partial_tiles(age)");
3015 :
3016 48 : if (m_poParentDS != nullptr)
3017 : {
3018 4 : m_poParentDS->m_osTempDBFilename = m_osTempDBFilename;
3019 4 : m_poParentDS->m_hTempDB = m_hTempDB;
3020 : }
3021 : }
3022 :
3023 4556 : if (m_poParentDS != nullptr)
3024 4020 : m_hTempDB = m_poParentDS->m_hTempDB;
3025 :
3026 : int nBlockXSize, nBlockYSize;
3027 4556 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
3028 4556 : const int nBands = IGetRasterCount();
3029 4556 : const size_t nBandBlockSize =
3030 4556 : static_cast<size_t>(nBlockXSize) * nBlockYSize * m_nDTSize;
3031 :
3032 4556 : int iQuadrantFlag = 0;
3033 4556 : if (nDstXOffset == 0 && nDstYOffset == 0)
3034 1031 : iQuadrantFlag |= (1 << 0);
3035 4556 : if ((nDstXOffset != 0 || nDstXOffset + nDstXSize == nBlockXSize) &&
3036 : nDstYOffset == 0)
3037 1125 : iQuadrantFlag |= (1 << 1);
3038 4556 : if (nDstXOffset == 0 &&
3039 1031 : (nDstYOffset != 0 || nDstYOffset + nDstYSize == nBlockYSize))
3040 1135 : iQuadrantFlag |= (1 << 2);
3041 4556 : if ((nDstXOffset != 0 || nDstXOffset + nDstXSize == nBlockXSize) &&
3042 1125 : (nDstYOffset != 0 || nDstYOffset + nDstYSize == nBlockYSize))
3043 1323 : iQuadrantFlag |= (1 << 3);
3044 4556 : int l_nFlags = iQuadrantFlag << (4 * (nBand - 1));
3045 4556 : int nFullFlags = (1 << (4 * nBands)) - 1;
3046 4556 : int nOldFlags = 0;
3047 :
3048 18224 : for (int i = 1; i <= 3; i++)
3049 : {
3050 13668 : m_asCachedTilesDesc[i].nRow = -1;
3051 13668 : m_asCachedTilesDesc[i].nCol = -1;
3052 13668 : m_asCachedTilesDesc[i].nIdxWithinTileData = -1;
3053 : }
3054 :
3055 4556 : int nExistingId = 0;
3056 4556 : const char *pszSQL = CPLSPrintf(
3057 : "SELECT id, partial_flag, tile_data_band_%d FROM partial_tiles WHERE "
3058 : "zoom_level = %d AND tile_row = %d AND tile_column = %d",
3059 : nBand, m_nZoomLevel, nRow, nCol);
3060 : #ifdef DEBUG_VERBOSE
3061 : CPLDebug("GPKG", "%s", pszSQL);
3062 : #endif
3063 4556 : sqlite3_stmt *hStmt = nullptr;
3064 4556 : int rc = SQLPrepareWithError(m_hTempDB, pszSQL, -1, &hStmt, nullptr);
3065 4556 : if (rc != SQLITE_OK)
3066 : {
3067 0 : return CE_Failure;
3068 : }
3069 :
3070 4556 : rc = sqlite3_step(hStmt);
3071 4556 : const int nTileBands = m_eDT == GDT_Byte ? 4 : 1;
3072 4556 : GByte *pabyTemp = m_pabyCachedTiles + nTileBands * nBandBlockSize;
3073 4556 : if (rc == SQLITE_ROW)
3074 : {
3075 4134 : nExistingId = sqlite3_column_int(hStmt, 0);
3076 : #ifdef DEBUG_VERBOSE
3077 : CPLDebug("GPKG", "Using partial_tile id=%d", nExistingId);
3078 : #endif
3079 4134 : nOldFlags = sqlite3_column_int(hStmt, 1);
3080 4134 : CPLAssert(nOldFlags != 0);
3081 4134 : if ((nOldFlags & (((1 << 4) - 1) << (4 * (nBand - 1)))) == 0)
3082 : {
3083 1021 : FillEmptyTileSingleBand(pabyTemp + (nBand - 1) * nBandBlockSize);
3084 : }
3085 : else
3086 : {
3087 3113 : CPLAssert(sqlite3_column_bytes(hStmt, 2) ==
3088 : static_cast<int>(nBandBlockSize));
3089 3113 : memcpy(pabyTemp + (nBand - 1) * nBandBlockSize,
3090 : sqlite3_column_blob(hStmt, 2), nBandBlockSize);
3091 : }
3092 : }
3093 : else
3094 : {
3095 422 : FillEmptyTileSingleBand(pabyTemp + (nBand - 1) * nBandBlockSize);
3096 : }
3097 4556 : sqlite3_finalize(hStmt);
3098 4556 : hStmt = nullptr;
3099 :
3100 : /* Copy the updated rectangle into the full tile */
3101 556108 : for (int iY = nDstYOffset; iY < nDstYOffset + nDstYSize; iY++)
3102 : {
3103 551552 : memcpy(pabyTemp +
3104 551552 : (static_cast<size_t>(nBand - 1) * nBlockXSize * nBlockYSize +
3105 551552 : static_cast<size_t>(iY) * nBlockXSize + nDstXOffset) *
3106 551552 : m_nDTSize,
3107 551552 : m_pabyCachedTiles +
3108 551552 : (static_cast<size_t>(nBand - 1) * nBlockXSize * nBlockYSize +
3109 551552 : static_cast<size_t>(iY) * nBlockXSize + nDstXOffset) *
3110 551552 : m_nDTSize,
3111 551552 : static_cast<size_t>(nDstXSize) * m_nDTSize);
3112 : }
3113 :
3114 : #ifdef notdef
3115 : static int nCounter = 1;
3116 : GDALDataset *poLogDS =
3117 : ((GDALDriver *)GDALGetDriverByName("GTiff"))
3118 : ->Create(CPLSPrintf("/tmp/partial_band_%d_%d.tif", 1, nCounter++),
3119 : nBlockXSize, nBlockYSize, nBands, m_eDT, NULL);
3120 : poLogDS->RasterIO(GF_Write, 0, 0, nBlockXSize, nBlockYSize,
3121 : pabyTemp + (nBand - 1) * nBandBlockSize, nBlockXSize,
3122 : nBlockYSize, m_eDT, 1, NULL, 0, 0, 0, NULL);
3123 : GDALClose(poLogDS);
3124 : #endif
3125 :
3126 4556 : if ((nOldFlags & l_nFlags) != 0)
3127 : {
3128 0 : CPLDebug("GPKG",
3129 : "Rewriting quadrant %d of band %d of tile (row=%d,col=%d)",
3130 : iQuadrantFlag, nBand, nRow, nCol);
3131 : }
3132 :
3133 4556 : l_nFlags |= nOldFlags;
3134 4556 : if (l_nFlags == nFullFlags)
3135 : {
3136 : #ifdef DEBUG_VERBOSE
3137 : CPLDebug("GPKG", "Got all quadrants for that tile");
3138 : #endif
3139 1192 : for (int iBand = 1; iBand <= nBands; iBand++)
3140 : {
3141 945 : if (iBand != nBand)
3142 : {
3143 698 : pszSQL = CPLSPrintf(
3144 : "SELECT tile_data_band_%d FROM partial_tiles WHERE "
3145 : "id = %d",
3146 : iBand, nExistingId);
3147 : #ifdef DEBUG_VERBOSE
3148 : CPLDebug("GPKG", "%s", pszSQL);
3149 : #endif
3150 698 : hStmt = nullptr;
3151 : rc =
3152 698 : SQLPrepareWithError(m_hTempDB, pszSQL, -1, &hStmt, nullptr);
3153 698 : if (rc != SQLITE_OK)
3154 : {
3155 0 : return CE_Failure;
3156 : }
3157 :
3158 698 : rc = sqlite3_step(hStmt);
3159 698 : if (rc == SQLITE_ROW)
3160 : {
3161 698 : CPLAssert(sqlite3_column_bytes(hStmt, 0) ==
3162 : static_cast<int>(nBandBlockSize));
3163 698 : memcpy(m_pabyCachedTiles + (iBand - 1) * nBandBlockSize,
3164 : sqlite3_column_blob(hStmt, 0), nBandBlockSize);
3165 : }
3166 698 : sqlite3_finalize(hStmt);
3167 698 : hStmt = nullptr;
3168 : }
3169 : else
3170 : {
3171 247 : memcpy(m_pabyCachedTiles + (iBand - 1) * nBandBlockSize,
3172 247 : pabyTemp + (iBand - 1) * nBandBlockSize, nBandBlockSize);
3173 : }
3174 : }
3175 :
3176 247 : m_asCachedTilesDesc[0].nRow = nRow;
3177 247 : m_asCachedTilesDesc[0].nCol = nCol;
3178 247 : m_asCachedTilesDesc[0].nIdxWithinTileData = 0;
3179 247 : m_asCachedTilesDesc[0].abBandDirty[0] = true;
3180 247 : m_asCachedTilesDesc[0].abBandDirty[1] = true;
3181 247 : m_asCachedTilesDesc[0].abBandDirty[2] = true;
3182 247 : m_asCachedTilesDesc[0].abBandDirty[3] = true;
3183 :
3184 494 : pszSQL = CPLSPrintf("UPDATE partial_tiles SET zoom_level = %d, "
3185 : "partial_flag = 0, age = -1 WHERE id = %d",
3186 247 : -1 - m_nZoomLevel, nExistingId);
3187 247 : SQLCommand(m_hTempDB, pszSQL);
3188 : #ifdef DEBUG_VERBOSE
3189 : CPLDebug("GPKG", "%s", pszSQL);
3190 : #endif
3191 247 : CPLErr eErr = WriteTile();
3192 :
3193 : // Call DoPartialFlushOfPartialTilesIfNecessary() after using
3194 : // m_pabyCachedTiles as it is going to mess with it.
3195 247 : if (DoPartialFlushOfPartialTilesIfNecessary() != CE_None)
3196 0 : eErr = CE_None;
3197 :
3198 247 : return eErr;
3199 : }
3200 :
3201 4309 : if (nExistingId == 0)
3202 : {
3203 : OGRErr err;
3204 844 : pszSQL = CPLSPrintf("SELECT id FROM partial_tiles WHERE "
3205 : "partial_flag = 0 AND zoom_level = %d "
3206 : "AND tile_row = %d AND tile_column = %d",
3207 422 : -1 - m_nZoomLevel, nRow, nCol);
3208 : #ifdef DEBUG_VERBOSE
3209 : CPLDebug("GPKG", "%s", pszSQL);
3210 : #endif
3211 422 : nExistingId = SQLGetInteger(m_hTempDB, pszSQL, &err);
3212 422 : if (nExistingId == 0)
3213 : {
3214 418 : pszSQL =
3215 : "SELECT id FROM partial_tiles WHERE partial_flag = 0 LIMIT 1";
3216 : #ifdef DEBUG_VERBOSE
3217 : CPLDebug("GPKG", "%s", pszSQL);
3218 : #endif
3219 418 : nExistingId = SQLGetInteger(m_hTempDB, pszSQL, &err);
3220 : }
3221 : }
3222 :
3223 4309 : const GIntBig nAge = (m_poParentDS) ? m_poParentDS->m_nAge : m_nAge;
3224 4309 : if (nExistingId == 0)
3225 : {
3226 338 : pszSQL = CPLSPrintf(
3227 : "INSERT INTO partial_tiles "
3228 : "(zoom_level, tile_row, tile_column, tile_data_band_%d, "
3229 : "partial_flag, age) VALUES (%d, %d, %d, ?, %d, " CPL_FRMT_GIB ")",
3230 : nBand, m_nZoomLevel, nRow, nCol, l_nFlags, nAge);
3231 : }
3232 : else
3233 : {
3234 3971 : pszSQL = CPLSPrintf(
3235 : "UPDATE partial_tiles SET zoom_level = %d, "
3236 : "tile_row = %d, tile_column = %d, "
3237 : "tile_data_band_%d = ?, partial_flag = %d, age = " CPL_FRMT_GIB
3238 : " WHERE id = %d",
3239 : m_nZoomLevel, nRow, nCol, nBand, l_nFlags, nAge, nExistingId);
3240 : }
3241 4309 : if (m_poParentDS)
3242 3801 : m_poParentDS->m_nAge++;
3243 : else
3244 508 : m_nAge++;
3245 :
3246 : #ifdef DEBUG_VERBOSE
3247 : CPLDebug("GPKG", "%s", pszSQL);
3248 : #endif
3249 :
3250 4309 : hStmt = nullptr;
3251 4309 : rc = SQLPrepareWithError(m_hTempDB, pszSQL, -1, &hStmt, nullptr);
3252 4309 : if (rc != SQLITE_OK)
3253 : {
3254 0 : return CE_Failure;
3255 : }
3256 :
3257 4309 : sqlite3_bind_blob(hStmt, 1, pabyTemp + (nBand - 1) * nBandBlockSize,
3258 : static_cast<int>(nBandBlockSize), SQLITE_TRANSIENT);
3259 4309 : rc = sqlite3_step(hStmt);
3260 4309 : CPLErr eErr = CE_Failure;
3261 4309 : if (rc == SQLITE_DONE)
3262 4309 : eErr = CE_None;
3263 : else
3264 : {
3265 0 : CPLError(CE_Failure, CPLE_AppDefined,
3266 : "Failure when inserting partial tile (row=%d,col=%d) at "
3267 : "zoom_level=%d : %s",
3268 : nRow, nCol, m_nZoomLevel, sqlite3_errmsg(m_hTempDB));
3269 : }
3270 :
3271 4309 : sqlite3_finalize(hStmt);
3272 :
3273 : // Call DoPartialFlushOfPartialTilesIfNecessary() after using
3274 : // m_pabyCachedTiles as it is going to mess with it.
3275 4309 : if (DoPartialFlushOfPartialTilesIfNecessary() != CE_None)
3276 0 : eErr = CE_None;
3277 :
3278 4309 : return eErr;
3279 : }
3280 :
3281 : /************************************************************************/
3282 : /* IWriteBlock() */
3283 : /************************************************************************/
3284 :
3285 5849 : CPLErr GDALGPKGMBTilesLikeRasterBand::IWriteBlock(int nBlockXOff,
3286 : int nBlockYOff, void *pData)
3287 : {
3288 : #ifdef DEBUG_VERBOSE
3289 : CPLDebug(
3290 : "GPKG",
3291 : "IWriteBlock(nBand=%d,nBlockXOff=%d,nBlockYOff=%d,m_nZoomLevel=%d)",
3292 : nBand, nBlockXOff, nBlockYOff, m_poTPD->m_nZoomLevel);
3293 : #endif
3294 5849 : if (!m_poTPD->ICanIWriteBlock())
3295 : {
3296 0 : return CE_Failure;
3297 : }
3298 5849 : if (m_poTPD->m_poParentDS)
3299 1139 : m_poTPD->m_poParentDS->m_bHasModifiedTiles = true;
3300 : else
3301 4710 : m_poTPD->m_bHasModifiedTiles = true;
3302 :
3303 5849 : int nRow = nBlockYOff + m_poTPD->m_nShiftYTiles;
3304 5849 : int nCol = nBlockXOff + m_poTPD->m_nShiftXTiles;
3305 :
3306 5849 : int nRowMin = nRow;
3307 5849 : int nRowMax = nRowMin;
3308 5849 : if (m_poTPD->m_nShiftYPixelsMod)
3309 1360 : nRowMax++;
3310 :
3311 5849 : int nColMin = nCol;
3312 5849 : int nColMax = nColMin;
3313 5849 : if (m_poTPD->m_nShiftXPixelsMod)
3314 1322 : nColMax++;
3315 :
3316 5849 : CPLErr eErr = CE_None;
3317 :
3318 13058 : for (nRow = nRowMin; eErr == CE_None && nRow <= nRowMax; nRow++)
3319 : {
3320 17062 : for (nCol = nColMin; eErr == CE_None && nCol <= nColMax; nCol++)
3321 : {
3322 9853 : if (nRow < 0 || nCol < 0 || nRow >= m_poTPD->m_nTileMatrixHeight ||
3323 9713 : nCol >= m_poTPD->m_nTileMatrixWidth)
3324 : {
3325 167 : continue;
3326 : }
3327 :
3328 9686 : if (m_poTPD->m_nShiftXPixelsMod == 0 &&
3329 4515 : m_poTPD->m_nShiftYPixelsMod == 0)
3330 : {
3331 4457 : if (!(nRow == m_poTPD->m_asCachedTilesDesc[0].nRow &&
3332 250 : nCol == m_poTPD->m_asCachedTilesDesc[0].nCol &&
3333 5 : m_poTPD->m_asCachedTilesDesc[0].nIdxWithinTileData == 0))
3334 : {
3335 4452 : eErr = m_poTPD->WriteTile();
3336 :
3337 4452 : m_poTPD->m_asCachedTilesDesc[0].nRow = nRow;
3338 4452 : m_poTPD->m_asCachedTilesDesc[0].nCol = nCol;
3339 4452 : m_poTPD->m_asCachedTilesDesc[0].nIdxWithinTileData = 0;
3340 : }
3341 : }
3342 :
3343 : // Composite block data into tile, and check if all bands for this
3344 : // block are dirty, and if so write the tile
3345 9686 : bool bAllDirty = true;
3346 42233 : for (int iBand = 1; iBand <= poDS->GetRasterCount(); iBand++)
3347 : {
3348 32547 : GDALRasterBlock *poBlock = nullptr;
3349 32547 : GByte *pabySrc = nullptr;
3350 32547 : if (iBand == nBand)
3351 : {
3352 9686 : pabySrc = static_cast<GByte *>(pData);
3353 : }
3354 : else
3355 : {
3356 22861 : if (!(m_poTPD->m_nShiftXPixelsMod == 0 &&
3357 8359 : m_poTPD->m_nShiftYPixelsMod == 0))
3358 14646 : continue;
3359 :
3360 : // If the block for this band is not dirty, it might be
3361 : // dirty in cache
3362 8215 : if (m_poTPD->m_asCachedTilesDesc[0].abBandDirty[iBand - 1])
3363 505 : continue;
3364 : else
3365 : {
3366 : poBlock =
3367 7710 : cpl::down_cast<GDALGPKGMBTilesLikeRasterBand *>(
3368 7710 : poDS->GetRasterBand(iBand))
3369 7710 : ->TryGetLockedBlockRef(nBlockXOff, nBlockYOff);
3370 7710 : if (poBlock && poBlock->GetDirty())
3371 : {
3372 : pabySrc =
3373 7670 : static_cast<GByte *>(poBlock->GetDataRef());
3374 7670 : poBlock->MarkClean();
3375 : }
3376 : else
3377 : {
3378 40 : if (poBlock)
3379 2 : poBlock->DropLock();
3380 40 : bAllDirty = false;
3381 40 : continue;
3382 : }
3383 : }
3384 : }
3385 :
3386 17356 : if (m_poTPD->m_nShiftXPixelsMod == 0 &&
3387 12185 : m_poTPD->m_nShiftYPixelsMod == 0)
3388 12127 : m_poTPD->m_asCachedTilesDesc[0].abBandDirty[iBand - 1] =
3389 : true;
3390 :
3391 17356 : int nDstXOffset = 0;
3392 17356 : int nDstXSize = nBlockXSize;
3393 17356 : int nDstYOffset = 0;
3394 17356 : int nDstYSize = nBlockYSize;
3395 : // Composite block data into tile data
3396 17356 : if (m_poTPD->m_nShiftXPixelsMod == 0 &&
3397 12185 : m_poTPD->m_nShiftYPixelsMod == 0)
3398 : {
3399 :
3400 : #ifdef DEBUG_VERBOSE
3401 : if (eDataType == GDT_Byte &&
3402 : nBlockXOff * nBlockXSize <=
3403 : nRasterXSize - nBlockXSize &&
3404 : nBlockYOff * nBlockYSize > nRasterYSize - nBlockYSize)
3405 : {
3406 : bool bFoundNonZero = false;
3407 : for (int y = nRasterYSize - nBlockYOff * nBlockYSize;
3408 : y < nBlockYSize; y++)
3409 : {
3410 : for (int x = 0; x < nBlockXSize; x++)
3411 : {
3412 : if (pabySrc[static_cast<GPtrDiff_t>(y) *
3413 : nBlockXSize +
3414 : x] != 0 &&
3415 : !bFoundNonZero)
3416 : {
3417 : CPLDebug("GPKG",
3418 : "IWriteBlock(): Found non-zero "
3419 : "content in ghost part of "
3420 : "tile(nBand=%d,nBlockXOff=%d,"
3421 : "nBlockYOff=%d,m_nZoomLevel=%d)\n",
3422 : iBand, nBlockXOff, nBlockYOff,
3423 : m_poTPD->m_nZoomLevel);
3424 : bFoundNonZero = true;
3425 : }
3426 : }
3427 : }
3428 : }
3429 : #endif
3430 :
3431 12127 : const size_t nBandBlockSize =
3432 12127 : static_cast<size_t>(nBlockXSize) * nBlockYSize *
3433 12127 : m_nDTSize;
3434 12127 : memcpy(m_poTPD->m_pabyCachedTiles +
3435 12127 : (iBand - 1) * nBandBlockSize,
3436 : pabySrc, nBandBlockSize);
3437 :
3438 : // Make sure partial blocks are zero'ed outside of the
3439 : // validity area but do that only when know that JPEG will
3440 : // not be used so as to avoid edge effects (although we
3441 : // should probably repeat last pixels if we really want to
3442 : // do that, but that only makes sense if readers only clip
3443 : // to the gpkg_contents extent). Well, ere on the safe side
3444 : // for now
3445 12127 : if (m_poTPD->m_eTF != GPKG_TF_JPEG &&
3446 11879 : (nBlockXOff * nBlockXSize >=
3447 11879 : nRasterXSize - nBlockXSize ||
3448 11338 : nBlockYOff * nBlockYSize >=
3449 11338 : nRasterYSize - nBlockYSize))
3450 : {
3451 1037 : int nXEndValidity =
3452 1037 : nRasterXSize - nBlockXOff * nBlockXSize;
3453 1037 : if (nXEndValidity > nBlockXSize)
3454 496 : nXEndValidity = nBlockXSize;
3455 1037 : int nYEndValidity =
3456 1037 : nRasterYSize - nBlockYOff * nBlockYSize;
3457 1037 : if (nYEndValidity > nBlockYSize)
3458 294 : nYEndValidity = nBlockYSize;
3459 1037 : if (nXEndValidity < nBlockXSize)
3460 : {
3461 9778 : for (int iY = 0; iY < nYEndValidity; iY++)
3462 : {
3463 9583 : m_poTPD->FillBuffer(
3464 9583 : m_poTPD->m_pabyCachedTiles +
3465 9583 : ((static_cast<size_t>(iBand - 1) *
3466 9583 : nBlockYSize +
3467 9583 : iY) *
3468 9583 : nBlockXSize +
3469 9583 : nXEndValidity) *
3470 9583 : m_nDTSize,
3471 9583 : nBlockXSize - nXEndValidity);
3472 : }
3473 : }
3474 1037 : if (nYEndValidity < nBlockYSize)
3475 : {
3476 235 : m_poTPD->FillBuffer(
3477 235 : m_poTPD->m_pabyCachedTiles +
3478 235 : (static_cast<size_t>(iBand - 1) *
3479 235 : nBlockYSize +
3480 235 : nYEndValidity) *
3481 235 : nBlockXSize * m_nDTSize,
3482 235 : static_cast<size_t>(nBlockYSize -
3483 : nYEndValidity) *
3484 235 : nBlockXSize);
3485 : }
3486 12127 : }
3487 : }
3488 : else
3489 : {
3490 5229 : const int nXValid =
3491 5229 : (nBlockXOff * nBlockXSize > nRasterXSize - nBlockXSize)
3492 5229 : ? (nRasterXSize - nBlockXOff * nBlockXSize)
3493 : : nBlockXSize;
3494 5229 : const int nYValid =
3495 5229 : (nBlockYOff * nBlockYSize > nRasterYSize - nBlockYSize)
3496 5229 : ? (nRasterYSize - nBlockYOff * nBlockYSize)
3497 : : nBlockYSize;
3498 :
3499 5229 : int nSrcXOffset = 0;
3500 5229 : if (nCol == nColMin)
3501 : {
3502 2643 : nDstXOffset = m_poTPD->m_nShiftXPixelsMod;
3503 2643 : nDstXSize = std::min(
3504 2643 : nXValid, nBlockXSize - m_poTPD->m_nShiftXPixelsMod);
3505 : }
3506 : else
3507 : {
3508 2586 : nDstXOffset = 0;
3509 2586 : if (nXValid > nBlockXSize - m_poTPD->m_nShiftXPixelsMod)
3510 : {
3511 2212 : nDstXSize = nXValid - (nBlockXSize -
3512 2212 : m_poTPD->m_nShiftXPixelsMod);
3513 : }
3514 : else
3515 374 : nDstXSize = 0;
3516 2586 : nSrcXOffset = nBlockXSize - m_poTPD->m_nShiftXPixelsMod;
3517 : }
3518 :
3519 5229 : int nSrcYOffset = 0;
3520 5229 : if (nRow == nRowMin)
3521 : {
3522 2607 : nDstYOffset = m_poTPD->m_nShiftYPixelsMod;
3523 2607 : nDstYSize = std::min(
3524 2607 : nYValid, nBlockYSize - m_poTPD->m_nShiftYPixelsMod);
3525 : }
3526 : else
3527 : {
3528 2622 : nDstYOffset = 0;
3529 2622 : if (nYValid > nBlockYSize - m_poTPD->m_nShiftYPixelsMod)
3530 : {
3531 2232 : nDstYSize = nYValid - (nBlockYSize -
3532 2232 : m_poTPD->m_nShiftYPixelsMod);
3533 : }
3534 : else
3535 390 : nDstYSize = 0;
3536 2622 : nSrcYOffset = nBlockYSize - m_poTPD->m_nShiftYPixelsMod;
3537 : }
3538 :
3539 : #ifdef DEBUG_VERBOSE
3540 : CPLDebug("GPKG",
3541 : "Copy source tile x=%d,w=%d,y=%d,h=%d into "
3542 : "buffer at x=%d,y=%d",
3543 : nDstXOffset, nDstXSize, nDstYOffset, nDstYSize,
3544 : nSrcXOffset, nSrcYOffset);
3545 : #endif
3546 :
3547 5229 : if (nDstXSize > 0 && nDstYSize > 0)
3548 : {
3549 556108 : for (GPtrDiff_t y = 0; y < nDstYSize; y++)
3550 : {
3551 551552 : GByte *pDst = m_poTPD->m_pabyCachedTiles +
3552 551552 : (static_cast<size_t>(iBand - 1) *
3553 551552 : nBlockXSize * nBlockYSize +
3554 551552 : (y + nDstYOffset) * nBlockXSize +
3555 551552 : nDstXOffset) *
3556 551552 : m_nDTSize;
3557 551552 : GByte *pSrc =
3558 551552 : pabySrc + ((y + nSrcYOffset) * nBlockXSize +
3559 551552 : nSrcXOffset) *
3560 551552 : m_nDTSize;
3561 551552 : GDALCopyWords(pSrc, eDataType, m_nDTSize, pDst,
3562 : eDataType, m_nDTSize, nDstXSize);
3563 : }
3564 : }
3565 : }
3566 :
3567 17356 : if (poBlock)
3568 7670 : poBlock->DropLock();
3569 :
3570 17356 : if (!(m_poTPD->m_nShiftXPixelsMod == 0 &&
3571 12185 : m_poTPD->m_nShiftYPixelsMod == 0))
3572 : {
3573 5229 : m_poTPD->m_asCachedTilesDesc[0].nRow = -1;
3574 5229 : m_poTPD->m_asCachedTilesDesc[0].nCol = -1;
3575 5229 : m_poTPD->m_asCachedTilesDesc[0].nIdxWithinTileData = -1;
3576 5229 : if (nDstXSize > 0 && nDstYSize > 0)
3577 : {
3578 4556 : eErr = m_poTPD->WriteShiftedTile(
3579 : nRow, nCol, iBand, nDstXOffset, nDstYOffset,
3580 : nDstXSize, nDstYSize);
3581 : }
3582 : }
3583 : }
3584 :
3585 9686 : if (m_poTPD->m_nShiftXPixelsMod == 0 &&
3586 4515 : m_poTPD->m_nShiftYPixelsMod == 0)
3587 : {
3588 4457 : if (bAllDirty)
3589 : {
3590 4435 : eErr = m_poTPD->WriteTile();
3591 : }
3592 : }
3593 : }
3594 : }
3595 :
3596 5849 : return eErr;
3597 : }
3598 :
3599 : /************************************************************************/
3600 : /* GetNoDataValue() */
3601 : /************************************************************************/
3602 :
3603 19907 : double GDALGPKGMBTilesLikeRasterBand::GetNoDataValue(int *pbSuccess)
3604 : {
3605 19907 : if (m_bHasNoData)
3606 : {
3607 401 : if (pbSuccess)
3608 400 : *pbSuccess = TRUE;
3609 401 : return m_dfNoDataValue;
3610 : }
3611 19506 : return GDALPamRasterBand::GetNoDataValue(pbSuccess);
3612 : }
3613 :
3614 : /************************************************************************/
3615 : /* SetNoDataValueInternal() */
3616 : /************************************************************************/
3617 :
3618 75 : void GDALGPKGMBTilesLikeRasterBand::SetNoDataValueInternal(double dfNoDataValue)
3619 : {
3620 75 : m_bHasNoData = true;
3621 75 : m_dfNoDataValue = dfNoDataValue;
3622 75 : }
3623 :
3624 : /************************************************************************/
3625 : /* GDALGeoPackageRasterBand() */
3626 : /************************************************************************/
3627 :
3628 1817 : GDALGeoPackageRasterBand::GDALGeoPackageRasterBand(
3629 1817 : GDALGeoPackageDataset *poDSIn, int nTileWidth, int nTileHeight)
3630 1817 : : GDALGPKGMBTilesLikeRasterBand(poDSIn, nTileWidth, nTileHeight)
3631 : {
3632 1817 : poDS = poDSIn;
3633 1817 : }
3634 :
3635 : /************************************************************************/
3636 : /* GetOverviewCount() */
3637 : /************************************************************************/
3638 :
3639 8 : int GDALGeoPackageRasterBand::GetOverviewCount()
3640 : {
3641 : GDALGeoPackageDataset *poGDS =
3642 8 : cpl::down_cast<GDALGeoPackageDataset *>(poDS);
3643 8 : return static_cast<int>(poGDS->m_apoOverviewDS.size());
3644 : }
3645 :
3646 : /************************************************************************/
3647 : /* GetOverviewCount() */
3648 : /************************************************************************/
3649 :
3650 41 : GDALRasterBand *GDALGeoPackageRasterBand::GetOverview(int nIdx)
3651 : {
3652 41 : GDALGeoPackageDataset *poGDS =
3653 : reinterpret_cast<GDALGeoPackageDataset *>(poDS);
3654 41 : if (nIdx < 0 || nIdx >= static_cast<int>(poGDS->m_apoOverviewDS.size()))
3655 1 : return nullptr;
3656 40 : return poGDS->m_apoOverviewDS[nIdx]->GetRasterBand(nBand);
3657 : }
3658 :
3659 : /************************************************************************/
3660 : /* SetNoDataValue() */
3661 : /************************************************************************/
3662 :
3663 23 : CPLErr GDALGeoPackageRasterBand::SetNoDataValue(double dfNoDataValue)
3664 : {
3665 : GDALGeoPackageDataset *poGDS =
3666 23 : cpl::down_cast<GDALGeoPackageDataset *>(poDS);
3667 :
3668 23 : if (eDataType == GDT_Byte)
3669 : {
3670 6 : if (!(dfNoDataValue >= 0 && dfNoDataValue <= 255 &&
3671 4 : static_cast<int>(dfNoDataValue) == dfNoDataValue))
3672 : {
3673 2 : CPLError(CE_Failure, CPLE_NotSupported,
3674 : "Invalid nodata value for a Byte band: %.17g",
3675 : dfNoDataValue);
3676 2 : return CE_Failure;
3677 : }
3678 :
3679 8 : for (int i = 1; i <= poGDS->nBands; ++i)
3680 : {
3681 5 : if (i != nBand)
3682 : {
3683 2 : int bHasNoData = FALSE;
3684 : double dfOtherNoDataValue =
3685 2 : poGDS->GetRasterBand(i)->GetNoDataValue(&bHasNoData);
3686 2 : if (bHasNoData && dfOtherNoDataValue != dfNoDataValue)
3687 : {
3688 1 : CPLError(
3689 : CE_Failure, CPLE_NotSupported,
3690 : "Only the same nodata value can be set on all bands");
3691 1 : return CE_Failure;
3692 : }
3693 : }
3694 : }
3695 :
3696 3 : SetNoDataValueInternal(dfNoDataValue);
3697 3 : poGDS->m_bMetadataDirty = true;
3698 3 : return CE_None;
3699 : }
3700 :
3701 17 : if (std::isnan(dfNoDataValue))
3702 : {
3703 0 : CPLError(CE_Warning, CPLE_NotSupported,
3704 : "A NaN nodata value cannot be recorded in "
3705 : "gpkg_2d_gridded_coverage_ancillary table");
3706 : }
3707 :
3708 17 : SetNoDataValueInternal(dfNoDataValue);
3709 :
3710 17 : char *pszSQL = sqlite3_mprintf(
3711 : "UPDATE gpkg_2d_gridded_coverage_ancillary SET data_null = ? "
3712 : "WHERE tile_matrix_set_name = '%q'",
3713 : poGDS->m_osRasterTable.c_str());
3714 17 : sqlite3_stmt *hStmt = nullptr;
3715 17 : int rc = SQLPrepareWithError(poGDS->IGetDB(), pszSQL, -1, &hStmt, nullptr);
3716 17 : if (rc == SQLITE_OK)
3717 : {
3718 17 : if (poGDS->m_eTF == GPKG_TF_PNG_16BIT)
3719 : {
3720 15 : if (eDataType == GDT_UInt16 && poGDS->m_dfOffset == 0.0 &&
3721 6 : poGDS->m_dfScale == 1.0 && dfNoDataValue >= 0 &&
3722 5 : dfNoDataValue <= 65535 &&
3723 5 : static_cast<GUInt16>(dfNoDataValue) == dfNoDataValue)
3724 : {
3725 5 : poGDS->m_usGPKGNull = static_cast<GUInt16>(dfNoDataValue);
3726 : }
3727 : else
3728 : {
3729 10 : poGDS->m_usGPKGNull = 65535;
3730 : }
3731 15 : sqlite3_bind_double(hStmt, 1, poGDS->m_usGPKGNull);
3732 : }
3733 : else
3734 : {
3735 2 : sqlite3_bind_double(hStmt, 1, static_cast<float>(dfNoDataValue));
3736 : }
3737 17 : rc = sqlite3_step(hStmt);
3738 17 : sqlite3_finalize(hStmt);
3739 : }
3740 17 : sqlite3_free(pszSQL);
3741 :
3742 17 : return (rc == SQLITE_OK) ? CE_None : CE_Failure;
3743 : }
3744 :
3745 : /************************************************************************/
3746 : /* LoadBandMetadata() */
3747 : /************************************************************************/
3748 :
3749 1239 : void GDALGeoPackageRasterBand::LoadBandMetadata()
3750 : {
3751 : GDALGeoPackageDataset *poGDS =
3752 1239 : cpl::down_cast<GDALGeoPackageDataset *>(poDS);
3753 :
3754 1239 : if (m_bHasReadMetadataFromStorage)
3755 1137 : return;
3756 :
3757 475 : m_bHasReadMetadataFromStorage = true;
3758 :
3759 475 : poGDS->TryLoadXML();
3760 :
3761 475 : if (!poGDS->HasMetadataTables())
3762 373 : return;
3763 :
3764 102 : char *pszSQL = sqlite3_mprintf(
3765 : "SELECT md.metadata, md.md_standard_uri, md.mime_type "
3766 : "FROM gpkg_metadata md "
3767 : "JOIN gpkg_metadata_reference mdr ON (md.id = mdr.md_file_id ) "
3768 : "WHERE "
3769 : "mdr.reference_scope = 'table' AND lower(mdr.table_name) = "
3770 : "lower('%q') ORDER BY md.id "
3771 : "LIMIT 1000", // to avoid denial of service
3772 : poGDS->m_osRasterTable.c_str());
3773 :
3774 102 : auto oResult = SQLQuery(poGDS->hDB, pszSQL);
3775 102 : sqlite3_free(pszSQL);
3776 102 : if (!oResult)
3777 : {
3778 0 : return;
3779 : }
3780 :
3781 : /* GDAL metadata */
3782 200 : for (int i = 0; i < oResult->RowCount(); i++)
3783 : {
3784 98 : const char *pszMetadata = oResult->GetValue(0, i);
3785 98 : const char *pszMDStandardURI = oResult->GetValue(1, i);
3786 98 : const char *pszMimeType = oResult->GetValue(2, i);
3787 98 : if (pszMetadata && pszMDStandardURI && pszMimeType &&
3788 98 : EQUAL(pszMDStandardURI, "http://gdal.org") &&
3789 74 : EQUAL(pszMimeType, "text/xml"))
3790 : {
3791 74 : CPLXMLNode *psXMLNode = CPLParseXMLString(pszMetadata);
3792 74 : if (psXMLNode)
3793 : {
3794 148 : GDALMultiDomainMetadata oLocalMDMD;
3795 74 : oLocalMDMD.XMLInit(psXMLNode, FALSE);
3796 :
3797 74 : CSLConstList papszDomainList = oLocalMDMD.GetDomainList();
3798 208 : for (CSLConstList papszIter = papszDomainList;
3799 208 : papszIter && *papszIter; ++papszIter)
3800 : {
3801 134 : if (STARTS_WITH(*papszIter, "BAND_"))
3802 : {
3803 5 : int l_nBand = atoi(*papszIter + strlen("BAND_"));
3804 5 : if (l_nBand >= 1 && l_nBand <= poGDS->GetRasterCount())
3805 : {
3806 : auto l_poBand =
3807 5 : cpl::down_cast<GDALGeoPackageRasterBand *>(
3808 : poGDS->GetRasterBand(l_nBand));
3809 5 : l_poBand->m_bHasReadMetadataFromStorage = true;
3810 :
3811 10 : char **papszMD = CSLDuplicate(
3812 5 : oLocalMDMD.GetMetadata(*papszIter));
3813 10 : papszMD = CSLMerge(
3814 : papszMD,
3815 5 : GDALGPKGMBTilesLikeRasterBand::GetMetadata(""));
3816 5 : l_poBand->GDALPamRasterBand::SetMetadata(papszMD);
3817 5 : CSLDestroy(papszMD);
3818 : }
3819 : }
3820 : }
3821 :
3822 74 : CPLDestroyXMLNode(psXMLNode);
3823 : }
3824 : }
3825 : }
3826 : }
3827 :
3828 : /************************************************************************/
3829 : /* GetMetadata() */
3830 : /************************************************************************/
3831 :
3832 952 : char **GDALGeoPackageRasterBand::GetMetadata(const char *pszDomain)
3833 : {
3834 952 : GDALGeoPackageDataset *poGDS =
3835 : reinterpret_cast<GDALGeoPackageDataset *>(poDS);
3836 952 : LoadBandMetadata(); /* force loading from storage if needed */
3837 :
3838 15 : if (poGDS->eAccess == GA_ReadOnly && eDataType != GDT_Byte &&
3839 11 : (pszDomain == nullptr || EQUAL(pszDomain, "")) &&
3840 11 : !m_bMinMaxComputedFromTileAncillary &&
3841 975 : !GDALGPKGMBTilesLikeRasterBand::GetMetadataItem("STATISTICS_MINIMUM") &&
3842 8 : !GDALGPKGMBTilesLikeRasterBand::GetMetadataItem("STATISTICS_MAXIMUM"))
3843 : {
3844 8 : m_bMinMaxComputedFromTileAncillary = true;
3845 :
3846 8 : const int nColMin = poGDS->m_nShiftXTiles;
3847 8 : const int nColMax =
3848 8 : (nRasterXSize - 1 + poGDS->m_nShiftXPixelsMod) / nBlockXSize +
3849 8 : poGDS->m_nShiftXTiles;
3850 8 : const int nRowMin = poGDS->m_nShiftYTiles;
3851 8 : const int nRowMax =
3852 8 : (nRasterYSize - 1 + poGDS->m_nShiftYPixelsMod) / nBlockYSize +
3853 8 : poGDS->m_nShiftYTiles;
3854 :
3855 8 : bool bOK = false;
3856 8 : if (poGDS->m_nShiftXPixelsMod == 0 && poGDS->m_nShiftYPixelsMod == 0 &&
3857 8 : (nRasterXSize % nBlockXSize) == 0 &&
3858 2 : (nRasterYSize % nBlockYSize) == 0)
3859 : {
3860 : // If the area of interest matches entire tiles, then we can
3861 : // use tile statistics
3862 2 : bOK = true;
3863 : }
3864 6 : else if (m_bHasNoData)
3865 : {
3866 : // Otherwise, in the case where we have nodata, we assume that
3867 : // if the area of interest is at least larger than the existing
3868 : // tiles, the tile statistics will be reliable.
3869 2 : char *pszSQL = sqlite3_mprintf(
3870 : "SELECT MIN(tile_column), MAX(tile_column), "
3871 : "MIN(tile_row), MAX(tile_row) FROM \"%w\" "
3872 : "WHERE zoom_level = %d",
3873 : poGDS->m_osRasterTable.c_str(), poGDS->m_nZoomLevel);
3874 4 : auto sResult = SQLQuery(poGDS->IGetDB(), pszSQL);
3875 2 : if (sResult && sResult->RowCount() == 1)
3876 : {
3877 2 : const char *pszMinX = sResult->GetValue(0, 0);
3878 2 : const char *pszMaxX = sResult->GetValue(1, 0);
3879 2 : const char *pszMinY = sResult->GetValue(2, 0);
3880 2 : const char *pszMaxY = sResult->GetValue(3, 0);
3881 2 : if (pszMinX && pszMaxX && pszMinY && pszMaxY)
3882 : {
3883 2 : bOK = atoi(pszMinX) >= nColMin &&
3884 2 : atoi(pszMaxX) <= nColMax &&
3885 4 : atoi(pszMinY) >= nRowMin && atoi(pszMaxY) <= nRowMax;
3886 : }
3887 : }
3888 2 : sqlite3_free(pszSQL);
3889 : }
3890 :
3891 8 : if (bOK)
3892 : {
3893 4 : char *pszSQL = sqlite3_mprintf(
3894 : "SELECT MIN(min), MAX(max) FROM "
3895 : "gpkg_2d_gridded_tile_ancillary WHERE tpudt_id "
3896 : "IN (SELECT id FROM \"%w\" WHERE "
3897 : "zoom_level = %d AND "
3898 : "tile_column >= %d AND tile_column <= %d AND "
3899 : "tile_row >= %d AND tile_row <= %d)",
3900 : poGDS->m_osRasterTable.c_str(), poGDS->m_nZoomLevel, nColMin,
3901 : nColMax, nRowMin, nRowMax);
3902 8 : auto sResult = SQLQuery(poGDS->IGetDB(), pszSQL);
3903 4 : CPLDebug("GPKG", "%s", pszSQL);
3904 4 : if (sResult && sResult->RowCount() == 1)
3905 : {
3906 4 : const char *pszMin = sResult->GetValue(0, 0);
3907 4 : const char *pszMax = sResult->GetValue(1, 0);
3908 4 : if (pszMin)
3909 : {
3910 4 : m_dfStatsMinFromTileAncillary = CPLAtof(pszMin);
3911 : }
3912 4 : if (pszMax)
3913 : {
3914 4 : m_dfStatsMaxFromTileAncillary = CPLAtof(pszMax);
3915 : }
3916 : }
3917 4 : sqlite3_free(pszSQL);
3918 : }
3919 : }
3920 :
3921 705 : if (m_bAddImplicitStatistics && m_bMinMaxComputedFromTileAncillary &&
3922 8 : (pszDomain == nullptr || EQUAL(pszDomain, "")) &&
3923 1665 : !GDALGPKGMBTilesLikeRasterBand::GetMetadataItem("STATISTICS_MINIMUM") &&
3924 8 : !GDALGPKGMBTilesLikeRasterBand::GetMetadataItem("STATISTICS_MAXIMUM"))
3925 : {
3926 : m_aosMD.Assign(CSLDuplicate(
3927 8 : GDALGPKGMBTilesLikeRasterBand::GetMetadata(pszDomain)));
3928 8 : if (!std::isnan(m_dfStatsMinFromTileAncillary))
3929 : {
3930 : m_aosMD.SetNameValue(
3931 : "STATISTICS_MINIMUM",
3932 4 : CPLSPrintf("%.14g", m_dfStatsMinFromTileAncillary));
3933 : }
3934 8 : if (!std::isnan(m_dfStatsMaxFromTileAncillary))
3935 : {
3936 : m_aosMD.SetNameValue(
3937 : "STATISTICS_MAXIMUM",
3938 4 : CPLSPrintf("%.14g", m_dfStatsMaxFromTileAncillary));
3939 : }
3940 8 : return m_aosMD.List();
3941 : }
3942 :
3943 944 : return GDALGPKGMBTilesLikeRasterBand::GetMetadata(pszDomain);
3944 : }
3945 :
3946 : /************************************************************************/
3947 : /* GetMetadataItem() */
3948 : /************************************************************************/
3949 :
3950 251 : const char *GDALGeoPackageRasterBand::GetMetadataItem(const char *pszName,
3951 : const char *pszDomain)
3952 : {
3953 251 : LoadBandMetadata(); /* force loading from storage if needed */
3954 :
3955 251 : if (m_bAddImplicitStatistics && eDataType != GDT_Byte &&
3956 7 : (pszDomain == nullptr || EQUAL(pszDomain, "")) &&
3957 4 : (EQUAL(pszName, "STATISTICS_MINIMUM") ||
3958 2 : EQUAL(pszName, "STATISTICS_MAXIMUM")))
3959 : {
3960 2 : return CSLFetchNameValue(GetMetadata(), pszName);
3961 : }
3962 :
3963 249 : return GDALGPKGMBTilesLikeRasterBand::GetMetadataItem(pszName, pszDomain);
3964 : }
3965 :
3966 : /************************************************************************/
3967 : /* SetMetadata() */
3968 : /************************************************************************/
3969 :
3970 16 : CPLErr GDALGeoPackageRasterBand::SetMetadata(char **papszMetadata,
3971 : const char *pszDomain)
3972 : {
3973 : GDALGeoPackageDataset *poGDS =
3974 16 : cpl::down_cast<GDALGeoPackageDataset *>(poDS);
3975 16 : LoadBandMetadata(); /* force loading from storage if needed */
3976 16 : poGDS->m_bMetadataDirty = true;
3977 72 : for (CSLConstList papszIter = papszMetadata; papszIter && *papszIter;
3978 : ++papszIter)
3979 : {
3980 56 : if (STARTS_WITH(*papszIter, "STATISTICS_"))
3981 51 : m_bStatsMetadataSetInThisSession = true;
3982 : }
3983 16 : return GDALPamRasterBand::SetMetadata(papszMetadata, pszDomain);
3984 : }
3985 :
3986 : /************************************************************************/
3987 : /* SetMetadataItem() */
3988 : /************************************************************************/
3989 :
3990 20 : CPLErr GDALGeoPackageRasterBand::SetMetadataItem(const char *pszName,
3991 : const char *pszValue,
3992 : const char *pszDomain)
3993 : {
3994 : GDALGeoPackageDataset *poGDS =
3995 20 : cpl::down_cast<GDALGeoPackageDataset *>(poDS);
3996 20 : LoadBandMetadata(); /* force loading from storage if needed */
3997 20 : poGDS->m_bMetadataDirty = true;
3998 20 : if (STARTS_WITH(pszName, "STATISTICS_"))
3999 20 : m_bStatsMetadataSetInThisSession = true;
4000 20 : return GDALPamRasterBand::SetMetadataItem(pszName, pszValue, pszDomain);
4001 : }
4002 :
4003 : /************************************************************************/
4004 : /* InvalidateStatistics() */
4005 : /************************************************************************/
4006 :
4007 344 : void GDALGeoPackageRasterBand::InvalidateStatistics()
4008 : {
4009 344 : bool bModified = false;
4010 688 : CPLStringList aosMD(CSLDuplicate(GetMetadata()));
4011 365 : for (CSLConstList papszIter = GetMetadata(); papszIter && *papszIter;
4012 : ++papszIter)
4013 : {
4014 21 : if (STARTS_WITH(*papszIter, "STATISTICS_"))
4015 : {
4016 20 : char *pszKey = nullptr;
4017 20 : CPLParseNameValue(*papszIter, &pszKey);
4018 20 : CPLAssert(pszKey);
4019 20 : aosMD.SetNameValue(pszKey, nullptr);
4020 20 : CPLFree(pszKey);
4021 20 : bModified = true;
4022 : }
4023 : }
4024 344 : if (bModified)
4025 4 : SetMetadata(aosMD.List());
4026 344 : }
|