8000 Fixes related to WKT import/export of DerivedGeodetic/GeographicCRS by rouault · Pull Request #4536 · OSGeo/PROJ · GitHub
[go: up one dir, main page]
More Web Proxy on the site http://driver.im/
Skip to content

Fixes related to WKT import/export of DerivedGeodetic/GeographicCRS #4536

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

Merged
merged 2 commits into from
Jul 4, 2025
Merged

Fixes related to WKT import/export of DerivedGeodetic/GeographicCRS #4536

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
6aff559
WKT import: add heuristics to figure out the baseCRS of a DerivedGeod…
rouault Jul 2, 2025
7ed4f9d
WKT2:2019 export: export ID of BaseCRS of DerivedGeographic/GeodeticC…
rouault Jul 2, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
8000 Loading
Diff view
Diff view
8 changes: 8 additions & 0 deletions src/iso19111/crs.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -6384,6 +6384,10 @@ void DerivedGeodeticCRS::_exportToWKT(io::WKTFormatter *formatter) const {
l_datumEnsemble->_exportToWKT(formatter);
}
l_baseCRS->primeMeridian()->_exportToWKT(formatter);
if (formatter->use2019Keywords() &&
!(formatter->idOnTopLevelOnly() && formatter->topLevelHasId())) {
l_baseCRS->formatID(formatter);
}
formatter->endNode();

