Merge remote-tracking branch 'osgeo/master' into doc-improvements

5 files changed, 182 insertions, 45 deletions

diff --git a/docs/source/operations/transformations/deformation.rst b/docs/source/operations/transformations/deformation.rst
index 4ae56f54..f10562dd 100644
--- a/docs/source/operations/transformations/deformation.rst
+++ b/docs/source/operations/transformations/deformation.rst
@@ -7,21 +7,21 @@ Kinematic datum shifting utilizing a deformation model
 Perform datum shifts means of a deformation/velocity model.
 
 +-----------------+--------------------------------------------------------------------+
-| **Input type**  | Cartesian coordinates.                                             |
+| **Input type**  | Cartesian coordinates (spatial), decimalyears (temporal).          |
 +-----------------+--------------------------------------------------------------------+
-| **Output type** | Cartesian coordinates.                                             |
+| **Output type** | Cartesian coordinates (spatial), decimalyears (temporal).          |
 +-----------------+--------------------------------------------------------------------+
 | **Options**                                                                          |
 +-----------------+--------------------------------------------------------------------+
-| `xy_grids`      | Comma-separated list of grids to load.                             |
+| `xy_grids`      | Comma-separated list of grids to load. *Required*.                 |
 +-----------------+--------------------------------------------------------------------+
-| `z_grids`       | Comma-separated list of grids to load.                             |
+| `z_grids`       | Comma-separated list of grids to load. *Required*.                 |
 +-----------------+--------------------------------------------------------------------+
-| `t_epoch`       | Central epoch of transformation. [decimalyear]. Only used in       |
-|                 | spatiotemporal transformations.                                    |
+| `t_epoch`       | Central epoch of transformation. [decimalyear]. *Required*.        |
 +-----------------+--------------------------------------------------------------------+
-| `t_obs`         | Observation time of coordinate(s). Mostly useful in 2D and 3D      |
-|                 | transformations. [decimalyear]. *Optional*.                        |
+| `t_obs`         | Observation time of coordinate(s). [decimalyear]. *Optional*.      |
+|                 | If not specified, will be get from the t component of 4D input     |
+|                 | points.                                                            |
 +-----------------+--------------------------------------------------------------------+
 
 The deformation operation is used to adjust coordinates for intraplate deformations.
diff --git a/docs/source/operations/transformations/helmert.rst b/docs/source/operations/transformations/helmert.rst
index 8c862866..941fc676 100644
--- a/docs/source/operations/transformations/helmert.rst
+++ b/docs/source/operations/transformations/helmert.rst
@@ -50,6 +50,8 @@ anoether by means of 3-, 4-and 7-parameter shifts, or one of their 6-, 8- and
 +----------------+--------------------------------------------------------------------+
 | `t_obs`        | Observation time of coordinate(s). Mostly useful in 2D and 3D      |
 |                | transformations. [decimalyear]. *Optional*.                        |
+|                | If not specified, will be get from the t component of 4D input     |
+|                | points.                                                            |
 +----------------+--------------------------------------------------------------------+
 | `exact`        | Use exact transformation equations. *Optional*.                    |
 +----------------+--------------------------------------------------------------------+
diff --git a/src/PJ_helmert.c b/src/PJ_helmert.c
index 229e30c2..5ff8c5d2 100644
--- a/src/PJ_helmert.c
+++ b/src/PJ_helmert.c
@@ -134,7 +134,8 @@ static void update_parameters(PJ *P) {
 
     /* debugging output */
     if (proj_log_level(P->ctx, PJ_LOG_TELL) >= PJ_LOG_TRACE) {
-        proj_log_trace(P, "Transformation parameters for observation t_epoch %g:", Q->t_obs);
+        proj_log_trace(P, "Transformation parameters for observation "
+                       "t_obs=%g (t_epoch=%g):", Q->t_obs, Q->t_epoch);
         proj_log_trace(P, "x: %g", Q->xyz.x);
         proj_log_trace(P, "y: %g", Q->xyz.y);
         proj_log_trace(P, "z: %g", Q->xyz.z);
diff --git a/src/PJ_ortho.c b/src/PJ_ortho.c
index 5b3770d0..7c27cd56 100644
--- a/src/PJ_ortho.c
+++ b/src/PJ_ortho.c
@@ -1,6 +1,7 @@
 #define PJ_LIB__
 #include <errno.h>
 #include "proj.h"
+#include "proj_internal.h"
 #include "projects.h"
 
 PROJ_HEAD(ortho, "Orthographic") "\n\tAzi, Sph.";
@@ -20,37 +21,39 @@ struct pj_opaque {
 
 #define EPS10 1.e-10
 
+static XY forward_error(PJ *P, LP lp, XY xy) {
+    proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION);
+    proj_log_trace(P, "Coordinate (%.3f, %.3f) is on the unprojected hemisphere",
+                   proj_todeg(lp.lam), proj_todeg(lp.phi));
+    return xy;
+}
 
 static XY s_forward (LP lp, PJ *P) {           /* Spheroidal, forward */
-    XY xy = {0.0,0.0};
+    XY xy;
     struct pj_opaque *Q = P->opaque;
     double  coslam, cosphi, sinphi;
 
+    xy.x = HUGE_VAL; xy.y = HUGE_VAL;
+
     cosphi = cos(lp.phi);
     coslam = cos(lp.lam);
     switch (Q->mode) {
     case EQUIT:
-        if (cosphi * coslam < - EPS10) {
-            proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION);
-            return xy;
-        }
+        if (cosphi * coslam < - EPS10)
+            return forward_error(P, lp, xy);
         xy.y = sin(lp.phi);
         break;
     case OBLIQ:
-        if (Q->sinph0 * (sinphi = sin(lp.phi)) + Q->cosph0 * cosphi * coslam < - EPS10) {
-            proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION);
-            return xy;
-        }
+        if (Q->sinph0 * (sinphi = sin(lp.phi)) + Q->cosph0 * cosphi * coslam < - EPS10)
+            return forward_error(P, lp, xy);
         xy.y = Q->cosph0 * sinphi - Q->sinph0 * cosphi * coslam;
         break;
     case N_POLE:
         coslam = - coslam;
                 /*-fallthrough*/
     case S_POLE:
-        if (fabs(lp.phi - P->phi0) - EPS10 > M_HALFPI) {
-            proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION);
-            return xy;
-        }
+        if (fabs(lp.phi - P->phi0) - EPS10 > M_HALFPI)
+            return forward_error(P, lp, xy);
         xy.y = cosphi * coslam;
         break;
     }
