A top-down assessment using OMI NO₂ suggests an underestimate in the NO_x emissions inventory in Seoul, South Korea, during KORUS-AQ

In this work, we investigate the NO_x emissions inventory in Seoul, South Korea, using a regional ozone monitoring instrument (OMI) NO₂ product derived from the standard NASA product. We first develop a regional OMI NO₂ product by recalculating the air mass factors using a high-resolution (4 km × 4 km) WRF-Chem model simulation, which better captures the NO₂ profile shapes in urban regions. We then apply a model-derived spatial averaging kernel to further downscale the retrieval and account for the subpixel variability. These two modifications yield OMI NO₂ values in the regional product that are 1.37 times larger in the Seoul metropolitan region and >2 times larger near substantial point sources. These two modifications also yield an OMI NO₂ product that is in better agreement with the Pandora NO₂ spectrometer measurements acquired during the South Korea–United States Air Quality (KORUS-AQ) field campaign. NO_x emissions are then derived for the Seoul metropolitan area during the KORUS-AQ field campaign using a top-down approach with the standard and regional NASA OMI NO₂ products. We first apply the top-down approach to a model simulation to ensure that the method is appropriate: the WRF-Chem simulation utilizing the bottom-up emissions inventory yields a NO_x emissions rate of 227±94 kt yr⁻¹, while the bottom-up inventory itself within a 40 km radius of Seoul yields a NO_x emissions rate of 198 kt yr⁻¹. Using the top-down approach on the regional OMI NO₂ product, we derive the NO_x emissions rate from Seoul to be 484±201 kt yr⁻¹, and a 353±146 kt yr⁻¹ NO_x emissions rate using the standard NASA OMI NO₂ product. This suggests an underestimate of 53 % and 36 % in the bottom-up inventory using the regional and standard NASA OMI NO₂ products respectively. To supplement this finding, we compare the NO₂ and NO_y simulated by WRF-Chem to observations of the same quantity acquired by aircraft and find a model underestimate. When NO_x emissions in the WRF-Chem model are increased by a factor of 2.13 in the Seoul metropolitan area, there is better agreement with KORUS-AQ aircraft observations and the recalculated OMI NO₂ tropospheric columns. Finally, we show that by using a WRF-Chem simulation with an updated emissions inventory to recalculate the air mass factor (AMF), there are small differences (∼8 %) in OMI NO₂ compared to using the original WRF-Chem simulation to derive the AMF. This suggests that changes in model resolution have a larger effect on the AMF calculation than modifications to the South Korean emissions inventory. Although the current work is focused on South Korea using OMI, the methodology developed in this work can be applied to other world regions using TROPOMI and future satellite datasets (e.g., GEMS and TEMPO) to produce high-quality region-specific top-down NO_x emissions estimates.

1 Introduction

Nitrogen oxides (NO_x≡NO+NO₂) are a group of reactive trace gases that are toxic to human health and can be converted in the atmosphere into other noxious chemical species. In the presence of abundant volatile organic compounds (VOCs) and strong sunlight, NO_x can participate in a series of chemical reactions to generate a net accumulation of O₃, another toxic air pollutant with a longer atmospheric lifetime. NO_x also participates in a series of reactions to create HNO₃, a key contributor to acid rain, and particulate nitrate (${\mathrm{NO}}_{\mathrm{3}}^{-}$), a component of fine particulate matter (PM_2.5), an additional health hazard. There are some natural emissions of NO_x (e.g., lightning), but the majority of the NO_x emissions are from anthropogenic sources (van Vuuren et al., 2011).

There is a rich legacy of NO₂ measurements by remote sensing instruments (Burrows et al., 1999). One of these instruments is the Dutch–Finnish ozone monitoring instrument (OMI), which measures the absorption of solar backscatter in the UV–visible spectral range. NO₂ can be observed from space because it has strong absorption features within the 400–465 nm wavelength region (Vandaele et al., 1998). By comparing observed spectra with a reference spectrum, the amount of NO₂ in the atmosphere between the instrument in low-Earth orbit and the surface can be derived; this technique is called differential optical absorption spectroscopy (DOAS) (Platt, 1994).

Tropospheric NO₂ column contents from OMI have been used to estimate NO_x emissions from various areas around the globe (Streets et al., 2013; Miyazaki et al., 2017) including North America (Boersma et al., 2008; Lu et al., 2015), Asia (Zhang et al., 2008; Han et al., 2015; Kuhlmann et al., 2015; Liu et al., 2017), the Middle East (Beirle et al., 2011), and Europe (Huijnen et al., 2010; Curier et al., 2014). It has also been used to produce and validate NO_x emissions estimates from sources such as soil (Hudman et al., 2010; Vinken et al., 2014a; Rasool et al., 2016), lightning (Allen et al., 2012; Liaskos et al., 2015; Pickering et al., 2016; Nault et al., 2017), power plants (de Foy et al., 2015), aircraft (Pujadas et al., 2011), marine vessels (Vinken et al., 2014b; Boersma et al., 2015), and urban centers (Lu et al., 2015; Canty et al., 2015; Souri et al., 2016).

With a pixel resolution varying from 13 km × 24 km to 26 km × 128 km, the OMI sensor was developed for global- to regional-scale studies rather than for individual urban areas. Even at the highest spatial resolution of 13 km × 24 km, the sensor has difficulty observing the fine structure of NO₂ plumes at or near the surface (e.g., highways, power plants, factories) (Chen et al., 2009; Ma et al., 2013; Flynn et al., 2014), which are often less than 10 km in width (Heue et al., 2008). This can lead to a spatial averaging of pollution (Hilboll et al., 2013). A temporary remedy, until higher spatial resolution satellite instruments are operational, is to use a regional air quality simulation to estimate the subpixel variability of OMI pixels. Kim et al. (2016) utilize the spatial variability in a regional air quality model to spatially downscale OMI NO₂ measurements using a spatial averaging kernel. The spatial averaging kernel technique has been shown to increase the OMI NO₂ signal within urban areas, which is in better agreement with observations in these regions (Goldberg et al., 2017).

