Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

gml:CircleByCenterPoint: correctly take into account radius.uom for projected CRS #11599

Merged
merged 1 commit into from
Jan 8, 2025
+71 −29
Merged

gml:CircleByCenterPoint: correctly take into account radius.uom for projected CRS #11599

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
30 changes: 30 additions & 0 deletions autotest/ogr/ogr_gml_geom.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2295,6 +2295,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
70 changes: 41 additions & 29 deletions ogr/gml2ogrgeometry.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -811,26 +811,26 @@ static bool GML2OGRGeometry_AddToMultiSurface(
}

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

static double GetDistanceInMetre(double dfDistance, const char *pszUnits,
const char *pszAttribute, const char *pszId)
static double GetUOMInMetre(const char *pszUnits, const char *pszAttribute,
const char *pszId)
{
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);

if (pszId)
{
Expand Down Expand Up @@ -1211,8 +1211,12 @@ GML2OGRGeometry_XMLNode_Internal(const CPLXMLNode *psNode, const char *pszId,
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, "radius", pszId);
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 @@ -1231,9 +1235,7 @@ GML2OGRGeometry_XMLNode_Internal(const CPLXMLNode *psNode, const char *pszId,
}
}
}
if (bSRSUnitIsDegree && pszUnits != nullptr &&
(dfRadius = GetDistanceInMetre(dfRadius, pszUnits, "radius",
pszId)) > 0)
if (bSRSUnitIsDegree && dfUOMConv > 0)
{
bIsApproximateArc = true;
dfLastCurveApproximateArcRadius = dfRadius;
Expand Down Expand Up @@ -1808,8 +1810,11 @@ GML2OGRGeometry_XMLNode_Internal(const CPLXMLNode *psNode, const char *pszId,
ReportFailure("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, "radius", pszId);
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 @@ -1856,17 +1861,19 @@ GML2OGRGeometry_XMLNode_Internal(const CPLXMLNode *psNode, const char *pszId,
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, "radius", pszId)) > 0)
if (bSRSUnitIsDegree && dfUOMConv > 0)
{
auto poLS = std::make_unique<OGRLineString>();
// coverity[tainted_data]
Expand All @@ -1882,7 +1889,7 @@ GML2OGRGeometry_XMLNode_Internal(const CPLXMLNode *psNode, const char *pszId,
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 @@ -1892,7 +1899,7 @@ GML2OGRGeometry_XMLNode_Internal(const CPLXMLNode *psNode, const char *pszId,
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 @@ -1904,15 +1911,15 @@ GML2OGRGeometry_XMLNode_Internal(const CPLXMLNode *psNode, const char *pszId,
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 @@ -1955,8 +1962,11 @@ GML2OGRGeometry_XMLNode_Internal(const CPLXMLNode *psNode, const char *pszId,
ReportFailure("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, "radius", pszId);
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 @@ -1985,17 +1995,19 @@ GML2OGRGeometry_XMLNode_Internal(const CPLXMLNode *psNode, const char *pszId,
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, "radius", pszId)) > 0)
if (bSRSUnitIsDegree && dfUOMConv > 0)
{
auto poLS = std::make_unique<OGRLineString>();
const double dfStep =
Expand All @@ -2007,15 +2019,15 @@ GML2OGRGeometry_XMLNode_Internal(const CPLXMLNode *psNode, const char *pszId,
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
Loading