Skip to content

Commit

Permalink
gml:CircleByCenterPoint: correctly take into account radius.uom for p…
Browse files Browse the repository at this point in the history 
…rojected CRS

Fixes #11597
  • Loading branch information
@rouault
rouault committed Jan 7, 2025
1 parent 2116577 commit 93144c2
Show file tree
Hide file tree
Showing 2 changed files with 71 additions and 25 deletions.
30 changes: 30 additions & 0 deletions autotest/ogr/ogr_gml_geom.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2296,6 +2296,36 @@ def test_gml_ArcByCenterPoint():
geom, ogr.CreateGeometryFromWkt("CIRCULARSTRING (1 4,-1 2,1 0)")
)

gml = "<gml:ArcByCenterPoint><gml:pos>1 2</gml:pos><gml:radius uom='km'>0.002</gml:radius><gml:startAngle>90</gml:startAngle><gml:endAngle>270</gml:endAngle></gml:ArcByCenterPoint>"
geom = ogr.CreateGeometryFromGML(gml)

ogrtest.check_feature_geometry(
geom, ogr.CreateGeometryFromWkt("CIRCULARSTRING (1 4,-1 2,1 0)")
)

gml = "<gml:ArcByCenterPoint srsName='http://www.opengis.net/def/crs/EPSG/0/25832'><gml:pos>1 2</gml:pos><gml:radius uom='km'>0.002</gml:radius><gml:startAngle>90</gml:startAngle><gml:endAngle>270</gml:endAngle></gml:ArcByCenterPoint>"
geom = ogr.CreateGeometryFromGML(gml)

ogrtest.check_feature_geometry(
geom, ogr.CreateGeometryFromWkt("CIRCULARSTRING (1 4,-1 2,1 0)")
)

# EPSG:2222 : unit is ft
gml = "<gml:ArcByCenterPoint srsName='http://www.opengis.net/def/crs/EPSG/0/2222'><gml:pos>1 2</gml:pos><gml:radius uom='m'>0.6096</gml:radius><gml:startAngle>90</gml:startAngle><gml:endAngle>270</gml:endAngle></gml:ArcByCenterPoint>"
geom = ogr.CreateGeometryFromGML(gml)

ogrtest.check_feature_geometry(
geom, ogr.CreateGeometryFromWkt("CIRCULARSTRING (1 4,-1 2,1 0)")
)

# EPSG:2222 : unit is ft
gml = "<gml:ArcByCenterPoint srsName='http://www.opengis.net/def/crs/EPSG/0/2222'><gml:pos>1 2</gml:pos><gml:radius uom='ft'>2</gml:radius><gml:startAngle>90</gml:startAngle><gml:endAngle>270</gml:endAngle></gml:ArcByCenterPoint>"
geom = ogr.CreateGeometryFromGML(gml)

ogrtest.check_feature_geometry(
geom, ogr.CreateGeometryFromWkt("CIRCULARSTRING (1 4,-1 2,1 0)")
)


###############################################################################
# Test compound curve of ArcByCenterPoint whose ends don't exactly match
Expand Down
66 changes: 41 additions & 25 deletions ogr/gml2ogrgeometry.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -811,25 +811,25 @@ static bool GML2OGRGeometry_AddToMultiSurface(
}

/************************************************************************/
/* GetDistanceInMetre() */
/* GetUOMInMetre() */
/************************************************************************/

static double GetDistanceInMetre(double dfDistance, const char *pszUnits)
static double GetUOMInMetre(const char *pszUnits)
{
if (EQUAL(pszUnits, "m"))
return dfDistance;
if (!pszUnits || EQUAL(pszUnits, "m"))
return 1.0;

if (EQUAL(pszUnits, "km"))
return dfDistance * 1000;
return 1000.0;

if (EQUAL(pszUnits, "nm") || EQUAL(pszUnits, "[nmi_i]"))
return dfDistance * CPLAtof(SRS_UL_INTL_NAUT_MILE_CONV);
return CPLAtof(SRS_UL_INTL_NAUT_MILE_CONV);

if (EQUAL(pszUnits, "mi"))
return dfDistance * CPLAtof(SRS_UL_INTL_STAT_MILE_CONV);
return CPLAtof(SRS_UL_INTL_STAT_MILE_CONV);

if (EQUAL(pszUnits, "ft"))
return dfDistance * CPLAtof(SRS_UL_INTL_FOOT_CONV);
return CPLAtof(SRS_UL_INTL_FOOT_CONV);

CPLDebug("GML2OGRGeometry", "Unhandled unit: %s", pszUnits);
return -1;
Expand Down Expand Up @@ -1159,8 +1159,12 @@ static std::unique_ptr<OGRGeometry> GML2OGRGeometry_XMLNode_Internal(
const CPLXMLNode *psRadius = FindBareXMLChild(psChild, "radius");
if (psRadius && psRadius->eType == CXT_Element)
{
double dfRadius = CPLAtof(CPLGetXMLValue(psRadius, nullptr, "0"));
const char *pszUnits = CPLGetXMLValue(psRadius, "uom", nullptr);
const double dfUOMConv = GetUOMInMetre(pszUnits);
const double dfRadiusRaw =
CPLAtof(CPLGetXMLValue(psRadius, nullptr, "0"));
const double dfRadius =
dfUOMConv > 0 ? dfRadiusRaw * dfUOMConv : dfRadiusRaw;
bool bSRSUnitIsDegree = false;
bool bInvertedAxisOrder = false;
if (l_pszSRSName != nullptr)
Expand All @@ -1179,8 +1183,8 @@ static std::unique_ptr<OGRGeometry> GML2OGRGeometry_XMLNode_Internal(
}
}
}
if (bSRSUnitIsDegree && pszUnits != nullptr &&
(dfRadius = GetDistanceInMetre(dfRadius, pszUnits)) > 0)

if (bSRSUnitIsDegree && dfUOMConv > 0)
{
bIsApproximateArc = true;
dfLastCurveApproximateArcRadius = dfRadius;
Expand Down Expand Up @@ -1759,8 +1763,11 @@ static std::unique_ptr<OGRGeometry> GML2OGRGeometry_XMLNode_Internal(
CPLError(CE_Failure, CPLE_AppDefined, "Missing radius element.");
return nullptr;
}
const double dfRadius = CPLAtof(CPLGetXMLValue(psChild, nullptr, "0"));
const char *pszUnits = CPLGetXMLValue(psChild, "uom", nullptr);
const double dfUOMConv = GetUOMInMetre(pszUnits);
const double dfRadiusRaw =
CPLAtof(CPLGetXMLValue(psChild, nullptr, "0"));
double dfRadius = dfUOMConv > 0 ? dfRadiusRaw * dfUOMConv : dfRadiusRaw;

psChild = FindBareXMLChild(psNode, "startAngle");
if (psChild == nullptr || psChild->eType != CXT_Element)
Expand Down Expand Up @@ -1808,16 +1815,19 @@ static std::unique_ptr<OGRGeometry> GML2OGRGeometry_XMLNode_Internal(
dfSemiMajor = GetSemiMajor(&oSRS);
bInvertedAxisOrder =
CPL_TO_BOOL(oSRS.EPSGTreatsAsNorthingEasting());

const double dfSRSUnitsToMetre =
oSRS.GetLinearUnits(nullptr);
if (dfSRSUnitsToMetre > 0)
dfRadius /= dfSRSUnitsToMetre;
}
}
}

double dfCenterX = p.getX();
double dfCenterY = p.getY();

