source: anuga_work/production/dampier_2006/build_dampier.py @ 3877

Last change on this file since 3877 was 3877, checked in by nick, 16 years ago

build_dampier works... run might needs some changes

File size: 7.5 KB
Line 
1"""Script for running tsunami inundation scenario for Dampier, WA, Australia.
2
3Source data such as elevation and boundary data is assumed to be available in
4directories specified by project.py
5The output sww file is stored in project.output_time_dir
6
7The scenario is defined by a triangular mesh created from project.polygon,
8the elevation data and a simulated submarine landslide.
9
10Ole Nielsen and Duncan Gray, GA - 2005 and Jane Sexton, Nick Bartzis, GA - 2006
11"""
12
13#------------------------------------------------------------------------------
14# Import necessary modules
15#------------------------------------------------------------------------------
16
17# Standard modules
18from os import sep
19from os.path import dirname, basename
20from os import mkdir, access, F_OK
21from shutil import copy
22import time
23import sys
24
25
26# Related major packages
27from anuga.shallow_water import Domain
28from anuga.shallow_water import Dirichlet_boundary
29from anuga.shallow_water import File_boundary
30from anuga.shallow_water import Reflective_boundary
31from anuga.shallow_water.data_manager import convert_dem_from_ascii2netcdf, dem2pts
32from anuga.pmesh.mesh_interface import create_mesh_from_regions
33from anuga.geospatial_data.geospatial_data import *
34from anuga.abstract_2d_finite_volumes.util import Screen_Catcher
35
36# Application specific imports
37import project   # Definition of file names and polygons
38
39#------------------------------------------------------------------------------
40# Copy scripts to time stamped output directory and capture screen
41# output to file
42#------------------------------------------------------------------------------
43
44if access(project.output_build_time_dir,F_OK) == 0:
45    mkdir (project.output_build_time_dir)
46copy (dirname(project.__file__) +sep+ project.__name__+'.py',
47      project.output_build_time_dir + project.__name__+'.py')     #copies project.py
48copy (__file__, project.output_build_time_dir + basename(__file__))
49print 'files '+ project.__name__+'.py and '+ basename(__file__)+' copied to '+ project.output_build_time_dir       #copies this file
50#import sys; sys.exit()
51
52
53#normal screen output is stored in
54screen_output_name = project.output_build_time_dir + "screen_output.txt"
55screen_error_name = project.output_build_time_dir + "screen_error.txt"
56
57#used to catch screen output to file
58sys.stdout = Screen_Catcher(screen_output_name)
59sys.stderr = Screen_Catcher(screen_error_name)
60
61print 'USER:    ', project.user
62
63#-------------------------------------------------------------------------------
64# Preparation of topographic data
65#
66# Convert ASC 2 DEM 2 PTS using source data and store result in source data
67# Do for coarse and fine data
68# Fine pts file to be clipped to area of interest
69#-------------------------------------------------------------------------------
70
71# topography directory filenames
72onshore_dir_name = project.onshore_dir_name
73coast_dir_name = project.coast_dir_name
74islands_dir_name = project.islands_dir_name
75offshore_dir_name = project.offshore_dir_name
76offshore_dir_name1 = project.offshore_dir_name1
77offshore_dir_name2 = project.offshore_dir_name2
78offshore_dir_name3 = project.offshore_dir_name3
79offshore_dir_name4 = project.offshore_dir_name4
80offshore_dir_name5 = project.offshore_dir_name5
81offshore_dir_name6 = project.offshore_dir_name6
82offshore_dir_name7 = project.offshore_dir_name7
83offshore_dir_name8 = project.offshore_dir_name8
84offshore_dir_name9 = project.offshore_dir_name9
85offshore_dir_name10 = project.offshore_dir_name10
86offshore_dir_name11 = project.offshore_dir_name11
87offshore_dir_name12 = project.offshore_dir_name12
88offshore_dir_name13 = project.offshore_dir_name13
89offshore_dir_name14 = project.offshore_dir_name14
90
91# creates DEM from asc data
92convert_dem_from_ascii2netcdf(onshore_dir_name, use_cache=True, verbose=True)
93convert_dem_from_ascii2netcdf(islands_dir_name, use_cache=True, verbose=True)
94
95#creates pts file for onshore DEM
96dem2pts(onshore_dir_name,
97#        easting_min=project.eastingmin,
98#        easting_max=project.eastingmax,
99#        northing_min=project.northingmin,
100#        northing_max= project.northingmax,
101        use_cache=True, 
102        verbose=True)
103
104#creates pts file for islands DEM
105dem2pts(islands_dir_name, use_cache=True, verbose=True)
106
107print'create Geospatial data objects from topographies'
108G1 = Geospatial_data(file_name = project.onshore_dir_name + '.pts')
109G2 = Geospatial_data(file_name = project.coast_dir_name + '.xya')
110G3 = Geospatial_data(file_name = project.islands_dir_name + '.pts')
111G_off = Geospatial_data(file_name = project.offshore_dir_name + '.xya')
112G_off1 = Geospatial_data(file_name = project.offshore_dir_name1 + '.xya')
113G_off2 = Geospatial_data(file_name = project.offshore_dir_name2 + '.xya')
114G_off3 = Geospatial_data(file_name = project.offshore_dir_name3 + '.xya')
115G_off4 = Geospatial_data(file_name = project.offshore_dir_name4 + '.xya')
116G_off5 = Geospatial_data(file_name = project.offshore_dir_name5 + '.xya')
117G_off6 = Geospatial_data(file_name = project.offshore_dir_name6 + '.xya')
118G_off7 = Geospatial_data(file_name = project.offshore_dir_name7 + '.xya')
119G_off8 = Geospatial_data(file_name = project.offshore_dir_name8 + '.xya')
120G_off9 = Geospatial_data(file_name = project.offshore_dir_name9 + '.xya')
121G_off10 = Geospatial_data(file_name = project.offshore_dir_name10 + '.xya')
122G_off11 = Geospatial_data(file_name = project.offshore_dir_name11 + '.xya')
123G_off12 = Geospatial_data(file_name = project.offshore_dir_name12 + '.xya')
124G_off13 = Geospatial_data(file_name = project.offshore_dir_name13 + '.xya')
125G_off14 = Geospatial_data(file_name = project.offshore_dir_name14 + '.xya')
126
127print'add all geospatial objects'
128G = G1 + G2 + G3 + G_off + G_off1 + G_off2 + G_off3 + G_off4 + G_off5 \
129    + G_off6 + G_off7 + G_off8 + G_off9 + G_off10 + G_off11 + G_off12 \
130    + G_off13 + G_off14
131
132print'clip combined geospatial object by bounding polygon'
133G.clip(project.bounding_polygon)
134#FIXME: add a clip function to pts
135
136print'export combined DEM file'
137if access(project.topographies_time_dir,F_OK) == 0:
138    mkdir (project.topographies_time_dir)
139G.export_points_file(project.combined_time_dir_name + '.pts')
140
141#-------------------------------------------------------------------------
142# Convert URS to SWW file for boundary conditions
143#-------------------------------------------------------------------------
144print 'starting to create boundary conditions'
145boundaries_in_dir_name = project.boundaries_in_dir_name
146
147from anuga.shallow_water.data_manager import urs2sww
148
149print 'minlat=project.north_boundary, maxlat=project.south_boundary',project.north_boundary, project.south_boundary
150print 'minlon= project.west_boundary, maxlon=project.east_boundary',project.west_boundary, project.east_boundary
151
152if access(project.boundaries_time_dir,F_OK) == 0:
153    mkdir (project.boundaries_time_dir)
154urs2sww(boundaries_in_dir_name,basename_out= project.boundaries_time_dir_name,
155        minlat=project.south_boundary, maxlat=project.north_boundary,
156        minlon= project.west_boundary, maxlon=project.east_boundary,
157        mint=0, maxt= 35000,
158        verbose='true')
159'''
160from caching import cache
161cache(urs2sww,
162      (boundaries_in_dir_name,
163       project.boundaries_time_dir_name),
164      {'verbose': True,
165       'minlat': project.south_boundary,
166       'maxlat': project.north_boundary,
167       'minlon': project.west_boundary,
168       'maxlon': project.east_boundary,
169       'mint': 0, 'maxt': 35000,
170#       'origin': domain.geo_reference.get_origin(),
171       'mean_stage': project.tide,
172#       'zscale': 1,                 #Enhance tsunami
173       'fail_on_NaN': False},
174       verbose = True,
175       )
176#       dependencies = source_dir + project.boundary_basename + '.sww')
177'''       
178
179
180
181
182
183
184
185
Note: See TracBrowser for help on using the repository browser.