High Resolution Coherent Doppler Wind Lidar Incorporating Phase-Shift Keying

A coherent Doppler wind lidar (CDWL) with spatial/temporal resolution of meter and sub-second scale is reported. Phase-shift keying (PSK) is introduced into the pulse pair to enhance the measurement efficiency, and the pulse shaping and perfect binary PSK modulation are both realized by a single Mach-Zehnder modulator (MZM). The identical intensity envelope of the paired pulses reduces the influence of amplifier's gain profile and provides a monitoring parameter for bias control, which improves the system stability in practice. Theoretical analysis and numerical simulation are conducted to illustrate the principle of the proposed method. In experiments, the bias drift of the MZM is monitored in real time. Continuous radial wind profile measurement over 800 m is demonstrated with spatial and temporal resolution of 3 m and 0.1 s, respectively.


I. INTRODUCTION
A S A powerful remote sensing technique, Doppler wind lidar has been widely used in scientific researches and engineering applications, such as gravity waves monitoring [1], [2], boundary layer detection [3], [4], wind energy exploitation [5], [6], weather forecast [7], and air pollution analysis [8], [9]. Due to the inherent properties of coherent detection, a coherent Doppler wind lidar (CDWL) has the advantage of acquiring the spectral characteristics of the atmospheric backscattered signals, which extends its scope to spectrum analysis based applications like precipitation observation [10], [11], cloud phase identification [12] and so on.
Despite all these successful examples, the performance of CDWLs still needs to be further improved in some aspects.
In areas like aerodynamics design and aviation safety, the CDWL has been used to detect the aircraft wake vortices [13], [14]. Since the wingspan of the aircraft is on the order of ten meters, meterscale detection of the wind field raises higher requirements for the CDWL. With higher spatial/temporal resolution, more detailed data can be got in the field experiment, thus improving the estimation accuracy of the wake vortex parameters [15], [16].
Nevertheless, it's quite challenging to enhance the spatial resolution of a conventional pulsed CDWL to meter scale, which depends on the pulse duration [17]. When the pulse duration is shortened to tens of nanoseconds, the estimation accuracy of the Doppler frequency shift will be greatly degraded due to the spectrum broadening of the backscattered signals [18]. Moreover, the limited peak power of fiber lasers further aggravates the carrier-to-noise ratio (CNR) of the backscattered signal under short pulse duration [19]. To mitigate these negative effects, complicated hardware and algorithms have to be adopted under conventional pulsed method, realizing the spatial resolutions from 3 m to 60 m [20], [21], [22]. Recently, several works are reported to enhance the spatial resolution without shortening the pulse duration. In the pseudo-random phase coding (PRPC) method [23], a 4.5-m spatial resolution is realized by decoupling the spatial resolution and pulse duration with the pseudo-random phase modulation of the lightwave. However, to achieve the fast phase transition between bits, a very high modulation bandwidth is required. Another attempt is the differential correlation pair (DCP) method [24], where the spatial resolution is determined by the duration of the short probing pulse and the window function. Under this scheme, the spatial resolution is enhanced to 3.3 m. However, due to the different intensity envelopes of the paired pulses, the pulse shape before amplification needs to be particularly designed to adapt to the gain characteristics of the amplifier, so as to ensure the cancellation of the common pulses.
In this work, we propose and demonstrate a novel CDWL incorporating the DCP technique and phase-shift keying (PSK) modulation. Inspired by the phase-shift pulse technique in distributed fiber sensors [25], the proposed PSK-DCP method adopts probing pulses with phase difference of π, thus enhancing the measurement efficiency compared with the original DCP scheme. Since the phase modulation does not change the intensity envelope, the influence of the amplifier's gain profile can be greatly reduced when the pulses are equally spaced, which makes the system more robust to environmental disturbance. To facilitate engineering applications, the PSK modulation together with the pulse shaping is realized by a single Mach-Zehnder This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ modulator (MZM) [26]. Meanwhile, a real-time bias monitoring scheme for the MZM is proposed.

