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