Index: trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/anuga/boundaries/Patong_bounding_polygon/Patong_bounding_polygon.prj
===================================================================
--- trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/anuga/boundaries/Patong_bounding_polygon/Patong_bounding_polygon.prj	(revision 9332)
+++ trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/anuga/boundaries/Patong_bounding_polygon/Patong_bounding_polygon.prj	(revision 9332)
@@ -0,0 +1,1 @@
+PROJCS["WGS_1984_UTM_Zone_47N",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",99],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",1]]
Index: trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/anuga/boundaries/Patong_bounding_polygon/Patong_bounding_polygon.qpj
===================================================================
--- trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/anuga/boundaries/Patong_bounding_polygon/Patong_bounding_polygon.qpj	(revision 9332)
+++ trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/anuga/boundaries/Patong_bounding_polygon/Patong_bounding_polygon.qpj	(revision 9332)
@@ -0,0 +1,1 @@
+PROJCS["WGS 84 / UTM zone 47N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",99],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32647"]]
Index: trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/scripts/project.py
===================================================================
--- trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/scripts/project.py	(revision 9331)
+++ trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/scripts/project.py	(revision 9332)
@@ -14,5 +14,6 @@
 from anuga_parallel import myid, barrier, send, receive, numprocs
 # We need the boundary tags info to correctly cache the mesh
-from setup_boundary_conditions import boundary_info
+#from setup_boundary_conditions import boundary_info
+import user_functions
 
 
@@ -34,5 +35,5 @@
 ## Total model runtime in seconds (starts at time=0.)
 #
-finaltime=0.01*3600.0
+finaltime=0.05*3600.0
 
 ## Model time between output file timesteps in seconds
@@ -54,7 +55,18 @@
 ## Bounding polygon with default mesh resolution (as a max triangle area in m^2)
 #
-# Filename can be .shp or .csv
-#
-bounding_polygon_file='../anuga/polygons/bounding_polygon.csv'
+# Filename MUST be a .shp file, and contain a set of lines making up the boundary polygon
+# The lines each have an attribute giving the boundary tag, in the attribute table
+# column named (the value of boundary_tags_attribute_name)
+#
+# The starting/ending points on each line MUST EXACTLY AGREE so the whole
+# thing makes a polygon.
+#
+# Make it in GIS by first making the bounding polygon, then splitting it at
+# the appropriate locations, then adding the boundary tag attribute.
+#
+#bounding_polygon_file='../anuga/polygons/bounding_polygon.csv'
+bounding_polygon_and_tags_file='../anuga/boundaries/Patong_bounding_polygon/Patong_bounding_polygon.shp'
+boundary_tags_attribute_name='bndryTag'
+
 default_res=0.5*500**2  
 
@@ -188,5 +200,10 @@
 
 
-bounding_polygon=su.read_polygon(bounding_polygon_file)
+#bounding_polygon=su.read_polygon(bounding_polygon_file)
+bounding_polygon_and_tags =\
+    user_functions.read_boundary_tags_line_shapefile(bounding_polygon_and_tags_file, boundary_tags_attribute_name)
+bounding_polygon = bounding_polygon_and_tags[0]
+boundary_tags = bounding_polygon_and_tags[1]
+
 breaklines=su.readListOfBreakLines(breakline_files)
 riverwalls, riverwall_par=su.readListOfRiverWalls(riverwall_csv_files) #{} # Dictionary
@@ -218,5 +235,5 @@
 mesh_id_hash=hashlib.md5(json.dumps(\
             [bounding_polygon, interior_regions, riverwalls, \
-             breaklines, region_point_areas, default_res, boundary_info]) ).hexdigest()
+             breaklines, region_point_areas, default_res, boundary_tags]) ).hexdigest()
 
 # Fix the output tif bounding polygon
Index: trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/scripts/setup_boundary_conditions.py
===================================================================
--- trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/scripts/setup_boundary_conditions.py	(revision 9331)
+++ trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/scripts/setup_boundary_conditions.py	(revision 9332)
@@ -2,9 +2,4 @@
 
  THE USER MUST EDIT THIS FILE TO SET THE BOUNDARY CONDITIONS
