diff options
Diffstat (limited to 'src/iso19111/coordinateoperation.cpp')
| -rw-r--r-- | src/iso19111/coordinateoperation.cpp | 180 |
1 files changed, 112 insertions, 68 deletions
diff --git a/src/iso19111/coordinateoperation.cpp b/src/iso19111/coordinateoperation.cpp index 0f5e4c4d..3ffd208f 100644 --- a/src/iso19111/coordinateoperation.cpp +++ b/src/iso19111/coordinateoperation.cpp @@ -2043,15 +2043,19 @@ SingleOperation::gridsNeeded(const io::DatabaseContextPtr &databaseContext, if (opParamvalue) { const auto &value = opParamvalue->parameterValue(); if (value->type() == ParameterValue::Type::FILENAME) { - GridDescription desc; - desc.shortName = value->valueFile(); - if (databaseContext) { - databaseContext->lookForGridInfo( - desc.shortName, considerKnownGridsAsAvailable, - desc.fullName, desc.packageName, desc.url, - desc.directDownload, desc.openLicense, desc.available); + const auto gridNames = split(value->valueFile(), ","); + for (const auto &gridName : gridNames) { + GridDescription desc; + desc.shortName = gridName; + if (databaseContext) { + databaseContext->lookForGridInfo( + desc.shortName, considerKnownGridsAsAvailable, + desc.fullName, desc.packageName, desc.url, + desc.directDownload, desc.openLicense, + desc.available); + } + res.insert(desc); } - res.insert(desc); } } } @@ -6527,7 +6531,7 @@ struct Transformation::Private { TransformationPtr forwardOperation_{}; - TransformationNNPtr registerInv(util::BaseObjectNNPtr thisIn, + TransformationNNPtr registerInv(const Transformation *thisIn, TransformationNNPtr invTransform); }; //! @endcond @@ -7996,12 +8000,11 @@ createApproximateInverseIfPossible(const Transformation *op) { //! @cond Doxygen_Suppress TransformationNNPtr -Transformation::Private::registerInv(util::BaseObjectNNPtr thisIn, +Transformation::Private::registerInv(const Transformation *thisIn, TransformationNNPtr invTransform) { - invTransform->d->forwardOperation_ = - util::nn_dynamic_pointer_cast<Transformation>(thisIn); + invTransform->d->forwardOperation_ = thisIn->shallowClone().as_nullable(); invTransform->setHasBallparkTransformation( - invTransform->d->forwardOperation_->hasBallparkTransformation()); + thisIn->hasBallparkTransformation()); return invTransform; } //! @endcond @@ -8039,12 +8042,24 @@ TransformationNNPtr Transformation::inverseAsTransformation() const { parameterValueNumericAsSI(EPSG_CODE_PARAMETER_Y_AXIS_TRANSLATION); double z = parameterValueNumericAsSI(EPSG_CODE_PARAMETER_Z_AXIS_TRANSLATION); + auto properties = createPropertiesForInverse(this, false, false); return d->registerInv( - shared_from_this(), - createGeocentricTranslations( - createPropertiesForInverse(this, false, false), l_targetCRS, - l_sourceCRS, negate(x), negate(y), negate(z), - coordinateOperationAccuracies())); + this, create(properties, l_targetCRS, l_sourceCRS, nullptr, + createMethodMapNameEPSGCode( + useOperationMethodEPSGCodeIfPresent( + properties, methodEPSGCode)), + VectorOfParameters{ + createOpParamNameEPSGCode( + EPSG_CODE_PARAMETER_X_AXIS_TRANSLATION), + createOpParamNameEPSGCode( + EPSG_CODE_PARAMETER_Y_AXIS_TRANSLATION), + createOpParamNameEPSGCode( + EPSG_CODE_PARAMETER_Z_AXIS_TRANSLATION), + }, + createParams(common::Length(negate(x)), + common::Length(negate(y)), + common::Length(negate(z))), + coordinateOperationAccuracies())); } if (methodEPSGCode == EPSG_CODE_METHOD_MOLODENSKY || @@ -8062,14 +8077,14 @@ TransformationNNPtr Transformation::inverseAsTransformation() const { if (methodEPSGCode == EPSG_CODE_METHOD_ABRIDGED_MOLODENSKY) { return d->registerInv( - shared_from_this(), + this, createAbridgedMolodensky( createPropertiesForInverse(this, false, false), l_targetCRS, l_sourceCRS, negate(x), negate(y), negate(z), negate(da), negate(df), coordinateOperationAccuracies())); } else { return d->registerInv( - shared_from_this(), + this, createMolodensky(createPropertiesForInverse(this, false, false), l_targetCRS, l_sourceCRS, negate(x), negate(y), negate(z), negate(da), negate(df), @@ -8082,10 +8097,9 @@ TransformationNNPtr Transformation::inverseAsTransformation() const { parameterValueMeasure(EPSG_CODE_PARAMETER_LONGITUDE_OFFSET); const common::Angle newOffset(negate(offset.value()), offset.unit()); return d->registerInv( - shared_from_this(), - createLongitudeRotation( - createPropertiesForInverse(this, false, false), l_targetCRS, - l_sourceCRS, newOffset)); + this, createLongitudeRotation( + createPropertiesForInverse(this, false, false), + l_targetCRS, l_sourceCRS, newOffset)); } if (methodEPSGCode == EPSG_CODE_METHOD_GEOGRAPHIC2D_OFFSETS) { @@ -8100,11 +8114,10 @@ TransformationNNPtr Transformation::inverseAsTransformation() const { offsetLong.unit()); return d->registerInv( - shared_from_this(), - createGeographic2DOffsets( - createPropertiesForInverse(this, false, false), l_targetCRS, - l_sourceCRS, newOffsetLat, newOffsetLong, - coordinateOperationAccuracies())); + this, createGeographic2DOffsets( + createPropertiesForInverse(this, false, false), + l_targetCRS, l_sourceCRS, newOffsetLat, newOffsetLong, + coordinateOperationAccuracies())); } if (methodEPSGCode == EPSG_CODE_METHOD_GEOGRAPHIC3D_OFFSETS) { @@ -8124,11 +8137,10 @@ TransformationNNPtr Transformation::inverseAsTransformation() const { offsetHeight.unit()); return d->registerInv( - shared_from_this(), - createGeographic3DOffsets( - createPropertiesForInverse(this, false, false), l_targetCRS, - l_sourceCRS, newOffsetLat, newOffsetLong, newOffsetHeight, - coordinateOperationAccuracies())); + this, createGeographic3DOffsets( + createPropertiesForInverse(this, false, false), + l_targetCRS, l_sourceCRS, newOffsetLat, newOffsetLong, + newOffsetHeight, coordinateOperationAccuracies())); } if (methodEPSGCode == EPSG_CODE_METHOD_GEOGRAPHIC2D_WITH_HEIGHT_OFFSETS) { @@ -8148,11 +8160,10 @@ TransformationNNPtr Transformation::inverseAsTransformation() const { offsetHeight.unit()); return d->registerInv( - shared_from_this(), - createGeographic2DWithHeightOffsets( - createPropertiesForInverse(this, false, false), l_targetCRS, - l_sourceCRS, newOffsetLat, newOffsetLong, newOffsetHeight, - coordinateOperationAccuracies())); + this, createGeographic2DWithHeightOffsets( + createPropertiesForInverse(this, false, false), + l_targetCRS, l_sourceCRS, newOffsetLat, newOffsetLong, + newOffsetHeight, coordinateOperationAccuracies())); } if (methodEPSGCode == EPSG_CODE_METHOD_VERTICAL_OFFSET) { @@ -8163,7 +8174,7 @@ TransformationNNPtr Transformation::inverseAsTransformation() const { offsetHeight.unit()); return d->registerInv( - shared_from_this(), + this, createVerticalOffset(createPropertiesForInverse(this, false, false), l_targetCRS, l_sourceCRS, newOffsetHeight, coordinateOperationAccuracies())); @@ -8173,18 +8184,17 @@ TransformationNNPtr Transformation::inverseAsTransformation() const { const double convFactor = parameterValueNumericAsSI( EPSG_CODE_PARAMETER_UNIT_CONVERSION_SCALAR); return d->registerInv( - shared_from_this(), - createChangeVerticalUnit( - createPropertiesForInverse(this, false, false), l_targetCRS, - l_sourceCRS, common::Scale(1.0 / convFactor), - coordinateOperationAccuracies())); + this, createChangeVerticalUnit( + createPropertiesForInverse(this, false, false), + l_targetCRS, l_sourceCRS, common::Scale(1.0 / convFactor), + coordinateOperationAccuracies())); } #ifdef notdef // We don't need that currently, but we might... if (methodEPSGCode == EPSG_CODE_METHOD_HEIGHT_DEPTH_REVERSAL) { return d->registerInv( - shared_from_this(), + this, createHeightDepthReversal( createPropertiesForInverse(this, false, false), l_targetCRS, l_sourceCRS, coordinateOperationAccuracies())); @@ -8989,20 +8999,41 @@ static void ThrowExceptionNotGeodeticGeographic(const char *trfrm_name) { // --------------------------------------------------------------------------- -static void setupPROJGeodeticSourceCRS(io::PROJStringFormatter *formatter, - const crs::CRSNNPtr &crs, bool addPushV3, - const char *trfrm_name) { - auto sourceCRSGeog = dynamic_cast<const crs::GeographicCRS *>(crs.get()); - if (!sourceCRSGeog) { - auto compoundCRS = dynamic_cast<const crs::CompoundCRS *>(crs.get()); +// If crs is a geographic CRS, or a compound CRS of a geographic CRS, +// or a compoundCRS of a bound CRS of a geographic CRS, return that +// geographic CRS +static crs::GeographicCRSPtr +extractGeographicCRSIfGeographicCRSOrEquivalent(const crs::CRSNNPtr &crs) { + auto geogCRS = util::nn_dynamic_pointer_cast<crs::GeographicCRS>(crs); + if (!geogCRS) { + auto compoundCRS = util::nn_dynamic_pointer_cast<crs::CompoundCRS>(crs); if (compoundCRS) { const auto &components = compoundCRS->componentReferenceSystems(); if (!components.empty()) { - sourceCRSGeog = dynamic_cast<const crs::GeographicCRS *>( - components[0].get()); + geogCRS = util::nn_dynamic_pointer_cast<crs::GeographicCRS>( + components[0]); + if (!geogCRS) { + auto boundCRS = + util::nn_dynamic_pointer_cast<crs::BoundCRS>( + components[0]); + if (boundCRS) { + geogCRS = + util::nn_dynamic_pointer_cast<crs::GeographicCRS>( + boundCRS->baseCRS()); + } + } } } } + return geogCRS; +} + +// --------------------------------------------------------------------------- + +static void setupPROJGeodeticSourceCRS(io::PROJStringFormatter *formatter, + const crs::CRSNNPtr &crs, bool addPushV3, + const char *trfrm_name) { + auto sourceCRSGeog = extractGeographicCRSIfGeographicCRSOrEquivalent(crs); if (sourceCRSGeog) { formatter->startInversion(); sourceCRSGeog->_exportToPROJString(formatter); @@ -9030,17 +9061,7 @@ static void setupPROJGeodeticSourceCRS(io::PROJStringFormatter *formatter, static void setupPROJGeodeticTargetCRS(io::PROJStringFormatter *formatter, const crs::CRSNNPtr &crs, bool addPopV3, const char *trfrm_name) { - auto targetCRSGeog = dynamic_cast<const crs::GeographicCRS *>(crs.get()); - if (!targetCRSGeog) { - auto compoundCRS = dynamic_cast<const crs::CompoundCRS *>(crs.get()); - if (compoundCRS) { - const auto &components = compoundCRS->componentReferenceSystems(); - if (!components.empty()) { - targetCRSGeog = dynamic_cast<const crs::GeographicCRS *>( - components[0].get()); - } - } - } + auto targetCRSGeog = extractGeographicCRSIfGeographicCRSOrEquivalent(crs); if (targetCRSGeog) { formatter->addStep("cart"); formatter->setCurrentStepInverted(true); @@ -9540,14 +9561,14 @@ void Transformation::_exportToPROJString( : CTABLE2Filename; if (!hGridShiftFilename.empty()) { auto sourceCRSGeog = - dynamic_cast<const crs::GeographicCRS *>(sourceCRS().get()); + extractGeographicCRSIfGeographicCRSOrEquivalent(sourceCRS()); if (!sourceCRSGeog) { throw io::FormattingException( concat("Can apply ", methodName, " only to GeographicCRS")); } auto targetCRSGeog = - dynamic_cast<const crs::GeographicCRS *>(targetCRS().get()); + extractGeographicCRSIfGeographicCRSOrEquivalent(targetCRS()); if (!targetCRSGeog) { throw io::FormattingException( concat("Can apply ", methodName, " only to GeographicCRS")); @@ -13191,7 +13212,7 @@ CoordinateOperationFactory::Private::createOperations( return applyInverse(createOperations(targetCRS, sourceCRS, context)); } - // boundCRS to boundCRS using the same geographic hubCRS + // boundCRS to boundCRS if (boundSrc && boundDst) { createOperationsBoundToBound(sourceCRS, targetCRS, context, boundSrc, boundDst, res); @@ -14327,6 +14348,7 @@ void CoordinateOperationFactory::Private::createOperationsBoundToBound( ENTER_FUNCTION(); + // BoundCRS to BoundCRS of horizontal CRS using the same (geographic) hub const auto &hubSrc = boundSrc->hubCRS(); auto hubSrcGeog = dynamic_cast<const crs::GeographicCRS *>(hubSrc.get()); const auto &hubDst = boundDst->hubCRS(); @@ -14375,6 +14397,28 @@ void CoordinateOperationFactory::Private::createOperationsBoundToBound( } } + // BoundCRS to BoundCRS of vertical CRS using the same vertical datum + // ==> ignore the bound transformation + auto baseOfBoundSrcAsVertCRS = + dynamic_cast<crs::VerticalCRS *>(boundSrc->baseCRS().get()); + auto baseOfBoundDstAsVertCRS = + dynamic_cast<crs::VerticalCRS *>(boundDst->baseCRS().get()); + if (baseOfBoundSrcAsVertCRS && baseOfBoundDstAsVertCRS) { + const auto datumSrc = baseOfBoundSrcAsVertCRS->datum(); + const auto datumDst = baseOfBoundDstAsVertCRS->datum(); + if (datumSrc && datumDst && + datumSrc->nameStr() == datumDst->nameStr() && + (datumSrc->nameStr() != "unknown" || + boundSrc->transformation()->_isEquivalentTo( + boundDst->transformation().get(), + util::IComparable::Criterion::EQUIVALENT))) { + res = createOperations(boundSrc->baseCRS(), boundDst->baseCRS(), + context); + return; + } + } + + // BoundCRS to BoundCRS of vertical CRS auto vertCRSOfBaseOfBoundSrc = boundSrc->baseCRS()->extractVerticalCRS(); auto vertCRSOfBaseOfBoundDst = boundDst->baseCRS()->extractVerticalCRS(); if (hubSrcGeog && hubDstGeog && |