Precipitation Type Classification of Micro Rain Radar Data Using an Improved Doppler Spectral Processing Methodology

Abstract

This paper describes a methodology for processing spectral raw data from Micro Rain Radar (MRR), a K-band vertically pointing Doppler radar designed to observe precipitation profiles. The objective is to provide a set of radar integral parameters and derived variables, including a precipitation type classification. The methodology first includes an improved noise level determination, peak signal detection and Doppler dealiasing, allowing us to consider the upward movements of precipitation particles. A second step computes for each of the height bin radar moments, such as equivalent reflectivity (Z_e), average Doppler vertical speed (W), spectral width (σ), the skewness and kurtosis. A third step performs a precipitation type classification for each bin height, considering snow, drizzle, rain, hail, and mixed (rain and snow or graupel). For liquid precipitation types, additional variables are computed, such as liquid water content (LWC), rain rate (RR), or gamma distribution parameters, such as the liquid water content normalized intercept (N_w) or the mean mass-weighted raindrop diameter (D_m) to classify stratiform or convective rainfall regimes. The methodology is applied to data recorded at the Eastern Pyrenees mountains (NE Spain), first with a detailed case study where results are compared with different instruments and, finally, with a 32-day analysis where the hydrometeor classification is compared with co-located Parsivel disdrometer precipitation-type present weather observations. The hydrometeor classification is evaluated with contingency table scores, including Probability of Detection (POD), False Alarm Rate (FAR), and Odds Ratio Skill Score (ORSS). The results indicate a very good capacity of Method3 to distinguish rainfall and snow (PODs equal or greater than 0.97), satisfactory results for mixed and drizzle (PODs of 0.79 and 0.69) and acceptable for a reduced number of hail cases (0.55), with relatively low rate of false alarms and good skill compared to random chance in all cases (FAR < 0.30, ORSS > 0.70). The methodology is available as a Python language program called RaProM at the public github repository.

1. Introduction

Precipitation is a key component of the hydrological cycle and a precise knowledge of the precipitating hydrometeor type is essential for remote quantitative precipitation estimates either from scanning or from vertically pointing ground-based or spaceborne radars. Precipitation observations of vertically pointing Doppler radars allow us to estimate the fall speed of hydrometeor particles which, in general, are the sum of their terminal fall speed and vertical air velocity (Atlas et al. [1], Hauser and Amayenc [2]). As it is well known, radar sensitivity to smaller particle detection increases with shorter wavelengths but, on the other hand, attenuation by intense precipitation, particularly rainfall, increases. Thus, a compromise exists between sensitivity and attenuation effects regarding the choice of operating frequency, a crucial aspect in radar design for different applications and platforms (Battaglia et al. [3], Kollias et al. [4]).

Unlike usual ground-based scanning precipitation weather radars operating at cm frequencies (S, C or X-band), Doppler radar profilers do not provide the precipitation field over a wide area, but a high spatial and temporal resolution vertical profile over the radar location. This approach has been employed for decades to study fine-scale vertical precipitation characteristics, for instance with the NOAA Aeronomy Laboratory S-band Doppler profiler (Ecklund et al. [5]), the X-band Precipitation Occurrence Sensor System (POSS, Sheppard [6]), the K-band Micro Rain Radar (MRR, Löffler-Mang et al. [7], Peters et al. [8]), and more recently, with shorter wavelength radars traditionally used for cloud studies, such as the Ka-band ARM zenith radar (Chandra et al. [9]) or the Milešovka observatory Ka-band cloud radar (Sokol et al. [10,11]), used to derive a hydrometeor classification, including four precipitation types. A related application has been the use of lidar observations from the NASA MPLNET network to resolve weak precipitation profiles (Lolli et al. [12,13]), particularly suited for light rain, drizzle and virga, as if more intense precipitation exists, attenuation becomes too important.

Among the above-mentioned Doppler radar profilers, the MRR stands out. This has been extensively used for a wide range of applications, including microphysical analysis of rainfall characteristics using collocated ground disdrometers (Adirosi et al. [14,15], Chang et al. [16], Gonzalez et al. [17], Jass et al. [18], Luo et al. [19], Tokay et al. [20]), diurnal and precipitation characteristics at low-latitude mountains (Bendix et al. [21], Seidel et al. [22]), orographic effects and low-level seeder-feeder processes (Arulraj and Barros [23]), bright-band (BB) radar signatures (Cha et al. [24], Brast and Markmann [25]), or the monitoring of absolute calibration of C-band polarimetric weather radars (Frech et al. [26]). Many of these applications rely on the separation of the liquid to solid precipitation phase, for example above and below the BB, which is crucial for accurate quantitative precipitation estimates (Fabry and Zawadzki [27], Sanchez-Diezma et al. [28], Bordoy et al. [29]). A recent study by Makino et al. [30] described the use of MRR to predict the hydrometeor type at ground level. However, to our best knowledge, no methodology based on MRR data to partition precipitation profiles in different hydrometeor types has yet been described. Such methodologies have been developed and improved in recent years for polarimetric scanning precipitation radars (see, for instance, Ryzhkov et al. [31], Park et al. [32], Schuur et al. [33], Dolan et al. [34], Chandrasekar et al. [35], Besic et al. [36]).

