## Abstract

This study investigates the influence of microscale turbulent clustering of cloud droplets on the radar reflectivity factor and proposes a new parameterization to account for it. A three-dimensional direct numerical simulation of particle-laden isotropic turbulence is performed to obtain turbulent clustering data. The clustering data are then used to calculate the power spectra of droplet number density fluctuations, which show a dependence on the Taylor microscale-based Reynolds number (Re_{λ}) and the Stokes number (St). First, the Reynolds number dependency of the turbulent clustering influence is investigated for 127 < Re_{λ} < 531. The spectra for this wide range of Re_{λ} values reveal that Re_{λ} = 204 is sufficiently large to be representative of the whole wavenumber range relevant for radar observations of atmospheric clouds. The authors then investigate the Stokes number dependency for Re_{λ} = 204 and propose an empirical model for the turbulent clustering influence assuming power laws for the number density spectrum. For Stokes numbers less than 2, the proposed model can estimate the influence of turbulence on the spectrum with an RMS error less than 1 dB when calculated over the wavenumber range relevant for radar observations. For larger Stokes number droplets, the model estimate has larger errors, but the influence of turbulence is likely negligible in typical clouds. Applications of the proposed model to two idealized cloud observing scenarios reveal that microscale turbulent clustering can cause a significant error in estimating cloud droplet amounts from radar observations with microwave frequencies less than 13.8 GHz.

## 1. Introduction

Clouds play crucial roles in the heat and water systems of Earth. To improve our understanding of cloud physics, a large number of observational studies have been conducted to estimate the spatial distribution of cloud microphysical properties, such as the cloud water mixing ratio and the effective droplet radius. Radar is one of the most powerful tools since it can provide two- or three-dimensional estimates of cloud microphysical properties over a large domain (Okamoto et al. 2007; Stephens et al. 2008; Ellis and Vivekanandan 2011). In radar observations, microwave radiation is transmitted from an antenna toward a target cloud and the reflected microwaves received and analyzed. The relation between the transmitted power *P*_{t} and the received power *P*_{r} of the microwaves is given by the following radar equation:

where *G* is the antenna gain, *k*_{m} is the microwave wavenumber, *R* is the distance between the antenna and the cloud, *K* is the dielectric coefficient of a water droplet, *V* is the measurement volume, and *Z* is the radar reflectivity factor (mm^{6} m^{−3}). Crucially, *Z* is dependent on the cloud microphysical properties, implying that cloud properties can be estimated from *Z*.

The relation between *Z* and cloud microphysical properties is explained by two mechanisms: incoherent scattering and coherent scattering. Incoherent scattering occurs when the cloud droplets are dispersed randomly and uniformly (Bohren and Huffman 1983). The radar reflectivity factor for the incoherent scattering case is proportional to the sum of the Rayleigh scattering intensity from each droplet and independent of the microwave frequency *f*_{m}. On the other hand, coherent scattering—often referred to as Bragg scattering—occurs when the droplets are distributed nonuniformly. The nonuniform distribution causes the interference of scattered microwaves, which in turn increases the radar reflectivity factor obtained from Eq. (1). This coherent scattering by discrete particles is more specifically referred to as “particulate” Bragg scattering (Kostinski and Jameson 2000). Coherent scattering can also be caused by a nonuniform distribution of the refractive index of clear air—which may be referred to as “clear-air Bragg scattering.” Most studies assume that particulate Bragg scattering is insignificant in atmospheric clouds (Gossard and Strauch 1983). However, this assumption is contradicted by the observations of developing cumulus clouds by Knight and Miller (1993) and Knight and Miller (1998). They observed significant differences between the radar reflectivity factors for 10- and 3-cm microwaves, which are classified in the S and X bands, respectively. A similar wavelength dependency of the radar reflectivity factor was found for the case of smoke plumes from an intense industrial fire by Rogers and Brown (1997), who compared the data observed by a UHF wind profiler (wavelength 32.8 cm) and an X-band radar (3.2 cm). Knight and Miller (1998) explained that these differences resulted from coherent scattering by nonuniform cloud droplet concentrations created by the turbulent mixing of cloud with environmental clear air (i.e., turbulent entrainment). That is, they attributed the differences to the large-scale nonuniform distribution of cloud droplets. Erkelens et al. (2001) investigated the influence of turbulent entrainment on the observations of Knight and Miller (1998). They analyzed the observational data using an equation for clear-air Bragg scattering based on the − power law of scalar concentration spectra in turbulence and concluded that turbulent entrainment is not the only relevant factor for coherent scattering in cumulus clouds. Kostinski and Jameson (2000) pointed out that microscale turbulent droplet clustering is also a cause of coherent scattering in cumulus clouds. The turbulent clustering is caused by an inertial effect of particles within turbulent flows, which generates microscale nonuniform particle distributions, often referred to as preferential concentration (Maxey 1987; Squires and Eaton 1991; Wang and Maxey 1993; Chen et al. 2006). Note that turbulent clustering can occur even without large-scale nonuniform particle distributions. Many authors have investigated the effect of turbulent clustering on collisions of cloud droplets (e.g., Sundaram and Collins 1997; Pinsky and Khain 1997; Reade and Collins 2000; Ayala et al. 2008b,a; Onishi et al. 2009; Woittiez et al. 2009; Wang et al. 2009; Jin et al. 2010; Onishi and Vassilicos 2014). The possible importance of turbulent clustering for the radar reflectivity factor, however, was first suggested by Kostinski and Jameson (2000). Recently, Dombrovsky and Zaichik (2010) estimated the influence of turbulent clustering based on the semianalytical clustering model of Zaichik and Alipchenkov (2007). Their estimate indicated that turbulent clustering considerably increases the radar reflectivity factor. These studies clearly suggest that the influence of turbulence should be carefully considered to obtain reliable estimates of cloud microphysical properties from radar observations. However, until now there has been no reliable way to estimate this influence. One recent approach is that of Dombrovsky and Zaichik (2010), but their estimate relied on a highly simplified clustering model that adopted a simple extrapolation for large scales.

