Parameter-retrieval of dry-air wake vortices with a scanning Doppler Lidar

Hang Gao; Jianbing Li; P. W. Chan; K. K. Hon; Xuesong Wang

doi:10.1364/OE.26.016377

1. Introduction

Wake turbulence is an inevitable by-product of lift for a flying aircraft. It is generally composed of two strongly counter-rotating trailing vortices and is regarded as one of the most dangerous hazards in aviation [1]. The wake vortices generated by a large aircraft may cause a following aircraft to roll out of control, particularly during the departing and landing phases. To avoid such accidents, ICAO (International Civil Aviation Organization) released a series of regulations for the minimal temporal–and spatial–separations between flights in the 1970s. These regulations can ensure the safety of flights most of the time, but they are conservative and believed to have limited the capacity of airports to a certain degree [2, 3]. With the rapid increase of air traffic, the detection and parameter-retrieval of wake vortices are greatly demanded by ATM (air traffic management) to adjust the departure and landing of flights in real-time for higher capacity. Programs such as SESAR (Single European Sky Air Traffic Management Research) [4] in Europe and NextGen (Next Generation Air Transportation System) [5, 6] in the US have also launched many projects on this issue.

The existing wake vortex detection sensors include Radar and Lidar. They are applicable for different weather conditions due to their different scattering mechanisms. Under wet weather conditions, Radar is advisable because the attenuation of Radar waves through precipitation is relatively weak [7–10]. But under dry-air condition, the Radar echo (caused by the fluctuation of reflectivity index) is extremely weak and Lidar becomes more favorable due to its sensitivity of the aerosols involved with the wake. In recent decades, Lidar detection of wake vortices has been well investigated [11, 12], and the representative parameter-retrieval methods can be classified into two categories.

The first group of parameter-retrieval algorithms by Lidar is based on velocity envelope extraction, and the basic procedures include: 1) locate vortex cores by velocity envelopes, 2) obtain the circulations from a given wake velocity model. For example, Rahm et al. [13, 14] obtained the velocity envelopes and vortex-core positions by applying a fixed and SNR-related (Signal-to-Noise Ratio) threshold to the Doppler spectra, while Wassaf et al. [15] proposed using a floating threshold instead; the circulations were then calculated according to the existing wake vortex velocity profile models. Ramsey et al. [16] estimated circulations and core positions with the method proposed by Rahm et al., but the process was specified to the region of interest inhabited by the vortex cores. For this group of algorithms, locating the vortex cores is crucial since it can greatly affect their robustness and accuracy.

The second group of algorithms is based on template matching. In these algorithms, mathematical models related to the vortices’ parameters were established, and the target parameters were obtained by fitting the mathematic models with the measured data. WakeMod algorithm and Lockheed Martin algorithm [17–19] both estimated vortex-core positions and circulations based on spectrum-matched filters and maximum likelihood (ML) estimators. Smalikho et al. [20, 21] evaluated the vortex-core positions from the measured velocities, and the circulations were then estimated by fitting the measured velocities to a theoretical model. Frehlich et al. [22] established two maximum likelihood (ML) algorithms (a three-parameter model for symmetric wake vortices and a six-parameter model for asymmetric wake cortices) based on the comparison of vortex model and spectrum data. Later, Hallermeyer et al. [23] improved Frehlich’s algorithm. They estimated the vortex-core positions from the velocity envelopes, and the ML algorithm was only applied to estimate the circulation; this can reduce the computational load and enhance the robustness to a certain degree. Yoshikawa et al. [24] obtained the target parameters by matching the Doppler spectra by ML method.

This paper proposes an algorithm to retrieve the characteristic parameters of wake vortices (circulations, vortex-core positions, and vortex-core radii) by sensing the kinematics of aerosols driven by the ambient air flow (wake vortices and background wind) based on a scanning Doppler Lidar. A governing equation is established to describe the relationship between the measured data and target parameters, and an optimization method is used to estimate the target parameters from this equation. In this method, the inhomogeneity of background wind is considered, and the distortions of Doppler velocities caused by the scanning of Lidar beam are adjusted by the Doppler acceleration. This algorithm is specialized for the wake phase after fully roll-up and before rebound. As a whole, this algorithm belongs to the second category of parameter-retrieval algorithms mentioned above.

This paper is organized as follows. In Section 2, a governing equation related to the target parameters and Doppler velocities is established from which the target parameters are then solved with an optimization method. Section 3 quantitatively assesses the performance of the algorithm by simulations. In Section 4, the algorithm is applied to the field data recorded in Hong Kong international airport (HKIA) by a Windcube Lidar in 2014, and the results are analyzed.

2. Methodology

2.1. Geometry setup

The geometric relationship between the Lidar, runway, and wake vortices is shown in Fig. 1. The Lidar is installed on one side of the runway and its beam scans up and down alternately in an elevation angle interval [α_min, α_max] on the plane perpendicular to the runway. The aircraft, which is supposed to depart or land along the centerline of the runway, generates a pair of wake vortices along the runway, and the circulations of them are Γ_L and Γ_R, respectively. The two vortex cores are located at O_L and O_R, and the spacing between them is b₀. For a measurement cell centered at (x_i, y_i), the distances from it to the two vortex cores are R_L and R_R.

Fig. 1 Geometry setup of the Lidar, runway and wake vortices.

Download Full Size | PDF

2.2. Governing equation

The main principle of this algorithm is to retrieve the characteristic parameters of wake vortices by sensing the kinematic of aerosols driven by the ambient airflow (wake vortices and background wind). To do this, it is necessary to work out the relationship between the kinematics of aerosols and that of wake vortices. Typically, the Stokes number, defined as S_t = t_s/t_flow, is used to measure the impact of ambient flow on particles. Here t_flow is the characteristic time of wake t_flow ∼ $O (1)$ [9], t_s ≈ 2ρ_wr²/9µ is the aerosols’ relaxation time, with ρ_w, r and µ being the aerosols’ density, radius, and kinematic viscosity coefficient respectively. Normally, if S_t ≪ 1, the particles are with weak inertia, and their velocities are good approximations of the ambient flow velocity. It is known that the radii of aerosols vary from 1 nm to 100 µm [25], so we have t_s ∈ [2.80 × 10⁻¹¹, 0.280]. As a result, the Stokes number of aerosols satisfies S_t ≪ 1, and the aerosols’ kinematics are assumed to be identical to that of the ambient flow.

