Index: trunk/MESOSCALE_DEV/PLOT/PYTHON/scripts/makefig.sh
===================================================================
--- trunk/MESOSCALE_DEV/PLOT/PYTHON/scripts/makefig.sh	(revision 232)
+++ trunk/MESOSCALE_DEV/PLOT/PYTHON/scripts/makefig.sh	(revision 233)
@@ -5,4 +5,5 @@
  file2="$fold/POLAR_THOMAS_ls20/wrfout_d03_2024-01-42_06:00:00"
  file3="$fold/POLAR_THOMAS_ls0/wrfout_d03_2024-01-06_06:00:00"
+ file4="$fold/POLAR_THOMAS_ls10/wrfout_d03_2024-01-19_06:00:00"
  dest="/u/aslmd/WWW/antichambre/thomas/"
 
@@ -12,3 +13,4 @@
            -v USTM -m   0. -M 0.8 \
            -v HGT  -m   0  -M 0 \
-           -f $file1 -f $file2 -f $file3 
+           -f $file1 -f $file2 -f $file3 -f $file4 \
+           -d
Index: trunk/MESOSCALE_DEV/PLOT/PYTHON/scripts/winds.py
===================================================================
--- trunk/MESOSCALE_DEV/PLOT/PYTHON/scripts/winds.py	(revision 232)
+++ trunk/MESOSCALE_DEV/PLOT/PYTHON/scripts/winds.py	(revision 233)
@@ -23,5 +23,7 @@
            vmax=None,\
            tile=False,\
-           zoom=None):
+           zoom=None,\
+           display=True,\
+           itstep=None):
 
     ####################################################################################################################
@@ -32,5 +34,6 @@
     from netCDF4 import Dataset
     from myplot import getcoord2d,define_proj,makeplotpng,simplinterv,vectorfield,ptitle,latinterv,getproj,wrfinterv,dumpbdy,\
-                       fmtvar,definecolorvec,getwinds,defcolorb,getprefix,putpoints
+                       fmtvar,definecolorvec,defcolorb,getprefix,putpoints,calculate_bounds,errormess,definesubplot,\
+                       zoomset,getcoorddef,getwinddef,whatkindfile,reducefield,bounds,getstralt,getfield
     from mymath import deg,max,min,mean
     from matplotlib.pyplot import contour,contourf, subplot, figure, rcParams, savefig, colorbar, pcolor
@@ -39,111 +42,23 @@
     from numpy.core.defchararray import find
 
-    ###
-    #rcParams['text.usetex'] = True
-    #rcParams['cairo.format'] = 'svg'
-
     ######################
     ### Load NETCDF object
     nc  = Dataset(namefile)
