%matplotlib inline
from pylab import *
from clawpack.geoclaw import topotools
from clawpack.visclaw import geoplot, plottools
import nc_tools
import os,sys,copy
loc = 'Flattery'
event = 'CSZ_L1'
fgmax_file1 = '../Runs/%s/%s_v570/results_%s_%s_8hr_20s.nc' % (loc,event,loc,event)
fgmax_file2 = '../Runs/%s/%s_v570/results_%s_%s_2hr_3s.nc' % (loc,event,loc,event)
fgmax1 = nc_tools.read_nc(fgmax_file1)
fgmax2 = nc_tools.read_nc(fgmax_file2)
print('max difference in B = %g' % abs(fgmax1.B - fgmax2.B).max())
print('max difference in dz = %g' % abs(fgmax1.dz - fgmax2.dz).max())
print('minimum change in h = %.2f m' % (fgmax1.h - fgmax2.h).min())
print('maximum change in h = %.2f m' % (fgmax1.h - fgmax2.h).max())
print('minimum change in s = %.2f m/s' % (fgmax1.s - fgmax2.s).min())
print('maximum change in s = %.2f m/s' % (fgmax1.s - fgmax2.s).max())
figure(figsize=(13,13))
plottools.pcolorcells(fgmax1.X, fgmax1.Y, fgmax1.s - fgmax2.s, cmap=cm.get_cmap('seismic'))
clim(-5,5)
colorbar(shrink=0.7, extend='both')
contour(fgmax1.X, fgmax1.Y, fgmax1.B, [0], colors='g', linewidths=0.5)
gca().set_aspect(1./cos(48*pi/180))
title('Difference in s (old - new)');
figure(figsize=(13,13))
plottools.pcolorcells(fgmax1.X, fgmax1.Y, fgmax1.h - fgmax2.h, cmap=cm.get_cmap('seismic'))
clim(-1,1)
colorbar(shrink=0.7, extend='both')
contour(fgmax1.X, fgmax1.Y, fgmax1.B, [0], colors='g', linewidths=0.5)
gca().set_aspect(1./cos(48*pi/180))
title('Difference in h (old - new)');
fgmax_file3 = '../Runs/%s/%s_v570/results_%s_%s_merged.nc' % (loc,event,loc,event)
os.system('cp %s %s' % (fgmax_file2, fgmax_file3))
fgmax3 = copy.copy(fgmax2)
fgmax3.h = maximum(fgmax1.h, fgmax2.h)
fgmax3.s = maximum(fgmax1.s, fgmax2.s)
fgmax3.hss = maximum(fgmax1.hss, fgmax2.hss)
fgmax3.hmin = minimum(fgmax1.hmin, fgmax2.hmin)
nc_tools.write_nc_output(fgmax_file3, fgmax3, new=False, force=True)