WRF-Chem simulations of CO₂ over Belgium and surrounding countries assessed by ground-based measurements

The Weather Research and Forecasting model coupled with Chemistry (WRF-Chem), in its passive tracer option (WRF-GHG), was used to simulate CO₂ concentrations over Western Europe during summer 2018. The model performance was evaluated against ground-based observations. Due to the large variety of anthropogenic emissions, we conducted five sensitivity tests using a combination of three different inventories (CAMS-REG-ANT, EDGAR, and TNO) and source-specific vertical emission profiles. Compared with observations from five Integrated Carbon Observation System (ICOS) atmospheric stations, the model captures diurnal CO₂ variations at different heights. At the ICOS site in Karlsruhe, Germany, simulated near-surface CO₂ mole fractions are highly sensitive to the choice of inventory, with discrepancies up to 14.99 ± 31.98 ppm, due to large nearby emission sources. Furthermore, incorporating source-specific vertical profiles notably improves accuracy, increasing the correlation coefficient from 0.53 to 0.78 when using EDGAR. The column-averaged dry-air mole fractions of CO₂ (XCO₂) from the Total Column Carbon Observing Network (TCCON) are well simulated by WRF-GHG. However, an overestimation of approximately 1.2 ppm was found at the Paris site, likely due to uncertainties in background fields and anthropogenic emissions. In addition, a negative bias was found in early June at most ICOS and TCCON sites, may be attributed to errors in simulated fluxes during the growing season. However, due to the lack of co-located flux observations, the exact cause remains uncertain. Overall, this study demonstrates the capability of WRF-GHG in simulating CO₂ over Belgium and surrounding countries and provides insights into the regional-scale characteristics of atmospheric CO₂.

1 Introduction

Intergovernmental Panel on Climate Change Sixth Assessment Report (IPCC AR6) points out that human-induced climate change has significantly influenced the frequency and intensity of extreme events such as heatwaves, heavy precipitation, droughts, and tropical cyclones (Pörtner et al., 2022). In recent years, extreme heat events have become increasingly frequent in western Europe, with prolonged durations of heatwaves (Della-Marta et al., 2007; Sousa et al., 2020; Sánchez-Benítez et al., 2022). The Paris Agreement proposes to limit the global temperature increase to within 2 °C above pre-industrial levels. To achieve this long-term goal, governments must implement measures to reduce carbon emissions. The atmospheric mole fractions of carbon dioxide (CO₂), a major greenhouse gas (GHG), have steadily risen due to human activities over the last centuries. By March 2025, the global mean mole fractions of CO₂ had increased to 426.40 ppm (Lan et al., 2025). Accurate estimation of carbon emissions is a crucial prerequisite for formulating scientifically sound and effective emission reduction strategies.

The 2019 refinement of the 2006 IPCC Guidelines for National Greenhouse Gas Inventories (Maksyutov et al., 2019) explicitly states that the top-down method, based on atmospheric inverse modelling, can serve as a potential way to support and verify national greenhouse gas inventories. However, several studies identified uncertainties associated with atmospheric transport models as one of the main sources of error in this approach (Díaz Isaac et al., 2014; Feng et al., 2016). Therefore, reducing transport errors and improving simulation accuracy are essential for improving the reliability of inversion results.

Regional atmospheric transport models have been widely applied to simulate CO₂ mole fractions, mainly focusing on national or urban scale (Zhao et al., 2019; Zhao et al., 2023; Thilakan et al., 2024; Yang et al., 2025), but their simulation results still present some drawbacks and limitations. Previous studies have shown that the quality of model simulations is highly dependent on the boundary conditions and emission inventories accuracy (Callewaert et al., 2022; Karbasi et al., 2025). Additionally, Brunner et al. (2019) pointed out that using CO₂ emissions only at the surface causes a significant overestimation of the simulated near-surface CO₂ mole fractions, highlighting the need for properly allocating the vertical signature of emissions.

In order to support accurate assessments of greenhouse gas budgets, various ground-based observation networks have been established to provide consistent and high-precision data, such as Integrated Carbon Observation System (ICOS) Atmosphere in Europe. In situ observations provide near-surface CO₂ mole fractions data, which can be used to constrain carbon sources and sinks at local scales. Using the ICOS atmospheric measurements, Ramonet et al. (2020) reported that a severe drought event in Europe in 2018 led to an atmospheric CO₂ signal of +1 to +2 ppm at most stations in summer. Ground-based remote sensing observations, on the other hand, offer information on total column abundances. The global Total Carbon Column Observing Network (TCCON, Wunch et al., 2011) provides long-term, high-precision column-averaged mole fraction measurements that commonly serve as validation data for satellite remote sensing observations (Velazco et al., 2019; Yang et al., 2020; Zhou et al., 2022) and for model verification (Saito et al., 2012; Turner et al., 2015; Ostler et al., 2016). A signal of 0.8 ppm was observed at the Sodankylä TCCON site during the 2018 drought (Ramonet et al., 2020). Since the magnitude of these signals is less than 2 ppm, it is still difficult for satellites to detect such variations (Connor et al., 2016).

Belgium is located in the central part of Western Europe, serving as an important hub connecting France, Germany, the Netherlands, and Luxembourg, and positioned at the intersection of continental transportation. To support climate change mitigation commitments and policy development, the project Towards a greenhouse gas emission monitoring and VERification system for BElgium (VERBE), led by the Royal Belgian Institute for Space Aeronomy (BIRA-IASB), proposes to establish an independent, top-down, temporally and spatially resolved Monitoring and Verification Support (MVS) capacity for greenhouse gas (GHG) emissions in Belgium. To our knowledge, there are no ground-based CO₂ observations, and there aren't any CO₂ model-based studies specifically for the Belgium region currently, while a few studies have been conducted in neighbouring Western European countries (Lian et al., 2021; van der Woude et al., 2023; Zhao et al., 2023).

The study presented in this paper has been performed in the framework of the VERBE project. It employs the Weather Research and Forecasting Greenhouse Gas model (WRF-GHG; Beck et al., 2011) and multi-source observational data to simulate and evaluate CO₂ mole fractions over Western Europe, with a focus on Belgium and surrounding regions during summer 2018. The aim is to evaluate the model performance in this region and analyse the spatial and temporal variations of CO₂ mole fractions. Different sensitivity tests were used to investigate the impact of different anthropogenic emission inventories on the model simulations performance and to optimize the WRF-GHG configuration, in order to improve the regional greenhouse gas simulation accuracy and carbon budget inversion in the near future.

