Stochastic simulation experiment to assess radar rainfall retrieval uncertainties associated with attenuation and its correction

As rainfall constitutes the main source of water for the terrestrial hydrological processes, accurate and reliable measurement and prediction of its spatial and temporal distribution over a wide range of scales is an important goal for hydrology. We investigate the potential of ground-based weather radar to provide such measurements through a theoretical analysis of some of the associated observation uncertainties. A stochastic model of range profiles of raindrop size distributions is employed in a Monte Carlo simulation experiment to investigate the rainfall retrieval uncertainties associated with weather radars operating at X-, C-, and S-band. We focus in particular on the errors and uncertainties associated with rain-induced signal attenuation and its correction for incoherent, non-polarimetric, single-frequency, operational weather radars. The performance of two attenuation correction schemes, the (forward) Hitschfeld-Bordan algorithm and the (backward) Marzoug-Amayenc algorithm, is analyzed for both moderate (assuming a 50 km path length) and intense Mediterranean rainfall (for a 30 km path). A comparison shows that the backward correction algorithm is more stable and accurate than the forward algorithm (with a bias in the order of a few percent for the former, compared to tens of percent for the latter), provided reliable estimates of the total path-integrated attenuation are available. Moreover, the bias and root mean square error associated with each algorithm are quantified as a function of path-averaged rain rate and distance from the radar in order to provide a plausible order of magnitude for the uncertainty in radar-retrieved rain rates for hydrological applications.