formatter->setUseDerivingConversion(true);
Expand Down Expand Up @@ -6522,6 +6526,10 @@ void DerivedGeographicCRS::_exportToWKT(io::WKTFormatter *formatter) const {
formatter->addQuotedString(l_baseCRS->nameStr());
l_baseCRS->exportDatumOrDatumEnsembleToWkt(formatter);
l_baseCRS->primeMeridian()->_exportToWKT(formatter);
if (formatter->use2019Keywords() &&
!(formatter->idOnTopLevelOnly() && formatter->topLevelHasId())) {
l_baseCRS->formatID(formatter);
}
formatter->endNode();

formatter->setUseDerivingConversion(true);
Expand Down
40 changes: 31 additions & 9 deletions src/iso19111/io.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1345,7 +1345,8 @@ struct WKTParser::Private {
const WKTNodeNNPtr &parentNode,
const UnitOfMeasure &defaultAngularUnit);

GeodeticCRSNNPtr buildGeodeticCRS(const WKTNodeNNPtr &node);
GeodeticCRSNNPtr buildGeodeticCRS(const WKTNodeNNPtr &node,
bool forceGeocentricIfNoCs = false);

CRSNNPtr buildDerivedGeodeticCRS(const WKTNodeNNPtr &node);

Expand Down Expand Up @@ -3155,7 +3156,8 @@ void WKTParser::Private::addExtensionProj4ToProp(const WKTNode::Private *nodeP,
// ---------------------------------------------------------------------------

GeodeticCRSNNPtr
WKTParser::Private::buildGeodeticCRS(const WKTNodeNNPtr &node) {
WKTParser::Private::buildGeodeticCRS(const WKTNodeNNPtr &node,
bool forceGeocentricIfNoCs) {
const auto *nodeP = node->GP();
auto &datumNode = nodeP->lookForChild(
WKTConstants::DATUM, WKTConstants::GEODETICDATUM, WKTConstants::TRF);
Expand Down Expand Up @@ -3262,6 +3264,10 @@ WKTParser::Private::buildGeodeticCRS(const WKTNodeNNPtr &node) {
}
}
}
if (forceGeocentricIfNoCs && isNull(csNode) &&
ci_equal(nodeName, WKTConstants::BASEGEODCRS)) {
cs = cs::CartesianCS::createGeocentric(UnitOfMeasure::METRE);
}

auto ellipsoidalCS = nn_dynamic_pointer_cast<EllipsoidalCS>(cs);
if (ellipsoidalCS) {
Expand Down Expand Up @@ -3378,8 +3384,6 @@ CRSNNPtr WKTParser::Private::buildDerivedGeodeticCRS(const WKTNodeNNPtr &node) {
// given the constraints enforced on calling code path
assert(!isNull(baseGeodCRSNode));

auto baseGeodCRS = buildGeodeticCRS(baseGeodCRSNode);

auto &derivingConversionNode =
nodeP->lookForChild(WKTConstants::DERIVINGCONVERSION);
if (isNull(derivingConversionNode)) {
Expand All @@ -3394,6 +3398,29 @@ CRSNNPtr WKTParser::Private::buildDerivedGeodeticCRS(const WKTNodeNNPtr &node) {
}
auto cs = buildCS(csNode, node, UnitOfMeasure::NONE);

bool forceGeocentricIfNoCs = false;
auto cartesianCS = nn_dynamic_pointer_cast<CartesianCS>(cs);
if (cartesianCS) {
if (cartesianCS->axisList().size() != 3) {
throw ParsingException(
"Cartesian CS for a GeodeticCRS should have 3 axis");
}
const int methodCode = derivingConversion->method()->getEPSGCode();
if ((methodCode == EPSG_CODE_METHOD_COORDINATE_FRAME_GEOCENTRIC ||
methodCode ==
EPSG_CODE_METHOD_COORDINATE_FRAME_FULL_MATRIX_GEOCENTRIC ||
methodCode == EPSG_CODE_METHOD_POSITION_VECTOR_GEOCENTRIC ||
methodCode == EPSG_CODE_METHOD_GEOCENTRIC_TRANSLATION_GEOCENTRIC ||
methodCode ==
EPSG_CODE_METHOD_TIME_DEPENDENT_POSITION_VECTOR_GEOCENTRIC ||
methodCode ==
EPSG_CODE_METHOD_TIME_DEPENDENT_COORDINATE_FRAME_GEOCENTRIC) &&
nodeP->lookForChild(WKTConstants::BASEGEODCRS) != nullptr) {
forceGeocentricIfNoCs = true;
}
}
auto baseGeodCRS = buildGeodeticCRS(baseGeodCRSNode, forceGeocentricIfNoCs);

auto ellipsoidalCS = nn_dynamic_pointer_cast<EllipsoidalCS>(cs);
if (ellipsoidalCS) {

Expand All @@ -3413,12 +3440,7 @@ CRSNNPtr WKTParser::Private::buildDerivedGeodeticCRS(const WKTNodeNNPtr &node) {
cs->getWKT2Type(true)));
}

auto cartesianCS = nn_dynamic_pointer_cast<CartesianCS>(cs);
if (cartesianCS) {
if (cartesianCS->axisList().size() != 3) {
throw ParsingException(
"Cartesian CS for a GeodeticCRS should have 3 axis");
}
return DerivedGeodeticCRS::create(buildProperties(node), baseGeodCRS,
derivingConversion,
NN_NO_CHECK(cartesianCS));
Expand Down
6 changes: 4 additions & 2 deletions test/unit/test_c_api.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5621,7 +5621,8 @@ TEST_F(CApi, proj_create_derived_geographic_crs) {
" LENGTHUNIT[\"metre\",1]],\n"
" ENSEMBLEACCURACY[2.0]],\n"
" PRIMEM[\"Greenwich\",0,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433]]],\n"
" ANGLEUNIT[\"degree\",0.0174532925199433]],\n"
" ID[\"EPSG\",4326]],\n"
" DERIVINGCONVERSION[\"Pole rotation (GRIB convention)\",\n"
" METHOD[\"Pole rotation (GRIB convention)\"],\n"
" PARAMETER[\"Latitude of the southern pole (GRIB "
Expand Down Expand Up @@ -5685,7 +5686,8 @@ TEST_F(CApi, proj_create_derived_geographic_crs_netcdf_cf) {
" ELLIPSOID[\"GRS 1980\",6378137,298.257222101,\n"
" LENGTHUNIT[\"metre\",1]]],\n"
" PRIMEM[\"Greenwich\",0,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433]]],\n"
" ANGLEUNIT[\"degree\",0.0174532925199433]],\n"
" ID[\"EPSG\",4019]],\n"
" DERIVINGCONVERSION[\"Pole rotation (netCDF CF convention)\",\n"
" METHOD[\"Pole rotation (netCDF CF convention)\"],\n"
" PARAMETER[\"Grid north pole latitude (netCDF CF "
Expand Down
6 changes: 4 additions & 2 deletions test/unit/test_crs.cpp
View file
Open in desktop
O EDBE riginal file line number Diff line number Diff line change
Expand Up @@ -5711,7 +5711,8 @@ TEST(crs, derivedGeographicCRS_WKT2_2019) {
" ELLIPSOID[\"WGS 84\",6378137,298.257223563,\n"
" LENGTHUNIT[\"metre\",1]]],\n"
" PRIMEM[\"Greenwich\",0,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433]]],\n"
" ANGLEUNIT[\"degree\",0.0174532925199433]],\n"
" ID[\"EPSG\",4326]],\n"
" DERIVINGCONVERSION[\"Atlantic pole\",\n"
" METHOD[\"Pole rotation\"],\n"
" PARAMETER[\"Latitude of rotated pole\",52,\n"
Expand Down Expand Up @@ -5918,7 +5919,8 @@ TEST(crs, derivedGeodeticCRS_WKT2_2019) {
" ELLIPSOID[\"WGS 84\",6378137,298.257223563,\n"
" LENGTHUNIT[\"metre\",1]]],\n"
" PRIMEM[\"Greenwich\",0,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433]]],\n"
" ANGLEUNIT[\"degree\",0.0174532925199433]],\n"
" ID[\"EPSG\",4326]],\n"
" DERIVINGCONVERSION[\"Some conversion\",\n"
" METHOD[\"Some method\"]],\n"
" CS[Cartesian,3],\n"
Expand Down
69 changes: 69 additions & 0 deletions test/unit/test_io.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5683,6 +5683,75 @@ TEST(wkt_parse, DerivedGeodeticCRS) {

// ---------------------------------------------------------------------------

TEST(wkt_parse, DerivedGeodeticCRS_where_base_is_geocentric) {
auto wkt = "GEODCRS[\"Local CRS derived from WGS-84\",\n"
" BASEGEODCRS[\"WGS 84\",\n"
" ENSEMBLE[\"World Geodetic System 1984 ensemble\",\n"
" MEMBER[\"World Geodetic System 1984 (Transit)\"],\n"
" MEMBER[\"World Geodetic System 1984 (G730)\"],\n"
" MEMBER[\"World Geodetic System 1984 (G873)\"],\n"
" MEMBER[\"World Geodetic System 1984 (G1150)\"],\n"
" MEMBER[\"World Geodetic System 1984 (G1674)\"],\n"
" MEMBER[\"World Geodetic System 1984 (G1762)\"],\n"
" MEMBER[\"World Geodetic System 1984 (G2139)\"],\n"
" MEMBER[\"World Geodetic System 1984 (G2296)\"],\n"
" ELLIPSOID[\"WGS 84\",6378137,298.257223563,\n"
" LENGTHUNIT[\"metre\",1]],\n"
" ENSEMBLEACCURACY[2.0]],\n"
" PRIMEM[\"Greenwich\",0,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433]]],\n"
" DERIVINGCONVERSION[\"Local origin shift\",\n"
" METHOD[\"Position Vector transformation (geocentric "
"domain)\",\n"
" ID[\"EPSG\",1033]],\n"
" PARAMETER[\"X-axis translation\",10,\n"
" LENGTHUNIT[\"metre\",1],\n"
" ID[\"EPSG\",8605]],\n"
" PARAMETER[\"Y-axis translation\",20,\n"
" LENGTHUNIT[\"metre\",1],\n"
" ID[\"EPSG\",8606]],\n"
" PARAMETER[\"Z-axis translation\",1,\n"
" LENGTHUNIT[\"metre\",1],\n"
" ID[\"EPSG\",8607]],\n"
" PARAMETER[\"X-axis rotation\",0,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433],\n"
" ID[\"EPSG\",8608]],\n"
" PARAMETER[\"Y-axis rotation\",0,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433],\n"
" ID[\"EPSG\",8609]],\n"
" PARAMETER[\"Z-axis rotation\",0,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433],\n"
" ID[\"EPSG\",86 1E0A 10]]],\n"
" CS[Cartesian,3],\n"
" AXIS[\"(X)\",geocentricX,\n"
" ORDER[1],\n"
" LENGTHUNIT[\"metre\",1,\n"
" ID[\"EPSG\",9001]]],\n"
" AXIS[\"(Y)\",geocentricY,\n"
" ORDER[2],\n"
" LENGTHUNIT[\"metre\",1,\n"
" ID[\"EPSG\",9001]]],\n"
" AXIS[\"(Z)\",geocentricZ,\n"
" ORDER[3],\n"
" LENGTHUNIT[\"metre\",1,\n"
" ID[\"EPSG\",9001]]]]";

auto obj = WKTParser().createFromWKT(wkt);
auto crs = nn_dynamic_pointer_cast<DerivedGeodeticCRS>(obj);
ASSERT_TRUE(crs != nullptr);

auto baseCRS = nn_dynamic_pointer_cast<GeodeticCRS>(crs->baseCRS());
ASSERT_TRUE(baseCRS != nullptr);

EXPECT_TRUE(baseCRS->isGeocentric());

auto exportedWKT = crs->exportToWKT(
WKTFormatter::create(WKTFormatter::Convention::WKT2_2019).get());
EXPECT_STREQ(exportedWKT.c_str(), wkt);
}

// ---------------------------------------------------------------------------

TEST(wkt_parse, DerivedGeographicCRS_GDAL_PROJ4_EXSTENSION_hack) {
// Note the lack of UNIT[] node
auto wkt =
Expand Down
0