Skip to content

Commit

Permalink
findsOpsInRegistryWithIntermediate(): restrict to using known source/…
Browse files Browse the repository at this point in the history
…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 #4244
(one could argue which of the results is the best, given that NAD83(86)
as a 3D geographic CRS has no solid foundation)
  • Loading branch information
rouault committed Sep 17, 2024
1 parent ed0a29c commit 78d5d84
Show file tree
Hide file tree
Showing 4 changed files with 98 additions and 2 deletions.
3 changes: 3 additions & 0 deletions include/proj/crs.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -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:
Expand Down
40 changes: 39 additions & 1 deletion src/iso19111/crs.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down Expand Up @@ -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.
*
Expand Down Expand Up @@ -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();
Expand Down
12 changes: 11 additions & 1 deletion src/iso19111/operation/coordinateoperationfactory.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -2001,12 +2001,21 @@ CoordinateOperationFactory::Private::findsOpsInRegistryWithIntermediate(
auto geodSrc = dynamic_cast<crs::GeodeticCRS *>(sourceCRS.get());
auto geodDst = dynamic_cast<crs::GeodeticCRS *>(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<CoordinateOperationNNPtr> 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<crs::GeographicCRS *>(sourceCRS.get()) !=
nullptr) ==
Expand All @@ -2022,6 +2031,7 @@ CoordinateOperationFactory::Private::findsOpsInRegistryWithIntermediate(
continue;
}
}

const auto opsWithIntermediate =
findsOpsInRegistryWithIntermediate(
candidateSrcGeod, targetCRS, context,
Expand Down
45 changes: 45 additions & 0 deletions test/unit/test_operationfactory.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -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")
Expand Down

0 comments on commit 78d5d84

Please sign in to comment.