 function c = svel80(s,t,p)
% function c = svel80(s,t,p)
%
%-----------------------------------------------------------------------------
%     Speed of Sound
%
%     Input:
%            s  -  Salinity in PSU-78
%            t  -  Temperature in IPTS-68
%            p  -  Pressure in decibars
%
%     Output:
%            c - Speed of sound in m/s
%
%
%     References:
%                Chen & Millero (1977),jasa(62),1129-1135.
%
%
%     Checkvalues:
%            svel80(40,40,10000) = 1731.9954 m/s
%---------------------------------------------------------------------
%
%      code by n.p. fofonoff
%     implemented at uri/gso by JLL 3/1/82.
%
%    matlab version S. chiswell 1991

% convert p to bars and sqrt salinity
      p=p/10.0                                                  ;
      sr=sqrt(abs(s))                                            ;
%
% s.**2 term
      d=1.727e-3-7.9836e-6.*p                                     ;
% s.**3/2 term
      b1=7.3637e-5+1.7945e-7.*t                                   ;
      b0=-1.922e-2-4.42e-5.*t                                     ;
      b=b0+b1.*p                                                  ;
% s.**1 term
      a3=(-3.389e-13.*t+6.649e-12).*t+1.100e-10                    ;
      a2=((7.988e-12.*t-1.6002e-10).*t+9.1041e-9).*t-3.9064e-7      ;
      a1=(((-2.0122e-10.*t+1.0507e-8).*t-6.4885e-8).*t-1.2580e-5).*t  ...
       +9.4742e-5                                                 ;
      a0=(((-3.21e-8.*t+2.006e-6).*t+7.164e-5).*t-1.262e-2).*t+1.389  ;
      a=((a3.*p+a2).*p+a1).*p+a0                                     ;
% s.**0 term
      c3=(-2.3643e-12.*t+3.8504e-10).*t-9.7729e-9                    ;
      c2=(((1.0405e-12.*t-2.5335e-10).*t+2.5974e-8).*t-1.7107e-6).*t    ...
       +3.1260e-5                                                    ;
      c1=(((-6.1185e-10.*t+1.3621e-7).*t-8.1788e-6).*t+6.8982e-4).*t     ...
       +0.153563                                                       ;
      c0=((((3.1464e-9.*t-1.47800e-6).*t+3.342e-4).*t-5.80852e-2).*t     ...
       +5.03711).*t+1402.388                                          ;
      c=((c3.*p+c2).*p+c1).*p+c0                                        ;
% sound speed return
%
      c=c+(a+b.*sr+d.*s).*s                                           ;
%


