European Air Quality Index Calculation
Contents
European Air Quality Index Calculation¶
About¶
This notebook provides you a practical introduction to the calculation and interpretation of the Air Quality Index (AQI) in Europe. Air pollution is the single largest environmental health risk in Europe, causing cardiovascular and respiratory diseases, that, in the most serious cases, lead to premature deaths.
The European Environment Agency’s European Air Quality Index allows users to understand more about air quality where they live. Even though the air quality in Europe has improved over recent decades, the levels of air pollutants still exceed EU standards and the most stringent World Health Organization (WHO) guidelines.
How is the European Air Quality index defined?¶
The European Air Quality index is computed for five main pollutants regulated in the European legislation:
Ozone (O3)Nitrogen Dioxide (NO2)Sulphur Dioxide (SO2)Fine particulate matter with a diameter smaller than 2.5 micrometers (PM2.5)Fine particulate matter with a diameter smaller than 10 micrometers (PM10)
The index ranges from 1 (good) to 6 (extremely poor). For each pollutant, the index is calculated separately according to the concentrations; the higher the concentrations, the higher the index (see Table below for index levels). The overall hourly European Air Quality index is simply defined as the highest value of the five individual pollutants indexes computed for the same hour. For instance, if the indices relative to O3, NO2, SO2, PM2.5 and PM10 are 1,3,1,2,2 respectively, the overall index will be 3. The overall daily European Air Quality index is the highest value of the overall hourly European Air Quality index in the corresponding day.
The notebook has the following outline:
1 - Compute the European Air Quality Index for Europe for one day in July 2021
1.1 - Request data from the ADS programmatically with the CDS API
1.2 - Load and prepare the European air quality forecast data
1.3 - Classify daily maximum pollutant values into European Air Quality Index levels
1.4 - Visualize a map of Air Quality Index levels in Europe for 15 July 2021
2 - Calculate the daily European Air Quality Index for London in December 2020
2.1 - Request data from the ADS programmatically with the CDS API
2.2 - Load and prepare the European air quality forecast data
2.3 - Select time-series information for London, UK and convert to a pandas dataframe
2.4 - Classify daily maximum values of key pollutants into European Air Quality Index levels
2.5 - Visualize daily European Air Quality Index for London in December 2020 as heatmap
Data¶
This notebook introduces you to the CAMS European air quality forecasts and analysis dataset. The data has the following specifications:
Data:
CAMS European air quality forecasts
Temporal coverage:three-year rolling archive
Temporal resolution:hourly
Spatial coverage:Europe
Spatial resolution:0.1° x 0.1°
Format:NetCDF
Install CDSAPI via pip¶
!pip install cdsapi
Requirement already satisfied: cdsapi in /Volumes/gdrive/cop_consultancy/venv310/lib/python3.10/site-packages (0.5.1)
Requirement already satisfied: requests>=2.5.0 in /Volumes/gdrive/cop_consultancy/venv310/lib/python3.10/site-packages (from cdsapi) (2.26.0)
Requirement already satisfied: tqdm in /Volumes/gdrive/cop_consultancy/venv310/lib/python3.10/site-packages (from cdsapi) (4.62.3)
Requirement already satisfied: urllib3<1.27,>=1.21.1 in /Volumes/gdrive/cop_consultancy/venv310/lib/python3.10/site-packages (from requests>=2.5.0->cdsapi) (1.26.7)
Requirement already satisfied: idna<4,>=2.5 in /Volumes/gdrive/cop_consultancy/venv310/lib/python3.10/site-packages (from requests>=2.5.0->cdsapi) (3.3)
Requirement already satisfied: certifi>=2017.4.17 in /Volumes/gdrive/cop_consultancy/venv310/lib/python3.10/site-packages (from requests>=2.5.0->cdsapi) (2021.10.8)
Requirement already satisfied: charset-normalizer~=2.0.0 in /Volumes/gdrive/cop_consultancy/venv310/lib/python3.10/site-packages (from requests>=2.5.0->cdsapi) (2.0.7)
Load libraries¶
# CDS API
import cdsapi
import os
# Libraries for working with multi-dimensional arrays
import numpy as np
import xarray as xr
import pandas as pd
# Libraries for plotting and visualising data
import matplotlib.path as mpath
import matplotlib.pyplot as plt
from matplotlib import animation
from IPython.display import HTML
from matplotlib.colors import ListedColormap
from datetime import datetime
import cartopy.crs as ccrs
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
import cartopy.feature as cfeature
import seaborn as sns
1. Compute the European Air Quality Index for Europe for one day in July 2021¶
As a first example, we want to compute and visualize the European Air Quality Index for Europe for one day in July 2021.
1.1 Request data from the ADS programmatically with the CDS API¶
Let us request data from the Atmosphere Data Store programmatically with the help of the CDS API. Let us set manually the CDS API credentials. First, you have to define two variables: URL and KEY which build together your CDS API key. Below, you have to replace the ######### with your personal ADS key. Please find here your personal ADS key.
URL = 'https://ads.atmosphere.copernicus.eu/api/v2'
KEY = '######################################'
The next step is then to request the data with the help of the CDS API. Below, we request analysis data from the CAMS European air quality forecasts dataset. We request hourly data for 15 July 2021 for the five main pollutants for which the European Air Quality Index is calculated: nitrogen_dioxide, ozone, particulate_matter_10um, particulate_matter_2.5um and sulphur_dioxide.
Let us store the dataset under the name 20210715_eaqi_europe.nc.
import cdsapi
c = cdsapi.Client(url=URL, key=KEY)
c.retrieve(
'cams-europe-air-quality-forecasts',
{
'model': 'euradim',
'date': '2021-07-15/2021-07-15',
'format': 'netcdf',
'type': 'analysis',
'time': [
'00:00', '01:00', '02:00',
'03:00', '04:00', '05:00',
'06:00', '07:00', '08:00',
'09:00', '10:00', '11:00',
'12:00', '13:00', '14:00',
'15:00', '16:00', '17:00',
'18:00', '19:00', '20:00',
'21:00', '22:00', '23:00',
],
'leadtime_hour': '0',
'variable': [
'nitrogen_dioxide', 'ozone', 'particulate_matter_10um',
'particulate_matter_2.5um', 'sulphur_dioxide',
],
'level': '0',
},
'./20210715_eaqi_europe.nc')
1.2 Load and prepare the European air quality forecast data¶
The European air quality forecast data are available either in GRIB or in NetCDF format. Above, we requested the data in the NetCDF format. You can use the Python library xarray and the function open_dataset() to load a NetCDF file as xarray.Dataset. A xarray Dataset is a collection of one or more variables that share the same dimensions. Below, you see that the Dataset has four dimensions, latitude, level, longitude and time, and five variables, no2_conc, o3_conc, pm10_conc, pm2p5_conc and so2_conc.
eaqi = xr.open_dataset('./20210715_eaqi_europe.nc')
eaqi
<xarray.Dataset>
Dimensions: (latitude: 420, level: 1, longitude: 700, time: 24)
Coordinates:
* longitude (longitude) float32 335.0 335.1 335.2 ... 44.75 44.85 44.95
* latitude (latitude) float32 71.95 71.85 71.75 71.65 ... 30.25 30.15 30.05
* level (level) float32 0.0
* time (time) timedelta64[ns] 00:00:00 01:00:00 ... 22:00:00 23:00:00
Data variables:
no2_conc (time, level, latitude, longitude) float32 ...
o3_conc (time, level, latitude, longitude) float32 ...
pm10_conc (time, level, latitude, longitude) float32 ...
pm2p5_conc (time, level, latitude, longitude) float32 ...
so2_conc (time, level, latitude, longitude) float32 ...
Attributes:
title: NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
institution: Data produced by RIUUK
source: Data from EURAD-IM model
history: Model EURAD ANALYSIS
ANALYSIS: Europe, 20210715+[0H_23H]
summary: EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concent...
project: MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)- latitude: 420
- level: 1
- longitude: 700
- time: 24
- longitude(longitude)float32335.0 335.1 335.2 ... 44.85 44.95
- long_name :
- longitude
- units :
- degrees_east
array([335.05, 335.15, 335.25, ..., 44.75, 44.85, 44.95], dtype=float32)
- latitude(latitude)float3271.95 71.85 71.75 ... 30.15 30.05
- long_name :
- latitude
- units :
- degrees_north
array([71.95, 71.85, 71.75, ..., 30.25, 30.15, 30.05], dtype=float32)
- level(level)float320.0
- long_name :
- level
- units :
- m
array([0.], dtype=float32)
- time(time)timedelta64[ns]00:00:00 01:00:00 ... 23:00:00
- long_name :
- ANALYSIS time from 20210715
array([ 0, 3600000000000, 7200000000000, 10800000000000, 14400000000000, 18000000000000, 21600000000000, 25200000000000, 28800000000000, 32400000000000, 36000000000000, 39600000000000, 43200000000000, 46800000000000, 50400000000000, 54000000000000, 57600000000000, 61200000000000, 64800000000000, 68400000000000, 72000000000000, 75600000000000, 79200000000000, 82800000000000], dtype='timedelta64[ns]')
- no2_conc(time, level, latitude, longitude)float32...
- species :
- Nitrogen Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_nitrogen_dioxide_in_air
[7056000 values with dtype=float32]
- o3_conc(time, level, latitude, longitude)float32...
- species :
- Ozone
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_ozone_in_air
[7056000 values with dtype=float32]
- pm10_conc(time, level, latitude, longitude)float32...
- species :
- PM10 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm10_ambient_aerosol_in_air
[7056000 values with dtype=float32]
- pm2p5_conc(time, level, latitude, longitude)float32...
- species :
- PM2.5 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm2p5_ambient_aerosol_in_air
[7056000 values with dtype=float32]
- so2_conc(time, level, latitude, longitude)float32...
- species :
- Sulphur Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_sulfur_dioxide_in_air
[7056000 values with dtype=float32]
- title :
- NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
- institution :
- Data produced by RIUUK
- source :
- Data from EURAD-IM model
- history :
- Model EURAD ANALYSIS
- ANALYSIS :
- Europe, 20210715+[0H_23H]
- summary :
- EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concentration at the Surface from 20210715+[0H_23H] on Europe
- project :
- MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)
Let us have a closer look at the coordinates of the dataset. You see above that the dataset consists of 24 time steps and that the longitude values are on a [0, 360] grid.
However, if you inspect the time dimension more in detail, you see that the time is given in nanoseconds [ns]. Hence, as a next step, let us convert the time information into a human-readable time format and bring the longitude coordinates to a [-180, 180] grid.
eaqi.time
<xarray.DataArray 'time' (time: 24)>
array([ 0, 3600000000000, 7200000000000, 10800000000000,
14400000000000, 18000000000000, 21600000000000, 25200000000000,
28800000000000, 32400000000000, 36000000000000, 39600000000000,
43200000000000, 46800000000000, 50400000000000, 54000000000000,
57600000000000, 61200000000000, 64800000000000, 68400000000000,
72000000000000, 75600000000000, 79200000000000, 82800000000000],
dtype='timedelta64[ns]')
Coordinates:
* time (time) timedelta64[ns] 00:00:00 01:00:00 ... 22:00:00 23:00:00
Attributes:
long_name: ANALYSIS time from 20210715- time: 24
- 00:00:00 01:00:00 02:00:00 03:00:00 ... 21:00:00 22:00:00 23:00:00
array([ 0, 3600000000000, 7200000000000, 10800000000000, 14400000000000, 18000000000000, 21600000000000, 25200000000000, 28800000000000, 32400000000000, 36000000000000, 39600000000000, 43200000000000, 46800000000000, 50400000000000, 54000000000000, 57600000000000, 61200000000000, 64800000000000, 68400000000000, 72000000000000, 75600000000000, 79200000000000, 82800000000000], dtype='timedelta64[ns]') - time(time)timedelta64[ns]00:00:00 01:00:00 ... 23:00:00
- long_name :
- ANALYSIS time from 20210715
array([ 0, 3600000000000, 7200000000000, 10800000000000, 14400000000000, 18000000000000, 21600000000000, 25200000000000, 28800000000000, 32400000000000, 36000000000000, 39600000000000, 43200000000000, 46800000000000, 50400000000000, 54000000000000, 57600000000000, 61200000000000, 64800000000000, 68400000000000, 72000000000000, 75600000000000, 79200000000000, 82800000000000], dtype='timedelta64[ns]')
- long_name :
- ANALYSIS time from 20210715
First, from the long_name information of the time dimension (see above), we can retrieve the initial timestamp. With the function strptime() from Python’s datetime library, we can convert it into a datetime.datetime object.
timestamp = eaqi.time.long_name[19:27]
timestamp_init=datetime.strptime(timestamp,'%Y%m%d' )
timestamp_init
datetime.datetime(2021, 7, 15, 0, 0)
In a next step, we then build a DateTimeIndex object with the help of Panda’s date_range() function, making use of the length of the time dimension. The result is a DateTimeIndex object, which can be used to newly assign the time coordinate information.
time_coords = pd.date_range(timestamp_init, periods=len(eaqi.time), freq='1h').strftime("%Y-%m-%d %H:%M:%S").astype('datetime64[ns]')
time_coords
DatetimeIndex(['2021-07-15 00:00:00', '2021-07-15 01:00:00',
'2021-07-15 02:00:00', '2021-07-15 03:00:00',
'2021-07-15 04:00:00', '2021-07-15 05:00:00',
'2021-07-15 06:00:00', '2021-07-15 07:00:00',
'2021-07-15 08:00:00', '2021-07-15 09:00:00',
'2021-07-15 10:00:00', '2021-07-15 11:00:00',
'2021-07-15 12:00:00', '2021-07-15 13:00:00',
'2021-07-15 14:00:00', '2021-07-15 15:00:00',
'2021-07-15 16:00:00', '2021-07-15 17:00:00',
'2021-07-15 18:00:00', '2021-07-15 19:00:00',
'2021-07-15 20:00:00', '2021-07-15 21:00:00',
'2021-07-15 22:00:00', '2021-07-15 23:00:00'],
dtype='datetime64[ns]', freq=None)
And the last step is to assign the converted time information to the loaded xarray.Dataset eaqi, with the function assign_coords().
eaqi_assigned = eaqi.assign_coords(time=time_coords)
eaqi_assigned
<xarray.Dataset>
Dimensions: (latitude: 420, level: 1, longitude: 700, time: 24)
Coordinates:
* longitude (longitude) float32 335.0 335.1 335.2 ... 44.75 44.85 44.95
* latitude (latitude) float32 71.95 71.85 71.75 71.65 ... 30.25 30.15 30.05
* level (level) float32 0.0
* time (time) datetime64[ns] 2021-07-15 ... 2021-07-15T23:00:00
Data variables:
no2_conc (time, level, latitude, longitude) float32 ...
o3_conc (time, level, latitude, longitude) float32 ...
pm10_conc (time, level, latitude, longitude) float32 ...
pm2p5_conc (time, level, latitude, longitude) float32 ...
so2_conc (time, level, latitude, longitude) float32 ...
Attributes:
title: NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
institution: Data produced by RIUUK
source: Data from EURAD-IM model
history: Model EURAD ANALYSIS
ANALYSIS: Europe, 20210715+[0H_23H]
summary: EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concent...
project: MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)- latitude: 420
- level: 1
- longitude: 700
- time: 24
- longitude(longitude)float32335.0 335.1 335.2 ... 44.85 44.95
- long_name :
- longitude
- units :
- degrees_east
array([335.05, 335.15, 335.25, ..., 44.75, 44.85, 44.95], dtype=float32)
- latitude(latitude)float3271.95 71.85 71.75 ... 30.15 30.05
- long_name :
- latitude
- units :
- degrees_north
array([71.95, 71.85, 71.75, ..., 30.25, 30.15, 30.05], dtype=float32)
- level(level)float320.0
- long_name :
- level
- units :
- m
array([0.], dtype=float32)
- time(time)datetime64[ns]2021-07-15 ... 2021-07-15T23:00:00
array(['2021-07-15T00:00:00.000000000', '2021-07-15T01:00:00.000000000', '2021-07-15T02:00:00.000000000', '2021-07-15T03:00:00.000000000', '2021-07-15T04:00:00.000000000', '2021-07-15T05:00:00.000000000', '2021-07-15T06:00:00.000000000', '2021-07-15T07:00:00.000000000', '2021-07-15T08:00:00.000000000', '2021-07-15T09:00:00.000000000', '2021-07-15T10:00:00.000000000', '2021-07-15T11:00:00.000000000', '2021-07-15T12:00:00.000000000', '2021-07-15T13:00:00.000000000', '2021-07-15T14:00:00.000000000', '2021-07-15T15:00:00.000000000', '2021-07-15T16:00:00.000000000', '2021-07-15T17:00:00.000000000', '2021-07-15T18:00:00.000000000', '2021-07-15T19:00:00.000000000', '2021-07-15T20:00:00.000000000', '2021-07-15T21:00:00.000000000', '2021-07-15T22:00:00.000000000', '2021-07-15T23:00:00.000000000'], dtype='datetime64[ns]')
- no2_conc(time, level, latitude, longitude)float32...
- species :
- Nitrogen Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_nitrogen_dioxide_in_air
[7056000 values with dtype=float32]
- o3_conc(time, level, latitude, longitude)float32...
- species :
- Ozone
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_ozone_in_air
[7056000 values with dtype=float32]
- pm10_conc(time, level, latitude, longitude)float32...
- species :
- PM10 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm10_ambient_aerosol_in_air
[7056000 values with dtype=float32]
- pm2p5_conc(time, level, latitude, longitude)float32...
- species :
- PM2.5 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm2p5_ambient_aerosol_in_air
[7056000 values with dtype=float32]
- so2_conc(time, level, latitude, longitude)float32...
- species :
- Sulphur Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_sulfur_dioxide_in_air
[7056000 values with dtype=float32]
- title :
- NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
- institution :
- Data produced by RIUUK
- source :
- Data from EURAD-IM model
- history :
- Model EURAD ANALYSIS
- ANALYSIS :
- Europe, 20210715+[0H_23H]
- summary :
- EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concentration at the Surface from 20210715+[0H_23H] on Europe
- project :
- MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)
Now, let us also re-assign the longitude values and shift the grid from [0,360] to [-180,180]. At the end, you might want to sort the longitude values in an ascending order. You can re-assign coordinate values with the xarray function assign_coords().
eaqi_assigned = eaqi_assigned.assign_coords(longitude=(((eaqi_assigned.longitude + 180) % 360) - 180)).sortby('longitude')
eaqi_assigned
<xarray.Dataset>
Dimensions: (latitude: 420, level: 1, longitude: 700, time: 24)
Coordinates:
* longitude (longitude) float32 -24.95 -24.85 -24.75 ... 44.75 44.85 44.95
* latitude (latitude) float32 71.95 71.85 71.75 71.65 ... 30.25 30.15 30.05
* level (level) float32 0.0
* time (time) datetime64[ns] 2021-07-15 ... 2021-07-15T23:00:00
Data variables:
no2_conc (time, level, latitude, longitude) float32 ...
o3_conc (time, level, latitude, longitude) float32 ...
pm10_conc (time, level, latitude, longitude) float32 ...
pm2p5_conc (time, level, latitude, longitude) float32 ...
so2_conc (time, level, latitude, longitude) float32 ...
Attributes:
title: NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
institution: Data produced by RIUUK
source: Data from EURAD-IM model
history: Model EURAD ANALYSIS
ANALYSIS: Europe, 20210715+[0H_23H]
summary: EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concent...
project: MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)- latitude: 420
- level: 1
- longitude: 700
- time: 24
- longitude(longitude)float32-24.95 -24.85 ... 44.85 44.95
array([-24.950012, -24.849976, -24.75 , ..., 44.75 , 44.850006, 44.949997], dtype=float32) - latitude(latitude)float3271.95 71.85 71.75 ... 30.15 30.05
- long_name :
- latitude
- units :
- degrees_north
array([71.95, 71.85, 71.75, ..., 30.25, 30.15, 30.05], dtype=float32)
- level(level)float320.0
- long_name :
- level
- units :
- m
array([0.], dtype=float32)
- time(time)datetime64[ns]2021-07-15 ... 2021-07-15T23:00:00
array(['2021-07-15T00:00:00.000000000', '2021-07-15T01:00:00.000000000', '2021-07-15T02:00:00.000000000', '2021-07-15T03:00:00.000000000', '2021-07-15T04:00:00.000000000', '2021-07-15T05:00:00.000000000', '2021-07-15T06:00:00.000000000', '2021-07-15T07:00:00.000000000', '2021-07-15T08:00:00.000000000', '2021-07-15T09:00:00.000000000', '2021-07-15T10:00:00.000000000', '2021-07-15T11:00:00.000000000', '2021-07-15T12:00:00.000000000', '2021-07-15T13:00:00.000000000', '2021-07-15T14:00:00.000000000', '2021-07-15T15:00:00.000000000', '2021-07-15T16:00:00.000000000', '2021-07-15T17:00:00.000000000', '2021-07-15T18:00:00.000000000', '2021-07-15T19:00:00.000000000', '2021-07-15T20:00:00.000000000', '2021-07-15T21:00:00.000000000', '2021-07-15T22:00:00.000000000', '2021-07-15T23:00:00.000000000'], dtype='datetime64[ns]')
- no2_conc(time, level, latitude, longitude)float32...
- species :
- Nitrogen Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_nitrogen_dioxide_in_air
[7056000 values with dtype=float32]
- o3_conc(time, level, latitude, longitude)float32...
- species :
- Ozone
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_ozone_in_air
[7056000 values with dtype=float32]
- pm10_conc(time, level, latitude, longitude)float32...
- species :
- PM10 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm10_ambient_aerosol_in_air
[7056000 values with dtype=float32]
- pm2p5_conc(time, level, latitude, longitude)float32...
- species :
- PM2.5 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm2p5_ambient_aerosol_in_air
[7056000 values with dtype=float32]
- so2_conc(time, level, latitude, longitude)float32...
- species :
- Sulphur Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_sulfur_dioxide_in_air
[7056000 values with dtype=float32]
- title :
- NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
- institution :
- Data produced by RIUUK
- source :
- Data from EURAD-IM model
- history :
- Model EURAD ANALYSIS
- ANALYSIS :
- Europe, 20210715+[0H_23H]
- summary :
- EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concentration at the Surface from 20210715+[0H_23H] on Europe
- project :
- MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)
Let us now also resample the hourly concentrations of the five pollutants and extract for each pollutant the maximum value for each day. The result is that the Dataset remains with only one entry of the time dimension.
eaqi_daily = eaqi_assigned.resample(time='1D').max()
eaqi_daily
<xarray.Dataset>
Dimensions: (latitude: 420, level: 1, longitude: 700, time: 1)
Coordinates:
* time (time) datetime64[ns] 2021-07-15
* longitude (longitude) float32 -24.95 -24.85 -24.75 ... 44.75 44.85 44.95
* latitude (latitude) float32 71.95 71.85 71.75 71.65 ... 30.25 30.15 30.05
* level (level) float32 0.0
Data variables:
no2_conc (time, level, latitude, longitude) float32 0.02787 ... 3.931
o3_conc (time, level, latitude, longitude) float32 61.97 61.97 ... 125.3
pm10_conc (time, level, latitude, longitude) float32 0.7773 ... 52.61
pm2p5_conc (time, level, latitude, longitude) float32 0.5984 ... 30.76
so2_conc (time, level, latitude, longitude) float32 0.001265 ... 22.91- latitude: 420
- level: 1
- longitude: 700
- time: 1
- time(time)datetime64[ns]2021-07-15
array(['2021-07-15T00:00:00.000000000'], dtype='datetime64[ns]')
- longitude(longitude)float32-24.95 -24.85 ... 44.85 44.95
array([-24.950012, -24.849976, -24.75 , ..., 44.75 , 44.850006, 44.949997], dtype=float32) - latitude(latitude)float3271.95 71.85 71.75 ... 30.15 30.05
- long_name :
- latitude
- units :
- degrees_north
array([71.95, 71.85, 71.75, ..., 30.25, 30.15, 30.05], dtype=float32)
- level(level)float320.0
- long_name :
- level
- units :
- m
array([0.], dtype=float32)
- no2_conc(time, level, latitude, longitude)float320.02787 0.02787 ... 4.235 3.931
array([[[[0.02787461, 0.02787461, 0.02771118, ..., 0.15140247, 0.14035352, 0.13047698], [0.02787461, 0.02787461, 0.02787461, ..., 0.17789151, 0.17789151, 0.17789151], [0.02787461, 0.02787461, 0.02787461, ..., 0.17789151, 0.17789151, 0.17789151], ..., [0.04484161, 0.04464976, 0.04462844, ..., 4.309593 , 4.309593 , 4.309593 ], [0.04396053, 0.04391233, 0.04391944, ..., 4.309593 , 4.309593 , 4.266726 ], [0.04391944, 0.04391944, 0.04391944, ..., 3.0567997 , 4.2351356 , 3.9313145 ]]]], dtype=float32) - o3_conc(time, level, latitude, longitude)float3261.97 61.97 61.81 ... 124.9 125.3
array([[[[ 61.965622, 61.965393, 61.80987 , ..., 50.012226, 49.533863, 49.057087], [ 61.965622, 61.965622, 61.965622, ..., 51.123573, 51.123573, 51.123573], [ 61.965622, 61.965622, 61.965622, ..., 51.123573, 51.123573, 51.123573], ..., [ 58.006176, 57.504444, 56.047306, ..., 125.32788 , 125.32788 , 125.32788 ], [ 56.17458 , 55.964684, 55.96062 , ..., 125.32788 , 125.32788 , 125.32767 ], [ 55.96062 , 55.96062 , 55.96062 , ..., 118.43542 , 124.91819 , 125.32606 ]]]], dtype=float32) - pm10_conc(time, level, latitude, longitude)float320.7773 0.7773 ... 52.86 52.61
array([[[[ 0.77728105, 0.77728105, 0.79217404, ..., 2.511054 , 2.5265799 , 2.546816 ], [ 0.77728105, 0.77728105, 0.77728105, ..., 2.4791365 , 2.4791365 , 2.4791365 ], [ 0.77728105, 0.77728105, 0.77728105, ..., 2.4791365 , 2.4791365 , 2.4791365 ], ..., [ 9.749049 , 9.740864 , 9.613364 , ..., 53.43701 , 53.43701 , 53.43701 ], [ 9.998307 , 10.912123 , 10.929857 , ..., 53.43701 , 53.43701 , 53.342907 ], [10.929857 , 10.929857 , 10.929857 , ..., 43.721504 , 52.85957 , 52.60667 ]]]], dtype=float32) - pm2p5_conc(time, level, latitude, longitude)float320.5984 0.5984 ... 30.56 30.76
array([[[[ 0.5983948, 0.5983948, 0.6123783, ..., 1.8387041, 1.8370837, 1.8356059], [ 0.5983948, 0.5983948, 0.5983948, ..., 1.8662741, 1.8662741, 1.8662741], [ 0.5983948, 0.5983948, 0.5983948, ..., 1.8662741, 1.8662741, 1.8662741], ..., [ 4.6096487, 4.6148496, 4.574519 , ..., 30.862682 , 30.862682 , 30.862682 ], [ 4.63375 , 4.6709533, 4.671948 , ..., 30.862682 , 30.862682 , 30.849081 ], [ 4.671948 , 4.671948 , 4.671948 , ..., 25.94779 , 30.555159 , 30.761532 ]]]], dtype=float32) - so2_conc(time, level, latitude, longitude)float320.001265 0.001265 ... 23.63 22.91
array([[[[1.2649833e-03, 1.2649833e-03, 8.8108622e-04, ..., 3.6806125e-02, 3.7403006e-02, 3.8014073e-02], [1.2649833e-03, 1.2649833e-03, 1.2649833e-03, ..., 4.0259570e-02, 4.0259570e-02, 4.0259570e-02], [1.2649833e-03, 1.2649833e-03, 1.2649833e-03, ..., 4.0259570e-02, 4.0259570e-02, 4.0259570e-02], ..., [1.3102422e-02, 1.2903470e-02, 1.2022397e-02, ..., 2.4014425e+01, 2.4014425e+01, 2.4014425e+01], [1.2306614e-02, 1.1993974e-02, 1.1993974e-02, ..., 2.4014425e+01, 2.4014425e+01, 2.3863991e+01], [1.1993974e-02, 1.1993974e-02, 1.1993974e-02, ..., 1.7468096e+01, 2.3625347e+01, 2.2913668e+01]]]], dtype=float32)
Let us now also drop the dimensions time and level, which only have one record each. The xarray function squeeze() drops all dimensions with only one entry. The result is a dataset with two remaining dimensions: latitude and longitude.
eaqi_daily = eaqi_daily.squeeze(drop=True)
eaqi_daily
<xarray.Dataset>
Dimensions: (latitude: 420, longitude: 700)
Coordinates:
* longitude (longitude) float32 -24.95 -24.85 -24.75 ... 44.75 44.85 44.95
* latitude (latitude) float32 71.95 71.85 71.75 71.65 ... 30.25 30.15 30.05
Data variables:
no2_conc (latitude, longitude) float32 0.02787 0.02787 ... 4.235 3.931
o3_conc (latitude, longitude) float32 61.97 61.97 61.81 ... 124.9 125.3
pm10_conc (latitude, longitude) float32 0.7773 0.7773 ... 52.86 52.61
pm2p5_conc (latitude, longitude) float32 0.5984 0.5984 ... 30.56 30.76
so2_conc (latitude, longitude) float32 0.001265 0.001265 ... 23.63 22.91- latitude: 420
- longitude: 700
- longitude(longitude)float32-24.95 -24.85 ... 44.85 44.95
array([-24.950012, -24.849976, -24.75 , ..., 44.75 , 44.850006, 44.949997], dtype=float32) - latitude(latitude)float3271.95 71.85 71.75 ... 30.15 30.05
- long_name :
- latitude
- units :
- degrees_north
array([71.95, 71.85, 71.75, ..., 30.25, 30.15, 30.05], dtype=float32)
- no2_conc(latitude, longitude)float320.02787 0.02787 ... 4.235 3.931
array([[0.02787461, 0.02787461, 0.02771118, ..., 0.15140247, 0.14035352, 0.13047698], [0.02787461, 0.02787461, 0.02787461, ..., 0.17789151, 0.17789151, 0.17789151], [0.02787461, 0.02787461, 0.02787461, ..., 0.17789151, 0.17789151, 0.17789151], ..., [0.04484161, 0.04464976, 0.04462844, ..., 4.309593 , 4.309593 , 4.309593 ], [0.04396053, 0.04391233, 0.04391944, ..., 4.309593 , 4.309593 , 4.266726 ], [0.04391944, 0.04391944, 0.04391944, ..., 3.0567997 , 4.2351356 , 3.9313145 ]], dtype=float32) - o3_conc(latitude, longitude)float3261.97 61.97 61.81 ... 124.9 125.3
array([[ 61.965622, 61.965393, 61.80987 , ..., 50.012226, 49.533863, 49.057087], [ 61.965622, 61.965622, 61.965622, ..., 51.123573, 51.123573, 51.123573], [ 61.965622, 61.965622, 61.965622, ..., 51.123573, 51.123573, 51.123573], ..., [ 58.006176, 57.504444, 56.047306, ..., 125.32788 , 125.32788 , 125.32788 ], [ 56.17458 , 55.964684, 55.96062 , ..., 125.32788 , 125.32788 , 125.32767 ], [ 55.96062 , 55.96062 , 55.96062 , ..., 118.43542 , 124.91819 , 125.32606 ]], dtype=float32) - pm10_conc(latitude, longitude)float320.7773 0.7773 ... 52.86 52.61
array([[ 0.77728105, 0.77728105, 0.79217404, ..., 2.511054 , 2.5265799 , 2.546816 ], [ 0.77728105, 0.77728105, 0.77728105, ..., 2.4791365 , 2.4791365 , 2.4791365 ], [ 0.77728105, 0.77728105, 0.77728105, ..., 2.4791365 , 2.4791365 , 2.4791365 ], ..., [ 9.749049 , 9.740864 , 9.613364 , ..., 53.43701 , 53.43701 , 53.43701 ], [ 9.998307 , 10.912123 , 10.929857 , ..., 53.43701 , 53.43701 , 53.342907 ], [10.929857 , 10.929857 , 10.929857 , ..., 43.721504 , 52.85957 , 52.60667 ]], dtype=float32) - pm2p5_conc(latitude, longitude)float320.5984 0.5984 ... 30.56 30.76
array([[ 0.5983948, 0.5983948, 0.6123783, ..., 1.8387041, 1.8370837, 1.8356059], [ 0.5983948, 0.5983948, 0.5983948, ..., 1.8662741, 1.8662741, 1.8662741], [ 0.5983948, 0.5983948, 0.5983948, ..., 1.8662741, 1.8662741, 1.8662741], ..., [ 4.6096487, 4.6148496, 4.574519 , ..., 30.862682 , 30.862682 , 30.862682 ], [ 4.63375 , 4.6709533, 4.671948 , ..., 30.862682 , 30.862682 , 30.849081 ], [ 4.671948 , 4.671948 , 4.671948 , ..., 25.94779 , 30.555159 , 30.761532 ]], dtype=float32) - so2_conc(latitude, longitude)float320.001265 0.001265 ... 23.63 22.91
array([[1.2649833e-03, 1.2649833e-03, 8.8108622e-04, ..., 3.6806125e-02, 3.7403006e-02, 3.8014073e-02], [1.2649833e-03, 1.2649833e-03, 1.2649833e-03, ..., 4.0259570e-02, 4.0259570e-02, 4.0259570e-02], [1.2649833e-03, 1.2649833e-03, 1.2649833e-03, ..., 4.0259570e-02, 4.0259570e-02, 4.0259570e-02], ..., [1.3102422e-02, 1.2903470e-02, 1.2022397e-02, ..., 2.4014425e+01, 2.4014425e+01, 2.4014425e+01], [1.2306614e-02, 1.1993974e-02, 1.1993974e-02, ..., 2.4014425e+01, 2.4014425e+01, 2.3863991e+01], [1.1993974e-02, 1.1993974e-02, 1.1993974e-02, ..., 1.7468096e+01, 2.3625347e+01, 2.2913668e+01]], dtype=float32)
1.3 Classify daily maximum concentration values into European Air Quality Index levels¶
Now, we can classify the daily maximum concentration values into European Air Quality index levels. For this, we first define for each key pollutant the thresholds and save these as a list of threshold limits. In a next step, we combine the threshold lists into a list of lists, which we call bin_list.
no2_limits = [0, 40, 90, 120, 230, 340, 1000]
ozone_limits = [0, 50, 100, 130, 240, 380, 800]
pm10_limits = [0, 10, 20, 25, 50, 75, 800]
pm25_limits = [0, 20, 40, 50, 100, 150, 1200]
so2_limits = [0, 100, 200, 350, 500, 750, 1250]
bin_list = [no2_limits, ozone_limits, pm10_limits, pm25_limits, so2_limits]
bin_list
[[0, 40, 90, 120, 230, 340, 1000],
[0, 50, 100, 130, 240, 380, 800],
[0, 10, 20, 25, 50, 75, 800],
[0, 20, 40, 50, 100, 150, 1200],
[0, 100, 200, 350, 500, 750, 1250]]
The next step is then to store all variables available in the xarray.Dataset eaqi_daily in a list with the name variable_list.
variable_list = list(eaqi_daily.variables)[2:]
variable_list
['no2_conc', 'o3_conc', 'pm10_conc', 'pm2p5_conc', 'so2_conc']
Now, we can loop over all data variables in the dataset and apply a customized function with the apply_ufunc() functionality of xarray:
function(): apply the numpy functiondigitize(), which returns the indices of the bins (thresholds) to which each value in the input array belongs toinput array: provide input data array, whose values shall be classifiedlist of thresholds (bin): provide list of thresholds, groups of bins
At the end, we create a list of all classified data arrays, which has the name classified_arrays.
classified_arrays = []
for i in range(0,len(eaqi_daily)):
variable = variable_list[i]
temp = xr.apply_ufunc(np.digitize,
eaqi_daily[variable],
bin_list[i])
classified_arrays.append(temp)
The function merge() allows us to merge the list of xarray.DataArray items into a xarray.Dataset. The resulting dataset has again five variables, but each grid point of the variable’s array is now classified according to the European Air Quality Index levels.
combined = xr.merge(classified_arrays)
combined
<xarray.Dataset>
Dimensions: (latitude: 420, longitude: 700)
Coordinates:
* longitude (longitude) float32 -24.95 -24.85 -24.75 ... 44.75 44.85 44.95
* latitude (latitude) float32 71.95 71.85 71.75 71.65 ... 30.25 30.15 30.05
Data variables:
no2_conc (latitude, longitude) int64 1 1 1 1 1 1 1 1 ... 1 1 1 1 1 1 1 1
o3_conc (latitude, longitude) int64 2 2 2 2 2 2 2 2 ... 3 3 3 3 3 3 3 3
pm10_conc (latitude, longitude) int64 1 1 1 1 1 1 1 1 ... 4 4 4 4 4 4 5 5
pm2p5_conc (latitude, longitude) int64 1 1 1 1 1 1 1 1 ... 2 2 2 2 2 2 2 2
so2_conc (latitude, longitude) int64 1 1 1 1 1 1 1 1 ... 1 1 1 1 1 1 1 1- latitude: 420
- longitude: 700
- longitude(longitude)float32-24.95 -24.85 ... 44.85 44.95
array([-24.950012, -24.849976, -24.75 , ..., 44.75 , 44.850006, 44.949997], dtype=float32) - latitude(latitude)float3271.95 71.85 71.75 ... 30.15 30.05
- long_name :
- latitude
- units :
- degrees_north
array([71.95, 71.85, 71.75, ..., 30.25, 30.15, 30.05], dtype=float32)
- no2_conc(latitude, longitude)int641 1 1 1 1 1 1 1 ... 1 1 1 1 1 1 1 1
array([[1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1], ..., [1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1]]) - o3_conc(latitude, longitude)int642 2 2 2 2 2 2 2 ... 3 3 3 3 3 3 3 3
array([[2, 2, 2, ..., 2, 1, 1], [2, 2, 2, ..., 2, 2, 2], [2, 2, 2, ..., 2, 2, 2], ..., [2, 2, 2, ..., 3, 3, 3], [2, 2, 2, ..., 3, 3, 3], [2, 2, 2, ..., 3, 3, 3]]) - pm10_conc(latitude, longitude)int641 1 1 1 1 1 1 1 ... 4 4 4 4 4 4 5 5
array([[1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1], ..., [1, 1, 1, ..., 5, 5, 5], [1, 2, 2, ..., 5, 5, 5], [2, 2, 2, ..., 4, 5, 5]]) - pm2p5_conc(latitude, longitude)int641 1 1 1 1 1 1 1 ... 2 2 2 2 2 2 2 2
array([[1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1], ..., [1, 1, 1, ..., 2, 2, 2], [1, 1, 1, ..., 2, 2, 2], [1, 1, 1, ..., 2, 2, 2]]) - so2_conc(latitude, longitude)int641 1 1 1 1 1 1 1 ... 1 1 1 1 1 1 1 1
array([[1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1], ..., [1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1]])
The last step before we can visualize the Air Quality Index in Europe for 15 July 2021 is to create the overall daily European Air Quality index based on the maximum index level among the five key pollutants of the respective day. For this, we can tranform the xarray.Dataset above into a xarray.DataArray with the function to_array(). The result is a three-dimensional data array and the five data variables will be added as an additional dimension with the name variable. From this data array, we can retrieve for each grid point the maximum value across all variables with the function max().
The result is a two-dimensional data array, which holds for each grid point the overall Air Quality Index level.
final_index = combined.to_array().max('variable')
final_index
<xarray.DataArray (latitude: 420, longitude: 700)>
array([[2, 2, 2, ..., 2, 1, 1],
[2, 2, 2, ..., 2, 2, 2],
[2, 2, 2, ..., 2, 2, 2],
...,
[2, 2, 2, ..., 5, 5, 5],
[2, 2, 2, ..., 5, 5, 5],
[2, 2, 2, ..., 4, 5, 5]])
Coordinates:
* longitude (longitude) float32 -24.95 -24.85 -24.75 ... 44.75 44.85 44.95
* latitude (latitude) float32 71.95 71.85 71.75 71.65 ... 30.25 30.15 30.05- latitude: 420
- longitude: 700
- 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 ... 4 4 4 4 4 4 4 4 4 4 4 4 4 4 5 5
array([[2, 2, 2, ..., 2, 1, 1], [2, 2, 2, ..., 2, 2, 2], [2, 2, 2, ..., 2, 2, 2], ..., [2, 2, 2, ..., 5, 5, 5], [2, 2, 2, ..., 5, 5, 5], [2, 2, 2, ..., 4, 5, 5]]) - longitude(longitude)float32-24.95 -24.85 ... 44.85 44.95
array([-24.950012, -24.849976, -24.75 , ..., 44.75 , 44.850006, 44.949997], dtype=float32) - latitude(latitude)float3271.95 71.85 71.75 ... 30.15 30.05
- long_name :
- latitude
- units :
- degrees_north
array([71.95, 71.85, 71.75, ..., 30.25, 30.15, 30.05], dtype=float32)
1.4 Visualize a map of Air Quality Index levels in Europe for 15 July 2021¶
The last step is now to visualize the overall daily Air Quality Index levels in Europe for 15 July 2021. The visualisation function consists of five main parts:
1. Initiate the matplotlib figure: Initiate a figure and axes objects and add a colorbar
2. Plotting function: Plot the xarray.DataArray with the pcolormesh() function from matplotlib and provide the customized colormap as color range
3. Set title of the map: Set a title to the resulting map
4. Customize colorbar entries: Customize the labels of the colorbar and indicate the index levels from ‘very good’ to ‘extremely poor’
5. Add additional features to the map: Add additional features such as coastlines and a geographic grid
Note: before plotting, we additionally define a customized colormap based on the official EAQI index colors. We also add additionally the labels to the six EAQI index levels, which range from ‘Very good’ to ‘Extremely poor’.
cmap = ListedColormap(['#5AAA5F', '#A7D25C', '#ECD347', '#F5BE41', '#F09235', '#D93322'])
labels = ['Very good', 'Good', 'Medium', 'Poor', 'Very Poor', 'Extremely Poor']
# Initiate a matplotlib figure
fig = plt.figure(figsize=(20,15))
ax = plt.subplot(1,1,1, projection=ccrs.PlateCarree())
# Plotting function
img = plt.pcolormesh(final_index.longitude,
final_index.latitude,
final_index,
cmap = cmap,
label = labels)
# Set the title of the plot
ax.set_title('European Air Quality Index levels - 15 July 2021', fontsize=16, pad=20.0)
# Customize colorbar
cbar = plt.colorbar(img,fraction=0.028, pad=0.05)
cbar.set_label(None)
cbar.ax.tick_params(labelsize=12)
cbar.set_ticks([1.5,2.5,3.5,4.5,5.5,6.5])
cbar.set_ticklabels(labels)
# Add additional features such as coastlines and a geographic grid
ax.coastlines(color='black')
gl = ax.gridlines(draw_labels=True, linewidth=1, color='gray', alpha=0.5, linestyle='--')
gl.xlabel_style={'size':12}
gl.ylabel_style={'size':12}
Text(0.5, 1.0, 'European Air Quality Index levels - 15 July 2021')
2. Calculate the daily European Air Quality Index for London in December 2020¶
Next, we want to calculate the daily European Air Quality Index for London, UK in December 2020.
2.1 Request data from the ADS programmatically with the CDS API¶
Let us request data from the Atmosphere Data Store programmatically with the help of the CDS API. We again set manually the CDS API credentials. First, you have to define two variables: URL and KEY which build together your CDS API key. Below, you have to replace the ######### with your personal ADS key. Please find here your personal ADS key.
URL = 'https://ads.atmosphere.copernicus.eu/api/v2'
KEY = '######################################'
The next step is then to request the data with the help of the CDS API. Below, we request analysis data from the CAMS European air quality forecasts dataset. We request hourly data for December 2020 for the five main pollutants for which the European Air Quality Index is calculated: nitrogen_dioxide, ozone, particulate_matter_10um, particulate_matter_2.5um and sulphur_dioxide. In the data request below, we additionally set the area key, as we want to tailor the data request to only a small geographical area around London, United Kingdom.
Let us store the dataset under the name 202012_eaqi_london.nc.
c = cdsapi.Client(url=URL, key=KEY)
c.retrieve(
'cams-europe-air-quality-forecasts',
{
'model': 'euradim',
'date': '2020-12-01/2020-12-31',
'format': 'netcdf',
'type': 'analysis',
'time': [
'00:00', '01:00', '02:00',
'03:00', '04:00', '05:00',
'06:00', '07:00', '08:00',
'09:00', '10:00', '11:00',
'12:00', '13:00', '14:00',
'15:00', '16:00', '17:00',
'18:00', '19:00', '20:00',
'21:00', '22:00', '23:00',
],
'leadtime_hour': '0',
'variable': [
'nitrogen_dioxide', 'ozone', 'particulate_matter_10um',
'particulate_matter_2.5um', 'sulphur_dioxide',
],
'area': [51.6, -0.3, 51.3, 0],
'level': '0',
},
'./202012_eaqi_london.nc')
2.2 Load and prepare the European air quality forecast data¶
The European air quality forecasts data are available either in GRIB or in NetCDF format. Above, we requested the data in the NetCDF format. You can use the Python library xarray and the function open_dataset() to load a NetCDF file as xarray.Dataset. A xarray Dataset is a collection of one or more variables that share the same dimensions. Below, you see that the Dataset has four dimensions, latitude, level, longitude and time, and five variables, no2_conc, o3_conc, pm10_conc, pm2p5_conc and so2_conc.
eaqi_london = xr.open_dataset('./202012_eaqi_london.nc')
eaqi_london
<xarray.Dataset>
Dimensions: (latitude: 3, level: 1, longitude: 3, time: 744)
Coordinates:
* longitude (longitude) float32 359.8 359.9 360.0
* latitude (latitude) float32 51.55 51.45 51.35
* level (level) float32 0.0
* time (time) timedelta64[ns] 00:00:00 01:00:00 ... 30 days 23:00:00
Data variables:
no2_conc (time, level, latitude, longitude) float32 15.51 10.58 ... 46.72
o3_conc (time, level, latitude, longitude) float32 37.22 36.29 ... 1.924
pm10_conc (time, level, latitude, longitude) float32 9.667 10.13 ... 19.78
pm2p5_conc (time, level, latitude, longitude) float32 8.624 9.112 ... 18.0
so2_conc (time, level, latitude, longitude) float32 2.99 3.44 ... 11.16
Attributes:
title: NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
institution: Data produced by RIUUK
source: Data from EURAD-IM model
history: Model EURAD ANALYSIS
ANALYSIS: Europe, 20201201-20201231+[0H_23H]
summary: EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concent...
project: MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)- latitude: 3
- level: 1
- longitude: 3
- time: 744
- longitude(longitude)float32359.8 359.9 360.0
- long_name :
- longitude
- units :
- degrees_east
array([359.75, 359.85, 359.95], dtype=float32)
- latitude(latitude)float3251.55 51.45 51.35
- long_name :
- latitude
- units :
- degrees_north
array([51.55, 51.45, 51.35], dtype=float32)
- level(level)float320.0
- long_name :
- level
- units :
- m
array([0.], dtype=float32)
- time(time)timedelta64[ns]00:00:00 ... 30 days 23:00:00
- long_name :
- ANALYSIS time from 20201201
array([ 0, 3600000000000, 7200000000000, ..., 2667600000000000, 2671200000000000, 2674800000000000], dtype='timedelta64[ns]')
- no2_conc(time, level, latitude, longitude)float32...
- species :
- Nitrogen Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_nitrogen_dioxide_in_air
array([[[[15.510345, ..., 18.174795], ..., [22.29408 , ..., 20.352665]]], ..., [[[36.479404, ..., 35.76146 ], ..., [55.333344, ..., 46.718143]]]], dtype=float32) - o3_conc(time, level, latitude, longitude)float32...
- species :
- Ozone
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_ozone_in_air
array([[[[37.22002 , ..., 37.39592 ], ..., [32.99458 , ..., 32.783318]]], ..., [[[17.871344, ..., 16.65876 ], ..., [ 2.548177, ..., 1.923809]]]], dtype=float32) - pm10_conc(time, level, latitude, longitude)float32...
- species :
- PM10 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm10_ambient_aerosol_in_air
array([[[[ 9.666528, ..., 10.224049], ..., [ 9.492814, ..., 10.637783]]], ..., [[[20.009241, ..., 21.280088], ..., [18.24215 , ..., 19.779423]]]], dtype=float32) - pm2p5_conc(time, level, latitude, longitude)float32...
- species :
- PM2.5 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm2p5_ambient_aerosol_in_air
array([[[[ 8.624495, ..., 9.175108], ..., [ 8.614178, ..., 9.701087]]], ..., [[[16.95187 , ..., 18.536068], ..., [16.331594, ..., 17.998186]]]], dtype=float32) - so2_conc(time, level, latitude, longitude)float32...
- species :
- Sulphur Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_sulfur_dioxide_in_air
array([[[[ 2.990305, ..., 3.666514], ..., [ 3.343928, ..., 4.098411]]], ..., [[[ 7.277549, ..., 8.871666], ..., [ 8.752181, ..., 11.158818]]]], dtype=float32)
- title :
- NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
- institution :
- Data produced by RIUUK
- source :
- Data from EURAD-IM model
- history :
- Model EURAD ANALYSIS
- ANALYSIS :
- Europe, 20201201-20201231+[0H_23H]
- summary :
- EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concentration at the Surface from 20201201-20201231+[0H_23H] on Europe
- project :
- MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)
Let us inspect the coordinates of the dataset more in detail. You see above that the dataset consists of 744 time steps and that the longitude values are on a [0, 360] grid.
However, if you inspect the time dimension more in detail, you see that the time is given in nanoseconds [ns]. As a next step, let us convert the time information into a human-readable time format and bring the longitude coordinates to a [-180, 180] grid.
eaqi_london.time
<xarray.DataArray 'time' (time: 744)>
array([ 0, 3600000000000, 7200000000000, ...,
2667600000000000, 2671200000000000, 2674800000000000],
dtype='timedelta64[ns]')
Coordinates:
* time (time) timedelta64[ns] 00:00:00 01:00:00 ... 30 days 23:00:00
Attributes:
long_name: ANALYSIS time from 20201201- time: 744
- 00:00:00 01:00:00 02:00:00 ... 30 days 22:00:00 30 days 23:00:00
array([ 0, 3600000000000, 7200000000000, ..., 2667600000000000, 2671200000000000, 2674800000000000], dtype='timedelta64[ns]') - time(time)timedelta64[ns]00:00:00 ... 30 days 23:00:00
- long_name :
- ANALYSIS time from 20201201
array([ 0, 3600000000000, 7200000000000, ..., 2667600000000000, 2671200000000000, 2674800000000000], dtype='timedelta64[ns]')
- long_name :
- ANALYSIS time from 20201201
First, from the long_name information of the time dimension (see above), we can retrieve the initial timestamp. With the function strptime() from Python’s datetime library, we can convert it into a datetime.datetime object.
timestamp = eaqi_london.time.long_name[19:27]
timestamp_init=datetime.strptime(timestamp,'%Y%m%d' )
timestamp_init
datetime.datetime(2020, 12, 1, 0, 0)
In a next step, we then build a DateTimeIndex object with the help of Panda’s date_range() function, making use of the length of the time dimension. The result is a DateTimeIndex object, which can be used to newly assign the time coordinate information.
time_coords = pd.date_range(timestamp_init, periods=len(eaqi_london.time), freq='1h').strftime("%Y-%m-%d %H:%M:%S").astype('datetime64[ns]')
time_coords
DatetimeIndex(['2020-12-01 00:00:00', '2020-12-01 01:00:00',
'2020-12-01 02:00:00', '2020-12-01 03:00:00',
'2020-12-01 04:00:00', '2020-12-01 05:00:00',
'2020-12-01 06:00:00', '2020-12-01 07:00:00',
'2020-12-01 08:00:00', '2020-12-01 09:00:00',
...
'2020-12-31 14:00:00', '2020-12-31 15:00:00',
'2020-12-31 16:00:00', '2020-12-31 17:00:00',
'2020-12-31 18:00:00', '2020-12-31 19:00:00',
'2020-12-31 20:00:00', '2020-12-31 21:00:00',
'2020-12-31 22:00:00', '2020-12-31 23:00:00'],
dtype='datetime64[ns]', length=744, freq=None)
And the last step is to assign the converted time information to the loaded xarray.Dataset eaqi_ldn_ts, with the function assign_coords().
eaqi_london = eaqi_london.assign_coords(time=time_coords)
eaqi_london
<xarray.Dataset>
Dimensions: (latitude: 3, level: 1, longitude: 3, time: 744)
Coordinates:
* longitude (longitude) float32 359.8 359.9 360.0
* latitude (latitude) float32 51.55 51.45 51.35
* level (level) float32 0.0
* time (time) datetime64[ns] 2020-12-01 ... 2020-12-31T23:00:00
Data variables:
no2_conc (time, level, latitude, longitude) float32 15.51 10.58 ... 46.72
o3_conc (time, level, latitude, longitude) float32 37.22 36.29 ... 1.924
pm10_conc (time, level, latitude, longitude) float32 9.667 10.13 ... 19.78
pm2p5_conc (time, level, latitude, longitude) float32 8.624 9.112 ... 18.0
so2_conc (time, level, latitude, longitude) float32 2.99 3.44 ... 11.16
Attributes:
title: NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
institution: Data produced by RIUUK
source: Data from EURAD-IM model
history: Model EURAD ANALYSIS
ANALYSIS: Europe, 20201201-20201231+[0H_23H]
summary: EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concent...
project: MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)- latitude: 3
- level: 1
- longitude: 3
- time: 744
- longitude(longitude)float32359.8 359.9 360.0
- long_name :
- longitude
- units :
- degrees_east
array([359.75, 359.85, 359.95], dtype=float32)
- latitude(latitude)float3251.55 51.45 51.35
- long_name :
- latitude
- units :
- degrees_north
array([51.55, 51.45, 51.35], dtype=float32)
- level(level)float320.0
- long_name :
- level
- units :
- m
array([0.], dtype=float32)
- time(time)datetime64[ns]2020-12-01 ... 2020-12-31T23:00:00
array(['2020-12-01T00:00:00.000000000', '2020-12-01T01:00:00.000000000', '2020-12-01T02:00:00.000000000', ..., '2020-12-31T21:00:00.000000000', '2020-12-31T22:00:00.000000000', '2020-12-31T23:00:00.000000000'], dtype='datetime64[ns]')
- no2_conc(time, level, latitude, longitude)float3215.51 10.58 18.17 ... 49.02 46.72
- species :
- Nitrogen Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_nitrogen_dioxide_in_air
array([[[[15.510345, ..., 18.174795], ..., [22.29408 , ..., 20.352665]]], ..., [[[36.479404, ..., 35.76146 ], ..., [55.333344, ..., 46.718143]]]], dtype=float32) - o3_conc(time, level, latitude, longitude)float3237.22 36.29 37.4 ... 2.409 1.924
- species :
- Ozone
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_ozone_in_air
array([[[[37.22002 , ..., 37.39592 ], ..., [32.99458 , ..., 32.783318]]], ..., [[[17.871344, ..., 16.65876 ], ..., [ 2.548177, ..., 1.923809]]]], dtype=float32) - pm10_conc(time, level, latitude, longitude)float329.667 10.13 10.22 ... 19.35 19.78
- species :
- PM10 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm10_ambient_aerosol_in_air
array([[[[ 9.666528, ..., 10.224049], ..., [ 9.492814, ..., 10.637783]]], ..., [[[20.009241, ..., 21.280088], ..., [18.24215 , ..., 19.779423]]]], dtype=float32) - pm2p5_conc(time, level, latitude, longitude)float328.624 9.112 9.175 ... 17.46 18.0
- species :
- PM2.5 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm2p5_ambient_aerosol_in_air
array([[[[ 8.624495, ..., 9.175108], ..., [ 8.614178, ..., 9.701087]]], ..., [[[16.95187 , ..., 18.536068], ..., [16.331594, ..., 17.998186]]]], dtype=float32) - so2_conc(time, level, latitude, longitude)float322.99 3.44 3.667 ... 10.33 11.16
- species :
- Sulphur Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_sulfur_dioxide_in_air
array([[[[ 2.990305, ..., 3.666514], ..., [ 3.343928, ..., 4.098411]]], ..., [[[ 7.277549, ..., 8.871666], ..., [ 8.752181, ..., 11.158818]]]], dtype=float32)
- title :
- NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
- institution :
- Data produced by RIUUK
- source :
- Data from EURAD-IM model
- history :
- Model EURAD ANALYSIS
- ANALYSIS :
- Europe, 20201201-20201231+[0H_23H]
- summary :
- EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concentration at the Surface from 20201201-20201231+[0H_23H] on Europe
- project :
- MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)
Now, let us also re-assign the longitude values and shift the grid from [0,360] to [-180,180]. At the end, you might want to sort the longitude values in an ascending order. You can re-assign coordinate values with the xarray function assign_coords().
eaqi_london_assigned = eaqi_london.assign_coords(longitude=(((eaqi_london.longitude + 180) % 360) - 180)).sortby('longitude')
eaqi_london_assigned
<xarray.Dataset>
Dimensions: (latitude: 3, level: 1, longitude: 3, time: 744)
Coordinates:
* longitude (longitude) float32 -0.25 -0.15 -0.04999
* latitude (latitude) float32 51.55 51.45 51.35
* level (level) float32 0.0
* time (time) datetime64[ns] 2020-12-01 ... 2020-12-31T23:00:00
Data variables:
no2_conc (time, level, latitude, longitude) float32 15.51 10.58 ... 46.72
o3_conc (time, level, latitude, longitude) float32 37.22 36.29 ... 1.924
pm10_conc (time, level, latitude, longitude) float32 9.667 10.13 ... 19.78
pm2p5_conc (time, level, latitude, longitude) float32 8.624 9.112 ... 18.0
so2_conc (time, level, latitude, longitude) float32 2.99 3.44 ... 11.16
Attributes:
title: NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
institution: Data produced by RIUUK
source: Data from EURAD-IM model
history: Model EURAD ANALYSIS
ANALYSIS: Europe, 20201201-20201231+[0H_23H]
summary: EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concent...
project: MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)- latitude: 3
- level: 1
- longitude: 3
- time: 744
- longitude(longitude)float32-0.25 -0.15 -0.04999
array([-0.25 , -0.150024, -0.049988], dtype=float32)
- latitude(latitude)float3251.55 51.45 51.35
- long_name :
- latitude
- units :
- degrees_north
array([51.55, 51.45, 51.35], dtype=float32)
- level(level)float320.0
- long_name :
- level
- units :
- m
array([0.], dtype=float32)
- time(time)datetime64[ns]2020-12-01 ... 2020-12-31T23:00:00
array(['2020-12-01T00:00:00.000000000', '2020-12-01T01:00:00.000000000', '2020-12-01T02:00:00.000000000', ..., '2020-12-31T21:00:00.000000000', '2020-12-31T22:00:00.000000000', '2020-12-31T23:00:00.000000000'], dtype='datetime64[ns]')
- no2_conc(time, level, latitude, longitude)float3215.51 10.58 18.17 ... 49.02 46.72
- species :
- Nitrogen Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_nitrogen_dioxide_in_air
array([[[[15.510345, ..., 18.174795], ..., [22.29408 , ..., 20.352665]]], ..., [[[36.479404, ..., 35.76146 ], ..., [55.333344, ..., 46.718143]]]], dtype=float32) - o3_conc(time, level, latitude, longitude)float3237.22 36.29 37.4 ... 2.409 1.924
- species :
- Ozone
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_ozone_in_air
array([[[[37.22002 , ..., 37.39592 ], ..., [32.99458 , ..., 32.783318]]], ..., [[[17.871344, ..., 16.65876 ], ..., [ 2.548177, ..., 1.923809]]]], dtype=float32) - pm10_conc(time, level, latitude, longitude)float329.667 10.13 10.22 ... 19.35 19.78
- species :
- PM10 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm10_ambient_aerosol_in_air
array([[[[ 9.666528, ..., 10.224049], ..., [ 9.492814, ..., 10.637783]]], ..., [[[20.009241, ..., 21.280088], ..., [18.24215 , ..., 19.779423]]]], dtype=float32) - pm2p5_conc(time, level, latitude, longitude)float328.624 9.112 9.175 ... 17.46 18.0
- species :
- PM2.5 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm2p5_ambient_aerosol_in_air
array([[[[ 8.624495, ..., 9.175108], ..., [ 8.614178, ..., 9.701087]]], ..., [[[16.95187 , ..., 18.536068], ..., [16.331594, ..., 17.998186]]]], dtype=float32) - so2_conc(time, level, latitude, longitude)float322.99 3.44 3.667 ... 10.33 11.16
- species :
- Sulphur Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_sulfur_dioxide_in_air
array([[[[ 2.990305, ..., 3.666514], ..., [ 3.343928, ..., 4.098411]]], ..., [[[ 7.277549, ..., 8.871666], ..., [ 8.752181, ..., 11.158818]]]], dtype=float32)
- title :
- NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
- institution :
- Data produced by RIUUK
- source :
- Data from EURAD-IM model
- history :
- Model EURAD ANALYSIS
- ANALYSIS :
- Europe, 20201201-20201231+[0H_23H]
- summary :
- EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concentration at the Surface from 20201201-20201231+[0H_23H] on Europe
- project :
- MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)
Let us now also drop the dimension level, which only has one record, as we requested the data of the five key pollutants on the surface. The xarray function squeeze() drops all dimensions with only one entry. The result is a dataset with three dimensions: latitude, longitude and time.
eaqi_london_assigned = eaqi_london_assigned.squeeze(drop=True)
eaqi_london_assigned
<xarray.Dataset>
Dimensions: (latitude: 3, longitude: 3, time: 744)
Coordinates:
* longitude (longitude) float32 -0.25 -0.15 -0.04999
* latitude (latitude) float32 51.55 51.45 51.35
* time (time) datetime64[ns] 2020-12-01 ... 2020-12-31T23:00:00
Data variables:
no2_conc (time, latitude, longitude) float32 15.51 10.58 ... 49.02 46.72
o3_conc (time, latitude, longitude) float32 37.22 36.29 ... 2.409 1.924
pm10_conc (time, latitude, longitude) float32 9.667 10.13 ... 19.35 19.78
pm2p5_conc (time, latitude, longitude) float32 8.624 9.112 ... 17.46 18.0
so2_conc (time, latitude, longitude) float32 2.99 3.44 ... 10.33 11.16
Attributes:
title: NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
institution: Data produced by RIUUK
source: Data from EURAD-IM model
history: Model EURAD ANALYSIS
ANALYSIS: Europe, 20201201-20201231+[0H_23H]
summary: EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concent...
project: MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)- latitude: 3
- longitude: 3
- time: 744
- longitude(longitude)float32-0.25 -0.15 -0.04999
array([-0.25 , -0.150024, -0.049988], dtype=float32)
- latitude(latitude)float3251.55 51.45 51.35
- long_name :
- latitude
- units :
- degrees_north
array([51.55, 51.45, 51.35], dtype=float32)
- time(time)datetime64[ns]2020-12-01 ... 2020-12-31T23:00:00
array(['2020-12-01T00:00:00.000000000', '2020-12-01T01:00:00.000000000', '2020-12-01T02:00:00.000000000', ..., '2020-12-31T21:00:00.000000000', '2020-12-31T22:00:00.000000000', '2020-12-31T23:00:00.000000000'], dtype='datetime64[ns]')
- no2_conc(time, latitude, longitude)float3215.51 10.58 18.17 ... 49.02 46.72
- species :
- Nitrogen Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_nitrogen_dioxide_in_air
array([[[15.510345, 10.579534, 18.174795], [20.994917, 18.346903, 22.862679], [22.29408 , 21.535889, 20.352665]], [[13.2516 , 8.708938, 13.546945], [18.658539, 14.884748, 17.36678 ], [20.586605, 18.495684, 17.110401]], ..., [[33.60502 , 18.792007, 38.690742], [50.68169 , 36.13545 , 40.703144], [57.09093 , 50.720074, 47.777363]], [[36.479404, 25.34729 , 35.76146 ], [49.230564, 35.99857 , 39.039207], [55.333344, 49.02131 , 46.718143]]], dtype=float32) - o3_conc(time, latitude, longitude)float3237.22 36.29 37.4 ... 2.409 1.924
- species :
- Ozone
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_ozone_in_air
array([[[37.22002 , 36.288414, 37.39592 ], [32.87107 , 30.231646, 29.552027], [32.99458 , 32.758266, 32.783318]], [[45.106533, 46.573887, 47.846783], [39.874947, 40.614536, 41.431847], [37.11564 , 38.48072 , 39.64118 ]], ..., [[21.514835, 31.547274, 14.364104], [ 6.591847, 11.323124, 7.030365], [ 1.501007, 2.097181, 2.34354 ]], [[17.871344, 27.412085, 16.65876 ], [ 7.608634, 12.21619 , 8.43275 ], [ 2.548177, 2.408655, 1.923809]]], dtype=float32) - pm10_conc(time, latitude, longitude)float329.667 10.13 10.22 ... 19.35 19.78
- species :
- PM10 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm10_ambient_aerosol_in_air
array([[[ 9.666528, 10.134321, 10.224049], [ 9.93193 , 11.177739, 11.461103], [ 9.492814, 10.298428, 10.637783]], [[ 8.71655 , 8.762962, 8.572779], [ 9.382599, 10.342698, 10.400551], [ 9.995328, 10.667473, 10.90721 ]], ..., [[19.123522, 20.140621, 20.78949 ], [18.499098, 19.941158, 20.598553], [16.874968, 17.37354 , 17.563284]], [[20.009241, 20.410782, 21.280088], [19.186546, 21.101486, 22.031332], [18.24215 , 19.353893, 19.779423]]], dtype=float32) - pm2p5_conc(time, latitude, longitude)float328.624 9.112 9.175 ... 17.46 18.0
- species :
- PM2.5 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm2p5_ambient_aerosol_in_air
array([[[ 8.624495, 9.111983, 9.175108], [ 8.910218, 10.065255, 10.328034], [ 8.614178, 9.357618, 9.701087]], [[ 7.747547, 7.855393, 7.650344], [ 8.371894, 9.294804, 9.34192 ], [ 9.043769, 9.675157, 9.929397]], ..., [[16.371088, 17.581823, 18.21415 ], [16.062485, 17.593449, 18.26474 ], [15.012928, 15.503884, 15.781394]], [[16.95187 , 17.51425 , 18.536068], [16.613737, 18.597458, 19.570078], [16.331594, 17.455702, 17.998186]]], dtype=float32) - so2_conc(time, latitude, longitude)float322.99 3.44 3.667 ... 10.33 11.16
- species :
- Sulphur Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_sulfur_dioxide_in_air
array([[[ 2.990305, 3.440334, 3.666514], [ 3.569994, 4.487958, 4.805997], [ 3.343928, 3.8807 , 4.098411]], [[ 2.548347, 2.805848, 3.030266], [ 3.150262, 3.59762 , 3.752746], [ 3.091941, 3.374168, 3.476543]], ..., [[ 7.829386, 9.040606, 9.668385], [ 9.034239, 11.045587, 11.954286], [ 8.767416, 10.257509, 10.950203]], [[ 7.277549, 8.221832, 8.871666], [ 8.627126, 10.639724, 11.633574], [ 8.752181, 10.331519, 11.158818]]], dtype=float32)
- title :
- NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
- institution :
- Data produced by RIUUK
- source :
- Data from EURAD-IM model
- history :
- Model EURAD ANALYSIS
- ANALYSIS :
- Europe, 20201201-20201231+[0H_23H]
- summary :
- EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concentration at the Surface from 20201201-20201231+[0H_23H] on Europe
- project :
- MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)
2.3 Select time-series information for London, UK and convert to a pandas dataframe¶
The xarray.Dataset above is now in a format that allows us to select the values of the pollutants for one single grid point. Let us first define variables for the latitude and longitude coordinates for London, United Kingdom. These two variables can then be used to select the values for one single grid point. You can use xarray’s function sel() to make a coordinate-based selection. Let us set the keyword argument method='nearest', which selects the information that belongs to the closest grid point to the provided latitude and longitude inputs.
The resulting dataset remains with time as the only dimension.
london_lat = 51.5
london_lon = -0.1
eaqi_ldn_ts = eaqi_london_assigned.sel(latitude=london_lat, longitude=london_lon, method='nearest')
eaqi_ldn_ts
<xarray.Dataset>
Dimensions: (time: 744)
Coordinates:
longitude float32 -0.04999
latitude float32 51.55
* time (time) datetime64[ns] 2020-12-01 ... 2020-12-31T23:00:00
Data variables:
no2_conc (time) float32 18.17 13.55 10.95 11.06 ... 40.85 38.69 35.76
o3_conc (time) float32 37.4 47.85 51.77 51.07 ... 14.49 14.36 16.66
pm10_conc (time) float32 10.22 8.573 7.056 6.757 ... 19.02 20.79 21.28
pm2p5_conc (time) float32 9.175 7.65 6.211 5.926 ... 16.62 18.21 18.54
so2_conc (time) float32 3.667 3.03 2.836 2.631 ... 9.585 9.668 8.872
Attributes:
title: NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
institution: Data produced by RIUUK
source: Data from EURAD-IM model
history: Model EURAD ANALYSIS
ANALYSIS: Europe, 20201201-20201231+[0H_23H]
summary: EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concent...
project: MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)- time: 744
- longitude()float32-0.04999
array(-0.04998779, dtype=float32)
- latitude()float3251.55
- long_name :
- latitude
- units :
- degrees_north
array(51.55, dtype=float32)
- time(time)datetime64[ns]2020-12-01 ... 2020-12-31T23:00:00
array(['2020-12-01T00:00:00.000000000', '2020-12-01T01:00:00.000000000', '2020-12-01T02:00:00.000000000', ..., '2020-12-31T21:00:00.000000000', '2020-12-31T22:00:00.000000000', '2020-12-31T23:00:00.000000000'], dtype='datetime64[ns]')
- no2_conc(time)float3218.17 13.55 10.95 ... 38.69 35.76
- species :
- Nitrogen Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_nitrogen_dioxide_in_air
array([18.174795, 13.546945, 10.945883, ..., 40.848576, 38.690742, 35.76146 ], dtype=float32) - o3_conc(time)float3237.4 47.85 51.77 ... 14.36 16.66
- species :
- Ozone
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_ozone_in_air
array([37.39592 , 47.846783, 51.769108, ..., 14.485551, 14.364104, 16.65876 ], dtype=float32) - pm10_conc(time)float3210.22 8.573 7.056 ... 20.79 21.28
- species :
- PM10 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm10_ambient_aerosol_in_air
array([10.224049, 8.572779, 7.05563 , ..., 19.01986 , 20.78949 , 21.280088], dtype=float32) - pm2p5_conc(time)float329.175 7.65 6.211 ... 18.21 18.54
- species :
- PM2.5 Aerosol
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_pm2p5_ambient_aerosol_in_air
array([ 9.175108, 7.650344, 6.210651, ..., 16.623259, 18.21415 , 18.536068], dtype=float32) - so2_conc(time)float323.667 3.03 2.836 ... 9.668 8.872
- species :
- Sulphur Dioxide
- units :
- µg/m3
- value :
- hourly values
- standard_name :
- mass_concentration_of_sulfur_dioxide_in_air
array([3.666514, 3.030266, 2.836089, ..., 9.584825, 9.668385, 8.871666], dtype=float32)
- title :
- NO2/O3/PM10/PM25/SO2 Air Pollutant ANALYSIS at the Surface
- institution :
- Data produced by RIUUK
- source :
- Data from EURAD-IM model
- history :
- Model EURAD ANALYSIS
- ANALYSIS :
- Europe, 20201201-20201231+[0H_23H]
- summary :
- EURAD model hourly ANALYSIS of NO2/O3/PM10/PM25/SO2 concentration at the Surface from 20201201-20201231+[0H_23H] on Europe
- project :
- MACC-RAQ (http://macc-raq.gmes-atmosphere.eu)
Our aim is to compute the daily European Air Quality Index for London for every day in December 2020. Hence, as a next step we want to resample the hourly data to a daily resolution. As resample function, we use max(), which selects for every key pollutant the maximum value per day. You can use the function resample() to aggregated the hourly information to daily maximum values. The result is a decrease of the time dimension entries from 744 to 31, one entry for every day in December 2020.
eaqi_lnd_ts_daily = eaqi_ldn_ts.resample(time='1D').max()
eaqi_lnd_ts_daily
<xarray.Dataset>
Dimensions: (time: 31)
Coordinates:
* time (time) datetime64[ns] 2020-12-01 2020-12-02 ... 2020-12-31
longitude float32 -0.04999
latitude float32 51.55
Data variables:
no2_conc (time) float64 36.3 65.69 70.01 66.55 ... 65.24 84.48 52.93
o3_conc (time) float64 51.77 24.67 49.51 77.38 ... 32.98 23.57 20.49
pm10_conc (time) float64 22.57 50.78 28.98 17.73 ... 41.05 36.28 33.25
pm2p5_conc (time) float64 20.82 47.19 26.39 16.21 ... 37.79 32.32 29.52
so2_conc (time) float64 8.053 41.21 48.89 18.53 ... 24.58 39.45 15.31- time: 31
- time(time)datetime64[ns]2020-12-01 ... 2020-12-31
array(['2020-12-01T00:00:00.000000000', '2020-12-02T00:00:00.000000000', '2020-12-03T00:00:00.000000000', '2020-12-04T00:00:00.000000000', '2020-12-05T00:00:00.000000000', '2020-12-06T00:00:00.000000000', '2020-12-07T00:00:00.000000000', '2020-12-08T00:00:00.000000000', '2020-12-09T00:00:00.000000000', '2020-12-10T00:00:00.000000000', '2020-12-11T00:00:00.000000000', '2020-12-12T00:00:00.000000000', '2020-12-13T00:00:00.000000000', '2020-12-14T00:00:00.000000000', '2020-12-15T00:00:00.000000000', '2020-12-16T00:00:00.000000000', '2020-12-17T00:00:00.000000000', '2020-12-18T00:00:00.000000000', '2020-12-19T00:00:00.000000000', '2020-12-20T00:00:00.000000000', '2020-12-21T00:00:00.000000000', '2020-12-22T00:00:00.000000000', '2020-12-23T00:00:00.000000000', '2020-12-24T00:00:00.000000000', '2020-12-25T00:00:00.000000000', '2020-12-26T00:00:00.000000000', '2020-12-27T00:00:00.000000000', '2020-12-28T00:00:00.000000000', '2020-12-29T00:00:00.000000000', '2020-12-30T00:00:00.000000000', '2020-12-31T00:00:00.000000000'], dtype='datetime64[ns]') - longitude()float32-0.04999
array(-0.04998779, dtype=float32)
- latitude()float3251.55
- long_name :
- latitude
- units :
- degrees_north
array(51.55, dtype=float32)
- no2_conc(time)float6436.3 65.69 70.01 ... 84.48 52.93
array([ 36.30299377, 65.68988037, 70.01146698, 66.55126953, 47.04573822, 36.37680435, 57.13370514, 57.63101196, 64.7775116 , 43.43681717, 113.46591949, 57.88596344, 56.97869873, 41.86167526, 54.94801712, 63.9841423 , 48.93317413, 21.91687584, 19.73365974, 25.81886673, 26.26212692, 69.27537537, 51.012043 , 12.58273888, 51.90653992, 46.18640137, 71.50334167, 67.38331604, 65.24152374, 84.48013306, 52.9269371 ]) - o3_conc(time)float6451.77 24.67 49.51 ... 23.57 20.49
array([51.76910782, 24.66534042, 49.51294708, 77.37693787, 55.20477676, 31.80633736, 19.87851334, 7.21585321, 32.25561905, 44.31040955, 44.1403656 , 45.81469727, 55.79682922, 75.31567383, 60.60375214, 56.93139648, 49.17032623, 64.97071838, 81.63644409, 68.82138062, 55.71593094, 74.09149933, 73.12327576, 80.72418976, 64.10220337, 69.87267303, 83.72303009, 41.2680397 , 32.98464584, 23.56501579, 20.48646545]) - pm10_conc(time)float6422.57 50.78 28.98 ... 36.28 33.25
array([22.57060051, 50.78494644, 28.97666359, 17.72842407, 19.83301544, 27.09029961, 42.65241241, 43.27971649, 32.22806168, 19.20415497, 46.63470078, 30.87269211, 19.65514374, 10.62569714, 12.21476555, 16.3072834 , 12.91017246, nan, 12.67010784, 8.51046944, 9.83042145, 38.75607681, 16.13843155, 5.15066624, 21.04909515, 17.39035606, 14.51557732, 35.5093956 , 41.04904938, 36.27640152, 33.24543762]) - pm2p5_conc(time)float6420.82 47.19 26.39 ... 32.32 29.52
array([20.82224083, 47.19441223, 26.39168358, 16.21479034, 18.2527771 , 24.14494705, 40.13755417, 40.59964752, 29.39426804, 17.02283287, 40.18713379, 25.86927795, 17.98646927, 8.36485767, 9.67144585, 14.95213985, 10.86158085, nan, 11.53662968, 6.05193615, 7.51580477, 33.04491425, 14.70893383, 4.40539551, 19.382658 , 16.17809486, 12.41596222, 29.79351234, 37.79423141, 32.32485199, 29.51881409]) - so2_conc(time)float648.053 41.21 48.89 ... 39.45 15.31
array([ 8.05340481, 41.21045685, 48.89346695, 18.53280258, 9.04989815, 8.82772636, 21.6830349 , 25.11410522, 31.41758537, 13.8028326 , 49.34526062, 14.57840443, 17.1786499 , 13.11939049, 12.9295187 , 25.2906456 , 10.98476315, nan, 4.22994614, 6.06269169, 6.2261734 , 25.56100845, 11.81103992, 3.35580802, 15.9163847 , 7.16892147, 13.61837578, 18.92266273, 24.57670975, 39.44808197, 15.31316185])
Above, the xarray.Dataset holds the daily time-series information for London for the five main pollutants. The Python library pandas is very effective in handling time-series data and offers an interface to xarray. Let us convert the xarray.Dataset above to a pandas.DataFrame with the function to_dataframe(). The result is a dataframe with 31 row entries and seven columns, including longitude and latitude information as well as the maximum daily values for each main pollutant.
eaqi_ts_daily_df = eaqi_lnd_ts_daily.to_dataframe()
eaqi_ts_daily_df
| longitude | latitude | no2_conc | o3_conc | pm10_conc | pm2p5_conc | so2_conc | |
|---|---|---|---|---|---|---|---|
| time | |||||||
| 2020-12-01 | -0.049988 | 51.549999 | 36.302994 | 51.769108 | 22.570601 | 20.822241 | 8.053405 |
| 2020-12-02 | -0.049988 | 51.549999 | 65.689880 | 24.665340 | 50.784946 | 47.194412 | 41.210457 |
| 2020-12-03 | -0.049988 | 51.549999 | 70.011467 | 49.512947 | 28.976664 | 26.391684 | 48.893467 |
| 2020-12-04 | -0.049988 | 51.549999 | 66.551270 | 77.376938 | 17.728424 | 16.214790 | 18.532803 |
| 2020-12-05 | -0.049988 | 51.549999 | 47.045738 | 55.204777 | 19.833015 | 18.252777 | 9.049898 |
| 2020-12-06 | -0.049988 | 51.549999 | 36.376804 | 31.806337 | 27.090300 | 24.144947 | 8.827726 |
| 2020-12-07 | -0.049988 | 51.549999 | 57.133705 | 19.878513 | 42.652412 | 40.137554 | 21.683035 |
| 2020-12-08 | -0.049988 | 51.549999 | 57.631012 | 7.215853 | 43.279716 | 40.599648 | 25.114105 |
| 2020-12-09 | -0.049988 | 51.549999 | 64.777512 | 32.255619 | 32.228062 | 29.394268 | 31.417585 |
| 2020-12-10 | -0.049988 | 51.549999 | 43.436817 | 44.310410 | 19.204155 | 17.022833 | 13.802833 |
| 2020-12-11 | -0.049988 | 51.549999 | 113.465919 | 44.140366 | 46.634701 | 40.187134 | 49.345261 |
| 2020-12-12 | -0.049988 | 51.549999 | 57.885963 | 45.814697 | 30.872692 | 25.869278 | 14.578404 |
| 2020-12-13 | -0.049988 | 51.549999 | 56.978699 | 55.796829 | 19.655144 | 17.986469 | 17.178650 |
| 2020-12-14 | -0.049988 | 51.549999 | 41.861675 | 75.315674 | 10.625697 | 8.364858 | 13.119390 |
| 2020-12-15 | -0.049988 | 51.549999 | 54.948017 | 60.603752 | 12.214766 | 9.671446 | 12.929519 |
| 2020-12-16 | -0.049988 | 51.549999 | 63.984142 | 56.931396 | 16.307283 | 14.952140 | 25.290646 |
| 2020-12-17 | -0.049988 | 51.549999 | 48.933174 | 49.170326 | 12.910172 | 10.861581 | 10.984763 |
| 2020-12-18 | -0.049988 | 51.549999 | 21.916876 | 64.970718 | NaN | NaN | NaN |
| 2020-12-19 | -0.049988 | 51.549999 | 19.733660 | 81.636444 | 12.670108 | 11.536630 | 4.229946 |
| 2020-12-20 | -0.049988 | 51.549999 | 25.818867 | 68.821381 | 8.510469 | 6.051936 | 6.062692 |
| 2020-12-21 | -0.049988 | 51.549999 | 26.262127 | 55.715931 | 9.830421 | 7.515805 | 6.226173 |
| 2020-12-22 | -0.049988 | 51.549999 | 69.275375 | 74.091499 | 38.756077 | 33.044914 | 25.561008 |
| 2020-12-23 | -0.049988 | 51.549999 | 51.012043 | 73.123276 | 16.138432 | 14.708934 | 11.811040 |
| 2020-12-24 | -0.049988 | 51.549999 | 12.582739 | 80.724190 | 5.150666 | 4.405396 | 3.355808 |
| 2020-12-25 | -0.049988 | 51.549999 | 51.906540 | 64.102203 | 21.049095 | 19.382658 | 15.916385 |
| 2020-12-26 | -0.049988 | 51.549999 | 46.186401 | 69.872673 | 17.390356 | 16.178095 | 7.168921 |
| 2020-12-27 | -0.049988 | 51.549999 | 71.503342 | 83.723030 | 14.515577 | 12.415962 | 13.618376 |
| 2020-12-28 | -0.049988 | 51.549999 | 67.383316 | 41.268040 | 35.509396 | 29.793512 | 18.922663 |
| 2020-12-29 | -0.049988 | 51.549999 | 65.241524 | 32.984646 | 41.049049 | 37.794231 | 24.576710 |
| 2020-12-30 | -0.049988 | 51.549999 | 84.480133 | 23.565016 | 36.276402 | 32.324852 | 39.448082 |
| 2020-12-31 | -0.049988 | 51.549999 | 52.926937 | 20.486465 | 33.245438 | 29.518814 | 15.313162 |
2.4 Classify daily maximum values of key pollutants into European Air Quality Index levels¶
The next step is now to classify the daily maximum values into the respective European Air Quality Index levels. The EAQI has six index levels and different thresholds for each of the five key pollutants. First, we define the limits for each pollutant and define additionally a list with the index class labels 1 to 6.
ozone_limits = [0, 50, 100, 130, 240, 380, 800]
no2_limits = [0, 40, 90, 120, 230, 340, 1000]
so2_limits = [0, 100, 200, 350, 500, 750, 1250]
pm25_limits = [0, 20, 40, 50, 100, 150, 1200]
pm10_limits = [0, 10, 20, 25, 50, 75, 800]
index_levels = [1, 2, 3, 4, 5, 6]
Based on the thresholds above, we now classify the daily maximum values of each pollutant into the respective European Air Quality Index level. The pandas function cut() allows us to categorize the daily maximum values into index levels. The function takes the following keyword arguments:
pandas.Series for each pollutant, which shall be categorized
list of thresholds for each pollutant
list with index labels
As a result, for each pollutant, the daily maximum value has been classified into one of the six index level classes, based on the pollutant thresholds.
ozone = pd.cut(eaqi_ts_daily_df['o3_conc'], ozone_limits , labels = index_levels)
no2 = pd.cut(eaqi_ts_daily_df['no2_conc'], no2_limits , labels = index_levels)
so2 = pd.cut(eaqi_ts_daily_df['so2_conc'], so2_limits , labels = index_levels)
pm25 = pd.cut(eaqi_ts_daily_df['pm2p5_conc'], pm25_limits , labels = index_levels)
pm10 = pd.cut(eaqi_ts_daily_df['pm10_conc'], pm10_limits , labels = index_levels)
In a next step, we want to bring together the categorized pandas.Series objects above into one pandas.DataFrame. For this, we can use the constructor pd.DataFrame() and combine the five pandas.Series objects into one DataFrame.
df = pd.DataFrame(dict(o3 = ozone, no2 = no2, so2 = so2, pm25=pm25, pm10=pm10))
| o3 | no2 | so2 | pm25 | pm10 | |
|---|---|---|---|---|---|
| time | |||||
| 2020-12-01 | 2 | 1 | 1 | 2 | 3 |
| 2020-12-02 | 1 | 2 | 1 | 3 | 5 |
| 2020-12-03 | 1 | 2 | 1 | 2 | 4 |
| 2020-12-04 | 2 | 2 | 1 | 1 | 2 |
| 2020-12-05 | 2 | 2 | 1 | 1 | 2 |
| 2020-12-06 | 1 | 1 | 1 | 2 | 4 |
| 2020-12-07 | 1 | 2 | 1 | 3 | 4 |
| 2020-12-08 | 1 | 2 | 1 | 3 | 4 |
| 2020-12-09 | 1 | 2 | 1 | 2 | 4 |
| 2020-12-10 | 1 | 2 | 1 | 1 | 2 |
| 2020-12-11 | 1 | 3 | 1 | 3 | 4 |
| 2020-12-12 | 1 | 2 | 1 | 2 | 4 |
| 2020-12-13 | 2 | 2 | 1 | 1 | 2 |
| 2020-12-14 | 2 | 2 | 1 | 1 | 2 |
| 2020-12-15 | 2 | 2 | 1 | 1 | 2 |
| 2020-12-16 | 2 | 2 | 1 | 1 | 2 |
| 2020-12-17 | 1 | 2 | 1 | 1 | 2 |
| 2020-12-18 | 2 | 1 | NaN | NaN | NaN |
| 2020-12-19 | 2 | 1 | 1 | 1 | 2 |
| 2020-12-20 | 2 | 1 | 1 | 1 | 1 |
| 2020-12-21 | 2 | 1 | 1 | 1 | 1 |
| 2020-12-22 | 2 | 2 | 1 | 2 | 4 |
| 2020-12-23 | 2 | 2 | 1 | 1 | 2 |
| 2020-12-24 | 2 | 1 | 1 | 1 | 1 |
| 2020-12-25 | 2 | 2 | 1 | 1 | 3 |
| 2020-12-26 | 2 | 2 | 1 | 1 | 2 |
| 2020-12-27 | 2 | 2 | 1 | 1 | 2 |
| 2020-12-28 | 1 | 2 | 1 | 2 | 4 |
| 2020-12-29 | 1 | 2 | 1 | 2 | 4 |
| 2020-12-30 | 1 | 2 | 1 | 2 | 4 |
| 2020-12-31 | 1 | 2 | 1 | 2 | 4 |
The last step before we can visualize the European Air Quality Index levels is to compute the overall index level based on the five pollutants. The overall index level for each day is defined by the pollutant with the highest level. Below, with the help of the function max(), we select the maximum index level for each day and define a new column with the name level.
df['level'] = df.max(axis=1).astype(int)
df
| o3 | no2 | so2 | pm25 | pm10 | level | |
|---|---|---|---|---|---|---|
| index | ||||||
| 2020-12-01 | 2 | 1 | 1 | 2 | 3 | 3 |
| 2020-12-02 | 1 | 2 | 1 | 3 | 5 | 5 |
| 2020-12-03 | 1 | 2 | 1 | 2 | 4 | 4 |
| 2020-12-04 | 2 | 2 | 1 | 1 | 2 | 2 |
| 2020-12-05 | 2 | 2 | 1 | 1 | 2 | 2 |
| 2020-12-06 | 1 | 1 | 1 | 2 | 4 | 4 |
| 2020-12-07 | 1 | 2 | 1 | 3 | 4 | 4 |
| 2020-12-08 | 1 | 2 | 1 | 3 | 4 | 4 |
| 2020-12-09 | 1 | 2 | 1 | 2 | 4 | 4 |
| 2020-12-10 | 1 | 2 | 1 | 1 | 2 | 2 |
| 2020-12-11 | 1 | 3 | 1 | 3 | 4 | 4 |
| 2020-12-12 | 1 | 2 | 1 | 2 | 4 | 4 |
| 2020-12-13 | 2 | 2 | 1 | 1 | 2 | 2 |
| 2020-12-14 | 2 | 2 | 1 | 1 | 2 | 2 |
| 2020-12-15 | 2 | 2 | 1 | 1 | 2 | 2 |
| 2020-12-16 | 2 | 2 | 1 | 1 | 2 | 2 |
| 2020-12-17 | 1 | 2 | 1 | 1 | 2 | 2 |
| 2020-12-18 | 2 | 1 | NaN | NaN | NaN | 2 |
| 2020-12-19 | 2 | 1 | 1 | 1 | 2 | 2 |
| 2020-12-20 | 2 | 1 | 1 | 1 | 1 | 2 |
| 2020-12-21 | 2 | 1 | 1 | 1 | 1 | 2 |
| 2020-12-22 | 2 | 2 | 1 | 2 | 4 | 4 |
| 2020-12-23 | 2 | 2 | 1 | 1 | 2 | 2 |
| 2020-12-24 | 2 | 1 | 1 | 1 | 1 | 2 |
| 2020-12-25 | 2 | 2 | 1 | 1 | 3 | 3 |
| 2020-12-26 | 2 | 2 | 1 | 1 | 2 | 2 |
| 2020-12-27 | 2 | 2 | 1 | 1 | 2 | 2 |
| 2020-12-28 | 1 | 2 | 1 | 2 | 4 | 4 |
| 2020-12-29 | 1 | 2 | 1 | 2 | 4 | 4 |
| 2020-12-30 | 1 | 2 | 1 | 2 | 4 | 4 |
| 2020-12-31 | 1 | 2 | 1 | 2 | 4 | 4 |
2.5 Visualize daily European Air Quality Index for London in December 2020 as heatmap¶
The last step is now to visualize the daily European Air Quality Index for London in December 2020 as heatmap. We can use the function heatmap() from the seaborn library to create a heatmap with days on the horizontal axis and the pollutants as well as the overall index level on the vertical axis. This representation allows for a quick interpretation on the overall index level for every day. It further allows you to identify the determining pollutant for every day.
For the heatmap() function, we need to transpose the dataframe, so that the key pollutants and overall index level are given as row entries and for every day in December 2020 a column entry. Additionally, with iloc[::-1], we want to revert the sequence of the pollutants, as we want to visualize the overall index level on top.
df = (df.T).iloc[::-1]
df
| index | 2020-12-01 | 2020-12-02 | 2020-12-03 | 2020-12-04 | 2020-12-05 | 2020-12-06 | 2020-12-07 | 2020-12-08 | 2020-12-09 | 2020-12-10 | ... | 2020-12-22 | 2020-12-23 | 2020-12-24 | 2020-12-25 | 2020-12-26 | 2020-12-27 | 2020-12-28 | 2020-12-29 | 2020-12-30 | 2020-12-31 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| level | 3 | 5 | 4 | 2 | 2 | 4 | 4 | 4 | 4 | 2 | ... | 4 | 2 | 2 | 3 | 2 | 2 | 4 | 4 | 4 | 4 |
| pm10 | 3 | 5 | 4 | 2 | 2 | 4 | 4 | 4 | 4 | 2 | ... | 4 | 2 | 1 | 3 | 2 | 2 | 4 | 4 | 4 | 4 |
| pm25 | 2 | 3 | 2 | 1 | 1 | 2 | 3 | 3 | 2 | 1 | ... | 2 | 1 | 1 | 1 | 1 | 1 | 2 | 2 | 2 | 2 |
| so2 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | ... | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 |
| no2 | 1 | 2 | 2 | 2 | 2 | 1 | 2 | 2 | 2 | 2 | ... | 2 | 2 | 1 | 2 | 2 | 2 | 2 | 2 | 2 | 2 |
| o3 | 2 | 1 | 1 | 2 | 2 | 1 | 1 | 1 | 1 | 1 | ... | 2 | 2 | 2 | 2 | 2 | 2 | 1 | 1 | 1 | 1 |
6 rows × 31 columns
The transpose function converts NaN values to a large negative number. For this reason, we have to set negative entries back to NaN (not a number).
df[df<0] = np.nan
df
| index | 2020-12-01 | 2020-12-02 | 2020-12-03 | 2020-12-04 | 2020-12-05 | 2020-12-06 | 2020-12-07 | 2020-12-08 | 2020-12-09 | 2020-12-10 | ... | 2020-12-22 | 2020-12-23 | 2020-12-24 | 2020-12-25 | 2020-12-26 | 2020-12-27 | 2020-12-28 | 2020-12-29 | 2020-12-30 | 2020-12-31 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| level | 3 | 5 | 4 | 2 | 2 | 4 | 4 | 4 | 4 | 2 | ... | 4 | 2 | 2 | 3 | 2 | 2 | 4 | 4 | 4 | 4 |
| pm10 | 3 | 5 | 4 | 2 | 2 | 4 | 4 | 4 | 4 | 2 | ... | 4 | 2 | 1 | 3 | 2 | 2 | 4 | 4 | 4 | 4 |
| pm25 | 2 | 3 | 2 | 1 | 1 | 2 | 3 | 3 | 2 | 1 | ... | 2 | 1 | 1 | 1 | 1 | 1 | 2 | 2 | 2 | 2 |
| so2 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | ... | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 |
| no2 | 1 | 2 | 2 | 2 | 2 | 1 | 2 | 2 | 2 | 2 | ... | 2 | 2 | 1 | 2 | 2 | 2 | 2 | 2 | 2 | 2 |
| o3 | 2 | 1 | 1 | 2 | 2 | 1 | 1 | 1 | 1 | 1 | ... | 2 | 2 | 2 | 2 | 2 | 2 | 1 | 1 | 1 | 1 |
6 rows × 31 columns
And now, we can finally visualize the heatmap. The code below has five main parts:
1. Initiate the matplotlib figure: Initiate a figure and axes objects and add a colorbar
2. Plotting function: Plot the pandas dataframe with the heatmap() function from the seaborn library
3. Set title of the heatmap: Set a title to the resulting heatmap
4. Customize x- and y-axis ticks and labels: Customize the x- and y-labels for correct representation of the data
5. Customize colorbar entries: Customize the labels of the colorbar and indicate the index levels from ‘very good’ to ‘extremely poor’
Note: before plotting, we additionally define a customized colormap based on the official EAQI index colors. We also add additionally the labels to the six EAQI index levels, which range from ‘Very good’ to ‘Extremely poor’.
cmap = ListedColormap(['#5AAA5F', '#A7D25C', '#ECD347', '#F5BE41', '#F09235', '#D93322'])
cmap
labels = ['Very good', 'Good', 'Medium', 'Poor', 'Very Poor', 'Extremely Poor']
# Initiate the matplotlib figure
fig, ax = plt.subplots(1,1,figsize=(35,5))
cbar_ax = fig.add_axes([.82, .13, .01, .75])
# Plotting function
g = sns.heatmap(df, cmap=cmap, linewidth=1, linecolor='w', square=True, cbar_ax = cbar_ax,vmin=1, vmax=6, ax=ax)
# Set title of the heatmap
g.set_title("\nDaily European Air Quality Index levels for London - December 2020\n", fontsize=20)
g.set(xlabel=None)
# Customize x- and y-axis ticks and labels
ylabels=['Index level', 'PM 10', 'PM 2.5', 'Sulphur dioxide', 'Nitrogen dioxide', 'Ozone']
g.set_yticklabels(ylabels,fontsize=14, rotation=0)
xlabels=df.columns.format('%Y-%m-%d')[1::]
g.set_xticklabels(xlabels, fontsize=14)
# Customize colorbar entries
cbar = ax.collections[0].colorbar
cbar.set_ticks([1.4,2.2,3.1,3.9,4.8,5.6])
cbar.set_ticklabels(['Very good', 'Good', 'Medium', 'Poor', 'Very Poor', 'Extremely Poor'])
cbar.ax.tick_params(labelsize=14)
This project is licensed under APACHE License 2.0. | View on GitHub