baltrad_pyart_rain_rate_example.ipynb

{"metadata":{"language_info":{"name":"python","version":"3.9.13","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kernelspec":{"name":"python3","display_name":"Python 3","language":"python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# In this notebook, an ODIM_H5 file is read using BALTRAD.  Then the rain rate is determined from the calculated specific attenuation using Py-ART.\n## This is a severe flooding case from July 8, 2013 in Toronto, Canada, with radar data from the King City, Ontario, radar.","metadata":{"tags":[]}},{"cell_type":"code","source":"%matplotlib inline","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Import the necessary modules.","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pyart\nimport baltrad_pyart_bridge as bridge  # routines to pass data from Py-ART and BALTRAD\nimport _raveio   # BALTRAD's input/output module","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Read in the data using RAVE (a component of BALTRAD)","metadata":{}},{"cell_type":"markdown","source":"# Rain rate retrieval using specific attenuation using BALTRAD and Py-ART","metadata":{}},{"cell_type":"code","source":"rio = _raveio.open('data/WKR_201307082030.h5')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Convert the data to a Py-ART Radar object.","metadata":{}},{"cell_type":"code","source":"radar = bridge.raveio2radar(rio)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Examine some of the radar moments.","metadata":{}},{"cell_type":"code","source":"display = pyart.graph.RadarDisplay(radar)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display.plot_ppi('DBZH', 0, vmin=-15, vmax=60)\ndisplay.plot_range_rings([50, 100, 150])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display.plot_ppi('PHIDP', 0, vmin=0, vmax=180)\ndisplay.plot_range_rings([50, 100, 150])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display.plot_ppi('RHOHV', 0, vmin=0, vmax=1.0, mask_outside=False)\ndisplay.plot_range_rings([50, 100, 150])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display.plot_ppi('SQIH', 0, vmin=0, vmax=1, mask_outside=False)\ndisplay.plot_range_rings([50, 100, 150])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Calculate the specific attenuation and attenuation corrected reflectivity using Py-ART, add these field to the radar object.","metadata":{}},{"cell_type":"code","source":"spec_at, cor_z = pyart.correct.calculate_attenuation(\n    radar, 0, doc=0, refl_field='DBZH', ncp_field='SQIH', \n    rhv_field='RHOHV', phidp_field='PHIDP', \n    fzl=8000,)\n# use the parameter below for a more 'cleanup up' attenuation field\n#ncp_min=-1, rhv_min=-1)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"radar.add_field('specific_attenuation', spec_at)\nradar.add_field('corrected_reflectivity', cor_z)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Examine these two new fields.","metadata":{}},{"cell_type":"code","source":"display.plot_ppi('specific_attenuation', 0, vmin=0, vmax=0.1)\ndisplay.plot_range_rings([50, 100, 150])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display.plot_ppi('corrected_reflectivity', 0, vmin=-15, vmax=60)\ndisplay.plot_range_rings([50, 100, 150])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Calculate the rain rate from the specific attenuation using a power law determined from the ARM Southern Great Plains site.  Mask values where the attenuation is not valid (when the cross correlation ratio or signal quality is low). Add this field to the radar object.","metadata":{}},{"cell_type":"code","source":"R = 300.0 * (radar.fields['specific_attenuation']['data']) ** 0.89\nrain_rate_dic = pyart.config.get_metadata('rain_rate')\nrain_rate_dic['units'] = 'mm/hr'\nrate_not_valid = np.logical_or(\n    (radar.fields['SQIH']['data'] < 0.4),\n    (radar.fields['RHOHV']['data'] < 0.8))\nrain_rate_dic['data'] = np.ma.masked_where(rate_not_valid, R)\n# fill the missing values with 0 for a nicer plot\nrain_rate_dic['data'] = np.ma.filled(rain_rate_dic['data'], 0)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"radar.add_field('RATE', rain_rate_dic)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Examine the rain rate","metadata":{}},{"cell_type":"code","source":"display.plot_ppi('RATE', 0, vmin=0, vmax=50.0)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Create a new RaveIO object from the Py-ART radar object and write this out using Rave","metadata":{}},{"cell_type":"code","source":"rio_out = bridge.radar2raveio(radar)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"container = _raveio.new()\ncontainer.object = rio_out.object\ncontainer.save(\"data/WKR_201307082030_with_rain_rate.h5\")\n\nimport os\nprint(\"ODIM_H5 file is %i bytes large\" % os.path.getsize(\"data/WKR_201307082030_with_rain_rate.h5\"))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Publish a time series of Cartesian products of corrected reflectivity to BALTRAD's GoogleMapsPlugin","metadata":{}},{"cell_type":"markdown","source":"### Using your Browser, preferably anything except Microsoft Internet Explorer, view a pre-loaded product:  http://localhost:8080  Use the small Calendar icon in the control panel to select 2013-07-08 20:30. The dropdown box under the date/time field should read \"King City, ON\".","metadata":{}},{"cell_type":"markdown","source":"### Fire up the RAVE Product Generation Framework's server. This is normally done on the command line.","metadata":{}},{"cell_type":"code","source":"import os\nos.system(\"rave_pgf start\")","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Connect to this XML-RPC server and feed it file strings of pre-generated products","metadata":{}},{"cell_type":"code","source":"import glob, xmlrpc.client\n\nipath = \"/home/vagrant/pyart2baltrad/data/cawkr\"\nopath = \"/home/vagrant/miniconda/envs/openradar/rave_gmap/web/data/cawkr_gmaps\"\n\nserver = xmlrpc.client.ServerProxy(\"http://localhost:8085/RAVE\")\n\nfstrs = glob.glob(ipath + \"/*.h5\")\n\nfor ifstr in fstrs:\n    # Output file name must only be date/time string with format: YYYYMMDDHHmm.png\n    dt = os.path.split(ifstr)[1].split('_')[2].split('.')[0]\n    \n    ofstr = opath + \"/%s/%s/%s/%s.png\" % (dt[:4], dt[4:6], dt[6:8], dt)\n    response = server.generate(\"se.smhi.rave.creategmapimage\", [ifstr], [\"outfile\",ofstr])\n\nprint(\"Generated %i PNG images for Google Maps\" % len(fstrs))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Go back to your browser and load a sequence of images. Stop the PGF server.","metadata":{}},{"cell_type":"code","source":"os.system(\"rave_pgf stop\")","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}