[Matplotlib-users] making multi panel figures

[Michael R. <raw...@ya...>]

I have about 140 lines of code that makes a map. I'd like to turn it into a program which makes a multiple panel (map) figure. I understand that subplot will help to do this. Ideally I would like the 140 lines to be like a subroutine called in a loop. In the current code there are two variable which would be passed to the subroutine, thetitle and ncfile. These are only two things different for each panel. So something like:

do irow = 1, 3
   do icolumn = 1, 3

      call mapping code (thetitle,ncfile)


    enddo
enddo

Here is the code. If the above method is not possible, I assume I'll need to repeat the 140 lines N times, where N is the number of panels. 


TIA
Mike
############################################################################


verbose=0 #verbose=2 says a bit more

import sys,getopt

from mpl_toolkits.basemap import Basemap, shiftgrid, cm 
#from netCDF3 import Dataset as NetCDFFile 
from mpl_toolkits.basemap import  NetCDFFile
from pylab import *
#from matplotlib.mlab import csv2rec

alloptions, otherargs= getopt.getopt(sys.argv[1:],'ro:p:X:Y:v:t:l:u:n:') # note the : after o and p
proj='lam'
#plotfile=None
#plotfile='testmap2.png'
usejetrev=False
colorbounds=[None,None]
extratext=""
xvar=None
yvar=None
thevar=None


#  Here set map title and the file containing gridded data to plot
thetitle='Map Title'
ncfile = NetCDFFile('simple_xy.nc', 'r')   # Here's filename


therec=None
thelev=None
cbot=None
ctop=None
startlon=-180 #default assumption for starting longitude
for theopt,thearg in alloptions:
    print theopt,thearg
    if theopt=='-o': # -o needs filename after it, which is now thearg
        plotfile=thearg    
    elif theopt=='-p': 
        proj=thearg
    elif theopt=='-X': 
        xvar=thearg
    elif theopt=='-Y': 
        yvar=thearg
    elif theopt=='-v': 
        thevar=thearg
    elif theopt=='-t': 
        thetitle=thearg
    elif theopt=='-l': 
        cbot=thearg
    elif theopt=='-u': 
        ctop=thearg
    elif theopt=='-n': 
        therec=thearg
    elif theopt=='-m': 
        thelev=thearg
    elif theopt=='-r': 
        usejetrev=True
    else: #something went wrong
        print "hmm, what are these??? ", theopt, thearg
        sys.exit()

print "\nPlotting, please wait...maybe more than 10 seconds"
if proj=='lam': #Lambert Conformal
    m = Basemap(llcrnrlon=-80.6,llcrnrlat=38.4,urcrnrlon=-66.0,urcrnrlat=47.7,\
            resolution='l',area_thresh=1000.,projection='lcc',\
            lat_1=65.,lon_0=-73.3)
    xtxt=200000. #offset for text
    ytxt=200000.
    parallels = arange(38.,48.,2.)
    meridians = arange(-80.,-64.,2.)
else: #cylindrical is default
#    m = Basemap(llcrnrlon=-180.,llcrnrlat=-90,urcrnrlon=180.,urcrnrlat=90.,\
#            resolution='c',area_thresh=10000.,projection='cyl')
    m = Basemap(llcrnrlon=startlon,llcrnrlat=-90,urcrnrlon=startlon+360.,urcrnrlat=90.,\
            resolution='c',area_thresh=10000.,projection='cyl')
    xtxt=1.
    ytxt=0.
    parallels = arange(-90.,90.,30.)
    if startlon==-180:
        meridians = arange(-180.,180.,60.)
    else:
        meridians = arange(0.,360.,60.)

if verbose>1: print m.__doc__ 
xsize = rcParams['figure.figsize'][0]
fig=figure(figsize=(xsize,m.aspect*xsize))
#ax = fig.add_axes([0.08,0.1,0.7,0.7],axisbg='white')
ax = fig.add_axes([0.07,0.00,0.86,1.0],axisbg='white')
# make a pcolor plot.
#x, y = m(lons, lats)
#p = m.pcolor(x,y,maskdat,shading='flat',cmap=cmap)
#clim(*colorbounds)

# axes units units are left, bottom, width, height
#cax = axes([0.85, 0.1, 0.05, 0.7])  #  colorbar axes for map w/ no graticule
#cax = axes([0.88, 0.1, 0.06, 0.81])  #  colorbar axes for map w/ graticule

axes(ax)  # make the original axes current again

#########    Plot symbol at station locations    #################

# draw coastlines and political boundaries.
m.drawcoastlines()
m.drawcountries()
m.drawstates()
# draw parallels and meridians.
# label on left, right and bottom of map.
#m.drawparallels(parallels,labels=[1,0,0,0])
#m.drawmeridians(meridians,labels=[1,1,0,1])

if not thetitle:
    title(thevar+extratext)
else:
    title(thetitle)

#data = csv2rec('latlon_GHCN.txt',delimiter=' ',names=['lat','lon'])
#for i in range(len(data)):
#    x,y=m(data['lon'][i],data['lat'][i]) # Translate to basemap (Lambert) coordinate space
##    ax.text(x,y,'.')
#    plot(x,y,color='black',marker='.',markersize=6.0)

xpt,ypt = m(-75.0,43.0)
ax.text(xpt,ypt,'*')

    
#if plotfile:
#    savefig(plotfile, dpi=72, facecolor='w', bbox_inches='tight', edgecolor='w', orientation='portrait')
#else:
#    show()

#plt.savefig('map.png')
plt.savefig('map.eps')
show()
#  comment show to mass produce