Changeset 5689


Ignore:
Timestamp:
Aug 26, 2008, 4:25:37 PM (15 years ago)
Author:
duncan
Message:

Getting rmsd calc going

Location:
anuga_validation/Hinwood_2008
Files:
3 edited
8 moved

Legend:

Unmodified
Added
Removed
  • anuga_validation/Hinwood_2008/calc_norm.py

    r5687 r5689  
    2525from anuga.utilities.numerical_tools import ensure_numeric, err, norm
    2626
    27 from slope import load_sensors
    2827from anuga.utilities.interp import interp
    2928
     
    5655        id = run_data['scenario_id']
    5756        outputdir_name = id + outputdir_tag
    58         end = id + ".csv"       
    59         file_sim = quantity + "_" + end
    60         file_exp = id + '_pressfilt_exp_' + quantity + '.csv'
    61         file_err = pro_instance.rmsd_dir + sep + outputdir_name + "_" \
    62                    + quantity + "_err.csv"
     57        file_sim = outputdir_name + '_' +  quantity + ".csv"
     58        file_exp = id + '_exp_' + quantity + '.csv'
     59        file_err = outputdir_name + "_" + quantity + "_err.csv"
    6360       
    6461
     
    121118
    122119
    123                  
     120
     121def load_sensors(quantity_file):
     122    """
     123    Load a csv file, where the first row is the column header and
     124    the first colum explains the rows.
     125    """
     126    #slope, _ = csv2dict(file_sim)
     127   
     128    # Read the depth file
     129    dfid = open(quantity_file)
     130    lines = dfid.readlines()
     131    dfid.close()
     132
     133    title = lines.pop(0)
     134    n_time = len(lines)
     135    n_sensors = len(lines[0].split(','))-1  # -1 to remove time
     136    times = zeros(n_time, Float)  #Time
     137    depths = zeros(n_time, Float)  #
     138    sensors = zeros((n_time,n_sensors), Float)
     139    quantity_locations = title.split(',') #(',')
     140    quantity_locations.pop(0) # remove 'time'
     141
     142    # Doing j.split(':')[0] drops the y location
     143    locations = [float(j.split(':')[0]) for j in quantity_locations]
     144   
     145    for i, line in enumerate(lines):
     146        fields = line.split(',') #(',')
     147        fields = [float(j) for j in fields]
     148        times[i] = fields[0]
     149        sensors[i] = fields[1:] # 1: to remove time
     150
     151    #print "times",times
     152    #print "locations", locations
     153    #print "sensors", sensors
     154    return times, locations, sensors                 
    124155   
    125156       
     
    128159# Return a bunch of lists
    129160# The err files, for all scenarios
    130 def err_files(scenarios, outputdir_tag, rmsd_dir, quantity):
     161def err_files(scenarios, outputdir_tag, quantity='stage'):
    131162    """
    132163    The err files, for a list of scenarios
     
    136167        id = scenario['scenario_id']
    137168        outputdir_name = id + outputdir_tag
    138         file_err = rmsd_dir + sep + outputdir_name + "_" \
    139                    + quantity + "_err.csv"
     169        file_err =  outputdir_name + "_" + quantity + "_err.csv"
    140170        file_errs.append(file_err)
    141171    return file_errs
    142172   
    143173
    144 def compare_different_settings(outputdir_tags, scenarios, quantity,
    145                    save_as=None,
    146                    is_interactive=False):
    147 
    148     # A bit hacky.  Getting a pro_instance to get the rmsd_dir.
    149     outputdir_name = scenarios[0]['scenario_id'] + outputdir_tags[0]
    150     pro_instance = project.Project(['data','flumes','Hinwood_2008'],
    151                                    outputdir_name=outputdir_name)
    152    
    153      
    154     settings = {} # keys are different settings.
    155     # For each setting there will be err and amount
    156        
    157     for outputdir_tag in outputdir_tags:
    158         files = err_files(scenarios, outputdir_tag, pro_instance.rmsd_dir,
    159                           quantity)
    160         sim = {}
    161         for run_data, file in map(None, scenarios, files):
    162            
    163             simulation, _ = csv2dict(file)
    164             #locations = [float(x) for x in simulation['x location']]
    165             err_list = [float(x) for x in simulation['err']]
    166             amount_list = [float(x) for x in simulation['Number_of_samples']]
    167            
    168             if sim.has_key('err'):
    169                 err_list.append(sim['err'])
    170                 amount_list.append(sim['amount'])
    171             sim['err'], sim['amount'] = err_addition(err_list,
    172                                                      amount_list)
    173         sim['rmsd'] = sim['err']/sqrt(sim['amount'])
    174         settings[outputdir_tag] = sim
    175     #print "settings", settings
    176 
    177     aux = [(settings[k]['rmsd'], k) for k in settings.keys()]
    178     aux.sort()
    179     for val, key in aux:
    180         print key + "   " + str(val)
    181    
    182    
    183    
    184 def err_addition(err_list, amount_list):
     174def compare_different_settings(outputdir_tag, scenarios, quantity='stage'):
     175
     176    files = err_files(scenarios, outputdir_tag, quantity=quantity)
     177    err = 0.0
     178    number_of_samples = 0
     179    for run_data, file in map(None, scenarios, files):
     180       
     181        simulation, _ = csv2dict(file)
     182        err_list = [float(x) for x in simulation['err']]
     183        number_of_samples_list = [float(x) for x in \
     184                                  simulation['Number_of_samples']]
     185       
     186        if number_of_samples is not 0:
     187            err_list.append(err)
     188            number_of_samples_list.append(number_of_samples)
     189        err, number_of_samples = err_addition(err_list, number_of_samples_list)
     190    rmsd = err/sqrt(number_of_samples)
     191    print outputdir_tag + "   " + str(rmsd)
     192   
     193   
     194   
     195def err_addition(err_list, number_of_samples_list):
    185196    """
    186197    err1 is the err value (sqrt(sum_over_x&y((xi - yi)^2))) for a set of values
    187     amount1 is the number of values associated with the err.
     198    number_of_samples1 is the number of values associated with the err.
    188199   
    189200    If this function gets used alot, maybe pull this out and make it an object
    190201    """
    191202    err = norm(ensure_numeric(err_list))
    192     amount = sum(ensure_numeric(amount_list))
    193 
    194     return err, amount
    195 
    196 
    197 def calc_max_rmsd(scenarios, outputdir_tags, quantity):
    198    
    199     max_rmsd = 0
    200    
    201     for run_data in scenarios:
    202         id = run_data['scenario_id']
    203         for outputdir_tag in outputdir_tags:
    204            
    205             outputdir_name = id + outputdir_tag
    206             pro_instance = project.Project(['data','flumes','Hinwood_2008'],
    207                                            outputdir_name=outputdir_name)
    208            
    209             file_err = pro_instance.rmsd_dir + sep + outputdir_name + "_" \
    210                        + quantity + "_err.csv"
    211            
    212             simulation, _ = csv2dict(file_err)
    213             rmsd_list = [float(x) for x in simulation['rmsd']]
    214             max_rmsd = max(max(rmsd_list), max_rmsd)
    215     return max_rmsd
     203    number_of_samples = sum(ensure_numeric(number_of_samples_list))
     204
     205    return err, number_of_samples
    216206
    217207                 
     
    222212    from scenarios import scenarios
    223213   
    224     outputdir_tags = []
    225     outputdir_tags.append("_lmts_wdth_0.1_z_0.0_ys_0.01_mta_0.01_I")
    226     outputdir_tags.append("_no_velocity_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.01_I")
    227     outputdir_tags.append("_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.0001_I")
    228     outputdir_tags.append("_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.001_I")
    229     outputdir_tags.append("_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.01_I")
    230     outputdir_tags.append("_nolmts_wdth_0.1_z_0.012_ys_0.01_mta_0.01_I")
    231     outputdir_tags.append("_nolmts_wdth_1.0_z_0.0_ys_0.01_mta_0.01_I")
    232     outputdir_tags.append("_nolmts_wdth_0.1_z_0.012_ys_0.01_mta_0.0001_I")
    233     outputdir_tags.append("_lmts_wdth_0.1_z_0.0_ys_0.01_mta_0.0001_I")
    234     outputdir_tags.append("_no_velocity_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.0001_I")
    235    
    236     #outputdir_tag = "_test_limiterC"
    237     #scenarios = [scenarios[0]] # !!!!!!!!!!!!!!!!!!!!!!
     214   
     215    scenarios = [scenarios[0]] # !!!!!!!!!!!!!!!!!!!!!!
    238216    #scenarios = scenarios[4:] # !!!!!!!!!!!!!!!!!!!!!!
     217
     218   
     219    outputdir_tag = "_nolmts_wdth_0.1_z_0.0_ys_0.5_mta_0.01"
    239220    calc_norms = True
    240221    #calc_norms = False
    241222    if calc_norms:
    242         for outputdir_tag in outputdir_tags:
    243             auto_rrms(outputdir_tag, scenarios, "stage", y_location_tag=':0.0')
    244    
    245     #scenarios = [scenarios[0]] # !!!!!!!!!!!!!!!!!!!!!!
    246     #auto_plot_rrms_sensor_settings(outputdir_tags, scenarios, "stage")
    247     compare_different_settings(outputdir_tags, scenarios, "stage")
    248 
    249     outputdir_tag = "_nolmts_wdth_0.1_z_0.0_ys_0.01_mta_0.01_I"
    250     #auto_plot_test_rmsd(scenarios, outputdir_tag)
     223        auto_rrms(outputdir_tag, scenarios, "stage", y_location_tag=':0.0')
     224    compare_different_settings(outputdir_tag, scenarios, "stage")
     225   
  • anuga_validation/Hinwood_2008/run_dam.py

    r5686 r5689  
    134134            points.append([gauge_x, 0.0])
    135135
    136         id = metadata_dic['scenario_id'] + ".csv"
     136        id = ".csv"
    137137        interpolate_sww2csv(basename +".sww",
    138138                            points,
    139                             basename + "_depth_" + id,
    140                             basename + "_velocity_x_" + id,
    141                             basename + "_velocity_y_" + id,
    142                             basename + "_stage_" + id)
     139                            basename + "_depth" + id,
     140                            basename + "_velocity_x" + id,
     141                            basename + "_velocity_y" + id,
     142                            basename + "_stage" + id)
    143143 
    144144
     
    196196    from scenarios import scenarios
    197197
    198     scenarios = [scenarios[0]]
     198    #scenarios = [scenarios[0]]
    199199   
    200200    width = 0.1
    201201    maximum_triangle_area=0.01
    202202    yieldstep = 0.01
    203     yieldstep = 0.5
     203    #yieldstep = 0.5
    204204    friction=0.0
    205205    isTest=True
    206     #isTest=False
     206    isTest=False
    207207   
    208208    for run_data in scenarios:
  • anuga_validation/Hinwood_2008/scenarios.py

    r5686 r5689  
    44# A list of scenario dictionaries
    55scenarios = []
    6 
    7 
    8 data = {'xleft':[-3.106,0.0],  # Av' of ADV and Gauge A
    9         'xtoe':[0.0,0.0],
    10         'xbeach':[1.285,0.090],
    11         'xright':[8,.4843334],
    12         'offshore_water_depth':.4,
    13         'scenario_id':'S1R2',
    14         'gauge_names':['B','1','2','3','4','5','7'],
    15         'gauge_x':[-0.68, 1.569, 2.568, 3.566, 4.564, 5.562, 7.559],
    16         'gauge_bed_elevation':[-0.400000, -0.293158,
    17                                -0.234473, -0.175788, -0.117104,
    18                                -0.058419, 0.058950],
    19         'start_slope_x':4,
    20         'finish_slope_x':6,
    21         'break_xs':[5.45, 5.45, 5.28, 5.28, 5.29, 5.29],
    22         'break_type':['none', 'plunge', 'plunge', 'plunge', 'plunge', \
    23                       'plunge'],
    24         'axis':[0,80,-0.040,0.040],
    25         'axis_maximum_x':6.0,
    26         'ANUGA_start_time':1.83,
    27         'band_offset':-1.0,
    28         'wave_times':[23.0,59.0]  # this is in Anuga time
    29                  }
    30 scenarios.append(data)
    316
    327
     
    5631
    5732
     33data = {'xleft':[-3.106,0.0],  # Av' of ADV and Gauge A
     34        'xtoe':[0.0,0.0],
     35        'xbeach':[1.285,0.090],
     36        'xright':[8,.4843334],
     37        'offshore_water_depth':.4,
     38        'scenario_id':'S1R2',
     39        'gauge_names':['B','1','2','3','4','5','7'],
     40        'gauge_x':[-0.68, 1.569, 2.568, 3.566, 4.564, 5.562, 7.559],
     41        'gauge_bed_elevation':[-0.400000, -0.293158,
     42                               -0.234473, -0.175788, -0.117104,
     43                               -0.058419, 0.058950],
     44        'start_slope_x':4,
     45        'finish_slope_x':6,
     46        'break_xs':[5.45, 5.45, 5.28, 5.28, 5.29, 5.29],
     47        'break_type':['none', 'plunge', 'plunge', 'plunge', 'plunge', \
     48                      'plunge'],
     49        'axis':[0,80,-0.040,0.040],
     50        'axis_maximum_x':6.0,
     51        'ANUGA_start_time':1.83,
     52        'band_offset':-1.0,
     53        'wave_times':[23.0,59.0]  # this is in Anuga time
     54                 }
     55scenarios.append(data)
     56
     57
    5858
    5959data = {'xleft':[-4.586,0.0],  # Av' of ADV and Gauge A
     
    109109        'xright':[14.,.3903881],
    110110        'offshore_water_depth':.336,
     111        'scenario_id':'S3R1',
     112        'gauge_names':['B','1','2','3','5','7',
     113                       '9','10','11','12'],
     114
     115        'gauge_x':[-0.325, 1.572, 2.571, 3.571, 5.57,
     116                   7.57, 9.569, 10.569, 11.569,
     117                   12.569],
     118
     119        # remove
     120        'gauge_bed_elevation':[-0.336000, -0.237263, -0.213789, -0.190315,
     121                               -0.143368,
     122                               -0.096420, -0.049472,
     123                               -0.025998, -0.002524, 0.020949],
     124           
     125        'start_slope_x':8.,
     126        'finish_slope_x':9.5,
     127        'break_xs':[9.063,9.063,9.043,9.043],
     128        'break_type':['none','collapse', 'collapse', 'collapse'],
     129        'axis':[0,80,-0.040,0.040],
     130        'axis_maximum_x':12.0,
     131        'ANUGA_start_time':10.48,
     132        'band_offset':-2.0,
     133        'wave_times':[30.0,74.0] # this is in Anuga time
     134                 }
     135scenarios.append(data)
     136
     137data = {'xleft':[-3.875,0.0],  # Av' of ADV and Gauge A
     138        'xtoe':[0.0,0.0],
     139        'xbeach':[1.285,0.090],
     140        'xright':[14.,.3903881],
     141        'offshore_water_depth':.336,
    111142        'scenario_id':'S3R2',
    112143        'gauge_names':['B','1','2','3','5','7',
     
    132163        'band_offset':-0.5,
    133164        'wave_times':[30.0,85.0] # this is in Anuga time
    134         }
    135 scenarios.append(data)
    136 
    137 data = {'xleft':[-3.875,0.0],  # Av' of ADV and Gauge A
    138         'xtoe':[0.0,0.0],
    139         'xbeach':[1.285,0.090],
    140         'xright':[14.,.3903881],
    141         'offshore_water_depth':.336,
    142         'scenario_id':'S3R1',
    143         'gauge_names':['B','1','2','3','5','7',
    144                        '9','10','11','12'],
    145 
    146         'gauge_x':[-0.325, 1.572, 2.571, 3.571, 5.57,
    147                    7.57, 9.569, 10.569, 11.569,
    148                    12.569],
    149 
    150         # remove
    151         'gauge_bed_elevation':[-0.336000, -0.237263, -0.213789, -0.190315,
    152                                -0.143368,
    153                                -0.096420, -0.049472,
    154                                -0.025998, -0.002524, 0.020949],
    155            
    156         'start_slope_x':8.,
    157         'finish_slope_x':9.5,
    158         'break_xs':[9.063,9.063,9.043,9.043],
    159         'break_type':['none','collapse', 'collapse', 'collapse'],
    160         'axis':[0,80,-0.040,0.040],
    161         'axis_maximum_x':12.0,
    162         'ANUGA_start_time':10.48,
    163         'band_offset':-2.0,
    164         'wave_times':[30.0,74.0] # this is in Anuga time
    165                  }
    166 scenarios.append(data)
    167 
    168 data = {'xleft':[-2.43,0.0],  # Av' of ADV and Gauge A
    169         'xtoe':[0.0,0.0],
    170         'xbeach':[1.285,0.090],
    171         'xright':[14.,.3903881],
    172         'offshore_water_depth':.336,
    173         'scenario_id':'S4R2',     
    174         'gauge_names':['B','1','2','3','5','7',
    175                        '9','10','11','12'],
    176 
    177         'gauge_x':[-0.325, 1.572, 2.571, 3.571, 5.57,
    178                    7.57, 9.569, 10.569, 11.569,
    179                    12.569],
    180 
    181         # remove
    182         'gauge_bed_elevation':[-0.336000, -0.237263, -0.213789, -0.190315,
    183                                -0.143368,
    184                                -0.096420, -0.049472,
    185                                -0.025998, -0.002524, 0.020949],
    186        
    187         'start_slope_x':6,
    188         'finish_slope_x':9,
    189         'break_xs':[8.53,7.46,7.492,7.492,7.444],
    190         'break_type':['collapse', 'spill', 'spill', 'spill', 'spill'],
    191         'axis':[0,80,-0.040,0.040],
    192         'axis_maximum_x':12.0,
    193         'ANUGA_start_time':11.63,
    194         'band_offset':-0.5,
    195         'wave_times':[34.0,75.0] # this is in Anuga time
    196165        }
    197166scenarios.append(data)
     
    230199
    231200
     201data = {'xleft':[-2.43,0.0],  # Av' of ADV and Gauge A
     202        'xtoe':[0.0,0.0],
     203        'xbeach':[1.285,0.090],
     204        'xright':[14.,.3903881],
     205        'offshore_water_depth':.336,
     206        'scenario_id':'S4R2',     
     207        'gauge_names':['B','1','2','3','5','7',
     208                       '9','10','11','12'],
     209
     210        'gauge_x':[-0.325, 1.572, 2.571, 3.571, 5.57,
     211                   7.57, 9.569, 10.569, 11.569,
     212                   12.569],
     213
     214        # remove
     215        'gauge_bed_elevation':[-0.336000, -0.237263, -0.213789, -0.190315,
     216                               -0.143368,
     217                               -0.096420, -0.049472,
     218                               -0.025998, -0.002524, 0.020949],
     219       
     220        'start_slope_x':6,
     221        'finish_slope_x':9,
     222        'break_xs':[8.53,7.46,7.492,7.492,7.444],
     223        'break_type':['collapse', 'spill', 'spill', 'spill', 'spill'],
     224        'axis':[0,80,-0.040,0.040],
     225        'axis_maximum_x':12.0,
     226        'ANUGA_start_time':11.63,
     227        'band_offset':-0.5,
     228        'wave_times':[34.0,75.0] # this is in Anuga time
     229        }
     230scenarios.append(data)
     231
Note: See TracChangeset for help on using the changeset viewer.