Coverage Report

Created: 2025-06-22 06:59

/src/gdal/frmts/vrt/vrtsources.cpp
Line
Count
Source (jump to first uncovered line)
1
/******************************************************************************
2
 *
3
 * Project:  Virtual GDAL Datasets
4
 * Purpose:  Implementation of VRTSimpleSource, VRTFuncSource and
5
 *           VRTAveragedSource.
6
 * Author:   Frank Warmerdam <warmerdam@pobox.com>
7
 *
8
 ******************************************************************************
9
 * Copyright (c) 2001, Frank Warmerdam <warmerdam@pobox.com>
10
 * Copyright (c) 2008-2013, Even Rouault <even dot rouault at spatialys.com>
11
 *
12
 * SPDX-License-Identifier: MIT
13
 ****************************************************************************/
14
15
#include "gdal_vrt.h"
16
#include "vrtdataset.h"
17
18
#include <cassert>
19
#include <climits>
20
#include <cmath>
21
#include <cstddef>
22
#include <cstdio>
23
#include <cstdlib>
24
#include <cstring>
25
#include <algorithm>
26
#include <limits>
27
#include <string>
28
29
#include "cpl_conv.h"
30
#include "cpl_error.h"
31
#include "cpl_hash_set.h"
32
#include "cpl_minixml.h"
33
#include "cpl_progress.h"
34
#include "cpl_string.h"
35
#include "cpl_vsi.h"
36
#include "gdal.h"
37
#include "gdal_priv.h"
38
#include "gdal_proxy.h"
39
#include "gdal_priv_templates.hpp"
40
41
/*! @cond Doxygen_Suppress */
42
43
// #define DEBUG_VERBOSE 1
44
45
/************************************************************************/
46
/* ==================================================================== */
47
/*                             VRTSource                                */
48
/* ==================================================================== */
49
/************************************************************************/
50
51
VRTSource::~VRTSource()
52
0
{
53
0
}
54
55
/************************************************************************/
56
/*                             GetFileList()                            */
57
/************************************************************************/
58
59
void VRTSource::GetFileList(char *** /* ppapszFileList */, int * /* pnSize */,
60
                            int * /* pnMaxSize */, CPLHashSet * /* hSetFiles */)
61
0
{
62
0
}
63
64
/************************************************************************/
65
/* ==================================================================== */
66
/*                          VRTSimpleSource                             */
67
/* ==================================================================== */
68
/************************************************************************/
69
70
/************************************************************************/
71
/*                          VRTSimpleSource()                           */
72
/************************************************************************/
73
74
0
VRTSimpleSource::VRTSimpleSource() = default;
75
76
/************************************************************************/
77
/*                          VRTSimpleSource()                           */
78
/************************************************************************/
79
80
VRTSimpleSource::VRTSimpleSource(const VRTSimpleSource *poSrcSource,
81
                                 double dfXDstRatio, double dfYDstRatio)
82
0
    : m_poMapSharedSources(poSrcSource->m_poMapSharedSources),
83
0
      m_poRasterBand(poSrcSource->m_poRasterBand),
84
0
      m_poMaskBandMainBand(poSrcSource->m_poMaskBandMainBand),
85
0
      m_aosOpenOptionsOri(poSrcSource->m_aosOpenOptionsOri),
86
0
      m_aosOpenOptions(poSrcSource->m_aosOpenOptions),
87
0
      m_bSrcDSNameFromVRT(poSrcSource->m_bSrcDSNameFromVRT),
88
0
      m_nBand(poSrcSource->m_nBand),
89
0
      m_bGetMaskBand(poSrcSource->m_bGetMaskBand),
90
0
      m_dfSrcXOff(poSrcSource->m_dfSrcXOff),
91
0
      m_dfSrcYOff(poSrcSource->m_dfSrcYOff),
92
0
      m_dfSrcXSize(poSrcSource->m_dfSrcXSize),
93
0
      m_dfSrcYSize(poSrcSource->m_dfSrcYSize),
94
0
      m_nMaxValue(poSrcSource->m_nMaxValue), m_bRelativeToVRTOri(-1),
95
0
      m_nExplicitSharedStatus(poSrcSource->m_nExplicitSharedStatus),
96
0
      m_osSrcDSName(poSrcSource->m_osSrcDSName),
97
0
      m_bDropRefOnSrcBand(poSrcSource->m_bDropRefOnSrcBand)
98
0
{
99
0
    if (!poSrcSource->IsSrcWinSet() && !poSrcSource->IsDstWinSet() &&
100
0
        (dfXDstRatio != 1.0 || dfYDstRatio != 1.0))
101
0
    {
102
0
        auto l_band = GetRasterBand();
103
0
        if (l_band)
104
0
        {
105
0
            m_dfSrcXOff = 0;
106
0
            m_dfSrcYOff = 0;
107
0
            m_dfSrcXSize = l_band->GetXSize();
108
0
            m_dfSrcYSize = l_band->GetYSize();
109
0
            m_dfDstXOff = 0;
110
0
            m_dfDstYOff = 0;
111
0
            m_dfDstXSize = l_band->GetXSize() * dfXDstRatio;
112
0
            m_dfDstYSize = l_band->GetYSize() * dfYDstRatio;
113
0
        }
114
0
    }
115
0
    else if (poSrcSource->IsDstWinSet())
116
0
    {
117
0
        m_dfDstXOff = poSrcSource->m_dfDstXOff * dfXDstRatio;
118
0
        m_dfDstYOff = poSrcSource->m_dfDstYOff * dfYDstRatio;
119
0
        m_dfDstXSize = poSrcSource->m_dfDstXSize * dfXDstRatio;
120
0
        m_dfDstYSize = poSrcSource->m_dfDstYSize * dfYDstRatio;
121
0
    }
122
0
}
123
124
/************************************************************************/
125
/*                          ~VRTSimpleSource()                          */
126
/************************************************************************/
127
128
VRTSimpleSource::~VRTSimpleSource()
129
130
0
{
131
0
    if (!m_bDropRefOnSrcBand)
132
0
        return;
133
134
0
    if (m_poMaskBandMainBand != nullptr)
135
0
    {
136
0
        if (m_poMaskBandMainBand->GetDataset() != nullptr)
137
0
        {
138
0
            m_poMaskBandMainBand->GetDataset()->ReleaseRef();
139
0
        }
140
0
    }
141
0
    else if (m_poRasterBand != nullptr &&
142
0
             m_poRasterBand->GetDataset() != nullptr)
143
0
    {
144
0
        m_poRasterBand->GetDataset()->ReleaseRef();
145
0
    }
146
0
}
147
148
/************************************************************************/
149
/*                           GetTypeStatic()                            */
150
/************************************************************************/
151
152
const char *VRTSimpleSource::GetTypeStatic()
153
0
{
154
0
    static const char *TYPE = "SimpleSource";
155
0
    return TYPE;
156
0
}
157
158
/************************************************************************/
159
/*                            GetType()                                 */
160
/************************************************************************/
161
162
const char *VRTSimpleSource::GetType() const
163
0
{
164
0
    return GetTypeStatic();
165
0
}
166
167
/************************************************************************/
168
/*                           FlushCache()                               */
169
/************************************************************************/
170
171
CPLErr VRTSimpleSource::FlushCache(bool bAtClosing)
172
173
0
{
174
0
    if (m_poMaskBandMainBand != nullptr)
175
0
    {
176
0
        return m_poMaskBandMainBand->FlushCache(bAtClosing);
177
0
    }
178
0
    else if (m_poRasterBand != nullptr)
179
0
    {
180
0
        return m_poRasterBand->FlushCache(bAtClosing);
181
0
    }
182
0
    return CE_None;
183
0
}
184
185
/************************************************************************/
186
/*                    UnsetPreservedRelativeFilenames()                 */
187
/************************************************************************/
188
189
void VRTSimpleSource::UnsetPreservedRelativeFilenames()
190
0
{
191
0
    if (m_bRelativeToVRTOri &&
192
0
        !STARTS_WITH(m_osSourceFileNameOri.c_str(), "http://") &&
193
0
        !STARTS_WITH(m_osSourceFileNameOri.c_str(), "https://"))
194
0
    {
195
0
        m_bRelativeToVRTOri = -1;
196
0
        m_osSourceFileNameOri = "";
197
0
    }
198
0
}
199
200
/************************************************************************/
201
/*                             SetSrcBand()                             */
202
/************************************************************************/
203
204
void VRTSimpleSource::SetSrcBand(const char *pszFilename, int nBand)
205
206
0
{
207
0
    m_nBand = nBand;
208
0
    m_osSrcDSName = pszFilename;
209
0
}
210
211
/************************************************************************/
212
/*                             SetSrcBand()                             */
213
/************************************************************************/
214
215
void VRTSimpleSource::SetSrcBand(GDALRasterBand *poNewSrcBand)
216
217
0
{
218
0
    m_poRasterBand = poNewSrcBand;
219
0
    m_nBand = m_poRasterBand->GetBand();
220
0
    auto poDS = poNewSrcBand->GetDataset();
221
0
    if (poDS != nullptr)
222
0
    {
223
0
        m_osSrcDSName = poDS->GetDescription();
224
0
        m_aosOpenOptions = CSLDuplicate(poDS->GetOpenOptions());
225
0
        m_aosOpenOptionsOri = m_aosOpenOptions;
226
0
    }
227
0
}
228
229
/************************************************************************/
230
/*                      SetSourceDatasetName()                          */
231
/************************************************************************/
232
233
void VRTSimpleSource::SetSourceDatasetName(const char *pszFilename,
234
                                           bool bRelativeToVRT)
235
0
{
236
0
    CPLAssert(m_nBand >= 0);
237
0
    m_osSrcDSName = pszFilename;
238
0
    m_osSourceFileNameOri = pszFilename;
239
0
    m_bRelativeToVRTOri = bRelativeToVRT;
240
0
}
241
242
/************************************************************************/
243
/*                          SetSrcMaskBand()                            */
244
/************************************************************************/
245
246
// poSrcBand is not the mask band, but the band from which the mask band is
247
// taken.
248
void VRTSimpleSource::SetSrcMaskBand(GDALRasterBand *poNewSrcBand)
249
250
0
{
251
0
    m_poRasterBand = poNewSrcBand->GetMaskBand();
252
0
    m_poMaskBandMainBand = poNewSrcBand;
253
0
    m_nBand = poNewSrcBand->GetBand();
254
0
    auto poDS = poNewSrcBand->GetDataset();
255
0
    if (poDS != nullptr)
256
0
    {
257
0
        m_osSrcDSName = poDS->GetDescription();
258
0
        m_aosOpenOptions = CSLDuplicate(poDS->GetOpenOptions());
259
0
        m_aosOpenOptionsOri = m_aosOpenOptions;
260
0
    }
261
0
    m_bGetMaskBand = true;
262
0
}
263
264
/************************************************************************/
265
/*                         RoundIfCloseToInt()                          */
266
/************************************************************************/
267
268
static double RoundIfCloseToInt(double dfValue)
269
0
{
270
0
    double dfClosestInt = floor(dfValue + 0.5);
271
0
    return (fabs(dfValue - dfClosestInt) < 1e-3) ? dfClosestInt : dfValue;
272
0
}
273
274
/************************************************************************/
275
/*                            SetSrcWindow()                            */
276
/************************************************************************/
277
278
void VRTSimpleSource::SetSrcWindow(double dfNewXOff, double dfNewYOff,
279
                                   double dfNewXSize, double dfNewYSize)
280
281
0
{
282
0
    m_dfSrcXOff = RoundIfCloseToInt(dfNewXOff);
283
0
    m_dfSrcYOff = RoundIfCloseToInt(dfNewYOff);
284
0
    m_dfSrcXSize = RoundIfCloseToInt(dfNewXSize);
285
0
    m_dfSrcYSize = RoundIfCloseToInt(dfNewYSize);
286
0
}
287
288
/************************************************************************/
289
/*                            SetDstWindow()                            */
290
/************************************************************************/
291
292
void VRTSimpleSource::SetDstWindow(double dfNewXOff, double dfNewYOff,
293
                                   double dfNewXSize, double dfNewYSize)
294
295
0
{
296
0
    m_dfDstXOff = RoundIfCloseToInt(dfNewXOff);
297
0
    m_dfDstYOff = RoundIfCloseToInt(dfNewYOff);
298
0
    m_dfDstXSize = RoundIfCloseToInt(dfNewXSize);
299
0
    m_dfDstYSize = RoundIfCloseToInt(dfNewYSize);
300
0
}
301
302
/************************************************************************/
303
/*                            GetDstWindow()                            */
304
/************************************************************************/
305
306
void VRTSimpleSource::GetDstWindow(double &dfDstXOff, double &dfDstYOff,
307
                                   double &dfDstXSize, double &dfDstYSize) const
308
0
{
309
0
    dfDstXOff = m_dfDstXOff;
310
0
    dfDstYOff = m_dfDstYOff;
311
0
    dfDstXSize = m_dfDstXSize;
312
0
    dfDstYSize = m_dfDstYSize;
313
0
}
314
315
/************************************************************************/
316
/*                        DstWindowIntersects()                         */
317
/************************************************************************/
318
319
bool VRTSimpleSource::DstWindowIntersects(double dfXOff, double dfYOff,
320
                                          double dfXSize, double dfYSize) const
321
0
{
322
0
    return IsDstWinSet() && m_dfDstXOff + m_dfDstXSize > dfXOff &&
323
0
           m_dfDstYOff + m_dfDstYSize > dfYOff &&
324
0
           m_dfDstXOff < dfXOff + dfXSize && m_dfDstYOff < dfYOff + dfYSize;
325
0
}
326
327
/************************************************************************/
328
/*                            IsSlowSource()                            */
329
/************************************************************************/
330
331
static bool IsSlowSource(const char *pszSrcName)
332
0
{
333
0
    return strstr(pszSrcName, "/vsicurl/http") != nullptr ||
334
0
           strstr(pszSrcName, "/vsicurl/ftp") != nullptr ||
335
0
           (strstr(pszSrcName, "/vsicurl?") != nullptr &&
336
0
            strstr(pszSrcName, "&url=http") != nullptr);
337
0
}
338
339
/************************************************************************/
340
/*                         AddSourceFilenameNode()                      */
341
/************************************************************************/
342
343
void VRTSimpleSource::AddSourceFilenameNode(const char *pszVRTPath,
344
                                            CPLXMLNode *psSrc)
345
0
{
346
347
0
    VSIStatBufL sStat;
348
0
    int bRelativeToVRT = FALSE;  // TODO(schwehr): Make this a bool?
349
0
    std::string osSourceFilename;
350
351
0
    if (m_bRelativeToVRTOri >= 0)
352
0
    {
353
0
        osSourceFilename = m_osSourceFileNameOri;
354
0
        bRelativeToVRT = m_bRelativeToVRTOri;
355
0
    }
356
0
    else if (IsSlowSource(m_osSrcDSName))
357
0
    {
358
        // Testing the existence of remote resources can be excruciating
359
        // slow, so let's just suppose they exist.
360
0
        osSourceFilename = m_osSrcDSName;
361
0
        bRelativeToVRT = FALSE;
362
0
    }
363
    // If this isn't actually a file, don't even try to know if it is a
364
    // relative path. It can't be !, and unfortunately CPLIsFilenameRelative()
365
    // can only work with strings that are filenames To be clear
366
    // NITF_TOC_ENTRY:CADRG_JOG-A_250K_1_0:some_path isn't a relative file
367
    // path.
368
0
    else if (VSIStatExL(m_osSrcDSName, &sStat, VSI_STAT_EXISTS_FLAG) != 0)
369
0
    {
370
0
        osSourceFilename = m_osSrcDSName;
371
0
        bRelativeToVRT = FALSE;
372
373
        // Try subdatasetinfo API first
374
        // Note: this will become the only branch when subdatasetinfo will become
375
        //       available for NITF_IM, RASTERLITE and TILEDB
376
0
        const auto oSubDSInfo{GDALGetSubdatasetInfo(osSourceFilename.c_str())};
377
0
        if (oSubDSInfo && !oSubDSInfo->GetPathComponent().empty())
378
0
        {
379
0
            auto path{oSubDSInfo->GetPathComponent()};
380
0
            std::string relPath{CPLExtractRelativePath(pszVRTPath, path.c_str(),
381
0
                                                       &bRelativeToVRT)};
382
0
            osSourceFilename = oSubDSInfo->ModifyPathComponent(relPath);
383
0
            GDALDestroySubdatasetInfo(oSubDSInfo);
384
0
        }
385
0
        else
386
0
        {
387
0
            for (const char *pszSyntax :
388
0
                 GDALDataset::apszSpecialSubDatasetSyntax)
389
0
            {
390
0
                CPLString osPrefix(pszSyntax);
391
0
                osPrefix.resize(strchr(pszSyntax, ':') - pszSyntax + 1);
392
0
                if (pszSyntax[osPrefix.size()] == '"')
393
0
                    osPrefix += '"';
394
0
                if (EQUALN(osSourceFilename.c_str(), osPrefix, osPrefix.size()))
395
0
                {
396
0
                    if (STARTS_WITH_CI(pszSyntax + osPrefix.size(), "{ANY}"))
397
0
                    {
398
0
                        const char *pszLastPart =
399
0
                            strrchr(osSourceFilename.c_str(), ':') + 1;
400
                        // CSV:z:/foo.xyz
401
0
                        if ((pszLastPart[0] == '/' || pszLastPart[0] == '\\') &&
402
0
                            pszLastPart - osSourceFilename.c_str() >= 3 &&
403
0
                            pszLastPart[-3] == ':')
404
0
                            pszLastPart -= 2;
405
0
                        CPLString osPrefixFilename(osSourceFilename);
406
0
                        osPrefixFilename.resize(pszLastPart -
407
0
                                                osSourceFilename.c_str());
408
0
                        osSourceFilename = CPLExtractRelativePath(
409
0
                            pszVRTPath, pszLastPart, &bRelativeToVRT);
410
0
                        osSourceFilename = osPrefixFilename + osSourceFilename;
411
0
                    }
412
0
                    else if (STARTS_WITH_CI(pszSyntax + osPrefix.size(),
413
0
                                            "{FILENAME}"))
414
0
                    {
415
0
                        CPLString osFilename(osSourceFilename.c_str() +
416
0
                                             osPrefix.size());
417
0
                        size_t nPos = 0;
418
0
                        if (osFilename.size() >= 3 && osFilename[1] == ':' &&
419
0
                            (osFilename[2] == '\\' || osFilename[2] == '/'))
420
0
                            nPos = 2;
421
0
                        nPos = osFilename.find(
422
0
                            pszSyntax[osPrefix.size() + strlen("{FILENAME}")],
423
0
                            nPos);
424
0
                        if (nPos != std::string::npos)
425
0
                        {
426
0
                            const CPLString osSuffix = osFilename.substr(nPos);
427
0
                            osFilename.resize(nPos);
428
0
                            osSourceFilename = CPLExtractRelativePath(
429
0
                                pszVRTPath, osFilename, &bRelativeToVRT);
430
0
                            osSourceFilename =
431
0
                                osPrefix + osSourceFilename + osSuffix;
432
0
                        }
433
0
                    }
434
0
                    break;
435
0
                }
436
0
            }
437
0
        }
438
0
    }
439
0
    else
440
0
    {
441
0
        std::string osVRTFilename = pszVRTPath;
442
0
        std::string osSourceDataset = m_osSrcDSName;
443
0
        char *pszCurDir = CPLGetCurrentDir();
444
0
        if (CPLIsFilenameRelative(osSourceDataset.c_str()) &&
445
0
            !CPLIsFilenameRelative(osVRTFilename.c_str()) &&
446
0
            pszCurDir != nullptr)
447
0
        {
448
0
            osSourceDataset = CPLFormFilenameSafe(
449
0
                pszCurDir, osSourceDataset.c_str(), nullptr);
450
0
        }
451
0
        else if (!CPLIsFilenameRelative(osSourceDataset.c_str()) &&
452
0
                 CPLIsFilenameRelative(osVRTFilename.c_str()) &&
453
0
                 pszCurDir != nullptr)
454
0
        {
455
0
            osVRTFilename =
456
0
                CPLFormFilenameSafe(pszCurDir, osVRTFilename.c_str(), nullptr);
457
0
        }
458
0
        CPLFree(pszCurDir);
459
0
        osSourceFilename = CPLExtractRelativePath(
460
0
            osVRTFilename.c_str(), osSourceDataset.c_str(), &bRelativeToVRT);
461
0
    }
462
463
0
    CPLSetXMLValue(psSrc, "SourceFilename", osSourceFilename.c_str());
464
465
0
    CPLCreateXMLNode(CPLCreateXMLNode(CPLGetXMLNode(psSrc, "SourceFilename"),
466
0
                                      CXT_Attribute, "relativeToVRT"),
467
0
                     CXT_Text, bRelativeToVRT ? "1" : "0");
468
469
    // Determine if we must write the shared attribute. The config option
470
    // will override the m_nExplicitSharedStatus value
471
0
    const char *pszShared = CPLGetConfigOption("VRT_SHARED_SOURCE", nullptr);
472
0
    if ((pszShared == nullptr && m_nExplicitSharedStatus == 0) ||
473
0
        (pszShared != nullptr && !CPLTestBool(pszShared)))
474
0
    {
475
0
        CPLCreateXMLNode(
476
0
            CPLCreateXMLNode(CPLGetXMLNode(psSrc, "SourceFilename"),
477
0
                             CXT_Attribute, "shared"),
478
0
            CXT_Text, "0");
479
0
    }
480
0
}
481
482
/************************************************************************/
483
/*                           SerializeToXML()                           */
484
/************************************************************************/
485
486
CPLXMLNode *VRTSimpleSource::SerializeToXML(const char *pszVRTPath)
487
488
0
{
489
0
    CPLXMLNode *const psSrc =
490
0
        CPLCreateXMLNode(nullptr, CXT_Element, GetTypeStatic());
491
492
0
    if (!m_osResampling.empty())
493
0
    {
494
0
        CPLCreateXMLNode(CPLCreateXMLNode(psSrc, CXT_Attribute, "resampling"),
495
0
                         CXT_Text, m_osResampling.c_str());
496
0
    }
497
498
0
    if (!m_osName.empty())
499
0
    {
500
0
        CPLAddXMLAttributeAndValue(psSrc, "name", m_osName.c_str());
501
0
    }
502
503
0
    if (m_bSrcDSNameFromVRT)
504
0
    {
505
0
        CPLAddXMLChild(psSrc, CPLParseXMLString(m_osSrcDSName.c_str()));
506
0
    }
507
0
    else
508
0
    {
509
0
        bool bDone = false;
510
0
        if (m_osSrcDSName.empty() && m_poRasterBand)
511
0
        {
512
0
            auto poSrcDS = m_poRasterBand->GetDataset();
513
0
            if (poSrcDS)
514
0
            {
515
0
                VRTDataset *poSrcVRTDS = nullptr;
516
                // For GDALComputedDataset
517
0
                void *pHandle = poSrcDS->GetInternalHandle("VRT_DATASET");
518
0
                if (pHandle && poSrcDS->GetInternalHandle(nullptr) == nullptr)
519
0
                {
520
0
                    poSrcVRTDS = static_cast<VRTDataset *>(pHandle);
521
0
                }
522
0
                else
523
0
                {
524
0
                    poSrcVRTDS = dynamic_cast<VRTDataset *>(poSrcDS);
525
0
                }
526
0
                if (poSrcVRTDS)
527
0
                {
528
0
                    poSrcVRTDS->UnsetPreservedRelativeFilenames();
529
0
                    CPLAddXMLChild(psSrc,
530
0
                                   poSrcVRTDS->SerializeToXML(pszVRTPath));
531
0
                    bDone = true;
532
0
                }
533
0
            }
534
0
        }
535
0
        if (!bDone)
536
0
        {
537
0
            AddSourceFilenameNode(pszVRTPath, psSrc);
538
0
        }
539
0
    }
540
541
0
    GDALSerializeOpenOptionsToXML(psSrc, m_aosOpenOptionsOri.List());
542
543
0
    if (m_bGetMaskBand)
544
0
        CPLSetXMLValue(psSrc, "SourceBand", CPLSPrintf("mask,%d", m_nBand));
545
0
    else
546
0
        CPLSetXMLValue(psSrc, "SourceBand", CPLSPrintf("%d", m_nBand));
547
548
    // TODO: in a later version, no longer emit SourceProperties, which
549
    // is no longer used by GDAL 3.4
550
0
    if (m_poRasterBand)
551
0
    {
552
        /* Write a few additional useful properties of the dataset */
553
        /* so that we can use a proxy dataset when re-opening. See XMLInit() */
554
        /* below */
555
0
        CPLSetXMLValue(psSrc, "SourceProperties.#RasterXSize",
556
0
                       CPLSPrintf("%d", m_poRasterBand->GetXSize()));
557
0
        CPLSetXMLValue(psSrc, "SourceProperties.#RasterYSize",
558
0
                       CPLSPrintf("%d", m_poRasterBand->GetYSize()));
559
0
        CPLSetXMLValue(
560
0
            psSrc, "SourceProperties.#DataType",
561
0
            GDALGetDataTypeName(m_poRasterBand->GetRasterDataType()));
562
563
0
        int nBlockXSize = 0;
564
0
        int nBlockYSize = 0;
565
0
        m_poRasterBand->GetBlockSize(&nBlockXSize, &nBlockYSize);
566
567
0
        CPLSetXMLValue(psSrc, "SourceProperties.#BlockXSize",
568
0
                       CPLSPrintf("%d", nBlockXSize));
569
0
        CPLSetXMLValue(psSrc, "SourceProperties.#BlockYSize",
570
0
                       CPLSPrintf("%d", nBlockYSize));
571
0
    }
572
573
0
    if (IsSrcWinSet())
574
0
    {
575
0
        CPLSetXMLValue(psSrc, "SrcRect.#xOff",
576
0
                       CPLSPrintf("%.15g", m_dfSrcXOff));
577
0
        CPLSetXMLValue(psSrc, "SrcRect.#yOff",
578
0
                       CPLSPrintf("%.15g", m_dfSrcYOff));
579
0
        CPLSetXMLValue(psSrc, "SrcRect.#xSize",
580
0
                       CPLSPrintf("%.15g", m_dfSrcXSize));
581
0
        CPLSetXMLValue(psSrc, "SrcRect.#ySize",
582
0
                       CPLSPrintf("%.15g", m_dfSrcYSize));
583
0
    }
584
585
0
    if (IsDstWinSet())
586
0
    {
587
0
        CPLSetXMLValue(psSrc, "DstRect.#xOff",
588
0
                       CPLSPrintf("%.15g", m_dfDstXOff));
589
0
        CPLSetXMLValue(psSrc, "DstRect.#yOff",
590
0
                       CPLSPrintf("%.15g", m_dfDstYOff));
591
0
        CPLSetXMLValue(psSrc, "DstRect.#xSize",
592
0
                       CPLSPrintf("%.15g", m_dfDstXSize));
593
0
        CPLSetXMLValue(psSrc, "DstRect.#ySize",
594
0
                       CPLSPrintf("%.15g", m_dfDstYSize));
595
0
    }
596
597
0
    return psSrc;
598
0
}
599
600
/************************************************************************/
601
/*                              XMLInit()                               */
602
/************************************************************************/
603
604
CPLErr VRTSimpleSource::XMLInit(const CPLXMLNode *psSrc, const char *pszVRTPath,
605
                                VRTMapSharedResources &oMapSharedSources)
