BiasAdjust_model_run.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "a83b8060-cdb7-4a14-96ce-b94fbf2417f8",
   "metadata": {},
   "source": [
    "# Description"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6dd2611b-6797-4643-ac5b-689d2512ab92",
   "metadata": {},
   "source": [
    "This code is written to bias adjust maximum temperature values from a global climate model run, calibrating with global overland gridded observations from Berkeley Earth. For simplicity, in this notebook we focus on a single user-chosen example model run, bias adjusting the '**forced**' (i.e. concatenated historical+RCP8.5) and '**natural**' Climate Model Intercomparison Project 5 (CMIP5) runs available in \\*.zarr format (saved out via XXXXX.ipynb).\n",
    "\n",
    "Bias adjustment is performed based on the methodology of [Lange (2019)](doi.org/10.5194/gmd-12-3055-2019); functions in the [cmip5_bias_adjustment.py](./cmip5_bias_adjustment.py). Accordingly, this notebook is licenced following GNU [Affero General Public License v3.0](https://www.gnu.org/licenses/agpl-3.0.en.html).\n",
    "\n",
    "Bias-adjusted model runs are saved out into a set of data files that mimic the structure of the inputs. Example distributions from a single month are output to visualize the results in accomapnying [BiasAdjust_comparison.ipynb](./BiasAdjust_comparison.ipynb) notebook."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9e37a7ea-0ce0-4e44-935f-7f03383b70be",
   "metadata": {},
   "source": [
    "# Setup"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "ae725422-5cc7-46fb-b504-e36875f483b8",
   "metadata": {},
   "outputs": [],
   "source": [
    "# import the packages we need to run this code\n",
    "import xarray as xr\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import utilities as util\n",
    "import pickle\n",
    "\n",
    "# import the bias-adjustment module modified from Lange (2019)\n",
    "import cmip5_bias_adjustment as ba"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "423de42b-d459-45d9-9fff-62b89f02f852",
   "metadata": {},
   "source": [
    "# Key considerations\n",
    "\n",
    "## Bias-adjustment parameters"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2e386c6d-4873-414e-8c32-b086470d17dd",
   "metadata": {},
   "source": [
    "The user chooses a set of parameters used to guide the bias-adjustment process, including the calibration ($cal$) and adjustment ($adj$) periods."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "725fbc2c-0920-4548-a5d9-71db8f7f13e5",
   "metadata": {},
   "outputs": [],
   "source": [
    "# choose the environmental state variable, the metric, we are looking at\n",
    "metric='TMAX'\n",
    "\n",
    "# define the training/calibration period which will use the relationship\n",
    "# between the observations and the model to guide the calibration over the\n",
    "# the adjustment period\n",
    "calibration_period=['1985-01-01','2015-12-31']\n",
    "\n",
    "# define the adjustment period for the forced and natural simulations \n",
    "forced_adj_period=['1885-01-01','2050-12-31']\n",
    "natural_adj_period=['1975-01-01','2005-12-31']"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "dc85191e-e588-4bc9-96a9-d8953650b20d",
   "metadata": {},
   "source": [
    "## Example Month for Comparison"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "be964734-cced-4881-8dc9-a5c9163106a8",
   "metadata": {},
   "source": [
    "To facilitate a (visual and statistical) comparison analysis between the observations, raw climate model data, and bias-adjusted model data, we select a single example month to explore."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "754a72ee-5d7d-4deb-979c-2788b6d9a3d5",
   "metadata": {},
   "outputs": [],
   "source": [
    "# for our comparison, take July\n",
    "eg_month=7\n",
    "\n",
    "# initialize comparison dictionary to store the results in\n",
    "eg_storage = {}\n",
    "eg_storage['month']=eg_month\n",
    "\n",
    "# locate the directory for saving the example dictionary\n",
    "eg_savepath='./data/cmip5/ba_comparison/'"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6a9588d1-ddd7-478e-ad21-f0c0bd37ed5d",
   "metadata": {},
   "source": [
    "# Load Observations"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8e6ca26f-76f9-47bd-9065-7450b7090bdd",
   "metadata": {},
   "source": [
    "Gridded Berkeley Earth observations are used to calibrate the CMIP5 global climate model runs.\n",
    "\n",
    "Observations are provided [Berkeley Earth](http://berkeleyearth.org/data/) (BE), and are identical to those used in the observational analysis of this work (e.g. XXXXX.ipynb). BE daily maximum temeprature data were first downloaded by Karsten Haustein in March 2019 from the [Climate Explorer](https://climexp.knmi.nl/select.cgi?id=someone@somewhere&field=berkeley_tmax_daily) database provided by the World Meteorological Organization; regridding to the standard N96 raster was performed by K. Haustein using [CDO tools](https://code.mpimet.mpg.de/projects/cdo/) following: ```cdo remapbil,r192x145```."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "74ce58b3-fcf4-46af-aa8a-4399764f7d0f",
   "metadata": {},
   "outputs": [],
   "source": [
    "# load the daily observed data from its path\n",
    "obs_path='./data/berkeley/Berkeley_'+metric+'_land_188001_201712_daily_N96_365days_degC.nc'\n",
    "obsdat=xr.open_dataset(obs_path)\n",
    "\n",
    "# reset observed daily data time grid\n",
    "timegridobs=util.get_dt64(obsdat)\n",
    "obsdat['time']=timegridobs"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b6549279-f2a5-4844-985c-9748577183c6",
   "metadata": {},
   "source": [
    "### Define observation training data during the calibration period"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "88ac2942-827d-4db6-9249-f5cc30dc6168",
   "metadata": {},
   "source": [
    "The climate model data will be bias-adjusted over the calibration period ($cal$) following the relationship between each model, $x_{cal}^{sim}$, and the Berkeley Earth observations, $x_{cal}^{obs}$, during that time. Here we take the calibration period to be 1985&mdash;2015 (as defined above).\n",
    "\n",
    "The adjustment period ($adj$) is taken to be either the entire concatenated historic/future period of record 1885&mdash;2050 in the case of forced simulations, or the last 31-years of analysis in case of natural simulations (that will make up the counterfactual, see XXXX.ipynb)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "8e350f37-7956-4dfa-b5d0-27867526d082",
   "metadata": {},
   "outputs": [],
   "source": [
    "# subset the data between the define training period endpoints\n",
    "xobscal=obsdat.sel(time=slice(calibration_period[0],calibration_period[-1]))\n",
    "\n",
    "# find the month indices associated with the training period note these\n",
    "# would identical between the model and observations b/c they use the \n",
    "# same training period, but if loaded separately we have to find separately\n",
    "monidx_calperiod=xobscal.groupby('time.month').groups"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "efb2297c-09bd-4357-a8e7-8c0069a127f4",
   "metadata": {},
   "source": [
    "# Setup load of (to-be-adjusted) CMIP5 TMAX data"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e1b53de6-1a97-4151-973f-19d8f2ab7426",
   "metadata": {},
   "source": [
    "### Approach\n",
    "\n",
    "As an example of the bias-adjustment process we will load in a single paired model, forced and natural, with both scenarios to be bias-adjusted. In the whole workflow, this process is repeated for each CMIP5 model (paired and unparied) in this study. This could be done directly with this script (by replacing the input name of the model) or through some paralellization (e.g. through [dask](https://dask.org/)).\n",
    "\n",
    "For our purposes, we illustrate the methods (and comparison with the observations in XXXX.ipynb) with a single model: **GFDL-CM3**.\n",
    "\n",
    "### File provenance and organization\n",
    "\n",
    "In November 2018 Karsten Haustein downloaded the CMIP5 TMAX data from the [Centre for Environmental Analysis (CEDA)](https://help.ceda.ac.uk/article/4465-cmip5-data), regridded the files onto a common N96 grid, and concatenated historical+RCP8.5 timeseries before saving out the model runs into \\*.nc formats. To organize our analysis herein, Daniel Gilford added datetime64 information and reorganized the files into individual monthly \\*.zarr formats in July 2021."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "6bfc1670-90ce-4ec1-883b-c1bd58dbde33",
   "metadata": {},
   "outputs": [],
   "source": [
    "# define the example model's name\n",
    "mname='GFDL-CM3'\n",
    "# store it in the example dictionary\n",
    "eg_storage['model']=mname\n",
    "\n",
    "# define the model historical and natural filepaths and names\n",
    "cmip5_raw_filepath='./data/cmip5/raw_TMAX/'\n",
    "forced_zarrpath=mname+'_forced/'\n",
    "nat_zarrpath=mname+'_natural/'\n",
    "\n",
    "# define the filepaths and name to save the BA model data into\n",
    "cmip5_ba_filepath='./data/cmip5/ba_TMAX/'\n",
    "\n",
    "# define the month names and month integer array\n",
    "monzarr='mon{:02d}.zarr'\n",
    "months=np.arange(1,13)\n",
    "\n",
    "# define the path to save the monthly bias-adjusted data to\n",
    "savepath='./data/cmip5/ba_TMAX/'"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2917c15c-bbd1-461c-bcbe-25b6c877dba0",
   "metadata": {},
   "source": [
    "# Apply Bias Adjustment"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e555b08f-65db-4c26-85df-f2a6c6b3b5b0",
   "metadata": {},
   "source": [
    "Bias adjustment is applied month-by-month for each model, and takes about "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "6ba1b0bb-3414-4b5c-bb7a-46f45f2a12a9",
   "metadata": {},
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "object of type 'NoneType' has no len()",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mKeyboardInterrupt\u001b[0m                         Traceback (most recent call last)",
      "    \u001b[0;31m[... skipping hidden 1 frame]\u001b[0m\n",
      "\u001b[0;32m/tmp/ipykernel_11765/1877972015.py\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m     24\u001b[0m     \u001b[0;31m# using the geospatial xarray wrapper defined in cmip5_bias_adjustment.py\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 25\u001b[0;31m     \u001b[0mmi_xsimadj_BA\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mba\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mxr_bias_adjust_singlemonth\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmi_xobscal\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mmi_xsimcal\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mmi_xsimadj\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m     26\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/gilford22_attframework/cmip5_bias_adjustment.py\u001b[0m in \u001b[0;36mxr_bias_adjust_singlemonth\u001b[0;34m(xobscal, xsimcal, xsimadj)\u001b[0m\n\u001b[1;32m     42\u001b[0m     \u001b[0;31m# using the xr.apply_ufunc feature\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 43\u001b[0;31m     result = xr.apply_ufunc(\n\u001b[0m\u001b[1;32m     44\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/xarray/core/computation.py\u001b[0m in \u001b[0;36mapply_ufunc\u001b[0;34m(func, input_core_dims, output_core_dims, exclude_dims, vectorize, join, dataset_join, dataset_fill_value, keep_attrs, kwargs, dask, output_dtypes, output_sizes, meta, dask_gufunc_kwargs, *args)\u001b[0m\n\u001b[1;32m   1152\u001b[0m     \u001b[0;32melif\u001b[0m \u001b[0many\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mis_dict_like\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0ma\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0ma\u001b[0m \u001b[0;32min\u001b[0m \u001b[0margs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 1153\u001b[0;31m         return apply_dataset_vfunc(\n\u001b[0m\u001b[1;32m   1154\u001b[0m             \u001b[0mvariables_vfunc\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/xarray/core/computation.py\u001b[0m in \u001b[0;36mapply_dataset_vfunc\u001b[0;34m(func, signature, join, dataset_join, fill_value, exclude_dims, keep_attrs, *args)\u001b[0m\n\u001b[1;32m    446\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 447\u001b[0;31m     result_vars = apply_dict_of_variables_vfunc(\n\u001b[0m\u001b[1;32m    448\u001b[0m         \u001b[0mfunc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0msignature\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0msignature\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mjoin\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mdataset_join\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mfill_value\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mfill_value\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/xarray/core/computation.py\u001b[0m in \u001b[0;36mapply_dict_of_variables_vfunc\u001b[0;34m(func, signature, join, fill_value, *args)\u001b[0m\n\u001b[1;32m    390\u001b[0m     \u001b[0;32mfor\u001b[0m \u001b[0mname\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mvariable_args\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mzip\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mnames\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mgrouped_by_name\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 391\u001b[0;31m         \u001b[0mresult_vars\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mname\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfunc\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mvariable_args\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    392\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/xarray/core/computation.py\u001b[0m in \u001b[0;36mapply_variable_ufunc\u001b[0;34m(func, signature, exclude_dims, dask, output_dtypes, vectorize, keep_attrs, dask_gufunc_kwargs, *args)\u001b[0m\n\u001b[1;32m    732\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 733\u001b[0;31m     \u001b[0mresult_data\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfunc\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0minput_data\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    734\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/numpy/lib/function_base.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, *args, **kwargs)\u001b[0m\n\u001b[1;32m   2303\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 2304\u001b[0;31m         \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_vectorize_call\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfunc\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mfunc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0margs\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mvargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m   2305\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/numpy/lib/function_base.py\u001b[0m in \u001b[0;36m_vectorize_call\u001b[0;34m(self, func, args)\u001b[0m\n\u001b[1;32m   2377\u001b[0m         \u001b[0;32mif\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0msignature\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 2378\u001b[0;31m             \u001b[0mres\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_vectorize_call_with_signature\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfunc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0margs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m   2379\u001b[0m         \u001b[0;32melif\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0margs\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/numpy/lib/function_base.py\u001b[0m in \u001b[0;36m_vectorize_call_with_signature\u001b[0;34m(self, func, args)\u001b[0m\n\u001b[1;32m   2417\u001b[0m         \u001b[0;32mfor\u001b[0m \u001b[0mindex\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mndindex\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mbroadcast_shape\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 2418\u001b[0;31m             \u001b[0mresults\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfunc\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0marg\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mindex\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0marg\u001b[0m \u001b[0;32min\u001b[0m \u001b[0margs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m   2419\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/gilford22_attframework/cmip5_bias_adjustment.py\u001b[0m in \u001b[0;36mbias_adjust_singleloc\u001b[0;34m(xobscal_loc, xsimcal_loc, xsimadj_loc, cal_yrgrid, adj_yrgrid)\u001b[0m\n\u001b[1;32m    103\u001b[0m     \u001b[0;31m# detrend the timeseries at this location for input to the Lange BA\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 104\u001b[0;31m     \u001b[0m_\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mdetrend_ts\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mtrend_ts\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mdetrend_each_ts\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mxobscal_loc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mxsimcal_loc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mxsimadj_loc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mcal_yrgrid\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0madj_yrgrid\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    105\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/gilford22_attframework/cmip5_bias_adjustment.py\u001b[0m in \u001b[0;36mdetrend_each_ts\u001b[0;34m(x_obs_cal, x_sim_cal, x_sim_fut, cal_yrgrid, adj_yrgrid)\u001b[0m\n\u001b[1;32m    280\u001b[0m     \u001b[0mx_sim_cal_DETREND\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mx_sim_cal_TREND\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0msubtract_or_add_trend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mx_sim_cal\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mcal_yrgrid\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtrend\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;32mNone\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 281\u001b[0;31m     \u001b[0mx_sim_fut_DETREND\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mx_sim_fut_TREND\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0msubtract_or_add_trend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mx_sim_fut\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0madj_yrgrid\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtrend\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;32mNone\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    282\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/gilford22_attframework/cmip5_bias_adjustment.py\u001b[0m in \u001b[0;36msubtract_or_add_trend\u001b[0;34m(x, years, trend)\u001b[0m\n\u001b[1;32m    334\u001b[0m     \u001b[0;32mif\u001b[0m \u001b[0mtrend\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 335\u001b[0;31m         \u001b[0mannual_means\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0marray\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmean\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0myears\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0my\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0my\u001b[0m \u001b[0;32min\u001b[0m \u001b[0munique_years\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    336\u001b[0m         \u001b[0ma\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mb\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mlinlstsq\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0munique_years\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mannual_means\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/gilford22_attframework/cmip5_bias_adjustment.py\u001b[0m in \u001b[0;36m<listcomp>\u001b[0;34m(.0)\u001b[0m\n\u001b[1;32m    334\u001b[0m     \u001b[0;32mif\u001b[0m \u001b[0mtrend\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 335\u001b[0;31m         \u001b[0mannual_means\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0marray\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmean\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0myears\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0my\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0my\u001b[0m \u001b[0;32min\u001b[0m \u001b[0munique_years\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    336\u001b[0m         \u001b[0ma\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mb\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mlinlstsq\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0munique_years\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mannual_means\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/numpy/core/overrides.py\u001b[0m in \u001b[0;36mmean\u001b[0;34m(*args, **kwargs)\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/numpy/core/fromnumeric.py\u001b[0m in \u001b[0;36mmean\u001b[0;34m(a, axis, dtype, out, keepdims, where)\u001b[0m\n\u001b[1;32m   3473\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 3474\u001b[0;31m     return _methods._mean(a, axis=axis, dtype=dtype,\n\u001b[0m\u001b[1;32m   3475\u001b[0m                           out=out, **kwargs)\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/numpy/core/_methods.py\u001b[0m in \u001b[0;36m_mean\u001b[0;34m(a, axis, dtype, out, keepdims, where)\u001b[0m\n\u001b[1;32m    188\u001b[0m         \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 189\u001b[0;31m             \u001b[0mret\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mret\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdtype\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtype\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mret\u001b[0m \u001b[0;34m/\u001b[0m \u001b[0mrcount\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    190\u001b[0m     \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mKeyboardInterrupt\u001b[0m: ",
      "\nDuring handling of the above exception, another exception occurred:\n",
      "\u001b[0;31mAttributeError\u001b[0m                            Traceback (most recent call last)",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/IPython/core/interactiveshell.py\u001b[0m in \u001b[0;36mshowtraceback\u001b[0;34m(self, exc_tuple, filename, tb_offset, exception_only, running_compiled_code)\u001b[0m\n\u001b[1;32m   2076\u001b[0m                         \u001b[0;31m# in the engines. This should return a list of strings.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 2077\u001b[0;31m                         \u001b[0mstb\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mvalue\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_render_traceback_\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m   2078\u001b[0m                     \u001b[0;32mexcept\u001b[0m \u001b[0mException\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mAttributeError\u001b[0m: 'KeyboardInterrupt' object has no attribute '_render_traceback_'",
      "\nDuring handling of the above exception, another exception occurred:\n",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "    \u001b[0;31m[... skipping hidden 1 frame]\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/IPython/core/interactiveshell.py\u001b[0m in \u001b[0;36mshowtraceback\u001b[0;34m(self, exc_tuple, filename, tb_offset, exception_only, running_compiled_code)\u001b[0m\n\u001b[1;32m   2077\u001b[0m                         \u001b[0mstb\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mvalue\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_render_traceback_\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   2078\u001b[0m                     \u001b[0;32mexcept\u001b[0m \u001b[0mException\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 2079\u001b[0;31m                         stb = self.InteractiveTB.structured_traceback(etype,\n\u001b[0m\u001b[1;32m   2080\u001b[0m                                             value, tb, tb_offset=tb_offset)\n\u001b[1;32m   2081\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/IPython/core/ultratb.py\u001b[0m in \u001b[0;36mstructured_traceback\u001b[0;34m(self, etype, value, tb, tb_offset, number_of_lines_of_context)\u001b[0m\n\u001b[1;32m   1365\u001b[0m         \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   1366\u001b[0m             \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtb\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mtb\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 1367\u001b[0;31m         return FormattedTB.structured_traceback(\n\u001b[0m\u001b[1;32m   1368\u001b[0m             self, etype, value, tb, tb_offset, number_of_lines_of_context)\n\u001b[1;32m   1369\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/IPython/core/ultratb.py\u001b[0m in \u001b[0;36mstructured_traceback\u001b[0;34m(self, etype, value, tb, tb_offset, number_of_lines_of_context)\u001b[0m\n\u001b[1;32m   1265\u001b[0m         \u001b[0;32mif\u001b[0m \u001b[0mmode\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mverbose_modes\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   1266\u001b[0m             \u001b[0;31m# Verbose modes need a full traceback\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 1267\u001b[0;31m             return VerboseTB.structured_traceback(\n\u001b[0m\u001b[1;32m   1268\u001b[0m                 \u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0metype\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mvalue\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtb\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtb_offset\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mnumber_of_lines_of_context\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   1269\u001b[0m             )\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/IPython/core/ultratb.py\u001b[0m in \u001b[0;36mstructured_traceback\u001b[0;34m(self, etype, evalue, etb, tb_offset, number_of_lines_of_context)\u001b[0m\n\u001b[1;32m   1122\u001b[0m         \u001b[0;34m\"\"\"Return a nice text document describing the traceback.\"\"\"\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   1123\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 1124\u001b[0;31m         formatted_exception = self.format_exception_as_a_whole(etype, evalue, etb, number_of_lines_of_context,\n\u001b[0m\u001b[1;32m   1125\u001b[0m                                                                tb_offset)\n\u001b[1;32m   1126\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/IPython/core/ultratb.py\u001b[0m in \u001b[0;36mformat_exception_as_a_whole\u001b[0;34m(self, etype, evalue, etb, number_of_lines_of_context, tb_offset)\u001b[0m\n\u001b[1;32m   1080\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   1081\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 1082\u001b[0;31m         \u001b[0mlast_unique\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mrecursion_repeat\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfind_recursion\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0morig_etype\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mevalue\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mrecords\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m   1083\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   1084\u001b[0m         \u001b[0mframes\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mformat_records\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mrecords\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mlast_unique\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mrecursion_repeat\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m/home/ubuntu/miniconda3/envs/attribution_py39/lib/python3.9/site-packages/IPython/core/ultratb.py\u001b[0m in \u001b[0;36mfind_recursion\u001b[0;34m(etype, value, records)\u001b[0m\n\u001b[1;32m    380\u001b[0m     \u001b[0;31m# first frame (from in to out) that looks different.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    381\u001b[0m     \u001b[0;32mif\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0mis_recursion_error\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0metype\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mvalue\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mrecords\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 382\u001b[0;31m         \u001b[0;32mreturn\u001b[0m \u001b[0mlen\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mrecords\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m0\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    383\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    384\u001b[0m     \u001b[0;31m# Select filename, lineno, func_name to track frames with\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mTypeError\u001b[0m: object of type 'NoneType' has no len()"
     ]
    }
   ],
   "source": [
    "# loop over the months\n",
    "for mi in months:\n",
    "    \n",
    "    ## FORCED\n",
    "    \n",
    "    # load the model data by month\n",
    "    filepath_in=cmip5_raw_filepath+forced_zarrpath\n",
    "    mondat=util.load_zarr_local(filepath_in,monzarr.format(mi)).load()\n",
    "    \n",
    "    # convert to Celsius\n",
    "    mondat=util.KtoC(mondat)\n",
    "    \n",
    "    # subset the model data over the defined calibration period\n",
    "    mi_xsimcal=mondat.sel(time=slice(calibration_period[0],calibration_period[1]))\n",
    "    \n",
    "    # subselect the observational data (by calendar month) to enable calibration\n",
    "    mi_xobscal=xobscal.isel(time=monidx_calperiod[mi]).load()\n",
    "    \n",
    "    # subset the model data over the defined adjustment period\n",
    "    mi_xsimadj=mondat.sel(time=slice(forced_adj_period[0],forced_adj_period[1]))\n",
    "    del mondat\n",
    "    \n",
    "    # bias adjust this calendar month of forced data over the full adjustment period timeseries\n",
    "    # using the geospatial xarray wrapper defined in cmip5_bias_adjustment.py\n",
    "    mi_xsimadj_BA=ba.xr_bias_adjust_singlemonth(mi_xobscal,mi_xsimcal,mi_xsimadj)\n",
    "    \n",
    "    # for keep the result in memory if it is our example month\n",
    "    if mi==eg_month:\n",
    "        \n",
    "        # store the bias adjusted dataset\n",
    "        eg_storage['forced_xsimadj_BA']=mi_xsimadj_BA\n",
    "        # store the raw dataset\n",
    "        eg_storage['forced_xsimadj_raw']=mi_xsimadj\n",
    "        \n",
    "        # store the calibration distributions\n",
    "        eg_storage['xsimcal']=mi_xsimcal\n",
    "        # store the raw dataset\n",
    "        eg_storage['xobscal']=mi_xobscal\n",
    "        \n",
    "    # save the result for the month to a local zarr folder/path\n",
    "    filepath_out=cmip5_ba_filepath+forced_zarrpath\n",
    "    util.save_zarr_local(mi_xsimadj_BA,filepath_out,zarrname=monzarr.format(mi))\n",
    "    \n",
    "    # remove excess forced model data from memory before natural analysis\n",
    "    del mi_xsimadj, mi_xsimadj_BA\n",
    "    \n",
    "    \n",
    "    ## NATURAL\n",
    "    \n",
    "    # load the Natural model data by month\n",
    "    filepath_in=cmip5_raw_filepath+nat_zarrpath\n",
    "    mondat_nat=util.load_zarr_local(filepath_in,monzarr.format(mi)).load()\n",
    "    \n",
    "    # convert to Celsius\n",
    "    mondat_nat=util.KtoC(mondat_nat)\n",
    "\n",
    "    # subset the natural model data over the defined adjustment period\n",
    "    mi_xsimadj_nat=mondat_nat.sel(time=slice(natural_adj_period[0],natural_adj_period[1]))\n",
    "    del mondat_nat\n",
    "    \n",
    "    # bias adjust this calendar month of natural data over the full adjustment period timeseries\n",
    "    # using the geospatial xarray wrapper defined in cmip5_bias_adjustment.py\n",
    "    mi_xsimadj_nat_BA=ba.xr_bias_adjust_singlemonth(mi_xobscal,mi_xsimcal,mi_xsimadj_nat)\n",
    "    \n",
    "    # for keep the result in memory if it is our example month\n",
    "    if mi==eg_month:\n",
    "        \n",
    "        # store the bias adjusted dataset\n",
    "        eg_storage['natural_xsimadj_BA']=mi_xsimadj_nat_BA\n",
    "        # store the raw dataset\n",
    "        eg_storage['natural_xsimadj_raw']=mi_xsimadj_nat\n",
    "        \n",
    "        ## SAVE EXAMPLE MONTH OF DATA\n",
    "        # define the save name of the example file to save out\n",
    "        eg_month_savename='egdata_BA_'+mname+'_mon{:02d}.pkl'.format(eg_month)\n",
    "        \n",
    "        # save out the data with pickle\n",
    "        with open(eg_savepath+eg_month_savename, 'wb') as handle:\n",
    "            pickle.dump(eg_storage, handle, protocol=pickle.HIGHEST_PROTOCOL)\n",
    "    \n",
    "    # save the result for the month to a local zarr folder/path\n",
    "    filepath_out=cmip5_ba_filepath+nat_zarrpath\n",
    "    util.save_zarr_local(mi_xsimadj_nat_BA,filepath_out,zarrname=monzarr.format(mi))\n",
    "    \n",
    "    # remove excess forced model data from memory before next month of analysis\n",
    "    del mi_xsimadj_nat, mi_xsimadj_nat_BA\n",
    "        \n",
    "    # tell the user where we are in the loop\n",
    "    print(mi)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2390fb85-5c2d-4138-bfaa-98194bb7a7af",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3.9 Attribution",
   "language": "python",
   "name": "attribution_py39"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}