From 4ee2ed65f728985e53ba0c76590f07b253b03587 Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Tue, 19 May 2020 01:31:25 +0200 Subject: createOperations(): fix 'caching' bugs causing potential exception about Inconsistent chainging of CRS (fixes #2232) --- src/iso19111/coordinateoperation.cpp | 4 ++++ 1 file changed, 4 insertions(+) (limited to 'src/iso19111/coordinateoperation.cpp') diff --git a/src/iso19111/coordinateoperation.cpp b/src/iso19111/coordinateoperation.cpp index 9815a22a..d1bc8eb5 100644 --- a/src/iso19111/coordinateoperation.cpp +++ b/src/iso19111/coordinateoperation.cpp @@ -6587,6 +6587,10 @@ TransformationNNPtr Transformation::shallowClone() const { auto transf = Transformation::nn_make_shared(*this); transf->assignSelf(transf); transf->setCRSs(this, false); + if (transf->d->forwardOperation_) { + transf->d->forwardOperation_ = + transf->d->forwardOperation_->shallowClone().as_nullable(); + } return transf; } -- cgit v1.2.3 From 166a1c681bba659095b1e1296cbb66e9f1637bbd Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Tue, 19 May 2020 01:33:28 +0200 Subject: createOperationsWithDatumPivot(): remove useless code --- src/iso19111/coordinateoperation.cpp | 11 ++--------- 1 file changed, 2 insertions(+), 9 deletions(-) (limited to 'src/iso19111/coordinateoperation.cpp') diff --git a/src/iso19111/coordinateoperation.cpp b/src/iso19111/coordinateoperation.cpp index d1bc8eb5..4bf5c41e 100644 --- a/src/iso19111/coordinateoperation.cpp +++ b/src/iso19111/coordinateoperation.cpp @@ -12946,7 +12946,6 @@ void CoordinateOperationFactory::Private::createOperationsWithDatumPivot( (isNullFirst || isNullThird || sourceAndTargetAre3D) ? opSecond->shallowClone() : opSecond); - CoordinateOperation *invCOForward = nullptr; if (isNullFirst || isNullThird) { if (opSecondCloned->identifiers().size() == 1 && (*opSecondCloned->identifiers()[0]->codeSpace()) @@ -12960,7 +12959,7 @@ void CoordinateOperationFactory::Private::createOperationsWithDatumPivot( auto invCO = dynamic_cast( opSecondCloned.get()); if (invCO) { - invCOForward = invCO->forwardOperation().get(); + auto invCOForward = invCO->forwardOperation().get(); if (invCOForward->identifiers().size() == 1 && (*invCOForward->identifiers()[0]->codeSpace()) .find("DERIVED_FROM") == @@ -12978,25 +12977,19 @@ void CoordinateOperationFactory::Private::createOperationsWithDatumPivot( auto invCO = dynamic_cast( opSecondCloned.get()); if (invCO) { - invCOForward = invCO->forwardOperation().get(); + auto invCOForward = invCO->forwardOperation().get(); invCOForward->getPrivate()->use3DHelmert_ = true; } } if (isNullFirst) { auto oldTarget(NN_CHECK_ASSERT(opSecondCloned->targetCRS())); setCRSs(opSecondCloned.get(), sourceCRS, oldTarget); - if (invCOForward) { - setCRSs(invCOForward, oldTarget, sourceCRS); - } } else { subOps.emplace_back(opFirst); } if (isNullThird) { auto oldSource(NN_CHECK_ASSERT(opSecondCloned->sourceCRS())); setCRSs(opSecondCloned.get(), oldSource, targetCRS); - if (invCOForward) { - setCRSs(invCOForward, targetCRS, oldSource); - } subOps.emplace_back(opSecondCloned); } else { subOps.emplace_back(opSecondCloned); -- cgit v1.2.3 From 1a715234754146ebe224fb849a87ca6575fdc88f Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Tue, 19 May 2020 01:46:28 +0200 Subject: createOperations(): speed optimizations for transforming between a BoundCRS of a datum and the same datum (relates to #2232) --- src/iso19111/coordinateoperation.cpp | 22 ++++++++++++++++++---- 1 file changed, 18 insertions(+), 4 deletions(-) (limited to 'src/iso19111/coordinateoperation.cpp') diff --git a/src/iso19111/coordinateoperation.cpp b/src/iso19111/coordinateoperation.cpp index 4bf5c41e..231d31a0 100644 --- a/src/iso19111/coordinateoperation.cpp +++ b/src/iso19111/coordinateoperation.cpp @@ -13997,6 +13997,20 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog( const auto &hubSrc = boundSrc->hubCRS(); auto hubSrcGeog = dynamic_cast(hubSrc.get()); auto geogCRSOfBaseOfBoundSrc = boundSrc->baseCRS()->extractGeographicCRS(); + auto geogDstDatum = geogDst->datum(); + + // If the underlying datum of the source is the same as the target, do + // not consider the boundCRS at all, but just its base + if (geogCRSOfBaseOfBoundSrc && geogDstDatum) { + auto geogCRSOfBaseOfBoundSrcDatum = geogCRSOfBaseOfBoundSrc->datum(); + if (geogCRSOfBaseOfBoundSrcDatum && + geogCRSOfBaseOfBoundSrcDatum->_isEquivalentTo( + geogDstDatum.get(), util::IComparable::Criterion::EQUIVALENT)) { + res = createOperations(boundSrc->baseCRS(), targetCRS, context); + return; + } + } + bool triedBoundCrsToGeogCRSSameAsHubCRS = false; // Is it: boundCRS to a geogCRS that is the same as the hubCRS ? if (hubSrcGeog && geogCRSOfBaseOfBoundSrc && @@ -14043,9 +14057,9 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog( } // If the datum are equivalent, this is also fine } else if (geogCRSOfBaseOfBoundSrc && hubSrcGeog && hubSrcGeog->datum() && - geogDst->datum() && + geogDstDatum && hubSrcGeog->datum()->_isEquivalentTo( - geogDst->datum().get(), + geogDstDatum.get(), util::IComparable::Criterion::EQUIVALENT)) { auto opsFirst = createOperations( boundSrc->baseCRS(), NN_NO_CHECK(geogCRSOfBaseOfBoundSrc), context); @@ -14088,14 +14102,14 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog( // +nadgrids=ntv1_can.dat,conus" // to "+proj=latlong +datum=NAD83" } else if (geogCRSOfBaseOfBoundSrc && hubSrcGeog && hubSrcGeog->datum() && - geogDst->datum() && + geogDstDatum && geogCRSOfBaseOfBoundSrc->ellipsoid()->_isEquivalentTo( datum::Ellipsoid::CLARKE_1866.get(), util::IComparable::Criterion::EQUIVALENT) && hubSrcGeog->datum()->_isEquivalentTo( datum::GeodeticReferenceFrame::EPSG_6326.get(), util::IComparable::Criterion::EQUIVALENT) && - geogDst->datum()->_isEquivalentTo( + geogDstDatum->_isEquivalentTo( datum::GeodeticReferenceFrame::EPSG_6269.get(), util::IComparable::Criterion::EQUIVALENT)) { auto nnGeogCRSOfBaseOfBoundSrc = NN_NO_CHECK(geogCRSOfBaseOfBoundSrc); -- cgit v1.2.3 From 8dc8e1de06dcb71da11facaaff6af9e46294deb0 Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Tue, 19 May 2020 12:57:22 +0200 Subject: createOperations(): fix wrong optimization in some instances of transforming a BoundCRS to a GeographicCRS (contributes to fixes #2232) --- src/iso19111/coordinateoperation.cpp | 39 ++++++++++++++++++++++++------------ 1 file changed, 26 insertions(+), 13 deletions(-) (limited to 'src/iso19111/coordinateoperation.cpp') diff --git a/src/iso19111/coordinateoperation.cpp b/src/iso19111/coordinateoperation.cpp index 231d31a0..32a4b9c7 100644 --- a/src/iso19111/coordinateoperation.cpp +++ b/src/iso19111/coordinateoperation.cpp @@ -14018,25 +14018,38 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog( geogDst, util::IComparable::Criterion::EQUIVALENT) || hubSrcGeog->is2DPartOf3D(NN_NO_CHECK(geogDst)))) { triedBoundCrsToGeogCRSSameAsHubCRS = true; + + CoordinateOperationPtr opIntermediate; + if (!geogCRSOfBaseOfBoundSrc->_isEquivalentTo( + boundSrc->transformation()->sourceCRS().get(), + util::IComparable::Criterion::EQUIVALENT)) { + auto opsIntermediate = createOperations( + NN_NO_CHECK(geogCRSOfBaseOfBoundSrc), + boundSrc->transformation()->sourceCRS(), context); + assert(!opsIntermediate.empty()); + opIntermediate = opsIntermediate.front(); + } + if (boundSrc->baseCRS() == geogCRSOfBaseOfBoundSrc) { - // Optimization to avoid creating a useless concatenated - // operation - res.emplace_back(boundSrc->transformation()); + if (opIntermediate) { + try { + res.emplace_back( + ConcatenatedOperation::createComputeMetadata( + {NN_NO_CHECK(opIntermediate), + boundSrc->transformation()}, + !allowEmptyIntersection)); + } catch (const InvalidOperationEmptyIntersection &) { + } + } else { + // Optimization to avoid creating a useless concatenated + // operation + res.emplace_back(boundSrc->transformation()); + } return; } auto opsFirst = createOperations( boundSrc->baseCRS(), NN_NO_CHECK(geogCRSOfBaseOfBoundSrc), context); if (!opsFirst.empty()) { - CoordinateOperationPtr opIntermediate; - if (!geogCRSOfBaseOfBoundSrc->_isEquivalentTo( - boundSrc->transformation()->sourceCRS().get(), - util::IComparable::Criterion::EQUIVALENT)) { - auto opsIntermediate = createOperations( - NN_NO_CHECK(geogCRSOfBaseOfBoundSrc), - boundSrc->transformation()->sourceCRS(), context); - assert(!opsIntermediate.empty()); - opIntermediate = opsIntermediate.front(); - } for (const auto &opFirst : opsFirst) { try { std::vector subops; -- cgit v1.2.3 From 5aa99d497433377fc08cb569b6826473dc0c46d1 Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Tue, 19 May 2020 13:15:40 +0200 Subject: createOperations(): take correctly into account vertical unit change in a case of transformation between Compound of BoundVerticalCRS to GeographicCRS (contributes to fixes #2232) --- src/iso19111/coordinateoperation.cpp | 53 ++++++++++++++++++++++++++++++++++-- 1 file changed, 51 insertions(+), 2 deletions(-) (limited to 'src/iso19111/coordinateoperation.cpp') diff --git a/src/iso19111/coordinateoperation.cpp b/src/iso19111/coordinateoperation.cpp index 32a4b9c7..6dff57c1 100644 --- a/src/iso19111/coordinateoperation.cpp +++ b/src/iso19111/coordinateoperation.cpp @@ -14221,8 +14221,57 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog( if (vertCRSOfBaseOfBoundSrc && hubSrcGeog) { auto opsFirst = createOperations(sourceCRS, hubSrc, context); if (context.skipHorizontalTransformation) { - if (!opsFirst.empty()) - res = opsFirst; + if (!opsFirst.empty()) { + const auto &hubAxisList = + hubSrcGeog->coordinateSystem()->axisList(); + const auto &targetAxisList = + geogDst->coordinateSystem()->axisList(); + if (hubAxisList.size() == 3 && hubAxisList.size() == 3 && + !hubAxisList[2]->_isEquivalentTo( + targetAxisList[2].get(), + util::IComparable::Criterion::EQUIVALENT)) { + + const auto &srcAxis = hubAxisList[2]; + const double convSrc = srcAxis->unit().conversionToSI(); + const auto &dstAxis = targetAxisList[2]; + const double convDst = dstAxis->unit().conversionToSI(); + const bool srcIsUp = + srcAxis->direction() == cs::AxisDirection::UP; + const bool srcIsDown = + srcAxis->direction() == cs::AxisDirection::DOWN; + const bool dstIsUp = + dstAxis->direction() == cs::AxisDirection::UP; + const bool dstIsDown = + dstAxis->direction() == cs::AxisDirection::DOWN; + const bool heightDepthReversal = + ((srcIsUp && dstIsDown) || (srcIsDown && dstIsUp)); + + const double factor = convSrc / convDst; + auto conv = Conversion::createChangeVerticalUnit( + util::PropertyMap().set( + common::IdentifiedObject::NAME_KEY, + "Change of vertical unit"), + common::Scale(heightDepthReversal ? -factor : factor)); + auto dbContext = context.context->getAuthorityFactory() + ->databaseContext(); + conv->setCRSs( + hubSrc, + hubSrc->demoteTo2D(std::string(), dbContext) + ->promoteTo3D(std::string(), dbContext, dstAxis), + nullptr); + + for (const auto &op : opsFirst) { + try { + res.emplace_back( + ConcatenatedOperation::createComputeMetadata( + {op, conv}, !allowEmptyIntersection)); + } catch (const InvalidOperationEmptyIntersection &) { + } + } + } else { + res = opsFirst; + } + } return; } else { auto opsSecond = createOperations(hubSrc, targetCRS, context); -- cgit v1.2.3 From b0b33b8447972ac6e60d68213d6c24b0a4989421 Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Tue, 19 May 2020 13:38:23 +0200 Subject: createOperations(): fix bug when transforming between CompoundCRS and BoundCRS in the general case let to 0 result (contributes to fixes #2232) --- src/iso19111/coordinateoperation.cpp | 5 +++++ 1 file changed, 5 insertions(+) (limited to 'src/iso19111/coordinateoperation.cpp') diff --git a/src/iso19111/coordinateoperation.cpp b/src/iso19111/coordinateoperation.cpp index 6dff57c1..bdb2ad2e 100644 --- a/src/iso19111/coordinateoperation.cpp +++ b/src/iso19111/coordinateoperation.cpp @@ -15218,10 +15218,15 @@ void CoordinateOperationFactory::Private::createOperationsBoundToCompound( } } } + return; } } } } + + // There might be better things to do, but for now just ignore the + // transformation of the bound CRS + res = createOperations(boundSrc->baseCRS(), targetCRS, context); } //! @endcond -- cgit v1.2.3