source: anuga_work/development/Hinwood_2008/plot.py @ 5681

Last change on this file since 5681 was 5681, checked in by duncan, 16 years ago

moving Honwood files around

File size: 13.1 KB
Line 
1
2"""
3Plot up files from the Hinwood project.
4"""
5from os import sep
6import project
7from time import localtime, strftime
8
9def plot_compare_csv(location_sim, file_sim, location_exp, file_exp,
10                     y_label,
11                     save_as=None,
12                     plot_title="",
13                     x_label='Time (s)',
14                     legend_sim='ANUGA simulation',
15                     legend_exp='Measured flume result',
16                     is_interactive=False,
17                     x_axis=None,
18                     y_axis=None):
19    """
20    """
21    from pylab import ion, plot, xlabel, ylabel, close, legend, \
22         savefig, title, axis, setp
23    from anuga.shallow_water.data_manager import csv2dict
24
25
26
27    # Load in the csv files and convert info from strings to floats
28    simulation, _ = csv2dict(file_sim)
29    experiment, _ = csv2dict(file_exp)
30    time_sim = [float(x) for x in simulation['time']]
31    quantity_sim = [float(x) for x in simulation[location_sim]]
32    #print "quantity_sim", quantity_sim
33    time_exp = [float(x) for x in experiment['Time']]
34    quantity_exp = [float(x) for x in experiment[location_exp]]
35
36    if is_interactive:
37        ion()
38   
39    l_sim, l_exp = plot(time_sim, quantity_sim, time_exp, quantity_exp)
40    setp(l_sim, color='r')
41    setp(l_exp, color='b')
42
43    # Add axis stuff and legend
44    xlabel(x_label)
45    ylabel(y_label)
46    # The order defines the label
47    #legend((legend_exp, legend_sim),'upper left')
48    legend((legend_sim, legend_exp),'upper left')
49    title(plot_title)
50    if x_axis is not None:
51        axis(xmin=x_axis[0], xmax=x_axis[1])
52   
53    if y_axis is not None:
54        axis(ymin=y_axis[0], ymax=y_axis[1])
55       
56    if is_interactive:
57        # Wait for enter pressed
58        raw_input()
59
60    if save_as is not None:
61        savefig(save_as)
62   
63    #Need to close this plot
64    close()
65   
66def Hinwood_files_locations(run_data, outputdir_tag, plot_type,
67                            quantity = "depth",
68                            y_location_tag=':0.0'):
69    """
70    run_data is a dictionary of data describing a Hinwood experiment
71    outputdir_tag a string at the end of an output dir; '_good_tri_area_0.01_A'
72    plot_type the file extension of the plot, eg '.pdf'
73    """
74    id = run_data['scenario_id']
75    outputdir_name = id + outputdir_tag
76    pro_instance = project.Project(['data','flumes','Hinwood_2008'],
77                                   outputdir_name=outputdir_name)
78   
79   
80    file_sim = pro_instance.outputdir + quantity + "_" + id + ".csv"
81    #print "file_exp",file_exp
82    file_exp = pro_instance.raw_data_dir + sep + id + 'pressfilt_exp_' \
83               + quantity + '.csv'
84    #print "file_sim", file_sim
85    location_sims = []
86    location_exps = []
87    save_as_list = []
88    for gauge_x in run_data['gauge_x']:
89        gauge_x = str(gauge_x)
90        location_sims.append(gauge_x + y_location_tag)
91        location_exps.append(gauge_x)
92        save_as_list.append(pro_instance.plots_dir + sep + \
93                            outputdir_name + "_" + quantity + "_" + \
94                            gauge_x + plot_type)
95    return file_exp, file_sim, location_sims, location_exps, outputdir_name, \
96           save_as_list
97
98def plot_exp_sim_comparision(scenarios, outputdir_tag,
99                             quantity = "stage",is_interactive=False,
100         y_location_tag=':0.0'):
101    plot_type = ".pdf"
102    plot_files = {} # Index scenario (T1R3) value, list of files
103    for run_data in scenarios:
104       
105        temp = Hinwood_files_locations(run_data, outputdir_tag,
106                                       plot_type, quantity,
107                                       y_location_tag=y_location_tag)
108                                   
109        file_exp, file_sim, location_sims, location_exps, outputdir_name, \
110                  save_as_list = temp
111        plot_files[run_data['scenario_id']] = save_as_list
112        print "file_exp",file_exp
113        print "run_data['scenario_id']", run_data['scenario_id']
114        #location_sims = [location_sims[0]]
115        #location_exps = [location_exps[0]]
116        #save_as_list = [save_as_list[0]]
117        for loc_sim, loc_exp, save_as, gauge, gauge_elev in \
118                map(None, location_sims,
119                    location_exps,
120                    save_as_list,
121                    run_data['gauge_x'],
122                    run_data['gauge_bed_elevation']):
123            time_date = strftime('plot date: %d/%m/%Y Time: %H:%M:%S',
124                                      localtime())
125            plot_title = "Scenario: " + outputdir_name + "\n" + \
126                         "X Gauge (m):" + str(gauge) + "    " + time_date
127
128            # When the shore gauges out of the are used, don't
129            # use the standard y axis scale
130            x_axis = run_data['wave_times']
131            if gauge < run_data['axis_maximum_x']:
132                y_axis = [run_data['axis'][2],run_data['axis'][3]]
133       
134            else:
135                y_axis = [run_data['axis'][2] + gauge_elev,
136                          run_data['axis'][3] + gauge_elev]
137            print "Doing ", plot_title
138            plot_compare_csv(location_sim=loc_sim,
139                             file_sim=file_sim,
140                             location_exp=loc_exp,
141                             file_exp=file_exp,
142                             plot_title=plot_title,
143                             y_label='Water '+ quantity +' (m)',
144                             is_interactive=is_interactive,
145                             save_as=save_as,
146                             x_axis=x_axis,
147                             y_axis=y_axis)
148            if is_interactive is True:
149                break
150    return plot_files
151
152
153def auto_plot_compare_csv_subfigures(scenarios, outputdir_tag,
154                                     to_publish_indexes,
155                                     quantity = "stage",is_interactive=False,
156                                     y_location_tag=':0.0',
157                                     add_run_info=False):
158    """
159    Used by validation graphs to produce the final figures
160    """
161    plot_type = ".pdf"
162    save_as_list = []
163    for run_data in scenarios:
164       
165        id = run_data['scenario_id']
166        if not to_publish_indexes.has_key(id):
167            continue
168           
169        outputdir_name = id + outputdir_tag
170        pro_instance = project.Project(['data','flumes','Hinwood_2008'],
171                                       outputdir_name=outputdir_name)
172   
173   
174        file_sim = pro_instance.outputdir + quantity + "_" + id + ".csv"
175        file_exp = pro_instance.raw_data_dir + sep + id + 'pressfilt_exp_' \
176                   + quantity + '.csv'
177
178        if add_run_info is True:
179            plot_title = "Title will not be in final draft \n" + \
180                         id + outputdir_tag
181        else:
182            plot_title = ""
183        #save_as = pro_instance.plots_dir + sep + \
184        #          outputdir_name + "_" + quantity + "_" + plot_type
185        save_as = pro_instance.plots_dir + sep + 'S' + id[1:2] + \
186                  '-' + quantity + '-compare' + plot_type
187        save_as_list.append(save_as)
188
189               
190        plot_compare_csv_subfigures(file_sim,
191                                    file_exp,
192                                    to_publish_indexes[id],
193                                    run_data,
194                                    plot_title=plot_title,
195                                    save_as=save_as,
196                                    y_label='Water '+ quantity +' (m)',
197                                    is_interactive=is_interactive,
198                                    y_location_tag=y_location_tag)
199    return save_as_list
200
201def plot_compare_csv_subfigures(file_sim,
202                                file_exp,
203                                gauge_indexs,
204                                run_data,
205                                save_as=None,
206                                plot_title="",
207                                x_label='Time (s)',
208                                y_label=None,
209                                legend_sim='ANUGA simulation',
210                                legend_exp='Measured flume result',
211                                is_interactive=False,
212                                x_axis=None,
213                                y_axis=None,
214                                y_location_tag=':0.0'):
215    """
216    Used by validation graphs to produce the final figures
217    """
218    from pylab import ion, plot, xlabel, ylabel, close, legend, \
219         savefig, title, axis, setp, subplot, grid, figlegend, gca, \
220         text, sqrt
221    import pylab
222   
223    from anuga.shallow_water.data_manager import csv2dict
224   
225    print "file_sim", file_sim
226
227    if False:
228        fig_width_pt = 246.0  # Get this from LaTeX using \showthe\columnwidth
229        inches_per_pt = 1.0/72.27               # Convert pt to inch
230        golden_mean = (sqrt(5)-1.0)/2.0         # Aesthetic ratio
231        fig_width = fig_width_pt*inches_per_pt  # width in inches
232        fig_height = fig_width*golden_mean      # height in inches
233        fig_size =  [fig_width,fig_height]
234        print "fig_size", fig_size
235        params = {'backend': 'ps',
236                  'axes.labelsize': 10,
237                  'text.fontsize': 10,
238                  'legend.fontsize': 10,
239                  'xtick.labelsize': 8,
240                  'ytick.labelsize': 8,
241                  'text.usetex': True,
242                  'figure.figsize': [3.4, 2.1]}
243    if False:
244        params = {'backend': 'ps',
245                  'axes.labelsize': 10,
246                  'text.fontsize': 10,
247                  'legend.fontsize': 10,
248                  'xtick.labelsize': 8,
249                  'ytick.labelsize': 8,
250                  'figure.figsize': [3.4, 3.1]}
251        pylab.rcParams.update(params)
252
253
254   
255    # Load in the csv files and convert info from strings to floats
256    simulation, _ = csv2dict(file_sim)
257    experiment, _ = csv2dict(file_exp)
258    time_sim = [float(x) for x in simulation['time']]
259    time_exp = [float(x) for x in experiment['Time']]
260
261   
262    if is_interactive:
263        ion()
264       
265    for position, i in enumerate(gauge_indexs):
266        gauge_x = run_data['gauge_x'][i]
267       
268        grid_position = (len(gauge_indexs))*100 + 10 + position +1
269        subplot(grid_position)
270        location_sim = str(gauge_x) + y_location_tag
271        location_exp = str(gauge_x)
272       
273        quantity_sim = [float(x) for x in simulation[location_sim]]
274        quantity_exp = [float(x) for x in experiment[location_exp]]
275
276   
277        l_sim, l_exp = plot(time_sim, quantity_sim, time_exp, quantity_exp)
278        setp(l_sim, color='k', linestyle='--')
279        setp(l_exp, color='k')
280        grid(True)
281       
282        # When the shore gauges out of the are used, don't
283        # use the standard y axis scale
284        gauge_elev = run_data['gauge_bed_elevation'][i]
285        x_axis = run_data['wave_times']
286        if gauge_x < run_data['axis_maximum_x']:
287            y_axis = [run_data['axis'][2],run_data['axis'][3]]   
288            setp(gca(), yticks=[-0.04,-0.02, 0.0, 0.02, 0.04])
289            y_text = -0.035
290        else:
291            y_axis = [run_data['axis'][2] + gauge_elev,
292                      run_data['axis'][3] + gauge_elev]   
293            setp(gca(), yticks=[0.02, 0.04, 0.06, 0.08, 0.1])
294            y_text = 0.025
295        text(x_axis[0]+5,y_text, str(round(gauge_x,1)) + " m gauge")   
296        if position == 0:
297            title(plot_title)
298        if position == 1:
299            ylabel(y_label)
300        #legend((legend_sim, legend_exp),'lower center')
301        if y_axis is not None:
302            axis(ymin=y_axis[0]-0.001, ymax=y_axis[1]-0.001)
303            # the fudge -0.001 is to reduce the ytick's
304
305        # need to drop the axis name being printed
306        if x_axis is not None:
307            axis(xmin=x_axis[0], xmax=x_axis[1])
308
309        # Only do tick marks on the final graph
310        if not position == 2:
311            setp(gca(), xticklabels=[])
312
313    # Add axis stuff and legend
314    xlabel(x_label)
315       
316    #ylabel(y_label)
317    # The order defines the label
318    #legend((legend_exp, legend_sim),'upper left')
319
320    # I couldn't get this working
321    #figlegend((l_sim, l_exp),
322     #         (legend_sim, legend_exp),
323      #        'lower center', axespad = -0.05)
324   
325    if is_interactive:
326        # Wait for enter pressed
327        raw_input()
328
329    if save_as is not None:
330        savefig(save_as)
331   
332    #Need to close this plot
333    close()
334   
335#-------------------------------------------------------------
336if __name__ == "__main__":
337    """ Plot the stage graph for the ANUGA validation papar
338    """
339    from scenarios import scenarios
340    outputdir_tag = "_good_tri_area_0.01_limiterE"
341    outputdir_tag = "_nolmts_wdth_0.1_z_0.012_ys_0.01_mta_0.01_H"
342    outputdir_tag = "_good_lmts_wdth_0.1_z_0.0_ys_0.01_mta_0.01_F"
343    outputdir_tag = "_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.01_G"
344    #outputdir_tag = "_nolmts_wdth_0.01_z_0.012_ys_0.01_mta_1e-05_G"
345    #outputdir_tag = "_test_C"
346    #scenarios = scenarios[1:]
347    #scenarios = [scenarios[7]]
348    is_interactive = False
349    #is_interactive = True
350    plot_exp_sim_comparision(scenarios, outputdir_tag,
351                             is_interactive=is_interactive,
352                             y_location_tag=':0.0')
353
Note: See TracBrowser for help on using the repository browser.