1

I need to calculate monthly, seasonal and annualy anomalies of air temperature of netcdf monthly files of 44 years with some function that allows me to obtain anomalies in the period on a monthly, seasonal and annualy automatically and save the results in a folder. I only know how to do It for one year and not for several years with a function.

from netCDF4 import Dataset, num2date
import xarray as xr
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as cfeat
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from mpl_toolkits.basemap import Basemap

ds = Dataset('./interim_t2m_19792017.nc')
lats = ds.variables['latitude'][:]  # extract/copy the data
lons = ds.variables['longitude'][:]
time = ds.variables['time']
var = ds.variables['t2m'] 

lon, lat = np.meshgrid(lons, lats)
dates = num2date(time[:], time.units)
dates_pd = pd.to_datetime(dates)
periods = dates_pd.to_period(freq='M')

def plt_map(data):
    m = Basemap(projection='mill',llcrnrlat=-80,urcrnrlat=80,\
            llcrnrlon=0,urcrnrlon=360,lat_ts=20,resolution='c')
    x, y = m(lon, lat)
    plt.figure(figsize=(10,7))
    m.drawcoastlines()
    m.drawparallels(np.arange(-80.,81.,20.))
    m.drawmeridians(np.arange(-180.,181.,20.))
    m.drawmapboundary(fill_color='white')
    m.contourf(x,y,data, extend="both",cmap="jet");
    plt.colorbar(orientation='horizontal', pad=0.05)
plt_map(var[0,:,:])

mask_2016 = periods.year==2016
data = var[mask_2016,:,:].mean(axis=0)
plt_map(data)
tmsppc
  • 107
  • 3
  • 10
  • What do you consider as anomaly? Is this difference from mean e.g for annual anomaly difference from annual means and for monthly anomaly (January) difference from mean on particular month values (mean January)? – msi_gerva Sep 23 '18 at 12:28
  • Yes you're right. I need to do your comment in a function for several years (44). – tmsppc Sep 23 '18 at 13:16

3 Answers3

3

I know you are looking for a python answer, but this is bread and butter of CDO (climate data operators), that allows you to do these sort of calculations in one or two commands from the terminal window.

For example, to get the annual means of your Era Interim data you can do 

cdo yearmean interim_t2m_19792017.nc erai_yearmean.nc

and then to calculate the annual anomaly you need to do the long term average and subtract it 

cdo timmean interim_t2m_19792017.nc erai_timemean.nc
cdo sub erai_yearmean.nc erai_timemean.nc yearanom.nc

You can combine all these above 3 commands using "piping", but I keep them separate here as it is easier to see what is going on.

you can get the mean monthly seasonal cycle with:

cdo ymonmean interim_t2m_19792017.nc erai_ymonmean.nc

this gives you a file with the average of all the januarys, Feb etc (12 time slices). And then you can calculate the monthly anomaly, each with respect to its own monthly mean with 

cdo monmean interim_t2m_19792017.nc erai_monmean.nc
cdo sub erai_monmean.nc erai_ymonmean.nc erai_monanom.nc

There are also functions for seasonal averages.

See the online documentation for further details: https://code.mpimet.mpg.de/projects/cdo/

Lastly, msi_gerva is correct in the comment, it is not clear in the question what the anomalies are with respect to, as you could also calculate monthly anomalies with respect to the annual or long term mean. Moreover, you ask for annual anomalies and say you only know how to do it for one year, but I don't think that makes much sense, as the anomalies would be zero. It may be helpful to clarify the question more precisely.

ClimateUnboxed
  • 7,106
  • 3
  • 41
  • 86
  • Thanks Adrian for your answer but I don't know yet how to use CDO in Anaconda on Windows. I've installed python-cdo packages but it's not working reading your lines. – tmsppc Sep 26 '18 at 00:49
  • hi there, I think CDO recommends using cygwin for direct use in Windows, see https://code.mpimet.mpg.de/projects/cdo/wiki/Win32, but you should be able to do it also using the python package. What is the issue, the installation of the package or trying to run those lines? If you let me know what the error is perhaps I can try to help – ClimateUnboxed Sep 26 '18 at 07:10
1

I think you cannot calculate the anomalies for different periods with just one call of some function, you need to do that in several steps. For example, if you have all the dates in vector datevec and data in var, you should be able to do:

a) for monthly means:

for mon in range(1,13):
    kkmon = [ii for ii,val in enumerate(datevec) if val.month == mon]
    monmean = var[kkmon,:,:].mean(axis=0)

b) for seasonal means:

seasons = {'DJF':[12,1,2],'MAM':[3,4,5],'JJA':[6,7,8],'SON':[9,10,11]}
for seaskey in seasons.keys():
    kkseas = [ii for ii,val in enumerate(datevec) if val.month in seasons[seaskey]]
    seasmean = var[kkseas,:,:].mean(axis=0)

This is just an example how you could calculate the overall means for different months and seasons. You have to calculate the anomalies based on the data and resolution that you have.

I would also suggest looking at possibilities of cdo (Climate Data Operators) as most of the averaging methods are built-in there and as it is written in C or C++, it works much faster compared to the Python. Of course, you can combine cdo and Python for your tasks - use the first to find the means and the second to calculate and plot anomalies afterwards.

msi_gerva
  • 2,021
  • 3
  • 22
  • 28
1

On Linux and macOS this can be done using nctoolkit (https://nctoolkit.readthedocs.io/en/latest/). If you want an annual anomaly, you can do the following:

import nctoolkit as nc
data = nc.open_data("./interim_t2m_19792017.nc")
data.annual_anomaly(baseline = [1997, 2017])

and for monthly anomalies:

import nctoolkit as nc
data = nc.open_data("./interim_t2m_19792017.nc")
data.monthly_anomaly(baseline = [1997, 2017])

The package does not have seasonal anomalies, but I will add that in a future release. Under the hood, nctoolkit uses CDO so the system commands will be similar to those mentioned by Adrian Tompkins.

Robert Wilson
  • 3,192
  • 11
  • 19