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