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