The Use of the 2D and 3D Complex Wavelet and Ridgelet Transforms in Geophysical Prospecting: Case of Potential Fields Data

Abstract

The complex wavelet and ridgelet transforms are used in the potential field data interpretation for identifying the buried structures responsible for potential field anomalies. Its basis is the use of particular analyzing wavelets belonging to the Poisson semigroup that possess amazing properties regarding potential fields. In fact, the analyzed anomaly displays a conical signature in the wavelet domain and whose apex is pointing out at the causative homogeneous structure. Fundamentally, the interpretation is performed in the upward-continued domain where, the dilation of the wavelet transform is the upward-continuation altitude. This confers on the wavelet transform a considerable advantage: its robustness with respect to noise. The method is also developed to identify the depth, horizontal positions, size, strike direction, dips and shape of elongated 3D structures such as finite-dimensional dykes and faults. For this type of anomaly, the 2D wavelet transform corresponds to the ridgelet transform performed in the Radon domain, where elongated anomalies are recognized by high amplitude signatures. A reminder of the developed theory and applications in the 2D and 3D cases on real case studies are shown.

1. Introduction

One of the most important applications in geophysical prospecting is the identification, localization and characterization of bodies of geological interest, especially the causative sources of potential field anomalies (gravitational, magnetic, electrical, thermal and so on) measured at the surface or since high-resolution marine and airborne surveys are carried out. Thus, it continues to stimulate many important methodological developments in analysis and inversion techniques [1]. The purpose of inversion methods is to recover the source distribution using an integral equation that relates the measured potential field and the causative source distribution [2, 3, 4, 5]. These analysis techniques contribute to reducing the non-uniqueness of the inverse problem by adding some a priori geological assumptions [6, 7]. Another category of processing methods that are not part of the inverse methods and which does not seek after the distribution of the source but can give information on the depth of the top of the causative sources [8, 9] or data transformation techniques such as downward and upward continuations, horizontal and vertical derivatives or reduction to the pole and oblique derivatives which produce the transformed fields, where the features of the original field are enhanced, using the relationship between the measured field and the distribution of its sources, based on a sum of convolution products using several transformation operators and an appropriate Green function [10, 11, 12, 13, 14, 15]. However, these methods do not allow localizing sources along the z dimension, so another analyzing method based on the shape of the anomalies is the Euler deconvolution [16], which needs to be improved to eliminate noise effects since it is based on the calculation of gradients [17].

The continuous wavelet transform is an approach that can make easier the analysis of large amounts of data [18, 19, 20, 21]. This method utilizes the homogeneity properties of the potential field, to identify and localize the causative sources [22, 23, 24, 25]. Further works show the robustness of this method with respect to noise [26, 27], as revealed by many applications in geophysical prospecting, such as in aeromagnetics data [28, 29, 30, 31], spontaneous electrical potential [32, 33, 34]; gravity data [35, 36, 37] and electromagnetic data [38]. The 2D wavelet method is then developed [30, 39] in order to localize and identify the potential fields anomalies causative structures in the case of elongated structures such as dykes, faults, etc. Thus, in the present work, after a brief recall of the theory in the cases 2D and 3D, an application to a part of aeromagnetic survey in the NW of Algeria will be presented and discussed.

2. The wavelet theory

The application of wavelets is recent in signal and image processing, but its mathematical history is much older since it basically works linking the Littlewood-Paley decomposition (1930), the version given by J.O. Strömberg from the basis of Franklin (1927) and the identity of Calderön (1960) [40].

Wavelet theory has experienced great development since the 1980s, to name a few references such as the work of Grossmann and Morlet (1984), [41, 42]. An interesting compilation of work on the developments and applications of wavelet theory for non-stationary signals in geophysics can be found in [19, 43, 44].

Wavelet analysis methods are essentially based on a representation of signals at different scales [45]. This is very interesting in geophysics since the information carried in the signals is carried by scaling laws or by non-stationaries. This is the case of seismological signals or potential fields whose variations represent the effects of multiscale or even fractal bodies of geological interest.

Other important developments have been made in wavelet theory since the work of Grossmann and Morlet [46], such as orthogonal wavelets, multiresolution analysis as well as the development of fast numerical algorithms [47, 48, 49, 50]. The orthogonal wavelets are used to develop discrete wavelet transforms, unlike the continuous wavelet transform. An overview of how the discrete wavelet transforms can be used in the analysis of geophysical time series can be found in [51]. For instance, the continuous wavelet transform is used in the interpretation of potential field theory, and a review can be found in [28].

3. The continuous wavelet and ridgelet transforms

Here is a brief recall of the mathematical framework used in this work. The continuous wavelet transform [46] makes it possible to process large amounts of data. The wavelet transform allows for the detecting and characterizing homogeneous singularities in signals using the property of the homogeneity degree of the analyzed function [52, 53, 54].

3.1 The 1D continuous wavelet analysis of potential fields data