-
-    ###################################
-    ### Recognize predefined file types
-    if 'controle' in nc.variables:   typefile = 'gcm'
-    elif 'vert' in nc.variables:     typefile = 'mesoapi'
-    elif 'U' in nc.variables:        typefile = 'meso'
-    elif 'HGT_M' in nc.variables:    typefile = 'geo'
-    else:                            
-        print "typefile not supported."
-        print nc.variables
-        exit()
-
-    ##############################################################
-    ### Try to guess the projection from wrfout if not set by user 
-    if typefile in ['mesoapi','meso','geo']:
-        if proj == None:       proj = getproj(nc)
-                                    ### (il faudrait passer CEN_LON dans la projection ?)
-    elif typefile in ['gcm']:
-        if proj == None:       proj = "cyl"   
-                                    ## pb avec les autres (de trace derriere la sphere ?)
-
-    ############################################
-    #### Choose underlying topography by default 
-    if not back:
-        if not var:                                        back = "mola"    ## if no var:         draw mola
-        elif typefile in ['mesoapi','meso','geo'] \
-           and proj not in ['merc','lcc','nsper','laea']:  back = "molabw"  ## if var but meso:   draw molabw
-        else:                                              pass             ## else:              draw None
-
-    ####################################################
-    ### Get geographical coordinates and plot boundaries
-    if typefile in ['mesoapi','meso']:
-        [lon2d,lat2d] = getcoord2d(nc)
-        lon2d = dumpbdy(lon2d)
-        lat2d = dumpbdy(lat2d)
-    elif typefile in ['gcm']:
-        [lon2d,lat2d] = getcoord2d(nc,nlat="latitude",nlon="longitude",is1d=True)
-    elif typefile in ['geo']:
-        [lon2d,lat2d] = getcoord2d(nc,nlat='XLAT_M',nlon='XLONG_M')
+    
+    ##################################
+    ### Initial checks and definitions
+    typefile = whatkindfile(nc)                                  ## TYPEFILE
+    if var not in nc.variables: var = False                      ## VAR
+    if winds: [uchar,vchar,metwind] = getwinddef(nc)             ## WINDS
+    if uchar == 'not found': winds = False
+    [lon2d,lat2d] = getcoorddef(nc)                              ## COORDINATES, could be moved below
+    if proj == None:   proj = getproj(nc)                        ## PROJECTION
+
+    ##########################
+    ### Define plot boundaries
     if proj == "npstere":             [wlon,wlat] = latinterv("North_Pole")
     elif proj in ["lcc","laea"]:      [wlon,wlat] = wrfinterv(lon2d,lat2d)
     else:                             [wlon,wlat] = simplinterv(lon2d,lat2d)
-    print wlon
-    print wlat
-    if zoom:  
-        dlon = abs(wlon[1]-wlon[0])/2.
-        dlat = abs(wlat[1]-wlat[0])/2.
-        [wlon,wlat] = [ [wlon[0]+zoom*dlon/100.,wlon[1]-zoom*dlon/100.],\
-                        [wlat[0]+zoom*dlat/100.,wlat[1]-zoom*dlat/100.] ]
-        print "zoom %",zoom,wlon,wlat
-
-    ##############################################################################
-    ### Get winds and know if those are meteorological winds (ie. zon, mer) or not
-    if winds:
-        if typefile is 'mesoapi': 
-            [u,v] = getwinds(nc)
-            metwind = True  ## meteorological (zon/mer)
-        elif typefile is 'gcm':
-            [u,v] = getwinds(nc,charu='u',charv='v')
-            metwind = True  ## meteorological (zon/mer)
-        elif typefile is 'meso':
-            [u,v] = getwinds(nc,charu='U',charv='V')
-            metwind = False ## geometrical (wrt grid)
-            print "Beware ! Not using meteorological winds. You trust numerical grid as being (x,y)."
-        elif typefile is 'geo':
-            winds = None
-
-    #####################################################
-    ### Load the chosen variables, whether it is 2D or 3D
-    if var:
-        if var not in nc.variables: 
-            print "not found in file:",var
-            exit()
-        else:   
-            dimension = np.array(nc.variables[var]).ndim
-            if dimension == 2:     field = nc.variables[var][:,:]
-            elif dimension == 3:   field = nc.variables[var][:,:,:]
-            elif dimension == 4:   field = nc.variables[var][:,nvert,:,:]  
-        fieldcalc = field[ field < 9e+35 ]
-        dev = np.std(fieldcalc)*2.0
-        if vmin is None:  zevmin = mean(fieldcalc) - dev
-        else:             zevmin = vmin
-        if vmax is None:  zevmax = mean(fieldcalc) + dev
-        else:             zevmax = vmax
-        if vmin == vmax:   
-                          zevmin = mean(fieldcalc) - dev  ### for continuity
-                          zevmax = mean(fieldcalc) + dev  ### for continuity            
-        print "field ", min(fieldcalc), max(fieldcalc)
-        print "bounds ", zevmin, zevmax
-        ### some already defined colormaps
-        if colorb is True:    colorb = defcolorb(var)
-    else:
-        dimension = 0
- 
-    ###########################
-    ### Get length of time axis
-    if winds:               nt = len(u[:,0,0,0])
-    elif var:  
-        if dimension == 2:  nt = 1
-        else             :  nt = len(field[:,0,0])
+    if zoom:                          [wlon,wlat] = zoomset(wlon,wlat,zoom) 
 
     #########################################
