{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "### Pyart_gridplot_panels.py\n", "Author: djv 7/19\n", "This script will compare Raw level II reflectivity to gridded radar output \n", "created from *pyart.map.grid_from_radars.py*.
Will plot out maps of raw data at given\n", "sweep and gridded data at level closest to sweep elevation as well as cross sections of the gridded datatest using **PYART** routines:\n", " 1. *plot_latitude_slice*\n", " 2. *plot_longitude_slice* \n", " 3. *plot_ppi_map*\n", " 4. *plot_grid*\n", " \n", "**This script has several plotting sections:**\n", " 1. 4-panel plot of radars+cross sections\n", " 2. Widget plots:\n", " *- 2-panel plot of map of gridded radar with longitudinal cross section with adjustable slider widget\n", " *- 2-panel plot of map of gridded radar with latitudinal cross section with adjustable slider widget\n", " *- 2-panel plot of latitudinal and longitudinal cross section with adjustable slider widgets\n", " \n", "Program Needs several user input files:\n", " A Level II radar archive file\n", " A gridded version of this data created from pyart.map.grid_from_radars.py\n", " \n", " Will also run python script args_gridplt.py to read in input data params" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "import pyart\n", "import math\n", "import numpy as np\n", "import cartopy\n", "import matplotlib.gridspec as gridspec\n", "import matplotlib.pyplot as plt\n", "import matplotlib.colors as col\n", "import matplotlib.cm as cm\n", "from ipywidgets import interactive\n", "import ipywidgets as widgets\n", "from IPython.display import display as wdisplay\n", "from IPython.html.widgets import *\n", "import sys\n", "sys.path.insert(0, '/home/djv/python/djvlibs/') #DEFINES DJV LIB PATH\n", "from dvpylib_shape import *\n", "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Run interactive parameter program \n", "edit following line to reflect desired sweep plotting parameter and plot output\n", "plttyp=png,ps,eps,jpg,screen" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "################### MAIN PROGRAM ##################################%\n", "%run /home/djv/python/notebooks/radar/args_gridplt.py 1 3000 reflectivity png\n", "#name = raw_input(\"What's your name? \") #inputs as strings\n", "print(sweep,type(sweep))\n", "print(grd_elev,type(grd_elev))\n", "print(pltvar,type(pltvar))" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "############# USER DEF VARIABLES #################\n", "radmax=250. #radius of Lev II plot map area\n", "grid_line_delta = 1.0\n", "xsec_vlim = [0,15] #vertical limit of xsec plots(km)\n", "xsec_hlim = [-200.,200.] #horizontal limit of xsec plot(km)\n", "nexrad_site=\"KCRP\"\n", "xdate=\"20170826_002504\"\n", "\n", "datpth = \"/jm11/djv/radar/pyart/\"\n", "ifile=datpth+nexrad_site+xdate+\"_V06.ar2v\"\n", "ifile=datpth+nexrad_site+\"20170826_0025_velocityallsw.nc\" #dealiased velocity\n", "gfile=datpth+\"KCRP201708260025_all_250km_grid_41x501_dist.nc\"\n", "\n", "ofile0 = nexrad_site+xdate+pltvar+\"grdtst_4pan\" #plot output name" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def map_add_on(ax,map_res,lonvals,latvals):\n", " import matplotlib.ticker as mticker\n", " import cartopy.feature as cfeature\n", " from metpy.plots import USCOUNTIES\n", "################# CARTOPY ADD_ONS ######################\n", " ax.add_feature(USCOUNTIES.with_scale('20m'), linewidth=1,linestyle='--',edgecolor='darkgray')\n", " ax.add_feature(cfeature.LAND)\n", " ax.add_feature(cfeature.OCEAN)\n", " ax.add_feature(cfeature.LAKES.with_scale(map_res), alpha=0.5)\n", " ax.add_feature(cfeature.RIVERS.with_scale(map_res), edgecolor='tab:blue') \n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def find_nearest(array, value):\n", "# Will return index of array_val closest to value\n", "# Return: index of nearest value\n", " array = np.asarray(array)\n", " idx = (np.abs(array - value)).argmin()\n", " return idx" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Main Program Starts Here#" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "############## Creates REF colormap ###############\n", "#'pyart_NWSRef',\n", "v0=-40\n", "v1=75\n", "cb_pth=\"/home11/staff/djv/python/\"\n", "cfil=cb_pth+\"djv_radar_bin5.rgb\"\n", "if pltvar != \"reflectivity\":\n", " cfil=cb_pth+\"djv_vad_bin5.rgb\"\n", " v0=-65.\n", " v1=65.\n", "chan=\"RADAR\"\n", "dvrad,ncolor=dvColormap(cfil,chan,\"\") #RGB colorfile,channel,path\n", "dv_cmap=col.LinearSegmentedColormap('my_colormap',dvrad,N=ncolor,gamma=1.)\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "################## GRIDDED SECTION #######################\n", "rg = pyart.io.read_grid(gfile)\n", "rg.fields.keys()\n", "displayg = pyart.graph.GridMapDisplay(rg)\n", "clatg=rg.origin_latitude['data'][0]\n", "clong=rg.origin_longitude['data'][0]\n", "xs_lon=clong+0.9 #location of lat(NS) xsec wrt radar\n", "xs_lat=clatg+.1 #location of lon(WE) xsec\n", "xvals=rg.x['data'][:]\n", "grd_res=str((int((xvals[-1]-xvals[0])/(len(xvals)-1))/1000))+\"km\" #determines horiz res\n", "glev =find_nearest(rg.z['data'][:],grd_elev) #find index of desired plot elev" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "########### READ AND PLOT PREP of LEVEL II ARCHIVE DATA ##########\n", "r = pyart.io.read(ifile)\n", "r.fields.keys()\n", "display = pyart.graph.RadarMapDisplay(r) #set up for map plot\n", "clat=r.latitude['data'][0]\n", "clon=r.longitude['data'][0]\n", "lalo_dphi=(1.1*radmax)/111.11\n", "wlon=clon-lalo_dphi\n", "elon=clon+lalo_dphi\n", "blat=clat-lalo_dphi\n", "ulat=clat+lalo_dphi\n", "lab_dphi=int(abs(wlon-elon)/3.)\n", "lonvals=np.around(np.arange(wlon,elon+lab_dphi,lab_dphi),1)\n", "latvals=np.around(np.arange(blat,ulat+lab_dphi,lab_dphi),1)\n", "sdate=r.time['units'].split()\n", "plt_title=nexrad_site+\" \"+pltvar+\" swp=\"+str(sweep)+\" \"+sdate[2]\n", "map_res = '50m' #map resolution of cartopy features\n", "projection = cartopy.crs.PlateCarree()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 4-panel Plot " ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Will compare actual radar map to gridded map at level closest to sweep \n", "# elevation and plot 2 vertical cross sections\n", "map_panel_axes0 = [-0.05, 0.52, .43, .43] # sets up panels for plotting \n", "map_panel_axes1 = [0.5, 0.54, .39, .39]\n", "x_cut_panel_axes = [0.05, 0.06, .39 ,.38]\n", "y_cut_panel_axes = [0.55, 0.06, .39, .38]\n", "clat=r.latitude['data'][0]\n", "clon=r.longitude['data'][0]\n", "fig = plt.figure(figsize=[15,7])\n", "########### PLOT LEV II RADAR ON MAP ####################\n", "plt_title=nexrad_site+\" \"+pltvar+\" swp=\"+str(sweep)+\" \"+sdate[2]\n", "ax0 = fig.add_axes(map_panel_axes0, projection=projection)\n", "display.plot_ppi_map(pltvar,sweep, vmin=v0, vmax=v1,\n", " min_lon=wlon, max_lon=elon, min_lat=blat, max_lat=ulat,\n", " lon_lines=np.arange(int(wlon)-1.0,int(elon)+1.0,grid_line_delta), \n", " lat_lines=np.arange(int(blat)-1.0,math.floor(ulat)+2.,grid_line_delta), \n", " projection=projection,cmap=dv_cmap,resolution=map_res,colorbar_label=\"m/s\",\n", " lat_0=clat,lon_0=clon,title=plt_title)\n", "map_add_on(ax0,map_res,lonvals,latvals)\n", "display.plot_range_rings([50, 100]) #plots rings \n", "ax0.scatter(clon,clat,color='darkblue',edgecolor='k',linewidth='1',s=30,zorder=500)\n", "\n", "########### PLOT GRIDDED RADAR ON MAP ####################\n", "ax1 = fig.add_axes(map_panel_axes1, projection=projection)\n", "plt_gtitle=(nexrad_site+\" \"+grd_res+\" grd \"+\n", " pltvar+\" lev=\"+str( rg.z['data'][glev]/1000.)+\"kkm \\n\"+sdate[2])\n", "displayg.plot_grid(pltvar, glev, vmin=v0, vmax=v1,projection=projection,\n", " cmap=dv_cmap,title=plt_gtitle,ax=ax1,colorbar_label=\"m/s\",axislabels=(' ',' '),\n", " axislabels_flag=True)\n", "map_add_on(ax1,map_res,lonvals,latvals)\n", "### Plots station location ####\n", "ax1.scatter(clon,clat,color='darkblue',edgecolor='k',linewidth='1',s=30,zorder=500)\n", "displayg.plot_crosshairs(lon=xs_lon, lat=xs_lat,linestyle='--', color='k')\n", "\n", "################## PLOT VERTICAL CROSS SECTIONS #############\n", "#### longitude slice ####\n", "title=\"Latitude slice\"\n", "ax2 = fig.add_axes(x_cut_panel_axes)\n", "ax2.set_ylim(xsec_vlim) ### Set vert/horiz ranges of XSEC\n", "ax2.set_xlim(xsec_hlim)\n", "displayg.plot_longitude_slice(pltvar,lon=xs_lon,vmin=v0,vmax=v1,cmap=dv_cmap,\n", " colorbar_label=\"m/s\",title=title)\n", " \n", "#### latitude slice ####\n", "title=\"Longitude slice\" \n", "ax3 = fig.add_axes(y_cut_panel_axes)\n", "ax3.set_ylim(xsec_vlim) ### Set vert/horiz ranges of XSEC\n", "ax3.set_xlim(xsec_hlim)\n", "displayg.plot_latitude_slice(pltvar,lat=xs_lat,vmin=v0,vmax=v1,cmap=dv_cmap,\n", " colorbar_label=\"m/s\",title=title)\n", "#%%%%%%%%%" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Save output plot if desired" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "if plttyp != \"screen\" and plttyp != \"Screen\":\n", " ofile=ofile0+\".\"+plttyp\n", " print(\"Saving as... \"+ofile)\n", " fig.savefig(ofile,format=plttyp,orientation='portrait',quality=100)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Plot Gridded Map and Interactive choice of Cross section Lats" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "\n", "def xhair(lt,ln,var,glev,vmin,vmax,proj,cmap,tit,mres,lonvals,latvals,wl,bl):\n", "# Routine to plot grided radar map and longitudingal cross setions\n", "# Uses python slider to move latitudinally to see cross secction evolution\n", " gs = gridspec.GridSpec(1, 1)\n", " gs.update(left=0.0, right=0.45, wspace=0.005)\n", " ax1 = plt.subplot(gs[0],projection=projection)#121) \n", " ax1.clear()\n", " displayg.plot_grid(var, glev, vmin=vmin, vmax=vmax,projection=proj,\n", " cmap=cmap,title=tit,ax=ax1,colorbar_flag=False)\n", " map_add_on(ax1,mres,lonvals,latvals)\n", " ax1.scatter(clong,clatg,color='darkblue',edgecolor='k',linewidth='1',s=30,zorder=500)\n", " displayg.plot_crosshairs(lon=ln,lat=lt,linestyle='--', color='k')\n", " \n", "# Plot Cross section\n", " xsec_vlim = [0,15]\n", " xsec_hlim = [-150.,150.]\n", " gs1 = gridspec.GridSpec(1, 1)\n", " gs1.update(left=0.53 ,right=1.,bottom=0.3,top=0.7, wspace=0.005)\n", " ax2 = plt.subplot(gs1[0])\n", " ax2.clear()\n", " title=\"Longitude slice\"\n", " displayg.plot_latitude_slice(var,lat=lt,vmin=vmin,vmax=vmax,cmap=cmap,\n", " colorbar_label=\"m/s\",title=title)\n", "\n", "plt_gtitle=(nexrad_site+\" \"+grd_res+\" grd \"+\n", " pltvar+\" lev=\"+str( rg.z['data'][glev]/1000.)+\"km \\n\"+sdate[2])\n", "interact(xhair,lt=(blat,ulat,0.2),ln=fixed(xs_lon),var=fixed(pltvar),glev=fixed(glev),\n", " vmin=fixed(v0),vmax=fixed(v1),proj=fixed(projection),cmap=fixed(dv_cmap),\n", " tit=fixed(plt_gtitle),mres=fixed(map_res),lonvals=fixed(lonvals),\n", " latvals=fixed(latvals),wl=fixed(wlon),bl=fixed(blat))\n", "#%%%%%%%%%%%%%%%%" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Plot Gridded Map and Interactive choice of Cross section Lons" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def yhair(lt,ln,var,glev,vmin,vmax,proj,cmap,tit,mres,lonvals,latvals,wl,bl):\n", "# Routine to plot grided radar map and latitudingal cross setions\n", "# Uses python slider to move longitudinally to see cross secction evolution\n", " gs = gridspec.GridSpec(1, 1)\n", " gs.update(left=0.0, right=0.45, wspace=0.005)\n", " ax1 = plt.subplot(gs[0],projection=projection)#121) \n", " ax1.clear()\n", " displayg.plot_grid(var, glev, vmin=vmin, vmax=vmax,projection=proj,\n", " cmap=cmap,title=tit,ax=ax1,colorbar_flag=False)\n", " map_add_on(ax1,mres,lonvals,latvals)\n", " ax1.scatter(clong,clatg,color='darkblue',edgecolor='k',linewidth='1',s=30,zorder=500)\n", " displayg.plot_crosshairs(lon=ln,lat=lt,linestyle='--', color='k')\n", " \n", "# Plot Cross section\n", " xsec_vlim = [0,15]\n", " xsec_hlim = [-150.,150.]\n", " gs1 = gridspec.GridSpec(1, 1)\n", " gs1.update(left=0.53 ,right=1.,bottom=0.3,top=0.7, wspace=0.005)\n", " ax2 = plt.subplot(gs1[0])\n", " ax2.clear()\n", " title=\"Latitude slice\"\n", " displayg.plot_longitude_slice(var,lon=ln,vmin=vmin,vmax=vmax,cmap=cmap,\n", " colorbar_label=\"m/s\",title=title)\n", "\n", "plt_gtitle=(nexrad_site+\" \"+grd_res+\" grd \"+\n", " pltvar+\" lev=\"+str( rg.z['data'][glev]/1000.)+\"km \\n\"+sdate[2])\n", "interact(yhair,lt=fixed(xs_lat),ln=(wlon,elon,0.2),var=fixed(pltvar),glev=fixed(glev),\n", " vmin=fixed(v0),vmax=fixed(v1),proj=fixed(projection),cmap=fixed(dv_cmap),\n", " tit=fixed(plt_gtitle),mres=fixed(map_res),lonvals=fixed(lonvals),\n", " latvals=fixed(latvals),wl=fixed(wlon),bl=fixed(blat))\n", "#%%%%%%%%%%%%%%%%" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Plot Interactive choice of Cross sections Lats+ Lons" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "scrolled": true }, "outputs": [], "source": [ "########### PLOT VERTICAL CROSS SECTIONS ################\n", "def latxsec(ln,var,vmin,vmax,cmap):\n", " displayg.plot_longitude_slice(var,lon=ln,vmin=vmin,vmax=vmax,cmap=cmap,\n", " colorbar_label=\"m/s\",title=title)\n", "\n", "def lonxsec(lt,var,vmin,vmax,cmap):\n", " displayg.plot_latitude_slice(var,lat=lt,vmin=vmin,vmax=vmax,cmap=cmap,\n", " colorbar_label=\"m/s\",title=title)\n", "\n", "#### latitude slice ####\n", "x_cut_panel_axes = [0.05, 0.30, .4, .25]\n", "ax2 = fig.add_axes(x_cut_panel_axes)\n", "ax2.set_ylim(xsec_vlim) ### Set vert/horiz ranges of XSEC\n", "ax2.set_xlim(xsec_hlim)\n", "title=\"Latitude Slice\"\n", "interact(latxsec,ln=(wlon,elon,0.2),var=fixed(pltvar),vmin=fixed(v0),vmax=fixed(v1),\n", " cmap=fixed(dv_cmap))\n", "\n", "#### longitude slice ####\n", "y_cut_panel_axes = [0.55, 0.30, .4, .25]\n", "ax3 = fig.add_axes(y_cut_panel_axes)\n", "ax3.set_ylim(xsec_vlim) ### Set vert/horiz ranges of XSEC\n", "ax3.set_xlim(xsec_hlim)\n", "title=\"Longitude Slice\"\n", "interact(lonxsec,lt=(blat,ulat,0.2),var=fixed(pltvar),vmin=fixed(v0),vmax=fixed(v1),\n", " cmap=fixed(dv_cmap))\n", "#%%%%%%%%%%%" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Will compare actual radat map to gridded map at level closest to sweep \n", "# elevation\n", "map_panel_axes0 = [0.02, 0.03, .5, .5] # sets up panels for plotting \n", "map_panel_axes1 = [0.5, 0.03, .5, .50]\n", "\n", "plt_title=nexrad_site+\" \"+pltvar+\" swp=\"+str(sweep)+\" \"+sdate[2]\n", "fig = plt.figure(figsize=[15,7])\n", "ax0 = fig.add_axes(map_panel_axes0, projection=projection)\n", "display.plot_ppi_map(pltvar,sweep, vmin=v0, vmax=v1,\n", " min_lon=wlon, max_lon=elon, min_lat=blat, max_lat=ulat,\n", " lon_lines=np.arange(int(wlon)-1.0,int(elon)+1.0,grid_line_delta), \n", " lat_lines=np.arange(int(blat)-1.0,math.floor(ulat)+2.,grid_line_delta), \n", " projection=projection,cmap=dv_cmap,resolution=map_res,colorbar_label=\"m/s\",\n", " lat_0=clat,lon_0=clon,title=plt_title)\n", "map_add_on(ax0,map_res,lonvals,latvals)\n", "display.plot_range_rings([50, 100]) #plots rings \n", "ax0.scatter(clon,clat,color='darkblue',edgecolor='k',linewidth='1',s=30,zorder=500)\n", "\n", "########### PLOT GRIDDED RADAR ON MAP ####################\n", "ax1 = fig.add_axes(map_panel_axes1, projection=projection)\n", "plt_gtitle=(nexrad_site+\" \"+grd_res+\" grd \"+\n", " pltvar+\" lev=\"+str( rg.z['data'][glev]/1000.)+\"kkm \\n\"+sdate[2])\n", "displayg.plot_grid(pltvar, glev, vmin=v0, vmax=v1,projection=projection,\n", " cmap=dv_cmap,title=plt_gtitle,ax=ax1,colorbar_label=\"m/s\",axislabels=(' ',' '),\n", " axislabels_flag=True)\n", "map_add_on(ax0,map_res,lonvals,latvals)\n", "### Plots station location ####\n", "ax1.scatter(clon,clat,color='darkblue',edgecolor='k',linewidth='1',s=30,zorder=500)\n", "#%%%%%%%%%" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "##Program End" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "##Program End" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "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.7" } }, "nbformat": 4, "nbformat_minor": 2 }