MRR units produce so-called raw data files, with spectral Doppler density data, which may be processed with the MRR manufacturer software (Metek [37], hereafter Method1) to obtain derived products, such as radar reflectivity Z, Doppler fall speed w or spectral width σ. Alternatively, researchers may develop their own MRR raw data processing, such as the methodology by Maahn and Kollias [38] which improves sensitivity and Doppler dealiasing and is suited for snowfall precipitation. This methodology is freely available as a Python program, called ImproToo (hereafter Method2).

The objective of this article is threefold. First, to present a novel MRR processing methodology, called RaProM (hereafter Method3), which includes enhanced spectral processing and Doppler dealiasing, a simplified novel hydrometeor classification scheme based on the assumption that the air vertical velocity is negligible compared to the precipitation particle speed, including drizzle, rain, snow, and hail, plus additional variables depending on the precipitation type. Second, to illustrate the methodology with a detailed case study, where results are compared with those from Method1 and Method2 and other co-located instruments, such as a Particle Size and Velocity (Parsivel) disdrometer or a microwave radiometer. Third, an evaluation of the hydrometeor type classification at the lowest MRR height bin compared with Parsivel present weather precipitation type observations recorded in 32 days.

The rest of the article is organized as follows. Section 2 introduces the instruments and location considered in this study, on the Eastern Pyrenees (NE Spain). Section 3 describes the new methodology proposed, dealing with spectral data processing with a multi peak detection procedure and details the derivation of spectral moments, hydrometeor classification and derived parameters. Section 4 presents a case study comparing, in detail, previous MRR processing methodologies and data from additional instruments and also provides the evaluation of the hydrometeor precipitation classification using contingency table scores and, finally, Section 5 presents a discussion and conclusion of the achievements and limitations of the proposed methodology and ideas for further research.

2. Instruments and Site Description

A Micro Rain Radar (MRR) is a compact frequency modulated continuous wave (FMCW) vertically pointing Doppler radar operating at 24.23 GHz, manufactured by Meteorologische Messtechnik GmbH (Metek) [37], which recently is manufacturing a newer version, the MRR-PRO. A summary of technical features of the MRR used here is given in

Table 1. Main features from MRR.

Frequency (GHz)	24.23
Radar Type	FMCW
Number of range gates	32
Number of spectral bins	64
Range resolution (m)	10–200
Frequency sampling (kHz)	125

. In this study, an MRR2 model was used and the range gate resolution was set to 100 m so observations extended up to 3.1 km above ground level. The unit was equipped with a heated antenna which prevented the accumulation of snow. Due to the relatively high operating frequency, possible attenuation by precipitation had to be checked for the data set analyzed.

Other instruments used here were a laser-optical disdrometer OTT Parsivel, a microwave radiometer RPG HATPRO (MWR) and two Automatic Weather Stations (AWSs) of the Meteorological Service of Catalonia [39]—see Supplementary Materials Table S1. All these instruments were deployed for the Cerdanya-2017 field campaign at the Das aerodrome (OACI code: LECD) in the Eastern Pyrenees mountain massif from December 2016 to April 2017. Part of the study makes use of additional data collected at the same site during the period 2018 and 2019 by the MRR and Parsivel disdrometer. The location is a relatively wide valley oriented west to east, with limited radar coverage due to orographic beam blockage (Bech et al. [40], Trapero et al. [41]) so most of the MRR beam cannot be directly compared with existing ground-based weather radar observations. The Cerdanya-2017 field campaign aimed at studying various complex terrain phenomena, including cold pool formation, mountain waves and orographic precipitation—see Gonzalez et al. [17]; Udina et al. [42] for more details.

The disdrometer records hydrometeor fall speed and size spectra at ground level and other derived variables, such as hydrometeor type (e.g., rain or snow), radar reflectivity factor and precipitation intensity for liquid precipitation, among others. The microwave radiometer provides air temperature vertical profiles so the freezing level and other isotherm levels can be calculated. The two AWS are located close or near the MRR location (at the same aerodrome and nearby but a higher altitude, see Figure S1 in Supplementary Materials) and provide independent measurements of temperature, precipitation and snow depth level.

3. New Methodology Proposed

The methodology has two different sections. The spectral data processing is detailed in the first section, and the equations and hypotheses are treated in the second section.

3.1. Spectral Data Processing

The initial processing stage of Method3 consists in the MRR spectral data processing and follows the flow chart described in Figure 1. The first step is to transform the original signal backscattered by hydrometeors to spectral reflectivity (ƞ) following the equation proposed by the manufacturer (Metek [37]):

$$\eta \left(n,i\right)=f\left(n,i\right)·\frac{i^2}{TF\left(i\right)}·\frac{C·\Delta h}{{10}^{20}} ,$$

(1)

where i is the range gate number (i = 0, …, 31), n is Doppler bin number (n = 0, …, 63), f(n,i) is the original MRR signal saved in the so-called raw data files, TF(i) is a transfer function specific for each height, C is the radar calibration constant and Δh is the range resolution in m. The spectral reflectivity ƞ(n,i) has units of m⁻¹ and TF(i) and C are stored in the original raw data files.

