1. Introduction

Optical fibers based distributed sensing is gathering increasing attention due to its many practical applications in safety-related fields like geotechnical monitoring [1] and intrusion detection [2, 3], or harsh environments, like cryogenic temperature monitoring [4]. This interest is justified by the characteristics of the optical fibers, namely the resistance to extreme temperature and stresses, the easiness of distributed measurements and the reduced size and cost.

Among the three scattering phenomena occurring in a fiber that can be exploited to realize a distributed optical fiber sensor [5], the easiest to observe is the Rayleigh scattering. Rayleigh scattering is an elastic scattering sensible to a wide range of static and dynamic perturbations, which can be interrogated by implementing time or frequency domain schemes [6]. In phase-sensitive optical time domain reflectometry (ϕ-OTDR) [7], the coherent optical source is modulated into a sequence of pulses that consecutively interrogates long sections of fibers, generally with spatial resolution of several meters. In optical frequency domain reflectometry (OFDR) [8] the central frequency of a strongly coherent continuous wave laser is linearly swept in time over a wide band enabling the scan of tens of meters of fiber with sub-millimeter spatial resolution [9]. Originally these techniques were developed to monitor static parameters [10, 11], but algorithms for distributed acoustic sensing (DAS) were reported [12–14] achieving acoustic bandwidths of few kHz. To increase the performance of the DAS different approaches were considered, like realizing polarization-diverse schemes dealing with the issue of signal fading [15] or modulating the central frequency of the laser [16, 17]. Eventually the DAS performance can also be improved in post-processing by implementing proper denoising algorithms [18], but the nature of the noise affecting the measurement must be analized throughfully in order to correctly implement such techniques.

Recently DASs exploiting the phase information of the Rayleigh backscattering trace were reported, with algorithms already achieving performance similar to the intensity-based counterparts [19]. The analysis of the phase of the backscattered trace is attractive because any perturbation acting on the fiber causes a corresponding distributed phase modulation, linear to the dynamic of the stress applied [20]. It is then possible to extract informations about the stress intensity other than just the dynamics. Unfortunately, any perturbation occurring on a particular position affects also the backscattered light coming from further points of the fiber, making all the subsequent contemporaneous perturbations hard to identify [21].

In this paper we propose a data analysis algorithm based on the OFDR setup that solves the issues of overlapped perturbations enabling the distributed detection of multiple high frequency stresses. High acoustic bandwidth of tens of kHz can be observed with spatial resolutions of centimeters. The algorithm performs a specific comparison of the phase of a reference OFDR trace, recorded before the perturbation, with the phase of a trace recorded during the perturbation. Moreover, a key aspect of the technique is a novel active compensation of the laser phase noise, described in Sec. 4. The phase noise affecting the measurement is analyzed and the performance of the sensor are evaluated by computing the SNR of a signal in a progressively noisy environment.

2. Algorithm description

Let’s consider an ideal fiber at rest, with a strongly coherent light signal $s_{in}\left(t\right)$ propagating along $\widehat{z}$ direction, characterized by a propagation constant $\beta \left(\omega \right)\approx {\beta }_0+\omega {\beta }_1$. The fiber is considered short enough to neglect chromatic dispersion. The light Rayleigh-backscattered by the fiber can be modeled as the superposition of the reflections due to each scattering center [22]:

If a time-varying stress $\epsilon \left(t\right)$ is applied on a section of the fiber $I=\left[z_1,z_2\right]$, of length $\delta z=z_2-z_1$, the elasto-optic effect causes a local change $\mathrm{\Delta }{\beta }_0\left(t\right)\propto \epsilon \left(t\right)$ on the propagation constant:

The change in β₁ is negligible for any practical stress intensity, therefore the stress changes the phase term ϕ_n in (3) as follows:

Notice that the $b(\cdot )$ function appearing in the first and third contributions is exactly the same of the unperturbed fiber of Eq. (4), because the fiber span before section I is unperturbed, while the fiber span after section I is entirely affected by the same z-invariant phase perturbation $\mathrm{\Delta }\varphi \left(t\right)$. The Rayleigh signal generated by an arbitrary fiber section beyond section I can be easily obtained by applying a proper bandpass filter with impulse response $w\left(t\right)$ to the beat signal. Then, the perturbation can be calculated as:

If K localized dynamic perturbations act on the fiber in different positions $z_{s,k},k=1,\mathrm{...},K$, the phase variations $\mathrm{\Delta }{\varphi }_k\left(t\right)$ induced by each stress sum up, causing cumulative phase change terms that read as:

When the filter bandwidth is completely inside the spectra spanned by $b\left(t;z_{s,k},z_{s,k+1}\right)$, the extracted phase coincides exactly with the cumulative phase terms reported in Eq. (10), otherwise some undesired distortions arise. However, by subtracting consecutive exact cumulative phase terms the single dynamic can be recovered.

The bandwidth B_w of the bandpass filter, expressed as a spatial width $L_w=c_0{B}_w/\left(2n_g\right)$, defines the minimum spatial separation between any two consecutive perturbations to correctly extract at least one replica of each exact cumulative phase change ${\displaystyle \sum }_k\mathrm{\Delta }{\varphi }_k\left(t\right)$. If the distance between any couple of consecutive perturbations is greater than L_w then more replicas of such phase changes ${\displaystyle \sum }_k\mathrm{\Delta }{\varphi }_k\left(t\right)$ are extracted. The phases extracted at this point of the analysis are the cumulative effect of the perturbations acting along the fiber. The local information can be extracted by subtracting the phase terms calculated from two different windows separated by a distance $\mathrm{\Delta }L$. As a result, the parameter $SR=\mathrm{\Delta }L+L_w$ represents the spatial resolution of our algorithm, since it defines the distance required to efficiently discriminate two consecutive perturbations. An increase of $\mathrm{\Delta }L$, which worsen the spatial resolution, allows for a SNR increase since a greater number of measured perturbation dynamic replicas can be averaged. The acoustic bandwidth achievable with this algorithm is $B_{a,max}=B_w/2=n_g{L}_w/c_0$ and can be computed from the filter bandwidth according to the Nyquist theorem. A larger L_w guarantees an higher $B_{a,max}$ at the cost, however, of a worse spatial resolution.

If the perturbations acting on the fiber are not localized, the terms $\tilde{b}(\cdot )$, neglected in the previous analysis, becomes significant. The phase change terms $\mathrm{\Delta }{\varphi }_k\left(n,t\right)$ cannot be factorized, due to the dependence on position, resulting in an incorrect ratio operation in Eq. (12) which introduces distortions on the amplitude and dynamic of the extracted signals. When the stress intensities are not too high, such distortions are limited and the results can still be approximated with a sequence of phase terms with linearly increasing intensity, due to the progressively longer section of perturbed fiber considered. The dynamics and the positions of the perturbations can still be extimated with sufficient precisions.

3. Phase noise analysis

The performance of the algorithm can be improved by averaging among replicas of the same perturbation measured at different fiber locations. The number of replicas available for averaging is $N_R=\mathrm{\Delta }L/\delta L_w$ and depends on $\mathrm{\Delta }L$ and on the spatial step of the moving bandpass filter $\delta L_w$. Assuming the noise affecting the different traces to be white, a smaller step should provide greater improvements at the cost of an increased computational time. Unfortunately any stress acting on the fiber, including undesired ones, generates a phase change term which sums coherently with all the others along the fiber, causing a non negligible degree of correlation even among traces extracted from well separated positions. This is expecially true if the undesired perturbations are relatively large, since they will cause strong and well defined phase change terms that strongly correlate all the subsequent extracted phase traces and make averaging mostly useless. Small and uniformly distributed environmental stresses however can be assumed to be low and random enough that, along the fiber, the sum of all their contributions tends to cancel out. Without strong perturbations then the degree of correlation of traces extracted from increasingly separated positions will decrease and a proper step can be identified to perform a fast analysis with a meaningful averaging operation.

To verify experimentally how the degree of correlation of the noise changes, the reference and stressed OFDR traces of a bare fiber of length L = 3 m have been measured without applying any perturbation, i.e. with only environmental noise acting, and the analysis has been performed setting L_w = 15 cm, $\delta L_w=0.1$ cm and $\mathrm{\Delta }L=15$cm. The noise correlation is computed among traces extracted in five increasingly separated positions of the fiber, respectively after $d_n=\left[0.1,1,5,10,15\right]$cm, $n=1,\mathrm{...},5$, with respect to the initial one located at the beginning of the fiber. The autocorrelation of the initial trace was computed to normalize, with respect to it, all the other values. No more positions have been investigated since the differentiation operation in the algorithm will remove any correlation in traces extracted from positions separated by a distance greater or equal to $\mathrm{\Delta }L$. The process was repeated along the whole fiber, by shifting forward the initial trace with steps of $15$cm, allowing the analysis on twenty non-overlapped sections. The results are reported in Fig. 1 where it can be seen that the behavior of the correlation degrees remain the same along the whole fiber. Traces extracted from positions $1$mm away are strongly correlated, eventually making the averaging operation less efficient, while a separation of $1$cm exhibits a very low correlation. For greater distances the degree of correlation decreases to zero as expected, but the number of traces collected becomes too low to guarantee a satisfactory noise reduction. For all those reasons the value of $\delta L_w=1$ cm was selected as the optimum one.

