A.1.  Model requirements

To be conformant with this Abstract Specification a deformation model shall satisfy the requirements in Table A.1. The attributes of entities specified in the requirements shall be unambiguously identified in an encoding of a deformation model. A specific encoding may use different names or structure for the entities provided that the encoding specification defines how the deformation model attributes are identified.

Table A.1 — Data requirements

Table A.2 — Time function mandatory (M) and optional (O) attributes

A.2.  Implementation requirements

To be conformant with this Abstract Specification a software implementation shall satisfy the requirements of Table A.3. Note that, as described in Clause 6, this does not prescribe the specific formulas that a software should use, but does require that the calculations by the software should obtain results that are practically identical to those defined by the these formulas.

Table A.3 — Implementation requirements

B.1.  Decomposition into elements

This OGC Abstract Specification (AS) assumes that the deformation can be partitioned into one or more elements each of which is characterized by a time function. There are many other ways deformation can be modeled. For example, in geophysical research it is often modeled using uniform dislocation on rectangular fault surfaces in a uniform elastic half space. While such models may be used to research and develop deformation models, common practice in publishing them for coordinate operations is to convert the model to a set of grids defining the surface displacement.

Deformation is caused by a variety of geophysical phenomena, each of which may be represented by one or more elements in a deformation model.

Typically, there is a velocity element representing the long-term tectonic deformation — this is a single element with a linear time function.

Overlaid on this there may be deformation due to earthquakes. The coseismic deformation will be represented by an element with a step time function. Post-seismic deformation following the earthquake may be represented by several elements depending on the complexity of the deformation and the quality of information available. Commonly, these will use exponential and logBaseE time functions. Several versions of a deformation model may be published after an earthquake as the ongoing post-seismic deformation evolves and as measurements and modeling of deformation are refined.

Alternatively, the deformation model may be based directly on observations of displacement rather than modelling from geophysical phenomena. For example, in Japan the dense network of continuous GNSS (cGNSS) stations provides a continuous record of deformation. Each year a new version of the deformation model is generated based on these measurements. Each of these could be represented by an element with a step time function.

B.2.  Requirement for the model to be continuous and invertible

Deformation models described in this AS are intended for use in coordinate operations — either to define time dependent transformations between different CRSs or to track the motion of objects between different epochs within a single CRS.

Generally, the user expectation of such operations is that they are invertible. If a coordinate is transformed from the source CRS to the target CRS, and then the resulting coordinate is transformed back to the source CRS, the final coordinate will be the same as the original coordinate. Similarly, if the model is used as a point motion model, transforming a coordinate from one epoch to another, and then transforming the resulting coordinate back to the original epoch, will result in the original coordinate.

To support this expectation, this AS requires that the spatial function of each element of the deformation model is continuous. This does not guarantee that the resulting coordinate operation is invertible. It is possible to define deformation models that are continuous but are not invertible. However, these are not realistic models. Ensuring that a spatial function is continuous is much easier than ensuring that the resulting coordinate operation is invertible. Continuity is used as a practical substitute for invertibility.

Within a single grid the spatial function is continuous — that is a consequence of the formulas such as bilinear interpolation used to interpolate between the displacements at the grid nodes.

If a spatial function uses more than one grid, the requirement for continuity places two constraints on deformation model producers.

1. Where an element does not cover the full extent of the deformation model then the deformation at the edge of the element must be zero. For example, an element might represent coseismic deformation from an earthquake. The deformation caused by the earthquake may not be significant over the full extent of the model so the element would only cover a subset of that area. To ensure continuity, the displacements at the edge of the element should be set to zero, even though the actual displacements will be slightly different from zero.

2. Where an element is represented by multiple overlapping or nested grids there will be transitions across which the selected grid changes from one grid to another. At these transitions the values interpolated from each grid must be the same. For example, the displacement at the edge of a child grid will generally have to exactly match that of its parent grid. However, if the edge abuts the edge of another child grid, then its displacement only needs to match that of the abutting grid.

In practice the second of these constraints means that parent and child grids must be aligned. That is to say that any node of a parent grid that lies within the child grid must also be a node of the child grid. The child grid will generally have smaller grid cells and will include nodes not in the parent grid. Where grids are abutting then along the common edge at least one of the grids should include all the nodes of the other. If one of the grids has a smaller grid cells than the other it may include additional nodes along the common edge that are not nodes of the coarser grid.

