Monitoring and forecasting of seismic wavefields in the subsurface

Monitoring seismic wavefields caused by induced seismicity in the subsurface is a difficult process. Ideally, it requires physical receivers in the subsurface, which is unpractical. Frequently, only measurements at the surface of the Earth are available, which give a limited amount of information about the subsurface. One way to improve the monitoring of the subsurface is through the use of virtual sources and receivers, which are not physically present but are created from the measured reflection data at the surface. This can be achieved through the use of the classical homogeneous Green's representation, however, this method requires two Green's functions measured on an enclosing boundary, which is an unrealistic requirement. Instead, a single-sided representation of the homogeneous Green's function can be used, where a focusing function, which is a wavefield that focuses from a single-sided boundary to a focal position in the subsurface without artifacts related to the internal multiples, is employed together with a Green's function. To obtain the Green's function and focusing function that are needed for this representation, the Marchenko method is used. This method employs reflection data, without free-surface multiples, at the surface of the Earth and an estimation of the first arrival, which can be modeled in a macro velocity model.

To test whether induced seismicity in the real subsurface can be monitored using the single-sided representation, synthetic data are first considered, which include a synthetic reflection response and macro velocity model. The Marchenko method is used in combination with these data to obtain the focusing functions and Green's functions that are required for the homogeneous Green's function representations. The classical representation and the single-sided representation of the homogeneous Green's function employ the Green's functions and focusing functions to obtain the homogeneous Green's function of the medium. The homogeneous Green's function is visualized by creating snapshots of the homogeneous Green's function and these snapshots are compared to a directly modeled reference wavefield. This demonstrates that the classical representation, when applied to data at an open acquisition boundary, yields significant artifacts in the results, while the single-sided representation obtains accurate results. It is also shown that the radiation pattern of a double-couple source can be included in the retrieval of the homogeneous Green's function. The synthetic reflection data are truncated by limiting the offsets and sampling distance and applying attenuation to simulate field conditions. These truncations show that the single-sided homogeneous Green's function contains artifacts and lacks physical events if the reflection data are not ideal. 2D field reflection data and a macro velocity model from the V\o ring basin are considered and pre-processed to account for these truncations. The classical and the single-sided homogeneous Green's function representation are both applied to the field data and the results show that the retrieval of the homogeneous Green's function is possible for 2D field data using point sources while employing the single-sided representation. The results of the classical representation contain a large amount of errors. It is also shown that a homogeneous Green's function can be retrieved that has a virtual source with a double-couple radiation pattern.

Next, the application of the single-sided representation is considered in greater detail. The representation is used to forecast a wavefield in the subsurface as well as to monitor a wavefield in the subsurface. For the monitoring of the wavefield, it is assumed that a physical source in the subsurface causes a wavefield which is measured at the surface of the Earth. The Marchenko method is used to create virtual receivers inside the subsurface, which are used in combination with the physical measurement in the single-sided representation. This is a one-step process, because the Marchenko method is only used to create the virtual receivers. The single-sided representation of the homogeneous Green's function requires the source wavelet to be symmetric in time, which is unlikely for physical sources. Hence, a different single-sided representation can be used, which retrieves the causal Green's function and does not require a symmetric source wavelet. The single-sided representation of the causal Green's function can retrieve a majority of the correct events, however, the results contain anti-symmetric artifacts when the physical source is located above the virtual receiver. To forecast a wavefield in the subsurface, given a specific source configuration, the single-sided representation of the homogeneous Green's function can be used. In this case, a two-step process is applied, where both the source and the receiver in the subsurface are created by the Marchenko method and are therefore both virtual. After the homogeneous Green's function is obtained, it can be convolved with a non-symmetric wavelet. To demonstrate the difference between the one-step monitoring process and the two-step forecasting process, 2D synthetic reflection data are utilized. For the source configuration, a rupture plane is considered, which is modeled by superposing and time-shifting point sources, which contain a double-couple radiation pattern and are all scaled differently to simulate the heterogeneity of the rupture plane. The total wavefield created by this rupture plane is monitored using the single-sided representation of the causal Green's function. There are anti-symmetric artifacts present in the result, related to each point source, however, the correct wavefield is retrieved above the shallowest source location and below this source location after the first arrivals of all sources. The single-sided representation of the homogeneous Green's function is applied to forecast a virtual rupture plane, by retrieving the homogeneous Green's function for each source separately. The retrieved homogeneous Green's functions are transformed to causal Green's functions, shifted in time and superposed to forecast the total wavefield, which is free of the anti-symmetric artifacts at any depth. Both the monitoring approach and the forecasting approach are tested on 2D field data and the retrieved wavefields show similar results as were seen when the synthetic data were used. When the total wavefield is forecasted, there are no anti-symmetric artifacts present and when the wavefield is monitored, there are artifacts, however, they are only present in part of the result, below the sources before and during the first arrival of each source.

