June 1990
Agronomy and Soils
Departmental Series No. 145
Alabama Agricultural Experiment Station
Lowell T. Frobish, Director
Auburn University

CALIBRATION AND AUTOMATION OF A DUAL-ENERGY GAMMA
SYSTEM FOR APPLICATIONS IN SOIL SCIENCE
M. Oostrom and J.H. Dane
Graduate Research Assistant and Professor
of Agronomy and Soils
Agronomy and Soils Departmental Series No. 145
June 1990
Alabama Agricultural Experiment Station
Auburn University
Auburn University, Alabama
Lowell T. Frobish, Director
ii
ACKNOWLEDGMENTS
This report was, in part, prepared for and funded by the
Electric Power Research Institute (EPRI) Land and Water
Quality Studies Program, Research Project 2938-1 (Evaluation
of Mixing and Plume Penetration Depth at the Groundwater
Table) and the United States Environmental Protection Agency
(USEPA) Research Project CR-814946-01-0 (Multiphase Chemical
Transport in Porous Media).
Principal investigators for the EPRI project are O. Guven
(Civil Engineering Dept., Auburn University), J.H. Dane
(Agronomy and Soils Dept., Auburn University), W.E. Hill
(Chemistry Department, Auburn University), and J.G. Melville
(Civil Engineering Dept., Auburn University). This project is
administered through the Water Resources Research Institute,
Auburn University, which is directed by J. Judkins. The EPRI
project manager is D.A. McIntosh, EPRI, Palo Alto, California.
Principal investigators for the USEPA project are G.F.
Pinder (College of Engineering, University of Vermont), P.R.
Jaffe and M.A. Celia (Dept. of Civil Engineering and
Operations Research, Princeton University), and J.H. Dane
(Agronomy and Soils Dept., Auburn University). This project is
administered through the Office of Research & Project
Administration, Princeton .. University, which is directed by
A.J. Sinisgalli. The USEPA project manager is J. Weaver,
USEPA, Ada, Oklahoma. This report does not necessarily reflect
iii
the views of the United States Environmental Protection Agency
and no official endorsement should be inferred.
iv
CONTENTS
ACKNOWLEDGMENTS ..................... iii
LIST OF FIGURES . . . . . . . . . . . . . .*.*.*.e....vil
LIST OF TABLES . . . . . . . . . . . . . . ....... ix
ABSTRACT .. . . .. . . . . . . . . . . . ......... x
INTRODUCTION .. . .. . . . . . . . . . . ........ 1
THEORY OF GAMMA RADIATION ATTENUATION . ........ 4
Simultaneous Determination of Dry Bulk Density
and Volumetric Water Content . . ............ 6
Simultaneous Determination of Volumetric NAPL Content
and Volumetric Water Content . . . . . . . . ... 8
Simultaneous Determination of the Volumetric Water
Content and the Salt Concentration .......... 9
GAMMA SYSTEM INSTRUMENTATION .............. 11
THE CALIBRATION PROCESS E.......... .. .*00 16
Determination of 137Cs and.
241
Am Spectra . . . . . . . 16
Dead Time Corrections. . . . . . . .. . . . . . . . 19
Correction for Compton Scattered 
137
Cs Photons
in the 
2 4 1
Am Window. *....... ... .... .25
Determination of Attenuation Coefficients . . . . . . 29
Volumetric Attenuation Coefficients of Liquids . 29
Mass Attenuation Coefficients of Soils . . . . . 35
Path Length Determination ............ ... 37
Determination of the Incident Count Rates 10 . . . . 38
ERROR ANALYSIS . . . . . . . . . . . . . . . . . .. . 40
Error Analysis for the Simultaneous Determination
of the Dry Soil Bulk Density, p,, and the Volumetric
Water Content, GO ........ .. ...... 41
V
Error Analysis for the Simultaneous Determination
of the Volumetric Water Content, O , and
the Volumetric NAPL Content, On . . ..........
Error Analysis for the Simultaneous Determination of the
Volumetric Water Content, GO, and the Salt
Concentration,C ............. . .. .
LITERATURE CITED . . . . . . . . . . . . . . . . . .
50
56
63
APPENDIX A. Asyst Software. . . . . . . . . . . . . . . . 66
APPENDIX B. Listings of Calibration Programs.. . . . . .. 76
FIRST PRINTING, JUNE 1990
Information contained herein is available
to all without regard to race, color, sex,
or national origin.
vi
LIST OF FIGURES
FIG. 1. Function diagram of gamma system components . . 12
FIG. 2. Energy spectrum of
137
Cs ... . . . . . 18
FIG. 3. The 241Am spectra . . . . . . . . . . . . . . . 20
FIG. 4. Relation between the integral 241Am count rates
of the NC-22 and the TC-450 single channel
analyzers . . . . . . . . . . . . . . . . . . . 26
FIG. 5. Compton scattering in the 
241
Am window versus
the high energy 1Cs intensity. All count rates are
dead time corrected . . . . . . . . . . . . . . 30
FIG. 6. Relation between path length through a solution
and the natural log of the corrected count rates
or deionized water and a 0.2 n Nal solution. The
absolute value of the slope is the volumetric
liquid attenuation coefficient . . . . . . . . 32
FIG. 7. Experimentally determined relations between the
NaI concentration of a solution and the
volumetric attenuation coefficient for 241Am and
137
Cs. Each data point is obtained with
calibration program 4 ........ ..... 34
FIG. 8. Relation between path length through TCE and the
natural log of the corrected count rates. The
absolute value of the slope is the TCE
volumetric attenuation coefficient . . . . . . 35
FIG. 9A. Relative standard deviation in p. and O versus
path length x for ps=1.6 g/cm
3
, =O0.3, and tm=
6
0
sec . . . . . . . . . . . . . . . . . . . . 47
FIG. 9B. Relative standard deviation in p. and 0 versus
path length x for p,=1.6 g/cm
3
, 9,=0.1, and t =60
sec . . . . . . . . . . . . . . . . . . . . . 47
FIG. 9C. Relative standard deviation in p. and 8, versus
path length x for p,=1l.2 g/cm
3
, Wm,=0.3, and tm=60
sec . . . . . .. . . . . . . . . . .. ... 48
FIG. 9D. Relative standard deviation in p and 0 versus
path length x for p
5
=1..6 g/cm, =0.1, Wand...
t,=300 sec . . . . . . . . . . . . . . . . . . 48
vii
10A. Relative standard deviation 
in Ow and 8 versus
path length x for Ow=O*l5, Gn=Oel5l 
ana tm=60
FIG. lOB.
FIG. 10C.
Relative standard deviation 
in Ow and On versus
path length x for Gw=Oe25, Gn=O.O5, 
and tm
6
O
sec
Relative standard deviation 
in Ow and On versus
path length x for OW=O.O5, G 
=O.25, and tm=6O
sec 0. 0. 0. 0. 0. 0. 0 0 0 0 0*0*0* 
0 0 0.. . . . . . . .
FIG.0D. Relative standard deviation 
in wand versus
path length x for gw=0.b15, 
U=Oe b5, ana tm=
3 0 0
54
55
55
FIG. hlA.
FIG. biB.
Relative standard deviation in OU 
and C versus
path length x for Uw=O.35, C=O9Ob 
meq/ml, and
tm=
6
Osec . e.9.9 .. 0.0 0 0. 
0 . 0. 0
Relative standard deviation in 0 and 
C versus
path length x for W=.03 5, C=O91 
meq/ml, and
t =60 sec . . . . . . . . . . * .
. . .
FIG. 11C. Relative standard deviation 
in Ow and C versus
-path length x for Uw=Oeb, C=0*b 
meq/ml, and
tn=6 0 sec a . . . 0. 0. 0. 0. . a .
. 0.0.0.0.
60
60
61
FIG. biD. Relative standard deviation 
in 0 and C versus
path length x for ,=.35, C=O.01 
meq/ml, and
t.=30 sec .*. . 0 0 . * 0 0 . .
. * . . *
0vai0
FIG.
54
61
LIST OF TABLES
Table 1'.
Table 2.
Table 3.
Table 4.
Dead Times (jssec) of SCA's and NA-22 * 0 0 . . 0 23
Expressions for Partial Derivatives aOIaXi
and a0
1
,/ax, 0. 0 0. 0. . 0. 0. 0 0. 0 0 0 0 0 0 .0 . 0. 0. 42
Expressions for Partial Derivatives aOjaX
and 8O.I,/8X 0.0.0. * 0 0 0 .0.0.0.........*....51
Expressions for Partial Derivatives aO.ax
1
and W4&, 0 . 0. 0. 9. 0. 0. 0. 0'.0 .* .0 .0 .0 .0 .0 .0 . 57
ix
ABSTRACT
Dual-energy gamma ray attenuation methods have been
traditionally used in soil science to determine water content
and bulk density values. Recently, interest in those
techniques revived with two new applications: (1) 
the
simultaneous determination of the volumetric non-aquous phase
liquid content and water content values, and (2) the
simultaneous determination of volumetric water content and
salt concentration values. Based on Beer's law, the
attenuation equations for all three applications were derived.
The dual-energy gamma system of the soil physics
laboratory at Auburn University has been recently adapted to
improve the calibration process and the data acquisition.
The calibration process consists of the determination of: (1)
the 
241
Am and 
137
Cs spectra, (2) the dead times of the single
channel analyzers, 
(3) the Compton scattering 
in the 
241
Am
window, (4) the attenuation coefficients of the liquids and
soils, (5) the path lengths through the porous medium, and (6).
the incident count rates. The programs used in the
calibration process are written in ASYST, a language
especially designed for data acquisition and analysis.
A probable error analysis, based on the error propagation
formula, is presented ?for the three applications. This
analysis gives an indication, of the attainable precision and
is a tool to determine the optimum path length through a
porous medium.
x
INTRODUCTION
Gamma ray attenuation methods were initially used in soil
science to determine water content or bulk density values.
The monoenergetic or single gamma ray attenuation method was
used by Van Bavel et al. (29) to determine attenuation
coefficients of soils and subsequently their dry bulk density
values. Assuming a constant dry bulk density, several
researchers (4, 12, 3, 23, 6) were able to determine water
content values. The assumption of a constant dry bulk density
implied that the monoenergetic gamma ray technique to
determine volumetric water contents could only be used for
non-swelling soils. 
Most investigators 
used the isotope 
137
Cs,
which emits gamma radiation photons at 662 KeV. King (15)
reported promising results measuring soil water content values
with the 241Am isotope, which emits gamma radiation primarily
at 60 KeV.
Soane (24) and Gardner and Calissendorff (8) were the
first to use a dual-energy gamma radiation system to measure
dry bulk density and volumetric water content simultaneously.
They used 
137
Cs and 241Am as the gamma radiation sources.
Theoretically, any two gamma sources with different gamma
energies, such that their mass attenuation coefficients are
sufficiently different, can be used to do this. Gardner and
Calissendorff (8) and Corey et al. (2) Showed that important
factors such as half-life of the isotopes, cost, shielding,
and self-absorption, make 
137
Cs and 241m the most logical*
1
2
choices for dual-energy gamma measurements.
Soane (24) placed his samples between a 
137
Cs source and
the detector and then between an 24Am source and the detector.
He did not combine the sources in a single collimator, because
the higher energy 
137
Cs gamma photons, through interaction with
matter, can produce lower energy photons (Compton scattering)
which may interfere with the 241Am gamma rays. The derivation
of the equations used to calculate the dry bulk density and
the volumetric water content is based on the simultaneous-
determination of 
1
37Cs and 241Am count rates. The use of these
equations, when using Soane's experimental setup, can only be
justified if the soil-water system is changing slowly.
Stroosnijder and de Swart (27) used two perpendicular beams
with two separate detectors to avoid difficulties with 
241
Am
interference. This configuration restricts the geometry of
soil-water systems that can be studied. Corey et al. (2) and
Mansell et al. (21) used a common collimator for the two
sources and a single detector. Corey et al. (2) used
multichannel analyzers to correct for the Compton scattering.
At least four measurements with different multichannel
analyzer settings were necessary to correct for the Compton
interference. This made his technique less suitable for
rapidly changing soil water contents. Nofziger and
Swartzendruber (22) expressed the correction due to Compton
scattering as a function of the 
137
Cs count rate. With the
24 1
Am removed and only the 
137
Cs source remaining, they found,
3
by varying the attenuation, that a third degree polynomial
expressed the Compton scattering the best. This correction
method was later improved by Stillwater and 
Klute (26).
Recently, the use of dual-energy gamma radiation
techniques revived with two new applications. Ferrand et al.
(5) and Lenhard et al. (18) experimented with the
nondestructive simultaneous determination of the volumetric
nonaqueous phase liquid (NAPL) and the volumetric water
content, while Grismer et al. (11) determined the volumetric
water content and the salt concentration in soils with low
water contents. In all these applications it was assumed that
the dry bulk density did not change with time. Ferrand et al.
(5) reported good results with halogenated hydrocarbons such
as trichloroethylene (TCE), tetrachloride and chloroform. The
relatively heavy Cl atoms caused the attenuation coefficients
of .the halogenated hydrocarbons to differ substantially from
that of water. Grismer et al. (11) experienced a linearly
increasing attenuation coefficient with salt concentration.
The best results were obtained with sodium iodide (NaI), lead
nitrate (Pb(N0
3
)
2
), and strontium chloride (SrCl
2
). The heavy
ions in these salts (I', Pb
2+
, Sr
2
)-influenced the attenuation
coefficient.
In the next section of this departmental series, a brief
overview is given of the theory of gamma attenuation. The
equations are derived for the simultaneous determination of
1. volumetric water content and dry bulk density,
4
2. volumetric 
nonaqueous phase 
liquid content 
and
volumetric water content, and
3. salt concentration and volumetric water content.
Subsequently, the instrumentation of the dual-energy gamma
radiation system of the soil physics laboratory at Auburn
University is discussed, followed by a description of the
calibration process. Finally, an error analysis for each of
the three applications is given. In Appendix A, a brief
overview of the computer language ASYST, used to automate the
calibration and data acquisition, is presented, followed by
listings of the programs used in the calibration process in
Appendix B.'
THEORY OF GAMMA RADIATION ATTENUATION
Gamma rays are electromagnetic photons which travel with
the speed of light and are uncharged and unaffected by
electrical or magnetic fields. When parent radionuclides
decay, gamma radiation is emitted in their transition to lower
lying levels. In most cases, the gamma radiation follows a or
B decay. The two sources used in the soil physics laboratory
at Auburn University are 
137
Cs and 141Am. The 
137
Cs source emits
gamma radiation primarily at 662 KeV following IB decay, while
the 241m gamma rays, emitted after a decay, have an energy of
60 KeV (17).
The three most important interactions of attenuation
processes of gamma rays with matter are photoelectric
5
absorption, Compton scattering and pair production (16). In
the photoelectric absorption process, the gamma ray collides
with an atom. The energy of the photon disappears completely
while a photoelectron is produced along with one or more x-ray
photons. For gamma rays with relatively low energy, this mode
of interaction is the most important. The photoelectric
absorption depends strongly on the atomic number Z of the
absorbing material. The higher the Z-value, the higher the
probability of interaction. In the process of Compton
scattering, an incoming photon with a certain energy collides
with an electron and is deflected. The deflected gamma ray
has an energy that is lower than the energy of the incoming
photon. Part of the energy of the incoming photon is
transferred to the electron, known as the recoil electron. The
third process of interaction is pair production. For this
process of combining an electron and a positron, a minimum
energy of 1020 KeV is required, which is impossible when using
137
Cs and 
241
Am as gamma sources. The total attenuation of 
137
Cs
and 
241
Am is, therefore, the sum of the attenuation due to
Compton scattering and photoelectric absorption. Below 200
KeV the photoelectric absorption is the predominant effect,
while higher energy gamma rays are more affected by Compton
scattering.
The attenuation of gamma rays can be described by Beers'
law:
6
I-I
oe x p  
Pc iP ix 
(1)
where I
o 
is the incident count rate (counts per second (cps)),
I is the emergent count rate (cps), I is the mass attenuation
coefficient of component i (L
2 
M-
1
) , x
i 
is the path length of
the gamma beam through i (L), and pi is the density of i
(ML'
3
). This equation forms the starting point for the
derivation of the equations for the simultaneous determination
of two components with a dual-energy gamma radiation system.
Simultaneous Determination of Dry Bulk Density
and Volumetric Water Content
The attenuation of gamma radiation through a soil with
dry bulk density p (ML-
3
) and volumetric water content O can
be written as
I- I(-ASPgX -A Pwx -P gPgogX) 
(2)
where I
0 
is the count rate through the empty column (cps), pi
is the density of phase i (ML'
3
), i is the mass attenuation
coefficient of phase i (L
2
M'
1
), and 0 is volumetric gas phase
content. The subscripts s, w, and g denote soil, water, and
gas phase, respectively. During the calibration process, the
attenuation coefficient and density of a liquid are not
7
determined separately. Instead, a value U (L
I
) representing
the product of the mass attenuation coefficient and the
density of the liquid, jp,, is obtained. Replacing the product
pw by the volumetric attenuation coefficient, Uw, and
ignoring the contribution of the gas phase to the attenuation
of the photons, eq. (2) becomes
I I- oexp(-A,px - Uww
x
) (3)
Introducing the subscripts c and a for the 
137
Cs and 
241
Am gamma
radiation sources, respectively, gives
I, " Ioaexp(- saPx 
x 
-Uwa X) 
(4)
I
C 
- Ioexp(-js
0
p x -UwcOx) (5)
Eq. (4) and (5) can be solved for ps and O to yield
P .Uw 1 Uwal()(6)
S a( al{) P~))
- ~(X PsaUwc - p scUwa) (7)
8
Simultaneous Determination of Volumetric NAPL Content
and Volumetric Water Content
The gamma ray attenuation through a. porous medium with
water (w), NAPL (n), and gas (g) can be written as
I - Ioexp(-p,pSx 
-pP,,8wx -nPnX? 
PgP OgX) 
(8 )
where I
o
' is the count rate through the empty column (cps).
The NAPL has attenuation coefficient n (L
2
M
"
'
) 
, density pn
(ML'
3
) and volumetric NAPL content On. Assuming the
contribution of the gas phase to be negligible and replacing
A~F and An n by U, and Un, respectively, eq. (8) can 
be written
as
I - Ioexp(-p,px -Uwx -Un6x) (9)
When the bulk density of the porous medium is constant in
time, eq. (9) becomes
I - IOexp(-Uw0wx -UnOnx) (10)
where I
o 
= I
0
o'exp(-spsx). For a dual-energy gamma radiation
system, with 
137
Cs (subscript c) and 
241
Am (subscript a) as the.
sources, the two attenuation equations are
I
a
-Ioaexp(-U,6 wx -UnaOX) (11)
I
c
- oexp(-Uwcjx -UncOnx)
(12)
9
Solving eq. (11) and (12) for 8, and 8, yields
Unal - Uncl(13)
S x(wcUna - UUnc)
i 
I
C - II o Ul l0
u- l I - Un ) 
(14)
x(UcUn - UaUnc)
Simultaneous Determination of the Volumetric Water Content
and the Salt Concentration.
The attenuation of gamma radiation through a porous
medium with a salt solution (subscript L) can be written as
I - Ioexp(-p.px -PL L
x 
-gPg ggx) (15)
Assuming the gas attenuation to be negligible and the dry bulk
density to be constant in time, eq. (15) becomes
I - Ioexp(- p LLx) (16)
where I. = 1
0
'exp(-spx). Replacing -LPL by UL gives
I- Ioexp(-UOLx) (17)
The attenuation equations for the sources 
241
Am and 1
37
Cs are,
respectively,
I
a
- IoaexP(-U aO) 
(18)
10
IC - Iocexp(-ULALX) (19)
When the salt solution consists of a single salt dissolved 
in
water, its volumetric attenuation coefficient can be written
as
ULa waCw + raCr (20)
ULc- IwcCw + P rcCr (21)
where C, and Cr are the water and salt concentrations 
(ML'),
respectively, and Mra and src are the mass attenuation
coefficients of the pure salt. At relatively low salt
concentrations, it can be assumed that 
9
L=
0
w and and that the
product CjC
w 
is a constant. ULj is then a linear function of
the salt concentration only. Replacing Cr by C, .waC. by a (L'
'), jwC by r (L'
1
), Ira by b (LM'
"1
) and Mrc by s (L
2
W), eq.
(20) and (21) become, respectively,
U- a + bC (22)
ULC - r + sC (23)
Grismer et al. (11) experimentally determined linear relations
between the volumetric attenuation coefficients and the salt
molality of various salts. Assuming UL=Ow, and substituting
eq. (22) in (18) and (23) in (19) yields
Ia M Ioaexp (-[a + bC]Owx) (24)
Ic Iocexp(-[r + sC]Ow
x
) (25)
Solving eq. (24) and (25) for C and Ow results in
11
sl I)- bl {)
ia 
IC
sl - bl 
(27)
-w x(as - br)
GAMMA SYSTEM INSTRUMENTATION
The dual-energy gamma 
system of the soil physics
laboratory has recently 
been adapted to improve 
the
calibration process and the 
data acquisition and, therefore,
differs from the system 
previously described by 
Hopmans and
Dane (13, 14). The instrumentation 
of the dual-energy gamma
system is presented in a 
block diagram in figure 1.
The two sources, a 200 mCi 
1
37
Cs- and a 200 mCi 
41
Am-
source, are mounted coaxially 
in the center of a 0.17-m-long
and 0.14-m-diameter lead shield. 
The Cs-source is located in
the center, while the Am-source 
is placed at the front end 
of
the lead shield and is contained 
in a brass holder, 10 mm
thick and 40 mm in diameter.
-- 
A thallium-activated sodium
iodide (NaI(TI)) scintillation 
crystal with a photomultiplier
tube form the gamma ray detector 
(Harshaw Chemical, Solon,
Ohio). The detector is enclosed 
in another cylindrical lead
shield with a length of 0.11 
m and a diameter of 0.14 m.
Collimation is provided by 
6.35-mm cylindrical holes in 
the
STEP MOTOR
CONTROL
ONTAL STEP MOTOR ASTON 800 
METRABYTE TC 589
MOTOR CONTROL INTERFACE 
MBC-488 INTERFACE
ZENITH
COMPUTER
ASYST
SOFTWARE
FIG. 1. Function diagram of gamma system components.
oCs .Am
13
lead shields. Thallium-activated sodium iodide is widely used
in the scintillation spectroscopy of gamma rays because of the
high Z-value of iodide and the high density of the Nal-
crystal. The interaction of gamma radiation and matter is a
strong function of the Z-value and density of the material.
The principle of this detection method is to create visible
light pulses upon interaction of the incident gamma photons
and the NaI(Tl)-crystal, followed by a conversion of the light
pulses into measurable electrical pulses. The thallium is
added to the NaI-crystal to increase the visible photon
emission. In a pure NaI-crystal, electrons have only two
discrete bands of energy available, the lower valence band and
the higher conduction band. The intermediate band of
energies, called the forbidden band, has never been shown to
contain electrons. A gamma photon can excite an electron from
the valence band to the conduction band. After this, the
electron has a tendency to return to the valence band, but
this process is not very effective and the photon emitted when
returning to the valence band has an energy that is too high
to lie in the visible range. The impurity thallium, however,
creates energy sites in the forbidden gap preferred by the
excited electrons. When these electrons drop back to the
valence band, the likelihood of visible photons being 
emitted
has strongly increased. The half-life or decay constant of 
the
excited state of thallium electrons is .23 js, which is 
much
shorter than that of the pure Nal excited state.
14
A photomultiplier (PM) is necessary to convert the very
weak light output of a scintillation pulse into a
corresponding electrical signal. The two major elements of a
PM are the photocathode and the electron multiplier. The
photocathode converts the 
incident light photons 
into low-
energy photoelectrons. The electron multiplier collects the
photoelectrons by focusing electrodes and multiplies them by
a factor 105 - 108 by means of a system of dynodes. This
process, known as charge amplification, is usually enough to
create a charge signal to be collected at the anode of the
electron multiplier.
A stable voltage source NV 25A (Harshaw) is required for
an efficient charge amplification in the PM. The first dynode
of the electron multiplier has to be positive with respect to
the photocathode to attract the electrons, and every
succeeding dynode must have a higher voltage than the previous
dynode. The charges Q, by the detector, are so small that
they have to be amplified. -The amplification is performed by
a preamplifier, NB-28 (Harshaw), and an amplifier, NA-17
(Harshaw). The preamplifier is located as close as possible
to the detector to keep the noise, occurring in
interconnecting cables, as low as possible. The gain of the
amplifier is adjustable through a combination of a course and
fine gain control.
Once amplified, the electrical pulses are linearly
transmitted through an automatic gain control (AGC) unit, NA-
15
22 (Harshaw) to two single channel analyzers (SCA), the NC-22
(Harshaw) and the TC 450 (Tennelec, Oak Ridge, Tennessee).
The AGC-unit is provided with a 'limit' display, which will
turn on when the product of input amplitude and the gain set
on the NA-17 exceeds the maximum amplitude detectable. The
NA-22 also has a built in SCA. A SCA converts the shaped
linear pulses into logic, countable, pulses by using its
integral or differential discriminator. The integral
discriminator produces a logic output pulse if the linear
input pulse exceeds a set discrimination level. No output
appears when the linear input pulse amplitude is below the
discrimination level. A differential discriminator produces
a logic output if the linear input pulse amplitude lies
between the lower level discriminator (LLD) and the upper
level discriminator (ULD). The range between the LLD and ULD
is known as the SCA "window." In the calibration process,
discussed in the next section, the settings of the TC 450,
NC-22 and NA-22 are determined.
The accumulation of the counts over a fixed time period
is recorded by the TC 535P timer/multiscaler (Tennelec). The
timer/multiscaler consists of four counters, three of them are
used to count the logic output pulses of the SCA's and the
NA-22 and the fourth is usecto &isplay the elapsed time of a
measurement. The timer/multiscaler is connected with a MBC
IEEE-488 data acquisition board (MetraByte Corporation,
Taunton, Massachusetts) via a TC 489 buffered interface
16
(Tennelec). The computer needs the data acquisition board to
exchange General Purpose Interface Bus (GPIB) commands with
GPIB devices such 'as the TC 489 buffered interface. The
computer can instruct the timer/multiscaler to start and stop
counting via the TC 489 buffered interface. Data can be
stored temporarily in the buffer of the TC 489. The software
package ASYST (ASYST Software Technologies Inc., Rochester,
New Jersey) is especially designed to communicate with GPIB
devices via data acquisition boards. GPIB interfacing is
discussed in more detail in Appendix A.
The interconnected sources and detector are moved
simultaneously, either horizontally or vertically, by stepper
motors. Measurements can be taken in a 2-dimensional area of
2.10 x 1.05 m. The SPI55A horizontal and vertical stepper
motor control units (Superior Electric, Bristol, Connecticut)
can be activated manually or by computer software. When the
movement of the sources and detector is computer controlled,
the Aston 800 actuator (Aston Company, Atlanta, Georgia)
serves as an intermediate GPIB device, interfacing the stepper
motor control units and the data acquisition board in the
computer.
THE CALIBRATION PROCESS
Determination of 
137
Cs and 
24 1
Am~ Spectra
The amplified output pulses have a broad energy spectrum.
To discriminate these pulses, the energy spectra of both
17
sources have to be known. First, the 
137
Cs spectrum is
determined using the pulse height analyzer part of the NA-22
.and calibration program 1, spectrum.txt (APPENDIX B).
Depending on the high voltage supply of the NV 25A and the
coarse gain setting of the amplifier NA-17, the ratio of the
NA-22 base line voltage to the energy of the processed pulses
is a certain constant. Knowing that 
3
7Cs emits gamma
radiation primarily of .66 MeV and that the E range of the NA-
22 is 0-10 V, it is convenient to manipulate the high voltage
setting and the coarse gain in such a way that the NA-22 base
line voltage to 
137
Cs gamma radiation energy (in MeV) ratio is
10:1. With the NA-22 in differential mode and an E setting of
6.5 V and AE of 0.2 V, a high voltage of 706.0 V and a coarse
gain of 80 result in the highest count rate and the desired
voltage to energy ratio. The high voltage should not be
disturbed, because the photon detection relies on a stable
voltage.
To obtain the 
137
Cs spectrum, the toggle switch of the NA-
22 should be in the disabled mode. The initial setting of E
is 3.9 V and AE is 0.2 V. The output mode should be
differential so that the LLD of the window is 3.9 V and the
ULD is 4.1 V. For subsequent measurements, E is usually
increased by 0.1 V. A 
137
Cs spectrum obtained with calibration
program 1, spectrum.txt, is presented in figure 2.
Based on this graph, the LLD, corresponding to E, is 5.5 V and
the ULD, corresponding to E + AE, is 7.5 V. This results in
18
2500
11/11/89
analyzer NA-22
S1500 
-
O
500.
4.00 4.50 5.00 5.50 6.00 6.50 7.00 7.50 
8.00
Voltage 
(V)
FIG. 2. Energy spectrum of 
137
Cs.
E and AE settings of 5.5 and 2.0 V, respectively. After the
137
Cs spectrum is determined, 
the toggle switch of the 
NA-22 is
set to track to prevent a shift of the Cs-window. Since the
241
Am count rates are obtained using both the TC 450 and the
NC-22 SCA, two 
241
Am energy spectra have to be determined. An
initial E of 0.05 V is recommended while 0.1 V is a convenient
window width (AE). For each subsequent measurement, E has to
be increased by 0.03 or 0.04 V.
Operation of the TC 450 is more complicated than the
NC-22, which is basically the same as the analyzer part of
19
the NA-22. 
During the determination 
of the TC 
450 
241
Am
spectrum, it is convenient to use the differential mode which
corresponds with the 'window/asym' setting on the TC 450. The
upper knob then represents 10*AE, while E can be manipulated
with the lower knob. Once the LLD and the ULD of the spectrum
are obtained, the operation is switched to the 'norm' mode,
where the lower knob represents the LLD and the upper knob
represents the ULD of the spectrum. The 241Am spectra
determined with calibration program 1 are shown in figure 3.
For the NC-22, an E setting of 0.42 V and a AE of 0.42 V were
obtained, while the TC 450 has a LLD of 0.47 and an ULD of
0.93 V. The 241Am spectra contain counts due to Compton
scattering of 
137
Cs photons; a method to correct for this
interference will be discussed in a later section.
Dead Time Corrections
The dead time or resolving time of a counting system is
defined as the minimum amount of time which must separate two
events so that they are recorded as two separate pulses
(16). When an event occurs in the detector, the system needs
a certain amount of time,. 1 (sec/count), to process the
incoming gamma photon. All events occurring during a dead
time r after a recorded event will be lost. In our system, a
lost event does not further extend the dead time by another
period r. The system is said to have a nonparalysable
response as opposed to a system with a paralysable response
20
7000 -
11/14/89
6000-
analyzer NC-22
5000-
0 4000-
4J)
j 3000-
0
- 2000-
1000 -
sooo
.100 .300 .500 .700 
.900
Voltage (V)
20000 -
t11/15/9
analyzer TC 450
15000 -
0
4 10000 -
4.)
M
0
S 5000 -
0.
.10 .300 
.500 .700 
.900
Voltage 
(V)
FIG. 3. The 
241
Am spectra
21
where any lost event during the dead time period extends this
period. When IM is the measured count rate (counts 
per
second, cps), the prodct IMT is the fraction per second 
the
system is dead. The true count rate I can then be calculated
with
I - Im(1 Im') 
(28)
In most counting systems, SCA's are responsible for the
dead time, since these are the only components with discrete
processing periods. Most often, the SCA dead time is 
unknown
or an approximate value is provided by the manufacturer. For
correction purposes accurate dead time values are needed.
Before Fritton (6) described a method to determine dead times,
most researchers did not correct their measured count rates
for system dead time. Only Gurr (12) mentioned the use 
of a
measured dead time of 9 gsec for his system, but he did not
outline his dead time determination procedure. Fritton varied
the path length through the materials water, soil, lead, and
Plexiglas and corrected the measured count rates according to
eq. (28) for dead times varying from 0.1 to 7.0 Msec. 
He
assumed that the linearized equation describing the
attenuation process of a single material was
ln (I) - -ppx + ln (I
0
) (29)
where I and I
o
are the corrected count rates through a
22
material with thickness x and 0, respectively. After plotting
ln(I) versus x for the different dead times and calculating
the slopes with linear regression, the best value of 
the
resolving time minimizes the weighted variance 
(w.v.). The
expression of this weighted variance is
n
w.v. - Ay
2
/I. 
(30)
i-1
where n is the number of points used for the line and Ay are
the deviations from the calculated regression line. For all
the materials used in his study, Fritton found a dead time 
of
5.0 Msec. He concluded that, when a correction for dead time
is made, the mass absorption coefficients are much closer to
the theoretical values. Nofziger and Swartzendruber (22) 
used
the same method, except they minimized the weighted sum 
of
variances in I compared to the weighted sum of variances 
in
In(I) used by Fritton.
In another effort to determine dead time, Stroosnijder
and de Swart (27) measured count rates with and without 
an
absorber for three different distances between source 
and
detector. The ratios between the corrected (eq. 28) count
rates with and without the absorber were calculated and
plotted against dead time. The three ratio-dead time curves
intersected in one point where the three ratios and the dead
time are the same. This dead time was then assumed the system
23
dead time. Hopmans and Dane (13, 14) modified 
this method by
using 3 different absorber thicknesses 
instead of three
source-detector distances.
A more precise, direct way to determine 
the dead time of
a SCA was presented by Stillwater and Klute 
(26). They used
a square wave generator to simulate 
input pulses to the SCA
and a dual channel oscilloscope 
to view input and output
pulses.. At low square wave frequencies, 
the ratio of input to
output pulses is 1:1, which means that 
all the input pulses
are processed, and at high frequencies 
this ratio is 2:1,
which means that only 50% of the 
incoming pulses are
processed. When starting at a low frequency 
and increasing
this frequency slowly, it should 
be possible to find the
transition point or a narrow. transition 
band where the input
to output ratio changes from 1:1 to 2:1. 
The reciprocal of
the transition frequency is then the 
dead time of the SCA.
This method was used 
to determine the dead 
time of our SCA's
and the SCA part of the NA-22. The 
results are shown in table
1.
Table 1. Dead Times
(Msec) of SCA's and 
NA-22.
Device 
Dead Time-
TC 450 0.60
NC-22 1.67
NA-22 1.25
Hopmans and Dane (13, 14) reported indirectly 
determined dead
24
times of 2.5 Asec for both the NA-22 and NC-22.
Stillwater and Klute (26) discussed the validity of eq.
(28), which has been used by most researchers to correct for
dead time. They argued that, when a SCA is operated in the
integral mode, all pulses exceeding the LLD are processed and
counted, while a SCA in differential mode also processes all
pulses exceeding the LLD, but only counts those pulses between
the preset LLD and ULD. The fraction of a second the SCA is
dead is therefore determined by the count rate of the
processed pulses (I ) and not by the count rate of the
measured pulses (Im). To obtain the dead time corrected count
rate, I, the measured count rate has to be divided by the
fraction of a second, -I 7, the system is not dead, i.e.,
1= I /(l-I P) (31)
The 137Cs radiation count rates measured in differential and
integral mode with the NA-22 are practically equal, because
only a few incoming photons have energies greater than the ULD
set. This means that eq. (28) can still be used for 
137
Cs
radiation count rates. For 
241
Am, however, the integral count
rate, which equals Ip, is usually much larger than the
measured rate in differential mode, I*, ...and eq. ...(31i) must be
used to correct for dead, time. Replacing Im by the
differential 241m count rate I and Ip by the integral count
rate Iat' eq. (31) becomes
25
Ia = Iad/(l-IaiTa) (32)
where I, and Ta are the dead time corrected count rate and dead
time of the 
241
Am SCA, respectively.
Stillwater and Klute (26) proposed the use of a SCA which
is able to supply both the integral and differential counts
simultaneously to two different scalers. With our system, we
use a different setup to obtain the same result. Since none
of our SCA's is able to produce both 'ai and I , we have to
use two SCA's supplying Iai and Iad. Since the TC 450 SCA has
the smallest ra, it was chosen to provide I. and Ta
. 
To use
eq. (32), however, the Iai produced by the TC 450 has to be
known also. An indirect way to obtain this Iai is to determine
the relation between the Iai of the NC-22 and the TC 450 during
calibration of the system. This relation can then be used to
transform a measured integral count rate supplied by the NC-22
into a calculated Iai of the TC 450. Calibration program 2,
rel.txt (APPENDIX B), is used to determine this relation,
which appears to be highly linear (fig. 4.). The data points
shown in figure 4 are obtained by moving the sources and
detector with respect to a teflon prism and taking subsequent
measurements.
Correction for Compton Scattered 
137
Cs Photons
in the 241Am Window.
Due to Compton scattering, 
13
Cs photons release part of
26
11/15/89
40
30
4 y - -99.3139 + 
1.0311 x
0 
r- .999952
Q4)
0
0 10 20 30 40
Count rate via the NC-22 (cps, thousands)
FIG. 4. Relation between the integral 241Am 
count rates of the
NC-22 and the TC-450 single channel analyzers.
their energy. When one detector is used to simultaneously
determine 
241
Ar and 137Cs photons, part of the Compton scattered
137
Cs photons cause interference in the 
241
Am window set on the
TC 450 SCA. The production of these lower energy 
137
Cs photons
originates in the NaI(TI) crystal and not in the 
material
placed between the sources and detector (22). 
Compton
scattering will cause the photon to deflect from its 
original
path while giving up some of its energy. These
characteristics of Compton scattering imply that the
27
probability for the attenuation material to reduce the energy
of a 13
7
Cs photon by 90%, without changing its path, is very
small. According to Nofziger 
and Schwartzendruber (22), 
a 
137
Cs
photon needs four collisions at an angle of 1800 to have its
energy reduced from 662 KeV to 60 KeV and to remain in the
collimated radiation beam. The results of experiments they
performed with different materials (glass, granular quartz,
water, brass, Salkum soil, Ruskell soil) showed the material
independence of the low 137Cs energy photon count rate in the
241Am window. Loch and Kay (19) 
also reported this independence
of the Compton scattering on the attenuation materials
Plexiglas, brass, aluminum, wood, lead, and water.
The Compton interference of the low energy 
137
Cs photon
count rate in the 
241
Am window (Ic, Low) as a function of the high
energy 
137
Cs count rate (Ic) is described by either a cubic
polynomial (22, 13, 19),
Iclow - a + bI
c 
+ cI + dIC (33)
or by a linear polynomial (30, 28, 26),
Ic,low- a + bI
c  
(34)
Goit et al. (10) did not provide an equation for the nonlinear
relation they observed. All the relations described by cubic
polynomials show a slight curvature at higher count rates and
the coefficients c and d are very small. For example, Hopmans
28
and Dane (13) reported a value of -2.315935*10
5 
sec/count for
c and 5.5831407*10
.9 
sec
2
/count
2 
for d. Swartzendruber et al.
(28) showed that a reduction of the high energy 
137
Cs count rate
range results in a linearization of the Compton scattering
relation.
Theoretically, since the gamma system geometry determines
the probability a photon will cause a pulse of a given energy,
the ratio of the low energy 
137
Cs photon count rate in the 
241
Am
window to the high energy 
137
Cs photon count rate is a constant
(26). The relation between the two count rates should be a
line passing through the origin. Stillwater and Klute (26)
attributed the nonlinearity at high count rates found by some
researchers to the lack of dead time correction of the count
rates. It has been common practice to do the Compton
correction before the dead time correction, because Compton-
Compton-corrected count rates were needed in the indirect
methods to determine the dead time. The direct dead time
determination proposed by Stillwater and Klute made it
possible to do the dead time correction before the Compton
correction, resulting in a more accurate determination of the
true 
241
Am count rate. In order to obtain the true 
241
Am count
rate (Ia,tru), the dead time corrected low energy 
137
Cs photon
count rate in the 
241
Am window (IcLt) has to be subtracted from
the total dead time corrected 
241
Am count rate (Ia,totat):
Ia, true " 
Ia,total 
Ic, low
(35)
29
With calibration program 3, compton.txt (APPENDIX B), the
relation between 
the dead time 
corrected 
high energy 
137
Cs
photon count 
rate and the 
dead time corrected 
low energy 
137
Cs
photon count rate in the 241AM window is determined. Two brass
plates were placed in front the detector to block the emitted
241
Am photons and all the photons detected in the 
241
AM window
were then assumed to be Compton scattered 137Cs photons. A
teflon prism placed between the sources and detector provided
variable radiation attenuation. The Compton relation obtained
is shown in figure 5. A linear polynomial provides the best
fit for the data points. Hopmans and Dane (13) reported a
cubic Compton relation for approximately the same high energy
137
Cs photon count rate, using uncorrected count rates.
Apparently, the curvilinear character of the Compton relation
disappears when a dead time correction precedes the Compton
correction in the calibration procedure.
Determination of Attenuation Coefficients
Volumetric Attenuation Coefficients of Liquids
The volumetric attenuation coefficient of a liquid, Uij
where i refers to ...the liquid and ...j to the _gamma source,
consists of the product of the mass attenuation coefficient,
iand the density of the liquid, p
1
. To obtain U a
Plexiglas cell with four parallel compartments is placed
30
1000 -
0
11/15/89
750
4)
xY .1930 X+ 4.48849
. I r - .99988
Q4 500
O 250
0
0 0
0 
1000 
2000 
3000 
4000
High energy 
137
Cs intensity 
(cps)
FIG. 5. Compton scattering in the 241Amwindow versus
t~he high energy 
137
Cs intensity. All count rates are
dead time corrected.
between the sources and detector'such that the gamma beam is
perpendicular to the compartmental walls and the radiation
traverses through all the compartments. The path lengths
through the compartments are determined with a caliper.
First, count readings are obtained at five different locations
on the empty cell, and subsequently on the cell with 1-4
compartments filled with the liquid of interest. This
procedure results in a total of 25 count rates for each of the
two gamma sources. Since the liquid is the only material
31
causing attenuation, the linearized attenuation equation for
the compartmental cell, with n filled compartments, can be
written as
1n (Inj) - ln(Ioj) -Uijy
n  
(36)
where IOj is the count rate (cps) through the empty cell, Inj
is the count rate (cps) through the cell with n filled
compartments, and yn is the path length through n filled
compartments, i.e.
n
Yn" 1Yk n- 0,1,2,3,4 
(37)
k-0
where k refers to each individual compartment. A linear
regression of the 25 Yn-ln(Inj) data pairs yields the slope of
eq. (36). The absolute value of this slope equals the Uij of
the liquid. Calibration program 4, atten.txt (APPENDIX B),
guides the user through the process of data acquisition and
the linear regression. The program produces a plot with the
241Am and 
137
Cs regression lines. Examples of such plots are
presented in figure 6.
For the simultaneous determination of salt concentration
and volumetric water content, .. the linear ...relations between
liquid salt concentration, C, and its volumetric liquid
attenuation coefficient, UiLl
, 
have to be determined. With
calibration program 4, the ULj-values of a series of salt
32
10.0
Am (solid line):
y-m-.i9496 x + 9.55309
deionized 
waterr 
9
Ca (dashed line):
9.0 yw-.08550 x + 8.36401
9 . 0 -
r a w .9 9 9 8 9 2
.1-
0
7.0
0.0 2.0 
4.0 600
Path length through 
filled compartments 
(cm)
10.0j
0.2 n Nal
Am (solid line):
ya-.39059 x + 9,55557
r -99959
Cs (dashed line):
ya-.08704 x + 8.36327
r -S.99829
W..
7.4
0.0 2.0
Path length 
through
*I I I
4.0 6.0
filled compartments 
(cm)
FIG. 6. Relation between path length 
through a solution and
the natural log of the corrected 
count rates for deionized
water and a 0.2 n NaI solution. The 
absolute value of the
slope is the volumetric liquid attenuation 
coefficient.
?9.0.
4-)
0
U
'~8.0.
33
solutions over the concentration range of interest can be
obtained. Linear regressions of ULc and ULa versus C result in
the coefficients a and b of eq. (22) and r and s of eq. (23),
respectively. Figure 7 shows C-ULj relations of the salt NaI
and their calculated regression coefficients. The high
correlation coefficients of both lines show that the
assumption of linearity is acceptable.
Since some NAPL's dissolve Plexiglas, the compartmental
cell cannot be used to determine their UW. An example of such
a NAPL is trichloroethylene (TCE), a solvent frequently used
in immiscible flow studies. To overcome this problem, an
alternative compartemental cell was created for these liquids
by glueing four glass French square bottles in series to a
wooden base. This setup allows the use of basically the same
method to determine volumetric attenuation coefficients
described previously for liquids with the Plexiglas
compartmental cell. While running program 4, the user is
prompted to choose either the Plexiglas or glass compartmental
cell. The path length through each of the four French square
bottles is determined with a method outlined in the next
section of this chapter. The 
241
Am and 
137
Cs regression lines
for TCE, determined with the glass compartmental cell, are
shown in figure 8.
34
AMERICIUM
*-4
U
*4
0
0
0OH
*4J
0
*4
4o
0
0
U
0
*4
4J
Z
.1-I
0
0
U
0
0
rH
.41
oq
.-
.1-I
.1u
40
0o
. 10iO .15
Nal concentration
.20 .25
(eq/1)
CESIUM
' 3
w
y- .08520 + .00780 x
r - .93473
.15 .20
concentration (eq/1)
FIG. 7. Experimentally determined relations between the Nal
concentration of a solution and the volumetric attenuation
coefficient for 
241
Am and 137Cs. Each data point is obtained
with calibration program 4.
.45
.35
.25
. 15
y- .19527 + .97533 x
r .99995
r = .99995
.05
N
.05
NalI
.10
35
Tri chloroethy lane
10.0
Am (solid line):
y-.54746 x + 9.07054
a
9.0 
r .999951
Ca (dashed line):
yo-.10783 x + 8.33286
- r - .999932
0
o 7.0-
6.0-
5.0
0.0 2.0 4.0 6.0 11.0
Path length through filled compartments (cm)
FIG. 8. Relation between path length through TCE and the
natural log of the corrected count rates. The absolute value
of the slope is the TCE volumetric attenuation coefficient.
Mass Attenuation Coefficients of Soils.
The Plexiglas compartmental cell is employed to
determine mass attenuation coefficients of soils, sj (L
2
M
1
').
When (air dry) soil is used to fill the compartments, it is
likely to have different bulk density values, Pk, and
different volumetric water content values, 
6
k, for the
compartments k, k=l,..,4. With the soil and water as the
36
attenuating materials, the linearized attenuation equation 
for
the compartmental cell with only compartment 1 filled can 
be
written as
ln(Ij) + Uj6ly
1 
-ln(Ioj) - uSjply, (38)
In generalized form, with n compartments filled, eq. (38)
becomes
n n
ln(I ) + Uvj kyk - ln(Ij) -IL jEPkyk (39)
k-0 k-0
The count rates I, and Ioj are determined with calibration
program 5a, sat5a.txt (APPENDIX B). First, count 
rate
readings are obtained at five different locations, with 
all
four compartments filled with soil. Then, starting 
with
compartment 4, the compartments are subsequently emptied 
in
moisture cans and count rates are obtained at exactly the same
locations. The bulk density and volumetric water 
content
values of each compartment are determined gravimetrically.
When, finally, Pkyk is plotted against the left hand side 
of
k-G
eq. (39), the absolute value of the slope of the regression
line is the mass attenuation coefficient of the soil. The
linear regression and the calculation of the mass attenuation
coefficients are performed by calibration program 5b,
37
sat5b.txt (APPENDIX B).
Path Length Determination
At each position of a column or two-dimensional flow
container where measurements will be taken, the path length
through the porous medium needs to be determined. To
accomplish this, count measurements have to be obtained at
each location when the column or container is empty, resulting
in Ioj, and when the column or container is filled with
deionized water, resulting in I. Solving the linearized
attenuation equation, with deionized water as the only
attenuating material, for path length x,, gives:
x = n(Ioj/Ij)/Uwj (40)
With calibration program 6, x.txt (APPENDIX B), the path
lengths x at each position of interest are determined. It
has been our experience that the xa and xc-values at the same
position vary by a maximum of 0.85%. Normally, the average
path length, (xa+xc)/2
, 
is accepted as the true path length
through the porous medium at a certain position. Program 6
creates a data file where the count rates IoJ and path lengths
x are stored inseparate .subfiles This information can be
read from the data file in future applications.
38
Determination of the Incident Count Rates
Besides the values of the attenuation coefficients and
the path lengths, the parameters Ioj are needed to allow
application of eq. (6) .and (7), eq. (13) and (14), or eq. (26)
and (27).
For the simultaneous determination of the dry bulk
density and volumetric liquid content, the incident count
rates Ioj represent count rates through an empty column or a
2-dimensional container. Due to variation in the wall
thickness, the Ioj values are generally not constant in
position and have to be determined at each location where
future measurements will be taken. It is not necessary,
however, to include a separate calibration step to obtain
these data since values of Ioj at each location were already
obtained during the determination of the path lengths with
program 6. The Ioa and IO data were stored in the first and
second subfile, respectively, of the data file created with
program 6.
For the simultaneous determination of the volumetric
water content and the volumetric NAPL content or the
volumetric liquid content and the salt concentration, the
incident count rates 1oj are the count rates through a column
or 2-dimensional container filled with an oven-dry, non
expanding soil, and can be expressed as
i0 oj= I 
0
'exp (-,spAx)
(41)
39
where Io' is the count rate through the empty column or
container. To obtain the Ioj values directly at each position
of interest, count measurements can be taken with oven dry
soil. Assuming the dry bulk density remains constant in time,
the determined Ij data can be used in the appropriate
equations. Normally, however, as a result of wetting and
drying, some consolidation and therefore changes in dry bulk
density may occur and the initially determined Ij values will
not represent the actual situation accurately. In our
calibration process, a method is used to determine I that
takes into account the possibility of soil consolidation.
First, the column or container containing the oven-dry soil,
is saturated, followed by a bulk density determination using
eq. (6). The column or container is then drained and again
saturated, followed by another bulk density determination.
This drainage-wetting cycle is repeated until the dry bulk
density remains constant. Once the dry bulk density is known,
the Ij values can be calculated using eq. (41), since all the
other variables, appearing on the right hand side of this
equation, were previously obtained.
Calibration program 7a, bulk.txt (APPENDIX B), guides the
user through this part of the calibration process. After
measurements are taken- at each location, a printout with the
bulk density and volumetric water content values is produced.
The user decides to quit the program when the values remain
constant, taking into account a small variation, due to the
40
random nature of the gamma radiation, of approximately 0.01
g/cm
3
. With calibration program 7b, bulkin.txt (APPENDIX B),
the bulk density values that will be used to calculate Ij at
each' position are entered in the data file created with
program 6, containing the I ' and path length values.
Once the I0j values are determined, the calibration
process has been completed.
ERROR ANALYSIS
A probable error analysis for quantities determined
utilizing the dual gamma method gives an indication of the
achievable precision and is a tool to determine the optimum
path length through a porous medium. According to Bevington
(I), the expected variance of any quantity y, which is
calculated based on n independent observations xi, can be
expressed with the error propagation formula:
To evaluate the variance of y, the variance of each set of
independent observations xP and the functional relationship
between y and x
1
must be known.
41
Error Analysis for the Simultaneous 
Determination of the Dry
Soil Bulk Density. -p, and the Volumetric 
Water Content. ,
Eq. (6) and (7) show that 
ps and Ow are both functions 
of
the variables I
a
' Ioa, c
, 
I , x, Asa Asc, Uwa, and 
UwC.
Replacing x
i 
in eq. (42) by these variables 
results in
2 / ay)2 
/ay2 
(43)
Y I a a 08 
(1 rz
a oaC
S&
I
1 2
I
2
0l 
+  
"oX 
+ 
U zI
a(oc81 wa
(SY2 2 2
WC 04sltSCo) so ao so
where y denotes ps or O
w
.
Expressions for the partial derivatives 
appearing in eq. (43)
are obtained using eq. (6) and (7) 
for ps and O
w, 
respectively,
and are listed in 
table 2.
" Gardner et al. (9) stated 
that the variance of a quantity
y, determined with the dual gamma 
method, is the sum of the
variance due to the randomness 
in the radioactive emission,
a2 
a.
ay (re) and the variance 
due to the measurement 
or calibration
error, aY(Ca). The expressions 
for a y(re) and ay(ca) are,
respectively,
2re a oa 
2
c 2 oc(44)
42
ay 2 ,2) y2 4 / r
O2+ I10 
u
y\axC&)at 
Uwa 
wa lAUWC)WC
( 2 / 2
apsa) ~sa 
a 
5
C) c
Table 2. Expressions for Partial
Derivatives ac lax, and ap ,/ax
1
.
xi ao~lax! 
apg/axi
'a ILsc/xkI
8  
Uwc/xk~a
Ioa igsc/xkloa 
Uwc/xkIoa
x e-G/x 
-PS/x
Uwa ILCe/k -Uwcew/k
Uwc iPsa~w/k UwaGw/k
ALsa !LSCOS /k 
u wps/k
!SC'sc sapS/k 
*UwaPs/k
(45)
Several investigators (9, 5, 11) showed that the variance due
to random emission is the dominant error source in dual gamma
measurements.
The 
241
Am~ and 
137
CS gamma emissions follow a Poisson
43
distribution for which the best estimate of the variance from
the true mean for a single count measurement, N, is the count
measurement itself (16). The best estimate of the standard
deviation is JN. A common error is the association of the
standard deviation with the square root of any quantity which
is not a directly measured number of counts. The standard
deviation of the count rate I is not equal to /I since I is
not measured directly but obtained by dividing N by the
counting time t. Using the error propagation formula, the
variance of a dead time corrected count rate I can be written
as
i "aN) N(46)
Replacing a
2
N by It and the partial derivative by l/t, eq. (46)
becomes
a- I/t (47)
The counting time, to, used to determine the initial count
rates Ioa and IO is not necessarily the same as tm, the
counting time used to determine I, and I
c
. Eq. (47) can be
used to evaluate the variances of the 
137
Cs count rates:
ep 
1 
Ic/tm
(48)
44
/ - oC/t
o  
(49)
Eq. (47), however, is not directly applicable 
to the 
241
Am
count rates, because these count rates are corrected for
Compton scattering (26). The corrected count rate Ia can be
written as
Ia - Iat - (pI
c 
+ q) (50)
where It is the uncorrected count rate and the term in
parentheses is the correction due to the Compton scattering in
the 241Am window. The constants p and q of the linear Compton
relation are determined in the calibration process with
program 3, compton.txt. Using the error propagation formula,
the variance of Ia can be expressed as
1 a 21 aJ (51)
a 
(3Iat 
at+ 
C(5
Expressions for the partial derivatives in eq. (51) are
obtained using eq. (50) i resulting in (ala/alat) = 1 and
(aI/a
c
) =-p. With (2 I at/tm and 2C = Ic/tm, eq. (50)
can be rewritten to yield
(52)
2 m (at + 
p21)/t.
45
Replacing Iat by Ia + (pI
c
+ -q) gives
(53)
a- (+ 
(p
2
p) I q)/t
3
For o r a similar expression can be derived.
At this point, expressions for the partial derivatives
and the variances of the count rates appearing in eq. (44) are
known, so that the a
2
y(re) can be written as a function of the
counting times tm, and to, the 
attenuation coefficients bsa'
Msc, Uwa' and Uwc, the path 
length x, and the count rates 
a,
Ioa Ic, and Ioc. When Ia is replaced with eq. (4) and I
c 
with
eq. (5), this results in the following 
relations forooe(re)
and aP ,(ze), respectively:
1 r sc exp (x saP s 
+ 
XUwa w)
m xk 
Ioa
(p
2 
+ p) Iocexp (-xp p -xUw w) + q +
Ioaexp (-2x, saP -2xUwa w) )
pi exp 
(xpscPs+xUWW) 
+
xk 
oc
1 s s2 Ioa + (P
2 
+ P)Ioc sa 1
t
o
xk 2 
xk Ioc
(54)
(720
46
4172 Uw2 exp (xp saP s + XUwaOw) (55)
p,(re) 
Io
tm ..o+
(p
2 
+ p) IoCexp(-xscp, 
-xUcw
) 
+ q
Ioaexp (-2xp, sa s -2xUwaw) )
(Uwa)exp (xp,
5 
CP + xUWCO,)
xk 
Ioc
Uw 2oa + (p
2 
+ P) Io
c 
+ q + 1
to I2a 
xk Ioc
These formulas should be used to determine the optimum path
length through a porous medium for a certain experiment. 
The
values of the attenuation coefficients and the coefficients of
the Compton relation (p, q) are determined in the calibration
process and are assumed to be constant. Once a decision is
made on the wall material and its thickness S, values for the
incident count rates I. and IC can be obtained by observing
the count rates through the material with a thickness 
of 2S.
With known values for #sa, so, U.a., Uwc, Ioa, and Ioc, the
variance of 0, and ps and functions of Ps, Ow, tm, to, and x
only. For a fixed set of quantities Ps, Ow, tm, and to, the
variances or standard deviations as a function of path length
can be calculated and plotted. Examples of such plots are
presented in figure 9 for four different sets of ps, O, and
tm values. The counting time to, used during calibration,
should be maximized. Most often, the limiting factor is the
overflow of the scaler that is used to determine the integral
241
Am count rates. We used to=200 sec. for all applications.
47
0 -
0
o230-
~1o.
_____ 
w
- -Ps
I
/
/
/
/
/
0-
0 48 
12 16 
2(
path length (cm)
FIG. 9A. Relative standard 
deviation in p,, and 0
length x for p,=1.6 g/cm
3
, Ow-O.3 and tm=6O sec.m
50-
0-
0
~20-
~10-.diI
-I
o 4* 
12 16
p ath length (cm)
FIG. 9B. Relative standard deviation in p. and 8 versus 
path
length x for ps=lo6 g/cm
3
, GW-O~l and t*=60 sec.
?versus path
-PS
/
/
20
II I 
,
m
48
_____ 
w
-- PS
0 --V C II
0 4 
5 12 
16 20
path length (cm)
FIG.* 9C.. Relative standard deviation 
in ps and Ow versus path
length x for p,=1,2 g/c
3
, GO=O.*3 and tm=
6
O sec.
4.
path
12
length (cm)
16
-- Ps
20
FIG. 9D. Relative standard deviation in ps and 
O versus path
length x for pr=l.6 g/cm
3
, Gw=Oml and t=300 sec.
0
r-40-
o2 30
'~20-
~10-
02
II
I'
50-
*40
o0-
0
~10-
02
0
I1
0
49
The attenuation coefficients used have the following values:
Msa=
0
.
2 5 2 6 
cm
2
/g, sc=0.0817 cm
2
/g (Flint Shot 2.8 
sand),
Ua=0.1950 cm
"1
, U =0.0854 cm
"I 
(deionized water). For the
incident count rates, values of I =17000 cps and Im=810 cps
were used (based on measurements), which correspond with 0.25"
Plexiglas walls. The errors in ps are less sensitive than the
errors in 8,, while the accuracy in both Ps and GO improves
when their magnitudes increase. The optimum range for dual
gamma measurements appears to be 6-10 cm for O,, and 4-12 cm
for ps. The probable error can be minimized by using large
counting times. However; for some experiments this may
conflict with rapidly changing water content values.
The variance in p, and 8, due to calibration errors, as
described with eq. (45), depends strongly on the methods used
to determine the attenuation coefficients and the path length.
Our attenuation coefficients are absolute values of slopes
determined with linear regression. According to Steel and
Torrie (25), for n data pairs (XiYi), an unbiased estimate of
the true variance of the slope is given by
Sn (X
i  
-X)(Yi
2  
(56)Y)
2%i-1
50
The order of magnitude of the variances of the attenuation
coefficients, calculated using 
eq. (56), varies between 
10
and l0 9, indicating an insignificant contribution to the total
variance.
Based on the calibration method to determine the path
length x, the variance in the path length x is estimated to be
0.02 cm,. The contribution of the estimated variance of the
path length x to the total variance in p. and Ow is calculated
by evaluating the first term of eq,. (45), which can be
rewritten as 0.02y
2
/x
2
, where x is the path length in cm. The
relative standard deviation due to the variance in path length
x is then 14.1/x%, for both p. and Owj, independent of their
values and the counting time. For p. determinations, this
additional error has the same order of magnitude as the error
due to random emissions. 
For the more important 
OW
determinations, however, this error is in. most cases
considerably smaller than the error due to random emissions
and it does not affect the optimum path length range observed.
Error Analysis for the Simultaneous Determination
of the Volumetric Water Content.,,..and
51
to express the variance of Ow and On,
2 2 2
,,2 & a2 ay 0,2 
+
I 
X a + 
O'y 
oa + 
X c
aia) aica ) 
( (3ic
2 My 2 2
Cly. 0.2 0 or2 cftj 2
ai 
X oc 
+ 
x + 
au 
a Uwa 
+
oc ) 
wa )
O'y 
O*Y or2 
&J 2
wc + ^T a U na + aaui 
Or Unc
OUWC) 
COVJna 
nc
(57)
where y denotes 0. or On. The 
partial derivatives in eq. (57)
are evaluated using eq. (13) and (14) for OW and On
respectively, and are listed in table 
3.
Table 3. Expressions for Partial
Derivatives ae~/ ax
1
and ae c/x
1
.
xi aojlax
1
o/ x
Ia U /XKIa Unc/XK~a
oa uwcxia 
n/K
Ic -UVXKI Un/xKIC
Uwa UbwcGw/K Uncew/K
Una -- Uce/K
Uric Uwaor/K 1-Unar/K
K = UcUna -Uw
4
aUfl
52
The total variance a
2
Y is the sum of the variance due to
random emission, a
2
y(re), consisting of the first four terms of
eq. (57), and the variance due to the 
calibration error,
a2y(a), consisting of the last five terms 
of this equation. By
replacing the partial derivatives by the appropriate
expressions in table 3 and the count rates by their
expressions obtained in the previous section, a
2
y(e can be
written as a function of the counting 
times, the path length,
the attenuation coefficients of the two liquids, and the count
rates. When Ia is substituted 
with eq. (11) and I
c 
with eq.
(12), this relation becomes, 
for the variance of , 
and One
respectively,
1 
Un( expW(xUwa, + xUn) +(58)
,(ze) t \ / _e
(p
2 
+ p) Ioexp(-xUwc 
, -xUn
0
O,) + q
Iaexp (-2xUwa6 
-2xUna )
SUna
2 
exp (xUwc 0 + xUncn) 
1+
xK 
Ioc
1 Unc 2 Ioa + (P 
+ 
P) Ioc 
+ 
q Una2
to- xK 
xK Ioc
a
9
~ I1 Uwc 2 exp (xUwa0v 
+ xUna n) 
(59)
tnxK 
IO
0% 1, - A
53
With these formulas, the optimum range of the thickness of the
porous medium for a certain experiment can be determined. The
values of the volumetric attenuation coefficients of water and
the NAPL are determined in the calibration process. Recalling
that the initial count rates I. are calculated from I, =
I' exp(- ssX), assumed or known values for Io' and ps reduce
the variables that determine a
2
y(re) to the counting times tin
to, the path length x, U, and on. For a fixed set of values
for tm, e, and on, a
2
y(re) can be calculated and plotted as a
function of the path length. Four examples are presented in
figure 11 for a column with 0.25" Plexiglas walls. Based on
measurements, it was assumed that I '=17000 cps, Io'=8100
cps, and to=200 sec. 
Trichloroethylene, with 
Un=0.54746 cm'I
and Un=0.10783 cm
"I
, is used as the NAPL, while the dry bulk
density of the soil is assumed to be 1.6 g/cm
3
. It appears
that the optimum ranges for both liquids are approximately
between 4 and 8 cm.
Since the variances of the volumetric attenuation
coefficients are very small and range from 108 to 10
9
, the
path length term is the only significant term of a
2  
.For
both 0, and n, the relative standard deviation, assuming a
path length variance of 0.02 cm, is 14.14/x%. This additional
error is too small to change the optimum ranges observed.
50-
0
o0-
230-
~10-
54
I.
I
I i
I
I
/
u - tI1 1 1 1 1 1
0 4< 5 12 16
path length (cm)
FIG. 1bA. Relative standard 
deviation in Ow
length x for Gw=O .15, 
6
On=O*15 and t,=60 sec.
8
20
and On~ versus path
~50Q
0
0
rd20-
~10.
02
I
I
I
w
--
n
~1I
.1
S/
/
/
/
0 i IIII
0 48 12 
16 20
path length (cm)
FIG. lOB. Relative standard deviation in Ow and O versus path
length x for Gw=O .25, Gn=O.O5 and t,=60 sec.
Aa50-
*0 -
0
o3.,0-
2 30
cz10-
L.j2
A
55
I
I
/
I
I
/
/
I
I '
-i
--
/
0 4' 12 16
path length (cm)
FIG. 1bC. Relative standard 
deviation in 8w
length x for Uw=0.05, On=0.25 and tm=
6
O sec.a
50-
*0-
o0-
20-
10- 10
20-
and On, versus path
Ow
O n
0 4~ 8 12 16 20
path length (cm)
FIG. 1OD. Relative standard deviation in Ow and O versus path
length x for 0 =0.l5, 8 =0.15 and tm=300 sec.
I
56
Error Analysis for the Simultaneous Determination of the
Volumetric Liquid Content, 8., and the Salt Concentration, C.
According to eq. (26) and (27), both 
OW and C are
functions of Ia, I , I
c, 
Io, a, b, r, s, and x. Using the
error propagation formula, the variance of 
W 
and C can be.
expressed as
( 2 2 (_
d oC 81x oa
( 2 +++
where y denotes O8 or C.
The partial derivatives in eq. (60) are evaluated using
eq. (26) and eq. (27) for C and 8
,
I respectively, and are
listed in table 4.
The first four terms of eq. (60) make up the variance 
due
to random emission, aZ while the remainder is the variance
due to the calibration error. Replacing the partial
derivatives by their expressions in table 4, and the variances
of the count rates by their expressions derived previously,
a
2
yre) can be written as a function of ..the counting time, the
path length, the count rates, and the coefficients of the
linear relations between salt concentration and the volumetric
attenuation coefficient. When I
a 
is replaced by eq. (24) and
I
c
by eq. (25), this relation becomes, for e
W
and C,
57
Table 4. Expressions for Partial Derivatives
clw/xiand a'c/axj; n =xGw(s~a + bC) - b(r + sC]).
Ia -s/(x(as- br)Ia (r + sC) /Ian
'oa s/{(x (as br) Ioa -(r + sC) /1
08
n
Icb/(x(as br)IC -(a + bC) /I~n
i~c -b/ (x (as br) IOC (a + bC) /Iocn
X Ow/X 0
a XGwULc 
ULc (SULa bULC)
b xeWuLa 
CLc (SULa bULC)
r O (r- 
ULC) (as- br) -ULa 
(SULa- bU
1
~)
s Ow (a + UL.) 
(as br) -Cka (SUia 
bULC),
respectively,
a8e(ze)a
1 s ~ ~2(exP (xO.w(a+ C
t[x(s-br) 
+ +sC) O))
(p 
2 
+ P) I Cexp 
(-xO,,(zr 
sC + q+
ga exp (-2xO,,(a +bC))
( b 2 eXp (XO-,(-3 + -SC))]
x (as -br-)) 
bc ---
br) 2( o + (P
2 
+ P)I
0
o + q'
(x (as b br) 
) ( ic)]
(61)
) +
58
O7
2 
-1 r + sC2( exp(xO,(a + bC) (62)
C 
t(re) t e [( 
n  
Ioa
(p
2 
+ p) Iocexp(-xOw(r + sC)) + q +
I exp(-2xO,(a + bC) 
) )
a +bC 2 exp(xO,(r + sC)) +
1 r + sC2( Ioa 2 + (P
2
+P)Ioc + q +a +bC 21
to n 
2 n Ioc
These formulas. are the basic tools to calculate the
optimum thickness of the porous medium for a certain
experiment. Following the same reasoning as in the previous
section, it can be shown that a
2
y(re) is a function of x, tm, to,
8
L 
and C only, when the count rates Ioj' and the dry bulk
density are known. For each fixed set of t*, Uw, and C values,
the relation between the variance (or standard deviation) and
the path length x can be calculated as will be demonstrated by
four examples, with Nal as the attenuating salt. It was
assumed that to=200 sec, IO'=17000 cps, I '=8000 cps and
ps=l.6 g/cm
3
. The coefficients a, b, r, and s for 
the Nal
concentration-volumetric attenuation coefficient relations
have values of 0.19527 cm
"I
, 0.97533 cm
2
/eq, 0.08520 cm'
I
, and
0.00780 cm
2
/eq, respectively (see figure .7) ..
Figure 11 shows that the optimum path length range for C
is 5-9 cm and that the error in 0, is rather insensitive to
variations in x over a range from 4 to 20 cm.
When Nal is used, it can be shown that the contribution
59
When Nal is used, it can be shown that the contribution
to the total variance of the a, b, r and s terms in eq. (60)
is negligible. Since b and s are the slopes of the salt
concentration-volumetric attenuation coefficient relations of
241
Am and 
137
Cs, respectively, their 
variances can be calculated
with eq. (49), resulting in a
2
b=-3.0*10-6 and a
2
s=4.8*10
"7
, for
the lines shown in figure 7. The variances of the intercepts
with the attenuation axes of the plots, a and r, are
calculated using Steel and Torrie (25):
S62
o2(- _______ 1 ' - KX?X)(i -i '
n
2x n2 n 
'
(
X)( - - ) (63)
i-i +- + 
(yi 2 2i- 
in
(n-2)
This results in a
2
a=-6.8*10
.7 
and U2r=1. 1*10
7
. All the
coefficient variances are very small compared to the variances
in the count rates and can, therefore, be ignored. Another
probable error analysis for the simultaneous determination of
salt concentrations and volumetric liquid contents was done by
Grismer et al. (11). They based their analysis on the
assumption that azY=
2
1N, where N is the observed count (Grismer
et al. used the symbol I for both the count rates and observed
counts, which is especially confusing in the probable error
analysis and should be avoided. For clarity, we use N for
observed counts and I for count rates). Generally, the
100 -60
* I
o 
1
0 .
80- \C
60 
/
CO
S4 _ 12 1
20
'Ct
0.1
o 0- /
40 -
0 III11I I1
0 4 8 12 
16 20
path length (cm)
FIG. 11A. Relative standard deviation in Ow and C versus path
length x for G =0.35, C=0.1 meg/mi, and tm=60 sec.
10o.0
II ow
~~'
8 0
i-1 
- C
$ 01
.261
0
14
.4.
%%wol
20
61
1
I
I
I
w
-- C
100
80,
60-
40,
I
/
04 8 12 16 20
path length (cm)
FIG. liC. Relative standard deviation in Ow and C versus path
l ength x for Ow=Ool, C=Ool mecT/ml, and tm=60 sec.
100-
00-
00.
40-
S20-
0
I
I I
I
I.
I ..
/
/
/
8 12
path length (cm)
FIG. liD. Relative standard deviation in Ow and C versus path
length x for OW=O.35, C=O.Oi meg/mi and tm=300 sec.
/
62
assumption that a
2
,=/N is used since radio-active 
decay
processes can be approximated by a Poisson 
distribution (9, 5).
Grismer et al. (11) referenced Knoll (16) to back up their
assumption. Knoll, however, clearly states 
that for
radioactive decay processes, the variance of an 
observed count
is the count itself. The errors due to random 
emission in the
paper by Grismer et al. (the first terms of their eqs. (11)
and (12), p. 166), should therefore be divided by a factor 4.
Furthermore, the form in which their eqs. (11) and (12) 
are
written allows easy calculation of the probable error of
actual measurements only, whereas our eqs. (61) and (62) can
also be used to predict probable errors for all theoretical
combinations of the dependent variables.
LITERATURE CITED
(1) Bevington, P.R. 1969. Data Reduction and Error Analysis
for the Physical Sciences. McGraw-Hill, New York.
(2) Corey, J.C., S.F. Peterson, and M.A. Wakat. 1971.
Measurement of Attenuation of 
137
Cs and 
24 1
Am Gamma Ray for
Soil Density and Water Content Determinations. Soil Sci.
Soc. Am. Proc. 35:215-219.
(3) Davidson, J.M., J.W. Biggar, and D.R. Nielsen. 1963.
Gamma-radiation Attenuation for Measuring Bulk Density
and Transient Water Flow in Porous Materials. J. Geophys.
Res. 68:4777-4783
(4) Ferguson, A.H., and W.H. Gardner. 1962. Water Content
Measurement in Soil Columns by Gamma Ray Absorption. Soil
Sci. Am. Proc. 26:11-14.
(5) Ferrand, L.A., P.C.D. Milly, and G.F. Pinder. 1986. Dual-
gamma attenuation for the determination of porous medium
saturation with respect to three fluids. Water Resour.
Res. 22:1657-1663.
(6) Fritton, D.D. 1969. Resolving Time, Mass Absorption
Coefficient and Water Content with Gamma-ray Attenuation.
Soil Sci. Soc. Am. Proc. 33:651-655.
(7) Gardner, W.H. 1986. Water 
content. In A.Klute (Ed.)
Methods of Soil Analysis, vol. 1, Physical and
Mineralogical Methods. American Society of Agronomy, Soil
Science Society of America, Madison, Wis.
(8) Gardner, W.H. and C. Calissendorff. 1967. Gamma-ray and
Neutron Attenuation in Measurement of Soil Bulk Density
and Water Content. In Proceedings of a Symposium on the
Use of Isotope and Radiation Techniques in Soil Physics
and Irrigation Studies, International Atomic Energy
Agency, Vienna, pp. 101-113.
(9) Gardner, W.H., G.S. Campbell, and C. Calissendorff. 1972.
Systematic and Random Errors in Dual Gamma Energy
Soil Bulk Density and Water Content Measurements. Soil
Sci. Soc. Am. Proc. 36:393-398.
(10) Goit, J.B., P.H. Groenevelt, B.D. Kay, and J.G.P. Loch.
1978. The Applicability of Dual Gamma Scanning to
Freezing Soils and the Problem of Stratification. Soil
Sci. Soc. Am. J. 42:858-863
63
64
(11) Grismer, M.E., D.B. McWhorter, and A. Klute. 1986.
Monitoring Water and Salt Movement in Soils at Low
Solution Contents. Soil Sci. 141:163-171
(12) Gurr, C.G. 1962. Use of Gamma Rays in Measuring Water
Content and Permeability in Unsaturated Columns 
of Soil.
Soil Sci. 94:224-229.
(13) Hopmans, J.W. and J.H. Dane. 1985. Calibration and Use
of a Dual-energy Gamma System, Pressure Transducers, and
Thermocouples to Determine Volumetric Water Content, Dry
Bulk Density, Soil Water Pressure Head, and Temperature
in Soil Columns. Agron. and Soils Dept. Ser. 101, Ala.
Agr. Exp. Sta., Auburn Univ.
(14) Hopmans, J.W. and J.H. Dane. 1986. Calibration of a Dual
Gamma Radiation System for Multiple Point Measurements in
a Soil. Water Resour. Res. 23:625-631.
(15) King, L.G. 1967. Gamma-ray Attenuation for Soil-
water-content Measurements Using 
241
Am. In Proceedings of
a Symposium on the Use of Isotope and Radiation
Techniques in Soil Physics and Irrigation Studies,
International Atomic Energy Agency, Vienna, pp. 101-113 
.
(16) Knoll, G.F. 1979.. Radiation Detection and Measurement.
Wiley, New York.
(17) Lederer, C.M., J.M. Hollander and I. Perlman. 1967. Table
of Isotopes, 
6
th edition. John Wiley and Sons, inc., New
York.
(18) Lenhard, R.J., J.H. Dane, J.C. Parker, and J.J.
Kaluarachchi. 1988. Measurement and Simulation of One-
dimensional Transient Three-phase Flow for Monotonic
Liquid Drainage. Water Resour. Res. 24:853-863.
(19) Loch, J.P.G. and B.D. Kay. 1978. Water Redistribution in
Partially Frozen, Saturated Silt under Several
Temperature Gradients and Overburden Loads. Soil Sci.
Soc. Am. J. 42:400-406. ..
(20) Loch, J.P.G., B.D Kay, and P.H. Groenevelt. 1980.
Comments on "Moisture Content and Bulk Density
Measurements Using a Dual-energy Beam of Gamma
Radiation." Soil Sci. Soc. Am. J. 44:663.
(21) Mansell, R.S., L.C. Hammond, and R.M. McCurdy. 1973.
Coincidence and Interference Corrections for Dual-energy
Gamma Ray Measurements of Soil Density and Water content.
Soil Sci. Soc. Am. Proc. 37:500-504.
65
(22) Nofziger, D.L. and D. Swartzendruber. 1974. Material
Content of Binary Physical Mixtures as Measured with a
Dual-energy Beam of Gamma Rays. J. Appl. Phys. 45:5443-
5449.
(23) Reginato, R.J., and C.H.M. van Bavel. 1964. Soil-water
Measurement with Gamma Attenuation. Soil Sci. Soc. Am.
Proc. 28:721-724.
(24) Soane, B.D. 1967. Dual Energy Gamma-ray Transmission for
Coincident Measurement of Water Content and Dry Bulk
Density of Soil. Nature 214:1273-1274.
(25) Steel, R.G.D. and J.H. Torrie. 1980. Principles and
Procedures of Statistics. A biometrical Approach. McGraw-
Hill, New York.
(26) Stillwater R. and A. Klute. 1988. Improved Methodology
for a Collinear Dual-energy Gamma Radiation System. Water
Resour. Res. 24:1411-1422.
(27) Stroosnijder, L., and J.G. de Swart. 1974. Column
Scanning with Simultaneous Use of 
241
Am and 
137
Cs gamma
Radiation. Soil Sci. 118:61-69.
(28) Swartzendruber, D., R.L. Warren, and D.L. Nofziger. 1980.
Reply to "Moisture Content and Bulk Density Measurements
Using a Dual-energy Beam of Gamma Radiation." Soil Sci.
Soc. Am J. 44:663
(29) Van Bavel, C.H.M., N. Underwood, and S.R. Ragar. 1957.
Transmission of Gamma Radiation by Soils and Soil
Densitimetry. Soil Sci. Soc. Am. J., 21:588-591.
(30) Wood, B.S.C. and N. Collis-George. 1980. Moisture Content
and Bulk Density Measurements Using a Dual-energy
Beam of gamma radiation. Soil Sci. Soc. Am. J. 44:662.
APPENDIX A
ASYST Software
1. Introduction
The software package ASYST 
3.0 is an important tool 
in
the automation of the dual-energy 
gamma radiation system.
ASYST is especially designed 
for data acquisition, analysis
and graphics. The package contains 
a pre-written acquisition
driver for the Metrabyte MBC-488 
interface board. This driver
makes it possible to communicate directly 
to General Purpose
Interface Bus 
(GPIB) devices. 
The GPIB devices 
in our dual-
energy gamma radiation system are:
- The Aston 800 interface (controls the stepper motors)
- The Tennelec TC 589 buffered interface 
(controls the
timer/multiscaler)
- The Hewlett Packard HP 3497A data acquisition 
control
unit (interfaces pressure transducers 
and
thermocouples).
ASYST can perform analysis. functions 
varying from curve
fitting to matrix algebra. The results 
can be plotted with
the graphics part of the package. This 
part of the software
is not as developed as the GPIB interfacing 
and the analysis
parts.
2. The proqramminj lancjuacje Asyst .. .........
ASYST is a programming language, allowing the entry 
of
simple commands or complete programs. The language 
and its
66
67
commands, called words, are extensively described 
in the ASYST
manuals. In this section, a brief overview of the language is
presented.
Numbers are manipulated in a special 
area, called the
number stack, which holds numbers for later 
operation. The
word stack.display shows the contents of 
the stack on the
screen. The last entered element is at the top 
of the stack.
After entering 25 and 1.6, the word stack.display gives 
as a
result 1.6 25. To erase the contents of the stack, 
stack.clear
has to be entered. There are many words available to
manipulate the number stack contents, e.g. the word 
swap will
cause the two numbers on top of the stack to switch 
places.
ASYST uses the Postfix Notation. Two numbers are
required on the stack before a function can be performed,
e.g., an addition. As in other high level languages, 
integer,
real, and complex numbers are supported, all with single 
or
double precision. In addition to numbers, ASYST supports
named variables (scalars) and arrays. Like numbers, a scalar
and an array element can be either a single or double
precision integer, real or complex number. For example, to
create an integer with the name alpha, one has to enter integer
scalar alpha. To assign the value 5 to this variable, 5 alpha :=
must be entered. The same words that manipulate numbers and
scalars can be used to manipulate arrays and their elements.
Arithmetic operations can be applied to arrays regardless of
the size of the array or the dimension. The maximum array
68
size is 64 K.
A very important feature of ASYST is its ability to
create personal commands. These commands are called colon
definitions. An example of such a colon definition is
? a
10 ramp sqrt ;
This colon definition is named a and it leaves the square root
of a 10-element ramped (each element contains its index
number) array on top of the number stack. Every colon
definition starts with a colon and ends with a semicolon. The
words inside a colon definition are not executed when created.
Only when the name is entered (in this case a) are the
individual words that make up the colon definition executed.
All personal words inside a colon definition have to be
defined previously. Colon definitions form the basis of ASYST
programming and groups of newly created words can be included
in subsequent new definitions.
ASYST also supports strings which are stored temporarily
on the symbol stack. An example of a string in ASYST is " this
is a string" . All strings have a space after the first ". The
space is called a delimiter between the " and the actual
string while the second " denotes the end of the string.
String variables can be defined and text can be assigned to
them. Strings are used extensively in the communication with
GPIB devices.
Execution of some ASYST words (e. g., comparisons) results
69
in a TRUE/FALSE value. These values are also placed on 
the
symbol stack. Besides the common mathematical operations (+,-
,*,:,**,exp, etc.), 
ASYST allows the 
use of trigonometric,
hyperbolic, and special functions. For data analysis
purposes, ASYST has words to carry out statistical functions.
For our applications, the words related to the analysis of
variance (ANOVA) are very useful. To create programs, ASYST
has a simple word processor. An array editor can be used to
view the array, change values, and perform calculations with
elements.
As for any other higher level computer language, the use
of conditional statements such as if...then and if...else...then
is relatively straightforward. The conditional statements are
frequently used in the calibration programs listed in Appendix
B. Although the statements resemble the Fortran and Pascal
statements, they work completely different. For example,
if... else...then checks the top of the symbol stack for a TRUE
or a FALSE condition. If this condition is TRUE, the series
of words between if and else are executed, if the condition is
FALSE, the words between else and then are executed. The
TRUE/FALSE condition on top of the symbol stack is removed by
the if word. Quite often a conditional statement is executed
after a comparison. Comparisons leave a TRUE or FALSE on top
of the symbol stack. The simplest comparison operator is the
word -, which compares the two top numbers on the number stack
and determines whether or not the two numbers are equal. If
70
they are equal, a logical TRUE is placed on the symbol stack,
if not, a logical FALSE.
To repeat a series of statements several times, loops can
be used. Frequently the statements that have to be repeated
are placed between the words do and loop. Other useful loops
are begin... until and begin.., while.. .until.
3. Asyst overlays and configuration
ASYST's possibilities are extensive. Most applications
use only a few of them simultaneously. To avoid the
occupation of computer memory by unused parts of ASYST, it is
possible to select the capabilities you want to have
permanently resident. ASYST can be configured into a custom
version that meets the requirements of a specific application
and the hardware being used. Such a system is smaller than a
complete system and is more memory efficient. Words (commands)
are grouped according to their function. These groups of
words are called overlays. An overlay is defined as a' file
containing pre-compiled word definitions. These word
definitions are. basically the same as colon definitions.
Examples of system overlays are:
gpib. sov contains all the GPIB words;
datafile.sov provides all the words -to read and write
-
numerical arrays from and to a disk.
The user can also create overlays containing his own
colon definitions. These overlays are called application
71
overlays. The ASYST configuration menu 
allows the user to
permanently load system and application overlays 
into a custom
version. This menu also helps to configure the hardware, 
the
memory, and the GPIB board. All configurations 
have a default
setup. These defaults can be changed by following 
the menu
commands. For example, to obtain a custom version 
that
contains all graphics -overlays available, one has to
permanently load the overlays charts and HP plotter driver 
from
the menu.
To work with the dual-energy gamma radiation system and
the HP plotter, a custom version is set up with the following
permanent overlays, chosen from the configuration menu:
GPIB master, type I nec GPIB driver, data files, HP plotter driver,
statistics and special functions, polynomial operations, matrix and curve
fit functions. With these overlays permanently loaded, we are
able to control the system, to acquire data, and to perform
data analysis and graphics.
4. Interfacina the computer with GPIB devices
The General Purpose Interface Bus, or GPIB, was developed
by Hewlett Packard to interface computers with laboratory
instrumentation. The GPIB was adopted in 1975 by the
Institute for Electrical and Electronic Engineers as the IEEE-..
488 standard. This standardization was necessary to improve
interfaces with highly specialized software and hardware.
Most of these interfaces were designed for only one specific
72
task and were very inflexible. The main advantage 
of the GPIB
is the reduction of the hardware requirements of a computer to
one interface for the entire system compared to one interface
for each instrument. Furthermore, GPIB interfacing allows for
board independent generalized software. ASYST is a software
package that makes use of this feature. During a configuration
procedure of a certain custom version in ASYST, the used board
has to be specified. From then on, the software package
internally takes care of the conversion of generalized
commands to Metrabyte specific commands.
A primary address has to be assigned to all of the GPIB
devices. This address can be set with DIP switches on the
back of, or inside, each device. For example, to assign the
name motor to the Aston 800 actuator with a DIP switch setting
of 2, the statement 2 gpib.device motor has to be entered.
Communication between the computer and the GPIB devices
on the bus is established by specifying one device which sends
information and at least one other device which receives this
information. The device receiving the commands is called
listener, while the device sending the information is called
talker. For most applications, the computer is the controller
of the system, which chooses the talker and listeners and
speCifies the commands to be sent. The-controller can also be
the talker or a listener. Frequently, the computer is the
talker. To set up communication with a GPIB device, one has
to enter the following commands:
73
send. interface. clear
remote. enab 7l e. on
With the first command the computer is made the controller and
it receives automatically the name me. The second command
allows all the GPIB devices connected to the bus to respond to
commands from the controller. To make the Aston 800 actuator
the listener and the computer the talker, the following
commands have to be entered:
untalk unlisten
motor listener me talker
It is important to distinguish between GPIB commands and
device commands. A GPIB command is a specific command from
the IEEE standard list, while the device commands are strings
of data defined by the manufacturer of the devices. The
commands differ from device to device. Commands like talk,
talker, unlisten are all GPIB commands. The string
"smO,ma+000000" is a device command which will cause the Aston
800 actuator to move sources and detector vertically to its
original position.
Sending a device command from the talker to the listener
is not enough for a command to be executed. End Of
Transmission delimiters are needed to tell the device
receiving the data that there are no more data to be sent.
Our three GPIB devices use End Of String (EOS) characters as
End Of Transmission delimiters. Every string sent from or to
the Aston 800 actuator and the HP 3497A has to end with a
74
carriage return-line feed as the EOS command. 
For example, in
order to let the sources and detector 
move to their original
position, the talk device command 
to be sent over is
smOma+000000" "cr. if" "cat. The "cat 
command catenates the two
substrings. The Tennelec 
TC589 uses ASCII character 13
(carriage return) as an EOS character. 
Before a device
command can be sent over to the Tennelec TC589, 
one has to
enter:
gamma eos.on
13 eos. character
One very important feature of 
the IEEE-bus is the
possibility to do service requests. A device 
can ask for
service from the controller. With our system 
this feature is
used to poll the Aston 800 actuator. It takes some 
time for
the sources and detector to move from one location 
to the
next. Since the next line in the program should 
only be
executed after the sources and detector have 
arrived at the
new location, the computer can poll the Aston 800 
actuator for
its status. The actuator returns what is called 
a "serial
poll byte," which determines the status of the device 
and the
condition which is responsible for this. It changes 
when the
stepper motor has completed the movement of the sources and
detector. The command sequence--of-a serial ..poll of the
actuator is:
begin
stack. clear
75
serial7.pol 1.enable
untalk unlisten
me listener motor talker
stack. 7listen
serial7.pol.disab7le
untalk unlisten
16 =
until
Once the sources and detector start to move, a serial 
poll
loop is done until the number 16 is returned. This number
represents the arrival of the sources and detector at the new
location. The next line in the program can then be executed.
APPENDIX B
Listing of Calibration Programs
This Appendix contains 
the ASYST programs for 
the
calibration of the gamma radiation 
system. The format is
similar for all the programs and 
is self explanatory.
SPECTRUM. TXT
REL. TXT
COMPTON.TXT
ATTEN.TXT
SAT5A. TXT
SAT5B. TXT
X.TXT
BULK. TXT
: Determination of the energy 
spectra.
: Determination of the-relation 
between the
integral 
241
Am count rate via the TC 450 and
the NC-22.
: Determination of the Compton relation.
: Determination of volumetric 
attenuation
coefficients of liquids.
: Count rate data acquisition for 
mass
attenuation coefficient determination 
of
soils.
: Calculationof mass attenuation 
coefficients.
: Determination of path lengths.
: Determination of dry bulk density 
values.
When the dry bulk density values 
are known,
the incident count rates can be 
calculated.
76
\ +++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ PROGRAM SPECTRUM.TXT
\ Calibration program 1
\ This program determines 
the spectrum of Americium-241 
or
\ Cesium-137.
\ Data pairs (E, count rate) 
are stored in a specified
\ datafile.
\ A plot is made of E vs. count rate.
\ Minimum ASYST overlay requirements 
on custom version:
\ gpib master,
\ type 1 nec gpib drive,
\ data files
\ HP plotter driver
\ charts
\ At the c:\ prompt enter c:\asyst\total 
in order to
\ get into a correct ASYST environment. The 
overlays
\ mentioned above are all part of the custom 
version
\ 'total'.
\ When in the correct ASYST environment (or 
custom version)
\ do the following:
\ 1. Load the program by entering:
\ load c:\asyst\spectrum.txt
\ 2. Enter at the subsequent OK prompt:
\ spectrum.determination
\ 3. Follow the instructions carefully.
07/19/89
M.Oostrom
\ ++++++++++++++++++++++++++++++++++++++++++* +
\ All scalars, arrays, and strings have to be 
declared in the
\ beginning of each program.
\ ************************************************** 
*******
real scalar counting.time
\ Counting time set on the timer/multiscaler
real scalar count.time..
\ Counting time transferred to buffer
real scalar rate
\ Calculated count rate
real scalar count
\ The total number of measurements done
77
real scalar deltaE
\ The width of the window (volts)
real dim[ 2 ] array dumx
\ A dummy array used in the graphics part
real dim( 2 ] array dumy
\ A dummy array used in the graphics part
real dim[ 51 , 2 ] array res
\ The data are temporarely stored in this 
array
\ the first column contains the E values, the
\ second the count rates.
\ Position (51,1) is reserved for 
the number of
\ measurements done. Position (51,2) 
is reserved
\ for the delta E used.
real dim[ 50 ] array x
\ For programming reasons the first column 
of array
\ res (except for position 51) is copied into 
array x
real dim( 50 ] array y
\ The second column of array res is copied into 
y
0. res := 0. x := 0. y :=
\ The arrays are initialized
real scalar dummy
\ While converting 
the string (transferred 
from
\ the buffer to the computer) into numbers, 
some
\ of the acquired numbers are useless.
\ These useless numbers are assigned 
to the variable
\ dummy. This saves memory space.
72 string sl
\ The string that is transferred to and 
from the
\ buffered interface. This string 
contains the infor-
\ mation obtained from the three scalers 
and the
\ timer of the timer/multiscaler TC 535P.
72 string namefile
\ The name of the output file.
2 gpib.device gamma
\ The buffered interface is named 'gamma'.
\ The number 2 refers to the way the dip switches
\ are set inside the bufferred interface.
\ The colon definition 's' is used after the 
prompt 'Press
\ any key when ready'. After pressing a key, a 
number value,
\ representing the key, and a T/F value are left 
on the
\ number and symbol stack, respectively. The words 
drop and
\ ?drop delete the expressions on both stacks. 
After this,
\ the screen is cleared and the next prompt will 
appear
\ in the middle of the screen. The 'p' macro also
\ clears the screen and moves the cursor to the 
middle of
\ the next screen.
78
: s
pokey drop ?drop screen.clear 10 1 do cr loop ;
: p
screen.clear 10 1 do or loop ;
\ This macro takes care of the GPIB connections. In this
\ program, the only GPIB device is the buffered interface.
\ All commands to the buffer need a carriage return as an
\ 'end of string'. The ASCII representation of a carriage
\ return is 13.
: interface
send. interface, clear
remote.enable.on
gamma eos.on
13 eos.character ;
\ This colon definition lets the user set up the gamma
\ system properly in order to do the spectrum determination.
\ After this, the computer is made the controller
\ (me talker) and the buffer the listener (gamma listener).
\ A 'control R' catenated with a 'control C'
\ is sent to the buffer in order to set the scalers to 0.
\ ******************************************************
: pre.comments
cr
." Turn power to NH85 off." cr
." The NH85 is to the right of the buffered interface." cr
cr ." Press any key when ready." s
. Set stop/count button of TC 535P on 'count' and" cr
" set .01 min/.1 sec button on '.1 sec'." cr cr
.Press any key when ready." s
" Preset timer." cr
" (Timer formula: timebase*n*10**m.}" cr
." Set CH A PRESET to '96'" cr cr
." Press any key when ready." s
." Disconnect all three analyzer-timer/multiscaler" cr
. connections." cr
." NA-22 from channel A, pos" cr
." TC 450 from channel B, pos" cr
." NC-22 from channel C, pos (on the back)" cr cr
." Press any key when ready." s
79
." The analyzer of interest should be in differential" cr
0" mode. Check this." cr
*" (For the TC 450 use 'window/asym'. The top knob then" 
cr
." represents delta E & 10. E can be manipulated" cr
." with the 'lower level' knob.)" cr cr
" Press any key when ready." s
." Connect the analyzer of interest to channel 
B, pos." cr
cr ." Press any key when ready." s
." Make sure the readout mode of the buffer TC 589 is" cr
." on auto recycle and that the output format is" cr
." on normal." cr cr
." Press any key when ready." s
." Turn power of NH85 on and wait till the counting" cr
." stops." cr cr
. Press any key when counting is finished." s
interface
untalk unlisten
gamma listener me talker
18 ASCII" 3 ASCII" "CAT TALK ;
\ *k*********************************************************
\ All data and some additional information must be stored in
\ a data file. The data file contains 7 comments and 1
\ subfile for the actual data. The exact organization of 
a
\ data file has to be specified after the command
\ 'file.template'. The 'file.create' word creates a file
\ according the latest 'file.template' definition.
\ The name of the data file is entered 
and is assigned
\ to a string named 'namefile'.
\ The user is prompted to enter the 7 comments:
\ 1. date
\ 2. source
\ 3. analyzer name
\ 4. delta E
\ 5. high voltage
\ 6. coarse gain
\ 7. fine gain.
\ The entered delta E is also converted to a number and is
\ stored in the (51,2) position of array res.
: specifications. file
file. template
7 comments
real dim[ 51 , 2 ] subfile
end
." Enter name of output file (including path)." cr
"input namefile ":=
namefile defer> file.create
80
namefile defer> file.open p
." Enter date (format: xx/xx/xx) ." cr
"input 1 >comment p
." Enter source (format: Am or Cs)." cr
." The TC 450 and the NC-22 are normally used 
for Am." cr
." The NA-22 is used for Cs."
"input 2 >comment p
." When 'the chosen source is Am, the 
toggle switch of" cr
."' the NA-22 should be set to track. 
When the Cs is" cr
." the source, the switch should be 
on disable." cr cr
." Press any key when ready." s
." Enter name of analyzer." cr
." Format: NA-22, NC-22 or TC 450." cr
"input 3 >comment p
." Set delta E on the analyzer." cr
." For the Cs spectrum a delta E of 0.2 volts 
and for" cr
." the Am spectrum a delta E of 0.1 volts is" cr
." recommended." cr
." When in the 'window asym' mode on the TC 
450, the" cr
." window (delta E) has a range of 0 - 1 volts." cr
." A setting of 1.00 volts on the top knob means 
a" cr
." delta E of 0.1 volts." cr cr
." Press any key when ready." s
." Enter delta E (format x.xx)." cr
"input "dup 4 >comment p
32 "number ?drop "drop
res [ 51 , 2 ] :=
res [ 51 , 2 ] deltaE :=
." Enter voltage on NV-25A" cr
"input 5 >comment p
." Enter coarse gain on amplifier NA-17" cr
"input 6 >comment p
." Enter fine gain on amplifier NA-17" cr
"input 7 >comment p
file.close ;
\ The counting time is entered and assigned to the scalar
\ counting.time.
: time.count
." Enter counting time in seconds. (format: 10.)" cr
." A counting time of 10 seconds is enough for a" cr
." spectrum determination." cr
#input 1000 * counting.time := p
stack.clear ;
\ **********************************************************
81
\ The array res is copied into the 
subfile of the datafile.
.\ The dimensions of the array have 
to be the same as the
\ subfile.
\ **********************************************************
: create. arrays
namefile defer> file.open
1 subfile
res file>array
file.close
\ This macro initializes the plotter.
: plotter cr
." Put paper in the HP 7475." cr
." Press any key when HP is ready to plot data." 
s
hp7475 device.init
7.8 10.15 plotter.size
normal.coords
\ The plot of the-Cs-137 spectrum is made by colon
\ definition 'GRAPHCS'. The plot of the Am-241 
spectrum,
\ obtained with the NC-22 analyzer, is prepared 
with the
\ macro 'GRAPHAMNC22'. The graph of the Am-241 spectrum,
\ obtained with the TC 450, is made by definition
\ 'GRAPHAMTC450'.
: graphCs
4. dumx [ 1 ] :=
7.99 dumx-[(-2 ] :=
.0 dumy [ 1 ] :=
2499 dumy [- 2 ] =
horizontal grid.off 0 1 label.points label.scale.off
vertical grid.off label.scale.off
.2 .25 axis.orig
.7 .65 axis.size normal.coords
.2 .25 axis.point
1. 1. tick.just
.01 .01 tick.size
8 5 axis.divisions
82
count 1.1 + 1 do
res [ i , 1 ] x [ i ] :=
x [ i ] deltaE 2. / + x i :=
res [ i , 2 ] y [ i ] :=
loop
" " symbol dumx dumy xy.auto.plot
" *" solid&symbol
x sub[ 1 , count ]
y sub[ 1 , count ]
xy.data.plot
normal.coords
.55 .14 position " Voltage (V)" centered.label
.55 .07 position " ENERGY SPECTRUM OF Cs-137" centered.label
90 label.dir normal.coords
.08 .573 position " Counts per second" centered.label
0 label.dir normal.coords
.26 .84 position 1 comment> label
.26 .78 position " analyzer " 3 comment> "cat label ;
: graphAmNC22
.1 dumx [ 1 ] :=
.99 dumx [ 2 ] :=
.0 dumy [ 1 ] :=
6999 dumy [ 2 ] :=
horizontal grid.off 0 2 label.points label.scale.off
vertical grid.off 1 1 label.points label.scale.off
.2 .25 axis.orig
.7 .65 axis.size normal.coords
.2 .25 axis.point
1. 1. tick.just
.01 .01 tick.size
9 7 axis.divisions
count 1.1 + 1 do
res i , 1 ] x [ i ] :=
x [ i ] deltaE 2. / + x[ i ] :=
res [ i , 2 ] y [ i ] =
loop
" " symbol dumx dumy 
xy.auto.plot
" * solid&symbol
x sub[ 1 , count ]
y sub[ 1 , count ]
xy.data.plot
normal.coords
.55 .14 position " Voltage (V)" centered.label
.55 .07 position " ENERGY SPECTRUM OF Am-241" centered.label
90 label.dir normal.coords
.08 .573 position " Counts per second" centered.label
0 label.dir normal.coords
.26 .84 position 1 comment> label
.26 .78 position " analyzer " 3 comment> "cat label ;
: graphAmTC450
83
.1 dumx [ 1 ] :=
.99 dumx [ 2 ] :=
.0 dumy [ 1 ] :=
19999 dumy [ 2 ] :=
horizontal grid.off 0 2 label.points label.scale.off
vertical grid.off no.labels label.scale.off
.2 .25 axis.orig
.7 .65 axis.size normal.coords
.2 .25 axis.point
i. 1. tick.just
.01 .01 tick.size
9 7 axis.divisions
count 1.1 + 1 do
res [ i , 1 ] x [ i ] :=
x [ i ] deltaE 2. / + x [ i ] :=
res [ i , 2 ] y [ ] :=
loop
" " symbol dumx dumy xy.auto.plot
" *" solid&symbol
x sub[ 1 , count ]
y sub[ 1 , count ]
xy.data.plot
normal. coords
.181 .257 position " 0" label
.141 .415 position " 5000" label
.125 .572 position " 10000" label
.125 .74 position " 15000" label
.125 .9 position " 20000" label
.55 .14 position " Voltage (V)" centered.label
.55 .07 position " ENERGY SPECTRUM OF Am-241" centered.label
90 label.dir normal.coords
.08 .573 position " Counts per second" centered.label
0 label.dir normal.coords
.26 .84 position 1 comment> label
.26 .78 position " analyzer " 3 comment> "cat label ;
\ Depending on the information stored in comment 3
\ (analyzer name) of the data file, a plot is constructed.
\ The comment is compared with an analyzer name.
\ When the two are the same, a True value is left on the
\ symbol stack and the commands between the 'if' and
\ 'then' words are carried out. If the comment is
\ not the same, a False value is returned and the program
\ skips to the first word after 'then'.
: graph
create. arrays
plotter
84
namefile defer> file.open
3 comment> " nc-22" "=
if graphAmNC22 then
3 comment> " NC-22" '"=
if graphAMNC22 then
3 comment> " na-22" "=
if graphCs then
3 comment> " NA-22" "=
if graphCs then
3 comment> " TC 450" "=
if graphAmTC450 then
3 comment> " tc 450" "=
if graphAmTC450 then
file.close
\ The actual data collection is done in this colon
\ definition. The 'P' command stores the string with all
\ the information in the buffer.
\ The 'gpib.read' word reads the string into the memory of
\ the computer. The string is analyzed and the
\ useful parts are stored in the scalars
\ 'count.time' and 'count.rate'. The maximum number of
\ measurements is 50.
\ After each measurement the user 
can quit.
: spectrum. loop
51 1 do
." Measurement " i . cr cr
" Set E on chosen analyzer" cr cr
i 1=
if
." When determining the spectrum with the NA-22 the" cr
." first E setting should be 3.9. For the increments" cr
." a value of 0.1 is advised." cr
." A spectrum determination with the TC 450 or the" cr
." NC-22 the first E setting should be 0.05 and the" cr
" increments 0.02." cr cr
then
." Press any key when ready." s
." Enter E (format x.xx)." cr
#input res [ i , 1 ] := ?drop p
." Press any key when system is ready to do measurement." cr
." Wait for the next prompt." s
gamma listener me talker
27 ASCII" TALK
counting.time 1000 + msec.delay
gamma listener me talker
85
" P" 3 ASCII" "CAT TALK
4 1 do
gamma si gpib.read
loop
sl
5 1 do
32 "number ?drop
loop
dummy :=
rate :=
dummy : =
count.time :=
"drop
rate count.time / 10. * res [ i , 2 ] :=
." Stored data" cr cr
." E: " res [ i , 1 ] . cr
." Count rate (cps): " res [ i , 2 ] . cr cr
." Enter 'q' if you want to stop taking measurements." cr
." Press enter to continue." cr
"input
" q" "= not p
if
else
i res [ 51 , 1 ] :=
i count :=
namefile defer> file.open cr cr
1 subfile
res array>file
file.close
." If you don't want to determine other energy spectra," cr
." reconnect the analyzers to the multiscalar in the" cr
." following way:" cr
." NA-22 to channel A, pos" cr
." TC 450 to channel B, pos" cr
." NC-22 to channel C, pos (on the back)." cr cr
" Press any key when ready." s
. The disable/track switch of the NA-22 has to be set" cr
." on track." cr cr
" Press any key when ready." s
graph
." Determine lower level E (LL E) and delta E from the" cr
" plot." cr
." Set these E and delta E values on the analyzer." cr
." The mode of the TC 450 should be in 'norm'. The" cr
." 'lower level' knob corresponds with E on the other" cr
." analyzers. The 'upper level or 'window' knob should" cr
." be set to be the sum of LL E and delta E." cr cr
." Both the NC-22 and the NA-22 have to be set in" cr
." integral mode." cr cr
." Press any key when ready." pckey drop ?drop screen.clear
quit
then
86
loop;
\The main program.
\Entering 'spectrum.determination' will execute the
\program.
*spectrum.determination p
pre. comments
specifications.file
time. count
spectrum. loop
87
\ ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ PROGRAM REL.TXT
\ Calibration program 2
\ This program determines the linear relation between
\ integral Am-241 count rates obtained via the NC-22
\ single channel analyzer and the TC 450 single channel
\ analyzer.
\ The vertical stepper motor is used to send the gamma
\ beam through different positions of a teflon prism.
\ With each step the path length through the teflon
\ increases.
\ Minimum ASYST overlay requirements on custom version:
\ gpib master,
\ type 1 nec gpib drive,
\ data files
\ HP plotter driver
\ statistics and special functions
\ polynomial operations
\ matrix and curve fit functions.
\ At the c:\prompt enter c:\asyst\total in order to get
\ into a correct ASYST environment. All the overlays
\ mentioned above are part of the custom version 'total'.
\ Then enter load c:\asyst\rel.txt To run the
\ program enter relation.determination.
\ 07/26/89
\ M.Oostrom
\ ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ ******************************************* *********
\ Variable declaration.
\ ********************************************************
dp.integer scalar counting.time
\ Counting time set on the timer/multiscaler
real scalar t
\ Counting time transferred to buffer
real scalar ratetc
\ Calculated count rate via the TC 450
real scalar ratenc
\ Calculated count rate via the NC-22
real dim[ 33 , 2 ] array res
\ The measured count rates are stored in this array.
88
\ The rates via the NC-22 are placed the first column;
\ the rates via the TC 450 in the second.
0. res :=
real scalar dummy
\ This dummy variable is used when converting the
\ data string coming from the buffer into numbers.
72 string sl
\ Datastring that contains count rates and measured time.
72 string namefile
\ The output file (including the path)
2 gpib.device gamma
\ The name 'gamma' is assigned to the buffered inter-
\ face.
3 gpib.device motor
\ The name 'motor' is assigned to the Aston 800 actuator.
real dim[ 33 ] array x
\ Array used in the graphics part of the program. It
\ contains the NC-22 count rates.
real dim[ 33 ] array y
\ Array used in the graphics part of the program. It
\ contains the TC 450 count rates.
real dim( 2 ] array dumx
real dim[ 2 ] array dumy
\ The use of the arrays 'dumx' and 'dumy' will keep
\ the plot within reasonable boundaries.
real dim( 2 ] array lin.coef
\ The coefficients of the regression line
9 gpib.device hp3497
\ The HP 3497 is named hp3497.
integer scalar number.days
\ The number of calender days the program is running.
integer scalar day
integer scalar days
\ Day of the month.
dp.integer scalar extra
integer scalar sec
integer scalar min
integer scalar hours
integer scalar month
dp.integer scalar time0
72 string s2
\ The time string 'month:day:hour:min:sec' is read from
\ the hp3497.
\ The colon definition 's' is used after the prompt 'Press
\ any key when ready'.
\ The 'p' macro clears the screen and the cursor will move
\ to the middle of the screen. Macro 's' does the same
\ thing after the user has pressed a key.
89
: s
pckey drop ?drop screen.clear 10 1 do cr loop ;
Sp
screen.clear 10 1 do cr loop ;
\ *******************************w***************************
\ The colon definition 'interface' takes care of the GPIB
\ connections.
\ In this program, the three GPIB devices are the buffered
\ interface, defined as 'gamma' ,the Aston 800, defined
\ as 'motor', and the HP3497, defined as 'hp3497'.
\ The buffered interface requires a 'carriage return'
\ (ASCII character 13) as an EOS (end of string) character.
\ The Aston 800 and the HP3497 require a "cr.lf",
\ catenated with the rest of the string that is sent over
\ as an EOS. The HP3497 accepts capitals only.
: interface
bus.init
send.interface.clear
remote.enable.on
gamma eos.on
13 eos.character ;
\ The Aston800 actuator is polled frequently when the
\ system moves from one location to the other. The moment
\ the sources and detector arrive at a new location, the
\ actuator should return the value 16 to the computer. The
\ next statement of a program can then be executed. For no
\ apparent reason, the actuator does not return the value 16
\ after the system has travelled. To avoid a complete crash
\ of a program, a routine should be inserted that takes over
\ the task of the actuator by producing a 16 after a certain
\ time lapse. The colon definitions 'day.determination' and
\ 'time.fetch' are part of this special routine. This
\ routine is part of the 'aston800.status.check' colon
\ definition.
: day.determination
1 number.days :=
0 extra :=
hp3497 listener me talker
90
" TD" "CR.LF" "CAT TALK
hp3497
s2 gpib.read
s2
3 1 do
ascii : "number ?drop
loop
"drop
day :=
month :=
untalk unlisten
: time. fetch
hp3497 listener me talker
" TD" "CR.LF" "CAT TALK
hp3497
s2 gpib.read
s2
5 1 do
ascii : "number ?drop
loop
32 "number
?drop "drop
sec :=
min :=
hours :=
days :=
month :=
days day =
if
else
number.days 1 + number.days :=
days day :=
number.days 1 - 86400 *
extra :=
then
hours 3600 * min 60 * + sec + extra +
untalk unlisten ;
\ The user has to set up the necessary devices properly in
\ order to determine the relation between the two integral
\ count rates. When this is done, the computer is made the
\ controller of the data acquisition process (me talker).
\ The buffer is the device that accepts command ( gamma
\ listener). A 'control R' (ASCII 18) and a 'control C'
\ (ASCII 3) are sent to the buffer to set the scalers to
\ zero.
91
: pre.comments
cr
." Turn power to NH85 off." cr
." (to the right of the buffered interface)" cr
cr ." Press any key when ready." s
." Set stop/count button of TC 535P on 'count'." cr cr
." Press any key when ready." s
0" Set .01 min/.l1 sec button of TC 535P on '.1 sec'." cr
." This is the 'timebase' of the counting time." cr cr
." Press any key when ready." s
." Preset timer." cr
." (Timer formula: timebase*n*10**m.)" cr
" A counting time of 60 seconds is recommended." cr
." Set timer to '62' if you want a 60 s counting time."
cr cr
." Press any key when ready." s
." Set CH A PRESET to '96'" cr cr
." Press any key when ready." s
." Disconnect the following connections:" cr
*1" 1: NA-22 from channel A of the TC 535P " cr
." 2: TC 450 from channel B of the TC 535P " cr
." 3: NC-22 from channel C on the back of the TC 535P."
cr cr
." Press any key when ready." s
." Connect the NC-22 to channel B (pos) and the TC 450" cr
." to channel A (pos) of the TC 535P." cr cr
." Press any key when ready." s
." Set the TC 450 and the NC-22 on integral mode." cr cr
." Press any key when ready." s
." Check the lower level setting on the TC 450 and" cr
." the E-value setting of the NC-22." cr
." The spectrum determinations (calibration program 1)" cr
." should provide the necessary information." cr cr
." Press any key when ready." s
.Make sure the readout mode of the buffer TC 589" cr
. is on auto recycle and that the output format" cr
" is on normal." cr cr
" Press any key when ready." s
.Make sure power of stepper motor indices is on." cr
. System must be able to move freely. Check this." cr cr
" Press any key when ready." s
.Turn power of NH85 on." cr
." Wait till 'Fn' appears on display actuator." cr
." Press any key when ready." s
interface
gamma listener me talker
18 ASCII" 3 ASCII" "CAT TALK ;
\ The following colon definition performs a status test of the
\ Aston 800. The serial poll check gives information on the
92
\ status of the actuator and the connected stepper motors.
\ When the actuator returns the value 16, the system 
has
\ arrived at its destination. The next command of 
the
\ program can then be executed.
: aston800.status.check
untalk unlisten
time.fetch time0 :=
begin
stack.clear
time. fetch
time0 -
60 >
if
16 16
else
1
then
16=
if
16=
else
serial.poll.enable
me listener motor talker
stack. listen
serial.poll.disable
16 =
then
until
untalk unlisten
\ A template is made up for file maintenance. The data
\ output file to be used in this program will contain 3
\ comments and a 33 by 2 array.
\ The first comment is the date. The second comment is
\ the equation of the regression line, while the value
\ of the correlation coefficient is stored in
\ comment 3.
\ The coefficients of the regression equation are stored in
\ position (33,1) and (33,2).
\ The name of the file is assigned to the string 'namefile'
\ The 'file.create' word creates a file according to the
\ latest 'file.template'.
: specifications.file
93
file.template
3 comments
real dim[ 33 , 2 ] subfile
end
." Enter name of output file (including 
path)." cr
"input namefile ":=
namefile defer> file.create
namefile defer> file.open p
." Enter date (format: xx/xx/xx) ." 
cr
"input 1 >comment
file.close
p ;
\ The counting time as set on 
the timer has to be entered.
: time.count
." Enter counting time in seconds." cr
." (enter this number as an integer)" 
cr
#input 1000 * counting.time := p
stack. clear
\ First, the coefficients of the regression 
line are
\ determined. Then the equation that 
represents the line
\ is stored as comment 2. An analysis 
of variance is
\ performed by the 'regress.anova' 
word. The correlation
\ coefficient is stored as comment 
3.
: regression.anova
-1 6 fix.format
x sub[ 1 , 32 ] y sub[ 1 , 32 ] 1 leastsq.poly.fit
[ 1 ] lin.coef [ 1 ] :=
lin.coef [ 1 ] res [ 33 , 1 ] :=
x sub[ 1 , 32 ] y sub[ 1 , 32 ] 1 leastsq.poly.fit
[ 2 ] lin.coef [ 2 ] :=
lin.coef [ 2 ] res [ 33 , 2 ] :=
namefile defer> file.open
-l 4 fix.format
" y = " lin.coef [ 2 ] "." "cat " + " "cat
lin.coef [ 1 ] "." "cat " x" "cat
2 >comment
-1 6 fix.format
y sub[ 1 , 32 ]
94
x sub[ 1 , 32 ] lin.coef poly[x]
2 regress.anova
[ 1 , 5 ] "."
-1 6 fix.format
3 >comment
file.close
." Coefficients of the regression line y=ax + b," cr
." where x is the integral Am count rate via the" cr
." NC-22 s.c.a. and y is the integral Am" cr
." count rate via the TC 450 s.c.a." cr cr
." a = " res [ 33 , 1 ] . cr
." b = " res [ 33 , 2 ] . cr
." Write these numbers down." cr
" You'll need them in the next part of the " cr
." calibration process." cr cr
" Press any key when ready."
pckey drop ?drop screen.clear
\ Initialization of the plotter.
: plotter p
." Press any key when HP is ready to plot data." cr
." Use quality plotter paper only and make sure the" cr
" error indicator is off." s
hp7475 device.init
7.8 10.15 plotter.size
normal.coords
\ The following colon definition contains all the commands
\ necessary to draw the plot.
: graph
horizontal grid.off no.labels label.scale.off
vertical grid.off no.labels label.scale.off
.25 .25 axis.orig
.54 .7 axis.size
normal .coords
.25 .25 axis.point
1. 1. tick.just
.01 .01 tick.size
9 9 axis.divisions
.0 dumx [ 1 ] :=
95
44999 dumx [ 2 ] :=
.0 dumy [ 1 ] :=
44999 dumy [ 2 ] :=
" " symbol dumx dumy xy.auto.plot
" *" symbol
x sub( 1 , 32 ] y sub[ 1 , 32 ] xy.data.plot
solid
45 ramp 1 - 1000 * dup
lin.coef poly[x] xy.data.plot
normal.coords
.25 .223 position " 0" label
.37 .223 position " 10" label
.486 .223 position " 20" label
.61 .223 position " 30" label
.73 .223 position " 40" label
.217 .25 position " 0" label
.203 .405 position " 10" label
.203 .56 position " 20" label
.203 .715 position " 30" label
.203 .875 position " 40" label
90 label.dir normal.coords
.152 .55 position
" count rate via TC 450 (cps, thousands)"
centered. label
0 label.dir normal.coords
.52 .18 position
" count rate via NC-22 (cps, 
thousands)"
centered. label
.25 .125 position
" RELATION OF THE INTEGRAL AM-241 COUNT "
" RATES MEASURED" "cat label
.25 .1 position
" VIA THE NC-22 SINGLE CHANNEL ANALYZER AND "
" VIA THE" "cat label
.25 .075 position
" TC 450 SINGLE CHANNEL ANALYZER " label
normal. coords
namefile defer> file.open
.3 .9 position 1 comment> label
.69 .62 position " r = " 3 comment> "cat label
.69 .66 position 2 comment> label
file.close
.008 .015 char.size
.703 .635 position " 2" label
\ The data acquisition is done by this macro.
\ '27 ASCII" ' will cause the scalers to count. The 'P'
\ command stores the string with the measured time and
\ the data from the three scalers.
96
\ The word 'gpib.read' reads the string into the computer.
\ A teflon prism is placed between sources and detector.
\ With each measurement the pathlength through the prism,
\ and therefore the attenuation, increases.
Srel.loop
33 1 do
i 1 =
if
" There should be nothing between sources and" cr
" detector during the first measurement." cr cr
" Press any key when system is ready to do measurement." cr
" Wait for the next prompt." s
else
i 2 =
if
." Place teflon prism on a jack on the platform." cr cr
." Press any key when ready." s
. Move sources and detector manually with the stepper" cr
. motor indices or manipulate the jack such that the" cr
. beam travels through the spot on the prism." cr cr
" Press any key when ready." s
. Before continuing make sure to set the stepper" cr
. motor indices back to external mode." cr cr
. Press any key when system is ready to do the 31" cr
.measurements through the teflons prism." s
motor listener me talker
" sm0,za" "cr.lf" "cat talk
50 msec.delay
" sm0,zp" "cr.lf" "cat talk
2500 msec.delay
else
motor listener me talker
" sm0,mr+000100" "cr.lf" "cat talk
50 msec.delay
aston800.status.check
2500 msec.delay
then
then
.1" Measurement " i . ." (out of 32)" cr
gamma listener me talker
27 ASCII" TALK
counting.time 2000 + msec.delay
gamma listener me talker
" P" 3 ASCII" "CAT TALK
4 1 do
gamma sl gpib.read
loop
sl
5 1 do
97
32 "number ?drop
loop
dummy :=
ratenc :=
ratetc :=
t :=
"drop
ratenc t / 10. * dup res [ i , 1 := x ( i ]
ratetc t / 10. * dup res i , 2 : y i ] :=
p ." Stored data" cr cr
." Count rate NC-22: " res [ i , 1 ] . cr
." Count rate TC-450: " res [ i , 2 ] . cr
5000 msec.delay p
loop
motor listener me talker
" sm0,ma+000000" "cr.lf" "cat talk
50 msec.delay
aston800. status. check
2500 msec.delay
." Reconnect the analyzers to the multiscaler in the" cr
0." following way:" cr
." NA-22 to channel A, pos" cr
." TC 450 to channel B, pos" cr
." NC-22 to channel C, pos (on the back)." cr cr
." Press any key when ready." s
." The NA-22 and the NC-22 have to be in integral and" cr
." the TC 450 in differential mode (norm). The" cr
" disable/track switch of the NA-22 has to" cr
." be set on track." cr cr
" Press any key when ready." pckey drop ?drop screen.clear
regression. anova
namefile defer> file.open cr cr
1 subfile
res array>file
file.close
stack.clear
plotter
graph
\ The main program. By entering 'relation.determination'
\ the program will run.
: relation.determination p
interface
day. determination
pre. comments
98
specifications. file
time. count
rel. loop
99
\ ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ PROGRAM COMPTON.TXT
\ Calibration program 3
\ This program determines the relation between the
\ high-energy Cs intensity and the Compton scattering in
\ the Am window. Two brass plates are placed between sources
\ and detector to absorb the radiation emitted by the Am
\ source.
\ The vertical stepper motor is used to send the gamma
\ beam through different positions of a teflon prism.
\ With each step the path length through the prism,
\ and therefore the attenuation, increases.
\ Minimum ASYST overlay requirements on custom version:
\ gpib master,
\ type 1 nec gpib drive,
\ data files
\ HP plotter driver
\ statistics and special functions
\ polynomial operations
\ matrix and curve fit functions.
\ At the c:\ prompt enter c:\asyst\total. The custom version
\ 'total' contains the overlays specified above.
\ Load this program by entering load c:\asyst\compton.txt.
\ To run the program enter compton.determination.
\ 08/01/89
\ M.Oostrom
\ ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ **********************************************************
\ Variable declaration.
\ **********************************************************
dp.integer scalar counting.time
\ Counting time set on the timer/multiscaler
real scalar a
real scalar b
\ Scalars a and b are the coefficients of the linear
\ relation of the integral Am count rates via the
\ NC-22 and the TC 450 single channel analyzers.
\ Relation: y = ax + b, where x is the count rate via
\ the NC-22 and y the count rate via the TC 450.
real scalar t
100
\ Counting time transferred to buffer
real scalar ratetc
\ Counts via TC 450 (Am)
real scalar ratena
\ Counts via NA-22 (Cs)
real scalar ratenc
\ Counts via NC-22 (Am)
real dim[ 33 , 2 ] array res
\ The measured count rates are stored in this array. The
\ rates via the NC-22 are placed in the first column;
\ the rates via the TC 450 in the second.
0. res :=
\ initialization of array res
real scalar dummy
\ This dummy variable is used when converting the
\ data string coming from the buffer into numbers.
72 string sl
\ Data string that contains count rates and measured time.
72 string namefile
\ The outputfile (including the path)
2 gpib.device gamma
\ The name 'gamma' is assigned to the buffered interface
3 gpib.device motor
\ The name 'motor' is assigned to the Aston 800 actuator.
real dim[ 33 ] array x
\ Array used in the graphics part of the program.
\ It contains the NC-22 count rates.
real dim[ 33 ] array y
\ Array used in the graphics part of the program.
\ It contains the TC 450 count rates.
real dim[ 2 ] array dumx
\ The use of the arrays 'dumx' and 'dumy' will keep
\ the plot within reasonable boundaries.
real dim[ 2 ] array dumy
real dim[ 2 ] array lin.coef
\ The coefficients of the regression line.
9 gpib.device hp3497
\ The HP 3497 is named hp3497.
integer scalar number.days
\ The number of calender days the program is running.
integer scalar day
integer scalar days
\ Day of the month.
dp.integer scalar extra
integer scalar sec
integer scalar min
integer scalar hours
integer scalar month
dp.integer scalar time0
72 string s2
\ The time string 'month:day:hour:min:sec' is read from
\ the hp3497.
101
\ The colon definition 's' is used after the prompt 'Press
\ any key when ready'. The 'p' macro clears the screen and
\ the cursor will move to the middle of the screen.
\ Macro 's' will do the same after the user has pressed a
\ key.
: s
pckey drop ?drop screen.clear 10 1 do cr loop ;
:p
screen.clear 10 1 do cr loop ;
\ The colon definition 'interface' takes care of the GPIB
\ connections.
\ In this program, the three GPIB devices are the buffered
\ interface, defined as 'gamma' ,the Aston 800, defined
\ as 'motor,' and the HP3497, defined as 'hp3497'.
\ The buffered interface requires a 'carriage return'
\ (ASCII character 13) as an EOS (end of string) character.
\ The Aston 800 and the HP3497 require a "cr.lf",
\ catenated with the rest of the string that is sent over
\ as an EOS. The HP3497 accepts capitals only.
: interface
bus.init
send.interface.clear
remote.enable.on
gamma eos.on
13 eos.character ;
\ The Aston800 actuator is polled frequently when the
\ system moves from one location to the other. The moment
\ the sources and detector arrive at a new location, the
\ actuator should return the value 16 to the computer. The
\ next statement of a program can then be executed. For no
\ apparent reason, the actuator does not return the value 16
\ after the system has travelled. To avoid a complete crash
\ of a program, a routine should be inserted that takes over
\ the task of the actuator by producing a 16 after a certain
\ time lapse. The colon definitions 'day.determination' and
\ 'time.fetch' are part of this special routine. This
\ routine is part of the 'aston800.status.check' colon
102
\ definition.
: day.determination
1 number.days :=
0 extra :=
hp3497 listener me talker
" TD" "CR.LF" "CAT TALK
hp3497
s2 gpib.read
s2
3 1 do
ascii : "number ?drop
loop
"drop
day :=
month :=
untalk unlisten
: time. fetch
hp3497 listener me talker
" TD" "CR.LF" "CAT TALK
hp3497
s2 gpib.read
s2
5 1 do
ascii : "number ?drop
loop
32 "number
?drop "drop
sec :=
min :=
hours :=
days :=
month :=
days day =
if
else
number.days 1 + number.days :=
days day :=
number.days 1 - 86400 *
extra :=
then
hours 3600 * min 60 * + sec + extra +
untalk unlisten ;
\ The user has to set up the necessary devices properly in
\ order to determine the relation between the high-energy
\ Cs and the compton scattering in the Am window.
103
\ When this is done, the computer is made the controller of
\ the data acquisition process (me talker). The buffer is
\ the device that accepts command ( gamma listener). A
\ 'control R' (ASCII 18) and a 'control C' (ASCII 3) are
\ sent to the buffer to set the scalers to zero.
: pre.comments cr
." Check the following connections:" cr
." NA-22 integral output to channel A (pos) of the TC 535P" cr
0" multiscaler;" cr
." TC 450 norm mode output to channel B (pos) of the" cr
." multiscaler;" cr
." NC-22 integral output to channel C (pos) on the" cr
." back of the multiscalar." cr cr
." Press any key when ready." s
." Turn power to NH85 off." cr
." (to the right of buffered interface) "
cr cr ." Press any key when ready." s
." Set stop/count button of TC 535P on 'count'." cr
." Press any key when ready." s
." Set .01 min/.l sec button of the TC 535P on '.1 sec'." cr
." This is the 'timebase' of the counting time."
cr cr
." Press any key when ready." s
." Preset timer." cr
" (Timer formula: timebase*n*10**m.)" cr
" A counting time of 60 s is recommended." cr
." Set timer to '62' if you want a 60 s counting time." cr cr
." Press any key when ready." s
." Set CH A PRESET to '96'" cr cr
." Press any key when ready." s
. Make sure the readout mode of the buffer TC 589 is"
cr ." on auto recycle and that the output format is"
cr ." on normal." cr cr
." Press any key when ready." s
." Make sure power of stepper motor indices is on." cr
." System must be able to move freely. Check this." cr cr
." Press any key when ready." s
." The coefficients of the linear relationship"
cr ." ( y = ax + b) between the integral Am-241"
cr ." count rate via the NC-22 (x) and the integral" cr
." count rate via the TC 450 (y) have to be entered." cr
" Enter a : " cr
#input a := cr
." Enter b : " cr
#input b := p
." Place two brass plates between sources and detector."
cr ." Be careful!!!!!!!!" cr cr
." Press any key when ready." s
." Turn power of NH85 on and wait till the " cr
104
." Aston800 actuator has performed its self test." cr
." The letters 'Fn' appear when the self test has" cr
." been completed." cr cr
." Press any key when counting is finished." s
interface
gamma listener me talker
18 ASCII" 3 ASCII" "CAT TALK ;
\ The following colon definition performs a status test of the
\ Aston 800. The serial poll check gives information on the
\ status of the actuator and the connected stepper motors.
\When the actuator returns the value 16, the system has
\ arrived at its destination. The next command of the
\ program can then be executed.
\ **********************************************************
: aston800.status.check
untalk unlisten
time.fetch time0 :=
begin
stack.clear
time.fetch
time0 -
60 >
if
16 16
else
1
then
16=
if
16=
else
serial.poll.enable
me listener motor talker
stack.listen
serial.poll.disable
16=
then
until
untalk unlisten
\
\ A template is made up for file maintenance. The data
\ output file to be used in this program will contain 3
\ comments and a 33 by 2 array. The first comment is
\ the date. The second comment is the equation of the
105
\ regression line, while the value of the correlation
\ coefficient is stored in comment 3.
\ The coefficients of the regression equation are
\ stored in position (33,1) and (33,2) of the 
array.
\ The name of the file is assigned to the 
string 'namefile'
\ The 'file.create' word creates a file according 
to the
\ latest 'file.template'.
: specifications. file
file.template
3 comments
real dim[ 33 , 2 ] subfile
end
." Enter name of output file (including path)." cr
"input namefile ":=
namefile defer> file.create
namefile defer> file.open p
." Enter date (format: xx/xx/xx) ." cr
"input 1 >comment
file.close
p ;
\ The counting time as set on the timer has to be entered.
: time.count
." Enter counting time in seconds." cr
." (enter this number as an integer)" cr cr
#input 1000 * counting.time := p
stack.clear ;
\ ******************************************************
\ The coefficients of the regression line are determined.
\ Subsequently, the equation that represents the line
\ is stored as comment 2. An analysis of variance is
\ performed by the 'regress.anova' word. The correlation
\ coefficient is stored as comment 3.
: regression. anova
-l 6 fix.format
x sub[ 1 , 32 ] y sub[ 1 , 32 ] 1 leastsq.poly.fit
[ 1 ] lin.coef [ 1 ] :=
lin.coef [ 1 ] res [ 33 , 1 ] :=
x sub[ 1 , 32 ] y sub[ 1 , 32 ] 1 leastsq.poly.fit
106
2 ] lin.coef [ 2 ] :=
lin.coef [ 2 ] res [ 33 , 2 ] :=
." The coefficients of the linear relationship"
cr ." (y = px + q) between the high energy Cs"
cr ." intensity (x) and the Compton scattering in the"
cr ." Am-window (y):" cr
0" p = " lin.coef [ 1 ] . cr
." q = " lin.coef [ 2 ] . cr
." Write these numbers down. You'll need them in later"
cr ." calibration programs and measurements." cr cr
." Press any key when ready." s
namefile defer> file.open
-1 5 fix.format
" y =" lin.coef [ 1 ] "." "cat
" x" "cat " +" "cat lin.coef [ 2 ">comment
y sub[ 1 , 32 ]
x sub[ 1 , 32 ] lin.coef poly[x]
2 regress.anova
[ 1 , 5 ] "."
3 >comment
file.close
\ Initialization of the plotter.
: hp p
." Press any key when HP is ready to plot data." cr
." Use quality plotted paper only and make sure the"
cr ." error indicator is off." s
hp7475 device.init
7.8 10.15 plotter.size
normal.coords ;
\ The following colon definition contains all the commands
\ necessary to draw the plot.
: graph
horizontal grid.off no. labels label. scale.off
vertical grid.off no. labels label.scale.off
.25 .25 axis.orig
.54 .7 axis.size
normal .coords
.25 .25 axis.point
1. 1. tick.just
.01 .01 tick.size
107
8 4 axis.divisions
.0 dumx [ 1 ] :=
4000 dumx [ 2 ] :=
.0 dumy C 1 ] :=
1000 dumy [ 2 ] :=
" " symbol dumx dumy xy.auto.plot
" *" 
symbol
x sub[ 1 , 32 ] y sub[ 1 , 32 ] xy.data.plot
solid
41 ramp 1 - 100 * dup
lin.coef poly[x] xy.data.plot
normal. coords
90 label.dir normal.coords
.125 .55 position
" Compton scattering in the Am-window, "
centered. label
.152 .55 position " counts per second" centered.label
0 label.dir normal.coords
.25 .21 position " 0" centered.label
.385 .22 position " 1000" centered.label
.52 .22 position " 2000" centered.label
.655 .22 position " 3000" centered.label
.79 .22 position " 4000" centered.label
.22 .25 position " 0" label
.193 .425 position " 250" label
.193 .6 position " 500" label
.193 .775 position " 750" label
.179 .95 position " 1000" label
.52 .17 position
" High-energy Cs intensity, counts per second"
centered. label
.25 .125 position
" COMPTON SCATTERING IN THE Am-WINDOW VERSUS" label
.25 .1 position " HIGH-ENERGY Cs INTENSITY." label
.25 .075 position
" (ALL COUNT RATES ARE CORRECTED FOR DEAD TIME) " label
normal.coords
namefile defer> file.open
.3 .9 position 1 comment> label
.69 .62 position " r =" 3 comment> "cat label
.69 .66 position 2 comment> label
file.close
.008 .015 char.size
.703 .635 position " 2" label
\ The data acquisition is done by this macro.
\ '27 ASCII" ' will cause the scalers to count. The 'P'
\ command stores the string with the measured time and
\ the data from the three scalers.
108
\ The word 'gpib.read' reads the string into the computer.
\ A teflon prism is placed between sources and detector.
: compton. loop
33 1 do
i =
if
." There should be nothing between sources and" cr
0" detector." cr cr
." Press any key when system is ready to do" cr
*" the first measurement." cr
." Wait for the next prompt." s
else
i 2
if
." Place teflon prism on a jack on the platform." cr cr
." Press any key when ready." s
." Move sources and detector manually with the stepper" cr
." motor indices or manipulate the jack in such a way" cr
." that the beam travels through the spot on the" cr
." prism." cr cr
." Press any key when ready." s
." Before continuing make sure to set the stepper" cr
." motor indices back to the external mode." cr cr
." Press any key when system is ready to do measurements." s
motor listener me talker
" sm0,za" "cr.lf" "cat talk 50 msec.delay
" sm0,zp" "cr.lf" "cat talk 50 msec.delay
2500 msec.delay
else
motor listener me talker
" sm0,mr+000100" "cr.lf" "cat talk
50 msec.delay
aston800.status. check
2500 msec.delay
then
then
." Measurement " i . ." (out of 32)"
gamma listener me talker
27 ASCII" TALK
counting.time 2000 + msec.delay
gamma listener me talker
" P" 3 ASCII" "CAT TALK
4 1 do
gamma si gpib.read
loop
si
5 1 do
32 "number ?drop
loop
109
ratenc :=
ratetc :=
ratena :=
t :=
"drop
ratetc t / 10. *
ratenc t / 10. * a * b + .0000006 * 1. - / abs
dup res [ i , 2 ] := y [ i ] :=
ratena t / 10. * dup
.00000167 * 1. - abs /
dup res [ i , 1 ] := x [ i ] :=
p ." Stored data" cr cr
." Dead time corrected Cs count rate:" res [ i , 1 ] . cr
." Dead time corrected Compton scatter in Am-window "
res [ i , 2 ] . cr
5000 msec.delay p
loop
motor listener me talker
" sm0,ma+000000" "cr.lf" "cat talk
50 msec.delay
aston800. status. check
2500 msec.delay
regression. anova p
namefile defer> file.open cr cr
1 subfile
res array>file
file.close
stack. clear
hp
graph
\ The main program. By entering 'compton.determination',
\ the program will run.
: compton.determination p
interface
day. determination
pre. comments
specifications.file
time. count
compton. loop
110
\ +++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ PROGRAM ATTEN.TXT
\ Calibration program 4
\ This program determines mass attenuation coefficients
\ of solutions.
\ A Plexiglas cell with four compartments is
\ used in this part of the calibration process.
\ Count rates will subsequently be obtained on the cell
\ with 0, i, 2, 3 and 4 compartments filled with the
\ solution. In each case, the count rates are determined
\ at five different locations. A total of 25 measure-
\ ments is done. With a linear regression analysis,
\ the products of the mass attenuation coefficient and
\ the density of the solution are calculated. The
\ variances of these quantities are also computed.
\ For liquids that dissolve Plexiglas, a glass
\ compartmental cell consisting of four French squared
\ bottles, glued to a wooden base, is used.
\ Minimum ASYST overlay requirements on custom version:
\ gpib master
\ type 1 nec gpib drive
\ HP plotter driver
\ statistics and special functions
\ polynomial operations
\ matrix and curve fit functions
\ At the C:\ prompt, enter c:\asyst\total in order to get
\ into a correct ASYST environment. Then enter
\ load c:\asyst\atten.txt. To run the program, enter
\ determination.attenuation.coef.
\ 08/18/89
\ M.Oostrom
\ +++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ variable declaration
integer scalar m
integer scalar 1
real scalar a
real scalar b
real scalar r
111
real scalar q
\ The coefficients of the linear relation
\ between the integral count rates via the NC-22 and
\ the TC 450 are a and b.
\ (calibration program 2)
\ The coefficients r and q determine the linear relation
\ between the high energy Cs count rate and the Cs
\ Compton scattering in the Am-window.
\ (calibration program 3)
real scalar ratenc
\ Raw integral Am count rate measured via the NC-22.
real scalar ratetc
\ Raw differential Am count rate measured via the TC 450.
real scalar ratena
\ Raw integral Cs count rate measured via the NA-22.
real dim[ 25 ] array yam
\ Natural log values of the corrected Am count rates of
\ the twenty five measurements.
real scalar yamav
\ Average value of the array yam
real dim[ 25 ] array ycs
\ Natural log values of the corrected Cs count rates of
\ the twenty five measurements.
real scalar ycsav
real dim[ 25 ] array x
\ Cumulative travel distances through the cells.
: init
6 1 do 0. x [ i ] := loop
11 6 do 1.505 x [ i ] := loop
16 11 do 2.975 x [ i ] := loop
21 16 do 4.445 x [ i ] := loop
26 21 do 5.95 x [ i ] := loop ;
\ Initialization of the x array.
real scalar xav
\ The average value of array x.
real scalar counting.time
\ Counting time entered by the user.
real scalar count.time
\ Counting time send over from buffer to computer.
real scalar rateam
\ Dead time and compton scattering corrected Am count rate.
real scalar ratecs
\ Dead time corrected Cs count rate.
72 string sl
\ String with the counting time and total counts stored in
\ the buffered interface and send over to the computer.
2 gpib.device gamma
\ The buffered interface is named gamma.
3 gpib.device motor
\ The buffered interface TC 589 is called motor.
real dim[ 2 ] array am.coef
\ Coefficients of the linear relation between the path
112
\ length through the compartmental 
cell and the natural
\ log of count rate of Am.
real dim[ 2 ] array cs.coef
\ Coefficients of the linear relation 
between the path
\ length through the compartmental 
cell and the natural
\ log of count rate of Cs.
real dim[ 2 ] array dumx
\ Dummy array used in the graphics part 
of the program.
real dim[ 2 ] array dumy
\ Dummy array used in the graphics 
part of the program.
10 string rsquaredam
10 string rsquaredcs
\ String containing the r
2 
value of the linear
\ relation mentioned above for Am 
and Cs, respectively.
real dim[ 2 ] array cs.coef
20 string information
\ Information on the solution that has to 
appear
\ on the plot.
real scalar varuwam
real scalar varuwcs
\ The variance of the product of the attenuation
\ coefficient and the density of the 
solution.
\ (Am and Cs, respectively)
real scalar pp
real scalar qq
real scalar rr
real scalar ss
real scalar tt
\ The temporary variables pp,qq,rr,ss, 
and tt are used
\ during the calculation of the variances.
9 gpib.device hp3497
\ The HP 3497 is named hp3497.
integer scalar number.days
\ The number of calender days the program 
is running.
integer scalar day
integer scalar days
\ Day of the month.
dp.integer scalar extra
integer scalar min
integer scalar sec
integer scalar hours
integer scalar month
dp.integer scalar time0
72 string s2
\ The time string 'month:day:hour:min:sec' is read from
\ the hp 3497
\ This colon definition takes care of 
the GPIB connections.
\ In this program, the three GPIB devices 
are the buffered
\ interface, defined as 'gamma', the Aston 
800, defined
113
\ as 'motor', and the HP 3497, defined as 'hp3497'.
\ The buffered interface requires a 'carriage 
return'
\ (ASCII character 13) as an EOS (end of string) 
character.
\ The Aston 800 and the HP 3497 require a "cr.lf",
\ catenated with the rest of the string that is 
sent over
\ as an EOS. The HP 3497 accepts capitals only.
\ **********************************************************
: interface
bus.init
send. interface.clear
remote.enable.on
gamma eos.on
13 eos.character ;
\ The Aston800 actuator is polled frequently when the
\ system moves from one location to the other. The moment
\ the sources and detector arrive at a new location, the
\ actuator should return the value 16 to the computer. The
\ next statement of a program can then be executed. For no
\ apparent reason, the actuator does not always return the
\ value 16 after the system has travelled. To avoid 
a
\ complete crash of a program, a routine should be inserted
\ that takes over the task of the actuator by producing a
\ 16 after a certain time lapse. The colon definitions
\ 'day.determination' and 'time.fetch' are part of this
\ special routine. This routine is part of
\ the 'aston800.status.check' colon definition.
: day.determination
1 number.days :=
0 extra :=
hp3497 listener me talker
" TD" "CR.LF" "CAT TALK
hp3497
s2 gpib.read
s2
3 1 do
ascii : "number ?drop
loop
"drop
day :=
month :=
untalk unlisten ;
: time.fetch
hp3497 listener me talker
114
" TD" "CR.LF" "CAT TALK
hp3497
s2 gpib.read
s2
5 1 do
ascii : "number ?drop
loop
32 "number
?drop "drop
sec :=
min :=
hours :=
days :=
month :=
days day =
if
else
number.days 1 + number.days :=
days day :=
number.days 1 - 86400 *
extra :=
then
hours 3600 * min 60 * + sec + extra +
untalk unlisten ;
\ The colon definition 's' is used after the prompt
\ 'Press any key when ready'. The 'p' macro clears the
\ screen and moves the cursor to the middle of the screen.
\ Macro 's' does the same after the user has pressed a key.
: s
pckey drop ?drop screen.clear 10 1 do cr loop ;
: p
screen.clear 10 1 do cr loop ;
\ Before the actual measurements can take place, the user has
\ to check a number of connections and settings. He also has
\ to enter coefficients of relations determined in previous
\ calibration programs.
: instructions
." Preset timer TC 535P to '1,1'" cr cr
115
." Press any key when ready." s
." Check the following connections with the multiscaler:" cr
." NA-22 integral output to channel A,pos;" cr
." TC 450 norm mode output to channel B,pos;" cr
." NC-22 integral output to channel C,pos." cr cr
." Press any key when ready." s
." Check E and delta E settings of the " cr
." single channel analyzers and the NA-22." cr cr
." Press any key when ready." s
." Turn power to NH85 off." cr
." (to the right of the buffered interface)" cr
cr ." Press any key when ready." s
." Set stop/count button of TC 535P on 'count' and" cr
." set .01 min/.l sec button on '.1 sec'." cr
I" This is the 'timebase' of the counting time." cr cr
." Press any key when ready." s
." Make sure the readout mode of the buffered inter-" cr
." face TC 589 is on auto recycle and that the output " cr
." format is on normal." cr cr
." Press any key when ready." s
." The coefficients of the linear relation (y=ax+b)" cr
." between the integral Am count rate via the NC-22" cr
." (x) and via the TC 450 (y) must be entered." cr cr
-1 6 fix.format
." Enter a:" cr
#input a := cr
." Enter b:" cr
#input b := p
." The coefficients of the linear relation (y=rx+q)" cr
" between the high energy Cs count rate (x) and the" cr
." Compton scattering in the Am-window (y) must be" cr
" entered." cr cr
." Enter r:" cr
#input r := cr
." Enter q:" cr
#input q := p
." Make sure power of stepper motors is on." cr
." System must be able to move freely. Check this." cr cr
." Press any key when ready." s
." Turn power of NH85 on." cr
." Wait till 'Fn' appears on the display of" cr
" the Aston actuator." cr cr
" Press any key when ready." s
interface
gamma listener me talker
18 ASCII" 3 ASCII" "CAT TALK
." Preset timer (minimum 200 sec)." cr
." Timer formula: timebase*n*10**m." cr
." A counting time of 200 s is recommended." cr
." Set the timer to '2,3' if you want a counting time" cr
." of 200 s.
' 
cr cr
." Press any key when ready." s
116
." Set CH A PRESET to '96'. " cr cr
" Press any key when ready." s
." Enter counting time in seconds." cr
" (enter as a real)" cr
#input 1000. * counting.time := p
stack.clear
\ ** ** 
************ **** 
* ******
\ The colon definition 'aston800.status.check' performs a
\ status test of the Aston 800. The serial poll check
\ gives information on the status of the actuator and the
\ connected stepper motors. When the actuator returns
\ the value 16, the system has arrived at its destination.
\ The next command of the program can then be executed.
: aston800.status.check
untalk unlisten
time.fetch time0 :=
begin
stack.clear
time.fetch
time0 -
60 >
if
16 16
else
1
then
16 -
if
16 =
else
serial.poll.enable
me listener motor talker
stack.listen
serial.poll.disable
16 =
then
until
untalk unlisten
\ ********************************************************
\ The determination of the attenuation coefficient is done
\ with a four-compartment Plexiglas calibration cell.
\ First, the count rate through the empty cell is determined.
\ Then the user is prompted four times to fill up a
\ compartment. In each case, a measurement is done at five
117
\ different locations.
\ By fitting a line through all 25 data points,
\ the product of the attenuation coefficient and the
\ density of the solution can be determined. This
\ product is the absolute value of the slope of the line.
\ (see next colon definition)
\ ******************************************************
: measurement
0 1 :=
6 1 do
i i =
if
." Place the calibration cell between sources and detector."
cr
." The radiation should go through compartment 1 first."
cr
." Make sure the long side of the calibration cell" cr
." is parallel to the radiation beam." cr cr
." Press any key when ready for measurement." s
else
." Some action is required after the tune!" cr
10000 1000 tune
p
." Fill up compartment " i 1 - . " with the solution."
cr cr
." Press any key when ready for next measurement." s
then
6 1 do
1 1 + 1 :=
im :=
m 1 =
if
motor listener me talker
" sm0,za" "cr.lf" "cat talk 50 msec.delay
" smO,zp" "cr.lf" "cat talk 50 msec.delay
" sml,za" "cr.lf" "cat talk 50 msec.delay
" sml,zp" "cr.lf" "cat talk
2500 msec.delay
then
m 2 =
if
motor listener me talker
" sml,mr+000500" "cr.lf" "cat talk
50 msec.delay
aston800. status. check
2500 msec.delay
then
if
motor listener me talker
118
" sml,mr-001000" "cr.lf" "cat talk
50 msec.delay
aston800.status.check
2500 msec.delay
then
m 4=
if
motor listener me talker
" sml,mr+000600" "cr.lf" "cat talk
50 msec.delay
aston800.status.check
2500 msec.delay
motor listener me talker
" sml,mr-000100" "cr.lf" "cat talk
50 msec.delay
aston800.status.check
2500 msec.delay
motor listener me talker
" sm0,mr+000500" "cr.lf" "cat talk
50 msec.delay
aston800.status. check
2500 msec.delay
then
m 5=
if
motor listener me talker
" smo,mr-001000" "cr.lf" "cat talk
50 msec.delay
aston800.status.check
2500 msec.delay
then
." Measurement " i . ." (out of 5) for this set up" cr
gamma listener me talker
27 ASCII" TALK
\ This command starts the scalers and the timer.
counting.time 2000 + msec.delay
gamma listener me talker
" P" 3 ASCII" "CAT TALK
\ This command transfers the data string with the information
\ on counts and time from the multiscaler to the buffered
\ interface.
4 1 do
gamma sl gpib.read
loop
\ The previous loop transfers the data string
\ from the buffer to the computer.
\ The next loop converts the data string into the counts and
\ counting time.
sl
5 1 do
32 "number ?drop
loop
119
ratenc =
ratetc :=
ratena :=
count.time :=
"drop
\ The raw data are corrected for dead time. The Am count
\ rate is also corrected for low energy Cs in the Am window.
ratetc count.time / 10. *
ratenc count.time / 10. * a * b + .00006 * 1. - / abs
rateam :=
ratena count.time / 10. * dup
.00000167 * 1. abs /
dup ratecs :=
r * q + rateam swap - rateam :=
ratecs In ycs [ 1 ] :=
rateam In yam [ 1 ] :=
p
m 5 =
if
motor listener me talker
" sm0,mr+000600" "cr.lf" "cat talk
50 msec.delay
aston800.status. check
2500 msec.delay
motor listener me talker
" sm0,mr-000100" "cr.1f" "cat talk
50 msec.delay
aston800.status.check
2500 msec.delay
then
loop
loop ;
\ The determination of the TCE attenuation coefficient is done
\ with a setup consisting of four French square glass
\ bottles glued on a wooden base.
\ First, the count rate through the empty setup is determined.
\ Then the user is prompted four times to fill up a
\ bottle. In each case, a measurement is done five times.
\ By fitting a line through all 25 data points,
\ the product of the attenuation coefficient and the
\ density of the solution can be determined. This
\ product is the absolute value of the slope of the line.
\ (see next colon definition)
: measurement.tce
6 1 do
120
i i =
if
" Place the calibration cell 
between sources and detector."
cr
." The radiation should go through 
bottle 1 first."
cr cr ." Press any key when 
ready for measurement." s
else
." Some action is required after 
the tune!" cr
10000 1000 tune
p
." Fill up bottle " i 1 -
. " with the NAPL."
cr cr
." Press any key when ready 
for next measurement." s
then
6 1 do
1 1 + 1 :=
." Measurement " i . ." (out of 5) 
for this set up" cr
gamma listener me talker
27 ASCII" TALK
\ This command starts the scalers 
and the timer.
counting.time 2000 + msec.delay
gamma listener me talker
" P" 3 ASCII" "CAT TALK
\ This command transfers 
the data string with the information
\ on counts and time from the 
multiscaler to the buffered
\ interface.
4 1 do
gamma sl gpib.read
loop
\ The previous loop transfers 
the data string
\ from the buffer to the computer.
\ The next loop converts the data 
string into the counts and
\ counting time.
sl
5 1 do
32 "number ?drop
loop
ratenc :=
ratetc :=
ratena :=
count.time :=
"drop
\ The raw data are corrected for 
dead time. The Am count
\ rate is also corrected for 
low energy Cs in the Am window.
ratetc count.time / 10. *
ratenc count.time / 10. * a * b + 
.0000006 * 1. - / abs
rateam :=
ratena count.time / 10. * dup
.00000167 * 1. - abs /
dup ratecs :=
r * q + rateam swap - rateam :=
ratecs in ycs [ 1 ] :=
121
rateam In yam [ 1 ] :=
p
loop
loop ;
\ This macro calculates the coefficients 
of the linear
\ relation between path length and the natural log of
\ the corrected count rates of Am and Cs.
\ The r
2 
of the relations and the variance of the
\ product of the mass attenuation coefficient 
and the
\ solution density are also determined.
: regression
x yam 1 leastsq.poly.fit
[ 1 ] am.coef [ 1 ] :=
x yam 1 leastsq.poly.fit
[ 2 ] am.coef [ 2 ] :=
x ycs 1 leastsq.poly.fit
[ 1 ] cs.coef [ 1 ] :=
x ycs 1 leastsq.poly.fit
[ 2 ] cs.coef [ 2 ] :=
yam
x am.coef poly[x]
2 regress.anova
1 , 5 ] "." rsquaredam ":=
ycs
x cs.coef poly[x]
2 regress.anova
[ 1 , 5 ] "." rsquaredcs ":=
\ Calculation of the variances.
0. yamav :=
0. ycsav :=
0. xav :=
26 1 do
yamav yam [ i ] + yamav :=
ycsav ycs [ i ] + ycsav :=
xav x [ i ] + xav :=
loop
yamav 25. / yamav :=
ycsav 25. / ycsav :=
xav 25. / xav :=
0. pp :=
0. qq :=
0. rr :=
0. ss :=
0. tt :=
26 1 do
yam [ i ] yamav- 2 ** tt + tt :=
122
ycs [ i ] ycsav - 2 ** ss + ss :=
x [ i ] xav - 2 ** rr + rr :=
yam [ i ] yamav -
x [ i ] xav- * qq + qq :=
ycs [ i ] ycsav -
x [ i ] xav - * pp + pp :=
loop
qq 2 ** qq :=
pp 2 ** pp :=
tt qq rr / - 23 / rr / varuwam :=
ss pp rr / - 23 / rr / varuwcs :=
screen. clear
." Linear relations:" cr cr
" Ln( Am count rate ) = " am.coef [ 1 ] . ." x + "
am.coef [ 2 ] . cr
." Volumetric attenuation coefficient: "
am.coef [ 1 ] abs . " /cm" cr
-1 10 fix.format
." Variance: " varuwam . cr cr
-1 6 fix.format
." Ln( Cs count rate ) = " cs.coef [ 1 ] . ." x + "
cs.coef [ 2 ] . cr
." Volumetric attenuation coefficient: "
cs.coef [ 1 ] abs . " /cm" cr
-1 10 fix.format
." Variance: " varuwcs cr cr
" Write this information down." cr cr
." Press any key when ready to proceed." s ;
\ *******************************************************
\ Initialization of the plotter
: plotter
." Use only quality paper 
to plot graphs." cr
." Make sure error indicator is off." cr
." Press any key when plotter is ready to plot data." s
hp7475
device. init
7.8 10.15 plotter.size
normal.coords ;
\ This colon definition contains the statements
\ necessary to plot the two linear relations in one graph.
: graph
123
." Enter information on the tested solution" cr
" that has to appear on the graph." cr
." ( maximum 20 characters )" cr cr
"input information ":= p
horizontal grid.off no.labels label.scale.off
vertical grid.off no.labels label.scale.off
.25 .25 axis.orig
.54 .7 axis.size
normal, .coords
.25 .25 axis.point
1. 1. tick.just
.01 .01 tick.size
6 3 axis.divisions
.0 dumx [ 1 ] :=
5.99 dumx [ 2 ] :=
7. dumy [ 1 ] :=
9.99 dumy [ 2 ] :=
" " symbol dumx dumy xy.auto.plot
" +" symbol
x yam xy.data.plot
solid
62 ramp 1 - 10 / dup
am.coef poly(x] xy.data.plot
" *" symbol
x ycs xy.data.plot
.02 .02 .02 .02 dashed
62 ramp 1 - 10 / dup
cs.coef poly[x] xy.data.plot
normal, .coords
.217 .25 position " 7.0" centered.label
.217 .483 position " 8.0" centered.label
.217 .717 position " 9.0" centered.label
.203 .95 position " 10.0" centered.label
.25 .223 position " 0.0" centered.label
.43 .223 position " 2.0" centered.label
.61 .223 position " 4.0" centered.label
.79 .223 position " 6.0" centered.label
normal.coords
.52 .97 position information centered.label
.7 .9 position " Am (solid line):" label
-1 5 fix.format
.7 .86 position
" y=" am.coef [ 1 ] "." "cat " x +" "cat
am.coef [ 2 ] "." "cat label
.7 .82 position
" r =" rsquaredam "cat label
.7 .78 position " Cs (dashed line):" label
.7 .74 position
" y=" cs.coef [ 1 ] "." "cat " x +" "cat
cs.coef [ 2 ] "." "cat label
.7 .70 position
" r =" rsquaredcs "cat label
124
90 label.dir
normal. coords
.152 .6 position " LN ( count rate )" centered.label
0 label.dir
normal.coords
.52 .18 position " PATH LENGTH THROUGH FILLED"
" COMPARTMENTS (CM)" "cat centered.label
.25 .125 position
" RELATION BETWEEN PATH LENGTH THROUGH A SOLUTION"
label
.25 .1 position
" AND THE NATURAL LOG OF THE CORRECTED COUNT RATES."
label
.25 .075 position
" THE ABSOLUTE VALUE OF THE SLOPE IS THE"
label
.25 .05 position
" VOLUMETRIC ATTENUATION COEFFICIENT."
label
.008 .015 char.size
normal. coords
.713 .835 position " 2" label
.713 .715 position " 2" label ;
: graph.tce
." Enter information on the tested solution" cr
." that has to appear on the graph." cr
." ( maximum 20 characters )" cr cr
"input information ":= p
horizontal grid.off no.labels label.scale.off
vertical grid.off no.labels label.scale.off
.25 .25 axis.orig
.54 .7 axis.size
normal.coords
.25 .25 axis.point
1. 1. tick.just
.01 .01 tick.size
4 5 axis.divisions
.0 dumx [ 1 ] :=
7.99 dumx [ 2 :=
5. dumy [ 1 ] :=
9.99 dumy [ 2 ] :=
" " symbol dumx dumy 
xy.auto.plot
" +" symbol
x yam xy.data.plot
solid
82 ramp 1- 10 / dup
am.coef poly[x] xy.data.plot
" *" symbol
x ycs xy.data.plot
.02 .02 .02 .02 dashed
82 ramp 1 - 10 / dup
125
cs.coef poly[x] xy.data.plot
normal.coords
.217 .25 position " 5.0" centered.label
.217 .39 position " 6.0" centered.label
.217 .53 position " 7.0" centered.label
.217 .67 position " 8.0" centered.label
.217 .81 position " 9.0" centered.label
.203 .95 position " 10.0" centered.label
.25 .223 position " 0.0" centered.label
.385 .223 position " 2.0" centered.label
.52 .223 position " 4.0" centered.label
.655 .223 position " 6.0" centered.label
.79 .223 position " 8.0" centered.label
normal. coords
.52 .97 position information centered.label
.7 .9 position " Am (solid line):" 
label
-1 5 fix.format
.7 .86 position
" y=" am.coef [ 1 ] "." "cat 
" x +" "cat
am.coef [ 2 ] "." "cat label
.7 .82 position
" r =" rsquaredam "cat label
.7 .78 position " Cs (dashed line):" 
label
.7 .74 position
" y=" cs.coef [ 1 ] "." "cat 
" x +" "cat
cs.coef [ 2 ] "." "cat label
.7 .70 position
" r =" rsquaredcs "cat label
90 label.dir
normal. coords
.152 .6 position " LN ( COUNT 
RATE )" centered.label
0 label.dir
normal.coords
.52 .18 position " PATH LENGTH THROUGH 
FILLED"
" COMPARTMENTS (CM)" "cat centered. label
.25 .125 position
" RELATION BETWEEN PATH LENGTH THROUGH 
A NAPL"
label
.25 .1 position
" AND THE NATURAL LOG OF THE CORRECTED 
COUNT RATES."
label
.25 .075 position
" THE ABSOLUTE VALUE OF THE SLOPE IS 
THE"
label
.25 .05 position
" VOLUMETRIC ATTENUATION COEFFICIENT."
label
.008 .015 char.size
normal .coords
.713 .835 position " 2" label
.713 .715 position " 2" label ;
126
\ The main program.
\ After this program is loaded into the right Asyst
\ environment, entering 'determination.attenuation.coef'
\ will run the program.
\ * ************* 
* ** ***************** * * *
: determination.attenuation.coef p
init
day. determination
instructions p
." Enter 'y' if you want to use the Plexiglass" cr
." compartmental cell." cr cr
"input p
" y,, ,,_
if
measurement
regression p
." Enter 'y' if you want to plot the relation" cr
." Enter 'n' if you want to quit" cr
"input p
of y" 
I,
if
plotter
graph
then
II yg II I
if
plotter
graph
else
quit
then
else
measurement.tce
regression p
." Enter 'y' if you want to plot the relation" cr
." Enter 'n' if you want to quit" cr
"input p
" y" ,=
if
plotter
graph. tce
then
if
plotter
graph. tce
else
quit
127
then
then;
128
\ +++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ PROGRAM SAT5A.TXT
\ Calibration program 5a
\ This program determines, together with program 5b,
\ mass attenuation coefficients of soils by using a cell
\ with 4 compartments.
\ Program 5a takes care of the data acquisition and
\ correction of the raw data for dead time and
\ Compton scattering. Before executing program 5b,
\ which carries out the linear regression anlysis,
\ the volumetric water content and the bulk density
\ of each of the four compartments must be determined
\ gravimetrically.
\ Before starting the determination, all compartments
\ of the cell must be filled with soil. Count rates
\ on 5 different positions will subsequently be obtained
\ on the cell with 4, 3, 2, 1, and 0 compartments filled.
\ This results in 25 data pairs for each source (Am and
\ Cs).
\ Minimum ASYST overlay requirements on custom version:
\ gpib master
\ type 1 nec gpib drive
\ data files
\ At the C:\ prompt, enter c:\asyst\total in order to get
\ into a correct ASYST environment. Then enter
\ load c:\asyst\sat5a.txt. To run the program, enter
\ count.rate.acquisition.
\ 09/18/89
\ M.Oostrom
\ ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ variable declaration
dp.integer scalar 1
real scalar a
real scalar b
real scalar r
real scalar q
\ The coefficients of the linear relation
129
\ between the integral count rates via the NC-22 and
\ the TC 450 are a and b.
\ (calibration program 2)
\ The coefficients r and q determine the linear relation
\ between the high energy Cs count rate and the Cs
\ Compton scattering in the Am-window.
\ (calibration program 3)
real scalar ratenc
\ Raw integral Am count rate measured via the NC-22
real scalar ratetc
\ Raw differential Am count rate measured via the TC 450
real scalar ratena
\ Raw integral Cs count rate measured via the NA-22
real dim[ 25 ] array yam
\ Natural log values of the corrected Am count rates of
\ the twenty five measurements
real dim[ 25 ] array ycs
\ Natural log values of the corrected Cs count rates of
\ the twenty five measurements
real scalar counting.time
\ Counting time entered by the user
real scalar t
\ Counting time send over from buffer to computer
real scalar rateam
\ Dead time and compton scattering corrected Am count rate
real scalar ratecs
\ Dead time corrected Cs count rate
72 string sl
\ String with the counting time and total counts stored in
\ the buffered interface and send over to the computer
2 gpib.device gamma
\ The buffered interface is named gamma
real dim[ 2 ] array am.coef
\ Coefficients of the linear relation to be determined
\ for Am
real dim[ 2 ] array cs.coef
\ Coefficients of the linear relation to be determined
\ for Cs
72 string namefile
\ String containing the path and name of data file
\ with natural log values of corrected Am and Cs count
\ rates.
3 gpib.device motor
\ The Aston 800 is named motor
9 gpib.device hp3497
\ The HP 3497 is named hp3497.
integer scalar number.days
\ The number of calender days the program is running.
integer scalar day
integer scalar days
\ Day of the month.
dp.integer scalar extra
130
integer scalar sec
integer scalar min
integer scalar hours
integer scalar month
dp.integer time0
72 string s2
\ The time string 'month:day:hour:min:sec' is read from
\ the hp3497.
\ The colon definition 's' is used after the prompt
\ 'Press any key when ready'. The 'p' macro clears the
\ screen and moves the cursor to the middle of the screen.
\ Macro 's' does the same after the user has pressed a key.
: s
pckey drop ?drop screen.clear 10 1 do cr loop ;
:p
screen.clear 10 1 do cr loop ;
\ The colon definition 'interface' takes care of the GPIB
\ connections.
\ In this program, the three GPIB devices are the buffered
\ interface, defined as 'gamma' ,the Aston 800, defined
\ as 'motor', and the HP3497, defined as 'hp3497'.
\ The buffered interface requires a 'carriage return'
\ (ASCII character 13) as an EOS (end of string) character.
\ The Aston 800 and the HP3497 require a "cr.lf",
\ catenated with the rest of the string that is sent over
\ as an EOS. The HP3497 accepts capitals only.
: interface
bus.init
send.interface.clear
remote.enable.on
gamma eos.on
13 eos.character ;
\ **********************************************************
\ The Aston800 actuator is polled frequently when the
\ system moves from one location to the other. The moment
\ the sources and detector arrive at a new location, the
\ actuator should return the value 16 to the computer. The
131
\ next statement of a program can then 
be executed. For no
\ apparent reason, the actuator does not 
return the value 16
\ after the system has travelled. To avoid 
a complete crash
\ of a program, a routine should 
be inserted that takes over
\ the task of the actuator by producing 
a 16 after a certain
\ time lapse. The colon definitions 
'day.determination' and
\ 'time.fetch' are part of this special 
routine. This
\ routine is part of the 'aston800.status.check' 
colon
\ definition.
: day.determination
1 number.days :=
0 extra :=
hp3497 listener me talker
" TD" "CR.LF" "CAT TALK
hp3497
s2 gpib.read
s2
3 1 do
ascii : "number ?drop
loop
"drop
day :=
month : =
untalk unlisten
: time. fetch
hp3497 listener me talker
" TD" "CR.LF" "CAT TALK
hp3497
s2 gpib.read
s2
5 1 do
ascii : "number ?drop
loop
32 "number
?drop "drop
sec :=
min :=
hours :=
days :=
month :=
days day =
if
else
number.days 1 + number.days :=
days day :=
number.days 1 - 86400 *
extra :=
then
132
hours 3600 * min 60 * + sec + extra +
untalk unlisten
\ Before the actual measurements can take place, the user has
\ to check a number of connections and settings. He also has
\ to enter coefficients of relations determined in previous
\ calibration programs.
: instructions
." Preset timer to '1,1'" cr cr
." Press any key when ready." s
" Check the following connections with the multiscaler:" cr
." NA-22 integral output to channel A,pos;" cr
." TC 450 norm mode output to channel B,pos;" cr
." NC-22 integral output to channel C,pos." cr cr
." Press any key when ready." s
" Check E and delta E settings of the " cr
." single channel analyzers and the NA-22." cr cr
." Press any key when ready." s
" Turn power to NH85 off." cr
." (to the right of buffered interface TC 589)" cr
cr ." Press any key when ready." s
." Set stop/count button of TC 535P on 'count' and" cr
." set .01 min/.l sec button on '.1 sec'." cr
." The .1 sec is the 'timebase' of the counting time."
cr cr ." Press any key when ready." s
." Make sure the readout mode of the buffered inter-" cr
." face is on auto recycle and that the output " cr
." format is on normal." cr cr
." Press any key when ready." s
." The coefficients of the linear relation (y=ax+b)" cr
." between the integral Am count rate via the NC-22" cr
." (x) and via the TC 450 (y) must be entered." cr cr
-1 6 fix.format
." Enter a:" cr
#input a := cr
." Enter b:" cr
#input b := p
." Make sure power to the stepper motor indices is on." cr
." System must be able to move freely. Check this." cr cr
." Press any key when ready." s
." The coefficients of the linear relation (y=rx+q)" cr
." between the high energy Cs count rate (x) and the" cr
." Compton scattering in the Am-window (y) must be" cr
." entered." cr cr
." Enter r:" cr
#input r := cr
." Enter q:" cr
133
#input q := p
." Turn power of NH85 on." cr
." Wait till 'Fn' appears on display 
Aston 800." cr cr
." Press any key when ready." s
interface
untalk unlisten
gamma listener me talker
18 ASCII" 3 ASCII" "CAT TALK
." Preset timer (minimum 300 sec)." cr
." Timer formula: timebase*n*10**m." 
cr
." If you want a counting time of 
300 sec, set the timer" cr
" to '3,3'." cr cr
." Press any key when ready." s
" Set CH A PRESET to '96'. " cr cr
." Press any key when ready." s
" Enter counting time in seconds (enter 
as a real)." cr
#input 1000 * counting.time := p
stack. clear
\ The following colon definition performs 
a status test of the
\ Aston 800. The serial poll check gives 
information on the
\ status of the actuator and the connected 
stepper motors.
\ When the actuator returns the value 
16, the system has
\ arrived at its destination. The 
next command of the
\ program can then be executed.
: aston800. status. check
untalk unlisten
time.fetch time0 :=
begin
stack.clear
time.fetch
time0 -
60 >
if
16 16
else
1
then
16 =
if
16 =
else
serial .poll. enable
me listener motor talker
stack. listen
serial .poll. disable
16 =
134
then
until
untalk unlisten
\ A file template is made up for file maintenance. 
The data
\ file will contain two subfiles, each with a storage
\ capacity of a 25 element array.
\ The name of the file is assigned to the 
string 'namefile'.
\ The 'file.create' word creates a file according 
to the
\ latest 'file.template'.
: specifications.file
file.template
real dim[ 25 ] subfile
2 times
end
." Enter name of data file (including path)." cr
"input namefile ":=
namefile defer> file.create
p
\ The determination of the attenuation coefficient 
is done
\ with a four-compartment Plexiglas calibration 
cell.
\ First, the count rate through the cell with 
all
\ compartments filled with soil is determined. These
\ count rates are obtained five times at different positions.
\ Source and detector are moved to the new positions 
under
\ computer control. This procedure is repeated four more
\ times, each time the user is prompted to empty one of 
the
\ compartments into a moisture can.
\ The corrected count rates are stored in a data
\ file. This data file will be used in program 5b 
to
\ determine the attenuation coefficient of the soil.
: measurement
26 1 :=
6 1 do
if
." Place the calibration cell between sources and detector."
cr
." All 4 compartments have to be filled with soil." 
cr
135
." Have four moisture cans ready." cr cr
" Press any key when ready." s
." The radiation should go through compartment 1 first."
cr
." Make sure the long side of the calibration cell" 
cr
" is parallel to the radiation beam." cr cr
." Press any key when ready for measurement." s
else
." Some action is required after the tune!" 
cr
10000 1000 tune
p
" Empty compartment " 6 i - " 
into a moisture can,"
cr ." and return Plexiglas cell." cr cr
." Press any key when ready for next measurement." 
s
then
6 1 do
1 1 - 1 :=
i 1 =
if
motor listener me talker
" sm0,za" "cr.lf" "cat talk 50 msec.delay
" sm0,zp" "cr.lf" "cat talk 50 msec.delay
" sml,za" "cr.lf" "cat talk 50 msec.delay
" sml,zp" "cr.lf" "cat talk 50 msec.delay
2500 msec.delay
then
i 2
if
motor listener me talker
" sml,mr+000500" "cr.lf" "cat talk
50 msec.delay
aston800. status, check
2500 msec.delay
then
i 3 =
if
motor listener me talker
" sml,mr-001000" "cr.lf" "cat talk
50 msec.delay
aston800.status.check
2500 msec.delay
then
i 4 =
if
motor listener me talker
" sml,mr+000600" "cr.lf" "cat talk
50 msec.delay
aston800. status. check
motor listener me talker
2500 msec.delay
" sml,mr-000100" "cr.lf" "cat talk
50 msec.delay
136
aston800.status.check
motor listener me talker
2500 msec.delay
" sm0,mr+000500" "cr.lf" "cat talk
50 msec.delay
aston800.status. check
motor listener me talker
2500 msec.delay
then
i 5 =
if
motor listener me talker
" sm0,mr-001000" "cr.lf" "cat talk
50 msec.delay
aston8 00.status. check
2500 msec.delay
then
." Measurement " i . ." (out of 5) for this set up" cr
gamma listener me talker
27 ASCII" TALK
\ This command starts the scalers and the timer.
counting.time 2000 + msec.delay
gamma listener me talker
" P" 3 ASCII" "CAT TALK
\ This command transfers the data string with the
\ information on counts and time from the multiscalar
\ to the buffered interface.
4 1 do
gamma sl gpib.read
loop
\ The previous loop transfers the data string
\ from the buffer to the computer.
\ The next loop converts the data string into the counts and
\ counting time.
sl
5 1 do
32 "number ?drop
loop
ratenc :=
ratetc :=
ratena :=
t :=
"drop
\ The raw data are corrected for dead time. The Am count
\ rate is also corrected for low energy Cs in the Am window.
ratetc t / 10. *
ratenc t / 10. * a * b + .0000006 * 1. - / abs
rateam :=
ratena t / 10. * dup
.00000167 * 1. - abs /
dup ratecs :=
r * q + rateam swap - rateam :=
137
ratecs In ycs [ 1 
rateam In yam [ 1 ] :=
p
i 5
if
motor listener me talker
" sm0,mr+000600" "cr.lf" "cat talk
50 msec.delay
aston800, status. check
2500 msec.delay
motor listener me talker
" sm0,mr-000100" "cr.lf" "cat talk
50 msec.delay
aston800. status. check
2500 msec.delay
motor listener me talker
then
loop
loop
namefile defer> file.open
1 subfile
ycs array>file
2 subfile
yam array>file
file.close
\ The main program.
\ After this program is loaded into the right ASYST
\ environment, entering 'count. rate. acquisition'
\ will run the program.
\ *****************************************************
: count.rate.acquisition p
interface
day.determination
instructions
specifications.file
measurement p
." Array of corrected Cs count rates:" cr
." (count rates through empty cell first)" cr cr
ycs exp . cr
." Press any key to continue." s
." Array of corrected Am count rates:" cr
." (count rates through empty cell first)" cr cr
yam exp .
138
\ +++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ PROGRAM SAT5B.TXT
\ Calibration program 5b
\ This program determines, 
together with program 
5a,
\ mass attenuation 
coefficients of soils.
\ In program 5a, 
the corrected count 
rates through the
\ compartmental 
cell are acquired. In 
program 5b, the
\ values of the 
volumetric water contents 
and dry bulk
\ densities of the 
four compartments 
are entered.
\ These data are necessary 
to calculate the 
mass
\ attenuation coefficient 
of the soil. The 
absolute
\ value of the slope, 
obtained with a 
linear regression,
\ equals the mass 
attenuation coefficient 
of the soil.
\ Finally, the variances 
of the mass attenuation
\ coefficients are calculated.
\ Minimum ASYST overlay 
requirements on custom 
version:
\ data files
\ statistics and special 
functions
\ polynomial functions
\ matrix and curve fit 
functions.
\ At the c:\ prompt 
enter c:\asyst\total. 
The custom version
\ 'total' contains 
the overlays specified 
above.
\ Load this program 
by entering load c:\asyst\sat5b.txt.
\ To run the program, 
enter soil.att.calculation.
\ 09/20/89
\ M.Oostrom
\ +++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ variable declaration
dp.real scalar uwam
dp.real scalar uwcs
\ The products of the 
mass attenuation coefficients 
and
\ the density of deionized 
water.
\ ( Am and Cs, respectively)
dp.real scalar v
\ The sum of the products 
of the dry bulk density 
and the
\ travel distance through 
the soil filled compartments.
dp.real scalar z
\ The sum of the products 
of the volumetric water 
content
139
\ and the travel distance through 
the soil filled
\ compartments.
real dim[ 25 ] array x
\ Array containing the sum 
of the products of the dry bulk
\ density and the travel 
distance through the 
soil filled
\ compartments.
real dim[ 25 ] array ycs
real dim[ 25 ] array yam
\ Arrays containing 
the natural log values 
of the corrected
\ count rates, obtained 
in program 5a. In this 
program the
\ product of the mass attenuation 
coefficient, density of
\ water and z are added to 
do the linear regression.
\ ( Cs and Am, respectively ).
real scalar yamav
real scalar ycsav
\ The average of the arrays 
yam and ycs, respectively.
real scalar xav
\ The average of the array x.
real scalar varusam
real scalar varuscs
\ The variance of the mass 
attenuation coefficients.
\ (Am and Cs, respectively)
real scalar p
real scalar q
real scalar r
real scalar s
real scalar t
\ The temporary variables p,q,r,s 
and t are used during the
\ calculation of the variances.
real dim[ 4 ] array y
\ Travel distance through a compartment.
real dim( 4 ] array bulk
\ The dry bulk densities of 
the compartments.
real dim[ 4 ] array theta
\ The volumetric water contents 
of the compartments.
72 string namefile
\ Name of data file (including 
path) used in 5a.
real dim[ 2 ] array am.coef
real dim[ 2 ] array cs.coef
\ Coefficients of the regression 
relation for Am and
\ Cs, respectively.
\ Initialization of the variables.
\ ********************************************************
: initialization
." Enter the product of the mass 
attenuation coefficient"
cr ." and the density of deionized 
water, for Am:" cr
#input uwam := screen.clear 10 
1 do cr loop
140
." Enter the product of 
the mass attenuation 
coefficient"
cr ." and the density 
of deionized water, for 
Cs:" cr
#input uwcs := screen.clear 
10 1 do cr loop
5 1 do
." Enter dry bulk density 
of compartment "i . cr
#input bulk [ i ] := screen.clear 
10 1 do cr loop
." Enter volumetric water 
content of compartment " 
i . cr
#input theta [ i ] := screen.clear 
10 1 do cr loop
loop
0. x :=
0. z :=
0. v :=
1.505 y [ 1 ] :=
1.47 y 2 1:
=
1.47 y [ 3 J :=
1.505 y [ 4 ] :=
\ In the following colon 
definition, the data file 
with the
\ corrected count rates 
(program 5a) is retrieved. 
The 2
\ subfiles are assigned 
to the arrays ycs and yam,
\ respectively.
********************************************************
: retrieve. file
." Enter name of data file 
(including path)," cr
." containing the corrected 
count rates." cr
." This data file was 
created with program 5a."
cr "input namefile ":= screen.clear 
10 1 do cr loop
namefile defer> file.open
1 subfile
ycs file>array
2 subfile
yam file>array
file.close
\ The arrays ycs, yam, 
and x are adapted according 
to the
\ volumetric water contents 
and dry bulk densities
\ of the soil in the four compartments.
\ The changed arrays are stored 
in the same data file.
********************************************************
: change.arrays
y [ 1 ] bulk [ 1 ] * v + v 
:=
theta [ 1 ] y [ 1 ] * z + 
z :=
141
11 6 do
vx C i] :=
uwcs z * ycs [ i ] + ycs 
[ i ] :=
uwam z * yam [ i ] + yam [ i ] 
:=
loop
y [ 2 ] bulk [ 2 ] * v + 
v :=
theta [ 2 ] y [ 2 ] * z + 
z :=
16 11 do
vx [ i :=
uwcs z * ycs [ i ] + ycs [ i 
] :=
uwam z * yam [ i ] + yam [ i ] :=
loop
y [ 3 ] bulk [ 3 ] * v + 
v :=
theta [ 3 ] y [ 3 ] * z 
+ z :=
21 16 do
vx[i] :=
uwcs z * yes [ i ] + ycs [ i ] :=
uwam z * yam [ i J + yam [ i ] :=
loop
y [ 4 ] bulk [ 4 ] * v 
+ v :=
theta [ 4 ] y [ 4 ] * z 
+ z :=
26 21 do
vx[i] :=
uwcs z * ycs [ i ] + ycs [ i :=
uwam z * yam [ i ] + yam [ i ] :=
loop
namefile defer> file.open
1 subfile
ycs array>file
2 subfile
yam array>file
file.close
\ The following colon 
definition calculates 
the mass
\ attenuation coefficients. 
The absolute values of 
the
\ slope of the linear 
relation between x and 
ycs or yam
\ are the mass attenuation 
coefficients of the soil
\ for Cs and Am, respectively.
\ The r
2 
of the relations and 
the variance of the 
mass
\ attenuation coefficients 
are also determined.
: regression
\ Determination of the 
regression coefficients.
-1 6 fix.format
x yam 1 leastsq.poly.fit
[ 1 ] am.coef [ 1 ] :=
x yam 1 leastsq.poly.fit
142
[ 2 ] am.coef [ 2 ] :=
x ycs 1 leastsq.poly.fit
[ 1 ] cs.coef [ 1 ] :=
x ycs 1 leastsq.poly.fit
[ 2 ] cs.coef [ 2 ] :=
\ Calculation of the variances.
0. yamav :=
0. ycsav :=
0. xav :=
26 1 do
yamav yam [ i ] + yamav :=
ycsav ycs [ i ] + ycsav :=
xav x [ i ] + xav :=
loop
yamav 25. / yamav :=
ycsav 25. / ycsav :=
xav 25. / xav :=
0. p :=
0. q :=
0. r :=
0. s :=
0. t :=
26 1 do
yam [ i ] yamav - 2 ** t + t :=
ycs [ i ] ycsav - 2 ** s + s :=
x [ i ] xav - 2 ** r + r :=
yam [ i ] yamav -
x [ i ] xav - * q + q :=
ycs [ i ] ycsav -
x [ i ] xav - * p + p :=
loop
q 2 ** q :=
p 2 ** p :=
t q r / - 23 / r / varusam :=
s p r / - 23 / r / varuscs :=
\ Determination of r squared and print out on screen.
am.coef [ 1 ] abs
screen.clear
." Am mass attenuation coefficient: " .
cr
yam
x am.coef poly[x]
2 regress.anova
." r squared of linear relation: " . cr
-1 10 fix.format varusam
." Variance of mass attenuation coefficient: " . cr cr
-1 6 fix.format cs.coef [ 1 ] abs
." Cs mass attenuation coefficient: " .
cr
ycs
x cs.coef poly[x]
143
2 regress.anova
[1,5]
." r squared of linear relation: 
" . cr
-1 10 fix.format varuscs
." Variance of mass attenuation 
coeffient: " . cr
." Write this information 
down."
-1 6 fix.format
\ The main program.
\ After this program is 
loaded into the right ASYST
\ environment, entering 
'soil.att.calculation' 
will run
\ the program.
: soil.att.calculation
screen.clear 10 1 do cr loop
initialization
retrieve.file
change.arrays
regression
144
\ +++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ PROGRAM X.TXT
\ Calibration program 6
\ This program determines the path lengths (x) at
\ positions along a vertical column where measurements
\ will be taken during an experiment.
\ Count rate determinations are done five times 
at a
\ maximum of 15 locations along the empty and water-filled
\ column. The locations are evenly spaced.
\ The average count rates through the empty column 
are
\ stored for future use.
\ An overshooting technique is used to overcome the
\ hysteresis of the stepper motors and to assure 
that the
\ locations where the measurements are taken are the
\ same for each loop.
\ Minimum ASYST overlay requirements on custom version:
\ gpib master
\ type 1 nec gpib drive
\ data files
\ At the C:\ prompt, enter c:\asyst\total in order 
to get
\ into a correct ASYST environment. Then enter
\ load c:\asyst\x.txt. To run the program, enter
\ path.length.determination.
\ 09/22/89
\ M.Oostrom
\ ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ variable declaration
integer scalar 1
integer scalar k
integer scalar m
integer scalar n
\ 1,k and m are counters in loops
real scalar a
real scalar b
real scalar r
real scalar q
\ The coefficients of the linear relation
\ between the integral count rates via the NC-22 and
145
\ the TC 450 are a and b.
\ (calibration program 2)
\ The coefficients r and 
q determine the linear 
relation
\ between the high energy 
Cs count rate and the Cs
\ Compton scattering in the 
Am-window.
\ (calibration program 3)
real scalar ratenc
\ Raw integral Am count 
rate measured via the NC-22
real scalar ratetc
\ Raw differential Am 
count rate measured 
via the TC 450
real scalar ratena
\ Raw integral Cs count 
rate measured via the 
NA-22
real dim[ 75 ] array yam.empty
real dim[ 75 ] array yam.filled
\ Natural log values of 
the corrected Am count rates
\ (with a maximum of 75) for 
the empty and filled
\ column, respectively.
real dim( 75 ] array ycs.empty
real dim[ 75 ] array ycs.filled
\ Natural log values of the 
corrected Cs count rates
\ (with a maximum of 75) for 
the empty and filled
\ column, respectively.
dp.integer scalar counting.time
\ Counting time entered by the user
real scalar t
\ Counting time sent over 
from buffer to computer
real scalar rateam
\ Dead time and Compton 
scattering corrected Am count 
rate
real scalar ratecs
\ Dead time corrected Cs count 
rate
72 string sl
\ String with the counting 
time and total counts stored 
in
\ the buffered interface and 
send over to the computer
2 gpib.device gamma
\ The buffered interface is 
named gamma
72 string namefile
\ String consisting of the path 
and name of a data file.
\ This data file will contain 
the natural log values
\ of the corrected count rates 
through the empty column
\ at the specified locations, 
the calculated path lengths
\ and the dry bulk densities. 
The dry bulk densities will
\ be determined in calibration 
program 7.
real dim[ 15 ] array path.length.am
real dim[ 15 ] array path.length.cs
\ Arrays with the calculated path lengths at a 
maximum
\ of fifteen locations along 
the column.
real dim[ 15 ] array am.empty
real dim[ 15 ] array cs.empty
\ Arrays containing the average 
corrected count rates
\ through the empty column. 
(Am and Cs, respectively)
real dim[ 15 ] array am.filled
real dim[ 15 ] array cs.filled
146
\ Arrays containing the average corrected 
count rates
\ through the water filled column. 
(Am and Cs, resp.)
3 gpib.device motor
\ The Aston 800 is named motor.
9 gpib.device hp3497
\ The HP 3497 is named hp3497.
integer scalar number.days
\ The number of calender days the 
program is running.
integer scalar day
integer scalar days
\ Day of the month.
dp.integer scalar extra
integer scalar sec
integer scalar min
integer scalar hours
integer scalar month
dp.integer time0
72 string s2
\ The time string 'month:day:hour:min:sec' 
is read from
\ the hp3497.
dp.real scalar uwam
dp.real scalar uwcs
\ Product of mass attenuation coefficient 
and density of
\ deionized water. (Am and Cs, respectively)
integer scalar step.size
\ Step size between two locations. 
All locations are
\ evenly spaced.
6 string step.size2
\ Step between two locations expressed 
as a six letter
\ string
6 string step.size3
\ Step used in overshooting procedure
integer scalar number.locations
\ Total number of locations.
\ The Aston800 actuator is polled frequently 
when the
\ system moves from one location to 
the other. The moment
\ the sources and detector arrive 
at a new location, the
\ actuator should return the value 16 
to the computer. The
\ next statement of a program can 
then be executed. For no
\ apparent reason, the actuator 
does not return the value 16
\ after the system has travelled. 
To avoid a complete crash
\ of a program, a routine should be inserted that takes over
\ the task of the actuator by producing 
a 16 after a certain
\ time lapse. The colon definitions 
'day.determination' and
\ 'time.fetch' are part of this special 
routine. This
\ routine is part of the 'aston800.status.check' 
colon
\ definition.
\*****************************
147
: day. determination
1 number.days :=
0 extra :=
hp3497 listener me talker
" TD" "CR.LF" "CAT TALK
hp3497
s2 gpib.read
s2
3 1 do
ascii : "number ?drop
loop
"drop
day :=
month :=
untalk unlisten ;
: time. fetch
hp3497 listener me talker
" TD" "CR.LF" "CAT TALK
hp3497
s2 gpib.read
s2
5 1 do
ascii : "number ?drop
loop
32 "number
?drop "drop
sec :=
min :=
hours :=
days :=
month :=
days day =
if
else
number.days 1 + number.days :=
days day :=
number.days 1 - 86400 *
extra :=
then
hours 3600 * min 60 * + sec + extra +
untalk unlisten ;
\ **********************************************************
\ The colon definition 's' is 
used after the prompt
\ 'Press any key when ready'. The 
'p' macro clears the
\ screen and moves the cursor 
to the middle of the screen.
\ Macro 's' does the same after 
the user has pressed a key.
\ **********************************************************
148
: S
pckey drop ?drop screen.clear 10 1 
do cr loop ;
: p
screen.clear 10 1 do cr loop 
;
\ The colon definition 'interface' takes 
care of the GPIB
\ connections.
\ In this program, the three GPIB devices 
are the buffered
\ interface, defined as 'gamma' 
,the Aston 800, defined
\ as 'motor', and the HP3497, defined 
as 'hp3497'.
\ The buffered interface requires 
a 'carriage return'
\ (ASCII character 13) as an EOS 
(end of string) character.
\ The Aston 800 and the HP3497 require a "cr.lf",
\ catenated with the rest of the string 
that is send over
\ as an EOS. The HP3497 accepts capitals 
only.
: interface
bus. init
send. interface. clear
remote. enable. on
gamma eos.on
13 eos.character ;
\ Before the actual measurements 
can take place, the user has
\ to check a number of connections 
and settings. He also has
\ to enter coefficients of relations 
determined in previous
\ calibration programs.
: instructions
." Preset timer to '1,1'" cr cr
" Press any key when ready." s
.Check the following connections with the 
multiscaler:" cr
" NA-22 integral output to channel A,pos;" cr
.TC 450 norm mode output to channel B,pos;" 
cr
.NC-22 integral output to channel C,pos." 
cr cr
." Press any key when ready." s
." Check E and delta E settings of the " cr
." single channel analyzers and the NA-22." cr cr
." Press any key when ready." s
." Turn power to NH85 (next to buffered interface) off." 
cr
cr ." Press any key when ready." s
" Move sources and detector manually to the first 
location"
cr
149
." where you want to know the path length. This point" cr
." has to be the highest of all locations." cr cr
." Press any key when ready." s
." Write down the coordinates of this location." cr 
cr
." Press any key when ready." s
" Determine the step size between two locations." 
cr cr
." Press any key when ready." s
" The stepper motor indexers should be on external 
mode" cr
." before you can continue. " cr
." Both stepper motor indexers should read zero." 
cr cr
." Press any key when ready." s
." Set stop/count button of TC 535P on 'count' 
and" cr
." set .01 min/.l sec button on '.1 sec'." cr 
cr
." Press any key when ready." s
." Make sure the readout mode of the buffered inter-" cr
." face is on auto recycle and that the output 
" cr
" format is on normal." cr cr
." Press any key when ready." s
." The coefficients of the linear relation (y=ax+b)" cr
." between the integral Am count rate via the NC-22" 
cr
." (x) and via the TC 450 (y) must be entered." cr 
cr
-1 6 fix.format
." Enter a:" cr
#input a := cr
." Enter b:" cr
#input b := p
." The coefficients of the linear relation (y=rx+q)" cr
." between the high energy Cs count rate (x) 
and the" cr
." Compton scattering in the Am-window (y) must 
be" cr
." entered." cr cr
." Enter r:" cr
#input r := cr
." Enter q:" cr
#input q := p
." Enter product of mass attenuation coefficient and" cr
." density of deionized water, for Am:" cr
#input uwam := p
." Enter product of mass attenuation coefficient and" cr
." density of deionized water, for Cs:" cr
#input uwcs := p
." Enter the # of steps between two locations:" cr
#input step.size := p
." Enter the number of steps again, but now as a six digit"
cr
." number with leading zero's. For instance, 1000" cr
." steps has to be entered as 001000." cr
"input step.size2 ":= p
." Multiply the # of steps between two locations" cr
." by (number of locations - 1) and add 1000." cr
." Enter this number as a six digit number" cr
." with leading zero's" cr cr
"input step.size3 ":= p
150
." Enter the total number of 
locations:" cr
." (maximum of 15)" cr
#input number.locations := p
." Place the empty column at exactly 
the same position" cr
." as during the actual experiment." cr
cr ." Press any key when ready." 
s
." System must be able to move 
freely. Check this." cr cr
." Press any key when ready." s
." Turn power of NH85 on." cr
." Wait till the Aston 800 
actuator has performed" cr
." its self test. The letters 'Fn' 
appear when" cr
." when the self test has been completed." 
cr cr
" Press any key when ready." s
interface
untalk unlisten
gamma listener me talker
18 ASCII" 3 ASCII" "CAT TALK
." Preset timer (minimum 200 sec)." 
cr
" Timer formula: timebase*n*10**m." 
cr
." A counting time of 200 sec. is recommended." 
cr
." Set timer to '2,3' if you want 
a 200 sec." cr
" counting period." cr cr
." Press any key when ready." s
." Set CH A PRESET to '96'. " cr cr
." Press any key when ready." s
" Enter counting time in- seconds (200 
sec. minimum)." cr
#input 1000 * counting.time := p
stack. clear
\ **********************************************************
\ The following colon definition performs 
a status test of the
\ Aston 800. The serial poll check 
gives information on the
\ status of the actuator and the connected 
stepper motors.
\ When the actuator returns the value 
16, the system has
\ arrived at its destination. The next 
command of the
\ program can then be executed.
\ ********************************************************
: aston800.status.check
untalk unlisten
time.fetch time0 :=
begin
stack. clear
time. fetch
time0 -
60 >
if
16 16
else
1
151
then
16 =
if
16 =
else
serial.poll.enable
me listener motor talker
stack. listen
serial.poll.disable
16 =
then
until
untalk unlisten
\ If the program execution is delayed (by using the word
\ msec.delay) past midnight, the program will crash. The
\ following colon definition will prevent this.
: no.crashes
86000000. rel.time <
if
begin
rel.time 500. <
until
then ;
\ **********************************************************
\ A file template is made up for file maintenance. The data
\ file will contain six subfiles, each with a storage
\ capacity of a 15-element array.
\ The first subfile contains the natural log values of the
\ count rates through the empty column at the specified
\ locations for Am. The second subfile contains these 
data
\ for Cs.
\ The third subfile is used for path length storage,
\ determined with the Am source.
\ Subfile 4 contains the path lengths determined with 
Cs.
\ Subfile 5 will be used to store the dry bulk densities,
\ while subfile 6 will contain the volumetric water
\ contents, to be determined with calibration program 7.
\ The name of the file is assigned to the string 'namefile'.
\ The 'file.create' word creates a file according to the
\ latest 'file.template'.
\ *********************************************************
152
: specifications.file
file. template
real dim[ 15 ] subfile
6 times
end
." Enter name of data file (including path)." cr
"input namefile ":=
namefile defer> file.create
p
\ Five measurements are done at each location for the
\ empty and deionized water filled column, respectively.
: measurement
0. ycs.filled :=
0. ycs.empty :=
0. yam.filled :=
0. yam.empty :=
0 k := 0 m := 0 n :=
3 1 do
0 1 :=
im :=
m 2=
if
p ." Some action required after this tune!" cr
10000 1000 tune
p
." Fill column up with deionized water." cr
. It is better not to remove the column." cr
. If this is unavoidable, be sure the column is" cr
.replaced at exactly the same position." cr cr
. Press any key when ready." s
then
6 1 do
i k :=
number.locations 1 + 1 do
in :=
1 1 + 1 :=
n 1 =
if
no. crashes
motor listener me talker
" sm0,za" "cr.lf" "cat talk
" sm0,zp" "cr.lf" "cat talk
2500 msec.delay
else
no. crashes
153
motor listener me talker
" sm0,mr+" step.size2 "cat 
"cr.lf" "cat talk
aston800 .status .check
2500 msec.delay
then
." Loop " k . ." (out of 5)." 
cr
." Location " i . ." (out 
of " number.locations
.l )" cr
gamma listener me talker
27 ASCII" TALK
\ This command starts the scalers 
and the timer.
counting.time 2000 + msec.delay
gamma listener me talker
" P" 3 ASCII" "CAT TALK
\ This command transfers the 
data string with the
\ information on counts and 
time from the multiscaler
\ to the buffered interface.
4 1 do
gamma sl gpib.read
loop
\ The previous loop transfers 
the data string
\ from the buffer to the computer.
\ The next loop converts 
the data string into the counts 
and
\ counting time.
sl
5 1 do
32 "number ?drop
loop
ratenc :=
ratetc :=
ratena :=
t :=
"drop
\ The raw data are corrected 
for dead time. The Am count
\ rate is also corrected for 
low energy Cs in the Am window.
ratetc t / 10. *
ratenc t / 10. * a * b + 
.0000006 * 1. - / abs
rateam :=
ratena t / 10. * dup
.00000167 * 1. - abs /
dup ratecs :=
r * q + rateam swap - rateam
m 1 =
if
ratecs ln ycs.empty [ 1 ] :=
rateam ln yam.empty [ 1 ] :=
p
else
ratecs In ycs.filled [ 1 ] :=
rateam ln yam.filled [ 1 ] :=
p
then
154
:=
n number.locations =
if
motor listener me talker
" smO,mr-" step.size3 "cat "cr.lf" "cat talk
aston800. status. check
2500 msec.delay
motor listener me talker
" sm0,mr+001000" "cat "cr.lf" "cat talk
aston800, status, check
2500 msec.delay
then
loop
loop
loop
\ In the colon definition 'calculations', the average count
\ rates through the empty column and the path lengths
\ at the specified locations are calculated. The data are
\ stored in a data file created with the colon definition
\ 'specifications.file'.
: calculations
0. k := 0. m :=
0. am.empty :=
0. am.filled :=
0. cs.empty :=
0. cs.filled :=
0. path.length.am :=
0. path.length.cs :=
number.locations m :=
i 1 :=
m 1 + 1 do
yam.empty [ i ] i m + k := yam.empty [ k ] +
k m + k := yam.empty [ k ] +
k m + k := yam.empty [ k ] +
k m + k := yam.empty [ k ] + 5. /
am.empty [ i ] :=
ycs.empty [ i ] i m + k := ycs.empty [ k ] +
k m + k := ycs.empty [ k ] +
k m + k := ycs.empty [ k ] +
k m + k := ycs.empty [ k ] + 5. /
cs.empty [ i ] :=
yam.filled [ i ] i m + k := yam.filled [ k ] +
k m + k := yam.filled [ k ] +
k m + k := yam.filled [ k ] +
k m + k := yam.filled [ k ] + 5. /
am.filled [ i ] :=
155
ycs.filled [ i ] i m + k := ycs.filled 
[ k ] +
k m + k := ycs.filled 
[ k ] +
k m + k := ycs.filled [ 
k ] +
k m + k := ycs.filled [ 
k ] + 5. /
cs.filled [ i ] :=
loop
m 1 + 1 do
am.empty [ i ] am.filled 
[ i ] - uwam /
path.length.am [ i ] 
:=
cs.empty [ i ] cs.filled 
[ i ] - uwcs
path.length.cs [ i ] :=
loop
screen.clear
." Path lengths determined 
with the Am source: 
" cr
" (upper location first)" 
cr cr
m 1 + 1 do
path.length.am [ 
i ] . cr
loop cr cr
." Write this information 
down!" cr cr
." Press any key if 
you want to continue." 
s screen.clear
" Path lengths determined 
with the Cs source: 
" cr
." (upper location first)" 
cr cr
m 1 + 1 do
path.length.cs [ i ] 
. cr
loop cr cr
namefile defer> file.open
1 subfile
am.empty array>file
2 subfile
cs.empty array>file
3 subfile
path.length.am array>file
4 subfile
path.length.cs array>file
file.close
\ The main program.
\ After this program 
is loaded into the right 
ASYST
\ environment, entering 
'path.length.determination'
\ will run the program.
: path.length.determination 
p
interface
day.determination
instructions
specifications.file
measurement
156
calculations
157
\ +++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ PROGRAM BULK.TXT
\ Calibration program 7a
\ This program determines the 
dry bulk densities at the
\ same locations where the gamma 
measurements will be taken
\ during the experiments. The dry 
bulk density values
\ are determined on a water saturated 
soil.
\ An overshooting technique 
is used to overcome the
\ hysteresis of the stepper motors 
and to assure that the
\ locations where the measurements 
are taken are the
\ same for each loop.
\ Minimum ASYST overlay requirements 
on custom version:
\ gpib master
\ type 1 nec gpib drive
\ data files
\ At the C:\ prompt, enter c:\asyst\total 
to create
\ a correct ASYST environment. Then enter
\ load o:\asyst\bulk.txt. To run 
the program, enter
\ bulk.density.determination.
\ 09/26/89
\ M.Oostrom
\ ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
\ variable declaration
integer scalar 1
integer scalar m
integer scalar k
integer scalar n
\ 1,k,m and n are counters in loops
real scalar a
real scalar b
real scalar r
real scalar q
\ The coefficients of the linear relation
\ between the integral count rates via 
the NC-22 and
\ the TC 450 are a and b.
\ (calibration program 2)
\ The coefficients r and q determine 
the linear relation
\ between the high energy Cs count rate 
and the Cs
158
\ Compton scattering in the Am-window.
\ (calibration program 3)
real scalar ratenc
\ Raw integral Am count rate measured via the NC-22
real scalar ratetc
\ Raw differential Am count rate measured via the TC 450
real scalar ratena
\ Raw integral Cs count rate measured via the NA-22
dp.integer scalar counting.time
\ Counting time entered by the user
real scalar t
\ Counting time send over from buffer to computer
real scalar rateam
\ Dead time and compton scattering corrected Am count rate
real scalar ratecs
\ Dead time corrected Cs count rate
72 string sl
\ String with the counting time and total counts stored in
\ the buffered interface and send over to the computer
2 gpib.device gamma
\ The buffered interface is named gamma
72 string namefile
\ String consisting of the path and name of a data file.
\ This data file will contain the natural log values
\ of the corrected count rates through the empty column
\ at the specified locations, the calculated path lengths
\ and the dry bulk densities. The dry bulk densities are
\ determined in this program.
real dim[ 15 ] array bulk.density
\ The dry bulk density values at a maximum of 15 locations.
real dim[ 15 ] array bd
\ Average bulk density
real dim[ 15 ] array theta
\ The volumetric water contents at a maximum of 15
\ locations.
real dim[ 15 ] array th
\ Average theta
real dim[ 15 ] array path.length.am
real dim[ 15 ] array path.length.cs
\ Arrays with the calculated path lengths at a maximum
\ of 15 locations along the column. The path lengths are
\ determined in program 6 and stored in a data file.
real dim[ 15 ] array am.empty
real dim[ 15 ] array cs.empty
\ Arrays containing the average corrected count rates
\ through the empty column. (Am and Cs, respectively)
\ These data are also determined in program 6.
real dim[ 15 ] array am.filled
real dim[ 15 ] array cs.filled
\ Arrays containing the average corrected count rates
\ through the water saturated soil column. (Am and Cs, resp.)
3 gpib.device motor
159
\ The Aston 800 is named motor
dp.real scalar uwam
dp.real scalar uwcs
\ Volumetric attenuation coefficient of
\ deionized water. (Am and Cs, respectively)
dp.real scalar usam
dp.real scalar uscs
\ Soil mass attenuation coefficients. (Am and Cs,
\ respectively)
integer scalar step.size
\ Step size between two locations. All locations are
\ evenly spaced.
6 string step.size2
\ Number of steps between locations expressed as a
\ six letter string
6 string step.size3
\ String used in overshooting process
integer scalar number.locations
\ Total number of locations.
9 gpib.device hp3497
\ The HP 3497 is named hp3497.
integer scalar number.days
\ The number of calender days the program is running.
integer scalar day
integer scalar days
\ Day of the month.
dp.integer scalar extra
integer scalar sec
integer scalar min
integer scalar hours
integer scalar month
dp.integer scalar time0
72 string s2
\ The time string 'month:day:hour:min:sec' is read from
\ the hp3497.
\ The colon definition 's' is used after the prompt
\ 'Press any key when ready'. The 'p' macro clears the
\ screen and moves the cursor to the middle of the screen.
\ Macro 's' does the same after the user has pressed a key.
: s
pckey drop ?drop screen.clear 10 1 do cr loop ;
:p
screen.clear 10 1 do cr loop ;
\ ******************************* 
* * * ** ****
160
\ The colon definition 'interface' takes care 
of the GPIB
\ connections.
\ In this program the three GPIB 
devices are the buffered
\ interface, defined as 'gamma' 
,the Aston 800, defined
\ as 'motor', and the HP3497, 
defined as 'hp3497'.
\ The buffered interface requires 
a 'carriage return'
\ (ASCII character 13) as an 
EOS (end of string) character.
\ The Aston 800 and the HP3497 require 
a "cr.lf",
\ catenated with the rest of the 
string that is send over
\ as an EOS. The HP3497 accepts 
capitals only.
: interface
bus. init
send. interface, clear
remote.enable.on
gamma eos.on
13 eos.character ;
\ The Aston800 actuator is polled 
frequently when the
\ system moves from one location 
to the other. The moment
\ the sources and detector 
arrive at a new location, the
\ actuator should return the 
value 16 to the computer. The
\ next statement of a program 
can then be executed. For no
\ apparent reason, the actuator 
does not return the value 16
\ after the system has travelled. 
To avoid a complete crash
\ of a program, a routine should 
be inserted that takes over
\ the task of the actuator 
by producing a 16 after a certain
\ time lapse. The colon definitions 
'day.determination' and
\ 'time.fetch' are part of this special 
routine. This
\ routine is part of the 'aston800.status.check' 
colon
\ definition.
\ ********************************************************
: day.determination
1 number.days :=
0 extra :=
hp3497 listener me talker
" TD" "CR.LF" "CAT TALK
hp3497
s2 gpib.read
s2
3 1 do
ascii : "number ?drop
loop
"drop
day :=
161
month :=
untalk unlisten
: time. fetch
hp3497 listener me talker
" TD" "CR.LF" "CAT TALK
hp3497
s2 gpib.read
s2
5 1 do
ascii : "number ?drop
loop
32 "number
?drop "drop
sec :=
min :
hours :=
days :=
month :=
days day =
if
else
number.days 1 + number.days :=
days day :=
number.days 1 - 86400 *
extra :=
then
hours 3600 * min 60 * + sec + extra +
untalk unlisten
\ * * * * * * * * * * * * * * 
* * * * * * * * * * * * * * * * 
* * * * * * * * * * * * * * * * * 
* * * * ** * * * * *
\ Before the actual measurements can take place, the user has
\ to check a number of connections and settings. He also has
\ to enter coefficients of relations determined in previous
\ calibration programs.
\ *******************************************www* * *** *
: instructions
." Turn power to NH85 (next to buffered interface) off." cr
cr ." Press any key when ready." s
." Preset timer to '1,1'" cr cr
." Press any key when ready." s
." Check the following connections with the multiscaler:" cr
." NA-22 integral output to channel A,pos;" cr
." TC 450 norm mode output to channel B,pos;" cr
." NC-22 integral output to channel C,pos." cr cr
." Press any key when ready." s
." Check E and delta E settings of the " cr
." single channel analyzers and the NA-22." cr cr
." Press any key when ready." s
162
." Set stop/count button of TC 535P on 'count' and" cr
" set .01 min/.l sec button on '.1 sec'." cr cr
." Press any key when ready." s
." Make sure the readout mode of the buffered inter-" cr
" face is on auto recycle and that the output " cr
." format is on normal." cr cr
" Press any key when ready." s
" The coefficients of the linear relation (y=ax+b)" cr
." between the integral Am count rate via the NC-22" cr
." (x) and via the TC 450 (y) must be entered." cr cr
-1 6 fix.format
." Enter a:" cr
#input a := cr
." Enter b:" cr
#input b := p
." The coefficients of the linear relation (y=rx+q)" cr
." between the high energy Cs count rate (x) and the" cr
" Compton scattering in the Am-window (y) must be" cr
." entered." cr cr
" Enter r:" cr
#input r := cr
." Enter q:" cr
#input q := p
." Enter volumetric attenuation coefficient" cr
." of deionized water, for Am:" cr
#input uwam := p
." Enter volumetric attenuation coefficient" cr
." of deionized water, for Cs:" cr
#input uwcs := p
." Enter Am soil mass attenuation coefficient:" cr
#input usam := p
." Enter Cs soil mass attenuation coefficient:" cr
#input uscs := p
." Enter the step size between two locations:" cr
#input step.size := p
." Enter the step size again, but now as a six digit" cr
." number with leading zero's. For instance, a step " cr
." size of 1000 has to be entered as 001000." cr
"input step.size2 ":= p
." Multiply step size between locations by" cr
." (number of locations - 1) and add 1000." cr
." Enter this number as a six letter digit with." cr
" leading zero's." cr cr
"input step.size3 ":= p
." Enter the total number of locations:" cr
." (maximum of 15)" cr
#input number.locations := p
." Place column exactly the way it was during" cr
." the executing of program 6." cr cr
." Press any key when ready." s
." System must be able to move freely. Check this." cr cr
." Press any key when ready." s
163
." Move stepper motors manually to the location where" 
cr
" you want to do the first measurement." cr
." The coordinates of the stepper motors for 
this" cr
." position were written down during program 6." 
cr cr
" Press any key when ready." s
." Both stepper motor indexers should be 
on external mode" cr
." before you can continue. All the counters on" 
cr
" both stepper motor indexers should read zero." 
cr cr
." Press any key when ready." s
." Turn power of NH85 on." cr
." Wait till 'Fn' appears on the display 
of the actuator." cr
cr
." Press any key when ready." s
interface
untalk unlisten
gamma listener me talker
18 ASCII" 3 ASCII" "CAT TALK
" Preset timer (minimum 200 sec)." cr
." Timer formula: timebase*n*10**m." cr
" If you want to use a counting time of 200 sec.," cr
." set timer to '2,3'." cr cr
." Press any key when ready." s
." Set CH A PRESET to '96'. " cr cr
." Press any key when ready." s
" Enter counting time in seconds (200 sec. minimum)." 
cr
#input 1000 * counting.time := p
stack. clear
\ The following colon definition performs a status 
test of the
\ Aston 800. The serial poll check gives information 
on the
\ status of the actuator and the connected stepper 
motors.
\ When the actuator returns the value 16, the system has
\ arrived at its destination. The next command of 
the
\ program can then be executed.
: aston800.status.check
untalk unlisten
time.fetch time0 :=
begin
stack.clear
time. fetch
time0 -
60 >
if
16 16
else
1
then
164
16 =
if
16=
else
serial.poll.enable
me listener motor talker
stack. listen
serial.poll.disable
16 =
then
until
untalk unlisten
\ If the program execution is delayed (by using 
the word
\ msec.delay) past midnight the program will 
crash. The
\ following colon definition will prevent 
this.
: no.crashes
86000000. rel.time <
if
begin
rel.time 500. <
until
then ;
\ Data collected with program x.txt are retrieved
\ from the data file. This data file consists of six
\ subfiles. The contents of each subfile are assigned 
to an
\ array.
\ The first subfile contains natural log values 
of the count
\ rates through the empty column at the specified 
locations
\ for Am. The second subfile contains these 
data for Cs.
\ The third subfile is used for path length storage,
\ determined with the Am source. Subfile 4 contains 
the
\ path lengths determined with Cs. Subfile 5 will be 
used
\ to store the dry bulk density values and subfile 6 
the
\ volumetric water content values, to be determined 
with
\ this program.
\ The name of the file is assigned to the string 'namefile'.
\ **********************************************************
: retrieve.data
." Enter name of data file (including path)." cr
"input namefile ":=
165
namefile defer> file.open
1 subfile
am.empty file>array
2 subfile
cs.empty file>array
3 subfile
path.length. am file>array
4 subfile
path.length.cs file>array
file.close
p
\ Count measurements are taken at all the locations. This
\ loop can be repeated 50 times.
\ The data are stored temporarely in the arrays 'yam.filled'
\ and 'ycs.filled'.
: measurement
0 bd :=
0 th :=
0 k := 0 1 := 0 n :=
51 1 do
i k :=
number.locations 1 + 1 do
in:=
n =
if
no. crashes
motor listener me talker
" sm0,za" "cr.lf" "cat talk
" sm0,zp" "cr.lf" "cat talk
2500 msec.delay
else
no.crashes
motor listener me talker
" sm0,mr+" step.size2 "cat "cr.lf" "cat talk
aston800.status.check
2500 msec.delay
then
." Location " i . ." (out of " number.locations .
." )" cr
gamma listener me talker
27 ASCII" TALK
\ This command starts the scalers and the timer.
counting.time 2000 + msec.delay
gamma listener me talker
" P" 3 ASCII" "CAT TALK
166
\ This command transfers 
the data string with 
the information
\ on counts and time 
from the multiscalar 
to the buffered
\ interface.
4 1 do
gamma sl gpib.read
loop
\ The previous loop transfers 
the data string
\ from the buffer to the 
computer.
\ The next loop 
converts the data 
string into the counts 
and
\ counting time.
sl
5 1 do
32 "number ?drop
loop
ratenc :=
ratetc :=
ratena :=
t :=
"drop
\ The raw data are corrected 
for dead time. The Am 
count
\ rate is also corrected 
for low energy Cs in 
the Am window.
ratetc t / 10. *
ratenc t / 10. * a * b + 
.0000006 * 1. - / abs
rateam :=
ratena t / 10. * dup
.00000167 * i. - abs /
dup ratecs :=
r * q + rateam swap - rateam 
:=
ratecs In cs.filled [ n 
] 
:=
rateam In am.filled [ n ] 
:=
p
n number.locations =
if
motor listener me talker
" sm0,mr-" step.size3 
"cat "cr.lf" "cat talk
100000 msec.delay
motor listener me talker
" sm0,mr+001000" "cr.lf" 
"cat talk
aston800.status. check
2500 msec.delay
number.locations 1 + 1 do
i 1 :=
uwcs am.empty [ 1 ] am.filled 
[ 1 ] - *
uwam cs.empty [ 
1 ] cs.filled [ 
1 ] - * -
usam uwcs * uscs uwam * - path.length.am [ 1 ] * 
/
bulk.density [ 1 ] :=
bulk.density [ 1 ] bd [ 1 
] + bd [ 1 ] :=
usam cs.empty [ 1 ] cs.filled 
[ 1 ] - *
uscs am.empty 
[ 1 ] am.filled 
[ 1 ] - * -
usam uwcs * uscs uwam * -
path.length.am [ 1 ] * /
theta i 1 
] 
:=
theta [ 1 3 th [ 1 ] + th 
[ 1 ] :=
167
loop
." Loop" k . cr cr
." Bulk density values: " cr
." (upper location first)" cr cr
number.locations 1 + 1 do
bulk.density [ i ] . cr
loop cr cr
screen. print
." Press any key when ready." s screen.clear
." Volumetric water content values: " cr
." (upper location first)" cr cr
number.locations 1 + 1 do
theta [ i ] . cr
loop cr cr
screen. print
." Enter 'y' if you want to do another run" cr
." Enter 'n' if you want to quit" cr
"input " y" "= p
if
else
bd k / bd :=
th k / th :=
." Average bulk density values: " cr
." (upper location first)" cr cr
number.locations 1 + 1 do
bd [ i ] . cr
loop cr cr
screen.print
." Press any key when ready." s screen.clear
." Average volumetric water content values: " cr
." (upper location first)" cr cr
number.locations 1 + 1 do
th [ i ] . cr
loop cr cr
screen. print
." Press any key when ready." s screen.clear
quit
then
then
loop
loop
namefile defer> file.open
5 subfile
bd array>file
6 subfile
th array>file
file.close
screen. clear
168
\ The main program.
\ After this program is loaded into the right Asyst
\ environment, entering 'bulk.density .determination'
\ will run the program.
bulk.density.determination p
interface
day.determination
instructions
retrieve.data
measurement
169