A. CDWL With Phase-Shift Keying
In the DCP method, the transmitted pulses' amplitude envelopes are composed of a long common pulse A L (t) and a followed short probing pulse A S (t)with truncation lengths of T L and T S , respectively [24]. When binary PSK is employed, the probing pulse in the even pulse is phase shifted by π. It is noted that the phase shift is defined by the phase difference between the common pulse and probing pulse in the even pulse, rather than the odd and even pulses. As illustrated in Fig. 1, the amplitude envelopes of the ith (i = 1 or 2) pulse in the pair can be expressed by: According to the classical layered model of atmospheric scattering [27], the backscattered signal can be modeled as the superposition of the random atmospheric reflections from different slices. Provided the pulse duration is less than the relaxation time of the scattering particles, only little phase distortion would be introduced in the scattering [27]. After heterodyne detection, the received signal is: where To retrieve the wind velocity at distance z 0 = ct 0 /2, a long window function W L (t − t 0 ) and a short window function W S (t − t 0 ) are used to truncate the received IF signal. The truncation lengths of the two window functions are the same as the common pulse and probing pulse, respectively. The windowed IF signals are denoted by r where the operation refers to averaging over pulses and refers to complex conjugate. Due to the random nature of aerosol scatters, the auto-correlation of κ(z) is β(z)δ(z − z ), with β(z) being the mean reflectivity at distance z and δ(z − z ) being the Dirac delta function. By defining z = z 0 + δz and t = t − t 0 with z 0 = ct 0 /2 and δz = cδt/2 , the R (i) (τ, t 0 )can be simplified as: in the range of "L" and "S". The contribution of F ab can be explicitly analyzed with the help of Fig. 1. Firstly, F SL is zero for all the cases. When δt < 0( Fig. 1 which makes their product always zero. Secondly, since the F LL and F SS terms exist in both R (1) (τ, t 0 ) and R (2) (τ, t 0 ), they will be eliminated in a differential operation. Thus, only the F LS term remains in the differential cross-correlation ( Fig. 1(c), (f)), leading to: To further understand the characteristics of the F LS , a numerical simulation is carried out. The window functions and the pulses' amplitude envelopes are set to be truncated Gaussian functions, where the half width at exp(−0.5) and truncated length of both W L and A L are 52 ns and 200 ns, while those parameters of W S and A S are 8 ns and 30 ns. The parameters are chosen for imitating the behavior of the experiment system described in Section III.
The contour map of the F LS is shown in Fig. 2(a) with its sections cut at τ = 115 ns and δt = 0shown in Fig. 2(b) and (c). The variations of the F LS in the τ and δt axes indicate the amplitude changes of the signal scattered by a specific spatial slice and the different contributions of the signals originate from different spatial slices, respectively. It can be observed that the non-zero part along the δt axis is determined by the duration of the short probing pulse while the non-zero part along the τ axis is mainly determined by the long common pulse. Note that δt = 2δz/c , the (5) indicates that the spatial resolution and the effective duration of the R d (τ, t 0 ) for a specific t 0 depend on the probing pulse and common pulse, respectively. Thus, by substituting the parameters of the probing pulse and the short window function, the recognized formula of the pulsed CDWL's spatial resolution can still hold, specifically [17]: where σ p and σ w are the half width at exp(-0.5) of the probing pulse's amplitude envelope and the Gaussian window function, respectively. More detailed information about the simulation can be found in the Appendix.

