Sensitive seismic sensors based on microwave frequency fiber interferometry in commercially deployed cables

The use of fiber infrastructures for environmental sensing is attracting global interest, as optical fibers emerge as low cost and easily accessible platforms exhibiting a large terrestrial deployment. Moreover, optical fiber networks offer the unique advantage of providing observations of submarine areas, where the sparse existence of permanent seismic instrumentation due to cost and difficulties in deployment limits the availability of high-resolution subsea information on natural hazards in both time and space. The use of optical techniques that leverage pre-existing fiber infrastructure can efficiently provide higher resolution coverage and pave the way for the identification of the detailed structure of the Earth especially on seismogenic submarine faults. The prevailing optical technique for use in earthquake detection and structural analysis is distributed acoustic sensing (DAS) which offers high spatial resolution and sensitivity, however is limited in range (< 100 km). In this work, we present a novel technique which relies on the dissemination of a stable microwave frequency along optical fibers in a closed loop configuration, thereby forming an interferometer that is sensitive to deformation. We call the proposed technique Microwave Frequency Fiber Interferometer (MFFI) and demonstrate its sensitivity to deformation induced by moderate-to-large earthquakes from either local or regional epicenters. MFFI signals are compared to signals recorded by accelerometers of the National Observatory of Athens, Institute of Geodynamics National Seismic Network and by a commercially available DAS interrogator operating in parallel at the same location. Remarkable agreement in dynamical behavior and strain rate estimation is achieved and demonstrated. Thus, MFFI emerges as a novel technique in the field of fiber seismometers offering critical advantages with respect to implementation cost, maximum range and simplicity.

shows the MFFI in a conceptual form. One can easily see that this is about a Mach-Zehnder delay interferometer where the two branches have a path difference of Lrtr where Lrtr is the total length of the fiber that is interrogated and it is a round-trip, closed-loop path starting and ending at the same point.  The only difference with a typical Mach-Zehnder is that the phase detection is accomplished after photodetection in the RF domain as depicted above and will explained mathematically in the next paragraphs.
The modulation signal VM applied on a Mach Zehnder modulator (MZM) is a sinusoidal voltage with angular frequency ωRF and amplitude Vm superimposed on a DC voltage VB: VM=VB+Vmsin(ωRFt) The output field as a function of the input voltage is given by 1 : Where A is the optical insertion loss, B=cos(πVB/Vπ), C=πVm/Vπ, D=sin(πVB/Vπ) and Vπ is the voltage required to change the phase of the modulated MZM arm by π, Jn(x) is the n-th order Bessel function of the first kind.
Depending on C, D values, one can either promote odd or even harmonics and preserve or eliminate the DC component. The latter -if removed -will make the system immune to dispersion as it will be explained later on.
If C is properly selected, then without loss of generality, we may deduce that the output field is mostly carrying the first harmonic ωRF of the sinusoidal signal.
The term m is the modulation depth, φnoise is the phase noise of the microwave source.
After transmission, and ignoring the chromatic dispersion effects at this point, the signal will acquire a phase shift in the microwave domain that is equal to φRF=ωRFTprop and a phase shift in the optical domain that is equal to φopt=ωoptTprop, where Tprop=Lrtr/vg Tprop is the time required for a round-trip propagation (from the transmitter back to the receiver) which depends on the round trip distance Lrtr and the group velocity vg=c/ng. The signal will also experience losses equal to Afiber=exp(aLrtr) and random polarization fluctuations.
If the signal is directly detected on a photodiode, then the photodiode will only reproduce the microwave power envelope of the modulated optical carrier, that is 56 = 7$8/9 • cosF )* + %!$./ F − :9!: H − )* :9!: H Where ΙPD is the photocurrent, after removing the DC component, and R is the responsivity of the photodiode. The photodiode output is sent to a microwave mixer and is combined with the local RF carrier that modulates the MZM. The output of the mixer will be the following: !"# = cosK )* + %!$./ F − :9!: H − )* :9!: Lcos ( )* + %!$./ ( )+ .;$7#/9 ) (4) Using the well-known formula cos( ) cos( ) = <=>(?-8)2<=> (?28) ' we can easily deduce that the term of interest that remains after low pass filtering the output of the mixer is the term of cos(a-b), hence: !"# = cos ( )* :9!: − %!$./ F − :9!: H + %!$./ ( )+ .;$7#/9 ) If the microwave source has a high coherence length (small linewidth), then one may safely assume that %!$./ F − :9!: H ≈ %!$./ ( ). This is easily obtained for microwave frequencies of high spectral purity that can preserve their coherence for distances > 5000 km. Using for instance the well-known formula of coherence length A!; = B $ CD7 , where Δf is the linewidth of the microwave source and assuming that its value is close to 1 Hz or even below (actually tens of mHz for a high stability 10 GHz carrier), then Lcoh > 63.000 km which proves that the system is immune to source phase noise even for transoceanic links whose length is usually below 15.000 km in the round-trip.
The signal is imprinted on Tprop which is disturbed by seismic events as explained in the manuscript. The amplitude V depends on the total amplification of the detected signal, both in RF (prior to the mixer) and low frequency regime and the RF power injection to the Local Oscillator port of the microwave mixer.

