Black Carbon in Bulgaria—Observed and Modelled Concentrations in Two Cities for Two Months

Abstract

Black carbon (BC) is one of the particulate matter (PM) components that both affects human health and contributes to climate change. In this study, we present the preliminary results of the investigation of BC concentrations in PM_2.5 for two Bulgarian cites—Sofia and Burgas. The parallel PM_2.5 samplings were organized in October 2020 and January 2021. The Multi-Wavelength Absorption Black carbon Instrument (MABI) was used for the evaluation of light-absorbing carbon. In addition, we compared the observed BC and PM_2.5 values to modelled ones and analyzed the spatial distribution over the country, using data from advanced operational chemical transport models (CTM)—the European (regional) air quality system established at the Copernicus Atmosphere Monitoring Service (CAMS). Generally, the observed BC and PM_2.5 values were higher in January than in October for both cities. In October, the model underestimated the observed BC concentrations (Sofia—2.44 μg.m⁻³, Burgas—1.63 μg.m⁻³) by 17% and 51%. In January 2021, the observed monthly BC concentrations were higher (Sofia—3.62 μg.m⁻³, Burgas—1.75 μg.m⁻³), and the bias of the model was less than that in October, with an overestimation of 22% for Sofia. The relative bias for PM_2.5 in October (17% for Sofia and −6% for Burgas) was less than the relative bias in January when the model underestimated PM_2.5 monthly mean concentrations by 20% (Sofia) and 42% (Burgas). In addition, we also elaborate on two episodes with high observed BC concentrations in view of the meteorological conditions.

1. Introduction

Air pollution with particulate matter (PM) is still a big problem in the large urban agglomerations in the world as well as in most of the cities of Bulgaria [1]. The urban air quality is determined by many factors: air masses, characteristic local flow and meteorology, low dispersion ability in built-up environments, the concentration of emission sources of different types, and various chemical processes. The PM chemical composition, thus, reflects site-specific complex processes. In the last decade, the studies on the chemical composition of PM in urban areas were widely developed and were associated mainly with the investigation of its effects on human health and understanding of its sources [2,3,4,5].

Studies on the chemical composition of PM in Bulgaria are extremely important. On one side, regular air quality stations of the national monitoring network often report high PM concentrations, whereas, on the other side, there is a lack of highly specialized monitoring sites like the ones of the EMEP (European Monitoring and Evaluation Programme) network [6]. Studies on PM speciation were conducted during field campaigns in the past years for Sofia only [7,8,9], however, without analysis for a significant class of pollutants such as carbon and carbonaceous compounds [10].

The carbonaceous fraction (elemental carbon (EC) and organic carbon (OC)) is a prime constituent of atmospheric aerosol. It generally contributes 20–45% of PM_2.5 (fine particles with an aerodynamic diameter of ≤2.5 µm) and 20–30% of PM₁₀ (particles with an aerodynamic diameter of ≤10 µm) across Europe [5,11]. The examination of carbon particles in urban areas is a matter of great importance for determining the contribution of transportation and domestic heating to air pollution, and thus providing a basis for mitigation measures towards air quality improvement.

There is no universally accepted definition of the term “black carbon” (BC) [12]. The difference between BC and EC is in the measuring method used—BC concentrations are obtained via the optical method, whereas EC—via the thermo-optical method. Not only the shape but also the chemical and physical properties of BC particles can affect the measurements of the optical instruments. Despite that, the average values of BC and EC in urban conditions obtained with different methods appear to be similar [13,14]. Therefore, as we go through the text, we will refer to them using the term BC.

Black carbon, a component of particulate matter, has significant effects on climate and human health [11,12,13,14,15,16,17,18]. BC is typically formed through the incomplete combustion of fossil fuels, biofuel, and biomass and is emitted by anthropogenic and natural sources. It consists of pure carbon in several forms, and the relevant particle size fraction can include known carcinogens and other toxic species [19]. BC stays in the atmosphere for periods varying from days to weeks [11]. Sources of BC can differ widely across regions. The main sources are on-road diesel vehicles, non-road transport, and residential combustion [14,15,16,17,18,19,20,21,22,23,24,25,26,27]. In general, higher concentrations of BC are reported during winter, likely due to higher energy consumption, unfavourable meteorological conditions, and the occasional but significant influence of fires at specific locations during the spring and the summer [21,22,23,24,25,26,27].

Countries around the world have not yet defined an air quality standard for black carbon concentrations and mandatory monitoring. This might explain the relatively low number of urban observational sites, despite some environmental agencies and institutions that monitor BC concentration levels [20]. In Europe, BC measurements were carried out during intensive measurement campaigns at rural sites of the EMEP network [6]. Such observational data represent snapshot information from a set number of locations and for limited periods, but they provide valuable information for evaluating modelling results.

