• Home
  • Features
  • Pricing
  • Docs
  • Announcements
  • Sign In

OSGeo / gdal / 13836648005

13 Mar 2025 02:09PM UTC coverage: 70.436% (-0.01%) from 70.446%
13836648005

push

github

web-flow
New Transform type: Homography (#11949)

Add new transform type, Homography.
Add functions to compute homography from a list of GCPs.
Add functions to serialize and deserialize a homography
Automatically select homography transfrom when there are 4 or 5 GCPs present.

Fixes #11940

231 of 274 new or added lines in 2 files covered. (84.31%)

16257 existing lines in 42 files now uncovered.

553736 of 786159 relevant lines covered (70.44%)

221595.72 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

85.92
/frmts/vrt/vrtsources.cpp
1
/******************************************************************************
2
 *
3
 * Project:  Virtual GDAL Datasets
4
 * Purpose:  Implementation of VRTSimpleSource, VRTFuncSource and
5
 *           VRTAveragedSource.
6
 * Author:   Frank Warmerdam <warmerdam@pobox.com>
7
 *
8
 ******************************************************************************
9
 * Copyright (c) 2001, Frank Warmerdam <warmerdam@pobox.com>
10
 * Copyright (c) 2008-2013, Even Rouault <even dot rouault at spatialys.com>
11
 *
12
 * SPDX-License-Identifier: MIT
13
 ****************************************************************************/
14

15
#include "gdal_vrt.h"
16
#include "vrtdataset.h"
17

18
#include <cassert>
19
#include <climits>
20
#include <cmath>
21
#include <cstddef>
22
#include <cstdio>
23
#include <cstdlib>
24
#include <cstring>
25
#include <algorithm>
26
#include <limits>
27
#include <string>
28

29
#include "cpl_conv.h"
30
#include "cpl_error.h"
31
#include "cpl_hash_set.h"
32
#include "cpl_minixml.h"
33
#include "cpl_progress.h"
34
#include "cpl_string.h"
35
#include "cpl_vsi.h"
36
#include "gdal.h"
37
#include "gdal_priv.h"
38
#include "gdal_proxy.h"
39
#include "gdal_priv_templates.hpp"
40

41
/*! @cond Doxygen_Suppress */
42

43
// #define DEBUG_VERBOSE 1
44

45
// See #5459
46
#ifdef isnan
47
#define HAS_ISNAN_MACRO
48
#endif
49
#include <algorithm>
50
#if defined(HAS_ISNAN_MACRO) && !defined(isnan)
51
#define isnan std::isnan
52
#endif
53

54
/************************************************************************/
55
/* ==================================================================== */
56
/*                             VRTSource                                */
57
/* ==================================================================== */
58
/************************************************************************/
59

60
VRTSource::~VRTSource()
13,424✔
61
{
62
}
13,424✔
63

64
/************************************************************************/
65
/*                             GetFileList()                            */
66
/************************************************************************/
67

68
void VRTSource::GetFileList(char *** /* ppapszFileList */, int * /* pnSize */,
×
69
                            int * /* pnMaxSize */, CPLHashSet * /* hSetFiles */)
70
{
71
}
×
72

73
/************************************************************************/
74
/* ==================================================================== */
75
/*                          VRTSimpleSource                             */
76
/* ==================================================================== */
77
/************************************************************************/
78

79
/************************************************************************/
80
/*                          VRTSimpleSource()                           */
81
/************************************************************************/
82

83
VRTSimpleSource::VRTSimpleSource() = default;
84

85
/************************************************************************/
86
/*                          VRTSimpleSource()                           */
87
/************************************************************************/
88

89
VRTSimpleSource::VRTSimpleSource(const VRTSimpleSource *poSrcSource,
76✔
90
                                 double dfXDstRatio, double dfYDstRatio)
76✔
91
    : m_poMapSharedSources(poSrcSource->m_poMapSharedSources),
76✔
92
      m_poRasterBand(poSrcSource->m_poRasterBand),
76✔
93
      m_poMaskBandMainBand(poSrcSource->m_poMaskBandMainBand),
76✔
94
      m_aosOpenOptionsOri(poSrcSource->m_aosOpenOptionsOri),
76✔
95
      m_aosOpenOptions(poSrcSource->m_aosOpenOptions),
76✔
96
      m_bSrcDSNameFromVRT(poSrcSource->m_bSrcDSNameFromVRT),
76✔
97
      m_nBand(poSrcSource->m_nBand),
76✔
98
      m_bGetMaskBand(poSrcSource->m_bGetMaskBand),
76✔
99
      m_dfSrcXOff(poSrcSource->m_dfSrcXOff),
76✔
100
      m_dfSrcYOff(poSrcSource->m_dfSrcYOff),
76✔
101
      m_dfSrcXSize(poSrcSource->m_dfSrcXSize),
76✔
102
      m_dfSrcYSize(poSrcSource->m_dfSrcYSize),
76✔
103
      m_nMaxValue(poSrcSource->m_nMaxValue), m_bRelativeToVRTOri(-1),
76✔
104
      m_nExplicitSharedStatus(poSrcSource->m_nExplicitSharedStatus),
76✔
105
      m_osSrcDSName(poSrcSource->m_osSrcDSName),
76✔
106
      m_bDropRefOnSrcBand(poSrcSource->m_bDropRefOnSrcBand)
76✔
107
{
108
    if (!poSrcSource->IsSrcWinSet() && !poSrcSource->IsDstWinSet() &&
76✔
109
        (dfXDstRatio != 1.0 || dfYDstRatio != 1.0))
×
110
    {
111
        auto l_band = GetRasterBand();
1✔
112
        if (l_band)
1✔
113
        {
114
            m_dfSrcXOff = 0;
1✔
115
            m_dfSrcYOff = 0;
1✔
116
            m_dfSrcXSize = l_band->GetXSize();
1✔
117
            m_dfSrcYSize = l_band->GetYSize();
1✔
118
            m_dfDstXOff = 0;
1✔
119
            m_dfDstYOff = 0;
1✔
120
            m_dfDstXSize = l_band->GetXSize() * dfXDstRatio;
1✔
121
            m_dfDstYSize = l_band->GetYSize() * dfYDstRatio;
1✔
122
        }
123
    }
124
    else if (poSrcSource->IsDstWinSet())
75✔
125
    {
126
        m_dfDstXOff = poSrcSource->m_dfDstXOff * dfXDstRatio;
75✔
127
        m_dfDstYOff = poSrcSource->m_dfDstYOff * dfYDstRatio;
75✔
128
        m_dfDstXSize = poSrcSource->m_dfDstXSize * dfXDstRatio;
75✔
129
        m_dfDstYSize = poSrcSource->m_dfDstYSize * dfYDstRatio;
75✔
130
    }
131
}
76✔
132

133
/************************************************************************/
134
/*                          ~VRTSimpleSource()                          */
135
/************************************************************************/
136

137
VRTSimpleSource::~VRTSimpleSource()
28,126✔
138

139
{
140
    if (!m_bDropRefOnSrcBand)
13,380✔
141
        return;
396✔
142

143
    if (m_poMaskBandMainBand != nullptr)
12,984✔
144
    {
145
        if (m_poMaskBandMainBand->GetDataset() != nullptr)
51✔
146
        {
147
            m_poMaskBandMainBand->GetDataset()->ReleaseRef();
51✔
148
        }
149
    }
150
    else if (m_poRasterBand != nullptr &&
23,328✔
151
             m_poRasterBand->GetDataset() != nullptr)
10,395✔
152
    {
153
        m_poRasterBand->GetDataset()->ReleaseRef();
10,395✔
154
    }
155
}
26,146✔
156

157
/************************************************************************/
158
/*                           GetTypeStatic()                            */
159
/************************************************************************/
160

161
const char *VRTSimpleSource::GetTypeStatic()
35,723✔
162
{
163
    static const char *TYPE = "SimpleSource";
164
    return TYPE;
35,723✔
165
}
166

167
/************************************************************************/
168
/*                            GetType()                                 */
169
/************************************************************************/
170

171
const char *VRTSimpleSource::GetType() const
16,199✔
172
{
173
    return GetTypeStatic();
16,199✔
174
}
175

176
/************************************************************************/
177
/*                           FlushCache()                               */
178
/************************************************************************/
179

180
CPLErr VRTSimpleSource::FlushCache(bool bAtClosing)
13,427✔
181

182
{
183
    if (m_poMaskBandMainBand != nullptr)
13,427✔
184
    {
185
        return m_poMaskBandMainBand->FlushCache(bAtClosing);
5✔
186
    }
187
    else if (m_poRasterBand != nullptr)
13,422✔
188
    {
189
        return m_poRasterBand->FlushCache(bAtClosing);
10,860✔
190
    }
191
    return CE_None;
2,562✔
192
}
193

194
/************************************************************************/
195
/*                    UnsetPreservedRelativeFilenames()                 */
196
/************************************************************************/
197

198
void VRTSimpleSource::UnsetPreservedRelativeFilenames()
18✔
199
{
200
    if (!STARTS_WITH(m_osSourceFileNameOri.c_str(), "http://") &&
35✔
201
        !STARTS_WITH(m_osSourceFileNameOri.c_str(), "https://"))
17✔
202
    {
203
        m_bRelativeToVRTOri = -1;
17✔
204
        m_osSourceFileNameOri = "";
17✔
205
    }
206
}
18✔
207

208
/************************************************************************/
209
/*                             SetSrcBand()                             */
210
/************************************************************************/
211

212
void VRTSimpleSource::SetSrcBand(const char *pszFilename, int nBand)
447✔
213

214
{
215
    m_nBand = nBand;
447✔
216
    m_osSrcDSName = pszFilename;
447✔
217
}
447✔
218

219
/************************************************************************/
220
/*                             SetSrcBand()                             */
221
/************************************************************************/
222

223
void VRTSimpleSource::SetSrcBand(GDALRasterBand *poNewSrcBand)
9,310✔
224

225
{
226
    m_poRasterBand = poNewSrcBand;
9,310✔
227
    m_nBand = m_poRasterBand->GetBand();
9,310✔
228
    auto poDS = poNewSrcBand->GetDataset();
9,310✔
229
    if (poDS != nullptr)
9,310✔
230
    {
231
        m_osSrcDSName = poDS->GetDescription();
9,310✔
232
        m_aosOpenOptions = CSLDuplicate(poDS->GetOpenOptions());
9,310✔
233
        m_aosOpenOptionsOri = m_aosOpenOptions;
9,310✔
234
    }
235
}
9,310✔
236

237
/************************************************************************/
238
/*                          SetSrcMaskBand()                            */
239
/************************************************************************/
240

241
// poSrcBand is not the mask band, but the band from which the mask band is
242
// taken.
243
void VRTSimpleSource::SetSrcMaskBand(GDALRasterBand *poNewSrcBand)
32✔
244

245
{
246
    m_poRasterBand = poNewSrcBand->GetMaskBand();
32✔
247
    m_poMaskBandMainBand = poNewSrcBand;
32✔
248
    m_nBand = poNewSrcBand->GetBand();
32✔
249
    auto poDS = poNewSrcBand->GetDataset();
32✔
250
    if (poDS != nullptr)
32✔
251
    {
252
        m_osSrcDSName = poDS->GetDescription();
32✔
253
        m_aosOpenOptions = CSLDuplicate(poDS->GetOpenOptions());
32✔
254
        m_aosOpenOptionsOri = m_aosOpenOptions;
32✔
255
    }
256
    m_bGetMaskBand = true;
32✔
257
}
32✔
258

259
/************************************************************************/
260
/*                         RoundIfCloseToInt()                          */
261
/************************************************************************/
262

263
static double RoundIfCloseToInt(double dfValue)
301,472✔
264
{
265
    double dfClosestInt = floor(dfValue + 0.5);
301,472✔
266
    return (fabs(dfValue - dfClosestInt) < 1e-3) ? dfClosestInt : dfValue;
301,472✔
267
}
268

269
/************************************************************************/
270
/*                            SetSrcWindow()                            */
271
/************************************************************************/
272

273
void VRTSimpleSource::SetSrcWindow(double dfNewXOff, double dfNewYOff,
12,755✔
274
                                   double dfNewXSize, double dfNewYSize)
275

276
{
277
    m_dfSrcXOff = RoundIfCloseToInt(dfNewXOff);
12,755✔
278
    m_dfSrcYOff = RoundIfCloseToInt(dfNewYOff);
12,755✔
279
    m_dfSrcXSize = RoundIfCloseToInt(dfNewXSize);
12,755✔
280
    m_dfSrcYSize = RoundIfCloseToInt(dfNewYSize);
12,755✔
281
}
12,755✔
282

283
/************************************************************************/
284
/*                            SetDstWindow()                            */
285
/************************************************************************/
286

287
void VRTSimpleSource::SetDstWindow(double dfNewXOff, double dfNewYOff,
12,754✔
288
                                   double dfNewXSize, double dfNewYSize)
289

290
{
291
    m_dfDstXOff = RoundIfCloseToInt(dfNewXOff);
12,754✔
292
    m_dfDstYOff = RoundIfCloseToInt(dfNewYOff);
12,754✔
293
    m_dfDstXSize = RoundIfCloseToInt(dfNewXSize);
12,754✔
294
    m_dfDstYSize = RoundIfCloseToInt(dfNewYSize);
12,754✔
295
}
12,754✔
296

297
/************************************************************************/
298
/*                            GetDstWindow()                            */
299
/************************************************************************/
300

301
void VRTSimpleSource::GetDstWindow(double &dfDstXOff, double &dfDstYOff,
516✔
302
                                   double &dfDstXSize, double &dfDstYSize) const
303
{
304
    dfDstXOff = m_dfDstXOff;
516✔
305
    dfDstYOff = m_dfDstYOff;
516✔
306
    dfDstXSize = m_dfDstXSize;
516✔
307
    dfDstYSize = m_dfDstYSize;
516✔
308
}
516✔
309

310
/************************************************************************/
311
/*                        DstWindowIntersects()                         */
312
/************************************************************************/
313

314
bool VRTSimpleSource::DstWindowIntersects(double dfXOff, double dfYOff,
1,072✔
315
                                          double dfXSize, double dfYSize) const
316
{
317
    return IsDstWinSet() && m_dfDstXOff + m_dfDstXSize > dfXOff &&
2,135✔
318
           m_dfDstYOff + m_dfDstYSize > dfYOff &&
1,063✔
319
           m_dfDstXOff < dfXOff + dfXSize && m_dfDstYOff < dfYOff + dfYSize;
2,135✔
320
}
321

322
/************************************************************************/
323
/*                            IsSlowSource()                            */
324
/************************************************************************/
325

326
static bool IsSlowSource(const char *pszSrcName)
2,501✔
327
{
328
    return strstr(pszSrcName, "/vsicurl/http") != nullptr ||
5,002✔
329
           strstr(pszSrcName, "/vsicurl/ftp") != nullptr ||
5,002✔
330
           (strstr(pszSrcName, "/vsicurl?") != nullptr &&
2,501✔
331
            strstr(pszSrcName, "&url=http") != nullptr);
2,501✔
332
}
333

334
/************************************************************************/
335
/*                         AddSourceFilenameNode()                      */
336
/************************************************************************/
337

338
void VRTSimpleSource::AddSourceFilenameNode(const char *pszVRTPath,
2,534✔
339
                                            CPLXMLNode *psSrc)
340
{
341

342
    VSIStatBufL sStat;
343
    int bRelativeToVRT = FALSE;  // TODO(schwehr): Make this a bool?
2,534✔
344
    std::string osSourceFilename;
5,068✔
345

346
    if (m_bRelativeToVRTOri >= 0)
2,534✔
347
    {
348
        osSourceFilename = m_osSourceFileNameOri;
33✔
349
        bRelativeToVRT = m_bRelativeToVRTOri;
33✔
350
    }
351
    else if (IsSlowSource(m_osSrcDSName))
2,501✔
352
    {
353
        // Testing the existence of remote resources can be excruciating
354
        // slow, so let's just suppose they exist.
355
        osSourceFilename = m_osSrcDSName;
×
356
        bRelativeToVRT = FALSE;
×
357
    }
358
    // If this isn't actually a file, don't even try to know if it is a
359
    // relative path. It can't be !, and unfortunately CPLIsFilenameRelative()
360
    // can only work with strings that are filenames To be clear
361
    // NITF_TOC_ENTRY:CADRG_JOG-A_250K_1_0:some_path isn't a relative file
362
    // path.
363
    else if (VSIStatExL(m_osSrcDSName, &sStat, VSI_STAT_EXISTS_FLAG) != 0)
2,501✔
364
    {
365
        osSourceFilename = m_osSrcDSName;
58✔
366
        bRelativeToVRT = FALSE;
58✔
367

368
        // Try subdatasetinfo API first
369
        // Note: this will become the only branch when subdatasetinfo will become
370
        //       available for NITF_IM, RASTERLITE and TILEDB
371
        const auto oSubDSInfo{GDALGetSubdatasetInfo(osSourceFilename.c_str())};
58✔
372
        if (oSubDSInfo && !oSubDSInfo->GetPathComponent().empty())
58✔
373
        {
374
            auto path{oSubDSInfo->GetPathComponent()};
10✔
375
            std::string relPath{CPLExtractRelativePath(pszVRTPath, path.c_str(),
376
                                                       &bRelativeToVRT)};
10✔
377
            osSourceFilename = oSubDSInfo->ModifyPathComponent(relPath);
5✔
378
            GDALDestroySubdatasetInfo(oSubDSInfo);
5✔
379
        }
380
        else
381
        {
382
            for (const char *pszSyntax : VRTDataset::apszSpecialSyntax)
303✔
383
            {
384
                CPLString osPrefix(pszSyntax);
253✔
385
                osPrefix.resize(strchr(pszSyntax, ':') - pszSyntax + 1);
253✔
386
                if (pszSyntax[osPrefix.size()] == '"')
253✔
387
                    osPrefix += '"';
50✔
388
                if (EQUALN(osSourceFilename.c_str(), osPrefix, osPrefix.size()))
253✔
389
                {
390
                    if (STARTS_WITH_CI(pszSyntax + osPrefix.size(), "{ANY}"))
3✔
391
                    {
392
                        const char *pszLastPart =
393
                            strrchr(osSourceFilename.c_str(), ':') + 1;
3✔
394
                        // CSV:z:/foo.xyz
395
                        if ((pszLastPart[0] == '/' || pszLastPart[0] == '\\') &&
1✔
396
                            pszLastPart - osSourceFilename.c_str() >= 3 &&
6✔
397
                            pszLastPart[-3] == ':')
2✔
398
                            pszLastPart -= 2;
×
399
                        CPLString osPrefixFilename(osSourceFilename);
3✔
400
                        osPrefixFilename.resize(pszLastPart -
3✔
401
                                                osSourceFilename.c_str());
3✔
402
                        osSourceFilename = CPLExtractRelativePath(
403
                            pszVRTPath, pszLastPart, &bRelativeToVRT);
3✔
404
                        osSourceFilename = osPrefixFilename + osSourceFilename;
3✔
405
                    }
406
                    else if (STARTS_WITH_CI(pszSyntax + osPrefix.size(),
×
407
                                            "{FILENAME}"))
408
                    {
409
                        CPLString osFilename(osSourceFilename.c_str() +
×
410
                                             osPrefix.size());
×
411
                        size_t nPos = 0;
×
412
                        if (osFilename.size() >= 3 && osFilename[1] == ':' &&
×
413
                            (osFilename[2] == '\\' || osFilename[2] == '/'))
×
414
                            nPos = 2;
×
415
                        nPos = osFilename.find(
×
416
                            pszSyntax[osPrefix.size() + strlen("{FILENAME}")],
×
417
                            nPos);
418
                        if (nPos != std::string::npos)
×
419
                        {
420
                            const CPLString osSuffix = osFilename.substr(nPos);
×
421
                            osFilename.resize(nPos);
×
422
                            osSourceFilename = CPLExtractRelativePath(
423
                                pszVRTPath, osFilename, &bRelativeToVRT);
×
424
                            osSourceFilename =
425
                                osPrefix + osSourceFilename + osSuffix;
×
426
                        }
427
                    }
428
                    break;
3✔
429
                }
430
            }
431
        }
432
    }
433
    else
434
    {
435
        std::string osVRTFilename = pszVRTPath;
4,886✔
436
        std::string osSourceDataset = m_osSrcDSName;
4,886✔
437
        char *pszCurDir = CPLGetCurrentDir();
2,443✔
438
        if (CPLIsFilenameRelative(osSourceDataset.c_str()) &&
2,443✔
439
            !CPLIsFilenameRelative(osVRTFilename.c_str()) &&
2,443✔
440
            pszCurDir != nullptr)
441
        {
442
            osSourceDataset = CPLFormFilenameSafe(
128✔
443
                pszCurDir, osSourceDataset.c_str(), nullptr);
64✔
444
        }
445
        else if (!CPLIsFilenameRelative(osSourceDataset.c_str()) &&
2,379✔
446
                 CPLIsFilenameRelative(osVRTFilename.c_str()) &&
2,379✔
447
                 pszCurDir != nullptr)
448
        {
449
            osVRTFilename =
450
                CPLFormFilenameSafe(pszCurDir, osVRTFilename.c_str(), nullptr);
72✔
451
        }
452
        CPLFree(pszCurDir);
2,443✔
453
        osSourceFilename = CPLExtractRelativePath(
454
            osVRTFilename.c_str(), osSourceDataset.c_str(), &bRelativeToVRT);
2,443✔
455
    }
456

457
    CPLSetXMLValue(psSrc, "SourceFilename", osSourceFilename.c_str());
2,534✔
458

459
    CPLCreateXMLNode(CPLCreateXMLNode(CPLGetXMLNode(psSrc, "SourceFilename"),
2,534✔
460
                                      CXT_Attribute, "relativeToVRT"),
461
                     CXT_Text, bRelativeToVRT ? "1" : "0");
2,534✔
462

463
    // Determine if we must write the shared attribute. The config option
464
    // will override the m_nExplicitSharedStatus value
465
    const char *pszShared = CPLGetConfigOption("VRT_SHARED_SOURCE", nullptr);
2,534✔
466
    if ((pszShared == nullptr && m_nExplicitSharedStatus == 0) ||
2,534✔
467
        (pszShared != nullptr && !CPLTestBool(pszShared)))
×
468
    {
469
        CPLCreateXMLNode(
×
470
            CPLCreateXMLNode(CPLGetXMLNode(psSrc, "SourceFilename"),
471
                             CXT_Attribute, "shared"),
472
            CXT_Text, "0");
473
    }
474
}
2,534✔
475

476
/************************************************************************/
477
/*                           SerializeToXML()                           */
478
/************************************************************************/
479

480
CPLXMLNode *VRTSimpleSource::SerializeToXML(const char *pszVRTPath)
2,535✔
481

482
{
483
    CPLXMLNode *const psSrc =
484
        CPLCreateXMLNode(nullptr, CXT_Element, GetTypeStatic());
2,535✔
485

486
    if (!m_osResampling.empty())
2,535✔
487
    {
488
        CPLCreateXMLNode(CPLCreateXMLNode(psSrc, CXT_Attribute, "resampling"),
24✔
489
                         CXT_Text, m_osResampling.c_str());
490
    }
491

492
    if (!m_osName.empty())
2,535✔
493
    {
494
        CPLAddXMLAttributeAndValue(psSrc, "name", m_osName);
13✔
495
    }
496

497
    if (m_bSrcDSNameFromVRT)
2,535✔
498
    {
499
        CPLAddXMLChild(psSrc, CPLParseXMLString(m_osSrcDSName.c_str()));
1✔
500
    }
501
    else
502
    {
503
        AddSourceFilenameNode(pszVRTPath, psSrc);
2,534✔
504
    }
505

506
    GDALSerializeOpenOptionsToXML(psSrc, m_aosOpenOptionsOri.List());
2,535✔
507

508
    if (m_bGetMaskBand)
2,535✔
509
        CPLSetXMLValue(psSrc, "SourceBand", CPLSPrintf("mask,%d", m_nBand));
18✔
510
    else
511
        CPLSetXMLValue(psSrc, "SourceBand", CPLSPrintf("%d", m_nBand));
2,517✔
512

513
    // TODO: in a later version, no longer emit SourceProperties, which
514
    // is no longer used by GDAL 3.4
515
    if (m_poRasterBand)
2,535✔
516
    {
517
        /* Write a few additional useful properties of the dataset */
518
        /* so that we can use a proxy dataset when re-opening. See XMLInit() */
519
        /* below */
520
        CPLSetXMLValue(psSrc, "SourceProperties.#RasterXSize",
2,423✔
521
                       CPLSPrintf("%d", m_poRasterBand->GetXSize()));
2,423✔
522
        CPLSetXMLValue(psSrc, "SourceProperties.#RasterYSize",
2,423✔
523
                       CPLSPrintf("%d", m_poRasterBand->GetYSize()));
2,423✔
524
        CPLSetXMLValue(
2,423✔
525
            psSrc, "SourceProperties.#DataType",
526
            GDALGetDataTypeName(m_poRasterBand->GetRasterDataType()));
2,423✔
527

528
        int nBlockXSize = 0;
2,423✔
529
        int nBlockYSize = 0;
2,423✔
530
        m_poRasterBand->GetBlockSize(&nBlockXSize, &nBlockYSize);
2,423✔
531

532
        CPLSetXMLValue(psSrc, "SourceProperties.#BlockXSize",
2,423✔
533
                       CPLSPrintf("%d", nBlockXSize));
534
        CPLSetXMLValue(psSrc, "SourceProperties.#BlockYSize",
2,423✔
535
                       CPLSPrintf("%d", nBlockYSize));
536
    }
537

538
    if (IsSrcWinSet())
2,535✔
539
    {
540
        CPLSetXMLValue(psSrc, "SrcRect.#xOff",
2,522✔
541
                       CPLSPrintf("%.15g", m_dfSrcXOff));
542
        CPLSetXMLValue(psSrc, "SrcRect.#yOff",
2,522✔
543
                       CPLSPrintf("%.15g", m_dfSrcYOff));
544
        CPLSetXMLValue(psSrc, "SrcRect.#xSize",
2,522✔
545
                       CPLSPrintf("%.15g", m_dfSrcXSize));
546
        CPLSetXMLValue(psSrc, "SrcRect.#ySize",
2,522✔
547
                       CPLSPrintf("%.15g", m_dfSrcYSize));
548
    }
549

550
    if (IsDstWinSet())
2,535✔
551
    {
552
        CPLSetXMLValue(psSrc, "DstRect.#xOff",
2,522✔
553
                       CPLSPrintf("%.15g", m_dfDstXOff));
554
        CPLSetXMLValue(psSrc, "DstRect.#yOff",
2,522✔
555
                       CPLSPrintf("%.15g", m_dfDstYOff));
556
        CPLSetXMLValue(psSrc, "DstRect.#xSize",
2,522✔
557
                       CPLSPrintf("%.15g", m_dfDstXSize));
558
        CPLSetXMLValue(psSrc, "DstRect.#ySize",
2,522✔
559
                       CPLSPrintf("%.15g", m_dfDstYSize));
560
    }
561

562
    return psSrc;
2,535✔
563
}
564

565
/************************************************************************/
566
/*                              XMLInit()                               */
567
/************************************************************************/
568

569
CPLErr VRTSimpleSource::XMLInit(const CPLXMLNode *psSrc, const char *pszVRTPath,
3,169✔
570
                                VRTMapSharedResources &oMapSharedSources)
571

572
{
573
    m_poMapSharedSources = &oMapSharedSources;
3,169✔
574

575
    m_osResampling = CPLGetXMLValue(psSrc, "resampling", "");
3,169✔
576
    m_osName = CPLGetXMLValue(psSrc, "name", "");
3,169✔
577

578
    /* -------------------------------------------------------------------- */
579
    /*      Prepare filename.                                               */
580
    /* -------------------------------------------------------------------- */
581
    const CPLXMLNode *psSourceFileNameNode =
582
        CPLGetXMLNode(psSrc, "SourceFilename");
3,169✔
583
    const CPLXMLNode *psSourceVRTDataset = CPLGetXMLNode(psSrc, "VRTDataset");
3,169✔
584
    const char *pszFilename =
585
        psSourceFileNameNode ? CPLGetXMLValue(psSourceFileNameNode, nullptr, "")
3,169✔
586
                             : "";
3,169✔
587

588
    if (pszFilename[0] == '\0' && !psSourceVRTDataset)
3,169✔
589
    {
590
        CPLError(CE_Warning, CPLE_AppDefined,
1✔
591
                 "Missing <SourceFilename> or <VRTDataset> element in <%s>.",
592
                 psSrc->pszValue);
1✔
593
        return CE_Failure;
1✔
594
    }
595

596
    // Backup original filename and relativeToVRT so as to be able to
597
    // serialize them identically again (#5985)
598
    m_osSourceFileNameOri = pszFilename;
3,168✔
599
    if (pszFilename[0])
3,168✔
600
    {
601
        m_bRelativeToVRTOri =
3,165✔
602
            atoi(CPLGetXMLValue(psSourceFileNameNode, "relativetoVRT", "0"));
3,165✔
603
        const char *pszShared =
604
            CPLGetXMLValue(psSourceFileNameNode, "shared", nullptr);
3,165✔
605
        if (pszShared == nullptr)
3,165✔
606
        {
607
            pszShared = CPLGetConfigOption("VRT_SHARED_SOURCE", nullptr);
3,163✔
608
        }
609
        if (pszShared != nullptr)
3,165✔
610
        {
611
            m_nExplicitSharedStatus = CPLTestBool(pszShared);
8✔
612
        }
613

614
        m_osSrcDSName = VRTDataset::BuildSourceFilename(
6,330✔
615
            pszFilename, pszVRTPath, CPL_TO_BOOL(m_bRelativeToVRTOri));
6,330✔
616
    }
617
    else if (psSourceVRTDataset)
3✔
618
    {
619
        CPLXMLNode sNode;
620
        sNode.eType = psSourceVRTDataset->eType;
3✔
621
        sNode.pszValue = psSourceVRTDataset->pszValue;
3✔
622
        sNode.psNext = nullptr;
3✔
623
        sNode.psChild = psSourceVRTDataset->psChild;
3✔
624
        char *pszXML = CPLSerializeXMLTree(&sNode);
3✔
625
        if (pszXML)
3✔
626
        {
627
            m_bSrcDSNameFromVRT = true;
3✔
628
            m_osSrcDSName = pszXML;
3✔
629
            CPLFree(pszXML);
3✔
630
        }
631
    }
632

633
    const char *pszSourceBand = CPLGetXMLValue(psSrc, "SourceBand", "1");
3,168✔
634
    m_bGetMaskBand = false;
3,168✔
635
    if (STARTS_WITH_CI(pszSourceBand, "mask"))
3,168✔
636
    {
637
        m_bGetMaskBand = true;
22✔
638
        if (pszSourceBand[4] == ',')
22✔
639
            m_nBand = atoi(pszSourceBand + 5);
22✔
640
        else
641
            m_nBand = 1;
×
642
    }
643
    else
644
    {
645
        m_nBand = atoi(pszSourceBand);
3,146✔
646
    }
647
    if (!GDALCheckBandCount(m_nBand, 0))
3,168✔
648
    {
649
        CPLError(CE_Warning, CPLE_AppDefined,
×
650
                 "Invalid <SourceBand> element in VRTRasterBand.");
651
        return CE_Failure;
×
652
    }
653

654
    m_aosOpenOptions = GDALDeserializeOpenOptionsFromXML(psSrc);
3,168✔
655
    m_aosOpenOptionsOri = m_aosOpenOptions;
3,168✔
656
    if (strstr(m_osSrcDSName.c_str(), "<VRTDataset") != nullptr)
3,168✔
657
        m_aosOpenOptions.SetNameValue("ROOT_PATH", pszVRTPath);
5✔
658

659
    return ParseSrcRectAndDstRect(psSrc);
3,168✔
660
}
661

662
/************************************************************************/
663
/*                        ParseSrcRectAndDstRect()                      */
664
/************************************************************************/
665

666
CPLErr VRTSimpleSource::ParseSrcRectAndDstRect(const CPLXMLNode *psSrc)
3,182✔
667
{
668
    const auto GetAttrValue = [](const CPLXMLNode *psNode,
23,740✔
669
                                 const char *pszAttrName, double dfDefaultVal)
670
    {
671
        if (const char *pszVal = CPLGetXMLValue(psNode, pszAttrName, nullptr))
23,740✔
672
            return CPLAtof(pszVal);
23,740✔
673
        else
674
            return dfDefaultVal;
×
675
    };
676

677
    /* -------------------------------------------------------------------- */
678
    /*      Set characteristics.                                            */
679
    /* -------------------------------------------------------------------- */
680
    const CPLXMLNode *const psSrcRect = CPLGetXMLNode(psSrc, "SrcRect");
3,182✔
681
    if (psSrcRect)
3,182✔
682
    {
683
        double xOff = GetAttrValue(psSrcRect, "xOff", UNINIT_WINDOW);
2,968✔
684
        double yOff = GetAttrValue(psSrcRect, "yOff", UNINIT_WINDOW);
2,968✔
685
        double xSize = GetAttrValue(psSrcRect, "xSize", UNINIT_WINDOW);
2,968✔
686
        double ySize = GetAttrValue(psSrcRect, "ySize", UNINIT_WINDOW);
2,968✔
687
        // Test written that way to catch NaN values
688
        if (!(xOff >= INT_MIN && xOff <= INT_MAX) ||
2,968✔
689
            !(yOff >= INT_MIN && yOff <= INT_MAX) ||
2,968✔
690
            !(xSize > 0 || xSize == UNINIT_WINDOW) || xSize > INT_MAX ||
2,968✔
691
            !(ySize > 0 || ySize == UNINIT_WINDOW) || ySize > INT_MAX)
2,967✔
692
        {
693
            CPLError(CE_Failure, CPLE_AppDefined, "Wrong values in SrcRect");
1✔
694
            return CE_Failure;
1✔
695
        }
696
        SetSrcWindow(xOff, yOff, xSize, ySize);
2,967✔
697
    }
698
    else
699
    {
700
        m_dfSrcXOff = UNINIT_WINDOW;
214✔
701
        m_dfSrcYOff = UNINIT_WINDOW;
214✔
702
        m_dfSrcXSize = UNINIT_WINDOW;
214✔
703
        m_dfSrcYSize = UNINIT_WINDOW;
214✔
704
    }
705

706
    const CPLXMLNode *const psDstRect = CPLGetXMLNode(psSrc, "DstRect");
3,181✔
707
    if (psDstRect)
3,181✔
708
    {
709
        double xOff = GetAttrValue(psDstRect, "xOff", UNINIT_WINDOW);
2,967✔
710
        ;
711
        double yOff = GetAttrValue(psDstRect, "yOff", UNINIT_WINDOW);
2,967✔
712
        double xSize = GetAttrValue(psDstRect, "xSize", UNINIT_WINDOW);
2,967✔
713
        ;
714
        double ySize = GetAttrValue(psDstRect, "ySize", UNINIT_WINDOW);
2,967✔
715
        // Test written that way to catch NaN values
716
        if (!(xOff >= INT_MIN && xOff <= INT_MAX) ||
2,967✔
717
            !(yOff >= INT_MIN && yOff <= INT_MAX) ||
2,967✔
718
            !(xSize > 0 || xSize == UNINIT_WINDOW) || xSize > INT_MAX ||
2,967✔
719
            !(ySize > 0 || ySize == UNINIT_WINDOW) || ySize > INT_MAX)
2,967✔
720
        {
721
            CPLError(CE_Failure, CPLE_AppDefined, "Wrong values in DstRect");
1✔
722
            return CE_Failure;
1✔
723
        }
724
        SetDstWindow(xOff, yOff, xSize, ySize);
2,966✔
725
    }
726
    else
727
    {
728
        m_dfDstXOff = UNINIT_WINDOW;
214✔
729
        m_dfDstYOff = UNINIT_WINDOW;
214✔
730
        m_dfDstXSize = UNINIT_WINDOW;
214✔
731
        m_dfDstYSize = UNINIT_WINDOW;
214✔
732
    }
733

734
    return CE_None;
3,180✔
735
}
736

737
/************************************************************************/
738
/*                             GetFileList()                            */
739
/************************************************************************/
740

741
void VRTSimpleSource::GetFileList(char ***ppapszFileList, int *pnSize,
112✔
742
                                  int *pnMaxSize, CPLHashSet *hSetFiles)
743
{
744
    if (!m_osSrcDSName.empty() && !m_bSrcDSNameFromVRT)
112✔
745
    {
746
        const char *pszFilename = m_osSrcDSName.c_str();
112✔
747

748
        /* --------------------------------------------------------------------
749
         */
750
        /*      Is it already in the list ? */
751
        /* --------------------------------------------------------------------
752
         */
753
        if (CPLHashSetLookup(hSetFiles, pszFilename) != nullptr)
112✔
754
            return;
15✔
755

756
        /* --------------------------------------------------------------------
757
         */
758
        /*      Grow array if necessary */
759
        /* --------------------------------------------------------------------
760
         */
761
        if (*pnSize + 1 >= *pnMaxSize)
97✔
762
        {
763
            *pnMaxSize = std::max(*pnSize + 2, 2 + 2 * (*pnMaxSize));
66✔
764
            *ppapszFileList = static_cast<char **>(
66✔
765
                CPLRealloc(*ppapszFileList, sizeof(char *) * (*pnMaxSize)));
66✔
766
        }
767

768
        /* --------------------------------------------------------------------
769
         */
770
        /*      Add the string to the list */
771
        /* --------------------------------------------------------------------
772
         */
773
        (*ppapszFileList)[*pnSize] = CPLStrdup(pszFilename);
97✔
774
        (*ppapszFileList)[(*pnSize + 1)] = nullptr;
97✔
775
        CPLHashSetInsert(hSetFiles, (*ppapszFileList)[*pnSize]);
97✔
776

777
        (*pnSize)++;
97✔
778
    }
779
}
780

781
/************************************************************************/
782
/*                           OpenSource()                               */
783
/************************************************************************/
784

785
void VRTSimpleSource::OpenSource() const
1,263✔
786
{
787
    CPLAssert(m_poRasterBand == nullptr);
1,263✔
788

789
    /* ----------------------------------------------------------------- */
790
    /*      Create a proxy dataset                                       */
791
    /* ----------------------------------------------------------------- */
792
    GDALProxyPoolDataset *proxyDS = nullptr;
1,263✔
793
    std::string osKeyMapSharedSources;
1,263✔
794
    if (m_poMapSharedSources)
1,263✔
795
    {
796
        osKeyMapSharedSources = m_osSrcDSName;
1,215✔
797
        for (int i = 0; i < m_aosOpenOptions.size(); ++i)
1,224✔
798
        {
799
            osKeyMapSharedSources += "||";
9✔
800
            osKeyMapSharedSources += m_aosOpenOptions[i];
9✔
801
        }
802

803
        proxyDS = cpl::down_cast<GDALProxyPoolDataset *>(
1,215✔
804
            m_poMapSharedSources->Get(osKeyMapSharedSources));
1,215✔
805
    }
806

807
    if (proxyDS == nullptr)
1,263✔
808
    {
809
        int bShared = true;
1,064✔
810
        if (m_nExplicitSharedStatus != -1)
1,064✔
811
            bShared = m_nExplicitSharedStatus;
8✔
812

813
        const CPLString osUniqueHandle(CPLSPrintf("%p", m_poMapSharedSources));
1,064✔
814
        proxyDS = GDALProxyPoolDataset::Create(
1,064✔
815
            m_osSrcDSName, m_aosOpenOptions.List(), GA_ReadOnly, bShared,
816
            osUniqueHandle.c_str());
817
        if (proxyDS == nullptr)
1,064✔
818
            return;
172✔
819
    }
820
    else
821
    {
822
        proxyDS->Reference();
199✔
823
    }
824

825
    if (m_bGetMaskBand)
1,091✔
826
    {
827
        GDALProxyPoolRasterBand *poMaskBand =
828
            cpl::down_cast<GDALProxyPoolRasterBand *>(
15✔
829
                proxyDS->GetRasterBand(m_nBand));
15✔
830
        poMaskBand->AddSrcMaskBandDescriptionFromUnderlying();
15✔
831
    }
832

833
    /* -------------------------------------------------------------------- */
834
    /*      Get the raster band.                                            */
835
    /* -------------------------------------------------------------------- */
836

837
    m_poRasterBand = proxyDS->GetRasterBand(m_nBand);
1,091✔
838
    if (m_poRasterBand == nullptr || !ValidateOpenedBand(m_poRasterBand))
1,091✔
839
    {
840
        proxyDS->ReleaseRef();
2✔
841
        return;
2✔
842
    }
843

844
    if (m_bGetMaskBand)
1,089✔
845
    {
846
        m_poRasterBand = m_poRasterBand->GetMaskBand();
15✔
847
        if (m_poRasterBand == nullptr)
15✔
848
        {
849
            proxyDS->ReleaseRef();
×
850
            return;
×
851
        }
852
        m_poMaskBandMainBand = m_poRasterBand;
15✔
853
    }
854

855
    if (m_poMapSharedSources)
1,089✔
856
    {
857
        m_poMapSharedSources->Insert(osKeyMapSharedSources, proxyDS);
1,041✔
858
    }
859
}
860

861
/************************************************************************/
862
/*                         GetRasterBand()                              */
863
/************************************************************************/
864

865
GDALRasterBand *VRTSimpleSource::GetRasterBand() const
105,982✔
866
{
867
    if (m_poRasterBand == nullptr)
105,982✔
868
        OpenSource();
1,261✔
869
    return m_poRasterBand;
105,982✔
870
}
871

872
/************************************************************************/
873
/*                        GetMaskBandMainBand()                         */
874
/************************************************************************/
875

876
GDALRasterBand *VRTSimpleSource::GetMaskBandMainBand()
1,301✔
877
{
878
    if (m_poRasterBand == nullptr)
1,301✔
879
        OpenSource();
2✔
880
    return m_poMaskBandMainBand;
1,301✔
881
}
882

883
/************************************************************************/
884
/*                       IsSameExceptBandNumber()                       */
885
/************************************************************************/
886

887
bool VRTSimpleSource::IsSameExceptBandNumber(
979✔
888
    const VRTSimpleSource *poOtherSource) const
889
{
890
    return m_dfSrcXOff == poOtherSource->m_dfSrcXOff &&
1,958✔
891
           m_dfSrcYOff == poOtherSource->m_dfSrcYOff &&
979✔
892
           m_dfSrcXSize == poOtherSource->m_dfSrcXSize &&
979✔
893
           m_dfSrcYSize == poOtherSource->m_dfSrcYSize &&
979✔
894
           m_dfDstXOff == poOtherSource->m_dfDstXOff &&
979✔
895
           m_dfDstYOff == poOtherSource->m_dfDstYOff &&
979✔
896
           m_dfDstXSize == poOtherSource->m_dfDstXSize &&
979✔
897
           m_dfDstYSize == poOtherSource->m_dfDstYSize &&
979✔
898
           !m_osSrcDSName.empty() &&
2,937✔
899
           m_osSrcDSName == poOtherSource->m_osSrcDSName;
1,958✔
900
}
901

902
/************************************************************************/
903
/*                              SrcToDst()                              */
904
/*                                                                      */
905
/*      Note: this is a no-op if the both src and dst windows are unset */
906
/************************************************************************/
907

908
void VRTSimpleSource::SrcToDst(double dfX, double dfY, double &dfXOut,
6,532✔
909
                               double &dfYOut) const
910

911
{
912
    dfXOut = ((dfX - m_dfSrcXOff) / m_dfSrcXSize) * m_dfDstXSize + m_dfDstXOff;
6,532✔
913
    dfYOut = ((dfY - m_dfSrcYOff) / m_dfSrcYSize) * m_dfDstYSize + m_dfDstYOff;
6,532✔
914
}
6,532✔
915

916
/************************************************************************/
917
/*                              DstToSrc()                              */
918
/*                                                                      */
919
/*      Note: this is a no-op if the both src and dst windows are unset */
920
/************************************************************************/
921

922
void VRTSimpleSource::DstToSrc(double dfX, double dfY, double &dfXOut,
4,022,310✔
923
                               double &dfYOut) const
924

925
{
926
    dfXOut = ((dfX - m_dfDstXOff) / m_dfDstXSize) * m_dfSrcXSize + m_dfSrcXOff;
4,022,310✔
927
    dfYOut = ((dfY - m_dfDstYOff) / m_dfDstYSize) * m_dfSrcYSize + m_dfSrcYOff;
4,022,310✔
928
}
4,022,310✔
929

930
/************************************************************************/
931
/*                          GetSrcDstWindow()                           */
932
/************************************************************************/
933

934
int VRTSimpleSource::GetSrcDstWindow(
70,264✔
935
    double dfXOff, double dfYOff, double dfXSize, double dfYSize, int nBufXSize,
936
    int nBufYSize, double *pdfReqXOff, double *pdfReqYOff, double *pdfReqXSize,
937
    double *pdfReqYSize, int *pnReqXOff, int *pnReqYOff, int *pnReqXSize,
938
    int *pnReqYSize, int *pnOutXOff, int *pnOutYOff, int *pnOutXSize,
939
    int *pnOutYSize, bool &bErrorOut)
940

941
{
942
    bErrorOut = false;
70,264✔
943

944
    if (m_dfSrcXSize == 0.0 || m_dfSrcYSize == 0.0 || m_dfDstXSize == 0.0 ||
70,264✔
945
        m_dfDstYSize == 0.0)
70,265✔
946
    {
947
        return FALSE;
×
948
    }
949

950
    const bool bDstWinSet = IsDstWinSet();
70,265✔
951

952
#ifdef DEBUG
953
    const bool bSrcWinSet = IsSrcWinSet();
70,265✔
954

955
    if (bSrcWinSet != bDstWinSet)
70,265✔
956
    {
957
        return FALSE;
×
958
    }
959
#endif
960

961
    /* -------------------------------------------------------------------- */
962
    /*      If the input window completely misses the portion of the        */
963
    /*      virtual dataset provided by this source we have nothing to do.  */
964
    /* -------------------------------------------------------------------- */
965
    if (bDstWinSet)
70,265✔
966
    {
967
        if (dfXOff >= m_dfDstXOff + m_dfDstXSize ||
69,991✔
968
            dfYOff >= m_dfDstYOff + m_dfDstYSize ||
59,563✔
969
            dfXOff + dfXSize <= m_dfDstXOff || dfYOff + dfYSize <= m_dfDstYOff)
59,027✔
970
            return FALSE;
20,311✔
971
    }
972

973
    /* -------------------------------------------------------------------- */
974
    /*      This request window corresponds to the whole output buffer.     */
975
    /* -------------------------------------------------------------------- */
976
    *pnOutXOff = 0;
49,954✔
977
    *pnOutYOff = 0;
49,954✔
978
    *pnOutXSize = nBufXSize;
49,954✔
979
    *pnOutYSize = nBufYSize;
49,954✔
980

981
    /* -------------------------------------------------------------------- */
982
    /*      If the input window extents outside the portion of the on       */
983
    /*      the virtual file that this source can set, then clip down       */
984
    /*      the requested window.                                           */
985
    /* -------------------------------------------------------------------- */
986
    bool bModifiedX = false;
49,954✔
987
    bool bModifiedY = false;
49,954✔
988
    double dfRXOff = dfXOff;
49,954✔
989
    double dfRYOff = dfYOff;
49,954✔
990
    double dfRXSize = dfXSize;
49,954✔
991
    double dfRYSize = dfYSize;
49,954✔
992

993
    if (bDstWinSet)
49,954✔
994
    {
995
        if (dfRXOff < m_dfDstXOff)
49,679✔
996
        {
997
            dfRXSize = dfRXSize + dfRXOff - m_dfDstXOff;
1,862✔
998
            dfRXOff = m_dfDstXOff;
1,862✔
999
            bModifiedX = true;
1,862✔
1000
        }
1001

1002
        if (dfRYOff < m_dfDstYOff)
49,679✔
1003
        {
1004
            dfRYSize = dfRYSize + dfRYOff - m_dfDstYOff;
166✔
1005
            dfRYOff = m_dfDstYOff;
166✔
1006
            bModifiedY = true;
166✔
1007
        }
1008

1009
        if (dfRXOff + dfRXSize > m_dfDstXOff + m_dfDstXSize)
49,679✔
1010
        {
1011
            dfRXSize = m_dfDstXOff + m_dfDstXSize - dfRXOff;
1,938✔
1012
            bModifiedX = true;
1,938✔
1013
        }
1014

1015
        if (dfRYOff + dfRYSize > m_dfDstYOff + m_dfDstYSize)
49,679✔
1016
        {
1017
            dfRYSize = m_dfDstYOff + m_dfDstYSize - dfRYOff;
267✔
1018
            bModifiedY = true;
267✔
1019
        }
1020
    }
1021

1022
    /* -------------------------------------------------------------------- */
1023
    /*      Translate requested region in virtual file into the source      */
1024
    /*      band coordinates.                                               */
1025
    /* -------------------------------------------------------------------- */
1026
    const double dfScaleX = m_dfSrcXSize / m_dfDstXSize;
49,954✔
1027
    const double dfScaleY = m_dfSrcYSize / m_dfDstYSize;
49,954✔
1028

1029
    *pdfReqXOff = (dfRXOff - m_dfDstXOff) * dfScaleX + m_dfSrcXOff;
49,954✔
1030
    *pdfReqYOff = (dfRYOff - m_dfDstYOff) * dfScaleY + m_dfSrcYOff;
49,954✔
1031
    *pdfReqXSize = dfRXSize * dfScaleX;
49,954✔
1032
    *pdfReqYSize = dfRYSize * dfScaleY;
49,954✔
1033

1034
    if (!std::isfinite(*pdfReqXOff) || !std::isfinite(*pdfReqYOff) ||
99,908✔
1035
        !std::isfinite(*pdfReqXSize) || !std::isfinite(*pdfReqYSize) ||
49,954✔
1036
        *pdfReqXOff > INT_MAX || *pdfReqYOff > INT_MAX || *pdfReqXSize < 0 ||
149,862✔
1037
        *pdfReqYSize < 0)
49,954✔
1038
    {
1039
        return FALSE;
×
1040
    }
1041

1042
    /* -------------------------------------------------------------------- */
1043
    /*      Clamp within the bounds of the available source data.           */
1044
    /* -------------------------------------------------------------------- */
1045
    if (*pdfReqXOff < 0)
49,954✔
1046
    {
1047
        *pdfReqXSize += *pdfReqXOff;
6✔
1048
        *pdfReqXOff = 0;
6✔
1049
        bModifiedX = true;
6✔
1050
    }
1051
    if (*pdfReqYOff < 0)
49,954✔
1052
    {
1053
        *pdfReqYSize += *pdfReqYOff;
6✔
1054
        *pdfReqYOff = 0;
6✔
1055
        bModifiedY = true;
6✔
1056
    }
1057

1058
    *pnReqXOff = static_cast<int>(floor(*pdfReqXOff));
49,954✔
1059
    *pnReqYOff = static_cast<int>(floor(*pdfReqYOff));
49,954✔
1060

1061
    constexpr double EPS = 1e-3;
49,954✔
1062
    constexpr double ONE_MINUS_EPS = 1.0 - EPS;
49,954✔
1063
    if (*pdfReqXOff - *pnReqXOff > ONE_MINUS_EPS)
49,954✔
1064
    {
1065
        (*pnReqXOff)++;
2✔
1066
        *pdfReqXOff = *pnReqXOff;
2✔
1067
    }
1068
    if (*pdfReqYOff - *pnReqYOff > ONE_MINUS_EPS)
49,954✔
1069
    {
1070
        (*pnReqYOff)++;
18✔
1071
        *pdfReqYOff = *pnReqYOff;
18✔
1072
    }
1073

1074
    if (*pdfReqXSize > INT_MAX)
49,954✔
1075
        *pnReqXSize = INT_MAX;
×
1076
    else
1077
        *pnReqXSize = static_cast<int>(floor(*pdfReqXSize + 0.5));
49,954✔
1078

1079
    if (*pdfReqYSize > INT_MAX)
49,954✔
1080
        *pnReqYSize = INT_MAX;
×
1081
    else
1082
        *pnReqYSize = static_cast<int>(floor(*pdfReqYSize + 0.5));
49,954✔
1083

1084
    /* -------------------------------------------------------------------- */
1085
    /*      Clamp within the bounds of the available source data.           */
1086
    /* -------------------------------------------------------------------- */
1087

1088
    if (*pnReqXSize == 0)
49,954✔
1089
        *pnReqXSize = 1;
675✔
1090
    if (*pnReqYSize == 0)
49,954✔
1091
        *pnReqYSize = 1;
22,678✔
1092

1093
    auto l_band = GetRasterBand();
49,954✔
1094
    if (!l_band)
49,954✔
1095
    {
1096
        bErrorOut = true;
86✔
1097
        return FALSE;
86✔
1098
    }
1099
    if (*pnReqXSize > INT_MAX - *pnReqXOff ||
99,736✔
1100
        *pnReqXOff + *pnReqXSize > l_band->GetXSize())
49,868✔
1101
    {
1102
        *pnReqXSize = l_band->GetXSize() - *pnReqXOff;
44✔
1103
        bModifiedX = true;
44✔
1104
    }
1105
    if (*pdfReqXOff + *pdfReqXSize > l_band->GetXSize())
49,868✔
1106
    {
1107
        *pdfReqXSize = l_band->GetXSize() - *pdfReqXOff;
44✔
1108
        bModifiedX = true;
44✔
1109
    }
1110

1111
    if (*pnReqYSize > INT_MAX - *pnReqYOff ||
99,736✔
1112
        *pnReqYOff + *pnReqYSize > l_band->GetYSize())
49,868✔
1113
    {
1114
        *pnReqYSize = l_band->GetYSize() - *pnReqYOff;
44✔
1115
        bModifiedY = true;
44✔
1116
    }
1117
    if (*pdfReqYOff + *pdfReqYSize > l_band->GetYSize())
49,868✔
1118
    {
1119
        *pdfReqYSize = l_band->GetYSize() - *pdfReqYOff;
62✔
1120
        bModifiedY = true;
62✔
1121
    }
1122

1123
    /* -------------------------------------------------------------------- */
1124
    /*      Don't do anything if the requesting region is completely off    */
1125
    /*      the source image.                                               */
1126
    /* -------------------------------------------------------------------- */
1127
    if (*pnReqXOff >= l_band->GetXSize() || *pnReqYOff >= l_band->GetYSize() ||
99,732✔
1128
        *pnReqXSize <= 0 || *pnReqYSize <= 0)
99,732✔
1129
    {
1130
        return FALSE;
9✔
1131
    }
1132

1133
    /* -------------------------------------------------------------------- */
1134
    /*      If we haven't had to modify the source rectangle, then the      */
1135
    /*      destination rectangle must be the whole region.                 */
1136
    /* -------------------------------------------------------------------- */
1137
    if (bModifiedX || bModifiedY)
49,859✔
1138
    {
1139
        /* --------------------------------------------------------------------
1140
         */
1141
        /*      Now transform this possibly reduced request back into the */
1142
        /*      destination buffer coordinates in case the output region is */
1143
        /*      less than the whole buffer. */
1144
        /* --------------------------------------------------------------------
1145
         */
1146
        double dfDstULX = 0.0;
3,266✔
1147
        double dfDstULY = 0.0;
3,266✔
1148
        double dfDstLRX = 0.0;
3,266✔
1149
        double dfDstLRY = 0.0;
3,266✔
1150

1151
        SrcToDst(*pdfReqXOff, *pdfReqYOff, dfDstULX, dfDstULY);
3,266✔
1152
        SrcToDst(*pdfReqXOff + *pdfReqXSize, *pdfReqYOff + *pdfReqYSize,
3,266✔
1153
                 dfDstLRX, dfDstLRY);
1154
#if DEBUG_VERBOSE
1155
        CPLDebug("VRT", "dfDstULX=%g dfDstULY=%g dfDstLRX=%g dfDstLRY=%g",
1156
                 dfDstULX, dfDstULY, dfDstLRX, dfDstLRY);
1157
#endif
1158

1159
        if (bModifiedX)
3,266✔
1160
        {
1161
            const double dfScaleWinToBufX = nBufXSize / dfXSize;
3,216✔
1162

1163
            const double dfOutXOff = (dfDstULX - dfXOff) * dfScaleWinToBufX;
3,216✔
1164
            if (dfOutXOff <= 0)
3,216✔
1165
                *pnOutXOff = 0;
1,356✔
1166
            else if (dfOutXOff > INT_MAX)
1,860✔
1167
                *pnOutXOff = INT_MAX;
×
1168
            else
1169
                *pnOutXOff = static_cast<int>(dfOutXOff + EPS);
1,860✔
1170

1171
            // Apply correction on floating-point source window
1172
            {
1173
                double dfDstDeltaX =
3,216✔
1174
                    (dfOutXOff - *pnOutXOff) / dfScaleWinToBufX;
3,216✔
1175
                double dfSrcDeltaX = dfDstDeltaX / m_dfDstXSize * m_dfSrcXSize;
3,216✔
1176
                *pdfReqXOff -= dfSrcDeltaX;
3,216✔
1177
                *pdfReqXSize = std::min(*pdfReqXSize + dfSrcDeltaX,
6,432✔
1178
                                        static_cast<double>(INT_MAX));
3,216✔
1179
            }
1180

1181
            double dfOutRightXOff = (dfDstLRX - dfXOff) * dfScaleWinToBufX;
3,216✔
1182
            if (dfOutRightXOff < dfOutXOff)
3,216✔
1183
                return FALSE;
×
1184
            if (dfOutRightXOff > INT_MAX)
3,216✔
1185
                dfOutRightXOff = INT_MAX;
×
1186
            const int nOutRightXOff =
3,216✔
1187
                static_cast<int>(ceil(dfOutRightXOff - EPS));
3,216✔
1188
            *pnOutXSize = nOutRightXOff - *pnOutXOff;
3,216✔
1189

1190
            if (*pnOutXSize > INT_MAX - *pnOutXOff ||
3,216✔
1191
                *pnOutXOff + *pnOutXSize > nBufXSize)
3,216✔
1192
                *pnOutXSize = nBufXSize - *pnOutXOff;
×
1193

1194
            // Apply correction on floating-point source window
1195
            {
1196
                double dfDstDeltaX =
3,216✔
1197
                    (nOutRightXOff - dfOutRightXOff) / dfScaleWinToBufX;
3,216✔
1198
                double dfSrcDeltaX = dfDstDeltaX / m_dfDstXSize * m_dfSrcXSize;
3,216✔
1199
                *pdfReqXSize = std::min(*pdfReqXSize + dfSrcDeltaX,
6,432✔
1200
                                        static_cast<double>(INT_MAX));
3,216✔
1201
            }
1202
        }
1203

1204
        if (bModifiedY)
3,266✔
1205
        {
1206
            const double dfScaleWinToBufY = nBufYSize / dfYSize;
377✔
1207

1208
            const double dfOutYOff = (dfDstULY - dfYOff) * dfScaleWinToBufY;
377✔
1209
            if (dfOutYOff <= 0)
377✔
1210
                *pnOutYOff = 0;
194✔
1211
            else if (dfOutYOff > INT_MAX)
183✔
1212
                *pnOutYOff = INT_MAX;
×
1213
            else
1214
                *pnOutYOff = static_cast<int>(dfOutYOff + EPS);
183✔
1215

1216
            // Apply correction on floating-point source window
1217
            {
1218
                double dfDstDeltaY =
377✔
1219
                    (dfOutYOff - *pnOutYOff) / dfScaleWinToBufY;
377✔
1220
                double dfSrcDeltaY = dfDstDeltaY / m_dfDstYSize * m_dfSrcYSize;
377✔
1221
                *pdfReqYOff -= dfSrcDeltaY;
377✔
1222
                *pdfReqYSize = std::min(*pdfReqYSize + dfSrcDeltaY,
754✔
1223
                                        static_cast<double>(INT_MAX));
377✔
1224
            }
1225

1226
            double dfOutTopYOff = (dfDstLRY - dfYOff) * dfScaleWinToBufY;
377✔
1227
            if (dfOutTopYOff < dfOutYOff)
377✔
1228
                return FALSE;
×
1229
            if (dfOutTopYOff > INT_MAX)
377✔
1230
                dfOutTopYOff = INT_MAX;
×
1231
            const int nOutTopYOff = static_cast<int>(ceil(dfOutTopYOff - EPS));
377✔
1232
            *pnOutYSize = nOutTopYOff - *pnOutYOff;
377✔
1233

1234
            if (*pnOutYSize > INT_MAX - *pnOutYOff ||
377✔
1235
                *pnOutYOff + *pnOutYSize > nBufYSize)
377✔
1236
                *pnOutYSize = nBufYSize - *pnOutYOff;
×
1237

1238
            // Apply correction on floating-point source window
1239
            {
1240
                double dfDstDeltaY =
377✔
1241
                    (nOutTopYOff - dfOutTopYOff) / dfScaleWinToBufY;
377✔
1242
                double dfSrcDeltaY = dfDstDeltaY / m_dfDstYSize * m_dfSrcYSize;
377✔
1243
                *pdfReqYSize = std::min(*pdfReqYSize + dfSrcDeltaY,
754✔
1244
                                        static_cast<double>(INT_MAX));
377✔
1245
            }
1246
        }
1247

1248
        if (*pnOutXSize < 1 || *pnOutYSize < 1)
3,266✔
1249
            return FALSE;
×
1250
    }
1251

1252
    *pdfReqXOff = RoundIfCloseToInt(*pdfReqXOff);
49,859✔
1253
    *pdfReqYOff = RoundIfCloseToInt(*pdfReqYOff);
49,859✔
1254
    *pdfReqXSize = RoundIfCloseToInt(*pdfReqXSize);
49,859✔
1255
    *pdfReqYSize = RoundIfCloseToInt(*pdfReqYSize);
49,859✔
1256

1257
    return TRUE;
49,859✔
1258
}
1259

1260
/************************************************************************/
1261
/*                          NeedMaxValAdjustment()                      */
1262
/************************************************************************/
1263

1264
int VRTSimpleSource::NeedMaxValAdjustment() const
42,469✔
1265
{
1266
    if (!m_nMaxValue)
42,469✔
1267
        return FALSE;
42,447✔
1268

1269
    auto l_band = GetRasterBand();
22✔
1270
    if (!l_band)
22✔
1271
        return FALSE;
×
1272
    const char *pszNBITS = l_band->GetMetadataItem("NBITS", "IMAGE_STRUCTURE");
22✔
1273
    const int nBits = (pszNBITS) ? atoi(pszNBITS) : 0;
22✔
1274
    if (nBits >= 1 && nBits <= 31)
22✔
1275
    {
1276
        const int nBandMaxValue = static_cast<int>((1U << nBits) - 1);
×
1277
        return nBandMaxValue > m_nMaxValue;
×
1278
    }
1279
    return TRUE;
22✔
1280
}
1281

1282
/************************************************************************/
1283
/*                              RasterIO()                              */
1284
/************************************************************************/
1285

1286
CPLErr VRTSimpleSource::RasterIO(GDALDataType eVRTBandDataType, int nXOff,
50,687✔
1287
                                 int nYOff, int nXSize, int nYSize, void *pData,
1288
                                 int nBufXSize, int nBufYSize,
1289
                                 GDALDataType eBufType, GSpacing nPixelSpace,
1290
                                 GSpacing nLineSpace,
1291
                                 GDALRasterIOExtraArg *psExtraArgIn,
1292
                                 WorkingState & /*oWorkingState*/)
1293

1294
{
1295
    GDALRasterIOExtraArg sExtraArg;
1296
    INIT_RASTERIO_EXTRA_ARG(sExtraArg);
50,687✔
1297
    GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
50,687✔
1298

1299
    double dfXOff = nXOff;
50,687✔
1300
    double dfYOff = nYOff;
50,687✔
1301
    double dfXSize = nXSize;
50,687✔
1302
    double dfYSize = nYSize;
50,687✔
1303
    if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
50,687✔
1304
    {
1305
        dfXOff = psExtraArgIn->dfXOff;
62✔
1306
        dfYOff = psExtraArgIn->dfYOff;
62✔
1307
        dfXSize = psExtraArgIn->dfXSize;
62✔
1308
        dfYSize = psExtraArgIn->dfYSize;
62✔
1309
    }
1310

1311
    // The window we will actually request from the source raster band.
1312
    double dfReqXOff = 0.0;
50,687✔
1313
    double dfReqYOff = 0.0;
50,687✔
1314
    double dfReqXSize = 0.0;
50,687✔
1315
    double dfReqYSize = 0.0;
50,687✔
1316
    int nReqXOff = 0;
50,687✔
1317
    int nReqYOff = 0;
50,687✔
1318
    int nReqXSize = 0;
50,687✔
1319
    int nReqYSize = 0;
50,687✔
1320

1321
    // The window we will actual set _within_ the pData buffer.
1322
    int nOutXOff = 0;
50,687✔
1323
    int nOutYOff = 0;
50,687✔
1324
    int nOutXSize = 0;
50,687✔
1325
    int nOutYSize = 0;
50,687✔
1326

1327
    bool bError = false;
50,687✔
1328
    if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
50,687✔
1329
                         &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
1330
                         &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
1331
                         &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
1332
    {
1333
        return bError ? CE_Failure : CE_None;
9,314✔
1334
    }
1335
#if DEBUG_VERBOSE
1336
    CPLDebug("VRT",
1337
             "nXOff=%d, nYOff=%d, nXSize=%d, nYSize=%d, nBufXSize=%d, "
1338
             "nBufYSize=%d,\n"
1339
             "dfReqXOff=%g, dfReqYOff=%g, dfReqXSize=%g, dfReqYSize=%g,\n"
1340
             "nReqXOff=%d, nReqYOff=%d, nReqXSize=%d, nReqYSize=%d,\n"
1341
             "nOutXOff=%d, nOutYOff=%d, nOutXSize=%d, nOutYSize=%d",
1342
             nXOff, nYOff, nXSize, nYSize, nBufXSize, nBufYSize, dfReqXOff,
1343
             dfReqYOff, dfReqXSize, dfReqYSize, nReqXOff, nReqYOff, nReqXSize,
1344
             nReqYSize, nOutXOff, nOutYOff, nOutXSize, nOutYSize);
1345
#endif
1346

1347
    /* -------------------------------------------------------------------- */
1348
    /*      Actually perform the IO request.                                */
1349
    /* -------------------------------------------------------------------- */
1350
    if (!m_osResampling.empty())
41,373✔
1351
    {
1352
        psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
2,052✔
1353
    }
1354
    else if (psExtraArgIn != nullptr)
39,321✔
1355
    {
1356
        psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
39,321✔
1357
    }
1358
    psExtraArg->bFloatingPointWindowValidity = TRUE;
41,373✔
1359
    psExtraArg->dfXOff = dfReqXOff;
41,373✔
1360
    psExtraArg->dfYOff = dfReqYOff;
41,373✔
1361
    psExtraArg->dfXSize = dfReqXSize;
41,373✔
1362
    psExtraArg->dfYSize = dfReqYSize;
41,373✔
1363
    if (psExtraArgIn)
41,373✔
1364
    {
1365
        psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
41,373✔
1366
        psExtraArg->pProgressData = psExtraArgIn->pProgressData;
41,373✔
1367
    }
1368

1369
    GByte *pabyOut = static_cast<unsigned char *>(pData) +
41,373✔
1370
                     nOutXOff * nPixelSpace +
41,373✔
1371
                     static_cast<GPtrDiff_t>(nOutYOff) * nLineSpace;
41,373✔
1372

1373
    auto l_band = GetRasterBand();
41,373✔
1374
    if (!l_band)
41,373✔
1375
        return CE_Failure;
×
1376

1377
    CPLErr eErr = CE_Failure;
41,373✔
1378
    if (GDALDataTypeIsConversionLossy(l_band->GetRasterDataType(),
41,373✔
1379
                                      eVRTBandDataType))
41,373✔
1380
    {
1381
        const int nBandDTSize = GDALGetDataTypeSizeBytes(eVRTBandDataType);
49✔
1382
        void *pTemp = VSI_MALLOC3_VERBOSE(nOutXSize, nOutYSize, nBandDTSize);
49✔
1383
        if (pTemp)
49✔
1384
        {
1385
            eErr = l_band->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize,
49✔
1386
                                    nReqYSize, pTemp, nOutXSize, nOutYSize,
1387
                                    eVRTBandDataType, 0, 0, psExtraArg);
1388
            if (eErr == CE_None)
49✔
1389
            {
1390
                GByte *pabyTemp = static_cast<GByte *>(pTemp);
49✔
1391
                for (int iY = 0; iY < nOutYSize; iY++)
779✔
1392
                {
1393
                    GDALCopyWords(
730✔
1394
                        pabyTemp +
730✔
1395
                            static_cast<size_t>(iY) * nBandDTSize * nOutXSize,
730✔
1396
                        eVRTBandDataType, nBandDTSize,
1397
                        pabyOut + static_cast<GPtrDiff_t>(iY * nLineSpace),
730✔
1398
                        eBufType, static_cast<int>(nPixelSpace), nOutXSize);
1399
                }
1400
            }
1401
            VSIFree(pTemp);
49✔
1402
        }
1403
    }
1404
    else
1405
    {
1406
        eErr = l_band->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize,
41,324✔
1407
                                nReqYSize, pabyOut, nOutXSize, nOutYSize,
1408
                                eBufType, nPixelSpace, nLineSpace, psExtraArg);
1409
    }