The wavelet transform $w$ which is a convolution product of an analyzing wavelet $g$ and a function ${\phi }_0$, is defined as follows:

$$w\left[g,{\phi }_0\right]\left(b,a\right)\equiv \underset{\Re }{\int }\frac{1}{a}g\left(\frac{b-x}{a}\right){\phi }_0\left(x\right)\mathit{dx}=\left(D_{a}g\ast {\phi }_0\right)\left(b\right)$$E1

$$D_{a}g\left(x\right)\equiv \frac{1}{a}g\left(\frac{x}{a}\right)$$E2

where $D_a$ is the dilation operator, the dilation parameter $a>0$, and $b$ is the translation parameter.

To be admissible, the analyzing wavelet $g\left(x\right)$ should have a zero-mean oscillating behavior localized in a finite interval including the origin [19], which enables the wavelet transform to perform a local analysis of the signal [38]. In the case of complex continuous wavelet transform, this analyzing wavelet must be an oscillating complex function, localized on the real line [39]. The dilated wavelets $D_{a}g\left(x\right)$ are $a$ then band-pass filters with bandwidth proportional to dilation $a$. This oscillating property means that the wavelet has a vanishing integral and allows the reconstruction of the analyzed signal from its wavelet transform [28, 38, 39]. Also, the covariance of the wavelet transform is an interesting property with respect to homogeneous functions [39], where the geometrical meaning is that the wavelet of homogeneous singularity displays a cone-like appearance whose apex points onto the singularity for $a\downarrow 0+$, so according to potential field theory, the analyzed field is produced by a homogeneous source located at $\left(\mathit{xs},\mathit{zs}\right)$ and, the dilation may correspond to an upward-continuation offset of the analyzed potential fields [24, 25, 26, 27, 39]. So, the cone-like structure apex is then located at $\left(\mathit{xs},a=-\mathit{zs}\right)$ below the positive half-plane of the wavelet.

3.2 The 2D continuous ridgelet transform of potential fields data

The 1D continuous wavelet transform method is applied to analyze 2D potential field anomalies data, measured in the horizontal plane by generalizing the Eq. (2) and then obtain the ridgelet transform [55, 56],