Numerical models that represent the atmospheric dynamics and the chemical transformation of air pollutants are powerful tools to complement sparse data from air quality monitoring networks. The models allow one to analyze the spatial and temporal distribution of pollutants and estimate the effects of different mitigation scenarios for air quality improvement. Chemical transport models (CTM) simulating BC concentrations were applied both at a global scale [28,29] and at a regional scale (over Europe) [30,31]. The uncertainties in modelling results are mainly associated with those in BC emissions, estimated to be a factor of two to five on a regional scale [32]. Another important source of uncertainty is linked to parameterizations of chemical mechanisms describing the black carbon ageing process [33]. The performance of models for BC concentrations in Europe was studied in model intercomparison exercises, comparing results from several CTMs to observations at background sites of the EMEP network [34,35]. The underestimation of EC in PM_2.5 was approximately 20–60% for 2005 [34] and a maximum of 60% for the intensive measurement campaigns in different seasons in the years 2006–2009 [35]. Modelling studies for other parts of the world also showed underestimation of BC by CTMs, for example, in South Africa, mean monthly BC concentrations were underestimated by approximately 50% [36]; in Japan, seasonal BC concentrations observed at urban sites were underestimated by a factor from 0.3 to 0.6; and in South Asia, daily BC concentrations were underestimated by CTM on average by a factor of 4 at urban sites [37,38]

Studies on BC concentration in Bulgaria are limited in number and are not executed regularly. Only a few studies on monitoring BC concentration in PM_2.5 in Bulgaria are available, conducted within the framework of the project “CARBOAEROSOL” [39,40], with an application of methods for estimation of BC from different sources.

In this study, we present and discuss preliminary results on BC concentrations in fine particulate matter. They were observed for months in two main Bulgarian cities: Sofia, found in a valley, and Burgas, located on the Black Sea coast. Furthermore, we compare the observed BC and PM_2.5 values to modelled ones and analyze the spatial distribution over the country, using data from an advanced operational CTM, which is the European (regional) air quality system established at the Copernicus Atmosphere Monitoring Service (CAMS) [41]. Two episodes with high observed BC concentrations are discussed in view of the meteorological conditions.

2. Materials and Methods

2.1. Sampling Sites Description

The PM_2.5 sampling field campaigns were organized in two large cities—Sofia and Burgas (Figure 1)—in the framework of the project “CARBOAEROSOL”. Sofia, located in the western part of the country, is the capital and the largest city in Bulgaria, with a population of 1.3 million inhabitants [42]. It is located in the semi-closed Sofia valley at approximately 550 m a.s.l., surrounded by different mountains, where the topography prevents the dispersion of pollutants and determines unfavourable air quality conditions [10]. The climate is continental with a mean temperature of −0.6 °C in January and a maximum of 21.1 °C in July.

Burgas, located on the Black Sea coast in eastern Bulgaria, is a regional and municipal center. It is the fourth largest town in Bulgaria, housing 277,922 inhabitants [42]. The climate is moderately continental with a clearly expressed sea influence and frequent sea breeze, most notable during the warm period of the year. Mean monthly temperatures are positive all over the year; the mean wind speed is ≥4.0 ms⁻¹ with a maximum of 4.8 ms⁻¹ in March. The city is surrounded by four large water bodies with different salinity, which determines the high relative humidity throughout the year. The winter is mild, usually without any snow. During the winter season, the average air temperature is 4.6 °C, and that of the sea water is 7.4 °C. Due to the influence of the sea, autumn at this site is long and significantly warmer than in inland regions, whereas spring tends to be cold and arrives a month later. These more dynamic atmospheric conditions affect the processes of atmospheric chemistry. The fewer residential, transport, and industrial sources are expected to lead to lower PM_2.5, BC, and nitrogen dioxide (NO₂) concentrations [43,44].

2.2. PM_2.5 and BC Emissions

Annual emissions of PM_2.5 and BC in Bulgaria reported to the Convention on Long-Range Transboundary Air Pollution (CLRTAP) [45] allow the analysis of the trend both of national totals and the contribution by different sectors (Figure 2). The PM_2.5 emissions in 2019 were 30 kt, which is slightly lower than in the previous few years. The main contributors were the sectors “other stationary combustion” (residential heating), with a share of 77%, road transportation, with a share of 7%, and industry (8%). The BC emissions in 2019 were 3.63 kt, representing 12% of the annual PM_2.5 emissions. The annual BC emissions in the period 2008–2019 exhibited little variation, with a maximum of 4.13 kt in 2012. In 2019, 63% of the BC emissions were attributed to residential heating, 33% to road transport. The emissions are unevenly distributed throughout the year. Although the residential sources prevail during winter months (November to February), the transport emissions are almost constant.

There have been no reliable estimations of emissions of particulate matter at city levels so far. Some estimates were made at the two municipalities concerning dispersion modelling for air quality planning. For example, the contribution from residential heating to PM_2.5 emissions in the Sofia region (about 20 km) in 2018 was estimated to be 73% [46]. For Burgas, the share of residential emissions to PM_2.5 in 2014 was estimated to be 68% [47].

A recent study for Sofia, based on the Positive Matrix Factorization (PMF) receptor model for source contribution of different groups of emission sources to PM₁₀ concentrations [10], indicated 23% for biomass burning and 9% for transport on an annual basis for 2019. Another recent study on the determination of the local pollution of PM and PAH in the Sofia municipality by the Gaussian dispersion-based ADMS-Urban model showed that domestic heating contributed 59% of the PM_2.5 concentrations [48].

2.3. Sampling Procedure, Equipment, and Black Carbon Analysis