1410

1411
    if (NeedMaxValAdjustment())
41,373✔
1412
    {
1413
        for (int j = 0; j < nOutYSize; j++)
234✔
1414
        {
1415
            for (int i = 0; i < nOutXSize; i++)
4,620✔
1416
            {
1417
                int nVal = 0;
4,400✔
1418
                GDALCopyWords(pabyOut + j * nLineSpace + i * nPixelSpace,
4,400✔
1419
                              eBufType, 0, &nVal, GDT_Int32, 0, 1);
1420
                if (nVal > m_nMaxValue)
4,400✔
1421
                    nVal = m_nMaxValue;
800✔
1422
                GDALCopyWords(&nVal, GDT_Int32, 0,
4,400✔
1423
                              pabyOut + j * nLineSpace + i * nPixelSpace,
4,400✔
1424
                              eBufType, 0, 1);
1425
            }
1426
        }
1427
    }
1428

1429
    if (psExtraArg->pfnProgress)
41,373✔
1430
        psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
2,991✔
1431

1432
    return eErr;
41,373✔
1433
}
1434

1435
/************************************************************************/
1436
/*                             GetMinimum()                             */
1437
/************************************************************************/
1438

1439
double VRTSimpleSource::GetMinimum(int nXSize, int nYSize, int *pbSuccess)
52✔
1440
{
1441
    // The window we will actually request from the source raster band.
1442
    double dfReqXOff = 0.0;
52✔
1443
    double dfReqYOff = 0.0;
52✔
1444
    double dfReqXSize = 0.0;
52✔
1445
    double dfReqYSize = 0.0;
52✔
1446
    int nReqXOff = 0;
52✔
1447
    int nReqYOff = 0;
52✔
1448
    int nReqXSize = 0;
52✔
1449
    int nReqYSize = 0;
52✔
1450

1451
    // The window we will actual set _within_ the pData buffer.
1452
    int nOutXOff = 0;
52✔
1453
    int nOutYOff = 0;
52✔
1454
    int nOutXSize = 0;
52✔
1455
    int nOutYSize = 0;
52✔
1456

1457
    bool bError = false;
52✔
1458
    auto l_band = GetRasterBand();
52✔
1459
    if (!l_band ||
104✔
1460
        !GetSrcDstWindow(0, 0, nXSize, nYSize, nXSize, nYSize, &dfReqXOff,
52✔
1461
                         &dfReqYOff, &dfReqXSize, &dfReqYSize, &nReqXOff,
1462
                         &nReqYOff, &nReqXSize, &nReqYSize, &nOutXOff,
1463
                         &nOutYOff, &nOutXSize, &nOutYSize, bError) ||
52✔
1464
        nReqXOff != 0 || nReqYOff != 0 || nReqXSize != l_band->GetXSize() ||
156✔
1465
        nReqYSize != l_band->GetYSize())
52✔
1466
    {
1467
        *pbSuccess = FALSE;
×
1468
        return 0;
×
1469
    }
1470

1471
    const double dfVal = l_band->GetMinimum(pbSuccess);
52✔
1472
    if (NeedMaxValAdjustment() && dfVal > m_nMaxValue)
52✔
1473
        return m_nMaxValue;
2✔
1474
    return dfVal;
50✔
1475
}
1476

