Commit 3da824de authored by Fabian Wachsmann's avatar Fabian Wachsmann
Browse files

Delete use-case_climate-extremes-indices_cdo.ipynb

parent 4e1f1975
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Climate Extremes Indices with CDOs\n",
"### according to the methods recommended by the Expert Team On Climate Change Detection and Indices (ETCCDI)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Some of the climate indices operators implemented in the CDOs (ECAs) are unprecise and have disadvantages when the goal is to produce ETCCDI conform indices.\n",
"Therefore, new cdo operators have been developed allowing the user to calculate these indices without a long and complicated script."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"from IPython.display import HTML\n",
"HTML('<iframe src=\"https://slides.com/wachsylon/cdoetccdi/embed\" width=\"576\" height=\"420\" scrolling=\"no\" frameborder=\"0\" webkitallowfullscreen mozallowfullscreen allowfullscreen></iframe>')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The next paragraph of code contains preparations required for cdo usage in python. Note that in the environment file we include the cdo installation via the dev channel of conda forge while installing the python cdo binder for the installation via pip."
]
},
{
"cell_type": "code",
"execution_count": 32,
"metadata": {},
"outputs": [],
"source": [
"#preparation\n",
"from cdo import Cdo\n",
"cdo = Cdo()\n",
"cdo.debug=True\n",
"#This prohibits that existing files are created a second time\n",
"cdo.forceOutput = False"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 0.1 If you work on mistral"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Set cdo to the binary which allows etccdi indices\n",
"cdoWithETCCDI=\"/work/bm0021/cdo_climateindices/cdo_cei/bin/cdo\"\n",
"cdo.setCdo(cdoWithETCCDI)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Set Input directories\n",
"timeseriesDir=\"/work/bm0021/k204210/INDICES2/timeseries/\"\n",
"pr=timeseriesDir+\"/pr*historical*\"\n",
"prHamburg=\"prHamburg.nc\"\n",
"tasmin=timeseriesDir+\"/tasmin*historical*\"\n",
"tasminHamburg=\"tasminHamburg.nc\"\n",
"!ls /work/bm0021/k204210/INDICES2/timeseries/"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Select a subarea because memory might not be large enough\n",
"#and change units\n",
"cdo.mulc(86400,\n",
" input=\"-sellonlatbox,9,10,53,54 \" + pr,\n",
" output=prHamburg)\n",
"cdo.sellonlatbox(9,10,53,54,\n",
" input=tasmin,\n",
" output=tasminHamburg)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 0.2 If you work locally"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import os\n",
"cwd = os.getcwd()\n",
"timeseriesDir=cwd+\"/hamburg_timeseries/\"\n",
"prHamburg=timeseriesDir+\"pr_1961-1990.nc\"\n",
"tasminHamburg=timeseriesDir+\"tasmin_2071-2100.nc\""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 1. Absolute indices"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Absolute indices are usually easy to compute with basic CDOs. The minimum of daily minimum temperature can be calculated with\n",
"***\n",
"*cdo yearmin tasmin tnn.nc*\n",
"***\n",
"The precipitation related indices however have been changed. The original eca_* operators can be modified to give yearly and monthly output where its default is to give one output over the whole time series. The new etccdi operators automatically produce yearly output."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Highest 5day percipitation sum\n",
"!export CDO_TIMESTAT_DATE=\"last\"\n",
"#$cdo eca_rx5day,50,freq=year -runsum,5\n",
"rx5day=\"rx5day_hamburg.nc\"\n",
"rx5day_values = cdo.etccdi_rx5day(input=\"-runsum,5 \"+prHamburg,\n",
" output=\"rx5day_hamburg.nc\",\n",
" returnCdf=True).variables[\n",
" \"rx5dayETCCDI\"][:]\n",
"rx5day_values = rx5day_values.flatten()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"plt.plot(rx5day_values)\n",
"plt.grid()\n",
"plt.xlabel(\"Year\")\n",
"plt.ylabel(\"Precipitation sum over 5 days [mm]\")\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 2. Threshold exceedances"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Threshold based CEIs, in contrast to percentile based CEIs, use predefined values for temperature and precipitation thresholds.\n",
"Similar to the absolute indices, their output can now be yearly or monthly."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#frost days\n",
"# $cdo eca_fd,freq=year\n",
"fd_values = cdo.etccdi_fd(input=tasminHamburg,\n",
" output=\"fd_hamburg.nc\",\n",
" returnCdf=True).variables[\n",
" \"fdETCCDI\"][:]\n",
"fd_values = fd_values.flatten()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"plt.plot(fd_values)\n",
"plt.grid()\n",
"plt.xlabel(\"Year\")\n",
"plt.ylabel(\"Number of Frost days per year\")\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 3. Create Percentiles"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"For some duration indices, percentiles are required. Two new arguments can be specified: <br>\n",
"***rm*** : The ***read_method*** defines whether the first days and last days of the time series should be taken into account. Right now, it can only be set to \"circular\". <br>\n",
"***pm*** : Since a lot of methods exist to calculate a percentile, CDO will allow to set ***percentileMethod*** in the operator call. ETCCDI recommends a method implemented in the software language R as type8. This is right now the only valid value for the argument.<br><br>\n",
"The running minimum and running maximum are needed as input files for *ydrunpctl* so that the operator can create bins for a histogram. The bins are not used if there is enough memory to save all values. This is only the case, if the environment parameter *CDO_PCTL_NBINS* can be set high enough (depending on the system)."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"windowDays=5\n",
"readMethod=\"circular\"\n",
"percentileMethod=\"rtype8\""
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# $cdo ydrunmin,5,rm=c $tasminMerged $tasminrunmin\n",
"cdo.ydrunmin(windowDays,\"rm=\"+readMethod,\n",
" input=tasminHamburg,\n",
" output=\"tasmin_runmin.nc\")"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# $cdo ydrunmax,5,rm=c $tasminMerged $tasminrunmin\n",
"cdo.ydrunmax(windowDays,\"rm=\"+readMethod,\n",
" input=tasminHamburg,\n",
" output=\"tasmin_runmax.nc\")"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#If you set this environment parameter,\n",
"#histograms are ordered values instead of bins\n",
"!export CDO_PCTL_NBINS=302\n",
"# $cdo subc,273.15 -ydrunpctl,10,5,pm=r8,rm=c $tasminMerged ${tasminrunmin} ${tasminrunmax} ${tn10thresh}\n",
"cdo.ydrunpctl(10,windowDays,\"rm=\"+readMethod,\"pm=\"+percentileMethod,\n",
" input=tasminHamburg+\" tasmin_runmin.nc tasmin_runmax.nc\",\n",
" output=\"tn10thresh.nc\")\n",
"# $cdo subc,273.15 -ydrunpctl,90,5,pm=r8,rm=c $tasminMerged ${tasminrunmin} ${tasminrunmax} ${tn10thresh}\n",
"cdo.ydrunpctl(90,windowDays,\"rm=\"+readMethod,\"pm=\"+percentileMethod,\n",
" input=tasminHamburg+\" tasmin_runmin.nc tasmin_runmax.nc\",\n",
" output=\"tn90thresh.nc\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 4. Duration indices"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Next to the output frequency updates (see Absolute and Thresholds), there have been two major changes in the duration indices:\n",
"1. The days are counted over overlapping years and the final period will get the time stamp of the last contributing day.\n",
"2. The percentile threshold file do not need to have the same number of time steps as the time series file."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Cold spell duration index (cold wave index)\n",
"csdi=\"csdi_hamburg.nc\"\n",
"#!/home/dkrz/k204210/cdo-git/src/cdo eca_cwfi,6,freq=year tasminHamburg.nc tn10thresh.nc csdi_hamburg.nc\n",
"csdiValues = cdo.etccdi_csdi(6,\"freq=year\",\n",
" input=tasminHamburg+\" tn10thresh.nc\",\n",
" output=csdi,\n",
" returnCdf=True).variables['csdiETCCDI'][:]\n",
"csdiValues = csdiValues.flatten()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"plt.plot(csdiValues)\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Consecutive Wet Days\n",
"#Precipitation threshold [mm]:\n",
"pt=1\n",
"#Minimum number of days exceeded for a second variable:\n",
"md=5\n",
"#!cdo eca_cwd,1,5,freq=year prHamburg.nc cwd_hamburg.nc"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"cwd_values = cdo.etccdi_cwd(input=prHamburg,\n",
" output=\"cwdHamburg.nc\",\n",
" returnCdf=True).variables[\n",
" \"cwdETCCDI\"][:] \n",
"cwd_values = cwd_values.flatten()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"scrolled": true
},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"plt.hist(cwd_values,bins= [5.5,6.5,7.5,8.5,\n",
" 9.5,10.5,11.5,12.5,13.5])\n",
"plt.grid()\n",
"plt.xlabel(\"Largest number of consecutive\"\n",
" \"wet days per year\")\n",
"plt.ylabel(\"Frequency\")\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 5. Percentile based indices"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The temperature related percentile based indices (*tx90p, tx10p, tn90p, tn10p*) require a special percentile calculation method for years that lie inside the base period. For that years, bootstrapping must be applied where the base period is modified: The index year is taken from the base period and is replaced by another year. Then, the percentile as well as the index is calculated for the new 30 year base period. This is done 29 times so that each year from the base period will be accounted twice. In the end, the mean of 29 indices is taken.<br>\n",
"Therefore, the operators need input arguments:\n",
"1. The window days\n",
"2. The start year of the bootstrapping interval\n",
"3. The end year of the bootstrapping interval\n",
"4. The output frequency"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"!export CDO_PCTL_NBINS=302\n",
"tn10p_values = cdo.etccdi_tn10p(5,2071,2100,\n",
" input=tasminHamburg+\" tasmin_runmin.nc tasmin_runmax.nc\",\n",
" output=\"tn10p_hamburg.nc\",\n",
" returnCdf=True).variables[\"tn10pETCCDI\"][:]\n",
"tn10p_values = tn10p_values.flatten()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"plt.plot(tn10p_values)\n",
"plt.grid()\n",
"plt.xlabel(\"Year\")\n",
"plt.ylabel(\"Number of days with tmin < tmin90\")\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python [conda env:cdoCEI]",
"language": "python",
"name": "conda-env-cdoCEI-py"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.7.3"
}
},
"nbformat": 4,
"nbformat_minor": 4
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment