Reduction of Moving Target Time-of-Flight Measurement Uncertainty in Femtosecond Laser Ranging by Singular Spectrum Analysis Based Filtering

: Femtosecond laser ranging has drawn great interest in recent years, particularly based on an asynchronous optical sampling implementation where a pair of femtosecond lasers are used. High precision absolute ranging either relies on tightly-phase-locked optical frequency combs (a dual-comb setup) or multiple averaging of the measurements from two free-running femtosecond lasers. The former technique is too complicated for practical applications, while the latter technique does not apply to moving targets. In this report, we propose a new route to utilizing a powerful singular spectrum analysis (SSA) ﬁltering method to improve femtosecond laser ranging precision for moving targets with acceleration. The SSA method is capable of separating complex patterns in signals without a priori knowledge of the dynamical model. Here, we utilize the basic SSA ﬁlter to extract the target trajectory in the presence of measurement noise both in numerical simulation and in the absolute ranging experiment based on a pair of free-running femtosecond lasers. The experimentally-achieved absolute ranging uncertainty of a moving target is well below 110 nm at a 200-Hz update rate by applying the basic SSA ﬁlter. This method paves the way to the practical applications of femtosecond absolute ranging for dynamic objects.


Introduction
Optical frequency combs (OFCs), a high precision metrological technology based on femtosecond lasers, have been widely applied in absolute distance measurement [1]. The OFCs can simply provide optical frequency calibration for continuous-wave lasers in multi-wavelength interferometric ranging [2][3][4][5]. Otherwise, they have been directly used for time-of-flight (TOF) ranging based on the ultrashort pulse trains or multi-wavelength ranging based on the equally-spaced come-lines [6][7][8][9][10]. In particular, a dual-comb configuration [11][12][13][14] emerges as a competitive absolute ranging method due to the sub-micrometer precision, high update rate (kHz-level) and dead-zone free operation, showing potential applications in large-scale manufacturing, aerospace engineering, as well as remote sensing.
A dual-comb ranging system implements a pair of femtosecond lasers with an offset repetition rate and utilizes asynchronous optical sampling (ASOPS) for absolute ranging. The two lasers are always tightly phase-locked to a common narrow linewidth continuous-wave laser, enabling the transition from micrometer precision TOF measurement to nanometer precision interferometric measurement [11]. However, the unprecedented measurement precision is at the expense of increased complexity of the system setup. For the sake of industrial applications, absolute ranging with a pair of free running femtosecond lasers provides a significantly simplified setup, thus attracting much attention in recent years [15][16][17]. Nevertheless, the technical noise, as well as the quantum noises inevitably degrade the distance measurement precision [16,17].
Several methods have been developed to increase the measurement precision of absolute range finders based on free-running femtosecond lasers. A dual-wavelength mode-locked laser allows two pulse trains with an offset repetition rate to oscillate in a single laser cavity, and the shared cavity effectively avoids the technical noises that are the common mode to the two pulse trains [18]. Quantum-limited timing jitter is another major contribution to measurement uncertainty. Its effect can be significantly reduced by proper engineering of the mode-locking regime of the femtosecond laser sources [16,19]. Besides the above experimental approaches, measurement precision can also be improved by post-processing. The most straightforward approach is multiple averaging (MA). Sub-micrometer ranging precision is allowed with hundreds times of averaging [15]. However, MA causes a loss of the update rate since it integrates several distance data points into their arithmetic mean. Moreover, MA will be limited when the measurement target is moving. A good remedy is to extract the acquired distance information by using a digital filter instead of MA. Nevertheless, most digital filters need to know certain priori knowledge of input signals, such as the mathematical models of varying signals, the noise intensity of the data acquiring system, as well as the signal bandwidth.
Basic singular spectrum analysis (SSA), a data-driven analysis, established in 1986 [20], synthesizing the ingredients of classical time series analysis, multivariate geometry, multivariate statistics, dynamical systems and signal processing [21], was designed for making an extraction of qualitative information from a noisy signal. A relatively short and noisy signal can be decomposed by basic SSA into a sum of components, each having a certain interpretation, e.g., signal or noise. The basic SSA method has been extending its application to multiple fields such as hyperspectral remote sensing [22,23], GPS data processing [24] and absolute distance measurement of stationary targets [25]. As an analytical method, the basic SSA method needs to be manually analyzed by analysts after the result is obtained, which leads to an experience-oriented category of signals and noises. However, when the system only bears a Gaussian variable or pure Brownian noise, the Hurst exponent [26] can be applied as the analytical standard for distinguishing the Gaussian variable, white Gaussian noise and long-range dependence automatically [27]. To this end, the application fields of the SSA method have been significantly expanded. In this paper, we apply a filter based on the basic SSA method and the Hurst exponent to reduce the femtosecond laser absolute distance measurement uncertainty automatically in the case of a wide range of absolute distance changes with a series of speeds. Both numerical simulation and an absolute ranging experiment based on a pair of free-running femtosecond lasers show that the basic SSA filter is capable of retrieving a moving target trajectory from the measurement data that is degraded by noise. The distance measurement uncertainty after filtering reaches 110 nm at a 200-Hz update rate. This method shows great potential to extend the femtosecond laser ranging technique to a wide range of LiDAR and dynamic measurement applications.

