Low-level Polarimetric Radar Signatures in Enkf Analyses and Forecasts of the May 8, 2003 Oklahoma City Tornadic Supercell: Impact of Multimoment Microphysics and Comparisons with Observation

The impact of increasing the number of predicted moments in a multimoment bulk microphysics scheme is investigated using ensemble Kalman filter analyses and forecasts of the May 8, 2003 Oklahoma City tornadic supercell storm and the analyses are validated using dual-polarization radar observations. The triple-moment version of the microphysics scheme exhibits the best performance, relative to the single-and double-moment versions, in reproducing the low-í µí± DR hail core and high-í µí± DR arc, as well as an improved probabilistic track forecast of the mesocyclone. A comparison of the impact of the improved microphysical scheme on probabilistic forecasts of the mesocyclone track with the observed tornado track is also discussed.


Introduction
The assimilation of radar data into storm scale models using the Ensemble Kalman Filter (EnKF) [1] approach has proven to be an extremely useful tool for the analysis and prediction of convective storms in recent years.There have been many recent successful uses of this approach for both analyses [2][3][4][5][6][7] and short-range forecasts [8][9][10] based on these analyses.In general, these studies have focused on improving techniques for assimilation of radar data, on the design of the overall data assimilation system, or on the impact of initial and boundary conditions.High-resolution numerical weather prediction has progressed during the past decade such that prediction of the dynamics of individual convective storms is now routinely attempted.One substantial challenge is the improvement and validation of the microphysics parameterization and the associated impacts on storm structure and behavior (e.g., through the development of the cold pool).Errors from the model's microphysical parameterization can significantly impact forecasts of these storms.Polarimetric radar observations offer a rich source of data to validate the output of such schemes within this context.
Several storm-scale simulation studies have shown that the microphysics parameterization has a profound impact on simulated storm structure and behavior [11][12][13][14][15][16][17] and even on tornadic potential [18].Here, we restrict our discussion to bulk microphysics schemes, which assume a priori a certain functional form for the underlying drop or particle size distribution (DSD/PSD) for several hydrometeor categories.Typically, one or more moments of the PSD for a given category are explicitly predicted within a scheme, with single-moment schemes that predict the mass mixing ratio

