Kinematic datum shifting utilizing a deformation model

Perform datum shifts means of a deformation/velocity model.

Input type

Cartesian coordinates.

Output type

Cartesian coordinates.

Options

xy_grids

Comma-separated list of grids to load.

z_grids

Comma-separated list of grids to load.

t_epoch

Central epoch of transformation. [decimalyear]. Only used in spatiotemporal transformations.

t_obs

Observation time of coordinate(s). Mostly useful in 2D and 3D transformations. [decimalyear]. Optional.

The deformation operation is used to adjust coordinates for intraplate deformations. Usually the transformation parameters for regional plate-fixed reference frames such as the ETRS89 does not take intraplate deformation into account. It is assumed that tectonic plate of the region is rigid. Often times this is true, but near the plate boundary and in areas with post-glacial uplift the assumption breaks. Tntraplate deformations can be modelled and then applied to the coordinates so that they represent the physical world better. In PROJ this is done with the deformation 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.

Example

In [Häkli2016] coordinate transformation including a deformation model is described. The paper describes how coordinates from the global ITRFxx frames are transformed to the local Nordic realisations of ERTS89. Scandinavia is an area with significant post-glacial rebound. The deformations from the post-glacial uplift is not accounted for in the offcial ETRS89 transformations so in order to get accurate transformations in the Nordic countries it is necesarry to apply the deformation model. The transformation from ITRF2008 to the Danish realisation of ETRS89 is in PROJ described as:

proj =  pipeline ellps = GRS80
        # ITRF2008@t_obs -> ITRF2000@t_obs
step    init = ITRF2008:ITRF2000
        # ITRF2000@t_obs -> ETRF2000@t_obs
step    proj=helmert t_epoch = 2000.0 transpose
        x =  0.054  rx =  0.000891 drx =  8.1e-05
        y =  0.051  ry =  0.00539  dry =  0.00049
        z = -0.048  rz = -0.008712 drz = -0.000792
        # ETRF2000@t_obs -> NKG_ETRF00@2000.0
step    proj = deformation t_epoch = 2000.0
        xy_grids = ./nkgrf03vel_realigned_xy.ct2
        z_grids  = ./nkgrf03vel_realigned_z.gtx
        # NKG_ETRF@2000.0 -> ETRF92@2000.0
step    proj=helmert transpose 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
        xy_grids = ./nkgrf03vel_realigned_xy.ct2
        z_grids  = ./nkgrf03vel_realigned_z.gtx

From this we can see that the transformation from ITRF2008 to the Danish realisation of ETRS89 is a combination of Helmert transformations and adjustments with a deformation model. The first use of the deformation operation is:

proj = deformation t_epoch = 2000.0
xy_grids = ./nkgrf03vel_realigned_xy.ct2
z_grids  = ./nkgrf03vel_realigned_z.gtx

Here we set the central epoch of the transformation, 2000.0. The observation epoch is expected as part of the input coordinate tuple. The deformation model is described by two grids, specified with xy_grids and z_grids. The first is the horizontal part of the model and the second is the vertical component.

Mathematical description

Mathematically speaking, application of a deformation model is simple. The deformation model is represented as a grid of velocities in three dimensions. Coordinate corrections are applied in cartesian space. For a given coordinate, \((X, Y, Z)\), velocities \((V_X, V_Y, V_Z)\) can be interpolated from the gridded model. The time span between \(t_c\) and \(t_{obs}\) determine the magnitude of the coordinate correcton as seen in eq. (1) below.

(1)\[\begin{split}\begin{align} \begin{pmatrix} X \\ Y \\ Z \\ \end{pmatrix}_B = \begin{pmatrix} X \\ Y \\ Z \\ \end{pmatrix}_A + (t_c - t_{obs}) \begin{pmatrix} V_X \\ V_Y \\ V_Z \\ \end{pmatrix} \end{align}\end{split}\]

Corrections are done in cartesian space.

Coordinates of the gridded model are in ENU (east, north, up) space because it would otherwise require an enormous 3 dimensional grid to handle the corrections in cartesian space. Keeping the correction in lat/long space reduces the complexity of the grid significantly. Consequentyly though, the input coordinates needs to be converted to lat/long space when searching for corrections in the grid. This is done with cart operation. The converted grid corrections can then be applied to the input coordinates in cartesian space. The conversion from ENU space to cartesian space is done in the following way:

(2)\[\begin{split}\begin{align} \begin{pmatrix} X \\ Y \\ Z \\ \end{pmatrix} = \begin{pmatrix} -\sin\phi \cos\lambda N - \sin\lambda E + \cos\phi \cos\lambda U \\ -\sin\phi \sin\lambda N + \sin\lambda E + \cos\phi \sin\lambda U \\ \cos\phi N + \sin\phi U \\ \end{pmatrix} \end{align}\end{split}\]

where \(\phi\) and \(\lambda\) are the latitude and longitude of the coordinate that is searched for in the grid. \((E, N, U)\) are the grid values in ENU-space and \((X, Y, Z)\) are the corrections converted to cartesian space.