Index: /inundation/ga/storm_surge/pyvolution/data_manager.py
===================================================================
--- /inundation/ga/storm_surge/pyvolution/data_manager.py	(revision 1154)
+++ /inundation/ga/storm_surge/pyvolution/data_manager.py	(revision 1155)
@@ -1835,9 +1835,14 @@
 
 
-def sww2domain(filename,t=None,fail_if_NaN=True,NaN_filler=0,verbose = True):
+def sww2domain(filename,boundary=None,t=None,\
+               fail_if_NaN=True,NaN_filler=0\
+               ,verbose = True,very_verbose = False):
     """
     Usage: domain = sww2domain('file.sww',t=time (default = last time in file))
 
-    If the sww has stages, but not 
+    Boundary is not recommended if domian.smooth is not selected, as it       
+    uses unique coordinates, but not unique boundaries. This means that 
+    the boundary file will not be compatable with the coordiantes, and will
+    give a different final boundary, or crash.
     """
     NaN=9.969209968386869e+036
@@ -1845,6 +1850,6 @@
 
     from Scientific.IO.NetCDF import NetCDFFile
-    from domain import Domain
-    from Numeric import asarray, transpose
+    from shallow_water import Domain
+    from Numeric import asarray, transpose, resize
 
     if verbose: print 'Reading from ', filename
@@ -1891,10 +1896,27 @@
 
     if verbose: print '    building domain'
-    domain = Domain(coordinates, volumes,\
-                    conserved_quantities = conserved_quantities,\
-                    other_quantities = other_quantities,zone=zone,\
-                    xllcorner=xllcorner, yllcorner=yllcorner)
-    domain.starttime=starttime
-    domain.time=t
+#    From domain.Domain:
+#    domain = Domain(coordinates, volumes,\
+#                    conserved_quantities = conserved_quantities,\
+#                    other_quantities = other_quantities,zone=zone,\
+#                    xllcorner=xllcorner, yllcorner=yllcorner)
+
+#   From shallow_water.Domain:
+    coordinates=coordinates.tolist()
+    volumes=volumes.tolist()
+    #FIXME:should this be in mesh?(peter row)
+    unique,coordinates,volumes,boundary=weed(coordinates,volumes,boundary)
+
+    domain = Domain(coordinates, volumes, boundary)
+
+    if not boundary is None:
+        domain.boundary = boundary
+    domain.zone=zone
+    domain.xllcorner=float(xllcorner)
+    domain.yllcorner=float(yllcorner)
+
+    domain.starttime=float(starttime)+float(t)
+    domain.time=0.0
+
     for quantity in other_quantities:
         try:
@@ -1918,4 +1940,6 @@
                 data = (X<>NaN)
                 X = (X*data)+(data==0)*NaN_filler
+        if unique:
+            X = resize(X,(len(X)/3,3))
         domain.set_quantity(quantity,X)
 #
@@ -1941,4 +1965,6 @@
                 data = (X<>NaN)
                 X = (X*data)+(data==0)*NaN_filler
+        if unique:
+            X = resize(X,(X.shape[0]/3,3))
         domain.set_quantity(quantity,X)
     fid.close()
@@ -1985,4 +2011,61 @@
 
 
+def weed(coordinates,volumes,boundary = None):
+    if type(coordinates)=='array':
+        coordinates = coordinates.tolist()
+    if type(volumes)=='array':
+        volumes = volumes.tolist()
+
+    unique = False
+    point_dict = {}
+    same_point = {}
+    for i in range(len(coordinates)):
+        point = tuple(coordinates[i])
+        if point_dict.has_key(point):
+            unique = True
+            same_point[i]=point
+            #to change all point i references to point j
+        else: 
+            point_dict[point]=i
+            same_point[i]=point
+
+    coordinates = []
+    i = 0
+    for point in point_dict.keys():
+        point = tuple(point)
+        coordinates.append(list(point))
+        point_dict[point]=i
+        i+=1
+
+
+    for volume in volumes:
+        for i in range(len(volume)):
+            index = volume[i]
+            if index>-1:
+                volume[i]=point_dict[same_point[index]]
+
+    new_boundary = {}
+    if not boundary is None:
+        for segment in boundary.keys():
+            point0 = point_dict[same_point[segment[0]]]
+            point1 = point_dict[same_point[segment[1]]]
+            label = boundary[segment]
+            #FIXME should the bounday attributes be concaterated
+            #('exterior, pond') or replaced ('pond')(peter row)
+
+            if new_boundary.has_key((point0,point1)):
+                new_boundary[(point0,point1)]=new_boundary[(point0,point1)]#\
+                                              #+','+label
+
+            elif new_boundary.has_key((point1,point0)):
+                new_boundary[(point1,point0)]=new_boundary[(point1,point0)]#\
+                                              #+','+label
+            else: new_boundary[(point0,point1)]=label
+
+        boundary = new_boundary               
+
+    return unique,coordinates,volumes,boundary
+        
+        
 
 #OBSOLETE STUFF
