Skip to content

Commit

Permalink
WMTS: for geographic CRS with official lat,lon order, be robust to bo…
Browse files Browse the repository at this point in the history 
…unding box and TopLeftCorner being in the wrong axis order and emits a warning

Partly fixes #11387
  • Loading branch information
@rouault
rouault committed Nov 29, 2024
1 parent a4aa0b5 commit 2a89c3b
Show file tree
Hide file tree
Showing 3 changed files with 190 additions and 17 deletions.
102 changes: 102 additions & 0 deletions .../gdrivers/data/wmts/WMTSCapabilities_THEMIS_NightIR_ControlledMosaics_100m_v2_oct2018.xml
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,102 @@
<?xml version="1.0" encoding="UTF-8"?>
<Capabilities xmlns="http://www.opengis.net/wmts/1.0" xmlns:ows="http://www.opengis.net/ows/1.1" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xmlns:gml="http://www.opengis.net/gml" xsi:schemaLocation="http://www.opengis.net/wmts/1.0 http://schemas.opengis.net/wmts/1.0/wmtsGetCapabilities_response.xsd" version="1.0.0">
<!-- Service Identification -->
<ows:ServiceIdentification>
<ows:Title>THEMIS_NightIR_ControlledMosaics_100m_v2_oct2018</ows:Title>
<ows:ServiceType>OGC WMTS</ows:ServiceType>
<ows:ServiceTypeVersion>1.0.0</ows:ServiceTypeVersion>
</ows:ServiceIdentification> <!-- Operations Metadata -->
<ows:OperationsMetadata>
<ows:Operation name="GetCapabilities">
<ows:DCP>
<ows:HTTP>
<ows:Get xlink:href="https://trek.nasa.gov/tiles/Mars/EQ/THEMIS_NightIR_ControlledMosaics_100m_v2_oct2018/1.0.0/WMTSCapabilities.xml">
<ows:Constraint name="GetEncoding">
<ows:AllowedValues>
<ows:Value>RESTful</ows:Value>
</ows:AllowedValues>
</ows:Constraint>
</ows:Get>

</ows:HTTP>
</ows:DCP>
</ows:Operation>
<ows:Operation name="GetTile">
<ows:DCP>
<ows:HTTP>
<ows:Get xlink:href="https://trek.nasa.gov/tiles/Mars/EQ/THEMIS_NightIR_ControlledMosaics_100m_v2_oct2018/1.0.0/">
<ows:Constraint name="GetEncoding">
<ows:AllowedValues>
<ows:Value>RESTful</ows:Value>
</ows:AllowedValues>
</ows:Constraint>
</ows:Get>

</ows:HTTP>
</ows:DCP>
</ows:Operation>
</ows:OperationsMetadata>
<Contents>
<!--Layer-->
<Layer>
<ows:Title>THEMIS_NightIR_ControlledMosaics_100m_v2_oct2018</ows:Title>
<ows:Identifier>THEMIS_NightIR_ControlledMosaics_100m_v2_oct2018</ows:Identifier>
<ows:BoundingBox crs="urn:ogc:def:crs:EPSG::104905">
<ows:LowerCorner>-179.9999997 -65.0006576</ows:LowerCorner>
<ows:UpperCorner>179.9998849 65.0007535</ows:UpperCorner>
</ows:BoundingBox>
<ows:WGS84BoundingBox crs="urn:ogc:def:crs:OGC:2:84">
<ows:LowerCorner>-179.9999997 -65.0006576</ows:LowerCorner>
<ows:UpperCorner>179.9998849 65.0007535</ows:UpperCorner>
</ows:WGS84BoundingBox>

