Changeset 4202


Ignore:
Timestamp:
Jan 31, 2007, 3:20:39 PM (17 years ago)
Author:
duncan
Message:

sigh. the 'bug' that I fixed in 4199 wasn't a bug. It was correct. Reverting by doing svn merge -r4199:4198 https://datamining.anu.edu.au/svn/ga/ in inundation dir

Location:
anuga_core/source/anuga
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • anuga_core/source/anuga/abstract_2d_finite_volumes/quantity.py

    r4166 r4202  
    707707                                      compression=False)
    708708        else:
    709 
    710709            vertex_attributes = apply(fit_to_mesh,
    711710                                      args, kwargs)
     
    714713        self.set_values_from_array(vertex_attributes,
    715714                                   location, indices, verbose)
    716 
    717 
    718715
    719716    def set_values_from_file(self, filename, attribute_name, alpha,
     
    747744                                        verbose=verbose)
    748745                                #, max_read_lines=max_read_lines)
    749        
    750746        #Call underlying method using array values
    751747        self.set_values_from_array(vertex_attributes,
  • anuga_core/source/anuga/abstract_2d_finite_volumes/test_quantity.py

    r4198 r4202  
    4343        self.mesh4.check_integrity()
    4444
     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       
    4556    def tearDown(self):
    4657        pass
     
    598609        os.remove(ptsfile)
    599610
     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       
    600705    def test_set_values_from_file_with_georef1(self):
    601706
     
    17631868    suite = unittest.makeSuite(Test_Quantity, 'test')
    17641869
    1765     #suite = unittest.makeSuite(Test_Quantity, 'Cache_cache_test_set_values_from_file')
     1870    suite = unittest.makeSuite(Test_Quantity, 'test_set_values_from_UTM')
    17661871    #print "restricted test"
    17671872    #suite = unittest.makeSuite(Test_Quantity,'test_set_values_from_file_with_georef2')
  • anuga_core/source/anuga/coordinate_transforms/redfearn.py

    r4199 r4202  
    199199       
    200200    for point in points:
     201       
    201202        zone, easting, northing = redfearn(float(point[0]),
    202203                                           float(point[1]),
     
    205206        new_geo = Geo_reference(zone)
    206207        old_geo.reconcile_zones(new_geo)       
    207         utm_points.append([northing, easting])
     208        utm_points.append([easting, northing])
    208209
    209210    return utm_points, old_geo.get_zone()
  • anuga_core/source/anuga/coordinate_transforms/test_redfearn.py

    r4199 r4202  
    158158        points, zone = convert_from_latlon_to_utm(latitudes=lats, longitudes=longs)
    159159
    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)
    164164        self.failUnless(zone == 56,
    165165                        'Bad zone error!')
     
    230230        points, zone = convert_from_latlon_to_utm(points=points)
    231231        #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)
    237236        self.failUnless(zone == 56,
    238237                        'Bad zone error!')
  • anuga_core/source/anuga/fit_interpolate/test_fit.py

    r4185 r4202  
    241241        os.remove(fileName)
    242242
    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       
    244281    def test_fit_to_mesh_pts(self):
    245282        a = [-1.0, 0.0]
     
    692729if __name__ == "__main__":
    693730    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')
    695732    #suite = unittest.makeSuite(Test_Fit,'test_smooth_attributes_to_mesh_one_point')
    696733    runner = unittest.TextTestRunner(verbosity=1)
  • anuga_core/source/anuga/geospatial_data/geospatial_data.py

    r4199 r4202  
    1616   
    1717from 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
     18from anuga.coordinate_transforms.geo_reference import Geo_reference, TitleError
    2219from anuga.coordinate_transforms.redfearn import convert_from_latlon_to_utm
    2320from anuga.utilities.anuga_exceptions import ANUGAError
     
    335332        return self.geo_reference
    336333       
    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
    342336
    343337        If absolute is False returned coordinates are relative to the
     
    350344        Default: absolute is True.
    351345        """
    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
    371347        if absolute is True and geo_reference is None:
    372348            return self.geo_reference.get_absolute(self.data_points)
     
    377353        else:
    378354            return self.data_points
     355       
    379356   
    380357    def get_attributes(self, attribute_name=None):
     
    845822    data_points, zone  = convert_from_latlon_to_utm(latitudes=latitudes,
    846823                                                    longitudes=longitudes)
    847    
    848824    return data_points, Geo_reference(zone=zone)
    849825   
  • anuga_core/source/anuga/geospatial_data/test_geospatial_data.py

    r4199 r4202  
    162162        assert allclose(results, points_rel)
    163163
    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       
    195165    def test_set_geo_reference(self):
    196166        points_ab = [[12.5,34.7],[-4.5,-60.0]]
     
    710680                        [0.0, 10.0, 40.0])
    711681
    712     def test_load_csv(self):
    713        
    714         import os
    715         import tempfile
    716        
    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 fileName
    725         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        
    735682       
    736683    def test_load_xya(self):
     
    18981845        fileName = tempfile.mktemp(".csv")
    18991846        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\
     1848150.916666667,-34.50,452.688000, 10\n\
     1849150.0,-34,459.126000, 10\n")
    19031850        file.close()
    19041851        results = Geospatial_data(fileName, delimiter=',')
     
    19061853        points = results.get_data_points()
    19071854       
    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)
    19131859       
    19141860     
     
    19281874        points = results.get_data_points()
    19291875       
    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
    19351881         
    19361882    def test_load_csv_lat_long_bad(self):
     
    19561902       
    19571903    def test_lat_long(self):
    1958         #Zone:   56   
    1959         #Easting:  308728.009  Northing: 6180432.601
    1960         #Latitude:   -34  30 ' 0.00000 ''  Longitude: 150  55 ' 0.00000 ''
    1961 
    19621904        lat_gong = degminsec2decimal_degrees(-34,30,0.)
    19631905        lon_gong = degminsec2decimal_degrees(150,55,0.)
     
    19721914        points = gsd.get_data_points(absolute=True)
    19731915       
    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)
    19781920        self.failUnless(gsd.get_geo_reference().get_zone() == 56,
    19791921                        'Bad zone error!')
     
    20251967        points = gsd.get_data_points(absolute=True)
    20261968       
    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)
    20321973        self.failUnless(gsd.get_geo_reference().get_zone() == 56,
    20331974                        'Bad zone error!')
  • anuga_core/source/anuga/shallow_water/test_data_manager.py

    r4199 r4202  
    45324532        points = gsd.get_data_points(absolute=True)
    45334533       
    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)
    45404538        self.failUnless(gsd.get_geo_reference().get_zone() == 56,
    45414539                        'Bad zone error!')
Note: See TracChangeset for help on using the changeset viewer.