606
607
0
{
608
0
    m_poMapSharedSources = &oMapSharedSources;
609
610
0
    m_osResampling = CPLGetXMLValue(psSrc, "resampling", "");
611
0
    m_osName = CPLGetXMLValue(psSrc, "name", "");
612
613
    /* -------------------------------------------------------------------- */
614
    /*      Prepare filename.                                               */
615
    /* -------------------------------------------------------------------- */
616
0
    const CPLXMLNode *psSourceFileNameNode =
617
0
        CPLGetXMLNode(psSrc, "SourceFilename");
618
0
    const CPLXMLNode *psSourceVRTDataset = CPLGetXMLNode(psSrc, "VRTDataset");
619
0
    const char *pszFilename =
620
0
        psSourceFileNameNode ? CPLGetXMLValue(psSourceFileNameNode, nullptr, "")
621
0
                             : "";
622
623
0
    if (pszFilename[0] == '\0' && !psSourceVRTDataset)
624
0
    {
625
0
        CPLError(CE_Warning, CPLE_AppDefined,
626
0
                 "Missing <SourceFilename> or <VRTDataset> element in <%s>.",
627
0
                 psSrc->pszValue);
628
0
        return CE_Failure;
629
0
    }
630
631
    // Backup original filename and relativeToVRT so as to be able to
632
    // serialize them identically again (#5985)
633
0
    m_osSourceFileNameOri = pszFilename;
634
0
    if (pszFilename[0])
635
0
    {
636
0
        m_bRelativeToVRTOri =
637
0
            atoi(CPLGetXMLValue(psSourceFileNameNode, "relativetoVRT", "0"));
638
0
        const char *pszShared =
639
0
            CPLGetXMLValue(psSourceFileNameNode, "shared", nullptr);
640
0
        if (pszShared == nullptr)
641
0
        {
642
0
            pszShared = CPLGetConfigOption("VRT_SHARED_SOURCE", nullptr);
643
0
        }
644
0
        if (pszShared != nullptr)
645
0
        {
646
0
            m_nExplicitSharedStatus = CPLTestBool(pszShared);
647
0
        }
648
649
0
        m_osSrcDSName = GDALDataset::BuildFilename(
650
0
            pszFilename, pszVRTPath, CPL_TO_BOOL(m_bRelativeToVRTOri));
651
0
    }
652
0
    else if (psSourceVRTDataset)
653
0
    {
654
0
        CPLXMLNode sNode;
655
0
        sNode.eType = psSourceVRTDataset->eType;
656
0
        sNode.pszValue = psSourceVRTDataset->pszValue;
657
0
        sNode.psNext = nullptr;
658
0
        sNode.psChild = psSourceVRTDataset->psChild;
659
0
        char *pszXML = CPLSerializeXMLTree(&sNode);
660
0
        if (pszXML)
661
0
        {
662
0
            m_bSrcDSNameFromVRT = true;
663
0
            m_osSrcDSName = pszXML;
664
0
            CPLFree(pszXML);
665
0
        }
666
0
    }
667
668
0
    const char *pszSourceBand = CPLGetXMLValue(psSrc, "SourceBand", "1");
669
0
    m_bGetMaskBand = false;
670
0
    if (STARTS_WITH_CI(pszSourceBand, "mask"))
671
0
    {
672
0
        m_bGetMaskBand = true;
673
0
        if (pszSourceBand[4] == ',')
674
0
            m_nBand = atoi(pszSourceBand + 5);
675
0
        else
676
0
            m_nBand = 1;
677
0
    }
678
0
    else
679
0
    {
680
0
        m_nBand = atoi(pszSourceBand);
681
0
    }
682
0
    if (!GDALCheckBandCount(m_nBand, 0))
683
0
    {
684
0
        CPLError(CE_Warning, CPLE_AppDefined,
685
0
                 "Invalid <SourceBand> element in VRTRasterBand.");
686
0
        return CE_Failure;
687
0
    }
688
689
0
    m_aosOpenOptions = GDALDeserializeOpenOptionsFromXML(psSrc);
690
0
    m_aosOpenOptionsOri = m_aosOpenOptions;
691
0
    if (strstr(m_osSrcDSName.c_str(), "<VRTDataset") != nullptr)
692
0
        m_aosOpenOptions.SetNameValue("ROOT_PATH", pszVRTPath);
693
694
0
    return ParseSrcRectAndDstRect(psSrc);
695
0
}
696
697
/************************************************************************/
698
/*                        ParseSrcRectAndDstRect()                      */
699
/************************************************************************/
700
701
CPLErr VRTSimpleSource::ParseSrcRectAndDstRect(const CPLXMLNode *psSrc)
702
0
{
703
0
    const auto GetAttrValue = [](const CPLXMLNode *psNode,
704
0
                                 const char *pszAttrName, double dfDefaultVal)
705
0
    {
706
0
        if (const char *pszVal = CPLGetXMLValue(psNode, pszAttrName, nullptr))
707
0
            return CPLAtof(pszVal);
708
0
        else
709
0
            return dfDefaultVal;
710
0
    };
711
712
    /* -------------------------------------------------------------------- */
713
    /*      Set characteristics.                                            */
714
    /* -------------------------------------------------------------------- */
715
0
    const CPLXMLNode *const psSrcRect = CPLGetXMLNode(psSrc, "SrcRect");
716
0
    if (psSrcRect)
717
0
    {
718
0
        double xOff = GetAttrValue(psSrcRect, "xOff", UNINIT_WINDOW);
719
0
        double yOff = GetAttrValue(psSrcRect, "yOff", UNINIT_WINDOW);
720
0
        double xSize = GetAttrValue(psSrcRect, "xSize", UNINIT_WINDOW);
721
0
        double ySize = GetAttrValue(psSrcRect, "ySize", UNINIT_WINDOW);
722
        // Test written that way to catch NaN values
723
0
        if (!(xOff >= INT_MIN && xOff <= INT_MAX) ||
724
0
            !(yOff >= INT_MIN && yOff <= INT_MAX) ||
725
0
            !(xSize > 0 || xSize == UNINIT_WINDOW) || xSize > INT_MAX ||
726
0
            !(ySize > 0 || ySize == UNINIT_WINDOW) || ySize > INT_MAX)
727
0
        {
728
0
            CPLError(CE_Failure, CPLE_AppDefined, "Wrong values in SrcRect");
729
0
            return CE_Failure;
730
0
        }
731
0
        SetSrcWindow(xOff, yOff, xSize, ySize);
732
0
    }
733
0
    else
734
0
    {
735
0
        m_dfSrcXOff = UNINIT_WINDOW;
736
0
        m_dfSrcYOff = UNINIT_WINDOW;
737
0
        m_dfSrcXSize = UNINIT_WINDOW;
738
0
        m_dfSrcYSize = UNINIT_WINDOW;
739
0
    }
740
741
0
    const CPLXMLNode *const psDstRect = CPLGetXMLNode(psSrc, "DstRect");
742
0
    if (psDstRect)
743
0
    {
744
0
        double xOff = GetAttrValue(psDstRect, "xOff", UNINIT_WINDOW);
745
0
        ;
746
0
        double yOff = GetAttrValue(psDstRect, "yOff", UNINIT_WINDOW);
747
0
        double xSize = GetAttrValue(psDstRect, "xSize", UNINIT_WINDOW);
748
0
        ;
749
0
        double ySize = GetAttrValue(psDstRect, "ySize", UNINIT_WINDOW);
750
        // Test written that way to catch NaN values
751
0
        if (!(xOff >= INT_MIN && xOff <= INT_MAX) ||
752
0
            !(yOff >= INT_MIN && yOff <= INT_MAX) ||
753
0
            !(xSize > 0 || xSize == UNINIT_WINDOW) || xSize > INT_MAX ||
754
0
            !(ySize > 0 || ySize == UNINIT_WINDOW) || ySize > INT_MAX)
755
0
        {
756
0
            CPLError(CE_Failure, CPLE_AppDefined, "Wrong values in DstRect");
757
0
            return CE_Failure;
758
0
        }
759
0
        SetDstWindow(xOff, yOff, xSize, ySize);
760
0
    }
761
0
    else
762
0
    {
763
0
        m_dfDstXOff = UNINIT_WINDOW;
764
0
        m_dfDstYOff = UNINIT_WINDOW;
765
0
        m_dfDstXSize = UNINIT_WINDOW;
766
0
        m_dfDstYSize = UNINIT_WINDOW;
767
0
    }
768
769
0
    return CE_None;
770
0
}
771
772
/************************************************************************/
773
/*                             GetFileList()                            */
774
/************************************************************************/
775
776
void VRTSimpleSource::GetFileList(char ***ppapszFileList, int *pnSize,
777
                                  int *pnMaxSize, CPLHashSet *hSetFiles)
778
0
{
779
0
    if (!m_osSrcDSName.empty() && !m_bSrcDSNameFromVRT)
780
0
    {
781
0
        const char *pszFilename = m_osSrcDSName.c_str();
782
783
        /* --------------------------------------------------------------------
784
         */
785
        /*      Is it already in the list ? */
786
        /* --------------------------------------------------------------------
787
         */
788
0
        if (CPLHashSetLookup(hSetFiles, pszFilename) != nullptr)
789
0
            return;
790
791
        /* --------------------------------------------------------------------
792
         */
793
        /*      Grow array if necessary */
794
        /* --------------------------------------------------------------------
795
         */
796
0
        if (*pnSize + 1 >= *pnMaxSize)
797
0
        {
798
0
            *pnMaxSize = std::max(*pnSize + 2, 2 + 2 * (*pnMaxSize));
799
0
            *ppapszFileList = static_cast<char **>(
800
0
                CPLRealloc(*ppapszFileList, sizeof(char *) * (*pnMaxSize)));
801
0
        }
802
803
        /* --------------------------------------------------------------------
804
         */
805
        /*      Add the string to the list */
806
        /* --------------------------------------------------------------------
807
         */
808
0
        (*ppapszFileList)[*pnSize] = CPLStrdup(pszFilename);
809
0
        (*ppapszFileList)[(*pnSize + 1)] = nullptr;
810
0
        CPLHashSetInsert(hSetFiles, (*ppapszFileList)[*pnSize]);
811
812
0
        (*pnSize)++;
813
0
    }
814
0
}
815
816
/************************************************************************/
817
/*                           OpenSource()                               */
818
/************************************************************************/
819
820
void VRTSimpleSource::OpenSource() const
821
0
{
822
0
    CPLAssert(m_poRasterBand == nullptr);
823
824
    /* ----------------------------------------------------------------- */
825
    /*      Create a proxy dataset                                       */
826
    /* ----------------------------------------------------------------- */
827
0
    GDALProxyPoolDataset *proxyDS = nullptr;
828
0
    std::string osKeyMapSharedSources;
829
0
    if (m_poMapSharedSources)
830
0
    {
831
0
        osKeyMapSharedSources = m_osSrcDSName;
832
0
        for (int i = 0; i < m_aosOpenOptions.size(); ++i)
833
0
        {
834
0
            osKeyMapSharedSources += "||";
835
0
            osKeyMapSharedSources += m_aosOpenOptions[i];
836
0
        }
837
838
0
        proxyDS = cpl::down_cast<GDALProxyPoolDataset *>(
839
0
            m_poMapSharedSources->Get(osKeyMapSharedSources));
840
0
    }
841
842
0
    if (proxyDS == nullptr)
843
0
    {
844
0
        int bShared = true;
845
0
        if (m_nExplicitSharedStatus != -1)
846
0
            bShared = m_nExplicitSharedStatus;
847
848
0
        const CPLString osUniqueHandle(CPLSPrintf("%p", m_poMapSharedSources));
849
0
        proxyDS = GDALProxyPoolDataset::Create(
850
0
            m_osSrcDSName, m_aosOpenOptions.List(), GA_ReadOnly, bShared,
851
0
            osUniqueHandle.c_str());
852
0
        if (proxyDS == nullptr)
853
0
            return;
854
0
    }
855
0
    else
856
0
    {
857
0
        proxyDS->Reference();
858
0
    }
859
860
0
    if (m_bGetMaskBand)
861
0
    {
862
0
        GDALProxyPoolRasterBand *poMaskBand =
863
0
            cpl::down_cast<GDALProxyPoolRasterBand *>(
864
0
                proxyDS->GetRasterBand(m_nBand));
865
0
        poMaskBand->AddSrcMaskBandDescriptionFromUnderlying();
866
0
    }
867
868
    /* -------------------------------------------------------------------- */
869
    /*      Get the raster band.                                            */
870
    /* -------------------------------------------------------------------- */
871
872
0
    m_poRasterBand = proxyDS->GetRasterBand(m_nBand);
873
0
    if (m_poRasterBand == nullptr || !ValidateOpenedBand(m_poRasterBand))
874
0
    {
875
0
        proxyDS->ReleaseRef();
876
0
        return;
877
0
    }
878
879
0
    if (m_bGetMaskBand)
880
0
    {
881
0
        m_poRasterBand = m_poRasterBand->GetMaskBand();
882
0
        if (m_poRasterBand == nullptr)
883
0
        {
884
0
            proxyDS->ReleaseRef();
885
0
            return;
886
0
        }
887
0
        m_poMaskBandMainBand = m_poRasterBand;
888
0
    }
889
890
0
    if (m_poMapSharedSources)
891
0
    {
892
0
        m_poMapSharedSources->Insert(osKeyMapSharedSources, proxyDS);
893
0
    }
894
0
}
895
896
/************************************************************************/
897
/*                         GetRasterBand()                              */
898
/************************************************************************/
899
900
GDALRasterBand *VRTSimpleSource::GetRasterBand() const
901
0
{
902
0
    if (m_poRasterBand == nullptr)
903
0
        OpenSource();
904
0
    return m_poRasterBand;
905
0
}
906
907
/************************************************************************/
908
/*                        GetMaskBandMainBand()                         */
909
/************************************************************************/
910
911
GDALRasterBand *VRTSimpleSource::GetMaskBandMainBand()
912
0
{
913
0
    if (m_poRasterBand == nullptr)
914
0
        OpenSource();
915
0
    return m_poMaskBandMainBand;
916
0
}
917
918
/************************************************************************/
919
/*                       IsSameExceptBandNumber()                       */
920
/************************************************************************/
921
922
bool VRTSimpleSource::IsSameExceptBandNumber(
923
    const VRTSimpleSource *poOtherSource) const
924
0
{
925
0
    return m_dfSrcXOff == poOtherSource->m_dfSrcXOff &&
926
0
           m_dfSrcYOff == poOtherSource->m_dfSrcYOff &&
927
0
           m_dfSrcXSize == poOtherSource->m_dfSrcXSize &&
928
0
           m_dfSrcYSize == poOtherSource->m_dfSrcYSize &&
929
0
           m_dfDstXOff == poOtherSource->m_dfDstXOff &&
930
0
           m_dfDstYOff == poOtherSource->m_dfDstYOff &&
931
0
           m_dfDstXSize == poOtherSource->m_dfDstXSize &&
932
0
           m_dfDstYSize == poOtherSource->m_dfDstYSize &&
933
0
           m_osSrcDSName == poOtherSource->m_osSrcDSName;
934
0
}
935
936
/************************************************************************/
937
/*                              SrcToDst()                              */
938
/*                                                                      */
939
/*      Note: this is a no-op if the both src and dst windows are unset */
940
/************************************************************************/
941
942
void VRTSimpleSource::SrcToDst(double dfX, double dfY, double &dfXOut,
943
                               double &dfYOut) const
944
945
0
{
946
0
    dfXOut = ((dfX - m_dfSrcXOff) / m_dfSrcXSize) * m_dfDstXSize + m_dfDstXOff;
947
0
    dfYOut = ((dfY - m_dfSrcYOff) / m_dfSrcYSize) * m_dfDstYSize + m_dfDstYOff;
948
0
}
949
950
/************************************************************************/
951
/*                              DstToSrc()                              */
952
/*                                                                      */
953
/*      Note: this is a no-op if the both src and dst windows are unset */
954
/************************************************************************/
955
956
void VRTSimpleSource::DstToSrc(double dfX, double dfY, double &dfXOut,
957
                               double &dfYOut) const
958
959
0
{
960
0
    dfXOut = ((dfX - m_dfDstXOff) / m_dfDstXSize) * m_dfSrcXSize + m_dfSrcXOff;
961
0
    dfYOut = ((dfY - m_dfDstYOff) / m_dfDstYSize) * m_dfSrcYSize + m_dfSrcYOff;
962
0
}
963
964
/************************************************************************/
965
/*                          GetSrcDstWindow()                           */
966
/************************************************************************/
967
968
int VRTSimpleSource::GetSrcDstWindow(
969
    double dfXOff, double dfYOff, double dfXSize, double dfYSize, int nBufXSize,
970
    int nBufYSize, double *pdfReqXOff, double *pdfReqYOff, double *pdfReqXSize,
971
    double *pdfReqYSize, int *pnReqXOff, int *pnReqYOff, int *pnReqXSize,
972
    int *pnReqYSize, int *pnOutXOff, int *pnOutYOff, int *pnOutXSize,
973
    int *pnOutYSize, bool &bErrorOut)
