Skip to content

Commit

Permalink
Merge pull request #3102 from rouault/fix_3076
Browse files Browse the repository at this point in the history 
Fix issues with WKT of concatenated operations (#3076)
  • Loading branch information
@rouault
rouault committed Mar 9, 2022
2 parents f026f99 + 273558e commit 99570fb
Show file tree
Hide file tree
Showing 3 changed files with 220 additions and 15 deletions.
24 changes: 18 additions & 6 deletions src/iso19111/operation/concatenatedoperation.cpp
View file
Expand Up @@ -276,14 +276,28 @@ void ConcatenatedOperation::fixStepsDirection(
}
}

const auto extractDerivedCRS =
[](const crs::CRS *crs) -> const crs::DerivedCRS * {
auto derivedCRS = dynamic_cast<const crs::DerivedCRS *>(crs);
if (derivedCRS)
return derivedCRS;
auto compoundCRS = dynamic_cast<const crs::CompoundCRS *>(crs);
if (compoundCRS) {
derivedCRS = dynamic_cast<const crs::DerivedCRS *>(
compoundCRS->componentReferenceSystems().front().get());
if (derivedCRS)
return derivedCRS;
}
return nullptr;
};

for (size_t i = 0; i < operationsInOut.size(); ++i) {
auto &op = operationsInOut[i];
auto l_sourceCRS = op->sourceCRS();
auto l_targetCRS = op->targetCRS();
auto conv = dynamic_cast<const Conversion *>(op.get());
if (conv && i == 0 && !l_sourceCRS && !l_targetCRS) {
if (auto derivedCRS = dynamic_cast<const crs::DerivedCRS *>(
concatOpSourceCRS.get())) {
if (auto derivedCRS = extractDerivedCRS(concatOpSourceCRS.get())) {
if (i + 1 < operationsInOut.size()) {
// use the sourceCRS of the next operation as our target CRS
l_targetCRS = operationsInOut[i + 1]->sourceCRS();
Expand Down Expand Up @@ -323,8 +337,7 @@ void ConcatenatedOperation::fixStepsDirection(
}
} else if (conv && i + 1 == operationsInOut.size() && !l_sourceCRS &&
!l_targetCRS) {
auto derivedCRS =
dynamic_cast<const crs::DerivedCRS *>(concatOpTargetCRS.get());
auto derivedCRS = extractDerivedCRS(concatOpTargetCRS.get());
if (derivedCRS) {
if (i >= 1) {
// use the sourceCRS of the previous operation as our source
Expand All @@ -350,8 +363,7 @@ void ConcatenatedOperation::fixStepsDirection(
} else if (i >= 1) {
l_sourceCRS = operationsInOut[i - 1]->targetCRS();
if (l_sourceCRS) {
derivedCRS = dynamic_cast<const crs::DerivedCRS *>(
l_sourceCRS.get());
derivedCRS = extractDerivedCRS(l_sourceCRS.get());
if (derivedCRS &&
conv->isEquivalentTo(
derivedCRS->derivingConversion().get(),
Expand Down
16 changes: 7 additions & 9 deletions src/iso19111/operation/projbasedoperation.cpp
View file
Expand Up @@ -232,15 +232,13 @@ void PROJBasedOperation::_exportToJSON(
method()->_exportToJSON(formatter);

const auto &l_parameterValues = parameterValues();
if (!l_parameterValues.empty()) {
writer->AddObjKey("parameters");
{
auto parametersContext(writer->MakeArrayContext(false));
for (const auto &genOpParamvalue : l_parameterValues) {
formatter->setAllowIDInImmediateChild();
formatter->setOmitTypeInImmediateChild();
genOpParamvalue->_exportToJSON(formatter);
}
writer->AddObjKey("parameters");
{
auto parametersContext(writer->MakeArrayContext(false));
for (const auto &genOpParamvalue : l_parameterValues) {
formatter->setAllowIDInImmediateChild();
formatter->setOmitTypeInImmediateChild();
genOpParamvalue->_exportToJSON(formatter);
}
}
}
Expand Down
195 changes: 195 additions & 0 deletions test/unit/test_io.cpp
View file
Expand Up @@ -3827,6 +3827,201 @@ TEST(

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

TEST(wkt_parse, CONCATENATEDOPERATION_with_inverse_conversion_of_compound) {

auto wkt =
"CONCATENATEDOPERATION[\"Inverse of RD New + Amersfoort to ETRS89 (9) "
"+ Inverse of ETRS89 to NAP height (2) + ETRS89 to WGS 84 (1)\",\n"
" SOURCECRS[\n"
" COMPOUNDCRS[\"Amersfoort / RD New + NAP height\",\n"
" PROJCRS[\"Amersfoort / RD New\",\n"
" BASEGEOGCRS[\"Amersfoort\",\n"
" DATUM[\"Amersfoort\",\n"
" ELLIPSOID[\"Bessel "
"1841\",6377397.155,299.1528128,\n"
" LENGTHUNIT[\"metre\",1]]],\n"
" PRIMEM[\"Greenwich\",0,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433]],\n"
" ID[\"EPSG\",4289]],\n"
" CONVERSION[\"RD New\",\n"
" METHOD[\"Oblique Stereographic\",\n"
" ID[\"EPSG\",9809]],\n"
" PARAMETER[\"Latitude of natural "
"origin\",52.1561605555556,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433],\n"
" ID[\"EPSG\",8801]],\n"
" PARAMETER[\"Longitude of natural "
"origin\",5.38763888888889,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433],\n"
" ID[\"EPSG\",8802]],\n"
" PARAMETER[\"Scale factor at natural "
"origin\",0.9999079,\n"
" SCALEUNIT[\"unity\",1],\n"
" ID[\"EPSG\",8805]],\n"
" PARAMETER[\"False easting\",155000,\n"
" LENGTHUNIT[\"metre\",1],\n"
" ID[\"EPSG\",8806]],\n"
" PARAMETER[\"False northing\",463000,\n"
" LENGTHUNIT[\"metre\",1],\n"
" ID[\"EPSG\",8807]]],\n"
" CS[Cartesian,2],\n"
" AXIS[\"easting (X)\",east,\n"
" ORDER[1],\n"
" LENGTHUNIT[\"metre\",1]],\n"
" AXIS[\"northing (Y)\",north,\n"
" ORDER[2],\n"
" LENGTHUNIT[\"metre\",1]]],\n"
" VERTCRS[\"NAP height\",\n"
" VDATUM[\"Normaal Amsterdams Peil\"],\n"
" CS[vertical,1],\n"
" AXIS[\"gravity-related height (H)\",up,\n"
" LENGTHUNIT[\"metre\",1]]],\n"
" ID[\"EPSG\",7415]]],\n"
" TARGETCRS[\n"
" GEOGCRS[\"WGS 84 (3D)\",\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"
" ELLIPSOID[\"WGS 84\",6378137,298.257223563,\n"
" LENGTHUNIT[\"metre\",1]],\n"
" ENSEMBLEACCURACY[2.0]],\n"
" PRIMEM[\"Greenwich\",0,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433]],\n"
" CS[ellipsoidal,3],\n"
" AXIS[\"geodetic latitude (Lat)\",north,\n"
" ORDER[1],\n"
" ANGLEUNIT[\"degree minute second "
"hemisphere\",0.0174532925199433]],\n"
" AXIS[\"geodetic longitude (Long)\",east,\n"
" ORDER[2],\n"
" ANGLEUNIT[\"degree minute second "
"hemisphere\",0.0174532925199433]],\n"
" AXIS[\"ellipsoidal height (h)\",up,\n"
" ORDER[3],\n"
" LENGTHUNIT[\"metre\",1]],\n"
" ID[\"EPSG\",4329]]],\n"
" STEP[\n"
" CONVERSION[\"Inverse of RD New\",\n"
" METHOD[\"Inverse of Oblique Stereographic\",\n"
" ID[\"INVERSE(EPSG)\",9809]],\n"
" PARAMETER[\"Latitude of natural "
"origin\",52.1561605555556,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433],\n"
" ID[\"EPSG\",8801]],\n"
" PARAMETER[\"Longitude of natural "
"origin\",5.38763888888889,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433],\n"
" ID[\"EPSG\",8802]],\n"
" PARAMETER[\"Scale factor at natural origin\",0.9999079,\n"
" SCALEUNIT[\"unity\",1],\n"
" ID[\"EPSG\",8805]],\n"
" PARAMETER[\"False easting\",155000,\n"
" LENGTHUNIT[\"metre\",1],\n"
" ID[\"EPSG\",8806]],\n"
" PARAMETER[\"False northing\",463000,\n"
" LENGTHUNIT[\"metre\",1],\n"
" ID[\"EPSG\",8807]],\n"
" ID[\"INVERSE(EPSG)\",19914]]],\n"
" STEP[\n"
" COORDINATEOPERATION[\"PROJ-based coordinate operation\",\n"
" SOURCECRS[\n"
" COMPOUNDCRS[\"Amersfoort + NAP height\",\n"
" GEOGCRS[\"Amersfoort\",\n"
" DATUM[\"Amersfoort\",\n"
" ELLIPSOID[\"Bessel "
"1841\",6377397.155,299.1528128,\n"
" LENGTHUNIT[\"metre\",1]]],\n"
" PRIMEM[\"Greenwich\",0,\n"
" "
"ANGLEUNIT[\"degree\",0.0174532925199433]],\n"
" CS[ellipsoidal,2],\n"
" AXIS[\"geodetic latitude (Lat)\",north,\n"
" ORDER[1],\n"
" "
"ANGLEUNIT[\"degree\",0.0174532925199433]],\n"
" AXIS[\"geodetic longitude (Lon)\",east,\n"
" ORDER[2],\n"
" "
"ANGLEUNIT[\"degree\",0.0174532925199433]],\n"
" ID[\"EPSG\",4289]],\n"
" VERTCRS[\"NAP height\",\n"
" VDATUM[\"Normaal Amsterdams Peil\"],\n"
" CS[vertical,1],\n"
" AXIS[\"gravity-related height (H)\",up,\n"
" LENGTHUNIT[\"metre\",1]],\n"
" ID[\"EPSG\",5709]]]],\n"
" TARGETCRS[\n"
" GEOGCRS[\"WGS 84 (3D)\",\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"
" ELLIPSOID[\"WGS 84\",6378137,298.257223563,\n"
" LENGTHUNIT[\"metre\",1]],\n"
" ENSEMBLEACCURACY[2.0]],\n"
" PRIMEM[\"Greenwich\",0,\n"
" ANGLEUNIT[\"degree\",0.0174532925199433]],\n"
" CS[ellipsoidal,3],\n"
" AXIS[\"geodetic latitude (Lat)\",north,\n"
" ORDER[1],\n"
" ANGLEUNIT[\"degree minute second "
"hemisphere\",0.0174532925199433]],\n"
" AXIS[\"geodetic longitude (Long)\",east,\n"
" ORDER[2],\n"
" ANGLEUNIT[\"degree minute second "
"hemisphere\",0.0174532925199433]],\n"
" AXIS[\"ellipsoidal height (h)\",up,\n"
" ORDER[3],\n"
" LENGTHUNIT[\"metre\",1]],\n"
" ID[\"EPSG\",4329]]],\n"
" METHOD[\"PROJ-based operation method: +proj=pipeline "
"+step +proj=axisswap +order=2,1 +step +proj=unitconvert +xy_in=deg "
"+xy_out=rad +step +proj=hgridshift +grids=nl_nsgi_rdtrans2018.tif "
"+step +proj=vgridshift +grids=nl_nsgi_nlgeo2018.tif +multiplier=1 "
"+step +proj=unitconvert +xy_in=rad +xy_out=deg +step +proj=axisswap "
"+order=2,1\"],\n"
" OPERATIONACCURACY[1.002]]],\n"
" USAGE[\n"
" SCOPE[\"unknown\"],\n"
" AREA[\"Netherlands - onshore, including Waddenzee, Dutch "
"Wadden Islands and 12-mile offshore coastal zone.\"],\n"
" BBOX[50.75,3.2,53.7,7.22]]]";

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

EXPECT_EQ(concat->exportToPROJString(PROJStringFormatter::create().get()),
"+proj=pipeline "
"+step +inv +proj=sterea +lat_0=52.1561605555556 "
"+lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 "
"+ellps=bessel "
"+step +proj=hgridshift +grids=nl_nsgi_rdtrans2018.tif "
"+step +proj=vgridshift +grids=nl_nsgi_nlgeo2018.tif "
"+multiplier=1 "
"+step +proj=unitconvert +xy_in=rad +xy_out=deg "
"+step +proj=axisswap +order=2,1");
}

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

TEST(wkt_parse, BOUNDCRS_transformation_from_names) {

auto projcrs = ProjectedCRS::create(
Expand Down

0 comments on commit 99570fb

Please sign in to comment.