1477
/************************************************************************/
1478
/*                             GetMaximum()                             */
1479
/************************************************************************/
1480

1481
double VRTSimpleSource::GetMaximum(int nXSize, int nYSize, int *pbSuccess)
51✔
1482
{
1483
    // The window we will actually request from the source raster band.
1484
    double dfReqXOff = 0.0;
51✔
1485
    double dfReqYOff = 0.0;
51✔
1486
    double dfReqXSize = 0.0;
51✔
1487
    double dfReqYSize = 0.0;
51✔
1488
    int nReqXOff = 0;
51✔
1489
    int nReqYOff = 0;
51✔
1490
    int nReqXSize = 0;
51✔
1491
    int nReqYSize = 0;
51✔
1492

1493
    // The window we will actual set _within_ the pData buffer.
1494
    int nOutXOff = 0;
51✔
1495
    int nOutYOff = 0;
51✔
1496
    int nOutXSize = 0;
51✔
1497
    int nOutYSize = 0;
51✔
1498

1499
    bool bError = false;
51✔
1500
    auto l_band = GetRasterBand();
51✔
1501
    if (!l_band ||
102✔
1502
        !GetSrcDstWindow(0, 0, nXSize, nYSize, nXSize, nYSize, &dfReqXOff,
51✔
1503
                         &dfReqYOff, &dfReqXSize, &dfReqYSize, &nReqXOff,
1504
                         &nReqYOff, &nReqXSize, &nReqYSize, &nOutXOff,
1505
                         &nOutYOff, &nOutXSize, &nOutYSize, bError) ||
51✔
1506
        nReqXOff != 0 || nReqYOff != 0 || nReqXSize != l_band->GetXSize() ||
153✔
1507
        nReqYSize != l_band->GetYSize())
51✔
1508
    {
1509
        *pbSuccess = FALSE;
×
1510
        return 0;
×
1511
    }
1512

1513
    const double dfVal = l_band->GetMaximum(pbSuccess);
51✔
1514
    if (NeedMaxValAdjustment() && dfVal > m_nMaxValue)
51✔
1515
        return m_nMaxValue;
2✔
1516
    return dfVal;
49✔
1517
}
1518

