Changeset 6288

Timestamp:

Feb 6, 2009, 9:04:49 AM (16 years ago)

Author:

jgriffin

Message:

new format scripts added, old ones deleted

Location:

anuga_work/production/australia_ph2/sydney

Files:

• Arc_asc2raster_GDA94z56.py (modified) (2 diffs)
• build_boundary_7875.py (modified) (1 diff)
• build_model.py (modified) (1 diff)
• build_sydney_250m.py (deleted)
• export_results_max.py (modified) (3 diffs)
• run_sydney_250m.py (modified) (5 diffs)

anuga_work/production/australia_ph2/sydney/Arc_asc2raster_GDA94z56.py

@@ -31 +31 @@
 ##time_dir1 = '20081031_133353_run_final_0.6_68693_alpha0.1_kvanputt'
 ##time_dir2 = '20081031_133511_run_final_0_68693_alpha0.1_kvanputt'
-time_dir3 = '20090123_155045_run_final_0_7875_250m_large_jgriffin'
-time_dir4 = '20090123_154945_run_final_0_7875_250m_small_jgriffin'
+time_dir3 = '20090130_160955_run_final_0_7875_250m_vlarge_jgriffin'
+#time_dir4 = '20090123_154945_run_final_0_7875_250m_small_jgriffin'
 ##time_dir5 = '20081031_133841_run_final_0_27283_alpha0.1_kvanputt'
 ##time_dir6 = '20081031_133925_run_final_0.6_27283_alpha0.1_kvanputt'
 
 
-time_dirs = [time_dir3, time_dir4] #, time_dir3, time_dir4, time_dir5, time_dir6]
+time_dirs = [time_dir3]#, time_dir4] #, time_dir3, time_dir4, time_dir5, time_dir6]
 
 for time_dir in time_dirs:
@@ -62 +62 @@
 
     generate_filename = []
-    input_ascii = glob.glob(folder + 'sydney_All_stage.asc')
+    input_ascii = glob.glob(folder + '*stage_max.asc')
     print time_dir
 

anuga_work/production/australia_ph2/sydney/build_boundary_7875.py

@@ -64 +64 @@
 
     scenario_name=project_250m.scenario_name
-    order_filename=os.path.join(project_250m.order_filename_dir)
+    order_filename=os.path.join(project_250m.urs_boundary_name)
 
     print 'reading', order_filename

anuga_work/production/australia_ph2/sydney/build_model.py

@@ -89 +89 @@
 
 print 'Export combined DEM file'
-G_clip.export_points_file(project.combined_dir_name + '.pts')
+G_clip.export_points_file(project.combined_elevation + '.pts')
 print 'Do txt version too'
 # Use for comparision in ARC
-G_clip.export_points_file(project.combined_dir_name + '.txt')
+G_clip.export_points_file(project.combined_elevation + '.txt')
 

anuga_work/production/australia_ph2/sydney/export_results_max.py

@@ -20 +20 @@
 directory = project_250m.output_dir
 
-time_dir1 = '20081125_135917_run_final_0_27283_250m_none_lfountai'
-
-time_dirs = [time_dir1]#, time_dir2, time_dir3, time_dir4, time_dir5, time_dir6]
+time_dir1 = '20090123_155045_run_final_0_7875_250m_large_jgriffin'
+#time_dir2 = '20090123_155045_run_final_0_7875_250m_large_jgriffin'
+#time_dir3 = '20090123_154945_run_final_0_7875_250m_small_jgriffin'
+time_dirs = [time_dir1]#,time_dir2,time_dir3]#, time_dir2, time_dir3, time_dir4, time_dir5, time_dir6]
 
 #cellsize = 20
@@ -48 +49 @@
 
 # one or more key strings from var_equations above
-var = ['stage', 'speed']
+#var = ['stage', 'speed','depth']
+var = ['elevation']
 
 ######
@@ -78 +80 @@
 
             name1 = directory+time_dir+sep+project_250m.scenario_name
-            name2 = directory+time_dir+sep+project_250m.scenario_name+'_time_59580_0'
+            #name2 = directory+time_dir+sep+project_250m.scenario_name+'_time_52620_0'
            # name3 = directory+time_dir+sep+project_250m.scenario_name+'_time_78360_0'
 
-            names = [name1, name2] #, name3]
+            names = [name1]#, name2] #, name3]
       
             asc_name = []   

anuga_work/production/australia_ph2/sydney/run_sydney_250m.py

@@ -26 +26 @@
 from os import mkdir, access, F_OK
 from shutil import copy
+from math import pi, sin, exp
 import time
 import sys
@@ -31 +32 @@
 # Related major packages
 from anuga.shallow_water import Domain
+from anuga.shallow_water.shallow_water_domain import Transmissive_stage_zero_momentum_boundary
 from anuga.shallow_water import Dirichlet_boundary
 from anuga.shallow_water import File_boundary
 from anuga.shallow_water import Reflective_boundary
 from anuga.shallow_water import Field_boundary
+from anuga.shallow_water import Time_boundary
 from Numeric import allclose
 from anuga.shallow_water.data_manager import export_grid, create_sts_boundary
