Skip to content

Commit 7247f06

Browse files
committed
fix: chain engineering CRS transforms through intermediate projected CRS
When an EngineeringCRS has registered transformations to projected CRSs (e.g. via similarity transformation), PROJ could not compose a pipeline to a *different* projected CRS on the same datum. For example, if engineering CRS A has a transformation to UTM zone 31N, PROJ returned 0 operations when asked to transform A to UTM zone 32N. Add findOpsInRegistryDirectFrom() and use it in createOperationsFromDatabase() to look up all registered transformations from/to the engineering CRS and chain through the intermediate CRS to reach the target. Tested with EPSG:5800 (Astra Minas Grid) -> EPSG:22193 (Campo Inchauspe / Argentina 3) which chains through EPSG:22192 (Argentina 2).
1 parent fc4c131 commit 7247f06

File tree

2 files changed

+183
-0
lines changed

2 files changed

+183
-0
lines changed

src/iso19111/operation/coordinateoperationfactory.cpp

Lines changed: 134 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -608,6 +608,10 @@ struct CoordinateOperationFactory::Private {
608608
findOpsInRegistryDirectTo(const crs::CRSNNPtr &targetCRS,
609609
Private::Context &context);
610610

611+
static std::vector<CoordinateOperationNNPtr>
612+
findOpsInRegistryDirectFrom(const crs::CRSNNPtr &sourceCRS,
613+
Private::Context &context);
614+
611615
static std::vector<CoordinateOperationNNPtr>
612616
findsOpsInRegistryWithIntermediate(
613617
const crs::CRSNNPtr &sourceCRS, const crs::CRSNNPtr &targetCRS,
@@ -1995,6 +1999,71 @@ CoordinateOperationFactory::Private::findOpsInRegistryDirectTo(
19951999
return std::vector<CoordinateOperationNNPtr>();
19962000
}
19972001

2002+
// ---------------------------------------------------------------------------
2003+
2004+
// Look in the authority registry for operations from sourceCRS
2005+
std::vector<CoordinateOperationNNPtr>
2006+
CoordinateOperationFactory::Private::findOpsInRegistryDirectFrom(
2007+
const crs::CRSNNPtr &sourceCRS, Private::Context &context) {
2008+
#ifdef TRACE_CREATE_OPERATIONS
2009+
ENTER_BLOCK("findOpsInRegistryDirectFrom(" +
2010+
objectAsStr(sourceCRS.get()) + " --> {any})");
2011+
#endif
2012+
2013+
const auto &authFactory = context.context->getAuthorityFactory();
2014+
assert(authFactory);
2015+
2016+
std::list<std::pair<std::string, std::string>> ids;
2017+
buildCRSIds(sourceCRS, context, ids);
2018+
2019+
const auto gridAvailabilityUse = context.context->getGridAvailabilityUse();
2020+
for (const auto &id : ids) {
2021+
const auto &srcAuthName = id.first;
2022+
const auto &srcCode = id.second;
2023+
2024+
const auto authorities(getCandidateAuthorities(
2025+
authFactory, srcAuthName, srcAuthName));
2026+
std::vector<CoordinateOperationNNPtr> res;
2027+
for (const auto &authName : authorities) {
2028+
const auto tmpAuthFactory = io::AuthorityFactory::create(
2029+
authFactory->databaseContext(), authName);
2030+
auto resTmp =
2031+
tmpAuthFactory->createFromCoordinateReferenceSystemCodes(
2032+
srcAuthName, srcCode, std::string(), std::string(),
2033+
context.context->getUsePROJAlternativeGridNames(),
2034+
2035+
gridAvailabilityUse ==
2036+
CoordinateOperationContext::GridAvailabilityUse::
2037+
DISCARD_OPERATION_IF_MISSING_GRID ||
2038+
gridAvailabilityUse ==
2039+
CoordinateOperationContext::GridAvailabilityUse::
2040+
KNOWN_AVAILABLE,
2041+
gridAvailabilityUse ==
2042+
CoordinateOperationContext::GridAvailabilityUse::
2043+
KNOWN_AVAILABLE,
2044+
context.context->getDiscardSuperseded(), true, true,
2045+
context.extent1, context.extent2);
2046+
res.insert(res.end(), resTmp.begin(), resTmp.end());
2047+
if (authName == "PROJ") {
2048+
continue;
2049+
}
2050+
if (!res.empty()) {
2051+
auto resFiltered =
2052+
FilterResults(res, context.context, context.extent1,
2053+
context.extent2, false)
2054+
.getRes();
2055+
#ifdef TRACE_CREATE_OPERATIONS
2056+
logTrace("filtering reduced from " +
2057+
toString(static_cast<int>(res.size())) + " to " +
2058+
toString(static_cast<int>(resFiltered.size())));
2059+
#endif
2060+
return resFiltered;
2061+
}
2062+
}
2063+
}
2064+
return std::vector<CoordinateOperationNNPtr>();
2065+
}
2066+
19982067
//! @endcond
19992068

20002069
// ---------------------------------------------------------------------------
@@ -3940,6 +4009,71 @@ bool CoordinateOperationFactory::Private::createOperationsFromDatabase(
39404009
return true;
39414010
}
39424011

4012+
// Special case for EngineeringCRS: if no direct transformation found,
4013+
// look up all registered transformations from the engineering CRS and
4014+
// try to chain through them to reach the target CRS.
4015+
// This handles cases like engineering CRS → ProjectedCRS_A (registered)
4016+
// chained with ProjectedCRS_A → ProjectedCRS_B (same datum reprojection).
4017+
if (res.empty()) {
4018+
auto engSrc =
4019+
dynamic_cast<const crs::EngineeringCRS *>(sourceCRS.get());
4020+
auto engDst =
4021+
dynamic_cast<const crs::EngineeringCRS *>(targetCRS.get());
4022+
if (engSrc && !engDst) {
4023+
auto opsFromEng = findOpsInRegistryDirectFrom(sourceCRS, context);
4024+
for (const auto &opFromEng : opsFromEng) {
4025+
auto intermediateCRS = opFromEng->targetCRS();
4026+
if (intermediateCRS &&
4027+
!intermediateCRS->_isEquivalentTo(
4028+
targetCRS.get(),
4029+
util::IComparable::Criterion::EQUIVALENT)) {
4030+
auto opsToTarget = createOperations(
4031+
NN_NO_CHECK(intermediateCRS), sourceEpoch, targetCRS,
4032+
targetEpoch, context);
4033+
for (const auto &opToTarget : opsToTarget) {
4034+
try {
4035+
res.emplace_back(
4036+
ConcatenatedOperation::createComputeMetadata(
4037+
{opFromEng, opToTarget},
4038+
context.disallowEmptyIntersection()));
4039+
} catch (
4040+
const InvalidOperationEmptyIntersection &) {
4041+
}
4042+
}
4043+
}
4044+
}
4045+
if (!res.empty()) {
4046+
return true;
4047+
}
4048+
} else if (engDst && !engSrc) {
4049+
auto opsToEng = findOpsInRegistryDirectTo(targetCRS, context);
4050+
for (const auto &opToEng : opsToEng) {
4051+
auto intermediateCRS = opToEng->sourceCRS();
4052+
if (intermediateCRS &&
4053+
!intermediateCRS->_isEquivalentTo(
4054+
sourceCRS.get(),
4055+
util::IComparable::Criterion::EQUIVALENT)) {
4056+
auto opsFromSource = createOperations(
4057+
sourceCRS, sourceEpoch,
4058+
NN_NO_CHECK(intermediateCRS), targetEpoch, context);
4059+
for (const auto &opFromSource : opsFromSource) {
4060+
try {
4061+
res.emplace_back(
4062+
ConcatenatedOperation::createComputeMetadata(
4063+
{opFromSource, opToEng},
4064+
context.disallowEmptyIntersection()));
4065+
} catch (
4066+
const InvalidOperationEmptyIntersection &) {
4067+
}
4068+
}
4069+
}
4070+
}
4071+
if (!res.empty()) {
4072+
return true;
4073+
}
4074+
}
4075+
}
4076+
39434077
bool doFilterAndCheckPerfectOp = false;
39444078

39454079
bool sameGeodeticDatum = false;

test/unit/test_operationfactory.cpp

Lines changed: 49 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -12142,3 +12142,52 @@ TEST(operation, createOperation_defmodel_from_database) {
1214212142
"+step +proj=unitconvert +xy_in=rad +z_in=m +xy_out=deg +z_out=m "
1214312143
"+step +proj=axisswap +order=2,1");
1214412144
}
12145+
12146+
// ---------------------------------------------------------------------------
12147+
12148+
TEST(
12149+
operation,
12150+
engineeringCRS_to_projectedCRS_through_intermediate_same_datum) {
12151+
// Test that PROJ can compose a pipeline from an engineering CRS to a
12152+
// projected CRS via chaining through an intermediate projected CRS that
12153+
// has a registered transformation from the engineering CRS.
12154+
//
12155+
// EPSG:5800 (Astra Minas Grid) has a registered transformation
12156+
// (EPSG:1035) to EPSG:22192 (Campo Inchauspe / Argentina 2).
12157+
// EPSG:22193 (Campo Inchauspe / Argentina 3) shares the same datum
12158+
// (Campo Inchauspe) but a different zone.
12159+
// PROJ should chain: engineering → Argentina 2 → inv(zone 2) → zone 3.
12160+
auto authFactory =
12161+
AuthorityFactory::create(DatabaseContext::create(), "EPSG");
12162+
12163+
// Direct transformation should work (EPSG:1035)
12164+
{
12165+
auto ctxt =
12166+
CoordinateOperationContext::create(authFactory, nullptr, 0.0);
12167+
ctxt->setSourceAndTargetCRSExtentUse(
12168+
CoordinateOperationContext::SourceTargetCRSExtentUse::NONE);
12169+
auto list = CoordinateOperationFactory::create()->createOperations(
12170+
authFactory->createCoordinateReferenceSystem("5800"),
12171+
authFactory->createCoordinateReferenceSystem("22192"), ctxt);
12172+
ASSERT_GE(list.size(), 1U);
12173+
EXPECT_EQ(list[0]->nameStr(),
12174+
"Astra Minas to Campo Inchauspe / Argentina 2 (1)");
12175+
}
12176+
12177+
// Chained transformation to a different zone on the same datum
12178+
{
12179+
auto ctxt =
12180+
CoordinateOperationContext::create(authFactory, nullptr, 0.0);
12181+
ctxt->setSourceAndTargetCRSExtentUse(
12182+
CoordinateOperationContext::SourceTargetCRSExtentUse::NONE);
12183+
auto list = CoordinateOperationFactory::create()->createOperations(
12184+
authFactory->createCoordinateReferenceSystem("5800"),
12185+
authFactory->createCoordinateReferenceSystem("22193"), ctxt);
12186+
ASSERT_GE(list.size(), 1U);
12187+
// Should chain through the registered transformation to Argentina 2,
12188+
// then reproject from zone 2 to zone 3 (same datum)
12189+
EXPECT_EQ(list[0]->nameStr(),
12190+
"Astra Minas to Campo Inchauspe / Argentina 2 (1) + "
12191+
"Inverse of Argentina zone 2 + Argentina zone 3");
12192+
}
12193+
}

0 commit comments

Comments
 (0)