Index: anuga_work/development/Hinwood_2008/create_mesh.py
===================================================================
--- anuga_work/development/Hinwood_2008/create_mesh.py	(revision 5350)
+++ anuga_work/development/Hinwood_2008/create_mesh.py	(revision 5370)
@@ -1,3 +1,3 @@
-"""Create mesh for University of Queensland dam break flume
+"""Create mesh for Jon Hinwoods wave flume
 """
 
@@ -6,12 +6,13 @@
 from anuga.coordinate_transforms.geo_reference import Geo_reference
 
-xslope = 5.0 # Distance between the boundary ande the start of the slope
+xslope = 5.0 # Distance between the boundary and the start of the slope
 
 
-def generate(mesh_filename, maximum_triangle_area=0.01):
+def generate(mesh_filename, slope, maximum_triangle_area=0.01):
     """
-    Generate mesh for University of Aberbeen dam break flume.
-    The gate position is the distance from the right wall of the wave tank
-    to the gate. 
+    Generate mesh for Monash University wave generation flume.
+    
+    Origin of x coordinate is the toe of the beach, x measured
+    positive shorewards.
     
     """
@@ -19,8 +20,10 @@
     global xslope
     
-    xright  = 19.0
+    xleft = slope['xleft'][0]
+    xtoe = slope['xtoe'][0]
+    xbeach = slope['xbeach'][0]
+    xright = slope['xright'][0]
     ybottom = 0
     ytop    = 1.00
-    xleft = 0.0
     ###xslope = slope
 
@@ -32,8 +35,12 @@
 
 
-    # slope seperation (middle)
-    point_slope_top = [xslope, ytop]
-    point_slope_bottom = [xslope, ybottom]    
+    # 1st slope seperation (middle)
+    point_toe_top = [xtoe, ytop]
+    point_toe_bottom = [xtoe, ybottom]    
 
+    # 2nd slope seperation (middle)
+    point_beach_top = [xbeach, ytop]
+    point_beach_bottom = [xbeach, ybottom]
+    
     m = Mesh()
 
@@ -41,8 +48,10 @@
     points = [point_sw,   #se
               point_nw,
-              point_slope_top,
+              point_toe_top, # 2
+              point_beach_top, # 3
               point_ne, 
               point_se,
-              point_slope_bottom
+              point_beach_bottom, # 6
+              point_toe_bottom # 7
               ]
     
@@ -53,20 +62,23 @@
         [3,4],
         [4,5],
-        [5,0],  #The outer border
-        [2,5]       # slope separator
+        [5,6],
+        [6,7],
+        [7,0],  #The outer border
+        [2,7],       # toe separator
+        [3,6]       # _beach separator
         ]    
     
-    segment_tags = {'wall':[1,2,3,4,5],
+    segment_tags = {'wall':[1,2,3,4,5,6,7],
                     'wave':[0]} # '':[6]
         
     m.add_points_and_segments(points, segments, segment_tags)
     
-    dam = m.add_region(xslope - 0.0000001,(ytop - ybottom)/2)
-    # this is the location of the reservoir region.
-    dam.setTag("flat")
+#     dam = m.add_region(xslope - 0.0000001,(ytop - ybottom)/2)
+#     # this is the location of the reservoir region.
+#     dam.setTag("flat")
     
-    slope = m.add_region(xslope + 0.0000001,(ytop - ybottom)/2)
-    # this is the location of the slope region.
-    slope.setTag("slope")
+#     slope = m.add_region(xslope + 0.0000001,(ytop - ybottom)/2)
+#     # this is the location of the slope region.
+#     slope.setTag("slope")
     
     m.generate_mesh(maximum_triangle_area=maximum_triangle_area)
Index: anuga_work/development/Hinwood_2008/prepare_time_boundary.py
===================================================================
--- anuga_work/development/Hinwood_2008/prepare_time_boundary.py	(revision 5350)
+++ anuga_work/development/Hinwood_2008/prepare_time_boundary.py	(revision 5370)
@@ -25,6 +25,5 @@
 def prepare_time_boundary(filename):
     """Convert benchmark 2 time series to NetCDF tms file.
-    This is a 'throw-away' code taylor made for files like
-    'Benchmark_2_input.txt' from the LWRU2 benchmark
+   
     """
 
@@ -35,7 +34,4 @@
     # Read the ascii (.csv) version of this file
     fid = open(filename[:-4] + '.csv')
-
-    # Skip first line
-    line = fid.readline()
 
     # Read remaining lines
@@ -62,6 +58,6 @@
 
 
+    
     # Create tms NetCDF file
-
     fid = NetCDFFile(filename, 'w')
     fid.institution = 'Geoscience Australia'
@@ -87,7 +83,9 @@
     """
     
-    Convert the rawish velocity and depth values, which have values at
+    Convert the raw velocity and depth values, which have values at
     different times to a csv file, with values at the same time, with
     SI units.
+
+    Set the depth values to be at the same times as the velocity values.
 
     The format for the velocity file is;
@@ -98,5 +96,6 @@
     [time, sec], [mm above SWL for sensor A], many other sensors...
     """
-
+    missing = 1e+20
+    
     # Read velocity file
     vfid = open(velocity_file)
@@ -107,11 +106,13 @@
     n_velocity = len(lines)
     vtimes = zeros(n_velocity, Float)  #Time
-    velocities = zeros(n_velocity, Float)  #
+    x_velocities = zeros(n_velocity, Float)  #
+    y_velocities = zeros(n_velocity, Float)  #
 
     for i, line in enumerate(lines):
-        fields = line.split(',')
+        fields = line.split() #(',')
 
         vtimes[i] = float(fields[0])
-        velocities[i] = float(fields[1])
+        x_velocities[i] = float(fields[1])
+        y_velocities[i] = float(fields[2])
 
     # Read the depth file
@@ -120,36 +121,91 @@
     dfid.close()
 
-
+    
     n_depth = len(lines)
+    n_sensors = len(lines[0].split())
     dtimes = zeros(n_depth, Float)  #Time
     depths = zeros(n_depth, Float)  #
+    sensors = zeros((n_depth,n_sensors), Float)
     
     for i, line in enumerate(lines):
-        fields = line.split(',')
-
-        dtimes[i] = float(fields[0])
-        depths[i] = float(fields[1])
-
-    depths_at_vtimes = interp(dtimes, depths, vtimes, missing=1e+20)
+        fields = line.split() #(',')
+        fields = [float(j) for j in fields]
+        dtimes[i] = fields[0]
+        depths[i] = fields[1]
+        sensors[i] = fields
+    #print "dtimes", dtimes
+    #print "depths", depths
+    #print "vtimes", vtimes
+    depths_at_vtimes = interp( depths, dtimes,
+                               vtimes, missing=missing)
     depths_at_vtimes = ensure_numeric(depths_at_vtimes)
+
+    #print "len(dtimes)", len(vtimes)
+    #print "len(depths_at_vtimes)", len(depths_at_vtimes)
+#     for i in range(len(depths_at_vtimes)):
+#         print "i", i
+#         print "vtimes[i]", vtimes[i]
+#         print "depths_at_vtimes[i]", depths_at_vtimes[i]
+#     print "depths_at_vtimes",  depths_at_vtimes
     depths_at_vtimes = depths_at_vtimes/1000.00 # convert from mm to m
-    velocities = ensure_numeric(velocities)
-    velocities = velocities * -1.0 # Swap axis around
+    missing=missing/1000.00 # Do to missing what is done to depths_at_vtimes
+    x_velocities = ensure_numeric(x_velocities)
+    x_velocities = x_velocities * -1.0 # Swap axis around
+    y_velocities = ensure_numeric(y_velocities)
+    y_velocities = y_velocities * -1.0 # Swap axis around
     
     fid = open(out_file,'w')
 
     assert len(depths_at_vtimes) == len(vtimes)
-    
+    start_time = None
+    #start_time = 0.0
     #for vtime, depth_at_vtime, velocity in map(vtimes, depths_at_vtimes,
     #                                           velocities):
     for i in xrange(len(vtimes)):
-        fid.write(str(vtimes[i]) + ',' + str(depths_at_vtimes[i]) \
-                  + ',' + str(velocities[i])+'\n')
-
-    fid.close()
+        if not depths_at_vtimes[i] == missing:
+            # Make the times start at zero.
+            if start_time is None:
+                start_time = vtimes[i]            
+            fid.write(str(vtimes[i]-start_time) \
+                      + ',' + str(depths_at_vtimes[i]) \
+                      + ',' + str(x_velocities[i]) \
+                      + ',' + str(y_velocities[i])+'\n')
+
+    fid.close()
+    
+    # Since there is a new time reference save the depth info using this
+    # new reference.
+    fid = open(depth_file[:-4] + '_new_time'+depth_file[-4:],'w')
+    sensors[:,0] -= start_time
+    #print "depth_file", depth_file
+    #print "start_time", start_time
+    for i in xrange(len(dtimes)):
+        fid.write(str(sensors[i,0]))
+        for j in xrange(1,len(sensors[0])):
+            fid.write(' ' + str(sensors[i,j]))
+        fid.write('\n')
+    fid.close()
+
     
     
 #-------------------------------------------------------------------
 if __name__ == "__main__":
-    combine_velocity_depth('T2R7velfilt.csv','T2R7pressfilt.csv', 'cyeah.csv')
-
+    from os import getenv
+    from os.path import join
+    home = getenv('INUNDATIONHOME')
+    Hinwood_dir = join(home,'data','flumes','Hinwood_2008')
+    raw_data_dir = join(Hinwood_dir, 'raw_data')
+
+    # Test 1 Run 5
+    combine_velocity_depth(join(raw_data_dir,'T1R5velfilt.txt'),
+                           join(raw_data_dir,'T1R5pressfilt.txt'),
+                           join(Hinwood_dir, 'T1R5_boundary.csv'))    
+    # Create the tsm file
+    prepare_time_boundary(join(Hinwood_dir, 'T1R5_boundary.tsm'))
+    
+    # Test 2 Run 7
+    combine_velocity_depth(join(raw_data_dir,'T2R7velfilt.txt'),
+                           join(raw_data_dir,'T2R7pressfilt.txt'),
+                           join(Hinwood_dir, 'T2R7_boundary.csv'))
+    # Create the tsm file
+    prepare_time_boundary(join(Hinwood_dir, 'T2R7_boundary.tsm'))
Index: anuga_work/development/Hinwood_2008/project.py
===================================================================
--- anuga_work/development/Hinwood_2008/project.py	(revision 5350)
+++ anuga_work/development/Hinwood_2008/project.py	(revision 5370)
@@ -38,5 +38,5 @@
         self.meshdir = scenariodir+sep+'meshes'+sep
         self.scenariodir = scenariodir+sep
-        self.boundary_file = self.scenariodir + 'boundary.tsm'
+        self.boundarydir = self.scenariodir
 
         # creates copy of output dir structure, if it doesn't exist
Index: anuga_work/development/Hinwood_2008/run_dam.py
===================================================================
--- anuga_work/development/Hinwood_2008/run_dam.py	(revision 5350)
+++ anuga_work/development/Hinwood_2008/run_dam.py	(revision 5370)
@@ -22,6 +22,6 @@
 from shutil import copy
 from os import path, sep
-from os.path import dirname  #, basename
-
+from os.path import dirname, join  #, basename
+from Numeric import zeros, size, Float
 
 # Related major packages
@@ -40,22 +40,26 @@
 import create_mesh
 from prepare_time_boundary import prepare_time_boundary
-
-
-
-def elevation_function(x,y):
-    from Numeric import zeros, size, Float
-
-    xslope = create_mesh.xslope  #4 ## Bit of a magic Number
-    print "xslope",xslope
-    z = zeros(size(x), Float)
-    for i in range(len(x)):
-        if x[i] < xslope:
-            z[i] = 0.0  #WARNING: the code in prepare_time_boundary
-                        # that calc's momentum assumes this is 0.0
-        else:
-            z[i] = (x[i]-xslope)*(1./16.)
-    return z 
-
-def main(friction=0.01, outputdir_name=None, is_trial_run=False):
+from interp import interp
+
+
+class Elevation_function:
+    def __init__(self, slope):
+        self.xslope_position = [slope['xleft'][0],slope['xtoe'][0],
+                  slope['xbeach'][0],slope['xright'][0]]
+        self.yslope_height = [slope['xleft'][1],slope['xtoe'][1],
+                  slope['xbeach'][1],slope['xright'][1]]
+        
+    def __call__(self, x,y):
+        
+        z = interp(self.yslope_height, self.xslope_position, x)
+        return z 
+
+def main(boundary_file,
+         boundary_location,
+         slope,
+         boundary_path=None,
+         friction=0.01,
+         outputdir_name=None,
+         is_trial_run=False):
 
     
@@ -84,11 +88,14 @@
           pro_instance.outputdir + 'create_mesh.py')
 
-    
+    # Boundary file manipulation
+    if boundary_path is None:
+        boundary_path = pro_instance.boundarydir
+    boundary_file_path = join(boundary_path, boundary_file)
     # Convert the boundary file, .csv to .tsm
     try:
-        temp = open(pro_instance.boundary_file)
+        temp = open(boundary_file_path)
         temp.close()
     except IOError:
-        prepare_time_boundary(pro_instance.boundary_file)
+        prepare_time_boundary(boundary_file_path)
     
     mesh_filename = pro_instance.meshdir + basename + '.msh'
@@ -110,5 +117,5 @@
     # this creates the mesh
     #gate_position = 12.0
-    create_mesh.generate(mesh_filename,
+    create_mesh.generate(mesh_filename, slope,
                          maximum_triangle_area=maximum_triangle_area)
 
@@ -138,5 +145,6 @@
 
     domain.set_quantity('stage', 0.4)
-    domain.set_quantity('friction', friction) 
+    domain.set_quantity('friction', friction)
+    elevation_function = Elevation_function(slope)
     domain.set_quantity('elevation', elevation_function)
 
@@ -148,5 +156,5 @@
     #Bts = Transmissive_Momentum_Set_Stage_boundary(domain, function)
     try:
-        function = file_function(pro_instance.boundary_file, domain,
+        function = file_function(boundary_file_path, domain,
                                  verbose=True)
     except IOError:
@@ -156,7 +164,8 @@
         
     Br = Reflective_boundary(domain)
-    #Bd = Dirichlet_boundary([0.3,0,0]) 
+    Bd = Dirichlet_boundary([0.3,0,0]) 
     Bts = Time_boundary(domain, function)
     domain.set_boundary( {'wall': Br, 'wave': Bts} )
+    domain.set_boundary( {'wall': Br, 'wave': Bd} )
 
     #-------------------------------------------------------------------------
@@ -165,4 +174,6 @@
     t0 = time.time()
 
+    # It seems that ANUGA can't handle a starttime that is >0.
+    domain.starttime = 1.0
     for t in domain.evolve(yieldstep, finaltime):    
         domain.write_time()
@@ -196,4 +207,10 @@
 #-------------------------------------------------------------
 if __name__ == "__main__":
-    main( is_trial_run = True,
+    #slopes = [[-4.5,0.0],[0.0,0.0],[1.285,0.090],[16.1,.960]]
+    slope = {'xleft':[-4.5,0.0],
+              'xtoe':[0.0,0.0],
+              'xbeach':[1.285,0.090],
+              'xright':[16.1,.960]}
+    main( 'T1R5_boundary.tsm', 8, slope,
+          is_trial_run = True,
          outputdir_name='Hinwood_low_stage_low_velocity_draft')
