14.2.3 Gridding Methods for Randomly Spaced Data

When your worksheet data are not regularly spaced in the X and Y dimensions, then you should use one of Origin's 5 gridding methods for random or irregularly-spaced XY data:

Renka Cline

This gridding method is based on the algorithm of Renka and Cline (1984). In sequence, the primary operations are:

  1. Triangulation. A triangulation is performed on the XY data. Triangles are as nearly equiangular as possible.
  2. Gradient Estimation. Estimate a gradient in the XY directions for each node, as the partial derivatives of a quadratic function.
  3. Interpolation. For an arbitrary point P, compute the interpolated value using the data values and gradient estimates at each of the three vertices of the triangle which contains P.

In the case of 200~1000 uniformly distributed data points, the Renka-Cline method is a good choice.

Shepard

This method implements the modified Shepard gridding method as described by Franke and Nielson. This is a distance-based weighted gridding method which interpolates data by:

F(x,y)=\sum_{i=1}^nF_iW_i(x,y) ,

where F_i\,\! is the underlying function at nodes (x_i\,\!, y_i\,\!), and W_i(x, y)\,\! is the weights. To localize the calculation, F_i\,\! and W_i\,\! are computed from the data points lying within the circle with center (x_i\,\!, y_i\,\!) of radius R_q\,\! and R_w\,\!, respectively.

To begin with, the weights are defined as:

w_i(x,y)=\frac{w_i(x,y)}{\sum_{k=1}^nw_k(x,y)}.

For a given radius R_w\,\!, the relative weight w_k\,\! is:

w_k=[\frac{(R_w-d_k)_{+}}{R_wd_k}]^2 for

\begin{cases} R_w-d_k \, if \, d_k<R_w,\\ 0 \, if \, d_k\geq R_w, \end{cases}

and d_k\,\! is the Euclidean distance between (x, y) and (x_k\,\!, y_k\,\!):

d_k=||(x_k,y_k)-(x,y)||_2\,\!.

For any R_w\,\! > 0, we have:

W_i(x_j,y_j)\,\!= \begin{cases} 1,i=j\\ 0,otherwise \end{cases}

\sum_{k=1}^nW_k(x,y)=1.

Subsequently, the nodal function F_i\,\! is replaced by a local approximation function Q_k\,\!:

Q_k(x,y)=c_{k1}(x-x_k)^2+c_{k2}(x-x_k)(y-y_k)+c_{k3}(y-y_k)^2+c_{k4}(x-x_k)+c_{k5}(y-y_k)+F_k\,\!

Q_k\,\! is the weighted least-square quadratic fitted function to the data located within R_q\,\! of nodal points. So the coefficients minimize:

\varepsilon _k=\sum_{i=1,j\neq k}^N\omega _i(x_k,y_k)[c_{k1}(x_i-x_k)^2+\ldots +c_{k5}(y_i-y_k)+F_k-F_i]^2

for

\omega _i(x,y)=[\frac{(R_q-d_i)_{+}}{R_qd_i}]^2=[\frac{(R_q-||(x_i,y_i)-(x,y)||_2)_{+}}{R_q||(x_i,y_i)-(x,y)||_2}]^2.

It can be seen above that the interpolation function is a local approximation function and depends on the radius of influence about nodal points, R_q\,\! and R_w\,\!. Two integers, N_q\,\! and N_w\,\! are used to calculate R_q\,\! and R_w\,\! (these are the parameters q and w of the function, and are referred to as the Quadratic Interpolant Locality Factor and Weight Function Locality Factor, respectively): (for more information, please see "The XYZ Gridding Dialog Box" later in this section.)

R_q=\frac D2\sqrt{\frac{N_q}n} and R_w=\frac D2\sqrt{\frac{N_w}n},

where n is the number of data points and D is the maximum distance between any pair of data points. So N_q\,\! and N_w\,\! can be considered to be the average numbers of data points lying within distances R_q\,\! and R_w\,\! respectively, of each node.

Increasing the values of N_q\,\! and N_w\,\! will make the calculation more global; likewise, decreasing their values will make the calculation more local. Generally speaking, setting N_q \approx 2N_w works well. By default, N_q= 18\,\!, and N_w= 9\,\!. However, the following constraints must be satisfied: 0 < N_w \leq N_q.

The XYZ Gridding dialog box calls the NAG library to perform the Shepard gridding method. Origin also provides these X-Function-based Shepard gridding methods: xyz_shep_nag and xyz_shep, which implements this method using R_q\,\! and R_w\,\!.

Thin Plate Spline (TPS)

This function provides a method for random matrix conversion based on Thin Plate Spline (TPS) algorithm. Thin Plate Spline is a physical interpolation method. To generate gridding data, this method assumes that all the data points are distributed on a thin, elastic plate, or spline. The plate is constrained at the grid points and forms a 2-dimensional surface by spanning the grid points. The surface is deformed between the points to form a likely fit to the data. The best results are generally found by minimizing the so-called "bending energy function" of the spline.

Since this method involves minimizing the bending energy, the less deformation of the plate the better. This calculation is similar to minimum curvature computing. Surface plots produced by TPS gridding may exhibit a greater degree of smoothness than those produced by other methods, so this method is best suited to interpolation of locally flat surfaces.

The mathematical description of the TPS algorithm is, given the bending energy function:

I=\iint_{R^2}(f_{xx}^2+2f_{xy}^2+f_{yy}^2)dxdy,

and the minimizing function is:

f(x,y)=a_1+a_xx+a_yy+\sum_{i=1}^pw_iU(||(x_i,y_i)-(x,y)||),

whereU(r)=r^2\log (r^2)\,\!.

To perform TPS gridding, you need to specify the Smoothing parameter, which controls the smoothness of the interpolated surface and the Extrapolation parameter, which influences matrix cell values lying outside of the original data range.

For more detailed information on TPS, please see Donato and Belongie, Approximation Methods for Thin Plate Spline Mappings and Principal Warps.

Kriging Correlation

Kriging -- named for mining engineer D. G. Krige -- is an established geostatistical method for interpolating spatial data. This technique employs a weighted moving average interpolation (extrapolation) method that minimizes the estimated variance of a predicted point (grid nodes) from the weighted average of its neighbors. The weighted value is determined by the spatial correlation structure of the original data.

This algorithm requires a model of spatial continuity, or dependence. Typically, we separate the process into two steps:

  1. We generate semivariograms of the regionalized variables being estimated. Origin calculates the distance information and produces a table of semivariance values stored in a table, and weights the semivariance data. This semivariogram describes the increasing difference or decreasing correlation between samples in the input dataset.
  2. We use this dependence model to calculate the spatial point value of z.

When performing kriging, you can control the process by configuring the following parameters:

  • Minimum Points: sets the minimum number of control points needed for calculation of z.
  • Search Radius: determines the set of nearby control points used in calculation of z.

For more detailed information on kriging, please see Stein, Interpolation of Spatial Data.

Weighted Average

The weighted average method is a simple weighted average of the points with the weight equal to 1/r, where r is the distance of each point from the cell within the searching radius. If there is no value within the search radius, the radius is increased until at least one point is encountered. Increasing the search radius means that each point is more interrelated to neighboring points, producing a smoother surface that may lose fine details.

English | Deutsch | 日本語