1519
/************************************************************************/
1520
/*                            GetHistogram()                            */
1521
/************************************************************************/
1522

1523
CPLErr VRTSimpleSource::GetHistogram(int nXSize, int nYSize, double dfMin,
4✔
1524
                                     double dfMax, int nBuckets,
1525
                                     GUIntBig *panHistogram,
1526
                                     int bIncludeOutOfRange, int bApproxOK,
1527
                                     GDALProgressFunc pfnProgress,
1528
                                     void *pProgressData)
1529
{
1530
    // The window we will actually request from the source raster band.
1531
    double dfReqXOff = 0.0;
4✔
1532
    double dfReqYOff = 0.0;
4✔
1533
    double dfReqXSize = 0.0;
4✔
1534
    double dfReqYSize = 0.0;
4✔
1535
    int nReqXOff = 0;
4✔
1536
    int nReqYOff = 0;
4✔
1537
    int nReqXSize = 0;
4✔
1538
    int nReqYSize = 0;
4✔
1539

1540
    // The window we will actual set _within_ the pData buffer.
1541
    int nOutXOff = 0;
4✔
1542
    int nOutYOff = 0;
4✔
1543
    int nOutXSize = 0;
4✔
1544
    int nOutYSize = 0;
4✔
1545

1546
    bool bError = false;
4✔
1547
    auto l_band = GetRasterBand();
4✔
1548
    if (!l_band || NeedMaxValAdjustment() ||
4✔
1549
        !GetSrcDstWindow(0, 0, nXSize, nYSize, nXSize, nYSize, &dfReqXOff,
4✔
1550
                         &dfReqYOff, &dfReqXSize, &dfReqYSize, &nReqXOff,
1551
                         &nReqYOff, &nReqXSize, &nReqYSize, &nOutXOff,
1552
                         &nOutYOff, &nOutXSize, &nOutYSize, bError) ||
4✔
1553
        nReqXOff != 0 || nReqYOff != 0 || nReqXSize != l_band->GetXSize() ||
12✔
1554
        nReqYSize != l_band->GetYSize())
4✔
1555
    {
1556
        return CE_Failure;
×
1557
    }
1558

1559
    return l_band->GetHistogram(dfMin, dfMax, nBuckets, panHistogram,
4✔
1560
                                bIncludeOutOfRange, bApproxOK, pfnProgress,
1561
                                pProgressData);
4✔
1562
}
1563

1564
/************************************************************************/
1565
/*                          DatasetRasterIO()                           */
1566
/************************************************************************/
1567

1568
CPLErr VRTSimpleSource::DatasetRasterIO(
1,002✔
1569
    GDALDataType eVRTBandDataType, int nXOff, int nYOff, int nXSize, int nYSize,
1570
    void *pData, int nBufXSize, int nBufYSize, GDALDataType eBufType,
1571
    int nBandCount, const int *panBandMap, GSpacing nPixelSpace,
1572
    GSpacing nLineSpace, GSpacing nBandSpace,
1573
    GDALRasterIOExtraArg *psExtraArgIn)
1574
{
1575
    if (GetType() != VRTSimpleSource::GetTypeStatic())
1,002✔
1576
    {
1577
        CPLError(CE_Failure, CPLE_NotSupported,
×
1578
                 "DatasetRasterIO() not implemented for %s", GetType());
×
1579
        return CE_Failure;
×
1580
    }
1581

1582
    GDALRasterIOExtraArg sExtraArg;
1583
    INIT_RASTERIO_EXTRA_ARG(sExtraArg);
1,002✔
1584
    GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
1,002✔
1585

1586
    double dfXOff = nXOff;
1,002✔
1587
    double dfYOff = nYOff;
1,002✔
1588
    double dfXSize = nXSize;
1,002✔
1589
    double dfYSize = nYSize;
1,002✔
1590
    if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
1,002✔
1591
    {
1592
        dfXOff = psExtraArgIn->dfXOff;
11✔
1593
        dfYOff = psExtraArgIn->dfYOff;
11✔
1594
        dfXSize = psExtraArgIn->dfXSize;
11✔
1595
        dfYSize = psExtraArgIn->dfYSize;
11✔
1596
    }
1597

1598
    // The window we will actually request from the source raster band.
1599
    double dfReqXOff = 0.0;
1,002✔
1600
    double dfReqYOff = 0.0;
1,002✔
1601
    double dfReqXSize = 0.0;
1,002✔
1602
    double dfReqYSize = 0.0;
1,002✔
1603
    int nReqXOff = 0;
1,002✔
1604
    int nReqYOff = 0;
1,002✔
1605
    int nReqXSize = 0;
1,002✔
1606
    int nReqYSize = 0;
1,002✔
1607

1608
    // The window we will actual set _within_ the pData buffer.
1609
    int nOutXOff = 0;
1,002✔
1610
    int nOutYOff = 0;
1,002✔
1611
    int nOutXSize = 0;
1,002✔
1612
    int nOutYSize = 0;
1,002✔
1613

1614
    bool bError = false;
1,002✔
1615
    if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
1,002✔
1616
                         &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
1617
                         &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
1618
                         &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
1619
    {
1620
        return bError ? CE_Failure : CE_None;
71✔
1621
    }
1622

1623
    auto l_band = GetRasterBand();
931✔
1624
    if (!l_band)
931✔
1625
        return CE_Failure;
×
1626

1627
    GDALDataset *poDS = l_band->GetDataset();
931✔
1628
    if (poDS == nullptr)
931✔
1629
        return CE_Failure;
×
1630

1631
    if (!m_osResampling.empty())
931✔
1632
    {
1633
        psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
6✔
1634
    }
1635
    else if (psExtraArgIn != nullptr)
925✔
1636
    {
1637
        psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
925✔
1638
    }
1639
    psExtraArg->bFloatingPointWindowValidity = TRUE;
931✔
1640
    psExtraArg->dfXOff = dfReqXOff;
931✔
1641
    psExtraArg->dfYOff = dfReqYOff;
931✔
1642
    psExtraArg->dfXSize = dfReqXSize;
931✔
1643
    psExtraArg->dfYSize = dfReqYSize;
931✔
1644
    if (psExtraArgIn)
931✔
1645
    {
1646
        psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
931✔
1647
        psExtraArg->pProgressData = psExtraArgIn->pProgressData;
931✔
1648
    }
1649

1650
    GByte *pabyOut = static_cast<unsigned char *>(pData) +
931✔
1651
                     nOutXOff * nPixelSpace +
931✔
1652
                     static_cast<GPtrDiff_t>(nOutYOff) * nLineSpace;
931✔
1653

1654
    CPLErr eErr = CE_Failure;
931✔
1655

1656
    if (GDALDataTypeIsConversionLossy(l_band->GetRasterDataType(),
931✔
1657
                                      eVRTBandDataType))
931✔
1658
    {
1659
        const int nBandDTSize = GDALGetDataTypeSizeBytes(eVRTBandDataType);
1✔
1660
        void *pTemp = VSI_MALLOC3_VERBOSE(
1✔
1661
            nOutXSize, nOutYSize, cpl::fits_on<int>(nBandDTSize * nBandCount));
1662
        if (pTemp)
1✔
1663
        {
1664
            eErr = poDS->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize,
1✔
1665
                                  nReqYSize, pTemp, nOutXSize, nOutYSize,
1666
                                  eVRTBandDataType, nBandCount, panBandMap, 0,
1667
                                  0, 0, psExtraArg);
1668
            if (eErr == CE_None)
1✔
1669
            {
1670
                GByte *pabyTemp = static_cast<GByte *>(pTemp);
1✔
1671
                const size_t nSrcBandSpace =
1✔
1672
                    static_cast<size_t>(nOutYSize) * nOutXSize * nBandDTSize;
1✔
1673
                for (int iBand = 0; iBand < nBandCount; iBand++)
2✔
1674
                {
1675
                    for (int iY = 0; iY < nOutYSize; iY++)
21✔
1676
                    {
1677
                        GDALCopyWords(
20✔
1678
                            pabyTemp + iBand * nSrcBandSpace +
20✔
1679
                                static_cast<size_t>(iY) * nBandDTSize *
20✔
1680
                                    nOutXSize,
20✔
1681
                            eVRTBandDataType, nBandDTSize,
1682
                            pabyOut + static_cast<GPtrDiff_t>(
20✔
1683
                                          iY * nLineSpace + iBand * nBandSpace),
20✔
1684
                            eBufType, static_cast<int>(nPixelSpace), nOutXSize);
1685
                    }
1686
                }
1687
            }
1688
            VSIFree(pTemp);
1✔
1689
        }
1690
    }
1691
    else
1692
    {
1693
        eErr = poDS->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
930✔
1694
                              pabyOut, nOutXSize, nOutYSize, eBufType,
1695
                              nBandCount, panBandMap, nPixelSpace, nLineSpace,
1696
                              nBandSpace, psExtraArg);
1697
    }
1698

1699
    if (NeedMaxValAdjustment())
931✔
1700
    {
1701
        for (int k = 0; k < nBandCount; k++)
×
1702
        {
1703
            for (int j = 0; j < nOutYSize; j++)
×
1704
            {
1705
                for (int i = 0; i < nOutXSize; i++)
×
1706
                {
1707
                    int nVal = 0;
×
1708
                    GDALCopyWords(pabyOut + k * nBandSpace + j * nLineSpace +
×
1709
                                      i * nPixelSpace,
×
1710
                                  eBufType, 0, &nVal, GDT_Int32, 0, 1);
1711

1712
                    if (nVal > m_nMaxValue)
×
1713
                        nVal = m_nMaxValue;
×
1714

1715
                    GDALCopyWords(&nVal, GDT_Int32, 0,
×
1716
                                  pabyOut + k * nBandSpace + j * nLineSpace +
×
1717
                                      i * nPixelSpace,
×
1718
                                  eBufType, 0, 1);
1719
                }
1720
            }
1721
        }
1722
    }
1723

1724
    if (psExtraArg->pfnProgress)
931✔
1725
        psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
470✔
1726

1727
    return eErr;
931✔
1728
}
1729

1730
/************************************************************************/
1731
/*                          SetResampling()                             */
1732
/************************************************************************/
1733

1734
void VRTSimpleSource::SetResampling(const char *pszResampling)
4,101✔
1735
{
1736
    m_osResampling = (pszResampling) ? pszResampling : "";
4,101✔
1737
}
4,101✔
1738

1739
/************************************************************************/
1740
/* ==================================================================== */
1741
/*                         VRTAveragedSource                            */
1742
/* ==================================================================== */
1743
/************************************************************************/
1744

1745
/************************************************************************/
1746
/*                         VRTAveragedSource()                          */
1747
/************************************************************************/
1748

1749
VRTAveragedSource::VRTAveragedSource()
16✔
1750
{
1751
}
16✔
1752

1753
/************************************************************************/
1754
/*                           GetTypeStatic()                            */
1755
/************************************************************************/
1756

1757
const char *VRTAveragedSource::GetTypeStatic()
3,166✔
1758
{
1759
    static const char *TYPE = "AveragedSource";
1760
    return TYPE;
3,166✔
1761
}
1762

1763
/************************************************************************/
1764
/*                            GetType()                                 */
1765
/************************************************************************/
1766

1767
const char *VRTAveragedSource::GetType() const
11✔
1768
{
1769
    return GetTypeStatic();
11✔
1770
}
1771

1772
/************************************************************************/
1773
/*                           SerializeToXML()                           */
1774
/************************************************************************/
1775

1776
CPLXMLNode *VRTAveragedSource::SerializeToXML(const char *pszVRTPath)
×
1777

1778
{
1779
    CPLXMLNode *const psSrc = VRTSimpleSource::SerializeToXML(pszVRTPath);
×
1780

1781
    if (psSrc == nullptr)
×
1782
        return nullptr;
×
1783

1784
    CPLFree(psSrc->pszValue);
×
1785
    psSrc->pszValue = CPLStrdup(GetTypeStatic());
×
1786

1787
    return psSrc;
×
1788
}
1789

1790
/************************************************************************/
1791
/*                           SetNoDataValue()                           */
1792
/************************************************************************/
1793

1794
void VRTAveragedSource::SetNoDataValue(double dfNewNoDataValue)
×
1795

1796
{
1797
    if (dfNewNoDataValue == VRT_NODATA_UNSET)
×
1798
    {
1799
        m_bNoDataSet = FALSE;
×
1800
        m_dfNoDataValue = VRT_NODATA_UNSET;
×
1801
        return;
×
1802
    }
1803

1804
    m_bNoDataSet = TRUE;
×
1805
    m_dfNoDataValue = dfNewNoDataValue;
×
1806
}
1807

1808
/************************************************************************/
1809
/*                              RasterIO()                              */
1810
/************************************************************************/
1811

1812
CPLErr VRTAveragedSource::RasterIO(GDALDataType /*eVRTBandDataType*/, int nXOff,
33✔
1813
                                   int nYOff, int nXSize, int nYSize,
1814
                                   void *pData, int nBufXSize, int nBufYSize,
1815
                                   GDALDataType eBufType, GSpacing nPixelSpace,
1816
                                   GSpacing nLineSpace,
1817
                                   GDALRasterIOExtraArg *psExtraArgIn,
1818
                                   WorkingState & /*oWorkingState*/)
1819

1820
{
1821
    GDALRasterIOExtraArg sExtraArg;
1822
    INIT_RASTERIO_EXTRA_ARG(sExtraArg);
33✔
1823
    GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
33✔
1824

1825
    double dfXOff = nXOff;
33✔
1826
    double dfYOff = nYOff;
33✔
1827
    double dfXSize = nXSize;
33✔
1828
    double dfYSize = nYSize;
33✔
1829
    if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
33✔
1830
    {
1831
        dfXOff = psExtraArgIn->dfXOff;
×
1832
        dfYOff = psExtraArgIn->dfYOff;
×
1833
        dfXSize = psExtraArgIn->dfXSize;
×
1834
        dfYSize = psExtraArgIn->dfYSize;
×
1835
    }
1836

1837
    // The window we will actually request from the source raster band.
1838
    double dfReqXOff = 0.0;
33✔
1839
    double dfReqYOff = 0.0;
33✔
1840
    double dfReqXSize = 0.0;
33✔
1841
    double dfReqYSize = 0.0;
33✔
1842
    int nReqXOff = 0;
33✔
1843
    int nReqYOff = 0;
33✔
1844
    int nReqXSize = 0;
33✔
1845
    int nReqYSize = 0;
33✔
1846

1847
    // The window we will actual set _within_ the pData buffer.
1848
    int nOutXOff = 0;
33✔
1849
    int nOutYOff = 0;
33✔
1850
    int nOutXSize = 0;
33✔
1851
    int nOutYSize = 0;
33✔
1852

1853
    bool bError = false;
33✔
1854
    if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
33✔
1855
                         &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
1856
                         &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
1857
                         &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
1858
    {
1859
        return bError ? CE_Failure : CE_None;
×
1860
    }
1861

1862
    auto l_band = GetRasterBand();
33✔
1863
    if (!l_band)
33✔
1864
        return CE_Failure;
×
1865

1866
    /* -------------------------------------------------------------------- */
1867
    /*      Allocate a temporary buffer to whole the full resolution        */
1868
    /*      data from the area of interest.                                 */
1869
    /* -------------------------------------------------------------------- */
1870
    float *const pafSrc = static_cast<float *>(
1871
        VSI_MALLOC3_VERBOSE(sizeof(float), nReqXSize, nReqYSize));
33✔
1872
    if (pafSrc == nullptr)
33✔
1873
    {
1874
        return CE_Failure;
×
1875
    }
1876

1877
    /* -------------------------------------------------------------------- */
1878
    /*      Load it.                                                        */
1879
    /* -------------------------------------------------------------------- */
1880
    if (!m_osResampling.empty())
33✔
1881
    {
1882
        psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
28✔
1883
    }
1884
    else if (psExtraArgIn != nullptr)
5✔
1885
    {
1886
        psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
5✔
1887
    }
1888

1889
    psExtraArg->bFloatingPointWindowValidity = TRUE;
33✔
1890
    psExtraArg->dfXOff = dfReqXOff;
33✔
1891
    psExtraArg->dfYOff = dfReqYOff;
33✔
1892
    psExtraArg->dfXSize = dfReqXSize;
33✔
1893
    psExtraArg->dfYSize = dfReqYSize;
33✔
1894
    if (psExtraArgIn)
33✔
1895
    {
1896
        psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
33✔
1897
        psExtraArg->pProgressData = psExtraArgIn->pProgressData;
33✔
1898
    }
1899

1900
    const CPLErr eErr = l_band->RasterIO(
33✔
1901
        GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize, pafSrc, nReqXSize,
1902
        nReqYSize, GDT_Float32, 0, 0, psExtraArg);
1903

1904
    if (eErr != CE_None)
33✔
1905
    {
1906
        VSIFree(pafSrc);
×
1907
        return eErr;
×
1908
    }
1909

1910
    /* -------------------------------------------------------------------- */
1911
    /*      Do the averaging.                                               */
1912
    /* -------------------------------------------------------------------- */
1913
    for (int iBufLine = nOutYOff; iBufLine < nOutYOff + nOutYSize; iBufLine++)
5,956✔
1914
    {
1915
        const double dfYDst =
5,923✔
1916
            (iBufLine / static_cast<double>(nBufYSize)) * nYSize + nYOff;
5,923✔
1917

1918
        for (int iBufPixel = nOutXOff; iBufPixel < nOutXOff + nOutXSize;
2,017,080✔
1919
             iBufPixel++)
1920
        {
1921
            double dfXSrcStart, dfXSrcEnd, dfYSrcStart, dfYSrcEnd;
1922
            int iXSrcStart, iYSrcStart, iXSrcEnd, iYSrcEnd;
1923

1924
            const double dfXDst =
2,011,160✔
1925
                (iBufPixel / static_cast<double>(nBufXSize)) * nXSize + nXOff;
2,011,160✔
1926

1927
            // Compute the source image rectangle needed for this pixel.
1928
            DstToSrc(dfXDst, dfYDst, dfXSrcStart, dfYSrcStart);
2,011,160✔
1929
            DstToSrc(dfXDst + 1.0, dfYDst + 1.0, dfXSrcEnd, dfYSrcEnd);
2,011,160✔
1930

1931
            // Convert to integers, assuming that the center of the source
1932
            // pixel must be in our rect to get included.
1933
            if (dfXSrcEnd >= dfXSrcStart + 1)
2,011,160✔
1934
            {
1935
                iXSrcStart = static_cast<int>(floor(dfXSrcStart + 0.5));
1,049,560✔
1936
                iXSrcEnd = static_cast<int>(floor(dfXSrcEnd + 0.5));
1,049,560✔
1937
            }
1938
            else
1939
            {
1940
                /* If the resampling factor is less than 100%, the distance */
1941
                /* between the source pixel is < 1, so we stick to nearest */
1942
                /* neighbour */
1943
                iXSrcStart = static_cast<int>(floor(dfXSrcStart));
961,600✔
1944
                iXSrcEnd = iXSrcStart + 1;
961,600✔
1945
            }
1946
            if (dfYSrcEnd >= dfYSrcStart + 1)
2,011,160✔
1947
            {
1948
                iYSrcStart = static_cast<int>(floor(dfYSrcStart + 0.5));
1,049,560✔
1949
                iYSrcEnd = static_cast<int>(floor(dfYSrcEnd + 0.5));
1,049,560✔
1950
            }
1951
            else
1952
            {
1953
                iYSrcStart = static_cast<int>(floor(dfYSrcStart));
961,600✔
1954
                iYSrcEnd = iYSrcStart + 1;
961,600✔
1955
            }
1956

1957
            // Transform into the coordinate system of the source *buffer*
1958
            iXSrcStart -= nReqXOff;
2,011,160✔
1959
            iYSrcStart -= nReqYOff;
2,011,160✔
1960
            iXSrcEnd -= nReqXOff;
2,011,160✔
1961
            iYSrcEnd -= nReqYOff;
2,011,160✔
1962

1963
            double dfSum = 0.0;
2,011,160✔
1964
            int nPixelCount = 0;
2,011,160✔
1965

1966
            for (int iY = iYSrcStart; iY < iYSrcEnd; iY++)
4,022,510✔
1967
            {
1968
                if (iY < 0 || iY >= nReqYSize)
2,011,360✔
1969
                    continue;
×
1970

1971
                for (int iX = iXSrcStart; iX < iXSrcEnd; iX++)
4,023,130✔
1972
                {
1973
                    if (iX < 0 || iX >= nReqXSize)
2,011,780✔
1974
                        continue;
×
1975

1976
                    const float fSampledValue =
2,011,780✔
1977
                        pafSrc[iX + static_cast<size_t>(iY) * nReqXSize];
2,011,780✔
1978
                    if (std::isnan(fSampledValue))
2,011,780✔
1979
                        continue;
×
1980

1981
                    if (m_bNoDataSet &&
4,023,550✔
1982
                        GDALIsValueInRange<float>(m_dfNoDataValue) &&
2,011,780✔
1983
                        ARE_REAL_EQUAL(fSampledValue,
×
1984
                                       static_cast<float>(m_dfNoDataValue)))
×
1985
                        continue;
×
1986

1987
                    nPixelCount++;
2,011,780✔
1988
                    dfSum += pafSrc[iX + static_cast<size_t>(iY) * nReqXSize];
2,011,780✔
1989
                }
1990
            }
1991

1992
            if (nPixelCount == 0)
2,011,160✔
1993
                continue;
×
1994

1995
            // Compute output value.
1996
            const float dfOutputValue = static_cast<float>(dfSum / nPixelCount);
2,011,160✔
1997

1998
            // Put it in the output buffer.
1999
            GByte *pDstLocation =
2,011,160✔
2000
                static_cast<GByte *>(pData) + nPixelSpace * iBufPixel +
2,011,160✔
2001
                static_cast<GPtrDiff_t>(nLineSpace) * iBufLine;
2,011,160✔
2002

2003
            if (eBufType == GDT_Byte)
2,011,160✔
2004
                *pDstLocation = static_cast<GByte>(
2,008,660✔
2005
                    std::min(255.0, std::max(0.0, dfOutputValue + 0.5)));
2,008,660✔
2006
            else
2007
                GDALCopyWords(&dfOutputValue, GDT_Float32, 4, pDstLocation,
2,500✔
2008
                              eBufType, 8, 1);
2009
        }
2010
    }
2011

2012
    VSIFree(pafSrc);
33✔
2013

2014
    if (psExtraArg->pfnProgress)
33✔
2015
        psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
×
2016

2017
    return CE_None;
33✔
2018
}
2019

