Technical Note

Calibration and Data Processing for Airborne Gamma-Ray Spectrometry

The multichannel gamma-ray spectra recorded at airborne altitudes are influenced by many variables. The role of data processing is to correct the observed data for those influences not related to the geology, and to transform the observed spectra to equivalent concentrations of the radioelements in the ground. The purpose of this Technical Note is to describe the methodology used by GammaSpec to calibrate airborne spectrometers and process airborne gamma-ray data. 

Each of the corrections for airborne gamma-ray spectrometric data is described below in the same sequence in which they must be applied. The processing sequence is as follows:

1. Noise reduction on multichannel spectra using the NASVD method.
2. Energy calibration of the gamma-ray spectra
3. Live time correction and windowing
4. Aircraft and cosmic background removal
5. Radon background estimation and removal
6. Stripping (channel interaction correction)
7. Calculation of effective height
8. Altitude correction
9. Sensitivity correction (reduction to elemental concentrations)

The calibration requirements for each correction are described at the end of each section. Table 2, at the end of this Technical Note, summarises the calibration requirements for the gamma-ray corrections, as well as the GammaSpec tools used to estimate the calibration coefficients.

Noise Reduction Using the NASVD Method

Noise-adjusted singular value decomposition (NASVD) is a spectral component analysis procedure for the removal of noise from gamma-ray spectra. The procedure transforms observed spectra into orthogonal spectral components. The lower-order components represent the signal in the original observed spectra, and the higher-order components represent uncorrelated noise. The component spectra that represent the noise are rejected, and the original spectra are "smoothed" by reconstructing noise-reduced spectra from the lower-order components only. 

The method is similar to principal component analysis. The principal components of a set of m spectra A (mx256 where m>=256) are the eigenvectors of the covariance matrix, ATA. The principal components are mutually orthogonal and are sorted, by eigenvalue, into descending order. The eigenvalues are the variances associated with each principal component. The observed spectra can then be represented as a linear combination of the principal components, i.e.

A = CS,

where C (mx256) is a matrix of concentrations (or amplitudes) and S (256x256) are the principal components (spectral shapes). Since the lower-order components represent the signal in the original observed spectra, and the higher-order components represent uncorrelated noise, the noise can be removed by reconstructing "smoothed" spectra from lower-order components only. 

There are, however, two issues with the use of the standard principal component method for the analysis of observed multichannel spectra. First, the variance associated with each channel count rate for each spectrum must be the same. Hovgaard (1997) suggested a simple solution to this problem. Since counting errors in gamma-ray spectrometry are Poisson distributed, the variance of a particular channel count rate is the same as the mean count rate for that channel. Also, since changes in spectral shape in gamma-ray spectra are typically small, the best fit of the mean spectrum to each input spectrum provides a good estimate of the mean count rate (and hence variance) for each channel. Thus, by scaling each observed spectrum by the best fit of the mean spectrum to each observed spectrum, we normalise the variance in each channel of each spectrum to unity. The second issue is that the principal component method requires the data to be mean-centred. This centring problem is overcome by using the singular value decomposition (SVD) method to analyse the dispersion of the data around the origin rather than the sample mean. This produces eigenvectors which are used in the same way as the principal components in a principal component analysis.

Figure 1 below shows the first 16 (lower-order) eigenvectors from a typical flight of airborne spectra. The signal is concentrated in the lower-order components. 

Figure 1. Typical eigenvectors derived from the analysis of a flight of airborne gamma-ray spectra. The figure shows the lower-order 16 eigenvectors ordered left-to-right and top-to-bottom. 

In summary, the NASVD method is applied as follows:

• Spectra are grouped into groups of around 20,000 to 50,000 spectra for NASVD processing. Sometimes even larger groups are used. Alternatively, GammaSpec facilitates the sorting of spectra into clusters (see Minty and McFadden, 1998) based on spectral shape for the purpose of NASVD analysis.
• For each group of spectra, the average spectrum is used to scale the spectra to approximately unit variance in each channel.
• An SVD analysis is applied to the scaled spectra to obtain the eigenvectors and their associated amplitudes.
• The eigenvectors, or amplitudes, or both, can be scaled in such a way as to effectively “undo” the scaling of the input spectra.
• GammaSpec saves the eigenvectors and amplitudes for later use. The lower-order eigenvectors and amplitudes (typically 8) are used to reconstruct smooth (noise-reduced) spectra for further processing.

  

Energy Calibration of Spectra

All spectrometers are affected by energy drift in the measured spectra. It is not unusual for older airborne spectrometers to drift up to 2-3 channels during the course of a day (1 channel = 11.7 keV), and spectrometer drift is thus a significant source of error with these instruments. Modern spectrometers have in-built self-stabilising features that can automatically adjust the gain of the instrument to maintain a prominent photo peak (usually Th at 2.61 MeV) at its nominal channel position to within a fraction of a channel. However, even self-stabilising spectrometers suffer a slow energy drift over the course of many flights due to drift in the base level of the spectra (the energy at the bottom of channel 1). Spectrometer data can be corrected for energy drift as long as all of the crystal/photo-multiplier assemblages comprising the detector drift sympathetically.