The innovations introduced in Method3 are the modification of the integration time, the average of spectra and the determination of maxima (peaks) in the signal. The integration time, usually set to 60 s, is now selectable by the user. Method3 performs the averaging of the spectrum checking that at least 50% of spectra contain a minimum valid signal. The 50% threshold is a value modifiable by the user on the code of Method3. The signal is considered valid if it verifies the Hildebrand and Sekhon [43] criterion:

$$\frac{{\overline{\eta }}^2}{var\left(\eta \right)} = \frac{{\left(\frac{\sum \eta \left(n\right)}{N}\right)}^2}{\frac{\sum {\left(\eta -\overline{\eta }\right)}^2}{N}}<60 ,$$

(2)

where var(ƞ) is the spectral reflectivity variance and N the number of Doppler bins considered.

The detection of the maxima of the signal at a given height consists in using all Doppler bins of the spectrum except the first and the last one, therefore all possible fall speed values are used, except the lowest and the highest. This approach allows us to recover more than one single peak of the signal more easily and provides more detail about the hydrometeor fall speed distribution potentially affected by noise in the lowest range bin, the closest one to the ground, which is excluded, as in Method1. Note that Method2, instead, excludes part of the highest and lowest Doppler bins, and discards the 3 lowest-range bins. The noise level is calculated for each height following the scheme detailed by the manufacturer [37], which implies subtracting the noise level to the signal.

Figure 2 shows two examples of noise subtraction and peak detection using Method3 applied to the lowest processed range bin (i = 1, here from 100 to 200 m a.g.l.) and the ninth range bin (i = 10, from 900 to 1000 m a.g.l.). Figure 2a,b display a single peak detection and Figure 2c,d a two peak detection. Note that these detections could not be possible with Method2, as it discards the 100 m to 200 m bin height (panels a and b) and does not consider all Doppler bins processed by Method3.

The next step is the conversion of the reflectivity spectra from n Doppler bins to velocity v according to Equation (3). Note that this transformation depends on specific MRR2 data acquisition features, such as wavelength or sampling frequency:

$$\eta \left(v,i\right) = \eta \left(n,i\right)·{\left(\frac{f_{sampling}}{2·n_{max}·i_{max}}·\frac{\lambda }{2}\right)}^{-1},$$

(3)

where f_sampling is 125 kHz, n_max is 64, i_max is 32 and λ is the wavelength (~1.24 cm).

The last step in the spectral processing of Method3 implements a dealiasing scheme, partly based on Kneifel et al. [44] who noticed that in some snowfall cases, MRR2 provided unrealistically high values of hydrometeor fall speed, much higher than terminal snowflake fall speeds. This was due to two assumptions of the original manufacturer software: (i) precipitation observed was always in liquid form, so there was an inherent dependence between hydrometeor terminal fall speed versus particle diameter and (ii) only downward velocities were allowed, which is the most usual situation but is not always the case for snowflakes or convective rainfall. The snowfall case is discussed in Maahn and Kollias [38] and a dealiasing system to solve it is implemented in Method2. The dealiasing scheme proposed in Method3 consists in, for a given height, combining information from spectra of adjacent (upper and lower) height levels to determine the hydrometeors’ fall speeds to extend the original speed range from 0 to 12 m/s to a dealiased range of −12 to 24 m/s. The lower-level spectra are used to expand the speed range to −12 to 0 m/s and the upper one to 12 to 24 m/s. The vertical continuity of the speed profile is used to provide the dealiased speed spectra, as illustrated in the example displayed in Figure 3. The Method2 dealiasing scheme works well for snowfall and finds possible snowflake upward movements but does not work properly for cases of intense rainfall with hydrometeor falling velocities greater than 8 m/s, unlike the proposed scheme implemented in Method3.

3.2. Parameters Calculation

Once Method3 spectral processing is completed, data are ready to proceed with the calculation of subsequent parameters, which is divided in two parts. The first part computes the basic Doppler radar parameters: the Doppler velocity, which is assumed here to be the hydrometeor fall speed, and the radar equivalent reflectivity factor. The second part provides an estimation of hydrometeor type and based on this, the calculation of derived parameters, such as precipitation rates for different precipitation types. Figure 4 provides an overview of parameters calculations.

3.2.1. Basic Parameters

The term basic parameters are applied here to integral parameters which can be calculated independent of the type of hydrometeor as they only depend on the spectral reflectivity. The parameters are the radar equivalent reflectivity (Z_e), the Doppler velocity, which here corresponds to the mean fall speed of hydrometeors ($\overline{w}$), and higher order moments of the Doppler speed distributions: Doppler spectral width (σ), skewness and kurtosis as described in Equations (4)–(8):

$$Z_e = {10}^{18}·\frac{{\mathsf{\lambda }}^4}{{\mathsf{\pi }}^5}·\frac{1}{{\left|K\right|}^2}·\Delta v·\sum \eta \left(v\right),$$

(4)

$$\overline{w} = \frac{\sum \eta \left(v,i\right)·v\left(i\right)}{\sum \eta \left(v,i\right)}$$

(5)

