Earthquakes are one of the most destructive geological disasters in the natural world, and crustal deformation and focal medium changes in the process of earthquake preparation cause variations of the gravity field near the earthquake area (Chen Yuntai et al., 1980; Chen S. et al., 2016; Chen Shi et al., 2014). Usually, signals of gravity field changes within the scope of coverage areas are obtained by continuous gravity observations or regular mobile repeated gravity observations, by which precursor information related to focal changes are likely to be captured. As one of important physical characteristics of the earth, gravity field changes have distinct geophysical significance, and can be thought of as a manifestation of tectonic movement. According to previous research, the application of the analytical method to the temporal and spatial dynamic changes of regional gravity field in the eastern margin of the Qinghai-Tibetan Plateau, Sichuan-Yunnan region and Xinjiang region provides references for the understanding on strong earthquake preparation and occurrence rules, the capture of earthquake precursor information and the determination of earthquake-prone zones (Zhu Yiqing et al., 2013). How to make better use of gravity change signals to study anomaly sources has always been an issue that has attracted much attention, and is also the basis for the construction of a quantitative index system for earthquake prediction (Chen S. et al., 2016; Chen Shi et al., 2015).
Moreover, the inversion method is an important part of the research and interpretation of field source characteristics that cause gravity changes. Earthquake precursor research based on geophysical observations usually follows the idea of seeking a source from the field and combining the field and source as guidelines. After the acquirement of field change information, further quantitative explanation on field source characteristics should be made. However, for the inversion of shock-to-shock field source parameters, the difficulty lies in the deficiency of available prior information, thus, selecting reasonable inversion techniques to study field source parameters is one of the main problems.
In this article, the three-dimensional Euler deconvolution method (Reid et al., 1990; Reid, 1995) fully developed in the field of geophysical survey is applied to studying field source parameter characteristics of gravity field changes in North China. This method is suitable for application with little available prior information, to automatically or semi-automatically determine field source locations and geometrical parameters, which certainly help us to minimize human errors in the interpretation, effectively delineate the scope of anomaly sources and calculate specific locations of anomalous bodies in order to carry out inversion and interpretation of gravity change signals. This article focuses on gravity field differential changes in North China during 2009-2014 and discusses field source characteristics of gravity changes in the region within this time frame. Firstly, gravity field changes of theoretical model are calculated through simulations, which are later inverted with the Euler deconvolution method, and then by repeated calculations, structural index and inversion parameters suitable for gravity field inversion are obtained. Finally, the three-dimensional Euler deconvolution method is applied to gravity field changes in the North China region to get its field source characteristic parameters.
Based on a full grasp of Euler deconvolution theory, we employed this method to inverting for field source characteristic parameters of gravity change signals in the North China region and carried out model experiments and real data applications. The studies in this article advanced quantitative research from field to source, and on the basis of qualitative explanation, provided quantitative explanation results, which greatly facilitated the inversion and interpretation work of mobile gravity data and provided a practical approach and a new idea of earthquake monitoring.
1 METHOD AND PRINCIPALThe three-dimensional Euler deconvolution method is suitable for application with little available prior information to automatically or semi-automatically determine field source locations, delineate the scope of structure and interpret the cause of anomaly sources (Fan Meining, 2006). This method, first put forward by Reid et al.(1990, 1995), provides the Euler homogeneous equation (formula (1)) which is independent of field source density parameters, but only related to geometrical parameters for field source locations, which simplifies the inversion. Usually by solving the Euler equation, three-dimensional spatial locations for field source locations can be determined, and estimates of parameters for different geology types can be presented.
$ \left({x - {x_0}} \right)\frac{{\partial T}}{{\partial x}} + \left({y - {y_0}} \right)\frac{{\partial T}}{{\partial y}} + \left({z - {z_0}} \right)\frac{{\partial T}}{{\partial z}} = N\left({B - T} \right) $ | (1) |
Where, (x, y, z) denotes location of observation point, (x_{0}, y_{0}, z_{0}) the field source location,
The process of inversion using three-dimensional Euler deconvolution method is as follows:
(1) Meshing the discrete data of gravity field;
(2) Calculating derivatives for potential-field anomaly in x, y and z direction;
(3) Selecting a suitable structural index N and window size. Structural index N is selected according to properties of anomalous bodies. Window size is chosen based on the scale of anomalous bodies, which should cover the scope of anomalies, and allowable error is determined according to accuracy requirements;
(4) Solving a system of equations to get parameters for field source locations and background values;
(5) Finally, adjusting parameters according to converged locations and depths of all solutions. If the inversion results are satisfactory, then calculations are completed. If not, then reselect parameters and repeat step (2) and (3) until the results are satisfied.
2 MODEL STUDYBecause the Euler deconvolution method can be directly used to invert for the field source parameters without the existing prior field source information, it has been widely applied in exploration geophysics. During inversion, a system of linear equations can be formed by selecting a suitable structural index N and determining derivatives for anomalies in the X, Y and Z direction, and then calculations are made to get solutions. The structural index, related to field-source geometric construction, is the most important parameter in Euler deconvolution. However, in practical applications, due to complex geological conditions, field source shapes are undetermined, thus the selection of the structural index usually needs to be determined by repeated calculations using various models. In addition, the selecting principles of sliding window size generally take into account the space between observation points and anomaly characteristics etc., and a window should contain enough observation points that can satisfy the equations and should cover one field-source anomaly shape. Finally, the inversion results are evaluated by getting horizontal and vertical errors of field-source location parameters and measuring if the aggregation extent of field-source location solutions matches with geological knowledge. More concentrated field-source locations indicate better accuracy and convergence of solutions, otherwise the convergence is poor.
By mastering the method above, before applying it to interpreting actual gravity field data, we designed a theoretical research model which is more practical for current observation ability, with 3 combined models of rectangular regular hexahedrons as research models, assuming that they are located respectively in the upper crust, middle crust and lower crust. Parameters for model M1, located in the upper crust, include X(80, 110) Y(150, 155) Z(5, 6), parameters for model M2 X(50, 60) Y(100, 110) Z(10, 13), located in the middle crust, parameters for M3 X(120, 150) Y(40, 70) Z(20, 25), located in the lower crust. Density change for all 3 models is assumed to be 3.0×10^{-3}g·cm^{-3}. A theoretical anomaly value on the surface observation network is obtained by forward calculation. The theoretical anomaly range calculated for this group of models is 0-120×10^{-8}m·s^{-2}, which is used to simulate observed signals for gravity changes and select preferred inversion parameters. Model anomalies are shown in Fig. 1, where the spacing between surface observation points is 5km and the lengths of observation network are both 200km in the X and Y directions. Parameters for 3 field-source models in different depths are shown in Table 1. Red rectangular frames represent their projections on the earth's surface.
The above model test results indicate that when the structural index N equals 1 and the sliding window covers a distance of 10 measuring points, the aggregation extent of inversion results is more aligned with geometric structure of model. Error of depth inversion for the combined model is limited within 15%. Fig. 2 provides a set of best inversion results obtained from calculations. Black rectangles are inverted modeling results, altogether including 564 sets of solutions, and red rectangles are theoretical locations for the combined model, which clearly show a better aggregation extent and by which field-source characteristics can be clearly distinguished. Model M1 and M3 are comparatively ideal, while inversion results for model M2 have a certain error for its horizontal position. Therefore, based on the above scientific calculations of the combined model, we select the best inversion parameters and apply them to invert the measured data of gravity changes from actual mobile gravimetric networks in North China, which have certain theoretical and guiding significance.
Since 2009, China Earthquake Administration has adjusted, optimized and reconstructed regional seismic gravity networks dispersed along major active tectonic belts in North China. 10 absolute gravity observation points and more than 120 relative gravity observation points have been newly added, measuring networks at all provincial administrations have been effectively connected, and necessary absolute gravity control has been fulfilled for measuring networks, forming an entire gravity monitoring network in North China, which makes measured data more abundant and complete.This should contribute to the mid-term prediction of strong earthquakes and great earthquakes in North China. Fig. 3(a) is the map of tectonics in observation zone in North China, where major faults are based on the Active Tectonic Map of China (Deng Qidong et al., 2007).
Figs. 3(b)-3(f) are maps of differential variations of annual gravity changes in North China during 2009-2014. Because point-by-point difference is applied to annual gravity changes using seasonal gravity data, the impact of seasonal periodical changes can be minimized. Red solid lines in the map represent positive changes of gravity value, blue solid lines are negative changes of gravity value, and black lines are the zero line. It can be seen from Figs. 3(b)-3(f) that from September 2009 to September 2010 (Fig. 3(b)), gravity field changes mainly show significant regional negative anomalous changes centering on the Shanxi fault depression belt and positive changes in the east of Shijiazhuang, spreading along the piedmont fault zone of the Taihang Mountains, where a high gradient belt of gravity changes is formed, with the biggest difference of 110×10^{-8}m/s^{2}. During September 2010-September 2011 (Fig. 3(c)), the whole Shanxi-Hebei-Inner Mongolia border area in the west of northern section of Taihang Mountains sub-block presents significant negative anomalies of gravity changes. During September 2011-September 2012 (Fig. 3(d)), the overall gravity changes demonstrate a variation trend opposite to the previous period. The overall trend of gravity changes in the northern part of the observation zone is the positive gravity change along the Shanxi fault zone. From September 2012 to September 2013(Fig. 3(e)), the gravity field overall shows positive changes, and gravity isolines mainly run in the NWW direction. Local gravity anomalous areas appeared. From September 2013 to September 2014 (Fig. 3(f)), the general trend of gravity isolines in the observation zone is cotnsistent with that of the previous period.
By using the Euler deconvolution method to invert annual gravity changes in North China duringthe 4 periods, we get field-source characteristics for gravity changes as shown in Fig. 4. Field-source locations are indicated by colored circles, and the size corresponds to field-source depth. 5 colors including yellow, pink, green, blue and purple are used respectively to represent various differential time-reversal solutions. In the process of Euler deconvolution calculation, the data gridding interval is 20km and the window size covers 8 times of grid spacing. The inversion results of field-source depths remain within 10km-50km, which is basically located in the earth's crust. According to the map of distribution characteristics of the inverted field-source locations, locations of field-source bodies that cause gravity changes show better coincidence along the Hetao fault zone, with the depths concentrating within 20km-40km, and in the strip area perpendicular to the Taihang Mountain piedmont fault, where 80% of field-sources are concentrated at depths of 20km-50km.
In this article, by the use of Euler deconvolution method, and based on theoretical model analysis, we select the optimal parameters to invert for annual gravity changes with actual measured gravity data in North China, and get the geometric parameters related to field-source locations. The main conclusions are as follows:
(1) According to inversion results in this article, we believe field source information can reasonably explain observation data, and the Euler deconvolution method is suitable for the interpretation of gravity field changes. The Euler deconvolution method can be applied to the research of seeking sources from the field without much existing prior information, which is the unique advantage of this method.
(2) On the basis of theoretical model analysis, areas with higher aggregation degree of field-source locations from Euler inversions are closer to location of real solution. In this article, annual gravity changes in North China during 4 periods show better coincidence in the NE anomaly area in the west of Datong (Hetao fault zone), with field-source depths concentrated within 20km-40km. The location corresponds to the Shanxi-Hebei-Inner Mongolia danger zone that has been intensively monitored in recent years. The strip area, perpendicular to the Taihang Mountain piedmont fault, has a certain correlation with the fault structure.
(3) With respect to the sea areas, there are also some inversion results revealing the consistency between depths and locations, which we believe is due to the meshed extrapolation of original data, and there is no actual measurement point, thus, we get false solutions from inversion. For locations with poor concentration, we consider the changes of these field-source are random to a certain extent, which can be interpreted mainly as a certain extent of field-source (substances) motions existing in these locations in the short term, but the motions are not persistent, which should belong to a type of relatively short-term crustal motions, or is related to the uncertainty of observation error.
In conclusion, on the basis of theoretical models for reference, by inversion model design and selection of optimal parameters, calculations of field-source parameters for observed gravity field signals can be realized to get quantitative interpretation results. In this article, the application of the Euler deconvolution method in the research of field-source characteristics of gravity changes in North China has certain scientific significance for further understanding of potential seismic risk in North China. Moreover, it also has some reference value for further construction of a multidisciplinary quantitative index system for earthquake prediction and forecasting.
ACKNOWLEDGEMENTS: We extend our heartfelt thanks to the Special Project of Enhanced Monitoring and Tracing of Strong Earthquakes in North China, China Earthquake Administration, who provided gravity data. Anonymous reviewers are thanked for their helpful and constructive suggestions for this article.
This paper has been published in Chinese in the journal of Earthquake, Volume 36, Number 4, 2016.
Chen S, Xu Weimin, Jiang Changsheng. Relationship between gravity variation and seismic hazards in the western China[J]. Acta Seismologica Sinica, 2015, 37(4): 575–587. |
Chen Shi, Wang Qianshen, Zhu Yiqing, et al. Regional gravity variation before Wenchuan M_{S} 8.0 earthquake and epicentroid research[J]. Progress in Geophys., 2011, 26(4): 1147–1156. |
Chen Shi, Wang Qinghua, Wang Qianshen, et al. The 3-D density structure and gravity change of Ludian M_{S} 6.5 Yunnan epicenter and surrounding regions[J]. Chinese J. Geophys., 2014, 57(9): 3080–3090. DOI:10.6038/cjg201409333. |
Chen Yuntai, Gu Haoding, Lu Zhaoxun. Variations of gravity before and after the Haicheng earthquake, 1975 and the Tangshan earthquake, 1976[J]. Acta Seismologica Sinica, 1980, 2(1): 21–31. |
Chen S., Jiang C.S., Zhuang J.C. Statistical evaluation of efficiency and possibility of earthquake predictions with gravity field variations and its analytic signal in western China[J]. Pure and Applied Geophysics, 2016, 173(1): 305–319. DOI:10.1007/s00024-015-1114-x. |
Chen S., Liu M., Xing L., Xu W., Wang W., et al. Gravity increase before the 2015M_{W} 7.8 Nepal earthquake[J]. Geophysical Research Letters, 2016, 43. DOI:10.1002/2015GL066595. |
Deng Qidong, Ran Yongkang, Yang Xiaoping, et al. Active Tectonic Map of China[M]. Beijing: Seismological Press, 2007 |
Fan Meining. The Study and Application of Euler Deconvolution Method [D]. Doctoral thesis. Changchun: Jilin University, 2006. 19-23 (in Chinese with English abstract). http: //cdmd. cnki. com. cn/article/cdmd-10183-2006109772. htm |
Lu Baoliang, Fan Meining, Zhang Yuanqing. The calculation and optimization of structure index in Euler deconvolution[J]. Progress in Geophys., 2009, 24(3): 1027–1031. |
Reid A.B. Euler deconvolution: past, present and future, a review[J]. 65th Ann. Internat. Mtng., Soc. Expl. Geophys., Expanded Abstracts, 1995: 272–273. |
Reid A.B., Allsop J.M., Granser H., Millett A.J., Somerton I.W. Magnetic interpretation in three dimensions using Euler deconvolution[J]. Geophysics, 1990, 55: 80–91. DOI:10.1190/1.1442774. |
Zhang Peizhen, Deng Qidong, Zhang Guomin, et al. Active tectonic blocks and strong earthquake in the continent of China[J]. Science in China (Ser. D), 2003, 33(Suppl): 12–20. |
Zhu Yiqing, Wen Xueze, Zhang Jing. Dynamic variation of the gravity field in middle North China and its implication for seismic potential[J]. Chinese J. Geophys., 2013, 56(2): 531–541. |