Changeset 3721
- Timestamp:
- Oct 10, 2006, 11:36:20 AM (18 years ago)
- Location:
- anuga_work/production/hobart_2006
- Files:
-
- 16 added
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
anuga_work/production/hobart_2006/export_results.py
r3680 r3721 6 6 from os import sep 7 7 8 # OLE - this is the smallest sww file with the complex bounding polygon 9 time_dir = '20060928_073318' #alpha = 0.5 8 time_dir = '20061006_062319' # Mw 8-5 10 9 11 #time_dir = '20060928_064732' #alpha = 0.0112 #time_dir = '20060928_220905' #alpha = 1.013 #time_dir = '20060926_064750' #this worked for cellsize 30 and 10014 #time_dir = '20060929_033100' #alpha = 0.01 poly anticlockwise - DIDN'T WORK15 #time_dir = '20060929_075954' #alpha = 0.1 grid with MOST16 #time_dir = '20060929_075009' #alpha = 0.1 points with MOST17 10 directory = project.outputdir 18 11 name = directory + time_dir +sep + 'source' 19 12 20 13 print 'output dir:', name 21 #print 'Which variable do you want to export?'22 #which_var = int(raw_input('Stage = 0, Absolute Momentum = 1, Depth = 2, Speed = 3 ' ))23 14 which_var = 4 24 15 if which_var == 0: # Stage … … 31 22 32 23 if which_var == 2: # Depth 33 outname = name + '_depth '24 outname = name + '_depth_bruny' 34 25 quantityname = 'stage-elevation' #Depth 35 26 … … 39 30 40 31 if which_var == 4: # Elevation 41 outname = name + '_elevation _2'32 outname = name + '_elevation' 42 33 quantityname = 'elevation' #Elevation 43 34 44 35 print 'start sww2dem' 36 # for hobart 37 ##sww2dem(name, basename_out = outname, 38 ## quantity = quantityname, 39 ## cellsize = 25, # would prefer this at 25 40 ## # define region for viz purposes 41 ## easting_min = project.eastingmin25, 42 ## easting_max = project.eastingmax25, 43 ## northing_min = project.northingmin25, 44 ## northing_max = project.northingmax25, 45 ## reduction = max, #this is because we want max quantityname 46 ## verbose = True, 47 ## format = 'asc') 48 49 # for bruny 45 50 sww2dem(name, basename_out = outname, 46 51 quantity = quantityname, 47 cellsize = 500, # would prefer this at 2552 cellsize = 25, # would prefer this at 25 48 53 # define region for viz purposes 49 #easting_min = project.e_min_area,50 #easting_max = project.e_max_area,51 #northing_min = project.n_min_area,52 #northing_max = project.n_max_area,54 easting_min = project.eastingmin25_2, 55 easting_max = project.eastingmax25_2, 56 northing_min = project.northingmin25_2, 57 northing_max = project.northingmax25_2, 53 58 reduction = max, #this is because we want max quantityname 54 59 verbose = True, 55 60 format = 'asc') 56 -
anuga_work/production/hobart_2006/get_timeseries.py
r3687 r3721 23 23 mkdir (reportdir) 24 24 25 # data26 production_dirs = {'20061003_062209': 'test'} 27 #'20061003_062209': 'Mw 8.7'}#, 28 #' ': 'Mw 8.5'}25 # 2500 interior resolution for both events 26 # use both to compare the events 27 production_dirs = {'20061006_062319': 'Mw 8.5'} 28 #'20061008_234702': 'Mw 8.7'} 29 29 30 30 # Generate figures -
anuga_work/production/hobart_2006/project.py
r3701 r3721 7 7 import sys 8 8 from time import localtime, strftime, gmtime 9 from anuga.utilities.polygon import read_polygon, plot_polygons 9 from anuga.utilities.polygon import read_polygon, plot_polygons, polygon_area 10 10 from anuga.coordinate_transforms.redfearn import degminsec2decimal_degrees, convert_points_from_latlon_to_utm 11 11 … … 26 26 scenario_dir_name = 'hobart_tsunami_scenario_2006' 27 27 28 # data provided by Tas SES and checked by NM&I 29 #onshore_name = 'hob3_topo' # original 30 onshore_name_25 = 'hob5_topo_25m' # 25m grid and clipped to 100m elevation or 3000m from coast 31 onshore_name = 'hob5_topo' # 12.5m 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' 28 37 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' 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' 60 55 61 boundary_basename = 'puysegur' # Mw 8.7 62 #boundary_basename = 'puysegur_clip' # Mw 8.5 56 # developed by NM&I 57 coast_name = 'coastline_points' 63 58 59 #boundary_basename = 'puysegur' # Mw 8.7 60 boundary_basename = 'puysegur_clip' # Mw 8.5 64 61 65 62 #swollen/ all data output … … 81 78 82 79 gauge_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'86 87 # boundary source data88 #MOST_dir = 'f:'+sep+'3'+sep+'ehn'+sep+'users'+sep+'davidb'+sep+'tsunami'+sep+'WA_project'+sep+'SU-AU_90'+sep+'most_2'+sep+'detailed'+sep89 80 90 81 codedir = getcwd()+sep … … 93 84 94 85 # 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 86 onshore_dem_name = datadir + onshore_name 87 onshore_dem_name_25 = datadir + onshore_name_25 88 all_onshore_dem_name = datadir + 'combined_onshore' 89 offshore_dem_name_local1 = datadir + offshore_name_tas1 90 offshore_dem_name_local2 = datadir + offshore_name_tas2 91 offshore_dem_name_local3 = datadir + offshore_name_tas3 92 offshore_dem_name_local4 = datadir + offshore_name_tas4 93 offshore_dem_name_aho1 = datadir + offshore_name1 94 offshore_dem_name_aho2 = datadir + offshore_name2 95 offshore_dem_name_aho3 = datadir + offshore_name3 96 offshore_dem_name_aho4 = datadir + offshore_name4 97 offshore_dem_name_aho5 = datadir + offshore_name5 98 offshore_dem_name_aho6 = datadir + offshore_name6 99 offshore_dem_name_aho7 = datadir + offshore_name7 100 offshore_dem_name_aho8 = datadir + offshore_name8 101 offshore_dem_name_aho9 = datadir + offshore_name9 102 offshore_dem_name_aho10 = datadir + offshore_name10 103 offshore_dem_name_aho11 = datadir + offshore_name11 104 offshore_dem_name_aho12 = datadir + offshore_name12 105 offshore_dem_name_aho13 = datadir + offshore_name13 106 offshore_dem_name_aho14 = datadir + offshore_name14 107 offshore_dem_name_aho15 = datadir + offshore_name15 108 offshore_dem_name_aho16 = datadir + offshore_name16 109 coast_dem_name = datadir + coast_name 120 110 121 111 # addition once total grid delivered … … 128 118 combined_dem_name = datadir + 'hobart_combined_elevation' 129 119 combined_dem_name_2 = datadir + 'hobart_combined_elevation_2' 120 combined_dem_name_3 = datadir + 'hobart_combined_elevation_jane' 121 combined_dem_name_4 = datadir + 'hobart_combined_elevation_test' 130 122 131 123 #outputname = outputtimedir + basename #Used by post processing … … 150 142 refzone = zone 151 143 152 # Second run - bottom bright, topr, top, left153 #polyAll = [[520000, 5170000],[580000, 5170000],[580000, 5200000],[590000,5240000],[520000,5260000]]154 155 156 144 # Mark run - morning Fri 29 Sep; surrounds -100m and 20mish elevation 157 polyAll = read_polygon(polygondir+'new_extent .csv')145 polyAll = read_polygon(polygondir+'new_extent_2.csv') 158 146 plot_polygons([polyAll, polyAll2],'boundingpoly',verbose=False) 159 147 print 'Area of bounding polygon', polygon_area(polyAll) 160 148 161 149 ################################################################### … … 163 151 ################################################################### 164 152 153 # clipping 12.5m onshore data set 154 eastingmin = 520000 155 eastingmax = 536000 156 northingmin = 5245000 157 northingmax = 5260000 165 158 166 159 # clipping 25m data set - Hobart … … 176 169 northingmax25_2 = 5212052.309 177 170 178 179 180 171 ############################### 181 172 # Interior region definitions 182 173 ############################### 183 184 ###Interior region - Hobart city area + Glenorchy, Kingston185 ##i0 = [517000, 5267000]186 ##i1 = [517000, 5255000]187 ##i2 = [520000, 5250000]188 ##i3 = [522000, 5239000]189 ##i4 = [524000, 5238000]190 ##i5 = [526000, 5236000]191 ##i6 = [530000, 5244000]192 ##i7 = [530000, 5250000]193 ##i8 = [534000, 5254000]194 ##i9 = [520000, 5270000]195 ##196 ##poly_hobart = [i0, i1, i2, i3, i4, i5, i6, i7, i8, i9]197 ##198 ### Tasman Peninsula199 ##l0 = [550000, 5247000]200 ##l1 = [550000, 5211000]201 ##l2 = [583000, 5211000]202 ##l3 = [583000, 5247000]203 ##204 ##poly_tasman_peninsula = [l0, l1, l2, l3]205 ##206 ### Bruny Island207 ##poly_bruny = read_polygon(polygondir+'bruny.csv')208 209 174 210 175 # Hobart digitized polygons … … 212 177 poly_hobart2 = read_polygon(polygondir+'Hob_poly2.csv') 213 178 poly_hobart3 = read_polygon(polygondir+'Hob_poly3.csv') 179 poly_hobart4 = read_polygon(polygondir+'Hob_poly4.csv') 214 180 215 plot_polygons([polyAll, poly_hobart1,poly_hobart2,poly_hobart3 ],'boundingpoly2',verbose=False)181 plot_polygons([polyAll, poly_hobart1,poly_hobart2,poly_hobart3,poly_hobart4],'boundingpoly2',verbose=False) -
anuga_work/production/hobart_2006/run_hobart.py
r3701 r3721 73 73 # creates pts file for combined 50m DEM 74 74 dem2pts(onshore_offshore_dem_name, use_cache=True, verbose=True) 75 75 """ 76 76 # 25m data (clipping the around the Hobart area) 77 77 convert_dem_from_ascii2netcdf(onshore_offshore_dem_name_25, use_cache=True, verbose=True) … … 113 113 #G.export_points_file(project.combined_dem_name + '.pts') 114 114 115 115 """ 116 116 #---------------------------------------------------------------------------- 117 117 # Create the triangular mesh based on overall clipping polygon with a tagged … … 146 146 'e12': [12], 'e13': [13], 'e14': [14], 147 147 'e15': [15]}, 148 'maximum_triangle_area': 750000,149 'filename': meshname ,150 'interior_regions': interior_regions},148 'maximum_triangle_area': 2500000, 149 'filename': meshname}, 150 #'interior_regions': interior_regions}, 151 151 verbose = True, evaluate=False) 152 152 … … 176 176 domain.set_quantity('friction', 0.0) 177 177 178 domain.set_quantity('elevation', 179 filename = project.combined_dem_name_2 + '.pts', 180 use_cache = True, 181 verbose = True, 182 alpha = 0.1 183 ) 178 #domain.set_quantity('elevation', 179 # filename = project.combined_dem_name_2 + '.pts', 180 # filename = onshore_offshore_dem_name + '.pts', 181 # use_cache = True, 182 # verbose = True, 183 # alpha = 0.1 184 # ) 184 185 185 186 … … 188 189 #------------------------------------------------------------------------------- 189 190 print 'start ferret2sww' 191 print '', project.boundary_basename 190 192 from anuga.shallow_water.data_manager import ferret2sww 191 193 … … 195 197 east = project.east 196 198 197 # ##note only need to do when an SWW file for the MOST boundary doesn't exist198 ##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 tsunami210 ##'fail_on_NaN': False,211 ##'inverted_bathymetry': True},212 ###evaluate = True,213 ##verbose = True,214 ##dependencies = source_dir + project.boundary_basename + '.sww')199 #note only need to do when an SWW file for the MOST boundary doesn't exist 200 cache(ferret2sww, 201 (source_dir + project.boundary_basename, 202 source_dir + project.boundary_basename), 203 {'verbose': True, 204 'minlat': south, 205 'maxlat': north, 206 'minlon': west, 207 'maxlon': east, 208 # 'origin': project.mesh_origin, 209 'origin': domain.geo_reference.get_origin(), 210 'mean_stage': tide, 211 'zscale': 1, #Enhance tsunami 212 'fail_on_NaN': False, 213 'inverted_bathymetry': True}, 214 #evaluate = True, 215 verbose = True, 216 dependencies = source_dir + project.boundary_basename + '.sww') 215 217 216 218
Note: See TracChangeset
for help on using the changeset viewer.