Gravity Anomaly Interpretation Using the R-Parameter Imaging Technique over a Salt Dome

Abstract

Rapid imaging technique, so-called “R-parameter”, utilized for interpreting a gravity anomaly profile. The R-parameter based on calculating the correlation factor between the analytic signal of the real anomaly and the analytic signal of the forward anomaly of assumed buried source denoted by simple geometric shapes. The model parameters (amplitude, origin, depth, and shape factor) picked at the maximum value of the R-parameter. The technique has been proved on noise free and noisy numerical example, numerical example showing the impact of interfering sources. Furthermore, the introduced technique has been successfully applied to visualize a salt dome gravity anomaly profile, USA. The obtained results are in good agreement with those reported in the published studies and that with that obtained from drilling.

1. Introduction

Gravity data interpretation has been widely used to appraise the different types of subsurface structures and their locations [1, 2, 3, 4, 5, 6, 7, 8]. Gravity methods have been widely applied to ore and mineral exploration [9, 10, 11, 12, 13], hydrocarbon exploration [14, 15, 16], cave detection [17, 18], hydrogeology [19, 20], geothermal and volcanic activity [21, 22, 23], locating of unexploded military ordnance [24], environmental and engineering application [25, 26] and archaeological investigations [27, 28].

The quantitative interpretation of gravity data using simple models (spheres and cylinders) is common in exploratory geophysics and continues to be of interest [29, 30, 31, 32, 33]. In geologic contexts with a single gravity anomaly, it can be quite appropriate [34]. A single isolated causal body can invert this recorded gravity anomaly to establish its distinctive inverted parameters and fit the recorded data.

The simple geometric models can be matched with the subsurface structures encountered during application of several approaches for inversion [35, 36, 37, 38, 39]. These methods include graphical and numerical characteristic points approaches [40, 41, 42], ratio technique [43], Fourier transform method [44], the neural network algorithms [45], Mellin transform technique [36], and Werner deconvolution technique [46]. However, the drawbacks of these methods based on tending to generate high number of invalid solution due to few numbers of points and data used, noise or window size incompatibility. As a result, these approaches are subjective, which can lead to significant inaccuracies in calculating the buried anomalous body's characteristic inverse parameters [41, 47], which is to be expected. Gupta [48] and Essa [49] developed techniques depending on successive minimization approaches, which utilize the whole measured data to assess the depth parameter and then used some of characteristic points to continue in estimating the rest parameters such as amplitude coefficient. Shaw and Agarwal [37] used the Walsh transform scheme to determine the depth of buried bodies. Mehanee [47] used the regularized conjugate gradient method to construct an effective iterative method based on the use of logarithms of the model parameters for gravity inversion. The method inverts the residual gravity data acquired along profile for evaluating a depth and amplitude coefficient of buried bodies and suitable for subsurface imaging and mineral exploration.

Here, the study proposed an application of the robust R-parameter imaging method to interpret residual gravity data along a profile over idealized geometric bodies such as semi-infinite vertical cylinder, infinitely long horizontal cylinder, and sphere models. The goal is to establish the underlying approximative model by determining the body parameters, which include its origin, depth, amplitude coefficient, and shape. The R-parameter imaging method depends on the correlation coefficient amongst the analytic signal of the collected and calculated gravity data. The optimum solution occurs at the maximum R-parameter value.

The benefit behind the use of this method is fall in estimating the depth and body location with an acceptable value compared to the true ones and used the whole gravity data points of the profile, instead of just a few characteristic points. In addition to the method does not require priori information of the subsurface and directly interpret the anomaly from the given observed data. This chapter begins with a layout of the forward modeling, which contains a theoretical gravity formula, an R-parameter imaging approach description, numerical models test without and including noise, and a field data for slat dome investigation.

2. Forward modeling

A closed-form solution for the gravity anomaly caused by simple geometric structures at a measured point (x_j) along a profile (Figure 1) is given by [29, 50, 51].

Figure 1.

