Writes Spectra to a HDF5 file
subroutine write_spectra
!+ Writes [[libfida:spectra]] to a HDF5 file
integer(HID_T) :: fid
integer(HSIZE_T), dimension(3) :: dims
integer(HSIZE_T), dimension(1) :: d
integer :: error
character(charlim) :: filename
integer :: i
real(Float64), dimension(:), allocatable :: lambda_arr
allocate(lambda_arr(inputs%nlambda))
do i=1,inputs%nlambda
lambda_arr(i) = (i-0.5)*inputs%dlambda + inputs%lambdamin
enddo
!! convert [Ph/(s*wavel_bin*cm^2*all_directions)] to [Ph/(s*nm*sr*m^2)]!
spec%brems=spec%brems/(inputs%dlambda)/(4.d0*pi)*1.d4
spec%bes=spec%bes/(inputs%dlambda)/(4.d0*pi)*1.d4
spec%fida=spec%fida/(inputs%dlambda)/(4.d0*pi)*1.d4
!! write to file
filename=trim(adjustl(inputs%result_dir))//"/"//trim(adjustl(inputs%runid))//"_spectra.h5"
!Open HDF5 interface
call h5open_f(error)
!Create file overwriting any existing file
call h5fcreate_f(filename, H5F_ACC_TRUNC_F, fid, error)
!Write variables
d(1) = 1
call h5ltmake_dataset_int_f(fid, "/nchan", 0, d, [spec_chords%nchan], error)
call h5ltmake_dataset_int_f(fid, "/nlambda", 0, d, [inputs%nlambda], error)
dims(1) = inputs%nlambda
dims(2) = spec_chords%nchan
dims(3) = particles%nclass
call h5ltmake_compressed_dataset_double_f(fid, "/lambda", 1, dims(1:1), &
lambda_arr, error)
call h5ltmake_compressed_dataset_double_f(fid, "/radius", 1, dims(2:2), &
spec_chords%radius, error)
!Add attributes
call h5ltset_attribute_string_f(fid,"/nchan", "description", &
"Number of channels", error)
call h5ltset_attribute_string_f(fid,"/nlambda", "description", &
"Number of wavelengths", error)
call h5ltset_attribute_string_f(fid,"/lambda","description", &
"Wavelength array", error)
call h5ltset_attribute_string_f(fid,"/lambda","units","nm", error)
call h5ltset_attribute_string_f(fid,"/radius", "description", &
"Line of sight radius at midplane or tangency point", error)
call h5ltset_attribute_string_f(fid,"/radius","units","cm", error)
if(inputs%calc_brems.ge.1) then
!Write variables
call h5ltmake_compressed_dataset_double_f(fid, "/brems", 2, &
dims(1:2), spec%brems, error)
!Add attributes
call h5ltset_attribute_string_f(fid,"/brems","description", &
"Visible Bremsstrahlung: brems(lambda,chan)", error)
call h5ltset_attribute_string_f(fid,"/brems","units",&
"Ph/(s*nm*sr*m^2)",error )
endif
if(inputs%calc_bes.ge.1) then
!Write variables
call h5ltmake_compressed_dataset_double_f(fid, "/full", 2, dims(1:2), &
spec%bes(:,:,nbif_type), error)
call h5ltmake_compressed_dataset_double_f(fid, "/half", 2, dims(1:2), &
spec%bes(:,:,nbih_type), error)
call h5ltmake_compressed_dataset_double_f(fid, "/third", 2, dims(1:2),&
spec%bes(:,:,nbit_type), error)
call h5ltmake_compressed_dataset_double_f(fid, "/halo", 2, dims(1:2), &
spec%bes(:,:,halo_type), error)
!Add attributes
call h5ltset_attribute_string_f(fid,"/full","description", &
"Full energy component of the beam emmision: full(lambda,chan)", error)
call h5ltset_attribute_string_f(fid,"/full","units","Ph/(s*nm*sr*m^2)",error )
call h5ltset_attribute_string_f(fid,"/half","description", &
"Half energy component of the beam emmision: half(lambda,chan)", error)
call h5ltset_attribute_string_f(fid,"/half","units","Ph/(s*nm*sr*m^2)",error )
call h5ltset_attribute_string_f(fid,"/third","description", &
"Third energy component of the beam emmision: third(lambda,chan)", error)
call h5ltset_attribute_string_f(fid,"/third","units","Ph/(s*nm*sr*m^2)",error )
call h5ltset_attribute_string_f(fid,"/halo","description", &
"Halo component of the beam emmision (includes dcx): halo(lambda,chan)", error)
call h5ltset_attribute_string_f(fid,"/halo","units","Ph/(s*nm*sr*m^2)",error )
endif
if(inputs%calc_fida.ge.1) then
!Write variables
if(particles%nclass.le.1) then
call h5ltmake_compressed_dataset_double_f(fid, "/fida", 2, &
dims(1:2), spec%fida(:,:,1), error)
!Add attributes
call h5ltset_attribute_string_f(fid,"/fida","description", &
"Fast-ion D-alpha (FIDA) emmision: fida(lambda,chan)", error)
else
call h5ltmake_dataset_int_f(fid,"/nclass", 0, d, [particles%nclass], error)
call h5ltmake_compressed_dataset_double_f(fid, "/fida", 3, &
dims, spec%fida, error)
!Add attributes
call h5ltset_attribute_string_f(fid,"/fida","description", &
"Fast-ion D-alpha (FIDA) emmision: fida(lambda,chan,class)", error)
endif
call h5ltset_attribute_string_f(fid,"/fida","units","Ph/(s*nm*sr*m^2)",error )
endif
call h5ltset_attribute_string_f(fid, "/", "version", version, error)
call h5ltset_attribute_string_f(fid,"/","description",&
"Spectra calculated by FIDASIM", error)
!Close file
call h5fclose_f(fid, error)
!Close HDF5 interface
call h5close_f(error)
if(inputs%verbose.ge.1) then
write(*,'(T4,a,a)') 'Spectra written to: ', trim(filename)
endif
end subroutine write_spectra