$$\sigma = \sqrt{\frac{\sum \eta \left(v,i\right)·{\left(v\left(i\right)-w\right)}^2}{\sum \eta \left(v,i\right)}}$$

(6)

$$skewness = \frac{\sum \eta \left(v,i\right)·{\left(v\left(i\right)-w\right)}^3}{\sum \eta \left(v,i\right)·{\sigma }^3}$$

(7)

$$kurtosis = \frac{\sum \eta \left(v,i\right)·{\left(v\left(i\right)-w\right)}^4}{\sum \eta \left(v,i\right)·{\sigma }^4}$$

(8)

Note that the calculation of the radar equivalent reflectivity does not take into account possible attenuation effects which may be relevant for high precipitation rates, considering that MRR operates at the K-band. However, this can be handled in the case of liquid hydrometeors, where the path attenuation is calculated (see Section 3.2.4).

3.2.2. Hydrometeor Type Classification

As mentioned earlier, there are a number of hydrometeor classification algorithms developed mainly for scanning polarimetric weather radars, some considering up to 10 different precipitation species. Here, a simplified approach is adopted aiming to distinguish, for each height bin, 5 possible precipitation types: drizzle, rain, snow, mixed and hail. For the purpose of this paper, we consider either wet snow, a mixture of snow and rain, or graupel in the mixed category.

The classification is based on a decision tree, considering empirical relations between hydrometeor fall speed and equivalent radar reflectivity, size and particle diameter characteristics for different hydrometeors and the existence or absence of the bright band. As a starting point, the empirical relations reported by Atlas et al. [1], linking radar reflectivity and fall speeds of rain (v_Rain) and snow (v_Snow) in the absence of bright band, are considered:

$$v_{Rain} = 2.65·Z{e}^{0.114}$$

(9)

$$v_{Snow}=0.817·Z{e}^{0.063}$$

(10)

These relationships are used to compute, given Z_e, v_Rain and v_Snow, which are the average fall speed expected for each precipitation type on those two cases. Additional parameters considered are the mean Doppler fall speed $\overline{w}$, the Doppler spectral width σ and the Skewness Sk, calculated for each height bin. Moreover, a BB detection scheme is used (Cha et al. [24]). In case there is BB, then the existence and height of its top (BB_Top) and bottom (BB_Bottom) levels are computed, following the methodology described by Wang et al. [45]. From all the above, the decision tree can be grouped in three main branches with additional conditions (Figure 4):

• If v_Snow is within the interval $\overline{w}$ ± σ and v_Rain exceeds $\overline{w}$ + σ, then:
  • If the bin height is lower than the BB_Bottom, the hydrometeor is classified as:
    Drizzle/Rain—Hail.
  • If the bin height is equal or above the BB_Bottom or BB_Bottom is not present:
    Mixed: if Sk > −0.5 and $\overline{w}$v_Snow;
    Snow: otherwise.
• If v_Rain and v_Snow are within the interval $\overline{w}$ ± σ, then:
  • If the bin height is below the BB_Bottom or BB_Bottom is not present:
    Drizzle/Rain—Hail.
  • If the bin height is above the BB_Bottom:
    Mixed: if the Sk > −0.5 and the $\overline{w}$v_Snow;
    Snow: otherwise.
• If v_Rain is within the interval $\overline{w}$ ± σ and v_Snow is lower than $\overline{w}$ − σ, then:
  • If the bin height is below the BB_Top or BB_Top is not present:
    Drizzle/Rain—Hail.
  • If the bin height is above the BB_Top:
    Mixed: if the Sk > −0.5 and the $\overline{w}$v_Snow;
    Snow: otherwise.
• Cases not included in any of the previous categories are labelled as unknown.

The category “Drizzle/Rain—Hail” is further disaggregated considering additional conditions (Figure 5). The basic criteria stem from the precipitation hydrometeor definitions. Drizzle and rain are formed only by liquid particles AMS (2020) [46,47] and are here distinguished by their skewness Sk in their fall speed distribution and for Z_e differences (ΔZ_e) between levels. According to Acquistapace et al. [48], if skewness (Sk) is lower than or equal to −0.5 and ΔZ_e ≥ 1 dBZ, then the hydrometeor is classified as drizzle and otherwise as rain. Hail is defined in the function of the maximum diameter on the Doppler velocity spectrum, considering here a threshold of maximum diameters greater than 5 mm. Snow and mixed class are stratified when the height is below the BB_Bottom as described by Kalesse et al. [49].

3.2.3. Snowfall Rate

Despite solid precipitation presents a greater variability than liquid precipitation, Matrosov and Heymsfield [50] studied the relation between equivalent radar reflectivity and snowfall rate at different wavelengths and proposed empirical relations between those variables. In particular, snowfall rate (SR) can be estimated from Z_e by inverting the Z_e–SR power–law relationship:

$$SR={\left(\frac{Z_e}{a}\right)}^{1/b}$$

(11)

where SR is in mm/h, Z_e in mm⁶ m⁻³, a and b are the coefficients from the corresponding Z_e–SR relation and their values, for K band, are 56.00 and 1.20. It should be noted that the estimated SR might differ from actual values, given the high variability of the mass–size relation of different snow particles, as discussed in Souverijns et al. [51].

