HTTPS SSH
WoF Processing/Visualization/Verification Python Scripts: 
(Please refer to README.PDF for more details)

Purpose:  

  This collection of python scripts is intended to be used to 
  post-process, plot, and verify Weather Research and Forecasting (WRF) 
  ensemble forecasts.  The scripts were developed for use with the 
  NSSL Experimental Warn-on-Forecast System for ensembles (NEWS-e; 
  Wheatley et al. 2015; Jones et al. 2016), but should be applicable 
  to any WRF output files in netCDF format with minimal modification.  

Version History:

  o Initial version by Patrick Skinner on 11/28/2016.

**Please drop an email to patrick.skinner@noaa.gov for using this package. Thanks for your interest.
**

Dependencies: 

  All code included in the collection has been written using the 
  Anaconda python build (version 2.7.11; https://www.continuum.io/downloads).  
  Non-standard libraries used by the code include: 

  o netCDF4:  Library used for reading and writing netCDF files 
              (http://unidata.github.io/netcdf4-python/)

  o basemap:  Library used for projecting and plotting geospatial data
              (http://matplotlib.org/basemap/)

  o scikit-image:  Computer vision library used for object identification 
                   and verification (http://scikit-image.org/)

  The code has been run on both Linux and OSX systems, but has not been 
  tested for Windows. 

Script Descriptions: 

  The scripts can be divided into five different categories: 

  1. ‘Cookbooks’:  These are collections of subroutines and objects 
     used by subsequent post-processing and plotting scripts.  

    • news_e_post_cbook.py – A collection of subroutines for calculating 
      derived variables from WRF output files.  Subroutines are available 
      for calculation of both common environmental quantities (e.g. 
      CAPE/Shear/STP) and storm-scale quantities (e.g. Updraft helicity/
      layer-mean vertical vorticity). 

    • news_e_plotting_cbook_v2.py – A collection of objects and subroutines 
      used to create plots found on the NEWS-e realtime website 
      (http://www.nssl.noaa.gov/projects/wof/news-e/images.php).  
      Customized colors from http://colorbrewer2.org/ are included as an 
      object; and several variations for plotting geospatial forecast 
      data are included as subroutines. 

    • radar_info.py – A simple class that contains the identifier, latitude, 
      and longitude for all CONUS WSR-88d sites. 

    • ctables.py – A collection of customized color maps for plotting. 

  2. Base post-processing/plotting scripts:   These are scripts used to 
     post-process WRF forecast output into smaller summary netCDF files 
     and produce basic environmental/storm-scale plots using the summary 
     files.  These scripts are written to process/plot single timesteps 
     or variables to facilitate multiprocessing using the ‘wrapper’ scripts. 

    • news_e_post_retro.py – The primary post-processing script.  This 
      will read over directories of WRF output files and produce a single 
      summary file of environmental and storm-scale quantities for each 
      forecast (i.e. for an ensemble forecast with 18 members, it will 
      produce 18 summary files for plotting).  The summary files produced 
      by this script are used by all subsequent plotting/verification scripts. 

    • news_e_pmm_retro.py – This script creates a file of 
      probability-matched mean fields (currently set to composite 
      reflectivity) for a given forecast.  It is done separately from the 
      other post-processing to speed up processing for real-time experiments.  
      The ‘pmm’ file produced by this script is used by most subsequent 
      plotting/verification scripts (but dependency can be easily removed 
      if desired).

    • news_e_timestep_retro.py – This script will produce ensemble-mean 
      environmental plots (e.g. Temp/Td/CAPE) when provided with a 
      directory of summary files produced by news_e_post_retro.py. 

    • news_e_swath_retro.py – This script will produce storm-scale 
      ensemble product plots (e.g. probability of exceedance/90th 
      percentile plots of 0-2 km UH) when provided with a directory of 
      summary files produced by news_e_post_retro.py

  3. Secondary post-processing/plotting scripts:  These scripts are 
     primarily used to produce verification products for the WRF 
     forecast using MRMS (http://www.nssl.noaa.gov/projects/mrms/) 
     merged radar products. These scripts are written to process/plot 
     single timesteps or variables to facilitate multiprocessing 
     using the ‘wrapper’ scripts.

    • mrms_post_cressman_qc_retro.py – This script will interpolate 
      merged MRMS products produced by the WDSSii software 
      (http://www.wdssii.org/) to the grid of the WRF forecasts being 
      processed using a Cressman scheme.  Additionally, for rotation 
      fields (i.e. azimuthal wind shear) rotation track objects will 
      be identified and quality controlled.

    • news_e_rot_object_retro.py – This script will produce 
      quality-controlled rotation track objects for WRF rotation products 
      (i.e. vertical vorticity/UH) corresponding to the MRMS objects 
      produced by mrms_post_cressman_qc_retro.py. 

    • news_e_verif_lsr_timestep.py – This script will produce similar 
      plots to news_e_swath_retro.py, except with NWS hail/wind/tornado 
      Local Storm Reports (LSR) overlain.  The script expects an LSR 
      shapefile downloaded from: 
      https://mesonet.agron.iastate.edu/request/gis/lsrs.phtml

    • news_e_verif_paintqc_timestep.py – This script will produce 
      paintball plots of ensemble rotation/dBZ objects vs. corresponding 
      MRMS objects.  It requires output from both mrms_post_cressman_qc_retro.py 
      and news_e_rot_object_retro.py.  Additionally, NWS warning 
      products and LSRs may be overlain, which require shapefiles from: 
      https://mesonet.agron.iastate.edu/request/gis/watchwarn.phtml and 
      https://mesonet.agron.iastate.edu/request/gis/lsrs.phtml, respectively. 

    • news_e_verif_objqc_timestep.py – This script will produce 
      object-based verification plots of ensemble rotation/dBZ objects vs. 
      corresponding MRMS objects (Skinner et al. 2016).  It requires 
      output from both mrms_post_cressman_qc_retro.py and 
      news_e_rot_object_retro.py.

    • mrms_timestep_retro.py – Produces plots of interpolated MRMS 
      reflectivity and rotation track objects similar to those for 
      the WRF forecast for comparison. 

  4. ‘Wrapper’ scripts:  These scripts are used to apply multiprocessing 
     to most of the post-processing/plotting scripts.  They are not necessary 
     to run any of the post/plotting/verification scripts, but can be used 
     to make them run faster.

    • news_e_post_wrapper_retro.py – Runs news_e_post_retro.py for each 
      WRF output file for a given forecast (assuming all WRF output 
      files are in a single directory). 

    • news_e_pmm_wrapper_retro.py – Runs news_e_pmm_retro.py for summary 
      files of a given WRF forecast. 

    • news_e_plot_wrapper_retro.py – Runs both news_e_timestep_retro.py 
      and news_e_swath_retro.py for summary files of a given WRF forecast. 

    • news_e_verif_lsr_wrapper_retro.py – Runs news_e_verif_lsr_retro.py 
      for summary files of a given WRF forecast. 

    • news_e_verif_paintqc_wrapper_retro.py – Runs 
      news_e_verif_paintqc_retro.py for summary files of a given WRF forecast. 

    • news_e_verif_objqc_wrapper_retro.py – Runs 
      news_e_verif_objqc_retro.py for summary files of a given WRF forecast. 

    • mrms_wrapper_retro.py – Runs mrms_timestep_retro.py for 
      interpolated MRMS dBZ/rotation track files. 

  5. ‘Run’ scripts:  These CSH scripts are used to run post/plotting scripts 
     from the compute nodes of the Loki supercomputer at NSSL.  They are 
     not needed if running off of Loki. 

Sample Workflow – Post-process and produce base plots for an ensemble of 
  WRF forecasts:

  1. Produce summary files for each ensemble member forecast using 
     news_e_post_wrapper_retro.py. Ex:

     $> python news_e_post_wrapper_retro.py -d <output directory for summary files> -f <directories of WRF output files> -n <number of ensemble members> 

     The script will expect the ensemble forecast to have a similar 
     directory structure to NEWS-e forecasts, with separate directories 
     for each member (named ‘ENS_MEM_#’) containing netCDF WRF output 
     files (beginning with ‘wrfout’) for each timestep in the forecast.  

     The script will output netCDF summary files named ‘YYYY-MM-DD-hh:mm:ss_nn.nc’ 
     for each member in the directory specified by ‘-d’.  The date/time in 
     the filename will correspond to the initialization time of the 
     forecast and ‘nn’ refers to the ensemble member.  

  2. Produce a probability-matched mean file for composite reflectivity 
     using the summary files using news_e_pmm_wrapper_retro.py.  If you 
     do not wish to plot probability-matched mean reflectivity this step 
     can be skipped with minor modifications to news_e_timestep_retro.py 
     and news_e_swath_retro.py. Ex:

     $> python news_e_pmm_wrapper_retro.py -d <output directory for summary files> -n <number of ensemble members>

     The script will read in summary files produced in the previous 
     step and output a file named ‘pmm_dz.nc’ containing the 
     probability-matched mean composite reflectivity for each 
     timestep in the forecast.

  3. Create environmental and ensemble product plots for the summary 
     files using news_e_plot_wrapper_retro.py. Ex:

     $> python news_e_plot_wrapper_retro.py -d <input directory of summary files> -i <output directory of .png images> -n <number of ensemble members>
 
     The script will read in summary files produced in step 1 as well 
     as the probability-matched mean file produced in step 2, then output 
     .png images for each product defined in news_e_timestep_retro.py and 
     news_e_swath_retro.py to the directory provided in the ‘-i’ option.  

     Images will be produced for each timestep in the forecast and are 
     named ‘<product_name>_f###’, where the ‘###’ following ‘f’ 
     corresponds to either forecast time or percentile/probability 
     threshold value of the plot. 

Extended Workflow – Produce secondary and verification plots:

  4. Produce interpolated MRMS verification fields/objects using 
     mrms_post_cressman_qc_retro.py. Ex:

     $> python mrms_post_cressman_qc_retro.py -d <directory of MRMS az. Shear WDSS ii files> 
                -z <directory of MRMS reflectivity WDSS ii files> -o <path to output netCDF file> 
                -f <path to example summary file of forecast to verify> -v <MRMS variable to process> 
                -s <start time (HHMMSS UTC)> -e <end time (HHMMSS UTC)>

     This script will read in directories of merged azimuthal wind shear 
     and reflectivity MRMS products then interpolate the raw MRMS data to 
     the WRF forecast grid using a Cressman scheme.  Additionally, 
     rotation track objects will be produced for azimuthal shear products, 
     objects will be determined based on a specified time window, product 
     layer (e.g. 2-5 km), and quality control thresholds for intensity, area, 
     continuity, and distance from a reflectivity object that are intended 
     to minimize spurious objects in the output.  Finally, a ‘radmask’ 
     field containing regions to mask either too close or too far from the 
     nearest WSR-88d site (determined using the radar_info.py script) will 
     be created. 

     A single output file containing interpolated and object data for the 
     specified variable at 5-min timesteps between the start and end times 
     will be produced.  

  5. Produce corresponding forecast rotation objects using 
     news_e_rot_object_retro.py.  Ex: 

     $> python news_e_rot_object_retro.py -d <input directory of summary files>  
               -o <path to output rotation object file> -f <path to MRMS azimuthal wind shear object file> 

     Rotation swath objects will be created in an identical manner to step 4, 
     except for WRF forecasts rather than MRMS azimuthal wind shear.  The 
     ‘radmask’ created for MRMS data in step 4 will be applied to the output 
     rotation track objects to ensure they are comparable to the MRMS 
     rotation objects. 

  6. Create LSR verification products using news_e_verif_lsr_wrapper.py. 

     $> python news_e_verif_lsr_wrapper.py -d <input directory of summary files> -i <output directory for .png images> 
               -w <input shapefile of NWS warnings> -l <input shapefile of NWS LSRs> 
               -z <input file of MRMS reflectivity interpolated to forecast grid> -a <input file of MRMS az. shear interpolated to forecast grid> 

     Storm-scale ensemble products similar to those created in step 3 will 
     be produced, except that NWS LSRs or flash flood warnings (QPF products 
     only) will be overlain. 

  7. Create paintball verification products using news_e_verif_paintqc_wrapper.py. 
     Ex: 

     $> python news_e_verif_paintqc_wrapper.py -d <input directory of summary files> -i <output directory for .png images> -w <input shapefile of NWS warnings> -l <input shapefile of NWS LSRs> -z <input file of MRMS reflectivity interpolated to forecast grid> -a <input file of low-level MRMS az. shear interpolated to forecast grid> -m <input file of midlevel MRMS az. shear interpolated to forecast grid> 

     Paintball plots of forecast (output of step 5) vs. observed (output of 
     step 4) rotation track and reflectivity objects will be produced.  Objects 
     from each ensemble member will be plotted in varying colors and the 
     observed objects will be overlain in dark gray.  Additionally, NWS 
     warnings and LSRs will be overlain. 


  8. Create object-based verification products using news_e_verif_objqc_wrapper.py. 
     Ex: 

     $> python news_e_verif_objqc_wrapper.py -d <input directory of summary files> -i <output directory for .png images> -w <input shapefile of NWS warnings> -l <input shapefile of NWS LSRs> -z <input file of MRMS reflectivity interpolated to forecast grid> -a <input file of low-level MRMS az. shear interpolated to forecast grid> -m <input file of midlevel MRMS az. shear interpolated to forecast grid>

     Forecast objects will be classified as ‘matches’ or ‘false alarms’ 
     according to an interest score determined by spatiotemporal distance 
     from the nearest observed object.  Gridpoint probability plots of 
     being within a matched or false alarm object will be produced with 
     observed objects overlain in dark gray.  Additionally, a weighted 
     and binary object-based threat score (OTS; Johnson and Wang 2011) 
     will be calculated and annotated in the lower right corner of the 
     plots.  

  9. Create plots of MRMS reflectivity/rotation tracks for comparison using mrms_wrapper_retro.py. 
     Ex: 

     $> python mrms_wrapper_retro.py -d <input directory of summary files> -i <output directory for .png images> -w <input shapefile of NWS warnings> -l <input shapefile of NWS LSRs> -z <input file of MRMS reflectivity interpolated to forecast grid> -a <input file of low-level MRMS az. shear interpolated to forecast grid> -m <input file of midlevel MRMS az. shear interpolated to forecast grid>

     Plots of interpolated MRMS reflectivity with either low-level or 
     midlevel azimuthal shear rotation tracks overlain in dark gray 
     will be produced.  Additionally, NWS warning and LSR products 
     will be overlain.