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 |
---|
5 | import sys, os, string |
---|
6 | |
---|
7 | # Related major packages |
---|
8 | from Numeric import zeros, Float, allclose, arange |
---|
9 | |
---|
10 | # Application specific imports |
---|
11 | from utilities.numerical_tools import histogram, create_bins |
---|
12 | from read_realisations import read_realisations |
---|
13 | import project |
---|
14 | |
---|
15 | # Initialise |
---|
16 | # Read in all realisations-timeseries |
---|
17 | |
---|
18 | #study = 'nautilus1' #~70 realisations, blocks of 100, sequential, stddev= 0.0006 |
---|
19 | #study = 'nautilus3' #~100 realisations, blocks of 10, sequential, stddev= 0.0013 |
---|
20 | study = 'cyclone1' #~4000 realisations, blocks of 100, parallel, stddev= 0.0006 |
---|
21 | #study = 'cyclone3' #~120 realisations, blocks of 10, parallel, stddev= 0.0013 |
---|
22 | |
---|
23 | time, data, filenames = read_realisations(study, |
---|
24 | #max_realisations = 200, |
---|
25 | gauge_number=0, |
---|
26 | sorting='randomised', |
---|
27 | #sorting='numerical', |
---|
28 | verbose=True, |
---|
29 | use_cache=True) |
---|
30 | number_of_realisations = data.shape[1] |
---|
31 | number_of_timesteps = data.shape[0] |
---|
32 | print 'Read %d realisations' %number_of_realisations |
---|
33 | |
---|
34 | |
---|
35 | # Plot |
---|
36 | from pylab import ion, hold, plot, title, xlabel, ylabel, legend, savefig, show |
---|
37 | |
---|
38 | ion() |
---|
39 | hold(True) |
---|
40 | #hold(False) |
---|
41 | |
---|
42 | |
---|
43 | for j in range(1): |
---|
44 | #print j, filenames[j] |
---|
45 | plot(time, data[:,j], 'k-') |
---|
46 | |
---|
47 | xlabel('time(s)') |
---|
48 | ylabel('stage (m)') |
---|
49 | |
---|
50 | shortname = string.join(os.path.split(filenames[j])[-1].split('_')[-2:], '_') |
---|
51 | title('Study %s: timeseries for %s (realisation %d of %d)'\ |
---|
52 | %(study, shortname, j, number_of_realisations-1)) |
---|
53 | |
---|
54 | #raw_input('Next') |
---|
55 | |
---|
56 | |
---|
57 | raw_input('Stats') |
---|
58 | |
---|
59 | |
---|
60 | |
---|
61 | hold(False) |
---|
62 | # Plot spread of stage values for each timestep |
---|
63 | for i in range(0,project.number_of_timesteps): |
---|
64 | # Plot histogram# |
---|
65 | |
---|
66 | w = data[i,:] |
---|
67 | plot(w, 'k.') |
---|
68 | xlabel('realisations') |
---|
69 | ylabel('stage (m)') |
---|
70 | title('Study %s: spread at timestep %d of %d (t=%.2f)'\ |
---|
71 | %(study, i, project.number_of_timesteps-1, time[i])) |
---|
72 | raw_input('Next') |
---|
73 | |
---|
74 | show() |
---|
75 | import sys; sys.exit() |
---|
76 | |
---|
77 | hold(False) |
---|
78 | # Plot histogram of stage values for each timestep |
---|
79 | for i in range(300,320):#project.number_of_timesteps): |
---|
80 | # Plot histogram |
---|
81 | |
---|
82 | w = data[i,:] |
---|
83 | |
---|
84 | bins = create_bins(w, project.number_of_bins) |
---|
85 | print 'bins', bins |
---|
86 | hist = histogram(w, bins, relative=True) |
---|
87 | #print 'hist', hist |
---|
88 | |
---|
89 | #plot(w, 'k.') |
---|
90 | plot(bins, hist, 'k.') |
---|
91 | xlabel('stage (m)') |
---|
92 | ylabel('frequency') |
---|
93 | title('Study %s: histogram at timestep %d of %d (t=%.2f)'\ |
---|
94 | %(study, i, project.number_of_timesteps-1, time[i])) |
---|
95 | raw_input('Next') |
---|
96 | |
---|
97 | |
---|
98 | #title('Gauge %s' %name) |
---|
99 | #xlabel('time(s)') |
---|
100 | #ylabel('stage (m)') |
---|
101 | #legend(('Observed', 'Modelled'), shadow=True, loc='upper left') |
---|
102 | #savefig(name, dpi = 300) |
---|
103 | |
---|
104 | |
---|
105 | |
---|
106 | |
---|
107 | |
---|
108 | |
---|
109 | show() |
---|
110 | |
---|
111 | |
---|
112 | |
---|
113 | #from pylab import * |
---|
114 | #plot(time, stage) |
---|
115 | #show() |
---|