Geometry and parameters of sphere (top), semi-infinite vertical cylinder (middle and infinitely long horizontal cylinder (bottom) (re-drawn from [39]).

$$\mathrm{g}\left(x_j,x_o,z,z_o,\eta ,q,A\right)=A\frac{{\left(z_o-z\right)}^{\eta }}{{\left[{\left(x_j-x_o\right)}^2+{\left(z_o-z\right)}^2\right]}^q},j=1,2,3,\dots ,n$$E1

where x_j and x_o are the X-coordinates of the measured points and the buried body center, z and z_o are the Z-coordinates of the measured point and the subsurface source (the z axis is chosen positive downward) (Figure 1), q is the shape, and A is the amplitude coefficient. The appropriate explanations of A, q, and η for the above-mentioned simple bodies are shown in Table 1.

Case	A	η	q
Sphere	$\frac{4}{3}\mathit{\pi \gamma \sigma }{r}^3$	1	3/2
Horizontal cylinder	2πγσr^2	1	1
Vertical cylinder	πγσr^2	0	1/2

Table 1.

Definitions of A, η and q of the simple-geometric bodies.

3. Methodology

The gravitational anomaly’s analytic signal is written as follows [52, 53]:

$${\xi }_s\left(x\right)=\frac{\partial g}{\partial x}-i\frac{\partial g}{\partial z},i=\sqrt{-1},$$E2

where $\frac{\partial g}{\partial x}$ and $\frac{\partial g}{\partial z}$ are the horizontal and vertical gradients of the gravity anomaly.

The amplitude of the analytic signal $\left|{\xi }_s\left(x\right)\right|$ can be calculated as follows [44]:

$$\left|{\xi }_s\left(x\right)\right|=\sqrt{{\left(\frac{\partial g}{\partial x}\right)}^2+{\left(\frac{\partial g}{\partial z}\right)}^2}.$$E3

By an adapting the horizontal and vertical derivatives to Eq. (1), and putting the obtained outcomes into Eq. (3), we get the following:

$$\left|{\xi }_s\left(x\right)\right|=A{\left(z_o-z\right)}^{\left(\eta -1\right)}\frac{\sqrt{4q^2{\left(x_j-x_o\right)}^2{\left(z_o-z\right)}^2+{\left(\eta {\left(x_j-x_o\right)}^2-{\left(z_o-z\right)}^2\left(2q-\eta \right)\right)}^2}}{{\left({\left(x_j-x_o\right)}^2+{\left(z_o-z\right)}^2\right)}^{\left(q+1\right)}},$$E4

where j = 1, 2, 3, …., n. To calculate the horizontal location (x_o) and depth (z_o) of the buried target (Figure 1), the 2-D X-Z mosaic of the correlation coefficient (R) is constructed from the analytic signal amplitudes $\left|{\xi }_{\mathit{so}}\left(x\right)\right|$ of the measured data and $\left|{\xi }_{\mathit{st}}\left(x\right)\right|$ calculated from the theoretical generated data by a supposed simple-geometric source S(x_o, z_o), and is expressed as:

$$R\left(x_o,z_o\right)=\frac{{\sum }_j^n{\left|{\xi }_{\mathit{so}}\left(x\right)\right|}_j{\left|{\xi }_{\mathit{st}}\left(x\right)\right|}_j}{\sqrt{{\sum }_j^n{{\left|{\xi }_{\mathit{so}}\left(x\right)\right|}_j}^2{\sum }_j^n{{\left|{\xi }_{\mathit{st}}\left(x\right)\right|}_j}^2}}.$$E5

The analytic signal $\left|{\xi }_{\mathit{so}}\left(x\right)\right|$ is assessed numerically along the profile using Eq. (3). To map the relevant discrepancy of the R-parameter from which x_o and z_o are appraised, discretization in the X- and Z-directions is done around the anticipated spatial location of the supposed source. The R-parameters value (R-max) reaches the maximum when the depth and location of the assumed source match the true ones.

4. Numerical examples

To verify constancy in performance of the proposed method, numerical example without noise (noise-free) and with a 20% random noise (noisy) is tested. Another numerical example to evaluate the accuracy and stability in assessing the model parameters in case of interference/neighboring influence.

4.1 Model 1

4.1.1 Noise-free data

The R-parameter imaging method is applied to noise-free numerical gravity anomaly due to simple model consisting of a a horizontal cylinder model (q = 1) with z_o = 5 m, A = 100 mGal m, and x_o = 51 m along a 101-m profile length (Figure 2a). This anomaly is the observable (measured) data that needs to be interpreted. The suggested methodology was initiated with estimating the horizontal and vertical derivative anomalies of the residual anomaly in Figure 2a (Figure 2b), and then calculate the amplitude of the analytic signal (Figure 2c). The 2-D mosaic surface S is constructed with depending on calculating the R-parameter values and discretized into 1-m intervals in X- and Z-directions and covers an area of 101 × 11 m in the X- and Z directions. Based on a priori information, the range of the model parameters was selected. The R-parameter values were obtained from (Eq. (5)) by setting the shape to its true value (i.e., q = 1) and utilizing the abovementioned ranges of the model location parameters x_o and z_o. The R-max value is represented by a black dot and equal 1 (Figure 2d). This is demonstrating that the methodology fruitfully recovered the real values of the location of the inferred gravity profile’s origin point (x_o = 51 m) and the depth of the target (z_o = 5 m).

Figure 2.

Model 1: Noise-free data. (a) Horizontal cylinder gravity anomaly, (b) Horizontal and vertical gradients of (a), (c) Analytic signal anomaly using the data of (b), and (d) 2-D mosaic of the R-parameter and the R-max value.

Table 2 shows a different shape values that employed in the interpretation process. The results (Figure 3 and Table 2) reveal that at q = 1, the R-max = 1 and indicating that the method is stable and capable of capturing thel exact values of the model parameters.

Shape factor	Maximum R-parameter
(q)	(R-max )
0.5	0.6696
0.6	0.7768
0.7	0.8897
0.8	0.9620
0.9	0.9930
1	1.0000
1.1	0.9962
1.2	0.9883
1.3	0.9791
1.4	0.9700
1.5	0.9611

Table 2.

Model 1: Noise-free data. The R-parameter computed for the different shape factors.

Figure 3.

Model 1: Noise-free data. The R-parameter, depth and shape factor relationship.

We applied the same procedures (by utilizing Eq. (4) as the forward modelling formula in this case) to the analytic signal data presented in Figure 2c to explore the recital of the current scheme when used to the analytic signal data themselves instead of the residual gravity data. Figure 4 shows the outcomes, which are match with those derived from the above-mentioned elucidation of gravity data (Table 3).

Figure 4.

Model 1: Noise-free data. (a) Analytic signal anomaly (Figure 2c) subjected to the same interpretation, (b) Horizontal and vertical gradients of (a), (c) Analytic signal anomaly using the data of (b), and (d) 2-D mosaic of the R-parameter and the R-max value.

Estimated model parameters	Analytic signal data	Gravity anomaly data
A (mGal m)	100	100
zo (m)	5	5
q	1	1
xo (m)	51	51

Table 3.

Model 1: Noise-free data. Comparison between the model parameters estimated from the interpretation of using residual anomaly and analytic signal anomaly.

4.1.2 Noisy data

Given the lack of totally noise-free gravity field data, a 20% random noise (Figure 5a) has been introduced to the data in Figure 2a. The horizontal and vertical derivatives, besides the magnitude of the analytic signal of the measured gravity anomaly, are depicted in Figure 5b and c. The R-parameter values were evaluated utilizing Eq. (5) and created a 2-D mosaic surface (Figure 5d). The maximum R-parameter value is 0.94. The imaging-derived model parameters (z_o = 6 m and x_o = 51 m) are quite near to the real values, indicating that the established technique is stable.

Figure 5.

Model 1: noisy data. (a) Noisy gravity anomaly of Figure 2a after adding 20% noise level, (b) Horizontal and vertical gradients of (a), (c) Analytic signal anomaly using the data of (b), and (d) 2-D mosaic of the R-parameter and the R-max value.

The amplitude coefficient (A) is substantially overstated (Figure 6a–d) when the analytic signal data (Figure 5c) is interpreted. This is unsurprising given that the technique tries to fit the real data, and the anomalous body's inferred depth (z) is likewise overstated (Figure 6a). Table 4 shows a comparison of model parameters derived by the established technique from the elucidation of each analytic signal and residual anomaly data. This investigation shows that using the given technique to analyses gravity data produces more precise findings than using analytic signal data.

Figure 6.

Model 1: Noisy data. (a) Analytic signal anomaly (Figure 5c) subjected to the same interpretation, (b) Horizontal and vertical gradients of (a), (c) Analytic signal anomaly using the data of (b), and (d) 2-D mosaic of the R-parameter and the R-max value.

Estimated model parameters	Analytic signal data	Gravity anomaly data
A (mGal m)	212.90	100
zo (m)	7	5
q	1	1
xo (m)	51	51

Table 4.

Model 1: Noisy data. Comparison between the model parameters estimated from the interpretation of using residual anomaly and analytic signal anomaly.

4.2 Model 2: neighboring effect

The performance of the proposed inversion method with complicated field anomalies and the effect of interfering subsurface structures was investigated. To achieve this, we once again generate a synthetic model data from multiple source bodies as a horizontal cylinder model with A₁ = 100 mGal m, z₁ = 3 m, x₁ = 30 m, and q₁ = 1) and a sphere model with A₂ = 400 mGal m², z₂ = 4 m, x₂ = 80 m, and q₂ = 1.5) (Figure 7a).