Singular Spectrum Analysis
The strategy of basic SSA is to decompose the observed time series into a sum of independent components. The procedure of basic SSA is in two isolated steps: decomposition and reconstruction. Each consists of another two steps separately. The decomposition step is composed of time lagged embedding and singular value decomposition (SVD), while the reconstruction step is composed of diagonal averaging and grouping. A flowchart of basic SSA, consisting of the sub-steps of decomposition and reconstruction, is shown in Figure 1.
In the first step, an initial dataset Y N , consisting of N samples with an equal time interval is projected in a vector space via delayed coordinate mapping over a lagged window length L, which is an integer satisfying 2 ≤ L ≤ N/2. The variable L is called the embedding dimension of the lagged vector space and should be adequately large. Through delayed coordinate mapping, we take L consecutive samples to create a time delay vector, x i , (2) Considering the relationship between the number of samples N and the lagging window length L, we get N − L + 1 lagged vectors. Next, we get the trajectory matrix X by connecting all the lagged vectors. Obviously, since the elements on the vice diagonal lines are equal, X is a Hankel matrix, In the second step, SVD is applied to decompose the trajectory matrix X into the form as follows: where U and V are unitary matrices. Their corresponding column vectors u i and v i are the left-singular vectors and right-singular vectors of X, which are orthogonal eigenvectors of matrices X X T and X T X, respectively. The vectors u i and v i are defined as the factor empirical orthogonal functions (EOFs) and principal components (PCs) of X, respectively [28]. Σ is a rectangular diagonal matrix whose elements on the diagonal line are in descending order. The elements on the diagonal line σ are singular values of X, which are the positive square roots of the eigenvalues of matrix X X T . Due to the dimension of decomposition d being equal to the rank of the matrix X, X can be written as: where X i is one of the decomposition matrices of X. For every decomposition X i , the corresponding pairs (u i ,σ i ,v i ) are sorted by descending order of singular values. Note that the square of each singular value can describe the contribution of the corresponding component to the initial signal; thus, the ratio can be considered as the contribution ratio of the k-th component X k . After SVD operation of X in the diagonal average step, the decomposition matrices X i of the trajectory matrix X should be converted back to a set of reconstructed time series Y (i) by the diagonal average, an inverse operation of time lagged embedding. In order to obtain k-th element y (i) k of the reconstructed time series Y (i) from the decomposition matrix X i , we need to take an average value of the elements on the k-th vice diagonal line. Then, any element of a reconstructed signal y (i) k can be described by the following formula: where L represents the embedding dimension, N is the length of the original signal and a i,j is the element of the decomposition matrix X i on the i-th row, j-th column.
Finally, in the grouping step, reconstructed time series Y (i) are divided into several classes and summed according to the category. In general, the SSA method can be classified as trend, oscillation and noise. In practice, it is very common that the leading component corresponds to the trend. In our application, the distance signal can be grouped into the classes of signal or noise. Hurst exponent H can be used as the grouping basis of basic SSA, allowing signal extraction from white noise automatically. The case of H > 0.5 indicates that the reconstructed time series is a long-term memory series, a part of the reconstructed signal; while the case of H ≤ 0.5 implies that the reconstructed time series is white noise that needs to be eliminated. In the calculation, we used the aggregated variance method to calculate the Hurst exponent [29].

Absolute Ranging Based on a Pair of Femtosecond Lasers
The principle of absolute distance measurement based on a pair of femtosecond lasers is shown in Figure 2. Two femtosecond lasers with similar repetition frequencies ( f r and f r − ∆ f r ) are used as the signal laser and the local oscillator (LO), respectively. The signal laser is directed at the reference and the target so that the target TOF is simply determined by the temporal separation of the retroreflected laser pulses from the reference and the target. LO is used to sample the retroreflected signal laser pulse envelopes nonlinearly via sum-frequency generation. Owning to the slightly different repetition rate, the sampling pulses will walk through the signal pulses and equivalently "stretch" the time frame by N = f r /∆ f r -times. The so-called nonlinear ASOPS process [12] makes the timing of ultrafast signal laser pulses accessible to standard electronics. In order to retrieve the target distance, we only need to measure the TOF of the target and the reference in the "stretched" time scale, as well as f r and ∆ f r . In practice, we do not measure the difference of repetition frequencies ∆ f r directly; rather, we use the following equation, where t re f 1 , t re f 2 are the measured timing of two successive temporally-broadened reference pulses. As a result, one can calculate the absolute distance L by: where c is light velocity in a vacuum, n g is the group refractive index in air, t TOF is the real TOF of the femtosecond laser, while t broaden and t tar represent the target-reference TOF and the timing of target pulse in the "stretched" time scale, respectively.

