blob: 030d1964f7ebe0edfc03c15130db081069b375a3 (
plain) (
blame)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
|
# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
.help
Procedure iipol_terp
A polynomial interpolator with x and y arrays given.
This algorithm is based on the Newton form as described in
de Boor's book, A Practical Guide to Splines, 1978.
There is no error checking - this is meant to be used only by calls
from more complete routines that take care of such things.
Maximum number of terms is 6.
.endhelp
real procedure iipol_terp(x,y,n,x0)
real x[ARB],y[ARB] # x and y array
real x0 # desired x
int n # number of points in x and y = number of
# terms in polynomial = order + 1
int k,i
real d[6]
begin
# Fill in entries for divided difference table.
do i = 1,n
d[i] = y[i]
# Generate divided differences
do k = 1,n-1
do i = 1,n-k
d[i] = (d[i+1] - d[i])/(x[i+k] - x[i])
# Shift divided difference table to center on x0
do i = 2,n
d[i] = d[i] + d[i-1] * (x0 - x[i])
return(d[n])
end
|