<Style isDefault="true">
<ows:Title>Default Style</ows:Title>
<ows:Identifier>default</ows:Identifier>
</Style>
<Format>image/png</Format>
<TileMatrixSetLink>
<TileMatrixSet>default028mm</TileMatrixSet>
</TileMatrixSetLink>
<ResourceURL format="image/png" resourceType="tile" template="https://trek.nasa.gov/tiles/Mars/EQ/THEMIS_NightIR_ControlledMosaics_100m_v2_oct2018/1.0.0//{Style}/{TileMatrixSet}/{TileMatrix}/{TileRow}/{TileCol}.png"/>
</Layer> <!--TileMatrixSet-->
<TileMatrixSet>
<ows:Title>default</ows:Title>
<ows:Abstract>The tile matrix set that has scale values calculated based on the dpi defined by OGC specification (dpi assumes 0.28mm as the physical distance of a pixel).</ows:Abstract>
<ows:Identifier>default028mm</ows:Identifier>
<ows:SupportedCRS>urn:ogc:def:crs:EPSG::104905</ows:SupportedCRS>
<TileMatrix><ows:Identifier>0</ows:Identifier><ScaleDenominator>2.7922763629807472E+08</ScaleDenominator><TopLeftCorner>-180.0 90.0</TopLeftCorner><TileWidth>256</TileWidth><TileHeight>256</TileHeight><MatrixWidth>2.0</MatrixWidth><MatrixHeight>1.0</MatrixHeight></TileMatrix>
<TileMatrix><ows:Identifier>1</ows:Identifier><ScaleDenominator>1.3961381814903736E+08</ScaleDenominator><TopLeftCorner>-180.0 90.0</TopLeftCorner><TileWidth>256</TileWidth><TileHeight>256</TileHeight><MatrixWidth>4.0</MatrixWidth><MatrixHeight>2.0</MatrixHeight></TileMatrix>
<TileMatrix><ows:Identifier>2</ows:Identifier><ScaleDenominator>6.9806909074518681E+07</ScaleDenominator><TopLeftCorner>-180.0 90.0</TopLeftCorner><TileWidth>256</TileWidth><TileHeight>256</TileHeight><MatrixWidth>8.0</MatrixWidth><MatrixHeight>4.0</MatrixHeight></TileMatrix>
<TileMatrix><ows:Identifier>3</ows:Identifier><ScaleDenominator>3.4903454537259340E+07</ScaleDenominator><TopLeftCorner>-180.0 90.0</TopLeftCorner><TileWidth>256</TileWidth><TileHeight>256</TileHeight><MatrixWidth>16.0</MatrixWidth><MatrixHeight>8.0</MatrixHeight></TileMatrix>
<TileMatrix><ows:Identifier>4</ows:Identifier><ScaleDenominator>1.7451727268629670E+07</ScaleDenominator><TopLeftCorner>-180.0 90.0</TopLeftCorner><TileWidth>256</TileWidth><TileHeight>256</TileHeight><MatrixWidth>32.0</MatrixWidth><MatrixHeight>16.0</MatrixHeight></TileMatrix>
<TileMatrix><ows:Identifier>5</ows:Identifier><ScaleDenominator>8.7258636343148351E+06</ScaleDenominator><TopLeftCorner>-180.0 90.0</TopLeftCorner><TileWidth>256</TileWidth><TileHeight>256</TileHeight><MatrixWidth>64.0</MatrixWidth><MatrixHeight>32.0</MatrixHeight></TileMatrix>
<TileMatrix><ows:Identifier>6</ows:Identifier><ScaleDenominator>4.3629318171574175E+06</ScaleDenominator><TopLeftCorner>-180.0 90.0</TopLeftCorner><TileWidth>256</TileWidth><TileHeight>256</TileHeight><MatrixWidth>128.0</MatrixWidth><MatrixHeight>64.0</MatrixHeight></TileMatrix>
<TileMatrix><ows:Identifier>7</ows:Identifier><ScaleDenominator>2.1814659085787088E+06</ScaleDenominator><TopLeftCorner>-180.0 90.0</TopLeftCorner><TileWidth>256</TileWidth><TileHeight>256</TileHeight><MatrixWidth>256.0</MatrixWidth><MatrixHeight>128.0</MatrixHeight></TileMatrix>
<TileMatrix><ows:Identifier>8</ows:Identifier><ScaleDenominator>1.0907329542893544E+06</ScaleDenominator><TopLeftCorner>-180.0 90.0</TopLeftCorner><TileWidth>256</TileWidth><TileHeight>256</TileHeight><MatrixWidth>512.0</MatrixWidth><MatrixHeight>256.0</MatrixHeight></TileMatrix>
<TileMatrix><ows:Identifier>9</ows:Identifier><ScaleDenominator>5.4536647714467719E+05</ScaleDenominator><TopLeftCorner>-180.0 90.0</TopLeftCorner><TileWidth>256</TileWidth><TileHeight>256</TileHeight><MatrixWidth>1024.0</MatrixWidth><MatrixHeight>512.0</MatrixHeight></TileMatrix>