@@ -152,89 +67,45 @@
     elif var:             basename = var
     elif winds:           basename = 'UV'
-    else:                 exit()
-    ###
-    if dimension == 4 or winds:
-        if typefile is 'meso':                      stralt = "_lvl" + str(nvert)
-        elif typefile is 'mesoapi':  
-            zelevel = int(nc.variables['vert'][nvert])
-            if abs(zelevel) < 10000.:   strheight=str(zelevel)+"m"
-            else:                       strheight=str(int(zelevel/1000.))+"km"
-            if 'altitude'       in nc.dimensions:   stralt = "_"+strheight+"-AMR"
-            elif 'altitude_abg' in nc.dimensions:   stralt = "_"+strheight+"-ALS"
-            elif 'bottom_top'   in nc.dimensions:   stralt = "_"+strheight
-            elif 'pressure'     in nc.dimensions:   stralt = "_"+str(zelevel)+"Pa" 
-        else:                                       stralt = ""          
     else:
-        stralt = ""
-    ###
-    basename = basename + stralt
+        print nc.variables                 
+        errormess("please set at least winds or var")
+    basename = basename + getstralt(nc,nvert)  ## can be moved elsewhere for a more generic routine
 
     ##################################
     ### Open a figure and set subplots
     fig = figure()
-    rcParams['font.size'] = 12.
-    if   typefile in ['geo']:  numplot = 1
-    if   numplot > 0:    
-        if   numplot == 4: 
-            sub = 221
-            fig.subplots_adjust(wspace = 0.3, hspace = 0.3)
-            rcParams['font.size'] = int( rcParams['font.size'] * 2. / 3. )
-        elif numplot == 2: 
-            sub = 121
-            fig.subplots_adjust(wspace = 0.35)
-            rcParams['font.size'] = int( rcParams['font.size'] * 3. / 4. )
-        elif numplot == 3: 
-            sub = 131
-            fig.subplots_adjust(wspace = 0.5)
-            rcParams['font.size'] = int( rcParams['font.size'] * 1. / 2. )
-        elif numplot == 6: 
-            sub = 231
-            fig.subplots_adjust(wspace = 0.4, hspace = 0.0)
-            rcParams['font.size'] = int( rcParams['font.size'] * 1. / 2. )
-        elif numplot == 8: 
-            sub = 331 #241
-            fig.subplots_adjust(wspace = 0.3, hspace = 0.3)
-            rcParams['font.size'] = int( rcParams['font.size'] * 1. / 2. )
-        elif numplot == 9:
-            sub = 331
-            fig.subplots_adjust(wspace = 0.3, hspace = 0.3)
-            rcParams['font.size'] = int( rcParams['font.size'] * 1. / 2. )
-        elif numplot == 1:
-            sub = 99999
-        else:
-            print "supported: 1,2,3,4,6,8,9"
-            exit()
-        ### Prepare time loop
-        if nt <= numplot or numplot == 1:  
-            tabrange = [0]
-            numplot = 1
-        else:                          
-            tabrange = range(0,nt,int(nt/numplot))  #nt-1
-            tabrange = tabrange[0:numplot]
-    else: 
-        tabrange = range(0,nt,1)
-        sub = 99999
-    print tabrange
-
+    sub = definesubplot( numplot, fig ) 
+ 
     #################################
     ### Time loop for plotting device
     found_lct = False
-    for i in tabrange:
+    itime = 0  ## could be an argument
+    nplot = 1
+    error = False
+    if itstep is None: itstep = int(24./numplot)
+    while error is False: 
 
        ### Which local time ?
