Earthquake Research in China  2020, Vol. 34 Issue (3): 328-342     DOI: 10.19743/j.cnki.0891-4176.202003008
Quick Phase Identification for Dense Seismic Array with Aid from Long Term Phase Records of Co-located Sparse Permanent Stations
CHEN Yini1, LI Jun1, WANG Zhenjie1, ZHAO Mengqi1, YU Tiehong1, LI Danning3, YU Junyi1, WANG Weitao2
1. Zhejiang Earthquake Agency, Hangzhou 310013, China;
2. Key Laboratory of Earthquake Source Physics, Institute of Geophysics, China Earthquake Administration, Beijing 100086, China;
3. Yunnan Earthquake Agency, Kunming 650224, China
Abstract: The phase identification and travel time picking are critical for seismic tomography, yet it will be challenging when the numbers of stations and earthquakes are huge. We here present a method to quickly obtain P and S travel times of pre-determined earthquakes from mobile dense array with the aid from long term phase records from co-located permanent stations. The records for 1 768 M ≥ 2.0 events from 2011 to 2013 recorded by 350 ChinArray stations deployed in Yunnan Province are processed with an improved AR-AIC method utilizing cumulative envelope and rectilinearity. The reference arrivals are predicted based on phase records from 88 permanent stations with similar spatial coverage, which are further refined with AR-AIC. Totally, 718 573 P picks and 512 035 S picks are obtained from mobile stations, which are 28 and 22 times of those from permanent stations, respectively. By comparing the automatic picks with manual picks from 88 permanent stations, for M ≥ 3.0 events, 81.5% of the P-pick errors are smaller than 0.5 second and 70.5% of S-pick errors are smaller than 1 second. For events with a lower magnitude, 76.5% P-pick errors fall into 0.5 second and 69.5% S-pick errors are smaller than 1 second. Moreover, the Pn and Sn phases are easily discriminated from directly P/S, indicating the necessity of combining traditional auto picking and integrating machine learning method.
Key words: Phase picking     Travel time     Dense array     Spatial overlap

INTRODUCTION

Body wave tomography based on seismic wave travel time is an important tool for studying the interior structure of the earth. Teleseismic tomography can provide deep information within a depth of several hundred kilometers, while imaging studies based on regional and local seismic data can study the fine scale of the crust and mantle (Zhao Dapeng, 2015).

The density and spatial distribution of body wave rays from an earthquake epicenter to a station are the key factors that determine the quality of body wave tomography. In order to obtain more detailed underground structure, researchers often choose long-term records to increase the number of earthquakes, as well as installing dense observation stations so that the ray coverage will be improved for better imaging (Huang Jinli et al., 2006; Zigone D. et al., 2019). In recent years, with the deployment of dense arrays such as USArray, ChinArray and other large-scale broadband seismic arrays, the observation basis of body wave imaging has been greatly improved. Furthermore, some small-scale and ultra-dense observations also provide excellent data for shallow structure imaging (Li Yonghua et al., 2014; Lin Fanchi et al., 2013; Liu Zhen et al., 2017).

The increasing number of array observations also poses challenges to the traditional manual travel time picking work. The program ChinArray Phase I has been conducted in Yunnan Province and adjacent regions by China Earthquake Administration since the September 2010. From 2011 to 2013, 350 broadband portable stations have been deployed in the southern section of the North-South Seismic Belt (Ding Zhifeng et al., 2013). Numerous earthquakes have occurred in Yunnan, yet the identification of seismic phases on portable stations could not be conducted due to the limitation of manual picking efficiency with only conventional manual quick reports and cataloging processes in the permanent network of Yunnan. In the previous studies of body wave tomography based on ChinArray, teleseismic events are used to image the deep structure of the earth, because of the limited quantity of teleseismic body wave phases and the low workload of manual picking (Zhang Fengxue et al., 2018; Huang Zhouchuan et al., 2015). As for regional studies on ChinArray, major events with large magnitudes are analyzed. For example, Gao Jiayi et al. (2016) manually pick the travel times of 45 M≥3.0 earthquakes in the region, and construct a one-dimensional velocity model in Yunnan. For a large number of small-magnitude earthquakes, limited by the efficiency of manual recognition, the automated processing methods still need to be developed to achieve rapid phase recognition and travel time extraction.