<!--<TileMatrix>
<ows:Identifier>0</ows:Identifier>
<ScaleDenominator>2.3623511904761901E8</ScaleDenominator>
<TopLeftCorner>-180.0 90.0</TopLeftCorner>
<TileWidth>256</TileWidth>
<TileHeight>256</TileHeight>
<MatrixWidth>3</MatrixWidth>
<MatrixHeight>2</MatrixHeight>
</TileMatrix>
<TileMatrix>
<ows:Identifier>1</ows:Identifier>
<ScaleDenominator>1.1811755952380951E8</ScaleDenominator>
<TopLeftCorner>-180.0 90.0</TopLeftCorner>
<TileWidth>256</TileWidth>
<TileHeight>256</TileHeight>
<MatrixWidth>5</MatrixWidth>
<MatrixHeight>3</MatrixHeight>
</TileMatrix>
-->

</TileMatrixSet>
</Contents>
<ServiceMetadataURL xlink:href="https://trek.nasa.gov/tiles/Mars/EQ/THEMIS_NightIR_ControlledMosaics_100m_v2_oct2018/1.0.0/WMTSCapabilities.xml"/>
</Capabilities>
29 changes: 29 additions & 0 deletions autotest/gdrivers/wmts.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2001,3 +2001,32 @@ def test_wmts_force_opening_no_match():

drv = gdal.IdentifyDriverEx("data/byte.tif", allowed_drivers=["WMTS"])
assert drv is None


###############################################################################
# Test bug fix for https://github.com/OSGeo/gdal/issues/11387


@pytest.mark.require_proj(9)
@gdaltest.enable_exceptions()
def test_wmts_read_esri_code_disguised_as_epsg_and_wrong_axis_order():

with gdaltest.error_handler():
with gdal.Open(
"data/wmts/WMTSCapabilities_THEMIS_NightIR_ControlledMosaics_100m_v2_oct2018.xml"
) as ds:
assert gdal.GetLastErrorMsg().startswith(
"Auto-correcting wrongly swapped TileMatrix.TopLeftCorner coordinates"
)
assert ds.GetSpatialRef().GetAuthorityName(None) == "ESRI"
assert ds.GetSpatialRef().GetAuthorityCode(None) == "104905"
assert ds.GetGeoTransform() == pytest.approx(
(
-180.0,
0.0013717509172233527,
0.0,
65.00121128452162,
0.0,
-0.0013717509172233527,
)
)
76 changes: 59 additions & 17 deletions frmts/wmts/wmtsdataset.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -158,7 +158,8 @@ class WMTSDataset final : public GDALPamDataset
const char *pszOperation);
static int ReadTMS(CPLXMLNode *psContents, const CPLString &osIdentifier,
const CPLString &osMaxTileMatrixIdentifier,
int nMaxZoomLevel, WMTSTileMatrixSet &oTMS);
int nMaxZoomLevel, WMTSTileMatrixSet &oTMS,
bool &bHasWarnedAutoSwap);
static int ReadTMLimits(
CPLXMLNode *psTMSLimits,
std::map<CPLString, WMTSTileMatrixLimits> &aoMapTileMatrixLimits);
Expand Down Expand Up @@ -599,10 +600,9 @@ CPLString WMTSDataset::FixCRSName(const char *pszCRS)

