# Assimilation of satellite lightning data in a storm-scale numerical weather prediction system using a 3D-EnVar data assimilation scheme

Pauline Combarnous, M. Martet, Olivier Caumont, Éric Defer

► To cite this version:

Pauline Combarnous, M. Martet, Olivier Caumont, Éric Defer. Assimilation of satellite lightning data in a storm-scale numerical weather prediction system using a 3D-EnVar data assimilation scheme. Monthly Weather Review, 2024, 152 (2), pp.553-570. 10.1175/MWR-D-23-0100.1. hal-04434925

HAL Id: hal-04434925

https://hal.science/hal-04434925v1

Submitted on 2 Feb 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Assimilation of satellite lightning data in a storm-scale numerical weather prediction system using a 3D-EnVar data assimilation scheme

Pauline Combarnous<sup>a</sup>, Maud Martet<sup>a</sup>, Olivier Caumont<sup>a,b</sup>, Éric Defer<sup>c</sup>

<sup>a</sup> CNRM, Université de Toulouse, Météo-France, CNRS, Toulouse, France

<sup>b</sup> Météo-France, Direction des opérations pour la prévision, Toulouse, France

<sup>c</sup> LAERO, Université de Toulouse, UT3, CNRS, IRD, Toulouse, France

<footnote>Corresponding author: Pauline Combarnous, pauline.combarnous@meteo.fr</footnote>
ABSTRACT: The objective of this study is to prepare the assimilation of the Meteosat Third Generation Lightning Imager (MTG-LI) observations in the regional numerical weather prediction model AROME-France using a 3D-EnVar data assimilation algorithm. As opposed to the current operational configuration of AROME-France data assimilation system, the hydrometeor specific contents are added in the control variable and are thus updated during the assimilation process. Consequently, a lightning observation operator based on the specific contents of snow and graupel is used. As the first MTG satellite was launched in December 2022, the real MTG-LI observations are not available yet, and proxy space-borne lightning observations generated from ground-based lightning detection network are used. We assess the forecast skills of experiments assimilating lightning in a close to operational configuration by comparing the results to a reference experiment that does not assimilate lightning. In the studied cases, the frequency bias and fraction skill score are improved for brightness temperatures forecasts lower than 280 K up to 4 hours after the assimilation, implying a better description of the general cloud cover. The convective cores, identified by brightness temperatures lower than 220 K, are better captured in the first forecast hour but this improvement is not maintained through time and quickly fades.
# 1. Introduction

Convective-scale numerical weather prediction (NWP) systems are now widely used to produce high resolution short term forecasts (Benjamin et al. 2019). However, the location and timing of severe event forecasts still have a margin of progress and could benefit from a better description of convection-related initial conditions through the assimilation of high spatio-temporal resolution observations.

Total lightning, i.e. cloud-to-ground (CG) and intra-cloud (IC) lightning, presents the advantage of being strongly related to atmospheric convection (Price and Rind 1992; Deierling et al. 2008) and microphysical contents (Goodman et al. 1988; Deierling et al. 2005) while being measured continuously at a high spatio-temporal resolution. Since 2016, total lightning can be observed from a geostationary orbit, with the Geostationary Lightning Mapper (GLM) on board the GOES16/17 R-series satellites (Goodman et al. 2013) and with the Lightning Mapping Imager on board FengYun-4A (Cao et al. 2021). Space-borne lightning detection offers the opportunity to monitor thunderstorms in data-sparse areas such as oceans, mountains and regions without any coverage of ground-based radar network. This makes it a rich potential source for NWP systems initialization at convective scale.

Lightning data assimilation (LDA) has already been investigated in several studies, with different assimilation approaches at convection-permitting resolution. Fierro et al. (2012) developed a method consisting in adjusting water vapor mixing ratio or temperature through nudging techniques where lightning is observed. They improved the simulated precipitation despite introducing a wet bias in the forecasts. This nudging method has been widely used since, adjusting either thermodynamic or hydrometeor variables linked to lightning (Marchand and Fuelberg 2014; Qie et al. 2014; Chen et al. 2019; Federico et al. 2019). Some variational approaches were also studied, such as a three-dimensional variational (3D-Var) assimilation of lightning-derived water vapor mixing ratio (Fierro et al. 2016; Hu et al. 2020; Erdmann et al. 2023). Even though a general improvement was observed in composite reflectivity or rainfall accumulation for short term forecasts (up to 3 hours), a wet bias was still present. A 4D-Var assimilation method of lightning-derived vertical velocities was developed by Xiao et al. (2021) and led to similar conclusions regarding very-short term forecast improvements and a wet bias.
The chief limitation of conventional 3/4D-Var methods is that they assume the background error covariances to be static (constant over time), isotropic and homogeneous. However, it has been shown by Ménétrier et al. (2014) that forecast errors at convective scales are strongly flow-dependent. Ensemble-based data assimilation algorithms, such as the ensemble Kalman filter (EnKF) or the ensemble variational algorithm (EnVar) use flow-dependent covariances updated at each analysis and estimated from an ensemble of forecasts from a previous analysis, instead of climatological covariances. In the field of LDA, methods using an EnKF have been developed by Mansell (2014) and Allen et al. (2016). The results were promising, showing a modulation in the intensity of simulated convection and a reduction of spurious deep convection. More recently, Kong et al. (2020) also used the EnKF to assimilate lightning data within an operational framework, building on Allen et al. (2016)'s work, and no wet bias was seen.

The present study aims at assimilating the future Meteosat Third Generation Lightning Imager (MTG-LI; Kokou et al. 2018) observations within the new 3D-EnVar data assimilation system of the French regional NWP system AROME-France. The recent launch of the first MTG satellite carrying the LI instrument will give us the opportunity to monitor continuously total lightning activity above Europe, the Mediterranean Sea, Africa and the Atlantic Ocean at a horizontal resolution of a few kilometers. Our long-term objective is the assimilation of those observations, therefore the data used here are generated from ground-based lightning detection network to simulate space-borne lightning observations (Erdmann et al. 2022) since LI records are not yet available.

LDA is particularly challenging because of the discrete nature of lightning and the non-linearity of its observation operator when lightning is not a prognostic variable. In variational DA, a linear tangent version of the observation operator is required, but its linearization can induce convergence issues during the minimization process when the observation operator is non-linear. This issue is addressed in this article and a solution based on a variable transformation is introduced. Furthermore, because lightning is a discrete variable whose values can be classified in a binary way (lightning/no-lightning), four possible configurations emerge when observations are compared to model background. The first two configurations are symmetric conditions: lightning/no-lightning in both observations and the model background. The other two configurations are asymmetric conditions: lightning is present either in the observation or in the model background but not in
82

both. The two latter configurations are studied in more detail in this article to examine how the
assimilation system deals with this asymmetry. Particular attention is given to the configuration
where lightning is observed but none is forecasted in the model background, because of the "zero
gradient" problem raised in several cloud and precipitation DA studies (e.g. Errico et al. 2007;
Lopez 2011), and solved with the variable change mentioned above.

83

84

85

86

87

The paper is organized as follows. Section 2 introduces the 3D-EnVar data assimilation system and the lightning data. Section 3 presents the lightning observation operator and discusses the non-linearity issues described above. In section 4, the design of the data assimilation experiments is described. The impact of LDA on the analysis fields and the forecasting skills is presented in section 5. A summary of the main conclusions and a discussion on future studies is given in section 6.

88

89

90

91

92

93

## 2. Assimilation system and lightning data

94

In this section, the NWP model used in this study is introduced as well as its data assimilation
algorithm. This section also contains a description of the lightning data and how they are generated
from a ground-based lightning detection network.

95

96

97

### a. AROME-France 3D-EnVar data assimilation system

98

The AROME-France model has been used operationally since 2008 (Seity et al. 2011). It has a
horizontal resolution of 1.3 km over a limited geographical domain shown in Fig. 1a. It has 90
vertical levels with a model top at 10 hPa, but the levels are not spaced evenly and the majority
of them are situated within the troposphere. It uses the single-moment microphysical scheme
ICE3 that predicts five hydrometeor species: the specific contents of rain, snow, graupel, cloud
droplets and ice crystals. The remaining prognostic variables are the temperature, two components
of the horizontal wind, specific content of water vapor, surface pressure, turbulent kinetic energy,
and two non-hydrostatic variables. Currently, the AROME-France data assimilation system is
a sequential 3D-Var data assimilation system (Brousseau et al. 2016) but a 3D-EnVar system is
under preparation for operational purpose and a prototype is used in this study, operating in an
hourly cycle (Brousseau et al. to be submitted). A 1 h forecast from the previous assimilation
hour is used as background, i.e. the model state to be blended with observations during the

100

101

102

103

104

105

106

107

108

109
analysis process. Unlike the current operational 3D-Var system, the prototype used here includes the hydrometeors in the control variable, i.e. the ensemble of model variables updated during the data assimilation process (Destouches et al. 2023). Indeed, in the current operational system, for scientific and technical reasons, the hydrometeors are not updated during the assimilation process and the analysis values are the same as the background values.

The aim of the 3D-Var is to retrieve the optimal analysis state $\mathbf{x}^a$, given by

$$ \mathbf{x}^a = \mathbf{x}^b + \delta\mathbf{x}, \qquad (1) $$

where $\mathbf{x}^b$ is the background state and $\delta\mathbf{x}$ the analysis increment. The analysis increment minimizes the cost function $J$:

$$ J(\delta\mathbf{x}) = \frac{1}{2}\delta\mathbf{x}^T\mathbf{B}^{-1}\delta\mathbf{x} + \frac{1}{2}(\mathbf{d} - \mathbf{H}\delta\mathbf{x})^T\mathbf{R}^{-1}(\mathbf{d} - \mathbf{H}\delta\mathbf{x}) \qquad (2) $$

