Changeset 7684


Ignore:
Timestamp:
Apr 13, 2010, 4:31:20 PM (10 years ago)
Author:
fountain
Message:
 
File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/anuga_1_1/anuga_work/production/bunbury_storm_surge_2009/process_gems_grids.py

    r7683 r7684  
    139139number_of_points = ncols * nrows
    140140
     141# insert nodata test here
     142
    141143for i in xrange(number_of_timesteps):
    142144    j = number_of_points * i
    143145    k = j + number_of_points
    144146    d = stage[j:k] - elevation     # Assumes elevations below sea level are negative
    145     depth[j:k] = [0 if x<-500 else x for x in d]   # sets stage to zero for no-data values of stage
     147    depth[j:k] = [0 if x<-9000 else x for x in d]   # sets stage to zero for no-data values of stage
    146148   
    147149    # depth[j:k] = [x if x>=0 and elevation>0 else 0 for x in d] # This is to correct for GEMS data which assumes stage=0 onshore when dry but should = elevation
     
    162164        mx = momentum[i]*math.sin(math.radians(t)) #Assuming t is in the "to" direction and 0 degrees is north
    163165        my = momentum[i]*math.cos(math.radians(t))
    164         assert mx >= 0 and my >= 0
    165     if t > 90.0 and t < 180.0:
     166        assert mx > 0 and my > 0
     167    elif t > 90.0 and t < 180.0:
    166168        mx = momentum[i]*math.sin(math.radians(t))
    167169        my = momentum[i]*math.cos(math.radians(t))
    168         assert mx >= 0 and my <= 0
    169     if t > 180.0 and t < 270.0:
     170        assert mx > 0 and my < 0
     171    elif t > 180.0 and t < 270.0:
    170172        mx = momentum[i]*math.cos(math.radians(t))
    171173        my = momentum[i]*math.sin(math.radians(t))
    172         assert mx <= 0 and my <= 0
    173     if t > 270.0 and t < 360.0:
     174        assert mx < 0 and my < 0
     175    elif t > 270.0 and t < 360.0:
    174176        mx = momentum[i]*math.sin(math.radians(t))
    175177        my = momentum[i]*math.cos(math.radians(t))
    176         assert mx <= 0 and my >= 0
    177     if t == 0.0 or t == 360.0:
     178        assert mx < 0 and my > 0
     179    elif t == 0.0 or t == 360.0:
    178180        mx = 0
    179181        my = momentum[i]
    180182        assert my >= 0
    181     if t == 90.0:
     183    elif t == 90.0:
    182184        mx = momentum[i]
    183185        my = 0
    184186        assert mx >= 0
    185     if t == 180.0:
     187    elif t == 180.0:
    186188        mx = 0
    187189        my = -momentum[i]
    188190        assert my <= 0
    189     if t == 270.0:
     191    elif t == 270.0:
    190192        mx = -momentum[i]
    191193        my = 0
    192194        assert mx <= 0
    193     if t == -9999:
     195    elif t == -9999:
    194196        mx = 0
    195197        my = 0
     198    else:
     199        print "Unexpected value of theta"
     200        exit()
    196201       
    197202    xmomentum[i] = mx
Note: See TracChangeset for help on using the changeset viewer.