974
975
0
{
976
0
    bErrorOut = false;
977
978
0
    if (m_dfSrcXSize == 0.0 || m_dfSrcYSize == 0.0 || m_dfDstXSize == 0.0 ||
979
0
        m_dfDstYSize == 0.0)
980
0
    {
981
0
        return FALSE;
982
0
    }
983
984
0
    const bool bDstWinSet = IsDstWinSet();
985
986
0
#ifdef DEBUG
987
0
    const bool bSrcWinSet = IsSrcWinSet();
988
989
0
    if (bSrcWinSet != bDstWinSet)
990
0
    {
991
0
        return FALSE;
992
0
    }
993
0
#endif
994
995
    /* -------------------------------------------------------------------- */
996
    /*      If the input window completely misses the portion of the        */
997
    /*      virtual dataset provided by this source we have nothing to do.  */
998
    /* -------------------------------------------------------------------- */
999
0
    if (bDstWinSet)
1000
0
    {
1001
0
        if (dfXOff >= m_dfDstXOff + m_dfDstXSize ||
1002
0
            dfYOff >= m_dfDstYOff + m_dfDstYSize ||
1003
0
            dfXOff + dfXSize <= m_dfDstXOff || dfYOff + dfYSize <= m_dfDstYOff)
1004
0
            return FALSE;
1005
0
    }
1006
1007
    /* -------------------------------------------------------------------- */
1008
    /*      This request window corresponds to the whole output buffer.     */
1009
    /* -------------------------------------------------------------------- */
1010
0
    *pnOutXOff = 0;
1011
0
    *pnOutYOff = 0;
1012
0
    *pnOutXSize = nBufXSize;
1013
0
    *pnOutYSize = nBufYSize;
1014
1015
    /* -------------------------------------------------------------------- */
1016
    /*      If the input window extents outside the portion of the on       */
1017
    /*      the virtual file that this source can set, then clip down       */
1018
    /*      the requested window.                                           */
1019
    /* -------------------------------------------------------------------- */
1020
0
    bool bModifiedX = false;
1021
0
    bool bModifiedY = false;
1022
0
    double dfRXOff = dfXOff;
1023
0
    double dfRYOff = dfYOff;
1024
0
    double dfRXSize = dfXSize;
1025
0
    double dfRYSize = dfYSize;
1026
1027
0
    if (bDstWinSet)
1028
0
    {
1029
0
        if (dfRXOff < m_dfDstXOff)
1030
0
        {
1031
0
            dfRXSize = dfRXSize + dfRXOff - m_dfDstXOff;
1032
0
            dfRXOff = m_dfDstXOff;
1033
0
            bModifiedX = true;
1034
0
        }
1035
1036
0
        if (dfRYOff < m_dfDstYOff)
1037
0
        {
1038
0
            dfRYSize = dfRYSize + dfRYOff - m_dfDstYOff;
1039
0
            dfRYOff = m_dfDstYOff;
1040
0
            bModifiedY = true;
1041
0
        }
1042
1043
0
        if (dfRXOff + dfRXSize > m_dfDstXOff + m_dfDstXSize)
1044
0
        {
1045
0
            dfRXSize = m_dfDstXOff + m_dfDstXSize - dfRXOff;
1046
0
            bModifiedX = true;
1047
0
        }
1048
1049
0
        if (dfRYOff + dfRYSize > m_dfDstYOff + m_dfDstYSize)
1050
0
        {
1051
0
            dfRYSize = m_dfDstYOff + m_dfDstYSize - dfRYOff;
1052
0
            bModifiedY = true;
1053
0
        }
1054
0
    }
1055
1056
    /* -------------------------------------------------------------------- */
1057
    /*      Translate requested region in virtual file into the source      */
1058
    /*      band coordinates.                                               */
1059
    /* -------------------------------------------------------------------- */
1060
0
    const double dfScaleX = m_dfSrcXSize / m_dfDstXSize;
1061
0
    const double dfScaleY = m_dfSrcYSize / m_dfDstYSize;
1062
1063
0
    *pdfReqXOff = (dfRXOff - m_dfDstXOff) * dfScaleX + m_dfSrcXOff;
1064
0
    *pdfReqYOff = (dfRYOff - m_dfDstYOff) * dfScaleY + m_dfSrcYOff;
1065
0
    *pdfReqXSize = dfRXSize * dfScaleX;
1066
0
    *pdfReqYSize = dfRYSize * dfScaleY;
1067
1068
0
    if (!std::isfinite(*pdfReqXOff) || !std::isfinite(*pdfReqYOff) ||
1069
0
        !std::isfinite(*pdfReqXSize) || !std::isfinite(*pdfReqYSize) ||
1070
0
        *pdfReqXOff > INT_MAX || *pdfReqYOff > INT_MAX || *pdfReqXSize < 0 ||
1071
0
        *pdfReqYSize < 0)
1072
0
    {
1073
0
        return FALSE;
1074
0
    }
1075
1076
    /* -------------------------------------------------------------------- */
1077
    /*      Clamp within the bounds of the available source data.           */
1078
    /* -------------------------------------------------------------------- */
1079
0
    if (*pdfReqXOff < 0)
1080
0
    {
1081
0
        *pdfReqXSize += *pdfReqXOff;
1082
0
        *pdfReqXOff = 0;
1083
0
        bModifiedX = true;
1084
0
    }
1085
0
    if (*pdfReqYOff < 0)
1086
0
    {
1087
0
        *pdfReqYSize += *pdfReqYOff;
1088
0
        *pdfReqYOff = 0;
1089
0
        bModifiedY = true;
1090
0
    }
1091
1092
0
    *pnReqXOff = static_cast<int>(floor(*pdfReqXOff));
1093
0
    *pnReqYOff = static_cast<int>(floor(*pdfReqYOff));
1094
1095
0
    constexpr double EPS = 1e-3;
1096
0
    constexpr double ONE_MINUS_EPS = 1.0 - EPS;
1097
0
    if (*pdfReqXOff - *pnReqXOff > ONE_MINUS_EPS)
1098
0
    {
1099
0
        (*pnReqXOff)++;
1100
0
        *pdfReqXOff = *pnReqXOff;
1101
0
    }
1102
0
    if (*pdfReqYOff - *pnReqYOff > ONE_MINUS_EPS)
1103
0
    {
1104
0
        (*pnReqYOff)++;
1105
0
        *pdfReqYOff = *pnReqYOff;
1106
0
    }
1107
1108
0
    if (*pdfReqXSize > INT_MAX)
1109
0
        *pnReqXSize = INT_MAX;
1110
0
    else
1111
0
        *pnReqXSize = static_cast<int>(floor(*pdfReqXSize + 0.5));
1112
1113
0
    if (*pdfReqYSize > INT_MAX)
1114
0
        *pnReqYSize = INT_MAX;
1115
0
    else
1116
0
        *pnReqYSize = static_cast<int>(floor(*pdfReqYSize + 0.5));
1117
1118
    /* -------------------------------------------------------------------- */
1119
    /*      Clamp within the bounds of the available source data.           */
1120
    /* -------------------------------------------------------------------- */
1121
1122
0
    if (*pnReqXSize == 0)
1123
0
        *pnReqXSize = 1;
1124
0
    if (*pnReqYSize == 0)
1125
0
        *pnReqYSize = 1;
1126
1127
0
    auto l_band = GetRasterBand();
1128
0
    if (!l_band)
1129
0
    {
1130
0
        bErrorOut = true;
1131
0
        return FALSE;
1132
0
    }
1133
0
    if (*pnReqXSize > INT_MAX - *pnReqXOff ||
1134
0
        *pnReqXOff + *pnReqXSize > l_band->GetXSize())
1135
0
    {
1136
0
        *pnReqXSize = l_band->GetXSize() - *pnReqXOff;
1137
0
        bModifiedX = true;
1138
0
    }
1139
0
    if (*pdfReqXOff + *pdfReqXSize > l_band->GetXSize())
1140
0
    {
1141
0
        *pdfReqXSize = l_band->GetXSize() - *pdfReqXOff;
1142
0
        bModifiedX = true;
1143
0
    }
1144
1145
0
    if (*pnReqYSize > INT_MAX - *pnReqYOff ||
1146
0
        *pnReqYOff + *pnReqYSize > l_band->GetYSize())
1147
0
    {
1148
0
        *pnReqYSize = l_band->GetYSize() - *pnReqYOff;
1149
0
        bModifiedY = true;
1150
0
    }
1151
0
    if (*pdfReqYOff + *pdfReqYSize > l_band->GetYSize())
1152
0
    {
1153
0
        *pdfReqYSize = l_band->GetYSize() - *pdfReqYOff;
1154
0
        bModifiedY = true;
1155
0
    }
1156
1157
    /* -------------------------------------------------------------------- */
1158
    /*      Don't do anything if the requesting region is completely off    */
1159
    /*      the source image.                                               */
1160
    /* -------------------------------------------------------------------- */
1161
0
    if (*pnReqXOff >= l_band->GetXSize() || *pnReqYOff >= l_band->GetYSize() ||
1162
0
        *pnReqXSize <= 0 || *pnReqYSize <= 0)
1163
0
    {
1164
0
        return FALSE;
1165
0
    }
1166
1167
    /* -------------------------------------------------------------------- */
1168
    /*      If we haven't had to modify the source rectangle, then the      */
1169
    /*      destination rectangle must be the whole region.                 */
1170
    /* -------------------------------------------------------------------- */
1171
0
    if (bModifiedX || bModifiedY)
1172
0
    {
1173
        /* --------------------------------------------------------------------
1174
         */
1175
        /*      Now transform this possibly reduced request back into the */
1176
        /*      destination buffer coordinates in case the output region is */
1177
        /*      less than the whole buffer. */
1178
        /* --------------------------------------------------------------------
1179
         */
1180
0
        double dfDstULX = 0.0;
1181
0
        double dfDstULY = 0.0;
1182
0
        double dfDstLRX = 0.0;
1183
0
        double dfDstLRY = 0.0;
1184
1185
0
        SrcToDst(*pdfReqXOff, *pdfReqYOff, dfDstULX, dfDstULY);
1186
0
        SrcToDst(*pdfReqXOff + *pdfReqXSize, *pdfReqYOff + *pdfReqYSize,
1187
0
                 dfDstLRX, dfDstLRY);
1188
#if DEBUG_VERBOSE
1189
        CPLDebug("VRT", "dfDstULX=%g dfDstULY=%g dfDstLRX=%g dfDstLRY=%g",
1190
                 dfDstULX, dfDstULY, dfDstLRX, dfDstLRY);
1191
#endif
1192
1193
0
        if (bModifiedX)
1194
0
        {
1195
0
            const double dfScaleWinToBufX = nBufXSize / dfXSize;
1196
1197
0
            const double dfOutXOff = (dfDstULX - dfXOff) * dfScaleWinToBufX;
1198
0
            if (dfOutXOff <= 0)
1199
0
                *pnOutXOff = 0;
1200
0
            else if (dfOutXOff > INT_MAX)
1201
0
                *pnOutXOff = INT_MAX;
1202
0
            else
1203
0
                *pnOutXOff = static_cast<int>(dfOutXOff + EPS);
1204
1205
            // Apply correction on floating-point source window
1206
0
            {
1207
0
                double dfDstDeltaX =
1208
0
                    (dfOutXOff - *pnOutXOff) / dfScaleWinToBufX;
1209
0
                double dfSrcDeltaX = dfDstDeltaX / m_dfDstXSize * m_dfSrcXSize;
1210
0
                *pdfReqXOff -= dfSrcDeltaX;
1211
0
                *pdfReqXSize = std::min(*pdfReqXSize + dfSrcDeltaX,
1212
0
                                        static_cast<double>(INT_MAX));
1213
0
            }
1214
1215
0
            double dfOutRightXOff = (dfDstLRX - dfXOff) * dfScaleWinToBufX;
1216
0
            if (dfOutRightXOff < dfOutXOff)
1217
0
                return FALSE;
1218
0
            if (dfOutRightXOff > INT_MAX)
1219
0
                dfOutRightXOff = INT_MAX;
1220
0
            const int nOutRightXOff =
1221
0
                static_cast<int>(ceil(dfOutRightXOff - EPS));
1222
0
            *pnOutXSize = nOutRightXOff - *pnOutXOff;
1223
1224
0
            if (*pnOutXSize > INT_MAX - *pnOutXOff ||
1225
0
                *pnOutXOff + *pnOutXSize > nBufXSize)
1226
0
                *pnOutXSize = nBufXSize - *pnOutXOff;
1227
1228
            // Apply correction on floating-point source window
1229
0
            {
1230
0
                double dfDstDeltaX =
1231
0
                    (nOutRightXOff - dfOutRightXOff) / dfScaleWinToBufX;
1232
0
                double dfSrcDeltaX = dfDstDeltaX / m_dfDstXSize * m_dfSrcXSize;
1233
0
                *pdfReqXSize = std::min(*pdfReqXSize + dfSrcDeltaX,
1234
0
                                        static_cast<double>(INT_MAX));
1235
0
            }
1236
0
        }
1237
1238
0
        if (bModifiedY)
1239
0
        {
1240
0
            const double dfScaleWinToBufY = nBufYSize / dfYSize;
1241
1242
0
            const double dfOutYOff = (dfDstULY - dfYOff) * dfScaleWinToBufY;
1243
0
            if (dfOutYOff <= 0)
1244
0
                *pnOutYOff = 0;
1245
0
            else if (dfOutYOff > INT_MAX)
1246
0
                *pnOutYOff = INT_MAX;
1247
0
            else
1248
0
                *pnOutYOff = static_cast<int>(dfOutYOff + EPS);
1249
1250
            // Apply correction on floating-point source window
1251
0
            {
1252
0
                double dfDstDeltaY =
1253
0
                    (dfOutYOff - *pnOutYOff) / dfScaleWinToBufY;
1254
0
                double dfSrcDeltaY = dfDstDeltaY / m_dfDstYSize * m_dfSrcYSize;
1255
0
                *pdfReqYOff -= dfSrcDeltaY;
1256
0
                *pdfReqYSize = std::min(*pdfReqYSize + dfSrcDeltaY,
1257
0
                                        static_cast<double>(INT_MAX));
1258
0
            }
1259
1260
0
            double dfOutTopYOff = (dfDstLRY - dfYOff) * dfScaleWinToBufY;
1261
0
            if (dfOutTopYOff < dfOutYOff)
1262
0
                return FALSE;
1263
0
            if (dfOutTopYOff > INT_MAX)
1264
0
                dfOutTopYOff = INT_MAX;
1265
0
            const int nOutTopYOff = static_cast<int>(ceil(dfOutTopYOff - EPS));
1266
0
            *pnOutYSize = nOutTopYOff - *pnOutYOff;
1267
1268
0
            if (*pnOutYSize > INT_MAX - *pnOutYOff ||
1269
0
                *pnOutYOff + *pnOutYSize > nBufYSize)
1270
0
                *pnOutYSize = nBufYSize - *pnOutYOff;
1271
1272
            // Apply correction on floating-point source window
1273
0
            {
1274
0
                double dfDstDeltaY =
1275
0
                    (nOutTopYOff - dfOutTopYOff) / dfScaleWinToBufY;
1276
0
                double dfSrcDeltaY = dfDstDeltaY / m_dfDstYSize * m_dfSrcYSize;
1277
0
                *pdfReqYSize = std::min(*pdfReqYSize + dfSrcDeltaY,
1278
0
                                        static_cast<double>(INT_MAX));
1279
0
            }
1280
0
        }
1281
1282
0
        if (*pnOutXSize < 1 || *pnOutYSize < 1)
1283
0
            return FALSE;
1284
0
    }
1285
1286
0
    *pdfReqXOff = RoundIfCloseToInt(*pdfReqXOff);
1287
0
    *pdfReqYOff = RoundIfCloseToInt(*pdfReqYOff);
1288
0
    *pdfReqXSize = RoundIfCloseToInt(*pdfReqXSize);
1289
0
    *pdfReqYSize = RoundIfCloseToInt(*pdfReqYSize);
1290
1291
0
    return TRUE;
1292
0
}
1293
1294
/************************************************************************/
1295
/*                          NeedMaxValAdjustment()                      */
1296
/************************************************************************/
1297
1298
int VRTSimpleSource::NeedMaxValAdjustment() const
1299
0
{
1300
0
    if (!m_nMaxValue)
1301
0
        return FALSE;
1302
1303
0
    auto l_band = GetRasterBand();
1304
0
    if (!l_band)
1305
0
        return FALSE;
1306
0
    const char *pszNBITS = l_band->GetMetadataItem("NBITS", "IMAGE_STRUCTURE");
1307
0
    const int nBits = (pszNBITS) ? atoi(pszNBITS) : 0;
1308
0
    if (nBits >= 1 && nBits <= 31)
1309
0
    {
1310
0
        const int nBandMaxValue = static_cast<int>((1U << nBits) - 1);
1311
0
        return nBandMaxValue > m_nMaxValue;
1312
0
    }
1313
0
    return TRUE;
1314
0
}
1315
1316
/************************************************************************/
1317
/*                              RasterIO()                              */
1318
/************************************************************************/
1319
1320
CPLErr VRTSimpleSource::RasterIO(GDALDataType eVRTBandDataType, int nXOff,
1321
                                 int nYOff, int nXSize, int nYSize, void *pData,
1322
                                 int nBufXSize, int nBufYSize,
1323
                                 GDALDataType eBufType, GSpacing nPixelSpace,
1324
                                 GSpacing nLineSpace,
1325
                                 GDALRasterIOExtraArg *psExtraArgIn,
1326
                                 WorkingState & /*oWorkingState*/)
1327
1328
0
{
1329
0
    GDALRasterIOExtraArg sExtraArg;
1330
0
    INIT_RASTERIO_EXTRA_ARG(sExtraArg);
1331
0
    GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
1332
1333
0
    double dfXOff = nXOff;
1334
0
    double dfYOff = nYOff;
1335
0
    double dfXSize = nXSize;
1336
0
    double dfYSize = nYSize;
1337
0
    if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
1338
0
    {
1339
0
        dfXOff = psExtraArgIn->dfXOff;
1340
0
        dfYOff = psExtraArgIn->dfYOff;
1341
0
        dfXSize = psExtraArgIn->dfXSize;
1342
0
        dfYSize = psExtraArgIn->dfYSize;
1343
0
    }
1344
1345
    // The window we will actually request from the source raster band.
1346
0
    double dfReqXOff = 0.0;
1347
0
    double dfReqYOff = 0.0;
1348
0
    double dfReqXSize = 0.0;
1349
0
    double dfReqYSize = 0.0;
1350
0
    int nReqXOff = 0;
1351
0
    int nReqYOff = 0;
1352
0
    int nReqXSize = 0;
1353
0
    int nReqYSize = 0;
1354
1355
    // The window we will actual set _within_ the pData buffer.
1356
0
    int nOutXOff = 0;
1357
0
    int nOutYOff = 0;
1358
0
    int nOutXSize = 0;
1359
0
    int nOutYSize = 0;
1360
1361
0
    bool bError = false;
1362
0
    if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
1363
0
                         &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
1364
0
                         &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
1365
0
                         &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
1366
0
    {
1367
0
        return bError ? CE_Failure : CE_None;
1368
0
    }
1369
#if DEBUG_VERBOSE
1370
    CPLDebug("VRT",
1371
             "nXOff=%d, nYOff=%d, nXSize=%d, nYSize=%d, nBufXSize=%d, "
1372
             "nBufYSize=%d,\n"
1373
             "dfReqXOff=%g, dfReqYOff=%g, dfReqXSize=%g, dfReqYSize=%g,\n"
1374
             "nReqXOff=%d, nReqYOff=%d, nReqXSize=%d, nReqYSize=%d,\n"
1375
             "nOutXOff=%d, nOutYOff=%d, nOutXSize=%d, nOutYSize=%d",
1376
             nXOff, nYOff, nXSize, nYSize, nBufXSize, nBufYSize, dfReqXOff,
1377
             dfReqYOff, dfReqXSize, dfReqYSize, nReqXOff, nReqYOff, nReqXSize,
1378
             nReqYSize, nOutXOff, nOutYOff, nOutXSize, nOutYSize);
1379
#endif
1380
1381
    /* -------------------------------------------------------------------- */
1382
    /*      Actually perform the IO request.                                */
1383
    /* -------------------------------------------------------------------- */
1384
0
    if (!m_osResampling.empty())
1385
0
    {
1386
0
        psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
1387
0
    }
1388
0
    else if (psExtraArgIn != nullptr)
1389
0
    {
1390
0
        psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
1391
0
    }
1392
0
    psExtraArg->bFloatingPointWindowValidity = TRUE;
1393
0
    psExtraArg->dfXOff = dfReqXOff;
1394
0
    psExtraArg->dfYOff = dfReqYOff;
1395
0
    psExtraArg->dfXSize = dfReqXSize;
1396
0
    psExtraArg->dfYSize = dfReqYSize;
1397
0
    if (psExtraArgIn)
1398
0
    {
1399
0
        psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
1400
0
        psExtraArg->pProgressData = psExtraArgIn->pProgressData;
1401
0
        if (psExtraArgIn->nVersion >= 2)
1402
0
        {
1403
0
            psExtraArg->bUseOnlyThisScale = psExtraArgIn->bUseOnlyThisScale;
1404
0
        }
1405
0
    }
1406
1407
0
    GByte *pabyOut = static_cast<unsigned char *>(pData) +
1408
0
                     nOutXOff * nPixelSpace +
1409
0
                     static_cast<GPtrDiff_t>(nOutYOff) * nLineSpace;
1410
1411
0
    auto l_band = GetRasterBand();
1412
0
    if (!l_band)
1413
0
        return CE_Failure;
1414
1415
0
    CPLErr eErr = CE_Failure;
1416
0
    if (GDALDataTypeIsConversionLossy(l_band->GetRasterDataType(),
1417
0
                                      eVRTBandDataType))
1418
0
    {
1419
0
        const int nBandDTSize = GDALGetDataTypeSizeBytes(eVRTBandDataType);
1420
0
        void *pTemp = VSI_MALLOC3_VERBOSE(nOutXSize, nOutYSize, nBandDTSize);
1421
0
        if (pTemp)
1422
0
        {
1423
0
            eErr = l_band->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize,
1424
0
                                    nReqYSize, pTemp, nOutXSize, nOutYSize,
1425
0
                                    eVRTBandDataType, 0, 0, psExtraArg);
1426
0
            if (eErr == CE_None)
1427
0
            {
1428
0
                GByte *pabyTemp = static_cast<GByte *>(pTemp);
1429
0
                for (int iY = 0; iY < nOutYSize; iY++)
1430
0
                {
1431
0
                    GDALCopyWords(
1432
0
                        pabyTemp +
1433
0
                            static_cast<size_t>(iY) * nBandDTSize * nOutXSize,
1434
0
                        eVRTBandDataType, nBandDTSize,
1435
0
                        pabyOut + static_cast<GPtrDiff_t>(iY * nLineSpace),
1436
0
                        eBufType, static_cast<int>(nPixelSpace), nOutXSize);
1437
0
                }
1438
0
            }
1439
0
            VSIFree(pTemp);
1440
0
        }
1441
0
    }
1442
0
    else
1443
0
    {
1444
0
        eErr = l_band->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize,
1445
0
                                nReqYSize, pabyOut, nOutXSize, nOutYSize,
1446
0
                                eBufType, nPixelSpace, nLineSpace, psExtraArg);
1447
0
    }
1448
1449
0
    if (NeedMaxValAdjustment())
1450
0
    {
1451
0
        for (int j = 0; j < nOutYSize; j++)
1452
0
        {
1453
0
            for (int i = 0; i < nOutXSize; i++)
1454
0
            {
1455
0
                int nVal = 0;
1456
0
                GDALCopyWords(pabyOut + j * nLineSpace + i * nPixelSpace,
1457
0
                              eBufType, 0, &nVal, GDT_Int32, 0, 1);
1458
0
                if (nVal > m_nMaxValue)
1459
0
                    nVal = m_nMaxValue;
1460
0
                GDALCopyWords(&nVal, GDT_Int32, 0,
1461
0
                              pabyOut + j * nLineSpace + i * nPixelSpace,
1462
0
                              eBufType, 0, 1);
1463
0
            }
1464
0
        }
1465
0
    }
1466
1467
0
    if (psExtraArg->pfnProgress)
1468
0
        psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
1469
1470
0
    return eErr;