Introduction
According to Beven (2006), "the most important (problem in hydrology of the 21st Century) is providing the techniques to measure integrated fluxes and storages at useful scales". The hydrological flux of interest here is precipitation, in particular in its liquid form: rain. Accurate and reliable measurement and prediction of the spatial and temporal distribution of rainfall over a wide range of scales is an important goal for hydrology, because rainfall constitutes the main source of water for the terrestrial hydrological processes. Traditionally, rain gauges have been employed for that purpose. A fundamental shortcoming of rain gauges from a hydrological perspective, however, is the fact that they represent point measurements. This causes an inevitable trade-off between spatial representativeness and temporal resolution. The application of rain gauges in networks has long been considered a solution to the problem. However, from a hydrological point of view, many operational rain gauge networks are too sparse to provide rainfall information at a satisfactory spatial and temporal resolution (e.g., Wood et al., 2000;Berne et al., 2004). Denser networks, on the other hand, would generally be very impractical (and quite expensive as well). An additional issue is that most spatial interpolation procedures, e.g., geostatistical techniques such as kriging (e.g., Krajewski, 1987;Creutin et al., 1988;Schuurmans et al., 2007), generally lack the ability to capture the extreme rainfall variability found in nature. The interpolated rainfall fields are often much smoother than what is known concerning this variability from weather radar observations. Operational issues such as wind effects and (lack of) maintenance also strongly affect the performace of rain gauges (e.g., Neff, 1977;Sevruk, 1989;Habib et al., 1999;Steiner et al., 1999).
Ground-based weather radars are in principle well-suited to provide rainfall measurements for hydrological applications because: (1) they provide complete spatial and Published by Copernicus Publications on behalf of the European Geosciences Union. temporal coverage of extended areas from one single measurement site; (2) they allow rapid access for real-time hydrological applications, both concerning rainfall measurement and short-term forecasting; (3) their combined spatial and temporal resolution is generally higher than what can be obtained using rain gauge networks. However, radar is a remote sensing technique, which implies that weather radars measure the backscattered signal from rain in the air, rather than the distribution of rain rates at the ground needed for hydrological applications. The conversion of the radar reflectivities measured aloft to rain rates at the ground constitutes the observer's problem in radar hydrometeorology (e.g., Austin, 1987;Smith and Krajewski, 1993). Both the reflectivity measurements themselves and the radar reflectivity -rain rate conversion are prone to errors and uncertainties. Quantification of the observation uncertainties associated with rainfall retrievals from ground-based weather radars is a prerequisite for the assimilation of radar-retrieved rainfall fields in hydrological models.
The first attempts to use weather radar to estimate the spatial and temporal distribution of rainfall for hydrological applications date from the early 1970s (e.g., Battan, 1973). Although these early studies yielded promising results, the difficult (real-time) access to the radar data and the need to treat the data for error sources often little known to hydrologists, prevented radar to become a standard hydrological instrument in those early years. During the 1980s, geostatistical techniques were developed to combine the information from radars with that from networks of rain gauges, the type of information hydrologists were used to working with (e.g., Krajewski, 1987;Creutin et al., 1988). The basic idea of this approach is that rain gauges provide the "ground truth" at various points in the area of interest and that the radar data can be used to interpolate between the gauges. However, because rain gauges themselves are prone to several error sources, the concept of ground truth is questionable: "ground truth is the amount of rain that would have reached the ground if the rain gauge had not been there". Moreover, after adjustment of the radar data using rain gauge measurements (called "calibration" at the time (e.g., Collier, 1986)), several errors and inconsistencies remained which this approach was not able to resolve. As opposed to the largely statistical approach of the 1980s, the more physical approach to radar rainfall retrieval adopted since the 1990s considers the principle of radar measurements and the microstructure of rainfall in quite some detail (e.g., Smith et al., 1996;Andrieu et al., 1997;Creutin et al., 1997;Serrar et al., 2000;Sánchez-Diezma et al., 2000;Berne et al., 2005a,b;Delrieu et al., 2005;Berenguer et al., 2005). Another aspect is that rain gauges are no longer used to "calibrate" the radar images, but mainly for verification purposes. Recent developments in radar technology have also demonstrated the potential of doppler and polarimetric techniques for rainfall estimation (e.g., Bringi and Chandrasekar, 2001;Meischner, 2004). However, currently most operational weather radars still employ incoherent, non-polarimetric algorithms for rainfall retrieval.
Operational radar rainfall estimates can be affected by several sources of error and uncertainty: -The spatial and temporal variation of the rainfall microstructure (i.e. the properties of rain at scales smaller than the radar spatial resolution) and macrostructure (at scales larger than the radar resolution) is an important source of uncertainty in radar remote sensing of rainfall (e.g., Durden et al., 1998;Uijlenhoet et al., 2003a,b); -Rain-induced signal attenuation may cause significant underestimation of rainfall, whereas the uncertainties associated with its correction can also lead to overestimation (e.g., Hitschfeld and Bordan, 1954;Marzoug and Amayenc, 1994;Uijlenhoet, 2005a, 2006); -An accurate absolute (power) calibration of the radar is crucial for reliable rainfall retrieval; a faulty calibration may lead to biased reflectivity and rainfall estimates, particularly when combined with an attenuation correction scheme (e.g., Atlas, 2002); -A rain-induced film of water on the protective cover of a radar antenna, the radome, may cause attenuation of the radar signal and result in underestimation of radar reflectivities and rain rates (e.g., Germann, 1999); -Undesired echoes from mountains or buildings which are intercepted by the sidelobes or even by the main lobe of the radar beam -ground clutter -may erroneously be interpreted as rain and thus lead to local overestimation of rain rates; at the same time, partial shielding of the radar beam by such targets may lead to underestimation at ranges further away from the radar (e.g., Andrieu et al., 1997;Creutin et al., 1997); -The vertical profile of reflectivity, combined with a radar beam which climbs and expands as it moves away from the radar, may lead to a systematic rangedependent bias in radar rainfall estimates (e.g., Joss and Waldvogel, 1990;; -Vertical gradients in the refractive index of the atmosphere may cause the radar beam to bend away from or towards the earth's surface; anomalous propagation associated with temperature and humidity inversions may lead to extended areas with ground echoes, which may be falsely interpreted as regions of rainfall (e.g., Pamment and Conway, 1998).
A comprehensive treatment of these aspects is beyond the scope of this article (e.g., Zawadzki, 1984;Sánchez-Diezma et al., 2001). Here we focus on the second of the mentioned Hydrol. Earth Syst. Sci., 12, 587-601, 2008 www.hydrol-earth-syst-sci.net/12/587/2008/ error sources, namely rainfall retrieval uncertainties associated with rain-induced signal attenuation and its correction. Many operational radar networks across Europe operate at relatively short wavelengths (C-band, ∼5 cm), which may be severely attenuated in heavy rainfall (e.g., Delrieu et al., 1999). In addition, there has recently been an increased interest in high-resolution radars operating at even shorter wavelengths (X-band, ∼3 cm), in particular for urban hydrological applications. Such X-band radars are much less expensive than C-band radars, mainly due to the smaller antennas needed to achieve the same angular resolution. Hence, there is potential for the application of such radar systems in relatively dense networks (e.g., CASA, http://www.casa.umass. edu). It has been recognized for a long time (e.g., Atlas and Banks, 1951) that quantitative radar rainfall estimation at X-and C-band is seriously hampered by attenuation of the radar signal by precipitation along its path (as is illustrated by Fig. 1). Therefore, the adverse effects of attenuation on radar-retrieved rainfall fields need to be identified and corrected.
This article presents the results of a simulation experiment designed to investigate the rainfall retrieval uncertainties associated with weather radars operating in different widely used radio frequency bands. We focus in particular on the rainfall retrieval errors and uncertainties associated with raininduced signal attenuation and its correction for incoherent, non-polarimetric, single-frequency operational weather radars. This provides an extension and generalization of previous work (e.g., Uijlenhoet, 2005a, 2006), which concentrated on the radar meteorological aspects of attenuation correction. Section 2 presents the stochastic rainfall model used to perform the simulation experiment. The range profiles of rain rate, radar reflectivity, and specific attenuation generated using the simulated profiles of raindrop size distributions are presented in Sect. 3. In Sect. 4 the two rainfall retrieval algorithms employed to estimate rain rate profiles from the simulated range profiles of attenuated ("measured") radar reflectivity are discussed. A discussion of the resulting radar rainfall retrieval uncertainties, both as a function of the path-average rain rate and as a function of the distance from the radar, is presented in Sect. 5. Finally, Section 6 presents the conclusions of this work.

Stochastic rainfall range profile simulator
We have developed a stochastic simulator of range profiles of raindrop size distributions (DSD), which provides a controlled experiment framework to investigate the accuracy and robustness of various attenuation correction algorithms (Berne and Uijlenhoet, 2005a). This simulator was recently employed to quantify the influence of uncertainties concerning radar calibration, parameterization of the power-law relation between the radar reflectivity Z and specific attenuation k, and total path-integrated attenuation (PIA) estimates  (Berne and Uijlenhoet, 2006). Here we focus on the uncertainty in the retrieved rainfall profiles from simulated single frequency, incoherent and non-polarimetric radar systems operating at X-, C-and S-band (∼10 cm) associated with the spatial variability of rainfall (and the corresponding DSD) on scales between 25 m and 50 km. S-band is used as a reference against which to compare the other two frequencies, because the former is known to be virtually immune to attenuation. This work complements previous experimental results concerning the uncertainty associated with attenuation correction due to the spatial variability of the DSD along a range profile (e.g., Delrieu et al., 1999) by posing the problem in a Monte Carlo framework, allowing a quantitative analysis of this uncertainty.

Model formulation
The description of the stochastic model of range profiles of raindrop size distributions used for the controlled simulation Table 1. Mean, standard deviation, and scale of fluctuation [km] of N ′ = ln N t (with N t in m −3 ) and ′ = ln (with in mm −1 ) deduced from HIRE'98 data (07/09/1998 event) for moderate (4-s time step) and intense (2-s time step) rainfall parameterization.

Mean
Std experiments in this article largely follows that of Berne and Uijlenhoet (2006); it is summarized here for the sake of completeness. The model assumes that the local drop size distribution (DSD) can be described adequately by an exponential DSD with two parameters, N t (total drop concentration) and (inverse of a characteristic diameter), considered to be random variables: where N(D|N t , )dD denotes the drop concentration in the diameter interval [D, D+dD] given N t and . The latter are assumed to be jointly lognormally distributed. A plausible spatial correlation structure is introduced in the range profiles by assuming N ′ = ln N t and ′ = ln to follow a first order discrete vector auto-regressive process (e.g., Bras and Rodríguez-Iturbe, 1985): with where j is the distance index, ρ N ′ (1) the auto-correlation at lag 1 (idem for ′ ), ρ N ′ ′ (1) the cross-correlation at lag 1, and ǫ N ′ a Gaussian white noise process (idem for ′ ). Hence, C 0 and C 1 represent the covariance matrices at lags 0 and 1, respectively. The variances of the white noise processes ǫ N ′ and ǫ ′ are determined such that X is a second order stationary stochastic process. For a first-order vector autoregressive process, the auto-correlation functions are exponential: where r denotes the distance lag and θ the characteristic spatial scale, also known as the scale of fluctuation (Vanmarcke, 1983): According to Eq. (3), θ essentially represents the decorrelation distance, in this case defined as the distance lag where the autocorrelation of the process has decreased to e −2 .

Model parameterization
The stochastic model described above allows the repeated generation of range profiles of DSDs of equivolumetric spherical raindrops. The model is parameterized using measurements of DSD time series collected with an optical spectropluviometer during the HIRE'98 experiment in Marseille, France . We have determined two sets of model parameters, a 'moderate' rainfall parameterization for which we used a 3 h-period of the rain event that occurred on 7 September 1998, and an 'intense' rainfall parameterization that was fitted on a period of 45 min of high-intensity rainfall during the same event. Taylor's hypothesis of "frozen turbulence" with a constant velocity of 12.5 m s −1 , consistent with the wind speed estimate of Berne et al. (2004), is invoked to convert the measured DSD time series to DSD range profiles. This implies that we implicitly assume that the simulated range profiles are oriented parallel to the prevailing wind direction and that non-uniform beam filling in the transverse direction does not play a role.
Both for the moderate and for the intense rainfall parameterization the zero-lag cross-correlations between the fitted N ′ and ′ values are found to be negligible. Moreover, the scales of fluctuation θ for N ′ and ′ are very close and will be assumed equal in what follows (although this is not a requirement of the model). Therefore, the proposed stochastic rainfall model (Eq. 2) effectively reduces to a combination of two uncorrelated first-order vector auto-regressive processes (for N ′ and ′ ) with a common auto-correlation function. The total number of model parameters has now reduced to five: the mean and standard deviation of N ′ and ′ , and the scale of fluctuation θ. Their values are given in Table 1, for both parameterizations.
In order to simulate the radar rainfall retrieval process over hydrologically relevant scales, we generate DSD profiles with a total length of 50 km for the moderate rainfall parameterization and 30 km for the intense parameterization. The spatial resolution for the moderate rainfall parameterization is taken to be 50 m (corresponding to a 4-s time step) and that for the intense parameterization 25 m (i.e. a 2-s time step).

Profiles of bulk rainfall variables
The radar equation relates the received power to the properties of the radar, those of the target (i.e. raindrops) and the distance (range) between radar and target. At attenuating wavelengths (such as X-and C-band) the classical radar equation (e.g., Uijlenhoet, 2001) should be multiplied by an exponential factor accounting for the attenuation of the received signal due to rainfall present on the path between the radar antenna and the target (e.g., Battan, 1973): with where P r [W] is the mean power received from raindrops at range r, C is the so-called radar constant (which is a function of the employed wavelength and antenna size, among other factors), |K| 2 is a coefficient related to the dielectric constant of water (≈0.93), Z A [mm 6 m −3 ] is the attenuated radar reflectivity factor, Z [mm 6 m −3 ] is the actual radar reflectivity factor (simply called "radar reflectivity" from now on), k [dB km −1 ] is the specific (one-way) attenuation coefficient (called "specific attenuation" hereafter), and c=0.2 ln(10). All three bulk rainfall variables relevant for radar rainfall retrieval using incoherent, single frequency, non-polarimetric radar systems, namely Z, k and the rain rate R [mm h −1 ], are (weighted) integrals over the raindrop size distribution. The radar reflectivity Z [mm 6 m −3 ] is defined as where λ [cm] denotes the wavelength of the radar signal and σ B [cm 2 ] is the backscattering cross-section. Similarly, the specific one-way attenuation k [dB km −1 ] is defined as where σ E [cm 2 ] is the extinction cross-section. Finally, the rain rate R [mm h −1 ] is defined as where v [m s −1 ] is the raindrop terminal fall velocity in still air. Using the Mie scattering theory for spherical particles (van de Hulst, 1981) to calculate the scattering cross-sections σ B and σ E and Beard's parameterization (Beard, 1976) to calculate the drop terminal fall speeds, profiles of the bulk rainfall variables Z, k and R are easily derived from the DSD profiles generated using the stochastic simulator described above. The radar equation (Eq. 6) is subsequently employed to simulate the corresponding profiles of the attenuated ("measured") radar reflectivity Z A . Examples of generated radar reflectivity profiles for both rainfall parameterizations are shown in Fig. 1. The stochastic simulation model described above allows controlled experiments in a Monte Carlo framework to quantify the rainfall retrieval uncertainty associated with spatial rainfall variability for weather radar systems operating in different widely used frequency bands.

Rainfall retrieval algorithms
For incoherent, single frequency, non-polarimetric radar systems the observer's problem of radar hydrometeorology consists of solving the inverse problem posed by Eq. (6). This implies inverting Eq. (6), i.e. reconstructing the range profile of Z given that of Z A , and subsequently converting the retrieved Z-profile to a rain rate (R) profile. Clearly, this inverse problem is ill-posed as long as no constraints on the relations between the bulk rain variables Z, k and R are specified. In accordance with all previous investigations in this field (e.g., Marshall and Palmer, 1948;Smith and Krajewski, 1993;Haddad and Rosenfeld, 1997;Uijlenhoet, 2001), we postulate the power-law relations We study two widely used attenuation correction algorithms. The first (Hitschfeld and Bordan, 1954) is based on the assumption that the measured reflectivity in the first range bin (i.e. the one closest to the radar) is not affected by attenuation (or by a radar calibration error). Using an a priori powerlaw relation between radar reflectivity and specific attenuation (Eq. 10), the path-integrated attenuation affecting the second range bin is calculated. Subsequently, the measured reflectivity in the second range bin is corrected and, using the same power-law relation, the path-integrated attenuation up to the third range bin is calculated and corrected for. In this manner an iterative correction for attenuation is carried out in the direction from the radar antenna towards the region of interest. Therefore this type of algorithm is termed "forward". Hitschfeld and Bordan (1954) ("HB" hereafter) derived a closed-form analytical solution to this problem under the assumption that the apparent radar reflectivity Z A (r) is known along the entire range profile r (see Appendix A). Their solution is reformulated here to express the retrieved (attenuation-corrected) rain rate (R ′ ) in terms of the measured (attenuated) reflectivities (Z A ): www.hydrol-earth-syst-sci.net/12/587/2008/ Hydrol. Earth Syst. Sci., 12, 587-601, 2008 Appendix A provides a derivation of the HB equation from Eq. (6). The fact that the integral is between 0 and r shows that the HB algorithm is a forward algorithm. Note that the difference in the denominator of Eq. (11) can reach values close to 0, which renders the HB algorithm potentially highly unstable (Hitschfeld and Bordan, 1954). Such numerical instabilities will later be referred to as 'diverging' corrections. The second attenuation correction algorithm considered here has been developed to avoid such instability problems. It is based on the assumption that the path-integrated attenuation (PIA) to a certain fixed target (e.g., a building or a mountain) at a given range r 0 is known. In practice the attenuation to this target can be estimated for instance by comparing the reflectivity of the target before and during a rainfall event. In this case the same iterative attenuation correction procedure is employed but this time in the direction from the fixed target towards the radar antenna. Therefore this type of algorithm is often referred to as "backward". Marzoug and Amayenc (1994) ("MA" hereafter) presented the corresponding analytical solution, reformulated here to express the retrieved (attenuation-corrected) rain rate (R ′ ) in terms of the measured (attenuated) reflectivities (Z A ): where A 0 =A(r 0 ) equals the exponential factor in Eq. (6) evaluated at the range r=r 0 , accounting for the (two-way) PIA between the radar antenna and the reference target. Appendix A provides a derivation of the MA equation from Eq. (6) as well.
The fact that the integral in Eq. (12) goes from r to r 0 (with r 0 >r) shows that the MA algorithm is a backward al-gorithm. Also note that the minus sign in the denominator of Eq. (11) has now become a plus sign. Therefore, this type of algorithm is numerically stable by definition. The only disadvantage of backward algorithms with respect to forward algorithms is that they require reliable PIA estimates at ranges beyond the region of (hydrological) interest. In practice, such reference targets may not be available in all directions. Delrieu et al. (1997) were the first to apply this algorithm, which was originally developed for correcting (vertical) spaceborne radar rainfall profiles, to correct (horizontal) ground-based radar rainfall profiles (employing echoes from mountains to estimate the PIA).
As the rain-induced signal attenuation tends to zero (e.g., at S-band), the denominators of Eqs. (11) and (12) become negligible and both attenuation correction schemes reduce to which is simply the inverse of the power-law Z-R relation (Eq. 10). All three rainfall retrieval algorithms presented above are based on two important assumptions, namely (see Appendix A): (1) that the radar system to which they are applied is perfectly calibrated; (2) that the coefficients of the power-law relations between Z, k and R (Eq. 10) are constant over the region to which they are applied (implying a region of homogeneous rainfall, i.e. with the same type of rain).

Resulting uncertainties in radar rainfall retrievals
The uncertainty associated with radar rainfall retrievals based on the two attenuation correction algorithms presented above is studied in a Monte Carlo framework. We focus on three frequency bands that are widely used operationally: X-band (3.2 cm wavelength), C-band (5.6 cm), and S-band (10.0 cm). The latter is used as a reference, because it is known that attenuation is negligible at S-band for all but the most extreme rainfall. We generate one thousand profiles of N t and and calculate from those (using Eqs. 6-9) the corresponding profiles of the bulk rainfall variables Z, k, Z A , and R. To mimic the typical sampling resolutions of operational radar systems, the high spatial resolution (25 m, 50 m) profiles are averaged at a lower spatial resolution of 500 m. Table 2 lists some statistics of the generated profiles of the bulk rainfall variables Z, k, and R.
In previous attenuation correction sensitivity studies using the stochastic DSD range profile simulator Uijlenhoet, 2005a, 2006), we fitted a Z-k power-law relation on each profile separately using a non-linear regression technique. These relations necessarily constituted the best possible power-law relations for the generated profiles. This approach was adopted because we wanted to study the sensitivity of attenuation correction schemes to spatial rainfall variability (Berne and Uijlenhoet, 2005a) and other sources of uncertainty (Berne and Uijlenhoet, 2006) per se.

Z-R
Z-k X-band (3.2 cm) Z=233 R 1.59 Z=1.18×10 5 k 1.26 C-band (5.6 cm) Z=256 R 1.45 Z=6.57×10 5 k 1.11 S-band (10.0 cm) Z=311 R 1.40 Z=1.70×10 7 k 1.33 Here we approach the radar rainfall retrieval problem from an operational perspective. In practice, it would never be possible to have real-time estimates of the coefficients of the power-law Z-k and Z-R relations needed for radar rainfall retrieval at attenuating wavelengths, unless a network of instruments for measuring raindrop size distributions (disdrometers) would be deployed under the radar umbrella. However, this would not be feasible from an operational and financial perspective. Therefore, we employ climatological power-law Z-k and Z-R relations, whose coefficients α, β, γ , and δ (Eq. 10) are estimated following the procedure outlined by Delrieu et al. (1999), using a large dataset of DSD measurements collected in southern France (Table 3).
For the MA algorithm, we calculate the PIA value for each of the generated profiles (corresponding to A 0 in Eq. (12)) as the difference between the non-attenuated and the attenuated Z values at the final range bin. In other words, we assume the PIA estimates to be exact. The effect of an error in the PIA estimates on the accuracy of the MA algorithm was studied by Uijlenhoet (2005a, 2006).

Influence of the path-average rain rate
We have applied the two attenuation correction algorithms (Eqs. 11 and 12) to the 1000 Z A profiles using the climatological Z-k and Z-R relations. Because for hydrological applications the retrieved rain rate profiles are more relevant than the retrieved reflectivity profiles, we concentrate here on the former -the latter have been dealt with in previous work Uijlenhoet, 2005a, 2006). For each of the 1000 generated profiles, we have calculated two statistics quantifying the accuracy and uncertainty associated with the radar rainfall retrievals: the mean bias error (MBE) and the root mean square error (RMSE) between the retrieved and the actual rain rate profiles. Figures 2-5 show the 10%, 50% (median), and 90% quantiles of these statistics as a function of the profile-average rain rate for the three frequency bands and the two rainfall parameterizations considered.
For the moderate rainfall parameterization (Figs. 2 and  3), the path-average rain rates (averaged over profiles of 50 km length) are found to vary between a few and almost 20 mm h −1 . The (backward) MA algorithm significantly out-  performs the (forward) HB algorithm only at X-band frequencies for such moderate rain rates (Figs. 2 and 3, top panels). At C-band and even more prominently at S-band, the differences between the two attenuation correction algorithms are insignificant, given the appreciable amount of uncertainty associated with both error statistics caused by the statistical variability among the generated rainfall profiles within the moderate rainfall climatology. Moreover, the HB algorithm does not significantly diverge for any of the frequencies in case of moderate rain rates ("div"=0% on all occasions), not even at X-band, where the path-integrated attenuation is expected to be strongest.
Interestingly, the biases are almost always negative for the moderate rainfall parameterization. At X-band, the biases for the HB algorithm increase from about 20% of the path-average rain rate at 5 mm h −1 to more than 50% for path-average rain rates above 15 mm h −1 (Fig. 2, top panel). Therefore, even at moderate rain rates, where numerical instabilities do not seem to play a major role, the HB algorithm should be applied with great care at X-band. At C-and Sband, on the other hand, the biases tend to be limited to 15-20% of the path-average rain rate for both attenuation correction algorithms. The fact that there is a remaining negative bias at S-band (see bottom panel of Fig. 2), for which attenuation is negligible, indicates that the climatological Z-R relation used to retrieve the rain rate is not optimal for each individual profile considered, which results in the observed negative biases.
Although the general picture for the intense rainfall parameterization (Figs. 4 and 5) seems to be the same, the detailed results differ appreciably from those for the moderate rainfall parameterization. First of all, at X-band frequencies the HB attenuation correction algorithm now diverges in approximately one out of every five cases (18% of the profiles are numerically unstable). In addition, the bias remaining after attenuation correction using the HB algorithm exceeds 70% of the path-average rain rate for the most intense rainfall profiles, indicating a recovered fraction of the path-average rain rate of less than 30%. This clearly shows the complete failure of the HB algorithm for rain rate retrieval in intense rainfall at X-band, which is in accordance with previous results (e.g., Hitschfeld and Bordan, 1954;Delrieu et al., 1999;Uijlenhoet, 2005a, 2006).
The MA algorithm, on the other hand, is able to correct almost entirely for the suffered signal loss at X-band on average, perhaps even better than for the moderate rainfall parameterization (Fig. 4, top panel). One should bear in mind, however, that the total path length in this case is only 30 km, as opposed to 50 km for the moderate rainfall parameterization. Moreover, the uncertainty associated with the retrieved rain rate profiles, as quantified by the RMSE in Fig. 5 (top  panel) is appreciable for the most intense rainfall profiles, also for the MA algorithm. Interestingly, at C-band the MA algorithm seems to have a tendency to overcompensate for attenuation, which may be caused by the fact that the employed Hydrol. Earth Syst. Sci., 12, 587-601, 2008 www.hydrol-earth-syst-sci.net/12/587/2008/   climatological Z-k and Z-R relations are less appropriate at this frequency. At S-band, finally, both rainfall retrieval algorithms provide satisfactory results, although at this frequency the loss of power due to rain-induced attenuation is obviously not going to be a major source of error and uncertainty in the first place.
5.2 Influence of the distance from the radar In the previous section, the focus was on the influence of the rain rate (averaged along the considered profile) on the uncertainty associated with attenuation correction. For practical applications (e.g., rainfall-runoff modelling), it is also useful to investigate the dependence of this uncertainty on the distance from the radar. Similarly to Figs. 2-5, Figs. 6-9 present the MBE and RMSE values as a function of the distance from the radar, for the moderate and intense rainfall parameterizations, and for the three frequency bands considered.
Concerning the MBE, there is a median underestimation of a few mm h −1 for both attenuation correction algorithms at X-, C-and S-band, and for both rainfall parameterizations (Figs. 6 and 7). The underestimation is slightly larger for the HB than for the MA algorithm. Again, part of this underestimation is due to the inadequacy of the climatological Z-R relations, especially at S-band, for which rain-induced attenuation is negligible. In comparison with the MA algorithm, the spread between the 10% and 90% quantiles for the HB algorithm (i.e., including 80% of the values) is much larger at X-band, slightly larger at C-band and similar at S-band. This behavior is accentuated for the intense rainfall parameterization. It is consistent with Figs. 2 and 4, which indicate that HB performance degrades quickly when the mean rain rate along the profile increases. Because the specific attenuation is larger at X-than at C-or S-band, the uncertainty associated with the HB algorithm is larger at X-band, specifically for long distances from the radar.
With regard to the RMSE (Figs. 8 and 9), the results are comparable to those for the MBE. The median RMSE values for the two attenuation correction algorithms are close, except at long distances from the radar, where the pathintegrated attenuation is large. The interquantile range (a measure for the spread of the RMSE values) is much larger at X-than at C-or S-band. RMSE values are larger than zero even at S-band because of (1) the use of a climatological Z-R relation (which may induce a bias error), and (2) the inherent uncertainty associated with the use of a deterministic power law to describe a relation between stochastic variables (which induces a random error, e.g., Berne and Uijlenhoet, 2005a).
In summary, the distance from the radar is found to have a limited influence (in the order of a few mm h −1 ) on the median error associated with rain rate retrievals corrected for attenuation, but has a significant influence on the dispersion of this error, in particular for intense rainfall (the associated Hydrol. Earth Syst. Sci., 12, 587-601, 2008 www.hydrol-earth-syst-sci.net/12/587/2008/   errors can attain tens of mm h −1 ). This behaviour results from the higher probability to have larger attenuations and hence larger errors in the retrieved rain rates when the distance from the radar increases. It must be noted that the stochastic simulation approach employed in this study only allows an investigation of the errors associated with the studied attenuation correction algorithms and with the considered climatological Z-R and Z-k relations. The effects of beam broadening and of increasing beam altitude are beyond the scope of the present work.

