close all
clear all


a = 2;
d = 2*a;
rho = 1025;
g = 9.81;

V = (2*pi*a^3)/3;

m = rho*V;
S = pi*a^2;

filename = 'hemisphere_r2';

FD_coefs = Readwdra3dResults([filename '.out']);
    
fh1 = figure
plot(FD_coefs.k*d,real(FD_coefs.F{3,3})/(rho*V), ...
     FD_coefs.k*d,imag(FD_coefs.F{3,3})./(rho*V*FD_coefs.w), ...
     FD_coefs.k*d,real(FD_coefs.F{7,3})./(rho*V*(FD_coefs.w).^2), ...
     FD_coefs.k*d,imag(FD_coefs.F{7,3})./(rho*V*(FD_coefs.w).^2) ...
  )
xlabel(['$ k d $'],'interpreter','latex','fontsize',18);
legend('mu','lambda','real(F)','imag(F)');

    
fh1b = figure
plot(FD_coefs.T,real(FD_coefs.F{3,3})/(rho*V), ...
     FD_coefs.T,imag(FD_coefs.F{3,3})./(rho*V*FD_coefs.w), ...
     FD_coefs.T,real(FD_coefs.F{7,3})./(rho*V*(FD_coefs.w).^2), ...
     FD_coefs.T,imag(FD_coefs.F{7,3})./(rho*V*(FD_coefs.w).^2) ...
  )
xlabel(['$ k d $'],'interpreter','latex','fontsize',18);
legend('mu','lambda','real(F)','imag(F)');


title('Hemisphere Heave Hydrodynamics');
%print('-dpdf','HemisphereHeaveHydrodynamics');
  
w = FD_coefs.w;
mu = real(FD_coefs.F{3,3});
lambda = imag(FD_coefs.F{3,3});
X = FD_coefs.F{7,3}./(-(m+mu).*w.^2+i*w.*lambda+rho*g*S);

X2 = abs(FD_coefs.F{7,3})./(-(m+mu).*w.^2+i*w.*lambda+rho*g*S);

fh2 = figure
subplot(2,1,1)
plot(FD_coefs.k*d,real(X),FD_coefs.k*d,imag(X),FD_coefs.k*d,abs(X))
xlabel(['$ k d $'],'interpreter','latex','fontsize',18)
title('Hemisphere Heave Motions');
subplot(2,1,2)
plot(FD_coefs.k*d,angle(X)*180/pi)
xlabel(['$ k d $'],'interpreter','latex','fontsize',18)
ylabel('Phase angle (deg)');


fh2 = figure
subplot(2,1,1)
plot(FD_coefs.T,real(X),FD_coefs.T,imag(X),FD_coefs.T,abs(X))
xlabel(['$ k d $'],'interpreter','latex','fontsize',18)
title('Hemisphere Heave Motions');
subplot(2,1,2)
plot(FD_coefs.T,angle(X)*180/pi)
xlabel(['$ k d $'],'interpreter','latex','fontsize',18)
ylabel('Phase angle (deg)');