Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: COG Driver
4 : * Purpose: Cloud optimized GeoTIFF write support.
5 : * Author: Even Rouault <even dot rouault at spatialys dot com>
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2019, Even Rouault <even dot rouault at spatialys dot com>
9 : *
10 : * SPDX-License-Identifier: MIT
11 : ****************************************************************************/
12 :
13 : #include "cpl_port.h"
14 :
15 : #include "gdal_priv.h"
16 : #include "gtiff.h"
17 : #include "gt_overview.h"
18 : #include "gdal_utils.h"
19 : #include "gdalwarper.h"
20 : #include "cogdriver.h"
21 : #include "geotiff.h"
22 :
23 : #include "tilematrixset.hpp"
24 :
25 : #include <algorithm>
26 : #include <memory>
27 : #include <vector>
28 :
29 : static bool gbHasLZW = false;
30 :
31 : extern "C" CPL_DLL void GDALRegister_COG();
32 :
33 : /************************************************************************/
34 : /* HasZSTDCompression() */
35 : /************************************************************************/
36 :
37 160 : static bool HasZSTDCompression()
38 : {
39 160 : TIFFCodec *codecs = TIFFGetConfiguredCODECs();
40 160 : bool bHasZSTD = false;
41 3198 : for (TIFFCodec *c = codecs; c->name; ++c)
42 : {
43 3198 : if (c->scheme == COMPRESSION_ZSTD)
44 : {
45 160 : bHasZSTD = true;
46 160 : break;
47 : }
48 : }
49 160 : _TIFFfree(codecs);
50 160 : return bHasZSTD;
51 : }
52 :
53 : /************************************************************************/
54 : /* GetTmpFilename() */
55 : /************************************************************************/
56 :
57 66 : static CPLString GetTmpFilename(const char *pszFilename, const char *pszExt)
58 : {
59 : const bool bSupportsRandomWrite =
60 66 : VSISupportsRandomWrite(pszFilename, false);
61 66 : CPLString osTmpFilename;
62 130 : if (!bSupportsRandomWrite ||
63 64 : CPLGetConfigOption("CPL_TMPDIR", nullptr) != nullptr)
64 : {
65 2 : osTmpFilename = CPLGenerateTempFilenameSafe(
66 4 : CPLGetBasenameSafe(pszFilename).c_str());
67 : }
68 : else
69 64 : osTmpFilename = pszFilename;
70 66 : osTmpFilename += '.';
71 66 : osTmpFilename += pszExt;
72 66 : VSIUnlink(osTmpFilename);
73 66 : return osTmpFilename;
74 : }
75 :
76 : /************************************************************************/
77 : /* GetResampling() */
78 : /************************************************************************/
79 :
80 75 : static const char *GetResampling(GDALDataset *poSrcDS)
81 : {
82 75 : return poSrcDS->GetRasterBand(1)->GetColorTable() ? "NEAREST" : "CUBIC";
83 : }
84 :
85 : /************************************************************************/
86 : /* GetPredictor() */
87 : /************************************************************************/
88 266 : static const char *GetPredictor(GDALDataset *poSrcDS, const char *pszPredictor)
89 : {
90 266 : if (pszPredictor == nullptr)
91 0 : return nullptr;
92 :
93 266 : if (EQUAL(pszPredictor, "YES") || EQUAL(pszPredictor, "ON") ||
94 265 : EQUAL(pszPredictor, "TRUE"))
95 : {
96 1 : if (GDALDataTypeIsFloating(
97 1 : poSrcDS->GetRasterBand(1)->GetRasterDataType()))
98 0 : return "3";
99 : else
100 1 : return "2";
101 : }
102 265 : else if (EQUAL(pszPredictor, "STANDARD") || EQUAL(pszPredictor, "2"))
103 : {
104 1 : return "2";
105 : }
106 264 : else if (EQUAL(pszPredictor, "FLOATING_POINT") || EQUAL(pszPredictor, "3"))
107 : {
108 0 : return "3";
109 : }
110 264 : return nullptr;
111 : }
112 :
113 : /************************************************************************/
114 : /* COGGetTargetSRS() */
115 : /************************************************************************/
116 :
117 39 : static bool COGGetTargetSRS(const char *const *papszOptions,
118 : CPLString &osTargetSRS,
119 : std::unique_ptr<gdal::TileMatrixSet> &poTM)
120 : {
121 39 : osTargetSRS = CSLFetchNameValueDef(papszOptions, "TARGET_SRS", "");
122 : CPLString osTilingScheme(
123 78 : CSLFetchNameValueDef(papszOptions, "TILING_SCHEME", "CUSTOM"));
124 39 : if (EQUAL(osTargetSRS, "") && EQUAL(osTilingScheme, "CUSTOM"))
125 2 : return false;
126 :
127 37 : if (!EQUAL(osTilingScheme, "CUSTOM"))
128 : {
129 30 : poTM = gdal::TileMatrixSet::parse(osTilingScheme);
130 30 : if (poTM == nullptr)
131 0 : return false;
132 30 : if (!poTM->haveAllLevelsSameTopLeft())
133 : {
134 0 : CPLError(CE_Failure, CPLE_NotSupported,
135 : "Unsupported tiling scheme: not all zoom levels have same "
136 : "top left corner");
137 0 : return false;
138 : }
139 30 : if (!poTM->haveAllLevelsSameTileSize())
140 : {
141 0 : CPLError(CE_Failure, CPLE_NotSupported,
142 : "Unsupported tiling scheme: not all zoom levels have same "
143 : "tile size");
144 0 : return false;
145 : }
146 30 : if (poTM->hasVariableMatrixWidth())
147 : {
148 0 : CPLError(CE_Failure, CPLE_NotSupported,
149 : "Unsupported tiling scheme: some levels have variable "
150 : "matrix width");
151 0 : return false;
152 : }
153 30 : if (!osTargetSRS.empty())
154 : {
155 0 : CPLError(CE_Warning, CPLE_AppDefined, "Ignoring TARGET_SRS option");
156 : }
157 30 : osTargetSRS = poTM->crs();
158 :
159 : // "Normalize" SRS as AUTH:CODE
160 60 : OGRSpatialReference oTargetSRS;
161 30 : oTargetSRS.SetFromUserInput(
162 : osTargetSRS,
163 : OGRSpatialReference::SET_FROM_USER_INPUT_LIMITATIONS_get());
164 30 : const char *pszAuthCode = oTargetSRS.GetAuthorityCode(nullptr);
165 30 : const char *pszAuthName = oTargetSRS.GetAuthorityName(nullptr);
166 30 : if (pszAuthName && pszAuthCode)
167 : {
168 30 : osTargetSRS = pszAuthName;
169 30 : osTargetSRS += ':';
170 30 : osTargetSRS += pszAuthCode;
171 : }
172 : }
173 :
174 37 : return true;
175 : }
176 :
177 : // Used by gdalwarp
178 3 : bool COGGetTargetSRS(const char *const *papszOptions, CPLString &osTargetSRS)
179 : {
180 3 : std::unique_ptr<gdal::TileMatrixSet> poTM;
181 6 : return COGGetTargetSRS(papszOptions, osTargetSRS, poTM);
182 : }
183 :
184 : // Used by gdalwarp
185 34 : std::string COGGetResampling(GDALDataset *poSrcDS,
186 : const char *const *papszOptions)
187 : {
188 : return CSLFetchNameValueDef(papszOptions, "WARP_RESAMPLING",
189 : CSLFetchNameValueDef(papszOptions, "RESAMPLING",
190 34 : GetResampling(poSrcDS)));
191 : }
192 :
193 : /************************************************************************/
194 : /* COGGetWarpingCharacteristics() */
195 : /************************************************************************/
196 :
197 36 : static bool COGGetWarpingCharacteristics(
198 : GDALDataset *poSrcDS, const char *const *papszOptions,
199 : CPLString &osResampling, CPLString &osTargetSRS, int &nXSize, int &nYSize,
200 : double &dfMinX, double &dfMinY, double &dfMaxX, double &dfMaxY,
201 : double &dfRes, std::unique_ptr<gdal::TileMatrixSet> &poTM, int &nZoomLevel,
202 : int &nAlignedLevels)
203 : {
204 36 : if (!COGGetTargetSRS(papszOptions, osTargetSRS, poTM))
205 1 : return false;
206 :
207 70 : CPLStringList aosTO;
208 35 : aosTO.SetNameValue("DST_SRS", osTargetSRS);
209 35 : void *hTransformArg = nullptr;
210 :
211 70 : OGRSpatialReference oTargetSRS;
212 35 : oTargetSRS.SetFromUserInput(
213 : osTargetSRS,
214 : OGRSpatialReference::SET_FROM_USER_INPUT_LIMITATIONS_get());
215 35 : const char *pszAuthCode = oTargetSRS.GetAuthorityCode(nullptr);
216 35 : const int nEPSGCode = pszAuthCode ? atoi(pszAuthCode) : 0;
217 :
218 : // Hack to compensate for GDALSuggestedWarpOutput2() failure (or not
219 : // ideal suggestion with PROJ 8) when reprojecting latitude = +/- 90 to
220 : // EPSG:3857.
221 : double adfSrcGeoTransform[6];
222 35 : std::unique_ptr<GDALDataset> poTmpDS;
223 62 : if (nEPSGCode == 3857 &&
224 27 : poSrcDS->GetGeoTransform(adfSrcGeoTransform) == CE_None &&
225 89 : adfSrcGeoTransform[2] == 0 && adfSrcGeoTransform[4] == 0 &&
226 27 : adfSrcGeoTransform[5] < 0)
227 : {
228 27 : const auto poSrcSRS = poSrcDS->GetSpatialRef();
229 32 : if (poSrcSRS && poSrcSRS->IsGeographic() &&
230 5 : !poSrcSRS->IsDerivedGeographic())
231 : {
232 5 : double maxLat = adfSrcGeoTransform[3];
233 5 : double minLat = adfSrcGeoTransform[3] +
234 5 : poSrcDS->GetRasterYSize() * adfSrcGeoTransform[5];
235 : // Corresponds to the latitude of below MAX_GM
236 5 : constexpr double MAX_LAT = 85.0511287798066;
237 5 : bool bModified = false;
238 5 : if (maxLat > MAX_LAT)
239 : {
240 4 : maxLat = MAX_LAT;
241 4 : bModified = true;
242 : }
243 5 : if (minLat < -MAX_LAT)
244 : {
245 4 : minLat = -MAX_LAT;
246 4 : bModified = true;
247 : }
248 5 : if (bModified)
249 : {
250 4 : CPLStringList aosOptions;
251 4 : aosOptions.AddString("-of");
252 4 : aosOptions.AddString("VRT");
253 4 : aosOptions.AddString("-projwin");
254 : aosOptions.AddString(
255 4 : CPLSPrintf("%.17g", adfSrcGeoTransform[0]));
256 4 : aosOptions.AddString(CPLSPrintf("%.17g", maxLat));
257 : aosOptions.AddString(
258 4 : CPLSPrintf("%.17g", adfSrcGeoTransform[0] +
259 4 : poSrcDS->GetRasterXSize() *
260 4 : adfSrcGeoTransform[1]));
261 4 : aosOptions.AddString(CPLSPrintf("%.17g", minLat));
262 : auto psOptions =
263 4 : GDALTranslateOptionsNew(aosOptions.List(), nullptr);
264 4 : poTmpDS.reset(GDALDataset::FromHandle(GDALTranslate(
265 : "", GDALDataset::ToHandle(poSrcDS), psOptions, nullptr)));
266 4 : GDALTranslateOptionsFree(psOptions);
267 4 : if (poTmpDS)
268 : {
269 8 : hTransformArg = GDALCreateGenImgProjTransformer2(
270 4 : GDALDataset::FromHandle(poTmpDS.get()), nullptr,
271 : aosTO.List());
272 4 : if (hTransformArg == nullptr)
273 : {
274 0 : return false;
275 : }
276 : }
277 : }
278 : }
279 : }
280 35 : if (hTransformArg == nullptr)
281 : {
282 : hTransformArg =
283 31 : GDALCreateGenImgProjTransformer2(poSrcDS, nullptr, aosTO.List());
284 31 : if (hTransformArg == nullptr)
285 : {
286 0 : return false;
287 : }
288 : }
289 :
290 35 : GDALTransformerInfo *psInfo =
291 : static_cast<GDALTransformerInfo *>(hTransformArg);
292 : double adfGeoTransform[6];
293 : double adfExtent[4];
294 :
295 35 : if (GDALSuggestedWarpOutput2(poTmpDS ? poTmpDS.get() : poSrcDS,
296 : psInfo->pfnTransform, hTransformArg,
297 : adfGeoTransform, &nXSize, &nYSize, adfExtent,
298 35 : 0) != CE_None)
299 : {
300 0 : GDALDestroyGenImgProjTransformer(hTransformArg);
301 0 : return false;
302 : }
303 :
304 35 : GDALDestroyGenImgProjTransformer(hTransformArg);
305 35 : hTransformArg = nullptr;
306 35 : poTmpDS.reset();
307 :
308 35 : dfMinX = adfExtent[0];
309 35 : dfMinY = adfExtent[1];
310 35 : dfMaxX = adfExtent[2];
311 35 : dfMaxY = adfExtent[3];
312 35 : dfRes = adfGeoTransform[1];
313 :
314 70 : const CPLString osExtent(CSLFetchNameValueDef(papszOptions, "EXTENT", ""));
315 70 : const CPLString osRes(CSLFetchNameValueDef(papszOptions, "RES", ""));
316 35 : if (poTM)
317 : {
318 28 : if (!osExtent.empty())
319 : {
320 0 : CPLError(CE_Warning, CPLE_AppDefined, "Ignoring EXTENT option");
321 : }
322 28 : if (!osRes.empty())
323 : {
324 0 : CPLError(CE_Warning, CPLE_AppDefined, "Ignoring RES option");
325 : }
326 : const bool bInvertAxis =
327 56 : oTargetSRS.EPSGTreatsAsLatLong() != FALSE ||
328 28 : oTargetSRS.EPSGTreatsAsNorthingEasting() != FALSE;
329 :
330 28 : const auto &bbox = poTM->bbox();
331 28 : if (bbox.mCrs == poTM->crs())
332 : {
333 28 : if (dfMaxX <
334 56 : (bInvertAxis ? bbox.mLowerCornerY : bbox.mLowerCornerX) ||
335 28 : dfMinX >
336 56 : (bInvertAxis ? bbox.mUpperCornerY : bbox.mUpperCornerX) ||
337 28 : dfMaxY <
338 56 : (bInvertAxis ? bbox.mLowerCornerX : bbox.mLowerCornerY) ||
339 28 : dfMinY >
340 28 : (bInvertAxis ? bbox.mUpperCornerX : bbox.mUpperCornerY))
341 : {
342 0 : CPLError(CE_Failure, CPLE_AppDefined,
343 : "Raster extent completely outside of tile matrix set "
344 : "bounding box");
345 2 : return false;
346 : }
347 : }
348 :
349 28 : const auto &tmList = poTM->tileMatrixList();
350 28 : const int nBlockSize = atoi(CSLFetchNameValueDef(
351 28 : papszOptions, "BLOCKSIZE", CPLSPrintf("%d", tmList[0].mTileWidth)));
352 28 : dfRes = 0.0;
353 :
354 : const char *pszZoomLevel =
355 28 : CSLFetchNameValue(papszOptions, "ZOOM_LEVEL");
356 28 : if (pszZoomLevel)
357 : {
358 3 : nZoomLevel = atoi(pszZoomLevel);
359 3 : if (nZoomLevel < 0 || nZoomLevel >= static_cast<int>(tmList.size()))
360 : {
361 2 : CPLError(CE_Failure, CPLE_AppDefined,
362 : "Invalid zoom level: should be in [0,%d]",
363 2 : static_cast<int>(tmList.size()) - 1);
364 2 : return false;
365 : }
366 : }
367 : else
368 : {
369 25 : double dfComputedRes = adfGeoTransform[1];
370 25 : double dfPrevRes = 0.0;
371 280 : for (; nZoomLevel < static_cast<int>(tmList.size()); nZoomLevel++)
372 : {
373 280 : dfRes = tmList[nZoomLevel].mResX * tmList[0].mTileWidth /
374 : nBlockSize;
375 280 : if (dfComputedRes > dfRes ||
376 260 : fabs(dfComputedRes - dfRes) / dfRes <= 1e-8)
377 : break;
378 255 : dfPrevRes = dfRes;
379 : }
380 25 : if (nZoomLevel == static_cast<int>(tmList.size()))
381 : {
382 0 : CPLError(CE_Failure, CPLE_AppDefined,
383 : "Could not find an appropriate zoom level");
384 0 : return false;
385 : }
386 :
387 25 : if (nZoomLevel > 0 && fabs(dfComputedRes - dfRes) / dfRes > 1e-8)
388 : {
389 19 : const char *pszZoomLevelStrategy = CSLFetchNameValueDef(
390 : papszOptions, "ZOOM_LEVEL_STRATEGY", "AUTO");
391 19 : if (EQUAL(pszZoomLevelStrategy, "LOWER"))
392 : {
393 1 : nZoomLevel--;
394 : }
395 18 : else if (EQUAL(pszZoomLevelStrategy, "UPPER"))
396 : {
397 : /* do nothing */
398 : }
399 : else
400 : {
401 17 : if (dfPrevRes / dfComputedRes < dfComputedRes / dfRes)
402 15 : nZoomLevel--;
403 : }
404 : }
405 : }
406 26 : CPLDebug("COG", "Using ZOOM_LEVEL %d", nZoomLevel);
407 26 : dfRes = tmList[nZoomLevel].mResX * tmList[0].mTileWidth / nBlockSize;
408 :
409 : const double dfOriX =
410 26 : bInvertAxis ? tmList[0].mTopLeftY : tmList[0].mTopLeftX;
411 : const double dfOriY =
412 26 : bInvertAxis ? tmList[0].mTopLeftX : tmList[0].mTopLeftY;
413 26 : const double dfTileExtent = dfRes * nBlockSize;
414 26 : constexpr double TOLERANCE_IN_PIXEL = 0.499;
415 26 : const double dfEps = TOLERANCE_IN_PIXEL * dfRes;
416 26 : int nTLTileX = static_cast<int>(
417 26 : std::floor((dfMinX - dfOriX + dfEps) / dfTileExtent));
418 26 : int nTLTileY = static_cast<int>(
419 26 : std::floor((dfOriY - dfMaxY + dfEps) / dfTileExtent));
420 26 : int nBRTileX = static_cast<int>(
421 26 : std::ceil((dfMaxX - dfOriX - dfEps) / dfTileExtent));
422 26 : int nBRTileY = static_cast<int>(
423 26 : std::ceil((dfOriY - dfMinY - dfEps) / dfTileExtent));
424 :
425 26 : nAlignedLevels =
426 26 : std::min(std::min(10, atoi(CSLFetchNameValueDef(
427 : papszOptions, "ALIGNED_LEVELS", "0"))),
428 26 : nZoomLevel);
429 26 : int nAccDivisor = 1;
430 34 : for (int i = 0; i < nAlignedLevels - 1; i++)
431 : {
432 8 : const int nCurLevel = nZoomLevel - i;
433 : const double dfResRatio =
434 8 : tmList[nCurLevel - 1].mResX / tmList[nCurLevel].mResX;
435 : // Magical number that has a great number of divisors
436 : // For example if previous scale denom was 50K and current one
437 : // is 20K, then dfResRatio = 2.5 and dfScaledInvResRatio = 24
438 : // We must then simplify 60 / 24 as 5 / 2, and make sure to
439 : // align tile coordinates on multiple of the 5 numerator
440 8 : constexpr int MAGICAL = 60;
441 8 : const double dfScaledInvResRatio = MAGICAL / dfResRatio;
442 16 : if (dfScaledInvResRatio < 1 || dfScaledInvResRatio > 60 ||
443 8 : std::abs(std::round(dfScaledInvResRatio) -
444 : dfScaledInvResRatio) > 1e-10)
445 : {
446 0 : CPLError(CE_Failure, CPLE_AppDefined,
447 : "Unsupported ratio of resolution for "
448 : "ALIGNED_LEVELS between zoom level %d and %d = %g",
449 : nCurLevel - 1, nCurLevel, dfResRatio);
450 0 : return false;
451 : }
452 8 : const int nScaledInvResRatio =
453 8 : static_cast<int>(std::round(dfScaledInvResRatio));
454 8 : int nNumerator = 0;
455 32 : for (int nDivisor = nScaledInvResRatio; nDivisor >= 2; --nDivisor)
456 : {
457 32 : if ((MAGICAL % nDivisor) == 0 &&
458 12 : (nScaledInvResRatio % nDivisor) == 0)
459 : {
460 8 : nNumerator = MAGICAL / nDivisor;
461 8 : break;
462 : }
463 : }
464 8 : if (nNumerator == 0)
465 : {
466 0 : CPLError(CE_Failure, CPLE_AppDefined,
467 : "Unsupported ratio of resolution for "
468 : "ALIGNED_LEVELS between zoom level %d and %d = %g",
469 : nCurLevel - 1, nCurLevel, dfResRatio);
470 0 : return false;
471 : }
472 8 : nAccDivisor *= nNumerator;
473 : }
474 26 : if (nAccDivisor > 1)
475 : {
476 4 : nTLTileX = (nTLTileX / nAccDivisor) * nAccDivisor;
477 4 : nTLTileY = (nTLTileY / nAccDivisor) * nAccDivisor;
478 4 : nBRTileY =
479 4 : ((nBRTileY + nAccDivisor - 1) / nAccDivisor) * nAccDivisor;
480 4 : nBRTileX =
481 4 : ((nBRTileX + nAccDivisor - 1) / nAccDivisor) * nAccDivisor;
482 : }
483 :
484 26 : if (nTLTileX < 0 || nTLTileY < 0 ||
485 74 : nBRTileX > tmList[nZoomLevel].mMatrixWidth ||
486 22 : nBRTileY > tmList[nZoomLevel].mMatrixHeight)
487 : {
488 4 : CPLError(CE_Warning, CPLE_AppDefined,
489 : "Raster extent partially outside of tile matrix "
490 : "bounding box. Clamping it to it");
491 : }
492 26 : nTLTileX = std::max(0, nTLTileX);
493 26 : nTLTileY = std::max(0, nTLTileY);
494 26 : nBRTileX = std::min(tmList[nZoomLevel].mMatrixWidth, nBRTileX);
495 26 : nBRTileY = std::min(tmList[nZoomLevel].mMatrixHeight, nBRTileY);
496 :
497 26 : dfMinX = dfOriX + nTLTileX * dfTileExtent;
498 26 : dfMinY = dfOriY - nBRTileY * dfTileExtent;
499 26 : dfMaxX = dfOriX + nBRTileX * dfTileExtent;
500 26 : dfMaxY = dfOriY - nTLTileY * dfTileExtent;
501 : }
502 7 : else if (!osExtent.empty() || !osRes.empty())
503 : {
504 1 : CPLStringList aosTokens(CSLTokenizeString2(osExtent, ",", 0));
505 1 : if (aosTokens.size() != 4)
506 : {
507 0 : CPLError(CE_Failure, CPLE_AppDefined, "Invalid value for EXTENT");
508 0 : return false;
509 : }
510 1 : dfMinX = CPLAtof(aosTokens[0]);
511 1 : dfMinY = CPLAtof(aosTokens[1]);
512 1 : dfMaxX = CPLAtof(aosTokens[2]);
513 1 : dfMaxY = CPLAtof(aosTokens[3]);
514 1 : if (!osRes.empty())
515 1 : dfRes = CPLAtof(osRes);
516 : }
517 :
518 33 : nXSize = static_cast<int>(std::round((dfMaxX - dfMinX) / dfRes));
519 33 : nYSize = static_cast<int>(std::round((dfMaxY - dfMinY) / dfRes));
520 :
521 33 : osResampling = COGGetResampling(poSrcDS, papszOptions);
522 :
523 33 : return true;
524 : }
525 :
526 : // Used by gdalwarp
527 3 : bool COGGetWarpingCharacteristics(GDALDataset *poSrcDS,
528 : const char *const *papszOptions,
529 : CPLString &osResampling,
530 : CPLString &osTargetSRS, int &nXSize,
531 : int &nYSize, double &dfMinX, double &dfMinY,
532 : double &dfMaxX, double &dfMaxY)
533 : {
534 3 : std::unique_ptr<gdal::TileMatrixSet> poTM;
535 3 : int nZoomLevel = 0;
536 3 : int nAlignedLevels = 0;
537 : double dfRes;
538 3 : return COGGetWarpingCharacteristics(poSrcDS, papszOptions, osResampling,
539 : osTargetSRS, nXSize, nYSize, dfMinX,
540 : dfMinY, dfMaxX, dfMaxY, dfRes, poTM,
541 6 : nZoomLevel, nAlignedLevels);
542 : }
543 :
544 : /************************************************************************/
545 : /* COGHasWarpingOptions() */
546 : /************************************************************************/
547 :
548 140 : bool COGHasWarpingOptions(CSLConstList papszOptions)
549 : {
550 273 : return CSLFetchNameValue(papszOptions, "TARGET_SRS") != nullptr ||
551 273 : !EQUAL(CSLFetchNameValueDef(papszOptions, "TILING_SCHEME", "CUSTOM"),
552 : "CUSTOM");
553 : }
554 :
555 : /************************************************************************/
556 : /* COGRemoveWarpingOptions() */
557 : /************************************************************************/
558 :
559 2 : void COGRemoveWarpingOptions(CPLStringList &aosOptions)
560 : {
561 2 : aosOptions.SetNameValue("TARGET_SRS", nullptr);
562 2 : aosOptions.SetNameValue("TILING_SCHEME", nullptr);
563 2 : aosOptions.SetNameValue("EXTENT", nullptr);
564 2 : aosOptions.SetNameValue("RES", nullptr);
565 2 : aosOptions.SetNameValue("ALIGNED_LEVELS", nullptr);
566 2 : aosOptions.SetNameValue("ZOOM_LEVEL_STRATEGY", nullptr);
567 2 : }
568 :
569 : /************************************************************************/
570 : /* CreateReprojectedDS() */
571 : /************************************************************************/
572 :
573 26 : static std::unique_ptr<GDALDataset> CreateReprojectedDS(
574 : const char *pszDstFilename, GDALDataset *poSrcDS,
575 : const char *const *papszOptions, const CPLString &osResampling,
576 : const CPLString &osTargetSRS, const int nXSize, const int nYSize,
577 : const double dfMinX, const double dfMinY, const double dfMaxX,
578 : const double dfMaxY, const double dfRes, GDALProgressFunc pfnProgress,
579 : void *pProgressData, double &dfCurPixels, double &dfTotalPixelsToProcess)
580 : {
581 26 : char **papszArg = nullptr;
582 : // We could have done a warped VRT, but overview building on it might be
583 : // slow, so materialize as GTiff
584 26 : papszArg = CSLAddString(papszArg, "-of");
585 26 : papszArg = CSLAddString(papszArg, "GTiff");
586 26 : papszArg = CSLAddString(papszArg, "-co");
587 26 : papszArg = CSLAddString(papszArg, "TILED=YES");
588 26 : papszArg = CSLAddString(papszArg, "-co");
589 26 : papszArg = CSLAddString(papszArg, "SPARSE_OK=YES");
590 26 : const char *pszBIGTIFF = CSLFetchNameValue(papszOptions, "BIGTIFF");
591 26 : if (pszBIGTIFF)
592 : {
593 0 : papszArg = CSLAddString(papszArg, "-co");
594 0 : papszArg = CSLAddString(papszArg,
595 0 : (CPLString("BIGTIFF=") + pszBIGTIFF).c_str());
596 : }
597 26 : papszArg = CSLAddString(papszArg, "-co");
598 26 : papszArg = CSLAddString(papszArg, HasZSTDCompression() ? "COMPRESS=ZSTD"
599 : : "COMPRESS=LZW");
600 26 : papszArg = CSLAddString(papszArg, "-t_srs");
601 26 : papszArg = CSLAddString(papszArg, osTargetSRS);
602 26 : papszArg = CSLAddString(papszArg, "-te");
603 26 : papszArg = CSLAddString(papszArg, CPLSPrintf("%.17g", dfMinX));
604 26 : papszArg = CSLAddString(papszArg, CPLSPrintf("%.17g", dfMinY));
605 26 : papszArg = CSLAddString(papszArg, CPLSPrintf("%.17g", dfMaxX));
606 26 : papszArg = CSLAddString(papszArg, CPLSPrintf("%.17g", dfMaxY));
607 26 : papszArg = CSLAddString(papszArg, "-ts");
608 26 : papszArg = CSLAddString(papszArg, CPLSPrintf("%d", nXSize));
609 26 : papszArg = CSLAddString(papszArg, CPLSPrintf("%d", nYSize));
610 :
611 : // to be kept in sync with gdalwarp_lib.cpp
612 26 : constexpr double RELATIVE_ERROR_RES_SHARED_BY_COG_AND_GDALWARP = 1e-8;
613 26 : if (fabs((dfMaxX - dfMinX) / dfRes - nXSize) <=
614 26 : RELATIVE_ERROR_RES_SHARED_BY_COG_AND_GDALWARP &&
615 26 : fabs((dfMaxY - dfMinY) / dfRes - nYSize) <=
616 : RELATIVE_ERROR_RES_SHARED_BY_COG_AND_GDALWARP)
617 : {
618 : // Try to produce exactly square pixels
619 26 : papszArg = CSLAddString(papszArg, "-tr");
620 26 : papszArg = CSLAddString(papszArg, CPLSPrintf("%.17g", dfRes));
621 26 : papszArg = CSLAddString(papszArg, CPLSPrintf("%.17g", dfRes));
622 : }
623 : else
624 : {
625 0 : CPLDebug("COG", "Cannot pass -tr option to GDALWarp() due to extent, "
626 : "size and resolution not consistent enough");
627 : }
628 :
629 26 : int bHasNoData = FALSE;
630 26 : poSrcDS->GetRasterBand(1)->GetNoDataValue(&bHasNoData);
631 52 : if (!bHasNoData &&
632 26 : CPLTestBool(CSLFetchNameValueDef(papszOptions, "ADD_ALPHA", "YES")))
633 : {
634 26 : papszArg = CSLAddString(papszArg, "-dstalpha");
635 : }
636 26 : papszArg = CSLAddString(papszArg, "-r");
637 26 : papszArg = CSLAddString(papszArg, osResampling);
638 26 : papszArg = CSLAddString(papszArg, "-wo");
639 26 : papszArg = CSLAddString(papszArg, "SAMPLE_GRID=YES");
640 26 : const char *pszNumThreads = CSLFetchNameValue(papszOptions, "NUM_THREADS");
641 26 : if (pszNumThreads)
642 : {
643 0 : papszArg = CSLAddString(papszArg, "-wo");
644 0 : papszArg = CSLAddString(
645 0 : papszArg, (CPLString("NUM_THREADS=") + pszNumThreads).c_str());
646 : }
647 :
648 26 : const auto poFirstBand = poSrcDS->GetRasterBand(1);
649 26 : const bool bHasMask = poFirstBand->GetMaskFlags() == GMF_PER_DATASET;
650 :
651 26 : const int nBands = poSrcDS->GetRasterCount();
652 : const char *pszOverviews =
653 26 : CSLFetchNameValueDef(papszOptions, "OVERVIEWS", "AUTO");
654 52 : const bool bUseExistingOrNone = EQUAL(pszOverviews, "FORCE_USE_EXISTING") ||
655 26 : EQUAL(pszOverviews, "NONE");
656 26 : dfTotalPixelsToProcess =
657 26 : double(nXSize) * nYSize * (nBands + (bHasMask ? 1 : 0)) +
658 26 : ((bHasMask && !bUseExistingOrNone) ? double(nXSize) * nYSize / 3 : 0) +
659 26 : (!bUseExistingOrNone ? double(nXSize) * nYSize * nBands / 3 : 0) +
660 26 : double(nXSize) * nYSize * (nBands + (bHasMask ? 1 : 0)) * 4. / 3;
661 :
662 26 : auto psOptions = GDALWarpAppOptionsNew(papszArg, nullptr);
663 26 : CSLDestroy(papszArg);
664 26 : if (psOptions == nullptr)
665 1 : return nullptr;
666 :
667 25 : const double dfNextPixels =
668 25 : double(nXSize) * nYSize * (nBands + (bHasMask ? 1 : 0));
669 50 : void *pScaledProgress = GDALCreateScaledProgress(
670 25 : dfCurPixels / dfTotalPixelsToProcess,
671 25 : dfNextPixels / dfTotalPixelsToProcess, pfnProgress, pProgressData);
672 25 : dfCurPixels = dfNextPixels;
673 :
674 25 : CPLDebug("COG", "Reprojecting source dataset: start");
675 25 : GDALWarpAppOptionsSetProgress(psOptions, GDALScaledProgress,
676 : pScaledProgress);
677 50 : CPLString osTmpFile(GetTmpFilename(pszDstFilename, "warped.tif.tmp"));
678 25 : auto hSrcDS = GDALDataset::ToHandle(poSrcDS);
679 :
680 25 : std::unique_ptr<CPLConfigOptionSetter> poWarpThreadSetter;
681 25 : if (pszNumThreads)
682 : {
683 0 : poWarpThreadSetter.reset(new CPLConfigOptionSetter(
684 0 : "GDAL_NUM_THREADS", pszNumThreads, false));
685 : }
686 :
687 25 : auto hRet = GDALWarp(osTmpFile, nullptr, 1, &hSrcDS, psOptions, nullptr);
688 25 : GDALWarpAppOptionsFree(psOptions);
689 25 : CPLDebug("COG", "Reprojecting source dataset: end");
690 :
691 25 : GDALDestroyScaledProgress(pScaledProgress);
692 :
693 25 : return std::unique_ptr<GDALDataset>(GDALDataset::FromHandle(hRet));
694 : }
695 :
696 : /************************************************************************/
697 : /* GDALCOGCreator */
698 : /************************************************************************/
699 :
700 : struct GDALCOGCreator final
701 : {
702 : std::unique_ptr<GDALDataset> m_poReprojectedDS{};
703 : std::unique_ptr<GDALDataset> m_poRGBMaskDS{};
704 : std::unique_ptr<GDALDataset> m_poVRTWithOrWithoutStats{};
705 : CPLString m_osTmpOverviewFilename{};
706 : CPLString m_osTmpMskOverviewFilename{};
707 :
708 : ~GDALCOGCreator();
709 :
710 : GDALDataset *Create(const char *pszFilename, GDALDataset *const poSrcDS,
711 : char **papszOptions, GDALProgressFunc pfnProgress,
712 : void *pProgressData);
713 : };
714 :
715 : /************************************************************************/
716 : /* GDALCOGCreator::~GDALCOGCreator() */
717 : /************************************************************************/
718 :
719 138 : GDALCOGCreator::~GDALCOGCreator()
720 : {
721 : // Destroy m_poRGBMaskDS before m_poReprojectedDS since the former
722 : // may reference the later
723 138 : m_poRGBMaskDS.reset();
724 :
725 : // Config option just for testing purposes
726 : const bool bDeleteTempFiles =
727 138 : CPLTestBool(CPLGetConfigOption("COG_DELETE_TEMP_FILES", "YES"));
728 138 : if (bDeleteTempFiles)
729 : {
730 137 : if (m_poReprojectedDS)
731 : {
732 50 : CPLString osProjectedDSName(m_poReprojectedDS->GetDescription());
733 25 : m_poReprojectedDS.reset();
734 25 : VSIUnlink(osProjectedDSName);
735 : }
736 137 : if (!m_osTmpOverviewFilename.empty())
737 : {
738 32 : VSIUnlink(m_osTmpOverviewFilename);
739 : }
740 137 : if (!m_osTmpMskOverviewFilename.empty())
741 : {
742 7 : VSIUnlink(m_osTmpMskOverviewFilename);
743 : }
744 : }
745 138 : }
746 :
747 : /************************************************************************/
748 : /* GDALCOGCreator::Create() */
749 : /************************************************************************/
750 :
751 138 : GDALDataset *GDALCOGCreator::Create(const char *pszFilename,
752 : GDALDataset *const poSrcDS,
753 : char **papszOptions,
754 : GDALProgressFunc pfnProgress,
755 : void *pProgressData)
756 : {
757 138 : if (pfnProgress == nullptr)
758 0 : pfnProgress = GDALDummyProgress;
759 :
760 138 : if (poSrcDS->GetRasterCount() == 0)
761 : {
762 1 : CPLError(CE_Failure, CPLE_NotSupported,
763 : "COG driver does not support 0-band source raster");
764 1 : return nullptr;
765 : }
766 :
767 : CPLConfigOptionSetter oSetterReportDirtyBlockFlushing(
768 274 : "GDAL_REPORT_DIRTY_BLOCK_FLUSHING", "NO", true);
769 :
770 : const char *pszStatistics =
771 137 : CSLFetchNameValueDef(papszOptions, "STATISTICS", "AUTO");
772 137 : auto poSrcFirstBand = poSrcDS->GetRasterBand(1);
773 : const bool bSrcHasStatistics =
774 137 : poSrcFirstBand->GetMetadataItem("STATISTICS_MINIMUM") &&
775 8 : poSrcFirstBand->GetMetadataItem("STATISTICS_MAXIMUM") &&
776 153 : poSrcFirstBand->GetMetadataItem("STATISTICS_MEAN") &&
777 8 : poSrcFirstBand->GetMetadataItem("STATISTICS_STDDEV");
778 137 : bool bNeedStats = false;
779 137 : bool bRemoveStats = false;
780 137 : bool bWrkHasStatistics = bSrcHasStatistics;
781 137 : if (EQUAL(pszStatistics, "AUTO"))
782 : {
783 : // nothing
784 : }
785 10 : else if (CPLTestBool(pszStatistics))
786 : {
787 5 : bNeedStats = true;
788 : }
789 : else
790 : {
791 5 : bRemoveStats = true;
792 : }
793 :
794 137 : double dfCurPixels = 0;
795 137 : double dfTotalPixelsToProcess = 0;
796 137 : GDALDataset *poCurDS = poSrcDS;
797 :
798 137 : std::unique_ptr<gdal::TileMatrixSet> poTM;
799 137 : int nZoomLevel = 0;
800 137 : int nAlignedLevels = 0;
801 137 : if (COGHasWarpingOptions(papszOptions))
802 : {
803 33 : CPLString osTargetResampling;
804 33 : CPLString osTargetSRS;
805 33 : int nTargetXSize = 0;
806 33 : int nTargetYSize = 0;
807 33 : double dfTargetMinX = 0;
808 33 : double dfTargetMinY = 0;
809 33 : double dfTargetMaxX = 0;
810 33 : double dfTargetMaxY = 0;
811 33 : double dfRes = 0;
812 33 : if (!COGGetWarpingCharacteristics(
813 : poCurDS, papszOptions, osTargetResampling, osTargetSRS,
814 : nTargetXSize, nTargetYSize, dfTargetMinX, dfTargetMinY,
815 : dfTargetMaxX, dfTargetMaxY, dfRes, poTM, nZoomLevel,
816 : nAlignedLevels))
817 : {
818 2 : return nullptr;
819 : }
820 :
821 : // Collect information on source dataset to see if it already
822 : // matches the warping specifications
823 31 : CPLString osSrcSRS;
824 31 : const auto poSrcSRS = poCurDS->GetSpatialRef();
825 31 : if (poSrcSRS)
826 : {
827 31 : const char *pszAuthName = poSrcSRS->GetAuthorityName(nullptr);
828 31 : const char *pszAuthCode = poSrcSRS->GetAuthorityCode(nullptr);
829 31 : if (pszAuthName && pszAuthCode)
830 : {
831 31 : osSrcSRS = pszAuthName;
832 31 : osSrcSRS += ':';
833 31 : osSrcSRS += pszAuthCode;
834 : }
835 : }
836 31 : double dfSrcMinX = 0;
837 31 : double dfSrcMinY = 0;
838 31 : double dfSrcMaxX = 0;
839 31 : double dfSrcMaxY = 0;
840 : double adfSrcGT[6];
841 31 : const int nSrcXSize = poCurDS->GetRasterXSize();
842 31 : const int nSrcYSize = poCurDS->GetRasterYSize();
843 31 : if (poCurDS->GetGeoTransform(adfSrcGT) == CE_None)
844 : {
845 31 : dfSrcMinX = adfSrcGT[0];
846 31 : dfSrcMaxY = adfSrcGT[3];
847 31 : dfSrcMaxX = adfSrcGT[0] + nSrcXSize * adfSrcGT[1];
848 31 : dfSrcMinY = adfSrcGT[3] + nSrcYSize * adfSrcGT[5];
849 : }
850 :
851 24 : if (nTargetXSize == nSrcXSize && nTargetYSize == nSrcYSize &&
852 12 : osTargetSRS == osSrcSRS &&
853 6 : fabs(dfSrcMinX - dfTargetMinX) < 1e-10 * fabs(dfSrcMinX) &&
854 5 : fabs(dfSrcMinY - dfTargetMinY) < 1e-10 * fabs(dfSrcMinY) &&
855 48 : fabs(dfSrcMaxX - dfTargetMaxX) < 1e-10 * fabs(dfSrcMaxX) &&
856 5 : fabs(dfSrcMaxY - dfTargetMaxY) < 1e-10 * fabs(dfSrcMaxY))
857 : {
858 5 : CPLDebug("COG",
859 : "Skipping reprojection step: "
860 : "source dataset matches reprojection specifications");
861 : }
862 : else
863 : {
864 52 : m_poReprojectedDS = CreateReprojectedDS(
865 : pszFilename, poCurDS, papszOptions, osTargetResampling,
866 : osTargetSRS, nTargetXSize, nTargetYSize, dfTargetMinX,
867 : dfTargetMinY, dfTargetMaxX, dfTargetMaxY, dfRes, pfnProgress,
868 26 : pProgressData, dfCurPixels, dfTotalPixelsToProcess);
869 26 : if (!m_poReprojectedDS)
870 1 : return nullptr;
871 25 : poCurDS = m_poReprojectedDS.get();
872 :
873 25 : if (bSrcHasStatistics && !bNeedStats && !bRemoveStats)
874 : {
875 1 : bNeedStats = true;
876 : }
877 25 : bWrkHasStatistics = false;
878 : }
879 : }
880 :
881 : CPLString osCompress = CSLFetchNameValueDef(papszOptions, "COMPRESS",
882 268 : gbHasLZW ? "LZW" : "NONE");
883 134 : if (EQUAL(osCompress, "JPEG") &&
884 144 : (poCurDS->GetRasterCount() == 2 || poCurDS->GetRasterCount() == 4) &&
885 10 : poCurDS->GetRasterBand(poCurDS->GetRasterCount())
886 10 : ->GetColorInterpretation() == GCI_AlphaBand)
887 : {
888 10 : char **papszArg = nullptr;
889 10 : papszArg = CSLAddString(papszArg, "-of");
890 10 : papszArg = CSLAddString(papszArg, "VRT");
891 10 : papszArg = CSLAddString(papszArg, "-b");
892 10 : papszArg = CSLAddString(papszArg, "1");
893 10 : if (poCurDS->GetRasterCount() == 2)
894 : {
895 1 : papszArg = CSLAddString(papszArg, "-mask");
896 1 : papszArg = CSLAddString(papszArg, "2");
897 : }
898 : else
899 : {
900 9 : CPLAssert(poCurDS->GetRasterCount() == 4);
901 9 : papszArg = CSLAddString(papszArg, "-b");
902 9 : papszArg = CSLAddString(papszArg, "2");
903 9 : papszArg = CSLAddString(papszArg, "-b");
904 9 : papszArg = CSLAddString(papszArg, "3");
905 9 : papszArg = CSLAddString(papszArg, "-mask");
906 9 : papszArg = CSLAddString(papszArg, "4");
907 : }
908 : GDALTranslateOptions *psOptions =
909 10 : GDALTranslateOptionsNew(papszArg, nullptr);
910 10 : CSLDestroy(papszArg);
911 10 : GDALDatasetH hRGBMaskDS = GDALTranslate(
912 : "", GDALDataset::ToHandle(poCurDS), psOptions, nullptr);
913 10 : GDALTranslateOptionsFree(psOptions);
914 10 : if (!hRGBMaskDS)
915 : {
916 0 : return nullptr;
917 : }
918 10 : m_poRGBMaskDS.reset(GDALDataset::FromHandle(hRGBMaskDS));
919 10 : poCurDS = m_poRGBMaskDS.get();
920 :
921 10 : if (bSrcHasStatistics && !bNeedStats && !bRemoveStats)
922 : {
923 1 : bNeedStats = true;
924 : }
925 9 : else if (bRemoveStats && bWrkHasStatistics)
926 : {
927 1 : poCurDS->ClearStatistics();
928 1 : bRemoveStats = false;
929 : }
930 : }
931 :
932 134 : const int nBands = poCurDS->GetRasterCount();
933 134 : const int nXSize = poCurDS->GetRasterXSize();
934 134 : const int nYSize = poCurDS->GetRasterYSize();
935 :
936 8 : const auto CreateVRTWithOrWithoutStats = [this, &poCurDS]()
937 : {
938 2 : const char *const apszOptions[] = {"-of", "VRT", nullptr};
939 : GDALTranslateOptions *psOptions =
940 2 : GDALTranslateOptionsNew(const_cast<char **>(apszOptions), nullptr);
941 2 : GDALDatasetH hVRTDS = GDALTranslate("", GDALDataset::ToHandle(poCurDS),
942 : psOptions, nullptr);
943 2 : GDALTranslateOptionsFree(psOptions);
944 2 : if (!hVRTDS)
945 0 : return false;
946 2 : m_poVRTWithOrWithoutStats.reset(GDALDataset::FromHandle(hVRTDS));
947 2 : poCurDS = m_poVRTWithOrWithoutStats.get();
948 2 : return true;
949 134 : };
950 :
951 134 : if (bNeedStats && !bWrkHasStatistics)
952 : {
953 5 : if (poSrcDS == poCurDS && !CreateVRTWithOrWithoutStats())
954 : {
955 0 : return nullptr;
956 : }
957 :
958 : // Avoid source files to be modified
959 : CPLConfigOptionSetter enablePamDirtyDisabler(
960 10 : "GDAL_PAM_ENABLE_MARK_DIRTY", "NO", true);
961 :
962 15 : for (int i = 1; i <= nBands; ++i)
963 : {
964 10 : poCurDS->GetRasterBand(i)->ComputeStatistics(
965 : /*bApproxOK=*/FALSE, nullptr, nullptr, nullptr, nullptr,
966 10 : nullptr, nullptr);
967 5 : }
968 : }
969 129 : else if (bRemoveStats && bWrkHasStatistics)
970 : {
971 1 : if (!CreateVRTWithOrWithoutStats())
972 0 : return nullptr;
973 :
974 1 : m_poVRTWithOrWithoutStats->ClearStatistics();
975 : }
976 :
977 268 : CPLString osBlockSize(CSLFetchNameValueDef(papszOptions, "BLOCKSIZE", ""));
978 134 : if (osBlockSize.empty())
979 : {
980 114 : if (poTM)
981 : {
982 20 : osBlockSize.Printf("%d", poTM->tileMatrixList()[0].mTileWidth);
983 : }
984 : else
985 : {
986 94 : osBlockSize = "512";
987 : }
988 : }
989 :
990 134 : const int nOvrThresholdSize = atoi(osBlockSize);
991 :
992 134 : const auto poFirstBand = poCurDS->GetRasterBand(1);
993 134 : const bool bHasMask = poFirstBand->GetMaskFlags() == GMF_PER_DATASET;
994 :
995 : CPLString osOverviews =
996 268 : CSLFetchNameValueDef(papszOptions, "OVERVIEWS", "AUTO");
997 : const bool bUseExistingOrNone =
998 134 : EQUAL(osOverviews, "FORCE_USE_EXISTING") || EQUAL(osOverviews, "NONE");
999 :
1000 : const int nOverviewCount =
1001 134 : atoi(CSLFetchNameValueDef(papszOptions, "OVERVIEW_COUNT", "-1"));
1002 :
1003 : const bool bGenerateMskOvr =
1004 130 : !bUseExistingOrNone && bHasMask &&
1005 7 : (nXSize > nOvrThresholdSize || nYSize > nOvrThresholdSize ||
1006 264 : nOverviewCount > 0) &&
1007 8 : (EQUAL(osOverviews, "IGNORE_EXISTING") ||
1008 8 : poFirstBand->GetMaskBand()->GetOverviewCount() == 0);
1009 : const bool bGenerateOvr =
1010 130 : !bUseExistingOrNone &&
1011 96 : (nXSize > nOvrThresholdSize || nYSize > nOvrThresholdSize ||
1012 264 : nOverviewCount > 0) &&
1013 38 : (EQUAL(osOverviews, "IGNORE_EXISTING") ||
1014 36 : poFirstBand->GetOverviewCount() == 0);
1015 :
1016 268 : std::vector<std::pair<int, int>> asOverviewDims;
1017 134 : int nTmpXSize = nXSize;
1018 134 : int nTmpYSize = nYSize;
1019 134 : if (poTM)
1020 : {
1021 23 : const auto &tmList = poTM->tileMatrixList();
1022 23 : int nCurLevel = nZoomLevel;
1023 : while (true)
1024 : {
1025 47 : if (nOverviewCount < 0)
1026 : {
1027 33 : if (nTmpXSize <= nOvrThresholdSize &&
1028 21 : nTmpYSize <= nOvrThresholdSize)
1029 20 : break;
1030 : }
1031 14 : else if (static_cast<int>(asOverviewDims.size()) ==
1032 26 : nOverviewCount ||
1033 12 : (nTmpXSize == 1 && nTmpYSize == 1))
1034 : {
1035 3 : break;
1036 : }
1037 : const double dfResRatio =
1038 : (nCurLevel >= 1)
1039 24 : ? tmList[nCurLevel - 1].mResX / tmList[nCurLevel].mResX
1040 24 : : 2;
1041 24 : nTmpXSize = static_cast<int>(nTmpXSize / dfResRatio + 0.5);
1042 24 : nTmpYSize = static_cast<int>(nTmpYSize / dfResRatio + 0.5);
1043 24 : if (nTmpXSize == 0)
1044 0 : nTmpXSize = 1;
1045 24 : if (nTmpYSize == 0)
1046 0 : nTmpYSize = 1;
1047 24 : asOverviewDims.emplace_back(std::pair(nTmpXSize, nTmpYSize));
1048 24 : nCurLevel--;
1049 24 : }
1050 : }
1051 111 : else if (bGenerateMskOvr || bGenerateOvr)
1052 : {
1053 27 : if (!bGenerateOvr)
1054 : {
1055 : // If generating only .msk.ovr, use the exact overview size as
1056 : // the overviews of the imagery.
1057 1 : int nIters = poFirstBand->GetOverviewCount();
1058 1 : if (nOverviewCount >= 0 && nOverviewCount < nIters)
1059 0 : nIters = nOverviewCount;
1060 2 : for (int i = 0; i < nIters; i++)
1061 : {
1062 1 : auto poOvrBand = poFirstBand->GetOverview(i);
1063 : asOverviewDims.emplace_back(
1064 1 : std::pair(poOvrBand->GetXSize(), poOvrBand->GetYSize()));
1065 : }
1066 : }
1067 : else
1068 : {
1069 : while (true)
1070 : {
1071 76 : if (nOverviewCount < 0)
1072 : {
1073 58 : if (nTmpXSize <= nOvrThresholdSize &&
1074 22 : nTmpYSize <= nOvrThresholdSize)
1075 22 : break;
1076 : }
1077 18 : else if (static_cast<int>(asOverviewDims.size()) ==
1078 33 : nOverviewCount ||
1079 15 : (nTmpXSize == 1 && nTmpYSize == 1))
1080 : {
1081 4 : break;
1082 : }
1083 50 : nTmpXSize /= 2;
1084 50 : nTmpYSize /= 2;
1085 50 : if (nTmpXSize == 0)
1086 0 : nTmpXSize = 1;
1087 50 : if (nTmpYSize == 0)
1088 1 : nTmpYSize = 1;
1089 50 : asOverviewDims.emplace_back(std::pair(nTmpXSize, nTmpYSize));
1090 : }
1091 : }
1092 : }
1093 :
1094 134 : if (dfTotalPixelsToProcess == 0.0)
1095 : {
1096 109 : dfTotalPixelsToProcess =
1097 109 : (bGenerateMskOvr ? double(nXSize) * nYSize / 3 : 0) +
1098 109 : (bGenerateOvr ? double(nXSize) * nYSize * nBands / 3 : 0) +
1099 109 : double(nXSize) * nYSize * (nBands + (bHasMask ? 1 : 0)) * 4. / 3;
1100 : }
1101 :
1102 268 : CPLStringList aosOverviewOptions;
1103 : aosOverviewOptions.SetNameValue(
1104 : "COMPRESS",
1105 : CPLGetConfigOption("COG_TMP_COMPRESSION", // only for debug purposes
1106 134 : HasZSTDCompression() ? "ZSTD" : "LZW"));
1107 : aosOverviewOptions.SetNameValue(
1108 134 : "NUM_THREADS", CSLFetchNameValue(papszOptions, "NUM_THREADS"));
1109 134 : aosOverviewOptions.SetNameValue("BIGTIFF", "YES");
1110 134 : aosOverviewOptions.SetNameValue("SPARSE_OK", "YES");
1111 :
1112 134 : if (bGenerateMskOvr)
1113 : {
1114 8 : CPLDebug("COG", "Generating overviews of the mask: start");
1115 8 : m_osTmpMskOverviewFilename = GetTmpFilename(pszFilename, "msk.ovr.tmp");
1116 8 : GDALRasterBand *poSrcMask = poFirstBand->GetMaskBand();
1117 8 : const char *pszResampling = CSLFetchNameValueDef(
1118 : papszOptions, "OVERVIEW_RESAMPLING",
1119 : CSLFetchNameValueDef(papszOptions, "RESAMPLING",
1120 : GetResampling(poSrcDS)));
1121 :
1122 8 : double dfNextPixels = dfCurPixels + double(nXSize) * nYSize / 3;
1123 8 : void *pScaledProgress = GDALCreateScaledProgress(
1124 : dfCurPixels / dfTotalPixelsToProcess,
1125 : dfNextPixels / dfTotalPixelsToProcess, pfnProgress, pProgressData);
1126 8 : dfCurPixels = dfNextPixels;
1127 :
1128 8 : CPLErr eErr = GTIFFBuildOverviewsEx(
1129 : m_osTmpMskOverviewFilename, 1, &poSrcMask,
1130 8 : static_cast<int>(asOverviewDims.size()), nullptr,
1131 8 : asOverviewDims.data(), pszResampling, aosOverviewOptions.List(),
1132 : GDALScaledProgress, pScaledProgress);
1133 8 : CPLDebug("COG", "Generating overviews of the mask: end");
1134 :
1135 8 : GDALDestroyScaledProgress(pScaledProgress);
1136 8 : if (eErr != CE_None)
1137 : {
1138 0 : return nullptr;
1139 : }
1140 : }
1141 :
1142 134 : if (bGenerateOvr)
1143 : {
1144 33 : CPLDebug("COG", "Generating overviews of the imagery: start");
1145 33 : m_osTmpOverviewFilename = GetTmpFilename(pszFilename, "ovr.tmp");
1146 33 : std::vector<GDALRasterBand *> apoSrcBands;
1147 101 : for (int i = 0; i < nBands; i++)
1148 68 : apoSrcBands.push_back(poCurDS->GetRasterBand(i + 1));
1149 33 : const char *pszResampling = CSLFetchNameValueDef(
1150 : papszOptions, "OVERVIEW_RESAMPLING",
1151 : CSLFetchNameValueDef(papszOptions, "RESAMPLING",
1152 : GetResampling(poSrcDS)));
1153 :
1154 33 : double dfNextPixels =
1155 33 : dfCurPixels + double(nXSize) * nYSize * nBands / 3;
1156 33 : void *pScaledProgress = GDALCreateScaledProgress(
1157 : dfCurPixels / dfTotalPixelsToProcess,
1158 : dfNextPixels / dfTotalPixelsToProcess, pfnProgress, pProgressData);
1159 33 : dfCurPixels = dfNextPixels;
1160 :
1161 33 : if (nBands > 1)
1162 : {
1163 18 : aosOverviewOptions.SetNameValue("INTERLEAVE", "PIXEL");
1164 : }
1165 33 : if (!m_osTmpMskOverviewFilename.empty())
1166 : {
1167 : aosOverviewOptions.SetNameValue("MASK_OVERVIEW_DATASET",
1168 7 : m_osTmpMskOverviewFilename);
1169 : }
1170 33 : CPLErr eErr = GTIFFBuildOverviewsEx(
1171 33 : m_osTmpOverviewFilename, nBands, &apoSrcBands[0],
1172 33 : static_cast<int>(asOverviewDims.size()), nullptr,
1173 33 : asOverviewDims.data(), pszResampling, aosOverviewOptions.List(),
1174 : GDALScaledProgress, pScaledProgress);
1175 33 : CPLDebug("COG", "Generating overviews of the imagery: end");
1176 :
1177 33 : GDALDestroyScaledProgress(pScaledProgress);
1178 33 : if (eErr != CE_None)
1179 : {
1180 1 : return nullptr;
1181 : }
1182 : }
1183 :
1184 266 : CPLStringList aosOptions;
1185 133 : aosOptions.SetNameValue("COPY_SRC_OVERVIEWS", "YES");
1186 133 : aosOptions.SetNameValue("COMPRESS", osCompress);
1187 133 : aosOptions.SetNameValue("TILED", "YES");
1188 133 : aosOptions.SetNameValue("BLOCKXSIZE", osBlockSize);
1189 133 : aosOptions.SetNameValue("BLOCKYSIZE", osBlockSize);
1190 : const char *pszPredictor =
1191 133 : CSLFetchNameValueDef(papszOptions, "PREDICTOR", "FALSE");
1192 133 : const char *pszPredictorValue = GetPredictor(poSrcDS, pszPredictor);
1193 133 : if (pszPredictorValue != nullptr)
1194 : {
1195 2 : aosOptions.SetNameValue("PREDICTOR", pszPredictorValue);
1196 : }
1197 :
1198 133 : const char *pszQuality = CSLFetchNameValue(papszOptions, "QUALITY");
1199 133 : if (EQUAL(osCompress, "JPEG"))
1200 : {
1201 11 : aosOptions.SetNameValue("JPEG_QUALITY", pszQuality);
1202 11 : if (nBands == 3)
1203 10 : aosOptions.SetNameValue("PHOTOMETRIC", "YCBCR");
1204 : }
1205 122 : else if (EQUAL(osCompress, "WEBP"))
1206 : {
1207 3 : if (pszQuality && atoi(pszQuality) == 100)
1208 2 : aosOptions.SetNameValue("WEBP_LOSSLESS", "YES");
1209 3 : aosOptions.SetNameValue("WEBP_LEVEL", pszQuality);
1210 : }
1211 119 : else if (EQUAL(osCompress, "DEFLATE") || EQUAL(osCompress, "LERC_DEFLATE"))
1212 : {
1213 : aosOptions.SetNameValue("ZLEVEL",
1214 8 : CSLFetchNameValue(papszOptions, "LEVEL"));
1215 : }
1216 111 : else if (EQUAL(osCompress, "ZSTD") || EQUAL(osCompress, "LERC_ZSTD"))
1217 : {
1218 : aosOptions.SetNameValue("ZSTD_LEVEL",
1219 3 : CSLFetchNameValue(papszOptions, "LEVEL"));
1220 : }
1221 108 : else if (EQUAL(osCompress, "LZMA"))
1222 : {
1223 : aosOptions.SetNameValue("LZMA_PRESET",
1224 1 : CSLFetchNameValue(papszOptions, "LEVEL"));
1225 : }
1226 :
1227 133 : if (STARTS_WITH_CI(osCompress, "LERC"))
1228 : {
1229 : aosOptions.SetNameValue("MAX_Z_ERROR",
1230 8 : CSLFetchNameValue(papszOptions, "MAX_Z_ERROR"));
1231 : aosOptions.SetNameValue(
1232 : "MAX_Z_ERROR_OVERVIEW",
1233 8 : CSLFetchNameValue(papszOptions, "MAX_Z_ERROR_OVERVIEW"));
1234 : }
1235 :
1236 133 : if (STARTS_WITH_CI(osCompress, "JXL"))
1237 : {
1238 8 : for (const char *pszKey : {"JXL_LOSSLESS", "JXL_EFFORT", "JXL_DISTANCE",
1239 10 : "JXL_ALPHA_DISTANCE"})
1240 : {
1241 8 : const char *pszValue = CSLFetchNameValue(papszOptions, pszKey);
1242 8 : if (pszValue)
1243 3 : aosOptions.SetNameValue(pszKey, pszValue);
1244 : }
1245 : }
1246 :
1247 : aosOptions.SetNameValue("BIGTIFF",
1248 133 : CSLFetchNameValue(papszOptions, "BIGTIFF"));
1249 : aosOptions.SetNameValue("NUM_THREADS",
1250 133 : CSLFetchNameValue(papszOptions, "NUM_THREADS"));
1251 : aosOptions.SetNameValue("GEOTIFF_VERSION",
1252 133 : CSLFetchNameValue(papszOptions, "GEOTIFF_VERSION"));
1253 : aosOptions.SetNameValue("SPARSE_OK",
1254 133 : CSLFetchNameValue(papszOptions, "SPARSE_OK"));
1255 133 : aosOptions.SetNameValue("NBITS", CSLFetchNameValue(papszOptions, "NBITS"));
1256 :
1257 133 : if (EQUAL(osOverviews, "NONE"))
1258 : {
1259 2 : aosOptions.SetNameValue("@OVERVIEW_DATASET", "");
1260 : }
1261 : else
1262 : {
1263 131 : if (!m_osTmpOverviewFilename.empty())
1264 : {
1265 : aosOptions.SetNameValue("@OVERVIEW_DATASET",
1266 32 : m_osTmpOverviewFilename);
1267 : }
1268 131 : if (!m_osTmpMskOverviewFilename.empty())
1269 : {
1270 : aosOptions.SetNameValue("@MASK_OVERVIEW_DATASET",
1271 8 : m_osTmpMskOverviewFilename);
1272 : }
1273 : aosOptions.SetNameValue(
1274 : "@OVERVIEW_COUNT",
1275 131 : CSLFetchNameValue(papszOptions, "OVERVIEW_COUNT"));
1276 : }
1277 :
1278 : const CPLString osTilingScheme(
1279 266 : CSLFetchNameValueDef(papszOptions, "TILING_SCHEME", "CUSTOM"));
1280 133 : if (osTilingScheme != "CUSTOM")
1281 : {
1282 23 : aosOptions.SetNameValue("@TILING_SCHEME_NAME", osTilingScheme);
1283 : aosOptions.SetNameValue("@TILING_SCHEME_ZOOM_LEVEL",
1284 23 : CPLSPrintf("%d", nZoomLevel));
1285 23 : if (nAlignedLevels > 0)
1286 : {
1287 : aosOptions.SetNameValue("@TILING_SCHEME_ALIGNED_LEVELS",
1288 4 : CPLSPrintf("%d", nAlignedLevels));
1289 : }
1290 : }
1291 133 : const char *pszOverviewCompress = CSLFetchNameValueDef(
1292 : papszOptions, "OVERVIEW_COMPRESS", osCompress.c_str());
1293 :
1294 : CPLConfigOptionSetter ovrCompressSetter("COMPRESS_OVERVIEW",
1295 266 : pszOverviewCompress, true);
1296 : const char *pszOverviewQuality =
1297 133 : CSLFetchNameValue(papszOptions, "OVERVIEW_QUALITY");
1298 : CPLConfigOptionSetter ovrQualityJpegSetter("JPEG_QUALITY_OVERVIEW",
1299 266 : pszOverviewQuality, true);
1300 :
1301 133 : std::unique_ptr<CPLConfigOptionSetter> poWebpLosslessSetter;
1302 133 : std::unique_ptr<CPLConfigOptionSetter> poWebpLevelSetter;
1303 133 : if (EQUAL(pszOverviewCompress, "WEBP"))
1304 : {
1305 3 : if (pszOverviewQuality && CPLAtof(pszOverviewQuality) == 100.0)
1306 : {
1307 0 : poWebpLosslessSetter.reset(new CPLConfigOptionSetter(
1308 0 : "WEBP_LOSSLESS_OVERVIEW", "TRUE", true));
1309 : }
1310 : else
1311 : {
1312 3 : poWebpLosslessSetter.reset(new CPLConfigOptionSetter(
1313 3 : "WEBP_LOSSLESS_OVERVIEW", "FALSE", true));
1314 3 : poWebpLevelSetter.reset(new CPLConfigOptionSetter(
1315 3 : "WEBP_LEVEL_OVERVIEW", pszOverviewQuality, true));
1316 : }
1317 : }
1318 :
1319 133 : std::unique_ptr<CPLConfigOptionSetter> poPhotometricSetter;
1320 133 : if (nBands == 3 && EQUAL(pszOverviewCompress, "JPEG"))
1321 : {
1322 10 : poPhotometricSetter.reset(
1323 10 : new CPLConfigOptionSetter("PHOTOMETRIC_OVERVIEW", "YCBCR", true));
1324 : }
1325 :
1326 : const char *osOvrPredictor =
1327 133 : CSLFetchNameValueDef(papszOptions, "OVERVIEW_PREDICTOR", "FALSE");
1328 133 : const char *pszOvrPredictorValue = GetPredictor(poSrcDS, osOvrPredictor);
1329 : CPLConfigOptionSetter ovrPredictorSetter("PREDICTOR_OVERVIEW",
1330 266 : pszOvrPredictorValue, true);
1331 :
1332 : GDALDriver *poGTiffDrv =
1333 133 : GDALDriver::FromHandle(GDALGetDriverByName("GTiff"));
1334 133 : if (!poGTiffDrv)
1335 0 : return nullptr;
1336 133 : void *pScaledProgress = GDALCreateScaledProgress(
1337 : dfCurPixels / dfTotalPixelsToProcess, 1.0, pfnProgress, pProgressData);
1338 :
1339 : CPLConfigOptionSetter oSetterInternalMask("GDAL_TIFF_INTERNAL_MASK", "YES",
1340 133 : false);
1341 :
1342 133 : const char *pszCopySrcMDD = CSLFetchNameValue(papszOptions, "COPY_SRC_MDD");
1343 133 : if (pszCopySrcMDD)
1344 2 : aosOptions.SetNameValue("COPY_SRC_MDD", pszCopySrcMDD);
1345 133 : char **papszSrcMDD = CSLFetchNameValueMultiple(papszOptions, "SRC_MDD");
1346 136 : for (CSLConstList papszSrcMDDIter = papszSrcMDD;
1347 136 : papszSrcMDDIter && *papszSrcMDDIter; ++papszSrcMDDIter)
1348 3 : aosOptions.AddNameValue("SRC_MDD", *papszSrcMDDIter);
1349 133 : CSLDestroy(papszSrcMDD);
1350 :
1351 133 : CPLDebug("COG", "Generating final product: start");
1352 : auto poRet =
1353 133 : poGTiffDrv->CreateCopy(pszFilename, poCurDS, false, aosOptions.List(),
1354 : GDALScaledProgress, pScaledProgress);
1355 :
1356 133 : GDALDestroyScaledProgress(pScaledProgress);
1357 :
1358 133 : if (poRet)
1359 121 : poRet->FlushCache(false);
1360 :
1361 133 : CPLDebug("COG", "Generating final product: end");
1362 133 : return poRet;
1363 : }
1364 :
1365 : /************************************************************************/
1366 : /* COGCreateCopy() */
1367 : /************************************************************************/
1368 :
1369 138 : static GDALDataset *COGCreateCopy(const char *pszFilename, GDALDataset *poSrcDS,
1370 : int /*bStrict*/, char **papszOptions,
1371 : GDALProgressFunc pfnProgress,
1372 : void *pProgressData)
1373 : {
1374 276 : return GDALCOGCreator().Create(pszFilename, poSrcDS, papszOptions,
1375 276 : pfnProgress, pProgressData);
1376 : }
1377 :
1378 : /************************************************************************/
1379 : /* GDALRegister_COG() */
1380 : /************************************************************************/
1381 :
1382 : class GDALCOGDriver final : public GDALDriver
1383 : {
1384 : bool m_bInitialized = false;
1385 :
1386 : bool bHasLZW = false;
1387 : bool bHasDEFLATE = false;
1388 : bool bHasLZMA = false;
1389 : bool bHasZSTD = false;
1390 : bool bHasJPEG = false;
1391 : bool bHasWebP = false;
1392 : bool bHasLERC = false;
1393 : std::string osCompressValues{};
1394 :
1395 : void InitializeCreationOptionList();
1396 :
1397 : public:
1398 : GDALCOGDriver();
1399 :
1400 40107 : const char *GetMetadataItem(const char *pszName,
1401 : const char *pszDomain) override
1402 : {
1403 40107 : if (EQUAL(pszName, GDAL_DMD_CREATIONOPTIONLIST))
1404 : {
1405 163 : InitializeCreationOptionList();
1406 : }
1407 40107 : return GDALDriver::GetMetadataItem(pszName, pszDomain);
1408 : }
1409 :
1410 59 : char **GetMetadata(const char *pszDomain) override
1411 : {
1412 59 : InitializeCreationOptionList();
1413 59 : return GDALDriver::GetMetadata(pszDomain);
1414 : }
1415 : };
1416 :
1417 1381 : GDALCOGDriver::GDALCOGDriver()
1418 : {
1419 : // We could defer this in InitializeCreationOptionList() but with currently
1420 : // released libtiff versions where there was a bug (now fixed) in
1421 : // TIFFGetConfiguredCODECs(), this wouldn't work properly if the LERC codec
1422 : // had been registered in between
1423 1381 : osCompressValues = GTiffGetCompressValues(bHasLZW, bHasDEFLATE, bHasLZMA,
1424 1381 : bHasZSTD, bHasJPEG, bHasWebP,
1425 1381 : bHasLERC, true /* bForCOG */);
1426 1381 : gbHasLZW = bHasLZW;
1427 1381 : }
1428 :
1429 222 : void GDALCOGDriver::InitializeCreationOptionList()
1430 : {
1431 222 : if (m_bInitialized)
1432 210 : return;
1433 12 : m_bInitialized = true;
1434 :
1435 24 : CPLString osOptions;
1436 : osOptions = "<CreationOptionList>"
1437 12 : " <Option name='COMPRESS' type='string-select' default='";
1438 12 : osOptions += bHasLZW ? "LZW" : "NONE";
1439 12 : osOptions += "'>";
1440 12 : osOptions += osCompressValues;
1441 12 : osOptions += " </Option>";
1442 :
1443 : osOptions +=
1444 12 : " <Option name='OVERVIEW_COMPRESS' type='string-select' default='";
1445 12 : osOptions += bHasLZW ? "LZW" : "NONE";
1446 12 : osOptions += "'>";
1447 12 : osOptions += osCompressValues;
1448 12 : osOptions += " </Option>";
1449 :
1450 12 : if (bHasLZW || bHasDEFLATE || bHasZSTD || bHasLZMA)
1451 : {
1452 12 : const char *osPredictorOptions =
1453 : " <Value>YES</Value>"
1454 : " <Value>NO</Value>"
1455 : " <Value alias='2'>STANDARD</Value>"
1456 : " <Value alias='3'>FLOATING_POINT</Value>";
1457 :
1458 : osOptions +=
1459 : " <Option name='LEVEL' type='int' "
1460 12 : "description='DEFLATE/ZSTD/LZMA compression level: 1 (fastest)'/>";
1461 :
1462 : osOptions +=
1463 12 : " <Option name='PREDICTOR' type='string-select' default='FALSE'>";
1464 12 : osOptions += osPredictorOptions;
1465 : osOptions += " </Option>"
1466 : " <Option name='OVERVIEW_PREDICTOR' "
1467 12 : "type='string-select' default='FALSE'>";
1468 12 : osOptions += osPredictorOptions;
1469 12 : osOptions += " </Option>";
1470 : }
1471 12 : if (bHasJPEG || bHasWebP)
1472 : {
1473 12 : std::string osJPEG_WEBP;
1474 12 : if (bHasJPEG)
1475 12 : osJPEG_WEBP = "JPEG";
1476 12 : if (bHasWebP)
1477 : {
1478 12 : if (!osJPEG_WEBP.empty())
1479 12 : osJPEG_WEBP += '/';
1480 12 : osJPEG_WEBP += "WEBP";
1481 : }
1482 : osOptions += " <Option name='QUALITY' type='int' "
1483 24 : "description='" +
1484 24 : osJPEG_WEBP +
1485 : " quality 1-100' min='1' max='100' default='75'/>"
1486 : " <Option name='OVERVIEW_QUALITY' type='int' "
1487 24 : "description='Overview " +
1488 24 : osJPEG_WEBP +
1489 : " quality 1-100' min='1' max='100' "
1490 12 : "default='75'/>";
1491 : }
1492 12 : if (bHasLERC)
1493 : {
1494 : osOptions +=
1495 : ""
1496 : " <Option name='MAX_Z_ERROR' type='float' description='Maximum "
1497 : "error for LERC compression' default='0'/>"
1498 : " <Option name='MAX_Z_ERROR_OVERVIEW' type='float' "
1499 : "description='Maximum error for LERC compression in overviews' "
1500 12 : "default='0'/>";
1501 : }
1502 : #ifdef HAVE_JXL
1503 : osOptions +=
1504 : ""
1505 : " <Option name='JXL_LOSSLESS' type='boolean' description='Whether "
1506 : "JPEGXL compression should be lossless' default='YES'/>"
1507 : " <Option name='JXL_EFFORT' type='int' description='Level of effort "
1508 : "1(fast)-9(slow)' min='1' max='9' default='5'/>"
1509 : " <Option name='JXL_DISTANCE' type='float' description='Distance "
1510 : "level for lossy compression (0=mathematically lossless, 1.0=visually "
1511 12 : "lossless, usual range [0.5,3])' default='1.0' min='0.01' max='25.0'/>";
1512 : #ifdef HAVE_JxlEncoderSetExtraChannelDistance
1513 : osOptions += " <Option name='JXL_ALPHA_DISTANCE' type='float' "
1514 : "description='Distance level for alpha channel "
1515 : "(-1=same as non-alpha channels, "
1516 : "0=mathematically lossless, 1.0=visually lossless, "
1517 12 : "usual range [0.5,3])' default='-1' min='-1' max='25.0'/>";
1518 : #endif
1519 : #endif
1520 : osOptions +=
1521 : " <Option name='NUM_THREADS' type='string' "
1522 : "description='Number of worker threads for compression. "
1523 : "Can be set to ALL_CPUS' default='1'/>"
1524 : " <Option name='NBITS' type='int' description='BITS for sub-byte "
1525 : "files (1-7), sub-uint16_t (9-15), sub-uint32_t (17-31), or float32 "
1526 : "(16)'/>"
1527 : " <Option name='BLOCKSIZE' type='int' "
1528 : "description='Tile size in pixels' min='128' default='512'/>"
1529 : " <Option name='BIGTIFF' type='string-select' description='"
1530 : "Force creation of BigTIFF file'>"
1531 : " <Value>YES</Value>"
1532 : " <Value>NO</Value>"
1533 : " <Value>IF_NEEDED</Value>"
1534 : " <Value>IF_SAFER</Value>"
1535 : " </Option>"
1536 : " <Option name='RESAMPLING' type='string' "
1537 : "description='Resampling method for overviews or warping'/>"
1538 : " <Option name='OVERVIEW_RESAMPLING' type='string' "
1539 : "description='Resampling method for overviews'/>"
1540 : " <Option name='WARP_RESAMPLING' type='string' "
1541 : "description='Resampling method for warping'/>"
1542 : " <Option name='OVERVIEWS' type='string-select' description='"
1543 : "Behavior regarding overviews'>"
1544 : " <Value>AUTO</Value>"
1545 : " <Value>IGNORE_EXISTING</Value>"
1546 : " <Value>FORCE_USE_EXISTING</Value>"
1547 : " <Value>NONE</Value>"
1548 : " </Option>"
1549 : " <Option name='OVERVIEW_COUNT' type='int' min='0' "
1550 : "description='Number of overviews'/>"
1551 : " <Option name='TILING_SCHEME' type='string-select' description='"
1552 : "Which tiling scheme to use pre-defined value or custom inline/outline "
1553 : "JSON definition' default='CUSTOM'>"
1554 12 : " <Value>CUSTOM</Value>";
1555 :
1556 24 : const auto tmsList = gdal::TileMatrixSet::listPredefinedTileMatrixSets();
1557 84 : for (const auto &tmsName : tmsList)
1558 : {
1559 144 : const auto poTM = gdal::TileMatrixSet::parse(tmsName.c_str());
1560 144 : if (poTM && poTM->haveAllLevelsSameTopLeft() &&
1561 216 : poTM->haveAllLevelsSameTileSize() &&
1562 72 : !poTM->hasVariableMatrixWidth())
1563 : {
1564 72 : osOptions += " <Value>";
1565 72 : osOptions += tmsName;
1566 72 : osOptions += "</Value>";
1567 : }
1568 : }
1569 :
1570 : osOptions +=
1571 : " </Option>"
1572 : " <Option name='ZOOM_LEVEL' type='int' description='Target zoom "
1573 : "level. "
1574 : "Only used for TILING_SCHEME != CUSTOM'/>"
1575 : " <Option name='ZOOM_LEVEL_STRATEGY' type='string-select' "
1576 : "description='Strategy to determine zoom level. "
1577 : "Only used for TILING_SCHEME != CUSTOM' default='AUTO'>"
1578 : " <Value>AUTO</Value>"
1579 : " <Value>LOWER</Value>"
1580 : " <Value>UPPER</Value>"
1581 : " </Option>"
1582 : " <Option name='TARGET_SRS' type='string' "
1583 : "description='Target SRS as EPSG:XXXX, WKT or PROJ string for "
1584 : "reprojection'/>"
1585 : " <Option name='RES' type='float' description='"
1586 : "Target resolution for reprojection'/>"
1587 : " <Option name='EXTENT' type='string' description='"
1588 : "Target extent as minx,miny,maxx,maxy for reprojection'/>"
1589 : " <Option name='ALIGNED_LEVELS' type='int' description='"
1590 : "Number of resolution levels for which the tiles from GeoTIFF and the "
1591 : "specified tiling scheme match'/>"
1592 : " <Option name='ADD_ALPHA' type='boolean' description='Can be set to "
1593 : "NO to "
1594 : "disable the addition of an alpha band in case of reprojection' "
1595 : "default='YES'/>"
1596 : #if LIBGEOTIFF_VERSION >= 1600
1597 : " <Option name='GEOTIFF_VERSION' type='string-select' default='AUTO' "
1598 : "description='Which version of GeoTIFF must be used'>"
1599 : " <Value>AUTO</Value>"
1600 : " <Value>1.0</Value>"
1601 : " <Value>1.1</Value>"
1602 : " </Option>"
1603 : #endif
1604 : " <Option name='SPARSE_OK' type='boolean' description='Should empty "
1605 : "blocks be omitted on disk?' default='FALSE'/>"
1606 : " <Option name='STATISTICS' type='string-select' default='AUTO' "
1607 : "description='Which to add statistics to the output file'>"
1608 : " <Value>AUTO</Value>"
1609 : " <Value>YES</Value>"
1610 : " <Value>NO</Value>"
1611 : " </Option>"
1612 12 : "</CreationOptionList>";
1613 :
1614 12 : SetMetadataItem(GDAL_DMD_CREATIONOPTIONLIST, osOptions.c_str());
1615 : }
1616 :
1617 1682 : void GDALRegister_COG()
1618 :
1619 : {
1620 1682 : if (GDALGetDriverByName("COG") != nullptr)
1621 301 : return;
1622 :
1623 1381 : auto poDriver = new GDALCOGDriver();
1624 1381 : poDriver->SetDescription("COG");
1625 1381 : poDriver->SetMetadataItem(GDAL_DCAP_RASTER, "YES");
1626 1381 : poDriver->SetMetadataItem(GDAL_DMD_LONGNAME,
1627 1381 : "Cloud optimized GeoTIFF generator");
1628 1381 : poDriver->SetMetadataItem(GDAL_DMD_HELPTOPIC, "drivers/raster/cog.html");
1629 1381 : poDriver->SetMetadataItem(GDAL_DMD_EXTENSIONS, "tif tiff");
1630 :
1631 1381 : poDriver->SetMetadataItem(
1632 : GDAL_DMD_CREATIONDATATYPES,
1633 : "Byte Int8 UInt16 Int16 UInt32 Int32 UInt64 Int64 Float32 "
1634 1381 : "Float64 CInt16 CInt32 CFloat32 CFloat64");
1635 :
1636 1381 : poDriver->SetMetadataItem(GDAL_DCAP_VIRTUALIO, "YES");
1637 :
1638 1381 : poDriver->SetMetadataItem(GDAL_DCAP_COORDINATE_EPOCH, "YES");
1639 :
1640 1381 : poDriver->pfnCreateCopy = COGCreateCopy;
1641 :
1642 1381 : GetGDALDriverManager()->RegisterDriver(poDriver);
1643 : }
|
