Changeset 3900
- Timestamp:
- Nov 1, 2006, 10:58:50 AM (18 years ago)
- Location:
- anuga_core/source/anuga
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
anuga_core/source/anuga/abstract_2d_finite_volumes/util.py
r3850 r3900 13 13 14 14 def file_function(filename, 15 domain = None, 16 quantities = None, 17 interpolation_points = None, 18 verbose = False, 19 use_cache = False): 15 domain=None, 16 quantities=None, 17 interpolation_points=None, 18 time_thinning=1, 19 verbose=False, 20 use_cache=False): 20 21 """Read time history of spatial data from NetCDF file and return 21 22 a callable object. … … 58 59 59 60 61 # Build arguments and keyword arguments for use with caching or apply. 62 args = (filename,) 63 64 kwargs = {'domain': domain, 65 'quantities': quantities, 66 'interpolation_points': interpolation_points, 67 'time_thinning': time_thinning, 68 'verbose': verbose} 69 70 71 # Call underlying engine with or without caching 60 72 if use_cache is True: 61 73 try: … … 66 78 raise msg 67 79 68 69 80 f = cache(_file_function, 70 filename, 71 {'domain': domain, 72 'quantities': quantities, 73 'interpolation_points': interpolation_points, 74 'verbose': verbose}, 75 dependencies = [filename], 76 compression = False, 77 verbose = verbose) 78 #FIXME (Ole): Pass cache arguments, such as compression, in some sort of 79 #structure 81 args, kwargs, 82 dependencies=[filename], 83 compression=False, 84 verbose=verbose) 85 86 else: 87 f = apply(_file_function, 88 args, kwargs) 89 90 91 #FIXME (Ole): Pass cache arguments, such as compression, in some sort of 92 #structure 80 93 81 else:82 f = _file_function(filename,83 domain,84 quantities,85 interpolation_points,86 verbose)87 94 88 95 return f … … 91 98 92 99 def _file_function(filename, 93 domain = None, 94 quantities = None, 95 interpolation_points = None, 96 verbose = False): 100 domain=None, 101 quantities=None, 102 interpolation_points=None, 103 time_thinning=1, 104 verbose=False): 97 105 """Internal function 98 106 … … 133 141 return get_netcdf_file_function(filename, domain, quantities, 134 142 interpolation_points, 143 time_thinning=time_thinning, 135 144 verbose = verbose) 136 145 else: … … 143 152 quantity_names=None, 144 153 interpolation_points=None, 154 time_thinning=1, 145 155 verbose = False): 146 156 """Read time history of spatial data from NetCDF sww file and … … 238 248 time = fid.variables['time'][:] 239 249 250 # Get time independent stuff 240 251 if spatial: 241 252 #Get origin … … 296 307 297 308 298 #Produce values for desired data points at 299 #each timestep 300 309 # Produce values for desired data points at 310 # each timestep for each quantity 301 311 quantities = {} 302 312 for i, name in enumerate(quantity_names): … … 318 328 triangles, 319 329 interpolation_points, 330 time_thinning=time_thinning, 320 331 verbose=verbose) 321 332 -
anuga_core/source/anuga/fit_interpolate/interpolate.py
r3896 r3900 405 405 FIXME (Ole): Need to allow vertex coordinates and interpolation points to be 406 406 geospatial data objects 407 408 Time assumed to be relative to starttime 407 409 """ 408 410 … … 415 417 triangles=None, 416 418 interpolation_points=None, 419 time_thinning=1, 417 420 verbose=False): 418 421 """Initialise object and build spatial interpolation if required 422 423 Time_thinning_number controls how many timesteps to use. Only timesteps with 424 index%time_thinning_number == 0 will used, or in other words a value of 3, say, 425 will cause the algorithm to use every third time step. 419 426 """ 420 427 … … 458 465 self.spatial = True 459 466 467 # Thin timesteps if needed 468 # Note array() is used to make the thinned arrays contiguous in memory 469 self.time = array(time[::time_thinning]) 470 for name in quantity_names: 471 if len(quantities[name].shape) == 2: 472 quantities[name] = array(quantities[name][::time_thinning,:]) 460 473 461 474 # Save for use with statistics … … 468 481 self.vertex_coordinates = vertex_coordinates 469 482 self.interpolation_points = interpolation_points 470 self.time = time[:] # Time assumed to be relative to starttime 483 484 471 485 self.index = 0 # Initial time index 472 486 self.precomputed_values = {} … … 495 509 496 510 # Build interpolator 497 if verbose: print 'Building interpolation matrix' 511 if verbose: 512 msg = 'Building interpolation matrix from source mesh ' 513 msg += '(%d vertices, %d triangles)' %(vertex_coordinates.shape[0], 514 triangles.shape[0]) 515 print msg 516 517 498 518 interpol = Interpolate(vertex_coordinates, 499 519 triangles, 500 #point_coordinates = \501 #self.interpolation_points,502 #alpha = 0,503 520 verbose=verbose) 504 521 505 if verbose: print 'Interpolate' 522 if verbose: 523 print 'Interpolating (%d interpolation points, %d timesteps).'\ 524 %(self.interpolation_points.shape[0], self.time.shape[0]), 525 526 if time_thinning > 1: 527 print 'Timesteps were thinned by a factor of %d' %time_thinning 528 else: 529 print 530 506 531 for i, t in enumerate(self.time): 507 532 # Interpolate quantities at this timestep -
anuga_core/source/anuga/fit_interpolate/test_interpolate.py
r3846 r3900 1124 1124 else: 1125 1125 raise 'Should raise exception' 1126 1127 1128 1129 def test_interpolation_interface_with_time_thinning(self): 1130 # Test spatio-temporal interpolation 1131 # Test that spatio temporal function performs the correct 1132 # interpolations in both time and space 1133 1134 #Three timesteps 1135 time = [1.0, 2.0, 4.0, 5.0, 7.0, 8.0, 9.0, 10.0] 1136 1137 #Setup mesh used to represent fitted function 1138 a = [0.0, 0.0] 1139 b = [0.0, 2.0] 1140 c = [2.0, 0.0] 1141 d = [0.0, 4.0] 1142 e = [2.0, 2.0] 1143 f = [4.0, 0.0] 1144 1145 points = [a, b, c, d, e, f] 1146 #bac, bce, ecf, dbe 1147 triangles = [[1,0,2], [1,2,4], [4,2,5], [3,1,4]] 1148 1149 1150 #New datapoints where interpolated values are sought 1151 interpolation_points = [[ 0.0, 0.0], 1152 [ 0.5, 0.5], 1153 [ 0.7, 0.7], 1154 [ 1.0, 0.5], 1155 [ 2.0, 0.4], 1156 [ 2.8, 1.2]] 1157 1158 #One quantity 1159 Q = zeros( (8,6), Float ) 1160 1161 #Linear in time and space 1162 for i, t in enumerate(time): 1163 Q[i, :] = t*linear_function(points) 1164 1165 # Check interpolation of one quantity using interpolaton points) using default 1166 # time_thinning of 1 1167 I = Interpolation_function(time, Q, 1168 vertex_coordinates=points, 1169 triangles=triangles, 1170 interpolation_points=interpolation_points, 1171 verbose=False) 1172 1173 answer = linear_function(interpolation_points) 1174 1175 1176 t = time[0] 1177 for j in range(50): #t in [1, 6] 1178 for id in range(len(interpolation_points)): 1179 assert allclose(I(t, id), t*answer[id]) 1180 t += 0.1 1181 1182 1183 # Now check time_thinning 1184 I = Interpolation_function(time, Q, 1185 vertex_coordinates=points, 1186 triangles=triangles, 1187 interpolation_points=interpolation_points, 1188 time_thinning=2, 1189 verbose=False) 1190 1191 1192 assert len(I.time) == 4 1193 assert( allclose(I.time, [1.0, 4.0, 7.0, 9.0] )) 1194 1195 answer = linear_function(interpolation_points) 1196 1197 t = time[0] 1198 for j in range(50): #t in [1, 6] 1199 for id in range(len(interpolation_points)): 1200 assert allclose(I(t, id), t*answer[id]) 1201 t += 0.1 1202 1203 1126 1204 1127 1205
Note: See TracChangeset
for help on using the changeset viewer.