preprocess_netcdf_RNN_public.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import xarray as xr\n",
    "import cdo\n",
    "CDO = cdo.Cdo()\n",
    "import netCDF4\n",
    "import os\n",
    "import datetime"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Authors\n",
    "* Martin Wegmann (martin.wegmann@giub.unibe.ch)\n",
    "* Fernando Jaume Santero (fernando.jaume@unige.ch)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set paths to data folders\n",
    "pathTo20cr= '/Volumes/SPARK/20crv3/' \n",
    "pathToMPIGE_t2m = '/Volumes/SPARK/mpi_ge/t2m/atmos/'\n",
    "pathToEKF = '/Volumes/SPARK/ekf400v2/ensmean/' \n",
    "pathToEKFmember = '/Volumes/SPARK/ekf400v2/member/' \n",
    "pathToLMR = '/Volumes/SPARK/lmr/' \n",
    "pathToCSMLME_t2m = '/Volumes/SPARK/csm_lme/t2m/' "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can download the ensemble mean 20CRv3 data here: https://psl.noaa.gov/data/gridded/data.20thC_ReanV3.html\n",
    "\n",
    "You can download the MPIGE data here: https://mpimet.mpg.de/en/grand-ensemble/\n",
    "\n",
    "You can download the CSMLME data here: https://www.cesm.ucar.edu/projects/community-projects/LME/data-sets.html\n",
    "\n",
    "You can download the EKF400v2 data here: https://cera-www.dkrz.de/WDCC/ui/cerasearch/entry?acronym=EKF400_v2.0\n",
    "\n",
    "You can download the LMR data here: https://atmos.washington.edu/~hakim/lmr/LMRv2/"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Preprocessing 20CRv3 data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.chdir(pathTo20cr)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2m temperature data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# piping in python CDO is faulty and 20cr is weirdly packed, so we do this wrap around. We delete this data later on.\n",
    "CDO.selname(\"air\",input=\"air.2m.mon.mean.nc\",output=\"air.2m.mon.mean_var.nc\")\n",
    "# remap to EKF400V2 resolution\n",
    "CDO.remapbil(\"r192x96\",input=\"air.2m.mon.mean_var.nc\",output=\"air.2m.mon.mean_var_remap.nc\")\n",
    "# remap the spatial dimensions of 20CRv3\n",
    "CDO.sellonlatbox(\"-180,180,-90,90\",input=\"air.2m.mon.mean_var_remap.nc\",output=\"air.2m.mon.mean_data_remap.nc\")\n",
    "CDO.invertlat(input=\"air.2m.mon.mean_data_remap.nc\",output=\"air.2m.mon.mean_data_remap_inv.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# creating the 20CR monthly climatology 1951-1980\n",
    "CDO.ymonmean(input=\"-selyear,1951/1980 air.2m.mon.mean_data_remap_inv.nc\",output=\"air.2m.mon.mean_remap_climate.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# creating the 20CR absolute and anomaly time windows\n",
    "CDO.selyear(\"1851/2015\",input=\"air.2m.mon.mean_data_remap_inv.nc\",output=\"air.2m.mon.mean_18512015_remap.nc\")\n",
    "CDO.selyear(\"1836/1850\",input=\"air.2m.mon.mean_data_remap_inv.nc\",output=\"air.2m.mon.mean_18361850_remap.nc\")\n",
    "CDO.selyear(\"1836/1990\",input=\"air.2m.mon.mean_data_remap_inv.nc\",output=\"air.2m.mon.mean_18361990_remap.nc\")\n",
    "CDO.selyear(\"1991/2015\",input=\"air.2m.mon.mean_data_remap_inv.nc\",output=\"air.2m.mon.mean_19912015_remap.nc\")\n",
    "#\n",
    "CDO.sub(input=\"air.2m.mon.mean_18512015_remap.nc air.2m.mon.mean_remap_climate.nc\", output=\"air.2m.mon.mean_18512015_anoms_remap.nc\")\n",
    "CDO.sub(input=\"air.2m.mon.mean_18361850_remap.nc air.2m.mon.mean_remap_climate.nc\", output=\"air.2m.mon.mean_18361850_anoms_remap.nc\")\n",
    "CDO.sub(input=\"air.2m.mon.mean_18361990_remap.nc air.2m.mon.mean_remap_climate.nc\", output=\"air.2m.mon.mean_18361990_anoms_remap.nc\")\n",
    "CDO.sub(input=\"air.2m.mon.mean_19912015_remap.nc air.2m.mon.mean_remap_climate.nc\", output=\"air.2m.mon.mean_19912015_anoms_remap.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.remove(\"air.2m.mon.mean_var.nc\")\n",
    "os.remove(\"air.2m.mon.mean_data.nc\")\n",
    "os.remove(\"air.2m.mon.mean_var_remap.nc\")\n",
    "os.remove(\"air.2m.mon.mean_data_remap.nc\")\n",
    "os.remove(\"air.2m.mon.mean_data_remap_inv.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create more data to compare with RNN output later on\n",
    "CDO.sellonlatbox(\"-180,180,90,0\",input=\"air.2m.mon.mean_18512015_anoms_remap.nc\",output=\"air.2m.mon.mean_18512015_anoms_remap_NH.nc\")\n",
    "CDO.sellonlatbox(\"-180,180,90,0\",input=\"air.2m.mon.mean_18361850_anoms_remap.nc\",output=\"air.2m.mon.mean_18361850_anoms_remap_NH.nc\")\n",
    "CDO.yearmean(input=\"air.2m.mon.mean_18361850_anoms_remap.nc\",output=\"air.2m.mon.mean_18361850_anoms_remap_ym.nc\")\n",
    "CDO.yearmean(input=\"-selmon,6,7,8 \"+\"air.2m.mon.mean_18361850_anoms_remap.nc\",output=\"air.2m.mon.mean_18361850_anoms_remap_JJA.nc\")\n",
    "CDO.selmon(\"1,2,12\",input=\"air.2m.mon.mean_18361850_anoms_remap.nc\",output=\"temp.nc\")\n",
    "CDO.seltimestep(\"3/44\",input=\"temp.nc\",output=\"tempp.nc\")\n",
    "os.remove(\"temp.nc\")\n",
    "CDO.timselmean(\"3\",input=\"tempp.nc\",output=\"air.2m.mon.mean_18361850_anoms_remap_DJF.nc\")\n",
    "os.remove(\"tempp.nc\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Preprocessing MPIGE data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2m temperature data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.chdir(pathToMPIGE_t2m)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# merging all members and all months into one file\n",
    "CDO.cat(input=\"tas_Amon*.nc\",output=\"tas_Amon_all_datapoints.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# remap to EKF400V2 resolution\n",
    "CDO.remapbil(pathTo20cr+\"air.2m.mon.mean_remap_climate.nc\",input=\"tas_Amon_all_datapoints.nc\",output=\"tas_Amon_all_datapoints_remap.nc\")\n",
    "os.remove(\"tas_Amon_all_datapoints.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create climatology\n",
    "CDO.ymonmean(input=\"-selyear,1951/1980 tas_Amon_all_datapoints_remap.nc\",output=\"tas_Amon_climate.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create anomalies\n",
    "CDO.sub(input=\"tas_Amon_all_datapoints_remap.nc tas_Amon_climate.nc\",output=\"tas_Amon_all_datapoints_anoms.nc\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Preprocessing CSMLME data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2m temperature data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.chdir(pathToCSMLME_t2m)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "arr = os.listdir(pathToCSMLME_t2m)\n",
    "input_paths=[]\n",
    "for file in os.listdir(pathToCSMLME_t2m):\n",
    "    if file.endswith(\".nc\"):\n",
    "        pathToCSMLME_t2m_file=os.path.join(pathToCSMLME_t2m, file)\n",
    "        input_paths.append(pathToCSMLME_t2m_file)\n",
    "len(input_paths)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "amount_files=len(input_paths)\n",
    "for i in range(amount_files):\n",
    "    print(i)\n",
    "    CDO.selname(\"TREFHT\",input=input_paths[i], output=os.path.join(pathToCSMLME_t2m, \"TREFHT_\"+str(i)+\"_data.nc\"))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# merging all members and all months into one file\n",
    "# remove the X for your data\n",
    "CDO.cat(input=\"*_data.nc\",output=\"TREFHT_Amon_all_datapoints.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# the year 850 and 1850 do not have 12 months for some reason\n",
    "# the year 2006 does have one month\n",
    "CDO.selyear(\"851/1849,1851/2005\",input=\"TREFHT_Amon_all_datapoints.nc\",output=\"TREFHT_Amon_all_datapoints_no850.nc\")\n",
    "os.remove(\"TREFHT_Amon_all_datapoints.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# remap to EKF400V2 resolution\n",
    "CDO.remapbil(pathTo20cr+\"air.2m.mon.mean_remap_climate.nc\",input=\"TREFHT_Amon_all_datapoints_no850.nc\",output=\"TREFHT_Amon_all_datapoints_remap.nc\")\n",
    "os.remove(\"TREFHT_Amon_all_datapoints_no850.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# remove leftover files\n",
    "for file in os.listdir(pathToCSMLME_t2m):\n",
    "    if file.endswith(\"_data.nc\"):\n",
    "        os.remove(file)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create climatology\n",
    "CDO.ymonmean(input=\"-selyear,1951/1980 TREFHT_Amon_all_datapoints_remap.nc\",output=\"TREFHT_Amon_climate.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create anomalies\n",
    "CDO.sub(input=\"TREFHT_Amon_all_datapoints_remap.nc TREFHT_Amon_climate.nc\",output=\"TREFHT_Amon_all_datapoints_anoms.nc\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Preprocessing EKF400v2 data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.chdir(pathToEKF)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# if there are too many open files, do the mergetime in batches of 100 years\n",
    "CDO.mergetime(input=\"EKF400_ensmean_????_v2.0.nc\",output=\"EKF400_ensmean_v2.0.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# select 2m temperature\n",
    "CDO.selname(\"air_temperature\",input=\"EKF400_ensmean_v2.0.nc\",output=\"EKF400_ensmean_v2.0_t2m.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# remap to EKF400V2 resolution\n",
    "CDO.remapbil(pathTo20cr+\"air.2m.mon.mean_remap_climate.nc\",input=\"EKF400_ensmean_v2.0_t2m.nc\",output=\"EKF400_ensmean_v2.0_t2m_remap.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create climatology\n",
    "CDO.ymonmean(input=\"-selyear,1951/1980 EKF400_ensmean_v2.0_t2m_remap.nc\",output=\"EKF400_ensmean_v2.0_t2m_climate.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "# create anomalies\n",
    "CDO.sub(input=\"EKF400_ensmean_v2.0_t2m_remap.nc EKF400_ensmean_v2.0_t2m_climate.nc\",output=\"EKF400_ensmean_v2.0_t2m_anoms.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create more data to compare with RNN output later on\n",
    "CDO.selyear(\"1836/1850\",input=\"EKF400_ensmean_v2.0_t2m_anoms.nc\",output=\"EKF400_ensmean_v2.0_t2m_short_anoms.nc\")\n",
    "CDO.sellonlatbox(\"-180,180,90,0\",input=\"EKF400_ensmean_v2.0_t2m_anoms.nc\",output=\"EKF400_ensmean_v2.0_t2m_anoms_NH.nc\")\n",
    "CDO.selyear(\"1836/1850\",input=\"EKF400_ensmean_v2.0_t2m_anoms_NH.nc\",output=\"EKF400_ensmean_v2.0_t2m_short_anoms_NH.nc\")\n",
    "CDO.yearmean(input=\"EKF400_ensmean_v2.0_t2m_anoms.nc\",output=\"EKF400_ensmean_v2.0_t2m_anoms_ym.nc\")\n",
    "CDO.selyear(\"1836/1850\",input=\"EKF400_ensmean_v2.0_t2m_anoms_ym.nc\",output=\"EKF400_ensmean_v2.0_t2m_short_anoms_ym.nc\")\n",
    "CDO.selyear(\"1602/2000\",input=\"EKF400_ensmean_v2.0_t2m_anoms_ym.nc\",output=\"EKF400_ensmean_v2.0_t2m_lmr_anoms_ym.nc\")\n",
    "CDO.yearmean(input=\"-selmon,6,7,8 \"+\"EKF400_ensmean_v2.0_t2m_anoms.nc\",output=\"EKF400_ensmean_v2.0_t2m_anoms_JJA.nc\")\n",
    "CDO.selyear(\"1836/1850\",input=\"EKF400_ensmean_v2.0_t2m_anoms_JJA.nc\",output=\"EKF400_ensmean_v2.0_t2m_short_anoms_JJA.nc\")\n",
    "CDO.selmon(\"1,2,12\",input=\"EKF400_ensmean_v2.0_t2m_anoms.nc\",output=\"temp.nc\")\n",
    "CDO.seltimestep(\"3/1205\",input=\"temp.nc\",output=\"tempp.nc\")\n",
    "os.remove(\"temp.nc\")\n",
    "CDO.timselmean(\"3\",input=\"tempp.nc\",output=\"EKF400_ensmean_v2.0_t2m_anoms_DJF.nc\")\n",
    "CDO.selyear(\"1837/1850\",input=\"EKF400_ensmean_v2.0_t2m_anoms_DJF.nc\",output=\"EKF400_ensmean_v2.0_t2m_short_anoms_DJF.nc\")\n",
    "os.remove(\"tempp.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.remove(\"EKF400_ensmean_v2.0_t2m.nc\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Preprocessing EKF400v2 member data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.chdir(pathToEKFmember)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# lets see what data we have\n",
    "arr = os.listdir(pathToEKFmember)\n",
    "input_paths=[]\n",
    "for file in os.listdir(pathToEKFmember):\n",
    "    if file.endswith(\"t2m.nc\"):\n",
    "        print(os.path.join(pathToEKFmember, file))\n",
    "        pathToEKFmember_file=os.path.join(pathToEKFmember, file)\n",
    "        input_paths.append(pathToEKFmember_file)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# \n",
    "member_index=[f\"{i:01}\" for i in range(1,31)]\n",
    "member_index_length = len(member_index)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# remap to EKF400V2 resolution\n",
    "remapbil_paths=[]\n",
    "for i in range(member_index_length):\n",
    "    output=os.path.join(pathToEKFmember, \"EKF400v2_ens_mem_\" + str(member_index[i])+'_allyears_t2m_remap.nc')\n",
    "    remapbil_paths.append(output)\n",
    "\n",
    "for i in range(member_index_length):\n",
    "    CDO.remapbil(pathTo20cr+\"air.2m.mon.mean_remap_climate.nc\", input=input_paths[i], output=remapbil_paths[i])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create climatology\n",
    "climate_paths=[]\n",
    "for i in range(member_index_length):\n",
    "    output=os.path.join(pathToEKFmember, \"EKF400v2_ens_mem_\" + str(member_index[i])+'_allyears_t2m_climate.nc')\n",
    "    climate_paths.append(output)\n",
    "\n",
    "for i in range(member_index_length): \n",
    "    CDO.ymonmean(input=\"-selyear,1951/1980 \"+ remapbil_paths[i], output=climate_paths[i])    \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create anomalies\n",
    "anom_paths=[]\n",
    "for i in range(member_index_length):\n",
    "    output=os.path.join(pathToEKFmember, \"EKF400v2_ens_mem_\" + str(member_index[i])+'_allyears_t2m_anoms.nc')\n",
    "    anom_paths.append(output)\n",
    "\n",
    "for i in range(member_index_length): \n",
    "    CDO.sub(input=remapbil_paths[i]+\" \"+climate_paths[i],output=anom_paths[i])\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# time windows for later analysis\n",
    "short_anom_paths=[]\n",
    "for i in range(member_index_length):\n",
    "    output=os.path.join(pathToEKFmember, \"EKF400v2_ens_mem_\" + str(member_index[i])+'_allyears_t2m_short_anoms.nc')\n",
    "    short_anom_paths.append(output)\n",
    "\n",
    "ym_anom_paths=[]\n",
    "for i in range(member_index_length):\n",
    "    output=os.path.join(pathToEKFmember, \"EKF400v2_ens_mem_\" + str(member_index[i])+'_allyears_t2m_anoms_ym.nc')\n",
    "    ym_anom_paths.append(output)\n",
    "    \n",
    "lmr_ym_anom_paths=[]\n",
    "for i in range(member_index_length):\n",
    "    output=os.path.join(pathToEKFmember, \"EKF400v2_ens_mem_\" + str(member_index[i])+'_allyears_t2m_lmr_anoms_ym.nc')\n",
    "    lmr_ym_anom_paths.append(output)\n",
    "    \n",
    "JJA_anom_paths=[]\n",
    "for i in range(member_index_length):\n",
    "    output=os.path.join(pathToEKFmember, \"EKF400v2_ens_mem_\" + str(member_index[i])+'_allyears_t2m_anoms_JJA.nc')\n",
    "    JJA_anom_paths.append(output)\n",
    "    \n",
    "short_JJA_anom_paths=[]\n",
    "for i in range(member_index_length):\n",
    "    output=os.path.join(pathToEKFmember, \"EKF400v2_ens_mem_\" + str(member_index[i])+'_allyears_t2m_short_anoms_JJA.nc')\n",
    "    short_JJA_anom_paths.append(output)\n",
    "    \n",
    "short_DJF_anom_paths=[]\n",
    "for i in range(member_index_length):\n",
    "    output=os.path.join(pathToEKFmember, \"EKF400v2_ens_mem_\" + str(member_index[i])+'_allyears_t2m_short_anoms_DJF.nc')\n",
    "    short_DJF_anom_paths.append(output)\n",
    "\n",
    "DJF_anom_paths=[]\n",
    "for i in range(member_index_length):\n",
    "    output=os.path.join(pathToEKFmember, \"EKF400v2_ens_mem_\" + str(member_index[i])+'_allyears_t2m_anoms_DJF.nc')\n",
    "    DJF_anom_paths.append(output)\n",
    "\n",
    "    \n",
    "\n",
    "for i in range(member_index_length): \n",
    "    CDO.selyear(\"1836/1850\",input=anom_paths[i],output=short_anom_paths[i])\n",
    "\n",
    "for i in range(member_index_length): \n",
    "    CDO.yearmean(input=anom_paths[i],output=ym_anom_paths[i])\n",
    "\n",
    "for i in range(member_index_length): \n",
    "    CDO.selyear(\"1602/2000\",input=ym_anom_paths[i],output=lmr_ym_anom_paths[i])\n",
    "\n",
    "for i in range(member_index_length): \n",
    "    CDO.yearmean(input=\"-selmon,6,7,8 \"+anom_paths[i],output=JJA_anom_paths[i])\n",
    "\n",
    "for i in range(member_index_length): \n",
    "    CDO.selyear(\"1836/1850\",input=JJA_anom_paths[i],output=short_JJA_anom_paths[i])\n",
    "\n",
    "for i in range(member_index_length): \n",
    "    CDO.selmon(\"1,2,12\",input=anom_paths[i],output=\"temp.nc\")\n",
    "    CDO.seltimestep(\"3/1205\",input=\"temp.nc\",output=\"tempp.nc\")\n",
    "    os.remove(\"temp.nc\")\n",
    "    CDO.timselmean(\"3\",input=\"tempp.nc\",output=DJF_anom_paths[i])\n",
    "\n",
    "for i in range(member_index_length): \n",
    "    CDO.selyear(\"1837/1850\",input=DJF_anom_paths[i],output=short_DJF_anom_paths[i])\n",
    "\n",
    "    \n",
    "os.remove(\"tempp.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# merging all members and all months into one file\n",
    "CDO.cat(input=\"*_allyears_t2m_remap.nc\",output=\"EKF400v2_allyears_allmembers_t2m_remap.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create climatology\n",
    "CDO.ymonmean(input=\"-selyear,1951/1980 EKF400v2_allyears_allmembers_t2m_remap.nc\",output=\"EKF400v2_allyears_allmembers_t2m_climate.nc\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create anomalies\n",
    "CDO.sub(input=\"EKF400v2_allyears_allmembers_t2m_remap.nc EKF400v2_allyears_allmembers_t2m_climate.nc\",output=\"EKF400v2_allyears_allmembers_t2m_anoms.nc\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Preprocessing Last Millenium Reanalysis data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So the LMR fields are already anomalies with respect to the 1951-1980 climatology"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.chdir(pathToLMR)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "CDO.remapbil(pathTo20cr+\"air.2m.mon.mean_remap_climate.nc\",input=\"air_MCruns_ensemble_mean_LMRv2.1_ensmean.nc\",output=\"air_MCruns_ensemble_mean_LMRv2.1_ensmean_remap.nc\")\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "CDO.selyear(\"1836/1850\",input=\"air_MCruns_ensemble_mean_LMRv2.1_ensmean_remap.nc\",output=\"air_MCruns_ensemble_mean_LMRv2.1_short_anoms_ym.nc\")\n",
    "CDO.selyear(\"1602/2000\",input=\"air_MCruns_ensemble_mean_LMRv2.1_ensmean_remap.nc\",output=\"air_MCruns_ensemble_mean_LMRv2.1_ekf400_anoms_ym.nc\")\n",
    "\n",
    "\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python (temprec)",
   "language": "python",
   "name": "temprec"
  },
  "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.7.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}