double dfDistance;
if (bSRSUnitIsDegree && pszUnits != nullptr &&
(dfDistance = GetDistanceInMetre(dfRadius, pszUnits)) > 0)
if (bSRSUnitIsDegree && dfUOMConv > 0)
{
auto poLS = std::make_unique<OGRLineString>();
// coverity[tainted_data]
Expand All @@ -1833,7 +1843,7 @@ static std::unique_ptr<OGRGeometry> GML2OGRGeometry_XMLNode_Internal(
if (bInvertedAxisOrder)
{
OGR_GreatCircle_ExtendPosition(
dfCenterX, dfCenterY, dfDistance,
dfCenterX, dfCenterY, dfRadius,
// See
// https://ext.eurocontrol.int/aixm_confluence/display/ACG/ArcByCenterPoint+Interpretation+Summary
dfAngle, dfSemiMajor, &dfLat, &dfLong);
Expand All @@ -1843,7 +1853,7 @@ static std::unique_ptr<OGRGeometry> GML2OGRGeometry_XMLNode_Internal(
else
{
OGR_GreatCircle_ExtendPosition(
dfCenterY, dfCenterX, dfDistance, 90 - dfAngle,
dfCenterY, dfCenterX, dfRadius, 90 - dfAngle,
dfSemiMajor, &dfLat, &dfLong);
p.setX(dfLong);
p.setY(dfLat);
Expand All @@ -1855,15 +1865,15 @@ static std::unique_ptr<OGRGeometry> GML2OGRGeometry_XMLNode_Internal(
double dfLat = 0.0;
if (bInvertedAxisOrder)
{
OGR_GreatCircle_ExtendPosition(dfCenterX, dfCenterY, dfDistance,
OGR_GreatCircle_ExtendPosition(dfCenterX, dfCenterY, dfRadius,
dfEndAngle, dfSemiMajor, &dfLat,
&dfLong);
p.setX(dfLat); // yes, external code will do the swap later
p.setY(dfLong);
}
else
{
OGR_GreatCircle_ExtendPosition(dfCenterY, dfCenterX, dfDistance,
OGR_GreatCircle_ExtendPosition(dfCenterY, dfCenterX, dfRadius,
90 - dfEndAngle, dfSemiMajor,
&dfLat, &dfLong);
p.setX(dfLong);
Expand Down Expand Up @@ -1906,8 +1916,11 @@ static std::unique_ptr<OGRGeometry> GML2OGRGeometry_XMLNode_Internal(
CPLError(CE_Failure, CPLE_AppDefined, "Missing radius element.");
return nullptr;
}
const double dfRadius = CPLAtof(CPLGetXMLValue(psChild, nullptr, "0"));
const char *pszUnits = CPLGetXMLValue(psChild, "uom", nullptr);
const double dfUOMConv = GetUOMInMetre(pszUnits);
const double dfRadiusRaw =
CPLAtof(CPLGetXMLValue(psChild, nullptr, "0"));
double dfRadius = dfUOMConv > 0 ? dfRadiusRaw * dfUOMConv : dfRadiusRaw;

OGRPoint p;
if (!ParseGMLCoordinates(psNode, &p, nSRSDimension))
Expand Down Expand Up @@ -1936,16 +1949,19 @@ static std::unique_ptr<OGRGeometry> GML2OGRGeometry_XMLNode_Internal(
dfSemiMajor = GetSemiMajor(&oSRS);
bInvertedAxisOrder =
CPL_TO_BOOL(oSRS.EPSGTreatsAsNorthingEasting());

const double dfSRSUnitsToMetre =
oSRS.GetLinearUnits(nullptr);
if (dfSRSUnitsToMetre > 0)
dfRadius /= dfSRSUnitsToMetre;
}
}
}

double dfCenterX = p.getX();
double dfCenterY = p.getY();

double dfDistance;
if (bSRSUnitIsDegree && pszUnits != nullptr &&
(dfDistance = GetDistanceInMetre(dfRadius, pszUnits)) > 0)
if (bSRSUnitIsDegree && dfUOMConv > 0)
{
auto poLS = std::make_unique<OGRLineString>();
const double dfStep =
Expand All @@ -1957,15 +1973,15 @@ static std::unique_ptr<OGRGeometry> GML2OGRGeometry_XMLNode_Internal(
if (bInvertedAxisOrder)
{
OGR_GreatCircle_ExtendPosition(
dfCenterX, dfCenterY, dfDistance, dfAngle, dfSemiMajor,
dfCenterX, dfCenterY, dfRadius, dfAngle, dfSemiMajor,
&dfLat, &dfLong);
p.setX(dfLat); // yes, external code will do the swap later
p.setY(dfLong);
}
else
{
OGR_GreatCircle_ExtendPosition(
dfCenterY, dfCenterX, dfDistance, dfAngle, dfSemiMajor,
dfCenterY, dfCenterX, dfRadius, dfAngle, dfSemiMajor,
&dfLat, &dfLong);
p.setX(dfLong);
p.setY(dfLat);
Expand Down

0 comments on commit 93144c2

Please sign in to comment.