Strain to ground motion conversion of distributed acoustic sensing data for earthquake magnitude and stress drop determination

The use of distributed acoustic sensing (DAS) presents unique advantages for earthquake monitoring compared with standard seismic networks: spatially dense measurements adapted for harsh environments and designed for remote operation. However, the ability to determine earthquake source parameters using DAS is yet to be fully established. In particular, resolving the magnitude and stress drop is a fundamental objective for seismic monitoring and earthquake early warning. To apply existing methods for source parameter estimation to DAS signals, they must first be converted from strain to ground motions. This conversion can be achieved using the waves’ apparent phase velocity, which varies for different seismic phases ranging from fast body waves to slow surface and scattered waves. To facilitate this conversion and improve its reliability, an algorithm for slowness determination is presented, based on the local slant-stack transform. This approach yields a unique slowness value at each time instance of a DAS time series. The ability to convert strain-rate signals to ground accelerations is validated using simulated data and applied to several earthquakes recorded by dark fibers of three ocean-bottom telecommunication cables in the Mediterranean Sea. The conversion emphasizes fast body waves compared to slow scattered waves and ambient noise and is robust even in the presence of correlated noise and varying wave propagation directions. Good agreement is found between source parameters determined using converted DAS waveforms and on-land seismometers for both P and S wave records. The demonstrated ability to resolve source parameters using P waves on horizontal ocean-bottom fibers is key for the implementation of DAS-based earthquake early warning, which will significantly improve hazard mitigation capabilities for offshore earthquakes, including those capable of generating tsunami.


Introduction
In recent years, the implementation of distributed acoustic sensing (DAS) for seismological purposes is rapidly expanding, for both on-land (e.g., Zhan, 2020;Fang et al., 2020;Jousset et al., 2018;Yu et al., 2019;Ajo-Franklin et al., 2019;Walter et al., 2020) and ocean-bottom (e.g., Sladen et al., 2019;Lior et al., 2021;Lindsey et al., 2019;Williams et al., 2019;Spica et al., 2020) applications. However, the use of DAS for several fundamental seismological tasks is yet to be fully established. These include earthquake location and source parameter (magnitude and stress drop) determination, both essential for harnessing DAS for earthquake early warning (EEW). But while applications of source location have been investigated in recent studies (e.g., van den Ende and Ampuero, 2020; Lellouch et al., 2020;Lindsey et al., 2017), the ability to infer source parameters has not been investigated in detail (e.g., Lellouch et al., 2020).
One of the major hindrances for source parameter determination using DAS stems from the measurement type: DAS produces strain or strain-rate recordings while source models (e.g., Brune, 1970;Madariaga 1976;Sato and Hirasawa, 1973) rely on ground motions, i.e., displacements, velocities or accelerations. Thus, the ability to invert for the source properties using conventional methods depends on the ability to reliably convert strain (rate) measurements to ground motions. The conversion between strain (rate) and ground motion has been demonstrated by various previous studies (e.g., Daley et al., 2016;Lindsey et al., 2020;Paitz et al., 2020;Wang et al., 2018;Lior et al., 2021;van den Ende and Ampuero 2020). This conversion can be accurately done by spatial integration if one seismometer is co-located with the fiber, along straight portions that have uniform coupling (van den Ende and Ampuero, 2020; Wang et al., 2018). When these conditions are unavailable, the common conversion approach consists of estimating the apparent phase velocities along the fiber and converting strain (rate) to ground velocity (acceleration) in the time (e.g., Wang et al., 2018) or frequency-(e.g., Lindsey et al., 2020) domains, using the following basic relations: where ,˙ , A, V and p are the strain, strain-rate, acceleration, velocity and apparent phase slowness along the fiber, respectively. These equations are valid for a single plane wave. Since different seismic phases travel at different velocities, and frequently in different directions, the apparent velocity required to convert strain (rate) to ground motions may vary rapidly, and a single time-invariant value may bias the analysis. In addition, velocities may vary along a fiber due to local velocity structure and fiber orientation variations. The ability to robustly convert DAS records to ground motions for signals containing varying phase velocities is key for harnessing DAS for early-warning and hazard mitigation purposes. Determining phase velocities via frequency-wavenumber (FK) analysis (e.g., Lior et al., 2021;Lindsey et al., 2020;Paitz et al., 2020) can be challenging since sufficiently long cable segments and time intervals are required to achieve adequate temporal and spatial frequency resolutions, in addition to delicate interpretation, as shown in later sections. To overcome this issue, a better-suited approach to retrieve phase velocities as a function of time should be sought. In this study, we propose a method for continuous apparent phase velocity estimation using semblance-based local slant-stack transform (e.g., Neidell and Taner, 1971;Taner et al., 1979;Shi and Huo, 2019). This technique, commonly applied in exploration seismology (e.g., Tatham et al., 1983), is used here to estimate phase velocities as a function of time, facilitating a time-dependent conversion of DAS strain (rate) signals to ground motion records. We validate this conversion method using synthetic signals and apply it to ocean-bottom DAS earthquake records. This paper is organized as follows. First, the slant-stack algorithm is presented and validated using synthetic waveforms. Then, the approach is used to convert earthquakes recorded by ocean-bottom DAS to ground motions. Finally, source parameters are determined by fitting DAS observations with a source model and compared with those determined using nearby seismometers. A comparison between the use of semblance-derived and FK-derived apparent velocities is presented throughout the paper.

