deformation: Replace +t_obs with +dt

The +t_obs parameter was confusing for users since it effectively
overwrote the observation time in input coordinates. To make it more clear
what is the operation is doing, users are now required to directly specify
the time span for which they wish to apply a given deformation. The parameter
+dt has been added for that purpose. The new parameter is mutually
exclusive with +t_epoch. +dt is used when deformation
for a set amount of time is needed and +t_epoch is used (in
conjunction with the observation time of the input coordinate) when
deformation from a specific epoch to the observation time is needed.

8 files changed, 96 insertions, 53 deletions

diff --git a/docs/source/operations/transformations/deformation.rst b/docs/source/operations/transformations/deformation.rst
index 0a772547..3a9d025c 100644
--- a/docs/source/operations/transformations/deformation.rst
+++ b/docs/source/operations/transformations/deformation.rst
@@ -31,9 +31,9 @@ they represent the physical world better. In PROJ this is done with the deformat
 operation.
 
 The horizontal grid is stored in CTable2 format and the vertical grid is stored in the
-GTX format. Both grids are expected to contain grid-values in units of mm/year.
-Details about the formats can be found in the GDAL documentation. GDAL both reads and
-writes both file formats. Using GDAL for construction of new grids is recommended.
+GTX format. Both grids are expected to contain grid-values in units of
+mm/year. GDAL both reads and writes both file formats. Using GDAL for
+construction of new grids is recommended.
 
 Example
 -------------------------------------------------------------------------------
@@ -59,13 +59,14 @@ to the Danish realisation of ETRS89 is in PROJ described as::
     step    proj = deformation t_epoch = 2000.0
             xy_grids = ./nkgrf03vel_realigned_xy.ct2
             z_grids  = ./nkgrf03vel_realigned_z.gtx
+            inv
             # NKG_ETRF@2000.0 -> ETRF92@2000.0
     step    proj=helmert convention=position_vector s = -0.009420e
             x = 0.03863 rx = 0.00617753
             y = 0.147   ry = 5.064e-05
             z = 0.02776 rz = 4.729e-05
             # ETRF92@2000.0 -> ETRF92@1994.704
-    step    proj = deformation t_epoch = 1994.704 t_obs = 2000.0
+    step    proj = deformation dt = -5.296
             xy_grids = ./nkgrf03vel_realigned_xy.ct2
             z_grids  = ./nkgrf03vel_realigned_z.gtx
 
@@ -86,10 +87,6 @@ component.
 Parameters
 -------------------------------------------------------------------------------
 
-Required
-################################################################################
-
-
 .. option:: +xy_grids=<list>
 
     Comma-separated list of grids to load. If a grid is prefixed by an `@` the
@@ -110,15 +107,22 @@ Required
 
 .. option:: +t_epoch=<value>
 
-    Central epoch of transformation given in decimalyears.
+    Central epoch of transformation given in decimalyears. Will be used in
+    conjunction with the observation time from the input coordinate to
+    determine :math:`dt` as used in eq. :eq:`apply_velocity` below.
+
+    .. note:: :option:`+t_epoch` is mutually exclusive with :option:`+dt`
+
+.. option:: +dt=<value>
 
-Optional
-################################################################################
+    .. versionadded:: 6.0.0
 
-.. option:: +t_obs=<value>
+    :math:`dt` as used in eq. :eq:`apply_velocity` below. Is useful when
+    no observation time is available in the input coordinate or when
+    a deformation for a specific timespan needs to be applied in a
+    transformation. :math:`dt` is given in units of decimalyears.
 
-    Observation time of coordinate(s) given in decimalyears. If not specified,
-    the observation time from the temporal component of 4D input points is used.
+    .. note:: :option:`+dt` is mutually exclusive with :option:`+t_epoch`
 
 Mathematical description
 -------------------------------------------------------------------------------
@@ -182,3 +186,9 @@ conversion from ENU space to cartesian space is done in the following way:
 where :math:`\phi` and :math:`\lambda` are the latitude and longitude of the coordinate
 that is searched for in the grid. :math:`(E, N, U)` are the grid values in ENU-space and
 :math:`(X, Y, Z)` are the corrections converted to cartesian space.
+
+
+See also
+-----------------------------------------------------------------------------
+
+#. :ref:`Behavioural changes from version 5 to 6<differences_deformation>`
\ No newline at end of file
diff --git a/docs/source/usage/differences.rst b/docs/source/usage/differences.rst
index b7c79f53..d5eebe8e 100644
--- a/docs/source/usage/differences.rst
+++ b/docs/source/usage/differences.rst
@@ -72,8 +72,15 @@ default ellipsoid is GRS80 (was WGS84 in previous PROJ versions).
 
 Changes to :ref:`deformation<deformation>`
 ------------------------------------------------------------------
