Seismic Structure of Local Crustal Earthquakes beneath the Zipingpu Reservoir of Longmenshan Fault Zone

Three-dimensional P wave velocity models under the Zipingpu reservoir in Longmenshan fault zone are obtained with a resolution of 2 km in the horizontal direction and 1 km in depth. We used a total of 8589 P wave arrival times from 1014 local earthquakes recorded by both the Zipingpu reservoir network and temporary stations deployed in the area. The 3-D velocity images at shallow depth show the low-velocity regions have strong correlation with the surface trace of the Zipingpu reservoir. According to the extension of those low-velocity regions, the infiltration depth directly from the Zipingpu reservoir itself is limited to 3.5 km depth, while the infiltration depth downwards along the Beichuan-Yingxiu fault in the study area is about 5.5 km depth. Results show the low-velocity region in the east part of the study area is related to the Proterozoic sedimentary rocks. The Guanxian-Anxian fault is well delineated by obvious velocity contrast and may mark the border between the Tibetan Plateau in the west and the Sichuan basin in the east.


Introduction
The 12 May Wenchuan earthquake occurred beneath the Longmenshan thrust belt and was the most destructive earthquake in China in the last 30 years.Due to the 10 km distance of the Wenchuan epicenter to the Zipingpu reservoir and abnormal seismicity increases around the Zipingpu reservoir, the possibility that the Wenchuan earthquake was a reservoir-induced earthquake is obviously an important issue and has become a debating issue [1][2][3].Most recent publications [4][5][6][7][8] give possible scenarios obtained from the calculation of the Coulomb stress change under a set of assumptions.However, those results are critically dependent on assumptions such as diffusivity and the fault plane orientation.Ge et al. presented that the Zipingpu reservoir potentially hastened the occurrence of the Wenchuan earthquake by tens to hundreds of years [4].Zhou et al. [7] repeated Ge et al.'s work and found that an improper dip angle parameter might lead to a wrong conclusion.Their modeling results based on the 2-D model and 3-D model with proper fault parameters show the Coulomb stress changes alone neither were large enough nor had the correct orientation to affect the occurrence of the Wenchuan earthquake.
During recent years, different seismic studies have been performed to explore the relationship between the crustal velocity structure and the occurrence of the Wenchuan earthquake [9,10] or aftershocks of the Wenchuan earthquake [11,12].However, the large distance between stations in Zipingpu reservoir has hindered the achievement of enough resolution of the seismic images to correlate the velocity anomalies with geological structures and possible connections between velocity structure and abnormal seismicity increases since the filling of the Zipingpu reservoir.In this paper, the results of a seismic tomography study are shown based on local earthquakes near the Zipingpu reservoir.The resolution obtained is higher than in previous studies.The seismic images are compared with the geological setting and the submerged region of the Zipingpu reservoir.

Geological Setting
The Zipingpu reservoir lies in the Longmenshan fault zone which is characterized by a NE strike and composed of International Journal of Geophysics three major faults: the Wenchuan-Maoxian fault in the northwest, the Beichuan-Yingxiu fault in the central, and the Guanxian-Anxian fault in the southeast (Figure 1).Analyses of field measurements of surface ruptures and coseismic deformations suggest that the Beichuan-Yingxiu part of the Longmenshan fault zone close to the Zipingpu reservoir is a northeast-striking northwest-dipping thrust with small right-lateral strike-slip component [13,14].The outcrops between the Wenchuan-Maoxian and Beichuan-Yingxiu faults in the study area are Proterozoic igneous Pengguan massif composed of diorites and granites and the other areas are mainly covered by Proterozoic sedimentary rocks which are older towards the Pengguan massif.

Aftershock Observations.
To investigate the seismic structures around the submerged region of the Zipingpu reservoir, we refined the reservoir monitoring network which includes seven stations around the Zipingpu reservoir (1 destroyed by the Wenchuan earthquake) with 6 temporary seismic stations at the beginning of the year of 2009.The seismicity in the studied area is characterized by abnormal increases since the filling of Zipingpu reservoir, especially after the occurrence of the Wenchuan earthquake.Most of earthquakes recorded by above 12 stations have low magnitude (M < 3).An et al. [10], using digital data recorded by 26 seismic stations, show the results of P-velocity model around the main shock epicenter region.They did not emphasize their focus on the Zipingpu reservoir during their study or previous studies.In this paper, we have more appropriate station distribution around the Zipingpu reservoir.

Data Selection.
We have used P wave arrival times from digital data recorded in 2009 by our temporary 6 seismic stations and the Zipingpu reservoir seismic network which belongs to reservoir institute of Sichuan Earthquake Administration.Figure 1 shows the distribution of seismic stations used in this study.We have selected the local earthquakes which have been recorded at least five stations and whose standard errors of epicenters are not further than 8 km.This procedure helps in the relocating of local earthquakes.Figure 2 shows the P wave and S wave travel time with epicentral distance at station LYS, before and after removal of outliers.After removal of any arrival time that is far away from the distribution trend with the eventstation distances, the local earthquakes were relocated using International Journal of Geophysics Zhang and Thurber's [15] program which can produce not only the relative relocations but also absolute relocations.
We treat all the events as one cluster.Because travel time tomography was done in this method, relocations more reliable than those obtained just from relative relocation algorithm.The tomography grids were arranged at intervals of 10 km and 3 km in horizontal and vertical directions separately.Careful tomography results were obtained by Hole and Zelt's tomography method [16].Finally, we selected a total of 1014 from 1085 local earthquakes, with 8589 P wave travel time observations.The average differences between initial and relocated hypocenters are 100 m, 91 m, and 198 m in horizontal and vertical directions.Due to topography and station elevations, the parameter of air depth is set as −1 km.So, the depth of the selected earthquakes ranges from −0.8 km to 10.5 km.The distribution of the relocation in horizontal profile is shown in Figure 3, where the ellipses denote the hypocenter of earthquakes given by doubledifference tomography program.

Methodology and Resolution.
The tomography method of Hole and Zelt [16] was used to determine the 3-D P wave velocity structure.This method introduces special headwave operators to better deal with the presence of strong velocity contrast and permits placing the shots at free location of the velocity model.Twelve seismic stations were set as the shots and local earthquakes as the receivers in order to simplify the input files.The area selected for the tomographic study is located between 30 • 55 N and 31 • 10 N and from 103 • 20 E to 103 • 38 E, comprising 36 km (latitude) × 32 km (longitude).
We set up a 3-D grid in this study with a grid spacing of 1 km in the horizontal and vertical directions.Lateral and vertical model smoothing as regularization was introduced to result in a more physically reasonable structure model.After removal of outliers, the average RMS of travel time residuals is 0.67 s for the reference 1-D model, and it is reduced to 0.25 s for the tomography model after the fourth iteration.We applied a checkerboard resolution test to examine the resolution scale of the present data set.We assigned positive and negative velocity anomalies of 3% to the reference model.The synthetic data was calculated by 3-D finite-difference forward modeling.Random noise was not added to the synthetic data.The inverted model of the checkerboard tests retrieved checkers with a size of 2 km in lateral directions beneath the submerged region of the Zipingpu reservoir (Figure 4).

Results and Discussion
The P wave velocities from the seismic tomography are displayed at different depths (Figure 5).The most robust features imaged at shallow depths are the several pronounced low-velocity regions at the background of the velocity larger than 6 km/s.We can see that low-velocity regions in the central of the study area were placed close to the surface trace of the Zipingpu reservoir, and those low-velocity regions extend downwards and vanish at 3.5 km depth.The obvious features at the depth of 3.5 km are two low-velocity regions trending NE-SW which correspond to the Beichuan-Yingxiu strike-slip thrust fault and the Guanxian-Anxian fault.
It is widely recognized that low-velocity regions are a feature of some old or active strike-slip zones [17].Especially for the place where there is a significant fracturing with contained fluids, the fault zone will exhibit low velocity [18,19].Hence, these features of the low-velocity distribution can be related with the infiltration water of the Zipingpu reservoir.From 5.5 km depth downwards, the low-velocity region indicating the range of infiltration water along the Beichuan-Yingxiu strike-slip thrust fault vanished.It is worth noting that the zone of high pore pressure is not the same as the zone of water infiltration.Due to probable existence of high pore pressure deeper, caused by the reservoir loading, we are not sure whether the impoundment of the Zipingpu reservoir triggers the main shock of 2008 or not.If the Zipingpu reservoir triggered the 2008 main shock, the epicenter of it should be located on the Beichuan-Yingxiu fault which is the main passage of leakage or pressure transmission deduced from the low-velocity regions.
As shown in Figure 5, with increasing depth, the range of the high velocity part indicating the Pengguan massif [10] gradually diminishes, and the range of the relative low velocity in the east of the study area broadens westwards.Those results support the recognitions that the Longmenshan fault thrusts eastwards to the Sichuan basin.The distribution of the low velocity means that the Guanxian-Anxian fault marks the border between the Tibetan Plateau in the west and the Sichuan basin in the east even though the surface rupture of this fault lies in the Sichuan basin.

