Abstract
Seismic interferometry (SI) has been recently employed to retrieve the reflection response from natural earthquakes. We perform experimental study to apply SI to Ocean Bottom Seismogram (OBS) records in the Nankai Trough, southwest Japan in order to reveal the relatively shallow geological boundaries including surface of oceanic crust. Although the local earthquakes with short raypath we use to retrieve reflection response are expected to contain the higher-frequency components to detect fine-scale structures by SI, they cannot be assumed as plane waves and are inhomogeneously distributed. Since the condition of inhomogeneous source distribution violates the assumption of SI, the conventional processing yields to the deteriorated subsurface images. Here we adopt the raypath calculation for stationary phase evaluation of SI in order to overcome this problem. To find stationary phase, we estimate the raypaths of two reflections: (1) sea-surface P-wave reflection and (2) sea-surface multiple P-wave reflection. From the estimated raypath, we choose the crosscorrelation traces which are expected to produce objective reflections considering the stationary phase points. We use the numerical-modeling data and field data with 6 localized earthquakes and show that choosing the crosscorrelation traces by stationary phase evaluation improves the quality of the reflections of the oceanic crust surface.
1. Introduction
Among various seismic exploration methods using the body-wave of natural earthquakes [1, 2], seismic interferometry (SI) has been recently employed to retrieve the reflection response. Although the receiver function method [1] has been broadly used to image the Moho and mantle discontinuities, there is a study claiming that retrieving and migrating the reflection response using SI is superior to the receiver function method [3].
SI retrieves Green’s function between receivers by crosscorrelating wavefield [4, 5]. This theory requires the physical sources homogeneously distributed along the enclosed surface which surrounds the receivers [4]. There are several successful applications of SI to natural earthquakes. Abe et al. [3] and Tonegawa et al. [6] retrieved the crustal reflection response in central Japan using P coda and S coda, respectively. Ruigrok et al. [7] used P wave to retrieve reflection response using Laramie array in USA. Abe et al. [3] further showed the comparison of migrated images of SI and those using receiver function analysis.
These applications focused on the teleseismic wavefields in which the epicentral distance is much longer than the length of receiver array. Consequently the wavefields can be assumed as plane wave. Ruigrok et al. [7] replaced the integral for source position of SI into a one for ray parameter. These teleseismic events are suitable to SI processing since the most earthquakes are generated at the numerous plate boundaries in the world, and consequently the teleseismic records contain the earthquakes propagating from the various directions. This can moderately enable us to assume that their source distribution is homogeneous.
Here, we performed experimental study to apply SI to the localized earthquake records acquired by Ocean Bottom Seismogram (OBS) in the Nankai Trough, southwest Japan in order to reveal the relatively shallow geological boundaries including surface of oceanic crust (i.e., plate boundary). There is no application of SI to the natural earthquake recorded by OBS. The OBS we used in this study was deployed to observe the local earthquake for earthquake observation in subduction zone [11, 12]. Those local earthquakes cannot be assumed as a teleseismic wavefield since the source-receiver distance is smaller than the length of receiver array. Furthermore, those localized earthquakes are usually inhomogeneously distributed and violate the assumption of SI. Therefore, the conventional SI processing (crosscorrelation and summation) yields to the deteriorated subsurface images. However, focusing on the local earthquakes with shorter raypath can give us the advantages over using the teleseismic wavefields because the teleseismic events usually have long propagating path and lose their higher-frequency components, which leads to detect only large-scale structure such as Moho reflections. On the other hand, the local earthquakes which have shorter raypath are expected to contain higher-frequency components and to resolve more fine-scale structures by using SI. In this study, we discuss a method to retrieve the reflection response using localized natural earthquakes observed by OBS; we adopt the raypath calculation for stationary phase evaluation in SI analysis in order to overcome the noise originated from the violation of homogeneous source distribution.
The physical interpretation of the condition to be posed on the SI can be explained by the stationary phase approximation [5, 13, 14]. In this approximation, it explains that the dominant contribution of the retrieval of the Green’s function comes from the crosscorrelation of the records from the physical source located at the stationary phase position (stationary phase source). Therefore, in the localized source distribution, crosscorrelation pairs with stationary phase positions have physical meaning. The other crosscorrelation traces will produce noise and deteriorate the quality of the imaging results. In this study, we identify crosscorrelation pairs with stationary phase position using the estimated raypath information of two reflections: (1) sea-surface P-wave reflection and (2) sea-surface multiple P-wave reflection.
Note that Chaput and Bostock [15] successfully retrieved reflection response using the subsurface noise sources located from 10 km to 60 km depth which is similar source depth to our study. They found that the source illumination is imperfect from the discussion using the stationary phase approximation. Our study is different from them in the point that we focus on the natural earthquakes.
2. Estimation of Stationary Phase Records Using Raypath Calculation
2.1. Stationary Phase Approximation of Seismic Interferometry
Seismic interferometry by crosscorrelation (e.g., [4]) is written as where and are the observed wavefield from the sources along the closed surface . is a Green’s function between two receiver positions.
We assume that the primary reflection is retrieved by crosscorrelating a direct wave and a specular reflection from the physical sources (Figure 1(a)). where denotes a travel time of a direct wave from to . denotes a travel time of a specular reflection from the source position to the receiver position through specular reflection point (Figure 1(a)). Superscript or indicates that we only consider a direct wave for and reflected wave for . Note that the amplitudes of these wavefields are assumed to be normalized in (2). Substituting these waves into (1) and applying stationary phase approximation (e.g., [5, 16]) yields the following equation: where denotes a travel time of the specular reflection between two receiver positions through its specular reflection point (Figure 1(b)). denotes a coefficient for stationary phase approximation. Note that we removed in (1) since we only consider causal part of Green’s function. The acausal part is obtained by considering a direct wave for and reflected wave for . The integral in (3) has a stationary point at , and the objective primary reflection is retrieved [5]. Furthermore, in this stationary point, the following relation is satisfied: This relation states that, for the stationary phase source position, the two events (direct wave and reflected wave in this case) have same raypath from the source position to the receiver position (Figure 1(b)). This corresponds to the fact that the crosscorrelation processing subtracts the travel times, cancels the common raypath, and produces the traveltime of the objective reflection event between receivers.
(a)
(b)
Note that the amplitudes derived from the crosscorrelation of the nonstationary phase sources are cancelled after the summation of the homogeneously distributed sources along the enclosed surface. In the case of the localized source distribution, the cancellation is insufficient and the unwanted amplitudes remain.
2.2. Selection of Receiver Pairs by Stationary Phase Evaluation
When the physical sources are widely distributed, the stationary phase sources effectively produce the objective reflection events. On the other hand, when the source distribution is localized, only the reflection events with the stationary phase sources are retrieved. Therefore, we evaluate the crosscorrelation traces for the existence of the stationary phase source and exclude those without the stationary phase sources in the summation of the crosscorrelation traces. We assume that we can estimate the raypath propagating from the source position to the two receivers. When these two raypaths have common pathway, we define that the crosscorrelated trace using these receivers contains objective reflections.
In order to evaluate the existence of the stationary phase sources, the raypaths for the direct waves and the arbitrary multiple reflected waves are needed to be estimated. We adopted a method developed by Tamagawa et al. [17] for raypath calculation. This method geometrically calculates the raypath for arbitrary multiple reflections with given 3 dimensionally deviated structures assuming straight raypath (Figure 2). When we have two reflector planes, the method calculates the mirror point of the source position for the reflector planes ( in Figure 2). A mirror point is defined as a projection of the position with plane symmetry. The arbitrary multiple reflection raypath can be derived by connecting those mirror points to the receiver position ( in Figure 2). This method is simple and fast to calculate raypaths.
We evaluate the stationary phase using these raypaths. For example, we assume that the receiver observes direct wave, and two different receivers and observe reflected wave (Figure 3). In this model, the candidate for crosscorrelation processing is and . Stationary phase evaluation using raypath can give us the information that only the crosscorrelation of has a stationary phase source and contain the objective reflection event. By evaluating this procedure for all source position and receiver combination, we remove the crosscorrelation traces which do not have stationary sources for objective reflections. Note that we need some knowledge of the position of reflector to evaluate raypaths before imaging. However, a rough estimate of the position of reflector can be sufficient since we do not need precise information about the location of the stationary points.
3. Numerical Modeling Results
We numerically simulate the wavefield with localized source distribution in order to evaluate the processing considering stationary phase point. For the simulation study, we consider the local earthquakes recorded by Ocean Bottom Seismogram (OBS) on the seafloor. The objective reflection events are assumed to be those generated from the oceanic crust surface. We use only one physical source for the simulation study because we evaluate our proposed method for localized source distribution. The velocity model is a three-layer structure including sea water (Figure 4). The dipping structure in the velocity model simulates a oceanic crust surface (i.e., plate interface) in the Nankai Trough.
We define the two events which produce the objective reflection event for the stationary phase evaluation. In our source-receiver configuration, the source is located in the subsurface, and the sea-surface (modeled as a free-surface) expects to produce the strong downpropagating reflections. This is the difference of the OBS recorded wavefields with those by the surface receivers. Therefore, we assume that the sea-surface P-wave reflections (denoted as pP in Figure 4) and the sea-surface multiple P-wave reflections (pPp in Figure 4) produce the reflection events from the oceanic crust surface after crosscorrelation. Although other events may contribute to produce the objective reflection event (e.g., the crosscorrelation of the direct P-wave and the multiple P-wave reflections; Abe et al. [3]), we consider the crosscorrelation of these two events (pP and pPp) since they propagate with high energy and are easy to be separated from S-waves.
We numerically modeled the seismic wavefield with a one subsurface source using the 2-dimensional staggered grid method [18]. The subsurface source is installed at the 30 km depth (Figure 4). We simulate the earthquake by giving horizontal stress as a source and using a Ricker wavelet with the central frequency of 5 Hz. The length of receiver array is 50 km with 200 m intervals.
Figure 5 shows the calculated wavefield (vertical component) at the receiver array. One can see that the sea-surface reflections (pP) and the sea-surface multiple reflections (pPp) are observed as well as the direct P-wave and direct S-wave (Figure 5). Since we focus on only P-waves in this study, we muted the amplitude of the direct S-wave and subsequent wavefields. Because we have 251 receivers, the total number of the crosscorrelation traces is which is the number of 2 combinations from 251 elements. Figure 6 shows the result of subsurface image derived from all crosscorrelation traces using prestack depth migration (PSDM). The target-dipping structure was appeared on the image. However, the signal-to-noise ratio at the near offset is low, and the artifact-dipping events are appeared (arrows in Figure 6). These artifacts are caused by the crosscorrelation traces which do not have a stationary phase source for objective reflections as pointed out by, for example, Snieder et al. [19].
We evaluate the stationary phase source in order to remove the unwanted crosscorrelation traces. We estimate the raypath of pP and pPp and compare their raypaths. We fix the receiver which observe pP (denoted as ) and calculate the raypath of pPp for all receivers (. We calculate the horizontal distance of the receiver and the point where pPp passes through the seafloor to a downward direction (Figure 7). We refer to the distance as an interferometric distance (Figure 4). The interferometric distance of zero indicates that the two events have the common raypath between the buried source to the receiver . Due to the receiver spacing, the interferometric distance is not always zero (Figure 7). Therefore, we define the threshold for the interferometric distance. The receivers with the interferometric distance less than the threshold value (1 km in this case) are assumed to have the common raypath. This threshold corresponds to the determination of the size of the first Fresnel zone around the stationary points. Note that the interferometric distance projected from the first Fresnel zone is dependent on the position of reflector and the receiver geometry. We iterate this procedure by changing the fixed receiver and evaluate all combination of crosscorrelation traces.
We remove the crosscorrelation traces which do not have a stationary phase source for objective reflections and apply PSDM using 2489 traces. Figure 8 shows the result of subsurface image considering stationary phase. One can see that the signal-to-noise ratio at the near offset improved, and the artifact events are suppressed. We conclude that evaluating stationary phase by raypath calculation and removing unwanted crosscorrelation traces improve reflection images in SI.
4. Extracting Reflected Waves of Oceanic Crust Surface from OBS Records in the Nankai Trough
We apply this analysis to the field passive-seismic data. The field data consists of the OBS records deployed in the Nankai Trough to observe local earthquakes. This dataset was originally obtained by JAMSTEC [11, 12] for the earthquake observation as well as for aftershocks of the September 2004 intraplate earthquake ruptured around this survey area [20]. The 28 OBS are 3 dimensionally deployed at the approximately 100 km-squared area (Figure 9(a)). The recorded period was for approximately 3 months from March 2005, and we used the data during 12 days in March. We extract the reflection events generated at the oceanic crust surface using crosscorrelation. Within the 653 local earthquakes detected in this duration, we extract 6 earthquakes whose magnitudes are larger than 3.0 with >20 km in depth (Figure 9(b)) because these earthquakes contain strong energy. We choose deep earthquakes since the P-wave energy is dominated in the vertical component of the records, and the arrival time of S-waves is late enough to be separated from P-waves. These earthquakes are localized at south-west of the survey area (Figure 9(a)).
(a)
(b)
We show the 168 traces from the 6 events aligned by the epicentral distance and corrected origin time with the earthquake nucleation time (Figure 10). The dominant frequency was5 Hz. One can see two events appearing with different propagating velocity on the seismic record. Since the first arrival is dominated in the vertical component (Figure 10), we assumed them as the P-waves and the second arrivals as the direct S-waves. We muted the amplitude of the S-waves as well as subsequent records in order to analyze only P-wave events.
(a)
(b)
Figure 11(a) shows the reflection profile across the Nankai Trough (Moore et al. [9]; Line IL in Figure 9(a)). The oceanic crust surface in this area is located from 7 km to 10 km in depth and dipping landward direction (dashed line in Figure 11(a)). We extended this structure (landward-dipping structure with the angle of 5.5 degree) perpendicular to the 2D cross-line (Line IL in Figure 9(a)) and obtain the 3D dipping structure used for raypath calculation (Figure 11(b)).
(a)
(b)
We estimate the raypaths of the sea-surface reflection (pP) and the sea-surface multiple reflection (pPp) as in the case of simulation study. To account for the difference of each elevation of OBS, we define the interferometric distance as the distance between the receiver and the point where pPp passes through the horizontal plane constructed by the receiver to a downward direction. We defined the threshold of the interferometric distance as 5 km. We removed the crosscorrelation traces without the stationary phase sources and obtained 55 crosscorrelation traces. Since we calculate the raypath, we can estimate the reflection point at the given structure (Figure 12). Due to the source localization, only a part of combination of the receivers is selected.
We applied the 3D PSDM to the crosscorrelation traces. We extended the tomographic velocity estimated by Nakanishi et al. [8] perpendicular to its survey line (Line KR9806 in Figure 9(a)) and used it as a velocity model for PSDM. The expected reflection points are sparsely distributed in the 3D region (Figure 12). Therefore, we spatially stacked the 3D imaging result and obtained the pseudo-2D profile.
We show the pseudo-2D profile projected to Line KR9806 (Figure 13(a)) running perpendicular to the trough axis. The result shows that the oceanic crust surface is imaged at the same depth of the active-source reflection profile of Line IL (arrows in Figure 13(a)). For a comparison, we show the PSDM result using all crosscorrelation traces projected to Line KR9806 (Figure 13(b)). Since it is difficult to detect the oceanic crust surface in the profile (Figure 13(b)), stationary phase evaluation improves the quality of the imaging result (Figure 13(a)).
(a)
(b)
(c)
We show the pseudo-2D profile projected to Line S3 (Figure 13(c)) running parallel to the trough axis. We convert the depth axis of our imaging result to the time axis using the migration velocity model. The imaged oceanic crust surface (arrows in Figure 13(c)) is discontinuous around the horizontal distance of 20 km due to the migration aperture and the sparse distribution of the reflection points (Figure 12). However, we can observe the dominated amplitudes at the same two-way time of the oceanic crust surface in Line S3 [10] as well as the local bulge of the oceanic crust surface.
5. Conclusion
We use the stationary phase interpretation to obtain high-quality imaging results in the localized source distribution. We estimate the raypath of two reflection events which are sea-surface P-wave reflection and sea-surface multiple P-wave reflection. We choose the crosscorrelation traces which is expected to produce objective reflections due to the stationary phase sources using the estimated raypath. We show the numerical modeling result to check the validity of this method. Furthermore, we use Ocean Bottom Seismogram (OBS) records which observe localized earthquakes. We show that choosing the crosscorrelation traces by stationary phase evaluation improves the quality of the imaged reflection boundary of the oceanic crust surface. This processing technique has a possibility to monitor the Nankai seismogenic fault without the active sources and higher resolution than using teleseismic records.
Acknowledgments
The OBS data were acquired by JAMSTEC. S. Minato is grateful for a support by JSPS Fellows (212666). T. Tsuij is grateful for a support by Grant-in-Aid for Scientific Research on Innovative Areas (21107003). This work is partially supported by JST/JICA, SATREPS and Kyoto University GCOE.