source: anuga_work/production/sydney_2006/process_gauges.py @ 5001

Last change on this file since 5001 was 3514, checked in by duncan, 18 years ago

Hi all,
I'm doing a change in the anuga structure, moving the code to

\anuga_core\source\anuga

After you have done an svn update, the PYTHONPATH has to be changed to;
PYTHONPATH = anuga_core/source/

This is part of changes required to make installation of anuga quicker and reducing the size of our sandpits.

If any imports are broken, try fixing them. With adding anuga. to them for example. If this seems to have really broken things, email/phone me.

Cheers
Duncan

File size: 2.5 KB
Line 
1"""Read in sww file, interpolate at specified locations and plot time series
2Additionally, store augmented building database with max inundation depths
3"""
4
5import project
6from anuga.pyvolution.util import file_function
7#from anuga.pyvolution.coordinate_transforms.redfearn import degminsec2decimal_degrees, redfearn
8from pylab import *
9
10
11plot = False
12   
13swwfile = project.outputname4 + '.sww'
14
15def get_gauges_from_file(filename):
16    fid = open(filename)
17    lines = fid.readlines()
18    fid.close()
19
20    gauges = []
21    gaugelocation = []
22    for line in lines[1:]:
23        fields = line.split(',')
24        location = fields[0]
25        easting = float(fields[1])
26        northing = float(fields[2])
27     
28        gauges.append([easting, northing])
29        gaugelocation.append(location)
30
31    #Return gauges and raw data for subsequent storage   
32    return gauges, lines, gaugelocation
33
34gauges, lines, locations = get_gauges_from_file(project.gauge_filename)
35
36#Read model output
37quantities = ['stage', 'elevation', 'xmomentum', 'ymomentum']
38f = file_function(swwfile,
39                  quantities = quantities,
40                  interpolation_points = gauges,
41                  verbose = True,
42                  use_cache = True)
43
44
45lines[0] = lines[0].strip() +\
46               ',MAX INUNDATION DEPTH (m)' +\
47               ',MAX MOMENTUM (m^2/s)\n'           
48               
49from math import sqrt
50N = len(gauges)
51for k, g in enumerate(gauges):
52    if k%((N+10)/10)==0:
53        print 'Doing row %d of %d' %(k, N)
54
55    model_time = []
56    stages = []
57    elevations = []
58    momenta = []
59
60    max_depth = 0
61    max_momentum = 0   
62    for i, t in enumerate(f.get_time()):
63        w = f(t, point_id = k)[0]
64        z = f(t, point_id = k)[1]
65        uh = f(t, point_id = k)[2]
66        vh = f(t, point_id = k)[3]       
67
68        m = sqrt(uh*uh + vh*vh)   #Absolute momentum
69       
70        model_time.append(t)       
71        stages.append(w)
72        elevations.append(z)  #Should be constant over time
73        momenta.append(m) 
74
75        if w-z > max_depth:
76            max_depth = w-z
77        if m > max_momentum:           
78            max_momentum = m   
79
80    #Augment building data
81    lines[k+1] = lines[k+1].strip() +\
82                     ',%f' %max_depth +\
83                     ',%f\n' %max_momentum
84   
85    print lines[k]
86   
87#Store new building file with mak depths added
88
89#FN = 'inundation_augmented_' + project.gauge_filename
90FN = project.gauge_outname
91fid = open(FN, 'w')
92for line in lines:
93    fid.write(line)
94fid.close()   
95
Note: See TracBrowser for help on using the repository browser.