This paper is structured as follows: the introduction to the WRF-GHG model setup and input datasets is given in Sect. 2. The observational datasets and the statistical metrics used to evaluate the model performance are described in Sect. 3. The evaluation of the simulation results, focusing on meteorological fields and CO₂ concentration fields are presented in Sect. 4. The simulation errors looking especially at anthropogenic emissions and biogenic sources are discussed in Sect. 5. Finally, the conclusions are drawn in Sect. 6.

2 WRF-GHG model

The passive tracer option in the WRF model coupled with Chemistry (WRF-Chem), also known as WRF-GHG, can be used to simulate the spatiotemporal distribution of long-lived GHGs like CO₂ and CH₄. In this configuration, the long-lived gases are transported in a passive way without any chemical reactions (Beck et al., 2011). Here WRF-Chem model version 4.5.1 was used.

Figure 1 Terrain elevation map of the simulated domains, with horizontal resolutions of 9 km (d01) and 3 km (d02), showing ICOS (yellow dots), TCCON (orange diamond) and co-located (red stars) sites within d02 (a), and synoptic stations in Belgium for which data are available for our study period (b).

2.1 Model settings

The simulation area covers Western Europe and is centered on Belgium. As shown in Fig. 1a using the Lambert Conformal Conic (LCC) projection, two nested domains are configured, with horizontal resolutions of 9 km for the outer domain and 3 km for the inner domain. The vertical grid consists of 60 hybrid sigma-pressure levels extending from the surface up to 50 hPa. Due to the similarity of the simulation domains, the physical parameterization settings were the same as the settings of Poraicu et al. (2023), except for the urban surface. The Multi-layer Building Effect Parameterization (BEP) model (Martilli et al., 2002) was used instead of the Single-layer Urban Canopy Model (UCM) (Kusaka et al., 2001) because the latter is incompatible with the GHG option. For the Planetary Boundary Layer, we have chosen the Yonsei University scheme (YSU) (Hong et al., 2006).

In the summer of 2018, Europe experienced an intense and widespread heatwave that had a significant impact on ecosystem processes (Bastos et al., 2020; Thompson et al., 2020; Smith et al., 2020). According to records from the Royal Meteorological Institute of Belgium (RMI) at the Uccle station, the total precipitation for the summer (June to August) of 2018 was 134.7 mm, which was substantially lower than the climatological mean for the summer period 1991–2020 (234.2 mm). Meanwhile, this period occurred before the outbreak of the COVID-19 pandemic, during which the ground-based observation network operated normally and provided abundant observational data, laying a solid foundation for model validation and analysis. Therefore, our simulation period covers the summer from 1 June to 31 August 2018. In this study, the meteorological fields are re-initialized every 24 h, with a 6 h spin-up applied before each re-initialization. By doing this, we can constrain the meteorological fields. This approach has been applied in many studies and has proven to improve the simulations accuracy (Pillai et al., 2011; Zhao et al., 2019; Ho et al., 2024).

2.2 Input dataset

In WRF-GHG, the simulated CO₂ (CO_2,total) is the sum of several tracer contributions distinguishing different sources and sinks that are driven by specific emission inventories or flux models, and the so-called background concentrations that are driven by the initial and boundary conditions. Thus, the simulated CO₂ concentrations are given by:

$$\begin{array}{}\text{(1)}& \begin{array}{rl}{\mathrm{CO}}_{\mathrm{2},\mathrm{total}}& ={\mathrm{CO}}_{\mathrm{2},\mathrm{bck}}+C{O}_{\mathrm{2},\mathrm{ant}}+{\mathrm{CO}}_{\mathrm{2},\mathrm{bio}}\\ & +{\mathrm{CO}}_{\mathrm{2},\mathrm{bbu}}+{\mathrm{CO}}_{\mathrm{2},\mathrm{oce}},\end{array}\end{array}$$

where CO_2,bck represents the background concentration, and the remaining terms represent contributions from anthropogenic emissions (CO_2,ant), biogenic activities (CO_2,bio), biomass burning emissions (CO_2,bbu), and ocean-atmosphere exchange (CO_2,oce), respectively. Table 1 gives an overview of the input datasets employed in the simulation, apart from anthropogenic emissions, along with their temporal and spatial resolutions. For initial and lateral boundary conditions, the meteorological fields are provided by the European Centre for Medium-Range Weather Forecasts (ECMWF) global ERA5 hourly reanalysis dataset on model levels, which includes 137 vertical levels (Hersbach et al., 2020), and the chemical fields are provided by the 3-hourly Copernicus Atmosphere Monitoring Service (CAMS) global greenhouse gas reanalysis (EGG4) which includes 60 vertical model levels (Inness et al., 2019). For the emissions, the open biomass burning flux is obtained from the daily Fire INventory from NCAR (FINN v2.5; Wiedinmyer et al., 2023), and the CO₂ fluxes from the oceans are taken from the observation-based global monthly gridded sea surface pCO₂ climatology (NCEI; Landschützer et al., 2019). Given the significant influence of anthropogenic activities throughout much of the study area, the choice of anthropogenic emission inventory can greatly affect the simulation accuracy. The sensitivity tests using different inventories and setups will be discussed later (see Sect. 2.3).

Table 1 Overview of input datasets, excluding anthropogenic emissions.

Biogenic CO₂ flux from the vegetation also plays a significant role during the summer. Here, we use the Vegetation Photosynthesis and Respiration Model (VPRM) (Mahadevan et al., 2008), which is coupled online with WRF-GHG. In VPRM, the calculation of net ecosystem exchange (NEE) consists of two components: gross primary production (GPP) and respiration (R_res). Since vegetation photosynthesis acts as a sink for CO₂, GPP is represented as a negative flux in the WRF-GHG model calculations.

$$\begin{array}{}\text{(2)}& {\displaystyle }\mathrm{NEE}=\mathrm{GPP}+R_{\mathrm{res}},\text{(3)}& {\displaystyle }\mathrm{GPP}=-\mathit{\lambda }{T}_{\mathrm{scale}}\cdot W_{\mathrm{scale}}\cdot P_{\mathrm{scale}}\cdot \mathrm{EVI}\cdot {\displaystyle \frac{\mathrm{PAR}}{\mathrm{1}+\frac{\mathrm{PAR}}{{\mathrm{PAR}}_{\mathrm{0}}}}},\text{(4)}& {\displaystyle }{R}_{\mathrm{res}}=\mathit{\alpha }{T}_{\mathrm{s}}+\mathit{\beta }.\end{array}$$