2020
/************************************************************************/
2021
/*                             GetMinimum()                             */
2022
/************************************************************************/
2023

2024
double VRTAveragedSource::GetMinimum(int /* nXSize */, int /* nYSize */,
×
2025
                                     int *pbSuccess)
2026
{
2027
    *pbSuccess = FALSE;
×
2028
    return 0.0;
×
2029
}
2030

2031
/************************************************************************/
2032
/*                             GetMaximum()                             */
2033
/************************************************************************/
2034

2035
double VRTAveragedSource::GetMaximum(int /* nXSize */, int /* nYSize */,
×
2036
                                     int *pbSuccess)
2037
{
2038
    *pbSuccess = FALSE;
×
2039
    return 0.0;
×
2040
}
2041

2042
/************************************************************************/
2043
/*                            GetHistogram()                            */
2044
/************************************************************************/
2045

2046
CPLErr VRTAveragedSource::GetHistogram(
×
2047
    int /* nXSize */, int /* nYSize */, double /* dfMin */, double /* dfMax */,
2048
    int /* nBuckets */, GUIntBig * /* panHistogram */,
2049
    int /* bIncludeOutOfRange */, int /* bApproxOK */,
2050
    GDALProgressFunc /* pfnProgress */, void * /* pProgressData */)
2051
{
2052
    return CE_Failure;
×
2053
}
2054

2055
/************************************************************************/
2056
/* ==================================================================== */
2057
/*                     VRTNoDataFromMaskSource                          */
2058
/* ==================================================================== */
2059
/************************************************************************/
2060

2061
/************************************************************************/
2062
/*                     VRTNoDataFromMaskSource()                        */
2063
/************************************************************************/
2064

2065
VRTNoDataFromMaskSource::VRTNoDataFromMaskSource()
23✔
2066
{
2067
}
23✔
2068

2069
/************************************************************************/
2070
/*                              XMLInit()                               */
2071
/************************************************************************/
2072

2073
CPLErr
2074
VRTNoDataFromMaskSource::XMLInit(const CPLXMLNode *psSrc,
8✔
2075
                                 const char *pszVRTPath,
2076
                                 VRTMapSharedResources &oMapSharedSources)
2077

2078
{
2079
    /* -------------------------------------------------------------------- */
2080
    /*      Do base initialization.                                         */
2081
    /* -------------------------------------------------------------------- */
2082
    {
2083
        const CPLErr eErr =
2084
            VRTSimpleSource::XMLInit(psSrc, pszVRTPath, oMapSharedSources);
8✔
2085
        if (eErr != CE_None)
8✔
2086
            return eErr;
×
2087
    }
2088

2089
    if (const char *pszNODATA = CPLGetXMLValue(psSrc, "NODATA", nullptr))
8✔
2090
    {
2091
        m_bNoDataSet = true;
8✔
2092
        m_dfNoDataValue = CPLAtofM(pszNODATA);
8✔
2093
    }
2094

2095
    m_dfMaskValueThreshold =
8✔
2096
        CPLAtofM(CPLGetXMLValue(psSrc, "MaskValueThreshold", "0"));
8✔
2097

2098
    if (const char *pszRemappedValue =
8✔
2099
            CPLGetXMLValue(psSrc, "RemappedValue", nullptr))
8✔
2100
    {
2101
        m_bHasRemappedValue = true;
×
2102
        m_dfRemappedValue = CPLAtofM(pszRemappedValue);
×
2103
    }
2104

2105
    return CE_None;
8✔
2106
}
2107

2108
/************************************************************************/
2109
/*                           GetTypeStatic()                            */
2110
/************************************************************************/
2111

2112
const char *VRTNoDataFromMaskSource::GetTypeStatic()
27✔
2113
{
2114
    static const char *TYPE = "NoDataFromMaskSource";
2115
    return TYPE;
27✔
2116
}
2117

2118
/************************************************************************/
2119
/*                            GetType()                                 */
2120
/************************************************************************/
2121

2122
const char *VRTNoDataFromMaskSource::GetType() const
11✔
2123
{
2124
    return GetTypeStatic();
11✔
2125
}
2126

2127
/************************************************************************/
2128
/*                           SerializeToXML()                           */
2129
/************************************************************************/
2130

2131
CPLXMLNode *VRTNoDataFromMaskSource::SerializeToXML(const char *pszVRTPath)
8✔
2132

2133
{
2134
    CPLXMLNode *const psSrc = VRTSimpleSource::SerializeToXML(pszVRTPath);
8✔
2135

2136
    if (psSrc == nullptr)
8✔
2137
        return nullptr;
×
2138

2139
    CPLFree(psSrc->pszValue);
8✔
2140
    psSrc->pszValue = CPLStrdup(GetTypeStatic());
8✔
2141

2142
    if (m_bNoDataSet)
8✔
2143
    {
2144
        CPLSetXMLValue(psSrc, "MaskValueThreshold",
8✔
2145
                       CPLSPrintf("%.17g", m_dfMaskValueThreshold));
2146

2147
        GDALDataType eBandDT = GDT_Unknown;
8✔
2148
        double dfNoDataValue = m_dfNoDataValue;
8✔
2149
        const auto kMaxFloat = std::numeric_limits<float>::max();
8✔
2150
        if (std::fabs(std::fabs(m_dfNoDataValue) - kMaxFloat) <
8✔
2151
            1e-10 * kMaxFloat)
2152
        {
2153
            auto l_band = GetRasterBand();
×
2154
            if (l_band)
×
2155
            {
2156
                eBandDT = l_band->GetRasterDataType();
×
2157
                if (eBandDT == GDT_Float32)
×
2158
                {
2159
                    dfNoDataValue =
2160
                        GDALAdjustNoDataCloseToFloatMax(m_dfNoDataValue);
×
2161
                }
2162
            }
2163
        }
2164
        CPLSetXMLValue(psSrc, "NODATA",
8✔
2165
                       VRTSerializeNoData(dfNoDataValue, eBandDT, 18).c_str());
16✔
2166
    }
2167

2168
    if (m_bHasRemappedValue)
8✔
2169
    {
2170
        CPLSetXMLValue(psSrc, "RemappedValue",
×
2171
                       CPLSPrintf("%.17g", m_dfRemappedValue));
2172
    }
2173

2174
    return psSrc;
8✔
2175
}
2176

2177
/************************************************************************/
2178
/*                           SetParameters()                            */
2179
/************************************************************************/
2180

2181
void VRTNoDataFromMaskSource::SetParameters(double dfNoDataValue,
15✔
2182
                                            double dfMaskValueThreshold)
2183
{
2184
    m_bNoDataSet = true;
15✔
2185
    m_dfNoDataValue = dfNoDataValue;
15✔
2186
    m_dfMaskValueThreshold = dfMaskValueThreshold;
15✔
2187
    if (!m_bHasRemappedValue)
15✔
2188
        m_dfRemappedValue = m_dfNoDataValue;
15✔
2189
}
15✔
2190

2191
/************************************************************************/
2192
/*                           SetParameters()                            */
2193
/************************************************************************/
2194

2195
void VRTNoDataFromMaskSource::SetParameters(double dfNoDataValue,
×
2196
                                            double dfMaskValueThreshold,
2197
                                            double dfRemappedValue)
2198
{
2199
    SetParameters(dfNoDataValue, dfMaskValueThreshold);
×
2200
    m_bHasRemappedValue = true;
×
2201
    m_dfRemappedValue = dfRemappedValue;
×
2202
}
×
2203

2204
/************************************************************************/
2205
/*                              RasterIO()                              */
2206
/************************************************************************/
2207

2208
CPLErr VRTNoDataFromMaskSource::RasterIO(
13✔
2209
    GDALDataType eVRTBandDataType, int nXOff, int nYOff, int nXSize, int nYSize,
2210
    void *pData, int nBufXSize, int nBufYSize, GDALDataType eBufType,
2211
    GSpacing nPixelSpace, GSpacing nLineSpace,
2212
    GDALRasterIOExtraArg *psExtraArgIn, WorkingState &oWorkingState)
2213

2214
{
2215
    if (!m_bNoDataSet)
13✔
2216
    {
2217
        return VRTSimpleSource::RasterIO(eVRTBandDataType, nXOff, nYOff, nXSize,
×
2218
                                         nYSize, pData, nBufXSize, nBufYSize,
2219
                                         eBufType, nPixelSpace, nLineSpace,
2220
                                         psExtraArgIn, oWorkingState);
×
2221
    }
2222

2223
    GDALRasterIOExtraArg sExtraArg;
2224
    INIT_RASTERIO_EXTRA_ARG(sExtraArg);
13✔
2225
    GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
13✔
2226

2227
    double dfXOff = nXOff;
13✔
2228
    double dfYOff = nYOff;
13✔
2229
    double dfXSize = nXSize;
13✔
2230
    double dfYSize = nYSize;
13✔
2231
    if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
13✔
2232
    {
2233
        dfXOff = psExtraArgIn->dfXOff;
×
2234
        dfYOff = psExtraArgIn->dfYOff;
×
2235
        dfXSize = psExtraArgIn->dfXSize;
×
2236
        dfYSize = psExtraArgIn->dfYSize;
×
2237
    }
2238

2239
    // The window we will actually request from the source raster band.
2240
    double dfReqXOff = 0.0;
13✔
2241
    double dfReqYOff = 0.0;
13✔
2242
    double dfReqXSize = 0.0;
13✔
2243
    double dfReqYSize = 0.0;
13✔
2244
    int nReqXOff = 0;
13✔
2245
    int nReqYOff = 0;
13✔
2246
    int nReqXSize = 0;
13✔
2247
    int nReqYSize = 0;
13✔
2248

2249
    // The window we will actual set _within_ the pData buffer.
2250
    int nOutXOff = 0;
13✔
2251
    int nOutYOff = 0;
13✔
2252
    int nOutXSize = 0;
13✔
2253
    int nOutYSize = 0;
13✔
2254

2255
    bool bError = false;
13✔
2256
    if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
13✔
2257
                         &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
2258
                         &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
2259
                         &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
2260
    {
2261
        return bError ? CE_Failure : CE_None;
×
2262
    }
2263

2264
    auto l_band = GetRasterBand();
13✔
2265
    if (!l_band)
13✔
2266
        return CE_Failure;
×
2267

2268
    /* -------------------------------------------------------------------- */
2269
    /*      Allocate temporary buffer(s).                                   */
2270
    /* -------------------------------------------------------------------- */
2271
    const auto eSrcBandDT = l_band->GetRasterDataType();
13✔
2272
    const int nSrcBandDTSize = GDALGetDataTypeSizeBytes(eSrcBandDT);
13✔
2273
    const auto eSrcMaskBandDT = l_band->GetMaskBand()->GetRasterDataType();
13✔
2274
    const int nSrcMaskBandDTSize = GDALGetDataTypeSizeBytes(eSrcMaskBandDT);
13✔
2275
    double dfRemappedValue = m_dfRemappedValue;
13✔
2276
    if (!m_bHasRemappedValue)
13✔
2277
    {
2278
        if (eSrcBandDT == GDT_Byte &&
19✔
2279
            m_dfNoDataValue >= std::numeric_limits<GByte>::min() &&
6✔
2280
            m_dfNoDataValue <= std::numeric_limits<GByte>::max() &&
25✔
2281
            static_cast<int>(m_dfNoDataValue) == m_dfNoDataValue)
6✔
2282
        {
2283
            if (m_dfNoDataValue == std::numeric_limits<GByte>::max())
6✔
2284
                dfRemappedValue = m_dfNoDataValue - 1;
1✔
2285
            else
2286
                dfRemappedValue = m_dfNoDataValue + 1;
5✔
2287
        }
2288
        else if (eSrcBandDT == GDT_UInt16 &&
10✔
2289
                 m_dfNoDataValue >= std::numeric_limits<uint16_t>::min() &&
3✔
2290
                 m_dfNoDataValue <= std::numeric_limits<uint16_t>::max() &&
13✔
2291
                 static_cast<int>(m_dfNoDataValue) == m_dfNoDataValue)
3✔
2292
        {
2293
            if (m_dfNoDataValue == std::numeric_limits<uint16_t>::max())
3✔
2294
                dfRemappedValue = m_dfNoDataValue - 1;
1✔
2295
            else
2296
                dfRemappedValue = m_dfNoDataValue + 1;
2✔
2297
        }
2298
        else if (eSrcBandDT == GDT_Int16 &&
6✔
2299
                 m_dfNoDataValue >= std::numeric_limits<int16_t>::min() &&
2✔
2300
                 m_dfNoDataValue <= std::numeric_limits<int16_t>::max() &&
8✔
2301
                 static_cast<int>(m_dfNoDataValue) == m_dfNoDataValue)
2✔
2302
        {
2303
            if (m_dfNoDataValue == std::numeric_limits<int16_t>::max())
2✔
2304
                dfRemappedValue = m_dfNoDataValue - 1;
1✔
2305
            else
2306
                dfRemappedValue = m_dfNoDataValue + 1;
1✔
2307
        }
2308
        else
2309
        {
2310
            constexpr double EPS = 1e-3;
2✔
2311
            if (m_dfNoDataValue == 0)
2✔
2312
                dfRemappedValue = EPS;
1✔
2313
            else
2314
                dfRemappedValue = m_dfNoDataValue * (1 + EPS);
1✔
2315
        }
2316
    }
2317
    const bool bByteOptim =
13✔
2318
        (eSrcBandDT == GDT_Byte && eBufType == GDT_Byte &&
6✔
2319
         eSrcMaskBandDT == GDT_Byte && m_dfMaskValueThreshold >= 0 &&
5✔
2320
         m_dfMaskValueThreshold <= 255 &&
5✔
2321
         static_cast<int>(m_dfMaskValueThreshold) == m_dfMaskValueThreshold &&
5✔
2322
         m_dfNoDataValue >= 0 && m_dfNoDataValue <= 255 &&
4✔
2323
         static_cast<int>(m_dfNoDataValue) == m_dfNoDataValue &&
4✔
2324
         dfRemappedValue >= 0 && dfRemappedValue <= 255 &&
23✔
2325
         static_cast<int>(dfRemappedValue) == dfRemappedValue);
4✔
2326
    GByte *pabyWrkBuffer;
2327
    try
2328
    {
2329
        if (bByteOptim && nOutXOff == 0 && nOutYOff == 0 &&
13✔
2330
            nOutXSize == nBufXSize && nOutYSize == nBufYSize &&
4✔
2331
            eSrcBandDT == eBufType && nPixelSpace == nSrcBandDTSize &&
4✔
2332
            nLineSpace == nPixelSpace * nBufXSize)
4✔
2333
        {
2334
            pabyWrkBuffer = static_cast<GByte *>(pData);
4✔
2335
        }
2336
        else
2337
        {
2338
            oWorkingState.m_abyWrkBuffer.resize(static_cast<size_t>(nOutXSize) *
9✔
2339
                                                nOutYSize * nSrcBandDTSize);
9✔
2340
            pabyWrkBuffer =
2341
                reinterpret_cast<GByte *>(oWorkingState.m_abyWrkBuffer.data());
9✔
2342
        }
2343
        oWorkingState.m_abyWrkBufferMask.resize(static_cast<size_t>(nOutXSize) *
13✔
2344
                                                nOutYSize * nSrcMaskBandDTSize);
13✔
2345
    }
2346
    catch (const std::exception &)
×
2347
    {
2348
        CPLError(CE_Failure, CPLE_OutOfMemory,
×
2349
                 "Out of memory when allocating buffers");
2350
        return CE_Failure;
×
2351
    }
2352

2353
    /* -------------------------------------------------------------------- */
2354
    /*      Load data.                                                      */
2355
    /* -------------------------------------------------------------------- */
2356
    if (!m_osResampling.empty())
13✔
2357
    {
2358
        psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
×
2359
    }
2360
    else if (psExtraArgIn != nullptr)
13✔
2361
    {
2362
        psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
13✔
2363
    }
2364

2365
    psExtraArg->bFloatingPointWindowValidity = TRUE;
13✔
2366
    psExtraArg->dfXOff = dfReqXOff;
13✔
2367
    psExtraArg->dfYOff = dfReqYOff;
13✔
2368
    psExtraArg->dfXSize = dfReqXSize;
13✔
2369
    psExtraArg->dfYSize = dfReqYSize;
13✔
2370
    if (psExtraArgIn)
13✔
2371
    {
2372
        psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
13✔
2373
        psExtraArg->pProgressData = psExtraArgIn->pProgressData;
13✔
2374
    }
2375

2376
    if (l_band->RasterIO(GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
13✔
2377
                         pabyWrkBuffer, nOutXSize, nOutYSize, eSrcBandDT, 0, 0,
2378
                         psExtraArg) != CE_None)
13✔
2379
    {
2380
        return CE_Failure;
×
2381
    }
2382

2383
    if (l_band->GetMaskBand()->RasterIO(
26✔
2384
            GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
2385
            oWorkingState.m_abyWrkBufferMask.data(), nOutXSize, nOutYSize,
13✔
2386
            eSrcMaskBandDT, 0, 0, psExtraArg) != CE_None)
13✔
2387
    {
2388
        return CE_Failure;
×
2389
    }
2390

2391
    /* -------------------------------------------------------------------- */
2392
    /*      Do the processing.                                              */
2393
    /* -------------------------------------------------------------------- */
2394

2395
    GByte *const pabyOut = static_cast<GByte *>(pData) +
13✔
2396
                           nPixelSpace * nOutXOff +
13✔
2397
                           static_cast<GPtrDiff_t>(nLineSpace) * nOutYOff;
13✔
2398
    if (bByteOptim)
13✔
2399
    {
2400
        // Special case when everything fits on Byte
2401
        const GByte nMaskValueThreshold =
4✔
2402
            static_cast<GByte>(m_dfMaskValueThreshold);
4✔
2403
        const GByte nNoDataValue = static_cast<GByte>(m_dfNoDataValue);
4✔
2404
        const GByte nRemappedValue = static_cast<GByte>(dfRemappedValue);
4✔
2405
        size_t nSrcIdx = 0;
4✔
2406
        for (int iY = 0; iY < nOutYSize; iY++)
8✔
2407
        {
2408
            GSpacing nDstOffset = iY * nLineSpace;
4✔
2409
            for (int iX = 0; iX < nOutXSize; iX++)
12✔
2410
            {
2411
                const GByte nMaskVal =
2412
                    oWorkingState.m_abyWrkBufferMask[nSrcIdx];
8✔
2413
                if (nMaskVal <= nMaskValueThreshold)
8✔
2414
                {
2415
                    pabyOut[static_cast<GPtrDiff_t>(nDstOffset)] = nNoDataValue;
4✔
2416
                }
2417
                else
2418
                {
2419
                    if (pabyWrkBuffer[nSrcIdx] == nNoDataValue)
4✔
2420
                    {
2421
                        pabyOut[static_cast<GPtrDiff_t>(nDstOffset)] =
2✔
2422
                            nRemappedValue;
2423
                    }
2424
                    else
2425
                    {
2426
                        pabyOut[static_cast<GPtrDiff_t>(nDstOffset)] =
2✔
2427
                            pabyWrkBuffer[nSrcIdx];
2✔
2428
                    }
2429
                }
2430
                nDstOffset += nPixelSpace;
8✔
2431
                nSrcIdx++;
8✔
2432
            }
2433
        }
2434
    }
2435
    else
2436
    {
2437
        size_t nSrcIdx = 0;
9✔
2438
        double dfMaskVal = 0;
9✔
2439
        const int nBufDTSize = GDALGetDataTypeSizeBytes(eBufType);
9✔
2440
        std::vector<GByte> abyDstNoData(nBufDTSize);
18✔
2441
        GDALCopyWords(&m_dfNoDataValue, GDT_Float64, 0, abyDstNoData.data(),
9✔
2442
                      eBufType, 0, 1);
2443
        std::vector<GByte> abyRemappedValue(nBufDTSize);
18✔
2444
        GDALCopyWords(&dfRemappedValue, GDT_Float64, 0, abyRemappedValue.data(),
9✔
2445
                      eBufType, 0, 1);
2446
        for (int iY = 0; iY < nOutYSize; iY++)
18✔
2447
        {
2448
            GSpacing nDstOffset = iY * nLineSpace;
9✔
2449
            for (int iX = 0; iX < nOutXSize; iX++)
28✔
2450
            {
2451
                if (eSrcMaskBandDT == GDT_Byte)
19✔
2452
                {
2453
                    dfMaskVal = oWorkingState.m_abyWrkBufferMask[nSrcIdx];
19✔
2454
                }
2455
                else
2456
                {
2457
                    GDALCopyWords(oWorkingState.m_abyWrkBufferMask.data() +
×
2458
                                      nSrcIdx * nSrcMaskBandDTSize,
×
2459
                                  eSrcMaskBandDT, 0, &dfMaskVal, GDT_Float64, 0,
2460
                                  1);
2461
                }
2462
                void *const pDst =
19✔
2463
                    pabyOut + static_cast<GPtrDiff_t>(nDstOffset);
19✔
2464
                if (!(dfMaskVal > m_dfMaskValueThreshold))
19✔
2465
                {
2466
                    memcpy(pDst, abyDstNoData.data(), nBufDTSize);
9✔
2467
                }
2468
                else
2469
                {
2470
                    const void *const pSrc =
10✔
2471
                        pabyWrkBuffer + nSrcIdx * nSrcBandDTSize;
10✔
2472
                    if (eSrcBandDT == eBufType)
10✔
2473
                    {
2474
                        // coverity[overrun-buffer-arg]
2475
                        memcpy(pDst, pSrc, nBufDTSize);
8✔
2476
                    }
2477
                    else
2478
                    {
2479
                        GDALCopyWords(pSrc, eSrcBandDT, 0, pDst, eBufType, 0,
2✔
2480
                                      1);
2481
                    }
2482
                    if (memcmp(pDst, abyDstNoData.data(), nBufDTSize) == 0)
10✔
2483
                        memcpy(pDst, abyRemappedValue.data(), nBufDTSize);
9✔
2484
                }
2485
                nDstOffset += nPixelSpace;
19✔
2486
                nSrcIdx++;
19✔
2487
            }
2488
        }
2489
    }
2490

2491
    if (psExtraArg->pfnProgress)
13✔
2492
        psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
×
2493

2494
    return CE_None;
13✔
2495
}
2496