in which vector $\mathbf{d}$ is the innovation vector, i.e. the difference between the observation vector $\mathbf{y}^o$ and their background counterpart $\mathcal{H}(\mathbf{x}^b)$ where $\mathcal{H}$ represents the nonlinear observation operators and $\mathbf{H}$ their corresponding linearized versions. $\mathbf{R}$ and $\mathbf{B}$ are, respectively, the observation- and the background-error covariance matrices. These background-error covariances describe how an increment induced by the assimilation of an observation will propagate spatially and among all the variables, even those that do not appear in the observation operator.

The principle of a 3D-EnVar data assimilation system is to use covariances estimated from an ensemble of background states to produce one analysis from one background and these covariances are recalculated for each analysis. It is thus a deterministic method, coupled to an ensemble data assimilation (EDA) cycle to estimate the flow-dependent covariances. The climatological $\mathbf{B}$ matrix introduced in Eq. (2) and used in the 3D-Var algorithm consists of homogeneous and isotropic background-error covariances (Brousseau et al. 2011), whereas the ensemble background-error covariance $\tilde{\mathbf{B}}$ matrix varies in each grid point of the model, accordingly to the current atmospheric conditions. In our simulations, the ensemble is composed of 50 members, using a horizontal grid spacing of 3.2 km. The reader is invited to consult Montmerle et al. (2018) for a more detailed description of the 3D-EnVar formulation and Michel and Brousseau (2021) for a use of
<figure><img src="image_2.png" /><figcaption>FIG. 1. (a) AROME-France domain (red box) and the sub-domains where the events of 3-4 October (SE France, black dashed) and 5-17 November (BLAS, black solid) are studied. (b) Observed 24 h accumulated rainfall at 18:00 UTC on 4 October 2021 in the SE France domain. Locations mentioned in the text are also indicated in italics.</figcaption></figure>

the algorithm in the same configuration as the operational one (i.e. hourly cycle and a horizontal grid spacing of 1.3 km).

An important factor in ensemble data assimilation is localization. Localization is necessary to cut off spurious correlations due to sampling noise in the estimation of $\tilde{B}$. This sampling noise is introduced because the number of ensemble members (50) used to calculate $\tilde{B}$ is less than the size of the state vector ($10^9$), which produces spuriously large magnitude background error covariance estimates between distant grid points. Then, localization aims at mitigating these spurious correlations by setting to zero the sample covariances at long range. Here, the horizontal localization length was chosen to range from 25 km at the lowest level of the model and increases linearly with the altitude to reach 150 km at the highest level. The vertical localization length is set to 0.3 units of the natural logarithm of pressure, as in Montmerle et al. (2018) and Caron et al. (2019) and similar to the value 0.4 used by Destouches et al. (2023). These localization lengths were selected out of sensitivity studies performed on the new 3D-EnVar scheme for AROME-France (Brousseau et al. to be submitted).
## b. Lightning data

The LI instrument is an optical sensor composed of 4 cameras, on board a geostationary satellite, first of the MTG mission, whose field of view covers Europe, Africa, Atlantic Ocean, the Mediterranean Sea and a part of South-America. The instrument detects the illumination of the cloud by lightning flashes at 777.4 nm, with a horizontal resolution of 7 km at Western European latitudes and 4.5 km at nadir and provides the three same products as the GLM, i.e. events, groups and flashes. One of the LI products is a count of flashes, including both IC and CG, per grid cell over a certain accumulation period, designated as the flash extent accumulation (FEA), measured in flashes (fl).

The geostationary lightning pseudo-observation generator developed by Erdmann et al. (2022) is used here to build data mimicking the future MTG-LI data. More specifically, a machine learning regressor was trained with coincident National Lightning Detection Network (NLDN) and Geostationary Lightning Mapper (GLM) data to learn how to reproduce flash characteristics as recorded from space from ground-based measured flashes. The generator was applied to lightning measured by the French low-frequency ground-based network Météorage to obtain the synthetic MTG-LI data, generated with a horizontal resolution of 7 km. Because MTG-LI will measure lightning similarly to GLM, and Météorage to the NLDN (Erdmann 2020, chapter II.2.4), the generator can indeed be used in this context. For more details on the method and its evaluation, please see Erdmann et al. (2022).

# 3. The lightning observation operator

In the following, we introduce and discuss the observation operator developed in a previous study and the modifications brought to overcome linearity issues and the "zero gradient" problem which reflects the impossibility of increasing hydrometeor contents when none is forecasted, due to a lack of sensitivity in the observation operator.

## a. Description

Direct assimilation of lightning data requires an observation operator that links the lightning observations to the predicted model variables. This observation operator allows the calculation of simulated observations from the model state. Combarnous et al. (2022) investigated a total of eight
parameters to determine which ones are best related to synthetic space-borne lightning observations.
Data from 44 stormy days were used to find the relationship between FEA observations and the
parameters, representing the annual distribution of thunderstorms in the AROME-France domain.
Several regression techniques were tested but a cubic polynomial regression was the one presenting
the best compromise between performances and implementation feasibility as an observation
operator for data assimilation. A cubic polynomial regression was also tested by Kong et al. (2022)
and was found to outperform forecasts relative to a linear operator. As for the parameters, it was
shown in Combarnous et al. (2022) that the ones based on microphysical contents were the most
successful at reproducing FEA areal coverage and amplitude, versus the ones based on updraft
characteristics.

For the present study, the parameter selected to calculate the background model equivalent of
the FEA with the observation operator is the vertically integrated precipitation ice mass, referred
to as the ice water path (IWP,kg m⁻²). The ability of the IWP to reproduce lightning is in
accordance with the non-inductive charging theory that suggests a correlation between the amount
of precipitation-sized ice masses and lightning production (MacGorman and Rust 1998). Petersen
et al. (2005) studied on the global scale the relationship between TRMM lightning and IWP derived
from radar observations and demonstrated that this theory is robust. Additionally, the relationship
they found between lightning and IWP is relatively invariant above land and ocean. The choice of
the IWP is in line with Combarnous et al. (2022)'s recommendation: because the IWP is calculated
with two different hydrometeor types, graupel and snow, it increases the chances to have a non-
zero total microphysical content in the background. The IWP is calculated column-wise and is
derived from the surface density of the specific contents of precipitating ice masses, snow qₛ and
graupel qₛ (specific contents are both in kg kg⁻¹), summed over each level from the level where
the temperature T goes below -10°C to the roof of the model (10 hPa) as follows:

$$ \text{IWP} = \sum_{T<-10^{\circ}\text{C}} \rho \times (q_s + q_g) \times \Delta z, \qquad (3) $$
<table><caption>TABLE 1. Coefficients values of the observation operator.</caption><thead><tr><th>Coefficient</th><th>Value</th></tr></thead><tbody><tr><td>A</td><td>-3.6161 × 10<sup>-2</sup></td></tr><tr><td>B</td><td>1.7296</td></tr><tr><td>C</td><td>-1.1723</td></tr><tr><td>D</td><td>2.2360 × 10<sup>-2</sup></td></tr><tr><td>IWP<sub>TH</sub></td><td>1.1844</td></tr><tr><td>a</td><td>0.60193</td></tr></tbody></table>

where $\rho$ is the local air density (kg m$^{-3}$) and $\Delta z$ is the thickness of level (m). The forward lightning observation operator can be expressed as follows:

$$FEA = A \times IWP^3 + B \times IWP^2 + C \times IWP + D \quad (4)$$

Values of coefficients $A$, $B$, $C$ and $D$ are listed in Table 1. Equation (4) is plotted in black in Fig. 2a. Because this function has a positive intercept (value of $D$) and is negative for IWP values in [0.019; 0.67] kg m$^{-2}$, the polynomial is replaced by a basic third power function for IWP values below a threshold. This threshold, $IWP_{TH}$, is selected so that the corresponding FEA value is 1, i.e. $FEA(IWP_{TH}) = 1$ fl. The final observation operator, plotted in red in Fig. 2a, is the following piece-wise defined function:

$$FEA = \begin{cases} a \times IWP^3 & \text{if } IWP \le IWP_{TH} \\ \text{Eq. (4)} & \text{if } IWP > IWP_{TH} \end{cases} \quad (5)$$

For continuity at the value $IWP_{TH}$, the coefficient $a$ is equal to $\frac{1}{(IWP_{TH})^3}$.

## b. Issues with non-linearity of the observation operator

As mentioned in the introduction, the linearization of a non-linear observation operator during the minimization process can lead to convergence problems. To get a better estimate of the linearized operator, we applied a cubic-root transformation to the FEA: $\sqrt[3]{FEA}$, inspired by a similar variable change to radar reflectivity in Hawkness-Smith and Simonin (2021)'s work. The function after transformation is plotted in blue in Fig. 2b. After this transformation, the forward operator for
<figure><img src="image_3.png" /><figcaption>FIG. 2. Lightning observation operator functions: (a) before and (b) after the cubic-root transformation. Note that the IWP ranges between 0 and 2.5 kg m⁻² for a better visualization of the operator behavior for low IWP values, but typical values range between 0 and 25 kg m⁻².</figcaption></figure>

IWP values below $IWP_{TH}$ becomes a simple linear function of IWP and its linear tangent term is a constant value equals to $\sqrt[3]{a}$. The linearization of the FEA around the background state $\mathbf{x}^b$ for IWP values below $IWP_{TH}$ becomes:

$$ \delta \text{FEA} = \left( \frac{d\text{FEA}}{d\text{IWP}} \right)_{\mathbf{x}^b} \delta \text{IWP} = \sqrt[3]{a} \times \delta \text{IWP}, \qquad (6) $$

and the gradient is thus never null. This cubic-root transformation is also applied to the observations.

