Index: /anuga_work/production/hobart_2006/export_results.py
===================================================================
--- /anuga_work/production/hobart_2006/export_results.py	(revision 3670)
+++ /anuga_work/production/hobart_2006/export_results.py	(revision 3671)
@@ -3,9 +3,8 @@
 
 from anuga.shallow_water.data_manager import sww2dem
-#from anuga.pyvolution.ermapper_grids import read_ermapper_grid
 from anuga.abstract_2d_finite_volumes.util import Screen_Catcher
 from os import sep
 
-time_dir = '20060925_115139' #MSL
+time_dir = '20060927_054931' #MSL
 directory = project.outputdir
 name = directory + time_dir +sep + 'source'
@@ -15,5 +14,4 @@
 #which_var = int(raw_input('Stage = 0, Absolute Momentum = 1, Depth = 2, Speed = 3  '  ))
 which_var = 4
-#sys.stderr.write(sys.stdout.data)
 if which_var == 0:  # Stage
     outname = name + '_stage'
@@ -30,17 +28,14 @@
 if which_var == 3:  # Speed
     outname = name + '_speed'
-    #quantityname = '((xmomentum/(stage-elevation))**2 + (ymomentum/(stage-elevation))**2)**0.5'  #Speed
-    #quantityname = 'xmomentum/(stage-elevation)'  #Speed
     quantityname = '(xmomentum**2 + ymomentum**2)**0.5/(stage-elevation+1.e-30)'  #Speed
 
 if which_var == 4:  # Elevation
-    outname = name + '_elevation_node3'
+    outname = name + '_elevation'
     quantityname = 'elevation'  #Elevation
 
 print 'start sww2dem'
-#sys.stderr.write(sys.stdout.data)
 sww2dem(name, basename_out = outname,
             quantity = quantityname,
-            cellsize = 30,      # would prefer this at 25
+            cellsize = 25,      # would prefer this at 25
             # define region for viz purposes
             #easting_min = project.e_min_area,
@@ -52,8 +47,2 @@
             format = 'asc')
 
-#sys.stderr.write(sys.stdout.data)
-
-#Check
-
-#data = read_ermapper_grid(name)
-#print 'Values from %s are in [%f, %f]' %(name, min(data.flat), max(data.flat))
Index: /anuga_work/production/hobart_2006/project.py
===================================================================
--- /anuga_work/production/hobart_2006/project.py	(revision 3670)
+++ /anuga_work/production/hobart_2006/project.py	(revision 3671)
@@ -7,5 +7,19 @@
 import sys
 from time import localtime, strftime, gmtime
-                
+from anuga.utilities.polygon import read_polygon, plot_polygons
+from anuga.coordinate_transforms.redfearn import degminsec2decimal_degrees, convert_points_from_latlon_to_utm
+
+if sys.platform == 'win32':
+    home = getenv('INUNDATIONHOME')
+    user = getenv('USERPROFILE')
+
+else:    
+    home = getenv('INUNDATIONHOME', sep+'d'+sep+'xrd'+sep+'gem'+sep+'2'+sep+'ramp'+sep+'risk_assessment_methods_project'+sep+'inundation')     
+    user = getenv('LOGNAME')
+    print 'USER:', user
+
+# INUNDATIONHOME is the inundation directory, not the data directory.
+home += sep +'data'
+
 #Making assumptions about the location of scenario data
 state = 'tasmania'
@@ -50,20 +64,7 @@
 codename = 'project.py'
 
-if sys.platform == 'win32':
-    home = getenv('INUNDATIONHOME')
-    user = getenv('USERPROFILE')
-
-else:    
-    home = getenv('INUNDATIONHOME', sep+'d'+sep+'xrd'+sep+'gem'+sep+'2'+sep+'ramp'+sep+'risk_assessment_methods_project'+sep+'inundation')     
-    user = getenv('LOGNAME')
-    print 'USER:', user
-
-# INUNDATIONHOME is the inundation directory, not the data directory.
-home += sep +'data'
-
 #Derive subdirectories and filenames
 #time = strftime('%Y%m%d_%H%M%S',localtime()) #gets time for new dir
 local_time = strftime('%Y%m%d_%H%M%S',gmtime()) #gets time for new dir
