Fast information acquisition using spectra subtraction for Brillouin distributed fiber sensors

The traditional methods of extracting sensing information of a Brillouin distributed sensor is by curve fitting the Brillouin gain profiles along the fiber, we propose two concise and time-saving methods to process signals from only the frequency shifted section(s) of the fiber by subtracting the original spectrum from the sensing Brillouin spectrum. Experimental results validate that our methods can provide up to over 9 times faster information acquisition for a 10 km sensing fiber with 800 m frequency shifted section comparing with the traditional way. The proposed methods could be potentially attractive in getting information for the Brillouin distributed sensors as well as the Raman and Rayleigh distributed sensors especially when the processing speed is concerned. © 2019 Optical Society of America under the terms of the OSA Open Access Publishing Agreement


Introduction
With distinguish advantages of high resolution, easy deployment, temperature and strain monitoring capability, Brillouin distributed fiber sensor is one of the main-stream scientific research topics of fiber optic sensors, and it has become commercially available in the last decade [1,2].A great amount of research focuses on the front end technologies, i.e., optical technologies, as coherent BOTDA [3], dynamic Brillouin [4], Differential Pulse Pair [5], prepulsed pump [6], sweep free BOTDA [7].In the last couple of years, some work tackling the ever growing distributed sensing data were reported.Many common signal processing techniques are employed to denoise the spectra, including, wavelets denoising [8], non-local means denoising [9], adaptive linear prediction and cyclic coding [10], et al.In Nature Communications, Soto reported an image processing method to handle with the data, which claimed to have improved 100 fold sensor performance [11].Some papers are aimed to process the spectra to obtain the sensing information efficiently as well, such as, continuous wavelet transform for phase optical time-domain reflectometer (Φ-OTDR) [12], twodimensional edge detection method to improve Signal-to-Noise Ratio (SNR) and spatial resolution [13], iterative subdivision method to improve spatial resolution [14], online strain profile estimation method with two-stage adaptive algorithm [15], artificial neural network [16], principal component analysis based pattern recognition [17], as well as similarity comparison to obtain Brillouin frequency shift [18,19], et al.
A major obstacle to further reducing the processing time for the above-mentioned methods is that they all need to calculate full spectrum to acquire sensing information.In practice, the Brillouin gain spectrum is normally obtained by Lorentzian curve fitting point by point and then the sensing information is retrieved from the central Brillouin frequency [20].However, on the way to a higher solution, especially with 1 million points and beyond [21][22][23][24], to in-site real-time dynamic measurements [25][26][27], and to extra fine frequency resolution [28][29][30], the curve fitting process takes a very long time and it is neither target oriented nor effective, which makes those goals hard to achieve and some even unrealistic.Actually, not every single point along the fiber needs to be curve fitted, because a large range of the fiber's Brillouin frequency may only have jitters in a measurement, which is not necessary to be calculated, and only those points with frequency changes carries information and need be fitted and calculated.
In this manuscript, we propose two methods based on straight-forward spectrum subtraction which subtracts the measurement spectrum from the un-shifted spectrum before the spectrum is fitted.With this technique, a portion of the regular time with the traditional processing algorithm is needed to acquire Brillouin frequency shift information depending on the frequency shifted section length.Theoretical analysis on how to acquire Brillouin shift and to reduce the algorithm complexity is carried out before experimentally verification.This method could hopefully be adopted with the state-of-the-arts methods, such as, neural network, similarity method, as well as in other distributed sensing schemes, for example, Raman and Rayleigh sensors, helping them to retrieve distributed sensing information faster.
To figure out the Brillouin frequency shifts (BFS) ν along the fiber, that is to find out differences between the original unshifted spectrum and, the strain and/or temperature changed spectrum, or often from the latter solely using itself as a frequency reference.Since no frequency-shifted points can be tracked directly, the BFS is commonly obtained point by point along the fiber, including the points with no frequency shift.We propose to subtract the reference spectrum from the measurement one before fitting, and the differential action could help to distinguish the shifted and informative portion of spectrum from the unchanged.Then one could employ either the traditional fitting or calculate to obtain the BFS from the frequency shifted section(s) of the subtract spectrum, above two methods are named the Subtract and Fit Method (SFM) and the Subtract and Calculate Method (SCM), respectively.The two are both based on subtracted spectrum but uses different way to obtain the BFS, and them along with the traditional fitting algorithm are further discussed and compared with experimental data in later paragraphs.
A simulation is performed on a 10 km fiber whose central Brillouin frequency is 10.84 GHz.As an example, the fiber's 9-10 km section has 20 MHz frequency shift, and Gaussian white noise is added to the spectrum to simulate 25 dB SNR at the tail end.Those parameters are chosen to facilitate comparison with the experimental results.Due to the subtraction, a big portion depends on the measurement situation and the circumstance of the stimulated Brillouin spectrum (SBS) could be removed, and the remaining part, as shown in Fig. 1(c where ∆ν is the sensing induced frequency shift.Define ∆ν PD as the frequency difference of the peak and dim of the subtracted curve (See Fig. 2(a)) of a sensing point, the original and shifted gain profiles at that point are normalized and simulated as in blue and red respectively.∆ν PD is determined by ∆ν and the FWHM ∆ν B of the gain profile, and the relationship is depicted in Fig. 2(b) and in its inset.It shows that ∆ν PD is always larger than ∆ν, and this is especially significant when the shift is small.To find out ∆ν with Eq. ( 2), one need to find out the peak and dim positions by differentiating Eq. ( 2) and then let the result to be zero as Find the distance of the real roots of Eq. ( 3), which is ∆ν PD .Then the Brillouin frequency shift ∆ν can be further deduced as: ∆ν B is previously known and ∆ν PD is available from the subtracted spectrum.Then in the SCM, ∆ν can be figured out with the above equation, and the complexity of the algorithm is reduced to minimum as no fitting is needed.However, it could be inaccurate for a noisy spectrum, thus discrete Meyer wavelet de-noising is employed in our algorithm.The simulation results of the SCM are illustrated in Fig. 1(d), ∆ν by the Lorentz fitting is in blue while the subtraction method calculated frequency shift is in black.The insets show that the Lorentz method fits every point but actually gets noise for those un-shifted points, nevertheless, the SCM provides cleaner un-shifted region but slightly more noisy in the sensing area.
The SFM uses a more straightforward way, which is to use the subtracted spectrum as a criterion to obtain shifted points indices and carry out Lorentzian fitting on the frequency shifted points of the sensing spectrum as in the traditional way.This method comprises between processing speed and sensitivity, and the simulated SFM deducted Brillouin central frequency is drawn in red in Fig. 1(d).

