Exploration Geophysics - Volume 34, Issue 1-2, 2003
Volume 34, Issue 1-2, 2003
- Articles
-
-
-
3D fracture assessment using AVAz and a layer-stripping approach
More LessAuthors Mu Luo and Brain J. EvansSubtle variations in P-wave multi-azimuthal seismic reflections can be used to map lineations in a horizon, representing seismic amplitude changes resulting from the presence of vertical fractures. This paper studies the phenomenon of P-wave propagation through vertically fractured media under controlled laboratory conditions, and shows how resulting amplitude maps can be used to predict fractured layers. A scheme of layer-stripping, based on Ruger’s reflection coefficient approximations, can reduce the overburden effect on reflection amplitudes of underlying layers. By removal of overburden transmission effects, greater clarity is achieved in mapping any deeper fractured horizons. The effect of fractures is examined using physical model data, and the technique is demonstrated using field data.
-
-
-
-
Anisotropic semblance analysis and NMO corrections for long-offset data
More LessAuthors Fanmin Zhang, Norman F. Uren and Milovan UrosevicSedimentary rocks can be regarded as transversely isotropic (TI) media. P-wave reflections from horizontal reflectors in transversely isotropic media have non-hyperbolic moveout. One difficulty in seismic processing is how to flatten such reflection events for long-offset data. For a multi-layered model, the reflection moveout formula is usually expressed as a Taylor series with higher order terms ignored. Alkhalifah and Tsvankin (1995) developed a three-term Taylor-series formula to calculate reflection travel times from a horizontal reflector in TI media with a vertical symmetry axis (VTI). Using this formula, NMO correction works well for short spread lengths, but not so well for long spreads.
Zhang and Uren (2001a, b) developed an approximate explicit analytical P-wave ray velocity function for transversely isotropic (TI) media. From this ray velocity function, a reflection travel-time formula for a horizontal reflector in TI media was derived. This formula can be used for anisotropic NMO correction. It works well for both large offsets and small offsets. In order to obtain the unknown parameters required for seismic processing, a 3D-semblance analysis technique has been developed. We tested this method on numerical data from TI models with a single horizontal reflector and from multiple horizontal reflector models. The method was also tested on an isotropic model with multiple horizontal reflectors. The results show that the events can be completely flattened even for very large offsets, and that both multi-layered TI and isotropic models may appear to be anisotropic.
-
-
-
Application of the radial basis function neural network to the prediction of log properties from seismic attributes
More LessAuthors Brain H. Russell, Laurence R. Lines and Daniel P. HampsonIn this paper, the radial basis function neural network (RBFN) is used to predict reservoir log properties from seismic attributes, and is compared to the generalized regression neural network (GRNN), discussed by Hampson et al. (2001). Both of these methods are related to the probabilistic neural network (PNN), which uses a Gaussian kernel estimator based on the distance between points in seismic attribute space. Our goal was to see if there are situations in which the RBFN could improve on the results found using the GRNN. The methodology consisted of first training the neural network at each well location, to find the optimum set of seismic attributes and weighting coefficients, and then applying the trained network to a 3D seismic volume. Each neural network was applied to the Blackfoot 3D seismic volume, which was recorded over a Cretaceous channel sand in central Alberta. Our results showed that, although the training results were quite close for each method, the RBFN approach generally produced higher-frequency results, especially as the number of training values was reduced. By computing the least-squared error between the predicted samples and the known training samples, we were also able to demonstrate the improvement of the results of RBFN over GRNN as the number of training samples decreased.
-
-
-
Topographic effects in frequency-domain helicopter-borne electromagnetics
More LessAuthors Yutaka Sasaki and Hiroomi NakazatoMany helicopter-borne electromagnetic (HEM) surveys have been carried out in mountainous areas. However, there have been few studies of the effects of topography on the HEM response reported in the literature. We simulate the response to topography using a staggered-grid finite-difference method. Modelling shows that a hill produces a high-resistivity anomaly over its top, and a low-resistivity anomaly over its foot, when the magnetic field response is transformed into the apparent resistivity. In order to reduce these effects, a simple correction procedure is presented and tested on synthetic data. Results indicate that the corrected data do not reproduce the effects of the actual resistivity structure accurately enough to permit quantitative interpretation assuming a flat-earth model. The reason for this is that the geometrical relationship between the coil system and the subsurface structure changes. The most rigorous and accurate approach to interpreting HEM data with topographic effects is to incorporate a forward-solution scheme, capable of modelling topography, into inversions. A 3D inversion method is successfully tested on synthetic data.
-
-
-
Approximate 2D inversion of AEM data
More LessAuthors Peter Wolfgram, Daniel Sattel and Niels B. ChristensenAirborne electromagnetic (AEM) data are presently inverted with one-dimensional (1D) models, either as Conductivity Depth Images (CDI) or with full non-linear inversion, to build model sections from concatenated 1D models. If lateral conductivity changes are small, 1D models are justified. However, AEM investigations are often carried out specifically to find localized conductors, and in this case, 1D inversion is inadequate and will often produce artefacts in the model section.
We have developed an approximate two-dimensional (2D) inversion method that deals with laterally inhomogeneous sections. The method is based on the adaptive Born approximation previously applied by one of the authors (NBC) to the interpretation of central-loop ground EM profiles. The technique produces synthetic models with moderate conductivity contrasts and with some improvement over CDI sections. The computing speed is comparable to that of stitched 1D inversions.
An example of processing field data with the approximate 2D inversion method over a massive nickel sulphide deposit shows results that are promising for its routine application on large AEM data sets.
-
-
-
Cross-well electromagnetic imaging in three dimensions
More LessAuthors Michael S. Zhdanov and Ken YoshiokaIn this paper, we develop a new technique for 3D cross-well electromagnetic tomography, based on an EM borehole survey consisting of a moving vertical magnetic dipole transmitter, located in one or several boreholes, and a tri-axial induction receiver, located in the other boreholes. The method is based on the LQL approximation for forward modelling, which results in a fast inversion scheme. The method incorporates both a smooth regularized inversion, which generates a smooth image of the inverted resistivity, and a focusing regularized inversion, producing a sharp focused image of the geoelectrical target. The practical application of the method to synthetic data demonstrates its ability to recover the resistivity, location, and shape of resistive and conductive rock formations.
-
-
-
An assessment of the accuracy of boundaries picked by AEM sounding, with application to mineral exploration and salinity mapping
More LessAuthors James Macnae, Ashley Grant and Richard LaneThe vertical distribution of conductivity in the ground exhibits varied characteristics. Abrupt changes may occur at geological unconformities, and gradational changes are common with saturation and salinity gradients, or with clay content through the weathering profile.
This paper describes a rapid, automatic method of picking probable vertical locations for abrupt changes in conductivity using the conductance-depth curve derived in an approximate conductivity-depth transform of airborne EM data. The method is based on identifying the location of slope changes on a cumulative conductance vs. depth plot. The location of probable layer boundaries shows good agreement with synthetic data when the conductance contrast is adequate and boundaries are well separated. Correspondence of structures picked on separate x and z components increases the confidence of interpretation on a conductivity-depth image (CDI).
In field data, the quantitative usefulness of including probable layer boundaries on a CDI display is yet to be determined, qualitatively however it appears to provide very useful information to aid in the assessment of CDI section accuracy. It is also very useful where conductivity contrasts are poorly imaged by the limitations of colour bar or greyscale choice.
-
-
-
A comparison of shipborne and airborne electromagnetic methods for Antarctic sea ice thickness measurements
More LessAuthors James E. Reid, Julian Vrbancich and Anthony P. WorbyThe three-dimensional modelling program MARCO_AIR has been used to calculate the response of idealized sea-ice pressure ridge models to practical airborne and shipborne electromagnetic systems. The model results clearly show the superior lateral resolution of the horizontal coplanar shipborne system compared to higher-altitude airborne measurements. However, sea-ice keel thicknesses estimated via one-dimensional inversion of shipborne single-frequency electromagnetic data are strongly dependent on relatively small variations in survey altitude. One-dimensional inversion of synthetic helicopter electromagnetic data over three-dimensional pressure ridge models shows that the maximum ice keel thickness is consistently underestimated, although airborne EM methods yield reliable thickness estimates over level ice.
The vertical-coaxial coil survey geometry offers excellent lateral resolution of multiple targets, but the anomalies of typical Antarctic sea-ice pressure ridges would be too small to be reliably detected in practical surveys using an HEM system with a transmitter-receiver separation of 2-3 m. For an HEM system with a coil separation of 8 m, the vertical coaxial responses are larger, and lateral resolution of the vertical coaxial measurements at a flight height of 20 m is superior to a close-coupled horizontal coplanar system flown at an altitude of 10 m.
-
-
-
The application of fractional calculus to potential field data
More LessAuthors Gordon Cooper and Duncan CowanGradients of magnetic and gravity data are used routinely to enhance the edges of anomalies, or as input to interpretation techniques such as analytic signal analysis or Euler deconvolution. The most commonly used gradients are of first and second order, higher orders being used less frequently due to noise problems. This paper discusses the benefits of a generalized approach using fractional gradients, demonstrating their usefulness as an aid to interpretation, both on theoretical models and on aeromagnetic data from the Free State province, South Africa. Fractional horizontal gradients are suggested as a means of avoiding the instability problems present when magnetic data from low latitudes is reduced to the pole. They also allow the use of an improved sunshading algorithm that is less affected by noise than the standard method. Fractional vertical gradients may be used to generate both enhanced analytic signal data and enhanced Euler deconvolution solutions.
-
-
-
Some comments on potential field tensor data
More LessAuthors Philip Heath, Graham Heinson and Stewart GreenhalghThe measurement of potential field tensor (vector gradient) data is rapidly becoming a new tool for geophysical exploration. Data potentially provide significant improvements in resolution, noise suppression, data interpolation, and discrimination of subtle geological basement and regolith features, particularly in areas of high remanence, low latitudes, and in steep terrains. We have constructed algorithms in MATLAB for the three-dimensional inversion of potential field tensor data using Monte Carlo and Downhill Simplex approaches. We use these algorithms to invert simulated magnetic and gravity tensor data generated from simple geological structures, such as linear dykes and faults. The algorithms have a set target misfit (e.g., RMS misfit equal to one), and the final geological models are illustrated in three dimensions.
-
-
-
Wavelet analysis of detailed drillhole magnetic susceptibility data, Brockman Iron Formation, Hamersley Basin, Western Australia
More LessAuthors Duncan Cowan and Gordon CooperClosely spaced drillcore susceptibility measurements for the complete Dales Gorge Member of the Brockman Iron Formation were made as part of a study of oxide mineralogy and rock magnetic properties of unmineralized banded iron formations (BIFs) in the Hamersley Basin, Western Australia. Magnetic properties of BIFs are dominated by magnetite because of its much higher intrinsic magnetization even when the opaque oxides are volumetrically dominated by hematite. Bulk magnetic susceptibility varies considerably throughout the BIF units due to variations in thickness and composition of oxide/silicate mesobands, variations in the proportion of microbands, and changes in the hematite/magnetite ratio. This local-scale variability and the presence of layering effects at several scales means that the drill-hole susceptibility data are non-stationary, and conventional tools for analysis such as autocorrelation, spectral analysis, and semivariogram techniques have limited value.
Both continuous wavelet (CWT) and discrete wavelet (DWT) techniques proved to be very powerful for analysis of these non-stationary data. CWT analysis using a Morlet wavelet helped to quantify the periodicity of the layering and DWT analysis using a Haar wavelet provided an effective means of ‘blocking’ the log so that susceptibilities of individual BIF macrobands can be compared. The results suggest that wavelet analysis can provide definitive information on the scaling properties of magnetic susceptibility needed to determine ‘formation response’ or facies variation within the major BIF units. Rapid changes in magnetic susceptibility in BIFs are mirrored by the oxide mineralogy. Opaque-oxide petrography and rock magnetism results show that replacement of depositional or diagenetic hematite by magnetite within and between sedimentary band structures by reduction during prograde metamorphism was incomplete, and equilibrium in the magnetite stability field was not reached even at the mesoband scale.
-
-
-
Improvement to long-wavelength Australian gravity anomalies expected from the CHAMP, GRACE and GOCE dedicated satellite gravimetry missions
More LessAuthors W. E. FeatherstoneThe concepts and mission parameters of the CHAMP, GRACE and GOCE dedicated satellite gravity missions are summarised, followed by the improvements that they are likely to make upon previous methods of determining the Earth’s long-wavelength gravity field. An example of an error in the EGM96 global geopotential model, caused by the use of an incorrect digital elevation model over Australia, is used to exemplify the deficiencies in current global gravity models. Summarised results of a preliminary study to quantify long-wavelength errors in Australian gravity anomalies using data from the GRACE-derived EIGEN-2 global geopotential model indicate the presence of long-wavelength (>1113 km) errors of over 10 mGal in Australian terrestrial gravity anomalies. Finally, the likely prospects for improved geophysical studies are speculated upon.
-
-
-
The importance of accurate altimetry in AEM surveys for land management
More LessAuthors Ross Brodie and Richard LaneAirborne electromagnetic (AEM) systems are increasingly being used for mapping conductivity in areas susceptible to secondary salinity, with particular attention to near-surface predictions (conductivities in the top 5 or 10 metres). Because measured AEM response is strongly dependent on the height of both the transmitter loop and receiver coil above conductive material, errors in measurements of terrain clearance translate directly into significant errors in predicted near-surface conductivity.
Radar altimetry has been the standard in airborne geophysical systems for measuring terrain clearance. In areas of agricultural activity, artefacts up to five metres in magnitude can be present. One class of error, found to correlate with ploughed paddocks and hence termed the "paddock effect", results in overestimation of terrain clearance. It is hypothesised that this is due to surface roughness and soil moisture levels in these paddocks. A second class of error, related to dense vegetation and hence termed the "canopy effect", results in underestimation of terrain clearance.
A survey example where terrain clearance was measured using both a radar and a laser altimeter illustrates the consequences of the paddock and canopy effects on shallow conductivity predictions. This shows that the combination, of the dependence of AEM response on terrain clearance and systematic radar altimeter artefacts spatially coincident with areas of differing land-use, may falsely imply that land-use practices are the controlling influence on conductivity variations in the near surface.
A laser altimeter is recommended for AEM applications because this device is immune to the paddock effect. Careful processing is still required to minimise canopy effects.
-
-
-
Airborne gravimetry – a new gravimeter system and test results
More LessSuccess of airborne gravity surveys mainly depends on determining the three-dimensional (3D) position of the moving platform. Recent advances in technology, especially the Global Positioning System (GPS), have made it possible to determine the velocity and position of the moving platform more frequently and with greater accuracy. Taking advantage of these advances in GPS technology, and using a newly developed system, helicopter-borne gravity measurements were successfully carried out over the Kanto and Tokai districts of Japan. This new gravimeter system is composed of servo accelerometer sensors, a stabilised platform, an optical-fibre gyroscope to control the stabilised platform, GPS receivers, and a data processor.
The 3D position of the helicopter at every second was accurately determined by the interferometric GPS method. These GPS data were also used to compute various correction factors which are applied to the measured gravity acceleration. Real-time differential GPS positioning was also conducted using a separate receiver mounted on the helicopter. These real-time positioning data were used for controlling the optical-fibre gyroscope. The gravity acceleration data were processed and all necessary corrections were applied. Numerical filtering was carried out to remove high-frequency noise in the data. The observed free-air gravity anomalies were then compared with upward continuation of the ground gravity data to the flight altitude. We also compiled an airborne gravity anomaly map from the airborne data, which was compared with upward-continued ground gravity data.
-
-
-
Drape-related problems in aeromagnetic surveys: the need for tight-drape surveys
More LessAuthors Duncan Cowan and Gordon CooperThe trend towards higher flight height and loose-drape fixed-wing surveys has important implications for high-resolution aeromagnetic surveys for mineral exploration, especially for targets such as kimberlites. Recent trends regarding safety considerations, implemented by leading airborne geophysical contractors and some Governments, make it virtually impossible to acquire adequate data in areas of even moderate relief using conventional aircraft.
Variations in survey elevation give rise to changes in magnetic relief and magnetic texture that are not related to magnetic sources but are simply artefacts produced by varying depth-to-source. Important high frequency and low amplitude magnetic signals are absent in areas of high survey elevation and cannot be recovered by drape corrections, which involve intrinsically unstable downward continuation.
Geophysicists are well aware of the loss of horizontal spatial resolution as survey elevation increases and the loss of subtle information on depth extent may be equally serious in exploration for pipe-like bodies such as kimberlites. Aeromagnetic data are relatively insensitive to depth extent because of the rapid decay of the kernel function with depth. At shallow depths there is a clear difference in anomaly shape between a pipe source and a thin sheet but as depth-to-source increases these differences become much less obvious.
For cost-effective high-resolution aeromagnetic surveys for kimberlites we need to fly tight-drape surveys at low level, using aircraft designed to fly at low level such as the Pacific Aerospace Cresco 750 and Fletcher FU-24. The terrain-following ability of these aircraft is impressive.
Comparison of tight-drape and simulated loose-drape data shows a significant loss of high frequency content for loose-drape surveys over shallow magnetic sources. The loose-drape examples show significant variations in magnetic resolution due entirely to variations in survey elevation. Profile-based drape corrections removed most of these artefacts but fail to restore high frequency signals.
Analysis of magnetic signatures of kimberlite models shows a significant deterioration in detection limits as survey elevation increases from 80 m to 150 m. Model studies comparing pipe models with thin-sheet models shows that our ability to distinguish pipe and sheet sources deteriorates as source depth increases. At low level there are clear differences in profile shapes, with the sheet models showing rapid decay away from the source, and pipe models showing slower decay rates. However, as source depth increases these differences are much less obvious.
-
-
-
Airborne vector magnetics mapping of remanently magnetized banded iron formations at Rocklea, Western Australia
More LessAuthors Mark Dransfield, Asbjorn Christensen and Guimin LiuStrongly remanent magnetic sources such as those found in the banded iron formations (BIFs) of the Hamersley Basin impair the interpretation of standard total magnetic intensity survey data. The use of aircraft attitude information, provided by the FALCONTM airborne gravity gradiometer (AGG) system, makes it possible to reference geographically the vector magnetic information collected by the fluxgate triad on board the aircraft. A full aircraft compensation processing method, allowing for permanent, induced, and eddy current effects, provides usable vector magnetic data. These data supply valuable extra information for the interpretation of strongly remanent BIFs.
An example using data from a FALCON survey near the Rocklea Dome in the Hamersley Basin demonstrates improved mapping of the BIFs by the use of the vector-residual magnetic intensity data.
-
-
-
Quietening the noise: an evaluation of noise reduction methods applied to aerial gamma-ray survey data
More LessAuthors B. L. Dickson and G. M. TaylorSynthetic survey data have been used to evaluate five aerial gamma-ray survey noise-reduction methods, four of them in industry use (NASVD, MNF, eMNF and Mathis adaptive filtering), as well as a new method, CSM. The CSM (short for CSIRO Self Organizing Map) method is based on deriving a self-organizing map from the full spectral data set and replacing individual spectra by others that best match them in a vector distance sense.
Synthetic surveys of 6400 spectra each were prepared to represent:
Case A: uniform ground,
Case B: ground with varying K-U-Th but fixed Th/U ratio,
Case C: ground as for B but with an additional area of different Th/U ratio and,
Case D: ground with uncorrelated K-U-Th concentrations.
For Cases A and B, the CSM and eMNF methods gave similar noise cleaning and best revealed the known Th-U correlation in Case B.
For Case C, both MNF and eMNF showed ghosting in the Th/U ratio in some areas. NASVD did not return the correct Th/U ratio in the area of different Th/U ratio. CSM showed only a small area of ghosting and gave the overall best cleaning. Pre-zoning the data by Th/U ratio and separate processing of the zones is shown to enable the MNF and eMNF methods to avoid ghosting.
For Case D, CSM performed poorly for individual radioelements, whilst the other methods gave some cleaning. However, for the Th/U image, the CSM method gave the best cleaning. This case illustrates the importance of some correlations in the data, for any of the spectral noise cleaning methods to give good results.
-
-
-
Three-dimensional kinematic imaging of borehole radar data
More LessAuthors Naomi Osman, Carina Simmat, Jonathan Hargreaves and Iain MasonBroadband, 10-100 MHz, VHF borehole radar (BHR) can be used to map ore bodies, faults and marker horizons, and to identify hazards in advance of mining. Mine boreholes are often drilled in fans. BHR data acquired from these boreholes can be used to reconstruct targets in 3D. Significant progress has been made in using synthetic aperture radar interferometry (InSAR) to reconstruct 3D images from sparse arrays. However, automatic methods of projecting data into 3D image space, such as migration and InSAR, make stringent demands upon rock homogeneity, translucence, and the accuracy of borehole trajectories. These demands can be relaxed by kinematic mapping using geologically plausible 3D primitives such as cylinders, planes, and spheres. In this paper we show how interactive kinematic mapping can be applied in practice to BHR data acquired in a Western Australian mine. BHR data were shot from a fan of ~200 m long boreholes. The holes were drilled from a single station at relative angles of ~5°. Both cross-hole transmission and single-hole reflection surveys were performed. Reflection ranges of 50-60 m were achieved in the rocks hosting the deposit. Echo patterns correlated using kinematic primitives were used to infer the geometry of an interface that was interpreted to be approximately 20 m above the uppermost borehole in the array.
-
-
-
On the use of magnetics and gravity to discriminate between gabbro and iron-rich ore-forming systems
More LessAuthors James E. HannesonIf the major rock-forming minerals are divided into three categories, and, if the density, magnetic susceptibility, and the proportions of each category are known, then the density and susceptibility of a mixture of the three categories can be determined. Conversely, if the physical properties of a mixture are known along with the physical properties of each of the three components, then the proportions of each component can be calculated. This reasoning can be used to superimpose category-percentage contours onto a density-susceptibility scatter diagram. Such a diagram can be used to relate petrophysical measurements on a rock specimen to its mineral category content, the veracity of which can be appraised by a geologist with a hand lens. Several hundred qualitative tests suggest that the approach is valid. The analysis is independent of scale, and can be applied to the bodies of a density-susceptibility model developed to simulate exploration data. Gabbro with magnetite plots along a specific locus on the combined phase/scatter diagram and can often be distinguished from denser accumulations of hematite and/or sulphides. Case histories for a number of economic and non-economic density-susceptibility anomalies give support for the usefulness of the method.
-
-
-
Minerals exploration methods modified for environmental targets
More LessAuthors Norman R. Carlson and Kenneth L. ZongeIn recent years, several traditional mining exploration methods have been successfully modified for use in the environmental geophysics field. In these cases, "successfully modified" refers primarily to acquiring data fast enough, and therefore economically efficiently enough, to accommodate the relatively small budgets that are available in most environmental studies. Modifications also include adjustments for the smaller scales of environmental surveys, often in the size of the project area, size of targets themselves, and required resolution.
For example, transient electromagnetics (TEM) methods have been increasingly applied to environmental problems, particularly in unexploded ordnance (UXO), underground storage tank (UST), and utilities detection. A major research effort is now underway to use, among other methods, multi-component, multi-time-gate mobile TEM systems (measuring Hx and Hy, as well as the standard Hz) in order to discriminate targets of interest (UXO, for example) from other anomalies (such as metallic debris).
A second good example is the induced polarization (IP) method. Although resistivity has been used extensively in shallow environmental applications, IP data acquisition has always been too slow, and therefore too expensive, for most environmental targets. Multi-channel receivers, multiplexers, and laptop computers now allow us to acquire IP data at rates of 2500 to 3000 data points per day (in the dipole-dipole configuration, for example), providing low-cost, high-density data. IP data have been shown to be particularly useful in delineating buried waste, such as at old landfills.
-
Volumes & issues
-
Volume 56 (2025)
-
Volume 55 (2024)
-
Volume 54 (2023)
-
Volume 53 (2022)
-
Volume 52 (2021)
-
Volume 51 (2020)
-
Volume 50 (2019)
-
Volume 49 (2018)
-
Volume 48 (2017)
-
Volume 47 (2016)
-
Volume 46 (2015)
-
Volume 45 (2014)
-
Volume 44 (2013)
-
Volume 43 (2012)
-
Volume 42 (2011)
-
Volume 41 (2010)
-
Volume 40 (2009)
-
Volume 39 (2008)
-
Volume 38 (2007)
-
Volume 37 (2006)
-
Volume 36 (2005)
-
Volume 35 (2004)
-
Volume 34 (2003)
-
Volume 33 (2002)
-
Volume 32 (2001)
-
Volume 31 (2000)
-
Volume 30 (1999)
-
Volume 29 (1998)
-
Volume 28 (1997)
-
Volume 27 (1996)
-
Volume 26 (1995)
-
Volume 25 (1994)
-
Volume 24 (1993)
-
Volume 23 (1992)
-
Volume 22 (1991)
-
Volume 21 (1990)
-
Volume 20 (1989)
-
Volume 19 (1988)
-
Volume 18 (1987)
-
Volume 17 (1986)
-
Volume 16 (1985)
-
Volume 15 (1984)
-
Volume 14 (1983)
-
Volume 13 (1982)
-
Volume 12 (1981)
-
Volume 11 (1980)
-
Volume 10 (1979)
-
Volume 9 (1978)
-
Volume 8 (1977)
-
Volume 7 (1976)
-
Volume 6 (1975)
-
Volume 5 (1974)
-
Volume 4 (1973)
-
Volume 3 (1972)
-
Volume 2 (1971)
-
Volume 1 (1970)
Most Read This Month