[192] | 1 | import numpy as np |
---|
| 2 | import os |
---|
| 3 | from netCDF4 import Dataset as NetCDFFile |
---|
| 4 | import drawing_tools as drw |
---|
| 5 | from optparse import OptionParser |
---|
| 6 | import sys |
---|
| 7 | from cStringIO import StringIO |
---|
| 8 | |
---|
| 9 | ## g.e. # drawing.py -f /media/data2/etudes/WRF_LMDZ/WL_HyMeX/IIphase/medic950116/wlmdza/wrfout/wrfout_d01_1995-01-13_00:00:00 -o create_movie -S 'draw_2D_shad#Time@WRFTimes@10@95@191@1#tas:East_West|-1,North_South|-1,Time|2:longitude:latitude:Summer:270.,300.:tas|at|t=0:pdf:None:None' -v T2 |
---|
| 10 | ## g.e. # drawing.py -f wrfout_d01_1980-03-01_00\:00\:00_Time_B0-E48-I1.nc -o draw_2D_shad -S 'tas:East_West|-1,North_South|-1,Time|2:longitude:latitude:Summer:270.,300.:tas|at|t=0:pdf:None:cyl,i' -v T2 |
---|
| 11 | ## g.e. # drawing.py -f ~/etudes/domains/MEDCORDEX/geo_em.d01.nc -o draw_2D_shad_cont -S 'landmask,height:Time|0:Time|0:XLONG_M:XLAT_M:rainbow:fixc,k:%.2f:0,1:0.,3000.,10:landmask & height:pdf:False:lcc,i' -v LANDMASK,HGT_M |
---|
| 12 | ## g.e. # drawing.py -f ~/etudes/domains/MEDCORDEX/geo_em.d01.nc -o draw_2D_shad_cont -S 'height,landmask:Time|0:Time|0:XLONG_M:XLAT_M:terrain:fixc,k:None:0.,3000.:0,1,10:MEDCORDEX height & landmask:pdf:False:lcc,i' -v HGT_M,LANDMASK |
---|
| 13 | ## g.e. # drawing.py -o draw_2D_shad_line -f 'mean_dtcon-pluc-pres_lat.nc,mean_dtcon-pluc-pres_lat.nc' -S 'dtcon,prc:bottom_top|-1,south_north|-1:latmean:presmean:seismic,k:-5.,5.:monthly|dtcon|&|prc:pdf:flip@y:None:True' -v 'dtconmean,prcmean' |
---|
| 14 | ## g.e. # drawing.py -f 'geo_em.d02.nc' -o draw_topo_geogrid -S '0.,3000.:None:FF_3dom d02:png:cyl,i' |
---|
| 15 | ## g.e. # drawing.py -o draw_topo_geogrid_boxes -S '0.,3000.:None:FF domains:pdf:lcc,i:d01,d02,d03' -f geo_em.d01.nc,geo_em.d02.nc,geo_em.d03.nc |
---|
| 16 | ## g.e. # drawing.py -o draw_trajectories -f 'WRF/control/trajectory.dat@-1@/home/lluis/etudes/domains/WL_HyMeX_HighRes_C/geo_em.d03.nc#XLONG_M#XLAT_M,WRF_LMDZ/wlmdza/trajectory.dat@-1@/home/lluis/etudes/domains/WL_HyMeX_C/geo_em.d01.nc#XLONG_M#XLAT_M,WRF_LMDZ/wlmdzb/trajectory.dat@-1@/home/lluis/etudes/domains/WL_HyMeX_C/geo_em.d01.nc#XLONG_M#XLAT_M,WRF_LMDZ/wlmdzb_ii/trajectory.dat@-1@/home/lluis/etudes/domains/WL_HyMeX_C/geo_em.d01.nc#XLONG_M#XLAT_M' -S '$WRF_{CRM}$,$LMDZ_{AR4.0}$,$LMDZ_{NPv3.1}$,$LMDZ_{NPv3.1b}$|None|medicane trajectories|pdf|cyl,i' -v obs/trajectory.dat,satellite,-1 |
---|
| 17 | ## g.e. # drawing.py -o draw_vals_trajectories -f WRF_LMDZ/wlmdza/tevolboxtraj_T2.nc,WRF_LMDZ/wlmdzb/tevolboxtraj_T2.nc,WRF/control/tevolboxtraj_T2.nc -S 'mean:-1:$WRF_{CRM}$,$LMDZ_{AR4.0}$,$LMDZ_{NPv3.1}$@4:tas:time|($[DD]^[HH]$):exct,6,h:$%d^{%H}$:trajectory|following|mean:pdf' -v T2 |
---|
| 18 | ## g.e. # drawing.py -o draw_2D_shad_time -f 'netcdf_concatenated.nc' -S 'dtcon:Time|-1,bottom_top|-1:presmean:time:seismic:-3.e-6,3.e-6:monthly|dtcon:pdf:transpose:time|hours!since!1949-12-01|exct,2,d|%d|date!([DD])|x:True' -v 'dtconmean' |
---|
| 19 | ## g.e. # drawing.py -o variable_values -S PSFC |
---|
| 20 | ## g.e. # drawing.py -o draw_timeSeries -f wrfout_d01_1979-12-01_00:00:00_bottom_top_B6-E6-I1_south_north_B3-E3-I1_west_east_B26-E26-I1.nc -S 'dt_con:time|($[DD]^{[HH]}$):exct,12,h:$%d^{%H}$:time|evolution|at|-1|6|3|26:1:pdf' -v 'LDQCON,time' |
---|
| 21 | ## g.e. # drawing.py -f wrfout_d01_1979-12-01_00:00:00 -o draw_Neighbourghood_evol -S 'q:Time|-1|Times,bottom_top|6|ZNU,south_north|3|XLAT,west_east|26|XLONG:south_north,west_east:5:auto:time|($[DD]^{[HH]}$),time|($[DD]^{[HH]}$):exct,2,h|exct,1,d:$%d^{%H}$,$%d^{%H}$:5|pts|neighbourghood|temporal|evolution:0.0,0.004:BuPu:pdf:True' -v QVAPOR |
---|
| 22 | |
---|
| 23 | main = 'drawing.py' |
---|
| 24 | |
---|
| 25 | errormsg = 'ERROR -- error -- ERROR -- error' |
---|
| 26 | warnmsg = 'WARNING -- waring -- WARNING -- warning' |
---|
| 27 | fillValue=1.e20 |
---|
| 28 | |
---|
| 29 | namegraphics = ['create_movie', 'draw_2D_shad', 'draw_2D_shad_time', \ |
---|
| 30 | 'draw_2D_shad_cont', 'draw_2D_shad_cont_time', 'draw_2D_shad_line', \ |
---|
| 31 | 'draw_2D_shad_line_time', 'draw_timeSeries', 'draw_topo_geogrid', \ |
---|
| 32 | 'draw_topo_geogrid_boxes', 'draw_trajectories', 'draw_vals_trajectories', \ |
---|
| 33 | 'draw_lines', 'draw_Neighbourghood_evol', 'list_graphics', 'variable_values'] |
---|
| 34 | |
---|
| 35 | def draw_2D_shad(ncfile, values, varn): |
---|
| 36 | """ plotting a fields with shading |
---|
| 37 | draw_2D_shad(ncfile, values, varn) |
---|
| 38 | ncfile= file to use |
---|
| 39 | values=[vnamefs]:[dimvals]:[dimxvn]:[dimyvn]:[colorbar]:[sminv],[smaxv]:[figt]: |
---|
| 40 | [kindfig]:[reverse]:[mapv]:[close] |
---|
| 41 | [vnamefs]: Name in the figure of the variable to be shaded |
---|
| 42 | [dimvals]: ',' list of [dimname]|[value] telling at which dimension of the |
---|
| 43 | variable a given value is required (-1, all the length) |
---|
| 44 | [dimx/yvn]: name of the variables with the values of the final dimensions (x,y) |
---|
| 45 | [colorbar]: name of the color bar |
---|
| 46 | [smin/axv]: minimum and maximum value for the shading or: |
---|
| 47 | 'Srange': for full range |
---|
| 48 | 'Saroundmean@val': for mean-xtrm,mean+xtrm where xtrm = np.min(mean-min@val,max@val-mean) |
---|
| 49 | 'Saroundminmax@val': for min*val,max*val |
---|
| 50 | 'Saroundpercentile@val': for median-xtrm,median+xtrm where xtrm = np.min(median-percentile_(val), |
---|
| 51 | percentile_(100-val)-median) |
---|
| 52 | 'Smean@val': for -xtrm,xtrm where xtrm = np.min(mean-min*@val,max*@val-mean) |
---|
| 53 | 'Smedian@val': for -xtrm,xtrm where xtrm = np.min(median-min@val,max@val-median) |
---|
| 54 | 'Spercentile@val': for -xtrm,xtrm where xtrm = np.min(median-percentile_(val), |
---|
| 55 | percentile_(100-val)-median) |
---|
| 56 | [figt]: title of the figure ('|' for spaces) |
---|
| 57 | [kindfig]: kind of figure |
---|
| 58 | [reverse]: Transformation of the values |
---|
| 59 | * 'transpose': reverse the axes (x-->y, y-->x) |
---|
| 60 | * 'flip'@[x/y]: flip the axis x or y |
---|
| 61 | [mapv]: map characteristics: [proj],[res] |
---|
| 62 | see full documentation: http://matplotlib.org/basemap/ |
---|
| 63 | [proj]: projection |
---|
| 64 | * 'cyl', cilindric |
---|
| 65 | * 'lcc', lamvbert conformal |
---|
| 66 | [res]: resolution: |
---|
| 67 | * 'c', crude |
---|
| 68 | * 'l', low |
---|
| 69 | * 'i', intermediate |
---|
| 70 | * 'h', high |
---|
| 71 | * 'f', full |
---|
| 72 | valules= 'rh:z|-1,x|-1:z|-1,x|-1:lat:pressure:BuPu:0.,100.:rh:pdf:flip@y:None' |
---|
| 73 | varn= [varsn] name of the variable to plot with shading |
---|
| 74 | """ |
---|
| 75 | |
---|
| 76 | fname = 'draw_2D_shad' |
---|
| 77 | if values == 'h': |
---|
| 78 | print fname + '_____________________________________________________________' |
---|
| 79 | print draw_2D_shad.__doc__ |
---|
| 80 | quit() |
---|
| 81 | |
---|
| 82 | expectargs = ['[vnamefs]', '[dimvals]', '[dimxvn]', '[dimyvn]', '[colorbar]', \ |
---|
| 83 | '[sminv],[smaxv]', '[figt]', '[kindfig]', '[reverse]', '[mapv]', '[close]'] |
---|
| 84 | |
---|
| 85 | drw.check_arguments(fname,len(expectargs),values,':',expectargs) |
---|
| 86 | |
---|
| 87 | vnamesfig = values.split(':')[0] |
---|
| 88 | dimvals= values.split(':')[1].replace('|',':') |
---|
| 89 | vdimxn = values.split(':')[2] |
---|
| 90 | vdimyn = values.split(':')[3] |
---|
| 91 | colbarn = values.split(':')[4] |
---|
| 92 | shadminmax = values.split(':')[5] |
---|
| 93 | figtitle = values.split(':')[6].replace('|',' ') |
---|
| 94 | figkind = values.split(':')[7] |
---|
| 95 | revals = values.split(':')[8] |
---|
| 96 | mapvalue = values.split(':')[9] |
---|
| 97 | # varn = values.split(':')[10] |
---|
| 98 | |
---|
| 99 | ncfiles = ncfile |
---|
| 100 | |
---|
| 101 | if not os.path.isfile(ncfiles): |
---|
| 102 | print errormsg |
---|
| 103 | print ' ' + fname + ': shading file "' + ncfiles + '" does not exist !!' |
---|
| 104 | quit(-1) |
---|
| 105 | |
---|
| 106 | objsf = NetCDFFile(ncfiles, 'r') |
---|
| 107 | |
---|
| 108 | varns = varn.split(',')[0] |
---|
| 109 | |
---|
| 110 | if not objsf.variables.has_key(varns): |
---|
| 111 | print errormsg |
---|
| 112 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 113 | '" does not have variable "' + varns + '" !!' |
---|
| 114 | quit(-1) |
---|
| 115 | |
---|
| 116 | # Variables' values |
---|
| 117 | objvars = objsf.variables[varns] |
---|
| 118 | |
---|
| 119 | valshad, dimsshad = drw.slice_variable(objvars, dimvals.replace(',','|')) |
---|
| 120 | |
---|
| 121 | # Dimensions names |
---|
| 122 | ## print fname + ' obj dimnames: ', objvars.dimensions, dimvals, len(dimvals.split(',')) |
---|
| 123 | ## dimnamesv = [] |
---|
| 124 | ## for idd in range(len(objvars.dimensions)): |
---|
| 125 | ## cutdim = False |
---|
| 126 | ## for idc in range(len(dimvals.split(','))): |
---|
| 127 | ## dimcutn = dimvals.split(',')[idc].split(':')[0] |
---|
| 128 | ## print objvars.dimensions[idd], dimcutn |
---|
| 129 | ## if objvars.dimensions[idd] == dimcutn: |
---|
| 130 | ## cutdim = True |
---|
| 131 | ## break |
---|
| 132 | ## if not cutdim: dimnamesv.append(objvars.dimensions[idd]) |
---|
| 133 | dimnamesv = [vdimyn, vdimxn] |
---|
| 134 | |
---|
| 135 | if drw.searchInlist(objvars.ncattrs(),'units'): |
---|
| 136 | varunits = objvars.getncattr('units') |
---|
| 137 | else: |
---|
| 138 | print warnmsg |
---|
| 139 | print ' ' + fname + ": variable '" + varn + "' without units!!" |
---|
| 140 | varunits = '-' |
---|
| 141 | |
---|
| 142 | if not objsf.variables.has_key(vdimxn): |
---|
| 143 | print errormsg |
---|
| 144 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 145 | '" does not have dimension variable "' + vdimxn + '" !!' |
---|
| 146 | quit(-1) |
---|
| 147 | if not objsf.variables.has_key(vdimyn): |
---|
| 148 | print errormsg |
---|
| 149 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 150 | '" does not have dimension variable "' + vdimyn + '" !!' |
---|
| 151 | quit(-1) |
---|
| 152 | |
---|
| 153 | objdimx = objsf.variables[vdimxn] |
---|
| 154 | objdimy = objsf.variables[vdimyn] |
---|
| 155 | if drw.searchInlist(objdimx.ncattrs(),'units'): |
---|
| 156 | odimxu = objdimx.getncattr('units') |
---|
| 157 | else: |
---|
| 158 | print warnmsg |
---|
| 159 | print ' ' + fname + ": variable dimension '" + vdimxn + "' without units!!" |
---|
| 160 | odimxu = '-' |
---|
| 161 | |
---|
| 162 | if drw.searchInlist(objdimy.ncattrs(),'units'): |
---|
| 163 | odimyu = objdimy.getncattr('units') |
---|
| 164 | else: |
---|
| 165 | print warnmsg |
---|
| 166 | print ' ' + fname + ": variable dimension '" + vdimyn + "' without units!!" |
---|
| 167 | odimyu = '-' |
---|
| 168 | |
---|
| 169 | if len(objdimx.shape) <= 2: |
---|
| 170 | # odimxv = objdimx[valshad.shape] |
---|
| 171 | # odimyv = objdimy[valshad.shape] |
---|
| 172 | odimxv = objdimx[:] |
---|
| 173 | odimyv = objdimy[:] |
---|
| 174 | |
---|
| 175 | elif len(objdimx.shape) == 3: |
---|
| 176 | # dimcut = [0, slice(0,valshad.shape[0]), slice(0,valshad.shape[1])] |
---|
| 177 | # odimxv = objdimx[tuple(dimcut)] |
---|
| 178 | # odimyv = objdimy[tuple(dimcut)] |
---|
| 179 | odimxv = objdimx[0,:] |
---|
| 180 | odimyv = objdimy[0,:] |
---|
| 181 | else: |
---|
| 182 | print errormsg |
---|
| 183 | print ' ' + fname + ': shape of dimension variable:', objdimx.shape, \ |
---|
| 184 | ' not ready!!' |
---|
| 185 | quit(-1) |
---|
| 186 | |
---|
| 187 | shading_nx = [] |
---|
| 188 | if shadminmax.split(',')[0][0:1] != 'S': |
---|
| 189 | shading_nx.append(np.float(shadminmax.split(',')[0])) |
---|
| 190 | else: |
---|
| 191 | shading_nx.append(shadminmax.split(',')[0]) |
---|
| 192 | |
---|
| 193 | if shadminmax.split(',')[1][0:1] != 'S': |
---|
| 194 | shading_nx.append(np.float(shadminmax.split(',')[1])) |
---|
| 195 | else: |
---|
| 196 | shading_nx.append(shadminmax.split(',')[1]) |
---|
| 197 | |
---|
| 198 | if mapvalue == 'None': mapvalue = None |
---|
| 199 | |
---|
| 200 | drw.plot_2D_shadow(valshad, vnamesfig, odimxv, odimyv, odimxu, odimyu, dimnamesv,\ |
---|
| 201 | colbarn, shading_nx, varunits, figtitle, figkind, revals, mapvalue, True) |
---|
| 202 | |
---|
| 203 | return |
---|
| 204 | |
---|
| 205 | def draw_2D_shad_time(ncfile, values, varn): |
---|
| 206 | """ plotting a fields with shading with time values |
---|
| 207 | draw_2D_shad(ncfile, values, varn) |
---|
| 208 | ncfile= file to use |
---|
| 209 | values=[vnamefs]:[dimvals]:[dimxvn]:[dimyvn]:[colorbar]:[sminv],[smaxv]:[figt]: |
---|
| 210 | [kindfig]:[reverse]:[timevals]:[close] |
---|
| 211 | [vnamefs]: Name in the figure of the variable to be shaded |
---|
| 212 | [dimvals]: ',' list of [dimname]|[value] telling at which dimension of the |
---|
| 213 | variable a given value is required (-1, all the length) |
---|
| 214 | [dimx/yvn]: name of the variables with the values of the final dimensions (x,y) |
---|
| 215 | [colorbar]: name of the color bar |
---|
| 216 | [smin/axv]: minimum and maximum value for the shading or: |
---|
| 217 | 'Srange': for full range |
---|
| 218 | 'Saroundmean@val': for mean-xtrm,mean+xtrm where xtrm = np.min(mean-min@val,max@val-mean) |
---|
| 219 | 'Saroundminmax@val': for min*val,max*val |
---|
| 220 | 'Saroundpercentile@val': for median-xtrm,median+xtrm where xtrm = np.min(median-percentile_(val), |
---|
| 221 | percentile_(100-val)-median) |
---|
| 222 | 'Smean@val': for -xtrm,xtrm where xtrm = np.min(mean-min*@val,max*@val-mean) |
---|
| 223 | 'Smedian@val': for -xtrm,xtrm where xtrm = np.min(median-min@val,max@val-median) |
---|
| 224 | 'Spercentile@val': for -xtrm,xtrm where xtrm = np.min(median-percentile_(val), |
---|
| 225 | percentile_(100-val)-median) |
---|
| 226 | [figt]: title of the figure ('|' for spaces) |
---|
| 227 | [kindfig]: kind of figure |
---|
| 228 | [reverse]: Transformation of the values |
---|
| 229 | * 'transpose': reverse the axes (x-->y, y-->x) |
---|
| 230 | * 'flip'@[x/y]: flip the axis x or y |
---|
| 231 | [timevals]: [timen]|[units]|[kind]|[tfmt]|[label]|[timeaxis] time labels characteristics |
---|
| 232 | [timen]; name of the time variable |
---|
| 233 | [units]; units string according to CF conventions ([tunits] since |
---|
| 234 | [YYYY]-[MM]-[DD] [[HH]:[MI]:[SS]], '!' for spaces) |
---|
| 235 | [kind]; kind of output |
---|
| 236 | 'Nval': according to a given number of values as 'Nval',[Nval] |
---|
| 237 | 'exct': according to an exact time unit as 'exct',[tunit]; |
---|
| 238 | tunit= [Nunits],[tu]; [tu]= 'c': centuries, 'y': year, 'm': month, |
---|
| 239 | 'w': week, 'd': day, 'h': hour, 'i': minute, 's': second, |
---|
| 240 | 'l': milisecond |
---|
| 241 | [tfmt]; desired format |
---|
| 242 | [label]; label at the graph ('!' for spaces) |
---|
| 243 | [close]: should figure be closed (finished) |
---|
| 244 | values='dtcon:Time|-1,bottom_top|-1:presmean:time:seismic:-3.e-6,3.e-6:monthly|' |
---|
| 245 | 'dtcon:pdf:transpose:time|hours!since!1949-12-01|exct,2,d|%d|date!([DD])|x:True |
---|
| 246 | varn= [varsn] name of the variable to plot with shading |
---|
| 247 | """ |
---|
| 248 | fname = 'draw_2D_shad_time' |
---|
| 249 | |
---|
| 250 | if values == 'h': |
---|
| 251 | print fname + '_____________________________________________________________' |
---|
| 252 | print draw_2D_shad_time.__doc__ |
---|
| 253 | quit() |
---|
| 254 | |
---|
| 255 | farguments = ['[vnamefs]', '[dimvals]', '[dimxvn]', '[dimyvn]', '[colorbar]', \ |
---|
| 256 | '[sminv],[smaxv]', '[figt]', '[kindfig]', '[reverse]', '[timevals]', '[close]'] |
---|
| 257 | drw.check_arguments(fname,11,values,':',farguments) |
---|
| 258 | |
---|
| 259 | vnamesfig = values.split(':')[0] |
---|
| 260 | dimvals= values.split(':')[1].replace('|',':') |
---|
| 261 | vdimxn = values.split(':')[2] |
---|
| 262 | vdimyn = values.split(':')[3] |
---|
| 263 | colbarn = values.split(':')[4] |
---|
| 264 | shadminmax = values.split(':')[5] |
---|
| 265 | figtitle = values.split(':')[6].replace('|',' ') |
---|
| 266 | figkind = values.split(':')[7] |
---|
| 267 | revals = values.split(':')[8] |
---|
| 268 | timevals = values.split(':')[9] |
---|
| 269 | close = values.split(':')[10] |
---|
| 270 | |
---|
| 271 | ncfiles = ncfile |
---|
| 272 | |
---|
| 273 | if not os.path.isfile(ncfiles): |
---|
| 274 | print errormsg |
---|
| 275 | print ' ' + fname + ': shading file "' + ncfiles + '" does not exist !!' |
---|
| 276 | quit(-1) |
---|
| 277 | |
---|
| 278 | objsf = NetCDFFile(ncfiles, 'r') |
---|
| 279 | |
---|
| 280 | varns = varn.split(',')[0] |
---|
| 281 | |
---|
| 282 | if not objsf.variables.has_key(varns): |
---|
| 283 | print errormsg |
---|
| 284 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 285 | '" does not have variable "' + varns + '" !!' |
---|
| 286 | quit(-1) |
---|
| 287 | |
---|
| 288 | # Variables' values |
---|
| 289 | objvars = objsf.variables[varns] |
---|
| 290 | |
---|
| 291 | valshad, dimsshad = drw.slice_variable(objvars, dimvals.replace(',','|')) |
---|
| 292 | |
---|
| 293 | dimnamesv = [vdimyn, vdimxn] |
---|
| 294 | |
---|
| 295 | varunits = objvars.getncattr('units') |
---|
| 296 | |
---|
| 297 | if not objsf.variables.has_key(vdimxn): |
---|
| 298 | print errormsg |
---|
| 299 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 300 | '" does not have dimension variable "' + vdimxn + '" !!' |
---|
| 301 | quit(-1) |
---|
| 302 | if not objsf.variables.has_key(vdimyn): |
---|
| 303 | print errormsg |
---|
| 304 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 305 | '" does not have dimensino variable "' + vdimyn + '" !!' |
---|
| 306 | quit(-1) |
---|
| 307 | |
---|
| 308 | objdimx = objsf.variables[vdimxn] |
---|
| 309 | objdimy = objsf.variables[vdimyn] |
---|
| 310 | odimxu = objdimx.getncattr('units') |
---|
| 311 | odimyu = objdimy.getncattr('units') |
---|
| 312 | |
---|
| 313 | if len(objdimx.shape) <= 2: |
---|
| 314 | odimxv = objdimx[:] |
---|
| 315 | odimyv = objdimy[:] |
---|
| 316 | |
---|
| 317 | elif len(objdimx.shape) == 3: |
---|
| 318 | odimxv = objdimx[0,:] |
---|
| 319 | odimyv = objdimy[0,:] |
---|
| 320 | else: |
---|
| 321 | print errormsg |
---|
| 322 | print ' ' + fname + ': shape of dimension variable:', objdimx.shape, \ |
---|
| 323 | ' not ready!!' |
---|
| 324 | quit(-1) |
---|
| 325 | |
---|
| 326 | timename = timevals.split('|')[0] |
---|
| 327 | timeunit = timevals.split('|')[1].replace('!',' ') |
---|
| 328 | timekind = timevals.split('|')[2] |
---|
| 329 | timefmt = timevals.split('|')[3] |
---|
| 330 | timelabel = timevals.split('|')[4].replace('!',' ') |
---|
| 331 | |
---|
| 332 | if vdimxn == timename: |
---|
| 333 | odimxv = objsf.variables[vdimxn][:] |
---|
| 334 | odimxu = timelabel |
---|
| 335 | timeaxis = 'x' |
---|
| 336 | odimyv = objsf.variables[vdimyn] |
---|
| 337 | odimyu = odimyv.getncattr('units') |
---|
| 338 | timepos, timelabels = drw.CFtimes_plot(odimxv, timeunit, timekind, timefmt) |
---|
| 339 | elif vdimyn == timename: |
---|
| 340 | odimyv = objsf.variables[vdimyn][:] |
---|
| 341 | odimyu = timelabel |
---|
| 342 | timeaxis = 'y' |
---|
| 343 | odimxv = objsf.variables[vdimxn] |
---|
| 344 | odimxu = odimxv.getncattr('units') |
---|
| 345 | timepos, timelabels = drw.CFtimes_plot(odimyv, timeunit, timekind, timefmt) |
---|
| 346 | else: |
---|
| 347 | print errormsg |
---|
| 348 | print ' ' + fname + ": time variable '" + timename + "' not found!!" |
---|
| 349 | quit(-1) |
---|
| 350 | |
---|
| 351 | shading_nx = [] |
---|
| 352 | if shadminmax.split(',')[0][0:1] != 'S': |
---|
| 353 | shading_nx.append(np.float(shadminmax.split(',')[0])) |
---|
| 354 | else: |
---|
| 355 | shading_nx.append(shadminmax.split(',')[0]) |
---|
| 356 | |
---|
| 357 | if shadminmax.split(',')[1][0:1] != 'S': |
---|
| 358 | shading_nx.append(np.float(shadminmax.split(',')[1])) |
---|
| 359 | else: |
---|
| 360 | shading_nx.append(shadminmax.split(',')[1]) |
---|
| 361 | |
---|
| 362 | closeval = drw.Str_Bool(close) |
---|
| 363 | |
---|
| 364 | drw.plot_2D_shadow_time(valshad, vnamesfig, odimxv, odimyv, odimxu, odimyu, \ |
---|
| 365 | dimnamesv, colbarn, shading_nx, varunits, figtitle, figkind, revals, timeaxis, \ |
---|
| 366 | timepos, timelabels, closeval) |
---|
| 367 | |
---|
| 368 | return |
---|
| 369 | |
---|
| 370 | def draw_2D_shad_cont(ncfile, values, varn): |
---|
| 371 | """ plotting two fields, one with shading and the other with contour lines |
---|
| 372 | draw_2D_shad_cont(ncfile, values, varn) |
---|
| 373 | ncfile= [ncfilevars],[ncfilevarc] files to use (one value, same file) |
---|
| 374 | values=[vnamefs]:[dimvals]:[dimvalc]:[dimxvn]:[dimyvn]:[colorbar]:[ckind]:[clabfmt]:[sminv],[smaxv]:[sminc],[smaxv],[Nlev]:[figt]:[kindfig]:[reverse]:[mapv] |
---|
| 375 | [vnamefs],[vnamefc]: Name in the figure of the shaded and the contour variables |
---|
| 376 | [dimvals/c]: list of [dimname]|[value] telling at which dimension of the |
---|
| 377 | variable a given value is required (no dimension name, all the length) |
---|
| 378 | [dimx/yvn]: ',' list with the name of the variables with the values of the dimensions |
---|
| 379 | [colorbar]: name of the color bar |
---|
| 380 | [ckind]: kind of contours |
---|
| 381 | 'cmap': as it gets from colorbar |
---|
| 382 | 'fixc,[colname]': fixed color [colname], all stright lines |
---|
| 383 | 'fixsigc,[colname]': fixed color [colname], >0 stright, <0 dashed line |
---|
| 384 | [clabfmt]: format of the labels in the contour (None, also possible) |
---|
| 385 | [smin/axv]: minimum and maximum value for the shading |
---|
| 386 | [sminc]:[smaxv]:[Nlev]: minimum, maximum and number of values for the contour |
---|
| 387 | [figt]: title of the figure ('|' for spaces) |
---|
| 388 | [kindfig]: kind of figure |
---|
| 389 | [reverse]: does the values be transposed? 'True/False', |
---|
| 390 | [mapv]: map characteristics: [proj],[res] |
---|
| 391 | see full documentation: http://matplotlib.org/basemap/ |
---|
| 392 | [proj]: projection |
---|
| 393 | * 'cyl', cilindric |
---|
| 394 | * 'lcc', lamvbert conformal |
---|
| 395 | [res]: resolution: |
---|
| 396 | * 'c', crude |
---|
| 397 | * 'l', low |
---|
| 398 | * 'i', intermediate |
---|
| 399 | * 'h', high |
---|
| 400 | * 'f', full |
---|
| 401 | valules= 'rh,ta:z|-1,x|-1:z|-1,x|-1:lat:pressure:BuPu:fixsigc,black:%d:0.,100.:195.,305.,7:Meridonal|average|of|rh|&|ta:pdf:flip@y:None' |
---|
| 402 | varn= [varsn],[varcn] name of the variable to plot with shading variable with contour |
---|
| 403 | """ |
---|
| 404 | |
---|
| 405 | fname = 'draw_2D_shad_cont' |
---|
| 406 | if values == 'h': |
---|
| 407 | print fname + '_____________________________________________________________' |
---|
| 408 | print draw_2D_shad_cont.__doc__ |
---|
| 409 | quit() |
---|
| 410 | |
---|
| 411 | vnamesfig = values.split(':')[0].split(',') |
---|
| 412 | dimvals= values.split(':')[1].replace('|',':') |
---|
| 413 | dimvalc= values.split(':')[2].replace('|',':') |
---|
| 414 | vdimxn = values.split(':')[3] |
---|
| 415 | vdimyn = values.split(':')[4] |
---|
| 416 | colbarn = values.split(':')[5] |
---|
| 417 | countkind = values.split(':')[6] |
---|
| 418 | countlabelfmt = values.split(':')[7] |
---|
| 419 | shadminmax = values.split(':')[8] |
---|
| 420 | contlevels = values.split(':')[9] |
---|
| 421 | figtitle = values.split(':')[10].replace('|',' ') |
---|
| 422 | figkind = values.split(':')[11] |
---|
| 423 | revals = values.split(':')[12] |
---|
| 424 | mapvalue = values.split(':')[13] |
---|
| 425 | |
---|
| 426 | if2filenames = ncfile.find(',') |
---|
| 427 | |
---|
| 428 | if if2filenames != -1: |
---|
| 429 | ncfiles = ncfile.split(',')[0] |
---|
| 430 | ncfilec = ncfile.split(',')[1] |
---|
| 431 | else: |
---|
| 432 | ncfiles = ncfile |
---|
| 433 | ncfilec = ncfile |
---|
| 434 | |
---|
| 435 | if not os.path.isfile(ncfiles): |
---|
| 436 | print errormsg |
---|
| 437 | print ' ' + fname + ': shading file "' + ncfiles + '" does not exist !!' |
---|
| 438 | quit(-1) |
---|
| 439 | |
---|
| 440 | if not os.path.isfile(ncfilec): |
---|
| 441 | print errormsg |
---|
| 442 | print ' ' + fname + ': contour file "' + ncfilec + '" does not exist !!' |
---|
| 443 | quit(-1) |
---|
| 444 | |
---|
| 445 | objsf = NetCDFFile(ncfiles, 'r') |
---|
| 446 | objcf = NetCDFFile(ncfilec, 'r') |
---|
| 447 | |
---|
| 448 | varns = varn.split(',')[0] |
---|
| 449 | varnc = varn.split(',')[1] |
---|
| 450 | |
---|
| 451 | if not objsf.variables.has_key(varns): |
---|
| 452 | print errormsg |
---|
| 453 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 454 | '" does not have variable "' + varns + '" !!' |
---|
| 455 | quit(-1) |
---|
| 456 | |
---|
| 457 | if not objcf.variables.has_key(varnc): |
---|
| 458 | print errormsg |
---|
| 459 | print ' ' + fname + ': contour file "' + ncfilec + \ |
---|
| 460 | '" does not have variable "' + varnc + '" !!' |
---|
| 461 | quit(-1) |
---|
| 462 | |
---|
| 463 | # Variables' values |
---|
| 464 | objvars = objsf.variables[varns] |
---|
| 465 | objvarc = objcf.variables[varnc] |
---|
| 466 | |
---|
| 467 | if len(objvars.shape) != len(objvarc.shape): |
---|
| 468 | print errormsg |
---|
| 469 | print ' ' + fname + ': shading variable "' + varns + '" has a shape: ', \ |
---|
| 470 | objvars.shape, 'different than contour variable "' + varnc + '": ', \ |
---|
| 471 | objvarc.shape,' !!!' |
---|
| 472 | quit(-1) |
---|
| 473 | |
---|
| 474 | for idim in range(len(objvars.shape)): |
---|
| 475 | if objvars.shape[idim] != objvarc.shape[idim]: |
---|
| 476 | print errormsg |
---|
| 477 | print ' ' + fname + ': shading variable "' + varns + '" has a shape: ', \ |
---|
| 478 | objvars.shape, 'different than contour variable "' + varnc + '": ', \ |
---|
| 479 | objvarc.shape,' !!!' |
---|
| 480 | quit(-1) |
---|
| 481 | |
---|
| 482 | valshad, dimsshad = drw.slice_variable(objvars, dimvals.replace(',','|')) |
---|
| 483 | valcont, dimscont = drw.slice_variable(objvarc, dimvalc.replace(',','|')) |
---|
| 484 | |
---|
| 485 | # Dimensions names |
---|
| 486 | ## print fname + ' obj dimnames: ', objvars.dimensions, dimvals, len(dimvals.split(',')) |
---|
| 487 | ## dimnamesv = [] |
---|
| 488 | ## for idd in range(len(objvars.dimensions)): |
---|
| 489 | ## cutdim = False |
---|
| 490 | ## for idc in range(len(dimvals.split(','))): |
---|
| 491 | ## dimcutn = dimvals.split(',')[idc].split(':')[0] |
---|
| 492 | ## print objvars.dimensions[idd], dimcutn |
---|
| 493 | ## if objvars.dimensions[idd] == dimcutn: |
---|
| 494 | ## cutdim = True |
---|
| 495 | ## break |
---|
| 496 | ## if not cutdim: dimnamesv.append(objvars.dimensions[idd]) |
---|
| 497 | dimnamesv = [vdimyn, vdimxn] |
---|
| 498 | |
---|
| 499 | varunits = [] |
---|
| 500 | varunits.append(objvars.getncattr('units')) |
---|
| 501 | varunits.append(objvarc.getncattr('units')) |
---|
| 502 | |
---|
| 503 | if not objsf.variables.has_key(vdimxn): |
---|
| 504 | print errormsg |
---|
| 505 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 506 | '" does not have dimension variable "' + vdimxn + '" !!' |
---|
| 507 | quit(-1) |
---|
| 508 | if not objsf.variables.has_key(vdimyn): |
---|
| 509 | print errormsg |
---|
| 510 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 511 | '" does not have dimensino variable "' + vdimyn + '" !!' |
---|
| 512 | quit(-1) |
---|
| 513 | |
---|
| 514 | objdimx = objsf.variables[vdimxn] |
---|
| 515 | objdimy = objsf.variables[vdimyn] |
---|
| 516 | odimxu = objdimx.getncattr('units') |
---|
| 517 | odimyu = objdimy.getncattr('units') |
---|
| 518 | |
---|
| 519 | if len(objdimx.shape) <= 2: |
---|
| 520 | # odimxv = objdimx[valshad.shape] |
---|
| 521 | # odimyv = objdimy[valshad.shape] |
---|
| 522 | odimxv = objdimx[:] |
---|
| 523 | odimyv = objdimy[:] |
---|
| 524 | |
---|
| 525 | elif len(objdimx.shape) == 3: |
---|
| 526 | # dimcut = [0, slice(0,valshad.shape[0]), slice(0,valshad.shape[1])] |
---|
| 527 | # odimxv = objdimx[tuple(dimcut)] |
---|
| 528 | # odimyv = objdimy[tuple(dimcut)] |
---|
| 529 | odimxv = objdimx[0,:] |
---|
| 530 | odimyv = objdimy[0,:] |
---|
| 531 | else: |
---|
| 532 | print errormsg |
---|
| 533 | print ' ' + fname + ': shape of dimension variable:', objdimx.shape, \ |
---|
| 534 | ' not ready!!' |
---|
| 535 | quit(-1) |
---|
| 536 | |
---|
| 537 | if countlabelfmt == 'None': |
---|
| 538 | countlfmt = None |
---|
| 539 | else: |
---|
| 540 | countlfmt = countlabelfmt |
---|
| 541 | |
---|
| 542 | shading_nx = np.zeros((2), dtype=np.float) |
---|
| 543 | shading_nx[0] = np.float(shadminmax.split(',')[0]) |
---|
| 544 | shading_nx[1] = np.float(shadminmax.split(',')[1]) |
---|
| 545 | |
---|
| 546 | clevmin = np.float(contlevels.split(',')[0]) |
---|
| 547 | clevmax = np.float(contlevels.split(',')[1]) |
---|
| 548 | Nclevels = int(contlevels.split(',')[2]) |
---|
| 549 | |
---|
| 550 | levels_cont = drw.pretty_int(clevmin, clevmax, Nclevels) |
---|
| 551 | |
---|
| 552 | if len(levels_cont) <= 1: |
---|
| 553 | print warnmsg |
---|
| 554 | print ' ' + fname + ': wrong contour levels:', levels_cont,' !!' |
---|
| 555 | del(levels_cont) |
---|
| 556 | levels_cont = np.zeros((Nclevels), dtype=np.float) |
---|
| 557 | levels_cont = np.arange(7)*(clevmax - clevmin)/(Nclevels-1) |
---|
| 558 | print ' generating default ones: ',levels_cont |
---|
| 559 | |
---|
| 560 | if mapvalue == 'None': mapvalue = None |
---|
| 561 | |
---|
| 562 | drw.plot_2D_shadow_contour(valshad, valcont, vnamesfig, odimxv, odimyv, odimxu, \ |
---|
| 563 | odimyu, dimnamesv, colbarn, countkind, countlfmt, shading_nx, levels_cont, \ |
---|
| 564 | varunits, figtitle, figkind, revals, mapvalue) |
---|
| 565 | |
---|
| 566 | return |
---|
| 567 | |
---|
| 568 | def draw_2D_shad_cont_time(ncfile, values, varn): |
---|
| 569 | """ plotting two fields, one with shading and the other with contour lines being |
---|
| 570 | one of the dimensions of time characteristics |
---|
| 571 | draw_2D_shad_cont(ncfile, values, varn) |
---|
| 572 | ncfile= [ncfilevars],[ncfilevarc] files to use (one value, same file) |
---|
| 573 | values=[vnamefs]:[dimvals]:[dimvalc]:[dimxvn]:[dimyvn]:[colorbar]:[ckind]:[clabfmt]:[sminv],[smaxv]:[sminc],[smaxv],[Nlev]:[figt]:[kindfig]:[reverse]:[timevals]:[mapv] |
---|
| 574 | [vnamefs],[vnamefc]: Name in the figure of the shaded and the contour variables |
---|
| 575 | [dimvals/c]: list of [dimname]|[value] telling at which dimension of the |
---|
| 576 | variable a given value is required (no dimension name, all the length) |
---|
| 577 | [dimx/yvn]: ',' list with the name of the variables with the values of the dimensions |
---|
| 578 | [colorbar]: name of the color bar |
---|
| 579 | [ckind]: kind of contours |
---|
| 580 | 'cmap': as it gets from colorbar |
---|
| 581 | 'fixc,[colname]': fixed color [colname], all stright lines |
---|
| 582 | 'fixsigc,[colname]': fixed color [colname], >0 stright, <0 dashed line |
---|
| 583 | [clabfmt]: format of the labels in the contour (None, also possible) |
---|
| 584 | [smin/axv]: minimum and maximum value for the shading |
---|
| 585 | [sminc]:[smaxv]:[Nlev]: minimum, maximum and number of values for the contour |
---|
| 586 | [figt]: title of the figure ('|' for spaces) |
---|
| 587 | [kindfig]: kind of figure |
---|
| 588 | [reverse]: modification to the dimensions: |
---|
| 589 | 'transposed': transpose matrices |
---|
| 590 | 'flip',[x/y]: flip only the dimension [x] or [y] |
---|
| 591 | [timevals]: [timen]|[units]|[kind]|[tfmt]|[label] time labels characteristics |
---|
| 592 | [timen]; name of the time variable |
---|
| 593 | [units]; units string according to CF conventions ([tunits] since |
---|
| 594 | [YYYY]-[MM]-[DD] [[HH]:[MI]:[SS]], '!' for spaces) |
---|
| 595 | [kind]; kind of output |
---|
| 596 | 'Nval': according to a given number of values as 'Nval',[Nval] |
---|
| 597 | 'exct': according to an exact time unit as 'exct',[tunit]; |
---|
| 598 | tunit= [Nunits],[tu]; [tu]= 'c': centuries, 'y': year, 'm': month, |
---|
| 599 | 'w': week, 'd': day, 'h': hour, 'i': minute, 's': second, |
---|
| 600 | 'l': milisecond |
---|
| 601 | [tfmt]; desired format |
---|
| 602 | [label]; label at the graph ('!' for spaces) |
---|
| 603 | [mapv]: map characteristics: [proj],[res] |
---|
| 604 | see full documentation: http://matplotlib.org/basemap/ |
---|
| 605 | [proj]: projection |
---|
| 606 | * 'cyl', cilindric |
---|
| 607 | * 'lcc', lamvbert conformal |
---|
| 608 | [res]: resolution: |
---|
| 609 | * 'c', crude |
---|
| 610 | * 'l', low |
---|
| 611 | * 'i', intermediate |
---|
| 612 | * 'h', high |
---|
| 613 | * 'f', full |
---|
| 614 | valules= 'rh,ta:z|-1,x|-1:z|-1,x|-1:lat:pressure:BuPu:fixsigc,black:%d:0.,100.:195.,305.,7:Meridonal|average|of|rh|&|ta:pdf:flip@y:time!hours!since!1949/12/01|exct,5d|%d|date!([DD]):None' |
---|
| 615 | varn= [varsn],[varcn] name of the variable to plot with shading variable with contour |
---|
| 616 | """ |
---|
| 617 | |
---|
| 618 | fname = 'draw_2D_shad_cont_time' |
---|
| 619 | if values == 'h': |
---|
| 620 | print fname + '_____________________________________________________________' |
---|
| 621 | print draw_2D_shad_cont_time.__doc__ |
---|
| 622 | quit() |
---|
| 623 | |
---|
| 624 | expectargs = ['[vnamefs]', '[dimvals]', '[dimvalc]', '[dimxvn]', '[dimyvn]', \ |
---|
| 625 | '[colorbar]', '[ckind]', '[clabfmt]', '[sminv],[smaxv]', \ |
---|
| 626 | '[sminc],[smaxv],[Nlev]', '[figt]', '[kindfig]', '[reverse]', '[timevals]', \ |
---|
| 627 | '[mapv]'] |
---|
| 628 | |
---|
| 629 | drw.check_arguments(fname,len(expectargs),values,':',expectargs) |
---|
| 630 | |
---|
| 631 | vnamesfig = values.split(':')[0].split(',') |
---|
| 632 | dimvals= values.split(':')[1].replace('|',':') |
---|
| 633 | dimvalc= values.split(':')[2].replace('|',':') |
---|
| 634 | vdimxn = values.split(':')[3] |
---|
| 635 | vdimyn = values.split(':')[4] |
---|
| 636 | colbarn = values.split(':')[5] |
---|
| 637 | countkind = values.split(':')[6] |
---|
| 638 | countlabelfmt = values.split(':')[7] |
---|
| 639 | shadminmax = values.split(':')[8] |
---|
| 640 | contlevels = values.split(':')[9] |
---|
| 641 | figtitle = values.split(':')[10].replace('|',' ') |
---|
| 642 | figkind = values.split(':')[11] |
---|
| 643 | revals = values.split(':')[12] |
---|
| 644 | timevals = values.split(':')[13] |
---|
| 645 | mapvalue = values.split(':')[14] |
---|
| 646 | |
---|
| 647 | if2filenames = ncfile.find(',') |
---|
| 648 | |
---|
| 649 | if if2filenames != -1: |
---|
| 650 | ncfiles = ncfile.split(',')[0] |
---|
| 651 | ncfilec = ncfile.split(',')[1] |
---|
| 652 | else: |
---|
| 653 | ncfiles = ncfile |
---|
| 654 | ncfilec = ncfile |
---|
| 655 | |
---|
| 656 | if not os.path.isfile(ncfiles): |
---|
| 657 | print errormsg |
---|
| 658 | print ' ' + fname + ': shading file "' + ncfiles + '" does not exist !!' |
---|
| 659 | quit(-1) |
---|
| 660 | |
---|
| 661 | if not os.path.isfile(ncfilec): |
---|
| 662 | print errormsg |
---|
| 663 | print ' ' + fname + ': contour file "' + ncfilec + '" does not exist !!' |
---|
| 664 | quit(-1) |
---|
| 665 | |
---|
| 666 | objsf = NetCDFFile(ncfiles, 'r') |
---|
| 667 | objcf = NetCDFFile(ncfilec, 'r') |
---|
| 668 | |
---|
| 669 | varns = varn.split(',')[0] |
---|
| 670 | varnc = varn.split(',')[1] |
---|
| 671 | |
---|
| 672 | if not objsf.variables.has_key(varns): |
---|
| 673 | print errormsg |
---|
| 674 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 675 | '" does not have variable "' + varns + '" !!' |
---|
| 676 | quit(-1) |
---|
| 677 | |
---|
| 678 | if not objcf.variables.has_key(varnc): |
---|
| 679 | print errormsg |
---|
| 680 | print ' ' + fname + ': contour file "' + ncfilec + \ |
---|
| 681 | '" does not have variable "' + varnc + '" !!' |
---|
| 682 | quit(-1) |
---|
| 683 | |
---|
| 684 | # Variables' values |
---|
| 685 | objvars = objsf.variables[varns] |
---|
| 686 | objvarc = objcf.variables[varnc] |
---|
| 687 | |
---|
| 688 | if len(objvars.shape) != len(objvarc.shape): |
---|
| 689 | print errormsg |
---|
| 690 | print ' ' + fname + ': shading variable "' + varns + '" has a shape: ', \ |
---|
| 691 | objvars.shape, 'different than contour variable "' + varnc + '": ', \ |
---|
| 692 | objvarc.shape,' !!!' |
---|
| 693 | quit(-1) |
---|
| 694 | |
---|
| 695 | for idim in range(len(objvars.shape)): |
---|
| 696 | if objvars.shape[idim] != objvarc.shape[idim]: |
---|
| 697 | print errormsg |
---|
| 698 | print ' ' + fname + ': shading variable "' + varns + '" has a shape: ', \ |
---|
| 699 | objvars.shape, 'different than contour variable "' + varnc + '": ', \ |
---|
| 700 | objvarc.shape,' !!!' |
---|
| 701 | quit(-1) |
---|
| 702 | |
---|
| 703 | valshad, dimsshad = drw.slice_variable(objvars, dimvals.replace(',','|')) |
---|
| 704 | valcont, dimscont = drw.slice_variable(objvarc, dimvalc.replace(',','|')) |
---|
| 705 | |
---|
| 706 | # Dimensions names |
---|
| 707 | ## print fname + ' obj dimnames: ', objvars.dimensions, dimvals, len(dimvals.split(',')) |
---|
| 708 | ## dimnamesv = [] |
---|
| 709 | ## for idd in range(len(objvars.dimensions)): |
---|
| 710 | ## cutdim = False |
---|
| 711 | ## for idc in range(len(dimvals.split(','))): |
---|
| 712 | ## dimcutn = dimvals.split(',')[idc].split(':')[0] |
---|
| 713 | ## print objvars.dimensions[idd], dimcutn |
---|
| 714 | ## if objvars.dimensions[idd] == dimcutn: |
---|
| 715 | ## cutdim = True |
---|
| 716 | ## break |
---|
| 717 | ## if not cutdim: dimnamesv.append(objvars.dimensions[idd]) |
---|
| 718 | dimnamesv = [vdimyn, vdimxn] |
---|
| 719 | |
---|
| 720 | varunits = [] |
---|
| 721 | varunits.append(objvars.getncattr('units')) |
---|
| 722 | varunits.append(objvarc.getncattr('units')) |
---|
| 723 | |
---|
| 724 | if not objsf.variables.has_key(vdimxn): |
---|
| 725 | print errormsg |
---|
| 726 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 727 | '" does not have dimension variable "' + vdimxn + '" !!' |
---|
| 728 | quit(-1) |
---|
| 729 | if not objsf.variables.has_key(vdimyn): |
---|
| 730 | print errormsg |
---|
| 731 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 732 | '" does not have dimension variable "' + vdimyn + '" !!' |
---|
| 733 | quit(-1) |
---|
| 734 | |
---|
| 735 | timename = timevals.split('|')[0] |
---|
| 736 | timeunit = timevals.split('|')[1].replace('!',' ') |
---|
| 737 | timekind = timevals.split('|')[2] |
---|
| 738 | timefmt = timevals.split('|')[3] |
---|
| 739 | timelabel = timevals.split('|')[4].replace('!',' ') |
---|
| 740 | |
---|
| 741 | if vdimxn == timename: |
---|
| 742 | timevals = objsf.variables[vdimxn][:] |
---|
| 743 | dimt = 'x' |
---|
| 744 | ovalaxis = objsf.variables[vdimyn] |
---|
| 745 | ovalu = ovalaxis.getncattr('units') |
---|
| 746 | elif vdimyn == timename: |
---|
| 747 | timevals = objsf.variables[vdimyn][:] |
---|
| 748 | dimt = 'y' |
---|
| 749 | ovalaxis = objsf.variables[vdimxn] |
---|
| 750 | ovalu = ovalaxis.getncattr('units') |
---|
| 751 | else: |
---|
| 752 | print errormsg |
---|
| 753 | print ' ' + fname + ": time variable '" + timename + "' not found!!" |
---|
| 754 | quit(-1) |
---|
| 755 | |
---|
| 756 | timepos, timelabels = drw.CFtimes_plot(timevals, timeunit, timekind, timefmt) |
---|
| 757 | |
---|
| 758 | if len(ovalaxis.shape) <= 2: |
---|
| 759 | ovalaxisv = ovalaxis[:] |
---|
| 760 | |
---|
| 761 | elif len(ovalaxis.shape) == 3: |
---|
| 762 | ovalaxisv = ovalaxis[0,:] |
---|
| 763 | else: |
---|
| 764 | print errormsg |
---|
| 765 | print ' ' + fname + ': shape of dimension variable:', ovalaxis.shape, \ |
---|
| 766 | ' not ready!!' |
---|
| 767 | quit(-1) |
---|
| 768 | |
---|
| 769 | if countlabelfmt == 'None': |
---|
| 770 | countlfmt = None |
---|
| 771 | else: |
---|
| 772 | countlfmt = countlabelfmt |
---|
| 773 | |
---|
| 774 | shading_nx = np.zeros((2), dtype=np.float) |
---|
| 775 | shading_nx[0] = np.float(shadminmax.split(',')[0]) |
---|
| 776 | shading_nx[1] = np.float(shadminmax.split(',')[1]) |
---|
| 777 | |
---|
| 778 | clevmin = np.float(contlevels.split(',')[0]) |
---|
| 779 | clevmax = np.float(contlevels.split(',')[1]) |
---|
| 780 | Nclevels = int(contlevels.split(',')[2]) |
---|
| 781 | |
---|
| 782 | levels_cont = drw.pretty_int(clevmin, clevmax, Nclevels) |
---|
| 783 | |
---|
| 784 | if len(levels_cont) <= 1: |
---|
| 785 | print warnmsg |
---|
| 786 | print ' ' + fname + ': wrong contour levels:', levels_cont,' !!' |
---|
| 787 | del(levels_cont) |
---|
| 788 | levels_cont = np.zeros((Nclevels), dtype=np.float) |
---|
| 789 | levels_cont = np.arange(7)*(clevmax - clevmin)/(Nclevels-1) |
---|
| 790 | print ' generating default ones: ',levels_cont |
---|
| 791 | |
---|
| 792 | if mapvalue == 'None': mapvalue = None |
---|
| 793 | |
---|
| 794 | drw.plot_2D_shadow_contour_time(valshad, valcont, vnamesfig, ovalaxisv, \ |
---|
| 795 | timevals, timepos, timelabels, ovalu, timelabel, dimt, dimnamesv, colbarn, \ |
---|
| 796 | countkind, countlfmt, shading_nx, levels_cont, varunits, figtitle, figkind, \ |
---|
| 797 | revals, mapvalue) |
---|
| 798 | |
---|
| 799 | return |
---|
| 800 | |
---|
| 801 | def draw_2D_shad_line(ncfile, values, varn): |
---|
| 802 | """ plotting a fields with shading and another with line |
---|
| 803 | draw_2D_shad_line(ncfile, values, varn) |
---|
| 804 | ncfile= [ncfiles],[ncfilel] file to use for the shading and for the line |
---|
| 805 | values=[vnamefs],[vnamefl]:[dimvals]:[dimxvn]:[dimyvn]:[colorbar],[colline]:[sminv],[smaxv]:[figt]: |
---|
| 806 | [kindfig]:[reverse]:[mapv]:[close] |
---|
| 807 | [vnamefs]: Name in the figure of the variable to be shaded |
---|
| 808 | [vnamefl]: Name in the figure of the variable to be lined |
---|
| 809 | [dimvals]: ',' list of [dimname]|[value] telling at which dimension of the |
---|
| 810 | variable a given value is required (-1, all the length) |
---|
| 811 | [dimx/yvn]: name of the variables with the values of the final dimensions (x,y) |
---|
| 812 | [colorbar]: name of the color bar |
---|
| 813 | [colline]: name of the color for the line |
---|
| 814 | [smin/axv]: minimum and maximum value for the shading |
---|
| 815 | [figt]: title of the figure ('|' for spaces) |
---|
| 816 | [kindfig]: kind of figure |
---|
| 817 | [reverse]: Transformation of the values |
---|
| 818 | * 'transpose': reverse the axes (x-->y, y-->x) |
---|
| 819 | * 'flip'@[x/y]: flip the axis x or y |
---|
| 820 | [mapv]: map characteristics: [proj],[res] |
---|
| 821 | see full documentation: http://matplotlib.org/basemap/ |
---|
| 822 | [proj]: projection |
---|
| 823 | * 'cyl', cilindric |
---|
| 824 | * 'lcc', lamvbert conformal |
---|
| 825 | [res]: resolution: |
---|
| 826 | * 'c', crude |
---|
| 827 | * 'l', low |
---|
| 828 | * 'i', intermediate |
---|
| 829 | * 'h', high |
---|
| 830 | * 'f', full |
---|
| 831 | valules= 'rh:z|-1,x|-1:z|-1,x|-1:lat:pressure:BuPu:0.,100.:rh:pdf:flip@y:None' |
---|
| 832 | varn= [varsn],[varnl] name of the variable to plot with shading and with line |
---|
| 833 | """ |
---|
| 834 | |
---|
| 835 | fname = 'draw_2D_shad_line' |
---|
| 836 | if values == 'h': |
---|
| 837 | print fname + '_____________________________________________________________' |
---|
| 838 | print draw_2D_shad_line.__doc__ |
---|
| 839 | quit() |
---|
| 840 | |
---|
| 841 | farguments = ['[vnamefs],[vnamefl]', '[dimvals]', '[dimxvn]', '[dimyvn]', \ |
---|
| 842 | '[colorbar],[colline]', '[sminv],[smaxv]', '[figt]', '[kindfig]', '[reverse]', \ |
---|
| 843 | '[mapv]', '[close]'] |
---|
| 844 | drw.check_arguments(fname,11,values,':',farguments) |
---|
| 845 | |
---|
| 846 | vnamesfig = values.split(':')[0].split(',')[0] |
---|
| 847 | dimvals= values.split(':')[1].replace('|',':') |
---|
| 848 | vdimxn = values.split(':')[2] |
---|
| 849 | vdimyn = values.split(':')[3] |
---|
| 850 | colbarn = values.split(':')[4].split(',')[0] |
---|
| 851 | shadminmax = values.split(':')[5] |
---|
| 852 | figtitle = values.split(':')[6].replace('|',' ') |
---|
| 853 | figkind = values.split(':')[7] |
---|
| 854 | revals = values.split(':')[8] |
---|
| 855 | mapvalue = values.split(':')[9] |
---|
| 856 | # varn = values.split(':')[10] |
---|
| 857 | |
---|
| 858 | ncfiles = ncfile.split(',')[0] |
---|
| 859 | |
---|
| 860 | if not os.path.isfile(ncfiles): |
---|
| 861 | print errormsg |
---|
| 862 | print ' ' + fname + ': shading file "' + ncfiles + '" does not exist !!' |
---|
| 863 | quit(-1) |
---|
| 864 | |
---|
| 865 | objsf = NetCDFFile(ncfiles, 'r') |
---|
| 866 | |
---|
| 867 | varns = varn.split(',')[0] |
---|
| 868 | |
---|
| 869 | if not objsf.variables.has_key(varns): |
---|
| 870 | print errormsg |
---|
| 871 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 872 | '" does not have variable "' + varns + '" !!' |
---|
| 873 | quit(-1) |
---|
| 874 | |
---|
| 875 | # Variables' values |
---|
| 876 | objvars = objsf.variables[varns] |
---|
| 877 | |
---|
| 878 | valshad, dimsshad = drw.slice_variable(objvars, dimvals.replace(',','|')) |
---|
| 879 | |
---|
| 880 | # Dimensions names |
---|
| 881 | ## print fname + ' obj dimnames: ', objvars.dimensions, dimvals, len(dimvals.split(',')) |
---|
| 882 | ## dimnamesv = [] |
---|
| 883 | ## for idd in range(len(objvars.dimensions)): |
---|
| 884 | ## cutdim = False |
---|
| 885 | ## for idc in range(len(dimvals.split(','))): |
---|
| 886 | ## dimcutn = dimvals.split(',')[idc].split(':')[0] |
---|
| 887 | ## print objvars.dimensions[idd], dimcutn |
---|
| 888 | ## if objvars.dimensions[idd] == dimcutn: |
---|
| 889 | ## cutdim = True |
---|
| 890 | ## break |
---|
| 891 | ## if not cutdim: dimnamesv.append(objvars.dimensions[idd]) |
---|
| 892 | dimnamesv = [vdimyn, vdimxn] |
---|
| 893 | |
---|
| 894 | varunits = objvars.getncattr('units') |
---|
| 895 | |
---|
| 896 | if not objsf.variables.has_key(vdimxn): |
---|
| 897 | print errormsg |
---|
| 898 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 899 | '" does not have dimension variable "' + vdimxn + '" !!' |
---|
| 900 | quit(-1) |
---|
| 901 | if not objsf.variables.has_key(vdimyn): |
---|
| 902 | print errormsg |
---|
| 903 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 904 | '" does not have dimensino variable "' + vdimyn + '" !!' |
---|
| 905 | quit(-1) |
---|
| 906 | |
---|
| 907 | objdimx = objsf.variables[vdimxn] |
---|
| 908 | objdimy = objsf.variables[vdimyn] |
---|
| 909 | odimxu = objdimx.getncattr('units') |
---|
| 910 | odimyu = objdimy.getncattr('units') |
---|
| 911 | |
---|
| 912 | if len(objdimx.shape) <= 2: |
---|
| 913 | # odimxv = objdimx[valshad.shape] |
---|
| 914 | # odimyv = objdimy[valshad.shape] |
---|
| 915 | odimxv = objdimx[:] |
---|
| 916 | odimyv = objdimy[:] |
---|
| 917 | |
---|
| 918 | elif len(objdimx.shape) == 3: |
---|
| 919 | # dimcut = [0, slice(0,valshad.shape[0]), slice(0,valshad.shape[1])] |
---|
| 920 | # odimxv = objdimx[tuple(dimcut)] |
---|
| 921 | # odimyv = objdimy[tuple(dimcut)] |
---|
| 922 | odimxv = objdimx[0,:] |
---|
| 923 | odimyv = objdimy[0,:] |
---|
| 924 | else: |
---|
| 925 | print errormsg |
---|
| 926 | print ' ' + fname + ': shape of dimension variable:', objdimx.shape, \ |
---|
| 927 | ' not ready!!' |
---|
| 928 | quit(-1) |
---|
| 929 | |
---|
| 930 | shading_nx = np.zeros((2), dtype=np.float) |
---|
| 931 | shading_nx[0] = np.float(shadminmax.split(',')[0]) |
---|
| 932 | shading_nx[1] = np.float(shadminmax.split(',')[1]) |
---|
| 933 | |
---|
| 934 | if mapvalue == 'None': mapvalue = None |
---|
| 935 | |
---|
| 936 | # line plot |
---|
| 937 | ## |
---|
| 938 | ncfilel = ncfile.split(',')[1] |
---|
| 939 | vnamelfig = values.split(':')[0].split(',')[1] |
---|
| 940 | varnl = varn.split(',')[1] |
---|
| 941 | colline = values.split(':')[4].split(',')[1] |
---|
| 942 | |
---|
| 943 | objlf = NetCDFFile(ncfilel,'r') |
---|
| 944 | objlvar = objlf.variables[varnl] |
---|
| 945 | |
---|
| 946 | linevals = objlvar[:] |
---|
| 947 | varlunits = objlvar.units |
---|
| 948 | |
---|
| 949 | drw.plot_2D_shadow_line(valshad, linevals, vnamesfig, vnamelfig, odimxv, odimyv, \ |
---|
| 950 | odimxu, odimyu, dimnamesv, colbarn, colline, shading_nx, varunits, varlunits, \ |
---|
| 951 | figtitle, figkind, revals, mapvalue, True) |
---|
| 952 | |
---|
| 953 | objsf.close() |
---|
| 954 | objlf.close() |
---|
| 955 | |
---|
| 956 | return |
---|
| 957 | |
---|
| 958 | def draw_2D_shad_line_time(ncfile, values, varn): |
---|
| 959 | """ plotting a fields with shading and a line with time values |
---|
| 960 | draw_2D_shad_line(ncfile, values, varn) |
---|
| 961 | ncfile= [ncfiles],[ncfilel] files to use to draw with shading and the line |
---|
| 962 | values= [vnamefs],[vanemefl]:[dimvals]:[dimxvn]:[dimyvn]:[colorbar]:[sminv],[smaxv]:[figt]: |
---|
| 963 | [kindfig]:[reverse]:[timevals]:[close] |
---|
| 964 | [vnamefs]: Name in the figure of the variable to be shaded |
---|
| 965 | [vnamefl]: Name in the figure of the variable to be lined |
---|
| 966 | [dimvals]: ',' list of [dimname]|[value] telling at which dimension of the |
---|
| 967 | variable a given value is required (-1, all the length) |
---|
| 968 | [dimx/yvn]: name of the variables with the values of the final dimensions (x,y) |
---|
| 969 | [colorbar]: name of the color bar |
---|
| 970 | [smin/axv]: minimum and maximum value for the shading |
---|
| 971 | [figt]: title of the figure ('|' for spaces) |
---|
| 972 | [kindfig]: kind of figure |
---|
| 973 | [reverse]: Transformation of the values |
---|
| 974 | * 'transpose': reverse the axes (x-->y, y-->x) |
---|
| 975 | * 'flip'@[x/y]: flip the axis x or y |
---|
| 976 | [timevals]: [timen]|[units]|[kind]|[tfmt]|[label]|[timeaxis] time labels characteristics |
---|
| 977 | [timen]; name of the time variable |
---|
| 978 | [units]; units string according to CF conventions ([tunits] since |
---|
| 979 | [YYYY]-[MM]-[DD] [[HH]:[MI]:[SS]], '!' for spaces) |
---|
| 980 | [kind]; kind of output |
---|
| 981 | 'Nval': according to a given number of values as 'Nval',[Nval] |
---|
| 982 | 'exct': according to an exact time unit as 'exct',[tunit]; |
---|
| 983 | tunit= [Nunits],[tu]; [tu]= 'c': centuries, 'y': year, 'm': month, |
---|
| 984 | 'w': week, 'd': day, 'h': hour, 'i': minute, 's': second, |
---|
| 985 | 'l': milisecond |
---|
| 986 | [tfmt]; desired format |
---|
| 987 | [label]; label at the graph ('!' for spaces) |
---|
| 988 | [close]: should figure be closed (finished) |
---|
| 989 | values='dtcon,prc:Time|-1,bottom_top|-1:presmean:time:seismic:-3.e-6,3.e-6:monthly|' |
---|
| 990 | 'dtcon:pdf:transpose:time|hours!since!1949-12-01|exct,2,d|%d|date!([DD])|x:True |
---|
| 991 | varn= [varsn].[varln] name of the variable to plot with shading and to plot with line |
---|
| 992 | """ |
---|
| 993 | fname = 'draw_2D_shad_line_time' |
---|
| 994 | if values == 'h': |
---|
| 995 | print fname + '_____________________________________________________________' |
---|
| 996 | print draw_2D_shad__line_time.__doc__ |
---|
| 997 | quit() |
---|
| 998 | |
---|
| 999 | farguments = ['[vnamefs],[vanemefl]', '[dimvals]', '[dimxvn]', '[dimyvn]', \ |
---|
| 1000 | '[colorbar]', '[sminv],[smaxv]', '[figt]', '[kindfig]', '[reverse]', \ |
---|
| 1001 | '[timevals]', '[close]'] |
---|
| 1002 | drw.check_arguments(fname,11,values,':',farguments) |
---|
| 1003 | |
---|
| 1004 | vnamesfig = values.split(':')[0].split(',')[0] |
---|
| 1005 | dimvals= values.split(':')[1].replace('|',':') |
---|
| 1006 | vdimxn = values.split(':')[2] |
---|
| 1007 | vdimyn = values.split(':')[3] |
---|
| 1008 | colbarn = values.split(':')[4] |
---|
| 1009 | shadminmax = values.split(':')[5] |
---|
| 1010 | figtitle = values.split(':')[6].replace('|',' ') |
---|
| 1011 | figkind = values.split(':')[7] |
---|
| 1012 | revals = values.split(':')[8] |
---|
| 1013 | timevals = values.split(':')[9] |
---|
| 1014 | close = values.split(':')[10] |
---|
| 1015 | |
---|
| 1016 | ncfiles = ncfile.split(',')[0] |
---|
| 1017 | |
---|
| 1018 | if not os.path.isfile(ncfiles): |
---|
| 1019 | print errormsg |
---|
| 1020 | print ' ' + fname + ': shading file "' + ncfiles + '" does not exist !!' |
---|
| 1021 | quit(-1) |
---|
| 1022 | |
---|
| 1023 | objsf = NetCDFFile(ncfiles, 'r') |
---|
| 1024 | |
---|
| 1025 | varns = varn.split(',')[0] |
---|
| 1026 | |
---|
| 1027 | if not objsf.variables.has_key(varns): |
---|
| 1028 | print errormsg |
---|
| 1029 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 1030 | '" does not have variable "' + varns + '" !!' |
---|
| 1031 | quit(-1) |
---|
| 1032 | |
---|
| 1033 | # Variables' values |
---|
| 1034 | objvars = objsf.variables[varns] |
---|
| 1035 | |
---|
| 1036 | valshad, dimsshad = drw.slice_variable(objvars, dimvals.replace(',','|')) |
---|
| 1037 | |
---|
| 1038 | dimnamesv = [vdimyn, vdimxn] |
---|
| 1039 | |
---|
| 1040 | varunits = objvars.getncattr('units') |
---|
| 1041 | |
---|
| 1042 | if not objsf.variables.has_key(vdimxn): |
---|
| 1043 | print errormsg |
---|
| 1044 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 1045 | '" does not have dimension variable "' + vdimxn + '" !!' |
---|
| 1046 | quit(-1) |
---|
| 1047 | if not objsf.variables.has_key(vdimyn): |
---|
| 1048 | print errormsg |
---|
| 1049 | print ' ' + fname + ': shading file "' + ncfiles + \ |
---|
| 1050 | '" does not have dimensino variable "' + vdimyn + '" !!' |
---|
| 1051 | quit(-1) |
---|
| 1052 | |
---|
| 1053 | objdimx = objsf.variables[vdimxn] |
---|
| 1054 | objdimy = objsf.variables[vdimyn] |
---|
| 1055 | odimxu = objdimx.getncattr('units') |
---|
| 1056 | odimyu = objdimy.getncattr('units') |
---|
| 1057 | |
---|
| 1058 | if len(objdimx.shape) <= 2: |
---|
| 1059 | odimxv = objdimx[:] |
---|
| 1060 | odimyv = objdimy[:] |
---|
| 1061 | |
---|
| 1062 | elif len(objdimx.shape) == 3: |
---|
| 1063 | odimxv = objdimx[0,:] |
---|
| 1064 | odimyv = objdimy[0,:] |
---|
| 1065 | else: |
---|
| 1066 | print errormsg |
---|
| 1067 | print ' ' + fname + ': shape of dimension variable:', objdimx.shape, \ |
---|
| 1068 | ' not ready!!' |
---|
| 1069 | quit(-1) |
---|
| 1070 | |
---|
| 1071 | timename = timevals.split('|')[0] |
---|
| 1072 | timeunit = timevals.split('|')[1].replace('!',' ') |
---|
| 1073 | timekind = timevals.split('|')[2] |
---|
| 1074 | timefmt = timevals.split('|')[3] |
---|
| 1075 | timelabel = timevals.split('|')[4].replace('!',' ') |
---|
| 1076 | |
---|
| 1077 | if vdimxn == timename: |
---|
| 1078 | odimxv = objsf.variables[vdimxn][:] |
---|
| 1079 | odimxu = timelabel |
---|
| 1080 | timeaxis = 'x' |
---|
| 1081 | odimyv = objsf.variables[vdimyn] |
---|
| 1082 | odimyu = odimyv.getncattr('units') |
---|
| 1083 | timepos, timelabels = drw.CFtimes_plot(odimxv, timeunit, timekind, timefmt) |
---|
| 1084 | elif vdimyn == timename: |
---|
| 1085 | odimyv = objsf.variables[vdimyn][:] |
---|
| 1086 | odimyu = timelabel |
---|
| 1087 | timeaxis = 'y' |
---|
| 1088 | odimxv = objsf.variables[vdimxn] |
---|
| 1089 | odimxu = odimxv.getncattr('units') |
---|
| 1090 | timepos, timelabels = drw.CFtimes_plot(odimyv, timeunit, timekind, timefmt) |
---|
| 1091 | else: |
---|
| 1092 | print errormsg |
---|
| 1093 | print ' ' + fname + ": time variable '" + timename + "' not found!!" |
---|
| 1094 | quit(-1) |
---|
| 1095 | |
---|
| 1096 | shading_nx = np.zeros((2), dtype=np.float) |
---|
| 1097 | shading_nx[0] = np.float(shadminmax.split(',')[0]) |
---|
| 1098 | shading_nx[1] = np.float(shadminmax.split(',')[1]) |
---|
| 1099 | |
---|
| 1100 | closeval = drw.Str_Bool(close) |
---|
| 1101 | |
---|
| 1102 | drw.plot_2D_shadow_time(valshad, vnamesfig, odimxv, odimyv, odimxu, odimyu, \ |
---|
| 1103 | dimnamesv, colbarn, shading_nx, varunits, figtitle, figkind, revals, timeaxis, \ |
---|
| 1104 | timepos, timelabels, False) |
---|
| 1105 | |
---|
| 1106 | # Line values |
---|
| 1107 | ## |
---|
| 1108 | ncfilel = ncfile.split(',')[1] |
---|
| 1109 | |
---|
| 1110 | vnamelfig = values.split(':')[0].split(',')[1] |
---|
| 1111 | varnl = varn.split(',')[1] |
---|
| 1112 | |
---|
| 1113 | objlf = NetCDFFile(ncfilel,'r') |
---|
| 1114 | objlvar = objlf.variables[varnl] |
---|
| 1115 | |
---|
| 1116 | linevals = objlvar[:] |
---|
| 1117 | if reva0 == 'tranpose': |
---|
| 1118 | plt.plot (linevals, odimxv, '-', color='k') |
---|
| 1119 | else: |
---|
| 1120 | plt.plot (odimxv, linevals, '-', color='k') |
---|
| 1121 | |
---|
| 1122 | objsf.close() |
---|
| 1123 | objsl.close() |
---|
| 1124 | |
---|
| 1125 | return |
---|
| 1126 | |
---|
| 1127 | def draw_topo_geogrid(ncfile, values): |
---|
| 1128 | """ plotting geo_em.d[nn].nc topography from WPS files |
---|
| 1129 | draw_topo_geogrid(ncfile, values) |
---|
| 1130 | ncfile= geo_em.d[nn].nc file to use |
---|
| 1131 | values= [minTopo],[maxTopo]:[lonlatL]:[title]:[graphic_kind]:[mapvalues] |
---|
| 1132 | [min/max]Topo: minimum and maximum values of topography to draw |
---|
| 1133 | lonlatL: limits of longitudes and latitudes [lonmin, latmin, lonmax, latmax] or None |
---|
| 1134 | title: title of the graph |
---|
| 1135 | graphic_kind: kind of figure (jpg, pdf, png) |
---|
| 1136 | mapvalues: map characteristics [proj],[res] |
---|
| 1137 | see full documentation: http://matplotlib.org/basemap/ |
---|
| 1138 | [proj]: projection |
---|
| 1139 | * 'cyl', cilindric |
---|
| 1140 | * 'lcc', lambert conformal |
---|
| 1141 | [res]: resolution: |
---|
| 1142 | * 'c', crude |
---|
| 1143 | * 'l', low |
---|
| 1144 | * 'i', intermediate |
---|
| 1145 | * 'h', high |
---|
| 1146 | * 'f', full |
---|
| 1147 | """ |
---|
| 1148 | fname = 'draw_topo_geogrid' |
---|
| 1149 | |
---|
| 1150 | if values == 'h': |
---|
| 1151 | print fname + '_____________________________________________________________' |
---|
| 1152 | print draw_topo_geogrid.__doc__ |
---|
| 1153 | quit() |
---|
| 1154 | |
---|
| 1155 | expectargs = ['[minTopo]','[maxTopo]', '[lonlatL]', '[title]', '[graphic_kind]', \ |
---|
| 1156 | '[mapvalues]'] |
---|
| 1157 | |
---|
| 1158 | drw.check_arguments(fname,5,values,':',expectargs) |
---|
| 1159 | |
---|
| 1160 | mintopo = values.split(':')[0].split(',')[0] |
---|
| 1161 | maxtopo = values.split(':')[0].split(',')[1] |
---|
| 1162 | |
---|
| 1163 | lonlatLS = values.split(':')[1] |
---|
| 1164 | lonlatLv = lonlatLS.split(',')[0] |
---|
| 1165 | |
---|
| 1166 | if lonlatLv == 'None': |
---|
| 1167 | lonlatL = None |
---|
| 1168 | else: |
---|
| 1169 | lonlatL = np.zeros((4), dtype=np.float) |
---|
| 1170 | lonlatL[0] = np.float(lonlatLS.split(',')[0]) |
---|
| 1171 | lonlatL[1] = np.float(lonlatLS.split(',')[1]) |
---|
| 1172 | lonlatL[2] = np.float(lonlatLS.split(',')[2]) |
---|
| 1173 | lonlatL[3] = np.float(lonlatLS.split(',')[3]) |
---|
| 1174 | |
---|
| 1175 | grtit = values.split(':')[2] |
---|
| 1176 | kindfig = values.split(':')[3] |
---|
| 1177 | mapvalues = values.split(':')[4] |
---|
| 1178 | |
---|
| 1179 | if not os.path.isfile(ncfile): |
---|
| 1180 | print errormsg |
---|
| 1181 | print ' ' + fname + ': domain file "' + ncfile + '" does not exist !!' |
---|
| 1182 | quit(-1) |
---|
| 1183 | |
---|
| 1184 | objdomf = NetCDFFile(ncfile, 'r') |
---|
| 1185 | |
---|
| 1186 | objhgt = objdomf.variables['HGT_M'] |
---|
| 1187 | objlon = objdomf.variables['XLONG_M'] |
---|
| 1188 | objlat = objdomf.variables['XLAT_M'] |
---|
| 1189 | |
---|
| 1190 | topography = objhgt[0,:,:] |
---|
| 1191 | |
---|
| 1192 | drw.plot_topo_geogrid(topography, objlon, objlat, mintopo, maxtopo, lonlatL, \ |
---|
| 1193 | grtit, kindfig, mapvalues, True) |
---|
| 1194 | |
---|
| 1195 | objdomf.close() |
---|
| 1196 | |
---|
| 1197 | return |
---|
| 1198 | |
---|
| 1199 | def draw_topo_geogrid_boxes(ncfiles, values): |
---|
| 1200 | """ plotting different geo_em.d[nn].nc topography from WPS files |
---|
| 1201 | draw_topo_geogrid_boxes(ncfile, values) |
---|
| 1202 | ncfiles= ',' list of geo_em.d[nn].nc files to use (fisrt as topographyc reference) |
---|
| 1203 | values= [minTopo],[maxTopo]:[lonlatL]:[title]:[graphic_kind]:[mapvalues]:[labels] |
---|
| 1204 | [min/max]Topo: minimum and maximum values of topography to draw |
---|
| 1205 | lonlatL: limits of longitudes and latitudes [lonmin, latmin, lonmax, latmax] or None |
---|
| 1206 | title: title of the graph |
---|
| 1207 | graphic_kind: kind of figure (jpg, pdf, png) |
---|
| 1208 | mapvalues: map characteristics [proj],[res] |
---|
| 1209 | see full documentation: http://matplotlib.org/basemap/ |
---|
| 1210 | [proj]: projection |
---|
| 1211 | * 'cyl', cilindric |
---|
| 1212 | * 'lcc', lambert conformal |
---|
| 1213 | [res]: resolution: |
---|
| 1214 | * 'c', crude |
---|
| 1215 | * 'l', low |
---|
| 1216 | * 'i', intermediate |
---|
| 1217 | * 'h', high |
---|
| 1218 | * 'f', full |
---|
| 1219 | labels= labels to write in the graph |
---|
| 1220 | """ |
---|
| 1221 | # import matplotlib as mpl |
---|
| 1222 | # mpl.use('Agg') |
---|
| 1223 | import matplotlib.pyplot as plt |
---|
| 1224 | |
---|
| 1225 | fname = 'draw_topo_geogrid_boxes' |
---|
| 1226 | |
---|
| 1227 | if values == 'h': |
---|
| 1228 | print fname + '_____________________________________________________________' |
---|
| 1229 | print draw_topo_geogrid_boxes.__doc__ |
---|
| 1230 | quit() |
---|
| 1231 | |
---|
| 1232 | mintopo = values.split(':')[0].split(',')[0] |
---|
| 1233 | maxtopo = values.split(':')[0].split(',')[1] |
---|
| 1234 | |
---|
| 1235 | lonlatLS = values.split(':')[1] |
---|
| 1236 | lonlatLv = lonlatLS.split(',')[0] |
---|
| 1237 | |
---|
| 1238 | if lonlatLv == 'None': |
---|
| 1239 | lonlatL = None |
---|
| 1240 | else: |
---|
| 1241 | lonlatL = np.zeros((4), dtype=np.float) |
---|
| 1242 | lonlatL[0] = np.float(lonlatLS.split(',')[0]) |
---|
| 1243 | lonlatL[1] = np.float(lonlatLS.split(',')[1]) |
---|
| 1244 | lonlatL[2] = np.float(lonlatLS.split(',')[2]) |
---|
| 1245 | lonlatL[3] = np.float(lonlatLS.split(',')[3]) |
---|
| 1246 | |
---|
| 1247 | grtit = values.split(':')[2] |
---|
| 1248 | kindfig = values.split(':')[3] |
---|
| 1249 | mapvalues = values.split(':')[4] |
---|
| 1250 | labels = values.split(':')[5] |
---|
| 1251 | |
---|
| 1252 | ncfile = ncfiles.split(',')[0] |
---|
| 1253 | if not os.path.isfile(ncfile): |
---|
| 1254 | print errormsg |
---|
| 1255 | print ' ' + fname + ': domain file "' + ncfile + '" does not exist !!' |
---|
| 1256 | quit(-1) |
---|
| 1257 | |
---|
| 1258 | objdomf = NetCDFFile(ncfile, 'r') |
---|
| 1259 | |
---|
| 1260 | objhgt = objdomf.variables['HGT_M'] |
---|
| 1261 | objlon0 = objdomf.variables['XLONG_M'] |
---|
| 1262 | objlat0 = objdomf.variables['XLAT_M'] |
---|
| 1263 | |
---|
| 1264 | topography = objhgt[0,:,:] |
---|
| 1265 | |
---|
| 1266 | Nfiles = len(ncfiles.split(',')) |
---|
| 1267 | boxlabels = labels.split(',') |
---|
| 1268 | |
---|
| 1269 | Xboxlines = [] |
---|
| 1270 | Yboxlines = [] |
---|
| 1271 | |
---|
| 1272 | for ifile in range(Nfiles): |
---|
| 1273 | ncfile = ncfiles.split(',')[ifile] |
---|
| 1274 | # print ifile, ncfile |
---|
| 1275 | if not os.path.isfile(ncfile): |
---|
| 1276 | print errormsg |
---|
| 1277 | print ' ' + fname + ': domain file "' + ncfile + '" does not exist !!' |
---|
| 1278 | quit(-1) |
---|
| 1279 | |
---|
| 1280 | objdomfi = NetCDFFile(ncfile, 'r') |
---|
| 1281 | |
---|
| 1282 | objlon = objdomfi.variables['XLONG_M'] |
---|
| 1283 | objlat = objdomfi.variables['XLAT_M'] |
---|
| 1284 | |
---|
| 1285 | dx = objlon.shape[2] |
---|
| 1286 | dy = objlon.shape[1] |
---|
| 1287 | |
---|
| 1288 | Xboxlines.append(objlon[0,0,:]) |
---|
| 1289 | Yboxlines.append(objlat[0,0,:]) |
---|
| 1290 | Xboxlines.append(objlon[0,dy-1,:]) |
---|
| 1291 | Yboxlines.append(objlat[0,dy-1,:]) |
---|
| 1292 | Xboxlines.append(objlon[0,:,0]) |
---|
| 1293 | Yboxlines.append(objlat[0,:,0]) |
---|
| 1294 | Xboxlines.append(objlon[0,:,dx-1]) |
---|
| 1295 | Yboxlines.append(objlat[0,:,dx-1]) |
---|
| 1296 | |
---|
| 1297 | objdomfi.close() |
---|
| 1298 | |
---|
| 1299 | drw.plot_topo_geogrid_boxes(topography, Xboxlines, Yboxlines, boxlabels, \ |
---|
| 1300 | objlon0, objlat0, mintopo, maxtopo, lonlatL, grtit, kindfig, mapvalues, True) |
---|
| 1301 | |
---|
| 1302 | objdomf.close() |
---|
| 1303 | |
---|
| 1304 | return |
---|
| 1305 | |
---|
| 1306 | def movievalslice(origslice, dimmovien, framenum): |
---|
| 1307 | """ Function to provide variable slice according to a geneation of a movie |
---|
| 1308 | movievals(origslice, dimmovien, framenum) |
---|
| 1309 | [origslice]= slice original as [dimname1]|[val1],[...,[dimnameN]|[valN]] |
---|
| 1310 | ([val] = -1, full length) |
---|
| 1311 | [dimmovien]= name of the dimension to produce the movie |
---|
| 1312 | [framenum]= value of the frame to substitue in [origslice] as |
---|
| 1313 | [dimmovien]|[framenum] |
---|
| 1314 | >>> movievalslice('East_West|-1,North_South|-1,Time|2','Time',0) |
---|
| 1315 | East_West|-1,North_South|-1,Time|0 |
---|
| 1316 | """ |
---|
| 1317 | |
---|
| 1318 | fname = 'movievalslice' |
---|
| 1319 | |
---|
| 1320 | if origslice == 'h': |
---|
| 1321 | print fname + '_____________________________________________________________' |
---|
| 1322 | print movievalslice.__doc__ |
---|
| 1323 | quit() |
---|
| 1324 | |
---|
| 1325 | dims = origslice.split(',') |
---|
| 1326 | |
---|
| 1327 | movieslice = '' |
---|
| 1328 | idim = 0 |
---|
| 1329 | |
---|
| 1330 | for dimn in dims: |
---|
| 1331 | dn = dimn.split('|')[0] |
---|
| 1332 | if dn == dimmovien: |
---|
| 1333 | movieslice = movieslice + dn + '|' + str(framenum) |
---|
| 1334 | else: |
---|
| 1335 | movieslice = movieslice + dimn |
---|
| 1336 | if idim < len(dims)-1: movieslice = movieslice + ',' |
---|
| 1337 | |
---|
| 1338 | idim = idim + 1 |
---|
| 1339 | |
---|
| 1340 | return movieslice |
---|
| 1341 | |
---|
| 1342 | class Capturing(list): |
---|
| 1343 | """ Class to capture function output as a list |
---|
| 1344 | from: http://stackoverflow.com/questions/16571150/how-to-capture-stdout-output-from-a-python-function-call |
---|
| 1345 | """ |
---|
| 1346 | # from cStringIO import StringIO |
---|
| 1347 | |
---|
| 1348 | def __enter__(self): |
---|
| 1349 | self._stdout = sys.stdout |
---|
| 1350 | sys.stdout = self._stringio = StringIO() |
---|
| 1351 | return self |
---|
| 1352 | def __exit__(self, *args): |
---|
| 1353 | self.extend(self._stringio.getvalue().splitlines()) |
---|
| 1354 | sys.stdout = self._stdout |
---|
| 1355 | |
---|
| 1356 | def create_movie(netcdfile, values, variable): |
---|
| 1357 | """ Function to create a movie assuming ImageMagick installed! |
---|
| 1358 | values= [graph]#[movie_dimension]#[graph_values] |
---|
| 1359 | [graph]: which graphic |
---|
| 1360 | [movie_dimension]: [dimnmovie]@[dimvmovie]@[moviedelay]@[interval] |
---|
| 1361 | [dimnmovie]; name of the dimension from which make the movie |
---|
| 1362 | [dimvmovie]; name of the variable with the values of the dimension |
---|
| 1363 | [moviedelay]; delay between frames |
---|
| 1364 | [interval]; [beg]@[end]@[freq] or -1 (all) |
---|
| 1365 | [graph_values]: values to generate the graphic |
---|
| 1366 | netcdfile= netCDF file |
---|
| 1367 | variable= variable to use (when applicable) |
---|
| 1368 | """ |
---|
| 1369 | fname = 'create_movie' |
---|
| 1370 | |
---|
| 1371 | if values == 'h': |
---|
| 1372 | print fname + '_____________________________________________________________' |
---|
| 1373 | print create_movie.__doc__ |
---|
| 1374 | quit() |
---|
| 1375 | |
---|
| 1376 | graph = values.split('#')[0] |
---|
| 1377 | movie_dim = values.split('#')[1] |
---|
| 1378 | graph_vals = values.split('#')[2] |
---|
| 1379 | |
---|
| 1380 | ncobj = NetCDFFile(netcdfile, 'r') |
---|
| 1381 | |
---|
| 1382 | # Movie dimension |
---|
| 1383 | ## |
---|
| 1384 | dimnmovie = movie_dim.split('@')[0] |
---|
| 1385 | dimvmovie = movie_dim.split('@')[1] |
---|
| 1386 | moviedelay = movie_dim.split('@')[2] |
---|
| 1387 | moviebeg = int(movie_dim.split('@')[3]) |
---|
| 1388 | |
---|
| 1389 | if not drw.searchInlist(ncobj.dimensions.keys(),dimnmovie): |
---|
| 1390 | print errormsg |
---|
| 1391 | print ' ' + fname + ": file '" + netcdfile + "' has not dimension '" + \ |
---|
| 1392 | dimnmovie + "' !!!" |
---|
| 1393 | quit(-1) |
---|
| 1394 | |
---|
| 1395 | objdmovie = ncobj.dimensions[dimnmovie] |
---|
| 1396 | dmovie = len(objdmovie) |
---|
| 1397 | if moviebeg != -1: |
---|
| 1398 | moviend = int(movie_dim.split('@')[4]) |
---|
| 1399 | moviefreq = int(movie_dim.split('@')[5]) |
---|
| 1400 | else: |
---|
| 1401 | moviebeg = 0 |
---|
| 1402 | moviend = dmovie |
---|
| 1403 | moviefreq = 1 |
---|
| 1404 | |
---|
| 1405 | if dimvmovie == 'WRFTimes': |
---|
| 1406 | objvdmovie = ncobj.variables['Times'] |
---|
| 1407 | vdmovieunits = '' |
---|
| 1408 | valsdmovie = [] |
---|
| 1409 | for it in range(objvdmovie.shape[0]): |
---|
| 1410 | valsdmovie.append(drw.datetimeStr_conversion(objvdmovie[it,:], \ |
---|
| 1411 | 'WRFdatetime', 'Y/m/d H-M-S')) |
---|
| 1412 | elif dimvmovie == 'CFtime': |
---|
| 1413 | objvdmovie = ncobj.variables['time'] |
---|
| 1414 | vdmovieunits = '' |
---|
| 1415 | print objvdmovie.units |
---|
| 1416 | valsdmovie0 = drw.netCDFdatetime_realdatetime(objvdmovie.units, 'standard', \ |
---|
| 1417 | objvdmovie[:]) |
---|
| 1418 | valsdmovie = [] |
---|
| 1419 | for it in range(objvdmovie.shape[0]): |
---|
| 1420 | valsdmovie.append(drw.datetimeStr_conversion(valsdmovie0[it,:], \ |
---|
| 1421 | 'matYmdHMS', 'Y/m/d H-M-S')) |
---|
| 1422 | else: |
---|
| 1423 | if not drw.searchInlist(ncobj.variables.keys(),dimvmovie): |
---|
| 1424 | print errormsg |
---|
| 1425 | print ' ' + fname + ": file '" + netcdfile + "' has not variable '" + \ |
---|
| 1426 | dimvmovie + "' !!!" |
---|
| 1427 | quit(-1) |
---|
| 1428 | vdmovieunits = objvdmovie.getncattr('units') |
---|
| 1429 | objvdmovie = ncobj.variables[dimvmovie] |
---|
| 1430 | if len(objvdmovie.shape) == 1: |
---|
| 1431 | vasldmovie = objvdmovie[:] |
---|
| 1432 | else: |
---|
| 1433 | print errormsg |
---|
| 1434 | print ' ' + fname + ': shape', objvdmovie.shape, 'of variable with ' + \ |
---|
| 1435 | 'dimension movie values not ready!!!' |
---|
| 1436 | quit(-1) |
---|
| 1437 | |
---|
| 1438 | ncobj.close() |
---|
| 1439 | os.system('rm frame_*.png > /dev/null') |
---|
| 1440 | |
---|
| 1441 | # graphic |
---|
| 1442 | ## |
---|
| 1443 | if graph == 'draw_2D_shad': |
---|
| 1444 | graphvals = graph_vals.split(':') |
---|
| 1445 | |
---|
| 1446 | for iframe in range(moviebeg,moviend,moviefreq): |
---|
| 1447 | iframeS = str(iframe).zfill(4) |
---|
| 1448 | |
---|
| 1449 | drw.percendone((iframe-moviebeg)/moviefreq,(moviend-moviebeg)/moviefreq, \ |
---|
| 1450 | 5, 'frames') |
---|
| 1451 | titgraph = dimnmovie + '|=|' + str(valsdmovie[iframe]) + '|' + \ |
---|
| 1452 | vdmovieunits |
---|
| 1453 | |
---|
| 1454 | graphvals[1] = movievalslice(graphvals[1],dimnmovie,iframe) |
---|
| 1455 | graphvals[6] = titgraph |
---|
| 1456 | graphvals[7] = 'png' |
---|
| 1457 | |
---|
| 1458 | graphv = drw.numVector_String(graphvals, ":") |
---|
| 1459 | |
---|
| 1460 | with Capturing() as output: |
---|
| 1461 | draw_2D_shad(netcdfile, graphv, variable) |
---|
| 1462 | |
---|
| 1463 | os.system('mv 2Dfields_shadow.png frame_' + iframeS + '.png') |
---|
| 1464 | else: |
---|
| 1465 | print errormsg |
---|
| 1466 | print ' ' + fname + ": graphic '" + graph + "' not defined !!!" |
---|
| 1467 | quit(-1) |
---|
| 1468 | |
---|
| 1469 | os.system('convert -delay ' + moviedelay + ' -loop 0 frame_*.png create_movie.gif') |
---|
| 1470 | |
---|
| 1471 | print "Succesfuly creation of movie file 'create_movie.gif' !!!" |
---|
| 1472 | |
---|
| 1473 | return |
---|
| 1474 | |
---|
| 1475 | def draw_lines(ncfilens, values, varname): |
---|
| 1476 | """ Function to draw different lines at the same time from different files |
---|
| 1477 | draw_trajectories(ncfilens, values, varname): |
---|
| 1478 | ncfilens= [filen] ',' separated list of netCDF files |
---|
| 1479 | values= [dimvname]:[valuesaxis]:[dimtit]:[leglabels]:[vtit]:[title]:[graphk] |
---|
| 1480 | [dimvname]: name of the variable with he values of the common dimension |
---|
| 1481 | [valuesaxis]: which axis will be used for the values ('x', or 'y') |
---|
| 1482 | [dimtit]: title for the common dimension |
---|
| 1483 | [leglabels]: ',' separated list of names for the legend |
---|
| 1484 | [vartit]: name of the variable in the graph |
---|
| 1485 | [title]: title of the plot |
---|
| 1486 | [graphk]: kind of the graphic |
---|
| 1487 | varname= variable to plot |
---|
| 1488 | values= 'XLAT:x:latitude:32x32:$wss^{*}$:wss Taylor's turbulence term:pdf' |
---|
| 1489 | """ |
---|
| 1490 | |
---|
| 1491 | fname = 'draw_lines' |
---|
| 1492 | |
---|
| 1493 | if values == 'h': |
---|
| 1494 | print fname + '_____________________________________________________________' |
---|
| 1495 | print draw_lines.__doc__ |
---|
| 1496 | quit() |
---|
| 1497 | |
---|
| 1498 | expectargs = '[dimvname]:[valuesaxis]:[dimtit]:[leglabels]:[vtit]:[title]:[graphk]' |
---|
| 1499 | drw.check_arguments(fname,len(expectargs.split(':')),values,':',expectargs) |
---|
| 1500 | |
---|
| 1501 | ncfiles = ncfilens.split(',') |
---|
| 1502 | dimvname = values.split(':')[0] |
---|
| 1503 | valuesaxis = values.split(':')[1] |
---|
| 1504 | dimtit = values.split(':')[2] |
---|
| 1505 | leglabels = values.split(':')[3] |
---|
| 1506 | vartit = values.split(':')[4] |
---|
| 1507 | title = values.split(':')[5] |
---|
| 1508 | graphk = values.split(':')[6] |
---|
| 1509 | |
---|
| 1510 | Nfiles = len(ncfiles) |
---|
| 1511 | |
---|
| 1512 | # Getting trajectotries |
---|
| 1513 | ## |
---|
| 1514 | |
---|
| 1515 | varvalues = [] |
---|
| 1516 | dimvalues = [] |
---|
| 1517 | |
---|
| 1518 | print ' ' + fname |
---|
| 1519 | ifn = 0 |
---|
| 1520 | for ifile in ncfiles: |
---|
| 1521 | filen = ifile.split('@')[0] |
---|
| 1522 | |
---|
| 1523 | print ' filen:',filen |
---|
| 1524 | |
---|
| 1525 | if not os.path.isfile(filen): |
---|
| 1526 | print errormsg |
---|
| 1527 | print ' ' + fname + ": netCDF file '" + filen + "' does not exist !!" |
---|
| 1528 | quit(-1) |
---|
| 1529 | |
---|
| 1530 | objfile = NetCDFFile(filen, 'r') |
---|
| 1531 | |
---|
| 1532 | if not objfile.variables.has_key(dimvname): |
---|
| 1533 | print errormsg |
---|
| 1534 | print ' ' + fname + ": netCDF file '" + filen + \ |
---|
| 1535 | "' does not have variable '" + dimvname + "' !!" |
---|
| 1536 | quit(-1) |
---|
| 1537 | |
---|
| 1538 | if not objfile.variables.has_key(varname): |
---|
| 1539 | print errormsg |
---|
| 1540 | print ' ' + fname + ": netCDF file '" + filen + \ |
---|
| 1541 | "' does not have variable '" + varname + "' !!" |
---|
| 1542 | quit(-1) |
---|
| 1543 | |
---|
| 1544 | vvobj = objfile.variables[varname] |
---|
| 1545 | if len(vvobj.shape) != 1: |
---|
| 1546 | print errormsg |
---|
| 1547 | print ' ' + fname + ': wrong shape:',vvobj.shape," of variable '" + \ |
---|
| 1548 | varname + "' !!" |
---|
| 1549 | quit(-1) |
---|
| 1550 | |
---|
| 1551 | vdobj = objfile.variables[dimvname] |
---|
| 1552 | if len(vdobj.shape) != 1: |
---|
| 1553 | print errormsg |
---|
| 1554 | print ' ' + fname + ': wrong shape:',vdobj.shape," of variable '" + \ |
---|
| 1555 | dimvname + "' !!" |
---|
| 1556 | quit(-1) |
---|
| 1557 | |
---|
| 1558 | varvalues.append(vvobj[:]) |
---|
| 1559 | dimvalues.append(vdobj[:]) |
---|
| 1560 | |
---|
| 1561 | if ifn == 0: |
---|
| 1562 | varunits = vvobj.units |
---|
| 1563 | |
---|
| 1564 | objfile.close() |
---|
| 1565 | |
---|
| 1566 | ifn = ifn + 1 |
---|
| 1567 | |
---|
| 1568 | |
---|
| 1569 | drw.plot_lines(dimvalues, varvalues, valuesaxis, dimtit, leglabels.split(','), \ |
---|
| 1570 | vartit, varunits, title, graphk) |
---|
| 1571 | |
---|
| 1572 | return |
---|
| 1573 | |
---|
| 1574 | def draw_Neighbourghood_evol(filen, values, variable): |
---|
| 1575 | """ Function to draw the temporal evolution of a neighbourghood around a point |
---|
| 1576 | draw_Neighbourghood_evol(filen, values, variable) |
---|
| 1577 | filen= netCDF file name |
---|
| 1578 | values= [gvarname]:[dimsval]:[neigdims]:[Nneig]:[Ncol]:[timetits]:[tkinds]: |
---|
| 1579 | [timefmts]:[gtitle]:[shadxtrms]:[cbar]:[gkind]:[ofile] |
---|
| 1580 | [dimsval]: [dimn1]|[val1]|[dimv1],...,[dimnN]|[valN]|[dimvN] dimension names, values to get |
---|
| 1581 | (-1, for all; no name/value pair given full length) and variable with values of the dimension |
---|
| 1582 | NOTE: when dimsval[X,Y] == neigdims[X,Y], valX,valY --> valX,valY-Nneig/2, valX,valY+Nneig/2 |
---|
| 1583 | [neigdims]: [dimnX],[dimnY] dimensions mnames along which the neigbourghood should be defined |
---|
| 1584 | [Nneig]: Number of grid points of the full side of the box (odd value) |
---|
| 1585 | [Ncol]: Number of columns ('auto': square final plot) |
---|
| 1586 | [gvarname]: name of the variable to appear in the graph |
---|
| 1587 | [timetits]: [titX],[titY] titles of the axes ('|' for spaces) |
---|
| 1588 | [tkinds]: [tkindX]|[tkindY] kinds of time to appear in the graph |
---|
| 1589 | 'Nval': according to a given number of values as 'Nval',[Nval] |
---|
| 1590 | 'exct': according to an exact time unit as 'exct',[tunit]; |
---|
| 1591 | tunit= [Nunits],[tu]; [tu]= 'c': centuries, 'y': year, 'm': month, |
---|
| 1592 | 'w': week, 'd': day, 'h': hour, 'i': minute, 's': second, |
---|
| 1593 | 'l': milisecond |
---|
| 1594 | [timefmts]: [tfmtX],[tfmtY] format of the time labels |
---|
| 1595 | [gtitle]: title of the graphic ('|' for spaces) |
---|
| 1596 | [shadxtrms]: Extremes for the shading |
---|
| 1597 | [cbar]: colorbar to use |
---|
| 1598 | [gkind]: kind of graphical output |
---|
| 1599 | [ofile]: True/False whether the netcdf with data should be created or not |
---|
| 1600 | variable= name of the variable |
---|
| 1601 | values = 'q:Time|-1|Times,bottom_top|6|ZNU,south_north|3|XLAT,west_east|26|XLONG:south_north,west_east:5:auto:time|($[DD]^{[HH]}$),time|($[DD]^{[HH]}$):exct,2,h|exct,1,d:$%d^{%H}$,$%d^{%H}$:5|pts|neighbourghood|temporal|evolution:0.0,0.004:BuPu:pdf:True' |
---|
| 1602 | """ |
---|
| 1603 | |
---|
| 1604 | fname = 'draw_Neighbourghood_evol' |
---|
| 1605 | |
---|
| 1606 | if values == 'h': |
---|
| 1607 | print fname + '_____________________________________________________________' |
---|
| 1608 | print draw_Neighbourghood_evol.__doc__ |
---|
| 1609 | quit() |
---|
| 1610 | |
---|
| 1611 | expectargs = ['[gvarname]', '[dimsval]', '[neigdims]', '[Nneig]', '[Ncol]', \ |
---|
| 1612 | '[timetits]', '[tkinds]', '[timefmts]', '[gtitle]', '[shadxtrms]', '[cbar]', \ |
---|
| 1613 | '[gkind]', '[ofile]'] |
---|
| 1614 | |
---|
| 1615 | drw.check_arguments(fname,len(expectargs),values,':',expectargs) |
---|
| 1616 | |
---|
| 1617 | gvarname = values.split(':')[0] |
---|
| 1618 | dimsval = values.split(':')[1].split(',') |
---|
| 1619 | neigdims = values.split(':')[2].split(',') |
---|
| 1620 | Nneig = int(values.split(':')[3]) |
---|
| 1621 | Ncol0 = values.split(':')[4] |
---|
| 1622 | timetits = values.split(':')[5].split(',') |
---|
| 1623 | timekinds = values.split(':')[6].split('|') |
---|
| 1624 | timefmts = values.split(':')[7].split(',') |
---|
| 1625 | gtitle = values.split(':')[8].replace('|',' ') |
---|
| 1626 | shadxtrms = values.split(':')[9].split(',') |
---|
| 1627 | cbar = values.split(':')[10] |
---|
| 1628 | gkind = values.split(':')[11] |
---|
| 1629 | ofile = values.split(':')[12] |
---|
| 1630 | |
---|
| 1631 | if Ncol0 != 'auto': |
---|
| 1632 | Ncol = int(Ncol0) |
---|
| 1633 | else: |
---|
| 1634 | Ncol = Ncol0 |
---|
| 1635 | |
---|
| 1636 | timetits[0] = timetits[0].replace('|',' ') |
---|
| 1637 | timetits[1] = timetits[1].replace('|',' ') |
---|
| 1638 | |
---|
| 1639 | if np.mod(Nneig,2) == 0: |
---|
| 1640 | print errormsg |
---|
| 1641 | print ' ' + fname + ": an odd value for 'Nneig':",Nneig,'is required !!!' |
---|
| 1642 | quit(-1) |
---|
| 1643 | |
---|
| 1644 | Nneig2 = int(Nneig/2) |
---|
| 1645 | |
---|
| 1646 | # Values to slice the variable |
---|
| 1647 | dimvslice = {} |
---|
| 1648 | dimvvalues = {} |
---|
| 1649 | for dimvs in dimsval: |
---|
| 1650 | dimn = dimvs.split('|')[0] |
---|
| 1651 | dimv = int(dimvs.split('|')[1]) |
---|
| 1652 | dimnv = dimvs.split('|')[2] |
---|
| 1653 | |
---|
| 1654 | dimvvalues[dimn] = dimnv |
---|
| 1655 | dimvslice[dimn] = dimv |
---|
| 1656 | |
---|
| 1657 | ncobj = NetCDFFile(filen, 'r') |
---|
| 1658 | |
---|
| 1659 | varobj = ncobj.variables[variable] |
---|
| 1660 | |
---|
| 1661 | slicevar = [] |
---|
| 1662 | newdimn = [] |
---|
| 1663 | newdimsvar = {} |
---|
| 1664 | |
---|
| 1665 | for dimn in varobj.dimensions: |
---|
| 1666 | if not drw.searchInlist(dimvslice.keys(), dimn): |
---|
| 1667 | dimsize = len(ncobj.dimensions[dimn]) |
---|
| 1668 | slicevar.append(slice(0, dimsize+1)) |
---|
| 1669 | newdimn.append(dimn) |
---|
| 1670 | newdimsvar[dimn] = dimseize |
---|
| 1671 | |
---|
| 1672 | for dimslicen in dimvslice.keys(): |
---|
| 1673 | if dimn == dimslicen: |
---|
| 1674 | if dimvslice[dimn] != -1: |
---|
| 1675 | if drw.searchInlist(neigdims, dimn): |
---|
| 1676 | slicevar.append(slice(dimvslice[dimn]-Nneig2, \ |
---|
| 1677 | dimvslice[dimn]+Nneig2+1)) |
---|
| 1678 | newdimn.append(dimn) |
---|
| 1679 | newdimsvar[dimn] = Nneig |
---|
| 1680 | break |
---|
| 1681 | else: |
---|
| 1682 | slicevar.append(slice(dimvslice[dimn], dimvslice[dimn]+1)) |
---|
| 1683 | break |
---|
| 1684 | else: |
---|
| 1685 | dimsize = len(ncobj.dimensions[dimn]) |
---|
| 1686 | slicevar.append(slice(0, dimsize+1)) |
---|
| 1687 | newdimn.append(dimn) |
---|
| 1688 | newdimsvar[dimn] = dimsize |
---|
| 1689 | break |
---|
| 1690 | |
---|
| 1691 | varv = varobj[tuple(slicevar)] |
---|
| 1692 | |
---|
| 1693 | if len(newdimn) != 3: |
---|
| 1694 | print errormsg |
---|
| 1695 | print ' ' + fname + ': sliced variable with shape=', varv.shape, \ |
---|
| 1696 | ' must have three dimensions',len(varv.shape),'given !!' |
---|
| 1697 | quit(-1) |
---|
| 1698 | |
---|
| 1699 | newdims = [] |
---|
| 1700 | for nwdims in newdimn: |
---|
| 1701 | newdims.append(newdimsvar[nwdims]) |
---|
| 1702 | |
---|
| 1703 | # The dimension which is not in the neighbourhood dimensions must be time! |
---|
| 1704 | for dim1 in newdimn: |
---|
| 1705 | if not drw.searchInlist(neigdims, dim1): |
---|
| 1706 | dimt = newdimsvar[dim1] |
---|
| 1707 | dimtime = dim1 |
---|
| 1708 | |
---|
| 1709 | if Ncol == 'auto': |
---|
| 1710 | dimtsqx = int(np.sqrt(dimt)) + 1 |
---|
| 1711 | dimtsqy = int(np.sqrt(dimt)) + 1 |
---|
| 1712 | else: |
---|
| 1713 | dimtsqx = int(Ncol) |
---|
| 1714 | dimtsqy = dimt/dimtsqx + 1 |
---|
| 1715 | |
---|
| 1716 | neighbourghood = np.ones((dimtsqy*Nneig,dimtsqx*Nneig), dtype=np.float)*fillValue |
---|
| 1717 | |
---|
| 1718 | for it in range(dimt): |
---|
| 1719 | ity = int(it/dimtsqx) |
---|
| 1720 | itx = it-ity*dimtsqx |
---|
| 1721 | |
---|
| 1722 | itty = (dimtsqy - ity - 1)*Nneig + Nneig2 |
---|
| 1723 | ittx = itx*Nneig + Nneig2 |
---|
| 1724 | |
---|
| 1725 | neighbourghood[itty-Nneig2:itty+Nneig2+1,ittx-Nneig2:ittx+Nneig2+1]= \ |
---|
| 1726 | varv[it,::-1,:] |
---|
| 1727 | |
---|
| 1728 | variablevals = drw.variables_values(variable) |
---|
| 1729 | if drw.searchInlist(varobj.ncattrs(), 'units'): |
---|
| 1730 | vunits = varobj.units |
---|
| 1731 | else: |
---|
| 1732 | vunits = variablevals[5] |
---|
| 1733 | |
---|
| 1734 | # Time values at the X/Y axes |
---|
| 1735 | if ncobj.variables[dimvvalues[dimtime]].dtype == '|S1': |
---|
| 1736 | print ' ' + fname + ': WRF time variable!' |
---|
| 1737 | refdate = '19491201000000' |
---|
| 1738 | tunitsval = 'hours' |
---|
| 1739 | dimtvalues = np.zeros((dimt), dtype=np.float) |
---|
| 1740 | tvals = ncobj.variables[dimvvalues[dimtime]] |
---|
| 1741 | yrref=refdate[0:4] |
---|
| 1742 | monref=refdate[4:6] |
---|
| 1743 | dayref=refdate[6:8] |
---|
| 1744 | horref=refdate[8:10] |
---|
| 1745 | minref=refdate[10:12] |
---|
| 1746 | secref=refdate[12:14] |
---|
| 1747 | |
---|
| 1748 | refdateS = yrref + '/' + monref + '/' + dayref + '_' + horref + ':' + \ |
---|
| 1749 | minref + ':' + secref |
---|
| 1750 | tunits = tunitsval + ' since ' + refdateS |
---|
| 1751 | for it in range(dimt): |
---|
| 1752 | wrfdates = drw.datetimeStr_conversion(tvals[it,:],'WRFdatetime', 'matYmdHMS') |
---|
| 1753 | dimtvalues[it] = drw.realdatetime1_CFcompilant(wrfdates, refdate, tunitsval) |
---|
| 1754 | else: |
---|
| 1755 | dimtvalues = ncobj.variables[dimvvalues[dimtime]][:] |
---|
| 1756 | tunits = ncobj.variables[newdimsvar[dimtime]].units |
---|
| 1757 | |
---|
| 1758 | dimxv = dimtvalues[0:dimtsqx] |
---|
| 1759 | dimyv = dimtvalues[0:dimt:dimtsqx] |
---|
| 1760 | |
---|
| 1761 | print 'Lluis lens: dimxv, dimyv:',len(dimxv), len(dimyv) |
---|
| 1762 | |
---|
| 1763 | dimn = ['time','time'] |
---|
| 1764 | |
---|
| 1765 | if ofile == 'True': |
---|
| 1766 | ofilen = 'Neighbourghood_evol.nc' |
---|
| 1767 | newnc = NetCDFFile(ofilen, 'w') |
---|
| 1768 | # Dimensions |
---|
| 1769 | newdim = newnc.createDimension('time',None) |
---|
| 1770 | newdim = newnc.createDimension('y',dimtsqy*Nneig) |
---|
| 1771 | newdim = newnc.createDimension('x',dimtsqx*Nneig) |
---|
| 1772 | # Dimension values |
---|
| 1773 | newvar = newnc.createVariable('time','f8',('time')) |
---|
| 1774 | newvar[:] = np.arange(dimt) |
---|
| 1775 | newattr = drw.basicvardef(newvar, 'time','time',tunits) |
---|
| 1776 | # Neighbourhghood variable |
---|
| 1777 | newvar = newnc.createVariable(variable + 'neigevol', 'f4', ('y','x'), \ |
---|
| 1778 | fill_value=fillValue) |
---|
| 1779 | newvar[:] = neighbourghood |
---|
| 1780 | |
---|
| 1781 | newnc.sync() |
---|
| 1782 | newnc.close() |
---|
| 1783 | print fname + ": Successfull generation of file '" + ofilen + "' !!" |
---|
| 1784 | |
---|
| 1785 | # Time ticks |
---|
| 1786 | timeposX, timelabelsX = drw.CFtimes_plot(dimxv, tunits, timekinds[0], timefmts[0]) |
---|
| 1787 | timeposY, timelabelsY = drw.CFtimes_plot(dimyv, tunits, timekinds[1], timefmts[1]) |
---|
| 1788 | |
---|
| 1789 | timepos = [timeposX, timeposY[::-1]] |
---|
| 1790 | timelabels = [timelabelsX, timelabelsY[::-1]] |
---|
| 1791 | |
---|
| 1792 | for i in range(2): |
---|
| 1793 | if shadxtrms[i][0:1] != 'S': |
---|
| 1794 | shadxtrms[i] = np.float(shadxtrms[i]) |
---|
| 1795 | |
---|
| 1796 | drw.plot_Neighbourghood_evol(neighbourghood, dimxv, dimyv, gvarname, timetits, \ |
---|
| 1797 | timepos, timelabels, cbar, Nneig, shadxtrms, vunits, gtitle, gkind, True) |
---|
| 1798 | |
---|
| 1799 | def draw_timeSeries(filen, values, variables): |
---|
| 1800 | """ Function to draw a time-series |
---|
| 1801 | draw_timeSeries(filen, values, variable): |
---|
| 1802 | filen= name of the file |
---|
| 1803 | values= [gvarname]:[timetit]:[tkind]:[timefmt]:[title]:[locleg]:[gkind] |
---|
| 1804 | [gvarname]: name of the variable to appear in the graph |
---|
| 1805 | [timetit]: title of the time axis (assumed x-axis, '|' for spaces) |
---|
| 1806 | [tkind]: kind of time to appear in the graph (assumed x-axis) |
---|
| 1807 | 'Nval': according to a given number of values as 'Nval',[Nval] |
---|
| 1808 | 'exct': according to an exact time unit as 'exct',[tunit]; |
---|
| 1809 | tunit= [Nunits],[tu]; [tu]= 'c': centuries, 'y': year, 'm': month, |
---|
| 1810 | 'w': week, 'd': day, 'h': hour, 'i': minute, 's': second, |
---|
| 1811 | 'l': milisecond |
---|
| 1812 | [timefmt]: format of the time labels |
---|
| 1813 | [title]: title of the graphic ('|' for spaces) |
---|
| 1814 | [locleg]: location of the legend (-1, autmoatic) |
---|
| 1815 | 1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right', |
---|
| 1816 | 5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center', |
---|
| 1817 | 9: 'upper center', 10: 'center' |
---|
| 1818 | [gkind]: kind of graphical output |
---|
| 1819 | variables= [varname],[timename] names of variable and variable with times |
---|
| 1820 | draw_timeSeries('wrfout_d01_1979-12-01_00:00:00_bottom_top_B6-E6-I1_south_north_B3-E3-I1_west_east_B26-E26-I1.nc', 'dt_con:time|($[DD]^{[HH]}$):exct,12,h:$%d^{%H}$:time|evolution|at|-1|6|3|26:1:pdf', 'LDQCON,time') |
---|
| 1821 | """ |
---|
| 1822 | |
---|
| 1823 | fname = 'draw_timeSeries' |
---|
| 1824 | |
---|
| 1825 | if values == 'h': |
---|
| 1826 | print fname + '_____________________________________________________________' |
---|
| 1827 | print draw_timeSeries.__doc__ |
---|
| 1828 | quit() |
---|
| 1829 | |
---|
| 1830 | expectargs = ['[gvarname]', '[timetit]', '[tkind]', '[timefmt]', '[title]', \ |
---|
| 1831 | '[locleg]', '[gkind]'] |
---|
| 1832 | |
---|
| 1833 | drw.check_arguments(fname,len(expectargs),values,':',expectargs) |
---|
| 1834 | |
---|
| 1835 | gvarname = values.split(':')[0] |
---|
| 1836 | timetit = values.split(':')[1].replace('|',' ') |
---|
| 1837 | tkind = values.split(':')[2] |
---|
| 1838 | timefmt = values.split(':')[3] |
---|
| 1839 | title = values.split(':')[4].replace('|',' ') |
---|
| 1840 | locleg = int(values.split(':')[5]) |
---|
| 1841 | gkind = values.split(':')[6] |
---|
| 1842 | |
---|
| 1843 | ncobj = NetCDFFile(filen, 'r') |
---|
| 1844 | |
---|
| 1845 | variable = variables.split(',')[0] |
---|
| 1846 | timevar = variables.split(',')[1] |
---|
| 1847 | |
---|
| 1848 | if not ncobj.variables.has_key(variable): |
---|
| 1849 | print errormsg |
---|
| 1850 | print ' ' + fname + ": file '" + filen + "' does not have variable '" + \ |
---|
| 1851 | variable + "' !!" |
---|
| 1852 | quit(-1) |
---|
| 1853 | |
---|
| 1854 | if not ncobj.variables.has_key(timevar): |
---|
| 1855 | print errormsg |
---|
| 1856 | print ' ' + fname + ": file '" + filen + "' does not have variable time '" \ |
---|
| 1857 | + timevar + "' !!" |
---|
| 1858 | quit(-1) |
---|
| 1859 | |
---|
| 1860 | varobj = ncobj.variables[variable] |
---|
| 1861 | timeobj = ncobj.variables[timevar] |
---|
| 1862 | |
---|
| 1863 | dimt = len(timeobj[:]) |
---|
| 1864 | varvals = np.zeros((2,dimt), dtype=np.float) |
---|
| 1865 | |
---|
| 1866 | gunits = varobj.getncattr('units') |
---|
| 1867 | tunits = timeobj.getncattr('units') |
---|
| 1868 | |
---|
| 1869 | varvals[0,:], valpot, newgunits, Spot = drw.pot_values(varobj[:].flatten(), gunits) |
---|
| 1870 | varvals[1,:] = timeobj[:] |
---|
| 1871 | |
---|
| 1872 | tseriesvals = [] |
---|
| 1873 | tseriesvals.append(varvals) |
---|
| 1874 | |
---|
| 1875 | drw.plot_TimeSeries(tseriesvals, Spot + drw.units_lunits(gunits), tunits, \ |
---|
| 1876 | 'TimeSeries', gvarname, timetit, tkind, timefmt, title, \ |
---|
| 1877 | [gvarname.replace('_','\_')], locleg, gkind) |
---|
| 1878 | |
---|
| 1879 | return |
---|
| 1880 | |
---|
| 1881 | #draw_timeSeries('wrfout_d01_1979-12-01_00:00:00_bottom_top_B6-E6-I1_south_north_B3-E3-I1_west_east_B26-E26-I1.nc', 'dt_con:time|($[DD]^{[HH]}$):exct,12,h:$%d^{%H}$:time|evolution|at|-1|6|3|26:1:pdf', 'LDQCON,time') |
---|
| 1882 | |
---|
| 1883 | def draw_trajectories(trjfilens, values, observations): |
---|
| 1884 | """ Function to draw different trajectories at the same time |
---|
| 1885 | draw_trajectories(trjfilens, values, observations): |
---|
| 1886 | trjfilens= [filen]@[Tint]@[map] ',' separated list of files with trajectories, |
---|
| 1887 | time intervals and reference maps (first one will be used to plot) |
---|
| 1888 | [filen]: name of the file to use (lines with '#', not readed) as: |
---|
| 1889 | [t-step] [x] [y] |
---|
| 1890 | [Tint]: interval of time as [Tbeg]@[Tend] or -1 for all the interval |
---|
| 1891 | [map]: [file]#[lonname]#[latname] |
---|
| 1892 | [file]; with the [lon],[lat] matrices |
---|
| 1893 | [lonname],[latname]; names of the longitudes and latitudes variables |
---|
| 1894 | values=[leglabels]|[lonlatlims]|[title]|[graphk]|[mapkind] |
---|
| 1895 | [leglabels]: ',' separated list of names for the legend |
---|
| 1896 | [lonlatlims]: limits of the map [lonmin, latmin, lonmax, latmax] or None |
---|
| 1897 | [title]: title of the plot |
---|
| 1898 | [graphk]: kind of the graphic |
---|
| 1899 | [mapkind]: drawing coastaline ([proj],[res]) or None |
---|
| 1900 | [proj]: projection |
---|
| 1901 | * 'cyl', cilindric |
---|
| 1902 | * 'lcc', lambert conformal |
---|
| 1903 | [res]: resolution: |
---|
| 1904 | * 'c', crude |
---|
| 1905 | * 'l', low |
---|
| 1906 | * 'i', intermediate |
---|
| 1907 | * 'h', high |
---|
| 1908 | * 'f', full |
---|
| 1909 | obsevations= [obsfile],[obsname],[Tint],[null] |
---|
| 1910 | [obsfile]: name fo the File with the observations as [t-step] [lat] [lon] |
---|
| 1911 | [obsname]: name of the observations in the graph |
---|
| 1912 | [Tint]: interval of time as [Tbeg]@[Tend] or -1 for all the interval |
---|
| 1913 | [null]: null value for the observed trajectory |
---|
| 1914 | """ |
---|
| 1915 | |
---|
| 1916 | fname = 'draw_trajectories' |
---|
| 1917 | |
---|
| 1918 | if values == 'h': |
---|
| 1919 | print fname + '_____________________________________________________________' |
---|
| 1920 | print draw_trajectories.__doc__ |
---|
| 1921 | quit() |
---|
| 1922 | |
---|
| 1923 | trjfiles = trjfilens.split(',') |
---|
| 1924 | leglabels = values.split('|')[0] |
---|
| 1925 | lonlatlims = values.split('|')[1] |
---|
| 1926 | title = values.split('|')[2] |
---|
| 1927 | graphk = values.split('|')[3] |
---|
| 1928 | mapkind = values.split('|')[4] |
---|
| 1929 | |
---|
| 1930 | Nfiles = len(trjfiles) |
---|
| 1931 | |
---|
| 1932 | # Getting trajectotries |
---|
| 1933 | ## |
---|
| 1934 | |
---|
| 1935 | lontrjvalues = [] |
---|
| 1936 | lattrjvalues = [] |
---|
| 1937 | |
---|
| 1938 | print ' ' + fname |
---|
| 1939 | ifn = 0 |
---|
| 1940 | for ifile in trjfiles: |
---|
| 1941 | filen = ifile.split('@')[0] |
---|
| 1942 | Tint = ifile.split('@')[1] |
---|
| 1943 | |
---|
| 1944 | print ' trajectory:',filen |
---|
| 1945 | |
---|
| 1946 | if Tint != '-1': |
---|
| 1947 | Tbeg = Tint |
---|
| 1948 | Tend = ifile.split('@')[2] |
---|
| 1949 | mapv = ifile.split('@')[3] |
---|
| 1950 | else: |
---|
| 1951 | mapv = ifile.split('@')[2] |
---|
| 1952 | |
---|
| 1953 | if not os.path.isfile(filen): |
---|
| 1954 | print errormsg |
---|
| 1955 | print ' ' + fname + ": trajectory file '" + filen + "' does not exist !!" |
---|
| 1956 | quit(-1) |
---|
| 1957 | |
---|
| 1958 | # Charging longitude and latitude values |
---|
| 1959 | ## |
---|
| 1960 | lonvals, latvals = drw.lonlat_values(mapv.split('#')[0], mapv.split('#')[1], \ |
---|
| 1961 | mapv.split('#')[2]) |
---|
| 1962 | |
---|
| 1963 | if ifn == 0: mapref = mapv |
---|
| 1964 | ifn = ifn + 1 |
---|
| 1965 | |
---|
| 1966 | objfile = open(filen, 'r') |
---|
| 1967 | trjtimev = [] |
---|
| 1968 | trjxv = [] |
---|
| 1969 | trjyv = [] |
---|
| 1970 | |
---|
| 1971 | for line in objfile: |
---|
| 1972 | if line[0:1] != '#': |
---|
| 1973 | trjtimev.append(int(line.split(' ')[0])) |
---|
| 1974 | trjxv.append(int(line.split(' ')[1])) |
---|
| 1975 | trjyv.append(int(line.split(' ')[2])) |
---|
| 1976 | |
---|
| 1977 | objfile.close() |
---|
| 1978 | |
---|
| 1979 | if Tint != '-1': |
---|
| 1980 | lontrjvalues.append(lonvals[trjyv[Tint:Tend+1], trjxv[Tint:Tend+1]]) |
---|
| 1981 | lattrjvalues.append(latvals[trjyv[Tint:Tend+1], trjxv[Tint:Tend+1]]) |
---|
| 1982 | else: |
---|
| 1983 | lontrjvalues.append(lonvals[trjyv[:], trjxv[:]]) |
---|
| 1984 | lattrjvalues.append(latvals[trjyv[:], trjxv[:]]) |
---|
| 1985 | |
---|
| 1986 | # lonlatlimits |
---|
| 1987 | ## |
---|
| 1988 | |
---|
| 1989 | if lonlatlims == 'None': |
---|
| 1990 | lonlatlimsv = None |
---|
| 1991 | else: |
---|
| 1992 | lonlatlimsv = np.zeros((4), dtype=np.float) |
---|
| 1993 | lonlatlimsv[0] = np.float(lonlatlims.split(',')[0]) |
---|
| 1994 | lonlatlimsv[1] = np.float(lonlatlims.split(',')[1]) |
---|
| 1995 | lonlatlimsv[2] = np.float(lonlatlims.split(',')[2]) |
---|
| 1996 | lonlatlimsv[3] = np.float(lonlatlims.split(',')[3]) |
---|
| 1997 | |
---|
| 1998 | # lon/lat objects |
---|
| 1999 | ## |
---|
| 2000 | objnc = NetCDFFile(mapref.split('#')[0]) |
---|
| 2001 | lonobj = objnc.variables[mapref.split('#')[1]] |
---|
| 2002 | latobj = objnc.variables[mapref.split('#')[2]] |
---|
| 2003 | |
---|
| 2004 | # map |
---|
| 2005 | ## |
---|
| 2006 | if mapkind == 'None': |
---|
| 2007 | mapkindv = None |
---|
| 2008 | else: |
---|
| 2009 | mapkindv = mapkind |
---|
| 2010 | |
---|
| 2011 | if observations is None: |
---|
| 2012 | obsname = None |
---|
| 2013 | else: |
---|
| 2014 | obsfile = observations.split(',')[0] |
---|
| 2015 | obsname = observations.split(',')[1] |
---|
| 2016 | Tint = observations.split(',')[2] |
---|
| 2017 | null = np.float(observations.split(',')[3]) |
---|
| 2018 | print ' observational trajectory:',obsfile |
---|
| 2019 | |
---|
| 2020 | if not os.path.isfile(obsfile): |
---|
| 2021 | print errormsg |
---|
| 2022 | print ' ' + fname + ": observations trajectory file '" + obsfile + \ |
---|
| 2023 | "' does not exist !!" |
---|
| 2024 | quit(-1) |
---|
| 2025 | |
---|
| 2026 | objfile = open(obsfile, 'r') |
---|
| 2027 | obstrjtimev = [] |
---|
| 2028 | obstrjxv = [] |
---|
| 2029 | obstrjyv = [] |
---|
| 2030 | |
---|
| 2031 | for line in objfile: |
---|
| 2032 | if line[0:1] != '#': |
---|
| 2033 | lon = np.float(line.split(' ')[2]) |
---|
| 2034 | lat = np.float(line.split(' ')[1]) |
---|
| 2035 | if not lon == null and not lat == null: |
---|
| 2036 | obstrjtimev.append(int(line.split(' ')[0])) |
---|
| 2037 | obstrjxv.append(lon) |
---|
| 2038 | obstrjyv.append(lat) |
---|
| 2039 | else: |
---|
| 2040 | obstrjtimev.append(int(line.split(' ')[0])) |
---|
| 2041 | obstrjxv.append(None) |
---|
| 2042 | obstrjyv.append(None) |
---|
| 2043 | |
---|
| 2044 | objfile.close() |
---|
| 2045 | |
---|
| 2046 | if Tint != '-1': |
---|
| 2047 | Tint = int(observations.split(',')[2].split('@')[0]) |
---|
| 2048 | Tend = int(observations.split(',')[2].split('@')[1]) |
---|
| 2049 | lontrjvalues.append(obstrjxv[Tint:Tend+1]) |
---|
| 2050 | lattrjvalues.append(obstrjyv[Tint:Tend+1]) |
---|
| 2051 | else: |
---|
| 2052 | lontrjvalues.append(obstrjxv[:]) |
---|
| 2053 | lattrjvalues.append(obstrjyv[:]) |
---|
| 2054 | |
---|
| 2055 | drw.plot_Trajectories(lontrjvalues, lattrjvalues, leglabels.split(','), \ |
---|
| 2056 | lonobj, latobj, lonlatlimsv, title, graphk, mapkindv, obsname) |
---|
| 2057 | |
---|
| 2058 | objnc.close() |
---|
| 2059 | |
---|
| 2060 | return |
---|
| 2061 | |
---|
| 2062 | def draw_vals_trajectories(ncfile, values, variable): |
---|
| 2063 | """ Function to draw values from the outputs from 'compute_tevolboxtraj' |
---|
| 2064 | draw_vals_trajectories(ncfile, values, variable) |
---|
| 2065 | ncfile= [ncfile] ',' list of files to use |
---|
| 2066 | values= [statistics]:[Tint]:[labels]@[locleg]:[gvarname]:[timetit]:[tkind]:[timefmt]:[title]:[gkind] |
---|
| 2067 | [statistics]: which statistics to use, from: 'center', 'min', 'max', 'mean', |
---|
| 2068 | 'mean2', 'stdev' |
---|
| 2069 | [Tint]: [Tbeg]@[Tend] or None, interval of time to plot or -1 for all the times |
---|
| 2070 | [labels]: ',' separated list of labels for the legend |
---|
| 2071 | [locleg]: location of the legend (-1, autmoatic) |
---|
| 2072 | 1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right', |
---|
| 2073 | 5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center', |
---|
| 2074 | 9: 'upper center', 10: 'center' |
---|
| 2075 | [gvarname]: name of the variable to appear in the graph |
---|
| 2076 | [timetit]: title of the time axis (assumed x-axis, '|' for spaces) |
---|
| 2077 | [tkind]: kind of time to appear in the graph (assumed x-axis) |
---|
| 2078 | 'Nval': according to a given number of values as 'Nval',[Nval] |
---|
| 2079 | 'exct': according to an exact time unit as 'exct',[tunit]; |
---|
| 2080 | tunit= [Nunits],[tu]; [tu]= 'c': centuries, 'y': year, 'm': month, |
---|
| 2081 | 'w': week, 'd': day, 'h': hour, 'i': minute, 's': second, |
---|
| 2082 | 'l': milisecond |
---|
| 2083 | [timefmt]: format of the time labels |
---|
| 2084 | [title]: title of the graphic ('|' for spaces) |
---|
| 2085 | [gkind]: kind of graphical output |
---|
| 2086 | variable= variable to use |
---|
| 2087 | """ |
---|
| 2088 | |
---|
| 2089 | fname = 'draw_vals_trajectories' |
---|
| 2090 | |
---|
| 2091 | if values == 'h': |
---|
| 2092 | print fname + '_____________________________________________________________' |
---|
| 2093 | print draw_vals_trajectories.__doc__ |
---|
| 2094 | quit() |
---|
| 2095 | |
---|
| 2096 | sims = ncfile.split(',') |
---|
| 2097 | |
---|
| 2098 | if len(values.split(':')) != 9: |
---|
| 2099 | print errormsg |
---|
| 2100 | print ' ' + fname + ': wrong number of values!', len(values.split(':')), \ |
---|
| 2101 | 'given 9 needed!!' |
---|
| 2102 | print ' ',values.split(':') |
---|
| 2103 | quit(-1) |
---|
| 2104 | |
---|
| 2105 | statistics = values.split(':')[0] |
---|
| 2106 | Tint = values.split(':')[1] |
---|
| 2107 | labels = values.split(':')[2] |
---|
| 2108 | gvarname = values.split(':')[3] |
---|
| 2109 | timetit = values.split(':')[4].replace('|',' ') |
---|
| 2110 | tkind = values.split(':')[5] |
---|
| 2111 | timefmt = values.split(':')[6] |
---|
| 2112 | title = values.split(':')[7].replace('|',' ') |
---|
| 2113 | gkind = values.split(':')[8] |
---|
| 2114 | |
---|
| 2115 | leglabels = labels.split('@')[0].split(',') |
---|
| 2116 | locleg = int(labels.split('@')[1]) |
---|
| 2117 | |
---|
| 2118 | Nsims = len(sims) |
---|
| 2119 | |
---|
| 2120 | if Tint != '-1': |
---|
| 2121 | tini = np.float(Tint.split('@')[0]) |
---|
| 2122 | tend = np.float(Tint.split('@')[1]) |
---|
| 2123 | else: |
---|
| 2124 | tini = -1. |
---|
| 2125 | tend = -1. |
---|
| 2126 | |
---|
| 2127 | vartimetrjv = [] |
---|
| 2128 | |
---|
| 2129 | print ' ' + fname |
---|
| 2130 | for trjfile in sims: |
---|
| 2131 | print ' ' + trjfile + ' ...' |
---|
| 2132 | if not os.path.isfile(trjfile): |
---|
| 2133 | print errormsg |
---|
| 2134 | print ' ' + fname + ": trajectory file: '" + trjfile + \ |
---|
| 2135 | "' does not exist !!" |
---|
| 2136 | quit(-1) |
---|
| 2137 | |
---|
| 2138 | trjobj = NetCDFFile(trjfile, 'r') |
---|
| 2139 | otim = trjobj.variables['time'] |
---|
| 2140 | if not trjobj.variables.has_key(statistics + 'box_' + variable): |
---|
| 2141 | print errormsg |
---|
| 2142 | print ' ' + fname + ": file '" + trjfile + "' does not have variable '"+\ |
---|
| 2143 | statistics + 'box_' + variable + "' !!" |
---|
| 2144 | quit(-1) |
---|
| 2145 | ovar = trjobj.variables[statistics + 'box_' + variable] |
---|
| 2146 | dimt = otim.shape[0] |
---|
| 2147 | |
---|
| 2148 | if trjfile == sims[0]: |
---|
| 2149 | gunits = ovar.getncattr('units') |
---|
| 2150 | lname = ovar.getncattr('long_name') |
---|
| 2151 | tunits = otim.getncattr('units') |
---|
| 2152 | |
---|
| 2153 | if tini != -1: |
---|
| 2154 | tiniid = -1 |
---|
| 2155 | tendid = -1 |
---|
| 2156 | for itv in range(dimt): |
---|
| 2157 | if otim[itv] <= tini and otim[itv+1] >= tini: tiniid = itv |
---|
| 2158 | if otim[itv] <= tend and otim[itv+1] >= tend: tendid = itv |
---|
| 2159 | |
---|
| 2160 | if tiniid == -1 or tendid == -1: |
---|
| 2161 | print errormsg |
---|
| 2162 | print ' ' + main + ' time interval ', tini,',',tend,' not found: ', \ |
---|
| 2163 | tendid, ',', tiniid, ' !!' |
---|
| 2164 | print ' data interval [',otim[0], otim[dimt-1],']' |
---|
| 2165 | quit(-1) |
---|
| 2166 | dimt = tendid - tiniid + 1 |
---|
| 2167 | |
---|
| 2168 | else: |
---|
| 2169 | dimt = otim.shape[0] |
---|
| 2170 | |
---|
| 2171 | valsv = np.zeros((2,dimt), dtype=np.float) |
---|
| 2172 | |
---|
| 2173 | if tini == -1: |
---|
| 2174 | valsv[1,:] = otim[:] |
---|
| 2175 | valsv[0,:] = ovar[:] |
---|
| 2176 | else: |
---|
| 2177 | valsv[1,:] = otim[tiniid:tendid+1] |
---|
| 2178 | valsv[0,:] = ovar[tiniid:tendid+1] |
---|
| 2179 | |
---|
| 2180 | vartimetrjv.append(valsv) |
---|
| 2181 | trjobj.close() |
---|
| 2182 | |
---|
| 2183 | drw.plot_TimeSeries(vartimetrjv, drw.units_lunits(gunits), tunits, \ |
---|
| 2184 | 'val_trajectories_' + statistics, gvarname, timetit, tkind, timefmt, title, \ |
---|
| 2185 | leglabels, locleg, gkind) |
---|
| 2186 | |
---|
| 2187 | def variable_values(values): |
---|
| 2188 | """ Function to give back values for a given variable |
---|
| 2189 | values= [varname] name of the variable |
---|
| 2190 | """ |
---|
| 2191 | |
---|
| 2192 | fname = 'variable_values' |
---|
| 2193 | |
---|
| 2194 | values = drw.variables_values(values) |
---|
| 2195 | |
---|
| 2196 | print fname,'values:',values |
---|
| 2197 | print fname,'variable_name:',values[0] |
---|
| 2198 | print fname,'standard_name:',values[1] |
---|
| 2199 | print fname,'min,max:',str(values[2]) + ',' + str(values[3]) |
---|
| 2200 | print fname,'long_name:',values[4] |
---|
| 2201 | print fname,'units:',values[5] |
---|
| 2202 | print fname,'shad_colors:',values[6] |
---|
| 2203 | print fname,'all_values:',drw.numVector_String(values,',') |
---|
| 2204 | |
---|
| 2205 | return |
---|
| 2206 | |
---|
| 2207 | ####### ###### ##### #### ### ## # |
---|
| 2208 | |
---|
| 2209 | ngraphics = "'" + drw.numVector_String(namegraphics, "', '") + "'" |
---|
| 2210 | |
---|
| 2211 | ### Options |
---|
| 2212 | ##string_operation="operation to make: " + '\n' + " out, output values -S inidim1,[inidim2,...]:enddim1,[enddim2,...]" |
---|
| 2213 | string_operation="""operation to make: |
---|
| 2214 | draw_topo_geogrid, draws topography from a WPS geo_em.d[nn].nc: -S [minTopo],[maxTopo]:[SW_lon],[SW_lat],[NE_lon],[NE_lat]:[title]:[graphic_kind]:[projection],[res_coastline] |
---|
| 2215 | draw_2D_shad_cont, draws two 2D fields, first with shading second with contour lines: -v [varns],[varnc] -S [vnamefs],[vnamefc],[dimxvn],[dimyvn],[colorbar],[ckind],[clabfmt],[sminv]:[smaxv],[sminc]:[smaxv]:[Nlev],[figt],[kindfig],[reverse] |
---|
| 2216 | [ckind]: |
---|
| 2217 | 'cmap': as it gets from colorbar |
---|
| 2218 | 'fixc,[colname]': fixed color [colname], all stright lines |
---|
| 2219 | 'fixsignc,[colname]': fixed color [colname], >0 stright, <0 dashed line |
---|
| 2220 | """ |
---|
| 2221 | |
---|
| 2222 | #print string_operation |
---|
| 2223 | |
---|
| 2224 | parser = OptionParser() |
---|
| 2225 | parser.add_option("-f", "--netCDF_file", dest="ncfile", |
---|
| 2226 | help="file to use", metavar="FILE") |
---|
| 2227 | parser.add_option("-o", "--operation", type='choice', dest="operation", |
---|
| 2228 | choices=namegraphics, |
---|
| 2229 | help="operation to make: " + ngraphics, metavar="OPER") |
---|
| 2230 | parser.add_option("-S", "--valueS", dest="values", |
---|
| 2231 | help="[WHEN APPLICABLE] values to use according to the operation", metavar="VALUES") |
---|
| 2232 | parser.add_option("-v", "--variable", dest="varname", |
---|
| 2233 | help="[WHEN APPLICABLE] variable to check", metavar="VAR") |
---|
| 2234 | |
---|
| 2235 | (opts, args) = parser.parse_args() |
---|
| 2236 | |
---|
| 2237 | ####### ####### |
---|
| 2238 | ## MAIN |
---|
| 2239 | ####### |
---|
| 2240 | |
---|
| 2241 | # Not checking file operation |
---|
| 2242 | Notcheckingfile = ['draw_2D_shad_cont', 'draw_2D_shad_cont_time', \ |
---|
| 2243 | 'draw_2D_shad_line', 'draw_2D_shad_line_time', 'draw_topo_geogrid_boxes', \ |
---|
| 2244 | 'draw_trajectories', 'draw_vals_trajectories', 'variable_values'] |
---|
| 2245 | |
---|
| 2246 | ####### ###### ##### #### ### ## # |
---|
| 2247 | errormsg='ERROR -- error -- ERROR -- error' |
---|
| 2248 | |
---|
| 2249 | varn=opts.varname |
---|
| 2250 | oper=opts.operation |
---|
| 2251 | |
---|
| 2252 | if opts.ncfile is not None and not os.path.isfile(opts.ncfile) and \ |
---|
| 2253 | not drw.searchInlist(Notcheckingfile, oper): |
---|
| 2254 | print errormsg |
---|
| 2255 | print ' ' + main + ': File ' + opts.ncfile + ' does not exist !!' |
---|
| 2256 | quit(-1) |
---|
| 2257 | |
---|
| 2258 | if oper == 'create_movie': |
---|
| 2259 | create_movie(opts.ncfile, opts.values, opts.varname) |
---|
| 2260 | elif oper == 'draw_2D_shad': |
---|
| 2261 | draw_2D_shad(opts.ncfile, opts.values, opts.varname) |
---|
| 2262 | elif oper == 'draw_2D_shad_time': |
---|
| 2263 | draw_2D_shad_time(opts.ncfile, opts.values, opts.varname) |
---|
| 2264 | elif oper == 'draw_2D_shad_cont': |
---|
| 2265 | draw_2D_shad_cont(opts.ncfile, opts.values, opts.varname) |
---|
| 2266 | elif oper == 'draw_2D_shad_cont_time': |
---|
| 2267 | draw_2D_shad_cont_time(opts.ncfile, opts.values, opts.varname) |
---|
| 2268 | elif oper == 'draw_2D_shad_line': |
---|
| 2269 | draw_2D_shad_line(opts.ncfile, opts.values, opts.varname) |
---|
| 2270 | elif oper == 'draw_2D_shad_line_time': |
---|
| 2271 | draw_2D_shad_line_time(opts.ncfile, opts.values, opts.varname) |
---|
| 2272 | elif oper == 'draw_Neighbourghood_evol': |
---|
| 2273 | draw_Neighbourghood_evol(opts.ncfile, opts.values, opts.varname) |
---|
| 2274 | elif oper == 'draw_lines': |
---|
| 2275 | draw_lines(opts.ncfile, opts.values, opts.varname) |
---|
| 2276 | elif oper == 'draw_timeSeries': |
---|
| 2277 | draw_timeSeries(opts.ncfile, opts.values, opts.varname) |
---|
| 2278 | elif oper == 'draw_topo_geogrid': |
---|
| 2279 | draw_topo_geogrid(opts.ncfile, opts.values) |
---|
| 2280 | elif oper == 'draw_topo_geogrid_boxes': |
---|
| 2281 | draw_topo_geogrid_boxes(opts.ncfile, opts.values) |
---|
| 2282 | elif oper == 'draw_trajectories': |
---|
| 2283 | draw_trajectories(opts.ncfile, opts.values, opts.varname) |
---|
| 2284 | elif oper == 'draw_vals_trajectories': |
---|
| 2285 | draw_vals_trajectories(opts.ncfile, opts.values, opts.varname) |
---|
| 2286 | elif oper == 'list_graphics': |
---|
| 2287 | # From: http://www.diveintopython.net/power_of_introspection/all_together.html |
---|
| 2288 | import drawing as myself |
---|
| 2289 | object = myself |
---|
| 2290 | for opern in namegraphics: |
---|
| 2291 | if opern != 'list_graphics': |
---|
| 2292 | print opern + '_______ ______ _____ ____ ___ __ _' |
---|
| 2293 | print getattr(object, opern).__doc__ |
---|
| 2294 | elif oper == 'variable_values': |
---|
| 2295 | variable_values(opts.values) |
---|
| 2296 | else: |
---|
| 2297 | print errormsg |
---|
| 2298 | print ' ' + main + ": the graphic '" + oper + "' is not ready !!" |
---|
| 2299 | print errormsg |
---|
| 2300 | quit() |
---|