Here, Photosynthetically Active Radiation (PAR) is assumed to be approximately equal to the downward shortwave radiation (SW) in the WRF-GHG model (i.e., PAR ≈SW). The 2 m temperature (T_s) and PAR are provided by WRF simulations. T_scale represents the temperature sensitivity of photosynthesis, which is defined by a minimum, maximum and optimum temperature (T_min, T_max, T_opt) for photosynthesis for each vegetation class (Evergreen Forest, Deciduous Forest, Mixed Forest, Shrubland, Wetland, Cropland, Grassland and Other). λ, PAR₀, α, β are four parameters that depend on the vegetation class. Considering their importance, and following sensitivity tests (not shown here), these parameters (see Table A1) are adopted from Table 3 in the preprint of Glauch et al. (2025). Notably, in Glauch et al. (2025), they employ the formulation $\mathrm{PAR}\approx \mathrm{0.505}\times \mathrm{SW}$, which differs from the default setting PAR ≈SW in WRF-GHG used here. W_scale and P_scale represent the effect of water stress and leaf age on photosynthesis, respectively. They are both calculated using the Land Surface Water Index (LSWI) (Xiao et al., 2004). Here, the Enhanced Vegetation Index (EVI) and LSWI are derived from the surface reflectance values of 500-m-resolution Moderate Resolution Imaging Spectroradiometer (MODIS) (Huete et al., 2002; Gao, 1996). The Copernicus Dynamic Land Cover Collection 3 (Buchhorn et al., 2020) with a high spatial resolution of 100 m is used to calculate the fraction of each vegetation class in every continental grid cell. The VPRM Preprocessor class in pyVPRM was used to generate the input data needed in VPRM (Glauch, et al., 2025).

2.3 Anthropogenic emission settings

Considering the importance of anthropogenic emissions inventories and the availability of multiple options, we conducted a sensitivity analysis of anthropogenic CO₂ emissions using five different input configurations, based on three different emission inventories, as summarized in Table 2. (1–2) Monthly sector-specific gridmaps of EDGAR v2024, with a spatial resolution of 0.1° × 0.1° (Crippa et al., 2024). (3–4) Yearly CAMS-REG-ANT v8.0 sector-specific gridmaps, with a spatial resolution of 0.1° × 0.05° (Kuenen et al., 2022). (5) Yearly TNO_GHGco_v4.1, with a spatial resolution of $\mathrm{1}/\mathrm{60}$° × $\mathrm{1}/\mathrm{120}$° (∼ 1 × 1 km over central Europe) (Super et al., 2020). As the TNO inventory doesn't encompass the entire simulation domain, emissions in areas outside its coverage are supplemented using the CAMS-REG-ANT inventory. Considering temporal variation, for EDGAR v2024, due to the lack of corresponding sector-specific temporal factors, we assumed constant hourly values within each month. For CAMS-REG-ANT and TNO, we downscaled the yearly fluxes to hourly emissions data, using the sector-specific factors from CAMS-REG-TEMPO (Guevara et al., 2021) and temporal profile factors from Nassar et al. (2013), respectively. In addition to comparing different emission inventories, we also evaluated the impact of accounting for the height of anthropogenic emission point sources. As pointed out by Brunner et al. (2019), more than 50 % of CO₂ emissions in Europe are emitted by large point sources, primarily released through stacks and cooling towers, underscoring the importance of accurately representing the vertical distribution of anthropogenic emissions in model simulations. Assuming all anthropogenic emissions are released from point sources, a vertical disaggregation was applied to the sector-specific emission inventories from (2) EDGAR and (4) CAMS-REG-ANT, using the vertical emission profiles provided by Brunner et al. (2019). Anthropogenic sources are classified into 10 different categories according to the Selected Nomenclature for Air Pollutants (SNAP) in the study by Brunner et al. (2019) on vertical profiles, and the sector mapping between different inventories is detailed in Table A2.

Table 2 Overview of the five different anthropogenic emissions inputs. The type “S” represents all emissions released at the surface. “P” represents anthropogenic emissions released according to source-specific vertical profiles (Brunner et al., 2019).

3 Observations and methods

We collected observational data for both meteorological and chemical fields, consisting of three types: meteorological observations, in-situ and ground-based remote sensing CO₂ observations. Unfortunately, CO₂ concentration measurements are lacking within Belgium during this period, therefore, we gathered observational data from the surrounding regions of Belgium.

3.1 Synoptic observations in Belgium

A meteorological observation network is operated across Belgium by the RMI, Meteorological Wing of the Air Component of Defense (Meteo Wing), and the Belgian Authority of airways that is the Belgian air navigation and traffic service provider for the civil airspace (Skeyes). To ensure data quality, the synoptic data provided by RMI undergo a quality control (QC) procedure consisting of an automatic process followed by manual supervision, whereas the QC of data from stations belonging to Skeyes and Meteo Wing is conducted independently of RMI (Bertrand et al., 2013). Figure 1b shows the locations of the 21 stations from which observations are currently available, where red triangles represent 13 stations operated by RMI, blue diamonds represent 7 stations operated by Skeyes, and the orange dot represents 1 station operated by Meteo Wing. Their detailed coordinates are listed in Table A3. These sites are mainly situated on areas covered by short grass and provide hourly near-surface weather parameters, including temperature, wind speed and wind direction. To evaluate the WRF-GHG model performance, the simulated 2-meter temperature (T2) and 10 m wind speed (WS10) and wind direction (WD10) from the grid nearest to each station within the inner domain were compared with the observations. It's worth noting that there are eight stations operated by Skeyes and Meteo Wing where wind speed and wind direction observations are recorded only as integer values.

3.2 ICOS – Atmospheric stations

The ICOS (ICOS RI, 2023) atmospheric observation network covers the European region and provides standardized, high-precision scientific data on the carbon cycle and greenhouse gas budgets. It currently comprises 46 stations in 16 countries and there are 5 sites located with available observational data within our inner domain, including Observatoire Pérenne de l'Environnement (OPE), Saclay (SAC), Cabauw (CBW), Trainou (TRN), and Karlsruhe (KIT). Their locations are shown in Fig. 1a, with details on site coordinates and observation heights listed in Table 3. Each site is equipped with meteorological instruments and Picarro CRDS (Cavity Ring-Down Spectroscopy) GHG analyzers installed at multiple heights on tall towers, providing measurements of hourly meteorological parameters and in situ CO₂ mole fractions. For comparison with the simulation results, the simulation data from the inner domain nearest grid cell for each observation site are first extracted. Given that the observation heights at each site do not correspond directly to the model levels, the extracted model data are subsequently interpolated to the corresponding observation heights.

