From 8bc7c00f12e81337f6b9fe4418a613276db434cd Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Tue, 10 Jul 2018 18:01:12 +0200 Subject: [PATCH] vgridshift: add a +multiplier={value} parameter As mentionned in #1071, it is often unclear how the offset of a vertical grid is applied. --- .../operations/transformations/vgridshift.rst | 16 ++++++++++++++++ src/PJ_vgridshift.c | 12 ++++++++++-- test/gie/more_builtins.gie | 8 ++++++-- 3 files changed, 32 insertions(+), 4 deletions(-) diff --git a/docs/source/operations/transformations/vgridshift.rst b/docs/source/operations/transformations/vgridshift.rst index e671779bbb..484ad4bfb0 100644 --- a/docs/source/operations/transformations/vgridshift.rst +++ b/docs/source/operations/transformations/vgridshift.rst @@ -113,3 +113,19 @@ Optional .. include:: ../options/t_final.rst .. versionadded:: 5.1.0 +.. option:: +multiplier= + + Specify the multiplier to apply to the grid value in the forward transformation + direction, such that: + + .. math:: + :label: multiplier_formula + + Z_{target} = Z_{source} + multiplier \times gridvalue + + The multiplier can be used to control whether the gridvalue should be added + or sustracted, and if unit conversion must be done (the multiplied gridvalue + must be expressed in metre). + + Note that the default is `-1.0` for historical reasons. +.. versionadded:: 5.2.0 diff --git a/src/PJ_vgridshift.c b/src/PJ_vgridshift.c index 1cd0880e78..39011eb2da 100644 --- a/src/PJ_vgridshift.c +++ b/src/PJ_vgridshift.c @@ -13,16 +13,18 @@ PROJ_HEAD(vgridshift, "Vertical grid shift"); struct pj_opaque_vgridshift { double t_final; double t_epoch; + int forward_multiplier; }; static XYZ forward_3d(LPZ lpz, PJ *P) { + struct pj_opaque_vgridshift *Q = (struct pj_opaque_vgridshift *) P->opaque; PJ_COORD point = {{0,0,0,0}}; point.lpz = lpz; if (P->vgridlist_geoid != NULL) { /* Only try the gridshift if at least one grid is loaded, * otherwise just pass the coordinate through unchanged. */ - point.xyz.z -= proj_vgrid_value(P, point.lp); + point.xyz.z += Q->forward_multiplier * proj_vgrid_value(P, point.lp); } return point.xyz; @@ -30,13 +32,14 @@ static XYZ forward_3d(LPZ lpz, PJ *P) { static LPZ reverse_3d(XYZ xyz, PJ *P) { + struct pj_opaque_vgridshift *Q = (struct pj_opaque_vgridshift *) P->opaque; PJ_COORD point = {{0,0,0,0}}; point.xyz = xyz; if (P->vgridlist_geoid != NULL) { /* Only try the gridshift if at least one grid is loaded, * otherwise just pass the coordinate through unchanged. */ - point.xyz.z += proj_vgrid_value(P, point.lp); + point.xyz.z -= Q->forward_multiplier * proj_vgrid_value(P, point.lp); } return point.lpz; @@ -110,6 +113,11 @@ PJ *TRANSFORMATION(vgridshift,0) { if (pj_param(P->ctx, P->params, "tt_epoch").i) Q->t_epoch = pj_param (P->ctx, P->params, "dt_epoch").f; + /* historical: the forward direction substracts the grid offset. */ + Q->forward_multiplier = -1.0; + if (pj_param(P->ctx, P->params, "tmultiplier").i) { + Q->forward_multiplier = pj_param(P->ctx, P->params, "dmultiplier").f; + } /* Build gridlist. P->vgridlist_geoid can be empty if +grids only ask for optional grids. */ proj_vgrid_init(P, "grids"); diff --git a/test/gie/more_builtins.gie b/test/gie/more_builtins.gie index 7e205e0afd..f2074d64be 100644 --- a/test/gie/more_builtins.gie +++ b/test/gie/more_builtins.gie @@ -205,8 +205,12 @@ operation proj=vgridshift grids=nonexistinggrid.gtx expect failure errno failed_to_load_grid ------------------------------------------------------------------------------- - - +------------------------------------------------------------------------------- +operation proj=vgridshift grids=egm96_15.gtx ellps=GRS80 multiplier=1.0 +tolerance 15 cm +ignore pjd_err_failed_to_load_grid +accept 12.5 55.5 0 0 +expect 12.5 55.5 36.021305084228516 0 ------------------------------------------------------------------------------- Some tests from PJ_hgridshift.c