Comparison of tropospheric NO2 columns from MAX-DOAS retrievals and regional air quality model simulations

Abstract

Multi-axis differential optical absorption spectroscopy (MAX-DOAS) tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ column retrievals from four European measurement stations are compared to simulations from five regional air quality models which contribute to the European regional ensemble forecasts and reanalyses of the operational Copernicus Atmosphere Monitoring Service (CAMS). Compared to other observational data usually applied for regional model evaluation, MAX-DOAS data are closer to the regional model data in terms of horizontal and vertical resolution, and multiple measurements are available during daylight, so that, for example, diurnal cycles of trace gases can be investigated.

In general, there is good agreement between simulated and retrieved ${\mathrm{NO}}_{\mathrm{2}}$ column values for individual MAX-DOAS measurements with correlations between 35 % and 70 % for individual models and 45 % to 75 % for the ensemble median for tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ vertical column densities (VCDs), indicating that emissions, transport and tropospheric chemistry of ${\mathrm{NO}}_x$ are on average well simulated. However, large differences are found for individual pollution plumes observed by MAX-DOAS. Most of the models overestimate seasonal cycles for the majority of MAX-DOAS sites investigated. At the urban stations, weekly cycles are reproduced well, but the decrease towards the weekend is underestimated and diurnal cycles are overall not well represented. In particular, simulated morning rush hour peaks are not confirmed by MAX-DOAS retrievals, and models fail to reproduce observed changes in diurnal cycles for weekdays versus weekends. The results of this study show that future model development needs to concentrate on improving representation of diurnal cycles and associated temporal scalings.

1 Introduction

Given the influence of ${\mathrm{NO}}_x$ on air quality and climate through effects on radiation, it is of high environmental and scientific interest to accurately observe and simulate spatial distribution and time evolution of ${\mathrm{NO}}_{\mathrm{2}}$ concentrations in the troposphere. Simulating ${\mathrm{NO}}_{\mathrm{2}}$ is a challenge for numerical models as it is chemically very active and hence depends on many factors including, for example, cloud cover which affects photolysis of this trace gas. Moreover, representation of ${\mathrm{NO}}_x$ emissions adds a large uncertainty to the model output.

Multi-axis differential optical absorption spectroscopy (MAX-DOAS) (e.g., ) measurements have been used to investigate air pollution in many studies, including the FORMAT campaign in northern Italy (), the Cabauw Intercomparison campaign for Nitrogen Dioxide measuring Instruments (CINDI) in the Netherlands , campaigns in Canada (), China (e.g., ) and during ship-borne measurements ().

MAX-DOAS observations of atmospheric composition are performed by taking measurements of the scattered sunlight at different elevation and sometimes also azimuthal angles. Depending on the viewing angle and solar position, the light path through the atmosphere is different, with the observation in the zenith direction usually providing the shortest light path through the lower troposphere. Therefore, using observations in low-elevation angles as measurement intensity and zenith measurements as reference intensity, the total amount of molecules of a certain species along the light path difference (zenith subtracted from non-zenith measurement), so-called differential slant column densities, can be determined using Beer–Lambert law. These can be inverted to tropospheric columns and lower altitude tropospheric profiles by radiative transfer modeling and optimal estimation techniques.

A large number of studies applied MAX-DOAS data for satellite validation (e.g., ) but up to now, comparisons to regional air quality model simulations of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ have, to our knowledge, only been carried out by and . Several studies compared regional air quality model simulations to satellite data (e.g., ), although satellite data are usually only available at much coarser time steps compared to regional model data. In this respect, the advantage of MAX-DOAS retrievals compared to satellite retrievals is the finer resolution in time. Moreover, several studies compared in situ ${\mathrm{NO}}_{\mathrm{2}}$ data to regional model results (e.g., ), although in situ data usually refer to a specific location (point measurements), whereas regional model results are available for a specific horizontal grid resolution and area depending on the model setup. As MAX-DOAS data represent a larger volume of air, they are much better suited for investigating performance of regional models than in situ data. According to , the horizontal averaging volume of MAX-DOAS data depends on aerosol loading, wavelength and viewing direction, and ranges from a few kilometers in the polluted boundary layer up to 80 km from the top of a mountain under clean air conditions. Another advantage of MAX-DOAS measurements is their ability to observe several pollution-related species at the same time (e.g., ${\mathrm{NO}}_{\mathrm{2}}$, HCHO, CHOCHO, ${\mathrm{SO}}_{\mathrm{2}}$, aerosols and potentially also ${\mathrm{O}}_{\mathrm{3}}$) and to provide ${\mathrm{NO}}_{\mathrm{2}}$ data which are virtually free of interferences from other species or nitrogen compounds such as ${\mathrm{NO}}_y$ (${\mathrm{NO}}_x$ and other oxidized nitrogen species). In contrast to ${\mathrm{NO}}_{\mathrm{2}}$, ${\mathrm{NO}}_x$ cannot be retrieved from MAX-DOAS measurements directly, so these measurements are of more interest for air quality than for atmospheric chemistry studies. Vertical profiles of trace gases can be retrieved from MAX-DOAS measurements, which is another advantage for model comparison studies.

