% Soundabsorption.m
%
% Computes sound absorption coefficient at TAPS-6 frequencies
%
% User must enter:
%        mean water temperature in C
%        mean water salinity in ppt
%
% Outputs are absorption in dB/m

clear all
close all

fxv=[104 165 265 420 700 1100 1850 3000]*1000;

format bank
fsqkhz = (fxv ./ 1000) .^ 2;
temp = input('Mean water temp (C) = ');
sal  = input('Mean water salinity (ppt) = ');

%          Coefs from Urick,1967, pp.88, 89 and 92.
%                   alpha is in dB/m.

acoef = 1.86e-2;
bcoef = 2.68e-2;
tranf = 21.9 * (10 .^ (6-(1520 ./ (temp + 273))));
lfterm = (acoef .* sal .* tranf .* fsqkhz) ...
             ./ (tranf .^ 2 + fsqkhz);
hfterm = (bcoef .* fsqkhz) ./ tranf;
x=((0.1 * fsqkhz) ./ (1 + fsqkhz));
y=((40 * fsqkhz) ./ (4100 + fsqkhz));
vlowterm = x+y;
alpha1 = vlowterm+lfterm+hfterm;
alpha = (39.37007874 ./ (36 * 1000)) * alpha1;


loglog(fxv,alpha,'b',fxv,alpha,'bo')
grid
s=['Seawater absorption @ T = ' num2str(temp) ...
    'C, S = ' num2str(sal) ' o/oo'];
title (s)
xlabel('Frequency in Hz')
ylabel('Absorption in dB/m')


disp('');disp(s)
for k=1:length(fxv)
fprintf('Freq = %4.0f KHz -- alpha = %6.2f dB/m\n',fxv(k)/1000,alpha(k))
end