Automatic identification of seismic phases has always been an important part of seismology (Allen R., 1982). As a result, the method of automatic identification has developed rapidly with the increase of seismic data volume. The classical STA/LTA (Short Term Averaging/Long Term Averaging) algorithm for automatic detection of microseismic signals based on waveform features (Allen R., 1982; Baer M. et al., 1987; Saragiotis C.D. et al., 2002; Baillard C. et al., 2014), AR-AIC (Leonard M. et al., 1999; Sleeman R. et al., 1999; Zhang Haijiang et al., 2003) have been widely used. In detail, the STA/LTA model is sensitive to variations in amplitude. AR-AIC is not only based on changes in amplitude, but also related to changes in frequency and phase. Normally we use STA/LTA to detect events, which can remove interference from pulse signals, and apply AR-AIC to mark the arrivals of specific phases (Akazawa T., 2004). However, when it comes to dealing with complex observation waveforms, the ability to pick up S-wave phases needs to be improved. By polarizing the three-component data and suppressing P code wave interference before S-wave arrives, the success rate of S detection can be effectively improved (Jurkevics A., 1988; Ross Z.E. et al., 2014). In recent years, deep learning techniques have also been applied to seismic phase recognition research to identify P and S phases automatically, including support vector machines, convolutional neural networks, etc (Yu Ziye et al., 2018; Ross Z.E. et al., 2018; Zhu Lijun et al., 2019).

Seismic phase recognition is mainly used in two aspects. Firstly, it is effective to identify seismic events from continuous data, distinguish seismic phases and pick up travel times, and then realize phase association and location monitoring through joint analysis of multiple stations (Zhang Miao et al., 2019). Secondly, it is used to pick up the travel time of different phases in the seismic record with known original time and location of the seismic event, including the fast picking of the travel time from active source signals (Xu Zhen et al., 2019), and the travel time picking of seismic events at different stations. Li Zefeng and Peng Zhigang (2016) use the existing catalog and initial one-dimensional velocity model to estimate the arrival time of events from each station, and then use a signal-to-noise-ratio-based characteristic function to pick up, updating the velocity model. After several times of iteration, more accurate arrival time and 1-D model can be obtained at the same time.

There are 127 seismic stations in Yunnan and its adjacent regions where the ChinArray stations are deployed (Fig. 1). These stations have been running steadily since 2008, and seismic phase reports and earthquake catalogs are regularly produced. Due to the overlap of the space area between the permanent and ChinArray stations, the underground space sample is relatively consistent, as well as the travel time of seismic events. Therefore, using the existing observation reports and earthquake catalogs of fixed stations, the corresponding seismic phases on portable stations of ChinArray can be quickly identified. Based on this, an improved AR-AIC method is applied to quickly analyze the travel time of the 2006 M≥2.0 earthquakes of on ChinArray from 2011 to 2013.

 Fig. 1 Station configuration in Yunnan area Solid triangles represent ChinArray distribution, hollow triangles indicate permanent stations; the circles mark the local seismic events occurred from 2011 to 2013
1 DATA AND METHOD

The target area, Yunnan and its adjacent regions (96°-108°E, 20.5°-30°N), includes 127 permanent seismic stations. From March 2011 to November 2013, the China Earthquake Administration has deployed 350 mobile broadband stations in the southern section of the North-South Seismic Belt as China Earthquake Scientific Exploration Array (ChinArray) Phase I, with an aperture of about 1 000 km and an average station distance of 35 km (Ding Zhifeng et al., 2013; Wang Weitao et al., 2018). The distribution of stations is shown in Fig. 1, in which all the stations are recorded with three components with a sampling rate of 100 Hz.

