Sentinel-5P and xarray

I’ve written about Sentinel-5P several times in the past. You can check these posts out below:

Why is Sentinel-5P important?

Sentinel-5P and Python

The Sentinel-5P and Python post is an outline on how to read in a Sentinel-5P netcdf file and plot it. The notebook (code) is available here.

I wanted to expand on this notebook and see if I could read a Sentinel-5P dataset into xarray and then visualise it. I am pretty new to xarray, but it seems to be a nice fit into a Geospatial Python toolbox. I was fortunate enough to chat with Robin Wilson about xarray on the latest #scenefromabove podcast (S5E1). Talking with him really opened my eyes to the potential of xarray. Check out his slides and code for more detail about xarray.

Setting up

If you have Anaconda python distribution then you can install xarray at the anaconda prompt with the command

conda install -c anaconda xarray

and for this post we will plot using Cartopy (unlike previously where I used basemap – which is being discontinued). You can get Cartopy on Anaconda using this command

conda install -c conda-forge cartopy

As an aside I saw this guide about cartopy this week; you should check it out for a brilliant overview:

https://rabernat.github.io/research_computing_2018/maps-with-cartopy.html

Reading the data and plotting

You can get Sentinel-5P data from the Copernicus scihub here. In this case I have downloaded a level 2 NO2 dataset – it is supplied in netcdf format. Have a look at my previous posts if you want to know more.

To read the data into xarray use the following code:

import xarray as xr
data5p = xr.open_dataset('S5P_NRTI_L2__NO2_xxx.nc', group='PRODUCT')
print(data5p)

This will print the contents of your xarray. To get a specific NO2 dataset specify the data variable you want. For example:

no2 = (data5p['nitrogendioxide_tropospheric_column_precision'])
print(no2)

This gives me a data array of the NO2 tropospheric column precision variable. In my previous post we plotted this data using the code shown below by first converting it to a numpy array and then plotting it.

### Old netcdf to array code from 
### http://www.acgeospatial.co.uk/sentinel-5p-and-python/
no2_data = data5p['nitrogendioxide_tropospheric_column_precision'][0,:,:]
no2_data.plot()

This is ok I suppose, but it doesn’t give me the projection. In fact the axis are labelled along-track and across track dimension – which isn’t really what I want.

The xarray data structure is shown in the diagram below:

Image courtesy of http://xarray.pydata.org/en/stable/data-structures.html

In Sentinel-5P xarray swap ‘temperature’ and ‘precipitation’ in the diagram above with the data variables like ‘nitrogendioxide_tropospheric_column’ and ‘nitrogendioxide_tropospheric_column_precision’. So, where are the latitude and longitude values?

We can find them here:

data5p.latitude.attrs
data5p.longitude.attrs

And, we can use them to plot out data like this:

from matplotlib import pyplot as plt
plt.figure(figsize=(14,6))
ax = plt.axes()


data5p.nitrogendioxide_tropospheric_column_precision[0].plot.pcolormesh \
(ax=ax, x='longitude', y='latitude', add_colorbar=True, cmap='jet')

giving us:

We can use cartopy to give this data some context:

import cartopy.crs as ccrs

ax = plt.axes(projection=ccrs.PlateCarree())
data5p.nitrogendioxide_tropospheric_column_precision[0].plot.pcolormesh(ax=ax, x='longitude', y='latitude', add_colorbar=True, cmap='jet')
ax.set_global(); ax.coastlines()

giving us:

Or, on a globe:

ax = plt.axes(projection=ccrs.Orthographic(0, 30))


data5p.nitrogendioxide_tropospheric_column_precision[0].plot.pcolormesh(ax=ax, x='longitude', y='latitude', add_colorbar=True, cmap='jet', transform=ccrs.PlateCarree())
ax.set_global(); ax.coastlines()

giving us:

Which is a pretty nice display, I think.

This is just the beginning for me with xarray, but it does seem perfectly suited for handling Sentinel-5P data.

Recap

I this post I have looked at:

  • Installing xarray
  • Loading a Sentinel-5P dataset into xarray
  • Plotting the data using the lat and long
  • Plotting the data with cartopy

The associated notebook to this post is available here


I am a freelancer able to help you with your projects. I offer consultancy, training and writing. I’d be delighted to hear from you. Please check out the books I have written on QGIS 3.4

https://www.packtpub.com/application-development/learn-qgis-fourth-edition

https://www.packtpub.com/application-development/qgis-quick-start-guide


I am @map_andrew on twitter