Table 3 The relevant information for each station. “m a.s.l.” stands for meters above sea level, and “m a.g.l.” stands for meters above ground level.

3.3 TCCON

TCCON is a global ground-based observation network of Fourier Transform Spectrometers (FTS). It uses the GGG software to retrieve gases mole fractions with high precision and is currently widely used for the validation of satellite measurements (Zhou et al., 2016; Karbasi et al., 2022; Wu et al., 2018). There are three observation sites, Orléans (47.97° N 2.11° E), Karlsruhe (49.10° N 8.44° E), and Paris (48.85° N 2.36° E), located within the inner domain. Among them, the Orléans and Karlsruhe sites are co-located with the ICOS TRN and KIT sites, respectively, and the Paris site is located in an urban area. Each site is equipped with a Bruker IFS 125HR instrument to record shortwave infrared (SWIR) spectra and use the GGG2020 code to retrieve the column-averaged dry air mole fractions of CO₂ (XCO₂) (Laughner et al., 2024). TCCON observations are limited to daytime and clear-sky only. To ensure a meaningful comparison with WRF-GHG outputs, the observational data within a 30 min interval before and after the corresponding model time step are averaged. Additionally, a smoothing correction to account for the a priori profile and averaging kernels (AVKs) associated with the TCCON data was applied to the simulations data before comparison (Rodgers and Connor, 2003). Details on this correction can be found in Appendix B1 of Callewaert et al. (2022).

3.4 Evaluation metrics

To evaluate the performance of the WRF-GHG model, we employed several statistical metrics. The mean bias error (MBE) quantifies the systematic bias between simulations and observations, the standard deviation (SD) of the simulation-observation differences reflects the variability of the simulations relative to the observations, the root mean square error (RMSE) of the difference quantifies the overall magnitude of simulation uncertainties, and the Pearson correlation coefficient (R) reflects the strength of the expected linear relationship between simulated and observed values. These metrics have been widely used in the assessment of model simulations (e.g. Callewaert et al., 2023; Yarragunta et al., 2025). Their calculation formulas are as follows:

$$\begin{array}{}\text{(5)}& {\displaystyle }{\mathrm{Diff}}_i={\mathrm{mod}}_i-{\mathrm{obs}}_i,\text{(6)}& {\displaystyle }\mathrm{MBE}={\displaystyle \frac{{\sum }_{i=\mathrm{1}}^N{\mathrm{Diff}}_i}{N}},\text{(7)}& {\displaystyle }\mathrm{SD}=\sqrt{{\displaystyle \frac{\mathrm{1}}{N}}{\sum }_{i=\mathrm{1}}^N{\left({\mathrm{Diff}}_i-\mathrm{MBE}\right)}^{\mathrm{2}}},\end{array}$$

$$\begin{array}{}\text{(8)}& {\displaystyle }\mathrm{RMSE}=\sqrt{{\displaystyle \frac{\mathrm{1}}{N}}{\sum }_{i=\mathrm{1}}^N{\mathrm{Diff}}_i^{\mathrm{2}}},\text{(9)}& {\displaystyle }R={\displaystyle \frac{{\sum }_{i=\mathrm{1}}^N\left({\mathrm{mod}}_i-\stackrel{\mathrm{‾}}{\mathrm{mod}}\right)\left({\mathrm{obs}}_i-\stackrel{\mathrm{‾}}{\mathrm{obs}}\right)}{\sqrt{{\sum }_{i=\mathrm{1}}^N{\left({\mathrm{mod}}_i-\stackrel{\mathrm{‾}}{\mathrm{mod}}\right)}^{\mathrm{2}}}\cdot \sqrt{{\sum }_{i=\mathrm{1}}^N{\left({\mathrm{obs}}_i-\stackrel{\mathrm{‾}}{\mathrm{obs}}\right)}^{\mathrm{2}}}}},\end{array}$$

where mod_i represents the WRF-GHG simulated values, obs_i represents the observed values, and N is the number of data pairs.

4 Model performance

4.1 Meteorological fields

The overall evaluation metrics between simulated and observed near-surface temperature, wind speed, and wind direction across all 21 synoptic observation stations in Belgium are given in Table 4, and N represents the number of data pairs. The detailed values at each station can be found in Table A3, and the map is shown in Fig. A1. As mentioned in Sect. 3.1, the wind direction and wind speed observations at eight stations are recorded only as integer values, which limits the ability to capture fine-scale variations and primarily reflect overall trends. Despite the lack of high-precision observations, the WRF model reproduces these variation patterns reasonably well at the trend level. The near-surface temperature was well simulated by the model, with MBE of −0.07 K and R as high as 0.95. Regarding wind fields, the observed wind speed and direction exhibit more pronounced fluctuations than the simulations, as also reported by Zhao et al. (2019). This may result from the rapid changes in real atmospheric conditions, which are difficult for models to capture. Besides, as shown in Fig. A1, the model tends to underestimate wind speed values significantly along coastal regions, with an MBE of −1.11 m s⁻¹ across the four coastal stations, which is probably due to the coastal effects in the WRF simulation (Hahmann et al., 2015). In inland areas, the model tends to overestimate at most sites with dense vegetation cover. Such overestimation of wind speed has also been found in previous studies (Duan et al., 2018; Liu et al., 2022; Che et al., 2024). This bias may be attributed to the complex wind distribution in areas with rugged terrain, where the WRF model fails to adequately account for the additional resistance effects of vegetation on unresolved terrain, ultimately leading to an overestimation of wind speed. While an underestimation of wind speed is also observed at some inland stations, most of these stations are located at airports or in areas dominated by cropland. The land-use classification of the corresponding model grid cells is mostly identified as urban and built-up, which may lead to an overestimation of the surface roughness near the station and consequently an underestimation of wind speed. Similar findings were reported by Aylas et al. (2020), who found that WRF underestimated wind speed at airport stations by approximately −0.124 m s⁻¹, while the bias was reduced to −0.079 m s⁻¹ after updating the land-use data.

