# Co-Processing Parallel Computation for Distributed Optical Fiber Vibration Sensing
## 4 Experimental Results and Discussions
### FFT Analysis Processing  
Figure 10 shows the cooperation process of the CPU and GPU to realize the spectrum analysis, which could be important for the position location and frequency component analysis of the external disturbance. Firstly, the original \(\Phi\)-OTDR data acquired by FPGA is transferred into the RAM of the host computer via the USB bus, and can then be processed in the form of the \(\Phi\)-OTDR data matrix (\(N^{\star}M\)), with \(N\) time traces and \(M\) sampling points in each time trace. Secondly, the CPU copies time traces into graphics memory, which can be read and modified by the GPU through the memory copy function. The GPU then allocates simultaneously numerous threads to achieve the one-dimensional FFT computation. In addition, \(n\) represents the number of time traces in a thread on the GPU. Finally, the processed result is sent back to the host computer via the PCI-E bus and is displayed on the screen. Thus, the frequency spectrum of temporal signal can be obtained, and it is beneficial to identify the external vibration position. In conclusion, this parallel operating method largely shortens the data processing time compared with the traditional serial computing method based on the CPU. The process of FFT analysis based on the GPU mainly contains two parts, which are multiple curves parallel computation and the butterfly operation, which relies on the butterfly-shaped arithmetic unit. It can simultaneously provide the needed operands for many arithmetic units, owing to the advantages of its parallel structure and in-place calculation. Thus, this parallel operation structure improves the speed of FFT calculation and shortens the time noticeably. Because of the butterfly operation of FFT, when the number of calculation points is \(M\) in one curve, the algorithm complexity is O(\(M\)/2 \(\times\) log2(\(M\))). The process of the first part is shown in Figure 11, the RBS data is generated and stored in the CPU and transmitted into the GPU. Then, the \(N\) threads are allocated for the time traces, which means each time curve is processed by one thread. It is necessary to verify whether the length of the time  
Figure 10: Cooperative data processing with GPU and CPU. trace (\(M\)) is a power of two, owing to the existence of the butterfly operation. If the length is not a power of two, the data needs to be complemented by adding a zero value at the end of the data string. Each curve is processed by these steps on the GPU at the same time. The second part is the butterfly operation of each time trace by using the butterfly arithmetic unit. This FFT analysis is called decimation in time (DIT) FFT. The process is shown in the equation as follows:  
\[X(k)=DFT[x(n)]=\sum\limits_{n=0}^{M-1}x(n)W_{M}^{\text{kn}}\qquad 0\leq k\leq M-1 \tag{5}\]  
Equation (5) represents the computation process of conventional discrete Fourier transform (DFT). In Equation (5), \(W_{M}^{\text{kn}}\) is the rotation factor which has features of symmetry and periodicity. This characteristic can decrease the computation time of FFT. The signal of the temporal domain (\(x(n)\)) could be decomposed into even and odd sequences, so the FFT of the even data sequence is called as \(X_{1}(k)\) and the FFT of the odd data sequence is called as \(X_{2}(k)\). It is shown in the equation as follows:  
\[ \begin{array}{l}X_{1}(k)=\sum\limits_{n=0}^{M_{1}/2-1}x(2n)W_{\frac{M}{2}}^ {\text{nk}}\\ X_{2}(k)=\sum\limits_{n=0}^{M_{1}/2-1}x(2n+1)W_{\frac{M}{2}}^{\text{nk}}\\ \end{array}  \tag{6}\]  
Therefore, Equation (5) could be developed as Equation (7). \[ \begin{array}{l}X_{thread}(k)=X_{1}(k)+W_{M}^{k}X_{2}(k)\\ X_{thread}(k+\
icefrac{{M}}{{2}})=X_{1}(k)-W_{M}^{k}X_{2}(k)\qquad 0\leq k\leq \
icefrac{{M}}{{2}}-1\\ \end{array}  \tag{7}\]  
Equation (7) expresses the principle of the butterfly operation in one thread. \(X_{thread}(k)\) indicates the frequency value of the first half of all points and \(X_{thread}\) (\(k+M/2\)) indicates the frequency value of the second half of all points.