Photo-elastic effect
The photoelastic effect relates the change of refractive index to the mechanical strain, described by the two strain-optic coefficients p11 and p12. By measuring the phase change of light propagating through the fiber due to a longitudinal strain, we can determine the individual values of the two coefficients p11 and p12. Using first-order theory of elasticity and the photoelastic effect, and assuming that the fiber is elastic and mechanically homogeneous, the phase change Δφ induced by an elongation ΔL is given by 2 : Δng is the refractive index change and v is the Poisson's ratio of the fiber material. The parameter = 1 − is the strain coefficient of fiber due to the photoelastic effect and is approximately equal to 0.78 3 .

Dispersion effects
Chromatic dispersion is one of the most significant deteriorating effects in fiber transmission that can cause temporal broadening and power fading. A microwave frequency can be superimposed on an optical carrier in many different ways as depicted in fig. 2. When both sidebands of the microwave frequencies are transmitted as in fig. 2a then the two sidebands acquire different phases due to dispersion. When they interact on the photodiode, the different phases result in power fading 4 . The expression providing how power fading evolves with respect to transmission distance L, second-order dispersion parameter β2 and ωRF is the following Using this simple formula, one can easily find how optical power fading evolves as a function of RF frequency for a given dispersion value (β2=-21 ps 2 /km in Single Mode Fibers) and L=2500 km and L=500 km. Based on fig. 3 one can easily deduce that power fading has different periodicity at different lengths and things become more challenging in terms of preserving the maximum power as the RF frequency increases. However, with the use of a Phase locked loop synthesizer as the one utilized in our experiment, one can set the microwave frequency at the nearest exact maxima for any given L and preserve the stability of the RF frequency in the level of 100s Hz for a long time period (2 ppm PLL reference crystals). The periodicity of maxima for 500 km of transmission is approximately 375 MHz @ 20 GHz and for 2500 km of transmission is approximately 75 MHz @ 20 GHz. Therefore, a general purpose PLL synthesizer can easily lock the frequency at a point of maximum even for very long distances.
Things become much easier when SSB (Fig. 2b) or DSB with suppressed carrier (Fig. 2c)   modulator by properly adjusting bias and RF modulation depth as explained above. In both cases, only two waves interact at the photodiode, therefore fading is eliminated. The use of DSB with suppressed carrier has a two-fold advantage 1 : on one hand it eliminates power fading, on the other hand it results in microwave frequency doubling, which will enhance the phase resolution of the system according to (4).
The following experimental results show how DSB with suppressed carrier copes with power fading 1 Fig. 4: Power fading as a function of frequency which clearly shows that carrier suppressed DSB is immune to power fading in contrast to DSB with unsuppressed carrier (taken from 1 )

