Abstract
The focal mechanisms of Mw 6.3 aftershocks, Chiang Rai Province, Northern Thailand, were determined by using a multistation waveform inversion. Three aftershocks were selected and their waveforms were inverted for moment tensor calculation. Waveform inversions were derived from three broadband stations with three components and epicentral distances less than 250 km after all seismic stations were considered. The deviatoric moment tensor inversion was used for focal mechanism calculations. Band-pass filtering in the range of 0.03–0.15 Hz was selected for reducing low- and high-frequency noise. Source positions were created by using a single-source inversion and a grid-search method computed to optimize the waveform match. The results showed stable moment tensors and fault geometries with the southwest azimuth in the northern part of the Payao Fault Zone (PFZ) with depths shallower than 10 km. Left-lateral strike-slip with a reverse component was detected. The tectonics of the PFZ is constrained by fault-plane solutions of earthquakes. WSW directional strikes are observed in the northern part of the PFZ.
1. Introduction
An Mw 6.3 earthquake occurred onshore on 05 May 2014 at 11:08:42 UTC in Mae Lao District, Chiang Rai Province, Northern Thailand, which directly affected Northern Thailand. Its hypocentre was reported by the Seismological Bureau, under the Thai Meteorological Department (TMD) as latitude 19.748°N, longitude 99.687°E, and 7 km depth. Global CMT Catalogue showed focal mechanisms: strike 1 = 67, dip 1 = 81, rake 1 = 0, strike 2 = 337, dip 2 = 90, and rake 2 = 171. The earthquake was felt by many people in Northern Thailand due to several shakings and the energy of the main shock dispersed to Chiang Mai City and far away to Bangkok.
After the main shock had occurred, 941 aftershocks were generated during 5–26 May 2014. The aftershocks consisted of eight events of Mw 5.0–5.9, 32 events of Mw 4.0–4.9, 154 events of Mw 3.0–3.9, and more than 747 events of Mw lower than 3.0 [1]. The main shock caused one person’s death and more than 1,000 people to be injured. Many buildings were damaged in seven provinces, such as temples, schools, and houses. Several earthquake ruptures made new overburden environments: sinkholes, surface cracks, and hot water upwelling.
Due to the occurrence of an Mw 6.3 earthquake, seismic waves generated at the earthquake source were propagated into the Earth’s crust and recorded by seismic stations on the Earth’s surface. The characteristics of earthquake waveforms can be used to determine fault-plane solutions and earthquake focal mechanisms using the deviatoric moment tensor (DMT) inversion. The earthquake focal mechanisms are important keys providing information on the stress field orientation [2]. The study area is located in Northern Thailand, where earthquakes of low to moderate magnitude and seismicity are characterized by continuous activity and frequency of occurrence. Geologically the area is characterized by basins, mountain ranges, and active faults (Figure 1). More than 40 basins appear in the tertiary age, with some basins containing oil fields. Most of the basins lie in N-S trending, are perpendicular to strike-slip tectonics, and also are separated by mountain ranges [3]. TMD located the Mw 6.3 epicentre (Figure 2) in the Payao Fault Zone (PFZ) separated into northern and southern parts [4, 5]. The northern part lies in the NE-SW direction with left-lateral strike-slip. The southern part of the fault lies in the N-S direction with right-lateral strike-slip.
This paper aims to present focal mechanisms of three aftershocks above magnitude 4. Only three aftershocks showed stable results. The synthetic and observed waveforms fit very well and nodal lines of P-wave polarities indicated in the same directions. The focal mechanisms of other aftershocks are not stable and are expected to be complex and difficult to identify for exact solutions. These may be assumed as a problematic model of the focal mechanisms.
2. Data and Method
After the Sumatra–Andaman Earthquake on 26 December 2004 occurred, more than 40 digital seismic stations were installed throughout Thailand and controlled by TMD. Seismic stations in Thailand are named TM network and are under coordinated by TMD. Eighteen broadband stations are available for data download using the Incorporated Research Institutions for Seismology (IRIS) system, but short-period stations are only available by direct contact with TMD. All stations are three-component seismometers of various models. Example models are Trillium 120 sec, BB KS2000M sec, and SP-S13-HZ. More details are described in the TMD website, http://www.seismology.tmd.go.th/. In this paper, two broadband seismometers (CMMT and MHIT) of Trillium 120 sec model with a nominal sensitivity of 1201 V/(m s) and one broadband seismometer (MHMT) of Trillium 40 sec model with a nominal sensitivity of 1500 V/(m s) were selected for DMT inversion (Figure 1). PHRA broadband station was not used because of a different model (BB KS2000M). Although instrument correction has been done, the amplitude was not on the same scale and may be caused by different companies. PHRA was rejected for DMT inversion.
Six stations, consisting of two short-period stations (LAMP and PAYA) and four broadband stations (CMMT, MHIT, MHMT, and PHRA), were used for hypocentral calculations. The hypocentres using the HYPOINVERSE computer program [6] integrated into SEISAN software [7] are presented in Table 1. An iasp91 velocity model [8] was selected as a 1D velocity model because of the Moho depth corresponding to Northern Thailand [9] and showed better results than other models. The errors of hypocentres are expected to be less than a few kilometres. These hypocentres were used for DMT inversion with a condition; epicentres fixed, time shifts, and depths varied.
The single-point source solution and the DMT inversion, which was composed from a DC (double-couple) and CLVD (compensated linear vector dipole) with VOL = 0%, were selected for focal mechanism studies and processing was done with freely available ISOLA Fortran code [10]. The code uses inverse problem formulations [11] based on six element moment tensors, published by Kikuchi and Kanamori [12] for evaluating the correlation between observed and synthetic waveforms. For the single-point source solution, latitudes and longitudes from HYPOINVERSE were fixed and depths varied from 0.5 to 35 km with 0.5 km increments. A distance weight was not applied because hypocentral distances from the used stations were assumed small. Centroid depths and Green’s function were calculated by a 3D spatial grid search and by a frequency-wavenumber method [13], respectively. To calculating the Green’s function, an iasp91 velocity model was also used. The maximum frequency of the Green’s function was a limit to 0.15 Hz. Densities of crustal media were calculated using the following equation [14]:
Instrument corrections were carried out before DMT inversion began. The corrections included DC and trend removal. Synthetic and observed waveforms were band-pass filtered in a frequency range of 0.03 Hz to 0.15 Hz. Lower than 0.03 Hz was not expected due to long-period noise and a frequency limit of seismometers. Although high-frequency waves are more sensitive than low-frequency waves for small to medium magnitude earthquakes, higher than 0.15 Hz was already tested and waveforms did not fit. This may be caused by hypocentral distances not small enough (distance < 10 km). The selected frequency ranges were tested and expanded in the results and discussion section. After band-pass filtering and instrument corrections were done, the data were converted from count into displacement units in metres. Finally, the corrected data were cut from the hypocentral time to 250 seconds and resampled from a frequency of 100 Hz to 33 Hz. The 250-second range covers all the earthquake events.
The DMT inversion was processed by minimizing the difference between the observed and synthetic data in the form of displacements. A least-square sense was set at a group of trial origin time and trial source position. As the inversion was running, an optimum depth and optimum time were searched. The depth increment was set by following a Green’s function parameter as shown in a paragraph before. The optimum time was performed by predefined time steps. A time step is 0.02 s that starts from −5 to +5 s referred to as the hypocentral time that was calculated from HYPOINVERSE. The optimum depth and optimum time are called the centroid depth and centroid time as summarized in Table 2.
A grid-search method performs the best centroid positions (epicentres and depths) and also a time in terms of the absolute value of the correlation coefficient between the data and synthetic values. The match between the observed and synthetic waveform is identified by variance reduction.where , , is synthesis, and is original waveforms along with the summation of all collected data. The higher value of var.red. indicates the better fit between observed and synthesis waveforms. The three-component waveform inversions were computed using an iterative deconvolution method [12]. A waveform inversion approach was followed and inversed without any separation of body and surface waves. The waveform fit was optimized during a grid search of various trial positions.
3. Results and Discussion
All aftershocks from TMD were analysed and only three events were selected to be presented in this paper with a condition of high-quality and stability of moment tensor results and also assumed nonproblematic results. Other aftershocks were affected by overlapped events, noise for low-magnitude events, DC < 70%, and low var.red. After inversions, distances from the earthquakes to the recording stations were calculated and less than 250 km. These distances are also less than global earthquake agencies reported, such as IRIS, USGS, and GEOFON. It can be assumed that the resolutions with higher frequencies of our results are better than other agencies with lower frequencies because of less distances. The results are summarized in Table 2. The hypocentres are very constrained within the northern part of PFZ and also located in the south of the Mw 6.3 epicentre.
The uncertainty of earthquake locations is shown in Table 1 by the values of RMS (s) and is less than 1 s. The stable inversion of focal mechanisms is shown in Figure 3 in the nodal planes of P-wave polarities. Three frequency ranges were designed to consider the results depending on the frequency ranges. The 0.03–0.15 Hz was selected for waveform inversions because P-wave polarities are in good agreement with nodal lines in the same direction. These indicate the stability of the results. DC% is higher than 70 and event time in Table 2 not too much different compared with the time calculated from HYPOINVERSE software in Table 1. The high-frequency wave (0.03–0.15 Hz) is better than the low-frequency waves (0.03–0.08 Hz and 0.03–0.10 Hz) for detecting small-scale features [15]. The centroid depths were observed shallower than 10 km. This may suggest that the 6.3 Mw main shock is a shallow earthquake. The optimum depths and times estimated by the grid-search method are shown in Figure 5. The result shows that this method is suitable and gives a depth and time shift with high DC percentages and high correlation coefficients. Focal depths and time shifts with maximum correlation value called the best correlation are shown by the largest beach balls in Figure 5. The DC percentages and the correlation coefficient were drawn in beach ball colours and background contours, respectively. The maximum correlation values of three aftershocks are located in pink. The maximum correlation value results that the grid-search method gives a positive value of time shift. The time shifts indicate that hypocentral times calculated by HYPOINVERSE are earlier than the grid-search method calculated for only three aftershocks studied in this paper. The hypocentral times come later than the HYPOINVERSE time because HYPOINVERSE gives the time of rupture initiation while the ISOLA results in the time of the moment tensor release.
A good fit between observed and synthetic waveforms is shown in Figure 4. All included in the waveform time range used in the inversion are complete in a time range of 0–150 s. The good fitting shows that the exact solutions, noncomplex focal mechanisms, and iasp91 velocity model can be used in this study area. Not only the main features but also the first motion of P-wave polarities is fitted. The blue texts are var.red. of each recorded component as represented in (2). The good agreement of P-wave polarities (Figure 3) and the well-fitted waveforms (Figure 4) shows that fault geometries in Table 2 are reliable. However, low fits at later times were found in components of MHIT station (Figure 4). A fit of the waveform at early times is more important for waveform inversions and a poor fit at later times does not lead to failure of waveform inversion [16].
Normally, the DMT inversion creates two nodal planes. Only one nodal plane agrees with fault lines in a topographic map (Figure 2), which is specified and presented in Table 2 as the values of strike, dip, and rake. These are drawn as black beach balls in Figure 2. The strike directions of three aftershocks showed the same WSW direction. Consideration of strike directions indicates that the northern part of the PFZ is left-lateral strike-slip. The strike directions are also parallel to the fault lines that appear on a topological map as shown in Figure 2. Moreover, the results also showed that the dip angle of the northern part of the PFZ is close to a vertical fault with more than 60°, especially 86° in Ev.1. Strike-slip, normal, and reverse faulting can be identified by rake angles. The rake angles of three events are positive values. The rake angles of Ev.1 and Ev.3 are 31° and 38°, respectively. These rakes indicate reverse faulting. Only for Ev.2 is the rake angle close to zero. It is identified as strike-slip faulting. In addition, the rake angles of all inverse events denote left-lateral strike-slip.
4. Conclusions
The fault-plane solutions of three aftershocks were obtained with the DMT inversion. The selected aftershocks revealed similar focal mechanisms and showed the southwest azimuth. The rake angles indicated that the northern part of the PFZ is characterized by left-lateral strike-slip and reverse faulting. The high fitting between observed and synthetic waveforms shows that the iasp91 velocity model can be used for focal mechanism observations within the northern part of the PFZ. The stability of nodal lines from P-wave polarities is good and an important key for considering results from the DMT inversion.
Conflicts of Interest
The author declares that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
The author sincerely acknowledges Seismological Bureau, Thai Meteorological Department, for earthquake waveforms. More thanks go to Efthimios Sokos and Jiri Zahradnik for the ISOLA Fortran code.