ASEG Extended Abstracts - 25th International Conference and Exhibition – Interpreting the Past, Discovering the Future, 2016
25th International Conference and Exhibition – Interpreting the Past, Discovering the Future, 2016
- Articles
-
-
-
Lithological mapping via Random Forests: Information Entropy as a proxy for inaccuracy
More LessAuthors Stephen Kuhn, Matthew J. Cracknell and Anya M. ReadingMachine Learning Algorithms (MLA) can be an effective means of lithological classification. The Random Forests™ (RF) supervised classification approach allows prediction of lithology from disparate geophysical, geochemical and remote sensing data. In this study, we examine the relationship between prediction accuracy and information entropy (H). Data were processed in accordance with industry best practice and input selection was optimised using RF. Using a training set containing 1.4% of available pixels, we produced a classified lithology map with an overall accuracy of 76% with regards to mapped geology. In addition, we produced a class membership probability for each pixel, a precursor to defining the ultimate class designation at each pixel. H was calculated at each pixel from output class membership probabilities; and in this context provides a measure of the state of disorder for each. H was normalised with 0-1 representing the minimum to maximum possible H for each pixel.
H equal to 1 at a pixel represents an equal probability of all candidate classes occurring, whereas H equal to 0 describes a 100% probability of single class occurring. In this study, we demonstrate that there is a significant difference in the distribution of H between correctly and incorrectly classified pixels. The median H of incorrectly classified samples occurs above the 75% percentile of H for correctly classified samples. Conversely, both the mean and median H for correctly classified pixels occurs below the 25% percentile level for incorrectly classified samples.
This information can be used to determine the well-defined transition range in H, above which classification is likely to be inaccurate. Using this approach, a geoscientist can produce a lithological map, a quantifiable measure of uncertainty and a quantifiable transition range above which they are likely to encounter incorrect classification, avoiding wasted expense in targeting based on an incorrect model.
-
-
-
-
Analysis of Electromagnetic Depth Sounding Responses over a Layered Earth: Investigating Oil & Gas Seeps in the Petroleum Provinces
More LessAuthors Shastri L Nimmagadda and Andrew OchanThe present work embodies the results of theoretical and practical investigations of electromagnetic depth sounding using central frequency sounding (CFS) system over a layered earth. Failure of conventional electrical resistivity sounding in the study of geological conditions under resistive overburden calls for variable frequency sounding techniques. Electromagnetic depth sounding which involves the measurement of variation in conductivity with depth is used for solving various geoengineering, hydrogeological and shallow cases of oil & gas seeps associated with stratified earth. The CFS, which is one of the depth sounding techniques involving the measurement of vertical component of magnetic field induced at the centre of a circular or square loop, is considered in the present study for obtaining theoretical responses over a layered earth and its interpretation with shallow oil and gas seeps.
Because of some limitations of contour integration and numerical integration approaches, used earlier, a more rapid digital linear filter technique is adopted for evaluation of the integral involved in the CFS theoretical expressions. Theoretical expressions for frequency-domain soundings written for layered earth models are suitably transformed for evaluation through digital linear filter. Dimensionless normalized vertical magnetic field is computed for different frequencies and loop radii for layered earth models with different layer conductivities and thicknesses. The responses computed for these cases are analysed in terms of resolution characteristics and detectability effects. In frequency-domain sounding, amplitude response curves of layer-sequences show the effect of layer conductivity, layer thickness and loop radius. Separation between individual curves on the sets for amplitude responses normally gives sufficient indications for subsurface conductivity variations of the layered earth cases. The author explores the CFS applicability and feasibility in investigating shallow oil & gas seeps in oil & gas provinces, in particular on the flanks of the rifted grabens and basin margin areas, where sediment - basement contact areas are interpreted.
-
-
-
Characterising cover and exploring under cover with AEM
More LessAuthors Shane MuléCharacterisation of cover is an important aspect of understanding many geological systems. Through a better understanding of cover and improved sensitivity to deeper structures the footprint of mineral systems can be extracted. This is important in providing greater confidence to mineral explorers and will assist with further exploration success at depth. UNCOVER is a government initiative, focussed on improving the mineral prospects of Australia through better exploration of the subsurface under cover. AEM is well placed in supporting this initiative, by providing a rapid and efficient way to undertake near surface geophysical exploration. To accurately characterise the subsurface, AEM systems must have both accuracy and precision which are achieved through robust calibration and low noise levels. Tempest has a long history of exploring cover and under cover. The characteristics of the system are tuned to provide a broad bandwidth which maximises its sensitivity and applicability. Years of development has led to significant improvements in signal-to-noise along with deployment of system variants and platform diversity, including High Moment Tempest which provides increased transmitter moment along with reduced noise levels leading to greater depth of penetration and sensitivity to low conductivity contrasts.
-
-
-
Designing Adiabatic Pulses for Surface NMR
More LessAuthors Denys Grombacher and Esben AukenSurface nuclear magnetic resonance (NMR) is a powerful technique providing non-invasive imaging of groundwater. One challenge with the method is that it commonly suffers from low signal to noise ratios (SNR). Two methods to increase SNR are to either develop noise cancellation approaches to reduce the noise level, or to perform the experiment in a manner capable of increasing the signal amplitude. A recent study adopted the latter approach by employing a novel transmit strategy. An adiabatic pulse was employed and observed to produce significant signal improvements compared to the standard transmit method in surface NMR. The advantage of an adiabatic pulse is that it is capable of producing a uniform excitation in the presence of a heterogeneous magnetic field, which describes exactly the transmit conditions in surface NMR.
Given the great potential of adiabatic pulses for surface NMR, we explore several factors related to the design of adiabatic pulses intended for application in surface NMR conditions. We investigate how various adiabatic pulses perform in a heterogeneous magnetic field given the limitation that current instrumentation couples the modulation of the current amplitude during the pulse to the instantaneous transmit frequency. That is, only the duration of the adiabatic sweep, the bandwidth through which the pulse sweeps, and frequency modulation throughout the pulse can be directly controlled. A numerical sensitivity analysis of each of these parameters is performed to gain insight into how to design optimal adiabatic pulses for surface NMR. Additionally, a numerically optimized modulation (NOM) approach is implemented to optimize the frequency sweep. The spatial resolution and depth penetration provided by an example adiabatic pulse is also investigated. A trade off between signal amplitude and spatial resolution is observed to be present when employing adiabatic pulses.
-
-
-
A Statistical Approach to Assessing Depth Conversion Uncertainty on a Regional Dataset: Cooper-Eromanga Basin, Australia
More LessAuthors David Kulikowski, Catherine Hochwald, Dennis Cooke and Khalid AmrouchDeciding on the most accurate grid based depth conversion method can often be an arbitrary choice made by geophysicists, particularly if previous research is limited. The importance of accurate depth conversion is particularly crucial in the Cooper-Eromanga Basin, where the presence of oil rich, low relief structural traps are questionable depending on the method used. Previous depth conversion studies are limited to local scales, limited well control and a focus on select horizons. To investigate the depth conversion uncertainty on a regional scale, this research performs a comprehensive and regionally extensive depth conversion analysis utilising 13 3D seismic surveys with 73 interpreted TWT grids and 657 wells. Depth conversions were performed using the 4 most commonly used methods; average (pseudo) velocity, time-depth trend, kriging with external drift using TWT, and kriging with external drift to tie stacking velocities to average well velocities.
To manage the large volume of data, a looping script was written to automate the depth conversion process and utilise the cross-validation, or blind-well method (use n wells to predict the nth +1 well). Statistics on several variables were captured after each loop, with cluster analysis performed on the final data set to test variable significance on depth conversion accuracy. A database of approximately 10000 error calculations found that although the average velocity method is the most accurate at a high level (average absolute error ~24.9 feet), the best method and the expected error changes significantly (tens of feet) depending on the combination and value of the most significant variables. The variables which impacted uncertainty the most were; location (3d survey), formation, distance to the nearest well, and the spatial location of the predicted well relative to the existing well data envelope.
-
-
-
Structural Interpretation of seismic, geological realism and 3D thinking
More LessAuthors Pete Boult, Brett Freeman and Graham YieldingOur notion of reality in seismic interpretation and structural geology usually follows a series of careful observations and ideas that eventually crystallize into a best-case model. In most other branches of science the strength or reality of such models (or hypotheses) is increased by the number of robust tests that either refine or fail to disprove the original idea. However, geological models in the Exploration and Production (E&P) sector differ because the starting point for testing is usually an interpretation of seismic or other remote measurements, rather than the direct observation of an effect. This means that whatever tests we are able to apply have, in themselves, significant margins for error and are further compounded by the intellectual issue of constructing a 3D view or model of the perceived geologic reality. This model building is an early stage of the E&P process, but errors and uncertainty at this point propagate throughout the subsequent workflow and arguably, amount to the single biggest factor affecting perceived value and in particular drilling risk.
-
-
-
P- and PS-wave vector wavefields for anisotropic petrophysics
More LessAuthors James GaiserPredicting petrophysical properties of lithology, density and fractures using seismic data is an essential part of reservoir evaluation. In addition to lithology characterization, “frackability” has become a very important area of investigation for high-grading locations for drilling and hydraulic fracture stimulation. Most seismic studies that estimate P- and S-wave impedance, density and brittleness or formation strength use conventional P-wave data and isotropic elastic inversion methods. However, converted-wave (PS-wave) joint inversion and S-wave splitting methods have successfully been used to improve determination of seismic properties for shale plays as well as other unconventional resource plays.
Anisotropic behaviour related to layered media (VTI), fracture properties, stress direction and the geomechanics of shales are in-creasingly more important for seismic analysis, imaging and reservoir characterization. Vector wavefields are sensitive to these prop-erties and can help identify optimal drilling and stimulation locations. Also, it has been shown that use of conventional elastic pa-rameters for characterizing “brittleness” should include anisotropic corrections to obtain a more accurate response. Including PS-wave seismic data is beneficial for isotropic elastic inversion and should improve anisotropy estimates for identification of potential fracture locations.
Elastic inversion of azimuthally anisotropic amplitude variations (AVAz) is also becoming more important. When layered media are fractured, orthorhombic symmetry of P-wave amplitude depends on S-wave birefringence. PS-waves are ideal for determining this S-wave splitting information from layerstripping and their reflectivity provides additional equations for joint inversion with P-waves. Two coefficients, a radial RPSV and transverse RPSH reinforce anisotropic signatures similar to P-wave reflectivity RP. Vector wave-fields contain all the necessary information for S-wave anisotropy from short wavelength AVAz as well as from long wavelength moveout behaviour.
-
-
-
Geophysics for Ni-Cu — Where are we at and where are we going?
More LessAuthors W. S. PetersExploration for Ni-Cu sulphide deposits has been possibly the most significant factor in the development of modern minerals geophysical technology. The strong response of Ni-Cu sulphides to multiple geophysical methods has resulted in a heavy reliance on geophysics for exploration. This has, in turn, spurred researchers, instrument manufacturers and contractors to continually develop new and/or improved technology.
As discovery of outcropping and near-surface deposits is becoming rarer, the challenge is moving into covered areas and exploring deeper. Historical technology has proven effective in the 0-250m depths; current technology has extended this to more than 500m; new developments are pushing this towards 1000m.
Significant advances in acquisition have been facilitated by factors such as faster processors, large data storage, faster A-D converters, miniaturisation of sensors, UAVs, integrated GPS and communications. Key developments, particularly in electrical and EM technology include high power transmitters, full waveform recording, and distributed arrays. Interpretation advances include small and large scale 3D inversion, integrated 3D visualisation, and advanced flexible modelling.
-
-
-
Development of Rapid Scanning Surface-NMR for Wide Area Hydrogeologic Mapping
More LessAuthors Elliot Grunewald and David O. WalshSurface nuclear magnetic resonance (Surface-NMR) measurements hold the valuable capability to directly image groundwater and to characterize aquifer flow and storage properties. Historically, implementations of surface-NMR have been limited by long stacking times and slow survey deployment, restricting applications primarily to 1D soundings and short 2D profiles. Through advancements in acquisition schemes, hardware, and deployment platforms, we demonstrate the ability to deploy surface-NMR as an efficient wide-area mapping technique. To increase measurement efficiency, we have developed acquisition schemes to improve the inherently low signal-to-noise-ratio of Earth's field NMR. Adiabatic pulse sequences are used to increase the detected NMR signal amplitude and to reduce requirements to scan over a wide range of pulse moments. Newly developed wireless noise-reference coil stations are used to cancel environmental noise without increasing the size or footprint of the signal-detection array. Smaller footprint wired signal-detection arrays are transported efficiently using mobile platforms, including towed coil mats and elevated coil forms. The detection array can be moved along a profile line and left in a static position for short time intervals to acquire measurements before being moved to the next position. These newly developed rapid scanning NMR technologies are demonstrated at a collection of sites in the Western United States.
-
-
-
3D imaging of the Earth's lithosphere using noise from ocean waves
More LessAuthors Yingjie Yang, Jun Xie and Kaifeng ZhaoThe lithosphere is the rigid outer shell of the Earth, composed of the crust and the rigid uppermost mantle. The lateral variation of lithosphere is believed to be closely related to the foci of intraplate earthquakes and volcanism, the location of large ore deposits and diamond-bearing kimberlites, and the formation of oil/gas-bearing sedimentary basins. Understanding the fine-scale structure of the lithosphere is therefore one of the most fundamental tasks in Earth sciences, and has important implications for society in terms of economic prosperity and hazard mitigation.
Seismic tomography is the main technique available to image the subsurface structure of the Earth across a range of scales. Conventional seismic tomography exploits the seismic waves emitted by earthquakes. However, because most large earthquakes occur at plate boundaries and most continents, including Australia, are not seismically active, earthquake-based tomography suffers from uneven distributions of earthquakes and has difficulties in deciphering fine-scale structures of the Earth in seismically quiet continents.
In the past decade, the advent of ambient noise tomography (ANT) has revolutionized seismic tomography because it can overcome the limitations of conventional earthquake surface wave tomography. This technique uses diffuse background seismic energy, mostly comes from the interaction of ocean waves with the crust. Empirical Green's functions of surface waves passing between a pair of stations are extracted by cross-correlating continuous time series of ambient noise. Within a regional seismic array, all inter-station measurements of surface-wave dispersion can be measured and tomography can be performed to image the underling lithospheric structure. In this study, we demonstrate the applications of broadband ANT in mapping fine-scale lithosphere structures around the world using continuous seismic data.
-
-
-
Microseismic characterization of brittle fracture mechanism in highly stressed surrounding rock mass
More LessAuthors Chunchi Ma, Tianbin Li, Yupeng Jiang and Guoqing ChenBrittle fracturing of rock mass is a major problem for deep tunnelling or mining in highly stressed rock mass, which could evolve into rockburst hazard and severely undermine the safety of engineering project. In this paper, an advanced microseismic approach is proposed to analyse the mechanism clearly. The results suggest that brittle fracturing contains three energy development stages (i.e. energy accumulation, energy transfer and energy release). Therefore, based on the seismic energy moment and the apparent stress criteria, microseismic events can be classified into six categories that corresponding to the three major categories of stress-adjustment event, deformation-driven event and bursting event. The Energy index develops steadily at first, then followed by a drastic drop and finally ends with a large increment. Cumulative apparent volume grows slowly before a sudden increase. For brittle fracturing development, tensile cracks appear firstly, then the shear cracks and finally the mixed types of cracks.
-
-
-
Achieving accurate interpretation results from full-waveform streamed data AEM surveys
More LessAuthors Magdel Combrinck and Richard WrightSuccessful interpretation of airborne electromagnetic (AEM) data depends strongly, but not exclusively, on the quality of data. The ability to accurately describe the system, knowledge of the data processing procedures as well as an understanding of the geological targets all contribute to derive accurate models. Xcite™ is a recently developed helicopter time domain electromagnetic system featuring an inflatable frame and full-waveform recording of streamed data. Working with the streamed data allows for restructuring of traditional data processing routines to customize data sets for specific applications. It also allows for a better understanding of the corrections that are required to isolate the secondary field response from measurements acquired on a helicopter platform. An accurate system description as well as transparency in the data processing phase allows interpreters to get the most value from airborne electromagnetic data.
-
-
-
Geophysical Responses from Mineral System Components in the Deep Crust and Upper Mantle
More LessAuthors Michael DentithThe concept of a mineral deposit forming via a mineral system that operates across areas of perhaps 1000s of squares kilometres and to mantle depths has important implications for greenfields mineral exploration. Geographically widespread datasets and deep penetrating geophysical methods are required to map key mineral system elements such as fluid/metal source zones and migration paths. Developed primarily for academic studies of the deep crust, there are several established geophysical techniques that can potentially be used to identify elements of mineral systems in the deep crust and upper mantle. Although the seismic reflection method produces the highest quality images, it is prohibitively expensive and the recommended approach is a combination of MT surveys and receiver function recordings with CCP stacking. Mineral system elements that can be detected in this fashion include major structures and geological boundaries which are potential controls on fluid flow and also areas of crust and mantle that have been altered by one or both of fluid creation and migration
-
-
-
Application of Fullwaveform Tomography to VSP Walkaway Data
More LessAuthors Eric Takam Takougang and Youcef BouzidiVSP walkaway data were collected in an oil field in the United Arab Emirates. 2D frequency domain fullwaveform tomography was used to obtain high resolution rock properties (P-wave), and structural information near and away from the borehole, following a specific data preconditioning and inversion strategy. The field data were inverted between 4 and 50 Hz, and the starting model was obtained from traveltime tomography. The results of the inversion show zones with anomalous low velocities that correlate in places with known presence of hydrocarbons, and highlight their possible extensions. A comparison of the results at the borehole location with the sonic log shows a generally good match. However, some mismatches are evident and can be explained by possible anisotropy and out of the plane structures not taken into account during the inversion.
-
-
-
The emperor's new clothes - opportunities and limitations applying AEM to geotechnical design work
More LessAuthors Andi A Pfaffhuber, Helgard Anschütz and Kristoffer KåsinIn summer 2015, we acquired close to 6.000 km of Helicopter, time-domain airborne electromagnetic (AEM) data for regional geotechnical mapping for the Norwegian National Rail Administration. This survey and further experience from related Norwegian road planning projects demonstrated the unprecedented accuracy of modern AEM data. The extent of geotechnical site investigations can be drastically reduced, both in terms of time schedule, and costs if AEM derived bedrock models are included when soil investigations are planned. Geotechnical projects demand high resolution (meter scale) and AEM data is to some extent capable of delivering that. Some of our data matched the resolution of corresponding ground geophysics data. Here we present the way in which AEM can be used as bedrock models, sensitive clay delineation and to determine bedrock types. Our discussion leads us to the missing link between high vertical resolution in the first tens of meters for geotechnical work and the focus on simple, sub-vertical structures in exploration AEM. Ultimately, we should strive for the best of both worlds, shouldn't we?
-
-
-
A new source parameters estimation method of airborne gravity gradient tensor data
More LessAuthors Shuai Zhou and Jian JiaoOne of the important tasks of potential field interpretation is source parameters estimation of the geological structures, including the horizontal position and buried depth. A new method to interpret airborne gravity tensor data is proposed in this paper based on Normalized Downward Continuation (NDC) of the tensor data directional total horizontal derivatives. The NDC method was introduced for source depth estimation, which can be applied to analytical signal modulus and potential fields themselves. And the maxima of the NDC map mainly correspond to the centre location of the geologic source. We applied the NDC method to the directional total horizontal derivatives which can be used for delineating the sources' horizontal edges. The maximum values of the result indicate the source edges horizontal position and buried depth simultaneously. During the calculation, the iteration method of the downward continuation is used in the NDC calculation process to improve the stability. The new method was tested on synthetic models and obtained satisfactory results. Compared with previous work, this new approach has a better lateral resolution.
-
-
-
X-ray Computed Tomography of Structures in Opalinus Clay from Large Scale to Small Scale after Mechanical Testing
More LessAuthors Gerhard Zacher, Annette Kaufhold and Matthias HalischIn recent years use of X-ray Computed Tomography (CT) has become more common for geoscientific applications and is used from the um-scale (e.g. for investigations of micro-fossils or pore scale structures) up to the dm-scale (full drill cores or soil columns). In this paper we present results from CT imaging and mineralogical investigations of an Opalinus Clay core on different scales and different regions of interest, emphasizing especially the 3D evaluation and distribution of cracks and fractures and their impact upon mechanical testing of such material. Enhanced knowledge of the behaviour of the Opalinus Clay as a result of these tests is of great interest, especially since this material is considered for a long term radioactive waste disposal and storage facility in Switzerland. Hence, results are compared regarding the mineral (i.e. phase) contrast resolution, the spatial resolution, and the overall scanning speed.
With this extensive interdisciplinary top-down approach it has been possible to characterize the general fracture propagation in comparison to mineralogical and textual features of the Opalinus Clay. Additionally, and to the best of our knowledge, a so called mylonitic zone has been observed for the first time in an experimentally deformed Opalinus sample. The multi-scale results are in good accordance with data from naturally deformed Opalinus Clay samples, which allows systematic analysis under controlled laboratory conditions. Accompanying 3D imaging greatly enhances the capability of data interpretation and assessment of such material.
-
-
-
Performance of Hankel transform filters for marine controlled-source electromagnetic surveys: a comparative study
More LessAuthors Hangilro Jang, Hee Joon Kim and Myung Jin NamFor the numerical calculation of electromagnetic responses by an arbitrary source in a one-dimensional model, integral-equation based approaches can be the best method since they are semi-analytic and thus very fast and accurate. In the numerical computation of the integral-based approaches, digital linear filters (e.g., Anderson’s filter, Kong’s filter and Mizunaga’s filter) play a key role. Using a closed-form solution of the Hankel transform in transverse magnetic mode for a homogeneous half-space model, we can assess the accuracy of digital linear filters for evaluating the Hankel transform. In this paper, we conduct comparative performance tests on the linear filter with known integral transforms. We examine three kinds of filters developed by Anderson, Kong and Mizunaga, which are known to be suitable for marine controlled-source electromagnetic (CSEM) applications. Kong’s filters perform best in the three kinds of filters over a practical range of offset distances in marine CSEM surveys. While the relative error versus distance appears as a V-shaped curve in semi-log scale, Mizunaga’s filters are the shortest in length and have a performance comparable to Kong’s filters. Anderson’s filters have a quite similar performance between the J0 and J1 filters, although these are somewhat inferior to Kong’s and Mizunaga’s filters when computing marine CSEM fields.
-
-
-
The Balboa ZTEM Cu-Mo-Au porphyry discovery at Cobre Panama
More LessAuthors Jean M. Legault, Chris Wijns, Carlos Izarra and Geoffrey PlastowThis paper describes the ZTEM airborne EM and magnetic results over the porphyry copper-gold camp at Cobre Panama, with a focus on the discovery of the blind Balboa deposit in 2010, the first documented case attributed to the ZTEM system. Balboa is the westernmost of six porphyry copper-gold deposits that make up Cobre Panama and it escaped detection in 40 years of exploration that relied primarily on soil geochemistry, airborne magnetics and drilling. ZTEM was flown in summer 2010 to detect resistivity variations related to hidden porphyry systems in a region of dense jungle, difficult access and thick (20-30 m) conductive saprolite cover. The ZTEM survey detected all the known porphyry systems, including Balboa, based on anomalous conductive response. Our study presents the geophysical survey results at Cobre Panama and is supported by 2D-3D ZTEM and magnetic inversions that appear to validate the survey evidence. 2D synthetic modeling appears to confirm the detectability of the weakly conductive Balboa orebody below 30 m of saprolite cover.
-
-
-
An analysis of MASW responses for urban ground subsidence
More LessAuthors Bitnarae Kim, Myung Jin Nam, Dongwoo Kim, Hannuree Jang and Seo Young SongGround subsidence occurred by sink-hole and cavities in the city has become a serious problem in Korea so these days geophysics survey study for the near-surface become active. Multichannel analysis of surface waves (MASW) is useful survey method to detect anomaly under the subsurface and characterize the structure of S wave velocity of the medium. In this study, we set subsurface seismic model with reference to the real coring data near the sink-hole outbreak point. To define elastic properties of the medium from the core data, we apply ground physics model (GPM) and use Gassman’s equation, the typical model for GPM. For Gassman’s equation input, we consider and assume rock properties as porosity and bulk modulus and so on. The algorithm for this study is time-domain 2D FDM elastic wave modelling algorithm which is based on staggered grid and it damp the reflected wave in the convolutional perfectly matched layer (CPML) zone. We carry out synthetic modelling experiments with the composited medium model and analyse the sensitivity of the underground cavity effect. To set variables of the cavity, we adjust the scale and the depth of the cavity. MASW modelling results show the effect of the cavity and we will study further for the synthetic experiments.
-
Volumes & issues
-
Volume 2019 (2019)
-
Volume 2018 (2018)
-
Volume 2016 (2016)
-
Volume 2015 (2015)
-
Volume 2013 (2013)
-
Volume 2012 (2012)
-
Volume 2010 (2010)
-
Volume 2009 (2009)
-
Volume 2007 (2007)
-
Volume 2006 (2006)
-
Volume 2004 (2004)
-
Volume 2003 (2003)
-
Volume 2001 (2001)
-
Volume 1999 (1999)
-
Volume 1994 (1994)
-
Volume 1987 (1987)
Most Read This Month