{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Atmospheric River Events" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In this notebook, we want to explore [Atmospheric River](https://en.wikipedia.org/wiki/Atmospheric_river) events in the REMO model data. We show how to integrate a 3D REMO field along the vertical axis using the pressure coordinate. We will sum the specific humidity $Q_D$ along the vertical axis weighted by the pressure difference to show atmospheric river events in the model data.\n", "\n", "The integral looks like this: $ Q_{vi} = \\int_{bottom}^{top} Q_d \\,dp$" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's fetch some REMO example output." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Workflow" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import pyremo as pr\n", "\n", "tfile = pr.tutorial.load_dataset(\"remo_output_3d\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We define the formula for the computation of pressure levels here using the vertical hybrid sigma coefficients." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "def pressure(ak, bk, ps):\n", " \"\"\"Returns pressure levels\n", "\n", " Computes pressure levels from hybrid sigma coeffiecients\n", " using the surface pressure ps.\n", "\n", " \"\"\"\n", " return ak + bk * ps" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The hybrid sigma coefficients can be defined at level interfaces or mid levels depending on which locations we need. For this integration, we'll use the pressure differences between level interfaces while the specific humidity is defined on the mid levels. Consequently, we'll need to use `hyai` and `hybi` to compute the pressure field define on the level interfaces." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "aki = tfile.hyai # i: defined at level interfaces\n", "bki = tfile.hybi\n", "akm = tfile.hybm # m: defined at mid levels\n", "bkm = tfile.hybm" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "pi = pressure(aki, bki, tfile.PS) # pressure at level interfaces\n", "pm = pressure(akm, bkm, tfile.PS) # pressure at mid levels" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's check, that the bottom level of `pi` should be equal to the surface pressure." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
<xarray.DataArray ()>\n",
"array(True)array(True)
<xarray.Dataset>\n",
"Dimensions: (rlat: 121, rlon: 129, time: 1)\n",
"Coordinates:\n",
" * time (time) datetime64[ns] 2006-01-01\n",
" * rlon (rlon) float64 -31.73 -31.29 -30.85 ... 23.71 24.15 24.59\n",
" * rlat (rlat) float64 -26.73 -26.29 -25.85 ... 25.19 25.63 26.07\n",
"Data variables:\n",
" QVI (time, rlat, rlon) float64 29.37 29.85 29.65 ... 6.588 6.583\n",
" rotated_pole int32 1array(['2006-01-01T00:00:00.000000000'], dtype='datetime64[ns]')
array([-31.73, -31.29, -30.85, -30.41, -29.97, -29.53, -29.09, -28.65, -28.21,\n",
" -27.77, -27.33, -26.89, -26.45, -26.01, -25.57, -25.13, -24.69, -24.25,\n",
" -23.81, -23.37, -22.93, -22.49, -22.05, -21.61, -21.17, -20.73, -20.29,\n",
" -19.85, -19.41, -18.97, -18.53, -18.09, -17.65, -17.21, -16.77, -16.33,\n",
" -15.89, -15.45, -15.01, -14.57, -14.13, -13.69, -13.25, -12.81, -12.37,\n",
" -11.93, -11.49, -11.05, -10.61, -10.17, -9.73, -9.29, -8.85, -8.41,\n",
" -7.97, -7.53, -7.09, -6.65, -6.21, -5.77, -5.33, -4.89, -4.45,\n",
" -4.01, -3.57, -3.13, -2.69, -2.25, -1.81, -1.37, -0.93, -0.49,\n",
" -0.05, 0.39, 0.83, 1.27, 1.71, 2.15, 2.59, 3.03, 3.47,\n",
" 3.91, 4.35, 4.79, 5.23, 5.67, 6.11, 6.55, 6.99, 7.43,\n",
" 7.87, 8.31, 8.75, 9.19, 9.63, 10.07, 10.51, 10.95, 11.39,\n",
" 11.83, 12.27, 12.71, 13.15, 13.59, 14.03, 14.47, 14.91, 15.35,\n",
" 15.79, 16.23, 16.67, 17.11, 17.55, 17.99, 18.43, 18.87, 19.31,\n",
" 19.75, 20.19, 20.63, 21.07, 21.51, 21.95, 22.39, 22.83, 23.27,\n",
" 23.71, 24.15, 24.59])array([-26.73, -26.29, -25.85, -25.41, -24.97, -24.53, -24.09, -23.65, -23.21,\n",
" -22.77, -22.33, -21.89, -21.45, -21.01, -20.57, -20.13, -19.69, -19.25,\n",
" -18.81, -18.37, -17.93, -17.49, -17.05, -16.61, -16.17, -15.73, -15.29,\n",
" -14.85, -14.41, -13.97, -13.53, -13.09, -12.65, -12.21, -11.77, -11.33,\n",
" -10.89, -10.45, -10.01, -9.57, -9.13, -8.69, -8.25, -7.81, -7.37,\n",
" -6.93, -6.49, -6.05, -5.61, -5.17, -4.73, -4.29, -3.85, -3.41,\n",
" -2.97, -2.53, -2.09, -1.65, -1.21, -0.77, -0.33, 0.11, 0.55,\n",
" 0.99, 1.43, 1.87, 2.31, 2.75, 3.19, 3.63, 4.07, 4.51,\n",
" 4.95, 5.39, 5.83, 6.27, 6.71, 7.15, 7.59, 8.03, 8.47,\n",
" 8.91, 9.35, 9.79, 10.23, 10.67, 11.11, 11.55, 11.99, 12.43,\n",
" 12.87, 13.31, 13.75, 14.19, 14.63, 15.07, 15.51, 15.95, 16.39,\n",
" 16.83, 17.27, 17.71, 18.15, 18.59, 19.03, 19.47, 19.91, 20.35,\n",
" 20.79, 21.23, 21.67, 22.11, 22.55, 22.99, 23.43, 23.87, 24.31,\n",
" 24.75, 25.19, 25.63, 26.07])array([[[29.37055973, 29.84539326, 29.6528811 , ..., 14.46666023,\n",
" 15.50615418, 13.73233009],\n",
" [29.37838781, 29.20898637, 29.01380292, ..., 12.57590328,\n",
" 13.4666962 , 13.69290153],\n",
" [29.53339433, 29.23728129, 28.35914653, ..., 15.62453258,\n",
" 18.20956416, 18.37268447],\n",
" ...,\n",
" [83.93557843, 83.95831155, 84.23090402, ..., 6.69006996,\n",
" 6.61244881, 6.63206705],\n",
" [86.2788733 , 86.28759534, 85.04863121, ..., 6.6652309 ,\n",
" 6.57409011, 6.59765059],\n",
" [86.32654304, 86.32800913, 85.06982852, ..., 6.6745295 ,\n",
" 6.58848365, 6.5826197 ]]])array(1, dtype=int32)