The parallel 24 h PM_2.5 sampling was performed at the Central Meteorological Observatory, Sofia (42.655 N, 23.384 E, at 586 m a.s.l.) and at “Prof. Dr Assen Zlatarov” University, Burgas (42.314 N, 27.264 E, at 30 m a.s.l.) during autumn (October 2020) and winter (January 2021). Both sampling sites can be classified as ‘urban background’. A total of 69 daily PM_2.5 samples on polytetrafluoroethylene membrane filters (PTFE) Whatman^® (47 mm) were collected with low volume air samplers following the EN-12341 standard (38.3 l min⁻¹). Although the PTFE filters are optically thinner and less homogeneous than the fiber media, they avoid interference from adsorbed organic gases that is associated with quartz and glass fiber media. This provides flexibility to perform subsequent analysis of the same filters using other techniques such as X-ray fluorescence, ICP, IC, etc. The samplers were installed on a grass field, and the sampling heads were at 2 m. a.g.l. Before and after sampling, the filters were weighed by analytical balances (in Sofia with Mettler Toledo, AG135 and in Burgas with KERN & Sohn GmbH, ABT 100-5M). PM_2.5 mass concentrations were determined from the weight differences and sampled air volume.

All collected PM_2.5 samples were analyzed for light-absorbing carbon (LAC) with the Multi-Wavelength Absorption Black carbon Instrument (MABI) at the National Institute of Meteorology and Hydrology (NIMH). This instrument consists of an optical module with a multi-wavelength light source (7 LEDs), a sampling holder, and a photodetector. Moreover, opaque glass is used in the MABI units to scatter the scattered light back through the filter to the detector. The MABI instrument measures light absorption at seven wavelengths, spanning from ultraviolet to infrared (405 nm, 465 nm, 525 nm, 639 nm, 870 nm, 940 nm, and 1050 nm). This allows the differentiation of the contributions from biomass burning (BC_bb) and fossil fuels combustion (BC_ff) based on different light absorption from different type particles at different wavelengths. The values obtained by subtracting the BC (1050 nm) data from BC (450 nm) data represent mainly BC from biomass burning [49,50].

The MABI is a research instrument developed at the Australian Nuclear Science and Technology Organisation that does not perform LAC or BC calculation [49,50]. This is done by the user, which ensures full control, understanding, and interpretation of the obtained BC values. MABI scans are initially performed on unexposed “blank filters” to determine the absorption at each wavelength from the filter substrate. These data for blank filters are referred to as Io. The MABI scans are then repeated on the same filters after exposure and the collection of particulate matter, i.e., “exposed filters”, which represent the absorption at each wavelength from both the collected particles and the filter substrate. This transmission data are referred to as I [49,50]. Both of these values along with the filter area (10.75 cm²) and sampled air volume (~55 m³) are used to determine the black carbon light absorption coefficient (b_abs) in (Mm⁻¹]) at each wavelength (Equation (1)):

$$b_{abs}={10}^2\frac{A}{V}ln\left(\frac{I_o}{I}\right)$$

(1)

where

A—filter collection area in cm²;

V—volume of air sampled through the filter in m³;

Io—measured light transmission through blank (unexposed) filter;

I—measured light transmission through an exposed filter.

The obtained values for b_abs (λ = 639 nm) in this study are in the range 6.74–45.10 Mm⁻¹ for Sofia and 1.19–24.3 Mm⁻¹ for Burgas.

The LAC mass concentration in [ng.m⁻³] is determined using a mass absorption coefficient (ε) in (m².g⁻¹) at each wavelength (Equation (2)):

$$LAC \left({\mathrm{n}\mathrm{g}.\mathrm{m}}^{-3}\right)=\frac{{10}^{5}A}{\epsilon V}ln\left(\frac{I_o}{I}\right)$$

(2)

The mass absorption coefficient (ε) is a function of wavelength (λ): ε = a*λ^−b and was obtained following the procedure presented by [49,50]. A plot of the log of all wavelength data against the log of all λ = 639 nm data is made to determinate the values of ε for each wavelength (example in Figure S1, Supplementary). According to authors [49], the mass absorption coefficient ε at λ = 639 is standard baseline data with a known value—6.6 m².g⁻¹. The received values of ε as a function of λ for Sofia and Burgas are presented in Figure S2 in the Supplementary. These values are in the range of 4.9–8.4 m².g⁻¹, which is in correspondence with literature for fine particulate size fraction PM_2.5 [49,50].

The light source calibration was made automatically in two calibration cycles. During measurements, the calibration was checked at regular intervals.

Filter-based optical methods for estimation of BC concentration suffer from two major artifacts—the filter loading (or shadowing) effect and the multiple scattering effect [51,52,53,54,55,56,57,58,59]. The filter loading effect (FLE) is a bias that affects the apparent concentrations relative to the ambient ones. The multiple scattering from aerosols on the filter and from the filter material leads also to attenuation, which might be misinterpreted as due to absorption.

Different approaches have been developed over time for compensating for these effects, proposing various empirically estimated correction procedures [53,54,55,56,57,58,59,60,61,62]. The correction algorithms depend on the wavelength and make use of parameters that account for instrument details. As there is no worldwide consensus on the most representative correction procedure [60], inter-comparison workshops and studies applying different laboratory and field instruments propose algorithms for widely used instruments (e.g., for AE-31, PSAP).

