J. Meteor. Res.  2017, Vol. 31 Issue (4): 731-746 PDF
http://dx.doi.org/10.1007/s13351-017-6133-3
The Chinese Meteorological Society
0

#### Article Information

Rong ZHANG, Yijun ZHANG, Liangtao XU, Dong ZHENG, Wen YAO . 2017.
Assimilation of Total Lightning Data Using the Three-Dimensional Variational Method at Convection-Allowing Resolution. 2017.
J. Meteor. Res., 31(4): 731-746
http://dx.doi.org/10.1007/s13351-017-6133-3

### Article History

Received August 8, 2016
in final form February 21, 2017
Assimilation of Total Lightning Data Using the Three-Dimensional Variational Method at Convection-Allowing Resolution
Rong ZHANG1,2,3, Yijun ZHANG1,3, Liangtao XU1,3, Dong ZHENG1,3, Wen YAO1,3
1. State Key Laboratory of Severe Weather, Chinese Academy of Meteorological Sciences, Beijing 100081;
2. College of Earth Sciences, University of Chinese Academy of Sciences, Beijing 100049;
3. Laboratory of Lightning Physics and Protection Engineering, Chinese Academy of Meteorological Sciences, Beijing 100081
ABSTRACT: A large number of observational analyses have shown that lightning data can be used to indicate areas of deep convection. It is important to assimilate observed lightning data into numerical models, so that more small-scale information can be incorporated to improve the quality of the initial condition and the subsequent forecasts. In this study, the empirical relationship between flash rate, water vapor mixing ratio, and graupel mixing ratio was used to adjust the model relative humidity, which was then assimilated by using the three-dimensional variational data assimilation system of the Weather Research and Forecasting model in cycling mode at 10-min intervals. To find the appropriate assimilation time-window length that yielded significant improvement in both the initial conditions and subsequent forecasts, four experiments with different assimilation time-window lengths were conducted for a squall line case that occurred on 10 July 2007 in North China. It was found that 60 min was the appropriate assimilation time-window length for this case, and longer assimilation window length was unnecessary since no further improvement was present. Forecasts of 1-h accumulated precipitation during the assimilation period and the subsequent 3-h accumulated precipitation were significantly improved compared with the control experiment without lightning data assimilation. The simulated reflectivity was optimal after 30 min of the forecast, it remained optimal during the following 42 min, and the positive effect from lightning data assimilation began to diminish after 72 min of the forecast. Overall, the improvement from lightning data assimilation can be maintained for about 3 h.
Key words: lightning     data assimilation     three-dimensional variational (3DVAR) method     Wether Research and Forecasting (WRF) model
1 Introduction

Severe convective weather is one of the most dangerous meteorological phenomena, and heavily affects routine work and daily life. It is small-scaled and short-lived, often accompanied by thunderstorms and strong winds. Owing to the small spatial and temporal scales of severe convective weather, as well as the lack of small-scale information in the large-scale reanalysis data driving numerical weather prediction (NWP) models, the forecast accuracy for such weather is relatively lower with current NWP models. To improve the forecasts of severe convective weather, it is necessary to assimilate high-resolution observational data that can initially resolve convective storms into NWP models, so that more storm-scale information can be incorporated into the model at the beginning of the forecast period.

One of the most evident characteristics of severe convective weather is that it is usually accompanied by frequent lightning activities. Previous observational analyses have shown that lightning data can be used to indicate areas of deep convection (Petersen and Rutledge, 1998; Schultz et al., 2011; Rudlosky and Fuelberg, 2013). With the development of state-of-the-art lightning detection technology and expansion of lightning locating systems, an increasing amount of near real-time lightning data can now be readily obtained. Lightning data have the advantage of a wide detection range, high precision, and being less influenced by terrain. Therefore, lightning data have the potential to be applied in the monitoring, early warning, and forecasting of severe convective weather (Qie et al., 2014a), especially in regions where radar data are generally not available, such as oceans and terrain-blocked areas. Assimilating lightning data into numerical models is of great importance because it allows more small-scale information to be incorporated, thus improving the quality of the initial conditions and the subsequent forecasts.

Compared with other observational data such as those from radar and satellites, lightning data and their assimilation into NWP models have received much less attention. Only a limited number of studies have attempted to assimilate lightning data into forecast models. Alexander et al. (1999) were the early investigators, who established a relationship between lightning data and rain rates. They used the rain rates retrieved from lightning data to modify the latent heating rates, and consequently greatly improved the quality of a 12–24-h precipitation forecast of a super storm. This lightning data assimilation method was later verified by Chang et al. (2001) using an alternative rainfall–flash relationship. Pessi and Businger (2009) also used a similar method to assimilate lightning data of a winter storm over the North Pacific Ocean, and showed improvements in forecasts of pressure and wind. Papadopoulos et al. (2005, 2009) nudged the model humidity profile towards the empirical humidity profile derived from atmospheric soundings during thunderstorm days. Their results showed that the prediction accuracy of convective precipitation during both the assimilation period and the subsequent 12-h forecast was significantly improved. Mansell et al. (2007) developed a method using lightning data to control the triggering of the convective parameterization scheme (CPS). The CPS was forced to produce convection where lightning indicated storms, and to prevent the production of spurious convection where no lightning was observed. In a case study, the method was successful in generating cold pools that were present in the surface observations at initialization of the forecast. Lagouvardos et al. (2013) also used this method to assimilate lightning data, and reported a significant improvement in the precipitation forecast. The advantage of this approach is that it not only produces convection, but also suppresses spurious convection. However, it relies on the use of the CPS, and hence may not be suitable for simulations at convection-allowing or convection-resolving resolutions.

With the continuous increase in model resolution, a number of attempts have been made to assimilate lightning data at convection-allowing resolution. Fierro et al. (2012) established an empirical relationship between the water vapor mixing ratio, total flash rate, and graupel mixing ratio. By using this relationship, the model water vapor mixing ratio was nudged based on the observed total lightning data. This work represented the first attempt to assimilate lightning data at convection-allowing resolution, and the technique significantly improved the representation of convection at the analysis and 1-h forecast timescales for a tornado event. The performance of this technique was further verified in some other cases (Fierro et al., 2014, 2015; Lynn et al., 2015). Marchand and Fuelberg (2014) presented a new nudging method that warmed the most unstable low levels of the atmosphere at locations where lightning was observed but deep convection was not simulated, and this method was compared with the nudging method developed by Fierro et al. (2012). The simulation results showed that this method generally performs better in simulating isolated thunderstorms and other weakly forced deep convection, while the method of Fierro et al. (2012) performs better for cases of strong synoptic forcing.

A number of studies of lightning data assimilation have been carried out in recent years in China. Li et al. (2008) assimilated the convective rainfall data retrieved from the Tropical Rainfall Measuring Mission (TRMM) satellite into the Advanced Regional Prediction System (ARPS), and the results showed that forecasts of both the storm center and rainfall density were improved. Ran and Zhou (2011) conducted a nudging method that combined the methods of Papadopoulos et al. (2005) and Mansell et al. (2007). According to the observed lightning data from the TRMM satellite, they adjusted not only the content and spatial distribution of water vapor and cloud hydrometeors, but also the model CPS. The results showed that the short-term forecasts of precipitation and hydrometeors were effectively enhanced. Qie et al. (2014b) constructed an empirical relationship between the total lightning flash rate and the ice-phase particle (graupel, ice, and snow) mixing ratio, which was used to adjust the mixing ratio of ice-phase particles within the mixed-phase region. They found that this method could improve short-term precipitation forecasts of mesoscale convective systems with high, and even moderate, lightning flash rates. Wang Y. et al. (2014) converted cloud-to-ground (CG) lightning data into proxy radar reflectivity using an assumed relationship between flash density and reflectivity in the Gridded Statistical Interpolation (GSI) system, and the proxy reflectivity was assimilated by using a physical initialization method. The results showed that predictions of reflectivity and precipitation were improved.