B. PSK Modulation and Bias Monitoring Based on MZM
Typically, the intensity and phase modulation can be realized by an MZM and a phase modulator (PM), respectively. However, an additional phase modulator will increase the complexity and decrease the robustness of the system, as the half-wave voltages of both the MZM and PM need to be calibrated. To avoid that, we employ a single MZM biased at the null point to realize pulse shaping and binary PSK modulation here. Meanwhile, the bias of the MZM can be monitored in real time from the modulated intensity waveform. The transfer function of an MZM can be expressed as [26]: where V m (t) is the modulation voltage, V DC is the DC offset, V πRF and V πDC are the half-wave voltages of the RF and DC electrodes, and φ 0 is the common phase shift introduced by the modulator. For simplicity, the φ 0 can be set to π without losing generality. When the modulator is biased at the null point, V DC = V πDC , the transfer function can be simplified as: As the sine function is an odd function, the sign of the output lightwave will follow the sign of the modulation voltage. In other words, when the modulation voltage goes through a zero-crossing, a phase shift of π will be introduced into the modulated light field. Thus, the PSK modulation of the probing pulse can be realized by applying modulation signals with the same magnitude but the opposite sign. When the DC offset drift happens, We mark the modulation voltage at the peak of the probing pulse as V M and the intensity of the lightwave entering the MZM as I 0 . As the intensity transfer function is the square of the lightwave transfer function's magnitude in (7), the output intensity of the paired pulses at the original peak position of the probing pulses can be expressed as: When |V b /V πDC | + |V M /V πRF | ≤ 1, I p1 and I p2 represent the peak intensity of the probing pulses, as illustrated in Fig. 3(a) and (b).
One can define the normalized difference (ND) of the paired pulses and obtain its small signal model as: 10) where k d is the detection sensitivity of bias drift. Fig. 3(c) shows the ND in simulation with V M , V πRF and V πDC of 5V, 6V and 6V, respectively. The good agreement of the simulated curve with its linear fit indicates that the ND can act as a linear sensor of the bias voltage. With a feedback control, the bias voltage can be maintained around zero.

A. System Layout
The system layout of the proposed PSK-DCP lidar is shown in Fig. 4(a). A continuous-wave laser (CW, Amonics ALiDAR-150-Seed) is split into two parts, where the minor portion (5%) serves as the local oscillator (LO). The continuous seed is frequency shifted by 80 MHz and chopped into pulses with repetition rate of 25 kHz by an acousto-optic modulator (AOM). The pulse shaping and PSK modulation are both realized by an   Fig. 4(b) illustrates the equally spaced paired pulse train.
A small portion (1%) of the modulated signal is detected by a photodetector (PD) to monitor the bias voltage of the MZM. In experiments, the MZM is biased at the null point. An erbium-doped fiber amplifier (EDFA, Amonics ALiDAR-150-AMP) is employed to amplify the paired pulses. In the receiver, the backscattered signal is mixed with the LO by a 50/50 beam splitter and down converted into electrical signal in IF band by beating on a balanced photodetector (BPD). After digitized by a 14-bit analog-to-digital converter (ADC), the raw data are stored in a computer for offline processing with the algorithm depicted in Section II. The AOM, MZM, and ADC are synchronized on a common reference provided by the arbitrary waveform generator (AWG, Keysight 33622A). The key system parameters are summarized in Table I.

B. Instrument Characterization
To get the detection sensitivity k d of the bias drift under practical conditions, the normalized differences of probing pulses  are measured at different bias voltages. The DC voltage on the MZM is changed from 5V to 6V and the ND is calculated from the observed intensity shapes. To distinguish the odd and even pulses, an extra trigger at the same repetition rate of the pulse pair is used. As shown in Fig. 5(a), the k d is extracted from the slope of the linear fit of the measured data, with the R-squared parameter of this linear regression being 99.91%. Fig. 5(b) shows a continuous monitoring of the bias drift for 16 s. The refresh time interval of the bias detection is set to be 20 ms, averaging on 250 pulse pairs. Intentional disturbances are applied in a step of 0.1 V, which is well followed by the bias detector.
When the MZM is biased at the null point, the modulated lightwave is observed through homodyne detection with a 180°o ptical hybrid. As an example, Fig. 6(a) shows the electric field's amplitude envelope of the paired pulses before amplification, where the truncated lengths of A L and A S are 200 ns and 30 ns, respectively. The half widths at exp(−0.5) of the common pulse (σ c ) and probing pulse (σ p ) are 50 ns and 8 ns. The peak magnitude of the probing pulse is designed to be slightly higher than the common pulse to reduce their intensity difference after amplification [28]. As shown in Fig. 6(b), the intensity envelopes of the amplified odd and even pulses are quite close. The full width at half-maximum (FWHM) pulse durations of the common pulse (τ c ) and probing pulse (τ p ) are measured to be 85 ns and 13 ns. As to the receiver, since the power of the LO is as high as 1.5 mW, it is operating at the shot noise limit [29].

