From 78d5d84ce0e4c1b75847295f4009fe53f1d96dec Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Tue, 17 Sep 2024 16:17:49 +0200 Subject: [PATCH] findsOpsInRegistryWithIntermediate(): restrict to using known source/target CRS that have the same originating authority Fixes https://lists.osgeo.org/pipermail/proj/2024-September/011531.html i.e. ``` projinfo -s EPSG:4269 -t EPSG:6318 --3d --spatial-test intersects ``` The fix consists in making sure that we recognize that the 3d-promoted object EPSG:4269 (NAD83(86)) is still linked to EPSG, and thus discard ESRI 3D objects in findsOpsInRegistryWithIntermediate() Fixes a "regression" of https://github.com/OSGeo/PROJ/pull/4244 (one could argue which of the results is the best, given that NAD83(86) as a 3D geographic CRS has no solid foundation) --- include/proj/crs.hpp | 3 ++ src/iso19111/crs.cpp | 40 ++++++++++++++++- .../operation/coordinateoperationfactory.cpp | 12 ++++- test/unit/test_operationfactory.cpp | 45 +++++++++++++++++++ 4 files changed, 98 insertions(+), 2 deletions(-) diff --git a/include/proj/crs.hpp b/include/proj/crs.hpp index a6efcd4fe1..3d0d9a3179 100644 --- a/include/proj/crs.hpp +++ b/include/proj/crs.hpp @@ -165,6 +165,9 @@ class PROJ_GCC_DLL CRS : public common::ObjectUsage, getResolvedCRS(const CRSNNPtr &crs, const io::AuthorityFactoryPtr &authFactory, metadata::ExtentPtr &extentOut); + + PROJ_INTERNAL std::string getOriginatingAuthName() const; + //! @endcond protected: diff --git a/src/iso19111/crs.cpp b/src/iso19111/crs.cpp index df606a41c4..e4d2d98e55 100644 --- a/src/iso19111/crs.cpp +++ b/src/iso19111/crs.cpp @@ -92,6 +92,10 @@ NS_PROJ_START namespace crs { +//! @cond Doxygen_Suppress +constexpr const char *PROMOTED_TO_3D_PRELUDE = "Promoted to 3D from "; +//! @endcond + // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress @@ -1257,6 +1261,40 @@ CRS::getNonDeprecated(const io::DatabaseContextNNPtr &dbContext) const { // --------------------------------------------------------------------------- +//! @cond Doxygen_Suppress + +/** \brief Return the authority name to which this object is registered, or + * has an indirect provenance. + * + * Typically this method called on EPSG:4269 (NAD83) promoted to 3D will return + * "EPSG". + * + * Returns empty string if more than an authority or no originating authority is + * found. + */ +std::string CRS::getOriginatingAuthName() const { + const auto &ids = identifiers(); + if (ids.size() == 1) { + return *(ids[0]->codeSpace()); + } + if (ids.size() > 1) { + return std::string(); + } + const auto &l_remarks = remarks(); + if (starts_with(l_remarks, PROMOTED_TO_3D_PRELUDE)) { + const auto pos = l_remarks.find(':'); + if (pos != std::string::npos) { + return l_remarks.substr(strlen(PROMOTED_TO_3D_PRELUDE), + pos - strlen(PROMOTED_TO_3D_PRELUDE)); + } + } + return std::string(); +} + +//! @endcond + +// --------------------------------------------------------------------------- + /** \brief Return a variant of this CRS "promoted" to a 3D one, if not already * the case. * @@ -1314,7 +1352,7 @@ CRSNNPtr CRS::promoteTo3D(const std::string &newName, const auto &l_identifiers = identifiers(); const auto &l_remarks = remarks(); if (l_identifiers.size() == 1) { - std::string remarks("Promoted to 3D from "); + std::string remarks(PROMOTED_TO_3D_PRELUDE); remarks += *(l_identifiers[0]->codeSpace()); remarks += ':'; remarks += l_identifiers[0]->code(); diff --git a/src/iso19111/operation/coordinateoperationfactory.cpp b/src/iso19111/operation/coordinateoperationfactory.cpp index 8e73713eee..0e1f35c997 100644 --- a/src/iso19111/operation/coordinateoperationfactory.cpp +++ b/src/iso19111/operation/coordinateoperationfactory.cpp @@ -2001,12 +2001,21 @@ CoordinateOperationFactory::Private::findsOpsInRegistryWithIntermediate( auto geodSrc = dynamic_cast(sourceCRS.get()); auto geodDst = dynamic_cast(targetCRS.get()); if (geodSrc) { + const std::string originatingAuthSrc = + geodSrc->getOriginatingAuthName(); const auto dbContext = authFactory->databaseContext().as_nullable(); const auto candidatesSrcGeod(findCandidateGeodCRSForDatum( authFactory, geodSrc, geodSrc->datumNonNull(dbContext))); std::vector res; for (const auto &candidateSrcGeod : candidatesSrcGeod) { - if (candidateSrcGeod->coordinateSystem()->axisList().size() == + // Restrict to using only objects that have the same authority + // as geodSrc. + const auto &candidateIds = candidateSrcGeod->identifiers(); + if ((originatingAuthSrc.empty() || + (candidateIds.size() == 1 && + *(candidateIds[0]->codeSpace()) == originatingAuthSrc) || + candidateIds.size() > 1) && + candidateSrcGeod->coordinateSystem()->axisList().size() == geodSrc->coordinateSystem()->axisList().size() && ((dynamic_cast(sourceCRS.get()) != nullptr) == @@ -2022,6 +2031,7 @@ CoordinateOperationFactory::Private::findsOpsInRegistryWithIntermediate( continue; } } + const auto opsWithIntermediate = findsOpsInRegistryWithIntermediate( candidateSrcGeod, targetCRS, context, diff --git a/test/unit/test_operationfactory.cpp b/test/unit/test_operationfactory.cpp index 6bdf952de0..423fadb945 100644 --- a/test/unit/test_operationfactory.cpp +++ b/test/unit/test_operationfactory.cpp @@ -1533,6 +1533,51 @@ TEST(operation, geogCRS_without_id_to_geogCRS_3D_context) { // --------------------------------------------------------------------------- +TEST(operation, geogCRS_promoted_to_3D_to_geogCRS_3D_context) { + auto dbContext = DatabaseContext::create(); + auto authFactory = AuthorityFactory::create(dbContext, std::string()); + auto authFactoryEPSG = AuthorityFactory::create(dbContext, "EPSG"); + auto ctxt = CoordinateOperationContext::create(authFactory, nullptr, 0.0); + ctxt->setSpatialCriterion( + CoordinateOperationContext::SpatialCriterion::PARTIAL_INTERSECTION); + ctxt->setGridAvailabilityUse( + CoordinateOperationContext::GridAvailabilityUse:: + IGNORE_GRID_AVAILABILITY); + { + auto list = CoordinateOperationFactory::create()->createOperations( + authFactoryEPSG->createCoordinateReferenceSystem("4269") + ->promoteTo3D(std::string(), + dbContext), // NAD83 (86) promoted to 3D + authFactoryEPSG->createCoordinateReferenceSystem( + "6319"), // NAD83 (2011) 3D + ctxt); + ASSERT_GE(list.size(), 1U); + // Check we don't get an ESRI operation + EXPECT_EQ(list[0]->nameStr(), "NAD83 to NAD83(HARN) (47) + " + "NAD83(HARN) to NAD83(FBN) (1) + " + "NAD83(FBN) to NAD83(NSRS2007) (1) + " + "NAD83(NSRS2007) to NAD83(2011) (1)"); + } + { + auto list = CoordinateOperationFactory::create()->createOperations( + authFactoryEPSG->createCoordinateReferenceSystem( + "6319"), // NAD83 (2011) 3D + authFactoryEPSG->createCoordinateReferenceSystem("4269") + ->promoteTo3D(std::string(), + dbContext), // NAD83 (86) promoted to 3D + ctxt); + ASSERT_GE(list.size(), 1U); + // Check we don't get an ESRI operation + EXPECT_EQ(list[0]->nameStr(), + "Inverse of NAD83(NSRS2007) to NAD83(2011) (1) + " + "Inverse of NAD83(FBN) to NAD83(NSRS2007) (1) + " + "Inverse of NAD83(HARN) to NAD83(FBN) (1) + " + "Inverse of NAD83 to NAD83(HARN) (47)"); + } +} + +// --------------------------------------------------------------------------- + static GeodeticCRSNNPtr createGeocentricDatumWGS84() { PropertyMap propertiesCRS; propertiesCRS.set(Identifier::CODESPACE_KEY, "EPSG")