close all
clear all

g = 9.81;

T = 4:.25:16;

w = 2*pi./T;

h = [10000000 100 50 25 10];

for j = 1:length(h)
 for i = 1:length(T)  
    k(i,j) = fzero(@(x) ShallowWaterDisp(x,h(j),w(i),g),(w(i)^2)/g);
    Vp(i,j) = w(i)/k(i,j);
    Vg(i,j) = (0.5+k(i,j)*h(j)/(sinh(2*k(i,j)*h(j))))*Vp(i,j);
    Lambda(i,j) = 2*pi/k(i,j);
 end
end

leg_str{1} = 'infinite depth';
for j = 2:length(h)
  leg_str{j} = ['depth = ' num2str(h(j)) ' m'];
end

figure('Position', [705 72 1138 1426],'PaperSize',[11 17],'PaperOrientation','portrait','PaperPositionMode','auto','PaperPosition',[.5 .5 10.5 16.5])

ah = subplot(3,1,1)
set(ah,'FontSize',12);
plot(T,Lambda)
xlabel('Wave Period (s)');
ylabel('Wave Length (m)');
legend(leg_str,'Location','NorthWest');
title('Linear Wave Theory');

ah = subplot(3,1,2)
set(ah,'FontSize',12);
plot(T,Vp)
xlabel('Wave Period (s)');
ylabel('Phase Velocity (m/s)');
legend(leg_str,'Location','NorthWest');



ah = subplot(3,1,3)
set(ah,'FontSize',12);
plot(T,Vg)
xlabel('Wave Period (s)');
ylabel('Group Velocity (m/s)');
legend(leg_str,'Location','NorthWest');


print('-dpdf','OceanWaveDispersionRelations');