Scilab Function
Last update : 24/03/2004

splin - cubic spline interpolation

### Calling Sequence

  `d = splin(x, y [,spline_type [, der]])`

### Parameters

- `x` : a strictly increasing (row or column) vector (x must have at least 2 components)
- `y` : a vector of same format than `x`
- `spline_type` : (optionnal) a string selecting the kind of spline to compute
- `der` : (optionnal) a vector with 2 components, with the end points derivatives (to provide when spline_type="clamped")
- `d` : vector of the same format than `x` (`di` is the derivative of the spline at `xi`)

### Description

This function computes a cubic spline or sub-spline s which interpolates the (xi,yi) points, ie, we have s(xi)=yi for all i=1,..,n. The resulting spline s is completly defined by the triplet `(x,y,d)` where `d` is the vector with the derivatives at the `xi`: s'(xi)=di (this is called the Hermite form). The evaluation of the spline at some points must be done by the `interp` function. Several kind of splines may be computed by selecting the appropriate `spline_type` parameter:

  - "not_a_knot": this is the default case, the cubic spline is computed by using the following conditions (considering n points x1,...,xn):

            s'''(x2-) = s'''(x2+)
            s'''(x{n-1}-) = s'''(x{n-1}+)

  - "clamped": in this case the cubic spline is computed by using the end points derivatives which must be provided as the last argument `der`:

            s'(x1) = der(1)
            s'(xn) = der(2)

  - "natural": the cubic spline is computed by using the conditions:

            s''(x1) = 0
            s''(xn) = 0

  - "periodic": a periodic cubic spline is computed (`y` must verify y1=yn) by using the conditions:

            s'(x1) = s'(xn)
            s''(x1) = s''(xn)

  - "monotone": in this case a sub-spline (s is only one continuously differentiable) is computed by using a local scheme for the di such that s is monotone on each interval:

            if y(i) <= y(i+1)  s is increasing on [x(i), x(i+1)]
            if y(i) >= y(i+1)  s is decreasing on [x(i), x(i+1)]

  - "fast": in this case a sub-spline is also computed by using a simple local scheme for the di : d(i) is the derivative at x(i) of the interpolation polynomial of (x(i-1),y(i-1)), (x(i),y(i)),(x(i+1),y(i+1)), except for the end points (d1 being computed from the 3 left most points and dn from the 3 right most points).
  - "fast_periodic": same as before but use also a centered formula for d1 = s'(x1) = dn = s'(xn) by using the periodicity of the underlying function (`y` must verify y1=yn).

### Remarks

From an accuracy point of view use essentially the clamped type if you know the end point derivatives, else use not_a_knot. But if the underlying approximated function is periodic use the periodic type. Under the good assumptions these kind of splines got an O(h^4) asymptotic behavior of the error. Don't use the natural type unless the underlying function have zero second end points derivatives.

The monotone, fast (or fast_periodic) type may be useful in some cases, for instance to limit oscillations (these kind of sub-splines have an O(h^3) asymptotic behavior of the error).

If n=2 (and `spline_type` is not clamped) linear interpolation is used. If n=3 and `spline_type` is not_a_knot, then a fast sub-spline type is in fact computed.

### Examples

    // example 1
    deff("y=runge(x)","y=1 ./(1 + x.^2)")
    a = -5; b = 5; n = 11; m = 400;
    x = linspace(a, b, n)';
    y = runge(x);
    d = splin(x, y);
    xx = linspace(a, b, m)';
    yyi = interp(xx, x, y, d);
    yye = runge(xx);
    xbasc()
    plot2d(xx, [yyi yye], style=[2 5], leg="interpolation spline@exact function")
    plot2d(x, y, -9)
    xtitle("interpolation of the Runge function")

    // example 2 : show behavior of different splines on random datas
    a = 0; b = 1;        // interval of interpolation
    n = 10;              // nb of interpolation  points
    m = 800;             // discretisation for evaluation
    x = linspace(a,b,n)'; // abscissae of interpolation points
    y = rand(x);          // ordinates of interpolation points
    xx = linspace(a,b,m)';
    yk = interp(xx, x, y, splin(x,y,"not_a_knot"));
    yf = interp(xx, x, y, splin(x,y,"fast"));
    ym = interp(xx, x, y, splin(x,y,"monotone"));
    xbasc()
    plot2d(xx, [yf ym yk], style=[5 2 3], strf="121", ...
           leg="fast@monotone@not a knot spline")
    plot2d(x,y,-9, strf="000")  // to show interpolation points
    xtitle("Various spline and sub-splines on random datas")
    xselect()

### See Also

`interp`, `lsq_splin`,

### Authors

  B. Pincon
  F. N. Fritsch (pchim.f Slatec routine is used for monotone interpolation)