-       ltst = ( interv[0] + 0.5*(wlon[0]+wlon[1])/15.) + i*interv[1]
+       print interv[0], interv[1], itime
+       ltst = ( interv[0] + 0.5*(wlon[0]+wlon[1])/15.) + itime*interv[1]
        ltst = int (ltst * 10) / 10.
        ltst = ltst % 24
 
        ### General plot settings
-       if numplot > 1: 
-           subplot(sub)
+       #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
-       elif numplot == 1:
-           found_lct = True 
-            ### If only one local time is requested (numplot < 0)
+       ### If only one local time is requested (numplot < 0)
        elif numplot <= 0: 
-           if int(ltst) + numplot != 0:		continue
-           else:                  		found_lct = True
+           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
@@ -244,55 +115,53 @@
        #### Contour plot
        if var:
-           if typefile in ['mesoapi','meso']:    what_I_plot = dumpbdy(field[i,:,:])
-           elif typefile in ['geo']:             what_I_plot = field[0,:,:]
-           elif typefile in ['gcm']:             
-               if dimension == 2:                what_I_plot = field[:,:]
-               elif dimension == 3:              what_I_plot = field[i,:,:]
-           palette = get_cmap(name=colorb)
-           #palette.set_over('b', 1.0)
-           #print np.array(x).shape
-           #print np.array(y).shape
-           #print np.array(what_I_plot).shape
-           if not tile:
-               zelevels = np.linspace(zevmin,zevmax)
-               hole = True
-               if not hole:
-                   what_I_plot[ what_I_plot < zevmin ] = zevmin*(1. + 1.e-7)
-                   what_I_plot[ what_I_plot > 9e+35  ] = -9e+35
-                   what_I_plot[ what_I_plot > zevmax ] = zevmax*(1. - 1.e-7)
-               contourf( x, y, what_I_plot, 10, cmap = palette, levels = zelevels )
-           else:    
-               pcolor( x, y, what_I_plot, cmap = palette, vmin=zevmin, vmax=zevmax )
-           #putpoints(m,fig)
-           if var in ['HGT']: pass
-           elif colorb:              
-                              ndiv = 10
-                              colorbar(fraction=0.05,pad=0.1,format=fmtvar(var),\
-                                       ticks=np.linspace(zevmin,zevmax,ndiv+1),\
-                                       extend='max',spacing='proportional')
-                                       # both min max neither
+           what_I_plot, error = reducefield( getfield(nc,var), d4=itime, d3=nvert )
+           if not error: 
+               if typefile in ['mesoapi','meso']:    what_I_plot = dumpbdy(what_I_plot)
+               zevmin, zevmax = calculate_bounds(what_I_plot,vmin=vmin,vmax=vmax)
+               if colorb is True: colorb = defcolorb(var)
+               palette = get_cmap(name=colorb)
+               if not tile:
+                   hole = True
+                   hole = None
+                   if not hole:  what_I_plot = bounds(what_I_plot,zevmin,zevmax)
+                   zelevels = np.linspace(zevmin,zevmax)
+                   contourf( x, y, what_I_plot, zelevels, cmap = palette )
+               else:    
+                   pcolor( x, y, what_I_plot, cmap = palette, vmin=zevmin, vmax=zevmax )
+               #putpoints(m,fig)
+               if var in ['HGT']: pass
+               elif colorb:              
+                         ndiv = 10
+                         colorbar(fraction=0.05,pad=0.1,format=fmtvar(var),\
+                                           ticks=np.linspace(zevmin,zevmax,ndiv+1),\
+                                           extend='max',spacing='proportional')
+                                           # both min max neither
+                  
        ### Vector plot
        if winds:
-           if   typefile in ['mesoapi','meso']:    
-               [vecx,vecy] = [dumpbdy(u[i,nvert,:,:]), dumpbdy(v[i,nvert,:,:])]
-               key = True
-           elif typefile in ['gcm']:               
-               [vecx,vecy] = [        u[i,nvert,:,:] ,         v[i,nvert,:,:] ]
-               key = False
-           if metwind:  [vecx,vecy] = m.rotate_vector(vecx, vecy, lon2d, lat2d)
-           if var == None:        colorvec = definecolorvec(back)
-           else:                  colorvec = definecolorvec(colorb)
-           vectorfield(vecx, vecy,\
-                      x, y, stride=stride, csmooth=2,\
-                      scale=15., factor=300., color=colorvec, key=key)
-                                        #200.         ## or csmooth=stride
-       
+           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,stag=uchar), dumpbdy(vecy,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)
+                                            #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)
+            if addchar:  plottitle = plottitle + addchar + "_LT"+str(ltst)
+            else:        plottitle = plottitle + "_LT"+str(ltst)
        ptitle( plottitle )
-       sub += 1
+       itime += itstep
+       nplot += 1
 
     ##########################################################################
@@ -309,5 +178,5 @@
     else:               zeplot = target + "/" + zeplot  
     ###
-    if found_lct:     makeplotpng(zeplot,pad_inches_value=0.35)   
+    if found_lct:     makeplotpng(zeplot,pad_inches_value=0.35,disp=display)   
     else:             print "Local time not found"
 
@@ -315,4 +184,19 @@
     ### 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:
+        print zelen, len(tab)
+        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
 
 ###########################################################################################
@@ -325,25 +209,26 @@
     from netCDF4 import Dataset
     from myplot import getlschar
-    from numpy import ones
 
     #############################
     ### Get options and variables
     parser = OptionParser()
-    parser.add_option('-f', action='append', dest='namefile',    type="string",  default=None,  help='[NEEDED] name of WRF file (append)')
-    parser.add_option('-l', action='store', dest='nvert',       type="float",   default=0,     help='vertical level (def=0)(-i 2: p,mbar)(-i 3,4: z,km)')
-    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=None,  help='background')
-    parser.add_option('-t', action='store', dest='target',      type="string",  default=None,  help='destination folder')
-    parser.add_option('-s', action='store', dest='stride',      type="int",     default=3,     help='stride vectors (def=3)')
-    parser.add_option('-v', action='append', dest='var',         type="string",  default=None,  help='variable contoured (append)')
-    parser.add_option('-n', action='store', dest='numplot',     type="int",     default=2,     help='number of plots (def=1)(<0: 1 plot of LT -*numplot*)')
-    parser.add_option('-i', action='store', dest='interp',      type="int",     default=None,  help='interpolation method (2: press, 3: z-amr, 4:z-als)')
-    parser.add_option('-c', action='store', dest='colorb',      type="string",  default=True,  help='change colormap')
-    parser.add_option('-x', action='store_false', dest='winds',                 default=True,  help='no wind vectors')
-    parser.add_option('-m', action='append', dest='vmin',        type="float",   default=None,  help='bounding minimum value (append)')    
-    parser.add_option('-M', action='append', dest='vmax',        type="float",   default=None,  help='bounding maximum value (append)') 
-    parser.add_option('-T', action='store_true', dest='tile',                   default=False, help='draw a tiled plot (no blank zone)')
-    parser.add_option('-z', action='store', dest='zoom',        type="float",   default=None,  help='zoom factor in %')
-    parser.add_option('-N', action='store_true', dest='nocall',                 default=False, help='do not recreate api file')
+    parser.add_option('-f', action='append',dest='namefile',    type="string",  default=None,  help='[NEEDED] name of WRF file (append)')
+    parser.add_option('-l', action='store',dest='nvert',        type="float",   default=0,     help='vertical level (def=0)(-i 2: p,mbar)(-i 3,4: z,km)')
+    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=None,  help='background image (def: None)')
+    parser.add_option('-t', action='store',dest='target',       type="string",  default=None,  help='destination folder')
+    parser.add_option('-s', action='store',dest='stride',       type="int",     default=3,     help='stride vectors (def=3)')
+    parser.add_option('-v', action='append',dest='var',         type="string",  default=None,  help='variable contoured (append)')
+    parser.add_option('-n', action='store',dest='numplot',      type="int",     default=2,     help='number of plots (def=1)(<0: 1 plot of LT -*numplot*)')
+    parser.add_option('-i', action='store',dest='interp',       type="int",     default=None,  help='interpolation method (2: press, 3: z-amr, 4:z-als)')
+    parser.add_option('-c', action='store',dest='colorb',       type="string",  default=True,  help='change colormap')
+    parser.add_option('-x', action='store_false',dest='winds',                  default=True,  help='no wind vectors')
+    parser.add_option('-m', action='append',dest='vmin',        type="float",   default=None,  help='bounding minimum value (append)')    
+    parser.add_option('-M', action='append',dest='vmax',        type="float",   default=None,  help='bounding maximum value (append)') 
+    parser.add_option('-T', action='store_true',dest='tile',                    default=False, help='draw a tiled plot (no blank zone)')
+    parser.add_option('-z', action='store',dest='zoom',         type="float",   default=None,  help='zoom factor in %')
+    parser.add_option('-N', action='store_true',dest='nocall',                  default=False, help='do not recreate api file')
+    parser.add_option('-d', action='store_false',dest='display',                default=True,  help='do not pop up created images')
+    parser.add_option('-e', action='store',dest='itstep',       type="int",     default=None,  help='stride time (def=4)')
     #parser.add_option('-V', action='store', dest='comb',        type="float",   default=None,  help='a defined combination of variables')
     (opt,args) = parser.parse_args()
