Plot_Energy_flux_maps.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d57c59b7",
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import scipy.io as scio\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib as mpl\n",
    "import cartopy.feature as cfeature\n",
    "import cartopy.crs as ccrs\n",
    "import netCDF4 as nc\n",
    "import numpy as np\n",
    "from Load_data import Data_from_nc\n",
    "import xarray as xr\n",
    "import math\n",
    "import gc\n",
    "\n",
    "def get_data_from_mat_for_calcu(file, variable):\n",
    "    var_dict = scio.loadmat(file)\n",
    "    var = var_dict[variable]\n",
    "    var = var[:, 29:]\n",
    "    return var.T\n",
    "\n",
    "def get_data_from_nc(file, variable):\n",
    "    file_obj = nc.Dataset(file)\n",
    "    data = file_obj.variables[variable]\n",
    "    var_data = np.array(data)\n",
    "    var_data = var_data[:, 29:, :]\n",
    "    var_data[var_data > 1000000] = np.nan\n",
    "    var_data = np.squeeze(var_data)\n",
    "    return var_data\n",
    "\n",
    "def plot_fre_single(ax, data_xarray, title, metric):\n",
    "    bwr = mpl.cm.get_cmap('Blues')\n",
    "    vmin=0\n",
    "    vmax=0.05\n",
    "    alpha=0.05\n",
    "    colors = [bwr(0),bwr(0.2),bwr(0.3), bwr(0.4), bwr(0.5), bwr(0.6), bwr(0.7), bwr(0.8),bwr(0.9), bwr(1.0)]\n",
    "    levels = [1, 10, 25, 50, 75, 100, 150, 200, 250]\n",
    "    cmap=mpl.colors.ListedColormap(colors)\n",
    "    im = data_xarray.plot(ax=ax, cmap=cmap,vmin=vmin,vmax=vmax,levels=levels, extend='both', add_colorbar=False, add_labels=False)\n",
    "    \n",
    "    cb = plt.colorbar(im,fraction= 0.02, pad= 0.04, extend='both')\n",
    "    ax.coastlines(color='white', linewidth=0.5)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.set_title(title, loc='right', fontsize=18)\n",
    "    ax.set_title(metric, loc='left', fontsize=18)\n",
    "    \n",
    "def plot_fre_single_custom(ax, data_xarray, title, metric, cmap, levels):\n",
    "    bwr = mpl.cm.get_cmap(cmap)\n",
    "    vmin=0\n",
    "    vmax=0.05\n",
    "    alpha=0.05\n",
    "    colors = [bwr(0),bwr(0.2),bwr(0.3), bwr(0.4), bwr(0.5), bwr(0.6), bwr(0.7), bwr(0.8),bwr(0.9), bwr(1.0)]\n",
    "    levels = levels\n",
    "    cmap=mpl.colors.ListedColormap(colors)\n",
    "    im = data_xarray.plot(ax=ax, cmap=cmap,vmin=vmin,vmax=vmax,levels=levels, extend='both', add_colorbar=False, add_labels=False)\n",
    "    \n",
    "    cb = plt.colorbar(im,fraction= 0.02, pad= 0.04, extend='both')\n",
    "    ax.coastlines(color='white', linewidth=0.5)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.set_title(title, loc='right', fontsize=18)\n",
    "    ax.set_title(metric, loc='left', fontsize=18)\n",
    "    \n",
    "def plot_fre_single_custom_div(ax, data_xarray,data_signal, title, metric, cmap, levels):\n",
    "    bwr = mpl.cm.get_cmap(cmap)\n",
    "    vmin=0\n",
    "    vmax=0.05\n",
    "    alpha=0.05\n",
    "    colors = [bwr(0),bwr(0.2),bwr(0.3), bwr(0.4), 'lightgray', bwr(0.6), bwr(0.7),bwr(0.8), bwr(1.0)]\n",
    "    levels = levels\n",
    "    cmap=mpl.colors.ListedColormap(colors)\n",
    "    im = data_xarray.plot(ax=ax, cmap=cmap,vmin=vmin,vmax=vmax,levels=levels, extend='both', add_colorbar=False, add_labels=False)\n",
    "    \n",
    "    cb = plt.colorbar(im,fraction= 0.3, pad= 0.04, extend='both', orientation='horizontal')\n",
    "    cb.ax.tick_params(labelsize=12)\n",
    "    ax.coastlines(color='black', linewidth=0.5)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.set_title(title, loc='right', fontsize=14)\n",
    "    ax.set_title(metric, loc='left', fontsize=14)\n",
    "    ax.contourf(data_lon, data_lat, data_signal, levels=[0, 3.5], colors= 'lightgray')\n",
    "    \n",
    "def plot_fre_single_custom_div_no_colorbar(ax, data_xarray, title, metric, cmap, levels):\n",
    "    bwr = mpl.cm.get_cmap(cmap)\n",
    "    vmin=0\n",
    "    vmax=0.05\n",
    "    alpha=0.05\n",
    "    colors = [bwr(0),bwr(0.2),bwr(0.3), bwr(0.4), 'white', bwr(0.6), bwr(0.7),bwr(0.8), bwr(1.0)]\n",
    "    levels = levels\n",
    "    cmap=mpl.colors.ListedColormap(colors)\n",
    "    im = data_xarray.plot(ax=ax, cmap=cmap,vmin=vmin,vmax=vmax,levels=levels, extend='both', add_colorbar=False, add_labels=False)\n",
    "    \n",
    "    # cb = plt.colorbar(im,fraction= 0.02, pad= 0.04, extend='both')\n",
    "    ax.coastlines(color='white', linewidth=0.5)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.set_title(title, loc='right', fontsize=18)\n",
    "    ax.set_title(metric, loc='left', fontsize=18)\n",
    "    \n",
    "    \n",
    "def plot_fre_single_1(ax, data_xarray, title, metric):\n",
    "    bwr = mpl.cm.get_cmap('Reds')\n",
    "    vmin=0\n",
    "    vmax=0.05\n",
    "    alpha=0.05\n",
    "    colors = [bwr(0),'b','c', 'g', 'm', 'r', 'y', 'orange', 'pink', 'purple']\n",
    "    levels = [1, 4, 5, 7, 20, 33, 36, 38, 39]\n",
    "    cmap=mpl.colors.ListedColormap(colors)\n",
    "    im = data_xarray.plot(ax=ax, cmap=cmap,vmin=vmin,vmax=vmax,levels=levels, extend='both', add_colorbar=False, add_labels=False)\n",
    "    \n",
    "    cb = plt.colorbar(im,fraction= 0.02, pad= 0.04, extend='both')\n",
    "    ax.coastlines(color='white', linewidth=0.5)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.set_title(title, loc='right', fontsize=18)\n",
    "    ax.set_title(metric, loc='left', fontsize=18)\n",
    "    \n",
    "def plot_fre_single_2(ax, data_xarray, title, metric):\n",
    "    bwr = mpl.cm.get_cmap('Reds')\n",
    "    vmin=0\n",
    "    vmax=0.05\n",
    "    alpha=0.05\n",
    "    colors = [bwr(0),bwr(0.2),bwr(0.3), bwr(0.4), bwr(0.5), bwr(0.6), bwr(0.7), bwr(0.8),bwr(0.9), bwr(1.0)]\n",
    "    levels = [0.01, 0.05, 0.10, 0.15, 0.20, 0.25, 0.30, 0.35, 0.40]\n",
    "    cmap=mpl.colors.ListedColormap(colors)\n",
    "    im = data_xarray.plot(ax=ax, cmap=cmap,vmin=vmin,vmax=vmax,levels=levels, extend='both', add_colorbar=False, add_labels=False)\n",
    "    \n",
    "    cb = plt.colorbar(im,fraction= 0.02, pad= 0.04, extend='both')\n",
    "    ax.coastlines(color='white', linewidth=0.5)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.set_title(title, loc='right', fontsize=18)\n",
    "    ax.set_title(metric, loc='left', fontsize=18)\n",
    "    \n",
    "def plot_land_single(ax, data_xarray, title, metric):\n",
    "    bwr = mpl.cm.get_cmap('Greens')\n",
    "    vmin=0\n",
    "    vmax=0.05\n",
    "    alpha=0.05\n",
    "    colors = ['white',bwr(0.2),bwr(0.3), bwr(0.4), bwr(0.5), bwr(0.6), bwr(0.7), bwr(0.8),bwr(0.9), bwr(1.0)]\n",
    "    levels = [0.01, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8]\n",
    "    cmap=mpl.colors.ListedColormap(colors)\n",
    "    im = data_xarray.plot(ax=ax, cmap=cmap,vmin=vmin,vmax=vmax,levels=levels, extend='both', add_colorbar=False, add_labels=False)\n",
    "    \n",
    "    cb = plt.colorbar(im,fraction= 0.02, pad= 0.04, extend='both')\n",
    "    ax.coastlines(color='dimgray', linewidth=0.5)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.set_title(title, loc='right', fontsize=18)\n",
    "    ax.set_title(metric, loc='left', fontsize=18)\n",
    "\n",
    "def get_data_from_mat(file, variable):\n",
    "    var_dict = scio.loadmat(file)\n",
    "    var = var_dict[variable]\n",
    "    var = var[:, 29:] \n",
    "    return var\n",
    "\n",
    "\n",
    "def get_data_from_mat_1(file, variable):\n",
    "    var_dict = scio.loadmat(file)\n",
    "    var = var_dict[variable]\n",
    "    return var\n",
    "\n",
    "ar6_region = get_data_from_mat_for_calcu('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/ar6_region.mat','ar6_region')\n",
    "\n",
    "irr_diff = get_data_from_mat_for_calcu('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/irr_diff_out.mat', 'irr_diff_out')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3d7cf7f4",
   "metadata": {},
   "outputs": [],
   "source": [
    "def calcu_everything_ens(variable, variable2):\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_CESM2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_CESM2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_CESM2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    All_forcings_CESM2 = IWW_IRR_1985_2014_CESM2 - IWW_IRR_1901_1930_CESM2\n",
    "    Oth_forcings_CESM2 = IWW_NOI_1985_2014_CESM2 - IWW_NOI_1901_1930_CESM2\n",
    "    Irr_forcings_CESM2 = IWW_IRR_1985_2014_CESM2 - IWW_NOI_1985_2014_CESM2\n",
    "    \n",
    "    \n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_gw = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_CESM2_gw = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_CESM2_gw = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_CESM2_gw = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    All_forcings_CESM2_gw = IWW_IRR_1985_2014_CESM2_gw - IWW_IRR_1901_1930_CESM2_gw\n",
    "    Oth_forcings_CESM2_gw = IWW_NOI_1985_2014_CESM2_gw - IWW_NOI_1901_1930_CESM2_gw\n",
    "    Irr_forcings_CESM2_gw = IWW_IRR_1985_2014_CESM2_gw - IWW_NOI_1985_2014_CESM2_gw\n",
    "    \n",
    "    IWW_IRR_1901_1930_NorESM = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_NorESM = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_NorESM = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_NorESM = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    All_forcings_NorESM = IWW_IRR_1985_2014_NorESM - IWW_IRR_1901_1930_NorESM\n",
    "    Oth_forcings_NorESM = IWW_NOI_1985_2014_NorESM - IWW_NOI_1901_1930_NorESM\n",
    "    Irr_forcings_NorESM = IWW_IRR_1985_2014_NorESM - IWW_NOI_1985_2014_NorESM\n",
    "    \n",
    "    IWW_IRR_1901_1930_E3SM = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'IRR_'+ variable + '_timmean_1901_1930.nc_0.9x1.25', variable)\n",
    "    IWW_NOI_1901_1930_E3SM = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'NOI_'+ variable + '_timmean_1901_1930.nc_0.9x1.25', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_E3SM = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'IRR_'+ variable + '_timmean_1985_2014.nc_0.9x1.25', variable)\n",
    "    IWW_NOI_1985_2014_E3SM = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'NOI_'+ variable + '_timmean_1985_2014.nc_0.9x1.25', variable)\n",
    "    \n",
    "    All_forcings_E3SM = IWW_IRR_1985_2014_E3SM - IWW_IRR_1901_1930_E3SM\n",
    "    Oth_forcings_E3SM = IWW_NOI_1985_2014_E3SM - IWW_NOI_1901_1930_E3SM\n",
    "    Irr_forcings_E3SM = IWW_IRR_1985_2014_E3SM - IWW_NOI_1985_2014_E3SM\n",
    "    \n",
    "    IWW_IRR_1901_1930_IPSL = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/IRR01/after_processing//IRR01_1901_2014_'+ variable2 + '_Month.nc_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1901_1930_IPSL = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/NOI01/after_processing//NOI01_1901_2014_'+ variable2 + '_Month.nc_1901_1930_timmean_0.9x1.25', variable2)\n",
    "\n",
    "    IWW_IRR_1985_2014_IPSL = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/IRR01/after_processing//IRR01_1901_2014_'+ variable2 + '_Month.nc_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1985_2014_IPSL = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/NOI01/after_processing//NOI01_1901_2014_'+ variable2 + '_Month.nc_1985_2014_timmean_0.9x1.25', variable2)\n",
    "\n",
    "    All_forcings_IPSL = IWW_IRR_1985_2014_IPSL - IWW_IRR_1901_1930_IPSL\n",
    "    Oth_forcings_IPSL = IWW_NOI_1985_2014_IPSL - IWW_NOI_1901_1930_IPSL\n",
    "    Irr_forcings_IPSL = IWW_IRR_1985_2014_IPSL - IWW_NOI_1985_2014_IPSL\n",
    "\n",
    "#     IWW_IRR_1901_1930_MIROC = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_IRR_0.9x1.25_1901_1930_timmean', variable2)\n",
    "#     IWW_NOI_1901_1930_MIROC = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_NOI_0.9x1.25_1901_1930_timmean', variable2)\n",
    "\n",
    "#     IWW_IRR_1985_2014_MIROC = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_IRR_0.9x1.25_1985_2014_timmean', variable2)\n",
    "#     IWW_NOI_1985_2014_MIROC = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_NOI_0.9x1.25_1985_2014_timmean', variable2)\n",
    "\n",
    "#     All_forcings_MIROC = IWW_IRR_1985_2014_MIROC - IWW_IRR_1901_1930_MIROC\n",
    "#     Oth_forcings_MIROC = IWW_NOI_1985_2014_MIROC - IWW_NOI_1901_1930_MIROC\n",
    "#     Irr_forcings_MIROC = IWW_IRR_1985_2014_MIROC - IWW_NOI_1985_2014_MIROC\n",
    "    \n",
    "    IWW_IRR_1901_1930_CNRM = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_IRR_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1901_1930_CNRM = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_NOI_1901_1930_timmean_0.9x1.25', variable2)\n",
    "\n",
    "    IWW_IRR_1985_2014_CNRM = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_IRR_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1985_2014_CNRM = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_NOI_1985_2014_timmean_0.9x1.25', variable2)\n",
    "\n",
    "    \n",
    "    All_forcings_CNRM = IWW_IRR_1985_2014_CNRM - IWW_IRR_1901_1930_CNRM\n",
    "    Oth_forcings_CNRM = IWW_NOI_1985_2014_CNRM - IWW_NOI_1901_1930_CNRM\n",
    "    Irr_forcings_CNRM = IWW_IRR_1985_2014_CNRM - IWW_NOI_1985_2014_CNRM\n",
    "    \n",
    "    All_forcings_all = (All_forcings_CESM2 + All_forcings_CESM2_gw + All_forcings_NorESM + All_forcings_E3SM + All_forcings_CNRM + All_forcings_IPSL) / 7\n",
    "    Oth_forcings_all = (Oth_forcings_CESM2 + Oth_forcings_CESM2_gw + Oth_forcings_NorESM + Oth_forcings_E3SM + Oth_forcings_CNRM + Oth_forcings_IPSL) / 7\n",
    "    Irr_forcings_all = (Irr_forcings_CESM2 + Irr_forcings_CESM2_gw + Irr_forcings_NorESM + Irr_forcings_E3SM + Irr_forcings_CNRM + Irr_forcings_IPSL) / 7\n",
    "    \n",
    "    All_forcings_pos = np.zeros([163, 288])\n",
    "    All_forcings_neg = np.zeros([163, 288])\n",
    "    All_forcings_pos[All_forcings_CESM2>=0.5] = All_forcings_pos[All_forcings_CESM2>=0.5]+1\n",
    "    All_forcings_neg[All_forcings_CESM2<=-0.5] = All_forcings_neg[All_forcings_CESM2<=-0.5]+1\n",
    "    \n",
    "    All_forcings_pos[All_forcings_NorESM>=0.5] = All_forcings_pos[All_forcings_NorESM>=0.5]+1\n",
    "    All_forcings_neg[All_forcings_NorESM<=-0.5] = All_forcings_neg[All_forcings_NorESM<=-0.5]+1\n",
    "    \n",
    "    All_forcings_pos[All_forcings_CESM2>=0.5] = All_forcings_pos[All_forcings_CESM2>=0.5]+1\n",
    "    All_forcings_neg[All_forcings_CESM2<=-0.5] = All_forcings_neg[All_forcings_CESM2<=-0.5]+1\n",
    "    \n",
    "    All_forcings_pos[All_forcings_CESM2_gw>=0.5] = All_forcings_pos[All_forcings_CESM2_gw>=0.5]+1\n",
    "    All_forcings_neg[All_forcings_CESM2_gw<=-0.5] = All_forcings_neg[All_forcings_CESM2_gw<=-0.5]+1\n",
    "    \n",
    "    All_forcings_pos[All_forcings_IPSL>=0.5] = All_forcings_pos[All_forcings_IPSL>=0.5]+1\n",
    "    All_forcings_neg[All_forcings_IPSL<=-0.5] = All_forcings_neg[All_forcings_IPSL<=-0.5]+1\n",
    "    \n",
    "#     All_forcings_pos[All_forcings_MIROC>=0.5] = All_forcings_pos[All_forcings_MIROC>=0.5]+1\n",
    "#     All_forcings_neg[All_forcings_MIROC<=-0.5] = All_forcings_neg[All_forcings_MIROC<=-0.5]+1\n",
    "    \n",
    "    All_forcings_pos[All_forcings_E3SM>=0.5] = All_forcings_pos[All_forcings_E3SM>=0.5]+1\n",
    "    All_forcings_neg[All_forcings_E3SM<=-0.5] = All_forcings_neg[All_forcings_E3SM<=-0.5]+1\n",
    "    \n",
    "    All_forcings_consis = np.maximum(All_forcings_pos, All_forcings_neg)\n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    Oth_forcings_pos = np.zeros([163, 288])\n",
    "    Oth_forcings_neg = np.zeros([163, 288])\n",
    "    Oth_forcings_pos[Oth_forcings_CESM2>=0.5] = Oth_forcings_pos[Oth_forcings_CESM2>=0.5]+1\n",
    "    Oth_forcings_neg[Oth_forcings_CESM2<=-0.5] = Oth_forcings_neg[Oth_forcings_CESM2<=-0.5]+1\n",
    "    \n",
    "    Oth_forcings_pos[Oth_forcings_NorESM>=0.5] = Oth_forcings_pos[Oth_forcings_NorESM>=0.5]+1\n",
    "    Oth_forcings_neg[Oth_forcings_NorESM<=-0.5] = Oth_forcings_neg[Oth_forcings_NorESM<=-0.5]+1\n",
    "    \n",
    "    Oth_forcings_pos[Oth_forcings_CESM2>=0.5] = Oth_forcings_pos[Oth_forcings_CESM2>=0.5]+1\n",
    "    Oth_forcings_neg[Oth_forcings_CESM2<=-0.5] = Oth_forcings_neg[Oth_forcings_CESM2<=-0.5]+1\n",
    "    \n",
    "    Oth_forcings_pos[Oth_forcings_CESM2_gw>=0.5] = Oth_forcings_pos[Oth_forcings_CESM2_gw>=0.5]+1\n",
    "    Oth_forcings_neg[Oth_forcings_CESM2_gw<=-0.5] = Oth_forcings_neg[Oth_forcings_CESM2_gw<=-0.5]+1\n",
    "    \n",
    "    Oth_forcings_pos[Oth_forcings_IPSL>=0.5] = Oth_forcings_pos[Oth_forcings_IPSL>=0.5]+1\n",
    "    Oth_forcings_neg[Oth_forcings_IPSL<=-0.5] = Oth_forcings_neg[Oth_forcings_IPSL<=-0.5]+1\n",
    "    \n",
    "#     Oth_forcings_pos[Oth_forcings_MIROC>=0.5] = Oth_forcings_pos[Oth_forcings_MIROC>=0.5]+1\n",
    "#     Oth_forcings_neg[Oth_forcings_MIROC<=-0.5] = Oth_forcings_neg[Oth_forcings_MIROC<=-0.5]+1\n",
    "    \n",
    "    Oth_forcings_pos[Oth_forcings_E3SM>=0.5] = Oth_forcings_pos[Oth_forcings_E3SM>=0.5]+1\n",
    "    Oth_forcings_neg[Oth_forcings_E3SM<=-0.5] = Oth_forcings_neg[Oth_forcings_E3SM<=-0.5]+1\n",
    "    \n",
    "    Oth_forcings_consis = np.maximum(Oth_forcings_pos, Oth_forcings_neg)\n",
    "    \n",
    "    \n",
    "    \n",
    "    Irr_forcings_pos = np.zeros([163, 288])\n",
    "    Irr_forcings_neg = np.zeros([163, 288])\n",
    "    Irr_forcings_pos[Irr_forcings_CESM2>=0.5] = Irr_forcings_pos[Irr_forcings_CESM2>=0.5]+1\n",
    "    Irr_forcings_neg[Irr_forcings_CESM2<=-0.5] = Irr_forcings_neg[Irr_forcings_CESM2<=-0.5]+1\n",
    "    \n",
    "    Irr_forcings_pos[Irr_forcings_NorESM>=0.5] = Irr_forcings_pos[Irr_forcings_NorESM>=0.5]+1\n",
    "    Irr_forcings_neg[Irr_forcings_NorESM<=-0.5] = Irr_forcings_neg[Irr_forcings_NorESM<=-0.5]+1\n",
    "    \n",
    "    Irr_forcings_pos[Irr_forcings_CESM2>=0.5] = Irr_forcings_pos[Irr_forcings_CESM2>=0.5]+1\n",
    "    Irr_forcings_neg[Irr_forcings_CESM2<=-0.5] = Irr_forcings_neg[Irr_forcings_CESM2<=-0.5]+1\n",
    "    \n",
    "    Irr_forcings_pos[Irr_forcings_CESM2_gw>=0.5] = Irr_forcings_pos[Irr_forcings_CESM2_gw>=0.5]+1\n",
    "    Irr_forcings_neg[Irr_forcings_CESM2_gw<=-0.5] = Irr_forcings_neg[Irr_forcings_CESM2_gw<=-0.5]+1\n",
    "    \n",
    "    Irr_forcings_pos[Irr_forcings_IPSL>=0.5] = Irr_forcings_pos[Irr_forcings_IPSL>=0.5]+1\n",
    "    Irr_forcings_neg[Irr_forcings_IPSL<=-0.5] = Irr_forcings_neg[Irr_forcings_IPSL<=-0.5]+1\n",
    "    \n",
    "#     Irr_forcings_pos[Irr_forcings_MIROC>=0.5] = Irr_forcings_pos[Irr_forcings_MIROC>=0.5]+1\n",
    "#     Irr_forcings_neg[Irr_forcings_MIROC<=-0.5] = Irr_forcings_neg[Irr_forcings_MIROC<=-0.5]+1\n",
    "    \n",
    "    Irr_forcings_pos[Irr_forcings_E3SM>=0.5] = Irr_forcings_pos[Irr_forcings_E3SM>=0.5]+1\n",
    "    Irr_forcings_neg[Irr_forcings_E3SM<=-0.5] = Irr_forcings_neg[Irr_forcings_E3SM<=-0.5]+1\n",
    "    \n",
    "    Irr_forcings_consis = np.maximum(Irr_forcings_pos, Irr_forcings_neg)\n",
    "    \n",
    "    \n",
    "    return All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4c3b0ade",
   "metadata": {},
   "outputs": [],
   "source": [
    "def calcu_net_radiation_ens():\n",
    "    IWW_IRR_1901_1930_CESM2_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ 'LWup' + '_1901_1930_timmean.nc', 'LWup')\n",
    "    IWW_NOI_1901_1930_CESM2_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ 'LWup' + '_1901_1930_timmean.nc', 'LWup')\n",
    "    IWW_IRR_1985_2014_CESM2_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ 'LWup' + '_1985_2014_timmean.nc', 'LWup')\n",
    "    IWW_NOI_1985_2014_CESM2_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ 'LWup' + '_1985_2014_timmean.nc', 'LWup')\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ 'SWup' + '_1901_1930_timmean.nc', 'SWup')\n",
    "    IWW_NOI_1901_1930_CESM2_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ 'SWup' + '_1901_1930_timmean.nc', 'SWup')\n",
    "    IWW_IRR_1985_2014_CESM2_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ 'SWup' + '_1985_2014_timmean.nc', 'SWup')\n",
    "    IWW_NOI_1985_2014_CESM2_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ 'SWup' + '_1985_2014_timmean.nc', 'SWup')\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ 'LWdown' + '_1901_1930_timmean.nc', 'LWdown')\n",
    "    IWW_NOI_1901_1930_CESM2_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ 'LWdown' + '_1901_1930_timmean.nc', 'LWdown')\n",
    "    IWW_IRR_1985_2014_CESM2_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ 'LWdown' + '_1985_2014_timmean.nc', 'LWdown')\n",
    "    IWW_NOI_1985_2014_CESM2_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ 'LWdown' + '_1985_2014_timmean.nc', 'LWdown')\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ 'SWdown' + '_1901_1930_timmean.nc', 'SWdown')\n",
    "    IWW_NOI_1901_1930_CESM2_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ 'SWdown' + '_1901_1930_timmean.nc', 'SWdown')\n",
    "    IWW_IRR_1985_2014_CESM2_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ 'SWdown' + '_1985_2014_timmean.nc', 'SWdown')\n",
    "    IWW_NOI_1985_2014_CESM2_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ 'SWdown' + '_1985_2014_timmean.nc', 'SWdown')\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2 = IWW_IRR_1901_1930_CESM2_LWdown + IWW_IRR_1901_1930_CESM2_SWdown - IWW_IRR_1901_1930_CESM2_LWup - IWW_IRR_1901_1930_CESM2_SWup\n",
    "    IWW_NOI_1901_1930_CESM2 = IWW_NOI_1901_1930_CESM2_LWdown + IWW_NOI_1901_1930_CESM2_SWdown - IWW_NOI_1901_1930_CESM2_LWup - IWW_NOI_1901_1930_CESM2_SWup\n",
    "    IWW_IRR_1985_2014_CESM2 = IWW_IRR_1985_2014_CESM2_LWdown + IWW_IRR_1985_2014_CESM2_SWdown - IWW_IRR_1985_2014_CESM2_LWup - IWW_IRR_1985_2014_CESM2_SWup\n",
    "    IWW_NOI_1985_2014_CESM2 = IWW_NOI_1985_2014_CESM2_LWdown + IWW_NOI_1985_2014_CESM2_SWdown - IWW_NOI_1985_2014_CESM2_LWup - IWW_NOI_1985_2014_CESM2_SWup\n",
    "    \n",
    "    All_forcings_CESM2 = IWW_IRR_1985_2014_CESM2 - IWW_IRR_1901_1930_CESM2\n",
    "    Oth_forcings_CESM2 = IWW_NOI_1985_2014_CESM2 - IWW_NOI_1901_1930_CESM2\n",
    "    Irr_forcings_CESM2 = IWW_IRR_1985_2014_CESM2 - IWW_NOI_1985_2014_CESM2\n",
    "    \n",
    "    \n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_gw_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ 'LWup' + '_1901_1930_timmean.nc', 'LWup')\n",
    "    IWW_NOI_1901_1930_CESM2_gw_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ 'LWup' + '_1901_1930_timmean.nc', 'LWup')\n",
    "    IWW_IRR_1985_2014_CESM2_gw_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ 'LWup' + '_1985_2014_timmean.nc', 'LWup')\n",
    "    IWW_NOI_1985_2014_CESM2_gw_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ 'LWup' + '_1985_2014_timmean.nc', 'LWup')\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_gw_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ 'SWup' + '_1901_1930_timmean.nc', 'SWup')\n",
    "    IWW_NOI_1901_1930_CESM2_gw_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ 'SWup' + '_1901_1930_timmean.nc', 'SWup')\n",
    "    IWW_IRR_1985_2014_CESM2_gw_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ 'SWup' + '_1985_2014_timmean.nc', 'SWup')\n",
    "    IWW_NOI_1985_2014_CESM2_gw_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ 'SWup' + '_1985_2014_timmean.nc', 'SWup')\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_gw_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ 'LWdown' + '_1901_1930_timmean.nc', 'LWdown')\n",
    "    IWW_NOI_1901_1930_CESM2_gw_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ 'LWdown' + '_1901_1930_timmean.nc', 'LWdown')\n",
    "    IWW_IRR_1985_2014_CESM2_gw_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ 'LWdown' + '_1985_2014_timmean.nc', 'LWdown')\n",
    "    IWW_NOI_1985_2014_CESM2_gw_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ 'LWdown' + '_1985_2014_timmean.nc', 'LWdown')\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_gw_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ 'SWdown' + '_1901_1930_timmean.nc', 'SWdown')\n",
    "    IWW_NOI_1901_1930_CESM2_gw_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ 'SWdown' + '_1901_1930_timmean.nc', 'SWdown')\n",
    "    IWW_IRR_1985_2014_CESM2_gw_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ 'SWdown' + '_1985_2014_timmean.nc', 'SWdown')\n",
    "    IWW_NOI_1985_2014_CESM2_gw_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ 'SWdown' + '_1985_2014_timmean.nc', 'SWdown')\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_gw = IWW_IRR_1901_1930_CESM2_gw_LWdown + IWW_IRR_1901_1930_CESM2_gw_SWdown - IWW_IRR_1901_1930_CESM2_gw_LWup - IWW_IRR_1901_1930_CESM2_gw_SWup\n",
    "    IWW_NOI_1901_1930_CESM2_gw = IWW_NOI_1901_1930_CESM2_gw_LWdown + IWW_NOI_1901_1930_CESM2_gw_SWdown - IWW_NOI_1901_1930_CESM2_gw_LWup - IWW_NOI_1901_1930_CESM2_gw_SWup\n",
    "    IWW_IRR_1985_2014_CESM2_gw = IWW_IRR_1985_2014_CESM2_gw_LWdown + IWW_IRR_1985_2014_CESM2_gw_SWdown - IWW_IRR_1985_2014_CESM2_gw_LWup - IWW_IRR_1985_2014_CESM2_gw_SWup\n",
    "    IWW_NOI_1985_2014_CESM2_gw = IWW_NOI_1985_2014_CESM2_gw_LWdown + IWW_NOI_1985_2014_CESM2_gw_SWdown - IWW_NOI_1985_2014_CESM2_gw_LWup - IWW_NOI_1985_2014_CESM2_gw_SWup\n",
    "    \n",
    "    \n",
    "    All_forcings_CESM2_gw = IWW_IRR_1985_2014_CESM2_gw - IWW_IRR_1901_1930_CESM2_gw\n",
    "    Oth_forcings_CESM2_gw = IWW_NOI_1985_2014_CESM2_gw - IWW_NOI_1901_1930_CESM2_gw\n",
    "    Irr_forcings_CESM2_gw = IWW_IRR_1985_2014_CESM2_gw - IWW_NOI_1985_2014_CESM2_gw\n",
    "    \n",
    "    IWW_IRR_1901_1930_NorESM_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ 'LWup' + '_1901_1930_timmean.nc', 'LWup')\n",
    "    IWW_NOI_1901_1930_NorESM_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ 'LWup' + '_1901_1930_timmean.nc', 'LWup')\n",
    "    IWW_IRR_1985_2014_NorESM_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ 'LWup' + '_1985_2014_timmean.nc', 'LWup')\n",
    "    IWW_NOI_1985_2014_NorESM_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ 'LWup' + '_1985_2014_timmean.nc', 'LWup')\n",
    "    \n",
    "    IWW_IRR_1901_1930_NorESM_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ 'SWup' + '_1901_1930_timmean.nc', 'SWup')\n",
    "    IWW_NOI_1901_1930_NorESM_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ 'SWup' + '_1901_1930_timmean.nc', 'SWup')\n",
    "    IWW_IRR_1985_2014_NorESM_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ 'SWup' + '_1985_2014_timmean.nc', 'SWup')\n",
    "    IWW_NOI_1985_2014_NorESM_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ 'SWup' + '_1985_2014_timmean.nc', 'SWup')\n",
    "    \n",
    "    IWW_IRR_1901_1930_NorESM_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ 'LWdown' + '_1901_1930_timmean.nc', 'LWdown')\n",
    "    IWW_NOI_1901_1930_NorESM_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ 'LWdown' + '_1901_1930_timmean.nc', 'LWdown')\n",
    "    IWW_IRR_1985_2014_NorESM_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ 'LWdown' + '_1985_2014_timmean.nc', 'LWdown')\n",
    "    IWW_NOI_1985_2014_NorESM_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ 'LWdown' + '_1985_2014_timmean.nc', 'LWdown')\n",
    "    \n",
    "    IWW_IRR_1901_1930_NorESM_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ 'SWdown' + '_1901_1930_timmean.nc', 'SWdown')\n",
    "    IWW_NOI_1901_1930_NorESM_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ 'SWdown' + '_1901_1930_timmean.nc', 'SWdown')\n",
    "    IWW_IRR_1985_2014_NorESM_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ 'SWdown' + '_1985_2014_timmean.nc', 'SWdown')\n",
    "    IWW_NOI_1985_2014_NorESM_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ 'SWdown' + '_1985_2014_timmean.nc', 'SWdown')\n",
    "    \n",
    "    IWW_IRR_1901_1930_NorESM = IWW_IRR_1901_1930_NorESM_LWdown + IWW_IRR_1901_1930_NorESM_SWdown - IWW_IRR_1901_1930_NorESM_LWup - IWW_IRR_1901_1930_NorESM_SWup\n",
    "    IWW_NOI_1901_1930_NorESM = IWW_NOI_1901_1930_NorESM_LWdown + IWW_NOI_1901_1930_NorESM_SWdown - IWW_NOI_1901_1930_NorESM_LWup - IWW_NOI_1901_1930_NorESM_SWup\n",
    "    IWW_IRR_1985_2014_NorESM = IWW_IRR_1985_2014_NorESM_LWdown + IWW_IRR_1985_2014_NorESM_SWdown - IWW_IRR_1985_2014_NorESM_LWup - IWW_IRR_1985_2014_NorESM_SWup\n",
    "    IWW_NOI_1985_2014_NorESM = IWW_NOI_1985_2014_NorESM_LWdown + IWW_NOI_1985_2014_NorESM_SWdown - IWW_NOI_1985_2014_NorESM_LWup - IWW_NOI_1985_2014_NorESM_SWup\n",
    "    \n",
    "    All_forcings_NorESM = IWW_IRR_1985_2014_NorESM - IWW_IRR_1901_1930_NorESM\n",
    "    Oth_forcings_NorESM = IWW_NOI_1985_2014_NorESM - IWW_NOI_1901_1930_NorESM\n",
    "    Irr_forcings_NorESM = IWW_IRR_1985_2014_NorESM - IWW_NOI_1985_2014_NorESM\n",
    "    \n",
    "    variable = 'LWup'\n",
    "    IWW_IRR_1901_1930_E3SM_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'IRR_'+ variable + '_timmean_1901_1930.nc_0.9x1.25', variable)\n",
    "    IWW_NOI_1901_1930_E3SM_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'NOI_'+ variable + '_timmean_1901_1930.nc_0.9x1.25', variable)\n",
    "    IWW_IRR_1985_2014_E3SM_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'IRR_'+ variable + '_timmean_1985_2014.nc_0.9x1.25', variable)\n",
    "    IWW_NOI_1985_2014_E3SM_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'NOI_'+ variable + '_timmean_1985_2014.nc_0.9x1.25', variable)\n",
    "    variable = 'SWup'\n",
    "    IWW_IRR_1901_1930_E3SM_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'IRR_'+ variable + '_timmean_1901_1930.nc_0.9x1.25', variable)\n",
    "    IWW_NOI_1901_1930_E3SM_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'NOI_'+ variable + '_timmean_1901_1930.nc_0.9x1.25', variable)\n",
    "    IWW_IRR_1985_2014_E3SM_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'IRR_'+ variable + '_timmean_1985_2014.nc_0.9x1.25', variable)\n",
    "    IWW_NOI_1985_2014_E3SM_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'NOI_'+ variable + '_timmean_1985_2014.nc_0.9x1.25', variable)\n",
    "    variable = 'LWdown'\n",
    "    IWW_IRR_1901_1930_E3SM_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'IRR_'+ variable + '_timmean_1901_1930.nc_0.9x1.25', variable)\n",
    "    IWW_NOI_1901_1930_E3SM_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'NOI_'+ variable + '_timmean_1901_1930.nc_0.9x1.25', variable)\n",
    "    IWW_IRR_1985_2014_E3SM_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'IRR_'+ variable + '_timmean_1985_2014.nc_0.9x1.25', variable)\n",
    "    IWW_NOI_1985_2014_E3SM_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'NOI_'+ variable + '_timmean_1985_2014.nc_0.9x1.25', variable)\n",
    "    variable = 'SWdown'\n",
    "    IWW_IRR_1901_1930_E3SM_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'IRR_'+ variable + '_timmean_1901_1930.nc_0.9x1.25', variable)\n",
    "    IWW_NOI_1901_1930_E3SM_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'NOI_'+ variable + '_timmean_1901_1930.nc_0.9x1.25', variable)\n",
    "    IWW_IRR_1985_2014_E3SM_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'IRR_'+ variable + '_timmean_1985_2014.nc_0.9x1.25', variable)\n",
    "    IWW_NOI_1985_2014_E3SM_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/E3SM_cdo/Water_Energy//'+ 'NOI_'+ variable + '_timmean_1985_2014.nc_0.9x1.25', variable)\n",
    "    \n",
    "    IWW_IRR_1901_1930_E3SM = IWW_IRR_1901_1930_E3SM_LWdown + IWW_IRR_1901_1930_E3SM_SWdown - IWW_IRR_1901_1930_E3SM_LWup - IWW_IRR_1901_1930_E3SM_SWup\n",
    "    IWW_NOI_1901_1930_E3SM = IWW_NOI_1901_1930_E3SM_LWdown + IWW_NOI_1901_1930_E3SM_SWdown - IWW_NOI_1901_1930_E3SM_LWup - IWW_NOI_1901_1930_E3SM_SWup\n",
    "    IWW_IRR_1985_2014_E3SM = IWW_IRR_1985_2014_E3SM_LWdown + IWW_IRR_1985_2014_E3SM_SWdown - IWW_IRR_1985_2014_E3SM_LWup - IWW_IRR_1985_2014_E3SM_SWup\n",
    "    IWW_NOI_1985_2014_E3SM = IWW_NOI_1985_2014_E3SM_LWdown + IWW_NOI_1985_2014_E3SM_SWdown - IWW_NOI_1985_2014_E3SM_LWup - IWW_NOI_1985_2014_E3SM_SWup\n",
    "    \n",
    "    All_forcings_E3SM = IWW_IRR_1985_2014_E3SM - IWW_IRR_1901_1930_E3SM\n",
    "    Oth_forcings_E3SM = IWW_NOI_1985_2014_E3SM - IWW_NOI_1901_1930_E3SM\n",
    "    Irr_forcings_E3SM = IWW_IRR_1985_2014_E3SM - IWW_NOI_1985_2014_E3SM\n",
    "    \n",
    "    variable2 = 'rlus'\n",
    "    IWW_IRR_1901_1930_IPSL_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/IRR01/after_processing//IRR01_1901_2014_'+ variable2 + '_Month.nc_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1901_1930_IPSL_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/NOI01/after_processing//NOI01_1901_2014_'+ variable2 + '_Month.nc_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_IRR_1985_2014_IPSL_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/IRR01/after_processing//IRR01_1901_2014_'+ variable2 + '_Month.nc_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1985_2014_IPSL_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/NOI01/after_processing//NOI01_1901_2014_'+ variable2 + '_Month.nc_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    variable2 = 'rsus'\n",
    "    IWW_IRR_1901_1930_IPSL_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/IRR01/after_processing//IRR01_1901_2014_'+ variable2 + '_Month.nc_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1901_1930_IPSL_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/NOI01/after_processing//NOI01_1901_2014_'+ variable2 + '_Month.nc_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_IRR_1985_2014_IPSL_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/IRR01/after_processing//IRR01_1901_2014_'+ variable2 + '_Month.nc_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1985_2014_IPSL_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/NOI01/after_processing//NOI01_1901_2014_'+ variable2 + '_Month.nc_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    variable2 = 'rlds'\n",
    "    IWW_IRR_1901_1930_IPSL_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/IRR01/after_processing//IRR01_1901_2014_'+ variable2 + '_Month.nc_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1901_1930_IPSL_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/NOI01/after_processing//NOI01_1901_2014_'+ variable2 + '_Month.nc_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_IRR_1985_2014_IPSL_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/IRR01/after_processing//IRR01_1901_2014_'+ variable2 + '_Month.nc_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1985_2014_IPSL_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/NOI01/after_processing//NOI01_1901_2014_'+ variable2 + '_Month.nc_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    variable2 = 'rsds'\n",
    "    IWW_IRR_1901_1930_IPSL_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/IRR01/after_processing//IRR01_1901_2014_'+ variable2 + '_Month.nc_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1901_1930_IPSL_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/NOI01/after_processing//NOI01_1901_2014_'+ variable2 + '_Month.nc_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_IRR_1985_2014_IPSL_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/IRR01/after_processing//IRR01_1901_2014_'+ variable2 + '_Month.nc_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1985_2014_IPSL_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP//IPSL-CM6/NOI01/after_processing//NOI01_1901_2014_'+ variable2 + '_Month.nc_1985_2014_timmean_0.9x1.25', variable2)\n",
    "\n",
    "    IWW_IRR_1901_1930_IPSL = IWW_IRR_1901_1930_IPSL_LWdown + IWW_IRR_1901_1930_IPSL_SWdown - IWW_IRR_1901_1930_IPSL_LWup - IWW_IRR_1901_1930_IPSL_SWup\n",
    "    IWW_NOI_1901_1930_IPSL = IWW_NOI_1901_1930_IPSL_LWdown + IWW_NOI_1901_1930_IPSL_SWdown - IWW_NOI_1901_1930_IPSL_LWup - IWW_NOI_1901_1930_IPSL_SWup\n",
    "    IWW_IRR_1985_2014_IPSL = IWW_IRR_1985_2014_IPSL_LWdown + IWW_IRR_1985_2014_IPSL_SWdown - IWW_IRR_1985_2014_IPSL_LWup - IWW_IRR_1985_2014_IPSL_SWup\n",
    "    IWW_NOI_1985_2014_IPSL = IWW_NOI_1985_2014_IPSL_LWdown + IWW_NOI_1985_2014_IPSL_SWdown - IWW_NOI_1985_2014_IPSL_LWup - IWW_NOI_1985_2014_IPSL_SWup\n",
    "    \n",
    "    All_forcings_IPSL = IWW_IRR_1985_2014_IPSL - IWW_IRR_1901_1930_IPSL\n",
    "    Oth_forcings_IPSL = IWW_NOI_1985_2014_IPSL - IWW_NOI_1901_1930_IPSL\n",
    "    Irr_forcings_IPSL = IWW_IRR_1985_2014_IPSL - IWW_NOI_1985_2014_IPSL\n",
    "\n",
    "#     variable2 = 'rlus'\n",
    "#     IWW_IRR_1901_1930_MIROC_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_IRR_0.9x1.25_1901_1930_timmean', variable2)\n",
    "#     IWW_NOI_1901_1930_MIROC_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_NOI_0.9x1.25_1901_1930_timmean', variable2)\n",
    "#     IWW_IRR_1985_2014_MIROC_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_IRR_0.9x1.25_1985_2014_timmean', variable2)\n",
    "#     IWW_NOI_1985_2014_MIROC_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_NOI_0.9x1.25_1985_2014_timmean', variable2)\n",
    "#     variable2 = 'rsus'\n",
    "#     IWW_IRR_1901_1930_MIROC_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_IRR_0.9x1.25_1901_1930_timmean', variable2)\n",
    "#     IWW_NOI_1901_1930_MIROC_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_NOI_0.9x1.25_1901_1930_timmean', variable2)\n",
    "#     IWW_IRR_1985_2014_MIROC_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_IRR_0.9x1.25_1985_2014_timmean', variable2)\n",
    "#     IWW_NOI_1985_2014_MIROC_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_NOI_0.9x1.25_1985_2014_timmean', variable2)\n",
    "#     variable2 = 'rlds'\n",
    "#     IWW_IRR_1901_1930_MIROC_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_IRR_0.9x1.25_1901_1930_timmean', variable2)\n",
    "#     IWW_NOI_1901_1930_MIROC_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_NOI_0.9x1.25_1901_1930_timmean', variable2)\n",
    "#     IWW_IRR_1985_2014_MIROC_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_IRR_0.9x1.25_1985_2014_timmean', variable2)\n",
    "#     IWW_NOI_1985_2014_MIROC_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_NOI_0.9x1.25_1985_2014_timmean', variable2)\n",
    "#     variable2 = 'rsds'\n",
    "#     IWW_IRR_1901_1930_MIROC_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_IRR_0.9x1.25_1901_1930_timmean', variable2)\n",
    "#     IWW_NOI_1901_1930_MIROC_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_NOI_0.9x1.25_1901_1930_timmean', variable2)\n",
    "#     IWW_IRR_1985_2014_MIROC_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_IRR_0.9x1.25_1985_2014_timmean', variable2)\n",
    "#     IWW_NOI_1985_2014_MIROC_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/MIROC/Water_energy//'+ variable2 + '_NOI_0.9x1.25_1985_2014_timmean', variable2)\n",
    "    \n",
    "#     IWW_IRR_1901_1930_MIROC = IWW_IRR_1901_1930_MIROC_LWdown + IWW_IRR_1901_1930_MIROC_SWdown - IWW_IRR_1901_1930_MIROC_LWup - IWW_IRR_1901_1930_MIROC_SWup\n",
    "#     IWW_NOI_1901_1930_MIROC = IWW_NOI_1901_1930_MIROC_LWdown + IWW_NOI_1901_1930_MIROC_SWdown - IWW_NOI_1901_1930_MIROC_LWup - IWW_NOI_1901_1930_MIROC_SWup\n",
    "#     IWW_IRR_1985_2014_MIROC = IWW_IRR_1985_2014_MIROC_LWdown + IWW_IRR_1985_2014_MIROC_SWdown - IWW_IRR_1985_2014_MIROC_LWup - IWW_IRR_1985_2014_MIROC_SWup\n",
    "#     IWW_NOI_1985_2014_MIROC = IWW_NOI_1985_2014_MIROC_LWdown + IWW_NOI_1985_2014_MIROC_SWdown - IWW_NOI_1985_2014_MIROC_LWup - IWW_NOI_1985_2014_MIROC_SWup\n",
    "    \n",
    "    \n",
    "#     All_forcings_MIROC = IWW_IRR_1985_2014_MIROC - IWW_IRR_1901_1930_MIROC\n",
    "#     Oth_forcings_MIROC = IWW_NOI_1985_2014_MIROC - IWW_NOI_1901_1930_MIROC\n",
    "#     Irr_forcings_MIROC = IWW_IRR_1985_2014_MIROC - IWW_NOI_1985_2014_MIROC\n",
    "    \n",
    "    variable2 = 'rlus'\n",
    "    IWW_IRR_1901_1930_CNRM_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_IRR_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1901_1930_CNRM_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_NOI_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_IRR_1985_2014_CNRM_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_IRR_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1985_2014_CNRM_LWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_NOI_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    variable2 = 'rsus'\n",
    "    IWW_IRR_1901_1930_CNRM_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_IRR_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1901_1930_CNRM_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_NOI_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_IRR_1985_2014_CNRM_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_IRR_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1985_2014_CNRM_SWup = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_NOI_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    variable2 = 'rlds'\n",
    "    IWW_IRR_1901_1930_CNRM_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_IRR_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1901_1930_CNRM_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_NOI_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_IRR_1985_2014_CNRM_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_IRR_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1985_2014_CNRM_LWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_NOI_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    variable2 = 'rsds'\n",
    "    IWW_IRR_1901_1930_CNRM_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_IRR_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1901_1930_CNRM_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_NOI_1901_1930_timmean_0.9x1.25', variable2)\n",
    "    IWW_IRR_1985_2014_CNRM_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_IRR_1985_2014_timmean_0.9x1.25', variable2)\n",
    "    IWW_NOI_1985_2014_CNRM_SWdown = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/CNRM/Water_energy//'+ variable2 + '_NOI_1985_2014_timmean_0.9x1.25', variable2)\n",
    "\n",
    "    IWW_IRR_1901_1930_CNRM = IWW_IRR_1901_1930_CNRM_LWdown + IWW_IRR_1901_1930_CNRM_SWdown - IWW_IRR_1901_1930_CNRM_LWup - IWW_IRR_1901_1930_CNRM_SWup\n",
    "    IWW_NOI_1901_1930_CNRM = IWW_NOI_1901_1930_CNRM_LWdown + IWW_NOI_1901_1930_CNRM_SWdown - IWW_NOI_1901_1930_CNRM_LWup - IWW_NOI_1901_1930_CNRM_SWup\n",
    "    IWW_IRR_1985_2014_CNRM = IWW_IRR_1985_2014_CNRM_LWdown + IWW_IRR_1985_2014_CNRM_SWdown - IWW_IRR_1985_2014_CNRM_LWup - IWW_IRR_1985_2014_CNRM_SWup\n",
    "    IWW_NOI_1985_2014_CNRM = IWW_NOI_1985_2014_CNRM_LWdown + IWW_NOI_1985_2014_CNRM_SWdown - IWW_NOI_1985_2014_CNRM_LWup - IWW_NOI_1985_2014_CNRM_SWup\n",
    "    \n",
    "    All_forcings_CNRM = IWW_IRR_1985_2014_CNRM - IWW_IRR_1901_1930_CNRM\n",
    "    Oth_forcings_CNRM = IWW_NOI_1985_2014_CNRM - IWW_NOI_1901_1930_CNRM\n",
    "    Irr_forcings_CNRM = IWW_IRR_1985_2014_CNRM - IWW_NOI_1985_2014_CNRM\n",
    "    \n",
    "    All_forcings_all = (All_forcings_CESM2 + All_forcings_CESM2_gw + All_forcings_NorESM + All_forcings_E3SM + All_forcings_CNRM + All_forcings_IPSL) / 7\n",
    "    Oth_forcings_all = (Oth_forcings_CESM2 + Oth_forcings_CESM2_gw + Oth_forcings_NorESM + Oth_forcings_E3SM + Oth_forcings_CNRM + Oth_forcings_IPSL) / 7\n",
    "    Irr_forcings_all = (Irr_forcings_CESM2 + Irr_forcings_CESM2_gw + Irr_forcings_NorESM + Irr_forcings_E3SM + Irr_forcings_CNRM + Irr_forcings_IPSL) / 7\n",
    "    \n",
    "    All_forcings_pos = np.zeros([163, 288])\n",
    "    All_forcings_neg = np.zeros([163, 288])\n",
    "    All_forcings_pos[All_forcings_CESM2>=0.5] = All_forcings_pos[All_forcings_CESM2>=0.5]+1\n",
    "    All_forcings_neg[All_forcings_CESM2<=-0.5] = All_forcings_neg[All_forcings_CESM2<=-0.5]+1\n",
    "    \n",
    "    All_forcings_pos[All_forcings_NorESM>=0.5] = All_forcings_pos[All_forcings_NorESM>=0.5]+1\n",
    "    All_forcings_neg[All_forcings_NorESM<=-0.5] = All_forcings_neg[All_forcings_NorESM<=-0.5]+1\n",
    "    \n",
    "    All_forcings_pos[All_forcings_CESM2>=0.5] = All_forcings_pos[All_forcings_CESM2>=0.5]+1\n",
    "    All_forcings_neg[All_forcings_CESM2<=-0.5] = All_forcings_neg[All_forcings_CESM2<=-0.5]+1\n",
    "    \n",
    "    All_forcings_pos[All_forcings_CESM2_gw>=0.5] = All_forcings_pos[All_forcings_CESM2_gw>=0.5]+1\n",
    "    All_forcings_neg[All_forcings_CESM2_gw<=-0.5] = All_forcings_neg[All_forcings_CESM2_gw<=-0.5]+1\n",
    "    \n",
    "    All_forcings_pos[All_forcings_IPSL>=0.5] = All_forcings_pos[All_forcings_IPSL>=0.5]+1\n",
    "    All_forcings_neg[All_forcings_IPSL<=-0.5] = All_forcings_neg[All_forcings_IPSL<=-0.5]+1\n",
    "    \n",
    "#     All_forcings_pos[All_forcings_MIROC>=0.5] = All_forcings_pos[All_forcings_MIROC>=0.5]+1\n",
    "#     All_forcings_neg[All_forcings_MIROC<=-0.5] = All_forcings_neg[All_forcings_MIROC<=-0.5]+1\n",
    "    \n",
    "    All_forcings_pos[All_forcings_E3SM>=0.5] = All_forcings_pos[All_forcings_E3SM>=0.5]+1\n",
    "    All_forcings_neg[All_forcings_E3SM<=-0.5] = All_forcings_neg[All_forcings_E3SM<=-0.5]+1\n",
    "    \n",
    "    All_forcings_consis = np.maximum(All_forcings_pos, All_forcings_neg)\n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    Oth_forcings_pos = np.zeros([163, 288])\n",
    "    Oth_forcings_neg = np.zeros([163, 288])\n",
    "    Oth_forcings_pos[Oth_forcings_CESM2>=0.5] = Oth_forcings_pos[Oth_forcings_CESM2>=0.5]+1\n",
    "    Oth_forcings_neg[Oth_forcings_CESM2<=-0.5] = Oth_forcings_neg[Oth_forcings_CESM2<=-0.5]+1\n",
    "    \n",
    "    Oth_forcings_pos[Oth_forcings_NorESM>=0.5] = Oth_forcings_pos[Oth_forcings_NorESM>=0.5]+1\n",
    "    Oth_forcings_neg[Oth_forcings_NorESM<=-0.5] = Oth_forcings_neg[Oth_forcings_NorESM<=-0.5]+1\n",
    "    \n",
    "    Oth_forcings_pos[Oth_forcings_CESM2>=0.5] = Oth_forcings_pos[Oth_forcings_CESM2>=0.5]+1\n",
    "    Oth_forcings_neg[Oth_forcings_CESM2<=-0.5] = Oth_forcings_neg[Oth_forcings_CESM2<=-0.5]+1\n",
    "    \n",
    "    Oth_forcings_pos[Oth_forcings_CESM2_gw>=0.5] = Oth_forcings_pos[Oth_forcings_CESM2_gw>=0.5]+1\n",
    "    Oth_forcings_neg[Oth_forcings_CESM2_gw<=-0.5] = Oth_forcings_neg[Oth_forcings_CESM2_gw<=-0.5]+1\n",
    "    \n",
    "    Oth_forcings_pos[Oth_forcings_IPSL>=0.5] = Oth_forcings_pos[Oth_forcings_IPSL>=0.5]+1\n",
    "    Oth_forcings_neg[Oth_forcings_IPSL<=-0.5] = Oth_forcings_neg[Oth_forcings_IPSL<=-0.5]+1\n",
    "    \n",
    "#     Oth_forcings_pos[Oth_forcings_MIROC>=0.5] = Oth_forcings_pos[Oth_forcings_MIROC>=0.5]+1\n",
    "#     Oth_forcings_neg[Oth_forcings_MIROC<=-0.5] = Oth_forcings_neg[Oth_forcings_MIROC<=-0.5]+1\n",
    "    \n",
    "    Oth_forcings_pos[Oth_forcings_E3SM>=0.5] = Oth_forcings_pos[Oth_forcings_E3SM>=0.5]+1\n",
    "    Oth_forcings_neg[Oth_forcings_E3SM<=-0.5] = Oth_forcings_neg[Oth_forcings_E3SM<=-0.5]+1\n",
    "    \n",
    "    Oth_forcings_consis = np.maximum(Oth_forcings_pos, Oth_forcings_neg)\n",
    "    \n",
    "    \n",
    "    \n",
    "    Irr_forcings_pos = np.zeros([163, 288])\n",
    "    Irr_forcings_neg = np.zeros([163, 288])\n",
    "    Irr_forcings_pos[Irr_forcings_CESM2>=0.5] = Irr_forcings_pos[Irr_forcings_CESM2>=0.5]+1\n",
    "    Irr_forcings_neg[Irr_forcings_CESM2<=-0.5] = Irr_forcings_neg[Irr_forcings_CESM2<=-0.5]+1\n",
    "    \n",
    "    Irr_forcings_pos[Irr_forcings_NorESM>=0.5] = Irr_forcings_pos[Irr_forcings_NorESM>=0.5]+1\n",
    "    Irr_forcings_neg[Irr_forcings_NorESM<=-0.5] = Irr_forcings_neg[Irr_forcings_NorESM<=-0.5]+1\n",
    "    \n",
    "    Irr_forcings_pos[Irr_forcings_CESM2>=0.5] = Irr_forcings_pos[Irr_forcings_CESM2>=0.5]+1\n",
    "    Irr_forcings_neg[Irr_forcings_CESM2<=-0.5] = Irr_forcings_neg[Irr_forcings_CESM2<=-0.5]+1\n",
    "    \n",
    "    Irr_forcings_pos[Irr_forcings_CESM2_gw>=0.5] = Irr_forcings_pos[Irr_forcings_CESM2_gw>=0.5]+1\n",
    "    Irr_forcings_neg[Irr_forcings_CESM2_gw<=-0.5] = Irr_forcings_neg[Irr_forcings_CESM2_gw<=-0.5]+1\n",
    "    \n",
    "    Irr_forcings_pos[Irr_forcings_IPSL>=0.5] = Irr_forcings_pos[Irr_forcings_IPSL>=0.5]+1\n",
    "    Irr_forcings_neg[Irr_forcings_IPSL<=-0.5] = Irr_forcings_neg[Irr_forcings_IPSL<=-0.5]+1\n",
    "    \n",
    "#     Irr_forcings_pos[Irr_forcings_MIROC>=0.5] = Irr_forcings_pos[Irr_forcings_MIROC>=0.5]+1\n",
    "#     Irr_forcings_neg[Irr_forcings_MIROC<=-0.5] = Irr_forcings_neg[Irr_forcings_MIROC<=-0.5]+1\n",
    "    \n",
    "    Irr_forcings_pos[Irr_forcings_E3SM>=0.5] = Irr_forcings_pos[Irr_forcings_E3SM>=0.5]+1\n",
    "    Irr_forcings_neg[Irr_forcings_E3SM<=-0.5] = Irr_forcings_neg[Irr_forcings_E3SM<=-0.5]+1\n",
    "    \n",
    "    Irr_forcings_consis = np.maximum(Irr_forcings_pos, Irr_forcings_neg)\n",
    "    \n",
    "    return All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9d4eb228",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_everything_ens(index1, index2, index3, variable, variable2, variable_for_show, unit, cmap_seq, cmap_div, levels_seq, levels_div):\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_everything_ens(variable, variable2)\n",
    "    fig = plt.figure(figsize=(6,12),dpi=300)\n",
    "    fig.subplots_adjust(hspace=0.1, wspace=0.01, left = 0.05, right = 0.95, top = 0.95, bottom = 0.05)\n",
    "    \n",
    "    \n",
    "    proj=ccrs.PlateCarree()\n",
    "    \n",
    "    data_xr_1 = xr.DataArray(All_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(311, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, All_forcings_consis, 'All forcings', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(-0.05, 0.92, index1, color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Oth_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(312, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Oth_forcings_consis, 'All forcings - irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(-0.05, 0.92, index2, color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(313, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(-0.05, 0.92, index3, color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "640789a5",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_LHF_SHF_Rnet():\n",
    "    variable = 'EFLX_LH_TOT'\n",
    "    variable2 = 'hfls'\n",
    "    variable_for_show = 'LHF'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'RdBu'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_everything_ens(variable, variable2)\n",
    "    fig = plt.figure(figsize=(18,12),dpi=300)\n",
    "    fig.subplots_adjust(hspace=0.1, wspace=0.01, left = 0.05, right = 0.95, top = 0.95, bottom = 0.05)\n",
    "    \n",
    "    proj=ccrs.PlateCarree()\n",
    "    \n",
    "    data_xr_1 = xr.DataArray(All_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(331, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, All_forcings_consis, 'All forcings', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'a', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Oth_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(334, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Oth_forcings_consis, 'All forcings - irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'd', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(337, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'g', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    variable = 'FSH'\n",
    "    variable2 = 'hfss'\n",
    "    variable_for_show = 'SHF'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'coolwarm'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_everything_ens(variable, variable2)\n",
    "\n",
    "    \n",
    "    data_xr_1 = xr.DataArray(All_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(332, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, All_forcings_consis, 'All forcings', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'b', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Oth_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(335, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Oth_forcings_consis, 'All forcings - irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'e', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(338, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'h', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "\n",
    "    variable_for_show = 'Rnet'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'BrBG'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_net_radiation_ens()\n",
    "\n",
    "    \n",
    "    data_xr_1 = xr.DataArray(All_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(333, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, All_forcings_consis, 'All forcings', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'c', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Oth_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(336, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Oth_forcings_consis, 'All forcings - irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'f', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(339, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'i', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    \n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ef2d8c9e",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_LHF_SHF_Rnet_new():\n",
    "    variable = 'EFLX_LH_TOT'\n",
    "    variable2 = 'hfls'\n",
    "    variable_for_show = 'LHF'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'RdBu'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_everything_ens(variable, variable2)\n",
    "    fig = plt.figure(figsize=(6,12),dpi=300)\n",
    "    fig.subplots_adjust(hspace=0.1, wspace=0.01, left = 0.05, right = 0.95, top = 0.95, bottom = 0.05)\n",
    "    \n",
    "    proj=ccrs.PlateCarree()\n",
    "    \n",
    "\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(311, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'a', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    variable = 'FSH'\n",
    "    variable2 = 'hfss'\n",
    "    variable_for_show = 'SHF'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'coolwarm'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_everything_ens(variable, variable2)\n",
    "\n",
    "    \n",
    "    \n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(312, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'b', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "\n",
    "    variable_for_show = 'Rnet'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'BrBG'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_net_radiation_ens()\n",
    "\n",
    "    \n",
    "    \n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(313, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'c', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    \n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "801eed40",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_LHF_SHF_Rnet_new()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5a03bdb4",
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_LHF_SHF_Rnet()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "13623c55",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_radiation():\n",
    "    variable = 'SWdown'\n",
    "    variable2 = 'rsds'\n",
    "    variable_for_show = 'SWdown'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'bwr'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_everything_ens(variable, variable2)\n",
    "    fig = plt.figure(figsize=(24,12),dpi=300)\n",
    "    fig.subplots_adjust(hspace=0.1, wspace=0.01, left = 0.05, right = 0.95, top = 0.95, bottom = 0.05)\n",
    "    \n",
    "    proj=ccrs.PlateCarree()\n",
    "    \n",
    "    data_xr_1 = xr.DataArray(All_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(341, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, All_forcings_consis, 'All forcings', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'a', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Oth_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(345, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Oth_forcings_consis, 'All forcings - irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'e', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(349, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'i', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    variable = 'SWup'\n",
    "    variable2 = 'rsus'\n",
    "    variable_for_show = 'SWup'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'bwr'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_everything_ens(variable, variable2)\n",
    "    \n",
    "    \n",
    "    data_xr_1 = xr.DataArray(All_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(342, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, All_forcings_consis, 'All forcings', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'b', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Oth_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(346, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Oth_forcings_consis, 'All forcings - irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'f', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(3,4,10, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'j', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    \n",
    "    \n",
    "    \n",
    "    variable = 'LWdown'\n",
    "    variable2 = 'rlds'\n",
    "    variable_for_show = 'LWdown'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'RdBu_r'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_everything_ens(variable, variable2)\n",
    "    \n",
    "    \n",
    "    data_xr_1 = xr.DataArray(All_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(343, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, All_forcings_consis, 'All forcings', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'c', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Oth_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(347, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Oth_forcings_consis, 'All forcings - irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'g', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(3,4,11, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'k', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    \n",
    "    \n",
    "    \n",
    "    variable = 'LWup'\n",
    "    variable2 = 'rlus'\n",
    "    variable_for_show = 'LWup'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'RdBu_r'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_everything_ens(variable, variable2)\n",
    "    \n",
    "    \n",
    "    data_xr_1 = xr.DataArray(All_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(344, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, All_forcings_consis, 'All forcings', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'd', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Oth_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(348, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Oth_forcings_consis, 'All forcings - irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'h', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(3,4,12, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'l', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1ae7a289",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_radiation_new():\n",
    "    variable = 'SWdown'\n",
    "    variable2 = 'rsds'\n",
    "    variable_for_show = 'SWdown'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'bwr'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_everything_ens(variable, variable2)\n",
    "    fig = plt.figure(figsize=(12,8),dpi=300)\n",
    "    fig.subplots_adjust(hspace=0.1, wspace=0.01, left = 0.05, right = 0.95, top = 0.95, bottom = 0.05)\n",
    "    \n",
    "    proj=ccrs.PlateCarree()\n",
    "    \n",
    "    \n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(221, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'a', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    variable = 'SWup'\n",
    "    variable2 = 'rsus'\n",
    "    variable_for_show = 'SWup'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'bwr'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_everything_ens(variable, variable2)\n",
    "    \n",
    "    \n",
    "    \n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(222, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'b', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    \n",
    "    \n",
    "    \n",
    "    variable = 'LWdown'\n",
    "    variable2 = 'rlds'\n",
    "    variable_for_show = 'LWdown'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'RdBu_r'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_everything_ens(variable, variable2)\n",
    "    \n",
    "    \n",
    "    \n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(223, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'c', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    \n",
    "    \n",
    "    \n",
    "    variable = 'LWup'\n",
    "    variable2 = 'rlus'\n",
    "    variable_for_show = 'LWup'\n",
    "    unit = '$\\mathregular{W/m^2}$'\n",
    "    cmap_div = 'RdBu_r'\n",
    "    levels_div = [-10.0, -5.0, -2.0, -0.5, 0.5, 2.0, 5.0, 10.0]\n",
    "    \n",
    "    All_forcings_all, Oth_forcings_all, Irr_forcings_all, All_forcings_consis, Oth_forcings_consis, Irr_forcings_consis = calcu_everything_ens(variable, variable2)\n",
    "    \n",
    "    \n",
    "    \n",
    "\n",
    "    data_xr_1 = xr.DataArray(Irr_forcings_all, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(224, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, Irr_forcings_consis, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(0.01, 0.92, 'd', color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "97e8e8ef",
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_everything_ens_net(index1, index2, index3, variable_for_show, unit, cmap_seq, cmap_div, levels_seq, levels_div):\n",
    "    variable = 'SWdown'\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_1 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_CESM2_1 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_CESM2_1 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_CESM2_1 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    \n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_gw_1 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_CESM2_gw_1 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_CESM2_gw_1 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_CESM2_gw_1 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    \n",
    "    \n",
    "    IWW_IRR_1901_1930_NorESM_1 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_NorESM_1 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_NorESM_1 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_NorESM_1 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    variable = 'LWdown'\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_CESM2_2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_CESM2_2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_CESM2_2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    \n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_gw_2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_CESM2_gw_2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_CESM2_gw_2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_CESM2_gw_2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    \n",
    "    \n",
    "    IWW_IRR_1901_1930_NorESM_2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_NorESM_2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_NorESM_2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_NorESM_2 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    variable = 'SWup'\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_3 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_CESM2_3 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_CESM2_3 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_CESM2_3 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    \n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_gw_3 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_CESM2_gw_3 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_CESM2_gw_3 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_CESM2_gw_3 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    \n",
    "    \n",
    "    IWW_IRR_1901_1930_NorESM_3 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_NorESM_3 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_NorESM_3 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_NorESM_3 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    variable = 'LWup'\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_4 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_CESM2_4 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_CESM2_4 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_CESM2_4 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2' + '/Water_Energy//'+ 'CESM2' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    \n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_gw_4 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_CESM2_gw_4 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_CESM2_gw_4 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_CESM2_gw_4 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'CESM2_gw' + '/Water_Energy//'+ 'CESM2_gw' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    \n",
    "    \n",
    "    IWW_IRR_1901_1930_NorESM_4 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "    IWW_NOI_1901_1930_NorESM_4 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1901_1930_timmean.nc', variable)\n",
    "\n",
    "    IWW_IRR_1985_2014_NorESM_4 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_IRR_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    IWW_NOI_1985_2014_NorESM_4 = get_data_from_nc('/dodrio/scratch/projects/2022_200/project_output/cesm/yi_yao_IRRMIP/after_process//'+ 'NorESM' + '/Water_Energy//'+ 'NorESM' + '_NOI_'+ variable + '_1985_2014_timmean.nc', variable)\n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    fig = plt.figure(figsize=(6,12),dpi=300)\n",
    "    fig.subplots_adjust(hspace=0.1, wspace=0.01, left = 0.05, right = 0.95, top = 0.95, bottom = 0.05)\n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2=IWW_IRR_1901_1930_CESM2_1+IWW_IRR_1901_1930_CESM2_2-IWW_IRR_1901_1930_CESM2_3-IWW_IRR_1901_1930_CESM2_4\n",
    "    IWW_NOI_1901_1930_CESM2=IWW_NOI_1901_1930_CESM2_1+IWW_NOI_1901_1930_CESM2_2-IWW_NOI_1901_1930_CESM2_3-IWW_NOI_1901_1930_CESM2_4\n",
    "    IWW_IRR_1985_2014_CESM2=IWW_IRR_1985_2014_CESM2_1+IWW_IRR_1985_2014_CESM2_2-IWW_IRR_1985_2014_CESM2_3-IWW_IRR_1985_2014_CESM2_4\n",
    "    IWW_NOI_1985_2014_CESM2=IWW_NOI_1985_2014_CESM2_1+IWW_NOI_1985_2014_CESM2_2-IWW_NOI_1985_2014_CESM2_3-IWW_NOI_1985_2014_CESM2_4\n",
    "    \n",
    "    \n",
    "    IWW_IRR_1901_1930_CESM2_gw=IWW_IRR_1901_1930_CESM2_gw_1+IWW_IRR_1901_1930_CESM2_gw_2-IWW_IRR_1901_1930_CESM2_gw_3-IWW_IRR_1901_1930_CESM2_gw_4\n",
    "    IWW_NOI_1901_1930_CESM2_gw=IWW_NOI_1901_1930_CESM2_gw_1+IWW_NOI_1901_1930_CESM2_gw_2-IWW_NOI_1901_1930_CESM2_gw_3-IWW_NOI_1901_1930_CESM2_gw_4\n",
    "    IWW_IRR_1985_2014_CESM2_gw=IWW_IRR_1985_2014_CESM2_gw_1+IWW_IRR_1985_2014_CESM2_gw_2-IWW_IRR_1985_2014_CESM2_gw_3-IWW_IRR_1985_2014_CESM2_gw_4\n",
    "    IWW_NOI_1985_2014_CESM2_gw=IWW_NOI_1985_2014_CESM2_gw_1+IWW_NOI_1985_2014_CESM2_gw_2-IWW_NOI_1985_2014_CESM2_gw_3-IWW_NOI_1985_2014_CESM2_gw_4\n",
    "    \n",
    "    IWW_IRR_1901_1930_NorESM=IWW_IRR_1901_1930_NorESM_1+IWW_IRR_1901_1930_NorESM_2-IWW_IRR_1901_1930_NorESM_3-IWW_IRR_1901_1930_NorESM_4\n",
    "    IWW_NOI_1901_1930_NorESM=IWW_NOI_1901_1930_NorESM_1+IWW_NOI_1901_1930_NorESM_2-IWW_NOI_1901_1930_NorESM_3-IWW_NOI_1901_1930_NorESM_4\n",
    "    IWW_IRR_1985_2014_NorESM=IWW_IRR_1985_2014_NorESM_1+IWW_IRR_1985_2014_NorESM_2-IWW_IRR_1985_2014_NorESM_3-IWW_IRR_1985_2014_NorESM_4\n",
    "    IWW_NOI_1985_2014_NorESM=IWW_NOI_1985_2014_NorESM_1+IWW_NOI_1985_2014_NorESM_2-IWW_NOI_1985_2014_NorESM_3-IWW_NOI_1985_2014_NorESM_4\n",
    "    \n",
    "    \n",
    "    proj=ccrs.PlateCarree()\n",
    "    \n",
    "    data_xr_1 = xr.DataArray((IWW_IRR_1985_2014_CESM2+IWW_IRR_1985_2014_CESM2_gw+IWW_IRR_1985_2014_NorESM-IWW_IRR_1901_1930_CESM2-IWW_IRR_1901_1930_CESM2_gw-IWW_IRR_1901_1930_NorESM)/3, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(311, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, 'All forcings', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(-0.05, 0.92, index1, color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray((IWW_NOI_1985_2014_CESM2+IWW_NOI_1985_2014_CESM2_gw+IWW_NOI_1985_2014_NorESM-IWW_NOI_1901_1930_CESM2-IWW_NOI_1901_1930_CESM2_gw-IWW_NOI_1901_1930_NorESM)/3, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(312, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, 'All forcings - irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(-0.05, 0.92, index2, color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "\n",
    "    data_xr_1 = xr.DataArray((IWW_IRR_1985_2014_CESM2+IWW_IRR_1985_2014_CESM2_gw+IWW_IRR_1985_2014_NorESM-IWW_NOI_1985_2014_CESM2-IWW_NOI_1985_2014_CESM2_gw-IWW_NOI_1985_2014_NorESM)/3, coords={'y': data_lat, 'x': data_lon}, dims=[\"y\", \"x\"])\n",
    "    ax = plt.subplot(313, projection=proj, frameon=False)\n",
    "    plot_fre_single_custom_div(ax, data_xr_1, 'irrigation expansion', variable_for_show + '(' + unit + ')', cmap_div, levels_div)\n",
    "    ax.add_feature(cfeature.OCEAN, color='lightgrey')\n",
    "    ax.text(-0.05, 0.92, index3, color='dimgrey', fontsize=16, transform=ax.transAxes, weight='bold')\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4daee627",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}