Figure 7.

Model 2: Interference/neighboring effect. (a) Gravity anomaly generated by two different adjacent bodies, (b) Horizontal and vertical gradients of (a), (c) Analytic signal anomaly using the data of (b), and (d) 2-D mosaic of the R-parameter and the R-max values.

Figure 7b and c illustrate the horizontal and vertical gradients of the composite gravity anomaly, as well as the amplitude of the analytic signal. The R-parameter values were determined using Eq. (5) for each source location and a 2-D mosaic surface S of 101 × 11 m in the X- and Z-directions constructed and discretized into 1-m intervals in both directions. The 2-D mosaic (Figure 7d) indicates that the R-max value for each source is 0.8 and 0.62 at q equal 1 and 1.5, respectively. The obtained model parameters for horizontal cylinder and sphere are A₁ = 120.1 mGal m, z₁ = 3.6 m, and x₁ = 30 m and A₂ = 443.1 mGal m², z₂ = 4.2 m, and x₂ = 80 m, respectively, which the results shows that the method is stable and robust.

To better understand the procedure, we tainted the composite anomaly (Figure 7a) with a 20% noise level (Figure 8a). The horizontal and vertical derivatives, as well as the corresponding amplitude of the analytic signal, are shown in Figure 8b and c. The retrieved R-parameter image is shown in Figure 8d, with R-max values of 0.79 and 0.61 for a horizontal cylinder and a sphere, respectively. The drop in maximum parameter values compared to (Figure 7d) is attributable to the noise introduced into the data as well as the effect of the nearby objects. The model parameters for the first and second bodies revealed by imaging are: A₁ = 126.9 mGal m, z₁ = 3.8 m, and x₁ = 30 m and A₂ = 556.3 mGal m², z₂ = 4.7 m, and x₂ = 80 m, respectively, which are quite adjacent to the real values.