Threshold for the subtraction algorithms
For both the SFM and the SCM, the points with value not equal to zero should be fitted if this is an ideal case without noise.However, in order to pick out the points with BFS from a noisy background, a proper threshold T is needed to distinguish subtracted gain from noise.Following the noise definition in [17,19,32], the variance along the fiber is calculated.In this manuscript, the noise is assumed to be similar along the fiber, and for convenience, T is defined as four times (i.e.SNR = 6 dB) the noise variance of the pre-fiber section of the spectrum.

Complexity of the algorithm
The traditional curve fitting method has a computational complexity on the order of O(r 1 N 2 ) for a spectrum with N points along the fiber [33], while the complexity of the SCM for the same object is on the order of O(N), and the complexity of calculating part of the spectrum as done in the SFM is O(r 2 N 2 /r 3 ), where the r 1 and r 2 is the number of iterations of the method, r 3 is the portion factor, which is the ratio of the whole fiber length over the shifted section(s) length.
The complexities of these algorithms for different numbers of resolved spatial points are theoretically calculated and illustrated in Fig. 3(a) in terms of processing time for a 10 km fiber with 9-10 km of it under 20 MHz strain.To fit all the curves in one figure, the Lorentzian fitting time curve is divided by a factor of 5.The same fitting algorithm is used in both the SFM and traditional methods to keep the results comparable.As shown by Ratio 1 in Fig. 3(a) that the SFM algorithm is about 10 times faster than that of the traditional way, while the SCM is even 3 times more effective than the former.With smaller numbers of spatial points processed, Ratio 1 go down due to the overhead calculations, but Ratio 2 stays relatively unchanged due to the fact that wavelet filtering consumes most of the time.Simulation of a 10 km fiber (N was set to 3000) with different lengths having a 20 MHz SBS shift is shown Fig. 3(b), the processing time and complexity of the SFM go up linearly with increasing frequency shifted fiber length, which is proportional to the fitted point number N. We find out that when the frequency-shifted fiber length is 8 km, that is about 80% of the whole length in this case, the processing time is similar for the SFM and traditional fitting method.Furthermore, SFM consumes 18.9% more time for a fully shifted fiber due to the fitting decision mechanism.Therefore, it is less effective to use the SFM if there is a larg processing tim

