Changeset 3701


Ignore:
Timestamp:
Oct 5, 2006, 5:09:30 PM (17 years ago)
Author:
sexton
Message:

hobart script updates

Location:
anuga_work/production/hobart_2006
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • anuga_work/production/hobart_2006/project.py

    r3695 r3701  
    2626scenario_dir_name = 'hobart_tsunami_scenario_2006'
    2727
    28 # data provided by Tas SES and checked by NM&I
    29 #onshore_name = 'hob3_topo' # original
    30 onshore_name = 'hob5_topo' # 12.5m grid and clipped to 100m elevation or 3000m from coast
    31 #onshore_name_25 = 'hob5_topo_25m' # 25m grid and clipped to 100m elevation or 3000m from coast
    32 onshore_name_25 = 'hob6_topo_25m' # 25m grid NOT clipped
    33 offshore_name_tas1 = 'derwent_2m'
    34 offshore_name_tas2 = 'derwent_5m'
    35 offshore_name_tas3 = 'south_east_tas' #actually this is AHO
    36 offshore_name_tas4 = 'hobart_1m'
    37 
    38 # AHO data and checked by NM&I
    39 offshore_name1 = 'xy100003760'
    40 offshore_name2 = 'xy100003761'
    41 offshore_name3 = 'xy100003762'
    42 offshore_name4 = 'xy100003907'
    43 offshore_name5 = 'xy100003908'
    44 offshore_name6 = 'xy100003909'
    45 offshore_name7 = 'xy100003910'
    46 offshore_name8 = 'xy100003932'
    47 offshore_name9 = 'xy100003933'
    48 offshore_name10 = 'xy100003934'
    49 offshore_name11 = 'xy100003935'
    50 offshore_name12 = 'xy100003936'
    51 offshore_name13 = 'xy100003964'
    52 offshore_name14 = 'xy100014250'
    53 offshore_name15 = 'xy100014253'
    54 offshore_name16 = 'xy100016142'
    55 
    56 # developed by NM&I
    57 coast_name = 'coastline_points'
    58 
    59 boundary_basename = 'puysegur'
     28
     29##
     30### data provided by Tas SES and checked by NM&I
     31###onshore_name = 'hob3_topo' # original
     32###onshore_name_25 = 'hob5_topo_25m' # 25m grid and clipped to 100m elevation or 3000m from coast
     33##onshore_name = 'hob5_topo' # 12.5m grid and clipped to 100m elevation or 3000m from coast
     34##onshore_name_25 = 'hob6_topo_25m' # 25m grid NOT clipped
     35##offshore_name_tas1 = 'derwent_2m'
     36##offshore_name_tas2 = 'derwent_5m'
     37##offshore_name_tas3 = 'south_east_tas' #actually this is AHO
     38##offshore_name_tas4 = 'hobart_1m'
     39##
     40### AHO data and checked by NM&I
     41##offshore_name1 = 'xy100003760'
     42##offshore_name2 = 'xy100003761'
     43##offshore_name3 = 'xy100003762'
     44##offshore_name4 = 'xy100003907'
     45##offshore_name5 = 'xy100003908'
     46##offshore_name6 = 'xy100003909'
     47##offshore_name7 = 'xy100003910'
     48##offshore_name8 = 'xy100003932'
     49##offshore_name9 = 'xy100003933'
     50##offshore_name10 = 'xy100003934'
     51##offshore_name11 = 'xy100003935'
     52##offshore_name12 = 'xy100003936'
     53##offshore_name13 = 'xy100003964'
     54##offshore_name14 = 'xy100014250'
     55##offshore_name15 = 'xy100014253'
     56##offshore_name16 = 'xy100016142'
     57##
     58### developed by NM&I
     59##coast_name = 'coastline_points'
     60
     61boundary_basename = 'puysegur' # Mw 8.7
     62#boundary_basename = 'puysegur_clip' # Mw 8.5
     63
    6064
    6165#swollen/ all data output
     
    7680polygondir = home+sep+state+sep+scenario_dir_name+sep+'anuga'+sep+'polygons'+sep
    7781
    78 gauge_filename = gaugedir + 'hobart_gauges.csv'
    79 buildings_filename = gaugedir + 'hobart_res.csv'
    80 buildings_filename_damage_out = 'hobart_res_modified.csv'
    81 gaugetimeseries = gaugedir + 'hobart'
     82gauge_filename = gaugedir + 'hobart_gauges_final.csv'
     83#buildings_filename = gaugedir + 'hobart_res.csv'
     84#buildings_filename_damage_out = 'hobart_res_modified.csv'
     85#gaugetimeseries = gaugedir + 'hobart'
    8286
    8387# boundary source data
     
    8892meshname = outputtimedir + 'mesh_' + basename
    8993
    90 onshore_dem_name = datadir + onshore_name
    91 onshore_dem_name_25 = datadir + onshore_name_25
    92 all_onshore_dem_name = datadir + 'combined_onshore'
    93 offshore_dem_name_local1 = datadir + offshore_name_tas1
    94 offshore_dem_name_local2 = datadir + offshore_name_tas2
    95 offshore_dem_name_local3 = datadir + offshore_name_tas3
    96 offshore_dem_name_local4 = datadir + offshore_name_tas4
    97 offshore_dem_name_aho1 = datadir + offshore_name1
    98 offshore_dem_name_aho2 = datadir + offshore_name2
    99 offshore_dem_name_aho3 = datadir + offshore_name3
    100 offshore_dem_name_aho4 = datadir + offshore_name4
    101 offshore_dem_name_aho5 = datadir + offshore_name5
    102 offshore_dem_name_aho6 = datadir + offshore_name6
    103 offshore_dem_name_aho7 = datadir + offshore_name7
    104 offshore_dem_name_aho8 = datadir + offshore_name8
    105 offshore_dem_name_aho9 = datadir + offshore_name9
    106 offshore_dem_name_aho10 = datadir + offshore_name10
    107 offshore_dem_name_aho11 = datadir + offshore_name11
    108 offshore_dem_name_aho12 = datadir + offshore_name12
    109 offshore_dem_name_aho13 = datadir + offshore_name13
    110 offshore_dem_name_aho14 = datadir + offshore_name14
    111 offshore_dem_name_aho15 = datadir + offshore_name15
    112 offshore_dem_name_aho16 = datadir + offshore_name16
    113 coast_dem_name = datadir + coast_name
     94# Necessary if using point datasets, rather than grid
     95##onshore_dem_name = datadir + onshore_name
     96##onshore_dem_name_25 = datadir + onshore_name_25
     97##all_onshore_dem_name = datadir + 'combined_onshore'
     98##offshore_dem_name_local1 = datadir + offshore_name_tas1
     99##offshore_dem_name_local2 = datadir + offshore_name_tas2
     100##offshore_dem_name_local3 = datadir + offshore_name_tas3
     101##offshore_dem_name_local4 = datadir + offshore_name_tas4
     102##offshore_dem_name_aho1 = datadir + offshore_name1
     103##offshore_dem_name_aho2 = datadir + offshore_name2
     104##offshore_dem_name_aho3 = datadir + offshore_name3
     105##offshore_dem_name_aho4 = datadir + offshore_name4
     106##offshore_dem_name_aho5 = datadir + offshore_name5
     107##offshore_dem_name_aho6 = datadir + offshore_name6
     108##offshore_dem_name_aho7 = datadir + offshore_name7
     109##offshore_dem_name_aho8 = datadir + offshore_name8
     110##offshore_dem_name_aho9 = datadir + offshore_name9
     111##offshore_dem_name_aho10 = datadir + offshore_name10
     112##offshore_dem_name_aho11 = datadir + offshore_name11
     113##offshore_dem_name_aho12 = datadir + offshore_name12
     114##offshore_dem_name_aho13 = datadir + offshore_name13
     115##offshore_dem_name_aho14 = datadir + offshore_name14
     116##offshore_dem_name_aho15 = datadir + offshore_name15
     117##offshore_dem_name_aho16 = datadir + offshore_name16
     118##coast_dem_name = datadir + coast_name
     119
     120
    114121# addition once total grid delivered
    115122onshore_offshore_dem_name_25 = datadir + '25m_se_tas' #25m grid
    116 onshore_offshore_dem_name = datadir + '50m_se_tas' #50m grid
     123onshore_offshore_dem_name    = datadir + '50m_se_tas' #50m grid
     124
    117125# output names
    118 bruny_dem_name_25 = datadir + 'bruny_25_dem'
    119 hobart_dem_name_25 = datadir + 'hobart_25_dem'
    120 combined_dem_name = datadir + 'hobart_combined_elevation'
     126bruny_dem_name_25   = datadir + 'bruny_25_dem'
     127hobart_dem_name_25  = datadir + 'hobart_25_dem'
     128combined_dem_name   = datadir + 'hobart_combined_elevation'
    121129combined_dem_name_2 = datadir + 'hobart_combined_elevation_2'
    122130
    123 outputname = outputtimedir + basename  #Used by post processing
     131#outputname = outputtimedir + basename  #Used by post processing
     132
    124133
    125134###############################
     
    127136###############################
    128137
    129 # bounding box
     138# bounding box for clipping MOST output (much bigger than study area)
    130139south = degminsec2decimal_degrees(-44,45,0)
    131140north = degminsec2decimal_degrees(-42,0,0)
    132 west = degminsec2decimal_degrees(146,45,0)
    133 east = degminsec2decimal_degrees(148,25,0) #degminsec2decimal_degrees(148,15,0)
    134 
    135 #Main Domain of Hobart: first run JS 18/9/06
     141west  = degminsec2decimal_degrees(146,45,0)
     142east  = degminsec2decimal_degrees(148,25,0)
     143
     144###Main Domain of Hobart:
    136145d0 = [south, west]
    137146d1 = [south, east]
     
    139148d3 = [north, west]
    140149polyAll2, zone = convert_points_from_latlon_to_utm([d0, d1, d2, d3])
    141 print "polyAll2", polyAll2
    142150refzone = zone
    143151
     
    145153#polyAll = [[520000, 5170000],[580000, 5170000],[580000, 5200000],[590000,5240000],[520000,5260000]]
    146154
    147 # Third run - afternoon Wed 27 Sep; surrounds -100m and 20mish elevation
    148 #polyAll = read_polygon(polygondir+'new_extent.csv')
    149 #plot_polygons([polyAll],'boundingpoly',verbose=False)
    150155
    151156# Mark run - morning Fri 29 Sep; surrounds -100m and 20mish elevation
     
    158163###################################################################
    159164
    160 # region to export for inundation map
    161 e_min_area = 520000#490000
    162 e_max_area = 580000#580000
    163 n_min_area = 5180000#5160000
    164 n_max_area = 5260000#5275000
    165 
    166 # clipping 12.5m onshore data set
    167 eastingmin = 520000
    168 eastingmax = 536000
    169 northingmin = 5245000
    170 northingmax = 5265000
    171 
    172 ### clipping 25m data set - Hobart
    173 ##eastingmin25  = 522000
    174 ##eastingmax25  = 555000
    175 ##northingmin25 = 5234000
    176 ##northingmax25 = 5260000
    177165
    178166# clipping 25m data set - Hobart
    179 eastingmin25  = 524208.387
    180 eastingmax25  = 554867.24
     167eastingmin25 = 524208.387
     168eastingmax25 = 554867.24
    181169northingmin25 = 5229154.555
    182170northingmax25 = 5258511.857
    183 
    184171
    185172# clipping 25m data set - Bruny
     
    225212poly_hobart2 = read_polygon(polygondir+'Hob_poly2.csv')
    226213poly_hobart3 = read_polygon(polygondir+'Hob_poly3.csv')
    227 poly_hobart4 = read_polygon(polygondir+'Hob_poly4.csv')
    228 
    229 plot_polygons([polyAll, poly_hobart1,poly_hobart2,poly_hobart3, poly_hobart4],'boundingpoly2',verbose=False)
     214
     215plot_polygons([polyAll, poly_hobart1,poly_hobart2,poly_hobart3],'boundingpoly2',verbose=False)
  • anuga_work/production/hobart_2006/run_hobart.py

    r3695 r3701  
    1010Ole Nielsen and Duncan Gray, GA - 2005 and Nick Bartzis, GA - 2006
    1111"""
    12 #-------------------------------------------------------------------------------# Import necessary modules
     12#-------------------------------------------------------------------------------
     13# Import necessary modules
    1314#-------------------------------------------------------------------------------
    1415
     
    4950#used to catch screen output to file
    5051sys.stdout = Screen_Catcher(screen_output_name)
    51 #sys.stderr = Screen_Catcher(screen_output_name)
    5252sys.stderr = Screen_Catcher(screen_error_name)
    5353
     
    7676# 25m data (clipping the around the Hobart area)
    7777convert_dem_from_ascii2netcdf(onshore_offshore_dem_name_25, use_cache=True, verbose=True)
     78
    7879# creates pts file for 25m data around Hobart
    7980dem2pts(onshore_offshore_dem_name_25, project.hobart_dem_name_25,
     
    9192
    9293# 25m data (clipping the around site 24 on Bruny Island)
    93 convert_dem_from_ascii2netcdf(onshore_offshore_dem_name_25, use_cache=True, verbose=True)
     94#convert_dem_from_ascii2netcdf(onshore_offshore_dem_name_25, use_cache=True, verbose=True)
     95
    9496# creates pts file for 25m data around site 24 at Bruny Island
    9597dem2pts(onshore_offshore_dem_name_25, project.bruny_dem_name_25,
     
    106108                                 project.combined_dem_name_2 + '.pts')
    107109
     110# If working with points, then add separate datasets together here
    108111# create geospatial data set and export
    109112#G = Geospatial_data(file_name = project.onshore_offshore_dem_name + '.pts')
    110113#G.export_points_file(project.combined_dem_name + '.pts')
     114
    111115
    112116#----------------------------------------------------------------------------
     
    119123
    120124# use 75 for onshore components (12.5m DEM)
    121 island_res = 35000
    122 hobart_res = 5000
    123 peninsula_res = 35000
     125# Resolution refers to max triangle area in a region
     126#island_res    = 35000
     127hobart_res    =  7000
     128#peninsula_res = 35000
    124129interior_regions = [[project.poly_hobart1, hobart_res],
    125130                    [project.poly_hobart2, hobart_res],
    126                     [project.poly_hobart3, hobart_res],
    127                     [project.poly_hobart4, hobart_res]]
     131                    [project.poly_hobart3, hobart_res]]
    128132
    129133print 'number of interior regions', len(interior_regions)
    130134
     135
     136# Number tags correspond to number segments from bounding polygon
     137# maximum_triangle_area refers to max triangle area in a region
     138# Interion_regions either comment out (perhaps empty list)
    131139from caching import cache
    132140_ = cache(create_mesh_from_regions,
     
    138146                              'e12': [12], 'e13': [13], 'e14': [14],
    139147                              'e15': [15]},
    140            'maximum_triangle_area': 200000,
     148           'maximum_triangle_area': 750000,
    141149           'filename': meshname,
    142150           'interior_regions': interior_regions},
     
    159167domain.set_store_vertices_uniquely(False)  # for writting to sww
    160168
     169
    161170#-------------------------------------------------------------------------------                                 
    162171# Setup initial conditions
     
    164173
    165174tide = 0.0
    166 
    167175domain.set_quantity('stage', tide)
    168176domain.set_quantity('friction', 0.0)
    169177
    170178domain.set_quantity('elevation',
    171 #                    filename = project.onshore_dem_name + '.pts',
    172179                    filename = project.combined_dem_name_2 + '.pts',
    173 #                    filename = project.offshore_dem_name + '.pts',
    174180                    use_cache = True,
    175181                    verbose = True,
     
    177183                    )
    178184
     185
    179186#-------------------------------------------------------------------------------                                 
    180187# Setup boundary conditions
     
    185192south = project.south
    186193north = project.north
    187 west = project.west
    188 east = project.east
    189 
    190 #note only need to do when an SWW file for the MOST boundary doesn't exist
    191 cache(ferret2sww,
    192       (source_dir + project.boundary_basename,
    193        source_dir + project.boundary_basename),
    194       {'verbose': True,
    195        'minlat': south,
    196        'maxlat': north,
    197        'minlon': west,
    198        'maxlon': east,
    199 #       'origin': project.mesh_origin,
    200        'origin': domain.geo_reference.get_origin(),
    201        'mean_stage': tide,
    202        'zscale': 1,                 #Enhance tsunami
    203        'fail_on_NaN': False,
    204        'inverted_bathymetry': True},
    205       #evaluate = True,
    206        verbose = True,
    207        dependencies = source_dir + project.boundary_basename + '.sww')
     194west  = project.west
     195east  = project.east
     196
     197###note only need to do when an SWW file for the MOST boundary doesn't exist
     198##cache(ferret2sww,
     199##      (source_dir + project.boundary_basename,
     200##       source_dir + project.boundary_basename),
     201##      {'verbose': True,
     202##       'minlat': south,
     203##       'maxlat': north,
     204##       'minlon': west,
     205##       'maxlon': east,
     206###       'origin': project.mesh_origin,
     207##       'origin': domain.geo_reference.get_origin(),
     208##       'mean_stage': tide,
     209##       'zscale': 1,                 #Enhance tsunami
     210##       'fail_on_NaN': False,
     211##       'inverted_bathymetry': True},
     212##      #evaluate = True,
     213##       verbose = True,
     214##       dependencies = source_dir + project.boundary_basename + '.sww')
    208215
    209216
     
    220227                   f=lambda t: [(60<t<480)*10, 0, 0])
    221228
    222 # for MOST BC
    223 #domain.set_boundary( {'top': Bd, 'left': Bd,
    224 #                      'bottom': Bf, 'right': Bf} )
    225 
    226 # for testing
    227 #domain.set_boundary( {'topr': Bd, 'left': Bd, 'top': Bd,
    228 #                      'bottom': Bw, 'bright': Bd} )
    229 domain.set_boundary( {'e0': Bd, 'e1': Bd, 'e2': Bd, 'e3': Bd, 'e4': Bd,
    230                         'e5': Bd, 'e6': Bd, 'e7': Bd, 'e8': Bd, 'e9': Bd,
    231                         'e10': Bd, 'e11': Bd, 'e12': Bf, 'e13': Bf, 'e14': Bf,
    232                         'e15': Bf} )
     229domain.set_boundary( {'e0': Bd,  'e1': Bd, 'e2': Bd, 'e3': Bd, 'e4': Bd,
     230                      'e5': Bd,  'e6': Bd, 'e7': Bd, 'e8': Bd, 'e9': Bd,
     231                      'e10': Bd, 'e11': Bd, 'e12': Bf, 'e13': Bf, 'e14': Bf,
     232                      'e15': Bf} )
     233
    233234
    234235#-------------------------------------------------------------------------------                                 
     
    240241for t in domain.evolve(yieldstep = 240, finaltime = 6800):
    241242    domain.write_time()
    242     domain.write_boundary_statistics(tags = 'bottom')     
     243    domain.write_boundary_statistics(tags = 'e14')     
    243244
    244245for t in domain.evolve(yieldstep = 30, finaltime = 9000
    245246                       ,skip_initial_step = True):
    246247    domain.write_time()
    247     domain.write_boundary_statistics(tags = 'bottom')     
     248    domain.write_boundary_statistics(tags = 'e14')     
    248249
    249250for t in domain.evolve(yieldstep = 240, finaltime = 15000
    250251                       ,skip_initial_step = True):
    251252    domain.write_time()
    252     domain.write_boundary_statistics(tags = 'bottom')
     253    domain.write_boundary_statistics(tags = 'e14')
    253254   
    254255print 'That took %.2f seconds' %(time.time()-t0)
Note: See TracChangeset for help on using the changeset viewer.