Slant-stack transform for strain to ground motion conversion
The semblance-based local slant-stack transform (e.g., Taner et al., 1979) is used to resolve apparent phase velocities as a function of time. This array-based technique measures the coherency (semblance) of plane waves recorded by several adjacent sensors. At each instance in time, a range of slowness values is tested to identify that with the maximum semblance. For each slowness, semblance is calculated by aligning the time series recorded at different locations with respect to the middle station of a linear array (Shi and Huo, 2019): where 2L + 1 is the number of adjacent stations over which slowness is estimated, x j − x 0 is the distance between station j and the middle station, g(t) is the seismic trace, and h(t) is the Hilbert transform of g(t). The slowness value with the highest semblance represents that of the most locally coherent plane wave at the specific time t. Including the Hilbert transform in Eq.
(2) is equivalent to applying the conventional definition of semblance (Taner et al., 1979) on the analytical signal g(t) + ih(t). This approach has the key advantage of allowing for reliable semblance calculation at the zero-crossings of the original signal, owing to the property that the amplitude of the analytical signal (which is the signal envelope) does not have zero-crossings. For optimal slowness determination, fiber segment lengths should correspond to the longest wavelength of interest. However, when implementing the local slant-stack transform, segment lengths should be chosen to be upper bounded by the wavelength and the segment for which the signal remains coherent (coherency length) and lower bounded by the desired spatial and slowness resolutions (e.g., Ventosa et al., 2012). Spatial and slowness resolutions exhibit a tradeoff since increasing the segment length will increase the slowness resolution and decrease the spatial resolution. Coherency lengths are expected to be small for DAS strain (rate) records since they may be dominated by scattered waves (e.g., Lior et al., 2021) and are extremely sensitive to local media heterogeneities (e.g., van den Ende and Ampuero, Figure 1. Flow chart detailing the conversion procedure from DAS strain (rate) signals to ground velocity (acceleration). 2020; Singh et al., 2019) and fiber coupling (e.g., Sladen et al., 2019;Lior et al., 2021). Thus, in the current application, short cable segments are used, as further described.
The process of converting DAS strain (rate) signals to ground velocity (acceleration) is detailed here, summarized in Fig. 1, and demonstrated using simulated data in the next section. Since seismic signals are a combination of various sources (e.g., earthquake waves, ambient noise, random signals) and may include dispersive waves, the signals need to be filtered at the frequency band of interest. Filtering will reduce noise and limit the dispersive effect; however, simultaneous wave arrivals, complex propagation and dispersion effects are still expected. The slant-stack transform is applied per virtual seismometer along the fiber, using the 2L + 1 adjacent traces (L on each side). The range of examined slowness values is chosen to be between −0.01 and 0.01 s m −1 with 0.0002 s m −1 slowness intervals (i.e., 100 slowness values), and at each time instance, the slowness value is determined based on the maximum semblance. The slowness time series (derived from semblance) may often be characterized by abrupt variations of value and sign (i.e., propagation direction) owing to complex wave propagation, interference and dispersion effects. Thus, a moving average, with window size set to be equal to the signal's longest period of interest, is applied to the absolute value of the slowness (preventing the averaging of slowness values with different signs). The sign is then determined as the one that dominates each averaged window, i.e., the most recurrent sign. The filtered strain-rate signals are then converted to ground accelerations using Eq. (1). Slowness sign variations result in discontinuities in the converted strain-rate records, which are smoothed by additional filtering. The second filter may be chosen to be identical to the first, yet since it is applied to smooth discontinuities, different filters containing an upper frequency limit (i.e., low-pass or band-pass), may be applied. In the following sections the same four-pole zero-phase Butterworth filter is used for both filtering operations.
Throughout this paper, the semblance-based slowness determination method is compared to an FK-based method, applied as follows. Frequency-wavenumber transforms are calculated on the filtered strain (rate) signals in consecutive windows using the same number of adjacent traces used for semblance calculation. Combinations of high-amplitude temporal frequency (f ) and spatial frequency (υ) are identified as those whose spectral amplitudes are higher than the 99th percentile of all amplitudes in the FK domain. The spectral amplitudes are summed separately for the two FK quadrant (positive f and positive υ or positive f and negative υ), and slowness is estimated using data from the higher sum quadrant by fitting f = υ/p x . The slowness time series is then smoothed and used for strain (rate) to ground motion conversion in the same manner as previously described for the semblance analysis.
3 Validation using simulated earthquake records

