1
-
!! Attenuation module for converting sound speed and attenuation to complex values
2
-
!! Provides routines to convert sound speed and attenuation in user units to complex sound speed
4
-
!! Acoustic attenuation calculations including volume attenuation formulas and unit conversions
7
-
! Routines to convert a sound speed and attenuation in user units to a complex sound speed
8
-
! Includes a formula for volume attenuation
15
-
INTEGER, PRIVATE, PARAMETER :: PRTFile = 6
16
-
INTEGER, PARAMETER :: MaxBioLayers = 200
17
-
INTEGER :: iBio, NBioLayers
18
-
REAL (KIND=8) :: T = 20, Salinity = 35, pH = 8, z_bar = 0, FG ! Francois-Garrison volume attenuation; temperature, salinity, ...
21
-
REAL (KIND=8) :: Z1, Z2, f0, Q, a0
22
-
END TYPE bioStructure
24
-
TYPE( bioStructure ) :: bio( MaxBioLayers )
27
-
!! Converts real wave speed and attenuation to complex wave speed
28
83*
FUNCTION CRCI( z, c, alpha, freq, freq0, AttenUnit, beta, fT )
30
-
! Converts real wave speed and attenuation to a single
31
-
! complex wave speed (with positive imaginary part)
34
-
! 6 Cases: N Nepers/meter
35
-
! M dB/meter (M for Meters)
36
-
! m dB/meter with a power law
37
-
! F dB/m-kHz (F for frequency dependent)
38
-
! W dB/wavelength (W for Wavelength)
42
-
! second letter adds volume attenuation according to standard laws:
44
-
! F for Francois Garrison
47
-
! freq is the current frequency
48
-
! freq0 is the reference frequency for which the dB/meter was specified (used only for 'm')
51
-
! c real part of sound speed
52
-
! alpha imaginary part of sound speed
56
-
REAL (KIND=8), INTENT( IN ) :: freq, freq0, alpha, c, z, beta, fT
57
-
CHARACTER (LEN=2), INTENT( IN ) :: AttenUnit
58
-
REAL (KIND=8) :: f2, omega, alphaT, Thorp, a, FG
59
-
COMPLEX (KIND=8) :: CRCI
61
83
omega = 2.0 * pi * freq
63
-
! Convert to Nepers/m
65
#####
SELECT CASE ( AttenUnit( 1 : 1 ) )
69
#####
alphaT = alpha / 8.6858896D0
70
-
CASE ( 'm' ) ! dB/m with power law
71
#####
alphaT = alpha / 8.6858896D0
72
#####
IF ( freq < fT ) THEN ! frequency raised to the power beta
73
#####
alphaT = alphaT * ( freq / freq0 ) ** beta
74
-
ELSE ! linear in frequency
75
#####
alphaT = alphaT * ( freq / freq0 ) * ( fT / freq0 ) ** ( beta - 1 )
77
-
CASE ( 'F' ) ! dB/(m kHz)
78
52
alphaT = alpha * freq / 8685.8896D0
79
-
CASE ( 'W' ) ! dB/wavelength
80
31*
IF ( c /= 0.0 ) alphaT = alpha * freq / ( 8.6858896D0 * c )
81
-
! The following lines give f^1.25 frequency dependence
82
-
! FAC = SQRT( SQRT( freq / 50.0 ) )
83
-
! IF ( c /= 0.0 ) alphaT = FAC * alpha * freq / ( 8.6858896D0 * c )
84
-
CASE ( 'Q' ) ! Quality factor
85
#####
IF ( c * alpha /= 0.0 ) alphaT = omega / ( 2.0 * c * alpha )
86
-
CASE ( 'L' ) ! loss parameter
87
83*
IF ( c /= 0.0 ) alphaT = alpha * omega / c
90
-
! added volume attenuation
91
#####
SELECT CASE ( AttenUnit( 2 : 2 ) )
93
#####
f2 = ( freq / 1000.0 ) ** 2
95
-
! Original formula from Thorp 1967
96
-
! Thorp = 40.0 * f2 / ( 4100.0 + f2 ) + 0.1 * f2 / ( 1.0 + f2 ) ! dB/kyard
97
-
! Thorp = Thorp / 914.4D0 ! dB / m
98
-
! Thorp = Thorp / 8.6858896D0 ! Nepers / m
100
-
! Updated formula from JKPS Eq. 1.34
101
#####
Thorp = 3.3d-3 + 0.11 * f2 / ( 1.0 + f2 ) + 44.0 * f2 / ( 4100.0 + f2 ) + 3d-4 * f2 ! dB/km
102
#####
Thorp = Thorp / 8685.8896 ! Nepers / m
104
#####
alphaT = alphaT + Thorp
105
-
CASE ( 'F' ) ! Francois-Garrison
106
10
FG = Franc_Garr( freq / 1000 ); ! dB/km
107
10
FG = FG / 8685.8896; ! Nepers / m
108
10*
alphaT = alphaT + FG
109
-
CASE ( 'B' ) ! biological attenuation per Orest Diachok
110
83*
DO iBio = 1, NBioLayers
111
#####
IF ( z >= bio( iBio )%Z1 .AND. z <= bio( iBio )%Z2 ) THEN
112
#####
a = bio( iBio )%a0 / ( ( 1.0 - bio( iBio )%f0 ** 2 / freq ** 2 ) ** 2 + 1.0 / bio( iBio )%Q ** 2 ) ! dB/km
113
#####
a = a / 8685.8896 ! Nepers / m
114
#####
alphaT = alphaT + a
119
-
! Convert Nepers/m to equivalent imaginary sound speed
120
83
alphaT = alphaT * c * c / omega
121
83
CRCI = CMPLX( c, alphaT, KIND=8 )
123
83
IF ( alphaT > c ) THEN
124
#####
WRITE( PRTFile, * ) 'Complex sound speed: ', CRCI
125
#####
WRITE( PRTFile, * ) 'Usually this means you have an attenuation that is way too high'
126
#####
CALL ERROUT( 'AttenMod : CRCI ', 'The complex sound speed has an imaginary part > real part' )
132
-
! **********************************************************************!
134
-
!! Computes Francois-Garrison volume attenuation
135
10
FUNCTION Franc_Garr( f )
137
-
! Francois Garrison formulas for attenuation
138
-
! Based on a Matlab version by D. Jackson APL-UW
141
-
! Verified using F-G Table IV
143
-
! alpha = attenuation (dB/km)
144
-
! f = frequency (kHz)
145
-
! T = temperature (deg C)
146
-
! S = salinity (psu)
147
-
! pH = 7 for neutral water
148
-
! z_bar = depth (m)
151
-
! alpha = volume attenuation in dB/km
153
-
REAL (KIND=8), INTENT( IN ) :: f
154
-
REAL (KIND=8) :: Franc_Garr
155
-
REAL (KIND=8) :: c, A1, A2, A3, P1, P2, P3, f1, f2
156
-
! LP: Bug (at least technically): Single-precision and double-precision
157
-
! literals are all mixed together here, both including values which are
158
-
! not representable as floats, thus leading to different results.
159
-
! Practically, these are approximate, experimentally-derived constants, so
160
-
! errors at the 1e-7 scale are completely not meaningful.
162
10
c = 1412 + 3.21 * T + 1.19 * Salinity + 0.0167 * z_bar
164
-
! Boric acid contribution
165
10
A1 = 8.86 / c * 10 ** ( 0.78 * pH - 5 )
167
10
f1 = 2.8 * sqrt( Salinity / 35 ) * 10 ** ( 4 - 1245 / ( T + 273 ) )
169
-
! Magnesium sulfate contribution
170
10
A2 = 21.44 * Salinity / c * ( 1 + 0.025 * T )
171
10
P2 = 1 - 1.37D-4 * z_bar + 6.2D-9 * z_bar ** 2
172
10
f2 = 8.17 * 10 ** ( 8 - 1990 / ( T + 273 ) ) / ( 1 + 0.0018 * ( Salinity - 35 ) )
175
10
P3 = 1 - 3.83D-5 * z_bar + 4.9D-10 * z_bar ** 2
176
10
if ( T < 20 ) THEN
177
10
A3 = 4.937D-4 - 2.59D-5 * T + 9.11D-7 * T ** 2 - 1.5D-8 * T ** 3
179
#####
A3 = 3.964D-4 -1.146D-5 * T + 1.45D-7 * T ** 2 - 6.5D-10 * T ** 3
182
-
Franc_Garr = A1 * P1 * ( f1 * f ** 2 ) / ( f1 ** 2 + f ** 2 ) + A2 * P2 * ( f2 * f ** 2 ) / ( f2 ** 2 + f ** 2 ) + &
185
10
END FUNCTION Franc_Garr
187
#####
END MODULE AttenMod