1471
0
}
1472
1473
/************************************************************************/
1474
/*                             GetMinimum()                             */
1475
/************************************************************************/
1476
1477
double VRTSimpleSource::GetMinimum(int nXSize, int nYSize, int *pbSuccess)
1478
0
{
1479
    // The window we will actually request from the source raster band.
1480
0
    double dfReqXOff = 0.0;
1481
0
    double dfReqYOff = 0.0;
1482
0
    double dfReqXSize = 0.0;
1483
0
    double dfReqYSize = 0.0;
1484
0
    int nReqXOff = 0;
1485
0
    int nReqYOff = 0;
1486
0
    int nReqXSize = 0;
1487
0
    int nReqYSize = 0;
1488
1489
    // The window we will actual set _within_ the pData buffer.
1490
0
    int nOutXOff = 0;
1491
0
    int nOutYOff = 0;
1492
0
    int nOutXSize = 0;
1493
0
    int nOutYSize = 0;
1494
1495
0
    bool bError = false;
1496
0
    auto l_band = GetRasterBand();
1497
0
    if (!l_band ||
1498
0
        !GetSrcDstWindow(0, 0, nXSize, nYSize, nXSize, nYSize, &dfReqXOff,
1499
0
                         &dfReqYOff, &dfReqXSize, &dfReqYSize, &nReqXOff,
1500
0
                         &nReqYOff, &nReqXSize, &nReqYSize, &nOutXOff,
1501
0
                         &nOutYOff, &nOutXSize, &nOutYSize, bError) ||
1502
0
        nReqXOff != 0 || nReqYOff != 0 || nReqXSize != l_band->GetXSize() ||
1503
0
        nReqYSize != l_band->GetYSize())
1504
0
    {
1505
0
        *pbSuccess = FALSE;
1506
0
        return 0;
1507
0
    }
1508
1509
0
    const double dfVal = l_band->GetMinimum(pbSuccess);
1510
0
    if (NeedMaxValAdjustment() && dfVal > m_nMaxValue)
1511
0
        return m_nMaxValue;
1512
0
    return dfVal;
1513
0
}
1514
1515
/************************************************************************/
1516
/*                             GetMaximum()                             */
1517
/************************************************************************/
1518
1519
double VRTSimpleSource::GetMaximum(int nXSize, int nYSize, int *pbSuccess)
1520
0
{
1521
    // The window we will actually request from the source raster band.
1522
0
    double dfReqXOff = 0.0;
1523
0
    double dfReqYOff = 0.0;
1524
0
    double dfReqXSize = 0.0;
1525
0
    double dfReqYSize = 0.0;
1526
0
    int nReqXOff = 0;
1527
0
    int nReqYOff = 0;
1528
0
    int nReqXSize = 0;
1529
0
    int nReqYSize = 0;
1530
1531
    // The window we will actual set _within_ the pData buffer.
1532
0
    int nOutXOff = 0;
1533
0
    int nOutYOff = 0;
1534
0
    int nOutXSize = 0;
1535
0
    int nOutYSize = 0;
1536
1537
0
    bool bError = false;
1538
0
    auto l_band = GetRasterBand();
1539
0
    if (!l_band ||
1540
0
        !GetSrcDstWindow(0, 0, nXSize, nYSize, nXSize, nYSize, &dfReqXOff,
1541
0
                         &dfReqYOff, &dfReqXSize, &dfReqYSize, &nReqXOff,
1542
0
                         &nReqYOff, &nReqXSize, &nReqYSize, &nOutXOff,
1543
0
                         &nOutYOff, &nOutXSize, &nOutYSize, bError) ||
1544
0
        nReqXOff != 0 || nReqYOff != 0 || nReqXSize != l_band->GetXSize() ||
1545
0
        nReqYSize != l_band->GetYSize())
1546
0
    {
1547
0
        *pbSuccess = FALSE;
1548
0
        return 0;
1549
0
    }
1550
1551
0
    const double dfVal = l_band->GetMaximum(pbSuccess);
1552
0
    if (NeedMaxValAdjustment() && dfVal > m_nMaxValue)
1553
0
        return m_nMaxValue;
1554
0
    return dfVal;
1555
0
}
1556
1557
/************************************************************************/
1558
/*                            GetHistogram()                            */
1559
/************************************************************************/
1560
1561
CPLErr VRTSimpleSource::GetHistogram(int nXSize, int nYSize, double dfMin,
1562
                                     double dfMax, int nBuckets,
1563
                                     GUIntBig *panHistogram,
1564
                                     int bIncludeOutOfRange, int bApproxOK,
1565
                                     GDALProgressFunc pfnProgress,
1566
                                     void *pProgressData)
1567
0
{
1568
    // The window we will actually request from the source raster band.
1569
0
    double dfReqXOff = 0.0;
1570
0
    double dfReqYOff = 0.0;
1571
0
    double dfReqXSize = 0.0;
1572
0
    double dfReqYSize = 0.0;
1573
0
    int nReqXOff = 0;
1574
0
    int nReqYOff = 0;
1575
0
    int nReqXSize = 0;
1576
0
    int nReqYSize = 0;
1577
1578
    // The window we will actual set _within_ the pData buffer.
1579
0
    int nOutXOff = 0;
1580
0
    int nOutYOff = 0;
1581
0
    int nOutXSize = 0;
1582
0
    int nOutYSize = 0;
1583
1584
0
    bool bError = false;
1585
0
    auto l_band = GetRasterBand();
1586
0
    if (!l_band || NeedMaxValAdjustment() ||
1587
0
        !GetSrcDstWindow(0, 0, nXSize, nYSize, nXSize, nYSize, &dfReqXOff,
1588
0
                         &dfReqYOff, &dfReqXSize, &dfReqYSize, &nReqXOff,
1589
0
                         &nReqYOff, &nReqXSize, &nReqYSize, &nOutXOff,
1590
0
                         &nOutYOff, &nOutXSize, &nOutYSize, bError) ||
1591
0
        nReqXOff != 0 || nReqYOff != 0 || nReqXSize != l_band->GetXSize() ||
1592
0
        nReqYSize != l_band->GetYSize())
1593
0
    {
1594
0
        return CE_Failure;
1595
0
    }
1596
1597
0
    return l_band->GetHistogram(dfMin, dfMax, nBuckets, panHistogram,
1598
0
                                bIncludeOutOfRange, bApproxOK, pfnProgress,
1599
0
                                pProgressData);
1600
0
}
1601
1602
/************************************************************************/
1603
/*                          DatasetRasterIO()                           */
1604
/************************************************************************/
1605
1606
CPLErr VRTSimpleSource::DatasetRasterIO(
1607
    GDALDataType eVRTBandDataType, int nXOff, int nYOff, int nXSize, int nYSize,
1608
    void *pData, int nBufXSize, int nBufYSize, GDALDataType eBufType,
1609
    int nBandCount, const int *panBandMap, GSpacing nPixelSpace,
1610
    GSpacing nLineSpace, GSpacing nBandSpace,
1611
    GDALRasterIOExtraArg *psExtraArgIn)
1612
0
{
1613
0
    if (GetType() != VRTSimpleSource::GetTypeStatic())
1614
0
    {
1615
0
        CPLError(CE_Failure, CPLE_NotSupported,
1616
0
                 "DatasetRasterIO() not implemented for %s", GetType());
1617
0
        return CE_Failure;
1618
0
    }
1619
1620
0
    GDALRasterIOExtraArg sExtraArg;
1621
0
    INIT_RASTERIO_EXTRA_ARG(sExtraArg);
1622
0
    GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
1623
1624
0
    double dfXOff = nXOff;
1625
0
    double dfYOff = nYOff;
1626
0
    double dfXSize = nXSize;
1627
0
    double dfYSize = nYSize;
1628
0
    if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
1629
0
    {
1630
0
        dfXOff = psExtraArgIn->dfXOff;
1631
0
        dfYOff = psExtraArgIn->dfYOff;
1632
0
        dfXSize = psExtraArgIn->dfXSize;
1633
0
        dfYSize = psExtraArgIn->dfYSize;
1634
0
    }
1635
1636
    // The window we will actually request from the source raster band.
1637
0
    double dfReqXOff = 0.0;
1638
0
    double dfReqYOff = 0.0;
1639
0
    double dfReqXSize = 0.0;
1640
0
    double dfReqYSize = 0.0;
1641
0
    int nReqXOff = 0;
1642
0
    int nReqYOff = 0;
1643
0
    int nReqXSize = 0;
1644
0
    int nReqYSize = 0;
1645
1646
    // The window we will actual set _within_ the pData buffer.
1647
0
    int nOutXOff = 0;
1648
0
    int nOutYOff = 0;
1649
0
    int nOutXSize = 0;
1650
0
    int nOutYSize = 0;
1651
1652
0
    bool bError = false;
1653
0
    if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
1654
0
                         &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
1655
0
                         &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
1656
0
                         &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
1657
0
    {
1658
0
        return bError ? CE_Failure : CE_None;
1659
0
    }
1660
1661
0
    auto l_band = GetRasterBand();
1662
0
    if (!l_band)
1663
0
        return CE_Failure;
1664
1665
0
    GDALDataset *poDS = l_band->GetDataset();
1666
0
    if (poDS == nullptr)
1667
0
        return CE_Failure;
1668
1669
0
    if (!m_osResampling.empty())
1670
0
    {
1671
0
        psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
1672
0
    }
1673
0
    else if (psExtraArgIn != nullptr)
1674
0
    {
1675
0
        psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
1676
0
    }
1677
0
    psExtraArg->bFloatingPointWindowValidity = TRUE;
1678
0
    psExtraArg->dfXOff = dfReqXOff;
1679
0
    psExtraArg->dfYOff = dfReqYOff;
1680
0
    psExtraArg->dfXSize = dfReqXSize;
1681
0
    psExtraArg->dfYSize = dfReqYSize;
1682
0
    if (psExtraArgIn)
1683
0
    {
1684
0
        psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
1685
0
        psExtraArg->pProgressData = psExtraArgIn->pProgressData;
1686
0
    }
1687
1688
0
    GByte *pabyOut = static_cast<unsigned char *>(pData) +
1689
0
                     nOutXOff * nPixelSpace +
1690
0
                     static_cast<GPtrDiff_t>(nOutYOff) * nLineSpace;
1691
1692
0
    CPLErr eErr = CE_Failure;
1693
1694
0
    if (GDALDataTypeIsConversionLossy(l_band->GetRasterDataType(),
1695
0
                                      eVRTBandDataType))
1696
0
    {
1697
0
        const int nBandDTSize = GDALGetDataTypeSizeBytes(eVRTBandDataType);
1698
0
        void *pTemp = VSI_MALLOC3_VERBOSE(
1699
0
            nOutXSize, nOutYSize, cpl::fits_on<int>(nBandDTSize * nBandCount));
1700
0
        if (pTemp)
1701
0
        {
1702
0
            eErr = poDS->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize,
1703
0
                                  nReqYSize, pTemp, nOutXSize, nOutYSize,
1704
0
                                  eVRTBandDataType, nBandCount, panBandMap, 0,
1705
0
                                  0, 0, psExtraArg);
1706
0
            if (eErr == CE_None)
1707
0
            {
1708
0
                GByte *pabyTemp = static_cast<GByte *>(pTemp);
1709
0
                const size_t nSrcBandSpace =
1710
0
                    static_cast<size_t>(nOutYSize) * nOutXSize * nBandDTSize;
1711
0
                for (int iBand = 0; iBand < nBandCount; iBand++)
1712
0
                {
1713
0
                    for (int iY = 0; iY < nOutYSize; iY++)
1714
0
                    {
1715
0
                        GDALCopyWords(
1716
0
                            pabyTemp + iBand * nSrcBandSpace +
1717
0
                                static_cast<size_t>(iY) * nBandDTSize *
1718
0
                                    nOutXSize,
1719
0
                            eVRTBandDataType, nBandDTSize,
1720
0
                            pabyOut + static_cast<GPtrDiff_t>(
1721
0
                                          iY * nLineSpace + iBand * nBandSpace),
1722
0
                            eBufType, static_cast<int>(nPixelSpace), nOutXSize);
1723
0
                    }
1724
0
                }
1725
0
            }
1726
0
            VSIFree(pTemp);
1727
0
        }
1728
0
    }
1729
0
    else
1730
0
    {
1731
0
        eErr = poDS->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
1732
0
                              pabyOut, nOutXSize, nOutYSize, eBufType,
1733
0
                              nBandCount, panBandMap, nPixelSpace, nLineSpace,
1734
0
                              nBandSpace, psExtraArg);
1735
0
    }
1736
1737
0
    if (NeedMaxValAdjustment())
1738
0
    {
1739
0
        for (int k = 0; k < nBandCount; k++)
1740
0
        {
1741
0
            for (int j = 0; j < nOutYSize; j++)
1742
0
            {
1743
0
                for (int i = 0; i < nOutXSize; i++)
1744
0
                {
1745
0
                    int nVal = 0;
1746
0
                    GDALCopyWords(pabyOut + k * nBandSpace + j * nLineSpace +
1747
0
                                      i * nPixelSpace,
1748
0
                                  eBufType, 0, &nVal, GDT_Int32, 0, 1);
1749
1750
0
                    if (nVal > m_nMaxValue)
1751
0
                        nVal = m_nMaxValue;
1752
1753
0
                    GDALCopyWords(&nVal, GDT_Int32, 0,
1754
0
                                  pabyOut + k * nBandSpace + j * nLineSpace +
1755
0
                                      i * nPixelSpace,
1756
0
                                  eBufType, 0, 1);
1757
0
                }
1758
0
            }
1759
0
        }
1760
0
    }
1761
1762
0
    if (psExtraArg->pfnProgress)
1763
0
        psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
1764
1765
0
    return eErr;
1766
0
}
1767
1768
/************************************************************************/
1769
/*                          SetResampling()                             */
1770
/************************************************************************/
1771
1772
void VRTSimpleSource::SetResampling(const char *pszResampling)
1773
0
{
1774
0
    m_osResampling = (pszResampling) ? pszResampling : "";
1775
0
}
1776
1777
/************************************************************************/
1778
/* ==================================================================== */
1779
/*                         VRTAveragedSource                            */
1780
/* ==================================================================== */
1781
/************************************************************************/
1782
1783
/************************************************************************/
1784
/*                         VRTAveragedSource()                          */
1785
/************************************************************************/
1786
1787
VRTAveragedSource::VRTAveragedSource()
1788
0
{
1789
0
}
1790
1791
/************************************************************************/
1792
/*                           GetTypeStatic()                            */
1793
/************************************************************************/
1794
1795
const char *VRTAveragedSource::GetTypeStatic()
1796
0
{
1797
0
    static const char *TYPE = "AveragedSource";
1798
0
    return TYPE;
1799
0
}
1800
1801
/************************************************************************/
1802
/*                            GetType()                                 */
1803
/************************************************************************/
1804
1805
const char *VRTAveragedSource::GetType() const
1806
0
{
1807
0
    return GetTypeStatic();
1808
0
}
1809
1810
/************************************************************************/
1811
/*                           SerializeToXML()                           */
1812
/************************************************************************/
1813
1814
CPLXMLNode *VRTAveragedSource::SerializeToXML(const char *pszVRTPath)
1815
1816
0
{
1817
0
    CPLXMLNode *const psSrc = VRTSimpleSource::SerializeToXML(pszVRTPath);
1818
1819
0
    if (psSrc == nullptr)
1820
0
        return nullptr;
1821
1822
0
    CPLFree(psSrc->pszValue);
1823
0
    psSrc->pszValue = CPLStrdup(GetTypeStatic());
1824
1825
0
    return psSrc;
1826
0
}
1827
1828
/************************************************************************/
1829
/*                           SetNoDataValue()                           */
1830
/************************************************************************/
1831
1832
void VRTAveragedSource::SetNoDataValue(double dfNewNoDataValue)
1833
1834
0
{
1835
0
    if (dfNewNoDataValue == VRT_NODATA_UNSET)
1836
0
    {
1837
0
        m_bNoDataSet = FALSE;
1838
0
        m_dfNoDataValue = VRT_NODATA_UNSET;
1839
0
        return;
1840
0
    }
1841
1842
0
    m_bNoDataSet = TRUE;
1843
0
    m_dfNoDataValue = dfNewNoDataValue;
1844
0
}
1845
1846
/************************************************************************/
1847
/*                              RasterIO()                              */
1848
/************************************************************************/
1849
1850
CPLErr VRTAveragedSource::RasterIO(GDALDataType /*eVRTBandDataType*/, int nXOff,
1851
                                   int nYOff, int nXSize, int nYSize,
1852
                                   void *pData, int nBufXSize, int nBufYSize,
1853
                                   GDALDataType eBufType, GSpacing nPixelSpace,
1854
                                   GSpacing nLineSpace,
1855
                                   GDALRasterIOExtraArg *psExtraArgIn,
1856
                                   WorkingState & /*oWorkingState*/)
1857
1858
0
{
1859
0
    GDALRasterIOExtraArg sExtraArg;
1860
0
    INIT_RASTERIO_EXTRA_ARG(sExtraArg);
1861
0
    GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
1862
1863
0
    double dfXOff = nXOff;
1864
0
    double dfYOff = nYOff;
1865
0
    double dfXSize = nXSize;
1866
0
    double dfYSize = nYSize;
1867
0
    if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
1868
0
    {
1869
0
        dfXOff = psExtraArgIn->dfXOff;
1870
0
        dfYOff = psExtraArgIn->dfYOff;
1871
0
        dfXSize = psExtraArgIn->dfXSize;
1872
0
        dfYSize = psExtraArgIn->dfYSize;
1873
0
    }
1874
1875
    // The window we will actually request from the source raster band.
1876
0
    double dfReqXOff = 0.0;
1877
0
    double dfReqYOff = 0.0;
1878
0
    double dfReqXSize = 0.0;
1879
0
    double dfReqYSize = 0.0;
1880
0
    int nReqXOff = 0;
1881
0
    int nReqYOff = 0;
1882
0
    int nReqXSize = 0;
1883
0
    int nReqYSize = 0;
1884
1885
    // The window we will actual set _within_ the pData buffer.
1886
0
    int nOutXOff = 0;
1887
0
    int nOutYOff = 0;
1888
0
    int nOutXSize = 0;
1889
0
    int nOutYSize = 0;
1890
1891
0
    bool bError = false;
1892
0
    if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
1893
0
                         &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
1894
0
                         &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
1895
0
                         &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
1896
0
    {
1897
0
        return bError ? CE_Failure : CE_None;
1898
0
    }
1899
1900
0
    auto l_band = GetRasterBand();
1901
0
    if (!l_band)
1902
0
        return CE_Failure;
1903
1904
    /* -------------------------------------------------------------------- */
1905
    /*      Allocate a temporary buffer to whole the full resolution        */
1906
    /*      data from the area of interest.                                 */
1907
    /* -------------------------------------------------------------------- */
1908
0
    float *const pafSrc = static_cast<float *>(
1909
0
        VSI_MALLOC3_VERBOSE(sizeof(float), nReqXSize, nReqYSize));
1910
0
    if (pafSrc == nullptr)
1911
0
    {
1912
0
        return CE_Failure;
1913
0
    }
1914
1915
    /* -------------------------------------------------------------------- */
1916
    /*      Load it.                                                        */
1917
    /* -------------------------------------------------------------------- */
1918
0
    if (!m_osResampling.empty())
1919
0
    {
1920
0
        psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
1921
0
    }
1922
0
    else if (psExtraArgIn != nullptr)
1923
0
    {
1924
0
        psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
1925
0
    }
1926
1927
0
    psExtraArg->bFloatingPointWindowValidity = TRUE;
1928
0
    psExtraArg->dfXOff = dfReqXOff;
1929
0
    psExtraArg->dfYOff = dfReqYOff;
1930
0
    psExtraArg->dfXSize = dfReqXSize;
1931
0
    psExtraArg->dfYSize = dfReqYSize;
1932
0
    if (psExtraArgIn)
1933
0
    {
1934
0
        psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
1935
0
        psExtraArg->pProgressData = psExtraArgIn->pProgressData;
1936
0
    }
1937
1938
0
    const CPLErr eErr = l_band->RasterIO(
1939
0
        GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize, pafSrc, nReqXSize,
1940
0
        nReqYSize, GDT_Float32, 0, 0, psExtraArg);
1941
1942
0
    if (eErr != CE_None)
1943
0
    {
1944
0
        VSIFree(pafSrc);
1945
0
        return eErr;
1946
0
    }
1947
1948
    /* -------------------------------------------------------------------- */
1949
    /*      Do the averaging.                                               */
1950
    /* -------------------------------------------------------------------- */
1951
0
    for (int iBufLine = nOutYOff; iBufLine < nOutYOff + nOutYSize; iBufLine++)
1952
0
    {
1953
0
        const double dfYDst =
1954
0
            (iBufLine / static_cast<double>(nBufYSize)) * nYSize + nYOff;
1955
1956
0
        for (int iBufPixel = nOutXOff; iBufPixel < nOutXOff + nOutXSize;
1957
0
             iBufPixel++)
1958
0
        {
1959
0
            double dfXSrcStart, dfXSrcEnd, dfYSrcStart, dfYSrcEnd;
1960
0
            int iXSrcStart, iYSrcStart, iXSrcEnd, iYSrcEnd;
1961
1962
0
            const double dfXDst =
1963
0
                (iBufPixel / static_cast<double>(nBufXSize)) * nXSize + nXOff;
1964
1965
            // Compute the source image rectangle needed for this pixel.
1966
0
            DstToSrc(dfXDst, dfYDst, dfXSrcStart, dfYSrcStart);
1967
0
            DstToSrc(dfXDst + 1.0, dfYDst + 1.0, dfXSrcEnd, dfYSrcEnd);
1968
1969
            // Convert to integers, assuming that the center of the source
1970
            // pixel must be in our rect to get included.
1971
0
            if (dfXSrcEnd >= dfXSrcStart + 1)
1972
0
            {
1973
0
                iXSrcStart = static_cast<int>(floor(dfXSrcStart + 0.5));
1974
0
                iXSrcEnd = static_cast<int>(floor(dfXSrcEnd + 0.5));
1975
0
            }
1976
0
            else
1977
0
            {
1978
                /* If the resampling factor is less than 100%, the distance */
1979
                /* between the source pixel is < 1, so we stick to nearest */
1980
                /* neighbour */
1981
0
                iXSrcStart = static_cast<int>(floor(dfXSrcStart));
1982
0
                iXSrcEnd = iXSrcStart + 1;
1983
0
            }
1984
0
            if (dfYSrcEnd >= dfYSrcStart + 1)
1985
0
            {
1986
0
                iYSrcStart = static_cast<int>(floor(dfYSrcStart + 0.5));
1987
0
                iYSrcEnd = static_cast<int>(floor(dfYSrcEnd + 0.5));
1988
0
            }
1989
0
            else
1990
0
            {
1991
0
                iYSrcStart = static_cast<int>(floor(dfYSrcStart));
1992
0
                iYSrcEnd = iYSrcStart + 1;
1993
0
            }
1994
1995
            // Transform into the coordinate system of the source *buffer*
1996
0
            iXSrcStart -= nReqXOff;
1997
0
            iYSrcStart -= nReqYOff;
1998
0
            iXSrcEnd -= nReqXOff;
1999
0
            iYSrcEnd -= nReqYOff;
2000
2001
0
            double dfSum = 0.0;
2002
0
            int nPixelCount = 0;
2003
2004
0
            for (int iY = iYSrcStart; iY < iYSrcEnd; iY++)
2005
0
            {
2006
0
                if (iY < 0 || iY >= nReqYSize)
2007
0
                    continue;
2008
2009
0
                for (int iX = iXSrcStart; iX < iXSrcEnd; iX++)
2010
0
                {
2011
0
                    if (iX < 0 || iX >= nReqXSize)
2012
0
                        continue;
2013
2014
0
                    const float fSampledValue =
2015
0
                        pafSrc[iX + static_cast<size_t>(iY) * nReqXSize];
2016
0
                    if (std::isnan(fSampledValue))
2017
0
                        continue;
2018
2019
0
                    if (m_bNoDataSet &&
2020
0
                        GDALIsValueInRange<float>(m_dfNoDataValue) &&
2021
0
                        ARE_REAL_EQUAL(fSampledValue,
2022
0
                                       static_cast<float>(m_dfNoDataValue)))
2023
0
                        continue;
2024
2025
0
                    nPixelCount++;
2026
0
                    dfSum += pafSrc[iX + static_cast<size_t>(iY) * nReqXSize];
2027
0
                }
2028
0
            }
2029
2030
0
            if (nPixelCount == 0)
2031
0
                continue;
2032
2033
            // Compute output value.
2034
0
            const float dfOutputValue = static_cast<float>(dfSum / nPixelCount);
2035
2036
            // Put it in the output buffer.
2037
0
            GByte *pDstLocation =
2038
0
                static_cast<GByte *>(pData) + nPixelSpace * iBufPixel +
2039
0
                static_cast<GPtrDiff_t>(nLineSpace) * iBufLine;
2040
2041
0
            if (eBufType == GDT_Byte)
2042
0
                *pDstLocation = static_cast<GByte>(
2043
0
                    std::min(255.0, std::max(0.0, dfOutputValue + 0.5)));
2044
0
            else
2045
0
                GDALCopyWords(&dfOutputValue, GDT_Float32, 4, pDstLocation,
2046
0
                              eBufType, 8, 1);
2047
0
        }
2048
0
    }
2049
2050
0
    VSIFree(pafSrc);
2051
2052
0
    if (psExtraArg->pfnProgress)
2053
0
        psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
2054
2055
0
    return CE_None;
2056
0
}
2057
2058
/************************************************************************/
2059
/*                             GetMinimum()                             */
2060
/************************************************************************/
2061
2062
double VRTAveragedSource::GetMinimum(int /* nXSize */, int /* nYSize */,
2063
                                     int *pbSuccess)
