Changeset 4202
- Timestamp:
- Jan 31, 2007, 3:20:39 PM (18 years ago)
- Location:
- anuga_core/source/anuga
- Files:
-
- 8 edited
Legend:
- Unmodified
- Added
- Removed
-
anuga_core/source/anuga/abstract_2d_finite_volumes/quantity.py
r4166 r4202 707 707 compression=False) 708 708 else: 709 710 709 vertex_attributes = apply(fit_to_mesh, 711 710 args, kwargs) … … 714 713 self.set_values_from_array(vertex_attributes, 715 714 location, indices, verbose) 716 717 718 715 719 716 def set_values_from_file(self, filename, attribute_name, alpha, … … 747 744 verbose=verbose) 748 745 #, max_read_lines=max_read_lines) 749 750 746 #Call underlying method using array values 751 747 self.set_values_from_array(vertex_attributes, -
anuga_core/source/anuga/abstract_2d_finite_volumes/test_quantity.py
r4198 r4202 43 43 self.mesh4.check_integrity() 44 44 45 # UTM round Onslow 46 a = [240000, 7620000] 47 b = [240000, 7680000] 48 c = [300000, 7620000] 49 50 points = [a, b, c] 51 elements = [[0,2,1]] 52 53 self.mesh_onslow = Domain(points, elements) 54 self.mesh_onslow.check_integrity() 55 45 56 def tearDown(self): 46 57 pass … … 598 609 os.remove(ptsfile) 599 610 611 def test_set_values_from_lat_long(self): 612 quantity = Quantity(self.mesh_onslow) 613 614 #Get (enough) datapoints 615 data_points = [[-21.5, 114.5],[-21.4, 114.6],[-21.45,114.65], 616 [-21.35, 114.65],[-21.45, 114.55],[-21.45,114.6]] 617 618 data_geo_spatial = Geospatial_data(data_points, 619 points_are_lats_longs=True) 620 points_UTM = data_geo_spatial.get_data_points(absolute=True) 621 attributes = linear_function(points_UTM) 622 att = 'elevation' 623 624 #Create .txt file 625 txt_file = tempfile.mktemp(".txt") 626 file = open(txt_file,"w") 627 file.write(" lat,long," + att + " \n") 628 for data_point, attribute in map(None, data_points, attributes): 629 row = str(data_point[0]) + ',' + str(data_point[1]) \ 630 + ',' + str(attribute) 631 #print "row", row 632 file.write(row + "\n") 633 file.close() 634 635 636 #Check that values can be set from file 637 quantity.set_values(filename = txt_file, 638 attribute_name = att, alpha = 0) 639 answer = linear_function(quantity.domain.get_vertex_coordinates()) 640 641 #print "quantity.vertex_values.flat", quantity.vertex_values.flat 642 #print "answer",answer 643 644 assert allclose(quantity.vertex_values.flat, answer) 645 646 647 #Check that values can be set from file using default attribute 648 quantity.set_values(filename = txt_file, alpha = 0) 649 assert allclose(quantity.vertex_values.flat, answer) 650 651 #Cleanup 652 import os 653 os.remove(txt_file) 654 655 def test_set_values_from_UTM(self): 656 quantity = Quantity(self.mesh_onslow) 657 658 #Get (enough) datapoints 659 data_points = [[-21.5, 114.5],[-21.4, 114.6],[-21.45,114.65], 660 [-21.35, 114.65],[-21.45, 114.55],[-21.45,114.6]] 661 662 data_geo_spatial = Geospatial_data(data_points, 663 points_are_lats_longs=True) 664 points_UTM = data_geo_spatial.get_data_points(absolute=True) 665 attributes = linear_function(points_UTM) 666 att = 'elevation' 667 668 #Create .txt file 669 txt_file = tempfile.mktemp(".txt") 670 file = open(txt_file,"w") 671 file.write(" x,y," + att + " \n") 672 for data_point, attribute in map(None, points_UTM, attributes): 673 row = str(data_point[0]) + ',' + str(data_point[1]) \ 674 + ',' + str(attribute) 675 #print "row", row 676 file.write(row + "\n") 677 file.close() 678 679 680 #Check that values can be set from file 681 quantity.set_values_from_file(txt_file, att, 0, 682 'vertices', None) 683 answer = linear_function(quantity.domain.get_vertex_coordinates()) 684 #print "quantity.vertex_values.flat", quantity.vertex_values.flat 685 #print "answer",answer 686 assert allclose(quantity.vertex_values.flat, answer) 687 688 #Check that values can be set from file 689 quantity.set_values(filename = txt_file, 690 attribute_name = att, alpha = 0) 691 answer = linear_function(quantity.domain.get_vertex_coordinates()) 692 #print "quantity.vertex_values.flat", quantity.vertex_values.flat 693 #print "answer",answer 694 assert allclose(quantity.vertex_values.flat, answer) 695 696 697 #Check that values can be set from file using default attribute 698 quantity.set_values(filename = txt_file, alpha = 0) 699 assert allclose(quantity.vertex_values.flat, answer) 700 701 #Cleanup 702 import os 703 os.remove(txt_file) 704 600 705 def test_set_values_from_file_with_georef1(self): 601 706 … … 1763 1868 suite = unittest.makeSuite(Test_Quantity, 'test') 1764 1869 1765 #suite = unittest.makeSuite(Test_Quantity, 'Cache_cache_test_set_values_from_file')1870 suite = unittest.makeSuite(Test_Quantity, 'test_set_values_from_UTM') 1766 1871 #print "restricted test" 1767 1872 #suite = unittest.makeSuite(Test_Quantity,'test_set_values_from_file_with_georef2') -
anuga_core/source/anuga/coordinate_transforms/redfearn.py
r4199 r4202 199 199 200 200 for point in points: 201 201 202 zone, easting, northing = redfearn(float(point[0]), 202 203 float(point[1]), … … 205 206 new_geo = Geo_reference(zone) 206 207 old_geo.reconcile_zones(new_geo) 207 utm_points.append([ northing, easting])208 utm_points.append([easting, northing]) 208 209 209 210 return utm_points, old_geo.get_zone() -
anuga_core/source/anuga/coordinate_transforms/test_redfearn.py
r4199 r4202 158 158 points, zone = convert_from_latlon_to_utm(latitudes=lats, longitudes=longs) 159 159 160 assert allclose(points[0][0], 6180432.601)161 assert allclose(points[0][1], 308728.009)162 assert allclose(points[1][0], 6233785.284)163 assert allclose(points[1][1], 222908.705)160 assert allclose(points[0][0], 308728.009) 161 assert allclose(points[0][1], 6180432.601) 162 assert allclose(points[1][0], 222908.705) 163 assert allclose(points[1][1], 6233785.284) 164 164 self.failUnless(zone == 56, 165 165 'Bad zone error!') … … 230 230 points, zone = convert_from_latlon_to_utm(points=points) 231 231 #print "points",points 232 233 assert allclose(points[0][0], 6180432.601) 234 assert allclose(points[0][1], 308728.009) 235 assert allclose(points[1][0], 6233785.284) 236 assert allclose(points[1][1], 222908.705) 232 assert allclose(points[0][0], 308728.009) 233 assert allclose(points[0][1], 6180432.601) 234 assert allclose(points[1][0], 222908.705) 235 assert allclose(points[1][1], 6233785.284) 237 236 self.failUnless(zone == 56, 238 237 'Bad zone error!') -
anuga_core/source/anuga/fit_interpolate/test_fit.py
r4185 r4202 241 241 os.remove(fileName) 242 242 243 243 def test_fit_to_mesh_UTM_file(self): 244 #Get (enough) datapoints 245 data_points = [[-21.5, 114.5],[-21.4, 114.6],[-21.45,114.65], 246 [-21.35, 114.65],[-21.45, 114.55],[-21.45,114.6]] 247 248 data_geo_spatial = Geospatial_data(data_points, 249 points_are_lats_longs=True) 250 points_UTM = data_geo_spatial.get_data_points(absolute=True) 251 attributes = linear_function(points_UTM) 252 att = 'elevation' 253 254 #Create .txt file 255 txt_file = tempfile.mktemp(".txt") 256 file = open(txt_file,"w") 257 file.write(" x,y," + att + " \n") 258 for data_point, attribute in map(None, points_UTM, attributes): 259 row = str(data_point[0]) + ',' + str(data_point[1]) \ 260 + ',' + str(attribute) 261 #print "row", row 262 file.write(row + "\n") 263 file.close() 264 265 # setting up the mesh 266 a = [240000, 7620000] 267 b = [240000, 7680000] 268 c = [300000, 7620000] 269 points = [a, b, c] 270 elements = [[0,2,1]] 271 f = fit_to_mesh(points, elements, txt_file, 272 alpha=0.0, max_read_lines=2) 273 answer = linear_function(points) 274 #print "f",f 275 #print "answer",answer 276 assert allclose(f, answer) 277 278 # Delete file! 279 os.remove(txt_file) 280 244 281 def test_fit_to_mesh_pts(self): 245 282 a = [-1.0, 0.0] … … 692 729 if __name__ == "__main__": 693 730 suite = unittest.makeSuite(Test_Fit,'test') 694 #suite = unittest.makeSuite(Test_Fit,'test_smooth_att_to_mesh_with_excess_verts')731 suite = unittest.makeSuite(Test_Fit,'test_fit_to_mesh_UTM_file') 695 732 #suite = unittest.makeSuite(Test_Fit,'test_smooth_attributes_to_mesh_one_point') 696 733 runner = unittest.TextTestRunner(verbosity=1) -
anuga_core/source/anuga/geospatial_data/geospatial_data.py
r4199 r4202 16 16 17 17 from anuga.utilities.numerical_tools import ensure_numeric 18 19 from anuga.coordinate_transforms.lat_long_UTM_conversion import UTMtoLL 20 from anuga.coordinate_transforms.geo_reference import Geo_reference, \ 21 TitleError, DEFAULT_ZONE 18 from anuga.coordinate_transforms.geo_reference import Geo_reference, TitleError 22 19 from anuga.coordinate_transforms.redfearn import convert_from_latlon_to_utm 23 20 from anuga.utilities.anuga_exceptions import ANUGAError … … 335 332 return self.geo_reference 336 333 337 def get_data_points(self, absolute=True, geo_reference=None, 338 as_lat_long=False): 339 """Get coordinates for all data points as an Nx2 array. 340 Column 0 is x values 341 Column 1 is y values 334 def get_data_points(self, absolute=True, geo_reference=None): 335 """Get coordinates for all data points as an Nx2 array 342 336 343 337 If absolute is False returned coordinates are relative to the … … 350 344 Default: absolute is True. 351 345 """ 352 if as_lat_long is True: 353 msg = "Points need a zone to be converted into lats and longs" 354 assert self.geo_reference is not None, msg 355 zone = self.geo_reference.get_zone() 356 assert self.geo_reference.get_zone() is not DEFAULT_ZONE, msg 357 lats_longs = [] 358 for point in self.get_data_points(True): 359 ### UTMtoLL(northing, easting, zone, 360 print "point[0]",point[0] 361 print "point[1]",point[1] 362 print "zone",zone 363 results = UTMtoLL(point[0],point[1], zone) 364 print "results", results 365 lat_calced, long_calced = UTMtoLL(point[0],point[1], zone) 366 print "lat_calced", lat_calced 367 print "long_calced", long_calced 368 lats_longs.append(UTMtoLL(point[0],point[1], zone)) 369 370 346 371 347 if absolute is True and geo_reference is None: 372 348 return self.geo_reference.get_absolute(self.data_points) … … 377 353 else: 378 354 return self.data_points 355 379 356 380 357 def get_attributes(self, attribute_name=None): … … 845 822 data_points, zone = convert_from_latlon_to_utm(latitudes=latitudes, 846 823 longitudes=longitudes) 847 848 824 return data_points, Geo_reference(zone=zone) 849 825 -
anuga_core/source/anuga/geospatial_data/test_geospatial_data.py
r4199 r4202 162 162 assert allclose(results, points_rel) 163 163 164 165 def DSG_test_get_data_points_lat_long(self): 166 # lat long [-30.],[130] 167 #Zone: 52 168 #Easting: 596450.153 Northing: 6680793.777 169 # lat long [-32.],[131] 170 #Zone: 52 171 #Easting: 688927.638 Northing: 6457816.509 172 173 points_Lat_long = [[-30.,130], [-32,131]] 174 175 spatial = Geospatial_data(latitudes=[-30, -32.], 176 longitudes=[130, 131]) 177 178 results = spatial.get_data_points() 179 print "results UTM",results 180 181 results = spatial.get_data_points(as_lat_long=True) 182 print "results",results 183 assert allclose(results, points_rel) 184 185 x_p = -1770 186 y_p = 4.01 187 geo_ref = Geo_reference(56, x_p, y_p) 188 points_rel = geo_ref.change_points_geo_ref(points_ab) 189 results = spatial.get_data_points \ 190 ( geo_reference=geo_ref) 191 192 assert allclose(results, points_rel) 193 194 164 195 165 def test_set_geo_reference(self): 196 166 points_ab = [[12.5,34.7],[-4.5,-60.0]] … … 710 680 [0.0, 10.0, 40.0]) 711 681 712 def test_load_csv(self):713 714 import os715 import tempfile716 717 fileName = tempfile.mktemp(".csv")718 file = open(fileName,"w")719 file.write("lat,long,elevation speed \n\720 1.0 0.0 10.0 0.0\n\721 0.0 1.0 0.0 10.0\n\722 1.0 0.0 10.4 40.0\n")723 file.close()724 #print fileName725 results = Geospatial_data(fileName, delimiter=',')726 os.remove(fileName)727 # print 'data', results.get_data_points()728 assert allclose(results.get_data_points(), [[1.0, 0.0],[0.0, 1.0],729 [1.0, 0.0]])730 assert allclose(results.get_attributes(attribute_name='elevation'),731 [10.0, 0.0, 10.4])732 assert allclose(results.get_attributes(attribute_name='speed'),733 [0.0, 10.0, 40.0])734 735 682 736 683 def test_load_xya(self): … … 1898 1845 fileName = tempfile.mktemp(".csv") 1899 1846 file = open(fileName,"w") 1900 file.write("long,lat, z\n\1901 150.916666667,-34.50,452.688000 \n\1902 150.0,-34,459.126000 \n")1847 file.write("long,lat, elevation, yeah \n\ 1848 150.916666667,-34.50,452.688000, 10\n\ 1849 150.0,-34,459.126000, 10\n") 1903 1850 file.close() 1904 1851 results = Geospatial_data(fileName, delimiter=',') … … 1906 1853 points = results.get_data_points() 1907 1854 1908 1909 assert allclose(points[0][0], 6180432.601) 1910 assert allclose(points[0][1], 308728.009) 1911 assert allclose(points[1][0], 6233785.284) 1912 assert allclose(points[1][1], 222908.705) 1855 assert allclose(points[0][0], 308728.009) 1856 assert allclose(points[0][1], 6180432.601) 1857 assert allclose(points[1][0], 222908.705) 1858 assert allclose(points[1][1], 6233785.284) 1913 1859 1914 1860 … … 1928 1874 points = results.get_data_points() 1929 1875 1930 assert allclose(points[0][0], 6180432.601)1931 assert allclose(points[0][1], 308728.009)1932 assert allclose(points[1][0], 6233785.284)1933 assert allclose(points[1][1], 222908.705)1934 1876 assert allclose(points[0][0], 308728.009) 1877 assert allclose(points[0][1], 6180432.601) 1878 assert allclose(points[1][0], 222908.705) 1879 assert allclose(points[1][1], 6233785.284) 1880 1935 1881 1936 1882 def test_load_csv_lat_long_bad(self): … … 1956 1902 1957 1903 def test_lat_long(self): 1958 #Zone: 561959 #Easting: 308728.009 Northing: 6180432.6011960 #Latitude: -34 30 ' 0.00000 '' Longitude: 150 55 ' 0.00000 ''1961 1962 1904 lat_gong = degminsec2decimal_degrees(-34,30,0.) 1963 1905 lon_gong = degminsec2decimal_degrees(150,55,0.) … … 1972 1914 points = gsd.get_data_points(absolute=True) 1973 1915 1974 assert allclose(points[0][0], 6180432.601)1975 assert allclose(points[0][1], 308728.009)1976 assert allclose(points[1][0], 6233785.284)1977 assert allclose(points[1][1], 222908.705)1916 assert allclose(points[0][0], 308728.009) 1917 assert allclose(points[0][1], 6180432.601) 1918 assert allclose(points[1][0], 222908.705) 1919 assert allclose(points[1][1], 6233785.284) 1978 1920 self.failUnless(gsd.get_geo_reference().get_zone() == 56, 1979 1921 'Bad zone error!') … … 2025 1967 points = gsd.get_data_points(absolute=True) 2026 1968 2027 assert allclose(points[0][0], 6180432.601) 2028 assert allclose(points[0][1], 308728.009) 2029 assert allclose(points[1][0], 6233785.284) 2030 assert allclose(points[1][1], 222908.705) 2031 1969 assert allclose(points[0][0], 308728.009) 1970 assert allclose(points[0][1], 6180432.601) 1971 assert allclose(points[1][0], 222908.705) 1972 assert allclose(points[1][1], 6233785.284) 2032 1973 self.failUnless(gsd.get_geo_reference().get_zone() == 56, 2033 1974 'Bad zone error!') -
anuga_core/source/anuga/shallow_water/test_data_manager.py
r4199 r4202 4532 4532 points = gsd.get_data_points(absolute=True) 4533 4533 4534 4535 assert allclose(points[0][0], 6180432.601) 4536 assert allclose(points[0][1], 308728.009) 4537 assert allclose(points[1][0], 6233785.284) 4538 assert allclose(points[1][1], 222908.705) 4539 4534 assert allclose(points[0][0], 308728.009) 4535 assert allclose(points[0][1], 6180432.601) 4536 assert allclose(points[1][0], 222908.705) 4537 assert allclose(points[1][1], 6233785.284) 4540 4538 self.failUnless(gsd.get_geo_reference().get_zone() == 56, 4541 4539 'Bad zone error!')
Note: See TracChangeset
for help on using the changeset viewer.