HTML
--> --> --> -->2.1. Advanced Regional Prediction System cloud analyses and 3DVar
We used the Advanced Regional Prediction System (ARPS) for DA and forecasts. Radial velocity data were assimilated by 3DVar using the observation operator proposed by (Gao et al., 2004) with the anelastic mass conservation equation as a weak constraint. Reflectivity data were analyzed by ARPS cloud analyses following ARPS 3DVar assimilation. Within the cloud analysis framework, precipitation types were diagnosed based on observed reflectivity and background states, and then hydrometeor state variables were retrieved and updated by applying radar reflectivity operators.Three sets of observation operators were used to convert reflectivity to hydrometeor mixing ratios. The first set of operators combined the Kessler reflectivity equation (Kessler, 1995) for warm rain and the (Rogers and Yau, 1989) reflectivity formula for snow and hail, which we refer to as KRY. The second set of operators (SMO) were defined by (Smith et al., 1975); their equations have also been applied using EnKF techniques to successfully form storm-scale analyses (Tong and Xue, 2005; Xue et al., 2006; Jung et al., 2012). However, both KRY and SMO are limited to single-moment microphysics schemes. The third set of operators (N0D), which use mixing ratio retrieval equations with diagnosed intercept parameters (Zhang et al., 2008; Wainwright et al., 2014), were developed by (Pan et al., 2016). These operators can be used in forecasting with double-moment microphysics schemes. Detailed descriptions of these procedures and equations can be found in (Pan et al., 2016).
2
2.2. DA design and forecasts
The mature stage of a trailing stratiform squall line that occurred in southern China on 23 April 2004 is shown in Fig. 1a. As shown in the cross section along A-B (Fig. 1b), the convective tower with the most intense radar echo (marked by the thick red line) was located at the leading edge of the system, and followed by a second, larger area of enhanced radar echo known as a tailing stratiform region (marked by the thick blue line). The leading convective line and the stratiform region were separated by a weak reflectivity region referred to as a transition zone. The leading edge lines of the convective region of the squall line at 2200 UTC 23 April to 0400 UTC 24 April, with a 2-h interval, are shown in Fig. 2 to demonstrate the southeastward squall-line movement across Guangdong Province, China. Six China New Generation 1998 Doppler radars (CINRAD-98D), including those at Guilin, Shaoguan, Guangzhou, Xiamen, Fuzhou, and Jianyang, captured the evolution of the squall line.Figure1. Observed squall line, in its mature stage (at 0200 UTC 24 April), which occurred on 23 April 2007: (a) composite mosaic radar reflectivity (color scale; units: dBZ); (b) cross sections of mosaic reflectivity (color scale; units: dBZ) along A-B, with the convective and stratiform region denoted by red and blue thick lines.
The domain configuration is defined in Fig. 2. The outer domain consisted of 323×323 grid points with 9-km grid spacing, and the inner domain consisted of 579×579 grid points with 3-km grid spacing. The vertical grids of both domains were 53 levels and stretched using 400-m averaged grid spacing with near-surface vertical spacing of about 50 m. The model top was at 20.4 km.
Figure2. Map of the 9-km and 3-km horizontal resolution model domains. The thick solid line indicates the fine domain. Circles indicate the maximum ranges (460 km) of the following radars: Guilin (GLRD), Shaoguan (SGRD), Guangzhou (GZRD), Xiamen (XMRD), Fuzhou (FZRD), and Jianyang (JYRD), which are indicated by blue triangles. Solid orange lines indicate frontal edges of the convective region of the squall line from 2200 UTC 23 April to 0400 UTC 24 April at 2-h intervals.
National Centers for Environmental Prediction Global Forecast System analyses were used as the initial condition for the outer domain. Forecasts were carried out for 16 h. Two control experiments were conducted by employing the Lin single-moment microphysics scheme (Lin et al., 1983) (CtlL) and Milbrandt and Yau double-moment (MY2) microphysics scheme (Milbrandt and Yau, 2005) (CtlM2) initialized from interpolated fields from a 10-h forecast at 2200 UTC in the outer domain. For the radar DA experiments, 8-h forecasts of the outer domain at 2000 UTC in CtrlM2 were interpolated into the 3-km grid to initialize the inner domain to guarantee identical initial fields. The DA window of all DA experiments was from 2000 to 2200 UTC, to cover the initial stage of squall-line evolution. To test the DA frequency sensitivity, we conducted two sets of experiments assimilating radar data with 10-, 20-, 30-, and 60-min intervals for the 2-h DA window from 2000 to 2200 UTC 23 April, using the SMO and N0D schemes (Fig. 3). These experiments were named S10L, S20L, S30L, S60L, D10M2, D20M2, D30M2, and D60M2 (Table 1). D10M2, D20M2, D30M2, and D60M2 employed the N0D DA scheme and the MY2 scheme for forecasting; whereas, S10L, S20L, S30L and S60L employed SMO for DA and the Lin scheme for forecasting. Besides, in high-frequency DA cycles, updating the water vapor mixing ratio and cloud water can lead to unrealistic warming in the midlevel troposphere, especially for mesoscale convective systems (Schenkman et al., 2011a). Thus, following (Pan et al., 2016), the cloud water and water vapor mixing ratio were not updated during DA cycles. The temperature inside cloud was adjusted by a modified moist air-parcel ascent that accounted for environmental air entrainment (Hu et al., 2006a). All radial velocity and reflectivity observations from the six radars were assimilated using ARPS 3DVar and cloud analyses. Because the observed reflectivity is usually more reliable than that forecast by the model, the precipitable hydrometeor variables from background fields were replaced by those retrieved directly from cloud analyses. The precipitable variables included rainwater, snow, and hail mixing ratio for SMO; and rainwater, snow, graupel and hail mixing ratio, and total number concentration for N0D. After 2 h of DA, 6-h forecasts were produced, extending to 0400 UTC 24 April.
Figure3. Timelines of experimental forecasts and analyses.
The main model physics parameterizations were as follows: fourth-order computational mixing, fourth-order advection in the horizontal and vertical, a 1.5-order turbulence kinetic energy-based subgrid-scale turbulent mixing scheme with planetary boundary layer parameterization, a rigid top boundary combined with a wave-absorbing layer, a two-layer land surface model, and Goddard Space Flight Center short- and longwave radiation parameterization.
Under identical DA and forecasting con?gurations, comparing S10L, S20L, S30L, and S60L or D10M2, D20M2, D30M2, and D60M2 determines the impact of DA frequency on forecasting results. Conversely, under identical DA frequencies, comparing S10L and D10M2, S20L and D20M2, S30L and D30M2, or S60L and D60M2 determines forecasting differences due to the use of single- and double-moment schemes. We investigated initial noise and model adjustment based on experiments employing 10-, 20-, 30-, and 60-min DA intervals. Then, the equitable threat scores (ETSs) and bias scores (BSs) were used to evaluate precipitation forecasting accuracy in experiments using different DA intervals. Observed 1-h precipitation data were obtained from rain gauge measurements.
-->
3.1. Adjustment following DA
3.1.1. Initial noiseBetter meso- or convective-scale forecasts are typically expected after assimilating radar observation data. However, the analyzed state variables are usually inconsistent with the model dynamic frame (Hu and Xue, 2007), particularly for a variational DA framework and cloud analyses, because these DA methods lack flow-dependent covariance and cross correlations among prognostic state variables. During forecasting, state variables are adjusted to balance each other, depending on the dynamics of the model. In every DA cycle, the model balance is disturbed, and then slowly rebuilt during forecasting. Experiments with different DA frequencies and a fixed DA window could result in different degrees of model noise. Therefore, the interaction between DA and forecasting is quite complex. Hence, the mean absolute tendency of surface pressure was employed to investigate balance adjustment following DA.
The mean absolute tendency of surface pressure N (Lynch and Huang, 1992) is defined as \[ N=(\frac{1}{IJ})\sum^{I}_{i=1}\sum^{J}_{j=1}|\frac{\partial p_{s}}{\partial t}|_{i,j} , \] where p s is the surface pressure and i and j are the horizontal x and y grid points. The absolute values of the surface pressure tendency were averaged over the whole domain, where IJ represents the total grid points of the whole domain. N reflects the overall balance of the model states (Lynch and Huang, 1992; Chen and Huang, 2006), and its time variants indicate the length of time during which spurious high-frequency noise is dumped.
N values for all experiments during first-hour forecasts from 2200 UTC to 2300 UTC were calculated every minute; the results are shown in Fig. 4. Without radar DA, N values from CtlM2 and CtlL were lower than those from experiments with assimilated radar data. For DA experiments, N values were all higher than for control experiments, and sharply increased during the first minute, before then decreasing to about 5 Pa min-1 after 1 h of forecasting. The dramatic oscillation observed during the first 20 min, before 2220 UTC, implied strong adjustment of the model. The pressure tendencies in experiments using a 10-min DA interval had the highest N values among all experiments, which indicates that a higher DA frequency could provoke greater imbalance into the final analysis, also known as the initial condition. At the same DA frequency, the time variants of N were very similar between experiments S10L and D10M2, S20L and D20M2, and between S30L and D30M2, which suggests that the selection of a single- or double-moment scheme did not significantly affect the forecasts. However, when a longer DA interval was used, i.e., 60 min, the N values in experiment S60L were greater than those in D60M2. To some extent, this indicates that using a different microphysics scheme could significantly affect dynamic and microphysical variables during 60 min of forecasting, and thus produce significant differences in the final analyses between S60L and D60M2.
Figure4. The 1-h evolution of the mean absolute surface pressure tendency in forecasts starting at 2200 UTC 23 April 2004.
3.1.2. Model adjustment
Following (Hu and Xue, 2007), we used the maximum vertical velocity (Wmax) as an indicator of dynamic variable adjustment. Updraft mainly occurred within the convective region of the squall line. To guarantee the representativeness of Wmax in the evolution of the squall line, we selected Wmax from nine slides cutting through the main updraft in the convective region from 2200 to 2300 UTC (Fig. 5). As shown in Fig. 6, Wmax from CtlM2 increased after initialization, and became stable at around 12-14 m s-1 after 2210 UTC; whereas, Wmax from CtlL increased slower and reached 8-10 m s-1 after 2240 UTC. For all radar DA experiments, dramatic high-frequency oscillations in Wmax were noticeable during the first 20 min, similar to those observed in N values. Wmax then became relatively stable at around 12-14 m s-1 in experiments D60M2 and S60L, and 8-12 m s-1 in the other experiments. Since all the DA experiments were initialized from CtlM2 at 2000 UTC, the Wmax was close to CtlM2, rather than CtlL, after oscillations. Also, using lower DA frequency, the Wmax from S60L and D60M2 was closer to CtrlM2 from 2210 to 2230 UTC.
Figure5. Forecast composite reflectivity (color scale; units: dBZ), vertical velocity (contours; units: m s-1), and wind vectors (units: m s-1) at 1 km above mean sea level for experiment D30M2 at (a) 2200 and (b) 2300 UTC 24 April 2007. Blue contours are 2, 6, and 10 m s-1 and black lines represent slides used for statistical analyses (Fig. 6).
Figure6. Maximum vertical velocity (m s-1) per min of 1-h forecasting started at 2200 UTC 23 April 2004.
The interactions among the microphysics, kinematic and thermodynamic variables were quite complex. Assuming that the observed reflectivity data were more reliable than those obtained from forecasts, the hydrometeor variables, including mixing ratio and total concentration, were replaced by those retrieved from reflectivity observations. Using the same set of radar equations, we were able to obtain the same analyzed hydrometeor variables among experiments with different DA intervals. For example, the analyzed hydrometeor fields at 2200 UTC from experiment S10L were the same as those from S20L, S30L and S60L. Rainfall, as a final phenomenon of sophisticated microphysical processes, could represent the adjustment of microphysics variables due to kinematic and thermodynamic differences following DA. The averaged accumulated rainfall rate (AAR) was calculated every minute from 2200 to 2300 UTC in regions where rainfall rates exceeded 0.005 mm s-1 (Fig. 7). Though the AAR from CtlL and CtlM2 increased from 2200 to 2215 UTC, it was around 0.011-0.013 mm s-1 in 1-h forecasts from 2200 to 2300 UTC. The variations were quite a lot smaller than those that occurred in the radar DA experiments. The AAR from the radar DA experiments dropped sharply within the first 1-2 min of forecasting, and then gradually increased until 2220 UTC, after which it remained stable at around 0.01-0.012 mm s-1. This indicated that the precipitable hydrometers retrieved by the cloud analysis fell to the ground as precipitation in 1-2 min. Thereafter, rainfall was produced by the evolution of the microphysics scheme, and reached its balance in 20 min. The lowest average AAR value was obtained from experiments S10L and D10M2 using a 10-min DA interval.
We also examined the evolutions of Wmax and AAR from 1-h forecasts after the first cycle at 2000 UTC of DA for all experiments. A similar spin-up time of around 20-30 min was found. Overall, both dynamic and microphysical fields took about 20 min to spin up and recover their balance under different microphysics schemes, and did not significantly change as the cycles advanced.
Figure7. Domain-averaged accumulated rain (units: mm s-1) exceeding 0.005 mm s-1 per min of forecasting from 2200 to 2300 UTC 23 April 2004.
Radar DA can create initial fields resolving more mesoscale structures, even convective-scale information, and lead to more accurate forecasts. However, a higher frequency of radar DA produced larger initial noise and had negative impacts on forecasts, while a lower frequency of radar DA may generate only limited improvements for forecasts. In our study, the investigations on initial noise and model adjustment suggested 20 min could be an optimal DA interval. Thus, we further compared the impacts of DA interval on forecasts, as reported in the next section.
2
3.2. Impacts of DA interval on forecasts
In the mature stage of the trailing stratiform squall line, along its cross section, a front inflow originated in the boundary layer of the cold pool, extended up through the convective region, and sloped into the stratiform region. Beneath the front-to-rear inflow, a rear inflow was induced by booklet vortices (Trapp and Weisman, 2003; Atkins et al., 2004; Meng et al., 2012) or counter-rotating circulation patterns (Wakimoto et al., 2015) and microphysical processes (Yang and Houze, 1995; Grim et al., 2009). The strong rear inflow usually appeared at the bow-shaped line segment and was associated with damaging winds. To investigate the impacts of DA frequency on forecasts, the structures of the squall line were compared. Studies (Hu et al., 2006a, Hu et al., 2006b; Schenkman et al., 2011a, Schenkman et al., 2011b) have already shown that radar DA can significantly improve forecasts. Thus, we did not further compare the results from CtlL and CtlM2 to those from the radar DA experiments.Figure 8 shows the composite reflectivity of the 4-h forecasts and wind fields from all experiments. With a 10-min DA interval, the apexes of the bow echo near x=900.0 km and y=800.0 km were more pointed than those using 20-, 30- and 60-min intervals, especially for forecasts using 60-min intervals. Horizontal velocities exceeding 18 m s-1 at 500 m were mainly located at the apex of the bow and the southwestern tail of the squall line. This area was much wider in experiments with higher DA frequencies, particularly in those experiments using the SMO and Lin schemes, which might lead to wider surface wind damage.
Figure8. Forecast composite reflectivity (color scale; units: dBZ) and wind vector perturbations at 500 m above mean sea level for the radar DA experiments at 0200 UTC 24 April 2007. Blue contours indicate regions where horizontal wind speed exceeded 18 m s-1. The blue arrows indicate the apexes of the bow echo. The relative wind vectors were defined as the perturbation relative to the average of mean wind vector of the box in front of squall line. The cross sections along A-B were shown in Fig. 9.
Cross sections along A-B at 0200 UTC 24 April are shown in Fig. 9. To improve representativeness, all physical variables were averaged over an 18-km band along A-B. The reflectivity at the convective region in experiment D10M2 was more intense than that in D20M2, D30M2 and D60M2, and the same as that in S10L, S20L, S30L and S60L. The front-to-rear inflow started to rise and leaned upshear into the stratiform region. Updrafts mainly occurred above the cold pool in the convective region and above the melting level (0°C) in the stratiform region. The rear inflow entered from the rear of the stratiform region, intensified through and under the melting layer, and reached a maximum at the low level toward the back of the convective region. The cold pool was deeper in experiment D10M2 and S10L than in experiments using 20-, 30- and 60-min DA intervals under the stratiform region. For experiments S10L, S20L, S30L, S60L and D10M2, the cold pool at z=1-2 km in the convective region was wider, and associated with areas where horizontal velocity exceeded 18 m s-1. In experiments using 20-, 30- and 60-min DA intervals, the stratiform regions were wider and clearer than those in experiments using 10-min DA intervals. Also, clear subsidence occurred at the middle to high level behind the main updraft of the convective region for experiments using 30- and 60-min DA intervals, which was associated with rearward gravity-wave propagation (Fovell et al., 1992) and could be attributed to the formation of the transition zone (Biggerstaff and Houze, 1993).
Figure9. Cross sections of reflectivity (color scale; units: dBZ) and the cold pool (blue lines) as defined by a -4 K potential temperature perturbation relative to the mean potential temperature of the box in front of the squall line shown in Fig. 8 and wind vector perturbations along line A-B in Fig. 8. Gray lines indicate vertical velocity (dashed line: downdraft, -0.1, -0.4, -1.0 and -3.0 m s-1; solid line: updraft, 0.2, 1.0, 2.0 and 4.0 m s-1). White bold lines delineate the areas where the horizontal velocity exceeded 18 m s-1 under 3 km above ground level. FRJ indicates front-to-rear jet; RIJ indicates rear inflow jet.
In the previous section, a spin-up time of 20 min was found. Using 20-min DA intervals, D20M2 reproduced squall line structures close to observation (Fig. 1). However, its pattern was similar to D30M2, and it was also hard to distinguish which DA interval was best. The same situation occurred for experiments using the SMO and Lin schemes. Thus, the objective, as reported in the next section, verifications were conducted to better demonstrate the impacts of DA frequency on forecasts.
2
3.3. Objective verification of forecasts
Figures 10 and 11 show the objective ETSs and BSs verified against 1-h accumulated precipitation at thresholds of 0.5 and 10 mm h-1 from 2200 UTC 23 April to 0400 UTC 24 April using different DA frequencies.Figure10. ETSs and BSs of predicted hourly accumulated rainfall at thresholds of (a, b) 0.5 mm h-1 and (c, d) 10 mm h-1 for experiments using the N0D and MY2 schemes (see Table 1) from 2200 UTC 23 April to 0400 UTC 24 April.
Because stratiform and convective rainfall coexisted in the squall line, scores for a threshold of 10 mm h-1 were used to assess convective rainfall predictions (Chang et al., 2009, 2015), such that scores of 0.5 mm h-1 generally represented accurate overall rainfall prediction. Using the N0D and MY2 schemes, the ETS score of D10M2 (0.5 mm h-1) was higher during forecasts of 0-2 h than those of D30M2 and D60M2 (Fig. 10a). However, scores decreased quickly within the first hour of the forecast and then increased. The BSs from all experiments were greater than 1.0 at a threshold of 0.5 mm h-1, which suggests an overestimation of rainfall rate; D10M2 experienced the least overestimation during forecasts of 0-2 h. These results indicate that high-frequency DA will lead to a better fit with observations. As the rainfall threshold increased, the performance of D10M2 decreased, relative to the other experiments. At 10 mm h-1, D20M2 performed best in terms of both ETSs and BSs before 0000 UTC, and then slightly worse than D30M2 after 0000 UTC. D10M2 exhibited the worst performance (Figs. 10c and d). Unlike the experiments employing the N0D and MY2 schemes, those using the SMO and Lin schemes performed similarly. However with different thresholds, S20L produced the highest ETSs and best BSs (closest to 1.0), except for the first 2 h at a threshold of 0.5 mm h-1. A detailed discussion of the impacts of using N0D and SMO schemes can be found in (Pan et al., 2016). Overall, the experiments using the N0D and MY2 schemes performed better than those using the SMO and LIN schemes, except for the first 1 h of forecasting at a specific threshold. (Pan et al., 2016) also indicated a significantly improved performance using double-moment scheme settings (N0D and MY2) at a higher precipitation threshold (10 mm h-1), caused by including the graupel category in cloud analyses.
Figure11. As in Fig. 10 but for experiments employing the SMO and Lin schemes (see Table 1).
These objective verifications indicate that DA frequency has substantial impacts on precipitation forecasting accuracy. Of the experiments that employed single-moment scheme settings (SMO and LIN), S20L provided the best convective rain forecasts, and was comparable with S10L in light rain forecasts. Similarly, D20M2 was comparable with D30M2, both having double-moment scheme settings, and performed better than the other experiments. Generally, when the DA interval was set close to the model adjustment time, forecasts were more accurate.