2497
/************************************************************************/
2498
/*                             GetMinimum()                             */
2499
/************************************************************************/
2500

2501
double VRTNoDataFromMaskSource::GetMinimum(int /* nXSize */, int /* nYSize */,
×
2502
                                           int *pbSuccess)
2503
{
2504
    *pbSuccess = FALSE;
×
2505
    return 0.0;
×
2506
}
2507

2508
/************************************************************************/
2509
/*                             GetMaximum()                             */
2510
/************************************************************************/
2511

2512
double VRTNoDataFromMaskSource::GetMaximum(int /* nXSize */, int /* nYSize */,
×
2513
                                           int *pbSuccess)
2514
{
2515
    *pbSuccess = FALSE;
×
2516
    return 0.0;
×
2517
}
2518

2519
/************************************************************************/
2520
/*                            GetHistogram()                            */
2521
/************************************************************************/
2522

2523
CPLErr VRTNoDataFromMaskSource::GetHistogram(
×
2524
    int /* nXSize */, int /* nYSize */, double /* dfMin */, double /* dfMax */,
2525
    int /* nBuckets */, GUIntBig * /* panHistogram */,
2526
    int /* bIncludeOutOfRange */, int /* bApproxOK */,
2527
    GDALProgressFunc /* pfnProgress */, void * /* pProgressData */)
2528
{
2529
    return CE_Failure;
×
2530
}
2531

2532
/************************************************************************/
2533
/* ==================================================================== */
2534
/*                          VRTComplexSource                            */
2535
/* ==================================================================== */
2536
/************************************************************************/
2537

2538
/************************************************************************/
2539
/*                          VRTComplexSource()                          */
2540
/************************************************************************/
2541

2542
VRTComplexSource::VRTComplexSource(const VRTComplexSource *poSrcSource,
5✔
2543
                                   double dfXDstRatio, double dfYDstRatio)
5✔
2544
    : VRTSimpleSource(poSrcSource, dfXDstRatio, dfYDstRatio),
2545
      m_nProcessingFlags(poSrcSource->m_nProcessingFlags),
5✔
2546
      m_dfNoDataValue(poSrcSource->m_dfNoDataValue),
5✔
2547
      m_osNoDataValueOri(poSrcSource->m_osNoDataValueOri),
5✔
2548
      m_dfScaleOff(poSrcSource->m_dfScaleOff),
5✔
2549
      m_dfScaleRatio(poSrcSource->m_dfScaleRatio),
5✔
2550
      m_bSrcMinMaxDefined(poSrcSource->m_bSrcMinMaxDefined),
5✔
2551
      m_dfSrcMin(poSrcSource->m_dfSrcMin), m_dfSrcMax(poSrcSource->m_dfSrcMax),
5✔
2552
      m_dfDstMin(poSrcSource->m_dfDstMin), m_dfDstMax(poSrcSource->m_dfDstMax),
5✔
2553
      m_dfExponent(poSrcSource->m_dfExponent),
5✔
2554
      m_nColorTableComponent(poSrcSource->m_nColorTableComponent),
5✔
2555
      m_adfLUTInputs(poSrcSource->m_adfLUTInputs),
5✔
2556
      m_adfLUTOutputs(poSrcSource->m_adfLUTOutputs)
5✔
2557
{
2558
}
5✔
2559

2560
/************************************************************************/
2561
/*                           GetTypeStatic()                            */
2562
/************************************************************************/
2563

2564
const char *VRTComplexSource::GetTypeStatic()
13,519✔
2565
{
2566
    static const char *TYPE = "ComplexSource";
2567
    return TYPE;
13,519✔
2568
}
2569

2570
/************************************************************************/
2571
/*                            GetType()                                 */
2572
/************************************************************************/
2573

2574
const char *VRTComplexSource::GetType() const
3,839✔
2575
{
2576
    return GetTypeStatic();
3,839✔
2577
}
2578

2579
/************************************************************************/
2580
/*                           SetNoDataValue()                           */
2581
/************************************************************************/
2582

2583
void VRTComplexSource::SetNoDataValue(double dfNewNoDataValue)
65✔
2584

2585
{
2586
    if (dfNewNoDataValue == VRT_NODATA_UNSET)
65✔
2587
    {
2588
        m_nProcessingFlags &= ~PROCESSING_FLAG_NODATA;
×
2589
        m_dfNoDataValue = VRT_NODATA_UNSET;
×
2590
        return;
×
2591
    }
2592

2593
    m_nProcessingFlags |= PROCESSING_FLAG_NODATA;
65✔
2594
    m_dfNoDataValue = dfNewNoDataValue;
65✔
2595
}
2596

2597
/************************************************************************/
2598
/*                      GetAdjustedNoDataValue()                        */
2599
/************************************************************************/
2600

2601
double VRTComplexSource::GetAdjustedNoDataValue() const
4,726✔
2602
{
2603
    if ((m_nProcessingFlags & PROCESSING_FLAG_NODATA) != 0)
4,726✔
2604
    {
2605
        auto l_band = GetRasterBand();
36✔
2606
        if (l_band && l_band->GetRasterDataType() == GDT_Float32)
36✔
2607
        {
2608
            return GDALAdjustNoDataCloseToFloatMax(m_dfNoDataValue);
15✔
2609
        }
2610
    }
2611
    return m_dfNoDataValue;
4,711✔
2612
}
2613

2614
/************************************************************************/
2615
/*                           SerializeToXML()                           */
2616
/************************************************************************/
2617

2618
CPLXMLNode *VRTComplexSource::SerializeToXML(const char *pszVRTPath)
65✔
2619

2620
{
2621
    CPLXMLNode *psSrc = VRTSimpleSource::SerializeToXML(pszVRTPath);
65✔
2622

2623
    if (psSrc == nullptr)
65✔
2624
        return nullptr;
×
2625

2626
    CPLFree(psSrc->pszValue);
65✔
2627
    psSrc->pszValue = CPLStrdup(GetTypeStatic());
65✔
2628

2629
    if ((m_nProcessingFlags & PROCESSING_FLAG_USE_MASK_BAND) != 0)
65✔
2630
    {
2631
        CPLSetXMLValue(psSrc, "UseMaskBand", "true");
30✔
2632
    }
2633

2634
    if ((m_nProcessingFlags & PROCESSING_FLAG_NODATA) != 0)
65✔
2635
    {
2636
        if (!m_osNoDataValueOri.empty() && GetRasterBandNoOpen() == nullptr)
25✔
2637
        {
2638
            CPLSetXMLValue(psSrc, "NODATA", m_osNoDataValueOri.c_str());
2✔
2639
        }
2640
        else
2641
        {
2642
            GDALDataType eBandDT = GDT_Unknown;
23✔
2643
            double dfNoDataValue = m_dfNoDataValue;
23✔
2644
            const auto kMaxFloat = std::numeric_limits<float>::max();
23✔
2645
            if (std::fabs(std::fabs(m_dfNoDataValue) - kMaxFloat) <
23✔
2646
                1e-10 * kMaxFloat)
2647
            {
2648
                auto l_band = GetRasterBand();
1✔
2649
                if (l_band)
1✔
2650
                {
2651
                    dfNoDataValue = GetAdjustedNoDataValue();
1✔
2652
                    eBandDT = l_band->GetRasterDataType();
1✔
2653
                }
2654
            }
2655
            CPLSetXMLValue(
23✔
2656
                psSrc, "NODATA",
2657
                VRTSerializeNoData(dfNoDataValue, eBandDT, 18).c_str());
46✔
2658
        }
2659
    }
2660

2661
    if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_LINEAR) != 0)
65✔
2662
    {
2663
        CPLSetXMLValue(psSrc, "ScaleOffset", CPLSPrintf("%g", m_dfScaleOff));
1✔
2664
        CPLSetXMLValue(psSrc, "ScaleRatio", CPLSPrintf("%g", m_dfScaleRatio));
1✔
2665
    }
2666
    else if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_EXPONENTIAL) != 0)
64✔
2667
    {
2668
        CPLSetXMLValue(psSrc, "Exponent", CPLSPrintf("%g", m_dfExponent));
×
2669
        if (m_bSrcMinMaxDefined)
×
2670
        {
2671
            CPLSetXMLValue(psSrc, "SrcMin", CPLSPrintf("%g", m_dfSrcMin));
×
2672
            CPLSetXMLValue(psSrc, "SrcMax", CPLSPrintf("%g", m_dfSrcMax));
×
2673
        }
2674
        CPLSetXMLValue(psSrc, "DstMin", CPLSPrintf("%g", m_dfDstMin));
×
2675
        CPLSetXMLValue(psSrc, "DstMax", CPLSPrintf("%g", m_dfDstMax));
×
2676
    }
2677

2678
    if (!m_adfLUTInputs.empty())
65✔
2679
    {
2680
        // Make sure we print with sufficient precision to address really close
2681
        // entries (#6422).
2682
        CPLString osLUT;
×
2683
        if (m_adfLUTInputs.size() >= 2 &&
×
2684
            CPLString().Printf("%g", m_adfLUTInputs[0]) ==
×
2685
                CPLString().Printf("%g", m_adfLUTInputs[1]))
×
2686
        {
2687
            osLUT = CPLString().Printf("%.17g:%g", m_adfLUTInputs[0],
×
2688
                                       m_adfLUTOutputs[0]);
×
2689
        }
2690
        else
2691
        {
2692
            osLUT = CPLString().Printf("%g:%g", m_adfLUTInputs[0],
×
2693
                                       m_adfLUTOutputs[0]);
×
2694
        }
2695
        for (size_t i = 1; i < m_adfLUTInputs.size(); i++)
×
2696
        {
2697
            if (CPLString().Printf("%g", m_adfLUTInputs[i]) ==
×
2698
                    CPLString().Printf("%g", m_adfLUTInputs[i - 1]) ||
×
2699
                (i + 1 < m_adfLUTInputs.size() &&
×
2700
                 CPLString().Printf("%g", m_adfLUTInputs[i]) ==
×
2701
                     CPLString().Printf("%g", m_adfLUTInputs[i + 1])))
×
2702
            {
2703
                // TODO(schwehr): An explanation of the 18 would be helpful.
2704
                // Can someone distill the issue down to a quick comment?
2705
                // https://trac.osgeo.org/gdal/ticket/6422
2706
                osLUT += CPLString().Printf(",%.17g:%g", m_adfLUTInputs[i],
×
2707
                                            m_adfLUTOutputs[i]);
×
2708
            }
2709
            else
2710
            {
2711
                osLUT += CPLString().Printf(",%g:%g", m_adfLUTInputs[i],
×
2712
                                            m_adfLUTOutputs[i]);
×
2713
            }
2714
        }
2715
        CPLSetXMLValue(psSrc, "LUT", osLUT);
×
2716
    }
2717

2718
    if (m_nColorTableComponent)
65✔
2719
    {
2720
        CPLSetXMLValue(psSrc, "ColorTableComponent",
7✔
2721
                       CPLSPrintf("%d", m_nColorTableComponent));
2722
    }
2723

2724
    return psSrc;
65✔
2725
}
2726

2727
/************************************************************************/
2728
/*                              XMLInit()                               */
2729
/************************************************************************/
2730

2731
CPLErr VRTComplexSource::XMLInit(const CPLXMLNode *psSrc,
217✔
2732
                                 const char *pszVRTPath,
2733
                                 VRTMapSharedResources &oMapSharedSources)
2734

2735
{
2736
    /* -------------------------------------------------------------------- */
2737
    /*      Do base initialization.                                         */
2738
    /* -------------------------------------------------------------------- */
2739
    {
2740
        const CPLErr eErr =
2741
            VRTSimpleSource::XMLInit(psSrc, pszVRTPath, oMapSharedSources);
217✔
2742
        if (eErr != CE_None)
217✔
2743
            return eErr;
×
2744
    }
2745

2746
    /* -------------------------------------------------------------------- */
2747
    /*      Complex parameters.                                             */
2748
    /* -------------------------------------------------------------------- */
2749
    const char *pszScaleOffset = CPLGetXMLValue(psSrc, "ScaleOffset", nullptr);
217✔
2750
    const char *pszScaleRatio = CPLGetXMLValue(psSrc, "ScaleRatio", nullptr);
217✔
2751
    if (pszScaleOffset || pszScaleRatio)
217✔
2752
    {
2753
        m_nProcessingFlags |= PROCESSING_FLAG_SCALING_LINEAR;
18✔
2754
        if (pszScaleOffset)
18✔
2755
            m_dfScaleOff = CPLAtof(pszScaleOffset);
18✔
2756
        if (pszScaleRatio)
18✔
2757
            m_dfScaleRatio = CPLAtof(pszScaleRatio);
15✔
2758
    }
2759
    else if (CPLGetXMLValue(psSrc, "Exponent", nullptr) != nullptr &&
199✔
2760
             CPLGetXMLValue(psSrc, "DstMin", nullptr) != nullptr &&
200✔
2761
             CPLGetXMLValue(psSrc, "DstMax", nullptr) != nullptr)
1✔
2762
    {
2763
        m_nProcessingFlags |= PROCESSING_FLAG_SCALING_EXPONENTIAL;
1✔
2764
        m_dfExponent = CPLAtof(CPLGetXMLValue(psSrc, "Exponent", "1.0"));
1✔
2765

2766
        const char *pszSrcMin = CPLGetXMLValue(psSrc, "SrcMin", nullptr);
1✔
2767
        const char *pszSrcMax = CPLGetXMLValue(psSrc, "SrcMax", nullptr);
1✔
2768
        if (pszSrcMin && pszSrcMax)
1✔
2769
        {
2770
            m_dfSrcMin = CPLAtof(pszSrcMin);
×
2771
            m_dfSrcMax = CPLAtof(pszSrcMax);
×
2772
            m_bSrcMinMaxDefined = true;
×
2773
        }
2774

2775
        m_dfDstMin = CPLAtof(CPLGetXMLValue(psSrc, "DstMin", "0.0"));
1✔
2776
        m_dfDstMax = CPLAtof(CPLGetXMLValue(psSrc, "DstMax", "0.0"));
1✔
2777
    }
2778

2779
    if (const char *pszNODATA = CPLGetXMLValue(psSrc, "NODATA", nullptr))
217✔
2780
    {
2781
        m_nProcessingFlags |= PROCESSING_FLAG_NODATA;
50✔
2782
        m_osNoDataValueOri = pszNODATA;
50✔
2783
        m_dfNoDataValue = CPLAtofM(m_osNoDataValueOri.c_str());
50✔
2784
    }
2785

2786
    const char *pszUseMaskBand = CPLGetXMLValue(psSrc, "UseMaskBand", nullptr);
217✔
2787
    if (pszUseMaskBand && CPLTestBool(pszUseMaskBand))
217✔
2788
    {
2789
        m_nProcessingFlags |= PROCESSING_FLAG_USE_MASK_BAND;
38✔
2790
    }
2791

2792
    const char *pszLUT = CPLGetXMLValue(psSrc, "LUT", nullptr);
217✔
2793
    if (pszLUT)
217✔
2794
    {
2795
        const CPLStringList aosValues(
2796
            CSLTokenizeString2(pszLUT, ",:", CSLT_ALLOWEMPTYTOKENS));
59✔
2797

2798
        const int nLUTItemCount = aosValues.size() / 2;
59✔
2799
        try
2800
        {
2801
            m_adfLUTInputs.resize(nLUTItemCount);
59✔
2802
            m_adfLUTOutputs.resize(nLUTItemCount);
59✔
2803
        }
2804
        catch (const std::bad_alloc &e)
×
2805
        {
2806
            CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
×
2807
            m_adfLUTInputs.clear();
×
2808
            m_adfLUTOutputs.clear();
×
2809
            return CE_Failure;
×
2810
        }
2811

2812
        for (int nIndex = 0; nIndex < nLUTItemCount; nIndex++)
544✔
2813
        {
2814
            m_adfLUTInputs[nIndex] = CPLAtof(aosValues[nIndex * 2]);
485✔
2815
            m_adfLUTOutputs[nIndex] = CPLAtof(aosValues[nIndex * 2 + 1]);
485✔
2816

2817
            // Enforce the requirement that the LUT input array is
2818
            // monotonically non-decreasing.
2819
            if (std::isnan(m_adfLUTInputs[nIndex]) && nIndex != 0)
485✔
2820
            {
2821
                CPLError(CE_Failure, CPLE_AppDefined,
×
2822
                         "A Not-A-Number (NaN) source value should be the "
2823
                         "first one of the LUT.");
2824
                m_adfLUTInputs.clear();
×
2825
                m_adfLUTOutputs.clear();
×
2826
                return CE_Failure;
×
2827
            }
2828
            else if (nIndex > 0 &&
911✔
2829
                     m_adfLUTInputs[nIndex] < m_adfLUTInputs[nIndex - 1])
426✔
2830
            {
2831
                CPLError(CE_Failure, CPLE_AppDefined,
×
2832
                         "Source values of the LUT are not listed in a "
2833
                         "monotonically non-decreasing order");
2834
                m_adfLUTInputs.clear();
×
2835
                m_adfLUTOutputs.clear();
×
2836
                return CE_Failure;
×
2837
            }
2838
        }
2839

2840
        m_nProcessingFlags |= PROCESSING_FLAG_LUT;
59✔
2841
    }
2842

2843
    const char *pszColorTableComponent =
2844
        CPLGetXMLValue(psSrc, "ColorTableComponent", nullptr);
217✔
2845
    if (pszColorTableComponent)
217✔
2846
    {
2847
        m_nColorTableComponent = atoi(pszColorTableComponent);
15✔
2848
        m_nProcessingFlags |= PROCESSING_FLAG_COLOR_TABLE_EXPANSION;
15✔
2849
    }
2850

2851
    return CE_None;
217✔
2852
}
2853

2854
/************************************************************************/
2855
/*                              LookupValue()                           */
2856
/************************************************************************/
2857

2858
double VRTComplexSource::LookupValue(double dfInput)
583,336✔
2859
{
2860
    auto beginIter = m_adfLUTInputs.begin();
583,336✔
2861
    auto endIter = m_adfLUTInputs.end();
583,336✔
2862
    size_t offset = 0;
583,336✔
2863
    if (std::isnan(m_adfLUTInputs[0]))
583,336✔
2864
    {
2865
        if (std::isnan(dfInput) || m_adfLUTInputs.size() == 1)
6✔
2866
            return m_adfLUTOutputs[0];
1✔
2867
        ++beginIter;
5✔
2868
        offset = 1;
5✔
2869
    }
2870

2871
    // Find the index of the first element in the LUT input array that
2872
    // is not smaller than the input value.
2873
    const size_t i =
2874
        offset +
2875
        std::distance(beginIter, std::lower_bound(beginIter, endIter, dfInput));
583,335✔
2876

2877
    if (i == offset)
583,335✔
2878
        return m_adfLUTOutputs[offset];
129,039✔
2879

2880
    // If the index is beyond the end of the LUT input array, the input
2881
    // value is larger than all the values in the array.
2882
    if (i == m_adfLUTInputs.size())
454,296✔
2883
        return m_adfLUTOutputs.back();
8✔
2884

2885
    if (m_adfLUTInputs[i] == dfInput)
454,288✔
2886
        return m_adfLUTOutputs[i];
179,297✔
2887

2888
    // Otherwise, interpolate.
2889
    return m_adfLUTOutputs[i - 1] +
274,991✔
2890
           (dfInput - m_adfLUTInputs[i - 1]) *
274,991✔
2891
               ((m_adfLUTOutputs[i] - m_adfLUTOutputs[i - 1]) /
274,991✔
2892
                (m_adfLUTInputs[i] - m_adfLUTInputs[i - 1]));
274,991✔
2893
}
2894

2895
/************************************************************************/
2896
/*                         SetLinearScaling()                           */
2897
/************************************************************************/
2898

2899
void VRTComplexSource::SetLinearScaling(double dfOffset, double dfScale)
94✔
2900
{
2901
    m_nProcessingFlags &= ~PROCESSING_FLAG_SCALING_EXPONENTIAL;
94✔
2902
    m_nProcessingFlags |= PROCESSING_FLAG_SCALING_LINEAR;
94✔
2903
    m_dfScaleOff = dfOffset;
94✔
2904
    m_dfScaleRatio = dfScale;
94✔
2905
}
94✔
2906

2907
/************************************************************************/
2908
/*                         SetPowerScaling()                           */
2909
/************************************************************************/
2910

2911
void VRTComplexSource::SetPowerScaling(double dfExponentIn, double dfSrcMinIn,
8✔
2912
                                       double dfSrcMaxIn, double dfDstMinIn,
2913
                                       double dfDstMaxIn)
2914
{
2915
    m_nProcessingFlags &= ~PROCESSING_FLAG_SCALING_LINEAR;
8✔
2916
    m_nProcessingFlags |= PROCESSING_FLAG_SCALING_EXPONENTIAL;
8✔
2917
    m_dfExponent = dfExponentIn;
8✔
2918
    m_dfSrcMin = dfSrcMinIn;
8✔
2919
    m_dfSrcMax = dfSrcMaxIn;
8✔
2920
    m_dfDstMin = dfDstMinIn;
8✔
2921
    m_dfDstMax = dfDstMaxIn;
8✔
2922
    m_bSrcMinMaxDefined = true;
8✔
2923
}
8✔
2924

2925
/************************************************************************/
2926
/*                    SetColorTableComponent()                          */
2927
/************************************************************************/
2928

2929
void VRTComplexSource::SetColorTableComponent(int nComponent)
186✔
2930
{
2931
    m_nProcessingFlags |= PROCESSING_FLAG_COLOR_TABLE_EXPANSION;
186✔
2932
    m_nColorTableComponent = nComponent;
186✔
2933
}
186✔
2934

2935
/************************************************************************/
2936
/*                              RasterIO()                              */
2937
/************************************************************************/
2938

2939
CPLErr VRTComplexSource::RasterIO(GDALDataType eVRTBandDataType, int nXOff,
15,703✔
2940
                                  int nYOff, int nXSize, int nYSize,
2941
                                  void *pData, int nBufXSize, int nBufYSize,
2942
                                  GDALDataType eBufType, GSpacing nPixelSpace,
2943
                                  GSpacing nLineSpace,
2944
                                  GDALRasterIOExtraArg *psExtraArgIn,
2945
                                  WorkingState &oWorkingState)
2946

