Changeset 2420

Timestamp:

Feb 16, 2006, 7:20:00 PM (18 years ago)

Author:

sexton

Message:

updated scripts for Sydney scenario to include new interior regions (refine existing ones and add two high resolution areas)

Location:

production/sydney_2006

Files:

• project.py (modified) (4 diffs)
• run_sydney_smf.py (modified) (2 diffs)

production/sydney_2006/project.py

@@ -34 +34 @@
 nmaxviz = 6283000
 
-basename = 'slump'
+basename = 'slump_TEST'
+basename2 = 'slump_newfriction'
 
 if sys.platform == 'win32':
@@ -54 +55 @@
 combineddemname = datadir + 'sydneytopo'
 outputname = outputdir + basename  #Used by post processing
+outputname2 = outputdir + basename2  #Used by post processing
 
 #csv file of coastline 50m epsilon belt
-manly_polygonname = polygondir + 'manly_polygon_UTM56_coarse' 
-manly_polygon = read_polygon(manly_polygonname + '.csv')
+#manly_polygonname = polygondir + 'manly_polygon_UTM56_coarse' 
+#manly_polygon = read_polygon(manly_polygonname + '.csv')
 #print manly_polygon
 
 gauge_filename = outputdir + 'sydney_gauges.xya'
 gauge_outname = outputdir + 'gauges_max_output.xya'
+polygonptsfile = outputdir + 'poly'
 
 #Georeferencing
@@ -106 +109 @@
 #diffpolygonall = [dp0, dp1, dp2, dp3, dp4, dp8]
 
-
-
-#Interior regions - the Harbour - take 2
+# testing new interior regions 15 Feb 06
+# these worked OK it seemed, no warnings and resulting mesh looked fine.
+pp0 = [343965, 6273229]
+pp1 = [342984, 6270664]
+pp2 = [343950, 6270005]
+pp3 = [343853, 6270399]
+pp4 = [343383, 6270925]
+pp5 = [343756, 6271926]
+pp6 = [344037, 6272401]
+pp7 = [344405, 6272503]
+pp8 = [344226, 6273244]
+
+poly1 = [pp0, pp1, pp2, pp3, pp4, pp5, pp6, pp7, pp8]
+
+qp0 = [343494.0, 6270650.0]
+qp1 = [343337.0, 6270303.0]
+qp2 = [343466.0, 6270228.0]
+qp3 = [343139.0, 6269901.0]
+qp4 = [342472.0, 6268573.0]
+qp5 = [342111.0, 6267115.0]
+qp6 = [342479.0, 6266121.0]
+qp7 = [342860.0, 6266386.0]
+qp8 = [342635.0, 6267438.0]
+qp9 = [343105.0, 6269070.0]
+qp10 = [343548.0, 6269567.0]
+qp11 = [343487.0, 6269928.0]
+qp12 = [343991.0, 6270269.0]
+
+poly2 = [qp0, qp1, qp2, qp3, qp4, qp5, qp6, qp7, qp8, qp9, qp10, qp11, qp12]
+#poly2 = [qp5, qp6, qp7, qp8, qp9, qp10]
+
+# didn't like this one - poly2 from Ingo
+# warning generated about vertex blah doesn't belong to an element
+# is there a maximum number of vertices per polygon?
+p0 = [343488.9455,      6270644.956]
+p1 = [343681.208,       6270500.759]
+p2 = [343997.0678,      6270253.564]
+p3 = [343482.0789,      6269930.838]
+p4 = [343550.7441,      6269546.313]
+p5 = [343111.287,       6269099.989]
+p6 = [342630.6307,      6267534.423]
+p7 = [342870.9588,      6266380.848]
+p8 = [342479.5673,      6266126.787]
+p9 = [342174.4759,      6266317.303]
+p10 = [342248.668,      6266727.125]
+p11 = [342101.9088,     6267108.699]
+p12 = [340986.2201,     6267454.634]
+p13 = [340663.6936,     6266113.602]
+p14 = [340120.4909,     6265892.926]
+p15 = [339560.3132,     6266317.303]
+p16 = [339237.7867,     6266928.406]
+p17 = [339000.1355,     6267912.961]
+p18 = [339169.8863,     6268642.889]
+p19 = [340001.6654,     6269016.341]
+p20 = [340392.0922,     6269423.743]
+p21 = [340867.3945,     6268863.566]
+p22 = [341614.2982,     6268659.865]
+p23 = [341495.4726,     6269220.042]
+p24 = [341750.0988,     6269287.943]
+p25 = [341902.8746,     6269882.07]
+p26 = [342700.7034,     6269780.22]
+p27 = [343125.0805,     6269746.27]
+p28 = [343125.0805,     6270221.572]
+p29 = [343057.1801,     6270527.124]
+p30 = [343488.9455,     6270644.956]
+
+testpoly = [p0,p1,p2,p3,p4,p5,p6,p7,p8,p9,p10,p11,p12,p13,p14,p15,p16,p17,p18,p19,p20,p21,p22,p23,p24,p25,p26,p27,p28,p29,p30]
+
+# test from Jane drawing
+# north
+np1 = [318200, 6253000]
+np2 = [327000, 6250500]
+np3 = [333000, 6249000]
+np4 = [342000, 6249000]
+np5 = [342000, 6255000]
+np6 = [343000, 6256000]
+np7 = [344000, 6258000]
+np8 = [343000, 6260000]
+np9 = [343000, 6264000]
+np10 = [345000, 6265000]
+np11 = [343000, 6266000]
+np12 = [344000, 6269000]
+np13 = [344000, 6272000]
+np14 = [342000, 6272000]
+np15 = [339000, 6269000]
+np16 = [339000, 6264000]
+#np17 = [332000, 6264000]
+np17 = [332500, 6262000]
+np18 = [334000, 6254000]
+np19 = [329000, 6257000]
+np20 = [330000, 6259000]
+#np21 = [327000, 6262000]
+np21 = [327000, 6259000]
+np22 = [327000, 6257000]
+np23 = [318200, 6257000]
+np24 = [335000, 6261000]
+np25 = [336000, 6262000]
+np26 = [338000, 6262000]
+np27 = [340000, 6264000]
+np28 = [334000, 6250000]
+np29 = [336000, 6250000]
+np30 = [341500, 6250000]
+np31 = [327000, 6254000]
+np32 = [323000, 6257000]
+np33 = [335000, 6255000]
+np34 = [337000, 6256000]
+np35 = [338000, 6255000]
+np36 = [319000, 6254000]
+np37 = [322000, 6252000]
+np38 = [324000, 6253000]
+np39 = [325000, 6251000]
+
+#testpoly1 = [np1, \
+#             np36, np37, np38, np39, \
+#             np2, np3, \
+#             np28, np29, np30, \
+#             np5, np6, np7, np8, np9, np10, \
+#             np11, np12, np13, np14, np15, \
+#             np27, np26, np25, np24, \
+#             np17, \
+#             np33, np34, np35, \
+#             np18, np19, np20, \
+#             np21, np22, \
+#             np31, np32, \
+#             np23]
+newpoly1 = [np2, np3, \
+            np28, np29, np30, \
+            np5, np6, np7, np8, np9, np10, \
+            np11, np12, np13, np14, np15, \
+            np27, np26, np25, np24, np17, \
+            np18, np19, np20, np21, np31]
+
+# south
+sp1 = [ 328000, 6231000]
+sp2 = [335000, 6231000]
+sp3 = [338000, 6235000]
+sp4 = [340000, 6240000]
+sp5 = [340000, 6244000]
+sp6 = [342000, 6248800]
+sp7 = [340000, 6248800]
+sp8 = [338000, 6244000]
+sp9 = [338000, 6237000]
+sp10 = [334000, 6242000]
+sp11 = [331000, 6245000]
+sp12 = [326000, 6247000]
+sp13 = [325000, 6246000]
+sp14 = [329000, 6243000]
+sp15 = [328000, 6237000]
+sp16 = [337000, 6236000]
+sp17 = [330000, 6237000]
+sp18 = [330000, 6239000]
+sp19 = [332000, 6239000]
+sp20 = [334000, 6240000]
+sp21 = [334000, 6238000]
+sp22 = [337000, 6236500]
+sp23 = [339000, 6236000]
+
+newpoly2 = [sp1, sp2, sp3, \
+             sp16, sp17, sp18, sp19, sp20, sp21, sp22, sp23, \
+             sp4, sp5, sp6, sp7, sp8, sp9, sp10, \
+             sp11, sp12, sp13, sp14, sp15]
+
+m1 = [340000, 6256000]
+m2 = [342800, 6256000]
+m3 = [342800, 6261000]
+m4 = [340000, 6260000]
+    
+finepolymanly = [m1, m2, m3, m4]
+
+q1 = [333000, 6250000]
+q2 = [340000, 6250000]
+q3 = [340000, 6254000]
+q4 = [333000, 6255000]
+    
+finepolyquay = [q1, q2, q3, q4]
+
 #Interior regions - the Harbour
 harbour_1x = degminsec2decimal_degrees(-33,51,0)
