<jupyter_start><jupyter_text>**Create and validate 31-day running mean climatology for ERA-5 v1.1 tas, tmax, tmin.**<jupyter_code>! pip install git+https://github.com/ClimateImpactLab/xclim.git@qdm_add_year_selection
%matplotlib inline

import xarray as xr
import numpy as np
import os as os
import matplotlib.pyplot as plt
from cartopy import config
import cartopy.crs as ccrs
import matplotlib.lines as mlines
import matplotlib.patches as mpatches
import pandas as pd
import itertools
import datetime 
import toolz
from datetime import date
from pathlib import Path, PurePath
from xclim.core.calendar import convert_calendar

import sys

import dask.distributed as dd
import dask
import rhg_compute_tools.kubernetes as rhgk
client, cluster = rhgk.get_big_cluster()
cluster.scale(10)
client
cluster.close()
yrs = np.arange(1994,2016) 
def pull_ERA5_variable(variable):
    filenames = []
    
    for num_yrs in range(len(yrs)):
        filename = '/gcs/impactlab-data/climate/source_data/ERA-5/{}/daily/netcdf/v1.1/{}_daily_{}-{}.nc'.format(variable, variable, yrs[num_yrs], yrs[num_yrs])
        filenames.append(filename)
    
    era5_var = xr.open_mfdataset(filenames, 
                               concat_dim='time', combine='by_coords')
    
    var_all = era5_var[variable]
    
    return var_all

da = pull_ERA5_variable('tas')
da = da.persist()
# remove leap days and convert calendar to no-leap
da = convert_calendar(da, 'noleap')
def assign_coords(da):
    years = np.arange(da.time.dt.year.min(),da.time.dt.year.max()+1)
    
    da_wcoords = da.assign_coords(dayofyear=xr.DataArray(np.array([np.arange(1,366)]*len(years)).flatten(),
                                                dims=('time'),
                                                coords={'time':da.time})).persist()
    return da_wcoords
da = assign_coords(da)
da_dayofyear_mean = da.sel(time=slice('1994-12-17','2015-01-15')).groupby('dayofyear').mean().persist()
da_dayofyear_mean.sel(latitude=35.0,longitude=180.0).plot()
da_dayofyear_mean.nbytes / (1024**3)
da_dayofyear_mean = da_dayofyear_mean.load()
def rolling_31day_mean(da):
    roll = da.pad(dayofyear=15, mode='wrap').rolling(dayofyear=31, center=True).mean().dropna('dayofyear')
    return roll
da_rolling_mean = rolling_31day_mean(da_dayofyear_mean)
da_rolling_mean
da_rolling_mean.sel(latitude=35.0,longitude=180.0).plot()<jupyter_output><empty_output><jupyter_text>Before creating the netcdf file, run a quick validation of the climatology.<jupyter_code>def test_for_nans(da):
    # no nans
    assert da.isnull().sum() == 0, "there are nans!"

def test_temp_range(da):
    # make sure temp values are in a valid range
    # asserts if statement below is false
    assert (da.min() > 150) or (ds[var].max() < 350), "temperature values are invalid" 
test_for_nans(da_rolling_mean)
test_temp_range(da_rolling_mean)
da_rolling_mean.mean(dim='dayofyear').plot()
def create_netcdf(da, variable):
    '''
    create netcdf file.
    '''
    out_direc = ('/gcs/impactlab-data/climate/source_data/ERA-5/{}/climatology/'.format(variable))
    out_filename = ('{}_1995_2015_climo.nc'.format(variable))
    
    if not Path(out_direc).exists():
                    Path(out_direc).mkdir(parents=True, exist_ok=True)
            
    # add attributes to dataset and save file 
    attrsdt = {
        'author': '<NAME>',
        'contact': '<EMAIL>',
        'project': ('downscaling CMIP6'),
        'source': ('/gcs/climate/source_data/ERA-5/day/tmin/v1.1'),
        'method': ('31-rolling mean climatology, no leap days'),
        'created': str(date.today()),
        'units': 'K'}

    da.attrs.update(attrsdt)
    
    # save file
    da.to_netcdf(Path(PurePath(out_direc).joinpath(out_filename)))
create_netcdf(da_rolling_mean, 'tas')<jupyter_output><empty_output><jupyter_text>**Validation back for the `pad` function**<jupyter_code>validate_pad = da_dayofyear_mean.sel(latitude=35.0,longitude=180.0).pad(dayofyear=15, mode='wrap')
validate_pad[:15]
validate_pad[365:380]
validate_pad[15:30]
validate_pad[380:]<jupyter_output><empty_output><jupyter_text>the `pad` function should wrap the last values to the beginning and the beginning values to the end. Given the original dataset of 1-365 we then add 15 days from the end to the beginning, which match the 351-365 days. Additionally, the 1-15 days match the days at the end of the new 'padded' `DataArray`. `Pad` function is validated. **Validation for the `rolling` function**<jupyter_code>validate_roll = da_dayofyear_mean.sel(latitude=35.0,longitude=180.0).rolling(dayofyear=31, center=True).mean().dropna('dayofyear')
validate_roll
validate_roll.sel(dayofyear=20)
validate_roll_v2 = da_dayofyear_mean.sel(dayofyear=slice(5,35),latitude=35.0,longitude=180.0)
validate_roll_v2
validate_roll_v2.mean()<jupyter_output><empty_output>