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