These constraints, and indeed the use of a gridded model, mean that the model cannot represent the true deformation exactly. For example, where there is surface faulting, the actual deformation may not be continuous across a fault line. The deformation also may not be invertible (at least as a function of horizontal position) in areas of thrust faulting where points originally on opposite sides of the fault may be moved to the same horizontal position (though at different heights).

For model producers, this can create a tension between the scientific desire to represent deformation as accurately as possible, and the practical constraints of a deformation model suitable for coordinate operations. This AS is framed within the context of coordinate operations — it does not purport to model deformation for scientific use.

The challenge for model producers is to decide the extent to which the deformation can be usefully modeled to best support the geospatial community.

B.3.  no-data and zero values

There may be areas within the extent of the deformation model grids where the deformation is either not known or is not within the model producer’s jurisdiction. In these areas the model producer may prefer not to specify an unreliable or unauthoritative displacement.

Grids must be defined by rectangular regions, which are unlikely to match the extent of the known deformation. To avoid grids being used in areas where the deformation is undefined, it is desirable that an identifiable no-data value is used for displacements at nodes where this applies.

The implementation of no-data will be dictated by the format used to carry the deformation model. This specification requires that an implementation supports a representation of a no-data value which indicates that a value is undefined.

The no-data value for a displacement is different from a zero displacement. A value of zero is used where there is no significant displacement. A no-data value is used where the displacement is unknown.

Instead of using no-data values, a producer may choose to specify a value with a large uncertainty. This may be preferable to ensure a continuous and invertible model, as discussed in Annex B.2. However, this will not be recognized by software that does not use uncertainty information.

If any element of the model evaluates to no-data at a given location and time, then the displacement (or uncertainty) at that location is undefined. Figure B.1 illustrates the use of no-data values in the model. In this figure the gray area shows the region where displacement is undefined and the square marks the total extent of the model. Outside this area the model cannot be evaluated. Within the gridded spatial function there are a number of nodes at which the displacement is not defined — these are identified with a no-data value. Where these nodes are required to interpolate a displacement, which is any grid cell they are on the boundary of, the displacement cannot be calculated.

Figure B.1 — Influence of no-data values

The gray area, where displacement is undefined, might include coastal regions where the deformation of the seabed is not measured, or it may cross a jurisdictional boundary. Since the grid is rectangular, it may include regions where the deformation is not known, which are represented in the model by a no-data value.

B.4.  Alternative interpolation methods

This AS only specifies bilinear interpolation. A producer may choose to use, and an implementation may support, other methods, such as bicubic or biquadratic interpolation. However they are not defined in this AS. The calculation of displacements by a specific implementation may differ from that expected by the producer.

The reasons this AS only specifies bilinear interpolation are:

• it is easier to ensure continuity in a nested grid structure if bilinear interpolation is used, and

• there are few implementation choices, and hence less ambiguity, in using bilinear interpolation.

There are fewer implementation choices for bilinear interpolation because calculating the displacement at a point only uses the values from the four nodes of the grid cell in which the point lies. However, more sophisticated methods, such as biquadratic and bicubic interpolation, use values from the nodes of cells adjacent to the grid cell in which the calculation point lies. These nodes may not exist for a point within a cell on the edge of the grid. There is more than one way that implementation can be written to handle this situation, and there is no “right way” to do it. Different choices may give different results.

B.5.  Deformation models near poles

The formulas defined in this AS are not suitable for transforming geographic (latitude and longitude) coordinates close to the north and south poles. Both the means of interpolating within a grid cell, and the formulas for applying the displacement to a coordinate, may yield unintended results as described below. The simplest way to represent deformation in such cases is to use a suitable projected (easting and northing) CRS. Future versions of this AS may offer better support for geographic CRSs in polar regions by offering alternative methods for interpolating displacements and adding displacements to coordinates.

Near the poles, grid cells may span a large longitude range. Consequently, the east and north components of a vector at adjacent grid nodes may lie in quite different directions, as illustrated in Figure B.2. Interpolating east and north components independently as described in Clause 6.1.2 fails to account for this. For example, if the grid cell spans 1° of longitude and the displacement is 1m then bilinear interpolation may give rise to an error of up to about 2cm (the length of the displacement vector multiplied by the cosine of the longitude range of the grid cell). This can be mitigated by using a finer longitude grid spacing.

Figure B.2 — Different directions of east and north components at grid nodes and a calculation point