## 4. Case description and experiments design

This section first describes the two Mediterranean events studied in this paper for which the new assimilation scheme was applied. The second subsection then discusses the FEA accumulation time, the lightning observation errors and the thinning method applied to the observations. The four different configurations mentioned in the Introduction are detailed along with the processing
<figure><img src="image_4.png" /></figure>

applied to each of them. The experiments set up to assess the performances of LDA are described
at the end of this section.

<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>

a. *Studied events*

<figure><img src="image_7.png" /></figure>

The first case study is a heavy precipitation event in the south-east of France on 3-4 October 2021. Thunderstorms hit the city of Marseille during the night of the 3rd to the 4th, causing devastating floods, with measured cumulative precipitation beyond 174 mm overnight. A severe weather watch (*Vigilance rouge*) for heavy rain and flooding was issued by Météo-France for the Bouches-du-Rhône department on the morning of the 4th because another rainstorm was forecasted for the afternoon of the 4th. However, the operational AROME-France misplaced the precipitation and the Var department, to the east of Bouches-du-Rhône, was hit instead. The island of Corsica was also hit by severe thunderstorms on the afternoon of the 4th, leading to cumulative rainfall exceeding 200 mm in the eastern part of the island over the whole day. The observed precipitation accumulated over 24 h from the 3rd at 18:00 to the 4th at 18:00 UTC (20:00 local time) is plotted in Fig. 1b. These observations are generated by merging radar-derived rainfall and rain gauge observations with a horizontal grid spacing of 0.1° and are called the ANTILOPE product (ANalyse par spaTiaLisation hOraire des PrEcipitations, Laurantin 2013). Grayed areas in Fig. 1b are outside the ANTILOPE domain and do not contain precipitation data. All the locations mentioned above are indicated in italics in Fig. 1.

<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /></figure>

<figure><img src="image_13.png" /></figure>

<figure><img src="image_14.png" /></figure>

<figure><img src="image_15.png" /></figure>

<figure><img src="image_16.png" /></figure>

<figure><img src="image_17.png" /></figure>

<figure><img src="image_18.png" /></figure>

<figure><img src="image_19.png" /></figure>

<figure><img src="image_20.png" /></figure>

<figure><img src="image_21.png" /></figure>

<figure><img src="image_22.png" /></figure>

<figure><img src="image_23.png" /></figure>

<figure><img src="image_24.png" /></figure>

<figure><img src="image_25.png" /></figure>

<figure><img src="image_26.png" /></figure>

<figure><img src="image_27.png" /></figure>

<figure><img src="image_28.png" /></figure>

<figure><img src="image_29.png" /></figure>

<figure><img src="image_30.png" /></figure>

<figure><img src="image_31.png" /></figure>

The second studied event is the Mediterranean storm Blas, that hit the Balearic Islands (Mallorca, Menorca and Ibiza), Sardinia and Corsica on 5-17 November 2021. The port of Valencia was closed on 6 November because of storm-induced waves of 2.5-4 m high. Wind gusts up to 20 m s⁻¹ measured at the Menorca airport forced all the ports of the island to close on 13 November. Heavy rainfall was measured in Corsica, for instance at Ghisoni (middle Corsica) where a rain accumulation of 229.1 mm was recorded on 9-10 November. The Blas storm was chosen for this study because it had a daily lightning activity in a region out of range of the radar network used operationally for assimilation in AROME-France (Martet et al. 2022). Indeed, data from one radar in Mallorca and two radars in Corsica are assimilated, but their range is not sufficient to monitor the events in between these two islands.
## b. Experimental setup

To assess the impact of LDA on the analysis and forecast fields, results of a LDA experiment (Ref+LI) are compared against results from a reference experiment (Ref). Both are using the AROME-France model and are performed in a 3D-EnVar framework. These experiments run in a continuous data assimilation cycle where each 1-hr forecast serves as the background state for the next analysis. Observations are thus assimilated every hour. The first background state is provided by operational AROME-France for both experiments. The first assimilation step is performed at 01:00 UTC on 3 October and the last one at 23:00 UTC on 4 October for the first case study and 01:00 UTC on 5 November to 23:00 UTC on 17 November for the second case study. Long forecasts, up to 12 hours, are performed every 6 hours: those are the forecasts that will be evaluated in Section 5c. The background error covariances for both experiments are provided by the AROME EDA, composed of 50 members and running at a numerically affordable horizontal resolution of 3.2 km (versus 1.3 km for the deterministic model) and a 3 h DA cycle to reduce numerical cost. The EDA runs independently of the deterministic Ref and Ref+LI experiments and only provides the background errors covariances to them. The EDA generates 50 independent analyses from 50 backgrounds. Variability among the members of the EDA is achieved through variability in the background members (provided by the previous ensemble since it is cycled), variability of the boundary conditions and perturbations of the assimilated observations. The assimilated observations in the EDA background members are the same as in the Ref experiment (see next section) and assimilation is performed in a 3D-Var framework.

## c. Assimilated observations

The Ref experiment assimilates in 3D-EnVar all the types of observations that are assimilated operationally, meaning that their observation errors are tuned optimally and they all have undergone the verification processes to be used operationally. These observations include: conventional observations (mainly temperature and humidity) measured by weather stations and radiosondes, radial wind and reflectivities measured by radars, radiances from geostationary and polar-orbiting satellites in clear-sky conditions. Radar reflectivites are assimilated as pseudo-profiles of relative humidity retrieved through a Bayesian procedure (see Caumont et al. 2010; Wattrelot et al. 2014). Because radiances are assimilated in clear-sky conditions and radar reflectivities are transformed
into relative humidity, none of the assimilated observations in **Ref** have a direct impact on hydrometeor contents. They can however have an impact through cross-correlations contained in the flow-dependent background-error covariances matrix. The LDA experiment assimilates lightning in addition to the types of observations already assimilated in the **Ref** experiment, that is why it is labeled **Ref+LI**.

The FEA is accumulated during 10 min, starting 5 min prior to the analysis time. Even though a sensitivity study on the FEA accumulation period was performed in Combarnous et al. (2022), it did not allow the identification of an optimal accumulation time for the assimilation of the FEA. A 10 min accumulation period was then selected because this duration was already used in several LDA studies (e.g. Marchand and Fuelberg 2014; Fierro et al. 2016). Furthermore, a shorter period of lightning activity is not sufficient to capture the extent of a convective cloud, and an accumulation period larger than 10 min would introduce displacement errors in the convection when assimilated against a fixed analysis time and too large non-zero FEA footprints. Finally, it was shown in Hu et al. (2020) that accumulated GLM lightning activity over 10 min was the one yielding the best results in terms of precipitation and composite reflectivity when assimilated in 3D-Var.

The lightning observation error covariances, contained in the diagonal of the **R** matrix, are associated with the measurement and its representativeness. Therefore, it will depend on both the performances and the resolution of the LI instrument. The error we used in our study is equal to 0.1 fl¹/³ and was chosen arbitrarily for a proof-of-concept, since the satellite lightning data are synthetic. As satellite observation footprints may extend over many grid points of a high spatial resolution model, and as horizontal integration is not implemented in the AROME-France data assimilation system, the observation error should additionally consider this potential representativeness error. The tuning of the optimal lightning observation error will be the subject of future work and will be assessed using long duration experiments and sensitivity studies once the real MTG-LI data are available and will take the representativeness error into account.

Observation error cross-correlations are difficult to estimate and can create problems in the calculation of the analysis and quality control algorithms and thus it is assumed that observation errors are not correlated. In practice, to reduce cross-correlation in the lightning observations, a thinning method selecting only one out of two grid points is applied to the observations with a geometrical pattern shown in Fig. 3 for an example grid of 4 by 4 observation points. The FEA
<figure><img src="image_32.png" /><figcaption>FIG. 3. Thinning method.</figcaption></figure>

values are assimilated at the position of the pixel center of the FEA grid. Other thinning patterns are currently tested in long duration assimilation experiments to identify the one yielding the best scores, but their evaluation is beyond the scope of this study.

