Changeset 3615

Timestamp:

Sep 18, 2006, 5:19:25 PM (18 years ago)

Author:

sexton

Message:

updates to Hobart script(data not yet received)

Location:

anuga_work/production/hobart_2006

Files:

• export_results.py (added)
• project.py (modified) (4 diffs)
• run_hobart.py (modified) (9 diffs)

anuga_work/production/hobart_2006/project.py

@@ -2 +2 @@
 """Common filenames and locations for topographic data, meshes and outputs.
 """
-
-
 
 from os import sep, environ, getenv, getcwd
@@ -67 +65 @@
 
 onshore_dem_name = datadir + onshore_name
-offshore_dem_name = datadir + offshore_name
+offshore_dem_name_local = datadir + offshore_name_local
+offshore_dem_name_aho = datadir + offshore_name_fairsheets
 coast_dem_name = datadir + coast_name
 combined_dem_name = datadir + 'hobart_combined_elevation'
@@ -73 +72 @@
 outputname = outputtimedir + basename  #Used by post processing
 
-from anuga.coordinate_transforms.redfearn import degminsec2decimal_degrees
-from anuga.coordinate_transforms.redfearn import redfearn
+# region to export
+e_min_area = 500000
+e_max_area = 540000
+n_min_area = 5220000
+n_max_area = 5275000
 
+from anuga.coordinate_transforms.redfearn import degminsec2decimal_degrees,
+                        convert_points_from_latlon_to_utm
+                        
 # bounding box
 south = degminsec2decimal_degrees(-43,45,0)
@@ -82 +87 @@
 east = degminsec2decimal_degrees(148,15,0)
 
-'''
-# region for visualisation
-eminviz = 260000
-emaxviz = 320000
-nminviz = 7590000
-nmaxviz = 7630000
-'''
-# region to export
+#Main Domain of Hobart: first run JS 18/9/06
+d0 = [south, west]
+d1 = [south, east]
+d2 = [north, east]
+d3 = [north, west]
 
-e_min_area = 300000
-e_max_area = 310000
-n_min_area = 7600000
-n_max_area = 7610000
+polyAll, zone = convert_points_from_latlon_to_utm([d0, d1, d2, d3])
+refzone = zone
 
-
-#Main Domain of Hobart: first run JS 15/9/06
-zone, e0, n0 = redfearn(south, west)
-zone, e1, n1 = redfearn(south, east)
-zone, e2, n2 = redfearn(north, east)
-zone, e3, n3 = redfearn(north, west)
-refzone = zone
-d0 = [e0, n0]
-d1 = [e1, n1]
-d2 = [e2, n2]
-d3 = [e3, n3]
-
-polyAll = [d0, d1, d2, d3]
-
-#Interior region - Hobart city area
-#from anuga.utilities.polygon import read_polygon
-polygonptsfile0 = polygondir + 'poly1’
-polygonptsfile1 = polygondir + 'poly2’
-northern_polygon = read_polygon(polygonptsfile0 + '.csv')
-southern_polygon = read_polygon(polygonptsfile1 + '.csv')
-
-i0 = [304000, 7607000]
-i1 = [302000, 7605000]
-i2 = [304000, 7603000]
-i3 = [307000, 7602000]
-i4 = [309000, 7603000]
-i5 = [307000, 7606000]
+#Interior region - Hobart city area + glenorchy, Kingston
+i0 = [517000, 5267000]
+i1 = [517000, 5255000]
+i2 = [520000, 5250000]
+i3 = [522000, 5239000]
+i4 = [524000, 5238000]
+i5 = [526000, 5236000]
+i6 = [530000, 5244000]
+i7 = [530000, 5250000]
+i7 = [534000, 5254000]
+i7 = [520000, 5270000]
 
 poly_hobart = [i0, i1, i2, i3, i4, i5]
 
-#med res around Hobart
-l0 = [300000, 7610000]
-l1 = [285000, 7600000]
-l2 = [300000, 7597500]
-l3 = [310000, 7600000]
-l4 = [315000, 7610000]
+# Tasman Peninsula
+l0 = [550000, 5247000]
+l1 = [550000, 5211000]
+l2 = [583000, 5211000]
+l3 = [583000, 5247000]
 
-poly_coast = [l0, l1, l2, l3, l4]
-
-#general coast and local area to Hobart region
-m0 = [270000, 7581000]
-m1 = [300000, 7591000]
-m2 = [339000, 7610000]
-m3 = [330000, 7630000]
-m4 = [290000, 7640000]
-m5 = [260000, 7600000]
-
-poly_region = [m0, m1, m2, m3, m4, m5]
+poly_tasman_peninsula = [l0, l1, l2, l3]
 
 # Bruny Island
-poly_bruny = [b0, b1, b2, b3, b4]
+from anuga.utilities.polygon import read_polygon
+poly_bruny = read_polygon(polygondir+'bruny.csv')

anuga_work/production/hobart_2006/run_hobart.py

@@ -6 +6 @@
 
 The scenario is defined by a triangular mesh created from project.polygon,
