Associate Editor: Gabriele Pfister; Atmospheric Chemistry Division, National Center for Atmospheric Research, US
Analyses and re-analyses weather products based on global circulation models (GCMs) are generated through the assimilation of various observations at difference temporal and spatial scales. They are usually designed to reproduce the large-scale or mesoscale patterns of atmospheric circulations. With the continuing development of models at higher resolution, more detailed products become available. The European Center for Medium Range Weather Forecast (ECMWF) produces analyses and forecasts at various spatial resolutions up to 0.08° × 0.08° globally. Although state-of-the-art for global scale products, such spatial resolutions may be insufficient to represent small-scale and terrain-driven meteorological features (Carvalho et al., 2014). Thus, regional numerical models, such as the Weather Research and Forecasting Model (WRF, Skamarock et al., 2008) are used to describe the atmospheric state for small-scale regional meteorological fields or to account for specific dynamical processes, e.g. in urban areas.
Over the years, many studies have been carried out to investigate the impact of various dynamical processes within the WRF/Urban modeling system, over different time and meteorological conditions. The selection of physics schemes has been the research priority when evaluating the performance of WRF in representing the urban meteorological patterns. For instance, Sarmiento et al. (2017) used the WRF model to test the sensitivity of simulated meteorological fields to different model packages for the city of Indianapolis, Indiana. According to this study, the Building Effect Parameterization (BEP) urban canopy model and the Bougeault & Lacarrere (BouLac) PBL scheme generated the most accurate meteorological variables in the wintertime, while the magnitude of urban canopy impacts was much smaller in summer. In addition to the choice of the physics schemes, another sensitive setting is the four-dimensional data assimilation (FDDA), also referred to as nudging. This approach adjusts the modeled values toward analysis fields or observations through the introduction of artificial forcing terms in the prognostic equations (Stauffer and Seaman 1994). The WRF model includes several nudging methods such as 3-D analysis nudging (grid nudging and spectral nudging), surface analysis nudging, and observation nudging. Within the WRF model package, the objective analysis program (OBSGRID) generates three output products that provide initial and boundary conditions, as well as the inputs for surface analysis and observation nudging based on surface and upper air meteorological observation data. Previous studies demonstrate that nudging improves the accuracy of WRF results (e.g. Otte, 2008; Rogers et al., 2013). The use of all of three outputs from the OBSGRID allows the assimilation of observations, which may be needed for the dynamical downscaling. For instance, Li et al. (2016) found a clear improvement of correlations between observed and simulated surface zonal and meridional winds (6~11%) and temperatures (9%) when observation nudging was used. While nudging creates an artificial forcing term that reduces deviations from observed values, it may also constrain the dynamics of the model reducing its ability to generate small-scale features in locations where no observations are present. Whether or not strong nudging is needed depends on the objective, such as the analysis of atmospheric dynamic processes, the evaluation of external forcing impacts, or the generation of meteorological fields in support to air quality analyses.
Urban areas are significant sources of fossil fuel CO2 emissions (Duren and Miller, 2012). CO2 monitoring networks combined with an atmospheric inversion approach through the use of an atmospheric transport model can provide a “top-down” method for the quantification of city emissions (e.g. Boon et al., 2016; Lauvaux et al., 2016; Gurney et al., 2012). The uncertainties associated with this method are mainly attributed to measurement errors, prior emission biases and atmospheric transport modeling errors (Feng et al., 2016). Typical atmospheric in situ CO2 measurements in urban areas are collected on top of buildings or on towers that integrate the signal of fluxes and transport on scales expected to be larger than 1 km, hence the importance of advection in transport models. Meteorological parameters such as temperature, wind and vertical mixing have a critical impact on the accuracy of atmospheric transport and chemistry models. Arguably, land use and land cover change (LUCC), urban geometry and anthropogenic heat emissions in urban areas generate specific forcing and surface boundary condition that impact the surface energy budget and create heat islands, modify wind direction and generally induce lower winds (Lee et al., 2011; Bonacquisti et al., 2006). To which extent these complex processes must be included in the generation of atmospheric transport fields for city scale CO2 inversion research remains an open question and may vary from city to city. Relevant studies that modeled CO2 concentrations over large complex urban areas such as Los Angeles (Feng et al., 2016) and Salt Lake City (Nehrkorn et al., 2013) show that the accounting of urban effects in high-resolution meteorological fields provided by the WRF model may improve the emission estimates. McKain et al. (2012) also indicated that high-resolution data, including small-scale and terrain-driven features of the atmospheric flow, must be accounted for to derive robust estimates of the CO2 city emissions, up to a level suitable for detecting temporal changes.
Nevertheless, the city of Paris is located over a rather flat area where transport is dominated by synoptic scale events. Previous studies have attempted to estimate the Paris urban CO2 emission from atmospheric concentration measurements combined with atmospheric transport modeling (Bréon et al., 2015; Staufer et al., 2016). These studies used the CHIMERE transport model driven by operational forecasts from ECMWF without any urban surface effects at an approximate spatial resolution of 16 km × 16 km and a temporal resolution of 3 h that were linearly interpolated to the CHIMERE grid on an hourly basis. Although successful, these analyses show some mismatches between the measured and modeled CO2 concentration, even after the emissions have been optimized. The impact of the mega-city urban area on the atmospheric transport may be a cause for the observed mismatch. Indeed, as the impacts of the urban land cover on the atmospheric transport, and thus on atmospheric CO2 concentrations, are not considered, it may impact the link between emissions and concentrations. In this context, the overall goal of this study is to assess whether or not the WRF model, in its various configurations, provides a better representation of the meteorological variables over the Île-de-France region (IdF region) containing the Paris megacity than the ECMWF high-resolution operational forecasts (HRES). This assessment focuses on three variables, namely the air temperature, the wind speed and the PBL height, for which observation data are available within Paris and its peri-urban area as well as in rural areas across the whole IdF. The results are discussed in the context of the ability of these two models to provide the necessary data as input to transport models for CO2 and other tracers, with a focus on the impact of urban land cover on this transport.
This paper is structured as follows: Section 2 describes the WRF model, with a focus on the settings that are necessary for this work and the experiment design, together with the ECMWF HRES and the atmospheric measurements that are used for the evaluation. Section 3 examines in detail the performances of ECMWF forecasts and the results from varying configurations of WRF against meteorological measurements. Conclusions are summarized in Section 4.
2.1 Description of WRF settings
We use the version 3.7.1 of WRF-ARW, coupled with the Single-Layer Urban Canopy Model (UCM) (Tewari et al., 2007; Chen et al., 2011). The simulations are initialized at 00 UTC on 1 January 2015 to allow enough spin-up before February 2015, which is the one-month period of model evaluation. The ERA-Interim global re-analyses data from ECMWF (0.75° × 0.75°) with 6 h interval are used as the boundary and initial conditions for the model. Three one-way nested domains are employed with horizontal grid spacing of 27, 9, and 3 km, covering Europe (Domain 01), France (Domain 02) and the IdF region (Domain 03) respectively. Each domain has 34 vertical layers extending from the surface to 100 hPa wherein 15 layers are arranged below 1.5 km and the height of the first layer top is approximately 19 m above the ground level.
Land surface information such as topography, land use and urban parameterization are key drivers of this modeling effort. Land use defines many vital parameters that impact the energy budget and atmospheric transport model elements, such as surface albedo, emissivity and moisture availability. The two default land use data sets provided in the WRF model were derived from the satellite measurements collected by the United States Geological Survey (USGS) from April 1992 to March 1993 (Hansen and Reed, 2000) and the Moderate Resolution Imaging Spectroradiometer (MODIS) over the period spanning January to December 2001 respectively. Neither of them considers the recent land use changes in the IdF region. We thus replaced the USGS land use data with the 2006 Coordination Information on the Environment (CORINE) 250-meter data from the European Environment Agency (EEA) (http://www.eea.europa.eu/data-and-maps/data/corine-land-cover-2006-raster). The conversion of the land use categories from CORINE to USGS types was based on the equivalence relations provided by Pineda et al. (2004) (Figure 1). As for the topography, the default WRF data were derived from USGS with 30 arc seconds resolution (≈0.9 km at the equator). Also in this study, the default terrain data were replaced with the 1 arc-second (≈30 meters) ASTER Global Digital Elevation Model Version 2 (GDEM v2) (Tachikawa et al., 2011).
The UCM is coupled with the Noah land surface model (Chen and Dudhia, 2001) of WRF to improve its description of lower boundary conditions, and better represent the physical processes involved in the exchange of heat, momentum, and water vapor in urban areas (Tewari et al., 2007). The geometric and thermal parameters for the UCM model over Paris are based on the work of Kim et al. (2013). The values of building height, roof width, road width are 12 m, 3.75 m and 11.25 m respectively. The default urban and build-up land fractions for the three different urban types in WRF, namely the low density residential area, high density residential area and industrial & commercial/transportation land use, were 0.5, 0.9 and 0.95 respectively. Previous results show that the use of the original constant urban canopy parameters produces unrealistic gradients of some meteorological factors, such as temperatures, at the junctions between urban and rural areas, whereas the incorporation of a 2-D UCM yields more consistent transitions (Ching et al., 2009; Lin et al., 2016). In this study, the existing UCM model with fixed urban land cover types within each grid-cell (FRC_URB) were replaced by with those obtained from a 2-D map generated from the CORINE land cover dataset as shown in Figure 1.
2.2 Experiment designs
The following physics options were used in the control run denoted as REF in Table 1: WSM6 microphysics scheme (Hong and Lim, 2006), RRTM longwave radiation scheme (Mlawer et al., 1997), Dudhia shortwave radiation scheme (Dudhia, 1989), YSU PBL scheme (Hong et al., 2006), Revised MM5 Monin-Obukhov surface layer (Jiménez et al., 2012), Noah land surface model. The BMJ cumulus convection scheme (Janjic, 1994, 2000) was applied for Domain 01 and Domain 02 only. In this work, we used the combination of grid nudging, surface analysis nudging and observation nudging together with the objective analysis to maximize the benefit of assimilating observations. The National Centers for Environmental Prediction (NCEP) operational global observation surface data (ds461.0) and upper-air data (ds351.0) from the Research Data Archive at the National Center for Atmospheric Research (https://rda.ucar.edu/) were used as input of the nudging routines to constrain the model. Grid and surface analysis nudging was applied to the mass fields and wind components for Domain 01 as well as for Domain 02 but with a reduced strength to ensure that the model create its own structures within the high-resolution nest. For the same reason, grid nudging was applied to vertical levels above the PBL only. Observation nudging was used for all domains and all layers. Details regarding the nudging parameterizations for the control runs are described in Table 2.
|mp6_bl1_sf1_ral1_ras1_cu2_su2||REF (WRF_OA)||Microphysics (mp): WSM6 (Hong and Lim, 2006),
PBL (bl): YSU (Hong et al., 2006),
Surface layer (sf): MM5 (Jiménez et al., 2012),
Longwave radiation (ral): RRTM (Mlawer et al., 1997),
Shortwave radiation (ras): Duhia (Dudhia, 1989),
Cumulus convection (cu): BMJ (Janjić, 1994) only in outer domains,
Land surface (su): Noah land surface model (Chen and Dudhia, 2001),
With UCM (Chen and Dudhia, 2001),
With OBSGRID: Grid nudging + Surface Nudging + Observation Nudging + Objective analysis
|Spectral/Grid nudging||Surface Nudging||Observation Nudging|
|Wind field||Nudging above PBL||Nudging above PBL||N/A||Nudging||Nudging||N/A||Nudging all layers||Nudging all layers||Nudging all layers|
|Mass field||Nudging above PBL||Nudging above PBL||N/A||Nudging||Nudging||N/A||Nudging all layers||Nudging all layers||Nudging all layers|
The WRF model offers a variety of alternative physics schemes to represent dynamic and physical processes as well as several nudging methods. All of these schemes and methods influence the numerical weather prediction model output results. As part of this work, we selected several schemes commonly adopted in relevant references (e.g., Borge et al., 2008; Stegehuis et al., 2015) and test 32 combinations. A set of numerical experiments were then carried out in addition to REF to assess the sensitivity of the WRF model to the choice of different physics schemes with changing one scheme at a time (16 combinations). Another set of WRF simulations without the OBSGRID processing, namely only with the spectral nudging and different physics schemes were conducted to examine the impact of the nudging strategy on model performances (16 combinations). Table 3 summarizes the sensitivity experiments of physics schemes and nudging methods used in this study.
|Microphysics||mp1_bl1_sf1_ral1_ras1_cu2_su2||MP1||Kessler (Kessler, 1969)|
|mp2_bl1_sf1_ral1_ras1_cu2_su2||MP2||Lin (Purdue) (Lin et al., 1983)|
|mp3_bl1_sf1_ral1_ras1_cu2_su2||MP3||WSM3 (Hong et al., 2004)|
|mp4_bl1_sf1_ral1_ras1_cu2_su2||MP4||WSM5 (Hong et al., 2004)|
|mp5_bl1_sf1_ral1_ras1_cu2_su2||MP5||Eta (Ferrier) (Rogers et al., 2001)|
|Radiation (Longwave + Shortwave)||mp6_bl1_sf1_ral1_ras5_cu2_su2||RA5||RRTM + New Goddard (Mlawer et al., 1997; Chou and Suarez, 1999)|
|mp6_bl1_sf1_ral4_ras4_cu2_su2||RA4||RRTMG + RRTMG (Iacono et al., 2008)|
|Cumulus Convection||mp6_bl1_sf1_ral1_ras1_cu1_su2||CU1||Kain-Fritsch (Kain, 2004)|
|mp6_bl1_sf1_ral1_ras1_cu3_su2||CU3||Grell-Freitas (Grell et al., 2013)|
|mp6_bl1_sf1_ral1_ras1_cu5_su2||CU5||New Grell (Grell, 1993; Grell and Dévényi, 2002)|
|mp6_bl1_sf1_ral1_ras1_cu93_su2||CU93||Grell-Devenyi (Grell and Dévényi, 2002)|
|PBL + Surface Layer||mp6_bl6_sf5_ral1_ras1_cu2_su2||PBL6||MYNN3 + MYNN (Nakanishi and Niino, 2006, 2009)|
|mp6_bl7_sf1_ral1_ras1_cu2_su2||PBL7||ACM2 + MM5 (Pleim, 2007a, b; Jiménez et al., 2012)|
|mp6_bl8_sf1_ral1_ras1_cu2_su2||PBL8||BouLac + MM5 (Bougeault and Lacarrere, 1989; Jiménez et al., 2012)|
|OBSGRID||mp6_bl1_sf1_ral1_ras1_cu2_su2||WRF_noOA||Without OBSGRID (Spectral nudging)|
2.3 Meteorological measurements
The meteorological observation data used in this study were obtained from Météo France, the French weather service, and the Site Instrumental de Recherche par Télédétection Atmosphérique (SIRTA) research data platform (http://sirta.ipsl.polytechnique.fr). Observation sites were located as shown in Figure 2. The names of the stations contained in the NCEP operational global observation subsets to be used in the observation nudging are shown in black, while those stations that are not used for nudging are shown in red. The stations were additionally classified as either urban, suburban or rural according to the dominant land use type and the proportion of built-up areas within the corresponding model grid cell. All meteorological measurements were acquired in accordance with the WMO Surface Meteorological Observations: the hourly air temperature was measured at a height of 1.50 meters above ground level, and the hourly wind speed and direction was obtained at 10 meters above the surface. In addition to these surface measurements, radiosonde vertical profiles of temperature and wind vectors were also obtained twice daily (00 and 12 UTC) at the Trappes station (shown by a blue triangle with a circle, WSW of Paris in Figure 2).
2.4 ECMWF operational forecasts
The ECMWF high-resolution operational forecasts (HRES) with the Integrated Forecasting System (IFS) were used in this study. This product is available every 12 h with a forecast step up to 10 days based on a spectral general circulation model and a 4D-Var data assimilation system which contains satellite data as well as ground-based observations. It provides meteorological fields with a horizontal resolution of about 16 km (TL1279) and 137 vertical levels wherein about 24 layers are below a height of 1.5 km (ECMWF, 2015). This high vertical resolution was chosen with the aim of achieving a better representation of physical processes, clouds, inversions and vertically propagating gravity waves. The feedback data from the ECMWF meteorological archival and retrieval system (MARS) archive were downloaded so as to get the list of meteorological observation stations within the IdF region that are assimilated in each ground-based observation type or platform such as SYNOP, METAR and BURF (http://apps.ecmwf.int/mars-catalogue/?class=od&stream=oper&expver=1).
3.1 Surface meteorology statistics
The ECMWF data were archived at three-hour resolutions, while the WRF model output meteorological variables were saved with hourly resolution. We therefore linearly interpolated the ECMWF data on an hourly basis and compared the performance of hourly ECMWF-analyzed and WRF-modeled data with respect to the observed surface temperature and wind speed for February 2015 based on the correlation coefficient (R), root-mean-square error (RMSE) and mean bias error (MBE) metrics (Table 4). The difference (Δd) between each pair of simulated and observed values and RMSE are used in the assessment of the ability to reproduce the surface wind direction considering its circular nature (Jiménez et al., 2013). In order to fit the specific point of observation, a horizontal bilinear interpolation was performed to extract the WRF gridded data to the exact point of interest. As for ECMWF, we simply used the cell that contains the observation location as the interpolation using such a coarse resolution would be less meaningful.
|Correlation coefficient (R)|
|Mean bias error (MBE)|
|Root-mean-square error (RMSE)||or for wind direction|
The statistics for surface temperature and wind speed considering all weather stations shown in Figure 2 are summarized in Table 5. The statistics for the WRF_noOA results as compared to the observations are substantially inferior to those of the two other models, i.e. WRF with the objective analysis (WRF_OA) and using the ECMWF outputs. The WRF model with only the spectral nudging as an external constrain is substantially different from the reality so that the objective analysis and multi-nudging with surface observations is necessary. Indeed, the WRF_OA results are much better and in some respect superior to those of ECMWF as shown below. Assimilation and nudging strongly decreases the RMSE and improves the correlation to observations, although it has little positive impact on the biases.
|Temperature 2 m (°C)||Urbana||0.95||0.72||0.93||1.31||2.47||1.46||0.82||0.77||–0.86||4.27|
|Wind speed 10 m (m/s)||Urbana||0.85||0.70||0.84||1.59||1.78||1.14||–1.07||–0.99||0.36||3.47|
|Median of Δd||RMSE|
|Wind direction 10 m (°)||Urbana||5.76||11.21||5.39||38.67||46.62||41.30|
Air temperature at 2 m. This variable is well reproduced by both WRF_OA and ECMWF forecasts with correlations of (0.95, 0.96, 0.94) and (0.93, 0.94, 0.93) respectively for (urban, suburban and rural) stations (hourly data). The RMSE is within the range of 0.93 to 1.31°C for WRF_OA with similar values for ECMWF (Table 5). In terms of the MBE, the surface temperatures are in general overestimated by WRF and underestimated by ECMWF.
Wind speed at 10 m. Generally, the performance for the modeled wind speed is not as good as for temperature fields, confirming previous findings that winds are strongly affected by various factors like topographical conditions, near-surface thermal and dynamic effects and complex turbulent motions in atmospheric boundary layer. The correlations with observations for the two models are similar but the values of RMSE and MBE are mostly better for ECMWF. The maximum RMSE of the model-observation discrepancies is 1.64 m/s for WRF_OA and 1.55 m/s for ECMWF. The analysis of the MBE shows that the surface wind speeds are mainly underestimated by the WRF model, with a bias of up to 1.07 m/s for the WRF_OA model at urban sites, while it gives positive biases for ECMWF.
Wind direction at 10 m. While the statistics for the wind speed tend to favor ECMWF with respect to WRF_OA, the opposite is true for the wind direction. The median of Δd at suburban and rural sites of WRF_OA is much smaller than that of ECMWF, indicating a better representation of the frictional turning. The performance of ECMWF and WRF_OA are comparable at urban stations, and it also shows the best scores of RMSE among three types of observation stations, with a value of 38.67° for WRF_OA and 41.30° for ECMWF.
Ten-meter winds in models are generally diagnosed from model winds using strong assumptions on the vertical profile to which they can be fairly sensitive. As a consequence, the comparison of grid-scale model results with a point observation may be a source of discrepancy. Since larger errors and more variability are typically associated with lower wind speeds both for observations and the model, the statistics are recalculated by only using wind speeds higher than or equal to 3 m/s (Table S1). As expected, the RMSEs of wind direction decrease when the wind speed threshold is increased. This is true both for ECMWF forecasts and the WRF model. While the wind speed bias of the WRF model could not be improved due to the less appropriate urban land surface characteristics in the UCM that were not calibrated for the urban form of the Paris city.
One may argue that the statistics provided in Table 5 and the associated interpretation maybe biased as they are partly based on observations that are used in both the ECMWF and WRF_OA data assimilation process. In order to address these concerns, we also calculated similar statistics with a reduced set of stations observations that are not used in the assimilation process. These results are shown in Table 6. The data in this table show that WRF_OA provides obvious improvements in modeled surface temperature and wind speed over WRF_noOA, even if the analysis is restricted to the stations that are not included in the WRF_OA assimilation. This indicates that the implementation of the refined nudging method does not only improve the statistical agreements at assimilated sites, but also the accuracy of WRF_OA results over other cross-validation sites.
|Temperature 2 m (°C)||Urbana||0.92||0.66||0.93||1.63||2.80||1.55||0.85||0.77||–0.88||4.32|
|Wind speed 10 m (m/s)||Urbana||0.85||0.66||0.85||1.30||1.58||1.35||–0.56||–0.50||0.86||2.94|
|Median of Δd||RMSE|
|Wind direction 10 m (°)||Urbana||7.56||10.85||6.20||50.92||56.14||53.10|
The results presented so far were for a single month. We have also analysed whether the results hold for other months in the year. A one-year simulation was carried out with WRF-ARW Version 3.9.1 from October 2015 to November 2016, with the first 2 months being the model spin-up and not used in the analysis. This experiment used the same model configurations as the reference one (REF), except that it has a higher horizontal resolution than REF with grid spacing of 9, 3 and 1 km for domains d01, d02 and d03 respectively. The seasonal statistics of the agreement between the new one-year WRF model and observations for temperatures and winds during that period are shown in the supplementary (Table S2). The results are mostly consistent with the conclusions that are based on a one-month simulation. Note that the full-year simulation further reveals an underestimation of temperature for ECMWF and an overestimation by WRF_OA during spring (March–April–May or MAM) and summer (June–July–August or JJA). In WRF_OA, these warm biases may be due to the fact that anthropogenic heat release in UCM over the Paris urban area was set to a constant for each month, whereas in reality, more heat should be released in cold months from heating and less in warm months, given the lack of systematic air-conditioning in IdF buildings. The values of RMSEs for wind speeds always tend to slightly favor ECMWF with respect to WRF_OA and the correlations provided by ECMWF are also a bit closer to the observations than those provided by WRF except for summer. In contrast, WRF_OA outperforms ECMWF in matching the observed wind direction for most of the year. Additionally, larger wind direction misfits are associated with lower wind speeds in summer and autumn (September–October–November or SON) as compared to the other two seasons.
In order to answer the question that whether the higher temporal resolution of the WRF model better describes other temporal variabilities, we computed the differences between the hourly data and the 24-hour moving average for each station. This procedure retains the sub-daily variability and eliminates the temperature variation at synoptic temporal scales. The statistics are presented in Figure 3. It shows that ECMWF and WRF_OA match the observed values better than WRF_noOA. The ECMWF and WRF_OA have comparable performance for the sub-daily variability of temperatures while the WRF_OA is slightly better at some suburban sites with statistics that are fairly close to the observations. This is not only true from a correlation perspective, but also in terms of biases. It is worth noting that despite their superior performance in temperature, both ECMWF and WRF_OA fail to adequately capture the sub-daily variability of wind speeds. The correlations with observations fluctuate between 0.6 and 0.9 depending on the station, which is much less than those of temperatures. This is again primarily due to the fact that the variations in temperature fields are mostly dominated by large-scale and mesoscale weather conditions. Although the wind field variability is further modulated by small scale processes, including turbulence at all scale, the higher resolution brought by the WRF model does not bring obvious improvement to reproduce these processes in synchronicity with the reality. According to the RMSE values and standard deviation, ECMWF performs better than WRF for wind speed variability, even with assimilation and nudging in WRF_OA. The conclusion, therefore, is that the higher temporal resolution and higher spatial resolution of WRF output are not primary drivers for the improved description of the temporal variability, especially for the wind.
3.2 Spatial gradients
3.2.1 Horizontal gradients
The statistics discussed in section 3.1 do not demonstrate any clear advantage of the higher resolution WRF model with respect to the ECMWF forecasts. One may argue that this apparent absence of added value may hide the capability of the WRF model to reproduce fine scale variability as opposed to the mesoscale features that are well represented by ECMWF forecasts. Indeed, the horizontal gradients of surface temperature and wind speed between urban and rural areas are critical drivers of the atmospheric transport and must have a substantial impact on variations in CO2 concentrations. For a deeper analysis, we now focus on the horizontal gradients between pairs of observational sites using Taylor diagrams to assess the statistics of the hourly observation and modeled values for February 2015 (Figure 4). The gradients were computed as the difference between hourly observations taken at Montsouris (an urban station at the southern edge of Paris) and all other atmospheric weather station shown in Figure 2. Montsouris is chosen as a reference because it is the only meteorological station that is within the Paris city area. It is therefore the most likely to be significantly influenced by urban effects. The distance between Montsouris and the other sites used here varies between 8 km and 73 km.
Due to their coarse resolution, the ECMWF forecasts cannot reproduce the fine-scale horizontal gradients of temperature. The statistical results shown in Figure 4a clearly illustrate that the ECMWF model data have, as expected, little skills at reproducing the observed fine-scale horizontal gradients in temperature between the urban station and other stations. The same figure also indicates that the finer scale WRF model provide improvements over the values derived from ECMWF data. The correlations with observed gradients for both the ECMWF and WRF results are similar and relatively weak with the maximum value less than 0.7. There are more substantial differences between the corresponding RMSEs and standard deviations. Interestingly, the noOA version of the WRF model shows a large variability, although poor correlation with observed gradients. The variability is somewhat reduced by the OA version of the model with some positive impacts on correlation. As for the wind speed gradients, the WRF model again shows more variability than those derived from ECMWF data. This is likely related to the fact that the finer horizontal resolution of the WRF model is better at representing the local orography and associated circulation, and thus better captures small-scale features of the atmospheric flow and improves accuracy of the analysis data in proximity to the complex terrain. However, the associated variability is lower than that of the observations, and its correlation is weak at 0.5. These results indicate that the ability of WRF to reproduce small-scale wind gradients is no better than ECMWF.
Although the correlation remains weak in WRF_OA and ECMWF, the Taylor diagrams in Figure 4 illustrate the improvements provided by the WRF model. Indeed, the simulated variability using WRF increases realistically towards that of the observations. Note that the WRF vs. ECMWF improvement with regards to the horizontal gradients of wind speed and temperature is mainly associated with the differences between Montsouris and stations located in suburban or rural areas, whereas it is not noticeable in the differences between Montsouris and stations in urban areas. A plausible explanation for this is that the total area of Paris (approximately 105 km2) and the horizontal resolution of the innermost domain (3 km), limits the number of 3 km grids cells that cover the highly-urbanized area. As the distance from one urban station to another shrinks, and input parameters such as surface albedo, emissivity and moisture availability become homogeneous, the WRF meteorological model variables cannot vary substantially.
3.2.2 Vertical gradients
The vertical profile of temperature is the key parameter that controls the air stability and the vertical mixing (Crétat et al., 2012). This is particularly true for the lower atmospheric layers, where the vertical gradients control the mixing within the PBL. Together with the height of the PBL, the vertical gradient is therefore an essential parameter for the mixing of emitted CO2 within the low levels of the atmosphere. In this section, we evaluate the ability of the models to reproduce the temperature vertical gradients. We define a proxy for this gradient as the difference between a higher altitude level and the surface level derived from the WRF model and ECMWF forecasts. In this study, the observed data are provided by the twice daily (00 and 12 UTC) radiosonde launches at Trappes, with typical PBL heights of 1700 m and 650 m respectively. These typical PBL heights lead us to choose the temperatures at 500 m, 1000 m and 1500 m above ground level during the daytime, and 200 m, 300 m and 400 m above ground level at night as the levels to compute the temperature gradients through performing a vertical linear interpolation both for ECMWF forecasts and the WRF model.
Figure 5a shows that WRF_OA outperforms ECMWF at all analyzed levels, with good results in the representation of temperature gradients at 12 UTC. WRF_OA not only enhances the correlation (from 0.70, 0.94, 0.97 to 0.84, 0.98, 0.98 respectively) but also reduces bias and RMSE for the daytime data. The opposite was found to be true for WRF_noOA, especially at 500 m which exhibited the lowest correlation of 0.58. It is also worth noting that both the WRF model and ECMWF do a better job at reproducing the vertical temperature profile gradients with better accuracy at middle and upper levels than near the surface level. This is likely due to the fact that the surface temperature is deeply influenced by the surface topography and air mass stability. With the decrease of the surface friction along with heights above the ground, the atmosphere becomes more homogeneous, being easier to reproduce its behavior. In the nighttime case (00 UTC), WRF_OA agrees best with observations (Figure 5b). Although the Trappes radiosonde station is not located within the highly-urbanized area, it is affected by the combined effect of the urban heat island and local atmospheric diffusion. Hence, WRF_OA and, to a lesser extent, WRF_noOA, strongly reduces negative biases in vertical temperature profiles compared to ECMWF (0.52~0.86°C). They however provide products that are biased slightly warmer (0.48~0.72°C). These statistics reveal a weaker correlation for WRF_OA (0.85) than ECMWF (0.90) at 200 m. Nevertheless, the WRF_OA configuration does a better job at resolving the temperature gradients between the surface and both 300 m and 400 m.
3.3 Sensitivity of WRF results to model options
3.3.1 Impact on surface temperature and wind speed
The WRF model offers many options for the choice of physics schemes and nudging methods. In this section, we selected different schemes based on their prevalence in the scientific literature (see Tables 2 and 3) and assessed which combination leads to the best agreement with the observations. Two sets of sensitivity tests were carried out to investigate the behavior of the WRF model with different physics schemes and two nudging methods with or without the nudging to the observation data. Figure 6 shows the Taylor diagram of the simulation results as a function of the observations for the various WRF options. The times of sunrise and sunset in IdF for February is about 7:30~8:30 am (local time) and 17:30~18:30 pm (local time). As the local time is one hour ahead of UTC (UTC+1) during winter, we therefore distinguish the daytime (6–18 UTC) and nighttime (19–5 UTC) for a representative urban site at Montsouris, as well as a representative rural site at Nangis.
The results shown in Figure 6 indicate that different PBL schemes have substantial impacts on the wind speed with an obvious variation of the normalized standard deviation (e.g. 0.63~0.89 for daytime and 0.56~0.97 for nighttime of WRF_OA at the urban site), but little impact on the temperature. In general, the winds in the boundary layer can be decomposed into average wind speed, turbulence and fluctuation. Each PBL scheme uses different turbulence closure theories to depict the turbulent motions that have crucial influences on the fluctuations of wind fields. No clear differences in the temperature and wind speed arise from using different microphysics schemes or cumulus convection schemes according to their small amplitude of the variations of statistics in the Taylor plots. Note that in most cases, the phases of diurnal cycles of temperature and wind speed reproduced by the multi-physics ensembles are just being parallel shifted to higher or lower values, without obvious changes in amplitude (figure not shown).
Nudging also has a strong impact on the results. Including the multi-nudging (WRF_OA) makes the results much more in line with the observations than without (WRF_noOA). This is particularly true for temperature, but also for the wind speed over the urban sites. However, Figure 6 clearly shows that, when using the multi-nudging, the sensitivity to the physics schemes is obviously reduced. Strong nudging kills the model ability to develop its own dynamics features. As our ultimate purpose is to produce accurate meteorological data for atmospheric transport rather than focusing on the dynamic mechanism and flexibility of the WRF model itself, we nevertheless recommend the inclusion of the multi-nudging option for applications related to accurate transport of tracers.
Figure 7 shows the observed and simulated mean diurnal variations of the surface temperature and wind speed for sites at Montsouris and Nangis. It can be seen that the observed urban temperatures are higher than those in rural areas due to the impacts of anthropogenic activities producing heat as modeled by the UCM. The parameterization of the ECMWF model is too coarse to account for urban heat emissions and urban surface processes, so that this model produces temperatures in urban and rural areas that are almost the same. In contrast, the resolution of the WRF model, together with its urban scheme, makes it possible to represent the urban heat island effects. Over urban areas, the WRF_OA increases the modeled surface temperatures, in particular during nighttime, which result in a smaller amplitude of the modeled diurnal cycle.
The mean diurnal cycles of wind speeds are shown to be similar for WRF_OA and ECMWF in the rural site, however they are well reproduced by ECMWF in urban areas while clearly underestimated by WRF_OA. On the other hand, the diurnal cycles of ECMWF in urban and rural site are close to each other when compared to observations, since the wind field is rather smooth in the forecasts, which is a limitation when compared to the WRF model.
It is necessary to point out that the maximum temperature and wind speed by WRF_noOA are delayed by about 1 h and 1~3 h respectively when the UCM is used, which leads to obvious discrepancies between model simulations and observations, but the deviation become quite smaller for WRF_OA. A hypothesis for the delay of temperature is that the building shadow effect reduces daytime temperatures and radiations are kept reflected among buildings, lengthening the time for heat energy accumulating on the earth surface and then being transited into the atmosphere (Kim et al., 2013). Meanwhile, the simulated urban wind speeds decrease after adopting the UCM. This is because the increase of building heights enlarges the urban friction coefficient and the damping effect, therefore the atmospheric turbulence and momentum flux are being dragged in the urban canopy layer. Note that the UCM has little impact on the meteorology elements in rural areas, nevertheless.
3.3.2 Impact on PBL height
The PBL height and the realistic representation of boundary layer processes are key parameters to reproduce the vertical extent of mixing for parameters such as heat, moisture, and momentum, which are of vital importance for accurate simulations of the atmospheric transport and air quality within the PBL (Cohen et al., 2015). During nighttime, the atmospheric convection is minimal and there is no clear PBL. It is well known that atmospheric models have difficulties in reproducing mixing processes in such conditions (Shin and Hong, 2011). We therefore focus this work on examining the daytime PBL when the mixing layer is usually better developed and its height can be defined from the temperature profiles. Thus, an evaluation of PBL heights derived from the ECMWF forecasts and the WRF model with different physics schemes is performed against the land-based radiosonde data at Trappes at 12 UTC.
The calculation of the PBL height in the radiosonde observations, the ECMWF forecasts and the WRF PBL schemes are based on different methods. The PBL heights at Trappes are retrieved using radiosonde profiles of the virtual potential temperature with a critical value of bulk Richardson number (Rib), typically 0.21 (Kim et al., 2013). The ECMWF defines the top of the PBL as the level where the Rib reaches a different threshold of 0.25 (Troen and Mahrt, 1986). The PBL heights obtained directly from each of the four PBL schemes in WRF are also estimated by scheme specific formulations. A detailed description of the methods used to calculate the PBL height in the schemes can be found in Hong (2006), Nakanishi and Niino (2006, 2009), Pleim (2007a, b), Bougeault and Lacarrere (1989). For the purpose of a fair comparison, the 1.5-theta-increase method (Nielsen-Gammon et al., 2008) was commonly used as a criterion to diagnose PBL heights. This method defines PBL heights as the level at which the potential temperature first exceeds the minimum potential temperature within the boundary layer by 1.5 K.
The results presented in Figure 8 show that the WRF_noOA experiments have a strong sensitivity of modeled PBL height to the physical configuration chosen. For WRF_OA, the sensitivity to the physical scheme is smaller, as expected. The PBL height remains nevertheless sensitive to the PBL schemes, and much less based on the chosen radiation, cumulus convection and microphysics schemes. It is likely that each PBL scheme with its own mathematical depiction of the turbulence is responsible for the turbulent mixing in the lower troposphere (Crétat et al., 2012). Compared with the observation, both ECMWF forecasts and WRF_OA reproduce higher PBL heights due to the difficulty in simulating stable conditions (Seidel et al., 2012). ECMWF overestimates the PBL heights with a MBE of 131.68 m. The local MYNN3 scheme (PBL6) represents accurately the PBL height variations, with the highest correlation (0.87) and the smallest bias (84.81 m). The nonlocal YSU scheme (REF) has a similar bias (88.51 m) as the YSU scheme. The hybrid local-nonlocal ACM2 scheme (PBL7) and the local BouLac scheme (PBL8) give estimations by reproducing relatively deeper PBL heights with positive MBEs of 111.1 m and 124.07 m respectively during the daytime in WRF_OA. According to previous studies, this could be explained by the fact that the YSU scheme and the ACM2 scheme generally produces warmer and drier daytime PBLs (Hu et al., 2010). However, the ACM2 scheme may show a deeper mixing than the YSU scheme, since the eddy diffusivity of the ACM2 scheme is not zero at the PBL height due to the use of local wind shear and local Richardson number in the convective regime (Pleim, 2007b).
These analyses confirm that the PBL height is highly sensitive to the choice made in the WRF configuration, in particular when no data assimilation is used. Our results provide some guidance on the most appropriate schemes to reproduce the observations made at Trappes. However, its representativeness remains limited so that it would be inappropriate to draw firm conclusions at this point.
4 Summary and Conclusions
Comparisons of the ECMWF-analyzed and WRF-modeled meteorology with observations have been performed with focus on three atmospheric variables relevant to tracer transport and thus atmospheric inversions of CO2 fluxes. These atmospheric state variables were the air temperature, the wind and the PBL height over the IdF region. WRF with no assimilation and single nudging shows poor performance when compared to the both WRF with assimilation and multi-nudging, as well as ECMWF despite of the coarser spatial resolution and 3 h temporal cadence. As a consequence, the OA version of WRF is recommended for the interpretation of temperature, wind or atmospheric transport when an accurate description of reality is needed.
The analysis of the comparison results, that involve both ECMWF and WRF_OA simulated fields, can be used to address the paper main question regarding our ability to simulate the CO2 atmospheric transport in an urban environment. The wind speeds provided by ECMWF are, in general, closer to the observations than those provided by WRF, whereas the statistics for the wind directions tend to favor WRF_OA with respect to ECMWF. On the other hand, WRF provides improved temperature fields that are more representative of the observed values. This is particularly true for the spatial gradients of temperature between urban and rural/suburban sites induced by the heat island of the Paris city, which is not represented in the ECMWF outputs. The WRF model also provides a slightly better description of the station-to-station gradients of wind speed, which validates the expected ability to generate small-scale patterns. The ability of the ECMWF forecasts to reproduce the fine-scale horizontal gradients between urban and rural areas is insufficient for this application due to its coarse resolution, while appreciable improvements are achieved by WRF_OA showing higher standard deviation and lower RMSE. WRF_OA also outperforms ECMWF in the depiction of vertical temperature gradients by improving the statistics for all analyzed levels at 12 UTC. This work also shows that it provides a better representation of vertical temperature gradients by strongly reducing negative biases from ECMWF (0.52~0.86°C) to slightly warm biases (0.48~0.72°C) at 00 UTC.
The results of these sensitivity tests with different physics schemes provide an objective method to select the appropriate WRF model setup for the best representation of meteorology, and hence atmospheric transport over the IdF region. The choice of PBL schemes has a large impact on the wind speed and the PBL height. Moreover, the multi-nudging in WRF_OA could aid in eliminating biases, although it also reduces the spatial gradient and suppresses the physical sensitivity. This observation, suggest that the nudging strengths should be adjusted on the inner domain to optimize the balance between accuracy and variability according to the aims of desired application.
Although the ECMWF forecasts compare favorably with the WRF model for the surface wind speed, WRF results in this study appears to have more potential for reproducing the variability of wind speeds between urban and non-urban stations. Whether this slight advantage of WRF justifies the use of this complex model for CO2 transport, as opposed to the ECMWF fields which are readily available, remains to be tested. As long as CO2 is measured at peri-urban sites like in Paris (Bréon et al., 2015; Staufer et al., 2016), the use of WRF may not be advantageous. On the other hand, with more urban CO2 measurement sites being available in metropolitan areas (e.g., Feng et al., 2016), the ability to resolve urban meteorology may become more important. One may also argue that a high-resolution transport model coupled to a high-resolution inventory of emissions will outperform coarser resolution fields (Oda et al., 2017). We conclude that there is a marginal advantage of WRF over ECMWF for the desired application, which is nevertheless sufficient to motive additional testing of WRF vs. ECMWF with prescribed CO2 flux maps for comparing modeled CO2 concentrations with available observations. Further work should be devoted to the study of a long term (one year) modeling of the CO2 concentration in the Paris area with several parameterizations of the WRF model. The result of this modeling shall be compared to CO2 observations from seven surface stations located within Paris and the IdF region as well as line observations made by the laser GreenLITE system (Dobler et al., 2015).
Data Accessibility Statement
All model results from this study are available upon request from the corresponding author.
The supplemental files for this article can be found as follows: