II. RAIN DROP INERTIAL PARAMETERS

In this section we aim at finding rain drop inertial parameters, which quantify how much a particle behaves as a perfect tracer or, alternatively speaking, how much a scatterer responds to fluctuations in wind velocities throughout its trajectory. To obtain rain drop inertial parameters that can be used in radar-based retrieval techniques for the estimation of EDR from the DSW, it makes sense to start from fundamental physics, which are the equations of motion for a rain drop. However the problem is that the equations of motions cannot be applied directly in a retrieval algorithm because only limited information is available. Therefore it is necessary to simplify the equations of motion down to scalars, the inertial parameters that can represent the behavior of an ensemble of particles. The advantage of this approach is that it is possible to verify all model assumptions and better adapt this approach to specific applications. In addition the same procedure can be applied to other classes of hydrometeors such as ice crystals. Works with such a fundamental approach already exist. Khvorostyanov and Curry [Reference Khvorostyanov and Curry25] started from the equations of motion for drop/ice particles and successfully estimated their terminal fall speeds. One of the most difficult parts in such works is the determination of the drag force coefficient for rain drops and ice particles. Here we reuse existing work from Khvorostyanov and Curry [Reference Khvorostyanov and Curry25] but now with the aim of estimating inertial parameters instead of the terminal fall speed. The equations of motion using the Cartesian coordinate system for a particle in the x/y- and z-direction can be written as:

(1) $$\displaystyle{{d{v_{p,xy}}} \over {dt}} = \displaystyle{{ \pm {F_{d,xy}}} \over m} = \mp\! {\eta _{xy}}{({v_{a,xy}} - {v_{p,xy}})^2},$$

(2) $$\displaystyle{{d{v_{p,z}}} \over {dt}} = \displaystyle{{ - {F_g} + {F_b} \pm {F_{d,z}}} \over m} = - {\eta _z}v_t^2 \mp {\eta _z}{({v_{a,z}} - {v_{p,z}})^2},$$

with

(3) $${\eta _{xy}} = \displaystyle{{{\rho _F}{A_{xy}}{C_{d,xy}}} \over {2m}}, \quad {\eta _z} = \displaystyle{{{\rho _F}{A_z}{C_{d,z}}} \over {2m}},$$

where v _p,* is the particle velocity, v _a,* air velocity, v _t the terminal particle fall speed, F _d,* = 1/2ρ _F A _* C _d,*(v _a,* − v _p,*)² the drag force, ρ _F the air density, A _* the projected particle surface area in direction *, C _d,* the drag force coefficient, F _g  = mg gravity, m the particle mass, g the gravitational acceleration and F _b  = ρ _F V _b g buoyancy force and V _b the volume of the particle. Here η _z [m^–1] and η _xy [m^–1] are scalars to write the equations of motion in a brief way. The calculations for the drag force coefficients C _d,* are extensive and can be obtained from the work of Khvorostyanov and Curry [Reference Khvorostyanov and Curry25]. For rain drops in the air the buoyancy force can be neglected, as the gravitational force is much larger [Reference Khvorostyanov and Curry25]. However, for crystallized ice particles, which can have much lower mass densities, the buoyancy force should also be considered. The direction of the drag force is opposite to the particle motion with respect to its surrounding air motion, which determines the right sign in the equations. For the remainder of this paper the terminal fall velocities and scalars η _xy and η _z are obtained from the work of Khvorostyanov and Curry [Reference Khvorostyanov and Curry25]. As big rain drops have a spheroidal shape, also an axis ratio formula as a function of equivolumetric spherical drop diameter is used from Beard and Chuang [Reference Beard and Chuang26]. When there is a balance between the gravitational, drag and buoyancy forces, the drop terminal fall speed is obtained as [Reference Khvorostyanov and Curry25]:

(4) $${v_t}(D) = {\left[ {\displaystyle{{2( \vert mg - {F_b} \vert )} \over {{\rho _F}{A_z}{C_{D,z}}}}} \right]^{1/2}},$$

where D is the equivolumetric spherical drop diameter. The total particle velocity can be written as:

(5) $${\vec v_p} = {\vec v_t} + {\vec v_a} + \vec v_p^{\,\prime}, $$