The other sources of noise that affect the system, mainly receiver noise and laser phase noise, can be modeled as uncorrelated white processes, for which the effectiveness of averaging is largely independent of $\delta L_w$.

 

Fig. 1 Noise degree of correlation.

Download Full Size | PDF

4. Frequency sweep non-linearity correction

The standard OFDR setup is based on a highly coherent source, linearly swept in frequency over a wide band, which generates the signal of Eq. (2). The linearity of the frequency sweep is usually guaranteed by an auxiliary interferometer which output allows to precisely trigger the data acquisition. As a consequence of this correction, the time scale of the measured data is distorted and not consistent with the time scale of the perturbation. If the data of the linearizing interferometer were available, this inconsistency could be corrected. However, this is not the case with commercial OFDRs.

To remove this distortions we stressed a section of fiber at the beginning of the link with a PZT driven by a known sinusoidal perturbation $\xi \left(t\right)={\xi }_0\sin (2\pi f_{p}t)$. The algorithm proposed previously can then be used to extract the distorted spectrum of such perturbation $\mathcal{F}\left\{\mathrm{\Delta }{\varphi }_p\left(t\right)\right\}=\left[S_b\left(f-f_p\right)\right]$, which will be a shifted version of the spectral broadening due to the spectral properties of $\xi \left(t\right)$. An homodyne demodulation scheme with a proper lowpass cut-off frequency, that should be low enough to cut-off the noise, but large enough to include all of the stress components, can then successfully retrieve an estimate of $S_b\left(f\right)$ and, by anti-transforming, an estimate of the undesired time variations. By resampling all the phase traces extracted from the subsequent spans of fiber, the spectral broadening can be effectively removed everywhere.

To experimentally verify the effects of the frequency sweep non linearity and the results of the proposed solution, a sinusoidal perturbation with a frequency arbitrarily set to $f_p=3$kHz was applied through a PZT in position $z_p=2$m of a $3$m long fiber. The algorithm parameters used where $L_w=15$cm, $\delta L_w=1$ cm and $\mathrm{\Delta }L=15$cm. Figure 2a shows the raw spectrum of the retrieved perturbation without any correction and, as can be seen, instead of a clean peak there is an evident spectral broadening with a 20dB-bandwidth of about $80$Hz. Through the homodyne demodulation scheme the spectral broadening is shifted back to baseband and a proper lowpass filter is applied. Considering our setup and instrumentation, a cut-off frequency $f_{co}=160$Hz, which is twice as large as the bandwidth of the observed perturbation, was selected. An instance of the undesired time sampling variations as function of the measurement time is shown in Fig. 2b. It is important to highlight that due to its intrinsically random nature, this is just an example and that it is thus not possible to find a general correction term. Finally, the perturbation spectrum after the resampling operation is shown in Fig. 2c where a clean spectral peak at f_p is clearly visible.

 

Fig. 2 (a) Perturbation spectrum before resampling. (b) Undesired sampling period variations. (c) Perturbation spectrum after resampling.

Download Full Size | PDF

 

Fig. 3 (a) Interrogation unit of the sensor. (b) Fiber loops to realize multiple identical localized perturbations.

Download Full Size | PDF

5. Experimental results

To experimentally demonstrate the possibility of sensing high frequency perturbations, thirty meters of a standard telecommunication fiber were connected to a commercial OFDR (OBR 4600, Luna Innovation) with a sweep rate $\sigma \approx 25$THz/s. To correct the frequency scan non-linearities, as explained in the previous section, a known reference stress was applied through a PZT to the fiber section located immediately after the instrument. Such setup, shown in Fig. 3a, represents the interrogation unit of the sensors and it is left unchanged during all the following experiments.

