Modelling and assimilation of lidar signals over Greater Paris during the MEGAPOLI summer campaign

Abstract. In this study, we investigate the ability of the chemistry tr ansport model (CTM) P OLAIR 3D of the air quality modelling platform P OLYPHEMUS to simulate lidar backscattered profiles from model aerosol concentration outputs. This investigat ion is an important preprocessing stage of data assimilation (validation of the observation operat or). To do so, simulated lidar signals are compared to hourly lidar observations performed during the MEGAPOLI (Megacities: Emissions, 5 urban, regional and Global Atmospheric POLlution and clima te effects, and Integrated tools for assessment and mitigation) summer experiment in July 2009, wh ere a ground-based mobile lidar was deployed around Paris on-board a van. The comparison is perf ormed for six different measurement days, 01, 04, 16, 21, 26 and 29 July 2009, corresponding to dif ferent levels of pollution and different atmospheric conditions. Overall, P OLYPHEMUS reproduces well the vertical distribution of 10 lidar signals and their temporal variability, especially f or 01, 16, 26 and 29 July 2009. Discrepancies on 04 and 21 July 2009 are due to high-altitude aerosol layers , which are not well modelled. In the second part of this study, two new algorithms for assimil ating lidar observations based on the optimal interpolation method are presented. One algorithm analyses PM10 (particulate matter with diameter less than 10 μm) concentrations. Another analyses PM 2.5 (particulate matter with diameter 15 less than2.5 μm) and PM2.5−10 (particulate matter with a diameter higher than 2.5 μm and lo wer than 10 μm) concentrations separately. The aerosol simulat ions without and with lidar data assimilation are evaluated using the Airparif (a regional operati onal network in charge of air quality survey around the Paris area) data base to demonstrate the feasibil ity and the usefulness of assimilating lidar profiles for aerosol forecasts. The evaluation shows that da t assimilation (DA) is more efficient at 20 correcting PM10 than PM2.5, probably because PM 2.5 is better modelled than PM 10. Furthermore,


