{ "cells": [ { "cell_type": "code", "execution_count": 54, "id": "29027796-dac6-4d58-a445-c1ad0b8e7325", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "/nird/datalake/NS9039K/www/ingo/people/meng/inputdata/v4f_volcanic_v5_mem01.nc\n", "/nird/datalake/NS9039K/www/ingo/people/meng/inputdata/v4f_volcanic_v6_mem01.nc\n", "volcanoes4future_rcp45Volc/volcanoes4future_rcp45Volc_61/volcanoes4future_rcp45Volc_61.cam2.h0.AODVIS-PRECT-TREFHT.gmean.200501-209912.nc\n", "['AODVIS', 'K', 'k']\n", "['PRECT', 'm', 'b']\n", "['TREFHT', 'm/s', 'r']\n" ] } ], "source": [ "####################\n", "# import libraries #\n", "####################\n", "import sys ; sys.path.append('../../jupyter/Modules')\n", "import numpy as np \n", "import matplotlib.pyplot as plt\n", "from scipy import interpolate\n", "from netCDF4 import Dataset\n", "import os \n", "\n", "settings = {\n", " 'experiment': 'volcanoes4future_rcp45Volc',\n", " 'member': 61,\n", " 'fields': [['AODVIS','K','k'],\n", " ['PRECT','m','b'],\n", " ['TREFHT','m/s','r']\n", " ],\n", " 'year1': 2005,\n", " 'yearn': 2099\n", "}\n", "s = settings\n", "nfield = len(s['fields'])\n", "\n", "# read volcanic forcing \n", "def get_forcing(version,forcing_member):\n", " fpath = f'/nird/datalake/NS9039K/www/ingo/people/meng/inputdata/v4f_volcanic_{version}_mem{forcing_member:0>2d}.nc'\n", " print(fpath)\n", " nc = Dataset(fpath,'r')\n", " lat = nc['lat'][:]\n", " area = np.zeros(lat.shape) \n", " for ilat in range(len(lat)):\n", " if ilat == 0:\n", " lats = -90 \n", " latn = 0.5*(lat[ilat]+lat[ilat+1])\n", " elif ilat == len(lat)-1:\n", " lats = 0.5*(lat[ilat]+lat[ilat-1])\n", " latn = 90\n", " else:\n", " lats = 0.5*(lat[ilat]+lat[ilat-1])\n", " latn = 0.5*(lat[ilat]+lat[ilat+1])\n", " area[ilat] = 2*np.pi*6371000**2*(np.sin(latn/180*np.pi)-np.sin(lats/180*np.pi))\n", " colmass = nc['colmass'][1:(2099-2005+1)*12+1,:]\n", " data = np.sum(colmass*area,axis=1) \n", " data = data.reshape([int(len(data)/12),12]).mean(axis=1)\n", " nc.close()\n", " return data\n", "colmass_v5 = get_forcing('v5',(s['member']-1)%60 + 1) \n", "colmass_v6 = get_forcing('v6',(s['member']-1)%60 + 1) \n", "\n", "# open file with monthly global mean output \n", "tag = [] \n", "for ifield in range(nfield):\n", " if ifield == 0:\n", " tag = s['fields'][ifield][0]\n", " else:\n", " tag = tag + '-' + s['fields'][ifield][0] \n", "ifile = f'{s[\"experiment\"]}/{s[\"experiment\"]}_{s[\"member\"]}/{s[\"experiment\"]}_{s[\"member\"]}.cam2.h0.{tag}.gmean.{s[\"year1\"]}01-{s[\"yearn\"]}12.nc'\n", "print(ifile)\n", "nc = Dataset(ifile,'r')\n", "\n", "# prepare plot\n", "SMALL_SIZE = 18\n", "MEDIUM_SIZE = 20\n", "BIGGER_SIZE = 20\n", "plt.rc('font', size=SMALL_SIZE) # controls default text sizes\n", "plt.rc('axes', titlesize=SMALL_SIZE) # fontsize of the axes title\n", "plt.rc('axes', labelsize=SMALL_SIZE) # fontsize of the x and y labels\n", "plt.rc('xtick', labelsize=SMALL_SIZE) # fontsize of the tick labels\n", "plt.rc('ytick', labelsize=SMALL_SIZE) # fontsize of the tick labels\n", "plt.rc('legend', fontsize=SMALL_SIZE) # legend fontsize\n", "plt.rc('figure', titlesize=BIGGER_SIZE) # fontsize of the figure title\n", "figsize=(16,3*(nfield+2))\n", "figres=100\n", "fig,axes = plt.subplots(nrows=(nfield+2),ncols=1, sharex=True, sharey=False, figsize=figsize,dpi=figres)\n", "\n", "# plot forcing versions\n", "tvec = np.arange(s['year1'],s['yearn']+1) \n", "axes[0].plot(tvec,colmass_v5,color='k',label=f'mass forcing of member {(s[\"member\"]-1)%60 + 1}')\n", "axes[1].plot(tvec,colmass_v6,color='k',label=f'mass forcing of member {s[\"member\"]}')\n", "\n", "# loop over output fields\n", "for ifield in range(nfield): \n", " print(s['fields'][ifield])\n", "\n", " # read monthly data \n", " data = nc[s['fields'][ifield][0]][:]\n", "\n", " # convert to annual means \n", " data = data.reshape([int(len(data)/12),12]).mean(axis=1)\n", "\n", " # plot annual time series\n", " axes[ifield+2].plot(tvec,data,color=s['fields'][ifield][2],label=f'{s[\"fields\"][ifield][0]} of member {s[\"member\"]}')\n", "\n", "\n", "# add legends \n", "for ipanel in range(nfield+2):\n", " axes[ipanel].legend(frameon=False,loc=\"right\")\n", "\n", "# save plot\n", "figpath = f'figs/tseries_{s[\"experiment\"]}_{s[\"member\"]}.{tag}.gmean.{s[\"year1\"]}-{s[\"yearn\"]}.png'\n", "plt.savefig(figpath,format='png',dpi=figres,bbox_inches='tight')\n", "\n", "# close input file\n", "nc.close()\n" ] }, { "cell_type": "code", "execution_count": 38, "id": "9128d1e9-ddf5-4506-9cf5-66cdff779275", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "masked_array(data=[0., 0., 0., ..., 0., 0., 0.],\n", " mask=False,\n", " fill_value=1e+20)" ] }, "execution_count": 38, "metadata": {}, "output_type": "execute_result" } ], "source": [ "colmass_v4\n" ] }, { "cell_type": "code", "execution_count": 27, "id": "cbb13b09-42a5-45b9-b4f5-c2f5d07a464f", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "2" ] }, "execution_count": 27, "metadata": {}, "output_type": "execute_result" } ], "source": [ "62%60" ] }, { "cell_type": "code", "execution_count": null, "id": "9d287aad-52de-4f26-9cb1-38371f5a268d", "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.11.15" } }, "nbformat": 4, "nbformat_minor": 5 }