3.2.4. Rainfall Parameters from Drizzle/Rain

Rainfall parameters can be calculated if hydrometeors are in liquid phase. Section 3.2.2 details that the Drizzle/Rain types are liquid hydrometeors, thus on these types the rainfall parameters are calculated, but it is also necessary to introduce a dependence between the hydrometeor terminal fall speed and the diameter of the hydrometeor, which implies the hydrometeor particle size distribution (N) for each n Doppler bin and height level i:

$$N_{n,i}\left(D_n\right) = \frac{\eta \left(D\right)}{{\sigma }_{n,i}}$$

(12)

where ƞ(D) is the spectral reflectivity as a function of the diameter and σ_n,i is the Mie backscattering cross section for liquid spherical particles. Note that the initial spectral reflectivity is only a function of the fall speed using the relation between Doppler bins n. The spectral reflectivity as a function of the raindrop diameter is determined using the Gunn and Kinzer [52] expression, which relates the fall speed with the raindrop diameter plus a correction factor for the fall speed δv(h) that takes into account air density changes with height:

$$\eta \left(D,i\right)\left[{\mathrm{ms}}^{-1}{\mathrm{mm}}^{-1}\right] = \eta \left(v,i\right)·6.18·\delta v\left(i·\Delta h\right)·e^{-0.6 {\mathrm{mm}}^{-1}·D\left[\mathrm{mm}\right]}$$

(13)

where D is expressed in mm and ƞ(D,i) is in m s⁻¹ mm⁻¹. The correction factor δv(h) is computed assuming the US Standard Atmosphere and a second order approximation following Foote and Du Toit [53]

$$\delta v\left(h\right) = \left(1+3.68·{10}^{-5}·h+1.71·{10}^{-9}·h^2\right)$$

(14)

So, the corrected terminal fall speed as a function of drop diameter and height is:

$$v\left(D\right)\left[\mathrm{m}/\mathrm{s}\right] = \delta v\left(h\right)·\left(9.65-10.3·e^{-0.6 {\mathrm{mm}}^{-1}·D\left[\mathrm{mm}\right]}\right) \mathrm{for} 0.109 \mathrm{mm} \le \mathrm{D} \le 6 \mathrm{mm}$$

(15)

After the drop size distribution is determined, the Path Integrated Attenuation (PIA) is determined using the approach detailed in Metek [37], wherein the single particle extinction coefficient is necessary. The result is the drop size distribution with attenuation correction (Na(D,i)) which allows us to calculate the reflectivity (Z), the liquid water content (LWC), and the rain rate (RR):

$$Z = {\displaystyle \sum }N{a}_n·D_n{}^6·\Delta D$$

(16)

$$LWC = {\rho }_w·\frac{\pi }{6}·{\displaystyle \sum }N{a}_n·D_n{}^3·\Delta D$$

(17)

$$RR = \frac{\pi }{6}·{\displaystyle \sum }N{a}_n·D_n{}^3·v\left(D\right)·\Delta D$$

(18)

Following Thurai et al. [54], Method3 implements the calculation of the mean mass-weighted raindrop diameter (D_m) and the intercept parameter of the gamma distribution normalized to the liquid water content (N_w), where it is assumed that D₀ is equal to D_m:

$$D_m = \frac{\sum N{a}_n·D_n{}^4·\Delta D}{\sum N{a}_n·D_n{}^3·\Delta D}$$

(19)

$$N_w = \frac{256}{\pi ·{\rho }_w}·\frac{LWC}{D_m^4}$$

(20)

These parameters are useful for discriminating between convective and stratiform rainfall, as discussed later.

4. Results

The results are divided into two parts. The first part examines a case study to assess the characteristics of Method3, compared with Method1 and Method2. In the second part, an objective validation of the hydrometeor classification is performed to show the performance of Method3 to distinguish different types of precipitation.

4.1. Case Study

The performance of the new methodology proposed (Method3) is assessed during a precipitation event, mostly stratiform, that took place on 27 March 2017. The event produced rainfall at Das AWS (4.9 mm) and 6.3 mm of equivalent rainfall amount, which fell as snow, at Malniu AWS located 1100 m aloft. These precipitation amounts are relatively modest in terms of daily amounts for the season and region (Gonzalez and Bech [55]). As the freezing level was about 750 m above ground, a substantial part of the profile observed was snow, which allows us to illustrate different features of Method3. The results are compared with Method1, Method2 and with data from other instruments—see Supplementary Materials Table S1. Comparisons may include different data subsets; for example, Method2 vs. Method3 profiles or Method3 lowest bin gate vs. disdrometer estimates.

4.1.1. Fall Speed

Precipitation fall speed profiles estimated with Method3 are displayed in Figure 6, overlaid with isotherm heights (0 °C, −10 °C and −20 °C levels), plus Parsivel fall speed at ground, indicating the consistency between the observations derived from the three independent instruments. On the one hand, a sharp increase in Method3 fall speed is generally observed below the 0 °C level, as expected when solid precipitation changes to liquid precipitation. On the other hand, the lowest Method3 height bin (100 m a.g.l.) presents fall speed values and temporal trends comparable to the Parsivel ones.

