LCOV - code coverage report
Current view: top level - frmts/wmts - wmtsdataset.cpp (source / functions) Hit Total Coverage
Test: gdal_filtered.info Lines: 1141 1259 90.6 %
Date: 2025-09-10 17:48:50 Functions: 30 31 96.8 %

          Line data    Source code
       1             : /******************************************************************************
       2             :  *
       3             :  * Project:  GDAL WMTS driver
       4             :  * Purpose:  Implement GDAL WMTS support
       5             :  * Author:   Even Rouault, <even dot rouault at spatialys dot com>
       6             :  * Funded by Land Information New Zealand (LINZ)
       7             :  *
       8             :  **********************************************************************
       9             :  * Copyright (c) 2015, Even Rouault <even dot rouault at spatialys dot com>
      10             :  *
      11             :  * SPDX-License-Identifier: MIT
      12             :  ****************************************************************************/
      13             : 
      14             : #include "cpl_http.h"
      15             : #include "cpl_minixml.h"
      16             : #include "gdal_frmts.h"
      17             : #include "gdal_pam.h"
      18             : #include "ogr_spatialref.h"
      19             : #include "../vrt/gdal_vrt.h"
      20             : #include "wmtsdrivercore.h"
      21             : 
      22             : #include <algorithm>
      23             : #include <array>
      24             : #include <cmath>
      25             : #include <map>
      26             : #include <set>
      27             : #include <vector>
      28             : #include <limits>
      29             : 
      30             : extern "C" void GDALRegister_WMTS();
      31             : 
      32             : // g++ -g -Wall -fPIC frmts/wmts/wmtsdataset.cpp -shared -o gdal_WMTS.so -Iport
      33             : // -Igcore -Iogr -Iogr/ogrsf_frmts -L. -lgdal
      34             : 
      35             : /* Set in stone by WMTS spec. In pixel/meter */
      36             : #define WMTS_PITCH 0.00028
      37             : 
      38             : #define WMTS_WGS84_DEG_PER_METER (180 / M_PI / SRS_WGS84_SEMIMAJOR)
      39             : 
      40             : typedef enum
      41             : {
      42             :     AUTO,
      43             :     LAYER_BBOX,
      44             :     TILE_MATRIX_SET,
      45             :     MOST_PRECISE_TILE_MATRIX
      46             : } ExtentMethod;
      47             : 
      48             : /************************************************************************/
      49             : /* ==================================================================== */
      50             : /*                            WMTSTileMatrix                            */
      51             : /* ==================================================================== */
      52             : /************************************************************************/
      53             : 
      54             : class WMTSTileMatrix
      55             : {
      56             :   public:
      57             :     CPLString osIdentifier{};
      58             :     double dfScaleDenominator = 0;
      59             :     double dfPixelSize = 0;
      60             :     double dfTLX = 0;
      61             :     double dfTLY = 0;
      62             :     int nTileWidth = 0;
      63             :     int nTileHeight = 0;
      64             :     int nMatrixWidth = 0;
      65             :     int nMatrixHeight = 0;
      66             : 
      67         263 :     OGREnvelope GetExtent() const
      68             :     {
      69         263 :         OGREnvelope sExtent;
      70         263 :         sExtent.MinX = dfTLX;
      71         263 :         sExtent.MaxX = dfTLX + nMatrixWidth * dfPixelSize * nTileWidth;
      72         263 :         sExtent.MaxY = dfTLY;
      73         263 :         sExtent.MinY = dfTLY - nMatrixHeight * dfPixelSize * nTileHeight;
      74         263 :         return sExtent;
      75             :     }
      76             : };
      77             : 
      78             : /************************************************************************/
      79             : /* ==================================================================== */
      80             : /*                          WMTSTileMatrixLimits                        */
      81             : /* ==================================================================== */
      82             : /************************************************************************/
      83             : 
      84             : class WMTSTileMatrixLimits
      85             : {
      86             :   public:
      87             :     CPLString osIdentifier{};
      88             :     int nMinTileRow = 0;
      89             :     int nMaxTileRow = 0;
      90             :     int nMinTileCol = 0;
      91             :     int nMaxTileCol = 0;
      92             : 
      93           4 :     OGREnvelope GetExtent(const WMTSTileMatrix &oTM) const
      94             :     {
      95           4 :         OGREnvelope sExtent;
      96           4 :         const double dfTileWidthUnits = oTM.dfPixelSize * oTM.nTileWidth;
      97           4 :         const double dfTileHeightUnits = oTM.dfPixelSize * oTM.nTileHeight;
      98           4 :         sExtent.MinX = oTM.dfTLX + nMinTileCol * dfTileWidthUnits;
      99           4 :         sExtent.MaxY = oTM.dfTLY - nMinTileRow * dfTileHeightUnits;
     100           4 :         sExtent.MaxX = oTM.dfTLX + (nMaxTileCol + 1) * dfTileWidthUnits;
     101           4 :         sExtent.MinY = oTM.dfTLY - (nMaxTileRow + 1) * dfTileHeightUnits;
     102           4 :         return sExtent;
     103             :     }
     104             : };
     105             : 
     106             : /************************************************************************/
     107             : /* ==================================================================== */
     108             : /*                          WMTSTileMatrixSet                           */
     109             : /* ==================================================================== */
     110             : /************************************************************************/
     111             : 
     112             : class WMTSTileMatrixSet
     113             : {
     114             :   public:
     115             :     OGRSpatialReference oSRS{};
     116             :     CPLString osSRS{};
     117             :     bool bBoundingBoxValid = false;
     118             :     OGREnvelope sBoundingBox{}; /* expressed in TMS SRS */
     119             :     std::vector<WMTSTileMatrix> aoTM{};
     120             : 
     121         110 :     WMTSTileMatrixSet()
     122         110 :     {
     123         110 :         oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
     124         110 :     }
     125             : };
     126             : 
     127             : /************************************************************************/
     128             : /* ==================================================================== */
     129             : /*                              WMTSDataset                             */
     130             : /* ==================================================================== */
     131             : /************************************************************************/
     132             : 
     133             : class WMTSDataset final : public GDALPamDataset
     134             : {
     135             :     friend class WMTSBand;
     136             : 
     137             :     CPLString osLayer{};
     138             :     CPLString osTMS{};
     139             :     CPLString osXML{};
     140             :     CPLString osURLFeatureInfoTemplate{};
     141             :     WMTSTileMatrixSet oTMS{};
     142             : 
     143             :     CPLStringList m_aosHTTPOptions{};
     144             : 
     145             :     std::vector<GDALDataset *> apoDatasets{};
     146             :     OGRSpatialReference m_oSRS{};
     147             :     GDALGeoTransform m_gt{};
     148             : 
     149             :     CPLString osLastGetFeatureInfoURL{};
     150             :     CPLString osMetadataItemGetFeatureInfo{};
     151             : 
     152             :     static CPLStringList BuildHTTPRequestOpts(CPLString osOtherXML);
     153             :     static CPLXMLNode *GetCapabilitiesResponse(const CPLString &osFilename,
     154             :                                                CSLConstList papszHTTPOptions);
     155             :     static CPLString FixCRSName(const char *pszCRS);
     156             :     static CPLString Replace(const CPLString &osStr, const char *pszOld,
     157             :                              const char *pszNew);
     158             :     static CPLString GetOperationKVPURL(CPLXMLNode *psXML,
     159             :                                         const char *pszOperation);
     160             :     static int ReadTMS(CPLXMLNode *psContents, const CPLString &osIdentifier,
     161             :                        const CPLString &osMaxTileMatrixIdentifier,
     162             :                        int nMaxZoomLevel, WMTSTileMatrixSet &oTMS,
     163             :                        bool &bHasWarnedAutoSwap);
     164             :     static int ReadTMLimits(
     165             :         CPLXMLNode *psTMSLimits,
     166             :         std::map<CPLString, WMTSTileMatrixLimits> &aoMapTileMatrixLimits);
     167             : 
     168             :   public:
     169             :     WMTSDataset();
     170             :     ~WMTSDataset() override;
     171             : 
     172             :     CPLErr GetGeoTransform(GDALGeoTransform &gt) const override;
     173             :     const OGRSpatialReference *GetSpatialRef() const override;
     174             :     virtual const char *GetMetadataItem(const char *pszName,
     175             :                                         const char *pszDomain) override;
     176             : 
     177             :     static GDALDataset *Open(GDALOpenInfo *);
     178             :     static GDALDataset *CreateCopy(const char *pszFilename,
     179             :                                    GDALDataset *poSrcDS, CPL_UNUSED int bStrict,
     180             :                                    CPL_UNUSED char **papszOptions,
     181             :                                    CPL_UNUSED GDALProgressFunc pfnProgress,
     182             :                                    CPL_UNUSED void *pProgressData);
     183             : 
     184             :   protected:
     185             :     int CloseDependentDatasets() override;
     186             : 
     187             :     CPLErr IRasterIO(GDALRWFlag eRWFlag, int nXOff, int nYOff, int nXSize,
     188             :                      int nYSize, void *pData, int nBufXSize, int nBufYSize,
     189             :                      GDALDataType eBufType, int nBandCount,
     190             :                      BANDMAP_TYPE panBandMap, GSpacing nPixelSpace,
     191             :                      GSpacing nLineSpace, GSpacing nBandSpace,
     192             :                      GDALRasterIOExtraArg *psExtraArg) override;
     193             : };
     194             : 
     195             : /************************************************************************/
     196             : /* ==================================================================== */
     197             : /*                               WMTSBand                               */
     198             : /* ==================================================================== */
     199             : /************************************************************************/
     200             : 
     201             : class WMTSBand final : public GDALPamRasterBand
     202             : {
     203             :   public:
     204             :     WMTSBand(WMTSDataset *poDS, int nBand, GDALDataType eDataType);
     205             : 
     206             :     GDALRasterBand *GetOverview(int nLevel) override;
     207             :     int GetOverviewCount() override;
     208             :     GDALColorInterp GetColorInterpretation() override;
     209             :     virtual const char *GetMetadataItem(const char *pszName,
     210             :                                         const char *pszDomain) override;
     211             : 
     212             :   protected:
     213             :     CPLErr IReadBlock(int nBlockXOff, int nBlockYOff, void *pImage) override;
     214             :     CPLErr IRasterIO(GDALRWFlag, int, int, int, int, void *, int, int,
     215             :                      GDALDataType, GSpacing, GSpacing,
     216             :                      GDALRasterIOExtraArg *psExtraArg) override;
     217             : };
     218             : 
     219             : /************************************************************************/
     220             : /*                            WMTSBand()                                */
     221             : /************************************************************************/
     222             : 
     223         188 : WMTSBand::WMTSBand(WMTSDataset *poDSIn, int nBandIn, GDALDataType eDataTypeIn)
     224             : {
     225         188 :     poDS = poDSIn;
     226         188 :     nBand = nBandIn;
     227         188 :     eDataType = eDataTypeIn;
     228         188 :     poDSIn->apoDatasets[0]->GetRasterBand(1)->GetBlockSize(&nBlockXSize,
     229             :                                                            &nBlockYSize);
     230         188 : }
     231             : 
     232             : /************************************************************************/
     233             : /*                            IReadBlock()                              */
     234             : /************************************************************************/
     235             : 
     236           0 : CPLErr WMTSBand::IReadBlock(int nBlockXOff, int nBlockYOff, void *pImage)
     237             : {
     238           0 :     WMTSDataset *poGDS = cpl::down_cast<WMTSDataset *>(poDS);
     239           0 :     return poGDS->apoDatasets[0]->GetRasterBand(nBand)->ReadBlock(
     240           0 :         nBlockXOff, nBlockYOff, pImage);
     241             : }
     242             : 
     243             : /************************************************************************/
     244             : /*                             IRasterIO()                              */
     245             : /************************************************************************/
     246             : 
     247          35 : CPLErr WMTSBand::IRasterIO(GDALRWFlag eRWFlag, int nXOff, int nYOff, int nXSize,
     248             :                            int nYSize, void *pData, int nBufXSize,
     249             :                            int nBufYSize, GDALDataType eBufType,
     250             :                            GSpacing nPixelSpace, GSpacing nLineSpace,
     251             :                            GDALRasterIOExtraArg *psExtraArg)
     252             : {
     253          35 :     WMTSDataset *poGDS = cpl::down_cast<WMTSDataset *>(poDS);
     254             : 
     255          35 :     if ((nBufXSize < nXSize || nBufYSize < nYSize) &&
     256          70 :         poGDS->apoDatasets.size() > 1 && eRWFlag == GF_Read)
     257             :     {
     258             :         int bTried;
     259           1 :         CPLErr eErr = TryOverviewRasterIO(
     260             :             eRWFlag, nXOff, nYOff, nXSize, nYSize, pData, nBufXSize, nBufYSize,
     261             :             eBufType, nPixelSpace, nLineSpace, psExtraArg, &bTried);
     262           1 :         if (bTried)
     263           1 :             return eErr;
     264             :     }
     265             : 
     266          34 :     return poGDS->apoDatasets[0]->GetRasterBand(nBand)->RasterIO(
     267             :         eRWFlag, nXOff, nYOff, nXSize, nYSize, pData, nBufXSize, nBufYSize,
     268          34 :         eBufType, nPixelSpace, nLineSpace, psExtraArg);
     269             : }
     270             : 
     271             : /************************************************************************/
     272             : /*                         GetOverviewCount()                           */
     273             : /************************************************************************/
     274             : 
     275          15 : int WMTSBand::GetOverviewCount()
     276             : {
     277          15 :     WMTSDataset *poGDS = cpl::down_cast<WMTSDataset *>(poDS);
     278             : 
     279          15 :     if (poGDS->apoDatasets.size() > 1)
     280          13 :         return static_cast<int>(poGDS->apoDatasets.size()) - 1;
     281             :     else
     282           2 :         return 0;
     283             : }
     284             : 
     285             : /************************************************************************/
     286             : /*                              GetOverview()                           */
     287             : /************************************************************************/
     288             : 
     289          11 : GDALRasterBand *WMTSBand::GetOverview(int nLevel)
     290             : {
     291          11 :     WMTSDataset *poGDS = cpl::down_cast<WMTSDataset *>(poDS);
     292             : 
     293          11 :     if (nLevel < 0 || nLevel >= GetOverviewCount())
     294           1 :         return nullptr;
     295             : 
     296          10 :     GDALDataset *poOvrDS = poGDS->apoDatasets[nLevel + 1];
     297          10 :     if (poOvrDS)
     298          10 :         return poOvrDS->GetRasterBand(nBand);
     299             :     else
     300           0 :         return nullptr;
     301             : }
     302             : 
     303             : /************************************************************************/
     304             : /*                   GetColorInterpretation()                           */
     305             : /************************************************************************/
     306             : 
     307           4 : GDALColorInterp WMTSBand::GetColorInterpretation()
     308             : {
     309           4 :     WMTSDataset *poGDS = cpl::down_cast<WMTSDataset *>(poDS);
     310           4 :     if (poGDS->nBands == 1)
     311             :     {
     312           0 :         return GCI_GrayIndex;
     313             :     }
     314           4 :     else if (poGDS->nBands == 3 || poGDS->nBands == 4)
     315             :     {
     316           4 :         if (nBand == 1)
     317           1 :             return GCI_RedBand;
     318           3 :         else if (nBand == 2)
     319           1 :             return GCI_GreenBand;
     320           2 :         else if (nBand == 3)
     321           1 :             return GCI_BlueBand;
     322           1 :         else if (nBand == 4)
     323           1 :             return GCI_AlphaBand;
     324             :     }
     325             : 
     326           0 :     return GCI_Undefined;
     327             : }
     328             : 
     329             : /************************************************************************/
     330             : /*                         GetMetadataItem()                            */
     331             : /************************************************************************/
     332             : 
     333           8 : const char *WMTSBand::GetMetadataItem(const char *pszName,
     334             :                                       const char *pszDomain)
     335             : {
     336           8 :     WMTSDataset *poGDS = cpl::down_cast<WMTSDataset *>(poDS);
     337             : 
     338             :     /* ==================================================================== */
     339             :     /*      LocationInfo handling.                                          */
     340             :     /* ==================================================================== */
     341           8 :     if (pszDomain != nullptr && EQUAL(pszDomain, "LocationInfo") &&
     342           7 :         pszName != nullptr && STARTS_WITH_CI(pszName, "Pixel_") &&
     343          16 :         !poGDS->oTMS.aoTM.empty() && !poGDS->osURLFeatureInfoTemplate.empty())
     344             :     {
     345             :         int iPixel, iLine;
     346             : 
     347             :         /* --------------------------------------------------------------------
     348             :          */
     349             :         /*      What pixel are we aiming at? */
     350             :         /* --------------------------------------------------------------------
     351             :          */
     352           6 :         if (sscanf(pszName + 6, "%d_%d", &iPixel, &iLine) != 2)
     353           0 :             return nullptr;
     354             : 
     355           6 :         const WMTSTileMatrix &oTM = poGDS->oTMS.aoTM.back();
     356             : 
     357           6 :         iPixel += static_cast<int>(
     358           6 :             std::round((poGDS->m_gt[0] - oTM.dfTLX) / oTM.dfPixelSize));
     359           6 :         iLine += static_cast<int>(
     360           6 :             std::round((oTM.dfTLY - poGDS->m_gt[3]) / oTM.dfPixelSize));
     361             : 
     362          12 :         CPLString osURL(poGDS->osURLFeatureInfoTemplate);
     363           6 :         osURL = WMTSDataset::Replace(osURL, "{TileMatrixSet}", poGDS->osTMS);
     364           6 :         osURL = WMTSDataset::Replace(osURL, "{TileMatrix}", oTM.osIdentifier);
     365          12 :         osURL = WMTSDataset::Replace(osURL, "{TileCol}",
     366          12 :                                      CPLSPrintf("%d", iPixel / oTM.nTileWidth));
     367          12 :         osURL = WMTSDataset::Replace(osURL, "{TileRow}",
     368          12 :                                      CPLSPrintf("%d", iLine / oTM.nTileHeight));
     369          12 :         osURL = WMTSDataset::Replace(osURL, "{I}",
     370          12 :                                      CPLSPrintf("%d", iPixel % oTM.nTileWidth));
     371          12 :         osURL = WMTSDataset::Replace(osURL, "{J}",
     372          12 :                                      CPLSPrintf("%d", iLine % oTM.nTileHeight));
     373             : 
     374           6 :         if (poGDS->osLastGetFeatureInfoURL.compare(osURL) != 0)
     375             :         {
     376           5 :             poGDS->osLastGetFeatureInfoURL = osURL;
     377           5 :             poGDS->osMetadataItemGetFeatureInfo = "";
     378           5 :             char *pszRes = nullptr;
     379             :             CPLHTTPResult *psResult =
     380           5 :                 CPLHTTPFetch(osURL, poGDS->m_aosHTTPOptions.List());
     381           5 :             if (psResult && psResult->nStatus == 0 && psResult->pabyData)
     382           3 :                 pszRes = CPLStrdup(
     383           3 :                     reinterpret_cast<const char *>(psResult->pabyData));
     384           5 :             CPLHTTPDestroyResult(psResult);
     385             : 
     386           5 :             if (pszRes)
     387             :             {
     388           3 :                 poGDS->osMetadataItemGetFeatureInfo = "<LocationInfo>";
     389           3 :                 CPLPushErrorHandler(CPLQuietErrorHandler);
     390           3 :                 CPLXMLNode *psXML = CPLParseXMLString(pszRes);
     391           3 :                 CPLPopErrorHandler();
     392           3 :                 if (psXML != nullptr && psXML->eType == CXT_Element)
     393             :                 {
     394           1 :                     if (strcmp(psXML->pszValue, "?xml") == 0)
     395             :                     {
     396           1 :                         if (psXML->psNext)
     397             :                         {
     398           1 :                             char *pszXML = CPLSerializeXMLTree(psXML->psNext);
     399           1 :                             poGDS->osMetadataItemGetFeatureInfo += pszXML;
     400           1 :                             CPLFree(pszXML);
     401             :                         }
     402             :                     }
     403             :                     else
     404             :                     {
     405           0 :                         poGDS->osMetadataItemGetFeatureInfo += pszRes;
     406           1 :                     }
     407             :                 }
     408             :                 else
     409             :                 {
     410             :                     char *pszEscapedXML =
     411           2 :                         CPLEscapeString(pszRes, -1, CPLES_XML_BUT_QUOTES);
     412           2 :                     poGDS->osMetadataItemGetFeatureInfo += pszEscapedXML;
     413           2 :                     CPLFree(pszEscapedXML);
     414             :                 }
     415           3 :                 if (psXML != nullptr)
     416           3 :                     CPLDestroyXMLNode(psXML);
     417             : 
     418           3 :                 poGDS->osMetadataItemGetFeatureInfo += "</LocationInfo>";
     419           3 :                 CPLFree(pszRes);
     420             :             }
     421             :         }
     422           6 :         return poGDS->osMetadataItemGetFeatureInfo.c_str();
     423             :     }
     424             : 
     425           2 :     return GDALPamRasterBand::GetMetadataItem(pszName, pszDomain);
     426             : }
     427             : 
     428             : /************************************************************************/
     429             : /*                          WMTSDataset()                               */
     430             : /************************************************************************/
     431             : 
     432          55 : WMTSDataset::WMTSDataset()
     433             : {
     434          55 :     m_oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
     435          55 : }
     436             : 
     437             : /************************************************************************/
     438             : /*                        ~WMTSDataset()                                */
     439             : /************************************************************************/
     440             : 
     441         110 : WMTSDataset::~WMTSDataset()
     442             : {
     443          55 :     WMTSDataset::CloseDependentDatasets();
     444         110 : }
     445             : 
     446             : /************************************************************************/
     447             : /*                      CloseDependentDatasets()                        */
     448             : /************************************************************************/
     449             : 
     450          55 : int WMTSDataset::CloseDependentDatasets()
     451             : {
     452          55 :     int bRet = GDALPamDataset::CloseDependentDatasets();
     453          55 :     if (!apoDatasets.empty())
     454             :     {
     455         147 :         for (size_t i = 0; i < apoDatasets.size(); i++)
     456         100 :             delete apoDatasets[i];
     457          47 :         apoDatasets.resize(0);
     458          47 :         bRet = TRUE;
     459             :     }
     460          55 :     return bRet;
     461             : }
     462             : 
     463             : /************************************************************************/
     464             : /*                             IRasterIO()                              */
     465             : /************************************************************************/
     466             : 
     467           2 : CPLErr WMTSDataset::IRasterIO(GDALRWFlag eRWFlag, int nXOff, int nYOff,
     468             :                               int nXSize, int nYSize, void *pData,
     469             :                               int nBufXSize, int nBufYSize,
     470             :                               GDALDataType eBufType, int nBandCount,
     471             :                               BANDMAP_TYPE panBandMap, GSpacing nPixelSpace,
     472             :                               GSpacing nLineSpace, GSpacing nBandSpace,
     473             :                               GDALRasterIOExtraArg *psExtraArg)
     474             : {
     475           2 :     if ((nBufXSize < nXSize || nBufYSize < nYSize) && apoDatasets.size() > 1 &&
     476             :         eRWFlag == GF_Read)
     477             :     {
     478             :         int bTried;
     479           1 :         CPLErr eErr = TryOverviewRasterIO(
     480             :             eRWFlag, nXOff, nYOff, nXSize, nYSize, pData, nBufXSize, nBufYSize,
     481             :             eBufType, nBandCount, panBandMap, nPixelSpace, nLineSpace,
     482             :             nBandSpace, psExtraArg, &bTried);
     483           1 :         if (bTried)
     484           1 :             return eErr;
     485             :     }
     486             : 
     487           1 :     return apoDatasets[0]->RasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize,
     488             :                                     pData, nBufXSize, nBufYSize, eBufType,
     489             :                                     nBandCount, panBandMap, nPixelSpace,
     490           1 :                                     nLineSpace, nBandSpace, psExtraArg);
     491             : }
     492             : 
     493             : /************************************************************************/
     494             : /*                          GetGeoTransform()                           */
     495             : /************************************************************************/
     496             : 
     497           9 : CPLErr WMTSDataset::GetGeoTransform(GDALGeoTransform &gt) const
     498             : {
     499           9 :     gt = m_gt;
     500           9 :     return CE_None;
     501             : }
     502             : 
     503             : /************************************************************************/
     504             : /*                         GetSpatialRef()                              */
     505             : /************************************************************************/
     506             : 
     507           9 : const OGRSpatialReference *WMTSDataset::GetSpatialRef() const
     508             : {
     509           9 :     return m_oSRS.IsEmpty() ? nullptr : &m_oSRS;
     510             : }
     511             : 
     512             : /************************************************************************/
     513             : /*                          WMTSEscapeXML()                             */
     514             : /************************************************************************/
     515             : 
     516         266 : static CPLString WMTSEscapeXML(const char *pszUnescapedXML)
     517             : {
     518         266 :     CPLString osRet;
     519         266 :     char *pszTmp = CPLEscapeString(pszUnescapedXML, -1, CPLES_XML);
     520         266 :     osRet = pszTmp;
     521         266 :     CPLFree(pszTmp);
     522         266 :     return osRet;
     523             : }
     524             : 
     525             : /************************************************************************/
     526             : /*                         GetMetadataItem()                            */
     527             : /************************************************************************/
     528             : 
     529           4 : const char *WMTSDataset::GetMetadataItem(const char *pszName,
     530             :                                          const char *pszDomain)
     531             : {
     532           4 :     if (pszName != nullptr && EQUAL(pszName, "XML") && pszDomain != nullptr &&
     533           2 :         EQUAL(pszDomain, "WMTS"))
     534             :     {
     535           2 :         return osXML.c_str();
     536             :     }
     537             : 
     538           2 :     return GDALPamDataset::GetMetadataItem(pszName, pszDomain);
     539             : }
     540             : 
     541             : /************************************************************************/
     542             : /*                          QuoteIfNecessary()                          */
     543             : /************************************************************************/
     544             : 
     545         121 : static CPLString QuoteIfNecessary(const char *pszVal)
     546             : {
     547         121 :     if (strchr(pszVal, ' ') || strchr(pszVal, ',') || strchr(pszVal, '='))
     548             :     {
     549          48 :         CPLString osVal;
     550          24 :         osVal += "\"";
     551          24 :         osVal += pszVal;
     552          24 :         osVal += "\"";
     553          24 :         return osVal;
     554             :     }
     555             :     else
     556          97 :         return pszVal;
     557             : }
     558             : 
     559             : /************************************************************************/
     560             : /*                             FixCRSName()                             */
     561             : /************************************************************************/
     562             : 
     563         104 : CPLString WMTSDataset::FixCRSName(const char *pszCRS)
     564             : {
     565         104 :     while (*pszCRS == ' ' || *pszCRS == '\r' || *pszCRS == '\n')
     566           0 :         pszCRS++;
     567             : 
     568             :     /* http://maps.wien.gv.at/wmts/1.0.0/WMTSCapabilities.xml uses
     569             :      * urn:ogc:def:crs:EPSG:6.18:3:3857 */
     570             :     /* instead of urn:ogc:def:crs:EPSG:6.18.3:3857. Coming from an incorrect
     571             :      * example of URN in WMTS spec */
     572             :     /* https://portal.opengeospatial.org/files/?artifact_id=50398 */
     573         104 :     if (STARTS_WITH_CI(pszCRS, "urn:ogc:def:crs:EPSG:6.18:3:"))
     574             :     {
     575             :         return CPLSPrintf("urn:ogc:def:crs:EPSG::%s",
     576          43 :                           pszCRS + strlen("urn:ogc:def:crs:EPSG:6.18:3:"));
     577             :     }
     578             : 
     579          61 :     if (EQUAL(pszCRS, "urn:ogc:def:crs:EPSG::102100"))
     580           0 :         return "EPSG:3857";
     581             : 
     582         122 :     CPLString osRet(pszCRS);
     583         122 :     while (osRet.size() && (osRet.back() == ' ' || osRet.back() == '\r' ||
     584          61 :                             osRet.back() == '\n'))
     585             :     {
     586           0 :         osRet.pop_back();
     587             :     }
     588          61 :     return osRet;
     589             : }
     590             : 
     591             : /************************************************************************/
     592             : /*                              ReadTMS()                               */
     593             : /************************************************************************/
     594             : 
     595          55 : int WMTSDataset::ReadTMS(CPLXMLNode *psContents, const CPLString &osIdentifier,
     596             :                          const CPLString &osMaxTileMatrixIdentifier,
     597             :                          int nMaxZoomLevel, WMTSTileMatrixSet &oTMS,
     598             :                          bool &bHasWarnedAutoSwap)
     599             : {
     600         110 :     for (CPLXMLNode *psIter = psContents->psChild; psIter != nullptr;
     601          55 :          psIter = psIter->psNext)
     602             :     {
     603         110 :         if (psIter->eType != CXT_Element ||
     604         110 :             strcmp(psIter->pszValue, "TileMatrixSet") != 0)
     605          55 :             continue;
     606          55 :         const char *pszIdentifier = CPLGetXMLValue(psIter, "Identifier", "");
     607          55 :         if (!EQUAL(osIdentifier, pszIdentifier))
     608           0 :             continue;
     609             :         const char *pszSupportedCRS =
     610          55 :             CPLGetXMLValue(psIter, "SupportedCRS", nullptr);
     611          55 :         if (pszSupportedCRS == nullptr)
     612             :         {
     613           1 :             CPLError(CE_Failure, CPLE_AppDefined, "Missing SupportedCRS");
     614           1 :             return FALSE;
     615             :         }
     616          54 :         oTMS.osSRS = pszSupportedCRS;
     617         108 :         if (oTMS.oSRS.SetFromUserInput(
     618         108 :                 FixCRSName(pszSupportedCRS),
     619          54 :                 OGRSpatialReference::SET_FROM_USER_INPUT_LIMITATIONS_get()) !=
     620             :             OGRERR_NONE)
     621             :         {
     622           0 :             CPLError(CE_Failure, CPLE_AppDefined, "Cannot parse CRS '%s'",
     623             :                      pszSupportedCRS);
     624           0 :             return FALSE;
     625             :         }
     626             :         const bool bSwap =
     627         154 :             !STARTS_WITH_CI(pszSupportedCRS, "EPSG:") &&
     628         100 :             (CPL_TO_BOOL(oTMS.oSRS.EPSGTreatsAsLatLong()) ||
     629          47 :              CPL_TO_BOOL(oTMS.oSRS.EPSGTreatsAsNorthingEasting()));
     630          54 :         CPLXMLNode *psBB = CPLGetXMLNode(psIter, "BoundingBox");
     631          54 :         oTMS.bBoundingBoxValid = false;
     632          54 :         if (psBB != nullptr)
     633             :         {
     634          18 :             CPLString osCRS = CPLGetXMLValue(psBB, "crs", "");
     635           9 :             if (EQUAL(osCRS, "") || EQUAL(osCRS, pszSupportedCRS))
     636             :             {
     637             :                 CPLString osLowerCorner =
     638          18 :                     CPLGetXMLValue(psBB, "LowerCorner", "");
     639             :                 CPLString osUpperCorner =
     640          18 :                     CPLGetXMLValue(psBB, "UpperCorner", "");
     641           9 :                 if (!osLowerCorner.empty() && !osUpperCorner.empty())
     642             :                 {
     643           9 :                     char **papszLC = CSLTokenizeString(osLowerCorner);
     644           9 :                     char **papszUC = CSLTokenizeString(osUpperCorner);
     645           9 :                     if (CSLCount(papszLC) == 2 && CSLCount(papszUC) == 2)
     646             :                     {
     647           9 :                         oTMS.sBoundingBox.MinX =
     648           9 :                             CPLAtof(papszLC[(bSwap) ? 1 : 0]);
     649           9 :                         oTMS.sBoundingBox.MinY =
     650           9 :                             CPLAtof(papszLC[(bSwap) ? 0 : 1]);
     651           9 :                         oTMS.sBoundingBox.MaxX =
     652           9 :                             CPLAtof(papszUC[(bSwap) ? 1 : 0]);
     653           9 :                         oTMS.sBoundingBox.MaxY =
     654           9 :                             CPLAtof(papszUC[(bSwap) ? 0 : 1]);
     655           9 :                         oTMS.bBoundingBoxValid = true;
     656             :                     }
     657           9 :                     CSLDestroy(papszLC);
     658           9 :                     CSLDestroy(papszUC);
     659             :                 }
     660             :             }
     661             :         }
     662             :         else
     663             :         {
     664             :             const char *pszWellKnownScaleSet =
     665          45 :                 CPLGetXMLValue(psIter, "WellKnownScaleSet", "");
     666          45 :             if (EQUAL(pszIdentifier, "GoogleCRS84Quad") ||
     667          45 :                 EQUAL(pszWellKnownScaleSet,
     668          45 :                       "urn:ogc:def:wkss:OGC:1.0:GoogleCRS84Quad") ||
     669          45 :                 EQUAL(pszIdentifier, "GlobalCRS84Scale") ||
     670          45 :                 EQUAL(pszWellKnownScaleSet,
     671             :                       "urn:ogc:def:wkss:OGC:1.0:GlobalCRS84Scale"))
     672             :             {
     673           0 :                 oTMS.sBoundingBox.MinX = -180;
     674           0 :                 oTMS.sBoundingBox.MinY = -90;
     675           0 :                 oTMS.sBoundingBox.MaxX = 180;
     676           0 :                 oTMS.sBoundingBox.MaxY = 90;
     677           0 :                 oTMS.bBoundingBoxValid = true;
     678             :             }
     679             :         }
     680             : 
     681          54 :         bool bFoundTileMatrix = false;
     682         308 :         for (CPLXMLNode *psSubIter = psIter->psChild; psSubIter != nullptr;
     683         254 :              psSubIter = psSubIter->psNext)
     684             :         {
     685         262 :             if (psSubIter->eType != CXT_Element ||
     686         262 :                 strcmp(psSubIter->pszValue, "TileMatrix") != 0)
     687         129 :                 continue;
     688             :             const char *l_pszIdentifier =
     689         133 :                 CPLGetXMLValue(psSubIter, "Identifier", nullptr);
     690             :             const char *pszScaleDenominator =
     691         133 :                 CPLGetXMLValue(psSubIter, "ScaleDenominator", nullptr);
     692             :             const char *pszTopLeftCorner =
     693         133 :                 CPLGetXMLValue(psSubIter, "TopLeftCorner", nullptr);
     694             :             const char *pszTileWidth =
     695         133 :                 CPLGetXMLValue(psSubIter, "TileWidth", nullptr);
     696             :             const char *pszTileHeight =
     697         133 :                 CPLGetXMLValue(psSubIter, "TileHeight", nullptr);
     698             :             const char *pszMatrixWidth =
     699         133 :                 CPLGetXMLValue(psSubIter, "MatrixWidth", nullptr);
     700             :             const char *pszMatrixHeight =
     701         133 :                 CPLGetXMLValue(psSubIter, "MatrixHeight", nullptr);
     702         133 :             if (l_pszIdentifier == nullptr || pszScaleDenominator == nullptr ||
     703         132 :                 pszTopLeftCorner == nullptr ||
     704         132 :                 strchr(pszTopLeftCorner, ' ') == nullptr ||
     705         132 :                 pszTileWidth == nullptr || pszTileHeight == nullptr ||
     706         132 :                 pszMatrixWidth == nullptr || pszMatrixHeight == nullptr)
     707             :             {
     708           1 :                 CPLError(CE_Failure, CPLE_AppDefined,
     709             :                          "Missing required element in TileMatrix element");
     710           1 :                 return FALSE;
     711             :             }
     712         132 :             WMTSTileMatrix oTM;
     713         132 :             oTM.osIdentifier = l_pszIdentifier;
     714         132 :             oTM.dfScaleDenominator = CPLAtof(pszScaleDenominator);
     715         132 :             oTM.dfPixelSize = oTM.dfScaleDenominator * WMTS_PITCH;
     716         132 :             if (oTM.dfPixelSize <= 0.0)
     717             :             {
     718           0 :                 CPLError(CE_Failure, CPLE_AppDefined,
     719             :                          "Invalid ScaleDenominator");
     720           0 :                 return FALSE;
     721             :             }
     722         132 :             if (oTMS.oSRS.IsGeographic())
     723          18 :                 oTM.dfPixelSize *= WMTS_WGS84_DEG_PER_METER;
     724         132 :             double dfVal1 = CPLAtof(pszTopLeftCorner);
     725         132 :             double dfVal2 = CPLAtof(strchr(pszTopLeftCorner, ' ') + 1);
     726         132 :             if (!bSwap)
     727             :             {
     728         113 :                 oTM.dfTLX = dfVal1;
     729         113 :                 oTM.dfTLY = dfVal2;
     730             :             }
     731             :             else
     732             :             {
     733          19 :                 oTM.dfTLX = dfVal2;
     734          19 :                 oTM.dfTLY = dfVal1;
     735             :             }
     736             : 
     737             :             // Hack for http://osm.geobretagne.fr/gwc01/service/wmts?request=getcapabilities
     738             :             // or https://trek.nasa.gov/tiles/Mars/EQ/THEMIS_NightIR_ControlledMosaics_100m_v2_oct2018/1.0.0/WMTSCapabilities.xml
     739         132 :             if (oTM.dfTLY == -180.0 &&
     740           0 :                 (STARTS_WITH_CI(l_pszIdentifier, "EPSG:4326:") ||
     741           0 :                  (oTMS.oSRS.IsGeographic() && oTM.dfTLX == 90)))
     742             :             {
     743           0 :                 if (!bHasWarnedAutoSwap)
     744             :                 {
     745           0 :                     bHasWarnedAutoSwap = true;
     746           0 :                     CPLError(CE_Warning, CPLE_AppDefined,
     747             :                              "Auto-correcting wrongly swapped "
     748             :                              "TileMatrix.TopLeftCorner coordinates. "
     749             :                              "They should be in latitude, longitude order "
     750             :                              "but are presented in longitude, latitude order. "
     751             :                              "This should be reported to the server "
     752             :                              "administrator.");
     753             :                 }
     754           0 :                 std::swap(oTM.dfTLX, oTM.dfTLY);
     755             :             }
     756             : 
     757             :             // Hack for "https://t0.tianditu.gov.cn/img_w/wmts?SERVICE=WMTS&REQUEST=GetCapabilities&VERSION=1.0.0&tk=ec899a50c7830ea2416ca182285236f3"
     758             :             // which returns swapped coordinates for WebMercator
     759         132 :             if (std::fabs(oTM.dfTLX - 20037508.3427892) < 1e-4 &&
     760           0 :                 std::fabs(oTM.dfTLY - (-20037508.3427892)) < 1e-4)
     761             :             {
     762           0 :                 if (!bHasWarnedAutoSwap)
     763             :                 {
     764           0 :                     bHasWarnedAutoSwap = true;
     765           0 :                     CPLError(CE_Warning, CPLE_AppDefined,
     766             :                              "Auto-correcting wrongly swapped "
     767             :                              "TileMatrix.TopLeftCorner coordinates. This "
     768             :                              "should be reported to the server administrator.");
     769             :                 }
     770           0 :                 std::swap(oTM.dfTLX, oTM.dfTLY);
     771             :             }
     772             : 
     773         132 :             oTM.nTileWidth = atoi(pszTileWidth);
     774         132 :             oTM.nTileHeight = atoi(pszTileHeight);
     775         132 :             if (oTM.nTileWidth <= 0 || oTM.nTileWidth > 4096 ||
     776         132 :                 oTM.nTileHeight <= 0 || oTM.nTileHeight > 4096)
     777             :             {
     778           0 :                 CPLError(CE_Failure, CPLE_AppDefined,
     779             :                          "Invalid TileWidth/TileHeight element");
     780           0 :                 return FALSE;
     781             :             }
     782         132 :             oTM.nMatrixWidth = atoi(pszMatrixWidth);
     783         132 :             oTM.nMatrixHeight = atoi(pszMatrixHeight);
     784             :             // http://datacarto.geonormandie.fr/mapcache/wmts?SERVICE=WMTS&REQUEST=GetCapabilities
     785             :             // has a TileMatrix 0 with MatrixWidth = MatrixHeight = 0
     786         132 :             if (oTM.nMatrixWidth < 1 || oTM.nMatrixHeight < 1)
     787           0 :                 continue;
     788         132 :             oTMS.aoTM.push_back(std::move(oTM));
     789         140 :             if ((nMaxZoomLevel >= 0 &&
     790         260 :                  static_cast<int>(oTMS.aoTM.size()) - 1 == nMaxZoomLevel) ||
     791         128 :                 (!osMaxTileMatrixIdentifier.empty() &&
     792           7 :                  EQUAL(osMaxTileMatrixIdentifier, l_pszIdentifier)))
     793             :             {
     794           7 :                 bFoundTileMatrix = true;
     795           7 :                 break;
     796             :             }
     797             :         }
     798          53 :         if (nMaxZoomLevel >= 0 && !bFoundTileMatrix)
     799             :         {
     800           2 :             CPLError(
     801             :                 CE_Failure, CPLE_AppDefined,
     802             :                 "Cannot find TileMatrix of zoom level %d in TileMatrixSet '%s'",
     803             :                 nMaxZoomLevel, osIdentifier.c_str());
     804           2 :             return FALSE;
     805             :         }
     806          51 :         if (!osMaxTileMatrixIdentifier.empty() && !bFoundTileMatrix)
     807             :         {
     808           2 :             CPLError(CE_Failure, CPLE_AppDefined,
     809             :                      "Cannot find TileMatrix '%s' in TileMatrixSet '%s'",
     810             :                      osMaxTileMatrixIdentifier.c_str(), osIdentifier.c_str());
     811           2 :             return FALSE;
     812             :         }
     813          49 :         if (oTMS.aoTM.empty())
     814             :         {
     815           1 :             CPLError(CE_Failure, CPLE_AppDefined,
     816             :                      "Cannot find TileMatrix in TileMatrixSet '%s'",
     817             :                      osIdentifier.c_str());
     818           1 :             return FALSE;
     819             :         }
     820          48 :         return TRUE;
     821             :     }
     822           0 :     CPLError(CE_Failure, CPLE_AppDefined, "Cannot find TileMatrixSet '%s'",
     823             :              osIdentifier.c_str());
     824           0 :     return FALSE;
     825             : }
     826             : 
     827             : /************************************************************************/
     828             : /*                              ReadTMLimits()                          */
     829             : /************************************************************************/
     830             : 
     831           2 : int WMTSDataset::ReadTMLimits(
     832             :     CPLXMLNode *psTMSLimits,
     833             :     std::map<CPLString, WMTSTileMatrixLimits> &aoMapTileMatrixLimits)
     834             : {
     835           4 :     for (CPLXMLNode *psIter = psTMSLimits->psChild; psIter;
     836           2 :          psIter = psIter->psNext)
     837             :     {
     838           2 :         if (psIter->eType != CXT_Element ||
     839           2 :             strcmp(psIter->pszValue, "TileMatrixLimits") != 0)
     840           0 :             continue;
     841           2 :         WMTSTileMatrixLimits oTMLimits;
     842             :         const char *pszTileMatrix =
     843           2 :             CPLGetXMLValue(psIter, "TileMatrix", nullptr);
     844             :         const char *pszMinTileRow =
     845           2 :             CPLGetXMLValue(psIter, "MinTileRow", nullptr);
     846             :         const char *pszMaxTileRow =
     847           2 :             CPLGetXMLValue(psIter, "MaxTileRow", nullptr);
     848             :         const char *pszMinTileCol =
     849           2 :             CPLGetXMLValue(psIter, "MinTileCol", nullptr);
     850             :         const char *pszMaxTileCol =
     851           2 :             CPLGetXMLValue(psIter, "MaxTileCol", nullptr);
     852           2 :         if (pszTileMatrix == nullptr || pszMinTileRow == nullptr ||
     853           2 :             pszMaxTileRow == nullptr || pszMinTileCol == nullptr ||
     854             :             pszMaxTileCol == nullptr)
     855             :         {
     856           0 :             CPLError(CE_Failure, CPLE_AppDefined,
     857             :                      "Missing required element in TileMatrixLimits element");
     858           0 :             return FALSE;
     859             :         }
     860           2 :         oTMLimits.osIdentifier = pszTileMatrix;
     861           2 :         oTMLimits.nMinTileRow = atoi(pszMinTileRow);
     862           2 :         oTMLimits.nMaxTileRow = atoi(pszMaxTileRow);
     863           2 :         oTMLimits.nMinTileCol = atoi(pszMinTileCol);
     864           2 :         oTMLimits.nMaxTileCol = atoi(pszMaxTileCol);
     865           2 :         aoMapTileMatrixLimits[pszTileMatrix] = std::move(oTMLimits);
     866             :     }
     867           2 :     return TRUE;
     868             : }
     869             : 
     870             : /************************************************************************/
     871             : /*                               Replace()                              */
     872             : /************************************************************************/
     873             : 
     874         380 : CPLString WMTSDataset::Replace(const CPLString &osStr, const char *pszOld,
     875             :                                const char *pszNew)
     876             : {
     877         380 :     size_t nPos = osStr.ifind(pszOld);
     878         380 :     if (nPos == std::string::npos)
     879          59 :         return osStr;
     880         642 :     CPLString osRet(osStr.substr(0, nPos));
     881         321 :     osRet += pszNew;
     882         321 :     osRet += osStr.substr(nPos + strlen(pszOld));
     883         321 :     return osRet;
     884             : }
     885             : 
     886             : /************************************************************************/
     887             : /*                       GetCapabilitiesResponse()                      */
     888             : /************************************************************************/
     889             : 
     890          72 : CPLXMLNode *WMTSDataset::GetCapabilitiesResponse(const CPLString &osFilename,
     891             :                                                  CSLConstList papszHTTPOptions)
     892             : {
     893             :     CPLXMLNode *psXML;
     894             :     VSIStatBufL sStat;
     895          72 :     if (VSIStatL(osFilename, &sStat) == 0)
     896          66 :         psXML = CPLParseXMLFile(osFilename);
     897             :     else
     898             :     {
     899           6 :         CPLHTTPResult *psResult = CPLHTTPFetch(osFilename, papszHTTPOptions);
     900           6 :         if (psResult == nullptr)
     901           0 :             return nullptr;
     902           6 :         if (psResult->pabyData == nullptr)
     903             :         {
     904           3 :             CPLHTTPDestroyResult(psResult);
     905           3 :             return nullptr;
     906             :         }
     907           6 :         psXML = CPLParseXMLString(
     908           3 :             reinterpret_cast<const char *>(psResult->pabyData));
     909           3 :         CPLHTTPDestroyResult(psResult);
     910             :     }
     911          69 :     return psXML;
     912             : }
     913             : 
     914             : /************************************************************************/
     915             : /*                          WMTSAddOtherXML()                           */
     916             : /************************************************************************/
     917             : 
     918         165 : static void WMTSAddOtherXML(CPLXMLNode *psRoot, const char *pszElement,
     919             :                             CPLString &osOtherXML)
     920             : {
     921         165 :     CPLXMLNode *psElement = CPLGetXMLNode(psRoot, pszElement);
     922         165 :     if (psElement)
     923             :     {
     924          60 :         CPLXMLNode *psNext = psElement->psNext;
     925          60 :         psElement->psNext = nullptr;
     926          60 :         char *pszTmp = CPLSerializeXMLTree(psElement);
     927          60 :         osOtherXML += pszTmp;
     928          60 :         CPLFree(pszTmp);
     929          60 :         psElement->psNext = psNext;
     930             :     }
     931         165 : }
     932             : 
     933             : /************************************************************************/
     934             : /*                          GetOperationKVPURL()                        */
     935             : /************************************************************************/
     936             : 
     937          68 : CPLString WMTSDataset::GetOperationKVPURL(CPLXMLNode *psXML,
     938             :                                           const char *pszOperation)
     939             : {
     940          68 :     CPLString osRet;
     941          68 :     CPLXMLNode *psOM = CPLGetXMLNode(psXML, "=Capabilities.OperationsMetadata");
     942         106 :     for (CPLXMLNode *psIter = psOM ? psOM->psChild : nullptr; psIter != nullptr;
     943          38 :          psIter = psIter->psNext)
     944             :     {
     945         139 :         if (psIter->eType != CXT_Element ||
     946          76 :             strcmp(psIter->pszValue, "Operation") != 0 ||
     947          38 :             !EQUAL(CPLGetXMLValue(psIter, "name", ""), pszOperation))
     948             :         {
     949          25 :             continue;
     950             :         }
     951          13 :         CPLXMLNode *psHTTP = CPLGetXMLNode(psIter, "DCP.HTTP");
     952          13 :         for (CPLXMLNode *psGet = psHTTP ? psHTTP->psChild : nullptr;
     953          26 :              psGet != nullptr; psGet = psGet->psNext)
     954             :         {
     955          13 :             if (psGet->eType != CXT_Element ||
     956          13 :                 strcmp(psGet->pszValue, "Get") != 0)
     957             :             {
     958           0 :                 continue;
     959             :             }
     960          13 :             if (!EQUAL(CPLGetXMLValue(psGet, "Constraint.AllowedValues.Value",
     961             :                                       "KVP"),
     962             :                        "KVP"))
     963           1 :                 continue;
     964          12 :             osRet = CPLGetXMLValue(psGet, "href", "");
     965             :         }
     966             :     }
     967          68 :     return osRet;
     968             : }
     969             : 
     970             : /************************************************************************/
     971             : /*                           BuildHTTPRequestOpts()                     */
     972             : /************************************************************************/
     973             : 
     974         127 : CPLStringList WMTSDataset::BuildHTTPRequestOpts(CPLString osOtherXML)
     975             : {
     976         127 :     osOtherXML = "<Root>" + osOtherXML + "</Root>";
     977         127 :     CPLXMLNode *psXML = CPLParseXMLString(osOtherXML);
     978         127 :     CPLStringList opts;
     979         635 :     for (const char *pszOptionName :
     980         762 :          {"Timeout", "UserAgent", "Accept", "Referer", "UserPwd"})
     981             :     {
     982         635 :         if (const char *pszVal = CPLGetXMLValue(psXML, pszOptionName, nullptr))
     983             :         {
     984           8 :             opts.SetNameValue(CPLString(pszOptionName).toupper(), pszVal);
     985             :         }
     986             :     }
     987         127 :     if (CPLTestBool(CPLGetXMLValue(psXML, "UnsafeSSL", "false")))
     988             :     {
     989         125 :         opts.SetNameValue("UNSAFESSL", "1");
     990             :     }
     991         127 :     CPLDestroyXMLNode(psXML);
     992         127 :     return opts;
     993             : }
     994             : 
     995             : /************************************************************************/
     996             : /*                                Open()                                */
     997             : /************************************************************************/
     998             : 
     999          71 : GDALDataset *WMTSDataset::Open(GDALOpenInfo *poOpenInfo)
    1000             : {
    1001          71 :     if (!WMTSDriverIdentify(poOpenInfo))
    1002           0 :         return nullptr;
    1003             : 
    1004          71 :     CPLXMLNode *psXML = nullptr;
    1005         142 :     CPLString osTileFormat;
    1006         142 :     CPLString osInfoFormat;
    1007             : 
    1008             :     CPLString osGetCapabilitiesURL =
    1009         142 :         CSLFetchNameValueDef(poOpenInfo->papszOpenOptions, "URL", "");
    1010             :     CPLString osLayer =
    1011         142 :         CSLFetchNameValueDef(poOpenInfo->papszOpenOptions, "LAYER", "");
    1012             :     CPLString osTMS =
    1013         142 :         CSLFetchNameValueDef(poOpenInfo->papszOpenOptions, "TILEMATRIXSET", "");
    1014             :     CPLString osMaxTileMatrixIdentifier =
    1015         142 :         CSLFetchNameValueDef(poOpenInfo->papszOpenOptions, "TILEMATRIX", "");
    1016         213 :     int nUserMaxZoomLevel = atoi(CSLFetchNameValueDef(
    1017          71 :         poOpenInfo->papszOpenOptions, "ZOOM_LEVEL",
    1018          71 :         CSLFetchNameValueDef(poOpenInfo->papszOpenOptions, "ZOOMLEVEL", "-1")));
    1019             :     CPLString osStyle =
    1020         142 :         CSLFetchNameValueDef(poOpenInfo->papszOpenOptions, "STYLE", "");
    1021             : 
    1022          71 :     int bExtendBeyondDateLine = CPLFetchBool(poOpenInfo->papszOpenOptions,
    1023          71 :                                              "EXTENDBEYONDDATELINE", false);
    1024             : 
    1025             :     CPLString osOtherXML =
    1026             :         "<Cache />"
    1027             :         "<UnsafeSSL>true</UnsafeSSL>"
    1028             :         "<ZeroBlockHttpCodes>204,404</ZeroBlockHttpCodes>"
    1029         142 :         "<ZeroBlockOnServerException>true</ZeroBlockOnServerException>";
    1030             : 
    1031          71 :     if (STARTS_WITH_CI(poOpenInfo->pszFilename, "WMTS:"))
    1032             :     {
    1033          56 :         char **papszTokens = CSLTokenizeString2(poOpenInfo->pszFilename + 5,
    1034             :                                                 ",", CSLT_HONOURSTRINGS);
    1035          56 :         if (papszTokens && papszTokens[0])
    1036             :         {
    1037          50 :             osGetCapabilitiesURL = papszTokens[0];
    1038          67 :             for (char **papszIter = papszTokens + 1; *papszIter; papszIter++)
    1039             :             {
    1040          17 :                 char *pszKey = nullptr;
    1041          17 :                 const char *pszValue = CPLParseNameValue(*papszIter, &pszKey);
    1042          17 :                 if (pszKey && pszValue)
    1043             :                 {
    1044          17 :                     if (EQUAL(pszKey, "layer"))
    1045           3 :                         osLayer = pszValue;
    1046          14 :                     else if (EQUAL(pszKey, "tilematrixset"))
    1047           4 :                         osTMS = pszValue;
    1048          10 :                     else if (EQUAL(pszKey, "tilematrix"))
    1049           3 :                         osMaxTileMatrixIdentifier = pszValue;
    1050           7 :                     else if (EQUAL(pszKey, "zoom_level") ||
    1051           4 :                              EQUAL(pszKey, "zoomlevel"))
    1052           3 :                         nUserMaxZoomLevel = atoi(pszValue);
    1053           4 :                     else if (EQUAL(pszKey, "style"))
    1054           3 :                         osStyle = pszValue;
    1055           1 :                     else if (EQUAL(pszKey, "extendbeyonddateline"))
    1056           1 :                         bExtendBeyondDateLine = CPLTestBool(pszValue);
    1057             :                     else
    1058           0 :                         CPLError(CE_Warning, CPLE_AppDefined,
    1059             :                                  "Unknown parameter: %s'", pszKey);
    1060             :                 }
    1061          17 :                 CPLFree(pszKey);
    1062             :             }
    1063             :         }
    1064          56 :         CSLDestroy(papszTokens);
    1065             : 
    1066          56 :         const CPLStringList aosHTTPOptions(BuildHTTPRequestOpts(osOtherXML));
    1067          56 :         psXML = GetCapabilitiesResponse(osGetCapabilitiesURL,
    1068             :                                         aosHTTPOptions.List());
    1069             :     }
    1070          17 :     else if (poOpenInfo->IsSingleAllowedDriver("WMTS") &&
    1071           2 :              (STARTS_WITH(poOpenInfo->pszFilename, "http://") ||
    1072           1 :               STARTS_WITH(poOpenInfo->pszFilename, "https://")))
    1073             :     {
    1074           1 :         const CPLStringList aosHTTPOptions(BuildHTTPRequestOpts(osOtherXML));
    1075           1 :         psXML = GetCapabilitiesResponse(poOpenInfo->pszFilename,
    1076             :                                         aosHTTPOptions.List());
    1077             :     }
    1078             : 
    1079          71 :     int bHasAOI = FALSE;
    1080          71 :     OGREnvelope sAOI;
    1081          71 :     int nBands = 4;
    1082          71 :     GDALDataType eDataType = GDT_Byte;
    1083         142 :     CPLString osProjection;
    1084         142 :     CPLString osExtraQueryParameters;
    1085             : 
    1086          53 :     if ((psXML != nullptr && CPLGetXMLNode(psXML, "=GDAL_WMTS") != nullptr) ||
    1087         183 :         STARTS_WITH_CI(poOpenInfo->pszFilename, "<GDAL_WMTS") ||
    1088          59 :         (poOpenInfo->nHeaderBytes > 0 &&
    1089          10 :          strstr(reinterpret_cast<const char *>(poOpenInfo->pabyHeader),
    1090             :                 "<GDAL_WMTS")))
    1091             :     {
    1092             :         CPLXMLNode *psGDALWMTS;
    1093          18 :         if (psXML != nullptr && CPLGetXMLNode(psXML, "=GDAL_WMTS") != nullptr)
    1094           8 :             psGDALWMTS = CPLCloneXMLTree(psXML);
    1095          10 :         else if (STARTS_WITH_CI(poOpenInfo->pszFilename, "<GDAL_WMTS"))
    1096           4 :             psGDALWMTS = CPLParseXMLString(poOpenInfo->pszFilename);
    1097             :         else
    1098           6 :             psGDALWMTS = CPLParseXMLFile(poOpenInfo->pszFilename);
    1099          18 :         if (psGDALWMTS == nullptr)
    1100           3 :             return nullptr;
    1101          17 :         CPLXMLNode *psRoot = CPLGetXMLNode(psGDALWMTS, "=GDAL_WMTS");
    1102          17 :         if (psRoot == nullptr)
    1103             :         {
    1104           1 :             CPLError(CE_Failure, CPLE_AppDefined,
    1105             :                      "Cannot find root <GDAL_WMTS>");
    1106           1 :             CPLDestroyXMLNode(psGDALWMTS);
    1107           1 :             return nullptr;
    1108             :         }
    1109          16 :         osGetCapabilitiesURL = CPLGetXMLValue(psRoot, "GetCapabilitiesUrl", "");
    1110          16 :         if (osGetCapabilitiesURL.empty())
    1111             :         {
    1112           1 :             CPLError(CE_Failure, CPLE_AppDefined,
    1113             :                      "Missing <GetCapabilitiesUrl>");
    1114           1 :             CPLDestroyXMLNode(psGDALWMTS);
    1115           1 :             return nullptr;
    1116             :         }
    1117             :         osExtraQueryParameters =
    1118          15 :             CPLGetXMLValue(psRoot, "ExtraQueryParameters", "");
    1119          15 :         if (!osExtraQueryParameters.empty() && osExtraQueryParameters[0] != '&')
    1120           0 :             osExtraQueryParameters = '&' + osExtraQueryParameters;
    1121             : 
    1122          15 :         osGetCapabilitiesURL += osExtraQueryParameters;
    1123             : 
    1124          15 :         osLayer = CPLGetXMLValue(psRoot, "Layer", osLayer);
    1125          15 :         osTMS = CPLGetXMLValue(psRoot, "TileMatrixSet", osTMS);
    1126             :         osMaxTileMatrixIdentifier =
    1127          15 :             CPLGetXMLValue(psRoot, "TileMatrix", osMaxTileMatrixIdentifier);
    1128          15 :         nUserMaxZoomLevel = atoi(CPLGetXMLValue(
    1129             :             psRoot, "ZoomLevel", CPLSPrintf("%d", nUserMaxZoomLevel)));
    1130          15 :         osStyle = CPLGetXMLValue(psRoot, "Style", osStyle);
    1131          15 :         osTileFormat = CPLGetXMLValue(psRoot, "Format", osTileFormat);
    1132          15 :         osInfoFormat = CPLGetXMLValue(psRoot, "InfoFormat", osInfoFormat);
    1133          15 :         osProjection = CPLGetXMLValue(psRoot, "Projection", osProjection);
    1134          15 :         bExtendBeyondDateLine = CPLTestBool(
    1135             :             CPLGetXMLValue(psRoot, "ExtendBeyondDateLine",
    1136             :                            (bExtendBeyondDateLine) ? "true" : "false"));
    1137             : 
    1138          15 :         osOtherXML = "";
    1139         165 :         for (const char *pszXMLElement :
    1140             :              {"Cache", "MaxConnections", "Timeout", "OfflineMode", "UserAgent",
    1141             :               "Accept", "UserPwd", "UnsafeSSL", "Referer", "ZeroBlockHttpCodes",
    1142         180 :               "ZeroBlockOnServerException"})
    1143             :         {
    1144         165 :             WMTSAddOtherXML(psRoot, pszXMLElement, osOtherXML);
    1145             :         }
    1146             : 
    1147          15 :         nBands = atoi(CPLGetXMLValue(psRoot, "BandsCount", "4"));
    1148          15 :         const char *pszDataType = CPLGetXMLValue(psRoot, "DataType", "Byte");
    1149          15 :         eDataType = GDALGetDataTypeByName(pszDataType);
    1150          15 :         if ((eDataType == GDT_Unknown) || (eDataType >= GDT_TypeCount))
    1151             :         {
    1152           0 :             CPLError(CE_Failure, CPLE_AppDefined,
    1153             :                      "GDALWMTS: Invalid value in DataType. Data type \"%s\" is "
    1154             :                      "not supported.",
    1155             :                      pszDataType);
    1156           0 :             CPLDestroyXMLNode(psGDALWMTS);
    1157           0 :             return nullptr;
    1158             :         }
    1159             : 
    1160             :         const char *pszULX =
    1161          15 :             CPLGetXMLValue(psRoot, "DataWindow.UpperLeftX", nullptr);
    1162             :         const char *pszULY =
    1163          15 :             CPLGetXMLValue(psRoot, "DataWindow.UpperLeftY", nullptr);
    1164             :         const char *pszLRX =
    1165          15 :             CPLGetXMLValue(psRoot, "DataWindow.LowerRightX", nullptr);
    1166             :         const char *pszLRY =
    1167          15 :             CPLGetXMLValue(psRoot, "DataWindow.LowerRightY", nullptr);
    1168          15 :         if (pszULX && pszULY && pszLRX && pszLRY)
    1169             :         {
    1170          14 :             sAOI.MinX = CPLAtof(pszULX);
    1171          14 :             sAOI.MaxY = CPLAtof(pszULY);
    1172          14 :             sAOI.MaxX = CPLAtof(pszLRX);
    1173          14 :             sAOI.MinY = CPLAtof(pszLRY);
    1174          14 :             bHasAOI = TRUE;
    1175             :         }
    1176             : 
    1177          15 :         CPLDestroyXMLNode(psGDALWMTS);
    1178             : 
    1179          15 :         CPLDestroyXMLNode(psXML);
    1180          15 :         const CPLStringList aosHTTPOptions(BuildHTTPRequestOpts(osOtherXML));
    1181          15 :         psXML = GetCapabilitiesResponse(osGetCapabilitiesURL,
    1182             :                                         aosHTTPOptions.List());
    1183             :     }
    1184          53 :     else if (!STARTS_WITH_CI(poOpenInfo->pszFilename, "WMTS:") &&
    1185           5 :              !STARTS_WITH(poOpenInfo->pszFilename, "http://") &&
    1186           4 :              !STARTS_WITH(poOpenInfo->pszFilename, "https://"))
    1187             :     {
    1188           4 :         osGetCapabilitiesURL = poOpenInfo->pszFilename;
    1189           4 :         psXML = CPLParseXMLFile(poOpenInfo->pszFilename);
    1190             :     }
    1191          68 :     if (psXML == nullptr)
    1192           4 :         return nullptr;
    1193          64 :     CPLStripXMLNamespace(psXML, nullptr, TRUE);
    1194             : 
    1195          64 :     CPLXMLNode *psContents = CPLGetXMLNode(psXML, "=Capabilities.Contents");
    1196          64 :     if (psContents == nullptr)
    1197             :     {
    1198           1 :         CPLError(CE_Failure, CPLE_AppDefined,
    1199             :                  "Missing Capabilities.Contents element");
    1200           1 :         CPLDestroyXMLNode(psXML);
    1201           1 :         return nullptr;
    1202             :     }
    1203             : 
    1204          63 :     if (STARTS_WITH(osGetCapabilitiesURL, "/vsimem/"))
    1205             :     {
    1206          59 :         osGetCapabilitiesURL = GetOperationKVPURL(psXML, "GetCapabilities");
    1207          59 :         if (osGetCapabilitiesURL.empty())
    1208             :         {
    1209             :             // (ERO) I'm not even sure this is correct at all...
    1210          55 :             const char *pszHref = CPLGetXMLValue(
    1211             :                 psXML, "=Capabilities.ServiceMetadataURL.href", nullptr);
    1212          55 :             if (pszHref)
    1213          26 :                 osGetCapabilitiesURL = pszHref;
    1214             :         }
    1215             :         else
    1216             :         {
    1217             :             osGetCapabilitiesURL =
    1218           4 :                 CPLURLAddKVP(osGetCapabilitiesURL, "service", "WMTS");
    1219           8 :             osGetCapabilitiesURL = CPLURLAddKVP(osGetCapabilitiesURL, "request",
    1220           4 :                                                 "GetCapabilities");
    1221             :         }
    1222             :     }
    1223         126 :     CPLString osCapabilitiesFilename(osGetCapabilitiesURL);
    1224          63 :     if (!STARTS_WITH_CI(osCapabilitiesFilename, "WMTS:"))
    1225          63 :         osCapabilitiesFilename = "WMTS:" + osGetCapabilitiesURL;
    1226             : 
    1227          63 :     int nLayerCount = 0;
    1228         126 :     CPLStringList aosSubDatasets;
    1229         126 :     CPLString osSelectLayer(osLayer), osSelectTMS(osTMS),
    1230         126 :         osSelectStyle(osStyle);
    1231         126 :     CPLString osSelectLayerTitle, osSelectLayerAbstract;
    1232         126 :     CPLString osSelectTileFormat(osTileFormat),
    1233         126 :         osSelectInfoFormat(osInfoFormat);
    1234          63 :     int nCountTileFormat = 0;
    1235          63 :     int nCountInfoFormat = 0;
    1236         126 :     CPLString osURLTileTemplate;
    1237         126 :     CPLString osURLFeatureInfoTemplate;
    1238         126 :     std::set<CPLString> aoSetLayers;
    1239         126 :     std::map<CPLString, OGREnvelope> aoMapBoundingBox;
    1240         126 :     std::map<CPLString, WMTSTileMatrixLimits> aoMapTileMatrixLimits;
    1241         126 :     std::map<CPLString, CPLString> aoMapDimensions;
    1242          63 :     bool bHasWarnedAutoSwap = false;
    1243          63 :     bool bHasWarnedAutoSwapBoundingBox = false;
    1244             : 
    1245             :     // Collect TileMatrixSet identifiers
    1246         126 :     std::set<std::string> oSetTMSIdentifiers;
    1247         202 :     for (CPLXMLNode *psIter = psContents->psChild; psIter != nullptr;
    1248         139 :          psIter = psIter->psNext)
    1249             :     {
    1250         139 :         if (psIter->eType != CXT_Element ||
    1251         139 :             strcmp(psIter->pszValue, "TileMatrixSet") != 0)
    1252          62 :             continue;
    1253             :         const char *pszIdentifier =
    1254          77 :             CPLGetXMLValue(psIter, "Identifier", nullptr);
    1255          77 :         if (pszIdentifier)
    1256          77 :             oSetTMSIdentifiers.insert(pszIdentifier);
    1257             :     }
    1258             : 
    1259         202 :     for (CPLXMLNode *psIter = psContents->psChild; psIter != nullptr;
    1260         139 :          psIter = psIter->psNext)
    1261             :     {
    1262         139 :         if (psIter->eType != CXT_Element ||
    1263         139 :             strcmp(psIter->pszValue, "Layer") != 0)
    1264          78 :             continue;
    1265          62 :         const char *pszIdentifier = CPLGetXMLValue(psIter, "Identifier", "");
    1266          62 :         if (aoSetLayers.find(pszIdentifier) != aoSetLayers.end())
    1267             :         {
    1268           0 :             CPLError(CE_Warning, CPLE_AppDefined,
    1269             :                      "Several layers with identifier '%s'. Only first one kept",
    1270             :                      pszIdentifier);
    1271             :         }
    1272          62 :         aoSetLayers.insert(pszIdentifier);
    1273          62 :         if (!osLayer.empty() && strcmp(osLayer, pszIdentifier) != 0)
    1274           1 :             continue;
    1275          61 :         const char *pszTitle = CPLGetXMLValue(psIter, "Title", nullptr);
    1276          61 :         if (osSelectLayer.empty())
    1277             :         {
    1278          47 :             osSelectLayer = pszIdentifier;
    1279             :         }
    1280          61 :         if (strcmp(osSelectLayer, pszIdentifier) == 0)
    1281             :         {
    1282          61 :             if (pszTitle != nullptr)
    1283          36 :                 osSelectLayerTitle = pszTitle;
    1284             :             const char *pszAbstract =
    1285          61 :                 CPLGetXMLValue(psIter, "Abstract", nullptr);
    1286          61 :             if (pszAbstract != nullptr)
    1287          17 :                 osSelectLayerAbstract = pszAbstract;
    1288             :         }
    1289             : 
    1290         122 :         std::vector<CPLString> aosTMS;
    1291         122 :         std::vector<CPLString> aosStylesIdentifier;
    1292         122 :         std::vector<CPLString> aosStylesTitle;
    1293             : 
    1294          61 :         CPLXMLNode *psSubIter = psIter->psChild;
    1295         483 :         for (; psSubIter != nullptr; psSubIter = psSubIter->psNext)
    1296             :         {
    1297         422 :             if (psSubIter->eType != CXT_Element)
    1298           1 :                 continue;
    1299         842 :             if (strcmp(osSelectLayer, pszIdentifier) == 0 &&
    1300         421 :                 strcmp(psSubIter->pszValue, "Format") == 0)
    1301             :             {
    1302          15 :                 const char *pszValue = CPLGetXMLValue(psSubIter, "", "");
    1303          21 :                 if (!osTileFormat.empty() &&
    1304           6 :                     strcmp(osTileFormat, pszValue) != 0)
    1305           3 :                     continue;
    1306          12 :                 nCountTileFormat++;
    1307          12 :                 if (osSelectTileFormat.empty() || EQUAL(pszValue, "image/png"))
    1308             :                 {
    1309          12 :                     osSelectTileFormat = pszValue;
    1310             :                 }
    1311             :             }
    1312         812 :             else if (strcmp(osSelectLayer, pszIdentifier) == 0 &&
    1313         406 :                      strcmp(psSubIter->pszValue, "InfoFormat") == 0)
    1314             :             {
    1315           4 :                 const char *pszValue = CPLGetXMLValue(psSubIter, "", "");
    1316           4 :                 if (!osInfoFormat.empty() &&
    1317           0 :                     strcmp(osInfoFormat, pszValue) != 0)
    1318           0 :                     continue;
    1319           4 :                 nCountInfoFormat++;
    1320           4 :                 if (osSelectInfoFormat.empty() ||
    1321           0 :                     (EQUAL(pszValue, "application/vnd.ogc.gml") &&
    1322           0 :                      !EQUAL(osSelectInfoFormat,
    1323           4 :                             "application/vnd.ogc.gml/3.1.1")) ||
    1324           0 :                     EQUAL(pszValue, "application/vnd.ogc.gml/3.1.1"))
    1325             :                 {
    1326           4 :                     osSelectInfoFormat = pszValue;
    1327             :                 }
    1328             :             }
    1329         804 :             else if (strcmp(osSelectLayer, pszIdentifier) == 0 &&
    1330         402 :                      strcmp(psSubIter->pszValue, "Dimension") == 0)
    1331             :             {
    1332             :                 /* Cf http://wmts.geo.admin.ch/1.0.0/WMTSCapabilities.xml */
    1333             :                 const char *pszDimensionIdentifier =
    1334          21 :                     CPLGetXMLValue(psSubIter, "Identifier", nullptr);
    1335             :                 const char *pszDefault =
    1336          21 :                     CPLGetXMLValue(psSubIter, "Default", "");
    1337          21 :                 if (pszDimensionIdentifier != nullptr)
    1338          21 :                     aoMapDimensions[pszDimensionIdentifier] = pszDefault;
    1339             :             }
    1340         381 :             else if (strcmp(psSubIter->pszValue, "TileMatrixSetLink") == 0)
    1341             :             {
    1342             :                 const char *pszTMS =
    1343          77 :                     CPLGetXMLValue(psSubIter, "TileMatrixSet", "");
    1344          77 :                 if (oSetTMSIdentifiers.find(pszTMS) == oSetTMSIdentifiers.end())
    1345             :                 {
    1346           2 :                     CPLDebug("WMTS",
    1347             :                              "Layer %s has a TileMatrixSetLink to %s, "
    1348             :                              "but it is not defined as a TileMatrixSet",
    1349             :                              pszIdentifier, pszTMS);
    1350           2 :                     continue;
    1351             :                 }
    1352          75 :                 if (!osTMS.empty() && strcmp(osTMS, pszTMS) != 0)
    1353          13 :                     continue;
    1354         124 :                 if (strcmp(osSelectLayer, pszIdentifier) == 0 &&
    1355          62 :                     osSelectTMS.empty())
    1356             :                 {
    1357          40 :                     osSelectTMS = pszTMS;
    1358             :                 }
    1359         124 :                 if (strcmp(osSelectLayer, pszIdentifier) == 0 &&
    1360          62 :                     strcmp(osSelectTMS, pszTMS) == 0)
    1361             :                 {
    1362             :                     CPLXMLNode *psTMSLimits =
    1363          56 :                         CPLGetXMLNode(psSubIter, "TileMatrixSetLimits");
    1364          56 :                     if (psTMSLimits)
    1365           2 :                         ReadTMLimits(psTMSLimits, aoMapTileMatrixLimits);
    1366             :                 }
    1367          62 :                 aosTMS.push_back(pszTMS);
    1368             :             }
    1369         304 :             else if (strcmp(psSubIter->pszValue, "Style") == 0)
    1370             :             {
    1371          77 :                 int bIsDefault = CPLTestBool(
    1372          77 :                     CPLGetXMLValue(psSubIter, "isDefault", "false"));
    1373             :                 const char *l_pszIdentifier =
    1374          77 :                     CPLGetXMLValue(psSubIter, "Identifier", "");
    1375          77 :                 if (!osStyle.empty() && strcmp(osStyle, l_pszIdentifier) != 0)
    1376          15 :                     continue;
    1377             :                 const char *pszStyleTitle =
    1378          62 :                     CPLGetXMLValue(psSubIter, "Title", l_pszIdentifier);
    1379          62 :                 if (bIsDefault)
    1380             :                 {
    1381          32 :                     aosStylesIdentifier.insert(aosStylesIdentifier.begin(),
    1382          64 :                                                CPLString(l_pszIdentifier));
    1383          32 :                     aosStylesTitle.insert(aosStylesTitle.begin(),
    1384          64 :                                           CPLString(pszStyleTitle));
    1385          32 :                     if (strcmp(osSelectLayer, l_pszIdentifier) == 0 &&
    1386           0 :                         osSelectStyle.empty())
    1387             :                     {
    1388           0 :                         osSelectStyle = l_pszIdentifier;
    1389             :                     }
    1390             :                 }
    1391             :                 else
    1392             :                 {
    1393          30 :                     aosStylesIdentifier.push_back(l_pszIdentifier);
    1394          30 :                     aosStylesTitle.push_back(pszStyleTitle);
    1395             :                 }
    1396             :             }
    1397         454 :             else if (strcmp(osSelectLayer, pszIdentifier) == 0 &&
    1398         227 :                      (strcmp(psSubIter->pszValue, "BoundingBox") == 0 ||
    1399         219 :                       strcmp(psSubIter->pszValue, "WGS84BoundingBox") == 0))
    1400             :             {
    1401          70 :                 CPLString osCRS = CPLGetXMLValue(psSubIter, "crs", "");
    1402          35 :                 if (osCRS.empty())
    1403             :                 {
    1404          22 :                     if (strcmp(psSubIter->pszValue, "WGS84BoundingBox") == 0)
    1405             :                     {
    1406          22 :                         osCRS = "EPSG:4326";
    1407             :                     }
    1408             :                     else
    1409             :                     {
    1410           0 :                         int nCountTileMatrixSet = 0;
    1411           0 :                         CPLString osSingleTileMatrixSet;
    1412           0 :                         for (CPLXMLNode *psIter3 = psContents->psChild;
    1413           0 :                              psIter3 != nullptr; psIter3 = psIter3->psNext)
    1414             :                         {
    1415           0 :                             if (psIter3->eType != CXT_Element ||
    1416           0 :                                 strcmp(psIter3->pszValue, "TileMatrixSet") != 0)
    1417           0 :                                 continue;
    1418           0 :                             nCountTileMatrixSet++;
    1419           0 :                             if (nCountTileMatrixSet == 1)
    1420             :                                 osSingleTileMatrixSet =
    1421           0 :                                     CPLGetXMLValue(psIter3, "Identifier", "");
    1422             :                         }
    1423           0 :                         if (nCountTileMatrixSet == 1)
    1424             :                         {
    1425             :                             // For
    1426             :                             // 13-082_WMTS_Simple_Profile/schemas/wmts/1.0/profiles/WMTSSimple/examples/wmtsGetCapabilities_response_OSM.xml
    1427           0 :                             WMTSTileMatrixSet oTMS;
    1428           0 :                             if (ReadTMS(psContents, osSingleTileMatrixSet,
    1429           0 :                                         CPLString(), -1, oTMS,
    1430           0 :                                         bHasWarnedAutoSwap))
    1431             :                             {
    1432           0 :                                 osCRS = oTMS.osSRS;
    1433             :                             }
    1434             :                         }
    1435             :                     }
    1436             :                 }
    1437             :                 CPLString osLowerCorner =
    1438          70 :                     CPLGetXMLValue(psSubIter, "LowerCorner", "");
    1439             :                 CPLString osUpperCorner =
    1440          70 :                     CPLGetXMLValue(psSubIter, "UpperCorner", "");
    1441          70 :                 OGRSpatialReference oSRS;
    1442          35 :                 oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
    1443          35 :                 if (!osCRS.empty() && !osLowerCorner.empty() &&
    1444         105 :                     !osUpperCorner.empty() &&
    1445          70 :                     oSRS.SetFromUserInput(FixCRSName(osCRS)) == OGRERR_NONE)
    1446             :                 {
    1447             :                     const bool bSwap =
    1448          58 :                         !STARTS_WITH_CI(osCRS, "EPSG:") &&
    1449          23 :                         (CPL_TO_BOOL(oSRS.EPSGTreatsAsLatLong()) ||
    1450          10 :                          CPL_TO_BOOL(oSRS.EPSGTreatsAsNorthingEasting()));
    1451          35 :                     char **papszLC = CSLTokenizeString(osLowerCorner);
    1452          35 :                     char **papszUC = CSLTokenizeString(osUpperCorner);
    1453          35 :                     if (CSLCount(papszLC) == 2 && CSLCount(papszUC) == 2)
    1454             :                     {
    1455          35 :                         OGREnvelope sEnvelope;
    1456          35 :                         sEnvelope.MinX = CPLAtof(papszLC[(bSwap) ? 1 : 0]);
    1457          35 :                         sEnvelope.MinY = CPLAtof(papszLC[(bSwap) ? 0 : 1]);
    1458          35 :                         sEnvelope.MaxX = CPLAtof(papszUC[(bSwap) ? 1 : 0]);
    1459          35 :                         sEnvelope.MaxY = CPLAtof(papszUC[(bSwap) ? 0 : 1]);
    1460             : 
    1461          38 :                         if (bSwap && oSRS.IsGeographic() &&
    1462           3 :                             (std::fabs(sEnvelope.MinY) > 90 ||
    1463           3 :                              std::fabs(sEnvelope.MaxY) > 90))
    1464             :                         {
    1465           0 :                             if (!bHasWarnedAutoSwapBoundingBox)
    1466             :                             {
    1467           0 :                                 bHasWarnedAutoSwapBoundingBox = true;
    1468           0 :                                 CPLError(
    1469             :                                     CE_Warning, CPLE_AppDefined,
    1470             :                                     "Auto-correcting wrongly swapped "
    1471             :                                     "ows:%s coordinates. "
    1472             :                                     "They should be in latitude, longitude "
    1473             :                                     "order "
    1474             :                                     "but are presented in longitude, latitude "
    1475             :                                     "order. "
    1476             :                                     "This should be reported to the server "
    1477             :                                     "administrator.",
    1478             :                                     psSubIter->pszValue);
    1479             :                             }
    1480           0 :                             std::swap(sEnvelope.MinX, sEnvelope.MinY);
    1481           0 :                             std::swap(sEnvelope.MaxX, sEnvelope.MaxY);
    1482             :                         }
    1483             : 
    1484          35 :                         aoMapBoundingBox[osCRS] = sEnvelope;
    1485             :                     }
    1486          35 :                     CSLDestroy(papszLC);
    1487          35 :                     CSLDestroy(papszUC);
    1488             :                 }
    1489             :             }
    1490         384 :             else if (strcmp(osSelectLayer, pszIdentifier) == 0 &&
    1491         192 :                      strcmp(psSubIter->pszValue, "ResourceURL") == 0)
    1492             :             {
    1493          72 :                 if (EQUAL(CPLGetXMLValue(psSubIter, "resourceType", ""),
    1494             :                           "tile"))
    1495             :                 {
    1496             :                     const char *pszFormat =
    1497          55 :                         CPLGetXMLValue(psSubIter, "format", "");
    1498          55 :                     if (!osTileFormat.empty() &&
    1499           0 :                         strcmp(osTileFormat, pszFormat) != 0)
    1500           0 :                         continue;
    1501          55 :                     if (osURLTileTemplate.empty())
    1502             :                         osURLTileTemplate =
    1503          55 :                             CPLGetXMLValue(psSubIter, "template", "");
    1504             :                 }
    1505          17 :                 else if (EQUAL(CPLGetXMLValue(psSubIter, "resourceType", ""),
    1506             :                                "FeatureInfo"))
    1507             :                 {
    1508             :                     const char *pszFormat =
    1509          17 :                         CPLGetXMLValue(psSubIter, "format", "");
    1510          17 :                     if (!osInfoFormat.empty() &&
    1511           0 :                         strcmp(osInfoFormat, pszFormat) != 0)
    1512           0 :                         continue;
    1513          17 :                     if (osURLFeatureInfoTemplate.empty())
    1514             :                         osURLFeatureInfoTemplate =
    1515          17 :                             CPLGetXMLValue(psSubIter, "template", "");
    1516             :                 }
    1517             :             }
    1518             :         }
    1519         122 :         if (strcmp(osSelectLayer, pszIdentifier) == 0 &&
    1520         122 :             osSelectStyle.empty() && !aosStylesIdentifier.empty())
    1521             :         {
    1522          41 :             osSelectStyle = aosStylesIdentifier[0];
    1523             :         }
    1524         123 :         for (size_t i = 0; i < aosTMS.size(); i++)
    1525             :         {
    1526         131 :             for (size_t j = 0; j < aosStylesIdentifier.size(); j++)
    1527             :             {
    1528          69 :                 int nIdx = 1 + aosSubDatasets.size() / 2;
    1529         138 :                 CPLString osName(osCapabilitiesFilename);
    1530          69 :                 osName += ",layer=";
    1531          69 :                 osName += QuoteIfNecessary(pszIdentifier);
    1532          69 :                 if (aosTMS.size() > 1)
    1533             :                 {
    1534          20 :                     osName += ",tilematrixset=";
    1535          20 :                     osName += QuoteIfNecessary(aosTMS[i]);
    1536             :                 }
    1537          69 :                 if (aosStylesIdentifier.size() > 1)
    1538             :                 {
    1539          16 :                     osName += ",style=";
    1540          16 :                     osName += QuoteIfNecessary(aosStylesIdentifier[j]);
    1541             :                 }
    1542             :                 aosSubDatasets.AddNameValue(
    1543          69 :                     CPLSPrintf("SUBDATASET_%d_NAME", nIdx), osName);
    1544             : 
    1545         138 :                 CPLString osDesc("Layer ");
    1546          69 :                 osDesc += pszTitle ? pszTitle : pszIdentifier;
    1547          69 :                 if (aosTMS.size() > 1)
    1548             :                 {
    1549          20 :                     osDesc += ", tile matrix set ";
    1550          20 :                     osDesc += aosTMS[i];
    1551             :                 }
    1552          69 :                 if (aosStylesIdentifier.size() > 1)
    1553             :                 {
    1554          16 :                     osDesc += ", style ";
    1555          16 :                     osDesc += QuoteIfNecessary(aosStylesTitle[j]);
    1556             :                 }
    1557             :                 aosSubDatasets.AddNameValue(
    1558          69 :                     CPLSPrintf("SUBDATASET_%d_DESC", nIdx), osDesc);
    1559             :             }
    1560             :         }
    1561          61 :         if (!aosTMS.empty() && !aosStylesIdentifier.empty())
    1562          55 :             nLayerCount++;
    1563             :         else
    1564           6 :             CPLError(CE_Failure, CPLE_AppDefined,
    1565             :                      "Missing TileMatrixSetLink and/or Style");
    1566             :     }
    1567             : 
    1568          63 :     if (nLayerCount == 0)
    1569             :     {
    1570           8 :         CPLDestroyXMLNode(psXML);
    1571           8 :         return nullptr;
    1572             :     }
    1573             : 
    1574          55 :     WMTSDataset *poDS = new WMTSDataset();
    1575             : 
    1576          55 :     if (aosSubDatasets.size() > 2)
    1577           6 :         poDS->SetMetadata(aosSubDatasets.List(), "SUBDATASETS");
    1578             : 
    1579          55 :     if (nLayerCount == 1)
    1580             :     {
    1581          55 :         if (!osSelectLayerTitle.empty())
    1582          35 :             poDS->SetMetadataItem("TITLE", osSelectLayerTitle);
    1583          55 :         if (!osSelectLayerAbstract.empty())
    1584          16 :             poDS->SetMetadataItem("ABSTRACT", osSelectLayerAbstract);
    1585             : 
    1586          55 :         poDS->m_aosHTTPOptions = BuildHTTPRequestOpts(osOtherXML);
    1587          55 :         poDS->osLayer = osSelectLayer;
    1588          55 :         poDS->osTMS = osSelectTMS;
    1589             : 
    1590          55 :         WMTSTileMatrixSet oTMS;
    1591          55 :         if (!ReadTMS(psContents, osSelectTMS, osMaxTileMatrixIdentifier,
    1592             :                      nUserMaxZoomLevel, oTMS, bHasWarnedAutoSwap))
    1593             :         {
    1594           7 :             CPLDestroyXMLNode(psXML);
    1595           7 :             delete poDS;
    1596           7 :             return nullptr;
    1597             :         }
    1598             : 
    1599          96 :         const char *pszExtentMethod = CSLFetchNameValueDef(
    1600          48 :             poOpenInfo->papszOpenOptions, "EXTENT_METHOD", "AUTO");
    1601          48 :         ExtentMethod eExtentMethod = AUTO;
    1602          48 :         if (EQUAL(pszExtentMethod, "LAYER_BBOX"))
    1603           0 :             eExtentMethod = LAYER_BBOX;
    1604          48 :         else if (EQUAL(pszExtentMethod, "TILE_MATRIX_SET"))
    1605           0 :             eExtentMethod = TILE_MATRIX_SET;
    1606          48 :         else if (EQUAL(pszExtentMethod, "MOST_PRECISE_TILE_MATRIX"))
    1607           0 :             eExtentMethod = MOST_PRECISE_TILE_MATRIX;
    1608             : 
    1609          48 :         bool bAOIFromLayer = false;
    1610             : 
    1611             :         // Use in priority layer bounding box expressed in the SRS of the TMS
    1612          48 :         if ((!bHasAOI || bExtendBeyondDateLine) &&
    1613          96 :             (eExtentMethod == AUTO || eExtentMethod == LAYER_BBOX) &&
    1614          82 :             aoMapBoundingBox.find(oTMS.osSRS) != aoMapBoundingBox.end())
    1615             :         {
    1616           4 :             if (!bHasAOI)
    1617             :             {
    1618           4 :                 sAOI = aoMapBoundingBox[oTMS.osSRS];
    1619           4 :                 bAOIFromLayer = true;
    1620           4 :                 bHasAOI = TRUE;
    1621             :             }
    1622             : 
    1623           4 :             int bRecomputeAOI = FALSE;
    1624           4 :             if (bExtendBeyondDateLine)
    1625             :             {
    1626           1 :                 bExtendBeyondDateLine = FALSE;
    1627             : 
    1628           2 :                 OGRSpatialReference oWGS84;
    1629           1 :                 oWGS84.SetFromUserInput(SRS_WKT_WGS84_LAT_LONG);
    1630           1 :                 oWGS84.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
    1631             :                 OGRCoordinateTransformation *poCT =
    1632           1 :                     OGRCreateCoordinateTransformation(&oTMS.oSRS, &oWGS84);
    1633           1 :                 if (poCT != nullptr)
    1634             :                 {
    1635           1 :                     double dfX1 = sAOI.MinX;
    1636           1 :                     double dfY1 = sAOI.MinY;
    1637           1 :                     double dfX2 = sAOI.MaxX;
    1638           1 :                     double dfY2 = sAOI.MaxY;
    1639           2 :                     if (poCT->Transform(1, &dfX1, &dfY1) &&
    1640           1 :                         poCT->Transform(1, &dfX2, &dfY2))
    1641             :                     {
    1642           1 :                         if (fabs(dfX1 + 180) < 1e-8 && fabs(dfX2 - 180) < 1e-8)
    1643             :                         {
    1644           1 :                             bExtendBeyondDateLine = TRUE;
    1645           1 :                             bRecomputeAOI = TRUE;
    1646             :                         }
    1647           0 :                         else if (dfX2 < dfX1)
    1648             :                         {
    1649           0 :                             bExtendBeyondDateLine = TRUE;
    1650             :                         }
    1651             :                         else
    1652             :                         {
    1653           0 :                             CPLError(
    1654             :                                 CE_Warning, CPLE_AppDefined,
    1655             :                                 "ExtendBeyondDateLine disabled, since "
    1656             :                                 "longitudes of %s "
    1657             :                                 "BoundingBox do not span from -180 to 180 but "
    1658             :                                 "from %.16g to %.16g, "
    1659             :                                 "or longitude of upper right corner is not "
    1660             :                                 "lesser than the one of lower left corner",
    1661             :                                 oTMS.osSRS.c_str(), dfX1, dfX2);
    1662             :                         }
    1663             :                     }
    1664           1 :                     delete poCT;
    1665             :                 }
    1666             :             }
    1667           4 :             if (bExtendBeyondDateLine && bRecomputeAOI)
    1668             :             {
    1669           1 :                 bExtendBeyondDateLine = FALSE;
    1670             : 
    1671             :                 std::map<CPLString, OGREnvelope>::iterator oIter =
    1672           1 :                     aoMapBoundingBox.begin();
    1673           2 :                 for (; oIter != aoMapBoundingBox.end(); ++oIter)
    1674             :                 {
    1675           2 :                     OGRSpatialReference oSRS;
    1676           4 :                     if (oSRS.SetFromUserInput(
    1677           4 :                             FixCRSName(oIter->first),
    1678             :                             OGRSpatialReference::
    1679           2 :                                 SET_FROM_USER_INPUT_LIMITATIONS_get()) ==
    1680             :                         OGRERR_NONE)
    1681             :                     {
    1682           2 :                         OGRSpatialReference oWGS84;
    1683           2 :                         oWGS84.SetFromUserInput(SRS_WKT_WGS84_LAT_LONG);
    1684           2 :                         oWGS84.SetAxisMappingStrategy(
    1685             :                             OAMS_TRADITIONAL_GIS_ORDER);
    1686             :                         auto poCT =
    1687             :                             std::unique_ptr<OGRCoordinateTransformation>(
    1688             :                                 OGRCreateCoordinateTransformation(&oSRS,
    1689           2 :                                                                   &oWGS84));
    1690           2 :                         double dfX1 = oIter->second.MinX;
    1691           2 :                         double dfY1 = oIter->second.MinY;
    1692           2 :                         double dfX2 = oIter->second.MaxX;
    1693           2 :                         double dfY2 = oIter->second.MaxY;
    1694           2 :                         if (poCT != nullptr &&
    1695           2 :                             poCT->Transform(1, &dfX1, &dfY1) &&
    1696           4 :                             poCT->Transform(1, &dfX2, &dfY2) && dfX2 < dfX1)
    1697             :                         {
    1698           1 :                             dfX2 += 360;
    1699           2 :                             OGRSpatialReference oWGS84_with_over;
    1700           1 :                             oWGS84_with_over.SetFromUserInput(
    1701             :                                 "+proj=longlat +datum=WGS84 +over +wktext");
    1702           1 :                             char *pszProj4 = nullptr;
    1703           1 :                             oTMS.oSRS.exportToProj4(&pszProj4);
    1704           1 :                             oSRS.SetFromUserInput(
    1705             :                                 CPLSPrintf("%s +over +wktext", pszProj4));
    1706           1 :                             CPLFree(pszProj4);
    1707           1 :                             poCT.reset(OGRCreateCoordinateTransformation(
    1708             :                                 &oWGS84_with_over, &oSRS));
    1709           2 :                             if (poCT && poCT->Transform(1, &dfX1, &dfY1) &&
    1710           1 :                                 poCT->Transform(1, &dfX2, &dfY2))
    1711             :                             {
    1712           1 :                                 bExtendBeyondDateLine = TRUE;
    1713           1 :                                 sAOI.MinX = std::min(dfX1, dfX2);
    1714           1 :                                 sAOI.MinY = std::min(dfY1, dfY2);
    1715           1 :                                 sAOI.MaxX = std::max(dfX1, dfX2);
    1716           1 :                                 sAOI.MaxY = std::max(dfY1, dfY2);
    1717           1 :                                 CPLDebug("WMTS",
    1718             :                                          "ExtendBeyondDateLine using %s "
    1719             :                                          "bounding box",
    1720           1 :                                          oIter->first.c_str());
    1721             :                             }
    1722           1 :                             break;
    1723             :                         }
    1724             :                     }
    1725             :                 }
    1726             :             }
    1727             :         }
    1728             :         else
    1729             :         {
    1730          44 :             if (bExtendBeyondDateLine)
    1731             :             {
    1732           0 :                 CPLError(CE_Warning, CPLE_AppDefined,
    1733             :                          "ExtendBeyondDateLine disabled, since BoundingBox of "
    1734             :                          "%s is missing",
    1735             :                          oTMS.osSRS.c_str());
    1736           0 :                 bExtendBeyondDateLine = FALSE;
    1737             :             }
    1738             :         }
    1739             : 
    1740             :         // Otherwise default to reproject a layer bounding box expressed in
    1741             :         // another SRS
    1742          48 :         if (!bHasAOI && !aoMapBoundingBox.empty() &&
    1743           0 :             (eExtentMethod == AUTO || eExtentMethod == LAYER_BBOX))
    1744             :         {
    1745             :             std::map<CPLString, OGREnvelope>::iterator oIter =
    1746          13 :                 aoMapBoundingBox.begin();
    1747          13 :             for (; oIter != aoMapBoundingBox.end(); ++oIter)
    1748             :             {
    1749          13 :                 OGRSpatialReference oSRS;
    1750          13 :                 oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
    1751          26 :                 if (oSRS.SetFromUserInput(
    1752          26 :                         FixCRSName(oIter->first),
    1753             :                         OGRSpatialReference::
    1754          13 :                             SET_FROM_USER_INPUT_LIMITATIONS_get()) ==
    1755             :                     OGRERR_NONE)
    1756             :                 {
    1757             :                     // Check if this doesn't match the most precise tile matrix
    1758             :                     // by densifying its contour
    1759          13 :                     const WMTSTileMatrix &oTM = oTMS.aoTM.back();
    1760             : 
    1761          13 :                     bool bMatchFound = false;
    1762             :                     const char *pszProjectionTMS =
    1763          13 :                         oTMS.oSRS.GetAttrValue("PROJECTION");
    1764             :                     const char *pszProjectionBBOX =
    1765          13 :                         oSRS.GetAttrValue("PROJECTION");
    1766          13 :                     const bool bIsTMerc =
    1767          12 :                         (pszProjectionTMS != nullptr &&
    1768          13 :                          EQUAL(pszProjectionTMS, SRS_PT_TRANSVERSE_MERCATOR)) ||
    1769           0 :                         (pszProjectionBBOX != nullptr &&
    1770           0 :                          EQUAL(pszProjectionBBOX, SRS_PT_TRANSVERSE_MERCATOR));
    1771             :                     // If one of the 2 SRS is a TMerc, try with classical tmerc
    1772             :                     // or etmerc.
    1773          23 :                     for (int j = 0; j < (bIsTMerc ? 2 : 1); j++)
    1774             :                     {
    1775             :                         CPLString osOldVal = CPLGetThreadLocalConfigOption(
    1776          17 :                             "OSR_USE_APPROX_TMERC", "");
    1777          17 :                         if (bIsTMerc)
    1778             :                         {
    1779           8 :                             CPLSetThreadLocalConfigOption(
    1780             :                                 "OSR_USE_APPROX_TMERC",
    1781             :                                 (j == 0) ? "NO" : "YES");
    1782             :                         }
    1783             :                         OGRCoordinateTransformation *poRevCT =
    1784          17 :                             OGRCreateCoordinateTransformation(&oTMS.oSRS,
    1785             :                                                               &oSRS);
    1786          17 :                         if (bIsTMerc)
    1787             :                         {
    1788           8 :                             CPLSetThreadLocalConfigOption(
    1789             :                                 "OSR_USE_APPROX_TMERC",
    1790           8 :                                 osOldVal.empty() ? nullptr : osOldVal.c_str());
    1791             :                         }
    1792          17 :                         if (poRevCT != nullptr)
    1793             :                         {
    1794          17 :                             const auto sTMExtent = oTM.GetExtent();
    1795          17 :                             const double dfX0 = sTMExtent.MinX;
    1796          17 :                             const double dfY1 = sTMExtent.MaxY;
    1797          17 :                             const double dfX1 = sTMExtent.MaxX;
    1798          17 :                             const double dfY0 = sTMExtent.MinY;
    1799          17 :                             double dfXMin =
    1800             :                                 std::numeric_limits<double>::infinity();
    1801          17 :                             double dfYMin =
    1802             :                                 std::numeric_limits<double>::infinity();
    1803          17 :                             double dfXMax =
    1804             :                                 -std::numeric_limits<double>::infinity();
    1805          17 :                             double dfYMax =
    1806             :                                 -std::numeric_limits<double>::infinity();
    1807             : 
    1808          17 :                             const int NSTEPS = 20;
    1809         374 :                             for (int i = 0; i <= NSTEPS; i++)
    1810             :                             {
    1811         357 :                                 double dfX = dfX0 + (dfX1 - dfX0) * i / NSTEPS;
    1812         357 :                                 double dfY = dfY0;
    1813         357 :                                 if (poRevCT->Transform(1, &dfX, &dfY))
    1814             :                                 {
    1815         357 :                                     dfXMin = std::min(dfXMin, dfX);
    1816         357 :                                     dfYMin = std::min(dfYMin, dfY);
    1817         357 :                                     dfXMax = std::max(dfXMax, dfX);
    1818         357 :                                     dfYMax = std::max(dfYMax, dfY);
    1819             :                                 }
    1820             : 
    1821         357 :                                 dfX = dfX0 + (dfX1 - dfX0) * i / NSTEPS;
    1822         357 :                                 dfY = dfY1;
    1823         357 :                                 if (poRevCT->Transform(1, &dfX, &dfY))
    1824             :                                 {
    1825         357 :                                     dfXMin = std::min(dfXMin, dfX);
    1826         357 :                                     dfYMin = std::min(dfYMin, dfY);
    1827         357 :                                     dfXMax = std::max(dfXMax, dfX);
    1828         357 :                                     dfYMax = std::max(dfYMax, dfY);
    1829             :                                 }
    1830             : 
    1831         357 :                                 dfX = dfX0;
    1832         357 :                                 dfY = dfY0 + (dfY1 - dfY0) * i / NSTEPS;
    1833         357 :                                 if (poRevCT->Transform(1, &dfX, &dfY))
    1834             :                                 {
    1835         357 :                                     dfXMin = std::min(dfXMin, dfX);
    1836         357 :                                     dfYMin = std::min(dfYMin, dfY);
    1837         357 :                                     dfXMax = std::max(dfXMax, dfX);
    1838         357 :                                     dfYMax = std::max(dfYMax, dfY);
    1839             :                                 }
    1840             : 
    1841         357 :                                 dfX = dfX1;
    1842         357 :                                 dfY = dfY0 + (dfY1 - dfY0) * i / NSTEPS;
    1843         357 :                                 if (poRevCT->Transform(1, &dfX, &dfY))
    1844             :                                 {
    1845         357 :                                     dfXMin = std::min(dfXMin, dfX);
    1846         357 :                                     dfYMin = std::min(dfYMin, dfY);
    1847         357 :                                     dfXMax = std::max(dfXMax, dfX);
    1848         357 :                                     dfYMax = std::max(dfYMax, dfY);
    1849             :                                 }
    1850             :                             }
    1851             : 
    1852          17 :                             delete poRevCT;
    1853             : #ifdef DEBUG_VERBOSE
    1854             :                             CPLDebug(
    1855             :                                 "WMTS",
    1856             :                                 "Reprojected densified bbox of most "
    1857             :                                 "precise tile matrix in %s: %.8g %8g %8g %8g",
    1858             :                                 oIter->first.c_str(), dfXMin, dfYMin, dfXMax,
    1859             :                                 dfYMax);
    1860             : #endif
    1861          17 :                             if (fabs(oIter->second.MinX - dfXMin) <
    1862          68 :                                     1e-5 * std::max(fabs(oIter->second.MinX),
    1863          17 :                                                     fabs(dfXMin)) &&
    1864          10 :                                 fabs(oIter->second.MinY - dfYMin) <
    1865          10 :                                     1e-5 * std::max(fabs(oIter->second.MinY),
    1866          10 :                                                     fabs(dfYMin)) &&
    1867          10 :                                 fabs(oIter->second.MaxX - dfXMax) <
    1868          10 :                                     1e-5 * std::max(fabs(oIter->second.MaxX),
    1869          37 :                                                     fabs(dfXMax)) &&
    1870          10 :                                 fabs(oIter->second.MaxY - dfYMax) <
    1871          10 :                                     1e-5 * std::max(fabs(oIter->second.MaxY),
    1872          27 :                                                     fabs(dfYMax)))
    1873             :                             {
    1874           7 :                                 bMatchFound = true;
    1875             : #ifdef DEBUG_VERBOSE
    1876             :                                 CPLDebug("WMTS",
    1877             :                                          "Matches layer bounding box, so "
    1878             :                                          "that one is not significant");
    1879             : #endif
    1880           7 :                                 break;
    1881             :                             }
    1882             :                         }
    1883             :                     }
    1884             : 
    1885          13 :                     if (bMatchFound)
    1886             :                     {
    1887           7 :                         if (eExtentMethod == LAYER_BBOX)
    1888           0 :                             eExtentMethod = MOST_PRECISE_TILE_MATRIX;
    1889           7 :                         break;
    1890             :                     }
    1891             : 
    1892             :                     // Otherwise try to reproject the bounding box of the
    1893             :                     // layer from its SRS to the TMS SRS. Except in some cases
    1894             :                     // where this would result in non-sense. (this could be
    1895             :                     // improved !)
    1896           8 :                     if (!(bIsTMerc && oSRS.IsGeographic() &&
    1897           2 :                           fabs(oIter->second.MinX - -180) < 1e-8 &&
    1898           1 :                           fabs(oIter->second.MaxX - 180) < 1e-8))
    1899             :                     {
    1900             :                         OGRCoordinateTransformation *poCT =
    1901           5 :                             OGRCreateCoordinateTransformation(&oSRS,
    1902             :                                                               &oTMS.oSRS);
    1903           5 :                         if (poCT != nullptr)
    1904             :                         {
    1905           5 :                             double dfX1 = oIter->second.MinX;
    1906           5 :                             double dfY1 = oIter->second.MinY;
    1907           5 :                             double dfX2 = oIter->second.MaxX;
    1908           5 :                             double dfY2 = oIter->second.MinY;
    1909           5 :                             double dfX3 = oIter->second.MaxX;
    1910           5 :                             double dfY3 = oIter->second.MaxY;
    1911           5 :                             double dfX4 = oIter->second.MinX;
    1912           5 :                             double dfY4 = oIter->second.MaxY;
    1913           5 :                             if (poCT->Transform(1, &dfX1, &dfY1) &&
    1914           5 :                                 poCT->Transform(1, &dfX2, &dfY2) &&
    1915          15 :                                 poCT->Transform(1, &dfX3, &dfY3) &&
    1916           5 :                                 poCT->Transform(1, &dfX4, &dfY4))
    1917             :                             {
    1918           5 :                                 sAOI.MinX = std::min(std::min(dfX1, dfX2),
    1919           5 :                                                      std::min(dfX3, dfX4));
    1920           5 :                                 sAOI.MinY = std::min(std::min(dfY1, dfY2),
    1921           5 :                                                      std::min(dfY3, dfY4));
    1922           5 :                                 sAOI.MaxX = std::max(std::max(dfX1, dfX2),
    1923           5 :                                                      std::max(dfX3, dfX4));
    1924           5 :                                 sAOI.MaxY = std::max(std::max(dfY1, dfY2),
    1925           5 :                                                      std::max(dfY3, dfY4));
    1926           5 :                                 bHasAOI = TRUE;
    1927           5 :                                 bAOIFromLayer = true;
    1928             :                             }
    1929           5 :                             delete poCT;
    1930             :                         }
    1931             :                     }
    1932           6 :                     break;
    1933             :                 }
    1934             :             }
    1935             :         }
    1936             : 
    1937             :         // Clip the computed AOI with the union of the extent of the tile
    1938             :         // matrices
    1939          48 :         if (bHasAOI && !bExtendBeyondDateLine)
    1940             :         {
    1941          22 :             OGREnvelope sUnionTM;
    1942         104 :             for (const WMTSTileMatrix &oTM : oTMS.aoTM)
    1943             :             {
    1944          82 :                 if (!sUnionTM.IsInit())
    1945          22 :                     sUnionTM = oTM.GetExtent();
    1946             :                 else
    1947          60 :                     sUnionTM.Merge(oTM.GetExtent());
    1948             :             }
    1949          22 :             sAOI.Intersect(sUnionTM);
    1950             :         }
    1951             : 
    1952             :         // Otherwise default to BoundingBox of the TMS
    1953          48 :         if (!bHasAOI && oTMS.bBoundingBoxValid &&
    1954           0 :             (eExtentMethod == AUTO || eExtentMethod == TILE_MATRIX_SET))
    1955             :         {
    1956           1 :             CPLDebug("WMTS", "Using TMS bounding box as layer extent");
    1957           1 :             sAOI = oTMS.sBoundingBox;
    1958           1 :             bHasAOI = TRUE;
    1959             :         }
    1960             : 
    1961             :         // Otherwise default to implied BoundingBox of the most precise TM
    1962          48 :         if (!bHasAOI && (eExtentMethod == AUTO ||
    1963             :                          eExtentMethod == MOST_PRECISE_TILE_MATRIX))
    1964             :         {
    1965          24 :             const WMTSTileMatrix &oTM = oTMS.aoTM.back();
    1966          24 :             CPLDebug("WMTS", "Using TM level %s bounding box as layer extent",
    1967             :                      oTM.osIdentifier.c_str());
    1968             : 
    1969          24 :             sAOI = oTM.GetExtent();
    1970          24 :             bHasAOI = TRUE;
    1971             :         }
    1972             : 
    1973          48 :         if (!bHasAOI)
    1974             :         {
    1975           0 :             CPLError(CE_Failure, CPLE_AppDefined,
    1976             :                      "Could not determine raster extent");
    1977           0 :             CPLDestroyXMLNode(psXML);
    1978           0 :             delete poDS;
    1979           0 :             return nullptr;
    1980             :         }
    1981             : 
    1982          96 :         if (CPLTestBool(CSLFetchNameValueDef(
    1983          48 :                 poOpenInfo->papszOpenOptions,
    1984             :                 "CLIP_EXTENT_WITH_MOST_PRECISE_TILE_MATRIX",
    1985             :                 bAOIFromLayer ? "NO" : "YES")))
    1986             :         {
    1987             :             // Clip with implied BoundingBox of the most precise TM
    1988             :             // Useful for http://tileserver.maptiler.com/wmts
    1989          39 :             const WMTSTileMatrix &oTM = oTMS.aoTM.back();
    1990          39 :             const OGREnvelope sTMExtent = oTM.GetExtent();
    1991          39 :             OGREnvelope sAOINew(sAOI);
    1992             : 
    1993             :             // For
    1994             :             // https://data.linz.govt.nz/services;key=XXXXXXXX/wmts/1.0.0/set/69/WMTSCapabilities.xml
    1995             :             // only clip in Y since there's a warp over dateline.
    1996             :             // Update: it sems that the content of the server has changed since
    1997             :             // initial coding. So do X clipping in default mode.
    1998          39 :             if (!bExtendBeyondDateLine)
    1999             :             {
    2000          39 :                 sAOINew.MinX = std::max(sAOI.MinX, sTMExtent.MinX);
    2001          39 :                 sAOINew.MaxX = std::min(sAOI.MaxX, sTMExtent.MaxX);
    2002             :             }
    2003          39 :             sAOINew.MaxY = std::min(sAOI.MaxY, sTMExtent.MaxY);
    2004          39 :             sAOINew.MinY = std::max(sAOI.MinY, sTMExtent.MinY);
    2005          39 :             if (sAOI != sAOINew)
    2006             :             {
    2007           0 :                 CPLDebug(
    2008             :                     "WMTS",
    2009             :                     "Layer extent has been restricted from "
    2010             :                     "(%f,%f,%f,%f) to (%f,%f,%f,%f) using the "
    2011             :                     "implied bounding box of the most precise tile matrix. "
    2012             :                     "You may disable this by specifying the "
    2013             :                     "CLIP_EXTENT_WITH_MOST_PRECISE_TILE_MATRIX open option "
    2014             :                     "to NO.",
    2015             :                     sAOI.MinX, sAOI.MinY, sAOI.MaxX, sAOI.MaxY, sAOINew.MinX,
    2016             :                     sAOINew.MinY, sAOINew.MaxX, sAOINew.MaxY);
    2017             :             }
    2018          39 :             sAOI = sAOINew;
    2019             :         }
    2020             : 
    2021             :         // Clip with limits of most precise TM when available
    2022          96 :         if (CPLTestBool(CSLFetchNameValueDef(
    2023          48 :                 poOpenInfo->papszOpenOptions,
    2024             :                 "CLIP_EXTENT_WITH_MOST_PRECISE_TILE_MATRIX_LIMITS",
    2025             :                 bAOIFromLayer ? "NO" : "YES")))
    2026             :         {
    2027          40 :             const WMTSTileMatrix &oTM = oTMS.aoTM.back();
    2028          40 :             if (aoMapTileMatrixLimits.find(oTM.osIdentifier) !=
    2029          80 :                 aoMapTileMatrixLimits.end())
    2030             :             {
    2031           2 :                 OGREnvelope sAOINew(sAOI);
    2032             : 
    2033             :                 const WMTSTileMatrixLimits &oTMLimits =
    2034           2 :                     aoMapTileMatrixLimits[oTM.osIdentifier];
    2035           2 :                 const OGREnvelope sTMLimitsExtent = oTMLimits.GetExtent(oTM);
    2036           2 :                 sAOINew.Intersect(sTMLimitsExtent);
    2037             : 
    2038           2 :                 if (sAOI != sAOINew)
    2039             :                 {
    2040           2 :                     CPLDebug(
    2041             :                         "WMTS",
    2042             :                         "Layer extent has been restricted from "
    2043             :                         "(%f,%f,%f,%f) to (%f,%f,%f,%f) using the "
    2044             :                         "implied bounding box of the most precise tile matrix. "
    2045             :                         "You may disable this by specifying the "
    2046             :                         "CLIP_EXTENT_WITH_MOST_PRECISE_TILE_MATRIX_LIMITS open "
    2047             :                         "option "
    2048             :                         "to NO.",
    2049             :                         sAOI.MinX, sAOI.MinY, sAOI.MaxX, sAOI.MaxY,
    2050             :                         sAOINew.MinX, sAOINew.MinY, sAOINew.MaxX, sAOINew.MaxY);
    2051             :                 }
    2052           2 :                 sAOI = sAOINew;
    2053             :             }
    2054             :         }
    2055             : 
    2056          48 :         if (!osProjection.empty())
    2057             :         {
    2058           0 :             poDS->m_oSRS.SetFromUserInput(
    2059             :                 osProjection,
    2060             :                 OGRSpatialReference::SET_FROM_USER_INPUT_LIMITATIONS_get());
    2061             :         }
    2062          48 :         if (poDS->m_oSRS.IsEmpty())
    2063             :         {
    2064          48 :             poDS->m_oSRS = oTMS.oSRS;
    2065             :         }
    2066             : 
    2067          48 :         if (osURLTileTemplate.empty())
    2068             :         {
    2069           5 :             osURLTileTemplate = GetOperationKVPURL(psXML, "GetTile");
    2070           5 :             if (osURLTileTemplate.empty())
    2071             :             {
    2072           1 :                 CPLError(CE_Failure, CPLE_AppDefined,
    2073             :                          "No RESTful nor KVP GetTile operation found");
    2074           1 :                 CPLDestroyXMLNode(psXML);
    2075           1 :                 delete poDS;
    2076           1 :                 return nullptr;
    2077             :             }
    2078             :             osURLTileTemplate =
    2079           4 :                 CPLURLAddKVP(osURLTileTemplate, "service", "WMTS");
    2080             :             osURLTileTemplate =
    2081           4 :                 CPLURLAddKVP(osURLTileTemplate, "request", "GetTile");
    2082             :             osURLTileTemplate =
    2083           4 :                 CPLURLAddKVP(osURLTileTemplate, "version", "1.0.0");
    2084             :             osURLTileTemplate =
    2085           4 :                 CPLURLAddKVP(osURLTileTemplate, "layer", osSelectLayer);
    2086             :             osURLTileTemplate =
    2087           4 :                 CPLURLAddKVP(osURLTileTemplate, "style", osSelectStyle);
    2088             :             osURLTileTemplate =
    2089           4 :                 CPLURLAddKVP(osURLTileTemplate, "format", osSelectTileFormat);
    2090             :             osURLTileTemplate =
    2091           4 :                 CPLURLAddKVP(osURLTileTemplate, "TileMatrixSet", osSelectTMS);
    2092           4 :             osURLTileTemplate += "&TileMatrix={TileMatrix}";
    2093           4 :             osURLTileTemplate += "&TileRow=${y}";
    2094           4 :             osURLTileTemplate += "&TileCol=${x}";
    2095             : 
    2096             :             std::map<CPLString, CPLString>::iterator oIter =
    2097           4 :                 aoMapDimensions.begin();
    2098           8 :             for (; oIter != aoMapDimensions.end(); ++oIter)
    2099             :             {
    2100          12 :                 osURLTileTemplate = CPLURLAddKVP(osURLTileTemplate,
    2101          12 :                                                  oIter->first, oIter->second);
    2102             :             }
    2103             :             // CPLDebug("WMTS", "osURLTileTemplate = %s",
    2104             :             // osURLTileTemplate.c_str());
    2105             :         }
    2106             :         else
    2107             :         {
    2108             :             osURLTileTemplate =
    2109          43 :                 Replace(osURLTileTemplate, "{Style}", osSelectStyle);
    2110             :             osURLTileTemplate =
    2111          43 :                 Replace(osURLTileTemplate, "{TileMatrixSet}", osSelectTMS);
    2112          43 :             osURLTileTemplate = Replace(osURLTileTemplate, "{TileCol}", "${x}");
    2113          43 :             osURLTileTemplate = Replace(osURLTileTemplate, "{TileRow}", "${y}");
    2114             : 
    2115             :             std::map<CPLString, CPLString>::iterator oIter =
    2116          43 :                 aoMapDimensions.begin();
    2117          57 :             for (; oIter != aoMapDimensions.end(); ++oIter)
    2118             :             {
    2119          42 :                 osURLTileTemplate = Replace(
    2120          14 :                     osURLTileTemplate, CPLSPrintf("{%s}", oIter->first.c_str()),
    2121          28 :                     oIter->second);
    2122             :             }
    2123             :         }
    2124          47 :         osURLTileTemplate += osExtraQueryParameters;
    2125             : 
    2126          47 :         if (osURLFeatureInfoTemplate.empty() && !osSelectInfoFormat.empty())
    2127             :         {
    2128             :             osURLFeatureInfoTemplate =
    2129           4 :                 GetOperationKVPURL(psXML, "GetFeatureInfo");
    2130           4 :             if (!osURLFeatureInfoTemplate.empty())
    2131             :             {
    2132             :                 osURLFeatureInfoTemplate =
    2133           4 :                     CPLURLAddKVP(osURLFeatureInfoTemplate, "service", "WMTS");
    2134           8 :                 osURLFeatureInfoTemplate = CPLURLAddKVP(
    2135           4 :                     osURLFeatureInfoTemplate, "request", "GetFeatureInfo");
    2136             :                 osURLFeatureInfoTemplate =
    2137           4 :                     CPLURLAddKVP(osURLFeatureInfoTemplate, "version", "1.0.0");
    2138           8 :                 osURLFeatureInfoTemplate = CPLURLAddKVP(
    2139           4 :                     osURLFeatureInfoTemplate, "layer", osSelectLayer);
    2140           8 :                 osURLFeatureInfoTemplate = CPLURLAddKVP(
    2141           4 :                     osURLFeatureInfoTemplate, "style", osSelectStyle);
    2142             :                 // osURLFeatureInfoTemplate =
    2143             :                 // CPLURLAddKVP(osURLFeatureInfoTemplate, "format",
    2144             :                 // osSelectTileFormat);
    2145           8 :                 osURLFeatureInfoTemplate = CPLURLAddKVP(
    2146           4 :                     osURLFeatureInfoTemplate, "InfoFormat", osSelectInfoFormat);
    2147           4 :                 osURLFeatureInfoTemplate += "&TileMatrixSet={TileMatrixSet}";
    2148           4 :                 osURLFeatureInfoTemplate += "&TileMatrix={TileMatrix}";
    2149           4 :                 osURLFeatureInfoTemplate += "&TileRow={TileRow}";
    2150           4 :                 osURLFeatureInfoTemplate += "&TileCol={TileCol}";
    2151           4 :                 osURLFeatureInfoTemplate += "&J={J}";
    2152           4 :                 osURLFeatureInfoTemplate += "&I={I}";
    2153             : 
    2154             :                 std::map<CPLString, CPLString>::iterator oIter =
    2155           4 :                     aoMapDimensions.begin();
    2156           8 :                 for (; oIter != aoMapDimensions.end(); ++oIter)
    2157             :                 {
    2158          12 :                     osURLFeatureInfoTemplate = CPLURLAddKVP(
    2159          12 :                         osURLFeatureInfoTemplate, oIter->first, oIter->second);
    2160             :                 }
    2161             :                 // CPLDebug("WMTS", "osURLFeatureInfoTemplate = %s",
    2162             :                 // osURLFeatureInfoTemplate.c_str());
    2163             :             }
    2164             :         }
    2165             :         else
    2166             :         {
    2167             :             osURLFeatureInfoTemplate =
    2168          43 :                 Replace(osURLFeatureInfoTemplate, "{Style}", osSelectStyle);
    2169             : 
    2170             :             std::map<CPLString, CPLString>::iterator oIter =
    2171          43 :                 aoMapDimensions.begin();
    2172          57 :             for (; oIter != aoMapDimensions.end(); ++oIter)
    2173             :             {
    2174          42 :                 osURLFeatureInfoTemplate = Replace(
    2175             :                     osURLFeatureInfoTemplate,
    2176          42 :                     CPLSPrintf("{%s}", oIter->first.c_str()), oIter->second);
    2177             :             }
    2178             :         }
    2179          47 :         if (!osURLFeatureInfoTemplate.empty())
    2180          18 :             osURLFeatureInfoTemplate += osExtraQueryParameters;
    2181          47 :         poDS->osURLFeatureInfoTemplate = osURLFeatureInfoTemplate;
    2182          47 :         CPL_IGNORE_RET_VAL(osURLFeatureInfoTemplate);
    2183             : 
    2184             :         // Build all TMS datasets, wrapped in VRT datasets
    2185         160 :         for (int i = static_cast<int>(oTMS.aoTM.size() - 1); i >= 0; i--)
    2186             :         {
    2187         123 :             const WMTSTileMatrix &oTM = oTMS.aoTM[i];
    2188         123 :             double dfRasterXSize = (sAOI.MaxX - sAOI.MinX) / oTM.dfPixelSize;
    2189         123 :             double dfRasterYSize = (sAOI.MaxY - sAOI.MinY) / oTM.dfPixelSize;
    2190         123 :             if (dfRasterXSize > INT_MAX || dfRasterYSize > INT_MAX)
    2191             :             {
    2192          13 :                 continue;
    2193             :             }
    2194             : 
    2195         111 :             if (poDS->apoDatasets.empty())
    2196             :             {
    2197             :                 // Align AOI on pixel boundaries with respect to TopLeftCorner
    2198             :                 // of this tile matrix
    2199          96 :                 poDS->m_gt[0] =
    2200          48 :                     oTM.dfTLX +
    2201          48 :                     floor((sAOI.MinX - oTM.dfTLX) / oTM.dfPixelSize + 1e-10) *
    2202          48 :                         oTM.dfPixelSize;
    2203          48 :                 poDS->m_gt[1] = oTM.dfPixelSize;
    2204          48 :                 poDS->m_gt[2] = 0.0;
    2205          96 :                 poDS->m_gt[3] =
    2206          48 :                     oTM.dfTLY +
    2207          48 :                     ceil((sAOI.MaxY - oTM.dfTLY) / oTM.dfPixelSize - 1e-10) *
    2208          48 :                         oTM.dfPixelSize;
    2209          48 :                 poDS->m_gt[4] = 0.0;
    2210          48 :                 poDS->m_gt[5] = -oTM.dfPixelSize;
    2211          48 :                 poDS->nRasterXSize =
    2212          48 :                     int(0.5 + (sAOI.MaxX - poDS->m_gt[0]) / oTM.dfPixelSize);
    2213          48 :                 poDS->nRasterYSize =
    2214          48 :                     int(0.5 + (poDS->m_gt[3] - sAOI.MinY) / oTM.dfPixelSize);
    2215             :             }
    2216             : 
    2217             :             const int nRasterXSize =
    2218         111 :                 int(0.5 + poDS->nRasterXSize / oTM.dfPixelSize * poDS->m_gt[1]);
    2219             :             const int nRasterYSize =
    2220         111 :                 int(0.5 + poDS->nRasterYSize / oTM.dfPixelSize * poDS->m_gt[1]);
    2221         170 :             if (!poDS->apoDatasets.empty() &&
    2222          59 :                 (nRasterXSize < 128 || nRasterYSize < 128))
    2223             :             {
    2224          10 :                 break;
    2225             :             }
    2226             :             CPLString osURL(
    2227         101 :                 Replace(osURLTileTemplate, "{TileMatrix}", oTM.osIdentifier));
    2228             : 
    2229         101 :             const double dfTileWidthUnits = oTM.dfPixelSize * oTM.nTileWidth;
    2230         101 :             const double dfTileHeightUnits = oTM.dfPixelSize * oTM.nTileHeight;
    2231             : 
    2232             :             // Get bounds of this tile matrix / tile matrix limits
    2233         101 :             auto sTMExtent = oTM.GetExtent();
    2234         101 :             if (aoMapTileMatrixLimits.find(oTM.osIdentifier) !=
    2235         202 :                 aoMapTileMatrixLimits.end())
    2236             :             {
    2237             :                 const WMTSTileMatrixLimits &oTMLimits =
    2238           2 :                     aoMapTileMatrixLimits[oTM.osIdentifier];
    2239           2 :                 sTMExtent.Intersect(oTMLimits.GetExtent(oTM));
    2240             :             }
    2241             : 
    2242             :             // Compute the shift in terms of tiles between AOI and TM origin
    2243             :             const int nTileX =
    2244         101 :                 static_cast<int>(floor(std::max(sTMExtent.MinX, poDS->m_gt[0]) -
    2245         101 :                                        oTM.dfTLX + 1e-10) /
    2246         101 :                                  dfTileWidthUnits);
    2247             :             const int nTileY = static_cast<int>(
    2248         101 :                 floor(oTM.dfTLY - std::min(poDS->m_gt[3], sTMExtent.MaxY) +
    2249         101 :                       1e-10) /
    2250         101 :                 dfTileHeightUnits);
    2251             : 
    2252             :             // Compute extent of this zoom level slightly larger than the AOI
    2253             :             // and aligned on tile boundaries at this TM
    2254         101 :             double dfULX = oTM.dfTLX + nTileX * dfTileWidthUnits;
    2255         101 :             double dfULY = oTM.dfTLY - nTileY * dfTileHeightUnits;
    2256         101 :             double dfLRX = poDS->m_gt[0] + poDS->nRasterXSize * poDS->m_gt[1];
    2257         101 :             double dfLRY = poDS->m_gt[3] + poDS->nRasterYSize * poDS->m_gt[5];
    2258         101 :             dfLRX = dfULX + ceil((dfLRX - dfULX) / dfTileWidthUnits - 1e-10) *
    2259             :                                 dfTileWidthUnits;
    2260         101 :             dfLRY = dfULY + floor((dfLRY - dfULY) / dfTileHeightUnits + 1e-10) *
    2261             :                                 dfTileHeightUnits;
    2262             : 
    2263             :             // Clip TMS extent to the one of this TM
    2264         101 :             if (!bExtendBeyondDateLine)
    2265          99 :                 dfLRX = std::min(dfLRX, sTMExtent.MaxX);
    2266         101 :             dfLRY = std::max(dfLRY, sTMExtent.MinY);
    2267             : 
    2268         101 :             const double dfSizeX = 0.5 + (dfLRX - dfULX) / oTM.dfPixelSize;
    2269         101 :             const double dfSizeY = 0.5 + (dfULY - dfLRY) / oTM.dfPixelSize;
    2270         101 :             if (dfSizeX > INT_MAX || dfSizeY > INT_MAX)
    2271             :             {
    2272           1 :                 continue;
    2273             :             }
    2274         100 :             if (poDS->apoDatasets.empty())
    2275             :             {
    2276          47 :                 CPLDebug("WMTS", "Using tilematrix=%s (zoom level %d)",
    2277          47 :                          oTMS.aoTM[i].osIdentifier.c_str(), i);
    2278          47 :                 oTMS.aoTM.resize(1 + i);
    2279          47 :                 poDS->oTMS = oTMS;
    2280             :             }
    2281             : 
    2282         100 :             const int nSizeX = static_cast<int>(dfSizeX);
    2283         100 :             const int nSizeY = static_cast<int>(dfSizeY);
    2284             : 
    2285         100 :             const double dfDateLineX =
    2286         100 :                 oTM.dfTLX + oTM.nMatrixWidth * dfTileWidthUnits;
    2287         100 :             const int nSizeX1 =
    2288         100 :                 int(0.5 + (dfDateLineX - dfULX) / oTM.dfPixelSize);
    2289         100 :             const int nSizeX2 =
    2290         100 :                 int(0.5 + (dfLRX - dfDateLineX) / oTM.dfPixelSize);
    2291         100 :             if (bExtendBeyondDateLine && dfDateLineX > dfLRX)
    2292             :             {
    2293           0 :                 CPLDebug("WMTS", "ExtendBeyondDateLine ignored in that case");
    2294           0 :                 bExtendBeyondDateLine = FALSE;
    2295             :             }
    2296             : 
    2297             : #define WMS_TMS_TEMPLATE                                                       \
    2298             :     "<GDAL_WMS>"                                                               \
    2299             :     "<Service name=\"TMS\">"                                                   \
    2300             :     "    <ServerUrl>%s</ServerUrl>"                                            \
    2301             :     "</Service>"                                                               \
    2302             :     "<DataWindow>"                                                             \
    2303             :     "    <UpperLeftX>%.16g</UpperLeftX>"                                       \
    2304             :     "    <UpperLeftY>%.16g</UpperLeftY>"                                       \
    2305             :     "    <LowerRightX>%.16g</LowerRightX>"                                     \
    2306             :     "    <LowerRightY>%.16g</LowerRightY>"                                     \
    2307             :     "    <TileLevel>0</TileLevel>"                                             \
    2308             :     "    <TileX>%d</TileX>"                                                    \
    2309             :     "    <TileY>%d</TileY>"                                                    \
    2310             :     "    <SizeX>%d</SizeX>"                                                    \
    2311             :     "    <SizeY>%d</SizeY>"                                                    \
    2312             :     "    <YOrigin>top</YOrigin>"                                               \
    2313             :     "</DataWindow>"                                                            \
    2314             :     "<BlockSizeX>%d</BlockSizeX>"                                              \
    2315             :     "<BlockSizeY>%d</BlockSizeY>"                                              \
    2316             :     "<BandsCount>%d</BandsCount>"                                              \
    2317             :     "<DataType>%s</DataType>"                                                  \
    2318             :     "%s"                                                                       \
    2319             :     "</GDAL_WMS>"
    2320             : 
    2321             :             CPLString osStr(CPLSPrintf(
    2322         100 :                 WMS_TMS_TEMPLATE, WMTSEscapeXML(osURL).c_str(), dfULX, dfULY,
    2323             :                 (bExtendBeyondDateLine) ? dfDateLineX : dfLRX, dfLRY, nTileX,
    2324             :                 nTileY, (bExtendBeyondDateLine) ? nSizeX1 : nSizeX, nSizeY,
    2325         100 :                 oTM.nTileWidth, oTM.nTileHeight, nBands,
    2326         200 :                 GDALGetDataTypeName(eDataType), osOtherXML.c_str()));
    2327         100 :             const auto eLastErrorType = CPLGetLastErrorType();
    2328         100 :             const auto eLastErrorNum = CPLGetLastErrorNo();
    2329         100 :             const std::string osLastErrorMsg = CPLGetLastErrorMsg();
    2330         100 :             GDALDataset *poWMSDS = GDALDataset::Open(
    2331             :                 osStr, GDAL_OF_RASTER | GDAL_OF_SHARED | GDAL_OF_VERBOSE_ERROR,
    2332             :                 nullptr, nullptr, nullptr);
    2333         100 :             if (poWMSDS == nullptr)
    2334             :             {
    2335           0 :                 CPLDestroyXMLNode(psXML);
    2336           0 :                 delete poDS;
    2337           0 :                 return nullptr;
    2338             :             }
    2339             :             // Restore error state to what it was prior to WMS dataset opening
    2340             :             // if WMS dataset opening did not cause any new error to be emitted
    2341         100 :             if (CPLGetLastErrorType() == CE_None)
    2342         100 :                 CPLErrorSetState(eLastErrorType, eLastErrorNum,
    2343             :                                  osLastErrorMsg.c_str());
    2344             : 
    2345         100 :             VRTDatasetH hVRTDS = VRTCreate(nRasterXSize, nRasterYSize);
    2346         500 :             for (int iBand = 1; iBand <= nBands; iBand++)
    2347             :             {
    2348         400 :                 VRTAddBand(hVRTDS, eDataType, nullptr);
    2349             :             }
    2350             : 
    2351             :             int nSrcXOff, nSrcYOff, nDstXOff, nDstYOff;
    2352             : 
    2353         100 :             nSrcXOff = 0;
    2354         100 :             nDstXOff = static_cast<int>(
    2355         100 :                 std::round((dfULX - poDS->m_gt[0]) / oTM.dfPixelSize));
    2356             : 
    2357         100 :             nSrcYOff = 0;
    2358         100 :             nDstYOff = static_cast<int>(
    2359         100 :                 std::round((poDS->m_gt[3] - dfULY) / oTM.dfPixelSize));
    2360             : 
    2361         100 :             if (bExtendBeyondDateLine)
    2362             :             {
    2363             :                 int nSrcXOff2, nDstXOff2;
    2364             : 
    2365           2 :                 nSrcXOff2 = 0;
    2366           2 :                 nDstXOff2 = static_cast<int>(std::round(
    2367           2 :                     (dfDateLineX - poDS->m_gt[0]) / oTM.dfPixelSize));
    2368             : 
    2369             :                 osStr = CPLSPrintf(
    2370           2 :                     WMS_TMS_TEMPLATE, WMTSEscapeXML(osURL).c_str(),
    2371           2 :                     -dfDateLineX, dfULY, dfLRX - 2 * dfDateLineX, dfLRY, 0,
    2372           2 :                     nTileY, nSizeX2, nSizeY, oTM.nTileWidth, oTM.nTileHeight,
    2373           4 :                     nBands, GDALGetDataTypeName(eDataType), osOtherXML.c_str());
    2374             : 
    2375             :                 GDALDataset *poWMSDS2 =
    2376           2 :                     GDALDataset::Open(osStr, GDAL_OF_RASTER | GDAL_OF_SHARED,
    2377             :                                       nullptr, nullptr, nullptr);
    2378           2 :                 CPLAssert(poWMSDS2);
    2379             : 
    2380          10 :                 for (int iBand = 1; iBand <= nBands; iBand++)
    2381             :                 {
    2382             :                     VRTSourcedRasterBandH hVRTBand =
    2383           8 :                         reinterpret_cast<VRTSourcedRasterBandH>(
    2384             :                             GDALGetRasterBand(hVRTDS, iBand));
    2385           8 :                     VRTAddSimpleSource(
    2386             :                         hVRTBand, GDALGetRasterBand(poWMSDS, iBand), nSrcXOff,
    2387             :                         nSrcYOff, nSizeX1, nSizeY, nDstXOff, nDstYOff, nSizeX1,
    2388             :                         nSizeY, "NEAR", VRT_NODATA_UNSET);
    2389           8 :                     VRTAddSimpleSource(
    2390             :                         hVRTBand, GDALGetRasterBand(poWMSDS2, iBand), nSrcXOff2,
    2391             :                         nSrcYOff, nSizeX2, nSizeY, nDstXOff2, nDstYOff, nSizeX2,
    2392             :                         nSizeY, "NEAR", VRT_NODATA_UNSET);
    2393             :                 }
    2394             : 
    2395           2 :                 poWMSDS2->Dereference();
    2396             :             }
    2397             :             else
    2398             :             {
    2399         490 :                 for (int iBand = 1; iBand <= nBands; iBand++)
    2400             :                 {
    2401             :                     VRTSourcedRasterBandH hVRTBand =
    2402         392 :                         reinterpret_cast<VRTSourcedRasterBandH>(
    2403             :                             GDALGetRasterBand(hVRTDS, iBand));
    2404         392 :                     VRTAddSimpleSource(
    2405             :                         hVRTBand, GDALGetRasterBand(poWMSDS, iBand), nSrcXOff,
    2406             :                         nSrcYOff, nSizeX, nSizeY, nDstXOff, nDstYOff, nSizeX,
    2407             :                         nSizeY, "NEAR", VRT_NODATA_UNSET);
    2408             :                 }
    2409             :             }
    2410             : 
    2411         100 :             poWMSDS->Dereference();
    2412             : 
    2413         100 :             poDS->apoDatasets.push_back(GDALDataset::FromHandle(hVRTDS));
    2414             :         }
    2415             : 
    2416          47 :         if (poDS->apoDatasets.empty())
    2417             :         {
    2418           0 :             CPLError(CE_Failure, CPLE_AppDefined, "No zoom level found");
    2419           0 :             CPLDestroyXMLNode(psXML);
    2420           0 :             delete poDS;
    2421           0 :             return nullptr;
    2422             :         }
    2423             : 
    2424          47 :         poDS->SetMetadataItem("INTERLEAVE", "PIXEL", "IMAGE_STRUCTURE");
    2425         235 :         for (int i = 0; i < nBands; i++)
    2426         188 :             poDS->SetBand(i + 1, new WMTSBand(poDS, i + 1, eDataType));
    2427             : 
    2428          47 :         poDS->osXML = "<GDAL_WMTS>\n";
    2429          47 :         poDS->osXML += "  <GetCapabilitiesUrl>" +
    2430         141 :                        WMTSEscapeXML(osGetCapabilitiesURL) +
    2431          47 :                        "</GetCapabilitiesUrl>\n";
    2432          47 :         if (!osSelectLayer.empty())
    2433             :             poDS->osXML +=
    2434          33 :                 "  <Layer>" + WMTSEscapeXML(osSelectLayer) + "</Layer>\n";
    2435          47 :         if (!osSelectStyle.empty())
    2436             :             poDS->osXML +=
    2437          33 :                 "  <Style>" + WMTSEscapeXML(osSelectStyle) + "</Style>\n";
    2438          47 :         if (!osSelectTMS.empty())
    2439          94 :             poDS->osXML += "  <TileMatrixSet>" + WMTSEscapeXML(osSelectTMS) +
    2440          47 :                            "</TileMatrixSet>\n";
    2441          47 :         if (!osMaxTileMatrixIdentifier.empty())
    2442           3 :             poDS->osXML += "  <TileMatrix>" +
    2443           9 :                            WMTSEscapeXML(osMaxTileMatrixIdentifier) +
    2444           3 :                            "</TileMatrix>\n";
    2445          47 :         if (nUserMaxZoomLevel >= 0)
    2446           0 :             poDS->osXML += "  <ZoomLevel>" +
    2447           8 :                            CPLString().Printf("%d", nUserMaxZoomLevel) +
    2448           4 :                            "</ZoomLevel>\n";
    2449          47 :         if (nCountTileFormat > 1 && !osSelectTileFormat.empty())
    2450           2 :             poDS->osXML += "  <Format>" + WMTSEscapeXML(osSelectTileFormat) +
    2451           1 :                            "</Format>\n";
    2452          47 :         if (nCountInfoFormat > 1 && !osSelectInfoFormat.empty())
    2453           0 :             poDS->osXML += "  <InfoFormat>" +
    2454           0 :                            WMTSEscapeXML(osSelectInfoFormat) +
    2455           0 :                            "</InfoFormat>\n";
    2456          47 :         poDS->osXML += "  <DataWindow>\n";
    2457             :         poDS->osXML +=
    2458          47 :             CPLSPrintf("    <UpperLeftX>%.16g</UpperLeftX>\n", poDS->m_gt[0]);
    2459             :         poDS->osXML +=
    2460          47 :             CPLSPrintf("    <UpperLeftY>%.16g</UpperLeftY>\n", poDS->m_gt[3]);
    2461             :         poDS->osXML +=
    2462             :             CPLSPrintf("    <LowerRightX>%.16g</LowerRightX>\n",
    2463          47 :                        poDS->m_gt[0] + poDS->m_gt[1] * poDS->nRasterXSize);
    2464             :         poDS->osXML +=
    2465             :             CPLSPrintf("    <LowerRightY>%.16g</LowerRightY>\n",
    2466          47 :                        poDS->m_gt[3] + poDS->m_gt[5] * poDS->nRasterYSize);
    2467          47 :         poDS->osXML += "  </DataWindow>\n";
    2468          47 :         if (bExtendBeyondDateLine)
    2469             :             poDS->osXML +=
    2470           1 :                 "  <ExtendBeyondDateLine>true</ExtendBeyondDateLine>\n";
    2471          47 :         poDS->osXML += CPLSPrintf("  <BandsCount>%d</BandsCount>\n", nBands);
    2472             :         poDS->osXML += CPLSPrintf("  <DataType>%s</DataType>\n",
    2473          47 :                                   GDALGetDataTypeName(eDataType));
    2474          47 :         poDS->osXML += "  <Cache />\n";
    2475          47 :         poDS->osXML += "  <UnsafeSSL>true</UnsafeSSL>\n";
    2476          47 :         poDS->osXML += "  <ZeroBlockHttpCodes>204,404</ZeroBlockHttpCodes>\n";
    2477             :         poDS->osXML +=
    2478          47 :             "  <ZeroBlockOnServerException>true</ZeroBlockOnServerException>\n";
    2479          47 :         poDS->osXML += "</GDAL_WMTS>\n";
    2480             :     }
    2481             : 
    2482          47 :     CPLDestroyXMLNode(psXML);
    2483             : 
    2484          47 :     poDS->SetPamFlags(poDS->GetPamFlags() & ~GPF_DIRTY);
    2485          47 :     return poDS;
    2486             : }
    2487             : 
    2488             : /************************************************************************/
    2489             : /*                             CreateCopy()                             */
    2490             : /************************************************************************/
    2491             : 
    2492          21 : GDALDataset *WMTSDataset::CreateCopy(const char *pszFilename,
    2493             :                                      GDALDataset *poSrcDS,
    2494             :                                      CPL_UNUSED int bStrict,
    2495             :                                      CPL_UNUSED char **papszOptions,
    2496             :                                      CPL_UNUSED GDALProgressFunc pfnProgress,
    2497             :                                      CPL_UNUSED void *pProgressData)
    2498             : {
    2499          42 :     if (poSrcDS->GetDriver() == nullptr ||
    2500          21 :         poSrcDS->GetDriver() != GDALGetDriverByName("WMTS"))
    2501             :     {
    2502          19 :         CPLError(CE_Failure, CPLE_NotSupported,
    2503             :                  "Source dataset must be a WMTS dataset");
    2504          19 :         return nullptr;
    2505             :     }
    2506             : 
    2507           2 :     const char *pszXML = poSrcDS->GetMetadataItem("XML", "WMTS");
    2508           2 :     if (pszXML == nullptr)
    2509             :     {
    2510           0 :         CPLError(CE_Failure, CPLE_AppDefined,
    2511             :                  "Cannot get XML definition of source WMTS dataset");
    2512           0 :         return nullptr;
    2513             :     }
    2514             : 
    2515           2 :     VSILFILE *fp = VSIFOpenL(pszFilename, "wb");
    2516           2 :     if (fp == nullptr)
    2517           0 :         return nullptr;
    2518             : 
    2519           2 :     VSIFWriteL(pszXML, 1, strlen(pszXML), fp);
    2520           2 :     VSIFCloseL(fp);
    2521             : 
    2522           4 :     GDALOpenInfo oOpenInfo(pszFilename, GA_ReadOnly);
    2523           2 :     return Open(&oOpenInfo);
    2524             : }
    2525             : 
    2526             : /************************************************************************/
    2527             : /*                       GDALRegister_WMTS()                            */
    2528             : /************************************************************************/
    2529             : 
    2530        2024 : void GDALRegister_WMTS()
    2531             : 
    2532             : {
    2533        2024 :     if (!GDAL_CHECK_VERSION("WMTS driver"))
    2534           0 :         return;
    2535             : 
    2536        2024 :     if (GDALGetDriverByName(DRIVER_NAME) != nullptr)
    2537         283 :         return;
    2538             : 
    2539        1741 :     GDALDriver *poDriver = new GDALDriver();
    2540        1741 :     WMTSDriverSetCommonMetadata(poDriver);
    2541             : 
    2542        1741 :     poDriver->pfnOpen = WMTSDataset::Open;
    2543        1741 :     poDriver->pfnCreateCopy = WMTSDataset::CreateCopy;
    2544             : 
    2545        1741 :     GetGDALDriverManager()->RegisterDriver(poDriver);
    2546             : }

Generated by: LCOV version 1.14