Where the source and target CRSs are geographic, then adding east north offsets to the longitude may not be appropriate. This is illustrated in Figure B.3 where the gray vector shows the result of adding an east displacement to the longitude coordinate, and the black vector shows the result applying the same east displacement in the direction of the east vector component. Close to the pole, applying an east displacement by changing the longitude coordinate, as defined by the formulas in Clause 6.4, may give a significantly different result than displacing the coordinate in the direction of the east vector. At the pole itself, changing the longitude coordinate makes no difference to the location at all.

Figure B.3 — Comparison of vector and angular displacement near a pole

This is only a significant issue very close to a geographic pole. The impact of this issue diminishes rapidly away from the pole. For displacement d at a point at distance R from the pole, the difference is approximately d*(1-cos(d/R)), or approximately d³/2R². For example, applying 1 m east displacement 1 km from the pole by changing the longitude would incur an error of only 5.10^-7m.

B.6.  Time function epochs

The time functions used for deformation model elements can be modified by specifying a function reference epoch t₀, a start epoch t_s, and an end epoch t_e. The effect of these is illustrated in Figure B.4 showing a base hyperbolic tangent function without these parameters (graph A) and the same function with each parameter applied in turn. The base function has an event epoch t_v = 2013.8 and a time constant τ = 0.5 years.

Figure B.4 — The effects of adding a start epoch (t_s), end epoch (t_e), and function reference epoch (t₀) to a base hyperbolic tangent function

The base function shows an event building in magnitude from 0.0 to 1.0. This might represent, for example, a slow slip event. However, there is not epoch at which the value of the function is exactly 0 or 1, even though it appears to be in the graph. The value approaches 0 exponentially as the epoch moves further into the past, and approaches 1 asymptotically as the epoch moves further into the future. The rate of change becomes infinitesimally small more than a few years from the event date for typical models of slow slip events.

The start and end epochs can be used to remove any change before the start epoch or after the end epoch. Setting a start epoch (graph B) applies the base function value at that epoch to any time before the start epoch. Similarly, setting an end epoch (graph C) applies the base function value at that epoch to any time after the end epoch. A consequence of this is that the start and end values are no longer 0 and 1. Generally a deformation model producer would use a function reference epoch as well as a start epoch, as described below, so that the time function is zero before the start date.

Using start and end epochs has practical value for deformation model producers managing a CRS. Using the start epoch means that displacements calculated before that epoch do not change. This may be more realistic, as there is at time at which the slow slip event starts, and before which it is causing no deformation.

The function reference epoch t₀ is an epoch at which the function is forced to be zero, and therefore when the calculated displacement for the event will be zero. A constant offset is added to the function to force it to be zero at this epoch. In graph D the value -0.8 is added to the time functions so that it is zero at the function reference epoch 2014.2.

In a more realistic example, the start epoch and function reference epoch would both be used. In graph B, for example, the producer could set the function reference epoch to the same value as the start epoch — 2013.2. This would ensure that the time function is exactly zero for any epoch before this.

One usage of the function reference epoch is to describe deformation that occurred before the reference epoch of a datum. For example, an earthquake occurring before the datum epoch may be modeled by a step function with a function reference epoch set to the datum epoch. This is illustrated in Figure B.5 which shows the modified step time function where the event epoch is 2014.5 and the datum epoch is 2015.2. As the datum is referenced to 2015.2 the earthquake has no effect on coordinates observed at or after that date. However, to calculate the location of an object before the earthquake, the displacements caused by the earthquake must be subtracted from the object coordinates. By setting the function reference epoch to 2015.2, the step function becomes -1.0 for dates before the earthquake, and 0.0 for dates after it — a “reverse step” function.

Figure B.5 — A "reverse step" function modeling an earthquake occurring before the datum reference epoch (defined by the function reference epoch t₀)

C.1.  Discontinuities in the displacement field

One difficulty with representing deformation using gridded models is that they cannot represent discontinuities in the displacement field such as occurs in areas of surface faulting. There may be very large relative coseismic movement across a fault trace which may traverse one or more grids in the model. There may be breaks between segments of a fault trace between which the ground is distorted in unpredictable and often unmeasurable ways.

These displacements cannot be represented in a gridded model in which displacements are interpolated smoothly across each grid cell. This is illustrated in Figure C.1 showing a modeled horizontal displacement field. The left-hand diagram shows the actual displacement in two grid cells. A right lateral strike slip fault runs through the upper grid cell, shown by the gray line. The right-hand diagram shows the displacement interpolated from the values at the grid nodes (red vectors). In the cell containing the fault the interpolated displacement is very different from the actual displacement, whereas in the adjacent cell they are virtually indistinguishable.