-
- Both the boundary_info list and the function 'setup_boundary_conditions'
- require adjustment
-
- Specify boundary conditions
 
 """
@@ -13,34 +8,4 @@
 import anuga
 import anuga.utilities.spatialInputUtil as su
-
-#
-# For each boundary type, we have a list with a 
-#  1) Boundary tag
-#  3) Points -- for each point, the closest segment on the bounding_polygon
-#      will be assigned to the boundary type. The default boundary type has
-#     'Default' in place of points
-#
-
-#@@@@@@@@@@@@@@@@@#
-#START EDITING HERE
-#@@@@@@@@@@@@@@@@@#
-
-## OCEAN BOUNDARY
-# These points define the 'ocean_boundary': 1 point per edge on the bounding polygon
-#
-ocean_boundary_points_file='../anuga/boundaries/ocean_boundary_points/ocean_boundary_points.shp'
-ocean_boundary_points=su.read_points(ocean_boundary_points_file)
-
-
-
-## BOUNDARY TAGS
-# Polygon edges will be determined from the points
-# Any boundary polygon edge which is not the closest edge to a point will be
-# assigned to the 'Default' boundary
-#
-boundary_info=[
-              ['ocean_waterlevel_boundary', ocean_boundary_points],
-              ['vertices_without_tagpoints', 'Default']
-             ]
 
 
@@ -53,23 +18,27 @@
 
         2) Next make sure the 'boundary_tags_and_conditions' dictionary
-           associates the correct 'tags' (defined in boundary_info) to the
-           correct boundary condition type
+           associates the correct 'tags' (contained as attribute names in the
+           bounding_polygon_and_tags_file) with the correct boundary condition type
     """
+
+#@@@@@@@@@@@@@@@@@#
+#START EDITING HERE
+#@@@@@@@@@@@@@@@@@#
 
     # Example:
     # "Standard" boundary conditions
     Br = anuga.Reflective_boundary(domain)
-    Bs = anuga.Transmissive_stage_zero_momentum_boundary(domain) # Be careful with this one
+    #Bs = anuga.Transmissive_stage_zero_momentum_boundary(domain) # Be careful with this one
    
     # Example: 
     # Transmissive with stage being set by a function 'stage_function'
     def stage_function(t):
-        return min( max( t-120., 0.)/360., 1.)
+        return min( max( t-60., 0.)/60., 1.)
     B_t_stage_function = anuga.shallow_water.boundaries.Transmissive_n_momentum_zero_t_momentum_set_stage_boundary(domain, stage_function)
 
     ## Here we associate each tag in boundary_info with its boundary condition defined above
     boundary_tags_and_conditions={ 
-                                  'ocean_waterlevel_boundary': B_t_stage_function,
-                                  'vertices_without_tagpoints':Br 
+                                  'Ocean': B_t_stage_function,
+                                  'Br':Br 
                                   }
 
@@ -79,6 +48,6 @@
 
     ## Check that all tags in boundary_info have been set
-    for i in range(len(boundary_info)):
-        tag=boundary_info[i][0]
+    for i in range(len(project.boundary_tags)):
+        tag=project.boundary_tags.keys()[i]
         if not boundary_tags_and_conditions.has_key(tag):
             msg = 'Need to set boundary_tags_and_conditions for tag = ' + tag
Index: trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/scripts/setup_mesh.py
===================================================================
--- trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/scripts/setup_mesh.py	(revision 9331)
+++ trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/scripts/setup_mesh.py	(revision 9332)
@@ -39,45 +39,45 @@
 #####################################################################################
 
-def make_boundary_tags(bounding_polygon, setup_boundary_conditions):
-    """
-        Make the boundary tags 
-    """
-
-    # Close the bounding polygon
-    bounding_polygon_closed=bounding_polygon+[bounding_polygon[0]]
-
-    # Assign all tags to 'Default', but remove them from this list during
-    # assignment
-    default_tags=range(len(bounding_polygon))
-
-    # Dict for boundary tags
-    boundary_tags={}
-
-    # Information on the tags
-    boundary_info=setup_boundary_conditions.boundary_info
-
-    for i in range(len(boundary_info)):
-        bi=boundary_info[i]
-        # Shorthand for i'th boundary tag
-        tag_i=bi[0]
-        # Shorthand for i'th set of boundary points
-        PTS_i=bi[1]
-
-        # Add the boundary tag to the list
-        boundary_tags[tag_i]=[]
-
-        if(PTS_i=='Default' or i==(len(boundary_info)-1)):
-            # Assign default boundaries
-            assert (i==(len(boundary_info)-1) and PTS_i=='Default'),\
-                'Default boundary condition must be at the end of the boundary_info list'
-            boundary_tags[tag_i]=default_tags
-        else:
-            # Loop over all points in the boundary condition
-            for j in range(len(PTS_i)):
-                new_index=su.find_nearest_segment(PTS_i[j], bounding_polygon_closed)[1]
-                boundary_tags[tag_i]=boundary_tags[tag_i]+[new_index]     
-                default_tags.remove(new_index)
-
-    return boundary_tags    
+#def make_boundary_tags(bounding_polygon, setup_boundary_conditions):
+#    """
+#        Make the boundary tags 
+#    """
+#
+#    # Close the bounding polygon
+#    bounding_polygon_closed=bounding_polygon+[bounding_polygon[0]]
+#
+#    # Assign all tags to 'Default', but remove them from this list during
+#    # assignment
+#    default_tags=range(len(bounding_polygon))
+#
+#    # Dict for boundary tags
+#    boundary_tags={}
+#
+#    # Information on the tags
+#    boundary_info=setup_boundary_conditions.boundary_info
+#
+#    for i in range(len(boundary_info)):
+#        bi=boundary_info[i]
+#        # Shorthand for i'th boundary tag
+#        tag_i=bi[0]
+#        # Shorthand for i'th set of boundary points
+#        PTS_i=bi[1]
+#
+#        # Add the boundary tag to the list
+#        boundary_tags[tag_i]=[]
+#
+#        if(PTS_i=='Default' or i==(len(boundary_info)-1)):
+#            # Assign default boundaries
+#            assert (i==(len(boundary_info)-1) and PTS_i=='Default'),\
+#                'Default boundary condition must be at the end of the boundary_info list'
+#            boundary_tags[tag_i]=default_tags
+#        else:
+#            # Loop over all points in the boundary condition
+#            for j in range(len(PTS_i)):
+#                new_index=su.find_nearest_segment(PTS_i[j], bounding_polygon_closed)[1]
+#                boundary_tags[tag_i]=boundary_tags[tag_i]+[new_index]     
+#                default_tags.remove(new_index)
+#
+#    return boundary_tags    
 
 #################################################################################################
