Changeset 7645
- Timestamp:
- Mar 1, 2010, 5:47:50 PM (15 years ago)
- Location:
- anuga_validation/automated_validation_tests/patong_beach_validation
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
anuga_validation/automated_validation_tests/patong_beach_validation/compare_model_timeseries.py
r7642 r7645 11 11 12 12 from anuga.abstract_2d_finite_volumes.util import file_function 13 from anuga.utilities.numerical_tools import cov, get_machine_precision 14 15 16 import pylab 13 17 import numpy as num 14 18 import sys, os … … 43 47 gauge_names = gauges.keys() 44 48 45 print sww_filename46 print gauge_locations 47 os.system('cp %s tempfile.sww' % sww_filename) # Has end with sww48 f = file_function( 'tempfile.sww', #_filename,49 tempfile = 'xyz1234tempfile.sww' # Has to end with sww 50 51 os.system('cp %s %s' % (sww_filename, tempfile)) 52 f = file_function(tempfile, 49 53 quantities='stage', 50 54 interpolation_points=gauge_locations, … … 62 66 63 67 64 return timevector, timeseries 65 66 67 68 def compare_timeseries(swwfile1, swwfile2): 68 return num.array(timevector), timeseries 69 70 71 72 def compare_timeseries(timevector, 73 timeseries1, 74 timeseries2, 75 name='', 76 legend = ('Time series 1', 'Time series 2')): 77 """Compare and plot two timeseries 78 """ 79 80 81 timeseries1 = num.array(timeseries1) 82 timeseries2 = num.array(timeseries2) 83 N = timevector.shape[0] 84 assert timeseries1.shape[0] == N 85 assert timeseries2.shape[0] == N 86 87 eps = 1.0e-3 88 89 # Covariance measure 90 #res = cov(timeseries1, timeseries2) 91 #msg = 'Covariance between timeseries was too large: %e' % res 92 #assert res < 0.01, msg 93 94 # Difference measures 95 print timeseries1 96 print timeseries1-timeseries2 97 98 99 # Maximum norm 100 nominator = max(abs(timeseries1-timeseries2)) 101 denominator = max(abs(timeseries1)) 102 print nominator, denominator, N 103 if denominator > 0: 104 # Relative measure 105 res = nominator/denominator 106 else: 107 # Absolute measure 108 res = nominator/N 109 110 msg = '%s: Difference between timeseries was too large: %e' % (name, res) 111 #assert res < eps, msg 112 113 114 nominator = sum(abs(timeseries1-timeseries2)) 115 denominator = sum(abs(timeseries1)) 116 print nominator, denominator, N 117 if denominator > 0: 118 # Relative measure 119 res = nominator/denominator 120 else: 121 # Absolute measure 122 res = nominator/N 123 124 msg = '%s: Difference between timeseries was too large: %e' % (name, res) 125 #assert res < eps, msg 126 127 # Extrema 128 #res = max(model) 129 #report_difference('Maximum', res, expected_maximum[name], rtol, atol) 130 131 132 133 # Generate plots 134 pylab.ion() # No plotting on screen 135 pylab.hold(False) 136 137 pylab.plot(timevector, timeseries1, 'r-', 138 timevector, timeseries2, 'k-') 139 140 pylab.title('Gauge %s' % name) 141 pylab.xlabel('time(s)') 142 pylab.ylabel('stage (m)') 143 pylab.legend(legend, shadow=True, loc='upper left') 144 pylab.savefig(name, dpi = 300) 145 146 # Error vector 147 pylab.ion() # No plotting on screen 148 pylab.hold(False) 149 pylab.plot(timevector, timeseries1-timeseries2, 'b-') 150 pylab.title('Gauge %s (difference)' % name) 151 pylab.xlabel('time(s)') 152 pylab.ylabel('stage difference (m)') 153 pylab.savefig(name + '_diff', dpi = 300) 154 155 156 157 def compare_all_timeseries(swwfile1, swwfile2): 69 158 70 159 timevector1, timeseries1 = get_timeseries(swwfile1, gauges) … … 74 163 (swwfile1, swwfile2) 75 164 assert num.allclose(timevector1, timevector2), msg 76 165 timevector = timevector1 166 167 # Check that both timeseries exist for all gauges 168 for name in timeseries2: 169 assert name in timeseries1 170 171 for name in timeseries1: 172 assert name in timeseries2 173 174 # Compare all timeseries data individually 175 for name in timeseries1: 176 compare_timeseries(timevector, 177 timeseries1[name], 178 timeseries2[name], 179 name=name) 180 181 182 183 184 77 185 print 'OK' 78 186 79 187 80 188 81 189 if __name__ == '__main__': … … 85 193 import sys; sys.exit(-1) 86 194 87 res = compare_timeseries(sys.argv[1], sys.argv[2])195 #res = compare_all_timeseries(sys.argv[1], sys.argv[2]) 88 196 try: 89 res = compare_ timeseries(sys.argv[1], sys.argv[2])197 res = compare_all_timeseries(sys.argv[1], sys.argv[2]) 90 198 except: 91 199 sys.exit(-1) -
anuga_validation/automated_validation_tests/patong_beach_validation/project.py
r7640 r7645 36 36 37 37 use_buildings = True 38 setup = ' trial' # This can be one of three values38 setup = 'basic' # This can be one of three values 39 39 # trial - coarsest mesh, fast 40 40 # basic - coarse mesh -
anuga_validation/automated_validation_tests/patong_beach_validation/validate.py
r7643 r7645 53 53 'patong.sww.FINAL.tgz' 54 54 ) 55 #Optional_Data_Objects = ('patong.sww.TRIAL.tgz',)55 Optional_Data_Objects = ('patong.sww.BASIC.tgz',) 56 56 57 57 # path to the local data directory … … 349 349 #cmd = 'python cmpsww.py %s %s > cmpsww.stdout' % (local_sww, new_output_sww) 350 350 cmd = 'python compare_model_timeseries.py %s %s > compare_model_timeseries.stdout' % (local_sww, new_output_sww) 351 print '-------------------------------------' 352 print cmd 353 print '-------------------------------------' 354 351 355 log.debug("check_that_output_is_as_expected: doing '%s'" % cmd) 352 356 res = os.system(cmd)
Note: See TracChangeset
for help on using the changeset viewer.