Furthermore, the air mass factor and surface reflectance used in obtaining the global OMI NO₂ retrievals are at a coarse spatial resolution (Lorente et al., 2017; Kleipool et al., 2008). While appropriate for a global operational retrieval, this is known to cause an underestimate in the OMI NO₂ signal in urban regions (Russell et al., 2011). The air mass factors in the operational OMI NO₂ retrieval are calculated using NO₂ profile shapes that are provided at a 1.25^∘ × 1^∘ spatial resolution in the NASA product (Krotkov et al., 2017) and 2^∘ × 3^∘ spatial resolution in the DOMINO product (Boersma et al., 2011). Developers of the NASA product provide scattering weights and additional auxiliary information so that users can develop their own tropospheric vertical column product a posteriori (Lamsal et al., 2015). Several users have recalculated the air mass factor using a regional air quality model (Russell et al., 2011; Kuhlmann et al., 2015; Lin et al., 2015; Goldberg et al., 2017; Laughner et al., 2019), which can better capture the NO₂ profile shapes in urban regions. Other techniques to improve the air mass factor involve correcting for the surface pressure in mountainous terrain (Zhou et al., 2009) and accounting for small-scale heterogeneities in surface reflectance (Zhou et al., 2010; Vasilkov et al., 2017). These a posteriori products have better agreement with ground-based spectrometers measuring tropospheric vertical column contents (Goldberg et al., 2017). When available, observations from aircraft can constrain the NO₂ profile shapes used in the air mass factor calculation (Goldberg et al., 2017).

In this paper, we apply both techniques (the spatial averaging kernel and an air mass factor adjustment) to develop a regional OMI NO₂ product for South Korea. We then use the regional product with only the air mass factor adjustment to derive NO_x emissions estimates for the Seoul metropolitan area using a statistical fit to an exponentially modified Gaussian (EMG) function (Beirle et al., 2011; Valin et al., 2013; de Foy et al., 2014; Lu et al., 2015); the methodology is described in depth in Sect. 2.5.

2 Methods

2.1 OMI NO₂

OMI has been operational on NASA's Earth Observing System Aura satellite since October 2004 (Levelt et al., 2006). The satellite follows a sun-synchronous, low-Earth (705 km) orbit with an Equator overpass time of approximately 13:45 local time. OMI measures total column amounts in a 2600 km swath divided into 60 unequal area “fields of view”, or pixels. At nadir (center of the swath), pixel size is 13 km × 24 km, but at the swath edges, pixels can be as large as 26 km × 128 km. In a single orbit, OMI measures approximately 1650 swaths and achieves daily global coverage over 14–15 orbits (99 min per orbit). Since June 2007, there has been a partial blockage of the detector's full field of view, which has limited the number of valid measurements by blocking consistent rows of data; this is known in the community as the row anomaly (Dobber et al., 2008): http://projects.knmi.nl/omi/research/product/rowanomaly-background.php (last access: 1 February 2019).

OMI measures radiance data between the instrument's detector and the Earth's surface. Comparison of these measurements with a reference spectrum (i.e., DOAS technique) enables the calculation of the total slant column density (SCD), which represents an integrated trace gas abundance from the sun to the surface and back to the instrument's detector, passing through the atmosphere twice. For tropospheric air quality studies, vertical column density (VCD) data are more relevant. This is done by subtracting the stratospheric slant column from the total (tropospheric + stratospheric) slant column and dividing by the tropospheric air mass factor (AMF), which is defined as the ratio of the SCD to the VCD, as shown in Eq. (1):

$$\begin{array}{ll}{\displaystyle }{\mathrm{VCD}}_{\mathrm{trop}}& {\displaystyle }={\displaystyle \frac{{\mathrm{SCD}}_{\mathrm{total}}-{\mathrm{SCD}}_{\mathrm{strat}}}{{\mathrm{AMF}}_{\mathrm{trop}}}},\text{where}\\ \text{(1)}& {\displaystyle }{\mathrm{AMF}}_{\mathrm{trop}}& {\displaystyle }={\displaystyle \frac{{\mathrm{SCD}}_{\mathrm{trop}}}{{\mathrm{VCD}}_{\mathrm{trop}}}}.\end{array}$$

The tropospheric AMF has been derived to be a function of the optical atmospheric and surface properties (viewing and solar angles, surface reflectivity, cloud radiance fraction, and cloud height) and a priori profile shape (Palmer et al., 2001; Martin et al., 2002) and can be calculated as follows (Lamsal et al., 2014) in Eq. (2):

$$\begin{array}{}\text{(2)}& {\mathrm{AMF}}_{\mathrm{trop}}={\displaystyle \frac{\sum _{n=\mathrm{surface}}^{\mathrm{tropopause}}{\mathrm{SW}}_n\times x_n}{\sum _{n=\mathrm{surface}}^{\mathrm{tropopause}}{x}_n}},\end{array}$$

where x is the partial column. The optical atmospheric and surface properties in the NASA retrieval are characterized by the scattering weight, which is calculated by the TOMRAD forward radiative transfer model. The scattering weights are first output as a lookup table and then adjusted to real time depending on observed viewing angles, surface albedo, cloud radiance fraction, and cloud pressure.

We follow previous studies (e.g., Palmer et al., 2001; Martin et al., 2002; Boersma et al., 2011; Bucsela et al., 2013) and assume that scattering weights and trace gas profile shapes are independent. The a priori trace gas profile shapes (x_a) must be provided by a model simulation. In an operational setting, NASA uses a monthly averaged and year-specific Global Model Initiative (GMI) global simulation with a spatial resolution of 1.25^∘ long × 1^∘ lat (∼110 km × 110 km in the midlatitudes) to provide the a priori profile shapes. For this study, we derive tropospheric VCDs using a priori NO₂ profile shapes from a regional WRF-Chem simulation. A full description of this methodology can be found in Goldberg et al. (2017); it is also described in brief in Sect. 2.1.1. We filter the level 2 OMI NO₂ data to ensure only valid pixels are used. Daily pixels with solar zenith angles ≥80^∘, cloud radiance fractions ≥0.5, or surface albedo ≥0.3 are removed as well as the five largest pixels at the swath edges (i.e., pixel numbers 1–5 and 56–60). Finally, we remove any pixel flagged by NASA including pixels with missing values, “XTrackQualityFlags” ≠0 or 255 (row anomaly flag), or “VcdQualityFlags” >0 and least significant bit ≠0 (ground pixel flag).

2.1.1 OMI–WRF-Chem NO₂

