- Home
- A-Z Publications
- Geophysical Prospecting
- Previous Issues
- Volume 64, Issue 3, 2016
Geophysical Prospecting - Volume 64, Issue 3, 2016
Volume 64, Issue 3, 2016
-
-
Piecewise‐continuous velocity inversion
Authors Zvi Koren and Igor RavveABSTRACTWe suggest a new method to determine the piecewise‐continuous vertical distribution of instantaneous velocities within sediment layers, using different order time‐domain effective velocities on their top and bottom points. We demonstrate our method using a synthetic model that consists of different compacted sediment layers characterized by monotonously increasing velocity, combined with hard rock layers, such as salt or basalt, characterized by constant fast velocities, and low velocity layers, such as gas pockets. We first show that, by using only the root‐mean‐square velocities and the corresponding vertical travel times (computed from the original instantaneous velocity in depth) as input for a Dix‐type inversion, many different vertical distributions of the instantaneous velocities can be obtained (inverted). Some geological constraints, such as limiting the values of the inverted vertical velocity gradients, should be applied in order to obtain more geologically plausible velocity profiles. In order to limit the non‐uniqueness of the inverted velocities, additional information should be added. We have derived three different inversion solutions that yield the correct instantaneous velocity, avoiding any a priori geological constraints. The additional data at the interface points contain either the average velocities (or depths) or the fourth‐order average velocities, or both. Practically, average velocities can be obtained from nearby wells, whereas the fourth‐order average velocity can be estimated from the quartic moveout term during velocity analysis. Along with the three different types of input, we consider two types of vertical velocity models within each interval: distribution with a constant velocity gradient and an exponential asymptotically bounded velocity model, which is in particular important for modelling thick layers. It has been shown that, in the case of thin intervals, both models lead to similar results. The method allows us to establish the instantaneous velocities at the top and bottom interfaces, where the velocity profile inside the intervals is given by either the linear or the exponential asymptotically bounded velocity models. Since the velocity parameters of each interval are independently inverted, discontinuities of the instantaneous velocity at the interfaces occur naturally. The improved accuracy of the inverted instantaneous velocities is particularly important for accurate time‐to‐depth conversion.
-
-
-
Theoretical analysis of non‐geometric PS‐waves recorded at the water bottom interface
Authors Nihed Allouche and Guy DrijkoningenABSTRACTWhen a seismic source is placed in the water at a height less than a wavelength from the water–solid interface, a prominent S‐wave arrival can be observed. It travels kinematically as if it was excited at the projection point of the source on the interface. This non‐geometric S‐wave has been investigated before, mainly for a free‐surface configuration. However, as was shown in a field experiment, the non‐geometric S‐wave can also be excited at a fluid–solid configuration if the S‐wave speed in the solid is less than the sound speed in the water. The amplitude of this wave exponentially decreases when the source is moved away from the interface revealing its evanescent character in the fluid. In the solid, this particular converted mode is propagating as an ordinary S‐wave and can be transmitted and reflected as such. There is a specific region of horizontal slownesses where this non‐geometric wave exists, depending on the ratio of the S‐wave velocity and the sound speed of water. Only for ratios smaller than 1, this wave appears. Lower ratios result in a wider region of appearance. Due to this property, this particular P‐S converted mode can be identified and filtered from other events in the Radon domain.
-
-
-
Common‐shot Fresnel beam migration based on wave‐field approximation in effective vicinity under complex topographic conditions
Authors Jianping Huang, Jidong Yang, Wenyuan Liao, Xin Wang and Zhenchun LiABSTRACTGaussian beam migration is a versatile imaging method for geologically complex land areas, which overcomes the limitation of Kirchhoff migration in imaging multiple arrivals and has no steep‐dip limits of one‐way wave‐equation migration. However, its imaging accuracy depends on the geometry of Gaussian beam that is determined by the initial parameter of dynamic ray tracing. As a result, its applications in exploration areas with strong variations in topography and near‐surface velocity are limited. Combined with the concept of Fresnel zone and the theory of wave‐field approximation in effective vicinity, we present a more robust common‐shot Fresnel beam imaging method for complex topographic land areas in this paper. Compared with the conventional Gaussian beam migration for irregular topography, our method improves the beam geometry by limiting its effective half‐width with Fresnel zone radius. Moreover, through a quadratic travel‐time correction and an amplitude correction that is based on the wave‐field approximation in effective vicinity, it gives an accurate method for plane‐wave decomposition at complex topography, which produces good imaging results in both shallow and deep zones. Trials of two typical models and its application in field data demonstrated the validity and robustness of our method.
-
-
-
Seismic directional random noise suppression by radial‐trace time–frequency peak filtering using the Hurst exponent statistic
Authors Chao Zhang, Yue Li, Hongbo Lin and Baojun YangABSTRACTRadial‐trace time–frequency peak filtering filters a seismic record along the radial‐trace direction rather than the conventional channel direction. It takes the spatial correlation of the reflected events between adjacent channels into account. Thus, radial‐trace time–frequency peak filtering performs well in denoising and enhancing the continuity of reflected events. However, in the seismic record there is often random noise whose energy is concentrated in certain directions; the noise in these directions is correlative. We refer to this kind of random noise (that is distributed randomly in time but correlative in the space) as directional random noise. Under radial‐trace time–frequency peak filtering, the directional random noise will be treated as signal and enhanced when this noise has same direction as the signal. Therefore, we need to identify the directional random noise before the filtering. In this paper, we test the linearity of signal and directional random noise in time using the Hurst exponent. The time series of signals with high linearity lead to large Hurst exponent value; however, directional random noise is a random series in time without a fixed waveform and thus its linearity is low; therefore, we can differentiate the signal and directional random noise by the Hurst exponent values. The directional random noise can then be suppressed by using a long filtering window length during the radial‐trace time–frequency peak filtering. Synthetic and real data examples show that the proposed method can remove most directional random noise and can effectively recover the reflected events.
-
-
-
Understanding the reverse time migration backscattering: noise or signal?
Authors Esteban Díaz and Paul SavaABSTRACTReverse time migration backscattered events are produced by the cross‐correlation between waves reflected from sharp interfaces (e.g., salt bodies). These events, along with head waves and diving waves, produce the so‐called reverse time migration artefacts, which are visible as low wavenumber energy on migrated images. Commonly, these events are seen as a drawback for the reverse time migration method because they obstruct the image of the geologic structure, which is the real objective for the process. In this paper, we perform numeric and theoretical analysis to understand the reverse time migration backscattering energy in conventional and extended images. We show that the reverse time migration backscattering contains a measure of the synchronization and focusing information between the source and receiver wavefields. We show that this synchronization and focusing information is sensitive to velocity errors; this implies that a correct velocity model produces reverse time migration backscattering with maximum energy. Therefore, before filtering the reverse time migration backscattered energy, we should try to obtain a model that maximizes it.
-
-
-
Optimal implicit staggered‐grid finite‐difference schemes based on the sampling approximation method for seismic modelling
Authors Lei Yang, Hongyong Yan and Hong LiuABSTRACTWe propose new implicit staggered‐grid finite‐difference schemes with optimal coefficients based on the sampling approximation method to improve the numerical solution accuracy for seismic modelling. We first derive the optimized implicit staggered‐grid finite‐difference coefficients of arbitrary even‐order accuracy for the first‐order spatial derivatives using the plane‐wave theory and the direct sampling approximation method. Then, the implicit staggered‐grid finite‐difference coefficients based on sampling approximation, which can widen the range of wavenumber with great accuracy, are used to solve the first‐order spatial derivatives. By comparing the numerical dispersion of the implicit staggered‐grid finite‐difference schemes based on sampling approximation, Taylor series expansion, and least squares, we find that the optimal implicit staggered‐grid finite‐difference scheme based on sampling approximation achieves greater precision than that based on Taylor series expansion over a wider range of wavenumbers, although it has similar accuracy to that based on least squares. Finally, we apply the implicit staggered‐grid finite difference based on sampling approximation to numerical modelling. The modelling results demonstrate that the new optimal method can efficiently suppress numerical dispersion and lead to greater accuracy compared with the implicit staggered‐grid finite difference based on Taylor series expansion. In addition, the results also indicate the computational cost of the implicit staggered‐grid finite difference based on sampling approximation is almost the same as the implicit staggered‐grid finite difference based on Taylor series expansion.
-
-
-
Least‐squares reverse‐time migration in a matrix‐based formulation
Authors Gang Yao and Helmut JakubowiczABSTRACTThis paper describes least‐squares reverse‐time migration. The method provides the exact adjoint operator pair for solving the linear inverse problem, thereby enhancing the convergence of gradient‐based iterative linear inversion methods. In this formulation, modified source wavelets are used to correct the source signature imprint in the predicted data. Moreover, a roughness constraint is applied to stabilise the inversion and reduce high‐wavenumber artefacts. It is also shown that least‐squares migration implicitly applies a deconvolution imaging condition. Three numerical experiments illustrate that this method is able to produce seismic reflectivity images with higher resolution, more accurate amplitudes, and fewer artefacts than conventional reverse‐time migration. The methodology is currently feasible in 2‐D and can naturally be extended to 3‐D when computational resources become more powerful.
-
-
-
Extending illumination using all multiples: application to 3D acquisition geometry analysis
Authors Amarjeet Kumar, Gerrit Blacquière and Eric VerschuurRecent advances in survey design have led to conventional common‐midpoint‐based analysis being replaced by subsurface‐based seismic acquisition analysis, with emphasis on advanced techniques of illumination analysis. Among them is the so‐called focal beam method, which is a wave‐equation‐based seismic illumination analysis method. The objective of the focal beam method is to provide a quantitative insight into the combined influence of acquisition geometry, overburden structure, and migration operators on the resolution and angle‐dependent amplitude fidelity of the image. The method distinguishes between illumination and sensing capability of a particular acquisition geometry by computing the focal source beam and the focal detector beam, respectively. Sensing is related to the detection properties of a detector configuration, whereas illumination is related to the emission properties of a source configuration. The focal source beam analyses the incident wavefield at a specific subsurface grid point from all available sources, whereas the focal detector beam analyses the sensing wavefield reaching at the detector locations from the same subsurface grid point. In the past, this method could only address illumination by primary reflections. In this paper, we will extend the concept of the focal beam method to incorporate the illumination due to the surface and internal multiples. This in fact complies with the trend of including multiples in the imaging process. Multiple reflections can illuminate a target location from other angles compared with primary reflections, resulting in a higher resolution and an improved illumination. We demonstrate how an acquisition‐related footprint can be corrected using both the surface and the internal multiples.
-
-
-
Smart averaging interpolation algorithm comparative test
Authors D. Krilov, Ya. Vaniarho and D. BasaevABSTRACTIn many practical cases, it is necessary to characterize the explored area with a regular set of geodata. Regular matrix data (e.g., ordinary maps) are calculated via existing data interpolation and extrapolation. For low frequency (oversampled) data acquired within a dense profile net (e.g., seismic three‐dimensional structural or gravity mapping), this procedure is mathematically more or less stable and, to a certain extent, unique since we might neglect discrepancies resulting from different interpolations. The situation is quite different for high‐resolution and high‐frequency contaminated data (e.g., raw seismic attributes or geochemistry measurements) represented by sparse profiling. Considering the variety of exploration cases, the investigation of different interpolation algorithm efficiency seems very important. Since it is impossible to compare all algorithms by means of formal mathematics, we have designed a test program. A representative set of seismic attribute maps has been artificially destroyed by introducing blank values (from 20% up to 95%) and then restored by different interpolation algorithms— bicubic, bilinear, nearest neighbor, and “smart averaging.” Smart averaging interpolation is done in a “live” window. The position, form, and size of the window are determined by some mathematical criterion on a trial‐and‐error basis. Discrepancies between restored and initial (true) data have been assessed and analysed. It is shown that the total (absolute) efficiency and comparative (relative) efficiency of the algorithms depend mostly upon the initial interpolant data characteristics. Identifying the best interpolation algorithm for all interpretive cases seems impossible. Some aspects of data processing are discussed in connection with interpolation accuracy.
-
-
-
Implementation aspects of attenuation compensation in reverse‐time migration
By Tieyuan ZhuABSTRACTAttenuation compensation in reverse‐time migration has been shown to improve the resolution of the seismic image. In this paper, three essential aspects of implementing attenuation compensation in reverse‐time migration are studied: the physical justification of attenuation compensation, the choice of imaging condition, and the choice of a low‐pass filter. The physical illustration of attenuation compensation supports the mathematical implementation by reversing the sign of the absorption operator and leaving the sign of the dispersion operator unchanged in the decoupled viscoacoustic wave equation. Further theoretical analysis shows that attenuation compensation in reverse‐time migration using the two imaging conditions (cross‐correlation and source‐normalized cross‐correlation) is able to effectively mitigate attenuation effects. In numerical experiments using a simple‐layered model, the source‐normalized cross‐correlation imaging condition may be preferable based on the criteria of amplitude corrections. The amplitude and phase recovery to some degree depend on the choice of a low‐pass filter. In an application to a realistic Marmousi model with added Q, high‐resolution seismic images with correct amplitude and kinematic phase are obtained by compensating for both absorption and dispersion effects. Compensating for absorption only can amplify the image amplitude but with a shifted phase.
-
-
-
Interaction of helically wound fibre‐optic cables with plane seismic waves
More LessABSTRACTDistributed acoustic sensing is a novel technology for seismic acquisition. In this technology, strain changes induced by seismic waves impinging on an optical fibre are monitored. Due to the fact that glass is relatively rigid, straight glass fibres are not sensitive to broadside waves. We suggest using distributed acoustic sensing systems with fibres helically wound around cables. One increases the fibre sensitivity to broadside waves by decreasing the fibre wrapping angle (the angle between the fibre axis and the plane normal to the cable axis). The optimal wrapping angle is chosen to minimize the impact of Rayleigh waves on the signal measured. This angle depends on the cable Poisson ratio, and it is approximately equal to 30° for cables composed of plastic. For reliable detection of seismic waves, one needs a good mechanical contact between the cable and the surrounding medium. On the other hand, the sensitivity of distributed acoustic sensing systems to primary waves can be significantly reduced if the cable is placed in a cemented borehole.
-
-
-
Using wavelet‐domain adaptive filtering to improve signal‐to‐noise ratio of nuclear magnetic resonance log data from tight gas sands
Authors Ranhong Xie, Youbin Wu, Kang Liu, Mi Liu and Xiangning MengABSTRACTIn tight gas sands, the signal‐to‐noise ratio of nuclear magnetic resonance log data is usually low, which limits the application of nuclear magnetic resonance logs in this type of reservoir. This project uses the method of wavelet‐domain adaptive filtering to denoise the nuclear magnetic resonance log data from tight gas sands. The principles of the maximum correlation coefficient and the minimum root mean square error are used to decide on the optimal basis function for wavelet transformation. The feasibility and the effectiveness of this method are verified by analysing the numerical simulation results and core experimental data. Compared with the wavelet thresholding denoise method, this adaptive filtering method is more effective in noise filtering, which can improve the signal‐to‐noise ratio of nuclear magnetic resonance data and the inversion precision of transverse relaxation time T2 spectrum. The application of this method to nuclear magnetic resonance logs shows that this method not only can improve the accuracy of nuclear magnetic resonance porosity but also can enhance the recognition ability of tight gas sands in nuclear magnetic resonance logs.
-
-
-
Seismoelectric measurements in a porous quartz‐sand sample with anisotropic permeability
Authors Zhenya Zhu, M.N. Toksӧz and Xin ZhanABSTRACTSeismoelectric coupling coefficients are difficult to predict theoretically because they depend on a large numbers of rock properties, including porosity, permeability, tortuosity, etc. The dependence of the coupling coefficient on rock properties such as permeability requires experimental data. In this study, we carry out a set of laboratory measurements to determine the dependence of seismoelectric coupling coefficient on permeability. We use both an artificial porous “sandstone” sample, with cracks, built using quartz‐sand and Berea sandstone samples. The artificial sample is a cube with 39% porosity. Its permeability levels are anisotropic: 14.7 D, 13.8 D, and 8.3 D in the x‐, y‐, and z‐directions, respectively. Seismoelectric measurements are performed in a water tank in the frequency range of 20 kHz–90 kHz. A piezoelectric P‐wave source is used to generate an acoustic wave that propagates through the sample from the three different (x, y, and z) directions. The amplitudes of the seismoelectric signal induced by the acoustic waves vary with the direction. The highest signal is in the direction of the highest permeability, and the lowest signal is in the direction of the lowest permeability. Since the porosity of the sample is constant, the results directly show the dependence of seismoelectric coefficients on permeability. Seismoelectric measurements with natural rocks are performed using Berea sandstone 500 and 100 samples. Because the Berea samples are nearly isotropic in permeability, the amplitudes of the seismoelectric signals induced in the different directions are the same within the measurement error. Because the permeability of Berea 500 is higher than that of Berea 100, the amplitude of the seismoelectric signals induced in Berea 500 is higher than those in Berea 100. To determine the relative contributions of porosity and permeability on seismoelectric conversion, we carried out an analysis, using Pride (1994) formulation and Kozeny–Carman relationship; the normalized amplitudes of seismoelectric coupling coefficients in three directions are calculated and compared with the experimental results. The results show that the seismoelectric conversion is related to permeability in the frequency range of measurements. This is an encouraging result since it opens the possibility of determining the permeability of a formation from seismoelectric measurements.
-
-
-
Permeability dependence of streaming potential coefficient in porous media
Authors Luong Duy Thanh and Rudolf SprikABSTRACTIn theory, the streaming potential coefficient depends not only on the zeta potential but also on the permeability of the rocks that partially determines the surface conductivity of the rocks. However, in practice, it is hard to show the permeability dependence of streaming potential coefficients because of the variation of zeta potential from sample to sample. To study permeability dependence of streaming potential, including the effects of the variation of the zeta potential and surface conductance due to the difference in mineral compositions between samples, we perform measurements on 12 consolidated samples, including natural and artificial samples saturated with 7 different NaCl solutions to determine the streaming potential coefficients. The results have shown that the streaming potential coefficients strongly depend on the permeability of the samples for low fluid conductivity. When the fluid conductivity is larger than than 0.50 S/m for the natural samples or 0.25 S/m for the artificial ceramic samples, the streaming potential coefficient is independent of permeability. This behavior is quantitatively explained by a theoretical model.
-
-
-
Mapping basement structures in the northwestern offshore of Abu Dhabi from high‐resolution aeromagnetic data
Authors Ahmed Salem and Mohammed Y. AliABSTRACTThis paper presents a case study of mapping basement structures in the northwestern offshore of Abu Dhabi using high‐resolution aeromagnetic data. Lineament analysis was carried out on the derivatives of the reduced‐to‐the‐pole magnetic data, along with supporting information from published geologic data. The lineament analysis suggests three well‐defined basement trends in the north–south, northeast–southwest, and northwest–southeast directions. The reduced‐to‐the‐pole magnetic data reveal high positive magnetic anomalies hypothesized to be related to intra‐basement bodies in the deep seated Arabian Shield. Depth to basement was estimated using spectral analysis and Source Parameter Imaging techniques. The spectral analysis suggests that the intruded basement blocks are at the same average depth level (around 8.5 km). The estimated Source Parameter Imaging depths from gridded reduced‐to‐the‐pole data are ranged between 4 km and 12 km with a large depth variation within small distances. These estimated depths prevent a reliable interpretation of the nature of the basement relief. However, low‐pass filtering of the horizontal local wavenumber data across two profiles shows that the basement terrain is characterized by a basin‐like structure trending in the northeast–southwest direction with a maximum depth of 10 km. Two‐dimensional forward magnetic modelling across the two profiles suggests that the high positive magnetic anomalies over the basin could be produced by intrusion of mafic igneous rocks with high susceptibility values (0.008 to 0.016 SI.
-
-
-
Artificial neural networks for removal of couplings in airborne transient electromagnetic data
ABSTRACTModern airborne transient electromagnetic surveys typically produce datasets of thousands of line kilometres, requiring careful data processing in order to extract as much and as reliable information as possible. When surveys are flown in populated areas, data processing becomes particularly time consuming since the acquired data are contaminated by couplings to man‐made conductors (power lines, fences, pipes, etc.). Coupled soundings must be removed from the dataset prior to inversion, and this is a process that is difficult to automate. The signature of couplings can be both subtle and difficult to describe in mathematical terms, rendering removal of couplings mostly an expensive manual task for an experienced geophysicist.
Here, we try to automate the process of removing couplings by means of an artificial neural network. We train an artificial neural network to recognize coupled soundings in manually processed reference data, and we use this network to identify couplings in other data. The approach provides a significant reduction in the time required for data processing since one can directly apply the network to the raw data. We describe the neural network put to use and present the inputs and normalizations required for maximizing its effectiveness. We further demonstrate and assess the training state and performance of the network before finally comparing inversions based on unprocessed data, manually processed data, and artificial neural network automatically processed data. The results show that a well‐trained network can produce high‐quality processing of airborne transient electromagnetic data, which is either ready for inversion or in need of minimal manual processing. We conclude that the use of artificial neural network scan significantly reduce the processing time and its costs by as much as 50%.
-
-
-
Three‐dimensional modelling of magnetotelluric data to image Sehqanat hydrocarbon reservoir in southwestern Iran
Authors I. Mansoori, B. Oskooi, L.B. Pedersen and R. JavaheriABSTRACTA detailed magnetotelluric survey was conducted in 2013 in the Sehqanat oil field, southwestern Iran to map the geoelectrical structures of the sedimentary Zagros zone, particularly the boundary between the Gachsaran Formation acting as cap rock and the Asmari Formation as the reservoir. According to the electrical well logs, a large resistivity contrast exists between the two formations. The Gachsaran Formation is formed by tens to hundreds of metres of evaporites and it is highly conductive (ca. 1 Ωm–10 Ωm), and the Asmari Formation consists of dense carbonates, which are considerably more resistive (more than 100 Ωm). Broadband magnetotelluric data were collected along five southwest–northeast directed parallel lines with more than 600 stations crossing the main geological trend. Although dimensionality and strike analysis of the magnetotelluric transfer functions showed that overall they satisfied local 2D conditions, there were also strong 3D conditions found in some of the sites. Therefore, in order to obtain a more reliable image of the resistivity distribution in the Sehqanat oil field, in addition to standard 2D inversion, we investigated to what extent 3D inversion of the data was feasible and what improvements in the resistivity image could be obtained. The 2D inversion models using the determinant average of the impedance tensor depict the main resistivity structures well, whereas the estimated 3D model shows significantly more details although problems were encountered in fitting the data with the latter. Both approaches resolved the Gachsaran–Asmari transition from high conductivity to moderate conductivity. The well‐known Sehqanat anticline could also be delineated throughout the 2D and 3D resistivity models as a resistive dome‐shaped body in the middle parts of the magnetotelluric profiles.
-
-
-
2.5D direct‐current resistivity forward modelling and inversion by finite‐element–infinite‐element coupled method
Authors Yuan Yuan, Jianke Qiang, Jingtian Tang, Zhengyong Ren and Xiao XiaoABSTRACTTo reduce the numerical errors arising from the improper enforcement of the artificial boundary conditions on the distant surface that encloses the underground part of the subsurface, we present a finite‐element–infinite‐element coupled method to significantly reduce the computation time and memory cost in the 2.5D direct‐current resistivity inversion. We first present the boundary value problem of the secondary potential. Then, a new type of infinite element is analysed and applied to replace the conventionally used mixed boundary condition on the distant boundary. In the internal domain, a standard finite‐element method is used to derive the final system of linear equations. With a novel shape function for infinite elements at the subsurface boundary, the final system matrix is sparse, symmetric, and independent of source electrodes. Through lower upper decomposition, the multi‐pole potentials can be swiftly obtained by simple back‐substitutions. We embed the newly developed forward solution to the inversion procedure. To compute the sensitivity matrix, we adopt the efficient adjoint equation approach to further reduce the computation cost. Finally, several synthetic examples are tested to show the efficiency of inversion.
-
Volumes & issues
-
Volume 72 (2023 - 2024)
-
Volume 71 (2022 - 2023)
-
Volume 70 (2021 - 2022)
-
Volume 69 (2021)
-
Volume 68 (2020)
-
Volume 67 (2019)
-
Volume 66 (2018)
-
Volume 65 (2017)
-
Volume 64 (2015 - 2016)
-
Volume 63 (2015)
-
Volume 62 (2014)
-
Volume 61 (2013)
-
Volume 60 (2012)
-
Volume 59 (2011)
-
Volume 58 (2010)
-
Volume 57 (2009)
-
Volume 56 (2008)
-
Volume 55 (2007)
-
Volume 54 (2006)
-
Volume 53 (2005)
-
Volume 52 (2004)
-
Volume 51 (2003)
-
Volume 50 (2002)
-
Volume 49 (2001)
-
Volume 48 (2000)
-
Volume 47 (1999)
-
Volume 46 (1998)
-
Volume 45 (1997)
-
Volume 44 (1996)
-
Volume 43 (1995)
-
Volume 42 (1994)
-
Volume 41 (1993)
-
Volume 40 (1992)
-
Volume 39 (1991)
-
Volume 38 (1990)
-
Volume 37 (1989)
-
Volume 36 (1988)
-
Volume 35 (1987)
-
Volume 34 (1986)
-
Volume 33 (1985)
-
Volume 32 (1984)
-
Volume 31 (1983)
-
Volume 30 (1982)
-
Volume 29 (1981)
-
Volume 28 (1980)
-
Volume 27 (1979)
-
Volume 26 (1978)
-
Volume 25 (1977)
-
Volume 24 (1976)
-
Volume 23 (1975)
-
Volume 22 (1974)
-
Volume 21 (1973)
-
Volume 20 (1972)
-
Volume 19 (1971)
-
Volume 18 (1970)
-
Volume 17 (1969)
-
Volume 16 (1968)
-
Volume 15 (1967)
-
Volume 14 (1966)
-
Volume 13 (1965)
-
Volume 12 (1964)
-
Volume 11 (1963)
-
Volume 10 (1962)
-
Volume 9 (1961)
-
Volume 8 (1960)
-
Volume 7 (1959)
-
Volume 6 (1958)
-
Volume 5 (1957)
-
Volume 4 (1956)
-
Volume 3 (1955)
-
Volume 2 (1954)
-
Volume 1 (1953)