The observed FEA takes discrete values, where zero values correspond to no lightning detected, and the next possible value is $1 \text{ fl}^{1/3}$, meaning that 1 flash illuminated in a pixel in 10 min. Even though observed FEA values cannot range in $]0; 1[$, it is nonetheless possible for the background FEA values to be assigned such range values. However, we considered that a background value below $1 \text{ fl}^{1/3}$ corresponds to no lightning. This results in four possible configurations when the observations are compared to their model equivalent, summarized in Table 2. Listed below are the treatments applied for each of those configurations.

* **Configuration 1:** both model and observation do not contain lightning. Those data are irrelevant because their assimilation will not bring useful information and are then discarded.

* **Configuration 2:** the model background contains lightning but none is observed. To take into account the ambiguity between a true "no-lightning" observation or an "undetected" flash, those observations are assimilated with an observation error larger than the non-zero observations. This larger observation error also allows weighting the assimilation less heavily toward these zero observations and to avoid excessive drying in the model. We chose an observation error 10 times larger for this study as a proof of concept, but this value will be set with sensitivity tests in future studies with the real MTG-LI observations.

* **Configuration 3:** lightning is observed but none is modeled in the background. This problem has been addressed in section 3b with the introduction of a variable change resulting in a constant gradient for low/null lightning values.
<table><caption>TABLE 2. Studied configurations of observed FEA values versus their background model equivalent.</caption><thead><tr><th></th><th>obs = 0</th><th>obs &ge; 1</th></tr></thead><tbody><tr><td>background &lt; 1</td><td>Configuration 1 (discarded)</td><td>Configuration 3</td></tr><tr><td>background &ge; 1</td><td>Configuration 2</td><td>Configuration 4</td></tr></tbody></table>

* **Configuration 4:** both observation and model background contain lightning. The assimilation of those points does not raise particular issues and will result in a correction of the modeled FEA value toward the observed one.

Because FEA observations are dominated by zero values, configurations 1 and 2 are the most commonly encountered. Once the observations belonging to configuration 1 are discarded, 2/3 of the remaining observations belong to configuration 2 and 1/3 to configurations 3 and 4. This distribution varies slightly depending on analysis time.

To study the behavior of the assimilation system for each one of the configurations described above, we carry out single-observation experiments; that is, a single-observation point at a specific spatial coordinate is assimilated. No other type of observation else than lightning is assimilated. Only configurations 2 and 3 are examined because samples belonging to configuration 1 are discarded and the assimilation of a point data in configuration 4 does not raise particular issue. The observed FEA values that are assimilated in this context are not the values that were actually measured but they were chosen to correspond to the desired configuration. The analysis increments (analysis minus background) resulting from the assimilation of these two points are studied in Section 5a. The single-observation experiments were conducted independently of the **Ref** and **Ref+LI** experiments. The observations were assimilated only at 10:00 UTC on 4 October, with a background obtained from a 1 h forecast of the operational AROME-France. This analysis time in particular was arbitrarily chosen, and single-observation assimilation at any other analysis time gives similar results.

# 5. Results

## a. Single-observation increments

The positions of the two assimilated single-observations, labeled A and B, are displayed on the FEA model background field in Fig. 4a. Point A belongs to configuration 2 i.e. background FEA
<table><caption>TABLE 3. single-observations description.</caption><thead><tr><th></th><th>Lon (°E)</th><th>Lat (°N)</th><th>Obs (fl<sup>1/3</sup>)</th><th>Background (fl<sup>1/3</sup>)</th><th>Analysis (fl<sup>1/3</sup>)</th><th>Configuration</th></tr></thead><tbody><tr><td>A</td><td>5.76</td><td>40.81</td><td>0.00</td><td>4.21</td><td>3.93</td><td>Configuration 2</td></tr><tr><td>B</td><td>7.00</td><td>43.80</td><td>3.00</td><td>9.64 × 10<sup>-7</sup></td><td>2.98</td><td>Configuration 3</td></tr></tbody></table>

value higher than 1 fl¹/³ while point B belongs to configuration 3, i.e. background FEA value lower than 1 fl¹/³. The values of the observation, background and analysis at points A and B are summarized in Table 3, along with their geo-locations. The examination of the analysis increments of FEA plotted in Fig. 4b and 4c shows that the quantity of FEA decreased in the vicinity of point A and increased in the vicinity of point B. There is an order of magnitude difference between the increments of points A and B that results from the 10-fold observation error applied to zero observations. For point A, the largest absolute FEA increment value is found tens of kilometers to the north of the position of the single-observation, reaching a value of 0.57 fl¹/³. For point B, the largest increment value is closer to the single-observation position and reaches almost 3 fl¹/³. Thus, it has been possible to add microphysical contents where there was none in the background, since the background value in point B is ~ 10⁻⁷ fl¹/³. The geographical distribution of the increments is dictated by the background-error covariances contained in the $\tilde{B}$ matrix, which depends on the current 3D atmospheric conditions. The horizontal extension of the increment depends on the horizontal localization of the $\tilde{B}$ matrix, addressed in Section 2a. In the troposphere, considered roughly below 12 km in the French latitudes, the horizontal localization length ranges between 25 km (at 5 m above sea level, a.s.l.) to approximately 70 km (at 12 km a.s.l.). Hydrometeors being present in the troposphere, the size of the increments plotted in Fig. 4 is consistent with these horizontal localization lengths.

The FEA being calculated from an integrated quantity over the atmospheric column, we study the increments of snow and graupel specific contents to get a sense of the altitudes at which changes occur after the assimilation of a single-observation. Fig. 5 shows the increment cross sections along points A and B longitudes, at 5.76 ° E and 7 ° E, respectively. The vertical profiles of the specific contents of snow and graupel for the background and analysis at points A and B positions are also depicted in Fig. 6, although null for the background at position B since the objective is to add content where none was forecasted. Those cross section positions are represented by the
<figure><img src="image_33.png" /><figcaption>FIG. 4. FEA background field (a) with green and blue rectangles indicating the subdomains where the analysis increments are plotted for single-observation A (b) and single-observation B (c). The positions of the single-observations A and B are indicated by the black circles. Background and analysis are valid on 4 October 2021 at 10:00 UTC. The dashed black lines in panels (b) and (c) indicate the position of the transects plotted in Fig. 5. Note the order of magnitude difference in the color scale of the increments.</figcaption></figure>

dashed black lines in Fig. 4b and 4c. Along these cross sections, the positions of A and B are represented by black vertical lines in Fig. 5. For the specific content of snow, highest absolute increment values are found around 7000 m a.s.l., whereas it is slightly below 6000 m a.s.l. for the specific content of graupel, which corresponds to the altitudes where the maximum of each of these specific contents can be found in the background (Fig. 6). At the position of point A, the background has a specific content of snow and graupel reaching a maximum value of $1.6×10^{-3}$ kg/kg and $4.25×10^{-3}$ kg/kg, respectively, and those values decrease of roughly 13% and 7% at the
<figure><img src="image_34.png" /><figcaption>FIG. 5. Specific contents of snow (a,b) and graupel (c,d) analysis increments cross sections along point A longitude (a,c) and point B longitude (b,d) represented by the black dashed lines in Fig. 4b and 4c. Valid on 4 October 2021 at 10:00 UTC.</figcaption></figure>

analysis, respectively for snow and graupel. For both of those hydrometeor types, the altitudes where modifications of contents occur range between 4000 and 9000 m a.s.l. The specific content of graupel has a maximum absolute increment value higher than the specific content of snow, both at point A and point B, reaching $-3.2 \times 10^{-4}$ and $1.73 \times 10^{-3}$ kg/kg, respectively. Again, an order of magnitude difference can be observed between the increments of point A and that of point B, resulting from the observation error difference. The increments along point A appear to be stratified into more stable layers than along point B where they appear noisy, which is probably due to the orography.

Because of the cross-correlations among all the control variables in the background-error covariances $\tilde{B}$ matrix, LDA impacts all the variables and not only the ones included in the observation operator. Fig. 7 compares the cross sections of temperature and specific humidity increments along points A and B longitudes. Because the quantities of temperature and specific humidity increase in the vicinity of point B, one can deduce that they are positively correlated with specific contents of snow and graupel in this area. The maxima of increase are located around 6000 m a.s.l., which is roughly the same altitude as for snow and graupel and reach values in excess of 1.3 K and 2.4 kg/kg for temperature and specific humidity, respectively. In the vicinity of point A, the specific humidity is positively correlated with graupel and snow specific contents, with a minimum value
<figure><img src="image_35.png" /><figcaption>FIG. 6. Vertical profiles of the specific contents of snow (a,c) and graupel (b,d) at the position of points A (a,b) and B (c,d) for the background (solid red line) and the analysis (dotted red line). Valid on 4 October 2021 at 10:00 UTC.</figcaption></figure>

<figure><img src="image_36.png" /><figcaption>FIG. 7. Specific humidity (a,b) and temperature (c,d) analysis increment cross sections along point A longitude (a,c) and point B longitude (b,d) represented by the black dashed lines in Fig. 4b and 4c. Valid on 4 October 2021 at 10:00 UTC.</figcaption></figure>

of $-4.1 \times 10^{-4}$ kg/kg. However, the correlation with temperature seems noisier, with a positive correlation between 4000 and 8000 m high and between 40.6 and 40.95 °N, resulting in a decrease of temperature, and a negative correlation beyond.
<figure><img src="image_37.png" /><figcaption>FIG. 8. For each hour of the 3-4 October 2021: (a) histogram of the number of observations belonging to configuration 2 (blue bars) or to configurations 3 and 4 (orange bars), and (b) root-mean-square of the innovation (solid line) and analysis residual (dashed line) of the FEA.</figcaption></figure>

b. Impact of lightning data assimilation on analysis fields

To measure how successful the LDA is in the studied case, we examine the variation of the root-mean-square (rms) of the residuals between FEA innovations and FEA analysis for each hour of the 3-4 October 2021. This metric has already been used in several recent LDA studies (e.g. Fierro et al. 2019; Hu et al. 2020; Xiao et al. 2021) and it represents how close the analyses came to the observations after the minimization process. Fig. 8a shows the number of observations belonging either to configuration 2 (no lightning detected) or configurations 3 and 4 (lightning detected). Over the whole period, the rms decreases by 21% on average (Fig. 8b). However, most of the lightning activity occurred on 4 October 2021, between 9:00 and 19:00 UTC. During this specific sub-period, the rms decreases by 33.6%. The decrease in rms indicates a successful assimilation of lightning data and the difference of decrease between a lightning-active period compared to a no-lightning one is consistent with the 10 times larger observation error applied to the zero observations. With this larger observation error, the assimilation system is weighted less heavily toward the observations, resulting in an analysis residual rms relatively close to the innovation rms as seen for the first ten hours on 3 October 2021. Note that the rms decrease presents a similar amplitude for the case study of 5-17 November (not shown).
Next, we compare analysis increments from the **Ref** and **Ref+LI** experiments to evaluate the impact of LDA on the analysis fields. Because lightning is the only type of observation assimilated that directly modifies specific contents of snow and graupel, we expect a larger increment of these two contents when lightning data is assimilated. This is particularly expected in the Mediterranean Sea, where a lot of areas are out of range of the radar network. An example is thus given for the Blas storm, since it occurred in the Mediterranean Sea, in Fig. 9 for an analysis at 18:00 UTC on 9 November 2021. Because the experiments are cycled, the analysis obtained at each hour after the assimilation is used as background for the next hour and thus the backgrounds used to calculate the increments of **Ref** and **Ref+LI** plotted in Fig. 9 are not the same. The IWP increment, calculated from the specific contents of snow and graupel as given in Eq. 3, is first examined because it is directly linked to the FEA via the observation operator. Fig. 9b shows the IWP increment for the **Ref** experiment: the only modifications in IWP are around south Corsica and north of the Balearic Islands, likely due to the assimilation of relative humidity derived from radars in Mallorca and Corsica. Fig. 9c shows the IWP increment for the **Ref+LI** experiment, which is in good agreement with the FEA innovation (observation minus background) plotted in Fig. 9a. Indeed, because the FEA is a positive function of the IWP in the observation operator, a positive FEA innovation means the model background lacks IWP, as opposed to a negative innovation which means the model background contains too much IWP. The areas where the IWP increased or decreased after the assimilation are consistent with the innovation sign at the corresponding position, demonstrating again a successful assimilation. The examination of the analysis increment at other hours and days reveals the same characteristics: a low IWP increment for **Ref** compared to **Ref+LI** because no other observation type else than lightning directly impact the graupel and snow specific contents, and a successful assimilation of the FEA. This demonstrates the potential of LDA in the Mediterranean Sea, by producing increments of hydrometeors where radar assimilation cannot.

As stated earlier, LDA does not only impact the variables that appear in the observation operator, i.e. graupel and snow, but also all the other control variables via the cross-correlations described by $\tilde{B}$. Thus, Fig. 9d and 9e compares the increment of integrated water vapor (IWV, in kg m$^{-2}$),
again at 18:00 UTC on 9 November 2021, calculated over the whole atmospheric column as follows

$$ \text{IWV} = \sum \rho \times q_v \times \Delta z, \qquad (7) $$

where $q_v$ is the specific humidity in kg kg⁻¹ and the other quantities were already introduced in Eq. 3. Fig. 9e shows that LDA increased the quantity of IWV where lightning was observed, in between the islands of Menorca and Sardinia. The IWV increment over the rest of the domain presents little to no difference with the **Ref** experiment as it comes from the assimilation of all other types of observations.

### c. Impact of lightning data assimilation on forecasting skills

In this subsection, we compare the results of forecasts from the **Ref+LI** experiment to the **Ref** experiment for the 3-4 October 2021 and the Blas storm on 5-17 November 2021. Accumulated rainfall (mm) and brightness temperature (K) forecasts are examined to evaluate the performance of LDA. The accumulative rainfall is compared to the ANTILOPE product, already mentioned in Section 4a. The RTTOV radiative transfer model (Saunders et al. 2018) is applied to forecasted fields to simulate 10.8 µm brightness temperature to be compared with Meteosat-11 SEVIRI product at a 3 km resolution. To evaluate the forecasts of the LDA experiments, we used the fraction skill score (FSS; Roberts and Lean 2008) and the frequency bias. For both precipitation and brightness temperature forecasts, the scores are calculated using the model point closest to the observation. The mean scores are calculated using forecasts produced every 6 hours. The frequency bias is the fraction of the number of forecasted events over the number of observed events:

$$ \frac{a+b}{a+c} \qquad (8) $$

where a, b and c are described in the contingency table (Tab. 4). An unbiased forecast has a frequency bias equal to 1. The FSS is a score that compares the fraction of forecasted and observed events within a neighborhood, behaving as a sliding spatial window. It ranges from 0 to 1, with a value of 1 indicating a perfect score. The FSS is typically calculated for different window sizes and different thresholds. Here, the FSS is given in terms of relative improvements compared with
obs-background - 2021-11-09 18:00 UTC

<figure><img src="image_38.png" /></figure>

IWP analysis increment - 2021-11-09 18:00 UTC

<figure><img src="image_39.png" /></figure>

IWV analysis increment - 2021-11-09 18:00 UTC

<figure><img src="image_40.png" /><figcaption>FIG. 9. Comparison of the analysis increments at 18:00 UTC on 9 November 2021 (Blas storm in the Mediterranean Sea). (a) innovation (observation minus background) at position of observation points. Ice water path increment for (b) **Ref** and (c) **Ref+LI** experiments. Integrated water vapor increment for (d) **Ref** and (e) **Ref+LI** experiments.</figcaption></figure>

the **Ref** experiment:

$$ \Delta \text{FSS} = \frac{\text{FSS(Ref+LI)} - \text{FSS(Ref)}}{\text{FSS(Ref)}}. \qquad (9) $$
<table><caption>TABLE 4. Contingency table.</caption><thead><tr><th rowspan="2">Forecast</th><th colspan="2">Observation</th></tr><tr><th>Yes</th><th>No</th></tr></thead><tbody><tr><td>Yes</td><td>Hit: <i>a</i></td><td>False alarm: <i>b</i></td></tr><tr><td>No</td><td>Miss: <i>c</i></td><td>Correct negative: <i>d</i></td></tr></tbody></table>

Statistical significance of score differences between **Ref** and **Ref+LI** is assessed with a 95% confidence interval bootstrap test (Efron and Tibshirani 1993).

Both FSS and frequency bias are calculated in the south-east (SE) and BLAS domains shown in Fig. 1a for the 3-4 October and 5-17 November, respectively. However, the ANTILOPE domain southern limit is 41°N, meaning it only partially covers the BLAS domain. Additionally, storm-induced precipitation mainly occurred over the sea during the Blas storm, but the ANTILOPE product is not reliable over the sea because it is too far from the radar network and a filter is applied to select only data over land. Consequently, the remaining available precipitation observations for the Blas storm do not cover a region large enough and verification scores are too noisy to be studied.

We first present overall scores for both case studies, and then some examples to qualitatively highlight specific strengths and limitations of LDA.

## 1) OVERALL SCORES

### (i) Precipitation
The impact of LDA is first examined on precipitation forecasts, only for the 3-4 October 2021 case. Fig. 10 compares the mean frequency bias calculated for different accumulated precipitation thresholds over periods of 1 h and plotted as a function of forecast lead time. The different precipitation thresholds studied are 0.1 mm (occurrence of precipitation, Fig. 10a), 0.5 mm (low intensity precipitation, Fig. 10b), 2 mm (medium intensity precipitation, Fig. 10c) and 5.0 mm (heavy precipitation, Fig. 10d). Because forecasts are produced every 6 h (at base hours 06:00, 12:00, 18:00 and 00:00 UTC), the mean bias for the 3-4 October case plotted in Fig. 10 is calculated only with 7 values and results must be interpreted with care. Generally, the frequency bias decreases with forecast lead time and, for all the thresholds, the frequency bias is below 1, ranging between 0.75 and 0.95 for the 0.1 mm, 0.5 mm and 2.0 mm thresholds (Fig. 10a-c) and between 0.65 and 1 for the 0.5 mm threshold (Fig. 10d). It means that precipitation is
generally underestimated for both **Ref** and **Ref+LI** experiments. Nevertheless, LDA increased the frequency bias for the 0.1 mm, 0.5 mm and 2 mm thresholds, bringing it closer to unity for forecasts up to 3 h after the assimilation. The increase in bias when lightning data are assimilated compared to **Ref** is however quite small, at a maximum of 0.05 points (0-1 h accumulated precipitation, threshold 0.1 mm, Fig. 10a). After 3 h, the trend reverses, and the bias of the **Ref+LI** experiment becomes smaller than that of the **Ref** experiment. Conversely, for the higher threshold, 5 mm per hour (Fig. 10d), the bias is deteriorated when lightning data are added, especially for the 2-3 h and 3-4 h forecasts, decreasing by roughly 0.1 points. Overall, these scores indicate that LDA increases the general precipitation area (when precipitation is higher than 0.1 mm per hour), which is beneficial for the 3-4 October case since it generally underestimates the precipitation. The small LDA-induced increase in precipitation is not maintained after 3 h of forecast lead time. The differences in FSS for both experiments for the accumulated precipitation for the 3-4 October case are not statistically significant for the majority of the forecast hours and therefore are not shown here.

(ii) *Cloud-cover* Modifying the hydrometeor contents has a direct impact on the cloud development and structure, and consequently the altitude of the top of the cloud, which is an indicator of the storm severity. Consequently, we compare the performances of **Ref** and **Ref+LI** in reproducing the brightness temperatures observed by SEVIRI at 10.8 µm by calculating frequency bias and FSS for both experiments and both case studies. We first examine the mean frequency bias of brightness temperatures lower than four thresholds, 280, 260, 240 and 220 K at forecast times of 1, 2 and 3 h in Fig. 11a for the November case and in Fig. 11b for the October case. Again, mean values are calculated with 7 values for the October case because forecasts are performed every 6 h and with 51 values for the November case. For both case studies the bias is closer to unity for **Ref+LI** for the three higher temperature thresholds: 280, 260 and 240 K and at all forecast times. For the lower temperature threshold, 220 K, representing the most convective areas, a sudden drop can be observed for forecast times higher than 1 h. For instance, for the 1 h forecasts of 5-17 November (Fig. 11a), the mean bias at 220 K is 0.92 for **Ref+LI** versus 0.44 for **Ref**. In contrast, it drops to 0.31 and 0.26 for the 2 h and 3 h forecasts, respectively, which is even lower than the bias for **Ref**. For the 3-4 October (Fig. 11b) a similar decrease is observed for **Ref+LI** with a mean bias that varies from 0.89 for the 1 h forecast to 0.57 and 0.48 for the 2 h and 3 h forecasts at
<figure><img src="image_41.png" /><figcaption>FIG. 10. Mean frequency bias of accumulated precipitation over periods of 1 h as a function of forecast lead time. The precipitation thresholds are (a) 0.1 mm per hour, (b) 0.5 mm per hour, (c) 2 mm per hour and (d) 5 mm per hour. Results for **Ref** are plotted in black and for **Ref+LI** in red.</figcaption></figure>

220 K. These scores suggest that, overall, the cloud cover is better represented when lightning is assimilated because the frequency bias of the warmer brightness temperatures (280, 260, 240 K) are closer to unity. Nevertheless, for both the October and November cases, the bias is below 1 for all temperature thresholds, implying that the brightness temperatures are generally underestimated both for **Ref** and **Ref+LI**. In the present cases, LDA greatly improves the forecast of the most convective area of the clouds for the first forecast hour, represented by the temperatures below 220 K. However, these colder cloud top areas are not maintained through time and dissipate faster than in the **Ref** experiment and in the SEVIRI observations.
<figure><img src="image_42.png" /></figure>

FIG. 11. Mean frequency bias of brightness temperature forecasts as a function of 4 different temperature thresholds: 280, 260, 240 and 220 K, for (a) 3-4 October and (b) 5-17 November 2021. Results for **Ref+LI** are plotted in red and for **Ref** in black. Solid, dashed and dotted lines represent 1 h, 2 h and 3 h forecasts, respectively.

The FSS relative improvement (ΔFSS) of **Ref+LI** compared to **Ref** is shown in Fig. 12b and 12d for temperature thresholds ranging from 220-290 K in steps of 10 K and forecast lead times up to 4 h after the assimilation. In Fig. 12a and 12c, the FSS values for the **Ref** experiments are given for information only. A positive (negative) ΔFSS indicates that **Ref+LI** (**Ref**) has a higher FSS. The FSS shown in Fig. 12 are all calculated for a window size of 0.5°, which is roughly the same as in Destouches et al. (2023). Window sizes of 0.05, 0.1, 0.25° were also tested, leading to very similar FSS values (not shown). Statistical significance of ΔFSS is indicated with "+" and "-" in Fig. 12b and 12d. Overall, **Ref+LI** is more skillful than **Ref**: a larger amount of positive ΔFSS are obtained, but significant improvements, up to 3%, are mostly observed for the higher brightness temperatures (higher than 270 K for the November case and higher than 250 K for the October case) and first hour of the forecasts. The larger improvements, higher than 5%, are found for colder cloud tops (lower than 240 K), although the differences are not statistically significant. Negative impact is observed for the lower brightness temperatures after 2-3 h of forecast, especially for the November case.

This FSS improvement at 270-290 K demonstrates that, for the studied cases, DA experiments additionally assimilating lightning data better capture the cloud cover and this improvement is maintained within the first 4 forecast hours. The benefit of lightning data for predicting lower temperature occurrences, defining the most convective regions of the cloud, is more mixed. It
<figure><img src="image_43.png" /></figure>

FIG. 13. (a) observed and (b,c) predicted brightness temperature at 10.8 µm valid at 13:00 UTC on 4 October 2021. Results from the **Ref** experiments and the **Ref+LI** experiment are shown in (b) and (c), respectively, and are obtained from a 1 h forecast. The transition from gray shades to rainbow shades indicates the temperatures corresponding to convective areas.

the SEVIRI instrument at 13:00 UTC (Fig. 13a) to 1 h forecasts valid at 13:00 UTC on 4 October 2021 for the **Ref** (Fig. 13b) and **Ref+LI** (Fig. 13c) experiments. Overall, both experiments are able to coarsely reproduce the large system over south-eastern France, between 4 and 7° E and 43 and 46° N, even though the lowest temperatures are overestimated in the **Ref+LI** experiment (Fig. 13c) compared to observations (Fig. 13a). **Ref** failed to initiate the convection that started in southern Corsica whereas **Ref+LI** reproduced it with an areal coverage and temperature values similar to the observations. The spatial extent of the convective cells forecasted in northern Corsica in **Ref** is less than what is observed, suggesting that the forecast might be delayed compared to the observations. This delay is corrected in **Ref+LI**, which successfully captures the northward extension of the convective system above northern Corsica. The system over northern Italy, roughly at 9° E and 45° N exhibits lower brightness temperatures, below 210 K, when lightning is assimilated, which is closer to the observed values than what the **Ref** forecasted. This example shows that LDA can produce a burst of convection during the first hour of forecast by enlarging the areas with the coldest brightness temperatures.

This "burst" however does not allow the convective systems to last as long as the observed one and the convective systems tend to decay faster when lightning data are assimilated. To illustrate this fading of the lowest brightness temperatures areas, an example is displayed in Fig. 14, comparing results from forecast times of 1, 2 and 3 h, for the November case. Among the two cells observed between Corsica and the Balearic Islands, only the northern one is present in the 1 h forecast
<figure><img src="image_44.png" /><figcaption>FIG. 14. (a)-(c) Observed and (d)-(i) forecasted brightness temperatures valid at (left to right) 13:00, 14:00 and 15:00 UTC on 7 November 2021. Data were assimilated at 12:00 the same day. Results from **Ref** and **Ref+LI** are shown in (d)-(f) and (g)-(i), respectively.</figcaption></figure>

for both **Ref** and **Ref+LI**, although its extent is better represented when lightning is assimilated. However, at forecast times of 2 h and 3 h, **Ref** improves and successfully captures the two-cell structure, whereas the colder area in **Ref+LI** quickly fades and almost no cloud tops below 235 K are forecasted 3 h after LDA.

# 6. Conclusions

We introduced in this study a method to perform direct assimilation of lightning data in a 3D-EnVar DA scheme using an observation operator based on the specific contents of snow and graupel. This work was performed in the framework of the development of a new 3D-EnVar system for the regional French NWP model AROME-France, which will also see the addition of hydrometeors in the control variable. Data assimilation experiments were conducted using lightning observations mimicking the future MTG-LI data, since the long-term objective of this
study is the assimilation of space-borne lightning observations, which present advantages in terms
of geographical coverage and detection efficiency compared to ground-based LF/VLF lightning
detection networks. The design of the experiments was very close to the operational configuration,
including all the observations that are currently assimilated in AROME-France **(Ref)** and results
were compared to experiments assimilating lightning in addition to them **(Ref+LI)**. In line with
Combarnous et al. (2022)'s conclusions, the observation operator is a cubic polynomial regression
between FEA observations and IWP calculated from AROME-France backgrounds of 44 days of
the year 2018. Modifications have been made to the operator to enhance the assimilation: (i)
the polynomial is replaced by a power function for small values of IWP and (ii) a cubic-root
transformation is applied to the FEA to solve convergence issues during the minimization process.
An observation error variance of 0.1 fl¹/³ was selected, which is smaller than the diagnosed
value from the consistency diagnostic of Desroziers et al. (2005) to highlight the effects of LDA,
and a thinning method keeping one out of two observation grid boxes was applied to reduce
cross-correlations in the observations. To compensate for the effect of the large number of zero
observations and to take into account the ambiguity those zero observations represent (true no-
lightning/undetected), a 10 times larger error is associated to them.

First, single-observation experiments were conducted to demonstrate the ability of the assimilation system to produce an increment of graupel and snow where none was forecasted in the background but where lightning was observed. This issue was mentioned in several studies (e.g. Errico et al. 2007; Lopez 2011) and a variable change allowed here to obtain a constant linear tangent term for low/null IWP value and get rid of a "zero gradient" when there is no IWP content. The study of the analysis increments demonstrated that the assimilation of a single-observation has an impact approximately between 3000 and 9000 m high and up to 30 km in the vicinity of the observation point, depending on the field (temperature, graupel, snow or humidity). The increment absolute value is roughly 10 times smaller when the observation is zero, which is coherent with the error difference applied to these observations. The impact of LDA on some other prognostic variables (temperature and specific humidity) than the ones in the observation operator is also discussed to show how they are updated through the ensemble background-error covariances. For the studied examples, the correlations between graupel, snow, temperature and specific humidity were mostly positive, meaning that the signs of their increments were the same. The examples shown
demonstrated that the LDA method developed here is both capable of reducing and increasing hydrometeor contents in the analysis fields.

Secondly, the examination of the analysis fields of **Ref+LI** demonstrates a successful assimilation of lightning data: a decrease in rms between the innovations and analysis residuals indicates that the values of the FEA in the analyses got closer to the observed values and the increments of IWP brought by the assimilation of FEA is consistent (sign and position) with the FEA innovation. Furthermore, when compared to IWP increments from the **Ref** experiment, it is shown that lightning is the assimilated observation that has the greater impact on snow and graupel contents, essential for thunderstorms development. For now, it is the only observation in the 3D-EnVar DA system of AROME-France that directly modifies the graupel and snow contents via its observation operator. This demonstrates the unequivocal potential of LDA to improve thunderstorm forecast.

Thirdly, precipitation and brightness temperature forecasts were studied for two cases of Mediterranean events, and the differences between **Ref** and **Ref+LI** were quantified using scores from the contingency table and the FSS. For the precipitation, the case of 3-4 October 2021 showed an improvement in frequency bias for precipitation thresholds lower than 2 mm per hour for forecasts up to 3 h. However, this improvement remains fairly small, at a maximum of 0.05 points. A larger number of precipitation case studies is obviously necessary to draw conclusions regarding the impact of LDA on rainfall accumulation forecasts, even though the results of 3-4 October 2021 are promising.

As for brightness temperatures, LDA improved the simulated fields and especially the coldest cloud top areas (<220 K) in the first forecast hour for the two cases studied. Nevertheless, the effects of LDA on cold cloud tops quickly fade with time. The description of the general cloud cover, defined by warmer brightness temperatures (<280 K) is better captured when lightning is assimilated, and this lasts until forecast times up to 4 hours after DA. The quick fading of coldest temperature areas through time needs further investigations, but could be explained by the fact that the assimilation of lightning mostly brings modifications to the cell hydrometeor contents but not enough to remote or larger-scale features like the low-level jet that may drive convection. Then, the convective systems lack energy from a warm and moist air and quickly dissipate. The assimilation system could benefit from a synergy between lightning and another type of observation that carries this kind of information to better forecast and maintain convective cells, like radar data and visible
and IR imagery that would point the presence of clouds. The very short-term improvement (3-4 h after DA) is consistent with other storm-scale LDA studies that were mentioned in the introduction (e.g. Hu et al. 2020; Kong et al. 2020; Xiao et al. 2021), but also with studies dealing with direct radar DA in convective-scale NWP (e.g. Bick et al. 2016) or simply when hydrometeors are included in the control variable (Destouches et al. 2023). This leads us to consider the assimilation of lightning data in a nowcasting model.

Finally, several parameters will have to be (re)calibrated when the real MTG-LI observations are available and this will be the subject of future works: the lightning observation error needs to be tuned. An optimal thinning method also needs to be identified and the coefficients of the regression function used in the lightning observation operator should be re-calculated. Although a solution has been designed in this paper to the "zero gradient" problem, which occurs when the background does not contain any hydrometeor content, an impediment remains when no members of the ensemble data assimilation system used to estimate $\tilde{\mathbf{B}}$ matrix covariances forecasted hydrometeors at the location where lightning is observed: the ensemble spread is then zero and convection initiation through assimilation is impossible. Several solutions are currently investigated to alleviate this issue.
Acknowledgments. We thank Météo-France and the Région Occitanie for funding the PhD of Pauline Combarnous. We also wish to thank Isabelle Couasnon and Tom Nicolau for providing the Météorage data. The script to calculate the verification scores is a collaborative work from the Groupe de Modélisation et d'Assimilation pour la Prévision of CNRM and we thank the authors for sharing the codes. This research has been supported by the CNES (Centre National d'Études Spatiales; SOLID project).

*Data availability statement.* The numerical model AROME and the data assimilation system are being developed at Météo-France. The code sources are not available under open source license and cannot be shared. Météorage data are provided by and are the property of Météorage as a company.

## References

Allen, B. J., E. R. Mansell, D. C. Dowell, and W. Deierling, 2016: Assimilation of pseudo-GLM data using the ensemble Kalman filter. *Monthly Weather Review*, **144** (9), 3465–3486, https://doi.org/10.1175/MWR-D-16-0117.1.

Benjamin, S. G., J. M. Brown, G. Brunet, P. Lynch, K. Saito, and T. W. Schlatter, 2019: 100 years of progress in forecasting and nwp applications. *Meteorological Monographs*, **59**, 13–1, https://doi.org/10.1175/AMSMONOGRAPHS-D-18-0020.1.

Bick, T., and Coauthors, 2016: Assimilation of 3d radar reflectivities with an ensemble kalman filter on the convective scale. *Quarterly Journal of the Royal Meteorological Society*, **142** (696), 1490–1504, https://doi.org/10.1002/qj.2751.

Brousseau, P., L. Berre, F. Bouttier, and G. Desroziers, 2011: Background-error covariances for a convective-scale data-assimilation system: AROME-France 3D-Var. *Quarterly Journal of the Royal Meteorological Society*, **137** (655), 409–422, https://doi.org/10.1002/qj.750.

