Changeset 2615 for production/onslow_2006/run_onslow_new.py

Timestamp:

Mar 28, 2006, 11:00:20 AM (19 years ago)

Author:

nick

Message:

update onslow

File:

• production/onslow_2006/run_onslow_new.py (moved) (moved from production/onslow_2006/run_onslow.py) (14 diffs)

production/onslow_2006/run_onslow_new.py

@@ -21 +21 @@
 # Related major packages
 from pyvolution.shallow_water import Domain, Reflective_boundary, \
-                            Dirichlet_boundary, Time_boundary
+                            Dirichlet_boundary, Time_boundary, File_boundary
 from pyvolution.data_manager import convert_dem_from_ascii2netcdf, dem2pts
 from pyvolution.combine_pts import combine_rectangular_points_files 
 from pyvolution.pmesh2domain import pmesh_to_domain_instance
+from geospatial_data import add_points_files
 
 # Application specific imports
@@ -30 +31 @@
 from smf import slump_tsunami  # Function for submarine mudslide 
 
+from shutil import copy
+from os import mkdir, access, F_OK
+
+from geospatial_data import *
 
 #-------------------------------------------------------------------------------
@@ -41 +46 @@
 # filenames
 coarsedemname = project.coarsedemname
-'''
-#finedemname = project.finedemname
-'''
+
+onshore_dem_name = project.onshore_dem_name
+
+offshore_points = project.offshore_dem_name
+
 meshname = project.meshname+'.msh'
 
+source_dir = project.boundarydir
+
+# creates copy of code in output dir
+if access(project.outputdir,F_OK) == 0 :
+    mkdir (project.outputdir)
+copy (project.codedirname, project.outputdir + project.codename)
+copy (project.codedir + 'run_onslow.py', project.outputdir + 'run_onslow.py')
+
+
+'''
 # coarse data
 convert_dem_from_ascii2netcdf(coarsedemname, use_cache=True, verbose=True)
 dem2pts(coarsedemname, use_cache=True, verbose=True)
-'''
+
+
 # fine data (clipping the points file to smaller area)
-convert_dem_from_ascii2netcdf(finedemname, use_cache=True, verbose=True)
-dem2pts(finedemname,
+convert_dem_from_ascii2netcdf(onshore_dem_name, use_cache=True, verbose=True)
+dem2pts(onshore_dem_name,
         easting_min=project.eastingmin,
         easting_max=project.eastingmax,
@@ -60 +78 @@
         verbose=True)
 
-
-# combining the coarse and fine data
-combine_rectangular_points_files(project.finedemname + '.pts',
-                                 project.coarsedemname + '.pts',
-                                 project.combineddemname + '.pts')
-'''
-
+        '''
+print'create G1'
+G1 = Geospatial_data(file_name = project.offshore_dem_name + '.xya')
+
+print'create G2'
+G2 = Geospatial_data(file_name = project.onshore_dem_name + '.pts')
+
+print'add G1+G2'
+G = G1 + G2
+
+print'export G'
+G.new_export_points_file(project.combined_dem_name + '.pts')
+
+'''
+add_points_files(
+                  project.offshore_dem_name + '.xya',
+                  project.onshore_dem_name + '.pts',
+                  project.combined_dem_name + '.pts')
+'''
 #-------------------------------------------------------------------------------                                 
 # Create the triangular mesh based on overall clipping polygon with a tagged
@@ -73 +103 @@
 #-------------------------------------------------------------------------------
 
-from pmesh.create_mesh import create_mesh_from_regions
+from pmesh.mesh_interface import create_mesh_from_regions
 
 # original
@@ -80 +110 @@
                     [project.poly_thevenard, interior_res],
                     [project.poly_direction, interior_res]]
-
-#FIXME: Fix caching of this one once the mesh_interface is ready
+                    #[project.testpoly, interior_res]]
+print 'number of interior regions', len(interior_regions)
+
 from caching import cache
 _ = cache(create_mesh_from_regions,
@@ -88 +119 @@
                              'left': [2], 'bottom': [3],
                              'bottomright': [4], 'topright': [5]},
-           'resolution': 100000,
+           'maximum_triangle_area': 100000,
            'filename': meshname,           
-           'interior_regions': interior_regions,
-           'UTM': True,
-           'refzone': project.refzone},
+           'interior_regions': interior_regions},
           verbose = True)
 
@@ -106 +135 @@
 print 'Number of triangles = ', len(domain)
 print 'The extent is ', domain.get_extent()
+print domain.statistics()
 
 domain.set_name(project.basename)
 domain.set_datadir(project.outputdir)
-domain.set_quantities_to_be_stored(['stage', 'xmomentum', 'ymomentum'])
+domain.set_quantities_to_be_stored(['stage'])
 
 
@@ -136 +166 @@
 domain.set_quantity('stage', tide)
 domain.set_quantity('friction', 0.0) 
+print 'hi and file',project.combined_dem_name + '.pts'
 domain.set_quantity('elevation', 
-#                    filename = project.combineddemname + '.pts',
-                    filename = project.coarsedemname + '.pts',
-                    use_cache = True,
-                    verbose = True
+#                    0.
+#                    filename = project.onshore_dem_name + '.pts',
+                    filename = project.combined_dem_name + '.pts',
+#                    filename = project.offshore_dem_name + '.pts',
+                    use_cache = False,
+                    verbose = True,
+                    alpha = 0.1
                     )
-
+print 'hi1'
 
 #-------------------------------------------------------------------------------                                 
@@ -156 +190 @@
 
 cache(ferret2sww,
-      (source_dir+project.boundary_basename,
-       project.boundary_basename), 
+      (source_dir + project.boundary_basename,
+       source_dir + project.boundary_basename), 
       {'verbose': True,
-       'minlat': south-1,
-       'maxlat': north+1,
-       'minlon': west-1,
-       'maxlon': east+1,
-       'origin': project.mesh_origin,
+# note didn't work with the below
+#       'minlat': south - 1,
+#       'maxlat': north + 1,
+#       'minlon': west - 1,
+#       'maxlon': east + 1,
+       'minlat': south,
+       'maxlat': north,
+       'minlon': west,
+       'maxlon': east,
+#       'origin': project.mesh_origin,
+       'origin': domain.geo_reference.get_origin(),
        'mean_stage': tide,
        'zscale': 1,                 #Enhance tsunami
@@ -174 +214 @@
 print 'Available boundary tags', domain.get_boundary_tags()
 
+Bf = File_boundary(source_dir + project.boundary_basename + '.sww', 
+                    domain, verbose = True)
 Br = Reflective_boundary(domain)
 Bd = Dirichlet_boundary([tide,0,0])
@@ -182 +224 @@
                    f=lambda t: [(60<t<480)*6, 0, 0])
 
-domain.set_boundary( {'top': Bw, 'topleft': Br,
+domain.set_boundary( {'top': Bf, 'topleft': Bf,
                              'left': Br, 'bottom': Br,
-                             'bottomright': Br, 'topright': Br} )
+                             'bottomright': Br, 'topright': Bf} )
 
 
@@ -190 +232 @@
 # Evolve system through time
 #-------------------------------------------------------------------------------
-
 import time
 t0 = time.time()
@@ -199 +240 @@
     
 print 'That took %.2f seconds' %(time.time()-t0)
+
+print 'finished'