Introduction
Aerosols are key air quality species to monitor and model as they impact vegetation as well as human health; impacts on the latter result from aerosol's penetration of the respiratory system, leading to possible respiratory and cardiovascular diseases (Kelly et al., 2011;Lauwerys et al., 2007;Dockery and Pope, 1996). They also impact visibility (Wang et al., 2009) and represent an uncertain component of climate changes due to their effects on the Earth's radiative budget (IPCC, 2007). For air quality, in order to simulate and predict particle concentrations, modellers have developed various chemistry transport models (CTM), e.g. EMEP (European Monitoring and Evaluation Programme) (Simpson et al., 2003), LOTOS (Long Term Ozone Simulation) -EUROS (European Operational Smog) (Schaap et al., 2004), CHIMERE (Hodzic et al., 2006), DEHM (Danish Eulerean Hemispheric Model) (Brandt et al., 2012)   . However, the aerosol vertical distribution is poorly quantified, because of numerous uncertainties on their sources (direct emissions) and on processes affecting their formation (e.g. nucleation, condensation, evaporation and coagulation), as well as on meteorological conditions. Since aerosol lifetime ranges from 1 to 10 days (Seinfeld and Pandis, 1998), improvements in the representation of their vertical distribution may lead to improved surface concentrations (lower error and higher correlation against observations) .
Various measurement types have been used to evaluate these models. The most frequently used data are in situ surface measurements, e.g. AirBase (http://www.eea.europa. eu/) and EMEP over Europe, BDQA (Base de Données de la Qualité de l'Air) Konovalov et al., 2009). However, they do not provide direct information on vertical profiles.
Satellite passive remote sensors (e.g. the Moderate Resolution Imaging Spectroradiometers -MODIS) and sunphotometer surface stations (e.g. the AErosol RObotic NETwork -AERONET) have greatly enhanced our ability to evaluate such models. Comparisons between observed and simulated aerosol optical depth (AOD) have been performed for global models and regional models (Kinne et al., 2006;Tombette et al., 2008;Péré et al., 2010). However, instruments such as sun photometers can only retrieve columnintegrated aerosol properties and can only work during daytime.
Since accurate vertical profiles of aerosols can be measured by aerosol lidars, lidar measurements were used in several campaigns, for example to evaluate the transport of particles (Chazette et al., 2012). Moreover, aerosol lidar networks such as the European Aerosol Research Lidar Network (EARLINET) are being developed at in situ sites. In space, measurements are performed with the Cloud-Aerosol Lidar with Orthogonal Polarization (CALIOP) lidar (Winker et al., 2007). Lidar measurements have been used for the validation of aerosol models. For example, Hodzic et al. (2004) compared vertical profiles simulated by CHIMERE with those observed by lidars from EARLINET, and Stromatas et al. (2012) used observations from the CALIOP space-based lidar. Royer et al. (2011) used optical-to-mass relationships (urban, pre-urban and rural) to retrieve the PM 10 (particulate matter with diameter less than 10 µm) concentrations from lidar signals (Raut et al., 2009a, b). In Royer et al. (2011), lidar-derived PM 10 concentrations were compared with simulations from POLYPHEMUS and CHIMERE during the MEGAPOLI (Megacities: Emissions, urban, regional and Global Atmospheric POLlution and climate effects, and Integrated tools for assessment and mitigation) summer experiment in July 2009. Data assimilation (DA hereafter) can reduce the uncertainties in input data such as initial or boundary conditions by coupling models to observations (Bouttier and Courtier, 2002). In air quality, applications of DA to PM 10 forecast using in situ surface measurements have been performed by Denby et al. (2008) and Tombette et al. (2009) over Europe, andPagowski et al. (2010), Pagowski and Grell (2012) and Li et al. (2013) over the United States of America. Over Europe, the efficiency of assimilating lidar measurements to improve PM 10 forecast has been compared to the efficiency of assimilating in situ surface measurements by Wang et al. (2013). Using an observing system simulation experiment (OSSE), they suggested that the assimilation of lidar observations may be more efficient to improve PM 10 forecast, although it depends on the number of lidar stations used. However, Wang et al. (2013) did not directly assimilate the lidar signal, but they used a relation between mass concentration and optical properties of pollution aerosol. Although this kind of relation has been determined for pollution aerosols over Greater Paris (Raut et al., 2009a), it needs to be generalised to other measurement sites before operationally assimilating the mass concentration converted from the lidar signal. Moreover, the uncertainly linked to the estimation of mass concentrations may be about 25 % (Raut et al., 2009a), which is mostly due to uncertainties in estimating the specific cross sections. Because uncertainties in the lidar signal may be less than 5 %, it is more accurate to directly assimilate lidar signals.
This paper aims at evaluating the lidar signals simulated by POLYPHEMUS and at testing new DA algorithms for assimilating lidar signals. We used measurements performed during the MEGAPOLI summer experiment, when a groundbased mobile lidar (GBML) was deployed around Paris onboard a van. Measurements from a ground-based in situ lidar at Saclay were also performed on 1 July 2009. The evaluation of lidar signals can also be regarded as a preprocessing stage of DA (validation of the observation operator). This paper is organised as follows. Section 2 describes the experiment setup, i.e. the chemistry transport model used (POLAIR3D) and the observations. In Sect. 3, the lidar observation operator is presented. Section 4 describes the evaluation of the simulation with in situ surface measurements and AERONET data. Results of the comparisons between observed and simulated lidar signals are shown in Sect. 5. Two new algorithms for the assimilation of lidar observations and results are shown in Sect. 6. The findings are summarised and discussed in Sect. 7.

POLAIR3D model
In this study, the POLAIR3D air quality model  of the air quality platform POLYPHEMUS, available at http://cerea.enpc.fr/polyphemus/ and described in Mallet et al. (2007), is used to simulate air quality over the Greater Paris area. Aerosols are modelled using the SIze-REsolved Aerosol Model (SIREAM-SuperSorgam), which is described in Debry et al. (2007) and Kim et al. (2011). SIREAM-SuperSorgam includes 20 aerosol species: 3 primary species (mineral dust, black carbon and primary organic species), 5 inorganic species (ammonium, sulfate, nitrate, chloride and sodium) and 12 organic species. Five bins logarithmically distributed over the size range 0.01-10 µm are used. The chemical mechanism CB05 (Carbon Bond version 5) is used for the gas chemistry (Yarwood et al., 2005). POLAIR3D/SIREAM has been used for several applications. For example, it was compared to in situ surface measurements for gas and aerosols over Europe by Sartelet et al. (2007Sartelet et al. ( , 2012 and Couvidat et al. (2012), and over Greater Paris by Couvidat et al. (2013); it was compared to AERONET data over Europe by Tombette et al. (2008) and to satellite data by Zhang et al. (2013); and it was compared to lidar-derived PM 10 over Greater Paris during MEGAPOLI by Royer et al. (2011).

Modelling setup and observational data
The modelling domain is the same as the one used in Royer et al. (2011) and Couvidat et al. (2013). It covers the Greater Royer et al. (2011) show that limited vertical model resolution leads to much smoother vertical profiles than those deduced from lidar signals, a finer vertical resolution is used with 23 vertical levels from the ground to 12 000 m, instead of the 9 vertical levels in Royer et al. (2011). The simulations are carried out for 1 month, from 28 June to 30 July 2009. Meteorological inputs are the same as in Couvidat et al. (2013). They are simulated with the Weather Research and Forecasting (WRF) model (Skamarock et al., 2008) using an urban canopy model and an undated Corine land-use database  with the Yonsei University (YSU) parameterisation (Hong et al., 2006) for the planetary boundary layer (PBL) dynamics . Anthropogenic emissions of gases and aerosols are generated using the Airparif (the Paris air quality agency) inventory for the year 2005. Boundary conditions for gaseous and particulate species were obtained from nested simulations over Europe and France, presented by Couvidat et al. (2013).
The GBML used during the MEGAPOLI campaign is based on an ALS450 lidar commercialised by the company LEOSPHERE and initially developed by the Commissariat à l'Energie Atomique (CEA) and the Centre National de la Recherche Scientifique (CNRS) (Chazette et al., 2007). It provides lidar measurements at 355 nm. The main characteristics of this lidar are detailed in Royer et al. (2011). This system is particularly well adapted to air pollution and tropospheric aerosol studies thanks to its full overlap reached at about 150-200 m height and its increased vertical resolution of 1.5 m. Measurement days of 1, 4, 16, 21, 26 and 29 July 2009, which correspond to different levels of pollution from Airparif (low, moderate or high), are used for comparisons to the lidar signal. Moreover, ground-based in situ lidar measurements were performed at Saclay (48.7 • N, 2.14 • E; 30 m a.s.l.) on 1 July 2009 from 06:49 to 16:44 UTC (the blue square in Fig. 1). These measurements are used for both the comparison and the assimilation of lidar bservations.
Airparif is the regional operational network in charge of air quality survey around the Paris area (http://www.airparif. asso.fr/). It provides hourly gas and/or aerosol concentrations (PM 10 and PM 2.5 ) measurements. Figure 1 shows the location of the Airparif stations by means of red squares and/or black triangles. There are 17 stations at which PM 10 and/or PM 2.5 concentration measurements are performed.
The AERONET programme is a federation of groundbased remote sensing aerosol networks established by NASA and PHOTONS (Univ. of Lille 1, CNES and CNRS-INSU), which provides a long-term, continuous and readily accessible public-domain database of aerosol optical measurements performed by sun photometers (Holben et al., 1998). Sun photometers measure AOD at different wavelengths ranging from 340 to 1024 nm. AOD data are computed for three data quality levels: level 1.0 (unscreened), level 1.5 (cloud-screened), and level 2.0 (cloud-screened and qualityassured). The uncertainty of AOD measurements is less than 0.02 (Holben et al., 2001). For this study, there are two available stations in Greater Paris: Paris (urban station; 48.87 • N, 2.33 • E; 50 m a.s.l.) and Palaiseau (suburban station; 48.70 • N, 2.21 • E; 156 m a.s.l.) (the green circles in Fig. 1). In this paper, level 2.0 AOD data at 340 and 380 nm are used to derive AOD data at 355 nm following the Ångström law where α is the Ångström exponent, defined by 3 Methodology Figure 2 describes the methodology used for lidar signal modelling from the outputs of the air quality model and for comparisons to measurements (aerosol concentration measurements, AOD data and lidar vertical profiles). This section presents the methodology used in POLYPHEMUS to derive the lidar observation operator.

Modelling of lidar signals
The range-corrected lidar signal PR 2 measured at an altitude z is defined by Collis and Russell (1976): where β m (β a ) is the molecular (aerosol) backscatter coefficient, α m (α a ) is the molecular (aerosol) extinction coefficient, and C is the instrumental constant for each channel depending on the technical characteristics of the emitting and receiving optics. In order to eliminate the instrumental constant C (because it is unknown), PR 2 is normalised as follows: where z ref is taken at an altitude in the molecular zone. In Eq. (4), to estimate the normalised lidar signal H , four optical parameters -β m , β a , α m and α a -are needed. The molecular backscatter coefficient (β m ) at the wavelength λ of the incident light is calculated by Nicolet (1984) where P is the pressure, T is the temperature, and k B is the Boltzmann constant, and the Rayleigh scattering cross section s Ray is given by s Ray = 4.678 × 10 −29 · λ −(3.916+0.074·λ+0.05/λ) .
The molecular extinction coefficient (α m ) is given by Nicolet (1984): Aerosol extinction and backscatter coefficients (β a and α a ) are functions of the particle sizes, of the aerosol complex refractive index (ACRI) of particles m, and of the wavelength λ of the incident light. With a population of different-sized particles of identical refractive index m and with a number size distribution function n(D wet ) with D wet the particle wet diameter, the aerosol extinction and backscatter coefficients are given by the following formulas: and where D max wet is a wet diameter upper limit for the particle population, a wet = π D wet λ a dimensionless size parameter, and Q ext (m, a wet ) and Q bsca (m, a wet ) are extinction and backscatter efficiencies respectively. These efficiencies are computed through the Mie code from ftp://ftp.giss.nasa. gov/pub/crmim/spher.f (de Rooij and van der Stap, 1984;Mishchenko et al., 2002). The dry complex refractive index (CRI) is interpolated from the OPAC package (Hess et al., 1998) for each species at the desired wavelength λ (355 nm). The CRI and densities used for calculation of optical properties are shown in Table 1. The wet diameter D wet is computed from the mean dry diameter of each section of the aerosol sectional model SIREAM and from the aerosol water content. The aerosol water content is calculated from the thermodynamic model ISORROPIA (Nenes et al., 1998a, b), which models the phase state (i.e. gas, liquid, solid) of inorganic aerosol species (i.e. ammonium, sodium, chloride, nitrate, sulfate) in equilibrium with gaseous precursors. The inputs of ISORROPIA are temperature, relative humidity (RH), gaseous precursor concentrations and inorganic aerosol concentrations. Because of the large amount of water vapour in the atmosphere, the ambient RH is assumed to be unaffected by the deliquescence of aerosol particles in ISORROPIA (Nenes et al., 1998a) and equals the water activity (referred to as a w ). The aerosol water content is estimated by the Zdanovskii-Stokes-Robinson (ZSR) relationship (Robinson and Stokes, 2002), where a w = RH, W is the aerosol water content concentration, M i is the molar concentration of species i (mol m −3 ) and m oi (a w ) is the molality of an aqueous solution of species i (mol kg −1 ). Computing the ACRI requires to make an assumption on the mixing state of the aerosol chemical species. The current version of POLYPHEMUS is based on an assumption of aerosol internal mixing: all the particles of a given size section at a given grid point of the domain are supposed to have the same chemical composition. Within this framework, Tombette et al. (2008) compared aerosol optical properties using two different assumptions for the black carbon (BC) mixing state: internally homogeneous mixing and core-shell mixing. In the internally homogeneous mixing case, BC is treated like the other components and a volumeweighted ACRI is calculated from the CRI of pure species. In the core-shell mixing case, each particle is assumed to have a structure: the core (BC) and the shell (all the other components). The hypothesis of an internally homogeneous mixing state seems to be unphysical as BC cannot be well mixed in the particle because of its complex geometry and solid state (Katrinak et al., 1993;Sachdeva and Attri, 2007). Tombette et al. (2008) have shown that the use of these two mixing states leads to negligible differences in AOD, but non-negligible differences on single-scattering albedo and absorption process. According to illustrations of Jacobson (2000), the BC mixing state influences the absorption cross section at small wavelengths (lower than 1 µm) for aerosols with diameters higher than 1 µm. Thus, a core-shell mixing hypothesis is used in this study. The Maxwell-Garnett approximation is used to calculate ACRI from the core CRI (i.e. BC in this study) and the shell CRI (where all the other components are well mixed) (Tombette et al., 2008).

Estimation of z ref
The altitude used to normalise the lidar signal does not need to correspond exactly to the beginning altitude of the molecular zone, but it could be any arbitrary chosen altitude in the molecular zone, where there is almost no aerosol. However, it is better to use an estimation of the normalisation altitude as close as possible to the beginning of the molecular zone, because lidar signals are attached to higher uncertainties at high altitudes because of a higher signal-to-noise ratio. Although the molecular zone is often determined visually from lidar vertical profiles, this method is not efficient for the treatment of large amounts of lidar profiles. We therefore created a new algorithm which can automatically estimate the normalisation altitude z ref from the lidar vertical profile.
The normalisation altitude z ref is estimated from the lidar signal and the simulated molecular signal PR 2,Ray , as follows: -Define a weight for each vertical point of the lidar signal (the vertical resolution is 1.5 m). The weights should be larger for the points that are more likely to be in the molecular zone, i.e. at high altitudes. We where h is the altitude of the points, h max is the maximal altitude considered (e.g. 4 km) and the parameter L is taken equal to 200 m.
-Fit all lidar signal vertical points (noted as a vector y) with a weighted least-absolute-deviations (LAD) regression (DasGupta and Mishra, 2007); the LAD regression is employed here because we are interested in the linear regression of lidar signal points at higher altitudes, e.g. the points between 2 and 3 km above the ground. However, it is difficult to know the altitude below which lidar signal points could be cut off for the estimation of z ref . When considering all available lidar signal points, the disturbances are prominently non-normally distributed and contain sizeable outliers (i.e. points at lower altitudes). In such cases, the least-squares method fails and the LAD method performs well (DasGupta and Mishra, 2007). In detail, we minimise to find a and b (cyan lines in Fig. 3).
-Calibrate the simulated molecular signal PR 2,Ray with the LAD regression line at altitude h max , and calculate the difference between the calibrated PR 2,Ray and the LAD regression line at each vertical point of the lidar signal in a loop starting from high altitudes to low altitudes. The altitude at which the difference becomes larger than a pre-assigned value (1 % of the value corresponding to the LAD regression line) corresponds to z ref . Figure 3 shows comparisons between the lidar signal and the simulated molecular signal PR 2,Ray for different lidar measurement days during MEGAPOLI. The simulated molecular signal (red lines in Fig. 3) agrees well with the lidar observations (black lines in Fig. 3) at high altitudes in the molecular zone, leading to the determination of the molecular zone and z ref .

Model evaluation
To evaluate air quality models, Boylan and Russell (2006) recommended a PM model performance goal and a criterion that are based upon an analysis of numerous PM and visibility-modelling studies. The PM model performance goal corresponds to the level of accuracy that is considered to be close to the best a model can be expected to achieve. The PM model performance criterion corresponds to the level of accuracy that is considered to be acceptable for modelling applications. The mean fractional bias (MFB) and the mean fractional error (MFE) are proposed by Boylan and Russell (2006) in order to evaluate model performances against observations: if the MFB is in the range [−30, 30 %] and the MFE is in the range [0, 50 %], then the PM model performance goal is met; if the MFB is in the range [−60, 60 %] and the MFE is in the range [0, 75 %], then the PM model performance criterion is met. The root-mean-square square (RMSE) and correlation are also often used in the aerosol modelling community. The statistical indicators are defined in Appendix A. As shown in Table 2, the model simulates PM 2.5 concentrations well, but PM 10 concentrations are underestimated. In other words, coarse particles (particulate matter with a diameter higher than 2.5 µm and lower than 10 µm) are underestimated. This may be because emissions and boundary conditions of coarse particles are underestimated: for example,    Boylan and Russell (2006) is met, despite a slight underestimation of AOD in agreement with the underestimation of PM 10 in comparison to Airparif observations (see Sect. 4.1).

Comparisons with lidar vertical profiles
The simulated lidar signal is compared with GBML observations performed during the MEGAPOLI summer experiment on the different measurement days (1, 4, 16, 21, 26 and 29 July 2009). The purpose of this section is to validate the ability of POLYPHEMUS to simulate lidar backscattered profiles and then choose suitable measurement days to do assimilation tests. On 1 July 2009, GBML measurements are performed leeward inside the pollution plume in the southwest of Paris between Saclay and Chateaudun during 3 h (black track in Fig. 1); this is the most polluted day of the MEGAPOLI experiment. High levels of PM 10 , on average about 45 µg m −3 (see Table 2), are measured by the Airparif network. Figure 4 presents the comparison between lidar observations and the simulation at 11:00, 12:00 and 13:00 UTC. It shows that POLYPHEMUS underestimates the lidar signal at 11:00 UTC but that it overestimates it at 12:00 UTC and agrees well with observations at 13:00 UTC. While the PBL height increases from about 1.2 to 1.8 km from 11:00 to 13:00 UTC and the GBML leaves the pollution plume (Royer et al., 2011), both the observed and simulated lidar signals decrease. Figures of the comparison between the simulation and observations from a ground-based in situ lidar at Saclay are shown later in this paper. The pollution plume covers Saclay because of the northeast wind. Thus high lidar signal values in both the simulation and observations are seen after 10:00 UTC, although the simulated lidar signals are underestimated. DA will be performed for this day, as it is the most polluted day with observations from both the GBML and a ground-based in situ lidar.
On 4 July 2009, GBML measurements are performed around Paris with a circular pattern from 14:49 to 17:24 UTC. Particle AOD and concentrations are underestimated in the simulation. The daily averaged AOD from the AERONET network is about 0.25; in the simulation it is 0.14 (see Table 3). The daily averaged PM 10 concentration from the Airparif network is about 18.37 µg m −3 , whereas it is 11.11 µg m −3 in the simulation (see Table 2). Figure 5 shows the comparison between the GBML measurements and the simulation at 15:03 and 16:00 UTC. The simulated lidar signals are underestimated. Moreover, lidar measurements show an aerosol layer between 2.0 and 3.0 km (probably from long-range transport), which is not present in the simulation, but impacts the lidar signal until low altitudes; this is mostly because boundary conditions do not provide information about this aerosol layer due to the large-scale model uncertainties.  Table 2). The observed and simulated AOD are 0.26 and 0.18 respectively (see Table 3).
The simulated AOD has a good correlation with AERONET data (up to 80 %). As deduced from the comparisons of the modelled and observed lidar signals in Fig. 6, the PBL height is well modelled until 12:00 UTC, but it is underestimated afterwards; for example, the PBL height is about 2.1 km from the observed lidar signal but it is about 1.6 km in the simulation. These differences in PBL height explain that the simulated lidar signal agrees better with the observation until 12:00 UTC. On 21 July 2009, the GBML travels from Saclay to the north of Paris across the Paris city centre. As shown in Fig. 7, the lidar signal is overestimated for this measurement day. However, the surface PM 10 concentration is underestimated: 27.84 and 16.84 µg m −3 (low-moderate level of pollution; see Table 2) from the Airparif network and from the simulation respectively. The large simulated lidar signals originate in high aerosol concentration at high altitudes, i.e. between 2.0 and 2.5 km, which leads to higher backscatter and extinction coefficients. This high-altitude aerosol layer originates in boundary conditions (large-scale model uncertainties), but is not present in the observations, and it impacts the lidar signal until low altitudes. This is why surface PM 10 is underestimated while lidar signal is overestimated.
On 26 July 2009, the GBML followed two circular patterns (the yellow and cyan tracks in Fig. 1). One is performed from 12:40 to 15:30 UTC at a distance between 15 and 30 km from the city centre. Another one is performed from 16:44 to 18:18 UTC in the south-southwest of Paris. Low levels of pollution are observed and simulated. Surface PM 10 concentration and AOD are underestimated. The daily averaged PM 10 concentration from Airparif is 18.04 µg m −3 , compared to 10.12 µg m −3 in the simulation. The mean observed AOD value is 0.15, compared to 0.08 in the simulation.
Although the lidar signal is slightly underestimated in the simulation, simulated and observed lidar signals agree fairly well, as shown in Fig. 8. The pollution from Paris is transported by the south wind to the north. This is why the lidar signal is higher at 14:00 UTC in Fig. 8. Because as much as 5 h of lidar measurements are performed, which is longer than on 4, 16, 21 and 29 July 2009, we will perform DA for this day.
On 29 July 2009, GBML measurements are performed from 12:22 to 15:10 UTC in the north of Paris and in peri-urban and rural areas. While low levels of pollution (12.33 µg m −3 of the mean PM 10 concentration in Table 2) are simulated, moderate levels of pollution (29.25 µg m −3 of the mean PM 10 concentration in Table 2) is observed by the Airparif network. As deduced from Fig. 9, at the beginning of measurement period, the PBL height is about 1.5 km and the simulated lidar signal agrees well with the lidar observations. At 15:00 UTC, the observed lidar signal has increased, due to an aerosol layer between 2.0 and 3.5 km. This layer is not simulated and the simulated lidar signal is underestimated.
For all measurement days, we also computed the statistics MFB ranges from −38 to 8 % and the MFE ranges from 3 to 38 %. Currently, there is no criterion to evaluate the comparisons for lidar signals. The criterion of Boylan and Russell (2006) was designed for PM concentration and light extinction. Because the scores of the lidar signal comparisons are extremely good compared to the criterion of Boylan and Russell (2006) with low errors and bias, the criterion of Boylan and Russell (2006) may not be restrictive enough for lidar signals.

Assimilation test of lidar observations
DA of lidar observations is performed for two out of the six different measurement days. Only these two days are retained due to the other days being cloudy and our algorithms not allowing us to assimilate lidar data when there are clouds. There is 13 h of cloud-cleaned measurements on 1 July, 5 h of cloud-cleaned measurements on 26 July and less than 3 h of cloud-cleaned measurements on the other measurement days. Therefore, DA run is performed on 1 and 26 July 2009 because too few data are available during the other measurement days.
In air quality, the large number of state variables leads to high computational costs when implementing DA algorithms. Among the widely used DA algorithms, the optimal interpolation (OI) is used here, as it is the most computationally efficient (Denby et al., 2008;Tombette et al., 2008;Wu et al., 2008;Li et al., 2013). In applications of DA to aerosol forecast, Tombette et al. (2009) have used the OI over western Europe for assimilating observations from the BDQA network, which covers France. Denby et al. (2008) have used two different DA techniques, the OI and ensemble Kalman filter, to assimilate PM 10 concentrations over Europe. Pagowski et al. (2010) used the OI over the United States of America for data assimilation of PM 2.5 observations, Li et al. (2013) used the OI for multiple aerosol species and for prediction of PM 2.5 in the Los Angeles Basin, and Wang et al. (2013) used the OI over Europe to investigate the potential impact of future ground-based lidar networks on analysis and short-term forecasts of PM 10 .

Basic formulation
A simple formulation for DA of lidar signals with OI is now described. Particles are represented in the model by mass concentrations of different chemical species for the different particle size sections. The state vector x is defined by where x h i,j,k is the mass concentration of the aerosol species j in section i for the horizontal spatial grid k at the model vertical level h, N b is the number of size sections, N s is the number of chemical species, n is the number of horizontal grid points at each vertical level h and l is the total number of vertical levels.
The analysed state vector is a solution to the variational optimisation problem: where J is the cost function defined by which leads to

Construction of error covariances
Since the measurements at different levels originate from the same lidar, the matrix R should not be diagonal because of measurement error correlations. However, in order to simplify R in the first tests of DA of lidar observations, one takes R = r I as a diagonal matrix, where I is the identity matrix and r is an error variance. The value of the observation error variance r is determined by a χ 2 diagnosis (Ménard et al., 1999), in which the scalar should be equal, on average, to the number of observations (N) at each DA step. Specifically, B plays a role in determining how the corrections of the concentrations should be distributed over the domain during DA. In practice, however, it is impossible to accurately know all coefficients of B. In our simulation, the number of model grid points is of the order of 10 5 . Thus the number of coefficients in the matrix B is about 10 10 multiplied by the square of the number of analysis variables (about 100 variables for particles are used here). Therefore, B is too large to be handled numerically.
In order to reduce the size of the error covariance matrices for background, we model the matrix B as follows where D is the error covariance matrix for PM 10 , defined by the Balgovind approach (Balgovind et al., 1983) obtained by considering the RMSE and correlation of simulated PM 10 concentrations. Thus, the size of D is much less than the one of B. The matrix P is defined by where M is equal to the dimension of the domain (l · n) and v k is a vector of size N b · N s (the number of state variables). Each component of v k corresponds to the proportion of the mass of particles for a given species in a given size section in PM 10 mass concentrations at grid point k. Let S ′ = S P be the directional derivative of S along a given direction, and let c b and c a be PM 10 concentration states before and after analysis respectively. In order to convert x into the PM 10 state c, we multiply each side of Eq. (18)  After the analysis, the concentrations c a are redistributed over particle species and size sections following the initial chemical and size distributions.

DA setup
DA experiments are carried out for 1 and 26 July 2009. All DA experiments are performed with a time step of 600 s and from 200 to 1800 m above the ground (10 model levels), since the lidar measurements are not available below the altitude of full overlap (200 m above the ground) and since aerosol concentrations above the PBL have limited impact on surface PM 10 in the short term . In the Balgovind approach (Balgovind et al., 1983), the horizontal correlation length is set to 0.2 • , which is estimated from numerical DA tests. The error variances are separately set for each DA level, depending on the RMSE of PM concentrations and the variability of PM concentrations at each model level.
Two new algorithms are tested for the assimilation of lidar observations. In the first algorithm, we use the assimilation of lidar observations to analyse PM 10 concentrations and the analysed PM 10 concentrations are redistributed over particle species and size sections following the initial chemical and size distributions (see Sect. 6.2). The background error variances of PM 10 concentrations are estimated by the simulation without DA and Airparif observations. The value of the observation error variance r is determined by a χ 2 diagnosis, which yields r = 1 µg 2 m −6 and r = 0.006 µg 2 m −6 respectively for 1 and 26 July, depending on the level of uncertainties (see Sect. 5). Let N be the number of lidar observations at one DA step. Figure 10 shows the time evolution of χ 2 /N (blue lines) for DA runs on 1 and 26 July. The mean over DA window of χ 2 /N is 1.02 (1.02) for 1 (26) July.
In the second algorithm, we separately analyse PM 2.5 and PM 2.5−10 (particulate matter with a diameter higher than 2.5 µm and lower than 10 µm) in the assimilation of lidar observations. We modify the matrices used in Sect. 6.2 to obtain c 2.5 and c 2.5−10 , the mass concentrations of PM 2.5 and PM 2.5−10 respectively (see Appendix B for details). We separately set the error variances for PM 2.5 and PM 2.5−10 in matrix D. Because of the lack of PM 2.5−10 observations, we can not directly estimate the background error variances. They are determined by the χ 2 diagnosis using the observation error variance r found in the first algorithm.
In the following, we note the assimilation with the first (second) DA algorithm as "DA (PM 10 )" ("DA (PM 2.5 and PM 2.5−10 )").

Results and discussions
The purpose of these DA tests is to verify whether these new algorithms are functional. Because we work at small scale, the corrections of DA are transported out of the simulation domain very quickly. Thus we only compute the statistics for the DA window to validate the DA tests. Table 4 presents statistics of the simulation results both without and with DA. Statistics are computed for both PM 10 and PM 2.5 concentrations. Overall, both DA algorithms lead to better scores (lower RMSE, MFB and MFE, and higher correlation) than the simulation without DA for PM 10 concentrations. Comparing the two DA algorithms, the simulation with DA (PM 2.5 and PM 2.5−10 ) leads to better scores than the simulation with DA (PM 10 ) for PM 10 concentrations (see Table 4). The RMSE of PM 10 is 11.63 µg m −3 in the simulation with DA (PM 2.5 and PM 2.5−10 ), compared to 13.69 µg m −3 in the simulation with DA (PM 10 ) on 1 July. The RMSE of PM 10 is 4.73 µg m −3 in the simulation with DA (PM 2.5 and PM 2.5−10 ), compared to 6.08 µg m −3 in the simulation with DA (PM 10 ) on 26 July. It is because higher background error variances are set for the coarse sections in the simulation with DA (PM 2.5 and PM 2.5−10 ). However, the simulation with DA (PM 2.5 and PM 2.5−10 ) leads to similar scores to the simulation with DA (PM 10 ) for PM 2.5 concentrations (see Table 4). It is because similar background error variances for PM 2.5 in the simulation with DA (PM 2.5 Table 4. Statistics (see Appendix A) of the simulation results (PM 10 and PM 2.5 ) without DA and with DA for the Airparif network for 1 and 26 July 2009. "With DA (PM 10 )" stands for the assimilation of lidar observations correcting PM 10 directly. "With DA (PM 2.5 and PM 2.5−10 )" stands for the assimilation of lidar observations correcting PM 2.5 and PM 2.5−10 separately. and PM 2.5−10 ) to the simulation with DA (PM 2.5 ) are used in the χ 2 diagnosis, since fine particles contribute to more than 80 % of the lidar signal (Randriamiarisoa et al., 2006). In the following, we compare the simulation without DA and the simulation with DA (PM 2.5 and PM 2.5−10 ). On 1 July, the averaged RMSE of PM 10 is 11.63 µg m −3 with DA (PM 2.5 and PM 2.5−10 ), compared to 17.74 µg m −3 without DA. The decrease of the RMSE are explained by the correlation length in the matrix D, since no Airparif station performs measurements in the southwest of Paris (the northeast wind). At the Issy-Les-Moulineaux station (48.82 • N, 2.27 • E; 36 m a.s.l.), the closest station to Saclay, the RMSE of PM 10 is 14.72 µg m −3 with DA (PM 2.5 and PM 2.5−10 ), compared to 22.81 µg m −3 without DA. However, the averaged RMSE of PM 2.5 is about 10.4 µg m −3 with DA (PM 2.5 and PM 2.5−10 ), compared to 8.54 µg m −3 without DA. This is due to the larger horizontal correlation length (see Sect. 6.3). While DA runs increase PM concentrations in the lidar measurement grids, PM concentrations are increased at Airparif stations, where PM 2.5 concentrations is well simulated and coarse particles are underestimated. This problem can be solved by decreasing the horizontal correlation length. Figure 11 shows that the model underestimates the lidar signal at Saclay. The simulation with DA better simulates the lidar signal than the one without DA. It means that DA corrects the model aerosol concentrations well (the closer to the truth the model aerosol concentrations are, the better the lidar signals are simulated).
On 26 July, the averaged RMSE of PM 10 is 4.73 µg m −3 with DA (PM 2.5 and PM 2.5−10 ), compared to 6.67 µg m −3 without DA. Because two circular GBML travelling patterns were performed around Paris (see Fig. 1), most of the Airparif stations are either leeward (the south wind) or close to the patterns of GBML. These patterns could validate improvements of PM concentrations. At the Paris 1er Les Halles station (48.86 • N, 2.35 • E; 35 m a.s.l.), the RMSE of PM 10 is 1.96 µg m −3 in the simulation with DA (PM 2.5 and PM 2.5−10 ), compared to 4.71 µg m −3 in the simulation without DA. Moreover, DA runs lead to better scores than the simulation without DA for PM 2.5 . At the Creil Faiencerie leeward station (49.26 • N, 2.47 • E; 28 m a.s.l.), the RMSE of PM 2.5 is 4.1 µg m −3 in the simulation with DA (PM 2.5 and PM 2.5−10 ), compared to 4.9 µg m −3 in the simulation without DA.

Conclusions
In order to investigate the ability of the CTM POLAIR3D of the air quality modelling platform POLYPHEMUS to simulate lidar vertical profiles, we performed a simulation over the Greater Paris area for the summer month of July 2009. The results (PM 10 and PM 2.5 concentrations) are evaluated using Airparif data. We simulated aerosol optical properties and lidar signals from the model aerosol concentration outputs using the ACRI and the wet particle diameter. The AOD was evaluated using AERONET data: the RMSE ranges from 0.07 to 0.20, the MFB ranges from −58 to −21 % and the MFE ranges from 29 to 58 %. According to the criterion of Boylan and Russell (2006), the model performance criterion is met for AOD. Hourly comparisons between simulated lidar signals and lidar observations were described for six measurement days during the MEGAPOLI summer campaign. These comparisons showed a good agreement between GBML measurements and the simulation except for 4 July 2009, when an aerosol layer was not modelled at high altitudes but observed in lidar measurements, and for 21 July 2009, when an aerosol layer was modelled at high altitudes but not observed in lidar measurements. The statistics obtained for the lidar comparison are extremely good compared to the criterion of Boylan and Russell (2006), with low errors and bias: the MFB ranges from −38 to 8 % and the MFE ranges from 3 to 38 %. Because the criterion of Boylan and Russell (2006) was designed for PM concentration and light extinction, they may not be restrictive enough for lidar signals. A specific criterion would therefore need to be designed. Overall, the results show that the optical property module of POLYPHEMUS models lidar signals well. Two new algorithms for the assimilation of lidar observations based on the optimal interpolation method were presented: one algorithm analyses PM 10 concentrations, another analyses PM 2.5 and PM 2.5−10 concentrations separately. DA tests were performed for only 1 and 26 July 2009, because the other measurement days were cloudy and our algorithms do not allow us to assimilate lidar data when there are clouds. Both of these algorithms lead to better scores (lower RMSE, MFB and MFE, and higher correlation) for PM 10 and PM 2.5 on 26 July 2009. However, they did not improve PM 2.5 on 1 July 2009, because of the large horizontal correlation length. The simulation with DA (PM 2.5 and PM 2.5−10 ) leads to better scores than the simulation with DA (PM 10 ) because the error variances for backgrounds are set separately for fine (PM 2.5 ) and coarse (PM 2.5−10 ) particles. The results shown in this paper suggest that the assimilation of lidar observations that analyses PM 2.5 and PM 2.5−10 would perform better than the assimilation of lidar observations that analyses PM 10 , but it is computationally more costly.

Atmos
Comparing the simulation without DA and the simulation with DA (PM 2.5 and PM 2.5−10 ), the averaged RMSE of PM 10 is 11.63 µg m −3 with DA (PM 2.5 and PM 2.5−10 ), compared to 17.74 µg m −3 without DA on 1 July 2009. The averaged RMSE of PM 10 is 4.73 µg m −3 with DA (PM 2.5 and PM 2.5−10 ), compared to 6.67 µg m −3 without DA on 26 July 2009.
A forthcoming paper will present results about the assimilation of continuous measurements from the AC-TRIS/EARLINET network during a 72 h period of intensive observations.