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
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
|
.help surfit Aug84 "Math Package"
.ih
NAME
surfit -- set of routines for fitting gridded surface data
.ih
SYNOPSIS
The surface is defined in the following way. The x and y values are
assumed to be integers and lie in the range 0 <= x <= ncols + 1 and
0 <= y <= nlines + 1. The package prefix is for image surface fit.
Package routines with a prefix followed by l operate on individual
image lines only.
.nf
isinit (sf, surf_type, xorder, yorder, xterms, ncols, nlines)
iszero (sf)
islzero (sf, lineno)
islaccum (sf, cols, lineno, z, w, ncols, wtflag)
islsolve (sf, lineno, ier)
islfit (sf, cols, lineno, z, w, ncols, wtflag, ier)
islrefit (sf, cols, lineno, z, w, ier)
issolve (sf, lines, nlines, ier)
isresolve (sf, lines, ier)
z = iseval (sf, x, y)
isvector (sf, x, y, zfit, npts)
issave (sf, surface)
isrestore (sf, surface)
isfree (sf)
.fi
.ih
DESCRIPTION
The SURFIT package provides a set of routines for fitting surfaces to functions
linear in their coefficients using the tensor product method and least squares
techniques. The basic numerical
technique employed is the solution of the normal equations by the Cholesky
method.
.ih
NOTES
The following series of steps illustrates the use of the package.
.ls 4
.ls (1)
Include an include <math/surfit.h> statement in the calling program to make the
SURFIT package definitions available to the user program.
.le
.ls (2)
Call ISINIT to initialize the surface fitting parameters.
.le
.ls (3)
Call ISLACCUM to select a weighting function and accumulate data points
for each line into the appropriate arrays and vectors. Call ISLSOLVE
to compute the coefficients in x for each line. The coefficients are
stored inside SURFIT. Repeat this procedure for each line. ISLACCUM
and SFLSOLVE can be combined by a call to SFLFIT. If the x values
and weights remain the same from line to line ISLREFIT can be called.
.le
.ls (4)
Call ISSOLVE to solve for the surface coefficients.
.le
.ls (5)
Call ISEVAL or ISVECTOR to evaluate the fitted surface at the x and y
value(s) of interest.
.le
.ls (6)
Call ISFREE to release the space allocated for the fit.
.le
.le
.ih
EXAMPLES
.nf
Example 1: Fit a 2nd order Lengendre polynomial in x and y to an image
and output the fitted image
include <imhdr.h>
include <math/surfit.h>
old = immap (old_image, READ_ONLY, 0)
new = immap (new_image, NEW_COPY, 0)
ncols = IM_LEN(old, 1)
nlines = IM_LEN(old, 2)
# initialize surface fit
call isinit (sf, LEGENDRE, 3, 3, YES, ncols, nlines)
# allocate space for lines, columns and weight arrays
call malloc (cols, ncols, TY_INT)
call malloc (lines, nlines, TY_INT)
call malloc (weight, ncols, TY_REAL)
# initialize lines and columns arrays
do i = 1, ncols
Memi[cols - 1 + i] = i
do i = 1, nlines
Memi[lines - 1 + i] = i
# fit the surface in x line by line
call amovkl (long(1), v, IM_MAXDIM)
do i = 1, nlines {
if (imgnlr (old, inbuf, v) == EOF)
call error (0, "Error reading image.")
if (i == 1) {
call islfit (sf, Memi[cols], i, Memr[inbuf], Memr[weight],
ncols, SF_WTSUNIFORM, ier)
if (ier != OK)
...
} else
call islrefit (sf, Memi[cols], i, Memr[inbuf], Memr[weight],
ier)
}
# solve for surface coefficients
call issolve (sf, Memi[lines], nlines, ier)
# free space used in fitting arrays
call mfree (cols, TY_INT)
call mfree (lines, TY_INT)
call mfree (weight, TY_REAL)
# allocate space for x and y arrays
call malloc (x, ncols, TY_REAL)
call malloc (y, ncols, TY_REAL)
# intialize z array
do i = 1, ncols
Memr[x - 1 + i] = real (i)
# create fitted image
call amovkl (long(10, v, IM_MAXDIM)
do i = 1, nlines {
if (impnlr (new, outbuf, v) == EOF)
call error (0, "Error writing image.")
call amovkr (real (i), Memr[y], ncols)
call isvector (sf, Memr[x], Memr[y], Memr[outbuf], ncols)
}
# close files and cleanup
call mfree (x, TY_REAL)
call mfree (y, TY_REAL
call isfree (sf)
call imunmap (old)
call imunmap (new)
.fi
.endhelp
|