Most of the correction algorithms were developed specifically for quartz or glass fiber filters. Correcting for loading effect is treated as quite straightforward [57], and it has already been implemented in an aethalometer using a dual-spot technology [58].

In our study we used a relatively new research type instrument (MABI) with a PTFE filter. There are still a limited number of published studies with this device [49,50], and, to our knowledge, they do not have details as to whether it is necessary to apply a filter loading correction.

In a recent work on MABI measurements [49], comparison of data from MABI with data from a MAGEE AE-33 aethalometer (instrument with automatic dual-spot ‘loading compensation’) for almost one year was addressed. The correlation coefficient of R² = 0.96 was reasonable for these two quite different instruments measuring the fine PM_2.5 black carbon concentrations at the same site. On average, the AE_33 aethalometer estimates were 12% higher than the MABI estimates. Furthermore, they quoted each MABI LAC measurement as having an error of ±14% at 639 nm.

The precision of the measurements in our study was verified by conducting repeated measurements on 20 samples. The results showed the difference between measurements on the same filter was <0.1% for all readings of Io and I. The minimum detection limit of LAC ranges from 60 to 70 ng.m⁻³. Considering the presented information in [49] and the [60] approach, the raw data could be further modified if researchers agree upon a correction method that could be used in combination with the MABI instrument. Therefore, in the next sections we will present and discuss raw BC (uncorrected) data.

2.4. Modelling Systems

We analysed data provided by the European (regional) air quality system established at the Copernicus Atmosphere Monitoring Service [63]. The system, further denoted as CAMS-ENS, is based on an ensemble approach, including results from nine state-of-the-art chemistry transport models [63]. All models use common data regarding meteorology, emissions, and boundary conditions. The meteorological driver is the Integrated Forecast system of the European Centre for Medium-range Weather Forecasts (ECMWF-IFS), [64]. The individual models use common anthropogenic emission input, the CAMS-REG-AP emission inventory, developed under the Copernicus Global and Regional emissions service [65]. The reference emission input (CAMS-REG-v.4.2) is for the year 2017, with upgrades in version 4.2 for 2018 and 2019 [66]. There are ongoing activities on including emissions of EC in the PM_2.5 fraction, still facing difficulties due to differences in reporting wood-burning emissions by several countries [67]. All models use common wildfire emission input, provided by the Global Fire Assimilation System at CAMS. The individual models use various chemistry algorithms and parameterizations and assimilation of in situ observations; further details can be found in [67]. The CAMS-ENS results (analysis and 4-day forecast hourly values) are available on a grid of 0.1° × 0.1° at eight vertical levels (from the surface to 5000 m) and refer to more than 20 species. For this study, we downloaded and assessed 1-h analysis data for PM₁₀, PM_2.5, and total elementary carbon from the CAMS atmosphere data store [68].

To investigate the effects of long-range transport, we considered back-trajectories from the Hybrid Single-Particle Lagrangian Integrated Trajectories (HYSPLIT) trajectory model [69,70]. For analysis of selected episodes, back-trajectories arriving at the sampling sites at 500 m, 1500 m, and 3000 m a.g.l. were calculated for three days, using GFS0p25 meteorological input with resolution 0.25° × 0.25°.

2.5. Data Analyses

We analyzed data for PM_2.5 and BC concentrations for October 2020 and January 2021 in Sofia and Burgas. We selected these months to investigate the variability of concentrations from the two cities under different thermal and dynamic conditions. Moreover, October is the month before the start of the heating season, whereas January is in the middle of it. The data for observed black carbon concentrations presented in our study are a sum of fossil fuels (BC_ff) and biomass burning (BC_bb) and are marked as BC.

Furthermore, the obtained data for PM_2.5 and BC concentrations, together with basic meteorological parameters provided by NIMH for both cities, were processed with statistical software R and StatSoft v.6.0.

The observed BC and PM_2.5 mean concentrations were compared with modelled ones for October 2020 and January 2021. Additionally, they were compared for two selected episodes with elevated BC concentrations. For briefly commenting on the synoptic situation during the two selected episodes, we analysed the spatial distribution of some meteorological variables from the global reanalysis dataset of the National Centers for Environmental Prediction/National Center for Atmospheric Research (NCEP/NCAR) [71]. NCEP/NCAR reanalysis products have been widely used to investigate large-scale meteorological conditions and transport paths linked to air pollution episodes. This reanalysis dataset was selected as it provides a possibility for fast visualisation of composite maps. The composite (mean) maps can be produced for selected dates and hours, facilitating evidence of large scale differences in the meteorological fields for low and high pollution days, or for periods in different seasons. We produced composite maps of the geopotential height and the wind at 700 hPa for each of the two selected episodes using the web tool at the National Oceanic and Atmospheric Administration, Earth System Research Laboratories (NOAA/ESRL) [72].

3. Results

3.1. Observed PM_2.5 and BC Concentrations in Sofia and Burgas

The comparison of daily PM_2.5 mass concentrations in Sofia and Burgas for October 2020 and January 2021 is presented in Figure 3. The number of the PM_2.5 samples in Burgas (10) is fewer than that in Sofia (19) for January due to technical problems during the experimental campaign. The PM_2.5 concentrations in Sofia were in the range 4.9–25.5 µg.m⁻³ in October and from 5.1 to 46.9 µg.m⁻³ in January. For Burgas, the PM_2.5 concentrations varied from 3.3 to 27.1 µg.m⁻³ in October and from 2.2 to 30.4 µg.m⁻³ in January.

