Calculates doppler shift, stark splitting, and intensities
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=Float64), | intent(in), | dimension(3) | :: | vecp | Vector directing towards optical head |
|
| real(kind=Float64), | intent(in), | dimension(3) | :: | vi | Particle velocity |
|
| type(LocalEMFields), | intent(in) | :: | fields | Electro-magnetic fields |
||
| real(kind=Float64), | intent(in) | :: | lambda0 | Reference wavelength [nm] |
||
| real(kind=Float64), | intent(in) | :: | sigma_pi | Sigma-pi ratio |
||
| real(kind=Float64), | intent(in) | :: | photons | Photon density from colrad |
||
| real(kind=Float64), | intent(in) | :: | dlength | LOS intersection length with beam_grid cell particle is in |
||
| real(kind=Float64), | intent(out), | dimension(n_stark) | :: | lambda | Wavelengths [nm] |
|
| real(kind=Float64), | intent(out), | dimension(n_stark) | :: | intensity | ||
| real(kind=Float64), | intent(out), | dimension(n_stark,4) | :: | stokes | Spectra intensities [Ph/(s cm^2 starkline)] |
subroutine spectrum(vecp, vi, fields, lambda0, sigma_pi, photons, dlength, lambda, intensity, stokes)
!+ Calculates doppler shift, stark splitting, and intensities
real(Float64), dimension(3), intent(in) :: vecp
!+ Vector directing towards optical head
real(Float64), dimension(3), intent(in) :: vi
!+ Particle velocity
type(LocalEMFields), intent(in) :: fields
!+ Electro-magnetic fields
real(Float64), intent(in) :: lambda0
!+ Reference wavelength [nm]
real(Float64), intent(in) :: sigma_pi
!+ Sigma-pi ratio
real(Float64), intent(in) :: photons
!+ Photon density from [[libfida:colrad]]
real(Float64), intent(in) :: dlength
!+ LOS intersection length with [[libfida:beam_grid]] cell particle is in
real(Float64), dimension(n_stark), intent(out) :: lambda
!+ Wavelengths [nm]
real(Float64), dimension(n_stark), intent(out) :: intensity
real(Float64), dimension(n_stark,4), intent(out) :: stokes
!+ Spectra intensities [Ph/(s cm^2 starkline)]
real(Float64) :: m, h, normfactor
real(Float64), dimension(3) :: vp, vn
real(Float64), dimension(3) :: bfield, efield
real(Float64) :: E, B, cos_los_Efield, cos_los_Bfield, lambda_shifted, q0, q1, l0, szratio
integer, dimension(n_stark) :: stark_sign
real(Float64), dimension(n_stark) :: wavel
real(Float64), dimension(n_stark) :: circularity
stark_sign = +1*stark_sigma - 1*stark_pi
!! vector directing towards the optical head
vp=vecp/norm2(vecp)
! Calculate Doppler shift
vn=vi*0.01d0 ! [m/s]
lambda_shifted = lambda0*(1.d0 + dot_product(vn,vp)/c0)
!! Calculate Stark-Zeeman Splitting
! Calculate E-field
bfield = fields%b_norm*fields%b_abs
efield = fields%e_norm*fields%e_abs
efield(1) = efield(1) + vn(2)*bfield(3) - vn(3)*bfield(2)
efield(2) = efield(2) - (vn(1)*bfield(3) - vn(3)*bfield(1))
efield(3) = efield(3) + vn(1)*bfield(2) - vn(2)*bfield(1)
E = norm2(efield)
B = norm2(bfield)
!Stark-Zeeman Splitting
h=6.62607004d-34
! planck constant in SI units
m = 9.109384d-31
! mass of electron [kg]
l0 = lambda0*1d-9
! reference wavelength [m]
! stark-zeeman corrections to energy of n=2 states are -q0, 0 and q0
q0 = sqrt((e0*h*B/(4*pi*m))**2 + (3*a_0*e0*E)**2)
! stark-zeeman corrections to energy of n=3 states are -q1, -0.5*q1, 0, 0.5*q1, and q1/2
q1 = sqrt(4*(e0*h*B/(4*pi*m))**2 + 9*(3*a_0*e0*E)**2)
! szratio is gamma/epsilon factor. can be thought of as roughly zeeman energy shift/ stark energy shift.
if(E.eq.0.d0)then
szratio = 0.d0
else
szratio = (e0*h*B/(4*pi*m))/(3*a_0*e0*E)
endif
! wavelengths calculated from h*c0/lambda = E_i - E_j for transition from i to j energies
! order is small wavelengths to large wavelengths
if(n_stark.eq.15) then
wavel(1) = 2*c0*h*l0/(2*c0*h+(-2*q0*(-1)+q1*(2))*l0)
wavel(2) = 2*c0*h*l0/(2*c0*h+(-2*q0*(0)+q1*(2))*l0)
wavel(3) = 2*c0*h*l0/(2*c0*h+(-2*q0*(-1)+q1*(1))*l0)
wavel(4) = 2*c0*h*l0/(2*c0*h+(-2*q0*(1)+q1*(2))*l0)
wavel(5) = 2*c0*h*l0/(2*c0*h+(-2*q0*(0)+q1*(1))*l0)
wavel(6) = 2*c0*h*l0/(2*c0*h+(-2*q0*(-1)+q1*(0))*l0)
wavel(7) = 2*c0*h*l0/(2*c0*h+(-2*q0*(1)+q1*(1))*l0)
wavel(8) = 2*c0*h*l0/(2*c0*h+(-2*q0*(0)+q1*(0))*l0)
wavel(9) = 2*c0*h*l0/(2*c0*h+(-2*q0*(-1)+q1*(-1))*l0)
wavel(10) = 2*c0*h*l0/(2*c0*h+(-2*q0*(1)+q1*(0))*l0)
wavel(11) = 2*c0*h*l0/(2*c0*h+(-2*q0*(0)+q1*(-1))*l0)
wavel(12) = 2*c0*h*l0/(2*c0*h+(-2*q0*(-1)+q1*(-2))*l0)
wavel(13) = 2*c0*h*l0/(2*c0*h+(-2*q0*(1)+q1*(-1))*l0)
wavel(14) = 2*c0*h*l0/(2*c0*h+(-2*q0*(0)+q1*(-2))*l0)
wavel(15) = 2*c0*h*l0/(2*c0*h+(-2*q0*(1)+q1*(-2))*l0)
circularity = [1.0/3,4.0/3,5.0/3,-97.0/123,-1.0/2,-17.0/81,-1.0/3,0.0,1.0/3,17.0/81,1.0/2,97.0/123,-5.0/3,-4.0/3,-1.0/3]
circularity = circularity*szratio
else if(n_stark.eq.3) then
wavel(1) = c0*h*l0/(c0*h+q0*(-1)*l0)
wavel(2) = c0*h*l0/(c0*h+q0*(0)*l0)
wavel(3) = c0*h*l0/(c0*h+q0*(1)*l0)
circularity = [1.0,0.0,-1.0]
circularity = circularity*szratio
endif
!lambda = lambda_shifted + E * stark_wavel ![nm] <-old calculation for pure stark effect
wavel = (wavel-l0)*10**9 ! converting to [nm]
lambda = lambda_shifted + wavel
!Intensities of stark components
if (E .eq. 0.d0) then
cos_los_Efield = 0.d0
else
cos_los_Efield = dot_product(vp,efield) / E
endif
if (B .eq. 0.d0) then
cos_los_Bfield = 0.d0
else
cos_los_Bfield = dot_product(vp,bfield) / B
endif
!calculate stokes parameters
intensity = stark_intens
stokes(:,1) = intensity*(1.d0+ cos_los_Efield**2)*stark_sigma*sigma_pi + intensity*(1.d0- cos_los_Efield**2)*stark_pi
stokes(:,2) = intensity*(1.d0- cos_los_Efield**2)*stark_sigma*sigma_pi - intensity*(1.d0- cos_los_Efield**2)*stark_pi
stokes(:,3) = 0.d0
stokes(:,4) = 2*intensity*cos_los_Bfield*circularity
intensity = stark_intens*(1.d0+ stark_sign* cos_los_Efield**2)
!! E.g. mirrors may change the pi to sigma intensity ratio
where (stark_sigma .eq. 1)
intensity = intensity * sigma_pi
endwhere
!! normalize and multiply with photon density from colrad
normfactor = (1/sum(intensity))*photons*dlength
intensity = intensity/sum(intensity)*photons*dlength
stokes = stokes*normfactor
endsubroutine spectrum