Spectra are first grouped for the purpose of energy calibration into groups of at least 1000 consecutive spectra. The sum spectrum for each group is used to calculate a base level (energy at the bottom of channel 1) and gain (channel width) for that group, which is assigned to each spectrum in the group. Note that the raw gamma-ray spectra are never actually adjusted for drift. Rather, GammaSpec saves the base and gain for each spectrum as new database fields. The base level and gain are used to calculate channel positions for subsequent windowing (summing channels over designated energy windows) of spectra during processing. 

GammaSpec offers two methods for estimating the base level and gain of airborne spectra:

• The energies of at least two prominent photopeaks in the sum spectrum are determined as the maximum value of a quadratic fitted in the vicinity of each photopeak. A linear function is then fit to the photopeak positions (channel numbers vs energy) to obtain the energy at channel 1 and the gain (keV per channel) – see Figure 2. 
• An automatic method whereby the sum spectrum is iteratively shifted and scaled to obtain a best-fit with one of a number of standard spectra. This yields the base level and shift of the sum spectrum.

Figure 2. Energy calibration using the positions of prominent photo peaks in the spectrum.

 

Live-Time Correction and Widowing

Spectrometers require a finite time to process each pulse from the detector. Any pulse that arrives while another is being processed is automatically rejected, and the total counting time available is thus reduced by the time taken to process all pulses detected (the “dead time”). The “live time” may thus be smaller than the sample integration time; although for modern spectrometers the difference is usually small. A typical live time for a 1-s sample interval would be about 0.990 s. 

GammaSpec works with count rates (counts/sec), and the live time correction thus also corrects the raw data for sample intervals other than the typical 1 s. The correction is a simple scaling of the observed counts in a spectrum channel (or window) to convert them to counts per unit time.

                       
N = n/t

where N is the corrected count (counts/sec), n is the observed count and t the live time in seconds.

After NASVD and energy calibration, the spectra are summed over the conventional TC, K, U and Th windows recommended by the IAEA (see Table 3 in Tech Note 1). The base level and gain for each spectrum are used to calculate the channel positions of the windows for each spectrum. Adjacent channels are summed to produce the window counts. These are adjusted for live time and saved as new database fields. Where the window boundaries occur in the middle of a spectrum channel (which is usually the case), the channel counts are scaled downward in proportion to the fraction of the channel that appears within the window. 

Aircraft and Cosmic Background Correction

The "aircraft" spectrum is due to the radioactivity of the aircraft and its contents and is a constant (Figure 3). The "cosmic" spectrum (Figure 4) arises from the interaction of primary cosmic radiation with atoms and molecules in the upper atmosphere. The procedure for removing aircraft and cosmic background uses the fact that in the lower atmosphere the shape of the cosmic spectrum is constant, with its amplitude decreasing with decreasing altitude. Also, at energies above 3.0 MeV all radiation is cosmic in origin. All airborne spectrometers routinely monitor a "cosmic" window (typically 3.0-6.0 MeV) for estimating cosmic background. So, if we know the shape of the cosmic spectrum, we can scale it using the cosmic window count rate to calculate the cosmic contribution in any particular spectrum channel or window.

The aircraft and cosmic background spectra are estimated as follows:

Ni = ai + biNcos

where        Ni  =  aircraft plus cosmic background count rate in the i'th channel;

Ncos= cosmic window count rate;

ai   =  aircraft background in the i’th channel;

bi   =  cosmic background in the i’th channel normalised to unit counts in the cosmic window.

The aircraft and cosmic spectra are determined through suitable calibration. The cosmic window count rate can be lightly filtered prior to estimating the correction. GammaSpec facilitates both a full-spectrum method for aircraft and cosmic background correction, and a window method. The aircraft and cosmic background spectra (or window count rates) are subtracted from the livetime-corrected and energy-calibrated observed spectra (or window count rates).

Figure 3.  A typical “aircraft” spectrum.

Figure 4.  A typical “cosmic” spectrum

Atmospheric Radon Background Correction

GammaSpec offers two methods for removing radon background. The first is a spectral-ratio method (Minty 1992, Minty 1998) where the relative heights of uranium series photo peaks are used to determine the contributions to the spectrum of uranium in the ground and radon in the air. The second approach is through the use of upward-looking detectors (IAEA 1991). Minty (1998) also describes a full-spectrum method for estimating radon background, but this is not supported by GammaSpec.

 

Radon Removal (Spectral-Ratio Method)

