# Making Aircraft Vortices Visible to Radar by Spraying Water into the Wake

## 1. Introduction

### a. Motivation

Aircraft trailing vortices are a hazard to following airplanes during take-off and landing (e.g., Barbagallo 2014). For a review of their dynamics, see Spalart (1998). Due to their mutually induced velocity, a pair of trailing vortices generally descends from the altitude where it was generated. (For exceptions to this in strongly stratified conditions, see Spalart, 1996). Therefore, a wake encounter may occur when a following aircraft finds itself below^{1} the path of the leading aircraft. The possibility for this is increased during take-off and landing, and the present work addresses the latter situation. When cleared for a visual approach to landing, the pilot of the following aircraft can visually attempt to remain above the path of the leader by, for example, flying at a higher glide slope to the same touchdown point as the leader. Even in visual approaches, however, things can and do go awry (Barbagallo 2014, §2.3). This means that efforts to find an all-weather wake sensor should be continued.

When the ceiling is less than 1000 ft and the visibility less than 3 statute miles, operations must be conducted under Instrument Flight Regulations (IFR). In this case, air traffic controllers maintain separations according to the weight categories of the leading and following aircraft; see Table 1. These separations have started to limit capacity at some airports (Crouch et al. 2001b) and we refer the reader to a report (Broderick et al. 2008) by a committee of the National Research Council entitled “Wake Turbulence: An Obstacle to Increased Air Traffic Capacity.” About two decades ago, NASA initiated the AVOSS (Aircraft Vortex Spacing System) program whose aim was to make aircraft spacing dynamic through a combination of vortex sensing and real-time flowfield simulation. The present work is motivated by the vortex detection aspect of the AVOSS program. A current effort that has similar aims as the defunct AVOSS program is WakeNet-Europe (www.wakenet.eu). Finally, we mention the development of a system to hasten the breakup of trailing vortices by exciting vortex instabilities through periodic motion of control surfaces (Crouch et al. 2001a,b). Such a system could be implemented together with one for wake detection.

### TABLE 1

Follower | |||||
---|---|---|---|---|---|

Leader | Super | Heavy | B757 | Large | Small |

Super | 3 | 6 | 7 | 7 | 8 |

Heavy | 3 | 4 | 5 | 5 | 6 |

B757 | 3 | 4 | 4 | 4 | 5 |

Large | 3 | 3 | 3 | 3 | 4 |

Small | 3 | 3 | 3 | 3 | 3 |

### b. Previous work using radar to detect aircraft wakes

The main vortex detection technology tested by the AVOSS program was infrared lidar. A concern raised about lidar was that since water vapor strongly absorbs infrared, it would not be usable in IFR conditions. Another concern is that optical systems are more expensive and difficult to maintain than radar. These concerns motivate consideration of radar.

#### 1) Clear air reflectivity

The first observational and theoretical efforts on radar detection of wakes were in the context of *clear air*. Systematic tests by Gilson (1992) showed that the wake of a C-5A aircraft could be detected by radars having 2–7 MW of peak power at frequencies from 0.162 to 5.67 GHz, with no return at 35 GHz. Shariff and Wray (2002) analyzed the reflectivity in Gilson’s test using a model of a vortex pair descending in a stratified atmosphere, carrying with it an oval of atmospheric air from the altitude at which the wake formed. This leads to a gradient in refractive index between the oval and ambient air. Another mechanism investigated by Shariff and Wray (2002), which has peak reflectivity at a low frequency of 50 MHz, is the pressure (hence density) gradient in each vortex. Both of these mechanisms have some drawbacks for practical use: (a) the radar cross-section is small (−60 to −80 dB m^{); (b) the first mechanism depends on atmospheric stratification, which has seasonal, geographic, and diurnal variations; and (c) the frequency of the second mechanism is the same as that of Stratospheric-Tropospheric radars, which have an antenna array on the ground looking upward. For the present application, technology would have to be developed for aiming the radar by use of phasing. }Li et al. (2011) extended the work of Shariff and Wray (2002) with a better calculation of the compressibility-induced variation in air density in the wake. More importantly, they showed that the atmospheric gradient of water vapor is wound up into a spiral by each vortex, which allows scattering at high frequencies.

Babaresco et al. (2008) conducted tests using an X-band (9.6 GHz) radar for aircraft taking-off and flying straight and level at 1500 m altitude. For the take-off case, the range was about 700 m looking roughly sideways using peak pulse power of only 20 W. For the aircraft at an altitude of 1500 m, the radar was looking straight up and the peak pulse power was 75 W. A time doppler plot indicated a spiral structure in each vortex. However, in staring mode, a return is obtained for only five seconds, at most. Babaresco et al. (2008) did not report how far downstream of the aircraft the wake could be detected; it appears that this distance is very short.

In conclusion, clear air reflectivity is an interesting prospect, but more careful and better documented observational campaigns combined with theoretical efforts are needed.

#### 2) Exploiting natural precipitation

A different approach, which the author first heard about from Robert Neece (NASA Langley, personal communication) around 2000, exploits the fact that water droplets (in the form of fog or rain) are present in IFR weather conditions. They are strong radar reflectors, and, when an airplane wake sets these particles into motion, they can be separated from ambient droplets via their doppler signature. This strategy has the advantage that standard doppler weather radars could be used. Seliga and Mead (2009) demonstrated the feasibility of the approach using a W-band (94 GHz) radar in an opportunistic test with only 100 mW of peak power. An analysis of the reflectivity of this mechanism has been conducted by Liu et al. (2013), and further measurements are reported in Babaresco (2012). The main drawback of this approach is that there may be long stretches along the flight path where a sufficient number of natural droplets is unavailable.

### c. Present approach

The approach proposed here is related to the use of natural precipitation for wake detection discussed in the previous paragraph, except that, to provide a persistent radar target, water spray is injected into the wake. We envision one nozzle on each side of the airplane that injects the water spray near the wing trailing edge at a specified spanwise location. One possibility for nozzle placement is the aft tip of a flap track fairing. A pump and water tank could be located nearby within the wing.

The main constraint imposed by nature is droplet retention in the wake: a nozzle produces a distribution of droplet sizes; the smallest ones evaporate away while the largest ones fall out of the wake due to gravity. However, IFR conditions are correlated with high humidity (§2l). This reduces the rate of evaporation and makes the approach feasible. In a non-IFR case we consider that has moderate humidity, only the Ka and W band radars give sufficient reflectivity in a patch above each vortex. However, it is shown that spectral processing reduces noise and allows detection even when the signal to noise ratio in individual pulse returns is < 1. Finally, there is no reason why the present approach and that of using natural precipitation could not be combined using the same radar.

### a. Motivation

Aircraft trailing vortices are a hazard to following airplanes during take-off and landing (e.g., Barbagallo 2014). For a review of their dynamics, see Spalart (1998). Due to their mutually induced velocity, a pair of trailing vortices generally descends from the altitude where it was generated. (For exceptions to this in strongly stratified conditions, see Spalart, 1996). Therefore, a wake encounter may occur when a following aircraft finds itself below^{1} the path of the leading aircraft. The possibility for this is increased during take-off and landing, and the present work addresses the latter situation. When cleared for a visual approach to landing, the pilot of the following aircraft can visually attempt to remain above the path of the leader by, for example, flying at a higher glide slope to the same touchdown point as the leader. Even in visual approaches, however, things can and do go awry (Barbagallo 2014, §2.3). This means that efforts to find an all-weather wake sensor should be continued.

When the ceiling is less than 1000 ft and the visibility less than 3 statute miles, operations must be conducted under Instrument Flight Regulations (IFR). In this case, air traffic controllers maintain separations according to the weight categories of the leading and following aircraft; see Table 1. These separations have started to limit capacity at some airports (Crouch et al. 2001b) and we refer the reader to a report (Broderick et al. 2008) by a committee of the National Research Council entitled “Wake Turbulence: An Obstacle to Increased Air Traffic Capacity.” About two decades ago, NASA initiated the AVOSS (Aircraft Vortex Spacing System) program whose aim was to make aircraft spacing dynamic through a combination of vortex sensing and real-time flowfield simulation. The present work is motivated by the vortex detection aspect of the AVOSS program. A current effort that has similar aims as the defunct AVOSS program is WakeNet-Europe (www.wakenet.eu). Finally, we mention the development of a system to hasten the breakup of trailing vortices by exciting vortex instabilities through periodic motion of control surfaces (Crouch et al. 2001a,b). Such a system could be implemented together with one for wake detection.

### TABLE 1

Follower | |||||
---|---|---|---|---|---|

Leader | Super | Heavy | B757 | Large | Small |

Super | 3 | 6 | 7 | 7 | 8 |

Heavy | 3 | 4 | 5 | 5 | 6 |

B757 | 3 | 4 | 4 | 4 | 5 |

Large | 3 | 3 | 3 | 3 | 4 |

Small | 3 | 3 | 3 | 3 | 3 |

### b. Previous work using radar to detect aircraft wakes

The main vortex detection technology tested by the AVOSS program was infrared lidar. A concern raised about lidar was that since water vapor strongly absorbs infrared, it would not be usable in IFR conditions. Another concern is that optical systems are more expensive and difficult to maintain than radar. These concerns motivate consideration of radar.

#### 1) Clear air reflectivity

The first observational and theoretical efforts on radar detection of wakes were in the context of *clear air*. Systematic tests by Gilson (1992) showed that the wake of a C-5A aircraft could be detected by radars having 2–7 MW of peak power at frequencies from 0.162 to 5.67 GHz, with no return at 35 GHz. Shariff and Wray (2002) analyzed the reflectivity in Gilson’s test using a model of a vortex pair descending in a stratified atmosphere, carrying with it an oval of atmospheric air from the altitude at which the wake formed. This leads to a gradient in refractive index between the oval and ambient air. Another mechanism investigated by Shariff and Wray (2002), which has peak reflectivity at a low frequency of 50 MHz, is the pressure (hence density) gradient in each vortex. Both of these mechanisms have some drawbacks for practical use: (a) the radar cross-section is small (−60 to −80 dB m^{); (b) the first mechanism depends on atmospheric stratification, which has seasonal, geographic, and diurnal variations; and (c) the frequency of the second mechanism is the same as that of Stratospheric-Tropospheric radars, which have an antenna array on the ground looking upward. For the present application, technology would have to be developed for aiming the radar by use of phasing. }Li et al. (2011) extended the work of Shariff and Wray (2002) with a better calculation of the compressibility-induced variation in air density in the wake. More importantly, they showed that the atmospheric gradient of water vapor is wound up into a spiral by each vortex, which allows scattering at high frequencies.

Babaresco et al. (2008) conducted tests using an X-band (9.6 GHz) radar for aircraft taking-off and flying straight and level at 1500 m altitude. For the take-off case, the range was about 700 m looking roughly sideways using peak pulse power of only 20 W. For the aircraft at an altitude of 1500 m, the radar was looking straight up and the peak pulse power was 75 W. A time doppler plot indicated a spiral structure in each vortex. However, in staring mode, a return is obtained for only five seconds, at most. Babaresco et al. (2008) did not report how far downstream of the aircraft the wake could be detected; it appears that this distance is very short.

In conclusion, clear air reflectivity is an interesting prospect, but more careful and better documented observational campaigns combined with theoretical efforts are needed.

#### 2) Exploiting natural precipitation

A different approach, which the author first heard about from Robert Neece (NASA Langley, personal communication) around 2000, exploits the fact that water droplets (in the form of fog or rain) are present in IFR weather conditions. They are strong radar reflectors, and, when an airplane wake sets these particles into motion, they can be separated from ambient droplets via their doppler signature. This strategy has the advantage that standard doppler weather radars could be used. Seliga and Mead (2009) demonstrated the feasibility of the approach using a W-band (94 GHz) radar in an opportunistic test with only 100 mW of peak power. An analysis of the reflectivity of this mechanism has been conducted by Liu et al. (2013), and further measurements are reported in Babaresco (2012). The main drawback of this approach is that there may be long stretches along the flight path where a sufficient number of natural droplets is unavailable.

#### 1) Clear air reflectivity

The first observational and theoretical efforts on radar detection of wakes were in the context of *clear air*. Systematic tests by Gilson (1992) showed that the wake of a C-5A aircraft could be detected by radars having 2–7 MW of peak power at frequencies from 0.162 to 5.67 GHz, with no return at 35 GHz. Shariff and Wray (2002) analyzed the reflectivity in Gilson’s test using a model of a vortex pair descending in a stratified atmosphere, carrying with it an oval of atmospheric air from the altitude at which the wake formed. This leads to a gradient in refractive index between the oval and ambient air. Another mechanism investigated by Shariff and Wray (2002), which has peak reflectivity at a low frequency of 50 MHz, is the pressure (hence density) gradient in each vortex. Both of these mechanisms have some drawbacks for practical use: (a) the radar cross-section is small (−60 to −80 dB m^{); (b) the first mechanism depends on atmospheric stratification, which has seasonal, geographic, and diurnal variations; and (c) the frequency of the second mechanism is the same as that of Stratospheric-Tropospheric radars, which have an antenna array on the ground looking upward. For the present application, technology would have to be developed for aiming the radar by use of phasing. }Li et al. (2011) extended the work of Shariff and Wray (2002) with a better calculation of the compressibility-induced variation in air density in the wake. More importantly, they showed that the atmospheric gradient of water vapor is wound up into a spiral by each vortex, which allows scattering at high frequencies.

Babaresco et al. (2008) conducted tests using an X-band (9.6 GHz) radar for aircraft taking-off and flying straight and level at 1500 m altitude. For the take-off case, the range was about 700 m looking roughly sideways using peak pulse power of only 20 W. For the aircraft at an altitude of 1500 m, the radar was looking straight up and the peak pulse power was 75 W. A time doppler plot indicated a spiral structure in each vortex. However, in staring mode, a return is obtained for only five seconds, at most. Babaresco et al. (2008) did not report how far downstream of the aircraft the wake could be detected; it appears that this distance is very short.

In conclusion, clear air reflectivity is an interesting prospect, but more careful and better documented observational campaigns combined with theoretical efforts are needed.