Experime
To verify the system was se 10 km long centered at 10 was set to b temperature, i resolution) in collected by a Lase   The three methods were utilized to find out the BFS of the experimental data, and their results are illustrated respectively in Fig. 6(a)-6(c), the insets in those figures show that both the SFM and SCM results is less noisy than the traditional one in the unshifted section while the output of the SCM is more noisy than the other two in the shifted section.Figure 6(d) shows the variance of the data.Since the SFM and the SCM use a threshold to determine the fitting, they have much higher efficiency and largely reduce the amount of the points to be processed by finding out the frequency shifted section(s).The SFM keeps the same frequency shifted section details and accuracy comparing with the traditional way.However, not as expected, the original SCM doesn't work well with the presence of noises as it would be difficult to accurately obtain the peak and dim position.Therefore, it needs to be filtered and then fitted to find out the peak and dim values of the subtracted curve, which greatly slows it Then the processing time, as well as the accuracy of the methods were evaluated.Three data sets (40 -30, 50 -40, 50 -30 averaged for 1000 times, corresponding to Time 1, Time 2 and Time 3 in Table 1) were calculated for 20 times each to find out an average processing time and accuracy for the three methods.Figure 7(a)-7(c) shows a boxplot of the processing time, and Fig. 7(d) depicts the time ratios of the three.Table 1 illustrates the time consumption and accuracy figures.The total average time reduction for the SFM is 89.04% with no accuracy lose, while it is 86.25% for the SCM with 0.26% accuracy drop.

Discussion
The SFM we proposed can obtain the Brillouin frequency shift from the spectrum saves 89% of the processing time than the frequently used Lorentzian fitting method with the same accuracy, thanks to the idea of locating the frequency shifted sections.While the SCM can be used to get the results in a fast way but not as fast as expected in theory, it is because of filtering and fitting needed to fight the noises and the unwanted sub-peak of the fiber.
Comparing with the traditional way, our proposed methods are more vulnerable to noise, attributed to the lower signal power after deducting the common signal and the deteriorated SNR comparing with original signal.Therefore, in an extremely low SNR case, the traditional way is more preferable to retrieve the BFS.However, for a distributed sensor with a reasonable SNR, or if the induced frequency shift is significant, our proposed methods could save a large portion of the processing time, because instead of using all data fitting, the Brillouin peak shifted parts could be identified with a fitting threshold before being processed.

Conclusion
We have provided a fast subtraction technique to calculate the BFS of the fiber.Rather than fitting all the points along the fiber in the existing methods, only those frequency shifted points are fitted or calculated in the SFM and the SCM methods we proposed.In theoretical section, we obtained the relationship of the subtracted peak-dim distance and the actual Brillouin frequency shift, then discussed the algorithm complexity.The proposed methods are then verified by experimental data, which have shown that the SFM is about ten times faster than using the curve fitting way to process the Brillouin spectrum with the same accuracy, while the SCM is over 7 times faster with a small accuracy drop on average.The proposed methods could be easily combined with the front end techniques and they can be further applied to the Raman and Rayleigh sensors as well.

Funding
Beijing Jiaotong University Fundamental Research Fund Grant No. 2018JBM018; National Natural Science Foundation of China Grant No. 61805008, No. 61435006.
) which is the result of Fig. 1(b) minus Fig. 1(a), should theoretically be depicted in the form of

Fig. 1 .
Fig. 1.(a) Original spectrum of 10 km fiber, (b) 9-10 km section at the far end with + 20 MHz shift, (c) subtracted spectrum, and (d), frequency shift from the traditional Lorentzian fitting (blue), the SFM (red) and the SCM (black).Insets: Details of the curves.

Fig. 2 .
Fig. 2. Theoretical curves.(a) Demonstration of subtracting two Brillouin gain curves, (b) Frequency difference between peak-dim distance and sensing induced frequency shift, ∆ν.Insert: Peak dim distance relationship with the induced frequency shift.

Fig. 3 .
Fig. 3. Time consumption for the SCM, SFM, and traditional methods for (a), different spatial points number N, (b), different length of fiber with sensing induced frequency shift, ∆ν.

Fig. 6 .
Fig. 6.Frequency shift with (a) the traditional Lorentzian fitting, (b) the SFM, and (c) the SCM for 20 C temperature change, (d) is variance of the subtracted spectrum.The insets are details.

Fig. 7 .
Fig. 7. Time costs for the Lorentz fitting, the SFM and the SCM.
and SCM (a.u.) 40-30 L/SFM 40-30 L/SCM 50-40 L/SFM 50-40 L/SCM 50-30 L/SFM 50-30 L/SCMdown and its accuracy is also influenced by the filtering.The error also comes from the fiber's peak at ~10.92 MHz (as demonstrated in Figs.5(d) and 5(e)), which caused an unwanted peak in the subtracted curve in Fig.5(f), resulting in deviations from the subtraction curve in our model.