The method described here is the 2-component method described by Minty (1998). The method derives from the observation that the low energy 214Bi photopeak at 0.609 MeV from atmospheric radon suffers far less attenuation relative to the 214Bi peak at 1.76 MeV than is the case for radiation from the ground. So the ratio of the counts in each of these photopeaks is diagnostic of the relative contributions of atmospheric radon and uranium in the ground to the observed spectrum.

Minty (1992) also described a 4-component model that considers the contributions of atmospheric radon and terrestrial sources of potassium, uranium and thorium to the spectral regions of interest. But the calibration requirements for this model are stringent. Fortunately, the simple 2-component model based on an atmospheric radon spectrum and a composite potassium, uranium and thorium ‘ground’ spectrum can be used to estimate radon background almost as well as the complete 4-component model. This 2-component model works well in areas of relatively constant Th/U concentrations. It is also easily calibrated using normal survey data.

The 2-component method calls for the monitoring of two windows which we denote by L (low energy) and H (high energy) – see Figure 5. L represents counts above the Compton continuum in the 0.609 MeV photopeak, and H represents counts in the conventional uranium window (1.66 MeV–1.86 MeV).

Figure 5. Two-component spectral-ratio method: the radon and ground spectra relative to the two windows used by the method. The ‘observed’ spectrum has been corrected for aircraft and cosmic background.

Let Lob and Hob be the observed count rates in the low and high energy windows respectively, after correcting the spectrum for aircraft and cosmic background. Let Lr and Hr be the radon contributions to the low and high energy windows, and let Lg and Hg be the corresponding contribution to these windows due to radiation from the ground. Then

                       Lob = Lr + Lg                (1)

and 

                       Hob = Hr + Hg                (2)

Also, since the shape of the radon spectrum can be assumed constant for a particular altitude, and uranium and thorium concentrations correlate well over most lithological units, we have

                       Lr = c1Hr                (3)

and

                       Lg = c2Hg                (4)

where c1 and c2 are calibration coefficients that have to be determined. These quantities are shown graphically in Figure 5. Solving equations 1 to 4 gives

                                       
Lr = (Lob – c2Hob) / (1 – c2/c1)                (5)

The radon contribution in the conventional uranium window is then calculated from equation 5. The “radon stripping coefficients” relate the counts due to radon in the U window to the counts due to radon in each of the remaining windows as follows.

                       TCradon = radonStripTC x Uradon

                       Kradon   = radonStripK x Uradon

                       Thradon = radonStripTh x Uradon

The radon contributions to each window are subtracted from the window count rates. 

Radon Removal (Upward-looking Detector Method)

The upward-looking detector method uses an additional crystal pack which is partially shielded from radiation from below to give the system a directional sensitivity and the ability to discriminate between radiation from the atmosphere and from the ground. As with the spectral-ratio method, the radon is estimated by summing the window count rates (or spectra) over large time intervals – typically 200 s.

Grasty (1975) showed that the TC, K and Th backgrounds are linearly related to the uranium background, and the background in these channels can be derived from the background in the uranium channel by suitable calibration. After removal of aircraft and cosmic backgrounds from the upward U window and the downward TC, K, U and Th windows (IAEA, 1991): 

ur = auUr + bu

Kr = akUr + bk

Tr = atUr + bt

TCr = atcUr + btc

where ur is the radon contribution to the upward U window, and Kr, Tr and TCr are the radon contribution to the downward K, Th and TC windows. The ai and bi are coefficients that must be determined through calibration.

The measured count rates in the upward uranium window are also related to those in the downward uranium window for radiation due to uranium in the ground as follows:

ug = a1Ug + a2Tg

where ug, Ug and Tg are the ground components and a1 and a2 are coefficients that have to be determined through calibration. 

The radon contribution to the uranium window of the main detector package (i.e. the "downward" U window) is then given by 

       
Ur = (u – a1U – a2T + a2bt – bu) / (au – a1 – a2at)

       

where        Ur = radon background in the "downward" U window

       u   = count rate in the "upward" U window

       U  = count rate in the "downward" U window

       T  = count rate in the "downward" Th window

Stripping (Channel Interaction Correction)

This correction is used to correct the K, U and Th window count rates for those gamma-rays not originating from the radioelement or decay series being monitored by that particular window. For example, thorium series gamma-rays appear in both the uranium and potassium windows, and uranium series gamma-rays appear in the potassium window. The corrections are applied as follows:

NTh(corr) = NTh

NU(corr) = NU – αNTh(corr)

NK(corr) = NK – βNTh(corr) – γNU(corr)

where α are the counts in the U window per unit count in the Th window, β are the counts in the K window per unit count in the Th window, and γ are the counts in the K window per unit count in the U window. These quantities are shown graphically in Figure 6.

Figure 6. Observed, potassium, uranium and thorium spectra showing the stripping correction parameters.