EnKF Analysis and Forecast Experiments. All experiments use the NSSL Collaborative Model for Multiscale
Atmospheric Simulation (COMMAS) [20,40,41] and its associated EnKF radar data assimilation system [3,38].The model is run using a horizontally homogeneous background environment representative of the inflow conditions during the mature phase of the May 8, 2003 tornadic supercell storm (not shown).No surface or radiation physics are included in the simulations, and the bottom and top boundaries are free slip.We use a horizontal domain size of 100 km in each direction.The grid spacings are 1 km in the horizontal directions and utilize a stretched vertical grid with 60 total levels with the lowest 10 levels set to 150-m grid spacing, stretched thereafter to 600 m at the model top at 26 km AGL.This storm produced a long-track F4 tornado in the city of Moore, OK USA and the reader is referred to the work of Romine et al. [37] for a thorough overview of the storm.While some previous studies of this storm assimilate both radar reflectivity () and radial velocity (  ) data from nearby S-band radars [3,10], we choose to assimilate only   data from the KOUN polarimetric S-band radar.However, reflectivity data were used to define the regions for additive noise and thermal bubbles, in order to spin up the analysis and increase ensemble spread [38].These data were first quality controlled, dealiased, and objectively analyzed to a regular 2 km horizontal grid, but left on the original conical radar sweep surfaces [42].The   observations are then aggregated into 2 min bins and assimilated every two minutes from 2040 to 2208 UTC on May 8, 2003.This period covers the time from the beginning of the first echoes associated with the storm to just before tornadogenesis [38].Furthermore, we restrict the covariance updates of the model state variables by the assimilated   to only the three wind components (, V, ), potential temperature , and water vapor mixing ratio  V .The goal of this "minimalist" or "quasi-kinematic" EnKF update strategy is to allow the microphysics state variables as much freedom as possible to evolve according to the microphysics scheme itself, without any direct updating from the EnKF portion of the data assimilation system.In this regard, our study differs from the recent study of Jung et al. [6] who also verified their analyzed storm against polarimetric observations but directly assimilated reflectivity throughout the analysis period.
A total of four assimilation experiments with 30 members each are performed, each utilizing an increasing number of predicted moments using the NSSL multi-moment microphysics scheme [20,43], from one (1 M) to two (2 M and 2 MSC) to three (3 M).In 1 M, all fixed intercept parameter values for rain, snow, and hail are set as in Dawson et al. [12] for their "LINB" experiment, with the intercept parameter for rain  0 reduced from its default value of 8.0 × 10 6 m −4 to 4.0 × 10 5 m −4 , to reduce the overall strength of the cold pool.A reduction in  0 by a factor of 10 reduces the evaporation rate by roughly half.The single-moment graupel category has an intercept of 4.0 × 10 5 m −4 and particle density of 500 kg m −3 .For the double-moment configuration, two separate experiments are performed, one without correction to the sedimentation (2 M) and one using the sedimentation correction method (I+II, 2 MSC) of Mansell [44], which compensates for the excessive size sorting inherent in doublemoment sedimentation [39,[44][45][46].Finally, from the 2200 UTC analyses of each experiment, 1 hr forecasts are launched for each ensemble member, from which vorticity probabilities are computed (Section 3.2).
Excessive size sorting is characterized by increased reflectivity (unreasonably large particles) when the mass moment () falls too much faster than the number concentration (), particularly at the downward leading edge of a precipitation shaft, and is more pronounced for small fixed values of the shape parameter (e.g., [38]).Briefly, the correction scheme of Mansell [44] creates a temporary field of reflectivity () moments, and then all three moments sediment using the corresponding moment-weighted fall speeds.The sedimented mass and reflectivity moments (  and   ) are then used to generate a concentration number (  ), which is compared to the sedimented number concentration (  ).The final number concentration is the maximum of   and   at each point.The result is to artificially increase  to prevent  from increasing [43].In a 3-moment scheme, however, the size distribution shape parameter increases automatically in response to size sorting, which decreases the differences between the moment-weighted fall speeds and thus decelerates further size sorting.

Microphysics and Polarimetric Radar Emulator.
We use an updated version of the multi-moment microphysics scheme described in Mansell et al. [20] and originally based on Ziegler [43].The updated scheme predicts up to three moments (mixing ratio, total number concentration, and radar reflectivity factor-proportional to the 3rd, 0th, and 6th moments) of the assumed three-parameter gamma size distribution [47] for the hydrometeor categories of rain, graupel, and hail and up to two moments for cloud, cloud ice, and snow.The 6th moment closure mainly follows Milbrandt and Yau [23] and is described in detail in Dawson et al. [48].Additionally, as described in Mansell et al. [20] and Mansell and Ziegler [49], the bulk densities of graupel and hail are predicted for the double-and triple-moment versions of the scheme but are held fixed at 500 kg m −3 and 900 kg m −3 , respectively, for the single-moment scheme.
To derive polarimetric variables from the predicted model microphysics state variables, we utilize an updated version of the polarimetric emulator described in Jung et al. [34].The emulator makes use of the T-matrix method [50][51][52] to compute scattering amplitudes across size bins for each precipitating hydrometeor category, where the size bins have been discretized from the model predicted PSDs at each grid point.From these scattering amplitudes, the standard polarimetric variables of reflectivity at horizontal and vertical polarization (  ,   ), differential reflectivity  DR , specific differential phase  DP , and cross-correlation coefficient   are computed.
The updates primarily apply to how the water fraction on wet graupel and hail is diagnosed.Like Jung et al. [34], we diagnose a water fraction at each grid point where both ice and rain are present, where the ice can be any or all of the snow, graupel, and hail categories, by borrowing a portion of the available water from the rain field.However, whereas Jung et al. [34] assumed that the computed water fraction was applied evenly across the sizes in a given graupel or hail distribution, we make use of the empirical relationship of critical water fraction derived by Rasmussen et al. [53], which allows for a varying water fraction across the size distribution (possibly leading to completely melted hail or graupel that is smaller than 8 mm).This updated method provides more realistic polarimetric signatures in regions where partially melted hail and/or graupel are present.A complete description of the diagnostic water fraction technique is provided in Dawson et al. [48].

Low-Level Polarimetric Signatures.
From each of the four assimilation experiments, we choose the final analysis time (2208 UTC) and compute the polarimetric variables for the prior analysis (i.e., the 2 min mean forecast from the previous EnKF analysis update) using the emulator.The variables are then vertically interpolated to the 2.5 ∘ elevation sweep surface corresponding to the KOUN radar, to facilitate direct comparison with the objectively analyzed observations (Figure 1).For brevity, we focus only on the ,  DR , and   fields in this paper.The observed  DR signatures (Figure 1(b)) feature a well-defined  DR "arc" near the southeast edge of the forward flank and a low- DR (and high ) hail core between the hook echo region and the  DR arc [35,37].The  DR arc has been hypothesized to be a result of size sorting of large rain drops toward the upwind (in a storm-relative sense) side of the forward flank in the presence of low-level environmental wind shear by Kumjian and Ryzkhov [54].In   , the region of hail near the core is highlighted by relatively low values (∼0.90-0.96, Figure 2(a)), commensurate with the presence of mixed rain and tumbling hail.
Substantial differences in the simulated polarimetric fields are seen across each of the four experiments (Figures 1(c)-1(j) and 2(b)-2(e)).In general, the triple-moment experiment 3 M performs the best of all in qualitatively reproducing the  DR arc signature and performs reasonably well for the   pattern (Figures 1(j) and 2(e)), though it tends to over predict both the magnitude and coverage of high reflectivity in the core (Figure 1(i)).In contrast, the singlemoment experiment 1 M performs poorly in  DR and   (Figures 1(d) and 2(b)), as does the double-moment without sedimentation correction (2 M, Figures 1(f) and 2(c)), but on the other hand has better  in the hook region of the storm.In 1 M, the forward flank region is oriented differently than the observations (from southwest to northeast, Figure 1(c)) and the  DR arc and low- DR and   hail core are entirely absent (Figures 1(d previously, separate graupel and hail categories are predicted in the double-and triple-moment versions of the scheme.The graupel in the low levels, however, has hail-like characteristics due to a high predicted bulk density (not shown), and so is referred to as "hail" here.The separate hail category shows similar behavior (not shown), and so is omitted for the sake of brevity).In 2 M, the  DR arc is almost completely absent due to very large predicted mean hail diameters (∼20 mm) in the forward flank region where the  DR arc would otherwise be (Figure 3(d)), despite the presence of large ( mr ∼ 6 mm) rain here (Figure 3(c)).Additionally, the reflectivity values in the core are much too high (>70 dBZ, Figure 1(e)), and the   magnitudes are too much low (∼0.90-0.96)over a broad region in the forward flank (Figure 2(c)).The doublemoment experiment with sedimentation correction (2 MSC, Figures 1(g The patterns of mean volume diameter of both rain and hail in the double-and triple-moment experiments (Figures 3(c)-3(h)) reflect the action of size sorting in these experiments, with the largest hail and rain falling out closer to the updraft and south edge of the forward flank in each case [54,55], while the smaller rain and hail are advected further downstream in a storm-relative sense (i.e., to the north and northeast).The somewhat smaller  mr in the forward flank in 2 MSC relative to 3 M (compare Figures 3(e) and 3(g)) may be a result of the correction scheme increasing total number concentration  tr to prevent spurious growth in  via size sorting, whereas in 3 M, the shape parameter of the gamma distribution is also allowed to increase, allowing larger  mr for the same  [44].For rain in 3 M (Figure 3(g)), this leads to larger  DR values in the  DR arc region than elsewhere in the storm, due to the well-known increase in  DR with increasing drop oblateness and the fact that the hail sizes in this region are relatively small ( mh < 10mm) and relatively wet (Figure 4).The region of larger hail sizes in 3 M (Figure 3(h)) corresponds well with the low- DR hail core, as expected, due to the assumed tumbling characteristics for large hail in the emulator, which, combined with less diagnosed water fraction (Figure 4), tends to drive  DR towards lower values [34].
These results demonstrate that a triple-moment bulk microphysics scheme, in which hydrometeor size sorting is "fully" parameterized (i.e., all three parameters of the assumed gamma size distribution can vary independently), performs well in reproducing the salient polarimetric signatures in an analyzed supercell thunderstorm.In contrast, microphysics schemes that do not allow size sorting (singlemoment schemes) or parameterize it only partially (doublemoment schemes without some sort of correction mechanism to take into account the effects of the narrowing of the distribution during size sorting) exhibit degraded polarimetric signatures as compared with observations [39,44,45].The seemingly better  in the hook region in the 1 M case is not typical of pure simulation tests, and the 15 min forecast (below) shows a broadened reflectivity region more similar to the other cases.The 1 M ensemble members produce much less hail and have graupel that tends to be smaller and lower in density than the multimoment cases, both of which contribute to the lower reflectivity in the hook region.These results are consistent with those of Dawson et al. [48] who similarly investigated the impact of size sorting on the low-level polarimetric signatures in idealized simulations of a different supercell thunderstorm; they also found that the triple-moment scheme, by virtue of its complete parameterization of size sorting, performed better than the single-moment scheme as compared with observations.

Forecast Probabilistic Vorticity Swaths.
To further evaluate the impact of the microphysics scheme complexity, we investigate four ensemble forecast experiments, each of which is launched from the 2200 UTC prior ensemble state of the above four EnKF experiments and run out to 2300 UTC.Stensrud and Gao [56], Dawson et al. [8],and Yussouf et al.
[10] made use of probabilistic vorticity swaths to evaluate the ability of their ensemble supercell forecasts to predict the overall track of the mesocyclone.We follow the approach of these studies and compute the ensemble probability of vertical vorticity  at the surface and 1 km AGL exceeding a given threshold (0.01 s −1 ) for each of the four experiments, at any time over the 1 hr forecast period.The procedure is on a grid point-by-grid point basis for each level.For example, at a given grid point, if 3 out of the 30 ensemble members exceed the given threshold of , the ensemble probability at that point is 10%.In this manner, a probabilistic vorticity "swath" is derived.The 1 km level is used to indicate lowlevel mesocyclones, whereas the surface (lowest model level) yields an indication of how well the model can maintain the enhanced vorticity imposed by the data assimilation.The results (Figure 5) show that while each experiment shows an overall southward bias as compared with the observed tornado track (also seen in Yussouf et al. [10]), especially at the surface, nevertheless the double-and triple-moment experiments outperform the single-moment experiment, with the largest improvement seen for the surface mesocyclone (The presence of high probabilities of significant rotation at the surface is used as a proxy for the tornado track, since grid spacings of 1 km are far too coarse to explicitly resolve the tornadic circulations themselves [8,10].)Investigation of the surface cold pool at the 16 min forecast time for a select member for each experiment (Figure 6) reveals one possible reason for the differences in the mesocyclone track forecast at the surface: the cold pool strength in the double-and triple-moment schemes is somewhat larger than that in the corresponding single-moment scheme.This allows stronger rotation at the surface to persist, while it is suppressed in the single-moment scheme due possibly to lack of surface convergence.These results are overall consistent with the results of Dawson et al. [8] who also found an improved track forecast for the double-moment scheme they investigated.We hasten to point out, however, that whether the single-moment scheme will exhibit stronger or weaker cold pools than the multi-moment schemes depends strongly on the choice of tunable parameters in the scheme, particularly the intercept parameter for the various hydrometeor categories [8,14,18]; it may be difficult a priori to choose reasonable values, whereas the multi-moment schemes offer substantially more flexibility in allowing these parameters to vary in time and space.In any case, it appears that in the current study, the cause of the stronger cold pool in the multi-moment schemes relative to that of the single-moment scheme is tied to the presence of more hydrometeor mass in the low levels (not shown) which in turn results in more evaporation and melting.Further investigation of these differences is left for future work.

Conclusions
The goals of this study were to investigate the impact of varying the number of moments predicted in a multimoment bulk microphysics scheme on EnKF analyses and ensemble forecasts of a tornadic supercell thunderstorm.We focused on two main aspects: (1) the ability of the microphysics to qualitatively reproduce the low-level polarimetric (specifically ,  DR , and   ) signatures seen in the observed storm and (2) the impact on probabilistic vorticity swath forecasts (as a proxy for the tornado track).For the former aspect, the triple-moment scheme, at least in part by virtue of its complete parameterization of size sorting [48], was able to reproduce well the overall placement and magnitudes of the  DR arc and low- DR and   hail core.For the latter aspect, both the double-and triple-moment experiments performed much better than the single-moment scheme in predicting the overall track of the surface mesocyclone out to 1 hour, as compared with the observed tornado track.Though all experiments had a south-of-track bias, the single-moment scheme failed to maintain significant rotation near the surface during the duration of the forecast, while both the doubleand triple-moment schemes maintained rotation.The weaker cold pool in the single-moment experiment likely led to the suppression of significant near-ground rotation due to weaker surface convergence.We stress, however, that this reported sensitivity is not general to single-moment schemes, but rather depends substantially on the a priori choices of fixed parameters in such schemes (typically the intercept parameters of the assumed exponential size distributions).
Our overall results indicate that multi-moment schemes, in broad agreement with several recent studies, generally outperform their single-moment counterparts when various methods of analysis and forecast performance are assessed.In the current study, the simulated polarimetric signatures and ensemble vorticity forecasts were investigated.For future work, a more robust examination of the impact of microphysics might be undertaken across multiple cases.Additionally, sensitivity to various microphysical processes that become important in a multi-moment context such as the rain drop breakup parameterization [17] is needed.The results shown here and in related work demonstrate the large sensitivity of forecast output to the choice of microphysical parameters and/or parameterization that will need to be addressed prior to any operational implementation to the warn-on-forecast paradigm [57].
) and 1(h)) is intermediate between 2 M (Figures 1, 3(e), 3(f),and 2(d)), and 3 M (Figures 1, 3(g), 3(h), and 2(e)), and in fact 2 MSC appears to perform slightly better than 3 M in the magnitude and coverage of the low-  region (compare Figures 2(d) and 2(e) with Figure 2(a)), although this may come at the expense of a muted  DR arc.

Figure 2 :
Figure 2: As in Figure 1 but for cross-correlation coefficient   .

Figure 5 :
Figure 5: Ensemble probability of vertical vorticity  exceeding 0.01 s −1 (color fill) for the forecast period 2200-2300 UTC at (left column) the surface and (right column) 1 km AGL for each of the four ensemble forecast experiments: (a) and (b) 1 M, (c) and (d) 2 M, (e) and (f) 2 MSC, and (g) and (h) 3 M. Also overlaid in each panel is the observed May 8, 2003 F4 tornado track (bold line) and Oklahoma county boundaries (thin gray lines).