Changeset 5386


Ignore:
Timestamp:
May 30, 2008, 3:34:13 PM (16 years ago)
Author:
kristy
Message:

cleaned Perth python files

Location:
anuga_work/production/perth
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • anuga_work/production/perth/build_perth.py

    r5360 r5386  
    117117
    118118print'add all geospatial objects'
    119 G = G1 + G2 + G3 + G4 + G5 + G6 + G_off + G_off1
     119G = G1 + G2 + G3 + G4 + G5 + G6 + G_off
    120120
    121121print'clip combined geospatial object by bounding polygon'
    122 G_clipped = G.clip(project.bounding_polygon)
    123 #FIXME: add a clip function to pts
    124 #print'shape of clipped data', G_clipped.get_data_points().shape
     122G_clipped = G.clip(project.poly_all)
    125123
    126124print'export combined DEM file'
    … …  
    128126    mkdir (project.topographies_dir)
    129127G_clipped.export_points_file(project.combined_dir_name + '.pts')
    130 #G_clipped.export_points_file(project.combined_dir_name + '.xya')
    131 
    132 '''
    133 print'project.combined_dir_name + 1.xya',project.combined_dir_name + '1.xya'
    134 G_all=Geospatial_data(file_name = project.combined_dir_name + '1.xya')
    135 print'split'
    136 G_all_1, G_all_2 = G_all.split(.10)
    137 print'export 1'
    138 G_all_1.export_points_file(project.combined_dir_name+'_small1' + '.xya')
    139 print'export 2'
    140 G_all_2.export_points_file(project.combined_dir_name+'_other1' + '.xya')
    141 
     128#G_clipped.export_points_file(project.combined_dir_name + '.xya') - use for comparision in ARC
    142129
    143130#-------------------------------------------------------------------------
    … …  
    145132#-------------------------------------------------------------------------
    146133print 'starting to create boundary conditions'
     134from anuga.shallow_water.data_manager import urs2sww, urs_ungridded2sww
     135
    147136boundaries_in_dir_name = project.boundaries_in_dir_name
     137print 'boundaries_in_dir_name',project.boundaries_in_dir_name
    148138
    149 from anuga.shallow_water.data_manager import urs2sww
    150 
    151 print 'minlat=project.north_boundary, maxlat=project.south_boundary',project.north_boundary, project.south_boundary
    152 print 'minlon= project.west_boundary, maxlon=project.east_boundary',project.west_boundary, project.east_boundary
    153139
    154140#import sys; sys.exit()
    155141
    156 #if access(project.boundaries_dir,F_OK) == 0:
    157 #    mkdir (project.boundaries_dir)
    158 
    159 from caching import cache
    160 cache(urs2sww,
    161       (boundaries_in_dir_name,
    162        project.boundaries_dir_name1),
    163       {'verbose': True,
    164        'minlat': project.south_boundary,
    165        'maxlat': project.north_boundary,
    166        'minlon': project.west_boundary,
    167        'maxlon': project.east_boundary,
    168 #       'minlat': project.south,
    169 #       'maxlat': project.north,
    170 #       'minlon': project.west,
    171 #       'maxlon': project.east,
    172        'mint': 0, 'maxt': 40000,
    173 #       'origin': domain.geo_reference.get_origin(),
    174        'mean_stage': project.tide,
    175 #       'zscale': 1,                 #Enhance tsunami
    176        'fail_on_NaN': False},
    177        verbose = True,
    178        )
    179 #       dependencies = source_dir + project.boundary_basename + '.sww')
    180 
    181 '''
     142urs_ungridded2sww(project.boundaries_in_dir_name, project.boundaries_in_dir_name,
     143                  verbose=True, mint=4000, maxt=10000, zscale=1)
    182144
    183145
    … …  
    187149
    188150
     151
     152
     153

  • anuga_work/production/perth/project.py

    r5381 r5386  
    4242export_cellsize=50
    4343setup='final'
    44 source='test'
     44source='exmouth'
    4545
    4646if setup =='trial':
    … …  
    119119
    120120#boundaries_source = '1'
    121 
    122 if source =='dampier':
    123     boundaries_name = 'broome_3854_17042007' #Dampier gun
    124     boundaries_in_dir = anuga_dir+'boundaries'+sep+'urs'+sep+'dampier'+sep+'1_10000'+sep
    125 
    126 if source=='onslow':
    127     boundaries_name = 'broome_3859_16052007' #onslow_hedland_broome gun
    128     boundaries_in_dir = anuga_dir+'boundaries'+sep+'urs'+sep+'onslow_hedland_broome'+sep+'1_10000'+sep
    129121   
    130122if source=='exmouth':
    131     boundaries_name = 'broome_3103_18052007' #exmouth gun
     123    boundaries_name = 'perth_3103_28052008' #exmouth gun
    132124    boundaries_in_dir = anuga_dir+'boundaries'+sep+'urs'+sep+'exmouth'+sep+'1_10000'+sep
    133125
    … …  
    209201                     ,[poly_rottnest_ex, res_rottnest_ex]]
    210202
    211 boundary_tags={'back': [0,1,2], 'side': [3,6],'ocean': [4, 5]}
     203boundary_tags={'back': [0,1,2], 'side': [3,7],'ocean': [4, 5, 6]}
    212204
    213205trigs_min = number_mesh_triangles(interior_regions, poly_all, res_poly_all)
    214206
    215207print 'min number triangles', trigs_min
     208
     209poly_mainland = read_polygon(polygons_dir+'Initial_Condition.csv')
    216210
    217211###################################################################

  • anuga_work/production/perth/run_perth.py

    r5360 r5386  
    103103    print 'Setup computational domain'
    104104
    105     #domain = cache(Domain, (meshes_dir_name), {'use_cache':True, 'verbose':True}, verbose=True)
    106     #above don't work
    107105    domain = Domain(project.meshes_dir_name+'.msh', use_cache=False, verbose=True)
    108106    print 'memory usage before del domain',mem_usage()
    … …  
    123121        from polygon import Polygon_function
    124122        #following sets the stage/water to be offcoast only
    125 #        IC = Polygon_function( [(project.poly_mainland, -1.0)], default = kwargs['tide'],
    126 #                                 geo_reference = domain.geo_reference)
    127 #        domain.set_quantity('stage', IC)
    128         domain.set_quantity('stage',kwargs['tide'] )
    129 #        domain.set_quantity('stage', kwargs['tide'])
     123        IC = Polygon_function( [(project.poly_mainland, 0)], default = kwargs['tide'],
     124                                 geo_reference = domain.geo_reference)
     125        domain.set_quantity('stage', IC)
     126        #domain.set_quantity('stage',kwargs['tide'] )
    130127        domain.set_quantity('friction', kwargs['friction'])
    131128       
    132         print 'Start Set quantity',kwargs['bathy_file']
     129        print 'Start Set quantity',kwargs['elevation_file']
    133130
    134131        domain.set_quantity('elevation',
    135                             filename = kwargs['bathy_file'],
     132                            filename = kwargs['elevation_file'],
    136133                            use_cache = False,
    137134                            verbose = True,
    … …  
    172169    Br = Reflective_boundary(domain)
    173170    Bd = Dirichlet_boundary([kwargs['tide'],0,0])
    174     Bo = Dirichlet_boundary([kwargs['tide']+5.0,0,0])
    175 
    176     if project.source != 'test':
    177         print 'start reading boundary file'
    178 
    179         Bf = Field_boundary(kwargs['boundary_file'],
    180                     domain, time_thinning=kwargs['time_thinning'], mean_stage=kwargs['tide'],
    181                     use_cache=False, verbose=True)
    182         domain.set_boundary({'back': Br,
    183                              'side': Bd,
    184                              'ocean': Bf})
    185     else:
    186         domain.set_boundary({'back': Br,
    187                              'side': Bd,
    188                              'ocean': Bd})
     171
     172    print 'start reading boundary file'
     173
     174    Bf = Field_boundary(kwargs['boundary_file'],
     175                domain, time_thinning=kwargs['time_thinning'], mean_stage=kwargs['tide'],
     176                use_cache=False, verbose=True)
     177    domain.set_boundary({'back': Br,
     178                         'side': Bd,
     179                         'ocean': Bf})
     180
     181##    if project.source != 'test':
     182##        domain.set_boundary({'back': Br,
     183##                             'side': Bd,
     184##                             'ocean': Bd})
    189185
    190186    kwargs['input_start_time']=domain.starttime
    … …  
    197193    t0 = time.time()
    198194
    199     for t in domain.evolve(yieldstep = 240, finaltime = kwargs['starttime']):
     195    for t in domain.evolve(yieldstep = 240, finaltime = kwargs['finaltime']
     196                       ,skip_initial_step = True):
    200197        domain.write_time()
    201         domain.write_boundary_statistics(tags = 'ocean')     
    202 
    203         if allclose(t, 120):
    204             domain.set_boundary({'back': Br, 'side': Bd, 'ocean': Bo})
    205 
    206         if allclose(t, 720):
    207             domain.set_boundary({'back': Br, 'side': Bd, 'ocean': Bd})
    208 
    209 #    for t in domain.evolve(yieldstep = kwargs['yieldstep'], finaltime = kwargs['midtime']
    210 #                       ,skip_initial_step = True):
    211 #        domain.write_time()
    212 #        domain.write_boundary_statistics(tags = 'ocean')     
    213 #   
    214 #    for t in domain.evolve(yieldstep = 240, finaltime = kwargs['finaltime']
    215 #                       ,skip_initial_step = True):
    216 #        domain.write_time()
    217 #        domain.write_boundary_statistics(tags = 'ocean')   
     198        domain.write_boundary_statistics(tags = 'ocean')   
    218199
    219200    x, y = domain.get_maximum_inundation_location()
    … …  
    233214    print 'memory usage before del domain1',mem_usage()
    234215   
    235 def 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['aa_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 = -9999,
    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 = 'GDA94',
    259             format = 'asc')
    260            
    261     inundation_damage(swwfile+'.sww', project.buildings_filename,
    262                       project.buildings_filename_out)
    263216   
    264217#-------------------------------------------------------------
    … …  
    276229   
    277230    kwargs['output_dir']=project.output_run_time_dir
    278     kwargs['bathy_file']=project.combined_dir_name+'.txt'
    279 #    kwargs['bathy_file']=project.combined_small_dir_name + '.pts'
     231    kwargs['elevation_file']=project.combined_dir_name+'.txt'
     232#    kwargs['elevation_file']=project.combined_small_dir_name + '.pts'
    280233    kwargs['boundary_file']=project.boundaries_in_dir_name + '.sww'
    281234    kwargs['file_name']=project.home+'detail.csv'
Note: See TracChangeset for help on using the changeset viewer.