title: Calculate yearly anomalies for several years on netcdf files in python
body: <p>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.</p>
  
  <pre><code>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)
  </code></pre>