B. Processing of measured data
One of the most significant units at the reception state is the one undertaking system monitoring and realtime calibration and processing of the data in order to seamlessly maximize the phase signal. At the reception stage, the microcontroller continuously monitors the output signal and calculates its real-time deviation with respect to the reference voltage level that corresponds to approximately π/2 phase difference between the microwave generator and the received signal, which maximizes sensitivity. When the real-time error voltage overcomes a specific threshold, a suitable control byte is sent to the HMC642ALC5 microwave phase shifter and the resulting phase transition restores operation close to maximum sensitivity. In addition, a flag byte (outside the range of the output voltage signal) is written by the microcontroller in the output data, indicating the point where phase wrapping occurs. During post-processing, a suitable algorithm localizes the flag bytes in the acquired digital signal and unwraps phase by removing the voltage transition (the control byte) introduced during reception. Fig. 5a shows a snapshot of a raw signal recorded at the receiver where phase corrections are apparent as abrupt vertical transitions. In addition, specific bytes are used as flags to the recorded signal to ensure errorless phase unwrapping. Fig. 5b presents the phase unwrapped signal after post-processing. The arrow depicts the phase variation induced by the earthquake of Crete as an example (October 12, 2021, 09:24:03 UTC, (ML=6.3)) which is visible on the raw data before any post-processing. Step by step processing of phase raw data -the phase shifter provides a continuous correction of the phase shift so as to maximize sensitivity (2000 mV correspond to π/2). a. Raw data, b. data after phase unwrapping, c. strain rate, d. 1 st derivative of strain rate.
The output signal which is proportional to phase, is straightforwardly converted to strain by the equation: ξ@0.78 is the strain coefficient of the optical fiber due to the photoelastic effect. Further processing converts strain to strain rate ( JL J# ) or acceleration ( ' / ' ) for comparison with seismometer or accelerometer recordings respectively (see Fig.5c, Fig.5d respectively). It must be noted that such an active correction mechanism is more difficult to succeed when the phase of the link is interrogated using optical instead of microwave frequencies, due to the rapid changes of phase as a result of environmental effects (optical phase variations exceeding decades of π per second are expected only due to temperature variations). Therefore, despite the fact that microwave frequencies are characterized by lower resolution compared to their optical counterparts, they are also more robust to link noise. This is one of the underlying reasons of MFFI efficacy in providing real-time measurements and consequently event detection in ms scale, despite its use in a noisy terrestrial link positioned in a densely populated area. We anticipate that its performance will be better in less noisy submarine links.

Length of the sensor and its relation to measurement bandwidth
The length of the sensor is not directly related to the useful measurement bandwidth, although it is related to source phase noise (see in methods of the manuscript where fiber loop bandwidth, 4 kHz in our case, is investigated in terms of microwave frequency phase noise). The bandwidth of the final signal is mostly related to the ADC capabilities and the properties of the earthquake. For instance, a distant epicenter will provide signals that reside at lower frequencies (< 5 Hz), whilst a local epicenter will provide earthquakes that acquire higher frequency components (up to 10 Hz). Depending on the information content of the signal, we apply proper filtering in order to highlight the signal and suppress noise.

MFFI -DAS comparison
DAS, as a result of its distributed nature, can estimate strain rate at many locations along the fiber. On the contrary, MFFI provides measurements of the strain which are differentiated to provide strain rate accumulated over the entire link. Therefore, by averaging DAS values we see the agreement with MFFI differentiated signal. Despite the fact that the two systems have analog to digital converters (ADCs) of different sampling rates and resolution, we finally resample the signals at the same rate (50 Hz) and we filter them with the same band-pass filters (0.1 to 1.5 Hz for the case of Crete's earthquake). Their agreement is obvious.

C. Recording of significant earthquake events
We selected more than ten significant earthquakes that took place during the system testing period (July 2021 to February 2022). In this document, we present the recordings from earthquakes with epicenters in Marousi, Ikaria, Kefalonia, Florina, Chalkidiki with epicenters depicted on the map in Figure 2, (left side) in the manuscript. The comparison between MFFI and ATHP data for the selected events are shown below. In all cases a satisfactory agreement is observed, especially noted for the horizontal components. It must be mentioned that MFFI records the seismic wave effect across a 25 km long link, whilst ATHP accelerometric station records the seismic wave at a specific point location. Our next goal is to install accelerometers along the specific link and study the correlation of the accumulated strain provided by the fiber with the superposition of the responses of the in-line accelerometers. This study will provide a better understanding of the strain distribution along the line and will guide us in transforming the MFFI technique to a distributed measurement technique.

D. Phase noise of the terrestrial link and MFFI building blocks/stages
Here we present a study of the phase noise of the overall system including the fiber link, the optical and optoelectronic components (optical amplifiers, laser, and photodetector) and the noise of the electronic components (PLL, electrical amplifiers, mixers). We perform the analysis by measuring the one-sided spectral density of phase noise in dBc/Hz as depicted in Fig. 11. In all cases, the output signal is the one digitized by the 10-bit embedded ADC of a low-cost micro-controller (Arduino Uno). Green curve: The phase noise at the output of MFFI is measured when the entire link is included. The highlevel spectral density at low frequencies is mostly attributed to link noise (temperature fluctuations, human activity) Blue curve: The phase noise at the output of MFFI excluding the transmission link which has been replaced by optical attenuators in order to regulate the OSNR at the same level as in the green curve. It is obvious that phase noise has been reduced dramatically, although we observe strong components at low frequencies mostly due to the disturbance of optical components (patch-cords, EDFAs, etc) by mechanical/thermal vibrations inside OTE Academy building. Even at very low frequencies the spectral density of phase noise has been reduced by 20 dB. Black curve: Here the phase noise is attributed only to the electronic part of the system; it corresponds to the noise floor of the measurement system. Hence the Arduino sampler quantization noise sets the lower limit of phase noise power spectral density approximately to -95 to -100 dBc/Hz. When we measured the phase noise of the MFFI electronic subsystem with the use of an oscilloscope ( Fig.  11 -dashed line) we obtained an improvement of 15-20 dB in the noise floor as depicted, which proves that a better sampling system could further improve MFFI's sensitivity, especially in less noisy environments. This is an improvement that will be applied in the next weeks, by substituting the Arduino 10-bit ADC of 50 Hz with a 14-bit 700 Hz ADC.

E. Localization techniques -converting MFFI to a distributed sensor.
Localization of events is very significant for seismology applications. The following scheme shows the principle of operation already demonstrated in ref.
[26] and [51] of the manuscript.  If we rely only on a single MFFI interrogator similarly to [51], we could localize the events as follows: Let us assume that a vibration event affects the fiber at distance L1 from end A. Then this event will be also repeated after 2L2 km as depicted in the figure above (L1+L2=L of the link). This length difference will be detected as a time difference in the autocorrelator unit of end A. The time difference will be Δt=2L2/vg. This technique has been evaluated in ref. 51 of the manuscript and showed that hundreds of meters resolution is possible. It is true that DAS can offer resolution in the order of a meter, however this is quite necessary in other types of applications related to intrusion detection, structural health monitoring, etc.
In the field of seismology, the minimum wavelengths that we are dealing with are on the order of several hundred meters or several kilometers. Hence, with a resolution of around 100 m, we would still be able to sample the wavefield completely. Therefore, 100 m resolution is considered adequate.
An extension of this technique is to use simultaneously both A end B MFFI systems. At end B, the same event will provide an autocorrelation signal with Δt=2L1/vg. Thus, a cross correlation of the signals captured by end A and end B will also reveal the position of the event. Marra et al in ref. 26 utilize and demonstrate localization in this manner.
The spatial resolution depends on the ADC sampling rate, that is if we want to target ΔL resolution, then we need an ADC with temporal resolution lower than ΔL/vg. For 10 meters of spatial resolution we need an ADC of more than 20 MHz bandwidth. The accuracy of the measurement and thus spatial resolution can be enhanced through averaging [26].
The authors also study other techniques towards the goal of transforming MFFI to a distributed sensor (see ref. 52, 53 of the manuscript), which however tackle other mechanisms which are beyond the scope of this work.