The above review of literature demonstrates that previous lightning data assimilation efforts have mostly employed the nudging method. However, a rather limited amount of research has been carried out using the more sophisticated variational or ensemble Kalman filter (EnKF) methods. Mansell (2014) performed an observing system simulation experiment, in which simulated total lightning data were assimilated by using EnKF, and the results showed that the method modulated the strength of convection and suppressed spurious convection. This method was further verified by real-data cases (Allen et al., 2016). Wang Y. D. et al. (2014) assimilated the CG-converted rainfall rates within the Weather Research and Forecasting (WRF) four-dimensional variational framework, and the results showed that both the initial field and the 6-h forecast were improved. Very recently, Fierro et al. (2016) set the water vapor mixing ratio to its saturation value in the layer between the lifted condensation level and an assumed fixed height of 15 km at observed lightning locations, and then the pseudo-observations for water vapor were assimilated by the three-dimensional variational (3DVAR) system of ARPS. The results demonstrated the potential value of assimilating total lightning data with the 3DVAR method.

In the present study, total lightning data were assimilated by using the 3DVAR method in the WRF data assimilation system (Barker et al., 2003, 2004, 2012; hereafter referred to as WRFDA-3DVAR). Basically, the method of the present study is similar to that of Fierro et al. (2016), both of which assimilate lightning data in the 3DVAR framework. However, instead of assimilating the water vapor mixing ratio, we assimilate the relative humidity retrieved from the total lighting data as a sounding observation using WRFDA-3DVAR in cycling mode at 10-min intervals. Since relative humidity is a variable that can be observed directly, and can be directly assimilated into any operational assimilation system, assimilating lightning data retrieved relative humidity is generally easier for implementation in operational runs. To find the appropriate time-window length for cycled assimilation having significant improvements in both the initial condition and subsequent forecast, four experiments with different assimilation time-window lengths were conducted. After finding the appropriate time length of the assimilation window, the improvement in the subsequent forecast precipitation and reflectivity was further analyzed.

It should be noted that most current operational ground-based lightning location systems can only detect CG flashes, which account for only a small fraction (about one third) of total flashes. In fact, intracloud flashes, which account for the majority of flashes, have been shown to be better correlated with updraft in clouds (MacGorman et al., 1989; Schultz et al., 2011). Hence, assimilating all flashes outperforms assimilating only CG flashes. The lightning data used in this study were from SAFIR3000 (see Section 3.1), which detects both CG and intracloud flashes.

The remainder of this paper is organized as follows. In Section 2, the thunderstorm case is introduced. In Section 3, the data sources, lightning data assimilation procedure, and experimental design are described. In Section 4, a detailed analysis of the results is presented, followed by a summary and discussion in Section 5.

2 Thunderstorm case

The case in the present study was a squall line that occurred on 10 July 2007 in North China. The weather process spanned a period of more than 6 h, during which it was accompanied by frequent lightning activity. The squall line initially developed in northwestern Beijing and propagated southeastwards across Beijing, Tianjin, and Hebei Province, before finally dissipating in Bohai Bay. During 1300–1500 UTC, hail events were recorded in some areas of Beijing and Tianjin. Figure 1 shows the 500-hPa geopotential height, temperature, and wind vectors at 1200 UTC 10 July 2007. It can be seen that North China and Northeast China were controlled by a trough. The low-level convergence and high-level divergence (Xu et al., 2016) resulted in the upward movement of air. In addition, there was abundant water vapor in the lower atmosphere (Xu et al., 2016), which was favorable for convection.

 Figure 1 Horizontal distribution of geopotential height (blue line; gpm), temperature (red line; K), and wind vectors at 500 hPa from NCEP FNL (Final) analysis data at 1200 UTC 10 July 2007.

The sounding data for Beijing at 1200 UTC (Fig. 2) were analyzed to further investigate the thermodynamic and dynamic conditions. The results indicated highly unstable atmospheric conditions. In particular, the atmosphere was wet at low levels and dry at high levels. The convective available potential energy exceeded 2000 J kg–1, and the thickness of the unstable layers was more than 8 km. Moreover, there was obvious shear in wind direction and speed at low levels. From the above analy-sis, it is clear that the thermodynamic conditions were favorable for the occurrence of the squall line.

 Figure 2 Sounding observation of Beijing station at 1200 UTC 10 July 2007.
3 Data and methods 3.1 Data sources

The lightning data used in this study were from the SAFIR3000 lightning detection system operated by the Beijing Meteorological Bureau, which consists of three substations (shown as triangles in Fig. 3) and detects total lightning. SAFIR3000 simultaneously detects lightning in two different frequency bands [very high frequency (110–118 MHz) and low frequency (300 Hz to 3 MHz)], which can discriminate intracloud and CG lightning. Its detection efficiency can reach 90% in a 200-km radius with a less than 2-km location error (Zheng et al., 2009; Liu et al., 2013). Owing to the operational running requirements, the triggering threshold of SAFIR3000 was set relatively high; hence the number of radiation sources detected per flash was only a few (typically 2–3 sources). Because of the sparseness of the radiation sources, we did not attempt to define the flashes in this study—each radiation source was regarded as a flash.

The precipitation observations used in this study were from the hourly and three-hourly ground-based measurements of MICAPS (Meteorological Information Comprehensive Analysis and Process System). The radar echo data were from the mosaic of S-band Doppler weather radar in Beijing, Tianjin, and Qinhuangdao (shown as red dots in Fig. 3), which provided a reasonable coverage of the whole lifecycle of the storm.

3.2 Assimilation method

Data assimilation is the technique by which observations are combined with an NWP product (called the first guess or background forecast) and their respective error statistics, to provide an improved estimate (i.e., the analysis) of the atmospheric state. The 3DVAR assimilation technique achieves this through the iterative minimization of a prescribed cost function,

 $\begin{split}J({x}) =& \frac{1}{2}{({x} - {{x}_{\rm b}})^{\rm{T}}}{{B}^{ - 1}}({x} - {{x}_{\rm b}}) \\ & + \frac{1}{2}{(H({x}) - {{y}^{\rm o}})^{\rm{T}}}{{R}^{ - 1}}(H({x}) - {{y}^{\rm o}}),\end{split}$ (1)

where x is the analysis to be found that minimizes the cost function J (x), xb is the first guess of the NWP mo-del, yo is the assimilated observation, B is the background error covariance, R is the observation error covariance, H is the observation operator that performs the necessary interpolation and transformation from model variables to observation space, and the superscript "T" represents the transposition of the matrix.

