Changeset 3969

Timestamp:

Nov 13, 2006, 1:53:33 PM (18 years ago)

Author:

nick

Message:

added split and get sample modules to geospatial_data
also removed xxx_add_data_points this was not being used and is old

Location:

anuga_core/source/anuga/geospatial_data

Files:

• geospatial_data.py (modified) (7 diffs)
• test_geospatial_data.py (modified) (1 diff)

anuga_core/source/anuga/geospatial_data/geospatial_data.py

@@ -7 +7 @@
 from types import DictType
 
-from Numeric import concatenate, array, Float, shape, reshape, ravel, take
+from Numeric import concatenate, array, Float, shape, reshape, ravel, take, \
+                        size, shape
+from RandomArray import randint
+#from MA import tolist
 
 from anuga.utilities.numerical_tools import ensure_numeric
@@ -223 +226 @@
           regarded as belonging to the polygon (closed = True)
           or not (closed = False). Default is True.
-          
+           
         Output
           New geospatial data object representing points inside
@@ -239 +242 @@
         inside_indices = inside_polygon(points, polygon, closed)
 
-        clipped_points = take(points, inside_indices)
+        clipped_G = self.get_sample(inside_indices)
+#        clipped_points = take(points, inside_indices)
 
         # Clip all attributes
-        attributes = self.get_all_attributes()
-
-        clipped_attributes = {}
-        if attributes is not None:
-            for key, att in attributes.items():
-                clipped_attributes[key] = take(att, inside_indices)
-
-        return Geospatial_data(clipped_points, clipped_attributes)
+#        attributes = self.get_all_attributes()
+
+#        clipped_attributes = {}
+#        if attributes is not None:
+#            for key, att in attributes.items():
+#                clipped_attributes[key] = take(att, inside_indices)
+
+#        return Geospatial_data(clipped_points, clipped_attributes)
+        return clipped_G
         
         
@@ -276 +281 @@
         outside_indices = outside_polygon(points, polygon, closed)
 
-        clipped_points = take(points, outside_indices)
+        clipped_G = self.get_sample(outside_indices)
+
+#        clipped_points = take(points, outside_indices)
 
         # Clip all attributes
-        attributes = self.get_all_attributes()
-
-        clipped_attributes = {}
-        if attributes is not None:
-            for key, att in attributes.items():
-                clipped_attributes[key] = take(att, outside_indices)
-
-        return Geospatial_data(clipped_points, clipped_attributes)
+#        attributes = self.get_all_attributes()
+
+#        clipped_attributes = {}
+#        if attributes is not None:
+#            for key, att in attributes.items():
+#                clipped_attributes[key] = take(att, outside_indices)
+
+#        return Geospatial_data(clipped_points, clipped_attributes)
+        return clipped_G
 
     
@@ -565 +573 @@
             msg = 'Unknown file type %s ' %file_name
             raise IOError, msg 
-    
-
+        
+    def get_sample(self, indices):
+        """ Returns a object which is a subset of the original
+        and the data points and attributes in this new object refer to
+        the indices provided
+        
+        Input
+            indices- a list of integers that represent the new object
+        Output
+            New geospatial data object representing points specified by 
+            the indices 
+            """
+        #FIXME: add the geo_reference to this
+        
+        points = self.get_data_points()
+        sampled_points = take(points, indices)
+
+        attributes = self.get_all_attributes()
+
+        sampled_attributes = {}
+        if attributes is not None:
+            for key, att in attributes.items():
+                sampled_attributes[key] = take(att, indices)
+
+        return Geospatial_data(sampled_points, sampled_attributes)
+    
+    
+    def split(self, factor=0.5):
+        """Returns two geospatial_data object, first is size of the 'factor'
+        smaller the original one and the second is the remainer. The two new 
+        object are disjoin set of each other. 
+        
+        Points of the two new object have selected RANDOMLY. 
+        AND if factor is a decimal it will round (2.25 to 2 and 2.5 to 3)
+        
+        
+        Input - the factor which to split the object, if 0.1 then 10% of the
+            object will be returned
+        
+        Output - two geospatial_data objects that are disjoint sets of the 
+            original
+            """
+        
+        i=0
+        self_size = len(self)
+#        print 'size', self_size
+        random_list = []
+        remainder_list = []
+        new_size = round(factor*self_size)
+   #     print'Split original %s by %s' %(self_size, factor)
+   #     print'New samples are %s and %s in size' %(int(round(factor*self_size)),int(self_size-new_size))
+        
+        #find unique random numbers
+        while i < new_size:
+            random_num = randint(0,self_size)
+            if random_num not in random_list:
+                random_list.append(random_num)
+                i=i+1
+
+        #Make list of opposite to random_list
+        for i in range(0,self_size,1):
+            remainder_list.append(i)
+#        print 'start remainer',remainder_list
+
+        #need to sort and reverse so the pop() works correctly
+        random_list.sort()
+        random_list.reverse()
+        for i in random_list:
+            remainder_list.pop(i)
+
+        #get new samples
+        G1 = self.get_sample(random_list)
+        G2 = self.get_sample(remainder_list)
+
+        return G1, G2
 
 
@@ -830 +911 @@
     return numbers
             