2947
{
2948
    GDALRasterIOExtraArg sExtraArg;
2949
    INIT_RASTERIO_EXTRA_ARG(sExtraArg);
15,703✔
2950
    GDALRasterIOExtraArg *psExtraArg = &sExtraArg;
15,703✔
2951

2952
    double dfXOff = nXOff;
15,703✔
2953
    double dfYOff = nYOff;
15,703✔
2954
    double dfXSize = nXSize;
15,703✔
2955
    double dfYSize = nYSize;
15,703✔
2956
    if (psExtraArgIn != nullptr && psExtraArgIn->bFloatingPointWindowValidity)
15,703✔
2957
    {
2958
        dfXOff = psExtraArgIn->dfXOff;
87✔
2959
        dfYOff = psExtraArgIn->dfYOff;
87✔
2960
        dfXSize = psExtraArgIn->dfXSize;
87✔
2961
        dfYSize = psExtraArgIn->dfYSize;
87✔
2962
    }
2963

2964
    // The window we will actually request from the source raster band.
2965
    double dfReqXOff = 0.0;
15,703✔
2966
    double dfReqYOff = 0.0;
15,703✔
2967
    double dfReqXSize = 0.0;
15,703✔
2968
    double dfReqYSize = 0.0;
15,703✔
2969
    int nReqXOff = 0;
15,703✔
2970
    int nReqYOff = 0;
15,703✔
2971
    int nReqXSize = 0;
15,703✔
2972
    int nReqYSize = 0;
15,703✔
2973

2974
    // The window we will actual set _within_ the pData buffer.
2975
    int nOutXOff = 0;
15,703✔
2976
    int nOutYOff = 0;
15,703✔
2977
    int nOutXSize = 0;
15,703✔
2978
    int nOutYSize = 0;
15,703✔
2979

2980
    bool bError = false;
15,703✔
2981
    if (!GetSrcDstWindow(dfXOff, dfYOff, dfXSize, dfYSize, nBufXSize, nBufYSize,
15,703✔
2982
                         &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
2983
                         &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
2984
                         &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize, bError))
2985
    {
2986
        return bError ? CE_Failure : CE_None;
10,946✔
2987
    }
2988
#if DEBUG_VERBOSE
2989
    CPLDebug("VRT",
2990
             "nXOff=%d, nYOff=%d, nXSize=%d, nYSize=%d, nBufXSize=%d, "
2991
             "nBufYSize=%d,\n"
2992
             "dfReqXOff=%g, dfReqYOff=%g, dfReqXSize=%g, dfReqYSize=%g,\n"
2993
             "nReqXOff=%d, nReqYOff=%d, nReqXSize=%d, nReqYSize=%d,\n"
2994
             "nOutXOff=%d, nOutYOff=%d, nOutXSize=%d, nOutYSize=%d",
2995
             nXOff, nYOff, nXSize, nYSize, nBufXSize, nBufYSize, dfReqXOff,
2996
             dfReqYOff, dfReqXSize, dfReqYSize, nReqXOff, nReqYOff, nReqXSize,
2997
             nReqYSize, nOutXOff, nOutYOff, nOutXSize, nOutYSize);
2998
#endif
2999

3000
    auto poSourceBand = GetRasterBand();
4,757✔
3001
    if (!poSourceBand)
4,757✔
3002
        return CE_Failure;
×
3003

3004
    if (!m_osResampling.empty())
4,757✔
3005
    {
3006
        psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
28✔
3007
    }
3008
    else if (psExtraArgIn != nullptr)
4,729✔
3009
    {
3010
        psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
4,729✔
3011
    }
3012
    psExtraArg->bFloatingPointWindowValidity = TRUE;
4,757✔
3013
    psExtraArg->dfXOff = dfReqXOff;
4,757✔
3014
    psExtraArg->dfYOff = dfReqYOff;
4,757✔
3015
    psExtraArg->dfXSize = dfReqXSize;
4,757✔
3016
    psExtraArg->dfYSize = dfReqYSize;
4,757✔
3017
    if (psExtraArgIn)
4,757✔
3018
    {
3019
        psExtraArg->pfnProgress = psExtraArgIn->pfnProgress;
4,757✔
3020
        psExtraArg->pProgressData = psExtraArgIn->pProgressData;
4,757✔
3021
    }
3022

3023
    GByte *const pabyOut = static_cast<GByte *>(pData) +
4,757✔
3024
                           nPixelSpace * nOutXOff +
4,757✔
3025
                           static_cast<GPtrDiff_t>(nLineSpace) * nOutYOff;
4,757✔
3026
    const auto eSourceType = poSourceBand->GetRasterDataType();
4,757✔
3027
    if (m_nProcessingFlags == PROCESSING_FLAG_NODATA)
4,757✔
3028
    {
3029
        // Optimization if doing only nodata processing
3030
        if (eSourceType == GDT_Byte)
65✔
3031
        {
3032
            if (!GDALIsValueInRange<GByte>(m_dfNoDataValue))
36✔
3033
            {
3034
                return VRTSimpleSource::RasterIO(
1✔
3035
                    eVRTBandDataType, nXOff, nYOff, nXSize, nYSize, pData,
3036
                    nBufXSize, nBufYSize, eBufType, nPixelSpace, nLineSpace,
3037
                    psExtraArgIn, oWorkingState);
1✔
3038
            }
3039

3040
            return RasterIOProcessNoData<GByte, GDT_Byte>(
35✔
3041
                poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3042
                nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3043
                nLineSpace, psExtraArg, oWorkingState);
35✔
3044
        }
3045
        else if (eSourceType == GDT_Int16)
29✔
3046
        {
3047
            if (!GDALIsValueInRange<int16_t>(m_dfNoDataValue))
2✔
3048
            {
3049
                return VRTSimpleSource::RasterIO(
1✔
3050
                    eVRTBandDataType, nXOff, nYOff, nXSize, nYSize, pData,
3051
                    nBufXSize, nBufYSize, eBufType, nPixelSpace, nLineSpace,
3052
                    psExtraArgIn, oWorkingState);
1✔
3053
            }
3054

3055
            return RasterIOProcessNoData<int16_t, GDT_Int16>(
1✔
3056
                poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3057
                nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3058
                nLineSpace, psExtraArg, oWorkingState);
1✔
3059
        }
3060
        else if (eSourceType == GDT_UInt16)
27✔
3061
        {
3062
            if (!GDALIsValueInRange<uint16_t>(m_dfNoDataValue))
8✔
3063
            {
3064
                return VRTSimpleSource::RasterIO(
1✔
3065
                    eVRTBandDataType, nXOff, nYOff, nXSize, nYSize, pData,
3066
                    nBufXSize, nBufYSize, eBufType, nPixelSpace, nLineSpace,
3067
                    psExtraArgIn, oWorkingState);
1✔
3068
            }
3069

3070
            return RasterIOProcessNoData<uint16_t, GDT_UInt16>(
7✔
3071
                poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3072
                nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3073
                nLineSpace, psExtraArg, oWorkingState);
7✔
3074
        }
3075
    }
3076

3077
    const bool bIsComplex = CPL_TO_BOOL(GDALDataTypeIsComplex(eBufType));
4,711✔
3078
    CPLErr eErr;
3079
    // For Int32, float32 isn't sufficiently precise as working data type
3080
    if (eVRTBandDataType == GDT_CInt32 || eVRTBandDataType == GDT_CFloat64 ||
4,711✔
3081
        eVRTBandDataType == GDT_Int32 || eVRTBandDataType == GDT_UInt32 ||
4,709✔
3082
        eVRTBandDataType == GDT_Int64 || eVRTBandDataType == GDT_UInt64 ||
4,707✔
3083
        eVRTBandDataType == GDT_Float64 || eSourceType == GDT_Int32 ||
4,688✔
3084
        eSourceType == GDT_UInt32 || eSourceType == GDT_Int64 ||
4,687✔
3085
        eSourceType == GDT_UInt64)
3086
    {
3087
        eErr = RasterIOInternal<double>(
25✔
3088
            poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3089
            nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3090
            nLineSpace, psExtraArg, bIsComplex ? GDT_CFloat64 : GDT_Float64,
3091
            oWorkingState);
3092
    }
3093
    else
3094
    {
3095
        eErr = RasterIOInternal<float>(
4,686✔
3096
            poSourceBand, eVRTBandDataType, nReqXOff, nReqYOff, nReqXSize,
3097
            nReqYSize, pabyOut, nOutXSize, nOutYSize, eBufType, nPixelSpace,
3098
            nLineSpace, psExtraArg, bIsComplex ? GDT_CFloat32 : GDT_Float32,
3099
            oWorkingState);
3100
    }
3101

3102
    if (psExtraArg->pfnProgress)
4,711✔
3103
        psExtraArg->pfnProgress(1.0, "", psExtraArg->pProgressData);
43✔
3104

3105
    return eErr;
4,711✔
3106
}
3107

3108
/************************************************************************/
3109
/*                              hasZeroByte()                           */
3110
/************************************************************************/
3111

3112
CPL_NOSANITIZE_UNSIGNED_INT_OVERFLOW
3113
static inline bool hasZeroByte(uint32_t v)
3,351,380✔
3114
{
3115
    // Cf https://graphics.stanford.edu/~seander/bithacks.html#ZeroInWord
3116
    return (((v)-0x01010101U) & ~(v)&0x80808080U) != 0;
3,351,380✔
3117
}
3118

3119
/************************************************************************/
3120
/*                                CopyWord()                            */
3121
/************************************************************************/
3122

3123
template <class SrcType>
3124
static void CopyWord(const SrcType *pSrcVal, GDALDataType eSrcType,
6,418,544✔
3125
                     void *pDstVal, GDALDataType eDstType)
3126
{
3127
    switch (eDstType)
6,418,544✔
3128
    {
3129
        case GDT_Byte:
2,144,660✔
3130
            GDALCopyWord(*pSrcVal, *static_cast<uint8_t *>(pDstVal));
2,144,660✔
3131
            break;
2,144,660✔
3132
        case GDT_Int8:
80✔
3133
            GDALCopyWord(*pSrcVal, *static_cast<int8_t *>(pDstVal));
80✔
3134
            break;
80✔
3135
        case GDT_UInt16:
1,348,720✔
3136
            GDALCopyWord(*pSrcVal, *static_cast<uint16_t *>(pDstVal));
1,348,720✔
3137
            break;
1,348,720✔
3138
        case GDT_Int16:
2,920,180✔
3139
            GDALCopyWord(*pSrcVal, *static_cast<int16_t *>(pDstVal));
2,920,180✔
3140
            break;
2,920,180✔
3141
        case GDT_UInt32:
20✔
3142
            GDALCopyWord(*pSrcVal, *static_cast<uint32_t *>(pDstVal));
20✔
3143
            break;
20✔
3144
        case GDT_Int32:
638✔
3145
            GDALCopyWord(*pSrcVal, *static_cast<int32_t *>(pDstVal));
638✔
3146
            break;
638✔
3147
        case GDT_UInt64:
20✔
3148
            GDALCopyWord(*pSrcVal, *static_cast<uint64_t *>(pDstVal));
20✔
3149
            break;
20✔
3150
        case GDT_Int64:
40✔
3151
            GDALCopyWord(*pSrcVal, *static_cast<int64_t *>(pDstVal));
40✔
3152
            break;
40✔
3153
        case GDT_Float32:
1,650✔
3154
            GDALCopyWord(*pSrcVal, *static_cast<float *>(pDstVal));
1,650✔
3155
            break;
1,650✔
3156
        case GDT_Float64:
1,525✔
3157
            GDALCopyWord(*pSrcVal, *static_cast<double *>(pDstVal));
1,525✔
3158
            break;
1,525✔
3159
        default:
1,016✔
3160
            GDALCopyWords(pSrcVal, eSrcType, 0, pDstVal, eDstType, 0, 1);
1,016✔
3161
            break;
1,016✔
3162
    }
3163
}
6,418,544✔
3164

3165
/************************************************************************/
3166
/*                       RasterIOProcessNoData()                        */
3167
/************************************************************************/
3168

3169
// This method is an optimization of the generic RasterIOInternal()
3170
// that deals with a VRTComplexSource with only a NODATA value in it, and
3171
// no other processing flags.
3172

3173
// nReqXOff, nReqYOff, nReqXSize, nReqYSize are expressed in source band
3174
// referential.
3175
template <class SourceDT, GDALDataType eSourceType>
3176
CPLErr VRTComplexSource::RasterIOProcessNoData(
43✔
3177
    GDALRasterBand *poSourceBand, GDALDataType eVRTBandDataType, int nReqXOff,
3178
    int nReqYOff, int nReqXSize, int nReqYSize, void *pData, int nOutXSize,
3179
    int nOutYSize, GDALDataType eBufType, GSpacing nPixelSpace,
3180
    GSpacing nLineSpace, GDALRasterIOExtraArg *psExtraArg,
3181
    WorkingState &oWorkingState)
3182
{
3183
    CPLAssert(m_nProcessingFlags == PROCESSING_FLAG_NODATA);
43✔
3184
    CPLAssert(GDALIsValueInRange<SourceDT>(m_dfNoDataValue));
43✔
3185

3186
    /* -------------------------------------------------------------------- */
3187
    /*      Read into a temporary buffer.                                   */
3188
    /* -------------------------------------------------------------------- */
3189
    try
3190
    {
3191
        // Cannot overflow since pData should at least have that number of
3192
        // elements
3193
        const size_t nPixelCount = static_cast<size_t>(nOutXSize) * nOutYSize;
43✔
3194
        if (nPixelCount >
43✔
3195
            static_cast<size_t>(std::numeric_limits<ptrdiff_t>::max()) /
43✔
3196
                sizeof(SourceDT))
3197
        {
3198
            CPLError(CE_Failure, CPLE_OutOfMemory,
×
3199
                     "Too large temporary buffer");
3200
            return CE_Failure;
×
3201
        }
3202
        oWorkingState.m_abyWrkBuffer.resize(sizeof(SourceDT) * nPixelCount);
43✔
3203
    }
3204
    catch (const std::bad_alloc &e)
×
3205
    {
3206
        CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
×
3207
        return CE_Failure;
×
3208
    }
3209
    const auto paSrcData =
3210
        reinterpret_cast<const SourceDT *>(oWorkingState.m_abyWrkBuffer.data());
43✔
3211

3212
    const GDALRIOResampleAlg eResampleAlgBack = psExtraArg->eResampleAlg;
42✔
3213
    if (!m_osResampling.empty())
42✔
3214
    {
3215
        psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
×
3216
    }
3217

3218
    const CPLErr eErr = poSourceBand->RasterIO(
43✔
3219
        GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
3220
        oWorkingState.m_abyWrkBuffer.data(), nOutXSize, nOutYSize, eSourceType,
42✔
3221
        sizeof(SourceDT), sizeof(SourceDT) * static_cast<GSpacing>(nOutXSize),
8✔
3222
        psExtraArg);
3223
    if (!m_osResampling.empty())
43✔
3224
        psExtraArg->eResampleAlg = eResampleAlgBack;
×
3225

3226
    if (eErr != CE_None)
43✔
3227
    {
3228
        return eErr;
×
3229
    }
3230

3231
    const auto nNoDataValue = static_cast<SourceDT>(m_dfNoDataValue);
43✔
3232
    size_t idxBuffer = 0;
43✔
3233
    if (eSourceType == eBufType &&
74✔
3234
        !GDALDataTypeIsConversionLossy(eSourceType, eVRTBandDataType))
31✔
3235
    {
3236
        // Most optimized case: the output type is the same as the source type,
3237
        // and conversion from the source type to the VRT band data type is
3238
        // not lossy
3239
        for (int iY = 0; iY < nOutYSize; iY++)
14✔
3240
        {
3241
            GByte *pDstLocation = static_cast<GByte *>(pData) +
19,601✔
3242
                                  static_cast<GPtrDiff_t>(nLineSpace) * iY;
19,601✔
3243

3244
            int iX = 0;
19,601✔
3245
            if (sizeof(SourceDT) == 1 && nPixelSpace == 1)
19,589✔
3246
            {
3247
                // Optimization to detect more quickly if source pixels are
3248
                // at nodata.
3249
                const GByte byNoDataValue = static_cast<GByte>(nNoDataValue);
19,598✔
3250
                const uint32_t wordNoData =
19,598✔
3251
                    (static_cast<uint32_t>(byNoDataValue) << 24) |
19,598✔
3252
                    (byNoDataValue << 16) | (byNoDataValue << 8) |
19,598✔
3253
                    byNoDataValue;
19,598✔
3254

3255
                // Warning: hasZeroByte() assumes WORD_SIZE = 4
3256
                constexpr int WORD_SIZE = 4;
19,598✔
3257
                for (; iX < nOutXSize - (WORD_SIZE - 1); iX += WORD_SIZE)
3,393,620✔
3258
                {
3259
                    uint32_t v;
3260
                    static_assert(sizeof(v) == WORD_SIZE,
3261
                                  "sizeof(v) == WORD_SIZE");
3262
                    memcpy(&v, paSrcData + idxBuffer, sizeof(v));
3,579,280✔
3263
                    // Cf https://graphics.stanford.edu/~seander/bithacks.html#ValueInWord
3264
                    if (!hasZeroByte(v ^ wordNoData))
3,579,280✔
3265
                    {
3266
                        // No bytes are at nodata
3267
                        memcpy(pDstLocation, &v, WORD_SIZE);
3,398,020✔
3268
                        idxBuffer += WORD_SIZE;
3,398,020✔
3269
                        pDstLocation += WORD_SIZE;
3,398,020✔
3270
                    }
UNCOV
3271
                    else if (v == wordNoData)
×
3272
                    {
3273
                        // All bytes are at nodata
3274
                        idxBuffer += WORD_SIZE;
131,075✔
3275
                        pDstLocation += WORD_SIZE;
131,075✔
3276
                    }
3277
                    else
3278
                    {
3279
                        // There are both bytes at nodata and valid bytes
3280
                        for (int k = 0; k < WORD_SIZE; ++k)
×
3281
                        {
3282
                            if (paSrcData[idxBuffer] != nNoDataValue)
48✔
3283
                            {
3284
                                memcpy(pDstLocation, &paSrcData[idxBuffer],
36✔
3285
                                       sizeof(SourceDT));
3286
                            }
3287
                            idxBuffer++;
48✔
3288
                            pDstLocation += nPixelSpace;
48✔
3289
                        }
3290
                    }
3291
                }
3292
            }
3293

3294
            for (; iX < nOutXSize;
72✔
3295
                 iX++, pDstLocation += nPixelSpace, idxBuffer++)
3,709✔
3296
            {
3297
                if (paSrcData[idxBuffer] != nNoDataValue)
3,709✔
3298
                {
3299
                    memcpy(pDstLocation, &paSrcData[idxBuffer],
3,679✔
3300
                           sizeof(SourceDT));
3301
                }
3302
            }
3303
        }
3304
    }
3305
    else if (!GDALDataTypeIsConversionLossy(eSourceType, eVRTBandDataType))
13✔
3306
    {
3307
        // Conversion from the source type to the VRT band data type is
3308
        // not lossy, so we can directly convert from the source type to
3309
        // the the output type
3310
        for (int iY = 0; iY < nOutYSize; iY++)
106✔
3311
        {
3312
            GByte *pDstLocation = static_cast<GByte *>(pData) +
94✔
3313
                                  static_cast<GPtrDiff_t>(nLineSpace) * iY;
94✔
3314

3315
            for (int iX = 0; iX < nOutXSize;
1,556✔
3316
                 iX++, pDstLocation += nPixelSpace, idxBuffer++)
1,462✔
3317
            {
3318
                if (paSrcData[idxBuffer] != nNoDataValue)
1,462✔
3319
                {
3320
                    CopyWord(&paSrcData[idxBuffer], eSourceType, pDstLocation,
658✔
3321
                             eBufType);
3322
                }
3323
            }
3324
        }
3325
    }
3326
    else
3327
    {
3328
        GByte abyTemp[2 * sizeof(double)];
3329
        for (int iY = 0; iY < nOutYSize; iY++)
7✔
3330
        {
3331
            GByte *pDstLocation = static_cast<GByte *>(pData) +
6✔
3332
                                  static_cast<GPtrDiff_t>(nLineSpace) * iY;
6✔
3333

3334
            for (int iX = 0; iX < nOutXSize;
36✔
3335
                 iX++, pDstLocation += nPixelSpace, idxBuffer++)
30✔
3336
            {
3337
                if (paSrcData[idxBuffer] != nNoDataValue)
30✔
3338
                {
3339
                    // Convert first to the VRTRasterBand data type
3340
                    // to get its clamping, before outputting to buffer data type
3341
                    CopyWord(&paSrcData[idxBuffer], eSourceType, abyTemp,
20✔
3342
                             eVRTBandDataType);
3343
                    GDALCopyWords(abyTemp, eVRTBandDataType, 0, pDstLocation,
20✔
3344
                                  eBufType, 0, 1);
3345
                }
3346
            }
3347
        }
3348
    }
3349

3350
    return CE_None;
8✔
3351
}
3352

3353
/************************************************************************/
3354
/*                          RasterIOInternal()                          */
3355
/************************************************************************/
3356

3357
// nReqXOff, nReqYOff, nReqXSize, nReqYSize are expressed in source band
3358
// referential.
3359
template <class WorkingDT>
3360
CPLErr VRTComplexSource::RasterIOInternal(
4,725✔
3361
    GDALRasterBand *poSourceBand, GDALDataType eVRTBandDataType, int nReqXOff,
3362
    int nReqYOff, int nReqXSize, int nReqYSize, void *pData, int nOutXSize,
3363
    int nOutYSize, GDALDataType eBufType, GSpacing nPixelSpace,
3364
    GSpacing nLineSpace, GDALRasterIOExtraArg *psExtraArg,
3365
    GDALDataType eWrkDataType, WorkingState &oWorkingState)
3366
{
3367
    const GDALColorTable *poColorTable = nullptr;
4,725✔
3368
    const bool bIsComplex = CPL_TO_BOOL(GDALDataTypeIsComplex(eBufType));
4,725✔
3369
    const int nWordSize = GDALGetDataTypeSizeBytes(eWrkDataType);
4,725✔
3370
    assert(nWordSize != 0);
4,725✔
3371

3372
    // If no explicit <NODATA> is set, but UseMaskBand is set, and the band
3373
    // has a nodata value, then use it as if it was set as <NODATA>
3374
    int bNoDataSet = (m_nProcessingFlags & PROCESSING_FLAG_NODATA) != 0;
4,725✔
3375
    double dfNoDataValue = GetAdjustedNoDataValue();
4,725✔
3376

3377
    if ((m_nProcessingFlags & PROCESSING_FLAG_USE_MASK_BAND) != 0 &&
4,828✔
3378
        poSourceBand->GetMaskFlags() == GMF_NODATA)
103✔
3379
    {
3380
        dfNoDataValue = poSourceBand->GetNoDataValue(&bNoDataSet);
×
3381
    }
3382

3383
    const bool bNoDataSetIsNan = bNoDataSet && std::isnan(dfNoDataValue);
4,725✔
3384
    const bool bNoDataSetAndNotNan =
4,725✔
3385
        bNoDataSet && !std::isnan(dfNoDataValue) &&
4,752✔
3386
        GDALIsValueInRange<WorkingDT>(dfNoDataValue);
27✔
3387
    const auto fWorkingDataTypeNoData = static_cast<WorkingDT>(dfNoDataValue);
4,725✔
3388

3389
    const GByte *pabyMask = nullptr;
4,725✔
3390
    const WorkingDT *pafData = nullptr;
4,725✔
3391
    if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_LINEAR) != 0 &&
4,725✔
3392
        m_dfScaleRatio == 0 && bNoDataSet == FALSE &&
