INTRODUCTION
Air quality modeling is a growing interest in both urban and rural areas, particularly in the former due to the high density of population. Between 1950 and 2018 the world’s population underwent rapid urbanization, with the urban proportion rising from 30% in 1950 to 55% in 2018. By 2050, the urban fraction of the world’s population is projected to reach 68% (UN 2019). For this reason, urban modeling has become an important aspect not only for air quality but also for local meteorology assessments that require a better representation of the physical processes involved in urban environment mesoscale modeling. Thus, the Urban Canopy Model (UCM) with a single layer has been coupled to the Weather Research and Forecasting (WRF) model to improve the description of lower conditions of the urban boundary layer (UBL) and to provide more accurate forecasts for urban regions (Tewari et al. 2007). One key requirement for urban applications is that the WRF accurately represents city morphology influences on wind, temperature, and humidity in the atmospheric boundary layer and their collective influences on the atmospheric mesoscale motions (Chen et al. 2011).
Several urban parameterizations in the WRF have been used and coupled with the Noah Land Surface Model (LSM) (Chen and Dudhia 2001, Kong et al. 2021) through the urban percentage parameter (or urban fraction, Furb), which represents the proportion of impervious surfaces in the WRF sub-grid scale. Noah LSM provides surface-sensible, latent heat fluxes, and surface skin temperature as lower boundary conditions for coupled atmospheric models, whereas UCM provides the fluxes from anthropogenic surfaces. The relevance of urban-scale modeling is not only increasing with WRF but also with different modeling approaches; however, a small number of studies has been conducted in the field of environmental sciences (Wong et al. 2021). A bulk urban parameterization was included in WRF V2 (Chen et al. 2011). Then, the Single-Layer Urban Canopy Model (SLUCM) was developed, adopting infinitely long street canyons parameterized to represent urban geometry, but recognizing the three-dimensional nature of urban surfaces (Chen et al. 2011). Subsequently, the multi-layer Urban Canopy Model known as BEP for Building Effect Parameterization was developed by Martilli et al. (2002). Several modeling studies have been made using different kind of models and parameterizations (Ulpiani 2021). The study of Liao et al. (2014) evaluated diverse urban canopy schemes in the WRF model coupled with chemistry (WRF-Chem), including the SLUCM, the multi-layer UCM (BEP), and the multi-layer urban model with a Building Energy Model (BEM) that includes anthropogenic heat due to air conditioning (BEP + BEM). They determined that wind speed at 10 m decreases when urban canopy schemes are coupled to the model (Liao et al. 2014), which is important for air quality issues.
Several modeling studies have been carried out in the central region of Mexico; however, only a small number report the use of urban parameterizations. Ozone (O3) production and response under different meteorological conditions was studied in the Mexico City Metropolitan Area (MCMA) during a field campaign in 2003, precising that O3 production depends on pollutant emissions and meteorological conditions (Lei et al. 2008). Correspondingly, in the MILAGRO field campaign numerous studies were performed in the Mexico City region in order to characterize chemical and physical transformations, and to assess impacts of exported pollutants from urban areas on regional and global air quality, ecosystems and climate. The effects of spatial resolution over trace gases and aerosols suggest that simulations at a 3-km horizontal grid resolution adequately reproduce the overall transport and mixing of trace gases and aerosols downwind of Mexico City, whereas a large grid resolution (such as 75 km) is insufficient to represent both local emissions and the terrain impact on meteorological fields affecting dispersion and transport of trace gases and aerosols and their sub-grid variability (Qian et al. 2010).
Zhang et.al. (2009) assessed meteorological variables and criteria pollutants simulated with the WRF-Chem model against measurements from the Red Automática de Monitoreo Atmosférico (Automatic Air Quality Monitoring Network, RAMA), finding there is a cold bias in temperature surface between 1 and 2 ºC from daytime to nighttime, being probably one of the reasons that the 3-km grid resolution is insufficient to resolve small scale circulations in urban environments and the urban infrastructure effects, which were not included in this work (Zhang et al. 2009). This is supported in the study made by Tie et al. (2010) to illustrate the effects of horizontal spatial resolutions on concentrations of O3 and its precursors in a mega city, whose results show that model resolution has important effects on calculations of air pollutants dispersion in urban areas and photochemical O3 production, due to calculated meteorological conditions, spatial distribution of emissions and the non-linearity of photochemical O3 production (Tie et al. 2010). A previous study found that O3 concentrations in a city plume are very sensitive to meteorological and ambient chemical conditions, including extra-urban scale transport winds, wind field vertical structure, and mixing processes (Tie et al. 2009); therefore, we can expect that including an urban parameterization and increasing grid resolution in the model could lead to an improve in local meteorology and air quality simulations.
According to Zhang et al. (2009), characterizing the impacts of urban pollutants requires detailed modeling studies. Mexico City’s air quality and meteorological forecast (AQFS-CDMX), run by the Secretaría del Medio Ambiente de la Ciudad de México (Mexico City Environment Secretariat, SEDEMA), incorporates such characteristics by connecting a meteorological model (WRF-ARW) and a regional emissions model with the Community Multi-scale Air Quality (CMAQ) chemical transport model and the Noah LSM with single-layer Urban Canopy Model (UCM), in order to estimate the air quality for the next 24 h in the MCMA (SEDEMA 2017) and support the activation and suspension of atmospheric environmental contingencies (SEDEMA 2019). In addition, to characterize Mexico City’s urban meteorological conditions, the WRF model was coupled with the multi-layer UCM (BEP) and the WUDAPT level 0 (Martilli et al. 2016), and then configured to run four nested domains of 13.5, 4.5, 1.5, and 0.5-km grid resolutions. A comparison of the urban integration (WRFu) against the traditional (WRF) simulations and observations was made for both dry and wet seasons in 2016, showing that WRFu better captures the influence of the city over changes in temperature, wind speed, and planetary boundary layer (PBL) in comparison with the WRF simulations. WRFu also showed a greater urban heat island effect during daytime and simulations were closer to observations in both seasons (Fernández 2017). A similar study was developed in the city of Barcelona to evaluate the performance of urban schemes integrated in the WRF, including BEP. In average, highest correlations with observations for temperature and wind speed were obtained by BEM in comparison with WRF (Ribeiro et al. 2021).
The multi-layer UCM parameterization has already been implemented for the central region of Mexico; however, there are some gaps that still need to be addressed, such as the implementation of the Chem module in WRFu in order to assess air quality simulations, which is the main purpose of this study. In this work, we describe the implementation of the urban WRF, including a multi-layer urban canopy model (BEP) with chemistry to simulate air quality pollutants. In order to select an optimal configuration for this model, sensitivity assessments were conducted by increasing the grid resolution from 1 m to 0.5 km and by changing the initial and boundary conditions to evaluate the impact over local meteorology. In addition, a comparison between the urban WRF-Chem (WRFu-Chem) and WRF-Chem simulations without urban parameterization and local observations is presented for carbon monoxide and O3 as model performance evaluation.
MATERIALS AND METHODS
Study area description
The MCMA is defined by the conurbation of 76 municipalities from three entities: 16 of Mexico City (CDMX), 59 of the State of Mexico, and one of the state of Hidalgo (SEDATU 2018). The study area comprises some municipalities of the State of Mexico, Morelos and the entire CDMX, as shown in figure 1, where urban and rural regions of the latter, surrounded by urban areas of the State of Mexico (12 M inhabitants), are also displayed. The CDMX (19.0482-19.5928º N, 99.3649-98.9403º W, with elevations between 2240-3930 masl [INEGI 2022]) has a total population of 8 985 339 inhabitants. It sits in a closed basin surrounded by mountains, with the nearest mountain range located to the southeast, south and west, and the highest mountains found to the east about 40 km from the city center. The study area is displayed in figure 1. In 2014, the CDMX annual mean temperature was 17.7 ºC, with an annual mean maximum of 23.8 ºC and a minimum of 11.5 ºC. The mean precipitation for this year was 655.9 mm (CONAGUA 2014).
Air quality observations data
In order to evaluate the model performance, air quality data from the Sistema de Monitoreo Atmosférico (Atmospheric Monitoring System, SIMAT) of Mexico City were used. Criteria pollutants (O3, SO2, NO2, CO, PM10, PM2.5, are measured by RAMA in 34 stations throughout the city, and surface meteorological parameters (temperature, relative humidity, wind direction, and wind speed) are measured continuously in 26 different sites of the Red de Meteorología y Radiación Solar (Meteorology and Solar Radiation Network, REDMET). These data are available in hourly average format and can be accessed within the air quality webpage of Mexico City (GCDMX 2019). The air quality monitoring stations (AQMS) used to evaluate the model performance are indicated in figure 1. It should be noted that stations Santa Fe (SFE), Cuajimalpa (CUA), and Tlahuac (TAH) are in the edge of the mountain range. On the other hand, Atizapán (ATI), San Agustín (SAG), FES Acatlán (FAC) and Villa de las Flores (VIF) are located in semi-urban regions, mostly surrounded by rural areas, while Acolman (ACO) is located in a rural area. The rest of the AQMS are located in urban areas, as shown in table I.
Site ID | Latitude | Longitude | Height (m) | Site name | Observations |
ACO | 19.6355 | -98.9120 | 2198 | Acolman | Rural |
CUA | 19.3653 | -99.2917 | 2704 | Cuajimalpa | Urban |
FAC | 19.4824 | -99.2435 | 2299 | FES Acatlán | Semiurban |
MER | 19.4246 | -99.1196 | 2245 | Merced | Urban |
PED | 19.3251 | -99.2041 | 2326 | Pedregal | Urban |
SAG | 19.5329 | -99.0303 | 2241 | San Agustín | Urban |
ATI | 19.5769 | -99.2541 | 2341 | Atizapán de Zaragoza | Semiurban |
CAM | 19.4684 | -99.1698 | 2233 | Camarones | Urban |
XAL | 19.5259 | -99.0824 | 2160 | Xalostoc | Urban |
SFE | 19.3573 | -99.2628 | 2599 | Santa Fe | Semiurban |
TAH | 19.2464 | -99.0105 | 2297 | Tlahuac | Semiurban |
VIF | 19.6582 | -99.0966 | 2242 | Villa de las Flores | Semiurban |
Description of the urban WRF-Chem modeling system
The Weather Research and Forecasting (WRF) model is a next-generation mesoscale numerical weather prediction system designed both for atmospheric research and operational numerical weather predictions (NWP), and is suitable for applications on air quality modeling, among others (Skamarock et al. 2005). The WRF model incorporates several physics options. In this study the BEP model, which is considered the most high-level modeling in WRF for urban applications, has been included (Chen et al. 2011), alongside WRF v. 3.2. The core of the urban WRF-Chem is shown figure 2, in which red boxes correspond to the built-in modules in the WRFu-Chem modeling system and the main application. These modules will be described briefly in the following sections.
Building Effect Parameterization (BEP) model
The BEP model comprises a multi-layer UCM developed by Martilli et al. (2002), which allows direct interaction with the PBL and recognizes the 3D nature of urban surfaces considering vertical effects (walls) and horizontal surfaces (streets and roofs) on the momentum (drag-force approach), turbulent kinetic energy (TKE), and potential temperature that substantially impacts the thermodynamic structure of the urban roughness sub-layer, and hence the lower part of the UBL. The BEP scheme is operational with Noah LSM and has been coupled with two turbulence schemes: Bougeault and Lacarrère (1989) and Mellor-Yamada-Janjić (Janjić 1994), making it able to simulate some of the most observed features in urban atmosphere, such as the nocturnal urban heat island (UHI) and the elevated inversion layer above the city (Chen et al. 2011).
A key requirement for urban applications in the WRF is to accurately capture influences of cities on wind, temperature, and humidity in the atmospheric boundary layer, as well as their collective influences on the atmospheric mesoscale motions. Chen et al. (2011) established that to take full advantage of BEP, it is necessary to have high vertical resolution close to the ground. In order to consider influences of buildings over meteorological parameters we use 51 vertical levels in this application with a height of about 20 m for the first layer and a 0.5-km grid resolution for the finest domain.
Urban land use
To run the multi-layer urban canopy model version of BEP embedded in WRFV3.2, data from the World Urban Database and Access Portal Tools (WUDAPT) level 0 was used. WUDAPT level 0 is based on the LCZ classification scheme (Steward and Oke 2012, Stewart et al. 2014). Figure 3 shows a summary of the methodology followed for implementing the WRFu-Chem modeling system (Martilli et al. 2016).
To use WUDAPT level 0 as input for BEP in WRFV3.2 it is necessary to follow up the following steps:
Extension of the number of urban classes from 3 to 10 as stated in Steward and Oke (2012).
Using the LCZ map created by WUDAPT level 0 as a foundation, we conducted an update of the land use field (LU_INDEX).
Modification of the URBPARM.TBL table to define supplementary parameters such as urban fraction, building heights, heat capacity, etc., for each urban class (Martilli et al. 2016).
Various studies have demonstrated that changes in land use impact meteorological parameters and pollutants, and these changes vary according to the urban canopy schemes used (Liao et al. 2014). In this study, the default USGS land use configuration of LANDUSE.TBL, structured in 24 categories, was used in WPS. Then, urban categories were extended to 10 different categories according to the method used by Martilli et al. (2016) in WUDAPT level 0 based on the LCZ classification. In this way, the new categories of urban classification correspond to the categories from item 31 to 40.
WRF-Chem module
The WRF-Chem model simulates emission, transport, mixing, and chemical transformation of trace gases and aerosols simultaneously with meteorology. It is used for investigation of regional-scale air quality, field program analysis, and cloud-scale interactions between clouds and chemistry (Grell et al. 2005, NCAR 2020a). WRF-Chem is coupled with several modules for working with anthropogenic and biogenic emissions, allowing several choices for gas-phase chemical mechanisms like RADM2, RACM, SAPRC99, among others, and also for photolysis and aerosol schemes (Peckham 2010). Detailed information about the WRF-chem model is presented in Grell et al.( 2005) and Fast et al. (2006).
Emissions inventory
The emissions used for chemical simulation corresponds to the National Emissions Inventory, which is considered a hybrid inventory since it consists of different years emissions as follows: 2013 for point source emissions; 2014 for source area emissions; and 2015 for mobile emissions. The Inventory only comprises the emissions of Mexico. Regarding time distribution it considers hourly emissions data, as well as UTM and latitude-longitude in degrees coordinates (INECC 2017). The emissions information contains organic and inorganic species according to RADM2 mechanism used for modeling the atmospheric chemistry (Stockwell et al. 1990). To produce emission input data for WRFu-Chem simulations, these emissions were interpolated using a mass conservative method1 into the three domains throughout the initial condition files wrfinput_d01, wrfinput_d02 and wrfinput_d03 generated by the real.exe program. For details on the spatial and temporal distribution and speciation refers to (García-Reynoso et al. 2018).
WRFu-Chem model configuration
Meteorological input data
The meteorological input data used as boundary conditions to initialize the urban modeling system through the WRF Preprocessing System (WPS) are the National Center for Environmental Prediction (NCEP) final (FNL) Operational Global Analysis data, with 1-degree by 1-degree spatial resolution and prepared operationally every six hours of temporal resolution (NCAR 2020) and the NCEP North American Mesoscale (NAM) Forecast System data with 12 km spatial resolution every 6 hours (NCAR 2020).
Study case
During February 19th 2014 to February 27th 2014, measurements of O3, PM10 and PM2.5 with the nearest monitoring sites (CUA, SFE) were compared, showing a high correlation in tendency of these pollutants (Noyola 2014). The WRFu-Chem modeling system was set up to simulate this case. This modeling system will be helpful to study further episodes of high pollutant concentrations in the region, for network monitoring design and also to evaluate spatial representativeness of air quality monitoring stations that will be presented in a separate study. The simulations run continuously over 150 hours, starting from February-20-2014 00:00:00 UTC and ending at February-26-2014 06:00:00 UTC. Only the simulations from days 21-00:00:00 to 25-23:00:00 of local time (UTC-06:00) were used to evaluate the model performance against measurements of the monitoring stations of the RAMA.
Domain configuration
For this study, three nested domains were set up for simulations. The domains are shown in figure 4 as follows: the coarse domain (d01) has 79 x 79 cells with a 4.5 km x 4.5 km grid resolution; the second domain (d02) has 1.5 km x 1.5 km grid resolution; and the third and the finest domain (d03) has 151 x 163 cells with 0.5 km x 0.5 km grid resolution that comprise the entire CDMX and part of Mexico State where monitoring stations are also installed and used for the model evaluation. The map in figure 1 illustrates the finest domain d03 in detail. Due to computational reasons, usually the model resolution cannot be higher enough to consider the effects induced by buildings on wind and dispersion and thus, the turbulent flow around the urban obstacles cannot be resolved explicitly. In order to consider the atmospheric processes inside the Urban Canopy Layer induced by the urban morphology, a high resolution of 0.5 km was set for the finest domain (Santiago and Martilli 2010).
Physics
In order to simulate the PBL effects in WRFu-Chem model over wind speed and temperature, the option bl_pbl_physics = 8 was used. This PBL scheme was designed for using with BEP urban model integrated in WRF model. To select between Noah Land Surface Model with no urban physics (WRF) and the multi-layer urban canopy model urban parameterization embedded in BEP (WRFu), the parameterizations described in table II were used. Additional physic parameterization configured for running the model is described in table III.
Scenario | Parameterization | Observation |
WRF | sf_urban_physics = 0 | No urban physics is used |
WRFu | sf_urban_physics = 2 | BEP model is used |
Physics variable name | Option | Scheme description |
mp_physics | 4 | WRF Single-Moment 5-class scheme: A simple efficient scheme with ice and snow processes suitable for mesoscale grid sizes. Allows for mixed-phase processes and super-cooled water. |
ra_lw_physics | 1 | RRTM scheme: Rapid Radiative Transfer Model. An accurate scheme using lookup tables for efficiency. Accounts for multiple bands, trace gases and microphysics species (Mlawer et al. 1997). |
ra_sw_physics | 2 | Goddard shortwave: Two-stream multi-band scheme with ozone from climatology and cloud effects. |
sf_sfclay_physics | 2 | Eta similarity: Used in Eta model. Based in Monin-Obukhov with Zilitinkevich thermal roughness length and standard similarity functions from look-up tables. |
sf_surface_physics | 2 | Noah Land Surface Model: Unified NCEP/NCAR/AFWA scheme with soil temperature and moisture in four layers, fractional snow cover and frozen soil physics. |
bl_pbl_physics | 8 | BouLac: Bougeault-Lacarrère PBL. A TKE-prediction option (Bougeault and Lacarrère 1989) |
sf_urban_physics | 0 | Noah Land Surface Model: Unified NCEP/NCAR/AFWA scheme with soil temperature and moisture in four layers, fractional snow cover and frozen soil physics. No urban physics. |
2 | BEP: Building Environment Parameterization: Multi-layer urban canopy model that allows for buildings higher than the lowest model levels. | |
cu_physics | 5 | Grell 3d ensemble cumulus scheme. Scheme for higher resolution domains allowing for subsidence in neighboring columns. |
*For the planetary boundary layer scheme (BL_PBL_PHYSICS), option 8 (Bougeault-Lacarrère) designed to be used with Building Effect Parameterization (BEP) was selected.
RESULTS AND DISCUSSION
In order to identify the best configuration for WRFu-Chem, six different experiments were conducted, including two with input data for initial and boundary conditions, changes in grid resolution with two and three nested domains, and urban physics parameterization. The results obtained are presented in this section.
Sensitivity to initial boundary conditions
The experiments configured to assess the sensitivity of the system to initial boundary conditions are described in table IV.
No. | Experiment | Grid Resolution | Domains | Sf_urban _physics | Initial and Boundary conditions |
1 | WRF _ fnl | 0.5 km | 3 | 0 | FNL |
2 | WRFu_fnl | 2 | |||
3 | WRF _ nam | 0 | NAM | ||
4 | WRFu_nam | 2 |
*All experiments have three domains for simulation and a 0.5 × 0.5 km grid resolution for the finest domain (d03). Urban surface (sf_urban_physics) option is switched between 0 and 2 to select the Noah Land Surface Model (LSM) with no urban physics (WRF) and the Multi-Layer Urban Canopy Model (WRFU) with Building Effect Parameterization (BEP), respectively.
Temperature
To evaluate the initial boundary conditions impact on local weather, NCEP FNL (NCAR 2020b) and NAM (NCAR 2020c) input data were used for the same configuration of WRF and WRFu. Figure 5 shows a time series of 2-m temperature simulations with WRF and WRFu using NAM and FNL data, and their comparison with observations (in black) of different monitoring sites. The WRF model underestimates high values of observed temperature; in contrast, WRFu simulations are close to observed high values in most of the stations. Nonetheless, in rural areas (ACO) and in zones with low or moderate urbanization (such as FAC, PED, SFE, and TAH) temperature observations tend to be lower than in urban areas. In these cases, both WRF and WRFu fail to reproduce the lowest observed values, though WRFu is more accurate (Fig. 5). It means that urban parameterization enhances simulations, increasing the reproducibility of observations; however, a better LCZ classification is needed in the WRFu model to improve local weather simulations in these zones.
To illustrate the differences between NAM and FLN input data in temperature simulations, a statistical comparison for nine sites was made. Bias variability is reduced when using urban parameterization against no urban physic simulations. Also, a reduction of bias error was observed when using NCEP FNL data in both WRF and WRFu in comparison with NAM data. Table V shows the verification measures of simulations made with WRF, when no urban physics is used. In this case, 2-m temperature simulations showed a bias error from -0.87 to 0.3 ºC for NAM input data, while for FNL input data the bias error was from -0.95 to 0.21 ºC. The bias error for WRFu simulations (Table VI) was from -1.24 to 0.60 ºC for NAM data whereas for simulations with FNL data it was from -1.30 to 0.47 ºC. (In tables V-VIII and XII-XIII numbers in bold denote minimum and maximum values for each variable).
T2 WRF Site | BIAS (ºC) | RMSE (ºC) | Correlation | Ia | ||||
NAM | FNL | NAM | FNL | NAM | FNL | NAM | FNL | |
ACO | 0.30 | 0.21 | 1.71 | 1.83 | 0.96 | 0.95 | 0.96 | 0.96 |
CUA | -0.87 | -0.95 | 2.28 | 2.06 | 0.91 | 0.92 | 0.88 | 0.91 |
FAC | -0.50 | -0.71 | 3.24 | 2.91 | 0.92 | 0.94 | 0.86 | 0.90 |
MER | -0.61 | -0.83 | 1.42 | 1.36 | 0.96 | 0.96 | 0.96 | 0.96 |
PED | 0.05 | -0.19 | 2.60 | 2.29 | 0.90 | 0.91 | 0.88 | 0.92 |
SAG | -0.11 | -0.20 | 2.03 | 1.78 | 0.95 | 0.96 | 0.93 | 0.95 |
SFE | -0.64 | -0.77 | 2.30 | 2.16 | 0.89 | 0.90 | 0.89 | 0.91 |
TAH | 0.02 | -0.18 | 2.31 | 1.97 | 0.95 | 0.95 | 0.91 | 0.94 |
VIF | 0.03 | -0.13 | 1.71 | 1.51 | 0.96 | 0.97 | 0.95 | 0.97 |
T2 WRFu Site | BIAS (ºC) | RMSE (ºC) | Correlation | Ia | ||||
NAM | FNL | NAM | FNL | NAM | FNL | NAM | FNL | |
ACO | 0.11 | -0.01 | 1.80 | 1.89 | 0.95 | 0.94 | 0.96 | 0.96 |
CUA | 0.11 | 0.14 | 1.63 | 1.50 | 0.93 | 0.93 | 0.95 | 0.96 |
FAC | -0.37 | -0.43 | 2.38 | 2.12 | 0.94 | 0.95 | 0.94 | 0.96 |
MER | -1.24 | -1.30 | 1.79 | 1.95 | 0.94 | 0.93 | 0.94 | 0.94 |
PED | 0.41 | 0.29 | 1.83 | 1.59 | 0.94 | 0.95 | 0.95 | 0.97 |
SAG | -0.04 | 0.04 | 1.45 | 1.42 | 0.95 | 0.96 | 0.97 | 0.98 |
SFE | -0.01 | 0.01 | 1.75 | 1.69 | 0.92 | 0.92 | 0.94 | 0.95 |
TAH | -1.03 | -1.02 | 2.37 | 2.02 | 0.95 | 0.96 | 0.91 | 0.94 |
VIF | 0.60 | 0.47 | 1.92 | 1.97 | 0.92 | 0.92 | 0.96 | 0.96 |
High correlation values were obtained in temperature for all simulations. WRF correlations were between 0.89 and 0.97 with a variability of 0.08, whereas correlations in WRFu simulations were from 0.92 to 0.96 with a variability of 0.04, which means a reduction between sites. In this situation, WRFu simulations with NAM and FNL data are very similar. In addition, the RMSE is reduced in most sites when WRFu is used and the variability from WRF to WRFu is also reduced in 49% for NAM and 55% for FNL data. A graphical comparison of these metrics is shown in the Taylor Diagram on figure 6, where results of WRFu simulations are closer to the observed values. The index of agreement (Ia) also presented higher values, between 0.86 and 0.97 for WRF and above 0.91 for WRFu. In this last situation, the Ia variability is 0.06 for NAM and 0.04 for FNL, which shows in general better results for temperature.
Wind speed
Wind speed simulations with WRF and WRFu underestimate observations; however, a lower bias error for each site is observed when NCEP FNL data are used as initial boundary conditions. The bias error for WRF simulations was from -0.95 to -0.08 m/s for NAM data and from -0.78 to -0.02 m/s for FNL data, as shown in table VII.
WRF Site | BIAS (m/s) | RMSE (m/s) | Correlation | Ia | ||||
NAM | FNL | NAM | FNL | NAM | FNL | NAM | FNL | |
ACO | -0.08 | -0.02 | 1.45 | 1.56 | 0.41 | 0.45 | 0.61 | 0.62 |
CUA | -0.45 | -0.27 | 0.97 | 0.97 | 0.25 | 0.29 | 0.5 | 0.52 |
FAC | -0.53 | -0.45 | 1.12 | 1.25 | 0.5 | 0.42 | 0.67 | 0.63 |
MER | -0.77 | -0.6 | 1.21 | 1.22 | 0.54 | 0.59 | 0.61 | 0.65 |
PED | -0.82 | -0.76 | 1.16 | 1.19 | 0.45 | 0.46 | 0.53 | 0.55 |
SAG | -0.15 | -0.02 | 0.95 | 1.03 | 0.72 | 0.7 | 0.78 | 0.76 |
SFE | -0.95 | -0.78 | 1.33 | 1.23 | 0.24 | 0.38 | 0.44 | 0.53 |
TAH | -0.81 | -0.59 | 1.42 | 1.36 | 0.35 | 0.34 | 0.56 | 0.57 |
VIF | -0.56 | -0.39 | 1.21 | 1.16 | 0.34 | 0.39 | 0.58 | 0.61 |
On the other hand, the bias error for simulations with urban parameterization (WRFu) and NAM data was from -0.70 to -0.04 m/s, whereas for WRFu with FNL data it was from -0.52 to 0.25 m/s (Table VIII). The maximum bias errors for all simulations were in the Santa Fe (SFE) site, which is in the mountain range side.
WRFu Site | BIAS (m/s) | RMSE (m/s) | Correlation | Ia | ||||
NAM | FNL | NAM | FNL | NAM | FNL | NAM | FNL | |
ACO | -0.13 | 0.25 | 1.49 | 1.66 | 0.41 | 0.43 | 0.6 | 0.59 |
CUA | -0.29 | -0.13 | 0.82 | 0.89 | 0.52 | 0.43 | 0.65 | 0.6 |
FAC | -0.31 | -0.24 | 0.81 | 0.87 | 0.66 | 0.61 | 0.78 | 0.77 |
MER | -0.6 | -0.42 | 1.07 | 1.09 | 0.6 | 0.65 | 0.68 | 0.71 |
PED | -0.53 | -0.49 | 0.9 | 0.9 | 0.58 | 0.61 | 0.66 | 0.68 |
SAG | -0.04 | 0.07 | 0.82 | 0.89 | 0.67 | 0.68 | 0.8 | 0.78 |
SFE | -0.7 | -0.52 | 1.22 | 1.22 | 0.35 | 0.3 | 0.52 | 0.51 |
TAH | -0.66 | -0.29 | 1.26 | 1.2 | 0.49 | 0.46 | 0.66 | 0.67 |
VIF | -0.23 | -0.2 | 0.9 | 0.88 | 0.49 | 0.54 | 0.71 | 0.72 |
A wind-speed Taylor diagram for different sites is shown in figure 7. The WRF and WRFu simulations with NAM and FNL initial boundary conditions are graphically and statistically compared. WRFu simulations exhibited greater correlation and lower RMSE than WRF simulations in the stations compared. Furthermore, the index of agreement was also computed and higher values were observed for WRFu simulations compared to WRF. Correlation variability for WRF was 0.48 for NAM and 0.42 for FNL data. For WRFu, the variability of correlation was 0.32 for NAM and 0.38 for FNL data. Variability of RMSE for WRF was 0.50 m/s for NAM data and 0.59 m/s for FNL, whereas variability of RMSE for NAM and FNL simulations with WRFu was 0.68 and 0.79 m/s, respectively. The relatively large variability in RMSE is likely attributed to variations in urban classification as perceived by the WRFu model at each individual station and variations in measurement height.
The previous results suggest that the WRFu model improves simulations of local weather conditions using either NAM or FNL data. In general, when using FNL data in WRFu simulations the bias error is reduced in most sites where the model is compared against observations. In this way and considering the limitations of geographical coverage of NAM data, we decided to use FNL data for the next simulations.
Influence of grid resolution
We evaluated the spatial resolution influence of the model over the study area. For this, four experiments were configured as described in table IX. The experiments with ID 1 and 2 have three domains with a grid resolution of 0.5 by 0.5 km for the finest domain (d03); and experiments 3 and 4 have two domains and 1 by 1 km grid resolution for the finest domain (d02). All experiments were run using NCEP FNL input data. The flag sf_urban_physics = 0 was set in order to select between Noah LSF with no urban physics (WRF), while sf_urban_physics = 2 was set for using the multi-layer urban canopy model BEP (WRFu).
No. | Run Experiment | Spatial Resolution | Domains | sf_urban_physics | Initial and Boundary conditions |
1 | WRF - 0.5 km | 0.5 km | 3 | 0 | FNL |
2 | WRFu-0.5 km | 0.5 km | 3 | 2 | |
3 | WRF - 1 km | 1 km | 2 | 0 | |
4 | WRFu-1 km | 1 km | 2 | 2 |
*Experiments 1 and 2 have three domains for simulation and a 0.5 × 0.5 km grid resolution for the finest domain (d03). Experiments 3 and 4 have two domains for simulation and a 1 × 1 km grid resolution for the finest domain (d02). Urban surface (sf_urban_physics) option is switched between 0 and 2 to select the Noah Land Surface Model (LSM) with no urban physics (WRF) and the Multi-Layer Urban Canopy Model (WRFU) with Building Effect Parameterization (BEP), respectively.
In figure 8 a diurnal variation of 2 m temperature observations and simulations of 1 and 0.5 km spatial resolution are shown. It is observed that WRF simulations tend to underestimate high values of temperature, whereas simulations with WRFu increase temperature, reaching high observed values; however, simulations at 1-km resolution overestimate the observations.
A reduction in bias variability is observed for both WRF and WRFu simulations when grid resolution is increased from 1 to 0.5 km, as observed in table X. Variability for WRF simulations is 1.9 and 1.2 ºC for 1 and 0.5 km, respectively, and 2.6 and 1.8 ºC for 1 and 0.5 km, respectively, in the case of WRFu simulations, which correspond to a bias reduction of 39 and 33%, respectively. Even though RMSE has relatively high values, variability was also reduced about 17 and 39% for WRF and WRFu, respectively, when grid resolution was increased from 1 to 0.5 km. In terms of correlation and index of agreement, similar values were obtained, resulting in values above 0.9, but regarding urban parameterization the observed WRFu variability is lower. The correlation coefficient is lower in WRFu because in some stations (MER, SAG) temperature from WRFu is higher than WRF during early hours and lower during afternoon hours, which indicates that the heat capacity used in those areas has to be reviewed.
No. | Run experiment | BIAS (ºC) | RMSE (ºC) | Correlation | Ia | ||||
Max | Min | Max | Min | Max | Min | Max | Min | ||
1 | WRF-0.5 km | 0.22 | -0.95 | 2.91 | 1.36 | 0.97 | 0.90 | 0.97 | 0.90 |
2 | WRFu-0.5 km | 0.47 | -1.30 | 2.12 | 1.42 | 0.96 | 0.92 | 0.98 | 0.94 |
3 | WRF-1 km | 1.50 | -0.39 | 2.81 | 0.93 | 0.98 | 0.90 | 0.98 | 0.90 |
4 | WRFu-1 km | 1.93 | -0.71 | 2.52 | 1.37 | 0.97 | 0.92 | 0.98 | 0.94 |
*Values correspond only to the maximum and minimum obtained from the nine sites evaluated.
Wind speed simulations generally have lower values than observations; however, observed peak values are overestimated in all experiments. Higher peak values were obtained in simulations with 1-km grids, whereas simulations with 0.5-km grid resolutions are closer to the maximum value of observations, causing a 18% reduction of in WRFu bias variability and 8% in WRF. RMSE variability is greater in 0.5 km than 1 km simulations; nonetheless, specific results by site show that WRFu simulations have lower RMSE values than WRF. On the contrary, correlation and index of agreement are generally greater in WRFu than in WRF. The maximum and minimum values of nine sites compared in the experiments are shown in table XI.
No. | Run Experiment | BIAS (m/s) | RMSE (m/s) | Correlation | Ia | |||||
Max | Min | Max | Min | Max | Min | Max | Min | |||
1 | WRF-0.5 km | -0.02 | -0.78 | 1.56 | 0.97 | 0.70 | 0.29 | 0.76 | 0.52 | |
2 | WRFu-0.5 km | 0.25 | -0.52 | 1.66 | 0.87 | 0.68 | 0.30 | 0.78 | 0.51 | |
3 | WRF-1 km | -0.02 | -0.86 | 1.50 | 1.12 | 0.60 | 0.31 | 0.75 | 0.50 | |
4 | WRFu-1 km | 0.24 | -0.71 | 1.36 | 0.89 | 0.68 | 0.28 | 0.82 | 0.48 |
*Values correspond only to the maximum and minimum obtained from the nine sites evaluated.
In conclusion, increasing the resolution in simulations has a positive impact in model performance since bias and RMSE are reduced and the correlation coefficient and index of agreement increase. On the other hand, simulations of urban parameterization with WRFu + BEP show a better performance in 2-m temperature (T2) and wind speed. Thus, increasing resolution and including urban details in the WRFu model gives a better approach to observations; nevertheless, an improved definition of the LCZ classification must be included in order to increase the WRFu model performance and reduce bias variability, which deserves further investigation.
Effects of urban parameterization on air quality
Experiments 1 and 2 of table IX were used to evaluate the effects of urban parameterization on air quality modeling. Experiment 1 corresponds to WRF simulations without urban parameterization and experiment 2 to WRFu, including BEP. These simulations were run at a 0.5-km grid resolution to better capture influences of local meteorology and urban canopy on air quality. In this study, results of a primary pollutant and tracer, carbon monoxide, and a secondary photochemical pollutant, O3, are presented.
Ozone
Time series of O3 observations at different sites are shown in figure 9, alongside their comparison with WRF and WRFu simulations. In both cases the temporal profile is well simulated; however, WRFu has slightly lower values than WRF. The simulations consistently underestimate peak observed values during daylight at all sites. This is primarily because it is the first time that a high-resolution emissions inventory is used, and improved meteorology data are employed, which produces a change in ambient concentrations of primary pollutants such as carbon monoxide (CO) and nitric oxide (NO). An increase in CO indicates an increment in NO concentrations, then a titration of O3 occurs in the WRFu and a lower peak is obtained. A sensitivity analysis for nitrogen oxides (NOx) and volatile organic compounds (VOC) emissions is required to improve O3 concentrations and it deserves more research.
O3 simulations and observations from different stations are graphically and statistically compared in the Taylor diagram shown in figure 10. In this case, no large differences are observed; however, WRFu results show some improvement, according to metrics on table XII. Correlation coefficients vary from 0.76 to 0.89 for WRF and from 0.79 to 0.88 for WRFu. The index of agreement fluctuates between 0.67 and 0.85 in both experiments. Not only correlation values are high but also the index of agreement, even though RMSE and biases are high. Bias variability is 17.7 ppb for WRF and 17.2 ppb for WRFu; and RMSE variability is 11.4 and 13.7 ppb, respectively. The variability is obtained from the differences between maximum and minimum values highlighted in bold in each column.
O3 Site | BIAS (ppb) | RMSE (ppb) | Correlation | Ia | ||||
WRF | WRFu | WRF | WRFu | WRF | WRFu | WRF | WRFu | |
ACO | -5.66 | -8.41 | 19.57 | 19.49 | 0.77 | 0.82 | 0.72 | 0.75 |
CUA | -12.36 | -14.72 | 28.81 | 28.99 | 0.8 | 0.83 | 0.7 | 0.71 |
FAC | -11.51 | -15.42 | 29.2 | 30.11 | 0.81 | 0.83 | 0.74 | 0.74 |
ATI | -10.42 | -13.21 | 29.77 | 30.39 | 0.88 | 0.87 | 0.73 | 0.73 |
CAM | 2.33 | -2.29 | 21.71 | 21.57 | 0.89 | 0.88 | 0.85 | 0.85 |
XAL | 5.37 | 1.79 | 18.4 | 16.68 | 0.83 | 0.86 | 0.82 | 0.85 |
SFE | -10.54 | -13.05 | 26.2 | 26.32 | 0.83 | 0.86 | 0.77 | 0.78 |
TAH | -7.75 | -9.84 | 25.99 | 25.39 | 0.76 | 0.79 | 0.74 | 0.77 |
VIF | -7.41 | -11.55 | 28.86 | 29.06 | 0.79 | 0.84 | 0.67 | 0.69 |
Carbon monoxide
Highest observed CO values (above 2 ppm, with peak values up to 4 ppm) were obtained north of the city. In general, peak values are underestimated by simulations at all sites, except in CAM, where WRFu presented a higher peak than observations in day 23, suggesting an extreme episode of ambient CO concentrations north of the MCMA, since concentrations in WRFu and observations in stations ATI, FAC, VIF, and XAL also reproduced that behavior during the mentioned day, as observed in figure 11. It is also observed that CO ambient concentrations are larger in WRFu than in WRF at all sites. CO concentrations in the south region are lower in both observations and simulations, with values below 2 ppm, except for the rural station ACO located to the northeast, where observations are less than 1 ppm and simulations of WRF and WRFu present similar patterns with values near to 0 ppm, underestimating observations. This could suggest that the LCZ classification used for the simulation period better represent the local conditions of the urban area north of the MCMA, as compared to the south region and the rural site ACO. On the other hand, since local meteorological conditions change when using the urban parameterization, a revised emissions inventory must be used.
The diurnal variation in figure 12 shows measured and simulated CO concentrations as average daily profiles of the nine sites compared in figure 11. This figure shows that, in general, WRFu simulations had a better performance when compared to observations, since peak values are better represented than in WRF simulations.
CO metrics are shown in table XIII. Bias variability is similar in both simulations when compared to observations; however, bias errors are lower in WRFu. Sites located north of the CDMX present larger reductions in bias with respect to WRF simulations. The bias error varies from -0.45 to 0.02 in WRFu and from 0.59 to -0.11 in WRF. RMSE values are also lower in WRFu than in WRF, except for CAM station, where RMSE increased in 24%. The RMSE for WRF was between 0.19-0.78 and between 0.16-0.66 for WRFu, with variabilities of 0.59 and 0.50 for WRF and WRFu, respectively. On the other hand, the correlation coefficient for WRF was from 0.19 to 0.75 and from 0.19 to 0.70 for WRFu, with an index of agreement from 0.38 to 0.50 for WRF and from 0.4 to 0.72 for WRFu. In this case, values obtained by WRFu were higher than those obtained by WRF when compared against observations, showing evidence of a better performance of the WRFu simulation when using urban parameterization with BEP.
CO Site | BIAS (ppm) | RMSE (ppm) | Correlation | Ia | ||||
WRF | WRFu | WRF | WRFu | WRF | WRFu | WRF | WRFu | |
ACO | -0.11 | -0.08 | 0.19 | 0.16 | 0.6 | 0.68 | 0.47 | 0.58 |
CUA | -0.48 | -0.42 | 0.53 | 0.48 | 0.32 | 0.36 | 0.38 | 0.4 |
FAC | -0.52 | -0.29 | 0.78 | 0.66 | 0.52 | 0.43 | 0.44 | 0.58 |
ATI | -0.59 | -0.45 | 0.72 | 0.61 | 0.48 | 0.42 | 0.41 | 0.46 |
CAM | -0.22 | 0.02 | 0.41 | 0.51 | 0.52 | 0.38 | 0.5 | 0.62 |
XAL | -0.57 | -0.34 | 0.73 | 0.52 | 0.75 | 0.7 | 0.5 | 0.72 |
SFE | -0.31 | -0.24 | 0.39 | 0.36 | 0.19 | 0.19 | 0.44 | 0.46 |
TAH | -0.38 | -0.33 | 0.49 | 0.43 | 0.55 | 0.56 | 0.48 | 0.52 |
VIF | -0.53 | -0.39 | 0.66 | 0.55 | 0.6 | 0.52 | 0.46 | 0.52 |
The model performance evaluation for CO is represented by a Taylor diagram in figure 13, showing that WRFu simulations are closer to observations compared to WRF simulations.
For further explanation on the verification measures, please refer to section S1 in the supplementary material.
CONCLUSIONS
Using urban parameterization and FNL data reduces the bias error in temperature and bias variability. RMSE values are lower in most of the sites with WRFu simulations, reducing RMSE variability in 55% as compared to WRF simulations. Correspondingly, WRFu simulations of temperature show a better agreement with observations even when these are underestimated, and peak values are reached in most of the sites. Nonetheless, in rural areas (ACO) and zones with low or moderate urbanization such as FAC, PED, SFE and TAH, observations of temperature tend to be lower than in urban areas. Observations of temperature often tend to be lower in non-urban areas compared to urban areas. This highlights the role of urban parameterization in improving simulations and increasing the accuracy of observations.
Increasing the resolution to 0.5 km in simulations has a positive impact in local weather since BIAS and RMSE are reduced and the correlation coefficient and index of agreement increase. WRFu temperature simulation at a 1-km grid resolution overestimates observations; however, bias variability and RMSE variability are reduced in 33 and 39%, respectively, when grid resolution is increased to 0.5 km and peak observation values are also reached, while WRF simulations tend to underestimate high values of temperature in both cases.
The highest CO concentrations are measured north of the MCMA, while lower values are measured south of the city, with the exception of the ACO rural site located northeast, which presented the lowest values. These patterns were reproduced by the WRFu-Chem model. Although simulations underestimate peak observation values at all sites, CO ambient concentrations are larger in WRFu-Chem than in WRF-Chem, indicating a less dilution rate and showing a better agreement with peak observation values. This suggests that the LCZ classification used for the simulation period better represents the urban area local conditions north of the MCMA than in the south and in the ACO rural site.
On the other hand, ozone simulations have differences in peak values. WRFu-Chem has slightly lower values than WRF-Chem. Peak observation values occurring during daylight are underestimated by simulations at all sites due to an increase in wind speed that reduces accumulation of O3 precursors such as NO and VOC. In contrast, lower observation values during nighttime are overestimated, suggesting that simulated concentrations are influenced from O3 regional transport. Bias and RMSE are relatively high; however, index of agreement values between 0.67 and 0.85 were obtained at all sites.
Improving the model performance in meteorological variables using urban parameterization with WRFu leads to a change in environmental concentrations. This experience deserves further investigation on local urban physical parameters and evaluation of emission inventories. Moreover, the LCZ classification used in this study should be actualized to improve local weather and reducing bias and RMSE errors to increase the WRFu-Chem model performance over the study area.