plot_mrms-CONUS.py 3.98 KB
Newer Older
Jared Rennie's avatar
Jared Rennie committed
1
#!/usr/bin/python
Jared Rennie's avatar
Jared Rennie committed
2
3
4

# Written by Jared Rennie, NCICS

Jared Rennie's avatar
Jared Rennie committed
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
import matplotlib as mpl
mpl.use('Agg')
import matplotlib.pyplot as plt
import matplotlib.colors as mcolors

import numpy as np
import numpy.ma as ma
import sys
import time
import datetime
import os
import xarray as xr

import cartopy
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import warnings
warnings.filterwarnings("ignore")

from metpy.cbook import get_test_data
from metpy.io import Level2File
from metpy.plots import colortables
from metpy.plots import ctables
from metpy.plots import USCOUNTIES

# Define Directories
main_directory='/store/sfcnet/datasets/mrms'
in_directory=main_directory+'/source'
out_directory=main_directory+'/results_plots'

# Set Up for Plotting
dpi=300
plt.style.use('dark_background')

# Set Bounds
minLat = 22.5    
maxLat = 50.5   
minLon = -120 
maxLon = -73 

#################################################
# BEGIN PROGRAM
start=time.time()

# Read in Arguments 
Jared Rennie's avatar
Jared Rennie committed
50
# DATA SOURCE: https://mrms.ncep.noaa.gov/data/2D/MergedReflectivityQComposite/
Jared Rennie's avatar
Jared Rennie committed
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
if len(sys.argv) < 2:
    sys.exit("USAGE: <INFILE>")  
file = sys.argv[1]
print(file)

# Get Time
current_time=file[40:55]
if time.localtime().tm_isdst == 1: # DST in Affect (EDT)
    utcAdjust=4
    timeZone='EDT'
else:
    utcAdjust=5 # DST NOT in Affect (EST)
    timeZone='EST'

current_date=datetime.datetime(int(current_time[0:4]),int(current_time[4:6]),int(current_time[6:8]),int(current_time[9:11]),int(current_time[11:13])) - datetime.timedelta(hours=utcAdjust)  
current_day= "%04i%02i%02i" % (current_date.year,current_date.month,current_date.day)                                   
current_hour= "%02i:%02i" % (current_date.hour,current_date.minute) 

# Read in Data
inData=xr.load_dataset(in_directory+'/'+file,engine='cfgrib')
lats=inData.latitude.values
lons=inData.longitude.values-360
vals=inData.unknown.values
vals_nonmiss=ma.masked_where(vals <= 0, vals)

print("PLOT")
unit='dBz'
color_map='Greens' 
bounds=np.arange(-10,120,10)

ref_norm, ref_cmap = ctables.registry.get_with_steps('NWSReflectivity', 5, 5)

# Set Up Figure
#vmin=np.min(bounds)
#vmax=np.max(bounds)
extend='max'
#orm = mcolors.BoundaryNorm(boundaries=bounds, ncolors=256)
cm = plt.cm.ScalarMappable(cmap=ref_cmap)
#cm.set_array(np.array(vals_nonmiss))
#cm.set_clim(vmin, vmax)

fig= plt.figure(num=1, figsize=(8,5), dpi=dpi, facecolor='w', edgecolor='k')
ax = fig.add_axes([0, 0, 1, 1], projection=ccrs.LambertConformal())
ax.set_extent([minLon, maxLon, minLat, maxLat], crs=ccrs.Geodetic())
ax.set_facecolor('black')
ax.add_feature(cfeature.COASTLINE,linewidth=0.5,edgecolor='white',facecolor='black')
ax.add_feature(cfeature.OCEAN,linewidth=0.5,edgecolor='white',facecolor='black')
ax.add_feature(cfeature.BORDERS,linewidth=0.5,edgecolor='white',facecolor='black')
ax.add_feature(cfeature.STATES,linewidth=0.5,edgecolor='white',facecolor='#4a4a4a')

# Plot Data
ax.pcolormesh(lons, lats,vals_nonmiss,transform=ccrs.PlateCarree(),cmap=ref_cmap,norm=ref_norm,zorder=10)
value_max="%5.2f" % np.nanmax(vals_nonmiss)
value_min="%5.2f" % np.nanmin(vals_nonmiss)

# Add Colormap
#cax = fig.add_axes([0.1, -0.035, 0.8, 0.03])
#cbar=plt.colorbar(cm, cax=cax,boundaries=bounds,orientation='horizontal',extend=extend,spacing='uniform')
#cbar=plt.colorbar(cm, cax=cax,orientation='horizontal',extend=extend,spacing='uniform')
#cbar.ax.tick_params(labelsize=15)
#cbar.set_label(unit,size=15)

# Add Titles
plt.annotate('Time='+str(current_day)+' '+str(current_hour)+' '+str(timeZone)+'\nMade By Jared Rennie (@jjrennie)',xy=(0.02, 0.02), xycoords='axes fraction', fontsize=5,backgroundcolor='black',color='white',horizontalalignment='left', verticalalignment='bottom',zorder=4)
plt.suptitle('MRMS Composite Radar',size=17,color='white',y=1.05) 

# Save Figure
plt.savefig(out_directory+"/MRMS_CONUS_"+str(current_day)+"-"+str(current_hour.replace(':',''))+".png",bbox_inches='tight') 
plt.clf()
plt.close()

####################
# DONE
####################
print("DONE!")
end=time.time()
print ("Runtime: %8.1f seconds." % (end-start))
sys.exit()