createFromCRSCodesWithIntermediates(): do not consider intermediate CRS whose datum has a publication date older than the source and target datums

2 files changed, 83 insertions, 16 deletions

diff --git a/src/iso19111/factory.cpp b/src/iso19111/factory.cpp
index 98b05317..987136cb 100644
--- a/src/iso19111/factory.cpp
+++ b/src/iso19111/factory.cpp
@@ -3886,11 +3886,45 @@ AuthorityFactory::createFromCRSCodesWithIntermediates(
         joinArea +
         "WHERE v1.source_crs_auth_name = ? AND v1.source_crs_code = ? "
         "AND v2.target_crs_auth_name = ? AND v2.target_crs_code = ? ");
+    std::string minDate;
+    std::string criterionOnIntermediateCRS;
     if (allowedIntermediateObjectType == ObjectType::GEOGRAPHIC_CRS) {
-        sql += "AND EXISTS(SELECT 1 FROM geodetic_crs x WHERE "
-               "x.auth_name = v1.target_crs_auth_name AND "
-               "x.code = v1.target_crs_code AND "
-               "x.type IN ('geographic 2D', 'geographic 3D')) ";
+        auto sourceCRS = d->createFactory(sourceCRSAuthName)
+                             ->createGeodeticCRS(sourceCRSCode);
+        auto targetCRS = d->createFactory(targetCRSAuthName)
+                             ->createGeodeticCRS(targetCRSCode);
+        const auto &sourceDatum = sourceCRS->datum();
+        const auto &targetDatum = targetCRS->datum();
+        if (sourceDatum && sourceDatum->publicationDate().has_value() &&
+            targetDatum && targetDatum->publicationDate().has_value()) {
+            const auto sourceDate(sourceDatum->publicationDate()->toString());
+            const auto targetDate(targetDatum->publicationDate()->toString());
+            minDate = std::min(sourceDate, targetDate);
+            // Check that the datum of the intermediateCRS has a publication
+            // date most recent that the one of the source and the target CRS
+            // Except when using the usual WGS84 pivot which happens to have a
+            // NULL publication date.
+            criterionOnIntermediateCRS =
+                "AND EXISTS(SELECT 1 FROM geodetic_crs x "
+                "JOIN geodetic_datum y "
+                "ON "
+                "y.auth_name = x.datum_auth_name AND "
+                "y.code = x.datum_code "
+                "WHERE "
+                "x.auth_name = v1.target_crs_auth_name AND "
+                "x.code = v1.target_crs_code AND "
+                "x.type IN ('geographic 2D', 'geographic 3D') AND "
+                "(y.publication_date IS NULL OR "
+                "(y.publication_date >= '" +
+                minDate + "'))) ";
+        } else {
+            criterionOnIntermediateCRS =
+                "AND EXISTS(SELECT 1 FROM geodetic_crs x WHERE "
+                "x.auth_name = v1.target_crs_auth_name AND "
+                "x.code = v1.target_crs_code AND "
+                "x.type IN ('geographic 2D', 'geographic 3D')) ";
+        }
+        sql += criterionOnIntermediateCRS;
     }
     auto params = ListOfParams{sourceCRSAuthName, sourceCRSCode,
                                targetCRSAuthName, targetCRSCode};