+.. _differences_deformation:
 
-In the initial version of the of :ref:`deformation<deformation>` operation the time span between :math:`t_{obs}` and :math:`t_c` was determine by the expression
+
+Reversed order of operation
+...........................
+
+In the initial version of the of :ref:`deformation<deformation>` operation
+the time span between :math:`t_{obs}` and :math:`t_c` was determined by the
+expression
 
 .. math::
 
@@ -90,4 +97,20 @@ Effectively this means that the direction of the operation has been reversed,
 so that what in PROJ 5 was a forward operation is now an inverse operation and
 vice versa.
 
-Pipelines written for PROJ 5 can be migrated to PROJ 6 by adding :option:`+inv` to forward steps involving the deformation operation. Similarly :option:`+inv` should be removed from inverse steps to be compatible with PROJ 6.
+Pipelines written for PROJ 5 can be migrated to PROJ 6 by adding :option:`
++inv` to forward steps involving the deformation operation. Similarly
+:option:`+inv` should be removed from inverse steps to be compatible with
+PROJ 6.
+
+Removed ``+t_obs``  parameter
+.............................
+
+The ``+t_obs`` parameter was confusing for users since it effectively
+overwrote the observation time in input coordinates. To make it more clear
+what is the operation is doing, users are now required to directly specify
+the time span for which they wish to apply a given deformation. The parameter
+:option:`+dt` has been added for that purpose. The new parameter is mutually
+exclusive with :option:`+t_epoch`. :option:`+dt` is used when deformation
+for a set amount of time is needed and :option:`+t_epoch` is used (in
+conjunction with the observation time of the input coordinate) when
+deformation from a specific epoch to the observation time is needed.
\ No newline at end of file
diff --git a/src/apps/gie.cpp b/src/apps/gie.cpp
index c3622d52..abed11c0 100644
--- a/src/apps/gie.cpp
+++ b/src/apps/gie.cpp
@@ -1082,6 +1082,7 @@ static const struct errno_vs_err_const lookup[] = {
     {"pjd_err_too_many_inits"           ,  -57},
     {"pjd_err_invalid_arg"              ,  -58},
     {"pjd_err_inconsistent_unit"        ,  -59},
+    {"pjd_err_mutually_exclusive_args"  ,  -60},
     {"pjd_err_dont_skip"                ,  5555},
     {"pjd_err_unknown"                  ,  9999},
     {"pjd_err_enomem"                   ,  ENOMEM},
diff --git a/src/proj_internal.h b/src/proj_internal.h
index 9ffcc2b3..b7e3de9e 100644
--- a/src/proj_internal.h
+++ b/src/proj_internal.h
@@ -677,6 +677,7 @@ struct FACTORS {
 #define PJD_ERR_TOO_MANY_INITS          -57
 #define PJD_ERR_INVALID_ARG             -58
 #define PJD_ERR_INCONSISTENT_UNIT       -59
+#define PJD_ERR_MUTUALLY_EXCLUSIVE_ARGS -60
 /* NOTE: Remember to update src/strerrno.cpp, src/apps/gie.cpp and transient_error in */
 /* src/transform.cpp when adding new value */
 
diff --git a/src/strerrno.cpp b/src/strerrno.cpp
index 13e9c757..01097a42 100644
--- a/src/strerrno.cpp
+++ b/src/strerrno.cpp
@@ -68,6 +68,7 @@ pj_err_list[] = {
     "only one +init allowed for non-pipeline operations",              /* -57 */
     "argument not numerical or out of range",                          /* -58 */
     "inconsistent unit type between input and output",                 /* -59 */
+    "arguments are mutually exclusive",                                /* -60 */
 
     /* When adding error messages, remember to update ID defines in
        projects.h, and transient_error array in pj_transform                  */
diff --git a/src/transform.cpp b/src/transform.cpp
index 486178d0..d111d835 100644
--- a/src/transform.cpp
+++ b/src/transform.cpp
@@ -75,14 +75,15 @@ static int adjust_axis( projCtx ctx, const char *axis, int denormalize_flag,
 **
 */
 
-static const int transient_error[60] = {
+static const int transient_error[70] = {
     /*             0  1  2  3  4  5  6  7  8  9   */
     /* 0 to 9 */   0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
     /* 10 to 19 */ 0, 0, 0, 0, 1, 1, 0, 1, 1, 1,
     /* 20 to 29 */ 1, 0, 0, 0, 0, 0, 0, 1, 0, 0,
     /* 30 to 39 */ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
     /* 40 to 49 */ 0, 0, 0, 0, 0, 0, 0, 0, 1, 0,
-    /* 50 to 59 */ 1, 0, 1, 0, 1, 1, 1, 1, 0, 0 };
+    /* 50 to 59 */ 1, 0, 1, 0, 1, 1, 1, 1, 0, 0,
+    /* 60 to 69 */ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,};
 
 
 /* -------------------------------------------------------------------- */
diff --git a/src/transformations/deformation.cpp b/src/transformations/deformation.cpp
index 14f34e4e..c28e1489 100644
--- a/src/transformations/deformation.cpp
+++ b/src/transformations/deformation.cpp
@@ -65,7 +65,7 @@ PROJ_HEAD(deformation, "Kinematic grid shift");
 
 namespace { // anonymous namespace
 struct pj_opaque {
-    double t_obs;
+    double dt;
     double t_epoch;
     PJ *cart;
 };
@@ -167,23 +167,20 @@ static PJ_XYZ forward_3d(PJ_LPZ lpz, PJ *P) {
     struct pj_opaque *Q = (struct pj_opaque *) P->opaque;
     PJ_COORD out, in;
     PJ_XYZ shift;
-    double dt = 0.0;
     in.lpz = lpz;
     out = in;
 
-    if (Q->t_obs != HUGE_VAL) {
-        dt = Q->t_obs - Q->t_epoch;
-    } else {
+    if (Q->dt == HUGE_VAL) {
         out = proj_coord_error(); /* in the 3D case +t_obs must be specified */
-        proj_log_debug(P, "deformation: +t_obs must be specified");
+        proj_log_debug(P, "deformation: +dt must be specified");
         return out.xyz;
     }
 
     shift = get_grid_shift(P, in.xyz);
 
-    out.xyz.x += dt * shift.x;
-    out.xyz.y += dt * shift.y;
-    out.xyz.z += dt * shift.z;
+    out.xyz.x += Q->dt * shift.x;
+    out.xyz.y += Q->dt * shift.y;
+    out.xyz.z += Q->dt * shift.z;
 
     return out.xyz;
 }
@@ -195,10 +192,10 @@ static PJ_COORD forward_4d(PJ_COORD in, PJ *P) {
     PJ_XYZ shift;
     PJ_COORD out = in;
 
-    if (Q->t_obs != HUGE_VAL) {
-            dt = Q->t_obs - Q->t_epoch;
-        } else {
-            dt = in.xyzt.t - Q->t_epoch ;
+    if (Q->dt != HUGE_VAL) {
+        dt = Q->dt;
+    } else {
+        dt = in.xyzt.t - Q->t_epoch ;
     }
 
     shift = get_grid_shift(P, in.xyz);
@@ -215,18 +212,15 @@ static PJ_COORD forward_4d(PJ_COORD in, PJ *P) {
 static PJ_LPZ reverse_3d(PJ_XYZ in, PJ *P) {
     struct pj_opaque *Q = (struct pj_opaque *) P->opaque;
     PJ_COORD out;
-    double dt = 0.0;
     out.xyz = in;
 
-    if (Q->t_obs != HUGE_VAL) {
-        dt = Q->t_obs - Q->t_epoch;
-    } else {
+    if (Q->dt == HUGE_VAL) {
         out = proj_coord_error(); /* in the 3D case +t_obs must be specified */
-        proj_log_debug(P, "deformation: +t_obs must be specified");
+        proj_log_debug(P, "deformation: +dt must be specified");
         return out.lpz;
     }
 
-    out.xyz = reverse_shift(P, in, dt);
+    out.xyz = reverse_shift(P, in, Q->dt);
 
     return out.lpz;
 }
@@ -237,8 +231,8 @@ static PJ_COORD reverse_4d(PJ_COORD in, PJ *P) {
     double dt;
 
 
-    if (Q->t_obs != HUGE_VAL) {
-            dt = Q->t_obs - Q->t_epoch;
+    if (Q->dt != HUGE_VAL) {
+            dt = Q->dt;
         } else {
             dt = in.xyzt.t - Q->t_epoch;
     }
@@ -298,16 +292,29 @@ PJ *TRANSFORMATION(deformation,1) {
         return destructor(P, PJD_ERR_FAILED_TO_LOAD_GRID);
     }
 
-    Q->t_obs = HUGE_VAL;
-    if (pj_param(P->ctx, P->params, "tt_obs").i) {
-       Q->t_obs = pj_param(P->ctx, P->params, "dt_obs").f;
+    Q->dt = HUGE_VAL;
+    if (pj_param(P->ctx, P->params, "tdt").i) {
+       Q->dt = pj_param(P->ctx, P->params, "ddt").f;
+    }
+
+    if (pj_param_exists(P->params, "t_obs")) {
+        proj_log_error(P, "deformation: +t_obs parameter is deprecated. Use +dt instead.");
+        return destructor(P, PJD_ERR_MISSING_ARGS);
     }
 
+    Q->t_epoch = HUGE_VAL;
     if (pj_param(P->ctx, P->params, "tt_epoch").i) {
         Q->t_epoch = pj_param(P->ctx, P->params, "dt_epoch").f;
-    } else {
-       proj_log_error(P, "deformation: +t_epoch parameter missing.");
-       return destructor(P, PJD_ERR_MISSING_ARGS);
+    }
+
+    if (Q->dt == HUGE_VAL && Q->t_epoch == HUGE_VAL) {
+        proj_log_error(P, "deformation: either +dt or +t_epoch needs to be set.");
+        return destructor(P, PJD_ERR_MISSING_ARGS);
+    }
+
+    if (Q->dt != HUGE_VALL && Q->t_epoch != HUGE_VALL) {
+        proj_log_error(P, "deformation: +dt or +t_epoch are mutually exclusive.");
+        return destructor(P, PJD_ERR_MUTUALLY_EXCLUSIVE_ARGS);
     }
 
     P->fwd4d = forward_4d;
diff --git a/test/gie/deformation.gie b/test/gie/deformation.gie
index e4d87fe2..37f62e6e 100644
--- a/test/gie/deformation.gie
+++ b/test/gie/deformation.gie
@@ -12,13 +12,12 @@ The input coordinate is located at lon=60, lam=-160 - somewhere in Alaska.
 <gie>
 
 -------------------------------------------------------------------------------
-Test using both horizontal and vertical grids as well as the +tobs parameter
+Test the +dt parameter
 -------------------------------------------------------------------------------
 operation   +proj=deformation +xy_grids=alaska +z_grids=egm96_15.gtx
-            +t_epoch=2016.0 +t_obs=2000.0 +ellps=GRS80
+            +ellps=GRS80 +dt=16.0 # 2016.0 - 2000.0
 -------------------------------------------------------------------------------
 tolerance   0.1 mm
-direction   inverse
 ignore      pjd_err_failed_to_load_grid
 accept      -3004295.5882503074 -1093474.1690603832 5500477.1338251457
 expect      -3004295.7025       -1093474.2106       5500477.3444
@@ -28,7 +27,7 @@ roundtrip   5
 # Here we test 120W 40N which is well outside the alaska grid.
 accept      -2446353.8001  -4237209.0750  4077985.572
 expect      failure errno grid_area
-accept      -2446353.8001  -4237209.0750  4077985.572   2000
+accept      -2446353.8001  -4237209.0750  4077985.572
 expect      failure errno grid_area
 
 
@@ -46,20 +45,20 @@ expect      -3004295.7025       -1093474.2106       5500477.3444       2000.0
 roundtrip   5
 
 -------------------------------------------------------------------------------
-operation   proj=deformation xy_grids=alaska t_epoch=2016.0 t_obs=2017.0
+operation   proj=deformation xy_grids=alaska +dt=1.0
             ellps=GRS80
 expect      failure pjd_err_no_args
 
-operation   proj=deformation z_grids=egm96_15.gtx t_epoch=2016.0 t_obs=2017.0
+operation   proj=deformation z_grids=egm96_15.gtx +dt=1.0
             ellps=GRS80
 expect      failure pjd_err_no_args
 
 operation   proj=deformation xy_grids=nonexisting z_grids=egm96_15.gtx
-            t_epoch=2016.0 t_obs=2017.0 ellps=GRS80
+            +dt=1.0 ellps=GRS80
 expect      failure pjd_err_failed_to_load_grid
 
 operation   proj=deformation xy_grids=alaska z_grids=nonexisting
-            t_epoch=2016.0 t_obs=2017.0 ellps=GRS80
+            +dt=1.0 ellps=GRS80
 expect      failure pjd_err_failed_to_load_grid
 
 operation   proj=deformation xy_grids=alaska z_grids=nonexisting ellps=GRS80