In this manner, for aerosols in the i^th measurement cell, the negative Doppler velocity $(- V_{D}^{(i)})$ can represent the projection of ambient flow velocity (V^flow,(ⁱ⁾) on the LOS (line of sight). Considering that the ambient flow is composed of background wind and wake vortices, we have the following postulate:

- V_{D}^{(i)} = {V^{flow, (i)} |}_{r} {= V^{wind, (i)} |}_{r} + {V^{wake, (i)} |}_{r},

where V^wind,(ⁱ⁾ and V^wake,(ⁱ⁾ are respectively the background wind velocity and the wake vortex velocity at the i^th measurement cell, X|_r = k·X denotes the projection of a velocity X on the radial direction k = [ cos α sin α ] with α being the elevation angle. This postulate directly leads to the following equation (governing equation in this paper) of the three velocity components

{V^{wind, (i)} |}_{r}

,

{V^{wake, (i)} |}_{r}

, and

V_{D}^{(i)}

:

ℱ^{(i)} = {V^{wind, (i)} |}_{r} + {V^{wake, (i)} |}_{r} + V_{D}^{(i)} = 0 .

We first give detailed analysis of these three velocity components in the following subsection.

2.3. Analysis of the three velocity components

2.3.1. ${V^{wind, (i)} |}_{r}$ : Projection of background wind velocity on LOS

There are generally two ways to estimate the background wind. The first assumes static background wind, where the wind measured at a time can represent that for all scans. The second way takes into account the inhomogeneity of background wind, i.e., its temporal and spatial variations. We propose calculating the background wind in the second way.

For a single scan, the RHI (Range-Height Indicator) is divided into WV (wake vortex) region and non-WV region (the region away from the wake). In section 2.4, we show that the two vortex cores [O_L, O_R] can be roughly located from the gradient extrema, then the region $[O_{L}^{x} - Δ x, O_{R}^{x} + Δ x]$ is assumed to be the WV region and those beyond WV region are taken as non-WV regions. Here, we assign a proper value comparable to wingspan to Δx, for example Δx = 60 m. Normally, the non-WV region is assumed to be independent of wake vortices, and the Doppler velocity is directly the background wind velocity on the $LOS: - V_{D}^{(j)} = {V^{wind, (j)} |}_{r}$ , where V⁽^j⁾ represents the velocity of the j^th measurement cell in the non-WV region. Based on this assumption, we can obtain ${V^{wind, (i)} |}_{r}$ for the i^th measurement cell in the WV region by the bilinear interpolation of ${V^{wind, (i)} |}_{r}$ .

2.3.2. ${V^{wake, (i)} |}_{r}$ : Projection of wake vortex velocity on LOS

The radial velocity of wake vortices at a given point can be obtained from existing wake velocity profile models. For example, the Burnham-Hallock (B-H) velocity model for a single vortex follows [26]:

{V^{wake} |}_{θ} = Γ \cdot {{\tilde{V}}^{wake, (i)} |}_{θ} = \frac{Γ}{2 π r} \cdot \frac{r^{2}}{r^{2} + r_{c}^{2}},

where Γ and r_c are the circulation and core radius respectively,

{{\tilde{V}}^{wake, (i)} |}_{θ}

is the normalized tangential velocity at the i^th measurement cell, and r is the distance from the measurement cell to the vortex core.

Then, for the i^th measurement cell, the normalized velocity for the left and right vortices can be expressed as:

{\tilde{V}}_{Lw}^{(i)} = {[{{\tilde{V}}_{L}^{(i)} |}_{x} {{\tilde{V}}_{L}^{(i)} |}_{y}]}^{T} = \frac{{{\tilde{V}}_{L}^{(i)} |}_{θ}}{R_{L}} {[y_{i} - O_{L}^{y} - (x_{i} - O_{L}^{x})]}^{T}

{\tilde{V}}_{Rw}^{(i)} = {[{{\tilde{V}}_{R}^{(i)} |}_{x} {{\tilde{V}}_{R}^{(i)} |}_{y}]}^{T} = \frac{{{\tilde{V}}_{R}^{(i)} |}_{θ}}{R_{R}} {[- (y_{i} - O_{R}^{y}) x_{i} - O_{R}^{x}]}^{T}

where

{{\tilde{V}}_{L}^{(i)} |}_{θ}

and

{{\tilde{V}}_{R}^{(i)} |}_{θ}

are the normalized tangential velocities for the left and right vortices obtained from Eq. (3)R_L and R_R are the distances from the i^th measurement cell to the two vortex cores:

R_{L} = {[{(x_{i} - O_{L}^{x})}^{2} + {(y_{i} - O_{L}^{y})}^{2}]}^{1 / 2}

,

R_{R} = {[{(x_{i} - O_{R}^{x})}^{2} + {(y_{i} - O_{R}^{y})}^{2}]}^{1 / 2}

.

It is seen that ${\tilde{V}}_{Lw}^{(i)}$ and ${\tilde{V}}_{Rw}^{(i)}$ are strongly related to the vortex-core positions and the radius, so the radial velocity of wake vortices at the point (x_i, y_i) can be expressed as a function of Γ_L, Γ_R, O_L, O_R and r_c:

{V^{wake, (i)} |}_{r} = k \cdot [Γ_{L} \cdot {\tilde{V}}_{Lw}^{(i)} (O_{L}, O_{R}, r_{c}) + Γ_{R} \cdot {\tilde{V}}_{Rw}^{(i)} (O_{L}, O_{R}, r_{c})]

where k is the unit wave number vector indicating the line of sight (LOS), and the velocity of a point is assumed to be induced by the two vortices.

2.3.3. $V_{D}^{(i)}$ : Doppler velocity and the adjustment of its distortion

We try to retrieve the target parameters using the data over a whole RHI. That is to say, the parameters we get correspond to the center time of the RHI. For the detection scenario in this paper, the Lidar beam scans up and down alternately, and each RHI is composed of dozens of beams at different time samples. Compared to the velocity distribution at the center time sample, the measured Doppler velocities are distorted due to the movement of ambient flow, especially when a slow scan rate is used. Therefore, we adjust the distorted Doppler velocities for better estimations of the target parameters.

The Doppler velocity at the center time $\bar{t}$ for the i^th measurement cell can be adjusted by Doppler acceleration as follows [8]:

V_{D}^{(i)} (\bar{t}) = V_{D}^{(i)} (t_{i}) + \int_{t_{i}}^{\bar{t}} A_{D}^{(i)} (τ) d τ,

where t_i is the sample time of the measurement cell,

V_{D}^{(i)} (t_{i})

and

V_{D}^{(i)}

are the measured Doppler velocity and Doppler acceleration of the measurement cell respectively. Theoretically, the Doppler velocity of the measurement cell located at (x_i, y_i) can be derived from Eq. (4) and Eq. (5):

V_{D}^{(i)} (t) = - V_{r}^{(i)} (t) = - [\cos α_{i} \sin α_{i}] \cdot [(Γ_{L} \cdot {\tilde{V}}_{Lw}^{(i)} + Γ_{R} \cdot {\tilde{V}}_{Rw}^{(i)}) + V^{wind}],

Basically, the variations of Γ and V^wind in a RHI are very small, so the acceleration $A_{D}^{(i)}$ can be approximated as follows if these variations are neglected:

A_{D}^{(i)} = \frac{\partial V_{D}^{(i)} (t)}{\partial t} \approx - \frac{1}{2 π} [\cos α_{i} \sin α_{i}] \cdot [\begin{array}{l} Γ_{L} \cdot Ψ_{L}^{(i)} - Γ_{R} \cdot Ψ_{R}^{(i)} \\ - Γ_{L} \cdot Φ_{L}^{(i)} + Γ_{R} \cdot Φ_{R}^{(i)} \end{array}],

where

Ψ_{L}^{(i)} = {V_{L}^{down, (i)} |}_{y} \cdot \frac{{(y_{i} - O_{L}^{y})}^{2} - {(x_{i} - O_{L}^{x})}^{2} - r_{c}^{2}}{{(R_{L}^{2} + r_{c}^{2})}^{2}} + {V_{L}^{wind, (i)} |}_{x} \cdot \frac{2 (x_{i} - O_{L}^{x}) (y_{i} - O_{L}^{y})}{{(R_{L}^{2} + r_{c}^{2})}^{2}},

Ψ_{R}^{(i)} = {V_{R}^{down, (i)} |}_{y} \cdot \frac{{(y_{i} - O_{R}^{y})}^{2} - {(x_{i} - O_{R}^{x})}^{2} - r_{c}^{2}}{{(R_{R}^{2} + r_{c}^{2})}^{2}} + {V_{R}^{wind, (i)} |}_{x} \cdot \frac{2 (x_{i} - O_{R}^{x}) (y_{i} - O_{R}^{y})}{{(R_{R}^{2} + r_{c}^{2})}^{2}},

Φ_{L}^{(i)} = {V_{L}^{wind, (i)} |}_{x} \cdot \frac{{(x_{i} - O_{L}^{x})}^{2} - {(y_{i} - O_{L}^{y})}^{2} - r_{c}^{2}}{{(R_{L}^{2} + r_{c}^{2})}^{2}} + {V_{L}^{down, (i)} |}_{y} \cdot \frac{2 (x_{i} - O_{L}^{x}) (y_{i} - O_{L}^{y})}{{(R_{L}^{2} + r_{c}^{2})}^{2}},

Φ_{R}^{(i)} = {V_{R}^{wind, (i)} |}_{x} \cdot \frac{{(x_{i} - O_{R}^{x})}^{2} - {(y_{i} - O_{R}^{y})}^{2} - r_{c}^{2}}{{(R_{R}^{2} + r_{c}^{2})}^{2}} + {V_{R}^{down, (i)} |}_{y} \cdot \frac{2 (x_{i} - O_{R}^{x}) (y_{i} - O_{R}^{y})}{{(R_{R}^{2} + r_{c}^{2})}^{2}} .

In the above expressions, ${V_{L}^{wind, (i)} |}_{x}$ and ${V_{R}^{wind, (i)} |}_{x}$ are the horizontal component of the background wind velocity at the left and right vortex cores, ${V_{L}^{down, (i)} |}_{y}$ and ${V_{R}^{down, (i)} |}_{y}$ are the descending velocity of the left and right vortex cores:

\begin{array}{l} {V_{L}^{down, (i)} |}_{y} = - (O_{R}^{x} - O_{L}^{x}) \cdot \frac{Γ_{R}}{2 π (b_{0}^{2} + r_{c}^{2})}, \\ {V_{R}^{down, (i)} |}_{y} = - (O_{R}^{x} - O_{L}^{x}) \cdot \frac{Γ_{L}}{2 π (b_{0}^{2} + r_{c}^{2})}, \end{array}

with

b_{0} = {[{(O_{R}^{x} - O_{L}^{x})}^{2} + {(O_{R}^{y} - O_{L}^{y})}^{2}]}^{1 / 2}

being the spacing between the two vortex cores.

It has to be mentioned that the vertical component of background wind velocity ${V^{wind} |}_{y}$ and the horizontal component of descending velocity ${V^{down} |}_{x}$ caused by the height difference of the two vortex cores are generally very small, thus the corresponding displacements of vortex-core positions are neglected here. From Eq. (9) and Eq. (10), it is observed that $A_{D}^{(i)}$ is related to the target parameters [Γ_L, Γ_R, O_L, O_R, r_c] and the background wind velocity V^wind. Substituting $A_{D}^{(i)}$ into Eq. (7), the distorted Doppler velocity $V_{D}^{(i)}$ is adjusted.

2.4. Preliminary determination of the vortex-core positions

From the above analysis, it is seen that the vortex-core positions are very important to the velocity of wake vortices, thus we first give a preliminary estimation of the vortex-core positions.

As shown in Fig. 2(a), the velocity field of the aircraft wake usually has a crossover structure of positive and negative values. According to the velocity model, when an aircraft flies from the right side of the Lidar to the left, the velocity gradient in the wake region should have a negative peak on the left and a positive one on the right i.e., ${\bar{O}}_{L} = {[{\bar{O}}_{L}^{x}, {\bar{O}}_{L}^{y}]}^{T} = \underset{x, y}{\arg \min} [\partial V_{D} (x, y) / \partial y]$ and ${\bar{O}}_{R} = {[{\bar{O}}_{R}^{x}, {\bar{O}}_{R}^{y}]}^{T} = \underset{x, y}{\arg \max} [\partial V_{D} (x, y) / \partial y]$ . Theoretically, there is only a pair of negative-positive peaks in the wake region, but multiple gradient maxima could be observed in some non-ideal cases: 1) the vortices are not fully rolled up, 2) old vortices of previous aircraft still exist, 3) the presence of noise and measurement errors, and others. In order to enhance the robustness of the algorithm, this paper proposes to discriminate the vortex cores from two positive gradient peaks and two negative gradient peaks by the geometric constraints. As shown in Fig. 2(b), points A and B are the two negative gradient peaks, while points C and D are the two positive ones. Then four negative-positive pairs (AC, AD, BC, BD) are developed and their spacings are denoted as r_AC, r_AD, r_BC and r_BD, respectively.

Direction constraint: as mentioned above, the positive peak should be on the right side of the negative one if the aircraft flies from the right to the left side of the Lidar. In this way, the pair of B and C is ruled out.
Spacing constraints: the distance, horizontal, and vertical spacings between negative and positive peaks should within certain bounds. For example, 25 m < r_pq < 90 m, 25 m < |x_p−x_q| < 90 m and |y_p−y_q| < 30 m. Here, p represents the gradient’s minimum (A or B), while q represents the gradient’s maximum (C or D). The combination of B and D is ruled out because |x_B − x_D| = 21.9 m < 25 m.

Fig. 2 Illustration of vortex-core positions’ preliminary estimation: distribution of Doppler velocity (a) and the corresponding velocity gradient (b).

Download Full Size | PDF

Among the remaining pairs, the pair with the largest absolute gradient values is chosen as the initial vortex-core pair. Here, only the pair A and C, and the pair A and D satisfy the constraints mentioned above; the former pair is kept due to its larger absolute gradient values. The resulting vortex-core positions are only used to set up the boundaries of vortex-core positions for the optimization which is used in next subsection, thus the above preliminary estimation can fulfill this purpose well.

2.5. Solving the governing equation with an optimization method

When substituting the three velocity components in section 2.3 into Eq. (2), the governing equation for a point (x_i, y_i) can be further rewritten as:

{ℱ^{(i)} (Γ_{L}, Γ_{R}, O_{L}, O_{R}, r_{c}; V_{D}^{(i)}, x_{i}, y_{i}) = V^{wind, (i)} |}_{r} + {V^{wake, (i)} (Γ_{L}, Γ_{R}, O_{L}, O_{R}, r_{c}) |}_{r} + V_{D}^{(i)} = 0

When sufficient measurement data $[x_{i}, y_{i}, V_{D}^{(i)}]$ are used, the estimated target parameters $[{\hat{Γ}}_{L}, {\hat{Γ}}_{R}, {\hat{O}}_{L}, {\hat{r}}_{c}]$ can be obtained from the governing equation with an optimization method (such as the nonlinear least square method):

\begin{array}{l} [{\hat{Γ}}_{L}, {\hat{Γ}}_{R}, {\hat{O}}_{L}, {\hat{O}}_{R}, {\hat{r}}_{c}] = \underset{Γ_{L}, Γ_{R}, O_{L}, O_{R}, r_{c}}{\arg \min} \sum_{i = 1}^{M} {| ℱ^{(i)} (Γ_{L}, Γ_{R}, O_{L}, O_{R}, r_{c}; V_{D}^{(i)}, x_{i}, y_{i}) |}^{2} \\ s . t . {\hat{O}}_{L} \in ({\hat{O}}_{L} + [- 20, 20] \times [- 20, 20]), \\ {\hat{O}}_{R} \in ({\hat{O}}_{R} + [- 20, 20] \times [- 20, 20]), \\ {\hat{Γ}}_{L} \in [0, 800], {\hat{Γ}}_{R} \in [0, 800], {\hat{r}}_{c} \in [0, 6], \end{array}

where M is the number of measurement cells in the WV region, and the unknown variables are limited to certain boundaries for better convergence.

For each scan, we first solve the optimization problem Eq. (12) to get the global optimization values $[{\hat{Γ}}_{L}, {\hat{Γ}}_{R}, {\hat{O}}_{L}, {\hat{O}}_{R}, {\hat{r}}_{c}]$ . These optimization values are substituted into Eq. (9) in section 2.3.3 to obtain the acceleration $A_{D}^{(i)}$ , which is further used to adjust the distorted Doppler velocities in the scan according to Eq. (7). When the Doppler velocities of each RHI are adjusted, we can then substitute them into Eq. (12) to estimate the target parameters $[{\hat{Γ}}_{L}, {\hat{Γ}}_{R}, {\hat{O}}_{L}, {\hat{O}}_{R}, {\hat{r}}_{c}]$ .

During the optimization process for Eq. (12), we may have multiple extrema (global and local) of $\sum_{i = 1}^{M} {| ℱ^{(i)} |}^{2}$ which correspond to different groups of target parameter estimations [see Fig. 3(a)]. Generally, the global optimization values would give the best estimations of the target parameters, but the presence of non-ideal effects (noise, measurement error, etc.) may make the global optimization values deviate from the real parameters occasionally, and a certain local optimization result could be even more reasonable.

Fig. 3 Sketch of global and local optimization values (a), and illustration of determining the trajectories of vortices, where the local optimization results for point D are presented (b).

Download Full Size | PDF