-the elevation data and a simulated submarine landslide.
+the elevation data and a tsunami wave generated by MOST.
 
 Ole Nielsen and Duncan Gray, GA - 2005 and Nick Bartzis, GA - 2006
 """
-
-
 #-------------------------------------------------------------------------------# Import necessary modules
 #-------------------------------------------------------------------------------
@@ -67 +65 @@
 # filenames
 onshore_dem_name = project.onshore_dem_name
-coast_points = project.coast_dem_name
-offshore_points = project.offshore_dem_name
 meshname = project.meshname+'.msh'
 source_dir = project.boundarydir
@@ -81 +77 @@
 
 #creates pts file for onshore DEM
-dem2pts(onshore_dem_name,
-        easting_min=project.eastingmin,
-        easting_max=project.eastingmax,
-        northing_min=project.northingmin,
-        northing_max= project.northingmax,
-        use_cache=True, 
-        verbose=True)
-
-convert_dem_from_ascii2netcdf(islands_dem_name, use_cache=True, verbose=True)
-
-#creates pts file for islands DEM
-dem2pts(islands_dem_name, use_cache=True, verbose=True)
+dem2pts(onshore_dem_name, use_cache=True, verbose=True)
 
 print'create G1'
-G1 = Geospatial_data(file_name = project.offshore_dem_name + '.xya')
+G1 = Geospatial_data(file_name = project.offshore_dem_name_local + '.xya')
 print'create G2'
-G2 = Geospatial_data(file_name = project.onshore_dem_name + '.pts')
+G2 = Geospatial_data(file_name = project.offshore_dem_name_aho + '.xya')
 print'create G3'
-G3 = Geospatial_data(file_name = project.coast_dem_name + '.xya')
-print'add G1+G2+G3'
-G = G1 + G2 + G3
+G3 = Geospatial_data(file_name = project.onshore_dem_name + '.pts')
+print'create G4'
+G4 = Geospatial_data(file_name = project.coast_dem_name + '.xya')
+print'add G1+G2+G3+G4'
+G = G1 + G2 + G3 + G4
 print'export G'
 G.export_points_file(project.combined_dem_name + '.pts')
@@ -113 +100 @@
 from anuga.pmesh.mesh_interface import create_mesh_from_regions
 
-#new
-region_res = 200000
-coast_res = 25000
-hobart_res = 5000
+# use 75 for onshore components (12.5m DEM)
+island_res = 10000
+hobart_res = 10000
+peninsular_res = 10000
 interior_regions = [[project.poly_hobart, hobart_res],
-                    [project.poly_coast, coast_res],
-                    [project.poly_region, region_res]]
+                    [project.poly_tasman_peninsula, peninsula_res],
+                    [project.poly_bruny, island_res]]
 
 print 'number of interior regions', len(interior_regions)
@@ -126 +113 @@
 _ = cache(create_mesh_from_regions,
           project.polyAll,
-          {'boundary_tags': {'top': [0], 'topleft': [1],
-                             'topleft1': [2], 'bottomleft': [3],
-                             'bottom': [4], 'bottomright': [5],
-                             'topright':[6]},
+          {'boundary_tags': {'bottom': [0], 'right': [1],
+                             'top': [2], 'left': [3]},
            'maximum_triangle_area': 100000,
            'filename': meshname,           
@@ -177 +162 @@
 print 'start ferret2sww'
 from anuga.pyvolution.data_manager import ferret2sww
-
+'''
 south = project.south 
 north = project.north
@@ -209 +194 @@
        dependencies = source_dir + project.boundary_basename + '.sww')
 
-
+'''
 print 'Available boundary tags', domain.get_boundary_tags()
 
-Bf = File_boundary(source_dir + project.boundary_basename + '.sww', 
+#Bf = File_boundary(source_dir + project.boundary_basename + '.sww', 
                     domain, verbose = True)
 Br = Reflective_boundary(domain)
@@ -222 +207 @@
                    f=lambda t: [(60<t<480)*6, 0, 0])
 
-domain.set_boundary( {'top': Bf, 'topleft': Bf,
-                      'topleft1': Bf, 'bottomleft': Bd,
-                      'bottom': Br, 'bottomright': Br, 'topright': Bd} )
+# for MOST BC
+#domain.set_boundary( {'top': Bd, 'left': Bd,
+#                      'bottom': Bf, 'right': Bf} )
+
+# for testing
+domain.set_boundary( {'top': Bd, 'left': Bd,
+                      'bottom': Bd, 'right': Bw} )
 
 #-------------------------------------------------------------------------------                                 
@@ -234 +223 @@
 for t in domain.evolve(yieldstep = 240, finaltime = 7200): 
     domain.write_time()
-    domain.write_boundary_statistics(tags = 'top')     
+    domain.write_boundary_statistics(tags = 'bottom')     
 
 for t in domain.evolve(yieldstep = 120, finaltime = 12600
                        ,skip_initial_step = True): 
     domain.write_time()
-    domain.write_boundary_statistics(tags = 'top')     
-
-for t in domain.evolve(yieldstep = 60, finaltime = 19800
-                       ,skip_initial_step = True): 
-    domain.write_time()
-    domain.write_boundary_statistics(tags = 'top')     
-    
-for t in domain.evolve(yieldstep = 120, finaltime = 25200
-                       ,skip_initial_step = True): 
-    domain.write_time()
-    domain.write_boundary_statistics(tags = 'top')     
-
-for t in domain.evolve(yieldstep = 240, finaltime = 36000
-                       ,skip_initial_step = True): 
-    domain.write_time()
-    domain.write_boundary_statistics(tags = 'top')
+    domain.write_boundary_statistics(tags = 'bottom')     
     
 print 'That took %.2f seconds' %(time.time()-t0)