Coverage Report

Created: 2025-08-11 09:23

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