Pablo Romero wrote:
Jeff,
Here's a link to the website that's creating the basemap plots with
the "H"'s and "L"'s:
http://magicseaweed.com/msw-surf-charts2.php?chart=64&res=750&type=pressure&starttime=
As you can see, this website's "pressure chart" interface is based on
creating individual images for each point in a time series. And, if
you click on any one of the "region" links below the plot, you'll be
presented with the same pressure plot for a sub-region; obviously,
this needs to be based on an automated process, since I seriously
doubt someone is manually creating each plot for each time step for
each "sub-region" in their catalog...bottom line, there MUST be a
systematic way to add those "H" and "L" characters to closed contours
that are local min/max AND that are LARGE enough to hold an "H" or "L"
character (as you mentioned, some of the smaller closed contours dont
display the 'H'/'L' characters).
So the big question is: how could one go about setting something like
this up???
Do you think it would be possible to accomplish this with the
function/tutorial you provided???
If so, could you please help me understand a little more how I would
go about using that function? (sorry, I wasnt very clear on how it
should be applied to this problem).
I know for a fact that this website is using matplotlib/basemap to
create these plots. Theyre using a python interface to
matplolib/basemap named 'pygrads'; it provides basemap/mpl plotting
capabilities to an existing application named GrADS. I just dont know
how the heck they got those "H/L" characters there....
Thanks again for the help,
Pablo
Pablo: I've added a script (attached to this email) called
"plothighsandlows.py" in the basemap examples directory that shows how
to do this. It uses scipy.ndimage.filters, so you'll need to have scipy
installed to run it.
-Jeff
----------------------------------------
Date: Fri, 6 Mar 2009 06:07:10 -0700
From: jsw...@fastmail.fm
To: romero...@hotmail.com
CC: matplotlib-users@lists.sourceforge.net
Subject: Re: [Matplotlib-users] plotting air pressure data with
contour() function
Jeff Whitaker wrote:
Pablo Romero wrote:
Hello,
I have a question about plotting pressure contours with matplotlib.
I've seen other applications using matplotlib where the pressure
contours are drawn with annoted text "H" and "L" characters being
drawn in the centers of closed contours...
i.e., if there is a closed contour line, and its value is over a
certain threshold value, plot an "H" to indicate a "high pressure
zone", else plot a "L" to indicate a "low pressure zone."
This is standard weather map plotting stuff, and Ive definitely
seen other plots produced using matplotlib that inlude these
annoted "H" and "L" characters. In the application Ive seen, the
process was most definitely automated, since it was applied to a
large number of plots (i.e., the "H"'s and "L"'s werent added
manually, since its not feasible). Unfortunately, the creators of
these plots are not willing to share their technique.
I dont know how to set this up with matplotlib.
Does anyone have any experience with this???
Is there any way to identify a "closed contour" & its value from a
"contour class" that is returned from matplotlib's contour()
function???
Pablo: There is no easy way to do this in matplotlib. I guess you
would try to find closed contours with no other contours inside them,
the place the label at the center of that region. This might end up
being quite tricky. I see from your example that there are many such
regions that are not labelled (some are, some aren't).
Or, you might just try to find local minima and maxima in your gridded
data and plot H's and L's there. This ought to be easier.
-Jeff
Pablo: Regarding the latter method, here's a relevant thread from the
scipy list:
http://www.nabble.com/Finding-local-minima-of-greater-than-a-given-depth-td18988309.html
-Jeff
_________________________________________________________________
Windows Live™ Contacts: Organize your contact list.
http://windowslive.com/connect/post/marcusatmicrosoft.spaces.live.com-Blog-cns!503D1D86EBB2B53C!2285.entry?ocid=TXT_TAGLM_WL_UGC_Contacts_032009
"""
plot H's and L's on a sea-level pressure map
(uses scipy.ndimage.filters)
"""
import numpy as np
import matplotlib.pyplot as plt
import sys
from mpl_toolkits.basemap import Basemap, NetCDFFile, addcyclic
from scipy.ndimage.filters import minimum_filter, maximum_filter
def extrema(mat,mode='wrap',window=10):
mn = minimum_filter(mat, size=window, mode=mode)
mx = maximum_filter(mat, size=window, mode=mode)
# (mat == mx) true if pixel is equal to the local max
# The next computation suppresses responses where
# the function is flat.
local_maxima = ((mat == mx) & (mat != mn))
local_minima = ((mat == mn) & (mat != mx))
# Get the indices of the maxima, minima
return np.nonzero(local_minima), np.nonzero(local_maxima)
if len(sys.argv) < 2:
print 'enter date to plot (YYYYMMDDHH) on command line'
raise SystemExit
# get date from command line.
YYYYMMDDHH = sys.argv[1]
# open OpenDAP dataset.
try:
data=NetCDFFile("http://nomad1.ncep.noaa.gov:9090/dods/gdas/rotating/gdas"+YYYYMMDDHH)
except:
data=NetCDFFile("http://nomad1.ncep.noaa.gov:9090/dods/gdas/rotating/"+YYYYMMDDHH[0:6]+"/gdas"+YYYYMMDDHH)
# read lats,lons.
lats = data.variables['lat'][:]
lons1 = data.variables['lon'][:]
nlats = len(lats)
nlons = len(lons1)
# read prmsl, convert to hPa (mb).
prmsl = 0.01*data.variables['prmslmsl'][0]
# the window parameter controls the number of highs and lows detected.
# (higher value, fewer highs and lows)
local_min, local_max = extrema(prmsl, mode='wrap', window=25)
# create Basemap instance.
m =
Basemap(llcrnrlon=0,llcrnrlat=-65,urcrnrlon=360,urcrnrlat=65,projection='merc')
# add wrap-around point in longitude.
prmsl, lons = addcyclic(prmsl, lons1)
# contour levels
clevs = np.arange(900,1100.,5.)
# find x,y of map projection grid.
lons, lats = np.meshgrid(lons, lats)
x, y = m(lons, lats)
# create figure.
fig=plt.figure(figsize=(12,6))
cs = m.contour(x,y,prmsl,clevs,colors='k',linewidths=1.)
m.drawcoastlines(linewidth=1.25)
m.fillcontinents(color='0.8')
m.drawparallels(np.arange(-80,81,20))
m.drawmeridians(np.arange(0,360,60))
xlows = x[local_min]; xhighs = x[local_max]
ylows = y[local_min]; yhighs = y[local_max]
lowvals = prmsl[local_min]; highvals = prmsl[local_max]
# plot lows as blue L's, with min pressure value underneath.
xyplotted = []
# don't plot if there is already a L or H within dmin meters.
dmin = 500000
for x,y,p in zip(xlows, ylows, lowvals):
if x < m.xmax and x > m.xmin and y < m.ymax and y > m.ymin:
dist = [np.sqrt((x-x0)**2+(y-y0)**2) for x0,y0 in xyplotted]
if not dist or min(dist) > dmin:
plt.text(x,y,'L',fontsize=14,fontweight='bold',
horizontalalignment='center',
verticalalignment='center',color='blue')
plt.text(x,y-400000,repr(int(p)),fontsize=9,
horizontalalignment='center',
verticalalignment='top',color='blue')
xyplotted.append((x,y))
# plot highs as red H's, with max pressure value underneath.
xyplotted = []
for x,y,p in zip(xhighs, yhighs, highvals):
if x < m.xmax and x > m.xmin and y < m.ymax and y > m.ymin:
dist = [np.sqrt((x-x0)**2+(y-y0)**2) for x0,y0 in xyplotted]
if not dist or min(dist) > dmin:
plt.text(x,y,'H',fontsize=14,fontweight='bold',
horizontalalignment='center',
verticalalignment='center',color='red')
plt.text(x,y-400000,repr(int(p)),fontsize=9,
horizontalalignment='center',
verticalalignment='top',color='red')
xyplotted.append((x,y))
plt.title('Mean Sea-Level Pressure (with Highs and Lows) %s' % YYYYMMDDHH)
plt.show()
------------------------------------------------------------------------------
Open Source Business Conference (OSBC), March 24-25, 2009, San Francisco, CA
-OSBC tackles the biggest issue in open source: Open Sourcing the Enterprise
-Strategies to boost innovation and cut costs with open source participation
-Receive a $600 discount off the registration fee with the source code: SFAD
http://p.sf.net/sfu/XcvMzF8H
_______________________________________________
Matplotlib-users mailing list
Matplotlib-users@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/matplotlib-users
