Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: High Performance Image Reprojector
4 : * Purpose: Implementation of high level convenience APIs for warper.
5 : * Author: Frank Warmerdam, warmerdam@pobox.com
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2003, Frank Warmerdam <warmerdam@pobox.com>
9 : * Copyright (c) 2008-2012, Even Rouault <even dot rouault at spatialys.com>
10 : *
11 : * SPDX-License-Identifier: MIT
12 : ****************************************************************************/
13 :
14 : #include "cpl_port.h"
15 : #include "gdalwarper.h"
16 :
17 : #include <stdlib.h>
18 : #include <string.h>
19 :
20 : #include <algorithm>
21 : #include <cmath>
22 : #include <limits>
23 :
24 : #include "cpl_conv.h"
25 : #include "cpl_error.h"
26 : #include "cpl_float.h"
27 : #include "cpl_mask.h"
28 : #include "cpl_minixml.h"
29 : #include "cpl_progress.h"
30 : #include "cpl_string.h"
31 : #include "cpl_vsi.h"
32 : #include "gdal.h"
33 : #include "gdal_priv.h"
34 : #include "ogr_api.h"
35 : #include "ogr_core.h"
36 : #include "vrtdataset.h" // for VRTSerializeNoData
37 :
38 : #if (defined(__x86_64) || defined(_M_X64))
39 : #include <emmintrin.h>
40 : #endif
41 :
42 : /************************************************************************/
43 : /* GDALReprojectImage() */
44 : /************************************************************************/
45 :
46 : /**
47 : * Reproject image.
48 : *
49 : * This is a convenience function utilizing the GDALWarpOperation class to
50 : * reproject an image from a source to a destination. In particular, this
51 : * function takes care of establishing the transformation function to
52 : * implement the reprojection, and will default a variety of other
53 : * warp options.
54 : *
55 : * Nodata values set on destination dataset are taken into account.
56 : *
57 : * No metadata, projection info, or color tables are transferred
58 : * to the output file. Source overviews are not considered.
59 : *
60 : * For more advanced warping capabilities, consider using GDALWarp().
61 : *
62 : * @param hSrcDS the source image file.
63 : * @param pszSrcWKT the source projection. If NULL the source projection
64 : * is read from from hSrcDS.
65 : * @param hDstDS the destination image file.
66 : * @param pszDstWKT the destination projection. If NULL the destination
67 : * projection will be read from hDstDS.
68 : * @param eResampleAlg the type of resampling to use.
69 : * @param dfWarpMemoryLimit the amount of memory (in bytes) that the warp
70 : * API is allowed to use for caching. This is in addition to the memory
71 : * already allocated to the GDAL caching (as per GDALSetCacheMax()). May be
72 : * 0.0 to use default memory settings.
73 : * @param dfMaxError maximum error measured in input pixels that is allowed
74 : * in approximating the transformation (0.0 for exact calculations).
75 : * @param pfnProgress a GDALProgressFunc() compatible callback function for
76 : * reporting progress or NULL.
77 : * @param pProgressArg argument to be passed to pfnProgress. May be NULL.
78 : * @param psOptions warp options, normally NULL.
79 : *
80 : * @return CE_None on success or CE_Failure if something goes wrong.
81 : * @see GDALWarp()
82 : */
83 :
84 52 : CPLErr CPL_STDCALL GDALReprojectImage(
85 : GDALDatasetH hSrcDS, const char *pszSrcWKT, GDALDatasetH hDstDS,
86 : const char *pszDstWKT, GDALResampleAlg eResampleAlg,
87 : CPL_UNUSED double dfWarpMemoryLimit, double dfMaxError,
88 : GDALProgressFunc pfnProgress, void *pProgressArg,
89 : GDALWarpOptions *psOptions)
90 :
91 : {
92 : /* -------------------------------------------------------------------- */
93 : /* Setup a reprojection based transformer. */
94 : /* -------------------------------------------------------------------- */
95 52 : void *hTransformArg = GDALCreateGenImgProjTransformer(
96 : hSrcDS, pszSrcWKT, hDstDS, pszDstWKT, TRUE, 1000.0, 0);
97 :
98 52 : if (hTransformArg == nullptr)
99 0 : return CE_Failure;
100 :
101 : /* -------------------------------------------------------------------- */
102 : /* Create a copy of the user provided options, or a defaulted */
103 : /* options structure. */
104 : /* -------------------------------------------------------------------- */
105 : GDALWarpOptions *psWOptions = psOptions == nullptr
106 52 : ? GDALCreateWarpOptions()
107 1 : : GDALCloneWarpOptions(psOptions);
108 :
109 52 : psWOptions->eResampleAlg = eResampleAlg;
110 :
111 : /* -------------------------------------------------------------------- */
112 : /* Set transform. */
113 : /* -------------------------------------------------------------------- */
114 52 : if (dfMaxError > 0.0)
115 : {
116 3 : psWOptions->pTransformerArg = GDALCreateApproxTransformer(
117 : GDALGenImgProjTransform, hTransformArg, dfMaxError);
118 :
119 3 : psWOptions->pfnTransformer = GDALApproxTransform;
120 : }
121 : else
122 : {
123 49 : psWOptions->pfnTransformer = GDALGenImgProjTransform;
124 49 : psWOptions->pTransformerArg = hTransformArg;
125 : }
126 :
127 : /* -------------------------------------------------------------------- */
128 : /* Set file and band mapping. */
129 : /* -------------------------------------------------------------------- */
130 52 : psWOptions->hSrcDS = hSrcDS;
131 52 : psWOptions->hDstDS = hDstDS;
132 :
133 52 : int nSrcBands = GDALGetRasterCount(hSrcDS);
134 : {
135 52 : GDALRasterBandH hBand = GDALGetRasterBand(hSrcDS, nSrcBands);
136 52 : if (hBand && GDALGetRasterColorInterpretation(hBand) == GCI_AlphaBand)
137 : {
138 2 : psWOptions->nSrcAlphaBand = nSrcBands;
139 2 : nSrcBands--;
140 : }
141 : }
142 :
143 52 : int nDstBands = GDALGetRasterCount(hDstDS);
144 : {
145 52 : GDALRasterBandH hBand = GDALGetRasterBand(hDstDS, nDstBands);
146 52 : if (hBand && GDALGetRasterColorInterpretation(hBand) == GCI_AlphaBand)
147 : {
148 2 : psWOptions->nDstAlphaBand = nDstBands;
149 2 : nDstBands--;
150 : }
151 : }
152 :
153 52 : GDALWarpInitDefaultBandMapping(psWOptions, std::min(nSrcBands, nDstBands));
154 :
155 : /* -------------------------------------------------------------------- */
156 : /* Set source nodata values if the source dataset seems to have */
157 : /* any. Same for target nodata values */
158 : /* -------------------------------------------------------------------- */
159 108 : for (int iBand = 0; iBand < psWOptions->nBandCount; iBand++)
160 : {
161 56 : GDALRasterBandH hBand = GDALGetRasterBand(hSrcDS, iBand + 1);
162 :
163 56 : int bGotNoData = FALSE;
164 56 : double dfNoDataValue = GDALGetRasterNoDataValue(hBand, &bGotNoData);
165 56 : if (bGotNoData)
166 : {
167 2 : GDALWarpInitSrcNoDataReal(psWOptions, -1.1e20);
168 2 : psWOptions->padfSrcNoDataReal[iBand] = dfNoDataValue;
169 : }
170 :
171 : // Deal with target band.
172 56 : hBand = GDALGetRasterBand(hDstDS, iBand + 1);
173 :
174 56 : dfNoDataValue = GDALGetRasterNoDataValue(hBand, &bGotNoData);
175 56 : if (bGotNoData)
176 : {
177 2 : GDALWarpInitDstNoDataReal(psWOptions, -1.1e20);
178 2 : psWOptions->padfDstNoDataReal[iBand] = dfNoDataValue;
179 : }
180 : }
181 :
182 : /* -------------------------------------------------------------------- */
183 : /* Set the progress function. */
184 : /* -------------------------------------------------------------------- */
185 52 : if (pfnProgress != nullptr)
186 : {
187 3 : psWOptions->pfnProgress = pfnProgress;
188 3 : psWOptions->pProgressArg = pProgressArg;
189 : }
190 :
191 : /* -------------------------------------------------------------------- */
192 : /* Create a warp options based on the options. */
193 : /* -------------------------------------------------------------------- */
194 52 : GDALWarpOperation oWarper;
195 52 : CPLErr eErr = oWarper.Initialize(psWOptions);
196 :
197 52 : if (eErr == CE_None)
198 52 : eErr = oWarper.ChunkAndWarpImage(0, 0, GDALGetRasterXSize(hDstDS),
199 : GDALGetRasterYSize(hDstDS));
200 :
201 : /* -------------------------------------------------------------------- */
202 : /* Cleanup. */
203 : /* -------------------------------------------------------------------- */
204 52 : GDALDestroyGenImgProjTransformer(hTransformArg);
205 :
206 52 : if (dfMaxError > 0.0)
207 3 : GDALDestroyApproxTransformer(psWOptions->pTransformerArg);
208 :
209 52 : GDALDestroyWarpOptions(psWOptions);
210 :
211 52 : return eErr;
212 : }
213 :
214 : /************************************************************************/
215 : /* GDALCreateAndReprojectImage() */
216 : /* */
217 : /* This is a "quickie" reprojection API. */
218 : /************************************************************************/
219 :
220 : /** Reproject an image and create the target reprojected image */
221 0 : CPLErr CPL_STDCALL GDALCreateAndReprojectImage(
222 : GDALDatasetH hSrcDS, const char *pszSrcWKT, const char *pszDstFilename,
223 : const char *pszDstWKT, GDALDriverH hDstDriver, char **papszCreateOptions,
224 : GDALResampleAlg eResampleAlg, double dfWarpMemoryLimit, double dfMaxError,
225 : GDALProgressFunc pfnProgress, void *pProgressArg,
226 : GDALWarpOptions *psOptions)
227 :
228 : {
229 0 : VALIDATE_POINTER1(hSrcDS, "GDALCreateAndReprojectImage", CE_Failure);
230 :
231 : /* -------------------------------------------------------------------- */
232 : /* Default a few parameters. */
233 : /* -------------------------------------------------------------------- */
234 0 : if (hDstDriver == nullptr)
235 : {
236 0 : hDstDriver = GDALGetDriverByName("GTiff");
237 0 : if (hDstDriver == nullptr)
238 : {
239 0 : CPLError(CE_Failure, CPLE_AppDefined,
240 : "GDALCreateAndReprojectImage needs GTiff driver");
241 0 : return CE_Failure;
242 : }
243 : }
244 :
245 0 : if (pszSrcWKT == nullptr)
246 0 : pszSrcWKT = GDALGetProjectionRef(hSrcDS);
247 :
248 0 : if (pszDstWKT == nullptr)
249 0 : pszDstWKT = pszSrcWKT;
250 :
251 : /* -------------------------------------------------------------------- */
252 : /* Create a transformation object from the source to */
253 : /* destination coordinate system. */
254 : /* -------------------------------------------------------------------- */
255 0 : void *hTransformArg = GDALCreateGenImgProjTransformer(
256 : hSrcDS, pszSrcWKT, nullptr, pszDstWKT, TRUE, 1000.0, 0);
257 :
258 0 : if (hTransformArg == nullptr)
259 0 : return CE_Failure;
260 :
261 : /* -------------------------------------------------------------------- */
262 : /* Get approximate output definition. */
263 : /* -------------------------------------------------------------------- */
264 0 : double adfDstGeoTransform[6] = {};
265 0 : int nPixels = 0;
266 0 : int nLines = 0;
267 :
268 : CPLErr eErr =
269 0 : GDALSuggestedWarpOutput(hSrcDS, GDALGenImgProjTransform, hTransformArg,
270 : adfDstGeoTransform, &nPixels, &nLines);
271 :
272 0 : GDALDestroyGenImgProjTransformer(hTransformArg);
273 :
274 0 : if (eErr != CE_None)
275 0 : return eErr;
276 :
277 : /* -------------------------------------------------------------------- */
278 : /* Create the output file. */
279 : /* -------------------------------------------------------------------- */
280 0 : GDALDatasetH hDstDS = GDALCreate(
281 : hDstDriver, pszDstFilename, nPixels, nLines, GDALGetRasterCount(hSrcDS),
282 : GDALGetRasterDataType(GDALGetRasterBand(hSrcDS, 1)),
283 : papszCreateOptions);
284 :
285 0 : if (hDstDS == nullptr)
286 0 : return CE_Failure;
287 :
288 : /* -------------------------------------------------------------------- */
289 : /* Write out the projection definition. */
290 : /* -------------------------------------------------------------------- */
291 0 : GDALSetProjection(hDstDS, pszDstWKT);
292 0 : GDALSetGeoTransform(hDstDS, adfDstGeoTransform);
293 :
294 : /* -------------------------------------------------------------------- */
295 : /* Perform the reprojection. */
296 : /* -------------------------------------------------------------------- */
297 0 : eErr = GDALReprojectImage(hSrcDS, pszSrcWKT, hDstDS, pszDstWKT,
298 : eResampleAlg, dfWarpMemoryLimit, dfMaxError,
299 : pfnProgress, pProgressArg, psOptions);
300 :
301 0 : GDALClose(hDstDS);
302 :
303 0 : return eErr;
304 : }
305 :
306 : /************************************************************************/
307 : /* GDALWarpNoDataMaskerT() */
308 : /************************************************************************/
309 :
310 : template <class T>
311 760 : static CPLErr GDALWarpNoDataMaskerT(const double *padfNoData, size_t nPixels,
312 : const T *pData, GUInt32 *panValidityMask,
313 : int *pbOutAllValid)
314 : {
315 : // Nothing to do if value is out of range.
316 760 : if (padfNoData[0] < cpl::NumericLimits<T>::min() ||
317 1473 : padfNoData[0] > cpl::NumericLimits<T>::max() + 0.000001 ||
318 713 : padfNoData[1] != 0.0)
319 : {
320 47 : *pbOutAllValid = TRUE;
321 47 : return CE_None;
322 : }
323 :
324 713 : const int nNoData = static_cast<int>(floor(padfNoData[0] + 0.000001));
325 713 : int bAllValid = TRUE;
326 115341239 : for (size_t iOffset = 0; iOffset < nPixels; ++iOffset)
327 : {
328 115341089 : if (pData[iOffset] == nNoData)
329 : {
330 110875912 : bAllValid = FALSE;
331 110875912 : CPLMaskClear(panValidityMask, iOffset);
332 : }
333 : }
334 713 : *pbOutAllValid = bAllValid;
335 :
336 713 : return CE_None;
337 : }
338 :
339 : /************************************************************************/
340 : /* GDALWarpNoDataMasker() */
341 : /* */
342 : /* GDALMaskFunc for establishing a validity mask for a source */
343 : /* band based on a provided NODATA value. */
344 : /************************************************************************/
345 :
346 1025 : CPLErr GDALWarpNoDataMasker(void *pMaskFuncArg, int nBandCount,
347 : GDALDataType eType, int /* nXOff */,
348 : int /* nYOff */, int nXSize, int nYSize,
349 : GByte **ppImageData, int bMaskIsFloat,
350 : void *pValidityMask, int *pbOutAllValid)
351 :
352 : {
353 1025 : const double *padfNoData = static_cast<double *>(pMaskFuncArg);
354 1025 : GUInt32 *panValidityMask = static_cast<GUInt32 *>(pValidityMask);
355 1025 : const size_t nPixels = static_cast<size_t>(nXSize) * nYSize;
356 :
357 1025 : *pbOutAllValid = FALSE;
358 :
359 1025 : if (nBandCount != 1 || bMaskIsFloat)
360 : {
361 0 : CPLError(
362 : CE_Failure, CPLE_AppDefined,
363 : "Invalid nBandCount or bMaskIsFloat argument in SourceNoDataMask");
364 0 : return CE_Failure;
365 : }
366 :
367 1025 : CPLErr eErr = CE_None;
368 :
369 1025 : switch (eType)
370 : {
371 608 : case GDT_UInt8:
372 608 : return GDALWarpNoDataMaskerT(padfNoData, nPixels,
373 : *ppImageData, // Already a GByte *.
374 608 : panValidityMask, pbOutAllValid);
375 :
376 28 : case GDT_Int8:
377 28 : return GDALWarpNoDataMaskerT(
378 : padfNoData, nPixels, reinterpret_cast<int8_t *>(*ppImageData),
379 28 : panValidityMask, pbOutAllValid);
380 :
381 97 : case GDT_Int16:
382 97 : return GDALWarpNoDataMaskerT(
383 : padfNoData, nPixels, reinterpret_cast<GInt16 *>(*ppImageData),
384 97 : panValidityMask, pbOutAllValid);
385 :
386 27 : case GDT_UInt16:
387 27 : return GDALWarpNoDataMaskerT(
388 : padfNoData, nPixels, reinterpret_cast<GUInt16 *>(*ppImageData),
389 27 : panValidityMask, pbOutAllValid);
390 :
391 116 : case GDT_Float32:
392 : {
393 116 : const float fNoData = static_cast<float>(padfNoData[0]);
394 116 : const float *pafData = reinterpret_cast<float *>(*ppImageData);
395 116 : const bool bIsNoDataNan = std::isnan(fNoData);
396 :
397 : // Nothing to do if value is out of range.
398 116 : if (padfNoData[1] != 0.0)
399 : {
400 0 : *pbOutAllValid = TRUE;
401 0 : return CE_None;
402 : }
403 :
404 116 : int bAllValid = TRUE;
405 219554 : for (size_t iOffset = 0; iOffset < nPixels; ++iOffset)
406 : {
407 219438 : float fVal = pafData[iOffset];
408 438628 : if ((bIsNoDataNan && std::isnan(fVal)) ||
409 219190 : (!bIsNoDataNan && ARE_REAL_EQUAL(fVal, fNoData)))
410 : {
411 123734 : bAllValid = FALSE;
412 123734 : CPLMaskClear(panValidityMask, iOffset);
413 : }
414 : }
415 116 : *pbOutAllValid = bAllValid;
416 : }
417 116 : break;
418 :
419 41 : case GDT_Float64:
420 : {
421 41 : const double dfNoData = padfNoData[0];
422 41 : const double *padfData = reinterpret_cast<double *>(*ppImageData);
423 41 : const bool bIsNoDataNan = std::isnan(dfNoData);
424 :
425 : // Nothing to do if value is out of range.
426 41 : if (padfNoData[1] != 0.0)
427 : {
428 0 : *pbOutAllValid = TRUE;
429 0 : return CE_None;
430 : }
431 :
432 41 : int bAllValid = TRUE;
433 5912470 : for (size_t iOffset = 0; iOffset < nPixels; ++iOffset)
434 : {
435 5912430 : double dfVal = padfData[iOffset];
436 5913270 : if ((bIsNoDataNan && std::isnan(dfVal)) ||
437 840 : (!bIsNoDataNan && ARE_REAL_EQUAL(dfVal, dfNoData)))
438 : {
439 5912040 : bAllValid = FALSE;
440 5912040 : CPLMaskClear(panValidityMask, iOffset);
441 : }
442 : }
443 41 : *pbOutAllValid = bAllValid;
444 : }
445 41 : break;
446 :
447 108 : default:
448 : {
449 108 : const int nWordSize = GDALGetDataTypeSizeBytes(eType);
450 :
451 108 : const bool bIsNoDataRealNan = std::isnan(padfNoData[0]);
452 :
453 108 : eErr = CE_Failure;
454 : double *padfWrk = static_cast<double *>(
455 108 : VSI_MALLOC2_VERBOSE(nXSize, sizeof(double) * 2));
456 108 : if (padfWrk)
457 : {
458 108 : eErr = CE_None;
459 108 : bool bAllValid = true;
460 714 : for (size_t iLine = 0; iLine < static_cast<size_t>(nYSize);
461 : iLine++)
462 : {
463 606 : GDALCopyWords64((*ppImageData) + nWordSize * iLine * nXSize,
464 : eType, nWordSize, padfWrk, GDT_CFloat64, 16,
465 : nXSize);
466 :
467 606 : const size_t iOffsetLine = iLine * nXSize;
468 8122 : for (size_t iPixel = 0;
469 8122 : iPixel < static_cast<size_t>(nXSize); ++iPixel)
470 : {
471 0 : if (((bIsNoDataRealNan &&
472 15032 : std::isnan(padfWrk[iPixel * 2])) ||
473 15032 : (!bIsNoDataRealNan &&
474 7516 : ARE_REAL_EQUAL(padfWrk[iPixel * 2],
475 : padfNoData[0]))))
476 : {
477 2213 : const size_t iOffset = iOffsetLine + iPixel;
478 :
479 2213 : bAllValid = false;
480 2213 : CPLMaskClear(panValidityMask, iOffset);
481 : }
482 : }
483 : }
484 108 : *pbOutAllValid = bAllValid;
485 :
486 108 : VSIFree(padfWrk);
487 : }
488 : }
489 108 : break;
490 : }
491 :
492 265 : return eErr;
493 : }
494 :
495 : /************************************************************************/
496 : /* GDALWarpSrcAlphaMasker() */
497 : /* */
498 : /* GDALMaskFunc for reading source simple 8bit alpha mask */
499 : /* information and building a floating point density mask from */
500 : /* it. */
501 : /************************************************************************/
502 :
503 142 : CPLErr GDALWarpSrcAlphaMasker(void *pMaskFuncArg, int /* nBandCount */,
504 : GDALDataType /* eType */, int nXOff, int nYOff,
505 : int nXSize, int nYSize, GByte ** /*ppImageData */,
506 : int bMaskIsFloat, void *pValidityMask,
507 : int *pbOutAllOpaque)
508 :
509 : {
510 142 : GDALWarpOptions *psWO = static_cast<GDALWarpOptions *>(pMaskFuncArg);
511 142 : float *pafMask = static_cast<float *>(pValidityMask);
512 142 : *pbOutAllOpaque = FALSE;
513 142 : const size_t nPixels = static_cast<size_t>(nXSize) * nYSize;
514 :
515 : /* -------------------------------------------------------------------- */
516 : /* Do some minimal checking. */
517 : /* -------------------------------------------------------------------- */
518 142 : if (!bMaskIsFloat)
519 : {
520 0 : CPLAssert(false);
521 : return CE_Failure;
522 : }
523 :
524 142 : if (psWO == nullptr || psWO->nSrcAlphaBand < 1)
525 : {
526 0 : CPLAssert(false);
527 : return CE_Failure;
528 : }
529 :
530 : /* -------------------------------------------------------------------- */
531 : /* Read the alpha band. */
532 : /* -------------------------------------------------------------------- */
533 : GDALRasterBandH hAlphaBand =
534 142 : GDALGetRasterBand(psWO->hSrcDS, psWO->nSrcAlphaBand);
535 142 : if (hAlphaBand == nullptr)
536 0 : return CE_Failure;
537 :
538 : // Rescale.
539 142 : const float inv_alpha_max = static_cast<float>(
540 142 : 1.0 / CPLAtof(CSLFetchNameValueDef(psWO->papszWarpOptions,
541 142 : "SRC_ALPHA_MAX", "255")));
542 142 : bool bOutAllOpaque = true;
543 :
544 142 : size_t iPixel = 0;
545 : CPLErr eErr;
546 :
547 : #if (defined(__x86_64) || defined(_M_X64))
548 142 : GDALDataType eDT = GDALGetRasterDataType(hAlphaBand);
549 : // Make sure that pafMask is at least 8-byte aligned, which should
550 : // normally be always the case if being a ptr returned by malloc().
551 142 : if ((eDT == GDT_UInt8 || eDT == GDT_UInt16) && CPL_IS_ALIGNED(pafMask, 8))
552 : {
553 : // Read data.
554 262 : eErr = GDALRasterIOEx(
555 : hAlphaBand, GF_Read, nXOff, nYOff, nXSize, nYSize, pafMask, nXSize,
556 : nYSize, eDT, static_cast<GSpacing>(sizeof(int)),
557 131 : static_cast<GSpacing>(sizeof(int)) * nXSize, nullptr);
558 :
559 131 : if (eErr != CE_None)
560 0 : return eErr;
561 :
562 : // Make sure we have the correct alignment before doing SSE
563 : // On Linux x86_64, the alignment should be always correct due
564 : // the alignment of malloc() being 16 byte.
565 131 : const GUInt32 mask = (eDT == GDT_UInt8) ? 0xff : 0xffff;
566 131 : if (!CPL_IS_ALIGNED(pafMask, 16))
567 : {
568 0 : pafMask[iPixel] =
569 0 : (reinterpret_cast<GUInt32 *>(pafMask)[iPixel] & mask) *
570 0 : inv_alpha_max;
571 0 : if (pafMask[iPixel] >= 1.0f)
572 0 : pafMask[iPixel] = 1.0f;
573 : else
574 0 : bOutAllOpaque = false;
575 0 : iPixel++;
576 : }
577 131 : CPLAssert(CPL_IS_ALIGNED(pafMask + iPixel, 16));
578 131 : const __m128 xmm_inverse_alpha_max = _mm_load1_ps(&inv_alpha_max);
579 131 : const float one_single = 1.0f;
580 131 : const __m128 xmm_one = _mm_load1_ps(&one_single);
581 262 : const __m128i xmm_i_mask = _mm_set1_epi32(mask);
582 131 : __m128 xmmMaskNonOpaque0 = _mm_setzero_ps();
583 131 : __m128 xmmMaskNonOpaque1 = _mm_setzero_ps();
584 131 : __m128 xmmMaskNonOpaque2 = _mm_setzero_ps();
585 231086 : for (; iPixel + 6 * 4 - 1 < nPixels; iPixel += 6 * 4)
586 : {
587 461910 : __m128 xmm_mask0 = _mm_cvtepi32_ps(_mm_and_si128(
588 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
589 230955 : pafMask + iPixel + 4 * 0))));
590 461910 : __m128 xmm_mask1 = _mm_cvtepi32_ps(_mm_and_si128(
591 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
592 230955 : pafMask + iPixel + 4 * 1))));
593 461910 : __m128 xmm_mask2 = _mm_cvtepi32_ps(_mm_and_si128(
594 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
595 230955 : pafMask + iPixel + 4 * 2))));
596 461910 : __m128 xmm_mask3 = _mm_cvtepi32_ps(_mm_and_si128(
597 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
598 230955 : pafMask + iPixel + 4 * 3))));
599 461910 : __m128 xmm_mask4 = _mm_cvtepi32_ps(_mm_and_si128(
600 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
601 230955 : pafMask + iPixel + 4 * 4))));
602 692865 : __m128 xmm_mask5 = _mm_cvtepi32_ps(_mm_and_si128(
603 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
604 230955 : pafMask + iPixel + 4 * 5))));
605 230955 : xmm_mask0 = _mm_mul_ps(xmm_mask0, xmm_inverse_alpha_max);
606 230955 : xmm_mask1 = _mm_mul_ps(xmm_mask1, xmm_inverse_alpha_max);
607 230955 : xmm_mask2 = _mm_mul_ps(xmm_mask2, xmm_inverse_alpha_max);
608 230955 : xmm_mask3 = _mm_mul_ps(xmm_mask3, xmm_inverse_alpha_max);
609 230955 : xmm_mask4 = _mm_mul_ps(xmm_mask4, xmm_inverse_alpha_max);
610 230955 : xmm_mask5 = _mm_mul_ps(xmm_mask5, xmm_inverse_alpha_max);
611 : xmmMaskNonOpaque0 =
612 461910 : _mm_or_ps(xmmMaskNonOpaque0, _mm_cmplt_ps(xmm_mask0, xmm_one));
613 : xmmMaskNonOpaque1 =
614 461910 : _mm_or_ps(xmmMaskNonOpaque1, _mm_cmplt_ps(xmm_mask1, xmm_one));
615 : xmmMaskNonOpaque2 =
616 461910 : _mm_or_ps(xmmMaskNonOpaque2, _mm_cmplt_ps(xmm_mask2, xmm_one));
617 : xmmMaskNonOpaque0 =
618 461910 : _mm_or_ps(xmmMaskNonOpaque0, _mm_cmplt_ps(xmm_mask3, xmm_one));
619 : xmmMaskNonOpaque1 =
620 461910 : _mm_or_ps(xmmMaskNonOpaque1, _mm_cmplt_ps(xmm_mask4, xmm_one));
621 : xmmMaskNonOpaque2 =
622 461910 : _mm_or_ps(xmmMaskNonOpaque2, _mm_cmplt_ps(xmm_mask5, xmm_one));
623 230955 : xmm_mask0 = _mm_min_ps(xmm_mask0, xmm_one);
624 230955 : xmm_mask1 = _mm_min_ps(xmm_mask1, xmm_one);
625 230955 : xmm_mask2 = _mm_min_ps(xmm_mask2, xmm_one);
626 230955 : xmm_mask3 = _mm_min_ps(xmm_mask3, xmm_one);
627 230955 : xmm_mask4 = _mm_min_ps(xmm_mask4, xmm_one);
628 230955 : xmm_mask5 = _mm_min_ps(xmm_mask5, xmm_one);
629 230955 : _mm_store_ps(pafMask + iPixel + 4 * 0, xmm_mask0);
630 230955 : _mm_store_ps(pafMask + iPixel + 4 * 1, xmm_mask1);
631 230955 : _mm_store_ps(pafMask + iPixel + 4 * 2, xmm_mask2);
632 230955 : _mm_store_ps(pafMask + iPixel + 4 * 3, xmm_mask3);
633 230955 : _mm_store_ps(pafMask + iPixel + 4 * 4, xmm_mask4);
634 230955 : _mm_store_ps(pafMask + iPixel + 4 * 5, xmm_mask5);
635 : }
636 262 : if (_mm_movemask_ps(
637 : _mm_or_ps(_mm_or_ps(xmmMaskNonOpaque0, xmmMaskNonOpaque1),
638 131 : xmmMaskNonOpaque2)))
639 : {
640 93 : bOutAllOpaque = false;
641 : }
642 1519 : for (; iPixel < nPixels; iPixel++)
643 : {
644 1388 : pafMask[iPixel] =
645 1388 : (reinterpret_cast<GUInt32 *>(pafMask)[iPixel] & mask) *
646 1388 : inv_alpha_max;
647 1388 : if (pafMask[iPixel] >= 1.0f)
648 406 : pafMask[iPixel] = 1.0f;
649 : else
650 982 : bOutAllOpaque = false;
651 131 : }
652 : }
653 : else
654 : #endif
655 : {
656 : // Read data.
657 11 : eErr = GDALRasterIO(hAlphaBand, GF_Read, nXOff, nYOff, nXSize, nYSize,
658 : pafMask, nXSize, nYSize, GDT_Float32, 0, 0);
659 :
660 11 : if (eErr != CE_None)
661 0 : return eErr;
662 :
663 : // TODO(rouault): Is loop unrolling by hand (r34564) actually helpful?
664 12969 : for (; iPixel + 3 < nPixels; iPixel += 4)
665 : {
666 12958 : pafMask[iPixel] = pafMask[iPixel] * inv_alpha_max;
667 12958 : if (pafMask[iPixel] >= 1.0f)
668 3028 : pafMask[iPixel] = 1.0f;
669 : else
670 9930 : bOutAllOpaque = false;
671 12958 : pafMask[iPixel + 1] = pafMask[iPixel + 1] * inv_alpha_max;
672 12958 : if (pafMask[iPixel + 1] >= 1.0f)
673 3018 : pafMask[iPixel + 1] = 1.0f;
674 : else
675 9940 : bOutAllOpaque = false;
676 12958 : pafMask[iPixel + 2] = pafMask[iPixel + 2] * inv_alpha_max;
677 12958 : if (pafMask[iPixel + 2] >= 1.0f)
678 3066 : pafMask[iPixel + 2] = 1.0f;
679 : else
680 9892 : bOutAllOpaque = false;
681 12958 : pafMask[iPixel + 3] = pafMask[iPixel + 3] * inv_alpha_max;
682 12958 : if (pafMask[iPixel + 3] >= 1.0f)
683 3016 : pafMask[iPixel + 3] = 1.0f;
684 : else
685 9942 : bOutAllOpaque = false;
686 : }
687 :
688 12 : for (; iPixel < nPixels; iPixel++)
689 : {
690 1 : pafMask[iPixel] = pafMask[iPixel] * inv_alpha_max;
691 1 : if (pafMask[iPixel] >= 1.0f)
692 0 : pafMask[iPixel] = 1.0f;
693 : else
694 1 : bOutAllOpaque = false;
695 : }
696 : }
697 :
698 142 : *pbOutAllOpaque = bOutAllOpaque;
699 :
700 142 : return CE_None;
701 : }
702 :
703 : /************************************************************************/
704 : /* GDALWarpSrcMaskMasker() */
705 : /* */
706 : /* GDALMaskFunc for reading source simple 8bit validity mask */
707 : /* information and building a one bit validity mask. */
708 : /************************************************************************/
709 :
710 2 : CPLErr GDALWarpSrcMaskMasker(void *pMaskFuncArg, int /* nBandCount */,
711 : GDALDataType /* eType */, int nXOff, int nYOff,
712 : int nXSize, int nYSize, GByte ** /*ppImageData */,
713 : int bMaskIsFloat, void *pValidityMask)
714 :
715 : {
716 2 : GDALWarpOptions *psWO = static_cast<GDALWarpOptions *>(pMaskFuncArg);
717 2 : GUInt32 *panMask = static_cast<GUInt32 *>(pValidityMask);
718 :
719 : /* -------------------------------------------------------------------- */
720 : /* Do some minimal checking. */
721 : /* -------------------------------------------------------------------- */
722 2 : if (bMaskIsFloat)
723 : {
724 0 : CPLAssert(false);
725 : return CE_Failure;
726 : }
727 :
728 2 : if (psWO == nullptr)
729 : {
730 0 : CPLAssert(false);
731 : return CE_Failure;
732 : }
733 :
734 : /* -------------------------------------------------------------------- */
735 : /* Allocate a temporary buffer to read mask byte data into. */
736 : /* -------------------------------------------------------------------- */
737 : GByte *pabySrcMask =
738 2 : static_cast<GByte *>(VSI_MALLOC2_VERBOSE(nXSize, nYSize));
739 2 : if (pabySrcMask == nullptr)
740 : {
741 0 : return CE_Failure;
742 : }
743 :
744 : /* -------------------------------------------------------------------- */
745 : /* Fetch our mask band. */
746 : /* -------------------------------------------------------------------- */
747 2 : GDALRasterBandH hMaskBand = nullptr;
748 : GDALRasterBandH hSrcBand =
749 2 : GDALGetRasterBand(psWO->hSrcDS, psWO->panSrcBands[0]);
750 2 : if (hSrcBand != nullptr)
751 2 : hMaskBand = GDALGetMaskBand(hSrcBand);
752 :
753 2 : if (hMaskBand == nullptr)
754 : {
755 0 : CPLAssert(false);
756 : return CE_Failure;
757 : }
758 :
759 : /* -------------------------------------------------------------------- */
760 : /* Read the mask band. */
761 : /* -------------------------------------------------------------------- */
762 2 : CPLErr eErr = GDALRasterIO(hMaskBand, GF_Read, nXOff, nYOff, nXSize, nYSize,
763 : pabySrcMask, nXSize, nYSize, GDT_UInt8, 0, 0);
764 :
765 2 : if (eErr != CE_None)
766 : {
767 0 : CPLFree(pabySrcMask);
768 0 : return eErr;
769 : }
770 :
771 : /* -------------------------------------------------------------------- */
772 : /* Pack into 1 bit per pixel for validity. */
773 : /* -------------------------------------------------------------------- */
774 2 : const size_t nPixels = static_cast<size_t>(nXSize) * nYSize;
775 234579 : for (size_t iPixel = 0; iPixel < nPixels; iPixel++)
776 : {
777 234577 : if (pabySrcMask[iPixel] == 0)
778 31660 : CPLMaskClear(panMask, iPixel);
779 : }
780 :
781 2 : CPLFree(pabySrcMask);
782 :
783 2 : return CE_None;
784 : }
785 :
786 : /************************************************************************/
787 : /* GDALWarpDstAlphaMasker() */
788 : /* */
789 : /* GDALMaskFunc for reading or writing the destination simple */
790 : /* 8bit alpha mask information and building a floating point */
791 : /* density mask from it. Note, writing is distinguished */
792 : /* negative bandcount. */
793 : /************************************************************************/
794 :
795 3424 : CPLErr GDALWarpDstAlphaMasker(void *pMaskFuncArg, int nBandCount,
796 : CPL_UNUSED GDALDataType /* eType */, int nXOff,
797 : int nYOff, int nXSize, int nYSize,
798 : GByte ** /*ppImageData */, int bMaskIsFloat,
799 : void *pValidityMask)
800 : {
801 : /* -------------------------------------------------------------------- */
802 : /* Do some minimal checking. */
803 : /* -------------------------------------------------------------------- */
804 3424 : if (!bMaskIsFloat)
805 : {
806 0 : CPLAssert(false);
807 : return CE_Failure;
808 : }
809 :
810 3424 : GDALWarpOptions *psWO = static_cast<GDALWarpOptions *>(pMaskFuncArg);
811 3424 : if (psWO == nullptr || psWO->nDstAlphaBand < 1)
812 : {
813 0 : CPLAssert(false);
814 : return CE_Failure;
815 : }
816 :
817 3424 : float *pafMask = static_cast<float *>(pValidityMask);
818 3424 : const size_t nPixels = static_cast<size_t>(nXSize) * nYSize;
819 :
820 : GDALRasterBandH hAlphaBand =
821 3424 : GDALGetRasterBand(psWO->hDstDS, psWO->nDstAlphaBand);
822 3424 : if (hAlphaBand == nullptr)
823 0 : return CE_Failure;
824 :
825 3424 : size_t iPixel = 0;
826 :
827 : /* -------------------------------------------------------------------- */
828 : /* Read alpha case. */
829 : /* -------------------------------------------------------------------- */
830 3424 : if (nBandCount >= 0)
831 : {
832 : const char *pszInitDest =
833 1712 : CSLFetchNameValue(psWO->papszWarpOptions, "INIT_DEST");
834 :
835 : // Special logic for destinations being initialized on-the-fly.
836 1712 : if (pszInitDest != nullptr)
837 : {
838 883 : memset(pafMask, 0, nPixels * sizeof(float));
839 883 : return CE_None;
840 : }
841 :
842 : // Rescale.
843 829 : const float inv_alpha_max = static_cast<float>(
844 829 : 1.0 / CPLAtof(CSLFetchNameValueDef(psWO->papszWarpOptions,
845 829 : "DST_ALPHA_MAX", "255")));
846 :
847 : #if (defined(__x86_64) || defined(_M_X64))
848 829 : const GDALDataType eDT = GDALGetRasterDataType(hAlphaBand);
849 : // Make sure that pafMask is at least 8-byte aligned, which should
850 : // normally be always the case if being a ptr returned by malloc().
851 829 : if ((eDT == GDT_UInt8 || eDT == GDT_UInt16) &&
852 820 : CPL_IS_ALIGNED(pafMask, 8))
853 : {
854 : // Read data.
855 1640 : const CPLErr eErr = GDALRasterIOEx(
856 : hAlphaBand, GF_Read, nXOff, nYOff, nXSize, nYSize, pafMask,
857 : nXSize, nYSize, eDT, static_cast<GSpacing>(sizeof(int)),
858 820 : static_cast<GSpacing>(sizeof(int)) * nXSize, nullptr);
859 :
860 820 : if (eErr != CE_None)
861 0 : return eErr;
862 :
863 : // Make sure we have the correct alignment before doing SSE
864 : // On Linux x86_64, the alignment should be always correct due
865 : // the alignment of malloc() being 16 byte.
866 820 : const GUInt32 mask = (eDT == GDT_UInt8) ? 0xff : 0xffff;
867 820 : if (!CPL_IS_ALIGNED(pafMask, 16))
868 : {
869 0 : pafMask[iPixel] =
870 0 : (reinterpret_cast<GUInt32 *>(pafMask)[iPixel] & mask) *
871 0 : inv_alpha_max;
872 0 : pafMask[iPixel] = std::min(1.0f, pafMask[iPixel]);
873 0 : iPixel++;
874 : }
875 820 : CPLAssert(CPL_IS_ALIGNED(pafMask + iPixel, 16));
876 820 : const __m128 xmm_inverse_alpha_max = _mm_load1_ps(&inv_alpha_max);
877 820 : const float one_single = 1.0f;
878 820 : const __m128 xmm_one = _mm_load1_ps(&one_single);
879 1640 : const __m128i xmm_i_mask = _mm_set1_epi32(mask);
880 1677150 : for (; iPixel + 31 < nPixels; iPixel += 32)
881 : {
882 3352660 : __m128 xmm_mask0 = _mm_cvtepi32_ps(_mm_and_si128(
883 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
884 1676330 : pafMask + iPixel + 4 * 0))));
885 3352660 : __m128 xmm_mask1 = _mm_cvtepi32_ps(_mm_and_si128(
886 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
887 1676330 : pafMask + iPixel + 4 * 1))));
888 3352660 : __m128 xmm_mask2 = _mm_cvtepi32_ps(_mm_and_si128(
889 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
890 1676330 : pafMask + iPixel + 4 * 2))));
891 3352660 : __m128 xmm_mask3 = _mm_cvtepi32_ps(_mm_and_si128(
892 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
893 1676330 : pafMask + iPixel + 4 * 3))));
894 3352660 : __m128 xmm_mask4 = _mm_cvtepi32_ps(_mm_and_si128(
895 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
896 1676330 : pafMask + iPixel + 4 * 4))));
897 3352660 : __m128 xmm_mask5 = _mm_cvtepi32_ps(_mm_and_si128(
898 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
899 1676330 : pafMask + iPixel + 4 * 5))));
900 3352660 : __m128 xmm_mask6 = _mm_cvtepi32_ps(_mm_and_si128(
901 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
902 1676330 : pafMask + iPixel + 4 * 6))));
903 5028980 : __m128 xmm_mask7 = _mm_cvtepi32_ps(_mm_and_si128(
904 : xmm_i_mask, _mm_load_si128(reinterpret_cast<__m128i *>(
905 1676330 : pafMask + iPixel + 4 * 7))));
906 1676330 : xmm_mask0 = _mm_mul_ps(xmm_mask0, xmm_inverse_alpha_max);
907 1676330 : xmm_mask1 = _mm_mul_ps(xmm_mask1, xmm_inverse_alpha_max);
908 1676330 : xmm_mask2 = _mm_mul_ps(xmm_mask2, xmm_inverse_alpha_max);
909 1676330 : xmm_mask3 = _mm_mul_ps(xmm_mask3, xmm_inverse_alpha_max);
910 1676330 : xmm_mask4 = _mm_mul_ps(xmm_mask4, xmm_inverse_alpha_max);
911 1676330 : xmm_mask5 = _mm_mul_ps(xmm_mask5, xmm_inverse_alpha_max);
912 1676330 : xmm_mask6 = _mm_mul_ps(xmm_mask6, xmm_inverse_alpha_max);
913 1676330 : xmm_mask7 = _mm_mul_ps(xmm_mask7, xmm_inverse_alpha_max);
914 1676330 : xmm_mask0 = _mm_min_ps(xmm_mask0, xmm_one);
915 1676330 : xmm_mask1 = _mm_min_ps(xmm_mask1, xmm_one);
916 1676330 : xmm_mask2 = _mm_min_ps(xmm_mask2, xmm_one);
917 1676330 : xmm_mask3 = _mm_min_ps(xmm_mask3, xmm_one);
918 1676330 : xmm_mask4 = _mm_min_ps(xmm_mask4, xmm_one);
919 1676330 : xmm_mask5 = _mm_min_ps(xmm_mask5, xmm_one);
920 1676330 : xmm_mask6 = _mm_min_ps(xmm_mask6, xmm_one);
921 1676330 : xmm_mask7 = _mm_min_ps(xmm_mask7, xmm_one);
922 1676330 : _mm_store_ps(pafMask + iPixel + 4 * 0, xmm_mask0);
923 1676330 : _mm_store_ps(pafMask + iPixel + 4 * 1, xmm_mask1);
924 1676330 : _mm_store_ps(pafMask + iPixel + 4 * 2, xmm_mask2);
925 1676330 : _mm_store_ps(pafMask + iPixel + 4 * 3, xmm_mask3);
926 1676330 : _mm_store_ps(pafMask + iPixel + 4 * 4, xmm_mask4);
927 1676330 : _mm_store_ps(pafMask + iPixel + 4 * 5, xmm_mask5);
928 1676330 : _mm_store_ps(pafMask + iPixel + 4 * 6, xmm_mask6);
929 1676330 : _mm_store_ps(pafMask + iPixel + 4 * 7, xmm_mask7);
930 : }
931 1564 : for (; iPixel < nPixels; iPixel++)
932 : {
933 744 : pafMask[iPixel] =
934 744 : (reinterpret_cast<GUInt32 *>(pafMask)[iPixel] & mask) *
935 744 : inv_alpha_max;
936 744 : pafMask[iPixel] = std::min(1.0f, pafMask[iPixel]);
937 820 : }
938 : }
939 : else
940 : #endif
941 : {
942 : // Read data.
943 : const CPLErr eErr =
944 9 : GDALRasterIO(hAlphaBand, GF_Read, nXOff, nYOff, nXSize, nYSize,
945 : pafMask, nXSize, nYSize, GDT_Float32, 0, 0);
946 :
947 9 : if (eErr != CE_None)
948 0 : return eErr;
949 :
950 842 : for (; iPixel < nPixels; iPixel++)
951 : {
952 833 : pafMask[iPixel] = pafMask[iPixel] * inv_alpha_max;
953 833 : pafMask[iPixel] = std::min(1.0f, pafMask[iPixel]);
954 : }
955 : }
956 :
957 829 : return CE_None;
958 : }
959 :
960 : /* -------------------------------------------------------------------- */
961 : /* Write alpha case. */
962 : /* -------------------------------------------------------------------- */
963 : else
964 : {
965 1712 : GDALDataType eDT = GDALGetRasterDataType(hAlphaBand);
966 : const float cst_alpha_max =
967 1712 : static_cast<float>(CPLAtof(CSLFetchNameValueDef(
968 1712 : psWO->papszWarpOptions, "DST_ALPHA_MAX", "255"))) +
969 30 : ((eDT == GDT_UInt8 || eDT == GDT_Int16 || eDT == GDT_UInt16 ||
970 1 : eDT == GDT_Int32 || eDT == GDT_UInt32)
971 1742 : ? 0.1f
972 1712 : : 0.0f);
973 :
974 1712 : CPLErr eErr = CE_None;
975 :
976 : #if (defined(__x86_64) || defined(_M_X64))
977 : // Make sure that pafMask is at least 8-byte aligned, which should
978 : // normally be always the case if being a ptr returned by malloc()
979 1712 : if ((eDT == GDT_UInt8 || eDT == GDT_Int16 || eDT == GDT_UInt16) &&
980 1711 : CPL_IS_ALIGNED(pafMask, 8))
981 : {
982 : // Make sure we have the correct alignment before doing SSE
983 : // On Linux x86_64, the alignment should be always correct due
984 : // the alignment of malloc() being 16 byte
985 1711 : if (!CPL_IS_ALIGNED(pafMask, 16))
986 : {
987 0 : reinterpret_cast<int *>(pafMask)[iPixel] =
988 0 : static_cast<int>(pafMask[iPixel] * cst_alpha_max);
989 0 : iPixel++;
990 : }
991 1711 : CPLAssert(CPL_IS_ALIGNED(pafMask + iPixel, 16));
992 1711 : const __m128 xmm_alpha_max = _mm_load1_ps(&cst_alpha_max);
993 3347320 : for (; iPixel + 31 < nPixels; iPixel += 32)
994 : {
995 3345610 : __m128 xmm_mask0 = _mm_load_ps(pafMask + iPixel + 4 * 0);
996 3345610 : __m128 xmm_mask1 = _mm_load_ps(pafMask + iPixel + 4 * 1);
997 3345610 : __m128 xmm_mask2 = _mm_load_ps(pafMask + iPixel + 4 * 2);
998 3345610 : __m128 xmm_mask3 = _mm_load_ps(pafMask + iPixel + 4 * 3);
999 3345610 : __m128 xmm_mask4 = _mm_load_ps(pafMask + iPixel + 4 * 4);
1000 3345610 : __m128 xmm_mask5 = _mm_load_ps(pafMask + iPixel + 4 * 5);
1001 3345610 : __m128 xmm_mask6 = _mm_load_ps(pafMask + iPixel + 4 * 6);
1002 6691220 : __m128 xmm_mask7 = _mm_load_ps(pafMask + iPixel + 4 * 7);
1003 3345610 : xmm_mask0 = _mm_mul_ps(xmm_mask0, xmm_alpha_max);
1004 3345610 : xmm_mask1 = _mm_mul_ps(xmm_mask1, xmm_alpha_max);
1005 3345610 : xmm_mask2 = _mm_mul_ps(xmm_mask2, xmm_alpha_max);
1006 3345610 : xmm_mask3 = _mm_mul_ps(xmm_mask3, xmm_alpha_max);
1007 3345610 : xmm_mask4 = _mm_mul_ps(xmm_mask4, xmm_alpha_max);
1008 3345610 : xmm_mask5 = _mm_mul_ps(xmm_mask5, xmm_alpha_max);
1009 3345610 : xmm_mask6 = _mm_mul_ps(xmm_mask6, xmm_alpha_max);
1010 3345610 : xmm_mask7 = _mm_mul_ps(xmm_mask7, xmm_alpha_max);
1011 : // Truncate to int.
1012 3345610 : _mm_store_si128(
1013 3345610 : reinterpret_cast<__m128i *>(pafMask + iPixel + 4 * 0),
1014 : _mm_cvttps_epi32(xmm_mask0));
1015 3345610 : _mm_store_si128(
1016 3345610 : reinterpret_cast<__m128i *>(pafMask + iPixel + 4 * 1),
1017 : _mm_cvttps_epi32(xmm_mask1));
1018 3345610 : _mm_store_si128(
1019 3345610 : reinterpret_cast<__m128i *>(pafMask + iPixel + 4 * 2),
1020 : _mm_cvttps_epi32(xmm_mask2));
1021 3345610 : _mm_store_si128(
1022 3345610 : reinterpret_cast<__m128i *>(pafMask + iPixel + 4 * 3),
1023 : _mm_cvttps_epi32(xmm_mask3));
1024 3345610 : _mm_store_si128(
1025 3345610 : reinterpret_cast<__m128i *>(pafMask + iPixel + 4 * 4),
1026 : _mm_cvttps_epi32(xmm_mask4));
1027 3345610 : _mm_store_si128(
1028 3345610 : reinterpret_cast<__m128i *>(pafMask + iPixel + 4 * 5),
1029 : _mm_cvttps_epi32(xmm_mask5));
1030 3345610 : _mm_store_si128(
1031 3345610 : reinterpret_cast<__m128i *>(pafMask + iPixel + 4 * 6),
1032 : _mm_cvttps_epi32(xmm_mask6));
1033 3345610 : _mm_store_si128(
1034 3345610 : reinterpret_cast<__m128i *>(pafMask + iPixel + 4 * 7),
1035 : _mm_cvttps_epi32(xmm_mask7));
1036 : }
1037 3051 : for (; iPixel < nPixels; iPixel++)
1038 1340 : reinterpret_cast<int *>(pafMask)[iPixel] =
1039 1340 : static_cast<int>(pafMask[iPixel] * cst_alpha_max);
1040 :
1041 : // Write data.
1042 : // Assumes little endianness here.
1043 3422 : eErr = GDALRasterIOEx(
1044 : hAlphaBand, GF_Write, nXOff, nYOff, nXSize, nYSize, pafMask,
1045 : nXSize, nYSize, eDT, static_cast<GSpacing>(sizeof(int)),
1046 1711 : static_cast<GSpacing>(sizeof(int)) * nXSize, nullptr);
1047 : }
1048 : else
1049 : #endif
1050 : {
1051 9 : for (; iPixel + 3 < nPixels; iPixel += 4)
1052 : {
1053 8 : pafMask[iPixel + 0] = static_cast<float>(
1054 8 : static_cast<int>(pafMask[iPixel + 0] * cst_alpha_max));
1055 8 : pafMask[iPixel + 1] = static_cast<float>(
1056 8 : static_cast<int>(pafMask[iPixel + 1] * cst_alpha_max));
1057 8 : pafMask[iPixel + 2] = static_cast<float>(
1058 8 : static_cast<int>(pafMask[iPixel + 2] * cst_alpha_max));
1059 8 : pafMask[iPixel + 3] = static_cast<float>(
1060 8 : static_cast<int>(pafMask[iPixel + 3] * cst_alpha_max));
1061 : }
1062 2 : for (; iPixel < nPixels; iPixel++)
1063 1 : pafMask[iPixel] = static_cast<float>(
1064 1 : static_cast<int>(pafMask[iPixel] * cst_alpha_max));
1065 :
1066 : // Write data.
1067 :
1068 : eErr =
1069 1 : GDALRasterIO(hAlphaBand, GF_Write, nXOff, nYOff, nXSize, nYSize,
1070 : pafMask, nXSize, nYSize, GDT_Float32, 0, 0);
1071 : }
1072 1712 : return eErr;
1073 : }
1074 : }
1075 :
1076 : /************************************************************************/
1077 : /* GDALWarpGetOptionList() */
1078 : /************************************************************************/
1079 :
1080 : /** Return a XML string describing options accepted by
1081 : * GDALWarpOptions::papszWarpOptions.
1082 : *
1083 : * @since 3.11
1084 : */
1085 1787 : const char *GDALWarpGetOptionList(void)
1086 : {
1087 : return "<OptionList>"
1088 : "<Option name='INIT_DEST' type='string' description='"
1089 : "Numeric value or NO_DATA. This option forces the destination image "
1090 : "to be initialized to the indicated value (for all bands) "
1091 : "or indicates that it should be initialized to the NO_DATA value in "
1092 : "padfDstNoDataReal/padfDstNoDataImag. If this value is not set, the "
1093 : "destination image will be read and overlaid.'/>"
1094 : "<Option name='RESET_DEST_PIXELS' type='boolean' description='"
1095 : "Whether the whole destination image must be re-initialized to the "
1096 : "destination nodata value of padfDstNoDataReal/padfDstNoDataImag "
1097 : "if set, or 0 otherwise. The main difference with INIT_DEST is that "
1098 : "it also affects regions not related to the source dataset.' "
1099 : "default='NO'/>"
1100 : "<Option name='WRITE_FLUSH' type='boolean' description='"
1101 : "This option forces a flush to disk of data after "
1102 : "each chunk is processed. In some cases this helps ensure a serial "
1103 : " writing of the output data otherwise a block of data may be "
1104 : "written to disk each time a block of data is read for the input "
1105 : "buffer resulting in a lot of extra seeking around the disk, and "
1106 : "reduced IO throughput.' default='NO'/>"
1107 : "<Option name='SKIP_NOSOURCE' type='boolean' description='"
1108 : "Skip all processing for chunks for which there is no corresponding "
1109 : "input data. This will disable initializing the destination "
1110 : "(INIT_DEST) and all other processing, and so should be used "
1111 : "carefully. Mostly useful to short circuit a lot of extra work "
1112 : "in mosaicing situations. gdalwarp will automatically enable this "
1113 : "option when it is assumed to be safe to do so.' default='NO'/>"
1114 : #ifdef undocumented
1115 : "<Option name='ERROR_OUT_IF_EMPTY_SOURCE_WINDOW' type='boolean' "
1116 : "description='By default, if the source window corresponding to the "
1117 : "current target window fails to be determined due to reprojection "
1118 : "errors, the warping fails. Setting this option to NO prevent such "
1119 : "failure from happening. The warped VRT mechanism automatically "
1120 : "sets it to NO.'/>"
1121 : #endif
1122 : "<Option name='UNIFIED_SRC_NODATA' type='string-select' "
1123 : "description='"
1124 : "This setting determines how to take into account nodata values "
1125 : "when there are several input bands. Consult "
1126 : "GDALWarpOptions::papszWarpOptions documentation for more details.'>"
1127 : " <Value>AUTO</Value>"
1128 : " <Value>PARTIAL</Value>"
1129 : " <Value>YES</Value>"
1130 : " <Value>NO</Value>"
1131 : "</Option>"
1132 : "<Option name='CUTLINE' type='string' description='"
1133 : "This may contain the WKT geometry for a cutline. It will be "
1134 : "converted into a geometry by GDALWarpOperation::Initialize() and "
1135 : "assigned to the GDALWarpOptions hCutline field. The coordinates "
1136 : "must be expressed in source pixel/line coordinates. Note: this is "
1137 : "different from the assumptions made for the -cutline option "
1138 : "of the gdalwarp utility !'/>"
1139 : "<Option name='CUTLINE_BLEND_DIST' type='float' description='"
1140 : "This may be set with a distance in pixels which will be assigned "
1141 : "to the dfCutlineBlendDist field in the GDALWarpOptions.'/>"
1142 : "<Option name='CUTLINE_ALL_TOUCHED' type='boolean' description='"
1143 : "This may be set to TRUE to enable ALL_TOUCHED mode when "
1144 : "rasterizing cutline polygons. This is useful to ensure that that "
1145 : "all pixels overlapping the cutline polygon will be selected, not "
1146 : "just those whose center point falls within the polygon.' "
1147 : "default='NO'/>"
1148 : "<Option name='XSCALE' type='float' description='"
1149 : "Ratio expressing the resampling factor (number of destination "
1150 : "pixels per source pixel) along the target horizontal axis. The "
1151 : "scale is used to determine the number of source pixels along the "
1152 : "x-axis that are considered by the resampling algorithm. "
1153 : "Equals to one for no resampling, below one for downsampling "
1154 : "and above one for upsampling. This is automatically computed, "
1155 : "for each processing chunk, and may thus vary among them, depending "
1156 : "on the shape of output regions vs input regions. Such variations "
1157 : "can be undesired in some situations. If the resampling factor "
1158 : "can be considered as constant over the warped area, setting a "
1159 : "constant value can lead to more reproducible pixel output.'/>"
1160 : "<Option name='YSCALE' type='float' description='"
1161 : "Same as XSCALE, but along the horizontal axis.'/>"
1162 : "<Option name='OPTIMIZE_SIZE' type='boolean' description='"
1163 : "This defaults to FALSE, but may be set to TRUE typically when "
1164 : "writing to a compressed dataset (GeoTIFF with COMPRESS creation "
1165 : "option set for example) for achieving a smaller file size. This "
1166 : "is achieved by writing at once data aligned on full blocks of the "
1167 : "target dataset, which avoids partial writes of compressed blocks "
1168 : "and lost space when they are rewritten at the end of the file. "
1169 : "However sticking to target block size may cause major processing "
1170 : "slowdown for some particular reprojections. OPTIMIZE_SIZE mode "
1171 : "is automatically enabled when it is safe to do so. "
1172 : "As this parameter influences the shape of warping chunk, and by "
1173 : "default the XSCALE and YSCALE parameters are computed per warping "
1174 : "chunk, this parameter may influence the pixel output.' "
1175 : "default='NO'/>"
1176 : "<Option name='NUM_THREADS' type='string' description='"
1177 : "Can be set to a numeric value or ALL_CPUS to set the number of "
1178 : "threads to use to parallelize the computation part of the warping. "
1179 : "If not set, computation will be done in a single thread..'/>"
1180 : "<Option name='STREAMABLE_OUTPUT' type='boolean' description='"
1181 : "This defaults to FALSE, but may be set to TRUE typically when "
1182 : "writing to a streamed file. The gdalwarp utility automatically "
1183 : "sets this option when writing to /vsistdout/ or a named pipe "
1184 : "(on Unix). This option has performance impacts for some "
1185 : "reprojections. Note: band interleaved output is "
1186 : "not currently supported by the warping algorithm in a streamable "
1187 : "compatible way.' default='NO'/>"
1188 : "<Option name='SRC_COORD_PRECISION' type='float' description='"
1189 : "Advanced setting. This defaults to 0, to indicate that no rounding "
1190 : "of computing source image coordinates corresponding to the target "
1191 : "image must be done. If greater than 0 (and typically below 1), "
1192 : "this value, expressed in pixel, will be used to round computed "
1193 : "source image coordinates. The purpose of this option is to make "
1194 : "the results of warping with the approximated transformer more "
1195 : "reproducible and not sensitive to changes in warping memory size. "
1196 : "To achieve that, SRC_COORD_PRECISION must be at least 10 times "
1197 : "greater than the error threshold. The higher the "
1198 : "SRC_COORD_PRECISION/error_threshold ratio, the higher the "
1199 : "performance will be, since exact reprojections must statistically "
1200 : "be done with a frequency of "
1201 : "4*error_threshold/SRC_COORD_PRECISION.' default='0'/>"
1202 : "<Option name='SRC_ALPHA_MAX' type='float' description='"
1203 : "Maximum value for the alpha band of the source dataset. If the "
1204 : "value is not set and the alpha band has a NBITS metadata item, "
1205 : "it is used to set SRC_ALPHA_MAX = 2^NBITS-1. Otherwise, if the "
1206 : "value is not set and the alpha band is of type UInt16 "
1207 : "(resp Int16), 65535 (resp 32767) is used. "
1208 : "Otherwise, 255 is used.'/>"
1209 : "<Option name='DST_ALPHA_MAX' type='float' description='"
1210 : "Maximum value for the alpha band of the destination dataset. "
1211 : "If the value is not set and the alpha band has a NBITS metadata "
1212 : "item, it is used to set SRC_ALPHA_MAX = 2^NBITS-1. Otherwise, if "
1213 : "the value is not set and the alpha band is of type UInt16 "
1214 : "(resp Int16), 65535 (resp 32767) is used. "
1215 : "Otherwise, 255 is used.'/>"
1216 : "<Option name='SAMPLE_GRID' type='boolean' description='"
1217 : "Setting this option to YES will force the sampling to "
1218 : "include internal points as well as edge points which can be "
1219 : "important if the transformation is esoteric inside out, or if "
1220 : "large sections of the destination image are not transformable into "
1221 : "the source coordinate system.' default='NO'/>"
1222 : "<Option name='SAMPLE_STEPS' type='string' description='"
1223 : "Modifies the density of the sampling grid. Increasing this can "
1224 : "increase the computational cost, but improves the accuracy with "
1225 : "which the source region is computed. This can be set to ALL to "
1226 : "mean to sample along all edge points of the destination region "
1227 : "(if SAMPLE_GRID=NO or not specified), or all points of the "
1228 : "destination region if SAMPLE_GRID=YES.' default='21'/>"
1229 : "<Option name='SOURCE_EXTRA' type='int' description='"
1230 : "This is a number of extra pixels added around the source "
1231 : "window for a given request, and by default it is 1 to take care "
1232 : "of rounding error. Setting this larger will increase the amount of "
1233 : "data that needs to be read, but can avoid missing source data.' "
1234 : "default='1'/>"
1235 : "<Option name='APPLY_VERTICAL_SHIFT' type='boolean' description='"
1236 : "Force the use of vertical shift. This option is generally not "
1237 : "necessary, except when using an explicit coordinate transformation "
1238 : "(COORDINATE_OPERATION), and not specifying an explicit source and "
1239 : "target SRS.'/>"
1240 : "<Option name='MULT_FACTOR_VERTICAL_SHIFT' type='float' "
1241 : "description='"
1242 : "Multiplication factor for the vertical shift' default='1.0'/>"
1243 : "<Option name='EXCLUDED_VALUES' type='string' "
1244 : "description='"
1245 : "Comma-separated tuple of values (thus typically \"R,G,B\"), that "
1246 : "are ignored as contributing source pixels during resampling. "
1247 : "The number of values in the tuple must be the same as the number "
1248 : "of bands, excluding the alpha band. Several tuples of excluded "
1249 : "values may be specified using the \"(R1,G1,B2),(R2,G2,B2)\" syntax."
1250 : " Only taken into account by Average currently. This concept is a "
1251 : "bit similar to nodata/alpha, but the main difference is that "
1252 : "pixels matching one of the excluded value tuples are still "
1253 : "considered as valid, when determining the target pixel "
1254 : "validity/density.'/>"
1255 : "<Option name='EXCLUDED_VALUES_PCT_THRESHOLD' type='float' "
1256 : "min='0' max='100' description='"
1257 : "Minimum percentage of source pixels that must be set at one of "
1258 : "the EXCLUDED_VALUES to cause the excluded value, that is in "
1259 : "majority among source pixels, to be used as the target pixel "
1260 : "value. Only taken into account by Average currently.' "
1261 : "default='50'/>"
1262 : "<Option name='NODATA_VALUES_PCT_THRESHOLD' type='float' "
1263 : "min='0' max='100' description='"
1264 : "Minimum percentage of source pixels that must be at nodata (or "
1265 : "alpha=0 or any other way to express transparent pixel) to cause "
1266 : "the target pixel value to not be set. Default value is 100 (%), "
1267 : "which means that a target pixel is not set only if all "
1268 : "contributing source pixels are not set. Note that "
1269 : "NODATA_VALUES_PCT_THRESHOLD is taken into account before "
1270 : "EXCLUDED_VALUES_PCT_THRESHOLD. Only taken into account by Average "
1271 : "currently.' default='100'/>"
1272 : "<Option name='MODE_TIES' type='string-select' "
1273 : "description='"
1274 : "Strategy to use when breaking ties with MODE resampling. "
1275 : "By default, the first value encountered will be used. "
1276 : "Alternatively, the minimum or maximum value can be selected.' "
1277 : "default='FIRST'>"
1278 : " <Value>FIRST</Value>"
1279 : " <Value>MIN</Value>"
1280 : " <Value>MAX</Value>"
1281 : "</Option>"
1282 1787 : "</OptionList>";
1283 : }
1284 :
1285 : /************************************************************************/
1286 : /* ==================================================================== */
1287 : /* GDALWarpOptions */
1288 : /* ==================================================================== */
1289 : /************************************************************************/
1290 :
1291 : /**
1292 : * \var char **GDALWarpOptions::papszWarpOptions;
1293 : *
1294 : * A string list of additional options controlling the warp operation in
1295 : * name=value format. A suitable string list can be prepared with
1296 : * CSLSetNameValue().
1297 : *
1298 : * The available options can also be retrieved programmatically with
1299 : * GDALWarpGetOptionList().
1300 : *
1301 : * The following values are currently supported:
1302 : * <ul>
1303 : * <li>INIT_DEST=[value] or INIT_DEST=NO_DATA: This option forces the
1304 : * destination image to be initialized to the indicated value (for all bands)
1305 : * or indicates that it should be initialized to the NO_DATA value in
1306 : * padfDstNoDataReal/padfDstNoDataImag. If this value isn't set the
1307 : * destination image will be read and overlaid.</li>
1308 : *
1309 : * <li>RESET_DEST_PIXELS=YES/NO (since GDAL 3.13): Defaults to NO.
1310 : * Whether the whole destination image must be re-initialized to the
1311 : * destination nodata value of padfDstNoDataReal/padfDstNoDataImag if set,
1312 : * or 0 otherwise.
1313 : * The main difference with INIT_DEST is that it also affects regions
1314 : * not related to the source dataset.
1315 : * </li>
1316 : *
1317 : * <li>WRITE_FLUSH=YES/NO: This option forces a flush to disk of data after
1318 : * each chunk is processed. In some cases this helps ensure a serial
1319 : * writing of the output data otherwise a block of data may be written to disk
1320 : * each time a block of data is read for the input buffer resulting in a lot
1321 : * of extra seeking around the disk, and reduced IO throughput. The default
1322 : * is NO.</li>
1323 : *
1324 : * <li>SKIP_NOSOURCE=YES/NO: Skip all processing for chunks for which there
1325 : * is no corresponding input data. This will disable initializing the
1326 : * destination (INIT_DEST) and all other processing, and so should be used
1327 : * carefully. Mostly useful to short circuit a lot of extra work in mosaicing
1328 : * situations. gdalwarp will automatically enable this
1329 : * option when it is assumed to be safe to do so.</li>
1330 : *
1331 : * <li>UNIFIED_SRC_NODATA=YES/NO/PARTIAL: This setting determines
1332 : * how to take into account nodata values when there are several input bands.
1333 : * <ul>
1334 : * <li>When YES, all bands are considered as nodata if and only if, all bands
1335 : * match the corresponding nodata values.
1336 : * Note: UNIFIED_SRC_NODATA=YES is set by default, when called from gdalwarp
1337 : * / GDALWarp() with an explicit -srcnodata setting.
1338 : *
1339 : * Example with nodata values at (1, 2, 3) and target alpha band requested.
1340 : * <ul>
1341 : * <li>input pixel = (1, 2, 3) ==> output pixel = (0, 0, 0, 0)</li>
1342 : * <li>input pixel = (1, 2, 127) ==> output pixel = (1, 2, 127, 255)</li>
1343 : * </ul>
1344 : * </li>
1345 : * <li>When NO, nodata masking values is considered independently for each band.
1346 : * A potential target alpha band will always be valid if there are multiple
1347 : * bands.
1348 : *
1349 : * Example with nodata values at (1, 2, 3) and target alpha band requested.
1350 : * <ul>
1351 : * <li>input pixel = (1, 2, 3) ==> output pixel = (0, 0, 0, 255)</li>
1352 : * <li>input pixel = (1, 2, 127) ==> output pixel = (0, 0, 127, 255)</li>
1353 : * </ul>
1354 : *
1355 : * Note: NO was the default behavior before GDAL 3.3.2
1356 : * </li>
1357 : * <li>When PARTIAL, or not specified at all (default behavior),
1358 : * nodata masking values is considered independently for each band.
1359 : * But, and this is the difference with NO, if for a given pixel, it
1360 : * evaluates to the nodata value of each band, the target pixel is
1361 : * considered as globally invalid, which impacts the value of a potential
1362 : * target alpha band.
1363 : *
1364 : * Note: PARTIAL is new to GDAL 3.3.2 and should not be used with
1365 : * earlier versions. The default behavior of GDAL < 3.3.2 was NO.
1366 : *
1367 : * Example with nodata values at (1, 2, 3) and target alpha band requested.
1368 : * <ul>
1369 : * <li>input pixel = (1, 2, 3) ==> output pixel = (0, 0, 0, 0)</li>
1370 : * <li>input pixel = (1, 2, 127) ==> output pixel = (0, 0, 127, 255)</li>
1371 : * </ul>
1372 : * </li>
1373 : * </ul>
1374 : * </li>
1375 : *
1376 : * <li>CUTLINE: This may contain the WKT geometry for a cutline. It will
1377 : * be converted into a geometry by GDALWarpOperation::Initialize() and assigned
1378 : * to the GDALWarpOptions hCutline field. The coordinates must be expressed
1379 : * in source pixel/line coordinates. Note: this is different from the
1380 : * assumptions made for the -cutline option of the gdalwarp utility !</li>
1381 : *
1382 : * <li>CUTLINE_BLEND_DIST: This may be set with a distance in pixels which
1383 : * will be assigned to the dfCutlineBlendDist field in the GDALWarpOptions.</li>
1384 : *
1385 : * <li>CUTLINE_ALL_TOUCHED: This defaults to FALSE, but may be set to TRUE
1386 : * to enable ALL_TOUCHED mode when rasterizing cutline polygons. This is
1387 : * useful to ensure that that all pixels overlapping the cutline polygon
1388 : * will be selected, not just those whose center point falls within the
1389 : * polygon.</li>
1390 : *
1391 : * <li>XSCALE: Ratio expressing the resampling factor (number of destination
1392 : * pixels per source pixel) along the target horizontal axis.
1393 : * The scale is used to determine the number of source pixels along the x-axis
1394 : * that are considered by the resampling algorithm.
1395 : * Equals to one for no resampling, below one for downsampling
1396 : * and above one for upsampling. This is automatically computed, for each
1397 : * processing chunk, and may thus vary among them, depending on the
1398 : * shape of output regions vs input regions. Such variations can be undesired
1399 : * in some situations. If the resampling factor can be considered as constant
1400 : * over the warped area, setting a constant value can lead to more reproducible
1401 : * pixel output.</li>
1402 : *
1403 : * <li>YSCALE: Same as XSCALE, but along the horizontal axis.</li>
1404 : *
1405 : * <li>OPTIMIZE_SIZE: This defaults to FALSE, but may be set to TRUE
1406 : * typically when writing to a compressed dataset (GeoTIFF with
1407 : * COMPRESS creation option set for example) for achieving a smaller
1408 : * file size. This is achieved by writing at once data aligned on full
1409 : * blocks of the target dataset, which avoids partial writes of
1410 : * compressed blocks and lost space when they are rewritten at the end
1411 : * of the file. However sticking to target block size may cause major
1412 : * processing slowdown for some particular reprojections. Starting
1413 : * with GDAL 3.8, OPTIMIZE_SIZE mode is automatically enabled when it is safe
1414 : * to do so.
1415 : * As this parameter influences the shape of warping chunk, and by default the
1416 : * XSCALE and YSCALE parameters are computed per warping chunk, this parameter may
1417 : * influence the pixel output.
1418 : * </li>
1419 : *
1420 : * <li>NUM_THREADS: (GDAL >= 1.10) Can be set to a numeric value or ALL_CPUS to
1421 : * set the number of threads to use to parallelize the computation part of the
1422 : * warping. If not set, computation will be done in a single thread.</li>
1423 : *
1424 : * <li>STREAMABLE_OUTPUT: This defaults to FALSE, but may
1425 : * be set to TRUE typically when writing to a streamed file. The
1426 : * gdalwarp utility automatically sets this option when writing to
1427 : * /vsistdout/ or a named pipe (on Unix). This option has performance
1428 : * impacts for some reprojections. Note: band interleaved output is
1429 : * not currently supported by the warping algorithm in a streamable
1430 : * compatible way.</li>
1431 : *
1432 : * <li>SRC_COORD_PRECISION: Advanced setting. This
1433 : * defaults to 0, to indicate that no rounding of computing source
1434 : * image coordinates corresponding to the target image must be
1435 : * done. If greater than 0 (and typically below 1), this value,
1436 : * expressed in pixel, will be used to round computed source image
1437 : * coordinates. The purpose of this option is to make the results of
1438 : * warping with the approximated transformer more reproducible and not
1439 : * sensitive to changes in warping memory size. To achieve that,
1440 : * SRC_COORD_PRECISION must be at least 10 times greater than the
1441 : * error threshold. The higher the SRC_COORD_PRECISION/error_threshold
1442 : * ratio, the higher the performance will be, since exact
1443 : * reprojections must statistically be done with a frequency of
1444 : * 4*error_threshold/SRC_COORD_PRECISION.</li>
1445 : *
1446 : * <li>SRC_ALPHA_MAX: Maximum value for the alpha band of the
1447 : * source dataset. If the value is not set and the alpha band has a NBITS
1448 : * metadata item, it is used to set SRC_ALPHA_MAX = 2^NBITS-1. Otherwise, if the
1449 : * value is not set and the alpha band is of type UInt16 (resp Int16), 65535
1450 : * (resp 32767) is used. Otherwise, 255 is used.</li>
1451 : *
1452 : * <li>DST_ALPHA_MAX: Maximum value for the alpha band of the
1453 : * destination dataset. If the value is not set and the alpha band has a NBITS
1454 : * metadata item, it is used to set DST_ALPHA_MAX = 2^NBITS-1. Otherwise, if the
1455 : * value is not set and the alpha band is of type UInt16 (resp Int16), 65535
1456 : * (resp 32767) is used. Otherwise, 255 is used.</li>
1457 : * </ul>
1458 : *
1459 : * Normally when computing the source raster data to
1460 : * load to generate a particular output area, the warper samples transforms
1461 : * 21 points along each edge of the destination region back onto the source
1462 : * file, and uses this to compute a bounding window on the source image that
1463 : * is sufficient. Depending on the transformation in effect, the source
1464 : * window may be a bit too small, or even missing large areas. Problem
1465 : * situations are those where the transformation is very non-linear or
1466 : * "inside out". Examples are transforming from WGS84 to Polar Stereographic
1467 : * for areas around the pole, or transformations where some of the image is
1468 : * untransformable. The following options provide some additional control
1469 : * to deal with errors in computing the source window:
1470 : * <ul>
1471 : *
1472 : * <li>SAMPLE_GRID=YES/NO: Setting this option to YES will force the sampling to
1473 : * include internal points as well as edge points which can be important if
1474 : * the transformation is esoteric inside out, or if large sections of the
1475 : * destination image are not transformable into the source coordinate
1476 : * system.</li>
1477 : *
1478 : * <li>SAMPLE_STEPS: Modifies the density of the sampling grid. The default
1479 : * number of steps is 21. Increasing this can increase the computational
1480 : * cost, but improves the accuracy with which the source region is
1481 : * computed.
1482 : * Starting with GDAL 3.7, this can be set to ALL to mean to sample
1483 : * along all edge points of the destination region (if SAMPLE_GRID=NO or not
1484 : * specified), or all points of the destination region if SAMPLE_GRID=YES.</li>
1485 : *
1486 : * <li>SOURCE_EXTRA: This is a number of extra pixels added around the source
1487 : * window for a given request, and by default it is 1 to take care of rounding
1488 : * error. Setting this larger will increase the amount of data that needs to
1489 : * be read, but can avoid missing source data.</li>
1490 : * <li>APPLY_VERTICAL_SHIFT=YES/NO: Force the use of vertical shift.
1491 : * This option is generally not necessary, except when using an explicit
1492 : * coordinate transformation (COORDINATE_OPERATION), and not specifying
1493 : * an explicit source and target SRS.</li>
1494 : * <li>MULT_FACTOR_VERTICAL_SHIFT: Multiplication factor for the vertical
1495 : * shift. Default 1.0</li>
1496 : *
1497 : * <li>EXCLUDED_VALUES: (GDAL >= 3.9) Comma-separated tuple of values
1498 : * (thus typically "R,G,B"), that are ignored as contributing source
1499 : * pixels during resampling. The number of values in the tuple must be the same
1500 : * as the number of bands, excluding the alpha band.
1501 : * Several tuples of excluded values may be specified using the
1502 : * "(R1,G1,B2),(R2,G2,B2)" syntax.
1503 : * Only taken into account by Average currently.
1504 : * This concept is a bit similar to nodata/alpha, but the main difference is
1505 : * that pixels matching one of the excluded value tuples are still considered
1506 : * as valid, when determining the target pixel validity/density.
1507 : * </li>
1508 : *
1509 : * <li>EXCLUDED_VALUES_PCT_THRESHOLD=[0-100]: (GDAL >= 3.9) Minimum percentage
1510 : * of source pixels that must be set at one of the EXCLUDED_VALUES to cause
1511 : * the excluded value, that is in majority among source pixels, to be used as the
1512 : * target pixel value. Default value is 50 (%).
1513 : * Only taken into account by Average currently.</li>
1514 : *
1515 : * <li>NODATA_VALUES_PCT_THRESHOLD=[0-100]: (GDAL >= 3.9) Minimum percentage
1516 : * of source pixels that must be at nodata (or alpha=0 or any other way to express
1517 : * transparent pixel) to cause the target pixel value to not be set. Default
1518 : * value is 100 (%), which means that a target pixel is not set only if all
1519 : * contributing source pixels are not set.
1520 : * Note that NODATA_VALUES_PCT_THRESHOLD is taken into account before
1521 : * EXCLUDED_VALUES_PCT_THRESHOLD.
1522 : * Only taken into account by Average currently.</li>
1523 : *
1524 : * <li>MODE_TIES=FIRST/MIN/MAX: (GDAL >= 3.11) Strategy to use when breaking
1525 : * ties with MODE resampling. By default, the first value encountered will be used.
1526 : * Alternatively, the minimum or maximum value can be selected.</li>
1527 : *
1528 : * </ul>
1529 : */
1530 :
1531 : /************************************************************************/
1532 : /* GDALCreateWarpOptions() */
1533 : /************************************************************************/
1534 :
1535 : /** Create a warp options structure.
1536 : *
1537 : * Must be deallocated with GDALDestroyWarpOptions()
1538 : */
1539 3818 : GDALWarpOptions *CPL_STDCALL GDALCreateWarpOptions()
1540 :
1541 : {
1542 : GDALWarpOptions *psOptions =
1543 3818 : static_cast<GDALWarpOptions *>(CPLCalloc(sizeof(GDALWarpOptions), 1));
1544 :
1545 3818 : psOptions->nBandCount = 0;
1546 3818 : psOptions->eResampleAlg = GRA_NearestNeighbour;
1547 3818 : psOptions->pfnProgress = GDALDummyProgress;
1548 3818 : psOptions->eWorkingDataType = GDT_Unknown;
1549 3818 : psOptions->eTieStrategy = GWKTS_First;
1550 :
1551 3818 : return psOptions;
1552 : }
1553 :
1554 : /************************************************************************/
1555 : /* GDALDestroyWarpOptions() */
1556 : /************************************************************************/
1557 :
1558 : /** Destroy a warp options structure. */
1559 3818 : void CPL_STDCALL GDALDestroyWarpOptions(GDALWarpOptions *psOptions)
1560 :
1561 : {
1562 3818 : if (psOptions == nullptr)
1563 0 : return;
1564 :
1565 3818 : CSLDestroy(psOptions->papszWarpOptions);
1566 3818 : CPLFree(psOptions->panSrcBands);
1567 3818 : CPLFree(psOptions->panDstBands);
1568 3818 : CPLFree(psOptions->padfSrcNoDataReal);
1569 3818 : CPLFree(psOptions->padfSrcNoDataImag);
1570 3818 : CPLFree(psOptions->padfDstNoDataReal);
1571 3818 : CPLFree(psOptions->padfDstNoDataImag);
1572 3818 : CPLFree(psOptions->papfnSrcPerBandValidityMaskFunc);
1573 3818 : CPLFree(psOptions->papSrcPerBandValidityMaskFuncArg);
1574 :
1575 3818 : if (psOptions->hCutline != nullptr)
1576 54 : delete static_cast<OGRGeometry *>(psOptions->hCutline);
1577 :
1578 3818 : CPLFree(psOptions);
1579 : }
1580 :
1581 : #define COPY_MEM(target, type, count) \
1582 : do \
1583 : { \
1584 : if ((psSrcOptions->target) != nullptr && (count) != 0) \
1585 : { \
1586 : (psDstOptions->target) = \
1587 : static_cast<type *>(CPLMalloc(sizeof(type) * (count))); \
1588 : memcpy((psDstOptions->target), (psSrcOptions->target), \
1589 : sizeof(type) * (count)); \
1590 : } \
1591 : else \
1592 : (psDstOptions->target) = nullptr; \
1593 : } while (false)
1594 :
1595 : /************************************************************************/
1596 : /* GDALCloneWarpOptions() */
1597 : /************************************************************************/
1598 :
1599 : /** Clone a warp options structure.
1600 : *
1601 : * Must be deallocated with GDALDestroyWarpOptions()
1602 : */
1603 : GDALWarpOptions *CPL_STDCALL
1604 2067 : GDALCloneWarpOptions(const GDALWarpOptions *psSrcOptions)
1605 :
1606 : {
1607 2067 : GDALWarpOptions *psDstOptions = GDALCreateWarpOptions();
1608 :
1609 2067 : memcpy(psDstOptions, psSrcOptions, sizeof(GDALWarpOptions));
1610 :
1611 2067 : if (psSrcOptions->papszWarpOptions != nullptr)
1612 1846 : psDstOptions->papszWarpOptions =
1613 1846 : CSLDuplicate(psSrcOptions->papszWarpOptions);
1614 :
1615 2067 : COPY_MEM(panSrcBands, int, psSrcOptions->nBandCount);
1616 2067 : COPY_MEM(panDstBands, int, psSrcOptions->nBandCount);
1617 2067 : COPY_MEM(padfSrcNoDataReal, double, psSrcOptions->nBandCount);
1618 2067 : COPY_MEM(padfSrcNoDataImag, double, psSrcOptions->nBandCount);
1619 2067 : COPY_MEM(padfDstNoDataReal, double, psSrcOptions->nBandCount);
1620 2067 : COPY_MEM(padfDstNoDataImag, double, psSrcOptions->nBandCount);
1621 : // cppcheck-suppress pointerSize
1622 2067 : COPY_MEM(papfnSrcPerBandValidityMaskFunc, GDALMaskFunc,
1623 : psSrcOptions->nBandCount);
1624 2067 : psDstOptions->papSrcPerBandValidityMaskFuncArg = nullptr;
1625 :
1626 2067 : if (psSrcOptions->hCutline != nullptr)
1627 10 : psDstOptions->hCutline =
1628 10 : OGR_G_Clone(static_cast<OGRGeometryH>(psSrcOptions->hCutline));
1629 2067 : psDstOptions->dfCutlineBlendDist = psSrcOptions->dfCutlineBlendDist;
1630 :
1631 2067 : return psDstOptions;
1632 : }
1633 :
1634 : namespace
1635 : {
1636 138 : void InitNoData(int nBandCount, double **ppdNoDataReal, double dDataReal)
1637 : {
1638 138 : if (nBandCount <= 0)
1639 : {
1640 0 : return;
1641 : }
1642 138 : if (*ppdNoDataReal != nullptr)
1643 : {
1644 32 : return;
1645 : }
1646 :
1647 106 : *ppdNoDataReal =
1648 106 : static_cast<double *>(CPLMalloc(sizeof(double) * nBandCount));
1649 :
1650 244 : for (int i = 0; i < nBandCount; ++i)
1651 : {
1652 138 : (*ppdNoDataReal)[i] = dDataReal;
1653 : }
1654 : }
1655 : } // namespace
1656 :
1657 : /************************************************************************/
1658 : /* GDALWarpInitDstNoDataReal() */
1659 : /************************************************************************/
1660 :
1661 : /**
1662 : * \brief Initialize padfDstNoDataReal with specified value.
1663 : *
1664 : * @param psOptionsIn options to initialize.
1665 : * @param dNoDataReal value to initialize to.
1666 : *
1667 : */
1668 38 : void CPL_STDCALL GDALWarpInitDstNoDataReal(GDALWarpOptions *psOptionsIn,
1669 : double dNoDataReal)
1670 : {
1671 38 : VALIDATE_POINTER0(psOptionsIn, "GDALWarpInitDstNoDataReal");
1672 38 : InitNoData(psOptionsIn->nBandCount, &psOptionsIn->padfDstNoDataReal,
1673 : dNoDataReal);
1674 : }
1675 :
1676 : /************************************************************************/
1677 : /* GDALWarpInitSrcNoDataReal() */
1678 : /************************************************************************/
1679 :
1680 : /**
1681 : * \brief Initialize padfSrcNoDataReal with specified value.
1682 : *
1683 : * @param psOptionsIn options to initialize.
1684 : * @param dNoDataReal value to initialize to.
1685 : *
1686 : */
1687 34 : void CPL_STDCALL GDALWarpInitSrcNoDataReal(GDALWarpOptions *psOptionsIn,
1688 : double dNoDataReal)
1689 : {
1690 34 : VALIDATE_POINTER0(psOptionsIn, "GDALWarpInitSrcNoDataReal");
1691 34 : InitNoData(psOptionsIn->nBandCount, &psOptionsIn->padfSrcNoDataReal,
1692 : dNoDataReal);
1693 : }
1694 :
1695 : /************************************************************************/
1696 : /* GDALWarpInitNoDataReal() */
1697 : /************************************************************************/
1698 :
1699 : /**
1700 : * \brief Initialize padfSrcNoDataReal and padfDstNoDataReal with specified
1701 : * value.
1702 : *
1703 : * @param psOptionsIn options to initialize.
1704 : * @param dNoDataReal value to initialize to.
1705 : *
1706 : */
1707 1 : void CPL_STDCALL GDALWarpInitNoDataReal(GDALWarpOptions *psOptionsIn,
1708 : double dNoDataReal)
1709 : {
1710 1 : GDALWarpInitDstNoDataReal(psOptionsIn, dNoDataReal);
1711 1 : GDALWarpInitSrcNoDataReal(psOptionsIn, dNoDataReal);
1712 1 : }
1713 :
1714 : /************************************************************************/
1715 : /* GDALWarpInitDstNoDataImag() */
1716 : /************************************************************************/
1717 :
1718 : /**
1719 : * \brief Initialize padfDstNoDataImag with specified value.
1720 : *
1721 : * @param psOptionsIn options to initialize.
1722 : * @param dNoDataImag value to initialize to.
1723 : *
1724 : */
1725 35 : void CPL_STDCALL GDALWarpInitDstNoDataImag(GDALWarpOptions *psOptionsIn,
1726 : double dNoDataImag)
1727 : {
1728 35 : VALIDATE_POINTER0(psOptionsIn, "GDALWarpInitDstNoDataImag");
1729 35 : InitNoData(psOptionsIn->nBandCount, &psOptionsIn->padfDstNoDataImag,
1730 : dNoDataImag);
1731 : }
1732 :
1733 : /************************************************************************/
1734 : /* GDALWarpInitSrcNoDataImag() */
1735 : /************************************************************************/
1736 :
1737 : /**
1738 : * \brief Initialize padfSrcNoDataImag with specified value.
1739 : *
1740 : * @param psOptionsIn options to initialize.
1741 : * @param dNoDataImag value to initialize to.
1742 : *
1743 : */
1744 31 : void CPL_STDCALL GDALWarpInitSrcNoDataImag(GDALWarpOptions *psOptionsIn,
1745 : double dNoDataImag)
1746 : {
1747 31 : VALIDATE_POINTER0(psOptionsIn, "GDALWarpInitSrcNoDataImag");
1748 31 : InitNoData(psOptionsIn->nBandCount, &psOptionsIn->padfSrcNoDataImag,
1749 : dNoDataImag);
1750 : }
1751 :
1752 : /************************************************************************/
1753 : /* GDALWarpResolveWorkingDataType() */
1754 : /************************************************************************/
1755 :
1756 : /**
1757 : * \brief If the working data type is unknown, this method will determine
1758 : * a valid working data type to support the data in the src and dest
1759 : * data sets and any noData values.
1760 : *
1761 : * @param psOptions options to initialize.
1762 : *
1763 : */
1764 1855 : void CPL_STDCALL GDALWarpResolveWorkingDataType(GDALWarpOptions *psOptions)
1765 : {
1766 1855 : if (psOptions == nullptr)
1767 : {
1768 0 : return;
1769 : }
1770 : /* -------------------------------------------------------------------- */
1771 : /* If no working data type was provided, set one now. */
1772 : /* */
1773 : /* Ensure that the working data type can encapsulate any value */
1774 : /* in the target, source, and the no data for either. */
1775 : /* -------------------------------------------------------------------- */
1776 1855 : if (psOptions->eWorkingDataType != GDT_Unknown)
1777 : {
1778 471 : return;
1779 : }
1780 :
1781 : // If none of the provided input nodata values can be represented in the
1782 : // data type of the corresponding source band, ignore them.
1783 1384 : if (psOptions->hSrcDS && psOptions->padfSrcNoDataReal)
1784 : {
1785 140 : int nCountInvalidSrcNoDataReal = 0;
1786 349 : for (int iBand = 0; iBand < psOptions->nBandCount; iBand++)
1787 : {
1788 418 : GDALRasterBandH hSrcBand = GDALGetRasterBand(
1789 209 : psOptions->hSrcDS, psOptions->panSrcBands[iBand]);
1790 :
1791 418 : if (hSrcBand &&
1792 209 : !GDALIsValueExactAs(psOptions->padfSrcNoDataReal[iBand],
1793 : GDALGetRasterDataType(hSrcBand)))
1794 : {
1795 2 : nCountInvalidSrcNoDataReal++;
1796 : }
1797 : }
1798 140 : if (nCountInvalidSrcNoDataReal == psOptions->nBandCount)
1799 : {
1800 2 : CPLFree(psOptions->padfSrcNoDataReal);
1801 2 : psOptions->padfSrcNoDataReal = nullptr;
1802 2 : CPLFree(psOptions->padfSrcNoDataImag);
1803 2 : psOptions->padfSrcNoDataImag = nullptr;
1804 : }
1805 : }
1806 :
1807 3161 : for (int iBand = 0; iBand < psOptions->nBandCount; iBand++)
1808 : {
1809 1777 : if (psOptions->hDstDS != nullptr)
1810 : {
1811 3464 : GDALRasterBandH hDstBand = GDALGetRasterBand(
1812 1732 : psOptions->hDstDS, psOptions->panDstBands[iBand]);
1813 :
1814 1732 : if (hDstBand != nullptr)
1815 : {
1816 1732 : psOptions->eWorkingDataType =
1817 1732 : GDALDataTypeUnion(psOptions->eWorkingDataType,
1818 : GDALGetRasterDataType(hDstBand));
1819 : }
1820 : }
1821 :
1822 1777 : if (psOptions->hSrcDS != nullptr)
1823 : {
1824 3514 : GDALRasterBandH hSrcBand = GDALGetRasterBand(
1825 1757 : psOptions->hSrcDS, psOptions->panSrcBands[iBand]);
1826 :
1827 1757 : if (hSrcBand != nullptr)
1828 : {
1829 1757 : psOptions->eWorkingDataType =
1830 1757 : GDALDataTypeUnion(psOptions->eWorkingDataType,
1831 : GDALGetRasterDataType(hSrcBand));
1832 : }
1833 : }
1834 :
1835 1777 : if (psOptions->padfSrcNoDataReal != nullptr)
1836 : {
1837 217 : psOptions->eWorkingDataType = GDALDataTypeUnionWithValue(
1838 : psOptions->eWorkingDataType,
1839 217 : psOptions->padfSrcNoDataReal[iBand], false);
1840 : }
1841 :
1842 1777 : if (psOptions->padfSrcNoDataImag != nullptr &&
1843 4 : psOptions->padfSrcNoDataImag[iBand] != 0.0)
1844 : {
1845 3 : psOptions->eWorkingDataType = GDALDataTypeUnionWithValue(
1846 : psOptions->eWorkingDataType,
1847 3 : psOptions->padfSrcNoDataImag[iBand], true);
1848 : }
1849 :
1850 1777 : if (psOptions->padfDstNoDataReal != nullptr)
1851 : {
1852 341 : psOptions->eWorkingDataType = GDALDataTypeUnionWithValue(
1853 : psOptions->eWorkingDataType,
1854 341 : psOptions->padfDstNoDataReal[iBand], false);
1855 : }
1856 :
1857 1777 : if (psOptions->padfDstNoDataImag != nullptr &&
1858 231 : psOptions->padfDstNoDataImag[iBand] != 0.0)
1859 : {
1860 3 : psOptions->eWorkingDataType = GDALDataTypeUnionWithValue(
1861 : psOptions->eWorkingDataType,
1862 3 : psOptions->padfDstNoDataImag[iBand], true);
1863 : }
1864 : }
1865 :
1866 1384 : if (psOptions->eWorkingDataType == GDT_Unknown)
1867 1 : psOptions->eWorkingDataType = GDT_UInt8;
1868 :
1869 2768 : const bool bApplyVerticalShift = CPLFetchBool(
1870 1384 : psOptions->papszWarpOptions, "APPLY_VERTICAL_SHIFT", false);
1871 1403 : if (bApplyVerticalShift &&
1872 19 : GDALDataTypeIsInteger(psOptions->eWorkingDataType))
1873 : {
1874 16 : const double dfMultFactorVerticalShift = CPLAtof(CSLFetchNameValueDef(
1875 16 : psOptions->papszWarpOptions, "MULT_FACTOR_VERTICAL_SHIFT", "1.0"));
1876 16 : if (dfMultFactorVerticalShift != 1)
1877 : {
1878 0 : psOptions->eWorkingDataType =
1879 0 : GDALDataTypeUnion(psOptions->eWorkingDataType, GDT_Float32);
1880 : }
1881 : }
1882 : }
1883 :
1884 : /************************************************************************/
1885 : /* GDALWarpInitDefaultBandMapping() */
1886 : /************************************************************************/
1887 :
1888 : /**
1889 : * \brief Init src and dst band mappings such that Bands[i] = i+1
1890 : * for nBandCount
1891 : * Does nothing if psOptionsIn->nBandCount is non-zero.
1892 : *
1893 : * @param psOptionsIn options to initialize.
1894 : * @param nBandCount bands to initialize for.
1895 : *
1896 : */
1897 221 : void CPL_STDCALL GDALWarpInitDefaultBandMapping(GDALWarpOptions *psOptionsIn,
1898 : int nBandCount)
1899 : {
1900 221 : if (psOptionsIn->nBandCount != 0)
1901 : {
1902 0 : return;
1903 : }
1904 :
1905 221 : psOptionsIn->nBandCount = nBandCount;
1906 :
1907 221 : psOptionsIn->panSrcBands =
1908 221 : static_cast<int *>(CPLMalloc(sizeof(int) * psOptionsIn->nBandCount));
1909 221 : psOptionsIn->panDstBands =
1910 221 : static_cast<int *>(CPLMalloc(sizeof(int) * psOptionsIn->nBandCount));
1911 :
1912 523 : for (int i = 0; i < psOptionsIn->nBandCount; i++)
1913 : {
1914 302 : psOptionsIn->panSrcBands[i] = i + 1;
1915 302 : psOptionsIn->panDstBands[i] = i + 1;
1916 : }
1917 : }
1918 :
1919 : /************************************************************************/
1920 : /* GDALSerializeWarpOptions() */
1921 : /************************************************************************/
1922 :
1923 42 : CPLXMLNode *CPL_STDCALL GDALSerializeWarpOptions(const GDALWarpOptions *psWO)
1924 :
1925 : {
1926 : /* -------------------------------------------------------------------- */
1927 : /* Create root. */
1928 : /* -------------------------------------------------------------------- */
1929 : CPLXMLNode *psTree =
1930 42 : CPLCreateXMLNode(nullptr, CXT_Element, "GDALWarpOptions");
1931 :
1932 : /* -------------------------------------------------------------------- */
1933 : /* WarpMemoryLimit */
1934 : /* -------------------------------------------------------------------- */
1935 42 : CPLCreateXMLElementAndValue(
1936 : psTree, "WarpMemoryLimit",
1937 84 : CPLString().Printf("%g", psWO->dfWarpMemoryLimit));
1938 :
1939 : /* -------------------------------------------------------------------- */
1940 : /* ResampleAlg */
1941 : /* -------------------------------------------------------------------- */
1942 42 : const char *pszAlgName = nullptr;
1943 :
1944 42 : if (psWO->eResampleAlg == GRA_NearestNeighbour)
1945 42 : pszAlgName = "NearestNeighbour";
1946 0 : else if (psWO->eResampleAlg == GRA_Bilinear)
1947 0 : pszAlgName = "Bilinear";
1948 0 : else if (psWO->eResampleAlg == GRA_Cubic)
1949 0 : pszAlgName = "Cubic";
1950 0 : else if (psWO->eResampleAlg == GRA_CubicSpline)
1951 0 : pszAlgName = "CubicSpline";
1952 0 : else if (psWO->eResampleAlg == GRA_Lanczos)
1953 0 : pszAlgName = "Lanczos";
1954 0 : else if (psWO->eResampleAlg == GRA_Average)
1955 0 : pszAlgName = "Average";
1956 0 : else if (psWO->eResampleAlg == GRA_RMS)
1957 0 : pszAlgName = "RootMeanSquare";
1958 0 : else if (psWO->eResampleAlg == GRA_Mode)
1959 0 : pszAlgName = "Mode";
1960 0 : else if (psWO->eResampleAlg == GRA_Max)
1961 0 : pszAlgName = "Maximum";
1962 0 : else if (psWO->eResampleAlg == GRA_Min)
1963 0 : pszAlgName = "Minimum";
1964 0 : else if (psWO->eResampleAlg == GRA_Med)
1965 0 : pszAlgName = "Median";
1966 0 : else if (psWO->eResampleAlg == GRA_Q1)
1967 0 : pszAlgName = "Quartile1";
1968 0 : else if (psWO->eResampleAlg == GRA_Q3)
1969 0 : pszAlgName = "Quartile3";
1970 0 : else if (psWO->eResampleAlg == GRA_Sum)
1971 0 : pszAlgName = "Sum";
1972 : else
1973 0 : pszAlgName = "Unknown";
1974 :
1975 42 : CPLCreateXMLElementAndValue(psTree, "ResampleAlg", pszAlgName);
1976 :
1977 : /* -------------------------------------------------------------------- */
1978 : /* Working Data Type */
1979 : /* -------------------------------------------------------------------- */
1980 42 : CPLCreateXMLElementAndValue(psTree, "WorkingDataType",
1981 42 : GDALGetDataTypeName(psWO->eWorkingDataType));
1982 :
1983 : /* -------------------------------------------------------------------- */
1984 : /* Name/value warp options. */
1985 : /* -------------------------------------------------------------------- */
1986 174 : for (int iWO = 0; psWO->papszWarpOptions != nullptr &&
1987 174 : psWO->papszWarpOptions[iWO] != nullptr;
1988 : iWO++)
1989 : {
1990 132 : char *pszName = nullptr;
1991 : const char *pszValue =
1992 132 : CPLParseNameValue(psWO->papszWarpOptions[iWO], &pszName);
1993 :
1994 : // EXTRA_ELTS is an internal detail that we will recover
1995 : // no need to serialize it.
1996 : // And CUTLINE is also serialized in a special way
1997 132 : if (pszName != nullptr && !EQUAL(pszName, "EXTRA_ELTS") &&
1998 90 : !EQUAL(pszName, "CUTLINE"))
1999 : {
2000 : CPLXMLNode *psOption =
2001 90 : CPLCreateXMLElementAndValue(psTree, "Option", pszValue);
2002 :
2003 90 : CPLCreateXMLNode(CPLCreateXMLNode(psOption, CXT_Attribute, "name"),
2004 : CXT_Text, pszName);
2005 : }
2006 :
2007 132 : CPLFree(pszName);
2008 : }
2009 :
2010 : /* -------------------------------------------------------------------- */
2011 : /* Source and Destination Data Source */
2012 : /* -------------------------------------------------------------------- */
2013 42 : if (psWO->hSrcDS != nullptr)
2014 : {
2015 42 : CPLCreateXMLElementAndValue(psTree, "SourceDataset",
2016 42 : GDALGetDescription(psWO->hSrcDS));
2017 :
2018 : CSLConstList papszOpenOptions =
2019 42 : GDALDataset::FromHandle(psWO->hSrcDS)->GetOpenOptions();
2020 42 : GDALSerializeOpenOptionsToXML(psTree, papszOpenOptions);
2021 : }
2022 :
2023 84 : if (psWO->hDstDS != nullptr &&
2024 42 : strlen(GDALGetDescription(psWO->hDstDS)) != 0)
2025 : {
2026 0 : CPLCreateXMLElementAndValue(psTree, "DestinationDataset",
2027 0 : GDALGetDescription(psWO->hDstDS));
2028 : }
2029 :
2030 : /* -------------------------------------------------------------------- */
2031 : /* Serialize transformer. */
2032 : /* -------------------------------------------------------------------- */
2033 42 : if (psWO->pfnTransformer != nullptr)
2034 : {
2035 : CPLXMLNode *psTransformerContainer =
2036 42 : CPLCreateXMLNode(psTree, CXT_Element, "Transformer");
2037 :
2038 84 : CPLXMLNode *psTransformerTree = GDALSerializeTransformer(
2039 42 : psWO->pfnTransformer, psWO->pTransformerArg);
2040 :
2041 42 : if (psTransformerTree != nullptr)
2042 42 : CPLAddXMLChild(psTransformerContainer, psTransformerTree);
2043 : }
2044 :
2045 : /* -------------------------------------------------------------------- */
2046 : /* Band count and lists. */
2047 : /* -------------------------------------------------------------------- */
2048 42 : CPLXMLNode *psBandList = nullptr;
2049 :
2050 42 : if (psWO->nBandCount != 0)
2051 42 : psBandList = CPLCreateXMLNode(psTree, CXT_Element, "BandList");
2052 :
2053 107 : for (int i = 0; i < psWO->nBandCount; i++)
2054 : {
2055 : CPLXMLNode *psBand;
2056 :
2057 65 : psBand = CPLCreateXMLNode(psBandList, CXT_Element, "BandMapping");
2058 65 : if (psWO->panSrcBands != nullptr)
2059 65 : CPLCreateXMLNode(CPLCreateXMLNode(psBand, CXT_Attribute, "src"),
2060 : CXT_Text,
2061 130 : CPLString().Printf("%d", psWO->panSrcBands[i]));
2062 65 : if (psWO->panDstBands != nullptr)
2063 65 : CPLCreateXMLNode(CPLCreateXMLNode(psBand, CXT_Attribute, "dst"),
2064 : CXT_Text,
2065 130 : CPLString().Printf("%d", psWO->panDstBands[i]));
2066 :
2067 65 : if (psWO->padfSrcNoDataReal != nullptr)
2068 : {
2069 13 : CPLCreateXMLElementAndValue(
2070 : psBand, "SrcNoDataReal",
2071 13 : VRTSerializeNoData(psWO->padfSrcNoDataReal[i],
2072 13 : psWO->eWorkingDataType, 16)
2073 : .c_str());
2074 : }
2075 :
2076 65 : if (psWO->padfSrcNoDataImag != nullptr)
2077 : {
2078 2 : if (std::isnan(psWO->padfSrcNoDataImag[i]))
2079 0 : CPLCreateXMLElementAndValue(psBand, "SrcNoDataImag", "nan");
2080 : else
2081 2 : CPLCreateXMLElementAndValue(
2082 : psBand, "SrcNoDataImag",
2083 4 : CPLString().Printf("%.16g", psWO->padfSrcNoDataImag[i]));
2084 : }
2085 : // Compatibility with GDAL <= 2.2: if we serialize a SrcNoDataReal,
2086 : // it needs a SrcNoDataImag as well
2087 63 : else if (psWO->padfSrcNoDataReal != nullptr)
2088 : {
2089 11 : CPLCreateXMLElementAndValue(psBand, "SrcNoDataImag", "0");
2090 : }
2091 :
2092 65 : if (psWO->padfDstNoDataReal != nullptr)
2093 : {
2094 13 : CPLCreateXMLElementAndValue(
2095 : psBand, "DstNoDataReal",
2096 13 : VRTSerializeNoData(psWO->padfDstNoDataReal[i],
2097 13 : psWO->eWorkingDataType, 16)
2098 : .c_str());
2099 : }
2100 :
2101 65 : if (psWO->padfDstNoDataImag != nullptr)
2102 : {
2103 2 : if (std::isnan(psWO->padfDstNoDataImag[i]))
2104 0 : CPLCreateXMLElementAndValue(psBand, "DstNoDataImag", "nan");
2105 : else
2106 2 : CPLCreateXMLElementAndValue(
2107 : psBand, "DstNoDataImag",
2108 4 : CPLString().Printf("%.16g", psWO->padfDstNoDataImag[i]));
2109 : }
2110 : // Compatibility with GDAL <= 2.2: if we serialize a DstNoDataReal,
2111 : // it needs a SrcNoDataImag as well
2112 63 : else if (psWO->padfDstNoDataReal != nullptr)
2113 : {
2114 11 : CPLCreateXMLElementAndValue(psBand, "DstNoDataImag", "0");
2115 : }
2116 : }
2117 :
2118 : /* -------------------------------------------------------------------- */
2119 : /* Alpha bands. */
2120 : /* -------------------------------------------------------------------- */
2121 42 : if (psWO->nSrcAlphaBand > 0)
2122 0 : CPLCreateXMLElementAndValue(
2123 : psTree, "SrcAlphaBand",
2124 0 : CPLString().Printf("%d", psWO->nSrcAlphaBand));
2125 :
2126 42 : if (psWO->nDstAlphaBand > 0)
2127 6 : CPLCreateXMLElementAndValue(
2128 : psTree, "DstAlphaBand",
2129 12 : CPLString().Printf("%d", psWO->nDstAlphaBand));
2130 :
2131 : /* -------------------------------------------------------------------- */
2132 : /* Cutline. */
2133 : /* -------------------------------------------------------------------- */
2134 42 : if (psWO->hCutline != nullptr)
2135 : {
2136 0 : char *pszWKT = nullptr;
2137 0 : if (OGR_G_ExportToWkt(static_cast<OGRGeometryH>(psWO->hCutline),
2138 0 : &pszWKT) == OGRERR_NONE)
2139 : {
2140 0 : CPLCreateXMLElementAndValue(psTree, "Cutline", pszWKT);
2141 : }
2142 0 : CPLFree(pszWKT);
2143 : }
2144 :
2145 42 : if (psWO->dfCutlineBlendDist != 0.0)
2146 0 : CPLCreateXMLElementAndValue(
2147 : psTree, "CutlineBlendDist",
2148 0 : CPLString().Printf("%.5g", psWO->dfCutlineBlendDist));
2149 :
2150 42 : return psTree;
2151 : }
2152 :
2153 : /************************************************************************/
2154 : /* GDALDeserializeWarpOptions() */
2155 : /************************************************************************/
2156 :
2157 147 : GDALWarpOptions *CPL_STDCALL GDALDeserializeWarpOptions(CPLXMLNode *psTree)
2158 :
2159 : {
2160 147 : CPLErrorReset();
2161 :
2162 : /* -------------------------------------------------------------------- */
2163 : /* Verify this is the right kind of object. */
2164 : /* -------------------------------------------------------------------- */
2165 147 : if (psTree == nullptr || psTree->eType != CXT_Element ||
2166 147 : !EQUAL(psTree->pszValue, "GDALWarpOptions"))
2167 : {
2168 0 : CPLError(CE_Failure, CPLE_AppDefined,
2169 : "Wrong node, unable to deserialize GDALWarpOptions.");
2170 0 : return nullptr;
2171 : }
2172 :
2173 : /* -------------------------------------------------------------------- */
2174 : /* Create pre-initialized warp options. */
2175 : /* -------------------------------------------------------------------- */
2176 147 : GDALWarpOptions *psWO = GDALCreateWarpOptions();
2177 :
2178 : /* -------------------------------------------------------------------- */
2179 : /* Warp memory limit. */
2180 : /* -------------------------------------------------------------------- */
2181 147 : psWO->dfWarpMemoryLimit =
2182 147 : CPLAtof(CPLGetXMLValue(psTree, "WarpMemoryLimit", "0.0"));
2183 :
2184 : /* -------------------------------------------------------------------- */
2185 : /* resample algorithm */
2186 : /* -------------------------------------------------------------------- */
2187 147 : const char *pszValue = CPLGetXMLValue(psTree, "ResampleAlg", "Default");
2188 :
2189 147 : if (EQUAL(pszValue, "NearestNeighbour"))
2190 89 : psWO->eResampleAlg = GRA_NearestNeighbour;
2191 58 : else if (EQUAL(pszValue, "Bilinear"))
2192 8 : psWO->eResampleAlg = GRA_Bilinear;
2193 50 : else if (EQUAL(pszValue, "Cubic"))
2194 9 : psWO->eResampleAlg = GRA_Cubic;
2195 41 : else if (EQUAL(pszValue, "CubicSpline"))
2196 9 : psWO->eResampleAlg = GRA_CubicSpline;
2197 32 : else if (EQUAL(pszValue, "Lanczos"))
2198 4 : psWO->eResampleAlg = GRA_Lanczos;
2199 28 : else if (EQUAL(pszValue, "Average"))
2200 6 : psWO->eResampleAlg = GRA_Average;
2201 22 : else if (EQUAL(pszValue, "RootMeanSquare"))
2202 5 : psWO->eResampleAlg = GRA_RMS;
2203 17 : else if (EQUAL(pszValue, "Mode"))
2204 4 : psWO->eResampleAlg = GRA_Mode;
2205 13 : else if (EQUAL(pszValue, "Maximum"))
2206 3 : psWO->eResampleAlg = GRA_Max;
2207 10 : else if (EQUAL(pszValue, "Minimum"))
2208 2 : psWO->eResampleAlg = GRA_Min;
2209 8 : else if (EQUAL(pszValue, "Median"))
2210 3 : psWO->eResampleAlg = GRA_Med;
2211 5 : else if (EQUAL(pszValue, "Quartile1"))
2212 2 : psWO->eResampleAlg = GRA_Q1;
2213 3 : else if (EQUAL(pszValue, "Quartile3"))
2214 2 : psWO->eResampleAlg = GRA_Q3;
2215 1 : else if (EQUAL(pszValue, "Sum"))
2216 1 : psWO->eResampleAlg = GRA_Sum;
2217 0 : else if (EQUAL(pszValue, "Default"))
2218 : /* leave as is */;
2219 : else
2220 : {
2221 0 : CPLError(CE_Failure, CPLE_AppDefined,
2222 : "Unrecognised ResampleAlg value '%s'.", pszValue);
2223 : }
2224 :
2225 : /* -------------------------------------------------------------------- */
2226 : /* Working data type. */
2227 : /* -------------------------------------------------------------------- */
2228 147 : psWO->eWorkingDataType = GDALGetDataTypeByName(
2229 : CPLGetXMLValue(psTree, "WorkingDataType", "Unknown"));
2230 :
2231 : /* -------------------------------------------------------------------- */
2232 : /* Name/value warp options. */
2233 : /* -------------------------------------------------------------------- */
2234 1324 : for (CPLXMLNode *psItem = psTree->psChild; psItem != nullptr;
2235 1177 : psItem = psItem->psNext)
2236 : {
2237 1177 : if (psItem->eType == CXT_Element && EQUAL(psItem->pszValue, "Option"))
2238 : {
2239 227 : const char *pszName = CPLGetXMLValue(psItem, "Name", nullptr);
2240 227 : pszValue = CPLGetXMLValue(psItem, "", nullptr);
2241 :
2242 227 : if (pszName != nullptr && pszValue != nullptr)
2243 : {
2244 227 : psWO->papszWarpOptions =
2245 227 : CSLSetNameValue(psWO->papszWarpOptions, pszName, pszValue);
2246 : }
2247 : }
2248 : }
2249 :
2250 : /* -------------------------------------------------------------------- */
2251 : /* Source Dataset. */
2252 : /* -------------------------------------------------------------------- */
2253 147 : pszValue = CPLGetXMLValue(psTree, "SourceDataset", nullptr);
2254 :
2255 147 : if (pszValue != nullptr)
2256 : {
2257 : CPLXMLNode *psGeoLocNode =
2258 147 : CPLSearchXMLNode(psTree, "GeoLocTransformer");
2259 147 : if (psGeoLocNode)
2260 : {
2261 1 : CPLCreateXMLElementAndValue(psGeoLocNode, "SourceDataset",
2262 : pszValue);
2263 : }
2264 :
2265 294 : CPLConfigOptionSetter oSetter("CPL_ALLOW_VSISTDIN", "NO", true);
2266 :
2267 147 : char **papszOpenOptions = GDALDeserializeOpenOptionsFromXML(psTree);
2268 147 : psWO->hSrcDS =
2269 147 : GDALOpenEx(pszValue, GDAL_OF_RASTER | GDAL_OF_VERBOSE_ERROR,
2270 : nullptr, papszOpenOptions, nullptr);
2271 147 : CSLDestroy(papszOpenOptions);
2272 : }
2273 :
2274 : /* -------------------------------------------------------------------- */
2275 : /* Destination Dataset. */
2276 : /* -------------------------------------------------------------------- */
2277 147 : pszValue = CPLGetXMLValue(psTree, "DestinationDataset", nullptr);
2278 :
2279 147 : if (pszValue != nullptr)
2280 : {
2281 0 : psWO->hDstDS = GDALOpenShared(pszValue, GA_Update);
2282 : }
2283 :
2284 : /* -------------------------------------------------------------------- */
2285 : /* First, count band mappings so we can establish the bandcount. */
2286 : /* -------------------------------------------------------------------- */
2287 147 : CPLXMLNode *psBandTree = CPLGetXMLNode(psTree, "BandList");
2288 :
2289 147 : int nBandCount = 0;
2290 147 : CPLXMLNode *psBand = psBandTree ? psBandTree->psChild : nullptr;
2291 369 : for (; psBand != nullptr; psBand = psBand->psNext)
2292 : {
2293 222 : if (psBand->eType != CXT_Element ||
2294 222 : !EQUAL(psBand->pszValue, "BandMapping"))
2295 0 : continue;
2296 :
2297 222 : nBandCount++;
2298 : }
2299 :
2300 147 : GDALWarpInitDefaultBandMapping(psWO, nBandCount);
2301 :
2302 : /* ==================================================================== */
2303 : /* Now actually process each bandmapping. */
2304 : /* ==================================================================== */
2305 147 : int iBand = 0;
2306 :
2307 147 : psBand = psBandTree ? psBandTree->psChild : nullptr;
2308 :
2309 369 : for (; psBand != nullptr; psBand = psBand->psNext)
2310 : {
2311 222 : if (psBand->eType != CXT_Element ||
2312 222 : !EQUAL(psBand->pszValue, "BandMapping"))
2313 0 : continue;
2314 :
2315 : /* --------------------------------------------------------------------
2316 : */
2317 : /* Source band */
2318 : /* --------------------------------------------------------------------
2319 : */
2320 222 : pszValue = CPLGetXMLValue(psBand, "src", nullptr);
2321 222 : if (pszValue != nullptr)
2322 222 : psWO->panSrcBands[iBand] = atoi(pszValue);
2323 :
2324 : /* --------------------------------------------------------------------
2325 : */
2326 : /* Destination band. */
2327 : /* --------------------------------------------------------------------
2328 : */
2329 222 : pszValue = CPLGetXMLValue(psBand, "dst", nullptr);
2330 222 : if (pszValue != nullptr)
2331 222 : psWO->panDstBands[iBand] = atoi(pszValue);
2332 :
2333 66 : const auto NormalizeValue = [](const char *pszValueIn,
2334 : GDALDataType eDataType) -> double
2335 : {
2336 78 : if (eDataType == GDT_Float32 &&
2337 78 : CPLString().Printf("%.16g",
2338 : static_cast<double>(
2339 12 : std::numeric_limits<float>::lowest())) ==
2340 : pszValueIn)
2341 : {
2342 : return static_cast<double>(
2343 0 : std::numeric_limits<float>::lowest());
2344 : }
2345 78 : else if (eDataType == GDT_Float32 &&
2346 78 : CPLString().Printf(
2347 : "%.16g", static_cast<double>(
2348 12 : std::numeric_limits<float>::max())) ==
2349 : pszValueIn)
2350 : {
2351 0 : return static_cast<double>(std::numeric_limits<float>::max());
2352 : }
2353 : else
2354 : {
2355 66 : return CPLAtof(pszValueIn);
2356 : }
2357 : };
2358 :
2359 : /* --------------------------------------------------------------------
2360 : */
2361 : /* Source nodata. */
2362 : /* --------------------------------------------------------------------
2363 : */
2364 222 : pszValue = CPLGetXMLValue(psBand, "SrcNoDataReal", nullptr);
2365 222 : if (pszValue != nullptr)
2366 : {
2367 31 : GDALWarpInitSrcNoDataReal(psWO, -1.1e20);
2368 62 : psWO->padfSrcNoDataReal[iBand] =
2369 31 : NormalizeValue(pszValue, psWO->eWorkingDataType);
2370 : }
2371 :
2372 222 : pszValue = CPLGetXMLValue(psBand, "SrcNoDataImag", nullptr);
2373 222 : if (pszValue != nullptr)
2374 : {
2375 31 : GDALWarpInitSrcNoDataImag(psWO, 0);
2376 31 : psWO->padfSrcNoDataImag[iBand] = CPLAtof(pszValue);
2377 : }
2378 :
2379 : /* --------------------------------------------------------------------
2380 : */
2381 : /* Destination nodata. */
2382 : /* --------------------------------------------------------------------
2383 : */
2384 222 : pszValue = CPLGetXMLValue(psBand, "DstNoDataReal", nullptr);
2385 222 : if (pszValue != nullptr)
2386 : {
2387 35 : GDALWarpInitDstNoDataReal(psWO, -1.1e20);
2388 70 : psWO->padfDstNoDataReal[iBand] =
2389 35 : NormalizeValue(pszValue, psWO->eWorkingDataType);
2390 : }
2391 :
2392 222 : pszValue = CPLGetXMLValue(psBand, "DstNoDataImag", nullptr);
2393 222 : if (pszValue != nullptr)
2394 : {
2395 35 : GDALWarpInitDstNoDataImag(psWO, 0);
2396 35 : psWO->padfDstNoDataImag[iBand] = CPLAtof(pszValue);
2397 : }
2398 :
2399 222 : iBand++;
2400 : }
2401 :
2402 : /* -------------------------------------------------------------------- */
2403 : /* Alpha bands. */
2404 : /* -------------------------------------------------------------------- */
2405 147 : psWO->nSrcAlphaBand = atoi(CPLGetXMLValue(psTree, "SrcAlphaBand", "0"));
2406 147 : psWO->nDstAlphaBand = atoi(CPLGetXMLValue(psTree, "DstAlphaBand", "0"));
2407 :
2408 : /* -------------------------------------------------------------------- */
2409 : /* Cutline. */
2410 : /* -------------------------------------------------------------------- */
2411 147 : const char *pszWKT = CPLGetXMLValue(psTree, "Cutline", nullptr);
2412 147 : if (pszWKT)
2413 : {
2414 7 : char *pszWKTTemp = const_cast<char *>(pszWKT);
2415 7 : OGRGeometryH hCutline = nullptr;
2416 7 : OGR_G_CreateFromWkt(&pszWKTTemp, nullptr, &hCutline);
2417 7 : psWO->hCutline = hCutline;
2418 : }
2419 :
2420 147 : psWO->dfCutlineBlendDist =
2421 147 : CPLAtof(CPLGetXMLValue(psTree, "CutlineBlendDist", "0"));
2422 :
2423 : /* -------------------------------------------------------------------- */
2424 : /* Transformation. */
2425 : /* -------------------------------------------------------------------- */
2426 147 : CPLXMLNode *psTransformer = CPLGetXMLNode(psTree, "Transformer");
2427 :
2428 147 : if (psTransformer != nullptr && psTransformer->psChild != nullptr)
2429 : {
2430 147 : GDALDeserializeTransformer(psTransformer->psChild,
2431 : &(psWO->pfnTransformer),
2432 : &(psWO->pTransformerArg));
2433 : }
2434 :
2435 : /* -------------------------------------------------------------------- */
2436 : /* If any error has occurred, cleanup else return success. */
2437 : /* -------------------------------------------------------------------- */
2438 147 : if (CPLGetLastErrorType() != CE_None)
2439 : {
2440 0 : if (psWO->pTransformerArg)
2441 : {
2442 0 : GDALDestroyTransformer(psWO->pTransformerArg);
2443 0 : psWO->pTransformerArg = nullptr;
2444 : }
2445 0 : if (psWO->hSrcDS != nullptr)
2446 : {
2447 0 : GDALClose(psWO->hSrcDS);
2448 0 : psWO->hSrcDS = nullptr;
2449 : }
2450 0 : if (psWO->hDstDS != nullptr)
2451 : {
2452 0 : GDALClose(psWO->hDstDS);
2453 0 : psWO->hDstDS = nullptr;
2454 : }
2455 0 : GDALDestroyWarpOptions(psWO);
2456 0 : return nullptr;
2457 : }
2458 :
2459 147 : return psWO;
2460 : }
2461 :
2462 : /************************************************************************/
2463 : /* GDALGetWarpResampleAlg() */
2464 : /************************************************************************/
2465 :
2466 : /** Return a GDALResampleAlg from a string */
2467 871 : bool GDALGetWarpResampleAlg(const char *pszResampling,
2468 : GDALResampleAlg &eResampleAlg, bool bThrow)
2469 : {
2470 871 : if (STARTS_WITH_CI(pszResampling, "near"))
2471 183 : eResampleAlg = GRA_NearestNeighbour;
2472 688 : else if (EQUAL(pszResampling, "bilinear"))
2473 149 : eResampleAlg = GRA_Bilinear;
2474 539 : else if (EQUAL(pszResampling, "cubic"))
2475 255 : eResampleAlg = GRA_Cubic;
2476 284 : else if (EQUAL(pszResampling, "cubicspline"))
2477 54 : eResampleAlg = GRA_CubicSpline;
2478 230 : else if (EQUAL(pszResampling, "lanczos"))
2479 52 : eResampleAlg = GRA_Lanczos;
2480 178 : else if (EQUAL(pszResampling, "average"))
2481 92 : eResampleAlg = GRA_Average;
2482 86 : else if (EQUAL(pszResampling, "rms"))
2483 3 : eResampleAlg = GRA_RMS;
2484 83 : else if (EQUAL(pszResampling, "mode"))
2485 42 : eResampleAlg = GRA_Mode;
2486 41 : else if (EQUAL(pszResampling, "max"))
2487 3 : eResampleAlg = GRA_Max;
2488 38 : else if (EQUAL(pszResampling, "min"))
2489 3 : eResampleAlg = GRA_Min;
2490 35 : else if (EQUAL(pszResampling, "med"))
2491 3 : eResampleAlg = GRA_Med;
2492 32 : else if (EQUAL(pszResampling, "q1"))
2493 8 : eResampleAlg = GRA_Q1;
2494 24 : else if (EQUAL(pszResampling, "q3"))
2495 3 : eResampleAlg = GRA_Q3;
2496 21 : else if (EQUAL(pszResampling, "sum"))
2497 20 : eResampleAlg = GRA_Sum;
2498 : else
2499 : {
2500 1 : if (bThrow)
2501 : {
2502 1 : throw std::invalid_argument("Unknown resampling method");
2503 : }
2504 : else
2505 : {
2506 0 : CPLError(CE_Failure, CPLE_IllegalArg,
2507 : "Unknown resampling method: %s.", pszResampling);
2508 0 : return false;
2509 : }
2510 : }
2511 870 : return true;
2512 : }
|