Brousseau, P., Y. Seity, D. Ricard, and J. Léger, 2016: Improvement of the forecast of convective activity from the AROME-France system. *Quarterly Journal of the Royal Meteorological Society*, **142** (699), 2231–2243, https://doi.org/10.1002/qj.2822.
750 Brousseau, P., V. Vogt, E. Arbogast, and M. Martet, to be submitted: Configuration and evaluation
751 of a 3denvar data assimilation scheme for the météo-france operational convective model arome-
752 france. *Quarterly Journal of the Royal Meteorological Society*.

753 Cao, D., F. Lu, X. Zhang, and J. Yang, 2021: Lightning activity observed by the fengyun-4a
754 lightning mapping imager. *Remote Sensing, 13 (15)*, 3013, https://doi.org/10.3390/rs13153013.

755 Caron, J.-F., Y. Michel, T. Montmerle, and É. Arbogast, 2019: Improving background error
756 covariances in a 3d ensemble-variational data assimilation system for regional nwp. *Monthly
757 Weather Review, 147 (1)*, 135-151, https://doi.org/10.1175/MWR-D-18-0248.1.

758 Caumont, O., V. Ducrocq, É. Wattrelot, G. Jaubert, and S. Pradier-Vabre, 2010: 1d+3dvar as-
759 simulation of radar reflectivity data: a proof of concept. *Tellus A: Dynamic Meteorology and
760 Oceanography, 62 (2)*, 173-187, https://doi.org/10.1111/j.1600-0870.2009.00430.x.

