| 12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455565758596061626364656667686970717273747576777879 |
- function eos_nemo, t_in, psal
- ;------------------------------------------------------------------------------
- ;
- ; Purpose:
- ; To calculate the density using the equation of state used in NEMO.
- ; This is based on the Jackett and McDougall (1994) equation of state
- ; for calculating the in situ density based on potential temperature
- ; and salinity.
- ;
- ; Inputs:
- ; temperature => 1d array potential temperature (deg C)
- ; salinity => 1d array salinity (PSU)
- ;
- ; Outputs:
- ; rhop => 1d array pot density (kg/m**3)
- ;
- ; could add use a difference reference pressure
- ; observations are in-situ temperature?
- ;
- ; Author:
- ; D. J. Lea. Dec 2006.
- ;------------------------------------------------------------------------------
- zws = SQRT( ABS( psal ) )
- sz=size(t_in)
- ndim=sz(0)
- NO_LEVS=sz(ndim)
-
- prhop=psal*0.
- ;
- for jk = 0L, NO_LEVS-1 do begin
- zt = T_IN(jk)
- zs = psal(jk)
- ; * depth
- ; zh = O_DEP_LEVS(jk) ;used in calculating insitu density only
- zsr= zws(jk)
- ; * compute volumic mass pure water at atm pressure
- zr1= ( ( ( ( 6.536332e-9*zt-1.120083e-6 )*zt+1.001685e-4)*zt $
- -9.095290e-3 )*zt+6.793952e-2 )*zt+999.842594
- ; * seawater volumic mass atm pressure
- zr2= ( ( ( 5.3875e-9*zt-8.2467e-7 ) *zt+7.6438e-5 ) *zt $
- -4.0899e-3 ) *zt+0.824493
- zr3= ( -1.6546e-6*zt+1.0227e-4 ) *zt-5.72466e-3
- zr4= 4.8314e-4
- ; * potential volumic mass (reference to the surface)
- zrhop= ( zr4*zs + zr3*zsr + zr2 ) *zs + zr1
- ; * save potential volumic mass
- prhop(jk) = zrhop
- ; * add the compression terms
- ze = ( -3.508914e-8*zt-1.248266e-8 ) *zt-2.595994e-6
- zbw= ( 1.296821e-6*zt-5.782165e-9 ) *zt+1.045941e-4
- zb = zbw + ze * zs
- zd = -2.042967e-2
- zc = (-7.267926e-5*zt+2.598241e-3 ) *zt+0.1571896
- zaw= ( ( 5.939910e-6*zt+2.512549e-3 ) *zt-0.1028859 ) *zt - 4.721788
- za = ( zd*zsr + zc ) *zs + zaw
- zb1= (-0.1909078*zt+7.390729 ) *zt-55.87545
- za1= ( ( 2.326469e-3*zt+1.553190)*zt-65.00517 ) *zt+1044.077
- zkw= ( ( (-1.361629e-4*zt-1.852732e-2 ) *zt-30.41638 ) *zt + 2098.925 ) *zt+190925.6
- zk0= ( zb1*zsr + za1 )*zs + zkw
- ; ; in situ density anomaly
- ; prd(jk) = zrhop / ( 1.0 - zh / ( zk0 - zh * ( za - zh * zb ) ) )
- ;
- endfor
- return, prhop
- END
|