Index: /inundation/ga/storm_surge/pyvolution/test_data_manager.py
===================================================================
--- /inundation/ga/storm_surge/pyvolution/test_data_manager.py	(revision 1154)
+++ /inundation/ga/storm_surge/pyvolution/test_data_manager.py	(revision 1155)
@@ -1341,6 +1341,7 @@
 
 	#Create basic mesh
-	points, vertices, boundary = rectangular(2,2)
-
+
+        yiel=0.01
+	points, vertices, boundary = rectangular(10,10)
 	#Create shallow water domain
 	domain = Domain(points, vertices, boundary)
@@ -1358,8 +1359,9 @@
 	Bt = Transmissive_boundary(domain)
 	Bd = Dirichlet_boundary([0.2,0.,0.])
-	Bw = Time_boundary(domain=domain,
+ 	Bw = Time_boundary(domain=domain,
 			   f=lambda t: [(0.1*sin(t*2*pi)), 0.0, 0.0])
 
 	domain.set_boundary({'left': Bd, 'right': Br, 'top': Br, 'bottom': Br})
+
 	domain.quantities_to_be_stored.extend(['xmomentum','ymomentum'])
 	#Initial condition
@@ -1367,11 +1369,10 @@
 	elevation = domain.quantities['elevation'].vertex_values
 	domain.set_quantity('stage', elevation + h)
-	#elevation = domain.get_quantity('elevation',location='unique vertices')
-	#domain.set_quantity('stage', elevation + h,location='unique vertices')
+	#elevation = domain.get_quantity('elevation')
+	#domain.set_quantity('stage', elevation + h)
 
 	domain.check_integrity()
-	dir(domain)
 	#Evolution
-	for t in domain.evolve(yieldstep = 1, finaltime = 2.0):
+	for t in domain.evolve(yieldstep = yiel, finaltime = 0.05):
 	#    domain.write_time()
 	    pass
@@ -1383,9 +1384,65 @@
 	from data_manager import sww2domain
 	from Numeric import allclose
-        import os
-
-	filename = domain.datadir+os.sep+domain.filename+'.sww'
-
-	domain2 = sww2domain(filename,fail_if_NaN=False,verbose = False)
+
+	filename = domain.datadir+'\\'+domain.filename+'.sww'
+	domain2 = sww2domain(filename,None,fail_if_NaN=False,verbose = False)
+	#points, vertices, boundary = rectangular(15,15)
+        #domain2.boundary = boundary
+        ###################
+	##NOW TEST IT!!!
+        ###################
+
+	bits = ['xllcorner','yllcorner','vertex_coordinates']
+
+	for quantity in ['elevation']+domain.quantities_to_be_stored:
+	    bits.append('quantities["%s"].get_integral()'%quantity)
+	    bits.append('get_quantity("%s")'%quantity)
+
+	for bit in bits:
+	#    print 'testing that domain.'+bit+' has been restored'
+        #    print bit
+        #print 'done'
+	    assert allclose(eval('domain.'+bit),eval('domain2.'+bit))
+
+        ######################################
+        #Now evolve them both, just to be sure
+        visualise = False
+        #visualise = True
+	domain.visualise = visualise
+        domain.time = 0.
+        from time import sleep
+
+        final = .1
+        domain.set_quantity('friction', 0.1)
+	domain.store = False
+	domain.set_boundary({'left': Bd, 'right': Bd, 'top': Bd, 'bottom': Bd})
+
+
+	for t in domain.evolve(yieldstep = yiel, finaltime = final):
+            if visualise: sleep(1.)
+	    #domain.write_time()
+	    pass
+
+	domain2.smooth = False
+	domain2.visualise = visualise
+	domain2.store = False
+	domain2.default_order=2
+	domain2.set_quantity('friction', 0.1)
+	#Bed-slope and friction
+	# Boundary conditions
+        Bd2=Dirichlet_boundary([0.2,0.,0.])
+        domain2.boundary = domain.boundary
+        #print 'domain2.boundary'
+        #print domain2.boundary
+	domain2.set_boundary({'left': Bd, 'right': Bd, 'top': Bd, 'bottom': Bd})
+	#domain2.set_boundary({'exterior': Bd})
+
+	domain2.check_integrity()
+
+	for t in domain2.evolve(yieldstep = yiel, finaltime = final):
+            if visualise: sleep(1.)
+	    #domain2.write_time()
+	    pass
+
 
         ###################
@@ -1393,15 +1450,15 @@
         ##################
 
-	bits = ['xllcorner','yllcorner','vertex_coordinates','time','starttime']
-
-	for quantity in ['elevation']+domain.quantities_to_be_stored:
+	bits = ['xllcorner','yllcorner','vertex_coordinates']
+
+	for quantity in ['elevation','xmomentum','ymomentum']:#+domain.quantities_to_be_stored:
+	    bits.append('quantities["%s"].get_integral()'%quantity)
 	    bits.append('get_quantity("%s")'%quantity)
 
 	for bit in bits:
-	#    print 'testing that domain.'+bit+' has been restored'
+            #print bit
 	    assert allclose(eval('domain.'+bit),eval('domain2.'+bit))
 
-	#print 'passed'
-
+        
 
     def test_sww2domain2(self):
@@ -1443,6 +1500,4 @@
 	elevation = domain.quantities['elevation'].vertex_values
 	domain.set_quantity('stage', elevation + h)
-	#elevation = domain.get_quantity('elevation',location='unique vertices')
-	#domain.set_quantity('stage', elevation + h,location='unique vertices')
 
 	domain.check_integrity()
@@ -1451,4 +1506,6 @@
 	    pass
 	    #domain.write_time()
+
+	
 
         ##################################
@@ -1459,19 +1516,19 @@
         import os
 
-
 	filename = domain.datadir+os.sep+domain.filename+'.sww'
 
 	#Fail because NaNs are present
 	try:
-	    domain2 = sww2domain(filename,fail_if_NaN=True,verbose=False)
+	    domain2 = sww2domain(filename,boundary,fail_if_NaN=True,verbose=False)
 	    assert True == False
 	except:
 	    #Now import it, filling NaNs to be 0
 	    filler = 0
-	    domain2 = sww2domain(filename,fail_if_NaN=False,NaN_filler = filler,verbose=False)
-
-	bits = ['xllcorner','yllcorner','vertex_coordinates','time','starttime']
+	    domain2 = sww2domain(filename,None,fail_if_NaN=False,NaN_filler = filler,verbose=False)
+
+	bits = ['xllcorner','yllcorner','vertex_coordinates']
 
 	for quantity in ['elevation']+domain.quantities_to_be_stored:
+	    bits.append('quantities["%s"].get_integral()'%quantity)
 	    bits.append('get_quantity("%s")'%quantity)
 
@@ -1479,9 +1536,4 @@
 	#    print 'testing that domain.'+bit+' has been restored'
 	    assert allclose(eval('domain.'+bit),eval('domain2.'+bit))
-
-	#print max(max(domain2.get_quantity('xmomentum')))
-	#print min(min(domain2.get_quantity('xmomentum')))
-	#print max(max(domain2.get_quantity('ymomentum')))
-	#print min(min(domain2.get_quantity('ymomentum')))
 
 	assert max(max(domain2.get_quantity('xmomentum')))==filler
@@ -1496,8 +1548,24 @@
 	#os.remove(domain.datadir+'/'+domain.filename+'.sww')
 
+
+    #def test_weed(self):
+	#from data_manager import weed
+
+        #coordinates1 = [[0.,0.],[1.,0.],[1.,1.],[1.,0.],[2.,0.],[1.,1.]]
+        #volumes1 = [[0,1,2],[3,4,5]]
+        #boundary1= {(0,1): 'external',(1,2): 'not external',(2,0): 'external',(3,4): 'external',(4,5): 'external',(5,3): 'not external'}
+        #unique2,coordinates2,volumes2,boundary2=weed(coordinates1,volumes1,boundary1)
+        #print 'coordinates2'
+        #print coordinates2
+        #print 'volumes2'
+        #print volumes2
+        #print 'boundary2'
+        #print boundary2
+
+
 #-------------------------------------------------------------
 if __name__ == "__main__":
     suite = unittest.makeSuite(Test_Data_Manager,'test')    
-    #suite = unittest.makeSuite(Test_Data_Manager,'test_sww2asc_mis')
+    #suite = unittest.makeSuite(Test_Data_Manager,'test_sww2asc_stage')
     runner = unittest.TextTestRunner()
     runner.run(suite)
