source: development/stochastic_study/plot_spread.py @ 3010

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

beautification

File size: 2.4 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
17#study = 'cyclone1'
18study = 'nautilus3'
19
20time, data, filenames = read_realisations(study, 
21                                          max_realisations = 200,
22                                          use_cache=True)
23number_of_realisations = data.shape[1]
24
25
26# Plot
27from pylab import ion, hold, plot, title, xlabel, ylabel, legend, savefig, show
28
29ion()
30hold(True)
31#hold(False)
32
33
34for j in range(20):
35    #print j, filenames[j]
36    plot(time, data[:,j], 'k-')
37   
38    xlabel('time(s)')
39    ylabel('stage (m)')
40
41    shortname = string.join(os.path.split(filenames[j])[-1].split('_')[-2:], '_')
42    title('Study %s: timeseries for %s (realisation %d of %d)'\
43          %(study, shortname, j, number_of_realisations-1))
44   
45    #title('Study %s: timeseries for %s (realisation %d of %d)'\
46    #      %(study, filenames[j], j, number_of_realisations-1))   
47    raw_input('Next')
48 
49
50
51hold(False)
52# Plot spread of stage values for each timestep
53for i in range(300,320): #project.number_of_timesteps):
54    # Plot histogram#
55
56    w = data[i,:]
57    plot(w, 'k.')
58    xlabel('realisations')
59    ylabel('stage (m)')
60    title('Study %s: spread at timestep %d of %d (t=%.2f)'\
61          %(study, i, project.number_of_timesteps-1, time[i]))       
62    raw_input('Next')
63
64
65hold(False)
66# Plot histogram of stage values for each timestep
67for i in range(300,320):#project.number_of_timesteps):
68    # Plot histogram
69
70    w = data[i,:]
71     
72    bins = create_bins(w, project.number_of_bins)
73    print 'bins', bins
74    hist = histogram(w, bins, relative=True)
75    #print 'hist', hist
76   
77    #plot(w, 'k.')
78    plot(bins, hist, 'k.')
79    xlabel('stage (m)')
80    ylabel('frequency')
81    title('Study %s: histogram at timestep %d of %d (t=%.2f)'\
82          %(study, i, project.number_of_timesteps-1, time[i]))       
83    raw_input('Next')
84
85
86#title('Gauge %s' %name)
87#xlabel('time(s)')
88#ylabel('stage (m)')   
89#legend(('Observed', 'Modelled'), shadow=True, loc='upper left')
90#savefig(name, dpi = 300)
91
92
93
94
95
96
97show()
98
99
100
101#from pylab import *
102#plot(time, stage)
103#show()
Note: See TracBrowser for help on using the repository browser.