In the area of 20.5°-30.5°N and 96°-108°E, we select 1 796 events of ML≥2.0, among which 326 events are ML≥3.0. The distribution of seismic magnitude against time is shown in Fig. 2. The catalog information is provided by China Earthquake Networks Center (CENC).

 Fig. 2 Time-Magnitude distribution between March 2011 and December 2013

The theoretical travel times of the corresponding seismic phases are calculated using TauP (Crotwell HP et al., 1999) with the 1-D velocity model from Yunnan Station Network (Fig. 3). The results are roughly consistent with the manually recognized travel times. The crustal thickness in Yunnan is relatively thick with large lateral variation. We select the phase of the events with an epicenter larger than 200 km to analyze the Pn phases.

 Fig. 3 (a) 1-D velocity model from bulletin in Yunnan; (b) theoretical travel times corresponding to the 1-D model Black dots indicate the travel times from manual picking of events above ML2.0 between the year 2011 and 2013

According to the theoretical calculation, it is possible to estimate the travel time of various seismic phases of the local earthquakes, and conduct seismic phase identification. We firstly take event locations from CENC bulletin and a preliminary velocity model of Yunnan and its adjacent regions as input. Given the location of an event and the preliminary velocity model, the phase arrivals are predicted at several stations. We then apply an improved AR-AIC detector based on the characteristic function to search the actual arrival within a window around the predicted arrival. To determine the window for searching the real phase arrival, we assume that the actual velocity deviates within a percentage compared with the velocity model

 ${ \in = \frac{{\Delta {v_{{\rm{max}}}}}}{v}}$ (1)

in which v is a given velocity model and Δvmax is the maximal allowed difference between the true velocity structure and the given velocity model. By propagating of the velocity perturbation to arrival-time perturbation, the time window becomes

 ${\left[ {{t_0} + \frac{{{t_T}}}{{1 + \in }}, {t_0} + \frac{{{t_T}}}{{1 - \in }}} \right]}$ (2)

P phase is picked on the vertical component. 1) We first apply bandpass filter from 1 to 10 Hz in the time window; 2) and then calculate the cumulative envelope function of the filtered record (Akazawa T., 2004):

 $\begin{array}{l} {a_{v2}}t = {a_{v1}}\left({t - 1} \right)\left({if{a_{v1}}\left(t \right) \le {a_{v1}}\left({t - 1} \right)} \right){\rm{ }}\\ \;\;{a_{v2}}t = {a_{v1}}\left(t \right)\left({if{a_{v1}}\left(t \right) > {a_{v1}}\left({t - 1} \right)} \right) \end{array}$ (3)

in which av1 (t) is the filtered record from step 1), and av2 (t) is the achieved cumulative envelope function.3) Finally, we apply improved AR-AIC detection to the cumulative envelope (Fig. 4).

 Fig. 4 Picking procedure

After P phase arrival is obtained, the corresponding S phase arrival can be estimated based on the theoretical P-S time difference. And then the picking time window for S phase can be determined by Equa.(2) as well. For S phase identification, an S-wave polarization filter is applied to the horizontal components to suppress the energy of P coda as much as possible. After that, cumulative envelope function calculation and AR-AIC picker is applied as the same procedure in P phase identification (Fig. 4). The design of the S filter requires polarization analysis of the three-component waveform to calculate rectilinearity r and apparent vertical incidence angle ϕ. The details are introduced as follows (Jurkevics A., 1988; Ross Z.E. et al., 2014):