The first test was performed by placing a PZT, driven by a sinusoidal signal at $f_{s,1}=41$ kHz under the fiber in position $z_s=4.60$ m. The algorithm parameters chosen for this first measure were L_w = 40 cm, $\mathrm{\Delta }L=10$ cm, guaranteeing an acoustic bandwidth of $B_{a,max}=50$ kHz and an effective spatial resolution of$50$cm. The filter step was fixed at the optimum value of $\delta L_w=1$ cm, to guarantee a good number of traces for averaging while limiting the computational time. To assess the intrinsic performances of the technique, the reference and stressed traces were collected one right after the other, since the longer the time that elapses between the two measurements, the higher the environmental noise affecting the result. As can be seen in Fig. 4, which represents the spectrogram of the retrieved dynamic, the perturbation oscillating at $f_{s,1}$ is clearly visible around z_s. The extracted perturbation is spread over a length of 50 cm which corresponds to the effective spatial resolution achievable with the selected algorithm parameters.

The fiber has then been arranged in five loops with a common point placed above the PZT, as can be seen in Fig. 3b. The PZT was driven by an arbitrary signal generated as the sum of three sinusoidal waves oscillating at frequencies $f_{a,1}=2.5$kHz, $f_{a,2}=5$kHz and $f_{a,3}=12.5$kHz. The fiber layout ensures the application of the same perturbation to five different positions, and it allows to verify the capability of the proposed algorithm to discriminate between spatially distinguished, but spectrally overlapped stresses. To sharpen the spatial resolution, the filter spatial width was reduced to $L_w=25$cm, guaranteeing an acoustic bandwidth of $B_{a,max}\simeq 30$ kHz. The other algorithm parameters were left untouched resulting in a spatial resolution of $35$cm. After collecting the stressed and reference traces, the raw spectrogram was computed and reported in Fig. 5a, where the five perturbations oscillating at the three frequencies $f_{a,1}$, $f_{a,2}$ and $f_{a,3}$ are clearly discriminated. Even if the setup in Fig. 3b show that the same perturbation is applied to all fiber loops, the coupling coefficient between the PZT and the fiber is different for each loop, thus the experienced perturbation intensity will vary. This justifies the different amplitudes of the signals in Fig. 5a and is confirmed by measurements taken from the other end of the fiber. As can be seen in Fig. 5b, the amplitudes from the other end of the fiber appears almost identical to the ones in Fig. 5a, but with reversed positions.

 

Fig. 4 Spectrum measured for a localized high-frequency sinusoidal perturbation at 41 kHz.

Download Full Size | PDF

 

Fig. 5 Spectra of the phase measured when three different frequencies (2.5, 5 and 12.5 kHz) are applied simultaneously in five different positions along the fiber. Only the portions of the spectra around the main frequencies are shown for clarity. (a) Measurement taken in the forward direction. (b) Measurement taken in the reverse direction.

Download Full Size | PDF

To evaluate if the performance of the algorithm varies when multiple localized perturbations act on the fiber, two PZTs were placed under two different positions of a fiber $z_{s,1}$ and $z_{s,2}$, 1.5 meters away. The first PZT was used as an acoustic disturbance oscillating at the same frequency of the useful signal generated by the second PZT. The signal to noise ratio (SNR) was defined as the ratio between the power of the extracted useful signal with respect to the power of the noise and it was computed multiple times as a function of the acoustic signal-to-disturbance ratio (acoustic SDR), defined as the ratio between the power of the extracted useful signal with respect to the power of the disturbance, when the algorithm parameters $\mathrm{\Delta }L$ and L_w were varied. Our analysis started by fixing the filter spatial width to L_w = 15 cm and by computing the SNR variations for different values of $\mathrm{\Delta }L$. The curves were normalized by the value of SNR obtained for $\mathrm{\Delta }L=15$ cm and L_w = 15 cm, at maximum acoustic SDR. The result is shown in Fig. 6a, where it is possible to see that the SNR is almost constant with respect to the acoustic SDR. Since, as discussed in the previous section, all the traces in the spatial interval defined by $\mathrm{\Delta }L$ represent the same perturbation, increasing such parameter allows to average more terms, obtaining an increase in the SNR at the cost of a duller spatial resolution. Figure 6b shows the SNR variations for different values of L_w, fixing $\mathrm{\Delta }L=15$ cm and normalizing the results for the same value as in the previous case. Once again the curves appear almost constant with respect to the acoustic SDR, except for high values of L_w. The reason is that for a filter with a wide spatial width the environmental noise affecting the traces impact more than for smaller ones, limiting the SNR and increasing the measurement uncertainty. For small value of L_w however, the influence of the receiver noise affects significantly the measurement, resulting in a lower SNR.

 

Fig. 6 (a) SNR as function of the acoustic SDR for different values of ΔL; (b) SNR as function of the acoustic SDR for different values of L_w.