In the present study, MAX-DOAS measurements are compared to regional air quality model simulations to evaluate the quality of the simulations and to identify possible weaknesses of the models. Parts of this approach are already applied within scientific reports of the operational Copernicus Atmosphere Monitoring Service (CAMS; http://atmosphere.copernicus.eu/, last access: 7 February 2020) (see, e.g., ) but mainly to model results provided on eight output levels only, which introduces uncertainty to comparison results. CAMS is the operational follow-up of the former GEMS (Global and regional Earth-system Monitoring using Satellite and in-situ data) and three succeeding MACC (Monitoring Atmospheric Composition and Climate) projects. The global component of CAMS extends weather services of the ECMWF (European Centre for Medium-Range Weather Forecasts) with simulations of atmospheric trace gases and aerosols, while operational air quality forecasts and analyses for Europe are provided at much finer resolution through the regional component. Hourly ${\mathrm{NO}}_{\mathrm{2}}$ vertical column densities (VCDs) from six different regional model runs based on five models which are used within CAMS will be compared to MAX-DOAS measurements from one rural and three urban European stations: Bremen (operated by IUP-Bremen), De Bilt (operated by KNMI), Uccle and OHP (Observatoire de Haute-Provence) (the latter two operated by Royal Belgian Institute for Space Aeronomy – BIRA-IASB). The locations of the stations are plotted on top of mean ${\mathrm{NO}}_{\mathrm{2}}$ tropospheric columns from the Ozone Monitoring Instrument (OMI) satellite observations for February 2011 as well as on a map of anthropogenic ${\mathrm{NO}}_x$ emissions used by the models in Fig.  as an indicator of pollution levels in these and surrounding regions. The spatial distribution of ${\mathrm{NO}}_x$ emissions agrees well with pollution hotspots and cleaner areas identified by OMI.

Figure 1

Maps of (a) average tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs (molec cm${}^{-\mathrm{2}}$) observed by OMI for February 2011 and (b) TNO/MACC-II anthropogenic ${\mathrm{NO}}_x$ emissions (g m${}^{-\mathrm{2}}$ yr${}^{-\mathrm{1}}$) over Europe. The locations of MAX-DOAS measurement sites investigated in this study are marked by black crosses on the maps. The satellite data have been gridded to 0.1${}^{\circ }$ lat $\times$ 0.1${}^{\circ }$ long; the resolution of the emission database is 0.125${}^{\circ }$ lat $\times$ 0.0625${}^{\circ }$ long.

[Figure omitted. See PDF]

This paper starts with an overview of regional model and MAX-DOAS data (Sect. ), followed by a description of the comparison method (Sect. ). Results are described and discussed in Sect. . Finally, a summary and conclusions are given.

2 Data basis

2.1 Regional air quality model simulations

CHIMERE , LOTOS-EUROS (LOng Term Ozone Simulation – EURopean Operational Smog) , EMEP MSC-W (European Monitoring and Evaluation Programme Meteorological Synthesizing Centre – West) , MOCAGE (Model Of atmospheric Chemistry At larGE scale) () and SILAM (System for Integrated modeLling of Atmospheric coMposition) () contributed to the European regional ensemble forecasts and reanalyses of the former MACC projects and are currently used within CAMS. These models have been used in many studies for investigating atmospheric composition on a regional scale (e.g., ).

All of these models use ECMWF-IFS and MACC reanalysis data as meteorological and chemical input data and boundary conditions, respectively. Anthropogenic emissions are taken from the MACC emissions database , and the Global Fire Assimilation System (GFAS) is used to account for fire emissions. The input to these models is thus consistent and differences in model results are mainly due to different representation of chemical processes, advection, convection, turbulent mixing, wet and dry deposition .

The model runs investigated in the present study were performed by different European institutions and are based on different horizontal and vertical grid spacings and chemistry schemes (see Table  for further details). Apart from SILAM, the models were run without chemical data assimilation. The SILAM simulations included assimilation of surface observations of ${\mathrm{NO}}_{\mathrm{2}}$ as described in .

Table 1

Overview of regional air quality model simulations.

Model	Institution	Grid spacing	Number of vertical levels;	Chemistry scheme
		(zonal $\times$ meridional)	model top
CHIMERE	LISA-CNRS/UPEC/UPD	0.25${}^{\circ }$ $\times$ 0.25${}^{\circ }$	8;	MELCHIOR II
	INERIS	($\sim \mathrm{18}\times \mathrm{28}$ km${}^{\mathrm{2}}$)	500 hPa
EMEP-MACCEVA	MetNo	0.25${}^{\circ }$ $\times$ 0.125${}^{\circ }$	20;	EMEP-EmChem09soa
		($\sim \mathrm{18}\times \mathrm{14}$ km${}^{\mathrm{2}}$)	100 hPa	(;
				)
EMEP	MetNo	$\mathrm{50}\times \mathrm{50}$ km${}^{\mathrm{2}}$	20;	EMEP-EmChem09soa
			100 hPa	(;
				)
LOTOS-EUROS	TNO	0.125${}^{\circ }$ $\times$ 0.0625${}^{\circ }$	3;	TNO CBM-IV
		($\sim \mathrm{9}\times \mathrm{7}$ km${}^{\mathrm{2}}$)	$\sim \mathrm{3.5}$ km	(;
				)
MOCAGE	CNRS Météo-France	0.2${}^{\circ }$ $\times$ 0.2${}^{\circ }$	47;	Troposphere:
		($\sim \mathrm{15}\times \mathrm{22}$ km${}^{\mathrm{2}}$)	5 hPa	RACM
				Stratosphere:
				REPROBUS
SILAM	FMI	years 2010/2011:	9 (2010), 8 (2011–2012);	DMAT
		0.2${}^{\circ }$ $\times$ 0.2${}^{\circ }$	6.725 km (2010),
		($\sim \mathrm{15}\times \mathrm{22}$ km${}^{\mathrm{2}}$)	6.7 km (2011–2012)
		year 2012:
		0.15${}^{\circ }$ $\times$ 0.15${}^{\circ }$
		($\sim \mathrm{11}\times \mathrm{17}$ km${}^{\mathrm{2}}$)

Two different sets of EMEP model runs are investigated in this study. The first one uses the same setup as the other regional models described above and is termed EMEP-MACCEVA in the following. EVA (validated assessments for air quality in Europe) was a subproject of MACC dedicated to the development and implementation of operational yearly production of European air quality assessment reports (https://policy.atmosphere.copernicus.eu/Reports.php, last access: 10 February 2020). The second set of simulations (called EMEP in the following) uses the same setup as in the EMEP status reports (see http://www.emep.int, last access: 7 February 2020) for each year based on the EMEP subdomain, ECMWF-IFS as meteorological driver, EMEP emissions, Fire INventory from NCAR version 1.0 (FINNv1; ), initial conditions described by for the years 2010–2011 and for 2012 and climatological boundary conditions described by .

According to , chemistry transport models in general account for seasonal, daily and diurnal emission changes by applying average time profiles given for different energy sectors and regions to totals of annual emissions across the model domain. Temporal emission patterns used by the regional air quality models listed above are country and SNAP (Selected Nomenclature for Sources of Air Pollution) sector dependent and are based on . A list of the SNAP sectors is given by . Moreover, different vertical emission profiles are applied for each regional model. These are described in more detail by for EMEP and CHIMERE, for SILAM and for LOTOS-EUROS. For MOCAGE, emissions are injected into the five lowest model levels using a hyperbolic decay.

More details on specific model setups and scores with respect to surface observations can be found in and in the model specification/validation dossiers which are available online at https://atmosphere.copernicus.eu/regional-services (last access: 7 February 2020).

2.2 MAX-DOAS retrievals

This study makes use of MAX-DOAS measurements from four European stations: Bremen (Germany), De Bilt (the Netherlands), Uccle (Belgium) and OHP (France). Characteristics of the data available from the stations, such as exact location and time period of retrievals investigated here, are briefly summarized in Table  and will be described below.

Table 2

Overview of MAX-DOAS station data.

Station	Lat, long;	Institution	Time period	Type	Retrieved	Number of layers;	Additional data
location	height (m a.s.l.)				quantity	layer top (km)
De Bilt,	52.1${}^{\circ }$ N, 5.18${}^{\circ }$ E;	KNMI	Mar 2011–Dec 2012	urban	column	12;	wind data (in situ)
Netherlands	$\sim \mathrm{23}$ m					4.0 km
Bremen,	53.11${}^{\circ }$ N, 8.86${}^{\circ }$ E;	IUP-UB	Jan 2011–Dec 2011	urban	column	81;	wind data
Germany	21 m				profile	4.025 km	(in situ data from
							airport weather station
							$\sim \mathrm{9}$ km southwards
							at 53.05${}^{\circ }$ N, 8.79${}^{\circ }$ E)
Uccle,	50.8${}^{\circ }$ N, 4.32${}^{\circ }$ E;	BIRA-IASB	Jan 2012–Dec 2012	urban	column	13;	wind data (in situ)
Belgium	120 m				profile	3.5 km	clouds from MAX-DOAS
OHP,	43.92${}^{\circ }$ N, 5.7${}^{\circ }$ E;	BIRA-IASB	Jan 2010–Dec 2010	rural	column	13;	wind data (in situ)
France	650 m				profile	3.5 km

For Bremen, Uccle and OHP, ${\mathrm{NO}}_{\mathrm{2}}$ slant column densities (SCDs) are obtained by a DOAS analysis for a specific wavelength window using a series of low-elevation angles as measurement intensity and zenith measurements as reference intensity. Cross sections of different trace gases are accounted for in the retrieval. Resulting SCDs of ${\mathrm{NO}}_{\mathrm{2}}$ and ${\mathrm{O}}_{\mathrm{4}}$ are then used as input for a radiative transfer model which is a two-step approach. First, an aerosol extinction profile is retrieved by comparing the measured ${\mathrm{O}}_{\mathrm{4}}$ SCDs to ${\mathrm{O}}_{\mathrm{4}}$ SCDs simulated by the radiative transfer models SCIATRAN for Bremen and bePRO for Uccle and OHP. In the second step, the derived aerosol extinction, measured ${\mathrm{NO}}_{\mathrm{2}}$ SCDs and an a priori ${\mathrm{NO}}_{\mathrm{2}}$ profile are used to retrieve the ${\mathrm{NO}}_{\mathrm{2}}$ profile of interest. This is an inverse problem solved by the means of the optimal estimation method .

De Bilt (52.10${}^{\circ }$ N, 5.18${}^{\circ }$ E; see Fig. ) is the hometown of KNMI and located just outside the city of Utrecht. The De Bilt experimental research site is surrounded by local and regional roads, with a lot of traffic which can affect regional air quality significantly. According to , it can also be affected by pollution sources which are located more far away in the Rotterdam region to the southwest, Amsterdam to the northwest and the German Ruhr region to the southeast of De Bilt. The MAX-DOAS instrument operated at De Bilt is a commercial system obtained from Hoffmann Messtechnik GmbH. It has an Ocean Optics spectrograph, diffraction grating and a charge-coupled device (CCD) detector. It operates at a wavelength range of 400–600 nm. The pointing direction of the instrument is 80${}^{\circ }$ (east to northeast); the wavelength window of the DOAS fit for ${\mathrm{NO}}_{\mathrm{2}}$ is 425–490 nm. Wavelength calibration and slit-function width are determined using a high-resolution solar spectrum. Cross sections of ${\mathrm{O}}_{\mathrm{3}}$, ${\mathrm{NO}}_{\mathrm{2}}$, ${\mathrm{O}}_{\mathrm{4}}$, ${\mathrm{H}}_{\mathrm{2}}\mathrm{O}$ and a pseudo cross section accounting for the Ring effect are applied. The choice of fitting parameters complies with the standards agreed by the MAX-DOAS community, following from homogenization efforts within, e.g., CINDI, GEOMON, NORS and QA4ECV as much as possible. Air mass factor (AMF) calculations are performed by the DAK radiative transfer model (Doubling Adding KNMI; ). A priori profiles of ${\mathrm{NO}}_{\mathrm{2}}$ are based on a block profile with ${\mathrm{NO}}_{\mathrm{2}}$ present in the boundary layer; boundary layer heights were taken from a climatology based on ECMWF data. For De Bilt, averaging kernels refer to the altitude-dependent (or box-)differential AMFs divided by the total differential AMF. The differential AMF is derived at a specific altitude by simulating the radiance with and without an added partial column of ${\mathrm{NO}}_{\mathrm{2}}$ at this altitude with the DAK model. ${\mathrm{NO}}_{\mathrm{2}}$ tropospheric columns are retrieved from the 30${}^{\circ }$ elevation angles; ${\mathrm{NO}}_{\mathrm{2}}$ profiles were not available at the time of the study.

The IUP-Bremen MAX-DOAS instrument consists of an outdoor telescope unit collecting light in different directions, and an indoor grating spectrometer (Shamrock 163 equipped with an Andor LOT257U CCD with $\mathrm{2048}\times \mathrm{512}$ pixels) covering a wavelength interval from 430 to 516 nm at a resolution of approximately 0.7 nm. Both components are connected via an optical fiber bundle which simplifies handling and overcomes polarization effects. The telescope unit is installed at an altitude of approximately 20 m above ground level (m a.g.l.) at the roof of the Institute of Environmental Physics building (53.11${}^{\circ }$ N, 8.86${}^{\circ }$ E) at the University of Bremen which is located to the northeast of the city center. The azimuthal pointing direction is northwest, which means that some of the measured pollution peaks are due to the exhaust of an industrial area, predominantly a steel plant, as well as a nearby highway. However, over longer time periods, the retrievals are dominated by pollution from the city center. ${\mathrm{NO}}_{\mathrm{2}}$ SCDs are derived by DOAS analysis using a fitting window of 450–497 nm. Cross sections of ${\mathrm{O}}_{\mathrm{3}}$, ${\mathrm{NO}}_{\mathrm{2}}$, ${\mathrm{O}}_{\mathrm{4}}$, ${\mathrm{H}}_{\mathrm{2}}\mathrm{O}$ and a pseudo cross section accounting for the ring effect are accounted for in the fit. Profiles of ${\mathrm{NO}}_{\mathrm{2}}$ are derived from SCDs applying the BRemian Advanced MAX-DOAS retrieval algorithm (BREAM) which incorporates SCIATRAN radiative transfer simulations. An ${\mathrm{NO}}_{\mathrm{2}}$ a priori profile which is constant with altitude is assumed and iterated in the retrieval. Detailed information about the profile retrieval is given by and .

BIRA-IASB has been operating a MAX-DOAS instrument at OHP (43.92${}^{\circ }$ N, 5.7${}^{\circ }$ E) since 2005. OHP is a background remote site in the south of France located on a plateau at 650 m above sea level (m a.s.l.), temporarily affected by transport of pollution from regional sources (e.g., from the petrochemical plants of Étang de Berre close to Marseille in the southwest) and the Po Valley (Italy) to the northeast of the station. The MAX-DOAS instrument, which points towards the SSW direction, consists of a grating spectrometer Jobin-Yvon TRIAX 180 (1800 grooves mm${}^{-\mathrm{1}}$) covering the 330–390 nm wavelength range coupled to a thermoelectrically cooled ($-\mathrm{40}$ ${}^{\circ }$C) Hamamatsu CCD detector (1024 pixels). ${\mathrm{NO}}_{\mathrm{2}}$ SCDs are obtained by applying the DOAS technique to a 364–384 nm wavelength interval, taking into account spectral signatures of ${\mathrm{O}}_{\mathrm{3}}$, ${\mathrm{O}}_{\mathrm{4}}$, the ring effect and ${\mathrm{NO}}_{\mathrm{2}}$. In Uccle (50.8${}^{\circ }$ N, 4.32${}^{\circ }$ E), which is located southwest of the Brussels city center, a mini-MAX-DOAS from Hoffmann Messtechnik GmbH covering the 290–435 nm wavelength range has been operated by BIRA-IASB since 2011. The instrument points north to northeast towards the city center. ${\mathrm{NO}}_{\mathrm{2}}$ SCDs are retrieved in a 407–432 nm wavelength interval including the same spectral signatures as for OHP. It should be noted that a sequential zenith reference spectrum has been implemented in order to minimize the impact of changes in shift and resolution due to temperature instabilities. The DOAS fit for ${\mathrm{NO}}_{\mathrm{2}}$ has also been improved by introducing pseudo-absorber cross sections derived from principal component analysis of residuals on days affected by large thermal instabilities. This approach allows for a better correction of fast-changing slit-function variations, resulting in more stable residuals and therefore more realistic random uncertainty estimates. For ${\mathrm{NO}}_{\mathrm{2}}$ vertical profile retrievals at both stations, the bePRO radiative transfer code is used. ${\mathrm{NO}}_{\mathrm{2}}$ profiles are retrieved at 420 nm for Uccle and 372 nm for OHP. For ${\mathrm{NO}}_{\mathrm{2}}$ vertical profile retrievals, exponentially decreasing a priori profiles have been constructed, based on an estimation of ${\mathrm{NO}}_{\mathrm{2}}$ vertical column densities derived from the so-called geometrical approximation () and using scaling heights of 0.5 and 1 km for OHP and Uccle, respectively. A priori and measurement-uncertainty covariance matrices are constructed as by with adopted correlation lengths of 0.05 and covariance scaling values of 0.5 and 0.35 for Uccle and OHP, respectively. A more detailed description of the model and trace gas profile retrievals can be found in . Although there has not been formal side-by-side operation of both instruments for verification purpose, a good overall agreement has been obtained between the mini-DOAS and other BIRA research-grade spectrometers similar to the one operated at OHP, e.g., like during the CINDI campaign .

Previous studies (e.g., ) have shown that the typical error on MAX-DOAS retrieved VCDs is around 20 %, including uncertainties related to the optimal estimation method, trace gas cross sections and aerosol retrievals, and can be higher for sites with low trace gas concentrations like OHP or due to instrumental conditions. Moreover, the uncertainty of the retrieval is increased in cloudy conditions in particular as clouds are not included in the MAX-DOAS forward calculations.

For Uccle, information on cloud conditions was retrieved according to the method by , which is based on analysis of the MAX-DOAS retrievals but not applied for results shown in the present study. No cloud flags are available for Bremen, De Bilt and OHP. For Bremen, MAX-DOAS retrievals were filtered for broken cloud situations by comparing radiative forward calculations of ${\mathrm{O}}_{\mathrm{4}}$ to retrieved ${\mathrm{O}}_{\mathrm{4}}$ columns and removing scans with correlations $<$ 0.6. For the other stations, the following data were not considered: (De Bilt) retrievals with fitting residuals $>$ $\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ and ${\mathrm{NO}}_{\mathrm{2}}$ or ${\mathrm{O}}_{\mathrm{4}}$ fitting errors $>$ 50 %, (Uccle/OHP) retrievals with ${\mathrm{NO}}_{\mathrm{2}}$ fitting root mean square error (RMSE) $>$ 50 % and degrees of freedom for signal $<$ 1. Clouds may alter MAX-DOAS retrievals in several ways depending on the viewing direction at which they are present. This is described in detail by . In addition to the direct effect of clouds on the measurements, clouds also affect photolysis rates and hence ${\mathrm{NO}}_x$ chemistry and NO-to-${\mathrm{NO}}_{\mathrm{2}}$ partitioning, which can impact on tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ columns and profiles retrieved under cloudy weather conditions. The influence of clouds on comparison results is complex and regarded as a separate topic for future studies.

2.3 Wind measurements

In order to investigate the ability of the models to reproduce transport of ${\mathrm{NO}}_{\mathrm{2}}$ towards the stations, the MAX-DOAS data described above are complemented by meteorological in situ station data of wind speed and wind direction. Wind data for Bremen were provided by the German Weather Service/Deutscher Wetterdienst through their website at http://www.dwd.de (last access: 7 February 2020). The weather station in Bremen is located at the main airport, approximately 9 km southwards of the MAX-DOAS station. This may result in differences to the actual wind direction and wind speed at time and location of the MAX-DOAS retrievals. Wind data for OHP were taken from the weather station at the observatory and downloaded from the corresponding website at http://pc-meteo.obs-hp.fr/intervalle.php (last access: 10 February 2020). Wind speed and direction measurements at Uccle are performed using a commercial rugged wind sensor from Young (model 05103) and were provided by BIRA-IASB (Francois Hendrick, personal communication, 2016). For De Bilt, wind measurements (within 300 m from the MAX-DOAS instrument) carried out by KNMI were downloaded from https://www.knmi.nl/nederland-nu/klimatologie/uurgegevens (last access: 7 February 2020).

3 Methodology for regional model evaluation

The sensitivity of MAX-DOAS retrievals is largest in the boundary layer, which needs to be taken into account when comparing MAX-DOAS retrievals to model simulated values. This is achieved here by applying column-averaging kernels (AVKs) to the model data prior to comparison. The AVKs are part of the MAX-DOAS profiling output and represent the sensitivity of the retrieved column to the amount of ${\mathrm{NO}}_{\mathrm{2}}$ at different altitudes.

In this study, model VCDs are derived by two different methods in order to test the influence of AVKs on the data analysis. Non-AVK-weighted model VCDs are calculated by summing up ${\mathrm{NO}}_{\mathrm{2}}$ partial columns (VCD${}_i$) over all $N$ model levels in the vertical. AVK-weighted model VCDs are calculated by applying column AVKs of the retrievals to model ${\mathrm{NO}}_{\mathrm{2}}$ partial columns before summing up ${\mathrm{NO}}_{\mathrm{2}}$ partial columns in the vertical:

1$${\mathrm{VCD}}_{\text{AVK-weighted}}^{\mathrm{model}}=\sum _{i=\mathrm{1}}^{N_{\mathrm{obs}}}{\mathrm{AVK}}_i\cdot {\mathrm{VCD}}_i^{\mathrm{model}},$$ where $N_{\mathrm{obs}}$ is the number of measurement altitudes. ${\mathrm{NO}}_{\mathrm{2}}$ partial columns were derived at measurement altitudes from the model output prior to the application of AVKs.

Only those model values closest to the measurement time are used for the comparisons. As the model output is given in hourly time steps, the maximum possible time difference between measurements and simulations is 30 min. VCDs are calculated separately for each model and constitute the basis for calculating ensemble mean values.

Following studies by, e.g., , , and , the present paper focuses on results of the model ensemble, i.e., the median of individual model results of a given quantity. As an even number of six different model runs (based on five different models) constitutes the model ensemble in the present study, the median is calculated by ordering the six different model values (e.g., for seasonal cycles, these values refer to the average of individual model runs for each month) in terms of magnitude and taking the average of the two middle numbers. An exception is OHP, as MOCAGE data are not available for this station, so the median refers to the middle number here. In addition, results from separate models are briefly discussed and shown in the main part of the paper to understand characteristics of the model ensemble output. It is beyond the scope of this study to describe the performance of each individual model in detail. The reader is referred to the Appendix for additional comparison figures of individual model simulations and MAX-DOAS data.

While the calculation of an ensemble median is a common approach to reduce individual model outliers, it is mainly used here for the sake of simplicity, allowing easier overall evaluation of how the models compare to MAX-DOAS retrievals. The model ensemble is based on five of the seven models (though with partly different setups) which constitute the CAMS regional model ensemble (http://www.regional.atmosphere.copernicus.eu/, last access: 7 February 2020) for which have shown that, at least for ozone, the ensemble median performs on average best in terms of statistical indicators compared to the seven individual models and that the ensemble is also robust against reducing the ensemble size by one member. Statistical indicators for ${\mathrm{NO}}_{\mathrm{2}}$ (see Tables  to ) show that the ensemble median of the present study performs best in terms of overall correlation to individual MAX-DOAS measurements at each station. Compared to individual models for other statistical indicators and also comparisons for seasonal, diurnal and weekly cycles, reasonable results are achieved by the ensemble median.

Table 3

Statistics on how AVK-weighted tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$) from regional models compare to MAX-DOAS retrievals at the four MAX-DOAS stations. Each column entry shows (from left to right) the root mean squared error (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$), bias (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$) and Pearson correlation coefficient. MOCAGE data are not available for the measurement time period at OHP.

	De Bilt			Bremen			Uccle			OHP
	RMSE	Bias	$r$	RMSE	Bias	$r$	RMSE	Bias	$r$	RMSE	Bias	$r$
ENSEMBLE	7.955	0.833	0.465	4.142	0.898	0.583	6.862	$-\mathrm{2.220}$	0.719	1.369	$-\mathrm{0.341}$	0.683
LOTOS-EUROS	10.516	5.254	0.352	5.180	1.598	0.461	7.815	$-\mathrm{0.897}$	0.598	3.004	0.217	0.666
CHIMERE	8.221	0.573	0.402	3.960	0.198	0.533	7.059	$-\mathrm{1.950}$	0.686	1.269	$-\mathrm{0.563}$	0.627
EMEP	8.680	0.919	0.393	4.558	$-\mathrm{0.167}$	0.521	8.134	$-\mathrm{3.609}$	0.624	1.824	$-\mathrm{0.338}$	0.600
EMEP-MACCEVA	8.340	$-\mathrm{0.182}$	0.397	5.308	2.427	0.554	7.591	$-\mathrm{2.777}$	0.659	2.012	$-\mathrm{0.473}$	0.532
SILAM	8.516	0.115	0.408	4.506	0.570	0.550	7.985	$-\mathrm{3.385}$	0.633	2.577	1.195	0.482
MOCAGE	9.731	3.082	0.427	5.651	1.801	0.520	7.413	1.476	0.692

As the typical error on MAX-DOAS retrieved VCDs is around 20 % but can be higher for sites with low trace gas concentrations like OHP or due to instrumental conditions (see Sect. ), a conservative overall uncertainty of MAX-DOAS retrievals of 30 % is assumed for all stations within this paper and given along with the data plots, where appropriate. Data products with more detailed uncertainty information are currently in development, for example, in the framework of the FRM4DOAS project (http://frm4doas.aeronomie.be/, last access: 7 February 2020), and once available, these data and related uncertainty information can be used for other comparison studies.

The horizontal grid spacing (Table ) differs for the six model runs evaluated in the present study, with a resolution of approximately 9 $\times$ 7 km${}^{\mathrm{2}}$ for the finest resolution run (LOTOS-EUROS) and 50 $\times$ 50 km${}^{\mathrm{2}}$ for the coarsest one (EMEP). The resolution of the remaining model runs is approximately 20 $\times$ 20 km${}^{\mathrm{2}}$. As described in Sect. , the horizontal averaging volume of MAX-DOAS retrievals strongly depends on aerosol loading, viewing direction and wavelength . As a rough estimate, it ranges from 5 to 10 km for the stations used in the present study. Therefore, the horizontal averaging volume is (apart from the coarsest resolution run) expected to be either on the same spatial scale as the horizontal model resolution or by a factor of 1–4 smaller. From the latter (i.e., horizontal averaging volume of MAX-DOAS finer than model resolution), one would expect an underestimation of enhancements in tropospheric columns observed by MAX-DOAS in the event of horizontal changes in tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ columns below the model resolution and, similarly, an overestimation of local minima in tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ columns. However, in reality, the comparison between horizontal averaging volume of MAX-DOAS and horizontal resolution of the models is much more complicated, as MAX-DOAS instruments usually measure in one azimuthal pointing direction meaning that measurements, though integrated along the line of sight and over several elevation angles, do not necessarily fully represent the ${\mathrm{NO}}_{\mathrm{2}}$ amount in the three-dimensional grid boxes of the model simulations. This could, for example, mean that a pollution plume with a horizontal extent on the order of the model resolution and hence showing up in the simulations is not observed by the MAX-DOAS instrument.

4 Results

Figure  shows time series of AVK-weighted tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs from MAX-DOAS and model ensemble data. The magnitude of VCDs from the measurements for Bremen and OHP is reproduced by the model ensemble. For De Bilt, the spread of retrieved ${\mathrm{NO}}_{\mathrm{2}}$ VCDs is larger than for simulated ones. At all of the four stations, measurements and simulations show large deviations for some of the time steps. Larger ${\mathrm{NO}}_{\mathrm{2}}$ values inside individual pollution plumes are generally underestimated by the model ensemble, especially at Uccle and De Bilt. This is in agreement with , who compared car MAX-DOAS measurements and OMI retrievals with a regional model (CHIMERE) and found that values inside emission plumes are systematically underestimated. However, especially in Bremen and OHP, the model ensemble also shows peaks at time steps when there was no or much less enhancement in retrieved ${\mathrm{NO}}_{\mathrm{2}}$ VCDs. Out of the four MAX-DOAS sites investigated, the magnitude of ${\mathrm{NO}}_{\mathrm{2}}$ VCDs is lowest at the rural station (OHP), which is sometimes affected by nearby pollution plumes that show up in the time series. Further investigation shows that most of these peaks are associated with northeasterly wind directions and hence pollution sources to the northeast of the station such as the Po Valley (Italy). At OHP, retrieved tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ columns are generally a bit higher than simulated ones, resulting in a negative bias of $-\mathrm{0.34}\times {\mathrm{10}}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$ for the ensemble (Table ).

Figure 2

Time series of AVK-weighted tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$) from (black circles) MAX-DOAS and (colored circles) model ensemble hourly data for (a–d) De Bilt, Bremen, Uccle and OHP.

[Figure omitted. See PDF]

The evolution of time series of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs is largely determined by the evolution of surface partial columns (see Fig. ) which already account for about 25 % of the magnitude of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs. In the present study, surface partial columns refer to the partial column of the lowest measurement layer (Bremen 50 m, De Bilt 180 m, Uccle 180 m, OHP 150 m a.g.l.). As vertical profiles are not available from the MAX-DOAS output for De Bilt, comparisons of surface partial columns are not given for this station in the present paper. The same conclusions as for tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs described in the previous paragraph arise for surface partial columns when comparing the model ensemble to MAX-DOAS data. However, the negative bias found for OHP for tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs is not present when looking at the surface partial column time series for this station (see Tables  and ). This result is consistent with average vertical profiles for OHP shown in Fig. , with the regional models slightly overestimating partial columns at the surface, while the majority of the models underestimate partial columns at altitudes between 1 and 1.5 km (Fig. ).

Figure 3

As in Fig.  but for ${\mathrm{NO}}_{\mathrm{2}}$ surface partial columns (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$). Surface partial columns from MAX-DOAS are not available for De Bilt for the investigated time period.

[Figure omitted. See PDF]

Table 4

As in Table  but for ${\mathrm{NO}}_{\mathrm{2}}$ surface partial columns (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$). Surface partial columns from MAX-DOAS are not available for De Bilt for the investigated time period.

	Bremen			Uccle			OHP
	RMSE	Bias	$r$	RMSE	Bias	$r$	RMSE	Bias	$r$
ENSEMBLE	0.715	0.123	0.374	2.905	$-\mathrm{1.124}$	0.586	0.351	0.058	0.439
LOTOS-EUROS	0.783	0.181	0.336	3.309	$-\mathrm{1.659}$	0.509	0.517	0.048	0.393
CHIMERE	0.927	0.364	0.252	2.902	$-\mathrm{0.554}$	0.531	0.337	0.081	0.400
EMEP	0.723	$-\mathrm{0.133}$	0.318	3.330	$-\mathrm{1.819}$	0.533	0.443	0.068	0.417
EMEP-MACCEVA	0.869	0.414	0.320	3.229	$-\mathrm{1.658}$	0.548	0.428	0.014	0.344
SILAM	0.681	$-\mathrm{0.054}$	0.397	2.910	$-\mathrm{0.498}$	0.572	0.659	0.388	0.318
MOCAGE	0.750	0.101	0.372	2.886	0.272	0.596

Although there are larger differences between simulations and retrievals especially for individual pollution plumes, Fig.  shows that frequency distributions of AVK-weighted tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs are similar for ensemble simulations and observations. However, for OHP, the number of data values with tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs lower than $\mathrm{1}\times {\mathrm{10}}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$ is significantly larger for model simulated values (about 1400 model values compared to about 200 observed data counts), in agreement with the negative bias in tropospheric columns described above.

Figure 4

Frequency distributions of AVK-weighted tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$) from (blue) MAX-DOAS and (red) model ensemble data for (a) De Bilt, (b) Bremen, (c) Uccle and (d) OHP. The distance between vertical gray lines on the $x$ axis corresponds to the size of the bins used to calculate the number of values given on the $y$ axis.

[Figure omitted. See PDF]

Figure  shows model simulated and MAX-DOAS retrieved vertical profiles of ${\mathrm{NO}}_{\mathrm{2}}$ partial columns averaged over the whole time period of measurements together with a priori profiles and AVKs for completeness. Averages of vertical profiles over different seasons are given in Fig. , in order to investigate consistency between profiles throughout different times of the year. In general, differences between retrievals and simulations are largest for larger ${\mathrm{NO}}_{\mathrm{2}}$ partial columns, which means for the lower altitude layers and during the colder winter and autumn seasons. Many of the values simulated by individual models do not fall into the uncertainty range of MAX-DOAS retrievals assumed here. For example, SILAM largely overestimates ${\mathrm{NO}}_{\mathrm{2}}$ partial columns up to 1.5 km altitude at OHP, while MOCAGE (apart from the lowest observation layer) strongly overestimates values up to about 1 km altitude at Uccle. Although model ensemble profiles show some differences from the retrievals regarding the shape and magnitude of the profiles, they also show the largest partial columns close to the surface for all of the three stations investigated. This result also shows up throughout different seasons.

Figure 5

Average vertical profiles of ${\mathrm{NO}}_{\mathrm{2}}$ partial columns (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$) from (black) MAX-DOAS, (brown) a priori profiles used for MAX-DOAS retrievals, (gray) model ensemble median, (blue) LOTOS-EUROS, (yellow) CHIMERE, (green) EMEP, (orange) EMEP-MACCEVA, (pink) SILAM and (red) MOCAGE as well as (purple) column AVKs (unitless) for (a) Bremen, (b) Uccle and (c) OHP. Black error bars refer to the uncertainty associated with the MAX-DOAS retrievals (assumed to be 30 % for all stations). MAX-DOAS vertical profiles are not available for De Bilt for the investigated time period.

[Figure omitted. See PDF]

As the sensitivity of MAX-DOAS retrievals is largest in the boundary layer, a feature which is independent of the retrieval method, we initially expected the application of column AVKs from the measurements to model simulations to be of crucial importance for evaluation results. However, further analysis showed that applying column AVKs to model ${\mathrm{NO}}_{\mathrm{2}}$ partial columns before summing these up in the vertical does not have a big impact on derived tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs and therefore has a minor effect on the data analysis presented in this paper. Only AVK-weighted simulations of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs are therefore shown here. Statistical values (root mean squared error, bias, Pearson correlation coefficient) which will be described below are quite similar for AVK-weighted model ensemble VCDs and those from non-AVK-weighted ones. One of the reasons for this is that (as shown in Figs.  and ), AVKs are close to 1 around the boundary layer where MAX-DOAS instruments have the highest sensitivity (generally a bit larger than 1 close to the surface and smaller than 1 higher up, which has a balancing effect) and that the vertical shape of the column AVK curve is in principal agreement with the shape of simulated ${\mathrm{NO}}_{\mathrm{2}}$ partial columns. At altitudes above 1 km, AVKs are smaller than 1, but simulated ${\mathrm{NO}}_{\mathrm{2}}$ partial columns are also significantly smaller at these compared to lower altitudes, so that the contribution to the tropospheric column is limited. At higher altitudes, MAX-DOAS retrievals tend to follow the a priori profile, while retrievals in the boundary layer are not much influenced by the a priori profile in general. This is in contrast to the situation for satellite observations of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$, which usually have a minimum of the AVK in the boundary layer, i.e., where the largest fraction of ${\mathrm{NO}}_{\mathrm{2}}$ is usually located in polluted situations. A priori profiles used within the MAX-DOAS retrievals (see Sect. ) are in principal agreement with the ones simulated by the models. The vertical weighting caused by application of AVKs to partial columns therefore does not significantly impact derived tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs.

Scatter density plots of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs from MAX-DOAS against model values corresponding to the time series displayed in Fig.  are shown in Fig.  (see Fig.  for individual model results). Statistical values (root mean squared error, bias, Pearson correlation) and least squares regression lines are given along with the plots to draw further conclusions on the ability of the model ensemble to reproduce MAX-DOAS retrievals and are listed in Table  (together with statistics on ensemble members). Statistical values for surface partial columns are given in Table . The ensemble median performs best in terms of overall correlation with values between 45 % and 75 % (compared to 35 % to 70 % for individual models) for tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs for all stations, the highest correlation is found for Uccle. Note, however, that for other statistical indicators, some of the individual models perform better. Correlations are generally lower than the ones based on tropospheric columns for surface partial columns which are on the order of 40 % for Bremen and OHP but much higher again for Uccle ($\sim$ 60 %) for the ensemble. The model ensemble has a negative bias on the order of $-\mathrm{0.3}$ and $-\mathrm{2}\times {\mathrm{10}}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$ at OHP and Uccle, respectively, and a positive bias on the order of $\mathrm{1}\times {\mathrm{10}}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$ at De Bilt and Bremen for tropospheric columns. The largest RMSE and bias (10.5 and $\mathrm{5}\times {\mathrm{10}}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$, respectively) are found for LOTOS-EUROS at De Bilt. Considering that tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs for OHP are generally lower than for the three urban sites, SILAM also shows a high RMSE and bias (2.6 and $\mathrm{1.2}\times {\mathrm{10}}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$, respectively) at this station. Vertical profile comparisons described above show that the overestimation is largest at altitudes up to about 1.5 km. Our findings agree with , who found that SILAM tends to overestimate ${\mathrm{NO}}_{\mathrm{2}}$ at rural sites based on in situ data and concluded that this is due to an overestimation of the lifetime of ${\mathrm{NO}}_{\mathrm{2}}$, which is also consistent with findings by . For surface partial columns, biases are negligibly small for OHP and Bremen for the ensemble and most of the individual models, while the ensemble is negatively biased by $\mathrm{1.1}\times {\mathrm{10}}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$ at Uccle. The largest RMSE and bias in surface partial columns are found for EMEP at Uccle (3.3 and $-\mathrm{1.8}\times {\mathrm{10}}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$, respectively). The spread between models and observations is large for some individual data points. Regression lines show that the model ensemble tends to overestimate low and underestimate high tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs. The underestimation of larger tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs is most pronounced for De Bilt, followed by Uccle.

Figure 6

Scatter density plots of AVK-weighted tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$) from MAX-DOAS against model ensemble data for (a) De Bilt, (b) Bremen, (c) Uccle and (d) OHP. The data are shown for different bins with a size of 10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$ and are colored according to the number of data points per bin (%). The dashed line is the reference line ($f\left(x\right)=x$). The solid line is the regression line (see top left of each plot for $f\left(x\right)$ of this line). The root mean squared error (RMSE) (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$), bias (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$), Pearson correlation coefficient (cor, not squared) as well as the number of data points $N$ are given in the top left of each plot.

