aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--src/iso19111/coordinateoperation.cpp51
-rw-r--r--src/iso19111/io.cpp5
-rw-r--r--test/unit/test_operation.cpp47
3 files changed, 93 insertions, 10 deletions
diff --git a/src/iso19111/coordinateoperation.cpp b/src/iso19111/coordinateoperation.cpp
index c51b4d3d..287611f8 100644
--- a/src/iso19111/coordinateoperation.cpp
+++ b/src/iso19111/coordinateoperation.cpp
@@ -13100,10 +13100,16 @@ std::vector<CoordinateOperationNNPtr> CoordinateOperationFactory::Private::
// which is likely/hopefully the only one.
for (const auto &opFirst : resTmp) {
if (hasIdentifiers(opFirst)) {
- res.emplace_back(
- ConcatenatedOperation::createComputeMetadata(
- {opFirst, opsSecond.front()},
- !allowEmptyIntersection));
+ if (candidateVert->_isEquivalentTo(
+ targetCRS.get(),
+ util::IComparable::Criterion::EQUIVALENT)) {
+ res.emplace_back(opFirst);
+ } else {
+ res.emplace_back(
+ ConcatenatedOperation::createComputeMetadata(
+ {opFirst, opsSecond.front()},
+ !allowEmptyIntersection));
+ }
}
}
if (!res.empty())
@@ -13456,6 +13462,24 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog(
hubSrcGeog->_isEquivalentTo(geogDst,
util::IComparable::Criterion::EQUIVALENT) &&
dynamic_cast<const crs::VerticalCRS *>(boundSrc->baseCRS().get())) {
+ auto transfSrc = boundSrc->transformation()->sourceCRS();
+ if (dynamic_cast<const crs::VerticalCRS *>(transfSrc.get()) &&
+ !boundSrc->baseCRS()->_isEquivalentTo(
+ transfSrc.get(), util::IComparable::Criterion::EQUIVALENT)) {
+ auto opsFirst =
+ createOperations(boundSrc->baseCRS(), transfSrc, context);
+ for (const auto &opFirst : opsFirst) {
+ try {
+ res.emplace_back(
+ ConcatenatedOperation::createComputeMetadata(
+ {opFirst, boundSrc->transformation()},
+ !allowEmptyIntersection));
+ } catch (const InvalidOperationEmptyIntersection &) {
+ }
+ }
+ return;
+ }
+
res.emplace_back(boundSrc->transformation());
return;
}
@@ -13590,7 +13614,7 @@ void CoordinateOperationFactory::Private::createOperationsVertToVert(
common::Scale(heightDepthReversal ? -factor : factor), {});
conv->setHasBallparkTransformation(true);
res.push_back(conv);
- } else if (convSrc != convDst) {
+ } else if (convSrc != convDst || !heightDepthReversal) {
auto conv = Conversion::createChangeVerticalUnit(
util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, name),
// In case of a height depth reversal, we should probably have
@@ -13598,7 +13622,7 @@ void CoordinateOperationFactory::Private::createOperationsVertToVert(
common::Scale(heightDepthReversal ? -factor : factor));
conv->setCRSs(sourceCRS, targetCRS, nullptr);
res.push_back(conv);
- } else if (heightDepthReversal) {
+ } else {
auto conv = Conversion::createHeightDepthReversal(
util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, name));
conv->setCRSs(sourceCRS, targetCRS, nullptr);
@@ -14131,10 +14155,13 @@ void CoordinateOperationFactory::Private::createOperationsCompoundToCompound(
if (componentsSrc[0]->extractGeographicCRS() &&
componentsDst[0]->extractGeographicCRS()) {
+ bool verticalTransfIsNoOp = false;
std::vector<CoordinateOperationNNPtr> verticalTransforms;
if (componentsSrc.size() >= 2 &&
componentsSrc[1]->extractVerticalCRS() &&
componentsDst[1]->extractVerticalCRS()) {
+ verticalTransfIsNoOp =
+ componentsSrc[1]->_isEquivalentTo(componentsDst[1].get());
verticalTransforms = createOperations(
componentsSrc[1], componentsDst[1], context);
}
@@ -14181,9 +14208,15 @@ void CoordinateOperationFactory::Private::createOperationsCompoundToCompound(
for (const auto &opDst : opGeogCRStoDstCRS) {
try {
- auto op = createHorizVerticalHorizPROJBased(
- sourceCRS, targetCRS, opSrc, verticalTransform,
- opDst, interpolationGeogCRS, true);
+ auto op = verticalTransfIsNoOp
+ ? ConcatenatedOperation::
+ createComputeMetadata(
+ {opSrc, opDst},
+ !allowEmptyIntersection)
+ : createHorizVerticalHorizPROJBased(
+ sourceCRS, targetCRS, opSrc,
+ verticalTransform, opDst,
+ interpolationGeogCRS, true);
res.emplace_back(op);
} catch (const InvalidOperationEmptyIntersection &) {
}
diff --git a/src/iso19111/io.cpp b/src/iso19111/io.cpp
index 645bec0b..a0a87f65 100644
--- a/src/iso19111/io.cpp
+++ b/src/iso19111/io.cpp
@@ -8364,7 +8364,10 @@ PROJStringParser::Private::buildBoundOrCompoundCRSIfNeeded(int iStep,
Transformation::createGravityRelatedHeightToGeographic3D(
PropertyMap().set(IdentifiedObject::NAME_KEY,
"unknown to WGS84 ellipsoidal height"),
- crs, GeographicCRS::EPSG_4979, nullptr, geoidgrids,
+ VerticalCRS::create(createMapWithUnknownName(), vdatum,
+ VerticalCS::createGravityRelatedHeight(
+ common::UnitOfMeasure::METRE)),
+ GeographicCRS::EPSG_4979, nullptr, geoidgrids,
std::vector<PositionalAccuracyNNPtr>());
auto boundvcrs =
BoundCRS::create(vcrs, GeographicCRS::EPSG_4979, transformation);
diff --git a/test/unit/test_operation.cpp b/test/unit/test_operation.cpp
index 07b0daea..ce01f967 100644
--- a/test/unit/test_operation.cpp
+++ b/test/unit/test_operation.cpp
@@ -6638,6 +6638,53 @@ TEST(operation, compoundCRS_with_boundProjCRS_and_boundVerticalCRS_to_geogCRS) {
// ---------------------------------------------------------------------------
+TEST(operation,
+ compoundCRS_with_boundVerticalCRS_from_geoidgrids_with_m_to_geogCRS) {
+
+ auto objSrc = PROJStringParser().createFromPROJString(
+ "+proj=longlat +datum=WGS84 +geoidgrids=@foo.gtx +type=crs");
+ auto src = nn_dynamic_pointer_cast<CRS>(objSrc);
+ ASSERT_TRUE(src != nullptr);
+ auto op = CoordinateOperationFactory::create()->createOperation(
+ NN_NO_CHECK(src), GeographicCRS::EPSG_4979);
+ ASSERT_TRUE(op != nullptr);
+ EXPECT_EQ(op->nameStr(), "axis order change (2D) + "
+ "unknown to WGS84 ellipsoidal height");
+ EXPECT_EQ(op->exportToPROJString(PROJStringFormatter::create().get()),
+ "+proj=pipeline "
+ "+step +proj=unitconvert +xy_in=deg +xy_out=rad "
+ "+step +proj=vgridshift +grids=@foo.gtx +multiplier=1 "
+ "+step +proj=unitconvert +xy_in=rad +xy_out=deg "
+ "+step +proj=axisswap +order=2,1");
+}
+
+// ---------------------------------------------------------------------------
+
+TEST(operation,
+ compoundCRS_with_boundVerticalCRS_from_geoidgrids_with_ftus_to_geogCRS) {
+
+ auto objSrc = PROJStringParser().createFromPROJString(
+ "+proj=longlat +datum=WGS84 +geoidgrids=@foo.gtx +vunits=us-ft "
+ "+type=crs");
+ auto src = nn_dynamic_pointer_cast<CRS>(objSrc);
+ ASSERT_TRUE(src != nullptr);
+ auto op = CoordinateOperationFactory::create()->createOperation(
+ NN_NO_CHECK(src), GeographicCRS::EPSG_4979);
+ ASSERT_TRUE(op != nullptr);
+ EXPECT_EQ(op->nameStr(), "axis order change (2D) + "
+ "Transformation from unknown to unknown + "
+ "unknown to WGS84 ellipsoidal height");
+ EXPECT_EQ(
+ op->exportToPROJString(PROJStringFormatter::create().get()),
+ "+proj=pipeline "
+ "+step +proj=unitconvert +xy_in=deg +z_in=us-ft +xy_out=rad +z_out=m "
+ "+step +proj=vgridshift +grids=@foo.gtx +multiplier=1 "
+ "+step +proj=unitconvert +xy_in=rad +xy_out=deg "
+ "+step +proj=axisswap +order=2,1");
+}
+
+// ---------------------------------------------------------------------------
+
TEST(operation, geocent_to_compoundCRS) {
auto objSrc = PROJStringParser().createFromPROJString(
"+proj=geocent +datum=WGS84 +units=m +type=crs");