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
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
|
include <mach.h>
include <gset.h>
include <pkg/gtools.h>
include "hdicfit.h"
define MSIZE 2. # Mark size
define SZ_TPLOT 7 # String length for plot type
define EB_WTS 10
define EB_SDEV 11
define MAX_SZMARKER 4
# ICG_GRAPH -- Graph data and fit.
procedure icg_graphd (ic, gp, gt, cv, x, y, wts, ebw, npts)
pointer ic # ICFIT pointer
pointer gp # GIO pointer
pointer gt # GTOOLS pointers
pointer cv # Curfit pointer
double x[npts] # Independent variable
double y[npts] # Dependent variable
double wts[npts] # Weights
double ebw[npts] # Half the error bar width
int npts # Number of points
pointer xout, yout
int nvalues
errchk malloc
begin
call malloc (xout, npts, TY_DOUBLE)
call malloc (yout, npts, TY_DOUBLE)
call icg_axesd (ic, gt, cv, 1, x, y, Memd[xout], npts)
call icg_axesd (ic, gt, cv, 2, x, y, Memd[yout], npts)
call icg_paramsd (ic, cv, x, y, wts, npts, gt)
call icg_g1d (ic, gp, gt, Memd[xout], Memd[yout], wts, ebw, npts)
if (npts != IC_NFIT(ic)) {
if ((abs (IC_NAVERAGE(ic)) > 1) || (IC_NREJECT(ic) > 0)) {
call realloc (xout, IC_NFIT(ic), TY_DOUBLE)
call realloc (yout, IC_NFIT(ic), TY_DOUBLE)
call icg_axesd (ic, gt, cv, 1, Memd[IC_XFIT(ic)],
Memd[IC_YFIT(ic)], Memd[xout], IC_NFIT(ic))
call icg_axesd (ic, gt, cv, 2, Memd[IC_XFIT(ic)],
Memd[IC_YFIT(ic)], Memd[yout], IC_NFIT(ic))
call icg_g2d (ic, gp, gt, Memd[xout], Memd[yout],
IC_NFIT(ic))
}
} else if (IC_NREJECT(ic) > 0)
call icg_g2d (ic, gp, gt, Memd[xout], Memd[yout], npts)
nvalues = NVALS_FIT
call icg_gfd (ic, gp, gt, cv, nvalues)
# Mark the the sample regions.
call icg_sampled (ic, gp, gt, x, npts, 1)
call mfree (xout, TY_DOUBLE)
call mfree (yout, TY_DOUBLE)
end
# ICG_G1D -- Plot data vector.
procedure icg_g1d (ic, gp, gt, x, y, wts, ebw, npts)
pointer ic # Pointer to ic structure
pointer gp # Pointer to graphics stream
pointer gt # Pointer to gtools structure
double x[npts] # Array of independent variables
double y[npts] # Array of dependent variables
double wts[npts] # Array of weights
double ebw[npts] # Error bar half width in WCS (positive density)
int npts # Number of points
pointer sp, xr, yr, sz, szmk, gt1, gt2
char tplot[SZ_TPLOT]
int i, xaxis, yaxis, symbols[3], markplot
real size, rmin, rmax, big
bool fp_equald(), streq(), fp_equalr()
int strncmp()
data symbols/GM_PLUS, GM_BOX, GM_CIRCLE/
include "hdic.com"
begin
call smark (sp)
call salloc (xr, npts, TY_REAL)
call salloc (yr, npts, TY_REAL)
call salloc (sz, npts, TY_REAL)
call salloc (szmk, npts, TY_REAL)
call achtdr (x, Memr[xr], npts)
call achtdr (y, Memr[yr], npts)
xaxis = (IC_AXES (ic, IC_GKEY(ic), 1))
yaxis = (IC_AXES (ic, IC_GKEY(ic), 2))
# Set up gtools structure for deleted (gt1) and added (gt2) points
call gt_copy (gt, gt1)
call gt_setr (gt1, GTXSIZE, MSIZE)
call gt_setr (gt1, GTYSIZE, MSIZE)
call gt_sets (gt1, GTMARK, "cross")
call gt_copy (gt, gt2)
call gt_setr (gt2, GTXSIZE, MSIZE)
call gt_setr (gt2, GTYSIZE, MSIZE)
call gt_sets (gt2, GTMARK, "circle")
markplot = NO
call gt_gets (gt, GTTYPE, tplot, SZ_TPLOT)
if (strncmp (tplot, "mark", 4) == 0)
markplot = YES
if (IC_OVERPLOT(ic) == NO) {
# Start a new plot
call gclear (gp)
# Set the graph scale and axes
call gascale (gp, Memr[xr], npts, 1)
call gascale (gp, Memr[yr], npts, 2)
# If plotting HD curve, set wy2 to maxden, which may have
# been updated if a new endpoint was added.
if (xaxis == 'y' && yaxis == 'u')
call gswind (gp, INDEF, INDEF, INDEF, real (maxden))
call gt_swind (gp, gt)
call gt_labax (gp, gt)
}
# Calculate size of markers if error bars are being used. If the
# weights are being used as the marker size, they are first scaled
# between 1.0 and the maximum marker size.
call gt_gets (gt, GTMARK, tplot, SZ_TPLOT)
if (streq (tplot, "hebar") || streq (tplot, "vebar")) {
if (IC_EBARS(ic) == EB_WTS) {
call achtdr (wts, Memr[sz], npts)
call alimr (Memr[sz], npts, rmin, rmax)
if (fp_equalr (rmin, rmax))
call amovr (Memr[sz], Memr[szmk], npts)
else {
big = real (MAX_SZMARKER)
call amapr (Memr[sz], Memr[szmk], npts,rmin,rmax,1.0, big)
}
} else {
call achtdr (ebw, Memr[sz], npts)
call hd_szmk (gp, gt, xaxis, yaxis, tplot, Memr[sz],
Memr[szmk], npts)
}
} else
call amovkr (MSIZE, Memr[szmk], npts)
do i = 1, npts {
# Check for deleted point
if (fp_equald (wts[i], 0.0D0))
call gt_plot (gp, gt1, Memr[xr+i-1], Memr[yr+i-1], 1)
# Check for added point
else if (fp_equald (ebw[i], ADDED_PT))
call gt_plot (gp, gt2, Memr[xr+i-1], Memr[yr+i-1], 1)
else {
size = Memr[szmk+i-1]
call gt_setr (gt, GTXSIZE, size)
call gt_setr (gt, GTYSIZE, size)
call gt_plot (gp, gt, Memr[xr+i-1], Memr[yr+i-1], 1)
}
}
IC_OVERPLOT(ic) = NO
call sfree (sp)
end
# HD_SZMK -- Calculate size of error bar markers. This procedure is
# called when the marker type is hebar or vebar and the marker size is
# to be the error in the point, not its weight in the fit.
procedure hd_szmk (gp, gt, xaxis, yaxis, mark, insz, outsz, npts)
pointer gp # Pointer to gio structure
pointer gt # Pointer to gtools structure
int xaxis, yaxis # Codes for x and y axis types
char mark[SZ_TPLOT] # Type of marker to use
real insz[npts] # Standard deviations in WCS units
real outsz[npts] # Output size array in NDC units
int npts # Number of points in arrays
int gt_geti()
char tplot[SZ_TPLOT]
real dx, dy
bool streq()
begin
# Check validity of axis types
if (xaxis != 'x' && xaxis != 'u' && yaxis != 'x' && yaxis != 'u') {
call eprintf ("Choose graph type with axes 'u' or 'x'; ")
call eprintf ("Using marker size 2.0\n")
call amovkr (MSIZE, outsz, npts)
return
}
call gt_gets (gt, GTMARK, tplot, SZ_TPLOT)
if (streq (tplot, "hebar")) {
if (yaxis == 'x' || yaxis == 'u') {
call gt_sets (gt, GTMARK, "vebar")
call eprintf ("Marker switched to vebar\n")
# call flush (STDOUT)
}
} else if (streq (tplot, "vebar")) {
if (xaxis == 'x' || xaxis == 'u') {
call gt_sets (gt, GTMARK, "hebar")
call eprintf ("Marker switched to hebar\n")
# call flush (STDOUT)
}
}
# Need to scale standard deviation from density to NDC units
call ggscale (gp, 0.0, 0.0, dx, dy)
if (gt_geti (gt, GTTRANSPOSE) == NO)
call adivkr (insz, dx, outsz, npts)
else
call adivkr (insz, dy, outsz, npts)
end
# ICG_G2D -- Show sample range and rejected data on plot.
procedure icg_g2d (ic, gp, gt, x, y, npts)
pointer ic # ICFIT pointer
pointer gp # GIO pointer
pointer gt # GTOOLS pointer
double x[npts], y[npts] # Data points
int npts # Number of data points
int i
pointer sp, xr, yr, gt1
begin
call smark (sp)
call salloc (xr, npts, TY_REAL)
call salloc (yr, npts, TY_REAL)
call achtdr (x, Memr[xr], npts)
call achtdr (y, Memr[yr], npts)
call gt_copy (gt, gt1)
call gt_sets (gt1, GTTYPE, "mark")
# Mark the sample points.
if (abs (IC_NAVERAGE(ic)) > 1) {
call gt_sets (gt1, GTMARK, "plus")
call gt_setr (gt1, GTXSIZE, real (-abs (IC_NAVERAGE(ic))))
call gt_setr (gt1, GTYSIZE, 1.)
call gt_plot (gp, gt1, Memr[xr], Memr[yr], npts)
}
# Mark the rejected points.
if (IC_NREJECT(ic) > 0) {
call gt_sets (gt1, GTMARK, "diamond")
call gt_setr (gt1, GTXSIZE, MSIZE)
call gt_setr (gt1, GTYSIZE, MSIZE)
do i = 1, npts {
if (Memi[IC_REJPTS(ic)+i-1] == YES)
call gt_plot (gp, gt1, Memr[xr+i-1], Memr[yr+i-1], 1)
}
}
call gt_free (gt1)
call sfree (sp)
end
# ICG_GF[RD] -- Overplot the fit on a graph of the data points. A vector
# is constructed and evaluated, then plotted.
procedure icg_gfd (ic, gp, gt, cv, npts)
pointer ic # ICFIT pointer
pointer gp # GIO pointer
pointer gt # GTOOL pointer
pointer cv # CURFIT pointer
int npts # Number of points to plot
pointer sp, xr, yr, x, y, xo, yo, gt1
int gstati()
begin
call smark (sp)
if (IC_FITERROR(ic) == YES)
return
call salloc (xr, npts, TY_REAL)
call salloc (yr, npts, TY_REAL)
call salloc (x, npts, TY_DOUBLE)
call salloc (y, npts, TY_DOUBLE)
call salloc (xo, npts, TY_DOUBLE)
call salloc (yo, npts, TY_DOUBLE)
# Calculate big vector of independent variable values. Note value
# of npts can change with this call.
call hdic_gvec (ic, Memd[x], npts, IC_TRANSFORM(ic))
# Calculate vector of fit values.
call dcvvector (cv, Memd[x], Memd[y], npts)
# Convert to user function or transpose axes. Change type to reals
# for plotting.
call icg_axesd (ic, gt, cv, 1, Memd[x], Memd[y], Memd[xo], npts)
call icg_axesd (ic, gt, cv, 2, Memd[x], Memd[y], Memd[yo], npts)
call achtdr (Memd[xo], Memr[xr], npts)
call achtdr (Memd[yo], Memr[yr], npts)
call gt_copy (gt, gt1)
call gt_sets (gt1, GTTYPE, "line")
call gt_seti (gt1, GTLINE, gstati (gp, G_PLTYPE))
call gt_plot (gp, gt1, Memr[xr], Memr[yr], npts)
call gt_free (gt1)
call sfree (sp)
end
|