[Figure omitted. See PDF]

Figure  shows comparisons between MAX-DOAS and the model ensemble of wind directional distributions of average tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs based on wind measurements from station data (note that further analysis has shown a good agreement between measured wind speeds and wind directions and those of the simulations). Changes of ${\mathrm{NO}}_{\mathrm{2}}$ mean values from one wind direction bin to another are reproduced well by the model ensemble (and in general also by ensemble members; see Fig. ), with an overall slightly better agreement with retrievals for tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs compared to surface partial columns (not shown). Both MAX-DOAS and model ensemble show the highest ${\mathrm{NO}}_{\mathrm{2}}$ mean values for wind directions mainly where an influence from pollution sources is expected (i.e., the Ruhr area to the southeast of De Bilt, the Bremen city center to the southwest of the Bremen MAX-DOAS, Brussels city center to the northeast of Uccle, the Po Valley to the northeast of OHP; see Sect. ).

Figure 7

Average AVK-weighted tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$) in 45${}^{\circ }$ wide wind direction bins from (solid blue lines) MAX-DOAS and (dashed red lines) model ensemble data for (a) De Bilt, (b) Bremen, (c) Uccle and (d) OHP. Wind directions correspond to the direction towards the station and are taken from weather station measurements. The numbers close to the center of each plot refer to the number of data values used for calculating average values for each bin.

