Changeset 5012


Ignore:
Timestamp:
Feb 8, 2008, 5:13:08 PM (16 years ago)
Author:
nick
Message:

update geraldton

Location:
anuga_work/production/geraldton
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • anuga_work/production/geraldton/build_geraldton.py

    r5007 r5012  
    2525# Related major packages
    2626from anuga.shallow_water import Domain
    27 from anuga.shallow_water.data_manager import convert_dem_from_ascii2netcdf, dem2pts
     27from anuga.shallow_water.data_manager import convert_dem_from_ascii2netcdf, dem2pts,start_screen_catcher, copy_code_files
    2828from anuga.geospatial_data.geospatial_data import *
    29 from anuga.abstract_2d_finite_volumes.util import start_screen_catcher, copy_code_files
    3029
    3130# Application specific imports
     
    3736#------------------------------------------------------------------------------
    3837
    39 start_screen_catcher(project.output_build_time_dir)
    4038
    41 print 'time stamp: ',project.gtime
     39print 'time stamp: ',project.build_time
    4240print 'USER: ', project.user
    4341
     
    4644               dirname(project.__file__)+sep+ project.__name__+'.py' )
    4745
     46start_screen_catcher(project.output_build_time_dir)
    4847
    4948
     
    5958onshore_in_dir_name = project.onshore_in_dir_name
    6059coast_in_dir_name = project.coast_in_dir_name
    61 #island_in_dir_name = project.island_in_dir_name
     60island_in_dir_name = project.island_in_dir_name
    6261offshore_in_dir_name = project.offshore_in_dir_name
    6362
    6463onshore_dir_name = project.onshore_dir_name
    6564coast_dir_name = project.coast_dir_name
    66 #island_dir_name = project.island_dir_name
     65island_dir_name = project.island_dir_name
    6766offshore_dir_name = project.offshore_dir_name
    6867
     
    7069print "creates DEMs from ascii data"
    7170convert_dem_from_ascii2netcdf(onshore_in_dir_name, basename_out=onshore_dir_name, use_cache=True, verbose=True)
    72 #convert_dem_from_ascii2netcdf(island_in_dir_name, basename_out=island_dir_name, use_cache=True, verbose=True)
     71convert_dem_from_ascii2netcdf(island_in_dir_name, basename_out=island_dir_name, use_cache=True, verbose=True)
    7372
    7473#creates pts file for onshore DEM
    7574print "creates pts file for onshore DEM"
    7675dem2pts(onshore_dir_name,
    77 #        easting_min=project.eastingmin,
    78 #        easting_max=project.eastingmax,
    79 #        northing_min=project.northingmin,
    80 #        northing_max= project.northingmax,
    8176        use_cache=True,
    8277        verbose=True)
    8378
    8479#creates pts file for island DEM
    85 #dem2pts(island_dir_name, use_cache=True, verbose=True)
     80dem2pts(island_dir_name, use_cache=True, verbose=True)
    8681
    8782print'create Geospatial data1 objects from topographies'
    8883G1 = Geospatial_data(file_name = onshore_dir_name + '.pts')
    89 G2 = Geospatial_data(file_name = coast_in_dir_name + '.xya')
    90 #G3 = Geospatial_data(file_name = island_dir_name + '.pts')
    91 G_off = Geospatial_data(file_name = offshore_in_dir_name + '.xya')
     84print'finished reading in %s.pts' %onshore_dir_name
     85G2 = Geospatial_data(file_name = coast_in_dir_name + '.txt')
     86G3 = Geospatial_data(file_name = island_dir_name + '.pts')
     87G4 = Geospatial_data(file_name = offshore_in_dir_name + '.txt')
     88print'finished reading in files'
     89
     90G_onshore_clip = G1.clip(project.poly_all, verbose=True)
     91print 'finished clip'
     92G_onshore_clip_clipout=G_onshore_clip.clip_outside(project.poly_cbd, verbose=True)
     93print 'finished clipout'
     94# reduce resolution by 9 time (1/9) of the original file. This was determined
     95# as the resolution for this region is going to be 20000 and the asc grid cellsize is
     96# 30m there (900). The sqrt of 20000 is about 140 or 150, a 150m x150m grid would contain
     97# 36 points at a 30m grid spacing and we only need 4. So 4/36 is 1/9 or 0.11111
     98G_onshore_clip_clipout_reduced=G_onshore_clip_clipout.split(0.111, verbose=True)
     99
     100G_cbd = G1.clip(projec.poly_cbd,verbose=True)
     101G_cbd_reduced = G_cbd.split(0.25, verbose=True)
     102
     103G_bathy = G4.clip(project.poly_all)
     104G_island = G3.clip(project.poly_all)
     105G_coast = G2.clip(project.poly_all)
    92106
    93107print'add all geospatial objects'
    94 G = G1 + G2 + G_off
     108G = G_onshore_clip_clipout_reduced + G_cbd_reduced + G_bathy + G_island + G_coast
    95109
    96110print'clip combined geospatial object by bounding polygon'
     
    102116if access(project.topographies_dir,F_OK) == 0:
    103117    mkdir (project.topographies_dir)
    104 G.export_points_file(project.combined_dir_name + '.xya')
     118G.export_points_file(project.combined_dir_name + '.txt')
    105119#G_clipped.export_points_file(project.combined_dir_name + '.xya')
     120G_small, = G.split(0.1, verbose)
     121G_small.export_points_file(project.combined_dir_name_small + '.txt')
    106122
    107123'''
  • anuga_work/production/geraldton/project.py

    r5006 r5012  
    3131#Making assumptions about the location of scenario data
    3232state = 'western_australia'
    33 scenario_name = 'gerladton'
    34 scenario = 'busselton_tsunami_scenario_2006'
     33scenario_name = 'geraldton'
     34scenario = 'geraldton_tsunami_scenario'
    3535
    3636#Maybe will try to make project a class to allow these parameters to be passed in.
     
    4141finaltime=25000
    4242export_cellsize=50
    43 setup='final'
     43setup='trial'
    4444source='other'
    4545
     
    6565
    6666# onshore data provided by WA DLI
    67 onshore_name = 'DLI_orthophoto_DEM' # original
    68 #islands
     67onshore_name = 'landgate_dem_clip' # original
     68island_name = 'dted_srtm_islands'
    6969
    7070# AHO + DPI data
    71 coast_name = '100k_coastline'
    72 offshore_name = 'Bathymetry'
     71coast_name = 'coastline'
     72offshore_name = 'geraldton_bathy'
    7373
    7474#final topo name
    75 combined_name ='busselton_combined_elevation.pts'
    76 combined_smaller_name = 'busselton_combined_elevation_smaller'
     75combined_name ='geraldton_combined_elevation.pts'
     76combined_smaller_name = 'geraldton_combined_elevation_smaller'
    7777
    7878anuga_dir = home+state+sep+scenario+sep+'anuga'+sep
     
    8383#input topo file location
    8484onshore_in_dir_name = topographies_in_dir + onshore_name
    85 #island_in_dir_name = topographies_in_dir + island_name
     85island_in_dir_name = topographies_in_dir + island_name
    8686
    8787coast_in_dir_name = topographies_in_dir + coast_name
     
    8989
    9090onshore_dir_name = topographies_dir + onshore_name
    91 #island_dir_name = topographies_dir + island_name
     91island_dir_name = topographies_dir + island_name
    9292coast_dir_name = topographies_dir + coast_name
    9393offshore_dir_name = topographies_dir + offshore_name
     
    128128output_dir = anuga_dir+'outputs'+sep
    129129output_build_time_dir = output_dir+build_time+'_'+source+sep
    130 #output_run_time_dir = output_dir +run_time+dir_comment+sep
    131 output_run_time_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'outputs'+sep+run_time+sep
     130output_run_time_dir = output_dir +run_time+dir_comment+sep
     131#output_run_time_dir = output_dir+sep+run_time+sep
    132132output_run_time_dir_name = output_run_time_dir + scenario_name  #Used by post processing
    133133
     
    150150from anuga.utilities.polygon import read_polygon, plot_polygons, polygon_area, is_inside_polygon
    151151
    152 # bounding polygon for study area
    153 poly_all = read_polygon(polygons_dir+'poly_all.csv')
    154 res_poly_all = 100000*res_factor
    155152
    156153###################################################################
     
    169166
    170167#digitized polygons
    171 poly_busselton1 = read_polygon(polygons_dir+'neg20_pos10_polygon.csv')
    172 res_busselton1 = 10000*res_factor
    173 poly_busselton2 = read_polygon(polygons_dir+'neg5_pos5_poly_.csv')
    174 res_busselton2 = 500*res_factor
     168# bounding polygon for study area
     169poly_all = read_polygon(polygons_dir+'extent_new_points.csv')
     170res_poly_all = 100000*res_factor
     171
     172poly_neg20_pos20 = read_polygon(polygons_dir+'neg20_pos20_points.csv')
     173res_neg20_pos20 = 25000*res_factor
     174
     175poly_neg10_pos10= read_polygon(polygons_dir+'neg10_pos10_points.csv')
     176res_neg10_pos10 = 5000*res_factor
     177
     178poly_cbd = read_polygon(polygons_dir+'cbd_points.csv')
     179res_cbd = 500*res_factor
     180
     181poly_wallabi = read_polygon(polygons_dir+'island_wallabi_poly_points.csv')
     182res_wallabi = 25000*res_factor
     183poly_dingiville = read_polygon(polygons_dir+'island_dingiville_poly_points.csv')
     184res_dingiville = 25000*res_factor
     185poly_pelsaert = read_polygon(polygons_dir+'island_pelsaert_poly_points.csv')
     186res_pelsaert = 25000*res_factor
     187
    175188
    176189#plot_polygons([polyAll,poly_broome1,poly_broome2,poly_broome3],figname='boundingpoly2',verbose=False)
    177190
    178 interior_regions = [[poly_busselton1,res_busselton1],[poly_busselton2,res_busselton2]]
    179 
    180 boundary_tags={'back': [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 14],
    181                'side': [10,13], 'ocean': [11, 12]}
     191interior_regions = [[poly_neg20_pos20,res_neg20_pos20],[poly_neg10_pos10,res_neg10_pos10],
     192                    [poly_cbd,res_cbd],[poly_wallabi,res_wallabi],
     193                    [poly_dingiville,res_dingiville],[poly_pelsaert, res_pelsaert]]
     194
     195boundary_tags={'back': [2, 3, 4, 5],
     196               'side': [1, 6], 'ocean': [0, 7, 8, 9, 10]}
    182197
    183198trigs_min = number_mesh_triangles(interior_regions, poly_all, res_poly_all)
    184199
    185 poly_mainland=read_polygon(polygons_dir+'initial_condition_south.csv')
     200poly_mainland=read_polygon(polygons_dir+'initial_condition.csv')
    186201
    187202print 'min number triangles', trigs_min
Note: See TracChangeset for help on using the changeset viewer.