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
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
|
include <tbset.h>
# tu_getput -- do the interpolation
# This routine reads the independent and dependent variable values,
# does the interpolation, and writes the results to the output table.
#
# Phil Hodge, 24-April-2000 subroutine created
procedure tu_getput (itp, otp, iv_icp, iv_ocp,
icp, ocp, ncols,
xout, outrows, iv_step,
i_func, extrapolate, ext_value, padvalue, array, verbose)
pointer itp, otp # i: descr of input & output tables
pointer iv_icp # i: column descr for input indep var column
pointer iv_ocp # i: column descr for output indep var column
pointer icp[ncols] # i: column descriptors for input table
pointer ocp[ncols] # i: column descriptors for output table
int ncols # i: number of columns to copy
double xout[ARB] # i: output indep var values
int outrows # i: size of xout
double iv_step # i: increment in independent variable
int i_func # i: interpolation function code
bool extrapolate # i: true means don't use ext_value
double ext_value # i: value to assign to extrapolated points
double padvalue # i: value at end of input indep. var. to ignore
bool array # i: true if indep var column contains arrays
bool verbose # i: print info?
#--
pointer sp
pointer xin # scratch for input indep var values
pointer yin # scratch for input data values
pointer y2 # scratch for second derivatives (spline only)
pointer xa # scratch for non-indef indep var values
pointer ya # scratch for non-indef dep var values
pointer yout # array of interpolated values
double dbuf # one interpolated value
int inrows # number of rows in input table
int xnelem # number of elements in xin, INDEF & padvalue trimmed
int nelem # allocated number of elements in array
int nvals # number of array elements actually gotten
int n # size of xa, ya, y2 (counting only non-INDEF values)
int row, col # loop indices
int i
int tbpsta(), tbcigi(), tbagtd()
errchk tugcol, tbagtd, tbaptd, tbegtd, tbeptd, tbrcsc
begin
call smark (sp)
inrows = tbpsta (itp, TBL_NROWS)
if (array)
nelem = tbcigi (iv_icp, TBL_COL_LENDATA)
else
nelem = inrows
call salloc (xin, nelem, TY_DOUBLE)
call salloc (yin, nelem, TY_DOUBLE)
call salloc (y2, nelem, TY_DOUBLE)
call salloc (xa, nelem, TY_DOUBLE)
call salloc (ya, nelem, TY_DOUBLE)
if (array) {
call salloc (yout, outrows, TY_DOUBLE)
# for each row in the input table ...
do row = 1, inrows {
# Get input independent variable array from current row.
call tugcol (itp, iv_icp, row, Memd[xin], xnelem,
padvalue, array)
# Put output indep var values into current row.
call tbaptd (otp, iv_ocp, row, xout, 1, outrows)
# for each column to be interpolated ...
do col = 1, ncols {
nelem = tbcigi (icp[col], TBL_COL_LENDATA)
if (nelem == 1) {
# just copy scalar column
call tbrcsc (itp, otp, icp[col], ocp[col], row, row, 1)
} else {
nvals = tbagtd (itp, icp[col], row, Memd[yin], 1, nelem)
call tuifit (i_func, Memd[xin], Memd[yin], xnelem,
Memd[xa], Memd[ya], Memd[y2], n)
if (n > 0) {
do i = 1, outrows {
# interpolate
call tuival (i_func,
Memd[xa], Memd[ya], Memd[y2], n,
extrapolate, ext_value, xout,
i, outrows, iv_step, Memd[yout+i-1])
}
# put the result
call tbaptd (otp, ocp[col], row,
Memd[yout], 1, outrows)
} else if (verbose && row == 1) {
call printf ("not enough data for interpolation\n")
call flush (STDOUT)
}
}
}
}
} else {
# Get input independent variable column.
row = 1
call tugcol (itp, iv_icp, row, Memd[xin], xnelem, padvalue, array)
# Put output independent variable values into output table column.
do row = 1, outrows
call tbeptd (otp, iv_ocp, row, xout[row])
# for each column to be interpolated ...
do col = 1, ncols {
do row = 1, inrows
call tbegtd (itp, icp[col], row, Memd[yin+row-1])
# xnelem will be less than or equal to inrows.
call tuifit (i_func, Memd[xin], Memd[yin], xnelem,
Memd[xa], Memd[ya], Memd[y2], n)
if (n > 0) {
do row = 1, outrows {
# interpolate, and put the result
call tuival (i_func, Memd[xa], Memd[ya], Memd[y2], n,
extrapolate, ext_value,
xout, row, outrows, iv_step, dbuf)
call tbeptd (otp, ocp[col], row, dbuf)
}
} else if (verbose) {
call printf ("not enough data for interpolation\n")
call flush (STDOUT)
}
}
}
call sfree (sp)
end
|