Summary and Conclusions
Through the use of digital data from more appropriate station distribution around the Zipingpu reservoir, seismic images are obtained with high resolution.Results not only confirmed previous sharp velocity changes at both sides of the fault and higher velocity of the Pengguan massif, but also found the correlation between low velocity and possible infiltration zone of the Zipingpu water reservoir.Velocity distribution provided a sign of water diffusion, but the conclusions from the results may need to be examined thoroughly in the future by other observable geophysical properties, such as electric and magnetic properties.
The principal findings of this study are as follows.
(1) The low-velocity regions placed close to the surface trace of the Zipingpu reservoir are correlated with the submerged region of the Zipingpu reservoir.Those low-velocity regions under the Zipingpu reservoir extend downwards and vanish at 3.5 km depth, which means that water possibly infiltrates downwards to this depth.
(2) At 3.5 km depth, two pronounced NE-SW trending low-velocity regions can be associated with the Beichuan-Yingxiu fault and the Guanxian-Anxian fault.The low-velocity region along the Beichuan-Yingxiu fault can be traced up to 2.5 km depth, where it is connected with the other low-velocity regions below the western branch of the Zipingpu reservoir.At 5.5 km depth, the low-velocity region along the Beichuan-Yingxiu fault fades.So, possible downwards infiltration depth along the Beichuan-Yingxiu fault in the study area is deduced at about 5.5 km depth.
(3) The Guanxian-Anxian fault is well delineated by obvious velocity contrast.The distribution of the low velocity in the east of the study area indicates that the Guanxian-Anxian fault marks the border between the Tibetan Plateau in the west and the Sichuan basin in the east.

Figure 1 :
Figure 1: Geological sketch map of the study area.The areas filled with plus symbols are igneous rocks of the Pengguan massif, and the other areas are sedimentary rocks.Black open triangles are seismic stations of the Zipingpu reservoir network.Black closed triangles are the temporary seismic stations.Bold lines denote the faults.Fine lines denote the river.Gray region represents the Zipingpu reservoir.Squares represent the county.LYS represents the name of the station.

Figure 2 :
Figure 2: Travel time variations of the station LYS with change of the event-station distances before (a) and after (b) removal of abnormal travel times.

Figure 3 :
Figure 3: Distribution of the local earthquakes used in the present study.Ellipses denote earthquake epicenters given by the doubledifference relocation.Solid triangles are all stations (temporary + permanent) used in the relocation and tomography inversion.The plus symbols present the area of the tomography grids and that of velocity inversion as shown in Figure 5.The others are the same as in Figure 1.

Figure 4 :
Figure 4: Results of checkerboard resolution test.Open and closed circles denote slow and fast velocities, respectively.The depth of each layer is shown at the upper part of the map.The grid spacing is 2 km in the horizontal direction and 1 km in depth.The perturbation scale is shown on the right.The coordinates are in km.

5 DepthFigure 5 :
Figure 5: P wave velocity models at different depth.The depth of each layer is shown at the upper part of the map.Red triangles represent the stations.The P wave velocity scale is shown on the bottom.Color scale is absolute P wave velocity (km/s).Background velocity near the edges denotes the initial 1-D velocity model.The coordinates are in km.