write_spectra Subroutine

public subroutine write_spectra()

Writes Spectra to a HDF5 file

Arguments

None

Calls

proc~~write_spectra~~CallsGraph proc~write_spectra write_spectra interface~h5ltmake_compressed_dataset_double_f h5ltmake_compressed_dataset_double_f proc~write_spectra->interface~h5ltmake_compressed_dataset_double_f h5gcreate_f h5gcreate_f proc~write_spectra->h5gcreate_f h5fcreate_f h5fcreate_f proc~write_spectra->h5fcreate_f h5ltmake_dataset_int_f h5ltmake_dataset_int_f proc~write_spectra->h5ltmake_dataset_int_f h5open_f h5open_f proc~write_spectra->h5open_f h5fclose_f h5fclose_f proc~write_spectra->h5fclose_f h5gclose_f h5gclose_f proc~write_spectra->h5gclose_f h5close_f h5close_f proc~write_spectra->h5close_f h5ltset_attribute_string_f h5ltset_attribute_string_f proc~write_spectra->h5ltset_attribute_string_f proc~h5ltmake_compressed_dataset_double_f_4 h5ltmake_compressed_dataset_double_f_4 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_4 proc~h5ltmake_compressed_dataset_double_f_2 h5ltmake_compressed_dataset_double_f_2 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_2 proc~h5ltmake_compressed_dataset_double_f_1 h5ltmake_compressed_dataset_double_f_1 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_1 proc~h5ltmake_compressed_dataset_double_f_5 h5ltmake_compressed_dataset_double_f_5 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_5 proc~h5ltmake_compressed_dataset_double_f_3 h5ltmake_compressed_dataset_double_f_3 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_3 proc~h5ltmake_compressed_dataset_double_f_6 h5ltmake_compressed_dataset_double_f_6 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_6 proc~h5ltmake_compressed_dataset_double_f_7 h5ltmake_compressed_dataset_double_f_7 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_7 proc~chunk_size chunk_size proc~h5ltmake_compressed_dataset_double_f_4->proc~chunk_size h5sclose_f h5sclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5sclose_f h5ltmake_dataset_double_f h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_4->h5ltmake_dataset_double_f h5pset_chunk_f h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_4->h5pset_chunk_f h5pcreate_f h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_4->h5pcreate_f h5pset_shuffle_f h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_4->h5pset_shuffle_f h5pset_deflate_f h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_4->h5pset_deflate_f h5dcreate_f h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_4->h5dcreate_f h5dwrite_f h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_4->h5dwrite_f h5dclose_f h5dclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5dclose_f h5pclose_f h5pclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5pclose_f h5screate_simple_f h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_4->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_2->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_2->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_2->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_2->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_2->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_2->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_2->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_2->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_2->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_2->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_2->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_2->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_1->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_1->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_1->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_1->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_1->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_1->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_1->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_1->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_1->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_1->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_1->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_1->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_5->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_5->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_5->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_5->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_5->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_5->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_5->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_5->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_5->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_5->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_5->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_5->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_3->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_3->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_3->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_3->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_3->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_3->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_3->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_3->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_3->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_3->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_3->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_3->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_6->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_6->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_6->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_6->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_6->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_6->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_6->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_6->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_6->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_6->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_6->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_6->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_7->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_7->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_7->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_7->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_7->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_7->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_7->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_7->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_7->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_7->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_7->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_7->h5screate_simple_f

Called by

proc~~write_spectra~~CalledByGraph proc~write_spectra write_spectra program~fidasim fidasim program~fidasim->proc~write_spectra

Contents

Source Code


Source Code