2064
0
{
2065
0
    *pbSuccess = FALSE;
2066
0
    return 0.0;
2067
0
}
2068
2069
/************************************************************************/
2070
/*                             GetMaximum()                             */
2071
/************************************************************************/
2072
2073
double VRTAveragedSource::GetMaximum(int /* nXSize */, int /* nYSize */,
2074
                                     int *pbSuccess)
2075
0
{
2076
0
    *pbSuccess = FALSE;
2077
0
    return 0.0;
2078
0
}
2079
2080
/************************************************************************/
2081
/*                            GetHistogram()                            */
2082
/************************************************************************/
2083
2084
CPLErr VRTAveragedSource::GetHistogram(
2085
    int /* nXSize */, int /* nYSize */, double /* dfMin */, double /* dfMax */,
2086
    int /* nBuckets */, GUIntBig * /* panHistogram */,
2087
    int /* bIncludeOutOfRange */, int /* bApproxOK */,
2088
    GDALProgressFunc /* pfnProgress */, void * /* pProgressData */)
2089
0
{
2090
0
    return CE_Failure;
2091
0
}
2092
2093
/************************************************************************/
2094
/* ==================================================================== */
2095
/*                     VRTNoDataFromMaskSource                          */
2096
/* ==================================================================== */
2097
/************************************************************************/
2098
2099
/************************************************************************/
2100
/*                     VRTNoDataFromMaskSource()                        */
2101
/************************************************************************/
2102
2103
VRTNoDataFromMaskSource::VRTNoDataFromMaskSource()
2104
0
{
2105
0
}
2106
2107
/************************************************************************/
2108
/*                              XMLInit()                               */
2109
/************************************************************************/
2110
2111
CPLErr
2112
VRTNoDataFromMaskSource::XMLInit(const CPLXMLNode *psSrc,
2113
                                 const char *pszVRTPath,
2114
                                 VRTMapSharedResources &oMapSharedSources)
2115
2116
0
{
2117
    /* -------------------------------------------------------------------- */
2118
    /*      Do base initialization.                                         */
2119
    /* -------------------------------------------------------------------- */
2120
0
    {
2121
0
        const CPLErr eErr =
2122
0
            VRTSimpleSource::XMLInit(psSrc, pszVRTPath, oMapSharedSources);
2123
0
        if (eErr != CE_None)
2124
0
            return eErr;
2125
0
    }
2126
2127
0
    if (const char *pszNODATA = CPLGetXMLValue(psSrc, "NODATA", nullptr))
2128
0
    {
2129
0
        m_bNoDataSet = true;
2130
0
        m_dfNoDataValue = CPLAtofM(pszNODATA);
2131
0
    }
2132
2133
0
    m_dfMaskValueThreshold =
2134
0
        CPLAtofM(CPLGetXMLValue(psSrc, "MaskValueThreshold", "0"));
2135
2136
0
    if (const char *pszRemappedValue =
2137
0
            CPLGetXMLValue(psSrc, "RemappedValue", nullptr))
2138
0
    {
2139
0
        m_bHasRemappedValue = true;
2140
0
        m_dfRemappedValue = CPLAtofM(pszRemappedValue);
2141
0
    }
2142
2143
0
    return CE_None;
2144
0
}
2145
2146
/************************************************************************/
2147
/*                           GetTypeStatic()                            */
2148
/************************************************************************/
2149
2150
const char *VRTNoDataFromMaskSource::GetTypeStatic()
2151
0
{
2152
0
    static const char *TYPE = "NoDataFromMaskSource";
2153
0
    return TYPE;
2154
0
}
2155
2156
/************************************************************************/
2157
/*                            GetType()                                 */
2158
/************************************************************************/
2159
2160
const char *VRTNoDataFromMaskSource::GetType() const
2161
0
{
2162
0
    return GetTypeStatic();
2163
0
}
2164
2165
/************************************************************************/
2166
/*                           SerializeToXML()                           */
2167
/************************************************************************/
2168
2169
CPLXMLNode *VRTNoDataFromMaskSource::SerializeToXML(const char *pszVRTPath)
2170
2171
0
{
2172
0
    CPLXMLNode *const psSrc = VRTSimpleSource::SerializeToXML(pszVRTPath);
2173
2174
0
    if (psSrc == nullptr)
2175
0
        return nullptr;
2176
2177
0
    CPLFree(psSrc->pszValue);
2178
0
    psSrc->pszValue = CPLStrdup(GetTypeStatic());
2179
2180
0
    if (m_bNoDataSet)
2181
0
    {
2182
0
        CPLSetXMLValue(psSrc, "MaskValueThreshold",
2183
0
                       CPLSPrintf("%.17g", m_dfMaskValueThreshold));
2184
2185
0
        GDALDataType eBandDT = GDT_Unknown;
2186
0
        double dfNoDataValue = m_dfNoDataValue;
2187
0
        const auto kMaxFloat = std::numeric_limits<float>::max();
2188
0
        if (std::fabs(std::fabs(m_dfNoDataValue) - kMaxFloat) <
2189
0
            1e-10 * kMaxFloat)
2190
0
        {
2191
0
            auto l_band = GetRasterBand();
2192
0
            if (l_band)
2193
0
            {
2194
0
                eBandDT = l_band->GetRasterDataType();
2195
0
                if (eBandDT == GDT_Float32)
2196
0
                {
2197
0
                    dfNoDataValue =
2198
0
                        GDALAdjustNoDataCloseToFloatMax(m_dfNoDataValue);
2199
0
                }
2200
0
            }
2201
0
        }
2202
0
        CPLSetXMLValue(psSrc, "NODATA",
2203
0
                       VRTSerializeNoData(dfNoDataValue, eBandDT, 18).c_str());
2204
0
    }
2205
2206
0
    if (m_bHasRemappedValue)
2207
0
    {
2208
0
        CPLSetXMLValue(psSrc, "RemappedValue",
2209
0
                       CPLSPrintf("%.17g", m_dfRemappedValue));
2210
0
    }
2211
2212
0
    return psSrc;
2213
0
}
2214
2215
/************************************************************************/
2216
/*                           SetParameters()                            */
2217
/************************************************************************/
2218
2219
void VRTNoDataFromMaskSource::SetParameters(double dfNoDataValue,
2220
                                            double dfMaskValueThreshold)
2221
0
{
2222
0
    m_bNoDataSet = true;
2223
0
    m_dfNoDataValue = dfNoDataValue;
2224
0
    m_dfMaskValueThreshold = dfMaskValueThreshold;
2225
0
    if (!m_bHasRemappedValue)
2226
0
        m_dfRemappedValue = m_dfNoDataValue;
2227
0
}
2228
2229
/************************************************************************/
2230
/*                           SetParameters()                            */
2231
/************************************************************************/
2232
2233
void VRTNoDataFromMaskSource::SetParameters(double dfNoDataValue,
2234
                                            double dfMaskValueThreshold,
2235
                                            double dfRemappedValue)
2236
0
{
2237
0
    SetParameters(dfNoDataValue, dfMaskValueThreshold);
2238
0
    m_bHasRemappedValue = true;
2239
0
    m_dfRemappedValue = dfRemappedValue;
2240
0
}
2241
2242
/************************************************************************/
2243
/*                              RasterIO()                              */
2244
/************************************************************************/
2245
2246
CPLErr VRTNoDataFromMaskSource::RasterIO(
2247
    GDALDataType eVRTBandDataType, int nXOff, int nYOff, int nXSize, int nYSize,
2248
    void *pData, int nBufXSize, int nBufYSize, GDALDataType eBufType,
2249
    GSpacing nPixelSpace, GSpacing nLineSpace,
2250
    GDALRasterIOExtraArg *psExtraArgIn, WorkingState &oWorkingState)
2251
2252
0
{
2253
0
    if (!m_bNoDataSet)
2254
0
    {
2255
0
        return VRTSimpleSource::RasterIO(eVRTBandDataType, nXOff, nYOff, nXSize,
2256
0
                                         nYSize, pData, nBufXSize, nBufYSize,
2257
0
                                         eBufType, nPixelSpace, nLineSpace,
2258
0
                                         psExtraArgIn, oWorkingState);
2259
0
    }
2260
2261
0
    GDALRasterIOExtraArg sExtraArg;
2262
0
    INIT_RASTERIO_EXTRA_ARG(sExtraArg);
2263
0
    GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
2264
2265
0
    double dfXOff = nXOff;
2266
0
    double dfYOff = nYOff;
2267
0
    double dfXSize = nXSize;
2268
0
    double dfYSize = nYSize;
2269
0
    if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
2270
0
    {
2271
0
        dfXOff = psExtraArgIn->dfXOff;
2272
0
        dfYOff = psExtraArgIn->dfYOff;
2273
0
        dfXSize = psExtraArgIn->dfXSize;
2274
0
        dfYSize = psExtraArgIn->dfYSize;
2275
0
    }
2276
2277
    // The window we will actually request from the source raster band.
2278
0
    double dfReqXOff = 0.0;
2279
0
    double dfReqYOff = 0.0;
2280
0
    double dfReqXSize = 0.0;
2281
0
    double dfReqYSize = 0.0;
2282
0
    int nReqXOff = 0;
2283
0
    int nReqYOff = 0;
2284
0
    int nReqXSize = 0;
2285
0
    int nReqYSize = 0;
2286
2287
    // The window we will actual set _within_ the pData buffer.
2288
0
    int nOutXOff = 0;
2289
0
    int nOutYOff = 0;
2290
0
    int nOutXSize = 0;
2291
0
    int nOutYSize = 0;
2292
2293
0
    bool bError = false;
2294
0
    if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
2295
0
                         &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
2296
0
                         &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
2297
0
                         &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
2298
0
    {
2299
0
        return bError ? CE_Failure : CE_None;
2300
0
    }
2301
2302
0
    auto l_band = GetRasterBand();
2303
0
    if (!l_band)
2304
0
        return CE_Failure;
2305
2306
    /* -------------------------------------------------------------------- */
2307
    /*      Allocate temporary buffer(s).                                   */
2308
    /* -------------------------------------------------------------------- */
2309
0
    const auto eSrcBandDT = l_band->GetRasterDataType();
2310
0
    const int nSrcBandDTSize = GDALGetDataTypeSizeBytes(eSrcBandDT);
2311
0
    const auto eSrcMaskBandDT = l_band->GetMaskBand()->GetRasterDataType();
2312
0
    const int nSrcMaskBandDTSize = GDALGetDataTypeSizeBytes(eSrcMaskBandDT);
2313
0
    double dfRemappedValue = m_dfRemappedValue;
2314
0
    if (!m_bHasRemappedValue)
2315
0
    {
2316
0
        if (eSrcBandDT == GDT_Byte &&
2317
0
            m_dfNoDataValue >= std::numeric_limits<GByte>::min() &&
2318
0
            m_dfNoDataValue <= std::numeric_limits<GByte>::max() &&
2319
0
            static_cast<int>(m_dfNoDataValue) == m_dfNoDataValue)
2320
0
        {
2321
0
            if (m_dfNoDataValue == std::numeric_limits<GByte>::max())
2322
0
                dfRemappedValue = m_dfNoDataValue - 1;
2323
0
            else
2324
0
                dfRemappedValue = m_dfNoDataValue + 1;
2325
0
        }
2326
0
        else if (eSrcBandDT == GDT_UInt16 &&
2327
0
                 m_dfNoDataValue >= std::numeric_limits<uint16_t>::min() &&
2328
0
                 m_dfNoDataValue <= std::numeric_limits<uint16_t>::max() &&
2329
0
                 static_cast<int>(m_dfNoDataValue) == m_dfNoDataValue)
2330
0
        {
2331
0
            if (m_dfNoDataValue == std::numeric_limits<uint16_t>::max())
2332
0
                dfRemappedValue = m_dfNoDataValue - 1;
2333
0
            else
2334
0
                dfRemappedValue = m_dfNoDataValue + 1;
2335
0
        }
2336
0
        else if (eSrcBandDT == GDT_Int16 &&
2337
0
                 m_dfNoDataValue >= std::numeric_limits<int16_t>::min() &&
2338
0
                 m_dfNoDataValue <= std::numeric_limits<int16_t>::max() &&
2339
0
                 static_cast<int>(m_dfNoDataValue) == m_dfNoDataValue)
2340
0
        {
2341
0
            if (m_dfNoDataValue == std::numeric_limits<int16_t>::max())
2342
0
                dfRemappedValue = m_dfNoDataValue - 1;
2343
0
            else
2344
0
                dfRemappedValue = m_dfNoDataValue + 1;
2345
0
        }
2346
0
        else
2347
0
        {
2348
0
            constexpr double EPS = 1e-3;
2349
0
            if (m_dfNoDataValue == 0)
2350
0
                dfRemappedValue = EPS;
2351
0
            else
2352
0
                dfRemappedValue = m_dfNoDataValue * (1 + EPS);
2353
0
        }
2354
0
    }
2355
0
    const bool bByteOptim =
2356
0
        (eSrcBandDT == GDT_Byte && eBufType == GDT_Byte &&
2357
0
         eSrcMaskBandDT == GDT_Byte && m_dfMaskValueThreshold >= 0 &&
2358
0
         m_dfMaskValueThreshold <= 255 &&
2359
0
         static_cast<int>(m_dfMaskValueThreshold) == m_dfMaskValueThreshold &&
2360
0
         m_dfNoDataValue >= 0 && m_dfNoDataValue <= 255 &&
2361
0
         static_cast<int>(m_dfNoDataValue) == m_dfNoDataValue &&
2362
0
         dfRemappedValue >= 0 && dfRemappedValue <= 255 &&
2363
0
         static_cast<int>(dfRemappedValue) == dfRemappedValue);
2364
0
    GByte *pabyWrkBuffer;
2365
0
    try
2366
0
    {
2367
0
        if (bByteOptim && nOutXOff == 0 && nOutYOff == 0 &&
2368
0
            nOutXSize == nBufXSize && nOutYSize == nBufYSize &&
2369
0
            eSrcBandDT == eBufType && nPixelSpace == nSrcBandDTSize &&
2370
0
            nLineSpace == nPixelSpace * nBufXSize)
2371
0
        {
2372
0
            pabyWrkBuffer = static_cast<GByte *>(pData);
2373
0
        }
2374
0
        else
2375
0
        {
2376
0
            oWorkingState.m_abyWrkBuffer.resize(static_cast<size_t>(nOutXSize) *
2377
0
                                                nOutYSize * nSrcBandDTSize);
2378
0
            pabyWrkBuffer =
2379
0
                reinterpret_cast<GByte *>(oWorkingState.m_abyWrkBuffer.data());
2380
0
        }
2381
0
        oWorkingState.m_abyWrkBufferMask.resize(static_cast<size_t>(nOutXSize) *
2382
0
                                                nOutYSize * nSrcMaskBandDTSize);
2383
0
    }
2384
0
    catch (const std::exception &)
2385
0
    {
2386
0
        CPLError(CE_Failure, CPLE_OutOfMemory,
2387
0
                 "Out of memory when allocating buffers");
2388
0
        return CE_Failure;
2389
0
    }
2390
2391
    /* -------------------------------------------------------------------- */
2392
    /*      Load data.                                                      */
2393
    /* -------------------------------------------------------------------- */
2394
0
    if (!m_osResampling.empty())
2395
0
    {
2396
0
        psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
2397
0
    }
2398
0
    else if (psExtraArgIn != nullptr)
2399
0
    {
2400
0
        psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
2401
0
    }
2402
2403
0
    psExtraArg->bFloatingPointWindowValidity = TRUE;
2404
0
    psExtraArg->dfXOff = dfReqXOff;
2405
0
    psExtraArg->dfYOff = dfReqYOff;
2406
0
    psExtraArg->dfXSize = dfReqXSize;
2407
0
    psExtraArg->dfYSize = dfReqYSize;
2408
0
    if (psExtraArgIn)
2409
0
    {
2410
0
        psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
2411
0
        psExtraArg->pProgressData = psExtraArgIn->pProgressData;
2412
0
    }
2413
2414
0
    if (l_band->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
2415
0
                         pabyWrkBuffer, nOutXSize, nOutYSize, eSrcBandDT, 0, 0,
2416
0
                         psExtraArg) != CE_None)
2417
0
    {
2418
0
        return CE_Failure;
2419
0
    }
2420
2421
0
    if (l_band->GetMaskBand()->RasterIO(
2422
0
            GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
2423
0
            oWorkingState.m_abyWrkBufferMask.data(), nOutXSize, nOutYSize,
2424
0
            eSrcMaskBandDT, 0, 0, psExtraArg) != CE_None)
2425
0
    {
2426
0
        return CE_Failure;
2427
0
    }
2428
2429
    /* -------------------------------------------------------------------- */
2430
    /*      Do the processing.                                              */
2431
    /* -------------------------------------------------------------------- */
2432
2433
0
    GByte *const pabyOut = static_cast<GByte *>(pData) +
2434
0
                           nPixelSpace * nOutXOff +
2435
0
                           static_cast<GPtrDiff_t>(nLineSpace) * nOutYOff;
2436
0
    if (bByteOptim)
2437
0
    {
2438
        // Special case when everything fits on Byte
2439
0
        const GByte nMaskValueThreshold =
2440
0
            static_cast<GByte>(m_dfMaskValueThreshold);
2441
0
        const GByte nNoDataValue = static_cast<GByte>(m_dfNoDataValue);
2442
0
        const GByte nRemappedValue = static_cast<GByte>(dfRemappedValue);
2443
0
        size_t nSrcIdx = 0;
2444
0
        for (int iY = 0; iY < nOutYSize; iY++)
2445
0
        {
2446
0
            GSpacing nDstOffset = iY * nLineSpace;
2447
0
            for (int iX = 0; iX < nOutXSize; iX++)
2448
0
            {
2449
0
                const GByte nMaskVal =
2450
0
                    oWorkingState.m_abyWrkBufferMask[nSrcIdx];
2451
0
                if (nMaskVal <= nMaskValueThreshold)
2452
0
                {
2453
0
                    pabyOut[static_cast<GPtrDiff_t>(nDstOffset)] = nNoDataValue;
2454
0
                }
2455
0
                else
2456
0
                {
2457
0
                    if (pabyWrkBuffer[nSrcIdx] == nNoDataValue)
2458
0
                    {
2459
0
                        pabyOut[static_cast<GPtrDiff_t>(nDstOffset)] =
2460
0
                            nRemappedValue;
2461
0
                    }
2462
0
                    else
2463
0
                    {
2464
0
                        pabyOut[static_cast<GPtrDiff_t>(nDstOffset)] =
2465
0
                            pabyWrkBuffer[nSrcIdx];
2466
0
                    }
2467
0
                }
2468
0
                nDstOffset += nPixelSpace;
2469
0
                nSrcIdx++;
2470
0
            }
2471
0
        }
2472
0
    }
2473
0
    else
2474
0
    {
2475
0
        size_t nSrcIdx = 0;
2476
0
        double dfMaskVal = 0;
2477
0
        const int nBufDTSize = GDALGetDataTypeSizeBytes(eBufType);
2478
0
        std::vector<GByte> abyDstNoData(nBufDTSize);
2479
0
        GDALCopyWords(&m_dfNoDataValue, GDT_Float64, 0, abyDstNoData.data(),
2480
0
                      eBufType, 0, 1);
2481
0
        std::vector<GByte> abyRemappedValue(nBufDTSize);
2482
0
        GDALCopyWords(&dfRemappedValue, GDT_Float64, 0, abyRemappedValue.data(),
2483
0
                      eBufType, 0, 1);
2484
0
        for (int iY = 0; iY < nOutYSize; iY++)
2485
0
        {
2486
0
            GSpacing nDstOffset = iY * nLineSpace;
2487
0
            for (int iX = 0; iX < nOutXSize; iX++)
2488
0
            {
2489
0
                if (eSrcMaskBandDT == GDT_Byte)
2490
0
                {
2491
0
                    dfMaskVal = oWorkingState.m_abyWrkBufferMask[nSrcIdx];
2492
0
                }
2493
0
                else
2494
0
                {
2495
0
                    GDALCopyWords(oWorkingState.m_abyWrkBufferMask.data() +
2496
0
                                      nSrcIdx * nSrcMaskBandDTSize,
2497
0
                                  eSrcMaskBandDT, 0, &dfMaskVal, GDT_Float64, 0,
2498
0
                                  1);
2499
0
                }
2500
0
                void *const pDst =
2501
0
                    pabyOut + static_cast<GPtrDiff_t>(nDstOffset);
2502
0
                if (!(dfMaskVal > m_dfMaskValueThreshold))
2503
0
                {
2504
0
                    memcpy(pDst, abyDstNoData.data(), nBufDTSize);
2505
0
                }
2506
0
                else
2507
0
                {
2508
0
                    const void *const pSrc =
2509
0
                        pabyWrkBuffer + nSrcIdx * nSrcBandDTSize;
2510
0
                    if (eSrcBandDT == eBufType)
2511
0
                    {
2512
                        // coverity[overrun-buffer-arg]
2513
0
                        memcpy(pDst, pSrc, nBufDTSize);
2514
0
                    }
2515
0
                    else
2516
0
                    {
2517
0
                        GDALCopyWords(pSrc, eSrcBandDT, 0, pDst, eBufType, 0,
2518
0
                                      1);
2519
0
                    }
2520
0
                    if (memcmp(pDst, abyDstNoData.data(), nBufDTSize) == 0)
2521
0
                        memcpy(pDst, abyRemappedValue.data(), nBufDTSize);
2522
0
                }
2523
0
                nDstOffset += nPixelSpace;
2524
0
                nSrcIdx++;
2525
0
            }
2526
0
        }
2527
0
    }