In this study, statistics of the differences between daily 24-and 12-h forecasts valid at 0000 UTC and 1200 UTC for the whole month of July 2007 were paired (i.e., a total of 62 samples) and used to estimate background error covariance with the widely-used NMC [National Meteorological Center, now known as the National Centers for Environmental Prediction (NCEP)] method (Parrish and Derber, 1992) for our own forecasting domain (i.e., inner domain) with the utility packages in WRFDA. The observation error used in this study was from the default error file "obserr.txt" provided in WRFDA, in which observation errors for different meteorological parameters (e.g., pressure, temperature, height, wind, and moisture) and for various types of data sources (e.g., ship reports, radiosonde/rawinsonde, etc.) at standard pressure levels are given. These values are originally from the NCEP, but have been modified at the NCAR after comparisons against observation-minus-background statistics (Barker et al., 2003).

In fact, lightning data are neither a model variable nor a diagnosed variable in most current NWP models; hence, it is necessary to establish a relationship between observed lightning data and model variables. In Fierro et al. (2012), a smooth continuous equation [Eq. (2)], which is a function of gridded flash rate and simulated graupel mixing ratio, was constructed to adjust the model water vapor mixing ratio:

 ${Q_\rm{v}} = A{Q_\rm{sat}} + B{Q_\rm{sat}}\tanh (CX)[1 - \tanh (DQ_{\rm{g}}^\alpha )].$ (2)

Here, X is the gridded flash rate at 10-min intervals; Qsat is the water vapor saturation mixing ratio (g kg–1); Qv is the water vapor mixing ratio (g kg–1); Qg is the graupel mixing ratio (g kg–1); and A, B, C, D, and $\alpha$ are constant coefficients. In Fierro et al. (2012), the number of total lightning data per grid element from the Earth Networks Total Lightning Network was gridded onto 9-km grids so as to mimic the expected (8–12 km) resolution of the Geostationary Lightning Mapper instrument onboard the Geostationary Operational Environmental Satellite R series (Goodman et al., 2012, 2013) over the Americas. Through this function, the model water vapor mixing ratio was adjusted by using nudging, based on the observed total lightning data.

Instead of directly assimilating the water vapor mixing ratio, which cannot currently be assimilated in the variational framework of WRFDA (at least in version 3.5.1), the water vapor mixing ratio in Eq. (2) was converted to relative humidity with the following relationship (Wang et al., 2013):

 ${\rm {R\!H}} = {Q_\rm{v}}/{Q_\rm{{sat}}}.$ (3)

Dividing both sides of Eq. (2) by Qsat, we obtain

 ${\rm {R\!H}} = A + B\tanh (CX)[1 - \tanh (DQ_{\rm{g}}^\alpha )],$ (4)

where RH is relative humidity, which is a control variable and could be assimilated in the variational framework of WRFDA; and X is the gridded flash rate at 10-min intervals (in our experiments, total lightning data were gridded onto 3-km grid cells). Coefficient A sets the relative humidity threshold when Eq. (4) can be used, and the meaning of the other coefficients can be found in Fierro et al. (2012). To prevent too much water being added into the model, relative humidity was adjusted according to Eq. (4) only if the model relative humidity was below 85% and the graupel mixing ratio was below 3 g kg–1, and the upper bound of relative humidity was limited to 100%. Sensitivity tests with different values (A = 81%, 90%, 95%) showed that 85% was the optimum value of parameter A for this case, and the result was not sensitive to changes in other parameters. Therefore, in the present study, A = 85%, and the other coefficients were the same as in Fierro et al. (2012), i.e., B = 0.2, C = 0.02, D = 0.25, and α = 2.2. Next, we describe the lightning data assimilation procedure.

Firstly, the lightning data detected by SAFIR3000 were parsed into 10-min intervals with the center at the analysis time (i.e., ± 5 min from the analysis time), and gridded onto the 3-km grid cells of the WRF inner domain. Whenever a flash was located within a model grid cell, and the relative humidity and graupel mixing ratio in the mixed-phase region [defined as the layer between the 0℃ and –20℃ isotherms, as in Fierro et al. (2012) and Qie et al. (2014b)] were below the thresholds of 85% and 3 g kg–1, respectively, then the relative humidity in the mixed-phase region of that grid column was adjusted according to Eq. (4). The adjusted relative humidity was then assimilated into WRF-3DVAR as a sounding observation in cycling mode at 10-min intervals. A time interval of 10 min was selected for the cycled lightning data assimilation so as to capture the development of convection.

The basic rationale of lightning data assimilation is that lightning data can be used to indicate areas of deep convection within a model domain. Based on Eq. (4), whenever a flash occurs in a given grid cell, the relative humidity in the mixed-phase region of that grid column is increased, ensuring that the relative humidity in the mixed-phase region is no less than 85%. Then, the release of latent heat increases as the water vapor condenses, thus leading to the acceleration of atmospheric updraft and, ultimately, the development of convection in locations where lightning is observed. The height range of assimilation is confined between the temperature levels of 0 and –20℃, which represents the mixed-phase graupel-rich region, and hence is most likely associated with electrification and lightning activity according to the non-inductive charging mechanism (Fierro et al., 2012; Qie et al., 2014b).

3.3 Model and simulation setup

The numerical model used in this study was the three-dimensional compressible nonhydrostatic WRF model with the advanced research dynamic solver (WRF-ARW, version 3.5.1) and WRF data assimilation system (WRFDA, version 3.5.1). The experiments were performed on two nested domains (D01 and D02; see Fig. 3), with horizontal resolutions of 9 and 3 km, respectively. The dimensions of the outer domain (D01) were 160 × 160 and those of the inner domain (D02) were 181 × 181. Both domains had 43 uneven vertical levels [the same as in Zhang et al. (2011)], with the model top at 50 hPa. The time steps for the outer and inner domains were 30 and 10 s, respectively. The main physics schemes used were the single-moment microphysics scheme WSM6 (Hong and Lim, 2006), the Dudhia shortwave radiation scheme (Dudhia, 1989), the RRTM (Rapid Radiative Transfer Model) longwave radiation scheme (Mlawer et al., 1997), and the YSU (Yonsei University) planetary boundary layer scheme (Hong et al., 2006). The Kain–Fritsch cumulus parameterization scheme (Kain, 2004) was utilized in the outer domain, while no cumulus parameterization was utilized in the inner domain. The initial and boundary conditions were from the six-hourly NCEP FNL (Final) analysis data.

 Figure 3 Model domains superimposed with topography (m) for the outer domain (D01) and the inner domain (D02). The black triangles show the locations of SAFIR3000 substations. The red dots show the radar stations used to plot the composite reflectivity in Figs. 4 and 9. Beijing, Tianjin, and Hebei Province are indicated by "BJ", "TJ", and "HB", respectively.

A total of five experiments were conducted, including one control experiment (CTRL) that assimilated no observational data and four other experiments that assimilated lightning data in WRFDA-3DVAR (Table 1). One of the lightning data assimilation experiments was a single time assimilation, and the other three were performed in cycling mode with 10-min intervals but different assimilation time-window lengths. All five experiments used the same physics schemes as described above. The outer domain started at 0000 UTC 10 July, while the inner domain was introduced at 0600 UTC, thus allowing the parent 9-km grids a 6-h spin-up prior to the nested grid initialization. Since the radiation sources began to show up in northwestern Beijing at about 1300 UTC, the start time of the lightning data assimilation in this study was set at 1300 UTC. Limited by the detection range of SAFIR3000, lightning data assimilation was performed only for the inner domain.