[Figure omitted. See PDF]

Comparisons for seasonal cycles (i.e., monthly averages) of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs are given in Fig.  together with corresponding statistical values in Table . There is good agreement between MAX-DOAS and the model ensemble for Uccle regarding the magnitude of ${\mathrm{NO}}_{\mathrm{2}}$ VCDs and seasonality, with simulated ensemble median values within the estimated uncertainty interval of the retrievals. The same is true for De Bilt, apart from the strong overestimation of MAX-DOAS retrieved values for January, March and April. However, a low number of observations is available during these compared to other months, so the results are statistically less significant. The model ensemble overestimates seasonal cycles for Bremen and OHP. More explicitly, there is an overestimation of wintertime values, while summertime values are better reproduced by the model ensemble. compared an ensemble of regional and global models to satellite data over Europe and found an overestimation of seasonal cycles by the simulations, which is in agreement with results for Bremen and OHP shown in the present paper. In , however, model values were closer to satellite retrievals during winter, whereas for summer a strong underestimation was found, while comparisons to Dutch surface observations showed that this could be partly attributed to a high bias of satellite retrievals in summer at least over the Netherlands. Moreover, showed that the GEMS regional model ensemble median underestimates background values of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ columns compared to OMI satellite retrievals in summer and spring. This is in agreement with results shown in Fig. , where the majority of regional models underestimate summer and springtime values at OHP. The spread between individual models is quite large for seasonal cycles at OHP, meaning that some of the models perform better than others. Seasonal cycles are generally more pronounced compared to the other model runs and retrievals for LOTOS-EUROS and MOCAGE. Especially LOTOS-EUROS largely overestimates the observed seasonal cycle at OHP. Low to moderate correlations in seasonal cycles are found for De Bilt, followed by moderate ones for Bremen. All models perform well in terms of correlation at Uccle and OHP (values around 0.8).