#### 2) Exploiting natural precipitation

A different approach, which the author first heard about from Robert Neece (NASA Langley, personal communication) around 2000, exploits the fact that water droplets (in the form of fog or rain) are present in IFR weather conditions. They are strong radar reflectors, and, when an airplane wake sets these particles into motion, they can be separated from ambient droplets via their doppler signature. This strategy has the advantage that standard doppler weather radars could be used. Seliga and Mead (2009) demonstrated the feasibility of the approach using a W-band (94 GHz) radar in an opportunistic test with only 100 mW of peak power. An analysis of the reflectivity of this mechanism has been conducted by Liu et al. (2013), and further measurements are reported in Babaresco (2012). The main drawback of this approach is that there may be long stretches along the flight path where a sufficient number of natural droplets is unavailable.

### c. Present approach

The approach proposed here is related to the use of natural precipitation for wake detection discussed in the previous paragraph, except that, to provide a persistent radar target, water spray is injected into the wake. We envision one nozzle on each side of the airplane that injects the water spray near the wing trailing edge at a specified spanwise location. One possibility for nozzle placement is the aft tip of a flap track fairing. A pump and water tank could be located nearby within the wing.

The main constraint imposed by nature is droplet retention in the wake: a nozzle produces a distribution of droplet sizes; the smallest ones evaporate away while the largest ones fall out of the wake due to gravity. However, IFR conditions are correlated with high humidity (§2l). This reduces the rate of evaporation and makes the approach feasible. In a non-IFR case we consider that has moderate humidity, only the Ka and W band radars give sufficient reflectivity in a patch above each vortex. However, it is shown that spectral processing reduces noise and allows detection even when the signal to noise ratio in individual pulse returns is < 1. Finally, there is no reason why the present approach and that of using natural precipitation could not be combined using the same radar.

## 2. Calculation methods

### a. General procedure

We use aircraft centered coordinates (*x,y,z*), where *x* is streamwise, *y* is spanwise, and *z* is vertical. Air and water are denoted by subscripts ‘a’ and ‘w’, respectively.

The calculation has three steps. The first creates a sample of droplet radii *a*(*t*) from a size distribution pertinent to aerial spray nozzles. The second step generates a spray trail on the starboard side of the aircraft wake. This involves injecting droplets from the sample into the wake, tracking the position **X**(*t*) of each droplet and its radius *a*(*t*) as it evaporates. The flow-field model consists of two counter-rotating vortices whose height decreases with downstream distance *x* behind the aircraft. The third and final step mirrors the starboard trail to the port side and, for a given set of radar parameters, computes the reflectivity of both trails together.

Since it is prohibitive to track all of the actual droplets, the spray trail computed in step two consists of a certain number, *N*_{comp}, of computational droplets. In the reflectivity calculation, each computational droplet is taken to represent a multiplicity, *M*_{true}, of actual droplets, which is simply the ratio of the desired injected volume to the volume injected in the computation.

The equations for droplet motion and evaporation are integrated using the routine LSODE which is described in Radhakrishnan and Hindmarsh (1993) and available from NETLIB. The routine chooses a time step based on a specified error tolerance. Since the flow-field model we have adopted is steady in a reference frame moving with the aircraft, it is not necessary to inject droplets at every time step, which would continuously increase the number of computational droplets. Instead, droplets are injected only during a certain time interval Δ*t*_{0} at the beginning of the calculation. This set of droplets is then evolved for successive Δ*t*_{0} intervals and appended to a file at the end of every Δ*t*_{0} interval. At the end of the computation, the file contains a spray trail that is 7 nm miles long, whose radar reflectivity we subsequently analyze. In actual practice, an aircraft would likely not generate a trail of several nautical miles behind it. Rather, it would release the spray for a short period at pre-selected locations during its approach.

A detailed description of each part of the procedure is described in the following subsections. The evaporation calculation is described in Appendix A2. To avoid stiffness of the system of evolution equations, droplets are removed from the calculation when their radius becomes < 20 *μ*m; at this point their reflectivity is too small to significantly affect the received power.

### b. Droplet trajectory

The position **X**(*t*) and velocity **U**(*t*) of a droplet of mass *m*_{d} evolves according to

where *g*_{eff} = (1 − ρ_{a}*/ρ*_{w})*g* is the effective gravity accounting for buoyancy. The drag force **F**_{D} is given by

where

**u**

_{rel}=

**u**(

**X**) −

**U**,

is the velocity of the air flow, **u**(**X**), relative to the droplet. Evaluation of the drag coefficient, *C*_{D}, is described in Appendix A1.

### c. Droplet size distribution of aerial spray nozzles

When the water jet issues from the nozzle, it will encounter a blast of free-stream air with a speed of 77.2 m s^{−}^{. A comparable situation in the literature is that of a cylindrical liquid jet surrounded by an annulus of co-flowing air (}Lorenzetto and Lefebvre 1977;Varga et al. 2003) which shows that much smaller droplets are produced than for the case of still air (with the same velocity of the liquid jet). This is due to the occurrence of both the Kelvin-Helmholtz instability and Rayleigh-Taylor instabilities. The former is driven by the shear between the air and liquid flow and leads to a smaller instability wavelength. The latter arises due to the acceleration of liquid droplets by the drag force of the air stream. If the droplets are too small, they quickly evaporate. To produce larger droplets, both instabilities can be mitigated by reducing the relative velocity between the liquid jet and free-stream air. This is accomplished by increasing the driving pressure. However, if the droplets are too large, they fall out of the wake. Hence, there is an optimum droplet size.

We were fortunate that an experimental study (Fritz and Hoffmann 2015), which uses a wind-tunnel to mimic the free stream air flow, has recently been performed to characterize the droplet size distribution produced by nozzles used for aerial agricultural spraying. This study did indeed show that larger pressures produce larger droplets. Dr. B. Fritz kindly provided us with an Excel program, developed from that study, which gives parameters of the droplet size distribution for various nozzles, free-stream speeds, and driving pressures. Use of these parameters is now described.

Let *p*(*a*) be the probability density such that *p*(*a*)*da* is the probability that the droplet radius is in the interval [*a, a* + *da*]. The log-normal distribution

with parameters *a*_{0} and *σ*, is commonly used in the spray literature. Fritz’s Excel program provides information about the function *Q*(*a*), defined to be the fractional volume occupied by droplets of radius ≤ *a*. One can show from appropriately integrating (5) that

*Q*(

*a*) = ½(1 + erf

*ξ*),

where

The Excel program provides *a*_{0}_{.}_{5} and *a*_{0}_{.}_{9} defined such that *Q*(*a*_{0}_{.}_{5}) = 0.5 and *Q*(*a*_{0}_{.}_{9}) = 0.9. Using them, (6) can be numerically inverted to yield the parameters *a*_{0} and *σ* of the log-normal distribution (5). The Excel program also provides *a*_{0}_{.}_{1}, which we did not use because we wished to nail the size distribution for large droplets which contribute most to reflectivity.

### d. Choice of nozzle

What is the best drop size distribution? A set of *N*_{d} droplets over which the incident beam is assumed to have uniform intensity, has a reflectivity proportional to (e.g., Doviak and Zrnić 1984, p. 58)

assuming Rayleigh scattering. Maximizing this subject to fixed volume of water and fixed *N*_{d} gives the result that all droplets must be of the same size. Given this result, maximizing *ζ* with respect to the number *N*_{d} subject to fixed volume gives *N*_{d} = 1, i.e., all the volume must be in one droplet. However, such a droplet would likely fall too rapidly. To minimize droplet loss by sedimentation, droplets must not have a terminal velocity larger than the vortex descent speed, *W*_{descent} = 1.75 m s^{−}^{ in the present case. Consulting the terminal velocity plot in }Pruppacher and Klett (1997, p. 416) we conclude that *a* must be ≤ 200 *μ*m. Hence the best distribution is uniform with a drop radius *a* = 200 *μ*m. Since the rate of droplet evaporation is ∝ 1*/a*, i.e., small droplets evaporate faster than larger ones, the above conclusion is not altered by including evaporation.

The above considerations suggest that the following “rate of *ζ*” could be used to initially evaluate different nozzles without having to perform an evaporation and reflectivity calculation:

where the < superscript denotes the exclusion from the sum of droplets that fall away and *N*_{d}(Δ*t*) denotes the number of droplets produced by the nozzle in a period Δ*t*. This suggestion will be tested in §3h.

Table 2 lists the parameters for four candidate nozzles and operating conditions that were selected from B. Fritz’s Excel program based on having peak probability density at a large radius (which rarely exceeded 100 *μ*m) and a high flow-rate. The flow-rate for nozzle 1 was provided by Calvin Kroes (private communication) of CP Products. For nozzles 2 and 3, flow-rates were extrapolated from the values of 5.30 gpm and 2.45 gpm, respectively, at 60 psi reported on the manufacturer’s data sheets (www.translandllc.com/wp-content/uploads/2015/08/Aerial-Flow-Chart-20152.pdf and www.cpproductsinc.com/images/stories/downloads/Misc-Tables/A1-Web%20Aerial%20Tip%20Rate%20Chart.pdf). The extrapolations assumed a square-root dependence of flow-rate on pressure (Lefebvre 1989, p. 157) expected from Bernoulli’s principle, and which the manufacturer’s data follows well. For nozzle 4, the flow-rate versus pressure provided on the manufacturer’s website has a linear rather than square-root dependence. We are grateful to Dr. Brad Fritz (USDA) for measuring the actual flow-rate for us. It turned out to be much lower than the value provided by on the website.

### TABLE 2

Nozzle 1 | Nozzle 2 | Nozzle 3 | Nozzle 4 | |
---|---|---|---|---|

Model | CP-09 | CP-09 | CP11TT | Davidon-Triset |

Pressure (psi) | 90 | 90 | 90 | 90 |

Airspeed (mph) | 175 | 175 | 175 | 175 |

Deflection-plane/body angle | 0° | 0° | 0° | 0° |

Fan angle | Str.^{1} | Str. | Str. | Str. |

Orifice Code | 20 | |||

Orifice diameter (in) | 0.125 | 0.172 | 0.105 | 0.125 |

a_{0 5} (μm)^{Note 2} | 178.5 | 146.5 | 183 | 239.5 |

a_{0 9} (μm)^{Note 2} | 315 | 303 | 359 | 463 |

a_{0} of log-normal | 99.02 | 55.83 | 79.85 | 108.3 |

σ of log-normal | 0.443 | 0.567 | 0.526 | 0.514 |

Flow-rate (gpm)^{Note 3} | 3.70 | 6.49 | 3.00 | 3.06 |

U_{exit} (m s^{−}^{)} | 29.5 | 27.3 | 33.9 | 24.4 |

No. of nozzles per side | 1 | 1 | 1 | 1 |

Gallons per nm (two sides) | 2.96 | 5.19 | 2.4 | 2.45 |

### e. Droplet injection in the computation

