Changeset 7645


Ignore:
Timestamp:
Mar 1, 2010, 5:47:50 PM (15 years ago)
Author:
ole
Message:

Work on Patong validation regression test

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  
    1111
    1212from anuga.abstract_2d_finite_volumes.util import file_function
     13from anuga.utilities.numerical_tools import cov, get_machine_precision
     14
     15
     16import pylab
    1317import numpy as num
    1418import sys, os
     
    4347    gauge_names = gauges.keys()
    4448
    45     print sww_filename
    46     print gauge_locations
    47     os.system('cp %s tempfile.sww' % sww_filename) # Has end with sww
    48     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,
    4953                      quantities='stage',
    5054                      interpolation_points=gauge_locations,
     
    6266
    6367
    64     return timevector, timeseries
    65 
    66 
    67 
    68 def compare_timeseries(swwfile1, swwfile2):
     68    return num.array(timevector), timeseries
     69
     70
     71   
     72def 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
     157def compare_all_timeseries(swwfile1, swwfile2):
    69158   
    70159    timevector1, timeseries1 = get_timeseries(swwfile1, gauges)
     
    74163          (swwfile1, swwfile2)
    75164    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   
    77185    print 'OK'
    78186   
    79 
     187   
    80188
    81189if __name__ == '__main__':
     
    85193        import sys; sys.exit(-1)
    86194
    87     res = compare_timeseries(sys.argv[1], sys.argv[2])
     195    #res = compare_all_timeseries(sys.argv[1], sys.argv[2])
    88196    try:
    89         res = compare_timeseries(sys.argv[1], sys.argv[2])
     197        res = compare_all_timeseries(sys.argv[1], sys.argv[2])
    90198    except:
    91199        sys.exit(-1)
  • anuga_validation/automated_validation_tests/patong_beach_validation/project.py

    r7640 r7645  
    3636
    3737use_buildings = True
    38 setup = 'trial'         # This can be one of three values
     38setup = 'basic'         # This can be one of three values
    3939                        #    trial - coarsest mesh, fast
    4040                        #    basic - coarse mesh
  • anuga_validation/automated_validation_tests/patong_beach_validation/validate.py

    r7643 r7645  
    5353                         'patong.sww.FINAL.tgz'
    5454                        )
    55 #Optional_Data_Objects = ('patong.sww.TRIAL.tgz',)
     55Optional_Data_Objects = ('patong.sww.BASIC.tgz',)
    5656
    5757# path to the local data directory
     
    349349    #cmd = 'python cmpsww.py %s %s > cmpsww.stdout' % (local_sww, new_output_sww)
    350350    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   
    351355    log.debug("check_that_output_is_as_expected: doing '%s'" % cmd)
    352356    res = os.system(cmd)
Note: See TracChangeset for help on using the changeset viewer.