Changeset 5751


Ignore:
Timestamp:
Sep 9, 2008, 2:14:12 PM (14 years ago)
Author:
kristy
Message:
 
Location:
anuga_work/production/geraldton
Files:
4 edited

Legend:

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

    r5732 r5751  
    9191print'finished reading in files'
    9292
     93G_off1_clip = G_off1.clip(project.poly_CBD_1km, verbose=True)
     94
    9395print'add all geospatial objects'
    94 G = G1 + G2 + G3 + G_off + G_off1 + G_off2 + G_off3
     96G = G1 + G2 + G3 + G_off + G_off1_clip + G_off2 + G_off3
    9597
    9698print'clip combined geospatial object by bounding polygon'
  • anuga_work/production/geraldton/get_gauges.py

    r5703 r5751  
    1414        compress,zeros,fabs,take
    1515    fid = NetCDFFile(filename+'.sts', 'r')    #Open existing file for read
     16    permutation = fid.variables['permutation'][:]
    1617    x = fid.variables['x'][:]+fid.xllcorner   #x-coordinates of vertices
    1718    y = fid.variables['y'][:]+fid.yllcorner   #y-coordinates of vertices
     
    3132    fid_max.write(s)   
    3233    for j in range(len(x)):
    33         locx=int(x[j])
    34         locy=int(y[j])
     34        index= permutation[j]
    3535        stage = quantities['stage'][:,j]
    3636        xmomentum = quantities['xmomentum'][:,j]
    3737        ymomentum = quantities['ymomentum'][:,j]
    3838
    39         s = '%.6f, %.6f, %.6f\n' %(x[j], y[j], max(stage))
     39        s = '%d, %.6f, %.6f, %.6f\n' %(index, x[j], y[j], max(stage))
    4040        fid_max.write(s)
    4141       
    42         fid_sts = open(project.boundaries_dir+sep+basename+'_'+str(locx)+'_'+str(locy)+'.csv', 'w')
     42        fid_sts = open(project.boundaries_dir+sep+basename+'_'+ str(index)+'.csv', 'w')
    4343        s = 'time, stage, xmomentum, ymomentum \n'
    4444        fid_sts.write(s)
  • anuga_work/production/geraldton/project.py

    r5702 r5751  
    152152
    153153poly_neg20_pos20 = read_polygon(polygons_dir+'neg20_pos20.csv')
    154 res_neg20_pos20 = 20000*res_factor
     154res_neg20_pos20 = 25000*res_factor
    155155
    156156poly_CBD_1km= read_polygon(polygons_dir+'CBD_1km.csv')
     
    161161
    162162poly_wallabi = read_polygon(polygons_dir+'island_wallabi_poly.csv')
    163 res_wallabi = 1000*res_factor
     163res_wallabi = 5000*res_factor
    164164
    165165poly_dingiville = read_polygon(polygons_dir+'island_dingiville_poly.csv')
    166 res_dingiville = 1000*res_factor
     166res_dingiville = 5000*res_factor
    167167
    168168poly_pelsaert = read_polygon(polygons_dir+'island_pelsaert_poly.csv')
    169 res_pelsaert = 1000*res_factor
     169res_pelsaert = 5000*res_factor
    170170
    171171
     
    179179print 'min number triangles', trigs_min
    180180
    181 #For no input boundary file
     181##For no input boundary file
     182#boundary_tags={'back': [2,3,4,5], 'side': [0,1,6], 'ocean': [7,8,9]}
    182183
    183 boundary_tags={'back': [2,3,4,5], 'side': [0,1,6], 'ocean': [7,8,9]}
    184184poly_mainland=read_polygon(topographies_in_dir +'initial_condition.csv')
    185185
  • anuga_work/production/geraldton/run_geraldton.py

    r5669 r5751  
    7878    # Domain definitions
    7979    #-----------------------------------------------------------------------
    80 ##
    81 ##    # Read in boundary from ordered sts file
    82 ##    urs_bounding_polygon=create_sts_boundary(os.path.join(project.boundaries_dir,project.scenario_name))
    83 ##
    84 ##    # Reading the landward defined points, this incorporates the original clipping
    85 ##    # polygon minus the 100m contour
    86 ##    landward_bounding_polygon = read_polygon(project.boundaries_dir+'landward_boundary.txt')
    87 ##
    88 ##    # Combine sts polyline with landward points
    89 ##    bounding_polygon = urs_bounding_polygon + landward_bounding_polygon
    90 ##   
    91 ##    # counting segments
    92 ##    N = len(urs_bounding_polygon)-1
    93 ##    boundary_tags={'back': [N+2,N+3], 'side': [N,N+1,N+4],'ocean': range(N)}
    94 
    95     bounding_polygon = project.poly_all
    96    
    97 
    98    
     80
     81    # Read in boundary from ordered sts file
     82    urs_bounding_polygon=create_sts_boundary(os.path.join(project.boundaries_dir,project.scenario_name))
     83
     84    # Reading the landward defined points, this incorporates the original clipping
     85    # polygon minus the 100m contour
     86    landward_bounding_polygon = read_polygon(project.boundaries_dir+'landward_boundary_polygon.txt')
     87
     88    # Combine sts polyline with landward points
     89    bounding_polygon = urs_bounding_polygon + landward_bounding_polygon
     90   
     91    # counting segments
     92    N = len(urs_bounding_polygon)-1
     93    boundary_tags={'back': [N+2,N+3], 'side': [N,N+1,N+4],'ocean': range(N)}
     94
     95    print 'boundary tags',boundary_tags
     96       
    9997    #--------------------------------------------------------------------------
    10098    # Create the triangular mesh based on overall clipping polygon with a
     
    111109
    112110        create_mesh_from_regions(bounding_polygon,
    113                              boundary_tags=project.boundary_tags,
     111                             boundary_tags=boundary_tags,
    114112                             maximum_triangle_area=project.res_poly_all,
    115113                             interior_regions=project.interior_regions,
     
    208206   
    209207    boundary_urs_out=project.boundaries_dir_name
     208
     209    Br = Reflective_boundary(domain)
     210    Bd = Dirichlet_boundary([kwargs['tide'],0,0])
    210211   
    211212    print 'Available boundary tags', domain.get_boundary_tags()
     
    213214                   domain, mean_stage= project.tide,
    214215                   time_thinning=1,
     216                   default_boundary=Bd,
    215217                   use_cache=True,
    216218                   verbose = True,
    217219                   boundary_polygon=bounding_polygon)
    218220
    219    
    220     Br = Reflective_boundary(domain)
    221     Bd = Dirichlet_boundary([kwargs['tide'],0,0])
    222 
    223 ##    Bf = Field_boundary(kwargs['boundary_file'],
    224 ##                domain, time_thinning=kwargs['time_thinning'], mean_stage=kwargs['tide'],
    225 ##                use_cache=False, verbose=True)
    226 
    227221    domain.set_boundary({'back': Br,
    228222                         'side': Bd,
    229                          'ocean': Bd}) # change baxk to Bf when running properly
     223                         'ocean': Bf})
    230224
    231225    kwargs['input_start_time']=domain.starttime
     
    290284    run_model(**kwargs)
    291285     
    292     if myid==0:
    293         export_model(**kwargs)
    294     #barrier
     286       #barrier
Note: See TracChangeset for help on using the changeset viewer.