Typically, for wake vortices after roll-up and before rebound, the vortex cores’ trajectory should be relatively smooth, so the global optimization result, with the corresponding vortex-core position deviating significantly from the overall tendency of the vortex-core trajectory, is more likely to be an outlier; it should be replaced by a more proper local optimization result. For example, the red triangles (A_L, B_L, ⋯, G_L) and blue circles (A_R, B_R, ⋯, G_R) in Fig. 3(b) respectively indicate the left and right vortex-core positions of the global optimization results for seven consecutive RHIs. For the left vortex cores, points A_L, B_L, ⋯, G_L can be fitted to a smooth dark line after removing some aberrant points. This line shows the tendency of the left vortex core, and two dashed lines constrained by a user-defined threshold are used to indicate whether a point is far away from the smooth line. For the right vortex, we have similar smooth line and threshold lines. Basically, if a vortex-core position is far away from the smooth line [see D_L, D_R and F_L in Fig. 3(b)], the global optimization value is more likely to be an outlier for this scan, and we should find out whether there exists a local optimization value closer to the smooth lines. For example, there are four pairs of local optimization values for point D: {D_L, D_R}, ${D_{L}^{(1)}, D_{R}^{(1)}}$ , ${D_{L}^{(2)}, D_{R}^{(2)}}$ and ${D_{L}^{(3)}, D_{R}^{(3)}}$ and their derivations to the smooth lines are respectively $d_{L}^{D} + d_{R}^{D}$ , $d_{L}^{D, (1)} + d_{R}^{D, (1)}$ , $d_{L}^{D, (2)} + d_{R}^{D, (2)}$ and $d_{L}^{D, (3)} + d_{R}^{D, (3)}$ , with ${d_{L}^{D}, d_{R}^{D}}$ , ${d_{L}^{D, (1)}, d_{R}^{D, (1)}}$ , ${d_{L}^{D, (2)}, d_{R}^{D, (2)}}$ , ${d_{L}^{D, (3)}, d_{R}^{D, (3)}}$ being the distances of these points to the corresponding smooth lines. It is found that $d_{L}^{D, (2)} + d_{R}^{D, (2)}$ is the minimum among them, so the local optimization result for ${D_{L}^{(2)}, D_{R}^{(2)}}$ is used to replace the global result {D_L, D_R}. Of course, other outliers such as F can be dealt with in the same way.

3. Validation with simulations

To quantitatively verify the algorithm’s performance, a Doppler velocity RHI simulator has been developed. In the simulator, the parameters are set as follows.

Dynamics simulation parameters: the background wind field is composed of crosswind (with vertical wind shear) and turbulence (in Von Kalman model) [27]. The crosswind is [−1 − 0.03h, 0] m/s with h being the height from ground, and the eddy dissipation rate (EDR) of atmospheric turbulence is ε = 0.003 m²/s³. The evolution of wake vortices is assumed to obey the D2P model [28]. The initial circulations of the two vortices are both Γ₀ = 400 m²/s, and the initial positions for the left and right vortices are O_L(t₀) = (450, 67) m and O_R(t₀) = (510, 67) m respectively.
Lidar echo simulation parameters: the wavelength is λ = 1.54 µm, the sampling rate is B = 500 MHz, the pulse duration is τ = 120 ns, the signal-to-noise ratio of a single shot is SNR = −10 dB, and 500 shots are used for spectra accumulation. The Lidar beam scans up and down alternately at a scanning rate of ±1.99 °/s, and a full scan takes about five seconds. The 500 MHz sampling rate corresponds to a range gate of 0.3 m. Twenty consecutive samples are adopted to get the Doppler spectrum of these range gates with FFT, then the Doppler velocity is consequently estimated from the spectrum peak.

The first and last RHIs of 13 simulated scans are shown in Fig. 4.

Fig. 4 Simulated Doppler velocity distribution: the first scan (a) and the last scan (b).

Download Full Size | PDF

By applying the proposed algorithm to the simulated Doppler velocity distributions, the parameter-retrieval results are obtained as shown in Fig. 5. For comparison, three sets of results are presented.

Results of the new method when the distortions of Doppler velocities are NOT adjusted [see Figs. 5(a) and 5(b)].
Results of the new method when the distortions of Doppler velocities are adjusted by the method proposed in section 2.3.3 [see Figs. 5(c) and 5(d)].
Results of the Tangential Velocities (TV) method introduced in [29] [see Figs. 5(e) and 5(f)].

Fig. 5 Vortex-core trajectory and circulations retrieved from simulated Doppler velocity distribution. Results of the new method without adjusting the Doppler velocity: vortex-core positions (a), circulations (b); Results of the new method with the Doppler velocity adjusted: vortex-core positions (c), circulations (d); Results of the TV method: vortex-core positions (e), circulations (f).

Download Full Size | PDF

From Figs. 5(a), 5(c) and 5(e), special phenomenon of the vortex-core’s positions is observed. When the Doppler velocity’s distribution is not adjusted, the vortex-core positions retrieved by the new method gather in pairs alternately and deviate from the true values, i.e., the pairs marked by ellipses in Fig. 5(a), with Δ₁ and Δ₂ being the separations between two consecutive vortex cores. This pair-wise phenomenon [marked by ellipses in Fig. 5(e)] is also evident for the TV method (Δ₁ < Δ₂). But when the distortions of Doppler velocities are adjusted [see Fig. 5(c)], the separations (Δ₁, Δ₂, ⋯) are almost equal, and they agree with the true values much better. This shows that the Doppler velocity adjustment can effectively help to improve the performance of the new method.

For the circulations [Figs. 5(b), 5(d) and 5(f)], the overall tendency of the estimated results agree with the true values, but the new method’s results [Figs. 5(b) and 5(d)] are more accurate in comparison with those using the TV method [Fig. 5(f)].

To give statistically analysis of the algorithm, a dozen of realizations of the Doppler velocity distributions are obtained based on the random realizations of the background turbulence (with the same intensity as that shown above). The estimated vortex-core positions and circulations are gathered in Fig. 6.

Fig. 6 The estimated vortex-core trajectories (a) and circulations (b) for different background turbulence.

Download Full Size | PDF

The relative error and relative root mean square error (RMSE) are defined as follows:

\begin{array}{l} E_{r} = \frac{1}{N} \sum_{i = 1}^{N} \frac{1}{P_{T}^{(i)}} \cdot | 〈 P_{E}^{(i)} 〉 - P_{T}^{(i)} | \times 100 % (Relative error) \\ R M S E_{r} = \frac{1}{N} \sum_{i = 1}^{N} (\frac{1}{Γ_{0}} \sqrt{〈 {(P_{E}^{(i)} - P_{T}^{(i)})}^{2} 〉}) \times 100 % (Relative RMSE of circulation) \\ R M S E_{r} = \frac{1}{N} \sum_{i = 1}^{N} (\frac{1}{b_{0}} \sqrt{〈 {(P_{E}^{(i)} - P_{T}^{(i)})}^{2} 〉}) \times 100 % (Relative RMSE of vortex - core height) \end{array}