This study, therefore, aims to investigate the influence of microscale turbulent clustering on the radar reflectivity factor and construct a reliable model for estimating it. A three-dimensional direct numerical simulation (DNS) of particle-laden isotropic turbulence is performed in order to obtain turbulent clustering data, and then the influence of turbulence is analyzed and modeled. The model is then applied to two idealized radar observation scenarios to assess the influence quantitatively.

## 2. Computational method

### a. Air turbulence

The governing equations of turbulent airflow are the continuity and Navier–Stokes equations for three-dimensional incompressible flows:

where *u*_{i} is the fluid velocity in the *i*th direction, *ρ*_{a} is the air density, *p* is the pressure, *ν* is the kinematic viscosity, and *F*_{i} is the external forcing term.

The fourth-order central-difference scheme (Morinishi et al. 1998) was used for the advection term and the second-order Runge–Kutta scheme was used for time integration. The velocity and pressure were coupled by the highly simplified marker and cell (HSMAC) method (Hirt and Cook 1972). Statistically steady-state turbulence was formed by applying an external forcing using the reduced-communication forcing (RCF) method of Onishi et al. (2011), which maintains the intensity of large-scale eddies while keeping a high parallel efficiency.

It should be noted that atmospheric turbulence is typically neither homogeneous nor isotropic. However, the assumptions of homogeneity and isotropy are reasonable for the small scales corresponding to the wavenumber range relevant to radar observations (see section 4a). Although energy-containing large-scale eddies generate large-scale inhomogeneity and anisotropy, dissipative small-scale eddies work to flatten the inhomogeneities, leading to local homogeneity and isotropy. This local homogeneity assumption is the basis of most turbulence models.

### b. Droplet motions

Droplet motions are tracked by the Lagrangian method. The governing equation for droplet motion is