@@ -60,13 +63,16 @@ static XY s_forward (LP lp, PJ *P) {           /* Spheroidal, forward */
 
 
 static LP s_inverse (XY xy, PJ *P) {           /* Spheroidal, inverse */
-    LP lp = {0.0,0.0};
+    LP lp;
     struct pj_opaque *Q = P->opaque;
     double  rh, cosc, sinc;
 
+    lp.lam = HUGE_VAL; lp.phi = HUGE_VAL;
+
     if ((sinc = (rh = hypot(xy.x, xy.y))) > 1.) {
         if ((sinc - 1.) > EPS10) {
             proj_errno_set(P, PJD_ERR_TOLERANCE_CONDITION);
+            proj_log_trace(P, "Point (%.3f, %.3f) is outside the projection boundary");
             return lp;
         }
         sinc = 1.;
diff --git a/test/gie/builtins.gie b/test/gie/builtins.gie
index f2ee2e16..9cfb7cce 100644
--- a/test/gie/builtins.gie
+++ b/test/gie/builtins.gie
@@ -3078,28 +3078,156 @@ Orthographic
 ===============================================================================
 
 -------------------------------------------------------------------------------
-operation +proj=ortho   +a=6400000    +lat_1=0.5 +lat_2=2
--------------------------------------------------------------------------------
-tolerance 0.1 mm
-accept  2 1
-expect  223322.760576727 111695.401198614
-accept  2 -1
-expect  223322.760576727 -111695.401198614
-accept  -2 1
-expect  -223322.760576727 111695.401198614
-accept  -2 -1
-expect  -223322.760576727 -111695.401198614
-
-direction inverse
-accept  200 100
-expect  0.001790493 0.000895247
-accept  200 -100
-expect  0.001790493 -0.000895247
-accept  -200 100
-expect  -0.001790493 0.000895247
-accept  -200 -100
-expect  -0.001790493 -0.000895247
-
+Test the equatorial aspect of the Orthopgraphic projection.
+
+Test data from Snyder (1987), table 22, p. 151.
+-------------------------------------------------------------------------------
+operation   +proj=ortho +R=1 +lat_0=0 +lon_0=0
+-------------------------------------------------------------------------------
+tolerance   0.1 mm
+accept      0       0
+expect      0       0
+roundtrip   100
+accept      0       90
+expect      0       1
+roundtrip   100
+accept      0       30
+expect      0       0.5000
+roundtrip   100
+accept      10      50
+expect      0.1116  0.7660
+roundtrip   100
+accept      20      40
+expect      0.2620  0.6428
+roundtrip   100
+accept      30      80
+expect      0.0868  0.9848
+roundtrip   100
+accept      40      70
+expect      0.2198  0.9397
+roundtrip   100
+accept      50      60
+expext      0.3830  0.8660
+roundtrip   100
+accept      60      50
+expect      0.5567  0.7660
+roundtrip   100
+accept      70      20
+expect      0.8830  0.3420
+roundtrip   100
+accept      80      10
+expect      0.9698  0.1736
+roundtrip   100
+accept      90      90
+expect      0       1
+roundtrip   100
+accept      120     0
+expect      failure errno tolerance_condition
+
+direction   inverse
+accept      2   2
+expect      failure errno tolerance_condition
+
+
+-------------------------------------------------------------------------------
+Test the oblique aspect of the Orthopgraphic projection.
+
+Test data from Snyder (1987), table 23, pp. 152-153.
+-------------------------------------------------------------------------------
+operation   +proj=ortho +R=1 +lat_0=40 +lon_0=0
+-------------------------------------------------------------------------------
+tolerance   0.1 mm
+accept      0.0     90
+expect      0.0     0.7660
+roundtrip   100
+accept      20      60
+expect      0.1710  0.3614
+roundtrip   100
+accept      40     -30
+expect      0.5567  -0.8095
+roundtrip   100
+accept      100     70
+expect      0.3368  0.7580
+roundtrip   100
+accept      130     40
+expect      0.5868  0.8089
+roundtrip   100
+accept      170     60
+expect      0.0868  0.9799
+roundtrip   100
+accept      140     20
+expect      failure errno tolerance_condition
+
+direction   inverse
+accept      2   2
+expect      failure errno tolerance_condition
+
+
+-------------------------------------------------------------------------------
+Test the north polar aspect of the Orthopgraphic projection.
+-------------------------------------------------------------------------------
+operation   +proj=ortho +R=1 +lat_0=90 +lon_0=0
+-------------------------------------------------------------------------------
+tolerance 0.1 mm
+accept      0   0
+expect      0   -1
+roundtrip   100
+accept      180 0
+expect      0   1
+roundtrip   100
+accept      90  90
+expect      0   0
+roundtrip   100
+accept      0   90
+expect      0   0
+roundtrip   100
+accept      90  0
+expect      1   0
+roundtrip   100
+accept      180 -90
+expect      failure errno tolerance_condition
+accept      0   -45
+expect      failure errno tolerance_condition
+
+direction   inverse
+accept      2   2
+expect      failure errno tolerance_condition
+
+-------------------------------------------------------------------------------
+Test the south polar aspect of the Orthopgraphic projection.
+-------------------------------------------------------------------------------
+operation   +proj=ortho +R=1 +lat_0=-90 +lon_0=0
+-------------------------------------------------------------------------------
+tolerance 0.1 mm
+accept      0   0
+expect      0   1
+roundtrip   100
+accept      180 0
+expect      0   -1
+roundtrip   100
+accept      90  -90
+expect      0   0
+roundtrip   100
+accept      0   -90
+expect      0   0
+roundtrip   100
+accept      90  0
+expect      1   0
+roundtrip   100
+accept      180 90
+expect      failure errno tolerance_condition
+accept      0   45
+expect      failure errno tolerance_condition
+
+direction   inverse
+accept      2   2
+expect      failure errno tolerance_condition
+
+# Put a point a tiny tiny bit outside the radius of the sphere.
+# Since we are right at the numerical limit of floating point representation
+# it should still results in a correct coordinate.
+accept      0.70710678118 0.7071067812
+expect      45  0
 
 ===============================================================================
 Perspective Conic