- Email: [email protected]
Balanced horizontal derivative of potential field data to recognize the edges and estimate location parameters of the source
Journal of Applied Geophysics 108 (2014) 12–18
Contents lists available at ScienceDirect
Journal of Applied Geophysics journal homepage: www.elsevier.com/locate/jappgeo
Balanced horizontal derivative of potential field data to recognize the edges and estimate location parameters of the source Guoqing Ma ⁎, Cai Liu, Lili Li College of Geoexploration Science and Technology, Jilin University, Changchun 130061, China
a r t i c l e
i n f o
Article history: Received 15 May 2014 Accepted 11 June 2014 Available online 18 June 2014 Keywords: Balanced horizontal derivative Edge detection Linear equation Depth Normalized
a b s t r a c t Edge detection is a common task in the interpretation of potential field data, and many edge detection filters are presented to accomplish this task, which are the functions composed by the first-order horizontal derivative and vertical derivative of potential field data, but we find that the edges recognized by the existing edge detection filters are bigger than the true edges. In this paper, we present the balanced horizontal derivative (BHD) edge detection filter, which uses the ratio of the first-order horizontal derivative to the second-order horizontal derivative to recognize the edges of the source, and the recognized edges by the BHD edge detection filter are more correct and are more insensitive to noise. We derive a linear equation based on the derivatives of the BHD to estimate the depth of the source, and we also present a normalized total horizontal derivative (NTHD) method to image the location of the source. We demonstrate the presented methods on synthetic potential field data, and the results show that the presented methods can provide the edges and location parameters of the sources correctly, and the BHD filter can display the edges more clearly and correctly. At last, we apply the presented methods to real potential field data, and the inversion results computed by the presented methods are in accord with the geology information. © 2014 Elsevier B.V. All rights reserved.
1. Introduction Edge detection is a commonly requested task in the interpretation of potential field data, and the maxima of total horizontal derivative and the zero of the vertical derivative are corresponding to the edges of the source (Even, 1936; Hood and Tasked, 1989; Cordell, 1979), so the edge detection filters are all the functions based on the horizontal derivative and vertical derivatives of the data (Blakely, 1995; Cordell and Grauch, 1985; Roest et al., 1992; Hsu et al., 1996; Fedi and Florio, 2001). Initially, people directly used the total horizontal derivative or vertical derivative to recognize the source edges (Even, 1936; Cordell, 1979), but these methods cannot display the edges of deeper source clearly, so many people began to develop the balanced edge detection filters. Miller and Singh (1994) presented the first balanced edge detection filter-tilt angle, which is the ratio of the vertical derivative to total horizontal derivative, and this filter can effectively balance the amplitudes of the anomalies generated by the sources with different depths. Some authors also presented other format of balanced edge detection filters that can display the edges of the shallow and deep sources simultaneously (Miller and Singh, 1994; Rajagopalan and Milligan, 1995; Verduzco et al., 2004; Wijns et al., 2005; Cooper and Cowan, 2006; 2008; Ma and Li; 2012; Ma, 2013a,b).
We find that the horizontal coordinates of the maxima of total horizontal derivative and the zero of vertical derivative are both bigger than the true locations of the edges, and the coordinates of the maxima of total horizontal derivative are more close to the true edges. We prove that the ratio of the total horizontal derivative to the second-order vertical derivative can obtain more accurate edge detection results. We also introduce a Laplace equation that uses the sum of two second-order horizontal derivatives to compute the second-order vertical derivative, which cannot increase the interference of noise dramatically, and this method is called as balanced horizontal derivative (BHD) filter. The depth estimation is another important task in the interpretation of potential field data. We derive a linear equation based on the derivatives of the BHD, which can obtain the horizontal location and depth of the source without any a priori knowledge. We also present a normalized total horizontal derivative (NTHD) method, and the location and the depth of the source can be obtained depending on the closed maxima of the results. We demonstrate the presented methods on synthetic and real potential field data, and the inversion results are consistent with the real values.
2. Balanced horizontal derivate (BHD) filter ⁎ Corresponding author at. College of Geoexploration Science and Technology, Jilin University, Changchun, China. E-mail address: [email protected] (G. Ma).
http://dx.doi.org/10.1016/j.jappgeo.2014.06.005 0926-9851/© 2014 Elsevier B.V. All rights reserved.
Total horizontal derivative, Theta map and normalized horizontal tilt angle filter are three commonly used edge detection filters (Cooper
G. Ma et al. / Journal of Applied Geophysics 108 (2014) 12–18
and Cowan, 2006). Total horizontal derivative (THD) can be given by s ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 ∂f 2 ∂f : þ THD ¼ ∂x ∂y
The expression of second-order vertical derivative of the 2-D vertical prism is ð1Þ
"
V zzz
# xþa xþa x−a x−a : ¼ 2Gρ − − þ ðx þ aÞ2 þ h2 ðx þ aÞ2 þ H 2 ðx−aÞ2 þ h2 ðx−aÞ2 þ H 2
ð9Þ
The Theta map can be written by r ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2ffi ∂f 2 þ ∂∂yf ∂x ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Theta ¼ r 2 2 : ∂f 2 ∂f þ þ ∂∂zf ∂x ∂y
The coordinates that make the second-order vertical derivative get zero are ð2Þ x0;2
The normalized horizontal tilt angle (TDX) filter can be expressed as 0qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi1 2 2 −1 @ ð∂f =∂xÞ þ ð∂f =∂yÞ A TDX ¼ tan : ∂f =∂z
ð3Þ
The maxima of the above three edge detection filters automatically delineate the edges of the sources (Cooper and Cowan, 2008; Ma, 2013a,b). Depending on Eqs. (2) and (3) we can find that the Theta and TDX filters get the maxima when the first-order vertical derivative is equal to zero, which is not related to the value of horizontal derivative. The prism is a most common model in the interpretation of real data, and the expression of the first-order vertical derivative of 2-D vertical prism (Blakely, 1995) can be given by V zz
h −1 H −1 −1 H −1 h ¼ 2Gρ tan − tan − tan þ tan ð4Þ xþa xþa x−a x−a
where, 2a is the width of the prism, and h and H are the top and bottom depths of the prism, respectively. The coordinates that make the vertical derivative get zero are pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi x0;1 ¼ a2 þ Hh:
ð5Þ
We can find that the edges recognized by the Theta map and TDX filter are related to the depth of the source, and both are bigger than the true edges. The expressions of the horizontal vertical derivatives of 2-D vertical prism can be given by h
V zx
ih i 2 2 2 2 ðx þ aÞ þ H ðx−aÞ þ h : ¼ Gρ ln ðx þ aÞ2 þ h2 ðx−aÞ2 þ h2
ð6Þ
The coordinates that make the horizontal derivative getthe maxima are ffi sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 1 2 2 2 2 2 2 2 þ 12h H þ 2a − h þ H xm ¼ : 4a þ h þ H 6
ð7Þ
The relationship between the maxima of the horizontal derivative and the zero of vertical derivative is x Njx jNa: 0;1 m
13
ð8Þ
We can see that the coordinates of the maxima of horizontal derivative are bigger than the true edges, but are smaller than the coordinates of the zero of vertical derivative, so the edges recognized by the total horizontal derivative are more accurate, but the total horizontal derivative cannot give the edges of the deeper sources clearly.
ffi sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 1 2 2 2 2 2 2 2 ¼ þ 12h H þ 2a − h þ H 4a þ h þ H : ð10Þ 6
The coordinates of the zero of second-order vertical derivative are the same as the coordinates of the maxima of total horizontal derivative, so we suggest using the ratio of total horizontal derivative to the second-order vertical derivative to recognize the source edges. In order to lower the noise effect, we use the Laplace equation (Blakely, 1995) to compute the second order vertical derivative. 2
2
2
∂ f ∂ f ∂ f ¼− þ ∂z2 ∂x2 ∂y2
! ð11Þ
The second-order horizontal derivative can be computed in the space domain. This method is called as balanced horizontal derivative (BHD) edge detection filter, which can be given by qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 ð∂f =∂xÞ2 þ ð∂f =∂yÞ2 A k − ∂2 f =∂x2 þ ∂2 f =∂y2
0 −1 @
BHD ¼ tan
ð12Þ
where, ∂f/∂x, ∂f/∂y and ∂2f/∂z2 are the derivatives of the data f, and meanð∂ f =∂zÞ , and the maxima of the absolute value of the BHD are k ¼ mean ð∂2 f =∂z2 Þ corresponding to the edges of the source. We give a series of gravity profiles to show the distribution of different edge detection filters. Fig. 1a shows the synthetic gravity anomaly generated by two prisms with top depths of 15 and 20 m, respectively. Fig. 1b shows the total horizontal derivative (THD) of the data in Fig. 1a, and the THD can recognize the source edges, but cannot display the edges of the deep sources clearly. The results of the Theta map and TDX are similar, so we only give the Theta map of the data, and the Theta map of the data is shown in Fig. 1c. The Theta map can display the edges of shallow and deep sources simultaneously, but the recognized edges are significantly bigger than the true edges. Fig. 1d shows the absolute value of the BHD of the data in Fig. 1a. As can be seen from the results, the edges recognized by the THD and BHD are more accurate, and the BHD filter can display the edges of the shallow and deep sources simultaneously. Fig. 1e shows the absolute value of the BHD of the magnetic anomaly generated by the prism with an inclination of 60°, and the recognized edges have big difference with the real edges. Fig. 1f shows the absolute value of the BHD of the magnetic anomaly generated by the prism with an inclination of 90°, and the recognized edges are consistent with the true edges. We can see that the horizontal derivatives of magnetic anomaly are sensitive to magnetization direction, so we should reduce the magnetic anomaly to the North Pole when using the BHD to interpret the magnetic anomaly. In order to show the application effect of the BHD filter, we tested it on gridding gravity anomaly. Fig. 2a shows the synthetic gravity anomaly generated by two prisms with top depths of 15 and 20 m by adding random noise with mean square error (MSE) of 1 mGal, and the white dotted lines represent the true horizontal locations of the sources. Fig. 2b shows the THD of the data in Fig. 2a, and the THD cannot display the edges of deep source clearly. Fig. 2c and d shows the Theta map and TDX of the data in Fig. 2a respectively, and they can display the edges of
14
G. Ma et al. / Journal of Applied Geophysics 108 (2014) 12–18
Fig. 1. (a) Synthetic gravity anomaly generated by two prisms with top depths of 15 m and 20 m. (b) Total horizontal derivative of the data in (a). (c) Theta map of the data in (a). (d) Balanced horizontal derivative of the data in (a). (e) BHD of the magnetic anomaly generated by the model with an inclination of 60°.(f) BHD of the magnetic anomaly generated by the model with an inclination of 90°.
the sources clearly, but the recognized edges are diffused and bigger than the true edges. Fig. 2e shows the absolute value of the BHD of the data in Fig. 2a, and the BHD filter can display the source edges more clearly and precisely, and the edges recognized by BHD filter are not destroyed seriously because of the computation of the higher-order derivatives. We can see that the higher-order derivative can obtain more accurate results, and the second-order normalized horizontal tilt angle (TDX2) can be given by 0qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi1 2 2 −1 B ð∂f z =∂xÞ þ ð∂f z =∂yÞ C TDX 2 ¼ tan @ A 2 ∂ f =∂z2
ð13Þ
where, fz is the first-order vertical derivative of potential field f. Fig. 2f shows the second-order TDX of the data in Fig. 2a. There are lower values surrounding the recognized edges, and the results are sensitive to noise, which are unwanted and they lead to complicate the interpretation. The depth estimation is another main goal in the interpretation of potential field data. We derive a linear equation based on the derivatives of the BHD to estimate the horizontal location and the depth of the
source. The conventional Euler equation can be given by (Thompson, 1982; Reid et al., 1990) ðx−x0 Þ
∂f ∂f ∂f þ ðy−y0 Þ þ ðz−z0 Þ ¼ −Nð f −BÞ ∂x ∂y ∂z
ð14Þ
where, x, y and z are the observation coordinates,x0, y0 and z0 are the coordinates of the source,B is the background field, and N is the structural index of the source, which represents the nature of the source (Stavrev and Reid, 2007). Computing the derivatives of Eq. (14) in the x and y directions and assuming a constant background B, we obtain ðx−x0 Þ
∂2 f ∂2 f ∂2 f ∂f þ ðy−y0 Þ þ ðz−z0 Þ ¼ −ðN þ 1Þ 2 ∂x∂y ∂x∂z ∂x ∂x
ð15Þ
∂2 f ∂2 f ∂2 f ∂f þ ðy−y0 Þ 2 þ ðz−z0 Þ ¼ −ðN þ 1Þ : ∂x∂y ∂y∂z ∂y ∂y
ð16Þ
and
ðx−x0 Þ
G. Ma et al. / Journal of Applied Geophysics 108 (2014) 12–18
15
Fig. 2. (a)Synthetic gravity anomaly generated by two prisms with top depths of 10 m and 15 m by adding random noise. (b) Total horizontal derivative of the data in (a). (c) Theta map of the data in (a). (d) TDX of the data in (a). (e) Balanced horizontal derivative of the data in (a). (f) Second-order TDX of the data in (a).
Computing the second-order vertical derivative of Eq. (14), we obtain 3
ðx−x0 Þ
3
3
ð20Þ
2
∂ f ∂ f ∂ f ∂ f þ ðy−y0 Þ þ ðz−z0 Þ 3 ¼ −ðN þ 2Þ 2 : ∂x∂z2 ∂y∂z2 ∂z ∂z
Multiplying Eqs. (15), (16) and (17) by
! ! ∂BHD 1 ∂2 f ∂2 f ∂f ∂2 f ∂f ∂3 f −THD ¼ Ak þ THD ∂z2 ∂y ∂x∂y ∂x ∂y2 ∂y ∂y∂z2
∂ f ∂2 f ∂ f ∂2 f ∂x ∂z2 ∂y ∂z2
,
THD THD
ð17Þ
and THD respectively,
and adding the first two and then subtracting the third one, we get ! # 1 ∂2 f ∂2 f ∂f ∂2 f ∂f ∂3 f −THD þ THD ∂z2 ∂x2 ∂x ∂x∂y ∂y ∂x∂z2 ! " # 1 ∂2 f ∂2 f ∂f ∂2 f ∂f ∂3 f −THD þ þ ðy−y0 Þ THD ∂z2 ∂x∂y ∂x ∂y2 ∂y ∂y∂z2 ! " # 2 2 2 1 ∂ f ∂ f ∂f ∂ f ∂f ∂3 f −THD 3 þ ðz−z0 Þ þ THD ∂z2 ∂x∂z ∂x ∂y∂z ∂y ∂z "
and ! ! ∂BHD 1 ∂2 f ∂2 f ∂f ∂2 f ∂f ∂3 f −THD ¼ Ak þ THD ∂z2 ∂z ∂x∂z ∂x ∂y∂z ∂y ∂z3 where, Ak ¼
ðx−x0 Þ
k 2. ð∂ f =∂xÞ2 þð∂ f =∂yÞ2 þðk∂2 f =∂z2 Þ
ð21Þ
Substituting Eqs. (19), (20) and
(21) into Eq. (18), we can obtain
ð18Þ
ðx−x0 Þ
∂BHD ∂BHD ∂BHD þ ðy−y0 Þ þ ðz−z0 Þ ¼ −Ak THD ∂x ∂y ∂z
! ∂2 f ∂2 f : þ ∂x2 ∂y2
ð22Þ
2
¼ THD
∂ f : ∂z2
Computing the derivatives of the BHD in the x, y and z directions, we get ∂BHD 1 ∂2 f ¼ Ak THD ∂z2 ∂x
! ! ∂2 f ∂f ∂2 f ∂f ∂3 f −THD þ ∂x2 ∂x ∂x∂y ∂y ∂x∂z2
ð19Þ
The location parameters x0, y0, and z0 of the source can be obtained by solving Eq. (22). We also use a clustering method (Ma, 2013a,b) to get better results. The inversion results of conventional Euler deconvolution method depend on the difference between the true structural index and the given structural index (FitzGerald et al, 2004; Barbosa et al, 1999), but the structural indexes of an area are hard to know, and Eq. (22) can obtain the location parameters without any a priori information about the source, so the inversion results computed by Eq. (22) are more stable and correct.
16
G. Ma et al. / Journal of Applied Geophysics 108 (2014) 12–18
We also present the normalized total horizontal derivative (NTHD) method to obtain the depth of the source, which can be given by NTHDðx; y; zÞ ¼
THDðx; y; zÞ THDðx; y; zÞ ¼ M T ðzÞ 1X THDðx; y; zÞ M 1
ð23Þ
where, M is the number of total observation points. Depending on the theory of the normalized full gradient (NFG) method (Griffin, 1949; Berezkin, 1967; Aydin, 1997; Zeng et al, 2002; Aydin, 2005, 2007) we know that horizontal location and the depth of the source can be obtained depending on the closed maxima of the NTHD. Fig. 3a shows the gravity anomaly generated by two prisms with top depths of 10 and 13 m. Fig. 3b shows the absolute value of the BHD of the data in Fig. 3a, and we can see that the BHD filter can display the source edges clearly. Fig. 3c shows the estimated location and depth of the source by Eq. (22), and the black lines represent the true locations of the sources. We can see that the estimated depths are close to real depths. Fig. 3d shows the NTHD of the data computed by Eq. (23), and the maxima of the NTHD can identify the horizontal location and depth of the source correctly. We can see that the maxima of total horizontal derivative are totally consistent with the true edges of the source when the top depth of the source is zero depending on Eq. (6), so the horizontal location identified by the NTHD is most accurate. Fig. 4a shows the synthetic gravity anomaly generated by two prisms with top depths of 10 and 15 m, respectively. Fig. 4b shows the THD of the data in Fig. 4a, and the THD cannot display the edges of the deeper body clearly. Fig. 4c and d shows the Theta map and TDX of the data in Fig. 4a, and they can display the edges of shallow and deep bodies simultaneously, but the recognized edges are diffused. Fig. 4e shows the absolute value of the BHD of the data in Fig. 4a, and the recognized edges are more accurate and clear compared to other edge detection filters. Fig. 4f shows the NTHD of the data in Fig. 4a computed by Eq. (23), and we can see that the top depths of the sources are 10 and 15 m depending on the closed maxima of the results, respectively. Fig. 4g and h shows
the estimated horizontal locations and the depths of the sources computed by Eq. (22). The estimated horizontal locations are consistent with the true locations of the sources, and the estimated depths are 9.8 and 14.7 m, which are close to the true values. 3. Application to real magnetic data Fig. 5a shows the reduction to the pole of the magnetic data in Sanjiang area, Northeast China, and the interval of the point is 100 m. Fig. 5b shows the THD of the data in Fig. 5a, and Fig. 5c shows the absolute value of the BHD of the data in Fig. 5a. The BHD filter can display the edges of the stratums more clearly, and show more details. Fig. 5d shows the NTHD of the data in Fig. 5a, and we can deduce that the depths of the sources are about at the range of 130–200 m depending on the results. Fig. 5e shows the estimated location parameters of the sources by Eq. (22), and the depths of the sources are mostly at the range of 120–200 m, which is consistent with the results of the NTHD. Fig. 5f shows the geology map of Sanjiang area. The Neoproterozoic Supersequence consists predominantly of continental sediments, except for some iron-rich SE–NW trending sandstone dykes. As can be seen from the results the horizontal locations of the sources obtained by the proposed methods are in accord with the geology information, so the proposed methods have better practicality. 4. Conclusions We present the balanced horizontal derivative (BHD) filter, which uses the combination of the first-order horizontal derivative and the second-order vertical derivatives to recognize source edges. The BHD filter can recognize the edges more precisely and clearly compared to the existing edge detection filters. We derive a linear equation based on the derivatives of the BHD to estimate the horizontal location and depth of the source without any a priori information about the source, and we also present a normalized total horizontal derivative (NTHD) method to show the depth of the source. We demonstrate the proposed methods on synthetic potential field data, and the recognized edges are more precise and clear,
Fig. 3. (a) Synthetic gravity anomaly generated by two prisms with top depths of 10 and 13 m. (b) Balanced horizontal derivative of the data in (a). (c) Horizontal location and depth of the source estimated by Eq. (22). (d) Normalized total horizontal derivative (NTHD) of the data in (a).
G. Ma et al. / Journal of Applied Geophysics 108 (2014) 12–18
17
Fig. 4. (a)Synthetic gravity anomaly generated by two prisms with top depths of 10 m and 15 m. (b) Total horizontal derivative of the data in (a). (c) Theta map of the data in (a). (d) TDX of the data in (a). (e) Balanced horizontal derivative of the data in (a). (f) Normalized total horizontal derivative (NTHD) of the data in (a). (g) Estimated horizontal locations of the sources computed by Eq. (22). (h) Estimated depths of the sources computed by Eq. (22).
and the estimated depths by the presented methods are consistent with the true values. We also apply the proposed methods to real potential field data, and the inversion results are in accord with the geology information.
Acknowledgment This research is supported by the China Postdoctoral Science Foundation (2014M550173) and deep exploration technique and experiment
18
G. Ma et al. / Journal of Applied Geophysics 108 (2014) 12–18
Fig. 5. (a) Reduction to the pole of magnetic anomaly in Sanjiang area; (b) Total horizontal derivative of the data in 5a; (c) Balanced horizontal derivative of the data in 5a. (d) Normalized total horizontal derivative (NTHD) of the data in 5a. (e) Estimated location parameters by the Eq. (22). (f) Geology information map of Sanjiang area.
study SinoProbe-09-01(201011078). We sincerely acknowledge the assistance of the associate editors and the reviewers in improving the paper. References Aydin, A., 1997. Evaluation of Gravity Data in Terms of Hydrocarbon by Normalized Full Gradient, Variation and Statistic Methods, Model Studies and Application in Hasankale-Horasan Basin (Erzurum). (Ph.D. Thesis) Karadeniz Technical Univ., Natural and Applied Sciences Institute, Trabzon, Turkey. Aydin, A., 2005. Evaluation of gravity anomalies by direct interpretation techniques: an application from Hasankale-Horasan region. J. Eng. Sci. 11, 95–102. Aydin, A., 2007. Interpretation of gravity anomalies with the normalized full gradient (NFG) method and an example. Pure Appl. Geophys. 164, 2329–2344. Barbosa, V.C.F., Silva, J.B.C., Medeiros, W.E., 1999. Stability analysis and improvement of structural index estimation in Euler deconvolution. Geophysics 64, 48–60. Berezkin, V.M., 1967. Application of the total vertical gradient of gravity for determination of the depths to the sources of gravity anomalies. Explor. Geophys. 18, 69–79. Blakely, R.J., 1995. Potential Theory in Gravity and Magnetic Application. Cambridge University Press. Cooper, G.R.J., Cowan, D.R., 2006. Enhancing potential field data using filters based on the local phase. Comput. Geosci. 32, 1585–1591. Cooper, G.R.J., Cowan, D.R., 2008. Edge enhancement of potential-field data using normalized statistics. Geophysics 73, H1–H4. Cordell, L., 1979. Gravimetric expression of graben faulting in Santa Fe Country and the Espanola Basin, New Mexico. New Mexico Geol. Soc. Guidebook, 30th Field Conf, pp. 59–64. Cordell, L., Grauch, V.J.S., 1985. Mapping basement magnetization zones from aeromagnetic data in the San Juan Basin, New Mexico. In: Hinzc, W.J. (Ed.), The Utility of Regional Gravity and Magnetic Anomaly. Society of Exploration Geophysicists, pp. 181–197. Even, H.M., 1936. The place of the vertical gradient in gravitational interpretations. Geophysics 1, 127–136.
Fedi, M., Florio, G., 2001. Detection of potential fields source boundaries by enhanced horizontal derivative method. Geophys. Prospect. 49, 40–58. FitzGerald, D., Reid, A., McInerney, P., 2004. New discrimination techniques for Euler deconvolution. Comput. Geosci. 30, 461–469. Griffin, W.R., 1949. Residual gravity in theory and practice. Geophysics 14, 39–58. Hood, P.J., Tasked, D.J., 1989. Aeromagnetic gradiometer program of the Geological Survey of Canada. Geophysics 54, 1012–1022. Hsu, S., Sibuet, J.C., Shyu, C., 1996. High-resolution detection of geologic boundaries from potential field anomalies: an enhanced analytic signal technique. Geophysics 61, 373–386. Ma, G., 2013a. Edge detection of potential field data using improved local phase filter. Explor. Geophys. 44, 36–41. Ma, G., 2013b. Combination of horizontal gradient ratio and Euler (HGR-EUL) methods for the interpretation of potential field data. Geophysics 78, J53–J60. Ma, G., Li, L., 2012. Edge detection in potential fields with the normalized total horizontal derivative. Comput. Geosci. 41, 83–87. Miller, H.G., Singh, V., 1994. Potential field tilt—a new concept for location of potential field sources. J. Appl. Geophys. 32, 213–217. Rajagopalan, S., Milligan, P., 1995. Image enhancement of aeromagnetic data using automatic gain control. Explor. Geophys. 25, 173–178. Reid, A., Allsop, J.M., Granser, H., et al., 1990. Magnetic interpretation in three dimensions using Euler deconvolution. Geophysics 55, 80–91. Roest, W.R., Verhoef, J., Pilkington, M., 1992. Magnetic interpretation using the 3-D analytic signal. Geophysics 57, 116–125. Stavrev, P., Reid, A., 2007. Degrees of homogeneity of potential fields and structural indices of Euler deconvolution. Geophysics 72, L1–L12. Thompson, D.T., 1982. EULDPH: a new technique for making computer-assisted depth estimates from magnetic data. Geophysics 47, 31–37. Verduzco, B., Fairhead, J.D., Green, C.M., 2004. New insights into magnetic derivatives for structural mapping. Lead. Edge 23, 116–119. Wijns, C., Perez, C., Kowalczyk, P., 2005. Theta map: edge detection in magnetic data. Geophysics 70, 39–43. Zeng, H., Meng, X., Yao, C., et al., 2002. Detection of reservoirs from normalized full gradient of gravity anomalies and its application to Shengli Oil Field, East China. Geophysics 67, 1138–1147.
Contents lists available at ScienceDirect
Journal of Applied Geophysics journal homepage: www.elsevier.com/locate/jappgeo
Balanced horizontal derivative of potential field data to recognize the edges and estimate location parameters of the source Guoqing Ma ⁎, Cai Liu, Lili Li College of Geoexploration Science and Technology, Jilin University, Changchun 130061, China
a r t i c l e
i n f o
Article history: Received 15 May 2014 Accepted 11 June 2014 Available online 18 June 2014 Keywords: Balanced horizontal derivative Edge detection Linear equation Depth Normalized
a b s t r a c t Edge detection is a common task in the interpretation of potential field data, and many edge detection filters are presented to accomplish this task, which are the functions composed by the first-order horizontal derivative and vertical derivative of potential field data, but we find that the edges recognized by the existing edge detection filters are bigger than the true edges. In this paper, we present the balanced horizontal derivative (BHD) edge detection filter, which uses the ratio of the first-order horizontal derivative to the second-order horizontal derivative to recognize the edges of the source, and the recognized edges by the BHD edge detection filter are more correct and are more insensitive to noise. We derive a linear equation based on the derivatives of the BHD to estimate the depth of the source, and we also present a normalized total horizontal derivative (NTHD) method to image the location of the source. We demonstrate the presented methods on synthetic potential field data, and the results show that the presented methods can provide the edges and location parameters of the sources correctly, and the BHD filter can display the edges more clearly and correctly. At last, we apply the presented methods to real potential field data, and the inversion results computed by the presented methods are in accord with the geology information. © 2014 Elsevier B.V. All rights reserved.
1. Introduction Edge detection is a commonly requested task in the interpretation of potential field data, and the maxima of total horizontal derivative and the zero of the vertical derivative are corresponding to the edges of the source (Even, 1936; Hood and Tasked, 1989; Cordell, 1979), so the edge detection filters are all the functions based on the horizontal derivative and vertical derivatives of the data (Blakely, 1995; Cordell and Grauch, 1985; Roest et al., 1992; Hsu et al., 1996; Fedi and Florio, 2001). Initially, people directly used the total horizontal derivative or vertical derivative to recognize the source edges (Even, 1936; Cordell, 1979), but these methods cannot display the edges of deeper source clearly, so many people began to develop the balanced edge detection filters. Miller and Singh (1994) presented the first balanced edge detection filter-tilt angle, which is the ratio of the vertical derivative to total horizontal derivative, and this filter can effectively balance the amplitudes of the anomalies generated by the sources with different depths. Some authors also presented other format of balanced edge detection filters that can display the edges of the shallow and deep sources simultaneously (Miller and Singh, 1994; Rajagopalan and Milligan, 1995; Verduzco et al., 2004; Wijns et al., 2005; Cooper and Cowan, 2006; 2008; Ma and Li; 2012; Ma, 2013a,b).
We find that the horizontal coordinates of the maxima of total horizontal derivative and the zero of vertical derivative are both bigger than the true locations of the edges, and the coordinates of the maxima of total horizontal derivative are more close to the true edges. We prove that the ratio of the total horizontal derivative to the second-order vertical derivative can obtain more accurate edge detection results. We also introduce a Laplace equation that uses the sum of two second-order horizontal derivatives to compute the second-order vertical derivative, which cannot increase the interference of noise dramatically, and this method is called as balanced horizontal derivative (BHD) filter. The depth estimation is another important task in the interpretation of potential field data. We derive a linear equation based on the derivatives of the BHD, which can obtain the horizontal location and depth of the source without any a priori knowledge. We also present a normalized total horizontal derivative (NTHD) method, and the location and the depth of the source can be obtained depending on the closed maxima of the results. We demonstrate the presented methods on synthetic and real potential field data, and the inversion results are consistent with the real values.
2. Balanced horizontal derivate (BHD) filter ⁎ Corresponding author at. College of Geoexploration Science and Technology, Jilin University, Changchun, China. E-mail address: [email protected] (G. Ma).
http://dx.doi.org/10.1016/j.jappgeo.2014.06.005 0926-9851/© 2014 Elsevier B.V. All rights reserved.
Total horizontal derivative, Theta map and normalized horizontal tilt angle filter are three commonly used edge detection filters (Cooper
G. Ma et al. / Journal of Applied Geophysics 108 (2014) 12–18
and Cowan, 2006). Total horizontal derivative (THD) can be given by s ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 ∂f 2 ∂f : þ THD ¼ ∂x ∂y
The expression of second-order vertical derivative of the 2-D vertical prism is ð1Þ
"
V zzz
# xþa xþa x−a x−a : ¼ 2Gρ − − þ ðx þ aÞ2 þ h2 ðx þ aÞ2 þ H 2 ðx−aÞ2 þ h2 ðx−aÞ2 þ H 2
ð9Þ
The Theta map can be written by r ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2ffi ∂f 2 þ ∂∂yf ∂x ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Theta ¼ r 2 2 : ∂f 2 ∂f þ þ ∂∂zf ∂x ∂y
The coordinates that make the second-order vertical derivative get zero are ð2Þ x0;2
The normalized horizontal tilt angle (TDX) filter can be expressed as 0qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi1 2 2 −1 @ ð∂f =∂xÞ þ ð∂f =∂yÞ A TDX ¼ tan : ∂f =∂z
ð3Þ
The maxima of the above three edge detection filters automatically delineate the edges of the sources (Cooper and Cowan, 2008; Ma, 2013a,b). Depending on Eqs. (2) and (3) we can find that the Theta and TDX filters get the maxima when the first-order vertical derivative is equal to zero, which is not related to the value of horizontal derivative. The prism is a most common model in the interpretation of real data, and the expression of the first-order vertical derivative of 2-D vertical prism (Blakely, 1995) can be given by V zz
h −1 H −1 −1 H −1 h ¼ 2Gρ tan − tan − tan þ tan ð4Þ xþa xþa x−a x−a
where, 2a is the width of the prism, and h and H are the top and bottom depths of the prism, respectively. The coordinates that make the vertical derivative get zero are pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi x0;1 ¼ a2 þ Hh:
ð5Þ
We can find that the edges recognized by the Theta map and TDX filter are related to the depth of the source, and both are bigger than the true edges. The expressions of the horizontal vertical derivatives of 2-D vertical prism can be given by h
V zx
ih i 2 2 2 2 ðx þ aÞ þ H ðx−aÞ þ h : ¼ Gρ ln ðx þ aÞ2 þ h2 ðx−aÞ2 þ h2
ð6Þ
The coordinates that make the horizontal derivative getthe maxima are ffi sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 1 2 2 2 2 2 2 2 þ 12h H þ 2a − h þ H xm ¼ : 4a þ h þ H 6
ð7Þ
The relationship between the maxima of the horizontal derivative and the zero of vertical derivative is x Njx jNa: 0;1 m
13
ð8Þ
We can see that the coordinates of the maxima of horizontal derivative are bigger than the true edges, but are smaller than the coordinates of the zero of vertical derivative, so the edges recognized by the total horizontal derivative are more accurate, but the total horizontal derivative cannot give the edges of the deeper sources clearly.
ffi sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 1 2 2 2 2 2 2 2 ¼ þ 12h H þ 2a − h þ H 4a þ h þ H : ð10Þ 6
The coordinates of the zero of second-order vertical derivative are the same as the coordinates of the maxima of total horizontal derivative, so we suggest using the ratio of total horizontal derivative to the second-order vertical derivative to recognize the source edges. In order to lower the noise effect, we use the Laplace equation (Blakely, 1995) to compute the second order vertical derivative. 2
2
2
∂ f ∂ f ∂ f ¼− þ ∂z2 ∂x2 ∂y2
! ð11Þ
The second-order horizontal derivative can be computed in the space domain. This method is called as balanced horizontal derivative (BHD) edge detection filter, which can be given by qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 ð∂f =∂xÞ2 þ ð∂f =∂yÞ2 A k − ∂2 f =∂x2 þ ∂2 f =∂y2
0 −1 @
BHD ¼ tan
ð12Þ
where, ∂f/∂x, ∂f/∂y and ∂2f/∂z2 are the derivatives of the data f, and meanð∂ f =∂zÞ , and the maxima of the absolute value of the BHD are k ¼ mean ð∂2 f =∂z2 Þ corresponding to the edges of the source. We give a series of gravity profiles to show the distribution of different edge detection filters. Fig. 1a shows the synthetic gravity anomaly generated by two prisms with top depths of 15 and 20 m, respectively. Fig. 1b shows the total horizontal derivative (THD) of the data in Fig. 1a, and the THD can recognize the source edges, but cannot display the edges of the deep sources clearly. The results of the Theta map and TDX are similar, so we only give the Theta map of the data, and the Theta map of the data is shown in Fig. 1c. The Theta map can display the edges of shallow and deep sources simultaneously, but the recognized edges are significantly bigger than the true edges. Fig. 1d shows the absolute value of the BHD of the data in Fig. 1a. As can be seen from the results, the edges recognized by the THD and BHD are more accurate, and the BHD filter can display the edges of the shallow and deep sources simultaneously. Fig. 1e shows the absolute value of the BHD of the magnetic anomaly generated by the prism with an inclination of 60°, and the recognized edges have big difference with the real edges. Fig. 1f shows the absolute value of the BHD of the magnetic anomaly generated by the prism with an inclination of 90°, and the recognized edges are consistent with the true edges. We can see that the horizontal derivatives of magnetic anomaly are sensitive to magnetization direction, so we should reduce the magnetic anomaly to the North Pole when using the BHD to interpret the magnetic anomaly. In order to show the application effect of the BHD filter, we tested it on gridding gravity anomaly. Fig. 2a shows the synthetic gravity anomaly generated by two prisms with top depths of 15 and 20 m by adding random noise with mean square error (MSE) of 1 mGal, and the white dotted lines represent the true horizontal locations of the sources. Fig. 2b shows the THD of the data in Fig. 2a, and the THD cannot display the edges of deep source clearly. Fig. 2c and d shows the Theta map and TDX of the data in Fig. 2a respectively, and they can display the edges of
14
G. Ma et al. / Journal of Applied Geophysics 108 (2014) 12–18
Fig. 1. (a) Synthetic gravity anomaly generated by two prisms with top depths of 15 m and 20 m. (b) Total horizontal derivative of the data in (a). (c) Theta map of the data in (a). (d) Balanced horizontal derivative of the data in (a). (e) BHD of the magnetic anomaly generated by the model with an inclination of 60°.(f) BHD of the magnetic anomaly generated by the model with an inclination of 90°.
the sources clearly, but the recognized edges are diffused and bigger than the true edges. Fig. 2e shows the absolute value of the BHD of the data in Fig. 2a, and the BHD filter can display the source edges more clearly and precisely, and the edges recognized by BHD filter are not destroyed seriously because of the computation of the higher-order derivatives. We can see that the higher-order derivative can obtain more accurate results, and the second-order normalized horizontal tilt angle (TDX2) can be given by 0qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi1 2 2 −1 B ð∂f z =∂xÞ þ ð∂f z =∂yÞ C TDX 2 ¼ tan @ A 2 ∂ f =∂z2
ð13Þ
where, fz is the first-order vertical derivative of potential field f. Fig. 2f shows the second-order TDX of the data in Fig. 2a. There are lower values surrounding the recognized edges, and the results are sensitive to noise, which are unwanted and they lead to complicate the interpretation. The depth estimation is another main goal in the interpretation of potential field data. We derive a linear equation based on the derivatives of the BHD to estimate the horizontal location and the depth of the
source. The conventional Euler equation can be given by (Thompson, 1982; Reid et al., 1990) ðx−x0 Þ
∂f ∂f ∂f þ ðy−y0 Þ þ ðz−z0 Þ ¼ −Nð f −BÞ ∂x ∂y ∂z
ð14Þ
where, x, y and z are the observation coordinates,x0, y0 and z0 are the coordinates of the source,B is the background field, and N is the structural index of the source, which represents the nature of the source (Stavrev and Reid, 2007). Computing the derivatives of Eq. (14) in the x and y directions and assuming a constant background B, we obtain ðx−x0 Þ
∂2 f ∂2 f ∂2 f ∂f þ ðy−y0 Þ þ ðz−z0 Þ ¼ −ðN þ 1Þ 2 ∂x∂y ∂x∂z ∂x ∂x
ð15Þ
∂2 f ∂2 f ∂2 f ∂f þ ðy−y0 Þ 2 þ ðz−z0 Þ ¼ −ðN þ 1Þ : ∂x∂y ∂y∂z ∂y ∂y
ð16Þ
and
ðx−x0 Þ
G. Ma et al. / Journal of Applied Geophysics 108 (2014) 12–18
15
Fig. 2. (a)Synthetic gravity anomaly generated by two prisms with top depths of 10 m and 15 m by adding random noise. (b) Total horizontal derivative of the data in (a). (c) Theta map of the data in (a). (d) TDX of the data in (a). (e) Balanced horizontal derivative of the data in (a). (f) Second-order TDX of the data in (a).
Computing the second-order vertical derivative of Eq. (14), we obtain 3
ðx−x0 Þ
3
3
ð20Þ
2
∂ f ∂ f ∂ f ∂ f þ ðy−y0 Þ þ ðz−z0 Þ 3 ¼ −ðN þ 2Þ 2 : ∂x∂z2 ∂y∂z2 ∂z ∂z
Multiplying Eqs. (15), (16) and (17) by
! ! ∂BHD 1 ∂2 f ∂2 f ∂f ∂2 f ∂f ∂3 f −THD ¼ Ak þ THD ∂z2 ∂y ∂x∂y ∂x ∂y2 ∂y ∂y∂z2
∂ f ∂2 f ∂ f ∂2 f ∂x ∂z2 ∂y ∂z2
,
THD THD
ð17Þ
and THD respectively,
and adding the first two and then subtracting the third one, we get ! # 1 ∂2 f ∂2 f ∂f ∂2 f ∂f ∂3 f −THD þ THD ∂z2 ∂x2 ∂x ∂x∂y ∂y ∂x∂z2 ! " # 1 ∂2 f ∂2 f ∂f ∂2 f ∂f ∂3 f −THD þ þ ðy−y0 Þ THD ∂z2 ∂x∂y ∂x ∂y2 ∂y ∂y∂z2 ! " # 2 2 2 1 ∂ f ∂ f ∂f ∂ f ∂f ∂3 f −THD 3 þ ðz−z0 Þ þ THD ∂z2 ∂x∂z ∂x ∂y∂z ∂y ∂z "
and ! ! ∂BHD 1 ∂2 f ∂2 f ∂f ∂2 f ∂f ∂3 f −THD ¼ Ak þ THD ∂z2 ∂z ∂x∂z ∂x ∂y∂z ∂y ∂z3 where, Ak ¼
ðx−x0 Þ
k 2. ð∂ f =∂xÞ2 þð∂ f =∂yÞ2 þðk∂2 f =∂z2 Þ
ð21Þ
Substituting Eqs. (19), (20) and
(21) into Eq. (18), we can obtain
ð18Þ
ðx−x0 Þ
∂BHD ∂BHD ∂BHD þ ðy−y0 Þ þ ðz−z0 Þ ¼ −Ak THD ∂x ∂y ∂z
! ∂2 f ∂2 f : þ ∂x2 ∂y2
ð22Þ
2
¼ THD
∂ f : ∂z2
Computing the derivatives of the BHD in the x, y and z directions, we get ∂BHD 1 ∂2 f ¼ Ak THD ∂z2 ∂x
! ! ∂2 f ∂f ∂2 f ∂f ∂3 f −THD þ ∂x2 ∂x ∂x∂y ∂y ∂x∂z2
ð19Þ
The location parameters x0, y0, and z0 of the source can be obtained by solving Eq. (22). We also use a clustering method (Ma, 2013a,b) to get better results. The inversion results of conventional Euler deconvolution method depend on the difference between the true structural index and the given structural index (FitzGerald et al, 2004; Barbosa et al, 1999), but the structural indexes of an area are hard to know, and Eq. (22) can obtain the location parameters without any a priori information about the source, so the inversion results computed by Eq. (22) are more stable and correct.
16
G. Ma et al. / Journal of Applied Geophysics 108 (2014) 12–18
We also present the normalized total horizontal derivative (NTHD) method to obtain the depth of the source, which can be given by NTHDðx; y; zÞ ¼
THDðx; y; zÞ THDðx; y; zÞ ¼ M T ðzÞ 1X THDðx; y; zÞ M 1
ð23Þ
where, M is the number of total observation points. Depending on the theory of the normalized full gradient (NFG) method (Griffin, 1949; Berezkin, 1967; Aydin, 1997; Zeng et al, 2002; Aydin, 2005, 2007) we know that horizontal location and the depth of the source can be obtained depending on the closed maxima of the NTHD. Fig. 3a shows the gravity anomaly generated by two prisms with top depths of 10 and 13 m. Fig. 3b shows the absolute value of the BHD of the data in Fig. 3a, and we can see that the BHD filter can display the source edges clearly. Fig. 3c shows the estimated location and depth of the source by Eq. (22), and the black lines represent the true locations of the sources. We can see that the estimated depths are close to real depths. Fig. 3d shows the NTHD of the data computed by Eq. (23), and the maxima of the NTHD can identify the horizontal location and depth of the source correctly. We can see that the maxima of total horizontal derivative are totally consistent with the true edges of the source when the top depth of the source is zero depending on Eq. (6), so the horizontal location identified by the NTHD is most accurate. Fig. 4a shows the synthetic gravity anomaly generated by two prisms with top depths of 10 and 15 m, respectively. Fig. 4b shows the THD of the data in Fig. 4a, and the THD cannot display the edges of the deeper body clearly. Fig. 4c and d shows the Theta map and TDX of the data in Fig. 4a, and they can display the edges of shallow and deep bodies simultaneously, but the recognized edges are diffused. Fig. 4e shows the absolute value of the BHD of the data in Fig. 4a, and the recognized edges are more accurate and clear compared to other edge detection filters. Fig. 4f shows the NTHD of the data in Fig. 4a computed by Eq. (23), and we can see that the top depths of the sources are 10 and 15 m depending on the closed maxima of the results, respectively. Fig. 4g and h shows
the estimated horizontal locations and the depths of the sources computed by Eq. (22). The estimated horizontal locations are consistent with the true locations of the sources, and the estimated depths are 9.8 and 14.7 m, which are close to the true values. 3. Application to real magnetic data Fig. 5a shows the reduction to the pole of the magnetic data in Sanjiang area, Northeast China, and the interval of the point is 100 m. Fig. 5b shows the THD of the data in Fig. 5a, and Fig. 5c shows the absolute value of the BHD of the data in Fig. 5a. The BHD filter can display the edges of the stratums more clearly, and show more details. Fig. 5d shows the NTHD of the data in Fig. 5a, and we can deduce that the depths of the sources are about at the range of 130–200 m depending on the results. Fig. 5e shows the estimated location parameters of the sources by Eq. (22), and the depths of the sources are mostly at the range of 120–200 m, which is consistent with the results of the NTHD. Fig. 5f shows the geology map of Sanjiang area. The Neoproterozoic Supersequence consists predominantly of continental sediments, except for some iron-rich SE–NW trending sandstone dykes. As can be seen from the results the horizontal locations of the sources obtained by the proposed methods are in accord with the geology information, so the proposed methods have better practicality. 4. Conclusions We present the balanced horizontal derivative (BHD) filter, which uses the combination of the first-order horizontal derivative and the second-order vertical derivatives to recognize source edges. The BHD filter can recognize the edges more precisely and clearly compared to the existing edge detection filters. We derive a linear equation based on the derivatives of the BHD to estimate the horizontal location and depth of the source without any a priori information about the source, and we also present a normalized total horizontal derivative (NTHD) method to show the depth of the source. We demonstrate the proposed methods on synthetic potential field data, and the recognized edges are more precise and clear,
Fig. 3. (a) Synthetic gravity anomaly generated by two prisms with top depths of 10 and 13 m. (b) Balanced horizontal derivative of the data in (a). (c) Horizontal location and depth of the source estimated by Eq. (22). (d) Normalized total horizontal derivative (NTHD) of the data in (a).
G. Ma et al. / Journal of Applied Geophysics 108 (2014) 12–18
17
Fig. 4. (a)Synthetic gravity anomaly generated by two prisms with top depths of 10 m and 15 m. (b) Total horizontal derivative of the data in (a). (c) Theta map of the data in (a). (d) TDX of the data in (a). (e) Balanced horizontal derivative of the data in (a). (f) Normalized total horizontal derivative (NTHD) of the data in (a). (g) Estimated horizontal locations of the sources computed by Eq. (22). (h) Estimated depths of the sources computed by Eq. (22).
and the estimated depths by the presented methods are consistent with the true values. We also apply the proposed methods to real potential field data, and the inversion results are in accord with the geology information.
Acknowledgment This research is supported by the China Postdoctoral Science Foundation (2014M550173) and deep exploration technique and experiment
18
G. Ma et al. / Journal of Applied Geophysics 108 (2014) 12–18
Fig. 5. (a) Reduction to the pole of magnetic anomaly in Sanjiang area; (b) Total horizontal derivative of the data in 5a; (c) Balanced horizontal derivative of the data in 5a. (d) Normalized total horizontal derivative (NTHD) of the data in 5a. (e) Estimated location parameters by the Eq. (22). (f) Geology information map of Sanjiang area.
study SinoProbe-09-01(201011078). We sincerely acknowledge the assistance of the associate editors and the reviewers in improving the paper. References Aydin, A., 1997. Evaluation of Gravity Data in Terms of Hydrocarbon by Normalized Full Gradient, Variation and Statistic Methods, Model Studies and Application in Hasankale-Horasan Basin (Erzurum). (Ph.D. Thesis) Karadeniz Technical Univ., Natural and Applied Sciences Institute, Trabzon, Turkey. Aydin, A., 2005. Evaluation of gravity anomalies by direct interpretation techniques: an application from Hasankale-Horasan region. J. Eng. Sci. 11, 95–102. Aydin, A., 2007. Interpretation of gravity anomalies with the normalized full gradient (NFG) method and an example. Pure Appl. Geophys. 164, 2329–2344. Barbosa, V.C.F., Silva, J.B.C., Medeiros, W.E., 1999. Stability analysis and improvement of structural index estimation in Euler deconvolution. Geophysics 64, 48–60. Berezkin, V.M., 1967. Application of the total vertical gradient of gravity for determination of the depths to the sources of gravity anomalies. Explor. Geophys. 18, 69–79. Blakely, R.J., 1995. Potential Theory in Gravity and Magnetic Application. Cambridge University Press. Cooper, G.R.J., Cowan, D.R., 2006. Enhancing potential field data using filters based on the local phase. Comput. Geosci. 32, 1585–1591. Cooper, G.R.J., Cowan, D.R., 2008. Edge enhancement of potential-field data using normalized statistics. Geophysics 73, H1–H4. Cordell, L., 1979. Gravimetric expression of graben faulting in Santa Fe Country and the Espanola Basin, New Mexico. New Mexico Geol. Soc. Guidebook, 30th Field Conf, pp. 59–64. Cordell, L., Grauch, V.J.S., 1985. Mapping basement magnetization zones from aeromagnetic data in the San Juan Basin, New Mexico. In: Hinzc, W.J. (Ed.), The Utility of Regional Gravity and Magnetic Anomaly. Society of Exploration Geophysicists, pp. 181–197. Even, H.M., 1936. The place of the vertical gradient in gravitational interpretations. Geophysics 1, 127–136.
Fedi, M., Florio, G., 2001. Detection of potential fields source boundaries by enhanced horizontal derivative method. Geophys. Prospect. 49, 40–58. FitzGerald, D., Reid, A., McInerney, P., 2004. New discrimination techniques for Euler deconvolution. Comput. Geosci. 30, 461–469. Griffin, W.R., 1949. Residual gravity in theory and practice. Geophysics 14, 39–58. Hood, P.J., Tasked, D.J., 1989. Aeromagnetic gradiometer program of the Geological Survey of Canada. Geophysics 54, 1012–1022. Hsu, S., Sibuet, J.C., Shyu, C., 1996. High-resolution detection of geologic boundaries from potential field anomalies: an enhanced analytic signal technique. Geophysics 61, 373–386. Ma, G., 2013a. Edge detection of potential field data using improved local phase filter. Explor. Geophys. 44, 36–41. Ma, G., 2013b. Combination of horizontal gradient ratio and Euler (HGR-EUL) methods for the interpretation of potential field data. Geophysics 78, J53–J60. Ma, G., Li, L., 2012. Edge detection in potential fields with the normalized total horizontal derivative. Comput. Geosci. 41, 83–87. Miller, H.G., Singh, V., 1994. Potential field tilt—a new concept for location of potential field sources. J. Appl. Geophys. 32, 213–217. Rajagopalan, S., Milligan, P., 1995. Image enhancement of aeromagnetic data using automatic gain control. Explor. Geophys. 25, 173–178. Reid, A., Allsop, J.M., Granser, H., et al., 1990. Magnetic interpretation in three dimensions using Euler deconvolution. Geophysics 55, 80–91. Roest, W.R., Verhoef, J., Pilkington, M., 1992. Magnetic interpretation using the 3-D analytic signal. Geophysics 57, 116–125. Stavrev, P., Reid, A., 2007. Degrees of homogeneity of potential fields and structural indices of Euler deconvolution. Geophysics 72, L1–L12. Thompson, D.T., 1982. EULDPH: a new technique for making computer-assisted depth estimates from magnetic data. Geophysics 47, 31–37. Verduzco, B., Fairhead, J.D., Green, C.M., 2004. New insights into magnetic derivatives for structural mapping. Lead. Edge 23, 116–119. Wijns, C., Perez, C., Kowalczyk, P., 2005. Theta map: edge detection in magnetic data. Geophysics 70, 39–43. Zeng, H., Meng, X., Yao, C., et al., 2002. Detection of reservoirs from normalized full gradient of gravity anomalies and its application to Shengli Oil Field, East China. Geophysics 67, 1138–1147.