@@ -162 +338 @@
 
 #Interior region - Botany Bay
-
-#Interior region - Botany Bay - take 2
 bb_1x = degminsec2decimal_degrees(-34,3,0)
 bb_1y = degminsec2decimal_degrees(151,2,30)

production/sydney_2006/run_sydney_smf.py

@@ -68 +68 @@
 # boundary and interior regions defined in project.py along with
 # resolutions (maximal area of per triangle) for each polygon
-#-------------------------------------------------------------------------------                                 
+#-------------------------------------------------------------------------------
+
+#def get_polygon_from_file(filename):
+#    fid = open(filename)
+#    lines = fid.readlines()
+#    fid.close()
+    
+#    polygon = []
+#    for line in lines[1:]:
+#       fields = line.split(',')
+#        x = float(fields[3])
+#        y = float(fields[4])
+#        polygon.append([x, y])
+        
+#    return interior_polygon
 
 from pmesh.create_mesh import create_mesh_from_regions
 
-interior_res = 5000 
-interior_regions = [[project.harbour_polygon_2, interior_res],
-                    [project.botanybay_polygon_2, interior_res]]
+num_polygons = 1
+filelist = []
+fileext = '.csv'
+filename = project.polygonptsfile
+interior_res = 50000
 
+for p in range(1, num_polygons+1):
+    thefilename = filename + str(p) + fileext
+    interior_polygon = get_polygon_from_file(thefilename)
+    interior_regions.append([interior_polygon, interior_res])
+    
+# test
+#interior_regions = [[project.poly1, interior_res],
+#                    [project.poly2, interior_res]]
+# original
+#interior_regions = [[project.harbour_polygon_2, interior_res],
+#                    [project.botanybay_polygon_2, interior_res]]
 
 #FIXME: Fix caching of this one once the mesh_interface is ready
@@ -130 +157 @@
 
 domain.set_quantity('stage', tsunami_source)
-domain.set_quantity('friction', 0)
+domain.set_quantity('friction', 0.03) # supplied by Benfield
 domain.set_quantity('elevation',
                     filename = project.combineddemname + '.pts',