Downhole distributed acoustic sensing reveals the wavefield structure of the coastal microseisms

Ocean-generated seismic waves are omnipresent in passive seismic records around the world and present both a challenge for earthquakes observations and an input signal for interferometric methods for characterisation of the Earth’s interior. Understanding of these waves requires the knowledge of the depthdependence of the oceanic noise at the transition into continent. To this end, we examine 80 days of continuous acquisition with Distributed Acoustic Sensor (DAS) system deployed in two deep boreholes near the south-eastern coast of Australia. The data has excellent Signal-to-Noise Ratio (SNR) in a range from 0.03Hz to over 100Hz. By analysing the seismograms and correlation with wave climate, the DAS response are confidently decomposed into the microseisms generated by swell from remote storms (∼0.15Hz) and local winds (between 0.3Hz and 2Hz), and strong body waves energy from large surf break at the coast (from 2Hz to 20Hz). The depth dependence of the microseims allows for robust normal modes analysis of the Rayleigh waves with only one borehole. The results of this analysis agree with the data from conventional dense seismological arrays. Overall, we found that the link between the amplitudes at each channel ∗Corresponding author Email address: stanislav.glubokovskikh@curtin.edu.au (Stanislav Glubokovskikh) Preprint submitted to Earth and Planetary Science Letters September 24, 2020 along the borehole and wave climate is so strong and stable that with sufficient amount of training data, the passive seismic records on downhole DAS may be used for high-precision monitoring of both formations surrounding the borehole and remote storms in the ocean.


