write_cfpd_weights Subroutine

public subroutine write_cfpd_weights()

Writes cfpd to a HDF5 file

Arguments

None

Calls

proc~~write_cfpd_weights~~CallsGraph proc~write_cfpd_weights write_cfpd_weights interface~h5ltmake_compressed_dataset_double_f h5ltmake_compressed_dataset_double_f proc~write_cfpd_weights->interface~h5ltmake_compressed_dataset_double_f h5fcreate_f h5fcreate_f proc~write_cfpd_weights->h5fcreate_f h5ltmake_dataset_int_f h5ltmake_dataset_int_f proc~write_cfpd_weights->h5ltmake_dataset_int_f h5fclose_f h5fclose_f proc~write_cfpd_weights->h5fclose_f h5open_f h5open_f proc~write_cfpd_weights->h5open_f h5close_f h5close_f proc~write_cfpd_weights->h5close_f h5ltset_attribute_string_f h5ltset_attribute_string_f proc~write_cfpd_weights->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_cfpd_weights~~CalledByGraph proc~write_cfpd_weights write_cfpd_weights proc~cfpd_f cfpd_f proc~cfpd_f->proc~write_cfpd_weights program~fidasim fidasim program~fidasim->proc~cfpd_f

Contents

Source Code


Source Code

subroutine write_cfpd_weights
    !+ Writes [[libfida:cfpd]] to a HDF5 file
    integer(HID_T) :: fid
    integer(HSIZE_T), dimension(1) :: dim1
    integer(HSIZE_T), dimension(2) :: dim2
    integer(HSIZE_T), dimension(4) :: dim4
    integer :: error

    character(charlim) :: filename

    !! write to file
    filename=trim(adjustl(inputs%result_dir))//"/"//trim(adjustl(inputs%runid))//"_cfpd.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
    if(inputs%dist_type.eq.1) then
        dim1(1) = 1
        dim2 = [ctable%nenergy, ctable%nchan]
        dim4 = [ctable%nenergy, ctable%nchan, fbm%nenergy, fbm%npitch]

        call h5ltmake_compressed_dataset_double_f(fid, "/flux", 2, dim2, cfpd%flux, error)
        call h5ltmake_compressed_dataset_double_f(fid, "/prob", 2, dim2, cfpd%prob, error)
        call h5ltmake_compressed_dataset_double_f(fid, "/gam", 2, dim2, cfpd%gam, error)
        call h5ltmake_compressed_dataset_double_f(fid, "/weight", 4, dim4, cfpd%weight, error)

        call h5ltmake_dataset_int_f(fid,"/nenergy",0,dim1,[fbm%nenergy], error)
        call h5ltmake_dataset_int_f(fid,"/npitch",0,dim1,[fbm%npitch], error)
        call h5ltmake_compressed_dataset_double_f(fid,"/energy", 1, dim4(3:3), fbm%energy, error)
        call h5ltmake_compressed_dataset_double_f(fid,"/pitch", 1, dim4(4:4), fbm%pitch, error)
        call h5ltmake_compressed_dataset_double_f(fid,"/earray", 1, dim2(1:1), ctable%earray, error)

        call h5ltset_attribute_string_f(fid,"/flux", "description", &
             "CFPD flux: flux(energy,chan)", error)
        call h5ltset_attribute_string_f(fid,"/flux", "units", &
             "kHz", error)

        call h5ltset_attribute_string_f(fid,"/prob", "description", &
             "CFPD average nonzero probability: prob(energy,chan)", error)
        call h5ltset_attribute_string_f(fid,"/gam", "description", &
             "CFPD average nonzero gyroangle: gam(energy,chan)", error)
        call h5ltset_attribute_string_f(fid,"/gam", "units", "rad", error)
        call h5ltset_attribute_string_f(fid,"/weight", "description", &
             "CFPD Weight Function: weight(Ch,E3,E,p), rate = sum(f*weight)", error)
        call h5ltset_attribute_string_f(fid,"/weight", "units","products*cm^3*dE*dp/fast-ion*s", error)

        call h5ltset_attribute_string_f(fid,"/nenergy", "description", &
             "Number of distribution function energy values", error)
        call h5ltset_attribute_string_f(fid,"/npitch", "description", &
             "Number of distribution function pitch values", error)

        call h5ltset_attribute_string_f(fid,"/energy","description", &
             "Energy array", error)
        call h5ltset_attribute_string_f(fid,"/energy", "units","keV", error)
        call h5ltset_attribute_string_f(fid,"/pitch", "description", &
             "Pitch array: p = v_parallel/v  w.r.t. the magnetic field", error)
        call h5ltset_attribute_string_f(fid,"/earray","description", &
             "E3 energy array", error)
        call h5ltset_attribute_string_f(fid,"/earray", "units","keV", error)
    endif

    call h5ltset_attribute_string_f(fid, "/", "version", version, error)
    call h5ltset_attribute_string_f(fid,"/","description",&
         "CFPD signals 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)') 'Charged fusion products written to: ', trim(filename)
    endif

end subroutine write_cfpd_weights