The analysis showed that 30% of the measured daily PM_2.5 concentrations in Sofia and 43% in Burgas exceeded the most recent WHO air quality guideline level (15 µg.m⁻³) [19]. In October, a higher mean PM_2.5 concentration (14.5 µg.m⁻³) was observed in Burgas than in Sofia (11.8 µg.m⁻³). This is probably due to the fewer rainy days and monthly accumulated precipitation in Burgas (12 mm) in comparison to Sofia (46.7 mm). In January, a higher mean PM_2.5 concentration (18.6 µg.m⁻³) was observed in Sofia than in Burgas (16.2 µg.m⁻³). Although some daily PM_2.5 concentrations in Sofia were higher than the WHO guideline value, the monthly mean PM_2.5 concentration (18.6 µg.m⁻³) was lower than the mean value for January 2019 (22.8 µg.m⁻³), the period before the COVID-19 pandemic [27]. The highest daily PM_2.5 concentrations in Sofia were observed on 21 and 22 January, when the air temperature was around zero and the wind speed was very low (0.2–0.5 ms⁻¹). The highest PM_2.5 concentration in Burgas was measured on 20 January, when the air temperature was 1.9 °C, and the wind speed was 2.9 ms⁻¹.

Figure 4 presents the variation in BC concentrations obtained in October and in January, as well as the BC/PM_2.5 ratio.

Generally, the derived BC values are higher for January than those for October. The BC concentrations in Sofia samples varied from 1.25 to 4.36 µg.m⁻³ for October and from 1.67 to 8.58 µg.m⁻³ for January. The BC concentrations in Burgas samples ranged from 0.37 to 2.76 µg.m⁻³ in October and from 0.42 to 3.36 µg.m⁻³ in January. Based on the measured PM_2.5 and BC concentrations, the derived BC/PM_2.5 ratio for Sofia varies from 17% to 29% and for Burgas, from 1% to 17% in October. This ratio in January covers a wider range in Sofia (14–34%), reflecting the increased impact of combustion sources on the mass concentration of PM and more variable meteorological/synoptic situations. For Burgas the ratio is limited to a narrower range (9–11%) in January, probably because of the similar meteorological conditions in terms of dynamical mixing (high mean diurnal wind speed) for measurement days.

The monthly mean BC concentration for Sofia was lower in October (2.4 µg.m⁻³) than in January (3.6 µg.m⁻³). During the latter, the number of days with inversion was high, and there was an influence from biomass burning for heating. There was no significant difference between the mean BC concentration in Burgas in October (1.63 µg.m⁻³) and January (1.75 µg.m⁻³). This is most likely due to the more intense dispersion of air pollutants and the higher wind speed in Burgas during winter. Moreover, October 2020 was characterized by lower wind speed than the multi-annual (1981–2010) value (norm): in Burgas, it was 3.3 ms⁻¹ (compared to the norm of 4.6 ms⁻¹); and in Sofia, it was 1.0 ms⁻¹ (compared to the norm of 1.2 ms⁻¹). The monthly sum of precipitation was 30.2 mm in Burgas (less than the norm of 57 mm), whereas in Sofia monthly amount of rainfall was slightly higher—60.4 mm with a standard amount of 50 mm. The air temperatures in October 2020 in Burgas and Sofia were above the norm with 3.3 °C and 1.4 °C, respectively. In addition, January 2021 was warmer than usual: 2.8 °C above the norm in Burgas and 2 °C above the norm in Sofia. It was more humid, which is depicted by the precipitation sum—above the norm by 432% in Burgas and 335% in Sofia. The mean wind speed was near the norm for both cities. These specific meteorological conditions led to smaller than expected variations in daily and mean concentrations.

Because the BC and PM_2.5 concentrations depend on the meteorological conditions, a statistical analysis involving main meteorological parameters (temperature, wind speed, relative humidity) was carried out. The correlations between black carbon concentration from fossil fuel (BC_ff) and biomass burning (BC_bb), PM_2.5 mass concentration, and the noted meteorological parameters for both sites are illustrated in Figure 5. A high correlation between BC_ff and PM_2.5 in Sofia was observed, with a statistically significant correlation coefficient equal to 0.9. A moderately high correlation (0.7) between BC_ff and PM_2.5 was observed in Burgas. There were significant negative correlations (−0.5) between BC and the temperature (T) and between BC_ff and wind speed (WS) in Sofia, as expected as high concentrations of pollutants are observed in cold and calm conditions. There is a significant negative correlation between BC_bb concentration and temperature in Sofia (−0.5) and Burgas (−0.6) due to enhanced biomass burning during colder days.

The comparison between the obtained values for BC in Sofia and Burgas and those reported for other European countries is presented in

Table 1. Mean BC values for different European sities.