Table 5

As in Table  but for (upper rows) seasonal cycles, (middle rows) diurnal cycles and (lower rows) weekly cycles of AVK-weighted tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$). In addition, values for diurnal cycles are given based on data during weekdays only and during weekends only for the ensemble.

		De Bilt			Bremen			Uccle			OHP
		RMSE	Bias	$r$	RMSE	Bias	$r$	RMSE	Bias	$r$	RMSE	Bias	$r$
Seasonal	ENSEMBLE	3.750	2.494	0.511	2.218	1.347	0.710	2.857	$-\mathrm{2.074}$	0.863	1.055	0.018	0.755
	LOTOS-EUROS	6.982	5.920	0.413	2.808	1.807	0.708	2.648	$-\mathrm{0.892}$	0.777	1.975	0.752	0.859
	CHIMERE	4.001	2.592	0.343	1.291	0.519	0.627	3.337	$-\mathrm{1.853}$	0.810	0.563	$-\mathrm{0.507}$	0.878
	EMEP	3.237	1.837	0.575	1.960	0.337	0.743	4.000	$-\mathrm{3.071}$	0.820	1.187	$-\mathrm{0.020}$	0.733
	EMEP-MACCEVA	3.647	1.587	0.304	3.344	2.595	0.702	3.455	$-\mathrm{2.581}$	0.832	0.903	$-\mathrm{0.208}$	0.795
	SILAM	2.904	1.631	0.617	1.933	0.905	0.688	3.323	$-\mathrm{2.759}$	0.853	1.506	1.270	0.798
	MOCAGE	7.443	5.315	0.214	3.823	2.414	0.643	2.621	1.594	0.853