Table 4 Evaluation metrics between simulated and observed temperature, wind speed, and wind direction across all 21 stations.

In comparison with ICOS observations, there were no meteorological observations available at KIT and CBW stations, as well as at the 100 m height at TRN site during the simulation period. For the remaining sites, the diurnal cycles of temperature and wind speed, along with the time series of wind direction, are shown in Fig. 2. The values in each subplot represent the corresponding MBE ± SD values between the simulations and the observations. In the diurnal variation plots, the solid lines represent the mean values at the same time of day throughout the simulation period, while the shaded areas indicate the standard deviation. Similar to the RMI observation sites, the model simulates temperature well at different heights across all ICOS stations, with MBE of less than 0.22 K. As for wind speed, the model captures the diurnal variation at each site, showing higher wind speeds at night and lower speeds during the day. However, it tends to overestimate the wind speed values, especially as the height gets closer to the surface, such as at 50 m at TRN and OPE stations, and at 60 m at SAC station. This agrees with Tuccella et al. (2012), who pointed out that the model can capture the upper-level wind speed profile but tends to overestimate it in the lower layers. Besides the model tends to significantly overestimate wind speed during the night, as also shown in previous studies (Zhang and Zheng, 2004; Ngan et al., 2013; Dayal et al., 2020). These wind speed simulation errors may be due to limitations of some parameterization of physical processes during the daytime or nighttime, such as turbulence and surface roughness. Additionally, at the OPE and SAC stations, low-level jets (LLJs) frequently occur near the top of the nocturnal boundary layer (NBL), where excessive downward momentum transport may also lead to overestimated wind speeds (Zhang and Zheng, 2004). For wind direction, the model can simulate it well, except at the height of 50 m at OPE site. We additionally compared the wind direction observations at the 50 m height at the OPE site, obtained from an independent ICOS instrument through personal communication with the Principal Investigator (PI), with the model simulation results (not shown). Compared to the ICOS data, this independent observation shows better agreement with the model. There might be uncertainties or measurement errors in the wind direction data from the ICOS anemometer at the 50 m height, which could stem from the sensor itself or be related to its setup, possibly affected by disturbances from the tower structure. The wind speed measurements at the 50m level at OPE are also most probably affected by these limitations.

Overall, the model captures the near-surface variations in meteorological fields well, and exhibits high accuracy in the vertical profiles, which are similar to the performances displayed in previous studies (Tuccella et al., 2012; Mar et al., 2016; Zhao et al., 2023).

Figure 2 The diurnal cycles of temperature (a–f) and wind speed (g–l), and the time series of wind direction (m–r) at different heights at various ICOS stations. The time is local time (UTC+2).

4.2 Chemical fields

The simulated CO₂ mole fractions are compared with corresponding observations in this section.

4.2.1 Comparisons with ICOS – near surface CO₂ mole fractions

The statistical metrics for differences in near-surface CO₂ mole fractions between simulations using five different anthropogenic emission settings and observations at different heights above ground across various ICOS sites are given in Fig. 3. At most sites, the values of SD and RMSE between the simulations and observations tend to increase at lower heights. Additionally, when anthropogenic emissions are released only at the surface, the differences in simulation results between emission inventories can be significant. For example, using surface emission from either EDGAR or TNO can lead to differences up to 14.99 ± 31.98 ppm (MBE ± SD) at the KIT site at 30 m height. It reflects that the anthropogenic emission inventory has a significant impact on the simulation of near-surface CO₂ mole fractions, especially at lower heights above ground-level and for peri-urban stations.

Figure 3 The MBE (a), SD (b), RMSE (c), and R (d) of near surface CO₂ mole fractions between five different simulations and observations at different heights at various ICOS stations. The colors in (a) represent values from negative to positive, with a blue–white–red gradient: negative values in blue, positive values in red, and values near zero in white, whereas in (b)–(d), more intense colors represent larger values and less intense colors represent smaller values.

At the KIT site, taking into account the vertical distribution of anthropogenic emissions has a significant impact on the simulation results. Simulations using elevated emissions show much better agreement with observations, especially for the lower levels. At the height of 30 m, when using the CAMS-REG-ANT emission inventory and considering elevated emissions, the model-observation bias has a low MBE of −5.32 ppm, a SD of 14.37 ppm, a RMSE of 15.33 ppm, and high correlation R of 0.79. When only surface emissions are considered, MBE, SD, and RMSE increase significantly to 11.08, 37.07, and 38.69 ppm, respectively, and R drops to 0.57. The simulations using EDGAR v2024 display similar characteristics, with even larger differences between elevated emissions and surface emissions. As shown in Fig. 4i–l, compared to simulations accounting for elevated emissions, the simulations that only consider surface emissions tend to overestimate CO₂ mole fractions at 30, 60 and 100 m heights. Their diurnal cycles (Fig. 4a–h) indicate that all five simulations are able to capture the diurnal variation of CO₂ and reproduce the lower CO₂ mole fractions well in the afternoon, but the simulation considering only surface emissions significantly overestimates higher CO₂ mole fractions in the morning, leading to large discrepancies with the observations. This overestimation is especially pronounced when using the EDGAR inventory. However, at the other four observation sites (TRN, SAC, OPE, CBW), the differences among the five simulation results are relatively small, and the simulations considering only surface emissions do not show a significant overestimation (see Fig. A2). The discussion in Sect. 5.1 will provide a more detailed analysis of this phenomenon.

Figure 4 Diurnal cycles (local time) of simulations with different anthropogenic emissions and observations (a–h), where the values represent the MBE between each simulation and observations, along with scatterplots comparing each simulation to the observations (i–l) at different heights at the ICOS KIT site.

As shown in Fig. 5, the contributions from biomass burning and oceanic sources to the diurnal cycles are negligible at these five sites. At the KIT, SAC and CBW sites, the diurnal variability of near-surface CO₂ mole fractions is mainly driven by local anthropogenic emissions and biogenic processes, with the biogenic signal systematically stronger. The diurnal variation of anthropogenic CO₂ is attributed to planetary boundary layer (PBL) dynamics and regional transport. Between 15:00 and 17:00 in the afternoon, the PBL (see Fig. A3) reaches its maximum height, coinciding with the minimum anthropogenic CO₂ mole fractions. After sunset, the PBL gradually decreases and remains shallow from midnight until sunrise, leading to the accumulation of anthropogenic emissions. The mole fractions of CO₂ increase during this period and typically peak between 06:00 and 10:00 in the morning. As the PBL rises again after sunrise and turbulent mixing emerges, CO₂ mole fractions begin to decrease, forming a distinct diurnal cycle. The diurnal variation of biogenic processes is driven by vegetation photosynthesis during the day and respiration at night. Biogenic processes are the dominant drivers of the CO₂ mole fractions at the OPE and TRN sites. These results align with Storm et al. (2023), which found that in the summer of 2020, the biogenic flux signals at most ICOS sites were stronger than anthropogenic emissions, and the largest signal is associated with cropland.

