source: anuga_work/production/busselton/build_busselton.py @ 5525

Last change on this file since 5525 was 5525, checked in by kristy, 15 years ago

updated data input

File size: 5.1 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# Related major packages
26from anuga.shallow_water import Domain
27from anuga.shallow_water.data_manager import convert_dem_from_ascii2netcdf, dem2pts
28from anuga.geospatial_data.geospatial_data import *
29from anuga.shallow_water.data_manager import start_screen_catcher, copy_code_files,store_parameters
30
31# Application specific imports
32import project   # Definition of file names and polygons
33
34#------------------------------------------------------------------------------
35# Copy scripts to time stamped output directory and capture screen
36# output to file
37#------------------------------------------------------------------------------
38
39copy_code_files(project.output_build_time_dir,__file__, 
40               dirname(project.__file__)+sep+ project.__name__+'.py' )
41
42start_screen_catcher(project.output_build_time_dir)
43
44print 'USER: ', project.user
45
46#-------------------------------------------------------------------------------
47# Preparation of topographic data
48#
49# Convert ASC 2 DEM 2 PTS using source data and store result in source data
50# Do for coarse and fine data
51# Fine pts file to be clipped to area of interest
52#-------------------------------------------------------------------------------
53print"project.poly_all",project.poly_all
54print"project.combined_dir_name",project.combined_dir_name
55
56# topography directory filenames
57onshore_in_dir_name = project.onshore_in_dir_name
58coast_in_dir_name = project.coast_in_dir_name
59coast_in_dir_name1 = project.coast_in_dir_name1
60offshore_in_dir_name = project.offshore_in_dir_name
61offshore_in_dir_name1 = project.offshore_in_dir_name1
62offshore_in_dir_name2 = project.offshore_in_dir_name2
63offshore_in_dir_name3 = project.offshore_in_dir_name3
64
65onshore_dir_name = project.onshore_dir_name
66coast_dir_name = project.coast_dir_name
67coast_dir_name1 = project.coast_dir_name1
68offshore_dir_name = project.offshore_dir_name
69offshore_dir_name1 = project.offshore_dir_name1
70offshore_dir_name2 = project.offshore_dir_name2
71offshore_dir_name3 = project.offshore_dir_name3
72
73# creates DEM from asc data
74print "creates DEMs from ascii data"
75convert_dem_from_ascii2netcdf(onshore_in_dir_name, basename_out=onshore_dir_name, use_cache=True, verbose=True)
76
77#creates pts file for onshore DEM
78print "creates pts file for onshore DEM"
79dem2pts(onshore_dir_name,
80#        easting_min=project.eastingmin,
81#        easting_max=project.eastingmax,
82#        northing_min=project.northingmin,
83#        northing_max= project.northingmax,
84        use_cache=True, 
85        verbose=True)
86
87#creates pts file for island DEM
88#dem2pts(island_dir_name, use_cache=True, verbose=True)
89
90print'create Geospatial data1 objects from topographies'
91G1 = Geospatial_data(file_name = onshore_dir_name + '.pts',verbose=True)
92G2 = Geospatial_data(file_name = coast_in_dir_name + '.txt',verbose=True)
93G3 = Geospatial_data(file_name = coast_in_dir_name1 + '.txt',verbose=True)
94G_off = Geospatial_data(file_name = offshore_in_dir_name + '.txt',verbose=True)
95G_off1 = Geospatial_data(file_name = offshore_in_dir_name1 + '.txt',verbose=True)
96G_off2 = Geospatial_data(file_name = offshore_in_dir_name2 + '.txt',verbose=True)
97G_off3 = Geospatial_data(file_name = offshore_in_dir_name3 + '.txt',verbose=True)
98                         
99print'add all geospatial objects'
100G = G1 + G2 + G3 + G_off + G_off1 + G_off2 + G_off3
101
102print'clip combined geospatial object by bounding polygon'
103G_clipped = G.clip(project.poly_all)
104
105print'export combined DEM file'
106if access(project.topographies_dir,F_OK) == 0:
107    mkdir (project.topographies_dir)
108G_clipped.export_points_file(project.combined_dir_name + '.txt')
109
110print'project.combined_dir_name + .txt',project.combined_dir_name + '.txt'
111
112#-------------------------------------------------------------------------
113# Convert URS to SWW file for boundary conditions
114#-------------------------------------------------------------------------
115print 'starting to create boundary conditions'
116from anuga.shallow_water.data_manager import urs2sww, urs_ungridded2sww
117
118boundaries_in_dir_name = project.boundaries_in_dir_name
119print 'boundaries_in_dir_name',project.boundaries_in_dir_name
120
121
122#import sys; sys.exit()
123
124urs_ungridded2sww(project.boundaries_in_dir_name, project.boundaries_in_dir_name,
125                  verbose=True, mint=4000, maxt=80000, zscale=1)
126
Note: See TracBrowser for help on using the repository browser.