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
|
include <mach.h>
include "hdicfit.h"
# HDIC_INIT -- Initialize hdfit/icgfit interface.
procedure hdic_init (density, nvalues, dmax)
double density[nvalues] # Reference density above fog values
int nvalues # Number of values in sample
double dmax # Maximum possible density
int i
pointer sp, base
double xxmax, xxmin, delta_den
bool fp_equald()
include "hdic.com"
errchk malloc
begin
call smark (sp)
if (den == NULL || big_den == NULL) {
call malloc (den, nvalues, TY_DOUBLE)
call malloc (big_den, NVALS_FIT, TY_DOUBLE)
} else if (nvalues != nraw)
call realloc (den, nvalues, TY_DOUBLE)
nraw = nvalues
call salloc (base, NVALS_FIT, TY_DOUBLE)
# Copy density array to pointer location
call amovd (density, Memd[den], nraw)
# Calculate big vector of density values. The points are spaced
# linear in log space, to yield adequate spacing at low density values.
call alimd (density, nraw, xxmin, xxmax)
# Put user value for maximum density in common block if it is valid.
if (! fp_equald (dmax, 0.0D0))
maxden = dmax
# Make big_den go all the way to maxden, not just xxmax. Make sure
# the value of xxmin won't cause the log function to blow up.
if (xxmin > 0.0D0)
;
else
xxmin = 2.0 * EPSILOND
delta_den = (log10 (maxden) - log10 (xxmin)) / double (NVALS_FIT - 1)
do i = 1, NVALS_FIT
Memd[big_den+i-1] = log10 (xxmin) + double (i-1) * delta_den
call amovkd (10.0D0, Memd[base], NVALS_FIT)
call aexpd (Memd[base], Memd[big_den], Memd[big_den], NVALS_FIT)
call sfree (sp)
end
|