We modify the air mass factor in the OMI NO₂ retrieval based on the vertical profiles from a high spatial (4 km × 4 km) resolution WRF-Chem simulation. The vertical profiles are scaled based on a comparison with in situ aircraft observations; this accounts for any consistent biases in the model simulation. We use a campaign mean comparison over all land-based areas (34–38^∘ N, 126–130^∘ E) and scale all modeled profiles in this box by this ratio; there are not enough measurements in any one grid box to scale each individual model grid cell differently. For example, if the aircraft observations during the campaign show that mean NO₂ concentrations between 0 and 500 m are low by 50 %, then we scale all of the modeled NO₂ in this altitude bin by this same amount. To recalculate the air mass factor for each OMI pixel, we first compute subpixel air mass factors for each WRF-Chem model grid cell, using the same method as outlined in Goldberg et al. (2017). The subpixel air mass factor for each WRF-Chem grid cell is a function of the modeled NO₂ profile shape and the scattering weight calculated by a radiative transfer model. We then average all subpixel air mass factors within an OMI pixel (usually 10–100) to generate a single tropospheric air mass factor for each individual OMI pixel. This new air mass factor is used to convert the total slant column into a total vertical column using Eq. (1). Model outputs were sampled at the local time of OMI overpass. For May 2016, we used daily NO₂ profiles and terrain pressures (e.g., Zhou et al., 2009; Laughner et al., 2016) to recalculate the AMF. For other months and years, we used May 2016 monthly mean values of NO₂ and tropopause pressures for the a priori profiles, which are used in the calculation of the AMF.

Once the tropospheric vertical column of each OMI pixel was recalculated, the product was oversampled (de Foy et al., 2009; Russell et al., 2010) for April–June over a 3-year period (2015–2017; 9 months total). During this time frame, there are approximately 9 valid OMI NO₂ pixels per month over any given location on the Korean Peninsula. In the top-down emissions derivation, we use all 9 months of OMI data for the analysis.

2.2 NO₂ observations during KORUS-AQ

We use in situ NO₂ observations from the KORUS-AQ field campaign to test the regional satellite product. KORUS-AQ was a joint South Korea–United States field experiment designed to better understand the trace gas and aerosol composition above the Korean Peninsula using aircrafts, ground station networks, and satellites. The campaign took place between 1 May and 15 June 2016 and measurements were primarily focused in the Seoul metropolitan area. In this paper, we utilize data acquired by the ground-based Pandora spectrometer network, the thermally dissociated laser-induced fluorescence NO₂ instrument on DC-8 aircraft, and the chemiluminescence NO_y instrument on the DC-8 aircraft (${\mathrm{NO}}_y=\mathrm{NO}+{\mathrm{NO}}_{\mathrm{2}}+{\mathrm{HNO}}_{\mathrm{3}}+\mathrm{2}\times {\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}$ + peroxy nitrates + alkyl nitrates + …). KORUS-AQ observations were retrieved from the online data archive: http://www-air.larc.nasa.gov/cgi-bin/ArcView/korusaq (last access: 1 February 2019). A further description of this field campaign can be found in the KORUS-AQ White Paper (https://espo.nasa.gov/korus-aq/content/KORUS-AQ_Science_Overview_0, last access: 1 February 2019).

2.2.1 Pandora NO₂ data

Measurements of total column NO₂ from the Pandora instrument (Herman et al., 2009) are used to evaluate the OMI NO₂ satellite products. The Pandora instrument is a stationary, ground-based, sun-tracking spectrometer, which measures direct sunlight in the UV–visible spectral range (280–525 nm) with a sampling period of 90 s. The Pandora spectrometer measures total column NO₂ using a DOAS technique similar to OMI. A distinct advantage of the Pandora instrument is that it does not require complex assumptions for converting slant columns into vertical columns, compared to zenith sky measurements (e.g., MAX-DOAS). The spatial and temporal variation of NO₂ columns in South Korea as observed by Pandora has been documented elsewhere (Chong et al., 2018; Herman et al., 2018).

In our comparison, valid OMI NO₂ pixels are matched spatially and temporally to Pandora total column NO₂ observations. To smooth the data and eliminate brief small-scale plumes that would be undetectable by a satellite, we average the Pandora observations over a 2 h period (±1 h of the overpass time) before matching to the OMI NO₂ data (Goldberg et al., 2017). During May 2016, there were seven Pandora NO₂ spectrometers operating during the experiment (five instruments were situated within the Seoul metropolitan area and their locations are shown in Fig. 5); this corresponded to 50 instances in which valid Pandora NO₂ observations matched valid OMI NO₂ column data.

2.2.2 DC-8 aircraft data

We compare the model simulation to in situ NO₂ data gathered by the UC-Berkeley Cohen Group (Thornton et al., 2000; Day et al., 2002) on the DC-8 aircraft. The instrument quantifies NO₂ via laser-induced fluorescence at 585 nm. This instrument does not have the same positive bias as chemiluminescence NO₂ detectors, so there is no need to modify NO₂ concentrations by applying an empirical equation (e.g., Lamsal et al., 2008). We also compare the model simulation to chemiluminescence NO_y data gathered by the NCAR Weinheimer group (Ridley et al., 2004).

We utilize 1 min averaged DC-8 data from all 26 flights during May–June 2016. A typical flight path included several low-altitude spirals over the Seoul metropolitan area and a long-distance transect over the Korean Peninsula or the Yellow Sea. The 1 min averaged data are already pregenerated in the data archive. Hourly output from the model simulation is spatially and temporally matched to the observations. We then bin the data into different altitude ranges for our comparison.

2.3 WRF-Chem model simulation