Firstly, take a sliding window of 3s for the three-component waveform, and calculate the covariance matrix σ of the waveform in each window

 $\sigma = \left(\begin{array}{l} {\rm{Cov}}\left({N, N} \right)\;\;\;{\rm{Cov}}\left({N, E} \right)\;\;\;\;{\rm{Cov}}\left({N, Z} \right)\\ {\rm{Cov}}\left({E, N} \right)\;\;\;\;{\rm{Cov}}\left({E, E} \right)\;\;\;\;{\rm{Cov}}\left({E, Z} \right)\\ {\rm{Cov}}\left({Z, N} \right)\;\;\;\;{\rm{Cov}}\left({Z, E} \right)\;\;\;\;\;{\rm{Cov}}\left({Z, Z} \right) \end{array} \right)$ (4)

and then calculate the eigenvalues λ1, λ2, λ3(λ1λ2λ3) and eigenvectors u=(u1, u2, u3) of this covariance matrix. Then, the linearity r and the apparent normal incidence angle ϕ are calculated as

 ${r = 1 - \frac{{{\lambda _2} + {\lambda _3}}}{{2{\lambda _1}}}}$ (5)
 ${\phi = {\rm{co}}{{\rm{s}}^{ - 1}}\left({{u_{11}}} \right)}$ (6)

The S-wave polarization filter is defined as

 ${S = r\left({1 - cos\left(\phi \right)} \right)}$ (7)
2 RESULTS 2.1 The Number of Picked Phases

The statistical result shows that based on the improved AR-AIC method, P-wave and S-wave can be well identified. Fig. 5 shows the ray coverage distributions without and with ChinArray stations of a ML3.3 earthquake at 19 :54 on November 9th, 2012, and its waveform plotting with automatic P and S pickers marked. It can be seen from Fig. 5 that the automatically picked seismic phases are in good agreement with the actual seismic phases.

 Fig. 5 The picking result and ray path coverage map ofML3.3 event The black triangles mark permanent stations, and the blue triangles mark the mobile stations

The density of ChinArray stations is much higher than that of the permanent stations. Through automatic phase picking, the number of regional phases and seismic rays can be greatly increased. Fig. 6 is aML2.1 event compared with the event showed in Fig. 5. It also shows the comparison between the number of phases picked up at fixed stations and the number of phases picked up by the ChinArray stations and their ray path coverages respectively. The number of phases corresponding to the two earthquakes is shown in Table 1. The original seismic stations are relatively sparse, thus the number of phases P recorded for a single event is between 10 and 40, and the number of S phases is between 10 and 30. After adding dense array data, the number of P phases that can be effectively picked up is between 200 and 400, and the number of S phases is between 200 and 300, which is approximately 10 times more than that with only fixed stations.

 Fig. 6 Ray path coverage distribution of aML2.1 event (a) is with only fixed stations, and (b) is with additional mobile stations

Table 1 Comparison between the number of phases at permanent stations and the number of phases adding ChinArray

The integration of mobile array data can also expand the recording scale of vibrations triggered by seismic events. For theML3.3 earthquake, we compare the distribution of seismic phases with the epicentral distance. Fig. 7 shows that when only permanent stations are used, the results can identify the distribution of seismic phases as far as about 290 km. After combining the array data, the phases picked can reach as far as 600 km, and it increases greatly as the epicentral distance getting further.

 Fig. 7 Distance distribution as a function of phase time Black traces are recorded by permanent stations, while blue traces are from dense array stations
2.2 The Accuracy of Phase Picking

Manual phase picking is relatively precise, but it requires massive labor cost. On average, the time cost for one phase is approximately 10-30 seconds, including identification, confirmation, labeling, review, and archiving. By contrast, such complex process only takes no longer than 1 second in automatic picking. Phase picking not only needs to be efficient, but also requires high accuracy for subsequent applications. Therefore, P-wave phase picking error is required to be within 0.5 second and the S phase picking error is required to be within 1 second.