Figure 5 Diurnal cycles (local time) of simulated tracer contributions at 30 m at the KIT (a), 50 m at the TRN (b), 60 m at the SAC (c), 50 m at the OPE (d) and 207 m at the CBW (e) sites. Here, the anthropogenic emissions are based on EDGAR v2024, taking into account the vertical emission profiles.

4.2.2 Comparison with TCCON – XCO₂

Figure 6b shows the time series of observed and simulated XCO₂ at three TCCON sites, Fig. 6c is the corresponding scatterplot, with N representing the number of data pairs, and Fig. 6a shows the differences between simulations and observations at each site. Here the simulations are based on the TNO inventory. Due to the lack of observational data at the Paris and Karlsruhe sites in June, the number of valid data pairs for comparison with the simulations is only 123 and 243, respectively. It is obvious that at the Orléans site, the simulated XCO₂ values show a significant underestimation in early June, which will be discussed in detail in the following section.

The statistical metrics between the simulated XCO₂ using five different anthropogenic emission inputs and observed XCO₂ at the three TCCON sites are given in Fig. 7. Among the three sites, the choice of emission inventory has the largest impact on the simulations at the Paris site. At this site, the difference in simulated XCO₂ between using the EDGAR and TNO inventories reaches 0.50 ± 0.34 ppm. Among the five sensitivity tests, all the simulations show an overestimation, with MBEs around 1.2 ppm, and up to 1.65 ppm when using TNO. An inspection of the different tracer contributions to XCO₂ (see Fig. 8) indicates that despite uncertainties in the estimation of biogenic carbon fluxes, the overestimation of background XCO₂ from the CAMS and anthropogenic emissions plays a dominant role in causing biases in the simulated XCO₂. During the period from late July to early August, the simulations show an overestimation, consistent with the relatively large contribution from anthropogenic tracers during that period. These results indicate that the uncertainty in anthropogenic emission inventories remains relatively high in urban areas (Gately and Hutyra, 2017; Super et al., 2020).

Figure 6 Time series (local time) of simulated XCO₂ using TNO inventory and observed XCO₂ at three TCCON sites (b), their absolute differences (WRF-GHG – TCCON) (a), and scatterplots (c). “N” represents the number of data pairs, “a” represents the slope of the linear regression curve and “R” represents the correlation coefficient.

Additionally, as XCO₂ is less sensitive to vertical transport processes (Wunch et al., 2011), one expects less impact of considering elevated anthropogenic emissions heights on the simulation results for XCO₂. Indeed, at the Orléans and Paris sites, the impact is negligible for a given emission inventory, while at the Karlsruhe site, we observe a small improvement in the simulation of XCO₂ but not as notable as that observed for near-surface mole fractions. At the Karlsruhe site, when accounting for elevated emission heights in the CAMS-REG-ANT inventory, the SD and RMSE of differences between simulations and observations improves from 1.13 to 1.05 and 1.06 ppm, respectively while R increases from 0.44 to 0.47; for EDGAR, SD and RMSE decrease from 1.16 to 1.03 and 1.04 ppm, respectively, while R increases from 0.43 to 0.48.

Figure 7 The MBE (a), SD (b), RMSE (c), and R (d) of XCO₂ between five different simulations and observations at three TCCON stations.

Figure 8 Time series (local time) of (a) observed and simulated total, background, and (b) tracer-specific XCO₂ at Paris site. The simulated values here are all without AVK smoothing.

5 Discussion

We will first discuss the anthropogenic emissions sensitivity tests results and then focus on the biogenic fluxes impacts on model biases reported in June 2018.

5.1 Impact of taking into account the height of anthropogenic emissions

As previously noted, at the KIT site, whether anthropogenic emissions are released according to source-specific vertical profiles or only at the surface has a significant impact on the simulation of near-surface CO₂ mole fractions, and also shows a slight influence on XCO₂. Such large impacts are not observed at the other sites.

Figure 9 The maps of sector-specific CO₂ emissions for 2018 from CAMS-REG-ANT (a) public power, (b) industry and (c) fugitives sectors, and EDGAR v2024 (d) power industry, (e) industrial combustion and (f) fuel exploitation. All values are shown on a logarithmic scale (base 10, unit: Tg).

Figure 9 shows the 2018 emission maps of the major contributing sectors from two inventories, including the Public Power (SNAP 1), Industry (SNAP 3) and Fugitives (SNAP 5) sectors from CAMS-REG-ANT, and the Power Industry (SNAP 1), Industrial Combustion (SNAP 3) and Fuel Exploitation (SNAP 5) sectors from EDGAR. Although monthly EDGAR emission data for 2018 were used in the simulations, in these maps the monthly EDGAR inventories were aggregated by summing over all months to produce an annual emission inventory for 2018 for enabling comparison with the CAMS-REG-ANT inventory. It is clearly evident that, except for Industrial Combustion, there are large emission sources near the KIT site, as well as near the CBW and SAC sites, whereas there are almost no emission sources in the vicinity of the TRN and OPE sites. Figure A4 shows the diurnal cycles of the five near-surface CO₂ mole fractions at each ICOS site as simulated by WRF-GHG fractions at the lowest model level (approximately 25 m above ground level). It can be found that at the TRN and OPE sites, the five simulations show nearly identical patterns, consistent with the relatively weak anthropogenic emissions in their surrounding areas. In contrast, at the SAC and CBW sites, whether vertical emissions are taken into account does indeed affect the simulations. However, at the observation heights of these two sites, the impact is small and does not exhibit characteristics similar to those observed at the KIT site. For the SAC site, this may be due to the relatively weaker emission sources nearby, whereas for the CBW site, this may be related to the relatively high observation height used (207 m above ground level), as shown for the KIT site (Fig. 4), the discrepancies between sensitivity tests decrease with increasing observation height.