Here, $P_{E}^{(i)}$ and $P_{T}^{(i)}$ are respectively the estimated value and true value of a parameter for the i^th scan, and N is the total number of the considered scans. $〈 x^{(i)} 〉 = \frac{1}{M} \sum_{j = 1}^{M} x^{(i, j)}$ is the average value of x over the realizations at the i^th scan, and M is the total number of random realizations. The relative errors and relative RMSEs of the estimated parameters are listed in Table 1. When the distortions of Doppler velocities are adjusted, the relative errors of the estimated circulations and two vortex-core positions by the new method are 6.24%, 3.15%, and 2.39%, and the corresponding relative RMSEs are 7.91%, 3.94%, and 3.78% respectively.

Table 1. Average relative errors and relative RMSE of the estimated parameters

View Table | View all tables in this article

From the above error analysis, it can be seen that the new method can give very good estimations of the target parameters, especially when the Doppler velocities are adjusted.

4. Retrieval results of field data

This algorithm has been applied to the Lidar measurement data recorded at Hong Kong International Airport (HKIA) in 2014. This section presents some basic analysis of that. More detailed descriptions of the measurement campaigns can be found in [30].

In the campaign, a Windcube 200s Lidar, whose parameters are listed in Table 2, was temporarily set on the roof of Asia World-Expo (about 19 m above the ground) to observe the north runway (25 RA) of HKIA at a distance of 400 m from south (see Fig. 7). The Lidar scans up and down alternately on the plane perpendicular to the runway, and outputs the Doppler velocities in RHI form.

Table 2. Main parameters of the pulsed Lidar in HKIA

View Table | View all tables in this article

Fig. 7 Layout of Lidar and the north runway (25RA) in HKIA, where the vertical scanning plane is the solid line in red (a), and a detailed setup of the Lidar (b).

Download Full Size | PDF

From section 2.3.1 and section 2.4, a RHI shown in Fig. 8 can be divided into three regions: WV region (Region ①, between the two black lines) inhabited by the vortices, non-WV region (Region ②, between the red lines and black lines) where only the ambient wind is presented, and the rest region (Region ③) which is discarded in the parameter-retrieval process.

Fig. 8 Doppler velocity RHI of wake vortex obtained in HKIA at 07:17:46–07:17:51 (UTC) on August 3, 2014, where the origin is the Lidar position (at a height of 19 m to the ground).

Download Full Size | PDF

Figure 8 shows the distribution of measured Doppler velocity at 07:17:46–07:17:51 (UTC) on August 3, 2014. The distance between the centerline of Region ① and the Lidar is 375 m, and the widths of Region ① and Region ② are set as 225 m and 200 m respectively. The initial positions for vortex cores are estimated approximately with the method in section 2.4 (see symbols Δ and ◯ in Fig. 8). The background velocity in Region ① is calculated with the method introduced in section 2.3.1 with the measured data in Region ②. At the same time, wake vortices move along with the background wind, and the corresponding displacement should be considered for Region ①.

Based on the above configuration, the data in Region ① are adopted to solve the governing equation to get the target parameters, where the distortions of Doppler velocities are fixed by the method proposed in section 2.3.3. Four sets of data are used in this section to show the applicability of the new method. The aircraft types retroactively obtained from the departing and landing schedules are A300F4-200, B747-400, B777-300ER, and A330-300 respectively.

4.1. Retrieval results of vortex-core positions

The retrieved vortex-core positions from the four sets of field data are shown in Fig. 9, where we can see that all the left and right vortex cores descend along with time. The initial vortex-core separations for A300F4-200, B747-400, B777-300ER, and A330-300 are approximately obtained as 36 m, 50 m, 50 m, 50 m according to b₀ = πB/4 with B is the wingspan of an aircraft [31]. The separations between two estimated vortex cores in the initial period are approximately equal to the theoretic initial values accordingly.

Fig. 9 Vortex-core position retrieval results from four sets of data on August 3, 2014, where the range is the horizontal distance to the Lidar, and the height is relative to the ground, Airbus A300F4-200 at 00:38:22–00:39:36 (UTC) (a); Boeing 747-400 at 00:52:47–00:54:07 (UTC) (b); Boeing 777-300ER at 07:17:25–07:18:23 (UTC) (c); Airbus A330-300 at 08:35:46–08:36:39 (UTC) (d).

Download Full Size | PDF

It should be noted that the wake vortices were invisible to the Lidar when the vortex cores descended to a height below 24.8 m (relative to the ground). The main reason is that the Lidar was temporarily set on the roof of a high building (Asia World-Expo) and the lowest elevation angle is 0.83°, resulting in a lowest observable height of 19 + 400 · sin(0.83°) = 24.8 m. Normally, for wake vortex above this height (about b₀/2) the ground effect is not obvious, so in this paper we don’t consider the ground effect.

4.2. Retrieval results of circulations

The circulation retrieval results are shown in Fig. 10. It is seen that even though some fluctuations are observed, the circulation curves decrease overall. The initial circulations of vortices generated by the aircrafts are calculated theoretically with Γ₀ = M_Ag/ρ_ab₀V_A, where M_A, g, ρ_a, V_A are the aircraft mass, gravity acceleration, air density, and aircraft speed, respectively [31]. Since we did not record the exact flight information, we can only present some empirical values of V_A and M_A. In this paper, V_A is approximately set as 130 kn, M_A is a value between the maximum landing weight and no-load weight of the aircraft. These approximations result in the boundaries of initial circulations shown in Fig. 10 (the gray regions). It is seen that the circulations estimated in the initial period are within the gray region.

Fig. 10 Circulation retrieval results from four sets of data on August 3, 2014. Airbus A300F4-200 at 00:38:22–00:39:36 (UTC) (a); Boeing 747-400 at 00:52:47–00:54:07 (UTC) (b); Boeing 777-300ER at 07:17:25–07:18:23 (UTC) (c); Airbus A330-300 at 08:35:46–08:36:39 (UTC) (d).

Download Full Size | PDF

It has been mentioned that the Lidar was set on the roof of a 19 m-high building and the lowest elevation angle was 0.83°, so the wake vortices of longer evolution time (when the vortices descend below 24.8 m) were invisible and we are unable to analysis the further phase of wake vortices from these data.