Simulated data
To validate the proposed method, simulated earthquake waveforms are produced for a simple 2D velocity model representing an underwater sedimentary basin. A basin was simulated to generate the commonly observed interactions of waves propagating in opposite directions. A spectral-element-based numerical simulation is done using the SPECFEM2D 7.0.0 code published under the CECILL V2 License (Komatitsch et al., 2012). The simple, though adequate, physical model is composed of two linear elastic sub-domains: a trapezoidal sedimentary layer with maximum depth of 100 m, characterized by P-and S-wave velocities of 1600 and 400 m s −1 , respectively, and density of 2000 kg m −3 (yellow in Fig. 2a), overlying a bedrock with P-and S-wave velocities of 2500 and 1200 m s −1 , respectively, and density of 2200 kg m −3 (orange in Fig. 2a). These layers lie beneath a thin water layer at 20 m depth with P-and S-wave velocities of 1500 and 0 m s −1 , respectively, and density of 1000 kg m −3 (green in Fig. 2a). The unstructured numerical mesh is generated using the Gmsh software (Geuzaine and Remacle, 2009) distributed under the terms of the GNU General Public License (GPL). It contains 2344 quadrilateral elements and, using an interpolation polynomial order of 5, allows simulations up to 12 Hz maximum frequency. A double-couple point source is located at I. Lior et al.: Strain to Ground Motion Conversion X = 2500 m, Z = 500 m and is time-modulated by a Ricker wavelet with a central frequency of 4 Hz, corresponding to a M w = 1 earthquake. Receivers are regularly spaced at the bottom of the water layer from X = 400 m to X = 2600 m. Spatial and temporal sampling were set to be 5 m and 5 ms, respectively. These velocity waveforms represent 12.5 s of seismic recordings at 441 sensors, as seen in Fig. 2b. Since DAS can only record deformations along the optical fiber, simulated seismograms have a single component, oriented along the array axis. These waves exhibit complex propagation patterns, with visible P and S waves, as well as surface waves reflected from the edges of the basin, providing an excellent test case for the proposed algorithm.
Simulated velocity waveforms were differentiated in time and space to obtain ground accelerations and strain rates, respectively; noise was added to the latter, and the ability to convert the latter to the former via the proposed slant-stack approach is examined. Strain-rate records were calculated aṡ (x) = (V (x+GL/2)−V (x−GL/2))/GL and ground acceleration records were calculated as A(t) = (V (t +dt)−V (t − dt))/(2dt), where V , GL and dt are the simulated velocity, gauge length and temporal sampling, respectively. Simulated strain-rate signals are thus characterized by gauge length and spatial sampling of 10 and 5 m, respectively. To reliably generate ocean-bottom DAS records, ambient noise measurements were added to simulated strain rates, keeping the noise records' spatial correlation (Fig. 2c). In the simulated water depth (i.e., 20 m), ocean-bottom DAS records are typically dominated by surface gravity waves (e.g., Lior et al., 2021), which may be easily filtered due to their lower frequency content compared to the simulated earthquake. Thus, ambient noise recorded at a water depth of ∼ 800 m was used. These records are composed of instrumental noise and secondary microseisms at similar frequencies to those of the simulated earthquake (Fig. 2d). The added noise measurements were recorded on 22 July 2019 by an underwater cable deployed offshore from Toulon, southern France (Sect. 4). Noise records, sampled at 100 Hz and 10 m (spatial sampling is equal to the gauge length), were resampled to match the simulated data using a 2D interpolation function. Noise records were then differentiated to strain rate and rescaled to simulate challenging noise conditions, with an average signal-to-noise ratio (SNR) of 8.2 (Fig. 2d). Here, SNR was calculated as the root-mean-square (RMS) ratio of average signal and noise amplitude spectra between 2 and 12 Hz. In spite of the added noise, accelerations converted from strain rates are compared to simulated accelerations (derived from simulated velocities by a finite-difference time derivative). Noise was not added to the latter, constituting a stringent algorithm validation.

