Index: /trunk/UTIL/PYTHON/myplot.py =================================================================== --- /trunk/UTIL/PYTHON/myplot.py (revision 568) +++ /trunk/UTIL/PYTHON/myplot.py (revision 569) @@ -35,4 +35,16 @@ ltst = ltst % 24 return ltst + +## Author: AC +def check_localtime(time): + a=-1 + for i in range(len(time)-1): + if time[i] > time[i+1]: a=i + if a >= 0: + print "Sorry, time axis is not regular." + print "Contourf needs regular axis... recasting" + for i in range(a+1): + time[i]=time[i]-24. + return time ## Author: AS, AC, JL @@ -71,5 +83,11 @@ mask = np.ma.getmask(masked) if (True in np.array(mask)):out=masked - else:out=field + else: + # missing values from api are 1e36 + masked=np.ma.masked_where(field > 1e35,field) + masked.set_fill_value([np.NaN]) + mask = np.ma.getmask(masked) + if (True in np.array(mask)):out=masked + else:out=field return out Index: /trunk/UTIL/PYTHON/myscript.py =================================================================== --- /trunk/UTIL/PYTHON/myscript.py (revision 568) +++ /trunk/UTIL/PYTHON/myscript.py (revision 569) @@ -38,4 +38,5 @@ parser.add_option('--xlabel', action='store',dest='xlab', type="string", default=None, help='customize the x-axis label') parser.add_option('--ylabel', action='store',dest='ylab', type="string", default=None, help='customize the y-axis label') + parser.add_option('--labels', action='store',dest='labels', type="string", default=None, help='customize 1D curve labels. String separated by commas. [None]') ### SPECIFIC FOR MAPPING [MAPMODE 1] Index: /trunk/UTIL/PYTHON/planetoplot.py =================================================================== --- /trunk/UTIL/PYTHON/planetoplot.py (revision 568) +++ /trunk/UTIL/PYTHON/planetoplot.py (revision 569) @@ -62,5 +62,6 @@ seevar=False,\ xlab=None,\ - ylab=None): + ylab=None,\ + lbls=None): #################################################################################################################### @@ -73,5 +74,5 @@ fmtvar,definecolorvec,defcolorb,getprefix,putpoints,calculate_bounds,errormess,definesubplot,\ zoomset,getcoorddef,getwinddef,whatkindfile,reducefield,bounds,getstralt,getfield,smooth,nolow,\ - getname,localtime,polarinterv,getsindex,define_axis,determineplot,readslices,bidimfind,getlschar,hole_bounds + getname,localtime,check_localtime,polarinterv,getsindex,define_axis,determineplot,readslices,bidimfind,getlschar,hole_bounds from mymath import deg,max,min,mean,get_tsat,writeascii,fig2data,fig2img import matplotlib as mpl @@ -139,5 +140,5 @@ ### ... VAR varname=var[vvv] - if varname not in nc.variables: + if ((varname not in nc.variables) and not ((typefile in ['meso']) and (varname in ['UV','uv','uvmet']))): if len(varinfile) == 1: varname = varinfile[0] else: varname = False @@ -196,4 +197,7 @@ ###### principle: calculate correct indices then repopulate slon and slat if slon is not None or slat is not None: + if slon is not None: sslon = np.zeros([1,2]) + if slat is not None: sslat = np.zeros([1,2]) + if firstfile and save == 'png' and typefile == 'meso' and "HGT" in varinfile: iwantawhereplot = nc #show a topo map with a cross on the chosen point else: iwantawhereplot = None #do not show anything, just select indices @@ -208,9 +212,8 @@ indices[iii,jjj,:] = bidimfind(lon2d,lat2d,vlon,vlat,file=iwantawhereplot) lonp,latp = ( lon2d[indices[iii,jjj,0],indices[iii,jjj,1]] , lat2d[indices[iii,jjj,0],indices[iii,jjj,1]] ) - #print vlon, lonp, vlat, latp for iii in range(numlon): for jjj in range(numlat): - if slon is not None: slon[iii][0] = indices[iii,0,1] ; slon[iii][1] = indices[iii,0,1] #...this is idx - if slat is not None: slat[jjj][0] = indices[0,jjj,0] ; slat[jjj][1] = indices[0,jjj,0] #...this is idy + if slon is not None: sslon[iii][0] = indices[iii,0,1] ; sslon[iii][1] = indices[iii,0,1] #...this is idx + if slat is not None: sslat[jjj][0] = indices[0,jjj,0] ; sslat[jjj][1] = indices[0,jjj,0] #...this is idy lonp,latp = ( lon2d[indices[0,0,0],indices[0,0,1]] , lat2d[indices[0,0,0],indices[0,0,1]] ) ###### @@ -235,9 +238,16 @@ elif "Time" in nc.variables: time = count + np.arange(0,len(nc.variables["Time"]),1) else: time = count + np.arange(0,1,1) - if nnn > 0: count = time[-1] + 1 ## so that a cat is possible with simple subscripts - else: count = 0 + if ope in ["cat"]: + if nnn > 0: count = time[-1] + 1 ## so that a cat is possible with simple subscripts + else: count = 0 + else: count = 0 if axtime in ["lt"]: - for i in range(len(time)): time[i] = localtime ( interv[0]+time[i]*interv[1], 0.5*(wlon[0]+wlon[1]) ) + ftime = np.zeros(len(time)) + for i in range(len(time)): + ftime[i] = localtime ( interv[0]+time[i]*interv[1], 0.5*(wlon[0]+wlon[1]) ) + time=ftime + time=check_localtime(time) print "LOCAL TIMES.... ", time + ### if typefile in ['geo']: vert = [0.] ; stime = readslices(str(0)) @@ -271,6 +281,25 @@ else: ##### GENERAL STUFF HERE - all_var[k] = getfield(nc,varname) - + if ((varname in ['UV','uv','uvmet']) and (typefile in ['meso']) and (varname not in nc.variables)): + if "U" in varinfile: zu=getfield(nc,'U') + elif "Um" in varinfile: zu=getfield(nc,'Um') + if "V" in varinfile: zv=getfield(nc,'V') + elif "Vm" in varinfile: zv=getfield(nc,'Vm') + znt,znz,zny,znx = np.array(zu).shape + if "U" in varinfile:znx=znx-1 + zuint = np.zeros([znt,znz,zny,znx]) + zvint = np.zeros([znt,znz,zny,znx]) + if "U" in varinfile: + for xx in np.arange(znx): + zuint[:,:,:,xx] = (zu[:,:,:,xx] + zu[:,:,:,xx+1])/2. + for yy in np.arange(zny): + zvint[:,:,yy,:] = (zv[:,:,yy,:] + zv[:,:,yy+1,:])/2. + else: + zuint=zu + zvint=zv + all_var[k] = np.sqrt(zuint**2 + zvint**2) + else: + all_var[k] = getfield(nc,varname) + print "**** GOT SUBDATA:",k," NAMEFILE:",namefile," VAR:",varname, var2 ; k += 1 ; firstfile = False #### End of for namefile in namefiles @@ -283,6 +312,30 @@ if len(var) > 1: errormess("for this operation... please set only one var !") if ope in ["-","+","-%"]: - if fileref is not None: all_var[k] = getfield(Dataset(fileref),all_varname[k-1]) ; all_varname[k] = all_varname[k-1] ; all_time[k] = all_time[k-1] ; all_namefile[k] = all_namefile[k-1] - else: errormess("fileref is missing!") + if fileref is not None: + + if ((all_varname[k-1] in ['UV','uv','uvmet']) and (typefile in ['meso']) and (all_varname[k-1] not in Dataset(fileref).variables)): + if "U" in varinfile: zu=getfield(Dataset(fileref),'U') + elif "Um" in varinfile: zu=getfield(Dataset(fileref),'Um') + if "V" in varinfile: zv=getfield(Dataset(fileref),'V') + elif "Vm" in varinfile: zv=getfield(Dataset(fileref),'Vm') + znt,znz,zny,znx = np.array(zu).shape + if "U" in varinfile:znx=znx-1 + zuint = np.zeros([znt,znz,zny,znx]) + zvint = np.zeros([znt,znz,zny,znx]) + if "U" in varinfile: + for xx in np.arange(znx): + zuint[:,:,:,xx] = (zu[:,:,:,xx] + zu[:,:,:,xx+1])/2. + for yy in np.arange(zny): + zvint[:,:,yy,:] = (zv[:,:,yy,:] + zv[:,:,yy+1,:])/2. + else: + zuint=zu + zvint=zv + all_var[k] = np.sqrt(zuint**2 + zvint**2) + else: + all_var[k] = getfield(Dataset(fileref),all_varname[k-1]) + + all_varname[k] = all_varname[k-1] ; all_time[k] = all_time[k-1] ; all_namefile[k] = all_namefile[k-1] ; all_var2[k] = all_var2[k-1] + + else: errormess("fileref is missing!") if ope == "-": all_var[k+1]= all_var[k-1] - all_var[k] elif ope == "+": all_var[k+1]= all_var[k-1] + all_var[k] @@ -291,5 +344,5 @@ masked.set_fill_value([np.NaN]) all_var[k+1]= 100.*(all_var[k-1] - masked)/masked - all_varname[k+1] = all_varname[k] ; all_time[k+1] = all_time[k] ; all_namefile[k+1] = all_namefile[k] ; numplot = numplot+2 + all_varname[k+1] = all_varname[k] ; all_time[k+1] = all_time[k] ; all_namefile[k+1] = all_namefile[k] ; all_var2[k+1] = all_var2[k] ; numplot = numplot+2 elif ope in ["cat"]: tabtime = all_time[0];tab = all_var[0];k = 1 @@ -329,6 +382,6 @@ ## get all indexes to be taken into account for this subplot and then reduce field ## 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 - indexlon = getsindex(slon,(nplot-1)%nlon,lon) - indexlat = getsindex(slat,((nplot-1)//nlon)%nlat,lat) + indexlon = getsindex(sslon,(nplot-1)%nlon,lon) + indexlat = getsindex(sslat,((nplot-1)//nlon)%nlat,lat) indexvert = getsindex(svert,((nplot-1)//(nlon*nlat))%nvert,vert) if ope is not None: @@ -456,4 +509,5 @@ if which == "unidim": + #lbls = [" tau=0.5"," tau=0.5"," tau=1"," tau=1"] lbl = "" if indexlat is not None: lbl = lbl + " ix" + str(indexlat[0]) @@ -462,4 +516,7 @@ if indextime is not None: lbl = lbl + " it" + str(indextime[0]) if lbl == "": lbl = namefiles[index_f] + + if lbls is not None: lbl=lbls[index_f] + if mrate is not None: x = y ## because swapaxes... #what_I_plot_frame = np.diff(what_I_plot_frame, n=1) ; x = x[1:] @@ -473,5 +530,5 @@ if save == 'txt': writeascii(np.transpose(what_I_plot),'profile'+str(nplot*1000+imov)+'.txt') - elif which == "regular": + elif which == "regular": if hole: what_I_plot_frame = hole_bounds(what_I_plot_frame,zevmin,zevmax) else: what_I_plot_frame = bounds(what_I_plot_frame,zevmin,zevmax) @@ -570,6 +627,6 @@ if mrate is not None: basename = "movie_" + basename if typefile in ['meso']: - if slon is not None: basename = basename + "_lon_" + str(int(round(lonp))) - if slat is not None: basename = basename + "_lat_" + str(int(round(latp))) + if sslon is not None: basename = basename + "_lon_" + str(int(round(lonp))) + if sslat is not None: basename = basename + "_lat_" + str(int(round(latp))) plottitle = basename ### dans le nouveau systeme time=ls,sol,lt cette ligne pourrait ne servir a rien (ou deplacer au dessus) Index: /trunk/UTIL/PYTHON/pp.py =================================================================== --- /trunk/UTIL/PYTHON/pp.py (revision 568) +++ /trunk/UTIL/PYTHON/pp.py (revision 569) @@ -86,6 +86,4 @@ ze_interp_levels = np.linspace(float(inputnvert[0]),float(inputnvert[1]),float(inputnvert[2])) - ######################################################### - ### Call Fortran routines for vertical interpolations ### ######################################################### if opt.itp is not None: @@ -114,4 +112,12 @@ if fff == 0: zetab = newname else: zetab = np.append(zetab,newname) + if interpref: + reffile = api_onelevel ( path_to_input = '', \ + input_name = opt.fref, \ + fields = zefields, \ + interp_method = opt.itp, \ + interp_level = ze_interp_levels, \ + onelevel = zelevel, \ + nocall = opt.nocall ) zefiles = zetab #; print zefiles zelevel = 0 ## so that zelevel could play again the role of nvert @@ -158,5 +164,5 @@ blat=opt.blat,blon=opt.blon,tsat=opt.tsat,flagnolow=opt.nolow,\ mrate=opt.rate,mquality=opt.quality,trans=opt.trans,zarea=opt.area,axtime=opt.axtime,\ - redope=opt.redope,seevar=opt.seevar,xlab=opt.xlab,ylab=opt.ylab) + redope=opt.redope,seevar=opt.seevar,xlab=opt.xlab,ylab=opt.ylab,lbls=separatenames(opt.labels)) print 'DONE: '+name system("rm -f to_be_erased")