Table 1 Description of the experiments
 Experiment name Assimilation window length Number of assimilations Coverage (UTC) Assimilation interval CTRL – – – – aLDA_single – 1 Only once at 1300 – bLDA_30min 30 min 4 1300–1330 10 min cLDA_60min 60 min 7 1300–1400 10 min dLDA_90min 90 min 10 1300–1430 10 min aLDA_single conducted lightning data assimilation only once at 1300 UTC. bLDA_30min conducted lightning data assimilation at 1300, 1310, 1320, and 1330 UTC. cLDA_60min conducted lightning data assimilation at 1300, 1310, 1320, 1330, 1340, 1350, and 1400 UTC. dLDA_90min conducted lightning data assimilation at 1300, 1310, 1320, 1330, 1340, 1350, 1400, 1410, 1420, and 1430 UTC.
4 Results

In this section, we begin by discussing the results of the control experiment. The results from the four lightning data assimilation experiments with different assimilation window lengths are then compared and discussed. Finally, the forecasts based on the LDA_60min lightning data assimilation experiment are further discussed.

4.1 Impact of different assimilation lengths

Considering that convective weather is usually short-lived, it is very important for the subsequent forecast if the assimilation can be completed early. To find the appropriate assimilation time length that yielded significant improvement in both the initial condition and subsequent forecast, four experiments with different assimilation time-window lengths (see Table 1) were conducted. The composite reflectivity of the observation, CTRL, and LDA (lightning data assimilation) experiments were then compared and discussed. Because reflectivity can be adjusted only after model integrating for several minutes (Wang Y. et al., 2014), the reflectivity of 6-min forecasting, rather than the analysis of each LDA experiment, is presented (6-min forecasting was selected for consistency with the observation, since a radar volume scan took 6 min). Figure 4 shows the composite reflectivity of the observations, CTRL, and LDA experiments at 6 min after the analysis time of each LDA experiment, as well as the total flashes used in each LDA experiment. It was found that, at 1306 UTC, there was almost no reflectivity information in both CTRL and LDA_single, and LDA_single showed little or no improvement compared with CTRL. At 1336 UTC, the CTRL experiment failed to simulate the main body and convection core, while the reflectivity of LDA_30min showed improvement compared with CTRL; both the placement and morphology of LDA_30min were closer to observation. However, the area of reflectivity of LDA_30min was much smaller than observed, the placement was a little further south than observed, and the intensity was lower than observed. At 1406 UTC, CTRL simulated only the small part of convection near the adjacent region of Tianjin and Tangshan City, failing to simulate the main part of the observations. The reflectivity of LDA_60min showed significant improvement compared with CTRL, and was in good agreement with observations. In addition, LDA_60min simulated the peripheral stratiform precipitation region around the strong reflectivity well, although the area was somewhat smaller than observed. The convection area simulated by CTRL at 1436 UTC was smaller and situated further south than observed. Although the composite reflectivity of LDA_90min showed salient improvement compared with CTRL, there was no further improvement compared with LDA_60min.

From Fig. 4, we can see that LDA_single, which assimilated lightning data only once, showed little or no improvement compared with CTRL, while the other three experiments (LDA_30min, LDA_60min, and LDA_90min), which assimilated lightning data in cycling mode, showed salient improvement compared with CTRL. Also, the strong echo ( > 35 dBZ) area of the LDA experiments was generally in agreement with the flash distribution used in the assimilation period, indicating that the 3DVAR lightning data assimilation can effectively initiate convection at locations where flashes occur. As the assimilation time-window length increased, the improvement was more pronounced; and yet, after 1-h cycled assimilation (i.e., 7 assimilations), there was no further improvement detected as the assimilation time-window length increased. Moreover, although the stratiform precipitation region extending from the strong reflectivity of the three cycled assimilation experiments (LDA_ 30min, LDA_60min, and LDA_90min) showed salient improvement compared with CTRL, the area was smaller than observed. Two reasons may account for the deficiency in the stratiform precipitation region: firstly, it could be related to the inherent characteristic that flashes primarily occur in strong radar echo areas, yet occur infrequently in stratiform regions (Liu et al., 2013); and secondly, the microphysics scheme used in the present study was single-moment (not double-moment), and the model horizontal resolution was convection-allowing rather than convection-resolving, which may also lead to deficiency in the stratiform region (Bryan and Morrison, 2012; Fierro et al., 2012).

Since the improvement in initial conditions was simi-lar for LDA_60min and LDA_90min, forecasts from LDA_60min and LDA_90min were further compared. The results showed that forecasts from the analyses of LDA_60min and LDA_90min valid at the same time were also generally the same [e.g., 60-min forecast of LDA_60min was similar to 30-min forecast of LDA_90min at 1500 UTC; 90-min forecast of LDA_60min was similar to 60-min forecast of LDA_90min at 1530 UTC; and 120-min forecast of LDA_60min was similar to 90-min forecast of LDA_90min at 1600 UTC (figures omitted)]. The above analysis indicated that the improvements in both the initial condition and the forecast from LDA_60min and LDA_90min were basically the same; considering that LDA_60min needed only a 60-min assimilation window, LDA_60min was obviously preferred, and 60 min was the appropriate assimilation time-window length for this case.

 Figure 4 Composite reflectivity of observation, CTRL (control), and LDA (lightning data assimilation) experiments (6-min forecast), as well as the total flashes used in each LDA experiment. From top to bottom, the time is 1306, 1336, 1406, and 1436 UTC, respectively. Columns 1 to 4 are the composite reflectivity of observation, CTRL, and LDA, and the flashes used in each LDA experiment, respectively. The black line in subplot LDA_60min indicates the position of the cross-section in Fig. 7. In subplot CTRL_1306, Beijing, Tianjin, Hebei, and Tangshan are indicated by "BJ", "TJ", "HB", and "TS", respectively.
4.2 Analysis of increment

Figure 5 shows the horizontal distribution of water vapor mixing ratio (qv) analysis increment (analysis minus background field) of each LDA experiment (at the end of the assimilation window) at the 20th model level (approximately –4.4℃). The analysis increments of qv of all the LDA experiments were positive. The qv analysis increment of LDA_single was the most evident, with a maximum value of 2.4 g kg–1. As the assimilation time-window length increased, the analysis increment decreased. The maximum value of LDA_60min was 0.6 g kg–1, which was 1/4 that of the LDA_single value. The maximum value of LDA_90min was 0.4 g kg–1, which was only 1/6 that of the LDA_single value. This means that the impact of lightning data assimilation on the mo-del background field decreases as the length of the assi-milation period increases. This effect may explain why there was no further improvement in LDA_90min compared with LDA_60min, although the assimilation time window of LDA_90min was 30 min longer than LDA_60min. Since LDA_60min was preferred to LDA_90min, as concluded in Section 4.1, the following analysis focuses only on experiment LDA_ 60min.

 Figure 5 Analysis increment (analysis minus background) of the water vapor mixing ratio (qv) of each LDA (lightning data assimilation) experiment at the 20th model level: (a) LDA_single; (b) LDA_30min; (c) LDA_60min; and (d) LDA_90min. The black line in (c) indicates the position of the cross-section in Fig. 6.

Figure 6 shows the vertical distribution of qv analysis increment (analysis minus background field) of the LDA_60min (at the end of the assimilation window) through the cross-section indicated by the black line in Fig. 5c. The horizontal lines represent the isotherms of 0℃ (solid line) and –20℃ (dashed line). The vertical distribution of qv analysis increment was also positive, and the increment was primarily in the range of 0℃ to –20℃, which was an intentional constraint in this study to adjust the relative humidity. The increment of qv can also be seen beyond 1 km of the upper and lower limit of the designated range, which may be the result of propagation of the background error covariance.

 Figure 6 Vertical distribution of the water vapor mixing ratio (qv) analysis increment of the LDA_60min experiment at 1400 UTC through the cross-section shown in Fig. 5c. The horizontal lines represent the isotherms of 0℃ (solid line) and –20℃ (dashed line).