@@ -91,5 +91,5 @@
 
     # Create the boundary tags
-    boundary_tags=make_boundary_tags(project.bounding_polygon, setup_boundary_conditions)
+    #boundary_tags=make_boundary_tags(project.bounding_polygon, setup_boundary_conditions)
     
     # Ensure mesh_breaklines include riverwalls and breaklines
@@ -98,5 +98,5 @@
     ## Make the mesh
     anuga.create_mesh_from_regions(project.bounding_polygon,
-                             boundary_tags=boundary_tags,
+                             boundary_tags=project.boundary_tags,
                              maximum_triangle_area=project.default_res,
                              filename=project.meshname,
Index: trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/scripts/user_functions.py
===================================================================
--- trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/scripts/user_functions.py	(revision 9331)
+++ trunk/anuga_work/development/gareth/template_scenarios/tsunami/tsunami1/scripts/user_functions.py	(revision 9332)
@@ -39,2 +39,64 @@
     # Make a newline
     if(myid==0): print ''
+
+############################################################
+
+def read_boundary_tags_line_shapefile(shapefile_name, tag_attribute='bndryTag'):
+    """
+    Read in the boundary lines + tags from an appropriately structured line shapefile
+
+    Return the bounding polygon + dict with tags + indices
+    """
+
+    ## Step 1: Read the data
+    import ogr
+    
+    data_source=ogr.Open(shapefile_name)
+
+    assert data_source.GetLayerCount()==1, ' Bounding polygon + boundary tags file can only have 1 layer'
+
+    layer=data_source.GetLayer()
+
+    boundary_line_and_tags=[]
+
+    for feature in layer:
+        line = feature.GetGeometryRef().GetPoints()
+        line = [ list(l) for l in line]
+        tag = feature.GetField(tag_attribute)
+        boundary_line_and_tags.append([line, tag])
+
+
+    ## Step 2: Convert to bounding polygon + boundary tags in ANUGA format
+    boundary_segnum = len(boundary_line_and_tags)
+    # Initial values
+    bounding_polygon=boundary_line_and_tags[0][0]
+    boundary_tags={boundary_line_and_tags[0][1]: range(len(bounding_polygon)-1) }
+    # 
+    #import pdb
+    #pdb.set_trace()
+    for i in range(boundary_segnum-1):
+        # Find the line which starts at the end point of bounding_polygon
+        for j in range(boundary_segnum):
+            blj = boundary_line_and_tags[j][0]
+            if blj[0]==bounding_polygon[-1]:
+                # Append all but the first point of blj to bounding_polygon
+                old_lbp=len(bounding_polygon)
+                extra_polygon = [blj[k] for k in range(1, len(blj))]
+                bounding_polygon = bounding_polygon + extra_polygon
+                # If we are on the last segment, then drop the final point
+                # (for consistency with ANUGA format)
+                if i==(boundary_segnum-2):
+                    bounding_polygon = [bounding_polygon[k] for k in range(len(bounding_polygon)-1)]
+                    new_tags = range(old_lbp-1, len(bounding_polygon))
+                else:
+                    new_tags = range(old_lbp-1, len(bounding_polygon)-1)
+                
+                # Either the key already exists and we append edge indices, OR
+                # we add the key and the edge indices
+                key = boundary_line_and_tags[j][1]
+                if boundary_tags.has_key(key):
+                    boundary_tags[key] = boundary_tags[key]+ new_tags
+                else:
+                    boundary_tags[key] = new_tags
+
+    return bounding_polygon, boundary_tags 
