Coverage Report

Created: 2025-11-16 06:25

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