Changeset 392


Ignore:
Timestamp:
Nov 17, 2011, 2:17:48 AM (13 years ago)
Author:
aslmd
Message:

PYTHON: quite a lot of modifs and tests. extended multivar subplot capabilities + cosmetic changes + general cleaning. corrected some stuff that were not working with mesoscale. and there is now a common script for meso and gcm svn status -uq! it is named pp.py [stands for: planetoplot]. put some examples in README file. testing is probably a bit necessary for most complex options.

Location:
trunk/UTIL/PYTHON
Files:
1 deleted
4 edited
2 moved

Legend:

Unmodified
Added
Removed

  • trunk/UTIL/PYTHON/README.PP

    r391 r392  
     1**************************************
     2**************************************
     3**************************************
     4    PLANETOPLOT TUTORIAL EXAMPLES
     5**************************************
     6         Authors : AC + AS
     7**************************************
     8  DON'T FORGET YOUR BEST FRIEND IS
     9     pp.py -h [or] pp.py --help
     10**************************************
     11**************************************
     12**************************************
    113
    2 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    314
    4 PYTHON COMMAND LINE PLOTS EXAMPLES
     15*************************************
     16MAPMODE 1
     17MAPPING MODE
     18SIMPLE EXAMPLES on a SAMPLE GCM FILE
     19*************************************
    520
    6 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     21Goal: The simplest, most minimal example. Mapping near-surface wind vectors over topography.
     22pp.py -f diagfired.nc
    723
     24Goal: I would like finer contours.
     25pp.py -f diagfired.nc --div 30
     26
     27Goal: I would like more vectors [i.e. lower the stride].
     28pp.py -f diagfired.nc --div 30 -s 1
     29
     30Goal: I want to get rid of those damn vectors.
     31pp.py -f diagfired.nc -x
     32
     33Goal: I want to map a given field (surface temperature).
     34pp.py -f diagfired.nc -x -v tsurf
     35
     36Goal: I want to map two fields next to one another (topography and tauice).
     37pp.py -f diagfired.nc -x -v phisinit,tauice
     38
     39Goal: I want to map two fields, tauice shaded, topography contoured, same plot.
     40pp.py -f diagfired.nc -x -v tauice -w phisinit
     41
     42Goal: I want to map a field but projected on the sphere.
     43pp.py -f diagfired.nc -x -v tauice -p ortho
     44
     45Goal: I want to redefine the minimum and maximum values shown.
     46pp.py -f diagfired.nc -x -v tauice -m 0.2 -M 0.9
     47
     48Goal: I want to insert holes wherever values are lower than 0.2 and higher than 0.9
     49pp.py -f diagfired.nc -x -v tauice -m 0.2 -M 0.9 -H
     50
     51Goal: I want to fill holes with an background image of Mars [you have to be connected to Internet]
     52pp.py -f diagfired.nc -x -v tauice -m 0.2 -M 0.9 -H -b vishires
     53
     54Goal: I want the same map, but projected on the sphere
     55pp.py -f diagfired.nc -x -v tauice -m 0.2 -M 0.9 -H -b vishires -p ortho
     56
     57Goal: I want the same map, but projected with north polar stereographic view
     58pp.py -f diagfired.nc -x -v tauice -m 0.2 -M 0.9 -H -b vishires -p npstere
     59
     60Goal: I want to save this in PNG format
     61pp.py -f diagfired.nc -x -v tauice -m 0.2 -M 0.9 -H -b vishires -p ortho -S png
     62
     63Goal: I want to plot results from two simulation files next to one another
     64pp.py -f diagfired.nc,diagfired.nc -x -v tsurf
    865
    966***********************************************************************************

  • trunk/UTIL/PYTHON/api_wrapper.py

    r351 r392  
    3737        onelevel = -99999.
    3838
    39     #print input_name, output_name
     39    print input_name, output_name
    4040
    4141    if nocall:     pass

  • trunk/UTIL/PYTHON/myplot.py

    r391 r392  
    4343    elif 'HGT_M' in nc.variables:    typefile = 'geo'
    4444    #else:                            errormess("whatkindfile: typefile not supported.")
    45     else: typefile = 'gcm' # for lslin-ed files from gcm
     45    else:                            typefile = 'gcm' # for lslin-ed files from gcm
    4646    return typefile
    4747
    … …  
    5050    ## this allows to get much faster (than simply referring to nc.variables[var])
    5151    dimension = len(nc.variables[var].dimensions)
    52     print "   Opening variable with", dimension, "dimensions ..."
     52    #print "   Opening variable with", dimension, "dimensions ..."
    5353    if dimension == 2:    field = nc.variables[var][:,:]
    5454    elif dimension == 3:  field = nc.variables[var][:,:,:]
    … …  
    6565    shape = np.array(input).shape
    6666    #print 'd1,d2,d3,d4: ',d1,d2,d3,d4
    67     print 'dim,shape: ',dimension,shape
     67    print 'IN  REDUCEFIELD dim,shape: ',dimension,shape
    6868    output = input
    6969    error = False
    … …  
    9898        elif max(d1) >= shape[3]: error = True
    9999        elif d4 is not None and d3 is not None and d2 is not None and d1 is not None:
    100              output = mean(input[d4,:,:,:],axis=0); output = reduce_zaxis(output[d3,:,:],ax=0,yint=yint,vert=alt[d3]); output = mean(output[d2,:],axis=0); output = mean(output[d1],axis=0)
     100             output = mean(input[d4,:,:,:],axis=0)
     101             output = reduce_zaxis(output[d3,:,:],ax=0,yint=yint,vert=alt,indice=d3)
     102             output = mean(output[d2,:],axis=0)
     103             output = mean(output[d1],axis=0)
    101104        elif d4 is not None and d3 is not None and d2 is not None:
    102             output = mean(input[d4,:,:,:],axis=0); output = reduce_zaxis(output[d3,:,:],ax=0,yint=yint,vert=alt[d3]); output = mean(output[d2,:],axis=0)
     105             output = mean(input[d4,:,:,:],axis=0)
     106             output = reduce_zaxis(output[d3,:,:],ax=0,yint=yint,vert=alt,indice=d3)
     107             output = mean(output[d2,:],axis=0)
    103108        elif d4 is not None and d3 is not None and d1 is not None:
    104             output = mean(input[d4,:,:,:],axis=0); output = reduce_zaxis(output[d3,:,:],ax=0,yint=yint,vert=alt[d3]); output = mean(output[:,d1],axis=1)
     109             output = mean(input[d4,:,:,:],axis=0)
     110             output = reduce_zaxis(output[d3,:,:],ax=0,yint=yint,vert=alt,indice=d3)
     111             output = mean(output[:,d1],axis=1)
    105112        elif d4 is not None and d2 is not None and d1 is not None:
    106             output = mean(input[d4,:,:,:],axis=0); output = mean(output[:,d2,:],axis=1); output = mean(output[:,d1],axis=1)
     113             output = mean(input[d4,:,:,:],axis=0)
     114             output = mean(output[:,d2,:],axis=1)
     115             output = mean(output[:,d1],axis=1)
    107116        elif d3 is not None and d2 is not None and d1 is not None:
    108             output = reduce_zaxis(input[:,d3,:,:],ax=1,yint=yint,vert=alt[d3]); output = mean(output[:,d2,:],axis=1); output = mean(output[:,d1],axis=1)
    109         elif d4 is not None and d3 is not None:  output = mean(input[d4,:,:,:],axis=0); output = reduce_zaxis(output[d3,:,:],ax=0,yint=yint,vert=alt[d3])
    110         elif d4 is not None and d2 is not None:  output = mean(input[d4,:,:,:],axis=0); output = mean(output[:,d2,:],axis=1)
    111         elif d4 is not None and d1 is not None:  output = mean(input[d4,:,:,:],axis=0); output = mean(output[:,:,d1],axis=2)
    112         elif d3 is not None and d2 is not None:  output = reduce_zaxis(input[:,d3,:,:],ax=1,yint=yint,vert=alt[d3]); output = mean(output[:,d2,:],axis=1)
    113         elif d3 is not None and d1 is not None:  output = reduce_zaxis(input[:,d3,:,:],ax=1,yint=yint,vert=alt[d3]); output = mean(output[:,:,d1],axis=0)
    114         elif d2 is not None and d1 is not None:  output = mean(input[:,:,d2,:],axis=2); output = mean(output[:,:,d1],axis=2)
    115         elif d1 is not None:                     output = mean(input[:,:,:,d1],axis=3)
    116         elif d2 is not None:                     output = mean(input[:,:,d2,:],axis=2)
    117         elif d3 is not None:                     output = reduce_zaxis(input[:,d3,:,:],ax=0,yint=yint,vert=alt[d3])
    118         elif d4 is not None:                     output = mean(input[d4,:,:,:],axis=0)
     117             output = reduce_zaxis(input[:,d3,:,:],ax=1,yint=yint,vert=alt,indice=d3)
     118             output = mean(output[:,d2,:],axis=1)
     119             output = mean(output[:,d1],axis=1)
     120        elif d4 is not None and d3 is not None: 
     121             output = mean(input[d4,:,:,:],axis=0)
     122             output = reduce_zaxis(output[d3,:,:],ax=0,yint=yint,vert=alt,indice=d3)
     123        elif d4 is not None and d2 is not None: 
     124             output = mean(input[d4,:,:,:],axis=0)
     125             output = mean(output[:,d2,:],axis=1)
     126        elif d4 is not None and d1 is not None: 
     127             output = mean(input[d4,:,:,:],axis=0)
     128             output = mean(output[:,:,d1],axis=2)
     129        elif d3 is not None and d2 is not None: 
     130             output = reduce_zaxis(input[:,d3,:,:],ax=1,yint=yint,vert=alt,indice=d3)
     131             output = mean(output[:,d2,:],axis=1)
     132        elif d3 is not None and d1 is not None: 
     133             output = reduce_zaxis(input[:,d3,:,:],ax=1,yint=yint,vert=alt,indice=d3)
     134             output = mean(output[:,:,d1],axis=0)
     135        elif d2 is not None and d1 is not None: 
     136             output = mean(input[:,:,d2,:],axis=2)
     137             output = mean(output[:,:,d1],axis=2)
     138        elif d1 is not None:        output = mean(input[:,:,:,d1],axis=3)
     139        elif d2 is not None:        output = mean(input[:,:,d2,:],axis=2)
     140        elif d3 is not None:        output = reduce_zaxis(input[:,d3,:,:],ax=0,yint=yint,vert=alt,indice=d3)
     141        elif d4 is not None:        output = mean(input[d4,:,:,:],axis=0)
    119142    dimension = np.array(output).ndim
    120143    shape = np.array(output).shape
    121     print 'dim,shape: ',dimension,shape
     144    print 'OUT REDUCEFIELD dim,shape: ',dimension,shape
    122145    return output, error
    123146
    124 ## Author: AC
    125 
    126 def reduce_zaxis (input,ax=None,yint=False,vert=None):
     147## Author: AC + AS
     148def reduce_zaxis (input,ax=None,yint=False,vert=None,indice=None):
    127149    from mymath import max,mean
    128150    from scipy import integrate
    129     if yint:
     151    if yint and vert is not None and indice is not None:
    130152       if type(input).__name__=='MaskedArray':
    131153          input.set_fill_value([np.NaN])
    132           output = integrate.trapz(input.filled(),x=vert,axis=ax)
     154          output = integrate.trapz(input.filled(),x=vert[indice],axis=ax)
    133155       else:
    134           output = integrate.trapz(input.filled(),x=vert,axis=ax)
     156          output = integrate.trapz(input.filled(),x=vert[indice],axis=ax)
    135157    else:
    136158       output = mean(input,axis=ax)
    … …  
    257279            if cen_lat > 10.:   
    258280                proj="npstere"
    259                 print "NP stereographic polar domain"
     281                #print "NP stereographic polar domain"
    260282            else:           
    261283                proj="spstere"
    262                 print "SP stereographic polar domain"
     284                #print "SP stereographic polar domain"
    263285        elif map_proj == 1:
    264             print "lambert projection domain"
     286            #print "lambert projection domain"
    265287            proj="lcc"
    266288        elif map_proj == 3:
    267             print "mercator projection"
     289            #print "mercator projection"
    268290            proj="merc"
    269291        else:
    … …  
    460482        elif wlat[1] >= 0.:    blat = wlat[0]
    461483        elif wlat[0] <= 0.:    blat = wlat[1]
    462     print "blat ", blat
     484    #print "blat ", blat
    463485    h = 50.  ## en km
    464486    radius = 3397200.
    … …  
    484506    #    step = np.min([5.,step])
    485507    #    steplon = step
    486     print step
    487508    m.drawmeridians(np.r_[-180.:180.:steplon], labels=[0,0,0,1], color='grey', fontsize=fontsizemer)
    488509    m.drawparallels(np.r_[-90.:90.:step], labels=[1,0,0,0], color='grey', fontsize=fontsizemer)
    … …  
    534555    ###
    535556    if zevmin < 0. and min(fieldcalc) > 0.: zevmin = 0.
    536     print "field ", min(fieldcalc), max(fieldcalc)
    537     print "bounds ", zevmin, zevmax
     557    print "BOUNDS field ", min(fieldcalc), max(fieldcalc)
     558    print "BOUNDS adopted ", zevmin, zevmax
    538559    return zevmin, zevmax
    539560
    … …  
    545566    if zevmin < 0: what_I_plot[ what_I_plot < zevmin*(1. - 1.e-7) ] = zevmin*(1. - 1.e-7)
    546567    else:          what_I_plot[ what_I_plot < zevmin*(1. + 1.e-7) ] = zevmin*(1. + 1.e-7)
    547     print "new min ", min(what_I_plot)
     568    print "NEW MIN ", min(what_I_plot)
    548569    what_I_plot[ what_I_plot > 9e+35  ] = -9e+35
    549570    what_I_plot[ what_I_plot > zevmax ] = zevmax
    550     print "new max ", max(what_I_plot)
    551    
     571    print "NEW MAX ", max(what_I_plot)
    552572    return what_I_plot
    553573
    … …  
    556576    from mymath import max,min
    557577    lim = 0.15*0.5*(abs(max(what_I_plot))+abs(min(what_I_plot)))
    558     print "on vire en dessous de ", lim
     578    print "NO PLOT BELOW VALUE ", lim
    559579    what_I_plot [ abs(what_I_plot) < lim ] = 1.e40
    560580    return what_I_plot
    … …  
    566586    [wlon,wlat] = [ [wlon[0]+zoom*dlon/100.,wlon[1]-zoom*dlon/100.],\
    567587                    [wlat[0]+zoom*dlat/100.,wlat[1]-zoom*dlat/100.] ]
    568     print "zoom %",zoom,wlon,wlat
     588    print "ZOOM %",zoom,wlon,wlat
    569589    return wlon,wlat
    570590
    … …  
    642662#W --> spectral ou jet
    643663#spectral BrBG RdBu_r
    644     print "predefined colorbars"
     664    #print "predefined colorbars"
    645665    if whichone not in whichcolorb:
    646666        whichone = "def"
    … …  
    813833      if axis[imin] == -180 and axis[imax] == 180:
    814834         zeindex = zeindex[0:len(zeindex)-1]
    815          print "whole longitude averaging asked, so last point is not taken into account."
     835         print "INFO: whole longitude averaging asked, so last point is not taken into account."
    816836  return zeindex
    817837     
    … …  
    828848   shape = what_I_plot.shape
    829849   if indextime is None:
    830       print "axis is time"
     850      print "AXIS is time"
    831851      x = time
    832852      count = count+1
    833853   if indexlon is None:
    834       print "axis is lon"
     854      print "AXIS is lon"
    835855      if count == 0: x = lon
    836856      else: y = lon
    837857      count = count+1
    838858   if indexlat is None:
    839       print "axis is lat"
     859      print "AXIS is lat"
    840860      if count == 0: x = lat
    841861      else: y = lat
    842862      count = count+1
    843863   if indexvert is None and dim0 is 4:
    844       print "axis is vert"
     864      print "AXIS is vert"
    845865      if vertmode == 0: # vertical axis is as is (GCM grid)
    846866         if count == 0: x=range(len(vert))
    … …  
    853873   x = array(x)
    854874   y = array(y)
    855    print "what_I_plot.shape", what_I_plot.shape
    856    print "x.shape, y.shape", x.shape, y.shape
     875   print "CHECK: what_I_plot.shape", what_I_plot.shape
     876   print "CHECK: x.shape, y.shape", x.shape, y.shape
    857877   if len(shape) == 1:
    858878      print shape[0]
    … …  
    864884      if shape[1] == len(y) and shape[0] == len(x) and shape[0] != shape[1]:
    865885         what_I_plot = swapaxes(what_I_plot,0,1)
    866          print "swapaxes", what_I_plot.shape, shape
     886         print "INFO: swapaxes", what_I_plot.shape, shape
    867887         #x0 = x
    868888         #x = y
    … …  
    896916    if slon is None and slat is None:
    897917       mapmode = 1 # in this case we plot a map, with the given projection
    898     #elif proj is not None:
    899     #   print "WARNING: you specified a", proj,\
    900     #   "projection but asked for slices", slon,"in longitude and",slat,"in latitude"
    901     print "mapmode: ", mapmode
    902918
    903919    return nlon, nlat, nvert, ntime, mapmode, nslices

  • trunk/UTIL/PYTHON/myscript.py

    r388 r392  
    3434    parser.add_option('-s', '--stride', action='store',dest='ste',       type="int",     default=3,     help='stride vectors [3]')
    3535    parser.add_option('-x', '--no-vect',action='store_false',dest='winds',               default=True,  help='no wind vectors')
    36     parser.add_option('-n', '--num',    action='store',dest='num',       type="int",     default=None,  help='plot number (<0: plot LT -*numplot*) [2]')
     36    parser.add_option('-n', '--num',    action='store',dest='num',       type="int",     default=None,  help='plot number (<0: plot LT -*numplot*) [1]')
    3737    parser.add_option('-z', '--zoom',   action='store',dest='zoom',      type="float",   default=None,  help='zoom factor in %')
    3838    parser.add_option('-e', '--itstep', action='store',dest='it',        type="int",     default=None,  help='stride time [4]')
    … …  
    4444    parser.add_option('--lon',          action='append',dest='slon',   type="string",  default=None, help='slices along lon. 2 comma-separated values: averaging')
    4545    parser.add_option('--vert',         action='append',dest='svert',  type="string",  default=None, help='slices along vert. 2 comma-separated values: averaging')
    46     parser.add_option('--column',         action='store_true',dest='column',               default=False,help='changes the function of --vert z1,z2 from MEAN to INTEGRATE along z')
     46    parser.add_option('--column',         action='store_true',dest='column',           default=False,help='changes the function of --vert z1,z2 from MEAN to INTEGRATE along z')
    4747    parser.add_option('--time',         action='append',dest='stime',  type="string",  default=None, help='slices along time. 2 comma-separated values: averaging')
    4848    parser.add_option('--xmax',         action='store',dest='xmax',    type="float",   default=None, help='max value for x-axis in contour-plots [max(xaxis)]')

  • trunk/UTIL/PYTHON/planetoplot.py

    r391 r392  
    55### A. Spiga     -- LMD -- 06~09/2011 -- General building and mapping capabilities
    66### T. Navarro   -- LMD -- 10~11/2011 -- Improved use for GCM and added sections + 1Dplot capabilities
    7 ### A. Colaitis  -- LMD --            -- Mostly minor improvements and inter-plot operation capabilities + zrecast interpolation for gcm
     7### A. Colaitis  -- LMD --    11/2011 -- Mostly minor improvements and inter-plot operation capabilities + zrecast interpolation for gcm
     8### A. Spiga     -- LMD --    11/2011 -- Extended multivar subplot capabilities + cosmetic changes + general cleaning and tests
    89
    910def planetoplot (namefiles,\
    … …  
    6364                       fmtvar,definecolorvec,defcolorb,getprefix,putpoints,calculate_bounds,errormess,definesubplot,\
    6465                       zoomset,getcoorddef,getwinddef,whatkindfile,reducefield,bounds,getstralt,getfield,smooth,nolow,\
    65                        getname,localtime,polarinterv,getsindex,define_axis,determineplot
     66                       getname,localtime,polarinterv,getsindex,define_axis,determineplot,readslices
    6667    from mymath import deg,max,min,mean,get_tsat
    6768    import matplotlib as mpl
    68     from matplotlib.pyplot import contour,contourf, subplot, figure, rcParams, savefig, colorbar, pcolor, show, plot, clabel
     69    from matplotlib.pyplot import contour,contourf, subplot, figure, rcParams, savefig, colorbar, pcolor, show, plot, clabel, title
    6970    from matplotlib.cm import get_cmap
    7071    import numpy as np
    … …  
    7273
    7374    ################################
     75    ### Preliminary stuff
     76    ################################
     77    print "********************************************"
     78    print "********** WELCOME TO PLANETOPLOT **********"
     79    print "********************************************"
     80    if not isinstance(namefiles, np.ndarray): namefiles = [namefiles]
     81    if not isinstance(var, np.ndarray):       var = [var]
     82    if ope is not None and len(var) > 1:      errormess("you are using an operation... please set only one var !")
     83
     84    ################################
    7485    ### Which plot needs to be done?
     86    ################################
    7587    nlon, nlat, nvert, ntime, mapmode, nslices = determineplot(slon, slat, svert, stime)
    76     if mapmode == 0: winds=False
    77     if not isinstance(namefiles, np.ndarray): namefiles = [namefiles]
    78     zelen = len(namefiles)
    79     if numplot == None:  numplot = zelen*nslices
    80     print "len(namefiles), nslices, numplot: ", zelen, nslices, numplot
    81     if fileref is not None:
    82         all_var  = [[]]*(zelen+2)
    83     else:
    84         all_var   = [[]]*zelen
    85         all_var2  = [[]]*zelen
    86         all_title = [[]]*zelen
    87    
    88     #########################
    89     ### Loop over the files initially separated by comas to be plotted on the same figure
     88    if mapmode == 0:       winds=False
     89    elif mapmode == 1:     
     90        if svert is None:  svert = readslices(str(level)) ; nvert=1
     91    zelen = len(namefiles)*len(var)
     92    if numplot is None:  numplot = zelen*nslices
     93    print "********** FILES, SLICES, VARS, TOTAL PLOTS: ", len(namefiles), nslices, len(var), numplot
     94    print "********** MAPMODE: ", mapmode
     95    if fileref is not None:       all_var  = [[]]*(zelen+2) ;  all_varname = [[]]*(zelen+2)
     96    else:                         all_var  = [[]]*zelen     ;  all_var2  = [[]]*zelen ;  all_title = [[]]*zelen ;  all_varname = [[]]*zelen
     97 
     98    #################################################################################################
     99    ### Loop over the files + vars initially separated by commas to be plotted on the same figure ###
     100    #################################################################################################
    90101    k = 0
    91102    firstfile = True
    92     for namefile in namefiles:
    93       print namefile
     103    for nnn in range(len(namefiles)):
     104     for vvv in range(len(var)):
     105
     106      print "********** LOOP..... THIS IS SUBPLOT NUMBER.....",k
     107
    94108      ######################
    95109      ### Load NETCDF object
     110      namefile = namefiles[nnn] ; print "********** THE NAMEFILE IS....", namefile
    96111      nc  = Dataset(namefile)
    97112
    … …  
    100115      ### ... TYPEFILE
    101116      typefile = whatkindfile(nc)                                 
    102       if firstfile:
    103          typefile0 = typefile
    104       elif typefile != typefile0:
    105          errormess("Not the same kind of files !", [typefile0, typefile])
     117      if firstfile:                 typefile0 = typefile
     118      elif typefile != typefile0:   errormess("Not the same kind of files !", [typefile0, typefile])
    106119      ### ... VAR
    107       varname=var
    108       if var not in nc.variables: var = False
     120      varname=var[vvv]
     121      print "********** THE VAR IS....",varname, var2
     122      if varname not in nc.variables: varname = False
    109123      ### ... WINDS
    110124      if winds:                                                   
    111125         [uchar,vchar,metwind] = getwinddef(nc)             
    112126         if uchar == 'not found': winds = False
    113       if not var and not winds: errormess("please set at least winds or var",printvar=nc.variables)
     127      if not varname and not winds: errormess("please set at least winds or var",printvar=nc.variables)
    114128      ### ... COORDINATES, could be moved below
    115129      [lon2d,lat2d] = getcoorddef(nc)                       
    116130      ### ... PROJECTION
    117       if proj == None:   proj = getproj(nc)                  
     131      if proj == None:   proj = getproj(nc)                 
    118132
    119133##########################################################
    … …  
    131145      elif typefile in ['meso','mesoapi']:
    132146          if mapmode == 0:
    133               if var in ['PHTOT','W']:    vertdim='BOTTOM-TOP_PATCH_END_STAG'
    134               else:                       vertdim='BOTTOM-TOP_PATCH_END_UNSTAG'
    135               if var in ['V']:  latdim='SOUTH-NORTH_PATCH_END_STAG'
    136               else:             latdim='SOUTH-NORTH_PATCH_END_UNSTAG'
    137               if var in ['U']:  londim='WEST-EAST_PATCH_END_STAG'
    138               else:             londim='WEST-EAST_PATCH_END_UNSTAG'
     147              if varname in ['PHTOT','W']:    vertdim='BOTTOM-TOP_PATCH_END_STAG'
     148              else:                           vertdim='BOTTOM-TOP_PATCH_END_UNSTAG'
     149              if varname in ['V']:  latdim='SOUTH-NORTH_PATCH_END_STAG'
     150              else:                 latdim='SOUTH-NORTH_PATCH_END_UNSTAG'
     151              if varname in ['U']:  londim='WEST-EAST_PATCH_END_STAG'
     152              else:                 londim='WEST-EAST_PATCH_END_UNSTAG'
    139153              lon = np.arange(0,getattr(nc,londim),1)
    140154              lat = np.arange(0,getattr(nc,latdim),1)
    141               if vertmode is None: vertmode=0
    142               if vertmode == 0:   vert = np.arange(0,getattr(nc,vertdim),1)
    143               else:               vert = nc.variables["vert"][:]
     155              if vertmode is None:  vertmode=0
     156              if vertmode == 0:     vert = np.arange(0,getattr(nc,vertdim),1)
     157              else:                 vert = nc.variables["vert"][:]
    144158              time = np.arange(0,len(nc.variables["Times"]),1)
     159          else:
     160              lon=None ; lat=None ; vert=None ; time=None
    145161       #if firstfile:
    146162       #   lat0 = lat
    … …  
    160176         else:                             [wlon,wlat] = simplinterv(lon2d,lat2d)
    161177         if zoom:                          [wlon,wlat] = zoomset(wlon,wlat,zoom)
    162 
    163          #########################################
    164          ### Name for title and graphics save file
    165          basename = getname(var=var,winds=winds,anomaly=anomaly)
     178         ### ... NAME FOR TITLE and FILE (have to do something for multiple vars...)
     179         basename = getname(var=varname,winds=winds,anomaly=anomaly)
    166180         basename = basename + getstralt(nc,level)  ## can be moved elsewhere for a more generic routine
    167181
    168       print "var, var2: ", var, var2
    169       if varname in ["temp","t","T_nadir_nit","T_nadir_day"] and var and tsat:
    170           print "computing Tsat-T, I ASSUME Z-AXIS IS PRESSURE"
    171           tt=getfield(nc,var)
    172           if type(tt).__name__=='MaskedArray':
    173              tt.set_fill_value([np.NaN])
    174              tinput=tt.filled()
    175           else:
    176              tinput=tt
     182      ##### SPECIFIC
     183      if varname in ["temp","t","T_nadir_nit","T_nadir_day"] and tsat:
     184          tt=getfield(nc,varname) ; print "computing Tsat-T, I ASSUME Z-AXIS IS PRESSURE"
     185          if type(tt).__name__=='MaskedArray':  tt.set_fill_value([np.NaN]) ; tinput=tt.filled()
     186          else:                                 tinput=tt
    177187          all_var[k]=get_tsat(vert,tinput,zlon=lon,zlat=lat,zalt=vert,ztime=time)
    178188      else:
    179          if var: all_var[k] = getfield(nc,var)
     189      ##### GENERAL STUFF HERE
     190          all_var[k] = getfield(nc,varname)
     191          all_varname[k] = varname
    180192      if var2: all_var2[k] = getfield(nc,var2)
    181193     
    182       print "k", k
    183       print "all_var[k].shape", all_var[k].shape
     194      print "********** all_var[k].shape", all_var[k].shape
    184195      k += 1
    185196      firstfile = False
    … …  
    188199    ##################################
    189200    ### Operation on files
     201    ### ... for the moment only OK when 1 var only
    190202    if ope is not None:
    191        if var not in nc.variables: var = False
    192        if var:
    193203             print ope
    194204             if ope in ["-","+"]:
    195                 if fileref is not None:   all_var[k] = getfield(Dataset(fileref),var)
     205                if fileref is not None:   all_var[k] = getfield(Dataset(fileref),all_varname[k-1]) ; all_varname[k] = all_varname[k-1]
    196206                else:                     errormess("fileref is missing!")
    197207                if ope == "-":     all_var[k+1]= all_var[k-1] - all_var[k]
    198208                elif ope == "+":   all_var[k+1]= all_var[k-1] + all_var[k]
     209                all_varname[k+1] = all_varname[k]
    199210                numplot = numplot+2
    200211             elif ope in ["cat"]:
    … …  
    219230    elif numplot <= 0:                 itstep = 1
    220231   
    221     #for nplot in range(numplot):
     232    print "********************************************"
    222233    while error is False:
    223        print "nplot", nplot
    224        print error
     234       print "********** nplot", nplot, "itime",itime,"error",error
    225235       
    226236       ### Which local time ?
    … …  
    249259
    250260####################################################################
    251        if typefile in ['meso','mesoapi'] and mapmode == 1:
    252                indextime = itime
    253                indexlon = None
    254                indexlat = None
    255                indexvert = level  ## ou svert ???
    256                nlon = 1
    257                nlat = 1
    258                nvert = 1
    259                ntime = 1
     261       ## get all indexes to be taken into account for this subplot and then reduce field
     262       ## We plot 1) all lon slices 2) all lat slices 3) all vert slices 4) all time slices and then go to the next slice
     263       indexlon  = getsindex(slon,(nplot-1)%nlon,lon)
     264       indexlat  = getsindex(slat,((nplot-1)//nlon)%nlat,lat)
     265       indexvert = getsindex(svert,((nplot-1)//(nlon*nlat))%nvert,vert)
     266       if mapmode == 1 and stime is None:
     267           indextime = itime
     268           if len(var) != 1 or len(namefiles) != 1: indextime = first
    260269       else:
    261            ## get all indexes to be taken into account for this subplot and then reduce field
    262            ## We plot 1) all lon slices 2) all lat slices 3) all vert slices 4) all time slices and then go to the next slice
    263            indexlon  = getsindex(slon,(nplot-1)%nlon,lon)
    264            indexlat  = getsindex(slat,((nplot-1)//nlon)%nlat,lat)
    265            indexvert = getsindex(svert,((nplot-1)//(nlon*nlat))%nvert,vert) 
    266270           indextime = getsindex(stime,((nplot-1)//(nlon*nlat*nvert))%ntime,time) 
    267271
    268272       if fileref is not None:
    269            index_f = ((nplot-1)//(nlon*nlat*nvert*ntime))%(len(namefiles)+2)
     273           index_f = ((nplot-1)//(nlon*nlat*nvert*ntime))%(len(namefiles)+2)  ## OK only 1 var, see test in the beginning
    270274       else:
    271            index_f = ((nplot-1)//(nlon*nlat*nvert*ntime))%len(namefiles)
    272 
    273        print nlon, nlat, nvert, ntime 
    274        print slon, slat, svert, stime       
    275        print index_f
    276        print "index lon,lat,vert,time", max(indexlon), max(indexlat), max(indexvert), max(indextime)
     275           yeah = len(namefiles)*len(var)
     276           index_f = ((nplot-1)//(nlon*nlat*nvert*ntime))%yeah
     277       #print nlon, nlat, nvert, ntime 
     278       #print slon, slat, svert, stime       
     279       #print index_f
     280       #print "index lon,lat,vert,time", max(indexlon), max(indexlat), max(indexvert), max(indextime)
     281       #if varname != all_varname[index_f]:  indextime = first
    277282####################################################################
     283
     284       ticks = ndiv + 1
    278285
    279286       #### Contour plot
    … …  
    285292                  if typefile in ['mesoapi','meso']:    what_I_plot = dumpbdy(what_I_plot,6)
    286293              zevmin, zevmax = calculate_bounds(what_I_plot)
    287               zelevels = np.linspace(zevmin,zevmax,num=ndiv+1) #20)
     294              zelevels = np.linspace(zevmin,zevmax,ticks) #20)
    288295              if var2 == 'HGT':  zelevels = np.arange(-10000.,30000.,2000.)
    289296              if mapmode == 0:
    … …  
    293300              if len(what_I_plot.shape) is 2:
    294301                  cs = contour(x,y,what_I_plot, zelevels, colors='k', linewidths = 1 ) #0.33 colors='w' )# , alpha=0.5)
    295                   if typefile in ['gcm']: clabel(cs,fmt = '%1.2e')
     302                  #if typefile in ['gcm']: clabel(cs,fmt = '%1.2e')
    296303              ### If we plot a 1-D field
    297304              elif len(what_I_plot.shape) is 1:
    … …  
    301308
    302309       #### Shaded plot
    303        if var:
     310       varname = all_varname[index_f]
     311       if varname:
    304312           what_I_plot, error = reducefield(all_var[index_f], d4=indextime, d1=indexlon, d2=indexlat , d3=indexvert , yint=yintegral, alt=vert)
    305313           what_I_plot = what_I_plot*mult
    306314           if not error:
    307                fvar = var
     315               fvar = varname
    308316               ###
    309317               if anomaly:
    … …  
    329337                     if not hole: what_I_plot = bounds(what_I_plot,zevmin,zevmax)
    330338                     #zelevels = np.linspace(zevmin*(1. + 1.e-7),zevmax*(1. - 1.e-7)) #,num=20)
    331                      zelevels = np.linspace(zevmin,zevmax,num=2.*ndiv+1)
     339                     zelevels = np.linspace(zevmin,zevmax,num=ticks)
    332340                     #contourf(what_I_plot, zelevels, cmap = palette )
    333341
    … …  
    336344                     elif mapmode == 0:
    337345                       contourf( x, y, what_I_plot, zelevels, cmap = palette)
     346
    338347                     zxmin, zxmax = xaxis
    339348                     zymin, zymax = yaxis
    340                      if zxmin is not None:
    341                         mpl.pyplot.xlim(xmin=zxmin)
    342                      if zxmax is not None:
    343                         mpl.pyplot.xlim(xmax=zxmax)
    344                      if zymin is not None:
    345                         mpl.pyplot.ylim(ymin=zymin)
    346                      if zymax is not None:
    347                         mpl.pyplot.ylim(ymax=zymax)
    348 
    349                      if invert_y:
    350                         lima,limb = mpl.pyplot.ylim()
    351                         mpl.pyplot.ylim(limb,lima)
    352                      if ylog:
    353                         mpl.pyplot.semilogy()
     349                     if zxmin is not None: mpl.pyplot.xlim(xmin=zxmin)
     350                     if zxmax is not None: mpl.pyplot.xlim(xmax=zxmax)
     351                     if zymin is not None: mpl.pyplot.ylim(ymin=zymin)
     352                     if zymax is not None: mpl.pyplot.ylim(ymax=zymax)
     353                     if invert_y:     lima,limb = mpl.pyplot.ylim() ; mpl.pyplot.ylim(limb,lima)
     354                     if ylog:         mpl.pyplot.semilogy()
    354355
    355356                 else:
    356357                     if hole:  what_I_plot = nolow(what_I_plot)
    357358                     pcolor( x, y, what_I_plot, cmap = palette, vmin=zevmin, vmax=zevmax )
    358                  if colorb != 'nobar' and var != 'HGT' :       
     359                 if colorb != 'nobar' and varname != 'HGT' :       
    359360                     if (fileref is not None) and (index_f is numplot-1):
    360361                        colorbar(fraction=0.05,pad=0.03,format="%.2f",\
    361                                            ticks=np.linspace(zevmin,zevmax,min([ndiv+1,20])),\
     362                                           ticks=np.linspace(zevmin,zevmax,num=min([ticks/2+1,20])),\
    362363                                           extend='neither',spacing='proportional')
    363364                     else:
    364365                        colorbar(fraction=0.05,pad=0.03,format=fmtvar(fvar.upper()),\
    365                                            ticks=np.linspace(zevmin,zevmax,min([ndiv+1,10])),\
     366                                           ticks=np.linspace(zevmin,zevmax,num=min([ticks/2+1,20])),\
    366367                                           extend='neither',spacing='proportional')
    367368
    … …  
    391392                   key = False
    392393               if metwind:  [vecx,vecy] = m.rotate_vector(vecx, vecy, lon2d, lat2d)
    393                if var == False:       colorvec = definecolorvec(back)
     394               if varname == False:       colorvec = definecolorvec(back)
    394395               else:                  colorvec = definecolorvec(colorb)
    395396               vectorfield(vecx, vecy,\
    … …  
    400401   
    401402       ### Next subplot
     403       basename = getname(var=varname,winds=winds,anomaly=anomaly)
    402404       if typefile in ['mesoapi','meso']:
    403405            plottitle = basename
    … …  
    411413                   plottitle = basename+' '+fileref
    412414                else:
    413                    plottitle = basename+' '+namefiles[index_f]
     415                   plottitle = basename+' '+namefiles[0]#index_f]
    414416            else:
    415                 plottitle = basename+' '+namefiles[index_f]
     417                plottitle = basename+' '+namefiles[0]#index_f]
    416418       if mult != 1:           plottitle = str(mult) + "*" + plottitle
    417419       if zetitle != "fill":   
    … …  
    432434#       if indextime is not None:
    433435#         plottitle = plottitle + " time: " + str(min(time[indextime])) +" "+ str(max(time[indextime]))
    434        ptitle( plottitle )
     436       title( plottitle )
    435437       itime += itstep
    436438       if nplot >= numplot:
    … …  
    464466    if found_lct:     
    465467        pad_inches_value = 0.35
    466         print "save", save
     468        print "********** SAVE ", save
    467469        if save == 'png':
    468470            if display: makeplotres(zeplot,res=100.,pad_inches_value=pad_inches_value) #,erase=True)  ## a miniature
    … …  
    473475            show()
    474476        else:
    475             print "save mode not supported. using gui instead."
     477            print "INFO: save mode not supported. using gui instead."
    476478            show()
    477     else:   print "Local time not found"
     479    else:   print "!!! Local time not found"
    478480
    479481    ###############

  • trunk/UTIL/PYTHON/pp.py

    r391 r392  
    11#!/usr/bin/env python
    22
    3 ### A. Spiga
     3### A. Spiga + T. Navarro + A. Colaitis
    44
    55###########################################################################################
    … …  
    1010    from optparse import OptionParser    ### to be replaced by argparse
    1111    from api_wrapper import api_onelevel
     12    from zrecast_wrapper import call_zrecast
    1213    from netCDF4 import Dataset
    13     from myplot import getlschar, separatenames, readslices, adjust_length
     14    from myplot import getlschar, separatenames, readslices, adjust_length, whatkindfile, errormess
    1415    from os import system
    1516    from planetoplot import planetoplot
    … …  
    1920    #############################
    2021    ### Get options and variables
    21     parser = OptionParser()
    22     getparseroptions(parser)
    23     (opt,args) = parser.parse_args()
     22    parser = OptionParser() ; getparseroptions(parser) ; (opt,args) = parser.parse_args()
     23    if opt.file is None:                                errormess("I want to eat one file at least ! Use winds.py -f name_of_my_file. Or type winds.py -h")
     24    if opt.var is None and opt.anomaly is True:         errormess("Cannot ask to compute anomaly if no variable is set")
     25    if opt.fref is not None and opt.operat is None:     errormess("you must specify an operation when using a reference file")
     26    if opt.operat in ["+","-"] and opt.fref is None:    errormess("you must specifiy a reference file when using inter-file operations")
     27    if opt.fref is not None and opt.operat is not None and opt.itp is not None: interpref=True
     28    else:   interpref=False
    2429
    25     if opt.file is None:
    26         print "I want to eat one file at least ! Use winds.py -f name_of_my_file. Or type winds.py -h"
    27         exit()
    28     if opt.var is None and opt.anomaly is True:
    29         print "Cannot ask to compute anomaly if no variable is set"
    30         exit()
    31     print "Options:", opt
     30    #############################
     31    ### Get infos about slices
     32    zeslat  = readslices(opt.slat) ; zeslon  = readslices(opt.slon) ; zesvert = readslices(opt.svert) ; zestime = readslices(opt.stime)
     33    reffile = opt.fref
     34    zexaxis = [opt.xmin,opt.xmax] ; zeyaxis=[opt.ymin,opt.ymax]
    3235
    33     listvar = ''
    34     if opt.var is None:
    35         zerange = [-999999]
    36     else:
    37         zelen = len(opt.var)
    38         zerange = range(zelen)
    39         #if zelen == 1: listvar = opt.var[0] + ','
    40         #else         :
    41         for jjj in zerange: listvar += opt.var[jjj] + ','
    42         listvar = listvar[0:len(listvar)-1]
    43         vmintab = adjust_length (opt.vmin, zelen) 
    44         vmaxtab = adjust_length (opt.vmax, zelen)
    45 
    46 
    47     ################################
    48     ### General check
    49  
    50     if opt.fref is not None:
    51        if opt.operat is None:
    52           print "you must specify an operation when using a reference file"
    53           exit()
    54     if opt.operat in ["+","-"]:
    55        if opt.fref is None:
    56           print "you must specifiy a reference file when using inter-file operations"
    57           exit()
    58 
    59     interpref=False
    60     if opt.fref is not None:
    61        if opt.operat is not None:
    62           if opt.itp is not None:
    63              interpref=True
    64 
    65     ################################
    66 
    67 
    68     print "file, length", opt.file, len(opt.file)
    69 
    70     zeslat  = readslices(opt.slat)
    71     zeslon  = readslices(opt.slon)
    72     zesvert = readslices(opt.svert)
    73     zestime = readslices(opt.stime)
    74     print "slat,zeslat", opt.slat, zeslat
    75     print "slon,zeslon", opt.slon, zeslon
    76     print "svert,zesvert", opt.svert, zesvert
    77     print "stime,zestime", opt.stime, zestime
    78 
     36    #############################
     37    ### 1. LOOP ON FILE LISTS TO BE PUT IN DIFFERENT FIGURES
    7938    for i in range(len(opt.file)):
    8039
    81         zefiles = separatenames(opt.file[i])
    82         print "zefiles", zefiles
     40      zefiles = separatenames(opt.file[i])
    8341
    84         #zefile = opt.file[i]
    85         #print zefile   
    86         zefile = zefiles[0]
    87              
    88         #zelevel = opt.lvl   
    89         stralt = None
    90         [lschar,zehour,zehourin] = getlschar ( zefile )  ## getlschar from wrfout (or simply return "" if another file)
    91    
     42      typefile = whatkindfile(Dataset(zefiles[0]))
     43      stralt = None
     44      if typefile in ["meso","mesoapi"]:         
     45          [lschar,zehour,zehourin] = getlschar ( zefiles[0] )
     46          if opt.var is None:  opt.var = ["HGT"]
     47      else:                                       
     48          [lschar,zehour,zehourin] = ["",0,0]
     49          if opt.var is None:  opt.var = ["phisinit"]
     50
     51      if opt.vmin is not None : zevmin  = opt.vmin[min(i,len(opt.vmin)-1)]
     52      else:                     zevmin = None
     53      if opt.vmax is not None : zevmax  = opt.vmax[min(i,len(opt.vmax)-1)]
     54      else:                     zevmax = None
     55      #print "vmin, zevmin", opt.vmin, zevmin ; print "vmax, zevmax", opt.vmax, zevmax
     56
     57      #############################
     58      ### 2. LOOP ON VAR LISTS TO BE PUT IN DIFFERENT FIGURES
     59      for j in range(len(opt.var)):
     60
     61        zevars = separatenames(opt.var[j])
     62
    9263        inputnvert = separatenames(opt.lvl)
    9364        if np.array(inputnvert).size == 1:
    … …  
    9768            zelevel = -99.
    9869            ze_interp_levels = np.linspace(float(inputnvert[0]),float(inputnvert[1]),float(inputnvert[2]))
    99         print 'level: ', zelevel
    100         print 'interp_levels: ',ze_interp_levels
    10170
    102         #####################################################
    103         ### Call Fortran routines for vertical interpolations       
     71        #########################################################
     72        ### Call Fortran routines for vertical interpolations ###     
     73        #########################################################
    10474        if opt.itp is not None:
     75          #####
     76          ##### MESOSCALE : written by AS
     77          #####
     78          if typefile == "meso":
    10579            if zelevel == 0. and opt.itp == 4:  zelevel = 0.010
    10680            ### winds or no winds
    … …  
    10882            else                    :  zefields = ''
    10983            ### var or no var
    110             #if opt.var is None      :  pass
    111             if zefields == ''       :  zefields = listvar
    112             else                    :  zefields = zefields + "," + listvar
    113             if opt.var2 is not None : zefields = zefields + "," + opt.var2 
    114             print zefields
    115             zefile = api_onelevel (  path_to_input   = '', \
    116                                      input_name      = zefile, \
    117                                      fields          = zefields, \
    118                                      interp_method   = opt.itp, \
    119                                      interp_level    = ze_interp_levels, \
    120                                      onelevel        = zelevel, \
    121                                      nocall          = opt.nocall )
    122             print zefile
     84            if zefields == ''       :  zefields = opt.var[j]
     85            else                    :  zefields = zefields + "," + opt.var[j]
     86            if opt.var2 is not None :  zefields = zefields + "," + opt.var2 
     87            ### call fortran routines
     88            for fff in range(len(zefiles)):
     89                newname = api_onelevel (  path_to_input   = '', \
     90                                               input_name      = zefiles[fff], \
     91                                               fields          = zefields, \
     92                                               interp_method   = opt.itp, \
     93                                               interp_level    = ze_interp_levels, \
     94                                               onelevel        = zelevel, \
     95                                               nocall          = opt.nocall )
     96                if fff == 0: zetab = newname
     97                else:        zetab = np.append(zetab,newname)
     98            zefiles = zetab #; print zefiles
    12399            zelevel = 0 ## so that zelevel could play again the role of nvert
     100          #####
     101          ##### GCM : written by AC
     102          #####
     103          elif typefile == "gcm":
     104            interpolated_files=""
     105            interpolated_files=call_zrecast(interp_mode=opt.itp,\
     106                    input_name=zefiles,\
     107                    fields=zevars)
    124108
    125         if opt.var is None: zerange = [-999999]
    126         else:               zerange = range(zelen)
    127 
    128 
    129 
    130 
    131 
    132 
    133 
    134 
    135         if interpref:
    136            interpolated_ref=""
    137            interpolated_ref=call_zrecast(interp_mode=opt.itp,\
     109            zefiles=interpolated_files
     110            if interpref:
     111               interpolated_ref=""
     112               interpolated_ref=call_zrecast(interp_mode=opt.itp,\
    138113                    input_name=[opt.fref],\
    139114                    fields=zevars)
    140115
    141            reffile=interpolated_ref[0]
    142         else:
    143            reffile=opt.fref
    144 # Divers ####################################################
    145    
    146         zexaxis=[opt.xmin,opt.xmax]
    147         zeyaxis=[opt.ymin,opt.ymax]
     116               reffile=interpolated_ref[0]
     117          else:
     118            print "not supported"
     119            exit()
    148120
    149         for jjj in zerange:
    150             if jjj == -999999:
    151                 argvar  = None
    152                 zevmin = None
    153                 zevmax = None
    154             else:
    155                 argvar = opt.var[jjj]
    156                 if vmintab[jjj] != -999999:  zevmin = vmintab[jjj]
    157                 else:                        zevmin = None
    158                 if vmaxtab[jjj] != -999999:  zevmax = vmaxtab[jjj]
    159                 else:                        zevmax = None
    160 
    161             #############
    162             ### Main call
    163             name = planetoplot (zefiles,level=int(zelevel),vertmode=opt.itp,\
    164                 proj=opt.proj,back=opt.back,target=opt.tgt,stride=opt.ste,var=argvar,\
     121        #############
     122        ### Main call
     123        name = planetoplot (zefiles,level=int(zelevel),vertmode=opt.itp,\
     124                proj=opt.proj,back=opt.back,target=opt.tgt,stride=opt.ste,var=zevars,\
    165125                numplot=opt.num,colorb=opt.clb,winds=opt.winds,\
    166126                addchar=lschar,interv=[zehour,zehourin],vmin=zevmin,vmax=zevmax,\
    … …  
    173133                ope=opt.operat,fileref=reffile,minop=opt.vminope,maxop=opt.vmaxope,titleref=opt.titref,\
    174134                invert_y=opt.inverty,xaxis=zexaxis,yaxis=zeyaxis,ylog=opt.logy,yintegral=opt.column,\
    175                 blat=opt.blat,tsat=opt.tstat)
    176         print 'Done: '+name
     135                blat=opt.blat,tsat=opt.tsat)
     136        print 'DONE: '+name
    177137        system("rm -f to_be_erased")
    178138 
    179139    #########################################################
    180140    ### Generate a .sh file with the used command saved in it
    181     command = ""
     141    command = "" 
    182142    for arg in sys.argv: command = command + arg + ' '
     143    if typefile not in ["meso","mesoapi"]: name = 'pycommand'
    183144    f = open(name+'.sh', 'w')
    184145    f.write(command)
     146
     147    print "********** OPTIONS: ", opt
     148    print "********************************************************** END"
Note: See TracChangeset for help on using the changeset viewer.