Add new option to proj_create_crs_to_crs_from_pj method to force +over on transformation operations (#2914)

Fixes #2512

1 files changed, 184 insertions, 0 deletions

diff --git a/test/unit/gie_self_tests.cpp b/test/unit/gie_self_tests.cpp
index 5f96a93a..23684e68 100644
--- a/test/unit/gie_self_tests.cpp
+++ b/test/unit/gie_self_tests.cpp
@@ -890,4 +890,188 @@ TEST(gie, proj_trans_generic) {
     proj_destroy(P);
 }
 
+TEST(gie, proj_create_crs_to_crs_from_pj_force_over) {
+
+    PJ_CONTEXT* ctx;
+
+    ctx = proj_context_create();
+    ASSERT_TRUE(ctx != nullptr);
+
+    auto epsg27700 = proj_create(ctx, "EPSG:27700");
+    ASSERT_TRUE(epsg27700 != nullptr);
+
+    auto epsg4326 = proj_create(ctx, "EPSG:4326");
+    ASSERT_TRUE(epsg4326 != nullptr);
+    
+    auto epsg3857 = proj_create(ctx, "EPSG:3857");
+    ASSERT_TRUE(epsg3857 != nullptr);
+
+    {
+        const char* const options[] = { "FORCE_OVER=YES", nullptr };
+        auto P = proj_create_crs_to_crs_from_pj(ctx, epsg4326, epsg3857, nullptr, options);
+        ASSERT_TRUE(P != nullptr);
+        ASSERT_TRUE(P->over);
+        PJ_COORD input;
+        PJ_COORD input_over;
+
+        //Test a point along the equator. 
+        //The same point, but in two different representations.
+        input.xyz.x = 0; // Lat in deg
+        input.xyz.y = 140;  // Long in deg
+        input.xyz.z = 0;
+
+        input_over.xyz.x = 0; // Lat in deg
+        input_over.xyz.y = -220;  // Long in deg
+        input_over.xyz.z = 0;
+
+        auto output = proj_trans(P, PJ_FWD, input);
+        auto output_over = proj_trans(P, PJ_FWD, input_over);
+
+        auto input_inv = proj_trans(P, PJ_INV, output);
+        auto input_over_inv = proj_trans(P, PJ_INV, output_over);
+
+        //Web Mercator x's between 0 and 180 longitude come out positive.
+        //But when forcing the over flag, the -220 calculation makes it flip.
+        EXPECT_GT(output.xyz.x, 0);
+        EXPECT_LT(output_over.xyz.x, 0);
+
+        EXPECT_NEAR(output.xyz.x, 15584728.711058298, 1e-8);
+        EXPECT_NEAR(output_over.xyz.x, -24490287.974520184, 1e-8);
+
+        //The distance from 140 to 180 and -220 to -180 should be pretty much the same.
+        auto dx_o = fabs(output.xyz.x - 20037508.342789244);
+        auto dx_over = fabs(output_over.xyz.x + 20037508.342789244);
+        auto dx = fabs(dx_o - dx_over);
+
+        EXPECT_NEAR(dx, 0, 1e-8);
+
+        //Check the inverse operations get us back close to our original input values.
+        EXPECT_NEAR(input.xyz.x, input_inv.xyz.x, 1e-8);
+        EXPECT_NEAR(input.xyz.y, input_inv.xyz.y, 1e-8);
+        EXPECT_NEAR(input_over.xyz.x, input_over_inv.xyz.x, 1e-8);
+        EXPECT_NEAR(input_over.xyz.y, input_over_inv.xyz.y, 1e-8);
+
+        proj_destroy(P);
+    }
+
+    {
+        //Try again with force over set to anything but YES to verify it didn't do anything.
+        const char* const options[] = { "FORCE_OVER=NO", nullptr };
+        auto P = proj_create_crs_to_crs_from_pj(ctx, epsg4326, epsg3857, nullptr, options);
+        ASSERT_TRUE(P != nullptr);
+        ASSERT_FALSE(P->over);
+        PJ_COORD input;
+        PJ_COORD input_notOver;
+
+        input.xyz.x = 0; // Lat in deg
+        input.xyz.y = 140;  // Long in deg
+        input.xyz.z = 0;
+
+        input_notOver.xyz.x = 0; // Lat in deg
+        input_notOver.xyz.y = -220;  // Long in deg
+        input_notOver.xyz.z = 0;
+
+        auto output = proj_trans(P, PJ_FWD, input);
+        auto output_notOver = proj_trans(P, PJ_FWD, input_notOver);
+
+        EXPECT_GT(output.xyz.x, 0);
+        EXPECT_GT(output_notOver.xyz.x, 0);
+
+        EXPECT_NEAR(output.xyz.x, 15584728.711058298, 1e-8);
+        EXPECT_NEAR(output_notOver.xyz.x, 15584728.711058298, 1e-8);
+
+        proj_destroy(P);
+    }
+
+    {
+        //Try again with no options to verify it didn't do anything.
+        auto P = proj_create_crs_to_crs_from_pj(ctx, epsg4326, epsg3857, nullptr, nullptr);
+        ASSERT_TRUE(P != nullptr);
+        ASSERT_FALSE(P->over);
+        PJ_COORD input;
+        PJ_COORD input_notOver;
+
+        input.xyz.x = 0; // Lat in deg
+        input.xyz.y = 140;  // Long in deg
+        input.xyz.z = 0;
+
+        input_notOver.xyz.x = 0; // Lat in deg
+        input_notOver.xyz.y = -220;  // Long in deg
+        input_notOver.xyz.z = 0;
+
+        auto output = proj_trans(P, PJ_FWD, input);
+        auto output_notOver = proj_trans(P, PJ_FWD, input_notOver);
+
+        EXPECT_GT(output.xyz.x, 0);
+        EXPECT_GT(output_notOver.xyz.x, 0);
+
+        EXPECT_NEAR(output.xyz.x, 15584728.711058298, 1e-8);
+        EXPECT_NEAR(output_notOver.xyz.x, 15584728.711058298, 1e-8);
+
+        proj_destroy(P);
+    }
+
+    {
+        //EPSG:4326 -> EPSG:27700 has more than one coordinate operation candidate.
+        const char* const options[] = { "FORCE_OVER=YES", nullptr };
+        auto P = proj_create_crs_to_crs_from_pj(ctx, epsg4326, epsg27700, nullptr, options);
+        ASSERT_TRUE(P != nullptr);
+        ASSERT_TRUE(P->over);
+        PJ_COORD input;
+        PJ_COORD input_over;
+
+        input.xyz.x = 0; // Lat in deg
+        input.xyz.y = 140;  // Long in deg
+        input.xyz.z = 0;
+
+        input_over.xyz.x = 0; // Lat in deg
+        input_over.xyz.y = -220;  // Long in deg
+        input_over.xyz.z = 0;
+
+        auto output = proj_trans(P, PJ_FWD, input);
+        auto output_over = proj_trans(P, PJ_FWD, input_over);
+
+        //Doesn't actually change the result for this tmerc transformation.
+        EXPECT_NEAR(output.xyz.x, 4980122.749364435, 1e-8);
+        EXPECT_NEAR(output.xyz.y, 14467212.882603768, 1e-8);
+        EXPECT_NEAR(output_over.xyz.x, 4980122.749364435, 1e-8);
+        EXPECT_NEAR(output_over.xyz.y, 14467212.882603768, 1e-8);
+        
+        proj_destroy(P);
+    }
+
+    {
+        //Negative test for 27700.
+        const char* const options[] = { "FORCE_OVER=NO", nullptr };
+        auto P = proj_create_crs_to_crs_from_pj(ctx, epsg4326, epsg27700, nullptr, options);
+        ASSERT_TRUE(P != nullptr);
+        ASSERT_FALSE(P->over);
+        PJ_COORD input;
+        PJ_COORD input_over;
+
+        input.xyz.x = 0; // Lat in deg
+        input.xyz.y = 140;  // Long in deg
+        input.xyz.z = 0;
+
+        input_over.xyz.x = 0; // Lat in deg
+        input_over.xyz.y = -220;  // Long in deg
+        input_over.xyz.z = 0;
+
+        auto output = proj_trans(P, PJ_FWD, input);
+        auto output_over = proj_trans(P, PJ_FWD, input_over);
+
+        EXPECT_NEAR(output.xyz.x, 4980122.749364435, 1e-8);
+        EXPECT_NEAR(output.xyz.y, 14467212.882603768, 1e-8);
+        EXPECT_NEAR(output_over.xyz.x, 4980122.749364435, 1e-8);
+        EXPECT_NEAR(output_over.xyz.y, 14467212.882603768, 1e-8);
+
+        proj_destroy(P);
+    }
+
+    proj_destroy(epsg27700);
+    proj_destroy(epsg4326);
+    proj_destroy(epsg3857);
+    proj_context_destroy(ctx);
+}
+
 } // namespace