Since the manual picking results are not recorded in the mobile stations, we also apply automatic picking on the permanent local network records with the same parameters and procedures as those for mobile network, and then compare the results with the bulletin from the local network. Theoretically, the accuracy of automatic picking of mobile network and local network should be consistent. It is indicated that for the events with ML≥3.0, 81.5% deviation time of P-wave arrival is less than 0.5 s, and for distance within 200 km, the value is 83.8%; 70.5% deviation time of S-wave arrival is less than 1s, for distance within 200 km, the value reaches 80.7% (Fig. 8(a)). For the events with ML≥2.0 andML < 3.0, 76.5% deviation time of P-wave arrival is less than 0.5 s, and with distances within 200 km the percentage is 76.9%; 69.5% deviation time of S-wave arrival is less than 1s, and for distance within 200 km the value is 79.9% (Fig. 8(b)). Obviously, the picking error for S-waves significantly increases with distances greater than 200 km, while P-wave picking do not show much difference. It's mainly because that after 200 km, the S-wave phase is relatively blurring and with lower signal-to-noise ratio, which is also a limitation for manual picking. Computer picking is often at the place where the amplitude changes the most, as performed in characteristic curve like AIC, STA/LTAs. In these cases, the smoothness of S-wave is better than the P-wave's. Thus, the manually picked points are usually prior to the computer-recognized time. In actual work, seismic phase identification is often corrected based on manual experience according to the seismic phase cycle changes.

 Fig. 8 (a) Arrival time deviation of events with ML≥3.0; (b) Arrival time deviation of events with 2.0≤ ML < 3.0

In seismic phase recognition the signal to noise ratio is an important factor that determines the accuracy of travel time picking no matter the phase is picked manually or automatically. Using fixed stations to automatically identify the signal-to-noise ratio of seismic phases, the relationship between the pickup error and the signal-to-noise ratio is analyzed. In the signal-to-noise ratio calculation, we select the signal waveform within T seconds after the arrival time of the signal to calculate the sum of squared amplitudes, while the noise is selected as T seconds before the signal arrives. The waveform is filtered through the frequency of 1-20 Hz, and the signal-to-noise ratio is calculated as

 ${{\rm{SNR}} = \sum\limits_{t = {t_0}}^{{t_0} + T} {{x^2}\left(t \right)} /\sum\limits_{t = {t_0} - T}^{{t_0}} {{x^2}\left(t \right)}, }$ (8)

in which t0 is the signal arrival time, x(t) is the bandpass filtered waveform. The window time T is 0.5 second for P-wave, and 1 second for S-wave.

Fig. 9 shows the relationship between the SNR and the picking error of 8 819 P waves generated by 267 earthquakes with ML>3.0. It can be seen that when the signal-to-noise ratio is relatively low, the pick-up error is relatively high. As the SNR increases, the pickup error gradually decreases. Overall, most errors are distributed within 0.5s.

 Fig. 9 Relationship between picking error and signal to noise ratio (SNR) of events with ML≥3.0

We extract the travel times of 2 006 earthquakes with ML>2.0 recorded by ChinArray and local network. Finally, 718 573 P-wave phases and 512 035 S-wave phases are extracted from the dense array, which are 28 times and 22 times of those in the bulletin reported by local network respectively.

Fig. 10 shows the travel time distribution of various phases within the epicentral distance of 1 200 km. It can be seen that compared with the number of phases at fixed stations, the number of various phases increases greatly after the automatic picking of phases combined with the array, and the corresponding epicentral distance also expands. Among various types of seismic phases, the number of Sn and S phases increase most, but the travel time deviation of the S phase is higher than that of the P-wave phase, which is associated with the complex shape of S-wave and the difficulty of identification. For different distances, the automatic recognition of seismic phases and the artificial recognition of seismic phases in the range of 0-400 km are relatively consistent. In the range of 600-700 km, the recognition of Pg and Sg phases becomes inferior, while the number of P and S phases increases correspondingly. The first arrivals of Pn and Sn are difficult to be recognized in both manual picking and machine learning results. By the assistance of the predicted travel time, in this paper, we can also effectively identify the first wave signals, especially the Pn wave. The number of Pn phases greatly increases in the range of 150-500 km, providing excellent data for subsequent imaging studies on the top of the Moho (Li Zhiwei et al., 2012).

 Fig. 10 Plots showing T travel time against distance distribution (a) Phases from manual results on local permanent network; (b) Total phases after automatic identification; (c) Details in (b) zooming into the distance between 550 km and 750 km

The identification of seismic phase is closely related to the signal-to-noise ratio and the complexity of the waveform. Although some seismic events are difficult to be identified manually, the automatic identification method provides the possibility to identify Pg, Sg, Pn and Sn phases (Fig. 11(a), (b)). Meanwhile, there are some relatively clear phases on the waveform that are not automatically recognized (Fig. 11(c), (d)). In the next work, it is hopeful to combine the traditional automatic recognition method with the newly developed method based on machine learning so that the recognition rate may be further improved.

 Fig. 11 Results of phase picking at several stations
3 CONCLUSION

In this paper, the improved AR-AIC method is used to automatically pick the phase arrival time from ML>2.0 earthquakes recorded by dense arrays deployed in the Yunnan area. Credit to the spatial coincidence of permanent network and mobile stations, existing fixed stations can be used to obtain earthquake prediction travel time and rapid phase recognition. A variety of seismic phases such as Pg, Sg, P, S, Pn, Sn in the mobile array can be obtained by using related methods. The automatic picking accuracy for most of the P-wave phases is within 0.5 s, and the S phase picking accuracy is within 1 s. The number of P-wave and S-wave phases recognized on the obtained mobile dense array stations is 28 times and 22 times of that of the original fixed station, respectively.

The method proposed in this paper has fast processing speed, and much higher efficiency than manual picking, providing foundation for the following body wave tomography. In order to achieve fast recognition, this paper uses a one-dimensional velocity model to predict arrival times. Due to the fact that Yunnan area has a complex velocity structure, three-dimensional velocity structure information should be added in the future. It is expected to further improve the recognition accuracy and perform detailed analysis on the travel time of events with smaller magnitudes.