To better understand the impact of LDA_60min on the microphysical, thermodynamic, and dynamic fields, as well as the mechanism of lightning data assimilation, the vertical distribution of model variable differences between CTRL and LDA_60min (LDA_60min minus CTRL) at 1400 UTC were analyzed. Figure 7 shows the vertical distribution of model variable differences through the cross-section indicated by the black line in Fig. 4. After continuous cycled assimilation (a total of 7 cycles), the impact on graupel was most obvious, with the maximum difference being 4 g kg–1. The impact on ice was least obvious, with the maximum difference being only 0.18 g kg–1. The vertical distribution of the difference in hydrometeors was quite different. The water vapor mixing ratio was mainly in the temperature range of 0℃ to –20℃; the rain water mixing ratio was concentrated below the 0℃ isotherm line; the graupel was mainly in the range of 3–9 km; the cloud water was mainly in the range of 4–6 km; and the snow was mainly in the range of 6–9 km. The difference in virtual temperature between the 0℃ and –20℃ isotherms was mostly positive, with a maximum of 4 K. The difference in vertical velocity was mostly positive above the 0℃ isotherm, and was most obvious in the range of 6–9 km, with a maximum value exceeding 4.5 m s–1. From Fig. 7, the mechanism of lightning data assimilation proposed by Fierro et al. (2012) is clearly substantiated: the increase in relative humidity directly caused the increment of virtual temperature and promoted the condensation of water vapor, accompanied by an increase in latent heat, thus leading to the acceleration of atmospheric updraft and ultimately the development of convection in locations where lightning was observed.

 Figure 7 Vertical distribution of model variable differences between CTRL (control) and LDA_60min (LDA_60min minus CTRL; LDA: lightning data assimilation) at 1400 UTC (the position of the cross-section is indicated by the black line in Fig. 4): (a) water vapor mixing ratio (qv); (b) rain water mixing ratio (qr); (c) graupel mixing ratio (qg); (d) cloud water mixing ratio (qc); (e) ice mixing ratio (qi); (f) snow mixing ratio (qs); (g) virtual temperature (Tv); and (h) vertical velocity (w). The horizontal lines represent the isotherms of 0℃ (solid line) and –20℃ (dashed line).
4.3 Forecast

To examine the improvement in precipitation forecast during the assimilation period (i.e., 1300–1400 UTC), the 1-h accumulated precipitation from the observation, LDA_60min experiment (for convenience, the term LDA_60min and LDA are used interchangeably hereafter), and CTRL was compared (Fig. 8). The 1-h accumulated precipitation of the LDA experiment was obtained by adding the six 10-min accumulated precipitation forecasts using the analysis field at each assimilation time (i.e., 1300, 1310, 1320, 1330, 1340, and 1350 UTC) as the initial conditions. The observation showed an east–west zonal distribution, with two centers located over northeastern Beijing and eastern Hebei Province. In CTRL, only a small part of the rainfall area near the adjacent region of Tianjin and Tangshan was simulated. The LDA experiment showed a salient improvement compared with CTRL: the rainfall area increased significantly and was much closer to observed. The two precipitation centers were also present in the LDA experiment. Nonetheless, both the area and intensity of precipitation was slightly smaller than observed. Overall, the forecast precipitation in LDA during the assimilation period showed significant improvement compared with CTRL.

 Figure 8 One-hour accumulated precipitation during the assimilation period 1300–1400 UTC: (a) observation; (b) CTRL (control); and (c) LDA (lightning data assimilation) experiment. Subplot (c) was obtained by adding the six 10-min accumulated precipitation forecasts using the analy-sis field at each assimilation time (i.e., 1300, 1310, 1320, 1330, 1340, and 1350 UTC) as the initial conditions. Beijing, Tianjin, Hebei, Tangshan, and Chengde are indicated by "BJ", "TJ", "HB", "TS" and "CD", respectively.

To examine the improvement in the forecast after the lightning data assimilation period, the composite reflectivity and the 3-h accumulated precipitation of the observation, CTRL, and LDA experiment were compared. Figure 9 shows the composite reflectivity of the observation, CTRL, and LDA experiment at 1430, 1530, and 1630 UTC. From the observation, it can be seen that the reflectivity was strong and mainly in the vicinity of Beijing and Tianjin City at 1430 UTC, and showed a zonal distribution. Thereafter, the echo moved southeastwards. At 1530 UTC, the echo intensity was lower than at 1430 UTC and the echo split into two parts—one in the north and the other in the south, of which the northern part consisted of two small cells. Then, the echo moved further southeastwards. At 1630 UTC, the two cells in the northern part merged into one, both the echo area and intensity increased, and the southern part split into three small cells. The echo then moved further southeastwards into Bohai Bay and gradually dissipated (figure omitted). However, there were only a few scattered cells at all the three times in CTRL, which were situated more southeastwards than observed. The reflectivity of LDA at 1430 UTC was improved significantly. The morphology and position were basically consistent with observations, although the intensity was slightly weaker than observed. At 1530 UTC, the reflectivity in the LDA experiment split into two parts—one in the north and the other in the south, of which the northern part consisted of two small cells. However, the intensity of the northern part was stronger than the southern part, which was in contrast to observations. At 1630 UTC, the two small cells in the northern part merged into one, and the southern part split into three small cells, yet the intensity of the northern part was again stronger than the southern part. Overall, the forecast reflectivity obtained after cycled assimilation of lightning data improved significantly in position, morphology, and intensity, compared with CTRL. The deve-lopment trend was generally consistent with observations. As the forecast lead time increased, the forecasting ability for reflectivity decreased.

 Figure 9 Composite reflectivity of the observation (left column), CTRL (control; middle column), and LDA (lightning data assimilation) experiment (right column) at 1430 UTC (top row), 1530 UTC (middle row), and 1630 UTC (bottom row), as labeled over the upper left corner of each subplot. Beijing, Tianjin, and Hebei Province are indicated by "BJ", "TJ", and "HB", respectively.

Figure 10 shows the 3-h (1400–1700 UTC) accumulated precipitation from observation, CTRL, and LDA experiments. Precipitation was observed mainly in eastern Beijing, northern Tianjin, and Tangshan City of Hebei Province. Heavy precipitation was observed in eastern Beijing and northern Tianjin, with a maximum of 30 mm. CTRL only predicted the precipitation in the south of Tangshan and the regions adjacent to Tianjin. The rainfall area was obviously too small, and the intensity was weaker than observed. For LDA, the rainfall area was significantly improved compared with CTRL, and was generally in agreement with observation. However, the precipitation center was in the east of Chengde City of Hebei Province, which was not consistent with observation, and the precipitation intensity in eastern Beijing and northern Tianjin was weaker and the rainfall area in eastern Beijing was slightly smaller than observed.

 Figure 10 As in Fig. 8, but for 1400–1700 UTC. Subplot (c) shows the 3-h accumulated precipitation forecast using the analysis field at 1400 UTC as the initial condition.

To quantitatively examine the quality of the forecasts from LDA_60min, the threat score (TS) and bias score (Bias) were used to evaluate the predicted precipitation. To verify the simulated precipitation at the various observation sites, precipitation from the model grid was interpolated to each observation site by using bilinear interpolation. For composite reflectivity, we used the root-mean-square error (RMSE) and spatial correlation coefficient (SC) (Zhang et al., 2012) to evaluate the predicted reflectivity every 6 min. Since the resolution of the radar observations was higher than the model resolution, observed composite reflectivity was interpolated to the model grid by using inverse-distance-squared weighted interpolation.

TS and Bias were computed by using a 2 × 2 contingency table of possible forecast outcomes at individual grid points, wherein the table elements were Hits (correctly forecast events), Misses (observed but not forecast events), False alarms (forecast but not observed events), and Correct negatives (correctly forecast non-events). The following formulae were used to compute TS, Bias, RMSE, and SC:

 ${\rm{T\!S}} = {\rm{Hits}}/\left( {{\rm{Hits}} + {\rm{Misses}} + {\rm{False}}\;{\rm{alarms}}} \right),\quad\,\,\,\,$ (5)
 ${\rm{Bias}} = \left( {{\rm{Hits}} + {\rm{False}}\;{\rm{alarms}}} \right)/\left( {{\rm{Hits}} + {\rm{Misses}}} \right),$ (6)
 ${\rm{RMSE = }}\sqrt {\frac{1}{N}\sum\limits_{i = 1}^N {{{({F_i} - {O_i})}^2}} } , \qquad\qquad\qquad\qquad$ (7)
 ${\rm{SC = }}\frac{{\sum\limits_{i = 1}^N {({F_i} - \overline {{F_i}} )} ({O_i} - \overline {{O_i}} )}}{{\sqrt {\sum\limits_{i = 1}^N {{{({F_i} - \overline {{F_i}} )}^2}} + \sum\limits_{i = 1}^N {{{({O_i} - \overline {{O_i}} )}^2}} } }}.$ (8)

Here, N is the total number of grids; F and O represent the forecast and observed reflectivity, respectively; and the over bar in Eq. (8) represents the spatial average over all the grids. As indicated by Eqs. (5) and (6), TS measures the fraction of observed and forecast events that were correctly predicted, and Bias measures the ratio of the frequency of forecast events to the frequency of observed events.

Figure 11 shows the TS and Bias of the predicted 3-h (1400–1700 UTC) accumulated precipitation of the CTRL and LDA experiment. It can be seen that both the TS and Bias of CTRL were zero when the precipitation threshold was greater than 1 mm. When the threshold was less than or equal to 1 mm, both the Bias and TS were very low ( < 0.12), indicating that CTRL almost had no forecast skill for 3-h precipitation. In contrast, the precipitation prediction of the LDA experiment was improved significantly. When the precipitation threshold was 1 mm, the TS was 0.57, which was much higher than that of CTRL. With the threshold value increased, the TS decreased gradually. When the threshold was greater than 10 mm, the TS was less than 0.1. The Bias result was generally the same as the TS result. When the threshold was 1 mm, the Bias was 0.8, indicating that the predicted rainfall area was comparable to that observed. When the threshold was greater than 5 mm, the Bias decreased to less than 0.5, indicating a smaller predicted rainfall area. Overall, the predicted 3-h precipitation of the LDA experiment showed a notable improvement compared with CTRL. The prediction for smaller threshold values was reasonably credible. As the threshold increased, the prediction became worse. When the threshold was greater than 15 mm, there was no obvious improvement for the LDA experiment.

Figure 12 shows the change in RMSE and SC of the predicted composite reflectivity of the CTRL and LDA experiments during 1400–1700 UTC. It can be seen that the RMSE of CTRL was greater than 10 dBZ during the whole period, and the SC was near zero, indicating that there was almost no spatial correlation between the composite reflectivity of CTRL and the observations. The RMSE of the LDA experiment decreased from 8.8 to 7.8 dBZ during the first 30 min (1400–1430 UTC) of the forecast. During the following 42 min (1430–1512 UTC), the RMSE remained steady near 8 dBZ. The RMSE began to increase after 1512 UTC, reaching 11 dBZ at 1648 UTC, which was equal to that of CTRL. It then further increased, and exceeded the value in CTRL at 1700 UTC. The SC of the LDA experiment increased from 0.52 to 0.65 during the first 30 min (1400–1430 UTC) of the forecast. During the following 42 min (1430–1512 UTC), the SC remained near 0.65, and began to decrease after 1512 UTC, dropping to 0.28 at 1700 UTC. From the above analysis, we can see that the predicted reflectivity was optimal after 30 min from the start of the forecast, which was due to the fact that the model was not quite dynamically consistent with the analysis field at first, needing some time to spin up. During the following 42 min, the prediction remained steady and was reasonably credible. From 1512 UTC, i.e., 72 min after the start of the forecast, with the continuously accumulated model error as well as the pre-existing errors on synoptic scales, the positive effect obtained from the lightning data assimilation gradually diminished. At 1700 UTC, the SC was only 0.28, and the RMSE was roughly equivalent to that in CTRL. In view of the above analyses, both qualitatively and quantitatively, the improvement from assimilating lightning data could be maintained for about 3 h.

 Figure 11 The TS (threat score; top) and Bias (bias score; bottom) of the 3-h (1400–1700 UTC) accumulated precipitation of the LDA (lightning data assimilation) and CTRL (control) experiments.
 Figure 12 RMSE (root-mean-square error) and SC (spatial correlation coefficient) of the composite reflectivity of the LDA (lightning data assimilation) and CTRL (control) experiments during 1400–1700 UTC.
5 Summary and discussion

Lightning data have been widely recognized as a valuable proxy variable for identifying the occurrence of deep convection (Petersen and Rutledge, 1998; Schultz et al., 2011; Rudlosky and Fuelberg, 2013). Therefore, lightning data have the potential to be applied in the monitoring, early warning, and forecasting of severe convective weather, and are especially important over oce-anic and terrain-blocked areas where radar reflectivity data are generally not available. Assimilating lightning data is of great importance for numerical models because it enables more small-scale information to be incorporated, thus improving the quality of the initial conditions and the subsequent forecasts. The relationship between the water vapor mixing ratio, flash rate, and graupel mixing ratio (Fierro et al., 2012) was used in the present study to adjust the model relative humidity, which was then assimilated into WRF-3DVAR in a cycling mode at 10-min intervals. The method was applied to a squall line case that occurred on 10 July 2007 in North China. To find the appropriate time-window length for the cycled assimilation window, with significant improvement both in the initial field and subsequent forecast, four experiments with different assimilation window lengths (including one single time assimilation) were conducted. After finding the appropriate time length of the assimilation window, the improvement in forecast precipitation and reflectivity was further analyzed. The main conclusions of the study can be summarized as follows.

(1) It is feasible and effective to assimilate the relative humidity retrieved from lightning data by utilizing the 3DVAR method, with significant improvement in both the initial condition and subsequent forecasts.

(2) The experiment in which lightning data were assimilated only once (LDA_single) showed little or no improvement in the initial condition. The other three experiments, in which lightning data were assimilated in cycling mode, showed varying degrees of improvement. The improvement in the convective zone was particularly evident, while the region of simulated stratiform precipitation was smaller than observed. From the comparative analyses, it was concluded that 60 min was the appropriate assimilation window length for the studied case, and longer assimilation window length was probably unnecessary.

