GMTB Workflow Documentation
GMTB_graphics.py
Go to the documentation of this file.
1 #!/usr/bin/env python
2 ##\file GMTB_graphics.py
3 
4 import matplotlib as mpl
5 mpl.use('Agg')
6 import pygrib
7 from mpl_toolkits.basemap import Basemap
8 import numpy as np
9 import matplotlib.pyplot as plt
10 import datetime
11 import os,sys,math
12 import maps, regrib
13 import get_grib_field as ggf
14 import get_spd_field as gsf
15 
16 ##\ingroup graphics
17 # This function defines a variable table to plot.
18 #\param lev The value of grb.level
19 #\param pvar The value of grb.shortName
20 #\param barbs logical to plot barbs
21 #\param hgt
22 #\param vc The value for grb.typeOfLevel
23 #\param acc_field
25 
26  map_confs = {
27  '200_wind': {'lev': '200', 'pvar': 'wind', 'barbs': True, 'hgt': False, 'vc': 'isobaricInhPa'},
28  '200_temp': {'lev': '200', 'pvar': 't', 'barbs': True, 'hgt': False, 'vc': 'isobaricInhPa'},
29  '500_hgt': {'lev': '500', 'pvar': 'gh', 'barbs': True, 'hgt': True, 'vc': 'isobaricInhPa'},
30  '500_temp': {'lev': '500', 'pvar': 't', 'barbs': True, 'hgt': False, 'vc': 'isobaricInhPa'},
31  '500_vort': {'lev': '500', 'pvar': 'absv', 'barbs': False, 'hgt': False, 'vc': 'isobaricInhPa'},
32  '700_temp': {'lev': '700', 'pvar': 't', 'barbs': True, 'hgt': False, 'vc': 'isobaricInhPa'},
33  '700_vvel': {'lev': '700', 'pvar': 'w', 'barbs': False, 'hgt': False, 'vc': 'isobaricInhPa'},
34  '850_hgt': {'lev': '850', 'pvar': 'gh', 'barbs': True, 'hgt': True, 'vc': 'isobaricInhPa'},
35  '850_temp': {'lev': '850', 'pvar': 't', 'barbs': True, 'hgt': False, 'vc': 'isobaricInhPa'},
36  '850_rh': {'lev': '850', 'pvar': 'r', 'barbs': True, 'hgt': False, 'vc': 'isobaricInhPa'},
37  '2m_tmp': {'lev': '2', 'pvar': '2t', 'barbs': False, 'hgt': False, 'vc': 'heightAboveGround'},
38  '2m_dpt': {'lev': '2', 'pvar': '2d', 'barbs': False, 'hgt': False, 'vc': 'heightAboveGround'},
39  '10m_u': {'lev': '10', 'pvar': '10u', 'barbs': False, 'hgt': False, 'vc': 'heightAboveGround'},
40  '10m_v': {'lev': '10', 'pvar': '10v', 'barbs': False, 'hgt': False, 'vc': 'heightAboveGround'},
41  '10m_wind': {'lev': '10', 'pvar':'10wind', 'barbs': True, 'hgt': False, 'vc': 'heightAboveGround'},
42  'prmsl': {'lev': '0', 'pvar': 'prmsl', 'pline': True, 'barbs': True, 'hgt': False, 'vc': 'meanSea'},
43  'pwat': {'lev': '0', 'pvar': 'pwat', 'pline': False, 'barbs': False, 'hgt': False, 'vc': 'unknown'},
44  'prate': {'lev': '0', 'pvar': 'prate','pline': False, 'barbs': False, 'hgt': False, 'vc': 'surface', 'acc_field': True},
45  'cprat': {'lev': '0', 'pvar': 'cprat','pline': False, 'barbs': False, 'hgt': False, 'vc': 'surface', 'acc_field': True},
46  '6h_acpcp': {'lev': '0', 'pvar': 'acpcp', 'barbs': False, 'hgt': False, 'vc': 'surface', 'acc_field': True},
47  '6h_apcp': {'lev': '0', 'pvar': 'tp', 'barbs': False, 'hgt': False, 'vc': 'surface', 'acc_field': True}
48  }
49 
50  return map_confs
51 
52 
53 ##\defgroup graphics graphics
54 def main():
55 
56  ## Get the cycle information from env vars
57  moad_dataroot=os.getenv('MOAD_DATAROOT')
58  rotdir=os.getenv('ROTDIR')
59  date=os.getenv('START_TIME') # @Y@m@d@H
60  fhr=os.getenv('FCST_TIME') # 3 digit fcst hour
61  grid_list=os.getenv('GRID_VX_LIST')
62  vx_grid=grid_list.split(" ")
63  filename=''.join(['pgrbq', fhr, '.gfs.', date, '.grib2'])
64  gribin='/'.join([rotdir, filename])
65  user_def = os.getenv('USER_MAPS',False)
66 
67  ## Set the output file location (template set in maps.py).
68  outdir='/'.join([moad_dataroot, 'figprd'])
69 
70 
71  # Change projection if needed.
72  for grid in vx_grid:
73  outmap = grid
74  print grid
75 
76  # Wind rotation is a user choice, not necessarily a predefined choice, so is included here.
77  wind_rot = {'G104': 'grid', 'G218': 'grid', 'G3': 'grid'}
78 
79  # Assumes that NCEP pre-defined grids will be defined with a leading "G"
80  if grid[0] == 'G':
81  newgrib = ''.join(['pgrbq', fhr, '.gfs.', date,'_', grid, '.grib2'])
82  newgrib = '/'.join([outdir,newgrib])
83  if not os.path.exists(outdir):
84  os.makedirs(outdir)
85  vx = regrib.new_grid(grid,gribin,newgrib,wind_rot.get(grid,None))
86  vx.regrid()
87  grib = newgrib
88  else:
89  grib = gribin
90 
91  plot_wind=True
92 #z print gribin, outdir
93 
94 
95  # Set forecast hour
96  pvars=['gh', 't']
97  plevels=[ 850, 700, 500, 200 ]
98  coord = [ "isobaricInhPa", "isobaricInhPa", "isobaricInhPa", "isobaricInhPa" ]
99  # Variables include: t, r, u, v, w, gh, absv, clwmr,soilw, 2t, 2d, pwat,
100  # cwat,press,q,pt,tp, acpcp
101  if user_def:
102  for var in pvars:
103  for level in plevels:
104 
105  if plot_wind and var == 'wind':
106  u = ggf.get_ua_field(grib,level,coord,shortName='u')
107  v = ggf.get_ua_field(grib,level,coord,shortName='v')
108  winds=[u,v]
109  elif plot_wind and var == '10wind':
110  u = ggf.get_ua_field(grib,level,coord,shortName='10u')
111  v = ggf.get_ua_field(grib,level,coord,shortName='10v')
112  winds=[u,v]
113  else:
114  winds=None
115 
116  # Get the variable from file by passing key,value pairs from the grib file. If you do not want to
117  # use a particular field to search for the grib record, "None" is also acceptable.
118  if var != 'wind':
119  field=ggf.get_ua_field(grib,level,coord,shortname=var)
120  else:
121  field=u
122 
123 
124  mymap=maps.make_basemap(field,date,var,level,
125  winds=winds, fhr=fhr)
126  mymap.set_bm()
127  mymap.run()
128  mymap.save_figure(outdir)
129  else:
130  def_maps=set_def_maps()
131  for fig in def_maps:
132  # Optional if statement here can control evaluating one type of figure at a time.
133  #if fig == '6h_apcp' or fig == '6h_acpcp':
134  acc_field=def_maps[fig].get('acc_field',False)
135  if (acc_field and fhr != '000') or (not acc_field):
136  level = def_maps[fig].get('lev')
137  var = def_maps[fig].get('pvar')
138  hgt = def_maps[fig].get('hgt')
139  plot_wind = def_maps[fig].get('barbs',False)
140  plot_line = def_maps[fig].get('pline',False)
141  coord = def_maps[fig].get('vc')
142 #z print level,var, hgt,plot_wind, hgt
143 
144  # Grab winds from file if plotting the barbs
145  if plot_wind and var == 'wind':
146  u = ggf.get_ua_field(grib,level,coord,shortName='u')
147  v = ggf.get_ua_field(grib,level,coord,shortName='v')
148  winds=[u,v]
149  elif plot_wind and var == '10wind':
150  u = ggf.get_ua_field(grib,level,coord,shortName='10u')
151  v = ggf.get_ua_field(grib,level,coord,shortName='10v')
152  winds=[u,v]
153  else:
154  winds=None
155 
156  height = None
157  if hgt:
158  height = ggf.get_ua_field(grib,level,coord,shortName='gh')
159 
160  # Get the variable from file
161  if var != 'wind' and var != '10wind':
162  field=ggf.get_ua_field(grib,level,coord,shortName=var)
163  mymap=maps.make_basemap(field,date,var,level,area_flag=outmap, line=plot_line,
164  winds=winds, plot_height=height,def_maps=def_maps)
165  mymap.set_bm()
166  mymap.run()
167  mymap.save_figure(outdir)
168 
169  else:
170  # Calculate and plot wind speed on level
171  spd=gsf.get_spd_field(grib,level,coord,outdir,outmap,fhr,var)
172 
173 
174 main()
def set_def_maps()
This function defines a variable table to plot.