Changeset 5139


Ignore:
Timestamp:
Mar 7, 2008, 4:37:36 PM (17 years ago)
Author:
herve
Message:

deleting data files

Location:
anuga_work/production/tonga
Files:
18 deleted
2 edited

Legend:

Unmodified
Added
Removed
  • anuga_work/production/tonga/project.py

    r5135 r5139  
    22"""Common filenames and locations for topographic data, meshes and outputs.
    33"""
    4 print 'hello'
     4
    55from os import sep, environ, getenv, getcwd
    66from os.path import expanduser
    77import sys
    88from time import localtime, strftime, gmtime
    9 from anuga.utilities.polygon import read_polygon, plot_polygons, \
    10                                     polygon_area, is_inside_polygon
     9from anuga.utilities.polygon import read_polygon, plot_polygons, is_inside_polygon, number_mesh_triangles
     10#from anuga.coordinate_transforms.redfearn import degminsec2decimal_degrees, convert_points_from_latlon_to_utm
     11from anuga.utilities.system_tools import get_user_name, get_host_name
     12
     13# file and system info
     14#---------------------------------
     15#codename = 'project.py'
     16
     17home = getenv('INUNDATIONHOME') + sep +'data'+sep #Sandpit's parent diruser = get_user_name()
     18user = get_user_name()
     19host = get_host_name()
     20
     21# INUNDATIONHOME is the inundation directory, not the data directory.
     22
     23#time stuff
     24time = strftime('%Y%m%d_%H%M%S',localtime()) #gets time for new dir
     25gtime = strftime('%Y%m%d_%H%M%S',gmtime()) #gets time for new dir
     26build_time = time+'_build'
     27run_time = time+'_run'
     28print 'gtime: ', gtime
     29
     30#Making assumptions about the location of scenario data
     31state = 'sw_pacific'
     32scenario_name = 'Tonga'
     33scenario = 'fixed_wave'
     34
     35tide = 0
     36
     37alpha = 0.1
     38friction=0.01
     39starttime=10000
     40midtime=21600
     41finaltime=10000
     42export_cellsize=50
     43setup='trial'
     44
     45if setup =='trial':
     46    print'trial'
     47    res_factor=10
     48    time_thinning=48
     49    yieldstep=240
     50if setup =='basic':
     51    print'basic'
     52    res_factor=4
     53    time_thinning=12
     54    yieldstep=120
     55if setup =='final':
     56    print'final'
     57    res_factor=1
     58    time_thinning=4
     59    yieldstep=60
     60
     61dir_comment='_'+setup+'_'+str(tide)+'_'+str(source)+'_'+str(user)
     62
     63
     64topo_gridfile ='tongatapu_10mgrid'
     65
     66anuga_dir = home+state+sep+scenario+sep+'anuga'+sep
     67
     68topographies_in_dir = home+sep+state+sep+scenario+sep+'elevation_final'+sep+'points'+sep
     69topographies_dir = anuga_dir+'topographies'+sep
     70#topographies_time_dir = topographies_dir+build_time+sep
     71
     72meshes_dir = anuga_dir+'meshes'+sep
     73meshes_dir_name = meshes_dir + scenario_name
     74
     75polygons_dir = anuga_dir+'polygons'+sep
     76tide_dir = anuga_dir+'tide_data'+sep
     77
     78boundaries_in_dir_name = boundaries_in_dir + scenario_name
     79boundaries_dir = anuga_dir+'boundaries'+sep
     80boundaries_dir_name = boundaries_dir + scenario_name
     81
     82#output locations
     83output_dir = anuga_dir+'outputs'+sep
     84output_build_time_dir = anuga_dir+'outputs'+sep+build_time+dir_comment+sep
     85output_run_time_dir = anuga_dir+'outputs'+sep+run_time+dir_comment+sep
     86output_run_time_dir_name = output_run_time_dir + scenario_name  #Used by post pr
    1187
    1288###############################
    1389# Domain definitions
    1490###############################
    15 
    16 # bounding polygon for study area
    17 bounding_polygon = read_polygon('extent.txt')
    18 
    19 print 'Area of bounding polygon', polygon_area(bounding_polygon)/1000000.0
     91from anuga.utilities.polygon import read_polygon, plot_polygons, polygon_area, is_inside_polygon
     92
     93poly_all = read_polygon(polygons_dir+'extent.txt')
     94res_poly_all = 100000*res_factor
     95
    2096
    2197###############################
     
    39115poly_island13= read_polygon('islands13.txt')
    40116
    41 #plot_polygons([bounding_polygon,poly_cairns,poly_island0,poly_island1,\
    42 #               poly_island2,poly_island3,poly_shallow],\
    43 #              'boundingpoly',verbose=False)
     117
     118islands_res = 5000
     119cairns_res = 5000
     120interior_regions = [[project.poly_tongatapu, tongatapu_res],
     121                    [project.poly_island1, islands_res],
     122                    [project.poly_island2, islands_res],
     123                    [project.poly_island3, islands_res],
     124                    [project.poly_island4, islands_res],
     125                    [project.poly_island5, islands_res],
     126                    [project.poly_island6, islands_res],
     127                    [project.poly_island7, islands_res],
     128                    [project.poly_island8, islands_res],
     129                    [project.poly_island9, islands_res],
     130                    [project.poly_island10, islands_res],
     131                    [project.poly_island11, islands_res],
     132                    [project.poly_island12, islands_res],
     133                    [project.poly_island13, islands_res]]
     134
     135boundary_tags={'ocean_east': [0],'ocean_north': [1],'ocean_west': [2],'land': [3,4,5,6,7],'ocean_southeast':[8,9]}
     136
     137trigs_min = number_mesh_triangles(interior_regions, poly_all, res_poly_all)
     138
     139print 'min number triangles', trigs_min
    44140
    45141###################################################################
     
    57153slide_depth = 207.
    58154
    59 gauge_filename = 'gauges.csv'
     155
     156
     157
     158
     159
     160
     161
     162
     163
     164
     165# topo Filenames
     166dem_name = 'tongatapu_10mgrid'
     167meshname = 'tongatapu.msh'
     168
     169# Create DEM from asc data
     170convert_dem_from_ascii2netcdf(dem_name, use_cache=True, verbose=True)
     171
     172# Create pts file for onshore DEM
     173dem2pts(dem_name, use_cache=True, verbose=True)
     174
     175
     176
     177
     178
     179
     180# -*- coding: cp1252 -*-
     181"""Common filenames and locations for topographic data, meshes and outputs.
     182"""
     183
     184from os import sep, environ, getenv, getcwd
     185from os.path import expanduser
     186import sys
     187from time import localtime, strftime, gmtime
     188from anuga.utilities.polygon import read_polygon, plot_polygons, \
     189                                    polygon_area, is_inside_polygon
     190
     191###############################
     192# Domain definitions
     193###############################
     194
     195# bounding polygon for study area
     196bounding_polygon = read_polygon('Y:\data\sw_pacific\tonga\anuga\boundaries\extent.txt')
     197
     198print 'Area of bounding polygon', polygon_area(bounding_polygon)/1000000.0
     199
     200###############################
     201# Interior region definitions
     202###############################
     203
     204# interior polygons
     205poly_tongatapu = read_polygon('tongatapu_finemesharea.txt')
     206poly_island1 = read_polygon('island1.txt')
     207poly_island2 = read_polygon('island2.txt')
     208poly_island3 = read_polygon('islands3.txt')
     209poly_island4 = read_polygon('islands4.txt')
     210poly_island5= read_polygon('islands5.txt')
     211poly_island6= read_polygon('islands6.txt')
     212poly_island7= read_polygon('islands7.txt')
     213poly_island8= read_polygon('islands8.txt')
     214poly_island9= read_polygon('islands9.txt')
     215poly_island10= read_polygon('islands10.txt')
     216poly_island11= read_polygon('islands11.txt')
     217poly_island12= read_polygon('islands12.txt')
     218poly_island13= read_polygon('islands13.txt')
     219
     220#plot_polygons([bounding_polygon,poly_cairns,poly_island0,poly_island1,\
     221#               poly_island2,poly_island3,poly_shallow],\
     222#              'boundingpoly',verbose=False)
     223
     224###################################################################
     225# Clipping regions for export to asc and regions for clipping data
     226###################################################################
     227
     228# exporting asc grid
     229eastingmin = 670500
     230eastingmax = 712750
     231northingmin = 7646000
     232northingmax = 7677000
     233
     234
     235slide_origin = [701290, 7665750] # move onto the continental shelf, depth = 500
     236slide_depth = 207.
     237
     238#gauge_filename = 'gauges.csv'
  • anuga_work/production/tonga/tongatapu.py

    r5134 r5139  
    1 """Script for running a tsunami inundation scenario for Cairns, QLD Australia.
     1"""Script for running tsunami inundation scenario for Dampier, WA, Australia.
    22
    33Source data such as elevation and boundary data is assumed to be available in
    44directories specified by project.py
    5 The output sww file is stored in directory named after the scenario, i.e
    6 slide or fixed_wave.
     5The output sww file is stored in project.output_run_time_dir
    76
    87The scenario is defined by a triangular mesh created from project.polygon,
    9 the elevation data and a tsunami wave generated by a submarine mass failure.
    10 
    11 Ole Nielsen and Duncan Gray, GA - 2005 and Jane Sexton and
    12 Nick Bartzis, GA - 2006
     8the elevation data and a simulated tsunami generated with URS code.
     9
     10Ole Nielsen and Duncan Gray, GA - 2005 and Jane Sexton, Nick Bartzis, GA - 2006
    1311"""
    1412
     
    1816
    1917# Standard modules
    20 import os
     18from os import sep
     19from os.path import dirname, basename
     20from os import mkdir, access, F_OK
     21from shutil import copy
    2122import time
    2223import sys
     
    2425# Related major packages
    2526from anuga.shallow_water import Domain
     27from anuga.shallow_water import Dirichlet_boundary
     28from anuga.shallow_water import File_boundary
    2629from anuga.shallow_water import Reflective_boundary
    27 from anuga.shallow_water import Dirichlet_boundary
    28 from anuga.shallow_water import Time_boundary
    29 from anuga.shallow_water import File_boundary
     30from anuga.shallow_water import Field_boundary
     31from Numeric import allclose
     32from anuga.shallow_water.data_manager import export_grid
     33
    3034from anuga.pmesh.mesh_interface import create_mesh_from_regions
    31 from anuga.shallow_water.data_manager import convert_dem_from_ascii2netcdf
    32 from anuga.shallow_water.data_manager import dem2pts
     35from anuga.shallow_water.data_manager import start_screen_catcher, copy_code_files,store_parameters
     36from anuga_parallel.parallel_api import distribute, numprocs, myid, barrier
     37from anuga_parallel.parallel_abstraction import get_processor_name
     38from anuga.caching import myhash
     39from anuga.damage_modelling.inundation_damage import add_depth_and_momentum2csv, inundation_damage
     40from anuga.fit_interpolate.benchmark_least_squares import mem_usage
    3341
    3442# Application specific imports
    3543import project                 # Definition of file names and polygons
    3644
    37 
    38 #------------------------------------------------------------------------------
    39 # Define scenario as either slide or fixed_wave.
    40 #------------------------------------------------------------------------------
    41 scenario = 'slide'
    42 #scenario = 'fixed_wave'
    43 
    44 if os.access(scenario, os.F_OK) == 0:
    45     os.mkdir(scenario)
    46 basename = scenario + 'source'
    47 
    48 
    49 #------------------------------------------------------------------------------
    50 # Preparation of topographic data
    51 # Convert ASC 2 DEM 2 PTS using source data and store result in source data
    52 #------------------------------------------------------------------------------
    53 
    54 # Filenames
    55 dem_name = 'tongatapu_10mgrid'
    56 meshname = 'tongatapu.msh'
    57 
    58 # Create DEM from asc data
    59 convert_dem_from_ascii2netcdf(dem_name, use_cache=True, verbose=True)
    60 
    61 # Create pts file for onshore DEM
    62 dem2pts(dem_name, use_cache=True, verbose=True)
    63 
    64 
    65 #------------------------------------------------------------------------------
    66 # Create the triangular mesh based on overall clipping polygon with a tagged
    67 # boundary and interior regions defined in project.py along with
    68 # resolutions (maximal area of per triangle) for each polygon
    69 #------------------------------------------------------------------------------
    70 
    71 remainder_res = 1000000
    72 islands_res = 5000
    73 cairns_res = 5000
    74 interior_regions = [[project.poly_tongatapu, tongatapu_res],
    75                     [project.poly_island1, islands_res],
    76                     [project.poly_island2, islands_res],
    77                     [project.poly_island3, islands_res],
    78                     [project.poly_island4, islands_res],
    79                     [project.poly_island5, islands_res],
    80                     [project.poly_island6, islands_res],
    81                     [project.poly_island7, islands_res],
    82                     [project.poly_island8, islands_res],
    83                     [project.poly_island9, islands_res],
    84                     [project.poly_island10, islands_res],
    85                     [project.poly_island11, islands_res],
    86                     [project.poly_island12, islands_res],
    87                     [project.poly_island13, islands_res]]
    88 
    89 create_mesh_from_regions(project.bounding_polygon,
    90                          boundary_tags={'top': [0],
    91                                         'ocean_east': [1],
    92                                         'bottom': [2],
    93                                         'onshore': [3]},
    94                          maximum_triangle_area=remainder_res,
    95                          filename=meshname,
    96                          interior_regions=interior_regions,
    97                          use_cache=True,
    98                          verbose=True)
    99 
    100 
    101 #------------------------------------------------------------------------------
    102 # Setup computational domain
    103 #------------------------------------------------------------------------------
    104 
    105 domain = Domain(meshname, use_cache=True, verbose=True)
    106 
    107 print 'Number of triangles = ', len(domain)
    108 print 'The extent is ', domain.get_extent()
    109 print domain.statistics()
    110 
    111 domain.set_name(basename)
    112 domain.set_datadir(scenario)
    113 domain.set_quantities_to_be_stored(['stage', 'xmomentum', 'ymomentum'])
    114 domain.set_minimum_storable_height(0.01)
    115 
    116 #print 'domain.tight_slope_limiters', domain.tight_slope_limiters
    117 domain.tight_slope_limiters = 0
    118 print 'domain.tight_slope_limiters', domain.tight_slope_limiters
    119 
    120 
    121 
    122 #------------------------------------------------------------------------------
    123 # Setup initial conditions
    124 #------------------------------------------------------------------------------
    125 
    126 tide = 0.0
    127 domain.set_quantity('stage', tide)
    128 domain.set_quantity('friction', 0.0)
    129 domain.set_quantity('elevation',
    130                     filename=dem_name + '.pts',
    131                     use_cache=True,
    132                     verbose=True,
    133                     alpha=0.1)
    134 
    135 
    136 
    137 #------------------------------------------------------------------------------
    138 # Setup information for slide scenario (to be applied 1 min into simulation
    139 #------------------------------------------------------------------------------
    140 
    141 if scenario == 'slide':
    142     # Function for submarine slide
    143     from anuga.shallow_water.smf import slide_tsunami 
    144     tsunami_source = slide_tsunami(length=35000.0,
    145                                    depth=project.slide_depth,
    146                                    slope=6.0,
    147                                    thickness=500.0,
    148                                    x0=project.slide_origin[0],
    149                                    y0=project.slide_origin[1],
    150                                    alpha=0.0,
    151                                    domain=domain,
    152                                    verbose=True)
    153 
    154 
    155 #------------------------------------------------------------------------------
    156 # Setup boundary conditions
    157 #------------------------------------------------------------------------------
    158 
    159 print 'Available boundary tags', domain.get_boundary_tags()
    160 
    161 Br = Reflective_boundary(domain)
    162 Bd = Dirichlet_boundary([tide,0,0])
    163 
    164 # 60 min square wave starting at 1 min, 50m high
    165 if scenario == 'fixed_wave':
     45def run_model(**kwargs):
     46   
     47
     48    #------------------------------------------------------------------------------
     49    # Copy scripts to time stamped output directory and capture screen
     50    # output to file
     51    #------------------------------------------------------------------------------
     52    print "Processor Name:",get_processor_name()
     53
     54    #copy script must be before screen_catcher
     55    #print kwargs
     56
     57    print 'output_dir',kwargs['output_dir']
     58    if myid == 0:
     59        copy_code_files(kwargs['output_dir'],__file__,
     60                 dirname(project.__file__)+sep+ project.__name__+'.py' )
     61
     62        store_parameters(**kwargs)
     63
     64    barrier()
     65
     66    start_screen_catcher(kwargs['output_dir'], myid, numprocs)
     67
     68    print "Processor Name:",get_processor_name()
     69
     70    # filenames
     71#    meshes_dir_name = project.meshes_dir_name+'.msh'
     72
     73    # creates copy of code in output dir
     74    print 'min triangles', project.trigs_min,
     75    print 'Note: This is generally about 20% less than the final amount'
     76
     77    #--------------------------------------------------------------------------
     78    # Create the triangular mesh based on overall clipping polygon with a
     79    # tagged
     80    # boundary and interior regions defined in project.py along with
     81    # resolutions (maximal area of per triangle) for each polygon
     82    #--------------------------------------------------------------------------
     83
     84    #IMPORTANT don't cache create_mesh_from_region and Domain(mesh....) as it
     85    # causes problems with the ability to cache set quantity which takes alot of times
     86   
     87    if myid == 0:
     88   
     89        print 'start create mesh from regions'
     90
     91        create_mesh_from_regions(project.poly_all,
     92                             boundary_tags=project.boundary_tags,
     93                             maximum_triangle_area=project.res_poly_all,
     94                             interior_regions=project.interior_regions,
     95                             filename=project.meshes_dir_name+'.msh',
     96                             use_cache=False,
     97                             verbose=True)
     98    barrier()
     99
     100    #-------------------------------------------------------------------------
     101    # Setup computational domain
     102    #-------------------------------------------------------------------------
     103    print 'Setup computational domain'
     104
     105    #domain = cache(Domain, (meshes_dir_name), {'use_cache':True, 'verbose':True}, verbose=True)
     106    #above don't work
     107    domain = Domain(project.meshes_dir_name+'.msh', use_cache=False, verbose=True)
     108    print 'memory usage before del domain',mem_usage()
     109       
     110    print domain.statistics()
     111    print 'triangles',len(domain)
     112   
     113    kwargs['act_num_trigs']=len(domain)
     114   
     115    #-------------------------------------------------------------------------
     116    # Setup initial conditions
     117    #-------------------------------------------------------------------------
     118    if myid == 0:
     119
     120        print 'Setup initial conditions'
     121
     122        from polygon import Polygon_function
     123        #following sets the stage/water to be offcoast only
     124#        IC = Polygon_function( [(project.poly_mainland, -1.0)], default = kwargs['tide'],
     125#                                 geo_reference = domain.geo_reference)
     126#        domain.set_quantity('stage', IC)
     127        domain.set_quantity('stage',kwargs['tide'] )
     128#        domain.set_quantity('stage', kwargs['tide'])
     129        domain.set_quantity('friction', kwargs['friction'])
     130       
     131        print 'Start Set quantity',kwargs['bathy_file']
     132
     133        domain.set_quantity('elevation',
     134                            filename = kwargs['bathy_file'],
     135                            use_cache = True,
     136                            verbose = True,
     137                            alpha = kwargs['alpha'])
     138        print 'Finished Set quantity'
     139    barrier()
     140   
     141    #------------------------------------------------------
     142    # Distribute domain to implement parallelism !!!
     143    #------------------------------------------------------
     144
     145    if numprocs > 1:
     146        domain=distribute(domain)
     147
     148    #------------------------------------------------------
     149    # Set domain parameters
     150    #------------------------------------------------------
     151    print 'domain id', id(domain)
     152    domain.set_name(kwargs['aa_scenario_name'])
     153    domain.set_datadir(kwargs['output_dir'])
     154    domain.set_default_order(2) # Apply second order scheme
     155    domain.set_minimum_storable_height(0.01) # Don't store anything less than 1cm
     156    domain.set_store_vertices_uniquely(False)
     157    domain.set_quantities_to_be_stored(['stage', 'xmomentum', 'ymomentum'])
     158    domain.tight_slope_limiters = 1
     159    #domain.set_maximum_allowed_speed(0.1) # Allow a little runoff (0.1 is OK)
     160    print 'domain id', id(domain)
     161    domain.beta_h = 0
     162
     163    #------------------------------------------------------------------------------
     164    # Setup information for slide scenario (to be applied 1 min into simulation
     165    #------------------------------------------------------------------------------
     166
     167    if scenario == 'slide':
     168        # Function for submarine slide
     169        from anuga.shallow_water.smf import slide_tsunami 
     170        tsunami_source = slide_tsunami(length=35000.0,
     171                                       depth=project.slide_depth,
     172                                       slope=6.0,
     173                                       thickness=500.0,
     174                                       x0=project.slide_origin[0],
     175                                       y0=project.slide_origin[1],
     176                                       alpha=0.0,
     177                                       domain=domain,
     178                                       verbose=True)
     179   
     180    #-------------------------------------------------------------------------
     181    # Setup boundary conditions
     182    #-------------------------------------------------------------------------
     183    print 'Available boundary tags', domain.get_boundary_tags()
     184    print 'domain id', id(domain)
     185    #print 'Reading Boundary file',project.boundaries_dir_namea + '.sww'
     186
     187    Br = Reflective_boundary(domain)
     188    Bd = Dirichlet_boundary([kwargs['tide'],0,0])
     189    Bo = Dirichlet_boundary([kwargs['tide']+5.0,0,0])
     190
     191    if scenario == 'fixed_wave':
    166192    # FIXME (Ole): Change this to Transmissive Momentum as
    167193    # suggested by Rajaraman (ticket:208)
    168     Bw = Time_boundary(domain = domain,
    169                        f=lambda t: [(60<t<3660)*50, 0, 0])
    170     domain.set_boundary({'ocean_east': Bw,
    171                          'bottom': Bd,
    172                          'onshore': Bd,
    173                          'top': Bd})
    174 
    175 # boundary conditions for slide scenario
    176 if scenario == 'slide':
    177     domain.set_boundary({'ocean_east': Bd,
    178                          'bottom': Bd,
    179                          'onshore': Bd,
    180                          'top': Bd})
    181    
    182 
    183 #------------------------------------------------------------------------------
    184 # Evolve system through time
    185 #------------------------------------------------------------------------------
    186 
    187 import time
    188 t0 = time.time()
    189 
    190 from Numeric import allclose
    191 from anuga.abstract_2d_finite_volumes.quantity import Quantity
    192 
    193 if scenario == 'slide':
    194    
    195     for t in domain.evolve(yieldstep = 10, finaltime = 60):
     194        Bw = Time_boundary(domain = domain,
     195                           f=lambda t: [(60<t<3660)*100, 0, 0])
     196        domain.set_boundary({'ocean_east': Bw,
     197                             'ocean_southeast': Bd,
     198                             'land': Bd,
     199                             'ocean_north': Bd,
     200                             'ocean_west':})
     201
     202    # boundary conditions for slide scenario
     203    if scenario == 'slide':
     204        domain.set_boundary({'ocean_east': Bd,
     205                             'ocean_southeast': Bd,
     206                             'land': Bd,
     207                             'ocean_north': Bd,
     208                             'ocean_west':})
     209
     210    #----------------------------------------------------------------------------
     211    # Evolve system through time
     212    #--------------------------------------------------------------------
     213    t0 = time.time()
     214
     215    for t in domain.evolve(yieldstep = 240, finaltime = kwargs['finaltime']):
    196216        domain.write_time()
    197217        domain.write_boundary_statistics(tags = 'ocean_east')     
    198        
    199     # add slide
    200     thisstagestep = domain.get_quantity('stage')
    201     if allclose(t, 60):
    202         slide = Quantity(domain)
    203         slide.set_values(tsunami_source)
    204         domain.set_quantity('stage', slide + thisstagestep)
    205 
    206     for t in domain.evolve(yieldstep = 10, finaltime = 5000,
    207                            skip_initial_step = True):
    208         domain.write_time()
    209         domain.write_boundary_statistics(tags = 'ocean_east')
    210 
    211 if scenario == 'fixed_wave':
    212 
    213     # save every two mins leading up to wave approaching land
    214     for t in domain.evolve(yieldstep = 120, finaltime = 5000):
    215         domain.write_time()
    216         domain.write_boundary_statistics(tags = 'ocean_east')     
    217 
    218     # save every 30 secs as wave starts inundating ashore
    219     for t in domain.evolve(yieldstep = 10, finaltime = 10000,
    220                            skip_initial_step = True):
    221         domain.write_time()
    222         domain.write_boundary_statistics(tags = 'ocean_east')
    223            
    224 print 'That took %.2f seconds' %(time.time()-t0)
     218
     219    x, y = domain.get_maximum_inundation_location()
     220    q = domain.get_maximum_inundation_elevation()
     221
     222    print 'Maximum runup observed at (%.2f, %.2f) with elevation %.2f' %(x,y,q)
     223
     224    print 'That took %.2f seconds' %(time.time()-t0)
     225
     226    #kwargs 'completed' must be added to write the final parameters to file
     227    kwargs['completed']=str(time.time()-t0)
     228   
     229    if myid==0:
     230        store_parameters(**kwargs)
     231    barrier
     232   
     233    print 'memory usage before del domain1',mem_usage()
     234
     235def export_model(**kwargs):   
     236    #store_parameters(**kwargs)
     237   
     238#    print 'memory usage before del domain',mem_usage()
     239    #del domain
     240    print 'memory usage after del domain',mem_usage()
     241   
     242    swwfile = kwargs['output_dir']+kwargs['scenario_name']
     243    print'swwfile',swwfile
     244   
     245    export_grid(swwfile, extra_name_out = 'town',
     246            quantities = ['speed','depth','elevation','stage'], # '(xmomentum**2 + ymomentum**2)**0.5' defaults to elevation
     247            #quantities = ['speed','depth'], # '(xmomentum**2 + ymomentum**2)**0.5' defaults to elevation
     248            timestep = None,
     249            reduction = max,
     250            cellsize = kwargs['export_cellsize'],
     251            NODATA_value = -1E-030,
     252            easting_min = project.eastingmin,
     253            easting_max = project.eastingmax,
     254            northing_min = project.northingmin,
     255            northing_max = project.northingmax,
     256            verbose = False,
     257            origin = None,
     258            datum = 'WGS84',
     259            format = 'txt')
     260   
     261#-------------------------------------------------------------
     262if __name__ == "__main__":
     263   
     264    kwargs={}
     265    kwargs['est_num_trigs']=project.trigs_min
     266    kwargs['num_cpu']=numprocs
     267    kwargs['host']=project.host
     268    kwargs['res_factor']=project.res_factor
     269    kwargs['starttime']=project.starttime
     270    kwargs['yieldstep']=project.yieldstep
     271    kwargs['midtime']=project.midtime
     272    kwargs['finaltime']=project.finaltime
     273   
     274    kwargs['output_dir']=project.output_run_time_dir
     275    kwargs['bathy_file']=project.combined_dir_name+'.txt'
     276#    kwargs['bathy_file']=project.combined_small_dir_name + '.pts'
     277    kwargs['boundary_file']=project.boundaries_in_dir_name + '.sww'
     278    kwargs['file_name']=project.home+'detail.csv'
     279    kwargs['aa_scenario_name']=project.scenario_name
     280    kwargs['ab_time']=project.time
     281    kwargs['res_factor']= project.res_factor
     282    kwargs['tide']=project.tide
     283    kwargs['user']=project.user
     284    kwargs['alpha'] = project.alpha
     285    kwargs['friction']=project.friction
     286    kwargs['time_thinning'] = project.time_thinning
     287    kwargs['dir_comment']=project.dir_comment
     288    kwargs['export_cellsize']=project.export_cellsize
     289   
     290
     291    run_model(**kwargs)
     292     
     293    if myid==0:
     294        export_model(**kwargs)
     295    barrier
Note: See TracChangeset for help on using the changeset viewer.