Diurnal	ENSEMBLE	2.571	0.749	$-\mathrm{0.389}$	1.499	1.134	0.148	2.401	$-\mathrm{1.735}$	0.815	0.599	$-\mathrm{0.371}$	$-\mathrm{0.069}$
	ENSEMBLE (weekdays)	2.762	0.103	$-\mathrm{0.237}$	1.343	0.752	0.125	3.181	$-\mathrm{2.527}$	0.827	0.639	$-\mathrm{0.390}$	$-\mathrm{0.251}$
	ENSEMBLE (weekends)	3.438	3.006	0.559	2.664	2.463	0.405	1.294	0.341	0.806	0.461	$-\mathrm{0.253}$	0.567
	LOTOS-EUROS	6.283	4.959	$-\mathrm{0.525}$	1.777	1.362	0.413	2.055	$-\mathrm{0.448}$	0.715	0.756	0.106	$-\mathrm{0.255}$
	CHIMERE	2.560	0.813	$-\mathrm{0.456}$	0.614	0.250	0.335	1.927	$-\mathrm{1.379}$	0.919	0.703	$-\mathrm{0.641}$	0.585
	EMEP	2.690	0.482	$-\mathrm{0.219}$	1.339	$-\mathrm{0.251}$	0.052	3.568	$-\mathrm{2.739}$	0.567	0.711	$-\mathrm{0.435}$	$-\mathrm{0.298}$
	EMEP-MACCEVA	2.681	$-\mathrm{0.450}$	$-\mathrm{0.424}$	2.971	2.574	$-\mathrm{0.235}$	2.911	$-\mathrm{2.194}$	0.724	0.662	$-\mathrm{0.540}$	0.128
	SILAM	2.620	$-\mathrm{0.199}$	$-\mathrm{0.327}$	1.429	0.497	$-\mathrm{0.085}$	3.410	$-\mathrm{2.895}$	0.763	1.468	1.301	0.270
	MOCAGE	3.782	2.639	$-\mathrm{0.262}$	2.070	1.670	0.149	2.277	1.433	0.812
Weekly	ENSEMBLE	0.129	0.009	0.917	0.131	$-\mathrm{0.009}$	0.820	0.101	0.006	0.967	0.060	$-\mathrm{0.009}$	0.802
	LOTOS-EUROS	0.122	0.014	0.970	0.136	$-\mathrm{0.008}$	0.699	0.122	$-\mathrm{0.005}$	0.973	0.129	0.005	0.572
	CHIMERE	0.147	0.029	0.973	0.138	$-\mathrm{0.009}$	0.878	0.134	$-\mathrm{0.006}$	0.954	0.036	0.001	0.935
	EMEP	0.186	0.014	0.708	0.176	$-\mathrm{0.013}$	0.278	0.095	$-\mathrm{0.004}$	0.928	0.082	0.004	0.696
	EMEP-MACCEVA	0.146	0.014	0.865	0.141	$-\mathrm{0.010}$	0.727	0.126	$-\mathrm{0.006}$	0.828	0.130	0.005	0.335
	SILAM	0.097	0.007	0.930	0.083	$-\mathrm{0.006}$	0.898	0.028	0.000	0.993	0.102	$-\mathrm{0.001}$	0.557
	MOCAGE	0.139	$-\mathrm{0.009}$	0.850	0.117	$-\mathrm{0.008}$	0.899	0.120	$-\mathrm{0.006}$	0.872

Figure 8

Seasonal cycles (monthly averages) of AVK-weighted tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$) from (black) MAX-DOAS, (gray) model ensemble median, (blue) LOTOS-EUROS, (yellow) CHIMERE, (green) EMEP, (orange) EMEP-MACCEVA, (pink) SILAM and (red) MOCAGE for (a) De Bilt, (b) Bremen, (c) Uccle and (d) OHP. Black error bars refer to the uncertainty associated with the MAX-DOAS retrievals (assumed to be 30 % for all stations). The number of data values used for calculating average values is shown in the upper $x$ axis of each plot as an indicator of statistical significance.

[Figure omitted. See PDF]