$$\begin{array}{c}R\left[r,{\phi }_0\right]\left(b,a,s_{\perp }\right)=\underset{{\Re }^2}{\int }\frac{1}{a}r\left(\frac{b-x}{a},y,s_{\perp }\right){\phi }_0\left(x,y\right)\mathit{\text{dxdy}}=\underset{\Re }{\int }\frac{\mathit{dx}}{a}g\left(\frac{b-x}{a}\right)\underset{\Re }{\int }{\phi }_0\left(x,y\right)\mathit{\text{dxdy}}\\ =W\left[g,\mathit{RT}\left[{\phi }_0,s_{//}\right]\right]\left(b,a\right)\end{array}$$E3

where $s_{\perp }$ a unit vector is perpendicular to the anomaly strike, $s_{//}$ is a unit vector in the direction of the elongated anomaly. The analyzing ridgelet is obtained by steering a 1D Poisson wavelet $g\left(x\right)$ in the perpendicular direction $y$, $\mathit{RT}$ is the Radon transform of the potential field anomaly. The ridgelet transform of 2D anomalies, as defined in Eq. (3) is obtained by computing the wavelet transform for each value of the angular parameter in the Radon domain [28, 38, 39]. We use the complex ridgelet transform since we use complex analyzing wavelets [39]. The complex wavelet transform allows to analyze a signal using modulus and phase. The phase of the complex wavelet transform provides information about the inclination of the source. The imaginary part is obtained from the Hilbert transform of the real part. In this case, the complex wavelets correspond to analytical signal and their modulus and phase can be determined [27, 38].

Both the modulus and the phase of the ridgelet transform are used to localize the sources and, display conspicuous cone-like patterns associated with each analyzed anomaly. This cone is pointing to the source depth $\mathit{zs}$ [39]. These locations are obtained by testing the geometry of the cone pointing towards the depth of the source for each grid point, in the half-space $\left(\mathit{x},\mathit{z}\right)$ [38, 39]. Subsequently, a statistical method is introduced, and the likelihood for the occurrence of an apex at source location $\left({\mathit{x}}_{\mathit{s}},{\mathit{z}}_{\mathit{s}}\right)$ is evaluated by the maximum entropy criteria $\mathit{\rho }$ [38, 39, 57]:

$$\rho \left(x_s,z_s\right)=\frac{lnN+{\sum }_{i=1}^N{h}_{i}ln{h}_i}{lnN}$$E4

N corresponds to the number of grid points; hi correspond to the histogram values of the slopes along the modulus or phase lines forming the cone; $\mathit{\rho }$ takes its values in the range [0, 1]. The result is a tomography map of the sources.

4. Application to aeromagnetic data

In this section, we show some results obtained by applying the complex wavelet and ridgelet analysis, to identify and localize structures responsible for aeromagnetic field anomalies, in the seismogenic Cheliff basin (NW of Algeria). The geological and seismotectonics framework can be found in [58, 59]. This region is one of the most seismically active zones of the western Mediterranean Sea [58, 60], related to the collision between African and European plates since the Upper Cretaceous [61]. The kinematic models derived from the Atlantic Ocean magnetic anomalies study have shown that this convergence is linked to a counter-clockwise rotation of Africa relative to Eurasia [61, 62]. The seismic activity in this area is directly associated with the plate boundary between Europe and Africa. This region is known for having been the site of two destructive earthquakes: 9 September 1954, earthquake of Orleansville with a magnitude of 6.5 and 10 October 1980, El Asnam one, with a 7.3 magnitude.

The aeromagnetic data used in this work resulted from the digitization of aeromagnetic maps issued from the aeromagnetic survey of Algeria. This survey was carried out by Aero Service Corporation between 1970 and 1974. The maps, digitized by the shape recognition method, were interpolated at the nodes of a 325 × 325 regular grid and reduced to a pole (Figure 1). This map displays important magnetic anomalies in North, along the coast, in the Mediterranean Sea and in the South. The ridgelet analysis results are shown in Figure 2. The 3D image obtained attests to the geological and seismotectonic complexity of the area. The elongated structures identified are a juxtaposition of prismatic bodies at different depths. The structures bordering the Chelliff basin are elongated in the E-W direction. In the North, in the Mediterranean Sea, their depth reaches 31 km. At the coast, within the volcanic structures, these structures reach a depth of 20 km. To the south of the basin, at the Ouarsenis Mountains, the depth of these structures reaches 29 km. In the Cheliff basin, the structures are oriented in NE-SW, NW-SE and E-W directions and located between near-surface and 25 km, while those oriented in N-S direction are at depth ranging between 9 and 16 km. Elongated structures oriented N-S appear to the North, limiting offshore and coastal anomalies. These structures reach depths of 20 km.

Figure 1.

Aeromagnetic map of Chlef region: total field anomaly reduced to the pole (Modified from [31]).

Figure 2.

3D imaging of magnetized structures identified from the complex ridgelet transform. The North direction is given by the latitude axis. The color scale corresponds to the maximum entropy criteria to select the source location (modified from [31]).

In order to skecth out the topography of the magnetized substratum and identify the structures in depth, the complex wavelet transform is applied to a N-S magnetic anomalies profile, located at 1°40′ and crossing magnetic anomalies from Mediterranean Sea to the North, as far as the Ouarsenis Mountains in South. Figure 3 shows the intensity of magnetic anomaly (top of figure) varying from −40nT to +40nT. The middle of the figure corresponds to the modulus of the wavelet transform and the 2D image (bottom), where the magnetized substratum top depth, identified by the maximum entropy criteria, ranges between 6 and 30 km. The deepest is located in North (Med. Sea) and South (Ouarsenis Mounts); less deep in the sedimentary basin (thickness of basin). Many faults and contacts are identified along this profile. We can cite the Oued Fodda region (site of 1980 earthquake), where the O. Fodda fault is identified at the latitude of 4010 km at a depth of 6 km, with an inclination of 30° identified from the phase of the wavelet (Figure 4), these results are confirmed by the results of dislocation model of vertical movements for this area [59]. A good correlation is shown with N-S geological profile (Figure 5).

Figure 3.

The shape of the magnetized substratum and the identified structures in-depth, along a N-S profile (1°40′E) (modified from [31]).

Figure 4.

Phase map of the complex continuous wavelet transform (modified from [31]).

Figure 5.

Correlation with geological N-S profile (modified from [31]).

5. Conclusion

The use of the wavelet and ridgelet transform in geophysical exploration allows for the identification of sources of potential field anomalies in 2D and 3D cases. The wavelet transform possesses many interesting mathematical properties with respect to potential fields theory; studies show that when applied to potential fields, it can have a deep physical sense, since the idea in the use of a homogeneous source is that an elongated geological structure may be replaced by a small number of equivalent point sources. The homogeneity of these point sources depends on the shape of the geological structure. Also, in the range of dilations, the signal-to-noise ratio is much better and makes the method more robust.

The ridgelet transform helps in the automatic detection of elongated structures in 3D, and the information provided by the wavelet transform concerning the identified sources, such as dip, depth, and dimensions, can be used to reduce the non-uniqueness of the inverse problem considerably. The application of real aeromagnetic data, without any a priori assumptions, shows a good correlation with known geological structures and identifies many more unknown structures.

Acknowledgments

This work is supported by the Centre de Recherche en Astronomie, Astrophysique & Géophysique. Our gratitude to D. Gibert, F. Moreau, P.G. Saracco, Sailhac and H. Beldjoudi. We are grateful to the editor and reviewers and special thanks to L. Divic and Tonči Lučić.

Conflict of interest

The author declares no conflict of interest.