Changeset 7647


Ignore:
Timestamp:
Mar 2, 2010, 2:31:30 PM (15 years ago)
Author:
Leharne Fountain
Message:

updates to mandurah and bunbury storm surge models

Location:
anuga_work/production
Files:
1 added
1 deleted
11 edited

Legend:

Unmodified
Added
Removed
  • anuga_work/production/bunbury_storm_surge_2009/build_elevation.py

    r7613 r7647  
    8484keylist = project.ascii_grid_filenames + project.point_filenames
    8585
    86 G1 = geospatial_data['swwa_10m'].clip_outside(project.elevation_clip_box)
    87 G2 = geospatial_data['bunbury_5m'].clip(project.elevation_clip_box)
     86G1 = geospatial_data['swwa_10m_IC'].clip_outside(project.elevation_clip_box)
     87G2 = geospatial_data['bunbury_5m_IC'].clip(project.elevation_clip_box)
    8888#~ G2 = geospatial_data['bunbury_nth_a'].clip(project.elevation_clip_box)
    8989#~ G3 = geospatial_data['bunbury_nth_b'].clip(project.elevation_clip_box)
     
    9494G5 = geospatial_data['Busselton_NavyFinal_Clip_ss.txt']
    9595G6 = geospatial_data['Leschenault_TIN.txt']
    96 G7 = geospatial_data['DPI5U1A02_01a_edited.txt']
    97 G8 = geospatial_data['DPI5U1A02_01b_edited.txt']
    98 G9 = geospatial_data['DPI5U1A02_01c_edited.txt']
    99 G10 = geospatial_data['DPI5U1A02_01d_edited.txt']
    100 G11 = geospatial_data['DPI5U1A02_01e_edited.txt']
     96G7 = geospatial_data['Leschenault_Inlet_TIN.txt']
     97G8 = geospatial_data['DPI5U1A02_01a_edited.txt']
     98G9 = geospatial_data['DPI5U1A02_01b_edited.txt']
     99G10 = geospatial_data['DPI5U1A02_01c_edited.txt']
     100G11 = geospatial_data['DPI5U1A02_01d_edited.txt']
     101G12 = geospatial_data['DPI5U1A02_01e_edited.txt']
    101102
    102103#~ G_list = [G1, G2, G3, G4, G5, G6, G7, G8, G9, G10]
     
    113114#~ print 'G10', G10.attributes
    114115
    115 G = G1 + G2 + G3 + G4 + G5 + G6 + G7 + G8 + G9 + G10 + G11
     116G = G1 + G2 + G3 + G4 + G5 + G6 + G7 + G8 + G9 + G10 + G11 + G12
    116117
    117118# G = None
     
    145146G10.export_points_file(join(project.topographies_folder, keylist[9] +'_export.txt'))
    146147G11.export_points_file(join(project.topographies_folder, keylist[10] +'_export.txt'))
     148G12.export_points_file(join(project.topographies_folder, keylist[11] +'_export.txt'))
  • anuga_work/production/bunbury_storm_surge_2009/get_timeseries.py

    r7613 r7647  
    1919directory = project.output_folder
    2020
    21 time_dir1 = '20090508_150215_run_final_0_51469_lfountai'
    22 time_dir2 = '20090511_161526_run_final_1.1_51469_kvanputt'
    23 time_dir3 = '20090511_165539_run_final_0_50863_lfountai'
    24 time_dir4 = '20090521_220101_run_final_1.1_50863_kvanputt'
    25 time_dir5 = '20090522_164526_run_final_0_51392_lfountai'
    26 time_dir6 = '20090522_164640_run_final_1.1_51392_lfountai'
    27 time_dir7 = '20090518_154710_run_final_0_50994_lfountai'
    28 time_dir8 = '20090519_160510_run_final_1.1_50994_lfountai'
    29 time_dir9 = '20090522_164948_run_final_0_51423_lfountai'
    30 time_dir10 = '20090522_165600_run_final_1.1_51423_lfountai'
     21time_dir1 = '20100103_102505_run_storm_surge_final_0_alby_coarse_lfountai'
    3122
    32 time_dirs = [time_dir1, time_dir2, time_dir3, time_dir4, time_dir5,
    33              time_dir6, time_dir7, time_dir8, time_dir9, time_dir10]
     23
     24time_dirs = [time_dir1]
    3425
    3526for time_dir in time_dirs:
  • anuga_work/production/bunbury_storm_surge_2009/project.py

    r7613 r7647  
    4646                                                #    storm_surge_final - as 'final' but with a longer yieldstep (12 mins)
    4747                                               
    48 #------------------------------------------------------------------------------
    49 #Floodgate parameters
    50 #------------------------------------------------------------------------------
    51 
    52 #changeable parameters specific to the floodgates
    53 
    54 floodgate_height = 4.4          #the relief of the floodgates above adjacent topography
    55 floodgate_thickness = 1         #the thickness of the floodgates in meters
    56 floodgate_LHS_x = 373546.02     #the x coordinate of the left hand side of the floodgates
    57 floodgate_LHS_y = 6312332.62    #the y coordinate of the left hand side of the floodgates
    58 floodgate_RHS_x = 373553.97     #the x coordinate of the right hand side of the floodgates
    59 floodgate_RHS_y = 6312331.97    #the y coordinate of the right hand side of the floodgates
    60 
    61 floodgate_start = 86400 #the time in seconds at which the floodgates begin to close
    62 floodgate_duration = 10 #the number of seconds that it takes the floodgates to close
    63 
    64 #polygon around the floodgates used to define an interior region when building mesh
    65 floodgate_boundary = [[floodgate_LHS_x,(floodgate_LHS_y - floodgate_thickness)],
    66                       [floodgate_RHS_x,(floodgate_RHS_y - floodgate_thickness)],
    67                       [floodgate_RHS_x,(floodgate_RHS_y + floodgate_thickness)],
    68                       [floodgate_LHS_x,(floodgate_LHS_y + floodgate_thickness)]]       
    69 floodgate_resolution = 1.0 #max triangle area to use when building mesh around floodgates
    7048
    7149#-------------------------------------------------------------------------------
     
    9169# Used in build_elevation.py
    9270# Format for ascii grids, as produced in ArcGIS + a projection file
    93 ascii_grid_filenames = ['swwa_10m',                     # LiDAR mosaiced and resampled to 10m
    94                                                                 'bunbury_5m']           # Bunbury LiDar grid resampled to 5 m - 1m data caused memory errors
    95                                                                 #~ 'bunbury_nth_a',     # 1m LiDAR from WA Dot - nothern quarter (file split as too large for dem2pts)
    96                                                                 #~ 'bunbury_nth_b',     # 1m LiDAR from WA Dot - 2nd northmost quarter (file split as too large for dem2pts)
    97                                                                 #~ 'bunbury_sth_a',     # 1m LiDAR from WA Dot - 2nd southmost quarter (file split as too large for dem2pts)
    98                                                                 #~ 'bunbury_sth_b']     # 1m LiDAR from WA Dot - southern quarter (file split as too large for dem2pts)
    99                                                 # 'bunbury_aoi']                        # 1m LiDAR from WA DoT - clipped to area_of_interest.csv
     71ascii_grid_filenames = ['swwa_10m_IC',                  # LiDAR mosaiced and resampled to 10m, clipped to IC to remove values in water
     72                                                'bunbury_5m_IC']        # Bunbury LiDar grid resampled to 5 m (1m data caused memory errors)
     73                                            # and clipped to Initial Conditions (to remove values in water)
     74                                                # 'bunbury_nth_a',      # 1m LiDAR from WA Dot - nothern quarter (file split as too large for dem2pts)
     75                                                # 'bunbury_nth_b',      # 1m LiDAR from WA Dot - 2nd northmost quarter (file split as too large for dem2pts)
     76                                                # 'bunbury_sth_a',      # 1m LiDAR from WA Dot - 2nd southmost quarter (file split as too large for dem2pts)
     77                                                # 'bunbury_sth_b']      # 1m LiDAR from WA Dot - southern quarter (file split as too large for dem2pts)
     78                                                # 'bunbury_aoi']                # 1m LiDAR from WA DoT - clipped to area_of_interest.csv
    10079                                                # 'bunbury_nth',                # 1m LiDAR from WA Dot - nothern half (file split as too large for dem2pts)
    10180                                                # 'bunbury_sth']                # 1m LiDAR data from WA DoT - southern half
     
    10382               
    10483# Format for point is x,y,elevation (with header)
    105 point_filenames = ['DPI.txt',                                                                           # Bathymetry data from DPI
    106                                                 'Busselton_Chart_Clip_ss.txt',                  # Clipped from Busselton_Chart  - see Busselton Tsunami Scenario 2009
    107                                                 'Busselton_NavyFinal_Clip_ss.txt',      # Clipped from Busselton_NavyFinal - see Busselton Tsunami Scenario 2009
    108                                                 'Leschenault_TIN.txt',                                                  # TIN created over the Leschenault Estuary
    109                                                 'DPI5U1A02_01a_edited.txt', 'DPI5U1A02_01b_edited.txt',
    110                                                 'DPI5U1A02_01c_edited.txt','DPI5U1A02_01d_edited.txt',
    111                                                 'DPI5U1A02_01e_edited.txt']                             # Bathymetric LiDAR from DPI     - split into manageable pieces         
     84point_filenames = ['DPI.txt',                                                   # Bathymetry data from DPI
     85                   'Busselton_Chart_Clip_ss.txt',               # Clipped from Busselton_Chart  - see Busselton Tsunami Scenario 2009
     86                   'Busselton_NavyFinal_Clip_ss.txt',   # Clipped from Busselton_NavyFinal - see Busselton Tsunami Scenario 2009
     87                   'Leschenault_TIN.txt',               # TIN created over the Leschenault Estuary
     88                   'Leschenault_inlet_TIN.txt'          # TIN created over the Leschenault Inlet
     89                   'DPI5U1A02_01a_edited.txt',          # Bathymetric LiDAR from DPI - split into manageable pieces
     90                   'DPI5U1A02_01b_edited.txt', 
     91                   'DPI5U1A02_01c_edited.txt',
     92                   'DPI5U1A02_01d_edited.txt',
     93                   'DPI5U1A02_01e_edited.txt']                                 
    11294
    11395# BOUNDING POLYGON - for data clipping and estimate of triangles in mesh
     
    122104interior_regions_data = [['intermediate.csv', 2500],
    123105                         ['area_of_interest.csv', 100],
    124                          [floodgate_boundary, floodgate_resolution]]
     106                         ['storm_gate_area.csv', 1],
     107                         ['floodgates.csv', floodgate_resolution]]
    125108
    126109# LAND - used to set the initial stage/water to be offcoast only
  • anuga_work/production/bunbury_storm_surge_2009/run_model.py

    r7613 r7647  
    159159                     'ocean': Bf})
    160160
    161 #-----------------------------------------------------------------------------
    162 #  Setup dynamic topography
    163 #-----------------------------------------------------------------------------
    164 from project import floodgate_height, floodgate_thickness,floodgate_LHS_x,floodgate_LHS_y,floodgate_RHS_x,floodgate_RHS_y,floodgate_start,floodgate_duration
    165 
    166 floodgate_end = floodgate_start + floodgate_duration #the time in seconds at which the floodgates are fully closed
    167 
    168 def floodgate(x,y):
    169     z = 0*x
    170     N = len(x)
    171     for i in range(N):
    172         ymin = ((floodgate_RHS_y - floodgate_LHS_y)/(floodgate_RHS_x - floodgate_LHS_x))*(x[i]-floodgate_LHS_x)+ floodgate_LHS_y - 0.5*floodgate_thickness
    173         ymax = ((floodgate_RHS_y - floodgate_LHS_y)/(floodgate_RHS_x - floodgate_LHS_x))*(x[i]-floodgate_LHS_x)+ floodgate_LHS_y + 0.5*floodgate_thickness
    174         xmin = floodgate_LHS_x
    175         xmax = floodgate_RHS_x
    176         if ymin < y[i] < ymax and xmin < x[i] < xmax:
    177             z[i] += (floodgate_height/((floodgate_duration)/project.yieldstep))
    178 
    179     return z
    180    
    181161#-------------------------------------------------------------------------------
    182162# Evolve system through time
    183163#-------------------------------------------------------------------------------
    184 
    185 raising_floodgates = False
    186 lowering_floodgates = False
    187164
    188165t0 = time.time()
     
    193170    print domain.boundary_statistics(tags='ocean')
    194171
    195 #raise floodgate at a given time interval (to make floodgate lower again repeat this code below
    196 #for a later time interval and change rising to false and falling to true
    197     if floodgate_start < t < floodgate_end:
    198         rising = True
    199         falling = False
    200 
    201     if rising:
    202         domain.add_quantity('elevation',floodgate)
    203 
    204     if falling:
    205         domain.add_quantity('elevation', lambda x,y: -1*floodgate(x,y))
    206172
    207173print 'Simulation took %.2f seconds' % (time.time()-t0)
  • anuga_work/production/mandurah_storm_surge_2009/build_elevation.py

    r7627 r7647  
    8787G1 = geospatial_data['m_peel_aoi'].clip(project.elevation_clip_box)
    8888G2 = geospatial_data['ph10m_ss'].clip_outside(project.elevation_clip_box)
    89 G3 = geospatial_data['MA-46893-SNDS_AHD.csv']
    90 G4 = geospatial_data['MS0205HY_AHD.csv']
    91 G5 = geospatial_data['MS0404_AHD.csv']
    92 G6 = geospatial_data['YU0403HY_AHD.csv']
    93 G7 = geospatial_data['original_data_ss.csv']
     89G3 = geospatial_data['DpiU1A03_ss.csv']
     90G4 = geospatial_data['MA-46893-SNDS_AHD.csv']
     91G5 = geospatial_data['MS0205HY_AHD.csv']
     92G6 = geospatial_data['MS0404_AHD.csv']
     93G7 = geospatial_data['YU0403HY_AHD.csv']
     94G8 = geospatial_data['original_data_ss.csv']
    9495
    95 G = G1 + G2 + G3 + G4 + G5 + G6 + G7
     96G = G1 + G2 + G3 + G4 + G5 + G6 + G7 + G8
    9697       
    9798print 'Export combined DEM file'
  • anuga_work/production/mandurah_storm_surge_2009/gems_comparison/build_elevation_GEMS.py

    r7628 r7647  
    8282print 'project_GEMS.elevation_clip_box', project_GEMS.elevation_clip_box
    8383
    84 G1 = geospatial_data['Man_25m'].clip(project_GEMS.elevation_clip_box)
     84G1 = geospatial_data['man_25m_prj'].clip(project_GEMS.elevation_clip_box)
    8585G2 = geospatial_data['man10m_ss'].clip_outside(project_GEMS.elevation_clip_box)
    8686
  • anuga_work/production/mandurah_storm_surge_2009/gems_comparison/project_GEMS.py

    r7628 r7647  
    6262#-------------------------------------------------------------------------------
    6363
    64 output_comment = [setup, tide, event, '_GEMS_compare']    # event_number will have to
     64output_comment = [setup, tide, event, 'GEMS_compare']    # event_number will have to
    6565                                                                                                # change to something relevent
    6666                                                                                                # for storm surge
     
    7777##                        'm_peel_10m.asc'
    7878##                        'm_harvey_10m.asc']
    79 ascii_grid_filenames = ['Man_25m',              # this is the latest 25m DEM from GEMS
     79ascii_grid_filenames = ['man_25m_prj',          # this is the latest 25m DEM from GEMS
    8080                                                'man10m_ss']    # this is to fill in areas not covered by the
    8181                                        # GEMS grid and is derived from the 10m resampled LiDAR from DoT,
  • anuga_work/production/mandurah_storm_surge_2009/gems_comparison/run_model_GEMS.py

    r7628 r7647  
    110110print domain.statistics()
    111111
    112 domain.starttime(project_GEMS.starttime)
     112domain.set_starttime(project_GEMS.starttime)
    113113domain.set_name(project_GEMS.scenario_name)
    114114domain.set_datadir(project_GEMS.output_run)
  • anuga_work/production/mandurah_storm_surge_2009/get_timeseries.py

    r7614 r7647  
    1919directory = project.output_folder
    2020
    21 time_dir1 = '20090508_150215_run_final_0_51469_lfountai'
    22 time_dir2 = '20090511_161526_run_final_1.1_51469_kvanputt'
    23 time_dir3 = '20090511_165539_run_final_0_50863_lfountai'
    24 time_dir4 = '20090521_220101_run_final_1.1_50863_kvanputt'
    25 time_dir5 = '20090522_164526_run_final_0_51392_lfountai'
    26 time_dir6 = '20090522_164640_run_final_1.1_51392_lfountai'
    27 time_dir7 = '20090518_154710_run_final_0_50994_lfountai'
    28 time_dir8 = '20090519_160510_run_final_1.1_50994_lfountai'
    29 time_dir9 = '20090522_164948_run_final_0_51423_lfountai'
    30 time_dir10 = '20090522_165600_run_final_1.1_51423_lfountai'
     21time_dir1 = '20091229_100831_run_storm_surge_final_0_alby_waves_lfountai'
    3122
    32 time_dirs = [time_dir1, time_dir2, time_dir3, time_dir4, time_dir5,
    33              time_dir6, time_dir7, time_dir8, time_dir9, time_dir10]
     23time_dirs = [time_dir1]
    3424
    3525for time_dir in time_dirs:
  • anuga_work/production/mandurah_storm_surge_2009/project.py

    r7627 r7647  
    8686# Format for point is x,y,elevation (with header)
    8787# Don't use these for comparison with GEMS inundation
    88 point_filenames = ['MA-46893-SNDS_AHD.csv',     # These files contain inf-fill bathymetry for
     88point_filenames = ['DpiU1A03_ss.csv',           # Bathymetric LiDAR data from DPI, clipped to bounding polygon
     89                  'MA-46893-SNDS_AHD.csv',      # These files contain inf-fill bathymetry for
    8990                  'MS0205HY_AHD.csv',           # the Peel & Harvey estuaries, as well as around
    9091                  'MS0404_AHD.csv',             # the canals
  • anuga_work/production/mandurah_storm_surge_2009/run_model.py

    r7614 r7647  
    110110print domain.statistics()
    111111
     112domain.set_starttime(project.starttime)
    112113domain.set_name(project.scenario_name)
    113114domain.set_datadir(project.output_run)
Note: See TracChangeset for help on using the changeset viewer.