-print 'home', home
 meshdir = home+sep+state+sep+scenario_dir_name+sep+'anuga'+sep+'meshes'+sep
 datadir = home+sep+state+sep+scenario_dir_name+sep+'anuga'+sep+'topographies'+sep
@@ -71,5 +72,4 @@
 polygondir = home+sep+state+sep+scenario_dir_name+sep+'anuga'+sep+'polygons'+sep
 boundarydir = home+sep+sep+state+sep+scenario_dir_name+sep+'anuga'+sep+'boundaries'+sep
-#output dir without time
 outputdir = home+sep+state+sep+scenario_dir_name+sep+'anuga'+sep+'outputs'+sep
 outputtimedir = outputdir + local_time + sep
@@ -79,5 +79,4 @@
 buildings_filename = gaugedir + 'hobart_res.csv'
 buildings_filename_damage_out = 'hobart_res_modified.csv'
-
 gaugetimeseries = gaugedir + 'hobart'
 
@@ -85,8 +84,6 @@
 #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'+sep
 
-codedir = getcwd()+sep
-                                
+codedir = getcwd()+sep                            
 codedirname = codedir + 'project.py'
-
 meshname = outputtimedir + 'mesh_' + basename
 
@@ -115,7 +112,47 @@
 offshore_dem_name_aho16 = datadir + offshore_name16
 coast_dem_name = datadir + coast_name
+# addition once total grid delivered
+#onshore_offshore_dem_name = datadir + '25m_se_tas' #25m grid
+onshore_offshore_dem_name = datadir + '50m_se_tas' #50m grid
+
 combined_dem_name = datadir + 'hobart_combined_elevation'
+combined_dem_name2 = datadir + 'hobart_combined_elevation2' # for alpha checking
 
 outputname = outputtimedir + basename  #Used by post processing
+
+###############################
+# Domain definitions
+###############################
+
+# bounding box
+south = degminsec2decimal_degrees(-43,45,0)
+north = degminsec2decimal_degrees(-42,30,0)
+west = degminsec2decimal_degrees(146,45,0)
+east = degminsec2decimal_degrees(148,0,0)#degminsec2decimal_degrees(148,15,0)
+
+#Main Domain of Hobart: first run JS 18/9/06
+d0 = [south, west]
+d1 = [south, east]
+d2 = [north, east]
+d3 = [north, west]
+polyAll, zone = convert_points_from_latlon_to_utm([d0, d1, d2, d3])
+refzone = zone
+
+# Second run - bottom bright, topr, top, left
+polyAll = [[520000, 5170000],[580000, 5170000],[580000, 5200000],[590000,5240000],[520000,5260000]]
+
+# Third run - afternoon Wed 27 Sep; surrounds -100m and 20mish elevation
+polyAll = read_polygon(polygondir+'new_extent.csv')
+plot_polygons([polyAll],'boundingpoly',verbose=False)
+
+###################################################################
+# Clipping regions for export to asc and regions for clipping data
+###################################################################
+
+# region to export for inundation map
+e_min_area = 521000#490000
+e_max_area = 522000#580000
+n_min_area = 5190000#5160000
+n_max_area = 5200000#5275000
 
 # clipping 12.5m onshore data set
@@ -130,38 +167,7 @@
 northingmax25 = 5260000
 
-# region to export for inundation map
-e_min_area = 521000#490000
-e_max_area = 522000#580000
-n_min_area = 5190000#5160000
-n_max_area = 5200000#5275000
-
-from anuga.coordinate_transforms.redfearn import degminsec2decimal_degrees, convert_points_from_latlon_to_utm
-                        
-# bounding box
-south = degminsec2decimal_degrees(-43,45,0)
-north = degminsec2decimal_degrees(-42,30,0)
-west = degminsec2decimal_degrees(146,45,0)
-east = degminsec2decimal_degrees(148,0,0)#degminsec2decimal_degrees(148,15,0)
-
-#Main Domain of Hobart: first run JS 18/9/06
-d0 = [south, west]
-d1 = [south, east]
-d2 = [north, east]
-d3 = [north, west]
-
-polyAll, zone = convert_points_from_latlon_to_utm([d0, d1, d2, d3])
-refzone = zone
-#polyAll = [[500000, 5200000],[580000, 5200000],[580000, 5270000],[500000,5270000]] #073937
-#polyAll = [[500000, 5160000],[580000, 5160000],[580000, 5270000],[500000,5270000]] # try this is the above works
-#polyAll = [[510000, 5170000],[570000, 5170000],[570000, 5250000],[510000,5250000]] # try this is the above works
-#polyAll = [[520000, 5160000],[580000, 5170000],[580000, 5200000],[600000,5260000],[510000,5280000]] # try this is the above works
-polyAll = [[520000, 5170000],[580000, 5170000],[580000, 5200000],[590000,5240000],[520000,5260000]] # try this is the above works
-# bottom bright, topr, top, left
-
-# region to export for Alex to make bathymetry map
-e_min_area = 560000#500000 #490000
-e_max_area = 570000#560000 #580000
-n_min_area = 5200000#5240000#5160000
-n_max_area = 5230000#5260000#5270000
+###############################
+# Interior region definitions
+###############################
 
 #Interior region - Hobart city area + Glenorchy, Kingston
