Merge pull request #3750 from rouault/fix_3749

omerc: catch invalid value of gamma (when specified without alpha)
OSGeo · Jun 25, 2023 · fbe4a45 · fbe4a45
2 parents df8e98a + 09542ec
commit fbe4a45
Show file tree

Hide file tree

Showing 3 changed files with 48 additions and 4 deletions.
diff --git a/docs/source/operations/projections/omerc.rst b/docs/source/operations/projections/omerc.rst
@@ -128,6 +128,11 @@ Central point and azimuth method
     bearing of centre line. If :option:`+alpha` is not given, then
     :option:`+gamma` is used to determine :option:`+alpha`.
 
+    If specifying only :option:`+gamma` without :option:`+alpha`, the maximum
+    value of the absolute value of :option:`+gamma` is a function of the
+    absolute value of :option:`+lat_0`, equal to :math:`90° - |\phi_0|` on a
+    sphere and slightly above on a non-spherical ellipsoid.
+
 .. option:: +lonc=<value>
 
     Longitude of the projection centre. Note that this value is used to

diff --git a/src/projections/omerc.cpp b/src/projections/omerc.cpp
@@ -187,7 +187,7 @@ PJ *PROJECTION(omerc) {
 
         if (fabs(fabs(P->phi0) - M_HALFPI) <= TOL) {
             proj_log_error(
-                P, _("Invalid value for lat_01: |lat_0| should be < 90°"));
+                P, _("Invalid value for lat_0: |lat_0| should be < 90°"));
             return pj_default_destructor(P,
                                          PROJ_ERR_INVALID_OP_ILLEGAL_ARG_VALUE);
         }
@@ -225,12 +225,25 @@ PJ *PROJECTION(omerc) {
             gamma0 = aasin(P->ctx, sin(alpha_c) / D);
             if (!gam)
                 gamma = alpha_c;
-        } else
-            alpha_c = aasin(P->ctx, D * sin(gamma0 = gamma));
+        } else {
+            gamma0 = gamma;
+            alpha_c = aasin(P->ctx, D * sin(gamma0));
+            if (proj_errno(P) != 0) {
+                // For a sphere, |gamma| must be <= 90 - |lat_0|
+                // On an ellipsoid, this is very slightly above
+                proj_log_error(P,
+                               ("Invalid value for gamma: given lat_0 value, "
+                                "|gamma| should be <= " +
+                                std::to_string(asin(1. / D) / M_PI * 180) + "°")
+                                   .c_str());
+                return pj_default_destructor(
+                    P, PROJ_ERR_INVALID_OP_ILLEGAL_ARG_VALUE);
+            }
+        }
 
         if (fabs(fabs(P->phi0) - M_HALFPI) <= TOL) {
             proj_log_error(
-                P, _("Invalid value for lat_01: |lat_0| should be < 90°"));
+                P, _("Invalid value for lat_0: |lat_0| should be < 90°"));
             return pj_default_destructor(P,
                                          PROJ_ERR_INVALID_OP_ILLEGAL_ARG_VALUE);
         }

diff --git a/test/gie/builtins.gie b/test/gie/builtins.gie
@@ -4720,6 +4720,32 @@ operation +proj=omerc    +lat_2=91
 expect     failure errno invalid_op_illegal_arg_value
 
 
+===============================================================================
+# Test limitations on +gamma when specifying without +alpha
+# cf https://github.com/OSGeo/PROJ/issues/3749
+===============================================================================
+
+-------------------------------------------------------------------------------
+operation +proj=omerc    +lat_0=10 +R=6400000 +gamma=80
+-------------------------------------------------------------------------------
+# OK
+
+-------------------------------------------------------------------------------
+operation +proj=omerc    +lat_0=10 +R=6400000 +gamma=80.0000001
+-------------------------------------------------------------------------------
+expect     failure errno invalid_op_illegal_arg_value
+
+-------------------------------------------------------------------------------
+operation +proj=omerc    +lat_0=10 +R=6400000 +rf=300 +gamma=80.01
+-------------------------------------------------------------------------------
+# OK
+
+-------------------------------------------------------------------------------
+operation +proj=omerc    +lat_0=10 +a=6400000 +rf=300 +gamma=80.1
+-------------------------------------------------------------------------------
+expect     failure errno invalid_op_illegal_arg_value
+
+
 ===============================================================================
 # Ortelius Oval
 # 	Misc Sph, no inv.