Changeset 5569

Timestamp:

Jul 25, 2008, 8:53:51 AM (17 years ago)

Author:

kristy

Message:

updated scripts for polyline boundary

Location:

anuga_work/production/perth

Files:

• build_boundary.py (modified) (2 diffs)
• export_results.py (modified) (1 diff)
• get_timeseries.py (modified) (1 diff)
• project.py (modified) (4 diffs)
• run_perth.py (modified) (6 diffs)

anuga_work/production/perth/build_boundary.py

@@ -34 +34 @@
 urs_filenames = [os.path.join(dir, basename) for basename in urs_filenames.keys()]
 
-# Reduce this list slightly
-urs_filenames = urs_filenames[:30]
-for name in urs_filenames:
-    print name
+# AS per David Burbidge email on friday 4th July the mag 9.3 event
+# has 1m worth of slip on each sub fault therefore mutliple each unit
+# soucre by the slip (10.4544) and sum the 44 time series together to
+# get the time series for this event at the points on your boundary. 
 
 weights=10.4544*ones(len(urs_filenames),Float)
 
-base_name=project.scenario_name
-order_filename=os.path.join(project.boundaries_dir, 'boundary_ordering.txt')
+scenario_name=project.scenario_name
+order_filename=os.path.join(project.boundaries_dir, 'thinned_bound_order.txt')
 
 print 'reading', order_filename
 # Create ordered sts file
 print 'creating sts file'
-urs2sts(urs_filenames,basename_out=base_name,
+
+urs2sts(urs_filenames,basename_out=scenario_name,
         ordering_filename=order_filename,
         weights=weights,
@@ -54 +55 @@
 
 # Read in boundary from ordered sts file
-urs_boundary=create_sts_boundary(base_name)
+urs_boundary_polygon=create_sts_boundary(scenario_name)
 
+
+
+
+
+
+

anuga_work/production/perth/export_results.py

@@ -8 +8 @@
 #time_dir = '20080530_170833_run_final_0.6_exmouth_kvanputt'
 #time_dir = '20080619_115643_run_trial_0.6_exmouth_kvanputt'
-time_dir = '20080714_095738_run_final_0.6_exmouth_kvanputt'
-
+#time_dir = '20080714_095738_run_final_0.6_exmouth_kvanputt'
+#time_dir = '20080714_134935_run_final_0_exmouth_kvanputt'
+time_dir = '20080718_112228_run_final_0_exmouth_kvanputt'
 
 cellsize = 25

anuga_work/production/perth/get_timeseries.py

@@ -16 +16 @@
 
 
-timestamp='20080619_115643_run_trial_0.6_exmouth_kvanputt'
+timestamp='20080724_121200_run_trial_0.6_polyline_alpha0.1_kvanputt'
 
 filename=project.output_dir+timestamp+sep+project.scenario_name+'.sww'

anuga_work/production/perth/project.py

@@ -10 +10 @@
 #from anuga.coordinate_transforms.redfearn import degminsec2decimal_degrees, convert_points_from_latlon_to_utm
 from anuga.utilities.system_tools import get_user_name, get_host_name
+from anuga.shallow_water.data_manager import urs2sts,create_sts_boundary
+from anuga.utilities.polygon import read_polygon, plot_polygons, polygon_area, is_inside_polygon
 
 # file and system info
@@ -36 +38 @@
 tide = 0.6
 
-alpha = 0.2
+alpha = 0.1
 friction=0.01
 starttime=0
 finaltime=80000 
 export_cellsize=25
-setup='trial'
+setup='final'
 source='polyline'
 
@@ -160 +162 @@
 buildings_filename_out = 'Perth_res_Project_modified.csv'
 
-###############################
-# Domain definitions
-###############################
-from anuga.utilities.polygon import read_polygon, plot_polygons, polygon_area, is_inside_polygon
-
-poly_all = read_polygon(polygons_dir+'poly_all.csv')
-res_poly_all = 100000*res_factor
-
-#refzone = 50 
+
 
 ###############################
 # Interior region definitions
 ###############################
+
+#Initial bounding polygon for data clipping 
+poly_all = read_polygon(polygons_dir+'poly_all.csv')
+res_poly_all = 100000*res_factor
 
 #Polygon designed by 20m contours, or 3km from the coastline
@@ -217 +215 @@
                      ,[poly_rottnest_ex, res_rottnest_ex], [poly_garden, res_garden]]
 
-boundary_tags={'back': [0,1,2], 'side': [3,7],'ocean': [4, 5, 6]}
-
+    
 trigs_min = number_mesh_triangles(interior_regions, poly_all, res_poly_all)
-
 print 'min number triangles', trigs_min
+   
 
 poly_mainland = read_polygon(polygons_dir+'initial_condition.csv')

anuga_work/production/perth/run_perth.py

@@ -30 +30 @@
 from anuga.shallow_water import Field_boundary
 from Numeric import allclose
-from anuga.shallow_water.data_manager import export_grid
+from anuga.shallow_water.data_manager import export_grid, create_sts_boundary
 
 from anuga.pmesh.mesh_interface import create_mesh_from_regions
@@ -39 +39 @@
 from anuga.damage_modelling.inundation_damage import add_depth_and_momentum2csv, inundation_damage
 from anuga.fit_interpolate.benchmark_least_squares import mem_usage 
+from anuga.utilities.polygon import read_polygon, plot_polygons, polygon_area, is_inside_polygon
 
 # Application specific imports
@@ -70 +71 @@
     print "Processor Name:",get_processor_name()
 
-    # filenames
-#    meshes_dir_name = project.meshes_dir_name+'.msh'
-
-    # creates copy of code in output dir
-    print 'min triangles', project.trigs_min,
-    print 'Note: This is generally about 20% less than the final amount'
-
+    
+    #-----------------------------------------------------------------------
+    # Domain definitions
+    #-----------------------------------------------------------------------
+
+    # Read in boundary from ordered sts file
+    urs_bounding_polygon=create_sts_boundary(project.scenario_name)
+
+    # Reading the landward defined points, this incorporates the original clipping
+    # polygon minus the 100m contour
+    landward_bounding_polygon = read_polygon(project.polygons_dir+'landward_bounding_polygon.txt')
+
+    # Combine sts polyline with landward points
+    bounding_polygon = urs_bounding_polygon + landward_bounding_polygon
+    
+    # counting segments
+    N = len(urs_bounding_polygon)-1
+    boundary_tags={'back': [N+1,N+2,N+3], 'side': [N,N+4],'ocean': range(N)}
+
+    
     #--------------------------------------------------------------------------
     # Create the triangular mesh based on overall clipping polygon with a
@@ -90 +104 @@
         print 'start create mesh from regions'
 
-        create_mesh_from_regions(project.poly_all,
-                             boundary_tags=project.boundary_tags,
+        create_mesh_from_regions(bounding_polygon,
+                             boundary_tags=boundary_tags,
                              maximum_triangle_area=project.res_poly_all,
                              interior_regions=project.interior_regions,
@@ -168 +182 @@
     print 'domain id', id(domain)
    
-
+    boundary_urs_out=project.boundaries_dir_name
+    
     print 'Available boundary tags', domain.get_boundary_tags()
     Bf = File_boundary(boundary_urs_out+'.sts', 
@@ -174 +189 @@
                    use_cache=True,
                    verbose = True,
-                   boundary_polygon=project.bounding_polygon)
+                   boundary_polygon=bounding_polygon)
 
     Br = Reflective_boundary(domain)