Place	Period	Mean BC, μg.m−3	Source
Zabrze, Poland (urban)	2019	3.4	[21]
Zabrze, Poland (urban)	2020	2.9	[21]
Zabrze, Poland (urban)	October 2020	4.4	[21]
Zabrze, Poland (urban)	January 2020	6.7	[21]
Athens, Greece (suburban)	October 2013 January 2014	1.8 3.7	[22]
Helsinki, Filand (urban)	October 2015–May 2017	1.7	[23]
Paris, France (urban)	September 2009–September 2010	3.0	[24]
London, UK (urban)	2009	2.0	[24]
Madrid, Spain, (urban background)	Autumn 2015 Winter 2014–2015	3.8 2.4	[25]
Madrid, Spain (urban traffic)	Autumn 2015 Winter 2014–2015	4.9 4.2	[25]
urban station in Germany	2005–2014	2.1	[26]
Ostrava, Czech Republic (urban)	2012–2014	3.5	[27]
Sofia, Bulgaria (urban background)	October 2018	3.8	[40]
Sofia, Bulgaria	January 2019	3.7	[40]
Sofia, Bulgaria (urban background)	October 2020 January 2021	2.4 3.6	This study
Burgas, Bulgaria (urban background)	October 2020 January 2021	1.6 1.8	This study

. It was found that the obtained values for BC in Sofia were close to those measured in Madrid and Athens (January 2020) and lower than those obtained in Zabrze, Poland. The mean BC concentrations in Burgas were generally lower compared to the values recorded at most urban or urban background stations. In contrast, they were closer to those for Helsinki and Athens (October 2020). It should be noted that the applied method and instrument for estimation of BC concentrations in our study are relatively new. According to the authors who created this device, there is a possible underestimation of BC concentration by approximately 12% [49]. The possible underestimation of BC concentration might lead to a higher BC/PM_2.5 ratio in some particular days with high PM_2.5 concentrations in January in Sofia, not more than 5–6% of what was presented originally in Figure 4. New research and results from comparative measurements with different instruments are required for more accuracy in these estimates.

3.2. PM_2.5 and BC Concentrations: Observed vs. Model

In this section, we present the results from a comparison of observed and modelled BC and PM_2.5 concentrations for October 2020 and January 2021. Moreover, we elaborate on selected episodes with elevated values of observed BC concentration. As we noted in Section 2.5, the designation “BC” refers to observed BC and modelled total EC by CAMS-ENS. We must keep in mind that there is some uncertainty in the measurement of BC, estimated at ±14%, which we assume will have a lower effect in averaging the data, compared to single values.

3.2.1. Monthly Mean PM_2.5 and BC Concentrations: Observed vs. Model

The comparison for monthly mean BC and PM_2.5 concentrations at the two sampling sites for October 2020 and January 2021 is presented in Figure 6. The observed BC concentrations in October (Sofia—2.44 μg.m⁻³, Burgas—1.63 μg.m⁻³) are underestimated by the model, respectively, by 17% and 51%. In January 2021, the observed monthly BC concentrations are higher (Sofia—3.62 μg.m⁻³, Burgas—1.75 μg.m⁻³), and the bias of the model is less than in the autumn month. An overestimation of 22% for Sofia modelled BC in January compared to the observed might be partly a result of the higher PM_2.5 concentration and bias in measured BC concentrations because of the filter loading effect.

It has to be noted that model values are representative of the computational grid (about 10 × 10 km) and background concentrations in urban environments located in regions with complex topography, whereas monitored concentrations reflect local emission sources and meteorological and dispersion conditions. The relative bias for PM_2.5 in October (17% for Sofia and −6% for Burgas) was less than the relative bias in January when the model underestimated PM_2.5 monthly mean concentrations by 20% (Sofia) and 42% (Burgas). The PM_2.5 underestimation and, at the same time, BC overestimation in the winter month at the Sofia sampling site suggest that model BC emissions from residential combustion around the capital are high. The underestimation of BC at both locations during autumn, with biases higher than for PM_2.5 values, suggest that BC emitted from all types of transportation might be underestimated. The observations and the model results show that BC and PM_2.5 concentrations were higher in Sofia than in Burgas.

The spatial distribution of monthly mean BC and PM_2.5 concentrations based on CAMS-ENS data is shown over part of the Balkans in Figure 7 for October 2020 and Figure 8 for January 2021.

The regions with higher BC concentrations in both months were mainly outside of the country. The hot spots in the country were around Sofia and in the lowlands of southern Bulgaria, especially in the winter month. Elevated BC concentrations were also noted in the northern part of the country, in the Lower Danube plain, covering wider areas in October. The PM_2.5 pattern was similar to the one for BC, although some new features were noted. The model simulated higher PM_2.5 concentrations in October than in January in the region south-east of Bulgaria. This area includes the megalopolis of Istanbul and the industrialized area of Turkey along the northern Aegean coast. The higher PM_2.5 concentrations were most likely linked to the influence of emission sources in the region and the atmospheric dynamics, more easily noticed in October than under winter conditions.

Model concentrations depend on many factors—emission sources, chemistry mechanisms, and meteorological processes. The black carbon in atmospheric particulates is inert and directly linked to emission sources. The spatial distribution of PM_2.5 emissions used by CASM-ENS is shown in Figure 9. The emission data are available in the Emissions of Atmospheric Compounds and Compilation of Ancillary Data (ECCAD) at [73]. The maps in Figure 9, produced by the ECCAD data tools, refer to the 2019 and CAMS-REG-AP-v4.2 emission inventory on a 0.05° × 0.1° grid. The annual total emissions of PM_2.5 were mainly due to residential combustion, especially in the lowland regions in the northern and southern parts of the country.

