Requesting Grid Parameters and Levels with python-awips

The Python AWIPS Data Access Framework can be used to query available grid parameters and levels if given a known Grid name (as of AWIPS 15.1.3 we can not query derived parameters, only parameters which have been directly decoded).

This example requests the U and V wind components for the GFS 40km CONUS and plots the wind speed (total vector) as a gridded contour (color-filled isotachs, essentially):

from awips.dataaccess import DataAccessLayer 
# Select grid 
DataAccessLayer.changeEDEXHost('edex-cloud.unidata.ucar.edu') 
request = DataAccessLayer.newDataRequest() 
request.setDatatype('grid') 
request.setLocationNames('RAP13') 

grid is given as the data type, and RAP13 is given as the grid name. The function getAvailableParameters() returns a list of parameters:

# Print parm 
list available_parms = DataAccessLayer.getAvailableParameters(request) 
available_parms.sort() 
for parm in available_parms: 
    print parm 

printed as

 AV BLI CAPE CFRZR6hr CICEP6hr CIn CP6hr CRAIN6hr CSNOW6hr GH P P6hr PMSL PVV PW RH SLI T TP6hr VSS WEASD WGH uW vW 

List Available Levels for Parameter

Selecting the u wind component (uW) first and printing the return from getAvailableLevels():

request.setParameters('uW')
available_levels = DataAccessLayer.getAvailableLevels(request) 
available_levels.sort() 
for level in available_levels: 
    print level 

prints as

 1000.0MB 950.0MB 925.0MB 900.0MB 875.0MB 850.0MB 825.0MB 800.0MB 775.0MB 725.0MB 600.0MB 575.0MB 0.0_30.0BL 60.0_90.0BL 90.0_120.0BL 0.5PV 2.0PV 30.0_60.0BL 1.0PV 750.0MB 120.0_150.0BL 975.0MB 700.0MB 675.0MB 650.0MB 625.0MB 550.0MB 525.0MB 500.0MB 450.0MB 400.0MB 300.0MB 250.0MB 200.0MB 150.0MB 100.0MB 0.0TROP 1.5PV 150.0_180.0BL 350.0MB 10.0FHAG 0.0MAXW 

Construct Wind Field from U and V Components

Calling setLevels("10.0FHAG") and selecting the last available time t[-1]:

import numpy from metpy.units 
import units 
request.setLevels('10.0FHAG') 
t = DataAccessLayer.getAvailableTimes(request) 
# Select last time for u-wind 
response = DataAccessLayer.getGridData(request, [t[-1]]) 
data_uw = response[-1] 
lons,lats = data_uw.getLatLonCoords() 

Requesting the grid for v wind componeents (vW):

request.setParameters('vW') 
response = DataAccessLayer.getGridData(request, [t[-1]]) 
data_vw = response[-1] 
print 'Time :', t[-1] 
print 'Model:', data_vw.getLocationName() 
print 'Unit :', data_vw.getUnit() 
print 'Parms :', data_uw.getParameter(), data_vw.getParameter() 
print data_vw.getRawData().shape 

prints as

 Time : 2016-04-20 18:00:00 (240) 
Model: GFS40 
Unit : m*sec^-1 
arms : vW vW (185, 129) windArray = [[ 1.47078204 1.69705617 0.69296461 ..., 6.98621511 ..., 0.91923875 1.24450791 1.28693426]] 

Calculate absolute wind speed from U and V

The data objects data_uw.getRawData() and data_vw.getRawData()

spd = numpy.sqrt( data_uw.getRawData()**2 + data_vw.getRawData()**2 ) 
spd = spd * units.knot 
print 'windArray =', spd 

Plotting a Grid with Basemap

Using matplotlib, numpy, and basemap:

python 
%matplotlib inline 
import matplotlib.tri as mtri 
import matplotlib.pyplot as plt 
from matplotlib.transforms import offset_copy 
from mpl_toolkits.basemap import Basemap, cm 
import numpy as np 
from numpy import linspace, transpose 
from numpy import meshgrid 

plt.figure(figsize=(12, 12), dpi=100) 
map = Basemap(projection='cyl', resolution = 'c', llcrnrlon = lons.min(), llcrnrlat = lats.min(), urcrnrlon =lons.max(), urcrnrlat = lats.max() ) 
map.drawcoastlines() 
map.drawstates() 
map.drawcountries() 
#  # We have to reproject our grid, see https://stackoverflow.com/questions/31822553/m 
# x = linspace(0, map.urcrnrx, data_uw.getRawData().shape[1]) 

y = linspace(0, map.urcrnry, data_uw.getRawData().shape[0]) 
xx, yy = meshgrid(x, y) 
ngrid = len(x) 
rlons = np.repeat(np.linspace(np.min(lons), np.max(lons), ngrid), ngrid).reshape(ngrid, ngrid) 
rlats = np.repeat(np.linspace(np.min(lats), np.max(lats), ngrid), ngrid).reshape(ngrid, ngrid).T 
tli = mtri.LinearTriInterpolator(mtri.Triangulation(lons.flatten(), lats.flatten()), spd.flatten()) 
rdata = tli(rlons, rlats) 
#cs = map.contourf(rlons, rlats, rdata, latlon=True) 
cs = map.contourf(rlons, rlats, rdata, latlon=True, vmin=0, vmax=20, cmap='BuPu') 
# Add colorbar 
cbar = map.colorbar(cs,location='bottom',pad='5%') 
cbar.set_label('Wind Speed (Knots)') 
# Show plot 
plt.show() 

or use pcolormesh rather than contourf

python 
plt.figure(figsize=(12, 12), dpi=100) 
map = Basemap(projection='cyl', resolution = 'c', llcrnrlon = lons.min(), llcrnrlat = lats.min(), urcrnrlon =lons.max(), urcrnrlat = lats.max() ) 
map.drawcoastlines() 
map.drawstates() 
map.drawcountries() 
cs = map.pcolormesh(rlons, rlats, rdata, latlon=True, vmin=0, vmax=20, cmap='BuPu') 

Plotting a Grid with Cartopy

python 
import os 
import matplotlib.pyplot as plt 
import numpy as np 
import iris 
import cartopy.crs as ccrs 
from cartopy import config 

lon,lat = data.getLatLonCoords() 
plt.figure(figsize=(12, 12), dpi=100) 
ax = plt.axes(projection=ccrs.PlateCarree()) 
cs = plt.contourf(rlons, rlats, rdata, 60, transform=ccrs.PlateCarree(), vmin=0, vmax=20, cmap='BuPu') 
ax.coastlines() 
ax.gridlines() 

# add colorbar 
cbar = plt.colorbar(orientation='horizontal') 
cbar.set_label('Wind Speed (Knots)') 

plt.show()