The stripping ratio ‘a’ (uranium into thorium) is small, and is often neglected during processing. If this stripping ratio is included, then equation 15 becomes

NTh(corr) = (NTh – αNU) / (1 – aα)

NU(corr) = (NU – αNTh) / (1 – aα)

NK(corr) = NK – βNTh(corr) – γNU(corr)

The stripping ratios α, β and γ and ‘a’ can be calculated directly from pure spectra due to U and Th (Figure 6). The calculated stripping ratios are for data acquired at ground level. They are adjusted for the height of the detector (STP height) using correction factors (Table 1) published by the IAEA (1991). The altitude correction for stripping ratio ‘a’ is small and is neglected.

Table 1. Increase in stripping ratio with altitude (IAEA(1991).

Calculation of Effective Height

The measured radar (or laser) altimeter height cannot be used directly for correcting the window count rate deviations from the nominal survey height. The height must first be corrected for the ambient temperature and pressure, since both affect the density and thus the attenuating properties of the air. The equivalent height at standard temperature (273.15 °K) and pressure (101.325 kPa) – i.e. the STP height, is given by (IAEA, 1991)

hSTP = (273.15 x P x hobs) / (T + 273.15) x (1013.25)

where         hobs   = observed height above ground level (m);

       hSTP  = equivalent height at STP (m);

       T     = air temperature (°C);

       P     = barometric pressure (mbars).

The temperature, pressure and radio altimeter data can be lightly filtered prior to the calculation of the STP height.

Altitude Correction

The TC, K, U and Th window count rates can now be corrected to the nominal survey height to remove the effect of variations in the terrain clearance of the detector. The fall-off of radiation with height is approximately exponential for broad (effectively infinite) sources.

The altitude correction corrects the window counts for the difference in height between the measurement height and the nominal survey height as follows:

N = N0e-µ(H-h)

where        μ  = the empirically derived height attenuation coefficient for a broad source type;

       N0 = the observed count rate at the STP height, h, and

       N  = the corrected count rate for the nominal survey terrain clearance H.

Sensitivity Correction (Reduction to Elemental Count Rates)

The final window count rates are converted to concentrations by scaling the count rates by a window sensitivity factor as follows:

       C = N/S

where        N = background−corrected and stripped count rate at the nominal survey altitude;

       C  = ground concentration of the element;

       S  = sensitivity coefficient.

Sensitivity coefficients for each window are calculated from data acquired over a calibration range

Calibration Summary

Table 2. Calibration Tool Summary

References

Grasty, R.L., Holman, P.B., and Blanchard, Y.B., 1991. Transportable calibration pads for ground and airborne gamma-ray spectrometers. Geological Survey of Canada, Paper No. 90-23, 25p.

Grasty, R.L., 1979. Gamma-ray spectrometric methods in uranium exploration – theory and operational procedures. Paper 10B in Geophysics and geochemistry in the search for metallic ores. Geological Survey of Canada, Economic Geology Report 31.

IAEA, 1991.  Airborne gamma ray spectrometer surveying.  Technical Report Series, No. 323.  International Atomic Energy Agency, Vienna, 1991.

Minty, B., 1997.  Fundamentals of airborne gamma-ray spectrometry. AGSO Journal of Australian Geology and Geophysics, 17 (2), 39-50.

Recommended Reading

Grasty, R.L., and Minty, B.R.S., 1995.  The standardisation of airborne gamma-ray surveys in Australia. Exploration Geophysics, 26, 276-283.

Hovgaard, J., and Grasty, R.L., 1997.  Reducing statistical noise in airborne gamma-ray data through spectral component analysis.  In “Proceedings of Exploration 97: Fourth Decennial Conference on Mineral Exploration” edited by A.G.  Gubins, 1997, 753-764.

IAEA, 1991.  Airborne gamma ray spectrometer surveying.  Technical Report Series, No. 323.  International Atomic Energy Agency, Vienna, 1991.

IAEA, 2003.  Guidelines for radioelement mapping using gamma-ray spectrometry data. IAEA-TECDOC-1363, International Atomic Energy Agency, Vienna.

Minty, B.R.S. 1992.  Airborne gamma-ray spectrometric background estimation using full spectrum analysis. Geophysics, 57(2), 279-287. 

Minty, B.R.S. 1998.  Multichannel models for the estimation of radon background in airborne gamma-ray spectrometry. Geophysics, 63(6), 1986-1996.

Minty, B.R.S., and Hovgaard, J., 2002.  Reducing noise in gamma-ray spectrometry using spectral component analysis. Exploration Geophysics, 33, 172-176.

Minty, B., Luyendyk, A. and Brodie, R., 1997.  Calibration and data processing for airborne gamma-ray spectrometry. AGSO Journal of Australian Geology and Geophysics, 17 (2), 51-62.