@@ -188,4 +194,3 @@
 
 # Bruny Island
-from anuga.utilities.polygon import read_polygon
 poly_bruny = read_polygon(polygondir+'bruny.csv')
Index: /anuga_work/production/hobart_2006/run_hobart.py
===================================================================
--- /anuga_work/production/hobart_2006/run_hobart.py	(revision 3670)
+++ /anuga_work/production/hobart_2006/run_hobart.py	(revision 3671)
@@ -61,6 +61,5 @@
 
 # filenames
-onshore_dem_name = project.onshore_dem_name
-onshore_dem_name_25 = project.onshore_dem_name_25
+onshore_offshore_dem_name = project.onshore_offshore_dem_name
 meshname = project.meshname+'.msh'
 source_dir = project.boundarydir
@@ -68,54 +67,14 @@
 copied_files = False
 
-# creates DEM from asc data - 12.5m
-#convert_dem_from_ascii2netcdf(onshore_dem_name, use_cache=True, verbose=True)
-
-#creates pts file for onshore DEM - 12.5
-#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)
-
 # create DEM from asc data - 25m data
-convert_dem_from_ascii2netcdf(onshore_dem_name_25, use_cache=True, verbose=True)
-
-#creates pts file for onshore DEM - 25
-dem2pts(onshore_dem_name_25,
-        easting_min=project.eastingmin25,
-        easting_max=project.eastingmax25,
-        northing_min=project.northingmin25,
-        northing_max= project.northingmax25,
-        use_cache=True, verbose=True)
-
-#combine_rectangular_points_files(project.onshore_dem_name + '.pts',
-#                                 project.onshore_dem_name_25 + '.pts',
-#                                 project.all_onshore_dem_name + '.pts')
-
-print 'adding data sets'
-G = Geospatial_data(file_name = project.offshore_dem_name_local1 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_local2 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_local3 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_local4 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho1 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho2 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho3 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho4 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho5 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho6 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho7 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho8 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho9 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho10 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho11 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho12 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho13 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho14 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho15 + '.xya')+\
-    Geospatial_data(file_name = project.offshore_dem_name_aho16 + '.xya')+\
-    Geospatial_data(file_name = project.onshore_dem_name_25 + '.pts')
+convert_dem_from_ascii2netcdf(onshore_offshore_dem_name, use_cache=True, verbose=True)
+
+#creates pts file for combined DEM - 25
+dem2pts(onshore_offshore_dem_name, use_cache=True, verbose=True)
+
+# create geospatial data set and export 
+G = Geospatial_data(file_name = project.onshore_offshore_dem_name + '.pts')
 G.export_points_file(project.combined_dem_name + '.pts')
-#Geospatial_data(file_name = project.all_onshore_dem_name + '.pts')
+
 #----------------------------------------------------------------------------
 # Create the triangular mesh based on overall clipping polygon with a tagged
@@ -228,5 +187,5 @@
 # 7 min square wave starting at 1 min, 6m high
 Bw = Time_boundary(domain = domain,
-                   f=lambda t: [(60<t<480)*6, 0, 0])
+                   f=lambda t: [(60<t<480)*10, 0, 0])
 
 # for MOST BC
@@ -236,5 +195,5 @@
 # for testing
 domain.set_boundary( {'topr': Bd, 'left': Bd, 'top': Bd,
-                      'bottom': Bd, 'bright': Bw} )
+                      'bottom': Bw, 'bright': Bd} )
 
 #-------------------------------------------------------------------------------                                 