where ${\vec v_t}$ is the particle terminal fall velocity, ${\vec v_a}$ is the air velocity and $\vec v{^{\,\prime}_p}$ the additional velocity difference due to relaxation. The direction of the terminal fall velocity is always in the negative z-direction, pointing towards the earth. The direction and magnitude of $\vec v_p^{\,\prime} $ can be obtained by numerically integrating the equations of motions, equations (1) and (2).

From analytical solutions to simple cases it is possible to estimate typical numbers that characterize the motion of an ensemble of particles. One of such simple cases is a sudden jump in velocity, from which a response time can be estimated. For the x/y-direction, the “sudden jump” case can be formulated as:

(6) $${v_{a,xy}} = \left\{ {\matrix{ {{v_{a,xy}}( - \infty )} \hfill & {{\rm for}\;t \lt 0} \hfill \cr {{v_{a,xy}}(\infty )} \hfill & {{\rm for}\;t \gt 0} \hfill \cr}} \right.,$$

i.e. a step function at t = 0, for which the solution for the particle velocity at t > 0 is:

(7) $${v_{p,xy}}(t) = {v_{a,xy}}(\infty ) + \displaystyle{1 \over {1/[{v_{a,xy}}( - \infty ) - {v_{a,xy}}(\infty )] - {\eta _{xy}}t}}.$$

From the analysis of the “sudden jump” case we can define the inertial time for the x/y-direction τ _I,xy as the time that is required such that the relative difference in velocity has decreased with exp(1):

(8) $$\displaystyle{{{v_{p,xy}}(t) - {v_{a,xy}}(\infty )} \over {{v_{a,xy}}( - \infty ) - {v_{a,xy}}(\infty )}} = \exp ( - 1),$$

from which follows:

(9) $${\tau _{I,xy}} = \displaystyle{{1 - \exp ( - 1)} \over {\exp ( - 1){\eta _{xy}}[{v_{a,xy}}(\infty ) - {v_{a,xy}}( - \infty )]}}.$$

It should be mentioned that such an inertial time is only a figure of merit and the reality is more complex. For the z-direction no analytical solution was found for the “sudden jump” case. It is however possible to find analytical solutions for more specific cases where the terminal fall speed is large or small compared with the order of magnitude of velocity differences due to relaxation ${{\rm \cal O}}(v_p^{\prime} )$ . For the first case with a small terminal fall speed, i.e. ${v_t}{\rm\,\ll {\cal O}}(v_p^{\prime} )$ , the solution is similar for the x/y-direction and the inertial time τ _I,z1 is:

(10) $${\tau _{I,z1}} = \displaystyle{{1 - \exp ( - 1)} \over {\exp ( - 1){\eta _z}[{v_{a,z}}(\infty ) - {v_{a,z}}( - \infty )]}}.$$

In case the terminal fall speed is large, i.e. ${v_t}{\rm \,\gg {\cal O}}(v_p^{\prime} )$ , the equation of motion for the z-direction can be written as:

(11) $$\displaystyle{{d{v_{p,z}}} \over {dt}} = \displaystyle{{dv_{p,z}^{\prime}} \over {dt}} \approx \mp 2{\eta _z}{v_t}v_{p,z}^{\prime} = \mp \displaystyle{{v_{p,z}^{\prime}} \over {{\tau _{I,z2}}}},$$

with:

(12) $${\tau _{I,z2}} = \displaystyle{1 \over {2{\eta _z}{v_t}}},$$

where τ _I,z2 is the inertial time for the z-direction. Again the sign is such that the relaxation velocity, v′ _p,z , is reduced, e.g. negative derivative, dv′ _p,z /dt, when v′ _p,z is positive. For the “sudden jump” case the analytical solution is:

(13) $$\eqalign{{v_{p,z}}(t) &= - {v_t} + {v_{a,z}}(\infty ) \cr & \quad + [{v_{a,z}}( - \infty ) - {v_{a,z}}(\infty )]\exp ( - \displaystyle{t \over {{\tau _{I,z2}}}}).} $$

Given the inertial times, it is also possible now to approximate the inertial distances by assuming that the particles are moving with the terminal fall speed in the z-direction for ${v_t}\,{\rm \gg {\cal O}}(v_p^{\prime} )$ , or otherwise moving with the final speed v _a,z (∞) − v _a,z (− ∞). The inertial distances are given in combination with an overview of all inertial parameters for rain drops in Table 1 Calculated values for the inertial distance are shown in Fig. 1. The inertial distance increases rapidly with the drop size. For large terminal fall speeds, ${v_t}\,{\rm \gg {\cal O}}(v_p^{\prime} )$ , the inertia parameters are reduced. In Fig. 1 the calculated inertial parameters are also calculated for spherical drops, which shows that the spheroidal shape has only a minor impact. From the equations of motions for a rain drop inertial parameters have now been estimated by making several assumptions. Without any additional information, the derived inertial parameters give a best guess for the time (distance) that a particle needs to respond to a change in the wind field. These assumptions include the simplification of the equations of motion to analytical solutions to the “sudden jump” case. The analytical solutions are different for direction, either x/y or z, and in case of the z-direction the solutions deviate for relatively large or small terminal fall speeds. An additional assumption is made that the particle maintains its shape, i.e. no microphysics (evaporation/condensation) and no vibrations are considered. Although such simplifications are rather crude, they are necessary for the development of a radar forward model that can be used in retrievals, where information is the limiting factor.

Fig. 1. The inertial distance d _I for rain drops is plotted as function of drop equivolumetric diameter D. The calculations are done both for spherical as spheroid droplets, for which the axis ratio relation from Beard and Chuang [Reference Beard and Chuang26] is used. Note that “z-dir, spherical, ${v_t}\,{\rm \ll {\cal O}}(v_p^{\prime} )$ ” is not in the figure because it is the same solution as “x/y-dir, spherical”.

Table 1. Inertial parameters for rain drops.

III. A NEW MODEL FOR INERTIA CORRECTION

Here we propose a new way on how to account for the inertia of rain droplets, that can be used in radar forward models of the radar Doppler spectrum or the radar DSW. The proposed method provides a correction for the observed variance of Doppler velocities, given the equivolumetric drop size or DSD and the observed radar resolution volume parameters.

In the case of homogeneous isotropic turbulence, the Kolmogorov hypothesis [Reference Kolmogorov27] states that within the inertial subrange the statistical representation of the turbulent energy spectrum S(k) is given by:

(14) $$S(k) = C{{\rm \epsilon} ^{2/3}}{k^{ - 5/3}},$$

where C is a Kolmogorov constant, ε the EDR and k the wavenumber. The variance of velocities, σ ², is then obtained as:

(15) $${\sigma ^2} = \int_{{k_{min}}}^{{k_{max}}} \!S(k)dk = \displaystyle{3 \over 2}C{{\epsilon} ^{2/3}}[k_{min}^{ - 2/3} - k_{max}^{ - 2/3} ],$$

where k is the wavenumber that is related to a length scale L via k = 2π/L. This leads to the proportionality between the variance of velocities and length scales:

(16) $${\sigma ^2} \propto [L_{max}^{2/3} - L_{min}^{2/3} ].$$

A correction ζ to the variance due to turbulence intensity, equation (16), can now be defined, which can be directly applied to the radar DSW:

(17) $$\zeta = {\sigma _{inertia}}/{\sigma _{no,inertia}},$$

where σ _inertia is the standard deviation of velocities for the inertial particles and σ _no,inertia is the standard deviation of velocities as if there was no inertia. Typically in a radar retrieval algorithm σ _inertia is measured and one would like to have σ _no,inertia , which can be directly related to the turbulence intensity. The term σ _no,inertia is related to perfect tracers of the air, as if the radar was looking to particles that were following the air motion exactly.

Based on equation (16) and the derived inertial parameters (Table 1) a new inertia correction model is proposed. The concept of the correction is twofold. Given a radar resolution volume scale L _*, both the maximal and the minimal scale for which the inertial particles are sensors are adjusted. The maximum scale is adjusted, as an ensemble of inertial particles can have more variations of velocities that are representative for a larger volume then the radar resolution volume. To do this a typical distance is defined, d _I,xyz , which is a 3-D distance on which the particle responds to variations in the wind field. On the other side, the inertial particles need a certain distance to respond to variations in the wind velocity, and also the smallest scale is adjusted. The inertia correction model can then be formulated as:

(18) $${d_{I,xyz}} = \sqrt {2^{\ast}d_{I,xy}^2 + d_{I,z2}^2}, $$

(19) $$\zeta _{xy}^2 ({L_{xy}},D) = \displaystyle{{{{({L_{xy}} + {\alpha _{xy}}{d_{I,xyz}})}^{2/3}} - {{({\beta _{xy}}{d_{I,xy}})}^{2/3}}} \over {L_{xy}^{2/3}}}, $$

(20) $$\zeta _z^2 ({L_z},D) = \displaystyle{{{{({L_z} + {\alpha _z}{d_{I,xyz}})}^{2/3}} - {{({\beta _z}{d_{I,z2}})}^{2/3}}} \over {L_z^{2/3}}}, $$

where L _* is the radar resolution volume scale for direction * and d _I,xyz is 3-D inertial distance. Here we introduce α _* and β _*, which are dimensionless tuning parameters that can be used to adapt this model based on additional information such as simulations or experiments. Note that for (α _* = 0,  β _* = 0) the result is that ζ = 1 and there is no correction. Here we choose d _I,z2, i.e. the solution for ${v_t}\,{\rm \gg {\cal O}}(v_p^{\prime} )$ , which seems legitimate for small turbulence intensities in precipitation. The alternative, d _I,z1, has a similar dependency on D, which can be seen in Fig. 1. Estimated uncertainties in α _* and β _* from simulations or experiments should be able to cover for such an introduced uncertainty. In Fig. 2 the model for inertia correction is shown for different adaptation parameters α _* and β _* for a drop of 0.5 mm, which demonstrates both effects. The parameter α _* determines the importance of the correction for large scales. In Fig. 2 it can be seen that an increased value of α _* leads to an increased value of σ _inertia , which can be attributed to the “transport of velocity fluctuations” in the radar resolution volume. The parameter β _* determines the relative importance of the classical inertia effect, i.e. the limited response to the small scale velocity fluctuations. In Fig. 2 it can be seen that an increased value of β _* leads to a decrease in the value of σ _inertia , which can be interpreted as little (/no) response of inertial particles to small scale fluctuations.

Fig. 2. The inertia correction ζ is plotted as a function of the observation length scale L _* for the direction * which is either x/y or z for a droplet with a diameter of D = 0.5 mm. The plots are for different tuning parameters α _* and β _*. (a) D = 0.5 mm, x/y-direction (b) D = 0.5 mm, z-direction.

The correction for any radar line of sight, relevant for application, with radar elevation θ is:

(21) $$\eqalign{& {\zeta ^2}({L_{xy}},{L_z},D,\theta ) = \mathop {\cos} \nolimits^2 (\theta )\zeta _{xy}^2 ({L_{xy}},D) \cr & \quad + \mathop {\sin} \nolimits^2 (\theta )\zeta _z^2 ({L_z},D).} $$

IV. SIMULATIONS OF RAIN DROP INERTIA

In this section we show how wind field simulations can be used to adapt the model for inertia correction. A 3-D homogeneous isotropic turbulence model from Mann [Reference Mann28] is used to estimate the adaptation parameters α _* and β _* in the inertia correction model (equation (19) and (20)). The simulated wind field and the terminal fall speed of droplets can be used to calculate backward trajectories of rain drops. Given the backward trajectories, the additional velocity difference due to relaxation can then be calculated by implicit integration of the equations of motion. The additional velocity difference due to relaxation, i.e. the inertia effect, can then be used to estimate the difference in variance of Doppler velocities within a radar resolution volume. The simulations will also give an idea of the uncertainties of EDR estimation from the DSW due to the inertia effect.

Here we use the work of Mann [Reference Mann28] to obtain a 3-D homogeneous isotropic turbulence periodic field with spatial scales of 30 × 30 × 30 m³, a spatial resolution of 1.2 m and an EDR of 0.01 m² s^–3, which is a typical value of daytime atmospheric EDR. The scale of 30 m is chosen because this is a typical maximum distance where inertia becomes less important (see Fig. 1). The parameters in the wind field simulation are chosen such that the simulated wind field satisfies the Kolmogorov energy spectrum for the inertial subrange, equation (14). Details of such simulations can be found in Mann [Reference Mann28]. The simulated wind fields are shown in the background of Fig. 3. Given the terminal particle velocity, the wind field and the assumption of no inertia, $\vec v_p^{\,\prime} = 0$ in equation (5), the backward trajectories of particles can then be calculated. In Fig. 3 it can be seen that for small rain droplets (0.01 and 0.1 mm) the trajectories and the origins have a random nature, whereas for larger rain drops (0.5 and 4 mm) the trajectories have an orchestrated nature and the drops mainly originates from above.

Fig. 3. Backward trajectories of rain drops are shown for different equivolumetric drop sizes D for the same 3-D wind field. The black dots represent the final locations of the drops, which are used in the simulations. In the background wind vectors used for the calculations are shown on the xy, xz, and yz-planes of a 3-D wind field simulation [Reference Mann28]. (a) D = 0.01 mm (b) D = 0.1 mm, (c) D = 0.5 mm (d) D = 4 mm.

To estimate the effect of inertia on the radar DSW, the following method is applied. On a line, either parallel to the x-direction or the z-direction, at a random position in the simulation, 100 equidistant points are taken for which the radar Doppler velocities are calculated, see Fig. 3. For the Doppler velocities the projected velocities in the line of sight are taken, e.g. u for the x-direction. Consequently the standard deviation of the Doppler velocities is calculated two times, with and without the velocity difference due to relaxation v′ _p,*. The standard deviations can be used to calculate the inertia factor ζ _*,i for direction *, which is specific for the realization of this wind field and specific for the random location of the line spanning 100 equidistant drop locations for simulation number i:

(22) $${\zeta _{*,i}} = {\sigma _{inertia}}/{\sigma _{no,inertia}},$$

where σ _inertia is the standard deviation of line of sight velocities with v′ _p,* included. For σ _no,inertia the relaxation term v′ _p,* is not included. An ensemble of such calculations can then be used to calculate optimal parameters in the model for inertia correction and to estimate the uncertainty in EDR from the DSW due to rain drop inertia.

For the estimation of the rain drop inertia effect, the velocity difference due to relaxation v′ _p,* needs to be calculated. First a backward trajectory is calculated, under the assumption that v′ _p,* is not too large and has a negligible influence on the trajectory. An advantage of this assumption is that the final position is not altered. Consequently the equations of motions are integrated over the trajectory to find the particle velocity difference due to relaxation v′ _p,*. The equations of motion, equations (1) and (2), can be used to write the differential equations that describe the evolution of the relaxation term v′ _p,*:

(23) $$\displaystyle{{dv_{p,xy}^{\prime}} \over {dt}} = \left\{ {\matrix{ { + {\eta _{xy}}v_{p,xy}^{{\prime}2}} \hfill & {{\rm for}\;v_{p,xy}^{\prime} \lt 0} \hfill \cr { - {\eta _{xy}}v_{p,xy}^{{\prime}2}} \hfill & {{\rm for}\;v_{p,xy}^{\prime} \gt 0} \hfill \cr}} \right.,$$

(24) $$\displaystyle{{dv_{p,z}^{\prime}} \over {dt}} = \left\{ {\matrix{ { + {\eta _z}\left[ {v_{p,z}^{{\prime}2} - 2{v_t}v_{p,z}^{\prime}} \right]} \hfill & {{\rm for}\;v_{p,z}^{\prime} \lt 0} \hfill \cr { - {\eta _z}\left[ {v_{p,z}^{{\prime}2} - 2{v_t}v_{p,z}^{\prime}} \right]} \hfill & {{\rm for}\;0 \lt v_{p,z}^{\prime} \lt \displaystyle{{{v_t}} \over 2}} \hfill \cr { + {\eta _z}\left[ { - 2v_t^2 + 2{v_t}v_{p,z}^{\prime} - v_{p,z}^{{\prime}2}} \right]} \hfill & {{\rm for}\;v_{p,z}^{\prime} \gt \displaystyle{{{v_t}} \over 2}} \hfill \cr}} \right..$$

Implicit integration of the differential equation:

(25) $$\displaystyle{{dy} \over {dt}} = a{y^2} + by + c,$$

has the following solution:

(26) $${y_{i + 1}} = - \displaystyle{b \over {2a}} + \displaystyle{1 \over {2a\Delta t}} \pm \sqrt {{{\left( {\displaystyle{b \over {2a}} - \displaystyle{1 \over {2a\Delta t}}} \right)}^2} - \displaystyle{{{y_i} + c\Delta t} \over {a\Delta t}}}, $$

which can be applied to find the following solutions of equations (23) and (24):

(27) $$v_{p,xy,i + 1}^{\prime} = \left\{ {\matrix{ { + {Q_1} - \sqrt {Q_1^2 + {Q_2}}} \hfill & {{\rm for}\;{U_{xy,i}} \lt 0} \hfill \cr { - {Q_1} + \sqrt {Q_1^2 + {Q_2}}} \hfill & {{\rm for}\;{U_{xy,i}} \gt 0} \hfill \cr}} \right.,$$

(28) $$v_{p,z,i + 1}^{\prime} = \left\{ {\matrix{ { + {Q_3} - \sqrt {Q_3^2 + {Q_4}}} \hfill & {{\rm for}\;{U_{z,i}} \lt 0} \hfill \cr { - {Q_3} + \sqrt {Q_3^2 + {Q_4}}} \hfill & {{\rm for}\;0 \lt {U_{z,i}} \lt \displaystyle{{{v_t}} \over 2}} \hfill \cr { - {Q_3} + \sqrt {Q_3^2 + {Q_5}}} \hfill & {{\rm for}\;{U_{z,i}} \gt \displaystyle{{{v_t}} \over 2}} \hfill \cr}} \right.,$$

with

(29) $${U_{xy,i}} = {v_{a,xy,i + 1}} - (v_{p,xy,i}^{\prime} + {v_{a,xy,i}}),$$

(30) $${U_{z,i}} = {v_{a,z,i + 1}} - (v_{p,i,z}^{\prime} + {v_{a,z,i}}),$$

(31) $${Q_1} = \displaystyle{1 \over {2{\eta _{xy}}\Delta t}}{\kern 1pt} ,\;{\kern 1pt} {Q_2} = \left \vert {\displaystyle{{{U_{xy,i}}} \over {{\eta _{xy}}\Delta t}}} \right \vert, $$

(32) $${Q_3} = {v_t} + \displaystyle{1 \over {2{\eta _z}\Delta t}}{\kern 1pt},\;{\kern 1pt} {Q_4} = \left \vert {\displaystyle{{{U_{z,i}}} \over {{\eta _z}\Delta t}}} \right \vert {\kern 1pt},\;{\kern 1pt} {Q_5} = {Q_4} - 2v_t^2. $$

Given the terminal fall speed v _t , rain drop parameters to calculate η _* and a 3-D wind field simulation is it now possible to integrate the relaxation term v′ _p,* Typically it is necessary to integrate over a distance on the order of several d _I,xyz to obtain a stable solution.

The results of the inertia correction factors ζ _*,i , that are calculated from the simulations are shown in Fig. 4. On the x-axis the length scale L _* is varied, which defines the total distance of the line that spans the rain drop locations. On the y-axis the rain drop inertia correction factor is shown. Simulations are done for multiple drop sizes and for different looking directions (x and z). In Figs 4(e) and 4(f) it can be seen that for larger droplets the correction becomes more important as it deviates more from 1. The smaller the spatial scale L _*, i.e. shorter radar resolution volume scale/higher spatial radar resolution, the more important the inertia correction becomes. With a minimization of a cost function (Nelder-Mead Simplex) for the inertia correction ζ, it is possible to obtain optimal parameters for α _* and β _* which are summarized in Table 2. The minimization is executed separately for the x- and z-directions, and for different drop sizes. In addition the optimization of inertial parameters is done for all the drop sizes at the same time. The Root Mean Square Error (RMSE) in ζ increases for larger droplets, which can be expected as large rain droplets behave more as inertial particles. In Fig. 4 it can be seen that the model is sufficiently adaptable for the different cases. Especially for large droplets (e.g. 4 mm) and high resolution L < 10 m the uncertainties become very large.

Fig. 4. Calculation of inertia factor ζ _*,i for different simulations i, that have a random position in the 3-D wind field simulation. The calculations are done for different equivolumetric rain drop diameters D and for different looking directions. On the x-axis the length of the line L _* is varied, on which 100 equidistant points are put to calculate the Doppler velocity variance. (a) D = 0.1 mm, x-direction (b) D = 0.1 mm, z-direction, (c) D = 0.5 mm, x-direction (d) D = 0.5 mm, z-direction, (e) D = 4 mm, x-direction (f) D = 4.0 mm, z-direction.

Table 2. Optimized tuning parameters α _* and β _* for the inertia correction model for rain drops with various diameters.

V. RADAR FORWARD MODEL

In this section the mathematical equations are given that are related to the calculation of the radar observables that can be used in a retrieval of the turbulence intensity. Here it is accounted for the fact that the contributions to the radar spectral width due to the terminal fall speed of droplets and turbulence are not independent. This issue is discussed in e.g. Fang and Doviak [Reference Fang and Doviak29]. The Doppler mean fall velocity v _D and DSW σ _D can be calculated as a sum over rain drop sizes:

(33) $${v_D} = {{\sum\limits_i N({D_i})\sigma ({D_i}){v_f}({D_i})\sin \theta} / {\sum\limits_i N({D_i})\sigma ({D_i})}},$$

(34) $$\sigma _D^2 = \displaystyle{{\sum\limits_i N({D_i})\sigma ({D_i})\left[ {{\zeta ^2}({D_i})\sigma _T^2 + {{({v_f}({D_i})\sin \theta - {v_D})}^2}} \right]} \over {\sum\limits_i N({D_i})\sigma ({D_i})}},$$

where N(D _i ) is the DSD, σ(D _i ) the radar cross section, D _i the rain drop equivolumetric diameter, v _f (D _i ) the rain drop fall velocity, θ the elevation angle and ζ the inertia correction, which is equal to 1 in case there is no inertia correction. The turbulence spectral width without inertia correction, σ _T , can be related to the EDR as [Reference White, Lataitis and Lawrence30]:

(35) $$\sigma _T^2 = \displaystyle{{C{{\epsilon} ^{2/3}}} \over {4\pi}} I,$$

where C is a Kolmogorov constant, ε the EDR and I an elaborate integral for which we refer to White et al. [Reference White, Lataitis and Lawrence30]. The equivalent reflectivity factor Z(v ₁, v ₂) in the velocity interval [v ₁, v ₂], which can be used to calculate the power Doppler spectrum, is equal to:

(36) $$Z({v_1},{v_2}) = \sum\limits_i N({D_i})\sigma ({D_i})f({D_i},{v_1},{v_2})$$

where f(v ₁, v ₂, D _i ) is the fraction of Doppler velocities that are between v ₁ and v ₂. It can be calculated with the inverse cumulative normal distribution function F ⁻¹(p|μ,σ) with mean μ and standard deviation σ as:

(37) $$\eqalign{f({v_1},{v_2},{D_i}) &= {F^{ - 1}}({v_2} \vert {v_f}({D_i})\sin \theta, {\zeta _T}({D_i}){\sigma _T}) \cr & \quad - {F^{ - 1}}({v_1} \vert {v_f}({D_i})\sin \theta, {\zeta _T}({D_i}){\sigma _T}).} $$

The advantage of this formulation is that the fall velocity and cross section are still in symbolic form and the choice for cross sections and fall velocities calculations is unconstrained. It is general enough to extend it to polarimetric radar observables and/or use Mie cross sections, e.g. cross sections from Mishchenko [Reference Mishchenko31]. However these new formulations do not solve the problem as the radar observables typically do not carry sufficient information on the DSD and additional constraints for the retrieval are required. An example of such constraints are linked DSD parameters, which are used in radar DSD retrievals [Reference Brandes, Zhang and Vivekanandan32]. However, constraints that are typically used for weather radars are limited to radars with a low elevation angle and cannot be applied in a general context. The application of the new radar forward model, including finding appropriate constraints for a specific radar configuration, is part of future research.