Introduction
Motion of water in stormy seas induces omnipresent ambient seismic wavefield in the frequency band from ∼0.01Hz to ∼20Hz (Webb, 1998). These frequencies overlap with the frequency band of earthquake signals critical for seismological observations by Ocean Bottom Stations (OBS), and hence has often 5 been regarded as unwanted noise (e.g., Webb, 1998). Over the years, a link has been established between some components of the ocean-generated microseisms and wave climate in stormy seas (McCreery et al., 1993;Bromirski et al., 1999;Bromirski and Duennebier, 2002;Aucan et al., 2006), which allows the seismic observations to be used for monitoring the sea conditions and the rock proper-10 ties in the vicinity of the sensor. Furthermore, rapid development of interferometric methods expanded the research focus to the subsurface characterisation using passive seismic data (e.g., Nakata et al., 2019). However, inversion of the passive seismic observations often assumes stationary isotropic wavefield (Wape-peak frequency of the ocean waves and a much stronger Double Frequency (DF) microseisms at twice that frequency. The SF energy is related to the transfer 25 of water surface oscillations into elastic energy through a direct interaction of the water waves with ocean floor (Hasselmann, 1963), and hence the intensity of these signals drops exponentially with the water depth. Typically SF microseisms become undetectable a few kilometeres into the continent Haubrich and McCamy (1969); Bromirski and Duennebier (2002). In turn, DF microseims 30 originate from nonlinear interaction of ocean wave trains of a similar frequency and opposite direction (Longuet-Higgins and Jeffreys, 1950;Hasselmann, 1963;Kibblewhite and Wu, 1991). The resulting oscillations penetrate almost lossless to the sea bottom and couple into leaky Rayleigh wave modes, which dominate the low-frequency part of the seismic noise on land. Amplitudes of the seismic 35 waves depend on the storm parameters in a large source region and the seismic properties of the ocean bottom (Webb, 1992;Tanimoto, 2007;Ardhuin et al., 2013;Gimbert and Tsai, 2015).
Overall, the mechanisms generating the seismic waves at the ocean bottom are relatively well understood, although a quantitative prediction of the mi-40 croseisms remains a challenging problem. At the same time, transition of the microseisms into the continent is poorly understood. Conventional seismological arrays are too sparse for the ambient seismic wavefield (∼10km). As a result, we only have coarse scale distribution of the velocities and intensity of the propagating coherent signals from the kinematic and polarization analy-ent seismic wavefield at the coast remains uncertain. Numerical simulations of this process are computationally expensive and are based on poorly constrained models of the subsurface and seismic source, and are able to capture only wellestablished features of the microseisms such as dependence on water depth or presence of the soft sediments at the bottom (Levchenko et al., 2011;Ying et al., 2014). More detailed information may be obtained from measurement of a DAS system in a sufficiently deep borehole. With the rapid development of the borehole seismic monitoring using DAS (e.g., Correa et al., 2017;Egorov et al., 2017), the instrumented wells may be employed for ambient seismic monitoring.
of the DAS measurements in deep boreholes for monitoring both subsurface properties and wave climate.

The data set
Ambient seismic wavefield was recorded from 21 October 2018 to 07 December 2019 by DAS systems in CRC-2 and CRC-3 wells (see fig. 1b). These wells Several closely spaced wells ∼ 1600m deep have been drilled and instrumented 95 with modern DAS systems. Seismic properties in the vicinity of the wells are available from sonic logs and numerous 3D borehole and surface seismic surveys with active sources (Glubokovskikh et al., 2016;Egorov et al., 2017Egorov et al., , 2018. Our analysis of the oceanic noise relies on remarkable performance of the modern DAS systems in a broad frequency range, especially at low frequencies 100 (Lindsey et al., 2017). The seismic measurements by DAS estimate a rate of the axial deformation of the optical fibre zz induced by a seismic wave. To this end, a DAS interrogator measures the phase difference between laser pulses backscattered from adjacent segments of the optical fibre. Appendix Appendix A presents an original derivation of the measurement principle implemented 105 in iDASv3™ and iDASv2™ systems manufactured by Silixa (Elstree, Hertfordshire, UK, http://www.silixa.com), which are deployed in CRC-3 and CRC-2 wells, respectively. CRC-2 has a previous generation DAS system, which uses a standard single-mode optical communication fibre deployed on the production tubing. As a result, CRC-2 DAS records feature much higher instrumentation 110 noise than the data from CRC-3, which has the cable cemented behind the casing ( fig. 2). Note that the increased noise in the uppermost ∼270 m of both wells is due to surface casing, which causes intense mechanical reverberations of seismic signals. In the following, we focus on the CRC-3 records. The seismic signal on DAS deployed in CRC-3 may be converted into the particle displacement along the fibre u z using factor F DAS (k z ) derived in appendix Appendix A. The factor depends on the frequency ω of seismic signals and parameters of instrumentation: distance along the fibre between the reflection points of the backscattered pulse known as gauge length L 0 , duration of the pulse τ , and speed of light c 120 in the fibre material. For small k z , a projection of the seismic wave vector on the fibre, equation (A.10) for the DAS response reduces to Equation (1) shows that for most types of seismic waves, the DAS sensitivity drops with decreasing frequency relatively quickly, as ω 2 . However, the phase interferometry approach implemented in iDASv3™ provides sufficient SNR start-

125
ing from a few millihertz. Pevzner et al. (2020a) showed that eq. (1) holds true for the DAS system in CRC-3 and the data still provides clear records of teleseismic earthquakes. In the following sections we divide the spectral amplitude of the DAS response by ω, to make DAS data more consistent with the conventional seismological records, which measure either particle displacement or its 130 temporal derivatives, velocity and acceleration.

DAS response versus wave climate
The DAS record in CRC-3 ( fig. 2a) shows an abundance of ambient seismic energy at all depth levels. The Otway site is located in the area of active farming; previous analysis (Pevzner et al., 2020b) showed DAS records containing 135 distinct anthropogenic noise: monochromatic signals from various farming machinery, impulse sources, and even routine movement of cattle. Yet, the records are dominated by low-frequency surface waves with an apparent period ∼2s propagating horizontally and frequently repeated body waves and their scattered modes with an average period below 0.5s propagating at an angle to the 140 borehole receiver array. In the following we show that these signals are generated by the ocean.   5. from 2Hz to 20Hz: frequently recurring body waves; 6. >20Hz: vibrations of the recording unit caused by local wind.
In group 1, microseisms, if any, are weaker than the instrumental noise. Between 50mHz and 120mHz, although the correlation coefficients gradually increase, the intensity of the primary microseisms at ∼70mHz is relatively low and a peak at <0.1Hz only seldom appears in the data. Group 6 contains highfrequency clutter, which is unlikely to come from the ocean. Indeed, seismic Q in the area is about 100 or smaller (Pirogova et al., 2019), and hence these waves would attenuate substantially at distances on the order of several kilome-165 ters. Instead, signals at frequencies over 20Hz may be caused by wind-induced mechanical vibrations of the hut housing the recording unit. We omit a further discussion of groups 1, 2 and 6, because they are affected by characteristics of the acquisition system, rather than the ambient seismic field. In the remaining frequency range, 0.1Hz -20Hz, the DAS data stably provides high SNR as 170 evident from high correlation with features of the wave climate.
The correlation peak for H S at ∼150mHz corresponds to twice the peak ocean wave frequency and hence represents the classical DF microseisms. These seismic signals are studied extensively in the literature and our data agree with numerous reports from around the world (e.g., Webb, 2007 by signals with frequency equal to twice the peak frequency of ocean waves (Webb, 1998;Nakata et al., 2019). Commonly accepted models (e.g., Pierson Jr. and Moskowitz, 1964) predict that with increasing wind speed, the short-period ocean waves reach maximum intensity quickly, while progressively lower-frequency waves become stronger with increasing wind speed. Therefore, 220 the wave height energy and hence microseisms spectrum decrease rapidly with increasing frequency. It is believed that above 0.3Hz, microseisms have a worldwide constant Holu spectrum (McCreery et al., 1993;Bromirski et al., 1999).
But the nonlinear generation of the microseisms requires interaction of ocean waves of a similar frequency and opposite direction. Hence the strength of 225 DF microseisms at each frequency, expressed as the power spectral density of the vertical displacement I(ω), depends on the distribution of the intensity of ocean waves A(ω/2, θ) over azimuth θ. With a few reasonable assumptions, the relationship has the following simple form (see equation (7) in Webb, 1992) The buoy data in the area shows the tendency of the high-frequency ocean 230 waves to form pairs of sufficiently strong opposing wavetrains (see fig. 6). Predicted frequency spectra of the microseisms are in good agreement with the theoretical predictions of eq. (2). This suggests that the peak frequency in the DAS response may be a consequence of the peculiar wave climate at the Victorian coast. The frequency shift may be further increased due to factor k Z in 235 the DAS response (see eq. (1)), which increases with frequency.

Variation of strain with depth in the microseisms
Spectrogram in fig. 7a shows a typical intensity distribution of the DAS measurements along the CRC-3 borehole: first, amplitudes decrease towards a pronounced minimum at ∼700m depth; then the amplitudes increase steadily  tra. At the bottom, power spectral density in the DAS data is plotted against the theoretically predicted I(ω) using eq. (2), which is normalized to the measured DAS response at 0.6Hz to facilitate the comparison. The directional spectra correspond to WBAST1 buoy (see fig. 1a). modes controlled by the parameters of the source region: frequency of the ocean waves and bathymetry (e.g., Ewing et al., 1957). However, transition of the 250 excited normal modes into the continent is a much more complicated process, which depends on onshore geology. Thus, an accurate estimate of the magnitude of the microseims would require the knowledge of the location of the source region and the wave climate in that region along with the bathymetry and elastic properties of the sediments and crustal rocks (Hasselmann, 1963;Webb, 255 1992;Tanimoto, 2007;Gimbert and Tsai, 2015;Ardhuin et al., 2013). Such a comprehensive modelling is impractical with the available data and lies beyond the scope of this paper. However it is still possible to find a best-fit combination of the normal modes that may be supported by the Otway subsurface properties distribution.

260
To this end, we assume a 1D layered subsurface structure, where the upper >0.3Hz, while in our data, the peak occurs at ∼0.6Hz.
Another feature of the measured DAS intensity is a clear dependence on the stiffness of surrounding rocks. In fig. 7b we may clearly identify the peaks and 280 troughs corresponding to stiff and soft geological formations. Similar behaviour was observed in the analysis of amplitudes along CRC-3 of compressional waves from remote earthquakes (Pevzner et al., 2020a). Unlike the stress traction, vertical strain is discontinuous at the lithological boundaries, and hence its finescale variability reflects the variation of the stiffness of rocks along the borehole. Compared with traditional means for passive seismic monitoring, we see two lim-295 itations of our data set. First, we were unable to obtain confident characteristics of very low frequency signals <40mHz,including infragravity oscillations known as Earth's hum (e.g., Webb, 2007). Nevertheless, we believe that the acquisition settings may be adjusted to provide a stronger response in millehertz range, because our data set contains clear signal from numerous teleseismic earthquakes 300 and primary microseisms at these frequencies. Moreover, a similar DAS has been successfully used to monitor very slow strain due to changes of hydraulic head in boreholes (Becker et al., 2017). Second, DAS provides only axial strain, and thus precludes a polarization analysis of the recorded signals. Moreover, DAS has strong directivity pattern -its sensitivity deteriorates rapidly with the

Conclusions
This study demonstrated a substantial potential of passive seismic records obtained with a downhole distributed acoustic sensing system for studying the 340 ocean-generated seismic signals. The data acquired at the CO2CRC Otway Project site has high signal-to-noise ratio for frequencies from 10mHz to 100Hz. at variance with a commonly-accepted hypothesis of a worldwide constant Holu spectrum. However, this discrepancy is fully explained by the analysis of local winds, which induce high frequency wave trains moving in opposite directions, whose nonlinear interaction becomes an efficient generator of microseisms.
The microseisms have a very specific amplitude variation with depth at all 360 frequencies and for variety of the source regions in the ocean: the surface waves change polarity at the same depth where rock stiffness changes dramatically.
Such stability suggests a strong dependence of microseisms structure on the geological section in the immediate vicinity of the borehole. Normal modes analysis shows that observed amplitude variation with depth requires relatively strong 365 contribution from higher Rayleigh modes even at low frequencies (∼100mHz).
A more comprehensive analysis of the ocean-generated seismic signals re-quires a better control of the source contribution into the response on distributed acoustic sensing data. However, the link between the amplitudes at each channel along the borehole and wave climate appears so strong and stable 370 that with sufficient amount of training data, the passive seismic records may be used for high-precision monitoring of both formations surrounding the borehole and remote storms in the ocean.

Acknowledgements
The Otway Project received CO2CRC funding through its industry mem- where t is a two-way travel-time. This travel-time is related to the distance to a reflection point z = ct /2, where c is a speed of light in the fibre (typically, 1m corresponds to 10ns). With some reasonable assumptions, the backscattered 400 field may be modelled as a linear function of the fibre reflectivity r(z) where * denotes the convolution operator. We assume an idealised laser source that generates a monochromatic pulse of temporal frequency ω and width τ . In the presence of time-dependent displacements u(z, t), the reflectivity becomes a function of v(z, t) = ∂u(z,t) ∂t , a relative particle displacement velocity at different parts of optical fiber at a moment in time t. To retrieve v(z, t) induced by a seismic wave, DAS systems must accurately measure the temporal varia-410 tion of the OTDR intensity I(t ) = E 2 (t ). The main challenge for a practical implementation of the OTDR-like approach is generation of stable laser pulses, otherwise temporal fluctuations of s(t) would obscure any signal due to v(z, t).
In communication networks, OTDR uses incoherent pulses, and its output is insensitive to variations of the optical phase along fiber. But even with the 415 best existing laser sources, the coherent OTDR can provide only qualitative estimates of the particle velocity v (?).
In all the quantitative DAS approaches, the effect of fibre movement on the reflected optical field is estimated using interferometric techniques (Hartog, 2017). The main idea is to compare pulses reflected at adjacent points on the 420 fibre. Interference of these pulses is controlled by the reflection coefficients and their relative movement, but a specific form of the DAS response depends on the instrumentation configuration. Figure A.8 illustrates an implementation of this approach using an interferometer with a delay line of length 2L 0 . In the interferometer, the pulses are duplicated into the straight and delay lines, and 425 the delayed pulse interferes with a pulse that was reflected from a point offset by L 0 and then passed through the straight line. Thus, L 0 is often called a gauge length, a key parameter that controls the spatial resolution of the DAS measurements.
First, we derive an optical response to a uniform displacement velocity shows that the absolute value |v 0 | defines the beat frequency of the photon counts. For example, the magnitude of the oscillations I Ω of ∂I12(t) ∂t may be expressed as where the time dependence is converted into a distance along the fibre using a For slow variation of the displacement velocity field along the fibre ( L 0 ), 460 eq. (A.4) shows that DAS response is proportional to an axial strain rate of the fibre˙ zz = ∂v(z,t) ∂z . Equation (A.4) shows that the effect of the fibre deformation bears an overprint of the spatial variations of the reflectivity along the fibre. The approach to removal of the effect of r(z) constitutes the main difference between iDASv3™ 465 and iDASv2™ systems. The reflectivity in a conventional single-mode fibre is associated with Rayleigh scattering on low-contrast fluctuations of the optical refraction index, an inherently random process. Strictly speaking, elimination of this random effect would require an ensemble of the DAS measurements for the same velocity field v(z, t) but different fibres followed by analysis of the 470 ensemble average A(z) . For Rayleigh scattering, eq. (A.4) simplifies to To illustrate the effect of random fibre reflectivity, ? simulated numerically the DAS response corresponding to the movement of a 40m segment of fibre (no deformation occurs inside the segment), for τ = 10-100ns and the gauge length L 0 = 10m ( fig. A.9). If the pulse width is small, agreement with the 475 theory is good. But for a typical pulse width, this effect causes significant errors in estimates of the velocity field v(z).
This observation explains a general approach to choosing the pulse width. We can increase SNR and reduce the distortions simultaneously by using an engineered fiber with regularly spaced high reflectivity markers. With L 0 = 2τ , such a design prevents an overlap between the reflected pulses and minimises the effects of fluctuations of the fibre reflectivity. This idea is implemented in In Fourier domain with a spatial frequency k z , the spectrum of the DAS response from eqs. (A.7) and (A.8), F E (k z ), can be expressed via the spatial spectrum of the displacement velocity field along the fibre F v (k z ) as where we assumed that the pulse is a unit rectangular function. The DAS receiver function F DAS is F DAS (k z ) = ω · |sinc(k z τ /2) · sin(k z L 0 /2)| . (A.10) Because v = ∂u(t) ∂t , conversion of the DAS measurements to displacements will feature also a factor ω (see 1).
We can conclude that, unlike DAS with a conventional fibre, the spectrum of 505 the DAS response with an engineered fiber is subjected to aliasing similar to an array of geophones. Figure A.10 shows the spectra F DAS for typical parameters of DAS measurements. Clearly, an engineered fibre may be used until the cut-off frequency if an appropriate anti-aliasing filter is deployed, while a conventional fibre has no distortions in a wider frequency range, but its sensitivity drops with 510 increasing frequency.
Design of the engineered fiber dictates the gauge length. However, it is possible to synthesise a long optical gauge length by shifting and stacking the measurements obtained with a short gauge length, as follows from eq. (A.8). In