Figure 8.

Model 2: Interference/neighboring effect with noise. (a) Noisy gravity anomaly of Figure 7a after adding 20% noise level, (b) Horizontal and vertical gradients of (a), (c) Analytic signal anomaly using the data of (b), and (d) 2-D mosaic of the R-parameter and the R-max values.

Figure 9 depicts the results of the study of the noisy analytic signal data seen in Figure 8c. The amplitude coefficients and burial depths recovered from the elucidation are exaggerated (Figure 9a–d), as shown in Table 5, which coincides with and confirms the aforementioned results.

Figure 9.

Model 2: Interference/neighboring effect with noise. (a) Analytic signal data (Figure 8c) subjected to the same interpretation, (b) Horizontal and vertical gradients of (a), (c) Analytic signal anomaly using the data of (b), and (d) 2-D mosaic of the R-parameter and the R-max values.

Estimated model	Noisy contaminated interference/neighboring effect
	Analytic signal data		Gravity anomaly data
Parameters	First anomaly	Second anomaly	First anomaly	Second anomaly
A (mGal m^2q−η)	174.2 mGal m	1554.3 mGal m^2	126.9 mGal m	556.3 mGal m^2
zo (m)	4.2	7.6	3.8	4.7
q	1	1.5	1	1.5
xo (m)	30	80	30	80

