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 2526 : GDALGPKGMBTilesLikePseudoDataset::GDALGPKGMBTilesLikePseudoDataset()
34 : : m_bForceTempDBCompaction(
35 2526 : CPLTestBool(CPLGetConfigOption("GPKG_FORCE_TEMPDB_COMPACTION", "NO")))
36 : {
37 12630 : for (int i = 0; i < 4; i++)
38 : {
39 10104 : m_asCachedTilesDesc[i].nRow = -1;
40 10104 : m_asCachedTilesDesc[i].nCol = -1;
41 10104 : m_asCachedTilesDesc[i].nIdxWithinTileData = -1;
42 10104 : m_asCachedTilesDesc[i].abBandDirty[0] = FALSE;
43 10104 : m_asCachedTilesDesc[i].abBandDirty[1] = FALSE;
44 10104 : m_asCachedTilesDesc[i].abBandDirty[2] = FALSE;
45 10104 : m_asCachedTilesDesc[i].abBandDirty[3] = FALSE;
46 : }
47 2526 : }
48 :
49 : /************************************************************************/
50 : /* ~GDALGPKGMBTilesLikePseudoDataset() */
51 : /************************************************************************/
52 :
53 2526 : GDALGPKGMBTilesLikePseudoDataset::~GDALGPKGMBTilesLikePseudoDataset()
54 : {
55 2526 : 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 2526 : CPLFree(m_pabyCachedTiles);
68 2526 : delete m_poCT;
69 2526 : CPLFree(m_pabyHugeColorArray);
70 2526 : }
71 :
72 : /************************************************************************/
73 : /* SetDataType() */
74 : /************************************************************************/
75 :
76 306 : void GDALGPKGMBTilesLikePseudoDataset::SetDataType(GDALDataType eDT)
77 : {
78 306 : CPLAssert(eDT == GDT_Byte || eDT == GDT_Int16 || eDT == GDT_UInt16 ||
79 : eDT == GDT_Float32);
80 306 : m_eDT = eDT;
81 306 : m_nDTSize = GDALGetDataTypeSizeBytes(m_eDT);
82 306 : }
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 2526 : GDALGPKGMBTilesLikeRasterBand::GDALGPKGMBTilesLikeRasterBand(
107 2526 : GDALGPKGMBTilesLikePseudoDataset *poTPD, int nTileWidth, int nTileHeight)
108 2526 : : m_poTPD(poTPD)
109 : {
110 2526 : eDataType = m_poTPD->m_eDT;
111 2526 : m_nDTSize = m_poTPD->m_nDTSize;
112 2526 : nBlockXSize = nTileWidth;
113 2526 : nBlockYSize = nTileHeight;
114 2526 : }
115 :
116 : #ifdef __GNUC__
117 : #pragma GCC diagnostic pop
118 : #endif
119 :
120 : /************************************************************************/
121 : /* FlushCache() */
122 : /************************************************************************/
123 :
124 5451 : CPLErr GDALGPKGMBTilesLikeRasterBand::FlushCache(bool bAtClosing)
125 : {
126 5451 : m_poTPD->m_nLastSpaceCheckTimestamp = -1; // disable partial flushes
127 5451 : CPLErr eErr = GDALPamRasterBand::FlushCache(bAtClosing);
128 5451 : if (eErr == CE_None)
129 5448 : eErr = m_poTPD->IFlushCacheWithErrCode(bAtClosing);
130 5451 : m_poTPD->m_nLastSpaceCheckTimestamp = 0;
131 5451 : return eErr;
132 : }
133 :
134 : /************************************************************************/
135 : /* FlushTiles() */
136 : /************************************************************************/
137 :
138 3248 : CPLErr GDALGPKGMBTilesLikePseudoDataset::FlushTiles()
139 : {
140 3248 : CPLErr eErr = CE_None;
141 3248 : GDALGPKGMBTilesLikePseudoDataset *poMainDS =
142 3248 : m_poParentDS ? m_poParentDS : this;
143 3248 : if (poMainDS->m_nTileInsertionCount < 0)
144 7 : return CE_Failure;
145 :
146 3241 : if (IGetUpdate())
147 : {
148 1932 : if (m_nShiftXPixelsMod || m_nShiftYPixelsMod)
149 : {
150 607 : eErr = FlushRemainingShiftedTiles(false /* total flush*/);
151 : }
152 : else
153 : {
154 1325 : eErr = WriteTile();
155 : }
156 : }
157 :
158 3241 : if (poMainDS->m_nTileInsertionCount > 0)
159 : {
160 180 : if (poMainDS->ICommitTransaction() != OGRERR_NONE)
161 : {
162 0 : poMainDS->m_nTileInsertionCount = -1;
163 0 : eErr = CE_Failure;
164 : }
165 : else
166 : {
167 180 : poMainDS->m_nTileInsertionCount = 0;
168 : }
169 : }
170 3241 : return eErr;
171 : }
172 :
173 : /************************************************************************/
174 : /* GetColorTable() */
175 : /************************************************************************/
176 :
177 335 : GDALColorTable *GDALGPKGMBTilesLikeRasterBand::GetColorTable()
178 : {
179 335 : if (poDS->GetRasterCount() != 1)
180 127 : return nullptr;
181 :
182 208 : if (!m_poTPD->m_bTriedEstablishingCT)
183 : {
184 60 : m_poTPD->m_bTriedEstablishingCT = true;
185 60 : 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 62 : for (int i = 0; i < 2; i++)
195 : {
196 57 : bool bRetry = false;
197 57 : char *pszSQL = nullptr;
198 57 : if (i == 0)
199 : {
200 52 : pszSQL = sqlite3_mprintf("SELECT tile_data FROM \"%w\" "
201 : "WHERE zoom_level = %d LIMIT 1",
202 52 : m_poTPD->m_osRasterTable.c_str(),
203 52 : m_poTPD->m_nZoomLevel);
204 : }
205 : else
206 : {
207 : // Try a tile in the middle of the raster
208 5 : pszSQL = sqlite3_mprintf(
209 : "SELECT tile_data FROM \"%w\" "
210 : "WHERE zoom_level = %d AND tile_column = %d AND tile_row = "
211 : "%d",
212 5 : m_poTPD->m_osRasterTable.c_str(), m_poTPD->m_nZoomLevel,
213 5 : m_poTPD->m_nShiftXTiles + nRasterXSize / 2 / nBlockXSize,
214 5 : m_poTPD->GetRowFromIntoTopConvention(
215 5 : m_poTPD->m_nShiftYTiles +
216 5 : nRasterYSize / 2 / nBlockYSize));
217 : }
218 57 : sqlite3_stmt *hStmt = nullptr;
219 57 : int rc = sqlite3_prepare_v2(m_poTPD->IGetDB(), pszSQL, -1, &hStmt,
220 : nullptr);
221 57 : if (rc == SQLITE_OK)
222 : {
223 57 : rc = sqlite3_step(hStmt);
224 72 : if (rc == SQLITE_ROW &&
225 15 : sqlite3_column_type(hStmt, 0) == SQLITE_BLOB)
226 : {
227 15 : const int nBytes = sqlite3_column_bytes(hStmt, 0);
228 : GByte *pabyRawData = reinterpret_cast<GByte *>(
229 15 : const_cast<void *>(sqlite3_column_blob(hStmt, 0)));
230 : const CPLString osMemFileName(
231 30 : VSIMemGenerateHiddenFilename("gpkg_read_tile"));
232 15 : VSILFILE *fp = VSIFileFromMemBuffer(
233 : osMemFileName.c_str(), pabyRawData, nBytes, FALSE);
234 15 : VSIFCloseL(fp);
235 :
236 : // Only PNG can have color table.
237 15 : 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 30 : apszDrivers, nullptr, nullptr));
242 15 : if (poDSTile != nullptr)
243 : {
244 15 : 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 10 : bRetry = true;
254 : }
255 : }
256 : else
257 0 : bRetry = true;
258 :
259 15 : VSIUnlink(osMemFileName);
260 : }
261 : }
262 57 : sqlite3_free(pszSQL);
263 57 : sqlite3_finalize(hStmt);
264 57 : if (!bRetry)
265 47 : break;
266 : }
267 : }
268 :
269 200 : 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 208 : GDALColorInterp GDALGPKGMBTilesLikeRasterBand::GetColorInterpretation()
316 : {
317 208 : if (m_poTPD->m_eDT != GDT_Byte)
318 22 : return GCI_Undefined;
319 186 : if (poDS->GetRasterCount() == 1)
320 48 : 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 13842 : void GDALGPKGMBTilesLikePseudoDataset::FillBuffer(GByte *pabyData,
383 : size_t nPixels)
384 : {
385 13842 : int bHasNoData = FALSE;
386 13842 : const double dfNoDataValue = IGetRasterBand(1)->GetNoDataValue(&bHasNoData);
387 13842 : if (!bHasNoData || dfNoDataValue == 0.0)
388 : {
389 : // cppcheck-suppress nullPointer
390 13540 : 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 13842 : }
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 2581 : CPLErr GDALGPKGMBTilesLikePseudoDataset::ReadTile(
430 : const CPLString &osMemFileName, GByte *pabyTileData, double dfTileOffset,
431 : double dfTileScale, bool *pbIsLossyFormat)
432 : {
433 2581 : const char *const apszDriversByte[] = {"JPEG", "PNG", "WEBP", nullptr};
434 2581 : const char *const apszDriversInt[] = {"PNG", nullptr};
435 2581 : const char *const apszDriversFloat[] = {"GTiff", nullptr};
436 : int nBlockXSize, nBlockYSize;
437 2581 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
438 2581 : const int nBands = IGetRasterCount();
439 : auto poDSTile = std::unique_ptr<GDALDataset>(GDALDataset::Open(
440 : osMemFileName.c_str(), GDAL_OF_RASTER | GDAL_OF_INTERNAL,
441 2581 : (m_eDT == GDT_Byte) ? apszDriversByte
442 43 : : (m_eTF == GPKG_TF_TIFF_32BIT_FLOAT) ? apszDriversFloat
443 : : apszDriversInt,
444 5162 : nullptr, nullptr));
445 2581 : 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 2579 : const int nTileBandCount = poDSTile->GetRasterCount();
453 :
454 5158 : if (!(poDSTile->GetRasterXSize() == nBlockXSize &&
455 2579 : poDSTile->GetRasterYSize() == nBlockYSize &&
456 7737 : (nTileBandCount >= 1 && nTileBandCount <= 4)) ||
457 2579 : (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 2579 : GDALDataType eRequestDT = GDT_Byte;
466 2579 : 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 2539 : 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 2579 : if (poDSTile->RasterIO(GF_Read, 0, 0, nBlockXSize, nBlockYSize,
479 : pabyTileData, nBlockXSize, nBlockYSize, eRequestDT,
480 : poDSTile->GetRasterCount(), nullptr, 0, 0, 0,
481 2579 : nullptr) != CE_None)
482 : {
483 0 : FillEmptyTile(pabyTileData);
484 0 : return CE_Failure;
485 : }
486 :
487 2579 : 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 2536 : GDALColorTable *poCT = nullptr;
566 2536 : if (nBands == 1 || nTileBandCount == 1)
567 : {
568 152 : poCT = poDSTile->GetRasterBand(1)->GetColorTable();
569 152 : IGetRasterBand(1)->GetColorTable();
570 : }
571 :
572 2536 : 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 2536 : const GPtrDiff_t nBlockPixels =
579 2536 : static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize;
580 2536 : 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 30 : if (nBands == 1 && nTileBandCount == 1 && poCT != nullptr &&
626 2565 : m_poCT != nullptr && !poCT->IsSame(m_poCT))
627 : {
628 0 : CPLError(CE_Warning, CPLE_NotSupported,
629 : "Different color tables. Unhandled for now");
630 : }
631 2535 : else if ((nBands == 1 && nTileBandCount >= 3) ||
632 30 : (nBands == 1 && nTileBandCount == 1 && m_poCT != nullptr &&
633 2535 : poCT == nullptr) ||
634 2535 : ((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 2535 : 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 2228 : else if (nTileBandCount == 2)
657 : {
658 : /* Do Grey+Alpha -> RGBA */
659 361 : memcpy(pabyTileData + 3 * nBlockPixels, pabyTileData + 1 * nBlockPixels,
660 : nBlockPixels);
661 361 : memcpy(pabyTileData + 1 * nBlockPixels, pabyTileData, nBlockPixels);
662 361 : 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 2535 : return CE_None;
712 : }
713 :
714 : /************************************************************************/
715 : /* ReadTile() */
716 : /************************************************************************/
717 :
718 7541 : GByte *GDALGPKGMBTilesLikePseudoDataset::ReadTile(int nRow, int nCol)
719 : {
720 : int nBlockXSize, nBlockYSize;
721 7541 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
722 7541 : const int nBands = IGetRasterCount();
723 7541 : const size_t nBandBlockSize =
724 7541 : static_cast<size_t>(nBlockXSize) * nBlockYSize * m_nDTSize;
725 7541 : const int nTileBands = m_eDT == GDT_Byte ? 4 : 1;
726 7541 : 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 2738 : GByte *pabyDest = m_pabyCachedTiles + 2 * nTileBands * nBandBlockSize;
778 2738 : bool bAllNonDirty = true;
779 11225 : for (int i = 0; i < nBands; i++)
780 : {
781 8487 : if (m_asCachedTilesDesc[0].abBandDirty[i])
782 2 : bAllNonDirty = false;
783 : }
784 2738 : if (bAllNonDirty)
785 : {
786 2736 : 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 2581 : void GDALGPKGMBTilesLikePseudoDataset::GetTileOffsetAndScale(
812 : GIntBig nTileId, double &dfTileOffset, double &dfTileScale)
813 : {
814 2581 : dfTileOffset = 0.0;
815 2581 : dfTileScale = 1.0;
816 :
817 2581 : 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 2581 : }
841 :
842 : /************************************************************************/
843 : /* ReadTile() */
844 : /************************************************************************/
845 :
846 5519 : GByte *GDALGPKGMBTilesLikePseudoDataset::ReadTile(int nRow, int nCol,
847 : GByte *pabyData,
848 : bool *pbIsLossyFormat)
849 : {
850 5519 : int nBlockXSize = 0;
851 5519 : int nBlockYSize = 0;
852 5519 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
853 5519 : const int nBands = IGetRasterCount();
854 :
855 5519 : if (pbIsLossyFormat)
856 18 : *pbIsLossyFormat = false;
857 :
858 5519 : const size_t nBandBlockSize =
859 5519 : static_cast<size_t>(nBlockXSize) * nBlockYSize * m_nDTSize;
860 5519 : if (nRow < 0 || nCol < 0 || nRow >= m_nTileMatrixHeight ||
861 5462 : 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 16338 : 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 5446 : m_eDT != GDT_Byte ? ", id" : "", // MBTiles do not have an id
875 : m_osRasterTable.c_str(), m_nZoomLevel,
876 5446 : GetRowFromIntoTopConvention(nRow), nCol,
877 5446 : !m_osWHERE.empty() ? CPLSPrintf(" AND (%s)", m_osWHERE.c_str()) : "");
878 :
879 : #ifdef DEBUG_VERBOSE
880 : CPLDebug("GPKG", "%s", pszSQL);
881 : #endif
882 :
883 5446 : sqlite3_stmt *hStmt = nullptr;
884 5446 : int rc = sqlite3_prepare_v2(IGetDB(), pszSQL, -1, &hStmt, nullptr);
885 5446 : 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 5446 : sqlite3_free(pszSQL);
893 5446 : rc = sqlite3_step(hStmt);
894 :
895 5446 : if (rc == SQLITE_ROW && sqlite3_column_type(hStmt, 0) == SQLITE_BLOB)
896 : {
897 2564 : const int nBytes = sqlite3_column_bytes(hStmt, 0);
898 : GIntBig nTileId =
899 2564 : (m_eDT == GDT_Byte) ? 0 : sqlite3_column_int64(hStmt, 1);
900 : GByte *pabyRawData = static_cast<GByte *>(
901 2564 : const_cast<void *>(sqlite3_column_blob(hStmt, 0)));
902 : const CPLString osMemFileName(
903 5128 : VSIMemGenerateHiddenFilename("gpkg_read_tile"));
904 2564 : VSILFILE *fp = VSIFileFromMemBuffer(osMemFileName.c_str(), pabyRawData,
905 : nBytes, FALSE);
906 2564 : VSIFCloseL(fp);
907 :
908 2564 : double dfTileOffset = 0.0;
909 2564 : double dfTileScale = 1.0;
910 2564 : GetTileOffsetAndScale(nTileId, dfTileOffset, dfTileScale);
911 2564 : ReadTile(osMemFileName, pabyData, dfTileOffset, dfTileScale,
912 : pbIsLossyFormat);
913 2564 : VSIUnlink(osMemFileName);
914 2564 : 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 5446 : return pabyData;
986 : }
987 :
988 : /************************************************************************/
989 : /* IReadBlock() */
990 : /************************************************************************/
991 :
992 3949 : 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 3949 : if (m_poTPD->m_pabyCachedTiles == nullptr)
1002 0 : return CE_Failure;
1003 :
1004 3949 : const int nRowMin = nBlockYOff + m_poTPD->m_nShiftYTiles;
1005 3949 : int nRowMax = nRowMin;
1006 3949 : if (m_poTPD->m_nShiftYPixelsMod)
1007 1211 : nRowMax++;
1008 :
1009 3949 : const int nColMin = nBlockXOff + m_poTPD->m_nShiftXTiles;
1010 3949 : int nColMax = nColMin;
1011 3949 : if (m_poTPD->m_nShiftXPixelsMod)
1012 1189 : nColMax++;
1013 :
1014 2760 : retry:
1015 : /* Optimize for left to right reading at constant row */
1016 3952 : 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 9112 : for (int nRow = nRowMin; nRow <= nRowMax; nRow++)
1055 : {
1056 12701 : for (int nCol = nColMin; nCol <= nColMax; nCol++)
1057 : {
1058 7541 : if (m_poTPD->m_nShiftXPixelsMod == 0 &&
1059 2785 : m_poTPD->m_nShiftYPixelsMod == 0)
1060 : {
1061 2738 : 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 2736 : if (m_poTPD->WriteTile() != CE_None)
1066 0 : return CE_Failure;
1067 : }
1068 : }
1069 :
1070 7541 : GByte *pabyTileData = m_poTPD->ReadTile(nRow, nCol);
1071 7541 : if (pabyTileData == nullptr)
1072 0 : return CE_Failure;
1073 :
1074 34598 : for (int iBand = 1; iBand <= poDS->GetRasterCount(); iBand++)
1075 : {
1076 27060 : GDALRasterBlock *poBlock = nullptr;
1077 27060 : GByte *pabyDest = nullptr;
1078 27060 : if (iBand == nBand)
1079 : {
1080 7541 : 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 27055 : if (m_poTPD->m_nShiftXPixelsMod == 0 &&
1114 8639 : m_poTPD->m_nShiftYPixelsMod == 0)
1115 : {
1116 8485 : const size_t nBandBlockSize =
1117 8485 : static_cast<size_t>(nBlockXSize) * nBlockYSize *
1118 8485 : m_nDTSize;
1119 8485 : memcpy(pabyDest,
1120 8485 : pabyTileData + (iBand - 1) * nBandBlockSize,
1121 8485 : 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 27055 : if (poBlock)
1211 19514 : poBlock->DropLock();
1212 : }
1213 : }
1214 : }
1215 :
1216 3949 : 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 13350 : CPLErr GDALGPKGMBTilesLikePseudoDataset::WriteTile()
1561 : {
1562 13350 : GDALGPKGMBTilesLikePseudoDataset *poMainDS =
1563 13350 : m_poParentDS ? m_poParentDS : this;
1564 13350 : if (poMainDS->m_nTileInsertionCount < 0)
1565 501 : return CE_Failure;
1566 :
1567 12849 : 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 12849 : GDALRasterBlock::EnterDisableDirtyBlockFlush();
1577 12849 : m_bInWriteTile = true;
1578 12849 : CPLErr eErr = WriteTileInternal();
1579 12849 : m_bInWriteTile = false;
1580 12849 : GDALRasterBlock::LeaveDisableDirtyBlockFlush();
1581 12849 : return eErr;
1582 : }
1583 :
1584 : /* should only be called by WriteTile() */
1585 12849 : CPLErr GDALGPKGMBTilesLikePseudoDataset::WriteTileInternal()
1586 : {
1587 17463 : if (!(IGetUpdate() && m_asCachedTilesDesc[0].nRow >= 0 &&
1588 4614 : m_asCachedTilesDesc[0].nCol >= 0 &&
1589 4614 : m_asCachedTilesDesc[0].nIdxWithinTileData == 0))
1590 8235 : return CE_None;
1591 :
1592 4614 : int nRow = m_asCachedTilesDesc[0].nRow;
1593 4614 : int nCol = m_asCachedTilesDesc[0].nCol;
1594 :
1595 4614 : bool bAllDirty = true;
1596 4614 : bool bAllNonDirty = true;
1597 4614 : const int nBands = IGetRasterCount();
1598 17980 : for (int i = 0; i < nBands; i++)
1599 : {
1600 13366 : if (m_asCachedTilesDesc[0].abBandDirty[i])
1601 13333 : bAllNonDirty = false;
1602 : else
1603 33 : bAllDirty = false;
1604 : }
1605 4614 : if (bAllNonDirty)
1606 0 : return CE_None;
1607 :
1608 : int nBlockXSize, nBlockYSize;
1609 4614 : 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 4614 : bool bIsLossyFormat = false;
1614 4614 : const size_t nBandBlockSize =
1615 4614 : static_cast<size_t>(nBlockXSize) * nBlockYSize * m_nDTSize;
1616 4614 : 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 4614 : int nXOff = (nCol - m_nShiftXTiles) * nBlockXSize - m_nShiftXPixelsMod;
1639 4614 : int nYOff = (nRow - m_nShiftYTiles) * nBlockYSize - m_nShiftYPixelsMod;
1640 :
1641 : /* Assert that the tile at least intersects some of the GDAL raster space */
1642 4614 : CPLAssert(nXOff > -nBlockXSize);
1643 4614 : CPLAssert(nYOff > -nBlockYSize);
1644 : /* Can happen if the tile of the raster is less than the block size */
1645 4614 : const int nRasterXSize = IGetRasterBand(1)->GetXSize();
1646 4614 : const int nRasterYSize = IGetRasterBand(1)->GetYSize();
1647 4614 : 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 4614 : int iXOff = 0;
1688 4614 : GPtrDiff_t iYOff = 0;
1689 4614 : int iXCount = nBlockXSize;
1690 4614 : int iYCount = nBlockYSize;
1691 :
1692 4614 : bool bPartialTile = false;
1693 4614 : int nAlphaBand = (nBands == 2) ? 2 : (nBands == 4) ? 4 : 0;
1694 4614 : if (nAlphaBand == 0)
1695 : {
1696 3734 : if (nXOff < 0)
1697 : {
1698 4 : bPartialTile = true;
1699 4 : iXOff = -nXOff;
1700 4 : iXCount += nXOff;
1701 : }
1702 3734 : if (nXOff > nRasterXSize - nBlockXSize)
1703 : {
1704 85 : bPartialTile = true;
1705 85 : iXCount -= static_cast<int>(static_cast<GIntBig>(nXOff) +
1706 85 : nBlockXSize - nRasterXSize);
1707 : }
1708 3734 : if (nYOff < 0)
1709 : {
1710 8 : bPartialTile = true;
1711 8 : iYOff = -nYOff;
1712 8 : iYCount += nYOff;
1713 : }
1714 3734 : if (nYOff > nRasterYSize - nBlockYSize)
1715 : {
1716 99 : bPartialTile = true;
1717 99 : iYCount -= static_cast<int>(static_cast<GIntBig>(nYOff) +
1718 99 : nBlockYSize - nRasterYSize);
1719 : }
1720 3734 : CPLAssert(iXOff >= 0);
1721 3734 : CPLAssert(iYOff >= 0);
1722 3734 : CPLAssert(iXCount > 0);
1723 3734 : CPLAssert(iYCount > 0);
1724 3734 : CPLAssert(iXOff + iXCount <= nBlockXSize);
1725 3734 : CPLAssert(iYOff + iYCount <= nBlockYSize);
1726 : }
1727 :
1728 4614 : m_asCachedTilesDesc[0].nRow = -1;
1729 4614 : m_asCachedTilesDesc[0].nCol = -1;
1730 4614 : m_asCachedTilesDesc[0].nIdxWithinTileData = -1;
1731 4614 : m_asCachedTilesDesc[0].abBandDirty[0] = false;
1732 4614 : m_asCachedTilesDesc[0].abBandDirty[1] = false;
1733 4614 : m_asCachedTilesDesc[0].abBandDirty[2] = false;
1734 4614 : m_asCachedTilesDesc[0].abBandDirty[3] = false;
1735 :
1736 4614 : CPLErr eErr = CE_Failure;
1737 :
1738 4614 : int bHasNoData = FALSE;
1739 4614 : double dfNoDataValue = IGetRasterBand(1)->GetNoDataValue(&bHasNoData);
1740 4614 : const bool bHasNanNoData = bHasNoData && std::isnan(dfNoDataValue);
1741 :
1742 4614 : 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 4614 : 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 3744 : else if (bAllDirty && m_eDT == GDT_Byte && m_poCT == nullptr &&
1773 3649 : (!bHasNoData || dfNoDataValue == 0.0))
1774 : {
1775 3649 : bool bAllEmpty = true;
1776 3649 : const auto nPixels =
1777 3649 : static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize * nBands;
1778 1207860 : for (GPtrDiff_t i = 0; i < nPixels; i++)
1779 : {
1780 1207840 : if (m_pabyCachedTiles[i] != 0)
1781 : {
1782 3632 : bAllEmpty = false;
1783 3632 : break;
1784 : }
1785 : }
1786 3649 : 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 3632 : }
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 4157 : 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 8314 : VSIMemGenerateHiddenFilename("gpkg_write_tile"));
1836 4157 : const char *pszDriverName = "PNG";
1837 4157 : CPL_IGNORE_RET_VAL(pszDriverName); // Make CSA happy
1838 4157 : bool bTileDriverSupports1Band = false;
1839 4157 : bool bTileDriverSupports2Bands = false;
1840 4157 : bool bTileDriverSupports4Bands = false;
1841 4157 : bool bTileDriverSupportsCT = false;
1842 :
1843 4157 : if (nBands == 1 && m_eDT == GDT_Byte)
1844 72 : IGetRasterBand(1)->GetColorTable();
1845 :
1846 4157 : GDALDataType eTileDT = GDT_Byte;
1847 : // If not all bands are dirty, then (temporarily) use a lossless format
1848 4157 : if (m_eTF == GPKG_TF_PNG_JPEG || (!bAllDirty && m_eTF == GPKG_TF_JPEG))
1849 : {
1850 96 : bTileDriverSupports1Band = true;
1851 96 : if (bPartialTile || !bAllDirty || (nBands == 2 && !bAllOpaque) ||
1852 19 : (nBands == 4 && !bAllOpaque) || m_poCT != nullptr)
1853 : {
1854 88 : pszDriverName = "PNG";
1855 88 : bTileDriverSupports2Bands = m_bPNGSupports2Bands;
1856 88 : bTileDriverSupports4Bands = true;
1857 88 : bTileDriverSupportsCT = m_bPNGSupportsCT;
1858 : }
1859 : else
1860 8 : pszDriverName = "JPEG";
1861 : }
1862 4061 : else if (m_eTF == GPKG_TF_PNG || m_eTF == GPKG_TF_PNG8)
1863 : {
1864 3841 : pszDriverName = "PNG";
1865 3841 : bTileDriverSupports1Band = true;
1866 3841 : bTileDriverSupports2Bands = m_bPNGSupports2Bands;
1867 3841 : bTileDriverSupports4Bands = true;
1868 3841 : bTileDriverSupportsCT = m_bPNGSupportsCT;
1869 : }
1870 220 : else if (m_eTF == GPKG_TF_JPEG)
1871 : {
1872 82 : pszDriverName = "JPEG";
1873 82 : bTileDriverSupports1Band = true;
1874 : }
1875 138 : else if (m_eTF == GPKG_TF_WEBP)
1876 : {
1877 83 : pszDriverName = "WEBP";
1878 83 : bTileDriverSupports4Bands = WEBPSupports4Bands();
1879 : }
1880 55 : else if (m_eTF == GPKG_TF_PNG_16BIT)
1881 : {
1882 51 : pszDriverName = "PNG";
1883 51 : eTileDT = GDT_UInt16;
1884 51 : bTileDriverSupports1Band = true;
1885 : }
1886 4 : else if (m_eTF == GPKG_TF_TIFF_32BIT_FLOAT)
1887 : {
1888 4 : pszDriverName = "GTiff";
1889 4 : eTileDT = GDT_Float32;
1890 4 : bTileDriverSupports1Band = true;
1891 : }
1892 : else
1893 : {
1894 0 : CPLAssert(false);
1895 : }
1896 :
1897 : GDALDriver *l_poDriver =
1898 4157 : GDALDriver::FromHandle(GDALGetDriverByName(pszDriverName));
1899 4157 : if (l_poDriver != nullptr)
1900 : {
1901 4156 : auto poMEMDS = MEMDataset::Create("", nBlockXSize, nBlockYSize, 0,
1902 : eTileDT, nullptr);
1903 4156 : int nTileBands = nBands;
1904 4156 : if (bPartialTile && nBands == 1 && m_poCT == nullptr &&
1905 : bTileDriverSupports2Bands)
1906 27 : nTileBands = 2;
1907 4129 : else if (bPartialTile && bTileDriverSupports4Bands)
1908 39 : nTileBands = 4;
1909 : // only use (somewhat lossy) PNG8 if all bands are dirty
1910 4090 : else if (bAllDirty && m_eTF == GPKG_TF_PNG8 && nBands >= 3 &&
1911 10 : bAllOpaque && !bPartialTile)
1912 10 : nTileBands = 1;
1913 4080 : else if (nBands == 2)
1914 : {
1915 383 : if (bAllOpaque)
1916 : {
1917 51 : if (bTileDriverSupports2Bands)
1918 30 : nTileBands = 1;
1919 : else
1920 21 : nTileBands = 3;
1921 : }
1922 332 : else if (!bTileDriverSupports2Bands)
1923 : {
1924 135 : if (bTileDriverSupports4Bands)
1925 68 : nTileBands = 4;
1926 : else
1927 67 : nTileBands = 3;
1928 : }
1929 : }
1930 3697 : else if (nBands == 4 && (bAllOpaque || !bTileDriverSupports4Bands))
1931 21 : nTileBands = 3;
1932 3676 : else if (nBands == 1 && m_poCT != nullptr && !bTileDriverSupportsCT)
1933 : {
1934 1 : nTileBands = 3;
1935 1 : if (bTileDriverSupports4Bands)
1936 : {
1937 0 : for (int i = 0; i < m_poCT->GetColorEntryCount(); i++)
1938 : {
1939 0 : const GDALColorEntry *psEntry = m_poCT->GetColorEntry(i);
1940 0 : if (psEntry->c4 == 0)
1941 : {
1942 0 : nTileBands = 4;
1943 0 : break;
1944 : }
1945 : }
1946 1 : }
1947 : }
1948 3675 : else if (nBands == 1 && m_poCT == nullptr && !bTileDriverSupports1Band)
1949 1 : nTileBands = 3;
1950 :
1951 4156 : if (bPartialTile && (nTileBands == 2 || nTileBands == 4))
1952 : {
1953 66 : int nTargetAlphaBand = nTileBands;
1954 66 : memset(m_pabyCachedTiles + (nTargetAlphaBand - 1) * nBandBlockSize,
1955 : 0, nBandBlockSize);
1956 4065 : for (GPtrDiff_t iY = iYOff; iY < iYOff + iYCount; iY++)
1957 : {
1958 3999 : memset(m_pabyCachedTiles +
1959 3999 : (static_cast<size_t>(nTargetAlphaBand - 1) *
1960 3999 : nBlockYSize +
1961 3999 : iY) *
1962 3999 : nBlockXSize +
1963 3999 : iXOff,
1964 : 255, iXCount);
1965 : }
1966 : }
1967 :
1968 4156 : GUInt16 *pTempTileBuffer = nullptr;
1969 4156 : GPtrDiff_t nValidPixels = 0;
1970 4156 : double dfTileMin = 0.0;
1971 4156 : double dfTileMax = 0.0;
1972 4156 : double dfTileMean = 0.0;
1973 4156 : double dfTileStdDev = 0.0;
1974 4156 : double dfTileOffset = 0.0;
1975 4156 : double dfTileScale = 1.0;
1976 4156 : if (m_eTF == GPKG_TF_PNG_16BIT)
1977 : {
1978 : pTempTileBuffer = static_cast<GUInt16 *>(
1979 51 : VSI_MALLOC3_VERBOSE(2, nBlockXSize, nBlockYSize));
1980 :
1981 51 : if (m_eDT == GDT_Int16)
1982 : {
1983 72 : ProcessInt16UInt16Tile<GInt16>(
1984 36 : m_pabyCachedTiles,
1985 36 : static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize, true,
1986 36 : CPL_TO_BOOL(bHasNoData), dfNoDataValue, m_usGPKGNull,
1987 : m_dfOffset, m_dfScale, pTempTileBuffer, dfTileOffset,
1988 : dfTileScale, dfTileMin, dfTileMax, dfTileMean, dfTileStdDev,
1989 : nValidPixels);
1990 : }
1991 15 : else if (m_eDT == GDT_UInt16)
1992 : {
1993 20 : ProcessInt16UInt16Tile<GUInt16>(
1994 10 : m_pabyCachedTiles,
1995 10 : static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize, false,
1996 10 : CPL_TO_BOOL(bHasNoData), dfNoDataValue, m_usGPKGNull,
1997 : m_dfOffset, m_dfScale, pTempTileBuffer, dfTileOffset,
1998 : dfTileScale, dfTileMin, dfTileMax, dfTileMean, dfTileStdDev,
1999 : nValidPixels);
2000 : }
2001 5 : else if (m_eDT == GDT_Float32)
2002 : {
2003 5 : const float *pSrc =
2004 : reinterpret_cast<float *>(m_pabyCachedTiles);
2005 5 : float fMin = 0.0f;
2006 5 : float fMax = 0.0f;
2007 5 : double dfM2 = 0.0;
2008 327685 : for (GPtrDiff_t i = 0;
2009 327685 : i < static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize;
2010 : i++)
2011 : {
2012 327680 : const float fVal = pSrc[i];
2013 327680 : if (bHasNanNoData)
2014 : {
2015 0 : if (std::isnan(fVal))
2016 196206 : continue;
2017 : }
2018 720494 : else if (bHasNoData &&
2019 327680 : fVal == static_cast<float>(dfNoDataValue))
2020 : {
2021 196206 : continue;
2022 : }
2023 131474 : if (std::isinf(fVal))
2024 0 : continue;
2025 :
2026 131474 : if (nValidPixels == 0)
2027 : {
2028 5 : fMin = fVal;
2029 5 : fMax = fVal;
2030 : }
2031 : else
2032 : {
2033 131469 : fMin = std::min(fMin, fVal);
2034 131469 : fMax = std::max(fMax, fVal);
2035 : }
2036 131474 : nValidPixels++;
2037 131474 : const double dfDelta = fVal - dfTileMean;
2038 131474 : dfTileMean += dfDelta / nValidPixels;
2039 131474 : dfM2 += dfDelta * (fVal - dfTileMean);
2040 : }
2041 5 : dfTileMin = fMin;
2042 5 : dfTileMax = fMax;
2043 5 : if (nValidPixels)
2044 5 : dfTileStdDev = sqrt(dfM2 / nValidPixels);
2045 :
2046 5 : double dfGlobalMin = (fMin - m_dfOffset) / m_dfScale;
2047 5 : double dfGlobalMax = (fMax - m_dfOffset) / m_dfScale;
2048 5 : if (dfGlobalMax > dfGlobalMin)
2049 : {
2050 4 : if (bHasNoData && m_usGPKGNull == 65535 &&
2051 2 : dfGlobalMax - dfGlobalMin >= 65534.0)
2052 : {
2053 1 : dfTileOffset = dfGlobalMin;
2054 1 : dfTileScale = (dfGlobalMax - dfGlobalMin) / 65534.0;
2055 : }
2056 3 : else if (bHasNoData && m_usGPKGNull == 0 &&
2057 0 : (dfNoDataValue - m_dfOffset) / m_dfScale != 0)
2058 : {
2059 0 : dfTileOffset =
2060 0 : (65535.0 * dfGlobalMin - dfGlobalMax) / 65534.0;
2061 0 : dfTileScale = dfGlobalMin - dfTileOffset;
2062 : }
2063 : else
2064 : {
2065 3 : dfTileOffset = dfGlobalMin;
2066 3 : dfTileScale = (dfGlobalMax - dfGlobalMin) / 65535.0;
2067 : }
2068 : }
2069 : else
2070 : {
2071 1 : dfTileOffset = dfGlobalMin;
2072 : }
2073 :
2074 327685 : for (GPtrDiff_t i = 0;
2075 327685 : i < static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize;
2076 : i++)
2077 : {
2078 327680 : const float fVal = pSrc[i];
2079 327680 : if (bHasNanNoData)
2080 : {
2081 0 : if (std::isnan(fVal))
2082 : {
2083 0 : pTempTileBuffer[i] = m_usGPKGNull;
2084 0 : continue;
2085 : }
2086 : }
2087 327680 : else if (bHasNoData)
2088 : {
2089 196608 : if (fVal == static_cast<float>(dfNoDataValue))
2090 : {
2091 196206 : pTempTileBuffer[i] = m_usGPKGNull;
2092 196206 : continue;
2093 : }
2094 : }
2095 : double dfVal =
2096 131474 : std::isfinite(fVal)
2097 131474 : ? ((fVal - m_dfOffset) / m_dfScale - dfTileOffset) /
2098 : dfTileScale
2099 : : (fVal > 0) ? 65535
2100 131474 : : 0;
2101 131474 : CPLAssert(dfVal >= 0.0 && dfVal < 65535.5);
2102 131474 : pTempTileBuffer[i] = static_cast<GUInt16>(dfVal + 0.5);
2103 131474 : if (bHasNoData && pTempTileBuffer[i] == m_usGPKGNull)
2104 : {
2105 1 : if (m_usGPKGNull > 0)
2106 1 : pTempTileBuffer[i]--;
2107 : else
2108 0 : pTempTileBuffer[i]++;
2109 : }
2110 : }
2111 : }
2112 :
2113 51 : auto hBand = MEMCreateRasterBandEx(
2114 : poMEMDS, 1, reinterpret_cast<GByte *>(pTempTileBuffer),
2115 : GDT_UInt16, 0, 0, false);
2116 51 : poMEMDS->AddMEMBand(hBand);
2117 : }
2118 4105 : else if (m_eTF == GPKG_TF_TIFF_32BIT_FLOAT)
2119 : {
2120 4 : const float *pSrc = reinterpret_cast<float *>(m_pabyCachedTiles);
2121 4 : float fMin = 0.0f;
2122 4 : float fMax = 0.0f;
2123 4 : double dfM2 = 0.0;
2124 262148 : for (GPtrDiff_t i = 0;
2125 262148 : i < static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize; i++)
2126 : {
2127 262144 : const float fVal = pSrc[i];
2128 262144 : if (bHasNanNoData)
2129 : {
2130 0 : if (std::isnan(fVal))
2131 65136 : continue;
2132 : }
2133 392816 : else if (bHasNoData &&
2134 262144 : fVal == static_cast<float>(dfNoDataValue))
2135 : {
2136 65136 : continue;
2137 : }
2138 :
2139 197008 : if (nValidPixels == 0)
2140 : {
2141 4 : fMin = fVal;
2142 4 : fMax = fVal;
2143 : }
2144 : else
2145 : {
2146 197004 : fMin = std::min(fMin, fVal);
2147 197004 : fMax = std::max(fMax, fVal);
2148 : }
2149 197008 : nValidPixels++;
2150 197008 : const double dfDelta = fVal - dfTileMean;
2151 197008 : dfTileMean += dfDelta / nValidPixels;
2152 197008 : dfM2 += dfDelta * (fVal - dfTileMean);
2153 : }
2154 4 : dfTileMin = fMin;
2155 4 : dfTileMax = fMax;
2156 4 : if (nValidPixels)
2157 4 : dfTileStdDev = sqrt(dfM2 / nValidPixels);
2158 :
2159 4 : auto hBand = MEMCreateRasterBandEx(poMEMDS, 1, m_pabyCachedTiles,
2160 : GDT_Float32, 0, 0, false);
2161 4 : poMEMDS->AddMEMBand(hBand);
2162 : }
2163 : else
2164 : {
2165 4101 : CPLAssert(m_eDT == GDT_Byte);
2166 16173 : for (int i = 0; i < nTileBands; i++)
2167 : {
2168 12072 : int iSrc = i;
2169 12072 : if (nBands == 1 && m_poCT == nullptr && nTileBands == 3)
2170 3 : iSrc = 0;
2171 12069 : else if (nBands == 1 && m_poCT == nullptr && bPartialTile &&
2172 : nTileBands == 4)
2173 8 : iSrc = (i < 3) ? 0 : 3;
2174 12061 : else if (nBands == 2 && nTileBands >= 3)
2175 536 : iSrc = (i < 3) ? 0 : 1;
2176 :
2177 24144 : auto hBand = MEMCreateRasterBandEx(
2178 : poMEMDS, i + 1,
2179 12072 : m_pabyCachedTiles + iSrc * nBlockXSize * nBlockYSize,
2180 : GDT_Byte, 0, 0, false);
2181 12072 : poMEMDS->AddMEMBand(hBand);
2182 :
2183 12072 : if (i == 0 && nTileBands == 1 && m_poCT != nullptr)
2184 15 : poMEMDS->GetRasterBand(1)->SetColorTable(m_poCT);
2185 : }
2186 : }
2187 :
2188 4156 : if ((m_eTF == GPKG_TF_PNG_16BIT || m_eTF == GPKG_TF_TIFF_32BIT_FLOAT) &&
2189 55 : nValidPixels == 0)
2190 : {
2191 : // If tile is fully transparent, don't serialize it and remove
2192 : // it if it exists.
2193 1 : GIntBig nId = GetTileId(nRow, nCol);
2194 1 : if (nId > 0)
2195 : {
2196 1 : DeleteTile(nRow, nCol);
2197 :
2198 1 : DeleteFromGriddedTileAncillary(nId);
2199 : }
2200 :
2201 1 : CPLFree(pTempTileBuffer);
2202 1 : delete poMEMDS;
2203 2 : return CE_None;
2204 : }
2205 :
2206 4155 : if (m_eTF == GPKG_TF_PNG8 && nTileBands == 1 && nBands >= 3)
2207 : {
2208 10 : CPLAssert(bAllDirty);
2209 :
2210 10 : auto poMEM_RGB_DS = MEMDataset::Create("", nBlockXSize, nBlockYSize,
2211 : 0, GDT_Byte, nullptr);
2212 40 : for (int i = 0; i < 3; i++)
2213 : {
2214 60 : auto hBand = MEMCreateRasterBandEx(
2215 30 : poMEMDS, i + 1, m_pabyCachedTiles + i * nBandBlockSize,
2216 : GDT_Byte, 0, 0, false);
2217 30 : poMEM_RGB_DS->AddMEMBand(hBand);
2218 : }
2219 :
2220 10 : if (m_pabyHugeColorArray == nullptr)
2221 : {
2222 5 : if (nBlockXSize <= 65536 / nBlockYSize)
2223 4 : m_pabyHugeColorArray =
2224 4 : VSIMalloc(MEDIAN_CUT_AND_DITHER_BUFFER_SIZE_65536);
2225 : else
2226 1 : m_pabyHugeColorArray =
2227 1 : VSIMalloc2(256 * 256 * 256, sizeof(GUInt32));
2228 : }
2229 :
2230 10 : GDALColorTable *poCT = new GDALColorTable();
2231 20 : GDALComputeMedianCutPCTInternal(
2232 10 : poMEM_RGB_DS->GetRasterBand(1), poMEM_RGB_DS->GetRasterBand(2),
2233 10 : poMEM_RGB_DS->GetRasterBand(3),
2234 : /*NULL, NULL, NULL,*/
2235 10 : m_pabyCachedTiles, m_pabyCachedTiles + nBandBlockSize,
2236 10 : m_pabyCachedTiles + 2 * nBandBlockSize, nullptr,
2237 : 256, /* max colors */
2238 : 8, /* bit depth */
2239 : static_cast<GUInt32 *>(
2240 10 : m_pabyHugeColorArray), /* preallocated histogram */
2241 : poCT, nullptr, nullptr);
2242 :
2243 20 : GDALDitherRGB2PCTInternal(
2244 10 : poMEM_RGB_DS->GetRasterBand(1), poMEM_RGB_DS->GetRasterBand(2),
2245 10 : poMEM_RGB_DS->GetRasterBand(3), poMEMDS->GetRasterBand(1), poCT,
2246 : 8, /* bit depth */
2247 : static_cast<GInt16 *>(
2248 10 : m_pabyHugeColorArray), /* pasDynamicColorMap */
2249 10 : m_bDither, nullptr, nullptr);
2250 10 : poMEMDS->GetRasterBand(1)->SetColorTable(poCT);
2251 10 : delete poCT;
2252 10 : GDALClose(poMEM_RGB_DS);
2253 : }
2254 4145 : else if (nBands == 1 && m_poCT != nullptr && nTileBands > 1)
2255 : {
2256 : GByte abyCT[4 * 256];
2257 7 : const int nEntries = std::min(256, m_poCT->GetColorEntryCount());
2258 1799 : for (int i = 0; i < nEntries; i++)
2259 : {
2260 1792 : const GDALColorEntry *psEntry = m_poCT->GetColorEntry(i);
2261 1792 : abyCT[4 * i] = static_cast<GByte>(psEntry->c1);
2262 1792 : abyCT[4 * i + 1] = static_cast<GByte>(psEntry->c2);
2263 1792 : abyCT[4 * i + 2] = static_cast<GByte>(psEntry->c3);
2264 1792 : abyCT[4 * i + 3] = static_cast<GByte>(psEntry->c4);
2265 : }
2266 7 : for (int i = nEntries; i < 256; i++)
2267 : {
2268 0 : abyCT[4 * i] = 0;
2269 0 : abyCT[4 * i + 1] = 0;
2270 0 : abyCT[4 * i + 2] = 0;
2271 0 : abyCT[4 * i + 3] = 0;
2272 : }
2273 7 : if (iYOff > 0)
2274 : {
2275 0 : memset(m_pabyCachedTiles + 0 * nBandBlockSize, 0,
2276 0 : nBlockXSize * iYOff);
2277 0 : memset(m_pabyCachedTiles + 1 * nBandBlockSize, 0,
2278 0 : nBlockXSize * iYOff);
2279 0 : memset(m_pabyCachedTiles + 2 * nBandBlockSize, 0,
2280 0 : nBlockXSize * iYOff);
2281 0 : memset(m_pabyCachedTiles + 3 * nBandBlockSize, 0,
2282 0 : nBlockXSize * iYOff);
2283 : }
2284 1057 : for (GPtrDiff_t iY = iYOff; iY < iYOff + iYCount; iY++)
2285 : {
2286 1050 : if (iXOff > 0)
2287 : {
2288 0 : const GPtrDiff_t i = iY * nBlockXSize;
2289 0 : memset(m_pabyCachedTiles + 0 * nBandBlockSize + i, 0,
2290 : iXOff);
2291 0 : memset(m_pabyCachedTiles + 1 * nBandBlockSize + i, 0,
2292 : iXOff);
2293 0 : memset(m_pabyCachedTiles + 2 * nBandBlockSize + i, 0,
2294 : iXOff);
2295 0 : memset(m_pabyCachedTiles + 3 * nBandBlockSize + i, 0,
2296 : iXOff);
2297 : }
2298 183250 : for (int iX = iXOff; iX < iXOff + iXCount; iX++)
2299 : {
2300 182200 : const GPtrDiff_t i = iY * nBlockXSize + iX;
2301 182200 : GByte byVal = m_pabyCachedTiles[i];
2302 182200 : m_pabyCachedTiles[i] = abyCT[4 * byVal];
2303 182200 : m_pabyCachedTiles[i + 1 * nBandBlockSize] =
2304 182200 : abyCT[4 * byVal + 1];
2305 182200 : m_pabyCachedTiles[i + 2 * nBandBlockSize] =
2306 182200 : abyCT[4 * byVal + 2];
2307 182200 : m_pabyCachedTiles[i + 3 * nBandBlockSize] =
2308 182200 : abyCT[4 * byVal + 3];
2309 : }
2310 1050 : if (iXOff + iXCount < nBlockXSize)
2311 : {
2312 800 : const GPtrDiff_t i = iY * nBlockXSize + iXOff + iXCount;
2313 800 : memset(m_pabyCachedTiles + 0 * nBandBlockSize + i, 0,
2314 800 : nBlockXSize - (iXOff + iXCount));
2315 800 : memset(m_pabyCachedTiles + 1 * nBandBlockSize + i, 0,
2316 800 : nBlockXSize - (iXOff + iXCount));
2317 800 : memset(m_pabyCachedTiles + 2 * nBandBlockSize + i, 0,
2318 800 : nBlockXSize - (iXOff + iXCount));
2319 800 : memset(m_pabyCachedTiles + 3 * nBandBlockSize + i, 0,
2320 800 : nBlockXSize - (iXOff + iXCount));
2321 : }
2322 : }
2323 7 : if (iYOff + iYCount < nBlockYSize)
2324 : {
2325 7 : const GPtrDiff_t i = (iYOff + iYCount) * nBlockXSize;
2326 7 : memset(m_pabyCachedTiles + 0 * nBandBlockSize + i, 0,
2327 7 : nBlockXSize * (nBlockYSize - (iYOff + iYCount)));
2328 7 : memset(m_pabyCachedTiles + 1 * nBandBlockSize + i, 0,
2329 7 : nBlockXSize * (nBlockYSize - (iYOff + iYCount)));
2330 7 : memset(m_pabyCachedTiles + 2 * nBandBlockSize + i, 0,
2331 7 : nBlockXSize * (nBlockYSize - (iYOff + iYCount)));
2332 7 : memset(m_pabyCachedTiles + 3 * nBandBlockSize + i, 0,
2333 7 : nBlockXSize * (nBlockYSize - (iYOff + iYCount)));
2334 : }
2335 : }
2336 :
2337 : char **papszDriverOptions =
2338 4155 : CSLSetNameValue(nullptr, "_INTERNAL_DATASET", "YES");
2339 4155 : if (EQUAL(pszDriverName, "JPEG") || EQUAL(pszDriverName, "WEBP"))
2340 : {
2341 : // If not all bands are dirty, then use lossless WEBP
2342 172 : if (!bAllDirty && EQUAL(pszDriverName, "WEBP"))
2343 : {
2344 2 : papszDriverOptions =
2345 2 : CSLSetNameValue(papszDriverOptions, "LOSSLESS", "YES");
2346 : }
2347 : else
2348 : {
2349 : papszDriverOptions =
2350 170 : CSLSetNameValue(papszDriverOptions, "QUALITY",
2351 : CPLSPrintf("%d", m_nQuality));
2352 : }
2353 : }
2354 3983 : else if (EQUAL(pszDriverName, "PNG"))
2355 : {
2356 3979 : papszDriverOptions = CSLSetNameValue(papszDriverOptions, "ZLEVEL",
2357 : CPLSPrintf("%d", m_nZLevel));
2358 : }
2359 4 : else if (EQUAL(pszDriverName, "GTiff"))
2360 : {
2361 : papszDriverOptions =
2362 4 : CSLSetNameValue(papszDriverOptions, "COMPRESS", "LZW");
2363 4 : if (nBlockXSize * nBlockYSize <= 512 * 512)
2364 : {
2365 : // If tile is not too big, create it as single-strip TIFF
2366 : papszDriverOptions =
2367 4 : CSLSetNameValue(papszDriverOptions, "BLOCKYSIZE",
2368 : CPLSPrintf("%d", nBlockYSize));
2369 : }
2370 : }
2371 : #ifdef DEBUG
2372 : VSIStatBufL sStat;
2373 4155 : CPLAssert(VSIStatL(osMemFileName, &sStat) != 0);
2374 : #endif
2375 : GDALDataset *poOutDS =
2376 4155 : l_poDriver->CreateCopy(osMemFileName, poMEMDS, FALSE,
2377 : papszDriverOptions, nullptr, nullptr);
2378 4155 : CSLDestroy(papszDriverOptions);
2379 4155 : CPLFree(pTempTileBuffer);
2380 :
2381 4155 : if (poOutDS)
2382 : {
2383 4155 : GDALClose(poOutDS);
2384 4155 : vsi_l_offset nBlobSize = 0;
2385 : GByte *pabyBlob =
2386 4155 : VSIGetMemFileBuffer(osMemFileName, &nBlobSize, TRUE);
2387 :
2388 : /* Create or commit and recreate transaction */
2389 4155 : GDALGPKGMBTilesLikePseudoDataset *poMainDS =
2390 4155 : m_poParentDS ? m_poParentDS : this;
2391 4155 : if (poMainDS->m_nTileInsertionCount == 0)
2392 : {
2393 181 : poMainDS->IStartTransaction();
2394 : }
2395 3974 : else if (poMainDS->m_nTileInsertionCount == 1000)
2396 : {
2397 3 : if (poMainDS->ICommitTransaction() != OGRERR_NONE)
2398 : {
2399 1 : poMainDS->m_nTileInsertionCount = -1;
2400 1 : CPLFree(pabyBlob);
2401 1 : VSIUnlink(osMemFileName);
2402 1 : delete poMEMDS;
2403 1 : return CE_Failure;
2404 : }
2405 2 : poMainDS->IStartTransaction();
2406 2 : poMainDS->m_nTileInsertionCount = 0;
2407 : }
2408 4154 : poMainDS->m_nTileInsertionCount++;
2409 :
2410 : char *pszSQL =
2411 4154 : sqlite3_mprintf("INSERT OR REPLACE INTO \"%w\" "
2412 : "(zoom_level, tile_row, tile_column, "
2413 : "tile_data) VALUES (%d, %d, %d, ?)",
2414 : m_osRasterTable.c_str(), m_nZoomLevel,
2415 4154 : GetRowFromIntoTopConvention(nRow), nCol);
2416 : #ifdef DEBUG_VERBOSE
2417 : CPLDebug("GPKG", "%s", pszSQL);
2418 : #endif
2419 4154 : sqlite3_stmt *hStmt = nullptr;
2420 4154 : int rc = sqlite3_prepare_v2(IGetDB(), pszSQL, -1, &hStmt, nullptr);
2421 4154 : if (rc != SQLITE_OK)
2422 : {
2423 0 : CPLError(CE_Failure, CPLE_AppDefined,
2424 : "failed to prepare SQL %s: %s", pszSQL,
2425 0 : sqlite3_errmsg(IGetDB()));
2426 0 : CPLFree(pabyBlob);
2427 : }
2428 : else
2429 : {
2430 4154 : sqlite3_bind_blob(hStmt, 1, pabyBlob,
2431 : static_cast<int>(nBlobSize), CPLFree);
2432 4154 : rc = sqlite3_step(hStmt);
2433 4154 : if (rc == SQLITE_DONE)
2434 4154 : eErr = CE_None;
2435 : else
2436 : {
2437 0 : CPLError(CE_Failure, CPLE_AppDefined,
2438 : "Failure when inserting tile (row=%d,col=%d) at "
2439 : "zoom_level=%d : %s",
2440 0 : GetRowFromIntoTopConvention(nRow), nCol,
2441 0 : m_nZoomLevel, sqlite3_errmsg(IGetDB()));
2442 : }
2443 : }
2444 4154 : sqlite3_finalize(hStmt);
2445 4154 : sqlite3_free(pszSQL);
2446 :
2447 4154 : if (m_eTF == GPKG_TF_PNG_16BIT || m_eTF == GPKG_TF_TIFF_32BIT_FLOAT)
2448 : {
2449 54 : GIntBig nTileId = GetTileId(nRow, nCol);
2450 54 : if (nTileId == 0)
2451 0 : eErr = CE_Failure;
2452 : else
2453 : {
2454 54 : DeleteFromGriddedTileAncillary(nTileId);
2455 :
2456 54 : pszSQL = sqlite3_mprintf(
2457 : "INSERT INTO gpkg_2d_gridded_tile_ancillary "
2458 : "(tpudt_name, tpudt_id, scale, offset, min, max, "
2459 : "mean, std_dev) VALUES "
2460 : "('%q', ?, %.17g, %.17g, ?, ?, ?, ?)",
2461 : m_osRasterTable.c_str(), dfTileScale, dfTileOffset);
2462 : #ifdef DEBUG_VERBOSE
2463 : CPLDebug("GPKG", "%s", pszSQL);
2464 : #endif
2465 54 : hStmt = nullptr;
2466 54 : rc = sqlite3_prepare_v2(IGetDB(), pszSQL, -1, &hStmt,
2467 : nullptr);
2468 54 : if (rc != SQLITE_OK)
2469 : {
2470 0 : eErr = CE_Failure;
2471 0 : CPLError(CE_Failure, CPLE_AppDefined,
2472 : "failed to prepare SQL %s: %s", pszSQL,
2473 0 : sqlite3_errmsg(IGetDB()));
2474 : }
2475 : else
2476 : {
2477 54 : sqlite3_bind_int64(hStmt, 1, nTileId);
2478 54 : sqlite3_bind_double(hStmt, 2, dfTileMin);
2479 54 : sqlite3_bind_double(hStmt, 3, dfTileMax);
2480 54 : sqlite3_bind_double(hStmt, 4, dfTileMean);
2481 54 : sqlite3_bind_double(hStmt, 5, dfTileStdDev);
2482 54 : rc = sqlite3_step(hStmt);
2483 54 : if (rc == SQLITE_DONE)
2484 : {
2485 54 : eErr = CE_None;
2486 : }
2487 : else
2488 : {
2489 0 : CPLError(CE_Failure, CPLE_AppDefined,
2490 : "Cannot insert into "
2491 : "gpkg_2d_gridded_tile_ancillary");
2492 0 : eErr = CE_Failure;
2493 : }
2494 : }
2495 54 : sqlite3_finalize(hStmt);
2496 54 : sqlite3_free(pszSQL);
2497 : }
2498 : }
2499 : }
2500 :
2501 4154 : VSIUnlink(osMemFileName);
2502 4154 : delete poMEMDS;
2503 : }
2504 : else
2505 : {
2506 1 : CPLError(CE_Failure, CPLE_NotSupported, "Cannot find driver %s",
2507 : pszDriverName);
2508 : }
2509 :
2510 4155 : return eErr;
2511 : }
2512 :
2513 : /************************************************************************/
2514 : /* FlushRemainingShiftedTiles() */
2515 : /************************************************************************/
2516 :
2517 : CPLErr
2518 609 : GDALGPKGMBTilesLikePseudoDataset::FlushRemainingShiftedTiles(bool bPartialFlush)
2519 : {
2520 609 : if (m_hTempDB == nullptr)
2521 460 : return CE_None;
2522 :
2523 745 : for (int i = 0; i <= 3; i++)
2524 : {
2525 596 : m_asCachedTilesDesc[i].nRow = -1;
2526 596 : m_asCachedTilesDesc[i].nCol = -1;
2527 596 : m_asCachedTilesDesc[i].nIdxWithinTileData = -1;
2528 : }
2529 :
2530 : int nBlockXSize, nBlockYSize;
2531 149 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
2532 149 : const int nBands = IGetRasterCount();
2533 149 : const int nRasterXSize = IGetRasterBand(1)->GetXSize();
2534 149 : const int nRasterYSize = IGetRasterBand(1)->GetYSize();
2535 149 : const int nXBlocks = DIV_ROUND_UP(nRasterXSize, nBlockXSize);
2536 149 : const int nYBlocks = DIV_ROUND_UP(nRasterYSize, nBlockYSize);
2537 :
2538 149 : int nPartialActiveTiles = 0;
2539 149 : if (bPartialFlush)
2540 : {
2541 2 : sqlite3_stmt *hStmt = nullptr;
2542 4 : CPLString osSQL;
2543 : osSQL.Printf("SELECT COUNT(*) FROM partial_tiles WHERE zoom_level = %d "
2544 : "AND partial_flag != 0",
2545 2 : m_nZoomLevel);
2546 2 : if (sqlite3_prepare_v2(m_hTempDB, osSQL.c_str(), -1, &hStmt, nullptr) ==
2547 : SQLITE_OK)
2548 : {
2549 2 : if (sqlite3_step(hStmt) == SQLITE_ROW)
2550 : {
2551 2 : nPartialActiveTiles = sqlite3_column_int(hStmt, 0);
2552 2 : CPLDebug("GPKG", "Active partial tiles before flush: %d",
2553 : nPartialActiveTiles);
2554 : }
2555 2 : sqlite3_finalize(hStmt);
2556 : }
2557 : }
2558 :
2559 298 : CPLString osSQL = "SELECT tile_row, tile_column, partial_flag";
2560 616 : for (int nBand = 1; nBand <= nBands; nBand++)
2561 : {
2562 467 : osSQL += CPLSPrintf(", tile_data_band_%d", nBand);
2563 : }
2564 : osSQL += CPLSPrintf(" FROM partial_tiles WHERE "
2565 : "zoom_level = %d AND partial_flag != 0",
2566 149 : m_nZoomLevel);
2567 149 : if (bPartialFlush)
2568 : {
2569 2 : osSQL += " ORDER BY age";
2570 : }
2571 149 : const char *pszSQL = osSQL.c_str();
2572 :
2573 : #ifdef DEBUG_VERBOSE
2574 : CPLDebug("GPKG", "%s", pszSQL);
2575 : #endif
2576 149 : sqlite3_stmt *hStmt = nullptr;
2577 149 : int rc = sqlite3_prepare_v2(m_hTempDB, pszSQL, -1, &hStmt, nullptr);
2578 149 : if (rc != SQLITE_OK)
2579 : {
2580 0 : CPLError(CE_Failure, CPLE_AppDefined,
2581 : "sqlite3_prepare_v2(%s) failed: %s", pszSQL,
2582 : sqlite3_errmsg(m_hTempDB));
2583 0 : return CE_Failure;
2584 : }
2585 :
2586 149 : CPLErr eErr = CE_None;
2587 149 : bool bGotPartialTiles = false;
2588 149 : int nCountFlushedTiles = 0;
2589 149 : const size_t nBandBlockSize =
2590 149 : static_cast<size_t>(nBlockXSize) * nBlockYSize * m_nDTSize;
2591 175 : do
2592 : {
2593 324 : rc = sqlite3_step(hStmt);
2594 324 : if (rc == SQLITE_ROW)
2595 : {
2596 177 : bGotPartialTiles = true;
2597 :
2598 177 : int nRow = sqlite3_column_int(hStmt, 0);
2599 177 : int nCol = sqlite3_column_int(hStmt, 1);
2600 177 : int nPartialFlags = sqlite3_column_int(hStmt, 2);
2601 :
2602 177 : if (bPartialFlush)
2603 : {
2604 : // This method assumes that there are no dirty blocks alive
2605 : // so check this assumption.
2606 : // When called with bPartialFlush = false, FlushCache() has
2607 : // already been called, so no need to check.
2608 16 : bool bFoundDirtyBlock = false;
2609 16 : int nBlockXOff = nCol - m_nShiftXTiles;
2610 16 : int nBlockYOff = nRow - m_nShiftYTiles;
2611 96 : for (int iX = 0; !bFoundDirtyBlock &&
2612 48 : iX < ((m_nShiftXPixelsMod != 0) ? 2 : 1);
2613 : iX++)
2614 : {
2615 32 : if (nBlockXOff + iX < 0 || nBlockXOff + iX >= nXBlocks)
2616 12 : continue;
2617 120 : for (int iY = 0; !bFoundDirtyBlock &&
2618 60 : iY < ((m_nShiftYPixelsMod != 0) ? 2 : 1);
2619 : iY++)
2620 : {
2621 40 : if (nBlockYOff + iY < 0 || nBlockYOff + iY >= nYBlocks)
2622 17 : continue;
2623 69 : for (int iBand = 1;
2624 69 : !bFoundDirtyBlock && iBand <= nBands; iBand++)
2625 : {
2626 : GDALRasterBlock *poBlock =
2627 : cpl::down_cast<GDALGPKGMBTilesLikeRasterBand *>(
2628 46 : IGetRasterBand(iBand))
2629 46 : ->AccessibleTryGetLockedBlockRef(
2630 : nBlockXOff + iX, nBlockYOff + iY);
2631 46 : if (poBlock)
2632 : {
2633 0 : if (poBlock->GetDirty())
2634 0 : bFoundDirtyBlock = true;
2635 0 : poBlock->DropLock();
2636 : }
2637 : }
2638 : }
2639 : }
2640 16 : if (bFoundDirtyBlock)
2641 : {
2642 : #ifdef DEBUG_VERBOSE
2643 : CPLDebug("GPKG",
2644 : "Skipped flushing tile row = %d, column = %d "
2645 : "because it has dirty block(s) in GDAL cache",
2646 : nRow, nCol);
2647 : #endif
2648 0 : continue;
2649 : }
2650 : }
2651 :
2652 177 : nCountFlushedTiles++;
2653 177 : if (bPartialFlush && nCountFlushedTiles >= nPartialActiveTiles / 2)
2654 : {
2655 2 : CPLDebug("GPKG", "Flushed %d tiles", nCountFlushedTiles);
2656 2 : break;
2657 : }
2658 :
2659 705 : for (int nBand = 1; nBand <= nBands; nBand++)
2660 : {
2661 530 : if (nPartialFlags & (((1 << 4) - 1) << (4 * (nBand - 1))))
2662 : {
2663 498 : CPLAssert(sqlite3_column_bytes(hStmt, 2 + nBand) ==
2664 : static_cast<int>(nBandBlockSize));
2665 498 : memcpy(m_pabyCachedTiles + (nBand - 1) * nBandBlockSize,
2666 : sqlite3_column_blob(hStmt, 2 + nBand),
2667 : nBandBlockSize);
2668 : }
2669 : else
2670 : {
2671 32 : FillEmptyTileSingleBand(m_pabyCachedTiles +
2672 32 : (nBand - 1) * nBandBlockSize);
2673 : }
2674 : }
2675 :
2676 175 : int nFullFlags = (1 << (4 * nBands)) - 1;
2677 :
2678 : // In case the partial flags indicate that there's some quadrant
2679 : // missing, check in the main database if there is already a tile
2680 : // In which case, use the parts of that tile that aren't in the
2681 : // temporary database
2682 175 : if (nPartialFlags != nFullFlags)
2683 : {
2684 525 : char *pszNewSQL = sqlite3_mprintf(
2685 : "SELECT tile_data%s FROM \"%w\" "
2686 : "WHERE zoom_level = %d AND tile_row = %d AND tile_column = "
2687 : "%d%s",
2688 175 : m_eDT != GDT_Byte ? ", id"
2689 : : "", // MBTiles do not have an id
2690 : m_osRasterTable.c_str(), m_nZoomLevel,
2691 175 : GetRowFromIntoTopConvention(nRow), nCol,
2692 175 : !m_osWHERE.empty()
2693 0 : ? CPLSPrintf(" AND (%s)", m_osWHERE.c_str())
2694 : : "");
2695 : #ifdef DEBUG_VERBOSE
2696 : CPLDebug("GPKG", "%s", pszNewSQL);
2697 : #endif
2698 175 : sqlite3_stmt *hNewStmt = nullptr;
2699 175 : rc = sqlite3_prepare_v2(IGetDB(), pszNewSQL, -1, &hNewStmt,
2700 : nullptr);
2701 175 : if (rc == SQLITE_OK)
2702 : {
2703 175 : rc = sqlite3_step(hNewStmt);
2704 192 : if (rc == SQLITE_ROW &&
2705 17 : sqlite3_column_type(hNewStmt, 0) == SQLITE_BLOB)
2706 : {
2707 17 : const int nBytes = sqlite3_column_bytes(hNewStmt, 0);
2708 : GIntBig nTileId =
2709 17 : (m_eDT == GDT_Byte)
2710 17 : ? 0
2711 0 : : sqlite3_column_int64(hNewStmt, 1);
2712 : GByte *pabyRawData =
2713 : const_cast<GByte *>(static_cast<const GByte *>(
2714 17 : sqlite3_column_blob(hNewStmt, 0)));
2715 : const CPLString osMemFileName(
2716 34 : VSIMemGenerateHiddenFilename("gpkg_read_tile"));
2717 17 : VSILFILE *fp = VSIFileFromMemBuffer(
2718 : osMemFileName.c_str(), pabyRawData, nBytes, FALSE);
2719 17 : VSIFCloseL(fp);
2720 :
2721 17 : double dfTileOffset = 0.0;
2722 17 : double dfTileScale = 1.0;
2723 17 : GetTileOffsetAndScale(nTileId, dfTileOffset,
2724 : dfTileScale);
2725 17 : const int nTileBands = m_eDT == GDT_Byte ? 4 : 1;
2726 17 : GByte *pabyTemp =
2727 17 : m_pabyCachedTiles + nTileBands * nBandBlockSize;
2728 17 : ReadTile(osMemFileName, pabyTemp, dfTileOffset,
2729 : dfTileScale);
2730 17 : VSIUnlink(osMemFileName);
2731 :
2732 17 : int iYQuadrantMax = (m_nShiftYPixelsMod) ? 1 : 0;
2733 17 : int iXQuadrantMax = (m_nShiftXPixelsMod) ? 1 : 0;
2734 51 : for (int iYQuadrant = 0; iYQuadrant <= iYQuadrantMax;
2735 : iYQuadrant++)
2736 : {
2737 100 : for (int iXQuadrant = 0;
2738 100 : iXQuadrant <= iXQuadrantMax; iXQuadrant++)
2739 : {
2740 226 : for (int nBand = 1; nBand <= nBands; nBand++)
2741 : {
2742 160 : int iQuadrantFlag = 0;
2743 160 : if (iXQuadrant == 0 && iYQuadrant == 0)
2744 42 : iQuadrantFlag |= (1 << 0);
2745 160 : if (iXQuadrant == iXQuadrantMax &&
2746 : iYQuadrant == 0)
2747 42 : iQuadrantFlag |= (1 << 1);
2748 160 : if (iXQuadrant == 0 &&
2749 : iYQuadrant == iYQuadrantMax)
2750 42 : iQuadrantFlag |= (1 << 2);
2751 160 : if (iXQuadrant == iXQuadrantMax &&
2752 : iYQuadrant == iYQuadrantMax)
2753 42 : iQuadrantFlag |= (1 << 3);
2754 160 : int nLocalFlag = iQuadrantFlag
2755 160 : << (4 * (nBand - 1));
2756 160 : if (!(nPartialFlags & nLocalFlag))
2757 : {
2758 : int nXOff, nYOff, nXSize, nYSize;
2759 118 : if (iXQuadrant == 0 &&
2760 60 : m_nShiftXPixelsMod != 0)
2761 : {
2762 56 : nXOff = 0;
2763 56 : nXSize = m_nShiftXPixelsMod;
2764 : }
2765 : else
2766 : {
2767 62 : nXOff = m_nShiftXPixelsMod;
2768 62 : nXSize = nBlockXSize -
2769 62 : m_nShiftXPixelsMod;
2770 : }
2771 118 : if (iYQuadrant == 0 &&
2772 63 : m_nShiftYPixelsMod != 0)
2773 : {
2774 63 : nYOff = 0;
2775 63 : nYSize = m_nShiftYPixelsMod;
2776 : }
2777 : else
2778 : {
2779 55 : nYOff = m_nShiftYPixelsMod;
2780 55 : nYSize = nBlockYSize -
2781 55 : m_nShiftYPixelsMod;
2782 : }
2783 118 : for (int iY = nYOff;
2784 13888 : iY < nYOff + nYSize; iY++)
2785 : {
2786 13770 : memcpy(m_pabyCachedTiles +
2787 13770 : ((static_cast<size_t>(
2788 13770 : nBand - 1) *
2789 13770 : nBlockYSize +
2790 13770 : iY) *
2791 13770 : nBlockXSize +
2792 13770 : nXOff) *
2793 13770 : m_nDTSize,
2794 13770 : pabyTemp +
2795 13770 : ((static_cast<size_t>(
2796 13770 : nBand - 1) *
2797 13770 : nBlockYSize +
2798 13770 : iY) *
2799 13770 : nBlockXSize +
2800 13770 : nXOff) *
2801 13770 : m_nDTSize,
2802 13770 : static_cast<size_t>(nXSize) *
2803 13770 : m_nDTSize);
2804 : }
2805 : }
2806 : }
2807 : }
2808 : }
2809 : }
2810 158 : else if (rc != SQLITE_DONE)
2811 : {
2812 0 : CPLError(CE_Failure, CPLE_AppDefined,
2813 : "sqlite3_step(%s) failed: %s", pszNewSQL,
2814 : sqlite3_errmsg(m_hTempDB));
2815 : }
2816 175 : sqlite3_finalize(hNewStmt);
2817 : }
2818 : else
2819 : {
2820 0 : CPLError(CE_Failure, CPLE_AppDefined,
2821 : "sqlite3_prepare_v2(%s) failed: %s", pszNewSQL,
2822 : sqlite3_errmsg(m_hTempDB));
2823 : }
2824 175 : sqlite3_free(pszNewSQL);
2825 : }
2826 :
2827 175 : m_asCachedTilesDesc[0].nRow = nRow;
2828 175 : m_asCachedTilesDesc[0].nCol = nCol;
2829 175 : m_asCachedTilesDesc[0].nIdxWithinTileData = 0;
2830 175 : m_asCachedTilesDesc[0].abBandDirty[0] = true;
2831 175 : m_asCachedTilesDesc[0].abBandDirty[1] = true;
2832 175 : m_asCachedTilesDesc[0].abBandDirty[2] = true;
2833 175 : m_asCachedTilesDesc[0].abBandDirty[3] = true;
2834 :
2835 175 : eErr = WriteTile();
2836 :
2837 175 : if (eErr == CE_None && bPartialFlush)
2838 : {
2839 : pszSQL =
2840 14 : CPLSPrintf("DELETE FROM partial_tiles WHERE zoom_level = "
2841 : "%d AND tile_row = %d AND tile_column = %d",
2842 : m_nZoomLevel, nRow, nCol);
2843 : #ifdef DEBUG_VERBOSE
2844 : CPLDebug("GPKG", "%s", pszSQL);
2845 : #endif
2846 14 : if (SQLCommand(m_hTempDB, pszSQL) != OGRERR_NONE)
2847 0 : eErr = CE_None;
2848 : }
2849 : }
2850 : else
2851 : {
2852 147 : if (rc != SQLITE_DONE)
2853 : {
2854 0 : CPLError(CE_Failure, CPLE_AppDefined,
2855 : "sqlite3_step(%s) failed: %s", pszSQL,
2856 : sqlite3_errmsg(m_hTempDB));
2857 : }
2858 147 : break;
2859 : }
2860 175 : } while (eErr == CE_None);
2861 :
2862 149 : sqlite3_finalize(hStmt);
2863 :
2864 149 : if (bPartialFlush && nCountFlushedTiles < nPartialActiveTiles / 2)
2865 : {
2866 0 : CPLDebug("GPKG", "Flushed %d tiles. Target was %d", nCountFlushedTiles,
2867 : nPartialActiveTiles / 2);
2868 : }
2869 :
2870 149 : if (bGotPartialTiles && !bPartialFlush)
2871 : {
2872 : #ifdef DEBUG_VERBOSE
2873 : pszSQL = CPLSPrintf("SELECT p1.id, p1.tile_row, p1.tile_column FROM "
2874 : "partial_tiles p1, partial_tiles p2 "
2875 : "WHERE p1.zoom_level = %d AND p2.zoom_level = %d "
2876 : "AND p1.tile_row = p2.tile_row AND p1.tile_column "
2877 : "= p2.tile_column AND p2.partial_flag != 0",
2878 : -1 - m_nZoomLevel, m_nZoomLevel);
2879 : rc = sqlite3_prepare_v2(m_hTempDB, pszSQL, -1, &hStmt, nullptr);
2880 : CPLAssert(rc == SQLITE_OK);
2881 : while ((rc = sqlite3_step(hStmt)) == SQLITE_ROW)
2882 : {
2883 : CPLDebug("GPKG",
2884 : "Conflict: existing id = %d, tile_row = %d, tile_column = "
2885 : "%d, zoom_level = %d",
2886 : sqlite3_column_int(hStmt, 0), sqlite3_column_int(hStmt, 1),
2887 : sqlite3_column_int(hStmt, 2), m_nZoomLevel);
2888 : }
2889 : sqlite3_finalize(hStmt);
2890 : #endif
2891 :
2892 110 : pszSQL = CPLSPrintf("UPDATE partial_tiles SET zoom_level = %d, "
2893 : "partial_flag = 0, age = -1 WHERE zoom_level = %d "
2894 : "AND partial_flag != 0",
2895 55 : -1 - m_nZoomLevel, m_nZoomLevel);
2896 : #ifdef DEBUG_VERBOSE
2897 : CPLDebug("GPKG", "%s", pszSQL);
2898 : #endif
2899 55 : SQLCommand(m_hTempDB, pszSQL);
2900 : }
2901 :
2902 149 : return eErr;
2903 : }
2904 :
2905 : /************************************************************************/
2906 : /* DoPartialFlushOfPartialTilesIfNecessary() */
2907 : /************************************************************************/
2908 :
2909 : CPLErr
2910 4556 : GDALGPKGMBTilesLikePseudoDataset::DoPartialFlushOfPartialTilesIfNecessary()
2911 : {
2912 4556 : time_t nCurTimeStamp = time(nullptr);
2913 4556 : if (m_nLastSpaceCheckTimestamp == 0)
2914 3 : m_nLastSpaceCheckTimestamp = nCurTimeStamp;
2915 4556 : if (m_nLastSpaceCheckTimestamp > 0 &&
2916 81 : (m_bForceTempDBCompaction ||
2917 13 : nCurTimeStamp - m_nLastSpaceCheckTimestamp > 10))
2918 : {
2919 68 : m_nLastSpaceCheckTimestamp = nCurTimeStamp;
2920 : GIntBig nFreeSpace =
2921 68 : VSIGetDiskFreeSpace(CPLGetDirnameSafe(m_osTempDBFilename).c_str());
2922 68 : bool bTryFreeing = false;
2923 68 : if (nFreeSpace >= 0 && nFreeSpace < 1024 * 1024 * 1024)
2924 : {
2925 0 : CPLDebug("GPKG",
2926 : "Free space below 1GB. Flushing part of partial tiles");
2927 0 : bTryFreeing = true;
2928 : }
2929 : else
2930 : {
2931 : VSIStatBufL sStat;
2932 68 : if (VSIStatL(m_osTempDBFilename, &sStat) == 0)
2933 : {
2934 68 : GIntBig nTempSpace = sStat.st_size;
2935 136 : if (VSIStatL((m_osTempDBFilename + "-journal").c_str(),
2936 68 : &sStat) == 0)
2937 0 : nTempSpace += sStat.st_size;
2938 136 : else if (VSIStatL((m_osTempDBFilename + "-wal").c_str(),
2939 68 : &sStat) == 0)
2940 0 : nTempSpace += sStat.st_size;
2941 :
2942 : int nBlockXSize, nBlockYSize;
2943 68 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
2944 68 : const int nBands = IGetRasterCount();
2945 :
2946 68 : if (nTempSpace >
2947 68 : 4 * static_cast<GIntBig>(IGetRasterBand(1)->GetXSize()) *
2948 68 : nBlockYSize * nBands * m_nDTSize)
2949 : {
2950 2 : CPLDebug("GPKG",
2951 : "Partial tiles DB is " CPL_FRMT_GIB
2952 : " bytes. Flushing part of partial tiles",
2953 : nTempSpace);
2954 2 : bTryFreeing = true;
2955 : }
2956 : }
2957 : }
2958 68 : if (bTryFreeing)
2959 : {
2960 2 : if (FlushRemainingShiftedTiles(true /* partial flush*/) != CE_None)
2961 : {
2962 0 : return CE_Failure;
2963 : }
2964 2 : SQLCommand(m_hTempDB,
2965 : "DELETE FROM partial_tiles WHERE zoom_level < 0");
2966 2 : SQLCommand(m_hTempDB, "VACUUM");
2967 : }
2968 : }
2969 4556 : return CE_None;
2970 : }
2971 :
2972 : /************************************************************************/
2973 : /* WriteShiftedTile() */
2974 : /************************************************************************/
2975 :
2976 4556 : CPLErr GDALGPKGMBTilesLikePseudoDataset::WriteShiftedTile(
2977 : int nRow, int nCol, int nBand, int nDstXOffset, int nDstYOffset,
2978 : int nDstXSize, int nDstYSize)
2979 : {
2980 4556 : CPLAssert(m_nShiftXPixelsMod || m_nShiftYPixelsMod);
2981 4556 : CPLAssert(nRow >= 0);
2982 4556 : CPLAssert(nCol >= 0);
2983 4556 : CPLAssert(nRow < m_nTileMatrixHeight);
2984 4556 : CPLAssert(nCol < m_nTileMatrixWidth);
2985 :
2986 4556 : if (m_hTempDB == nullptr &&
2987 55 : (m_poParentDS == nullptr || m_poParentDS->m_hTempDB == nullptr))
2988 : {
2989 : const char *pszBaseFilename =
2990 48 : m_poParentDS ? m_poParentDS->IGetFilename() : IGetFilename();
2991 : m_osTempDBFilename =
2992 48 : CPLResetExtensionSafe(pszBaseFilename, "partial_tiles.db");
2993 48 : CPLPushErrorHandler(CPLQuietErrorHandler);
2994 48 : VSIUnlink(m_osTempDBFilename);
2995 48 : CPLPopErrorHandler();
2996 48 : m_hTempDB = nullptr;
2997 48 : int rc = 0;
2998 48 : if (STARTS_WITH(m_osTempDBFilename, "/vsi"))
2999 : {
3000 33 : m_pMyVFS = OGRSQLiteCreateVFS(nullptr, nullptr);
3001 33 : sqlite3_vfs_register(m_pMyVFS, 0);
3002 33 : rc = sqlite3_open_v2(m_osTempDBFilename, &m_hTempDB,
3003 : SQLITE_OPEN_READWRITE | SQLITE_OPEN_CREATE |
3004 : SQLITE_OPEN_NOMUTEX,
3005 33 : m_pMyVFS->zName);
3006 : }
3007 : else
3008 : {
3009 15 : rc = sqlite3_open(m_osTempDBFilename, &m_hTempDB);
3010 : }
3011 48 : if (rc != SQLITE_OK || m_hTempDB == nullptr)
3012 : {
3013 0 : CPLError(CE_Failure, CPLE_AppDefined,
3014 : "Cannot create temporary database %s",
3015 : m_osTempDBFilename.c_str());
3016 0 : return CE_Failure;
3017 : }
3018 48 : SQLCommand(m_hTempDB, "PRAGMA synchronous = OFF");
3019 48 : SQLCommand(m_hTempDB, "PRAGMA journal_mode = OFF");
3020 48 : SQLCommand(m_hTempDB, "CREATE TABLE partial_tiles("
3021 : "id INTEGER PRIMARY KEY AUTOINCREMENT,"
3022 : "zoom_level INTEGER NOT NULL,"
3023 : "tile_column INTEGER NOT NULL,"
3024 : "tile_row INTEGER NOT NULL,"
3025 : "tile_data_band_1 BLOB,"
3026 : "tile_data_band_2 BLOB,"
3027 : "tile_data_band_3 BLOB,"
3028 : "tile_data_band_4 BLOB,"
3029 : "partial_flag INTEGER NOT NULL,"
3030 : "age INTEGER NOT NULL,"
3031 : "UNIQUE (zoom_level, tile_column, tile_row))");
3032 48 : SQLCommand(m_hTempDB, "CREATE INDEX partial_tiles_partial_flag_idx "
3033 : "ON partial_tiles(partial_flag)");
3034 48 : SQLCommand(m_hTempDB, "CREATE INDEX partial_tiles_age_idx "
3035 : "ON partial_tiles(age)");
3036 :
3037 48 : if (m_poParentDS != nullptr)
3038 : {
3039 4 : m_poParentDS->m_osTempDBFilename = m_osTempDBFilename;
3040 4 : m_poParentDS->m_hTempDB = m_hTempDB;
3041 : }
3042 : }
3043 :
3044 4556 : if (m_poParentDS != nullptr)
3045 4020 : m_hTempDB = m_poParentDS->m_hTempDB;
3046 :
3047 : int nBlockXSize, nBlockYSize;
3048 4556 : IGetRasterBand(1)->GetBlockSize(&nBlockXSize, &nBlockYSize);
3049 4556 : const int nBands = IGetRasterCount();
3050 4556 : const size_t nBandBlockSize =
3051 4556 : static_cast<size_t>(nBlockXSize) * nBlockYSize * m_nDTSize;
3052 :
3053 4556 : int iQuadrantFlag = 0;
3054 4556 : if (nDstXOffset == 0 && nDstYOffset == 0)
3055 1031 : iQuadrantFlag |= (1 << 0);
3056 4556 : if ((nDstXOffset != 0 || nDstXOffset + nDstXSize == nBlockXSize) &&
3057 : nDstYOffset == 0)
3058 1125 : iQuadrantFlag |= (1 << 1);
3059 4556 : if (nDstXOffset == 0 &&
3060 1031 : (nDstYOffset != 0 || nDstYOffset + nDstYSize == nBlockYSize))
3061 1135 : iQuadrantFlag |= (1 << 2);
3062 4556 : if ((nDstXOffset != 0 || nDstXOffset + nDstXSize == nBlockXSize) &&
3063 1125 : (nDstYOffset != 0 || nDstYOffset + nDstYSize == nBlockYSize))
3064 1323 : iQuadrantFlag |= (1 << 3);
3065 4556 : int l_nFlags = iQuadrantFlag << (4 * (nBand - 1));
3066 4556 : int nFullFlags = (1 << (4 * nBands)) - 1;
3067 4556 : int nOldFlags = 0;
3068 :
3069 18224 : for (int i = 1; i <= 3; i++)
3070 : {
3071 13668 : m_asCachedTilesDesc[i].nRow = -1;
3072 13668 : m_asCachedTilesDesc[i].nCol = -1;
3073 13668 : m_asCachedTilesDesc[i].nIdxWithinTileData = -1;
3074 : }
3075 :
3076 4556 : int nExistingId = 0;
3077 4556 : const char *pszSQL = CPLSPrintf(
3078 : "SELECT id, partial_flag, tile_data_band_%d FROM partial_tiles WHERE "
3079 : "zoom_level = %d AND tile_row = %d AND tile_column = %d",
3080 : nBand, m_nZoomLevel, nRow, nCol);
3081 : #ifdef DEBUG_VERBOSE
3082 : CPLDebug("GPKG", "%s", pszSQL);
3083 : #endif
3084 4556 : sqlite3_stmt *hStmt = nullptr;
3085 4556 : int rc = sqlite3_prepare_v2(m_hTempDB, pszSQL, -1, &hStmt, nullptr);
3086 4556 : if (rc != SQLITE_OK)
3087 : {
3088 0 : CPLError(CE_Failure, CPLE_AppDefined,
3089 : "sqlite3_prepare_v2(%s) failed: %s", pszSQL,
3090 : sqlite3_errmsg(m_hTempDB));
3091 0 : return CE_Failure;
3092 : }
3093 :
3094 4556 : rc = sqlite3_step(hStmt);
3095 4556 : const int nTileBands = m_eDT == GDT_Byte ? 4 : 1;
3096 4556 : GByte *pabyTemp = m_pabyCachedTiles + nTileBands * nBandBlockSize;
3097 4556 : if (rc == SQLITE_ROW)
3098 : {
3099 4134 : nExistingId = sqlite3_column_int(hStmt, 0);
3100 : #ifdef DEBUG_VERBOSE
3101 : CPLDebug("GPKG", "Using partial_tile id=%d", nExistingId);
3102 : #endif
3103 4134 : nOldFlags = sqlite3_column_int(hStmt, 1);
3104 4134 : CPLAssert(nOldFlags != 0);
3105 4134 : if ((nOldFlags & (((1 << 4) - 1) << (4 * (nBand - 1)))) == 0)
3106 : {
3107 1021 : FillEmptyTileSingleBand(pabyTemp + (nBand - 1) * nBandBlockSize);
3108 : }
3109 : else
3110 : {
3111 3113 : CPLAssert(sqlite3_column_bytes(hStmt, 2) ==
3112 : static_cast<int>(nBandBlockSize));
3113 3113 : memcpy(pabyTemp + (nBand - 1) * nBandBlockSize,
3114 : sqlite3_column_blob(hStmt, 2), nBandBlockSize);
3115 : }
3116 : }
3117 : else
3118 : {
3119 422 : FillEmptyTileSingleBand(pabyTemp + (nBand - 1) * nBandBlockSize);
3120 : }
3121 4556 : sqlite3_finalize(hStmt);
3122 4556 : hStmt = nullptr;
3123 :
3124 : /* Copy the updated rectangle into the full tile */
3125 556108 : for (int iY = nDstYOffset; iY < nDstYOffset + nDstYSize; iY++)
3126 : {
3127 551552 : memcpy(pabyTemp +
3128 551552 : (static_cast<size_t>(nBand - 1) * nBlockXSize * nBlockYSize +
3129 551552 : static_cast<size_t>(iY) * nBlockXSize + nDstXOffset) *
3130 551552 : m_nDTSize,
3131 551552 : m_pabyCachedTiles +
3132 551552 : (static_cast<size_t>(nBand - 1) * nBlockXSize * nBlockYSize +
3133 551552 : static_cast<size_t>(iY) * nBlockXSize + nDstXOffset) *
3134 551552 : m_nDTSize,
3135 551552 : static_cast<size_t>(nDstXSize) * m_nDTSize);
3136 : }
3137 :
3138 : #ifdef notdef
3139 : static int nCounter = 1;
3140 : GDALDataset *poLogDS =
3141 : ((GDALDriver *)GDALGetDriverByName("GTiff"))
3142 : ->Create(CPLSPrintf("/tmp/partial_band_%d_%d.tif", 1, nCounter++),
3143 : nBlockXSize, nBlockYSize, nBands, m_eDT, NULL);
3144 : poLogDS->RasterIO(GF_Write, 0, 0, nBlockXSize, nBlockYSize,
3145 : pabyTemp + (nBand - 1) * nBandBlockSize, nBlockXSize,
3146 : nBlockYSize, m_eDT, 1, NULL, 0, 0, 0, NULL);
3147 : GDALClose(poLogDS);
3148 : #endif
3149 :
3150 4556 : if ((nOldFlags & l_nFlags) != 0)
3151 : {
3152 0 : CPLDebug("GPKG",
3153 : "Rewriting quadrant %d of band %d of tile (row=%d,col=%d)",
3154 : iQuadrantFlag, nBand, nRow, nCol);
3155 : }
3156 :
3157 4556 : l_nFlags |= nOldFlags;
3158 4556 : if (l_nFlags == nFullFlags)
3159 : {
3160 : #ifdef DEBUG_VERBOSE
3161 : CPLDebug("GPKG", "Got all quadrants for that tile");
3162 : #endif
3163 1192 : for (int iBand = 1; iBand <= nBands; iBand++)
3164 : {
3165 945 : if (iBand != nBand)
3166 : {
3167 698 : pszSQL = CPLSPrintf(
3168 : "SELECT tile_data_band_%d FROM partial_tiles WHERE "
3169 : "id = %d",
3170 : iBand, nExistingId);
3171 : #ifdef DEBUG_VERBOSE
3172 : CPLDebug("GPKG", "%s", pszSQL);
3173 : #endif
3174 698 : hStmt = nullptr;
3175 698 : rc = sqlite3_prepare_v2(m_hTempDB, pszSQL, -1, &hStmt, nullptr);
3176 698 : if (rc != SQLITE_OK)
3177 : {
3178 0 : CPLError(CE_Failure, CPLE_AppDefined,
3179 : "sqlite3_prepare_v2(%s) failed: %s", pszSQL,
3180 : sqlite3_errmsg(m_hTempDB));
3181 0 : return CE_Failure;
3182 : }
3183 :
3184 698 : rc = sqlite3_step(hStmt);
3185 698 : if (rc == SQLITE_ROW)
3186 : {
3187 698 : CPLAssert(sqlite3_column_bytes(hStmt, 0) ==
3188 : static_cast<int>(nBandBlockSize));
3189 698 : memcpy(m_pabyCachedTiles + (iBand - 1) * nBandBlockSize,
3190 : sqlite3_column_blob(hStmt, 0), nBandBlockSize);
3191 : }
3192 698 : sqlite3_finalize(hStmt);
3193 698 : hStmt = nullptr;
3194 : }
3195 : else
3196 : {
3197 247 : memcpy(m_pabyCachedTiles + (iBand - 1) * nBandBlockSize,
3198 247 : pabyTemp + (iBand - 1) * nBandBlockSize, nBandBlockSize);
3199 : }
3200 : }
3201 :
3202 247 : m_asCachedTilesDesc[0].nRow = nRow;
3203 247 : m_asCachedTilesDesc[0].nCol = nCol;
3204 247 : m_asCachedTilesDesc[0].nIdxWithinTileData = 0;
3205 247 : m_asCachedTilesDesc[0].abBandDirty[0] = true;
3206 247 : m_asCachedTilesDesc[0].abBandDirty[1] = true;
3207 247 : m_asCachedTilesDesc[0].abBandDirty[2] = true;
3208 247 : m_asCachedTilesDesc[0].abBandDirty[3] = true;
3209 :
3210 494 : pszSQL = CPLSPrintf("UPDATE partial_tiles SET zoom_level = %d, "
3211 : "partial_flag = 0, age = -1 WHERE id = %d",
3212 247 : -1 - m_nZoomLevel, nExistingId);
3213 247 : SQLCommand(m_hTempDB, pszSQL);
3214 : #ifdef DEBUG_VERBOSE
3215 : CPLDebug("GPKG", "%s", pszSQL);
3216 : #endif
3217 247 : CPLErr eErr = WriteTile();
3218 :
3219 : // Call DoPartialFlushOfPartialTilesIfNecessary() after using
3220 : // m_pabyCachedTiles as it is going to mess with it.
3221 247 : if (DoPartialFlushOfPartialTilesIfNecessary() != CE_None)
3222 0 : eErr = CE_None;
3223 :
3224 247 : return eErr;
3225 : }
3226 :
3227 4309 : if (nExistingId == 0)
3228 : {
3229 : OGRErr err;
3230 844 : pszSQL = CPLSPrintf("SELECT id FROM partial_tiles WHERE "
3231 : "partial_flag = 0 AND zoom_level = %d "
3232 : "AND tile_row = %d AND tile_column = %d",
3233 422 : -1 - m_nZoomLevel, nRow, nCol);
3234 : #ifdef DEBUG_VERBOSE
3235 : CPLDebug("GPKG", "%s", pszSQL);
3236 : #endif
3237 422 : nExistingId = SQLGetInteger(m_hTempDB, pszSQL, &err);
3238 422 : if (nExistingId == 0)
3239 : {
3240 418 : pszSQL =
3241 : "SELECT id FROM partial_tiles WHERE partial_flag = 0 LIMIT 1";
3242 : #ifdef DEBUG_VERBOSE
3243 : CPLDebug("GPKG", "%s", pszSQL);
3244 : #endif
3245 418 : nExistingId = SQLGetInteger(m_hTempDB, pszSQL, &err);
3246 : }
3247 : }
3248 :
3249 4309 : const GIntBig nAge = (m_poParentDS) ? m_poParentDS->m_nAge : m_nAge;
3250 4309 : if (nExistingId == 0)
3251 : {
3252 338 : pszSQL = CPLSPrintf(
3253 : "INSERT INTO partial_tiles "
3254 : "(zoom_level, tile_row, tile_column, tile_data_band_%d, "
3255 : "partial_flag, age) VALUES (%d, %d, %d, ?, %d, " CPL_FRMT_GIB ")",
3256 : nBand, m_nZoomLevel, nRow, nCol, l_nFlags, nAge);
3257 : }
3258 : else
3259 : {
3260 3971 : pszSQL = CPLSPrintf(
3261 : "UPDATE partial_tiles SET zoom_level = %d, "
3262 : "tile_row = %d, tile_column = %d, "
3263 : "tile_data_band_%d = ?, partial_flag = %d, age = " CPL_FRMT_GIB
3264 : " WHERE id = %d",
3265 : m_nZoomLevel, nRow, nCol, nBand, l_nFlags, nAge, nExistingId);
3266 : }
3267 4309 : if (m_poParentDS)
3268 3801 : m_poParentDS->m_nAge++;
3269 : else
3270 508 : m_nAge++;
3271 :
3272 : #ifdef DEBUG_VERBOSE
3273 : CPLDebug("GPKG", "%s", pszSQL);
3274 : #endif
3275 :
3276 4309 : hStmt = nullptr;
3277 4309 : rc = sqlite3_prepare_v2(m_hTempDB, pszSQL, -1, &hStmt, nullptr);
3278 4309 : if (rc != SQLITE_OK)
3279 : {
3280 0 : CPLError(CE_Failure, CPLE_AppDefined, "failed to prepare SQL %s: %s",
3281 : pszSQL, sqlite3_errmsg(m_hTempDB));
3282 0 : return CE_Failure;
3283 : }
3284 :
3285 4309 : sqlite3_bind_blob(hStmt, 1, pabyTemp + (nBand - 1) * nBandBlockSize,
3286 : static_cast<int>(nBandBlockSize), SQLITE_TRANSIENT);
3287 4309 : rc = sqlite3_step(hStmt);
3288 4309 : CPLErr eErr = CE_Failure;
3289 4309 : if (rc == SQLITE_DONE)
3290 4309 : eErr = CE_None;
3291 : else
3292 : {
3293 0 : CPLError(CE_Failure, CPLE_AppDefined,
3294 : "Failure when inserting partial tile (row=%d,col=%d) at "
3295 : "zoom_level=%d : %s",
3296 : nRow, nCol, m_nZoomLevel, sqlite3_errmsg(m_hTempDB));
3297 : }
3298 :
3299 4309 : sqlite3_finalize(hStmt);
3300 :
3301 : // Call DoPartialFlushOfPartialTilesIfNecessary() after using
3302 : // m_pabyCachedTiles as it is going to mess with it.
3303 4309 : if (DoPartialFlushOfPartialTilesIfNecessary() != CE_None)
3304 0 : eErr = CE_None;
3305 :
3306 4309 : return eErr;
3307 : }
3308 :
3309 : /************************************************************************/
3310 : /* IWriteBlock() */
3311 : /************************************************************************/
3312 :
3313 5839 : CPLErr GDALGPKGMBTilesLikeRasterBand::IWriteBlock(int nBlockXOff,
3314 : int nBlockYOff, void *pData)
3315 : {
3316 : #ifdef DEBUG_VERBOSE
3317 : CPLDebug(
3318 : "GPKG",
3319 : "IWriteBlock(nBand=%d,nBlockXOff=%d,nBlockYOff=%d,m_nZoomLevel=%d)",
3320 : nBand, nBlockXOff, nBlockYOff, m_poTPD->m_nZoomLevel);
3321 : #endif
3322 5839 : if (!m_poTPD->ICanIWriteBlock())
3323 : {
3324 0 : return CE_Failure;
3325 : }
3326 5839 : if (m_poTPD->m_poParentDS)
3327 1139 : m_poTPD->m_poParentDS->m_bHasModifiedTiles = true;
3328 : else
3329 4700 : m_poTPD->m_bHasModifiedTiles = true;
3330 :
3331 5839 : int nRow = nBlockYOff + m_poTPD->m_nShiftYTiles;
3332 5839 : int nCol = nBlockXOff + m_poTPD->m_nShiftXTiles;
3333 :
3334 5839 : int nRowMin = nRow;
3335 5839 : int nRowMax = nRowMin;
3336 5839 : if (m_poTPD->m_nShiftYPixelsMod)
3337 1360 : nRowMax++;
3338 :
3339 5839 : int nColMin = nCol;
3340 5839 : int nColMax = nColMin;
3341 5839 : if (m_poTPD->m_nShiftXPixelsMod)
3342 1322 : nColMax++;
3343 :
3344 5839 : CPLErr eErr = CE_None;
3345 :
3346 13038 : for (nRow = nRowMin; eErr == CE_None && nRow <= nRowMax; nRow++)
3347 : {
3348 17042 : for (nCol = nColMin; eErr == CE_None && nCol <= nColMax; nCol++)
3349 : {
3350 9843 : if (nRow < 0 || nCol < 0 || nRow >= m_poTPD->m_nTileMatrixHeight ||
3351 9703 : nCol >= m_poTPD->m_nTileMatrixWidth)
3352 : {
3353 167 : continue;
3354 : }
3355 :
3356 9676 : if (m_poTPD->m_nShiftXPixelsMod == 0 &&
3357 4505 : m_poTPD->m_nShiftYPixelsMod == 0)
3358 : {
3359 4447 : if (!(nRow == m_poTPD->m_asCachedTilesDesc[0].nRow &&
3360 250 : nCol == m_poTPD->m_asCachedTilesDesc[0].nCol &&
3361 5 : m_poTPD->m_asCachedTilesDesc[0].nIdxWithinTileData == 0))
3362 : {
3363 4442 : eErr = m_poTPD->WriteTile();
3364 :
3365 4442 : m_poTPD->m_asCachedTilesDesc[0].nRow = nRow;
3366 4442 : m_poTPD->m_asCachedTilesDesc[0].nCol = nCol;
3367 4442 : m_poTPD->m_asCachedTilesDesc[0].nIdxWithinTileData = 0;
3368 : }
3369 : }
3370 :
3371 : // Composite block data into tile, and check if all bands for this
3372 : // block are dirty, and if so write the tile
3373 9676 : bool bAllDirty = true;
3374 42207 : for (int iBand = 1; iBand <= poDS->GetRasterCount(); iBand++)
3375 : {
3376 32531 : GDALRasterBlock *poBlock = nullptr;
3377 32531 : GByte *pabySrc = nullptr;
3378 32531 : if (iBand == nBand)
3379 : {
3380 9676 : pabySrc = static_cast<GByte *>(pData);
3381 : }
3382 : else
3383 : {
3384 22855 : if (!(m_poTPD->m_nShiftXPixelsMod == 0 &&
3385 8353 : m_poTPD->m_nShiftYPixelsMod == 0))
3386 14646 : continue;
3387 :
3388 : // If the block for this band is not dirty, it might be
3389 : // dirty in cache
3390 8209 : if (m_poTPD->m_asCachedTilesDesc[0].abBandDirty[iBand - 1])
3391 505 : continue;
3392 : else
3393 : {
3394 : poBlock =
3395 7704 : cpl::down_cast<GDALGPKGMBTilesLikeRasterBand *>(
3396 7704 : poDS->GetRasterBand(iBand))
3397 7704 : ->TryGetLockedBlockRef(nBlockXOff, nBlockYOff);
3398 7704 : if (poBlock && poBlock->GetDirty())
3399 : {
3400 : pabySrc =
3401 7664 : static_cast<GByte *>(poBlock->GetDataRef());
3402 7664 : poBlock->MarkClean();
3403 : }
3404 : else
3405 : {
3406 40 : if (poBlock)
3407 2 : poBlock->DropLock();
3408 40 : bAllDirty = false;
3409 40 : continue;
3410 : }
3411 : }
3412 : }
3413 :
3414 17340 : if (m_poTPD->m_nShiftXPixelsMod == 0 &&
3415 12169 : m_poTPD->m_nShiftYPixelsMod == 0)
3416 12111 : m_poTPD->m_asCachedTilesDesc[0].abBandDirty[iBand - 1] =
3417 : true;
3418 :
3419 17340 : int nDstXOffset = 0;
3420 17340 : int nDstXSize = nBlockXSize;
3421 17340 : int nDstYOffset = 0;
3422 17340 : int nDstYSize = nBlockYSize;
3423 : // Composite block data into tile data
3424 17340 : if (m_poTPD->m_nShiftXPixelsMod == 0 &&
3425 12169 : m_poTPD->m_nShiftYPixelsMod == 0)
3426 : {
3427 :
3428 : #ifdef DEBUG_VERBOSE
3429 : if (eDataType == GDT_Byte &&
3430 : nBlockXOff * nBlockXSize <=
3431 : nRasterXSize - nBlockXSize &&
3432 : nBlockYOff * nBlockYSize > nRasterYSize - nBlockYSize)
3433 : {
3434 : bool bFoundNonZero = false;
3435 : for (int y = nRasterYSize - nBlockYOff * nBlockYSize;
3436 : y < nBlockYSize; y++)
3437 : {
3438 : for (int x = 0; x < nBlockXSize; x++)
3439 : {
3440 : if (pabySrc[static_cast<GPtrDiff_t>(y) *
3441 : nBlockXSize +
3442 : x] != 0 &&
3443 : !bFoundNonZero)
3444 : {
3445 : CPLDebug("GPKG",
3446 : "IWriteBlock(): Found non-zero "
3447 : "content in ghost part of "
3448 : "tile(nBand=%d,nBlockXOff=%d,"
3449 : "nBlockYOff=%d,m_nZoomLevel=%d)\n",
3450 : iBand, nBlockXOff, nBlockYOff,
3451 : m_poTPD->m_nZoomLevel);
3452 : bFoundNonZero = true;
3453 : }
3454 : }
3455 : }
3456 : }
3457 : #endif
3458 :
3459 12111 : const size_t nBandBlockSize =
3460 12111 : static_cast<size_t>(nBlockXSize) * nBlockYSize *
3461 12111 : m_nDTSize;
3462 12111 : memcpy(m_poTPD->m_pabyCachedTiles +
3463 12111 : (iBand - 1) * nBandBlockSize,
3464 : pabySrc, nBandBlockSize);
3465 :
3466 : // Make sure partial blocks are zero'ed outside of the
3467 : // validity area but do that only when know that JPEG will
3468 : // not be used so as to avoid edge effects (although we
3469 : // should probably repeat last pixels if we really want to
3470 : // do that, but that only makes sense if readers only clip
3471 : // to the gpkg_contents extent). Well, ere on the safe side
3472 : // for now
3473 12111 : if (m_poTPD->m_eTF != GPKG_TF_JPEG &&
3474 11863 : (nBlockXOff * nBlockXSize >=
3475 11863 : nRasterXSize - nBlockXSize ||
3476 11338 : nBlockYOff * nBlockYSize >=
3477 11338 : nRasterYSize - nBlockYSize))
3478 : {
3479 1021 : int nXEndValidity =
3480 1021 : nRasterXSize - nBlockXOff * nBlockXSize;
3481 1021 : if (nXEndValidity > nBlockXSize)
3482 496 : nXEndValidity = nBlockXSize;
3483 1021 : int nYEndValidity =
3484 1021 : nRasterYSize - nBlockYOff * nBlockYSize;
3485 1021 : if (nYEndValidity > nBlockYSize)
3486 294 : nYEndValidity = nBlockYSize;
3487 1021 : if (nXEndValidity < nBlockXSize)
3488 : {
3489 9366 : for (int iY = 0; iY < nYEndValidity; iY++)
3490 : {
3491 9187 : m_poTPD->FillBuffer(
3492 9187 : m_poTPD->m_pabyCachedTiles +
3493 9187 : ((static_cast<size_t>(iBand - 1) *
3494 9187 : nBlockYSize +
3495 9187 : iY) *
3496 9187 : nBlockXSize +
3497 9187 : nXEndValidity) *
3498 9187 : m_nDTSize,
3499 9187 : nBlockXSize - nXEndValidity);
3500 : }
3501 : }
3502 1021 : if (nYEndValidity < nBlockYSize)
3503 : {
3504 219 : m_poTPD->FillBuffer(
3505 219 : m_poTPD->m_pabyCachedTiles +
3506 219 : (static_cast<size_t>(iBand - 1) *
3507 219 : nBlockYSize +
3508 219 : nYEndValidity) *
3509 219 : nBlockXSize * m_nDTSize,
3510 219 : static_cast<size_t>(nBlockYSize -
3511 : nYEndValidity) *
3512 219 : nBlockXSize);
3513 : }
3514 12111 : }
3515 : }
3516 : else
3517 : {
3518 5229 : const int nXValid =
3519 5229 : (nBlockXOff * nBlockXSize > nRasterXSize - nBlockXSize)
3520 5229 : ? (nRasterXSize - nBlockXOff * nBlockXSize)
3521 : : nBlockXSize;
3522 5229 : const int nYValid =
3523 5229 : (nBlockYOff * nBlockYSize > nRasterYSize - nBlockYSize)
3524 5229 : ? (nRasterYSize - nBlockYOff * nBlockYSize)
3525 : : nBlockYSize;
3526 :
3527 5229 : int nSrcXOffset = 0;
3528 5229 : if (nCol == nColMin)
3529 : {
3530 2643 : nDstXOffset = m_poTPD->m_nShiftXPixelsMod;
3531 2643 : nDstXSize = std::min(
3532 2643 : nXValid, nBlockXSize - m_poTPD->m_nShiftXPixelsMod);
3533 : }
3534 : else
3535 : {
3536 2586 : nDstXOffset = 0;
3537 2586 : if (nXValid > nBlockXSize - m_poTPD->m_nShiftXPixelsMod)
3538 : {
3539 2212 : nDstXSize = nXValid - (nBlockXSize -
3540 2212 : m_poTPD->m_nShiftXPixelsMod);
3541 : }
3542 : else
3543 374 : nDstXSize = 0;
3544 2586 : nSrcXOffset = nBlockXSize - m_poTPD->m_nShiftXPixelsMod;
3545 : }
3546 :
3547 5229 : int nSrcYOffset = 0;
3548 5229 : if (nRow == nRowMin)
3549 : {
3550 2607 : nDstYOffset = m_poTPD->m_nShiftYPixelsMod;
3551 2607 : nDstYSize = std::min(
3552 2607 : nYValid, nBlockYSize - m_poTPD->m_nShiftYPixelsMod);
3553 : }
3554 : else
3555 : {
3556 2622 : nDstYOffset = 0;
3557 2622 : if (nYValid > nBlockYSize - m_poTPD->m_nShiftYPixelsMod)
3558 : {
3559 2232 : nDstYSize = nYValid - (nBlockYSize -
3560 2232 : m_poTPD->m_nShiftYPixelsMod);
3561 : }
3562 : else
3563 390 : nDstYSize = 0;
3564 2622 : nSrcYOffset = nBlockYSize - m_poTPD->m_nShiftYPixelsMod;
3565 : }
3566 :
3567 : #ifdef DEBUG_VERBOSE
3568 : CPLDebug("GPKG",
3569 : "Copy source tile x=%d,w=%d,y=%d,h=%d into "
3570 : "buffer at x=%d,y=%d",
3571 : nDstXOffset, nDstXSize, nDstYOffset, nDstYSize,
3572 : nSrcXOffset, nSrcYOffset);
3573 : #endif
3574 :
3575 5229 : if (nDstXSize > 0 && nDstYSize > 0)
3576 : {
3577 556108 : for (GPtrDiff_t y = 0; y < nDstYSize; y++)
3578 : {
3579 551552 : GByte *pDst = m_poTPD->m_pabyCachedTiles +
3580 551552 : (static_cast<size_t>(iBand - 1) *
3581 551552 : nBlockXSize * nBlockYSize +
3582 551552 : (y + nDstYOffset) * nBlockXSize +
3583 551552 : nDstXOffset) *
3584 551552 : m_nDTSize;
3585 551552 : GByte *pSrc =
3586 551552 : pabySrc + ((y + nSrcYOffset) * nBlockXSize +
3587 551552 : nSrcXOffset) *
3588 551552 : m_nDTSize;
3589 551552 : GDALCopyWords(pSrc, eDataType, m_nDTSize, pDst,
3590 : eDataType, m_nDTSize, nDstXSize);
3591 : }
3592 : }
3593 : }
3594 :
3595 17340 : if (poBlock)
3596 7664 : poBlock->DropLock();
3597 :
3598 17340 : if (!(m_poTPD->m_nShiftXPixelsMod == 0 &&
3599 12169 : m_poTPD->m_nShiftYPixelsMod == 0))
3600 : {
3601 5229 : m_poTPD->m_asCachedTilesDesc[0].nRow = -1;
3602 5229 : m_poTPD->m_asCachedTilesDesc[0].nCol = -1;
3603 5229 : m_poTPD->m_asCachedTilesDesc[0].nIdxWithinTileData = -1;
3604 5229 : if (nDstXSize > 0 && nDstYSize > 0)
3605 : {
3606 4556 : eErr = m_poTPD->WriteShiftedTile(
3607 : nRow, nCol, iBand, nDstXOffset, nDstYOffset,
3608 : nDstXSize, nDstYSize);
3609 : }
3610 : }
3611 : }
3612 :
3613 9676 : if (m_poTPD->m_nShiftXPixelsMod == 0 &&
3614 4505 : m_poTPD->m_nShiftYPixelsMod == 0)
3615 : {
3616 4447 : if (bAllDirty)
3617 : {
3618 4425 : eErr = m_poTPD->WriteTile();
3619 : }
3620 : }
3621 : }
3622 : }
3623 :
3624 5839 : return eErr;
3625 : }
3626 :
3627 : /************************************************************************/
3628 : /* GetNoDataValue() */
3629 : /************************************************************************/
3630 :
3631 19474 : double GDALGPKGMBTilesLikeRasterBand::GetNoDataValue(int *pbSuccess)
3632 : {
3633 19474 : if (m_bHasNoData)
3634 : {
3635 401 : if (pbSuccess)
3636 400 : *pbSuccess = TRUE;
3637 401 : return m_dfNoDataValue;
3638 : }
3639 19073 : return GDALPamRasterBand::GetNoDataValue(pbSuccess);
3640 : }
3641 :
3642 : /************************************************************************/
3643 : /* SetNoDataValueInternal() */
3644 : /************************************************************************/
3645 :
3646 75 : void GDALGPKGMBTilesLikeRasterBand::SetNoDataValueInternal(double dfNoDataValue)
3647 : {
3648 75 : m_bHasNoData = true;
3649 75 : m_dfNoDataValue = dfNoDataValue;
3650 75 : }
3651 :
3652 : /************************************************************************/
3653 : /* GDALGeoPackageRasterBand() */
3654 : /************************************************************************/
3655 :
3656 1792 : GDALGeoPackageRasterBand::GDALGeoPackageRasterBand(
3657 1792 : GDALGeoPackageDataset *poDSIn, int nTileWidth, int nTileHeight)
3658 1792 : : GDALGPKGMBTilesLikeRasterBand(poDSIn, nTileWidth, nTileHeight)
3659 : {
3660 1792 : poDS = poDSIn;
3661 1792 : }
3662 :
3663 : /************************************************************************/
3664 : /* GetOverviewCount() */
3665 : /************************************************************************/
3666 :
3667 8 : int GDALGeoPackageRasterBand::GetOverviewCount()
3668 : {
3669 : GDALGeoPackageDataset *poGDS =
3670 8 : cpl::down_cast<GDALGeoPackageDataset *>(poDS);
3671 8 : return poGDS->m_nOverviewCount;
3672 : }
3673 :
3674 : /************************************************************************/
3675 : /* GetOverviewCount() */
3676 : /************************************************************************/
3677 :
3678 41 : GDALRasterBand *GDALGeoPackageRasterBand::GetOverview(int nIdx)
3679 : {
3680 41 : GDALGeoPackageDataset *poGDS =
3681 : reinterpret_cast<GDALGeoPackageDataset *>(poDS);
3682 41 : if (nIdx < 0 || nIdx >= poGDS->m_nOverviewCount)
3683 1 : return nullptr;
3684 40 : return poGDS->m_papoOverviewDS[nIdx]->GetRasterBand(nBand);
3685 : }
3686 :
3687 : /************************************************************************/
3688 : /* SetNoDataValue() */
3689 : /************************************************************************/
3690 :
3691 23 : CPLErr GDALGeoPackageRasterBand::SetNoDataValue(double dfNoDataValue)
3692 : {
3693 : GDALGeoPackageDataset *poGDS =
3694 23 : cpl::down_cast<GDALGeoPackageDataset *>(poDS);
3695 :
3696 23 : if (eDataType == GDT_Byte)
3697 : {
3698 6 : if (!(dfNoDataValue >= 0 && dfNoDataValue <= 255 &&
3699 4 : static_cast<int>(dfNoDataValue) == dfNoDataValue))
3700 : {
3701 2 : CPLError(CE_Failure, CPLE_NotSupported,
3702 : "Invalid nodata value for a Byte band: %.17g",
3703 : dfNoDataValue);
3704 2 : return CE_Failure;
3705 : }
3706 :
3707 8 : for (int i = 1; i <= poGDS->nBands; ++i)
3708 : {
3709 5 : if (i != nBand)
3710 : {
3711 2 : int bHasNoData = FALSE;
3712 : double dfOtherNoDataValue =
3713 2 : poGDS->GetRasterBand(i)->GetNoDataValue(&bHasNoData);
3714 2 : if (bHasNoData && dfOtherNoDataValue != dfNoDataValue)
3715 : {
3716 1 : CPLError(
3717 : CE_Failure, CPLE_NotSupported,
3718 : "Only the same nodata value can be set on all bands");
3719 1 : return CE_Failure;
3720 : }
3721 : }
3722 : }
3723 :
3724 3 : SetNoDataValueInternal(dfNoDataValue);
3725 3 : poGDS->m_bMetadataDirty = true;
3726 3 : return CE_None;
3727 : }
3728 :
3729 17 : if (std::isnan(dfNoDataValue))
3730 : {
3731 0 : CPLError(CE_Warning, CPLE_NotSupported,
3732 : "A NaN nodata value cannot be recorded in "
3733 : "gpkg_2d_gridded_coverage_ancillary table");
3734 : }
3735 :
3736 17 : SetNoDataValueInternal(dfNoDataValue);
3737 :
3738 17 : char *pszSQL = sqlite3_mprintf(
3739 : "UPDATE gpkg_2d_gridded_coverage_ancillary SET data_null = ? "
3740 : "WHERE tile_matrix_set_name = '%q'",
3741 : poGDS->m_osRasterTable.c_str());
3742 17 : sqlite3_stmt *hStmt = nullptr;
3743 17 : int rc = sqlite3_prepare_v2(poGDS->IGetDB(), pszSQL, -1, &hStmt, nullptr);
3744 17 : if (rc == SQLITE_OK)
3745 : {
3746 17 : if (poGDS->m_eTF == GPKG_TF_PNG_16BIT)
3747 : {
3748 15 : if (eDataType == GDT_UInt16 && poGDS->m_dfOffset == 0.0 &&
3749 6 : poGDS->m_dfScale == 1.0 && dfNoDataValue >= 0 &&
3750 5 : dfNoDataValue <= 65535 &&
3751 5 : static_cast<GUInt16>(dfNoDataValue) == dfNoDataValue)
3752 : {
3753 5 : poGDS->m_usGPKGNull = static_cast<GUInt16>(dfNoDataValue);
3754 : }
3755 : else
3756 : {
3757 10 : poGDS->m_usGPKGNull = 65535;
3758 : }
3759 15 : sqlite3_bind_double(hStmt, 1, poGDS->m_usGPKGNull);
3760 : }
3761 : else
3762 : {
3763 2 : sqlite3_bind_double(hStmt, 1, static_cast<float>(dfNoDataValue));
3764 : }
3765 17 : rc = sqlite3_step(hStmt);
3766 17 : sqlite3_finalize(hStmt);
3767 : }
3768 17 : sqlite3_free(pszSQL);
3769 :
3770 17 : return (rc == SQLITE_OK) ? CE_None : CE_Failure;
3771 : }
3772 :
3773 : /************************************************************************/
3774 : /* LoadBandMetadata() */
3775 : /************************************************************************/
3776 :
3777 1206 : void GDALGeoPackageRasterBand::LoadBandMetadata()
3778 : {
3779 : GDALGeoPackageDataset *poGDS =
3780 1206 : cpl::down_cast<GDALGeoPackageDataset *>(poDS);
3781 :
3782 1206 : if (m_bHasReadMetadataFromStorage)
3783 1105 : return;
3784 :
3785 464 : m_bHasReadMetadataFromStorage = true;
3786 :
3787 464 : poGDS->TryLoadXML();
3788 :
3789 464 : if (!poGDS->HasMetadataTables())
3790 363 : return;
3791 :
3792 101 : char *pszSQL = sqlite3_mprintf(
3793 : "SELECT md.metadata, md.md_standard_uri, md.mime_type "
3794 : "FROM gpkg_metadata md "
3795 : "JOIN gpkg_metadata_reference mdr ON (md.id = mdr.md_file_id ) "
3796 : "WHERE "
3797 : "mdr.reference_scope = 'table' AND lower(mdr.table_name) = "
3798 : "lower('%q') ORDER BY md.id "
3799 : "LIMIT 1000", // to avoid denial of service
3800 : poGDS->m_osRasterTable.c_str());
3801 :
3802 101 : auto oResult = SQLQuery(poGDS->hDB, pszSQL);
3803 101 : sqlite3_free(pszSQL);
3804 101 : if (!oResult)
3805 : {
3806 0 : return;
3807 : }
3808 :
3809 : /* GDAL metadata */
3810 198 : for (int i = 0; i < oResult->RowCount(); i++)
3811 : {
3812 97 : const char *pszMetadata = oResult->GetValue(0, i);
3813 97 : const char *pszMDStandardURI = oResult->GetValue(1, i);
3814 97 : const char *pszMimeType = oResult->GetValue(2, i);
3815 97 : if (pszMetadata && pszMDStandardURI && pszMimeType &&
3816 97 : EQUAL(pszMDStandardURI, "http://gdal.org") &&
3817 73 : EQUAL(pszMimeType, "text/xml"))
3818 : {
3819 73 : CPLXMLNode *psXMLNode = CPLParseXMLString(pszMetadata);
3820 73 : if (psXMLNode)
3821 : {
3822 146 : GDALMultiDomainMetadata oLocalMDMD;
3823 73 : oLocalMDMD.XMLInit(psXMLNode, FALSE);
3824 :
3825 73 : CSLConstList papszDomainList = oLocalMDMD.GetDomainList();
3826 206 : for (CSLConstList papszIter = papszDomainList;
3827 206 : papszIter && *papszIter; ++papszIter)
3828 : {
3829 133 : if (STARTS_WITH(*papszIter, "BAND_"))
3830 : {
3831 5 : int l_nBand = atoi(*papszIter + strlen("BAND_"));
3832 5 : if (l_nBand >= 1 && l_nBand <= poGDS->GetRasterCount())
3833 : {
3834 : auto l_poBand =
3835 5 : cpl::down_cast<GDALGeoPackageRasterBand *>(
3836 : poGDS->GetRasterBand(l_nBand));
3837 5 : l_poBand->m_bHasReadMetadataFromStorage = true;
3838 :
3839 10 : char **papszMD = CSLDuplicate(
3840 5 : oLocalMDMD.GetMetadata(*papszIter));
3841 10 : papszMD = CSLMerge(
3842 : papszMD,
3843 5 : GDALGPKGMBTilesLikeRasterBand::GetMetadata(""));
3844 5 : l_poBand->GDALPamRasterBand::SetMetadata(papszMD);
3845 5 : CSLDestroy(papszMD);
3846 : }
3847 : }
3848 : }
3849 :
3850 73 : CPLDestroyXMLNode(psXMLNode);
3851 : }
3852 : }
3853 : }
3854 : }
3855 :
3856 : /************************************************************************/
3857 : /* GetMetadata() */
3858 : /************************************************************************/
3859 :
3860 921 : char **GDALGeoPackageRasterBand::GetMetadata(const char *pszDomain)
3861 : {
3862 921 : GDALGeoPackageDataset *poGDS =
3863 : reinterpret_cast<GDALGeoPackageDataset *>(poDS);
3864 921 : LoadBandMetadata(); /* force loading from storage if needed */
3865 :
3866 15 : if (poGDS->eAccess == GA_ReadOnly && eDataType != GDT_Byte &&
3867 11 : (pszDomain == nullptr || EQUAL(pszDomain, "")) &&
3868 11 : !m_bMinMaxComputedFromTileAncillary &&
3869 944 : !GDALGPKGMBTilesLikeRasterBand::GetMetadataItem("STATISTICS_MINIMUM") &&
3870 8 : !GDALGPKGMBTilesLikeRasterBand::GetMetadataItem("STATISTICS_MAXIMUM"))
3871 : {
3872 8 : m_bMinMaxComputedFromTileAncillary = true;
3873 :
3874 8 : const int nColMin = poGDS->m_nShiftXTiles;
3875 8 : const int nColMax =
3876 8 : (nRasterXSize - 1 + poGDS->m_nShiftXPixelsMod) / nBlockXSize +
3877 8 : poGDS->m_nShiftXTiles;
3878 8 : const int nRowMin = poGDS->m_nShiftYTiles;
3879 8 : const int nRowMax =
3880 8 : (nRasterYSize - 1 + poGDS->m_nShiftYPixelsMod) / nBlockYSize +
3881 8 : poGDS->m_nShiftYTiles;
3882 :
3883 8 : bool bOK = false;
3884 8 : if (poGDS->m_nShiftXPixelsMod == 0 && poGDS->m_nShiftYPixelsMod == 0 &&
3885 8 : (nRasterXSize % nBlockXSize) == 0 &&
3886 2 : (nRasterYSize % nBlockYSize) == 0)
3887 : {
3888 : // If the area of interest matches entire tiles, then we can
3889 : // use tile statistics
3890 2 : bOK = true;
3891 : }
3892 6 : else if (m_bHasNoData)
3893 : {
3894 : // Otherwise, in the case where we have nodata, we assume that
3895 : // if the area of interest is at least larger than the existing
3896 : // tiles, the tile statistics will be reliable.
3897 2 : char *pszSQL = sqlite3_mprintf(
3898 : "SELECT MIN(tile_column), MAX(tile_column), "
3899 : "MIN(tile_row), MAX(tile_row) FROM \"%w\" "
3900 : "WHERE zoom_level = %d",
3901 : poGDS->m_osRasterTable.c_str(), poGDS->m_nZoomLevel);
3902 4 : auto sResult = SQLQuery(poGDS->IGetDB(), pszSQL);
3903 2 : if (sResult && sResult->RowCount() == 1)
3904 : {
3905 2 : const char *pszMinX = sResult->GetValue(0, 0);
3906 2 : const char *pszMaxX = sResult->GetValue(1, 0);
3907 2 : const char *pszMinY = sResult->GetValue(2, 0);
3908 2 : const char *pszMaxY = sResult->GetValue(3, 0);
3909 2 : if (pszMinX && pszMaxX && pszMinY && pszMaxY)
3910 : {
3911 2 : bOK = atoi(pszMinX) >= nColMin &&
3912 2 : atoi(pszMaxX) <= nColMax &&
3913 4 : atoi(pszMinY) >= nRowMin && atoi(pszMaxY) <= nRowMax;
3914 : }
3915 : }
3916 2 : sqlite3_free(pszSQL);
3917 : }
3918 :
3919 8 : if (bOK)
3920 : {
3921 4 : char *pszSQL = sqlite3_mprintf(
3922 : "SELECT MIN(min), MAX(max) FROM "
3923 : "gpkg_2d_gridded_tile_ancillary WHERE tpudt_id "
3924 : "IN (SELECT id FROM \"%w\" WHERE "
3925 : "zoom_level = %d AND "
3926 : "tile_column >= %d AND tile_column <= %d AND "
3927 : "tile_row >= %d AND tile_row <= %d)",
3928 : poGDS->m_osRasterTable.c_str(), poGDS->m_nZoomLevel, nColMin,
3929 : nColMax, nRowMin, nRowMax);
3930 8 : auto sResult = SQLQuery(poGDS->IGetDB(), pszSQL);
3931 4 : CPLDebug("GPKG", "%s", pszSQL);
3932 4 : if (sResult && sResult->RowCount() == 1)
3933 : {
3934 4 : const char *pszMin = sResult->GetValue(0, 0);
3935 4 : const char *pszMax = sResult->GetValue(1, 0);
3936 4 : if (pszMin)
3937 : {
3938 4 : m_dfStatsMinFromTileAncillary = CPLAtof(pszMin);
3939 : }
3940 4 : if (pszMax)
3941 : {
3942 4 : m_dfStatsMaxFromTileAncillary = CPLAtof(pszMax);
3943 : }
3944 : }
3945 4 : sqlite3_free(pszSQL);
3946 : }
3947 : }
3948 :
3949 685 : if (m_bAddImplicitStatistics && m_bMinMaxComputedFromTileAncillary &&
3950 8 : (pszDomain == nullptr || EQUAL(pszDomain, "")) &&
3951 1614 : !GDALGPKGMBTilesLikeRasterBand::GetMetadataItem("STATISTICS_MINIMUM") &&
3952 8 : !GDALGPKGMBTilesLikeRasterBand::GetMetadataItem("STATISTICS_MAXIMUM"))
3953 : {
3954 : m_aosMD.Assign(CSLDuplicate(
3955 8 : GDALGPKGMBTilesLikeRasterBand::GetMetadata(pszDomain)));
3956 8 : if (!std::isnan(m_dfStatsMinFromTileAncillary))
3957 : {
3958 : m_aosMD.SetNameValue(
3959 : "STATISTICS_MINIMUM",
3960 4 : CPLSPrintf("%.14g", m_dfStatsMinFromTileAncillary));
3961 : }
3962 8 : if (!std::isnan(m_dfStatsMaxFromTileAncillary))
3963 : {
3964 : m_aosMD.SetNameValue(
3965 : "STATISTICS_MAXIMUM",
3966 4 : CPLSPrintf("%.14g", m_dfStatsMaxFromTileAncillary));
3967 : }
3968 8 : return m_aosMD.List();
3969 : }
3970 :
3971 913 : return GDALGPKGMBTilesLikeRasterBand::GetMetadata(pszDomain);
3972 : }
3973 :
3974 : /************************************************************************/
3975 : /* GetMetadataItem() */
3976 : /************************************************************************/
3977 :
3978 251 : const char *GDALGeoPackageRasterBand::GetMetadataItem(const char *pszName,
3979 : const char *pszDomain)
3980 : {
3981 251 : LoadBandMetadata(); /* force loading from storage if needed */
3982 :
3983 251 : if (m_bAddImplicitStatistics && eDataType != GDT_Byte &&
3984 7 : (pszDomain == nullptr || EQUAL(pszDomain, "")) &&
3985 4 : (EQUAL(pszName, "STATISTICS_MINIMUM") ||
3986 2 : EQUAL(pszName, "STATISTICS_MAXIMUM")))
3987 : {
3988 2 : return CSLFetchNameValue(GetMetadata(), pszName);
3989 : }
3990 :
3991 249 : return GDALGPKGMBTilesLikeRasterBand::GetMetadataItem(pszName, pszDomain);
3992 : }
3993 :
3994 : /************************************************************************/
3995 : /* SetMetadata() */
3996 : /************************************************************************/
3997 :
3998 14 : CPLErr GDALGeoPackageRasterBand::SetMetadata(char **papszMetadata,
3999 : const char *pszDomain)
4000 : {
4001 : GDALGeoPackageDataset *poGDS =
4002 14 : cpl::down_cast<GDALGeoPackageDataset *>(poDS);
4003 14 : LoadBandMetadata(); /* force loading from storage if needed */
4004 14 : poGDS->m_bMetadataDirty = true;
4005 64 : for (CSLConstList papszIter = papszMetadata; papszIter && *papszIter;
4006 : ++papszIter)
4007 : {
4008 50 : if (STARTS_WITH(*papszIter, "STATISTICS_"))
4009 46 : m_bStatsMetadataSetInThisSession = true;
4010 : }
4011 14 : return GDALPamRasterBand::SetMetadata(papszMetadata, pszDomain);
4012 : }
4013 :
4014 : /************************************************************************/
4015 : /* SetMetadataItem() */
4016 : /************************************************************************/
4017 :
4018 20 : CPLErr GDALGeoPackageRasterBand::SetMetadataItem(const char *pszName,
4019 : const char *pszValue,
4020 : const char *pszDomain)
4021 : {
4022 : GDALGeoPackageDataset *poGDS =
4023 20 : cpl::down_cast<GDALGeoPackageDataset *>(poDS);
4024 20 : LoadBandMetadata(); /* force loading from storage if needed */
4025 20 : poGDS->m_bMetadataDirty = true;
4026 20 : if (STARTS_WITH(pszName, "STATISTICS_"))
4027 20 : m_bStatsMetadataSetInThisSession = true;
4028 20 : return GDALPamRasterBand::SetMetadataItem(pszName, pszValue, pszDomain);
4029 : }
4030 :
4031 : /************************************************************************/
4032 : /* InvalidateStatistics() */
4033 : /************************************************************************/
4034 :
4035 334 : void GDALGeoPackageRasterBand::InvalidateStatistics()
4036 : {
4037 334 : bool bModified = false;
4038 668 : CPLStringList aosMD(CSLDuplicate(GetMetadata()));
4039 354 : for (CSLConstList papszIter = GetMetadata(); papszIter && *papszIter;
4040 : ++papszIter)
4041 : {
4042 20 : if (STARTS_WITH(*papszIter, "STATISTICS_"))
4043 : {
4044 20 : char *pszKey = nullptr;
4045 20 : CPLParseNameValue(*papszIter, &pszKey);
4046 20 : CPLAssert(pszKey);
4047 20 : aosMD.SetNameValue(pszKey, nullptr);
4048 20 : CPLFree(pszKey);
4049 20 : bModified = true;
4050 : }
4051 : }
4052 334 : if (bModified)
4053 4 : SetMetadata(aosMD.List());
4054 334 : }
|
