- Email: [email protected]
The semiautomatic interpretation of gravity profile data
Computers & Geosciences 37 (2011) 1102–1109
Contents lists available at ScienceDirect
Computers & Geosciences journal homepage: www.elsevier.com/locate/cageo
The semiautomatic interpretation of gravity profile data$ G.R.J. Cooper School of Geosciences, University of the Witwatersrand, Johannesburg, South Africa
a r t i c l e i n f o
a b s t r a c t
Article history: Received 2 March 2010 Received in revised form 1 February 2011 Accepted 22 February 2011 Available online 11 March 2011
Potential field semiautomatic interpretation techniques can significantly decrease the interpretation workload of a geophysicist and are widely used. This paper generalises the tilt-depth method to gravity data using the horizontal cylinder and the buried sphere models. The method is compared to the Euler deconvolution and is applied to gravity data from South Africa. Matlab source code can be downloaded from the IAMG server at www.iamg.org. & 2011 Elsevier Ltd. All rights reserved.
Keywords: Tilt-depth Euler deconvolution
1. Introduction
2. Application to gravity profile data
Semiautomatic interpretation techniques are widely used in potential field interpretation projects. There are several techniques available such as Werner deconvolution (Hartman et al., 1971), SLUTH (Thurston and Smith, 2007), and wavelet methods (Sailhac and Gibert, 2003). Euler deconvolution (Thomson, 1982; Reid et al., 1990) is possibly the technique in widest use currently. In 2007 Salem et al. (2007) introduced the tilt-depth method for the magnetic anomaly over a contact. Previously Miller and Singh had developed the tilt angle as a method of enhancing images of the vertical derivative of potential field data. The tilt angle T of the field f was defined as (Miller and Singh, 1994)
The tilt-depth method will be both generalised (by applying it to gravity models) and extended (by using all values of the ratio of the field gradients, not just a single value) here.
0
1
@f =@z B C T ¼ tan1 @qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiA: ðð@f =@xÞ2 þ ð@f =@yÞ2 Þ
ð1Þ
Assuming both contact and direction of magnetisation to be vertical, then when the analytic expressions for the field gradients are substituted into Eq. (1) the result is T ¼ tan1
h zc
2.1. Horizontal cylinder model The gravity anomaly over a horizontal cylinder is given by (Beck, 1981, p. 79) g¼
2pGa2 rz x2 þ z2
ð3Þ
where the cylinder radius is a and its density contrast is r. The first horizontal and vertical derivatives of g are then @g 4pGa2 rxz ¼ @x ðx2 þz2 Þ2
ð4Þ
@g 4pGa2 rðx2 z2 Þ ¼ : @z ðx2 þ z2 Þ2
ð5Þ
The ratio of the vertical to the horizontal derivative is ð2Þ
where h and zc are the horizontal and vertical distances, respectively, to the contact. Salem et al. pointed out that this indicated that the depth to the contact could be determined as half the distance between the 451 contours of the tilt angle of the magnetic field.
r¼
x2 z2 2xz
ð6Þ
which can be rearranged as x2 þ 2xzrz2 ¼ 0
ð7Þ
which has solutions x1 ¼ rz þ zOðr 2 þ 1Þ
$
Code available from: http://www.iamg.org/CGEditor/index.htm. E-mail address: [email protected]
0098-3004/$ - see front matter & 2011 Elsevier Ltd. All rights reserved. doi:10.1016/j.cageo.2011.02.016
x2 ¼ rzzOðr 2 þ 1Þ
ð8Þ
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
1103
Fig. 1. (a) (i) Gravity anomaly from two cylinders located at positions of 55 and 145 m. The cylinders have depths of 15 m and density contrasts of 1 kg m 3. The data sampling interval was 1 m. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. A minimum data point separation of 2 m and a maximum data point separation of 30 m were used. All solutions are shown. Solutions from Euler deconvolution (using a structural index of 1) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (b) (i) Gravity anomaly from two cylinders located at positions of 85 and 115 m. The cylinders have depths of 15 m and a density contrast of 1 kg m 3. The data sampling interval was 1 m. (ii) Gradient ratio r as defined in Eq. (6). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 30 m were used. Solutions from Euler deconvolution (using a structural index of 1) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (c) (i) Gravity anomaly from two cylinders located at positions of 92.5 and 107.5 m. The cylinders have depths of 15 m and a density contrast of 1 kg m 3. The data sampling interval was 1 m. (ii) Gradient ratio r as defined in Eq. (6). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 30 m were used. Solutions from Euler deconvolution (using a structural index of 1) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
1104
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
Fig. 2. (a) (i) Gravity anomaly from two cylinders located at positions of 55 and 145 m. The cylinders have depths of 15 m and density contrasts of 1 kg m 3. The data sampling interval was 1 m. Uniformly distributed random noise with an amplitude equal to 1% of the synthetic data amplitude has been added. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. A minimum data point separation of 2 m and a maximum data point separation of 30 m were used. All solutions are shown. Solutions from Euler deconvolution (using a structural index of 1) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (b) (i) Gravity anomaly from two cylinders located at positions of 55 and 145 m. The cylinders have depths of 15 m and density contrasts of 1 kg m 3. The data sampling interval was 1 m. Uniformly distributed random noise with an amplitude equal to 5% of the synthetic data amplitude has been added. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. A minimum data point separation of 2 m and a maximum data point separation of 30 m were used. All solutions are shown. Solutions from Euler deconvolution (using a structural index of 1) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (c) (i) Gravity anomaly from two cylinders located at positions of 55 and 145 m. The cylinders have depths of 15 m and density contrasts of 1 kg m 3. The data sampling interval was 1 m. Uniformly distributed random noise with an amplitude equal to 10% of the synthetic data amplitude has been added. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. A minimum data point separation of 2 m and a maximum data point separation of 30 m were used. All solutions are shown. Solutions from Euler deconvolution (using a structural index of 1) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
1105
Fig. 3. (a) Gravity data over the Trompsburg structure, South Africa, showing location of extracted profiles. (b) Extracted profile P1–P2. (c) Gradient ratio r as defined in Eq. (6). (d) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 5 km and a maximum data point separation of 40 km were used. Solutions from Euler deconvolution (using a structural index of 1 and a window size of 21 points) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
so that x1 x2 z ¼ pffiffiffiffiffiffiffiffiffiffiffiffi : 2 r2 þ 1
ð9Þ
Hence if the locations x1 and x2 of two points on the profile, which have the same gradient ratio r are known, then the cylinder depth can be obtained from Eq. (9) and subsequently its location from Eq. (8). The strategy is, therefore, for every value on the gradient ratio profile r find another point with the same value. In practise
the second point will lie between two adjacent values of r. To reduce problems with noise and interference only those data pairs with a separation lying between chosen minimum and maximum values are used. The maximum value should be approximately the width of the anomaly. Fig. 1a shows a synthetic dataset consisting of anomalies from two cylinders separated by a distance equal to six times their depth. The method correctly locates both cylinders very accurately, as can be seen. For comparison solutions obtained from the Euler
1106
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
Fig. 4. (a) (i) Gravity anomaly from two spheres located at positions of 55 and 145 m. There are 200 points on the profile. The spheres have depths of 15 m and masses of 1010 kg. The data sampling interval was 1 m. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 40 m were used. Solutions from Euler deconvolution (using a structural index of 2) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (b) (i) Gravity anomaly from two spheres located at positions of 92.5 and 112.5 m. There are 200 points on the profile. The spheres have depths of 15 m and masses of 1010 kg. The data sampling interval was 1 m. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 40 m were used. Solutions from Euler deconvolution (using a structural index of 2) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (c) (i) Gravity anomaly from two spheres located at positions of 92.5 and 112.5 m. There are 200 points on the profile. The spheres have depths of 15 m and masses of 1010 kg. The data sampling interval was 1 m. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as redþ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 40 m were used. Solutions from Euler deconvolution (using a structural index of 2) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
1107
Fig. 5. (a) (i) Gravity anomaly from two spheres located at positions of 55 and 145 m. There are 200 points on the profile. The spheres have depths of 15 m and masses of 1010 kg. The data sampling interval was 1 m. Uniformly distributed random noise with an amplitude equal to 1% of the synthetic data amplitude has been added. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 40 m were used. Solutions from Euler deconvolution (using a structural index of 2 and a window size of 11 points) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. (b) (i) Gravity anomaly from two spheres located at positions of 55 and 145 m. There are 200 points on the profile. The spheres have depths of 15 m and masses of 1010 kg. The data sampling interval was 1 m. Uniformly distributed random noise with an amplitude equal to 5% of the synthetic data amplitude has been added. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 40 m were used. Solutions from Euler deconvolution (using a structural index of 2 and a window size of 11 points) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. (c) (i) Gravity anomaly from two spheres located at positions of 55 and 145 m. There are 200 points on the profile. The spheres have depths of 15 m and masses of 1010 kg. The data sampling interval was 1 m. Uniformly distributed random noise with an amplitude equal to 10% of the synthetic data amplitude has been added. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 40 m were used. Solutions from Euler deconvolution (using a structural index of 2 and a window size of 11 points) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
1108
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
deconvolution are also shown. They also locate the cylinders, but with much poorer horizontal and vertical resolution. When the cylinders are moved closer together then interference becomes an issue. In Fig. 1b the cylinders are now separated by a distance equal to twice their depth, and this causes a scattering of the solutions. In Fig. 1c the separation is reduced to a distance equal to their depth. Both the gradient ratio method and the Euler deconvolution are now unable to resolve the individual sources and produce a set of solutions, which encompasses both. The sensitivity of the method to noise is investigated in Fig. 2. Fig. 2a shows the synthetic gravity datasets from Fig. 1 when uniformly distributed random noise equal to 1% of the data amplitude has been added. In Fig. 2b and c the amount of noise is increased to 5% and then 10%. The noise causes a scattering of the previously tightly focussed solutions, and their depth is reduced. In addition other solutions appear along the profile because the different peaks that the noise produces are indistinguishable from small anomalies. These solutions can be reduced by increasing the minimum distance between the data points that are used by the method. Fig. 3a shows a gravity dataset from South Africa. The circular feature at the centre of the image is the Trompsburg high. A profile across the linear feature in the South of the data was extracted (P1–P2) and analysed using the gradient method. Fig. 3b shows the extracted profile, Fig. 3c shows the gradient ratio r, and Fig. 3d shows the solutions obtained. Euler deconvolution gives a spray of solutions with depths from 0 to 25 km grouped at about 20 km in the vicinity of the anomaly. The gradient method produces far fewer solutions than the Euler deconvolution but they are much more tightly grouped and give a depth of 15–20 km. 2.2. Buried sphere model The analysis for the buried sphere model follows closely that of the horizontal cylinder. The vertical component of the gravity
anomaly along a profile over the centre of a buried sphere of mass m is (Beck, 1981, p.73) g¼
Gmz
:
ð10Þ
@g 3Gmzx ¼ @x ðx2 þ z2 Þ5=2
ð11Þ
ðx2 þ z2 Þ3=2
Hence
and @g Gmðx2 2z2 Þ ¼ @z ðx2 þz2 Þ5=2
ð12Þ
so the ratio of the vertical to the horizontal derivative is r¼
ðx2 2z2 Þ 3zx
x2 þ 3xzr2z2 ¼ 0 which has solutions pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 z ð9r 2 þ8Þ x1 ¼ rz 2 2 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 z ð9r 2 þ 8Þ x2 ¼ rz þ 2 2
ð13Þ
ð14Þ
so that ðx2 x1 Þ z ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi : ð9r 2 þ 8Þ
ð15Þ
The strategy is identical to that applied to the horizontal cylinder model; two points on the profile are found for each value of r, enabling the depth and then the location of the source to be established. Fig. 4a shows a synthetic dataset consisting of the gravity anomalies from two spheres. The gradient method correctly locates the spheres and the solutions cluster more
Fig. 6. (a) Gravity data over the Trompsburg structure, South Africa. The location of the profile is marked as P3–P4 in Fig. 3a. The data points were 1 km apart. (b) Gradient ratio r as defined in Eq. (13). (c) Solutions from the gradient ratio method are marked as red þsymbols. All solutions are shown. A minimum data point separation of 10 km and a maximum data point separation of 60 km were used. Solutions from Euler deconvolution (using a structural index of 2 and a window size of 21 data points) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
tightly than those from the Euler deconvolution. As the separation of the spheres is reduced, interference causes the solutions to scatter (Fig. 4b and c); until once their separation is equal to their depth the individual spheres can no longer be distinguished. Fig. 5 shows the effect of adding random noise to the gravity data. The dispersion of the solutions from the gradient method and from Euler deconvolution increases as the noise level is increased, and additional spurious solutions appear throughout the profile. Fig. 3a shows the location of a profile extracted over the Trompsburg structure gravity data (P3–P4). This structure is a large layered intrusion that was discovered in 1942 by B.D. Maree (Buchmann, 1960). It is a gravity high 100 mGal in amplitude, and little is known about its age or composition. Fig. 6a plots the profile, Fig. 6b shows the gradient ratio r, and Fig. 6c shows the extracted solutions. Euler deconvolution gave a spray of sources associated with the main anomaly, whose depths varied from zero to greater than 40 km, whereas the gradient method produced sources with depths lying mostly between 10 and 20 km. Previous work using different interpretation techniques (e.g., Cooper, 2004) gave similar depths for this structure when a buried sphere model was used. The Trompsburg gravity anomaly is clearly more complex than that from a single simple source such as a buried sphere.
3. Conclusions A new technique, based on the tilt-depth method, has been introduced for the semiautomatic interpretation of gravity data profiles. The method uses the horizontal cylinder and sphere models. Its sensitivity to random noise and to interference was tested on synthetic data. The method was also applied to gravity data from South Africa.
1109
Acknowledgments The Council for Geoscience, Pretoria, are thanked for permission to use the gravity data shown in Fig. 3a. The NRF (Pretoria) is thanked for funding this project. Dr. P. Keating is thanked for his constructive comments during the review process.
Appendix A. Supplementary information Supplementary data associated with this article can be found in the online version at doi:10.1016/j.cageo.2011.02.016.
References Buchmann, J.P., 1960. Exploration of a geophysical anomaly at Trompsburg, Orange Free State, South Africa. Transactions of the Geological Society of South Africa 63, 1–10. Cooper, G.R.J., 2004. A semi-automatic procedure for the interpretation of geophysical data. Exploration Geophysics 35 (3), 180–185. Hartman, R.R., Teskey, D.J., Friedberg, J.L., 1971. A system for rapid digital aeromagnetic interpretation. Geophysics 36, 891–918. Miller, H.G., Singh, V., 1994. Potential field tilt—a new concept for location of potential field sources. Journal of Applied Geophysics 32, 213–217. Reid, A.B., Allsop, J.M., Granser, H., Millet, A.J., Somerton, I.W., 1990. Magnetic interpretation in three dimensions using Euler deconvolution. Geophysics 55, 80–91. Sailhac, P., Gibert, D., 2003. Identification of sources of potential fields with the continuous wavelet transform: two dimensional wavelets and multipolar approximations. Journal of Geophysical Research 108, 2262. Salem, A., Williams, S., Fairhead, J.D., Ravat, D., 2007. Tilt-depth method: a simple depth estimation method using first-order magnetic derivatives. The Leading Edge, 1502–1505. Thomson, D.T., 1982. Euldph: a new technique for making computer assisted depth estimates from magnetic data. Geophysics 47 (1), 31–37. Thurston, J., Smith, R., 2007. Source location using total-field homogeneity: introducing the SLUTH method for depth estimation. The Leading Edge, 1272–1276.
Contents lists available at ScienceDirect
Computers & Geosciences journal homepage: www.elsevier.com/locate/cageo
The semiautomatic interpretation of gravity profile data$ G.R.J. Cooper School of Geosciences, University of the Witwatersrand, Johannesburg, South Africa
a r t i c l e i n f o
a b s t r a c t
Article history: Received 2 March 2010 Received in revised form 1 February 2011 Accepted 22 February 2011 Available online 11 March 2011
Potential field semiautomatic interpretation techniques can significantly decrease the interpretation workload of a geophysicist and are widely used. This paper generalises the tilt-depth method to gravity data using the horizontal cylinder and the buried sphere models. The method is compared to the Euler deconvolution and is applied to gravity data from South Africa. Matlab source code can be downloaded from the IAMG server at www.iamg.org. & 2011 Elsevier Ltd. All rights reserved.
Keywords: Tilt-depth Euler deconvolution
1. Introduction
2. Application to gravity profile data
Semiautomatic interpretation techniques are widely used in potential field interpretation projects. There are several techniques available such as Werner deconvolution (Hartman et al., 1971), SLUTH (Thurston and Smith, 2007), and wavelet methods (Sailhac and Gibert, 2003). Euler deconvolution (Thomson, 1982; Reid et al., 1990) is possibly the technique in widest use currently. In 2007 Salem et al. (2007) introduced the tilt-depth method for the magnetic anomaly over a contact. Previously Miller and Singh had developed the tilt angle as a method of enhancing images of the vertical derivative of potential field data. The tilt angle T of the field f was defined as (Miller and Singh, 1994)
The tilt-depth method will be both generalised (by applying it to gravity models) and extended (by using all values of the ratio of the field gradients, not just a single value) here.
0
1
@f =@z B C T ¼ tan1 @qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiA: ðð@f =@xÞ2 þ ð@f =@yÞ2 Þ
ð1Þ
Assuming both contact and direction of magnetisation to be vertical, then when the analytic expressions for the field gradients are substituted into Eq. (1) the result is T ¼ tan1
h zc
2.1. Horizontal cylinder model The gravity anomaly over a horizontal cylinder is given by (Beck, 1981, p. 79) g¼
2pGa2 rz x2 þ z2
ð3Þ
where the cylinder radius is a and its density contrast is r. The first horizontal and vertical derivatives of g are then @g 4pGa2 rxz ¼ @x ðx2 þz2 Þ2
ð4Þ
@g 4pGa2 rðx2 z2 Þ ¼ : @z ðx2 þ z2 Þ2
ð5Þ
The ratio of the vertical to the horizontal derivative is ð2Þ
where h and zc are the horizontal and vertical distances, respectively, to the contact. Salem et al. pointed out that this indicated that the depth to the contact could be determined as half the distance between the 451 contours of the tilt angle of the magnetic field.
r¼
x2 z2 2xz
ð6Þ
which can be rearranged as x2 þ 2xzrz2 ¼ 0
ð7Þ
which has solutions x1 ¼ rz þ zOðr 2 þ 1Þ
$
Code available from: http://www.iamg.org/CGEditor/index.htm. E-mail address: [email protected]
0098-3004/$ - see front matter & 2011 Elsevier Ltd. All rights reserved. doi:10.1016/j.cageo.2011.02.016
x2 ¼ rzzOðr 2 þ 1Þ
ð8Þ
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
1103
Fig. 1. (a) (i) Gravity anomaly from two cylinders located at positions of 55 and 145 m. The cylinders have depths of 15 m and density contrasts of 1 kg m 3. The data sampling interval was 1 m. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. A minimum data point separation of 2 m and a maximum data point separation of 30 m were used. All solutions are shown. Solutions from Euler deconvolution (using a structural index of 1) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (b) (i) Gravity anomaly from two cylinders located at positions of 85 and 115 m. The cylinders have depths of 15 m and a density contrast of 1 kg m 3. The data sampling interval was 1 m. (ii) Gradient ratio r as defined in Eq. (6). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 30 m were used. Solutions from Euler deconvolution (using a structural index of 1) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (c) (i) Gravity anomaly from two cylinders located at positions of 92.5 and 107.5 m. The cylinders have depths of 15 m and a density contrast of 1 kg m 3. The data sampling interval was 1 m. (ii) Gradient ratio r as defined in Eq. (6). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 30 m were used. Solutions from Euler deconvolution (using a structural index of 1) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
1104
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
Fig. 2. (a) (i) Gravity anomaly from two cylinders located at positions of 55 and 145 m. The cylinders have depths of 15 m and density contrasts of 1 kg m 3. The data sampling interval was 1 m. Uniformly distributed random noise with an amplitude equal to 1% of the synthetic data amplitude has been added. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. A minimum data point separation of 2 m and a maximum data point separation of 30 m were used. All solutions are shown. Solutions from Euler deconvolution (using a structural index of 1) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (b) (i) Gravity anomaly from two cylinders located at positions of 55 and 145 m. The cylinders have depths of 15 m and density contrasts of 1 kg m 3. The data sampling interval was 1 m. Uniformly distributed random noise with an amplitude equal to 5% of the synthetic data amplitude has been added. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. A minimum data point separation of 2 m and a maximum data point separation of 30 m were used. All solutions are shown. Solutions from Euler deconvolution (using a structural index of 1) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (c) (i) Gravity anomaly from two cylinders located at positions of 55 and 145 m. The cylinders have depths of 15 m and density contrasts of 1 kg m 3. The data sampling interval was 1 m. Uniformly distributed random noise with an amplitude equal to 10% of the synthetic data amplitude has been added. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. A minimum data point separation of 2 m and a maximum data point separation of 30 m were used. All solutions are shown. Solutions from Euler deconvolution (using a structural index of 1) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
1105
Fig. 3. (a) Gravity data over the Trompsburg structure, South Africa, showing location of extracted profiles. (b) Extracted profile P1–P2. (c) Gradient ratio r as defined in Eq. (6). (d) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 5 km and a maximum data point separation of 40 km were used. Solutions from Euler deconvolution (using a structural index of 1 and a window size of 21 points) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
so that x1 x2 z ¼ pffiffiffiffiffiffiffiffiffiffiffiffi : 2 r2 þ 1
ð9Þ
Hence if the locations x1 and x2 of two points on the profile, which have the same gradient ratio r are known, then the cylinder depth can be obtained from Eq. (9) and subsequently its location from Eq. (8). The strategy is, therefore, for every value on the gradient ratio profile r find another point with the same value. In practise
the second point will lie between two adjacent values of r. To reduce problems with noise and interference only those data pairs with a separation lying between chosen minimum and maximum values are used. The maximum value should be approximately the width of the anomaly. Fig. 1a shows a synthetic dataset consisting of anomalies from two cylinders separated by a distance equal to six times their depth. The method correctly locates both cylinders very accurately, as can be seen. For comparison solutions obtained from the Euler
1106
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
Fig. 4. (a) (i) Gravity anomaly from two spheres located at positions of 55 and 145 m. There are 200 points on the profile. The spheres have depths of 15 m and masses of 1010 kg. The data sampling interval was 1 m. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 40 m were used. Solutions from Euler deconvolution (using a structural index of 2) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (b) (i) Gravity anomaly from two spheres located at positions of 92.5 and 112.5 m. There are 200 points on the profile. The spheres have depths of 15 m and masses of 1010 kg. The data sampling interval was 1 m. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 40 m were used. Solutions from Euler deconvolution (using a structural index of 2) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (c) (i) Gravity anomaly from two spheres located at positions of 92.5 and 112.5 m. There are 200 points on the profile. The spheres have depths of 15 m and masses of 1010 kg. The data sampling interval was 1 m. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as redþ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 40 m were used. Solutions from Euler deconvolution (using a structural index of 2) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. The window size used was 11 points. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
1107
Fig. 5. (a) (i) Gravity anomaly from two spheres located at positions of 55 and 145 m. There are 200 points on the profile. The spheres have depths of 15 m and masses of 1010 kg. The data sampling interval was 1 m. Uniformly distributed random noise with an amplitude equal to 1% of the synthetic data amplitude has been added. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 40 m were used. Solutions from Euler deconvolution (using a structural index of 2 and a window size of 11 points) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. (b) (i) Gravity anomaly from two spheres located at positions of 55 and 145 m. There are 200 points on the profile. The spheres have depths of 15 m and masses of 1010 kg. The data sampling interval was 1 m. Uniformly distributed random noise with an amplitude equal to 5% of the synthetic data amplitude has been added. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 40 m were used. Solutions from Euler deconvolution (using a structural index of 2 and a window size of 11 points) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. (c) (i) Gravity anomaly from two spheres located at positions of 55 and 145 m. There are 200 points on the profile. The spheres have depths of 15 m and masses of 1010 kg. The data sampling interval was 1 m. Uniformly distributed random noise with an amplitude equal to 10% of the synthetic data amplitude has been added. (ii) Gradient ratio r as defined in Eq. (13). (iii) Solutions from the gradient ratio method are marked as red þ symbols. All solutions are shown. A minimum data point separation of 2 m and a maximum data point separation of 40 m were used. Solutions from Euler deconvolution (using a structural index of 2 and a window size of 11 points) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
1108
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
deconvolution are also shown. They also locate the cylinders, but with much poorer horizontal and vertical resolution. When the cylinders are moved closer together then interference becomes an issue. In Fig. 1b the cylinders are now separated by a distance equal to twice their depth, and this causes a scattering of the solutions. In Fig. 1c the separation is reduced to a distance equal to their depth. Both the gradient ratio method and the Euler deconvolution are now unable to resolve the individual sources and produce a set of solutions, which encompasses both. The sensitivity of the method to noise is investigated in Fig. 2. Fig. 2a shows the synthetic gravity datasets from Fig. 1 when uniformly distributed random noise equal to 1% of the data amplitude has been added. In Fig. 2b and c the amount of noise is increased to 5% and then 10%. The noise causes a scattering of the previously tightly focussed solutions, and their depth is reduced. In addition other solutions appear along the profile because the different peaks that the noise produces are indistinguishable from small anomalies. These solutions can be reduced by increasing the minimum distance between the data points that are used by the method. Fig. 3a shows a gravity dataset from South Africa. The circular feature at the centre of the image is the Trompsburg high. A profile across the linear feature in the South of the data was extracted (P1–P2) and analysed using the gradient method. Fig. 3b shows the extracted profile, Fig. 3c shows the gradient ratio r, and Fig. 3d shows the solutions obtained. Euler deconvolution gives a spray of solutions with depths from 0 to 25 km grouped at about 20 km in the vicinity of the anomaly. The gradient method produces far fewer solutions than the Euler deconvolution but they are much more tightly grouped and give a depth of 15–20 km. 2.2. Buried sphere model The analysis for the buried sphere model follows closely that of the horizontal cylinder. The vertical component of the gravity
anomaly along a profile over the centre of a buried sphere of mass m is (Beck, 1981, p.73) g¼
Gmz
:
ð10Þ
@g 3Gmzx ¼ @x ðx2 þ z2 Þ5=2
ð11Þ
ðx2 þ z2 Þ3=2
Hence
and @g Gmðx2 2z2 Þ ¼ @z ðx2 þz2 Þ5=2
ð12Þ
so the ratio of the vertical to the horizontal derivative is r¼
ðx2 2z2 Þ 3zx
x2 þ 3xzr2z2 ¼ 0 which has solutions pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 z ð9r 2 þ8Þ x1 ¼ rz 2 2 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 z ð9r 2 þ 8Þ x2 ¼ rz þ 2 2
ð13Þ
ð14Þ
so that ðx2 x1 Þ z ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi : ð9r 2 þ 8Þ
ð15Þ
The strategy is identical to that applied to the horizontal cylinder model; two points on the profile are found for each value of r, enabling the depth and then the location of the source to be established. Fig. 4a shows a synthetic dataset consisting of the gravity anomalies from two spheres. The gradient method correctly locates the spheres and the solutions cluster more
Fig. 6. (a) Gravity data over the Trompsburg structure, South Africa. The location of the profile is marked as P3–P4 in Fig. 3a. The data points were 1 km apart. (b) Gradient ratio r as defined in Eq. (13). (c) Solutions from the gradient ratio method are marked as red þsymbols. All solutions are shown. A minimum data point separation of 10 km and a maximum data point separation of 60 km were used. Solutions from Euler deconvolution (using a structural index of 2 and a window size of 21 data points) are plotted as blue symbols. All solutions that fell below the ground surface are plotted. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
G.R.J. Cooper / Computers & Geosciences 37 (2011) 1102–1109
tightly than those from the Euler deconvolution. As the separation of the spheres is reduced, interference causes the solutions to scatter (Fig. 4b and c); until once their separation is equal to their depth the individual spheres can no longer be distinguished. Fig. 5 shows the effect of adding random noise to the gravity data. The dispersion of the solutions from the gradient method and from Euler deconvolution increases as the noise level is increased, and additional spurious solutions appear throughout the profile. Fig. 3a shows the location of a profile extracted over the Trompsburg structure gravity data (P3–P4). This structure is a large layered intrusion that was discovered in 1942 by B.D. Maree (Buchmann, 1960). It is a gravity high 100 mGal in amplitude, and little is known about its age or composition. Fig. 6a plots the profile, Fig. 6b shows the gradient ratio r, and Fig. 6c shows the extracted solutions. Euler deconvolution gave a spray of sources associated with the main anomaly, whose depths varied from zero to greater than 40 km, whereas the gradient method produced sources with depths lying mostly between 10 and 20 km. Previous work using different interpretation techniques (e.g., Cooper, 2004) gave similar depths for this structure when a buried sphere model was used. The Trompsburg gravity anomaly is clearly more complex than that from a single simple source such as a buried sphere.
3. Conclusions A new technique, based on the tilt-depth method, has been introduced for the semiautomatic interpretation of gravity data profiles. The method uses the horizontal cylinder and sphere models. Its sensitivity to random noise and to interference was tested on synthetic data. The method was also applied to gravity data from South Africa.
1109
Acknowledgments The Council for Geoscience, Pretoria, are thanked for permission to use the gravity data shown in Fig. 3a. The NRF (Pretoria) is thanked for funding this project. Dr. P. Keating is thanked for his constructive comments during the review process.
Appendix A. Supplementary information Supplementary data associated with this article can be found in the online version at doi:10.1016/j.cageo.2011.02.016.
References Buchmann, J.P., 1960. Exploration of a geophysical anomaly at Trompsburg, Orange Free State, South Africa. Transactions of the Geological Society of South Africa 63, 1–10. Cooper, G.R.J., 2004. A semi-automatic procedure for the interpretation of geophysical data. Exploration Geophysics 35 (3), 180–185. Hartman, R.R., Teskey, D.J., Friedberg, J.L., 1971. A system for rapid digital aeromagnetic interpretation. Geophysics 36, 891–918. Miller, H.G., Singh, V., 1994. Potential field tilt—a new concept for location of potential field sources. Journal of Applied Geophysics 32, 213–217. Reid, A.B., Allsop, J.M., Granser, H., Millet, A.J., Somerton, I.W., 1990. Magnetic interpretation in three dimensions using Euler deconvolution. Geophysics 55, 80–91. Sailhac, P., Gibert, D., 2003. Identification of sources of potential fields with the continuous wavelet transform: two dimensional wavelets and multipolar approximations. Journal of Geophysical Research 108, 2262. Salem, A., Williams, S., Fairhead, J.D., Ravat, D., 2007. Tilt-depth method: a simple depth estimation method using first-order magnetic derivatives. The Leading Edge, 1502–1505. Thomson, D.T., 1982. Euldph: a new technique for making computer assisted depth estimates from magnetic data. Geophysics 47 (1), 31–37. Thurston, J., Smith, R., 2007. Source location using total-field homogeneity: introducing the SLUTH method for depth estimation. The Leading Edge, 1272–1276.