T1 = 10;
Hs = 5;

f = 0:.001:0.5;  %Hz
w = f* 2 * pi;

Swhoi = (2 * pi) * (Hs^2*T1) .* (0.11/(2*pi)) .* ((w.*T1./(2*pi)).^-5) .* exp(-0.44.*(w.*T1./(2*pi)).^-4);

A = 172.75 * Hs^2 / T1^4;
B = 691 / T1^4;

Sittc = (A./w.^5) .* exp(-B./w.^4) * (2 * pi);

a = .0081;
g = 980;
B = 0.74;

Vw =  1530;      % 30 kt wind speed in cm/sec at 19.5 m above the surface

Spm = ((a .* g^2 ./ w.^5) .* exp( (-B * g^4 / Vw^4)./ w.^4) ) * 2 * pi / 10000.0;

plot(f,Swhoi, f,Sittc, f,Spm)
title('Wave Spectral Density')
xlabel('Frequency (Hz)')
ylabel('Wave Spectral Density (m^2/Hz)')