int WMTSDataset::ReadTMS(CPLXMLNode *psContents, const CPLString &osIdentifier,
const CPLString &osMaxTileMatrixIdentifier,
int nMaxZoomLevel, WMTSTileMatrixSet &oTMS)
int nMaxZoomLevel, WMTSTileMatrixSet &oTMS,
bool &bHasWarnedAutoSwap)
{
bool bHasWarnedAutoSwap = false;

for (CPLXMLNode *psIter = psContents->psChild; psIter != nullptr;
psIter = psIter->psNext)
{
Expand All @@ -629,9 +629,10 @@ int WMTSDataset::ReadTMS(CPLXMLNode *psContents, const CPLString &osIdentifier,
pszSupportedCRS);
return FALSE;
}
int bSwap = !STARTS_WITH_CI(pszSupportedCRS, "EPSG:") &&
(oTMS.oSRS.EPSGTreatsAsLatLong() ||
oTMS.oSRS.EPSGTreatsAsNorthingEasting());
const bool bSwap =
!STARTS_WITH_CI(pszSupportedCRS, "EPSG:") &&
(CPL_TO_BOOL(oTMS.oSRS.EPSGTreatsAsLatLong()) ||
CPL_TO_BOOL(oTMS.oSRS.EPSGTreatsAsNorthingEasting()));
CPLXMLNode *psBB = CPLGetXMLNode(psIter, "BoundingBox");
oTMS.bBoundingBoxValid = false;
if (psBB != nullptr)
Expand Down Expand Up @@ -740,16 +741,21 @@ int WMTSDataset::ReadTMS(CPLXMLNode *psContents, const CPLString &osIdentifier,
}

// Hack for http://osm.geobretagne.fr/gwc01/service/wmts?request=getcapabilities
if (STARTS_WITH_CI(l_pszIdentifier, "EPSG:4326:") &&
oTM.dfTLY == -180.0)
// or https://trek.nasa.gov/tiles/Mars/EQ/THEMIS_NightIR_ControlledMosaics_100m_v2_oct2018/1.0.0/WMTSCapabilities.xml
if (oTM.dfTLY == -180.0 &&
(STARTS_WITH_CI(l_pszIdentifier, "EPSG:4326:") ||
(oTMS.oSRS.IsGeographic() && oTM.dfTLX == 90)))
{
if (!bHasWarnedAutoSwap)
{
bHasWarnedAutoSwap = true;
CPLError(CE_Warning, CPLE_AppDefined,
"Auto-correcting wrongly swapped "
"TileMatrix.TopLeftCorner coordinates. This "
"should be reported to the server administrator.");
"TileMatrix.TopLeftCorner coordinates. "
"They should be in latitude, longitude order "
"but are presented in longitude, latitude order. "
"This should be reported to the server "
"administrator.");
}
std::swap(oTM.dfTLX, oTM.dfTLY);
}
Expand Down Expand Up @@ -1246,6 +1252,8 @@ GDALDataset *WMTSDataset::Open(GDALOpenInfo *poOpenInfo)
std::map<CPLString, OGREnvelope> aoMapBoundingBox;
std::map<CPLString, WMTSTileMatrixLimits> aoMapTileMatrixLimits;
std::map<CPLString, CPLString> aoMapDimensions;
bool bHasWarnedAutoSwap = false;
bool bHasWarnedAutoSwapBoundingBox = false;