2528
2529
0
    if (psExtraArg->pfnProgress)
2530
0
        psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
2531
2532
0
    return CE_None;
2533
0
}
2534
2535
/************************************************************************/
2536
/*                             GetMinimum()                             */
2537
/************************************************************************/
2538
2539
double VRTNoDataFromMaskSource::GetMinimum(int /* nXSize */, int /* nYSize */,
2540
                                           int *pbSuccess)
2541
0
{
2542
0
    *pbSuccess = FALSE;
2543
0
    return 0.0;
2544
0
}
2545
2546
/************************************************************************/
2547
/*                             GetMaximum()                             */
2548
/************************************************************************/
2549
2550
double VRTNoDataFromMaskSource::GetMaximum(int /* nXSize */, int /* nYSize */,
2551
                                           int *pbSuccess)
2552
0
{
2553
0
    *pbSuccess = FALSE;
2554
0
    return 0.0;
2555
0
}
2556
2557
/************************************************************************/
2558
/*                            GetHistogram()                            */
2559
/************************************************************************/
2560
2561
CPLErr VRTNoDataFromMaskSource::GetHistogram(
2562
    int /* nXSize */, int /* nYSize */, double /* dfMin */, double /* dfMax */,
2563
    int /* nBuckets */, GUIntBig * /* panHistogram */,
2564
    int /* bIncludeOutOfRange */, int /* bApproxOK */,
2565
    GDALProgressFunc /* pfnProgress */, void * /* pProgressData */)
2566
0
{
2567
0
    return CE_Failure;
2568
0
}
2569
2570
/************************************************************************/
2571
/* ==================================================================== */
2572
/*                          VRTComplexSource                            */
2573
/* ==================================================================== */
2574
/************************************************************************/
2575
2576
/************************************************************************/
2577
/*                          VRTComplexSource()                          */
2578
/************************************************************************/
2579
2580
VRTComplexSource::VRTComplexSource(const VRTComplexSource *poSrcSource,
2581
                                   double dfXDstRatio, double dfYDstRatio)
2582
0
    : VRTSimpleSource(poSrcSource, dfXDstRatio, dfYDstRatio),
2583
0
      m_nProcessingFlags(poSrcSource->m_nProcessingFlags),
2584
0
      m_dfNoDataValue(poSrcSource->m_dfNoDataValue),
2585
0
      m_osNoDataValueOri(poSrcSource->m_osNoDataValueOri),
2586
0
      m_dfScaleOff(poSrcSource->m_dfScaleOff),
2587
0
      m_dfScaleRatio(poSrcSource->m_dfScaleRatio),
2588
0
      m_bSrcMinMaxDefined(poSrcSource->m_bSrcMinMaxDefined),
2589
0
      m_dfSrcMin(poSrcSource->m_dfSrcMin), m_dfSrcMax(poSrcSource->m_dfSrcMax),
2590
0
      m_dfDstMin(poSrcSource->m_dfDstMin), m_dfDstMax(poSrcSource->m_dfDstMax),
2591
0
      m_dfExponent(poSrcSource->m_dfExponent), m_bClip(poSrcSource->m_bClip),
2592
0
      m_nColorTableComponent(poSrcSource->m_nColorTableComponent),
2593
0
      m_adfLUTInputs(poSrcSource->m_adfLUTInputs),
2594
0
      m_adfLUTOutputs(poSrcSource->m_adfLUTOutputs)
2595
0
{
2596
0
}
2597
2598
/************************************************************************/
2599
/*                           GetTypeStatic()                            */
2600
/************************************************************************/
2601
2602
const char *VRTComplexSource::GetTypeStatic()
2603
0
{
2604
0
    static const char *TYPE = "ComplexSource";
2605
0
    return TYPE;
2606
0
}
2607
2608
/************************************************************************/
2609
/*                            GetType()                                 */
2610
/************************************************************************/
2611
2612
const char *VRTComplexSource::GetType() const
2613
0
{
2614
0
    return GetTypeStatic();
2615
0
}
2616
2617
/************************************************************************/
2618
/*                           SetNoDataValue()                           */
2619
/************************************************************************/
2620
2621
void VRTComplexSource::SetNoDataValue(double dfNewNoDataValue)
2622
2623
0
{
2624
0
    if (dfNewNoDataValue == VRT_NODATA_UNSET)
2625
0
    {
2626
0
        m_nProcessingFlags &= ~PROCESSING_FLAG_NODATA;
2627
0
        m_dfNoDataValue = VRT_NODATA_UNSET;
2628
0
        return;
2629
0
    }
2630
2631
0
    m_nProcessingFlags |= PROCESSING_FLAG_NODATA;
2632
0
    m_dfNoDataValue = dfNewNoDataValue;
2633
0
}
2634
2635
/************************************************************************/
2636
/*                      GetAdjustedNoDataValue()                        */
2637
/************************************************************************/
2638
2639
double VRTComplexSource::GetAdjustedNoDataValue() const
2640
0
{
2641
0
    if ((m_nProcessingFlags & PROCESSING_FLAG_NODATA) != 0)
2642
0
    {
2643
0
        auto l_band = GetRasterBand();
2644
0
        if (l_band && l_band->GetRasterDataType() == GDT_Float32)
2645
0
        {
2646
0
            return GDALAdjustNoDataCloseToFloatMax(m_dfNoDataValue);
2647
0
        }
2648
0
    }
2649
0
    return m_dfNoDataValue;
2650
0
}
2651
2652
/************************************************************************/
2653
/*                           SerializeToXML()                           */
2654
/************************************************************************/
2655
2656
CPLXMLNode *VRTComplexSource::SerializeToXML(const char *pszVRTPath)
2657
2658
0
{
2659
0
    CPLXMLNode *psSrc = VRTSimpleSource::SerializeToXML(pszVRTPath);
2660
2661
0
    if (psSrc == nullptr)
2662
0
        return nullptr;
2663
2664
0
    CPLFree(psSrc->pszValue);
2665
0
    psSrc->pszValue = CPLStrdup(GetTypeStatic());
2666
2667
0
    if ((m_nProcessingFlags & PROCESSING_FLAG_USE_MASK_BAND) != 0)
2668
0
    {
2669
0
        CPLSetXMLValue(psSrc, "UseMaskBand", "true");
2670
0
    }
2671
2672
0
    if ((m_nProcessingFlags & PROCESSING_FLAG_NODATA) != 0)
2673
0
    {
2674
0
        if (!m_osNoDataValueOri.empty() && GetRasterBandNoOpen() == nullptr)
2675
0
        {
2676
0
            CPLSetXMLValue(psSrc, "NODATA", m_osNoDataValueOri.c_str());
2677
0
        }
2678
0
        else
2679
0
        {
2680
0
            GDALDataType eBandDT = GDT_Unknown;
2681
0
            double dfNoDataValue = m_dfNoDataValue;
2682
0
            const auto kMaxFloat = std::numeric_limits<float>::max();
2683
0
            if (std::fabs(std::fabs(m_dfNoDataValue) - kMaxFloat) <
2684
0
                1e-10 * kMaxFloat)
2685
0
            {
2686
0
                auto l_band = GetRasterBand();
2687
0
                if (l_band)
2688
0
                {
2689
0
                    dfNoDataValue = GetAdjustedNoDataValue();
2690
0
                    eBandDT = l_band->GetRasterDataType();
2691
0
                }
2692
0
            }
2693
0
            CPLSetXMLValue(
2694
0
                psSrc, "NODATA",
2695
0
                VRTSerializeNoData(dfNoDataValue, eBandDT, 18).c_str());
2696
0
        }
2697
0
    }
2698
2699
0
    if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_LINEAR) != 0)
2700
0
    {
2701
0
        CPLSetXMLValue(psSrc, "ScaleOffset", CPLSPrintf("%g", m_dfScaleOff));
2702
0
        CPLSetXMLValue(psSrc, "ScaleRatio", CPLSPrintf("%g", m_dfScaleRatio));
2703
0
    }
2704
0
    else if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_EXPONENTIAL) != 0)
2705
0
    {
2706
0
        CPLSetXMLValue(psSrc, "Exponent", CPLSPrintf("%g", m_dfExponent));
2707
0
        if (m_bSrcMinMaxDefined)
2708
0
        {
2709
0
            CPLSetXMLValue(psSrc, "SrcMin", CPLSPrintf("%g", m_dfSrcMin));
2710
0
            CPLSetXMLValue(psSrc, "SrcMax", CPLSPrintf("%g", m_dfSrcMax));
2711
0
        }
2712
0
        CPLSetXMLValue(psSrc, "DstMin", CPLSPrintf("%g", m_dfDstMin));
2713
0
        CPLSetXMLValue(psSrc, "DstMax", CPLSPrintf("%g", m_dfDstMax));
2714
0
        CPLSetXMLValue(psSrc, "Clip", m_bClip ? "true" : "false");
2715
0
    }
2716
2717
0
    if (!m_adfLUTInputs.empty())
2718
0
    {
2719
        // Make sure we print with sufficient precision to address really close
2720
        // entries (#6422).
2721
0
        CPLString osLUT;
2722
0
        if (m_adfLUTInputs.size() >= 2 &&
2723
0
            CPLString().Printf("%g", m_adfLUTInputs[0]) ==
2724
0
                CPLString().Printf("%g", m_adfLUTInputs[1]))
2725
0
        {
2726
0
            osLUT = CPLString().Printf("%.17g:%g", m_adfLUTInputs[0],
2727
0
                                       m_adfLUTOutputs[0]);
2728
0
        }
2729
0
        else
2730
0
        {
2731
0
            osLUT = CPLString().Printf("%g:%g", m_adfLUTInputs[0],
2732
0
                                       m_adfLUTOutputs[0]);
2733
0
        }
2734
0
        for (size_t i = 1; i < m_adfLUTInputs.size(); i++)
2735
0
        {
2736
0
            if (CPLString().Printf("%g", m_adfLUTInputs[i]) ==
2737
0
                    CPLString().Printf("%g", m_adfLUTInputs[i - 1]) ||
2738
0
                (i + 1 < m_adfLUTInputs.size() &&
2739
0
                 CPLString().Printf("%g", m_adfLUTInputs[i]) ==
2740
0
                     CPLString().Printf("%g", m_adfLUTInputs[i + 1])))
2741
0
            {
2742
                // TODO(schwehr): An explanation of the 18 would be helpful.
2743
                // Can someone distill the issue down to a quick comment?
2744
                // https://trac.osgeo.org/gdal/ticket/6422
2745
0
                osLUT += CPLString().Printf(",%.17g:%g", m_adfLUTInputs[i],
2746
0
                                            m_adfLUTOutputs[i]);
2747
0
            }
2748
0
            else
2749
0
            {
2750
0
                osLUT += CPLString().Printf(",%g:%g", m_adfLUTInputs[i],
2751
0
                                            m_adfLUTOutputs[i]);
2752
0
            }
2753
0
        }
2754
0
        CPLSetXMLValue(psSrc, "LUT", osLUT);
2755
0
    }
2756
2757
0
    if (m_nColorTableComponent)
2758
0
    {
2759
0
        CPLSetXMLValue(psSrc, "ColorTableComponent",
2760
0
                       CPLSPrintf("%d", m_nColorTableComponent));
2761
0
    }
2762
2763
0
    return psSrc;
2764
0
}
2765
2766
/************************************************************************/
2767
/*                              XMLInit()                               */
2768
/************************************************************************/
2769
2770
CPLErr VRTComplexSource::XMLInit(const CPLXMLNode *psSrc,
2771
                                 const char *pszVRTPath,
2772
                                 VRTMapSharedResources &oMapSharedSources)
2773
2774
0
{
2775
    /* -------------------------------------------------------------------- */
2776
    /*      Do base initialization.                                         */
2777
    /* -------------------------------------------------------------------- */
2778
0
    {
2779
0
        const CPLErr eErr =
2780
0
            VRTSimpleSource::XMLInit(psSrc, pszVRTPath, oMapSharedSources);
2781
0
        if (eErr != CE_None)
2782
0
            return eErr;
2783
0
    }
2784
2785
    /* -------------------------------------------------------------------- */
2786
    /*      Complex parameters.                                             */
2787
    /* -------------------------------------------------------------------- */
2788
0
    const char *pszScaleOffset = CPLGetXMLValue(psSrc, "ScaleOffset", nullptr);
2789
0
    const char *pszScaleRatio = CPLGetXMLValue(psSrc, "ScaleRatio", nullptr);
2790
0
    if (pszScaleOffset || pszScaleRatio)
2791
0
    {
2792
0
        m_nProcessingFlags |= PROCESSING_FLAG_SCALING_LINEAR;
2793
0
        if (pszScaleOffset)
2794
0
            m_dfScaleOff = CPLAtof(pszScaleOffset);
2795
0
        if (pszScaleRatio)
2796
0
            m_dfScaleRatio = CPLAtof(pszScaleRatio);
2797
0
    }
2798
0
    else if (CPLGetXMLValue(psSrc, "Exponent", nullptr) != nullptr &&
2799
0
             CPLGetXMLValue(psSrc, "DstMin", nullptr) != nullptr &&
2800
0
             CPLGetXMLValue(psSrc, "DstMax", nullptr) != nullptr)
2801
0
    {
2802
0
        m_nProcessingFlags |= PROCESSING_FLAG_SCALING_EXPONENTIAL;
2803
0
        m_dfExponent = CPLAtof(CPLGetXMLValue(psSrc, "Exponent", "1.0"));
2804
2805
0
        const char *pszSrcMin = CPLGetXMLValue(psSrc, "SrcMin", nullptr);
2806
0
        const char *pszSrcMax = CPLGetXMLValue(psSrc, "SrcMax", nullptr);
2807
0
        if (pszSrcMin && pszSrcMax)
2808
0
        {
2809
0
            m_dfSrcMin = CPLAtof(pszSrcMin);
2810
0
            m_dfSrcMax = CPLAtof(pszSrcMax);
2811
0
            m_bSrcMinMaxDefined = true;
2812
0
        }
2813
2814
0
        m_dfDstMin = CPLAtof(CPLGetXMLValue(psSrc, "DstMin", "0.0"));
2815
0
        m_dfDstMax = CPLAtof(CPLGetXMLValue(psSrc, "DstMax", "0.0"));
2816
0
        m_bClip = CPLTestBool(CPLGetXMLValue(psSrc, "Clip", "true"));
2817
0
    }
2818
2819
0
    if (const char *pszNODATA = CPLGetXMLValue(psSrc, "NODATA", nullptr))
2820
0
    {
2821
0
        m_nProcessingFlags |= PROCESSING_FLAG_NODATA;
2822
0
        m_osNoDataValueOri = pszNODATA;
2823
0
        m_dfNoDataValue = CPLAtofM(m_osNoDataValueOri.c_str());
2824
0
    }
2825
2826
0
    const char *pszUseMaskBand = CPLGetXMLValue(psSrc, "UseMaskBand", nullptr);
2827
0
    if (pszUseMaskBand && CPLTestBool(pszUseMaskBand))
2828
0
    {
2829
0
        m_nProcessingFlags |= PROCESSING_FLAG_USE_MASK_BAND;
2830
0
    }
2831
2832
0
    const char *pszLUT = CPLGetXMLValue(psSrc, "LUT", nullptr);
2833
0
    if (pszLUT)
2834
0
    {
2835
0
        const CPLStringList aosValues(
2836
0
            CSLTokenizeString2(pszLUT, ",:", CSLT_ALLOWEMPTYTOKENS));
2837
2838
0
        const int nLUTItemCount = aosValues.size() / 2;
2839
0
        try
2840
0
        {
2841
0
            m_adfLUTInputs.resize(nLUTItemCount);
2842
0
            m_adfLUTOutputs.resize(nLUTItemCount);
2843
0
        }
2844
0
        catch (const std::bad_alloc &e)
2845
0
        {
2846
0
            CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
2847
0
            m_adfLUTInputs.clear();
2848
0
            m_adfLUTOutputs.clear();
2849
0
            return CE_Failure;
2850
0
        }
2851
2852
0
        for (int nIndex = 0; nIndex < nLUTItemCount; nIndex++)
2853
0
        {
2854
0
            m_adfLUTInputs[nIndex] = CPLAtof(aosValues[nIndex * 2]);
2855
0
            m_adfLUTOutputs[nIndex] = CPLAtof(aosValues[nIndex * 2 + 1]);
2856
2857
            // Enforce the requirement that the LUT input array is
2858
            // monotonically non-decreasing.
2859
0
            if (std::isnan(m_adfLUTInputs[nIndex]) && nIndex != 0)
2860
0
            {
2861
0
                CPLError(CE_Failure, CPLE_AppDefined,
2862
0
                         "A Not-A-Number (NaN) source value should be the "
2863
0
                         "first one of the LUT.");
2864
0
                m_adfLUTInputs.clear();
2865
0
                m_adfLUTOutputs.clear();
2866
0
                return CE_Failure;
2867
0
            }
2868
0
            else if (nIndex > 0 &&
2869
0
                     m_adfLUTInputs[nIndex] < m_adfLUTInputs[nIndex - 1])
2870
0
            {
2871
0
                CPLError(CE_Failure, CPLE_AppDefined,
2872
0
                         "Source values of the LUT are not listed in a "
2873
0
                         "monotonically non-decreasing order");
2874
0
                m_adfLUTInputs.clear();
2875
0
                m_adfLUTOutputs.clear();
2876
0
                return CE_Failure;
2877
0
            }
2878
0
        }
2879
2880
0
        m_nProcessingFlags |= PROCESSING_FLAG_LUT;
2881
0
    }
2882
2883
0
    const char *pszColorTableComponent =
2884
0
        CPLGetXMLValue(psSrc, "ColorTableComponent", nullptr);
2885
0
    if (pszColorTableComponent)
2886
0
    {
2887
0
        m_nColorTableComponent = atoi(pszColorTableComponent);
2888
0
        m_nProcessingFlags |= PROCESSING_FLAG_COLOR_TABLE_EXPANSION;
2889
0
    }
2890
2891
0
    return CE_None;
2892
0
}
2893
2894
/************************************************************************/
2895
/*                              LookupValue()                           */
2896
/************************************************************************/
2897
2898
double VRTComplexSource::LookupValue(double dfInput)
2899
0
{
2900
0
    auto beginIter = m_adfLUTInputs.begin();
2901
0
    auto endIter = m_adfLUTInputs.end();
2902
0
    size_t offset = 0;
2903
0
    if (std::isnan(m_adfLUTInputs[0]))
2904
0
    {
2905
0
        if (std::isnan(dfInput) || m_adfLUTInputs.size() == 1)
2906
0
            return m_adfLUTOutputs[0];
2907
0
        ++beginIter;
2908
0
        offset = 1;
2909
0
    }
2910
2911
    // Find the index of the first element in the LUT input array that
2912
    // is not smaller than the input value.
2913
0
    const size_t i =
2914
0
        offset +
2915
0
        std::distance(beginIter, std::lower_bound(beginIter, endIter, dfInput));
2916
2917
0
    if (i == offset)
2918
0
        return m_adfLUTOutputs[offset];
2919
2920
    // If the index is beyond the end of the LUT input array, the input
2921
    // value is larger than all the values in the array.
2922
0
    if (i == m_adfLUTInputs.size())
2923
0
        return m_adfLUTOutputs.back();
2924
2925
0
    if (m_adfLUTInputs[i] == dfInput)
2926
0
        return m_adfLUTOutputs[i];
2927
2928
    // Otherwise, interpolate.
2929
0
    return m_adfLUTOutputs[i - 1] +
2930
0
           (dfInput - m_adfLUTInputs[i - 1]) *
2931
0
               ((m_adfLUTOutputs[i] - m_adfLUTOutputs[i - 1]) /
2932
0
                (m_adfLUTInputs[i] - m_adfLUTInputs[i - 1]));
2933
0
}
2934
2935
/************************************************************************/
2936
/*                         SetLinearScaling()                           */
2937
/************************************************************************/
2938
2939
void VRTComplexSource::SetLinearScaling(double dfOffset, double dfScale)
2940
0
{
2941
0
    m_nProcessingFlags &= ~PROCESSING_FLAG_SCALING_EXPONENTIAL;
2942
0
    m_nProcessingFlags |= PROCESSING_FLAG_SCALING_LINEAR;
2943
0
    m_dfScaleOff = dfOffset;
2944
0
    m_dfScaleRatio = dfScale;
2945
0
}
2946
2947
/************************************************************************/
2948
/*                         SetPowerScaling()                           */
2949
/************************************************************************/
2950
2951
void VRTComplexSource::SetPowerScaling(double dfExponentIn, double dfSrcMinIn,
2952
                                       double dfSrcMaxIn, double dfDstMinIn,
2953
                                       double dfDstMaxIn, bool bClip)
2954
0
{
2955
0
    m_nProcessingFlags &= ~PROCESSING_FLAG_SCALING_LINEAR;
2956
0
    m_nProcessingFlags |= PROCESSING_FLAG_SCALING_EXPONENTIAL;
2957
0
    m_dfExponent = dfExponentIn;
2958
0
    m_dfSrcMin = dfSrcMinIn;
2959
0
    m_dfSrcMax = dfSrcMaxIn;
2960
0
    m_dfDstMin = dfDstMinIn;
2961
0
    m_dfDstMax = dfDstMaxIn;
2962
0
    m_bSrcMinMaxDefined = true;
2963
0
    m_bClip = bClip;
2964
0
}
2965
2966
/************************************************************************/
2967
/*                    SetColorTableComponent()                          */
2968
/************************************************************************/
2969
2970
void VRTComplexSource::SetColorTableComponent(int nComponent)
2971
0
{
2972
0
    m_nProcessingFlags |= PROCESSING_FLAG_COLOR_TABLE_EXPANSION;
2973
0
    m_nColorTableComponent = nComponent;
2974
0
}
2975
2976
/************************************************************************/
2977
/*                              RasterIO()                              */
2978
/************************************************************************/
2979
2980
CPLErr VRTComplexSource::RasterIO(GDALDataType eVRTBandDataType, int nXOff,
2981
                                  int nYOff, int nXSize, int nYSize,
2982
                                  void *pData, int nBufXSize, int nBufYSize,
2983
                                  GDALDataType eBufType, GSpacing nPixelSpace,
2984
                                  GSpacing nLineSpace,
2985
                                  GDALRasterIOExtraArg *psExtraArgIn,
2986
                                  WorkingState &oWorkingState)