The lowest range bin (100 m height) estimated with Method3 is compared with Parsivel measurements (Figure 7a). Some discrepancies are expected due to the different measurement principle of each instrument and also the different heights compared. Both data sets compare reasonably well for speeds up to 8 m/s, and particularly well for 3 to 6 m/s fall speeds (note the higher density of data for that range, as shown in Figure 7a). All Method3 speeds are below 8 m/s, whereas disdrometer data exceed 12 m/s, so major discrepancies occur for disdrometer speeds above 8 m/s. This discrepancy could be partly due to raindrop coalescence in the lowest 100 m above ground level.

An analysis of the first three height bins obtained with Method3 is performed with a comparison with Method1 (Figure 7b). Note that Method2 cannot be included in this comparison as it discards the two lowest height bins. All negative fall speeds detected by Method1 are discarded according to the manufacturer approach, which assumes only positive values. The agreement found is generally good, with a few cases where Method3 overestimates Method1.

A comparison of profiles except for the lowest two bins not processed by Method2 is performed between Method2 and Method3 (Figure 7c). As displayed by the data density the fall speed agreement is generally good (R² of 0.995); however, a few cases present discrepancies, typical of the order of 1 m/s. Figure 7d shows a scatter plot between Method2 and Method3 where data are labelled to hydrometeor type classified by Method3 and illustrates that the largest discrepancies shown in Figure 7c are originated by snow and mixed cases.

The fall speed difference between the two methods, expressed here as:

$$\Delta w = w_{Method3}-w_{Method2},$$

(21)

is further examined in terms of hydrometeor type according to the new hydrometeor classification methodology performed by Method3 and fall speed range, considering the different speed classes (

Table 2. Fall speed differences Δw between Method3 and Method2 stratified by speed values and hydrometeor type. Hail, mixed and unknown values were not found so only four hydrometeor classes (Rain, Drizzle, Mixed and Snow) are considered. Values listed correspond to number of height bins for each fall speed and hydrometeor class and mean error (ME) and root mean squared error (RMSE).

(m/s)	Rain	Drizzle	Mixed	Snow
$\left|\Delta w\right|<1$	1142	88	1057	7441
$1\le \left|\Delta w\right|\le 2$	0	0	0	4
$\left|\Delta w\right|>2$	0	0	2	1
ME	−0.01	−0.02	0.00	0.01
RMSE	0.06	0.03	0.16	0.08

). It can be seen that, for all classes, the absolute value of the mean error is equal or lower than 0.02 m/s. Moreover, a few snow and mixed cases present speed differences above 1 m/s. Root mean square errors are similar for all hydrometeor classes.

The distinct behavior of snow cases is also seen in Figure S2 in the Supplementary Materials, which shows the difference in distribution for the four hydrometeor types, shown in Table 2. Rain and mixed cases present similar quasi symmetric distribution patterns, while drizzle is much more leptokurtic and snow is platykurtic. The systematic differences between Method2 and Method3 may be due to differences in the spectral processing of the methods. However, note that these differences are very small in absolute value.

4.1.2. Equivalent Reflectivity

Equivalent reflectivity (Z_e) profiles obtained with Method3 are displayed with selected temperature levels (0 °C, –10 °C and –20 °C) retrieved from the microwave radiometer (MWR) and the reflectivity observed at ground-level by the disdrometer (Figure 8). It can be seen that the 0 °C level, around 750 m A.G.L., matches approximately with an abrupt increase in Z_e consistent with a bright band signature caused by the change from solid to liquid hydrometeors and an increase in the fall speed (shown previously). The disdrometer reflectivity is also consistent with the profiles and reproduces particularly well the timing of the local maxima (>25 dBZ).

A comparison of the reflectivities provided by Method3 at the lowest bin and the disdrometer is shown in Figure 9a, which indicates an overall agreement but a slight overestimation of the disdrometer compared to Method3 for values lower than 20 dBZ. The three lowest height bins provided by Method1 and Method3 present generally similar values (Figure 9b); however, some discrepancies are found in a few cases, mainly due to Method3 overestimating Method1. These variations appear because Method1 does not apply the dealiasing. Figure 9c shows a scatter plot comparing Method3 and Method2, which also show a good global agreement between the two methods (R² of 0.993), except in a few cases.

Those discrepancies are examined in more detail by considering the difference ΔZ_e between the two methods:

$$\Delta Z_e=Z_e{}_{Method3}-Z_e{}_{Method2}$$

(22)

Aside from a few isolated snow and mixed cases, both snow and rain provide similar differences in terms of RMSE (~1 dB) and ME (~–0.40 dB), so Method3 provides slightly lower reflectivity values than Method2—see

Table 3. As Table 2 but for reflectivity differences ΔZ_e between Method3 and Method2.

(dBZ)	Rain	Drizzle	Mixed	Snow
$\left|\Delta Z_e\right|<1$	1003	88	1023	6518
$1\le \left|\Delta Z_e\right|\le 5$	135	0	32	914
$\left|\Delta Z_e\right|>5$	4	0	4	14
ME	−0.38	−0.01	−0.14	−0.45
RMSE	1.28	0.04	0.75	0.80

