Changeset 5670


Ignore:
Timestamp:
Aug 20, 2008, 2:56:51 PM (15 years ago)
Author:
duncan
Message:

Current Hinwood - working on fig's for paper

Location:
anuga_work/development/Hinwood_2008
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • anuga_work/development/Hinwood_2008/calc_norm.py

    r5664 r5670  
    181181        axvspan(break_x-0.001,break_x+0.001, facecolor='g')
    182182       
    183     figlegend(lines, outputdir_tags,'upper left')
     183    legend(lines, outputdir_tags,'upper left')
    184184   
    185185    if is_interactive:
    … …  
    197197def auto_plot_test_rmsd(scenarios, outputdir_tag, quantity="stage",
    198198              save_as=True,
    199               is_interactive=False, max_rmsd=1.0):
     199              is_interactive=False, max_rmsd=0.010):
     200    """
     201    Produce the RMSD graphs for the anuga article.
     202
     203    to do:
     204    labels
     205    """
    200206   
    201207    tests = []
     208    file_names = []
    202209   
    203210    for i in range(0,len(scenarios),2):
    … …  
    205212
    206213    for test in tests:
    207         plot_test_rmsd(test, outputdir_tag, quantity,
     214        file_name = plot_test_rmsd(test, outputdir_tag, quantity,
    208215              save_as=save_as,
    209216              is_interactive=is_interactive, max_rmsd=max_rmsd)
     217        file_names.append(file_name)
     218    return file_names
    210219       
    211220def plot_test_rmsd(test, outputdir_tag, quantity,
    … …  
    238247
    239248    lines = []
     249    legend_name = []
     250    i = 0
    240251    for run in test:
    241            
     252        i += 1
     253        legend_name.append(run['scenario_id'][:2]+"T"+str(i))   
    242254        outputdir_name = run['scenario_id'] + outputdir_tag
    243255        pro_instance = project.Project(['data','flumes','Hinwood_2008'],
    … …  
    250262        locations = [float(x) for x in simulation['x location']]
    251263        rmsd_list = [float(x) for x in simulation['rmsd']]
    252         lines.append(plot(locations, rmsd_list))
     264        lines.append(plot(locations, rmsd_list, '-o'))
     265        #plot(locations, rmsd_list, 'o')
    253266       
    254267        for break_x in run['break_xs']:
    … …  
    258271        #print "setting axis"
    259272        axis(ymin=0, ymax=max_rmsd)
    260        
    261     figlegend(lines, outputdir_tags,'upper left')
     273        
     274    legend(lines, legend_name,'upper left')
    262275   
    263276    if is_interactive:
    … …  
    266279    save_as = pro_instance.plots_dir + sep + \
    267280              id + "_rmsd" + plot_type
     281    print "save_as", save_as
    268282    if save_as is not None:
    269283        savefig(save_as)
    270284   
    271285    close()
     286
     287    return save_as
    272288
    273289# Return a bunch of lists
    … …  
    377393   
    378394    outputdir_tags = []
    379     #outputdir_tags.append("_nolmts_wdth_0.1_z_0.012_ys_0.01_mta_0.01_H")
    380     outputdir_tags.append("_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.01_G")
    381     #outputdir_tag = "_good_tri_area_0.01_D"
    382     outputdir_tags.append("_nolmts_wdth_0.1_z_0.012_ys_0.01_mta_0.01_G")
    383     outputdir_tags.append("_lmts_wdth_0.1_z_0.012_ys_0.01_mta_0.01_G")
    384     outputdir_tags.append("_nolmts_wdth_0.1_z_0.012_ys_0.05_mta_1e-05_G")
    385     outputdir_tags.append("_nolmts_wdth_0.1_z_0.012_ys_0.01_mta_0.001_G")
    386     outputdir_tags.append("_nolmts_wdth_0.01_z_0.012_ys_0.01_mta_1e-05_G")
    387     outputdir_tags = []
    388     outputdir_tags.append("_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.01_G")
    389    
    390     #outputdir_tag = "_good_lmts_wdth_1.0_z_0.012_ys_0.01_mta_0.001_H"
     395    outputdir_tags.append("_lmts_wdth_0.1_z_0.0_ys_0.01_mta_0.01_I")
     396    outputdir_tags.append("_no_velocity_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.01_I")
     397    outputdir_tags.append("_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.0001_I")
     398    outputdir_tags.append("_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.001_I")
     399    outputdir_tags.append("_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.01_I")
     400    outputdir_tags.append("_nolmts_wdth_0.1_z_0.012_ys_0.01_mta_0.01_I")
     401    outputdir_tags.append("_nolmts_wdth_1.0_z_0.0_ys_0.01_mta_0.01_I")
     402   
    391403    #outputdir_tag = "_test_limiterC"
    392404    #scenarios = [scenarios[0]] # !!!!!!!!!!!!!!!!!!!!!!
    393405    #scenarios = scenarios[4:] # !!!!!!!!!!!!!!!!!!!!!!
    394     if False:
     406    calc_norms = True
     407    calc_norms = False
     408    if calc_norms:
    395409        for outputdir_tag in outputdir_tags:
    396410            auto_rrms(outputdir_tag, scenarios, "stage", y_location_tag=':0.0')
    … …  
    400414    #compare_different_settings(outputdir_tags, scenarios, "stage")
    401415
    402     outputdir_tag = "_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.01_G"
     416    outputdir_tag = "_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.01_I"
    403417    auto_plot_test_rmsd(scenarios, outputdir_tag)

  • anuga_work/development/Hinwood_2008/plot.py

    r5664 r5670  
    9696           save_as_list
    9797
    98 def plot(scenarios, outputdir_tag, quantity = "stage",is_interactive=False,
     98def plot_exp_sim_comparision(scenarios, outputdir_tag,
     99                             quantity = "stage",is_interactive=False,
    99100         y_location_tag=':0.0'):
    100101    plot_type = ".pdf"
    101    
     102    plot_files = {} # Index scenario (T1R3) value, list of files
    102103    for run_data in scenarios:
    103104       
    … …  
    108109        file_exp, file_sim, location_sims, location_exps, outputdir_name, \
    109110                  save_as_list = temp
     111        plot_files[run_data['scenario_id']] = save_as_list
    110112        print "file_exp",file_exp
    111113        print "run_data['scenario_id']", run_data['scenario_id']
    … …  
    146148            if is_interactive is True:
    147149                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    plot_type = ".pdf"
     159    save_as_list = []
     160    for run_data in scenarios:
     161       
     162        id = run_data['scenario_id']
     163        if not to_publish_indexes.has_key(id):
     164            continue
     165           
     166        outputdir_name = id + outputdir_tag
     167        pro_instance = project.Project(['data','flumes','Hinwood_2008'],
     168                                       outputdir_name=outputdir_name)
     169   
     170   
     171        file_sim = pro_instance.outputdir + quantity + "_" + id + ".csv"
     172        file_exp = pro_instance.raw_data_dir + sep + id + 'pressfilt_exp_' \
     173                   + quantity + '.csv'
     174
     175        if add_run_info is True:
     176            plot_title = "Title will not be in final draft \n" + \
     177                         id + outputdir_tag
     178        else:
     179            plot_title = ""
     180        #save_as = pro_instance.plots_dir + sep + \
     181        #          outputdir_name + "_" + quantity + "_" + plot_type
     182        save_as = pro_instance.plots_dir + sep + id + quantity + \
     183                  '_compare' + plot_type
     184        save_as_list.append(save_as)
     185
     186               
     187        plot_compare_csv_subfigures(file_sim,
     188                                    file_exp,
     189                                    to_publish_indexes[id],
     190                                    run_data,
     191                                    plot_title=plot_title,
     192                                    save_as=save_as,
     193                                    y_label='Water '+ quantity +' (m)',
     194                                    is_interactive=is_interactive,
     195                                    y_location_tag=y_location_tag)
     196    return save_as_list
     197
     198def plot_compare_csv_subfigures(file_sim,
     199                                file_exp,
     200                                gauge_indexs,
     201                                run_data,
     202                                save_as=None,
     203                                plot_title="",
     204                                x_label='Time (s)',
     205                                y_label=None,
     206                                legend_sim='ANUGA simulation',
     207                                legend_exp='Measured flume result',
     208                                is_interactive=False,
     209                                x_axis=None,
     210                                y_axis=None,
     211                                y_location_tag=':0.0'):
     212    """
     213    """
     214    from pylab import ion, plot, xlabel, ylabel, close, legend, \
     215         savefig, title, axis, setp, subplot, grid, figlegend, gca, \
     216         text
     217   
     218    from anuga.shallow_water.data_manager import csv2dict
     219   
     220    print "file_sim", file_sim
     221    # Load in the csv files and convert info from strings to floats
     222    simulation, _ = csv2dict(file_sim)
     223    experiment, _ = csv2dict(file_exp)
     224    time_sim = [float(x) for x in simulation['time']]
     225    time_exp = [float(x) for x in experiment['Time']]
     226
     227   
     228    if is_interactive:
     229        ion()
     230       
     231    for position, i in enumerate(gauge_indexs):
     232        gauge_x = run_data['gauge_x'][i]
     233        grid_position = (len(gauge_indexs)+1)*100 + 10 + position +1
     234        subplot(grid_position)
     235        location_sim = str(gauge_x) + y_location_tag
     236        location_exp = str(gauge_x)
     237       
     238        quantity_sim = [float(x) for x in simulation[location_sim]]
     239        quantity_exp = [float(x) for x in experiment[location_exp]]
     240
     241   
     242        l_sim, l_exp = plot(time_sim, quantity_sim, time_exp, quantity_exp)
     243        setp(l_sim, color='r', linestyle='--')
     244        setp(l_exp, color='b')
     245        grid(True)
     246       
     247        # When the shore gauges out of the are used, don't
     248        # use the standard y axis scale
     249        gauge_elev = run_data['gauge_bed_elevation'][i]
     250        x_axis = run_data['wave_times']
     251        if gauge_x < run_data['axis_maximum_x']:
     252            y_axis = [run_data['axis'][2],run_data['axis'][3]]   
     253            setp(gca(), yticks=[-0.04,-0.02, 0.0, 0.02, 0.04])
     254            y_text = -0.035
     255        else:
     256            y_axis = [run_data['axis'][2] + gauge_elev,
     257                      run_data['axis'][3] + gauge_elev]   
     258            setp(gca(), yticks=[0.02, 0.04, 0.06, 0.08, 0.1])
     259            y_text = 0.025
     260        text(x_axis[0]+5,y_text, str(round(gauge_x,1)) + " m gauge")   
     261        if position == 0:
     262            title(plot_title)
     263        if position == 1:
     264            ylabel(y_label)
     265        #legend((legend_sim, legend_exp),'lower center')
     266        if y_axis is not None:
     267            axis(ymin=y_axis[0]-0.001, ymax=y_axis[1]-0.001)
     268            # the fudge -0.001 is to reduce the ytick's
     269
     270        # need to drop the axis name being printed
     271        if x_axis is not None:
     272            axis(xmin=x_axis[0], xmax=x_axis[1])
     273
     274        # Only do tick marks on the final graph
     275        if not position == 2:
     276            setp(gca(), xticklabels=[])
     277
     278    # Add axis stuff and legend
     279    xlabel(x_label)
     280       
     281    #ylabel(y_label)
     282    # The order defines the label
     283    #legend((legend_exp, legend_sim),'upper left')
     284    figlegend((l_sim, l_exp),
     285              (legend_sim, legend_exp),
     286              'lower center')
     287   
     288    if is_interactive:
     289        # Wait for enter pressed
     290        raw_input()
     291
     292    if save_as is not None:
     293        savefig(save_as)
     294   
     295    #Need to close this plot
     296    close()
    148297   
    149298#-------------------------------------------------------------
    … …  
    162311    is_interactive = False
    163312    #is_interactive = True
    164     plot(scenarios, outputdir_tag,is_interactive=is_interactive,
    165          y_location_tag=':0.0')
    166 
     313    plot_exp_sim_comparision(scenarios, outputdir_tag,
     314                             is_interactive=is_interactive,
     315                             y_location_tag=':0.0')
     316

  • anuga_work/development/Hinwood_2008/run_dam.py

    r5664 r5670  
    5858def main(boundary_file,
    5959         metadata_dic,
     60         maximum_triangle_area,
     61         yieldstep,
    6062         boundary_path=None,
    6163         friction=0.012,  # planed wood. http://www.lmnoeng.com/manningn.htm
    6264         outputdir_name=None,
    63          run_type=0,
     65         run_type=1,
    6466         width=1.0,
    6567         use_limits=True,
    … …  
    6971    basename = 'zz_' + metadata_dic['scenario_id']
    7072   
    71     if run_type == 1:
     73    if run_type == 0:
    7274        yieldstep = 1.0
    7375        finaltime = 15.
    7476        maximum_triangle_area=0.1
    7577        outputdir_name += '_test'
     78         
     79    elif run_type == 1:
     80        finaltime = None
    7681       
    7782    elif run_type == 2:
    … …  
    131136        maximum_triangle_area=0.00001
    132137        #outputdir_name += '_good'
    133        
    134     if use_limits is True:
    135         outputdir_name += '_lmts'
    136     else:
    137         outputdir_name += '_nolmts'
    138     outputdir_name += '_wdth_' + str(width)
    139     outputdir_name += '_z_' + str(friction)
    140     outputdir_name += '_ys_' + str(yieldstep)
    141     outputdir_name += '_mta_' + str(maximum_triangle_area)
    142     outputdir_name += end_tag
    143    
     138
     139    outputdir_name = paras2outputdir_tag(mta=maximum_triangle_area,
     140                        yieldstep=yieldstep,
     141                        width=width,
     142                        run_type=run_type,
     143                        use_limits=use_limits,
     144                        friction=friction,
     145                        end_tag=end_tag,
     146                        outputdir_name=outputdir_name
     147                        )
     148
    144149    metadata_dic = set_z_origin_to_water_depth(metadata_dic)   
    145150       
    … …  
    180185
    181186    # creates copy of code in output dir
    182     if run_type >= 2:
     187    if run_type >= 1:
    183188        #start_screen_catcher(pro_instance.outputdir, rank, pypar.size())
    184189        start_screen_catcher(pro_instance.outputdir)
    … …  
    270275    for gauge_x in metadata_dic['gauge_x']:
    271276        points.append([gauge_x, flume_y_middle])
    272     print "points",points
    273277
    274278
    … …  
    289293    return pro_instance
    290294
     295
     296
     297def paras2outputdir_tag(mta,
     298                        yieldstep,
     299                        width,
     300                        use_limits,
     301                        friction,
     302                        end_tag,
     303                        outputdir_name=None
     304                        ):
     305
     306    if outputdir_name is None:
     307        outputdir_name = ''
     308       
     309    if use_limits is True:
     310        outputdir_name += '_lmts'
     311    else:
     312        outputdir_name += '_nolmts'
     313    outputdir_name += '_wdth_' + str(width)
     314    outputdir_name += '_z_' + str(friction)
     315    outputdir_name += '_ys_' + str(yieldstep)
     316    outputdir_name += '_mta_' + str(mta)
     317    outputdir_name += end_tag
     318
     319    return outputdir_name
     320
     321   
    291322def set_z_origin_to_water_depth(seabed_coords):
    292323    offset = seabed_coords['offshore_water_depth']
    … …  
    310341   
    311342    run_type = 1
    312     run_type = 4
     343    run_type = 0
     344   
    313345    #for run_data in [scenarios[5]]:
    314346    #scenarios = scenarios[3:]
    315347    #scenarios = [scenarios[0]]
     348   
    316349    width = 1.0
    317350    width = 0.1
    318351    #width = 0.01
     352
     353    maximum_triangle_area=0.01
     354    yieldstep = 0.01
     355    friction=0.0
     356   
    319357    for run_data in scenarios:
    320358        pro_instance = main( run_data['scenario_id'] + '_boundary.tsm'  ,
    321359                             run_data,
     360                             maximum_triangle_area=maximum_triangle_area,
     361                             yieldstep=yieldstep,
    322362                             width=width,
    323363                             run_type=run_type,
    324364                             outputdir_name=run_data['scenario_id'],
    325365                             use_limits=False,
    326                              friction=0.0,
     366                             friction=friction,
    327367                             end_tag='_I')
    328368        #gauges_for_slope(pro_instance.outputdir,[run_data])

  • anuga_work/development/Hinwood_2008/scenarios.py

    r5664 r5670  
    255255        }
    256256scenarios.append(data)
     257
     258# These are the indexes from 'gauge_x' of the gauges for the journal article.
     259to_publish_indexes = {'T1R5':[0, 5, 6],
     260              'T2R7':[1, 3, 5],
     261              'T3R29':[1, 6, 8],
     262              'T4R32':[0, 5, 8]
     263              }
     264
     265
     266
Note: See TracChangeset for help on using the changeset viewer.