source: anuga_work/development/Hinwood_2008/slope.py @ 5459

"""
Plot up files from the Hinwood project.
"""
from os import sep
import project
#from scipy import arange

from Numeric import arange, array, zeros, Float

from anuga.fit_interpolate.interpolate import interpolate_sww2csv
from anuga.shallow_water.data_manager import csv2dict

def load_depths(slope_file):
    #slope, _ = csv2dict(file_sim)
    
    # Read the depth file
    dfid = open(slope_file)
    lines = dfid.readlines()
    dfid.close()

    title = lines.pop(0)
    n_time = len(lines)
    n_sensors = len(lines[0].split(','))-1  # -1 to remove time
    dtimes = zeros(n_time, Float)  #Time
    depths = zeros(n_time, Float)  #
    sensors = zeros((n_time,n_sensors), Float)
    depth_locations = title.split(',') #(',')
    depth_locations.pop(0) # remove 'time'
    # !!!! -3 assumes a 
    depths = [float(j.split(':')[0]) for j in depth_locations]
    
    for i, line in enumerate(lines):
        fields = line.split(',') #(',')
        fields = [float(j) for j in fields]
        dtimes[i] = fields[0]
        sensors[i] = fields[1:] # 1: to remove time

    #print "dtimes",dtimes
    #print "depths", depths
    #print "sensors", sensors
    return dtimes, depths, sensors
   
def load_slopes(slope_file):
    times, depths, sensors = load_depths(slope_file)
    n_slope_locations = len(depths)-1
    n_time = len(times)
    slope_locations = zeros(n_slope_locations, Float)  #
    slopes = zeros((n_time,n_slope_locations), Float)

    # An array of the sensor spacing values
    delta_depths = zeros(n_slope_locations, Float)
    
    for i in arange(n_slope_locations):
        slope_locations[i] = (depths[i+1] + depths[i+1])/2.
        delta_depths[i] = (depths[i+1] - depths[i])
    
    for j in arange(n_time):
        for i in arange(n_slope_locations):
            slopes[j,i] = (sensors[j,i+1] - sensors[j,i])/delta_depths[i]

    return times, slope_locations, slopes

def graph_slopes(slope_file, break_xs=None, break_times=None):
    from pylab import meshgrid, cm, contourf, contour, ion, plot, xlabel, \
         ylabel, close, legend, savefig, title, figure ,colorbar, show , axis
    origin = 'lower'
    times, slope_locations, slopes = load_slopes(slope_file)
    #print "times", times
    #print "slope_locations", slope_locations
    #print "slopes", slopes

    # Can't seem to reshape this info once it is in the function
    CS = contourf(slope_locations, times, slopes, 10, # [-1, -0.1, 0, 0.1],
                  #alpha=0.5,
                  cmap=cm.bone,
                  origin=origin)
    
    #CS2 = contour(slope_locations, times, slopes, CS.levels[::2],
       #           colors = 'r',
        #          origin=origin,
         #         hold='on')

    title('slope')
    xlabel('x location')
    ylabel('Time, seconds')
    #axis([5.0, 5.5, 30, 60])


    if break_times is not None and break_xs is not None:
        plot(break_xs, break_times, 'ro')
    
    # Make a colorbar for the ContourSet returned by the contourf call.
    cbar = colorbar(CS)
    cbar.ax.set_ylabel('slope')
    # Add the contour line levels to the colorbar
    #cbar.add_lines(CS2)

    figure()
    show()

def auto_graph_slopes():
    from scenarios import scenarios
    
    outputdir_tag = "_good_tri_area_0.01_A"
    outputdir_tag = "_good_tri_area_0.01_old"
    #outputdir_tag = "_test"
    scenarios = [scenarios[1]] # !!!!!!!!!!!!!!!!!!!!!!
    for run_data in scenarios:
        id = run_data['scenario_id']
        outputdir_name = id + outputdir_tag
        pro_instance = project.Project(['data','flumes','Hinwood_2008'],
                                       outputdir_name=outputdir_name)
        end = id + ".csv"
        slope_file = pro_instance.outputdir + "bslope_stage_" + end
        graph_slopes(slope_file, run_data['break_xs'],
                     run_data['break_times']) 
    
def gauges_for_slope(outputdir_tag, scenarios):
    
    outputdir_tag = "_good_tri_area_0.01_A"
    #outputdir_tag = "_test"
    #scenarios = [scenarios[0]] # !!!!!!!!!!!!!!!!!!!!!!
    dx = 0.05
    for run_data in scenarios:
        point_x = arange(run_data['start_slope_x'],
                        run_data['finish_slope_x'],
                        dx).tolist()
        flume_y_middle = 0.5
        points = []
        for gauge_x in point_x:
            points.append([gauge_x, flume_y_middle])
        id = run_data['scenario_id']
    
        basename = 'zz_' + run_data['scenario_id']
        outputdir_name = id + outputdir_tag
        pro_instance = project.Project(['data','flumes','Hinwood_2008'],
                                       outputdir_name=outputdir_name)
        end = id + ".csv"
        interpolate_sww2csv(pro_instance.outputdir + basename +".sww",
                            points,
                            pro_instance.outputdir + "fslope_depth_" + end,
                            pro_instance.outputdir + "fslope_velocity_x_" + end,
                            pro_instance.outputdir + "fslope_velocity_y_" + end,
                            pro_instance.outputdir + "fslope_stage_" + end,
                            time_thinning=1)
  
    
    

#-------------------------------------------------------------
if __name__ == "__main__":
    """ 
    """
    from scenarios import scenarios
    scenarios = [scenarios[0]]
    #gauges_for_slope("_good_tri_area_0.01_A", scenarios)
    #auto_graph_slopes()