Table 5.

Model 2: Interference/neighboring effect with noise. Comparison between the model parameters estimated from the interpretation of using residual anomaly and analytic signal anomaly.

On the basis of the theoretical models presented above, it can be inferred that the technique described here is stable and robust.

5. Field data

A published field example over a salt dome anomaly is examined in order to thoroughly test the applicability of the established methodology. For a variety of reasons, this case was chosen. First, the residual gravity profile was created by a simple body that may be truthfully inferred. Second, the drilling information helps in estimating the density contrast of the underlying body. Knowing the density contrast, the radius can be calculated and the depth to the top also can be inferred by using the definition of the amplitude coefficient (Table 1). Moreover, the depth of the vertical cylinder model is measured to the top but the depth of a horizontal cylinder and sphere model is measured to the center of the body (Figure 1). Third, the gravity data was taken from an area with recognized drilling information, allowing the results obtained from the technique proposed here to be cross-validated against those received via drilling.

5.1 Humble dome anomaly, Houston, USA

Gravity map was acquired over the Humble Dome, Houston, Texas. At the earth's surface, the measurements of this salt dome structure reveal a negative circular contoured Bouguer anomaly ([41], Figures 8–16). A Bouguer gravity profile is taken across the center of the Humble salt dome gravity map in Houston ([41], Figures 8–16). The Bouguer gravity profile was subject to a suitable separation method to remove the regional anomaly and obtain the residual gravity anomaly. The residual gravity anomaly profile of about 26 km long was digitized at an interval of 0.26 km (Figure 10a).

Figure 10.

The Humble dome anomaly, USA. (a) Gravity anomaly profile (red dotted lines) and the optimum-fitting model (solid black line), (b) Horizontal and vertical gradients of (a), (c) Analytic signal anomaly using the data of (b), and (d) 2-D mosaic of the R-parameter and the R-max value.

The R-parameter method procedures were applied to the residual gravity anomaly profile for the available shape parameters (Table 6). Figure 10b–d express the horizontal and vertical derivative anomalies, the amplitude analytic signal, and the R-parameter 2-D mosaic. It is found that the R-max value is 0.99 corresponds to a spherical shape (q = 1.5). The best estimated model parameters are: q = 1.5, z_o = 4.90 m, x_o = 0 m, and A = −269.39 mGal m².

q	R-max
0.5	0.999462
0.6	0.997341
0.7	0.996274
0.8	0.992509
0.9	0.986895
1	0.984110
1.1	0.987549
1.2	0.992912
1.3	0.996863
1.4	0.998792
1.5	0.999501

Table 6.

The Humble dome gravity anomaly, USA. The R-parameter calculated from different shape factors.

The humble dome anomaly has been interpreted by several authors assuming a spherical source to decide the depth of the salt body. The obtained results agree well with those depths to the center that obtained by the published literatures of [36, 41, 47, 54, 55, 56] (Table 7).

Model parameters	Approaches and techniques of					Present study
	[28]	[33]	[39]	[46]	[47]
A mGal km^2	—	—	−292.54	—	−279.81	−269.39
zo (km)	4.96	4.97	4.62	4.81	4.58	4.90
q	1.5	1.5	1.5	1.5	1.5	1.5
xo (km)	—	—	—	—	—	0

Table 7.

The Humble dome gravity anomaly, USA. The Estimated parameters.

By using a density contrast of −0.13 gm/cm³ of [41], then the depth to the top of the spherical body of the humble dome obtained from the proposed technique is 315 m, which in excellent covenant with the true depth (305 m) confirmed by drilling and seismic information [41]. Table 8 shows that several other researchers utilizing the same density contrast found some differences in the depths to the top of this spherical source. The use of simple geometrically bodies in the constrained class of spheres, horizontal cylinders, and vertical cylinders is thus suggested as a way to accurately apply the current methodology to extract depth information. As a result, if exact density contrasts are used, the related radii can be correctly computed as well.

Approaches and techniques of	Depth to the top (m)
[33]	426
[39]	299
[46]	326
[47]	326
Present study	315

