Regional exploration of continental margin

from ultra-long offset OBN data by slope tomography and FWI

In 2001, the JAMSTEC Institute carried out the first dense OBS (Ocean Bottom Seismometers) survey in the frame of the French-Japanese collaborative SFJ (Seize France Japan) project with the aim to assess FWI for deep crustal exploration. JAMSTEC deployed 100 OBSs along a 100-km long profile across the Tokai segment of the eastern Nankai subduction zone (Figure 1). After some preliminary applications of FWI on this data set ( Dessa et al., 2004; Operto et al., 2006), A. Gorszczyk recently revisited this dataset with a more careful quality control of the initial velocity model built by first-arrival traveltime tomography and a more complex Laplace-Fourier FWI workflow combining data-driven frequency, traveltime and offset continuation strategies ( Gorszczyk et al., 2017). The inversion was performed in the 1.5-8.5 Hz frequency band, where virtual low frequencies have been generated through aggressive time damping of the data. The FWI velocity models reveals unprecedented details of the deep structure of the Tokai segment of the eastern-Nankai subduction zone (Figure 2). The reliability of the velocity model was checked by seismic modeling of the OBS data (Figure 3) and Kirchhoff depth migration of co-incident towed-streamer data (Figure 4 and Gorszczyk et al., 2019).

Figure 1: The left inset shows the geodynamical context and the rupture zones of the big earthquakes in the Nankai zone which occurred during last century. The most oriental segment in the Tokai area remains un-ruptured. Therefore, a soon-to-come big earthquake is expected in this area. Right inset: map of the SFJ-OBS survey across the Tokai segment of the eastern Nankai subduction zone. The red dots denote the position of the OBSs. Note the Zenisu ridge. Some analogous topographic reliefs probably entered the subduction zone. The figure shows the anatomy of ultra-long offset data (offset max. 120km) recorded by one OBS located nearby the northern end of the profile. The recorded wavefield is dominated by diving waves and post-critical reflections. Some weak first arrivals at 30-40km offset are head waves trapped along thrusts affecting the backstop. The acronyms PmP and Pn denote the reflection from the base of the subducting slab and the refracted waves from the Moho. FWI of these data in the 1.5-8.5 Hz frequency band produced the velocity model shown in Figure 2.

Figure 2: (Top) P-wave velocity model obtained by frequency-domain acoustic Laplace-Fourier FWI of the SFJ-OBS dataset. The derivative of the velocity model along the structural dips is superimposed to highlight the short-scale structure mapped by FWI. (Bottom) Same figure with superimposed prestack depth migrated image from coincident towed-streamer data. From right to left, the subduction front (85km distance), the accretionary wedge, the Tokai thrust (65km distance), the backstop with underthrusting of crustal sheets. See Gorszczyk et al., 2017 and Gorszczyk et al., 2019 for more details. The figure is borrowed from the poster presented at the 2018 AGU fall meeting by Sambolian et al. (2018) where first-arrival slope tomography was applied to build the initial model of FWI.

Figure 3: Fitting ultra-long offset OBS data. (a) Synthetic seismograms computed in the initial velocity model built by first-arrival traveltime tomography. The insets show interleaved recorded and simulated seismograms. Simulated amplitudes are inaccurate because energy partitioning at interfaces has not been modeled in the smooth velocity model built by first-arrival traveltime tomography: amplitudes of the simulated first-arrivals are overestimated, while those of the simulated late arrivals (i.e., reflections) are underestimated. These amplitude mismatches are corrected as FWI reconstructs short-scale features of the subsurface. The red and green curves show the first-arrival traveltime picks and the simulated counterparts. (b) Same as (a) for the seismograms computed in the FWI model. Both traveltimes and amplitudes are matched as short-scale features have been injected by FWI in the subsurface medium. (c) Recorded data. (d) Comparison between the AVO curves of the three sets of seismograms extracted from a narrow window centered on the early arrivals. Figure from Gorszczyk et al., 2017.

Figure 4: FWI velocity model as background model for prestack depth migration of towed-streamer data. Common-Image gathers generated by Kirchhoff prestack-depth migration of vintage towed-streamed data (4.5km long streamer) using as background velocity model: (a) A basic velocity model built by velocity analysis. (b) Same as (a) after refinement of the velocity model by slope tomography. (c) the FWI model developed by Gorszczyk et al., 2017. (d) Same as (c) after refinement of the FWI velocity model by slope tomography. Although the FWI velocity model was inferred from wide-angle arrivals recorded by the OBSs, it is accurate enough to flatten reasonably well the common-image gathers generated from short-spread reflections (see Gorszczyk et al. (2019) for the joint analysis of the towed-streamer and OBS data of the SFJ experiment).

We are going to benchmark FWI with the GO_3D_OBS geomodel representative of subduction zone. A 102 km x 20 km x 28.3km target of this geomodel is shown below in (a). Our goal is to design optimal FWI workflow for crustal scale imaging and study the footprint of sparse OBN acquisition with this benchmark. We will perform illumination analysis in the theoretical framework of diffraction tomography to check this footprint (b).

Large-scale finite-difference frequency-domain simulation in a 102 km x 20 km x 28.3 km of the GO_3D_OBS subduction-zone velocity model. Velocities range between 1.5km/s and 8.6km/s. Frequency is 3.75 Hz. Number of degrees of freedom if 67.5 millions. Eighty nodes of the occigen supercomputer of CINES (https://www.cines.fr/calcul/materiels/occigen) are used. The linear system is solved with the massively parallel MUMPS multifrontal solver (https://mumps-consortium.org) using the block low rank function with a compression threshold parameter of epsilon=1e-5. Elapsed time for LU factorization is 869 s. Elapses time to compute 130 wavefields is 34 s. I am still reading in papers that sparse direct solvers cannot be used for 3D applications. Please check our papers.