@@ -4049,10 +4083,7 @@ AuthorityFactory::createFromCRSCodesWithIntermediates(
           "WHERE v1.source_crs_auth_name = ? AND v1.source_crs_code = ? "
           "AND v2.source_crs_auth_name = ? AND v2.source_crs_code = ? ";
     if (allowedIntermediateObjectType == ObjectType::GEOGRAPHIC_CRS) {
-        sql += "AND EXISTS(SELECT 1 FROM geodetic_crs x WHERE "
-               "x.auth_name = v1.target_crs_auth_name AND "
-               "x.code = v1.target_crs_code AND "
-               "x.type IN ('geographic 2D', 'geographic 3D')) ";
+        sql += criterionOnIntermediateCRS;
     }
     intermediateWhere =
         buildIntermediateWhere(intermediateCRSAuthCodes, "target", "target");
@@ -4096,10 +4127,28 @@ AuthorityFactory::createFromCRSCodesWithIntermediates(
           "WHERE v1.target_crs_auth_name = ? AND v1.target_crs_code = ? "
           "AND v2.target_crs_auth_name = ? AND v2.target_crs_code = ? ";
     if (allowedIntermediateObjectType == ObjectType::GEOGRAPHIC_CRS) {
-        sql += "AND EXISTS(SELECT 1 FROM geodetic_crs x WHERE "
-               "x.auth_name = v1.source_crs_auth_name AND "
-               "x.code = v1.source_crs_code AND "
-               "x.type IN ('geographic 2D', 'geographic 3D')) ";
+        if (!minDate.empty()) {
+            criterionOnIntermediateCRS =
+                "AND EXISTS(SELECT 1 FROM geodetic_crs x "
+                "JOIN geodetic_datum y "
+                "ON "
+                "y.auth_name = x.datum_auth_name AND "
+                "y.code = x.datum_code "
+                "WHERE "
+                "x.auth_name = v1.source_crs_auth_name AND "
+                "x.code = v1.source_crs_code AND "
+                "x.type IN ('geographic 2D', 'geographic 3D') AND "
+                "(y.publication_date IS NULL OR "
+                "(y.publication_date >= '" +
+                minDate + "'))) ";
+        } else {
+            criterionOnIntermediateCRS =
+                "AND EXISTS(SELECT 1 FROM geodetic_crs x WHERE "
+                "x.auth_name = v1.source_crs_auth_name AND "
+                "x.code = v1.source_crs_code AND "
+                "x.type IN ('geographic 2D', 'geographic 3D')) ";
+        }
+        sql += criterionOnIntermediateCRS;
     }
     intermediateWhere =
         buildIntermediateWhere(intermediateCRSAuthCodes, "source", "source");
@@ -4143,10 +4192,7 @@ AuthorityFactory::createFromCRSCodesWithIntermediates(
           "WHERE v1.target_crs_auth_name = ? AND v1.target_crs_code = ? "
           "AND v2.source_crs_auth_name = ? AND v2.source_crs_code = ? ";
     if (allowedIntermediateObjectType == ObjectType::GEOGRAPHIC_CRS) {
-        sql += "AND EXISTS(SELECT 1 FROM geodetic_crs x WHERE "
-               "x.auth_name = v1.source_crs_auth_name AND "
-               "x.code = v1.source_crs_code AND "
-               "x.type IN ('geographic 2D', 'geographic 3D')) ";
+        sql += criterionOnIntermediateCRS;
     }
     intermediateWhere =
         buildIntermediateWhere(intermediateCRSAuthCodes, "source", "target");
diff --git a/test/unit/test_operation.cpp b/test/unit/test_operation.cpp
index 73f5dbdc..01ebf219 100644
--- a/test/unit/test_operation.cpp
+++ b/test/unit/test_operation.cpp
@@ -4783,6 +4783,27 @@ TEST(operation, geogCRS_to_geogCRS_context_concatenated_operation) {
 
 // ---------------------------------------------------------------------------
 
+TEST(operation, geogCRS_to_geogCRS_context_ED50_to_WGS72_no_NTF_intermediate) {
+    auto authFactory =
+        AuthorityFactory::create(DatabaseContext::create(), "EPSG");
+    auto ctxt = CoordinateOperationContext::create(authFactory, nullptr, 0.0);
+    ctxt->setSpatialCriterion(
+        CoordinateOperationContext::SpatialCriterion::PARTIAL_INTERSECTION);
+    auto list = CoordinateOperationFactory::create()->createOperations(
+        authFactory->createCoordinateReferenceSystem("4230"), // ED50
+        authFactory->createCoordinateReferenceSystem("4322"), // WGS 72
+        ctxt);
+    ASSERT_GE(list.size(), 2U);
+    // We should not use the ancient NTF as an intermediate when looking for
+    // ED50 -> WGS 72 operations.
+    for (const auto &op : list) {
+        EXPECT_TRUE(op->nameStr().find("NTF") == std::string::npos)
+            << op->nameStr();
+    }
+}
+
+// ---------------------------------------------------------------------------
+
 TEST(operation, geogCRS_to_geogCRS_context_same_grid_name) {
     auto authFactory =
         AuthorityFactory::create(DatabaseContext::create(), "EPSG");