Figure  shows comparisons of diurnal cycles for the whole time series. Overall, the model ensemble fails to reproduce diurnal cycles, reflected by generally low correlations (Table ) for all models in De Bilt, Bremen and OHP. All models show negative correlations in De Bilt, while some of the models only reach negative correlations in Bremen as well. MAX-DOAS retrieved values increase from the morning towards the afternoon, while simulated values in general decrease from the morning towards the afternoon. In Uccle, however, high or at least moderate correlations are achieved. CHIMERE performs best in terms of correlation in Uccle and OHP (0.92 and 0.6, respectively). For this model, diurnal scaling factors of traffic emissions have been developed by analyzing measurements of ${\mathrm{NO}}_{\mathrm{2}}$ in European countries (). Although most of the model values fall within the estimated uncertainty interval of MAX-DOAS retrievals, the shape of diurnal cycles differs between observations and simulations. The ensemble shows a strong peak during the morning rush hour around 08:00 LT for Bremen, which is not confirmed by MAX-DOAS retrievals. In contrast to this, measurements show a maximum around 14:00 LT in the afternoon which coincides with a very weak local maximum simulated by the model ensemble. Looking at diurnal cycles for different seasons shown in Figs.  and reveals that these are in general much better reproduced for spring and summer compared to autumn and winter for all stations. This is in agreement with results for seasonal cycles described in the previous paragraph. Weak morning rush hour peaks are also simulated for the rural station (OHP) (possibly related to pollution transport towards the station from surrounding areas), which is not in agreement with the measurements. The morning rush hour peaks for Bremen and OHP occur for all models with the exception of SILAM for OHP, which however strongly overestimates values (by a factor of 1.5–2 for diurnal cycle values averaged over the whole time series) for this station, resulting in a bias of $\sim \mathrm{1.3}\times {\mathrm{10}}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$ (see Table ). The peak at 08:00 LT for Bremen is most pronounced for EMEP-MACCEVA, MOCAGE and LOTOS-EUROS. Individual model runs show the same shape of the diurnal cycle for Bremen, while the shape of diurnal cycles differs for OHP. Moreover, large differences regarding the magnitude of simulated values occur for both stations between the models. As described in Sect. , all models use the same emission inventory as a basis, except the EMEP run. There is a strong difference between the magnitude of the values simulated by EMEP and EMEP-MACCEVA specifically for the diurnal cycle at Bremen (while their shape agrees). EMEP and EMEP-MACCEVA runs were performed with different resolution and emission inventories which likely contributes to the differences found.

Figure 9

As in Fig.  but for diurnal cycles (averages over hourly bins) of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$).

[Figure omitted. See PDF]

The simulated morning rush hour peaks for Bremen coincide with maxima of diurnal scaling factors for the road transport and combustion pollutant sectors (see Fig. 1 of ). recently derived diurnal profiles from traffic counts for weekdays. The latter show, for several cities, less pronounced rush hours and much flatter time profiles during the day. This is in principal agreement with the MAX-DOAS observations of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs. The comparison to the MAX-DOAS retrievals in Bremen shows that the diurnal scaling factors for road transport and potentially combustion currently used by the models should be adjusted to match the observed temporal profiles.

Figure  shows comparisons of diurnal cycles for weekends (Saturdays and Sundays) only. A figure of diurnal cycles for weekdays only shows very similar results to Fig.  (and is therefore not shown here), meaning that overall diurnal cycles are mainly driven by weekday emissions. At the three urban stations, MAX-DOAS retrieved diurnal cycles show a different shape for weekends compared to diurnal cycles for the whole week (and hence weekdays only). This is in contrast to model simulated diurnal cycles, which do not change much going from cycles for the whole week to cycles for weekends only, apart from a general decrease in values towards weekends for both retrieved and simulated tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs. As expected, MAX-DOAS retrieved diurnal cycles are rather flat for weekends only at the urban stations, as emissions from traffic and industry are reduced during weekends compared to weekdays (e.g., ). As the shape of simulated diurnal cycles is similar for weekdays versus weekends, the difference between retrieved and simulated trends in tropospheric columns from morning to afternoon hours is reduced for weekends only resulting in significantly higher and positive correlations for diurnal cycles during weekends compared to weekdays for the ensemble at these stations (see Table ). In Uccle, correlations are equally high (about 80 %) for weekdays and weekends, which is due to the fact that the change in the shape of retrieved diurnal cycles from weekdays to weekends is comparably small at this urban station. Correlations are also significantly higher for the background station (OHP) for weekends for the ensemble, mainly due to a better agreement in the development from the afternoon towards the evening during weekends. The results described above show that models fail to reproduce observed changes in diurnal cycles from weekdays to weekends at urban stations. This is because the diurnal scalings incorporated in the runs stay constant throughout the week (see ). As observed diurnal cycles are rather flat at weekends, switching off diurnal scalings during weekends will likely lead to an improvement in model performance compared to MAX-DOAS.

Figure 10

As in Fig.  but for diurnal cycles (averages over hourly bins) of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$) during weekends only.

[Figure omitted. See PDF]

Weekly cycle comparisons are presented in Fig.  (see Figs.  and for different seasons). In contrast to diurnal cycles, weekly cycles and their seasonal variation measured by MAX-DOAS are much better simulated, reflected by high correlations (Table ) for the ensemble at all stations. Both MAX-DOAS and the model ensemble show a decrease in tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs towards the weekend when there is less traffic especially for the urban De Bilt, Bremen and Uccle stations. However, this observed weekly cycle is stronger than the simulated one, a feature which is most pronounced for Bremen. This can be attributed to the fact that the diurnal scalings used in the model runs do not change from weekdays towards weekends (see previous paragraph), so that the simulations are positively biased by emission peaks at the weekend. The underestimation of weekly cycles is in agreement with , who compared LOTOS-EUROS simulations to MAX-DOAS retrievals for De Bilt. As expected, only a very weak weekly cycle is observed by MAX-DOAS and simulated by the models for the rural station (OHP). Note that maxima of weekly cycles for specific days may just be a coincidence due to data sampling times. investigated weekly cycles of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ based on Global Ozone Monitoring Experiment (GOME) satellite observations and found a decrease in values of up to about 50 % towards Sundays over polluted regions and cities in Europe. This is in principal agreement with results of the present study, although the choice of the cities is different.

Figure 11

As in Fig.  but for weekly cycles (averages over daily bins divided by the mean over the whole week; unitless values) of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs.

[Figure omitted. See PDF]

Table shows that the overall correlations (Table ) reached at all stations are mainly driven by seasonal and weekly cycles, while significantly lower and in many cases negative correlations are found for diurnal cycles which decreases overall correlations. An exception for the latter is Uccle, where good correlations are also found for diurnal cycles.

A pollution plume and related increase in the time series of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs observed by MAX-DOAS is expected to be reproduced better by model runs with finer horizontal resolution compared to coarser resolution runs. Moreover, increasing the amount of ${\mathrm{NO}}_{\mathrm{2}}$ at relatively high ${\mathrm{NO}}_{\mathrm{2}}$ concentrations results in a decrease of OH and hence increasing lifetime of ${\mathrm{NO}}_{\mathrm{2}}$ . compared regional model simulations with differing horizontal resolution and found that an increase in resolution leads to a better agreement with ${\mathrm{NO}}_{\mathrm{2}}$ in situ data. In the present study, statistical values show that the LOTOS-EUROS run with significantly finer horizontal resolution than the other runs (see Table ) in general did not perform better than coarser horizontal resolution runs, but the LOTOS-EUROS run also has the lowest number of vertical levels. Similarly, the EMEP run with significantly coarser horizontal resolution did not perform worse than finer horizontal resolution runs but also has the largest number of vertical levels. Overall, no relation between the horizontal and vertical resolution of the models and their statistical performance is found.

In addition to the MAX-DOAS comparisons described above, we also carried out a comparison between the regional models and OMI satellite retrievals looking at maps of monthly means for a winter and summer month (February and August 2011, respectively) falling into the time period investigated by the present study. We found similar results to , which are therefore not shown here, i.e., an underestimation of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ columns over background regions during summer (in agreement with the general underestimation of means over summer months compared to MAX-DOAS shown by seasonal cycles for OHP for all models except SILAM) and a generally better agreement between satellite retrievals and models over pollution hotspots around Benelux countries, an underestimation however of values over large parts of Germany and over the Po Valley in many of the model runs. Some of the models also overestimated values to the south and southeast of OHP (roughly between Marseille and Genua along the southern coast of France) compared to OMI. However, due to the generally short lifetime of ${\mathrm{NO}}_{\mathrm{2}}$, to properly relate uncertainties in the simulations over emission hotspots indicated by the OMI-based comparisons to the ones derived from MAX-DOAS-based comparisons would generally require investigating transport patterns of individual model runs with much finer time resolution around the MAX-DOAS sites, which is not provided by the satellite data (only one OMI orbit per day over the stations).

5 Summary and conclusions

In this study, comparisons between ${\mathrm{NO}}_{\mathrm{2}}$ columns simulated by five regional models and retrieved from MAX-DOAS measurements for four European MAX-DOAS stations have been presented. According to , the main uncertainties of the regional model air quality forecasts are generally emissions, meteorological and chemical initial and boundary conditions used to drive the simulations as well as representation of physical and chemical processes. This also explains parts of the differences between model results and observations found in the present study. Uncertainties in tropospheric MAX-DOAS retrievals further contribute to differences described here (a conservative uncertainty of 30 % was assumed for all stations).

In general, the model ensemble does well represent tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ amounts observed by MAX-DOAS. Moderate correlations around 60 % are found for tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs at each station for the ensemble. Time series comparisons and corresponding scatterplots show that MAX-DOAS retrieved pollution peaks are often underestimated by the model ensemble, while models also show peaks at times when the instruments did not observe similar enhancements. Comparisons of wind directional distributions of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs and surface partial columns show good agreement between simulations and measurements. This indicates that transport of pollution towards the stations is, on average, well represented by the models.

Comparison of vertical profiles shows that the main differences between measurements and simulations are also not due to incorrect representation of the vertical ${\mathrm{NO}}_{\mathrm{2}}$ distribution. Hence, it is not necessarily needed to apply AVKs to model simulations in MAX-DOAS-based comparison studies of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs, a simple summation of simulated partial columns in the vertical can already be sufficient. This result was not expected, as the sensitivity of the MAX-DOAS profile retrievals is much larger close to the surface than at altitudes larger than approximately 1 km.