where *υ*_{i} is the particle velocity in the *i*th direction, *τ*_{p} is the droplet relaxation time, and *g*_{i} is the gravitational acceleration in the *i*th direction (Onishi et al. 2009; Woittiez et al. 2009; Wang et al. 2009; Jin et al. 2010; Onishi et al. 2013). Equation (4) is based on the following two assumptions: (i) the droplets are Stokes particles [i.e., spherical with small particle Reynolds numbers (Re_{p} ≡ 2*r*_{p}|**u** − **v**|/*ν*)] and (ii) the ratio of the density of droplets to that of the surrounding air is much larger than unity (Maxey and Riley 1983; Kim et al. 1998). The relaxation time for the Stokes particle is given by

In clouds, turbulent modulation and the frequency of droplet collisions likely remain small since the volume fraction *ϕ* is smaller than 10^{−6}. Thus, these effects were neglected for simplicity (Matsuda et al. 2012).

### c. Computational conditions

The computational domain was set to a cube with edges of length 2*πL*_{0}, where *L*_{0} is the representative length scale. Periodic boundary conditions were applied in all three directions. The domain was discretized uniformly into grid points, giving a grid spacing of Δ = 2*πL*_{0}/*N*_{g}. The DNS was performed for four turbulent flows, each with a different value of the Taylor-microscale-based turbulent Reynolds number, defined as Re_{λ} = *l*_{λ}*u*_{rms}/*ν*, where *u*_{rms} is the RMS value of the velocity fluctuations and *l*_{λ} is the Taylor microscale. Table 1 shows the computational parameters for the air turbulence simulations and the statistical results obtained. The kinematic viscosity was set to 1.5 × 10^{−5} m^{2} s^{−1}. Note that the flow conditions are the same as those of Onishi et al. (2011), who showed that turbulent flows are well resolved under all the chosen conditions. Since the resolutions were chosen to satisfy *k*_{max}*l*_{η} ≈ 2, where *k*_{max} is the maximum wavenumber given by *k*_{max} = *N*_{g}(2*L*_{0})^{−1}, in our DNS experiments the nondimensional energy dissipation rate was essentially the same for all the flows.

The droplet radius *r*_{p} was varied so that the Stokes number, defined as St = *τ*_{p}/*τ*_{η}, where *τ*_{η} = (*ν*/*ϵ*)^{1/2} is the Kolmogorov time scale, took values of 0.05, 0.1, 0.2, 0.5, 1.0, 2.0, and 5.0. The droplet radii for St = 1.0 were 22.9, 23.1, 20.2, and 23.4 *μ*m for Re_{λ} = 127, 204, 322, and 531, respectively. The number of droplets was set to 8 × 10^{6}, 1.5 × 10^{7}, 5 × 10^{7}, and 5 × 10^{7} for Re_{λ} = 127, 204, 322, and 531, respectively. For most of the simulations, the gravitational accelerations *g*_{i} were set to zero in order to focus on the Re_{λ} and St dependencies of turbulent clustering. However, we have also performed DNS experiments with (*g*_{1}, *g*_{2}, *g*_{3}) = (0, *g*, 0), where *g* = 9.8 m s^{−2}, to investigate the influence of gravitational droplet settling. Details of the numerical conditions are described in section 4d.

The code is fully parallelized for a three-dimensional domain decomposition using a Message Passing Interface (MPI) library (Onishi et al. 2013). The largest simulation (i.e., the case Re_{λ} = 531) was performed on 32 nodes of the Earth Simulator 2 supercomputer operated by the Japan Agency for Marine-Earth Science and Technology (JAMSTEC).

## 3. Radar reflectivity factor

The intensity of reflected microwaves is determined by the scattering intensity of each droplet and the interference between scattered microwaves. Since the radii of cloud droplets are much smaller than the wavelength of the microwaves, the scattering is classified as Rayleigh scattering, which gives intensities proportional to . In the case where droplets are randomly and uniformly dispersed—implying zero spatial correlations between droplets—the effects of interference cancel and become zero. Thus, the radar reflectivity factor for randomly and uniformly located monodispersed droplets *Z*_{random} is given by

where *n*_{p} is the droplet number density. Note that *Z*_{random} is independent of *k*_{m}. In the alternative case, where droplets form clusters, the effect of interference appears as an additional term and the radar reflectivity factor becomes dependent on *k*_{m}. The radar reflectivity factor for monodispersed clustering droplets *Z*_{cluster} is given by

where *κ* is the absolute value of the difference between the incident and scattered wavenumber vectors **k**_{inc} and **k**_{sca}; that is, *κ* = |**k**_{inc} − **k**_{sca}| (Gossard and Strauch 1983; Erkelens et al. 2001). Because the antenna receives backward scattering, *κ* becomes 2*k*_{m}, providing the Doppler effect is small enough. The power spectrum of droplet number density fluctuations *E*_{np}(*k*) represents the intensity of clustering for wavenumber *k*. It should be noted that Eq. (7) assumes isotropic turbulent clustering. Unfortunately, there is no widely accepted analytical model for *E*_{np}(*k*). Jeffery (2000, 2001a,b) derived theoretical power spectrum models for *kl*_{η} ≥ 0.1 based on a *δ*-correlated closure. His models, however, were obtained by assuming nonzero diffusivity of the particle number density. Recently, power spectra of the number density fluctuation have been obtained using DNS (Rani and Balachandar 2003; Shotorban and Balachandar 2007; Jin et al. 2010). For example, Jin et al. (2010) obtained the power spectra for several values of St. However, they did not discuss the Re_{λ} dependency or propose any model to predict the spectra. For this study, *E*_{np}(*k*) was calculated from the DNS data as

where is the spectral density function of droplet number density, given by

where the angle brackets represent the ensemble average. The variable is the Fourier coefficient of the spatial droplet number density distribution *n*_{p}(**x**), given by the following discrete Fourier transform

where *V*_{c} is the cubic domain with edge of length 2*πL*_{0}, and *n*_{p}(**x**) is given by

where **x**_{p,j} is the position vector of the *j*th droplet inside a target domain, *N*_{p} is the total number of droplets, and *δ*(**x**) is the Dirac delta function. The Fourier coefficients of Eq. (11) are then given by

Note that the transform implies periodicity of the droplet distribution. Finally, substitution of Eq. (12) into Eq. (9) yields

Note that terms for particle pairs with *j* = *j′* are removed from Eq. (13) in order to eliminate white noise from . For efficient computation of , Eq. (13) was transformed to

Equation (14) still requires *N*_{p} × *N*_{k} calculations, where *N*_{k} is the number of discrete wavenumber vectors **k** = (*m*_{1}/*L*_{0}, *m*_{2}/*L*_{0}, *m*_{3}/*L*_{0}), where *m*_{1}, *m*_{2}, and *m*_{3} are arbitrary integers. For this study, we chose 19 representative wavenumbers, giving *kL*_{0} = |**k**|*L*_{0} values of 1, 2, 3, 4, 6, 8, 12, 16, 24, 32, 48, 64, 96, 128, 192, 256, 384, 512, and 768. These 19 values of *k* cover the wavenumber range more or less uniformly on a log scale. We calculated for the discrete wavenumbers located between the spherical surfaces with radii of *k* − Δ*k*/2 and *k* + Δ*k*/2, where Δ*k* was set to 1/*L*_{0}, and summed to calculate *E*_{np}(*k*) from Eq. (8). We have checked that the results are insensitive to an increased number of representative wavenumbers. The ensemble average in Eq. (14) was obtained by averaging temporal slices of the droplet distributions. For Re_{λ} ≤ 322, *E*_{np}(*k*) was obtained by averaging 10 temporal slices, while for Re_{λ} = 531, just a single temporal slice was used as it provided a sufficiently large data volume to obtain reliable statistics. For Re_{λ} ≤ 322, the temporal slices of the droplet distributions were sampled for large intervals of *T*_{0} = *L*_{0}/*U*_{0} to eliminate the temporal correlations between the distributions.

## 4. Results and discussion

### a. Droplet distribution in turbulence

Figure 1 shows the spatial distributions of droplets within the range 0 < *z* < 4*l*_{η}, where *l*_{η} = (*ν*^{3}/*ϵ*)^{1/4} is the Kolmogorov scale, for St = 0.05, 0.2, 1.0, and 5.0 at Re_{λ} = 204. The number of the particles in each figure is similar; about 3.5 × 10^{4}. Void areas due to turbulent clustering are clearly observed for St = 1.0. For St < 1.0, the void areas are less clear. For St > 1.0, small void areas with dimensions less than 40*l*_{η} are less clear than for St = 1.0, but void areas larger than 40*l*_{η} are more prominent.

Figure 2 shows *E*_{np}(*k*) for St = 1.0, for different values of Re_{λ}. The arrow indicates the range of the nondimensional wavenumber relevant to actual radar observations, which we can estimate from the range of *f*_{m} used, and the typical *l*_{η} that apply in atmospheric clouds. The microwave frequencies used for radar observations of clouds or precipitation range from the S band (*f*_{m} ~ 2 GHz) to the W band (*f*_{m} ~ 100 GHz). The typical *l*_{η} in atmospheric clouds ranges from 5 × 10^{−4} to 1 × 10^{−3} m, which we estimate based on the energy dissipation rate *ϵ* ~ 10^{−3}–10^{−2} m^{2} s^{−3} and *ν* ~ 10^{−5} m^{2} s^{−1}. Since *E*_{np}(*k* = 2*k*_{m}) is used for estimating *Z*_{cluster} for *f*_{m} = *k*_{m}*c*_{m}/2*π*, where *c*_{m} is the speed of light, the relevant wavenumber range for radar observations is estimated to be 0.05 < *kl*_{η} < 4.0.

In atmospheric clouds, Re_{λ} ranges from 10^{3} to 10^{4}, higher than the maximum Re_{λ} value (=531) used within our simulations. However, for the wavenumber range 0.05 < *kl*_{η} < 4.0, the maximum difference between *E*_{np}(*k*) values for Re_{λ} = 204 and 531 in Fig. 2 is 11%, while for Re_{λ} = 127 and 531, the maximum difference is 22%. These differences correspond to differences of 0.47 and 1.1 dB in the increment to *Z* given by Eq. (7), respectively, where a value in units of decibels is defined as *A*^{dB} = 10 log_{10}*A* for a given value of *A*. Since errors of around 1 dB are unavoidable in radar observations (Bringi et al. 1990; Carey et al. 2000), the dependency of *E*_{np}(*k*) on Re_{λ} is sufficiently small for Re_{λ} > 200 and thus for the wavenumber range relevant for radar observations. Thus, this study uses *E*_{np}(*k*) at Re_{λ} = 204 to estimate *Z* for radar observations of atmospheric clouds.

Figure 3 shows *E*_{np}(*k*) of droplet number density fluctuations for different values of St at Re_{λ} = 204. The horizontal and vertical axes are normalized using *l*_{η} and the average number density 〈*n*_{p}〉. It is clear that *E*_{np}(*k*) depends strongly on St. For St ≤ 1.0, the peak values of *E*_{np}(*k*) are located around *kl*_{η} = 0.2 [i.e., (*kl*_{η})_{peak} ≈ 0.2] and become higher as St becomes closer to 1. This indicates that the representative void scale is almost constant, but that the number density difference between sparse (void) and dense (cluster) areas increases as St increases. This is because the number density of inertial particles tends to concentrate more in high-strain-rate and low-vorticity regions as *τ*_{p} increases (Maxey 1987). Since the Kolmogorov-scale eddies have the largest effect on the motions of St < 1 droplets, (*kl*_{η})_{peak} is almost fixed at about 0.2. On the other hand, for St ≥ 1.0 the peak location moves toward lower wavenumbers as St increases, indicating that the representative void scale becomes larger as St increases. This is because large-scale eddies preferentially concentrate large St droplets, and small-scale eddies tend to destroy this preferential concentration by uncorrelated stirring. This scale dependent clustering mechanism is explained by Goto and Vassilicos (2006) and Yoshimoto and Goto (2007). These features for St ≤ 1.0 and St ≥ 1.0 are consistent with what is observed in Fig. 1. Jin et al. (2010) examined the St dependency of *E*_{np}(*k*). Their power spectra show generally good agreement with ours, confirming the reliability of our simulation. It should be noted that Jin et al. (2010) used Re_{λ} = 102, which is too small for us to use their spectra to estimate the influence of turbulent clustering on radar observations.

### b. Influence of turbulent clustering on the radar reflectivity factor

In this section the influence of turbulence on the radar reflectivity factor *Z* is estimated from the *E*_{np}(*k*) curves shown in Fig. 3 and compared with the estimate by Dombrovsky and Zaichik (2010). The influence of turbulence is evaluated using the clustering coefficient *ζ* defined by

We estimate *ζ* from *E*_{np}(*k*) using the equation

which is obtained by substituting Eqs. (6) and (7) into Eq. (15). Dombrovsky and Zaichik (2010) semianalytically estimated *ζ* using the following equation (Kostinski and Jameson 2000):

where *g*(*r*) is the radial distribution function (RDF) defined as *g*(*r*) = 〈*n*_{p}(**x**)*n*_{p}(**x** + **r**)〉/〈*n*_{p}〉^{2}, where *r* = |**r**|. [Note that Eq. (17) can be considered as the Fourier transform of Eq. (16) under isotropic conditions.] Dombrovsky and Zaichik (2010) adopted an RDF model based on the probability density function (PDF) approach (Zaichik and Alipchenkov 2007). The RDF model is

where the model parameters *c* and Γ are given by

This RDF model is applicable to the case St < 0.6, Re_{λ} > 30 and *r* < *l*_{η}. For performing the integration from 0 to infinity in Eq. (17), Dombrovsky and Zaichik (2010) extrapolated the RDF model to the separation range *l*_{η} < *r* < *c*^{1/Γ}*l*_{η} and assumed *g*(*r*) − 1 = 0 for *r* > *c*^{1/Γ}*l*_{η}.

Figure 4 shows clustering coefficients *ζ* for Re_{λ} = 204. The horizontal axis is the microwave wavenumber difference normalized by *l*_{η}. The vertical axis is normalized by to eliminate the effect of droplet number density. The horizontal arrow indicates the typical range of the nondimensional wavenumber in actual radar observations; 0.05 < *κl*_{η} < 4.0, corresponding to the arrow in Fig. 2. In this wavenumber range, the *ζ* values obtained from the *E*_{np}(*k*) data show a strong dependency on St and a monotonically decreasing trend against *kl*_{η}. The St dependency was also analyzed in Dombrovsky and Zaichik (2010). However, they did not find the monotonically decreasing trends seen in our results: the *ζ* values of Dombrovsky and Zaichik (2010) are almost constant in the low-wavenumber region and decrease with wavy oscillations as the wavenumber increases. These two characteristics exist because the extrapolation of their RDF model to *r* > *l*_{η} is physically unrealistic, and the RDF for *r* > *l*_{η} has a large influence on *ζ*. Since the number density correlation function 〈*n*_{p}(**x**)*n*_{p}(**x** + **r**)〉 for *r* > *l*_{η} is considered in the calculation of *E*_{np}(*k*), the estimate of *ζ* based on *E*_{np}(*k*) is more reliable than that based on their RDF model.

### c. Modeling of the turbulent clustering influence on radar reflectivity factor

In this section we develop a new empirical model of the radar reflectivity factor as a function of *kl*_{η} and St by fitting curves to those of *E*_{np}(*k*) in Fig. 3. The *E*_{np}(*k*) curves in Fig. 3 show power-law-like slopes for *kl*_{η} < (*kl*_{η})_{peak} and *kl*_{η} > (*kl*_{η})_{peak}, so we assume that the power spectra of number density fluctuations asymptotically approach *E*_{np}(*k*)/(〈*n*_{p}〉^{2}*l*_{η}) ≈ *c*_{1}(*kl*_{η})^{α} and *E*_{np}(*k*)/(〈*n*_{p}〉^{2}*l*_{η}) ≈ *c*_{1}(*kl*_{η})^{β} in the small- and large-wavenumber regions, respectively. For intermediate wavenumbers, we connect the two asymptotic regimes with a function based on the power spectrum model of scalar concentration fluctuations suggested by Hill (1978) (model 2). That is, we assume

where *S*(*ξ*) is the nondimensional spectrum defined as *S*(*ξ*) = *S*(*kl*_{η}) = *E*_{np,model}(*k*)/(〈*n*_{p}〉^{2}*l*_{η}), *ξ* is the nondimensional wavenumber defined as *ξ* = *kl*_{η}, and *γ* is a positive-valued parameter enabling us to adjust the peak value in the transition region. The term *ξ*^{†} is defined as *ξ*^{†} = ln(*ξ*/*ξ*_{t}), where *ξ*_{t} is the nondimensional transition wavenumber. Under the conditions *S*(*ξ*) → *c*_{1}*ξ*^{α} for *ξ* → 0 and *S*(*ξ*) → *c*_{2}*ξ*^{β} for *ξ* → ∞, *ξ*_{t} becomes (*c*_{2}/*c*_{1})^{1/(α−β)} so that *S*(*ξ*) becomes

Figure 5 shows the values of the parameters *c*_{1}, *α*, *c*_{2}, and *β* for each of the Stokes numbers. The parameters *c*_{1} and *α* were obtained by finding a least squares fit within the wavenumber range *kl*_{η} < 0.1, while *c*_{2} and *β* were obtained by finding the best fit for *kl*_{η} > 0.7. The solid lines are the best-fit curves, given by

and

Note that the reference data for *c*_{1} and *α* for St = 5.0 were fewer and less reliable than those for the other St values. This is because the peak of *E*_{np}(*k*) is located in the region *kl*_{η} < 0.1. Thus, the *α* value for St = 5.0 was not considered in obtaining Eq. (23). Since the diffusion coefficient of the droplet number density is much smaller than *ν*, the value of *β* for St ≪ 1 should be −1, which is the power index of the power spectrum of scalar concentration fluctuation in the viscous-convective range (Batchelor 1959; Grant et al. 1968; Goto and Kida 1999). The behavior of *β* for St ≫ 1 is unknown. For this study, we simply assume that *β* approaches −1 for St ≫ 1. This simple assumption affects only large St values, where the influence of turbulence is negligibly small in radar observations (as will be seen in Fig. 9, described in section 4e).

Figure 6 shows the parameter *γ* which appears in Eq. (22). We see that *γ* is nearly constant for St ≤ 1.0 but becomes larger for St > 1.0. For this study we have ignored the *γ* values for St > 1.0 and averaged those for St ≤ 1.0, giving

This is justified by the fact that adjustment of *γ* for St > 1.0 resulted in only small improvements to the fit.

To summarize, our model of the influence of the microscale turbulent clustering *S*(*ξ*) is estimated from Eq. (22) together with Eqs. (23)–(25). The clustering coefficient (i.e., the increment due to turbulent clustering) is then obtained from Eq. (16) as

where *κ* = 2*k*_{m}.

Figure 7 shows the RMS error *e*_{rms} of the proposed model *E*_{np,model}(*k*), evaluated in units of decibels, in which *e*_{rms} is calculated for the wavenumber range relevant for radar observations—0.05 < *kl*_{η} < 4.0—using the equation

where the superscript dB denotes a value in units of decibels, *ξ*′ is defined as *ξ*′ = ln(*kl*_{η}), and and are set to ln(0.05) and ln(4.0), respectively. Also *e*_{rms} has its minimum value at St = 0.2. For St < 0.2, for which turbulent clustering is less pronounced than for St > 0.2, the increase of *e*_{rms} as St reduces is due to an increase of statistical error in the reference *E*_{np}(*k*) data, which is observed as fluctuations of *E*_{np}(*k*) for St = 0.05 and 0.1 and *kl*_{η} > 1 in Fig. 3. In the range St > 0.2, we see significant increases in error from St = 2.0, mainly caused by the simple assumptions for *β* when St ≫ 1. However, as will be discussed later in section 4e, the clustering effect for St ≥ 5.0 is irrelevant for actual radar observations. Except for St = 5.0, *e*_{rms} is well below 1 dB. As mentioned in section 4a, we consider this level of error to be acceptable.

### d. Influence of gravitational settling on the power spectra of number density fluctuations

The influence of gravitational settling on the power spectra of number density fluctuations has been investigated by performing additional DNSs with gravity included. Nondimensional parameters relevant for gravitational effects are *S*_{υ} = *υ*_{T}/*u*_{η}—where *υ*_{T} is the terminal velocity given by *τ*_{p}*g* and *u*_{η} is the Kolmogorov velocity (Wang and Maxey 1993; Grabowski and Vaillancourt 1999)—and the Froude number (Fr = *υ*_{T}/*u*_{rms}); *S*_{υ} measures the settling influence on small scales and Fr measures it on large scales. Strictly speaking, we need multiple parameters covering the wide range of clustering scales. However, we consider these two parameters—covering the two ends of the scale range—to be sufficient for our analysis. Table 2 shows the values of *S*_{υ} and Fr in the additional DNS runs. Re_{λ} was set to 204 and St to unity.

Figure 8 shows the settling influence on *E*_{np}(*k*). As the particle settling becomes stronger, *E*_{np}(*k*) decreases at small scales and increases at large scales. The decrease at small scales corresponds to the increase of *S*_{υ} and indicates that settling weakens small-scale clustering (Ayala et al. 2008a,b; Woittiez et al. 2009). The increase at large scales, on the other hand, corresponds to the increase of Fr and indicates that anisotropies generated by settling lead to large-scale clustering. [Woittiez et al. (2009) observed a nearly-two-dimensional “curtain shape” clustering.] However, the increase at large scales is outside of the wavenumber range relevant for radar observations, so we need only consider *S*_{υ}. The maximum differences between *E*_{np}(*k*) for *S*_{υ} > 0 and *S*_{υ} = 0 are 0.35, 0.68, 1.4, and 2.2 dB for *S*_{υ} = 1.37, 2.71, 6.88, and 11.1, respectively. That is, the errors of the proposed model are smaller than 1 dB for *S*_{υ} ≤ 2.7 ≈ 3. Thus the proposed model is reliable for *S*_{υ} < 3.

### e. Turbulent clustering influence in radar observations estimated by the proposed model

Recent radar observations of clouds and precipitation have been conducted using microwaves in six frequency bands: the S, C, X, Ku, Ka, and W bands, with typically used frequencies of 2.8, 5.3, 9.4, 13.8, 35, and 94 GHz, respectively. The S-, C-, X-, and Ku-band radars are often used for observing precipitation, while Ka- and W-band radars are used only for clouds. This is because the Rayleigh scattering approximation is invalid when *d*_{p}/*λ*_{m} is larger than about , where *d*_{p} is the droplet diameter and *λ*_{m} the microwave wavelength. That is, the Rayleigh scattering approximation is invalid for *r*_{p} larger than 430 and 160 *μ*m for the Ka and W bands, respectively. The S and X bands are also used for cloud observations, often using dual frequencies to retrieve liquid water content (LWC)—the liquid water mass contained in a unit volume of air. In dual-wavelength radar observations, LWC is estimated from the dual-wavelength ratio (DWR), which is defined as the difference of *Z*^{dB} (dB*Z*) [=10 log_{10}*Z* (mm^{6} m^{−3})] for two frequencies (Knight and Miller 1998; Vivekanandan et al. 1999; Wang et al. 2005; Ellis and Vivekanandan 2011) and analyzed by considering the frequency dependency of microwave attenuation. To take a recent example, Ellis and Vivekanandan (2011) proposed and tested a technique for estimating cloud LWC using the National Center for Atmospheric Research (NCAR) simultaneous S-band–Ka-band dual-polarimetric (S-PolKa) radar system.

We now estimate increments to *Z*^{dB} due to microscale turbulent clustering [i.e., ] under two idealized cloud scenarios: (i) a stratocumulus case, where *l*_{η} = 1 × 10^{−3} m and thus *ϵ* ≈ 3 × 10^{−3} m^{2} s^{−3}, and (ii) a cumulus case, where *l*_{η} = 5 × 10^{−4} m and thus *ϵ* ≈ 5 × 10^{−2} m^{2} s^{−3} (Pinsky et al. 2008). For both cases, the volume fraction was kept at *ϕ* = 10^{−6}. For simplicity, the cloud droplets are assumed to be monodispersed, though in real clouds, droplets have various kinds of size distribution. Since *ϕ* is fixed to a constant value, 〈*n*_{p}〉 varies depending on the droplet size as follows:

where and *ρ*_{p}/*ρ*_{a} = 840. For atmospheric clouds, *ϕ* is in the range 10^{−7} < *ϕ* < 10^{−6} (Kokhanovsky 2004), so the chosen condition *ϕ* = 10^{−6} corresponds to dense clouds.

It should be noted that turbulent entrainment can also be a cause of droplet number density fluctuations in real clouds. However, we can separate the influence of turbulent entrainment from that of turbulent clustering by using the well-known scalar concentration spectrum *E*_{θ}(*k*). In the inertial-convective range (*kl*_{η} < 0.1), *E*_{θ}(*k*) is given by

where *χ* is the scalar dissipation rate and *C*_{c} is the Obukhov–Corrsin constant (Sreenivasan 1996; Goto and Kida 1999). In the viscous-convective range (*kl*_{η} > 0.1), *E*_{θ}(*k*) is given by

where *C*_{b} is the Batchelor constant (Batchelor 1959; Grant et al. 1968; Oakey 1982; Goto and Kida 1999). Note that Eq. (30) is valid when the scalar diffusive coefficient *D* is much smaller than *υ*. Since the scales of clustering and entrainment are typically well separated, the correlation between the number density fluctuations due to clustering and entrainment should be negligible. Thus, *E*_{np}(*k*) for both clustering and entrainment should be given by , where and are the power spectra for clustering and entrainment, respectively. Here we focus on the influence of .

Figure 9 shows values of estimated from the proposed model together with the DNS results. For *r*_{p} smaller than 100 *μ*m, the increment for each value of *f*_{m} is larger in the cumulus case than in the stratocumulus case. The increment is larger than the observation error level (1 dB) for *f*_{m} ≤ 9.4 GHz in the stratocumulus case and for *f*_{m} ≤ 13.8 GHz in the cumulus case. However, the influence is not significant for larger drops. This is because, for St > 1, *E*_{np}(*k*)/(〈*n*_{p}〉^{2}*l*_{η}) becomes smaller as St increases, and 〈*n*_{p}〉 decreases as *r*_{p} increases, as indicated in Eq. (28). The peak of for each *f*_{m} is located at St ≈ 0.1, which corresponds to droplet radii of a few tens of micrometers—close to the typical cloud droplet size. Thus, the influence of turbulence can cause a significant error in retrieving the LWC of clouds in radar observations with microwave frequencies less than 13.8 GHz (S, C, X, and Ku bands).

For *S*_{υ} > 3 however, the influence of turbulent clustering would be overestimated significantly owing to the influence of gravitational settling, as discussed in section 4d. The threshold *S*_{υ} = 3 corresponds to *r*_{p} ≈ 20 *μ*m for the stratocumulus case and to *r*_{p} ≈ 30 *μ*m for the cumulus case. The possible overestimates for large *r*_{p} cannot be ignored, but they do not affect our main argument that the influence of turbulence can cause a significant error in radar observations using the S, C, X, and Ku bands.

Figure 9 indicates that the radar reflectivity factor becomes larger as *f*_{m} becomes lower and that the maximum difference between the S and X bands is approximately 8 dB. These characteristics are in good agreement with the observations of developing cumulus clouds by Knight and Miller (1998), in which the reflectivity factor for the S band is about 10 dB larger than for the X band. As mentioned in the introduction, a similar frequency dependency was found for the case of smoke plumes of an industrial fire by Rogers and Brown (1997). Although the constituents and sizes of smoke particles are different from cloud droplets, turbulence could also influence the radar reflectivity factor. Thus, as speculated by Kostinski and Jameson (2000) based on the theory of particulate Bragg scattering, turbulent clustering may influence the radar observation of cumulus clouds and smoke plumes.

## 5. Conclusions

This study has investigated the influence of microscale turbulent clustering of cloud droplets on the radar reflectivity factor and proposed an empirical parameterization to account for it. Three-dimensional direct numerical simulations (DNS) of particle-laden isotropic turbulence were performed in order to obtain turbulent clustering data, from which power spectra of droplet number density fluctuations were calculated. The calculated power spectra show dependencies on the Taylor microscale-based Reynolds number (Re_{λ}) and the Stokes number (St). To begin, we investigated the dependency of the turbulent clustering influence on Re_{λ}. Results for a wide range of Re_{λ} values (up to 531) reveal that Re_{λ} = 204 is large enough to be representative of the whole wavenumber range relevant to radar observations of atmospheric clouds (0.05 < *kl*_{η} < 4, where *k* is the wavenumber and *l*_{η} is the Kolmogorov scale). (Smaller Re_{λ} values were found to be unable to represent the power spectrum for low wavenumbers.) Setting Re_{λ} = 204, we then investigated the dependency on St. We observed that for St < 1 the peak of the power spectrum is located at around *kl*_{η} = 0.2 with the peak value increasing as the Stokes number increases toward unity. For St > 1, the peak location moves to lower wavenumbers as St increases. Based on these observations, and assuming that the power spectrum follows distinct power laws in the small- and large-wavenumber regions, we proposed an empirical model that approximately fits the power spectrum of number density fluctuations *E*_{np}(*k*). From this model, it was then possible to calculate the clustering coefficient *ζ* (i.e., the influence of turbulence on the radar reflectivity factor). A comparison between the model estimates and the DNS results for *E*_{np}(*k*) confirms the reliability of the model for droplets with Stokes number smaller than 2. For larger Stokes number droplets, the model estimate has larger errors, but the influence of turbulence of such large droplets is likely negligible in typical clouds. The proposed model has been applied to two idealized radar-observation scenarios: (i) a stratocumulus case, where *l*_{η} = 1 × 10^{−3} m, and (ii) a cumulus case, where *l*_{η} = 5 × 10^{−4} m. In both cases, the droplet volume fraction was 10^{−6} and the microwave frequency *f*_{m} ranged from 2.8 to 94 GHz. The results show that the influence of microscale turbulent clustering on the radar reflectivity factor is significant for droplets with radius smaller than 100 *μ*m for *f*_{m} ≤ 9.4 GHz in the stratocumulus case and for *f*_{m} ≤ 13.8 GHz in the cumulus case. That is, the influence of turbulent clustering can cause a significant error in retrieving cloud liquid water content from radar observations with microwave frequencies less than 13.8 GHz (S, C, X, and Ku bands). Additional DNSs with gravitational effects included reveal that the influence of gravitational settling causes significant errors in the model estimates when the nondimensional terminal velocity *S*_{υ} is larger than 3. These errors for large particles cannot be ignored but do not alter our main conclusions.

## Acknowledgments

The authors thank Dr. Masaki Katsumata of the Japan Agency for Marine-Earth Science and Technology for his helpful comments on radar observations. The numerical simulations presented here were carried out on the Earth Simulator 2 and ICE X supercomputer systems operated by the Japan Agency for Marine-Earth Science and Technology. This study was supported by Grant-in-Aid for Young Scientists (A) (20686015) and by Grant-in-Aid for JSPS Fellows (21·241).

## REFERENCES

*New J. Phys.,*

**10,**075015, doi:.

*New J. Phys.,*

**10,**075016, doi:.

*Absorption and Scattering of Light by Small Particles.*Wiley, 530 pp.

_{a}band radar measurements

*Radar Observation of Clear Air and Clouds.*Developments in Atmospheric Science, Vol. 14, Elsevier, 292 pp.

*Phys. Fluids,*

**18,**115103, doi:.

*Phys. Fluids,*

**22,**055106, doi:.

*Phys. Rev. Lett.,*

**108,**224502, doi:.

*J. Geophys. Res.,*

**112,**D08216, doi:.

*Phys. Fluids,*

**21,**125108, doi:.

*Phys. Fluids,*

**19,**118107, doi:.

*J. Geophys. Res.,*

**110,**D19201, doi:.

*Phys. Fluids,*

**19,**113308, doi:.