To test the application of the single-sided representation in 3D, a 3D implementation of the Marchenko method is required. The implementation is straightforward from a theoretical standpoint, as the surface integrals are performed over two dimensions instead of just one. The practical implementation is more difficult, however. The Marchenko method requires that the reflection data are well sampled in both space and time for sources and receivers, hence, the 3D reflection data are of a large size. As a result, not only a large amount of storage space is required, but the loading time of the reflection data is high, both of which are unpractical for efficient computation. We limit these problems by pre-transforming the reflection data to the frequency domain and compressing the data using floating point arrays, which reduces the storage space and loading time. Two datasets are considered, one modeled in a simple four layer model and the other in a subsection of the complex 3D Overthrust model. For both models, a Green's function inside the medium is retrieved, using a first arrival in the Marchenko method that was modeled in the exact medium, and compared to a reference Green's function that was directly modeled. The results for both models are accurate for the single Green's function. Next, imaging is performed for the models, however, instead of modeling the first arrivals, they are estimated using an Eikonal solver, because the modeling time of all the first arrivals is too high. The results of the imaging using the Marchenko method are compared to the results of conventional imaging, which demonstrates that artifacts, related to the internal multiples, are attenuated.

The 3D implementation of the Marchenko method is used to retrieve the Green's functions and focusing functions in 3D using 3D synthetic reflection data modeled in the Overhtrust model. The classical homogeneous Green's function representation and the single-sided representation of the causal Green's function and the homogeneous Green's function are all applied using these data, for three different combinations of a virtual source and a virtual receiver. The results are compared to a directly modeled wavefield, which shows that the result obtained by using the classical representation is contaminated by artifacts and lacks physical events. The result of the single-sided representation of the causal Green's function contains anti-symmetric artifacts related to the focusing function when the virtual receiver is located below the virtual source. The result of the single-sided representation of the homogeneous Green's function shows a good match to the reference result. The single-sided representation of the homogeneous Green's function is also applied using an Eikonal solver to obtain the first arrival that is required for the Marchenko method. The homogeneous Green's function that is obtained in this way shows a small decrease in quality for the result, however, this approach is more computationally feasible. The single-sided representation is used in combination with the Eikonal solver to retrieve a large amount of virtual receivers, so that the propagation of the wavefield in the subsurface can be visualized in time through the use of snapshots. This reveals that the part of the wavefield that is traveling at angles that are close to the normal of the surface is retrieved properly, while the part of the wavefield that is traveling at greater angles to the normal is reconstructed with less accuracy. This lack of proper retrieval is caused by the limited aperture of the reflection data. A rupture plane in 3D is considered and constructed in a similar way as is done for the 2D synthetic data. Point sources are used to model wavefields, which are time-shifted and superposed, however, to further represent the heterogeneity of the rupture plane, each wavefield is modeled using an unique causal wavelet. Both monitoring, using the single-sided causal Green's function representation, and forecasting, using the single-sided homogeneous Green's function representation, are performed on the rupture plane configuration. The two-step forecasting approach yields accurate results, for a given distribution of sources. The one-step monitoring approach retrieves accurate results above the shallowest source location, however, the result contains artifacts at the locations below the shallowest source, before and during the first arrival of each source.