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
|
include <mach.h>
include "hdicfit.h"
# HDIC_GVEC -- Get a vector of data to be plotted. The raw density
# values are stored in common. From these values, all possible transforms
# can be generated.
procedure hdic_gvec (ic, xout, npts, transform)
pointer ic # Pointer to ic structure
double xout[npts] # Output vector
int npts # Desired npoints to output - possibly changed on output
int transform # Integer code for desired type of transform
pointer sp, bigdenaf
int i, j, npts_desired
double fog, dval
real ic_getr()
include "hdic.com"
begin
npts_desired = npts
if (npts_desired != NVALS_FIT)
call error (0, "hdicgvec: nvals != NVALS_FIT")
call smark (sp)
call salloc (bigdenaf, npts, TY_DOUBLE)
j = 1
fog = double (ic_getr (ic, "fog"))
call asubkd (Memd[big_den], fog, Memd[bigdenaf], npts)
switch (transform) {
case HD_NONE:
call amovd (Memd[bigdenaf], xout, NVALS_FIT)
case HD_LOGO:
do i = 1, npts_desired {
dval = Memd[bigdenaf+i-1]
if (dval < 0.0D0)
npts = npts - 1
else {
xout[j] = log10 ((10.**dval) - 1.0)
j = j + 1
}
}
case HD_K75:
do i = 1, npts_desired {
dval = Memd[bigdenaf+i-1]
if (dval < 0.0D0)
npts = npts - 1
else {
xout[j] = dval + 0.75 * log10 (1.0 - (10.** (-dval)))
j = j + 1
}
}
case HD_K50:
do i = 1, npts_desired {
dval = Memd[bigdenaf+i-1]
if (dval < 0.0D0)
npts = npts - 1
else {
xout[j] = dval + 0.50 * log10 (1.0 - (10.** (-dval)))
j = j + 1
}
}
}
call sfree (sp)
end
|