Strain rate to ground acceleration conversion
For the simulated data, waveforms were low-pass-filtered at 12 Hz and short array segments of 100 m (L = 10) were used to calculate the semblance and to convert each strainrate signal into an acceleration seismogram. These segments are 25 % shorter than the longest apparent wavelength observed for P waves and provide a sufficient compromise between slowness and spatial resolution, allowing for reliable strain (rate) to ground motion conversion, as further shown. FK-based analysis was conducted on consecutive 1 s windows. An example of the slant-stack conversion, applied to trace no. 220 of the simulated waveforms ( Fig. 2b), is shown in Fig. 3. Panel a shows 7 s of strain-rate signals composed of P waves (0.6-1.2 s), S waves (1.2-2.5 s) and surface waves (2.5-7 s), including reflected waves with opposite propagation direction (4.5-7 s). The semblance analysis in panel b shows apparent phase slowness, corresponding to the maximum semblance and smoothed slowness in orange circles and the red curve, respectively. Red dots correspond to the ratio between strain rates and accelerations (both differentiated from the simulated velocity waveforms) plotted at acceleration maxima and minima, i.e., the optimal slowness values for strain-rate conversion. The similarity between the smoothed slowness (red curve) and optimal slowness for conversion (red dots) suggests that the resolved slowness may be reliably used for strain-rate to acceleration conversion. As expected, semblance and FK analysis are able to resolve different velocities for different phases. FKderived velocities are significantly slower than semblancederived velocities: semblance analysis resolves average apparent phase slowness of 0.77, 0.83 and 1.67 s km −1 (velocities of 1.3, 1.2 and 0.6 km s −1 , respectively) for P waves, S waves and surface waves, respectively, while FK analysis resolves average apparent phase slowness of 1.67, 1.67 and 2.5 s km −1 (velocities of 0.6, 0.6 and 0.4 km s −1 ) for P waves, S waves and surface waves, respectively (in the previously indicated intervals). When no noise is added, the same velocities are determined, with the exception of a lower P-wave apparent slowness of 0.43 s km −1 (velocity of 2.3 km s −1 ) for the semblance-based method. Figure 3c compares the acceleration signal (blue curve, differentiated from simulated velocity) with strain-rate-converted accelerations using semblance-derived slowness (red curve) and FK-derived slowness (black curve).
Excellent agreement is observed between acceleration (blue curve in Fig. 3c) and converted strain rate, when the latter is obtained using semblance-based slowness (red curve in Fig. 3c). This agreement persists for both fast (P waves and S waves) and slow (surface waves) phases, demonstrating the ability to reliably resolve phase velocities even for short-duration waves and wavelengths longer than the fiber segments used for slant-stack analysis, e.g., P waves. The low apparent velocities derived from FK analysis result in lower converted strain-rate signals (black curve in Fig. 3c) Figure 2. Velocity model and simulated velocity, strain rate and noise data. Panel (a) shows the unstructured 2D numerical mesh: the double-couple source location is indicated by a black star and sensors are indicated by black triangles, equally spaced between X = 400 m to X = 2600 m (first and last sensor channel numbers are indicated). Panel (b) shows simulated velocity waveforms, with visible P waves, S waves and Scholte waves. Panel (c) shows strain-rate signals added with ambient noise. Panel (d) shows simulated strain rate and noise signals. Noise curves for channel no. 220 and averaged for all channels are indicated by semi-transparent and solid blue curves, respectively. Signals added with noise for channel no. 220 and averaged for all channels are indicated by semi-transparent and solid orange curves, respectively. SNR values are reported in the legend. compared with acceleration records. Since reliable FK slowness estimation requires sufficiently long data segments, this method does not provide sufficient slowness resolution to resolve high apparent velocities in short spatial and temporal intervals. When different waves interfere, the single plane wave assumption in Eq. (1) does not hold, the phases of ve-locity and strain-rate signals may be misaligned (4.5-5.1 s), and phase velocities may not be reliably resolved, resulting in lower quality conversion. When the intensity of such effects fade, excellent agreement is again observed (5.1-7 s).
The use of time-dependent slowness is found to be particularly advantageous when velocities abruptly vary and when Red dots correspond to the ratio between strain rates and accelerations, both differentiated from the simulated velocity waveforms. (c) Accelerations that are time-differentiated from simulated velocities (blue curve) compared with accelerations that are space-differentiated from simulated velocities and converted using semblancederived slowness (red curve) and FK-derived slowness (black curve). Standard deviations of the residuals between converted strain rates (red and black curves) and time-differentiated velocities (blue curve) are indicated in the panel legend. The signals in the gray regions have been amplified by a factor of 6 for easy comparison of low-amplitude signals.
propagation direction changes. Specifically, the semblance approach is substantially more robust compared with the FK scheme when implemented on short spatial and temporal intervals. This is shown in Fig. 3c by comparing the FK-based converted strain rate (black curve) and semblance-based converted strain rate (red curve) and is further illustrated in Fig. 4, where strain-rate-converted accelerations, using both semblance and FK slowness, are compared with differentiated velocity for all simulated data. Comparing the residuals between strain-rate-converted accelerations (panels c and d) and differentiated velocity (panel a) using semblancederived slowness (panel e) and FK-derived slowness (panel f) demonstrates the benefit of using the former. The use of FK slowness produces larger residuals (panel f) than the use of semblance slowness (panel e), especially for the direct P and S waves (0-3 s). Standard deviations are thus systematically higher when using FK slowness, as seen in panel b.
To examine the performance of the semblance-based strain-rate conversion for dispersive waves (e.g., Scholte waves), this technique is applied to strain-rate records filtered at different frequency bands. For this test, noise was not added to the simulated data. An example for the same trace shown in Fig. 3 (no. 220) is shown in Fig. 5, where 6 different frequency bands are tested: filtered strain-rate data are shown in panel a, and a comparison between converted strain rate (red curves) and accelerations (blue curves) is shown in panel b. Different phase velocities are measured for dispersive waves (2.5-7 s) at different frequency bands, while the same velocity is measured for body waves (1-2.5 s) (panel c). The agreement between accelerations and converted strain rate in narrow frequency bands is only slightly better than that obtained for broadband signals, low-passfiltered at 12 Hz (top line in panels a and b), as indicated by the standard deviations of the residuals (top right corners of panels b). The small difference between bandpass filtered (e.g., 4.8-7.6 Hz in Fig. 5) and broadband (0-12 Hz in Fig. 5) standard deviations of the residuals suggests that the dispersive nature of these waves has a small effect on the conversion quality. Thus, applying time-specific and frequencyspecific slowness to convert broadband seismic signals is unlikely to result in a significant improvement in conversion robustness. Also, such an approach may require intricate processing, introduce artifacts, and is highly subjective and challenging to implement. Thus, this analysis is focused on the more generic case of a broadband signal, and a bandpasslimited conversion is not developed and implemented.
The conducted analysis and excellent agreement between simulated strain-rate-converted seismograms and acceleration signals demonstrate the advantages of using the proposed slant-stack-based approach for strain rate to ground motion conversion. Next, the technique is implemented on