REFERENCES
 Akazawa T. A technique for automatic detection of onset time of P-and S-phases in strong motion records: Proceedings of the 13th World Conference on Earthquake Engineering. International Association for Earthquake Engineering, paper. 2004(786). Aki K., Lee W.H.K.Aki K., Lee W.H.K. Determination of three-dimensional velocity anomalies under a seismic array using first P arrival times from local earthquakes:1.A homogeneous initial model[J]. Journal of Geophysical Research, 1976, 81(23): 4381-4399. DOI:10.1029/JB081i023p04381 Allen R.Allen R. Automatic phase pickers:their present use and future prospects[J]. Bulletin of the Seismological Society of America, 1982, 72(6B): S225-S242. Baer M., Kradolfer U.Baer M., Kradolfer U. An Automatic phase picker for local and teleseismic events[J]. Bulletin of the Seismological Society of America, 1987, 77(4): 1437-1445. Baillard C., Crawford W.C., Ballu V., Hibert C., Mangeney A.Baillard C., Crawford W.C., Ballu V., Hibert C., Mangeney A. An automatic kurtosis-based P- and S-phase picker designed for local seismic networks[J]. Bulletin of the Seismological Society of America, 2014, 104(1): 394-409. DOI:10.1785/0120120347 Ding Zhifeng, Wu Zhongliang, ChinArray Team. Advances of ChinArray program. In: Proceedings of AGU Fall Meeting. San Francisco: American Geophysical Union, 2013. Gao Jiayi, Li Yonghua, Xu Xiaoming, Zhang FengxueGao Jiayi, Li Yonghua, Xu Xiaoming, Zhang Fengxue. P-wave crustal velocity structure and seismicity in Yunnan area[J]. Progress in Geophysics, 2016, 31(2): 508-516 (in Chinese with English abstract). Huang Jinli, Zhao DapengHuang Jinli, Zhao Dapeng. High-resolution mantle tomography of China and surrounding regions[J]. Journal of Geophysical Research, 2006, 111(B9): B09305. Huang Zhouchuan, Wang Pan, Xu Mingjie, Wang Liangshu, Ding Zhifeng, Wu Yan, Xu Mijian, Mi Ning, Yu Dayong, Li HuaHuang Zhouchuan, Wang Pan, Xu Mingjie, Wang Liangshu, Ding Zhifeng, Wu Yan, Xu Mijian, Mi Ning, Yu Dayong, Li Hua. Mantle structure and dynamics beneath SE Tibet revealed by new seismic images[J]. Earth and Planetary Science Letters, 2015, 411: 100-111. DOI:10.1016/j.epsl.2014.11.040 Jurkevics A.Jurkevics A. Polarization analysis of three-component array data[J]. Bulletin of the Seismological Society of America, 1988, 78(5): 1725-1743. Leonard M., Kennett B.L.N.Leonard M., Kennett B.L.N. Multi-component autoregressive techniques for the analysis of seismograms[J]. Physics of the Earth and Planetary Interiors, 1999, 113(1/4): 247-263. Li Yonghua, Pan Jiatie, Wu Qingju, Ding ZhifengLi Yonghua, Pan Jiatie, Wu Qingju, Ding Zhifeng. Crustal and uppermost mantle structure of SE Tibetan plateau from Rayleigh-wave group-velocity measurements[J]. Earthquake Science, 2014, 27(4): 411-419. DOI:10.1007/s11589-014-0090-z Li Zefeng, Peng ZhigangLi Zefeng, Peng Zhigang. An automatic phase picker for local earthquakes with predetermined locations:combining a signal-to-noise ratio detector with 1D velocity model inversion[J]. Seismological Research Letters, 2016, 87(6): 1397-1405. DOI:10.1785/0220160027 Li Zhiwei, Ni Sidao, Hao Tianyao, Xu Yi, Roecker S.Li Zhiwei, Ni Sidao, Hao Tianyao, Xu Yi, Roecker S. Uppermost mantle structure of the eastern margin of the Tibetan plateau from interstation Pn traveltime difference tomography[J]. Earth and Planetary Science Letters, 2012, 335-336: 195-205. DOI:10.1016/j.epsl.2012.05.005 Lin Fanchi, Li Dunzhu, Clayton R.W., Hollis D.Lin Fanchi, Li Dunzhu, Clayton R.W., Hollis D. High-resolution 3D shallow crustal structure in Long Beach, California:application of ambient noise tomography on a dense seismic array[J]. Geophysics, 2013, 78(4): Q45-Q56. DOI:10.1190/geo2012-0453.1 Liu Zhen, Tian Xiaobo, Gao Rui, Wang Gaochun, Wu Zhenbo, Zhou Beibei, Tan Ping, Nie Shitan, Yu Guiping, Zhu Gaohua, Xu XiaoLiu Zhen, Tian Xiaobo, Gao Rui, Wang Gaochun, Wu Zhenbo, Zhou Beibei, Tan Ping, Nie Shitan, Yu Guiping, Zhu Gaohua, Xu Xiao. New images of the crustal structure beneath eastern Tibet from a high-density seismic array[J]. Earth and Planetary Science Letters, 2017, 480: 33-41. DOI:10.1016/j.epsl.2017.09.048 Ross Z.E., Ben-Zion Y.Ross Z.E., Ben-Zion Y. Automatic picking of direct P, S seismic phases and fault zone head waves[J]. Geophysical Journal International, 2014, 199(1): 368-381. DOI:10.1093/gji/ggu267 Ross Z.E., Meier M.A., Hauksson E.Ross Z.E., Meier M.A., Hauksson E. P wave arrival picking and first-motion polarity determination with deep learning[J]. Journal of Geophysical Research, 2018, 123(6): 5120-5129. Saragiotis C.D., Hadjileontiadis L.J., Panas S.M.Saragiotis C.D., Hadjileontiadis L.J., Panas S.M. PAI-S/K:a robust automatic seismic P phase arrival identification scheme[J]. IEEE Transactions on Geoscience and Remote Sensing, 2002, 40(6): 1395-1404. DOI:10.1109/TGRS.2002.800438 Sleeman R., Van Eck T.Sleeman R., Van Eck T. Robust automatic P-phase picking:an on-line implementation in the analysis of broadband seismogram recordings[J]. Physics of the Earth and Planetary Interiors, 1999, 113(1/4): 265-275. Wang Weitao, Gerstoft P., Wang BaoshanWang Weitao, Gerstoft P., Wang Baoshan. Seasonality of P wave microseisms from NCF-based beamforming using ChinArray[J]. Geophysical Journal International, 2018, 213(3): 1832-1848. DOI:10.1093/gji/ggy081 Xu Zhen, Wang Tao, Xu Shanhui, Wang Baoshan, Feng Xuping, Shi Jing, Yang MinghanXu Zhen, Wang Tao, Xu Shanhui, Wang Baoshan, Feng Xuping, Shi Jing, Yang Minghan. Active source seismic identification and automatic picking of the P-wave first arrival using a convolutional neural network[J]. Earthquake Research in China, 2019, 33(2): 288-304. Yu Ziye, Chu Risheng, Sheng MinhanYu Ziye, Chu Risheng, Sheng Minhan. Pick onset time of P and S phase by deep neural network[J]. Chinese Journal of Geophysics, 2018, 61(12): 4873-4886 (in Chinese with English abstract). Zhang Fengxue, Wu Qingju, Li Yonghua, Zhang Ruiqing, Sun Lian, Pan Jiatie, Ding ZhifengZhang Fengxue, Wu Qingju, Li Yonghua, Zhang Ruiqing, Sun Lian, Pan Jiatie, Ding Zhifeng. Seismic tomography of eastern Tibet:implications for the Tibetan Plateau growth[J]. Tectonics, 2018, 37(9): 2833-2847. DOI:10.1029/2018TC004977 Zhang Haijiang, Thurber C., Rowe C.Zhang Haijiang, Thurber C., Rowe C. Automatic P-wave arrival detection and picking with multiscale wavelet analysis for single-component recordings[J]. Bulletin of the Seismological Society of America, 2003, 93(5): 1904-1912. DOI:10.1785/0120020241 Zhang Miao, Ellsworth W.L., Beroza G.C.Zhang Miao, Ellsworth W.L., Beroza G.C. Rapid earthquake association and location[J]. Seismological Research Letters, 2019, 90(6): 2276-2284. DOI:10.1785/0220190052 Zhao Dapeng. Multiscale Seismic Tomography[M]. Tokyo: Springer, 2015. Zhu Lijun, Peng Zhigang, McClellan J., Li Chenyu, Yao Dongdong, Li Zefeng, Fang LihuaZhu Lijun, Peng Zhigang, McClellan J., Li Chenyu, Yao Dongdong, Li Zefeng, Fang Lihua. Deep learning for seismic phase detection and picking in the aftershock zone of 2008 MW7.9 Wenchuan Earthquake[J]. Physics of the Earth and Planetary Interiors, 2019, 293: 106261. DOI:10.1016/j.pepi.2019.05.004 Crotwell H.P., Owens T.J., Ritsema J.Crotwell H.P., Owens T.J., Ritsema J. The TauP Toolkit:Flexible seismic travel-time and ray-path utilities[J]. Seismological Research Letters., 1999, 70(2): 154-60. DOI:10.1785/gssrl.70.2.154 Zigone D., Ben-Zion Y., Lehujeur M., Campillo M., Hillers G., Vernon F.L.Zigone D., Ben-Zion Y., Lehujeur M., Campillo M., Hillers G., Vernon F.L. Imaging subsurface structures in the San Jacinto fault zone with high-frequency noise recorded by dense linear arrays[J]. Geophysical Journal International, 2019, 217(2): 879-893. DOI:10.1093/gji/ggz069