本文整理汇总了Python中mpl_toolkits.basemap.Basemap.makegrid方法的典型用法代码示例。如果您正苦于以下问题：Python Basemap.makegrid方法的具体用法？Python Basemap.makegrid怎么用？Python Basemap.makegrid使用的例子？那么恭喜您, 这里精选的方法代码示例或许可以为您提供帮助。您也可以进一步了解该方法所在类mpl_toolkits.basemap.Basemap的用法示例。
在下文中一共展示了Basemap.makegrid方法的15个代码示例，这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞，您的评价将有助于系统推荐出更棒的Python代码示例。
示例1: ray_density
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
def ray_density(lat1, lon1, lat2, lon2,
dt=1, gr_x=360, gr_y=180,
npts=180, projection='robin',
ray_coverage=False):
'''
Create the DATA array which contains the
info for ray density
'''
global long_0
mymap = Basemap(projection=projection, lon_0=long_0, lat_0=0)
#npts=max(gr_x, gr_y)
# grd[2]: longitude
# grd[3]: latitude
grd = mymap.makegrid(gr_x, gr_y, returnxy=True)
lons, lats = mymap.gcpoints(lon1, lat1, lon2, lat2, npts)
dist = locations2degrees(lat1, lon1, lat2, lon2)
bap = int((dist - 97.0)*npts/dist)/2
midlon = len(lons)/2
midlat = len(lats)/2
lons = lons[midlon-bap:midlon+1+bap]
lats = lats[midlat-bap:midlat+1+bap]
data = np.zeros([len(grd[2]), len(grd[3])])
for i in range(len(lons)):
xi, yi = point_finder(lons[i], lats[i], grd)
# first one is latitude and second longitude
try:
#data[yi][xi] = dt/float(dist-97.0)
data[yi][xi] += dt/len(lons)
except Exception, e:
print e示例2: plot_hmap_ortho
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
def plot_hmap_ortho(h, cmap='jet', mode='log', mx=None, drng=None,
res=0.25, verbose=False):
m = Basemap(projection='ortho',lat_0=90,lon_0=180,rsphere=1.)
if verbose:
print 'SCHEME:', h.scheme()
print 'NSIDE:', h.nside()
lons,lats,x,y = m.makegrid(360/res,180/res, returnxy=True)
lons = 360 - lons
lats *= a.img.deg2rad; lons *= a.img.deg2rad
y,x,z = a.coord.radec2eq(n.array([lons.flatten(), lats.flatten()]))
ax,ay,az = a.coord.latlong2xyz(n.array([0,0]))
data = h[x,y,z]
data.shape = lats.shape
data /= h[0,0,1]
#data = data**2 # only if a voltage beam
data = data_mode(data, mode)
m.drawmapboundary()
m.drawmeridians(n.arange(0, 360, 30))
m.drawparallels(n.arange(0, 90, 10))
if mx is None: mx = data.max()
if drng is None:
mn = data.min()
# if min < (max - 10): min = max-10
else: mn = mx - drng
step = (mx - mn) / 10
levels = n.arange(mn-step, mx+step, step)
return m.imshow(data, vmax=mx, vmin=mn, cmap=cmap)示例3: plot
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
def plot(mp):
from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
map = Basemap(projection='mill', llcrnrlon=0,llcrnrlat=-90,urcrnrlon=360,urcrnrlat=90)
map.drawcoastlines(linewidth=0.25)
map.drawmeridians(numpy.arange(0,360,30))
map.drawparallels(numpy.arange(-90,90,30))
data = mp.data
lons, lats = map.makegrid(mp.width, mp.height)
x, y = map(*numpy.meshgrid(mp.longitudes, mp.latitudes))
#clevs = range(200, 325, 5)
clevs = list(frange(0, 8, 0.25))
cs = map.contourf(x, y, data, clevs, cmap=plt.cm.jet)
cbar = map.colorbar(cs, location='bottom', pad="5%")
cbar.set_label('K')
lon, lat = 174.7772, -41.2889
xpt,ypt = map(lon,lat)
map.plot(xpt,ypt,'bo')
plt.title(mp.name)
plt.gcf().set_size_inches(10,10)
plt.savefig(mp.name + '.png',dpi=100)示例4: getMDAreas
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
def getMDAreas(self,dx=20000,dy=20000,
llcrnrlon=-119.2,
llcrnrlat=23.15,
urcrnrlon=-65.68,
urcrnrlat=48.7):
bmap = Basemap(projection="lcc",
llcrnrlon=llcrnrlon,
llcrnrlat=llcrnrlat,
urcrnrlon=urcrnrlon,
urcrnrlat=urcrnrlat,
resolution='l',
lat_0=38.5,
lat_1=38.5,
lon_0=-97.0)
from matplotlib.nxutils import points_inside_poly
xs = np.arange(bmap.llcrnrx,bmap.urcrnrx + dx,dx)
ys = np.arange(bmap.llcrnry,bmap.urcrnry + dy,dy)
lon,lat,x_grid,y_grid = bmap.makegrid(xs.shape[0],ys.shape[0],returnxy=True)
x, y = x_grid.flatten(), y_grid.flatten()
points = np.vstack((x,y)).T
nx = xs.shape[0]
ny = ys.shape[0]
areas = np.zeros((self.data.shape[0],))
for i in xrange(self.data.shape[0]):
md_x,md_y = bmap(self.data['Lon'][i],self.data['Lat'][i])
poly_xy = np.vstack((md_x,md_y)).T
areas[i] = np.nonzero(points_inside_poly(points,poly_xy))[0].shape[0] * dx * dy / 1000**2
return areas示例5: calculator
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
def calculator(DATA, passed_staev, gr_x, npts, start, end,
projection='robin', ray_coverage=False):
global long_0
mymap = Basemap(projection=projection, lon_0=long_0, lat_0=0)
nonzero = []
gr_y = gr_x
grd = mymap.makegrid(gr_x, gr_y, returnxy=True)
for i in range(start, end):
print i,
sys.stdout.flush()
data, exist_flag = ray_density(passed_staev[i][4], passed_staev[i][5],
passed_staev[i][0], passed_staev[i][1],
dt=passed_staev[i][2],
gr_x=gr_x, gr_y=gr_y, npts=npts,
projection=projection,
ray_coverage=ray_coverage)
if not i == end-1:
if not exist_flag:
continue
if DATA is None: DATA = data.copy()
else: DATA += data
nonzero_tmp = np.nonzero(data)
for j in range(len(nonzero_tmp[0])):
nonzero.append((nonzero_tmp[0][j], nonzero_tmp[1][j]))
fi = open('MAP_OUTPUT/DATA-' + str(start), 'w')
pickle.dump(DATA, fi)
fi.close()
fi = open('MAP_OUTPUT/nonzero-' + str(start), 'w')
pickle.dump(nonzero, fi)
fi.close()示例6: GridPlot
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
class GridPlot(Argument, Settings):
"""ROC绘图"""
parameters = dict()
# 命令行参数,参考 https://docs.python.org/2/howto/argparse.html
command_args = [
[('npy',), {'help': u'指定一个csv2npy后生成的.npy文件', 'type': str, 'nargs': 1}],
[('-n','--name',), {'help': u'输出结果文件名', 'type': str, 'nargs': 1}],
[('-t','--title',), {'help': 'plot image title', 'type': str, 'nargs': 1}],
[('-v','--verbose',), {'help': u'输出详细信息', 'action':"store_true"}]
]
def __init__(self):
"""init"""
Argument.__init__(self)
Settings.__init__(self)
self.parse_args()
def parse_args(self):
"""parse arguments"""
if self.args.npy and os.path.isfile(self.args.npy[0]):
self.parameters['npy'] = self.args.npy[0]
else:
print('{} not found!'.format(self.args.npy[0]))
exit(-1)
self.parameters['name'] = self.args.name[0] if self.args.name else tf.mktemp(suffix='.png',dir='.')
self.parameters['title'] = self.args.title[0] if self.args.title else 'Grid Plot'
if self.args.verbose:
print("============parameters===========")
for key in self.parameters:
print(key,self.parameters[key])
self.process()
def process(self):
""" process """
##load data
data = np.load(self.parameters['npy'])
##create basemap
self.bm = Basemap(projection='cyl',resolution='l',lon_0=120)
self.bm.drawcoastlines(linewidth=0.25)
self.bm.drawcountries(linewidth=0.25)
#self.bm.fillcontinents(color='grey')
lons,lats = self.bm.makegrid(360,181)
x,y = self.bm(lons,lats)
self.bm.contourf(x,y,data)
##add colorbar
self.bm.colorbar(location='bottom',size='5%',label="mm")
##add plot title
plt.title(self.parameters['title'])
##save plot
plt.savefig(self.parameters['name'])示例7: setMap
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
def setMap(lats, lons, data, showstations):
x = np.array(lons)
y = np.array(lats)
z = np.array(data)
# corners of the map
lon_min = np.min(x) - 1
lon_max = np.max(x) + 1
lat_min = np.min(y) - 1
lat_max = np.max(y) + 1
if individual:
data_min = np.min(z)
data_max = np.max(z)
else:
data_min = data_min_fix
data_max = data_max_fix
# general map
m = Basemap(
projection="merc",
llcrnrlat=lat_min,
urcrnrlat=lat_max,
llcrnrlon=lon_min,
urcrnrlon=lon_max,
resolution=map_res,
area_thresh=100,
)
m.drawcoastlines()
m.drawstates()
m.drawcountries()
if not showstations:
y_inc = (lat_max - lat_min) / grd_res
x_inc = (lon_max - lon_min) / grd_res
y_steps = np.linspace(lat_min, lat_max + grd_res, y_inc)
x_steps = np.linspace(lon_min, lon_max + grd_res, x_inc)
x_steps, y_steps = np.meshgrid(x_steps, y_steps)
zgrd = matplotlib.mlab.griddata(x, y, z, x_steps, y_steps)
xgrd, ygrd = m.makegrid(zgrd.shape[1], zgrd.shape[0])
# xgrd,ygrd= np.meshgrid(x_steps,y_steps)
# zgrd = scipy.interpolate.griddata((x, y), z, (xgrd, ygrd), method='nearest')
xgrd, ygrd = m(xgrd, ygrd)
m.contourf(xgrd, ygrd, zgrd, cmap=plt.cm.hot_r, vmin=data_min, vmax=data_max)
else:
# measurement stations
xpts, ypts = m(lons, lats)
m.plot(xpts, ypts, "bo")示例8: get_SAL_native_grid
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
def get_SAL_native_grid():
m = Basemap(projection='merc',llcrnrlat=34.0,llcrnrlon=-110.0,
urcrnrlat=45.0,urcrnrlon=-85.0,lat_ts=39.5,)
# m = Basemap(projection='lcc',llcrnrlat=34.0,llcrnrlon=-110.0,
# urcrnrlat=45.0,urcrnrlon=-85.0,lat_0=30.0,lat_1=50.0,lon_0=-95.0)
lons, lats, xx, yy = m.makegrid(537,308,returnxy=True)
# print(xx[5,-1] - xx[5,-2])
# print(xx[5,1] - xx[5,0])
# print(yy[-1,5] - yy[-2,5])
# print(yy[1,5] - yy[0,5])
# import pdb; pdb.set_trace()
# 4 km spacing
return m, lons, lats, xx[0,:], yy[:,0]示例9: ray_density
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
def ray_density(lat1, lon1, lat2, lon2, dt=1, gr_x=360, gr_y=180, npts=180, projection='robin', ray_coverage=False):
"""
Create the DATA array which contains the info for ray density
Procedure:
1. make a grid based on the inputs (grd)
grd: lon, lat, x, y
2. find the great circle points:
note that lon , lat are actually x and y!
3. calculate the distance and find the middle point
4. subtracting 97 degrees from the distance and find all the points on that section
5. data ---> zero array with x*y elements
"""
global long_0
mymap = Basemap(projection=projection, lon_0=long_0, lat_0=0)
#npts=max(gr_x, gr_y)
# grd[2]: longitude
# grd[3]: latitude
grd = mymap.makegrid(gr_x, gr_y, returnxy=True)
lons, lats = mymap.gcpoints(lon1, lat1, lon2, lat2, npts)
dist = locations2degrees(lat1, lon1, lat2, lon2)
# npts points on dist...how many on (dist-97)!: (dist-97)*npts/dist....but we also need to make it half!
bap = int((dist - 97.0)*npts/dist)/2
midlon = len(lons)/2
midlat = len(lats)/2
lons = lons[midlon-bap:midlon+1+bap]
lats = lats[midlat-bap:midlat+1+bap]
data = np.zeros([len(grd[2]), len(grd[3])])
if not len(lons) == len(lats):
sys.exit('ERROR: Lengths longitudes and latitudes are not the same! %s and %s' % (len(lons), len(lats)))
for i in range(len(lons)):
xi, yi = point_finder(lons[i], lats[i], grd)
# first one is latitude and second longitude
try:
#data[yi][xi] = dt/float(dist-97.0)
data[yi][xi] += dt/len(lons)
except Exception, e:
print '\nException: %s' % e示例10: basemap_lonlat_grid
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
def basemap_lonlat_grid(x,y):
"""Insert array or list of x and y from slab data, will return
proper grids for basemap ('cylindrical' projection) and projection
INPUT: x, y are arrays of coordinates
OUTPUT: x, y, m where x, y are grids in the projection 'm'"""
urcrnrlon, urcrnrlat = max(x), max(y)
llcrnrlon, llcrnrlat = min(x), min(y)
m = Basemap(projection='cyl', llcrnrlon=llcrnrlon, llcrnrlat=llcrnrlat, urcrnrlon=urcrnrlon,\
urcrnrlat=urcrnrlat, resolution='h')
ny = y.size
nx = x.size
lons, lats = m.makegrid(nx,ny)
x, y = m(lons, lats)
return x, y, m示例11: plot_observation_frequency
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
def plot_observation_frequency(locations,SEASONS,config):
for year in range(config['START_YEAR'],config['END_YEAR']+1):
for season in SEASONS:
wanted=SEASONS[season]
latitude = np.asarray(locations['LATITUDE'])
longitude = np.asarray(locations['LONGITUDE'])
Yearly_Data=(locations.loc[locations['YEAR']==year])
Seasonal_Data=(Yearly_Data.loc[Yearly_Data['MONTH'].isin(wanted)])
lats = np.asarray(Seasonal_Data['LATITUDE'])
lons = np.asarray(Seasonal_Data['LONGITUDE'])
Species_count=np.asarray(Seasonal_Data[config['SPECIES']])
Species_count=np.reshape(Species_count,len(Species_count))
lat_min = min(lats)
lat_max = max(lats)
lon_min = min(lons)
lon_max = max(lons)
spatial_resolution = 1
fig = plt.figure()
x = np.array(lons)
y = np.array(lats)
z = np.array(Species_count)
xinum = (lon_max - lon_min) / spatial_resolution
yinum = (lat_max - lat_min) / spatial_resolution
xi = np.linspace(lon_min, lon_max + spatial_resolution, xinum)
yi = np.linspace(lat_min, lat_max + spatial_resolution, yinum)
xi, yi = np.meshgrid(xi, yi)
zi = griddata(x, y, z, xi, yi, interp='linear')
m = Basemap(projection = 'merc',llcrnrlat=lat_min, urcrnrlat=lat_max,llcrnrlon=lon_min, urcrnrlon=lon_max,rsphere=6371200., resolution='l', area_thresh=10000)
m.drawcoastlines()
m.drawstates()
m.drawcountries()
m.drawparallels(np.arange(lat_min,lat_max,config['GRID_SIZE']),labels=[False,True,True,False])
m.drawmeridians(np.arange(lon_min,lon_max,config['GRID_SIZE']),labels=[True,False,False,True])
lat, lon = m.makegrid(zi.shape[1], zi.shape[0])
x,y = m(lat, lon)
z=zi.reshape(xi.shape)
levels=np.linspace(0,z.max(),25)
cm=plt.contourf(x, y, zi,levels=levels,cmap=plt.cm.Greys)
plt.colorbar()
plt.title(config['SPECIES']+"-"+str(year)+"-"+str(season))
#plt.show()
plt.savefig(config['SPECIES']+"-"+str(year)+"-"+str(season)+".png")
plt.close()
return示例12: PPI
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
def PPI(self, rs, data, elev, moment, figname):
figname = self.dirname+figname
data = rs.cartesian(data[:,:,elev], elev)
fig = plt.figure(figsize=(8,6))
ax = fig.add_axes([0.1,0.1,0.8,0.8])
latmin, latmax, lonmin, lonmax = rs.radar_limits(elev)
m = Basemap(projection='cyl', lon_0=rs.radar_lon, lat_0=rs.radar_lat,
llcrnrlat=latmin, urcrnrlat=latmax,
llcrnrlon=lonmin, urcrnrlon=lonmax,
resolution='h', suppress_ticks=True)
#m.drawcoastlines(color='0', linewidth=1)
#m.drawstates(color='0', linewidth=1)
#m.drawcountries(color='0', linewidth=1)
#m.drawrivers(color='blue')
#m.fillcontinents(color='#cc9955', lake_color='aqua', zorder = 0) #cor do continente
#m.drawmapboundary(fill_color='aqua') #cor do mar
parallels = np.arange(-90,0,1.)
m.drawparallels(parallels,labels=[1,0,0,0],fontsize=10, linewidth=0.0)
meridians = np.arange(180.,360.,1.)
m.drawmeridians(meridians,labels=[0,0,0,1],fontsize=10, linewidth=0.0)
numcols, numrows = data.shape
lons, lats = m.makegrid(numcols, numrows) # get lat/lons of ny by nx evenly space grid.
x, y = m(lons, lats) # compute map proj coordinates.
cmap, clevs, unit, title = cm.get_info(moment)
norm = mpl.colors.BoundaryNorm(clevs, cmap.N)
contour = m.contourf(x, y, data, clevs, cmap=cmap, norm=norm)
cbar = m.colorbar(contour, cmap=cmap, norm=norm, spacing='uniform', location='right', pad='2%', size='5%')
cbar.set_label(unit, rotation='horizontal')
plt.title(title+u' (%.1f°)\n%s' % (rs.fixed_angle[elev], rs.date))
self.draw_circle(rs, m, elev) #desenha o círculo do raio do radar
plt.savefig(figname)示例13: soda_temp_plot
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
def soda_temp_plot(file_name,t,d):
file = '/Users/yuewang/Documents/study/courses/OCVN689/project/'+ str(file_name)
nc = netCDF4.Dataset(file)
dates=netCDF4.num2date(nc.variables['TIME1'][:],'seconds since 1916-01-02 12:00:00')
fig = plt.figure(figsize=(20,10))
ax = fig.add_axes([0.1,0.1,0.8,0.8])
temp = nc.variables['TEMP']
temp_0 = temp[t,d,:,:]
lon = nc.variables['LON241_580'][:]
lat = nc.variables['LAT142_161'][:]
DEPTH = nc.variables['DEPTH1_4']
m = Basemap(llcrnrlat=lat[0],urcrnrlat=lat[-1],\
llcrnrlon=lon[0],urcrnrlon=lon[-1],\
projection='mill',resolution = 'h',ax=ax)
parallels = np.arange(-5.,5.,2.)
m.drawparallels(parallels,labels=[1,0,0,0],fontsize=10)
meridians = np.arange(120.,300.,30.)
m.drawmeridians(meridians,labels=[0,0,0,1],fontsize=10)
ny = temp_0.shape[0]; nx =temp_0.shape[1]
lons, lats = m.makegrid(nx, ny)
x, y = m(lons, lats)
cs = m.contourf(x,y,temp_0,cmap=cm.sstanom)
m.drawcoastlines()
cbar = m.colorbar(cs,location='bottom', size="15%", pad='35%')
cbar.set_label('Temperature(deg.C)')
ax.set_title('Sea Temperature at ' +str(dates[t])+ 'at depth of'+ str(DEPTH[d])+' m')
plt.show()示例14: interpolated_color_map
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
def interpolated_color_map(lats, lons, values, spatial_resolution=0.1, interp='nn', cmap=None):#cm.s3pcpn):
lat_0 = 51
lat_min = 47
lat_max = 55
lon_0 = 10
lon_min = 5
lon_max = 16
m = Basemap(projection='tmerc', lat_0=lat_0, lon_0=lon_0,
llcrnrlat=lat_min, llcrnrlon=lon_min, urcrnrlat=lat_max,
urcrnrlon=lon_max, resolution='i')
m.drawcoastlines()
m.drawcountries()
m.drawmapboundary()
lats = np.array(lats)
lons = np.array(lons)
values = np.array(values)
lat_inum = (lat_max - lat_min) / spatial_resolution
lon_inum = (lon_max - lon_min) / spatial_resolution
# xi = np.linspace(lat_min, lat_max + spatial_resolution, xinum)
# yi = np.linspace(lon_min, lon_max + spatial_resolution, yinum)
lat_i = np.linspace(lat_min, lat_max, lat_inum)
lon_i = np.linspace(lon_min, lon_max, lon_inum)
lat_i, lon_i = np.meshgrid(lat_i, lon_i)
value_i = griddata(lats, lons, values, lat_i, lon_i, interp=interp)
lat_grid, lon_grid = m.makegrid(value_i.shape[1], value_i.shape[0])
x_grid, y_grid = m(lat_grid, lon_grid)
m.contourf(x_grid, y_grid, value_i, cmap=cmap)
m.scatter(lats, lons, color='k', s=50, latlon=True)
plt.show()示例15: m
# 需要导入模块: from mpl_toolkits.basemap import Basemap [as 别名]
# 或者: from mpl_toolkits.basemap.Basemap import makegrid [as 别名]
llcrnrlat=latcorners[0],urcrnrlat=latcorners[2],\
llcrnrlon=loncorners[0],urcrnrlon=loncorners[2],\
rsphere=6371200.,resolution='l',area_thresh=1000)
# create figure
fig = plt.figure(figsize=(6,5))
ax = plt.gca()
# draw coastlines, state and country boundaries, edge of map.
m.drawcoastlines()
m.drawstates()
m.drawcountries()
# project data
ny = data.shape[0]; nx = data.shape[1]
lons, lats = m.makegrid(nx, ny) # get lat/lons of ny by nx evenly space grid.
x, y = m(lons, lats) # compute map proj coordinates.
# draw filled contours.
clevs = np.array([0,1,2.5,5,7.5,10,15,20,30,40,50,70,100,150,200,250,300,400,500,600,750])
cs = m.contourf(x,y,data,clevs,cmap=cm.s3pcpn)
# new axis for colorbar
pos = ax.get_position()
l, b, w, h = pos.bounds
cax = plt.axes([l+w+0.025, b, 0.025, h]) # setup colorbar axes
# draw colorbar.
plt.colorbar(cs, cax, format='%g', ticks=clevs, drawedges=False)
plt.axes(ax) # make the original axes current again