(3) The prediction of 1-h accumulated precipitation during the lightning data assimilation period and the subsequent 3-h accumulated precipitation showed significant improvement. As forecast lead time increased, the forecast for reflectivity became worse. Based on subjective comparisons and objective statistical scores, the improvement from lightning data assimilation could be maintained for about 3 h.

As mentioned in the introduction, previous efforts related to lightning data assimilation have mostly employed the nudging method, which substitutes the model variables directly with observations. The variational me-thod, which considers the system as a whole and combines the model and observations in a statistically opti-mal sense, is a more sophisticated and robust approach, and is currently widely used in research communities and operational centers. In the present study, lightning data assimilated by using the 3DVAR method showed significant improvement, with the assimilation window length being only 1 h, which may be an advantage over the nudging method.

However, some deficiencies nevertheless exist. Although the prediction of the rainfall area was generally in agreement with observations, the prediction of the precipitation centers and the intensity of reflectivity were not satisfactory. One possible reason for this is that the present study only assimilated the relative humidity retrieved from lightning data, which has limited ability in terms of improving the flow field of the model. Combining the assimilation of lightning data with radar radial velocity data can be considered in the future to improve the simulation of the flow field, and consequently improve the prediction of precipitation centers. Moreover, it could be seen from both the analysis fields and the forecasts that the stratiform region was generally smaller than observed. The deficiency in the stratiform region may be related to the inherent characteristic of lightning data that flashes primarily occur in strong radar echo areas, yet occur infrequently in stratiform regions (Liu et al., 2013). In future studies, we hope that the use of double-moment microphysics schemes and higher model resolution can partially relieve the overall deficiency in the simulated stratiform region. Lastly, but importantly, one inherent defect of the 3DVAR method is that it assumes the background error covariance is stationary, nearly homogeneous, and isotropic, while in fact the background error covariance may vary substantially with the flow of the atmosphere. A more sophisticated and advanced data assimilation method, such as EnKF (Mansell, 2014; Allen et al., 2016) or the hybrid variational–ensemble approach (Wang et al., 2008a, b), which can approximately obtain the flow-dependent background error covariance from an ensemble of background states, will hopefully further improve the effect of lightning data assimilation and prolong the forecast lead time.

As total lightning data can reflect convection activities more comprehensively than CG lightning data, the assimilation of total lightning data is strongly preferred relative to the assimilation of only CG lightning data. However, most current operational ground-based lightning location systems can only detect CG flashes, which may impede the wide use of lightning data assimilation. Fortunately, continuous total lightning data will soon be available through the lightning mapping imager onboard the geostationary FY-4 satellite (Cao et al, 2014), which was successfully launched at the end of 2016. By that time, the data assimilation method of this study will be further tested and refined. Since the present study focused on only one squall line case, the results may not be completely applicable to cases of other weather systems. In the future, other types of convective systems, such as tropical cyclones and mesoscale convective complexes, should also be investigated.

Acknowledgments. We thank Dr. Xinghua Bao of the Chinese Academy of Meteorological Sciences for her constructive suggestions.

