source: anuga_work/production/test_dampier_2006/project_old.py @ 7614

Last change on this file since 7614 was 4288, checked in by duncan, 17 years ago

A test for the dampier boundary problem.

File size: 15.9 KB
Line 
1"""Common filenames and locations for topographic data, meshes and outputs.
2"""
3
4from os import sep, environ, getenv, getcwd
5from os.path import expanduser
6import sys
7from time import localtime, strftime, gmtime
8from anuga.utilities.polygon import read_polygon, plot_polygons, polygon_area, is_inside_polygon, number_mesh_triangles
9from anuga.coordinate_transforms.redfearn import degminsec2decimal_degrees, convert_from_latlon_to_utm
10from anuga.utilities.system_tools import get_user_name
11
12# file and system info
13#---------------------------------
14codename = 'project.py'
15
16home = getenv('INUNDATIONHOME') #Sandpit's parent dir   
17user = get_user_name()
18
19# INUNDATIONHOME is the inundation directory, not the data directory.
20home += sep +'data'
21
22#time stuff
23time = strftime('%Y%m%d_%H%M%S',localtime()) #gets time for new dir
24gtime = strftime('%Y%m%d_%H%M%S',gmtime()) #gets time for new dir
25build_time = time+'_build'
26run_time = time+'_run'
27temp_time = time+'_temp'
28print 'gtime: ', gtime
29
30tide = 2.4
31
32 #Making assumptions about the location of scenario data
33 state = 'western_australia'
34 scenario_name = 'dampier'
35 scenario = 'dampier_tsunami_scenario_2006'
36
37# # onshore data provided by WA DLI
38# onshore_name = 'DLI_DTED_raster_clipped' # original
39
40# # AHO + DPI data + colin French coastline
41# coast_name = 'coastline_edited_w_DEM'
42# offshore_name = 'clipped_bathy'
43# offshore1_name = 'elev_501'
44# offshore2_name = 'inferrec_e'
45
46# #final topo name
47# combined_name ='dampier_combined_elevation'
48# combined_smaller_name = 'dampier_combined_elevation_small'
49# combined_smallest_name = 'dampier_combined_elevation_smallest'
50
51# topographies_in_dir = home+sep+state+sep+scenario+sep+'elevation_final'+sep+'070112'+sep+'points'+sep
52# topographies_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'topographies'+sep
53# topographies_time_dir = topographies_dir+build_time+sep
54
55# temp_dir = home+sep+state+sep+scenario+sep+'anuga'+sep
56# temp_dir_name = temp_dir + temp_time+sep
57
58# # input topo file location
59# onshore_in_dir_name = topographies_in_dir + onshore_name
60# coast_in_dir_name = topographies_in_dir + coast_name
61# offshore_in_dir_name = topographies_in_dir + offshore_name
62# offshore1_in_dir_name = topographies_in_dir + offshore1_name
63# offshore2_in_dir_name = topographies_in_dir + offshore2_name
64
65# onshore_dir_name = topographies_dir + onshore_name
66# coast_dir_name = topographies_dir + coast_name
67# offshore_dir_name = topographies_dir + offshore_name
68# offshore1_dir_name = topographies_dir + offshore1_name
69# offshore2_dir_name = topographies_dir + offshore2_name
70
71# #final topo files
72# combined_dir_name = topographies_dir + combined_name
73# combined_time_dir_name = topographies_time_dir + combined_name
74# combined_smaller_dir_name = topographies_dir + combined_smaller_name
75# combined_smallest_dir_name = topographies_dir + combined_smallest_name
76# #combined_time_dir_final_name = topographies_time_dir + combined_final_name
77
78# meshes_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'meshes'+sep
79# meshes_dir_name = meshes_dir + scenario_name
80
81# polygons_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'polygons'+sep+'2007polys'+sep
82# tide_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'tide_data'+sep
83
84# boundaries_source = ''
85 boundaries_name = 'o'
86# boundaries_name1 = 'o_new1'
87# boundaries_name2 = 'o_test'
88
89# #boundaries locations
90# #boundaries_in_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'boundaries'+sep
91# boundaries_in_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'boundaries'+sep+'urs'+sep+'1_10000'+sep
92# boundaries_in_dir_name = boundaries_in_dir + boundaries_name
93 boundaries_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'boundaries'+sep
94 boundaries_dir_name = boundaries_dir + boundaries_name
95# boundaries_dir_name1 = boundaries_dir + boundaries_name1
96# boundaries_dir_name2 = boundaries_dir + boundaries_name2
97# boundaries_dir_name3 = boundaries_dir + boundaries_name+'_test_8500_12000'
98# boundaries_dir_name4 = boundaries_dir + boundaries_name+'_8500_12000_no_zone'
99# boundaries_dir_name5 = boundaries_dir + boundaries_name+'_not-gridded'
100
101# #boundaries_time_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'boundaries'+sep+build_time+sep
102# #boundaries_time_dir_name = boundaries_time_dir + boundaries_name  #Used by post processing
103
104# #output locations
105# output_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'outputs'+sep
106# output_build_time_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'outputs'+sep+build_time+sep
107# output_run_time_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'outputs'+sep+run_time+sep
108# output_run_time_dir_name = output_run_time_dir + scenario_name  #Used by post processing
109
110# #gauges
111# gauge_name = 'dampier_gauges_up2.csv' #'dampier.csv'
112# gauge_name_test = 'dampier_gauges_up_test.csv' #'dampier.csv'
113# gauges_dir = home+sep+state+sep+scenario+sep+'anuga'+sep+'gauges'+sep
114# gauges_dir_name = gauges_dir + gauge_name
115# gauges_dir_name_test = gauges_dir + gauge_name_test
116
117
118# #buildings_filename = gauges_dir + 'dampier_res_nexis.csv'
119# buildings_filename_damage_out = 'dampier_res_nexis_modified.csv'
120# ###############################
121# # Domain definitions
122# ###############################
123
124# refzone = 50
125# south_boundary = degminsec2decimal_degrees(-20,58,0)
126# #north_boundary = degminsec2decimal_degrees(-20,13,0)
127# north_boundary = degminsec2decimal_degrees(-19,30,0)
128# west_boundary = degminsec2decimal_degrees(116,15,0)
129# east_boundary = degminsec2decimal_degrees(117,11,0)
130# ##
131# ##p0 = [south, degminsec2decimal_degrees(116,32,0)]
132# ##p1 = [south, west]
133# ##p2 = [degminsec2decimal_degrees(-20,23,0), west]
134# ##p3 = [north, degminsec2decimal_degrees(116,45,0)]
135# ##p4 = [north, degminsec2decimal_degrees(117,0,0)]
136# ##p5 = [p2[0], degminsec2decimal_degrees(117,8,0)]
137# ##p6 = [degminsec2decimal_degrees(-20,30,0), east]
138# ##p7 = [degminsec2decimal_degrees(-20,38,0), east]
139# ##p8 = [south, east]
140# ##
141# ##poly_all, zone = convert_from_latlon_to_utm([p0, p1, p2, p3, p4, p5, p6, p7, p8])
142# ##refzone = zone
143# poly_all = read_polygon(polygons_dir+'extent.csv')
144# print 'Area of bounding polygon', polygon_area(poly_all)/1000000.0
145
146# res_poly_all = 400000
147# #res_poly_all = 500000
148
149# ###############################
150# # Interior region definitions
151# ###############################
152
153# poly_region = read_polygon(polygons_dir+'region.csv')
154# res_region = 50000
155
156# poly_dampier = read_polygon(polygons_dir+'dampier_town.csv')
157# res_dampier = 1000
158# #res_dampier = 500
159# #res_dampier = 5000
160
161# poly_karratha = read_polygon(polygons_dir+'karrathav2.csv')
162# res_karratha = 15000
163
164# poly_karratha_town = read_polygon(polygons_dir+'karratha_townv2.csv')
165# res_karratha_town = 1000
166# #res_karratha_town = 500
167# #res_karratha_town = 5000
168
169# poly_facility = read_polygon(polygons_dir+'facility.csv')
170# res_facility = 1000
171
172# poly_delambre = read_polygon(polygons_dir+'delambre.csv')
173# res_delambre = 5000
174
175# poly_coast = read_polygon(polygons_dir+'coastpoly.csv')
176# res_coast = 1000
177# #res_coast = 10000
178
179# poly_NWislands = read_polygon(polygons_dir+'nw_islands_area.csv')
180# res_NWislands = 50000
181
182# poly_island0 = read_polygon(polygons_dir+'island0.csv')
183# res_island0 = res_poly_all
184
185# poly_island1 = read_polygon(polygons_dir+'island1.csv')
186# res_island0 = res_poly_all
187
188# poly_island2 = read_polygon(polygons_dir+'island2.csv')
189# res_island0 = res_poly_all
190
191# poly_island3 = read_polygon(polygons_dir+'island3.csv')
192# res_island0 = res_poly_all
193
194# res_islands = 5000
195# #res_islands = 15000
196
197# poly_ref_nw4 = read_polygon(polygons_dir+'ref_nw4.csv')
198# res_ref_nw4 = res_islands
199
200# poly_island4 = read_polygon(polygons_dir+'island4.csv')
201# res_island0 = res_poly_all
202
203# poly_ref_nw5 = read_polygon(polygons_dir+'ref_nw5.csv')
204# res_ref_nw5 = res_islands
205
206# poly_island5 = read_polygon(polygons_dir+'island5.csv')
207# res_island0 = res_poly_all
208
209# poly_ref_nw6 = read_polygon(polygons_dir+'ref_nw6.csv')
210# res_ref_nw6 = res_islands
211
212# poly_island6 = read_polygon(polygons_dir+'island6.csv')
213# res_island0 = res_poly_all
214
215# poly_ref_nw7 = read_polygon(polygons_dir+'ref_nw7.csv')
216# res_ref_nw7 = res_islands
217
218# poly_island7 = read_polygon(polygons_dir+'island7.csv')
219# res_island0 = res_poly_all
220
221# poly_ref_nw8 = read_polygon(polygons_dir+'ref_nw8.csv')
222# res_ref_nw8 = res_islands
223
224# poly_island8 = read_polygon(polygons_dir+'island8.csv')
225# res_island0 = res_poly_all
226
227
228# ##plot_polygons([poly_dampier,poly_karratha,poly_karratha_town,poly_delambre,
229# ##                poly_coast,poly_NWislands,poly_island0,poly_island1,poly_island2,
230# ##                poly_island3,poly_island4,poly_island5,poly_island6,
231# ##                poly_island7,poly_island8,poly_ref_nw4,poly_ref_nw5,
232# ##                poly_ref_nw6,poly_ref_nw7,poly_ref_nw8,poly_all],'poly_pic')
233
234# interior_regions = [[poly_dampier,res_dampier],
235#                     [poly_karratha,res_karratha],[poly_karratha_town,res_karratha_town],
236#                     [poly_delambre,res_delambre],[poly_coast,res_coast],
237#                     [poly_facility,res_facility],
238#                     #[poly_NWislands,res_NWislands],
239#                     [poly_island0,res_island0],[poly_island1,res_island0],
240#                     [poly_island2,res_island0],[poly_island3,res_island0],
241#                     [poly_island4,res_island0],[poly_island5,res_island0],
242#                     [poly_island6,res_island0],[poly_island7,res_island0],
243#                     [poly_island8,res_island0],[poly_ref_nw4,res_ref_nw4],
244#                     [poly_ref_nw5,res_ref_nw5],[poly_ref_nw6,res_ref_nw6],
245#                     [poly_ref_nw7,res_ref_nw7],[poly_ref_nw8,res_ref_nw8]]
246
247# interior_regions_test = [[poly_dampier,res_dampier],
248#                     [poly_karratha,res_karratha],[poly_karratha_town,res_karratha_town]]
249
250# trigs_min = number_mesh_triangles(interior_regions_test, poly_all, res_poly_all)
251
252# print 'min number triangles', trigs_min
253
254# ###################################################################
255# # Clipping regions for export to asc and regions for clipping data
256# ###################################################################
257
258# poly_bathy = read_polygon(polygons_dir+'mainland_only.csv')
259
260# # exporting asc grid - Dampier
261# e_min_area = 474000
262# e_max_area = 480000
263# n_min_area = 7719000
264# n_max_area = 7725000
265
266# # exporting asc grid - Karratha
267# ##e_min_area =
268# ##e_max_area =
269# ##n_min_area =
270# ##n_max_area =
271
272# """
273# # used in the CIPMA 2006 scenario
274# # CIPMA point of interest
275# cipma_latitude = -20.588456
276# cipma_longitude = 116.771527
277
278# k0 = [cipma_latitude-0.02, cipma_longitude-0.02]
279# k1 = [cipma_latitude-0.02, cipma_longitude+0.02]
280# k2 = [cipma_latitude+0.02, cipma_longitude+0.02]
281# k3 = [cipma_latitude+0.02, cipma_longitude-0.02]
282
283# cipma_polygon, zone = convert_from_latlon_to_utm([k0, k1, k2, k3])
284# assert zone == refzone
285
286# poly_facility = read_polygon(polygons_dir+'facility.csv')
287# poly_pipeline = read_polygon(polygons_dir+'pipeline2.csv')
288# poly_interior = read_polygon(polygons_dir+'interior.csv')
289# poly_coast = read_polygon(polygons_dir+'coast_final.csv')
290# clip_poly_e = read_polygon(polygons_dir+'gap_e.csv')
291# clip_poly_nw = read_polygon(polygons_dir+'gap_nw.csv')
292# clip_poly_mid_w = read_polygon(polygons_dir+'gap_mid_w.cvs')
293# clip_poly_mid_e = read_polygon(polygons_dir+'gap_mid_e.cvs')
294
295# # exporting asc grid
296# e_min_area = 474000
297# e_max_area = 480000
298# n_min_area = 7719000
299# n_max_area = 7725000
300
301# # used in the original 2005 scenario
302# #Interior regions
303# karratha_south = degminsec2decimal_degrees(-20,44,0)
304# karratha_north = degminsec2decimal_degrees(-20,42,0)
305# karratha_west = degminsec2decimal_degrees(116,48,0)
306# karratha_east = degminsec2decimal_degrees(116,53,30)
307
308# k0 = [karratha_south, karratha_west]
309# k1 = [karratha_south, karratha_east]
310# k2 = [karratha_north, karratha_east]
311# k3 = [karratha_north, karratha_west]   
312
313# karratha_polygon, zone = convert_from_latlon_to_utm([k0, k1, k2, k3])
314# assert zone == refzone
315
316# #Interior regions
317# dampier_south = degminsec2decimal_degrees(-20,40,0)
318# dampier_north = degminsec2decimal_degrees(-20,38,10)
319# dampier_west = degminsec2decimal_degrees(116,43,0)
320# dampier_east = degminsec2decimal_degrees(116,45,0)
321
322# d0 = [dampier_south, dampier_west]
323# d1 = [dampier_south, dampier_east]
324# d2 = [dampier_north, dampier_east]
325# d3 = [dampier_north, dampier_west]   
326
327# dampier_polygon, zone = convert_from_latlon_to_utm([d0, d1, d2, d3])
328# assert zone == refzone
329
330# #Interior regions
331# refinery_south = degminsec2decimal_degrees(-20,37,50)
332# refinery_north = degminsec2decimal_degrees(-20,36,0)
333# refinery_west = degminsec2decimal_degrees(116,44,0)
334# refinery_east = degminsec2decimal_degrees(116,46,10)
335
336# d0 = [refinery_south, refinery_west]
337# d1 = [refinery_south, refinery_east]
338# d2 = [refinery_north, refinery_east]
339# d3 = [refinery_north, refinery_west]   
340
341# refinery_polygon, zone = convert_from_latlon_to_utm([d0, d1, d2, d3])
342# assert zone == refzone
343
344# #Interior region around 468899, 7715177:
345# #lat (-20, 39, 44.93753), lon (116, 42, 5.09106)
346
347# point_south = degminsec2decimal_degrees(-20,39,46)
348# point_north = degminsec2decimal_degrees(-20,39,42)
349# point_west = degminsec2decimal_degrees(116,42,0)
350# point_east = degminsec2decimal_degrees(116,42,10)
351
352# d0 = [point_south, point_west]
353# d1 = [point_south, point_east]
354# d2 = [point_north, point_east]
355# d3 = [point_north, point_west]   
356
357# point_polygon, zone = convert_from_latlon_to_utm([d0, d1, d2, d3])
358# assert zone == refzone
359
360# #Neils areas around interesting points
361# neil1_point1 = [degminsec2decimal_degrees(-20,35,34),
362#                 degminsec2decimal_degrees(116,45,18)]
363# neil1_point2 = [degminsec2decimal_degrees(-20,36,15),
364#                 degminsec2decimal_degrees(116,46,18)]
365# neil1_point3 = [degminsec2decimal_degrees(-20,35,9),
366#                 degminsec2decimal_degrees(116,47,17)]
367# neil1_point4 = [degminsec2decimal_degrees(-20,34,26),
368#                 degminsec2decimal_degrees(116,46,17)]
369
370# neil1_polygon, zone = convert_from_latlon_to_utm([neil1_point1,
371#                                                   neil1_point2,
372#                                                   neil1_point3,
373#                                                   neil1_point4])
374# assert zone == refzone
375
376# neil2_point1 = [degminsec2decimal_degrees(-20,39,36),
377#                 degminsec2decimal_degrees(116,41,33)]
378# neil2_point2 = [degminsec2decimal_degrees(-20,40,10),
379#                 degminsec2decimal_degrees(116,42,13)]
380# neil2_point3 = [degminsec2decimal_degrees(-20,38,39),
381#                 degminsec2decimal_degrees(116,43,49)]
382# neil2_point4 = [degminsec2decimal_degrees(-20,38,5),
383#                 degminsec2decimal_degrees(116,43,9)]
384
385# neil2_polygon, zone = convert_from_latlon_to_utm([neil2_point1,
386#                                                   neil2_point2,
387#                                                   neil2_point3,
388#                                                   neil2_point4])
389# assert zone == refzone
390
391# #Withnell bay
392# wb_point1 = [degminsec2decimal_degrees(-20,35,34),
393#                 degminsec2decimal_degrees(116,45,18)]
394# wb_point2 = [degminsec2decimal_degrees(-20,36,15),
395#                 degminsec2decimal_degrees(116,46,18)]
396# wb_point3 = [degminsec2decimal_degrees(-20,35,9),
397#                 degminsec2decimal_degrees(116,47,17)]
398# wb_point4 = [degminsec2decimal_degrees(-20,34,26),
399#                 degminsec2decimal_degrees(116,46,17)]
400
401# wb_polygon, zone = convert_from_latlon_to_utm([wb_point1, wb_point2,
402#                                                wb_point3, wb_point4])
403# assert zone == refzone
404
405# #Larger Withnell bay
406# lwb_point1 = [degminsec2decimal_degrees(-20,35,59),
407#                 degminsec2decimal_degrees(116,42,00)]
408# lwb_point2 = [degminsec2decimal_degrees(-20,36,50),
409#                 degminsec2decimal_degrees(116,46,50)]
410# lwb_point3 = [degminsec2decimal_degrees(-20,34,00),
411#                 degminsec2decimal_degrees(116,47,39)]
412# lwb_point4 = [degminsec2decimal_degrees(-20,33,00),
413#                 degminsec2decimal_degrees(116,42,50)]
414
415# lwb_polygon, zone = convert_from_latlon_to_utm([lwb_point1, lwb_point2,
416#                                                 lwb_point3, lwb_point4])
417                                                     
418# assert zone == refzone
419# """
Note: See TracBrowser for help on using the repository browser.