@@ -76 +79 @@
     # Domain definitions
     #-----------------------------------------------------------------------
-
-    # Read in boundary from ordered sts file
-    urs_bounding_polygon=create_sts_boundary(os.path.join(project_250m.boundaries_dir_event,project_250m.scenario_name))
-
-    # Reading the landward defined points, this incorporates the original clipping
-    # polygon minus the 100m contour
-    landward_bounding_polygon = read_polygon(project_250m.landward_dir)
-
-    # Combine sts polyline with landward points
-    bounding_polygon = urs_bounding_polygon + landward_bounding_polygon
-    
-    # counting segments
-    N = len(urs_bounding_polygon)-1
-
-    # boundary tags refer to project_250m.landward 4 points equals 5 segments start at N
-    boundary_tags={'back': [N+1,N+2,N+3,N+4,N+5,N+6,N+7,N+8,N+9,N+10,N+11,N+12,N+13,N+14,N+15,N+16,N+17], 'side': [N,N+18], 'ocean': range(N)}
-
+##
+##    # Read in boundary from ordered sts file
+##    urs_bounding_polygon=create_sts_boundary(os.path.join(project_250m.boundaries_dir_event,project_250m.scenario_name))
+##
+##    # Reading the landward defined points, this incorporates the original clipping
+##    # polygon minus the 100m contour
+##    landward_bounding_polygon = read_polygon(project_250m.landward_dir)
+##
+##    # Combine sts polyline with landward points
+##    bounding_polygon = urs_bounding_polygon + landward_bounding_polygon
+##    
+##    # counting segments
+##    N = len(urs_bounding_polygon)-1
+##    print 'N'+str(N)
+##    # boundary tags refer to project_250m.landward 4 points equals 5 segments start at N
+##    boundary_tags={'back': [N+1,N+2,N+3,N+4,N+5,N+6,N+7,N+8,N+9,N+10,N+11,N+12,N+13,N+14,N+15,N+16,N+17,N+18,N+19], 'side': [N,N+20], 'ocean': range(N)}
+    if project_250m.area =='large':
+        land = [2,3,4,]
+        sea = [0]
+        side1 = [1]
+        side2 = [5]
+        bounding_polygon = project_250m.poly_all
+        boundary_tags={'back': land, 'ocean': sea, 'side1':side1,'side2':side2}
+    elif project_250m.area =='small':
+        land = [2,3,4,]
+        sea = [0]
+        side1 = [1]
+        side2 = [5]
+        bounding_polygon = project_250m.poly_all
+        boundary_tags={'back': land, 'ocean': sea, 'side1':side1,'side2':side2}
+    elif project_250m.area =='medium':
+        land = [2,3,4]
+        sea = [0]
+        side1 = [1]
+        side2 = [5]
+        bounding_polygon = project_250m.poly_all
+        boundary_tags={'back': land, 'ocean': sea, 'side1':side1,'side2':side2}
+    elif project_250m.area =='vlarge':
+        land = [9,10,11,12,13,14]
+        sea = [0,1,2,3,4,5,6,7]
+        side1 = [8]
+        side2 = [15]
+        bounding_polygon = project_250m.poly_all
+        boundary_tags={'back': land, 'ocean': sea, 'side1':side1,'side2':side2} 
+    else:
+        print 'area not defined therefore bondary_tags not defined'
     #--------------------------------------------------------------------------
     # Create the triangular mesh based on overall clipping polygon with a tagged
@@ -109 +141 @@
                          interior_regions=project_250m.interior_regions,
                          filename=project_250m.meshes_dir_name,
-                         use_cache=True,
+                         use_cache=False,
                          verbose=True)
    
@@ -177 +209 @@
     Br = Reflective_boundary(domain)
     Bd = Dirichlet_boundary([kwargs['tide'],0,0])
-    Btr = Transmissive_stage_zero_momentum_boundary(domain)
-    
-    print 'Available boundary tags', domain.get_boundary_tags()
-    Bf = Field_boundary(boundary_urs_out+'.sts',  # Change from file_boundary
-                   domain, mean_stage= project_250m.tide,
-                   time_thinning=1,
-                   default_boundary=Bd,
-                   use_cache=True,
-                   verbose = True,
-                   boundary_polygon=bounding_polygon)
+    Bt = Transmissive_stage_zero_momentum_boundary(domain)
+    Rb = Time_boundary (domain, f=lambda t: [(sin(t*pi/600) * exp(-(t/3600)**2)), 0.0, 0.0]) 
+##    print 'Available boundary tags', domain.get_boundary_tags()
+##    Bf = Field_boundary(boundary_urs_out+'.sts',  # Change from file_boundary
+##                   domain, mean_stage= project_250m.tide,
+##                   time_thinning=1,
+##                   default_boundary=Bd,
+##                   use_cache=False,
+##                   verbose = True,
+##                   boundary_polygon=bounding_polygon)
 
     domain.set_boundary({'back': Br,
-                         'side': Btr,
-                         'ocean': Bf}) 
+                         'side1': Bt,
+                         'side2':Bt,
+                         'ocean':Rb}) 
 
     kwargs['input_start_time']=domain.starttime