"""Read in timeseries from N simulations and plot spread for each timestep at specified locations (ch 5,7,9) """ # Standard modules import sys, os, string # Related major packages from Numeric import zeros, Float, allclose, arange # Application specific imports from utilities.numerical_tools import histogram, create_bins from read_realisations import read_realisations import project # Initialise # Read in all realisations-timeseries #study = 'nautilus1' #~70 realisations, blocks of 100, sequential, stddev= 0.0006 #study = 'nautilus3' #~100 realisations, blocks of 10, sequential, stddev= 0.0013 study = 'cyclone1' #~4000 realisations, blocks of 100, parallel, stddev= 0.0006 #study = 'cyclone3' #~120 realisations, blocks of 10, parallel, stddev= 0.0013 time, data, filenames = read_realisations(study, #max_realisations = 200, gauge_number=0, sorting='randomised', #sorting='numerical', verbose=True, use_cache=True) number_of_realisations = data.shape[1] number_of_timesteps = data.shape[0] print 'Read %d realisations' %number_of_realisations # Plot from pylab import ion, hold, plot, title, xlabel, ylabel, legend, savefig, show ion() hold(True) #hold(False) for j in range(1): #print j, filenames[j] plot(time, data[:,j], 'k-') xlabel('time(s)') ylabel('stage (m)') shortname = string.join(os.path.split(filenames[j])[-1].split('_')[-2:], '_') title('Study %s: timeseries for %s (realisation %d of %d)'\ %(study, shortname, j, number_of_realisations-1)) #raw_input('Next') raw_input('Stats') hold(False) # Plot spread of stage values for each timestep for i in range(0,project.number_of_timesteps): # Plot histogram# w = data[i,:] plot(w, 'k.') xlabel('realisations') ylabel('stage (m)') title('Study %s: spread at timestep %d of %d (t=%.2f)'\ %(study, i, project.number_of_timesteps-1, time[i])) raw_input('Next') show() import sys; sys.exit() hold(False) # Plot histogram of stage values for each timestep for i in range(300,320):#project.number_of_timesteps): # Plot histogram w = data[i,:] bins = create_bins(w, project.number_of_bins) print 'bins', bins hist = histogram(w, bins, relative=True) #print 'hist', hist #plot(w, 'k.') plot(bins, hist, 'k.') xlabel('stage (m)') ylabel('frequency') title('Study %s: histogram at timestep %d of %d (t=%.2f)'\ %(study, i, project.number_of_timesteps-1, time[i])) raw_input('Next') #title('Gauge %s' %name) #xlabel('time(s)') #ylabel('stage (m)') #legend(('Observed', 'Modelled'), shadow=True, loc='upper left') #savefig(name, dpi = 300) show() #from pylab import * #plot(time, stage) #show()