{ "cells": [ { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\nANVIL nowcast\n=============\n\nThis example demonstrates how to use ANVIL and the advantages compared to\nextrapolation nowcast and S-PROG.\n\nLoad the libraries.\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "from datetime import datetime, timedelta\nimport warnings\nwarnings.simplefilter(\"ignore\")\nimport matplotlib.pyplot as plt\nimport numpy as np\nfrom pysteps import motion, io, rcparams, utils\nfrom pysteps.nowcasts import anvil, extrapolation, sprog\nfrom pysteps.visualization import plot_precip_field" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Read the input data\n-------------------\n\nANVIL was originally developed to use vertically integrated liquid (VIL) as\nthe input data, but the model allows using any two-dimensional input fields.\nHere we use a composite of rain rates.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "date = datetime.strptime(\"201505151620\", \"%Y%m%d%H%M\")\n\n# Read the data source information from rcparams\ndata_source = rcparams.data_sources[\"mch\"]\n\nroot_path = data_source[\"root_path\"]\npath_fmt = data_source[\"path_fmt\"]\nfn_pattern = data_source[\"fn_pattern\"]\nfn_ext = data_source[\"fn_ext\"]\nimporter_name = data_source[\"importer\"]\nimporter_kwargs = data_source[\"importer_kwargs\"]\n\n# Find the input files in the archive. Use history length of 5 timesteps\nfilenames = io.archive.find_by_date(date, root_path, path_fmt, fn_pattern,\n fn_ext, timestep=5, num_prev_files=5)\n\n# Read the input time series\nimporter = io.get_method(importer_name, \"importer\")\nrainrate_field, quality, metadata = io.read_timeseries(filenames, importer,\n **importer_kwargs)\n\n# Convert to rain rate (mm/h)\nrainrate_field, metadata = utils.to_rainrate(rainrate_field, metadata)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Compute the advection field\n---------------------------\n\nApply the Lucas-Kanade method with the parameters given in Pulkkinen et al.\n(2020) to compute the advection field.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "fd_kwargs = {}\nfd_kwargs[\"max_corners\"] = 1000\nfd_kwargs[\"quality_level\"] = 0.01\nfd_kwargs[\"min_distance\"] = 2\nfd_kwargs[\"block_size\"] = 8\n\nlk_kwargs = {}\nlk_kwargs[\"winsize\"] = (15, 15)\n\noflow_kwargs = {}\noflow_kwargs[\"fd_kwargs\"] = fd_kwargs\noflow_kwargs[\"lk_kwargs\"] = lk_kwargs\noflow_kwargs[\"decl_scale\"] = 10\n\noflow = motion.get_method(\"lucaskanade\")\n\n# transform the input data to logarithmic scale\nrainrate_field_log, _ = utils.transformation.dB_transform(rainrate_field,\n metadata=metadata)\nvelocity = oflow(rainrate_field_log, **oflow_kwargs)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Compute the nowcasts and threshold rain rates below 0.5 mm/h\n~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "forecast_extrap = extrapolation.forecast(rainrate_field[-1], velocity, 3,\n extrap_kwargs={\"allow_nonfinite_values\": True})\nforecast_extrap[forecast_extrap < 0.5] = 0.0\n\nrainrate_field_finite = rainrate_field.copy()\nrainrate_field_finite[~np.isfinite(rainrate_field_finite)] = 0.0\nforecast_sprog = sprog.forecast(rainrate_field_finite[-3:], velocity, 3,\n n_cascade_levels=8, R_thr=0.5)\nforecast_sprog[~np.isfinite(forecast_extrap)] = np.nan\nforecast_sprog[forecast_sprog < 0.5] = 0.0\n\nforecast_anvil = anvil.forecast(rainrate_field[-4:], velocity, 3,\n ar_window_radius=25, ar_order=2)\nforecast_anvil[forecast_anvil < 0.5] = 0.0" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Read the reference observation field and threshold rain rates below 0.5 mm/h\n~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "filenames = io.archive.find_by_date(date, root_path, path_fmt, fn_pattern,\n fn_ext, timestep=5, num_next_files=3)\n\nrefobs_field, quality, metadata = io.read_timeseries(filenames, importer,\n **importer_kwargs)\n\nrefobs_field, metadata = utils.to_rainrate(refobs_field[-1], metadata)\nrefobs_field[refobs_field < 0.5] = 0.0" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Plot the extrapolation, S-PROG and ANVIL nowcasts.\n--------------------------------------------------\n\nFor comparison, the observed rain rate fields are also plotted. Growth and\ndecay areas are marked with red and blue circles, respectively.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "def plot_growth_decay_circles(ax):\n circle = plt.Circle((360, 300), 25, color=\"b\", clip_on=False, fill=False,\n zorder=1e9)\n ax.add_artist(circle)\n circle = plt.Circle((420, 350), 30, color=\"b\", clip_on=False, fill=False,\n zorder=1e9)\n ax.add_artist(circle)\n circle = plt.Circle((405, 380), 30, color=\"b\", clip_on=False, fill=False,\n zorder=1e9)\n ax.add_artist(circle)\n circle = plt.Circle((420, 500), 25, color=\"b\", clip_on=False, fill=False,\n zorder=1e9)\n ax.add_artist(circle)\n circle = plt.Circle((480, 535), 30, color=\"b\", clip_on=False, fill=False,\n zorder=1e9)\n ax.add_artist(circle)\n circle = plt.Circle((330, 470), 35, color=\"b\", clip_on=False, fill=False,\n zorder=1e9)\n ax.add_artist(circle)\n circle = plt.Circle((505, 205), 30, color=\"b\", clip_on=False, fill=False,\n zorder=1e9)\n ax.add_artist(circle)\n circle = plt.Circle((440, 180), 30, color=\"r\", clip_on=False, fill=False,\n zorder=1e9)\n ax.add_artist(circle)\n circle = plt.Circle((590, 240), 30, color=\"r\", clip_on=False, fill=False,\n zorder=1e9)\n ax.add_artist(circle)\n circle = plt.Circle((585, 160), 15, color=\"r\", clip_on=False, fill=False,\n zorder=1e9)\n ax.add_artist(circle)\n\nfig = plt.figure(figsize=(10, 13))\n\nax = fig.add_subplot(321)\nrainrate_field[-1][rainrate_field[-1] < 0.5] = 0.0\nplot_precip_field(rainrate_field[-1])\nplot_growth_decay_circles(ax)\nax.set_title(\"Obs. %s\" % str(date))\n\nax = fig.add_subplot(322)\nplot_precip_field(refobs_field)\nplot_growth_decay_circles(ax)\nax.set_title(\"Obs. %s\" % str(date + timedelta(minutes=15)))\n\nax = fig.add_subplot(323)\nplot_precip_field(forecast_extrap[-1])\nplot_growth_decay_circles(ax)\nax.set_title(\"Extrapolation +15 minutes\")\n\nax = fig.add_subplot(324)\nplot_precip_field(forecast_sprog[-1])\nplot_growth_decay_circles(ax)\nax.set_title(\"S-PROG (with post-processing)\\n +15 minutes\")\n\nax = fig.add_subplot(325)\nplot_precip_field(forecast_anvil[-1])\nplot_growth_decay_circles(ax)\nax.set_title(\"ANVIL +15 minutes\")\n\nplt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Remarks\n-------\n\nThe extrapolation nowcast is static, i.e. it does not predict any growth or\ndecay. While S-PROG is to some extent able to predict growth and decay, this\nthis comes with loss of small-scale features. In addition, statistical\npost-processing needs to be applied to correct the bias and incorrect wet-area\nratio introduced by the autoregressive process. ANVIL is able to do both:\npredict growth and decay and preserve the small-scale structure in a way that\npost-processing is not necessary.\n" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "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.8.0" } }, "nbformat": 4, "nbformat_minor": 0 }