References
 Alexer G. D., Weinman J. A., Karyampudi V. M., et al., 1999: The effect of assimilating rain rates derived from satellites and lightning on forecasts of the 1993 superstorm. Mon. Wea. Rev., 127, 1433–1457. DOI:10.1175/1520-0493(1999)127 < 1433:TEOARR > 2.0.CO; 2 Allen B. J., Mansell E. R., Dowell D. C., et al., 2016: Assimilation of pseudo-GLM data using the ensemble Kalman filter. Mon. Wea. Rev., 144, 3465–3486. DOI:10.1175/MWR-D-16-0117.1 Barker, D. , W. Huang, Y. -R. Guo, et al. , 2003: A Three-dimensional Variational (3DVAR) Data Assimilation System for Use with MM5. NCAR Tech Note, Boulder, CO, NCAR, 1-68. Barker D. M., Huang W., Guo Y.-R., et al., 2004: A three-dimensional variational data assimilation system for MM5: Implementation and initial results. Mon. Wea. Rev., 132, 897–914. DOI:10.1175/1520-0493(2004)132 < 0897:ATVDAS > 2.0.CO; 2 Barker D., Huang X.-Y., Liu Z. Q., et al., 2012: The weather research and forecasting model's community variational/ensemble data assimilation system: WRFDA. Bull. Amer. Me-teor. Soc., 93, 831–843. DOI:10.1175/BAMS-D-11-00167.1 Bryan H., G. Morrison, 2012: Sensitivity of a simulated squall line to horizontal resolution and parameterization of microphysics. Mon. Wea. Rev., 140, 202–225. DOI:10.1175/MWR-D-11-00046.1 Cao, D. J. , F. X. Huang, and X. S. Qie, 2014: Development and evaluation of detection algorithm for FY-4 Geostationary Lightning Imager (GLI) measurement. XV International Conference on Atmospheric Electricity, Norman, Oklahoma, 15-20 June, International Commission on Atmospheric Electricity, 1-6. Chang D.-E., Weinman J. A., Morales C. A., et al., 2001: The effect of spaceborne microwave and ground-based continuous lightning measurements on forecasts of the 1998 Groundhog Day storm. Mon. Wea. Rev., 129, 1809–1833. DOI:10.1175/1520-0493(2001)129 < 1809:TEOSMA > 2.0.CO; 2 Dudhia J., 1989: Numerical study of convection observed during the winter monsoon experiment using a mesoscale two-dimensional model. J. Atmos. Sci., 46, 3077–3107. DOI:10.1175/1520-0469(1989)046 < 3077:NSOCOD > 2.0.CO; 2 Fierro A. O., Clark A. J., Mansell E. R., et al., 2015: Impact of storm-scale lightning data assimilation on WRF-ARW precipitation forecasts during the 2013 warm season over the contiguous United States. Mon. Wea. Rev., 143, 757–777. DOI:10.1175/MWR-D-14-00183.1 Fierro A. O., Gao J. D., Ziegler C. L., et al., 2014: Evaluation of a cloud-scale lightning data assimilation technique and a 3DVAR method for the analysis and short-term forecast of the 29 June 2012 derecho event. Mon. Wea. Rev., 142, 183–202. DOI:10.1175/MWR-D-13-00142.1 Fierro A. O., Gao J. D., Ziegler C. L., et al., 2016: Assimilation of flash extent data in the variational framework at convection-allowing scales: Proof-of-concept and evaluation for the short-term forecast of the 24 May 2011 tornado outbreak. Mon. Wea. Rev., 144, 4373–4393. DOI:10.1175/MWR-D-16-0053.1 Fierro A. O., Mansell E. R., Ziegler C. L., et al., 2012: Application of a lightning data assimilation technique in the WRF-ARW model at cloud-resolving scales for the tornado outbreak of 24 May 2011. Mon. Wea. Rev., 140, 2609–2627. DOI:10.1175/MWR-D-11-00299.1 Goodman S. J., Blakeslee R. J., Koshak W. J., et al., 2013: The GOES-R geostationary lightning mapper (GLM). Atmos. Res., 125–126, 34–49. DOI:10.1016/j.atmosres.2013.01.006 Goodman S. J., Gurka J., DeMaria M., et al., 2012: The GOES-R proving ground: Accelerating user readiness for the next-generation geostationary environmental satellite system. Bull. Amer. Meteor. Soc., 93, 1029–1040. DOI:10.1175/BAMS-D-11-00175.1 Hong J. Lim, 2006: The WRF single-moment 6-class microphsics scheme (WSM6). J. Korean Meteor. Soc., 42, 192–151. Hong Noh, S.-Y. Dudhia, 2006: A new vertical diffusion package with an explicit treatment of entrainment processes. Mon. Wea. Rev., 134, 2318–2341. DOI:10.1175/MWR3199.1 Kain J. S., 2004: The Kain–Fritsch convective parameterization: An update. J. Appl. Meteor., 43, 170–181. DOI:10.1175/1520-0450(2004)043 < 0170:TKCPAU > 2.0.CO; 2 Lagouvardos K., Kotroni V., Defer E., et al., 2013: Study of a heavy precipitation event over southern France, in the frame of HYMEX project: Observational analysis and model results using assimilation of lightning. Atmos. Res., 134, 45–55. DOI:10.1016/j.atmosres.2013.07.003 Li, W. B. , G. Q. Song, and K. Tong, 2008: Applying TRMM-LIS lightning data to numerical model. Acta Sci. Natur. Univ. Pekinensis, 44, 399-406. (in Chinese) Liu D. X., Qie X. S., Pan L. X., et al., 2013: Some characteristics of lightning activity and radiation source distribution in a squall line over North China. Atmos. Res., 132–133, 423–433. DOI:10.1016/j.atmosres.2013.06.010 Lynn H., B. Kelman, G. Ellrod, 2015: An evaluation of the efficacy of using observed lightning to improve convective lightning forecasts. Wea. Forecasting, 30, 405–423. DOI:10.1175/WAF-D-13-00028.1 MacGorman D. R., Burgess D. W., Mazur V., et al., 1989: Lightning rates relative to tornadic storm evolution on 22 May 1981. J. Atmos. Sci., 46, 221–250. DOI:10.1175/1520-0469(1989)046 < 0221:LRRTTS > 2.0.CO; 2 Mansell E. R., 2014: Storm-scale ensemble Kalman filter assimilation of total lightning flash-extent data. Mon. Wea. Rev., 142, 3683–3695. DOI:10.1175/MWR-D-14-00061.1 Mansell R., E. L. Ziegler, C. R. MacGorman, 2007: A lightning data assimilation technique for mesoscale forecast models. Mon. Wea. Rev., 135, 1732–1748. DOI:10.1175/MWR3387.1 Marchand R., M. E. Fuelberg, 2014: Assimilation of lightning data using a nudging method involving low-level warming. Mon. Wea. Rev., 142, 4850–4871. DOI:10.1175/MWR-D-14-00076.1 Mlawer E. J., Taubman S. J., Brown P. D., et al., 1997: Radiative transfer for inhomogeneous atmospheres: RRTM, a validated correlated-k model for the longwave. J. Geophys. Res., 102, 16663–16682. DOI:10.1029/97JD00237 Papadopoulos G. Chronis, A. N. Anagnostou, 2005: Improving convective precipitation forecasting through assimilation of regional lightning measurements in a mesoscale mo-del. Mon. Wea. Rev., 133, 1961–1977. DOI:10.1175/MWR2957.1 Papadopoulos Serpetzoglou, A. N. Anagnostou, 2009: Evaluating the impact of lightning data assimilation on mesoscale model simulations of a flash flood inducing storm. Atmos. Res., 94, 715–725. DOI:10.1016/j.atmosres.2009.05.008 Parrish F., D. C. Derber, 1992: The National Meteorologi-cal Center's spectral statistical-interpolation analysis system. Mon. Wea. Rev., 120, 1747–1763. DOI:10.1175/1520-0493(1992)120 < 1747:TNMCSS > 2.0.CO; 2 Pessi T., A. Businger, 2009: The impact of lightning data assimilation on a winter storm simulation over the North Pacific Ocean. Mon. Wea. Rev., 137, 3177–3195. DOI:10.1175/2009MWR2765.1 Petersen A., W. A. Rutledge, 1998: On the relationship between cloud-to-ground lightning and convective rainfall. J. Geophys. Res., 103, 14025–14040. DOI:10.1029/97JD02064 Qie S., X. X. Liu, D. L. Sun, 2014a: Recent advances in research of lightning meteorology. J. Meteor. Res., 28, 983–1002. DOI:10.1007/s13351-014-3295-0 Qie X. S., Zhu R. P., Yuan T., et al., 2014b: Application of total-lightning data assimilation in a mesoscale convective system based on the WRF model. Atmos. Res., 145–146, 255–266. DOI:10.1016/j.atmosres.2014.04.012 Ran, L. K. , and Y. S. Zhou, 2011: Application of lightning observations of TRMM satellite to the mesoscale numerical model by a new nudging assimilation adjustment technique. Chinese J. Atmos. Sci. , 35, 1145-1158. (in Chinese) Rudlosky D., S. E. Fuelberg, 2013: Documenting storm severity in the mid-Atlantic region using lightning and radar information. Mon. Wea. Rev., 141, 3186–3202. DOI:10.1175/MWR-D-12-00287.1 Schultz J., C. A. Petersen, W. D. Carey, 2011: Lightning and severe weather: A comparison between total and cloud-to-ground lightning trends. Wea. Forecasting, 26, 744–755. DOI:10.1175/WAF-D-10-05026.1 Wang H. L., Sun J. Z., Fan S. Y., et al., 2013: Indirect assimilation of radar reflectivity with WRF 3D-Var and its impact on prediction of four summertime convective events. J. Appl. Meteor. Climatol., 52, 889–902. DOI:10.1175/JAMC-D-12-0120.1 Wang X. G., Barker D. M., Snyder C., et al., 2008a: A hybrid ETKF–3DVAR data assimilation scheme for the WRF model. Part Ⅰ: Observing system simulation experiment. Mon. Wea. Rev., 136, 5116–5131. DOI:10.1175/2008MWR2444.1 Wang X. G., Barker D. M., Snyder C., et al., 2008b: A ETKF–3DVAR data assimilation scheme for the WRF model. Part Ⅱ: Real observation experiments. Mon. Wea. Rev., 136, 5132–5147. DOI:10.1175/2008MWR2445.1 Wang Yang, Y. H. Wang, 2014: Improving forecasting of strong convection by assimilating cloud-to-ground lightning data using the physical initialization method. Atmos. Res., 150, 31–41. DOI:10.1016/j.atmosres.2014.06.017 Wang, Y. D. , Y. J. Zhou, X. Y. Wang, et al. , 2014: A study on the assimilation method of lightning data with mesoscale model WRF. J. Trop. Meteor. , 30, 281-292. (in Chinese) Xu S., Zheng D., Wang Y. Q., et al., 2016: Characteristics of the two active stages of lightning activity in two hailstorms. J. Meteor. Res., 30, 265–281. DOI:10.1007/s13351-016-5074-6 Zhang Q. Tian, D.-L. Yang, 2011: Genesis of Typhoon Nari (2001) from a mesoscale convective system. J. Geophys. Res., 116, D23104. DOI:10.1029/2011JD016640 Zhang C. X., Wang Y. Q., Lauer A., et al., 2012: Configuration and evaluation of the WRF model for the study of Hawaiian regional climate. Mon. Wea. Rev., 140, 3259–3277. DOI:10.1175/MWR-D-11-00260.1 Zheng D., Zhang Y. J., Meng Q., et al., 2009: Total lightning characteristics and electric structure evolution in a hailstorm. Acta Meteor. Sinica, 23, 233–249.