text: /******************************************************************************
   *
   * Project:  Virtual GDAL Datasets
   * Purpose:  Implementation of VRTSimpleSource, VRTFuncSource and
   *           VRTAveragedSource.
   * Author:   Frank Warmerdam <warmerdam@pobox.com>
   *
   ******************************************************************************
   * Copyright (c) 2001, Frank Warmerdam <warmerdam@pobox.com>
   * Copyright (c) 2008-2013, Even Rouault <even dot rouault at mines-paris dot org>
   *
   * Permission is hereby granted, free of charge, to any person obtaining a
   * copy of this software and associated documentation files (the "Software"),
   * to deal in the Software without restriction, including without limitation
   * the rights to use, copy, modify, merge, publish, distribute, sublicense,
   * and/or sell copies of the Software, and to permit persons to whom the
   * Software is furnished to do so, subject to the following conditions:
   *
   * The above copyright notice and this permission notice shall be included
   * in all copies or substantial portions of the Software.
   *
   * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
   * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
   * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
   * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
   * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
   * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
   * DEALINGS IN THE SOFTWARE.
   ****************************************************************************/
  
  #include "gdal_vrt.h"
  #include "vrtdataset.h"
  
  #include <climits>
  #include <cmath>
  #include <cstddef>
  #include <cstdio>
  #include <cstdlib>
  #include <cstring>
  #include <algorithm>
  #include <string>
  
  #include "cpl_conv.h"
  #include "cpl_error.h"
  #include "cpl_hash_set.h"
  #include "cpl_minixml.h"
  #include "cpl_progress.h"
  #include "cpl_string.h"
  #include "cpl_vsi.h"
  #include "gdal.h"
  #include "gdal_priv.h"
  #include "gdal_proxy.h"
  #include "gdal_priv_templates.hpp"
  
  /*! @cond Doxygen_Suppress */
  
  // #define DEBUG_VERBOSE 1
  
  // See #5459
  #ifdef isnan
  #define HAS_ISNAN_MACRO
  #endif
  #include <algorithm>
  #if defined(HAS_ISNAN_MACRO) && !defined(isnan)
  #define isnan std::isnan
  #endif
  
  CPL_CVSID("$Id$")
  
  /************************************************************************/
  /* ==================================================================== */
  /*                             VRTSource                                */
  /* ==================================================================== */
  /************************************************************************/
  
  VRTSource::~VRTSource() {}
  
  /************************************************************************/
  /*                             GetFileList()                            */
  /************************************************************************/
  
  void VRTSource::GetFileList(char*** /* ppapszFileList */,
                              int * /* pnSize */,
                              int * /* pnMaxSize */,
                              CPLHashSet * /* hSetFiles */)
  {}
  
  /************************************************************************/
  /* ==================================================================== */
  /*                          VRTSimpleSource                             */
  /* ==================================================================== */
  /************************************************************************/
  
  /************************************************************************/
  /*                          VRTSimpleSource()                           */
  /************************************************************************/
  
  VRTSimpleSource::VRTSimpleSource() :
      m_poRasterBand(nullptr),
      m_poMaskBandMainBand(nullptr),
      m_dfSrcXOff(0.0),
      m_dfSrcYOff(0.0),
      m_dfSrcXSize(0.0),
      m_dfSrcYSize(0.0),
      m_dfDstXOff(0.0),
      m_dfDstYOff(0.0),
      m_dfDstXSize(0.0),
      m_dfDstYSize(0.0),
      m_bNoDataSet(FALSE),
      m_dfNoDataValue(VRT_NODATA_UNSET),
      m_nMaxValue(0),
      m_bRelativeToVRTOri(-1),
      m_nExplicitSharedStatus(-1)
  {}
  
  /************************************************************************/
  /*                          VRTSimpleSource()                           */
  /************************************************************************/
  
  VRTSimpleSource::VRTSimpleSource( const VRTSimpleSource* poSrcSource,
                                    double dfXDstRatio, double dfYDstRatio ) :
      m_poRasterBand(poSrcSource->m_poRasterBand),
      m_poMaskBandMainBand(poSrcSource->m_poMaskBandMainBand),
      m_dfSrcXOff(poSrcSource->m_dfSrcXOff),
      m_dfSrcYOff(poSrcSource->m_dfSrcYOff),
      m_dfSrcXSize(poSrcSource->m_dfSrcXSize),
      m_dfSrcYSize(poSrcSource->m_dfSrcYSize),
      m_dfDstXOff(poSrcSource->m_dfDstXOff * dfXDstRatio),
      m_dfDstYOff(poSrcSource->m_dfDstYOff * dfYDstRatio),
      m_dfDstXSize(poSrcSource->m_dfDstXSize * dfXDstRatio),
      m_dfDstYSize(poSrcSource->m_dfDstYSize * dfYDstRatio),
      m_bNoDataSet(poSrcSource->m_bNoDataSet),
      m_dfNoDataValue(poSrcSource->m_dfNoDataValue),
      m_nMaxValue(poSrcSource->m_nMaxValue),
      m_bRelativeToVRTOri(-1),
      m_nExplicitSharedStatus(poSrcSource->m_nExplicitSharedStatus)
  {}
  
  /************************************************************************/
  /*                          ~VRTSimpleSource()                          */
  /************************************************************************/
  
  VRTSimpleSource::~VRTSimpleSource()
  
  {
      if( m_poMaskBandMainBand != nullptr )
      {
          if( m_poMaskBandMainBand->GetDataset() != nullptr )
          {
              m_poMaskBandMainBand->GetDataset()->ReleaseRef();
          }
      }
      else if( m_poRasterBand != nullptr && m_poRasterBand->GetDataset() != nullptr )
      {
          m_poRasterBand->GetDataset()->ReleaseRef();
      }
  }
  
  /************************************************************************/
  /*                           FlushCache()                               */
  /************************************************************************/
  
  CPLErr VRTSimpleSource::FlushCache()
  
  {
      if( m_poMaskBandMainBand != nullptr )
      {
          return m_poMaskBandMainBand->FlushCache();
      }
      else if( m_poRasterBand != nullptr )
      {
          return m_poRasterBand->FlushCache();
      }
      return CE_None;
  }
  
  /************************************************************************/
  /*                    UnsetPreservedRelativeFilenames()                 */
  /************************************************************************/
  
  void VRTSimpleSource::UnsetPreservedRelativeFilenames()
  {
      m_bRelativeToVRTOri = -1;
      m_osSourceFileNameOri = "";
  }
  
  /************************************************************************/
  /*                             SetSrcBand()                             */
  /************************************************************************/
  
  void VRTSimpleSource::SetSrcBand( GDALRasterBand *poNewSrcBand )
  
  {
      m_poRasterBand = poNewSrcBand;
  }
  
  /************************************************************************/
  /*                          SetSrcMaskBand()                            */
  /************************************************************************/
  
  // poSrcBand is not the mask band, but the band from which the mask band is
  // taken.
  void VRTSimpleSource::SetSrcMaskBand( GDALRasterBand *poNewSrcBand )
  
  {
      m_poRasterBand = poNewSrcBand->GetMaskBand();
      m_poMaskBandMainBand = poNewSrcBand;
  }
  
  /************************************************************************/
  /*                         RoundIfCloseToInt()                          */
  /************************************************************************/
  
  static double RoundIfCloseToInt(double dfValue)
  {
      double dfClosestInt = floor(dfValue + 0.5);
      return (fabs( dfValue - dfClosestInt ) < 1e-5) ? dfClosestInt : dfValue;
  }
  
  /************************************************************************/
  /*                            SetSrcWindow()                            */
  /************************************************************************/
  
  void VRTSimpleSource::SetSrcWindow( double dfNewXOff, double dfNewYOff,
                                      double dfNewXSize, double dfNewYSize )
  
  {
      m_dfSrcXOff = RoundIfCloseToInt(dfNewXOff);
      m_dfSrcYOff = RoundIfCloseToInt(dfNewYOff);
      m_dfSrcXSize = RoundIfCloseToInt(dfNewXSize);
      m_dfSrcYSize = RoundIfCloseToInt(dfNewYSize);
  }
  
  /************************************************************************/
  /*                            SetDstWindow()                            */
  /************************************************************************/
  
  void VRTSimpleSource::SetDstWindow( double dfNewXOff, double dfNewYOff,
                                      double dfNewXSize, double dfNewYSize )
  
  {
      m_dfDstXOff = RoundIfCloseToInt(dfNewXOff);
      m_dfDstYOff = RoundIfCloseToInt(dfNewYOff);
      m_dfDstXSize = RoundIfCloseToInt(dfNewXSize);
      m_dfDstYSize = RoundIfCloseToInt(dfNewYSize);
  }
  
  /************************************************************************/
  /*                           SetNoDataValue()                           */
  /************************************************************************/
  
  void VRTSimpleSource::SetNoDataValue( double dfNewNoDataValue )
  
  {
      if( dfNewNoDataValue == VRT_NODATA_UNSET )
      {
          m_bNoDataSet = FALSE;
          m_dfNoDataValue = VRT_NODATA_UNSET;
          return;
      }
  
      m_bNoDataSet = TRUE;
      m_dfNoDataValue = dfNewNoDataValue;
  }
  
  /************************************************************************/
  /*                           SerializeToXML()                           */
  /************************************************************************/
  
  static const char* const apszSpecialSyntax[] = {
      "HDF5:\"{FILENAME}\":{ANY}",
      "HDF5:{FILENAME}:{ANY}",
      "NETCDF:\"{FILENAME}\":{ANY}",
      "NETCDF:{FILENAME}:{ANY}",
      "NITF_IM:{ANY}:{FILENAME}",
      "PDF:{ANY}:{FILENAME}",
      "RASTERLITE:{FILENAME},{ANY}" };
  
  CPLXMLNode *VRTSimpleSource::SerializeToXML( const char *pszVRTPath )
  
  {
      if( m_poRasterBand == nullptr )
          return nullptr;
  
      GDALDataset *poDS = nullptr;
  
      if( m_poMaskBandMainBand )
      {
          poDS = m_poMaskBandMainBand->GetDataset();
          if( poDS == nullptr || m_poMaskBandMainBand->GetBand() < 1 )
              return nullptr;
      }
      else
      {
          poDS = m_poRasterBand->GetDataset();
          if( poDS == nullptr || m_poRasterBand->GetBand() < 1 )
              return nullptr;
      }
  
      CPLXMLNode * const psSrc =
          CPLCreateXMLNode( nullptr, CXT_Element, "SimpleSource" );
  
      if( !m_osResampling.empty() )
      {
          CPLCreateXMLNode(
              CPLCreateXMLNode( psSrc, CXT_Attribute, "resampling" ),
                  CXT_Text, m_osResampling.c_str() );
      }
  
      VSIStatBufL sStat;
      CPLString osTmp;
      int bRelativeToVRT = FALSE;  // TODO(schwehr): Make this a bool?
      const char *pszRelativePath = nullptr;
  
      if( m_bRelativeToVRTOri >= 0 )
      {
          pszRelativePath = m_osSourceFileNameOri;
          bRelativeToVRT = m_bRelativeToVRTOri;
      }
      else if( strstr(poDS->GetDescription(), "/vsicurl/http") != nullptr ||
               strstr(poDS->GetDescription(), "/vsicurl/ftp") != nullptr )
      {
          // Testing the existence of remote resources can be excruciating
          // slow, so let's just suppose they exist.
          pszRelativePath = poDS->GetDescription();
          bRelativeToVRT = FALSE;
      }
      // If this isn't actually a file, don't even try to know if it is a
      // relative path. It can't be !, and unfortunately CPLIsFilenameRelative()
      // can only work with strings that are filenames To be clear
      // NITF_TOC_ENTRY:CADRG_JOG-A_250K_1_0:some_path isn't a relative file
      // path.
      else if( VSIStatExL( poDS->GetDescription(), &sStat,
                           VSI_STAT_EXISTS_FLAG ) != 0 )
      {
          pszRelativePath = poDS->GetDescription();
          bRelativeToVRT = FALSE;
          for( size_t i = 0;
               i < sizeof(apszSpecialSyntax) / sizeof(apszSpecialSyntax[0]);
               ++i )
          {
              const char* const pszSyntax = apszSpecialSyntax[i];
              CPLString osPrefix(pszSyntax);
              osPrefix.resize(strchr(pszSyntax, ':') - pszSyntax + 1);
              if( pszSyntax[osPrefix.size()] == '"' )
                  osPrefix += '"';
              if( EQUALN(pszRelativePath, osPrefix, osPrefix.size()) )
              {
                  if( STARTS_WITH_CI(pszSyntax + osPrefix.size(), "{ANY}") )
                  {
                      const char* pszLastPart = strrchr(pszRelativePath, ':') + 1;
                      // CSV:z:/foo.xyz
                      if( (pszLastPart[0] == '/' || pszLastPart[0] == '\\') &&
                          pszLastPart - pszRelativePath >= 3 &&
                          pszLastPart[-3] == ':' )
                          pszLastPart -= 2;
                      CPLString osPrefixFilename(pszRelativePath);
                      osPrefixFilename.resize(pszLastPart - pszRelativePath);
                      pszRelativePath =
                          CPLExtractRelativePath( pszVRTPath, pszLastPart,
                                                  &bRelativeToVRT );
                      osTmp = osPrefixFilename + pszRelativePath;
                      pszRelativePath = osTmp.c_str();
                  }
                  else if( STARTS_WITH_CI(pszSyntax + osPrefix.size(),
                                          "{FILENAME}") )
                  {
                      CPLString osFilename(pszRelativePath + osPrefix.size());
                      size_t nPos = 0;
                      if( osFilename.size() >= 3 && osFilename[1] == ':' &&
                          (osFilename[2] == '\\' || osFilename[2] == '/') )
                          nPos = 2;
                      nPos =
                          osFilename.find(
                              pszSyntax[osPrefix.size() + strlen("{FILENAME}")],
                              nPos );
                      if( nPos != std::string::npos )
                      {
                          const CPLString osSuffix = osFilename.substr(nPos);
                          osFilename.resize(nPos);
                          pszRelativePath =
                              CPLExtractRelativePath( pszVRTPath, osFilename,
                                          &bRelativeToVRT );
                          osTmp = osPrefix + pszRelativePath + osSuffix;
                          pszRelativePath = osTmp.c_str();
                      }
                  }
                  break;
              }
          }
      }
      else
      {
          pszRelativePath =
              CPLExtractRelativePath( pszVRTPath, poDS->GetDescription(),
                                      &bRelativeToVRT );
      }
  
      CPLSetXMLValue( psSrc, "SourceFilename", pszRelativePath );
  
      CPLCreateXMLNode(
          CPLCreateXMLNode( CPLGetXMLNode( psSrc, "SourceFilename" ),
                            CXT_Attribute, "relativeToVRT" ),
          CXT_Text, bRelativeToVRT ? "1" : "0" );
  
      // Determine if we must write the shared attribute. The config option
      // will override the m_nExplicitSharedStatus value
      const char* pszShared = CPLGetConfigOption("VRT_SHARED_SOURCE", nullptr);
      if( (pszShared == nullptr && m_nExplicitSharedStatus == 0) ||
          (pszShared != nullptr && !CPLTestBool(pszShared)) )
      {
          CPLCreateXMLNode(
              CPLCreateXMLNode( CPLGetXMLNode( psSrc, "SourceFilename" ),
                                CXT_Attribute, "shared" ),
                                CXT_Text, "0" );
      }
  
      char** papszOpenOptions = poDS->GetOpenOptions();
      GDALSerializeOpenOptionsToXML(psSrc, papszOpenOptions);
  
      if( m_poMaskBandMainBand )
          CPLSetXMLValue( psSrc, "SourceBand",
                          CPLSPrintf("mask,%d",m_poMaskBandMainBand->GetBand()) );
      else
          CPLSetXMLValue( psSrc, "SourceBand",
                          CPLSPrintf("%d",m_poRasterBand->GetBand()) );
  
      /* Write a few additional useful properties of the dataset */
      /* so that we can use a proxy dataset when re-opening. See XMLInit() */
      /* below */
      CPLSetXMLValue( psSrc, "SourceProperties.#RasterXSize",
                      CPLSPrintf("%d",m_poRasterBand->GetXSize()) );
      CPLSetXMLValue( psSrc, "SourceProperties.#RasterYSize",
                      CPLSPrintf("%d",m_poRasterBand->GetYSize()) );
      CPLSetXMLValue( psSrc, "SourceProperties.#DataType",
                  GDALGetDataTypeName( m_poRasterBand->GetRasterDataType() ) );
  
      int nBlockXSize = 0;
      int nBlockYSize = 0;
      m_poRasterBand->GetBlockSize(&nBlockXSize, &nBlockYSize);
  
      CPLSetXMLValue( psSrc, "SourceProperties.#BlockXSize",
                      CPLSPrintf("%d",nBlockXSize) );
      CPLSetXMLValue( psSrc, "SourceProperties.#BlockYSize",
                      CPLSPrintf("%d",nBlockYSize) );
  
      if( m_dfSrcXOff != -1
          || m_dfSrcYOff != -1
          || m_dfSrcXSize != -1
          || m_dfSrcYSize != -1 )
      {
          CPLSetXMLValue( psSrc, "SrcRect.#xOff",
                          CPLSPrintf( "%.15g", m_dfSrcXOff ) );
          CPLSetXMLValue( psSrc, "SrcRect.#yOff",
                          CPLSPrintf( "%.15g", m_dfSrcYOff ) );
          CPLSetXMLValue( psSrc, "SrcRect.#xSize",
                          CPLSPrintf( "%.15g", m_dfSrcXSize ) );
          CPLSetXMLValue( psSrc, "SrcRect.#ySize",
                          CPLSPrintf( "%.15g", m_dfSrcYSize ) );
      }
  
      if( m_dfDstXOff != -1
          || m_dfDstYOff != -1
          || m_dfDstXSize != -1
          || m_dfDstYSize != -1 )
      {
          CPLSetXMLValue( psSrc, "DstRect.#xOff", CPLSPrintf( "%.15g", m_dfDstXOff ) );
          CPLSetXMLValue( psSrc, "DstRect.#yOff", CPLSPrintf( "%.15g", m_dfDstYOff ) );
          CPLSetXMLValue( psSrc, "DstRect.#xSize",CPLSPrintf( "%.15g", m_dfDstXSize ));
          CPLSetXMLValue( psSrc, "DstRect.#ySize",CPLSPrintf( "%.15g", m_dfDstYSize ));
      }
  
      return psSrc;
  }
  
  /************************************************************************/
  /*                              XMLInit()                               */
  /************************************************************************/
  
  CPLErr VRTSimpleSource::XMLInit( CPLXMLNode *psSrc, const char *pszVRTPath,
                                   void* pUniqueHandle )
  
  {
      m_osResampling = CPLGetXMLValue( psSrc, "resampling", "");
  
  /* -------------------------------------------------------------------- */
  /*      Prepare filename.                                               */
  /* -------------------------------------------------------------------- */
      CPLXMLNode* psSourceFileNameNode = CPLGetXMLNode(psSrc,"SourceFilename");
      const char *pszFilename =
          psSourceFileNameNode ?
          CPLGetXMLValue(psSourceFileNameNode, nullptr, nullptr) : nullptr;
  
      if( pszFilename == nullptr )
      {
          CPLError( CE_Warning, CPLE_AppDefined,
                    "Missing <SourceFilename> element in VRTRasterBand." );
          return CE_Failure;
      }
  
      // Backup original filename and relativeToVRT so as to be able to
      // serialize them identically again (#5985)
      m_osSourceFileNameOri = pszFilename;
      m_bRelativeToVRTOri =
          atoi( CPLGetXMLValue( psSourceFileNameNode, "relativetoVRT", "0") );
      const char* pszShared = CPLGetXMLValue( psSourceFileNameNode,
                                              "shared", nullptr );
      if( pszShared == nullptr )
      {
          pszShared = CPLGetConfigOption("VRT_SHARED_SOURCE", nullptr );
      }
      bool bShared = false;
      if( pszShared != nullptr )
      {
          bShared = CPLTestBool(pszShared);
          m_nExplicitSharedStatus = bShared;
      }
      else
      {
          bShared = true;
      }
  
      char *pszSrcDSName = nullptr;
      if( pszVRTPath != nullptr && m_bRelativeToVRTOri )
      {
          bool bDone = false;
          for( size_t i = 0;
               i < sizeof(apszSpecialSyntax) / sizeof(apszSpecialSyntax[0]);
               ++i )
          {
              const char* pszSyntax = apszSpecialSyntax[i];
              CPLString osPrefix(pszSyntax);
              osPrefix.resize(strchr(pszSyntax, ':') - pszSyntax + 1);
              if( pszSyntax[osPrefix.size()] == '"' )
                  osPrefix += '"';
              if( EQUALN(pszFilename, osPrefix, osPrefix.size()) )
              {
                  if( STARTS_WITH_CI(pszSyntax + osPrefix.size(), "{ANY}") )
                  {
                      const char * pszLastPart = strrchr(pszFilename, ':') + 1;
                      // CSV:z:/foo.xyz
                      if( ( pszLastPart[0] == '/' ||
                            pszLastPart[0] == '\\') &&
                          pszLastPart - pszFilename >= 3 &&
                          pszLastPart[-3] == ':' )
                      {
                          pszLastPart -= 2;
                      }
                      CPLString osPrefixFilename = pszFilename;
                      osPrefixFilename.resize(pszLastPart - pszFilename);
                      pszSrcDSName = CPLStrdup( (osPrefixFilename +
                          CPLProjectRelativeFilename( pszVRTPath, pszLastPart )).c_str() );
                      bDone = true;
                  }
                  else if( STARTS_WITH_CI(pszSyntax + osPrefix.size(), "{FILENAME}") )
                  {
                      CPLString osFilename(pszFilename + osPrefix.size());
                      size_t nPos = 0;
                      if( osFilename.size() >= 3 && osFilename[1] == ':' &&
                          (osFilename[2] == '\\' || osFilename[2] == '/') )
                          nPos = 2;
                      nPos = osFilename.find(
                          pszSyntax[osPrefix.size() + strlen("{FILENAME}")],
                          nPos);
                      if( nPos != std::string::npos )
                      {
                          const CPLString osSuffix = osFilename.substr(nPos);
                          osFilename.resize(nPos);
                          pszSrcDSName = CPLStrdup(
                              (osPrefix + CPLProjectRelativeFilename(
                                  pszVRTPath, osFilename ) + osSuffix).c_str() );
                          bDone = true;
                      }
                  }
                  break;
              }
          }
          if( !bDone )
          {
              pszSrcDSName = CPLStrdup(
                  CPLProjectRelativeFilename( pszVRTPath, pszFilename ) );
          }
      }
      else
      {
          pszSrcDSName = CPLStrdup( pszFilename );
      }
  
      const char* pszSourceBand = CPLGetXMLValue(psSrc,"SourceBand","1");
      int nSrcBand = 0;
      bool bGetMaskBand = false;
      if( STARTS_WITH_CI(pszSourceBand, "mask") )
      {
          bGetMaskBand = true;
          if( pszSourceBand[4] == ',' )
              nSrcBand = atoi(pszSourceBand + 5);
          else
              nSrcBand = 1;
      }
      else
      {
          nSrcBand = atoi(pszSourceBand);
      }
      if( !GDALCheckBandCount(nSrcBand, 0) )
      {
          CPLError( CE_Warning, CPLE_AppDefined,
                    "Invalid <SourceBand> element in VRTRasterBand." );
          CPLFree( pszSrcDSName );
          return CE_Failure;
      }
  
      // Newly generated VRT will have RasterXSize, RasterYSize, DataType,
      // BlockXSize, BlockYSize tags, so that we don't have actually to
      // open the real dataset immediately, but we can use a proxy dataset
      // instead. This is particularly useful when dealing with huge VRT
      // For example, a VRT with the world coverage of DTED0 (25594 files).
      CPLXMLNode* psSrcProperties = CPLGetXMLNode(psSrc,"SourceProperties");
      int nRasterXSize = 0;
      int nRasterYSize = 0;
      // TODO(schwehr): What is the difference between 0 (GDT_Unknown) and -1?
      // Does there need to be a GDT_Uninitialized?
      GDALDataType eDataType = static_cast<GDALDataType>(-1);
      int nBlockXSize = 0;
      int nBlockYSize = 0;
      if( psSrcProperties )
      {
          nRasterXSize =
              atoi(CPLGetXMLValue(psSrcProperties, "RasterXSize", "0"));
          nRasterYSize =
              atoi(CPLGetXMLValue(psSrcProperties, "RasterYSize", "0"));
          const char *pszDataType =
              CPLGetXMLValue(psSrcProperties, "DataType", nullptr);
          if( pszDataType != nullptr )
          {
              for( int iType = 0; iType < GDT_TypeCount; iType++ )
              {
                  const char *pszThisName =
                      GDALGetDataTypeName(static_cast<GDALDataType>(iType));
  
                  if( pszThisName != nullptr && EQUAL(pszDataType, pszThisName) )
                  {
                      eDataType = static_cast<GDALDataType>(iType);
                      break;
                  }
              }
          }
          nBlockXSize = atoi(CPLGetXMLValue(psSrcProperties, "BlockXSize", "0"));
          nBlockYSize = atoi(CPLGetXMLValue(psSrcProperties, "BlockYSize", "0"));
          if( nRasterXSize < 0 || nRasterYSize < 0 ||
              nBlockXSize < 0 || nBlockYSize < 0 )
          {
              CPLError( CE_Warning, CPLE_AppDefined,
                        "Invalid <SourceProperties> element in VRTRasterBand." );
              CPLFree( pszSrcDSName );
              return CE_Failure;
          }
      }
  
      char** papszOpenOptions = GDALDeserializeOpenOptionsFromXML(psSrc);
      if( strstr(pszSrcDSName,"<VRTDataset") != nullptr )
          papszOpenOptions =
              CSLSetNameValue(papszOpenOptions, "ROOT_PATH", pszVRTPath);
  
      GDALDataset *poSrcDS = nullptr;
      if( nRasterXSize == 0 || nRasterYSize == 0 ||
          eDataType == static_cast<GDALDataType>(-1) ||
          nBlockXSize == 0 || nBlockYSize == 0 )
      {
          /* ----------------------------------------------------------------- */
          /*      Open the file (shared).                                      */
          /* ----------------------------------------------------------------- */
          int nOpenFlags = GDAL_OF_RASTER | GDAL_OF_VERBOSE_ERROR;
          if( bShared )
              nOpenFlags |= GDAL_OF_SHARED;
          poSrcDS = static_cast<GDALDataset *>( GDALOpenEx(
                      pszSrcDSName, nOpenFlags, nullptr,
                      (const char* const* )papszOpenOptions, nullptr ) );
      }
      else
      {
          /* ----------------------------------------------------------------- */
          /*      Create a proxy dataset                                       */
          /* ----------------------------------------------------------------- */
          CPLString osUniqueHandle( CPLSPrintf("%p", pUniqueHandle) );
          GDALProxyPoolDataset * const proxyDS =
              new GDALProxyPoolDataset( pszSrcDSName, nRasterXSize, nRasterYSize,
                                        GA_ReadOnly, bShared, nullptr, nullptr,
                                        osUniqueHandle.c_str() );
          proxyDS->SetOpenOptions(papszOpenOptions);
          poSrcDS = proxyDS;
  
          // Only the information of rasterBand nSrcBand will be accurate
          // but that's OK since we only use that band afterwards.
          for( int i = 1; i <= nSrcBand; i++ )
              proxyDS->AddSrcBandDescription(eDataType, nBlockXSize, nBlockYSize);
  
          if( bGetMaskBand )
          {
            GDALProxyPoolRasterBand *poMaskBand =
                dynamic_cast<GDALProxyPoolRasterBand *>(
                proxyDS->GetRasterBand(nSrcBand) );
            if( poMaskBand == nullptr )
            {
                CPLError(
                    CE_Fatal, CPLE_AssertionFailed, "dynamic_cast failed." );
            }
            else
            {
                poMaskBand->AddSrcMaskBandDescription(
                    eDataType, nBlockXSize, nBlockYSize );
            }
          }
      }
  
      CSLDestroy(papszOpenOptions);
  
      CPLFree( pszSrcDSName );
  
      if( poSrcDS == nullptr )
          return CE_Failure;
  
  /* -------------------------------------------------------------------- */
  /*      Get the raster band.                                            */
  /* -------------------------------------------------------------------- */
  
      m_poRasterBand = poSrcDS->GetRasterBand(nSrcBand);
      if( m_poRasterBand == nullptr )
      {
          if( poSrcDS->GetShared() )
              GDALClose( poSrcDS );
          return CE_Failure;
      }
      if( bGetMaskBand )
      {
          m_poMaskBandMainBand = m_poRasterBand;
          m_poRasterBand = m_poRasterBand->GetMaskBand();
          if( m_poRasterBand == nullptr )
              return CE_Failure;
      }
  
  /* -------------------------------------------------------------------- */
  /*      Set characteristics.                                            */
  /* -------------------------------------------------------------------- */
      CPLXMLNode * const psSrcRect = CPLGetXMLNode(psSrc,"SrcRect");
      if( psSrcRect )
      {
          SetSrcWindow( CPLAtof(CPLGetXMLValue(psSrcRect,"xOff","-1")),
                        CPLAtof(CPLGetXMLValue(psSrcRect,"yOff","-1")),
                        CPLAtof(CPLGetXMLValue(psSrcRect,"xSize","-1")),
                        CPLAtof(CPLGetXMLValue(psSrcRect,"ySize","-1")) );
      }
      else
      {
          m_dfSrcXOff = -1;
          m_dfSrcYOff = -1;
          m_dfSrcXSize = -1;
          m_dfSrcYSize = -1;
      }
  
      CPLXMLNode * const psDstRect = CPLGetXMLNode(psSrc,"DstRect");
      if( psDstRect )
      {
          SetDstWindow( CPLAtof(CPLGetXMLValue(psDstRect,"xOff","-1")),
                        CPLAtof(CPLGetXMLValue(psDstRect,"yOff","-1")),
                        CPLAtof(CPLGetXMLValue(psDstRect,"xSize","-1")),
                        CPLAtof(CPLGetXMLValue(psDstRect,"ySize","-1")) );
      }
      else
      {
        m_dfDstXOff = -1;
        m_dfDstYOff = -1;
        m_dfDstXSize = -1;
        m_dfDstYSize = -1;
      }
  
      return CE_None;
  }
  
  /************************************************************************/
  /*                             GetFileList()                            */
  /************************************************************************/
  
  void VRTSimpleSource::GetFileList( char*** ppapszFileList, int *pnSize,
                                     int *pnMaxSize, CPLHashSet* hSetFiles )
  {
      const char* pszFilename = nullptr;
      if( m_poRasterBand != nullptr && m_poRasterBand->GetDataset() != nullptr &&
          (pszFilename = m_poRasterBand->GetDataset()->GetDescription()) != nullptr )
      {
  /* -------------------------------------------------------------------- */
  /*      Is the filename even a real filesystem object?                  */
  /* -------------------------------------------------------------------- */
          if( strstr(pszFilename, "/vsicurl/http") != nullptr ||
              strstr(pszFilename, "/vsicurl/ftp") != nullptr )
          {
              // Testing the existence of remote resources can be excruciating
              // slow, so just suppose they exist.
          }
          else
          {
              VSIStatBufL  sStat;
              if( VSIStatExL( pszFilename, &sStat, VSI_STAT_EXISTS_FLAG ) != 0 )
                  return;
          }
  
  /* -------------------------------------------------------------------- */
  /*      Is it already in the list ?                                     */
  /* -------------------------------------------------------------------- */
          if( CPLHashSetLookup(hSetFiles, pszFilename) != nullptr )
              return;
  
  /* -------------------------------------------------------------------- */
  /*      Grow array if necessary                                         */
  /* -------------------------------------------------------------------- */
          if( *pnSize + 1 >= *pnMaxSize )
          {
              *pnMaxSize = std::max(*pnSize + 2, 2 + 2 * (*pnMaxSize));
              *ppapszFileList = static_cast<char **>( CPLRealloc(
                  *ppapszFileList, sizeof(char*) * (*pnMaxSize) ) );
          }
  
  /* -------------------------------------------------------------------- */
  /*      Add the string to the list                                      */
  /* -------------------------------------------------------------------- */
          (*ppapszFileList)[*pnSize] = CPLStrdup(pszFilename);
          (*ppapszFileList)[(*pnSize + 1)] = nullptr;
          CPLHashSetInsert(hSetFiles, (*ppapszFileList)[*pnSize]);
  
          (*pnSize) ++;
      }
  }
  
  /************************************************************************/
  /*                             GetBand()                                */
  /************************************************************************/
  
  GDALRasterBand* VRTSimpleSource::GetBand()
  {
      return m_poMaskBandMainBand ? nullptr : m_poRasterBand;
  }
  
  /************************************************************************/
  /*                       IsSameExceptBandNumber()                       */
  /************************************************************************/
  
  int VRTSimpleSource::IsSameExceptBandNumber( VRTSimpleSource* poOtherSource )
  {
      return m_dfSrcXOff == poOtherSource->m_dfSrcXOff &&
             m_dfSrcYOff == poOtherSource->m_dfSrcYOff &&
             m_dfSrcXSize == poOtherSource->m_dfSrcXSize &&
             m_dfSrcYSize == poOtherSource->m_dfSrcYSize &&
             m_dfDstXOff == poOtherSource->m_dfDstXOff &&
             m_dfDstYOff == poOtherSource->m_dfDstYOff &&
             m_dfDstXSize == poOtherSource->m_dfDstXSize &&
             m_dfDstYSize == poOtherSource->m_dfDstYSize &&
             m_bNoDataSet == poOtherSource->m_bNoDataSet &&
             m_dfNoDataValue == poOtherSource->m_dfNoDataValue &&
             GetBand() != nullptr && poOtherSource->GetBand() != nullptr &&
             GetBand()->GetDataset() != nullptr &&
             poOtherSource->GetBand()->GetDataset() != nullptr &&
             EQUAL(GetBand()->GetDataset()->GetDescription(),
                   poOtherSource->GetBand()->GetDataset()->GetDescription());
  }
  
  /************************************************************************/
  /*                              SrcToDst()                              */
  /*                                                                      */
  /*      Note: this is a no-op if the dst window is -1,-1,-1,-1.         */
  /************************************************************************/
  
  void VRTSimpleSource::SrcToDst( double dfX, double dfY,
                                  double &dfXOut, double &dfYOut ) const
  
  {
      dfXOut = ((dfX - m_dfSrcXOff) / m_dfSrcXSize) * m_dfDstXSize + m_dfDstXOff;
      dfYOut = ((dfY - m_dfSrcYOff) / m_dfSrcYSize) * m_dfDstYSize + m_dfDstYOff;
  }
  
  /************************************************************************/
  /*                              DstToSrc()                              */
  /*                                                                      */
  /*      Note: this is a no-op if the dst window is -1,-1,-1,-1.         */
  /************************************************************************/
  
  void VRTSimpleSource::DstToSrc( double dfX, double dfY,
                                  double &dfXOut, double &dfYOut ) const
  
  {
      dfXOut = ((dfX - m_dfDstXOff) / m_dfDstXSize) * m_dfSrcXSize + m_dfSrcXOff;
      dfYOut = ((dfY - m_dfDstYOff) / m_dfDstYSize) * m_dfSrcYSize + m_dfSrcYOff;
  }
  
  /************************************************************************/
  /*                          GetSrcDstWindow()                           */
  /************************************************************************/
  
  int
  VRTSimpleSource::GetSrcDstWindow( int nXOff, int nYOff, int nXSize, int nYSize,
                                    int nBufXSize, int nBufYSize,
                                    double *pdfReqXOff, double *pdfReqYOff,
                                    double *pdfReqXSize, double *pdfReqYSize,
                                    int *pnReqXOff, int *pnReqYOff,
                                    int *pnReqXSize, int *pnReqYSize,
                                    int *pnOutXOff, int *pnOutYOff,
                                    int *pnOutXSize, int *pnOutYSize )
  
  {
      if( m_dfSrcXSize == 0.0 || m_dfSrcYSize == 0.0 ||
          m_dfDstXSize == 0.0 || m_dfDstYSize == 0.0 )
      {
          return FALSE;
      }
  
      const bool bDstWinSet = m_dfDstXOff != -1 || m_dfDstXSize != -1
          || m_dfDstYOff != -1 || m_dfDstYSize != -1;
  
  #ifdef DEBUG
      const bool bSrcWinSet = m_dfSrcXOff != -1 || m_dfSrcXSize != -1
          || m_dfSrcYOff != -1 || m_dfSrcYSize != -1;
  
      if( bSrcWinSet != bDstWinSet )
      {
          return FALSE;
      }
  #endif
  
  /* -------------------------------------------------------------------- */
  /*      If the input window completely misses the portion of the        */
  /*      virtual dataset provided by this source we have nothing to do.  */
  /* -------------------------------------------------------------------- */
      if( bDstWinSet )
      {
          if( nXOff >= m_dfDstXOff + m_dfDstXSize
              || nYOff >= m_dfDstYOff + m_dfDstYSize
              || nXOff + nXSize < m_dfDstXOff
              || nYOff + nYSize < m_dfDstYOff )
              return FALSE;
      }
  
  /* -------------------------------------------------------------------- */
  /*      This request window corresponds to the whole output buffer.     */
  /* -------------------------------------------------------------------- */
      *pnOutXOff = 0;
      *pnOutYOff = 0;
      *pnOutXSize = nBufXSize;
      *pnOutYSize = nBufYSize;
  
  /* -------------------------------------------------------------------- */
  /*      If the input window extents outside the portion of the on       */
  /*      the virtual file that this source can set, then clip down       */
  /*      the requested window.                                           */
  /* -------------------------------------------------------------------- */
      bool bModifiedX = false;
      bool bModifiedY = false;
      double dfRXOff = nXOff;
      double dfRYOff = nYOff;
      double dfRXSize = nXSize;
      double dfRYSize = nYSize;
  
      if( bDstWinSet )
      {
          if( dfRXOff < m_dfDstXOff )
          {
              dfRXSize = dfRXSize + dfRXOff - m_dfDstXOff;
              dfRXOff = m_dfDstXOff;
              bModifiedX = true;
          }
  
          if( dfRYOff < m_dfDstYOff )
          {
              dfRYSize = dfRYSize + dfRYOff - m_dfDstYOff;
              dfRYOff = m_dfDstYOff;
              bModifiedY = true;
          }
  
          if( dfRXOff + dfRXSize > m_dfDstXOff + m_dfDstXSize )
          {
              dfRXSize = m_dfDstXOff + m_dfDstXSize - dfRXOff;
              bModifiedX = true;
          }
  
          if( dfRYOff + dfRYSize > m_dfDstYOff + m_dfDstYSize )
          {
              dfRYSize = m_dfDstYOff + m_dfDstYSize - dfRYOff;
              bModifiedY = true;
          }
      }
  
  /* -------------------------------------------------------------------- */
  /*      Translate requested region in virtual file into the source      */
  /*      band coordinates.                                               */
  /* -------------------------------------------------------------------- */
      const double dfScaleX = m_dfSrcXSize / m_dfDstXSize;
      const double dfScaleY = m_dfSrcYSize / m_dfDstYSize;
  
      *pdfReqXOff = (dfRXOff - m_dfDstXOff) * dfScaleX + m_dfSrcXOff;
      *pdfReqYOff = (dfRYOff - m_dfDstYOff) * dfScaleY + m_dfSrcYOff;
      *pdfReqXSize = dfRXSize * dfScaleX;
      *pdfReqYSize = dfRYSize * dfScaleY;
  
      if( !CPLIsFinite(*pdfReqXOff) ||
          !CPLIsFinite(*pdfReqYOff) ||
          !CPLIsFinite(*pdfReqXSize) ||
          !CPLIsFinite(*pdfReqYSize) ||
          *pdfReqXOff > INT_MAX ||
          *pdfReqYOff > INT_MAX ||
          *pdfReqXSize < 0 ||
          *pdfReqYSize < 0 )
      {
          return FALSE;
      }
  
  /* -------------------------------------------------------------------- */
  /*      Clamp within the bounds of the available source data.           */
  /* -------------------------------------------------------------------- */
      if( *pdfReqXOff < 0 )
      {
          *pdfReqXSize += *pdfReqXOff;
          *pdfReqXOff = 0;
          bModifiedX = true;
      }
      if( *pdfReqYOff < 0 )
      {
          *pdfReqYSize += *pdfReqYOff;
          *pdfReqYOff = 0;
          bModifiedY = true;
      }
  
      *pnReqXOff = static_cast<int>( floor(*pdfReqXOff) );
      *pnReqYOff = static_cast<int>( floor(*pdfReqYOff) );
  
      if( *pdfReqXSize > INT_MAX )
          *pnReqXSize = INT_MAX;
      else
          *pnReqXSize = static_cast<int>( floor(*pdfReqXSize + 0.5) );
  
      if( *pdfReqYSize > INT_MAX )
          *pnReqYSize = INT_MAX;
      else
          *pnReqYSize = static_cast<int>( floor(*pdfReqYSize + 0.5) );
  
  /* -------------------------------------------------------------------- */
  /*      Clamp within the bounds of the available source data.           */
  /* -------------------------------------------------------------------- */
  
      if( *pnReqXSize == 0 )
          *pnReqXSize = 1;
      if( *pnReqYSize == 0 )
          *pnReqYSize = 1;
  
      if( *pnReqXSize > INT_MAX - *pnReqXOff ||
          *pnReqXOff + *pnReqXSize > m_poRasterBand->GetXSize() )
      {
          *pnReqXSize = m_poRasterBand->GetXSize() - *pnReqXOff;
          bModifiedX = true;
      }
      if( *pdfReqXOff + *pdfReqXSize > m_poRasterBand->GetXSize() )
      {
          *pdfReqXSize = m_poRasterBand->GetXSize() - *pdfReqXOff;
          bModifiedX = true;
      }
  
      if( *pnReqYSize > INT_MAX - *pnReqYOff ||
          *pnReqYOff + *pnReqYSize > m_poRasterBand->GetYSize() )
      {
          *pnReqYSize = m_poRasterBand->GetYSize() - *pnReqYOff;
          bModifiedY = true;
      }
      if( *pdfReqYOff + *pdfReqYSize > m_poRasterBand->GetYSize() )
      {
          *pdfReqYSize = m_poRasterBand->GetYSize() - *pdfReqYOff;
          bModifiedY = true;
      }
  
  /* -------------------------------------------------------------------- */
  /*      Don't do anything if the requesting region is completely off    */
  /*      the source image.                                               */
  /* -------------------------------------------------------------------- */
      if( *pnReqXOff >= m_poRasterBand->GetXSize()
          || *pnReqYOff >= m_poRasterBand->GetYSize()
          || *pnReqXSize <= 0 || *pnReqYSize <= 0 )
      {
          return FALSE;
      }
  
  /* -------------------------------------------------------------------- */
  /*      If we haven't had to modify the source rectangle, then the      */
  /*      destination rectangle must be the whole region.                 */
  /* -------------------------------------------------------------------- */
      if( !bModifiedX && !bModifiedY )
          return TRUE;
  
  /* -------------------------------------------------------------------- */
  /*      Now transform this possibly reduced request back into the       */
  /*      destination buffer coordinates in case the output region is     */
  /*      less than the whole buffer.                                     */
  /* -------------------------------------------------------------------- */
      double dfDstULX = 0.0;
      double dfDstULY = 0.0;
      double dfDstLRX = 0.0;
      double dfDstLRY = 0.0;
  
      SrcToDst( *pdfReqXOff, *pdfReqYOff, dfDstULX, dfDstULY );
      SrcToDst( *pdfReqXOff + *pdfReqXSize, *pdfReqYOff + *pdfReqYSize,
                dfDstLRX, dfDstLRY );
  #if DEBUG_VERBOSE
      CPLDebug( "VRT", "dfDstULX=%g dfDstULY=%g dfDstLRX=%g dfDstLRY=%g",
                dfDstULX, dfDstULY, dfDstLRX, dfDstLRY );
  #endif
  
      if( bModifiedX )
      {
          const double dfScaleWinToBufX =
              nBufXSize / static_cast<double>( nXSize );
  
          const double dfOutXOff = (dfDstULX - nXOff) * dfScaleWinToBufX+0.001;
          if( dfOutXOff <= 0 )
              *pnOutXOff = 0;
          else if( dfOutXOff > INT_MAX )
              *pnOutXOff = INT_MAX;
          else
              *pnOutXOff = (int) (dfOutXOff);
          double dfOutRightXOff = (dfDstLRX - nXOff) * dfScaleWinToBufX-0.001;
          if( dfOutRightXOff < dfOutXOff )
              return FALSE;
          if( dfOutRightXOff > INT_MAX )
              dfOutRightXOff = INT_MAX;
          *pnOutXSize = (int) ceil(dfOutRightXOff) - *pnOutXOff;
  
          *pnOutXOff = std::max(0, *pnOutXOff);
          if( *pnOutXSize > INT_MAX - *pnOutXOff ||
              *pnOutXOff + *pnOutXSize > nBufXSize )
              *pnOutXSize = nBufXSize - *pnOutXOff;
      }
  
      if( bModifiedY )
      {
          const double dfScaleWinToBufY =
              nBufYSize / static_cast<double>( nYSize );
  
          const double dfOutYOff = (dfDstULY - nYOff) * dfScaleWinToBufY+0.001;
          if( dfOutYOff <= 0 )
              *pnOutYOff = 0;
          else if( dfOutYOff > INT_MAX )
              *pnOutYOff = INT_MAX;
          else
              *pnOutYOff = static_cast<int>(dfOutYOff);
          *pnOutYOff = (int) (dfOutYOff);
          double dfOutTopYOff = (dfDstLRY - nYOff) * dfScaleWinToBufY-0.001;
          if( dfOutTopYOff < dfOutYOff )
              return FALSE;
          if( dfOutTopYOff > INT_MAX )
              dfOutTopYOff = INT_MAX;
          *pnOutYSize = static_cast<int>( ceil(dfOutTopYOff) ) - *pnOutYOff;
  
          *pnOutYOff = std::max(0, *pnOutYOff);
          if( *pnOutYSize > INT_MAX - *pnOutYOff ||
              *pnOutYOff + *pnOutYSize > nBufYSize )
              *pnOutYSize = nBufYSize - *pnOutYOff;
      }
  
      if( *pnOutXSize < 1 || *pnOutYSize < 1 )
          return FALSE;
  
      return TRUE;
  }
  
  /************************************************************************/
  /*                          NeedMaxValAdjustment()                      */
  /************************************************************************/
  
  int VRTSimpleSource::NeedMaxValAdjustment() const
  {
      if( !m_nMaxValue )
          return FALSE;
  
      const char* pszNBITS =
          m_poRasterBand->GetMetadataItem("NBITS", "IMAGE_STRUCTURE");
      const int nBits = (pszNBITS) ? atoi(pszNBITS) : 0;
      if( nBits >= 1 && nBits <= 31 )
      {
          const int nBandMaxValue = static_cast<int>((1U << nBits) - 1);
          return nBandMaxValue > m_nMaxValue;
      }
      return TRUE;
  }
  
  /************************************************************************/
  /*                              RasterIO()                              */
  /************************************************************************/
  
  CPLErr
  VRTSimpleSource::RasterIO( GDALDataType eBandDataType,
                             int nXOff, int nYOff, int nXSize, int nYSize,
                             void *pData, int nBufXSize, int nBufYSize,
                             GDALDataType eBufType,
                             GSpacing nPixelSpace,
                             GSpacing nLineSpace,
                             GDALRasterIOExtraArg* psExtraArgIn )
  
  {
      GDALRasterIOExtraArg sExtraArg;
      INIT_RASTERIO_EXTRA_ARG(sExtraArg);
      GDALRasterIOExtraArg* psExtraArg = &sExtraArg;
  
      // The window we will actually request from the source raster band.
      double dfReqXOff = 0.0;
      double dfReqYOff = 0.0;
      double dfReqXSize = 0.0;
      double dfReqYSize = 0.0;
      int nReqXOff = 0;
      int nReqYOff = 0;
      int nReqXSize = 0;
      int nReqYSize = 0;
  
      // The window we will actual set _within_ the pData buffer.
      int nOutXOff = 0;
      int nOutYOff = 0;
      int nOutXSize = 0;
      int nOutYSize = 0;
  
      if( !GetSrcDstWindow( nXOff, nYOff, nXSize, nYSize,
                            nBufXSize, nBufYSize,
                            &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
                            &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
                            &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize ) )
      {
          return CE_None;
      }
  #if DEBUG_VERBOSE
      CPLDebug(
          "VRT",
          "nXOff=%d, nYOff=%d, nXSize=%d, nYSize=%d, nBufXSize=%d, nBufYSize=%d,\n"
          "dfReqXOff=%g, dfReqYOff=%g, dfReqXSize=%g, dfReqYSize=%g,\n"
          "nReqXOff=%d, nReqYOff=%d, nReqXSize=%d, nReqYSize=%d,\n"
          "nOutXOff=%d, nOutYOff=%d, nOutXSize=%d, nOutYSize=%d",
          nXOff, nYOff, nXSize, nYSize,
          nBufXSize, nBufYSize,
          dfReqXOff, dfReqYOff, dfReqXSize, dfReqYSize,
          nReqXOff, nReqYOff, nReqXSize, nReqYSize,
          nOutXOff, nOutYOff, nOutXSize, nOutYSize );
  #endif
  
  /* -------------------------------------------------------------------- */
  /*      Actually perform the IO request.                                */
  /* -------------------------------------------------------------------- */
      if( !m_osResampling.empty() )
      {
          psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
      }
      else if( psExtraArgIn != nullptr )
      {
          psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
      }
      psExtraArg->bFloatingPointWindowValidity = TRUE;
      psExtraArg->dfXOff = dfReqXOff;
      psExtraArg->dfYOff = dfReqYOff;
      psExtraArg->dfXSize = dfReqXSize;
      psExtraArg->dfYSize = dfReqYSize;
  
      GByte* pabyOut =
          static_cast<unsigned char *>(pData)
          + nOutXOff * nPixelSpace
          + static_cast<GPtrDiff_t>(nOutYOff) * nLineSpace;
  
      CPLErr eErr = CE_Failure;
  
      if( GDALDataTypeIsConversionLossy(m_poRasterBand->GetRasterDataType(),
                                        eBandDataType) )
      {
          const int nBandDTSize = GDALGetDataTypeSizeBytes(eBandDataType);
          void* pTemp = VSI_MALLOC3_VERBOSE(nOutXSize, nOutYSize, nBandDTSize);
          if( pTemp )
          {
              eErr =
                  m_poRasterBand->RasterIO(
                      GF_Read,
                      nReqXOff, nReqYOff, nReqXSize, nReqYSize,
                      pTemp,
                      nOutXSize, nOutYSize,
                      eBandDataType, 0, 0, psExtraArg );
              if( eErr == CE_None )
              {
                  GByte* pabyTemp = reinterpret_cast<GByte*>(pTemp);
                  for( int iY = 0; iY < nOutYSize; iY++ )
                  {
                      GDALCopyWords(pabyTemp + static_cast<size_t>(iY) *
                                                      nBandDTSize * nOutXSize,
                                    eBandDataType, nBandDTSize,
                                    pabyOut +
                                      static_cast<size_t>(iY) * nLineSpace,
                                    eBufType,
                                    static_cast<int>(nPixelSpace),
                                    nOutXSize);
                  }
              }
              VSIFree(pTemp);
          }
      }
      else
      {
          eErr =
              m_poRasterBand->RasterIO(
                  GF_Read,
                  nReqXOff, nReqYOff, nReqXSize, nReqYSize,
                  pabyOut,
                  nOutXSize, nOutYSize,
                  eBufType, nPixelSpace, nLineSpace, psExtraArg );
      }
  
      if( NeedMaxValAdjustment() )
      {
          for( int j = 0; j < nOutYSize; j++ )
          {
              for( int i = 0; i < nOutXSize; i++ )
              {
                  int nVal = 0;
                  GDALCopyWords( pabyOut + j * nLineSpace + i * nPixelSpace,
                                 eBufType, 0,
                                 &nVal, GDT_Int32, 0,
                                 1 );
                  if( nVal > m_nMaxValue )
                      nVal = m_nMaxValue;
                  GDALCopyWords( &nVal, GDT_Int32, 0,
                                 pabyOut + j * nLineSpace + i * nPixelSpace,
                                 eBufType, 0,
                                 1);
              }
          }
      }
  
      return eErr;
  }
  
  /************************************************************************/
  /*                             GetMinimum()                             */
  /************************************************************************/
  
  double VRTSimpleSource::GetMinimum( int nXSize, int nYSize, int *pbSuccess )
  {
      // The window we will actually request from the source raster band.
      double dfReqXOff = 0.0;
      double dfReqYOff = 0.0;
      double dfReqXSize = 0.0;
      double dfReqYSize = 0.0;
      int nReqXOff = 0;
      int nReqYOff = 0;
      int nReqXSize = 0;
      int nReqYSize = 0;
  
      // The window we will actual set _within_ the pData buffer.
      int nOutXOff = 0;
      int nOutYOff = 0;
      int nOutXSize = 0;
      int nOutYSize = 0;
  
      if( !GetSrcDstWindow( 0, 0, nXSize, nYSize,
                            nXSize, nYSize,
                            &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
                            &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
                            &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize ) ||
          nReqXOff != 0 || nReqYOff != 0 ||
          nReqXSize != m_poRasterBand->GetXSize() ||
          nReqYSize != m_poRasterBand->GetYSize())
      {
          *pbSuccess = FALSE;
          return 0;
      }
  
      const double dfVal = m_poRasterBand->GetMinimum(pbSuccess);
      if( NeedMaxValAdjustment() && dfVal > m_nMaxValue )
          return m_nMaxValue;
      return dfVal;
  }
  
  /************************************************************************/
  /*                             GetMaximum()                             */
  /************************************************************************/
  
  double VRTSimpleSource::GetMaximum( int nXSize, int nYSize, int *pbSuccess )
  {
      // The window we will actually request from the source raster band.
      double dfReqXOff = 0.0;
      double dfReqYOff = 0.0;
      double dfReqXSize = 0.0;
      double dfReqYSize = 0.0;
      int nReqXOff = 0;
      int nReqYOff = 0;
      int nReqXSize = 0;
      int nReqYSize = 0;
  
      // The window we will actual set _within_ the pData buffer.
      int nOutXOff = 0;
      int nOutYOff = 0;
      int nOutXSize = 0;
      int nOutYSize = 0;
  
      if( !GetSrcDstWindow( 0, 0, nXSize, nYSize,
                            nXSize, nYSize,
                            &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
                            &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
                            &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize ) ||
          nReqXOff != 0 || nReqYOff != 0 ||
          nReqXSize != m_poRasterBand->GetXSize() ||
          nReqYSize != m_poRasterBand->GetYSize())
      {
          *pbSuccess = FALSE;
          return 0;
      }
  
      const double dfVal = m_poRasterBand->GetMaximum(pbSuccess);
      if( NeedMaxValAdjustment() && dfVal > m_nMaxValue )
          return m_nMaxValue;
      return dfVal;
  }
  
  /************************************************************************/
  /*                       ComputeRasterMinMax()                          */
  /************************************************************************/
  
  CPLErr VRTSimpleSource::ComputeRasterMinMax( int nXSize, int nYSize,
                                               int bApproxOK, double* adfMinMax )
  {
      // The window we will actually request from the source raster band.
      double dfReqXOff = 0.0;
      double dfReqYOff = 0.0;
      double dfReqXSize = 0.0;
      double dfReqYSize = 0.0;
      int nReqXOff = 0;
      int nReqYOff = 0;
      int nReqXSize = 0;
      int nReqYSize = 0;
  
      // The window we will actual set _within_ the pData buffer.
      int nOutXOff = 0;
      int nOutYOff = 0;
      int nOutXSize = 0;
      int nOutYSize = 0;
  
      if( !GetSrcDstWindow( 0, 0, nXSize, nYSize,
                            nXSize, nYSize,
                            &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
                            &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
                            &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize ) ||
          nReqXOff != 0 || nReqYOff != 0 ||
          nReqXSize != m_poRasterBand->GetXSize() ||
          nReqYSize != m_poRasterBand->GetYSize())
      {
          return CE_Failure;
      }
  
      const CPLErr eErr =
          m_poRasterBand->ComputeRasterMinMax( bApproxOK, adfMinMax );
      if( NeedMaxValAdjustment() )
      {
          if( adfMinMax[0] > m_nMaxValue )
              adfMinMax[0] = m_nMaxValue;
          if( adfMinMax[1] > m_nMaxValue )
              adfMinMax[1] = m_nMaxValue;
      }
      return eErr;
  }
  
  /************************************************************************/
  /*                         ComputeStatistics()                          */
  /************************************************************************/
  
  CPLErr VRTSimpleSource::ComputeStatistics(
      int nXSize, int nYSize,
      int bApproxOK,
      double *pdfMin, double *pdfMax,
      double *pdfMean, double *pdfStdDev,
      GDALProgressFunc pfnProgress, void *pProgressData )
  {
      // The window we will actually request from the source raster band.
      double dfReqXOff = 0.0;
      double dfReqYOff = 0.0;
      double dfReqXSize = 0.0;
      double dfReqYSize = 0.0;
      int nReqXOff = 0;
      int nReqYOff = 0;
      int nReqXSize = 0;
      int nReqYSize = 0;
  
      // The window we will actual set _within_ the pData buffer.
      int nOutXOff = 0;
      int nOutYOff = 0;
      int nOutXSize = 0;
      int nOutYSize = 0;
  
      if( NeedMaxValAdjustment() ||
          !GetSrcDstWindow( 0, 0, nXSize, nYSize,
                            nXSize, nYSize,
                            &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
                            &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
                            &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize ) ||
          nReqXOff != 0 || nReqYOff != 0 ||
          nReqXSize != m_poRasterBand->GetXSize() ||
          nReqYSize != m_poRasterBand->GetYSize())
      {
          return CE_Failure;
      }
  
      return m_poRasterBand->ComputeStatistics( bApproxOK, pdfMin, pdfMax,
                                                pdfMean, pdfStdDev,
                                                pfnProgress, pProgressData );
  }
  
  /************************************************************************/
  /*                            GetHistogram()                            */
  /************************************************************************/
  
  CPLErr VRTSimpleSource::GetHistogram(
      int nXSize, int nYSize,
      double dfMin, double dfMax,
      int nBuckets, GUIntBig * panHistogram,
      int bIncludeOutOfRange, int bApproxOK,
      GDALProgressFunc pfnProgress, void *pProgressData )
  {
      // The window we will actually request from the source raster band.
      double dfReqXOff = 0.0;
      double dfReqYOff = 0.0;
      double dfReqXSize = 0.0;
      double dfReqYSize = 0.0;
      int nReqXOff = 0;
      int nReqYOff = 0;
      int nReqXSize = 0;
      int nReqYSize = 0;
  
      // The window we will actual set _within_ the pData buffer.
      int nOutXOff = 0;
      int nOutYOff = 0;
      int nOutXSize = 0;
      int nOutYSize = 0;
  
      if( NeedMaxValAdjustment() ||
          !GetSrcDstWindow( 0, 0, nXSize, nYSize,
                            nXSize, nYSize,
                            &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
                            &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
                            &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize ) ||
          nReqXOff != 0 || nReqYOff != 0 ||
          nReqXSize != m_poRasterBand->GetXSize() ||
          nReqYSize != m_poRasterBand->GetYSize())
      {
          return CE_Failure;
      }
  
      return m_poRasterBand->GetHistogram( dfMin, dfMax, nBuckets,
                                           panHistogram,
                                           bIncludeOutOfRange, bApproxOK,
                                           pfnProgress, pProgressData );
  }
  
  /************************************************************************/
  /*                          DatasetRasterIO()                           */
  /************************************************************************/
  
  CPLErr VRTSimpleSource::DatasetRasterIO(
      GDALDataType eBandDataType,
      int nXOff, int nYOff, int nXSize, int nYSize,
      void * pData, int nBufXSize, int nBufYSize,
      GDALDataType eBufType,
      int nBandCount, int *panBandMap,
      GSpacing nPixelSpace, GSpacing nLineSpace,
      GSpacing nBandSpace,
      GDALRasterIOExtraArg* psExtraArgIn )
  {
      if( !EQUAL(GetType(), "SimpleSource") )
      {
          CPLError(CE_Failure, CPLE_NotSupported,
                   "DatasetRasterIO() not implemented for %s", GetType());
          return CE_Failure;
      }
  
      GDALRasterIOExtraArg sExtraArg;
      INIT_RASTERIO_EXTRA_ARG(sExtraArg);
      GDALRasterIOExtraArg* psExtraArg = &sExtraArg;
  
      // The window we will actually request from the source raster band.
      double dfReqXOff = 0.0;
      double dfReqYOff = 0.0;
      double dfReqXSize = 0.0;
      double dfReqYSize = 0.0;
      int nReqXOff = 0;
      int nReqYOff = 0;
      int nReqXSize = 0;
      int nReqYSize = 0;
  
      // The window we will actual set _within_ the pData buffer.
      int nOutXOff = 0;
      int nOutYOff = 0;
      int nOutXSize = 0;
      int nOutYSize = 0;
  
      if( !GetSrcDstWindow( nXOff, nYOff, nXSize, nYSize,
                            nBufXSize, nBufYSize,
                            &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
                            &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
                            &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize ) )
      {
          return CE_None;
      }
  
      GDALDataset* poDS = m_poRasterBand->GetDataset();
      if( poDS == nullptr )
          return CE_Failure;
  
      if( !m_osResampling.empty() )
      {
          psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
      }
      else if( psExtraArgIn != nullptr )
      {
          psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
      }
      psExtraArg->bFloatingPointWindowValidity = TRUE;
      psExtraArg->dfXOff = dfReqXOff;
      psExtraArg->dfYOff = dfReqYOff;
      psExtraArg->dfXSize = dfReqXSize;
      psExtraArg->dfYSize = dfReqYSize;
  
      GByte* pabyOut =
          static_cast<unsigned char *>(pData)
          + nOutXOff * nPixelSpace
          + (GPtrDiff_t)nOutYOff * nLineSpace;
  
      CPLErr eErr = CE_Failure;
  
      if( GDALDataTypeIsConversionLossy(m_poRasterBand->GetRasterDataType(),
                                        eBandDataType) )
      {
          const int nBandDTSize = GDALGetDataTypeSizeBytes(eBandDataType);
          void* pTemp = VSI_MALLOC3_VERBOSE(nOutXSize, nOutYSize,
                                            nBandDTSize * nBandCount);
          if( pTemp )
          {
              eErr = poDS->RasterIO(
                  GF_Read,
                  nReqXOff, nReqYOff, nReqXSize, nReqYSize,
                  pTemp,
                  nOutXSize, nOutYSize,
                  eBandDataType, nBandCount, panBandMap,
                  0, 0, 0, psExtraArg );
              if( eErr == CE_None )
              {
                  GByte* pabyTemp = reinterpret_cast<GByte*>(pTemp);
                  const size_t nSrcBandSpace = static_cast<size_t>(nOutYSize) *
                                                  nOutXSize * nBandDTSize;
                  for( int iBand = 0; iBand < nBandCount; iBand ++ )
                  {
                      for( int iY = 0; iY < nOutYSize; iY++ )
                      {
                          GDALCopyWords(pabyTemp + iBand * nSrcBandSpace +
                                          static_cast<size_t>(iY) * nBandDTSize * nOutXSize,
                                      eBandDataType, nBandDTSize,
                                      pabyOut +
                                          static_cast<size_t>(iY) * nLineSpace +
                                          static_cast<size_t>(iBand) * nBandSpace,
                                      eBufType, static_cast<int>(nPixelSpace),
                                      nOutXSize);
                      }
                  }
              }
              VSIFree(pTemp);
          }
      }
      else
      {
          eErr = poDS->RasterIO(
              GF_Read,
              nReqXOff, nReqYOff, nReqXSize, nReqYSize,
              pabyOut,
              nOutXSize, nOutYSize,
              eBufType, nBandCount, panBandMap,
              nPixelSpace, nLineSpace, nBandSpace, psExtraArg );
      }
  
      if( NeedMaxValAdjustment() )
      {
          for( int k = 0; k < nBandCount; k++ )
          {
              for( int j = 0; j < nOutYSize; j++ )
              {
                  for( int i = 0; i < nOutXSize; i++ )
                  {
                      int nVal = 0;
                      GDALCopyWords(
                          pabyOut + k * nBandSpace + j * nLineSpace +
                              i * nPixelSpace,
                          eBufType, 0,
                          &nVal, GDT_Int32, 0,
                          1 );
  
                      if( nVal > m_nMaxValue )
                          nVal = m_nMaxValue;
  
                      GDALCopyWords(
                          &nVal, GDT_Int32, 0,
                          pabyOut + k * nBandSpace + j * nLineSpace +
                              i * nPixelSpace,
                          eBufType, 0,
                          1 );
                  }
              }
          }
      }
  
      return eErr;
  }
  
  /************************************************************************/
  /*                          SetResampling()                             */
  /************************************************************************/
  
  void VRTSimpleSource::SetResampling( const char* pszResampling )
  {
      m_osResampling = (pszResampling) ? pszResampling : "";
  }
  
  /************************************************************************/
  /* ==================================================================== */
  /*                         VRTAveragedSource                            */
  /* ==================================================================== */
  /************************************************************************/
  
  /************************************************************************/
  /*                         VRTAveragedSource()                          */
  /************************************************************************/
  
  VRTAveragedSource::VRTAveragedSource() {}
  
  /************************************************************************/
  /*                           SerializeToXML()                           */
  /************************************************************************/
  
  CPLXMLNode *VRTAveragedSource::SerializeToXML( const char *pszVRTPath )
  
  {
      CPLXMLNode * const psSrc = VRTSimpleSource::SerializeToXML( pszVRTPath );
  
      if( psSrc == nullptr )
          return nullptr;
  
      CPLFree( psSrc->pszValue );
      psSrc->pszValue = CPLStrdup( "AveragedSource" );
  
      return psSrc;
  }
  
  /************************************************************************/
  /*                              RasterIO()                              */
  /************************************************************************/
  
  CPLErr
  VRTAveragedSource::RasterIO( GDALDataType /*eBandDataType*/,
                               int nXOff, int nYOff, int nXSize, int nYSize,
                               void *pData, int nBufXSize, int nBufYSize,
                               GDALDataType eBufType,
                               GSpacing nPixelSpace,
                               GSpacing nLineSpace,
                               GDALRasterIOExtraArg* psExtraArgIn )
  
  {
      GDALRasterIOExtraArg sExtraArg;
      INIT_RASTERIO_EXTRA_ARG(sExtraArg);
      GDALRasterIOExtraArg* psExtraArg = &sExtraArg;
  
      // The window we will actually request from the source raster band.
      double dfReqXOff = 0.0;
      double dfReqYOff = 0.0;
      double dfReqXSize = 0.0;
      double dfReqYSize = 0.0;
      int nReqXOff = 0;
      int nReqYOff = 0;
      int nReqXSize = 0;
      int nReqYSize = 0;
  
      // The window we will actual set _within_ the pData buffer.
      int nOutXOff = 0;
      int nOutYOff = 0;
      int nOutXSize = 0;
      int nOutYSize = 0;
  
      if( !GetSrcDstWindow( nXOff, nYOff, nXSize, nYSize, nBufXSize, nBufYSize,
                            &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
                            &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
                            &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize ) )
          return CE_None;
  
  /* -------------------------------------------------------------------- */
  /*      Allocate a temporary buffer to whole the full resolution        */
  /*      data from the area of interest.                                 */
  /* -------------------------------------------------------------------- */
      float * const pafSrc = static_cast<float *>(
          VSI_MALLOC3_VERBOSE(sizeof(float), nReqXSize, nReqYSize) );
      if( pafSrc == nullptr )
      {
          return CE_Failure;
      }
  
  /* -------------------------------------------------------------------- */
  /*      Load it.                                                        */
  /* -------------------------------------------------------------------- */
      if( !m_osResampling.empty() )
      {
          psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
      }
      else if( psExtraArgIn != nullptr )
      {
          psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
      }
  
      psExtraArg->bFloatingPointWindowValidity = TRUE;
      psExtraArg->dfXOff = dfReqXOff;
      psExtraArg->dfYOff = dfReqYOff;
      psExtraArg->dfXSize = dfReqXSize;
      psExtraArg->dfYSize = dfReqYSize;
  
      const CPLErr eErr =
          m_poRasterBand->RasterIO( GF_Read,
                                    nReqXOff, nReqYOff, nReqXSize, nReqYSize,
                                    pafSrc, nReqXSize, nReqYSize, GDT_Float32,
                                    0, 0, psExtraArg );
  
      if( eErr != CE_None )
      {
          VSIFree( pafSrc );
          return eErr;
      }
  
  /* -------------------------------------------------------------------- */
  /*      Do the averaging.                                               */
  /* -------------------------------------------------------------------- */
      for( int iBufLine = nOutYOff; iBufLine < nOutYOff + nOutYSize; iBufLine++ )
      {
          const double dfYDst =
              (iBufLine / static_cast<double>(nBufYSize) ) * nYSize + nYOff;
  
          for( int iBufPixel = nOutXOff;
               iBufPixel < nOutXOff + nOutXSize;
               iBufPixel++ )
          {
              double dfXSrcStart, dfXSrcEnd, dfYSrcStart, dfYSrcEnd;
              int    iXSrcStart, iYSrcStart, iXSrcEnd, iYSrcEnd;
  
              const double dfXDst =
                  (iBufPixel / static_cast<double>(nBufXSize)) * nXSize + nXOff;
  
              // Compute the source image rectangle needed for this pixel.
              DstToSrc( dfXDst, dfYDst, dfXSrcStart, dfYSrcStart );
              DstToSrc( dfXDst+1.0, dfYDst+1.0, dfXSrcEnd, dfYSrcEnd );
  
              // Convert to integers, assuming that the center of the source
              // pixel must be in our rect to get included.
              if( dfXSrcEnd >= dfXSrcStart + 1 )
              {
                  iXSrcStart = static_cast<int>(floor(dfXSrcStart+0.5));
                  iXSrcEnd = static_cast<int>(floor(dfXSrcEnd+0.5));
              }
              else
              {
                  /* If the resampling factor is less than 100%, the distance */
                  /* between the source pixel is < 1, so we stick to nearest */
                  /* neighbour */
                  iXSrcStart = static_cast<int>(floor(dfXSrcStart));
                  iXSrcEnd = iXSrcStart + 1;
              }
              if( dfYSrcEnd >= dfYSrcStart + 1 )
              {
                  iYSrcStart = static_cast<int>(floor(dfYSrcStart+0.5));
                  iYSrcEnd = static_cast<int>(floor(dfYSrcEnd+0.5));
              }
              else
              {
                  iYSrcStart = static_cast<int>(floor(dfYSrcStart));
                  iYSrcEnd = iYSrcStart + 1;
              }
  
              // Transform into the coordinate system of the source *buffer*
              iXSrcStart -= nReqXOff;
              iYSrcStart -= nReqYOff;
              iXSrcEnd -= nReqXOff;
              iYSrcEnd -= nReqYOff;
  
              double dfSum = 0.0;
              int nPixelCount = 0;
  
              for( int iY = iYSrcStart; iY < iYSrcEnd; iY++ )
              {
                  if( iY < 0 || iY >= nReqYSize )
                      continue;
  
                  for( int iX = iXSrcStart; iX < iXSrcEnd; iX++ )
                  {
                      if( iX < 0 || iX >= nReqXSize )
                          continue;
  
                      const float fSampledValue = pafSrc[iX + iY * nReqXSize];
                      if( CPLIsNan(fSampledValue) )
                          continue;
  
                      if( m_bNoDataSet &&
                          GDALIsValueInRange<float>(m_dfNoDataValue) &&
                          ARE_REAL_EQUAL(fSampledValue,
                                         static_cast<float>(m_dfNoDataValue)))
                          continue;
  
                      nPixelCount++;
                      dfSum += pafSrc[iX + iY * nReqXSize];
                  }
              }
  
              if( nPixelCount == 0 )
                  continue;
  
              // Compute output value.
              const float dfOutputValue = static_cast<float>(dfSum / nPixelCount);
  
              // Put it in the output buffer.
              GByte *pDstLocation =
                  static_cast<GByte *>(pData)
                  + nPixelSpace * iBufPixel
                  + static_cast<GPtrDiff_t>(nLineSpace) * iBufLine;
  
              if( eBufType == GDT_Byte )
                  *pDstLocation =
                      static_cast<GByte>(
                          std::min(255.0, std::max(0.0, dfOutputValue + 0.5)) );
              else
                  GDALCopyWords( &dfOutputValue, GDT_Float32, 4,
                                 pDstLocation, eBufType, 8, 1 );
          }
      }
  
      VSIFree( pafSrc );
  
      return CE_None;
  }
  
  /************************************************************************/
  /*                             GetMinimum()                             */
  /************************************************************************/
  
  double VRTAveragedSource::GetMinimum( int /* nXSize */,
                                        int /* nYSize */,
                                        int *pbSuccess )
  {
      *pbSuccess = FALSE;
      return 0.0;
  }
  
  /************************************************************************/
  /*                             GetMaximum()                             */
  /************************************************************************/
  
  double VRTAveragedSource::GetMaximum( int /* nXSize */,
                                        int /* nYSize */,
                                        int *pbSuccess )
  {
      *pbSuccess = FALSE;
      return 0.0;
  }
  
  /************************************************************************/
  /*                       ComputeRasterMinMax()                          */
  /************************************************************************/
  
  CPLErr VRTAveragedSource::ComputeRasterMinMax( int /* nXSize */,
                                                 int /* nYSize */,
                                                 int /* bApproxOK */,
                                                 double* /* adfMinMax */)
  {
      return CE_Failure;
  }
  
  /************************************************************************/
  /*                         ComputeStatistics()                          */
  /************************************************************************/
  
  CPLErr VRTAveragedSource::ComputeStatistics( int /* nXSize */,
                                               int /* nYSize */,
                                               int /* bApproxOK */,
                                               double * /* pdfMin */,
                                               double * /* pdfMax */,
                                               double * /* pdfMean */,
                                               double * /* pdfStdDev */,
                                               GDALProgressFunc /* pfnProgress */,
                                               void * /* pProgressData */ )
  {
      return CE_Failure;
  }
  
  /************************************************************************/
  /*                            GetHistogram()                            */
  /************************************************************************/
  
  CPLErr VRTAveragedSource::GetHistogram( int /* nXSize */,
                                          int /* nYSize */,
                                          double /* dfMin */,
                                          double /* dfMax */,
                                          int /* nBuckets */,
                                          GUIntBig * /* panHistogram */,
                                          int /* bIncludeOutOfRange */,
                                          int /* bApproxOK */,
                                          GDALProgressFunc /* pfnProgress */,
                                          void * /* pProgressData */)
  {
      return CE_Failure;
  }
  
  /************************************************************************/
  /* ==================================================================== */
  /*                          VRTComplexSource                            */
  /* ==================================================================== */
  /************************************************************************/
  
  /************************************************************************/
  /*                          VRTComplexSource()                          */
  /************************************************************************/
  
  VRTComplexSource::VRTComplexSource() :
      m_eScalingType(VRT_SCALING_NONE),
      m_dfScaleOff(0.0),
      m_dfScaleRatio(1.0),
      m_bSrcMinMaxDefined(FALSE),
      m_dfSrcMin(0.0),
      m_dfSrcMax(0.0),
      m_dfDstMin(0.0),
      m_dfDstMax(0.0),
      m_dfExponent(1.0),
      m_nColorTableComponent(0),
      m_padfLUTInputs(nullptr),
      m_padfLUTOutputs(nullptr),
      m_nLUTItemCount(0)
  {}
  
  VRTComplexSource::VRTComplexSource( const VRTComplexSource* poSrcSource,
                                      double dfXDstRatio, double dfYDstRatio ) :
      VRTSimpleSource(poSrcSource, dfXDstRatio, dfYDstRatio),
      m_eScalingType(poSrcSource->m_eScalingType),
      m_dfScaleOff(poSrcSource->m_dfScaleOff),
      m_dfScaleRatio(poSrcSource->m_dfScaleRatio),
      m_bSrcMinMaxDefined(poSrcSource->m_bSrcMinMaxDefined),
      m_dfSrcMin(poSrcSource->m_dfSrcMin),
      m_dfSrcMax(poSrcSource->m_dfSrcMax),
      m_dfDstMin(poSrcSource->m_dfDstMin),
      m_dfDstMax(poSrcSource->m_dfDstMax),
      m_dfExponent(poSrcSource->m_dfExponent),
      m_nColorTableComponent(poSrcSource->m_nColorTableComponent),
      m_padfLUTInputs(nullptr),
      m_padfLUTOutputs(nullptr),
      m_nLUTItemCount(poSrcSource->m_nLUTItemCount)
  {
      if( m_nLUTItemCount )
      {
          m_padfLUTInputs = static_cast<double *>(
              CPLMalloc(sizeof(double) * m_nLUTItemCount) );
          memcpy( m_padfLUTInputs, poSrcSource->m_padfLUTInputs,
                  sizeof(double) * m_nLUTItemCount );
  
          m_padfLUTOutputs = static_cast<double *>(
              CPLMalloc(sizeof(double) * m_nLUTItemCount) );
          memcpy( m_padfLUTOutputs, poSrcSource->m_padfLUTOutputs,
                  sizeof(double) * m_nLUTItemCount);
      }
  }
  
  VRTComplexSource::~VRTComplexSource()
  {
      VSIFree( m_padfLUTInputs );
      VSIFree( m_padfLUTOutputs );
  }
  
  /************************************************************************/
  /*                           SerializeToXML()                           */
  /************************************************************************/
  
  CPLXMLNode *VRTComplexSource::SerializeToXML( const char *pszVRTPath )
  
  {
      CPLXMLNode *psSrc = VRTSimpleSource::SerializeToXML( pszVRTPath );
  
      if( psSrc == nullptr )
          return nullptr;
  
      CPLFree( psSrc->pszValue );
      psSrc->pszValue = CPLStrdup( "ComplexSource" );
  
      if( m_bNoDataSet )
      {
          if( CPLIsNan(m_dfNoDataValue) )
              CPLSetXMLValue( psSrc, "NODATA", "nan");
          else
              CPLSetXMLValue( psSrc, "NODATA",
                              CPLSPrintf("%.16g", m_dfNoDataValue) );
      }
  
      switch( m_eScalingType )
      {
          case VRT_SCALING_NONE:
              break;
  
          case VRT_SCALING_LINEAR:
          {
              CPLSetXMLValue( psSrc, "ScaleOffset",
                              CPLSPrintf("%g", m_dfScaleOff) );
              CPLSetXMLValue( psSrc, "ScaleRatio",
                              CPLSPrintf("%g", m_dfScaleRatio) );
              break;
          }
  
          case VRT_SCALING_EXPONENTIAL:
          {
              CPLSetXMLValue( psSrc, "Exponent",
                              CPLSPrintf("%g", m_dfExponent) );
              CPLSetXMLValue( psSrc, "SrcMin",
                              CPLSPrintf("%g", m_dfSrcMin) );
              CPLSetXMLValue( psSrc, "SrcMax",
                              CPLSPrintf("%g", m_dfSrcMax) );
              CPLSetXMLValue( psSrc, "DstMin",
                              CPLSPrintf("%g", m_dfDstMin) );
              CPLSetXMLValue( psSrc, "DstMax",
                              CPLSPrintf("%g", m_dfDstMax) );
              break;
          }
      }
  
      if( m_nLUTItemCount )
      {
          // Make sure we print with sufficient precision to address really close
          // entries (#6422).
          CPLString osLUT;
          // TODO(schwehr): How is this not a read past the end of the array if
          // m_nLUTItemCount is 0 or 1?  Added in
          // https://trac.osgeo.org/gdal/changeset/33779
          if( m_nLUTItemCount > 0 &&
              CPLString().Printf("%g", m_padfLUTInputs[0]) ==
              CPLString().Printf("%g", m_padfLUTInputs[1]) )
          {
              osLUT = CPLString().Printf(
                  "%.18g:%g", m_padfLUTInputs[0], m_padfLUTOutputs[0]);
          }
          else
          {
              osLUT = CPLString().Printf(
                  "%g:%g", m_padfLUTInputs[0], m_padfLUTOutputs[0]);
          }
          for ( int i = 1; i < m_nLUTItemCount; i++ )
          {
              if( CPLString().Printf("%g", m_padfLUTInputs[i]) ==
                  CPLString().Printf("%g", m_padfLUTInputs[i-1]) ||
                  (i + 1 < m_nLUTItemCount &&
                   CPLString().Printf("%g", m_padfLUTInputs[i]) ==
                   CPLString().Printf("%g", m_padfLUTInputs[i+1])) )
              {
                  // TODO(schwehr): An explanation of the 18 would be helpful.
                  // Can someone distill the issue down to a quick comment?
                  // https://trac.osgeo.org/gdal/ticket/6422
                  osLUT += CPLString().Printf(
                      ",%.18g:%g", m_padfLUTInputs[i], m_padfLUTOutputs[i]);
              }
              else
              {
                  osLUT += CPLString().Printf(
                      ",%g:%g", m_padfLUTInputs[i], m_padfLUTOutputs[i]);
              }
          }
          CPLSetXMLValue( psSrc, "LUT", osLUT );
      }
  
      if( m_nColorTableComponent )
      {
          CPLSetXMLValue( psSrc, "ColorTableComponent",
                          CPLSPrintf("%d", m_nColorTableComponent) );
      }
  
      return psSrc;
  }
  
  /************************************************************************/
  /*                              XMLInit()                               */
  /************************************************************************/
  
  CPLErr VRTComplexSource::XMLInit( CPLXMLNode *psSrc, const char *pszVRTPath,
                                    void* pUniqueHandle )
  
  {
  /* -------------------------------------------------------------------- */
  /*      Do base initialization.                                         */
  /* -------------------------------------------------------------------- */
      {
          const CPLErr eErr = VRTSimpleSource::XMLInit( psSrc, pszVRTPath,
                                                        pUniqueHandle );
          if( eErr != CE_None )
              return eErr;
      }
  
  /* -------------------------------------------------------------------- */
  /*      Complex parameters.                                             */
  /* -------------------------------------------------------------------- */
      if( CPLGetXMLValue(psSrc, "ScaleOffset", nullptr) != nullptr
          || CPLGetXMLValue(psSrc, "ScaleRatio", nullptr) != nullptr )
      {
          m_eScalingType = VRT_SCALING_LINEAR;
          m_dfScaleOff = CPLAtof( CPLGetXMLValue(psSrc, "ScaleOffset", "0") );
          m_dfScaleRatio = CPLAtof( CPLGetXMLValue(psSrc, "ScaleRatio", "1") );
      }
      else if( CPLGetXMLValue(psSrc, "Exponent", nullptr) != nullptr &&
               CPLGetXMLValue(psSrc, "DstMin", nullptr) != nullptr &&
               CPLGetXMLValue(psSrc, "DstMax", nullptr) != nullptr )
      {
          m_eScalingType = VRT_SCALING_EXPONENTIAL;
          m_dfExponent = CPLAtof( CPLGetXMLValue(psSrc, "Exponent", "1.0") );
  
          if( CPLGetXMLValue(psSrc, "SrcMin", nullptr) != nullptr
           && CPLGetXMLValue(psSrc, "SrcMax", nullptr) != nullptr )
          {
              m_dfSrcMin = CPLAtof( CPLGetXMLValue(psSrc, "SrcMin", "0.0") );
              m_dfSrcMax = CPLAtof( CPLGetXMLValue(psSrc, "SrcMax", "0.0") );
              m_bSrcMinMaxDefined = TRUE;
          }
  
          m_dfDstMin = CPLAtof( CPLGetXMLValue(psSrc, "DstMin", "0.0") );
          m_dfDstMax = CPLAtof( CPLGetXMLValue(psSrc, "DstMax", "0.0") );
      }
  
      if( CPLGetXMLValue(psSrc, "NODATA", nullptr) != nullptr )
      {
          m_bNoDataSet = TRUE;
          m_dfNoDataValue = CPLAtofM( CPLGetXMLValue(psSrc, "NODATA", "0") );
      }
  
      if( CPLGetXMLValue(psSrc, "LUT", nullptr) != nullptr )
      {
          char **papszValues =
              CSLTokenizeString2(
                  CPLGetXMLValue(psSrc, "LUT", ""), ",:", CSLT_ALLOWEMPTYTOKENS);
  
          if( m_nLUTItemCount )
          {
              if( m_padfLUTInputs )
              {
                  VSIFree( m_padfLUTInputs );
                  m_padfLUTInputs = nullptr;
              }
              if( m_padfLUTOutputs )
              {
                  VSIFree( m_padfLUTOutputs );
                  m_padfLUTOutputs = nullptr;
              }
              m_nLUTItemCount = 0;
          }
  
          m_nLUTItemCount = CSLCount(papszValues) / 2;
  
          m_padfLUTInputs = static_cast<double *>(
              VSIMalloc2(m_nLUTItemCount, sizeof(double)) );
          if( !m_padfLUTInputs )
          {
              CSLDestroy(papszValues);
              m_nLUTItemCount = 0;
              return CE_Failure;
          }
  
          m_padfLUTOutputs = static_cast<double *>(
              VSIMalloc2(m_nLUTItemCount, sizeof(double)) );
          if( !m_padfLUTOutputs )
          {
              CSLDestroy(papszValues);
              VSIFree( m_padfLUTInputs );
              m_padfLUTInputs = nullptr;
              m_nLUTItemCount = 0;
              return CE_Failure;
          }
  
          for ( int nIndex = 0; nIndex < m_nLUTItemCount; nIndex++ )
          {
              m_padfLUTInputs[nIndex] = CPLAtof( papszValues[nIndex * 2] );
              m_padfLUTOutputs[nIndex] = CPLAtof( papszValues[nIndex * 2 + 1] );
  
              // Enforce the requirement that the LUT input array is
              // monotonically non-decreasing.
              if( nIndex > 0 &&
                  m_padfLUTInputs[nIndex] < m_padfLUTInputs[nIndex - 1] )
              {
                  CSLDestroy(papszValues);
                  VSIFree( m_padfLUTInputs );
                  VSIFree( m_padfLUTOutputs );
                  m_padfLUTInputs = nullptr;
                  m_padfLUTOutputs = nullptr;
                  m_nLUTItemCount = 0;
                  return CE_Failure;
              }
          }
  
          CSLDestroy(papszValues);
      }
  
      if( CPLGetXMLValue(psSrc, "ColorTableComponent", nullptr) != nullptr )
      {
          m_nColorTableComponent =
              atoi( CPLGetXMLValue(psSrc, "ColorTableComponent", "0") );
      }
  
      return CE_None;
  }
  
  /************************************************************************/
  /*                              LookupValue()                           */
  /************************************************************************/
  
  double
  VRTComplexSource::LookupValue( double dfInput )
  {
      // Find the index of the first element in the LUT input array that
      // is not smaller than the input value.
      int i = static_cast<int>(
          std::lower_bound( m_padfLUTInputs,
                            m_padfLUTInputs + m_nLUTItemCount,
                            dfInput)
          - m_padfLUTInputs );
  
      if( i == 0 )
          return m_padfLUTOutputs[0];
  
      // If the index is beyond the end of the LUT input array, the input
      // value is larger than all the values in the array.
      if( i == m_nLUTItemCount )
          return m_padfLUTOutputs[m_nLUTItemCount - 1];
  
      if( m_padfLUTInputs[i] == dfInput )
          return m_padfLUTOutputs[i];
  
      // Otherwise, interpolate.
      return
          m_padfLUTOutputs[i - 1] + (dfInput - m_padfLUTInputs[i - 1]) *
          ( (m_padfLUTOutputs[i] - m_padfLUTOutputs[i - 1]) /
            (m_padfLUTInputs[i] - m_padfLUTInputs[i - 1]) );
  }
  
  /************************************************************************/
  /*                         SetLinearScaling()                           */
  /************************************************************************/
  
  void VRTComplexSource::SetLinearScaling( double dfOffset, double dfScale )
  {
      m_eScalingType = VRT_SCALING_LINEAR;
      m_dfScaleOff = dfOffset;
      m_dfScaleRatio = dfScale;
  }
  
  /************************************************************************/
  /*                         SetPowerScaling()                           */
  /************************************************************************/
  
  void VRTComplexSource::SetPowerScaling( double dfExponentIn,
                                          double dfSrcMinIn,
                                          double dfSrcMaxIn,
                                          double dfDstMinIn,
                                          double dfDstMaxIn)
  {
      m_eScalingType = VRT_SCALING_EXPONENTIAL;
      m_dfExponent = dfExponentIn;
      m_dfSrcMin = dfSrcMinIn;
      m_dfSrcMax = dfSrcMaxIn;
      m_dfDstMin = dfDstMinIn;
      m_dfDstMax = dfDstMaxIn;
      m_bSrcMinMaxDefined = TRUE;
  }
  
  /************************************************************************/
  /*                    SetColorTableComponent()                          */
  /************************************************************************/
  
  void VRTComplexSource::SetColorTableComponent( int nComponent )
  {
      m_nColorTableComponent = nComponent;
  }
  
  /************************************************************************/
  /*                              RasterIO()                              */
  /************************************************************************/
  
  CPLErr
  VRTComplexSource::RasterIO( GDALDataType /*eBandDataType*/,
                              int nXOff, int nYOff, int nXSize, int nYSize,
                              void *pData, int nBufXSize, int nBufYSize,
                              GDALDataType eBufType,
                              GSpacing nPixelSpace,
                              GSpacing nLineSpace,
                              GDALRasterIOExtraArg* psExtraArgIn )
  
  {
      GDALRasterIOExtraArg sExtraArg;
      INIT_RASTERIO_EXTRA_ARG(sExtraArg);
      GDALRasterIOExtraArg* psExtraArg = &sExtraArg;
  
      // The window we will actually request from the source raster band.
      double dfReqXOff = 0.0;
      double dfReqYOff = 0.0;
      double dfReqXSize = 0.0;
      double dfReqYSize = 0.0;
      int nReqXOff = 0;
      int nReqYOff = 0;
      int nReqXSize = 0;
      int nReqYSize = 0;
  
      // The window we will actual set _within_ the pData buffer.
      int nOutXOff = 0;
      int nOutYOff = 0;
      int nOutXSize = 0;
      int nOutYSize = 0;
  
      if( !GetSrcDstWindow( nXOff, nYOff, nXSize, nYSize, nBufXSize, nBufYSize,
                            &dfReqXOff, &dfReqYOff, &dfReqXSize, &dfReqYSize,
                            &nReqXOff, &nReqYOff, &nReqXSize, &nReqYSize,
                            &nOutXOff, &nOutYOff, &nOutXSize, &nOutYSize ) )
          return CE_None;
  #if DEBUG_VERBOSE
      CPLDebug(
          "VRT",
          "nXOff=%d, nYOff=%d, nXSize=%d, nYSize=%d, nBufXSize=%d, nBufYSize=%d,\n"
          "dfReqXOff=%g, dfReqYOff=%g, dfReqXSize=%g, dfReqYSize=%g,\n"
          "nReqXOff=%d, nReqYOff=%d, nReqXSize=%d, nReqYSize=%d,\n"
          "nOutXOff=%d, nOutYOff=%d, nOutXSize=%d, nOutYSize=%d",
          nXOff, nYOff, nXSize, nYSize,
          nBufXSize, nBufYSize,
          dfReqXOff, dfReqYOff, dfReqXSize, dfReqYSize,
          nReqXOff, nReqYOff, nReqXSize, nReqYSize,
          nOutXOff, nOutYOff, nOutXSize, nOutYSize );
  #endif
  
      if( !m_osResampling.empty() )
      {
          psExtraArg->eResampleAlg = GDALRasterIOGetResampleAlg(m_osResampling);
      }
      else if( psExtraArgIn != nullptr )
      {
          psExtraArg->eResampleAlg = psExtraArgIn->eResampleAlg;
      }
      psExtraArg->bFloatingPointWindowValidity = TRUE;
      psExtraArg->dfXOff = dfReqXOff;
      psExtraArg->dfYOff = dfReqYOff;
      psExtraArg->dfXSize = dfReqXSize;
      psExtraArg->dfYSize = dfReqYSize;
  
      const bool bIsComplex = CPL_TO_BOOL( GDALDataTypeIsComplex(eBufType) );
      CPLErr eErr;
      // For Int32, float32 isn't sufficiently precise as working data type
      if( eBufType == GDT_CInt32 || eBufType == GDT_CFloat64 ||
          eBufType == GDT_Int32 || eBufType == GDT_UInt32 ||
          eBufType == GDT_Float64 )
      {
          eErr = RasterIOInternal<double>(
                  nReqXOff, nReqYOff, nReqXSize, nReqYSize,
                  static_cast<GByte *>(pData) + nPixelSpace * nOutXOff
                      + static_cast<GPtrDiff_t>(nLineSpace) * nOutYOff,
                  nOutXSize, nOutYSize,
                  eBufType,
                  nPixelSpace, nLineSpace, psExtraArg,
                  bIsComplex ? GDT_CFloat64 : GDT_Float64);
      }
      else
      {
          eErr = RasterIOInternal<float>(
                  nReqXOff, nReqYOff, nReqXSize, nReqYSize,
                  static_cast<GByte *>(pData) + nPixelSpace * nOutXOff
                      + static_cast<GPtrDiff_t>(nLineSpace) * nOutYOff,
                  nOutXSize, nOutYSize,
                  eBufType,
                  nPixelSpace, nLineSpace, psExtraArg,
                  bIsComplex ? GDT_CFloat32 : GDT_Float32);
      }
  
      return eErr;
  }
  
  /************************************************************************/
  /*                          RasterIOInternal()                          */
  /************************************************************************/
  
  // nReqXOff, nReqYOff, nReqXSize, nReqYSize are expressed in source band
  // referential.
  template <class WorkingDT>
  CPLErr VRTComplexSource::RasterIOInternal( int nReqXOff, int nReqYOff,
                                             int nReqXSize, int nReqYSize,
                                             void *pData,
                                             int nOutXSize, int nOutYSize,
                                             GDALDataType eBufType,
                                             GSpacing nPixelSpace,
                                             GSpacing nLineSpace,
                                             GDALRasterIOExtraArg* psExtraArg,
                                             GDALDataType eWrkDataType )
  {
  /* -------------------------------------------------------------------- */
  /*      Read into a temporary buffer.                                   */
  /* -------------------------------------------------------------------- */
      GDALColorTable* poColorTable = nullptr;
      const bool bIsComplex = CPL_TO_BOOL( GDALDataTypeIsComplex(eBufType) );
      const int nWordSize = GDALGetDataTypeSizeBytes(eWrkDataType);
      const bool bNoDataSetIsNan = m_bNoDataSet && CPLIsNan(m_dfNoDataValue);
      const bool bNoDataSetAndNotNan = m_bNoDataSet && !CPLIsNan(m_dfNoDataValue);
  
      WorkingDT *pafData = nullptr;
      if( m_eScalingType == VRT_SCALING_LINEAR &&
          m_bNoDataSet == FALSE &&
          m_dfScaleRatio == 0 )
      {
  /* -------------------------------------------------------------------- */
  /*      Optimization when writing a constant value                      */
  /*      (used by the -addalpha option of gdalbuildvrt)                  */
  /* -------------------------------------------------------------------- */
          // Already set to NULL when defined.
          // pafData = NULL;
      }
      else
      {
          pafData = static_cast<WorkingDT *>(
              VSI_MALLOC3_VERBOSE(nOutXSize,nOutYSize,nWordSize) );
          if( pafData == nullptr )
          {
              return CE_Failure;
          }
  
          const GDALRIOResampleAlg eResampleAlgBack = psExtraArg->eResampleAlg;
          if( !m_osResampling.empty() )
          {
              psExtraArg->eResampleAlg =
                  GDALRasterIOGetResampleAlg(m_osResampling);
          }
  
          const CPLErr eErr =
              m_poRasterBand->RasterIO( GF_Read,
                                        nReqXOff, nReqYOff,
                                        nReqXSize, nReqYSize,
                                        pafData,
                                        nOutXSize, nOutYSize,
                                        eWrkDataType,
                                        nWordSize,
                                        nWordSize *
                                        static_cast<GSpacing>(nOutXSize),
                                        psExtraArg );
          if( !m_osResampling.empty() )
              psExtraArg->eResampleAlg = eResampleAlgBack;
  
          if( eErr != CE_None )
          {
              CPLFree( pafData );
              return eErr;
          }
  
          if( m_nColorTableComponent != 0 )
          {
              poColorTable = m_poRasterBand->GetColorTable();
              if( poColorTable == nullptr )
              {
                  CPLError( CE_Failure, CPLE_AppDefined,
                            "Source band has no color table." );
                  CPLFree( pafData );
                  return CE_Failure;
              }
          }
      }
  
  /* -------------------------------------------------------------------- */
  /*      Selectively copy into output buffer with nodata masking,        */
  /*      and/or scaling.                                                 */
  /* -------------------------------------------------------------------- */
      for( int iY = 0; iY < nOutYSize; iY++ )
      {
          for( int iX = 0; iX < nOutXSize; iX++ )
          {
              GByte *pDstLocation =
                  static_cast<GByte *>(pData)
                  + nPixelSpace * iX
                  + static_cast<GPtrDiff_t>(nLineSpace) * iY;
  
              if( pafData && !bIsComplex )
              {
                  WorkingDT fResult = pafData[iX + iY * nOutXSize];
                  if( bNoDataSetIsNan && CPLIsNan(fResult) )
                      continue;
                  if( bNoDataSetAndNotNan &&
                      GDALIsValueInRange<WorkingDT>(m_dfNoDataValue) &&
                      ARE_REAL_EQUAL(fResult, static_cast<WorkingDT>(m_dfNoDataValue)) )
                      continue;
  
                  if( m_nColorTableComponent )
                  {
                      const GDALColorEntry* poEntry =
                          poColorTable->GetColorEntry(static_cast<int>(fResult));
                      if( poEntry )
                      {
                          if( m_nColorTableComponent == 1 )
                              fResult = poEntry->c1;
                          else if( m_nColorTableComponent == 2 )
                              fResult = poEntry->c2;
                          else if( m_nColorTableComponent == 3 )
                              fResult = poEntry->c3;
                          else if( m_nColorTableComponent == 4 )
                              fResult = poEntry->c4;
                      }
                      else
                      {
                          static bool bHasWarned = false;
                          if( !bHasWarned )
                          {
                              bHasWarned = true;
                              CPLError(
                                  CE_Failure, CPLE_AppDefined,
                                  "No entry %d.", static_cast<int>(fResult) );
                          }
                          continue;
                      }
                  }
  
                  if( m_eScalingType == VRT_SCALING_LINEAR )
                  {
                      fResult = static_cast<WorkingDT>(
                          fResult * m_dfScaleRatio + m_dfScaleOff );
                  }
                  else if( m_eScalingType == VRT_SCALING_EXPONENTIAL )
                  {
                      if( !m_bSrcMinMaxDefined )
                      {
                          int bSuccessMin = FALSE;
                          int bSuccessMax = FALSE;
                          double adfMinMax[2] = {
                              m_poRasterBand->GetMinimum(&bSuccessMin),
                              m_poRasterBand->GetMaximum(&bSuccessMax) };
                          if( (bSuccessMin && bSuccessMax) ||
                              m_poRasterBand->ComputeRasterMinMax( TRUE,
                                                                   adfMinMax )
                              == CE_None )
                          {
                              m_dfSrcMin = adfMinMax[0];
                              m_dfSrcMax = adfMinMax[1];
                              m_bSrcMinMaxDefined = TRUE;
                          }
                          else
                          {
                              CPLError(
                                  CE_Failure, CPLE_AppDefined,
                                  "Cannot determine source min/max value" );
                              return CE_Failure;
                          }
                      }
  
                      double dfPowVal =
                          (fResult - m_dfSrcMin) / (m_dfSrcMax - m_dfSrcMin);
                      if( dfPowVal < 0.0 )
                          dfPowVal = 0.0;
                      else if( dfPowVal > 1.0 )
                          dfPowVal = 1.0;
                      fResult = static_cast<WorkingDT>(
                          (m_dfDstMax - m_dfDstMin) *
                          pow( dfPowVal, m_dfExponent ) +
                          m_dfDstMin);
                  }
  
                  if( m_nLUTItemCount )
                      fResult = static_cast<WorkingDT>(LookupValue( fResult ));
  
                  if( m_nMaxValue != 0 && fResult > m_nMaxValue )
                      fResult = static_cast<WorkingDT>(m_nMaxValue);
  
                  if( eBufType == GDT_Byte )
                      *pDstLocation = static_cast<GByte>(
                          std::min(255.0f,
                                   std::max(0.0f,
                                            static_cast<float>(fResult) + 0.5f)));
                  else
                      GDALCopyWords( &fResult, eWrkDataType, 0,
                                  pDstLocation, eBufType, 0, 1 );
              }
              else if( pafData && bIsComplex )
              {
                  WorkingDT afResult[2] = {
                      pafData[2 * (iX + iY * nOutXSize)],
                      pafData[2 * (iX + iY * nOutXSize) + 1] };
  
                  // Do not use color table.
                  if( m_eScalingType == VRT_SCALING_LINEAR )
                  {
                      afResult[0] = static_cast<WorkingDT>(
                          afResult[0] * m_dfScaleRatio + m_dfScaleOff );
                      afResult[1] = static_cast<WorkingDT>(
                          afResult[1] * m_dfScaleRatio + m_dfScaleOff );
                  }
  
                  /* Do not use LUT */
  
                  if( eBufType == GDT_Byte )
                      *pDstLocation = static_cast<GByte>(
                          std::min(255.0, std::max(0.0, afResult[0] + 0.5) ) );
                  else
                      GDALCopyWords( afResult, eWrkDataType, 0,
                                     pDstLocation, eBufType, 0, 1 );
              }
              else
              {
                  WorkingDT afResult[2] = {
                      static_cast<WorkingDT>(m_dfScaleOff),
                      0 };
  
                  if( m_nLUTItemCount )
                      afResult[0] = static_cast<WorkingDT>(LookupValue( afResult[0] ));
  
                  if( m_nMaxValue != 0 && afResult[0] > m_nMaxValue )
                      afResult[0] = static_cast<WorkingDT>(m_nMaxValue);
  
                  if( eBufType == GDT_Byte )
                      *pDstLocation = static_cast<GByte>(
                          std::min(255.0, std::max(0.0, afResult[0] + 0.5)) );
                  else
                      GDALCopyWords( afResult, eWrkDataType, 0,
                                     pDstLocation, eBufType, 0, 1 );
              }
          }
      }
  
      CPLFree( pafData );
  
      return CE_None;
  }
  
  // Explicitly instantiate template method, as it is used in another file.
  template
  CPLErr VRTComplexSource::RasterIOInternal<float>( int nReqXOff, int nReqYOff,
                                      int nReqXSize, int nReqYSize,
                                      void *pData, int nOutXSize, int nOutYSize,
                                      GDALDataType eBufType,
                                      GSpacing nPixelSpace, GSpacing nLineSpace,
                                      GDALRasterIOExtraArg* psExtraArg,
                                      GDALDataType eWrkDataType );
  
  /************************************************************************/
  /*                             GetMinimum()                             */
  /************************************************************************/
  
  double VRTComplexSource::GetMinimum( int nXSize, int nYSize, int *pbSuccess )
  {
      if( m_dfScaleOff == 0.0 && m_dfScaleRatio == 1.0 &&
          m_nLUTItemCount == 0 && m_nColorTableComponent == 0 )
      {
          return VRTSimpleSource::GetMinimum(nXSize, nYSize, pbSuccess);
      }
  
      *pbSuccess = FALSE;
      return 0;
  }
  
  /************************************************************************/
  /*                             GetMaximum()                             */
  /************************************************************************/
  
  double VRTComplexSource::GetMaximum( int nXSize, int nYSize, int *pbSuccess )
  {
      if( m_dfScaleOff == 0.0 && m_dfScaleRatio == 1.0 &&
          m_nLUTItemCount == 0 && m_nColorTableComponent == 0 )
      {
          return VRTSimpleSource::GetMaximum(nXSize, nYSize, pbSuccess);
      }
  
      *pbSuccess = FALSE;
      return 0;
  }
  
  /************************************************************************/
  /*                       ComputeRasterMinMax()                          */
  /************************************************************************/
  
  CPLErr VRTComplexSource::ComputeRasterMinMax( int nXSize, int nYSize, int bApproxOK, double* adfMinMax )
  {
      if( m_dfScaleOff == 0.0 && m_dfScaleRatio == 1.0 &&
          m_nLUTItemCount == 0 && m_nColorTableComponent == 0 )
      {
          return VRTSimpleSource::ComputeRasterMinMax( nXSize, nYSize, bApproxOK,
                                                       adfMinMax);
      }
  
      return CE_Failure;
  }
  
  /************************************************************************/
  /*                            GetHistogram()                            */
  /************************************************************************/
  
  CPLErr VRTComplexSource::GetHistogram( int nXSize, int nYSize,
                                         double dfMin, double dfMax,
                                         int nBuckets, GUIntBig * panHistogram,
                                         int bIncludeOutOfRange, int bApproxOK,
                                         GDALProgressFunc pfnProgress,
                                         void *pProgressData )
  {
      if( m_dfScaleOff == 0.0 && m_dfScaleRatio == 1.0 &&
          m_nLUTItemCount == 0 && m_nColorTableComponent == 0 )
      {
          return VRTSimpleSource::GetHistogram( nXSize, nYSize,
                                                dfMin, dfMax, nBuckets,
                                                panHistogram,
                                                bIncludeOutOfRange, bApproxOK,
                                                pfnProgress, pProgressData );
      }
  
      return CE_Failure;
  }
  
  /************************************************************************/
  /*                         ComputeStatistics()                          */
  /************************************************************************/
  
  CPLErr VRTComplexSource::ComputeStatistics( int nXSize, int nYSize,
                                              int bApproxOK,
                                              double *pdfMin, double *pdfMax,
                                              double *pdfMean, double *pdfStdDev,
                                              GDALProgressFunc pfnProgress,
                                              void *pProgressData )
  {
      if( m_dfScaleOff == 0.0 && m_dfScaleRatio == 1.0 &&
          m_nLUTItemCount == 0 && m_nColorTableComponent == 0 )
      {
          return VRTSimpleSource::ComputeStatistics( nXSize, nYSize, bApproxOK,
                                                     pdfMin, pdfMax,
                                                     pdfMean, pdfStdDev,
                                                     pfnProgress, pProgressData );
      }
  
      return CE_Failure;
  }
  
  /************************************************************************/
  /* ==================================================================== */
  /*                          VRTFuncSource                               */
  /* ==================================================================== */
  /************************************************************************/
  
  /************************************************************************/
  /*                           VRTFuncSource()                            */
  /************************************************************************/
  
  VRTFuncSource::VRTFuncSource() :
      pfnReadFunc(nullptr),
      pCBData(nullptr),
      eType(GDT_Byte),
      fNoDataValue(static_cast<float>(VRT_NODATA_UNSET))
  {}
  
  /************************************************************************/
  /*                           ~VRTFuncSource()                           */
  /************************************************************************/
  
  VRTFuncSource::~VRTFuncSource() {}
  
  /************************************************************************/
  /*                           SerializeToXML()                           */
  /************************************************************************/
  
  CPLXMLNode *VRTFuncSource::SerializeToXML( CPL_UNUSED const char * pszVRTPath )
  {
      return nullptr;
  }
  
  /************************************************************************/
  /*                              RasterIO()                              */
  /************************************************************************/
  
  CPLErr
  VRTFuncSource::RasterIO( GDALDataType /*eBandDataType*/,
                           int nXOff, int nYOff, int nXSize, int nYSize,
                           void *pData, int nBufXSize, int nBufYSize,
                           GDALDataType eBufType,
                           GSpacing nPixelSpace,
                           GSpacing nLineSpace,
                           GDALRasterIOExtraArg* /* psExtraArg */ )
  {
      if( nPixelSpace*8 == GDALGetDataTypeSize( eBufType )
          && nLineSpace == nPixelSpace * nXSize
          && nBufXSize == nXSize && nBufYSize == nYSize
          && eBufType == eType )
      {
          return pfnReadFunc( pCBData,
                              nXOff, nYOff, nXSize, nYSize,
                              pData );
      }
      else
      {
          CPLError( CE_Failure, CPLE_AppDefined,
                    "VRTFuncSource::RasterIO() - Irregular request." );
          CPLDebug("VRT", "Irregular request: %d,%d  %d,%d, %d,%d %d,%d %d,%d",
                  static_cast<int>(nPixelSpace)*8,
                  GDALGetDataTypeSize(eBufType),
                  static_cast<int>(nLineSpace),
                  static_cast<int>(nPixelSpace) * nXSize,
                  nBufXSize, nXSize,
                  nBufYSize, nYSize,
                  static_cast<int>(eBufType),
                  static_cast<int>(eType) );
  
          return CE_Failure;
      }
  }
  
  /************************************************************************/
  /*                             GetMinimum()                             */
  /************************************************************************/
  
  double VRTFuncSource::GetMinimum( int /* nXSize */,
                                    int /* nYSize */,
                                    int *pbSuccess )
  {
      *pbSuccess = FALSE;
      return 0;
  }
  
  /************************************************************************/
  /*                             GetMaximum()                             */
  /************************************************************************/
  
  double VRTFuncSource::GetMaximum( int /* nXSize */,
                                    int /* nYSize */,
                                    int *pbSuccess )
  {
      *pbSuccess = FALSE;
      return 0;
  }
  
  /************************************************************************/
  /*                       ComputeRasterMinMax()                          */
  /************************************************************************/
  
  CPLErr VRTFuncSource::ComputeRasterMinMax( int /* nXSize */,
                                             int /* nYSize */,
                                             int /* bApproxOK */,
                                             double* /* adfMinMax */ )
  {
      return CE_Failure;
  }
  
  /************************************************************************/
  /*                         ComputeStatistics()                          */
  /************************************************************************/
  
  CPLErr VRTFuncSource::ComputeStatistics( int /* nXSize */,
                                           int /* nYSize */,
                                           int /* bApproxOK */,
                                           double * /* pdfMin */,
                                           double * /* pdfMax */,
                                           double * /* pdfMean */,
                                           double * /* pdfStdDev */,
                                           GDALProgressFunc /* pfnProgress */,
                                           void * /* pProgressData */ )
  {
      return CE_Failure;
  }
  
  /************************************************************************/
  /*                            GetHistogram()                            */
  /************************************************************************/
  
  CPLErr VRTFuncSource::GetHistogram( int /* nXSize */,
                                      int /* nYSize */,
                                      double /* dfMin */,
                                      double /* dfMax */,
                                      int /* nBuckets */,
                                      GUIntBig * /* panHistogram */,
                                      int /* bIncludeOutOfRange */,
                                      int /* bApproxOK */,
                                      GDALProgressFunc /* pfnProgress */,
                                      void * /* pProgressData */)
  {
      return CE_Failure;
  }
  
  /************************************************************************/
  /*                        VRTParseCoreSources()                         */
  /************************************************************************/
  
  VRTSource *VRTParseCoreSources( CPLXMLNode *psChild, const char *pszVRTPath,
                                  void* pUniqueHandle )
  
  {
      VRTSource * poSource = nullptr;
  
      if( EQUAL(psChild->pszValue,"AveragedSource")
          || (EQUAL(psChild->pszValue,"SimpleSource")
              && STARTS_WITH_CI(CPLGetXMLValue(psChild, "Resampling", "Nearest"),
                                "Aver")) )
      {
          poSource = new VRTAveragedSource();
      }
      else if( EQUAL(psChild->pszValue,"SimpleSource") )
      {
          poSource = new VRTSimpleSource();
      }
      else if( EQUAL(psChild->pszValue,"ComplexSource") )
      {
          poSource = new VRTComplexSource();
      }
      else
      {
          CPLError( CE_Failure, CPLE_AppDefined,
                    "VRTParseCoreSources() - Unknown source : %s", psChild->pszValue );
          return nullptr;
      }
  
      if( poSource->XMLInit( psChild, pszVRTPath, pUniqueHandle ) == CE_None )
          return poSource;
  
      delete poSource;
      return nullptr;
  }
  /*! @endcond */
  