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
|
include <math/nlfit.h>
include <pkg/inlfit.h>
# ING_ERRORS -- Compute error diagnostic information and print it on the
# screen.
procedure ing_errors$t (in, file, nl, x, y, wts, npts, nvars)
pointer in # INLFIT pointer
char file[ARB] # Output file name
pointer nl # NLFIT pointer
PIXEL x[ARB] # Ordinates (npts * nvars)
PIXEL y[ARB] # Abscissas
PIXEL wts[ARB] # Weights
int npts # Number of data points
int nvars # Number of variables
bool isfit
int i, j, deleted, rejected, nparams, fd
PIXEL chisqr, variance, rms
pointer sp, fit, wts1, params, errors, rejpts, plist
pointer name, pvnames, labels
int open(), nlstati(), inlstrwrd(), in_geti()
pointer in_getp()
PIXEL in_rms$t(), nlstat$t()
errchk open()
begin
# Open the output file.
if (file[1] == EOS)
return
fd = open (file, APPEND, TEXT_FILE)
# Determine the number of coefficients.
nparams = nlstati (nl, NLNPARAMS)
# Allocate memory for parameters, errors, and parameter list.
call smark (sp)
call salloc (params, nparams, TY_PIXEL)
call salloc (errors, nparams, TY_PIXEL)
call salloc (labels, SZ_LINE + 1, TY_CHAR)
# Allocate memory for the fit and strings.
call salloc (fit, npts, TY_PIXEL)
call salloc (wts1, npts, TY_PIXEL)
call salloc (name, SZ_LINE + 1, TY_CHAR)
call salloc (pvnames, SZ_LINE + 1, TY_CHAR)
# Get number of rejected points and rejected point list.
rejected = in_geti (in, INLNREJPTS)
rejpts = in_getp (in, INLREJPTS)
# Count deleted points.
deleted = 0
do i = 1, npts {
if (wts[i] == PIXEL (0.0))
deleted = deleted + 1
}
# Assign a zero weight to rejected points.
call amov$t (wts, Mem$t[wts1], npts)
if (rejected > 0) {
do i = 1, npts {
if (Memi[rejpts+i-1] == YES)
Mem$t[wts1+i-1] = PIXEL (0.0)
}
}
# Get the parameter values and errors.
call nlvector$t (nl, x, Mem$t[fit], npts, nvars)
call nlpget$t (nl, Mem$t[params], nparams)
call nlerrors$t (nl, y, Mem$t[fit], Mem$t[wts1], npts, variance,
chisqr, Mem$t[errors])
# Compute the RMS.
rms = in_rms$t (y, Mem$t[fit], Mem$t[wts1], npts)
# Print the error analysis.
call fprintf (fd, "\nniterations %d\n")
call pargi (nlstati (nl, NLITER))
call fprintf (fd, "total_points %d\n")
call pargi (npts)
call fprintf (fd, "rejected %d\n")
call pargi (in_geti (in, INLNREJPTS))
call fprintf (fd, "deleted %d\n")
call pargi (deleted)
call fprintf (fd, "standard deviation %10.7g\n")
call parg$t (sqrt (variance))
call fprintf (fd, "reduced chi %10.7g\n")
call parg$t (sqrt (chisqr))
call fprintf (fd, "average error %10.7g\n")
if (chisqr <= 0)
call parg$t (PIXEL(0.0))
else
call parg$t (sqrt (max (variance, PIXEL (0.0)) / chisqr))
call fprintf (fd, "average scatter %10.7g\n")
call parg$t (sqrt (nlstat$t (nl, NLSCATTER)))
call fprintf (fd, "RMS %10.7g\n")
call parg$t (rms)
# Print parameter values and errors.
call in_gstr (in, INLPLABELS, Memc[labels], SZ_LINE)
call strcpy (Memc[labels], Memc[pvnames], SZ_LINE)
call fprintf (fd, "\n%-10.10s %14.14s %14.14s\n")
call pargstr ("parameter")
call pargstr ("value")
call pargstr ("error")
plist = in_getp (in, INLPLIST)
do i = 1, nparams {
if (inlstrwrd (i, Memc[name], SZ_LINE, Memc[pvnames]) != 0) {
call fprintf (fd, "%-10.10s ")
call pargstr (Memc[name])
} else {
call fprintf (fd, "%-10.2d ")
call pargi (i)
}
call fprintf (fd, "%14.7f %14.7f (%s)\n")
call parg$t (Mem$t[params+i-1])
call parg$t (Mem$t[errors+i-1])
isfit = false
do j = 1, nparams {
if (Memi[plist+j-1] == i) {
isfit = true
break
}
}
if (isfit)
call pargstr ("fit")
else
call pargstr ("constant")
}
call fprintf (fd, "\n")
# Free allocated memory.
call sfree (sp)
call close (fd)
end
|
