subroutine integrate_spectra_datresol(igrid,wL,dwL,spectra,numWL,numXS,&
dir_loc,fl)
use mod_constants
integer, intent(in) :: igrid,numWL,numXS,dir_loc
type(te_fluid), intent(in) :: fl
double precision, intent(in) :: wL(numWL),dwL(numWL)
double precision, intent(inout) :: spectra(numWL,numXS)
integer :: direction_LOS
integer :: ixOmin1,ixOmin2,ixOmin3,ixOmax1,ixOmax2,ixOmax3,ixImin1,&
ixImin2,ixImin3,ixImax1,ixImax2,ixImax3,ix1,ix2,ix3,ixOnew
double precision, allocatable :: flux(:,:,:),v(:,:,:),pth(:,:,:),Te(:,:,&
:),rho(:,:,:)
double precision :: wlc,wlwd
integer :: mass
double precision :: logTe,lineCent
character (30) :: ion
double precision :: spaceRsl,wlRsl,sigma_PSF,wslit
integer :: levelg,rft,ixSmin,ixSmax,iwL
double precision :: flux_pix,dL
call get_line_info(spectrum_wl,ion,mass,logTe,lineCent,spaceRsl,wlRsl,&
sigma_PSF,wslit)
if (LOS_phi==0 .and. LOS_theta==90) then
direction_LOS=1
else if (LOS_phi==90 .and. LOS_theta==90) then
direction_LOS=2
else
direction_LOS=3
endif
ixOmin1=ixmlo1
ixOmin2=ixmlo2
ixOmin3=ixmlo3
ixOmax1=ixmhi1
ixOmax2=ixmhi2
ixOmax3=ixmhi3
ixImin1=ixglo1
ixImin2=ixglo2
ixImin3=ixglo3
ixImax1=ixghi1
ixImax2=ixghi2
ixImax3=ixghi3
allocate(flux(ixImin1:ixImax1,ixImin2:ixImax2,ixImin3:ixImax3),&
v(ixImin1:ixImax1,ixImin2:ixImax2,ixImin3:ixImax3),&
pth(ixImin1:ixImax1,ixImin2:ixImax2,ixImin3:ixImax3),&
Te(ixImin1:ixImax1,ixImin2:ixImax2,ixImin3:ixImax3),&
rho(ixImin1:ixImax1,ixImin2:ixImax2,ixImin3:ixImax3))
ix1=ixOmin1;ix2=ixOmin2;ix3=ixOmin3;
if (dir_loc==1) then
do ix1=ixOmin1,ixOmax1
if (location_slit>=(ps(igrid)%x(ix1,ix2,ix3,1)-half*ps(igrid)%dx(ix1,&
ix2,ix3,1)) .and. location_slit<(ps(igrid)%x(ix1,ix2,ix3,&
1)+half*ps(igrid)%dx(ix1,ix2,ix3,1))) then
ixOnew=ix1
endif
enddo
ixOmin1=ixOnew
ixOmax1=ixOnew
else if (dir_loc==2) then
do ix2=ixOmin2,ixOmax2
if (location_slit>=(ps(igrid)%x(ix1,ix2,ix3,2)-half*ps(igrid)%dx(ix1,&
ix2,ix3,2)) .and. location_slit<(ps(igrid)%x(ix1,ix2,ix3,&
2)+half*ps(igrid)%dx(ix1,ix2,ix3,2))) then
ixOnew=ix2
endif
enddo
ixOmin2=ixOnew
ixOmax2=ixOnew
else
do ix3=ixOmin3,ixOmax3
if (location_slit>=(ps(igrid)%x(ix1,ix2,ix3,3)-half*ps(igrid)%dx(ix1,&
ix2,ix3,3)) .and. location_slit<(ps(igrid)%x(ix1,ix2,ix3,&
3)+half*ps(igrid)%dx(ix1,ix2,ix3,3))) then
ixOnew=ix3
endif
enddo
ixOmin3=ixOnew
ixOmax3=ixOnew
endif
call get_EUV(spectrum_wl,ixImin1,ixImin2,ixImin3,ixImax1,ixImax2,ixImax3,&
ixOmin1,ixOmin2,ixOmin3,ixOmax1,ixOmax2,ixOmax3,ps(igrid)%w,&
ps(igrid)%x,fl,flux)
flux(ixOmin1:ixOmax1,ixOmin2:ixOmax2,&
ixOmin3:ixOmax3)=flux(ixOmin1:ixOmax1,ixOmin2:ixOmax2,&
ixOmin3:ixOmax3)/instrument_resolution_factor**2 !adjust flux due to artifical change of resolution
call fl%get_rho(ps(igrid)%w,ps(igrid)%x,ixImin1,ixImin2,ixImin3,ixImax1,&
ixImax2,ixImax3,ixOmin1,ixOmin2,ixOmin3,ixOmax1,ixOmax2,ixOmax3,rho)
v(ixOmin1:ixOmax1,ixOmin2:ixOmax2,ixOmin3:ixOmax3)=-&
ps(igrid)%w(ixOmin1:ixOmax1,ixOmin2:ixOmax2,ixOmin3:ixOmax3,&
iw_mom(direction_LOS))/rho(ixOmin1:ixOmax1,ixOmin2:ixOmax2,&
ixOmin3:ixOmax3)
call fl%get_pthermal(ps(igrid)%w,ps(igrid)%x,ixImin1,ixImin2,ixImin3,&
ixImax1,ixImax2,ixImax3,ixOmin1,ixOmin2,ixOmin3,ixOmax1,ixOmax2,&
ixOmax3,pth)
call fl%get_var_Rfactor(ps(igrid)%w,ps(igrid)%x,ixImin1,ixImin2,ixImin3,&
ixImax1,ixImax2,ixImax3,ixOmin1,ixOmin2,ixOmin3,ixOmax1,ixOmax2,&
ixOmax3,Te)
Te(ixOmin1:ixOmax1,ixOmin2:ixOmax2,ixOmin3:ixOmax3)=pth(ixOmin1:ixOmax1,&
ixOmin2:ixOmax2,ixOmin3:ixOmax3)/(Te(ixOmin1:ixOmax1,ixOmin2:ixOmax2,&
ixOmin3:ixOmax3)*rho(ixOmin1:ixOmax1,ixOmin2:ixOmax2,&
ixOmin3:ixOmax3))
! grid parameters
levelg=ps(igrid)%level
rft=2**(refine_max_level-levelg)
do ix1=ixOmin1,ixOmax1
do ix2=ixOmin2,ixOmax2
do ix3=ixOmin3,ixOmax3
if (flux(ix1,ix2,ix3)>smalldouble) then
if (SI_unit) then
wlc=lineCent*(1.d0+v(ix1,ix2,ix3)*unit_velocity*1.d2/const_c)
else
wlc=lineCent*(1.d0+v(ix1,ix2,ix3)*unit_velocity/const_c)
endif
wlwd=sqrt(kb_cgs*Te(ix1,ix2,ix3)*unit_temperature/(mass*mp_cgs))
wlwd=wlwd*lineCent/const_c
! involved pixel
select case(direction_slit)
case(1)
ixSmin=(block_nx1*(node(pig1_,igrid)-1)+(ix1-ixOmin1))*rft+1
ixSmax=(block_nx1*(node(pig1_,igrid)-1)+(ix1-ixOmin1+1))*rft
case(2)
ixSmin=(block_nx2*(node(pig2_,igrid)-1)+(ix2-ixOmin2))*rft+1
ixSmax=(block_nx2*(node(pig2_,igrid)-1)+(ix2-ixOmin2+1))*rft
case(3)
ixSmin=(block_nx3*(node(pig3_,igrid)-1)+(ix3-ixOmin3))*rft+1
ixSmax=(block_nx3*(node(pig3_,igrid)-1)+(ix3-ixOmin3+1))*rft
end select
! LOS depth
select case(direction_LOS)
case(1)
dL=ps(igrid)%dx(ix1,ix2,ix3,1)*unit_length
case(2)
dL=ps(igrid)%dx(ix1,ix2,ix3,2)*unit_length
case default
dL=ps(igrid)%dx(ix1,ix2,ix3,3)*unit_length
end select
if (SI_unit) dL=dL*1.d2
! integral pixel flux
do iwL=1,numWL
flux_pix=flux(ix1,ix2,ix3)*wlRsl*dL*exp(-(wL(iwL)-&
wlc)**2/(2*wlwd**2))/(sqrt(2*dpi)*wlwd)
if (flux_pix>smalldouble) then
flux_pix=flux_pix*wslit/spaceRsl
spectra(iwL,ixSmin:ixSmax)=spectra(iwL,ixSmin:ixSmax)+flux_pix
endif
enddo
endif
enddo
enddo
enddo
deallocate(flux,v,pth,Te,rho)
end subroutine integrate_spectra_datresol