Database: register grids for New Caledonia (added per https://github.com/OSGeo/PROJ-data/pull/99) and tune createOperations()

data/sql/grid_alternatives.sql

@@ -207,9 +207,11 @@ VALUES
 
 -- nc_dittt - Gouvernement de Nouvelle Calédonie - DITTT
 ('Ranc08_Circe.mnt','nc_dittt_Ranc08_Circe.tif',NULL,'GTiff','geoid_like',0,NULL,'https://cdn.proj.org/nc_dittt_Ranc08_Circe.tif',1,1,NULL),
+('RANC15.tac','nc_dittt_RANC15.tif',NULL,'GTiff','geoid_like',0,NULL,'https://cdn.proj.org/nc_dittt_RANC15.tif',1,1,NULL),
 ('gr3dnc01b.mnt','nc_dittt_gr3dnc01b.tif',NULL,'GTiff','geocentricoffset',0,NULL,'https://cdn.proj.org/nc_dittt_gr3dnc01b.tif',1,1,NULL),
 ('gr3dnc02b.mnt','nc_dittt_gr3dnc02b.tif',NULL,'GTiff','geocentricoffset',0,NULL,'https://cdn.proj.org/nc_dittt_gr3dnc02b.tif',1,1,NULL),
 ('gr3dnc03a.mnt','nc_dittt_gr3dnc03a.tif',NULL,'GTiff','geocentricoffset',0,NULL,'https://cdn.proj.org/nc_dittt_gr3dnc03a.tif',1,1,NULL),
+('gr3dncl08.tac','nc_dittt_gr3dncI08.tif',NULL,'GTiff','geocentricoffset',0,NULL,'https://cdn.proj.org/nc_dittt_gr3dncI08.tif',1,1,NULL),
 
 -- Netherlands / RDNAP (non-free grids). See https://salsa.debian.org/debian-gis-team/proj-rdnap/raw/master/debian/copyright
 -- Netherlands / RDNAP 2018

src/iso19111/operation/coordinateoperationfactory.cpp

@@ -2978,9 +2978,14 @@ void CoordinateOperationFactory::Private::createOperationsWithDatumPivot(
     auto createTransformations = [&](const crs::CRSNNPtr &candidateSrcGeod,
                                      const crs::CRSNNPtr &candidateDstGeod,
                                      const CoordinateOperationNNPtr &opFirst,
-                                     bool isNullFirst) {
+                                     bool isNullFirst,
+                                     bool useOnlyDirectRegistryOp) {
+        bool resNonEmptyBeforeFiltering;
         const auto opsSecond =
-            createOperations(candidateSrcGeod, candidateDstGeod, context);
+            useOnlyDirectRegistryOp
+                ? findOpsInRegistryDirect(candidateSrcGeod, candidateDstGeod,
+                                          context, resNonEmptyBeforeFiltering)
+                : createOperations(candidateSrcGeod, candidateDstGeod, context);
         const auto opsThird = createOperations(
             sourceAndTargetAre3D
                 ? candidateDstGeod->promoteTo3D(std::string(), dbContext)
@@ -3099,9 +3104,13 @@ void CoordinateOperationFactory::Private::createOperationsWithDatumPivot(
         }
     };
 
+    // The below logic is thus quite fragile, and attempts at changing it
+    // result in degraded results for other use cases...
+    //
     // Start in priority with candidates that have exactly the same name as
-    // the sourcCRS and targetCRS. Typically for the case of init=IGNF:XXXX
-
+    // the sourceCRS and targetCRS (Typically for the case of init=IGNF:XXXX);
+    // and then attempt candidate geodetic CRS with different names
+    //
     // Transformation from IGNF:NTFP to IGNF:RGF93G,
     // using
     // NTF geographiques Paris (gr) vers NTF GEOGRAPHIQUES GREENWICH (DMS) +
@@ -3110,72 +3119,98 @@ void CoordinateOperationFactory::Private::createOperationsWithDatumPivot(
     // of IGNF:RGF93G being returned before IGNF:RGF93GEO in candidatesDstGeod.
     // If RGF93GEO is returned before then we go through WGS84 and use
     // instead a Helmert transformation.
-    // The below logic is thus quite fragile, and attempts at changing it
-    // result in degraded results for other use cases...
-
-    for (const auto &candidateSrcGeod : candidatesSrcGeod) {
-        if (candidateSrcGeod->nameStr() == sourceCRS->nameStr()) {
-            auto sourceSrcGeodModified(
-                sourceAndTargetAre3D
-                    ? candidateSrcGeod->promoteTo3D(std::string(), dbContext)
-                    : candidateSrcGeod);
-            for (const auto &candidateDstGeod : candidatesDstGeod) {
-                if (candidateDstGeod->nameStr() == targetCRS->nameStr()) {
+    //
+    // Actually, in the general case, we do the lookup in 2 passes with the 2
+    // above steps in each pass:
+    // - one first pass where we only consider direct transformations (no
+    //   other intermediate CRS)
+    // - a second pass where we allow transformation through another
+    //   intermediate CRS.
+    // ... but when transforming between 2 IGNF CRS, we do just one single pass
+    // by allowing directly all transformation. There is no strong reason for
+    // that particular case, except that otherwise we'd get different results
+    // for thest test/cli/testIGNF script when transforming a point outside
+    // the area of validity... Not totally sure the behaviour we try to preserve
+    // here with the particular case is fundamentally better than the general
+    // case. The general case is needed typically for the RGNC91-93 -> RGNC15
+    // transformation where we we need to actually use a transformation between
+    // RGNC91-93 (lon-lat) -> RGNC15 (lon-lat), and not chain 2 no-op
+    // transformations RGNC91-93 -> WGS 84 and RGNC15 -> WGS84.
+
+    const auto isIGNF = [](const crs::CRSNNPtr &crs) {
+        const auto &ids = crs->identifiers();
+        return !ids.empty() && *(ids.front()->codeSpace()) == "IGNF";
+    };
+    const int nIters = (isIGNF(sourceCRS) && isIGNF(targetCRS)) ? 1 : 2;
+    for (int iter = 0; iter < nIters; ++iter) {
+        const bool useOnlyDirectRegistryOp = (iter == 0 && nIters == 2);
+        for (const auto &candidateSrcGeod : candidatesSrcGeod) {
+            if (candidateSrcGeod->nameStr() == sourceCRS->nameStr()) {
+                auto sourceSrcGeodModified(sourceAndTargetAre3D
+                                               ? candidateSrcGeod->promoteTo3D(
+                                                     std::string(), dbContext)
+                                               : candidateSrcGeod);
+                for (const auto &candidateDstGeod : candidatesDstGeod) {
+                    if (candidateDstGeod->nameStr() == targetCRS->nameStr()) {
 #ifdef TRACE_CREATE_OPERATIONS
-                    ENTER_BLOCK("try " + objectAsStr(sourceCRS.get()) + "->" +
-                                objectAsStr(candidateSrcGeod.get()) + "->" +
-                                objectAsStr(candidateDstGeod.get()) + "->" +
-                                objectAsStr(targetCRS.get()) + ")");
+                        ENTER_BLOCK("try " + objectAsStr(sourceCRS.get()) +
+                                    "->" + objectAsStr(candidateSrcGeod.get()) +
+                                    "->" + objectAsStr(candidateDstGeod.get()) +
+                                    "->" + objectAsStr(targetCRS.get()) + ")");
 #endif
-                    const auto opsFirst = createOperations(
-                        sourceCRS, sourceSrcGeodModified, context);
-                    assert(!opsFirst.empty());
-                    const bool isNullFirst =
-                        isNullTransformation(opsFirst[0]->nameStr());
-                    createTransformations(candidateSrcGeod, candidateDstGeod,
-                                          opsFirst[0], isNullFirst);
-                    if (!res.empty()) {
-                        if (hasResultSetOnlyResultsWithPROJStep(res)) {
-                            continue;
+                        const auto opsFirst = createOperations(
+                            sourceCRS, sourceSrcGeodModified, context);
+                        assert(!opsFirst.empty());
+                        const bool isNullFirst =
+                            isNullTransformation(opsFirst[0]->nameStr());
+                        createTransformations(
+                            candidateSrcGeod, candidateDstGeod, opsFirst[0],
+                            isNullFirst, useOnlyDirectRegistryOp);
+                        if (!res.empty()) {
+                            if (hasResultSetOnlyResultsWithPROJStep(res)) {
+                                continue;
+                            }
+                            return;
                         }
-                        return;
                     }
                 }
             }
         }
-    }
 
-    for (const auto &candidateSrcGeod : candidatesSrcGeod) {
-        const bool bSameSrcName =
-            candidateSrcGeod->nameStr() == sourceCRS->nameStr();
+        for (const auto &candidateSrcGeod : candidatesSrcGeod) {
+            const bool bSameSrcName =
+                candidateSrcGeod->nameStr() == sourceCRS->nameStr();
 #ifdef TRACE_CREATE_OPERATIONS
-        ENTER_BLOCK("");
+            ENTER_BLOCK("");
 #endif
-        auto sourceSrcGeodModified(
-            sourceAndTargetAre3D
-                ? candidateSrcGeod->promoteTo3D(std::string(), dbContext)
-                : candidateSrcGeod);
-        const auto opsFirst =
-            createOperations(sourceCRS, sourceSrcGeodModified, context);
-        assert(!opsFirst.empty());
-        const bool isNullFirst = isNullTransformation(opsFirst[0]->nameStr());
-
-        for (const auto &candidateDstGeod : candidatesDstGeod) {
-            if (bSameSrcName &&
-                candidateDstGeod->nameStr() == targetCRS->nameStr()) {
-                continue;
-            }
+            auto sourceSrcGeodModified(
+                sourceAndTargetAre3D
+                    ? candidateSrcGeod->promoteTo3D(std::string(), dbContext)
+                    : candidateSrcGeod);
+            const auto opsFirst =
+                createOperations(sourceCRS, sourceSrcGeodModified, context);
+            assert(!opsFirst.empty());
+            const bool isNullFirst =
+                isNullTransformation(opsFirst[0]->nameStr());
+
+            for (const auto &candidateDstGeod : candidatesDstGeod) {
+                if (bSameSrcName &&
+                    candidateDstGeod->nameStr() == targetCRS->nameStr()) {
+                    continue;
+                }
 
 #ifdef TRACE_CREATE_OPERATIONS
-            ENTER_BLOCK("try " + objectAsStr(sourceCRS.get()) + "->" +
-                        objectAsStr(candidateSrcGeod.get()) + "->" +
-                        objectAsStr(candidateDstGeod.get()) + "->" +
-                        objectAsStr(targetCRS.get()) + ")");
+                ENTER_BLOCK("try " + objectAsStr(sourceCRS.get()) + "->" +
+                            objectAsStr(candidateSrcGeod.get()) + "->" +
+                            objectAsStr(candidateDstGeod.get()) + "->" +
+                            objectAsStr(targetCRS.get()) + ")");
 #endif
-            createTransformations(candidateSrcGeod, candidateDstGeod,
-                                  opsFirst[0], isNullFirst);
-            if (!res.empty() && !hasResultSetOnlyResultsWithPROJStep(res)) {
-                return;
+                createTransformations(candidateSrcGeod, candidateDstGeod,
+                                      opsFirst[0], isNullFirst,
+                                      useOnlyDirectRegistryOp);
+                if (!res.empty() && !hasResultSetOnlyResultsWithPROJStep(res)) {
+                    return;
+                }
             }
         }
     }

test/unit/test_operationfactory.cpp

@@ -1046,6 +1046,31 @@ TEST(operation, geogCRS_to_geogCRS_longitude_rotation_context) {
 
 // ---------------------------------------------------------------------------
 
+TEST(operation, geogCRS_to_geogCRS_context_lonlat_vs_latlon_crs) {
+    auto authFactory =
+        AuthorityFactory::create(DatabaseContext::create(), "EPSG");
+    auto ctxt = CoordinateOperationContext::create(authFactory, nullptr, 0.0);
+    ctxt->setGridAvailabilityUse(
+        CoordinateOperationContext::GridAvailabilityUse::
+            IGNORE_GRID_AVAILABILITY);
+    ctxt->setSpatialCriterion(
+        CoordinateOperationContext::SpatialCriterion::PARTIAL_INTERSECTION);
+    auto list = CoordinateOperationFactory::create()->createOperations(
+        authFactory->createCoordinateReferenceSystem("4749"),  // RGNC91-93
+        authFactory->createCoordinateReferenceSystem("10310"), // RGNC15
+        ctxt);
+    ASSERT_EQ(list.size(), 3U);
+    // Check that we get direct transformation, and not through WGS 84
+    // The difficulty here is that the transformation is registered between
+    // "RGNC91-93 (lon-lat)" et "RGNC15 (lon-lat)"
+    EXPECT_EQ(list[0]->nameStr(), "axis order change (2D) + RGNC91-93 to "
+                                  "RGNC15 (2) + axis order change (2D)");
+    EXPECT_EQ(list[1]->nameStr(), "axis order change (2D) + RGNC91-93 to "
+                                  "RGNC15 (1) + axis order change (2D)");
+}
+
+// ---------------------------------------------------------------------------
+
 TEST(operation, geogCRS_to_geogCRS_context_concatenated_operation) {
     auto authFactory =
         AuthorityFactory::create(DatabaseContext::create(), "EPSG");