3.2.2. Daily Mean PM_2.5 and BC Concentrations: Observed vs. Model

Further, we discuss the comparison between daily BC and PM_2.5 concentrations at the two sites. Figure 10 shows the daily variability of BC and PM_2.5 in October, whereas

Table 2. Statistical parameters for BC and PM_2.5 in October 2020.

October 2020	n	Mean_OBS	Mean_MOD	MBE	RMSE	Corr	FGE	NMB %
BC_Sofia	14	2.44	2.02	−0.42	1.09	0.42	0.43	−17.21
PM2.5_Sofia	14	11.76	13.77	2.00	4.77	0.65	0.31	17.09
BC_Burgas	14	1.63	0.80	−0.83	1.04	0.51	0.67	−50.92
PM2.5_Burgas	14	14.47	13.6	−0.87	3.81	0.87	0.26	−6.01

Mean observed and modelled values: MB—mean bias error, RMSE—root mean square error, Corr—correlation coefficient, FGE—fractional gross error, NMB—normalized mean bias. Units: µg.m⁻³.

and

Table 3. Statistical parameters for BC and PM_2.5 in January 2021.

January 2021	n	Mean_OBS	Mean_MOD	MBE	RMSE	Corr	FGE	NMB %
BC_Sofia	19	3.62	4.67	1.05	1.85	0.46	0.39	29.01
PM2.5_Sofia	19	18.58	14.71	−3.87	8.80	0.65	0.28	−20.83
BC_Burgas	10	1.75	1.70	−0.04	0.55	0.80	0.25	−2.86
PM2.5_Burgas	10	16.16	9.32	−6.84	10.07	0.70	0.47	−42.33

Mean observed and modelled values: MBE—mean bias error, RMSE—root mean square error, Corr—correlation coefficient, FGE—fractional gross error, NMB—normalized mean bias. Units: µg.m⁻³.

present some main statistical parameters.

In Sofia, the observed BC daily concentrations ranged from 1.3 to 4.4 μg.m⁻³, whereas the simulated BC ranged from 0.8 to 3.9 μg.m⁻³, with a correlation coefficient equal to 0.42. For PM_2.5, the observed values were between 4.9 and 25.1 μg.m⁻³; the modelled ones between 7.1 and 26.2 μg.m⁻³, with a correlation coefficient equal to 0.65. In comparison, the correlation coefficients in Burgas were higher—0.51 for BC and 0.87 for PM_2.5.

The daily values in January (Figure 11) were, in general, higher at both sites. In Sofia, the observed BC and PM_2.5 were in the interval 1.7–8.6 μg.m⁻³ and 5.1–46.9 μg.m⁻³, whereas the model ones were in the range 2.6–6.7 for BC and 7.5–21 μg.m⁻³ for PM_2.5. BC was overestimated for almost all days in Sofia and half of the days in Burgas. It is interesting to note that on the days with the highest observed concentration of BC and PM_2.5, the modelled values are lower. The correlation coefficients in Sofia were almost the same as in October: −0.46 for BC and 0.65 for PM_2.5. The correlation coefficients in Burgas were higher (0.78 for BC and 0.7 for PM_2.5). However, one must interpret the numbers for Burgas with caution due to the limited number of daily values.

3.2.3. Daily Mean PM_2.5 and BC Concentrations: Observed vs. Model for Selected Episodes

Episode 1: 23–25 October 2020

Relatively high concentrations of PM_2.5 and BC were observed at both sites on 23 and 25 October 2020 and captured in the model in Figure 12. The CAMS-ENS model simulated an increase in the concentrations on 25 October, especially for PM_2.5 in Sofia.

The synoptic situation in the period 19–24 October 2020 was characterized by high pressure over southern and eastern Europe, favouring the accumulation of pollutants [74].

On 25 October, a perturbation, associated with the passage of a low-pressure center from the region of central Europe towards the Balkans, led to an upper air flow from the south-west to Bulgaria, as shown on the reanalysis maps for the geopotential height and wind at 700 hPa (Figure 13).

HYSPLIT back trajectories (Figure 14) indicated that the upper airflow (3000 m a.g.l.) was from the southern direction at both sites. However, at lower altitudes, due to the persistent high-pressure conditions, the winds were low, and temperature inversions on the surface and aloft were observed in Sofia in the morning and dissolved by midday.

PM_2.5 concentrations simulated by CAMS-ENS in the morning of 25 October 2020 (at 06UTC) were evident at a height of about 1000 m a.g.l (Figure 15a). Along with the plain areas, high PM_2.5 values were noted south-west and south-east (Black Sea) of the country.

Wildfires in the period 23–25 October 2020 were located south and west of Bulgaria, as detected by satellites and mapped by the NASA FIRMS (Fire Information for Resource Management System) service [73] (Figure 15b).

The analyses above suggest that the increased concentrations on 25 October 2010 were most likely due local emission sources in anti-cyclone atmospheric conditions over the country and an approaching perturbation from the west and the south-west, which led to a flow from regions affected by wildfires.