Simulation
A numerical simulation is conducted in the first place in order to testify the effectiveness of basic the SSA filter on the ASOPS-based femtosecond laser ranging system. To this end, we constructed an arbitrary target trajectory consisting of a motionless state, an accelerated motion and a uniform motion by computer. The motion parameters have been listed in Table 1. Here, a positive velocity is defined as an increase of the distance, while a negative velocity indicates a decrease of distance. The initial distance value of the simulation is set to 200,000 µm. Suppose the femtosecond laser ranging ideally reproduces the target trajectory; the retrieved target trajectory is shown as the red dashed line in Figure 3a. Here, the time span is 1.495 s, and the measurement sampling rate is 200 Hz (typical for ASOPS-based femtosecond laser ranging), containing 300 distance data points. The inset in Figure 3a shows the enlarged view for State ii, i.e., the transit from State i to State iii. In order to simulate the distance measurement uncertainty, a sequence of 300 Gaussian distributed random numbers with a standard deviation of 9.98 µm have been added to the ideal distance array. The simulated target trajectory considering measurement uncertainty is shown as the black curve in Figure 3a.
In the following, the basic SSA filter has been applied to the distance sequence with noise so as to retrieve the real target trajectory. As described in the Methods Section, only the lagged window length matters for effective implementation of an SSA filter. In order to choose an appropriate window length, we process the simulated distance signal under nine different window lengths, which ranges from 10 to 100 pts, corresponding to a time span from 0.05 to 0.5 s. The filtered distance sequence in each simulation is compared with the ideal target trajectory, and the effect of the SSA filter is evaluated by their correlation coefficient and the RMS residual that is simply the standard deviation of the errors obtained by subtracting the filtered distance sequence from the ideal one. The obtained correlation coefficient, as well as the RMS residual at different SSA filter window lengths are shown in Figure 3b. Obviously, the window length of 50 pts (0.25 s) provides the highest correlation coefficient, as well as the minimum RMS residual of 1.1 µm. The retrieved target distance sequence obtained by SSA filtering at a 50-pts window length is plotted as the blue dotted line in Figure 3a. The result explicitly demonstrates the effectiveness of the basic SSA filter for retrieving the target trajectory from Gaussian white noise. Note that a priori knowledge on the target kinetic behavior is not required.

Experimental Setup
An absolute ranging experiment based on a pair of free-running femtosecond lasers is conducted by using the basic SSA filter for ranging precision improvement in practice. The experimental setup is shown in Figure 4. The ranging system uses a pair of free running carbon nanotube (CNT) mode-locked Er-fiber lasers working at 1550 nm, serving as the signal laser and LO, respectively. The average output powers after the Er-doped fiber amplifier (EDFA) are 10 mW and 13 mW, respectively. Their repetition rates are~74 MHz, having a differential repetition rate of~2 kHz. The laser pulses emitted by the signal laser are separated by a polarization beam splitter (PBS) and are sent to the reference mirror and target retro-reflector, respectively. The target is a retro-reflector mounted on a motorized linear translator that allows one-dimensional movement at various speed. The retro-reflected pulses from the target and the reference are beam-combined and are nonlinearly sampled by the local oscillator with an~2-kHz repetition rate offset in a periodically-poled potassium titanyl phosphate (PPKTP)-based optical cross-correlator, thus accomplishing nonlinear ASOPS. The cross-correlation traces are detected by a high speed avalanche photo-detector, and the output voltage signals are digitized and stored by a 14-bit 100-MHz digitizer (PXIe-5122, National Instrument, Austin, TX, USA). The repetition rate of the signal laser is measured by a commercial frequency counter (53220A, Agilent, Santa Clara, CA, USA) referenced to a rubidium atomic clock. Detailed laser ranger design, as well as data gathering and processing can be found in our recent publications [16,18].