-def xxx_add_points_files(add_file1, add_file2, results_file):
-    """ adds the points and attruibutes of 2 pts or xya files and
-    writes it to a pts file
-    
-    NOTE will add the function to check and remove points from one set
-    that are shared. This will require some work and maybe __subtract__ function 
-    """
-    
-    G1 = Geospatial_data(file_name = add_file1)
-    G2 = Geospatial_data(file_name = add_file2)
-    new_add_file2 = add_file2[:-4] + '.pts' 
-
-    G = G1 + G2
-    
-    #FIXME remove dependance on points to dict in export only!
-#    G_points_dict = geospatial_data2points_dictionary(G)
-#    export_points_file(results_file, G_points_dict)
-
-#    G_points_dict = geospatial_data2points_dictionary(G)
-
-    G.export_points_file(results_file)
-
-# '
 def ensure_absolute(points, geo_reference = None):
     """
@@ -924 +982 @@
     points = Geospatial_data(data_points=points, geo_reference=geo)        
     return points
+
+#def file2xya(filename):
+    
+#    G = Geospatial_data(filename)
+#    G.export_points_file()
+
              

anuga_core/source/anuga/geospatial_data/test_geospatial_data.py

@@ -1627 +1627 @@
         G = Geospatial_data(points)
         self.failUnless(4 ==len(G),  'Len error!')
+        
+    def test_split(self):
+        """test if the results from spilt are disjoin sets"""
+        
+        points = [[1.0, 1.0], [1.0, 2.0],[1.0, 3.0], [1.0, 4.0], [1.0, 5.0],
+                  [2.0, 1.0], [2.0, 2.0],[2.0, 3.0], [2.0, 4.0], [2.0, 5.0],
+                  [3.0, 1.0], [3.0, 2.0],[3.0, 3.0], [3.0, 4.0], [3.0, 5.0],
+                  [4.0, 1.0], [4.0, 2.0],[4.0, 3.0], [4.0, 4.0], [4.0, 5.0],
+                  [5.0, 1.0], [5.0, 2.0],[5.0, 3.0], [5.0, 4.0], [5.0, 5.0]]
+        attributes = {'depth':[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 
+                      14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25],'speed':[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 
+                      14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25]}
+        G = Geospatial_data(points, attributes)
+#        print G.get_data_points()
+#        print G.get_attributes()
+
+        factor = 0.22
+
+        G1, G2  = G.split(factor) #will return G1 with 10% for points and G2 with 90%
+        
+#        print 'len(G): %s  len(G1): %s len(G2): %s' %(len(G), len(G1), len(G2))
+#        print 'G: ', len(G),'G1: ', len(G1), 'G2: ', len(G2)
+
+        assert allclose(len(G), len(G1)+len(G2))
+        assert allclose(round(len(G)*factor), len(G1))
+        
+#        assert allclose(G == G1 + G2) must implentent __equal__
+        
          
 if __name__ == "__main__":
 
-    #suite = unittest.makeSuite(Test_Geospatial_data, 'test_ensure_geospatial')
+#    suite = unittest.makeSuite(Test_Geospatial_data, 'test_split')
+#    suite = unittest.makeSuite(Test_Geospatial_data, 'test_clip0')
     suite = unittest.makeSuite(Test_Geospatial_data, 'test')
     runner = unittest.TextTestRunner()