Changeset 2622
- Timestamp:
- Mar 28, 2006, 2:36:03 PM (19 years ago)
- Files:
-
- 51 added
- 8 edited
Legend:
- Unmodified
- Added
- Removed
-
inundation/pyvolution/view_tsh.py
r1363 r2622 9 9 from shallow_water import Domain 10 10 from pmesh2domain import pmesh_to_domain_instance 11 from util import file_function, Polygon_function, read_polygon 11 from pyvolution.util import file_function 12 from utilities.polygon import Polygon_function, read_polygon 12 13 from Numeric import zeros, Float, maximum, minimum 13 14 from realtime_visualisation_new import * -
inundation/utilities/polygon.py
r2531 r2622 10 10 print 'Could not find scipy - using Numeric' 11 11 from Numeric import Float, Int, zeros, ones, array, concatenate, reshape, dot 12 12 13 13 14 14 from math import sqrt … … 45 45 def inside_polygon(points, polygon, closed = True, verbose = False): 46 46 """Determine points inside a polygon 47 47 48 48 Functions inside_polygon and outside_polygon have been defined in 49 terms af separate_by_polygon which will put all inside indices in 49 terms af separate_by_polygon which will put all inside indices in 50 50 the first part of the indices array and outside indices in the last 51 52 See separate_points_by_polygon for documentation 51 52 See separate_points_by_polygon for documentation 53 53 """ 54 54 55 55 if verbose: print 'Checking input to inside_polygon' 56 56 57 57 try: 58 58 points = ensure_numeric(points, Float) … … 66 66 polygon = ensure_numeric(polygon, Float) 67 67 except NameError, e: 68 raise NameError, e 68 raise NameError, e 69 69 except: 70 70 msg = 'Polygon %s could not be converted to Numeric array' %(str(polygon)) … … 90 90 def outside_polygon(points, polygon, closed = True, verbose = False): 91 91 """Determine points outside a polygon 92 92 93 93 Functions inside_polygon and outside_polygon have been defined in 94 terms af separate_by_polygon which will put all inside indices in 94 terms af separate_by_polygon which will put all inside indices in 95 95 the first part of the indices array and outside indices in the last 96 97 See separate_points_by_polygon for documentation 96 97 See separate_points_by_polygon for documentation 98 98 """ 99 99 … … 102 102 points = ensure_numeric(points, Float) 103 103 except NameError, e: 104 raise NameError, e 104 raise NameError, e 105 105 except: 106 106 msg = 'Points could not be converted to Numeric array' … … 110 110 polygon = ensure_numeric(polygon, Float) 111 111 except NameError, e: 112 raise NameError, e 112 raise NameError, e 113 113 except: 114 114 msg = 'Polygon could not be converted to Numeric array' … … 133 133 # No points are outside 134 134 return [] 135 else: 135 else: 136 136 return indices[count:][::-1] #return reversed 137 137 … … 180 180 #Input checks 181 181 182 182 183 183 try: 184 184 points = ensure_numeric(points, Float) 185 185 except NameError, e: 186 raise NameError, e 186 raise NameError, e 187 187 except: 188 188 msg = 'Points could not be converted to Numeric array' … … 192 192 polygon = ensure_numeric(polygon, Float) 193 193 except NameError, e: 194 raise NameError, e 194 raise NameError, e 195 195 except: 196 196 msg = 'Polygon could not be converted to Numeric array' … … 276 276 if verbose: print 'Checking input to separate_points_by_polygon' 277 277 278 278 279 279 #Input checks 280 280 281 281 assert isinstance(closed, bool), 'Keyword argument "closed" must be boolean' 282 assert isinstance(verbose, bool), 'Keyword argument "verbose" must be boolean' 283 284 282 assert isinstance(verbose, bool), 'Keyword argument "verbose" must be boolean' 283 284 285 285 try: 286 286 points = ensure_numeric(points, Float) 287 287 except NameError, e: 288 raise NameError, e 288 raise NameError, e 289 289 except: 290 290 msg = 'Points could not be converted to Numeric array' … … 295 295 polygon = ensure_numeric(polygon, Float) 296 296 except NameError, e: 297 raise NameError, e 297 raise NameError, e 298 298 except: 299 299 msg = 'Polygon could not be converted to Numeric array' … … 391 391 392 392 #Make points in polygons relative to geo_reference 393 self.regions = [] 393 self.regions = [] 394 394 for polygon, value in regions: 395 P = geo_reference.change_points_geo_ref(polygon) 395 P = geo_reference.change_points_geo_ref(polygon) 396 396 self.regions.append( (P, value) ) 397 397 … … 429 429 return z 430 430 431 def read_polygon(filename ):431 def read_polygon(filename,split=','): 432 432 """Read points assumed to form a polygon. 433 433 There must be exactly two numbers in each line separated by a comma. … … 441 441 polygon = [] 442 442 for line in lines: 443 fields = line.split( ',')443 fields = line.split(split) 444 444 polygon.append( [float(fields[0]), float(fields[1])] ) 445 445 … … 453 453 number_of_points - (optional) number of points 454 454 seed - seed for random number generator (default=None) 455 exclude - list of polygons (inside main polygon) from where points should be excluded 455 exclude - list of polygons (inside main polygon) from where points should be excluded 456 456 457 457 Output: … … 489 489 490 490 append = True 491 491 492 492 #Check exclusions 493 493 if exclude is not None: … … 495 495 if inside_polygon( [x,y], ex_poly ): 496 496 append = False 497 498 499 if append is True: 500 points.append([x,y]) 497 498 499 if append is True: 500 points.append([x,y]) 501 501 502 502 return points … … 530 530 else: 531 531 x_delta = x+x_mult*x*delta 532 532 533 533 if y == 0: 534 534 y_delta = y_mult*delta 535 535 else: 536 536 y_delta = y+y_mult*y*delta 537 537 538 538 point = [x_delta, y_delta] 539 539 #print "point",point 540 540 if inside_polygon( point, polygon, closed = False ): 541 541 raise Found 542 except Found: 542 except Found: 543 543 point_in = True 544 544 else: 545 delta = delta*0.1 545 delta = delta*0.1 546 546 return point 547 547 -
inundation/zeus/Merimbula.zpi
r2255 r2622 73 73 <ReleaseProjectSaveAll>Off</ReleaseProjectSaveAll> 74 74 <ReleaseProjectReload>Off</ReleaseProjectReload> 75 <file>..\merimbula\prepare.py</file> 76 <file>..\merimbula\project.py</file> 77 <file>..\merimbula\run_cairns.py</file> 78 <file>..\merimbula\run_least_squares_merimbula.py</file> 79 <file>..\merimbula\run_merimbula_lake.py</file> 80 <file>..\merimbula\run_merimbula_lake_checkpoint.py</file> 81 <file>..\merimbula\run_new_meribula.py</file> 82 <file>..\merimbula\view_sww.py</file> 75 <file>..\..\production\merimbula_2005\prepare.py</file> 76 <file>..\..\production\merimbula_2005\project.py</file> 77 <file>..\..\production\merimbula_2005\run_merimbula_lake.py</file> 78 <file>..\..\production\merimbula_2005\run_merimbula_lake_weed.py</file> 79 <file>..\..\production\merimbula_2005\run_new_meribula.py</file> 80 <file>..\..\production\merimbula_2005\run_new_meribula_weed.py</file> 83 81 <folder name="Header Files" /> 84 82 <folder name="Resource Files" /> -
inundation/zeus/anuga-workspace.zwi
r2386 r2622 2 2 <workspace name="anuga-workspace"> 3 3 <mode>Debug</mode> 4 <active> pyvolution</active>4 <active>Merimbula</active> 5 5 <project name="analytic_solutions">analytic_solutions.zpi</project> 6 6 <project name="euler">euler.zpi</project> -
inundation/zeus/pyvolution.zpi
r2255 r2622 76 76 <file>..\pyvolution\bed_w_eden_boundary.py</file> 77 77 <file>..\pyvolution\bed_w_file_boundary.py</file> 78 <file>..\pyvolution\cg_solve.py</file>79 78 <file>..\pyvolution\check_sww_tsh.py</file> 80 79 <file>..\pyvolution\combine_pts.py</file> … … 90 89 <file>..\pyvolution\island.py</file> 91 90 <file>..\pyvolution\least_squares.py</file> 92 <file>..\pyvolution\load_mesh\loadASCII.py</file>93 91 <file>..\pyvolution\log.ini</file> 94 92 <file>..\pyvolution\mesh.py</file> … … 112 110 <file>..\pyvolution\shallow_water_vtk.py</file> 113 111 <file>..\pyvolution\show_balanced_limiters.py</file> 114 <file>..\pyvolution\sparse.py</file>115 <file>..\pyvolution\sparse_ext.c</file>116 112 <file>..\pyvolution\test_advection.py</file> 117 113 <file>..\pyvolution\test_all.py</file> 118 <file>..\pyvolution\test_cg_solve.py</file>119 114 <file>..\pyvolution\test_combine_pts.py</file> 120 115 <file>..\pyvolution\test_data_manager.py</file> … … 131 126 <file>..\pyvolution\test_region.py</file> 132 127 <file>..\pyvolution\test_shallow_water.py</file> 133 <file>..\pyvolution\test_sparse.py</file>134 128 <file>..\pyvolution\test_util.py</file> 135 129 <file>..\pyvolution\treenode.py</file> … … 137 131 <file>..\pyvolution\twolevels.py</file> 138 132 <file>..\pyvolution\util.py</file> 139 <file>..\pyvolution\util_ext.c</file>140 <file>..\pyvolution\util_ext.h</file>141 133 <file>..\pyvolution\view_tsh.py</file> 142 134 <file>..\pyvolution\vtk_realtime_visualiser.py</file> -
production/merimbula_2005/prepare.py
r2225 r2622 1 1 2 2 3 def prepare_ timeboundary(filename):4 """Converting tide timeseries to NetCDF tms file.3 def prepare_wind_stress(filename): 4 """Converting wind timeseries to NetCDF tms file. 5 5 This is a 'throw-away' code taylor made for files like 6 6 'Benchmark_2_input.txt' from the LWRU2 benchmark 7 7 """ 8 8 9 print 'Preparing time boundary from%s' %filename9 print 'Preparing wind timeseries %s' %filename 10 10 from Numeric import array, zeros, Float, asarray 11 11 … … 13 13 14 14 #Skip first line 15 line = fid.readline()15 #line = fid.readline() 16 16 17 17 #Read remaining lines … … 22 22 N = len(lines) 23 23 T = zeros(N, Float) #Time 24 Q = zeros(N, Float) #Values 24 S = zeros(N, Float) #Speed 25 B = zeros(N, Float) #Bearing 25 26 26 27 Told = 0.0 27 28 Sold = ' ' 28 29 Lold = ' ' 30 29 31 for i, line in enumerate(lines): 30 32 fields = line.split() … … 46 48 T[i] = T[i] - Tstart 47 49 #this is specific to this data set. deals with daylight saving 48 if i>3270:49 T[i] = T[i]+360050 50 # if i>3270: 51 # T[i] = T[i]+3600 52 # 51 53 if T[i]<Told : 52 54 print Lold … … 58 60 print i, T[i]-Told 59 61 60 Q[i] = float(fields[2]) 62 S[i] = float(fields[2]) 63 B[i] = float(fields[3]) 61 64 62 65 Told = T[i] … … 65 68 66 69 70 #print T 71 #Create tms file 72 from Scientific.IO.NetCDF import NetCDFFile 73 74 outfile = filename[:-4] + '.tms' 75 print 'Writing to', outfile 76 fid = NetCDFFile(outfile, 'w') 77 78 fid.institution = 'Australian National University' 79 fid.description = 'Input wind for Merimbula' 80 fid.starttime = 0.0 81 fid.createDimension('number_of_timesteps', len(T)) 82 fid.createVariable('time', Float, ('number_of_timesteps',)) 83 fid.variables['time'][:] = T 84 85 fid.createVariable('speed', Float, ('number_of_timesteps',)) 86 fid.variables['speed'][:] = S[:] 87 88 fid.createVariable('bearing', Float, ('number_of_timesteps',)) 89 fid.variables['bearing'][:] = B[:] 90 91 92 fid.close() 93 94 95 def prepare_timeboundary(filename): 96 """Converting tide time series to NetCDF tms file. 97 This is a 'throw-away' code taylor made for files like 98 'Benchmark_2_input.txt' from the LWRU2 benchmark 99 """ 100 101 print 'Preparing time boundary from %s' %filename 102 from Numeric import array, zeros, Float, asarray 103 104 fid = open(filename) 105 106 #Skip first line 107 line = fid.readline() 108 109 #Read remaining lines 110 lines = fid.readlines() 111 fid.close() 112 113 114 N = len(lines) 115 T = zeros(N, Float) #Time 116 Q = zeros(N, Float) #Values 117 118 Told = 0.0 119 Sold = ' ' 120 Lold = ' ' 121 for i, line in enumerate(lines): 122 fields = line.split() 123 124 #print fields 125 126 l_time = (fields[0]+' '+fields[1])[0:-1] 127 from time import strptime, mktime 128 129 s_time = strptime(l_time,'%d/%m/%y %H:%M:%S') 130 131 #print s_time 132 133 T[i] = float(mktime(s_time)) 134 135 if i==0: 136 Tstart = T[0] 137 138 T[i] = T[i] - Tstart 139 #this is specific to this data set. deals with daylight saving 140 if i>3270: 141 T[i] = T[i]+3600 142 143 if T[i]<Told : 144 print Lold 145 print l_time 146 print Sold 147 print s_time 148 print Told 149 print T[i] 150 print i, T[i]-Told 151 152 Q[i] = float(fields[2]) 153 154 Told = T[i] 155 Sold = s_time 156 Lold = l_time 157 158 159 #print T 67 160 #Create tms file 68 161 from Scientific.IO.NetCDF import NetCDFFile … … 90 183 fid.close() 91 184 92 93 185 #------------------------------------------------------------- 94 186 if __name__ == "__main__": … … 96 188 print 'Prepare Open sea boundary condition from ',project.original_boundary_filename 97 189 prepare_timeboundary(project.original_boundary_filename ) 190 191 192 print 'Prepare wind from ',project.original_wind_filename 193 prepare_wind_stress(project.original_wind_filename ) 98 194 99 195 #Preparing points -
production/merimbula_2005/project.py
r2225 r2622 3 3 4 4 original_boundary_filename = 'Eden_tide_Sept03.dat' 5 original_wind_filename = 'merimbula_wind_sept_2003.dat' 5 6 boundary_filename = 'Eden_tide_Sept03.tms' 6 7 bathymetry_filename = 'merimbula_bathymetry.xya' 7 8 mesh_filename = 'merimbula_10785.tsh' 8
Note: See TracChangeset
for help on using the changeset viewer.