2987
2988
0
{
2989
0
    GDALRasterIOExtraArg sExtraArg;
2990
0
    INIT_RASTERIO_EXTRA_ARG(sExtraArg);
2991
0
    GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
2992
2993
0
    double dfXOff = nXOff;
2994
0
    double dfYOff = nYOff;
2995
0
    double dfXSize = nXSize;
2996
0
    double dfYSize = nYSize;
2997
0
    if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
2998
0
    {
2999
0
        dfXOff = psExtraArgIn->dfXOff;
3000
0
        dfYOff = psExtraArgIn->dfYOff;
3001
0
        dfXSize = psExtraArgIn->dfXSize;
3002
0
        dfYSize = psExtraArgIn->dfYSize;
3003
0
    }
3004
3005
    // The window we will actually request from the source raster band.
3006
0
    double dfReqXOff = 0.0;
3007
0
    double dfReqYOff = 0.0;
3008
0
    double dfReqXSize = 0.0;
3009
0
    double dfReqYSize = 0.0;
3010
0
    int nReqXOff = 0;
3011
0
    int nReqYOff = 0;
3012
0
    int nReqXSize = 0;
3013
0
    int nReqYSize = 0;
3014
3015
    // The window we will actual set _within_ the pData buffer.
3016
0
    int nOutXOff = 0;
3017
0
    int nOutYOff = 0;
3018
0
    int nOutXSize = 0;
3019
0
    int nOutYSize = 0;
3020
3021
0
    bool bError = false;
3022
0
    if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
3023
0
                         &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
3024
0
                         &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
3025
0
                         &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
3026
0
    {
3027
0
        return bError ? CE_Failure : CE_None;
3028
0
    }
3029
#if DEBUG_VERBOSE
3030
    CPLDebug("VRT",
3031
             "nXOff=%d, nYOff=%d, nXSize=%d, nYSize=%d, nBufXSize=%d, "
3032
             "nBufYSize=%d,\n"
3033
             "dfReqXOff=%g, dfReqYOff=%g, dfReqXSize=%g, dfReqYSize=%g,\n"
3034
             "nReqXOff=%d, nReqYOff=%d, nReqXSize=%d, nReqYSize=%d,\n"
3035
             "nOutXOff=%d, nOutYOff=%d, nOutXSize=%d, nOutYSize=%d",
3036
             nXOff, nYOff, nXSize, nYSize, nBufXSize, nBufYSize, dfReqXOff,
3037
             dfReqYOff, dfReqXSize, dfReqYSize, nReqXOff, nReqYOff, nReqXSize,
3038
             nReqYSize, nOutXOff, nOutYOff, nOutXSize, nOutYSize);
3039
#endif
3040
3041
0
    auto poSourceBand = GetRasterBand();
3042
0
    if (!poSourceBand)
3043
0
        return CE_Failure;
3044
3045
0
    if (!m_osResampling.empty())
3046
0
    {
3047
0
        psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
3048
0
    }
3049
0
    else if (psExtraArgIn != nullptr)
3050
0
    {
3051
0
        psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
3052
0
    }
3053
0
    psExtraArg->bFloatingPointWindowValidity = TRUE;
3054
0
    psExtraArg->dfXOff = dfReqXOff;
3055
0
    psExtraArg->dfYOff = dfReqYOff;
3056
0
    psExtraArg->dfXSize = dfReqXSize;
3057
0
    psExtraArg->dfYSize = dfReqYSize;
3058
0
    if (psExtraArgIn)
3059
0
    {
3060
0
        psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
3061
0
        psExtraArg->pProgressData = psExtraArgIn->pProgressData;
3062
0
    }
3063
3064
0
    GByte *const pabyOut = static_cast<GByte *>(pData) +
3065
0
                           nPixelSpace * nOutXOff +
3066
0
                           static_cast<GPtrDiff_t>(nLineSpace) * nOutYOff;
3067
0
    const auto eSourceType = poSourceBand->GetRasterDataType();
3068
0
    if (m_nProcessingFlags == PROCESSING_FLAG_NODATA)
3069
0
    {
3070
        // Optimization if doing only nodata processing
3071
0
        if (eSourceType == GDT_Byte)
3072
0
        {
3073
0
            if (!GDALIsValueInRange<GByte>(m_dfNoDataValue))
3074
0
            {
3075
0
                return VRTSimpleSource::RasterIO(
3076
0
                    eVRTBandDataType, nXOff, nYOff, nXSize, nYSize, pData,
3077
0
                    nBufXSize, nBufYSize, eBufType, nPixelSpace, nLineSpace,
3078
0
                    psExtraArgIn, oWorkingState);
3079
0
            }
3080
3081
0
            return RasterIOProcessNoData<GByte, GDT_Byte>(
3082
0
                poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3083
0
                nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3084
0
                nLineSpace, psExtraArg, oWorkingState);
3085
0
        }
3086
0
        else if (eSourceType == GDT_Int16)
3087
0
        {
3088
0
            if (!GDALIsValueInRange<int16_t>(m_dfNoDataValue))
3089
0
            {
3090
0
                return VRTSimpleSource::RasterIO(
3091
0
                    eVRTBandDataType, nXOff, nYOff, nXSize, nYSize, pData,
3092
0
                    nBufXSize, nBufYSize, eBufType, nPixelSpace, nLineSpace,
3093
0
                    psExtraArgIn, oWorkingState);
3094
0
            }
3095
3096
0
            return RasterIOProcessNoData<int16_t, GDT_Int16>(
3097
0
                poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3098
0
                nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3099
0
                nLineSpace, psExtraArg, oWorkingState);
3100
0
        }
3101
0
        else if (eSourceType == GDT_UInt16)
3102
0
        {
3103
0
            if (!GDALIsValueInRange<uint16_t>(m_dfNoDataValue))
3104
0
            {
3105
0
                return VRTSimpleSource::RasterIO(
3106
0
                    eVRTBandDataType, nXOff, nYOff, nXSize, nYSize, pData,
3107
0
                    nBufXSize, nBufYSize, eBufType, nPixelSpace, nLineSpace,
3108
0
                    psExtraArgIn, oWorkingState);
3109
0
            }
3110
3111
0
            return RasterIOProcessNoData<uint16_t, GDT_UInt16>(
3112
0
                poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3113
0
                nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3114
0
                nLineSpace, psExtraArg, oWorkingState);
3115
0
        }
3116
0
    }
3117
3118
0
    const bool bIsComplex = CPL_TO_BOOL(GDALDataTypeIsComplex(eBufType));
3119
0
    CPLErr eErr;
3120
    // For Int32, float32 isn't sufficiently precise as working data type
3121
0
    if (eVRTBandDataType == GDT_CInt32 || eVRTBandDataType == GDT_CFloat64 ||
3122
0
        eVRTBandDataType == GDT_Int32 || eVRTBandDataType == GDT_UInt32 ||
3123
0
        eVRTBandDataType == GDT_Int64 || eVRTBandDataType == GDT_UInt64 ||
3124
0
        eVRTBandDataType == GDT_Float64 || eSourceType == GDT_Int32 ||
3125
0
        eSourceType == GDT_UInt32 || eSourceType == GDT_Int64 ||
3126
0
        eSourceType == GDT_UInt64)
3127
0
    {
3128
0
        eErr = RasterIOInternal<double>(
3129
0
            poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3130
0
            nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3131
0
            nLineSpace, psExtraArg, bIsComplex ? GDT_CFloat64 : GDT_Float64,
3132
0
            oWorkingState);
3133
0
    }
3134
0
    else
3135
0
    {
3136
0
        eErr = RasterIOInternal<float>(
3137
0
            poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3138
0
            nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3139
0
            nLineSpace, psExtraArg, bIsComplex ? GDT_CFloat32 : GDT_Float32,
3140
0
            oWorkingState);
3141
0
    }
3142
3143
0
    if (psExtraArg->pfnProgress)
3144
0
        psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
3145
3146
0
    return eErr;
3147
0
}
3148
3149
/************************************************************************/
3150
/*                              hasZeroByte()                           */
3151
/************************************************************************/
3152
3153
CPL_NOSANITIZE_UNSIGNED_INT_OVERFLOW
3154
static inline bool hasZeroByte(uint32_t v)
3155
0
{
3156
    // Cf https://graphics.stanford.edu/~seander/bithacks.html#ZeroInWord
3157
0
    return (((v)-0x01010101U) & ~(v)&0x80808080U) != 0;
3158
0
}
3159
3160
/************************************************************************/
3161
/*                                CopyWord()                            */
3162
/************************************************************************/
3163
3164
template <class SrcType>
3165
static void CopyWord(const SrcType *pSrcVal, GDALDataType eSrcType,
3166
                     void *pDstVal, GDALDataType eDstType)
3167
0
{
3168
0
    switch (eDstType)
3169
0
    {
3170
0
        case GDT_Byte:
3171
0
            GDALCopyWord(*pSrcVal, *static_cast<uint8_t *>(pDstVal));
3172
0
            break;
3173
0
        case GDT_Int8:
3174
0
            GDALCopyWord(*pSrcVal, *static_cast<int8_t *>(pDstVal));
3175
0
            break;
3176
0
        case GDT_UInt16:
3177
0
            GDALCopyWord(*pSrcVal, *static_cast<uint16_t *>(pDstVal));
3178
0
            break;
3179
0
        case GDT_Int16:
3180
0
            GDALCopyWord(*pSrcVal, *static_cast<int16_t *>(pDstVal));
3181
0
            break;
3182
0
        case GDT_UInt32:
3183
0
            GDALCopyWord(*pSrcVal, *static_cast<uint32_t *>(pDstVal));
3184
0
            break;
3185
0
        case GDT_Int32:
3186
0
            GDALCopyWord(*pSrcVal, *static_cast<int32_t *>(pDstVal));
3187
0
            break;
3188
0
        case GDT_UInt64:
3189
0
            GDALCopyWord(*pSrcVal, *static_cast<uint64_t *>(pDstVal));
3190
0
            break;
3191
0
        case GDT_Int64:
3192
0
            GDALCopyWord(*pSrcVal, *static_cast<int64_t *>(pDstVal));
3193
0
            break;
3194
0
        case GDT_Float32:
3195
0
            GDALCopyWord(*pSrcVal, *static_cast<float *>(pDstVal));
3196
0
            break;
3197
0
        case GDT_Float64:
3198
0
            GDALCopyWord(*pSrcVal, *static_cast<double *>(pDstVal));
3199
0
            break;
3200
0
        default:
3201
0
            GDALCopyWords(pSrcVal, eSrcType, 0, pDstVal, eDstType, 0, 1);
3202
0
            break;
3203
0
    }
3204
0
}
Unexecuted instantiation: vrtsources.cpp:void CopyWord<unsigned char>(unsigned char const*, GDALDataType, void*, GDALDataType)
Unexecuted instantiation: vrtsources.cpp:void CopyWord<short>(short const*, GDALDataType, void*, GDALDataType)
Unexecuted instantiation: vrtsources.cpp:void CopyWord<unsigned short>(unsigned short const*, GDALDataType, void*, GDALDataType)
Unexecuted instantiation: vrtsources.cpp:void CopyWord<double>(double const*, GDALDataType, void*, GDALDataType)
Unexecuted instantiation: vrtsources.cpp:void CopyWord<float>(float const*, GDALDataType, void*, GDALDataType)
3205
3206
/************************************************************************/
3207
/*                       RasterIOProcessNoData()                        */
3208
/************************************************************************/
3209
3210
// This method is an optimization of the generic RasterIOInternal()
3211
// that deals with a VRTComplexSource with only a NODATA value in it, and
3212
// no other processing flags.
3213
3214
// nReqXOff, nReqYOff, nReqXSize, nReqYSize are expressed in source band
3215
// referential.
3216
template <class SourceDT, GDALDataType eSourceType>
3217
CPLErr VRTComplexSource::RasterIOProcessNoData(
3218
    GDALRasterBand *poSourceBand, GDALDataType eVRTBandDataType, int nReqXOff,
3219
    int nReqYOff, int nReqXSize, int nReqYSize, void *pData, int nOutXSize,
3220
    int nOutYSize, GDALDataType eBufType, GSpacing nPixelSpace,
3221
    GSpacing nLineSpace, GDALRasterIOExtraArg *psExtraArg,
3222
    WorkingState &oWorkingState)
3223
0
{
3224
0
    CPLAssert(m_nProcessingFlags == PROCESSING_FLAG_NODATA);
3225
0
    CPLAssert(GDALIsValueInRange<SourceDT>(m_dfNoDataValue));
3226
3227
    /* -------------------------------------------------------------------- */
3228
    /*      Read into a temporary buffer.                                   */
3229
    /* -------------------------------------------------------------------- */
3230
0
    try
3231
0
    {
3232
        // Cannot overflow since pData should at least have that number of
3233
        // elements
3234
0
        const size_t nPixelCount = static_cast<size_t>(nOutXSize) * nOutYSize;
3235
0
        if (nPixelCount >
3236
0
            static_cast<size_t>(std::numeric_limits<ptrdiff_t>::max()) /
3237
0
                sizeof(SourceDT))
3238
0
        {
3239
0
            CPLError(CE_Failure, CPLE_OutOfMemory,
3240
0
                     "Too large temporary buffer");
3241
0
            return CE_Failure;
3242
0
        }
3243
0
        oWorkingState.m_abyWrkBuffer.resize(sizeof(SourceDT) * nPixelCount);
3244
0
    }
3245
0
    catch (const std::bad_alloc &e)
3246
0
    {
3247
0
        CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
3248
0
        return CE_Failure;
3249
0
    }
3250
0
    const auto paSrcData =
3251
0
        reinterpret_cast<const SourceDT *>(oWorkingState.m_abyWrkBuffer.data());
3252
3253
0
    const GDALRIOResampleAlg eResampleAlgBack = psExtraArg->eResampleAlg;
3254
0
    if (!m_osResampling.empty())
3255
0
    {
3256
0
        psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
3257
0
    }
3258
3259
0
    const CPLErr eErr = poSourceBand->RasterIO(
3260
0
        GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
3261
0
        oWorkingState.m_abyWrkBuffer.data(), nOutXSize, nOutYSize, eSourceType,
3262
0
        sizeof(SourceDT), sizeof(SourceDT) * static_cast<GSpacing>(nOutXSize),
3263
0
        psExtraArg);
3264
0
    if (!m_osResampling.empty())
3265
0
        psExtraArg->eResampleAlg = eResampleAlgBack;
3266
3267
0
    if (eErr != CE_None)
3268
0
    {
3269
0
        return eErr;
3270
0
    }
3271
3272
0
    const auto nNoDataValue = static_cast<SourceDT>(m_dfNoDataValue);
3273
0
    size_t idxBuffer = 0;
3274
0
    if (eSourceType == eBufType &&
3275
0
        !GDALDataTypeIsConversionLossy(eSourceType, eVRTBandDataType))
3276
0
    {
3277
        // Most optimized case: the output type is the same as the source type,
3278
        // and conversion from the source type to the VRT band data type is
3279
        // not lossy
3280
0
        for (int iY = 0; iY < nOutYSize; iY++)
3281
0
        {
3282
0
            GByte *pDstLocation = static_cast<GByte *>(pData) +
3283
0
                                  static_cast<GPtrDiff_t>(nLineSpace) * iY;
3284
3285
0
            int iX = 0;
3286
0
            if (sizeof(SourceDT) == 1 && nPixelSpace == 1)
3287
0
            {
3288
                // Optimization to detect more quickly if source pixels are
3289
                // at nodata.
3290
0
                const GByte byNoDataValue = static_cast<GByte>(nNoDataValue);
3291
0
                const uint32_t wordNoData =
3292
0
                    (static_cast<uint32_t>(byNoDataValue) << 24) |
3293
0
                    (byNoDataValue << 16) | (byNoDataValue << 8) |
3294
0
                    byNoDataValue;
3295
3296
                // Warning: hasZeroByte() assumes WORD_SIZE = 4
3297
0
                constexpr int WORD_SIZE = 4;
3298
0
                for (; iX < nOutXSize - (WORD_SIZE - 1); iX += WORD_SIZE)
3299
0
                {
3300
0
                    uint32_t v;
3301
0
                    static_assert(sizeof(v) == WORD_SIZE,
3302
0
                                  "sizeof(v) == WORD_SIZE");
3303
0
                    memcpy(&v, paSrcData + idxBuffer, sizeof(v));
3304
                    // Cf https://graphics.stanford.edu/~seander/bithacks.html#ValueInWord
3305
0
                    if (!hasZeroByte(v ^ wordNoData))
3306
0
                    {
3307
                        // No bytes are at nodata
3308
0
                        memcpy(pDstLocation, &v, WORD_SIZE);
3309
0
                        idxBuffer += WORD_SIZE;
3310
0
                        pDstLocation += WORD_SIZE;
3311
0
                    }
3312
0
                    else if (v == wordNoData)
3313
0
                    {
3314
                        // All bytes are at nodata
3315
0
                        idxBuffer += WORD_SIZE;
3316
0
                        pDstLocation += WORD_SIZE;
3317
0
                    }
3318
0
                    else
3319
0
                    {
3320
                        // There are both bytes at nodata and valid bytes
3321
0
                        for (int k = 0; k < WORD_SIZE; ++k)
3322
0
                        {
3323
0
                            if (paSrcData[idxBuffer] != nNoDataValue)
3324
0
                            {
3325
0
                                memcpy(pDstLocation, &paSrcData[idxBuffer],
3326
0
                                       sizeof(SourceDT));
3327
0
                            }
3328
0
                            idxBuffer++;
3329
0
                            pDstLocation += nPixelSpace;
3330
0
                        }
3331
0
                    }
3332
0
                }
3333
0
            }
3334
3335
0
            for (; iX < nOutXSize;
3336
0
                 iX++, pDstLocation += nPixelSpace, idxBuffer++)
3337
0
            {
3338
0
                if (paSrcData[idxBuffer] != nNoDataValue)
3339
0
                {
3340
0
                    memcpy(pDstLocation, &paSrcData[idxBuffer],
3341
0
                           sizeof(SourceDT));
3342
0
                }
3343
0
            }
3344
0
        }
3345
0
    }
3346
0
    else if (!GDALDataTypeIsConversionLossy(eSourceType, eVRTBandDataType))
3347
0
    {
3348
        // Conversion from the source type to the VRT band data type is
3349
        // not lossy, so we can directly convert from the source type to
3350
        // the the output type
3351
0
        for (int iY = 0; iY < nOutYSize; iY++)
3352
0
        {
3353
0
            GByte *pDstLocation = static_cast<GByte *>(pData) +
3354
0
                                  static_cast<GPtrDiff_t>(nLineSpace) * iY;
3355
3356
0
            for (int iX = 0; iX < nOutXSize;
3357
0
                 iX++, pDstLocation += nPixelSpace, idxBuffer++)
3358
0
            {
3359
0
                if (paSrcData[idxBuffer] != nNoDataValue)
3360
0
                {
3361
0
                    CopyWord(&paSrcData[idxBuffer], eSourceType, pDstLocation,
3362
0
                             eBufType);
3363
0
                }
3364
0
            }
3365
0
        }
3366
0
    }
3367
0
    else
3368
0
    {
3369
0
        GByte abyTemp[2 * sizeof(double)];
3370
0
        for (int iY = 0; iY < nOutYSize; iY++)
3371
0
        {
3372
0
            GByte *pDstLocation = static_cast<GByte *>(pData) +
3373
0
                                  static_cast<GPtrDiff_t>(nLineSpace) * iY;
3374
3375
0
            for (int iX = 0; iX < nOutXSize;
3376
0
                 iX++, pDstLocation += nPixelSpace, idxBuffer++)
3377
0
            {
3378
0
                if (paSrcData[idxBuffer] != nNoDataValue)
3379
0
                {
3380
                    // Convert first to the VRTRasterBand data type
3381
                    // to get its clamping, before outputting to buffer data type
3382
0
                    CopyWord(&paSrcData[idxBuffer], eSourceType, abyTemp,
3383
0
                             eVRTBandDataType);
3384
0
                    GDALCopyWords(abyTemp, eVRTBandDataType, 0, pDstLocation,
3385
0
                                  eBufType, 0, 1);
3386
0
                }
3387
0
            }
3388
0
        }
3389
0
    }
3390
3391
0
    return CE_None;
3392
0
}
Unexecuted instantiation: CPLErr VRTComplexSource::RasterIOProcessNoData<unsigned char, (GDALDataType)1>(GDALRasterBand*, GDALDataType, int, int, int, int, void*, int, int, GDALDataType, long long, long long, GDALRasterIOExtraArg*, VRTSource::WorkingState&)
Unexecuted instantiation: CPLErr VRTComplexSource::RasterIOProcessNoData<short, (GDALDataType)3>(GDALRasterBand*, GDALDataType, int, int, int, int, void*, int, int, GDALDataType, long long, long long, GDALRasterIOExtraArg*, VRTSource::WorkingState&)
Unexecuted instantiation: CPLErr VRTComplexSource::RasterIOProcessNoData<unsigned short, (GDALDataType)2>(GDALRasterBand*, GDALDataType, int, int, int, int, void*, int, int, GDALDataType, long long, long long, GDALRasterIOExtraArg*, VRTSource::WorkingState&)
3393
3394
/************************************************************************/
3395
/*                          RasterIOInternal()                          */
3396
/************************************************************************/
3397
3398
// nReqXOff, nReqYOff, nReqXSize, nReqYSize are expressed in source band
3399
// referential.
3400
template <class WorkingDT>
3401
CPLErr VRTComplexSource::RasterIOInternal(
3402
    GDALRasterBand *poSourceBand, GDALDataType eVRTBandDataType, int nReqXOff,
3403
    int nReqYOff, int nReqXSize, int nReqYSize, void *pData, int nOutXSize,
3404
    int nOutYSize, GDALDataType eBufType, GSpacing nPixelSpace,
3405
    GSpacing nLineSpace, GDALRasterIOExtraArg *psExtraArg,
3406
    GDALDataType eWrkDataType, WorkingState &oWorkingState)
3407
0
{
3408
0
    const GDALColorTable *poColorTable = nullptr;
3409
0
    const bool bIsComplex = CPL_TO_BOOL(GDALDataTypeIsComplex(eBufType));
3410
0
    const int nWordSize = GDALGetDataTypeSizeBytes(eWrkDataType);
3411
0
    assert(nWordSize != 0);
3412
3413
    // If no explicit <NODATA> is set, but UseMaskBand is set, and the band
3414
    // has a nodata value, then use it as if it was set as <NODATA>
3415
0
    int bNoDataSet = (m_nProcessingFlags & PROCESSING_FLAG_NODATA) != 0;
3416
0
    double dfNoDataValue = GetAdjustedNoDataValue();
3417
3418
0
    if ((m_nProcessingFlags & PROCESSING_FLAG_USE_MASK_BAND) != 0 &&
3419
0
        poSourceBand->GetMaskFlags() == GMF_NODATA)
3420
0
    {
3421
0
        dfNoDataValue = poSourceBand->GetNoDataValue(&bNoDataSet);
3422
0
    }
3423
3424
0
    const bool bNoDataSetIsNan = bNoDataSet && std::isnan(dfNoDataValue);
3425
0
    const bool bNoDataSetAndNotNan =
3426
0
        bNoDataSet && !std::isnan(dfNoDataValue) &&
3427
0
        GDALIsValueInRange<WorkingDT>(dfNoDataValue);
3428
0
    const auto fWorkingDataTypeNoData = static_cast<WorkingDT>(dfNoDataValue);
3429
3430
0
    const GByte *pabyMask = nullptr;
3431
0
    const WorkingDT *pafData = nullptr;
3432
0
    if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_LINEAR) != 0 &&
3433
0
        m_dfScaleRatio == 0 && bNoDataSet == FALSE &&
3434
0
        (m_nProcessingFlags & PROCESSING_FLAG_USE_MASK_BAND) == 0)