// Collect TileMatrixSet identifiers
std::set<std::string> oSetTMSIdentifiers;
Expand Down Expand Up @@ -1431,7 +1439,8 @@ GDALDataset *WMTSDataset::Open(GDALOpenInfo *poOpenInfo)
// 13-082_WMTS_Simple_Profile/schemas/wmts/1.0/profiles/WMTSSimple/examples/wmtsGetCapabilities_response_OSM.xml
WMTSTileMatrixSet oTMS;
if (ReadTMS(psContents, osSingleTileMatrixSet,
CPLString(), -1, oTMS))
CPLString(), -1, oTMS,
bHasWarnedAutoSwap))
{
osCRS = oTMS.osSRS;
}
Expand All @@ -1448,9 +1457,10 @@ GDALDataset *WMTSDataset::Open(GDALOpenInfo *poOpenInfo)
!osUpperCorner.empty() &&
oSRS.SetFromUserInput(FixCRSName(osCRS)) == OGRERR_NONE)
{
int bSwap = !STARTS_WITH_CI(osCRS, "EPSG:") &&
(oSRS.EPSGTreatsAsLatLong() ||
oSRS.EPSGTreatsAsNorthingEasting());
const bool bSwap =
!STARTS_WITH_CI(osCRS, "EPSG:") &&
(CPL_TO_BOOL(oSRS.EPSGTreatsAsLatLong()) ||
CPL_TO_BOOL(oSRS.EPSGTreatsAsNorthingEasting()));
char **papszLC = CSLTokenizeString(osLowerCorner);
char **papszUC = CSLTokenizeString(osUpperCorner);
if (CSLCount(papszLC) == 2 && CSLCount(papszUC) == 2)
Expand All @@ -1460,6 +1470,30 @@ GDALDataset *WMTSDataset::Open(GDALOpenInfo *poOpenInfo)
sEnvelope.MinY = CPLAtof(papszLC[(bSwap) ? 0 : 1]);
sEnvelope.MaxX = CPLAtof(papszUC[(bSwap) ? 1 : 0]);
sEnvelope.MaxY = CPLAtof(papszUC[(bSwap) ? 0 : 1]);

if (bSwap && oSRS.IsGeographic() &&
(std::fabs(sEnvelope.MinY) > 90 ||
std::fabs(sEnvelope.MaxY) > 90))
{
if (!bHasWarnedAutoSwapBoundingBox)
{
bHasWarnedAutoSwapBoundingBox = true;
CPLError(
CE_Warning, CPLE_AppDefined,
"Auto-correcting wrongly swapped "
"ows:%s coordinates. "
"They should be in latitude, longitude "
"order "
"but are presented in longitude, latitude "
"order. "
"This should be reported to the server "
"administrator.",
psSubIter->pszValue);
}
std::swap(sEnvelope.MinX, sEnvelope.MinY);
std::swap(sEnvelope.MaxX, sEnvelope.MaxY);
}

aoMapBoundingBox[osCRS] = sEnvelope;
}
CSLDestroy(papszLC);
Expand Down Expand Up @@ -1568,7 +1602,7 @@ GDALDataset *WMTSDataset::Open(GDALOpenInfo *poOpenInfo)

WMTSTileMatrixSet oTMS;
if (!ReadTMS(psContents, osSelectTMS, osMaxTileMatrixIdentifier,
nUserMaxZoomLevel, oTMS))
nUserMaxZoomLevel, oTMS, bHasWarnedAutoSwap))
{
CPLDestroyXMLNode(psXML);
delete poDS;
Expand Down Expand Up @@ -2303,7 +2337,10 @@ GDALDataset *WMTSDataset::Open(GDALOpenInfo *poOpenInfo)
nTileY, (bExtendBeyondDateLine) ? nSizeX1 : nSizeX, nSizeY,
oTM.nTileWidth, oTM.nTileHeight, nBands,
GDALGetDataTypeName(eDataType), osOtherXML.c_str()));
GDALDataset *poWMSDS = (GDALDataset *)GDALOpenEx(
const auto eLastErrorType = CPLGetLastErrorType();
const auto eLastErrorNum = CPLGetLastErrorNo();
const std::string osLastErrorMsg = CPLGetLastErrorMsg();
GDALDataset *poWMSDS = GDALDataset::Open(
osStr, GDAL_OF_RASTER | GDAL_OF_SHARED | GDAL_OF_VERBOSE_ERROR,
nullptr, nullptr, nullptr);
if (poWMSDS == nullptr)
Expand All @@ -2312,6 +2349,11 @@ GDALDataset *WMTSDataset::Open(GDALOpenInfo *poOpenInfo)
delete poDS;
return nullptr;
}
// Restore error state to what it was prior to WMS dataset opening
// if WMS dataset opening did not cause any new error to be emitted
if (CPLGetLastErrorType() == CE_None)
CPLErrorSetState(eLastErrorType, eLastErrorNum,
osLastErrorMsg.c_str());

VRTDatasetH hVRTDS = VRTCreate(nRasterXSize, nRasterYSize);
for (int iBand = 1; iBand <= nBands; iBand++)
Expand Down

0 comments on commit 2a89c3b

Please sign in to comment.