source: trunk/anuga_work/development/convergence_okushiri_2008/export_results.py @ 7884

Last change on this file since 7884 was 5720, checked in by Leharne, 16 years ago

Updates to okushiri convergence study

File size: 5.7 KB
Line 
1import project_truescale, os
2import sys
3import time
4from anuga.shallow_water.data_manager import convert_dem_from_ascii2netcdf, dem2pts
5from anuga.shallow_water.data_manager import sww2dem
6from anuga.geospatial_data.geospatial_data import *
7from os import sep
8
9res_dir1 = '20080725_051221_run_octuple_contour_polygons_lfountai'
10res_dir2 = '20080715_023223_run_quadruple_contour_polygons_lfountai'
11res_dir3 = '20080708_064841_run_double_contour_polygons_lfountai'
12res_dir4 = '20080708_064714_run_original_contour_polygons_lfountai'
13res_dir5 = '20080708_064815_run_half_contour_polygons_lfountai'
14res_dir6 = '20080714_231759_run_quarter_contour_polygons_lfountai'
15res_dir7 = '20080714_234934_run_eighth_contour_polygons_lfountai'
16res_dir8 = '20080714_235634_run_sixteenth_contour_polygons_lfountai'
17res_dir9 = '20080715_015727_run_1-32_contour_polygons_lfountai'
18res_dir10 = '20080721_015521_run_1-64_contour_polygons_lfountai'
19res_dirs = [res_dir1, #res_dir2, res_dir3, res_dir4, res_dir5, res_dir6, res_dir7, res_dir8,res_dir9,
20            res_dir10]
21#timesteps = [200, 250, 300, 350, 400, 450]
22timesteps = [327] #this is time of maximum inundation
23
24
25for res_dir in res_dirs:
26   
27    for timestep in timesteps:
28               
29        directory = project_truescale.output_dir+res_dir+sep
30        name = directory+'okushiri_truescale'
31
32        is_parallel = False
33        #is_parallel = True
34
35        if is_parallel == True: nodes = 10
36        print 'output directory:', directory
37
38        area = ['Deep', 'Mid', 'Shallow']
39           
40        for which_area in area:
41            if which_area == 'Deep':
42                cellsize = 50
43                easting_min = project_truescale.xminDeep
44                easting_max = project_truescale.xmaxDeep
45                northing_min = project_truescale.yminDeep
46                northing_max = project_truescale.ymaxDeep
47
48            if which_area == 'Mid':
49                cellsize = 25
50                easting_min = project_truescale.xminMid
51                easting_max = project_truescale.xmaxMid
52                northing_min = project_truescale.yminMid
53                northing_max = project_truescale.ymaxMid
54
55            if which_area == 'Shallow':
56                cellsize = 10
57                easting_min = project_truescale.xminShallow
58                easting_max = project_truescale.xmaxShallow
59                northing_min = project_truescale.yminShallow
60                northing_max = project_truescale.ymaxShallow
61
62    ##    which_area='all'
63    ##    easting_min=0
64    ##    easting_max=2179.2
65    ##    northing_min=0
66    ##    northing_max=1360.8
67    ##    cellsize=10
68       
69           
70            var = [2]
71               
72            for which_var in var:
73                if which_var == 0:  # Stage
74                    outname = directory + which_area + '_stage_' + str(timestep)
75                    quantityname = 'stage'
76
77                if which_var == 1:  # Absolute Momentum
78                    outname = directory + which_area + '_momentum_' + str(timestep)
79                    quantityname = '(xmomentum**2 + ymomentum**2)**0.5' 
80
81                if which_var == 2:  # Depth
82                    outname = directory + which_area + '_depth_' + str(timestep)
83                    quantityname = 'stage-elevation' 
84
85                if which_var == 3:  # Speed
86                    outname = directory + which_area + '_speed_' + str(timestep)
87                    quantityname = '(xmomentum**2 + ymomentum**2)**0.5/(stage-elevation+1.e-6/(stage-elevation))'  #Speed
88
89                if which_var == 4:  # Elevation
90                    outname = directory + which_area + '_elevation_' + str(timestep)
91                    quantityname = 'elevation' 
92
93    ##            if is_parallel == True:
94    ##            #    print 'is_parallel',is_parallel
95    ##                for i in range(0,nodes):
96    ##                    namei = name + '_P%d_%d' %(i,nodes)
97    ##                    outnamei = outname + '_P%d_%d' %(i,nodes)
98    ##                    print 'start sww2dem for sww file %d' %(i)
99    ##                    sww2dem(namei, basename_out = outnamei,
100    ##                                quantity = quantityname,
101    ##                                timestep = timestep,
102    ##                                cellsize = cellsize,     
103    ##                                easting_min = project_grad.e_min_area,
104    ##                                easting_max = project_grad.e_max_area,
105    ##                                northing_min = project_grad.n_min_area,
106    ##                                northing_max = project_grad.n_max_area,       
107    ##                                reduction = max,
108    ##                                verbose = True,
109    ##                                format = 'asc')
110    ##            else:
111
112               
113                print 'start sww2dem', which_area
114                sww2dem(name, basename_out = outname,
115                            quantity = quantityname,
116                            timestep = timestep,
117                            cellsize = cellsize,     
118                            easting_min = easting_min,
119                            easting_max = easting_max,
120                            northing_min = northing_min,
121                            northing_max = northing_max,       
122                            reduction = max, 
123                            verbose = True,
124                            format = 'ers')
125
126
127##convert_dem_from_ascii2netcdf(outname, use_cache=True, verbose=True)
128##dem2pts(outname, use_cache=True, verbose=True)
129##G = Geospatial_data(file_name = outname + '.pts')
130##G.export_points_file('test' + '.txt')
Note: See TracBrowser for help on using the repository browser.