4,429✔
3393
        (m_nProcessingFlags & PROCESSING_FLAG_USE_MASK_BAND) == 0)
4,040✔
3394
    {
3395
        /* ------------------------------------------------------------------ */
3396
        /*      Optimization when writing a constant value */
3397
        /*      (used by the -addalpha option of gdalbuildvrt) */
3398
        /* ------------------------------------------------------------------ */
3399
        // Already set to NULL when defined.
3400
        // pafData = NULL;
3401
    }
3402
    else
3403
    {
3404
        /* ---------------------------------------------------------------- */
3405
        /*      Read into a temporary buffer.                               */
3406
        /* ---------------------------------------------------------------- */
3407
        const size_t nPixelCount = static_cast<size_t>(nOutXSize) * nOutYSize;
685✔
3408
        try
3409
        {
3410
            // Cannot overflow since pData should at least have that number of
3411
            // elements
3412
            if (nPixelCount >
685✔
3413
                static_cast<size_t>(std::numeric_limits<ptrdiff_t>::max()) /
685✔
3414
                    static_cast<size_t>(nWordSize))
685✔
3415
            {
3416
                CPLError(CE_Failure, CPLE_OutOfMemory,
×
3417
                         "Too large temporary buffer");
3418
                return CE_Failure;
×
3419
            }
3420
            oWorkingState.m_abyWrkBuffer.resize(nWordSize * nPixelCount);
685✔
3421
        }
3422
        catch (const std::bad_alloc &e)
×
3423
        {
3424
            CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
×
3425
            return CE_Failure;
×
3426
        }
3427
        pafData = reinterpret_cast<const WorkingDT *>(
3428
            oWorkingState.m_abyWrkBuffer.data());
685✔
3429

3430
        const GDALRIOResampleAlg eResampleAlgBack = psExtraArg->eResampleAlg;
685✔
3431
        if (!m_osResampling.empty())
685✔
3432
        {
3433
            psExtraArg->eResampleAlg =
28✔
3434
                GDALRasterIOGetResampleAlg(m_osResampling);
28✔
3435
        }
3436

3437
        const CPLErr eErr = poSourceBand->RasterIO(
685✔
3438
            GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
3439
            oWorkingState.m_abyWrkBuffer.data(), nOutXSize, nOutYSize,
685✔
3440
            eWrkDataType, nWordSize,
3441
            nWordSize * static_cast<GSpacing>(nOutXSize), psExtraArg);
685✔
3442
        if (!m_osResampling.empty())
685✔
3443
            psExtraArg->eResampleAlg = eResampleAlgBack;
28✔
3444

3445
        if (eErr != CE_None)
685✔
3446
        {
3447
            return eErr;
×
3448
        }
3449

3450
        // Allocate and read mask band if needed
3451
        if (!bNoDataSet &&
2,020✔
3452
            (m_nProcessingFlags & PROCESSING_FLAG_USE_MASK_BAND) != 0 &&
788✔
3453
            (poSourceBand->GetMaskFlags() != GMF_ALL_VALID ||
103✔
3454
             poSourceBand->GetColorInterpretation() == GCI_AlphaBand ||
51✔
3455
             GetMaskBandMainBand() != nullptr))
31✔
3456
        {
3457
            try
3458
            {
3459
                oWorkingState.m_abyWrkBufferMask.resize(nPixelCount);
103✔
3460
            }
3461
            catch (const std::exception &)
×
3462
            {
3463
                CPLError(CE_Failure, CPLE_OutOfMemory,
×
3464
                         "Out of memory when allocating mask buffer");
3465
                return CE_Failure;
×
3466
            }
3467
            pabyMask = reinterpret_cast<const GByte *>(
3468
                oWorkingState.m_abyWrkBufferMask.data());
103✔
3469
            auto poMaskBand =
51✔
3470
                (poSourceBand->GetColorInterpretation() == GCI_AlphaBand ||
103✔
3471
                 GetMaskBandMainBand() != nullptr)
83✔
3472
                    ? poSourceBand
3473
                    : poSourceBand->GetMaskBand();
52✔
3474
            if (poMaskBand->RasterIO(
103✔
3475
                    GF_Read, nReqXOff, nReqYOff, nReqXSize, nReqYSize,
3476
                    oWorkingState.m_abyWrkBufferMask.data(), nOutXSize,
103✔
3477
                    nOutYSize, GDT_Byte, 1, static_cast<GSpacing>(nOutXSize),
3478
                    psExtraArg) != CE_None)
103✔
3479
            {
3480
                return CE_Failure;
×
3481
            }
3482
        }
3483

3484
        if (m_nColorTableComponent != 0)
685✔
3485
        {
3486
            poColorTable = poSourceBand->GetColorTable();
64✔
3487
            if (poColorTable == nullptr)
64✔
3488
            {
3489
                CPLError(CE_Failure, CPLE_AppDefined,
×
3490
                         "Source band has no color table.");
3491
                return CE_Failure;
×
3492
            }
3493
        }
3494
    }
3495

3496
    /* -------------------------------------------------------------------- */
3497
    /*      Selectively copy into output buffer with nodata masking,        */
3498
    /*      and/or scaling.                                                 */
3499
    /* -------------------------------------------------------------------- */
3500

3501
    const bool bTwoStepDataTypeConversion =
4,725✔
3502
        CPL_TO_BOOL(
4,725✔
3503
            GDALDataTypeIsConversionLossy(eWrkDataType, eVRTBandDataType)) &&
9,393✔
3504
        !CPL_TO_BOOL(GDALDataTypeIsConversionLossy(eVRTBandDataType, eBufType));
4,668✔
3505

3506
    size_t idxBuffer = 0;
4,725✔
3507
    for (int iY = 0; iY < nOutYSize; iY++)
75,142✔
3508
    {
3509
        GByte *pDstLocation = static_cast<GByte *>(pData) +
70,417✔
3510
                              static_cast<GPtrDiff_t>(nLineSpace) * iY;
70,417✔
3511

3512
        for (int iX = 0; iX < nOutXSize;
11,976,442✔
3513
             iX++, pDstLocation += nPixelSpace, idxBuffer++)
11,906,048✔
3514
        {
3515
            WorkingDT afResult[2];
3516
            if (pafData && !bIsComplex)
11,906,048✔
3517
            {
3518
                WorkingDT fResult = pafData[idxBuffer];
6,296,746✔
3519
                if (bNoDataSetIsNan && std::isnan(fResult))
6,296,746✔
3520
                    continue;
1,188,632✔
3521
                if (bNoDataSetAndNotNan &&
6,297,986✔
3522
                    ARE_REAL_EQUAL(fResult, fWorkingDataTypeNoData))
1,266✔
3523
                    continue;
260✔
3524
                if (pabyMask && pabyMask[idxBuffer] == 0)
6,296,454✔
3525
                    continue;
1,188,350✔
3526

3527
                if (poColorTable)
5,108,114✔
3528
                {
3529
                    const GDALColorEntry *poEntry =
3530
                        poColorTable->GetColorEntry(static_cast<int>(fResult));
1,925,900✔
3531
                    if (poEntry)
1,925,900✔
3532
                    {
3533
                        if (m_nColorTableComponent == 1)
1,925,900✔
3534
                            fResult = poEntry->c1;
683,585✔
3535
                        else if (m_nColorTableComponent == 2)
1,242,320✔
3536
                            fResult = poEntry->c2;
529,008✔
3537
                        else if (m_nColorTableComponent == 3)
713,309✔
3538
                            fResult = poEntry->c3;
529,008✔
3539
                        else if (m_nColorTableComponent == 4)
184,301✔
3540
                            fResult = poEntry->c4;
184,301✔
3541
                    }
3542
                    else
3543
                    {
3544
                        CPLErrorOnce(CE_Failure, CPLE_AppDefined,
×
3545
                                     "No entry %d.", static_cast<int>(fResult));
3546
                        continue;
×
3547
                    }
3548
                }
3549

3550
                if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_LINEAR) != 0)
5,108,114✔
3551
                {
3552
                    fResult = static_cast<WorkingDT>(fResult * m_dfScaleRatio +
1,649,610✔
3553
                                                     m_dfScaleOff);
1,649,610✔
3554
                }
3555
                else if ((m_nProcessingFlags &
3,458,504✔
3556
                          PROCESSING_FLAG_SCALING_EXPONENTIAL) != 0)
3557
                {
3558
                    if (!m_bSrcMinMaxDefined)
3,543✔
3559
                    {
3560
                        int bSuccessMin = FALSE;
1✔
3561
                        int bSuccessMax = FALSE;
1✔
3562
                        double adfMinMax[2] = {
2✔
3563
                            poSourceBand->GetMinimum(&bSuccessMin),
1✔
3564
                            poSourceBand->GetMaximum(&bSuccessMax)};
1✔
3565
                        if ((bSuccessMin && bSuccessMax) ||
2✔
3566
                            poSourceBand->ComputeRasterMinMax(
1✔
3567
                                TRUE, adfMinMax) == CE_None)
3568
                        {
3569
                            m_dfSrcMin = adfMinMax[0];
1✔
3570
                            m_dfSrcMax = adfMinMax[1];
1✔
3571
                            m_bSrcMinMaxDefined = true;
1✔
3572
                        }
3573
                        else
3574
                        {
3575
                            CPLError(CE_Failure, CPLE_AppDefined,
×
3576
                                     "Cannot determine source min/max value");
3577
                            return CE_Failure;
×
3578
                        }
3579
                    }
3580

3581
                    double dfPowVal =
3,543✔
3582
                        (fResult - m_dfSrcMin) / (m_dfSrcMax - m_dfSrcMin);
3,543✔
3583
                    if (dfPowVal < 0.0)
3,543✔
3584
                        dfPowVal = 0.0;
×
3585
                    else if (dfPowVal > 1.0)
3,543✔
3586
                        dfPowVal = 1.0;
699✔
3587
                    fResult =
3,543✔
3588
                        static_cast<WorkingDT>((m_dfDstMax - m_dfDstMin) *
3,543✔
3589
                                                   pow(dfPowVal, m_dfExponent) +
3,543✔
3590
                                               m_dfDstMin);
3,543✔
3591
                }
3592

3593
                if (!m_adfLUTInputs.empty())
5,108,114✔
3594
                    fResult = static_cast<WorkingDT>(LookupValue(fResult));
583,336✔
3595

3596
                if (m_nMaxValue != 0 && fResult > m_nMaxValue)
5,108,114✔
3597
                    fResult = static_cast<WorkingDT>(m_nMaxValue);
800✔
3598

3599
                afResult[0] = fResult;
5,108,114✔
3600
                afResult[1] = 0;
5,108,114✔
3601
            }
3602
            else if (pafData && bIsComplex)
5,609,252✔
3603
            {
3604
                afResult[0] = pafData[2 * idxBuffer];
1,015✔
3605
                afResult[1] = pafData[2 * idxBuffer + 1];
1,015✔
3606

3607
                // Do not use color table.
3608
                if ((m_nProcessingFlags & PROCESSING_FLAG_SCALING_LINEAR) != 0)
1,015✔
3609
                {
3610
                    afResult[0] = static_cast<WorkingDT>(
1✔
3611
                        afResult[0] * m_dfScaleRatio + m_dfScaleOff);
1✔
3612
                    afResult[1] = static_cast<WorkingDT>(
1✔
3613
                        afResult[1] * m_dfScaleRatio + m_dfScaleOff);
1✔
3614
                }
3615

3616
                /* Do not use LUT */
3617
            }
3618
            else
3619
            {
3620
                afResult[0] = static_cast<WorkingDT>(m_dfScaleOff);
5,608,242✔
3621
                afResult[1] = 0;
5,608,242✔
3622

3623
                if (!m_adfLUTInputs.empty())
5,608,242✔
3624
                    afResult[0] =
×
3625
                        static_cast<WorkingDT>(LookupValue(afResult[0]));
×
3626

3627
                if (m_nMaxValue != 0 && afResult[0] > m_nMaxValue)
5,608,242✔
3628
                    afResult[0] = static_cast<WorkingDT>(m_nMaxValue);
×
3629
            }
3630

3631
            if (eBufType == GDT_Byte && eVRTBandDataType == GDT_Byte)
10,717,346✔
3632
            {
3633
                *pDstLocation = static_cast<GByte>(std::min(
4,299,480✔
3634
                    255.0f,
8,598,960✔
3635
                    std::max(0.0f, static_cast<float>(afResult[0]) + 0.5f)));
4,299,480✔
3636
            }
3637
            else if (!bTwoStepDataTypeConversion)
6,417,886✔
3638
            {
3639
                CopyWord<WorkingDT>(afResult, eWrkDataType, pDstLocation,
4,021✔
3640
                                    eBufType);
3641
            }
3642
            else if (eVRTBandDataType == GDT_Float32 && eBufType == GDT_Float64)
6,413,870✔
3643
            {
3644
                // Particular case of the below 2-step conversion.
3645
                // Helps a bit for some geolocation based warping with Sentinel3
3646
                // data where the longitude/latitude arrays are Int32 bands,
3647
                // rescaled in VRT as Float32 and requested as Float64
3648
                float fVal;
3649
                GDALCopyWord(afResult[0], fVal);
20✔
3650
                *reinterpret_cast<double *>(pDstLocation) = fVal;
20✔
3651
            }
3652
            else
3653
            {
3654
                GByte abyTemp[2 * sizeof(double)];
3655
                // Convert first to the VRTRasterBand data type
3656
                // to get its clamping, before outputting to buffer data type
3657
                CopyWord<WorkingDT>(afResult, eWrkDataType, abyTemp,
6,413,850✔
3658
                                    eVRTBandDataType);
3659
                GDALCopyWords(abyTemp, eVRTBandDataType, 0, pDstLocation,
6,413,850✔
3660
                              eBufType, 0, 1);
3661
            }
3662
        }
3663
    }
3664

3665
    return CE_None;
4,725✔
3666
}
3667

3668
// Explicitly instantiate template method, as it is used in another file.
3669
template CPLErr VRTComplexSource::RasterIOInternal<float>(
3670
    GDALRasterBand *poSourceBand, GDALDataType eVRTBandDataType, int nReqXOff,
3671
    int nReqYOff, int nReqXSize, int nReqYSize, void *pData, int nOutXSize,
3672
    int nOutYSize, GDALDataType eBufType, GSpacing nPixelSpace,
3673
    GSpacing nLineSpace, GDALRasterIOExtraArg *psExtraArg,
3674
    GDALDataType eWrkDataType, WorkingState &oWorkingState);
3675

3676
/************************************************************************/
3677
/*                        AreValuesUnchanged()                          */
3678
/************************************************************************/
3679

3680
bool VRTComplexSource::AreValuesUnchanged() const
9✔
3681
{
3682
    return m_dfScaleOff == 0.0 && m_dfScaleRatio == 1.0 &&
10✔
3683
           m_adfLUTInputs.empty() && m_nColorTableComponent == 0 &&
24✔
3684
           (m_nProcessingFlags & PROCESSING_FLAG_SCALING_EXPONENTIAL) == 0;
14✔
3685
}
3686

3687
/************************************************************************/
3688
/*                             GetMinimum()                             */
3689
/************************************************************************/
3690

3691
double VRTComplexSource::GetMinimum(int nXSize, int nYSize, int *pbSuccess)
2✔
3692
{
3693
    if (AreValuesUnchanged())
2✔
3694
    {
3695
        return VRTSimpleSource::GetMinimum(nXSize, nYSize, pbSuccess);
1✔
3696
    }
3697

3698
    *pbSuccess = FALSE;
1✔
3699
    return 0;
1✔
3700
}
3701

3702
/************************************************************************/
3703
/*                             GetMaximum()                             */
3704
/************************************************************************/
3705

3706
double VRTComplexSource::GetMaximum(int nXSize, int nYSize, int *pbSuccess)
2✔
3707
{
3708
    if (AreValuesUnchanged())
2✔
3709
    {
3710
        return VRTSimpleSource::GetMaximum(nXSize, nYSize, pbSuccess);
1✔
3711
    }
3712

3713
    *pbSuccess = FALSE;
1✔
3714
    return 0;
1✔
3715
}
3716

3717
/************************************************************************/
3718
/*                            GetHistogram()                            */
3719
/************************************************************************/
3720

3721
CPLErr VRTComplexSource::GetHistogram(int nXSize, int nYSize, double dfMin,
×
3722
                                      double dfMax, int nBuckets,
3723
                                      GUIntBig *panHistogram,
3724
                                      int bIncludeOutOfRange, int bApproxOK,
3725
                                      GDALProgressFunc pfnProgress,
3726
                                      void *pProgressData)
3727
{
3728
    if (AreValuesUnchanged())
×
3729
    {
3730
        return VRTSimpleSource::GetHistogram(
×
3731
            nXSize, nYSize, dfMin, dfMax, nBuckets, panHistogram,
3732
            bIncludeOutOfRange, bApproxOK, pfnProgress, pProgressData);
×
3733
    }
3734

3735
    return CE_Failure;
×
3736
}
3737

3738
/************************************************************************/
3739
/* ==================================================================== */
3740
/*                          VRTFuncSource                               */
3741
/* ==================================================================== */
3742
/************************************************************************/
3743

3744
/************************************************************************/
3745
/*                           VRTFuncSource()                            */
3746
/************************************************************************/
3747

3748
VRTFuncSource::VRTFuncSource()
12✔
3749
    : pfnReadFunc(nullptr), pCBData(nullptr), eType(GDT_Byte),
3750
      fNoDataValue(static_cast<float>(VRT_NODATA_UNSET))
12✔
3751
{
3752
}
12✔
3753

3754
/************************************************************************/
3755
/*                           ~VRTFuncSource()                           */
3756
/************************************************************************/
3757

3758
VRTFuncSource::~VRTFuncSource()
24✔
3759
{
3760
}
24✔
3761

3762
/************************************************************************/
3763
/*                            GetType()                                 */
3764
/************************************************************************/
3765

3766
const char *VRTFuncSource::GetType() const
×
3767
{
3768
    static const char *TYPE = "FuncSource";
3769
    return TYPE;
×
3770
}
3771

3772
/************************************************************************/
3773
/*                           SerializeToXML()                           */
3774
/************************************************************************/
3775

3776
CPLXMLNode *VRTFuncSource::SerializeToXML(CPL_UNUSED const char *pszVRTPath)
×
3777
{
3778
    return nullptr;
×
3779
}
3780

3781
/************************************************************************/
3782
/*                              RasterIO()                              */
3783
/************************************************************************/
3784

3785
CPLErr VRTFuncSource::RasterIO(GDALDataType /*eVRTBandDataType*/, int nXOff,
10✔
3786
                               int nYOff, int nXSize, int nYSize, void *pData,
3787
                               int nBufXSize, int nBufYSize,
3788
                               GDALDataType eBufType, GSpacing nPixelSpace,
3789
                               GSpacing nLineSpace,
3790
                               GDALRasterIOExtraArg * /* psExtraArg */,
3791
                               WorkingState & /* oWorkingState */)
3792
{
3793
    if (nPixelSpace * 8 == GDALGetDataTypeSize(eBufType) &&
10✔
3794
        nLineSpace == nPixelSpace * nXSize && nBufXSize == nXSize &&
10✔
3795
        nBufYSize == nYSize && eBufType == eType)
20✔
3796
    {
3797
        return pfnReadFunc(pCBData, nXOff, nYOff, nXSize, nYSize, pData);
10✔
3798
    }
3799
    else
3800
    {
3801
        CPLError(CE_Failure, CPLE_AppDefined,
×
3802
                 "VRTFuncSource::RasterIO() - Irregular request.");
3803
        CPLDebug("VRT", "Irregular request: %d,%d  %d,%d, %d,%d %d,%d %d,%d",
×
3804
                 static_cast<int>(nPixelSpace) * 8,
3805
                 GDALGetDataTypeSize(eBufType), static_cast<int>(nLineSpace),
3806
                 static_cast<int>(nPixelSpace) * nXSize, nBufXSize, nXSize,
3807
                 nBufYSize, nYSize, static_cast<int>(eBufType),
3808
                 static_cast<int>(eType));
×
3809

3810
        return CE_Failure;
×
3811
    }
3812
}
3813

3814
/************************************************************************/
3815
/*                             GetMinimum()                             */
3816
/************************************************************************/
3817

3818
double VRTFuncSource::GetMinimum(int /* nXSize */, int /* nYSize */,
×
3819
                                 int *pbSuccess)
3820
{
3821
    *pbSuccess = FALSE;
×
3822
    return 0;
×
3823
}
3824

3825
/************************************************************************/
3826
/*                             GetMaximum()                             */
3827
/************************************************************************/
3828

3829
double VRTFuncSource::GetMaximum(int /* nXSize */, int /* nYSize */,
×
3830
                                 int *pbSuccess)
3831
{
3832
    *pbSuccess = FALSE;
×
3833
    return 0;
×
3834
}
3835

3836
/************************************************************************/
3837
/*                            GetHistogram()                            */
3838
/************************************************************************/
3839

3840
CPLErr VRTFuncSource::GetHistogram(
×
3841
    int /* nXSize */, int /* nYSize */, double /* dfMin */, double /* dfMax */,
3842
    int /* nBuckets */, GUIntBig * /* panHistogram */,
3843
    int /* bIncludeOutOfRange */, int /* bApproxOK */,
3844
    GDALProgressFunc /* pfnProgress */, void * /* pProgressData */)
3845
{
3846
    return CE_Failure;
×
3847
}
3848

3849
/************************************************************************/
3850
/*                        VRTParseCoreSources()                         */
3851
/************************************************************************/
3852

3853
VRTSource *VRTParseCoreSources(const CPLXMLNode *psChild,
3,155✔
3854
                               const char *pszVRTPath,
3855
                               VRTMapSharedResources &oMapSharedSources)
3856

3857
{
3858
    VRTSource *poSource = nullptr;
3,155✔
3859

3860
    if (EQUAL(psChild->pszValue, VRTAveragedSource::GetTypeStatic()) ||
6,303✔
3861
        (EQUAL(psChild->pszValue, VRTSimpleSource::GetTypeStatic()) &&
3,148✔
3862
         STARTS_WITH_CI(CPLGetXMLValue(psChild, "Resampling", "Nearest"),
2,937✔
3863
                        "Aver")))
3864
    {
3865
        poSource = new VRTAveragedSource();
16✔
3866
    }
3867
    else if (EQUAL(psChild->pszValue, VRTSimpleSource::GetTypeStatic()))
3,139✔
3868
    {
3869
        poSource = new VRTSimpleSource();
2,928✔
3870
    }
3871
    else if (EQUAL(psChild->pszValue, VRTComplexSource::GetTypeStatic()))
211✔
3872
    {
3873
        poSource = new VRTComplexSource();
203✔
3874
    }
3875
    else if (EQUAL(psChild->pszValue, VRTNoDataFromMaskSource::GetTypeStatic()))
8✔
3876
    {
3877
        poSource = new VRTNoDataFromMaskSource();
8✔
3878
    }
3879
    else
3880
    {
3881
        CPLError(CE_Failure, CPLE_AppDefined,
×
3882
                 "VRTParseCoreSources() - Unknown source : %s",
3883
                 psChild->pszValue);
×
3884
        return nullptr;
×
3885
    }
3886

3887
    if (poSource->XMLInit(psChild, pszVRTPath, oMapSharedSources) == CE_None)
3,155✔
3888
        return poSource;
3,152✔
3889

3890
    delete poSource;
3✔
3891
    return nullptr;
3✔
3892
}
3893

3894
/*! @endcond */
STATUS · Troubleshooting · Open an Issue · Sales · Support · CAREERS · ENTERPRISE · START FREE · SCHEDULE DEMO
ANNOUNCEMENTS · TWITTER · TOS & SLA · Supported CI Services · What's a CI service? · Automated Testing

© 2025 Coveralls, Inc