761 Chen, Z., X. Qie, D. Liu, and Y. Xiong, 2019: Lightning data assimilation with comprehensively
762 nudging water contents at cloud-resolving scale using WRF model. *Atmospheric Research, 221*,
763 72-87, https://doi.org/10.1016/j.atmosres.2019.02.001.

764 Combarnous, P., F. Erdmann, O. Caumont, É. Defer, and M. Martet, 2022: A satellite lightning
765 observation operator for storm-scale numerical weather prediction. *Natural Hazards and Earth
766 System Sciences, 22 (9)*, 2943-2962, https://doi.org/10.5194/nhess-22-2943-2022.

767 Deierling, W., J. Latham, W. A. Petersen, S. M. Ellis, and H. J. Christian, 2005: On the relationship
768 of thunderstorm ice hydrometeor characteristics and total lightning measurements. *Atmospheric
769 Research, 76 (1-4)*, 114-126, https://doi.org/10.1016/j.atmosres.2004.11.023.

770 Deierling, W., W. A. Petersen, J. Latham, S. Ellis, and H. J. Christian, 2008: The relationship
771 between lightning activity and ice fluxes in thunderstorms. *Journal of Geophysical Research*,
772 **113** (*D15*), D15 210, https://doi.org/10.1029/2007JD009700.

773 Desroziers, G., L. Berre, B. Chapnik, and P. Poli, 2005: Diagnosis of observation, background
774 and analysis-error statistics in observation space. *Quarterly Journal of the Royal Meteorological
775 Society, 131 (613)*, 3385-3396, https://doi.org/10.1256/qj.05.108.
776 Destouches, M., T. Montmerle, Y. Michel, and J. Caron, 2023: Impact of hydrometeor control variables in a convective scale 3DEnVar data assimilation scheme. *Quarterly Journal of the Royal Meteorological Society*, https://doi.org/10.1002/qj.4426.

779 Efron, B., and R. J. Tibshirani, 1993: *An Introduction to the Bootstrap*. Chapman & Hall, 436 pp.

780 Erdmann, F., 2020: Préparation à l'utilisation des observations de l'imageur d'éclairs de météosat troisième génération pour la prévision numérique à courte échéance (preparation for the use of meteosat third generation lightning imager observations in short-term numerical weather prediction). Ph.D. thesis, Université Toulouse 3 - Paul Sabatier, 279 pp.

784 Erdmann, F., O. Caumont, and E. Defer, 2022: A Geostationary Lightning Pseudo-Observation Generator Utilizing Low-Frequency Ground-Based Lightning Observations. *Journal of Atmospheric and Oceanic Technology*, **39** (1), 3–30, https://doi.org/10.1175/JTECH-D-20-0160.1.

787 Erdmann, F., O. Caumont, and E. Defer, 2023: Assimilation of meteosat third generation (mtg) lightning imager (li) pseudo-observations in arome-france-proof of concept. *Natural Hazards and Earth System Sciences*, **23** (8), 2821–2840, https://doi.org/10.5194/nhess-23-2821-2023.

790 Errico, R. M., P. Bauer, and J.-F. Mahfouf, 2007: Issues Regarding the Assimilation of Cloud and Precipitation Data. *Journal of the Atmospheric Sciences*, **64** (11), 3785–3798, https://doi.org/10.1175/2006JAS2044.1.

793 Federico, S., R. C. Torcasio, E. Avolio, O. Caumont, M. Montopoli, L. Baldini, G. Vulpiani, and S. Dietrich, 2019: The impact of lightning and radar reflectivity factor data assimilation on the very short-term rainfall forecasts of RAMS@ISAC: application to two case studies in Italy. *Natural Hazards and Earth System Sciences*, **19** (8), 1839–1864, https://doi.org/10.5194/nhess-19-1839-2019.

798 Fierro, A. O., J. Gao, C. L. Ziegler, K. M. Calhoun, E. R. Mansell, and D. R. MacGorman, 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. *Monthly Weather Review*, **144** (11), 4373–4393, https://doi.org/10.1175/MWR-D-16-0053.1.
803 Fierro, A. O., E. R. Mansell, C. L. Ziegler, and D. R. MacGorman, 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. *Monthly Weather Review*, **140** (8), 2609–2627, https://doi.org/10.1175/MWR-D-11-00299.1.

803
804
805
806

Fierro, A. O., E. R. Mansell, C. L. Ziegler, and D. R. MacGorman, 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. *Monthly Weather Review*, **140** (8), 2609–2627, https://doi.org/10.1175/MWR-D-11-00299.1.

807 Fierro, A. O., Y. Wang, J. Gao, and E. R. Mansell, 2019: Variational Assimilation of Radar Data and GLM Lightning-Derived Water Vapor for the Short-Term Forecasts of High-Impact Convective Events. *Monthly Weather Review*, **147** (11), 4045–4069, https://doi.org/10.1175/MWR-D-18-0421.1.

807
808
809
810

