Merge pull request #1231 from rouault/fix_1229

proj_create_crs_to_crs(): defer selection of actual coordinate operation until proj_trans() is called (fixes #1229)

3 files changed, 66 insertions, 0 deletions

diff --git a/test/cli/td_out.dist b/test/cli/td_out.dist
index f6b2a219..cdbe0fa9 100644
--- a/test/cli/td_out.dist
+++ b/test/cli/td_out.dist
@@ -23,3 +23,7 @@ edge or even a wee bit outside (#141).
 -5.5000000000001 52.0000000000001	-5.498893534472	52.000109529717 0.000000000000
 -5.4999 51.9999	-5.498793541695	52.000009529743 0.000000000000
 -5.5001 52.0	-5.500100000000	52.000000000000 0.000000000000
+##############################################################
+NAD27 -> NAD83: 1st through ntv1, 2nd through conus
+44d00'00.000"N 111d00'00.000"W 0.0	43d59'59.732"N	111d0'3.208"W 0.000
+39d00'00.000"N 111d00'00.000"W 0.0	38d59'59.912"N	111d0'2.604"W 0.000
diff --git a/test/cli/testdatumfile b/test/cli/testdatumfile
index d048d8e6..e8995150 100755
--- a/test/cli/testdatumfile
+++ b/test/cli/testdatumfile
@@ -97,7 +97,16 @@ $EXE +proj=latlong +datum=WGS84 \
 -5.4999 51.9999
 -5.5001 52.0
 EOF
+#
+echo "##############################################################" >> ${OUT}
+echo "NAD27 -> NAD83: 1st through ntv1, 2nd through conus" >> ${OUT}
+#
+$EXE NAD27 NAD83 -E >>${OUT} <<EOF
+44d00'00.000"N 111d00'00.000"W 0.0
+39d00'00.000"N 111d00'00.000"W 0.0
+EOF
 
+# Cleanup
 rm -rf "dir with \" space"
 
 #
diff --git a/test/unit/gie_self_tests.cpp b/test/unit/gie_self_tests.cpp
index 1b5ce83b..2d298ea3 100644
--- a/test/unit/gie_self_tests.cpp
+++ b/test/unit/gie_self_tests.cpp
@@ -659,4 +659,57 @@ TEST(gie, horner_selftest) {
     proj_destroy(P);
 }
 
+// ---------------------------------------------------------------------------
+
+TEST(gie, proj_create_crs_to_crs_PULKOVO42_ETRS89) {
+    auto P = proj_create_crs_to_crs(PJ_DEFAULT_CTX, "EPSG:4179", "EPSG:4258",
+                                    nullptr);
+    ASSERT_TRUE(P != nullptr);
+    PJ_COORD c;
+
+    EXPECT_EQ(std::string(proj_pj_info(P).definition),
+              "unavailable until proj_trans is called");
+    EXPECT_EQ(proj_get_name(P), nullptr);
+    EXPECT_EQ(P->fwd, nullptr);
+    EXPECT_EQ(P->fwd3d, nullptr);
+    EXPECT_EQ(P->fwd4d, nullptr);
+
+    // Romania
+    c.xyz.x = 45; // Lat
+    c.xyz.y = 25; // Long
+    c.xyz.z = 0;
+    c = proj_trans(P, PJ_FWD, c);
+    EXPECT_NEAR(c.xy.x, 44.999701238, 1e-9);
+    EXPECT_NEAR(c.xy.y, 24.998474948, 1e-9);
+    EXPECT_EQ(std::string(proj_pj_info(P).definition),
+              "proj=pipeline step proj=axisswap order=2,1 step "
+              "proj=unitconvert xy_in=deg xy_out=rad step proj=cart "
+              "ellps=krass step proj=helmert x=2.3287 y=-147.0425 z=-92.0802 "
+              "rx=0.3092483 ry=-0.32482185 rz=-0.49729934 s=5.68906266 "
+              "convention=coordinate_frame step inv proj=cart ellps=GRS80 step "
+              "proj=unitconvert xy_in=rad xy_out=deg step proj=axisswap "
+              "order=2,1");
+
+    c = proj_trans(P, PJ_INV, c);
+    EXPECT_NEAR(c.xy.x, 45, 1e-9);
+    EXPECT_NEAR(c.xy.y, 25, 1e-9);
+
+    // Poland
+    c.xyz.x = 52; // Lat
+    c.xyz.y = 20; // Long
+    c.xyz.z = 0;
+    c = proj_trans(P, PJ_FWD, c);
+    EXPECT_NEAR(c.xy.x, 51.999714150, 1e-9);
+    EXPECT_NEAR(c.xy.y, 19.998187811, 1e-9);
+    EXPECT_EQ(std::string(proj_pj_info(P).definition),
+              "proj=pipeline step proj=axisswap order=2,1 step "
+              "proj=unitconvert xy_in=deg xy_out=rad step proj=cart "
+              "ellps=krass step proj=helmert x=33.4 y=-146.6 z=-76.3 rx=-0.359 "
+              "ry=-0.053 rz=0.844 s=-0.84 convention=position_vector step inv "
+              "proj=cart ellps=GRS80 step proj=unitconvert xy_in=rad "
+              "xy_out=deg step proj=axisswap order=2,1");
+
+    proj_destroy(P);
+}
+
 } // namespace