{ "cells": [ { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n# PQPF (Hong Kong)\nThis example demonstrates how to perform\noperational PQPF up to three hours from\nraingauge and radar data, using data from Hong Kong.\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Definitions\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import os\nimport numpy as np\nimport pandas as pd\nfrom pyresample import utils\nimport xarray\nimport cartopy.feature as cfeature\nimport cartopy.crs as ccrs\nimport cartopy.io.shapereader as shpreader\nimport matplotlib\nimport imageio\nfrom collections import OrderedDict\nimport matplotlib.pyplot as plt\nfrom matplotlib.colors import BoundaryNorm, LinearSegmentedColormap\n\nfrom swirlspy.core.resample import grid_resample\nfrom swirlspy.qpe.utils import timestamps_ending, locate_file\nfrom swirlspy.rad.iris import read_iris_grid\nfrom swirlspy.obs import Rain\nfrom swirlspy.qpf import rover\nfrom swirlspy.qpf import sla\nfrom swirlspy.utils import standardize_attr, FrameType\nfrom swirlspy.utils.conversion import to_rainfall_rate, to_rainfall_depth, acc_rainfall_depth, temporal_interpolate\n\nplt.switch_backend('agg')\nTHIS_DIR = os.getcwd()\nos.chdir(THIS_DIR)\n\nstart_time = pd.Timestamp.now()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Initialising\n\nThis section demonstrates extracting\nradar reflectivity data.\n\nStep 1: Define the basetime.\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# Define basetime\nbasetime = pd.Timestamp('201902190800')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Step 2: Using basetime, generate timestamps of desired radar files\ntimestamps_ending() and locate files using locate_file().\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# Obtain radar files\ndir = THIS_DIR + '/../tests/samples/iris/ppi'\nlocated_files = []\nradar_ts = timestamps_ending(\n basetime,\n duration=pd.Timedelta(60, 'm')\n)\n\nfor timestamp in radar_ts:\n located_files.append(locate_file(dir, timestamp))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Step 3: Read data from radar files into xarray.DataArray\nusing read_iris_grid().\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "reflectivity_list = [] # stores reflec from read_iris_grid()\nfor filename in located_files:\n reflec = read_iris_grid(\n filename\n )\n reflectivity_list.append(reflec)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Step 4: Define the target grid as a pyresample AreaDefinition.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# Defining target grid\narea_id = \"hk1980_250km\"\ndescription = (\"A 1km resolution rectangular grid \"\n \"centred at HKO and extending to 250 km \"\n \"in each direction in HK1980 easting/northing coordinates\")\nproj_id = 'hk1980'\nprojection = ('+proj=tmerc +lat_0=22.31213333333334 '\n '+lon_0=114.1785555555556 +k=1 +x_0=836694.05 '\n '+y_0=819069.8 +ellps=intl +towgs84=-162.619,-276.959,'\n '-161.764,0.067753,-2.24365,-1.15883,-1.09425 +units=m '\n '+no_defs')\nx_size = 500\ny_size = 500\narea_extent = (587000, 569000, 1087000, 1069000)\narea_def_tgt = utils.get_area_def(\n area_id, description, proj_id, projection, x_size, y_size, area_extent\n)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Step 5: Reproject the radar data from read_iris_grid() from Centered\nAzimuthal (source) projection to HK 1980 (target) projection.\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# Extracting the AreaDefinition of the source projection\narea_def_src = reflectivity_list[0].attrs['area_def']\n\n# Reprojecting\nreproj_reflectivity_list = []\nfor reflec in reflectivity_list:\n reproj_reflec = grid_resample(\n reflec, area_def_src, area_def_tgt,\n coord_label=['easting', 'northing']\n )\n reproj_reflectivity_list.append(reproj_reflec)\n\ninitialising_time = pd.Timestamp.now()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Step 6: Concatenate the reflectivity xarrays in the list along the time\ndimension.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# Concatenating all the (observed) reflectivity\n# xarray.DataArrays in reproj_reflectivity_list\nreflectivityObs = xarray.concat(reproj_reflectivity_list, dim='time')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Running ROVER and Semi-Lagrangian Advection\n\nThis section demonstrates how to run ROVER and\nSemi-Lagrangian Advection for multiple members.\n\nFirst, write a function to run these steps:\n\n#. Selecting the two relevant xarrays to be used as the ROVER input.\n#. Generate motion field using ROVER, with the selected xarray as the input.\n#. Perform Semi-Lagrangian Advection using the motion fields from ROVER.\n\nWriting the function makes it easier to implement multiprocessing if\nrequired.\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "def ensemble_rover_sla(\n reflectivityObs,\n basetime,\n interval,\n duration,\n member,\n start_level,\n max_level,\n rho, alpha, sigma,\n track_interval\n):\n \"\"\"\n A function to run ROVER and SLA.\n\n Parameters\n -------------\n reflectivityObs: xarray.DataArray\n The xarray containing observed reflectivity.\n basetime: pandas.Timestamp\n The basetime of the forecast.\n Equivalent to the end time of the latest radar scan.\n interval: pandas.Timedelta\n The interval between the radar scans.\n duration: pandas.Timedelta\n The duration of the forecast.\n member: str\n The name of the member.\n start_level, max_level, rho, alpha, sigma, track_interval: floats\n Parameters of the members.\n\n Returns\n ------------\n reflectivityFcst: xarray.DataArray\n Contains foreacsted reflectivity.\n \"\"\"\n\n # Generate timestamps of relevant radar scans\n latestRadarScanTime = basetime - interval\n\n roverTimestrings = timestamps_ending(\n latestRadarScanTime, duration=track_interval,\n interval=track_interval,\n format='%Y%m%d%H%M',\n exclude_end=False\n )\n\n roverTimestamps = [pd.Timestamp(t) for t in roverTimestrings]\n\n qpf_xarray = reflectivityObs.sel(time=roverTimestamps)\n standardize_attr(qpf_xarray, frame_type=FrameType.dBZ, zero_value=9999.0)\n\n # Generate motion field using ROVER\n motion = rover(\n qpf_xarray,\n start_level=start_level,\n max_level=max_level,\n rho=rho,\n sigma=sigma,\n alpha=alpha\n )\n\n # Perform SLA\n steps = int(duration / interval)\n reflectivityFcst = sla(qpf_xarray, motion, steps)\n\n reflectivityFcst = reflectivityFcst.expand_dims(\n dim=OrderedDict({'member': [member]})\n ).copy()\n\n return reflectivityFcst" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Then, run the different members by calling the function above.\nIn this example. only four members are run.\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# Define the member parameters\n# Ensure that values at the same position across\n# the tuples correspond to the same member\nstartLevelTuple = (1, 2, 1, 1)\nmaxLevelTuple = (7, 7, 7, 7)\nsigmaTuple = (2.5, 2.5, 1.5, 1.5)\nrhoTuple = (9., 9., 9., 9.)\nalphaTuple = (2000., 2000., 2000., 2000.)\nmemberTuple = ('Mem-1', 'Mem-2', 'Mem-3', 'Mem-4')\ntrackIntervalList = [pd.Timedelta(i, 'm') for i in [6, 12, 6, 12]]\n\nprocesses = 4\n\nreflectivityFcstList = []\nfor i in range(processes):\n args = (\n reflectivityObs,\n basetime,\n pd.Timedelta(6, 'm'),\n pd.Timedelta(3, 'h'),\n memberTuple[i],\n startLevelTuple[i],\n maxLevelTuple[i],\n rhoTuple[i],\n alphaTuple[i],\n sigmaTuple[i],\n trackIntervalList[i]\n )\n reflectivityFcst = ensemble_rover_sla(*args)\n reflectivityFcstList.append(reflectivityFcst)\n\nrover_sla_time = pd.Timestamp.now()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Concatenating observed and forecasted reflectivities\n\nFrom the results from all members:\n\n#. Change the timestamps of observed reflectivity\n and forecasted reflectivity of all members\n from start time to end time.\n#. Add forecasted reflectivity to reproj_reflectivity_list.\n#. Concatenate observed and forecasted reflectivity\n xarray.DataArrays along the time dimension.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# List to store combined reflectivity arrays of each member\nreflectivityConcatList = []\n\n# Changing the time coordinates from start time to end time\nreflectivityObs.coords['time'] = [\n pd.Timestamp(t) + pd.Timedelta(6, 'm')\n for t in reflectivityObs.time.values\n]\n\nfor reflectivityFcst in reflectivityFcstList:\n # Identify the member\n member = reflectivityFcst.member.values[0]\n # Add member dimension to observed reflectivity\n reflectivityObsExpanded = reflectivityObs.expand_dims(\n dim=OrderedDict({'member': [member]})\n ).copy()\n\n # Checking for the track interval of the forecasted\n # array\n timeIntervalsSet = set(np.diff(reflectivityFcst.time.values))\n if pd.Timedelta(12, 'm').to_timedelta64() in timeIntervalsSet:\n reflectivityFcst.coords['time'] = [\n pd.Timestamp(t) + pd.Timedelta(12, 'm')\n for t in reflectivityFcst.time.values\n ]\n else:\n reflectivityFcst.coords['time'] = [\n pd.Timestamp(t) + pd.Timedelta(6, 'm')\n for t in reflectivityFcst.time.values\n ]\n\n # Combining\n reflectivityList = [\n reflectivityObsExpanded,\n reflectivityFcst.isel(time=slice(1, None))\n ]\n reflectivity = xarray.concat(reflectivityList, dim='time')\n standardize_attr(reflectivity, frame_type=FrameType.dBZ,\n zero_value=reflectivityFcst.attrs['zero_value'])\n reflectivityConcatList.append(reflectivity)\n\nconcat_time = pd.Timestamp.now()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Accumulating hourly rainfall for 3-hour forecast\n\nHourly accumulated rainfall is calculated every 30 minutes, the first\nendtime is the basetime i.e. T+0min.\n\n#. Convert reflectivity in dBZ to rainrates in mm/h with dbz2rr().\n#. For rainrates every 12 minutes, interpolate rainfall\n temporally to once every 6 minutes to ensure that rainrates are\n equally spaced along the time axis.\n#. Convert rainrates to rainfalls with rr2rf().\n#. Compute hourly accumulated rainfall every 30 minutes.\n#. Concatenate accumulated rainfall of different members along the\n member dimension.\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "accList = []\nfor reflectivity in reflectivityConcatList:\n # Converting reflectivity to rainrates\n rainrates = to_rainfall_rate(\n reflectivity, logarithmic=False, a=58.53, b=1.56)\n\n # Interpolate rainfall to every 12 minutes if\n # the member has a track interval of 12 minutes\n if pd.Timedelta(12, 'm').to_timedelta64() in timeIntervalsSet:\n rainrates = temporal_interpolate(rainrates,\n basetime - pd.Timedelta(1, 'h'),\n basetime + pd.Timedelta(3, 'h'),\n pd.Timedelta(6, 'm'),\n 'linear')\n # Convert rainrates to rainfalls every 6 minutes\n rainfalls = to_rainfall_depth(rainrates)\n\n # Accumulate rainfalls every 30 minutes\n accr = acc_rainfall_depth(rainfalls,\n basetime,\n basetime + pd.Timedelta(3, 'h'),\n pd.Timedelta(30, 'm'),\n pd.Timedelta(1, 'h'))\n accList.append(accr)\n\n# Concatenating along the mmeber dimension\nacc_rf = xarray.concat(accList, dim='member')\n\nacc_time = pd.Timestamp.now()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Calculating Probability Exceeding Threshold\n\nIn this example, thresholds are 0.5mm, 5mm,\n10mm, 30mm, 50mm and 70 mm.\n\nProbabilities are 0%, 25%, 50%, 75% and 100%, as there are\nfour members.\n\nSteps:\n\n#. Define threshold.\n#. Use a loop to calculate the probability of exceeding threshold\n for each gridcell and store results in list.\n#. Concatenate the contents of the list. Result is an xarray with\n dimensions (threshold, time, y, x).\n#. Plot the results using xarray.plot(). In this example,\n only the probability exceeding 0.5mm rainfall\n every hour will be plotted.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# Define threshold\nthreshold = [0.5, 5., 10., 30., 50., 70.]\n\n# Define list to store probabilities of exceeding rainfall thresholds\n# List corresponds to threshold\nprob_list = []\n\n# Calculate probability\nfor th in threshold:\n bool_forecast = acc_rf >= th\n count = bool_forecast.sum(dim='member')\n prob = count/len(bool_forecast.coords['member'].values) * 100\n prob_list.append(prob)\n\n# Generate coordinate xarray for threshold\nth_xarray = xarray.IndexVariable('threshold', threshold)\n\n# concatenate\nprob_rainfall = xarray.concat(prob_list, dim=th_xarray)\n\nprob_rainfall.attrs['long_name'] = \"Probability of Exceeding Threshold\"\nprob_rainfall.attrs['units'] = \"%\"\n\n# Plot the results\n\n# Extracting the crs\ncrs = area_def_tgt.to_cartopy_crs()\n\n# Defining cmap\ncmap = LinearSegmentedColormap.from_list(\n 'custom blue', ['#FFFFFF', '#000099']\n)\n\n# Obtaining xarray slices to be plotted\nthreshold_list = [0.5]\ntimelist = [basetime + pd.Timedelta(hours=i) for i in range(4)]\nda_plot = prob_rainfall.sel(threshold=threshold_list, time=timelist)\n\n# Defining coastlines\nmap_shape_file = os.path.join(THIS_DIR, \"./../tests/samples/shape/rsmc\")\nocean_color = np.array([[[178, 208, 254]]], dtype=np.uint8)\nland_color = cfeature.COLORS['land']\ncoastline = cfeature.ShapelyFeature(\n list(shpreader.Reader(map_shape_file).geometries()),\n ccrs.PlateCarree()\n)\n\nfor th in threshold_list:\n # Plotting\n p = da_plot.sel(threshold=th).plot(\n col='time', col_wrap=2,\n subplot_kws={'projection': crs},\n cbar_kwargs={'ticks': [0, 25, 50, 75, 100],\n 'format': '%.3g'},\n cmap=cmap\n )\n\n for idx, ax in enumerate(p.axes.flat):\n # ocean\n ax.imshow(np.tile(ocean_color, [2, 2, 1]),\n origin='upper',\n transform=ccrs.PlateCarree(),\n extent=[-180, 180, -180, 180],\n zorder=-1)\n # coastline, color\n ax.add_feature(coastline,\n facecolor=land_color, edgecolor='none', zorder=0)\n # overlay coastline without color\n ax.add_feature(coastline, facecolor='none',\n edgecolor='gray', linewidth=0.5, zorder=3)\n ax.gridlines() # gridlines\n ax.set_title(\n f\"% Exceeding {th}mm rainfall\\n\"\n f\"Based @ {basetime.strftime('%H:%MH')}\",\n loc='left',\n fontsize=8\n )\n ax.set_title(\n ''\n )\n ax.set_title(\n f\"{basetime.strftime('%Y-%m-%d')} \\n\"\n f\"Valid @ {timelist[idx].strftime('%H:%MH')} \",\n loc='right',\n fontsize=8\n )\n plt.savefig(\n THIS_DIR +\n f\"/../tests/outputs/p_{th}.png\",\n dpi=300\n )\n\nprob_time = pd.Timestamp.now()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Rainfall percentiles\n\n#. Using xarray.DataArray.mean(),\n calculate the mean rainfall of all gridpoints.\n#. Using xarray.DataArray.min(),\n find the minimum rainfall of all gridpoints.\n#. Using xarray.DataArray.max(),\n find the maximum rainfall of all gridpoints.\n#. Using xarray.DataArray.quantile()\n find the 25th, 50th and 75th percentile rainfall\n of all gridpoints.\n#. Concatenate rainfall along percentile dimension.\n#. Plot results using xarray.plot(). In this example\n only the maximum percentile rainfall\n every hour will be plotted.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# mean\nmean_rainfall = acc_rf.mean(dim='member', keep_attrs=True)\n# max\nmax_rainfall = acc_rf.max(dim='member', keep_attrs=True)\n# min\nmin_rainfall = acc_rf.min(dim='member', keep_attrs=True)\n# quartiles\nq_rainfall = acc_rf.quantile(\n [.25, .5, .75], dim='member',\n interpolation='nearest', keep_attrs=True)\n\n# generate index\npercentile = xarray.IndexVariable(\n 'percentile',\n ['Minimum',\n '25th Percentile',\n '50th Percentile',\n 'Mean',\n '75th Percentile',\n 'Maximum']\n)\n\n# concatenating\n\np_rainfall = xarray.concat(\n [min_rainfall,\n q_rainfall.sel(quantile=.25).squeeze().drop('quantile'),\n q_rainfall.sel(quantile=.50).squeeze().drop('quantile'),\n mean_rainfall,\n q_rainfall.sel(quantile=.75).squeeze().drop('quantile'),\n max_rainfall],\n dim=percentile\n)\n\np_rainfall.attrs['long_name'] = \"Hourly Accumulated Rainfall\"\np_rainfall.attrs['units'] = \"mm\"\n\n# Defining levels\n\n# Defining the colour scheme\nlevels = [\n 0, 0.5, 2, 5, 10, 20,\n 30, 40, 50, 70, 100, 150,\n 200, 300, 400, 500, 600, 700\n]\n\ncmap = matplotlib.colors.ListedColormap([\n '#ffffff', '#9bf7f7', '#00ffff', '#00d5cc', '#00bd3d', '#2fd646',\n '#9de843', '#ffdd41', '#ffac33', '#ff621e', '#d23211', '#9d0063',\n '#e300ae', '#ff00ce', '#ff57da', '#ff8de6', '#ffe4fd'\n])\n\nnorm = BoundaryNorm(levels, ncolors=cmap.N, clip=True)\n\n# Obtaining xarray slices to be plotted\nplot_positions = ['Maximum']\ntimelist = [basetime + pd.Timedelta(hours=i) for i in range(4)]\nda_plot = p_rainfall.sel(percentile=plot_positions, time=timelist)\n\nfor pos in plot_positions:\n # Plotting\n p = da_plot.sel(percentile=pos).plot(\n col='time', col_wrap=2,\n subplot_kws={'projection': crs},\n cbar_kwargs={'ticks': levels,\n 'format': '%.3g'},\n cmap=cmap,\n norm=norm\n )\n\n for idx, ax in enumerate(p.axes.flat):\n # ocean\n ax.imshow(np.tile(ocean_color, [2, 2, 1]),\n origin='upper',\n transform=ccrs.PlateCarree(),\n extent=[-180, 180, -180, 180],\n zorder=-1)\n # coastline, color\n ax.add_feature(coastline,\n facecolor=land_color, edgecolor='none', zorder=0)\n # overlay coastline without color\n ax.add_feature(coastline, facecolor='none',\n edgecolor='gray', linewidth=0.5, zorder=3)\n ax.gridlines()\n ax.set_title(\n f\"{pos} Radar-Based Rainfall\\n\"\n f\"Based @ {basetime.strftime('%H:%MH')}\",\n loc='left',\n fontsize=8\n )\n ax.set_title(\n ''\n )\n ax.set_title(\n f\"{basetime.strftime('%Y-%m-%d')} \\n\"\n f\"Valid @ {timelist[idx].strftime('%H:%MH')} \",\n loc='right',\n fontsize=8\n )\n position = pos.split(\" \")[0]\n plt.savefig(\n THIS_DIR +\n f\"/../tests/outputs/rainfall_{position}.png\",\n dpi=300\n )\nptime = pd.Timestamp.now()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Extract the rainfall values at a specified location\n\nIn this example, the rainfall values at the location is assumed\nto be the same as the nearest gridpoint.\n\n#. Read information regarding the radar stations into a pandas.DataFrame.\n#. Extract the rainfall values at the nearest gridpoint to location\n for given timesteps (in this example, 30 minute intervals).\n#. Store rainfall values over time in an xarray.DataArray.\n#. Plot the time series of rainfall with boxplots at desired station.\n In this case, the 15th percentile member is plotted.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# Getting radar station coordinates\ndf = pd.read_csv(\n os.path.join(THIS_DIR, \"../tests/samples/hk_raingauge.csv\"),\n usecols=[0, 1, 2, 3, 4]\n)\n\n# Extract rainfall values at gridpoint closest to the\n# location specified for given timesteps and storing it\n# in xarray.DataArray.\nstation_rf_list = []\nstation_name = []\nfor index, row in df.iterrows():\n station_rf_list.append(p_rainfall.sel(\n northing=row[1], easting=row[2],\n method='nearest'\n ).drop('northing').drop('easting'))\n station_name.append(row[0])\n\nstation_name_index = xarray.IndexVariable(\n 'ID', station_name\n)\nstation_rf = xarray.concat(\n station_rf_list,\n dim=station_name_index\n)\n\n# Extracting the 15th ranked station\nxr_15_percentile = station_rf.quantile(\n .15, dim='ID', interpolation='nearest').drop('quantile')\n\n# Plotting\n_, tax = plt.subplots(figsize=(20, 15))\nplt.plot(\n np.arange(1, len(xr_15_percentile.coords['time'].values) + 1),\n xr_15_percentile.loc['Mean'].values,\n 'ko-'\n) # plot line\n\n# Storing percentiles as dictionary to call ax.bxp for boxplot\nstats_list = []\nfor i in range(len(xr_15_percentile.coords['time'].values)):\n stats = {\n 'med': xr_15_percentile.loc['50th Percentile'].values[i],\n 'q1': xr_15_percentile.loc['25th Percentile'].values[i],\n 'q3': xr_15_percentile.loc['75th Percentile'].values[i],\n 'whislo': xr_15_percentile.loc['Maximum'].values[i],\n 'whishi': xr_15_percentile.loc['Minimum'].values[i]\n }\n stats_list.append(stats)\n\n# Plot boxplot\ntax.bxp(\n stats_list, showfliers=False\n)\n\n# Labels\nxcoords = xr_15_percentile.coords['time'].values\nxticklabels = [pd.to_datetime(str(t)).strftime(\"%H:%M\") for t in xcoords]\ntax.set_xticklabels(xticklabels)\ntax.xaxis.set_tick_params(labelsize=20)\ntax.yaxis.set_tick_params(labelsize=20)\nplt.title('Time Series of Hourly Accumulated Rainfall', fontsize=25)\nplt.ylabel(\"Hourly Accumulated Rainfall [mm]\", fontsize=22)\nplt.xlabel(\"Time\", fontsize=18)\nplt.savefig(THIS_DIR+\"/../tests/outputs/pqpf_time_series.png\")\n\nextract_time = pd.Timestamp.now()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Checking run time of each component\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "print(f\"Start time: {start_time}\")\nprint(f\"Initialising time: {initialising_time}\")\nprint(f\"Rover and SLA time: {rover_sla_time}\")\nprint(f\"Concatenating time: {concat_time}\")\nprint(f\"Accumulating rainfall time: {acc_time}\")\nprint(\"Calculate and plot probability exceeding threshold: \"\n f\"{prob_time}\")\nprint(f\"Plotting rainfall maps: {ptime}\")\nprint(f\"Extracting and plotting time series time: {extract_time}\")\n\nprint(f\"Time to initialise: {initialising_time-start_time}\")\nprint(f\"Time to run rover and SLA: {rover_sla_time-initialising_time}\")\nprint(f\"Time to concatenate xarrays: {concat_time - rover_sla_time}\")\nprint(f\"Time to accumulate rainfall: {acc_time - concat_time}\")\nprint(\"Time to calculate and plot probability exceeding threshold: \"\n f\"{prob_time-acc_time}\")\nprint(f\"Time to plot rainfall maps: {ptime-prob_time}\")\nprint(f\"Time to extract station data and plot time series: \"\n f\"{extract_time-ptime}\")" ] } ], "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.6.15" } }, "nbformat": 4, "nbformat_minor": 0 }