Changeset 5517

Timestamp:

Jul 17, 2008, 7:48:55 PM (16 years ago)

Author:

ole

Message:

Work on test_quantity towards general ability to use a polygon.

File:

• anuga_core/source/anuga/abstract_2d_finite_volumes/test_quantity.py (modified) (14 diffs)

anuga_core/source/anuga/abstract_2d_finite_volumes/test_quantity.py

@@ -14 +14 @@
 from anuga.geospatial_data.geospatial_data import Geospatial_data
 from anuga.coordinate_transforms.geo_reference import Geo_reference
+from anuga.utilities.polygon import *
 
 #Aux for fit_interpolate.fit example
@@ -417 +418 @@
         quantity.set_values(3.14, indices=[0,2], location='vertices')
 
-        # Indices refer to triangle numbers here - not vertices (why?)
+        # Indices refer to triangle numbers
         assert allclose(quantity.vertex_values,
                         [[3.14,3.14,3.14], [0,0,0],
@@ -434 +435 @@
 
 
-        # Another polygon (pick triangle 1 and 2 (rightmost triangles)
+        # Another polygon (pick triangle 1 and 2 (rightmost triangles) using centroids
         polygon = [[2.1, 0.0], [3.5,0.1], [2,2.2], [0.2,2]]
         quantity.set_values(0.0)
-        quantity.set_values(3.14, polygon=polygon)
-
+        quantity.set_values(3.14, location='centroids', polygon=polygon)
         assert allclose(quantity.vertex_values,
                         [[0,0,0],
@@ -446 +446 @@
 
 
+        # Same polygon now use vertices (default)
+        polygon = [[2.1, 0.0], [3.5,0.1], [2,2.2], [0.2,2]]
+        quantity.set_values(0.0)
+        #print 'Here 2'
+        quantity.set_values(3.14, polygon=polygon)
+        assert allclose(quantity.vertex_values,
+                        [[0,0,0],
+                         [3.14,3.14,3.14],
+                         [3.14,3.14,3.14],                         
+                         [0,0,0]])                
+        
+        #print 'done'
 
         # Test input checking
@@ -488 +500 @@
         
         quantity.set_values(0.0)
-        quantity.set_values(3.14, polygon=polygon)
+        quantity.set_values(3.14, polygon=polygon, location='centroids')
         
         assert allclose(quantity.vertex_values,
@@ -680 +692 @@
         os.remove(ptsfile)
 
-    def Cache_cache_test_set_values_from_file(self):
+
+
+    def XXXXtest_set_values_from_file_using_polygon(self):
         quantity = Quantity(self.mesh4)
 
@@ -715 +729 @@
         file.close()
 
+        # Create restricting polygon (containing vertex #4 (2,2) or centroid of triangle #1 (bce)
+        polygon = [[1.0, 1.0], [4.0, 1.0],
+                   [4.0, 4.0], [1.0, 4.0]]
+
+        #print self.mesh4.nodes
+        #print inside_polygon(self.mesh4.nodes, polygon)
+        assert allclose(inside_polygon(self.mesh4.nodes, polygon), 4)
+
+        #print quantity.domain.get_vertex_coordinates()
+        #print quantity.domain.get_nodes()        
+        
+        # Check that values can be set from file
+        quantity.set_values(filename=ptsfile,
+                            polygon=polygon,
+                            location='unique vertices',
+                            alpha=0)
+
+        # Get indices for vertex coordinates in polygon
+        indices = inside_polygon(quantity.domain.get_vertex_coordinates(), polygon)
+        points = take(quantity.domain.get_vertex_coordinates(), indices)
+        
+        answer = linear_function(points)
+
+        #print quantity.vertex_values.flat
+        #print answer
+
+        # Check vertices in polygon have been set
+        assert allclose(take(quantity.vertex_values.flat, indices),
+                             answer)
+
+        # Check vertices outside polygon are zero
+        indices = outside_polygon(quantity.domain.get_vertex_coordinates(), polygon)        
+        assert allclose(take(quantity.vertex_values.flat, indices),
+                             0.0)        
+
+        #Cleanup
+        import os
+        os.remove(ptsfile)
+
+
+        
+
+    def Cache_cache_test_set_values_from_file(self):
+        quantity = Quantity(self.mesh4)
+
+        #Get (enough) datapoints
+        data_points = [[ 0.66666667, 0.66666667],
+                       [ 1.33333333, 1.33333333],
+                       [ 2.66666667, 0.66666667],
+                       [ 0.66666667, 2.66666667],
+                       [ 0.0, 1.0],
+                       [ 0.0, 3.0],
+                       [ 1.0, 0.0],
+                       [ 1.0, 1.0],
+                       [ 1.0, 2.0],
+                       [ 1.0, 3.0],
+                       [ 2.0, 1.0],
+                       [ 3.0, 0.0],
+                       [ 3.0, 1.0]]
+
+        data_geo_spatial = Geospatial_data(data_points,
+                         geo_reference = Geo_reference(56, 0, 0))
+        data_points_absolute = data_geo_spatial.get_data_points(absolute=True)
+        attributes = linear_function(data_points_absolute)
+        att = 'spam_and_eggs'
+        
+        #Create .txt file
+        ptsfile = tempfile.mktemp(".txt")
+        file = open(ptsfile,"w")
+        file.write(" x,y," + att + " \n")
+        for data_point, attribute in map(None, data_points_absolute
+                                         ,attributes):
+            row = str(data_point[0]) + ',' + str(data_point[1]) \
+                  + ',' + str(attribute)
+            file.write(row + "\n")
+        file.close()
+
 
         #Check that values can be set from file
@@ -884 +975 @@
         os.remove(pts_file)
         
-    def test_set_values_from_UTM_pts(self):
+    def verbose_test_set_values_from_UTM_pts(self):
         quantity = Quantity(self.mesh_onslow)
 
-        VERBOSE = False # Change this to True to see output
         #Get (enough) datapoints
         data_points = [[-21.5, 114.5],[-21.4, 114.6],[-21.45,114.65],
@@ -942 +1032 @@
         #Check that values can be set from file
         quantity.set_values_from_file(pts_file, att, 0,
-                                      'vertices', None, verbose = VERBOSE,
+                                      'vertices', None, verbose = True,
                                       max_read_lines=2)
         answer = linear_function(quantity.domain.get_vertex_coordinates())
@@ -950 +1040 @@
 
         #Check that values can be set from file
-        quantity.set_values(filename = pts_file, verbose = VERBOSE,
+        quantity.set_values(filename = pts_file,
                             attribute_name = att, alpha = 0)
         answer = linear_function(quantity.domain.get_vertex_coordinates())
@@ -962 +1052 @@
         assert allclose(quantity.vertex_values.flat, answer)
 
-        #Check that values can be block read from a text file
-        quantity.set_values_from_file(txt_file, att, 0,
-                                      'vertices', None, verbose = VERBOSE,
-                                      max_read_lines=2)
-        
         #Cleanup
         import os
@@ -1184 +1269 @@
 
 
-        #Negation
+        # Negation
         Q = -quantity1
         assert allclose(Q.vertex_values, -quantity1.vertex_values)
@@ -1190 +1275 @@
         assert allclose(Q.edge_values, -quantity1.edge_values)
 
-        #Addition
+        # Addition
         Q = quantity1 + 7
         assert allclose(Q.vertex_values, quantity1.vertex_values + 7)
@@ -2349 +2434 @@
 #-------------------------------------------------------------
 if __name__ == "__main__":
-    suite = unittest.makeSuite(Test_Quantity, 'test')
-
-    #suite = unittest.makeSuite(Test_Quantity, 'test_set_values_from_UTM')
+    suite = unittest.makeSuite(Test_Quantity, 'test')    
+    #suite = unittest.makeSuite(Test_Quantity, 'test_set_values_from_file_using_polygon')
+
+    #suite = unittest.makeSuite(Test_Quantity, 'test_set_vertex_values_using_general_interface_with_subset')
     #print "restricted test"
     #suite = unittest.makeSuite(Test_Quantity,'verbose_test_set_values_from_UTM_pts')
-    runner = unittest.TextTestRunner() #verbosity=2)
+    runner = unittest.TextTestRunner()
     runner.run(suite)