3435
0
    {
3436
        /* ------------------------------------------------------------------ */
3437
        /*      Optimization when writing a constant value */
3438
        /*      (used by the -addalpha option of gdalbuildvrt) */
3439
        /* ------------------------------------------------------------------ */
3440
        // Already set to NULL when defined.
3441
        // pafData = NULL;
3442
0
    }
3443
0
    else
3444
0
    {
3445
        /* ---------------------------------------------------------------- */
3446
        /*      Read into a temporary buffer.                               */
3447
        /* ---------------------------------------------------------------- */
3448
0
        const size_t nPixelCount = static_cast<size_t>(nOutXSize) * nOutYSize;
3449
0
        try
3450
0
        {
3451
            // Cannot overflow since pData should at least have that number of
3452
            // elements
3453
0
            if (nPixelCount >
3454
0
                static_cast<size_t>(std::numeric_limits<ptrdiff_t>::max()) /
3455
0
                    static_cast<size_t>(nWordSize))
3456
0
            {
3457
0
                CPLError(CE_Failure, CPLE_OutOfMemory,
3458
0
                         "Too large temporary buffer");
3459
0
                return CE_Failure;
3460
0
            }
3461
0
            oWorkingState.m_abyWrkBuffer.resize(nWordSize * nPixelCount);
3462
0
        }
3463
0
        catch (const std::bad_alloc &e)
3464
0
        {
3465
0
            CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
3466
0
            return CE_Failure;
3467
0
        }
3468
0
        pafData = reinterpret_cast<const WorkingDT *>(
3469
0
            oWorkingState.m_abyWrkBuffer.data());
3470
3471
0
        const GDALRIOResampleAlg eResampleAlgBack = psExtraArg->eResampleAlg;
3472
0
        if (!m_osResampling.empty())
3473
0
        {
3474
0
            psExtraArg->eResampleAlg =
3475
0
                GDALRasterIOGetResampleAlg(m_osResampling);
3476
0
        }
3477
3478
0
        const CPLErr eErr = poSourceBand->RasterIO(
3479
0
            GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
3480
0
            oWorkingState.m_abyWrkBuffer.data(), nOutXSize, nOutYSize,
3481
0
            eWrkDataType, nWordSize,
3482
0
            nWordSize * static_cast<GSpacing>(nOutXSize), psExtraArg);
3483
0
        if (!m_osResampling.empty())
3484
0
            psExtraArg->eResampleAlg = eResampleAlgBack;
3485
3486
0
        if (eErr != CE_None)
3487
0
        {
3488
0
            return eErr;
3489
0
        }
3490
3491
        // Allocate and read mask band if needed
3492
0
        if (!bNoDataSet &&
3493
0
            (m_nProcessingFlags & PROCESSING_FLAG_USE_MASK_BAND) != 0 &&
3494
0
            (poSourceBand->GetMaskFlags() != GMF_ALL_VALID ||
3495
0
             poSourceBand->GetColorInterpretation() == GCI_AlphaBand ||
3496
0
             GetMaskBandMainBand() != nullptr))
3497
0
        {
3498
0
            try
3499
0
            {
3500
0
                oWorkingState.m_abyWrkBufferMask.resize(nPixelCount);
3501
0
            }
3502
0
            catch (const std::exception &)
3503
0
            {
3504
0
                CPLError(CE_Failure, CPLE_OutOfMemory,
3505
0
                         "Out of memory when allocating mask buffer");
3506
0
                return CE_Failure;
3507
0
            }
3508
0
            pabyMask = reinterpret_cast<const GByte *>(
3509
0
                oWorkingState.m_abyWrkBufferMask.data());
3510
0
            auto poMaskBand =
3511
0
                (poSourceBand->GetColorInterpretation() == GCI_AlphaBand ||
3512
0
                 GetMaskBandMainBand() != nullptr)
3513
0
                    ? poSourceBand
3514
0
                    : poSourceBand->GetMaskBand();
3515
0
            if (poMaskBand->RasterIO(
3516
0
                    GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
3517
0
                    oWorkingState.m_abyWrkBufferMask.data(), nOutXSize,
3518
0
                    nOutYSize, GDT_Byte, 1, static_cast<GSpacing>(nOutXSize),
3519
0
                    psExtraArg) != CE_None)
3520
0
            {
3521
0
                return CE_Failure;
3522
0
            }
3523
0
        }
3524
3525
0
        if (m_nColorTableComponent != 0)
3526
0
        {
3527
0
            poColorTable = poSourceBand->GetColorTable();
3528
0
            if (poColorTable == nullptr)
3529
0
            {
3530
0
                CPLError(CE_Failure, CPLE_AppDefined,
3531
0
                         "Source band has no color table.");
3532
0
                return CE_Failure;
3533
0
            }
3534
0
        }
3535
0
    }
3536
3537
    /* -------------------------------------------------------------------- */
3538
    /*      Selectively copy into output buffer with nodata masking,        */
3539
    /*      and/or scaling.                                                 */
3540
    /* -------------------------------------------------------------------- */
3541
3542
0
    const bool bTwoStepDataTypeConversion =
3543
0
        CPL_TO_BOOL(
3544
0
            GDALDataTypeIsConversionLossy(eWrkDataType, eVRTBandDataType)) &&
3545
0
        !CPL_TO_BOOL(GDALDataTypeIsConversionLossy(eVRTBandDataType, eBufType));
3546
3547
0
    size_t idxBuffer = 0;
3548
0
    for (int iY = 0; iY < nOutYSize; iY++)
3549
0
    {
3550
0
        GByte *pDstLocation = static_cast<GByte *>(pData) +
3551
0
                              static_cast<GPtrDiff_t>(nLineSpace) * iY;
3552
3553
0
        for (int iX = 0; iX < nOutXSize;
3554
0
             iX++, pDstLocation += nPixelSpace, idxBuffer++)
3555
0
        {
3556
0
            WorkingDT afResult[2];
3557
0
            if (pafData && !bIsComplex)
3558
0
            {
3559
0
                WorkingDT fResult = pafData[idxBuffer];
3560
0
                if (bNoDataSetIsNan && std::isnan(fResult))
3561
0
                    continue;
3562
0
                if (bNoDataSetAndNotNan &&
3563
0
                    ARE_REAL_EQUAL(fResult, fWorkingDataTypeNoData))
3564
0
                    continue;
3565
0
                if (pabyMask && pabyMask[idxBuffer] == 0)
3566
0
                    continue;
3567
3568
0
                if (poColorTable)
3569
0
                {
3570
0
                    const GDALColorEntry *poEntry =
3571
0
                        poColorTable->GetColorEntry(static_cast<int>(fResult));
3572
0
                    if (poEntry)
3573
0
                    {
3574
0
                        if (m_nColorTableComponent == 1)
3575
0
                            fResult = poEntry->c1;
3576
0
                        else if (m_nColorTableComponent == 2)
3577
0
                            fResult = poEntry->c2;
3578
0
                        else if (m_nColorTableComponent == 3)
3579
0
                            fResult = poEntry->c3;
3580
0
                        else if (m_nColorTableComponent == 4)
3581
0
                            fResult = poEntry->c4;
3582
0
                    }
3583
0
                    else
3584
0
                    {
3585
0
                        CPLErrorOnce(CE_Failure, CPLE_AppDefined,
3586
0
                                     "No entry %d.", static_cast<int>(fResult));
3587
0
                        continue;
3588
0
                    }
3589
0
                }
3590
3591
0
                if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_LINEAR) != 0)
3592
0
                {
3593
0
                    fResult = static_cast<WorkingDT>(fResult * m_dfScaleRatio +
3594
0
                                                     m_dfScaleOff);
3595
0
                }
3596
0
                else if ((m_nProcessingFlags &
3597
0
                          PROCESSING_FLAG_SCALING_EXPONENTIAL) != 0)
3598
0
                {
3599
0
                    if (!m_bSrcMinMaxDefined)
3600
0
                    {
3601
0
                        int bSuccessMin = FALSE;
3602
0
                        int bSuccessMax = FALSE;
3603
0
                        double adfMinMax[2] = {
3604
0
                            poSourceBand->GetMinimum(&bSuccessMin),
3605
0
                            poSourceBand->GetMaximum(&bSuccessMax)};
3606
0
                        if ((bSuccessMin && bSuccessMax) ||
3607
0
                            poSourceBand->ComputeRasterMinMax(
3608
0
                                TRUE, adfMinMax) == CE_None)
3609
0
                        {
3610
0
                            m_dfSrcMin = adfMinMax[0];
3611
0
                            m_dfSrcMax = adfMinMax[1];
3612
0
                            m_bSrcMinMaxDefined = true;
3613
0
                        }
3614
0
                        else
3615
0
                        {
3616
0
                            CPLError(CE_Failure, CPLE_AppDefined,
3617
0
                                     "Cannot determine source min/max value");
3618
0
                            return CE_Failure;
3619
0
                        }
3620
0
                    }
3621
3622
0
                    double dfPowVal = (m_dfSrcMin == m_dfSrcMax)
3623
0
                                          ? 0
3624
0
                                          : (fResult - m_dfSrcMin) /
3625
0
                                                (m_dfSrcMax - m_dfSrcMin);
3626
0
                    if (m_bClip)
3627
0
                    {
3628
0
                        if (dfPowVal < 0.0)
3629
0
                            dfPowVal = 0.0;
3630
0
                        else if (dfPowVal > 1.0)
3631
0
                            dfPowVal = 1.0;
3632
0
                    }
3633
0
                    fResult =
3634
0
                        static_cast<WorkingDT>((m_dfDstMax - m_dfDstMin) *
3635
0
                                                   pow(dfPowVal, m_dfExponent) +
3636
0
                                               m_dfDstMin);
3637
0
                }
3638
3639
0
                if (!m_adfLUTInputs.empty())
3640
0
                    fResult = static_cast<WorkingDT>(LookupValue(fResult));
3641
3642
0
                if (m_nMaxValue != 0 && fResult > m_nMaxValue)
3643
0
                    fResult = static_cast<WorkingDT>(m_nMaxValue);
3644
3645
0
                afResult[0] = fResult;
3646
0
                afResult[1] = 0;
3647
0
            }
3648
0
            else if (pafData && bIsComplex)
3649
0
            {
3650
0
                afResult[0] = pafData[2 * idxBuffer];
3651
0
                afResult[1] = pafData[2 * idxBuffer + 1];
3652
3653
                // Do not use color table.
3654
0
                if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_LINEAR) != 0)
3655
0
                {
3656
0
                    afResult[0] = static_cast<WorkingDT>(
3657
0
                        afResult[0] * m_dfScaleRatio + m_dfScaleOff);
3658
0
                    afResult[1] = static_cast<WorkingDT>(
3659
0
                        afResult[1] * m_dfScaleRatio + m_dfScaleOff);
3660
0
                }
3661
3662
                /* Do not use LUT */
3663
0
            }
3664
0
            else
3665
0
            {
3666
0
                afResult[0] = static_cast<WorkingDT>(m_dfScaleOff);
3667
0
                afResult[1] = 0;
3668
3669
0
                if (!m_adfLUTInputs.empty())
3670
0
                    afResult[0] =
3671
0
                        static_cast<WorkingDT>(LookupValue(afResult[0]));
3672
3673
0
                if (m_nMaxValue != 0 && afResult[0] > m_nMaxValue)
3674
0
                    afResult[0] = static_cast<WorkingDT>(m_nMaxValue);
3675
0
            }
3676
3677
0
            if (eBufType == GDT_Byte && eVRTBandDataType == GDT_Byte)
3678
0
            {
3679
0
                *pDstLocation = static_cast<GByte>(std::min(
3680
0
                    255.0f,
3681
0
                    std::max(0.0f, static_cast<float>(afResult[0]) + 0.5f)));
3682
0
            }
3683
0
            else if (!bTwoStepDataTypeConversion)
3684
0
            {
3685
0
                CopyWord<WorkingDT>(afResult, eWrkDataType, pDstLocation,
3686
0
                                    eBufType);
3687
0
            }
3688
0
            else if (eVRTBandDataType == GDT_Float32 && eBufType == GDT_Float64)
3689
0
            {
3690
                // Particular case of the below 2-step conversion.
3691
                // Helps a bit for some geolocation based warping with Sentinel3
3692
                // data where the longitude/latitude arrays are Int32 bands,
3693
                // rescaled in VRT as Float32 and requested as Float64
3694
0
                float fVal;
3695
0
                GDALCopyWord(afResult[0], fVal);
3696
0
                *reinterpret_cast<double *>(pDstLocation) = fVal;
3697
0
            }
3698
0
            else
3699
0
            {
3700
0
                GByte abyTemp[2 * sizeof(double)];
3701
                // Convert first to the VRTRasterBand data type
3702
                // to get its clamping, before outputting to buffer data type
3703
0
                CopyWord<WorkingDT>(afResult, eWrkDataType, abyTemp,
3704
0
                                    eVRTBandDataType);
3705
0
                GDALCopyWords(abyTemp, eVRTBandDataType, 0, pDstLocation,
3706
0
                              eBufType, 0, 1);
3707
0
            }
3708
0
        }
3709
0
    }
3710
3711
0
    return CE_None;
3712
0
}
Unexecuted instantiation: CPLErr VRTComplexSource::RasterIOInternal<float>(GDALRasterBand*, GDALDataType, int, int, int, int, void*, int, int, GDALDataType, long long, long long, GDALRasterIOExtraArg*, GDALDataType, VRTSource::WorkingState&)
Unexecuted instantiation: CPLErr VRTComplexSource::RasterIOInternal<double>(GDALRasterBand*, GDALDataType, int, int, int, int, void*, int, int, GDALDataType, long long, long long, GDALRasterIOExtraArg*, GDALDataType, VRTSource::WorkingState&)
3713
3714
// Explicitly instantiate template method, as it is used in another file.
3715
template CPLErr VRTComplexSource::RasterIOInternal<float>(
3716
    GDALRasterBand *poSourceBand, GDALDataType eVRTBandDataType, int nReqXOff,
3717
    int nReqYOff, int nReqXSize, int nReqYSize, void *pData, int nOutXSize,
3718
    int nOutYSize, GDALDataType eBufType, GSpacing nPixelSpace,
3719
    GSpacing nLineSpace, GDALRasterIOExtraArg *psExtraArg,
3720
    GDALDataType eWrkDataType, WorkingState &oWorkingState);
3721
3722
/************************************************************************/
3723
/*                        AreValuesUnchanged()                          */
3724
/************************************************************************/
3725
3726
bool VRTComplexSource::AreValuesUnchanged() const
3727
0
{
3728
0
    return m_dfScaleOff == 0.0 && m_dfScaleRatio == 1.0 &&
3729
0
           m_adfLUTInputs.empty() && m_nColorTableComponent == 0 &&
3730
0
           (m_nProcessingFlags & PROCESSING_FLAG_SCALING_EXPONENTIAL) == 0;
3731
0
}
3732
3733
/************************************************************************/
3734
/*                             GetMinimum()                             */
3735
/************************************************************************/
3736
3737
double VRTComplexSource::GetMinimum(int nXSize, int nYSize, int *pbSuccess)
3738
0
{
3739
0
    if (AreValuesUnchanged())
3740
0
    {
3741
0
        return VRTSimpleSource::GetMinimum(nXSize, nYSize, pbSuccess);
3742
0
    }
3743
3744
0
    *pbSuccess = FALSE;
3745
0
    return 0;
3746
0
}
3747
3748
/************************************************************************/
3749
/*                             GetMaximum()                             */
3750
/************************************************************************/
3751
3752
double VRTComplexSource::GetMaximum(int nXSize, int nYSize, int *pbSuccess)
3753
0
{
3754
0
    if (AreValuesUnchanged())
3755
0
    {
3756
0
        return VRTSimpleSource::GetMaximum(nXSize, nYSize, pbSuccess);
3757
0
    }
3758
3759
0
    *pbSuccess = FALSE;
3760
0
    return 0;
3761
0
}
3762
3763
/************************************************************************/
3764
/*                            GetHistogram()                            */
3765
/************************************************************************/
3766
3767
CPLErr VRTComplexSource::GetHistogram(int nXSize, int nYSize, double dfMin,
3768
                                      double dfMax, int nBuckets,
3769
                                      GUIntBig *panHistogram,
3770
                                      int bIncludeOutOfRange, int bApproxOK,
3771
                                      GDALProgressFunc pfnProgress,
3772
                                      void *pProgressData)
3773
0
{
3774
0
    if (AreValuesUnchanged())
3775
0
    {
3776
0
        return VRTSimpleSource::GetHistogram(
3777
0
            nXSize, nYSize, dfMin, dfMax, nBuckets, panHistogram,
3778
0
            bIncludeOutOfRange, bApproxOK, pfnProgress, pProgressData);
3779
0
    }
3780
3781
0
    return CE_Failure;
3782
0
}
3783
3784
/************************************************************************/
3785
/* ==================================================================== */
3786
/*                          VRTFuncSource                               */
3787
/* ==================================================================== */
3788
/************************************************************************/
3789
3790
/************************************************************************/
3791
/*                           VRTFuncSource()                            */
3792
/************************************************************************/
3793
3794
VRTFuncSource::VRTFuncSource()
3795
0
    : pfnReadFunc(nullptr), pCBData(nullptr), eType(GDT_Byte),
3796
0
      fNoDataValue(static_cast<float>(VRT_NODATA_UNSET))
3797
0
{
3798
0
}
3799
3800
/************************************************************************/
3801
/*                           ~VRTFuncSource()                           */
3802
/************************************************************************/
3803
3804
VRTFuncSource::~VRTFuncSource()
3805
{
3806
}
3807
3808
/************************************************************************/
3809
/*                            GetType()                                 */
3810
/************************************************************************/
3811
3812
const char *VRTFuncSource::GetType() const
3813
0
{
3814
0
    static const char *TYPE = "FuncSource";
3815
0
    return TYPE;
3816
0
}
3817
3818
/************************************************************************/
3819
/*                           SerializeToXML()                           */
3820
/************************************************************************/
3821
3822
CPLXMLNode *VRTFuncSource::SerializeToXML(CPL_UNUSED const char *pszVRTPath)
3823
0
{
3824
0
    return nullptr;
3825
0
}
3826
3827
/************************************************************************/
3828
/*                              RasterIO()                              */
3829
/************************************************************************/
3830
3831
CPLErr VRTFuncSource::RasterIO(GDALDataType /*eVRTBandDataType*/, int nXOff,
3832
                               int nYOff, int nXSize, int nYSize, void *pData,
3833
                               int nBufXSize, int nBufYSize,
3834
                               GDALDataType eBufType, GSpacing nPixelSpace,
3835
                               GSpacing nLineSpace,
3836
                               GDALRasterIOExtraArg * /* psExtraArg */,
3837
                               WorkingState & /* oWorkingState */)
3838
0
{
3839
0
    if (nPixelSpace == GDALGetDataTypeSizeBytes(eBufType) &&
3840
0
        nLineSpace == nPixelSpace * nXSize && nBufXSize == nXSize &&
3841
0
        nBufYSize == nYSize && eBufType == eType)
3842
0
    {
3843
0
        return pfnReadFunc(pCBData, nXOff, nYOff, nXSize, nYSize, pData);
3844
0
    }
3845
0
    else
3846
0
    {
3847
0
        CPLError(CE_Failure, CPLE_AppDefined,
3848
0
                 "VRTFuncSource::RasterIO() - Irregular request.");
3849
0
        CPLDebug("VRT", "Irregular request: %d,%d  %d,%d, %d,%d %d,%d %d,%d",
3850
0
                 static_cast<int>(nPixelSpace),
3851
0
                 GDALGetDataTypeSizeBytes(eBufType),
3852
0
                 static_cast<int>(nLineSpace),
3853
0
                 static_cast<int>(nPixelSpace) * nXSize, nBufXSize, nXSize,
3854
0
                 nBufYSize, nYSize, static_cast<int>(eBufType),
3855
0
                 static_cast<int>(eType));
3856
3857
0
        return CE_Failure;
3858
0
    }
3859
0
}
3860
3861
/************************************************************************/
3862
/*                             GetMinimum()                             */
3863
/************************************************************************/
3864
3865
double VRTFuncSource::GetMinimum(int /* nXSize */, int /* nYSize */,
3866
                                 int *pbSuccess)
3867
0
{
3868
0
    *pbSuccess = FALSE;
3869
0
    return 0;
3870
0
}
3871
3872
/************************************************************************/
3873
/*                             GetMaximum()                             */
3874
/************************************************************************/
3875
3876
double VRTFuncSource::GetMaximum(int /* nXSize */, int /* nYSize */,
3877
                                 int *pbSuccess)
3878
0
{
3879
0
    *pbSuccess = FALSE;
3880
0
    return 0;
3881
0
}
3882
3883
/************************************************************************/
3884
/*                            GetHistogram()                            */
3885
/************************************************************************/
3886
3887
CPLErr VRTFuncSource::GetHistogram(
3888
    int /* nXSize */, int /* nYSize */, double /* dfMin */, double /* dfMax */,
3889
    int /* nBuckets */, GUIntBig * /* panHistogram */,
3890
    int /* bIncludeOutOfRange */, int /* bApproxOK */,
3891
    GDALProgressFunc /* pfnProgress */, void * /* pProgressData */)
3892
0
{
3893
0
    return CE_Failure;
3894
0
}
3895
3896
/************************************************************************/
3897
/*                        VRTParseCoreSources()                         */
3898
/************************************************************************/
3899
3900
VRTSource *VRTParseCoreSources(const CPLXMLNode *psChild,
3901
                               const char *pszVRTPath,
3902
                               VRTMapSharedResources &oMapSharedSources)
3903
3904
0
{
3905
0
    VRTSource *poSource = nullptr;
3906
3907
0
    if (EQUAL(psChild->pszValue, VRTAveragedSource::GetTypeStatic()) ||
3908
0
        (EQUAL(psChild->pszValue, VRTSimpleSource::GetTypeStatic()) &&
3909
0
         STARTS_WITH_CI(CPLGetXMLValue(psChild, "Resampling", "Nearest"),
3910
0
                        "Aver")))
3911
0
    {
3912
0
        poSource = new VRTAveragedSource();
3913
0
    }
3914
0
    else if (EQUAL(psChild->pszValue, VRTSimpleSource::GetTypeStatic()))
3915
0
    {
3916
0
        poSource = new VRTSimpleSource();
3917
0
    }
3918
0
    else if (EQUAL(psChild->pszValue, VRTComplexSource::GetTypeStatic()))
3919
0
    {
3920
0
        poSource = new VRTComplexSource();
3921
0
    }
3922
0
    else if (EQUAL(psChild->pszValue, VRTNoDataFromMaskSource::GetTypeStatic()))
3923
0
    {
3924
0
        poSource = new VRTNoDataFromMaskSource();
3925
0
    }
3926
0
    else
3927
0
    {
3928
0
        CPLError(CE_Failure, CPLE_AppDefined,
3929
0
                 "VRTParseCoreSources() - Unknown source : %s",
3930
0
                 psChild->pszValue);
3931
0
        return nullptr;
3932
0
    }
3933
3934
0
    if (poSource->XMLInit(psChild, pszVRTPath, oMapSharedSources) == CE_None)
3935
0
        return poSource;
3936
3937
0
    delete poSource;
3938
0
    return nullptr;
3939
0
}
3940
3941
/*! @endcond */