Data
For a real data application of the conversion technique, eight local earthquakes are used. These were recorded by three dark ocean-bottom telecommunication cables deployed in the Mediterranean Sea. This dataset is identical to the one used by Lior et al. (2021), and the cables' locations, bathy-metric profiles and layout are detailed there. This information, as well as earthquake data, is briefly repeated here. Data were acquired by Géoazur on two cables deployed offshore from Methoni, southwestern Greece, and one cable deployed offshore from Toulon, southern France. In addition to DAS data, several on-land seismometers, installed near the cables, were available during the measurements. These will be used later (Sect. 6) to compare DAS, and seismometerderived source parameters. The earthquake data are detailed in Table 1; earthquakes, cables and seismometer locations are shown in Appendix A; and magnitude and hypocentral  (Table 1). Hypocentral distances are measured to specific cable segments or broadband (BB) sensors. Solid and empty symbols correspond to seismometer and DAS data, respectively. Data for NESTOR, HCMR and MEUST are indicated by circles, triangles and diamonds, respectively, and the number of data points is indicated in parentheses in the legend. distance distribution are shown in Fig. 6. In the latter, variations of hypocentral distance correspond to the different analyzed fiber segments. The earthquake DAS records used in this study are dominated by scattered Scholte waves, with velocities between 100 and 400 m s −1 (Lior et al., 2021).
The two cables deployed offshore from Greece are used for the HCMR (Hellenic Centre for Marine Research) and NESTOR (Neutrino Extended Submarine Telescope with Oceanographic Research) projects. The acquisitions were conducted on 18 and 19 and 19-25 April 2019 on the HCMR and NESTOR cables, respectively. A Febus A1 DAS interrogator was used, measuring strain-rate signals. Data were sampled at 6 and 5 ms for HCMR and NESTOR, respectively, and gauge length and spatial sampling were both set to 19.2 m for both cables. These records amount to 688 and 1365 equally spaced channels for the 13.2 and 26.2 km long HCMR and NESTOR cables, respectively. In addition, two seismometers installed near the on-land end of the fibers were available during part of these campaigns, METN and METS.
The cable deployed offshore from Toulon is used for the MEUST-NUMerEnv project (Mediterranean Eurocentre for Underwater Sciences and Technologies -Neutrino Mer Environnement) (Lamare, 2016) and was previously harnessed for DAS measurements by Sladen et al. (2019). The acquisition was conducted on 11-31 July 2019 using a chirped-pulse hDAS interrogator developed by Aragon Photonics, measuring strain signals (Pastor-Graells et al., 2016;Fernandez-Ruiz et al., 2019;Williams et al., 2019). Data were sampled at 10 and 2 ms for the first and last 10 d of the campaign, respectively, and gauge length and spatial sampling were both set to 10 m. These records amount to 4480 equally spaced channels for the 44.8 km long cable. In addition, two seismometers were installed near the on-land end of the fiber, POSAN and POSAS.