C. Wind Profiling Results
During the experiments, the system is placed on the 9th floor of the laboratory building (32.20°N, 118.72°E) with the telescope pointing to the south at an elevation angle of 15°. In the data processing, window functions matched to the amplitude envelopes of the common and probing pulses are used to achieve optimum receiving [30]. As the half widths at exp(-0.5) of the probing pulse's amplitude envelope and the corresponding window function are both 8 ns, a spatial resolution of 3 m is realized according to (6). It is noted that the (6) has been verified in the experiment [24]. The temporal resolution is set as 0.1 s in this work, averaging the backscattered spectra of 1.25k pulse pairs. For each spectrum, the fast Fourier transform (FFT) length is extended to 4096 through zero padding.
As a reference, the atmospheric echoes of the odd pulses are regarded as the backscattered signals got by a conventional CDWL with low spatial resolution. Since the conventional pulsed CDWL accumulates the backscattered spectrum of each single pulse (SP) rather than pair, it would be called as SP Ref for brevity. The consistency with SP Ref is a necessary condition for the validity of the proposed lidar's measurement. As it's difficult to directly get the uncertainties of a wind lidar in such a high resolution due to the lack of calibrated references [28], the standard deviation (SD) of the velocity would be preliminarily estimated by the Cramer-Rao lower bound (CRLB) [31], [32]. While the center frequency is estimated by Gaussian fitting with trust-region algorithm here, the CRLB provides a lower bound of the estimation uncertainty [28].
In a first wind profiling experiment, measurement results are plotted in Fig. 7, where the blue solid line and black dash line stand for the PSK-DCP and SP Ref, respectively.
The velocity SD of the PSK-DCP lidar is derived from the narrowband carrier-to-noise ratio (CNR n ) [31], as shown by the    Fig. 7(a). The CNR n is defined in the spectral domain as the signal power divided by the noise power under the 1/e full width of the signal peak. The trend of CNR n versus distance is related to many factors, including the system configuration (like the focal length) and the weather condition (like the visibility). It can be observed that the SD of the velocity doesn't exceed 0.15 m/s within 800 m.
In a continuous measurement for 0.3 s, the wind velocity profiles of the two methods have the same trend versus distance, as shown in Fig. 7(b)-(d). The profile got by the high-spatialresolution PSK-DCP lidar changes continuously around that of the SP Ref. There are some differences between the two sets of profiles at some distances (like 420 m and 618 m), which can be attributed to the smoothing effect caused by the low spatial resolution of the SP Ref. The backscattered spectra at these distances are shown in Fig. 8. One can observe that the PSK-DCP method basically avoids the spectrum broadening while achieving high spatial resolution. Though the CNR n of the PSK-DCP is lower compared with the SP Ref due to the high spatial resolution, the signal peaks can still be discriminated from the noise effectively.
A second experiment was carried out on 16 March 2022. While a cold wave went through the East China, the temperature dropped by 20°C (from 29°C to 9°C) in one day. During the severe convective weather caused by the strong cold air, the wind field may change greatly in a relatively small temporal and spatial scale. As an initial prototype, the system cannot store long-term measurement results, limited by the huge amount of raw data. Therefore, a continuous observation for 8 s was conducted at 19:30 local time. As shown in Fig. 9, the evolution of the wind field is profiled with a spatial/temporal resolution of 3 m/0.1 s. In the experiment, continuous structures of the wind field can be observed. For instance, a band with positive wind velocity exists for over 6 s around 650 m away. One possible cause for the layering can be the gradient distribution of air temperature near the ground related to the sudden drop of temperature. This continuity conforms the time-varying characteristics of the wind field, which further verifies the reliability of the measurement.