According to Google Maps, we found that approximately 6.5 km southwest in a straight line from the KIT observation site lies the largest oil refinery of Germany (Junkermann et al., 2011), and about 45 km north of the site there is a gas-fired combined heat and power (CHP) plant located within the Badische Anilinund Sodafabrik (BASF) chemical production facility in Ludwigshafen, Germany. These two are likely the main contributors to the emissions in the inventory near KIT site.

In addition, based on footprint simulations from the regional Stochastic Time-Inverted Lagrangian Transport (STILT) model (Lin et al., 2003), we calculated the aggregated footprints for each afternoon hour (14:00 LT) over the simulation period (not shown here). The results show that, centered on the KIT site, at least 80 % of CO₂ enhancements from anthropogenic emissions are covered within a radius of approximately 1.5° . Therefore, we calculated the CO₂ emission fluxes from the sectors corresponding to those in Fig. 9 in 2018, within the area surrounding the KIT site, bounded by 47.6 to 50.6° N and 6.94 to 9.94° E (see Fig. 10a). For consistency, the SNAP sector classification is used here. The vertical emission profiles applied for the corresponding sectors are given in Fig. 10b, and the complete vertical profiles for all sectors can be found in Fig. 2a of Brunner et al. (2019). For the EDGAR inventory, emissions from sector SNAP 1 are the dominant source, with the emission flux in 2018 as high as 52.18 Tg, and for the CAMS-REG-ANT inventory, emissions from sectors SNAP 1 and SNAP 3 are the main sources, with fluxes of 24.39 and 32.35 Tg, respectively. According to their respective vertical emission profiles, emissions from SNAP 1 and SNAP 3 are primarily concentrated at higher altitudes. In particular, the emission profile of SNAP 1 begins at approximately 150 m above ground level, without emissions near the surface. This can explain the finding that in Sect. 4.2.1 at the KIT 30 m height, the difference between tests S and P is larger for EDGAR v2024 than for CAMS-REG-ANT.

In summary, at KIT, where anthropogenic emissions are significant, SNAP 1 and SNAP 3 are the dominant emissions sectors explaining the improvement of the simulation results by including vertical emissions profiles. Considering anthropogenic vertical profile emissions in the model setup significantly improve the simulations performance, especially in the vicinity of large sources.

Figure 10 (a) CO₂ emission fluxes from different sectors in 2018 within the area (47.6–50.6° N, 6.94–9.94° E) surrounding the KIT site, and (b) the vertical emission profiles applied for the corresponding sectors, the height refers to the altitude above ground level.

5.2 Underestimation of model CO₂ simulations in early June

A significant underestimation of the simulated XCO₂ was found at the Orléans site in early June. According to the contributions of each tracer to the total simulated XCO₂ (see Fig. 11), a similar pattern is found in the biogenic component. Additionally, a comparison between simulated and observed near-surface CO₂ mole fractions at various heights from the co-located ICOS TRN site also reveals a clear underestimation in early June (Fig. A5). This site pair is hereafter referred to as Orléans/TRN. Figure 11c shows the vertical distribution of the biogenic CO₂ tracer at Orléans/TRN over time, revealing a sink spanning a large vertical extent in early June. At the other two TCCON sites (Paris and Karlsruhe), due to the lack of observational data in June, it is difficult to determine whether a similar feature exists. However, the underestimation of CO₂ mole fractions in early June is not confined to the Orléans/TRN site. A consistent underestimation during this period is also observed at the other ICOS sites (OPE, SAC, CBW, see Fig. A5), except for the KIT site which includes significant anthropogenic emissions, indicating that this is likely a regional feature rather than a localized effect at a single site. In addition, the STILT simulation results published on the ICOS Carbon Portal also exhibit an underestimation at these ICOS sites (except for the KIT site) in early June. The biospheric fluxes are from the VPRM model (Gerbig and Koch, 2023), with VPRM parameters optimized for the year 2007 using 46 sites across Europe, and the land cover classification is based on SYNMAP.

Figure 11 Time series (local time) of (a) observed and simulated total, background, (b) tracer-specific XCO₂, and (c) variation of simulated biogenic CO₂ mole fractions over time (local time) and height above ground at Orléans/TRN site. The simulated values here are all without AVK smoothing.

At most ICOS and TCCON sites included in this study, GPP and daytime NEE in June exhibit a stronger carbon sink compared to July and August. Starting from July, LSWI and EVI values, which are used as input to VPRM, across all vegetation types declined significantly, reflecting an overall decrease in vegetation activity and greenness. This trend was observed at the Orléans/TRN site (see Fig. A6), and similar patterns were found at other sites. This could be attributed to a prolonged and intense drought that occurred in Europe during the summer of 2018. This event directly affected temperature, soil moisture, precipitation, and ecosystem functioning (Buras et al., 2020). Simulations by WRF-GHG and VPRM, driven by meteorological fields and vegetation indices, both responded to this drought event.

Regarding the mismatch between the model and observations in early June, the simulated NEE shows a pronounced carbon sink at most sites. Unfortunately, the lack of co-located flux measurements at these atmospheric observation sites limits the direct evaluation of flux accuracy. Alternatively, we evaluated the four ICOS flux stations within the model's inner domain: two sites in Germany, DE_RuR (50.621914° N, 6.304126° E; grasslands) and DE_RuW (50.50493° N, 6.330963° E; forest), and two in France, FR_EM2 (49.87211° N, 3.02065° E; crops) and FR_LGt (47.322918° N, 2.284102° E; acid fen) (not shown here). The simulated RES shows a significant underestimation at all sites, except at FR_LGt site, which is consistent with Hu et al. (2021), who pointed out that the VPRM model tends to underestimate ecosystem respiration during the growing season, which could introduce biases in the simulation of NEE. At the two French sites, the simulated GPP exhibits more gross uptake significantly.

Additionally, we tested four different sets of VPRM parameter settings (not show here), all of which resulted in a similar underestimation in early June, preliminarily ruling out parameter settings as the primary source of the bias. It could rather be the consequences of the limitations either of the model itself.

In summary, the underestimation of simulated CO₂ mole fractions observed at multiple ICOS and TCCON sites in early June might be attributed to an underestimation of ecosystem respiration fluxes or an overestimation of GPP. However, owing to the lack of co-located flux measurements, the exact source of this bias remains unclear and requires further investigation.

6 Conclusions