Conclusions
We have presented a theoretical analysis of the observation uncertainties associated with rainfall estimates from groundbased weather radar. Rainfall being the main source of water for the terrestrial hydrological processes, accurate and reliable measurement and prediction of its space-time distribution over a wide range of scales is an important goal for hydrology. The rainfall retrieval uncertainties associated with weather radars operating in different widely used frequency bands have been investigated using a recently developed stochastic simulation model of range profiles of rainfall microstructure. To better mimic an operational setting in which no real-time disdrometer observations are available, we have employed climatological power-law Z-k and Z-R relations instead of the optimal relations that were used in previous studies Uijlenhoet, 2005a, 2006). A detailed comparison between two different attenuation correction schemes, the (forward) Hitschfeld-Bordan (HB) algorithm and the (backward) Marzoug-Amayenc (MA) algorithm, both in moderate (assuming a 50 km path length) and in intense Mediterranean rainfall (for a 30 km path length), shows that the backward correction algorithm is more stable and accurate than the forward algorithm, provided reliable estimates of the total path-integrated attenuation are available. For moderate rain rates, the HB algorithm is numerically stable for all wavelengths considered, although the biases remaining after correction can be appreciable at X-band (up to about 10 mm h −1 ). For such rain rates, the MA algorithm outperforms the HB algorithm only at X-band, whereas at C-and S-band both algorithms yield comparable results. In intense Mediterranean rainfall, however, the HB attenuation correction algorithm diverges in approximately one out of every five cases. Moreover, the remaining biases of those profiles for which the correction does not diverge are enormous (up to about 30 mm h −1 ). The MA algorithm on the other hand still performs satisfactorilly for all radar wavelengths considered. Finally, the dependence of the median values of MBE and RMSE on the distance from the radar is found to be limited for both attenuation correction algorithms. The associated spread, however, increases significantly as a function of distance, in particular for the HB algorithm.
Hydrol. Earth Syst. Sci., 12, 587-601, 2008 www.hydrol-earth-syst-sci.net/12/587/2008/ The present investigation has focused entirely on incoherent, single-frequency, non-polarimetric radar rainfall retrieval algorithms, which are still widely used in operational applications in hydrology and meteorology. However, the presented stochastic simulator of rainfall microstructure also provides a potential test bed for multi-parameter attenuation correction techniques (e.g., Illingworth et al., 2000;Testud et al., 2000;Vulpiani et al., 2005Vulpiani et al., , 2006a. Moreover, a stochastic rainfall model which explicitly treats the spatial and/or temporal variations of raindrop size distributions could also be relevant for several other hydrological applications, e.g. in models of rainfall interception by vegetation canopies (Uijlenhoet and Stricker, 1999) or soil erosion by raindrop impact , as well as the study of sampling uncertainties in in situ rainfall observations using rain gauges or disdrometers . Z A to the actual radar reflectivity factor Z (or the two-way attenuation factor A) at the reference range r 0 is known exactly, and 4) the apparent radar reflectivity Z A (r) is known along the entire range profile r. This implies among others the absence of noise of any kind (either rain-induced or sensor-induced), the absence of a minimum detectable radar signal (which would limit the range up to which the inversion equations could be applied successfully) and moreover, an infinitely high range resolution of the radar system in question.