Episode 2: 20–21 January 2021

This episode was characterized by high observed BC and PM_2.5 concentrations in Sofia, with a maximum for the month recorded on 21 January 2021. The observations in Burgas showed the highest values on 20 January 2021 (Figure 16). Contrary to the previous case, the modelled PM_2.5 concentrations were much lower than the observed ones. The observed daily mean BC and PM_2.5 concentrations in Sofia on 21 January 2021 reached values by more than 50% higher than on 20 January 2021, whereas the model indicated less than 10% change for these two days.

The synoptic situation in the period 15–19 January 2021 was characterised by a deep upper trough, intensifying over northern and eastern Europe and leading to a deep surface cyclone with a centre north-east of Bulgaria (southern Ukraine) [75]. The cold spell led to daily mean temperatures in Sofia during this period between −2 °C and −7 °C. On 20 January 2021 and 21 January 2021, a ridge of high pressure from northern Africa extended towards the country; the upper air winds (3000 m a.g.l) were from the west and the north-west (Figure 17). This led to milder conditions and an increase in the daily temperatures of about 1 °C in Sofia. It also resulted in temperature inversions, affecting the air quality in the Sofia valley.

The sounding data for 21 January 2021 in Sofia indicate an inversion on 06 UTC, not only on the surface but also aloft (Figure 18). This elevated inversion layer remained through noon-time.

The bias between model and observed PM_2.5 concentrations is approximately −55% at both sites, indicating that the model has difficulties in capturing the increase in concentrations under low wind conditions.

Modelled PM_2.5 concentrations at 06 UTC on 21 January 2021 (Figure 19) showed higher values in Sofia and over the lower plains of central Bulgaria at approximately 50 m a.g.l. At a height of about 100 m, the accumulation of pollutants around Sofia was not evident.

4. Discussion and Conclusions

For the first time in Bulgaria, fine particles (PM_2.5) were analyzed for BC content. The reported results refer to data from field campaigns (October 2020 and January 2021) carried out in the cities of Sofia and Burgas. The mean contribution of BC to PM_2.5 was estimated to be about 20% in Sofia and 11% in Burgas. However, on some days, BC concentrations might have contributed up to 34% of PM_2.5 in Sofia and up to 19% of PM_2.5 in Burgas. The mean BC concentrations for the two months were 3.03 μg.m⁻³ in Sofia and 1.69 μg.m⁻³ in Burgas. The higher concentrations in Sofia, the largest city in the country, are linked both to more emissions (predominantly from road transport and biomass burning for heating) and to the particular topographic and meteorological conditions in the two cities. In general, wind speeds in Sofia are lower than at seaside Burgas. This affected the dispersion, especially in winter. At the same time, during the winter, the lower temperatures in Sofia, compared to Burgas, led to more emissions from household heating. The effects of the local and regional meteorological conditions on BC and PM_2.5 concentrations were investigated for two episodes with high daily mean BC concentrations (approximately 4.4–8.6 μg.m⁻³ in Sofia). The analysis showed that whereas the episode in October might have been influenced by regional airflow from areas with wildfires, the one in January was associated with temperature inversions.

The comparison of the obtained BC results in Sofia and those reported for other European cites shows values close to those measured in urban background stations in Madrid [25], in a suburban station of Athens [22], and in Ostrava [27]. These values are lower than those obtained in urban background stations in Zabrze, Poland [21]. The mean BC concentrations in Burgas were generally lower compared to the values recorded at most urban or urban background stations and closer to those for Helsinki [23] and Athens [22].

The maps for the spatial distribution of BC and PM_2.5 over the country, constructed using model data by CAMS-ENS [41], revealed common features for the two months in the areas with elevated concentrations around Sofia, the Thracian lowlands of Southern Bulgaria, and the Danube plain in Northern Bulgaria. These areas correspond to the emissions zones where combustion from stationary sources is a prevailing contributor to anthropogenic PM_2.5 emissions. It must be noted that the vast areas with elevated concentrations are outside of the country.

The comparison between observed and modelled concentrations showed a mean bias of approximately 21% (in absolute value) for PM_2.5, with the highest underestimation by—42%—occurring in Burgas during January. In general, BC was underestimated by the model, especially in autumn, by approximately—34%. Other studies on evaluation of model performance for carbonaceous aerosol in Europe [35] have revealed an underestimation of seasonally averaged BC concentrations in the PM_2.5 mass fraction by a maximum of 60%. The observations in [35] had a mean value of 1.5 μg.m⁻³ in the 2007 winter campaign and approximately 1.4 μg.m⁻³ in the 2008 autumn campaign. These values are lower than in our study, but they are representative of background conditions as obtained at sites of the EMEP network.

The reported results and analysis can be considered as the first step Bulgaria is taking towards gaining a better understanding of BC concentration in PM_2.5 in the urban areas of the country. We have applied a variety of tools, analyzing observational data, modelling results, and meteorological parameters. They proved to be helpful in the investigation of elevated BC and PM concentrations during particular episodes.

However, more observational data for different seasons and further studies are needed to draw firm conclusions. A comparative measurement for BC concentration with MABI and AE33 on the territory of Sofia is planned.

We believe that the presented data contribute to filling in the gaps relevant to BC observations in Southeastern Europe.