Index: trunk/MESOSCALE_DEV/PLOT/PYTHON/scripts/domain.py
===================================================================
--- trunk/MESOSCALE_DEV/PLOT/PYTHON/scripts/domain.py	(revision 254)
+++ 	(revision )
@@ -1,67 +1,0 @@
-#!/usr/bin/env python
-
-### A. Spiga -- LMD -- 30/06/2011 -- slight modif early 07/2011
-
-def domain (namefile,proj=None,back="vishires",target=None): 
-    from netCDF4 import Dataset
-    from myplot import getcoord2d,define_proj,makeplotres,simplinterv,getprefix,dumpbdy,getproj,latinterv,wrfinterv,simplinterv
-    from mymath import max,min
-    from matplotlib.pyplot import contourf,rcParams,pcolor
-    from numpy.core.defchararray import find
-    from numpy import arange
-    ###
-    nc  = Dataset(namefile)
-    ###
-    if proj == None:  proj = "ortho" #proj = getproj(nc)
-    ###
-    prefix = namefile[0] + namefile[1] + namefile[2]
-    if prefix == "geo":   
-        [lon2d,lat2d] = getcoord2d(nc,nlat='XLAT_M',nlon='XLONG_M')
-        var = 'HGT_M'
-        zeplot = "domain" 
-    else:                 
-        [lon2d,lat2d] = getcoord2d(nc)
-        var = "HGT"
-        zeplot = getprefix(nc) + "domain"
-    ###
-    lon2d = dumpbdy(lon2d,5)
-    lat2d = dumpbdy(lat2d,5)
-    if proj == "npstere":             [wlon,wlat] = latinterv("North_Pole")
-    elif proj in ["lcc","laea"]:      [wlon,wlat] = wrfinterv(lon2d,lat2d)
-    else:                             [wlon,wlat] = simplinterv(lon2d,lat2d)
-    ###
-    m = define_proj(proj,wlon,wlat,back=back)
-    x, y = m(lon2d, lat2d)
-    ###
-    what_I_plot = dumpbdy(nc.variables[var][0,:,:], 5)
-    #levinterv = 250.
-    #zelevels = arange(min(what_I_plot)-levinterv,max(what_I_plot)+levinterv,levinterv)
-    zelevels = 30
-    contourf(x, y, what_I_plot, zelevels)
-    #pcolor(x,y,what_I_plot)  ## on voit trop les lignes !
-    ###
-    if not target:   zeplot = namefile[0:find(namefile,'wrfout')] + zeplot
-    else:            zeplot = target + "/" + zeplot          
-    ###
-    pad_inches_value = 0.35
-    makeplotres(zeplot,res=100.,pad_inches_value=pad_inches_value) #,erase=True)  ## a miniature
-    makeplotres(zeplot,res=200.,pad_inches_value=pad_inches_value,disp=False)
-    #makeplotpng(zeplot,pad_inches_value=0.35)
-    #rcParams['savefig.facecolor'] = 'black'
-    #makeplotpng(zeplot+"b",pad_inches_value=0.35)
-
-if __name__ == "__main__":
-    import sys
-    ### to be replaced by argparse
-    from optparse import OptionParser
-    parser = OptionParser()
-    parser.add_option('-f', action='store', dest='namefile',    type="string",  default=None,       help='name of WRF file [NEEDED]')
-    parser.add_option('-p', action='store', dest='proj',        type="string",  default=None,       help='projection')
-    parser.add_option('-b', action='store', dest='back',        type="string",  default="vishires", help='background')
-    parser.add_option('-t', action='store', dest='target',      type="string",  default=None,       help='destination folder')
-    (opt,args) = parser.parse_args()
-    if opt.namefile is None:
-        print "I want to eat one file at least ! Use domain.py -f name_of_my_file. Or type domain.py -h"
-        exit()
-    print "Options:", opt
-    domain (opt.namefile,proj=opt.proj,back=opt.back,target=opt.target)
Index: trunk/MESOSCALE_DEV/PLOT/PYTHON/scripts/winds.py
===================================================================
--- trunk/MESOSCALE_DEV/PLOT/PYTHON/scripts/winds.py	(revision 254)
+++ 	(revision )
@@ -1,351 +1,0 @@
-#!/usr/bin/env python
-
-### A. Spiga -- LMD -- 30/06/2011 to 10/07/2011
-### Thanks to A. Colaitis for the parser trick
-
-
-####################################
-####################################
-### The main program to plot vectors
-def winds (namefile,\
-           nvert,\
-           proj=None,\
-           back=None,\
-           target=None,
-           stride=3,\
-           numplot=2,\
-           var=None,\
-           colorb="def",\
-           winds=True,\
-           addchar=None,\
-           interv=[0,1],\
-           vmin=None,\
-           vmax=None,\
-           tile=False,\
-           zoom=None,\
-           display=True,\
-           itstep=None,\
-           hole=False,\
-           save="gui",\
-           anomaly=False,\
-           var2=None,\
-           ndiv=10,\
-           first=1):
-
-    ####################################################################################################################
-    ### Colorbars http://www.scipy.org/Cookbook/Matplotlib/Show_colormaps?action=AttachFile&do=get&target=colormaps3.png
-
-    #################################
-    ### Load librairies and functions
-    from netCDF4 import Dataset
-    from myplot import getcoord2d,define_proj,makeplotres,simplinterv,vectorfield,ptitle,latinterv,getproj,wrfinterv,dumpbdy,\
-                       fmtvar,definecolorvec,defcolorb,getprefix,putpoints,calculate_bounds,errormess,definesubplot,\
-                       zoomset,getcoorddef,getwinddef,whatkindfile,reducefield,bounds,getstralt,getfield,smooth,nolow,\
-                       getname,localtime,polarinterv
-    from mymath import deg,max,min,mean
-    from matplotlib.pyplot import contour,contourf, subplot, figure, rcParams, savefig, colorbar, pcolor, show
-    from matplotlib.cm import get_cmap
-    import numpy as np
-    from numpy.core.defchararray import find
-
-    #if save == 'eps': rcParams['backend'] = 'PS'
-
-    ######################
-    ### Load NETCDF object
-    nc  = Dataset(namefile)  
-    
-    ##################################
-    ### Initial checks and definitions
-    typefile = whatkindfile(nc)                                  ## TYPEFILE
-    if var not in nc.variables: var = False                      ## VAR
-    if winds:                                                    ## WINDS
-        [uchar,vchar,metwind] = getwinddef(nc)             
-        if uchar == 'not found': winds = False
-    if not var and not winds: errormess("please set at least winds or var",printvar=nc.variables)
-    [lon2d,lat2d] = getcoorddef(nc)                              ## COORDINATES, could be moved below
-    if proj == None:   proj = getproj(nc)                        ## PROJECTION
-
-    ##########################
-    ### Define plot boundaries
-    ### todo: possible areas in latinterv in argument (ex: "Far_South_Pole")
-    if proj in ["npstere","spstere"]: [wlon,wlat] = polarinterv(lon2d,lat2d)
-    elif proj in ["lcc","laea"]:      [wlon,wlat] = wrfinterv(lon2d,lat2d)
-    else:                             [wlon,wlat] = simplinterv(lon2d,lat2d)
-    if zoom:                          [wlon,wlat] = zoomset(wlon,wlat,zoom) 
-
-    #########################################
-    ### Name for title and graphics save file
-    basename = getname(var=var,winds=winds,anomaly=anomaly)
-    basename = basename + getstralt(nc,nvert)  ## can be moved elsewhere for a more generic routine
-
-    ##################################
-    ### Open a figure and set subplots
-    fig = figure()
-    sub = definesubplot( numplot, fig ) 
- 
-    #################################
-    ### Time loop for plotting device
-    found_lct = False
-    itime = first
-    nplot = 1
-    error = False
-    if itstep is None and numplot > 0: itstep = int(24./numplot)
-    elif numplot <= 0:                 itstep = 1
-    while error is False: 
-
-       ### Which local time ?
-       ltst = localtime ( interv[0]+itime*interv[1], 0.5*(wlon[0]+wlon[1]) )
-
-       ### General plot settings
-       #print itime, int(ltst), numplot, nplot
-       if numplot >= 1: 
-           if nplot > numplot: break
-           if numplot > 1:     
-               if typefile not in ['geo']:  subplot(sub+nplot-1)
-           found_lct = True
-       ### If only one local time is requested (numplot < 0)
-       elif numplot <= 0: 
-           if int(ltst) + numplot != 0:	
-                 itime += 1 
-                 if found_lct is True: break     ## because it means LT was found at previous iteration
-                 else:                 continue  ## continue to iterate to find the correct LT
-           else: 
-                 found_lct = True
-
-       ### Map projection
-       m = define_proj(proj,wlon,wlat,back=back)
-       x, y = m(lon2d, lat2d)
-
-       #### Contour plot
-       if var2:
-           what_I_contour, error = reducefield( getfield(nc,var2), d4=itime, d3=nvert )
-           if not error:
-               if typefile in ['mesoapi','meso']:    what_I_contour = dumpbdy(what_I_contour,6)
-               zevmin, zevmax = calculate_bounds(what_I_contour)
-               zelevels = np.linspace(zevmin,zevmax,num=20)
-               if var2 == 'HGT':  zelevels = np.arange(-10000.,30000.,2000.)
-               contour( x, y, what_I_contour, zelevels, colors='k', linewidths = 0.33 ) #colors='w' )# , alpha=0.5)
-
-       #### Shaded plot
-       if var:
-           what_I_plot, error = reducefield( getfield(nc,var), d4=itime, d3=nvert )
-           if not error: 
-               fvar = var
-               ###
-               if anomaly:
-                   what_I_plot = 100. * ((what_I_plot / smooth(what_I_plot,12)) - 1.)
-                   fvar = 'anomaly'
-               ###
-               if typefile in ['mesoapi','meso']:    what_I_plot = dumpbdy(what_I_plot,6)
-               zevmin, zevmax = calculate_bounds(what_I_plot,vmin=vmin,vmax=vmax)
-               if colorb in ["def","nobar"]:   palette = get_cmap(name=defcolorb(fvar))
-               else:                           palette = get_cmap(name=colorb)
-               if not tile:
-                   if not hole: what_I_plot = bounds(what_I_plot,zevmin,zevmax)
-                   zelevels = np.linspace(zevmin*(1. + 1.e-7),zevmax*(1. - 1.e-7)) #,num=20)
-                   contourf( x, y, what_I_plot, zelevels, cmap = palette )
-               else:
-                   if hole:  what_I_plot = nolow(what_I_plot) 
-                   pcolor( x, y, what_I_plot, cmap = palette, \
-                           vmin=zevmin, vmax=zevmax )
-               if colorb != 'nobar' and var != 'HGT':              
-                         colorbar(fraction=0.05,pad=0.03,format=fmtvar(fvar),\
-                                           ticks=np.linspace(zevmin,zevmax,ndiv+1),\
-                                           extend='neither',spacing='proportional')
-                                           # both min max neither 
-           
-       ### Vector plot
-       if winds:
-           vecx, error = reducefield( getfield(nc,uchar), d4=itime, d3=nvert )
-           vecy, error = reducefield( getfield(nc,vchar), d4=itime, d3=nvert )
-           if not error:
-               if typefile in ['mesoapi','meso']:    
-                   [vecx,vecy] = [dumpbdy(vecx,6,stag=uchar), dumpbdy(vecy,6,stag=vchar)]
-                   key = True
-               elif typefile in ['gcm']:            
-                   key = False
-               if metwind:  [vecx,vecy] = m.rotate_vector(vecx, vecy, lon2d, lat2d)
-               if var == False:       colorvec = definecolorvec(back)
-               else:                  colorvec = definecolorvec(colorb)
-               vectorfield(vecx, vecy,\
-                          x, y, stride=stride, csmooth=2,\
-                          #scale=15., factor=300., color=colorvec, key=key)
-                          scale=20., factor=250., color=colorvec, key=key)
-                                            #200.         ## or csmooth=stride
-               
-       ### Next subplot
-       plottitle = basename
-       if typefile in ['mesoapi','meso']:
-            if addchar:  plottitle = plottitle + addchar + "_LT"+str(ltst)
-            else:        plottitle = plottitle + "_LT"+str(ltst)
-       ptitle( plottitle )
-       itime += itstep
-       nplot += 1
-
-    ##########################################################################
-    ### Save the figure in a file in the data folder or an user-defined folder
-    if typefile in ['meso','mesoapi']:   prefix = getprefix(nc)   
-    elif typefile in ['gcm','gcmex']:            prefix = 'LMD_GCM_'
-    else:                                prefix = ''
-    ###
-    zeplot = prefix + basename 
-    if addchar:         zeplot = zeplot + addchar
-    if numplot <= 0:    zeplot = zeplot + "_LT"+str(abs(numplot))
-    ###
-    if not target:      zeplot = namefile[0:find(namefile,'wrfout')] + zeplot
-    else:               zeplot = target + "/" + zeplot  
-    ###
-    if found_lct:     
-        pad_inches_value = 0.35
-        if save == 'png': 
-            makeplotres(zeplot,res=100.,pad_inches_value=pad_inches_value) #,erase=True)  ## a miniature
-            makeplotres(zeplot,res=200.,pad_inches_value=pad_inches_value,disp=False)
-        elif save in ['eps','svg','pdf']:
-            makeplotres(zeplot,         pad_inches_value=pad_inches_value,disp=False,ext=save)
-        elif save == 'gui':
-            show()
-        else: 
-            print "save mode not supported. using gui instead."
-            show()
-    else:             print "Local time not found"
-
-    ###############
-    ### Now the end
-    return zeplot
-
-##############################
-### A specific stuff for below
-def adjust_length (tab, zelen):
-    from numpy import ones
-    if tab is None:
-        outtab = ones(zelen) * -999999
-    else:
-        if zelen != len(tab):
-            print "not enough or too much values... setting same values all variables"
-            outtab = ones(zelen) * tab[0]
-        else:
-            outtab = tab
-    return outtab
-
-###########################################################################################
-###########################################################################################
-### What is below relate to running the file as a command line executable (very convenient)
-if __name__ == "__main__":
-    import sys
-    from optparse import OptionParser    ### to be replaced by argparse
-    from api_wrapper import api_onelevel
-    from netCDF4 import Dataset
-    from myplot import getlschar
-    from os import system
-
-    #############################
-    ### Get options and variables
-    parser = OptionParser()
-    parser.add_option('-f', '--file',   action='append',dest='namefile', type="string",  default=None,  help='[NEEDED] name of WRF file (append)')
-    parser.add_option('-l', '--level',  action='store',dest='nvert',     type="float",   default=0,     help='level (def=0)(-i 2: p,mbar)(-i 3,4: z,km)')
-    parser.add_option('-p', '--proj',   action='store',dest='proj',      type="string",  default=None,  help='projection')
-    parser.add_option('-b', '--back',   action='store',dest='back',      type="string",  default=None,  help='background image (def: None)')
-    parser.add_option('-t', '--target', action='store',dest='target',    type="string",  default=None,  help='destination folder')
-    parser.add_option('-s', '--stride', action='store',dest='stride',    type="int",     default=3,     help='stride vectors (def=3)')
-    parser.add_option('-v', '--var',    action='append',dest='var',      type="string",  default=None,  help='variable color-shaded (append)')
-    parser.add_option('-n', '--num',    action='store',dest='numplot',   type="int",     default=2,     help='plot number (def=2)(<0: plot LT -*numplot*)')
-    parser.add_option('-i', '--interp', action='store',dest='interp',    type="int",     default=None,  help='interpolation (2: p, 3: z-amr, 4:z-als)')
-    parser.add_option('-c', '--color',  action='store',dest='colorb',    type="string",  default="def", help='change colormap (nobar: no colorbar)')
-    parser.add_option('-x', '--no-vect',action='store_false',dest='winds',               default=True,  help='no wind vectors')
-    parser.add_option('-m', '--min',    action='append',dest='vmin',     type="float",   default=None,  help='bounding minimum value (append)')    
-    parser.add_option('-M', '--max',    action='append',dest='vmax',     type="float",   default=None,  help='bounding maximum value (append)') 
-    parser.add_option('-T', '--tiled',  action='store_true',dest='tile',                 default=False, help='draw a tiled plot (no blank zone)')
-    parser.add_option('-z', '--zoom',   action='store',dest='zoom',      type="float",   default=None,  help='zoom factor in %')
-    parser.add_option('-N', '--no-api', action='store_true',dest='nocall',               default=False, help='do not recreate api file')
-    parser.add_option('-d', '--display',action='store_false',dest='display',             default=True,  help='do not pop up created images')
-    parser.add_option('-e', '--itstep', action='store',dest='itstep',    type="int",     default=None,  help='stride time (def=4)')
-    parser.add_option('-H', '--hole',   action='store_true',dest='hole',                 default=False, help='holes above max and below min')
-    parser.add_option('-S', '--save',   action='store',dest='save',      type="string",  default="gui", help='save mode (png,eps,svg,pdf or gui)(def=gui)')
-    parser.add_option('-a', '--anomaly',action='store_true',dest='anomaly',              default=False, help='compute and plot relative anomaly in %')
-    parser.add_option('-w', '--with',   action='store',dest='var2',     type="string",   default=None,  help='variable contoured')
-    parser.add_option('--div',          action='store',dest='ndiv',     type="int",      default=10,    help='number of divisions in colorbar (def: 10)')
-    parser.add_option('-F', '--first',  action='store',dest='first',    type="int",      default=1,     help='first subscript to plot (def: 1)')
-    #parser.add_option('-V', action='store', dest='comb',        type="float",   default=None,  help='a defined combination of variables')
-    (opt,args) = parser.parse_args()
-    if opt.namefile is None: 
-        print "I want to eat one file at least ! Use winds.py -f name_of_my_file. Or type winds.py -h"
-        exit()
-    if opt.var is None and opt.anomaly is True:
-        print "Cannot ask to compute anomaly if no variable is set"
-        exit()    
-    print "Options:", opt
-
-    listvar = '' 
-    if opt.var is None:
-        zerange = [-999999]
-    else:
-        zelen = len(opt.var)
-        zerange = range(zelen)
-        #if zelen == 1: listvar = opt.var[0] + ','
-        #else         : 
-        for jjj in zerange: listvar += opt.var[jjj] + ','
-        listvar = listvar[0:len(listvar)-1]
-        vmintab = adjust_length (opt.vmin, zelen)  
-        vmaxtab = adjust_length (opt.vmax, zelen)
-
-    for i in range(len(opt.namefile)):
-
-        zefile = opt.namefile[i]
-        print zefile    
-        zelevel = opt.nvert   
-        stralt = None
-        [lschar,zehour,zehourin] = getlschar ( zefile )  ## getlschar from wrfout (or simply return "" if another file)
-    
-        #####################################################
-        ### Call Fortran routines for vertical interpolations
-        if opt.interp is not None:
-            if opt.nvert is 0 and opt.interp is 4:  zelevel = 0.010
-            ### winds or no winds
-            if opt.winds            :  zefields = 'uvmet'
-            else                    :  zefields = ''
-            ### var or no var
-            #if opt.var is None      :  pass
-            if zefields == ''       :  zefields = listvar 
-            else                    :  zefields = zefields + "," + listvar 
-            if opt.var2 is not None : zefields = zefields + "," + opt.var2  
-            print zefields
-            zefile = api_onelevel (  path_to_input   = '', \
-                                     input_name      = zefile, \
-                                     fields          = zefields, \
-                                     interp_method   = opt.interp, \
-                                     onelevel        = zelevel, \
-                                     nocall          = opt.nocall )
-            print zefile
-            zelevel = 0 ## so that zelevel could play again the role of nvert
-
-        if opt.var is None: zerange = [-999999]
-        else:               zerange = range(zelen) 
-        for jjj in zerange:
-            if jjj == -999999: 
-                argvar  = None
-                argvmin = None
-                argvmax = None
-            else:
-                argvar = opt.var[jjj]
-                if vmintab[jjj] != -999999:  argvmin = vmintab[jjj]
-                else:                        argvmin = None
-                if vmaxtab[jjj] != -999999:  argvmax = vmaxtab[jjj] 
-                else:                        argvmax = None
-            #############
-            ### Main call
-            name = winds (zefile,int(zelevel),\
-                proj=opt.proj,back=opt.back,target=opt.target,stride=opt.stride,var=argvar,\
-                numplot=opt.numplot,colorb=opt.colorb,winds=opt.winds,\
-                addchar=lschar,interv=[zehour,zehourin],vmin=argvmin,vmax=argvmax,\
-                tile=opt.tile,zoom=opt.zoom,display=opt.display,\
-                itstep=opt.itstep,hole=opt.hole,save=opt.save,\
-                anomaly=opt.anomaly,var2=opt.var2,ndiv=opt.ndiv,first=opt.first)
-            print 'Done: '+name
-            system("rm -f to_be_erased")
-  
-        #########################################################
-        ### Generate a .sh file with the used command saved in it
-        command = ""
-        for arg in sys.argv: command = command + arg + ' '
-        f = open(name+'.sh', 'w')
-        f.write(command)
