source: anuga_validation/automated_validation_tests/urs_mux_files_validation/run_Bf.py @ 4403

Last change on this file since 4403 was 4403, checked in by nick, 17 years ago

added Boundary_file test, this compares 2 time series from one location from 2 different sww files. One sww boundary file is created from MUX files and the other is from a simple evolve model which used the MUX sww file as boundary conditions.

This test was created to check URS data

NOTE this is still being developed and is not working yet

File size: 9.8 KB
Line 
1"""
2
3this test compares the gauge data from a boundary sww file which is precomputed and
4stored in this directory, and the same gauge location from a simple model that uses
5the same boundary sww file for its boundary conditions
6
7Basically tests that file_boundary and evolve is working correctly
8"""
9
10#------------------------------------------------------------------------------
11# Import necessary modules
12#------------------------------------------------------------------------------
13
14# Standard modules
15from os import sep,getcwd, access, F_OK, mkdir, getenv
16from os.path import dirname, basename,abspath
17from shutil import copy
18import time, sys, os, tempfile
19#import sys
20#import os
21#import tempfile
22
23# Related major packages
24from anuga.shallow_water import Domain,Dirichlet_boundary,File_boundary,Transmissive_boundary, Field_boundary
25#from anuga.shallow_water.shallow_water_domain import Field_boundary
26from Numeric import allclose, array
27from anuga.pmesh.mesh_interface import create_mesh_from_regions
28from anuga.abstract_2d_finite_volumes.util import start_screen_catcher, copy_code_files, sww2timeseries, get_data_from_file
29#from anuga.caching import myhash
30# Application specific imports
31
32#------------------------------------------------------------------------------
33# Copy scripts to time stamped output directory and capture screen
34# output to file
35#------------------------------------------------------------------------------
36out_dir = dirname(abspath(__file__))+sep
37#print'temp_name'
38fileName="temp"
39#import sys
40#sys.exit()
41meshes_dir_name = 'small.tsh' # this will be local
42
43tide = 2.4
44
45
46#EVOLVED MODEL
47
48#--------------------------------------------------------------------------
49# Create the triangular mesh based on overall clipping polygon with a
50# tagged
51# boundary and interior regions defined in project.py along with
52# resolutions (maximal area of per triangle) for each polygon
53#--------------------------------------------------------------------------
54
55#all = [[469000,7760000],[470000,7758000],[468000,7758000]]
56#all = [[470000,7760000],[469000,7758000],[468000,7760000]]
57
58#all_large = [[474000,7764000],[474000,7756000],[466000,7756000],[466000,7764000]]
59#all = [[470000,7764000],[470000,7756000],[466000,7756000],[466000,7764000]]
60all = [[465184,7764500],[470397,7764510],[470407,7758988],[465195,7758979]]
61# N, W, S, E
62
63create_mesh_from_regions(all,
64#                             boundary_tags={'ocean': [0,1,2]},
65#                             boundary_tags={'ocean': [ 2],'side': [0, 1]},
66                             boundary_tags={'ocean': [ 0],'side': [1, 3], 'back': [2]},
67                             maximum_triangle_area=50000,
68                             filename=meshes_dir_name,
69                             use_cache=False,
70#                             verbose=True
71                             )
72
73#-------------------------------------------------------------------------
74# Setup computational domain
75#-------------------------------------------------------------------------
76#print 'Setup computational domain'
77domain = Domain( meshes_dir_name, verbose=True)
78
79from anuga.shallow_water.data_manager import urs2sww
80boundaries_dir_name = 'o_test'
81
82# convert MUX urs files to an SWW file output
83#print 'boundary file is: ',boundaries_dir_name
84'''
85from caching import cache
86cache(urs2sww,
87          (boundaries_dir_name,
88           boundaries_dir_name),
89          {'verbose': False,
90           'mint': 9200, 'maxt': 11200,
91#           'origin': domain.geo_reference.get_origin(),
92           'fail_on_NaN': False},
93           verbose = True,
94           )
95'''
96urs2sww(boundaries_dir_name,boundaries_dir_name, verbose=True,
97           mint=9200, maxt= 11200,
98#           'origin': domain.geo_reference.get_origin(),
99           fail_on_NaN= False)
100   
101#-------------------------------------------------------------------------
102# Setup initial conditions
103#-------------------------------------------------------------------------
104
105#print 'Start Set quantity'
106
107domain.set_quantity('elevation', -42.3)
108
109#print 'Finished Set quantity'
110
111#------------------------------------------------------
112# Set domain parameters
113#------------------------------------------------------
114
115domain.set_quantity('stage', tide)
116domain.set_quantity('friction', 0.01)
117domain.set_name(fileName)
118domain.set_datadir(out_dir)
119
120#-------------------------------------------------------------------------
121# Setup boundary conditions
122#-------------------------------------------------------------------------
123#Bf = File_boundary(out_dir + 'o_test_8500_12000.sww',
124#                  domain, time_thinning=24, use_cache=True, verbose=True)
125Bf = Field_boundary(out_dir + 'o_test.sww',
126                  domain, mean_stage=tide, use_cache=True, verbose=False)
127
128#print 'finished reading boundary file'
129
130#Br = Reflective_boundary(domain)
131Bt = Transmissive_boundary(domain)
132Bd = Dirichlet_boundary([tide,0,0])
133
134#print'set_boundary'
135
136domain.set_boundary({'back': Bt,
137                     'side': Bd,
138                    'ocean': Bf
139                    })
140#print'finish set boundary'
141
142#----------------------------------------------------------------------------
143# Evolve system through time
144#----------------------------------------------------------------------------
145
146#remove mesh file!!! and o_test.sww
147t0 = time.time()
148
149for t in domain.evolve(yieldstep = 60, finaltime = 1920):
150    domain.write_time()
151    domain.write_boundary_statistics()
152
153
154
155#Gets timeseries from boundary sww and evolved sww
156home = getenv('INUNDATIONHOME') #Sandpit's parent dir   
157#user = get_user_name()
158data = 'data'
159state = 'western_australia'
160scenario_name = 'dampier.sww'
161
162scenario = 'dampier_tsunami_scenario_2006'
163#scenario = 'test_dampier'
164an = 'anuga'
165bo = 'boundaries'
166
167run_time = 'blank'
168#run_time = project.run_time
169production_dirs = {run_time: 'URS evolved data'#,
170                   #'boundaries': 'URS boundary condition'
171                   }
172
173topo = 'topographies'
174out = 'outputs'
175urs = 'urs'
176gridded = '1_10000'
177
178
179#gauge_boundary_filename = 'boundary_gauge_near_top.csv'
180gauge_boundary_filename = 'gauges_time_series_b_near_top.csv'
181gauge_evolved_filename = 'gauges_time_series_near_top.csv'
182
183boundary_dir_filename = os.path.join(out_dir,gauge_boundary_filename)
184#print'boundary_dir_filename',boundary_dir_filename
185
186evolved_dir_filename= os.path.join(out_dir,gauge_evolved_filename)
187
188#print'boundary_dir_filename',boundary_dir_filename
189#print'evolved_dir_filename',evolved_dir_filename
190                                     
191swwfiles = {}
192swwfile = out_dir + fileName + '.sww'
193swwfiles[swwfile] = run_time
194#print"swwfiles",swwfiles,"shallow_water"
195       
196texname, elev_output = sww2timeseries(swwfiles,
197                                      out_dir+sep+"gauges.csv",
198                                      production_dirs,
199                                      report = False,
200                                      plot_quantity = ['stage', 'xmomentum', 'ymomentum'],
201                                      surface = False,
202                                      time_min = None,
203                                      time_max = None,
204                                      title_on = False,
205                                      use_cache = True,
206                                      verbose = False)
207                                     
208swwfiles = {}
209swwfile = out_dir + boundaries_dir_name + '.sww'
210swwfiles[swwfile] = run_time
211#print"swwfiles",swwfiles,"shallow_water"
212       
213texname, elev_output = sww2timeseries(swwfiles,
214                                      out_dir+sep+"boundary_gauges.csv",
215                                      production_dirs,
216                                      report = False,
217                                      plot_quantity = ['stage', 'xmomentum', 'ymomentum'],
218                                      surface = False,
219                                      time_min = None,
220                                      time_max = None,
221                                      title_on = False,
222                                      use_cache = True,
223                                      verbose = False)
224
225#makes the csv files from the evolved model
226e_header, e_data = get_data_from_file(evolved_dir_filename)
227
228e_time = e_data[:,0]
229e_stage = e_data[:,1]
230e_momentum = e_data[:,2]
231e_speed = e_data[:,3]
232e_elevation = e_data[:,4]
233
234#print'boundary_dir_filename',boundary_dir_filename
235b_header, b_data = get_data_from_file(boundary_dir_filename)
236
237b_time = b_data[:,0]
238b_stage = b_data[:,1]
239b_momentum = b_data[:,2]
240b_speed = b_data[:,3]
241b_elevation = b_data[:,4]
242
243
244# compares the 2 models
245j=0
246k=0
247b_sample = []
248e_sample = []
249for i in range(len(b_time)):
250#    if j<(len(e_time)) and b_time[i] == e_time[j]:
251    if j<(len(e_time)) and k<(len(e_time)) and b_time[i] == e_time[j]:
252        b_sample.append(float(tide+b_stage[i]))
253        e_sample.append(float(e_stage[k]))
254#        if k <len(e_time): print 'time e equal b:', b_time[i],i, j, b_stage[i], b_sample[i], e_stage[k],(len(e_time)-1)
255        j = j +1
256        k = k +1
257
258
259#print len(b_sample), len(e_stage),e_stage,type(e_stage)
260e_stage = e_stage.tolist()
261e_stage.pop()
262#e_stage.pop()
263 
264#print len(b_sample), len(e_stage),e_stage
265
266
267#assert allclose (b_sample, e_stage, 0.2, 0.2)
268print "test successful" 
269print 'fileName',fileName+'.sww'
270os.remove(fileName+'.sww')
271#print 'evolved_dir_filename',evolved_dir_filename
272os.remove('gauges_time_series_near_top.csv')
273os.remove('gauges_time_series_b_near_top.csv')
274os.remove('gauges_t0.csv')
275os.remove('gauges_maxmins.csv')
276os.remove(meshes_dir_name)
277
278
279assert allclose (b_sample, e_sample, 0.5, 0.5)
280
281
282
283
284
285
286
287                             
288                             
289                                         
290
291
292
293
294
Note: See TracBrowser for help on using the repository browser.