For the high-resolution OMI NO₂ product, we use a regional simulation of the Weather Research and Forecasting model (Skamarock et al., 2008) coupled to Chemistry (WRF-Chem) (Grell et al., 2005) in forecast mode prepared for flight planning during the KORUS-AQ field campaign. The forecast simulations were performed daily and used National Centers for Environmental Prediction Global Forecast System (https://rda.ucar.edu/datasets/ds084.6/, last access: 1 February 2019) meteorological initial and boundary conditions from the 06:00 UTC cycle. Initial conditions for aerosols and gases were obtained from the previous forecasting cycle, while the Copernicus Atmosphere Monitoring Service (Inness et al., 2015) forecasts were used as boundary conditions. WRF-Chem was configured with two domains, with 20 and 4 km grid spacing. The 20 km domain included the major sources for transboundary pollution impacting the Korean Peninsula (deserts in China and Mongolia, wildfires in Siberia and anthropogenic sources from China). The 4 km domain provided a high-resolution simulation where detailed local sources could be modeled and where the KORUS-AQ flight tracks were contained. The inner domain was started 18 h after the outer domain and was simulated for 33 h (00:00 UTC from day 1 to 09:00 UTC of day 2 of the forecast); output was saved hourly. The last 24 h of each inner domain daily forecast over the course of KORUS-AQ was selected to allow spinup from the outer domain and was used in the analysis presented here.

WRF-Chem was configured with four-bin MOSAIC aerosols (Zaveri et al., 2008), a reduced hydrocarbon trace gas chemical mechanism (Pfister et al., 2014) including simplified secondary organic aerosol formation (Hodzic and Jimenez, 2011), and with capabilities to assimilate satellite aerosol optical depth both from low-Earth-orbiting and geostationary satellites (Saide et al., 2013, 2014).

2.4 Emissions inventory

The WRF-Chem simulation was driven by emissions developed by Konkuk University. Monthly emissions for South Korea were developed using the projected 2015 South Korean national emissions inventory, the Clean Air Policy Support System (CAPSS), provided by the National Institute of Environmental Research of South Korea and with enhancements by Konkuk University, which primarily include the addition of new power plants. The projected CAPSS 2015 emissions were estimated based on CAPSS 2012 and 3-year growth factors. Since the base year of the inventory is 2012, observed emissions from the post-2013 large point source inventory were not included. Emissions from China and North Korea were taken from the Comprehensive Regional Emissions for Atmospheric Transport Experiments (CREATE) v3.0 emissions inventory. In order to project the year 2010 emissions to 2015, the latest energy statistics from the International Energy Agency (http://www.iea.org/weo2017/, last access: 1 February 2019) and the China Statistical Yearbook 2016 (http://www.stats.gov.cn/tjsj/ndsj/2016/indexeh.htm, last access: 1 February 2019) were used to update the growth of fuel activities. In addition, the new emissions control policies in China, which were compiled by the International Institute for Applied Systems Analysis, were applied to consider efficiencies of emissions control (van der A et al., 2017).

Emissions were first processed to the monthly timescale at a spatial resolution of 3 km in South Korea and 0.1^∘ for the rest of Asia using SMOKE-Asia (Woo et al., 2012). Information from GIS, such as population, road network, and land cover, was applied to generate gridded emissions from the region-based (17 metropolitan and provincial boundaries of South Korea) emissions. The GIS-based population and regional boundary data compiled by the Ministry of Interior and Safety (https://www.mois.go.kr/, last access: 1 February 2019) and land cover data compiled by the Ministry of Environment (https://egis.me.go.kr/, last access: 1 February 2019) were used to generate population- and land-cover-based spatial surrogates. The Road and Railroad network data compiled by The Korea Transport Institute were used to generate spatial surrogates for on-road and non-road emissions (https://www.koti.re.kr/, last access: 1 February 2019). The emissions were downscaled temporally from monthly to hourly and spatially reallocated to 4 km over South Korea and 20 km over the rest of East Asia using the University of Iowa emissions preprocessor (EPRES).

Biogenic emissions are included using the on-line Model of Emissions of Gases and Aerosols from Nature (MEGAN) version 2; there are no NO_x emissions from MEGAN. For this simulation, the lightning NO_x parameterization was turned off. For wildfires we used the Quick Fire Emissions Dataset (QFED2), but there were only isolated, small fires in South Korea during this time frame.

2.5 Exponentially modified Gaussian fitting method

An EMG function is fit to a collection of NO₂ plumes observed from OMI in order to determine the NO₂ burden and lifetime from the Seoul metropolitan area. The original methodology, proposed by Beirle et al. (2011), involves the fitting of OMI NO₂ line densities to an EMG function. OMI NO₂ line densities are the integral of OMI NO₂ retrieval perpendicular to the path of the plume; the units are mass per distance. We define integration length scale as the across-plume width. The across-plume width is dependent on the NO₂ plume size and can vary between 10 km (for small point sources) and 240 km (for large urban areas). Visual inspection of the rotated oversampled OMI NO₂ plumes is the best way to determine the spatial extent of the emissions sources (Lu et al., 2015).

The EMG model is expressed as Eq. (3):

$$\begin{array}{ll}{\displaystyle }& {\displaystyle }\mathrm{OMI}{\mathrm{NO}}_{\mathrm{2}}\mathrm{line}\mathrm{density}=\\ \text{(3)}& {\displaystyle }& {\displaystyle }\mathit{\alpha }\left[{\displaystyle \frac{\mathrm{1}}{x_o}}\exp \left({\displaystyle \frac{\mathit{\mu }}{x_o}}+{\displaystyle \frac{{\mathit{\sigma }}^{\mathrm{2}}}{\mathrm{2}{x}_o^{\mathrm{2}}}}-{\displaystyle \frac{x}{x_o}}\right)\mathrm{\Phi }\left({\displaystyle \frac{x-\mathit{\mu }}{\mathit{\sigma }}}-{\displaystyle \frac{\mathit{\sigma }}{x_o}}\right)\right]+\mathit{\beta },\end{array}$$

where α is the total number of NO₂ molecules observed near the hotspot, excluding the effect of background NO₂, β; x_o is the e-folding distance downwind, representing the length scale of the NO₂ decay; μ is the location of the apparent source relative to the city center; σ is the standard deviation of the Gaussian function, representing the Gaussian smoothing length scale; and Φ is the cumulative distribution function. Using the “curvefit” function in IDL, we determine the five unknown parameters, α, x_o, σ, μ, and β, based on the independent (distance; x) and dependent (OMI NO₂ line density) variables.

Using the mean zonal wind speed, w, of the NO₂ line density domain, the mean effective NO₂ lifetime τ_effective and the mean NO_x emissions can be calculated from the fitted parameters x_o and α. The wind speed and direction are obtained from the ERA-Interim reanalysis project (Dee et al., 2011), instead of the WRF simulation because the WRF simulation is a forecast. We use the averaged wind fields of the bottom eight levels of the reanalysis (i.e., from the surface to ∼500 m). Only days in which the wind speeds are >3 m s⁻¹ are included in this analysis, because NO₂ decay under this condition is dominated by chemical removal, not variability in the winds (de Foy et al., 2014). The factor of 1.33 is the mean column-averaged NO_x∕NO₂ ratio in the WRF-Chem model simulation for the Seoul metropolitan area during the midafternoon. The NO_x∕NO₂ ratio is time-dependent, spatially varying and is primarily a function of the localized j(NO₂) and O₃ concentration.

$$\begin{array}{ll}{\displaystyle }{\mathrm{NO}}_x\text{emissions}& {\displaystyle }=\mathrm{1.33}\left({\displaystyle \frac{\mathit{\alpha }}{{\mathit{\tau }}_{\mathrm{effective}}}}\right),\text{where}\\ \text{(4)}& {\displaystyle }{\mathit{\tau }}_{\mathrm{effective}}& {\displaystyle }={\displaystyle \frac{x_o}{w}}\end{array}$$

The NO₂ plume concentration is a function of the emissions source strength, wind speed, and wind direction. Originally, the method separated all NO₂ plumes by wind direction and fit an EMG function to NO₂ in eight wind directions (Beirle et al., 2011; Ialongo et al., 2014; Liu et al., 2016). Newer methodologies rotate the plumes so that all plumes are in the same direction (Valin et al., 2013; de Foy et al., 2014; Lu et al., 2015). This process increases the signal-to-noise ratio and generates a more robust fit. In this work, we filter OMI NO₂ data and rotate the NO₂ plumes, as described in Lu et al. (2015), so that all plumes are decaying in the same direction. We rotate the retrieval based on the reanalyzed 0–500 m wind speed direction from the ERA-Interim reanalysis. In doing so, we develop a regridded satellite product in an x–y coordinate system, in which the urban plume is aligned along the x axis. Following de Foy et al. (2014) and Lu et al. (2015), we only use days in which the ERA-Interim wind speeds are >3 m s⁻¹ because there is more direct plume transport and less plume meandering on days with stronger winds; this yields more robust NO_x emissions estimates. We fit an EMG function to the line density as a function of the horizontal distance. This yields a single value at each point along the x direction.

3 Results

In this section, we describe the regional high-resolution satellite product and our validation efforts. All OMI NO₂ results presented here are vertical column densities. First, we show a continental snapshot of OMI NO₂ (OMI-standard) over East Asia using the standard NASA product. Then, we show a regional NASA OMI NO₂ satellite product (OMI-regional) using AMFs generated from the WRF-Chem a priori NO₂ profiles. We compare the OMI-regional product with NO₂ VCDs from the original WRF-Chem simulation. We evaluate the OMI-regional product by comparing to KORUS-AQ observations. Finally, we use the OMI-standard and OMI-regional products to estimate NO_x emissions from the Seoul metropolitan area.

3.1 OMI NO₂ in East Asia

Oversampled OMI NO₂ for May–September 2015–2017 (15 months total) in East Asia and the 4 km WRF-Chem model domain are shown in Fig. 1. The OMI NO₂ signals in East Asia over major metropolitan areas are 3 to 5 times larger than over similarly sized cities in the US (Krotkov et al., 2016). This is in spite of recent NO_x reductions in China since 2011 (de Foy et al., 2016; Souri et al., 2017; Zheng et al., 2018). OMI has observed a recent decrease in the NO₂ burden in the immediate Seoul metropolitan area in South Korea, but an increase in areas just outside the city center (Duncan et al., 2016). Oversampled values greater than 8×10¹⁵ molecules cm⁻² are still consistently seen in East Asia, while they are non-existent in the US during the warm season.

Figure 1 Warm-season averaged (May–September) NO₂ tropospheric vertical column content using the OMI-standard NO₂ product for the years 2015–2017 in East Asia. The 4 km × 4 km WRF-Chem domain is outlined over the Korean Peninsula.

Figure 2 (a) OMI-standard NO₂ product averaged over a 9-month period, April–June 2015–2017, (b) the OMI-regional NO₂ product with only the air mass factor adjustment averaged over the same time frame, and (c) the ratio between the two products: (b)∕(a). (d) Same as (a), (e) the OMI-regional NO₂ product with the air mass factor adjustment and spatial kernel averaged over the same time frame and (f) the ratio between the two products: (e)∕(d).

Figure 3 (a) The OMI-regional NO₂ product with the air mass factor adjustment and spatial kernel averaged during the month of May 2016, (b) the WRF-Chem model simulation showing only days with valid OMI measurements, and (c) the ratio between the two products. On average, there are only nine valid OMI pixels per month observed at any given location on the Korean Peninsula during May 2016.

3.2 Calculation of new OMI tropospheric column NO₂

In Fig. 2, we plot the OMI-standard and OMI-regional products over South Korea. The left panels are identical and show the OMI-standard product for April–June 2015–2017. Figure 2b shows a regional product in which only the air mass factor (AMF) correction is applied. Figure 2e shows a regional product in which the air mass factor correction and spatial averaging kernel (AMF + SK) are applied. The regional product yields larger OMI NO₂ values throughout the majority of the Korean Peninsula. Areas near major cities (e.g., Seoul), power plants, steel mills, and cement kilns have OMI NO₂ values that are >1.25 times larger in the regional AMF product and >2 times larger in the regional AMF + SK product. There are two reasons for the larger OMI NO₂ signals: the air mass factors in polluted regions are now smaller (Russell et al., 2011; Goldberg et al., 2017) and the spatial weighting kernel allocates a large portion of the OMI NO₂ signal into a smaller region (Kim et al., 2016).

Figure 4 The diurnal profile of NO_x emissions rates processed from the bottom-up inventory. (a) The diurnal profile of NO_x emissions rates during a weekday in New York City during July 2011 using SMOKE as the emissions preprocessor (Goldberg et al., 2016). (b) The diurnal profile of emissions rates during a weekday in Seoul during May 2016 using EPRES as the emissions preprocessor. Emissions profiles in the right panel were used in the WRF-Chem simulation.

3.3 OMI-regional vs. WRF-Chem

We now compare the OMI-regional product to tropospheric vertical columns from the WRF-Chem model simulation directly. In Fig. 3, we compare the regional satellite product (AMF + SK) to the WRF-Chem simulation over the Korean Peninsula. In most areas, the modeled tropospheric column NO₂ is of smaller magnitude than inferred by the satellite. In the area within 40 km of the Seoul city center, modeled tropospheric vertical columns are 44 % smaller than observed tropospheric vertical column in the regional AMF + SK product. We posit four reasons as to why the model simulation calculates columns that are consistently smaller. First, our WRF-Chem simulation uses a reduced hydrocarbon gas-phase chemical mechanism. This fast-calculating mechanism implemented in WRF-Chem for regional climate assessments (Pfister et al., 2014) and used during KORUS-AQ for forecasting does not quickly recycle alkyl nitrates back to NO₂; this will cause NO₂ to be too low. While an underestimate of the chemical conversion to NO₂ in WRF-Chem is a contributor to the underestimate, it likely does not account for the entire discrepancy; Canty et al. (2015) suggest that by shortening the lifetime of alkyl nitrates in the chemical mechanism, NO₂ will increase by roughly 3 % in urban areas and 18 % in rural areas. Second, an underestimate in VOC emissions would have an impact on peroxyacyl and alkyl nitrate formation, and should enhance the effective NO_x lifetime (Romer et al., 2016). Third, the temporal allocation of bottom-up emissions inventories can be a very significant source of uncertainty (Mues et al., 2014). The temporal allocation of the bottom-up South Korean NO_x emissions is such that the early afternoon rate during the OMI overpass time (between 12:00 and 14:00 local time) is approximately equal to the 24 h averaged rate (Fig. 4). For comparison, in the eastern US, the early afternoon emissions rate is 1.35 larger than the 24 h averaged emissions rate. Thus, there are scenarios in which the temporal allocation can be up to 35 % different in the midafternoon during the OMI overpass time. We are not suggesting that the South Korean emissions inventory should have the diurnal profile of the US or vice versa, but instead that there are scenarios in which the temporal allocation can vary widely. This substantiates the future use of geostationary satellites to better constrain this temporal allocation uncertainty. Lastly, the remaining difference will likely be due to an underestimate in the emissions inventory.

Figure 5 Measurements from the DC-8 aircraft binned by altitude in black. Co-located WRF-Chem observations within the same altitude bin as the aircraft observations are plotted above in red. Square dots represent the median values. Boxes represent the 25th and 75th percentiles, while whiskers represent the 5th and 95th percentiles. (a) Comparison of NO₂ in the Seoul plume (SW corner: 37.1^∘ N, 127.05^∘ E; NE corner: 37.75^∘ N, 127.85^∘ E), (b) comparison of NO_y in the Seoul plume, (c) comparison of NO₂ in areas outside of the Seoul metropolitan area on the Korean Peninsula (SW corner: 34.0^∘ N, 126.4^∘ E; NE corner: 37.1^∘ N, 130.0^∘ E), and (d) comparison of NO_y in areas outside of the Seoul metropolitan area on the Korean Peninsula.

3.4 Comparing WRF-Chem to aircraft measurements

When comparing the model simulation to in situ observations from the UC-Berkeley NO₂ instrument aboard the aircraft, we find that NO₂ concentrations are substantially larger than the model when spatially and temporally collocated in the immediate Seoul metropolitan area (Fig. 5). The comparison isolates the NO₂ within the lowermost boundary layer as the primary contributor to the tropospheric column underestimate. When comparing aircraft NO₂ to modeled NO₂ in other areas of the Korean Peninsula, the underestimate is smaller.

When comparing the model simulation of NO_y to observations of the same quantity observed from the aircraft, we find a similarly large underestimate. NO_y observed on the aircraft is roughly a factor of 2 larger at all altitudes below 2 km. This suggests that errors in NO₂ recycling (NO₂↔NO_y) are not the main cause of the NO₂ discrepancies seen in the satellite and aircraft comparison (also see Fig. 9). Instead, there must be errors in the NO_y production (i.e., NO_x emissions rates are too low) or removal rates (i.e., NO_y deposition rates are too slow).

3.5 Comparison of OMI NO₂ to Pandora NO₂

To quantify the skill of the regional OMI NO₂ product, we compare the new total NO₂ vertical columns from the satellite product to the same quantities observed by the Pandora instruments. In Fig. 6, monthly averaged observations during May 2016 from the Pandora instrument are overlaid onto the monthly average of the three OMI NO₂ satellite products. The two regional OMI NO₂ products capture the magnitude and spatial variability of monthly averaged NO₂ within the metropolitan region better.

Figure 6 (a) Total vertical column contents from the OMI-standard NO₂ product for May 2016, (b) same quantities from the OMI-regional product with only the air mass factor adjustment (AMF) during the same time frame, (c) same quantities from the OMI-regional product with the air mass factor adjustment and spatial kernel (AMF + SK) during the same time frame, and (d) a comparison between total column contents from the three OMI NO₂ products and Pandora NO₂ during May 2016. An average of Pandora 2 h means co-located to valid daily OMI overpasses are overlaid in the spatial plots.

Figure 7 (a) Bottom-up NO_x emissions inventory compiled for the KORUS-AQ field campaign, (b) the oversampled NO₂ plume rotated based on wind direction for Seoul, South Korea, from WRF-Chem (4 km × 4 km) for May 2016, and (c) NO₂ line densities integrating over the 240 km across-plume width (−120 to 120 km along the y axis) and the corresponding EMG fit. NO_x emissions estimates are shown in units of kt yr⁻¹ NO₂ equivalent and represent the midafternoon emissions rate.

We then compare daily Pandora observations to each daily OMI NO₂ value spatially and temporally collocated with the Pandora instrument (Fig. 6). The Pandora observation is a 2 h mean centered on the midafternoon OMI overpass. The slope of the linear best fit of the standard product is 0.58, indicating that there is a consistent low bias in the satellite product when the Pandora instrument observes large values. A similar result was also found by Herman et al. (2018). The best-fit slope of the OMI-regional product with only the air mass factor adjustment (AMF) is 0.76, and the OMI-regional product with the air mass factor adjustment and spatial kernel (AMF + SK) is 1.07, indicating that the regional products capture the polluted-to-clean spatial gradients best. The correlation of daily observations to the satellite retrievals does not improve between retrievals (OMI-standard: r²=0.57, OMI-regional (AMF): r²=0.57, and OMI-regional (AMF + SK): r²=0.58). The lack of improvement in the correlation suggests that the forecasted WRF-Chem simulation is unable to capture the daily variability of NO₂ plumes better than a longer-term average.

3.6 Estimating NO_x emissions from Seoul

To estimate NO_x emissions from the Seoul metropolitan area using a top-down satellite-based approach, we follow the exponentially modified Gaussian fitting methodology outlined in Sect. 2.5. When fit using the EMG method, the photochemical lifetime and OMI NO₂ burden can be derived. Using this information, a NO_x emissions rate can be inferred.

3.6.1 Validating the EMG method using WRF-Chem

The WRF-Chem simulation can serve as a test bed to assess the accuracy of the EMG method, since the bottom-up emissions used for the simulation are known. For this study, we find that for Seoul, an across-plume width of 160 km encompasses the entire NO₂ downwind plume. Using the NO₂ lifetime, NO₂ burden, and a 160 km across-plume width, we calculate the top-down NO_x emissions rate in the WRF-Chem simulation from the Seoul metropolitan area during the early afternoon (Fig. 7). We find the effective NO₂ photochemical lifetime to be 3.1±1.3 h and the emissions rate to be 227±94 kt yr⁻¹ NO₂ equivalent. Uncertainties in the top-down NO_x emissions are the square root of the sum of the squares of the NO_x∕NO₂ ratio (10 %), the OMI NO₂ vertical columns (25 %), the across-plume width (10 %), and the wind fields (30 %) (Lu et al., 2015). Only the latter three terms are used to calculate the uncertainty of the NO₂ lifetime (Lu et al., 2015).

Figure 8 Panels (a) and (b) represent the oversampled (4 km × 4 km) OMI NO₂ plume from Seoul rotated based on wind direction over a 9-month period, April–June 2015–2017, centered on May 2016. Panels (c) and (d) represent the OMI NO₂ line densities integrating over the 240 km across-plume width (−120 to 120 km along the y axis of a, b) and the corresponding EMG fit. Panels (a) and (c) are from the OMI-standard NO₂ product and (b) and (d) are from the OMI-regional NO₂ product. NO_x emissions estimates are shown in units of kt yr⁻¹ NO₂ equivalent and represent the midafternoon emissions rate.

The NO_x bottom-up emissions inventory calculated using a 40 km radius from the Seoul city center is 198 kt yr⁻¹ NO₂ equivalent. We use a 40 km radius in lieu of a larger radius because an assumption in the EMG method is that the emissions must be clustered around a single point (in this case, the city center). Therefore, the calculated emissions rate from the EMG fit is only measuring the magnitude of the perturbing emissions source, and not of smaller sources that are further from the city center. Previous studies (de Foy et al., 2014, 2015) suggest that the background level calculated by the EMG fit accounts for emissions outside the plume that are more regional and diffuse in nature. The agreement between the top-down (227 kt yr⁻¹) and bottom-up (198 kt yr⁻¹) approaches demonstrates the accuracy and effectiveness of the EMG method in estimating the emissions rate.

3.6.2 Deriving emissions using OMI NO₂

We now calculate the top-down NO_x emissions rate from the satellite data from the Seoul metropolitan area during the early afternoon (Fig. 8). Here we use the OMI standard product and the OMI NO₂ retrieval without the spatial averaging kernel; only the new air mass factor is applied to this retrieval. We do not use the retrieval with the spatial averaging kernel when calculating top-down NO_x emissions because the spatial averaging is strongly dependent on the wind fields in the WRF-Chem simulation, which are forecasted. Errors in the winds can greatly affect the estimate using this top-down approach (Valin et al., 2013; de Foy et al., 2014).

Figure 9 Measurements from the DC-8 aircraft binned by altitude in black. Co-located WRF-Chem within the same altitude bin as the aircraft observations are plotted above in red for the original and in orange for the 2.13×NO_x emissions simulation. Square dots represent the median values. Boxes represent the 25th and 75th percentiles, while whiskers represent the 5th and 95th percentiles. (a) Comparison of NO₂ in the Seoul plume (SW corner: 37.1^∘ N, 127.05^∘ E; NE corner: 37.75^∘ N, 127.85^∘ E), (b) comparison of NO_y in the Seoul plume, (c) comparison of the NO₂–NO_y ratio in the Seoul plume when coincident NO₂ and NO_y measurements are available.

Figure 10 Same as Fig. 3, but now showing the WRF-Chem simulation with NO_x emissions in the Seoul metropolitan area increased by a factor of 2.13 in panel (b).

Figure 11 (a) The OMI-standard product during the month of May 2016, (b) the OMI-regional NO₂ product with the WRF-Chem air mass factor adjustment and spatial kernel during the same period, (c) same as (b) but using WRF-Chem NO₂ profiles scaled based on the aircraft comparison, and (d) same as (b) but using the WRF-Chem simulation with NO_x in the Seoul metropolitan area emissions increased by a factor of 2.13.

For the standard product, the effective NO₂ photochemical lifetime is 4.2±1.7 h, while in the regional product, the effective lifetime is 3.4±1.4 h. In the standard product, we derive the NO_x emissions rate to be 353±146 kt yr⁻¹ NO₂ equivalent, while in the regional product it is 484±201 kt yr⁻¹ NO₂ equivalent. Emissions estimates using satellite products with coarse-resolution air mass factors will yield top-down emissions estimates that are lower than in reality. In this case, the regional satellite product yields NO_x emissions rates that are 37 % higher; we would expect similar results from other metropolitan regions. The top-down approach for the model simulation yielded a NO_x emissions rate of 227 kt yr⁻¹, while the top-down approach using the satellite data yielded a 484 kt yr⁻¹ NO_x emissions rate: a 53 % underestimate in the emissions inventory.

It should be noted that the NO₂ photochemical lifetime derived here is a fundamentally different quantity than the NO₂ lifetime observed by in situ measurements (de Foy et al., 2014; Lu et al., 2015) or derived by model simulations (Lamsal et al., 2010). This is because the lifetime calculation is extremely sensitive to the accuracy of the wind direction (de Foy et al., 2014) and spatial pattern of the emissions. Inaccuracies in the wind fields introduce noise that shortens the tail of the fit. As a result, NO₂ photochemical lifetimes derived here are considered “effective” photochemical lifetimes and are generally shorter than the tropospheric column NO₂ lifetimes derived by model simulations (Lamsal et al., 2010). NO_x sources at the outer portions of urban areas will lead to an artificially longer NO₂ lifetime. This partially compensates for the bias introduced by the wind direction. The heterogeneous topography and oscillating thermally driven wind flows (such as the Yellow Sea breeze) in the Seoul metropolitan area are effects that may bias the effective photochemical lifetime calculation. We partially account for this bias by only selecting days with strong winds (>3 m s⁻¹); on days with faster winds speeds, the sea and mountain breeze effects are secondary to the synoptic flow.

3.7 Model simulation with increased NO_x emissions

To test whether an increase in the NO_x emissions rate is appropriate for the Seoul metropolitan area, we conduct a simulation with NO_x emissions in the Seoul metropolitan area (within a 40 km radius of the city center), which increased by a factor of 2.13, and we analyze the results for May 2016. The 2.13 increase is representative of the change suggested by the top-down method (OMI-regional: 484 kt yr⁻¹ vs. WRF-Chem original: 227 kt yr⁻¹). This simulation was performed slightly differently than the original simulation in that it was a continuous month-long simulation and the outer domain was nudged to the reanalysis.

When comparing the new model simulation to in situ observations from the UC-Berkeley NO₂ and NCAR NO_y instruments aboard the DC-8 aircraft, we find that NO₂ concentrations are a bit high, but NO_y concentrations are in good agreement with WRF-Chem in the boundary layer when spatially and temporally collocated in the immediate Seoul metropolitan area (Fig. 9). The NO₂–NO_y partitioning is captured well by both model simulations, and there is no significant change in the NO₂–NO_y ratio when using increased NO_x emissions.

When comparing the new WRF-Chem simulation to the OMI-regional product for May 2016 (Fig. 10), we now find no significant biases in the Seoul metropolitan area. In the area within 40 km of the Seoul city center, NO₂ columns are now only 11 % smaller in the new model simulation. The better agreement in NO₂ and NO_y from a combination of aircraft and satellite data suggests that an increase in NO_x emissions by a factor of 2.13 is appropriate.

Finally, we reprocess the air mass factors for May 2016 using the newest WRF-Chem simulation. In Fig. 11, we show the OMI-standard product, the OMI-regional product with no scaling of the a priori profiles from the original WRF-Chem simulation, the OMI-regional product with scaling of the original a priori profiles, and the OMI-regional product with a priori profiles from the new WRF-Chem simulation. While using the new a priori profiles increases the OMI NO₂ retrieval further by 8 %, this change is much smaller than the 37 % increase associated with switching models and model resolution (i.e., standard vs. regional product).

4 Conclusions and discussion

In this work, we use a high-resolution (4 km × 4 km) WRF-Chem model simulation to recalculate satellite NO₂ air mass factors over South Korea. We also apply a spatial averaging kernel to better account for the subpixel variability that cannot be observed by OMI. The regional OMI NO₂ retrieval yields increased tropospheric columns in city centers and near large industrial areas. In the area within 40 km of the Seoul city center, OMI NO₂ values are 1.37 larger in the regional product. Areas near large industrial sources have OMI NO₂ values that are >2 times larger. The increase in remotely sensed tropospheric vertical column contents in the Seoul metropolitan area is in better agreement with the Pandora NO₂ spectrometer measurements acquired during the KORUS-AQ field campaign.

Using the regional OMI NO₂ product with only the air mass factor correction applied, we derive the NO_x emissions rate from the Seoul metropolitan area to be 484±201 kt yr⁻¹, while the standard NASA OMI NO₂ product gives an emissions rate of 353±146 kt yr⁻¹. The WRF-Chem simulation yields a midafternoon NO_x emissions rate of 227±94 kt yr⁻¹. This suggests an underestimate in the bottom-up NO_x emissions from Seoul metropolitan area by 53 %, when compared to the 484 kt yr⁻¹ emissions rate from our top-down method. When comparing observed OMI NO₂ to the WRF-Chem model simulation, we find similar underestimates of NO₂ in the Seoul metropolitan area. The effective photochemical lifetime derived in the Seoul plume is 4.2±1.7 h using the standard OMI NO₂ product and 3.4±1.4 h using the regional product. The regional product yields shorter NO₂ lifetimes, although it is not a statistically significant difference. Finally, we show that a WRF-Chem simulation with an increase in the NO_x emissions by a factor of 2.13 yields a better comparison with aircraft observations of NO₂ and NO_y and is in better agreement with the OMI-regional NO₂ product developed herein.

It should be noted that the Seoul metropolitan area has complex geographical features, which adds further uncertainty to this analysis. The area has large topographical changes over short distances, including many hills (>500 m) within the metropolitan area. Furthermore, the city is in close proximity to the Yellow Sea, which causes the area to be affected by sea breeze fronts, especially in the springtime, which is our period of focus. The localized mountain and sea breezes may not be fully captured by our 4 km × 4 km WRF-Chem simulation used to derive the OMI-regional product or the ERA-Interim dataset used to calculate top-down NO_x emissions. The effects of these features on local air quality have been documented elsewhere in the literature (Kim and Ghim, 2002; Lee et al., 2008; Ryu et al., 2013). Nevertheless, the 4 km × 4 km simulation will capture topography and mesoscale phenomena better than a coarse global model and further supports the benefits of WRF-Chem over a global model to derive NO₂ vertical column contents.

We hypothesize that the temporal allocation of NO_x emissions in the bottom-up emissions inventory is a large remaining uncertainty. The satellite-derived emissions rates are instantaneous rates at the time of the OMI overpass (∼ 13:45 local time). This is a different quantity than a bottom-up NO_x emissions inventory, which is often a daily averaged or monthly averaged emissions rate. For this study, we only attempt to derive a midafternoon NO_x emissions rate. Subsequently, we make sure to compare this to the midafternoon NO_x emissions rate from WRF-Chem. While bottom-up studies provide estimates of the diurnal variability of NO_x emissions, these are very difficult to confirm from top-down approaches. Due to a consistent midafternoon overpass time, OMI or TROPOMI cannot address this issue. Due to boundary layer dynamics, this is also very difficult to constrain from ground-based and aircraft measurements. In the future, observations from a geostationary satellite instruments such as the Geostationary Environment Monitoring Spectrometer (GEMS) and Tropospheric Emissions: Monitoring Pollution (TEMPO) will be integral in constraining the ratio of the midafternoon emissions rate to the 24 h averaged emissions rate.