In this study, we performed the first high-resolution simulation of the spatiotemporal distribution of CO₂ mole fractions over Western Europe during the summer (June to August) of 2018 using WRF-GHG. Given the significance of anthropogenic emissions and the diversity of emission inventories, we conducted several sensitivity tests and evaluated the model performance from both meteorological and chemical field perspectives by comparing with ground-based observations at multiple synoptic, ICOS and TCCON stations, on the purpose of optimizing the modelling set-up and – at a later stage – analyzing the CO₂ concentration variations over Belgium and surrounding countries systematically, ultimately offering a solid basis for related modelling and applied studies.

Overall, the WRF-GHG model reproduces the meteorological fields well, especially temperature, with high Pearson correlation coefficients ranging from 0.92 to 0.96 against observations from various synoptic sites. The model can capture the variations of wind speed at different heights, on both daily and monthly timescales. At inland stations, the model tends to overestimate wind speeds at most sites with dense vegetation cover, while it generally underestimates near airports or in areas dominated by cropland. These discrepancies highlight the importance of accurate and high-spatial-resolution land-cover maps for the reproduction of transport processes. Moreover, the diurnal variation of near-surface CO₂ mole fractions at different heights across the five ICOS observation sites was well captured by the model. During the summer 2018, variations in CO₂ mole fractions over Belgium and surrounding countries were mainly influenced by anthropogenic emissions and biogenic fluxes.

Sensitivity tests indicate that near large anthropogenic emission sources, the simulated near-surface CO₂ mole fractions are highly sensitive to the choice of anthropogenic emission inventory and the adoption of vertical emission profiles. At the KIT site in Germany, which is located near a very large oil refinery and power plant, differences between emission inventories can lead to discrepancies of up to 14.99 ± 31.98 ppm in simulated near surface CO₂ mole fractions. In addition, releasing anthropogenic emissions based on source-specific vertical profiles can significantly improve the accuracy of simulations, which is likely due to a more realistic representation of the real emissions. In contrast, at other observation sites where surrounding anthropogenic emissions are relatively low, the impact of vertical emission profiles on simulation results is much smaller.

Regarding XCO₂, the model seems to be less sensitive to the choice of emission inventories and anthropogenic emission heights, but certain effects are still evident. At the Paris urban site, all simulations overestimate XCO₂ by approximately 1.2–1.6 ppm. This bias is primarily attributed to the uncertainties introduced by anthropogenic emissions and using CAMS data as initial and boundary conditions. In addition, the differences caused by different anthropogenic emission inventories can reach up to 0.50 ± 0.34 ppm. Although it is less pronounced than for the near-surface mole fractions at the ICOS site, the simulation results confirm that considering vertical emission profiles leads to a modest improvement in model simulations at the Karlsruhe site.

Additionally, biogenic fluxes contribute significantly to CO₂ mole fractions during the growing season. The large negative bias observed in WRF-GHG simulations of CO₂ mole fractions in early June at most ICOS and TCCON sites may be attributed to the underestimation of RES or the overestimation of GPP by the VPRM model. However, due to the lack of flux observations, the exact cause remains uncertain.

This study systematically characterizes and evaluates the distribution of CO₂ concentrations across multiple observational sites, with a particular focus on regions surrounding Belgium, demonstrates the feasibility of using the WRF-GHG model to simulate CO₂ concentration variations in the core region of Western Europe, and provides a basis for applying this model to the study of other long-lived greenhouse gases like CH₄ and N₂O in this region. It further emphasizes the importance of optimizing initial and boundary conditions and refining the construction of source-specific vertical emission profiles to enhance simulation accuracy. Additionally, due to the relatively simplified parameterization of the current VPRM, its applicability across different climatic conditions remains limited. Nonlinear ecosystem responses under extreme temperature and moisture conditions may lead to biases in simulated respiration fluxes, which indirectly reflect the complexity of the biosphere-atmosphere system under extreme environmental stress and indicate that a modified VPRM model is necessary (Gourdji et al., 2022). At the same time, this study highlights the critical role of observations in evaluating regional simulations and supporting the understanding of carbon emissions and the distribution of CO₂ concentrations in the atmosphere. Overall, this work not only advances the understanding of the spatiotemporal variations of CO₂ in the core region of Western Europe, but also provides an important reference for developing more reliable regional-scale carbon-climate feedback modelling tools, thereby contributing to a deeper understanding of the interactions between extreme climate and the carbon cycle.

Appendix A

Table A1 Overview of the VPRM parameters. The temperature parameters are all in degrees Celsius. The abbreviations are as follows: EF – evergreen forest, DF – deciduous forest, MF – mixed forest, SHR – shrubland, WET – wetland, CRO – cropland, GRA – grassland.

Table A2 Sector mapping between different emission inventories (Granier et al., 2019).

Table A3 Evaluation metrics for temperature, wind speed, and wind direction between simulations and observations at the 21 synoptic stations. The eight stations marked with ^a are operated by Skeyes and Meteo Wing, where wind speed and wind direction observations are recorded only as integer values. The four stations marked with ^b represent coastal sites.

Figure A1 Map of evaluation metrics for temperature (a–d), wind speed (e–h), and wind direction (i–l) between simulated and observed values at 21 synoptic stations.

Figure A2 Diurnal cycles (local time) of simulations with different anthropogenic emissions and observations at 50 m at the TRN (a, b), 60 m at the SAC (c, d), 50 m at the OPE (e, f) and 207 m at the CBW (g, h) sites, where the values represent the MBE between each simulation and observations.

Figure A3 The diurnal cycles (local time) of planetary boundary layer height simulated by WRF-GHG at the KIT (a), the TRN (b), the SAC (c), the OPE (d) and the CBW (e) sites.

Figure A4 Diurnal cycles (local time) of simulations with different anthropogenic emissions at five ICOS sites at the model lowest layer (approximately 25 m above the ground).

Figure A5 Time series of the between the simulated near surface CO₂ mole fractions using TNO inventory and observed results at 30 m at the KIT (a), 50 m at the TRN (b), 60 m at the SAC (c), 50 m at the OPE (d) and 207 m at the CBW (e) sites, where the values represent the MBE ± SD between the simulations and observations. The red curve represents the 24 h moving average.

Figure A6 Time series of EVI (a) and LSWI (b) for different vegetation types from MODIS at Orléans/TRN site. According to the land cover data provided by Copernicus, there are no Shrubland (SHR) and Wetlands (WET) vegetation types near Orléans/TRN. Therefore, their corresponding EVI and LSWI values are 0 and are not shown here.