Changeset 475 for trunk/UTIL/PYTHON/planetoplot.py
- Timestamp:
- Dec 15, 2011, 4:05:49 PM (13 years ago)
- File:
-
- 1 edited
-
trunk/UTIL/PYTHON/planetoplot.py (modified) (11 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/UTIL/PYTHON/planetoplot.py
r468 r475 8 8 ### A. Spiga -- LMD -- 11~12/2011 -- Extended multivar subplot capabilities + cosmetic changes + general cleaning and tests 9 9 ### A. Colaitis -- LMD -- 12/2011 -- Added movie capability [mencoder must be installed] 10 ### A. Spiga -- LMD -- 12/2011 -- Added HTML animated page capability + general tests of consistency [winds, etc...] 10 ### A. Spiga -- LMD -- 12/2011 -- Added HTML animated page capability + general tests of consistency [winds, etc...] + consistent generic movie loop 11 11 12 12 def planetoplot (namefiles,\ … … 154 154 elif typefile in ['meso','mesoapi','geo','mesoideal']: 155 155 ### the following lines are kind of dirty... not possible to ask for several lats and lons 156 if vlon is not None or vlat is not None: indices = bidimfind(lon2d,lat2d,vlon,vlat) 157 if slon is not None: slon[0][0] = indices[0] ; slon[0][1] = indices[0] 158 if slat is not None: slat[0][0] = indices[1] ; slat[0][1] = indices[1] 156 if vlon is not None or vlat is not None: 157 indices = bidimfind(lon2d,lat2d,vlon,vlat) #,file=nc) #placer un point sur carte 158 lonp,latp = ( lon2d[indices[0],indices[1]] , lat2d[indices[0],indices[1]] ) 159 if typefile in ['mesoapi','meso'] and mapmode == 1: lon2d = dumpbdy(lon2d,6) ; lat2d = dumpbdy(lat2d,6) ### important to do that now and not before 160 if slon is not None: slon[0][0] = indices[1] ; slon[0][1] = indices[1] #...this is idx 161 if slat is not None: slat[0][0] = indices[0] ; slat[0][1] = indices[0] #...this is idy 159 162 if varname in ['PHTOT','W']: vertdim='BOTTOM-TOP_PATCH_END_STAG' 160 163 else: vertdim='BOTTOM-TOP_PATCH_END_UNSTAG' … … 236 239 tabtime = all_time[0];tab = all_var[0];k = 1 237 240 if var2: tab2 = all_var2[0] 238 while k != len(namefiles) :241 while k != len(namefiles) and len(all_time[k]) != 0: 239 242 if var2: tab2 = np.append(tab2,all_var2[k],axis=0) 240 243 tabtime = np.append(tabtime,all_time[k]) ; tab = np.append(tab,all_var[k],axis=0) ; k += 1 … … 292 295 print "********** INDEX LON:",indexlon," LAT:",indexlat," VERT:",indexvert," TIME:",indextime 293 296 #################################################################### 297 ########## REDUCE FIELDS 298 #################################################################### 294 299 error = False 295 300 varname = all_varname[index_f] … … 308 313 changesubplot = (numplot > 1) and (len(what_I_plot.shape) != 1) ## default for 1D plots: superimposed. to be reworked for better flexibility. 309 314 if changesubplot: subplot(subv,subh,nplot) 310 ### Map projection311 if mapmode == 1: m = define_proj(proj,wlon,wlat,back=back,blat=blat,blon=blon) ; x, y = m(lon2d, lat2d)312 elif mapmode ==0: m = None ; x = None ; y = None313 315 #################################################################### 314 315 if not error: 316 if error: 317 errormess("There is an error in reducing field !") 318 else: 316 319 ticks = ndiv + 1 317 320 fvar = varname 318 321 if anomaly: fvar = 'anomaly' 319 if mapmode == 0: 322 ### 323 if mapmode == 0: ### could this be moved inside imov loop ? 320 324 itime=indextime 321 325 if len(what_I_plot.shape) is 3:itime=[0] 326 m = None ; x = None ; y = None 322 327 what_I_plot, x, y = define_axis(lon,lat,vert,time,indexlon,indexlat,indexvert,\ 323 328 itime,what_I_plot, len(all_var[index_f].shape),vertmode) 324 zxmin, zxmax = xaxis ; zymin, zymax = yaxis 325 if zxmin is not None: mpl.pyplot.xlim(xmin=zxmin) 326 if zxmax is not None: mpl.pyplot.xlim(xmax=zxmax) 327 if zymin is not None: mpl.pyplot.ylim(ymin=zymin) 328 if zymax is not None: mpl.pyplot.ylim(ymax=zymax) 329 if invert_y: lima,limb = mpl.pyplot.ylim() ; mpl.pyplot.ylim(limb,lima) 330 if ylog: mpl.pyplot.semilogy() 331 332 if (fileref is not None) and (index_f is numplot-1): zevmin, zevmax = calculate_bounds(what_I_plot,vmin=minop,vmax=maxop) 329 ### 330 if (fileref is not None) and (index_f == numplot-1): zevmin, zevmax = calculate_bounds(what_I_plot,vmin=minop,vmax=maxop) 333 331 else: zevmin, zevmax = calculate_bounds(what_I_plot,vmin=vmin,vmax=vmax) 332 if (fileref is not None) and (index_f == numplot-1): colorb = "RdBu_r" 334 333 if colorb in ["def","nobar"]: palette = get_cmap(name=defcolorb(fvar.upper())) 335 elif (fileref is not None) and (index_f is numplot-1): palette = get_cmap(name="RdBu_r")336 334 else: palette = get_cmap(name=colorb) 337 338 ##### simple 2D field and movies of 2D fields 339 if len(what_I_plot.shape) >= 2: 340 if (len(what_I_plot.shape) is 3 and mrate is None): errormess("3D field. Use --rate RATE for movie or specify --time TIME. Exit.") 335 ##### 1. ELIMINATE >3D CASES 336 if len(what_I_plot.shape) >= 4: 337 print "WARNING!!! ",len(what_I_plot.shape),"-D PLOT NOT SUPPORTED !!! dimensions: ",what_I_plot.shape 338 errormess("Are you sure you did not forget to prescribe a dimension ?") 339 ##### 2. HANDLE simple 1D/2D field and movies of 1D/2D fields 340 else: 341 341 if mrate is not None: iend=len(time)-1 342 342 else: iend=0 343 343 imov = 0 344 if var2: which = "contour" ## have to start with contours rather than shading 345 else: which = "regular" 344 if len(what_I_plot.shape) == 3: 345 if var2: which = "contour" ## have to start with contours rather than shading 346 else: which = "regular" 347 if mrate is None: errormess("3D field. Use --rate RATE for movie or specify --time TIME. Exit.") 348 elif len(what_I_plot.shape) == 2: 349 if var2: which = "contour" ## have to start with contours rather than shading 350 else: which = "regular" 351 if mrate is not None: which = "unidim" 352 elif len(what_I_plot.shape) == 1: 353 which = "unidim" 354 ##### IMOV LOOP #### IMOV LOOP 346 355 while imov <= iend: 347 356 print "-> frame ",imov+1, which … … 355 364 if mrate is None or what_I_plot_contour.ndim < 3: what_I_plot_frame = what_I_plot_contour 356 365 else: what_I_plot_frame = what_I_plot_contour[imov,:,:] 357 if mrate is not None: 358 if mapmode == 1: 366 elif which == "unidim": 367 if mrate is None: what_I_plot_frame = what_I_plot 368 else: what_I_plot_frame = what_I_plot[:,imov] ## because swapaxes 369 #if mrate is not None: 370 if mapmode == 1: 359 371 m = define_proj(proj,wlon,wlat,back=back,blat=blat,blon=blon) ## this is dirty, defined above but out of imov loop 360 372 x, y = m(lon2d, lat2d) ## this is dirty, defined above but out of imov loop … … 362 374 # if typefile in ['mesoideal']: what_I_plot_frame = dumpbdy(what_I_plot_frame,0,stag='W',condition=dumped_vert_stag) 363 375 364 if imov >= 0: 365 # Renew axis directives for movie frames which are not the first one. 366 zxmin, zxmax = xaxis ; zymin, zymax = yaxis 367 if zxmin is not None: mpl.pyplot.xlim(xmin=zxmin) 368 if zxmax is not None: mpl.pyplot.xlim(xmax=zxmax) 369 if zymin is not None: mpl.pyplot.ylim(ymin=zymin) 370 if zymax is not None: mpl.pyplot.ylim(ymax=zymax) 371 if invert_y: lima,limb = mpl.pyplot.ylim() ; mpl.pyplot.ylim(limb,lima) 372 if ylog: mpl.pyplot.semilogy() 373 374 if which == "regular": 376 if which == "unidim": 377 lbl = "" 378 if indexlat is not None: lbl = lbl + " ix" + str(indexlat[0]) 379 if indexlon is not None: lbl = lbl + " iy" + str(indexlon[0]) 380 if indexvert is not None: lbl = lbl + " iz" + str(indexvert[0]) 381 if indextime is not None: lbl = lbl + " it" + str(indextime[0]) 382 if mrate is not None: x = y ## because swapaxes... 383 if indexvert is not None or indextime is None: plot(x,what_I_plot_frame,label=lbl) ## regular plot 384 else: plot(what_I_plot_frame,x,label=lbl) ## vertical profile 385 if nplot > 1: legend(loc='best') 386 if indextime is None and axtime is not None: xlabel(axtime.upper()) ## define the right label 387 if save == 'txt': writeascii(np.transpose(what_I_plot),'profile'+str(nplot*1000+imov)+'.txt') 388 389 elif which == "regular": 375 390 if hole: what_I_plot_frame = hole_bounds(what_I_plot_frame,zevmin,zevmax) 376 391 else: what_I_plot_frame = bounds(what_I_plot_frame,zevmin,zevmax) … … 384 399 if mapmode == 1: m.pcolor( x, y, what_I_plot_frame, cmap = palette, vmin=zevmin, vmax=zevmax, alpha=trans) 385 400 elif mapmode == 0: pcolor( x, y, what_I_plot_frame, cmap = palette, vmin=zevmin, vmax=zevmax, alpha=trans) 401 386 402 if colorb != 'nobar': 387 if (fileref is not None) and (index_f is numplot-1): daformat = "%.3f"403 if (fileref is not None) and (index_f == numplot-1): daformat = "%.3f" 388 404 elif mult != 1: daformat = "%.1f" 389 405 else: daformat = fmtvar(fvar.upper()) … … 413 429 elif mapmode == 1: cs = m.contour( x,y,what_I_plot_frame, zelevels, colors='k', linewidths = 1 ) #0.33 colors='w' )# , alpha=0.5) 414 430 415 416 if which == "regular": 431 if which in ["regular","unidim"]: 432 433 if nplot > 1 and which == "unidim": 434 pass ## because we superimpose nplot instances 435 else: 436 # Axis directives for movie frames [including the first one). 437 zxmin, zxmax = xaxis ; zymin, zymax = yaxis 438 if zxmin is not None: mpl.pyplot.xlim(xmin=zxmin) 439 if zxmax is not None: mpl.pyplot.xlim(xmax=zxmax) 440 if zymin is not None: mpl.pyplot.ylim(ymin=zymin) 441 if zymax is not None: mpl.pyplot.ylim(ymax=zymax) 442 if ylog: mpl.pyplot.semilogy() 443 if invert_y: ax = mpl.pyplot.gca() ; ax.set_ylim(ax.get_ylim()[::-1]) 444 417 445 if mrate is not None: 418 446 ### THIS IS A MENCODER MOVIE … … 437 465 myfile.write("last_image = "+str(iend)+";"+ '\n') 438 466 myfile.close() 439 if var2 : which = "contour"467 if var2 and which == "regular": which = "contour" 440 468 imov = imov+1 441 469 elif which == "contour": 442 470 which = "regular" 443 444 ##### 1D field445 elif len(what_I_plot.shape) is 1:446 lbl = ""447 if indexlat is not None: lbl = lbl + " ix" + str(indexlat)448 if indexlon is not None: lbl = lbl + " iy" + str(indexlon)449 if indexvert is not None: lbl = lbl + " iz" + str(indexvert)450 if indextime is not None: lbl = lbl + " it" + str(indextime)451 452 if indexvert is not None or indextime is None: plot(x,what_I_plot,label=lbl) ## regular plot453 else: plot(what_I_plot,x,label=lbl) ## vertical profile454 if nplot > 1: legend(loc='best')455 if indextime is None and axtime is not None: xlabel(axtime.upper()) ## define the right label456 if save == 'txt': writeascii(np.transpose(what_I_plot),'profile'+str(nplot)+'.txt')457 458 #### Other cases: (maybe plot 3-D field one day ??)459 else:460 print "WARNING!!! ",len(what_I_plot.shape),"-D PLOT NOT SUPPORTED !!! dimensions: ",what_I_plot.shape461 errormess("Are you sure you did not forget to prescribe a dimension ?")462 else:463 errormess("There is an error in reducing field !")464 471 465 472 ### Next subplot … … 469 476 if mrate is not None: basename = "movie_" + basename 470 477 if typefile in ['mesoapi','meso']: 471 if slon is not None: basename = basename + "_lon_" + str(int( lon2d[indices[1],indices[0]]))472 if slat is not None: basename = basename + "_lat_" + str(int( lat2d[indices[1],indices[0]]))478 if slon is not None: basename = basename + "_lon_" + str(int(round(lonp))) 479 if slat is not None: basename = basename + "_lat_" + str(int(round(latp))) 473 480 plottitle = basename 474 481 ### dans le nouveau systeme time=ls,sol,lt cette ligne pourrait ne servir a rien (ou deplacer au dessus)
Note: See TracChangeset
for help on using the changeset viewer.