Experimental Data Processing
An arbitrary target movement is generated by programming of the motorized translator. Specifically, the target maintains its initial position at about 214,777.80 µm in 5.1 s. In the following, the target moves away from the range finder at a speed of 0.34 mm/s in 2.3 s. After a 1.2-s pause, the target moves towards the same direction at a speed of 0.75 mm/s in 2.37 s. Then, the target moves towards the opposite direction at a speed of 0.5 mm/s for 0.47 s. The target remains stationary at the last 0.615 s, thus completing the movement. The detailed movements of the target are listed in Table 2. The measured target movement trajectory by the absolute range finder is shown in the black lines of Figure 5. The update rate of the absolute distance measurement is 200 Hz. Note that the ASOPS-based absolute range finder used here has an inherent measurement update rate of 2 kHz determined by the offset repetition rate. However, we intentionally lower the update rate to 200 Hz by taking only one datum out of every 10 successive measurements considering the limited storage depth of the digitizer (PXIe-5122), which only allows storage of 1-s data at a 2-kHz update rate.
The basic SSA filter is used to process the experimental data, and the procedure is identical to that in the simulation. The result after the basic SSA filter with a 0.25-s window length is shown in red lines in Figure 5a. The choice of a 0.25-s window length is on the basis that the motion behavior is similar to the simulation in the above section. The enlarged view of the transits are shown in Figure 5b-f. It can be seen that the basic SSA filter effectively improves the measurement signal-to-noise ratio and retrieved the trajectory of the moving target. In comparison, we also process data by moving averaging with a 64-pts sliding box, shown as green dots in Figure 5. Despite the same level of curve smoothness, the retrieved trace by moving average apparently deviates from the real trajectory, particularly at non-uniform motion stages.
It is worth mentioning that the basic SSA filter is able to extract the noise sources while filtering a distance signal. Figure 6 illustrates the noise sequence separated by the basic SSA filter and the corresponding power spectral density (PSD). The PSD indicates that the Gaussian white noise process dominates during ASOPS-based absolute laser ranging. The noise STD is 16.75 µm RMS. In addition, the PSD of the measured distance signal and reconstructed signal by the basic SSA filter are shown in Figure 7a, showing that the noises in the distance signal are well filtered out, while the spectral characteristics of the target movement at baseband are well preserved.  Allan deviation analysis is utilized to evaluate the effect of the basic SSA filter on the improvement of measurement uncertainty. To this end, we calculated the Allan deviation of a section of measured distance signal (from 0 to 1.25 s in Figure 5a when the target is stationary) after basic SSA filtering, and the result is shown in the red curve in Figure 7b. The minimum absolute measurement uncertainty after basic SSA filtering is 110 nm, obtained at the fundamental 200-Hz update rate (5-ms observation time). In comparison, the Allan deviation of the same section of measured distance signal without basic SSA filtering in shown as the black curve in Figure 7b. At the fundamental 200-Hz update rate (5-ms observation time), the measurement uncertainty is 18 µm, which is 163-times higher than the number obtained after basic SSA filtering. The measurement uncertainty deceases with observation time (multiple averaging) and reaches 2.4 µm at a 160-ms observation time. This number is still 21-times larger than the 110-nm measurement uncertainty, which was obtained without any loss of update rate based on basic SSA filtering. With even more times of averaging, the measurement uncertainty may approach that obtained by basic SSA filtering. At an expense, the significantly lowered update rate will not be suitable for capturing the motions of the moving targets.

Conclusions
To summarize, we introduce a data-based filter, namely the basic SSA filter, to femtosecond laser absolute range finders for significant reduction of measurement uncertainty in the conditions of irregular moving targets. We verify the effectiveness of the basic SSA filter both in numerical simulation and in an absolute ranging experiment based on a pair of free-running femtosecond lasers. The simulation shows that the basic SSA filter is capable of reproducing a target trajectory that is significantly degraded by a white Gaussian noise. The priori knowledge on the dynamic behavior and the noise characteristics is not required. In the experiment, we successfully separate the signal and noise components by applying the basic SSA filter. We show that Gaussian white noise is the dominant noise process in an absolute range finder based on two free-running femtosecond lasers. By eliminating the noise components, the range finder shows 110-nm measurement uncertainty at a 200-Hz update rate for an irregularly-moving state. This method can be potentially used in a number of ASOPS-based high precision measurement fields. Firstly, the capability of separating various noise and signal components allows basic research on the noise characters in the femtosecond laser ranging system and shines light on the improved femtosecond laser range finder design. Secondly, this method enables a variety of industrial dynamic and vibrational measurements based on femtosecond laser range finders. Finally, this method can also be used for LiDAR and remote sensing applications, where the basic SSA filter can be used for target identification from point clouds obtained by scanning femtosecond LiDAR.

Conflicts of Interest:
The authors declare no conflict of interest.