4.3. Retrieval results of vortex-core radii

The evolution of vortex-core radii along with time is shown in Fig. 11, where r_c = 0.052b₀ is plotted for comparison [31]. The results show that the estimated vortex-core radii for these four wake vortices are around 0.05b₀; this is very close to the theoretical prediction.

Fig. 11 Retrieval results of vortex-core radii from the four sets of data.

Download Full Size | PDF

In conclusion, for wake vortices after roll-up and before rebound, the following phenomena are observed from Fig. 9, Fig. 10, and Fig. 11.

The estimated vortex-core trajectories are smooth, and the estimated initial separations between two vortex cores are approximately equal to the corresponding theoretic initial values b₀.
The estimated circulations decrease along with vortex age, and the circulations estimated in the initial period are close to the theoretic initial values Γ₀.
The estimated vortex-core radii are almost constant, and their average value is approximately equal to 0.05b₀.

The above findings have roughly verified the applicability of the new method on the field data. For future research, we will set the Lidar on the ground to expand its lowest detection regions and take into account the ground effect.

5. Conclusion

This paper proposes an algorithm to retrieve the characteristic parameters of dry-air wake vortices by sensing the kinematic of aerosols driven by ambient airflow (wake vortices and background wind). Since the inertia of aerosols is very weak, the measured Doppler velocity can be directly regarded as the projected velocity of ambient airflow on the LOS. Based on this postulate, a governing equation related to the background velocity, wake vortex velocity, and Doppler velocity is established, from which the target parameters are obtained with an optimization method. During the process, the distortion of Doppler velocity caused by the scanning of Lidar beam is adjusted by the Doppler acceleration. Good performance of the new method has been verified by simulations and field Lidar detection data recorded at HKIA in 2014.

Compared with the existing parameter-retrieval methods for dry-air wake vortices, the newly proposed method has some advantages: 1) Distortions of Doppler velocities caused by Lidar beam scanning are adjusted, which can give better estimations of target parameters. For example, the pair-wise phenomenon of estimated vortex-core positions in existing methods is fixed well; 2) Better adaptability to complex background wind field since it can give a good estimation of the background wind field in both temporal (estimated for each scan) and spatial (interpolate from non-WV region to WV region) aspects; 3) For the wake phase after roll-up and before rebound, smoother vortex-core trajectories can be obtained by selecting the best vortex-core positions from all the optimization results.

Some improvements of the research are expected in the near future: 1) In this paper, we focus mainly on the wake vortex phase after roll-up and before rebound. Extension of the present work to other phases will be tried in the near future; 2) The Lidar will be placed on a ground site to expand its lowest detection region, and then the performance for longer evolution time can be verified.

Funding

National Natural Science Foundation of China (NSFC) (No. 61771479, No. 41375040, No. 61490649).

Acknowledgments

Hang Gao and Jianbing Li contributed equally to this work.

References and links

1. P. R. Veillette, “Data show that U. S. wake turbulence accidents are most frequent at low altitude and during approach and landing,” Flight Safty Digest 21(3–4), 1–47 (2002).

2. T. Gerz, F. Holzäpfel, W. Bryant, F. Köpp, M. Frech, A. Tafferner, and G. Winckelmans, “Research towards a wake-vortex advisory system for optimal aircraft spacing,” C. R. Phys. 6(4), 501–523 (2005). [CrossRef]

3. J. Li, H. Gao, T. Wang, and X. Wang, “A survey of the scattering characteristics and detection of aircraft wake vortices,” J. Radar 6(6), 653–672 (2017).

4. SESAR Consortium, European Air Traffic Management Master Plan, Edition 1 (2009).

5. J. Burden, P. A. Curry, D. Roby, and F. Love, “Introduction to the next generation automatic test system (NGATS),” in Proceedings of IEEE Autotestcon (IEEE, 2005), pp. 16–19.

6. P. A. Curry, J. Burden, and G. A. Lundy, “Next generation automatic test system (NGATS) update,” in Proceedings of IEEE Autotestcon (IEEE, 2006), pp. 318–322.

7. J. Li, X. Wang, T. Wang, J. Liu, H. Gao, and V. Chandrasekar, “Circulation retrieval of wake vortex under rainy condition with a vertically pointing radar,” IEEE Trans. Aerosp. Electron. Syst. 53(4), 1893–1906 (2017). [CrossRef]

8. J. Li, H. Gao, Y. Li, V. Chandrasekar, and X. Wang, “Circulation retrieval of simulated wake vortices under rainy condition with a side-looking scanning radar,” IEEE Trans. Aerosp. Electron. Syst. 54(2), 569–584 (2018). [CrossRef]

9. J. Li, T. Wang, L. Qu, and X. Wang, “Circulation retrieval of wake vortex in fog with an upward-looking monostatic radar,” IEEE Trans. Aerosp. Electron. Syst. 52(1), 169–180 (2016). [CrossRef]

10. J. Li, T. Wang, Z. Liu, and X. Wang, “Circulation retrieval of wake vortex in fog with a side-looking scanning radar,” IEEE Trans. Aerosp. Electron. Syst. 52(5), 2242–2254 (2016). [CrossRef]

11. F. Barbaresco, L. Thobois, A. Dolfi-Bouteyre, N. Jeannin, R. Wilson, M. Valla, A. Hallermeyer, P. Feneyrou, V. Brion, L. Besson, J. P. Cariou, L. Leviandier, G. Pillet, and D. Dolfi, “Monitoring wind, turbulence and aircraft wake vortices by high resolution radar and Lidar remote sensors in all weather conditions,” in Proceedings of URSI-France, wakenet.eu (2015), pp. 81–110.

12. L. Thobois, “Next generation scanning Lidar systems for optimizing wake turbulence separation minima,” J. Radar 6(6), 689–698 (2017).

13. S. Rahm, I. Smalikho, and F. Köpp, “Characterization of aircraft wake vortices by airborne coherent Doppler Lidar,” J. Aircr. 44(3), 799–805 (2007). [CrossRef]

14. S. Rahm and I. Smalikho, “Aircraft wake vortex measurement with airborne coherent Doppler Lidar,” J. Aircr. 45(4), 1148–1155 (2008). [CrossRef]

