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