Figure C.1 — Ground displacement (left diagram) and interpolated displacement (right diagram) in a grid spanning a fault trace (hashed line).

Even if the deformation model were perfect, it might not be useful for transforming geographic features near a fault trace. Following an earthquake there are many ways that spatial features may be changed — often not reflecting the exact deformation of the Earth’s surface. For example, the line of a cadastral boundary may be defined by legal considerations which override the physical movement of the ground. A road crossing a fault trace may be physically realigned to remain usable. These considerations cannot be accounted for by the deformation model. In many cases data must be recaptured in such areas to restore spatial integrity.

It may be valuable for the users of deformation models to know when a transformed features is likely to be affected by large errors which cannot be represented in the model. Conceptually, software would generate a warning to users. For example, this would happen when transforming a feature spanning a fault trace and over a time range that includes the time of the earthquake.

Four mechanisms have been suggested that could either alert users to such a situation when they do conversions, or allow them to investigate when it could affect their work. These mechanisms are as follows.

1. Add metadata that defines the affected location and time. This could be an area defined in a standard spatial format, such as well-known text (WKT), associated with a date or epoch attribute. The spatial definition could be either included in the deformation model metadata, or as an online resource referenced in the metadata. This AS does not describe any standard for presenting such metadata. It would be up to the user to identify and interpret it.

2. Use a large uncertainty to indicate where there is likely to be a large difference between the interpolated and actual displacements. While this can be used to indicate the approximate area affected by discontinuity, it is spatially imprecise. A large uncertainty value at a node affects all four grid cells surrounding the node. This differs from the interpolation error, shown in Figure C.1, which affects only the grid cells containing the fault trace — adjacent cells are not affected. Nonetheless this is an option available to producers, provided it is understood that it may overestimate the error in the vicinity of a fault trace.

3. Analyze the displacement vectors at the grid nodes in the implementation software. Grid cells likely to be affected by large interpolation errors could be identified by analyzing the displacement vectors themselves. A simple measure of this uncertainty within a grid cell is the maximum vector difference between any two of the displacement vectors at the four nodes of the cell. For a grid cell spanning the fault trace, the displacement vectors for nodes on one side of the fault are in the opposite direction to those of nodes on the other side of the fault so there is a large difference between these vectors. Conversely, for a grid cell to one side of the fault trace, the four displacement vectors are likely to have similar direction and magnitude, so the maximum difference between any two of these vectors is small. This method requires special purpose software to analyze the grids and determine the affected grid cells. The software would then have to alert users when these cells were used in a coordinate operation. A benefit of this approach is that it provides an indication of the magnitude of the potential error.

4. Use an additional grid node attribute — a “flag” value that specifically identifies affected cells. This is a Boolean value that is true or false. There are several options for how this could identify the affected cell. One option is to use a convention defining which node identifies a grid cell. For example, the node on row i and column j could identify the cell between grid nodes on rows i and i+1 and columns j and j+1. Another option is to set the flag as true on all nodes of an affected grid cell. A grid cell is identified as affected only if the flag value is true on all four of its nodes. Such an attribute would require special processing by software.

One difficulty with all these mechanisms is that coordinate operations are “black box” operations to most users — they often happen within software applications as data from different CRSs or different epochs are combined for analysis and presentation. Often this will be done on remote servers presenting information to users in web browsers. Users may not be aware that a coordinate operation is happening. The mechanisms above would need extensive implementation in the spatial software stack for a warning generated in a coordinate operation to be presented to users in a meaningful way. Given this practical difficulty, this AS does not explicitly specify how discontinuity should be handled.

Although handling discontinuities is not defined in this AS, the first two of these options, using metadata and using uncertainty, are available to model producers already.

The third method, analyzing the displacement vectors, can be used with any deformation model. An approach like this could be developed in software for analyzing deformation models independently of coordinate operation software. This is an area that may warrant further research once more deformation models, and particularly models incorporating coseismic displacements, become available for analysis.

C.2.  Vector interpolation in polar regions

As noted in Annex B.5, this AS does not recommend the use of geographic CRSs (i.e. latitude/longitude) for deformation models close to the geographic poles. In such cases, the horizontal displacement vector cannot be interpolated or applied correctly. Early versions of this AS did include an alternative interpolation method, and an alternative method of adding horizontal displacements to a geographic coordinate. In each case the method involved transforming the displacements to and from a geocentric Cartesian CRS.

These were removed from the AS for two reasons:

1. Currently there is no requirement for this capability, and

2. Further research should be undertaken before recommending specific methods.