Download Full Size | PDF

 

Fig. 7 Spectrum measured for a 12.3 meters long high-frequency sinusoidal perturbation at 50 kHz.

Download Full Size | PDF

To experimentally verify the possibility of monitoring a spread stress the setup was then changed again including a fiber stretcher with an effective length of 12 meters. The displacement generated by the fiber stretcher is uniformly spread over the whole length of the fiber and depends on both the frequency and the amplitude of the driving signal. A sinusoidal driving signal $s_d=V_d\sin (2\pi f_{d}t)$ at frequency $f_d=50$kHz was applied and the analysis was performed by setting L_w = 50 cm and $\mathrm{\Delta }{L}_w=50$ cm, guaranteeing an acoustic bandwidth of $62.5$kHz with and an effective spatial resolution of $SR=1$m. To assess the intrinsic performance of the technique, the rms amplitude V_d of the signal was set to just $0.15$ V, so to avoid phase variations larger than $2\pi$ and the related problems of phase unwrapping. As can be seen in Fig. 7 the extracted signal is centered at f_d and the result of the differentiation operation is a $12$m long signal. Differently from the localized perturbation reported in Fig. 4, in Fig. 7 it is possible to see that the extracted perturbation is mostly flat, due to the linearly increasing phase change induced by the spread perturbation, but slightly more distorted both in amplitude and frequency. These distortions are caused by both the non constant coupling of the perturbation of the fiber and by the above-mentioned non ideal ratio operation that affects the performance of the algorithm. The strong distortions visible around $z=23$ m in Fig. 7 are caused by the final connector of the fiber stretcher and are not significant.

The algorithm proposed in this paper allows also the analysis of the magnitude of the stress. In the case of the fiber stretcher, the generated displacement depends on both the amplitude and frequency of the driving signal. This relation is generally expressed through the modulation constant $\mathrm{\Delta }\gamma \left(f\right)\left[\frac{rad}{V}\right]$, which represents the phase change per Volt induced at a given frequency. If the differentiation operation of the algorithm proposed is avoided, it is possible to relate the total cumulated phase change induced by the fiber stretcher $\mathrm{\Delta }\gamma \left(f\right)V_d$ when driven by a signal s_d to the phase term $\mathrm{\Delta }\varphi \left(t\right)$, defined in Eq. (6), extracted just before the final connector of the fiber stretcher:

Experimental measurements have been performed applying a driving signal with variable voltage and frequency $s_d^{i,j}\left(t\right)=V_{d,i}\sin (2\pi f_{d,j}t)$, with $V_{d,i}\in \left[0.05,0.10,0.15,0.20,0.25\right]$V and $f_{d,j}\in \left[20,30,40,43,47,50\right]$kHz. We tested more frequencies between 40 and 50 kHz due to the closeness to the resonance peak of the piezoelectric component and the consequent sharper increase of the modulation constant. The setting values of $\mathrm{\Delta }\gamma \left(f_d\right)$ were computed by using the datasheet of the fiber stretcher and they are reported in the first row of Table 1. For each frequency $f_{d,j}$ the extracted values of $\mathrm{\Delta }\varphi \left(t\right)$ at different driving amplitudes were fitted with a first order polynomial function; the resulting slopes are reported on the second row of Table 1 and the R² value of the fitting was always around 0.95. As can be seen the data are correlated, even if there is a non negligible difference between them. The main cause of this variance can be found in the ratio operations described in Eqs. (9) and (12) which amplify the phase noise at those points where the signal amplitude is particularly low due to Rayleigh fading. A partial support to this argument is given by the fact that the ratio between the actual modulation constant (first row) and the measured one (second row) is on average about 1.02, i.e. very close to the ideal unity.

Table 1. Comparison between the Nominal and Measured Values of the Modulation Constant of the Fiber Stretcher.

View Table

6. Conclusions

In this paper a data analysis algorithm for OFDR-based distributed acoustic sensing has been proposed. The algorithm is based on the processing of the phase of the Rayleigh backscattering which allows linear measurement of the perturbations acting on the fiber. The setup non-idealities, namely the non-uniform sampling of the received signal and the phase noise, have been analyzed and methods to compensate them have been described and successfully implemented. It is proved that the algorithm guarantees acoustic bandwidths in the order of tens of kilohertz with spatial resolution in the orders of few tens of centimeters. Multiple simultaneous perturbations can be discriminated both in positions and frequency. Finally the SNR of the sensor has been studied as a function of the algorithm parameters for different noise intensities.