Calculates cold D-alpha emission
subroutine cold_spec
!+ Calculates cold D-alpha emission
integer :: ic, i, j, k, it, is
real(Float64), dimension(3) :: ri, vhalo, ri_cyl
integer,dimension(3) :: ind
!! Determination of the CX probability
type(LocalProfiles) :: plasma
real(Float64) :: cold_photons
integer :: n = 10000
!$OMP PARALLEL DO schedule(dynamic,1) private(i,j,k,ic,is,ind, &
!$OMP& cold_photons, ri, vhalo, plasma, ri_cyl)
loop_over_cells: do ic = istart, spec_chords%cyl_ncell, istep
call ind2sub(pass_grid%dims,spec_chords%cyl_cell(ic),ind)
i = ind(1) ; j = ind(2) ; k = ind(3)
ri_cyl(1) = pass_grid%r(i) + 0.5d0*pass_grid%dr
ri_cyl(2) = pass_grid%z(j) + 0.5d0*pass_grid%dz
ri_cyl(3) = pass_grid%phi(k) + 0.5d0*pass_grid%dphi
call cyl_to_uvw(ri_cyl, ri)
call get_plasma(plasma, pos=ri, input_coords=1)
do is = 1, n_thermal
cold_photons = plasma%denn(3,is)*tables%einstein(2,3)
if(cold_photons.le.0.0) cycle loop_over_cells
do it=1, n
!! Cold Spectra
call mc_halo(plasma, thermal_mass(is), vhalo)
call store_photons(ri, vhalo, thermal_lambda0(is), cold_photons/n, &
&spec%cold(:,:,:,is), spec%coldstokes(:,:,:,:,is), passive=.True.)
enddo
enddo
enddo loop_over_cells
!$OMP END PARALLEL DO
#ifdef _MPI
!! Combine Spectra
call parallel_sum(spec%cold)
#endif
end subroutine cold_spec