Table 8.

Depth to the top of the spherical body of the Humble dome anomaly, obtained by various approaches using an assumed density contrast (Δρ) of −0.13 gm/cm³ [41].

5.2 Louisiana dome anomaly, USA

A residual gravity map was acquired over a salt dome off the coast of Louisiana, USA ([41], Figures 8–20). The residual gravity anomaly profile [57] is redrawn across the center of the map, normal to the causal anomaly’s striking. The residual gravity anomaly profile of about 13,000 m long was digitized at sampling interval of 200 m (Figure 11a).

Figure 11.

The Louisiana dome anomaly, USA. (a) Gravity anomaly profile (red dotted lines) and the optimum-fitting model (solid black line), (b) Horizontal and vertical gradients of (a), (c) Analytic signal anomaly using the data of (b), and (d) 2-D mosaic of the R-parameter and the R-max value.

By applying the R-parameter method procedures mentioned before to the residual gravity anomaly profile of Louisiana we get the available shape parameters corresponding to the maximum R-parameter (R-max) as shown in Table 9. Figure 11b–d shows the horizontal and vertical derivative anomalies, the amplitude analytic signal, and the R-parameter 2-D mosaic of the Louisiana anomaly. It is found that the R-max value is 0.96 corresponds to q = 0.9 which is approximated by horizontal cylinder shape. The best estimated model parameters are: q = 0.9, z_o = 2950 m, x_o = 400 m, and A = −3282.61 mGal m.

q	R-max
0.5	0.943033
0.6	0.934839
0.7	0.926284
0.8	0.957334
0.9	0.963889
1	0.957410
1.1	0.942955
1.2	0.926284
1.3	0.909983
1.4	0.894917
1.5	0.881254

Table 9.

The Louisiana dome gravity anomaly, USA. The R-parameter calculated from different shape factors.

The Louisiana dome anomaly has been interpreted by different authors assuming a horizontal source to determine the depth to the center of the salt body. The obtained results have a good agreement with those depths to the center that obtained by the published literatures of [5, 41] (Table 10). In addition, the proposed method has the lowest misfit compared to the other method (Table 10).

Model parameters	Approaches and techniques of		Present study
	[47]	[5]
A mGal m	−16400	−16021	−3282.61
zo (m)	2899	2702.2	2950.00
q	1	1	0.9
xo (m)	—	506.5	400
Misfit	12.4%	2.9 × 10^−3	1.3 × 10^−3

Table 10.

The Louisiana dome gravity anomaly, USA. The estimated parameters.

6. Discussion

The proposed method of R-parameter imaging technique was deployed to visualizes the salt dome anomalies from the gravity data measured across a 2D profile. The method fitting the anomaly of the measured gravity profile by a single geometric shape body (sphere & cylinder). Such as the spherical source (q = 1.5) suggestion based on the maximum R-parameter (R-max = 0.99) of the Humble dome anomaly (Figure 10) and the approximated horizontal cylinder (q = 0.9) with maximum R-parameter (R-max = 0.96) of the Louisiana dome anomaly (Figure 11).

The obtained results by the R-parameter method of the Humble dome anomaly was compared with other results in the published literature (Table 7) and confirmed with drilling to insure the depth to the top of the buried anomaly (Table 8). For the Louisiana dome anomaly, the obtained results by R-parameter approach was compared with the pervious published literature and weighted by the misfit error between the observed and calculated anomaly for the different techniques used (Table 10) to increase the efficiency of the proposed method.

In over all the obtained results using the R-parameter method to investigate the salt dome anomalies is good and acceptable in the two given field examples.

7. Conclusions

In this study, we have introduced and investigated the applicability and the performance of the R-parameter imaging method in elucidating distinctive physical parameters (A, z_o, x_o, q) of simple geometrically-shaped geologic structures (spheres, horizontal cylinders and vertical cylinders) from gravity data along profiles. This inversion imaging method has been demonstrated successfully on numerically generated gravity anomalies corrupted with random noise, applied to cases of anomalies from multiple and interfering structures and finally experimented on a two different field cases for salt domes in USA. The approach presented here can be used to investigate salt domes and perform reconnaissance geological studies.