@@ -352,30 +237,17 @@
         exit()
     print "Options:", opt
-   
+
+    listvar = '' 
     if opt.var is None:
-        pass
+        zerange = [-999999]
     else:
-        listvar = ''
         zelen = len(opt.var)
-        if zelen == 1: listvar = opt.var[0] + ','
-        else         : 
-            for jjj in range(zelen): listvar += opt.var[jjj] + ','
+        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]
-        if opt.vmin:
-            if zelen != len(opt.vmin):
-                print "not enough or too much vmin values... setting same values all variables"
-                vmintab = ones(zelen) * opt.vmin[0]
-            else:
-                vmintab = opt.vmin
-        else:
-            vmintab = None
-        if opt.vmax:
-            if zelen != len(opt.vmax):
-                print "not enough or too much vmax values... setting same values all variables"
-                vmaxtab = ones(zelen) * opt.vmax[0]
-            else:
-                vmaxtab = opt.vmax
-        else:
-            vmaxtab = None
+        vmintab = adjust_length (opt.vmin, zelen)  
+        vmaxtab = adjust_length (opt.vmax, zelen)
 
     for i in range(len(opt.namefile)):
@@ -395,6 +267,6 @@
             else                    :  zefields = ''
             ### var or no var
-            if opt.var is None      :  pass
-            elif zefields == ''     :  zefields = listvar 
+            #if opt.var is None      :  pass
+            if zefields == ''       :  zefields = listvar 
             else                    :  zefields = zefields + "," + listvar 
             print zefields
@@ -409,23 +281,26 @@
             zelevel = 0 ## so that zelevel could play again the role of nvert
 
-        if opt.var is None:  
+        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=opt.var,numplot=opt.numplot,colorb=opt.colorb,winds=opt.winds,\
-               addchar=lschar,interv=[zehour,zehourin],vmin=opt.vmin,vmax=opt.vmax,tile=opt.tile,zoom=opt.zoom)
+                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)
             print 'Done: '+name
-        else:        
-            for jjj in range(len(opt.var)):  
-                if vmintab: argvmin = vmintab[jjj]
-                else:       argvmin = None
-                if vmaxtab: 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=opt.var[jjj],numplot=opt.numplot,colorb=opt.colorb,winds=opt.winds,\
-                   addchar=lschar,interv=[zehour,zehourin],vmin=argvmin,vmax=argvmax,tile=opt.tile,zoom=opt.zoom)
-                print 'Done: '+name
     
         #########################################################