e01tg:: Interpolation (NAG Toolbox)

nag_interp_3d_scat_shep (e01tg) constructs a smooth function

Q (x, y, z)

which interpolates a set of

m

scattered data points

(x_{r}, y_{r}, z_{r}, f_{r})

, for

r = 1, 2, \dots, m

, using a modification of Shepard's method. The surface is continuous and has continuous first partial derivatives.

The basic Shepard method, which is a generalization of the two-dimensional method described in Shepard (1968), interpolates the input data with the weighted mean

Q (x, y, z) = \frac{\sum_{r = 1}^{m} w_{r} (x, y, z) q_{r}}{\sum_{r = 1}^{m} w_{r} (x, y, z)},

where

q_{r} = f_{r} ​ and ​ w_{r} (x, y, z) = \frac{1}{d_{r}^{2}} ​ and ​ d_{r}^{2} = {(x - x_{r})}^{2} + {(y - y_{r})}^{2} + {(z - z_{r})}^{2} .

The basic method is global in that the interpolated value at any point depends on all the data, but this function uses a modification (see Franke and Nielson (1980) and Renka (1988a)), whereby the method becomes local by adjusting each

w_{r} (x, y, z)

to be zero outside a sphere with centre

(x_{r}, y_{r}, z_{r})

and some radius

R_{w}

. Also, to improve the performance of the basic method, each

q_{r}

above is replaced by a function

q_{r} (x, y, z)

, which is a quadratic fitted by weighted least squares to data local to

(x_{r}, y_{r}, z_{r})

and forced to interpolate

(x_{r}, y_{r}, z_{r}, f_{r})

. In this context, a point

(x, y, z)

is defined to be local to another point if it lies within some distance

R_{q}

of it. Computation of these quadratics constitutes the main work done by this function.

The radii

R_{w}

and

R_{q}

are chosen to be just large enough to include

N_{w}

and

N_{q}

data points, respectively, for user-supplied constants

N_{w}

and

N_{q}

. Default values of these arguments are provided by the function, and advice on alternatives is given in Choice of and .

References

Parameters

Compulsory Input Parameters

Optional Input Parameters

Output Parameters

Error Indicators and Warnings

Accuracy

Further Comments

Timing

The time taken for a call to nag_interp_3d_scat_shep (e01tg) will depend in general on the distribution of the data points. If x, y and z are uniformly randomly distributed, then the time taken should be O(m). At worst

O (m^{2})

time will be required.

Choice of Nw and Nq

If non-default values are required they must be supplied to nag_interp_3d_scat_shep (e01tg) through positive values of nw and nq. Increasing these arguments makes the method less local. This may increase the accuracy of the resulting interpolant at the expense of increased computational cost. The default values

nw = \min (32, m - 1)

and

nq = \min (17, m - 1)

have been chosen on the basis of experimental results reported in Renka (1988a). In these experiments the error norm was found to vary smoothly with

N_{w}

and

N_{q}

, generally increasing monotonically and slowly with distance from the optimal pair. The method is not therefore thought to be particularly sensitive to the argument values. For further advice on the choice of these arguments see Renka (1988a).

Example

This program reads in a set of

30

data points and calls nag_interp_3d_scat_shep (e01tg) to construct an interpolating function

Q (x, y, z)

. It then calls nag_interp_3d_scat_shep_eval (e01th) to evaluate the interpolant at a set of points.

function e01tg_example


fprintf('e01tg example results\n\n');

data = [ 0.80  0.23  0.37  0.51;
         0.23  0.88  0.05  1.80;
         0.18  0.43  0.04  0.11;
         0.58  0.95  0.62  2.65;
         0.64  0.69  0.20  0.93;
         0.88  0.35  0.49  0.72;
         0.30  0.10  0.78 -0.11;
         0.87  0.09  0.05  0.67;
         0.04  0.02  0.40  0.00;
         0.62  0.90  0.43  2.20;
         0.87  0.96  0.24  3.17;
         0.62  0.64  0.45  0.74;
         0.86  0.13  0.47  0.64;
         0.87  0.60  0.46  1.07;
         0.49  0.43  0.13  0.22;
         0.12  0.61  0.00  0.41;
         0.02  0.71  0.82  0.58;
         0.62  0.93  0.44  2.48;
         0.49  0.54  0.04  0.37;
         0.36  0.56  0.39  0.35;
         0.62  0.42  0.97 -0.20;
         0.01  0.72  0.45  0.78;
         0.41  0.36  0.52  0.11;
         0.17  0.99  0.65  2.82;
         0.51  0.29  0.59  0.14;
         0.85  0.05  0.04  0.61;
         0.20  0.20  0.87 -0.25;
         0.04  0.67  0.04  0.59;
         0.31  0.63  0.18  0.50;
         0.88  0.27  0.07  0.71];

x = data(:,1); y = data(:,2); z = data(:,3); f = data(:,4);
nw = int64(0);
nq = nw;
[iq, rq, ifail] = e01tg( ...
                         x, y, z, f, nw, nq);

% Evaluate at equispaced points on diagonal line
px = [0.1:0.1:0.6]'; py = px; pz = px;

[q, qx, qy, qz, ifail] = e01th( ...
                                x, y, z, f, iq, rq, px, py, pz);

fprintf('  Evaluation point       Q(x,y,z)\n');
fprintf(' (%4.1f, %4.1f, %4.1f)    %8.4f\n',[px py pz q]');

e01tg example results

  Evaluation point       Q(x,y,z)
 ( 0.1,  0.1,  0.1)      0.2630
 ( 0.2,  0.2,  0.2)      0.1182
 ( 0.3,  0.3,  0.3)      0.0811
 ( 0.4,  0.4,  0.4)      0.1552
 ( 0.5,  0.5,  0.5)      0.3019
 ( 0.6,  0.6,  0.6)      0.5712

On entry,	$m < 10$ ,
or	$0 < nq < 9$ ,
or	$nq > \min (40, m - 1)$ ,
or	$nw > \min (40, m - 1)$ ,
or	$liq < 2 \times m + 1$ ,
or	$lrq < 10 \times m + 7$ .

NAG Toolbox: nag_interp_3d_scat_shep (e01tg)

▸▿ Contents

Purpose

Syntax

Description