15. H. Wassaf, D. Burnham, and F. Wang, “Wake vortex tangential velocity adaptive spectral (TVAS) algorithm for pulsed Lidar systems,” presented at the 16th Coherent Laser Radar Conference, Long Beach, CA, 20–24 Jun. 2011.

16. D. J. Ramsey and C. Nguyen, “Characterizing aircraft wake vortices with ground-based pulsed coherent Lidar: effects of vortex circulation strength and Lidar signal-to-noise ratio on the spectral signature,” presented at the 3rd AIAA Atmospheric Space Environments Conference, Honolulu, Hawaii, 27–30 Jun. 2011.

17. D. Jacob, D. Y. Lai, D. P. Delisi, K. S. Barr, D. A. Hutton, S. Shald, M. H. Stephen, and G. Philip, “Assessment of Lockheed Martin’s aircraft wake vortex circulation estimation algorithms using simulated Lidar data,” presented at the 3rd AIAA Atmospheric Space Environments Conference, Honolulu, Hawaii, 27–30 Jun. 2011.

18. D. Jacob, D. Y. Lai, and D. P. Delisi, “Development of an improved pulsed Lidar circulation estimation algorithm and performance results for Denver OGE data,” presented at the 51st AIAA Aerospace Sciences Meeting including the New Horizons Forum and Aerospace Exposition, Grapevine, Texas, 07–10 Jan. 2013.

19. D. Jacob, M. J. Pruis, D. Y. Lai, and D. P. Delisi, “WakeMod 4: A new standalone wake vortex algorithm for estimating circulation strength and position,” presented at the 7th AIAA Atmospheric and Space Environments Conference, Dallas, Texas, 22–26 Jun. 2015.

20. I. N. Smalikho, V. A. Banakh, F. Holzäpfel, and S. Rahm, “Method of radial velocities for the estimation of aircraft wake vortex parameters from data measured by coherent Doppler Lidar,” Opt. Express 23(19), A1194–A1207 (2015). [CrossRef] [PubMed]

21. I. N. Smalikho and V. A. Banakh, “Estimation of aircraft wake vortex parameters from data measured by a stream line Lidar,” Proc. SPIE 9680, 968037 (2015). [CrossRef]

22. R. Frehlich and R. Sharman, “Maximum likelihood estimates of vortex parameters from simulated coherent Doppler Lidar data,” J. Atmos. Ocean. Technol. 22(2), 117–130 (2005). [CrossRef]

23. A. Hallermeyer, A. Dolfi-Bouteyre, M. Valla, L. Le Brusquet, G. Fleury, L. Thobois, J. P. Cariou, M. Duponcheel, and G. Winckelmans, “Development and assessment of a wake vortex characterization algorithm based on a hybrid Lidar signal processing,” presented at the 8th AIAA Atmospheric and Space Environments Conference, Washington, DC, 13–17 Jun. 2016.

24. E. Yoshikawa and N. Matayoshi, “Aircraft wake vortex retrieval method on Lidar lateral range-height indicator observation,” AIAA J. 55(7), 2269–2278 (2017). [CrossRef]

25. W. C. Hinds, Aerosol Technology: Properties, Behavior, and Measurement of Airborne Particles (Wiley, 1999).

26. D. C. Burnham and J. N. Hallock, “Chicago monostatic acoustic vortex sensing system,” U.S. Department of Transportation. Report No. DOT-TSC-FAA-79-103. 1982. 206 pp.

27. D. K. Wilson, V. E. Ostashev, G. H. Goedecke, and H. J. Auvermann, “Quasi-wavelet calculations of sound scattering behind barriers,” Appl. Acoust. 65(6), 605–627 (2004). [CrossRef]

28. F. Holzäpfel, “Probabilistic two-phase wake vortex decay and transport model,” J. Aircr. 40(2), 323–331 (2003). [CrossRef]

29. F. Holzäpfel, T. Gerz, F. Köpp, E. Stumpf, M. Harris, R. I. Young, and A. Dolfi-Bouteyre, “Strategies for circulation evaluation of aircraft wake vortices measured by Lidar,” J. Atmos. Ocean. Technol. 20, 1183–1195 (2003). [CrossRef]

30. K. K. Hon and P. W. Chan, “Aircraft wake vortex observations in Hong Kong,” J. Radar 6(6), 709–718 (2017).

31. T. Gerz, F. Holzäpfel, and D. Darracq, “Commercial aircraft wake vortices,” Prog. Aerospace Sci. 38, 181–208 (2002). [CrossRef]

Errors	New method without velocity adjustment	New method with velocity adjustment	TV method [29]
Γ_mean relative error%	7.32	6.24	23.58
Γ_mean relative RMSE %	9.73	7.91	22.44
$O_{L}^{y}$ relative error %	3.48	3.15	4.23
$O_{L}^{y}$ relative RMSE%	4.40	3.94	5.72
$O_{R}^{y}$ relative error%	2.93	2.39	8.47
$O_{R}^{y}$ relative RMSE%	5.14	3.78	7.45

Parameters	Default Value
Wavelength (μm)	1.54
Pulse Duration (ns)	200
Pulse Repetition Frequency (kHz)	20
Range (m)	[50, 6000]
Elevation Angle (°)	[0.83, 10.77]
Angular Resolution (°)	0.43
Scanning Rate (°/s)	1.87 (upward), −1.87 (downward)

Errors	New method without velocity adjustment	New method with velocity adjustment	TV method [29]
Γ_mean relative error%	7.32	6.24	23.58
Γ_mean relative RMSE %	9.73	7.91	22.44
$O_{L}^{y}$ relative error %	3.48	3.15	4.23
$O_{L}^{y}$ relative RMSE%	4.40	3.94	5.72
$O_{R}^{y}$ relative error%	2.93	2.39	8.47
$O_{R}^{y}$ relative RMSE%	5.14	3.78	7.45

Parameters	Default Value
Wavelength (μm)	1.54
Pulse Duration (ns)	200
Pulse Repetition Frequency (kHz)	20
Range (m)	[50, 6000]
Elevation Angle (°)	[0.83, 10.77]
Angular Resolution (°)	0.43
Scanning Rate (°/s)	1.87 (upward), −1.87 (downward)

Parameter-retrieval of dry-air wake vortices with a scanning Doppler Lidar

Abstract

1. Introduction