. Despite these similarities, the distribution of differences present distinct patterns, as highlighted in Figure S3 in the Supplementary Materials: they all present a mode value close to 0 dB but it is much less marked for snow than for others types.

4.1.3. Hydrometeor Classification

The evolution of the hydrometeor classification provided by Method3 is shown in Figure 10, overlaid with temperature levels (0, –10, –20 °C) from radiometer and Parsivel hydrometeor classification at ground level. Parsivel classification is derived from the World Meteorological Observations standard code 4677 used in surface synoptic observations (SYNOPs), grouped here as rain, drizzle, mixed, snow and hail (see details in Supplementary Materials Table S2).

Figure 10 shows clearly that the 0 °C level is slightly above the rain level, which may be explained by the local cooling caused by heat exchange with the environment due to snow melting. From 12 to 15 UTC the freezing level increases, as does the rain level, and later decreases, a trend also followed by the rain level. It can also be seen that the disdrometer hydrometeor classification detects rain and drizzle, consistently with Method3. About 17:30 UTC a short precipitation event with liquid precipitation over the freezing level is observed. In this case, no bright band was detected and both MRR Doppler fall speed and spectrum width did not change along the precipitation profile substantially (not shown), which is consistent with the fact that it was a brief shallow convective event.

The Malniu station, located about 1100 m above ground level from Das, recorded mostly air temperatures below 0 °C and an increase in snow depth at the time the precipitation occurred (Figure S4 in Supplementary Materials), confirming that, at that height, precipitation was falling as snow, as indicated by the Method3 hydrometeor classification. On the other hand, Das station—where, initially, there was no snow on the ground—did not record any increase in snow height, as expected for a rain event.

4.1.4. Rain Rate

Rain rate obtained from the lowest bin height (100 m a.g.l.) provided by Method3 is compared to ground level rain rate, calculated with disdrometer and AWS data located at Das. Disdrometer rain rate was calculated using the raindrop number concentration per unit volume of air following Friedrich et al. [56]. The smallest raindrops detected by the disdrometer had diameters about 0.312 mm, which limits the capacity of the instrument to measure weak precipitation formed by small raindrops. This is reflected in Figure 11, which shows a comparison of concurrent 1 min rain rates obtained with Method3 and Parsivel; the latter considerably underestimates Method3 values for rates below 0.1 mm/h, but overestimates them for some cases above 1 mm/h.

Another comparison is performed considering 30 min averages, which is the AWS time resolution, shown in Figure 12. It displays the first (lowest) and third height bin provided by Method3 (which is the first bin available for Method2) and the AWS and disdrometer rain rates. It shows a substantial agreement between the disdrometer and the AWS data and some variability when compared to Method3 rain rates.

4.1.5. Stratiform vs. Convective Rain

Following the criteria proposed by Thurai et al. [54], two parameters of the fitted gamma raindrop size distribution (D_m and N_w described in Section 4) are used by Method3 to classify bins identified as rainfall into three possible regimes: convective, stratiform or transition. Figure 13a illustrates the classification in the D_m and N_w space calculated from Method3 processing for all heights where rainfall is detected, confirming the predominantly stratiform character of the episode, with some periods of convective rain. Figure 13b shows a similar scatterplot comparing the lowest MRR height bin (100 m a.g.l.) and the values obtained with the disdrometer, where it is apparent that both instruments share a similar pattern but with differences that can be explained by the fact that the disdrometer has a detection limit on the smallest raindrops (~0.25 mm) and a much smaller sampling volume. This is particularly evident for log(N_w) values below 2 m⁻³ mm⁻¹. The agreement between these instruments is consistent with the recent results obtained by Adirosi et al. [15].

4.2. Hydrometeor Classification Verification

4.2.1. Verification Data and Methodology

The quality of Method3 hydrometeor classification is assessed by comparing the lowest height bin (from 100 to 200 m a.g.l.) with Parsivel present weather precipitation type at 1 min intervals. This high temporal resolution may easily introduce double penalty effects in the case of rapidly changing precipitation types, so a fuzzy verification approach is required considering neighborhoods either on the observations or on the forecast—see, for example, Ebert [57] or Trapero et al. [58] for two dimensional fuzzy verification procedures. Here, this is evaluated as a one-dimensional data set so the neighborhood is simply a time interval around the observation time.

To choose the interval length, two aspects are considered. First, the time required for hydrometeor particles to reach the ground due to their fall speed. Second, the precipitation drift caused by horizontal wind, which hampers matching the radar observation with ground records (Collier [59], Sandford [60]). To illustrate the first aspect, we may consider that the smallest raindrop detected by Parsivel (~0.25 mm) at 200 m a.g.l., which takes about 4.5 min to reach the ground so at least a 5 min window after the observation time should be considered—other precipitation particles as snowflakes may take even more. In the case of horizontal wind, the situation is more complex, as it may be impossible to observe aloft the same observation particle recorded on the ground, particularly for the long drifts possible for winter precipitation types (Thériault et al. [61]). Moreover, the relative position (upwind or downwind) from the ground record requires us to consider positive and negative time intervals (i.e., time windows centered on the observation time). Considering these aspects, a time window of ± 20 min was considered for evaluation.

