source: development/stochastic_study/plot_spread.py @ 3011

Last change on this file since 3011 was 3011, checked in by ole, 18 years ago

Stochastic plotting

File size: 2.6 KB
Line 
1"""Read in timeseries from N simulations and plot spread for each timestep at specified locations (ch 5,7,9)
2"""
3
4# Standard modules
5import sys, os, string
6
7# Related major packages
8from Numeric import zeros, Float, allclose, arange
9
10# Application specific imports
11from utilities.numerical_tools import histogram, create_bins
12from read_realisations import read_realisations
13import project
14
15# Initialise
16# Read in all realisations-timeseries
17study = 'cyclone1'
18study = 'cyclone3'
19#study = 'nautilus3'
20
21time, data, filenames = read_realisations(study, 
22                                          #max_realisations = 200,
23                                          gauge_number=2,
24                                          use_cache=False)
25number_of_realisations = data.shape[1]
26print 'Read %d realisations' %number_of_realisations
27
28# Plot
29from pylab import ion, hold, plot, title, xlabel, ylabel, legend, savefig, show
30
31ion()
32hold(True)
33#hold(False)
34
35
36for j in range(20):
37    #print j, filenames[j]
38    plot(time, data[:,j], 'k-')
39   
40    xlabel('time(s)')
41    ylabel('stage (m)')
42
43    shortname = string.join(os.path.split(filenames[j])[-1].split('_')[-2:], '_')
44    title('Study %s: timeseries for %s (realisation %d of %d)'\
45          %(study, shortname, j, number_of_realisations-1))
46   
47    #title('Study %s: timeseries for %s (realisation %d of %d)'\
48    #      %(study, filenames[j], j, number_of_realisations-1))   
49    raw_input('Next')
50 
51
52
53hold(False)
54# Plot spread of stage values for each timestep
55for i in range(300,320): #project.number_of_timesteps):
56    # Plot histogram#
57
58    w = data[i,:]
59    plot(w, 'k.')
60    xlabel('realisations')
61    ylabel('stage (m)')
62    title('Study %s: spread at timestep %d of %d (t=%.2f)'\
63          %(study, i, project.number_of_timesteps-1, time[i]))       
64    raw_input('Next')
65
66
67hold(False)
68# Plot histogram of stage values for each timestep
69for i in range(300,320):#project.number_of_timesteps):
70    # Plot histogram
71
72    w = data[i,:]
73     
74    bins = create_bins(w, project.number_of_bins)
75    print 'bins', bins
76    hist = histogram(w, bins, relative=True)
77    #print 'hist', hist
78   
79    #plot(w, 'k.')
80    plot(bins, hist, 'k.')
81    xlabel('stage (m)')
82    ylabel('frequency')
83    title('Study %s: histogram at timestep %d of %d (t=%.2f)'\
84          %(study, i, project.number_of_timesteps-1, time[i]))       
85    raw_input('Next')
86
87
88#title('Gauge %s' %name)
89#xlabel('time(s)')
90#ylabel('stage (m)')   
91#legend(('Observed', 'Modelled'), shadow=True, loc='upper left')
92#savefig(name, dpi = 300)
93
94
95
96
97
98
99show()
100
101
102
103#from pylab import *
104#plot(time, stage)
105#show()
Note: See TracBrowser for help on using the repository browser.