""" Script for building the elevation data to run a tsunami inundation scenario for dampier_ph2, WA, Australia. Input: elevation data from project_250m.py Output: pts file stored in project_250m.topographies_dir The run_dampier.py is reliant on the output of this script. """ #------------------------------------------------------------------------------ # Import necessary modules #------------------------------------------------------------------------------ # Standard modules from os import sep from os.path import dirname, basename from os import mkdir, access, F_OK from shutil import copy import time import sys # Related major packages from anuga.shallow_water import Domain from anuga.shallow_water.data_manager import convert_dem_from_ascii2netcdf, dem2pts from anuga.geospatial_data.geospatial_data import * from anuga.shallow_water.data_manager import start_screen_catcher, copy_code_files # Application specific imports import project_250m # Definition of file names and polygons #------------------------------------------------------------------------------ # Copy scripts to time stamped output directory and capture screen # output to file #------------------------------------------------------------------------------ copy_code_files(project_250m.output_build_time_dir,__file__, dirname(project_250m.__file__)+sep+ project_250m.__name__+'.py' ) start_screen_catcher(project_250m.output_build_time_dir) print 'USER: ', project_250m.user #------------------------------------------------------------------------------- # Preparation of topographic data # # Convert ASC 2 DEM 2 PTS using source data and store result in source data # Do for coarse and fine data # Fine pts file to be clipped to area of interest #------------------------------------------------------------------------------- print "project_250m.poly_all", project_250m.poly_all print "project_250m.combined_dir_name", project_250m.combined_dir_name # input elevation directory filenames onshore_in_dir_name = project_250m.onshore_in_dir_name # output elevation directory filenames onshore_dir_name = project_250m.onshore_dir_name # creates DEM from asc data print "creates DEMs from asc data" convert_dem_from_ascii2netcdf(onshore_in_dir_name, basename_out=onshore_dir_name, use_cache=True, verbose=True) # creates pts file for onshore DEM print "creates pts file for onshore DEM" dem2pts(onshore_dir_name ,use_cache=True,verbose=True) # create onshore pts files print'create Geospatial data1 objects from topographies' G = Geospatial_data(file_name = onshore_dir_name + '.pts') #------------------------------------------------------------------------------- # Combine, clip and export dataset #------------------------------------------------------------------------------- print'clip combined geospatial object by bounding polygon' G_clipped = G.clip(project_250m.poly_all) print'export combined DEM file' if access(project_250m.topographies_dir,F_OK) == 0: mkdir (project_250m.topographies_dir) G_clipped.export_points_file(project_250m.combined_dir_name + '.pts') G_clipped.export_points_file(project_250m.combined_dir_name + '.txt') #Use for comparision in ARC