A data set of 32 different days from January 2017 to October 2019 was selected (45,384 min), representing a wide variety of precipitation types and coverage of every season according to local climatology. The seasonal day distribution was: 5 winter days, 19 in spring, 4 in summer and 4 in autumn. Two examples of verification days are shown in Figure 14, a ground level transition from rain to snow, and a warm season convective event with some minutes of hail.

Verification scores based on a contingency table were calculated for each individual hydrometeor type, in particular the Probability of Detection (POD), the False Alarm Rate (FAR) and the Odds Ratio Skill Score (ORSS), which assess how good a forecast is compared to random chance (see Appendix A).

4.2.2. Verification Results

Table 4. Verification scores comparing Method3 precipitation type with Parsivel observations.

Class	POD	FAR	ORSS	Method3 (min)	Method3 (%)	Disdrometer (min)	Disdrometer (%)
Rain	0.99	0.29	0.99	7095	15.6	7173	15.8
Drizzle	0.69	0.26	0.72	3502	7.7	3108	6.8
Hail	0.55	0.01	0.98	49	0.1	88	0.2
Snow	0.97	0.14	0.99	3700	8.2	3897	8.6
Mixed	0.79	0.17	0.89	1001	2.2	933	2.1
No precipitation	0.94	0.05	0.99	30,037	66.2	30,185	66.5

shows for each precipitation type the value of POD, FAR and ORSS plus the total number of minutes of each type in the Method3 and disdrometer data sets, which share generally a similar proportion. Note that the number of hail minutes is rather limited; however, it is included to illustrate the relatively good results achieved. During the verification of the lowest bin, no cases of Method3 unclassified precipitation type arose, but they are marginally present in some bins aloft.

POD values indicate that both rain and snow, and also no precipitation, are the classes best detected (above 0.93), rain being the highest (0.99) and hail the lowest (0.55). F ranges from 0.29 (rain) to 0.01 (hail) and 0.05 (no precipitation). ORSS values indicate the substantial skill of Method3 hydrometeor classification, starting drizzle from 0.72 and yielding rain, snow, and no precipitation, the best results (0.99). Overall, these results illustrate that Method3 provides a reasonable classification of the precipitation types considered.

5. Discussion and Conclusions

A new methodology has been presented for processing K-band vertically pointing Doppler radar data recorded with Micro Rain Radar (MRR) systems. The methodology, referred to here as Method3, has been compared with two previously existing processing systems, Method1 (from the manufacturer) and Method2 (detailed in Maahn and Kollias [38]), using a vertical resolution of 100 m and time resolution of 1 min.

Method3 processes as input data spectral reflectivity (MRR raw data files) and produces as output data a number of fields. The first part of Method3 processing deals with spectral density processing and includes a new peak signal selection and noise treatment approach, considering all Doppler bins but the first and the last one. Then, a dealiasing method allowing for upward velocities, similar to the one included in Method2, is applied. With this new spectral processing, Method3 is able to extend the precipitation profile to the second lowest height level; in this case 100 to 200 m above ground. Moreover, the new methodology allows us to select the integration time (set here to 60 s), for example to improve sensitivity.

In Method3, the second processing part produces different variables, which include equivalent reflectivity (Z_e), Doppler fall speed and derived parameters, such as spectral width, skewness, and kurtosis, plus a simplified precipitation-type classification. The precipitation classes considered are drizzle, rain, snow, mixed, and hail. For liquid precipitation, Mie backscattering is assumed, which allows us to provide integral parameter reflectivity (Z), liquid water content (LWC), rainfall rate (RR) and a gamma drop size distribution fit, including the computation of the normalized intercept parameter with respect to the liquid water (N_w), and the mean mass-weighted raindrop diameter (D_m). Snow rate is also calculated for snow precipitation type. Compared to previously existing MRR processing methodologies, Method3 provides a comprehensive set of variables to study precipitation profiles to support precipitation microphysics analysis.

Method3 is illustrated with a case study comparing the results with Method2, Method3, microwave radiometer derived temperature profiles and ground data provided by a Parsivel disdrometer and two AWS, yielding consistent results. The comparison with Method2 denotes a high correlation for W and Z_e (R² of 0.995 and 0.993, respectively) with some exceptions for low reflectivity values, which may arise from differences in the signal detection and dealiasing. Comparisons with Method2 indicate that Method3 provides very similar patterns of fall speed and reflectivity.

Additionally, the Method3 precipitation type classification is compared with Parsivel present weather observations using contingency table scores. Results indicate a very good capacity of Method3 to distinguish rainfall and snow (PODs equal or greater than 0.97), satisfactory results for mixed and drizzle (PODs of 0.79 and 0.69) and acceptable for a reduced number of hail cases (0.55), with relatively low rate of false alarms and good skill compared to random chance in all cases (FAR < 0.30, ORSS > 0.70).

The methodology presented in this article has been implemented in Python and is freely available in the repository github as RaProM (https://github.com/AlbertGBena/RaProM). The parameters calculated in Method3 and the new hydrometeor classification proposed will facilitate the analysis of precipitation profiles for the MRR data-user community. Future work planned includes the extension of Method3 to process MRR-PRO data and the application to a larger data set to further verify the results presented here.