Changeset 2690
 Timestamp:
 Apr 11, 2006, 1:13:32 PM (18 years ago)
 Location:
 inundation/fit_interpolate
 Files:

 2 edited
Legend:
 Unmodified
 Added
 Removed

inundation/fit_interpolate/interpolate.py
r2684 r2690 39 39 vertex_coordinates, 40 40 triangles, 41 mesh_origin =None,41 mesh_origin=None, 42 42 verbose=False, 43 43 max_vertices_per_cell=30): … … 64 64 max_vertices_per_cell: Number of vertices in a quad tree cell 65 65 at which the cell is split into 4. 66 67 """66 """ 67 68 68 # Initialise variabels 69 69 self._A_can_be_reused = False … … 72 72 #Convert input to Numeric arrays 73 73 triangles = ensure_numeric(triangles, Int) 74 vertex_coordinates = ensure_numeric(vertex_coordinates, Float) 74 75 if isinstance(vertex_coordinates,Geospatial_data): 76 vertex_coordinates = vertex_coordinates.get_data_points( \ 77 absolute = True) 78 msg = "Use a Geospatial_data object or a meshorigin. Not both." 79 assert mesh_origin == None, msg 80 81 else: 82 vertex_coordinates = ensure_numeric(vertex_coordinates, Float) 75 83 #Build underlying mesh 76 84 if verbose: print 'Building mesh' … … 78 86 #FIXME: Trying the normal mesh while testing precrop, 79 87 # The functionality of boundary_polygon is needed for that 80 81 #FIXME: geo_ref can also be a geo_ref object82 #FIXME: move this to interpolate_block83 84 #Note: not so keen to use geospacial to represent the mesh,85 # since the data structure is86 # points, georef and triangles, rather than just points and georef87 # this info will usually be coming from a domain instance.88 88 89 89 if mesh_origin is None: 90 geo = None 90 geo = None #Geo_reference() 91 91 else: 92 geo = Geo_reference(mesh_origin[0],mesh_origin[1],mesh_origin[2]) 93 self.mesh = Mesh(vertex_coordinates, triangles, 94 geo_reference = geo) 95 92 if isinstance(mesh_origin, Geo_reference): 93 geo = mesh_origin 94 else: 95 geo = Geo_reference(mesh_origin[0], 96 mesh_origin[1], 97 mesh_origin[2]) 98 vertex_coordinates = geo.get_absolute(vertex_coordinates) 99 #Don't pass geo_reference to mesh. It doesn't work. 100 self.mesh = Mesh(vertex_coordinates, triangles) 96 101 self.mesh.check_integrity() 97 102 self.root = build_quadtree(self.mesh, … … 128 133 """ 129 134 #print "point_coordinates interpolate.interpolate",point_coordinates 130 # Can I interpolate, based on previous point_coordinates?131 135 if isinstance(point_coordinates,Geospatial_data): 132 136 point_coordinates = point_coordinates.get_data_points( \ 133 137 absolute = True) 134 138 139 # Can I interpolate, based on previous point_coordinates? 135 140 if point_coordinates is None: 136 141 if self._A_can_be_reused is True and \ 
inundation/fit_interpolate/test_interpolate.py
r2684 r2690 191 191 [0.5, 0.0, 0.5], #Affects vertex 0 and 2 192 192 [0.0, 0.5, 0.5]] 193 interp = Interpolate(points, vertices , data)193 interp = Interpolate(points, vertices) 194 194 195 195 assert allclose(interp._build_interpolation_matrix_A(data).todense(), … … 212 212 [0.0, 0.25, 0.75]] 213 213 214 interp = Interpolate(points, vertices , data)214 interp = Interpolate(points, vertices) 215 215 216 216 assert allclose(interp._build_interpolation_matrix_A(data).todense(), … … 232 232 data = [ [0.2, 1.5], [0.123, 1.768], [1.43, 0.44] ] 233 233 234 interp = Interpolate(points, vertices , data)234 interp = Interpolate(points, vertices) 235 235 #print "interp.get_A()", interp.get_A() 236 236 results = interp._build_interpolation_matrix_A(data).todense() … … 252 252 data = [ [0.2, 1.5], [0.123, 1.768], [1.43, 0.44], [5.0, 7.0]] 253 253 254 interp = Interpolate(points, vertices , data)254 interp = Interpolate(points, vertices) 255 255 results = interp._build_interpolation_matrix_A(data).todense() 256 256 assert allclose(sum(results, axis=1), [1,1,1,0]) … … 293 293 294 294 assert allclose(A, answer) 295 296 #FIXME  tests are hopefully failing due to points outsdie the 297 # mesh 295 296 def test_geo_ref(self): 297 v0 = [0.0, 0.0] 298 v1 = [0.0, 5.0] 299 v2 = [5.0, 0.0] 300 301 vertices_absolute = [v0, v1, v2] 302 triangles = [ [1,0,2] ] #bac 303 304 geo = Geo_reference(57,100, 500) 305 306 vertices = geo.change_points_geo_ref(vertices_absolute) 307 #print "vertices",vertices 308 309 d0 = [1.0, 1.0] 310 d1 = [1.0, 2.0] 311 d2 = [3.0, 1.0] 312 point_coords = [ d0, d1, d2] 313 314 interp = Interpolate(vertices, triangles, mesh_origin=geo) 315 f = linear_function(vertices_absolute) 316 z = interp.interpolate(f, point_coords) 317 answer = linear_function(point_coords) 318 319 #print "z",z 320 #print "answer",answer 321 assert allclose(z, answer) 322 323 def test_Geospatial_verts(self): 324 v0 = [0.0, 0.0] 325 v1 = [0.0, 5.0] 326 v2 = [5.0, 0.0] 327 328 vertices_absolute = [v0, v1, v2] 329 triangles = [ [1,0,2] ] #bac 330 331 geo = Geo_reference(57,100, 500) 332 vertices = geo.change_points_geo_ref(vertices_absolute) 333 geopoints = Geospatial_data(vertices,geo_reference = geo) 334 #print "vertices",vertices 335 336 d0 = [1.0, 1.0] 337 d1 = [1.0, 2.0] 338 d2 = [3.0, 1.0] 339 point_coords = [ d0, d1, d2] 340 341 interp = Interpolate(geopoints, triangles) 342 f = linear_function(vertices_absolute) 343 z = interp.interpolate(f, point_coords) 344 answer = linear_function(point_coords) 345 346 #print "z",z 347 #print "answer",answer 348 assert allclose(z, answer) 349 298 350 def test_interpolate_attributes_to_points(self): 299 351 v0 = [0.0, 0.0] … … 309 361 point_coords = [ d0, d1, d2] 310 362 311 interp = Interpolate(vertices, triangles , point_coords)363 interp = Interpolate(vertices, triangles) 312 364 f = linear_function(vertices) 313 365 z = interp.interpolate(f, point_coords) 314 366 answer = linear_function(point_coords) 315 367 316 368 #print "z",z 369 #print "answer",answer 317 370 assert allclose(z, answer) 318 371
Note: See TracChangeset
for help on using the changeset viewer.