source: branches/anuga_1_1/anuga_work/production/bunbury_storm_surge_2009/sts2csv.py @ 7729

"""
Function to read in sts file and convert to csv file for analysis purposes.
Adapted from build_urs_boundary.py from the tsunami modelling space.
May 2010.
"""

import os
import os.path
from time import localtime, strftime, gmtime

import project

from Scientific.IO.NetCDF import NetCDFFile
#import Numeric as num

home = os.path.join(os.getenv(ENV_INUNDATIONHOME), 'data') # Absolute path for data folder
anuga_folder = os.path.join(home, state, scenario_folder, 'anuga')
boundaries_folder = os.path.join(anuga_folder, 'boundaries')
event_folder = os.path.join(boundaries_folder, 'gcom_60min')

filename = join(event_folder,'gcom_60_min_boundary.sts')

#-------------------------------------------------------------------------------
# Get gauges (timeseries of index points)
#-------------------------------------------------------------------------------
def get_sts_gauge_data(filename, verbose=False):
    print 'get_sts_gauge_data: filename=%s' % filename
    fid = NetCDFFile(filename+'.sts', 'r')      #Open existing file for read
    #permutation = fid.variables['permutation'][:]
    x = fid.variables['x'][:] + fid.xllcorner   #x-coordinates of vertices
    y = fid.variables['y'][:] + fid.yllcorner   #y-coordinates of vertices
    #points = num.transpose(num.asarray([x.tolist(), y.tolist()]))
    time = fid.variables['time'][:] #+ fid.starttime
    elevation = fid.variables['elevation'][:]
       
    basename = 'sts_gauge'
    quantity_names = ['stage', 'xmomentum', 'ymomentum']
    quantities = {}
    for i, name in enumerate(quantity_names):
        quantities[name] = fid.variables[name][:]

    #---------------------------------------------------------------------------
    # Get timeseries values for wave height and components of momentum
    #---------------------------------------------------------------------------

    for j in range(len(x)):
        index = j # permutation[j]
        stage = quantities['stage'][:,j]
        xmomentum = quantities['xmomentum'][:,j]
        ymomentum = quantities['ymomentum'][:,j]

        out_file = os.path.join(project.event_folder,
                                basename+'_'+str(index)+'.csv')
        fid_sts = open(out_file, 'w')
        fid_sts.write('time, stage, xmomentum, ymomentum \n')

        #-----------------------------------------------------------------------
        # End of the get gauges
        #-----------------------------------------------------------------------
        for k in range(len(time)-1):
            fid_sts.write('%.6f, %.6f, %.6f, %.6f\n'
                          % (time[k], stage[k], xmomentum[k], ymomentum[k]))

        fid_sts.close()      

    fid.close()

    return quantities,elevation,time