Fierro, A. O., Y. Wang, J. Gao, and E. R. Mansell, 2019: Variational Assimilation of Radar Data and GLM Lightning-Derived Water Vapor for the Short-Term Forecasts of High-Impact Convective Events. *Monthly Weather Review*, **147** (11), 4045–4069, https://doi.org/10.1175/MWR-D-18-0421.1.

811 Goodman, S. J., D. E. Buechler, P. D. Wright, and W. D. Rust, 1988: Lightning and precipitation history of a microburst-producing storm. *Geophysical Research Letters*, **15** (11), 1185–1188, https://doi.org/10.1029/GL015i011p01185.

811
812
813

Goodman, S. J., D. E. Buechler, P. D. Wright, and W. D. Rust, 1988: Lightning and precipitation history of a microburst-producing storm. *Geophysical Research Letters*, **15** (11), 1185–1188, https://doi.org/10.1029/GL015i011p01185.

814 Goodman, S. J., and Coauthors, 2013: The goes-r geostationary lightning mapper (glm). *Atmospheric research*, **125**, 34–49, https://doi.org/10.1016/j.atmosres.2013.01.006.

814
815

Goodman, S. J., and Coauthors, 2013: The goes-r geostationary lightning mapper (glm). *Atmospheric research*, **125**, 34–49, https://doi.org/10.1016/j.atmosres.2013.01.006.

816 Hawkness-Smith, L. D., and D. Simonin, 2021: Radar reflectivity assimilation using hourly cycling 4D-Var in the Met Office Unified Model. *Quarterly Journal of the Royal Meteorological Society*, **147** (736), 1516–1538, https://doi.org/10.1002/qj.3977.

816
817
818

Hawkness-Smith, L. D., and D. Simonin, 2021: Radar reflectivity assimilation using hourly cycling 4D-Var in the Met Office Unified Model. *Quarterly Journal of the Royal Meteorological Society*, **147** (736), 1516–1538, https://doi.org/10.1002/qj.3977.

819 Hu, J., A. O. Fierro, Y. Wang, J. Gao, and E. R. Mansell, 2020: Exploring the Assimilation of GLM-Derived Water Vapor Mass in a Cycled 3DVAR Framework for the Short-Term Forecasts of High-Impact Convective Events. *Monthly Weather Review*, **148** (3), 1005–1028, https://doi.org/10.1175/MWR-D-19-0198.1.

819
820
821
822

Hu, J., A. O. Fierro, Y. Wang, J. Gao, and E. R. Mansell, 2020: Exploring the Assimilation of GLM-Derived Water Vapor Mass in a Cycled 3DVAR Framework for the Short-Term Forecasts of High-Impact Convective Events. *Monthly Weather Review*, **148** (3), 1005–1028, https://doi.org/10.1175/MWR-D-19-0198.1.

823 Kokou, P., P. Willemse, M. Lekouara, M. Arioua, A. Mora, P. Van den Braembussche, E. Neri, and D. M. A. Aminou, 2018: Algorithmic Chain for Lightning Detection and False Event Filtering Based on the MTG Lightning Imager. *IEEE Transactions on Geoscience and Remote Sensing*, **56** (9), 5115–5124, https://doi.org/10.1029/2007JD009598.

823
824
825
826

Kokou, P., P. Willemse, M. Lekouara, M. Arioua, A. Mora, P. Van den Braembussche, E. Neri, and D. M. A. Aminou, 2018: Algorithmic Chain for Lightning Detection and False Event Filtering Based on the MTG Lightning Imager. *IEEE Transactions on Geoscience and Remote Sensing*, **56** (9), 5115–5124, https://doi.org/10.1029/2007JD009598.

827 Kong, R., M. Xue, A. O. Fierro, Y. Jung, C. Liu, E. R. Mansell, and D. R. MacGorman, 2020: Assimilation of GOES-R Geostationary Lightning Mapper Flash Extent Density Data in GSI EnKF for the Analysis and Short-Term Forecast of a Mesoscale Convective System. *Monthly Weather Review*, **148** (5), 2111–2133, https://doi.org/10.1175/MWR-D-19-0192.1.

827
828
829
830

Kong, R., M. Xue, A. O. Fierro, Y. Jung, C. Liu, E. R. Mansell, and D. R. MacGorman, 2020: Assimilation of GOES-R Geostationary Lightning Mapper Flash Extent Density Data in GSI EnKF for the Analysis and Short-Term Forecast of a Mesoscale Convective System. *Monthly Weather Review*, **148** (5), 2111–2133, https://doi.org/10.1175/MWR-D-19-0192.1.
831 Kong, R., M. Xue, C. Liu, A. O. Fierro, and E. R. Mansell, 2022: Development of New Observation Operators for Assimilating GOES-R Geostationary Lightning Mapper Flash Extent Density Data Using GSI EnKF: Tests with Two Convective Events over the United States. *Monthly Weather Review*, **150 (8)**, 2091–2110, https://doi.org/10.1175/MWR-D-21-0326.1.

835 Laurantin, O., 2013: Antilope: hourly rainfall analysis over france merging radar and rain gauges data. *Proceedings of the 11th International Precipitation Conference*, H. Leijnse, and R. Uijlenhoet, Eds., KNMI, Ede-Wageningen, The Netherlands.

838 Lopez, P., 2011: Direct 4D-Var Assimilation of NCEP Stage IV Radar and Gauge Precipitation Data at ECMWF. *Monthly Weather Review*, **139 (7)**, 2098–2116, https://doi.org/10.1175/2010MWR3565.1.

841 MacGorman, D., and W. Rust, 1998: *The Electrical Nature of Storms*. Oxford University Press.

842 Mansell, E. R., 2014: Storm-scale ensemble Kalman filter assimilation of total lightning flash-extent data. *Monthly Weather Review*, **142 (10)**, 3683–3695, https://doi.org/10.1175/MWR-D-14-00061.1.

845 Marchand, M. R., and H. E. Fuelberg, 2014: Assimilation of Lightning Data Using a Nudging Method Involving Low-Level Warming. *Monthly Weather Review*, **142 (12)**, 4850–4871, https://doi.org/10.1175/MWR-D-14-00076.1.

848 Martet, M., P. Brousseau, E. Wattrelot, F. Guillaume, and J.-F. Mahfouf, 2022: Operational Assimilation of Radar Data from the European EUMETNET Programme OPERA in the Météo-France Convective-Scale Model AROME. *Data Assimilation for Atmospheric, Oceanic and Hydrologic Applications (Vol. IV)*, Springer International Publishing, Cham, 629–644, https://doi.org/10.1007/978-3-030-77722-7_24.

853 Ménétrier, B., T. Montmerle, L. Berre, and Y. Michel, 2014: Estimation and diagnosis of heterogeneous flow-dependent background-error covariances at the convective scale using either large or small ensembles. *Quarterly Journal of the Royal Meteorological Society*, **140 (683)**, 2050–2061, https://doi.org/10.1002/qj.2267.
857

Michel, Y., and P. Brousseau, 2021: A Square-Root, Dual-Resolution 3DEnVar for the AROME Model: Formulation and Evaluation on a Summertime Convective Period. *Monthly Weather Review*, **149** (9), 3135–3153, https://doi.org/10.1175/MWR-D-21-0026.1.

858

859

860

Montmerle, T., Y. Michel, E. Arbogast, B. Ménétrier, and P. Brousseau, 2018: A 3D ensemble variational data assimilation scheme for the limited-area AROME model: Formulation and preliminary results. *Quarterly Journal of the Royal Meteorological Society*, **144** (716), 2196–2215, https://doi.org/10.1002/qj.3334.

861

862

863

864

Petersen, W. A., H. J. Christian, and S. A. Rutledge, 2005: TRMM observations of the global relationship between ice water content and lightning. *Geophysical Research Letters*, **32** (14), 1–4, https://doi.org/10.1029/2005GL023236.

865

866

867

Price, C., and D. Rind, 1992: A simple lightning parameterization for calculating global lightning distributions. *Journal of Geophysical Research*, **97** (D9), 9919–9933, https://doi.org/10.1029/92JD00719.

868

869

870

Qie, X., R. Zhu, T. Yuan, X. Wu, W. Li, and D. Liu, 2014: Application of total-lightning data assimilation in a mesoscale convective system based on the WRF model. *Atmospheric Research*, **145-146**, 255–266, https://doi.org/10.1016/j.atmosres.2014.04.012.

871

872

873

Roberts, N. M., and H. W. Lean, 2008: Scale-selective verification of rainfall accumulations from high-resolution forecasts of convective events. *Monthly Weather Review*, **136** (1), 78–97, https://doi.org/10.1175/2007MWR2123.1.

874

875

876

Saunders, R., and Coauthors, 2018: An update on the RTTOV fast radiative transfer model (currently at version 12). *Geoscientific Model Development*, **11** (7), 2717–2737, https://doi.org/10.5194/gmd-11-2717-2018.

877

878

879

Seity, Y., P. Brousseau, S. Malardel, G. Hello, P. Bénard, F. Bouttier, C. Lac, and V. Masson, 2011: The AROME-France Convective-Scale Operational Model. *Monthly Weather Review*, **139** (3), 976–991, https://doi.org/10.1175/2010MWR3425.1.

880

881

882

Wattrelot, E., O. Caumont, and J.-F. Mahfouf, 2014: Operational implementation of the 1d+ 3d-var assimilation method of radar reflectivity data in the arome model. *Monthly Weather Review*, **142** (5), 1852–1873, https://doi.org/10.1175/MWR-D-13-00230.1.

883

884
Xiao, X., J. Sun, X. Qie, Z. Ying, L. Ji, M. Chen, and L. Zhang, 2021: Lightning Data Assimilation

Scheme in a 4DVAR System and Its Impact on Very Short-Term Convective Forecasting. Monthly

Weather Review, **149** (2), 353–373, https://doi.org/10.1175/MWR-D-19-0396.1.