Moderate to high correlations (50 %–97 %) are achieved by the ensemble for seasonal and weekly cycles, while correlations $<\mathrm{15}$ % are found for diurnal cycles at all stations except Uccle. This shows that simulation uncertainties in seasonal cycles (though overestimated by the model ensemble) and in reproducing the decrease in daily tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs towards the weekend (though underestimated by the models) do not explain the main part of the scatter between observations and simulations. The largest differences from the retrieved diurnal cycles generally occurred for LOTOS-EUROS and MOCAGE at Bremen and De Bilt, and for EMEP-MACCEVA at Bremen, LOTOS-EUROS and SILAM showed the largest differences to retrieved diurnal cycles for the background OHP station. As the diurnal scalings incorporated in the model runs stay constant throughout the week, models fail to reproduce observed changes in diurnal cycles from weekdays to weekends at urban stations, which is also consistent with the identified underestimation of weekly cycles. Switching off diurnal scalings during weekends is expected to already lead to an improvement in model performance. Moreover, the comparison to the MAX-DOAS retrievals at Bremen shows that the diurnal scaling factors for road transport and potentially combustion should be adjusted in order to match the observed temporal profiles. The results of this study clearly demonstrate that future model development needs to concentrate on improving representation of diurnal cycles and associated temporal scalings.

As described in Sect. , a pollution plume with a horizontal extent on the order of the model resolution and hence showing up in the simulations can be missed by the line of sight of the MAX-DOAS instrument. Multiple MAX-DOAS measurements should therefore be performed over a range of different azimuthal angles for each station and stricter quality filters could be applied on the MAX-DOAS data to reduce the impact of retrieval uncertainty in future MAX-DOAS comparison studies. Comparison results could also be complemented by further data sources where possible. It would be interesting to investigate the ability of the models to predict the scales of ${\mathrm{NO}}_{\mathrm{2}}$ spatial variations derived from timescales of ${\mathrm{NO}}_{\mathrm{2}}$ variations and wind speeds in the context of model resolution. Moreover, one could investigate the ability of the models to distribute ${\mathrm{NO}}_{\mathrm{2}}$ in the vertical in terms of characteristic layer height of ${\mathrm{NO}}_{\mathrm{2}}$, which is expected to be affected by vertical resolution of the models.

Similar comparisons to other regional models or other model setups (different chemistry schemes, resolutions, anthropogenic emission inventories, boundary layer physics, varying scalings of emissions in time and in the vertical) as well as for more MAX-DOAS sites should follow. As the stations investigated in the present study have, apart from the rural background OHP station, rather similar meteorological and pollution conditions, investigation of stations over a broader range of different conditions would be desirable. Further comparison studies could, for instance, include stations at pollution hotspots in the Mediterranean such as Athens with strong smog conditions especially during summer and clean mountain sites.

Appendix A Figure A1

Average vertical profiles of ${\mathrm{NO}}_{\mathrm{2}}$ partial columns (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$) from (black) MAX-DOAS, (brown) a priori profiles used for MAX-DOAS retrievals, (gray) model ensemble median, (blue) LOTOS-EUROS, (yellow) CHIMERE, (green) EMEP, (orange) EMEP-MACCEVA, (pink) SILAM and (red) MOCAGE as well as (purple) column AVKs (unitless). Black error bars refer to the uncertainty associated with the MAX-DOAS retrievals (assumed to be 30 % for all stations). Panel (a) shows profiles for data averaged over the whole time series, panel (b) shows profiles for DJF (December, January February), (c) MAM (March, April, May), (d) JJA (June, July, August) and (e) SON (September, October, November) months only. Panels (a) to (e) refer to (left) Bremen, (middle) Uccle and (right) OHP. MAX-DOAS vertical profiles are not available for De Bilt for the investigated time period.

[Figure omitted. See PDF]

Figure A2

Scatter density plots of AVK-weighted tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$) from MAX-DOAS against model data for (from left to right) De Bilt, Bremen, Uccle and OHP. The data are shown for different bins with a size of 10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$ and are colored according to the number of data points per bin (%). The different panels show different model runs: (a) LOTOS-EUROS, (b) CHIMERE, (c) EMEP, (d) EMEP-MACCEVA, (e) SILAM and (f) MOCAGE. MOCAGE data are not available for the measurement time period at OHP. The dashed line is the reference line ($f\left(x\right)=x$). The solid line is the regression line (see top left of each plot for $f\left(x\right)$ of this line). The root mean squared error (RMSE) (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$), bias (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$), Pearson correlation coefficient (cor, not squared) as well as the number of data points $N$ are given in the top left of each plot.

[Figure omitted. See PDF]

Figure A3

As in Fig.  but for average AVK-weighted tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$) in 45${}^{\circ }$ wide wind direction bins from (solid blue lines) MAX-DOAS and (dashed red lines) model data calculated. MOCAGE data are not available for the measurement time period at OHP.

[Figure omitted. See PDF]

Figure A4

Diurnal cycles (averages over hourly bins) of AVK-weighted tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs (10${}^{\mathrm{15}}$ molec cm${}^{-\mathrm{2}}$) from (black) MAX-DOAS, (gray) model ensemble median, (blue) LOTOS-EUROS, (yellow) CHIMERE, (green) EMEP, (orange) EMEP-MACCEVA, (pink) SILAM and (red) MOCAGE for (left) De Bilt and (right) Bremen. Model-based diurnal cycles were calculated from tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs. Black error bars refer to the uncertainty associated with the MAX-DOAS retrievals (assumed to be 30 % for all stations). Panel (a) shows cycles for the whole time series, (b) shows cycles for DJF, (c) MAM, (d) JJA and (e) SON months only. The number of data values used for calculating average values is shown in the top $x$ axis of each plot.

[Figure omitted. See PDF]

Figure A5

As in Fig.  but for (left) Uccle and (right) OHP.

[Figure omitted. See PDF]

Figure A6

As in Fig.  but for weekly cycles (averages over daily bins divided by the mean over the whole week; unitless values) of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs.

[Figure omitted. See PDF]

Figure A7

As in Fig.  but for weekly cycles (averages over daily bins divided by the mean over the whole week; unitless values) of tropospheric ${\mathrm{NO}}_{\mathrm{2}}$ VCDs.

[Figure omitted. See PDF]

Code availability

Source code and test data sets for the open-source EMEP/MSC-W model are available at https://github.com/metno/emep-ctm (last access: 7 February 2020) or by contacting EMEP/MSC-W ([email protected]). The SILAM code is available upon request from the authors ([email protected], [email protected]). The MOCAGE results in the present paper are based on source code which is presently incorporated in the MOCAGE model. The MOCAGE source code is the property of Météo-France and CERFACS, and it is based on libraries that belong to some other holders. The MOCAGE model is not open source, and routines from MOCAGE cannot be freely distributed. CHIMERE is an open-source code protected under the GNU General Public license. It can be found at http://www.lmd.polytechnique.fr/chimere/ (last access: 7 February 2020). LOTOS-EUROS is downloadable free of charge after signing a license agreement. All information concerning the LOTOS-EUROS code is available at http://lotos-euros.nl (last access: 7 February 2020). For further information, the reader can contact Astrid Manders ([email protected]).

Author contributions

AMB, AR and HE initialized and planned this study. AMB produced the comparison results, figures and tables, and wrote the main parts of the manuscript. JA, AC, LC, GF, VM, FM, JP, MP, AS, MS, AMV, JV and TV provided the model data and input to the text description about the models and on comparison results. CG, FH, EP, GP, AJMP, AR, MVR and TV provided the MAX-DOAS data and input to the text description about the measurements and on comparison results. AR provided the OMI retrievals for this study. All authors were involved in the discussion of the results.

Competing interests

The authors declare that they have no conflict of interest.

Acknowledgements

We thank the German Weather Service/Deutscher Wetterdienst for providing wind in situ data for Bremen through their website at http://www.dwd.de. We are also grateful to people behind the wind in situ data at Uccle (BIRA-IASB, Francois Hendrick, personal communication, 2016) and OHP (http://pc-meteo.obs-hp.fr/intervalle.php, last access: 7 February 2020) for providing these measurements.

Financial support

This study was funded by the European Commission under the EU Seventh Research Framework Programme (grant agreement no. 283576, MACC II), the EU Horizon 2020 Research and Innovation Programme (grant agreement no. 633080, MACC-III) and the Copernicus Atmosphere Monitoring Service (CAMS), implemented by the European Centre for Medium-Range Weather Forecasts (ECMWF) on behalf of the European Commission. It was also funded in part by the University of Bremen. The LOTOS-EUROS work was carried out within the ESA project, GLOB-EMISSION (grant no. AO/1-6721/11/I-NB). BIRA-IASB MAX-DOAS observations at Uccle and OHP were financially supported by the projects AGACC-II (BELSPO, Brussels) and NORS (EU FP7; contract 284421). The article processing charges for this open-access publication were covered by the University of Bremen.

Review statement

This paper was edited by Tim Butler and reviewed by three anonymous referees.