source: anuga_work/production/Broome_2008/build_broome.py @ 7159

Last change on this file since 7159 was 5813, checked in by kristy, 15 years ago
File size: 6.5 KB
Line 
1"""
2Script for building the elevation data to run tsunami inundation scenario
3for Broome, WA, Australia.
4
5Input: elevation data from project.py
6Output: pts file stored in project.topographies_dir
7The run_broome.py is reliant on the output of this script.
8
9Ole Nielsen and Duncan Gray, GA - 2005 and Jane Sexton, Nick Bartzis, GA - 2006
10"""
11
12#------------------------------------------------------------------------------
13# Import necessary modules
14#------------------------------------------------------------------------------
15
16# Standard modules
17from os import sep
18from os.path import dirname, basename
19from os import mkdir, access, F_OK
20from shutil import copy
21import time
22import sys
23
24# Related major packages
25from anuga.shallow_water import Domain
26from anuga.shallow_water.data_manager import convert_dem_from_ascii2netcdf, dem2pts
27from anuga.geospatial_data.geospatial_data import *
28from anuga.shallow_water.data_manager import start_screen_catcher, copy_code_files
29
30# Application specific imports
31import project   # Definition of file names and polygons
32
33#------------------------------------------------------------------------------
34# Copy scripts to time stamped output directory and capture screen
35# output to file
36#------------------------------------------------------------------------------
37
38copy_code_files(project.output_build_time_dir,__file__, 
39               dirname(project.__file__)+sep+ project.__name__+'.py' )
40
41start_screen_catcher(project.output_build_time_dir)
42
43print 'USER:    ', project.user
44
45#-------------------------------------------------------------------------------
46# Preparation of topographic data
47#
48# Convert ASC 2 DEM 2 PTS using source data and store result in source data
49# Do for coarse and fine data
50# Fine pts file to be clipped to area of interest
51#-------------------------------------------------------------------------------
52print"project.poly_all",project.poly_all
53print"project.combined_dir_name",project.combined_dir_name
54
55# input elevation directory filenames
56onshore_in_dir_name = project.onshore_in_dir_name
57
58coast_in_dir_name = project.coast_in_dir_name
59
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
64offshore_in_dir_name4 = project.offshore_in_dir_name4
65offshore_in_dir_name5 = project.offshore_in_dir_name5
66offshore_in_dir_name6 = project.offshore_in_dir_name6
67
68# output elevation directory filenames
69onshore_dir_name = project.onshore_dir_name
70
71coast_dir_name = project.coast_dir_name
72
73offshore_dir_name = project.offshore_dir_name
74offshore_dir_name1 = project.offshore_dir_name1
75offshore_dir_name2 = project.offshore_dir_name2
76offshore_dir_name3 = project.offshore_dir_name3
77offshore_dir_name4 = project.offshore_dir_name4
78offshore_dir_name5 = project.offshore_dir_name5
79offshore_dir_name6 = project.offshore_dir_name6
80
81# creates DEM from asc data
82print "creates DEMs from asc data"
83convert_dem_from_ascii2netcdf(onshore_in_dir_name, basename_out=onshore_dir_name, use_cache=True, verbose=True)
84
85convert_dem_from_ascii2netcdf(offshore_in_dir_name1, basename_out=offshore_dir_name1, use_cache=True, verbose=True)
86convert_dem_from_ascii2netcdf(offshore_in_dir_name2, basename_out=offshore_dir_name2, use_cache=True, verbose=True)
87convert_dem_from_ascii2netcdf(offshore_in_dir_name3, basename_out=offshore_dir_name3, use_cache=True, verbose=True)
88convert_dem_from_ascii2netcdf(offshore_in_dir_name4, basename_out=offshore_dir_name4, use_cache=True, verbose=True)
89convert_dem_from_ascii2netcdf(offshore_in_dir_name5, basename_out=offshore_dir_name5, use_cache=True, verbose=True)
90convert_dem_from_ascii2netcdf(offshore_in_dir_name6, basename_out=offshore_dir_name6, use_cache=True, verbose=True)
91
92#creates pts file for onshore DEM
93print "creates pts file for onshore DEM"
94dem2pts(onshore_dir_name ,use_cache=True,verbose=True)
95dem2pts(offshore_dir_name1 ,use_cache=True,verbose=True)
96dem2pts(offshore_dir_name2 ,use_cache=True,verbose=True)
97dem2pts(offshore_dir_name3 ,use_cache=True,verbose=True)
98dem2pts(offshore_dir_name4 ,use_cache=True,verbose=True)
99dem2pts(offshore_dir_name5 ,use_cache=True,verbose=True)
100dem2pts(offshore_dir_name6 ,use_cache=True,verbose=True)
101
102###creates pts file for island DEM
103##dem2pts(island_dir_name, use_cache=True, verbose=True)
104##dem2pts(island_dir_name1, use_cache=True, verbose=True)
105##dem2pts(island_dir_name2, use_cache=True, verbose=True)
106##dem2pts(island_dir_name3, use_cache=True, verbose=True)
107##dem2pts(island_dir_name4, use_cache=True, verbose=True)
108##dem2pts(island_dir_name5, use_cache=True, verbose=True)
109##dem2pts(island_dir_name6, use_cache=True, verbose=True)
110
111#-------------------------------------------------------------------------------
112# Combine datasets into project.combined_dir_name
113#-------------------------------------------------------------------------------
114
115print'create Geospatial data1 objects from topographies'
116G1 = Geospatial_data(file_name = onshore_dir_name + '.pts')
117print'create Geospatial data6 objects from topographies'
118G_coast = Geospatial_data(file_name = coast_in_dir_name)
119print'create Geospatial data7 objects from topographies'
120G_off = Geospatial_data(file_name = offshore_in_dir_name)
121print'create Geospatial data8 objects from topographies'
122G_off1 = Geospatial_data(file_name = offshore_dir_name1 + '.pts')
123print'create Geospatial data9 objects from topographies'
124G_off2 = Geospatial_data(file_name = offshore_dir_name2 + '.pts')
125print'create Geospatial data10 objects from topographies'
126G_off3 = Geospatial_data(file_name = offshore_dir_name3 + '.pts')
127print'create Geospatial data10 objects from topographies'
128G_off4 = Geospatial_data(file_name = offshore_dir_name4 + '.pts')
129print'create Geospatial data10 objects from topographies'
130G_off5 = Geospatial_data(file_name = offshore_dir_name5 + '.pts')
131print'create Geospatial data10 objects from topographies'
132G_off6 = Geospatial_data(file_name = offshore_dir_name6 + '.pts')
133print'add all geospatial objects'
134G = G1 + G_coast + G_off + G_off1 + G_off2 + G_off3  + G_off4 + G_off5 + G_off6
135print'clip combined geospatial object by bounding polygon'
136G_clipped = G.clip(project.poly_all)
137
138print'export combined DEM file'
139if access(project.topographies_dir,F_OK) == 0:
140    mkdir (project.topographies_dir)
141G_clipped.export_points_file(project.combined_dir_name + '.pts')
142#G_clipped.export_points_file(project.combined_dir_name + '.txt') #Use for comparision in ARC
143
144import sys
145sys.exit()
Note: See TracBrowser for help on using the repository browser.