Let the origin of coordinates be in the symmetry plane of the aircraft (which corresponds to *y* = 0) with same axial and vertical location as the droplet injector. Droplets are injected in a grid pattern within a square of width *w*_{square} in the *yz*-plane. The pattern consists of *n*_{square} × *n*_{square} droplets. The square is centered at (*x,y,z*) = (0*, fb/*2,0), where *f* represents the fractional spanwise distance from the aircraft center plane to the wing-tip, and was chosen to be *f* = 0.5. A square shape was chosen because the nozzles we have selected are not of the flat fan type. The streamwise extent of the computed mist trail was chosen to be *ℓ*_{trail} = 7 nm (168 seconds of elapsed time from injection) since we wish to detect the trail at 6 nm, the longest distance for which one would want to detect the wake of a heavy aircraft under current separation rules; see Table 1.

The square pattern is injected *n _{x}* times and the time interval between injections is Δ

*t*

_{inject}, whose value is is chosen so that the

*x*spacing between droplets is the same as in the cross-sectional (

*yz*) plane. After a time period Δ

*t*Δ

_{0}≡ n_{x}*t*

_{inject}, an

*n*

^{2}_{square}×

*n*slab of particles has been injected, which is then advanced for successive Δ

_{x}*t*

_{0}periods to form the entire trail. Inertial particles with a small Stokes number tend to an attractor (Haller and Sapsis 2008;Sapsis and Haller 2010) independent of injection location, and therefore where droplets end-up should be insensitive to where they are injected. A brief check on insensitivity to initial conditions will be presented in §3b.

The initial velocity of droplets is set equal to the air velocity, which is justified as follows. From the equations of droplet motion, (2)–(4), the characteristic relaxation time for a droplet to start following a new air speed, imposed at *t* = 0, say, is

where the subscript ‘0’ means that the quantity is calculated at *t* = 0. From this, the characteristic relaxation distance, *ℓ*_{relax} = *u*_{rel}(0)*τ*_{relax} can be evaluated. Note that the initial air speed relative to the water jet is given by *u*_{rel}(0) = *U*_{app} −*U*_{exit}, where *U*_{app} is the approach velocity of the aircraft given in Table 3, and *U*_{exit} is the exit velocity of the water jet given in Table 2. The experiment of Fritz and Hoffmann (2015) measured the size distribution 1.8 m downstream of the nozzle for all straight-stream nozzles. This value is from a private communication from B. Fritz and represents a correction from a value of 1.5 m reported in Fritz and Hoffmann (2015). Figure 1b plots *ℓ*_{relax} as a function of drop radius for nozzle 1. Inspecting it together with the size distribution in Figure 1a, we conclude that most of the droplets are following the air stream at the measurement station of the experiment. If this had not been the case and there had been a relative velocity large enough to give Weber numbers ≳ 10, then it would have been necessary to model further droplet break-up using a secondary break-up model (e.g., Apte et al. 2003).

### TABLE 3

Parameter | Value |
---|---|

Weight, W | 500,000 lb |

Wing span, b | 60 m |

Vortex spacing, b_{0} | 47.9 m |

Vortex circulation, Γ | 526 m^{ s}^{−}^{1} |

Approach speed, U_{app} | 150 knots |

### f. Flowfield of two counter-rotating wing-tip vortices

We consider an aircraft flying straight and level at altitude *z*_{0} = 0 in aircraft coordinates. The velocity field of the airflow in the wake is denoted by lower case **u**(**x**). This velocity field consists of the free stream, ${U}_{\text{app}}\widehat{\mathbf{x}}$(where *U*_{app} is the approach speed of the aircraft), superposed with the flow induced by a pair of counter-rotating vortices with circulations *±*Γ. The centerline of each vortex is at spanwise location *y*_{vort} = *±b*_{0}. Due do their mutually induced velocity, the height *z*_{vort}(*x*) of the vortex pair decreases with distance *x* behind the wing as follows:

*z*

_{vort}(

*x*) =

*z*

_{0}−

*W*

_{desct},

where *W*_{desc} = Γ*/*2*πb*_{0} is the descent speed of the vortex pair. The quantity *t* = *x/U*_{app} is time since the vortex at *x* was shed from the wing.

Each vortex induces a circumferential velocity *u _{θ}*(

*r*) in the cross-plane (

*yz*). For

*u*(

_{θ}*r*), a profile fit to flight data by P. Spalart (Private communication) of Boeing is used:

where *η ≡ r/b*_{0}.

For an elliptically loaded wing, lifting line theory (Batchelor 1967) gives the vortex spacing as

where *b* is the wingspan, and the vortex circulation as

where *W* is the aircraft weight. We use the parameters for a typical heavy aircraft given in Table 3.

### g. Calculation of the received signal and power

A pulse-doppler radar transmits a train of square wave pulses that modulate a carrier wave of frequency *f* = 2*π/ω*. The duration of each pulse is *τ* and the pulse repetition frequency is *f*_{PRF}. After each pulse is transmitted, the transmit-receive switch is set to the receive position and the incoming signal is sampled. Each sample at time *t* is said to come from the *range gate r* = *c*(*t − t*_{t})*/*2, where *t*_{t} is the transmit time of the pulse and the factor of two accounts for the round-trip. Throughout, we consider the case where the transmitting and receiving antenna are the same, the so-called mono-static case.

#### 1) Resolution shell

An important concept is that of the *resolution volume*, ℛ_{τ}(*t*) *R _{τ}*(

*t*), at time

*t*, associated with a single pulse of finite duration,

*τ*(e.g., Yuter 2003, p. 1836). It is defined as the volume from which a signal is received

*at a fixed time t*due to scattering by the pulse. Let

*t*

_{t}mark the beginning of the pulse at the transmitter. A signal from a scatterer at distance

*r*will be received in the time interval

*t*−

*t*

_{t}= 2

*r*/

*c*+

*ξ*

*τ*, 0 ≤

*ξ*≤ 1,

where *ξ* = 0 corresponds to leading edge of the pulse and *ξ* = 1 to its trailing edge. Solving (15) for *r* gives

*r*=

*c*(

*t*−

*t*

_{t}−

*ξ*

*τ*)/2, 0 ≤

*ξ*≤ 1.

Equation (16) defines a spherical shell (called the resolution volume) from which a signal is received at the fixed time *t*. The next subsection describes how the complex voltage received at a given time is evaluated by summing the complex voltages from each droplet in the resolution volume.

#### 2) Received power and signal-to-noise ratio

The material in this subsection is adapted from the texts Doviak and Zrnić (1984) and Ishimaru (1978). Let the “voltage” of the transmitted pulse be (the real part of)

where *V*_{0} is a complex amplitude, *τ* is the pulse width, and “voltage” is defined such that the instantaneous transmitted power is *P*_{t}(*t*) = *V*_{t}(*t*)*V*_{t}* ^{∗}*(

*t*). The voltage at the input terminals of the receiver is (the real part of) the following summation over droplets:

where *A _{m}* is a complex scattering amplitude and

*ω*

_{m}≡

*ω*(1 − 2

*u*

_{m}/

*c*)

is the twice doppler-shifted frequency (*u _{m}* being the radial velocity of the

*m*th droplet) and

*r*is the distance to each droplet. The summation in (18) is taken over the

_{m}*N*

_{d}droplets in the resolution volume associated with the pulse.

The complex scattering amplitude due to each droplet is

The first factor in (20), namely [.]^{, is copied from the square root of the radar equation (e.g., }Doviak and Zrnić 1984, p. 34) for power, where *λ* is the wavelength, and *σ*_{b}* _{m}* is the back-scattering cross-section of each droplet. The function

*G*(

*θ*) is the gain function of the antenna at the droplet, which we have assumed to depend on its angle

_{m}*θ*from the beam centerline. Two loss factors (< 1) have been included in (20):

_{m}*ℓ*

_{w}is the two-way waveguide loss and

*ℓ*

_{B}is the loss due to finite bandwidth of the receiver; values assumed for the present study are given in Table 4. The last factor in (20), exp(i

*ϕ*), accounts for the phase shift induced by back-scattering at the

_{m}*m*th droplet. The back-scattering cross-section and phase-shift (

*σ*

_{b}

*and*

_{m}*φ*, respectively) will be discussed further in §2k.

_{m}### TABLE 4

Manufacturer | EEC | EEC | EEC | Metek | ProSensing |
---|---|---|---|---|---|

Series | DWSR | DWSR | DWSR | ||

Model | 8501S | nnn1C^{Note 1} | 2001X | MIRA-35 | W-SACR |

Frequency, f (GHz) | 3 | 5.9 | 9.6 | 35.1 | 93.9 |

Peak power, P_{t} (kW) | 850 | 250–1000 | 200 | 30 | 1.7 |

Reflector diameter (m) | 4.2 | 4.2 | 2.4 | 1.2,2.0 | 0.9 |

½-power beam width, θ_{b} | 1.83° | 0.95° | 0.95° | 0.52°,0.31° | 0.30° |

Antenna gain, G (dB) | 39.5 | 45 | 45 | 50.4,53.5 | 54.5 |

Pulse width, τ (μs) | 0.4–2 | 0.2–3 | 0.2–2 | 0.1,0.2,0.4 | 0.05–2 |

Range resolution, cτ/2 (m) | 60–300 | 30–450 | 30–300 | 15,30,60 | 7.5–300 |

PRF (kHz) | 0.2–2.4 | 0.2–2.4 | 0.2–2.4 | 2.5,5,10 | ≤ 20 |

U_{max} = c PRF/4 f (m s^{−}^{)} | 5–60 | 2.5–31 | 1.6–19 | 5.3,11,21 | ≤ 16 |

Receiver noise figure (dB) | 2.0 | 2.0 | 2.0 | 6.2 | 6.0 |

2-way waveguide loss^{ (dB)} | 0.8 | 0.8 | 0.8 | 0.8 | 0.8 |

Finite bandwidth loss^{ (dB)} | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 |

Rain attenuation^{ (dB/km)} | 0.005 | 0.03 | 0.12 | 2 | 7 |

Minimum range (m) | 150^{Note 5} |

Radar receivers have electronics that can obtain the real and imaginary parts (denoted *I* and *Q*) of *V*(*t*) and compute the instantaneous received power

*P*

_{r}(

*t*) =

*I*(

*t*)

^{ + }

*Q*(

*t*)

^{ = }

*V*(

*t*)

*V*

^{(}

*t*).

We note in passing that *I*(*t*) and *Q*(*t*) are the components of the real part of *V*(*t*) that are in-phase and 90*°* out-of-phase with the transmitted carrier, respectively. Substituting (18) into (21) and splitting the sum into two parts, following Doviak and Zrnić (1984, §4,1), gives

where *k _{m} ≡ ω_{m}/c*. Arguments for using only the

*first term*in (23) in order to evaluate reflectivity, i.e., for summing the powers received from each droplet, are given by Rayleigh (1945, p. 37), Beckmann (1962), and Doviak and Zrnić (1984, §4,1). The important point, which was phrased eloquently by Rayleigh, is that it is not correct to say that the power in a single return from a random distribution of droplets is the sum of the powers scattered by each. Rather, the result is true only when a

*large ensemble of returns*from a statistically stationary target are averaged. This is most easily seen when we consider the case when all the

*A*are equal (to unity, say). Then, the first term

_{m}*T*

_{1}in (23) is

*T*

_{1}=

*N*

_{d}. If droplet distances are randomly distributed in the resolution shell (assumed to be several wavelengths wide), then the magnitude of the second term will be the average of the summands times the number of terms, i.e., ${T}_{2}\approx {N}_{t}^{-1/2}\times {N}_{t}={N}_{t}^{1/2}\approx {N}_{d}$, where

*N*

_{t}=

*N*

_{d}(

*N*

_{d}−1) is the number of terms in the double sum. Hence both terms in (23) are of

*similar magnitude*. To make the second term smaller than the first, an ensemble average must be taken over many pulses. If the ensemble has

*N*

_{s}phase-uncorrelated samples, then the second term will be ${N}_{\mathrm{S}}^{1/2}$

smaller than the first.

In the present application, the wake descends through the beam and so the target is not stationary, strictly speaking. In §3c we will explicitly verify that, for our case, an average of over a certain number of pulses does indeed yield the first term in (23). Note that a radar set does not have direct access to the first term in (23); only we as simulators do.

To provide a measure of detectability, we will present the signal-to-noise ratio

*P*

_{r1}/

*P*

_{noise},

where *P*_{r1} is the first term in (23) and the average noise power is

*P*

_{noise}=

*k*

_{B}

*T*

_{0}

*F*

_{N}/

*τ*,

where *k*_{B} = 1.381 × 10^{−}^{ J K}^{−}^{ is Boltzmann’s constant, }*T*_{0} = 290 K is a reference temperature set by convention, and *F*_{N} is the overall noise figure of the chain of components in the receiving cascade. Values for *F*_{N} and *τ* are listed in Table 4 for each radar set.

Finally, since each computational droplet represents a multiplicity *M*_{actual} of actual droplets, we have

where *N*_{comp} is the number of computational droplets in the resolution volume. Note that all *M*_{actual} copies of each computational droplet are assumed to be at the same location and therefore their scattered voltages at the receiver add constructively. This assumption does not bias SNR1 since its calculation involves summing individual scattered powers anyway. We claim that the statistics of individual pulse returns are also not affected by this assumption; this will be verified (§3d) in a computational test where the number of computational droplets is increased.

### h. Calculation of the doppler spectrum

Radars calculate a doppler spectrum for a given spatial observation location by performing a Fast Fourier Transform (FFT) of complex voltage returns (at the same range gate) from a sequence of pulses separated by Δ*t*_{pulses} = 1/PRF, where PRF is the pulse repetition frequency. We shall do the same for the simulations. From a series of pulse returns, *V _{n}*,

*n*= 0,1

*,…N*

_{FFT}− 1, the normalized transform

and then the power spectrum $S\left(k\right)\equiv \widehat{V}\left(k\right)\widehat{V}\ast \left(k\right)$is computed. Note that the frequency index *k* corresponds to an actual frequency

The frequency *k*_{actual} is in units of (the period of the sequence)^{−}^{ = (}*N*_{FFT}Δ*t*_{pulses})^{−}^{, so in units of s}^{−}^{}

*f*(

*k*

_{actual}) = PRF

*k*

_{actual}/

*N*

_{FFT}.

Equating this to −2(*u*_{doppler}*/c*) *f* gives the doppler velocity associated with each *k*_{actual}. Finally, we state that we use the Hamming window (e.g., Harris 1978).

### i. Antenna gain function

We assume a Gaussian beam with transmitted power flux (power per unit area) vector

where *A* is a coefficient, *θ* is the angle from the beam centerline, and

in terms of the half-power full-width, *θ*_{b}. The total power crossing a sphere of radius *r* is

for a narrow beam. Using the definition of the gain function, *G*(*θ*), we obtain

### j. Radars included in the study

Table 4 lists parameters of currently operational doppler weather/cloud radars considered in the present study. The S, C, and X-band radars chosen are the DWSR series manufactured by EEC (Enterprise Electronics Corporation, Enterprise, Alabama). ARC (Advanced Radar Corporation, Boulder, Co.) makes quite similar C and X-band radars, while Baron Services (Huntsville, Al.) makes similar S, C, and X-band radars. The power and beamwidth values of the C-Band TDWR (Terminal Doppler Weather Radar) deployed at many US airports is subsumed by the range of values provided by the EEC C-band radar, and is therefore not included here.

A number of descriptions of Ka-band (35 GHz) cloud radars have appeared in the literature (Hamazu et al. 2003; Görsdorf et al. 2015). In the present work, we use parameters of the MIRA-35 radar manufactured by Metek (Elmshorn, Germany) which is described in Görsdorf et al. (2015). This choice was motivated by its relatively high power (30 kW). Other Ka-band weather radars, operational at the time of writing are: (i) Scanning 2 kW radars operated by U.S. Department of Energy’s Atmospheric Radiation Measurement (ARM) Climate Research Facility (Widener et al. 2012). (ii) The Copernicus 1 kW radar at Chilbolton Observatory (UK). (iii) An airborne 25 kW multi-frequency (X, Ka, and W-band) radar developed by Prosensing that is being used by NASA’s Langley Research Center for research into the detection and avoidance of super-cooled water droplets.

The science and technology of W-band (94 GHz) radars for cloud and precipitation research is reviewed in Kollias et al. (2007). For the present work we chose the W-SACR radar, which has been developed by the U.S. Department of Energy’s Atmospheric Radiation Measurement (ARM) program (Widener et al. 2012;Kollias et al. 2014).

With increasing frequency, *f*, reflectivity increases as *f*^{ in Rayleigh’s formula (ignoring Mie-scattering corrections). Furthermore, the size of the antenna required to obtain the same beamwidth is reduced. The main drawback of high frequency is increased attenuation due to precipitation between the radar and target. For example, the last entry in }Table 4 gives the attenuation rate in medium rain at W-band as 7 dB/km. A compensating factor is that when there is precipitation, the ambient humidity is also very high and so there is minimal evaporation, and, if natural precipitation is present in the wake, it will also contribute to reflectivity.

### k. Mie cross-section and phase-shift

Since we have can rather large droplets in the present application and frequencies up to 94 GHz, the back-scattering cross-section and phase-shift are obtained using Mie’s formula instead of Rayleigh’s approximation. We used subroutine BHMIE, available from Prof. B.T. Draine’s website at Princeton University, and checked the results using subroutine MIEV0 developed by Dr. W.J. Wiscombe (NASA Goddard).

Some understanding of notation is required to properly use these routines. Let the incident field be of unit magnitude and polarized in the 2-direction (defined to be perpendicular to the plane containing the incident and observer directions). For a spherical target, the scattered field in the far-field is also polarized in the 2-direction and is given by

where *f*_{22} is complex, *k* = 2*π/λ*, *ϑ* is the angle of the observer relative to the direction of propagation of the incident wave, and *r* is distance from the center of the sphere. The backscattering cross-section and phase shift are obtained as

*ϕ*= arg(

*f*

_{22}(

*π*)).

At the start of the reflectivity calculation at a given frequency, we tabulate the ratio *σ*_{b}*/σ*_{b}_{,}_{Rayleigh} and the difference *φ −φ*_{Rayleigh} as a function of droplet radius *a*. Rayleigh’s formulas are (Ishimaru 1978, p. 19)

*σ*

_{b}, Rayleigh = 4|

*K*

_{ε}|

^{(}

*k*

*a*)

^{(}

*π*

*a*

^{),}

*ϕ*

_{Rayleigh}= arg(

*K*

_{ε}),

where *K _{ε}* (a complex number) is given by

The quantity *ε*(*f, T*) is the complex dielectric constant of water; our convention of e^{ωt} for the time dependence requires the imaginary part of *ε* be positive for an absorbing material. It is a function of frequency and temperature and was evaluated using the single Debye model of Liebe et al. (1991) as implemented in subroutines available from Prof. Chris O’Dell’s website at Colorado State University. Since the droplet temperature is almost the same for all drops (the spread was 4 C at most), the dielectric constant *ε* is evaluated at the average temperature of all the drops in the trail.

Figure 2 displays the Mie back-scattering cross-section *σ*_{b} (normalized by the Rayleigh value) and phase-shift *φ* as a function of droplet radius at the five frequencies considered in this work. At the largest radii in this work, *a ≈* 600 *μ*m, the error in using Rayleigh’s cross-section is about 20% at 35.1 GHz.

### l. Choice of ambient temperature and humidity

Droplet evaporation calculations require specification of the ambient temperature and humidity. For guidance on appropriate choices, METARs (Meteorological Aerodrome Reports) during 2000–2014 were downloaded from

and processed for the five busiest airports in the U.S. Figures 3a and b show the monthly-averaged temperature and relative humidity (RH), respectively, when IFR conditions prevailed. Figure 3c shows the percentage of reports that fall into the IFR category. One sees that the average RH is always above 90%. To further synthesize this data, yearly averages were taken (Table 5). Among the five airports, LAX has the highest rate of evaporation in IFR conditions on average since it is the warmest and driest on average. Our choice is the IFR average for LAX, namely, *T* = 15.2 C and RH = 92.7%. Looking at the monthly data for the other four airports, this appears to be a reasonable choice for them also: it is an approximate lower bound for their monthly RH and their temperature is higher only during the summer months when IFR reports are low.

### TABLE 5

Airport | T (C) | RH | % IFR Reports |
---|---|---|---|

ATL | 12.5 | 95.5% | 5.3% |

LAX | 15.2 | 92.7% | 3.5% |

DFW | 9.7 | 94.5% | 1.4% |

ORD | 5.3 | 94.2% | 3.7% |

JFK | 12.1 | 93.9% | 4.8% |

Since it is expensive for flight tests to wait for IFR conditions to occur, and it is desirable to have a wake sensor that can work in a wider variety of atmosphere conditions, we will also consider a case of lower humidity and higher temperature, namely, RH = 60% and *T* = 20 C.

### m. Radar placement with respect to the wake

Here, we choose the radar location (*x*_{rad}*,y*_{rad}*,z*_{rad}) in aircraft-centered coordinates. Based on current wake separations (Table 1), it should not be necessary to examine a wake more than 6 nm behind the aircraft. We therefore chose to present results for the reflectivity at *x* = 6 nm, the worst case for droplet loss by evaporation and sedimentation. The radar is also placed at *x*_{rad} = 6 nm so it can view the *x* = 6 nm wake cross-section at normal incidence. Next, we assume that the aircraft is at the touchdown point. At 6 nm from the touchdown point, the altitude of an aircraft flying a 3*°* glide slope is *H* = 582 m. Therefore, the vertical coordinate of the radar is *z*_{rad} = −582 m. For purposes of this study, we assume an aircraft flying straight and level at this altitude. Initial flight tests would also presumably have the aircraft fly straight and level. In this case, since the vortex pair descends at a speed of *W*_{desc} = 1.75 m s^{−}^{, its axis makes a downward angle of 1.7}*°* relative to the horizontal. For an aircraft on a 3*°* glideslope, the vortex axis would therefore be 1.3*°* upward from the wing. This difference in the angle of the wake axis is expected to have a very small effect on reflectivity.

To place the radar laterally with respect to the wake, we imagine several parallel approaches that are monitored by the same radar. The largest separation between parallel runways is about 5000 ft (Doyle and McGee 1998). At 6 nm from touchdown, the lateral width of the ILS (Instrument Landing System) approach is 3182 ft for a standard 5*°* splay, and we imagine an aircraft that has strayed to the outer edge of this zone. If the radar is placed in the middle of the two furthest runways we obtain a lateral distance of 0.67 nm. In the presence of a crosswind, we imagine that the wake would be monitored for as long as it remained between the outer edges of the left and right ILS zones. In conclusion, we select (*x*_{rad}*,y*_{rad}*,z*_{rad}) = (6 nm,−0.67 nm,−582 m) relative to the aircraft.

The elevation angle of the radar beam from this location varies between 10.6*°* and 17.1*°* as the scanned range of *z* on the wake center plane varies between *z* = −350 and −200 m (see Figure 4). Since the beamwidth of the radar likely to be used is ≤ 1*°*, ground reflection will be small. To significantly reduce ground and structure clutter, the radar can be placed directly under the flight path. This would require a separate radar for each parallel runway. Another issue is loss of radar sensitivity at smaller ranges; for MIRA-35 this happens for *r* < 360 m (Matthias Bauer-Pfundstein, private communication). However, this loss is probably offset by the increase in power from the *r*^{ factor in }Equation (20).

### a. General procedure

We use aircraft centered coordinates (*x,y,z*), where *x* is streamwise, *y* is spanwise, and *z* is vertical. Air and water are denoted by subscripts ‘a’ and ‘w’, respectively.

The calculation has three steps. The first creates a sample of droplet radii *a*(*t*) from a size distribution pertinent to aerial spray nozzles. The second step generates a spray trail on the starboard side of the aircraft wake. This involves injecting droplets from the sample into the wake, tracking the position **X**(*t*) of each droplet and its radius *a*(*t*) as it evaporates. The flow-field model consists of two counter-rotating vortices whose height decreases with downstream distance *x* behind the aircraft. The third and final step mirrors the starboard trail to the port side and, for a given set of radar parameters, computes the reflectivity of both trails together.

Since it is prohibitive to track all of the actual droplets, the spray trail computed in step two consists of a certain number, *N*_{comp}, of computational droplets. In the reflectivity calculation, each computational droplet is taken to represent a multiplicity, *M*_{true}, of actual droplets, which is simply the ratio of the desired injected volume to the volume injected in the computation.

The equations for droplet motion and evaporation are integrated using the routine LSODE which is described in Radhakrishnan and Hindmarsh (1993) and available from NETLIB. The routine chooses a time step based on a specified error tolerance. Since the flow-field model we have adopted is steady in a reference frame moving with the aircraft, it is not necessary to inject droplets at every time step, which would continuously increase the number of computational droplets. Instead, droplets are injected only during a certain time interval Δ*t*_{0} at the beginning of the calculation. This set of droplets is then evolved for successive Δ*t*_{0} intervals and appended to a file at the end of every Δ*t*_{0} interval. At the end of the computation, the file contains a spray trail that is 7 nm miles long, whose radar reflectivity we subsequently analyze. In actual practice, an aircraft would likely not generate a trail of several nautical miles behind it. Rather, it would release the spray for a short period at pre-selected locations during its approach.

A detailed description of each part of the procedure is described in the following subsections. The evaporation calculation is described in Appendix A2. To avoid stiffness of the system of evolution equations, droplets are removed from the calculation when their radius becomes < 20 *μ*m; at this point their reflectivity is too small to significantly affect the received power.

### b. Droplet trajectory

The position **X**(*t*) and velocity **U**(*t*) of a droplet of mass *m*_{d} evolves according to

where *g*_{eff} = (1 − ρ_{a}*/ρ*_{w})*g* is the effective gravity accounting for buoyancy. The drag force **F**_{D} is given by

where

**u**

_{rel}=

**u**(

**X**) −

**U**,

is the velocity of the air flow, **u**(**X**), relative to the droplet. Evaluation of the drag coefficient, *C*_{D}, is described in Appendix A1.

### c. Droplet size distribution of aerial spray nozzles

When the water jet issues from the nozzle, it will encounter a blast of free-stream air with a speed of 77.2 m s^{−}^{. A comparable situation in the literature is that of a cylindrical liquid jet surrounded by an annulus of co-flowing air (}Lorenzetto and Lefebvre 1977;Varga et al. 2003) which shows that much smaller droplets are produced than for the case of still air (with the same velocity of the liquid jet). This is due to the occurrence of both the Kelvin-Helmholtz instability and Rayleigh-Taylor instabilities. The former is driven by the shear between the air and liquid flow and leads to a smaller instability wavelength. The latter arises due to the acceleration of liquid droplets by the drag force of the air stream. If the droplets are too small, they quickly evaporate. To produce larger droplets, both instabilities can be mitigated by reducing the relative velocity between the liquid jet and free-stream air. This is accomplished by increasing the driving pressure. However, if the droplets are too large, they fall out of the wake. Hence, there is an optimum droplet size.

We were fortunate that an experimental study (Fritz and Hoffmann 2015), which uses a wind-tunnel to mimic the free stream air flow, has recently been performed to characterize the droplet size distribution produced by nozzles used for aerial agricultural spraying. This study did indeed show that larger pressures produce larger droplets. Dr. B. Fritz kindly provided us with an Excel program, developed from that study, which gives parameters of the droplet size distribution for various nozzles, free-stream speeds, and driving pressures. Use of these parameters is now described.

Let *p*(*a*) be the probability density such that *p*(*a*)*da* is the probability that the droplet radius is in the interval [*a, a* + *da*]. The log-normal distribution

with parameters *a*_{0} and *σ*, is commonly used in the spray literature. Fritz’s Excel program provides information about the function *Q*(*a*), defined to be the fractional volume occupied by droplets of radius ≤ *a*. One can show from appropriately integrating (5) that

*Q*(

*a*) = ½(1 + erf

*ξ*),

where

The Excel program provides *a*_{0}_{.}_{5} and *a*_{0}_{.}_{9} defined such that *Q*(*a*_{0}_{.}_{5}) = 0.5 and *Q*(*a*_{0}_{.}_{9}) = 0.9. Using them, (6) can be numerically inverted to yield the parameters *a*_{0} and *σ* of the log-normal distribution (5). The Excel program also provides *a*_{0}_{.}_{1}, which we did not use because we wished to nail the size distribution for large droplets which contribute most to reflectivity.

### d. Choice of nozzle

What is the best drop size distribution? A set of *N*_{d} droplets over which the incident beam is assumed to have uniform intensity, has a reflectivity proportional to (e.g., Doviak and Zrnić 1984, p. 58)

assuming Rayleigh scattering. Maximizing this subject to fixed volume of water and fixed *N*_{d} gives the result that all droplets must be of the same size. Given this result, maximizing *ζ* with respect to the number *N*_{d} subject to fixed volume gives *N*_{d} = 1, i.e., all the volume must be in one droplet. However, such a droplet would likely fall too rapidly. To minimize droplet loss by sedimentation, droplets must not have a terminal velocity larger than the vortex descent speed, *W*_{descent} = 1.75 m s^{−}^{ in the present case. Consulting the terminal velocity plot in }Pruppacher and Klett (1997, p. 416) we conclude that *a* must be ≤ 200 *μ*m. Hence the best distribution is uniform with a drop radius *a* = 200 *μ*m. Since the rate of droplet evaporation is ∝ 1*/a*, i.e., small droplets evaporate faster than larger ones, the above conclusion is not altered by including evaporation.

The above considerations suggest that the following “rate of *ζ*” could be used to initially evaluate different nozzles without having to perform an evaporation and reflectivity calculation:

where the < superscript denotes the exclusion from the sum of droplets that fall away and *N*_{d}(Δ*t*) denotes the number of droplets produced by the nozzle in a period Δ*t*. This suggestion will be tested in §3h.

Table 2 lists the parameters for four candidate nozzles and operating conditions that were selected from B. Fritz’s Excel program based on having peak probability density at a large radius (which rarely exceeded 100 *μ*m) and a high flow-rate. The flow-rate for nozzle 1 was provided by Calvin Kroes (private communication) of CP Products. For nozzles 2 and 3, flow-rates were extrapolated from the values of 5.30 gpm and 2.45 gpm, respectively, at 60 psi reported on the manufacturer’s data sheets (www.translandllc.com/wp-content/uploads/2015/08/Aerial-Flow-Chart-20152.pdf and www.cpproductsinc.com/images/stories/downloads/Misc-Tables/A1-Web%20Aerial%20Tip%20Rate%20Chart.pdf). The extrapolations assumed a square-root dependence of flow-rate on pressure (Lefebvre 1989, p. 157) expected from Bernoulli’s principle, and which the manufacturer’s data follows well. For nozzle 4, the flow-rate versus pressure provided on the manufacturer’s website has a linear rather than square-root dependence. We are grateful to Dr. Brad Fritz (USDA) for measuring the actual flow-rate for us. It turned out to be much lower than the value provided by on the website.

### TABLE 2

Nozzle 1 | Nozzle 2 | Nozzle 3 | Nozzle 4 | |
---|---|---|---|---|

Model | CP-09 | CP-09 | CP11TT | Davidon-Triset |

Pressure (psi) | 90 | 90 | 90 | 90 |

Airspeed (mph) | 175 | 175 | 175 | 175 |

Deflection-plane/body angle | 0° | 0° | 0° | 0° |

Fan angle | Str.^{1} | Str. | Str. | Str. |

Orifice Code | 20 | |||

Orifice diameter (in) | 0.125 | 0.172 | 0.105 | 0.125 |

a_{0 5} (μm)^{Note 2} | 178.5 | 146.5 | 183 | 239.5 |

a_{0 9} (μm)^{Note 2} | 315 | 303 | 359 | 463 |

a_{0} of log-normal | 99.02 | 55.83 | 79.85 | 108.3 |

σ of log-normal | 0.443 | 0.567 | 0.526 | 0.514 |

Flow-rate (gpm)^{Note 3} | 3.70 | 6.49 | 3.00 | 3.06 |

U_{exit} (m s^{−}^{)} | 29.5 | 27.3 | 33.9 | 24.4 |

No. of nozzles per side | 1 | 1 | 1 | 1 |

Gallons per nm (two sides) | 2.96 | 5.19 | 2.4 | 2.45 |

### e. Droplet injection in the computation

Let the origin of coordinates be in the symmetry plane of the aircraft (which corresponds to *y* = 0) with same axial and vertical location as the droplet injector. Droplets are injected in a grid pattern within a square of width *w*_{square} in the *yz*-plane. The pattern consists of *n*_{square} × *n*_{square} droplets. The square is centered at (*x,y,z*) = (0*, fb/*2,0), where *f* represents the fractional spanwise distance from the aircraft center plane to the wing-tip, and was chosen to be *f* = 0.5. A square shape was chosen because the nozzles we have selected are not of the flat fan type. The streamwise extent of the computed mist trail was chosen to be *ℓ*_{trail} = 7 nm (168 seconds of elapsed time from injection) since we wish to detect the trail at 6 nm, the longest distance for which one would want to detect the wake of a heavy aircraft under current separation rules; see Table 1.

The square pattern is injected *n _{x}* times and the time interval between injections is Δ

*t*

_{inject}, whose value is is chosen so that the

*x*spacing between droplets is the same as in the cross-sectional (

*yz*) plane. After a time period Δ

*t*Δ

_{0}≡ n_{x}*t*

_{inject}, an

*n*

^{2}_{square}×

*n*slab of particles has been injected, which is then advanced for successive Δ

_{x}*t*

_{0}periods to form the entire trail. Inertial particles with a small Stokes number tend to an attractor (Haller and Sapsis 2008;Sapsis and Haller 2010) independent of injection location, and therefore where droplets end-up should be insensitive to where they are injected. A brief check on insensitivity to initial conditions will be presented in §3b.

The initial velocity of droplets is set equal to the air velocity, which is justified as follows. From the equations of droplet motion, (2)–(4), the characteristic relaxation time for a droplet to start following a new air speed, imposed at *t* = 0, say, is

where the subscript ‘0’ means that the quantity is calculated at *t* = 0. From this, the characteristic relaxation distance, *ℓ*_{relax} = *u*_{rel}(0)*τ*_{relax} can be evaluated. Note that the initial air speed relative to the water jet is given by *u*_{rel}(0) = *U*_{app} −*U*_{exit}, where *U*_{app} is the approach velocity of the aircraft given in Table 3, and *U*_{exit} is the exit velocity of the water jet given in Table 2. The experiment of Fritz and Hoffmann (2015) measured the size distribution 1.8 m downstream of the nozzle for all straight-stream nozzles. This value is from a private communication from B. Fritz and represents a correction from a value of 1.5 m reported in Fritz and Hoffmann (2015). Figure 1b plots *ℓ*_{relax} as a function of drop radius for nozzle 1. Inspecting it together with the size distribution in Figure 1a, we conclude that most of the droplets are following the air stream at the measurement station of the experiment. If this had not been the case and there had been a relative velocity large enough to give Weber numbers ≳ 10, then it would have been necessary to model further droplet break-up using a secondary break-up model (e.g., Apte et al. 2003).

### TABLE 3

Parameter | Value |
---|---|

Weight, W | 500,000 lb |

Wing span, b | 60 m |

Vortex spacing, b_{0} | 47.9 m |

Vortex circulation, Γ | 526 m^{ s}^{−}^{1} |

Approach speed, U_{app} | 150 knots |

### f. Flowfield of two counter-rotating wing-tip vortices

We consider an aircraft flying straight and level at altitude *z*_{0} = 0 in aircraft coordinates. The velocity field of the airflow in the wake is denoted by lower case **u**(**x**). This velocity field consists of the free stream, ${U}_{\text{app}}\widehat{\mathbf{x}}$(where *U*_{app} is the approach speed of the aircraft), superposed with the flow induced by a pair of counter-rotating vortices with circulations *±*Γ. The centerline of each vortex is at spanwise location *y*_{vort} = *±b*_{0}. Due do their mutually induced velocity, the height *z*_{vort}(*x*) of the vortex pair decreases with distance *x* behind the wing as follows:

*z*

_{vort}(

*x*) =

*z*

_{0}−

*W*

_{desct},

where *W*_{desc} = Γ*/*2*πb*_{0} is the descent speed of the vortex pair. The quantity *t* = *x/U*_{app} is time since the vortex at *x* was shed from the wing.

Each vortex induces a circumferential velocity *u _{θ}*(

*r*) in the cross-plane (

*yz*). For

*u*(

_{θ}*r*), a profile fit to flight data by P. Spalart (Private communication) of Boeing is used:

where *η ≡ r/b*_{0}.

For an elliptically loaded wing, lifting line theory (Batchelor 1967) gives the vortex spacing as

where *b* is the wingspan, and the vortex circulation as

where *W* is the aircraft weight. We use the parameters for a typical heavy aircraft given in Table 3.

### g. Calculation of the received signal and power

A pulse-doppler radar transmits a train of square wave pulses that modulate a carrier wave of frequency *f* = 2*π/ω*. The duration of each pulse is *τ* and the pulse repetition frequency is *f*_{PRF}. After each pulse is transmitted, the transmit-receive switch is set to the receive position and the incoming signal is sampled. Each sample at time *t* is said to come from the *range gate r* = *c*(*t − t*_{t})*/*2, where *t*_{t} is the transmit time of the pulse and the factor of two accounts for the round-trip. Throughout, we consider the case where the transmitting and receiving antenna are the same, the so-called mono-static case.

#### 1) Resolution shell

An important concept is that of the *resolution volume*, ℛ_{τ}(*t*) *R _{τ}*(

*t*), at time

*t*, associated with a single pulse of finite duration,

*τ*(e.g., Yuter 2003, p. 1836). It is defined as the volume from which a signal is received

*at a fixed time t*due to scattering by the pulse. Let

*t*

_{t}mark the beginning of the pulse at the transmitter. A signal from a scatterer at distance

*r*will be received in the time interval

*t*−

*t*

_{t}= 2

*r*/

*c*+

*ξ*

*τ*, 0 ≤

*ξ*≤ 1,

where *ξ* = 0 corresponds to leading edge of the pulse and *ξ* = 1 to its trailing edge. Solving (15) for *r* gives

*r*=

*c*(

*t*−

*t*

_{t}−

*ξ*

*τ*)/2, 0 ≤

*ξ*≤ 1.

Equation (16) defines a spherical shell (called the resolution volume) from which a signal is received at the fixed time *t*. The next subsection describes how the complex voltage received at a given time is evaluated by summing the complex voltages from each droplet in the resolution volume.

#### 2) Received power and signal-to-noise ratio

The material in this subsection is adapted from the texts Doviak and Zrnić (1984) and Ishimaru (1978). Let the “voltage” of the transmitted pulse be (the real part of)

where *V*_{0} is a complex amplitude, *τ* is the pulse width, and “voltage” is defined such that the instantaneous transmitted power is *P*_{t}(*t*) = *V*_{t}(*t*)*V*_{t}* ^{∗}*(

*t*). The voltage at the input terminals of the receiver is (the real part of) the following summation over droplets:

where *A _{m}* is a complex scattering amplitude and

*ω*

_{m}≡

*ω*(1 − 2

*u*

_{m}/

*c*)

is the twice doppler-shifted frequency (*u _{m}* being the radial velocity of the

*m*th droplet) and

*r*is the distance to each droplet. The summation in (18) is taken over the

_{m}*N*

_{d}droplets in the resolution volume associated with the pulse.

The complex scattering amplitude due to each droplet is

The first factor in (20), namely [.]^{, is copied from the square root of the radar equation (e.g., }Doviak and Zrnić 1984, p. 34) for power, where *λ* is the wavelength, and *σ*_{b}* _{m}* is the back-scattering cross-section of each droplet. The function

*G*(

*θ*) is the gain function of the antenna at the droplet, which we have assumed to depend on its angle

_{m}*θ*from the beam centerline. Two loss factors (< 1) have been included in (20):

_{m}*ℓ*

_{w}is the two-way waveguide loss and

*ℓ*

_{B}is the loss due to finite bandwidth of the receiver; values assumed for the present study are given in Table 4. The last factor in (20), exp(i

*ϕ*), accounts for the phase shift induced by back-scattering at the

_{m}*m*th droplet. The back-scattering cross-section and phase-shift (

*σ*

_{b}

*and*

_{m}*φ*, respectively) will be discussed further in §2k.

_{m}### TABLE 4

Manufacturer | EEC | EEC | EEC | Metek | ProSensing |
---|---|---|---|---|---|

Series | DWSR | DWSR | DWSR | ||

Model | 8501S | nnn1C^{Note 1} | 2001X | MIRA-35 | W-SACR |

Frequency, f (GHz) | 3 | 5.9 | 9.6 | 35.1 | 93.9 |

Peak power, P_{t} (kW) | 850 | 250–1000 | 200 | 30 | 1.7 |

Reflector diameter (m) | 4.2 | 4.2 | 2.4 | 1.2,2.0 | 0.9 |

½-power beam width, θ_{b} | 1.83° | 0.95° | 0.95° | 0.52°,0.31° | 0.30° |

Antenna gain, G (dB) | 39.5 | 45 | 45 | 50.4,53.5 | 54.5 |

Pulse width, τ (μs) | 0.4–2 | 0.2–3 | 0.2–2 | 0.1,0.2,0.4 | 0.05–2 |

Range resolution, cτ/2 (m) | 60–300 | 30–450 | 30–300 | 15,30,60 | 7.5–300 |

PRF (kHz) | 0.2–2.4 | 0.2–2.4 | 0.2–2.4 | 2.5,5,10 | ≤ 20 |

U_{max} = c PRF/4 f (m s^{−}^{)} | 5–60 | 2.5–31 | 1.6–19 | 5.3,11,21 | ≤ 16 |

Receiver noise figure (dB) | 2.0 | 2.0 | 2.0 | 6.2 | 6.0 |

2-way waveguide loss^{ (dB)} | 0.8 | 0.8 | 0.8 | 0.8 | 0.8 |

Finite bandwidth loss^{ (dB)} | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 |

Rain attenuation^{ (dB/km)} | 0.005 | 0.03 | 0.12 | 2 | 7 |

Minimum range (m) | 150^{Note 5} |

Radar receivers have electronics that can obtain the real and imaginary parts (denoted *I* and *Q*) of *V*(*t*) and compute the instantaneous received power

*P*

_{r}(

*t*) =

*I*(

*t*)

^{ + }

*Q*(

*t*)

^{ = }

*V*(

*t*)

*V*

^{(}

*t*).

We note in passing that *I*(*t*) and *Q*(*t*) are the components of the real part of *V*(*t*) that are in-phase and 90*°* out-of-phase with the transmitted carrier, respectively. Substituting (18) into (21) and splitting the sum into two parts, following Doviak and Zrnić (1984, §4,1), gives

where *k _{m} ≡ ω_{m}/c*. Arguments for using only the

*first term*in (23) in order to evaluate reflectivity, i.e., for summing the powers received from each droplet, are given by Rayleigh (1945, p. 37), Beckmann (1962), and Doviak and Zrnić (1984, §4,1). The important point, which was phrased eloquently by Rayleigh, is that it is not correct to say that the power in a single return from a random distribution of droplets is the sum of the powers scattered by each. Rather, the result is true only when a

*large ensemble of returns*from a statistically stationary target are averaged. This is most easily seen when we consider the case when all the

*A*are equal (to unity, say). Then, the first term

_{m}*T*

_{1}in (23) is

*T*

_{1}=

*N*

_{d}. If droplet distances are randomly distributed in the resolution shell (assumed to be several wavelengths wide), then the magnitude of the second term will be the average of the summands times the number of terms, i.e., ${T}_{2}\approx {N}_{t}^{-1/2}\times {N}_{t}={N}_{t}^{1/2}\approx {N}_{d}$, where

*N*

_{t}=

*N*

_{d}(

*N*

_{d}−1) is the number of terms in the double sum. Hence both terms in (23) are of

*similar magnitude*. To make the second term smaller than the first, an ensemble average must be taken over many pulses. If the ensemble has

*N*

_{s}phase-uncorrelated samples, then the second term will be ${N}_{\mathrm{S}}^{1/2}$

smaller than the first.

In the present application, the wake descends through the beam and so the target is not stationary, strictly speaking. In §3c we will explicitly verify that, for our case, an average of over a certain number of pulses does indeed yield the first term in (23). Note that a radar set does not have direct access to the first term in (23); only we as simulators do.

To provide a measure of detectability, we will present the signal-to-noise ratio

*P*

_{r1}/

*P*

_{noise},

where *P*_{r1} is the first term in (23) and the average noise power is

*P*

_{noise}=

*k*

_{B}

*T*

_{0}

*F*

_{N}/

*τ*,

where *k*_{B} = 1.381 × 10^{−}^{ J K}^{−}^{ is Boltzmann’s constant, }*T*_{0} = 290 K is a reference temperature set by convention, and *F*_{N} is the overall noise figure of the chain of components in the receiving cascade. Values for *F*_{N} and *τ* are listed in Table 4 for each radar set.

Finally, since each computational droplet represents a multiplicity *M*_{actual} of actual droplets, we have

where *N*_{comp} is the number of computational droplets in the resolution volume. Note that all *M*_{actual} copies of each computational droplet are assumed to be at the same location and therefore their scattered voltages at the receiver add constructively. This assumption does not bias SNR1 since its calculation involves summing individual scattered powers anyway. We claim that the statistics of individual pulse returns are also not affected by this assumption; this will be verified (§3d) in a computational test where the number of computational droplets is increased.

#### 1) Resolution shell

An important concept is that of the *resolution volume*, ℛ_{τ}(*t*) *R _{τ}*(

*t*), at time

*t*, associated with a single pulse of finite duration,

*τ*(e.g., Yuter 2003, p. 1836). It is defined as the volume from which a signal is received

*at a fixed time t*due to scattering by the pulse. Let

*t*

_{t}mark the beginning of the pulse at the transmitter. A signal from a scatterer at distance

*r*will be received in the time interval

*t*−

*t*

_{t}= 2

*r*/

*c*+

*ξ*

*τ*, 0 ≤

*ξ*≤ 1,

where *ξ* = 0 corresponds to leading edge of the pulse and *ξ* = 1 to its trailing edge. Solving (15) for *r* gives

*r*=

*c*(

*t*−

*t*

_{t}−

*ξ*

*τ*)/2, 0 ≤

*ξ*≤ 1.

Equation (16) defines a spherical shell (called the resolution volume) from which a signal is received at the fixed time *t*. The next subsection describes how the complex voltage received at a given time is evaluated by summing the complex voltages from each droplet in the resolution volume.

#### 2) Received power and signal-to-noise ratio

The material in this subsection is adapted from the texts Doviak and Zrnić (1984) and Ishimaru (1978). Let the “voltage” of the transmitted pulse be (the real part of)

where *V*_{0} is a complex amplitude, *τ* is the pulse width, and “voltage” is defined such that the instantaneous transmitted power is *P*_{t}(*t*) = *V*_{t}(*t*)*V*_{t}* ^{∗}*(

*t*). The voltage at the input terminals of the receiver is (the real part of) the following summation over droplets:

where *A _{m}* is a complex scattering amplitude and

*ω*

_{m}≡

*ω*(1 − 2

*u*

_{m}/

*c*)

is the twice doppler-shifted frequency (*u _{m}* being the radial velocity of the

*m*th droplet) and

*r*is the distance to each droplet. The summation in (18) is taken over the

_{m}*N*

_{d}droplets in the resolution volume associated with the pulse.

The complex scattering amplitude due to each droplet is

The first factor in (20), namely [.]^{, is copied from the square root of the radar equation (e.g., }Doviak and Zrnić 1984, p. 34) for power, where *λ* is the wavelength, and *σ*_{b}* _{m}* is the back-scattering cross-section of each droplet. The function

*G*(

*θ*) is the gain function of the antenna at the droplet, which we have assumed to depend on its angle

_{m}*θ*from the beam centerline. Two loss factors (< 1) have been included in (20):

_{m}*ℓ*

_{w}is the two-way waveguide loss and

*ℓ*

_{B}is the loss due to finite bandwidth of the receiver; values assumed for the present study are given in Table 4. The last factor in (20), exp(i

*ϕ*), accounts for the phase shift induced by back-scattering at the

_{m}*m*th droplet. The back-scattering cross-section and phase-shift (

*σ*

_{b}

*and*

_{m}*φ*, respectively) will be discussed further in §2k.

_{m}### TABLE 4

Manufacturer | EEC | EEC | EEC | Metek | ProSensing |
---|---|---|---|---|---|

Series | DWSR | DWSR | DWSR | ||

Model | 8501S | nnn1C^{Note 1} | 2001X | MIRA-35 | W-SACR |

Frequency, f (GHz) | 3 | 5.9 | 9.6 | 35.1 | 93.9 |

Peak power, P_{t} (kW) | 850 | 250–1000 | 200 | 30 | 1.7 |

Reflector diameter (m) | 4.2 | 4.2 | 2.4 | 1.2,2.0 | 0.9 |

½-power beam width, θ_{b} | 1.83° | 0.95° | 0.95° | 0.52°,0.31° | 0.30° |

Antenna gain, G (dB) | 39.5 | 45 | 45 | 50.4,53.5 | 54.5 |

Pulse width, τ (μs) | 0.4–2 | 0.2–3 | 0.2–2 | 0.1,0.2,0.4 | 0.05–2 |

Range resolution, cτ/2 (m) | 60–300 | 30–450 | 30–300 | 15,30,60 | 7.5–300 |

PRF (kHz) | 0.2–2.4 | 0.2–2.4 | 0.2–2.4 | 2.5,5,10 | ≤ 20 |

U_{max} = c PRF/4 f (m s^{−}^{)} | 5–60 | 2.5–31 | 1.6–19 | 5.3,11,21 | ≤ 16 |

Receiver noise figure (dB) | 2.0 | 2.0 | 2.0 | 6.2 | 6.0 |

2-way waveguide loss^{ (dB)} | 0.8 | 0.8 | 0.8 | 0.8 | 0.8 |

Finite bandwidth loss^{ (dB)} | 1.8 | 1.8 | 1.8 | 1.8 | 1.8 |

Rain attenuation^{ (dB/km)} | 0.005 | 0.03 | 0.12 | 2 | 7 |

Minimum range (m) | 150^{Note 5} |

Radar receivers have electronics that can obtain the real and imaginary parts (denoted *I* and *Q*) of *V*(*t*) and compute the instantaneous received power

*P*

_{r}(

*t*) =

*I*(

*t*)

^{ + }

*Q*(

*t*)

^{ = }

*V*(

*t*)

*V*

^{(}

*t*).

We note in passing that *I*(*t*) and *Q*(*t*) are the components of the real part of *V*(*t*) that are in-phase and 90*°* out-of-phase with the transmitted carrier, respectively. Substituting (18) into (21) and splitting the sum into two parts, following Doviak and Zrnić (1984, §4,1), gives

where *k _{m} ≡ ω_{m}/c*. Arguments for using only the

*first term*in (23) in order to evaluate reflectivity, i.e., for summing the powers received from each droplet, are given by Rayleigh (1945, p. 37), Beckmann (1962), and Doviak and Zrnić (1984, §4,1). The important point, which was phrased eloquently by Rayleigh, is that it is not correct to say that the power in a single return from a random distribution of droplets is the sum of the powers scattered by each. Rather, the result is true only when a

*large ensemble of returns*from a statistically stationary target are averaged. This is most easily seen when we consider the case when all the

*A*are equal (to unity, say). Then, the first term

_{m}*T*

_{1}in (23) is

*T*

_{1}=

*N*

_{d}. If droplet distances are randomly distributed in the resolution shell (assumed to be several wavelengths wide), then the magnitude of the second term will be the average of the summands times the number of terms, i.e., ${T}_{2}\approx {N}_{t}^{-1/2}\times {N}_{t}={N}_{t}^{1/2}\approx {N}_{d}$, where

*N*

_{t}=

*N*

_{d}(

*N*

_{d}−1) is the number of terms in the double sum. Hence both terms in (23) are of

*similar magnitude*. To make the second term smaller than the first, an ensemble average must be taken over many pulses. If the ensemble has

*N*

_{s}phase-uncorrelated samples, then the second term will be ${N}_{\mathrm{S}}^{1/2}$

smaller than the first.

In the present application, the wake descends through the beam and so the target is not stationary, strictly speaking. In §3c we will explicitly verify that, for our case, an average of over a certain number of pulses does indeed yield the first term in (23). Note that a radar set does not have direct access to the first term in (23); only we as simulators do.

To provide a measure of detectability, we will present the signal-to-noise ratio

*P*

_{r1}/

*P*

_{noise},

where *P*_{r1} is the first term in (23) and the average noise power is

*P*

_{noise}=

*k*

_{B}

*T*

_{0}

*F*

_{N}/

*τ*,

where *k*_{B} = 1.381 × 10^{−}^{ J K}^{−}^{ is Boltzmann’s constant, }*T*_{0} = 290 K is a reference temperature set by convention, and *F*_{N} is the overall noise figure of the chain of components in the receiving cascade. Values for *F*_{N} and *τ* are listed in Table 4 for each radar set.

Finally, since each computational droplet represents a multiplicity *M*_{actual} of actual droplets, we have

where *N*_{comp} is the number of computational droplets in the resolution volume. Note that all *M*_{actual} copies of each computational droplet are assumed to be at the same location and therefore their scattered voltages at the receiver add constructively. This assumption does not bias SNR1 since its calculation involves summing individual scattered powers anyway. We claim that the statistics of individual pulse returns are also not affected by this assumption; this will be verified (§3d) in a computational test where the number of computational droplets is increased.

### h. Calculation of the doppler spectrum

Radars calculate a doppler spectrum for a given spatial observation location by performing a Fast Fourier Transform (FFT) of complex voltage returns (at the same range gate) from a sequence of pulses separated by Δ*t*_{pulses} = 1/PRF, where PRF is the pulse repetition frequency. We shall do the same for the simulations. From a series of pulse returns, *V _{n}*,

*n*= 0,1

*,…N*

_{FFT}− 1, the normalized transform

and then the power spectrum $S\left(k\right)\equiv \widehat{V}\left(k\right)\widehat{V}\ast \left(k\right)$is computed. Note that the frequency index *k* corresponds to an actual frequency

The frequency *k*_{actual} is in units of (the period of the sequence)^{−}^{ = (}*N*_{FFT}Δ*t*_{pulses})^{−}^{, so in units of s}^{−}^{}

*f*(

*k*

_{actual}) = PRF

*k*

_{actual}/

*N*

_{FFT}.

Equating this to −2(*u*_{doppler}*/c*) *f* gives the doppler velocity associated with each *k*_{actual}. Finally, we state that we use the Hamming window (e.g., Harris 1978).

### i. Antenna gain function

We assume a Gaussian beam with transmitted power flux (power per unit area) vector

where *A* is a coefficient, *θ* is the angle from the beam centerline, and

in terms of the half-power full-width, *θ*_{b}. The total power crossing a sphere of radius *r* is

for a narrow beam. Using the definition of the gain function, *G*(*θ*), we obtain

### j. Radars included in the study

Table 4 lists parameters of currently operational doppler weather/cloud radars considered in the present study. The S, C, and X-band radars chosen are the DWSR series manufactured by EEC (Enterprise Electronics Corporation, Enterprise, Alabama). ARC (Advanced Radar Corporation, Boulder, Co.) makes quite similar C and X-band radars, while Baron Services (Huntsville, Al.) makes similar S, C, and X-band radars. The power and beamwidth values of the C-Band TDWR (Terminal Doppler Weather Radar) deployed at many US airports is subsumed by the range of values provided by the EEC C-band radar, and is therefore not included here.

A number of descriptions of Ka-band (35 GHz) cloud radars have appeared in the literature (Hamazu et al. 2003; Görsdorf et al. 2015). In the present work, we use parameters of the MIRA-35 radar manufactured by Metek (Elmshorn, Germany) which is described in Görsdorf et al. (2015). This choice was motivated by its relatively high power (30 kW). Other Ka-band weather radars, operational at the time of writing are: (i) Scanning 2 kW radars operated by U.S. Department of Energy’s Atmospheric Radiation Measurement (ARM) Climate Research Facility (Widener et al. 2012). (ii) The Copernicus 1 kW radar at Chilbolton Observatory (UK). (iii) An airborne 25 kW multi-frequency (X, Ka, and W-band) radar developed by Prosensing that is being used by NASA’s Langley Research Center for research into the detection and avoidance of super-cooled water droplets.

The science and technology of W-band (94 GHz) radars for cloud and precipitation research is reviewed in Kollias et al. (2007). For the present work we chose the W-SACR radar, which has been developed by the U.S. Department of Energy’s Atmospheric Radiation Measurement (ARM) program (Widener et al. 2012;Kollias et al. 2014).

With increasing frequency, *f*, reflectivity increases as *f*^{ in Rayleigh’s formula (ignoring Mie-scattering corrections). Furthermore, the size of the antenna required to obtain the same beamwidth is reduced. The main drawback of high frequency is increased attenuation due to precipitation between the radar and target. For example, the last entry in }Table 4 gives the attenuation rate in medium rain at W-band as 7 dB/km. A compensating factor is that when there is precipitation, the ambient humidity is also very high and so there is minimal evaporation, and, if natural precipitation is present in the wake, it will also contribute to reflectivity.

### k. Mie cross-section and phase-shift

Since we have can rather large droplets in the present application and frequencies up to 94 GHz, the back-scattering cross-section and phase-shift are obtained using Mie’s formula instead of Rayleigh’s approximation. We used subroutine BHMIE, available from Prof. B.T. Draine’s website at Princeton University, and checked the results using subroutine MIEV0 developed by Dr. W.J. Wiscombe (NASA Goddard).

Some understanding of notation is required to properly use these routines. Let the incident field be of unit magnitude and polarized in the 2-direction (defined to be perpendicular to the plane containing the incident and observer directions). For a spherical target, the scattered field in the far-field is also polarized in the 2-direction and is given by

where *f*_{22} is complex, *k* = 2*π/λ*, *ϑ* is the angle of the observer relative to the direction of propagation of the incident wave, and *r* is distance from the center of the sphere. The backscattering cross-section and phase shift are obtained as

*ϕ*= arg(

*f*

_{22}(

*π*)).

At the start of the reflectivity calculation at a given frequency, we tabulate the ratio *σ*_{b}*/σ*_{b}_{,}_{Rayleigh} and the difference *φ −φ*_{Rayleigh} as a function of droplet radius *a*. Rayleigh’s formulas are (Ishimaru 1978, p. 19)

*σ*

_{b}, Rayleigh = 4|

*K*

_{ε}|

^{(}

*k*

*a*)

^{(}

*π*

*a*

^{),}

*ϕ*

_{Rayleigh}= arg(

*K*

_{ε}),

where *K _{ε}* (a complex number) is given by

The quantity *ε*(*f, T*) is the complex dielectric constant of water; our convention of e^{ωt} for the time dependence requires the imaginary part of *ε* be positive for an absorbing material. It is a function of frequency and temperature and was evaluated using the single Debye model of Liebe et al. (1991) as implemented in subroutines available from Prof. Chris O’Dell’s website at Colorado State University. Since the droplet temperature is almost the same for all drops (the spread was 4 C at most), the dielectric constant *ε* is evaluated at the average temperature of all the drops in the trail.

Figure 2 displays the Mie back-scattering cross-section *σ*_{b} (normalized by the Rayleigh value) and phase-shift *φ* as a function of droplet radius at the five frequencies considered in this work. At the largest radii in this work, *a ≈* 600 *μ*m, the error in using Rayleigh’s cross-section is about 20% at 35.1 GHz.

### l. Choice of ambient temperature and humidity

Droplet evaporation calculations require specification of the ambient temperature and humidity. For guidance on appropriate choices, METARs (Meteorological Aerodrome Reports) during 2000–2014 were downloaded from

and processed for the five busiest airports in the U.S. Figures 3a and b show the monthly-averaged temperature and relative humidity (RH), respectively, when IFR conditions prevailed. Figure 3c shows the percentage of reports that fall into the IFR category. One sees that the average RH is always above 90%. To further synthesize this data, yearly averages were taken (Table 5). Among the five airports, LAX has the highest rate of evaporation in IFR conditions on average since it is the warmest and driest on average. Our choice is the IFR average for LAX, namely, *T* = 15.2 C and RH = 92.7%. Looking at the monthly data for the other four airports, this appears to be a reasonable choice for them also: it is an approximate lower bound for their monthly RH and their temperature is higher only during the summer months when IFR reports are low.

### TABLE 5

Airport | T (C) | RH | % IFR Reports |
---|---|---|---|

ATL | 12.5 | 95.5% | 5.3% |

LAX | 15.2 | 92.7% | 3.5% |

DFW | 9.7 | 94.5% | 1.4% |

ORD | 5.3 | 94.2% | 3.7% |

JFK | 12.1 | 93.9% | 4.8% |

Since it is expensive for flight tests to wait for IFR conditions to occur, and it is desirable to have a wake sensor that can work in a wider variety of atmosphere conditions, we will also consider a case of lower humidity and higher temperature, namely, RH = 60% and *T* = 20 C.

### m. Radar placement with respect to the wake

Here, we choose the radar location (*x*_{rad}*,y*_{rad}*,z*_{rad}) in aircraft-centered coordinates. Based on current wake separations (Table 1), it should not be necessary to examine a wake more than 6 nm behind the aircraft. We therefore chose to present results for the reflectivity at *x* = 6 nm, the worst case for droplet loss by evaporation and sedimentation. The radar is also placed at *x*_{rad} = 6 nm so it can view the *x* = 6 nm wake cross-section at normal incidence. Next, we assume that the aircraft is at the touchdown point. At 6 nm from the touchdown point, the altitude of an aircraft flying a 3*°* glide slope is *H* = 582 m. Therefore, the vertical coordinate of the radar is *z*_{rad} = −582 m. For purposes of this study, we assume an aircraft flying straight and level at this altitude. Initial flight tests would also presumably have the aircraft fly straight and level. In this case, since the vortex pair descends at a speed of *W*_{desc} = 1.75 m s^{−}^{, its axis makes a downward angle of 1.7}*°* relative to the horizontal. For an aircraft on a 3*°* glideslope, the vortex axis would therefore be 1.3*°* upward from the wing. This difference in the angle of the wake axis is expected to have a very small effect on reflectivity.

To place the radar laterally with respect to the wake, we imagine several parallel approaches that are monitored by the same radar. The largest separation between parallel runways is about 5000 ft (Doyle and McGee 1998). At 6 nm from touchdown, the lateral width of the ILS (Instrument Landing System) approach is 3182 ft for a standard 5*°* splay, and we imagine an aircraft that has strayed to the outer edge of this zone. If the radar is placed in the middle of the two furthest runways we obtain a lateral distance of 0.67 nm. In the presence of a crosswind, we imagine that the wake would be monitored for as long as it remained between the outer edges of the left and right ILS zones. In conclusion, we select (*x*_{rad}*,y*_{rad}*,z*_{rad}) = (6 nm,−0.67 nm,−582 m) relative to the aircraft.

The elevation angle of the radar beam from this location varies between 10.6*°* and 17.1*°* as the scanned range of *z* on the wake center plane varies between *z* = −350 and −200 m (see Figure 4). Since the beamwidth of the radar likely to be used is ≤ 1*°*, ground reflection will be small. To significantly reduce ground and structure clutter, the radar can be placed directly under the flight path. This would require a separate radar for each parallel runway. Another issue is loss of radar sensitivity at smaller ranges; for MIRA-35 this happens for *r* < 360 m (Matthias Bauer-Pfundstein, private communication). However, this loss is probably offset by the increase in power from the *r*^{ factor in }Equation (20).

## 3. Results

### a. Signal-to-noise ratio for the IFR case

We begin by considering IFR ambient conditions (RH = 92.7%*,T* = 15.2 C) chosen as described in §2l. Nozzle 1 from Table 2 is used and parameters for the injected square of droplets, described in §2e, are *n*_{square} = 15, *n _{x}* = 120, and

*w*

_{square}= 1 m. Figure 4 shows simulated values of SNR1, calculated using (24), for the five radars listed is Table 4. Each SNR1 plot is an instantaneous range-elevation scan of the

*x*= 6 nm cross-section of the wake and each location on the plot corresponds to the mid-radius of a resolution shell along the beam centerline. The pulse width is chosen to be

*τ*= 0.2

*μ*s for all the radars except for DWSR-8501S, in which case the lowest available

*τ*of 0.4

*μ*s is used. Droplets in a 30 m thick axial slab centered at

*x*= 6 nm are shown in panel (a). Due to centrifugation, larger droplets lie at greater distances from the vortex center, which is devoid of droplets. The very large particles sediment due to gravity after being centrifuged. Except for the S-band radar, all radars give SNR1 > 10 dB at most points surrounding the vortices; the reflectivity is higher for the higher frequency radars. The W-SACR radar gives the highest reflectivity despite having the smallest pulse power.

For all radars, there is a drop in reflectivity near the 2 o’clock and 4 o’clock positions for the left vortex (8 and 10 o’clock positions for the right vortex). This manifests as a crescent-wrench shaped reflectivity pattern that is most prominent for the DWSR 2001X radar. How this feature is related to the droplet configuration, which in turn is related to the vortex flowfield, remains to be elucidated.

Figure 5 shows that with its lowest available pulse width of 0.05 *μ*s, the W-SACR radar is able to resolve some of the spiral structure of the droplet pattern at the expense of some loss in SNR1.

### b. Insensitivity to initial condition

To test sensitivity to initial conditions, instead of injecting droplets in a regular grid pattern on each square, droplets were randomly placed in the squares. Figures 6a and b show that both the droplet configuration at *x* = 6 nm as well SNR1 for DWSR-2001X are changed very little; compare with Figures 4a and d. We expect this to be true for all the radars as well. In another test (Figure 6c and d), the width of the square, *w*_{square}, was reduced from 1 m to 50 cm, keeping the number of droplets fixed. This increases the initial number density in the cross-plane (*yz*) and, to keep the spacing the same in the streamwise (*x*) direction, the injection interval Δ*t*_{inject} was also halved. One might think that this would increase the number density downstream. However, the flow tends to both reduce number density (where there is rotation) and increase it (where there is strain), and eventually, the number density tends to a distribution that is mostly independent of initial condition.

### c. Pulse to pulse fluctuation and averaging

The SNR1 results in §3a were based on using the first term in Equation (23), which sums the powers reflected by individual droplets. It was argued that for a statistically stationary target, this should equal the average power from many pulses. In the present case, the droplet configuration is not spatially homogeneous and is descending through a fixed beam. Hence, the question arises whether the powers returned from a sequence of pulses can be considered to be statistically stationary in a certain interval, and if so, how many pulses is sufficient to recover the SNR1 values presented.

To obtain complex voltage returns from a sequence of pulses one needs to evolve the wake in time, however, the method that was described in §2 gives a trail of droplets at a single instant of time, *t*. To evolve this configuration to time *t* + Δ*t*, the configuration at *t* is translated horizontally by Δ*x* = −*U*_{app}Δ*t*, i.e., the droplet trail is assumed to be invariant in a reference frame moving to the left with the airplane. This procedure does not correspond exactly to reality, but captures both the rotation of droplets around the vortices, and their vertical descent with time at a fixed location. The received complex voltage is evaluated using (18) at a sequence of times separated by the pulse repetition period, keeping the resolution volume centered at (*x,y,z*) = (6 nm,−50 m,−230 m). This location corresponds to the upper SNR1 peak of the crescent wrench in Figure 4. The value of the pulse repetition frequency (PRF) was chosen to be at or close to the highest value available for each radar.

Figure 7 shows the result for a time period during which a cluster of droplets enters and leaves the beam. The power in individual pulse returns is shown in gray. The total number of active pulses changes from radar to radar because their beam widths and PRFs are different. In particular, the period of activity was found to equal the time it would take the vortex pair to descend through roughly one-third of the vertical projection of the half-power beam width. The average of pulse powers is shown in green over an averaging segment whose length is 512 pulses. The red curve shows the value of SNR1. Our assumption was that SNR1 should equal the green level. This is seen to be true to a good degree. The fluctuations are due to statistical error and were found to decrease with increasing the averaging interval. It is worth remembering here that the radar has access to only the individual pulse returns and their average (for example the green values); only the simulation has access to the red curve (SNR1).

### d. Convergence of pulse statistics

Recall that each computational droplet at a single location represents *M*_{true} actual droplets located at different positions; in fact *M*_{true}*≈* 100 in the calculations presented. It was claimed (§2g) that this should not affect pulse statistics, provided the number of computational droplets is sufficiently large. To verify this, the number of computational droplets was increased by four. Random placement of droplets was employed in the injected squares. Four realizations of the droplet trail were generated using different random number seeds for the initial size distribution and droplet placement. The four realizations were then merged into one trail for the radar reflectivity calculation. Figure 8 shows that the probability density *p*(|*V*|) of the modulus |*V*| of complex voltage is unchanged by the resolution refinement. The probability densities are very well fit by the Rayleigh distribution (Beckmann 1962)

having the same mean as the data. The Rayleigh distribution results when the scattering amplitude is the same for all droplets and the phases uniformly distributed. The case selected is the same as that shown in Figure 7d (apart from the random placement of droplets in the injected squares). Pulses in the interval of stationarity were chosen, namely pulse number ∈ [−8000,8000].

### e. A non-IFR condition

It would be valuable to have the capability to detect wakes in non-IFR conditions. Furthermore, in a flight test study of the feasibility of the present proposal, it would be too costly to wait until IFR conditions occur before a test can be conducted. For this reason it is of interest to know what reflectivity is obtained at less humid and less cold conditions. We chose RH = 60% and *T* = 20 C.

With the previous choice of *n*_{square} = 15 and *n _{x}* = 120 as injection parameters, it was found that a high rate of evaporation resulted in a small number of computational droplets remaining near

*x*= 6 nm. This increased statistical error. To reduce sampling error, an ensemble of ten trails were computed with different random number seeds for the droplet size sample. The ensemble was then combined into a single trail for the reflectivity analysis. As a result, the total number of computational droplets is so large that each one presents only 9.8 true droplets in the reflectivity analysis.

Figure 9 shows that only the high-frequency radars, MIRA-35 and W-SACR, give positive values of SNR1 (dB) in the vicinity of the vortices, and even these values are marginal. To increase SNR1, the number of nozzles could be increased; for instance four nozzles on each side of the aircraft would increase SNR1 by 6 dB.

There is a powerful method that enables detection even when SNR1 (dB) < 0. It comes at the cost of increased dwell and processing time. We learnt about the method from notes on the sensitivity of the MIRA-35 radar given to us by Matthias Bauer-Pfundstein (Metek). It is also briefly described in Görsdorf et al. (2015, p 680). The idea is that in a discrete Fourier transform, the noise is spread equally to all the frequency bins, whereas the spectrum of the signal is confined to only a few of the bins. (The latter is true provided the probability distribution of droplet velocities in the resolution cell is narrow compared to 2*U*_{max} for the radar. For MIRA- 35, for example, at PRF = 10 kHz we have 2*U*_{max} = 42 m s^{−}^{ and so this is unlikely to be an issue.) Hence, an FFT effectively reduces the noise by a factor of }*N*_{FFT}.

To investigate this technique, complex white noise with a mean power equal to *P*_{noise} for the radar was added to complex voltages of pulse returns. Illustrative results are shown in Figure 10. Panels (a) and (b) are for a range cell centered at the left white dot in Figure 9e where SNR1 = −2.7 dB. Panels (c) and (d) are for the right white dot where SNR1 is even lower, namely, −7.2 dB. Consider panels (a) and (b). An averaging of pulse power returns by the radar would give values (the green line) only slightly above the noise, not enough for a positive detection. Averaging the doppler spectra from 10 segments gives panel (b) with a peak 40 dB above the noise. In the present example, this would require a dwell time of 0.5 s for each elevation angle. For the second location where SNR1 is weaker, the doppler spectrum has a peak that is about 25 dB above the noise (using the same dwell time).

### f. Power-weighted average radial velocity

It has been stated (Doviak and Zrnić 1984, §5.2) that the first moment of the doppler spectrum is the radial velocity of droplets in the resolution volume, weighted by their individual scattered powers. This is reasonable although we have neither encountered nor attempted a proof of it starting from Equation (2). Figure 11 displays the power-weighted average radial velocity corresponding to the cases previously shown in Figure 4. Only points where SNR (dB) > 0 are shown. The actual radial velocity (with respect to the radar) of the gas is shown in panel (f). The radar data appears as a filtered version of the actual velocity and, due to particle centrifugation, is unable to detect the maximum value of 19.8 m s^{−}^{ in the vortex core. Nevertheless, the radars give a good representation of the gas velocity where particles are present. To estimate vortex circulation, the observed velocity is multiplied by 2}*πr* where *r* is the distance from the vortex center and can be determined from the location of zero radial velocity.

### g. Effect of vortex core growth

The vortex model presented in §2f assumed that the vortex core radius does not grow with downstream distance. In reality, the core radius grows and the peak velocity diminishes. Experiments on aircraft wakes (see Govindaraju and Saffman 1971) have shown that the core radius, *r*_{1}, defined to be the radius where the tangential velocity peaks, diffuses as follows with distance *x* behind the plane:

where *b*_{1} is a constant that depends on the vortex Reynolds number, Γ*/ν*. In the present case Γ*/ν* = 3.5 × 10^{, or which Table I in }Govindaraju and Saffman (1971) gives *b*_{1} = 1.3 × 10^{−}^{ (using the corrected value listed in their table). The value of the parameter }*x*_{0} (called the virtual origin) can be obtained from the core radius at *x* = 0. Plotting the profile we find that *r*_{1} = 0.8 m for the present case which gives *x*_{0} = −558.8 m.

To include core growth in the vortex model of §2f, the variable *η* is redefined as

where *g* is the growth function

Figure 12 displays the decay of peak tangential velocity from *x* = 0 to *x* = 6 nm.

Figures 13a and b show the droplet configuration and SNR1, respectively, for the MIRA-35 radar when core growth is included; they should be compared with Figures 4a and e. This comparison shows that slightly less centrifugation takes places with core diffusion which makes the region of reflectivity slightly smaller.

### h. Comparison of four nozzles

Finally, Table 6 compares SNR1 reflectivity for the four nozzles listed in Table 2. The range cell for all cases is centered at (*x,y,z*) = (6nm,−50m,−230m), which corresponds to the top of the crescent wrench in Figure 4. IFR ambient conditions have been assumed. It is observed that the quantity ${\dot{\zeta}}_{\text{nozzle}}^{<}$(defined in Equation 9), which depends only on the droplet size distribution produced by a nozzle, provides an excellent indicator of the relative performance of different nozzles.

### TABLE 6

Nozzle no. | gpm | SNR1 (dB) | ˙ζ< (dB m^{ s}^{−}^{) nozzle} | ζ (dB m_{x}^{) near vortex at }x = 6 nm |
---|---|---|---|---|

1 | 3.70 | 19.15 | −157.4 | −178.8 |

2 | 6.49 | +1.03 | +1.0 | +0.6 |

3 | 3.00 | −1.04 | −1.2 | −1.4 |

4 | 3.06 | −0.5 | −0.3 | −1.0 |

To provide further insight we compute the quantity *ζ _{x}* which is defined to be

*ζ*(see Equation 8) per unit axial length of the wake. It is calculated as a diagnostic of the droplet trajectory and size evolution and is shown in Figure 14. The solid lines give the total value (over an entire cross-section) and diagnose total evaporative loss. The dashed lines give the contribution from droplets in a neighborhood (defined in the caption) of the vortex: these curves diagnose both evaporative loss and loss by sedimentation. The ordering of

*ζ*values (pertaining to the neighborhood of the vortex) at

_{x}*x*= 6 nm, which are also listed in the last column of Table 6, matches the ordering of SNR1 for the different nozzles.

All statements of comparison in the following are relative to nozzle 1 and make reference to Figure 14a. If the increased SNR1 reflectivity of nozzle 2 were due to increased volume alone, we would get a 2.4 dB increase in reflectivity. The actual increase is 1.03 dB. To understand this, we first observe that nozzle 2 (solid red line) initially has a 0.7 dB higher value of *ζ _{x}*, less than the 2.4 dB increase in its volume flow-rate. This is because nozzle 2 produces more small droplets. By

*x*= 6 nm the 0.7 dB increase has been reduced to 0.55 dB because the smaller droplets of nozzle 2 evaporate faster. The fact that

*ζ*in the vicinity of the vortex is 0.6 dB higher must arise from the fact that the smaller droplets of nozzle 2 have sedimented less.

_{x}Consider nozzle 3 (green curves). If its smaller reflectivity (relative to nozzle 1) were due to decrease in volume, then we would expect a −0.91 drop in SNR1 which is close to what is obtained. This is understandable given that its initial *ζ _{x}* is very nearly the same as for nozzle 1. This is surprising given that nozzle 3 has many more smaller droplets. However, close inspection of its size distribution (green line in Figure 14b) shows that it also has more droplets that are very large (specifically

*a*> 370

*μ*m). This fact also explains the more rapid loss of

*ζ*by sedimentation (dashed green curve in Figure 14a) and less rapid loss by evaporation (solid green curve). Overall, these two effects balance and the final effect that remains is that due to volume decrease.

_{x}Despite its smaller flow-rate, nozzle 4 has a higher initial value of *ζ _{x}*; see the solid blue curve. This is because it produces more large droplets. Unfortunately, they rapidly fall out of the wake (dashed blue curve).

### a. Signal-to-noise ratio for the IFR case

We begin by considering IFR ambient conditions (RH = 92.7%*,T* = 15.2 C) chosen as described in §2l. Nozzle 1 from Table 2 is used and parameters for the injected square of droplets, described in §2e, are *n*_{square} = 15, *n _{x}* = 120, and

*w*

_{square}= 1 m. Figure 4 shows simulated values of SNR1, calculated using (24), for the five radars listed is Table 4. Each SNR1 plot is an instantaneous range-elevation scan of the

*x*= 6 nm cross-section of the wake and each location on the plot corresponds to the mid-radius of a resolution shell along the beam centerline. The pulse width is chosen to be

*τ*= 0.2

*μ*s for all the radars except for DWSR-8501S, in which case the lowest available

*τ*of 0.4

*μ*s is used. Droplets in a 30 m thick axial slab centered at

*x*= 6 nm are shown in panel (a). Due to centrifugation, larger droplets lie at greater distances from the vortex center, which is devoid of droplets. The very large particles sediment due to gravity after being centrifuged. Except for the S-band radar, all radars give SNR1 > 10 dB at most points surrounding the vortices; the reflectivity is higher for the higher frequency radars. The W-SACR radar gives the highest reflectivity despite having the smallest pulse power.

For all radars, there is a drop in reflectivity near the 2 o’clock and 4 o’clock positions for the left vortex (8 and 10 o’clock positions for the right vortex). This manifests as a crescent-wrench shaped reflectivity pattern that is most prominent for the DWSR 2001X radar. How this feature is related to the droplet configuration, which in turn is related to the vortex flowfield, remains to be elucidated.

Figure 5 shows that with its lowest available pulse width of 0.05 *μ*s, the W-SACR radar is able to resolve some of the spiral structure of the droplet pattern at the expense of some loss in SNR1.

### b. Insensitivity to initial condition

To test sensitivity to initial conditions, instead of injecting droplets in a regular grid pattern on each square, droplets were randomly placed in the squares. Figures 6a and b show that both the droplet configuration at *x* = 6 nm as well SNR1 for DWSR-2001X are changed very little; compare with Figures 4a and d. We expect this to be true for all the radars as well. In another test (Figure 6c and d), the width of the square, *w*_{square}, was reduced from 1 m to 50 cm, keeping the number of droplets fixed. This increases the initial number density in the cross-plane (*yz*) and, to keep the spacing the same in the streamwise (*x*) direction, the injection interval Δ*t*_{inject} was also halved. One might think that this would increase the number density downstream. However, the flow tends to both reduce number density (where there is rotation) and increase it (where there is strain), and eventually, the number density tends to a distribution that is mostly independent of initial condition.

### c. Pulse to pulse fluctuation and averaging

The SNR1 results in §3a were based on using the first term in Equation (23), which sums the powers reflected by individual droplets. It was argued that for a statistically stationary target, this should equal the average power from many pulses. In the present case, the droplet configuration is not spatially homogeneous and is descending through a fixed beam. Hence, the question arises whether the powers returned from a sequence of pulses can be considered to be statistically stationary in a certain interval, and if so, how many pulses is sufficient to recover the SNR1 values presented.

To obtain complex voltage returns from a sequence of pulses one needs to evolve the wake in time, however, the method that was described in §2 gives a trail of droplets at a single instant of time, *t*. To evolve this configuration to time *t* + Δ*t*, the configuration at *t* is translated horizontally by Δ*x* = −*U*_{app}Δ*t*, i.e., the droplet trail is assumed to be invariant in a reference frame moving to the left with the airplane. This procedure does not correspond exactly to reality, but captures both the rotation of droplets around the vortices, and their vertical descent with time at a fixed location. The received complex voltage is evaluated using (18) at a sequence of times separated by the pulse repetition period, keeping the resolution volume centered at (*x,y,z*) = (6 nm,−50 m,−230 m). This location corresponds to the upper SNR1 peak of the crescent wrench in Figure 4. The value of the pulse repetition frequency (PRF) was chosen to be at or close to the highest value available for each radar.

Figure 7 shows the result for a time period during which a cluster of droplets enters and leaves the beam. The power in individual pulse returns is shown in gray. The total number of active pulses changes from radar to radar because their beam widths and PRFs are different. In particular, the period of activity was found to equal the time it would take the vortex pair to descend through roughly one-third of the vertical projection of the half-power beam width. The average of pulse powers is shown in green over an averaging segment whose length is 512 pulses. The red curve shows the value of SNR1. Our assumption was that SNR1 should equal the green level. This is seen to be true to a good degree. The fluctuations are due to statistical error and were found to decrease with increasing the averaging interval. It is worth remembering here that the radar has access to only the individual pulse returns and their average (for example the green values); only the simulation has access to the red curve (SNR1).

### d. Convergence of pulse statistics

Recall that each computational droplet at a single location represents *M*_{true} actual droplets located at different positions; in fact *M*_{true}*≈* 100 in the calculations presented. It was claimed (§2g) that this should not affect pulse statistics, provided the number of computational droplets is sufficiently large. To verify this, the number of computational droplets was increased by four. Random placement of droplets was employed in the injected squares. Four realizations of the droplet trail were generated using different random number seeds for the initial size distribution and droplet placement. The four realizations were then merged into one trail for the radar reflectivity calculation. Figure 8 shows that the probability density *p*(|*V*|) of the modulus |*V*| of complex voltage is unchanged by the resolution refinement. The probability densities are very well fit by the Rayleigh distribution (Beckmann 1962)

having the same mean as the data. The Rayleigh distribution results when the scattering amplitude is the same for all droplets and the phases uniformly distributed. The case selected is the same as that shown in Figure 7d (apart from the random placement of droplets in the injected squares). Pulses in the interval of stationarity were chosen, namely pulse number ∈ [−8000,8000].

### e. A non-IFR condition

It would be valuable to have the capability to detect wakes in non-IFR conditions. Furthermore, in a flight test study of the feasibility of the present proposal, it would be too costly to wait until IFR conditions occur before a test can be conducted. For this reason it is of interest to know what reflectivity is obtained at less humid and less cold conditions. We chose RH = 60% and *T* = 20 C.

With the previous choice of *n*_{square} = 15 and *n _{x}* = 120 as injection parameters, it was found that a high rate of evaporation resulted in a small number of computational droplets remaining near

*x*= 6 nm. This increased statistical error. To reduce sampling error, an ensemble of ten trails were computed with different random number seeds for the droplet size sample. The ensemble was then combined into a single trail for the reflectivity analysis. As a result, the total number of computational droplets is so large that each one presents only 9.8 true droplets in the reflectivity analysis.

Figure 9 shows that only the high-frequency radars, MIRA-35 and W-SACR, give positive values of SNR1 (dB) in the vicinity of the vortices, and even these values are marginal. To increase SNR1, the number of nozzles could be increased; for instance four nozzles on each side of the aircraft would increase SNR1 by 6 dB.

There is a powerful method that enables detection even when SNR1 (dB) < 0. It comes at the cost of increased dwell and processing time. We learnt about the method from notes on the sensitivity of the MIRA-35 radar given to us by Matthias Bauer-Pfundstein (Metek). It is also briefly described in Görsdorf et al. (2015, p 680). The idea is that in a discrete Fourier transform, the noise is spread equally to all the frequency bins, whereas the spectrum of the signal is confined to only a few of the bins. (The latter is true provided the probability distribution of droplet velocities in the resolution cell is narrow compared to 2*U*_{max} for the radar. For MIRA- 35, for example, at PRF = 10 kHz we have 2*U*_{max} = 42 m s^{−}^{ and so this is unlikely to be an issue.) Hence, an FFT effectively reduces the noise by a factor of }*N*_{FFT}.

To investigate this technique, complex white noise with a mean power equal to *P*_{noise} for the radar was added to complex voltages of pulse returns. Illustrative results are shown in Figure 10. Panels (a) and (b) are for a range cell centered at the left white dot in Figure 9e where SNR1 = −2.7 dB. Panels (c) and (d) are for the right white dot where SNR1 is even lower, namely, −7.2 dB. Consider panels (a) and (b). An averaging of pulse power returns by the radar would give values (the green line) only slightly above the noise, not enough for a positive detection. Averaging the doppler spectra from 10 segments gives panel (b) with a peak 40 dB above the noise. In the present example, this would require a dwell time of 0.5 s for each elevation angle. For the second location where SNR1 is weaker, the doppler spectrum has a peak that is about 25 dB above the noise (using the same dwell time).

### f. Power-weighted average radial velocity

It has been stated (Doviak and Zrnić 1984, §5.2) that the first moment of the doppler spectrum is the radial velocity of droplets in the resolution volume, weighted by their individual scattered powers. This is reasonable although we have neither encountered nor attempted a proof of it starting from Equation (2). Figure 11 displays the power-weighted average radial velocity corresponding to the cases previously shown in Figure 4. Only points where SNR (dB) > 0 are shown. The actual radial velocity (with respect to the radar) of the gas is shown in panel (f). The radar data appears as a filtered version of the actual velocity and, due to particle centrifugation, is unable to detect the maximum value of 19.8 m s^{−}^{ in the vortex core. Nevertheless, the radars give a good representation of the gas velocity where particles are present. To estimate vortex circul}