IV. DISCUSSION
Different from the previous work where only one probing pulse is included in the pulse pair [24], the PSK-DCP method employs two probing pulses with π-phase shift in pair. Therefore, the total energy of the probing pulses doubles, as well as the CNR n [19]. Another benefit from the PSK modulation is on the waveform control. When the EDFA is uniformly pumped, pulses with different intensity envelope will be amplified with different gains. To neutralize the time-dependent gain profile of the amplifier, the modulating signals need to be particularly designed to minimize the difference between the amplified common parts in one pulse pair [24]. Moreover, the design varies on different environment and amplifiers, which increases the difficulty of system deployment. In the PSK-DCP method, benefiting from the identical intensity envelope of the paired pulses ( Fig. 6(b)), the problem above can be avoided. By equally spacing the pulses with identical intensity envelope, the system is less sensitive to the gain profile of the EDFA, thus reducing the difficulty of waveform control and improving the robustness against environmental disturbance.
In the proposed method, while the duration of the probing pulse is directly related to the spatial resolution, the duration of the common pulse is flexible. Although a longer common pulse can reduce the spectral width by extending the effective duration of R d as explained in Section II-A, some negative effects will also be introduced. Firstly, the long specular reflection from the output optics in a monostatic configuration will increase the blind area [28]. Secondly, a longer common pulse may reduce the CNR n by increasing the variance of the R d . It is noted that the expectation value of the cross-correlation of backscattered signals from different distances is zero, but the variance always exists. To optimize the system performance, trade-offs need to be made.

V. CONCLUSION
A coherent Doppler wind lidar based on DCP technique incorporating PSK modulation is proposed and experimentally demonstrated. The principle of the PSK-DCP method is explained in detail through theoretical analysis and numerical simulation. In experiments, the binary PSK modulation together with the pulse shaping is realized by an MZM with real-time bias monitoring. At laser peak power of 400 W, continuous radial wind profiling over 800 m is demonstrated with the spatial and temporal resolution of 3 m and 0.1 s, respectively.
Considering that the high-speed signal processing has been realized by a field-programmable-gate-array (FPGA) in real time [33], the limitation of data storage on continuous observation duration could be greatly alleviated. To further demonstrate the role of the high resolution CDWL in aerodynamic structure observation, a measurement of the aircraft wake vortices is under preparing. The retrievals based on spectrum analysis [11] would be used to deal with the spectral deformation caused by the rapid dynamic process.

APPENDIX: NUMERICAL SIMULATION OF F LS
The truncated Gaussian function used in the simulation is formed by truncating a Gaussian function at the point where the intensity falls to 10% of its peak and shifting the DC level to the truncated point to avoid discontinuity. A truncated Gaussian function with truncated length of 200 ns is used for the amplitude envelope of the common pulse A L (t) and the long window function W L (t). Another truncated Gaussian function with truncated length of 30 ns is used for the amplitude envelope of the probing pulse A S (t) and the short window function W S (t). The F LS (δz, τ ) is defined as: (11) The variables are meshed by the parameters shown in Table II, and the values of the F LS (δz, τ ) at the 2D grid are calculated by numerical integral. The result has been shown in Fig. 2.
To clarify the spatial resolution of the DCP method, the contribution of the signals from different δz to the R d (τ, t 0 ) is analyzed base on (5). Assuming that the β(z 0 +δz) keeps unchanged when δz is in the range of a few meters, the power of each complex exponential function exp{−j2π[f IF + f d (z 0 + δz)]τ } is determined by the square of its amplitude envelope, which is F 2 LS (δz, τ ). Thus, the energy proportion (p E ) of the backscattered signal from a cell within range z c can be got by: As shown in Fig. 10, 87% energy of the differential cross correlation R d (τ, t 0 ) is contributed by signals from a range cell of 3 m. The result is consistent with the spatial resolution calculated by (6).