Application to DAS recorded earthquakes
To demonstrate the performance of the proposed conversion approach, DAS strain (rate) earthquake signals are converted to ground velocities (accelerations). For each of the three cables, short fiber segments that exhibit coherent and continuous waveform recordings are chosen. Each section contains 29 traces that are filtered between 1 and 5 Hz. For different applications, a different filter may be used, as demonstrated in the next section. Filtered signals are converted to ground velocities or accelerations using fiber segments of ∼ 380 m (L = 10 for NESTOR and HCMR and L = 19 for MEUST). Compared with simulated data, longer segments are used in order to resolve faster seismic phases and longer wavelengths. For the FK transform, 2 s data intervals were used. Applying the conversion to DAS-recorded earthquakes highlights body-wave arrivals, since these fast waves exhibit higher converted amplitudes compared with later-arriving scattered waves and presignal ambient noise. Figure 7 shows an example of strain conversion to ground velocities for an M L 2.6 earthquake recorded by the MEUST cable at a hypocentral distance of 185 km, and Fig. 8 shows an example of strain-rate conversion to ground accelerations for an M L 3.6 earthquake recorded by the NESTOR cable at a hypocentral distance of 140 km. In Fig. 7, mostly direct S waves are shown, exhibiting unilateral wave propagation (panels a-c). The apparent velocity of the direct S waves (1.2-2.2 s) is determined to be 2 and 1 km s −1 using semblance and FK-derived slowness, respectively, while later-arriving waves travel at ∼ 400 and ∼ 370 m s −1 using semblance-and FK-derived slowness, respectively (panel c). Thus, as observed for the simulated data (Sect. 3), direct waves (panels b and d) exhibit higher converted velocity amplitudes compared to later phases. This is visualized by comparing the color codes in panels a and b. Figure 8 shows both P and S waves, as well as presignal noise. In this example, scattered waves, and thus bilateral wave propagation, dominate the measurements (Lior et al., 2021), and several slowness sign flips are evident. The apparent velocity of the first arriving S waves (24-25 s) is 1.3 km s −1 and 800 m s −1 using semblance-and FK-derived slowness, respectively, while the average apparent velocity is 750 and 530 m s −1 , using semblance-and FK-derived slowness, respectively, resulting in higher converted acceleration amplitudes for the direct waves. Ocean-bottom presignal noise is dominated by instrument-related effects (e.g., Lior et al., 2021; and ambient noise. These signals are character- ized by low apparent velocities, which results in low acceleration amplitudes and facilitates easy identification of the initial P waves, subject to SNR conditions. Next, the ability to invert for the source parameters, i.e., seismic moment, corner frequency and stress drop, is examined by converting strain (rate) records for predefined P-and S-wave intervals and fitting their spectra with an earthquake source model.

Implications for source parameter inversion
Seismic moment, source corner frequency and stress drop are determined by fitting converted earthquake DAS signals with an earthquake source model. The source model chosen is the commonly used omega-squared model (e.g., Brune, 1970;Madariaga 1976;Sato and Hirasawa, 1973), which describes the far-field body wave radiation for ground displacements, velocities and accelerations. The model is fit to the data via the single-step inversion of Lior and Ziv (2018). This approach is advantageous as model fitting is done in the time domain, circumventing the time-to frequency-domain transformation and avoiding many spectral model fitting intricacies.

Source model
For ground displacements, velocities and accelerations, the omega-squared model subject to high-frequency attenuation (Anderson and Hough, 1984) is given by the following: and respectively, where 0 is the low-frequency displacement spectrum plateau, f is frequency, f 0 is the source corner frequency and κ is the attenuation parameter. The latter can be expressed as an attenuation corner frequency as f κ = 1/(π κ) (Eq. 4 of Lior and Ziv, 2018). The spectral parameters 0 and f 0 correspond to the seismic moment, M 0 , and stress drop, τ , as follows: where ρ is the density at the source, C is the wave velocity at the source (C P and C S for P and S waves, respectively), R is the hypocentral distance, U ϕθ is the radiation pattern, F s is the free-surface effect, and k is a constant which depends on the wave type and rupture speed (Madariaga, 1976). Equation (4b) applies to a circular crack (Eshelby, 1957) expanding isotropically at constant rupture speed. Parameter tuning is set as follows: ρ = 2600 kg m −3 ; C P = 5333 m s −1 ; C S = 3200 m s −1 (e.g., Lior and Ziv, 2020); U ϕθ equals 0.52 and 0.63 for P and S waves, respectively (Aki and Richards, 2002); F s = 2; and k equals 0.32 and 0.21 for P and S waves, respectively, corresponding to a rupture speed of 0.9C S (Madariaga, 1976). Using F s = 1.7 instead of F s = 2, as is sometimes used for ocean-bottom applications (e.g., Webb, 1998), will reduce magnitude estimates by ∼ 0.047, a minute difference compared to magnitude uncertainties, as shown in the next subsection.

From strain (rate) to source parameters
The spectral parameters are determined via the single-step inversion of Lior and Ziv (2018), which resolves 0 , f 0 and κ using the time-domain signals, circumventing the time-to frequency-domain transformation, required by most source parameter inversion methods. The approach is fully detailed in Lior and Ziv (2018) and briefly summarized in Appendix B. DAS strain (rate) data are converted to ground velocity (acceleration) for manually chosen P-and S-wave windows, and source parameters are resolved. The procedure of determining the frequency band of interest, strain (rate) to ground  Fig. 9 for the S waves of an M L 3.6 earthquake, recorded between 20.95 and 21.5 km along the NESTOR cable, at a hypocentral distance of 145 km. First, for the 29 traces composing each cable segment, strain (rate) signal and noise amplitude spectra (AS) are calculated, resampled (following the procedure described in McNamara and Buland, 2004) and stacked for signal and presignal time-windows of equal length (dashed black and solid gray curves in panel a, respectively). The analyzed presignal noise is that recorded 2 min before the signals. The bandwidth for which frequency-specific SNR (signal AS(f )/noise AS(f )) is larger than 2 is used for subsequent analysis (solid black curve in panel a). If less than three discrete frequencies have high SNR, the recording is disregarded. To fully preserve the frequency band of interest, the filter's lower and upper corner frequencies are slightly decreased and increased by factors of 10 −0.2 and 10 0.2 , respectively. Strain (rate) signals are then filtered (panel b), converted to ground velocities (accelerations) (panel c), and differentiated and/or integrated to obtain ground displacements, velocities and accelerations. Following each differentiation/integration, the aforementioned filter is applied. The signals' RMSs are calculated in the time domain, and source parameters are determined as detailed in Appendix B. An example of the single-step inversion's results is shown in Fig. 9d where the best-fitting model is plotted using Eq. (3c) and compared with observed stacked acceleration spectra, and Fig. 9e shows the best-fitting parameter combination ( 0 , f 0 and κ), indicated by a red star, in log(f 0 ) − log(f κ ) space. The color code in panel e corresponds to the bestfitting 0 (for each log(f 0 ) − log(f κ ) combination), and the contours correspond to the objective function's value. The solution exhibits a high degree of trade-off between the values of f 0 and κ (panel e), yet good agreement is found between observed (black curve) and modeled (blue curve) spectra (panel d).
To further compare the performance of FK-and semblance-based slowness, magnitude estimates are com- pared. Figure 10 plots the difference between magnitudes estimated following the FK-and semblance-based conversion schemes as a function of catalog magnitude. Magnitudes estimated for data converted using FK-derived slowness are generally lower than those resolved using semblance-derived slowness. This trend is expected given the lower apparent velocities determined via the FK-based approach, as previously shown for both simulated (Fig. 3) and observed (Figs. 7 and 8) earthquakes.
Moment magnitudes, stress drops and corner frequencies resolved using DAS semblance-based converted data are found to be in good agreement with those estimated on adjacent on-land broadband seismometers. Since the DAS fiber and the seismometers are not colocated, their estimated source parameters are not expected to be in perfect agreement. DAS P-and S-wave magnitude estimates are plotted as a function of average S-wave seismometer magnitudes in Fig. 11, and DAS S wave f 0 and τ estimates are plotted as a function of average S-wave seismometer-obtained parameters in Fig. 12. The low-SNR conditions observed for P waves, i.e., the narrow available frequency-band, did not allow for robust estimates of f 0 and τ , which are thus not shown. A similar comparison with catalog magnitude is shown in Fig. 13, noting that local and moment magnitudes may differ for small earthquakes (Deichmann, 2006). DAS parameter errors were calculated as the standard deviations of parameters determined for each individual seismogram in the analyzed fiber segment, while seismometer errors are the standard deviations of single seismometer estimates, when available. Data for the HCMR cable are not shown since seismometer gain was unavailable.

Discussion
The comparison between semblance-derived slowness and FK-derived slowness for strain (rate) to ground motion conversion reveals several advantages favoring the semblancebased approach. Unlike semblance analysis, whose implementation and interpretation is simple and objective, FK analysis introduces considerable subjectivity into the slowness determination procedure since interpreting an FK image may be done in various ways. The FK analysis generally yielded lower velocities owing to its poor slowness resolution, as observed for both simulated and recorded earthquakes. For these reasons, when comparing simulated accelerations and converted strain rates, a significantly better agreement was obtained using the semblance-based approach. The ground motion conversion difference between FK-and semblance-based approaches is most pronounced for direct arrivals, i.e., P and S waves. Thus, the semblancebased conversion is particularly advantageous for EEW, when short-duration, relatively fast propagating waves, are required for speedy and reliable source parameter estimation. A significant hindrance when using DAS is the mechanical coupling between the fiber and the solid Earth. This issue is particularly troublesome when standard telecommunication fibers are used, specifically those deployed underwater, as their coupling quality is often unknown and may prevent reliable seismic monitoring (e.g., Sladen et al., 2019;Lior et al., 2021). Here, specific segments for which recording quality is sufficiently uniform were manually identified. The effect of coupling along these limited segments is quantified by considering the signals' average absolute strain-rate amplitudes in decibels (dB), plotted in Fig. 14 for all traces, segments and cables. In Fig. 14, only earthquakes recorded at hypocentral distances longer than 80 km are plotted, to ensure slow propagation-related amplitude changes along the fiber. The variabilities are generally small, limited to 2-3 dB with several exceptions of ∼ 8 dB. These mostly minor deviations, along with the small DAS magnitude uncertainties (vertical error bars in Fig. 11 indicate standard deviations of magnitudes resolved independently for each DAS channel), indicate that these segments display sufficiently uniform coupling for ground motion conversion and source parameter estimation. Moreover, the fact that even for non-uniform coupled cables such as those used here, sufficiently uniform coupling is observed, even if in limited segments, demonstrates the potential of underwater fibers for reliable source parameter estimation.
The proposed slant-stack conversion approach relies on the ability to resolve the phase velocity of a single plane wave at every time instance. However, seismic records are often dominated by several waves, which may be dispersive (e.g., surface waves); characterized by different velocities and incidence angles; and exhibit complex propagation, scattering and interference patterns. The analysis on simulated data in Sect. 3 demonstrates that when a single plane wave is considered, converted strain rates are in excellent agreement with acceleration waveforms, while when two opposing plane waves interfere, the conversion's robustness is decreased (e.g., Fig. 3). Careful filtering of DAS signals needs to be applied to isolate specific plane waves from DAS earthquake records, as demonstrated in part in Fig. 5, yet conversion errors may still result from inadequate slowness resolution, incoherent plane waves and noise. The effect of noise is seen in Sect. 3 where synthetic P waves suffer from low SNR, which results in lower-than-expected apparent velocities and a slightly reduced conversion quality. The amplitudes of DAS P waves, which arrive at near-vertical incidence angles with respect to horizontal fibers, are especially low since they are reduced by a factor of squared cosine of incidence angle (e.g., Mateeva et al., 2014). Further consideration of these issues is beyond the scope of this paper. In spite of these complexities, converted DAS signals allow for reliable magnitude estimation, demonstrating the robustness of the conversion procedure.
Resolving source corner frequencies, and thus stress drops, for small and/or distant earthquakes in unfavorable SNR conditions, is a challenging task since high-frequency source effects, i.e., f 0 , are masked by high-frequency attenuation, i.e., f κ (e.g., Lior and Ziv, 2018). To address this issue, Lior and Ziv (2018) introduced an uncertainty parameter, δ, a quality control measure for the ambiguity between bestfitting f 0 and κ. The value of δ is proportional to the area enclosed within the 5 % contour in log(f 0 ) − log(f κ ) space (e.g., Fig. 9e). Lior and Ziv (2018) found that solutions with high δ values (typically > 6 %) usually exhibit a high degree of ambiguity between f 0 and κ. For such cases, they implemented a two-step inversion approach, which consists of determining a station, or in this case fiber segment, specific κ using low-δ solutions and repeating the inversion for 0 and f 0 . However, since for the data used in this study only a few low-δ solutions were obtained, this technique is not implemented here. Thus, corner frequencies and stress drops are not well constrained. The standard deviations to parameter residuals (Fig. 12) are within-event variabilities between the estimates of specific DAS segments and seismometers. That these values are only slightly higher than within-event variabilities reported by Lior and Ziv (2018) suggests that in spite of the inability to reliably determine these parameters, DAS and seismometer-derived f 0 and τ are found to be in good agreement (Fig. 12).
Even in cases where f 0 and κ (and thus τ ) may not be well constrained, 0 and thus moment magnitudes are reliably determined (Lior and Ziv, 2018). This is visualized in Fig. 9e, where 0 values (color code) are generally subparallel to the OF contours.
Magnitudes are reliably determined for both P and S waves (Fig. 11) in spite of the reduced sensitivity of horizontal fibers to transverse deformations, as expected for P waves. Recorded deformation amplitudes are modulated by cos 2 θ, where θ is the wave's incidence angle with respect to the fiber's axis (e.g., Ajo-Franklin et al., 2019;Kuvshinov, 2016;Mateeva et al., 2014;Papp et al., 2017;Yu et al., 2019); thus, DAS measurements are mostly sensitive to deformations along the fibers' axis, i.e., elongation and compression. Since direct P waves are expected to arrive at near-vertical incidence angles, they would not induce significant deformations: while direct S-wave arrivals are clearly identified for several fiber segments (e.g., Fig. 7), direct P waves are not. However, analyzed DAS records are dominated by low- Figure 11. Comparison between DAS and seismometer magnitude estimates. Moment magnitudes estimated using DAS recorded S waves (a) and P waves (b) are plotted as a function of average magnitudes obtained using seismometer records. DAS event averaged magnitudes are plotted in red squares. The black curve is a 1 : 1 line. The number of data points and standard deviations to magnitude residuals are indicated in the legends. DAS magnitude errors are the standard deviations for magnitudes determined for each individual seismogram in the analyzed fiber segment, while seismometer errors are the standard deviations for single seismometer estimates. Figure 12. Comparison between DAS and seismometer source corner frequencies and stress drops. log(f 0 ) (a) and log( τ ) (b) estimated using DAS recorded S waves are plotted as functions of average parameters obtained using seismometer records. DAS event averaged parameters are plotted in red squares. The black curve is a 1 : 1 line. The number of data points and standard deviations to parameter residuals are indicated in the legends. DAS parameter errors are the standard deviations for parameters determined for each individual seismogram in the analyzed fiber segment, while seismometer errors are the standard deviations for single seismometer estimates. velocity waves, following both direct P-and S-wave arrivals (e.g., Fig. 8). These are scattered P and S waves, propagating in a variety of horizontal directions (Lior et al., 2021), which are easily measured on horizontal fibers and used here to infer source parameters. That earthquake magnitudes, determined by DAS measurements of scattered waves, are in close agreement with both catalog magnitudes (Fig. 13) and seismometer-derived magnitudes (Fig. 11) indicates that these waves reliably represent the source characteristics and may be used for source parameter inversion.
The ability to infer source parameters using P waves recorded on horizontal fibers is key for harnessing DAS, specifically using underwater fibers, for EEW. To this end, the proposed algorithm will need to be adapted for real-time  Signal amplitude variability for all analyzed fiber segments for the three cables. Number of channels along the limited segments (29 traces) as a function of demeaned 10log 10 (X) with X being the average absolute strain-rate amplitudes. Data for NESTOR, HCMR and MEUST are plotted in red, green and blue, respectively.
performance. The goal of EEW systems is to robustly and rapidly predict ground shaking intensities, an objective that is typically achieved by estimating earthquake source properties in real time (e.g., Allen and Melgar, 2019; Lior and Ziv, 2020). To this end, and in order to issue ground shaking alerts as early as possible, seismic observations should be obtained at close proximity to earthquake epicenters, and source parameters should be estimated using both P and S waves. Since many of the most hazardous earthquakes on Earth occur at subduction zones, and therefore underwater, the ability to determine source parameters using both P and S waves recorded by ocean-bottom DAS will significantly improve the performance of EEW systems for underwater earthquakes and enhance hazard mitigation capabilities.

Conclusions
In this study, the ability to convert DAS strain (rate) signals to ground motion records and resolve earthquake source parameters is demonstrated. An algorithm for DAS data to ground motion conversion is presented: apparent phase slowness is determined at every time instance using a semblancebased local slant-stack transform and used to convert strain (rate) to ground velocities (accelerations). The algorithm is successful at resolving the apparent velocities of different seismic phases. Validation using simulated waveforms reveals excellent agreement between simulated accelerations and converted strain-rate signals even in the presence of correlated noise and propagation direction variations. Application of the algorithm to eight earthquakes recorded by ocean bottom DAS fibers in the Mediterranean Sea highlights fast waves (body waves) since they exhibit high converted ground motion amplitudes compared with low-velocity scattered waves and presignal ambient noise. Earthquake magnitudes and stress drops were determined for P and S waves using the single-step approach of Lior and Ziv (2018), circumventing the time-to frequency-domain transformation typically required for moment and corner frequency estimation. Close agreement is observed between source parameters determined using on-land broadband seismometers and ocean-bottom DAS, even when source corner frequencies and stress drops are not well constrained due to significant high-frequency attenuation. This ability to resolve earthquake magnitudes using P waves recorded by horizontal ocean-bottom fibers is key for implementing DAS for EEW. The algorithm for strain (rate) conversion may be adapted for real-time applications and used in conjunction with realtime source parameter determination schemes (e.g., Lior and Ziv, 2020) for a DAS-based EEW system. Harnessing DAS for EEW, specifically using ocean-bottom fibers, will significantly improve hazard mitigation capabilities for underwater earthquakes and tsunami earthquakes.   Code and data availability. Simulated and observed DAS earthquakes are available on https://osf.io/98cnk/ (last access: 15 January 2021) and https://osf.io/4bjph/ (last access: 27 December 2020), respectively. Broadband seismometer data were acquired by Géoazur; data for the POSAN and POSAS stations were downloaded from RESIF (http://seismology.resif.fr/, (last access: 1 May 2020).
Author contributions. IL designed the presented algorithms, performed the analysis and wrote the initial draft. DM designed and performed the simulations. SS helped to adapt and implement the slant-stack approach to DAS data. AS, JPA, DM, DR and SS contributed to the discussion, methodology, interpretation and presentation of the results.
Competing interests. The authors declare that they have no conflict of interest.
Special issue statement. This article is part of the special issue "Fibre-optic sensing in Earth sciences". It is not associated with a conference.