subroutine write_spectra
    !+ Writes [[libfida:spectra]] to a HDF5 file
    integer(HID_T) :: fid, gid
    integer(HSIZE_T), dimension(4) :: dims
    integer(HSIZE_T), dimension(5) :: dims_stokes
    integer(HSIZE_T), dimension(1) :: d
    integer :: error

    character(charlim) :: filename
    integer :: i, ir, ic
    real(Float64) :: factor
    real(Float64), dimension(:), allocatable :: lambda_arr

    real(Float64), dimension(inputs%nlambda,spec_chords%nchan) :: full, half, third
    real(Float64), dimension(4,inputs%nlambda,spec_chords%nchan) :: fullstokes, halfstokes, thirdstokes
    real(Float64), dimension(inputs%nlambda,spec_chords%nchan,n_thermal) :: dcx, halo, cold
    real(Float64), dimension(4,inputs%nlambda,spec_chords%nchan,n_thermal) :: dcxstokes, halostokes, coldstokes
    real(Float64), dimension(inputs%nlambda,spec_chords%nchan,particles%nclass) :: fida, pfida
    real(Float64), dimension(4,inputs%nlambda,spec_chords%nchan,particles%nclass) :: fidastokes, pfidastokes

    real(Float64), dimension(3,reservoir_size,spec_chords%nchan) :: dcx_spat, halo_spat, fida_spat, pfida_spat
    real(Float64), dimension(reservoir_size,spec_chords%nchan) :: dcx_photons, halo_photons, fida_photons, pfida_photons

    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)]!
    factor = 1.d0/(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)
    if(inputs%stark_components.ge.1) then
        call h5ltmake_dataset_int_f(fid, "/nstark", 0, d, [n_stark], error)
    endif
    dims(1) = n_stark
    dims(2) = inputs%nlambda
    dims(3) = spec_chords%nchan
    dims(4) = n_thermal
    dims_stokes(1) = n_stark
    dims_stokes(2) = 4
    dims_stokes(3) = inputs%nlambda
    dims_stokes(4) = spec_chords%nchan
    dims_stokes(5) = n_thermal
    call h5ltmake_compressed_dataset_double_f(fid, "/lambda", 1, dims(2:2), &
         lambda_arr, error)
    call h5ltmake_compressed_dataset_double_f(fid, "/radius", 1, dims(3:3), &
         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
        spec%brems = factor*spec%brems
        !Write variables
        call h5ltmake_compressed_dataset_double_f(fid, "/brems", 2, &
             dims(2:3), 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
        spec%full = factor*spec%full
        spec%fullstokes = factor*spec%fullstokes
        spec%half = factor*spec%half
        spec%halfstokes = factor*spec%halfstokes
        spec%third = factor*spec%third
        spec%thirdstokes = factor*spec%thirdstokes
        if (inputs%stark_components.eq.0) then
            full  = sum(spec%full, dim=1)
            half  = sum(spec%half, dim=1)
            third = sum(spec%third, dim=1)
            fullstokes = sum(spec%fullstokes,dim=1)
            halfstokes = sum(spec%halfstokes, dim=1)
            thirdstokes = sum(spec%thirdstokes, dim=1)
            !Write variables
            call h5ltmake_compressed_dataset_double_f(fid, "/full", 2, dims(2:3), &
                 full, error)
            call h5ltmake_compressed_dataset_double_f(fid, "/half", 2, dims(2:3), &
                 half, error)
            call h5ltmake_compressed_dataset_double_f(fid, "/third", 2, dims(2:3),&
                 third, error)
            call h5ltmake_compressed_dataset_double_f(fid, "/fullstokes", 3, dims_stokes(2:4), &
                 fullstokes, error)
            call h5ltmake_compressed_dataset_double_f(fid, "/halfstokes", 3, dims_stokes(2:4), &
                 halfstokes, error)
            call h5ltmake_compressed_dataset_double_f(fid, "/thirdstokes", 3, dims_stokes(2:4),&
                 thirdstokes, 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,"/fullstokes","description", &
                 "Full energy component of the beam emmision stokes parameters: &
                 &fullstokes(4,lambda,chan)", error)
            call h5ltset_attribute_string_f(fid,"/fullstokes","units","Ph/(s*nm*sr*m^2)",error )
            call h5ltset_attribute_string_f(fid,"/halfstokes","description", &
                 "Half energy component of the beam emmision stokes parameters: &
                 &halfstokes(4,lambda,chan)", error)
            call h5ltset_attribute_string_f(fid,"/halfstokes","units","Ph/(s*nm*sr*m^2)",error )
            call h5ltset_attribute_string_f(fid,"/thirdstokes","description", &
                 "Third energy component of the beam emmision stokes parameters: &
                 &thirdstokes(4,lambda,chan)", error)
            call h5ltset_attribute_string_f(fid,"/thirdstokes","units","Ph/(s*nm*sr*m^2)",error )
        else
            call h5ltmake_compressed_dataset_double_f(fid, "/full", 3, dims(1:3), &
                 spec%full, error)
            call h5ltmake_compressed_dataset_double_f(fid, "/half", 3, dims(1:3), &
                 spec%half, error)
            call h5ltmake_compressed_dataset_double_f(fid, "/third", 3, dims(1:3),&
                 spec%third, error)
            call h5ltmake_compressed_dataset_double_f(fid, "/fullstokes", 4, dims_stokes(1:4), &
                spec%fullstokes, error)
            call h5ltmake_compressed_dataset_double_f(fid, "/halfstokes", 4, dims_stokes(1:4), &
                spec%halfstokes, error)
            call h5ltmake_compressed_dataset_double_f(fid, "/thirdstokes", 4, dims_stokes(1:4),&
                spec%thirdstokes, error)
            call h5ltset_attribute_string_f(fid,"/full","description", &
                 "Full energy component of the beam emmision stark components: full(stark,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 stark components: half(stark,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 stark components: third(stark,lambda,chan)", &
                 error)
            call h5ltset_attribute_string_f(fid,"/third","units","Ph/(s*nm*sr*m^2)",error )
            call h5ltset_attribute_string_f(fid,"/fullstokes","description", &
                 "Full energy component of the beam emmision stark components stokes parameters: &
                 &fullstokes(stark,4,lambda,chan)", &
                 error)
            call h5ltset_attribute_string_f(fid,"/fullstokes","units","Ph/(s*nm*sr*m^2)",error )
            call h5ltset_attribute_string_f(fid,"/halfstokes","description", &
                 "Half energy component of the beam emmision stark components stokes parameters: &
                 &halfstokes(stark,4,lambda,chan)", &
                 error)
            call h5ltset_attribute_string_f(fid,"/halfstokes","units","Ph/(s*nm*sr*m^2)",error )
            call h5ltset_attribute_string_f(fid,"/thirdstokes","description", &
                 "Third energy component of the beam emmision stark components stokes parameters: &
                 &thirdstokes(stark,4,lambda,chan)", &
                 error)
            call h5ltset_attribute_string_f(fid,"/thirdstokes","units","Ph/(s*nm*sr*m^2)",error )
        endif
    endif

    if(inputs%calc_dcx.ge.1) then
        spec%dcx = factor*spec%dcx
        spec%dcxstokes = factor*spec%dcxstokes
        if (inputs%stark_components.eq.0) then
            dcx = sum(spec%dcx, dim=1)
            dcxstokes = sum(spec%dcxstokes, dim=1)
            call h5ltmake_compressed_dataset_double_f(fid, "/dcx", 3, dims(2:4), &
                 dcx, error)
            call h5ltset_attribute_string_f(fid,"/dcx","description", &
                 "Direct Charge Exchange (DCX) emission: dcx(lambda,chan,species)", error)
            call h5ltset_attribute_string_f(fid,"/dcx","units","Ph/(s*nm*sr*m^2)",error )
            call h5ltmake_compressed_dataset_double_f(fid,"/dcxstokes", 4, dims_stokes(2:5),&
                 dcxstokes, error)
            call h5ltset_attribute_string_f(fid,"/dcxstokes","description", &
                 "Direct Charge Exchange (DCX) emission stokes parameters: &
                 &dcx(4,lambda,chan,species)", error)
            call h5ltset_attribute_string_f(fid,"/dcxstokes","units","Ph/(s*nm*sr*m^2)",error )
        else
            call h5ltmake_compressed_dataset_double_f(fid, "/dcx", 4, dims, &
                 spec%dcx, error)
            call h5ltmake_compressed_dataset_double_f(fid, "/dcxstokes", 5, dims_stokes, &
                 spec%dcxstokes, error)
            call h5ltset_attribute_string_f(fid,"/dcx","description", &
                 "Direct Charge Exchange (DCX) emission stark components: dcx(stark,lambda,chan,species)", error)
            call h5ltset_attribute_string_f(fid,"/dcx","units","Ph/(s*nm*sr*m^2)",error )
            call h5ltset_attribute_string_f(fid,"/dcxstokes","description", &
                 "Direct Charge Exchange (DCX) emission stark components stokes parameters: &
                 &dcxstokes(stark,4,lambda,chan,species)", error)
            call h5ltset_attribute_string_f(fid,"/dcxstokes","units","Ph/(s*nm*sr*m^2)",error )
        endif
    endif

    if(inputs%calc_halo.ge.1) then
        spec%halo = factor*spec%halo
        spec%halostokes = factor*spec%halostokes
        if (inputs%stark_components.eq.0) then
            halo = sum(spec%halo, dim=1)
            halostokes = sum(spec%halostokes, dim=1)
            call h5ltmake_compressed_dataset_double_f(fid, "/halo", 3, dims(2:4), &
                 halo, error)
            call h5ltset_attribute_string_f(fid,"/halo","description", &
                 "Halo component of the beam emmision: halo(lambda,chan,species)", error)
            call h5ltset_attribute_string_f(fid,"/halo","units","Ph/(s*nm*sr*m^2)",error )
            call h5ltmake_compressed_dataset_double_f(fid, "/halostokes", 4, dims_stokes(2:5), &
                 halostokes, error)
            call h5ltset_attribute_string_f(fid,"/halostokes","description", &
                 "Halo component of the beam emmision stokes parameters: halostokes(4,lambda,chan,species)", error)
            call h5ltset_attribute_string_f(fid,"/halostokes","units","Ph/(s*nm*sr*m^2)",error )
        else
            call h5ltmake_compressed_dataset_double_f(fid, "/halo", 4, dims, &
                 spec%halo, error)
            call h5ltset_attribute_string_f(fid,"/halo","description", &
                 "Halo component of the beam emmision stark components: halo(stark,lambda,chan,species)", error)
            call h5ltset_attribute_string_f(fid,"/halo","units","Ph/(s*nm*sr*m^2)",error )
            call h5ltmake_compressed_dataset_double_f(fid, "/halostokes", 5, dims_stokes(1:5), &
                 spec%halostokes, error)
            call h5ltset_attribute_string_f(fid,"/halostokes","description", &
                 "Halo component of the beam emmision stark components stokes parameters: &
                 &halostokes(stark,4,lambda,chan,species)", error)
            call h5ltset_attribute_string_f(fid,"/halostokes","units","Ph/(s*nm*sr*m^2)",error )
        endif
    endif

    if((inputs%calc_halo.ge.1) .or. (inputs%calc_dcx.ge.1)) then
       call h5ltmake_compressed_dataset_double_f(fid, "/thermal_mass", 1, dims(4:4), thermal_mass(1:n_thermal), error)
       call h5ltset_attribute_string_f(fid,"/thermal_mass","description", &
            "Mass of each of the thermal ions: thermal_mass(species)", error)
       call h5ltset_attribute_string_f(fid,"/thermal_mass", "units", "amu", error)

       call h5ltmake_compressed_dataset_double_f(fid, "/thermal_lambda0", 1, dims(4:4), thermal_lambda0(1:n_thermal), error)
       call h5ltset_attribute_string_f(fid,"/thermal_lambda0","description", &
            "Rest wavelength of the thermal species lines: thermal_lambda0(species)", error)
       call h5ltset_attribute_string_f(fid,"/thermal_lambda0", "units", "nm", error)
    endif


    if(inputs%calc_cold.ge.1) then
        spec%cold = factor*spec%cold
        spec%coldstokes = factor*spec%coldstokes
        if (inputs%stark_components.eq.0) then
            cold = sum(spec%cold, dim=1)
            coldstokes = sum(spec%coldstokes, dim=1)
            call h5ltmake_compressed_dataset_double_f(fid, "/cold", 3, dims(2:4), &
                 cold, error)
            call h5ltset_attribute_string_f(fid,"/cold","description", &
                 "Cold D-alpha emission: cold(lambda,chan,species)", error)
            call h5ltset_attribute_string_f(fid,"/cold","units","Ph/(s*nm*sr*m^2)",error )
            call h5ltmake_compressed_dataset_double_f(fid, "/coldstokes", 4, dims_stokes(2:5), &
                 coldstokes, error)
            call h5ltset_attribute_string_f(fid,"/coldstokes","description", &
                 "Cold D-alpha emission stokes parameters: coldstokes(4,lambda,chan,species)", error)
            call h5ltset_attribute_string_f(fid,"/coldstokes","units","Ph/(s*nm*sr*m^2)",error )
        else
            call h5ltmake_compressed_dataset_double_f(fid, "/cold", 4, dims, &
                 spec%cold, error)
            call h5ltset_attribute_string_f(fid,"/cold","description", &
                 "Cold D-alpha emission stark components: cold(stark,lambda,chan,species)", &
                 error)
            call h5ltset_attribute_string_f(fid,"/cold","units","Ph/(s*nm*sr*m^2)",error )
            call h5ltmake_compressed_dataset_double_f(fid, "/coldstokes", 5, dims_stokes, &
                 spec%coldstokes, error)
            call h5ltset_attribute_string_f(fid,"/coldstokes","description", &
                 "Cold D-alpha emission stokes parameters: coldstokes(stark,4,lambda,chan,species)", error)
            call h5ltset_attribute_string_f(fid,"/coldstokes","units","Ph/(s*nm*sr*m^2)",error )
        endif
    endif

    dims(4) = particles%nclass
    if(inputs%calc_fida.ge.1) then
        spec%fida = factor*spec%fida
        spec%fidastokes = factor*spec%fidastokes
        if(inputs%stark_components.eq.0) then
            fida = sum(spec%fida, dim=1)
            fidastokes = sum(spec%fidastokes, dim=1)
        endif
        !Write variables
        if(particles%nclass.le.1) then
            if (inputs%stark_components.eq.0) then
                call h5ltmake_compressed_dataset_double_f(fid, "/fida", 2, &
                     dims(2:3), fida(:,:,1), error)
                call h5ltmake_compressed_dataset_double_f(fid, "/fidastokes", 3, &
                     dims_stokes(2:4), fidastokes(:,:,:,1), error)
                !Add attributes
                call h5ltset_attribute_string_f(fid,"/fida","description", &
                     "Active Fast-ion D-alpha (FIDA) emmision: fida(lambda,chan)", error)
                call h5ltset_attribute_string_f(fid,"/fidastokes","description", &
                     "Active Fast-ion D-alpha (FIDA) emmision stokes parameters: &
                     &fidastokes(4,lambda,chan)", error)
            else
                call h5ltmake_compressed_dataset_double_f(fid, "/fida", 3, &
                     dims(1:3), spec%fida(:,:,:,1), error)
                call h5ltmake_compressed_dataset_double_f(fid, "/fida", 4, &
                     dims_stokes(1:4), spec%fidastokes(:,:,:,:,1), error)
                !Add attributes
                call h5ltset_attribute_string_f(fid,"/fida","description", &
                     "Active Fast-ion D-alpha (FIDA) emmision stark components: fida(stark,lambda,chan)", &
                     error)
                call h5ltset_attribute_string_f(fid,"/fidastokes","description", &
                     "Active Fast-ion D-alpha (FIDA) emmision stark components stokes parameters: &
                     &fidastokes(stark,4,lambda,chan)", &
                     error)
            endif
        else
            if (inputs%stark_components.eq.0) then
                call h5ltmake_dataset_int_f(fid,"/nclass", 0, d, [particles%nclass], error)
                call h5ltmake_compressed_dataset_double_f(fid, "/fida", 3, &
                     dims(2:4), fida, error)
                call h5ltmake_compressed_dataset_double_f(fid, "/fidastokes", 4, &
                     dims_stokes(2:5), fidastokes, error)
                !Add attributes
                call h5ltset_attribute_string_f(fid,"/fida","description", &
                     "Active Fast-ion D-alpha (FIDA) emmision: fida(lambda,chan,class)", error)
                call h5ltset_attribute_string_f(fid,"/fidastokes","description", &
                     "Active Fast-ion D-alpha (FIDA) emmision stokes parameters: &
                     &fidastokes(4,lambda,chan,class)", error)
            else
                call h5ltmake_compressed_dataset_double_f(fid, "/fida", 4, &
                     dims, spec%fida, error)
                call h5ltmake_compressed_dataset_double_f(fid, "/fidastokes", 5, &
                     dims_stokes, spec%fidastokes, error)
                !Add attributes
                call h5ltset_attribute_string_f(fid,"/fida","description", &
                     "Active Fast-ion D-alpha (FIDA) emmision stark components: fida(stark,lambda,chan,class)", &
                     error)
                call h5ltset_attribute_string_f(fid,"/fidastokes","description", &
                     "Active Fast-ion D-alpha (FIDA) emmision stark components stokes parameters: &
                     &fidastokes(stark,4,lambda,chan,class)", &
                     error)
            endif
       endif
       call h5ltset_attribute_string_f(fid,"/fida","units","Ph/(s*nm*sr*m^2)",error )
       call h5ltset_attribute_string_f(fid,"/fidastokes","units","Ph/(s*nm*sr*m^2)",error )
    endif

    if(inputs%calc_pfida.ge.1) then
        spec%pfida = factor*spec%pfida
        spec%pfidastokes = factor*spec%pfidastokes
        if (inputs%stark_components.eq.0) then
            pfida = sum(spec%pfida, dim=1)
            pfidastokes = sum(spec%pfidastokes, dim=1)
        endif
        !Write variables
        if(particles%nclass.le.1) then
            if (inputs%stark_components.eq.0) then
                call h5ltmake_compressed_dataset_double_f(fid, "/pfida", 2, &
                     dims(2:3), pfida(:,:,1), error)
                call h5ltmake_compressed_dataset_double_f(fid, "/pfidastokes", 3, &
                     dims_stokes(2:4), pfidastokes(:,:,:,1), error)
                !Add attributes
                call h5ltset_attribute_string_f(fid,"/pfida","description", &
                     "Passive Fast-ion D-alpha (p-FIDA) emmision: pfida(lambda,chan)", error)
                call h5ltset_attribute_string_f(fid,"/pfidastokes","description", &
                     "Passive Fast-ion D-alpha (p-FIDA) emmision stokes parameters: &
                     &pfidastokes(4,lambda,chan)", error)
            else
                call h5ltmake_compressed_dataset_double_f(fid, "/pfida", 3, &
                     dims(1:3), spec%pfida(:,:,:,1), error)
                call h5ltmake_compressed_dataset_double_f(fid, "/pfidastokes", 4, &
                     dims_stokes(1:4), spec%pfidastokes(:,:,:,:,1), error)
                !Add attributes
                call h5ltset_attribute_string_f(fid,"/pfida","description", &
                     "Passive Fast-ion D-alpha (p-FIDA) emmision stark components: pfida(stark,lambda,chan)", &
                     error)
                call h5ltset_attribute_string_f(fid,"/pfidastokes","description", &
                     "Passive Fast-ion D-alpha (p-FIDA) emmision stokes parameters: &
                     &pfidastokes(stark,4,lambda,chan)", error)
            endif
        else
            if(inputs%calc_fida.le.0) then
                call h5ltmake_dataset_int_f(fid,"/nclass", 0, d, [particles%nclass], error)
            endif
            if (inputs%stark_components.eq.0) then
                call h5ltmake_compressed_dataset_double_f(fid, "/pfida", 3, &
                     dims(2:4), pfida, error)
                call h5ltmake_compressed_dataset_double_f(fid, "/pfidastokes", 4, &
                     dims_stokes(2:5), pfidastokes, error)
                !Add attributes
                call h5ltset_attribute_string_f(fid,"/pfida","description", &
                     "Passive Fast-ion D-alpha (p-FIDA) emmision: pfida(lambda,chan,class)", error)
                call h5ltset_attribute_string_f(fid,"/pfidastokes","description", &
                     "Passive Fast-ion D-alpha (p-FIDA) emmision stokes parameters: pfidastokes(4,lambda,chan, class)", error)
            else
                call h5ltmake_compressed_dataset_double_f(fid, "/pfida", 4, &
                     dims, spec%pfida, error)
                call h5ltmake_compressed_dataset_double_f(fid, "/pfidastokes", 5, &
                     dims_stokes, spec%pfidastokes, error)
                !Add attributes
                call h5ltset_attribute_string_f(fid,"/pfida","description", &
                     "Passive Fast-ion D-alpha (p-FIDA) emmision stark components: pfida(stark,lambda,chan,class)", &
                     error)
                call h5ltset_attribute_string_f(fid,"/pfidastokes","description", &
                     "Passive Fast-ion D-alpha (p-FIDA) emmision stokes parameters: pfidastokes(stark,4,lambda,chan, class)", error)
            endif
       endif
       call h5ltset_attribute_string_f(fid,"/pfida","units","Ph/(s*nm*sr*m^2)",error )
       call h5ltset_attribute_string_f(fid,"/pfidastokes","units","Ph/(s*nm*sr*m^2)",error )
    endif

    if(inputs%calc_res.ge.1) then 
       !Create spatial group
       call h5gcreate_f(fid, "spatial", gid, error)
       call h5ltset_attribute_string_f(fid,"spatial","description", &
            "Spatial Resolution",error)
       dims(1) = 3
       dims(2) = reservoir_size
       dims(3) = spec_chords%nchan
       if(inputs%calc_dcx.ge.1) then
           do ic=1,spec_chords%nchan
               do ir=1,reservoir_size
                   dcx_photons(ir,ic) = spatres%dcx(ic)%R(ir)%w
                   dcx_spat(:,ir,ic) = spatres%dcx(ic)%R(ir)%v
               enddo
           enddo
           call h5ltmake_compressed_dataset_double_f(gid, "dcx_spatial", 3, dims(1:3),&
                dcx_spat, error)
           call h5ltset_attribute_string_f(gid,"dcx_spatial","units","cm",error)
           call h5ltset_attribute_string_f(gid,"dcx_spatial","description",&
                "Birth position of neutral that produced a photon: dcx_spatial([x,y,z],sample,channel)",error)
           call h5ltmake_compressed_dataset_double_f(gid, "dcx_photons", 2, dims(2:3),&
                dcx_photons, error)
           call h5ltset_attribute_string_f(gid,"dcx_photons","units","Ph",error)
           call h5ltset_attribute_string_f(gid,"dcx_photons","description",&
                "Number of photons produced by neutral: dcx_photons(sample,channel)",error)
       endif
       if(inputs%calc_halo.ge.1) then
           do ic=1,spec_chords%nchan
               do ir=1,reservoir_size
                   halo_photons(ir,ic) = spatres%halo(ic)%R(ir)%w
                   halo_spat(:,ir,ic) = spatres%halo(ic)%R(ir)%v
               enddo
           enddo
           call h5ltmake_compressed_dataset_double_f(gid, "halo_spatial", 3, dims(1:3),&
                halo_spat, error)
           call h5ltset_attribute_string_f(gid,"halo_spatial","units","cm",error)
           call h5ltset_attribute_string_f(gid,"halo_spatial","description",&
                "Birth position of neutral that produced a photon: halo_spatial([x,y,z],sample,channel)",error)
           call h5ltmake_compressed_dataset_double_f(gid, "halo_photons", 2, dims(2:3),&
                halo_photons, error)
           call h5ltset_attribute_string_f(gid,"halo_photons","units","Ph",error)
           call h5ltset_attribute_string_f(gid,"halo_photons","description",&
                "Number of photons produced by neutral: halo_photons(sample,channel)",error)
       endif
       if(inputs%calc_fida.ge.1) then
           do ic=1,spec_chords%nchan
               do ir=1,reservoir_size
                   fida_photons(ir,ic) = spatres%fida(ic)%R(ir)%w
                   fida_spat(:,ir,ic) = spatres%fida(ic)%R(ir)%v
               enddo
           enddo
           call h5ltmake_compressed_dataset_double_f(gid, "fida_spatial", 3, dims(1:3),&
                fida_spat, error)
           call h5ltset_attribute_string_f(gid,"fida_spatial","units","cm",error)
           call h5ltset_attribute_string_f(gid,"fida_spatial","description",&
                "Birth position of neutral that produced a photon: fida_spatial([x,y,z],sample,channel)",error)
           call h5ltmake_compressed_dataset_double_f(gid, "fida_photons", 2, dims(2:3),&
                fida_photons, error)
           call h5ltset_attribute_string_f(gid,"fida_photons","units","Ph",error)
           call h5ltset_attribute_string_f(gid,"fida_photons","description",&
                "Number of photons produced by neutral: fida_photons(sample,channel)",error)
       endif
       if(inputs%calc_pfida.ge.1) then
           do ic=1,spec_chords%nchan
               do ir=1,reservoir_size
                   pfida_photons(ir,ic) = spatres%pfida(ic)%R(ir)%w
                   pfida_spat(:,ir,ic) = spatres%pfida(ic)%R(ir)%v
               enddo
           enddo
           call h5ltmake_compressed_dataset_double_f(gid, "pfida_spatial", 3, dims(1:3),&
                pfida_spat, error)
           call h5ltset_attribute_string_f(gid,"pfida_spatial","units","cm",error)
           call h5ltset_attribute_string_f(gid,"pfida_spatial","description",&
                "Birth position of neutral that produced a photon: pfida_spatial([x,y,z],sample,channel)",error)
           call h5ltmake_compressed_dataset_double_f(gid, "pfida_photons", 2, dims(2:3),&
                pfida_photons, error)
           call h5ltset_attribute_string_f(gid,"pfida_photons","units","Ph",error)
           call h5ltset_attribute_string_f(gid,"pfida_photons","description",&
                "Number of photons produced by neutral: pfida_photons(sample,channel)",error)
       endif

       !Close spatial group
       call h5gclose_f(gid, 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