LCOV - code coverage report
Current view: top level - frmts/wmts - wmtsdataset.cpp (source / functions) Hit Total Coverage
Test: gdal_filtered.info Lines: 1116 1242 89.9 %
Date: 2024-05-03 15:49:35 Functions: 28 29 96.6 %

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

Generated by: LCOV version 1.14