Welcome to the Harris Geospatial product documentation center. Here you will find reference guides, help documents, and product libraries.


IMSL_BSINTERP

IMSL_BSINTERP

The IMSL_BSINTERP function computes a one- or two-dimensional spline interpolant.



Syntax


Result = IMSL_BSINTERP(xdata, fdata [, /DOUBLE] [, XKNOTS=value] [, XORDER=value] [, YKNOTS=value] [, YORDER=value])

or

Result = IMSL_BSINTERP(xdata, ydata, fdata [, DOUBLE=value] [, XKNOTS=value] [, XORDER=value] [, YKNOTS=value] [, YORDER=value])

Return Value


A structure containing information that defines the one- or two-dimensional spline.

Arguments


If a one-dimensional spline is desired, then the arguments xdata and  fdata are required. If a two-dimensional, tensor-product spline is desired, then xdata, ydata, and fdata are required.

xdata

Array containing the abscissas in the x-direction of the interpolation problem.

ydata

Array containing the abscissas in the y-direction of the interpolation problem.

fdata

Array containing the ordinates of the interpolation problem. If a one-dimensional spline is being computed, then fdata (i) is the data value at xdata (i). If a two-dimensional spline is being computed, then fdata is a two-dimensional array, where fdata (i, j) is the data value at (xdata (i), ydata (i)).

Keywords


DOUBLE

If present and nonzero, double precision is used.

XKNOTS

Specifies the array of knots in the x-direction to be used when computing the definition of the spline. Default: knots are selected by the IMSL_BSKNOTS function using its defaults.

XORDER

Specifies the order of the spline in the x-direction. Default: XORDER = 4, i.e., cubic splines.

YKNOTS

Specifies the array of knots in the y-direction to be used when computing the definition of the spline. Default: knots are selected by the IMSL_BSKNOTS function using its defaults.

YORDER

Specifies the order of the spline in the y-direction. If a one-dimensional spline is being computed, then YORDER has no effect on the computations. Default: YORDER = 4, i.e., cubic splines.

Discussion


The IMSL_BSINTERP function is designed to compute either a one-dimensional spline interpolant or two-dimensional, tensor-product spline interpolant to input data. The decision of whether to compute the one- or two-dimensional spline is based on the number of arguments passed to the function. Keywords are provided to allow the user to specify the order of the spline and the knots used for the spline. When computing a one-dimensional spline, the available keywords are XORDER and XKNOTS. When computing a two-dimensional spline, the order and knots in x-direction and/or y-direction can be specified using the keywords XORDER, XKNOTS, YORDER, and YKNOTS.

Separate discussions on one- and two-dimensional splines follow.

One-dimensional B-splines

Given the data points x = xdata, f = fdata, and the number of elements (n) in xdata and fdata, the default action of IMSL_BSINTERP computes a cubic (k = 4) spline interpolant s to the data using the default knot sequence generated by the IMSL_BSKNOTS function.

Optional keyword XORDER allows the user to choose the order, k, of the spline interpolant; optional keyword XKNOTS allows user specification of knots.

The IMSL_BSINTERP function is based on the routine SPLINT by de Boor (1978, p. 204).

First, IMSL_BSINTERP sorts the xdata vector and stores the result in x. The elements of the fdata vector are permuted appropriately and stored in f, yielding the equivalent data (xi, fi) for i = 0 to n – 1.

The following preliminary checks are performed on the data:

xi < xi + 1       i = 0, ..., n – 2

ti < ti + k         i = 0, ..., n – 1

tt £ ti + 1              i = 0, ..., n + k – 2

The first test checks to see that the abscissas are distinct. The second and third inequalities verify that a valid knot sequence has been specified.

In order for the interpolation matrix to be nonsingular, tk – 1 £ xi £ tn is also checked for i = 0 to n – 1. This first inequality in the last check is necessary since the method used to generate the entries of the interpolation matrix requires that the k possibly nonzero B-splines at xi:

Bj – k + 1, ..., Bj              where j satisfies tj £ xi < tj + 1

be well-defined (that is, jk + 1 ³ 0).

General conditions are not known for the exact behavior of the error in spline interpolation; however, if t and x are selected properly and the data points arise from the values of a smooth (for example, Ck) function f, i.e., fi = f(xi), then the error behaves in a predictable fashion. The maximum absolute error satisfies:

where the following is true:

For more information on this, see de Boor (1978, Chapter 13) and the references therein. This function can be used in place of the IMSL_CSINTERP function.

The returned value for this function is a structure. This structure contains all the information to determine the spline (stored as a linear combination of B-splines) that is computed by this function. For example, the following code sequence evaluates this spline at x and returns the value in y:

y = IMSL_SPVALUE(x, spline)  

Two-dimensional, Tensor-product B-splines


If arguments xdata, ydata, and fdata are all included in the call to the IMSL_BSINTERP function, the function computes a two-dimensional, tensor-product spline interpolant. The tensor-product spline interpolant to data {(xiyjfij)}, where 0 £ i £ nx – 1 and 0 £ j £ ny – 1, has the form:

where kx and ky are the orders of the splines. These numbers are defaulted to 4 but can be set to any positive integer using keywords XORDER and YORDER. Likewise, tx and ty are the corresponding knot sequences (XKNOTS and YKNOTS). These values are defaulted to the knots returned by the IMSL_BSKNOTS function. The algorithm requires that the following is true:

tx (kx – 1) £ xi £ tx (nx)                 0 £ i £ nx – 1

ty (ky – 1) £ yj £ ty (ny)                 0 £ j £ ny – 1

Tensor-product spline interpolants in two dimensions can be computed quite efficiently by solving (repeatedly) two univariate interpolation problems. The computation is motivated by the following observations:

Setting:

note that for each fixed i from 0 to nx – 1, there are ny linear equations in the same number of unknowns as can be seen below:

The same matrix appears in the previous equation.

Thus, this matrix is factored only once, then the factorization to solve the nx right-hand sides is applied. Once this is done and hmi is computed, then the coefficients cnm are solved using the relation:

for m from 0 to ny – 1, which involves one factorization and ny solutions to the different right-hand sides. This ability of the IMSL_BSINTERP function is based on the SPLI2D routine by de Boor (1978, p. 347).

The returned value is a structure containing all the information to determine the spline (stored in B-spline format) that is computed by this function. For example, the following code sequence evaluates this spline at (xy) and returns the value in z:

z = IMSL_SPVALUE(x, y, spline)  

Examples


Example 1

In this example, a one-dimensional B-spline interpolant to function values is computed, as shown in Figure 6-5. Evaluations of the computed spline are then plotted along with the original data values. Since the default settings are being used, the interpolant is determined by the not-a-knot condition (see de Boor 1978).

x = FINDGEN(11)/10  
; Define data values.
f = SIN(15 * x)
bs = IMSL_BSINTERP(x, f)
; Compute interpolant.
bsval = IMSL_SPVALUE(FINDGEN(100)/99, bs)
PLOT, FINDGEN(100)/99, bsval
; Output results.
OPLOT, x, f, Psym = 6

 

Figure 6-5: B-Spline Interpolant

Example 2

In this example, a two-dimensional, tensor-product B-spline interpolant to gridded data is computed as shown in Figure 6-6.

x = FINDGEN(5)/4  
; Define the abscissas in the x-direction.
y = FINDGEN(5)/4
; Define the abscissas in the y-direction.
f = FLTARR(5, 5)
; Define the sample function values.
FOR i = 0, 4 DO $
   f(i, *) = SIN(2 * x(i)) - COS(5 * y)
bs = IMSL_BSINTERP(x, y, f)
; Compute the spline interpolant.
bsval = IMSL_SPVALUE(FINDGEN(20)/19, FINDGEN(20)/19, bs)
; Use IMSL_SPVALUE to evaluate the computed spline.
!P.Charsize = 1.5
!P.Multi = [0, 1, 2]
WINDOW, XSize = 400, YSize = 800
; Plot the original and computed surfaces in a tall window.
SURFACE, f, x, y
SURFACE, bsval, FINDGEN(20)/19, $
FINDGEN(20)/19

 

Figure 6-6: Two-Dimensional B-Spline

Errors


Warning Errors

MATH_ILL_COND_INTERP_PROB—Interpolation matrix is ill-conditioned. Solution might not be accurate.

Fatal Errors

MATH_DUPLICATE_XDATA_VALUES—The xdata values must be distinct.

MATH_YDATA_NOT_INCREASING—The ydata values must be strictly increasing.

MATH_KNOT_MULTIPLICITY—Multiplicity of the knots cannot exceed the order of the spline.

MATH_KNOT_NOT_INCREASING—Knots must be nondecreasing.

MATH_KNOT_XDATA_INTERLACING—The i-th smallest element of xdata (xi) must satisfy ti £ xi < ti + Order, where t is the knot sequence.

MATH_XDATA_TOO_LARGE—Array xdata must satisfy
xdatai £ tndata, for i = 1, ..., ndata.

MATH_XDATA_TOO_SMALL—Array xdata must satisfy
xdatai ³ tOrder – 1, for i = 1, ..., ndata.

MATH_KNOT_DATA_INTERLACING—The i-th smallest element of the data arrays xdata and ydata must satisfy ti £ datai + Order, where t is the knot sequence.

MATH_DATA_TOO_LARGE—Data arrays xdata and ydata must satisfy
datai £ tnum_data, for i = 1, ..., num_data.

MATH_DATA_TOO_SMALL—Data arrays xdata and ydata must satisfy
datai ³ tOrder – 1, for i = 1, ..., num_data.

Version History


6.4
Introduced



© 2017 Exelis Visual Information Solutions, Inc. |  Legal
My Account    |    Buy    |    Contact Us