Changeset 4465 for anuga_work


Ignore:
Timestamp:
May 21, 2007, 10:59:13 AM (18 years ago)
Author:
nick
Message:

updates to exmouth

Location:
anuga_work/production/exmouth_2006
Files:
2 added
1 edited

Legend:

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

    r4423 r4465  
    1 # -*- coding: cp1252 -*-
    21"""Common filenames and locations for topographic data, meshes and outputs.
     2Also includes origin for slump scenario.
    33"""
    44
    5 from os import sep, environ, getenv, getcwd
    6 from os.path import expanduser
     5from os import sep, environ, getenv, getcwd,umask
     6from os.path import expanduser, basename
     7from anuga.utilities.polygon import read_polygon, plot_polygons, polygon_area, is_inside_polygon, number_mesh_triangles
    78import sys
     9from anuga.coordinate_transforms.redfearn import degminsec2decimal_degrees
    810from time import localtime, strftime, gmtime
    9 from anuga.utilities.polygon import read_polygon, plot_polygons, polygon_area, is_inside_polygon
    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
    1112
    12 if sys.platform == 'win32':
    13     home = getenv('INUNDATIONHOME')
    14     user = getenv('USERPROFILE')
     13codename = 'project.py'
    1514
    16 else:   
    17     home = getenv('INUNDATIONHOME', sep+'d'+sep+'xrd'+sep+'gem'+sep+'2'+sep+'ramp'+sep+'risk_assessment_methods_project'+sep+'inundation')     
    18     user = getenv('LOGNAME')
    19     print 'USER:', user
    20 
    21 # INUNDATIONHOME is the inundation directory, not the data directory.
    22 home += sep +'data'
     15home = getenv('INUNDATIONHOME') + sep +'data'+sep #Sandpit's parent dir   
     16user = get_user_name()
     17host = get_host_name()
     18#needed when running using mpirun, mpirun doesn't inherit umask from .bashrc
     19umask(002)
    2320
    2421#Making assumptions about the location of scenario data
    2522state = 'western_australia'
    26 scenario_dir_name = 'exmouth_tsunami_scenario_2006'
     23scenario_name = 'exmouth'
     24scenario = 'exmouth_tsunami_scenario'
    2725
    28 # onshore data provided by WA DLI
    29 onshore_name =  'first' # original
     26#time stuff
     27time = strftime('%Y%m%d_%H%M%S',gmtime()) #gets time for new dir
     28build_time = time+'_build'
     29run_time = time+'_run'
    3030
    31 # AHO + DPI data
    32 offshore_name1 = 'XY100011610'
    33 offshore_name2 = 'XY100011611'
    34 offshore_name3 = 'XY100011613'
    35 offshore_name4 = 'XY100011614'
    36 offshore_name5 = 'XY100011616'
    37 offshore_name6 = 'XY100011617'
    38 offshore_name7 = 'XY100011618'
    39 offshore_name8 = 'XY100011621'
    40 offshore_name9 = 'XY100011623'
    41 offshore_name10 = 'XY100011745'
    42 offshore_name11 = 'XY100011746'
    43 offshore_name12 = 'XY100017530'
    44 offshore_name13 = 'XY100017532'
    45 offshore_name14 = 'XY100017538'
    46 offshore_name15 = 'XY100017540'
    47 offshore_name16 = 'XYBR66'
    48 offshore_name17 = 'XYBR70'
    49 offshore_name18 = 'XYBR80'
    50 offshore_name19 = 'XYBR88'
    51 offshore_name20 = 'XYBR93'
    52 offshore_name21 = 'XYBR0110'
    53 offshore_name22 = 'XYWADPI'
     31#tide = -1.4
     32#tide = 0
     33tide = 1.4
    5434
    55 # developed by NM&I
    56 coast_name = 'coastline'
     35#Maybe will try to make project a class to allow these parameters to be passed in.
     36alpha = 0.1
     37friction=0.01
     38starttime=3600
     39finaltime=25000
     40setup='trial'
    5741
    58 boundary_basename = 'SU-AU' # Mw ?
     42if setup =='trial':
     43    print'trial'
     44    res_factor=10
     45    time_thinning=48
     46    yieldstep=240
     47if setup =='basic':
     48    print'basic'
     49    res_factor=4
     50    time_thinning=12
     51    yieldstep=120
     52if setup =='final':
     53    print'final'
     54    res_factor=1
     55    time_thinning=1
     56    yieldstep=60
    5957
    60 #swollen/ all data output
    61 basename = 'source'
    62 codename = 'project.py'
     58dir_comment='_'+setup+'_'+str(tide)+'_'+str(user)
    6359
    64 #Derive subdirectories and filenames
    65 local_time = strftime('%Y%m%d_%H%M%S',gmtime())
    66 meshdir = home+sep+state+sep+scenario_dir_name+sep+'anuga'+sep+'meshes'+sep
    67 datadir = home+sep+state+sep+scenario_dir_name+sep+'anuga'+sep+'topographies'+sep
    68 gaugedir = home+sep+state+sep+scenario_dir_name+sep+'anuga'+sep+'gauges'+sep
    69 polygondir = home+sep+state+sep+scenario_dir_name+sep+'anuga'+sep+'polygons'+sep
    70 boundarydir = home+sep+sep+state+sep+scenario_dir_name+sep+'anuga'+sep+'boundaries'+sep
    71 outputdir = home+sep+state+sep+scenario_dir_name+sep+'anuga'+sep+'outputs'+sep
    72 outputtimedir = outputdir + local_time + sep
    73 polygondir = home+sep+state+sep+scenario_dir_name+sep+'anuga'+sep+'polygons'+sep
     60onshore_name = 'DLI'
     61onshore_name1 = 'DTED'
     62# offshore
     63offshore_name = 'Exmouth_bathymetry'
     64#offshore_name1 = 'inferred_north'
     65#offshore_name2 = 'inferred_south'
     66coast_name = 'Exmouth_coastline'
    7467
    75 gauge_filename = gaugedir + 'exmouth_gauges.csv'
    76 community_filename = gaugedir + 'CHINS_v2.csv'
    77 community_exmouth = gaugedir + 'community_exmouth.csv'
    78 codedir = getcwd()+sep                           
    79 codedirname = codedir + 'project.py'
    80 meshname = outputtimedir + 'mesh_' + basename
     68#final topo name
     69combined_name ='exmouth_combined_elevation'
     70#combined_name1 ='exmouth_combined_elevation1'
     71#combined_name_unclipped1 ='exmouth_combined_elevation_unclipped1'
     72combined_small_name = 'exmouth_combined_elevation_small'
    8173
    82 # Necessary if using point datasets, rather than grid
    83 onshore_dem_name = datadir + onshore_name
    84 offshore_dem_name1 = datadir + offshore_name1
    85 offshore_dem_name2 = datadir + offshore_name2
    86 offshore_dem_name3 = datadir + offshore_name3
    87 offshore_dem_name4 = datadir + offshore_name4
    88 offshore_dem_name5 = datadir + offshore_name5
    89 offshore_dem_name6 = datadir + offshore_name6
    90 offshore_dem_name7 = datadir + offshore_name7
    91 offshore_dem_name8 = datadir + offshore_name8
    92 offshore_dem_name9 = datadir + offshore_name9
    93 offshore_dem_name10 = datadir + offshore_name10
    94 offshore_dem_name11 = datadir + offshore_name11
    95 offshore_dem_name12 = datadir + offshore_name12
    96 offshore_dem_name13 = datadir + offshore_name13
    97 offshore_dem_name14 = datadir + offshore_name14
    98 offshore_dem_name15 = datadir + offshore_name15
    99 offshore_dem_name16 = datadir + offshore_name16
    100 offshore_dem_name17 = datadir + offshore_name17
    101 offshore_dem_name18 = datadir + offshore_name18
    102 offshore_dem_name19 = datadir + offshore_name19
    103 offshore_dem_name20 = datadir + offshore_name20
    104 offshore_dem_name21 = datadir + offshore_name21
    105 offshore_dem_name22 = datadir + offshore_name22
     74anuga_dir = home+state+sep+scenario+sep+'anuga'+sep
    10675
    107 coast_dem_name = datadir + coast_name
     76topographies_in_dir = home+state+sep+scenario+sep+'elevation_final'+sep+'points'+sep
     77topographies_dir = anuga_dir+'topographies'+sep
    10878
    109 combined_dem_name   = datadir + 'exmouth_combined_elevation'
     79# input topo file location
     80onshore_in_dir_name = topographies_in_dir + onshore_name
     81onshore1_in_dir_name = topographies_in_dir + onshore_name1
     82coast_in_dir_name = topographies_in_dir + coast_name
     83offshore_in_dir_name = topographies_in_dir + offshore_name
     84#offshore_in_dir_name1 = topographies_in_dir + offshore_name1
     85#offshore_in_dir_name2 = topographies_in_dir + offshore_name2
    11086
    111 #buildings_filename = gauges_dir + 'exmouth_res_nexis.csv'
     87onshore_dir_name = topographies_dir + onshore_name
     88onshore1_dir_name = topographies_dir + onshore_name1
     89coast_dir_name = topographies_dir + coast_name
     90offshore_dir_name = topographies_dir + offshore_name
     91#offshore_dir_name1 = topographies_dir + offshore_name1
     92#offshore_dir_name2 = topographies_dir + offshore_name2
     93
     94#final topo files
     95combined_dir_name = topographies_dir + combined_name
     96#combined_dir_name_unclipped1 = topographies_dir + combined_name_unclipped1
     97#combined_dir_name1 = topographies_dir + combined_name1
     98combined_small_name_dir = topographies_dir + combined_small_name
     99
     100meshes_dir = anuga_dir+'meshes'+sep
     101meshes_dir_name = meshes_dir + scenario_name
     102
     103polygons_dir = anuga_dir+'polygons'+sep
     104tide_dir = anuga_dir+'tide_data'+sep
     105
     106#boundaries locations
     107boundaries_name = 'exmouth_3854_17042007'
     108
     109boundaries_in_dir = anuga_dir+'boundaries'+sep+'urs'+sep+'1_10000'+sep
     110boundaries_in_dir_name = boundaries_in_dir + boundaries_name
     111boundaries_dir = anuga_dir+'boundaries'+sep
     112boundaries_dir_name = boundaries_dir + boundaries_name
     113
     114#output locations
     115output_dir = anuga_dir+'outputs'+sep
     116output_build_time_dir = output_dir+build_time+sep
     117output_run_time_dir = output_dir +run_time+dir_comment+sep
     118output_run_time_dir_name = output_run_time_dir + scenario_name  #Used by post processing
     119
     120#gauges
     121gauge_name = 'exmouth_gauges.csv'
     122gauges_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'gauges'+sep
     123gauges_dir_name = gauges_dir + gauge_name
     124
     125community_filename = gauges_dir + 'CHINS_v2.csv'
     126community_exmouth = gauges_dir + 'community_exmouth.csv'
     127
    112128buildings_filename_damage_out = 'exmouth_res_nexis_modified.csv'
    113 
    114 community_filename = gaugedir + 'CHINS_v2.csv'
    115 community_exmouth = gaugedir + 'community_exmouth.csv'
    116129
    117130###############################
    118131# Domain definitions
    119132###############################
     133from anuga.utilities.polygon import read_polygon, plot_polygons, polygon_area, is_inside_polygon
    120134
    121 # bounding box for clipping MOST/URS output (much bigger than study area)
    122 ##south = degminsec2decimal_degrees(-19,0,0)
    123 ##north = degminsec2decimal_degrees(-17,15,0)
    124 ##west  = degminsec2decimal_degrees(120,0,0)
    125 ##east  = degminsec2decimal_degrees(122,30,0)
    126 ##
    127 ##d0 = [south, west]
    128 ##d1 = [south, east]
    129 ##d2 = [north, east]
    130 ##d3 = [north, west]
    131 ##poly_bc, zone = convert_points_from_latlon_to_utm([d0, d1, d2, d3])
    132 ##refzone = zone
     135poly_all = read_polygon(polygons_dir+'extent_zone50a.csv')
    133136
    134 # bounding polygon for study area
    135 poly_all = read_polygon(polygondir+'extent_zone50a.csv')
    136 
    137 # plot bounding polygon and make sure BC info surrounds it
    138 #plot_polygons([poly_all, poly_bc],'boundingpoly',verbose=False)
    139137print 'Area of bounding polygon', polygon_area(poly_all)/1000000.0
    140138
     139res_poly_all = 150000*res_factor
     140
     141###############################
     142# Interior region definitions
     143###############################
     144'''
     145poly_0 = read_polygon(polygons_dir+'neg20_coast_contour_pts.csv')
     146res_0 = 20000*res_factor
     147
     148poly_1 = read_polygon(polygons_dir+'exmouth_north_coast_inside_extent.csv')
     149res_1 = 5000*res_factor
     150
     151poly_2 = read_polygon(polygons_dir+'exmouth_south_coast_inside_extent.csv')
     152res_2 = 5000*res_factor
     153
     154poly_3 = read_polygon(polygons_dir+'exmouth_town_pts.csv')
     155res_3 = 2000*res_factor
     156
     157poly_4 = read_polygon(polygons_dir+'exmouth_inner_town_pts.csv')
     158res_4 = 500*res_factor
     159
     160interior_regions = [[poly_0,res_0],[poly_1,res_1],[poly_2,res_2]
     161                     ,[poly_3,res_3],[poly_4,res_4]]
     162
     163trigs_min = number_mesh_triangles(interior_regions, poly_all, res_poly_all)
     164
     165print 'min number triangles', trigs_min
     166
     167poly_mainland = read_polygon(polygons_dir+'Initial_Condition.csv')
     168'''
    141169###################################################################
    142170# Clipping regions for export to asc and regions for clipping data
     
    148176northingmin = 7983427.73
    149177northingmax = 8032834.52
    150 
    151 ###############################
    152 # Interior region definitions
    153 ###############################
    154 '''
    155 # exmouth digitized polygons
    156 poly_exmouth1 = read_polygon(polygondir+'exmouth_Local_Polygon_update.csv')
    157 poly_exmouth2 = read_polygon(polygondir+'exmouth_Close2_update.csv')
    158 poly_exmouth3 = read_polygon(polygondir+'exmouth_Coast_update.csv')
    159 #poly_exmouth4 = read_polygon(polygondir+'Cable_Beach_revised.csv')
    160 
    161 plot_polygons([poly_all,poly_exmouth1,poly_exmouth2,poly_exmouth3],'boundingpoly2',verbose=False)
    162 print 'Area of local polygon', polygon_area(poly_exmouth1)/1000000.0
    163 print 'Area of close polygon', polygon_area(poly_exmouth2)/1000000.0
    164 print 'Area of coastal polygon', polygon_area(poly_exmouth3)/1000000.0
    165 #print 'Area of cable beach polygon', polygon_area(poly_exmouth4)/1000000.0
    166 
    167 for i in poly_exmouth3:
    168     v = is_inside_polygon(i,poly_exmouth1, verbose=False)
    169     if v == False: print v
    170 
    171 def number_mesh_triangles(interior_regions, bounding_poly, remainder_res):
    172     from anuga.utilities.polygon import polygon_area
    173    
    174     # TO DO check if any of the regions fall inside one another
    175     no_triangles = 0.0
    176     area = polygon_area(bounding_poly)
    177     for i,j in interior_regions:
    178         this_area = polygon_area(i)
    179         no_triangles += this_area/j
    180         area -= this_area
    181         print j, this_area/1000000., area/1000000.
    182     no_triangles += area/remainder_res
    183     return int(no_triangles/0.7)
    184 '''
Note: See TracChangeset for help on using the changeset viewer.