91 lines
2.1 KiB
Fortran
91 lines
2.1 KiB
Fortran
function voigte(a,vs)
|
|
c =====================
|
|
c
|
|
c computes a voigt function h = h(a,v)
|
|
c a=gamma/(4*pi*dnud) and v=(nu-nu0)/dnud. this is done after
|
|
c traving (landolt-b\rnstein, p. 449).
|
|
c
|
|
INCLUDE 'PARAMS.FOR'
|
|
dimension ak(19),a1(5)
|
|
data ak /-1.12470432, -0.15516677, 3.28867591, -2.34357915,
|
|
, 0.42139162, -4.48480194, 9.39456063, -6.61487486, 1.98919585,
|
|
, -0.22041650, 0.554153432, 0.278711796,-0.188325687, 0.042991293,
|
|
,-0.003278278, 0.979895023,-0.962846325, 0.532770573,-0.122727278/
|
|
data sqp/1.772453851/,sq2/1.414213562/
|
|
c
|
|
v = abs(vs)
|
|
u = a + v
|
|
v2 = v*v
|
|
if (a.eq.0.0) go to 140
|
|
if (a.gt.0.2) go to 120
|
|
if (v.ge.5.0) go to 121
|
|
c
|
|
ex=0.
|
|
if(v2.lt.100.)ex = exp(-v2)
|
|
k = 1
|
|
c
|
|
100 quo = 1.
|
|
if (v.lt.2.4) go to 101
|
|
quo = 1./(v2 - 1.5)
|
|
m = 11
|
|
go to 102
|
|
c
|
|
101 m = 6
|
|
if (v.lt.1.3) m = 1
|
|
102 do 103 i=1,5
|
|
a1(i) = ak(m)
|
|
m = m + 1
|
|
103 continue
|
|
h1 = quo*(a1(1) + v*(a1(2) + v*(a1(3) + v*(a1(4) + v*a1(5)))))
|
|
if (k.gt.1) go to 110
|
|
c
|
|
c a le 0.2 and v lt 5.
|
|
c
|
|
hh = h1*a + ex*(1. + a*a*(1. - 2.*v2))
|
|
voigte=hh
|
|
return
|
|
c
|
|
110 pqs = 2./sqp
|
|
h1p = h1 + pqs*ex
|
|
h2p = pqs*h1p - 2.*v2*ex
|
|
h3p = (pqs*(1. - ex*(1. - 2.*v2)) - 2.*v2*h1p)/3. + pqs*h2p
|
|
h4p = (2.*v2*v2*ex - pqs*h1p)/3. + pqs*h3p
|
|
psi = ak(16) + a*(ak(17) + a*(ak(18) + a*ak(19)))
|
|
c
|
|
c 0.2 lt a le 1.4 and a + v le 3.2
|
|
c
|
|
hh = psi*(ex + a*(h1p + a*(h2p + a*(h3p + a*h4p))))
|
|
voigte=hh
|
|
return
|
|
c
|
|
120 if (a.gt.1.4.or.u.gt.3.2) go to 130
|
|
ex=0.
|
|
if(v2.lt.100.)ex = exp(-v2)
|
|
k = 2
|
|
go to 100
|
|
c
|
|
c a le 0.2 and v ge 5.
|
|
c
|
|
121 hh = a*(15. + 6.*v2 + 4.*v2*v2)/(4.*v2*v2*v2*sqp)
|
|
voigte=hh
|
|
return
|
|
c
|
|
130 a2 = a*a
|
|
u = sq2*(a2 + v2)
|
|
u2 = 1./(u*u)
|
|
c
|
|
c a gt 1.4 or a + v gt 3.2
|
|
c
|
|
hh = sq2/sqp*a/u*(1. + u2*(3.*v2 - a2) +
|
|
, u2*u2*(15.*v2*v2 - 30.*v2*a2 + 3.*a2*a2))
|
|
voigte=hh
|
|
return
|
|
c
|
|
c a eq 0.
|
|
c
|
|
140 hh=0.
|
|
if(v2.lt.100.) hh=exp(-v2)
|
|
voigte=hh
|
|
return
|
|
end
|