Changeset 3969
- Timestamp:
- Nov 13, 2006, 1:53:33 PM (18 years ago)
- Location:
- anuga_core/source/anuga/geospatial_data
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
anuga_core/source/anuga/geospatial_data/geospatial_data.py
r3745 r3969 7 7 from types import DictType 8 8 9 from Numeric import concatenate, array, Float, shape, reshape, ravel, take 9 from Numeric import concatenate, array, Float, shape, reshape, ravel, take, \ 10 size, shape 11 from RandomArray import randint 12 #from MA import tolist 10 13 11 14 from anuga.utilities.numerical_tools import ensure_numeric … … 223 226 regarded as belonging to the polygon (closed = True) 224 227 or not (closed = False). Default is True. 225 228 226 229 Output 227 230 New geospatial data object representing points inside … … 239 242 inside_indices = inside_polygon(points, polygon, closed) 240 243 241 clipped_points = take(points, inside_indices) 244 clipped_G = self.get_sample(inside_indices) 245 # clipped_points = take(points, inside_indices) 242 246 243 247 # Clip all attributes 244 attributes = self.get_all_attributes() 245 246 clipped_attributes = {} 247 if attributes is not None: 248 for key, att in attributes.items(): 249 clipped_attributes[key] = take(att, inside_indices) 250 251 return Geospatial_data(clipped_points, clipped_attributes) 248 # attributes = self.get_all_attributes() 249 250 # clipped_attributes = {} 251 # if attributes is not None: 252 # for key, att in attributes.items(): 253 # clipped_attributes[key] = take(att, inside_indices) 254 255 # return Geospatial_data(clipped_points, clipped_attributes) 256 return clipped_G 252 257 253 258 … … 276 281 outside_indices = outside_polygon(points, polygon, closed) 277 282 278 clipped_points = take(points, outside_indices) 283 clipped_G = self.get_sample(outside_indices) 284 285 # clipped_points = take(points, outside_indices) 279 286 280 287 # Clip all attributes 281 attributes = self.get_all_attributes() 282 283 clipped_attributes = {} 284 if attributes is not None: 285 for key, att in attributes.items(): 286 clipped_attributes[key] = take(att, outside_indices) 287 288 return Geospatial_data(clipped_points, clipped_attributes) 288 # attributes = self.get_all_attributes() 289 290 # clipped_attributes = {} 291 # if attributes is not None: 292 # for key, att in attributes.items(): 293 # clipped_attributes[key] = take(att, outside_indices) 294 295 # return Geospatial_data(clipped_points, clipped_attributes) 296 return clipped_G 289 297 290 298 … … 565 573 msg = 'Unknown file type %s ' %file_name 566 574 raise IOError, msg 567 568 575 576 def get_sample(self, indices): 577 """ Returns a object which is a subset of the original 578 and the data points and attributes in this new object refer to 579 the indices provided 580 581 Input 582 indices- a list of integers that represent the new object 583 Output 584 New geospatial data object representing points specified by 585 the indices 586 """ 587 #FIXME: add the geo_reference to this 588 589 points = self.get_data_points() 590 sampled_points = take(points, indices) 591 592 attributes = self.get_all_attributes() 593 594 sampled_attributes = {} 595 if attributes is not None: 596 for key, att in attributes.items(): 597 sampled_attributes[key] = take(att, indices) 598 599 return Geospatial_data(sampled_points, sampled_attributes) 600 601 602 def split(self, factor=0.5): 603 """Returns two geospatial_data object, first is size of the 'factor' 604 smaller the original one and the second is the remainer. The two new 605 object are disjoin set of each other. 606 607 Points of the two new object have selected RANDOMLY. 608 AND if factor is a decimal it will round (2.25 to 2 and 2.5 to 3) 609 610 611 Input - the factor which to split the object, if 0.1 then 10% of the 612 object will be returned 613 614 Output - two geospatial_data objects that are disjoint sets of the 615 original 616 """ 617 618 i=0 619 self_size = len(self) 620 # print 'size', self_size 621 random_list = [] 622 remainder_list = [] 623 new_size = round(factor*self_size) 624 # print'Split original %s by %s' %(self_size, factor) 625 # print'New samples are %s and %s in size' %(int(round(factor*self_size)),int(self_size-new_size)) 626 627 #find unique random numbers 628 while i < new_size: 629 random_num = randint(0,self_size) 630 if random_num not in random_list: 631 random_list.append(random_num) 632 i=i+1 633 634 #Make list of opposite to random_list 635 for i in range(0,self_size,1): 636 remainder_list.append(i) 637 # print 'start remainer',remainder_list 638 639 #need to sort and reverse so the pop() works correctly 640 random_list.sort() 641 random_list.reverse() 642 for i in random_list: 643 remainder_list.pop(i) 644 645 #get new samples 646 G1 = self.get_sample(random_list) 647 G2 = self.get_sample(remainder_list) 648 649 return G1, G2 569 650 570 651 … … 830 911 return numbers 831 912 832 def xxx_add_points_files(add_file1, add_file2, results_file):833 """ adds the points and attruibutes of 2 pts or xya files and834 writes it to a pts file835 836 NOTE will add the function to check and remove points from one set837 that are shared. This will require some work and maybe __subtract__ function838 """839 840 G1 = Geospatial_data(file_name = add_file1)841 G2 = Geospatial_data(file_name = add_file2)842 new_add_file2 = add_file2[:-4] + '.pts'843 844 G = G1 + G2845 846 #FIXME remove dependance on points to dict in export only!847 # G_points_dict = geospatial_data2points_dictionary(G)848 # export_points_file(results_file, G_points_dict)849 850 # G_points_dict = geospatial_data2points_dictionary(G)851 852 G.export_points_file(results_file)853 854 # '855 913 def ensure_absolute(points, geo_reference = None): 856 914 """ … … 924 982 points = Geospatial_data(data_points=points, geo_reference=geo) 925 983 return points 984 985 #def file2xya(filename): 986 987 # G = Geospatial_data(filename) 988 # G.export_points_file() 989 926 990 -
anuga_core/source/anuga/geospatial_data/test_geospatial_data.py
r3797 r3969 1627 1627 G = Geospatial_data(points) 1628 1628 self.failUnless(4 ==len(G), 'Len error!') 1629 1630 def test_split(self): 1631 """test if the results from spilt are disjoin sets""" 1632 1633 points = [[1.0, 1.0], [1.0, 2.0],[1.0, 3.0], [1.0, 4.0], [1.0, 5.0], 1634 [2.0, 1.0], [2.0, 2.0],[2.0, 3.0], [2.0, 4.0], [2.0, 5.0], 1635 [3.0, 1.0], [3.0, 2.0],[3.0, 3.0], [3.0, 4.0], [3.0, 5.0], 1636 [4.0, 1.0], [4.0, 2.0],[4.0, 3.0], [4.0, 4.0], [4.0, 5.0], 1637 [5.0, 1.0], [5.0, 2.0],[5.0, 3.0], [5.0, 4.0], [5.0, 5.0]] 1638 attributes = {'depth':[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 1639 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, 1640 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25]} 1641 G = Geospatial_data(points, attributes) 1642 # print G.get_data_points() 1643 # print G.get_attributes() 1644 1645 factor = 0.22 1646 1647 G1, G2 = G.split(factor) #will return G1 with 10% for points and G2 with 90% 1648 1649 # print 'len(G): %s len(G1): %s len(G2): %s' %(len(G), len(G1), len(G2)) 1650 # print 'G: ', len(G),'G1: ', len(G1), 'G2: ', len(G2) 1651 1652 assert allclose(len(G), len(G1)+len(G2)) 1653 assert allclose(round(len(G)*factor), len(G1)) 1654 1655 # assert allclose(G == G1 + G2) must implentent __equal__ 1656 1629 1657 1630 1658 if __name__ == "__main__": 1631 1659 1632 #suite = unittest.makeSuite(Test_Geospatial_data, 'test_ensure_geospatial') 1660 # suite = unittest.makeSuite(Test_Geospatial_data, 'test_split') 1661 # suite = unittest.makeSuite(Test_Geospatial_data, 'test_clip0') 1633 1662 suite = unittest.makeSuite(Test_Geospatial_data, 'test') 1634 1663 runner = unittest.TextTestRunner()
Note: See TracChangeset
for help on using the changeset viewer.