1 | #!/usr/bin/env python |
---|
2 | |
---|
3 | ### T.Navarro |
---|
4 | |
---|
5 | ## |
---|
6 | # This routine reads sparse data, i.e. not binned in a lat/lon/vert/time grid and does plots |
---|
7 | # |
---|
8 | # There is nothing gnuplot could not be able to do for now, except command line and ease and all python plots possibilities, |
---|
9 | # but the goal is to provide a framework to compare with data from simulations on a constant grid. |
---|
10 | ### |
---|
11 | |
---|
12 | ########################################################################################### |
---|
13 | ########################################################################################### |
---|
14 | ### What is below relate to running the file as a command line executable (very convenient) |
---|
15 | if __name__ == "__main__": |
---|
16 | import sys |
---|
17 | from optparse import OptionParser ### to be replaced by argparse |
---|
18 | from os import system,path |
---|
19 | import numpy as np |
---|
20 | from sparse import sparse |
---|
21 | |
---|
22 | parser = OptionParser() |
---|
23 | |
---|
24 | ############################# |
---|
25 | ### Options |
---|
26 | parser.add_option('-f', '--file', action='store',dest='file', type="string", default=None, help='[NEEDED] filename, comma separated') |
---|
27 | parser.add_option('--ftype', action='store',dest='ftype', default=None, help=' Force data file type [\'sav\',\'txt\']') |
---|
28 | parser.add_option('-X', '--xvar', action='store',dest='xvar', type="string", default=None, help='[NEEDED] x variable: a name, a number, etc ...') |
---|
29 | parser.add_option('-Y', '--yvar', action='store',dest='yvar', type="string", default=None, help='y variable') |
---|
30 | parser.add_option('-Z', '--zvar', action='store',dest='zvar', type="string", default=None, help='z variable') |
---|
31 | parser.add_option('-C', '--cvar', action='store',dest='cvar', type="string", default=None, help='c variable for color') |
---|
32 | parser.add_option('--xmin', action='store',dest='xmin', type="float", default=None, help='min values for x var, comma separated') |
---|
33 | parser.add_option('--xmax', action='store',dest='xmax', type="float", default=None, help='max values for x var, comma separated') |
---|
34 | parser.add_option('--ymin', action='store',dest='ymin', type="float", default=None, help='min values for y var, comma separated') |
---|
35 | parser.add_option('--ymax', action='store',dest='ymax', type="float", default=None, help='max values for y var, comma separated') |
---|
36 | parser.add_option('--zmin', action='store',dest='zmin', type="float", default=None, help='min values for z var, comma separated') |
---|
37 | parser.add_option('--zmax', action='store',dest='zmax', type="float", default=None, help='max values for z var, comma separated') |
---|
38 | parser.add_option('--cmin', action='store',dest='cmin', type="string", default=None, help='min values for c var, comma separated') |
---|
39 | parser.add_option('--cmax', action='store',dest='cmax', type="string", default=None, help='max values for c var, comma separated') |
---|
40 | parser.add_option('-w','--with', action='append',dest='cond', type="string", default=None, help='conditions..') |
---|
41 | parser.add_option('--merge', action='store_true', dest='merge', default=False, help='merge datafiles in a single plot [False]') |
---|
42 | parser.add_option('-c','--color', action='store',dest='clb', type="string", default="def", help='change colormap (also: nobar,onebar)') |
---|
43 | parser.add_option('--trans', action='store',dest='trans', type="float", default=1., help='shaded plot transparency, 0 to 1 (=opaque) [1]') |
---|
44 | parser.add_option('-p', '--proj', action='store',dest='proj', type="string", default=None, help='projection') |
---|
45 | parser.add_option('--blat', action='store',dest='blat', type="int", default=None, help='reference lat (or bounding lat for stere) [computed]') |
---|
46 | parser.add_option('--blon', action='store',dest='blon', type="int", default=None, help='reference lon [computed]') |
---|
47 | parser.add_option('-b', '--back', action='store',dest='back', type="string", default=None, help='background image [None]') |
---|
48 | parser.add_option('-m', '--min', action='store',dest='vmin', type="float", default=None, help='bounding minimum value [min]') |
---|
49 | parser.add_option('-M', '--max', action='store',dest='vmax', type="float", default=None, help='bounding maximum value [max]') |
---|
50 | parser.add_option('--div', action='store',dest='ndiv', type="int", default=10, help='number of divisions in histogram [10]') |
---|
51 | parser.add_option('-S', '--save', action='store',dest='save', type="string", default="gui", help='save mode (gui,png,eps,svg,pdf,txt,html,avi) [gui]') |
---|
52 | |
---|
53 | |
---|
54 | ############################# |
---|
55 | ### Get options and variables |
---|
56 | (opt,args) = parser.parse_args() |
---|
57 | |
---|
58 | sparse(oplot=None,file=opt.file,ftype=opt.ftype,xvar=opt.xvar,yvar=opt.yvar,zvar=opt.zvar,cvar=opt.cvar,xmin=opt.xmin,xmax=opt.xmax,ymin=opt.ymin,ymax=opt.ymax,zmin=opt.zmin,zmax=opt.zmax,cmin=opt.cmin,cmax=opt.cmax,cond=opt.cond,merge=opt.merge,clb=opt.clb,trans=opt.trans,proj=opt.proj,blat=opt.blat,blon=opt.blon,back=opt.back,vmin=opt.vmin,vmax=opt.vmax,ndiv=opt.ndiv,save=opt.save) |
---|
59 | |
---|
60 | def sparse(oplot=None,file=None,ftype=None,xvar=None,yvar=None,zvar=None,cvar=None,xmin=None,xmax=None,ymin=None,ymax=None,zmin=None,zmax=None,cmin=None,cmax=None,cond=None,merge=False,clb="def",trans=1.,proj=None,blat=None,blon=None,back=None,vmin=None,vmax=None,ndiv=10,save="gui"): |
---|
61 | |
---|
62 | import sys |
---|
63 | from os import system,path |
---|
64 | from scipy.io.idl import readsav |
---|
65 | import numpy as np |
---|
66 | from myplot import separatenames, definesubplot, errormess, defcolorb, fmtvar, polarinterv, simplinterv, define_proj, readdata, makeplotres |
---|
67 | from mymath import min, max, writeascii |
---|
68 | import matplotlib as mpl |
---|
69 | from matplotlib.pyplot import subplot, figure, plot, scatter, colorbar, show, title, close, legend, xlabel, ylabel, axis, hist |
---|
70 | from matplotlib.cm import get_cmap |
---|
71 | from mpl_toolkits.mplot3d import Axes3D |
---|
72 | |
---|
73 | |
---|
74 | ############################# |
---|
75 | ### Load and check data |
---|
76 | |
---|
77 | ############################# |
---|
78 | ### Load and check data |
---|
79 | |
---|
80 | if file is None: |
---|
81 | print "You must specify at least a file to process with -f." |
---|
82 | exit() |
---|
83 | if xvar is None: |
---|
84 | print "You must specify at least a 1st field with -X." |
---|
85 | exit() |
---|
86 | if proj is not None and yvar is None: |
---|
87 | print "Why did you ask a projection with only one variable?" |
---|
88 | exit() |
---|
89 | |
---|
90 | filename=separatenames(file) |
---|
91 | |
---|
92 | #print 'conditions', cond |
---|
93 | |
---|
94 | print vmax |
---|
95 | print vmin |
---|
96 | |
---|
97 | if cond is None: nslices = 1 |
---|
98 | else: nslices = len(cond) |
---|
99 | numplot = nslices |
---|
100 | if merge is False: numplot = numplot*len(filename) |
---|
101 | |
---|
102 | print ' ' |
---|
103 | if merge is False: |
---|
104 | print nslices, 'condition(s) and', len(filename), 'file(s) without merging ---->', numplot, 'plot(s)' |
---|
105 | else: |
---|
106 | print nslices, 'condition(s) and', len(filename), 'file(s) with merging ---->', numplot, 'plot(s)' |
---|
107 | print ' ' |
---|
108 | |
---|
109 | all_type = [[]]*len(filename) |
---|
110 | all_x = [[]]*len(filename) |
---|
111 | all_y = [[]]*len(filename) |
---|
112 | all_z = [[]]*len(filename) |
---|
113 | all_c = [[]]*len(filename) |
---|
114 | all_data = [[]]*len(filename) |
---|
115 | #index = [[]]*len(filename) |
---|
116 | |
---|
117 | |
---|
118 | ############################## |
---|
119 | ################## READ DATA |
---|
120 | |
---|
121 | for k in range(len(filename)): |
---|
122 | |
---|
123 | ##### Find file type |
---|
124 | if ftype is None: |
---|
125 | if filename[k].find('.sav') is not -1: |
---|
126 | all_type[k] ='sav' |
---|
127 | print '.sav file detected for', filename[k],'!!' |
---|
128 | elif filename[k].find('.txt') is not -1: |
---|
129 | all_type[k] ='txt' |
---|
130 | #print '.txt file detected for', filename[k],'!!' |
---|
131 | else: |
---|
132 | all_type[k] ='txt' |
---|
133 | print 'no file type detected for', filename[k],'!!, default type is', all_type[k] |
---|
134 | |
---|
135 | if all_type[k] != all_type[0]: |
---|
136 | print 'Two different types were used: ', all_type[k], all_type[0] |
---|
137 | errormess('Not suported') |
---|
138 | else: |
---|
139 | all_type[k] = ftype |
---|
140 | |
---|
141 | ##### Read file |
---|
142 | if all_type[k] == 'sav': |
---|
143 | print 'reading .sav file', filename[k], '...' |
---|
144 | data = {} |
---|
145 | data = readsav(filename[k], idict=data, python_dict=False) #, verbose=True) |
---|
146 | all_data[k] = data |
---|
147 | elif all_type[k] == 'txt': |
---|
148 | print 'reading .text file', filename[k], '...' |
---|
149 | all_data[k] = np.loadtxt(filename[k]) |
---|
150 | |
---|
151 | all_x[k] = readdata(all_data,all_type[k],k,xvar) |
---|
152 | if yvar is not None: all_y[k] = readdata(all_data,all_type[k],k,yvar) |
---|
153 | else: all_y[k] = None |
---|
154 | if zvar is not None: all_z[k] = readdata(all_data,all_type[k],k,zvar) |
---|
155 | else: all_z[k] = None |
---|
156 | if cvar is not None: all_c[k] = readdata(all_data,all_type[k],k,cvar) |
---|
157 | else: all_c[k] = None |
---|
158 | |
---|
159 | if merge is True and k >=1 : |
---|
160 | all_x[0] = np.concatenate((all_x[0],all_x[k])) |
---|
161 | if yvar is not None: all_y[0] = np.concatenate((all_y[0],all_y[k])) |
---|
162 | if zvar is not None: all_z[0] = np.concatenate((all_z[0],all_z[k])) |
---|
163 | if cvar is not None: all_c[0] = np.concatenate((all_c[0],all_c[k])) |
---|
164 | |
---|
165 | |
---|
166 | |
---|
167 | ############################## |
---|
168 | ################## PLOT DATA |
---|
169 | |
---|
170 | ## Open a figure and set subplots |
---|
171 | if oplot is not None: fig=oplot |
---|
172 | else: fig = figure() |
---|
173 | subv,subh = definesubplot(numplot,fig) |
---|
174 | palette = get_cmap(name=clb) |
---|
175 | |
---|
176 | |
---|
177 | for nplot in range(numplot): |
---|
178 | |
---|
179 | print ' ' |
---|
180 | |
---|
181 | ###### Find which data and which file for plot nplot |
---|
182 | if merge is False: |
---|
183 | index_s = ((nplot)//len(filename))%nslices |
---|
184 | index_f = ((nplot-1)//nslices)%len(filename) |
---|
185 | else: |
---|
186 | index_s = ((nplot)//1)%nslices |
---|
187 | index_f = ((nplot-1)//nslices)%1 |
---|
188 | #print 'nplot,numplot,index_f,index_s', nplot,numplot,index_f,index_s |
---|
189 | |
---|
190 | |
---|
191 | ###### Select point under conditions defined in -w option |
---|
192 | index = np.isfinite(all_x[index_f]) |
---|
193 | if cond is not None: |
---|
194 | zecond = separatenames(cond[index_s]) |
---|
195 | #print 'hello', cond[index_s], zecond |
---|
196 | for i in range(len(zecond)): |
---|
197 | zecondi = zecond[i] |
---|
198 | if zecondi.find('>') != -1: |
---|
199 | variable = zecondi[0:zecondi.find('>')] |
---|
200 | value = float(zecondi[zecondi.find('>')+1:len(zecondi)]) |
---|
201 | if merge is True: # ultra moche de reconcatener a chaque fois, mais bon on s'en fout c'est du python |
---|
202 | zedata = [] |
---|
203 | for k in range(len(filename)): |
---|
204 | zedata = np.concatenate((zedata,readdata(all_data,all_type[k],k,variable))) |
---|
205 | else : zedata = readdata(all_data,all_type[k],index_f,variable) |
---|
206 | #print index.shape,zedata.shape, value |
---|
207 | index = index*(zedata>value) |
---|
208 | print 'All points such as', variable, '>', value, 'have been selected' |
---|
209 | elif zecondi.find('<') != -1: |
---|
210 | variable = zecondi[0:zecondi.find('<')] |
---|
211 | value = float(zecondi[zecondi.find('<')+1:len(zecondi)]) |
---|
212 | if merge is True: |
---|
213 | zedata = [] |
---|
214 | for k in range(len(filename)): |
---|
215 | zedata = np.concatenate((zedata,readdata(all_data,all_type[k],k,variable))) |
---|
216 | else : zedata = readdata(all_data,all_type[k],index_f,variable) |
---|
217 | #print index.shape,zedata.shape, value |
---|
218 | index = index*(zedata<value) |
---|
219 | print 'All points such as', variable, '<', value, 'have been selected' |
---|
220 | else: |
---|
221 | print '' |
---|
222 | print 'I do not understand that condition :', zecondi |
---|
223 | errormess('') |
---|
224 | |
---|
225 | |
---|
226 | if np.sum(index) == 0: |
---|
227 | print '*********** WARNING ***********' |
---|
228 | print '*********** NO DATA ***********' |
---|
229 | errormess('') |
---|
230 | else: |
---|
231 | print np.sum(index),'points have been selected among', len(all_x[index_f]), \ |
---|
232 | 'that is to say %2.2f' %(float(100*np.sum(index))/float(len(all_x[index_f]))), '%.' |
---|
233 | |
---|
234 | |
---|
235 | #print 'numplot', numplot |
---|
236 | changesubplot = (numplot > 1) #and (len(what_I_plot.shape) != 1) |
---|
237 | if changesubplot: subplot(subv,subh,nplot+1) |
---|
238 | #print 'subv,subh,nplot', subv,subh,nplot |
---|
239 | |
---|
240 | ###### Projection |
---|
241 | if proj is not None: # Nota : NEVER TRY TO DO A MESHGRID ON SPARSE DATA. WAY TOO MUCH POINTS. |
---|
242 | wlon = [min(all_x[index_f][index]),max(all_x[index_f][index])] |
---|
243 | wlat = [min(all_y[index_f][index]),max(all_y[index_f][index])] |
---|
244 | m = define_proj(proj,wlon,wlat,back=back,blat=blat,blon=blon) |
---|
245 | x, y = m(all_x[index_f][index],all_y[index_f][index]) |
---|
246 | else: |
---|
247 | x = all_x[index_f][index] |
---|
248 | if yvar is not None: y = all_y[index_f][index] |
---|
249 | |
---|
250 | ###### Plot: 1D histogram |
---|
251 | if yvar is None: |
---|
252 | hist(x,bins=ndiv,histtype='step',linewidth=2) |
---|
253 | if save == 'txt': |
---|
254 | print 'saving file profile'+str(nplot+1)+'.txt' |
---|
255 | writeascii(np.transpose(x),'profile'+str(nplot+1)+'.txt') |
---|
256 | ###### Plot: 2D cloud |
---|
257 | elif zvar is None and cvar is None : |
---|
258 | plot(x,y,'.b') |
---|
259 | if save == 'txt': |
---|
260 | print 'saving file profile'+str(nplot+1)+'.txt' |
---|
261 | writeascii(np.transpose(np.array([x,y])),'profile'+str(nplot+1)+'.txt') |
---|
262 | ###### Plot: 2D cloud with color |
---|
263 | elif zvar is None and cvar is not None : |
---|
264 | if save == 'txt': |
---|
265 | print 'saving file profile'+str(nplot+1)+'.txt' |
---|
266 | writeascii(np.transpose(np.array([x,y,all_c[index_f][index]])),'profile'+str(nplot+1)+'.txt') |
---|
267 | scatter(x,y,c=all_c[index_f][index],\ |
---|
268 | marker='o', edgecolor='None',cmap = palette, alpha=trans, vmin = vmin,vmax = vmax) |
---|
269 | ###### Plot: 3D cloud |
---|
270 | elif zvar is not None and cvar is None : |
---|
271 | ax = fig.add_subplot(subv,subh,nplot+1, projection='3d') |
---|
272 | ax.plot(x,y,all_z[index_f][index],'.b') |
---|
273 | ###### Plot: 3D cloud with color |
---|
274 | else : |
---|
275 | ax = fig.add_subplot(subv,subh,nplot+1, projection='3d') |
---|
276 | ax.scatter(x,y,all_z[index_f][index],c=all_c[index_f][index],\ |
---|
277 | marker='o', edgecolor='None', cmap = palette, alpha=trans,vmin = vmin,vmax = vmax) |
---|
278 | |
---|
279 | ###### Colorbar: http://www.scipy.org/Cookbook/Matplotlib/Show_colormaps?action=AttachFile&do=get&target=colormaps3.png |
---|
280 | if clb != 'nobar' and all_c[index_f] is not None and all_z[index_f] is None: # pourquoi la colorbar marche pas en 3d? |
---|
281 | colorbar( fraction=0.05,pad=0.03,format='%0.2f',\ |
---|
282 | extend='neither',spacing='proportional' ) |
---|
283 | |
---|
284 | ###### Plot limits (options --xmin, --xmax, etc ..) |
---|
285 | if proj is None: |
---|
286 | xlabel(xvar) |
---|
287 | ylabel(yvar) |
---|
288 | if zvar is None : |
---|
289 | if xmin is not None: mpl.pyplot.xlim(xmin=xmin) |
---|
290 | else: mpl.pyplot.xlim(xmin=min(all_x[index_f][index])) |
---|
291 | if xmax is not None: mpl.pyplot.xlim(xmax=xmax) |
---|
292 | else: mpl.pyplot.xlim(xmax=max(all_x[index_f][index])) |
---|
293 | if yvar is not None: |
---|
294 | if ymin is not None: mpl.pyplot.ylim(ymin=ymin) |
---|
295 | else: mpl.pyplot.ylim(ymin=min(all_y[index_f][index])) |
---|
296 | if ymax is not None: mpl.pyplot.ylim(ymax=ymax) |
---|
297 | else: mpl.pyplot.ylim(ymax=max(all_y[index_f][index])) |
---|
298 | |
---|
299 | |
---|
300 | if cond is not None: |
---|
301 | title(cond[index_s]) |
---|
302 | else: |
---|
303 | title('all point selected') |
---|
304 | |
---|
305 | zeplot = "output" |
---|
306 | if save in ['png','eps','svg','pdf']: makeplotres(zeplot,ext=save)#,res=resolution,pad_inches_value=pad_inches_value,disp=False,ext=save) |
---|
307 | elif save == 'gui': show() |
---|
308 | elif save == 'return': return mpl.pyplot |
---|
309 | else: print "INFO: save mode not supported. using gui instead." ; show() |
---|
310 | print '' |
---|
311 | |
---|
312 | command = "" |
---|
313 | for arg in sys.argv: command = command + arg + ' ' |
---|
314 | name = zeplot |
---|
315 | f = open(name+'.sh', 'w') |
---|
316 | f.write(command) |
---|
317 | |
---|
318 | |
---|
319 | |
---|
320 | |
---|
321 | |
---|