source: trunk/MESOSCALE_DEV/PLOT/PYTHON/mylib/myplot.py @ 240

Last change on this file since 240 was 240, checked in by aslmd, 14 years ago

MESOSCALE: corrected a bug in the previous commit in makemeso, works now fine. added anomaly calculations in python interface + better treatment for output; default is now GUI where one can zoom etc... and further explore the plotted field
File size: 23.3 KB
RevLine 
[233]1def errormess(text):
2    print text
3    exit()
4    return
5
6def whatkindfile (nc):
7    if 'controle' in nc.variables:   typefile = 'gcm'
8    elif 'vert' in nc.variables:     typefile = 'mesoapi'
9    elif 'U' in nc.variables:        typefile = 'meso'
10    elif 'HGT_M' in nc.variables:    typefile = 'geo'
11    else:                            errormess("whatkindfile: typefile not supported.")
12    return typefile
13
14def getfield (nc,var):
15    ## this allows to get much faster (than simply referring to nc.variables[var])
16    dimension = len(nc.variables[var].dimensions)
17    if dimension == 2:    field = nc.variables[var][:,:]
18    elif dimension == 3:  field = nc.variables[var][:,:,:]
19    elif dimension == 4:  field = nc.variables[var][:,:,:,:]
20    return field
21
22def reducefield (input,d4=None,d3=None,d2=None,d1=None):
23    ### do it the reverse way to be compliant with netcdf "t z y x" or "t y x" or "y x"
24    ### it would be actually better to name d4 d3 d2 d1 as t z y x
25    import numpy as np
26    dimension = np.array(input).ndim
27    shape = np.array(input).shape
28    print 'dim,shape: ',dimension,shape
29    output = input
30    error = False
31    if dimension == 2:
32        if   d2 >= shape[0]: error = True
33        elif d1 >= shape[1]: error = True
34        elif d1 is not None and d2 is not None:  output = input[d2,d1]
35        elif d1 is not None:         output = input[:,d1]
36        elif d2 is not None:         output = input[d2,:]
37    elif dimension == 3:
38        if   d4 >= shape[0]: error = True
39        elif d2 >= shape[1]: error = True
40        elif d1 >= shape[2]: error = True
41        elif d4 is not None and d2 is not None and d1 is not None:  output = input[d4,d2,d1]
42        elif d4 is not None and d2 is not None:         output = input[d4,d2,:]
43        elif d4 is not None and d1 is not None:         output = input[d4,:,d1]
44        elif d2 is not None and d1 is not None:         output = input[:,d2,d1]
45        elif d1 is not None:                output = input[:,:,d1]
46        elif d2 is not None:                output = input[:,d2,:]
47        elif d4 is not None:                output = input[d4,:,:]
48    elif dimension == 4:
49        if   d4 >= shape[0]: error = True
50        elif d3 >= shape[1]: error = True
51        elif d2 >= shape[2]: error = True
52        elif d1 >= shape[3]: error = True
53        elif d4 is not None and d3 is not None and d2 is not None and d1 is not None:  output = input[d4,d3,d2,d1]
54        elif d4 is not None and d3 is not None and d2 is not None:         output = input[d4,d3,d2,:]
55        elif d4 is not None and d3 is not None and d1 is not None:         output = input[d4,d3,:,d1]
56        elif d4 is not None and d2 is not None and d1 is not None:         output = input[d4,:,d2,d1]
57        elif d3 is not None and d2 is not None and d1 is not None:         output = input[:,d3,d2,d1]
58        elif d4 is not None and d3 is not None:                output = input[d4,d3,:,:]
59        elif d4 is not None and d2 is not None:                output = input[d4,:,d2,:]
60        elif d4 is not None and d1 is not None:                output = input[d4,:,:,d1]
61        elif d3 is not None and d2 is not None:                output = input[:,d3,d2,:]
62        elif d3 is not None and d1 is not None:                output = input[:,d3,:,d1]
63        elif d2 is not None and d1 is not None:                output = input[:,:,d2,d1]
64        elif d1 is not None:                       output = input[:,:,:,d1]
65        elif d2 is not None:                       output = input[:,:,d2,:]
66        elif d3 is not None:                       output = input[:,d3,:,:]
67        elif d4 is not None:                       output = input[d4,:,:,:]
68    dimension = np.array(output).ndim
69    shape = np.array(output).shape
70    print 'dim,shape: ',dimension,shape
71    return output, error
72
[180]73def latinterv (area):
74        if   area == "Europe": 
75                wlat = [20.,80.]
76                wlon = [-50.,50.]
77        elif area == "Central_America":
78                wlat = [-10.,40.]
79                wlon = [230.,300.]
80        elif area == "Africa":
81                wlat = [-20.,50.]
82                wlon = [-50.,50.]
83        elif area == "Whole":
84                wlat = [-90.,90.]
85                wlon = [-180.,180.]
86        elif area == "Southern_Hemisphere":
87                wlat = [-90.,60.]
88                wlon = [-180.,180.]
89        elif area == "Northern_Hemisphere":
90                wlat = [-60.,90.]
91                wlon = [-180.,180.]
92        elif area == "Tharsis":
93                wlat = [-30.,60.]
94                wlon = [-170.,-10.]
95        elif area == "Whole_No_High":
96                wlat = [-60.,60.]
97                wlon = [-180.,180.]
98        elif area == "Chryse":
99                wlat = [-60.,60.]
100                wlon = [-60.,60.]
101        elif area == "North_Pole":
[240]102                wlat = [50.,90.]
[180]103                wlon = [-180.,180.]
104        elif area == "Close_North_Pole":
105                wlat = [75.,90.]
106                wlon = [-180.,180.]
[240]107        elif area == "South_Pole":
108                wlat = [-90.,-50.]
109                wlon = [-180.,180.]
110        elif area == "Close_South_Pole":
111                wlat = [-90.,-75.]
112                wlon = [-180.,180.]
[180]113        return wlon,wlat
114
[233]115def definesubplot ( numplot, fig ):
116    from matplotlib.pyplot import rcParams
117    rcParams['font.size'] = 12. ## default (important for multiple calls)
118    if   numplot == 4:
119        sub = 221
120        fig.subplots_adjust(wspace = 0.3, hspace = 0.3)
121        rcParams['font.size'] = int( rcParams['font.size'] * 2. / 3. )
122    elif numplot == 2:
123        sub = 121
124        fig.subplots_adjust(wspace = 0.35)
125        rcParams['font.size'] = int( rcParams['font.size'] * 3. / 4. )
126    elif numplot == 3:
127        sub = 131
128        fig.subplots_adjust(wspace = 0.5)
129        rcParams['font.size'] = int( rcParams['font.size'] * 1. / 2. )
130    elif numplot == 6:
131        sub = 231
132        fig.subplots_adjust(wspace = 0.4, hspace = 0.0)
133        rcParams['font.size'] = int( rcParams['font.size'] * 1. / 2. )
134    elif numplot == 8:
135        sub = 331 #241
136        fig.subplots_adjust(wspace = 0.3, hspace = 0.3)
137        rcParams['font.size'] = int( rcParams['font.size'] * 1. / 2. )
138    elif numplot == 9:
139        sub = 331
140        fig.subplots_adjust(wspace = 0.3, hspace = 0.3)
141        rcParams['font.size'] = int( rcParams['font.size'] * 1. / 2. )
142    elif numplot == 1:
143        sub = 99999
144    elif numplot < 0: 
145        sub = 99999
146    else:
147        print "supported: 1,2,3,4,6,8,9"
148        exit()
149    return sub
150
151def getstralt(nc,nvert):
152    typefile = whatkindfile(nc)
153    if typefile is 'meso':                     
154        stralt = "_lvl" + str(nvert)
155    elif typefile is 'mesoapi':
156        zelevel = int(nc.variables['vert'][nvert])
157        if abs(zelevel) < 10000.:   strheight=str(zelevel)+"m"
158        else:                       strheight=str(int(zelevel/1000.))+"km"
159        if 'altitude'       in nc.dimensions:   stralt = "_"+strheight+"-AMR"
160        elif 'altitude_abg' in nc.dimensions:   stralt = "_"+strheight+"-ALS"
161        elif 'bottom_top'   in nc.dimensions:   stralt = "_"+strheight
162        elif 'pressure'     in nc.dimensions:   stralt = "_"+str(zelevel)+"Pa"
163        else:                                   stralt = ""
164    else:
165        stralt = ""
166    return stralt
167
[195]168def getlschar ( namefile ):
169    from netCDF4 import Dataset
170    from timestuff import sol2ls
[233]171    from numpy import array
[195]172    nc  = Dataset(namefile)
[237]173    zetime = None
[233]174    if 'Times' in nc.variables: 
175        zetime = nc.variables['Times'][0]
176        shape = array(nc.variables['Times']).shape
177        if shape[0] < 2: zetime = None
178    if zetime is not None \
[225]179       and 'vert' not in nc.variables:
[233]180       #### strangely enough this does not work for api or ncrcat results!
[195]181        zetimestart = getattr(nc, 'START_DATE')
182        zeday = int(zetime[8]+zetime[9]) - int(zetimestart[8]+zetimestart[9])
183        if zeday < 0:    lschar=""  ## might have crossed a month... fix soon
184        else:            lschar="_Ls"+str( int( sol2ls ( getattr( nc, 'JULDAY' ) + zeday ) ) )
[197]185        ###
186        zetime2 = nc.variables['Times'][1]
187        one  = int(zetime[11]+zetime[12]) + int(zetime[14]+zetime[15])/37.
188        next = int(zetime2[11]+zetime2[12]) + int(zetime2[14]+zetime2[15])/37. 
189        zehour    = one
190        zehourin  = abs ( next - one )
[195]191    else:
192        lschar=""
[197]193        zehour = 0
194        zehourin = 1 
195    return lschar, zehour, zehourin
[195]196
[202]197def getprefix (nc):
198    prefix = 'LMD_MMM_'
199    prefix = prefix + 'd'+str(getattr(nc,'GRID_ID'))+'_'
200    prefix = prefix + str(int(getattr(nc,'DX')/1000.))+'km_'
201    return prefix
202
[184]203def getproj (nc):
[233]204    typefile = whatkindfile(nc)
205    if typefile in ['mesoapi','meso','geo']:
206        ### (il faudrait passer CEN_LON dans la projection ?)
207        map_proj = getattr(nc, 'MAP_PROJ')
208        cen_lat  = getattr(nc, 'CEN_LAT')
209        if map_proj == 2:
210            if cen_lat > 10.:   
211                proj="npstere"
212                print "NP stereographic polar domain" 
213            else:           
214                proj="spstere"
215                print "SP stereographic polar domain"
216        elif map_proj == 1: 
217            print "lambert projection domain" 
218            proj="lcc"
219        elif map_proj == 3: 
220            print "mercator projection"
221            proj="merc"
222        else:
223            proj="merc"
224    elif typefile in ['gcm']:  proj="cyl"  ## pb avec les autres (de trace derriere la sphere ?)
225    else:                      proj="ortho"
[184]226    return proj   
227
[180]228def ptitle (name):
229    from matplotlib.pyplot import title
230    title(name)
231    print name
232
233def simplinterv (lon2d,lat2d):
234    import numpy as np
235    return [[np.min(lon2d),np.max(lon2d)],[np.min(lat2d),np.max(lat2d)]]
236
[184]237def wrfinterv (lon2d,lat2d):
238    nx = len(lon2d[0,:])-1
239    ny = len(lon2d[:,0])-1
[225]240    lon1 = lon2d[0,0] 
241    lon2 = lon2d[nx,ny] 
242    lat1 = lat2d[0,0] 
243    lat2 = lat2d[nx,ny] 
[233]244    if abs(0.5*(lat1+lat2)) > 60.:  wider = 0.5 * (abs(lon1)+abs(lon2)) * 0.1
245    else:                           wider = 0.
246    if lon1 < lon2:  wlon = [lon1, lon2 + wider] 
[225]247    else:            wlon = [lon2, lon1 + wider]
248    if lat1 < lat2:  wlat = [lat1, lat2]
249    else:            wlat = [lat2, lat1]
250    return [wlon,wlat]
[184]251
[240]252def makeplotres (filename,res=None,pad_inches_value=0.25,folder='',disp=True,ext='png',erase=False):
[180]253    import  matplotlib.pyplot as plt
[240]254    from os import system
255    addstr = ""
256    if res is not None:
257        res = int(res)
258        addstr = "_"+str(res)
259    name = filename+addstr+"."+ext
[186]260    if folder != '':      name = folder+'/'+name
[180]261    plt.savefig(name,dpi=res,bbox_inches='tight',pad_inches=pad_inches_value)
[240]262    if disp:              display(name)
263    if ext in ['eps','ps','svg']:   system("tar czvf "+name+".tar.gz "+name+" ; rm -f "+name)
264    if erase:   system("mv "+name+" to_be_erased")             
[180]265    return
266
[240]267def dumpbdy (field,n,stag=None):
[184]268    nx = len(field[0,:])-1
269    ny = len(field[:,0])-1
[233]270    if stag == 'U': nx = nx-1
271    if stag == 'V': ny = ny-1
[240]272    return field[n:ny-n,n:nx-n]
[180]273
[233]274def getcoorddef ( nc ):   
275    ## getcoord2d for predefined types
276    typefile = whatkindfile(nc)
277    if typefile in ['mesoapi','meso']:
278        [lon2d,lat2d] = getcoord2d(nc)
[240]279        lon2d = dumpbdy(lon2d,6)
280        lat2d = dumpbdy(lat2d,6)
[233]281    elif typefile in ['gcm']:
282        [lon2d,lat2d] = getcoord2d(nc,nlat="latitude",nlon="longitude",is1d=True)
283    elif typefile in ['geo']:
284        [lon2d,lat2d] = getcoord2d(nc,nlat='XLAT_M',nlon='XLONG_M')
285    return lon2d,lat2d   
286
[184]287def getcoord2d (nc,nlat='XLAT',nlon='XLONG',is1d=False):
288    import numpy as np
289    if is1d:
290        lat = nc.variables[nlat][:]
291        lon = nc.variables[nlon][:]
292        [lon2d,lat2d] = np.meshgrid(lon,lat)
293    else:
294        lat = nc.variables[nlat][0,:,:]
295        lon = nc.variables[nlon][0,:,:]
296        [lon2d,lat2d] = [lon,lat]
297    return lon2d,lat2d
298
[180]299def smooth (field, coeff):
300        ## actually blur_image could work with different coeff on x and y
301        if coeff > 1:   result = blur_image(field,int(coeff))
302        else:           result = field
303        return result
304
305def gauss_kern(size, sizey=None):
306        import numpy as np
307        ## FROM COOKBOOK http://www.scipy.org/Cookbook/SignalSmooth     
308        # Returns a normalized 2D gauss kernel array for convolutions
309        size = int(size)
310        if not sizey:
311                sizey = size
312        else:
313                sizey = int(sizey)
314        x, y = np.mgrid[-size:size+1, -sizey:sizey+1]
315        g = np.exp(-(x**2/float(size)+y**2/float(sizey)))
316        return g / g.sum()
317
318def blur_image(im, n, ny=None) :
319        from scipy.signal import convolve
320        ## FROM COOKBOOK http://www.scipy.org/Cookbook/SignalSmooth
321        # blurs the image by convolving with a gaussian kernel of typical size n.
322        # The optional keyword argument ny allows for a different size in the y direction.
323        g = gauss_kern(n, sizey=ny)
324        improc = convolve(im, g, mode='same')
325        return improc
326
[233]327def getwinddef (nc):   
328    ## getwinds for predefined types
329    typefile = whatkindfile(nc)
330    ###
331    if typefile is 'mesoapi':    [uchar,vchar] = ['Um','Vm']
332    elif typefile is 'gcm':      [uchar,vchar] = ['u','v']
333    elif typefile is 'meso':     [uchar,vchar] = ['U','V']
334    else:                        [uchar,vchar] = ['not found','not found']
335    ###
336    if typefile in ['meso']:     metwind = False ## geometrical (wrt grid)
337    else:                        metwind = True  ## meteorological (zon/mer)
338    if metwind is False:         print "Not using meteorological winds. You trust numerical grid as being (x,y)"
339    ###
340    return uchar,vchar,metwind
[202]341
[184]342def vectorfield (u, v, x, y, stride=3, scale=15., factor=250., color='black', csmooth=1, key=True):
343    ## scale regle la reference du vecteur
344    ## factor regle toutes les longueurs (dont la reference). l'AUGMENTER pour raccourcir les vecteurs.
345    import  matplotlib.pyplot               as plt
346    import  numpy                           as np
[187]347    posx = np.min(x) - np.std(x) / 10.
348    posy = np.min(y) - np.std(y) / 10.
[184]349    u = smooth(u,csmooth)
350    v = smooth(v,csmooth)
[188]351    widthvec = 0.003 #0.005 #0.003
[184]352    q = plt.quiver( x[::stride,::stride],\
353                    y[::stride,::stride],\
354                    u[::stride,::stride],\
355                    v[::stride,::stride],\
[228]356                    angles='xy',color=color,pivot='middle',\
[184]357                    scale=factor,width=widthvec )
[202]358    if color in ['white','yellow']:     kcolor='black'
359    else:                               kcolor=color
[184]360    if key: p = plt.quiverkey(q,posx,posy,scale,\
[194]361                   str(int(scale)),coordinates='data',color=kcolor,labelpos='S',labelsep = 0.03)
[184]362    return 
[180]363
364def display (name):
[184]365    from os import system
366    system("display "+name+" > /dev/null 2> /dev/null &")
367    return name
[180]368
369def findstep (wlon):
[184]370    steplon = int((wlon[1]-wlon[0])/4.)  #3
371    step = 120.
372    while step > steplon and step > 15. :       step = step / 2.
373    if step <= 15.:
374        while step > steplon and step > 5.  :   step = step - 5.
375    if step <= 5.:
376        while step > steplon and step > 1.  :   step = step - 1.
377    if step <= 1.:
378        step = 1. 
[180]379    return step
380
[233]381def define_proj (char,wlon,wlat,back=None):
[180]382    from    mpl_toolkits.basemap            import Basemap
383    import  numpy                           as np
384    import  matplotlib                      as mpl
[240]385    from mymath import max
[180]386    meanlon = 0.5*(wlon[0]+wlon[1])
387    meanlat = 0.5*(wlat[0]+wlat[1])
[184]388    if   wlat[0] >= 80.:   blat =  40. 
[240]389    elif wlat[1] <= -80.:  blat = -40. 
390    else:                  blat = max([wlat[0],wlat[1]])
391    print blat
[207]392    h = 50.  ## en km
[202]393    radius = 3397200.
[184]394    if   char == "cyl":     m = Basemap(rsphere=radius,projection='cyl',\
[180]395                              llcrnrlat=wlat[0],urcrnrlat=wlat[1],llcrnrlon=wlon[0],urcrnrlon=wlon[1])
[184]396    elif char == "moll":    m = Basemap(rsphere=radius,projection='moll',lon_0=meanlon)
397    elif char == "ortho":   m = Basemap(rsphere=radius,projection='ortho',lon_0=meanlon,lat_0=meanlat)
398    elif char == "lcc":     m = Basemap(rsphere=radius,projection='lcc',lat_1=meanlat,lat_0=meanlat,lon_0=meanlon,\
399                              llcrnrlat=wlat[0],urcrnrlat=wlat[1],llcrnrlon=wlon[0],urcrnrlon=wlon[1])
400    elif char == "npstere": m = Basemap(rsphere=radius,projection='npstere', boundinglat=blat, lon_0=0.)
401    elif char == "spstere": m = Basemap(rsphere=radius,projection='spstere', boundinglat=blat, lon_0=0.)
[207]402    elif char == "nplaea":  m = Basemap(rsphere=radius,projection='nplaea', boundinglat=wlat[0], lon_0=meanlon)
403    elif char == "laea":    m = Basemap(rsphere=radius,projection='laea',lon_0=meanlon,lat_0=meanlat,lat_ts=meanlat,\
404                              llcrnrlat=wlat[0],urcrnrlat=wlat[1],llcrnrlon=wlon[0],urcrnrlon=wlon[1])
[184]405    elif char == "nsper":   m = Basemap(rsphere=radius,projection='nsper',lon_0=meanlon,lat_0=meanlat,satellite_height=h*1000.)
406    elif char == "merc":    m = Basemap(rsphere=radius,projection='merc',lat_ts=0.,\
407                              llcrnrlat=wlat[0],urcrnrlat=wlat[1],llcrnrlon=wlon[0],urcrnrlon=wlon[1])
408    fontsizemer = int(mpl.rcParams['font.size']*3./4.)
[207]409    if char in ["cyl","lcc","merc","nsper","laea"]:   step = findstep(wlon)
410    else:                                             step = 10.
[238]411    steplon = step*2.
412    #if back in ["geolocal"]:                         
413    #    step = np.min([5.,step])
414    #    steplon = step
[225]415    print step
[238]416    m.drawmeridians(np.r_[-180.:180.:steplon], labels=[0,0,0,1], color='grey', fontsize=fontsizemer)
[180]417    m.drawparallels(np.r_[-90.:90.:step], labels=[1,0,0,0], color='grey', fontsize=fontsizemer)
[233]418    if back: m.warpimage(marsmap(back),scale=0.75)
419            #if not back:
420            #    if not var:                                        back = "mola"    ## if no var:         draw mola
421            #    elif typefile in ['mesoapi','meso','geo'] \
422            #       and proj not in ['merc','lcc','nsper','laea']:  back = "molabw"  ## if var but meso:   draw molabw
423            #    else:                                              pass             ## else:              draw None
[180]424    return m
425
[232]426#### test temporaire
427def putpoints (map,plot):
428    #### from http://www.scipy.org/Cookbook/Matplotlib/Maps
429    # lat/lon coordinates of five cities.
430    lats = [18.4]
431    lons = [-134.0]
432    points=['Olympus Mons']
433    # compute the native map projection coordinates for cities.
434    x,y = map(lons,lats)
435    # plot filled circles at the locations of the cities.
436    map.plot(x,y,'bo')
437    # plot the names of those five cities.
438    wherept = 0 #1000 #50000
439    for name,xpt,ypt in zip(points,x,y):
440       plot.text(xpt+wherept,ypt+wherept,name)
441    ## le nom ne s'affiche pas...
442    return
443
[233]444def calculate_bounds(field,vmin=None,vmax=None):
445    import numpy as np
446    from mymath import max,min,mean
447    ind = np.where(field < 9e+35)
448    fieldcalc = field[ ind ] # la syntaxe compacte ne marche si field est un tuple
449    ###
450    dev = np.std(fieldcalc)*3.0
451    ###
452    if vmin is None:
453        zevmin = mean(fieldcalc) - dev
454    else:             zevmin = vmin
455    ###
456    if vmax is None:  zevmax = mean(fieldcalc) + dev
457    else:             zevmax = vmax
458    if vmin == vmax:
459                      zevmin = mean(fieldcalc) - dev  ### for continuity
460                      zevmax = mean(fieldcalc) + dev  ### for continuity           
461    ###
462    if zevmin < 0. and min(fieldcalc) > 0.: zevmin = 0.
463    print "field ", min(fieldcalc), max(fieldcalc)
464    print "bounds ", zevmin, zevmax
465    return zevmin, zevmax
[232]466
[233]467def bounds(what_I_plot,zevmin,zevmax):
468    ### might be convenient to add the missing value in arguments
469    what_I_plot[ what_I_plot < zevmin ] = zevmin*(1. + 1.e-7)
470    what_I_plot[ what_I_plot > 9e+35  ] = -9e+35
471    what_I_plot[ what_I_plot > zevmax ] = zevmax*(1. - 1.e-7)
472    return what_I_plot
473
474def zoomset (wlon,wlat,zoom):
475    dlon = abs(wlon[1]-wlon[0])/2.
476    dlat = abs(wlat[1]-wlat[0])/2.
477    [wlon,wlat] = [ [wlon[0]+zoom*dlon/100.,wlon[1]-zoom*dlon/100.],\
478                    [wlat[0]+zoom*dlat/100.,wlat[1]-zoom*dlat/100.] ]
479    print "zoom %",zoom,wlon,wlat
480    return wlon,wlat
481
[201]482def fmtvar (whichvar="def"):
[204]483    fmtvar    =     { \
484             "tk":           "%.0f",\
485             "tpot":         "%.0f",\
486             "def":          "%.1e",\
487             "PTOT":         "%.0f",\
488             "HGT":          "%.1e",\
489             "USTM":         "%.2f",\
[225]490             "HFX":          "%.0f",\
[232]491             "ICETOT":       "%.1e",\
[237]492             "TAU_ICE":      "%.2f",\
[240]493             "anomaly":      "%.1f",\
[204]494                    }
495    if whichvar not in fmtvar:
496        whichvar = "def"
497    return fmtvar[whichvar]
[201]498
[233]499####################################################################################################################
500### Colorbars http://www.scipy.org/Cookbook/Matplotlib/Show_colormaps?action=AttachFile&do=get&target=colormaps3.png
[202]501def defcolorb (whichone="def"):
[204]502    whichcolorb =    { \
503             "def":          "spectral",\
504             "HGT":          "spectral",\
505             "tk":           "gist_heat",\
506             "QH2O":         "PuBu",\
507             "USTM":         "YlOrRd",\
[225]508             "HFX":          "RdYlBu",\
[232]509             "ICETOT":       "YlGnBu",\
[238]510             "TAU_ICE":      "Blues",\
[240]511             "anomaly":      "RdBu_r",\
[204]512                     }
[240]513#spectral BrBG RdBu_r
[204]514    if whichone not in whichcolorb:
515        whichone = "def"
516    return whichcolorb[whichone]
[202]517
518def definecolorvec (whichone="def"):
519        whichcolor =    { \
520                "def":          "black",\
521                "vis":          "yellow",\
522                "vishires":     "yellow",\
523                "molabw":       "yellow",\
524                "mola":         "black",\
525                "gist_heat":    "white",\
526                "hot":          "tk",\
527                "gist_rainbow": "black",\
528                "spectral":     "black",\
529                "gray":         "red",\
530                "PuBu":         "black",\
531                        }
532        if whichone not in whichcolor:
533                whichone = "def"
534        return whichcolor[whichone]
535
[180]536def marsmap (whichone="vishires"):
[233]537        from os import uname
538        mymachine = uname()[1]
539        ### not sure about speed-up with this method... looks the same
540        if "lmd.jussieu.fr" in mymachine: domain = "/u/aslmd/WWW/maps/"
541        else:                             domain = "http://www.lmd.jussieu.fr/~aslmd/maps/"
[180]542        whichlink =     { \
[233]543                #"vis":         "http://maps.jpl.nasa.gov/pix/mar0kuu2.jpg",\
544                #"vishires":    "http://www.lmd.jussieu.fr/~aslmd/maps/MarsMap_2500x1250.jpg",\
545                #"geolocal":    "http://dl.dropbox.com/u/11078310/geolocal.jpg",\
546                #"mola":        "http://www.lns.cornell.edu/~seb/celestia/mars-mola-2k.jpg",\
547                #"molabw":      "http://dl.dropbox.com/u/11078310/MarsElevation_2500x1250.jpg",\
548                "vis":         domain+"mar0kuu2.jpg",\
549                "vishires":    domain+"MarsMap_2500x1250.jpg",\
550                "geolocal":    domain+"geolocal.jpg",\
551                "mola":        domain+"mars-mola-2k.jpg",\
552                "molabw":      domain+"MarsElevation_2500x1250.jpg",\
[238]553                "clouds":      "http://www.johnstonsarchive.net/spaceart/marswcloudmap.jpg",\
554                "jupiter":     "http://www.mmedia.is/~bjj/data/jupiter_css/jupiter_css.jpg",\
555                "jupiter_voy": "http://www.mmedia.is/~bjj/data/jupiter/jupiter_vgr2.jpg",\
[180]556                        }
[238]557        ### see http://www.mmedia.is/~bjj/planetary_maps.html
[180]558        if whichone not in whichlink: 
559                print "marsmap: choice not defined... you'll get the default one... "
560                whichone = "vishires" 
561        return whichlink[whichone]
562
563def earthmap (whichone):
564        if   whichone == "contrast":    whichlink="http://users.info.unicaen.fr/~karczma/TEACH/InfoGeo/Images/Planets/EarthMapAtmos_2500x1250.jpg"
565        elif whichone == "bw":          whichlink="http://users.info.unicaen.fr/~karczma/TEACH/InfoGeo/Images/Planets/EarthElevation_2500x1250.jpg"
566        elif whichone == "nice":        whichlink="http://users.info.unicaen.fr/~karczma/TEACH/InfoGeo/Images/Planets/earthmap1k.jpg"
567        return whichlink
568
Note: See TracBrowser for help on using the repository browser.