Changeset 7574
- Timestamp:
- Dec 2, 2009, 5:35:02 PM (15 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
anuga_core/source/anuga/shallow_water/data_manager.py
r7569 r7574 2273 2273 quantity = 'elevation' 2274 2274 2275 # FIXME (Ole): I reckon we need a check so that 2276 # reduction and timestep cannot both be set. 2277 # if reduction is not None and timestep is not None: 2278 # msg = 'Either reduction or timestep must be set' 2279 # msg += ' but not both.' 2280 # raise Exception(msg) 2281 2282 2275 2283 if reduction is None: 2276 2284 reduction = max … … 2352 2360 % (min(times), max(times), len(times))) 2353 2361 log.critical(' Quantities [SI units]:') 2362 2354 2363 # Comment out for reduced memory consumption 2355 2364 for name in ['stage', 'xmomentum', 'ymomentum']: … … 2384 2393 2385 2394 for start_slice in xrange(0, number_of_points, block_size): 2386 # limit slice size to array end if at last block2395 # Limit slice size to array end if at last block 2387 2396 end_slice = min(start_slice + block_size, number_of_points) 2388 2397 2389 # get slices of all required variables2398 # Get slices of all required variables 2390 2399 q_dict = {} 2391 2400 for name in var_list: … … 2402 2411 new_res = num.zeros(res.shape[1], num.float) 2403 2412 for k in xrange(res.shape[1]): 2404 if reduction ==None:2413 if reduction is None: 2405 2414 new_res[k] = res[timestep,k] 2406 2415 else: … … 2410 2419 result[start_slice:end_slice] = res 2411 2420 2412 # Post condition: Now q has dimension: number_of_points2421 # Post condition: Now q has dimension: number_of_points 2413 2422 assert len(result.shape) == 1 2414 2423 assert result.shape[0] == number_of_points … … 2418 2427 % (quantity, min(result), max(result))) 2419 2428 2420 # Create grid and update xll/yll corner and x,y2421 # Relative extent2429 # Create grid and update xll/yll corner and x,y 2430 # Relative extent 2422 2431 if easting_min is None: 2423 2432 xmin = min(x) … … 2452 2461 nrows = int((ymax-ymin)/cellsize) + 1 2453 2462 2454 # New absolute reference and coordinates2463 # New absolute reference and coordinates 2455 2464 newxllcorner = xmin + xllcorner 2456 2465 newyllcorner = ymin + yllcorner … … 2468 2477 yg = i * cellsize 2469 2478 else: 2470 #this will flip the order of the y values for ers2479 # this will flip the order of the y values for ers 2471 2480 yg = (nrows-i) * cellsize 2472 2481 … … 2478 2487 grid_points[k, 1] = yg 2479 2488 2480 # Interpolate2489 # Interpolate 2481 2490 from anuga.fit_interpolate.interpolate import Interpolate 2482 2491 2483 2492 # Remove loners from vertex_points, volumes here 2484 2493 vertex_points, volumes = remove_lone_verts(vertex_points, volumes) 2485 #export_mesh_file('monkey.tsh',{'vertices':vertex_points, 'triangles':volumes}) 2486 #import sys; sys.exit() 2494 # export_mesh_file('monkey.tsh',{'vertices':vertex_points, 'triangles':volumes}) 2487 2495 interp = Interpolate(vertex_points, volumes, verbose = verbose) 2488 2496 2489 # Interpolate using quantity values2497 # Interpolate using quantity values 2490 2498 if verbose: log.critical('Interpolating') 2491 2499 grid_values = interp.interpolate(result, grid_points).flatten() … … 2495 2503 % (num.min(grid_values), num.max(grid_values))) 2496 2504 2497 # Assign NODATA_value to all points outside bounding polygon (from interpolation mesh)2505 # Assign NODATA_value to all points outside bounding polygon (from interpolation mesh) 2498 2506 P = interp.mesh.get_boundary_polygon() 2499 2507 outside_indices = outside_polygon(grid_points, P, closed=True) … … 2588 2596 2589 2597 ################################################################################ 2590 # Obsolete functions - main atined for backwards compatibility2598 # Obsolete functions - maintained for backwards compatibility 2591 2599 ################################################################################ 2592 2600
Note: See TracChangeset
for help on using the changeset viewer.