Changeset 202
- Timestamp:
- Jul 11, 2011, 12:02:13 AM (13 years ago)
- Location:
- trunk/MESOSCALE
- Files:
-
- 4 edited
-
LMD_MM_MARS/SRC/POSTPROC/api.F90 (modified) (3 diffs)
-
PLOT/PYTHON/mylib/myplot.py (modified) (7 diffs)
-
PLOT/PYTHON/scripts/domain.py (modified) (2 diffs)
-
PLOT/PYTHON/scripts/winds.py (modified) (13 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/MESOSCALE/LMD_MM_MARS/SRC/POSTPROC/api.F90
r186 r202 418 418 dvalj(j) = num_metgrid_levels 419 419 ELSE IF ( LINLOG .eq. 3 ) THEN 420 dnamej(j) = ' bottom_top' !'altitude'420 dnamej(j) = 'altitude' !'bottom_top' !'altitude' 421 421 dvalj(j) = num_metgrid_levels 422 422 ELSE 423 dnamej(j) = ' bottom_top' !'altitude_abg'423 dnamej(j) = 'altitude_abg' !'bottom_top' !'altitude_abg' 424 424 dvalj(j) = num_metgrid_levels 425 425 ENDIF … … 793 793 IF ( test_dim_name == 'bottom_top_stag' ) fix_meta_stag = .TRUE. 794 794 IF (LINLOG .le. 2) test_dim_name = 'pressure' 795 IF (LINLOG .eq. 3) test_dim_name = ' bottom_top' !'altitude'796 IF (LINLOG .eq. 4) test_dim_name = ' bottom_top' !'altitude_abg'795 IF (LINLOG .eq. 3) test_dim_name = 'altitude' !'bottom_top' !'altitude' 796 IF (LINLOG .eq. 4) test_dim_name = 'altitude_abg' !'bottom_top' !'altitude_abg' 797 797 interpolate = .TRUE. 798 798 ENDIF … … 970 970 IF ( ii == 2 ) test_dim_name = 'south_north' 971 971 IF (( ii == 3 ) .and. (LINLOG .le. 2)) test_dim_name = 'pressure' 972 IF (( ii == 3 ) .and. (LINLOG .eq. 3)) test_dim_name = ' bottom_top' !'altitude'973 IF (( ii == 3 ) .and. (LINLOG .eq. 4)) test_dim_name = ' bottom_top' !'altitude_abg'972 IF (( ii == 3 ) .and. (LINLOG .eq. 3)) test_dim_name = 'altitude' !'bottom_top' !'altitude' 973 IF (( ii == 3 ) .and. (LINLOG .eq. 4)) test_dim_name = 'altitude_abg' !'bottom_top' !'altitude_abg' 974 974 IF ( ii == 4 ) test_dim_name = 'Time' 975 975 DO jj = 1,j -
trunk/MESOSCALE/PLOT/PYTHON/mylib/myplot.py
r201 r202 62 62 return lschar, zehour, zehourin 63 63 64 def getprefix (nc): 65 prefix = 'LMD_MMM_' 66 prefix = prefix + 'd'+str(getattr(nc,'GRID_ID'))+'_' 67 prefix = prefix + str(int(getattr(nc,'DX')/1000.))+'km_' 68 return prefix 69 64 70 def api_onelevel ( path_to_input = None, \ 65 71 input_name = 'wrfout_d0?_????-??-??_??:00:00', \ … … 177 183 return improc 178 184 185 def getwinds (nc,charu='Um',charv='Vm'): 186 import numpy as np 187 u = nc.variables[charu] 188 v = nc.variables[charv] 189 if charu == 'U': u = u[:, :, :, 0:len(u[0,0,0,:])-1] 190 if charv == 'V': v = v[:, :, 0:len(v[0,0,:,0])-1, :] 191 ### ou alors prendre les coordonnees speciales 192 return u,v 193 179 194 def vectorfield (u, v, x, y, stride=3, scale=15., factor=250., color='black', csmooth=1, key=True): 180 195 ## scale regle la reference du vecteur … … 182 197 import matplotlib.pyplot as plt 183 198 import numpy as np 184 #posx = np.max(x) + np.std(x) / 3. ## pb pour les domaines globaux ...185 #posy = np.mean(y)186 #posx = np.min(x)187 #posy = np.max(x)188 #posx = np.max(x) - np.std(x) / 10.189 #posy = np.max(y) + np.std(y) / 10.190 199 posx = np.min(x) - np.std(x) / 10. 191 200 posy = np.min(y) - np.std(y) / 10. … … 199 208 angles='xy',color=color,\ 200 209 scale=factor,width=widthvec ) 201 if color=='white': kcolor='black' 202 elif color=='yellow': kcolor=color 203 else: kcolor=color 210 if color in ['white','yellow']: kcolor='black' 211 else: kcolor=color 204 212 if key: p = plt.quiverkey(q,posx,posy,scale,\ 205 213 str(int(scale)),coordinates='data',color=kcolor,labelpos='S',labelsep = 0.03) … … 233 241 else: blat = wlat[0] 234 242 h = 2000. 235 radius = 3397200 243 radius = 3397200. 236 244 if char == "cyl": m = Basemap(rsphere=radius,projection='cyl',\ 237 245 llcrnrlat=wlat[0],urcrnrlat=wlat[1],llcrnrlon=wlon[0],urcrnrlon=wlon[1]) … … 261 269 "def": "%.1e",\ 262 270 "PTOT": "%.0f",\ 271 "HGT": "%.1e",\ 263 272 "USTM": "%.2f",\ 264 273 } … … 266 275 whichvar = "def" 267 276 return fmtvar[whichvar] 277 278 def defcolorb (whichone="def"): 279 whichcolorb = { \ 280 "def": "spectral",\ 281 "HGT": "spectral",\ 282 "tk": "gist_heat",\ 283 "QH2O": "PuBu",\ 284 } 285 if whichone not in whichcolorb: 286 whichone = "def" 287 return whichcolorb[whichone] 288 289 def definecolorvec (whichone="def"): 290 whichcolor = { \ 291 "def": "black",\ 292 "vis": "yellow",\ 293 "vishires": "yellow",\ 294 "molabw": "yellow",\ 295 "mola": "black",\ 296 "gist_heat": "white",\ 297 "hot": "tk",\ 298 "gist_rainbow": "black",\ 299 "spectral": "black",\ 300 "gray": "red",\ 301 "PuBu": "black",\ 302 } 303 if whichone not in whichcolor: 304 whichone = "def" 305 return whichcolor[whichone] 268 306 269 307 def marsmap (whichone="vishires"): -
trunk/MESOSCALE/PLOT/PYTHON/scripts/domain.py
r184 r202 1 1 #!/usr/bin/env python 2 2 3 ### A. Spiga -- LMD -- 30/0 5/20113 ### A. Spiga -- LMD -- 30/06/2011 -- slight modif early 07/2011 4 4 5 def domain (namefile,proj="ortho",back=" molabw",target=None,var='HGT'):5 def domain (namefile,proj="ortho",back="vishires",target=None,var='HGT'): 6 6 from netCDF4 import Dataset 7 from myplot import getcoord2d,define_proj,makeplotpng,simplinterv 8 from matplotlib.pyplot import contourf 7 from myplot import getcoord2d,define_proj,makeplotpng,simplinterv,getprefix 8 from matplotlib.pyplot import contourf,rcParams 9 from numpy.core.defchararray import find 10 ### 9 11 nc = Dataset(namefile) 12 ### 10 13 [lon2d,lat2d] = getcoord2d(nc) 11 14 [wlon,wlat] = simplinterv(lon2d,lat2d) 15 ### 12 16 m = define_proj(proj,wlon,wlat,back=back) 13 17 x, y = m(lon2d, lat2d) 18 ### 14 19 contourf(x, y, nc.variables[var][0,:,:], 40) 15 if not target: zeplot = namefile+".domain" 16 else: zeplot = target+"/domain" 20 ### 21 zeplot = getprefix(nc) + "domain" 22 if not target: zeplot = namefile[0:find(namefile,'wrfout')] + zeplot 23 else: zeplot = target + "/" + zeplot 24 ### 17 25 makeplotpng(zeplot,pad_inches_value=0.35) 26 #rcParams['savefig.facecolor'] = 'black' 27 #makeplotpng(zeplot+"b",pad_inches_value=0.35) 18 28 19 29 if __name__ == "__main__": … … 22 32 from optparse import OptionParser 23 33 parser = OptionParser() 24 parser.add_option('-f', action='store', dest='namefile', type="string", default=None, help='name of WRF file [NEEDED]')25 parser.add_option('-p', action='store', dest='proj', type="string", default="ortho", help='projection')26 parser.add_option('-b', action='store', dest='back', type="string", default=" molabw", help='background')27 parser.add_option('-t', action='store', dest='target', type="string", default=None, help='destination folder')28 parser.add_option('-v', action='store', dest='var', type="string", default='HGT', help='variable contoured')34 parser.add_option('-f', action='store', dest='namefile', type="string", default=None, help='name of WRF file [NEEDED]') 35 parser.add_option('-p', action='store', dest='proj', type="string", default="ortho", help='projection') 36 parser.add_option('-b', action='store', dest='back', type="string", default="vishires", help='background') 37 parser.add_option('-t', action='store', dest='target', type="string", default=None, help='destination folder') 38 parser.add_option('-v', action='store', dest='var', type="string", default='HGT', help='variable contoured') 29 39 (opt,args) = parser.parse_args() 30 40 if opt.namefile is None: -
trunk/MESOSCALE/PLOT/PYTHON/scripts/winds.py
r201 r202 1 1 #!/usr/bin/env python 2 2 3 ### A. Spiga -- LMD -- 30/06/2011 to 04/07/20113 ### A. Spiga -- LMD -- 30/06/2011 to 10/07/2011 4 4 ### Thanks to A. Colaitis for the parser trick 5 5 … … 30 30 ### Load librairies and functions 31 31 from netCDF4 import Dataset 32 from myplot import getcoord2d,define_proj,makeplotpng,simplinterv,vectorfield,ptitle,latinterv,getproj,wrfinterv,dumpbdy,fmtvar 32 from myplot import getcoord2d,define_proj,makeplotpng,simplinterv,vectorfield,ptitle,latinterv,getproj,wrfinterv,dumpbdy,\ 33 fmtvar,definecolorvec,getwinds,defcolorb,getprefix 33 34 from mymath import deg,max,min,mean 34 from matplotlib.pyplot import contour f, subplot, figure, rcParams, savefig, colorbar, pcolor35 from matplotlib.pyplot import contour,contourf, subplot, figure, rcParams, savefig, colorbar, pcolor 35 36 from matplotlib.cm import get_cmap 36 37 import numpy as np 38 from numpy.core.defchararray import find 39 40 ### 41 #rcParams['text.usetex'] = True 42 #rcParams['cairo.format'] = 'svg' 37 43 38 44 ###################### … … 78 84 elif proj == "lcc": [wlon,wlat] = wrfinterv(lon2d,lat2d) 79 85 else: [wlon,wlat] = simplinterv(lon2d,lat2d) 80 81 ###################82 #### Get local time83 #if typefile in ['mesoapi','meso']:84 # ltst = interv[0]85 # #ltst = int(getattr(nc, 'GMT') + 0.5*(wlon[0]+wlon[1])/15.)86 #elif typefile in ['gcm']:87 # ltst = 088 86 89 87 ############################################################################## … … 118 116 else: zevmax = vmax 119 117 print "bounds ", zevmin, zevmax 118 ### some already defined colormaps 119 if colorb is True: colorb = defcolorb(var) 120 120 121 121 ########################### … … 128 128 ######################################### 129 129 ### Name for title and graphics save file 130 if winds: basename = " WINDS_"130 if winds: basename = "UV_" 131 131 else: basename = "" 132 if var: 133 basename = basename + var134 if typefile is 'meso':stralt = "_lvl" + str(nvert)132 if var: basename = basename + var 133 ### 134 if typefile is 'meso': stralt = "_lvl" + str(nvert) 135 135 elif typefile is 'mesoapi': 136 136 zelevel = int(nc.variables['vert'][nvert]) 137 if 'altitude_abg' in nc.dimensions: stralt = "_"+str(zelevel)+"m-ALS" 138 elif 'bottom_top' in nc.dimensions: stralt = "_"+str(zelevel)+"m-AMR" 137 if 'altitude' in nc.dimensions: stralt = "_"+str(zelevel)+"m-AMR" 138 elif 'altitude_abg' in nc.dimensions: stralt = "_"+str(zelevel)+"m-ALS" 139 elif 'bottom_top' in nc.dimensions: stralt = "_"+str(zelevel)+"m" 139 140 elif 'pressure' in nc.dimensions: stralt = "_"+str(zelevel)+"Pa" 140 141 else: stralt = "" 142 ### 141 143 basename = basename + stralt 142 144 … … 147 149 if numplot == 4: 148 150 sub = 221 149 #fig.subplots_adjust(wspace = 0.1, hspace = 0.3)150 151 fig.subplots_adjust(wspace = 0.3, hspace = 0.3) 151 152 rcParams['font.size'] = int( rcParams['font.size'] * 2. / 3. ) … … 156 157 elif numplot == 3: 157 158 sub = 131 158 fig.subplots_adjust(wspace = 0. 3)159 rcParams['font.size'] = int( rcParams['font.size'] * 2. / 3. )159 fig.subplots_adjust(wspace = 0.5) 160 rcParams['font.size'] = int( rcParams['font.size'] * 1. / 2. ) 160 161 elif numplot == 6: 161 162 sub = 231 162 163 fig.subplots_adjust(wspace = 0.4, hspace = 0.0) 163 rcParams['font.size'] = int( rcParams['font.size'] * 2. / 3. )164 rcParams['font.size'] = int( rcParams['font.size'] * 1. / 2. ) 164 165 elif numplot == 8: 165 166 sub = 331 #241 166 fig.subplots_adjust(wspace = 0. 1, hspace = 0.3)167 fig.subplots_adjust(wspace = 0.3, hspace = 0.3) 167 168 rcParams['font.size'] = int( rcParams['font.size'] * 1. / 2. ) 168 169 elif numplot == 9: 169 170 sub = 331 170 fig.subplots_adjust(wspace = 0. 1, hspace = 0.3)171 fig.subplots_adjust(wspace = 0.3, hspace = 0.3) 171 172 rcParams['font.size'] = int( rcParams['font.size'] * 1. / 2. ) 172 173 elif numplot == 1: … … 223 224 else: 224 225 pcolor( x, y, what_I_plot, cmap = get_cmap(name=colorb), vmin=zevmin, vmax=zevmax ) 225 if colorb: colorbar(fraction=0.05,pad=0.1,format=fmtvar(var)) 226 if var in ['HGT']: pass 227 elif colorb: colorbar(fraction=0.05,pad=0.1,format=fmtvar(var)) 226 228 227 229 ### Vector plot … … 234 236 key = False 235 237 if metwind: [vecx,vecy] = m.rotate_vector(vecx, vecy, lon2d, lat2d) 236 if var == None and back == "vishires": colorvec = 'w'237 else: colorvec = 'k'238 if var == None: colorvec = definecolorvec(back) 239 else: colorvec = definecolorvec(colorb) 238 240 vectorfield(vecx, vecy,\ 239 241 x, y, stride=stride, csmooth=stride,\ … … 242 244 243 245 ### Next subplot 244 plottitle = basename + "_LT"+str(ltst) 245 if addchar: plottitle = plottitle + addchar 246 plottitle = basename 247 if addchar: plottitle = plottitle + addchar + "_LT"+str(ltst) 248 else: plottitle = plottitle + "_LT"+str(ltst) 246 249 ptitle( plottitle ) 247 250 sub += 1 … … 249 252 ########################################################################## 250 253 ### Save the figure in a file in the data folder or an user-defined folder 251 if not target: zeplot = namefile+"_"+basename 252 else: zeplot = target+"/"+basename 253 if numplot <= 0: zeplot = zeplot + "_LT"+str(abs(numplot)) 254 if typefile in ['meso','mesoapi']: prefix = getprefix(nc) 255 elif typefile in ['gcm']: prefix = 'LMD_GCM_' 256 else: prefix = '' 257 ### 258 zeplot = prefix + basename 259 if addchar: zeplot = zeplot + addchar 260 if numplot <= 0: zeplot = zeplot + "_LT"+str(abs(numplot)) 261 ### 262 if not target: zeplot = namefile[0:find(namefile,'wrfout')] + zeplot 263 else: zeplot = target + "/" + zeplot 264 ### 254 265 if found_lct: makeplotpng(zeplot,pad_inches_value=0.35) 255 266 else: print "Local time not found" 267 268 ############### 269 ### Now the end 256 270 return zeplot 257 258 259 ####################################################260 ####################################################261 ### A simple program to get wind vectors' components262 def getwinds (nc,charu='Um',charv='Vm'):263 import numpy as np264 u = nc.variables[charu]265 v = nc.variables[charv]266 if charu == 'U': u = u[:, :, :, 0:len(u[0,0,0,:])-1]267 if charv == 'V': v = v[:, :, 0:len(v[0,0,:,0])-1, :]268 ### ou alors prendre les coordonnees speciales269 return u,v270 271 272 271 273 272 ########################################################################################### … … 293 292 parser.add_option('-n', action='store', dest='numplot', type="int", default=4, help='number of plots (def=1)(<0: 1 plot of LT -*numplot*)') 294 293 parser.add_option('-i', action='store', dest='interp', type="int", default=None, help='interpolation method (2: press, 3: z-amr, 4:z-als)') 295 parser.add_option('-c', action='store', dest='colorb', type="string", default=True, help='change colormap (and draw a colorbar)')294 parser.add_option('-c', action='store', dest='colorb', type="string", default=True, help='change colormap') 296 295 parser.add_option('-x', action='store_false', dest='winds', default=True, help='flag: no wind vectors') 297 296 parser.add_option('-m', action='store', dest='vmin', type="float", default=None, help='bounding minimum value for color plot') … … 334 333 proj=opt.proj,back=opt.back,target=opt.target,stride=opt.stride,var=opt.var,numplot=opt.numplot,colorb=opt.colorb,winds=opt.winds,\ 335 334 addchar=lschar,interv=[zehour,zehourin],vmin=opt.vmin,vmax=opt.vmax,tile=opt.tile) 335 print 'Done: '+name 336 336 337 337 #########################################################
Note: See TracChangeset
for help on using the changeset viewer.
