write_npa Subroutine

public subroutine write_npa()

Writes npa to a HDF5 file

Arguments

None

Calls

proc~~write_npa~~CallsGraph proc~write_npa write_npa proc~my_rank my_rank proc~write_npa->proc~my_rank interface~h5ltmake_compressed_dataset_double_f h5ltmake_compressed_dataset_double_f proc~write_npa->interface~h5ltmake_compressed_dataset_double_f proc~num_ranks num_ranks proc~write_npa->proc~num_ranks interface~h5ltmake_compressed_dataset_int_f h5ltmake_compressed_dataset_int_f proc~write_npa->interface~h5ltmake_compressed_dataset_int_f h5gcreate_f h5gcreate_f proc~write_npa->h5gcreate_f interface~parallel_sum parallel_sum proc~write_npa->interface~parallel_sum h5fcreate_f h5fcreate_f proc~write_npa->h5fcreate_f h5open_f h5open_f proc~write_npa->h5open_f h5ltmake_dataset_int_f h5ltmake_dataset_int_f proc~write_npa->h5ltmake_dataset_int_f h5fclose_f h5fclose_f proc~write_npa->h5fclose_f h5gclose_f h5gclose_f proc~write_npa->h5gclose_f h5close_f h5close_f proc~write_npa->h5close_f h5ltset_attribute_string_f h5ltset_attribute_string_f proc~write_npa->h5ltset_attribute_string_f 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_7 h5ltmake_compressed_dataset_double_f_7 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_7 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_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_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_6 h5ltmake_compressed_dataset_double_f_6 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_6 proc~h5ltmake_compressed_dataset_int_f_1 h5ltmake_compressed_dataset_int_f_1 interface~h5ltmake_compressed_dataset_int_f->proc~h5ltmake_compressed_dataset_int_f_1 proc~h5ltmake_compressed_dataset_int_f_2 h5ltmake_compressed_dataset_int_f_2 interface~h5ltmake_compressed_dataset_int_f->proc~h5ltmake_compressed_dataset_int_f_2 proc~h5ltmake_compressed_dataset_int_f_3 h5ltmake_compressed_dataset_int_f_3 interface~h5ltmake_compressed_dataset_int_f->proc~h5ltmake_compressed_dataset_int_f_3 proc~h5ltmake_compressed_dataset_int_f_4 h5ltmake_compressed_dataset_int_f_4 interface~h5ltmake_compressed_dataset_int_f->proc~h5ltmake_compressed_dataset_int_f_4 proc~h5ltmake_compressed_dataset_int_f_6 h5ltmake_compressed_dataset_int_f_6 interface~h5ltmake_compressed_dataset_int_f->proc~h5ltmake_compressed_dataset_int_f_6 proc~h5ltmake_compressed_dataset_int_f_7 h5ltmake_compressed_dataset_int_f_7 interface~h5ltmake_compressed_dataset_int_f->proc~h5ltmake_compressed_dataset_int_f_7 proc~h5ltmake_compressed_dataset_int_f_5 h5ltmake_compressed_dataset_int_f_5 interface~h5ltmake_compressed_dataset_int_f->proc~h5ltmake_compressed_dataset_int_f_5 proc~parallel_sum_i3 parallel_sum_i3 interface~parallel_sum->proc~parallel_sum_i3 proc~parallel_sum_d0 parallel_sum_d0 interface~parallel_sum->proc~parallel_sum_d0 proc~parallel_sum_d1 parallel_sum_d1 interface~parallel_sum->proc~parallel_sum_d1 proc~parallel_sum_d5 parallel_sum_d5 interface~parallel_sum->proc~parallel_sum_d5 proc~parallel_sum_i1 parallel_sum_i1 interface~parallel_sum->proc~parallel_sum_i1 proc~parallel_sum_i0 parallel_sum_i0 interface~parallel_sum->proc~parallel_sum_i0 proc~parallel_sum_i2 parallel_sum_i2 interface~parallel_sum->proc~parallel_sum_i2 proc~parallel_sum_d4 parallel_sum_d4 interface~parallel_sum->proc~parallel_sum_d4 proc~parallel_sum_d3 parallel_sum_d3 interface~parallel_sum->proc~parallel_sum_d3 proc~parallel_sum_d2 parallel_sum_d2 interface~parallel_sum->proc~parallel_sum_d2 proc~h5ltmake_compressed_dataset_int_f_1->h5ltmake_dataset_int_f proc~chunk_size chunk_size proc~h5ltmake_compressed_dataset_int_f_1->proc~chunk_size h5sclose_f h5sclose_f proc~h5ltmake_compressed_dataset_int_f_1->h5sclose_f h5pset_chunk_f h5pset_chunk_f proc~h5ltmake_compressed_dataset_int_f_1->h5pset_chunk_f h5pcreate_f h5pcreate_f proc~h5ltmake_compressed_dataset_int_f_1->h5pcreate_f h5pset_shuffle_f h5pset_shuffle_f proc~h5ltmake_compressed_dataset_int_f_1->h5pset_shuffle_f h5pset_deflate_f h5pset_deflate_f proc~h5ltmake_compressed_dataset_int_f_1->h5pset_deflate_f h5dcreate_f h5dcreate_f proc~h5ltmake_compressed_dataset_int_f_1->h5dcreate_f h5dwrite_f h5dwrite_f proc~h5ltmake_compressed_dataset_int_f_1->h5dwrite_f h5dclose_f h5dclose_f proc~h5ltmake_compressed_dataset_int_f_1->h5dclose_f h5pclose_f h5pclose_f proc~h5ltmake_compressed_dataset_int_f_1->h5pclose_f h5screate_simple_f h5screate_simple_f proc~h5ltmake_compressed_dataset_int_f_1->h5screate_simple_f proc~h5ltmake_compressed_dataset_int_f_2->h5ltmake_dataset_int_f proc~h5ltmake_compressed_dataset_int_f_2->proc~chunk_size proc~h5ltmake_compressed_dataset_int_f_2->h5sclose_f proc~h5ltmake_compressed_dataset_int_f_2->h5pset_chunk_f proc~h5ltmake_compressed_dataset_int_f_2->h5pcreate_f proc~h5ltmake_compressed_dataset_int_f_2->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_int_f_2->h5pset_deflate_f proc~h5ltmake_compressed_dataset_int_f_2->h5dcreate_f proc~h5ltmake_compressed_dataset_int_f_2->h5dwrite_f proc~h5ltmake_compressed_dataset_int_f_2->h5dclose_f proc~h5ltmake_compressed_dataset_int_f_2->h5pclose_f proc~h5ltmake_compressed_dataset_int_f_2->h5screate_simple_f proc~parallel_sum_i3->proc~parallel_sum_i1 mpi_allreduce mpi_allreduce proc~parallel_sum_i3->mpi_allreduce proc~h5ltmake_compressed_dataset_int_f_3->h5ltmake_dataset_int_f proc~h5ltmake_compressed_dataset_int_f_3->proc~chunk_size proc~h5ltmake_compressed_dataset_int_f_3->h5sclose_f proc~h5ltmake_compressed_dataset_int_f_3->h5pset_chunk_f proc~h5ltmake_compressed_dataset_int_f_3->h5pcreate_f proc~h5ltmake_compressed_dataset_int_f_3->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_int_f_3->h5pset_deflate_f proc~h5ltmake_compressed_dataset_int_f_3->h5dcreate_f proc~h5ltmake_compressed_dataset_int_f_3->h5dwrite_f proc~h5ltmake_compressed_dataset_int_f_3->h5dclose_f proc~h5ltmake_compressed_dataset_int_f_3->h5pclose_f proc~h5ltmake_compressed_dataset_int_f_3->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->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_1->h5pcreate_f h5ltmake_dataset_double_f h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_1->h5ltmake_dataset_double_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_int_f_4->h5ltmake_dataset_int_f proc~h5ltmake_compressed_dataset_int_f_4->proc~chunk_size proc~h5ltmake_compressed_dataset_int_f_4->h5sclose_f proc~h5ltmake_compressed_dataset_int_f_4->h5pset_chunk_f proc~h5ltmake_compressed_dataset_int_f_4->h5pcreate_f proc~h5ltmake_compressed_dataset_int_f_4->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_int_f_4->h5pset_deflate_f proc~h5ltmake_compressed_dataset_int_f_4->h5dcreate_f proc~h5ltmake_compressed_dataset_int_f_4->h5dwrite_f proc~h5ltmake_compressed_dataset_int_f_4->h5dclose_f proc~h5ltmake_compressed_dataset_int_f_4->h5pclose_f proc~h5ltmake_compressed_dataset_int_f_4->h5screate_simple_f proc~parallel_sum_d0->mpi_allreduce proc~parallel_sum_d1->proc~parallel_sum_d1 proc~parallel_sum_d1->mpi_allreduce proc~parallel_sum_d5->proc~parallel_sum_d1 proc~parallel_sum_d5->mpi_allreduce 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->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_5->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_5->h5ltmake_dataset_double_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_7->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_7->h5sclose_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->h5ltmake_dataset_double_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 proc~parallel_sum_i1->proc~parallel_sum_i1 proc~parallel_sum_i1->mpi_allreduce proc~parallel_sum_i0->mpi_allreduce proc~h5ltmake_compressed_dataset_double_f_4->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_4->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_4->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_4->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_4->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_4->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_4->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_4->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_4->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5screate_simple_f proc~parallel_sum_i2->proc~parallel_sum_i1 proc~parallel_sum_i2->mpi_allreduce 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->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_3->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_3->h5ltmake_dataset_double_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_int_f_6->h5ltmake_dataset_int_f proc~h5ltmake_compressed_dataset_int_f_6->proc~chunk_size proc~h5ltmake_compressed_dataset_int_f_6->h5sclose_f proc~h5ltmake_compressed_dataset_int_f_6->h5pset_chunk_f proc~h5ltmake_compressed_dataset_int_f_6->h5pcreate_f proc~h5ltmake_compressed_dataset_int_f_6->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_int_f_6->h5pset_deflate_f proc~h5ltmake_compressed_dataset_int_f_6->h5dcreate_f proc~h5ltmake_compressed_dataset_int_f_6->h5dwrite_f proc~h5ltmake_compressed_dataset_int_f_6->h5dclose_f proc~h5ltmake_compressed_dataset_int_f_6->h5pclose_f proc~h5ltmake_compressed_dataset_int_f_6->h5screate_simple_f proc~h5ltmake_compressed_dataset_int_f_7->h5ltmake_dataset_int_f proc~h5ltmake_compressed_dataset_int_f_7->proc~chunk_size proc~h5ltmake_compressed_dataset_int_f_7->h5sclose_f proc~h5ltmake_compressed_dataset_int_f_7->h5pset_chunk_f proc~h5ltmake_compressed_dataset_int_f_7->h5pcreate_f proc~h5ltmake_compressed_dataset_int_f_7->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_int_f_7->h5pset_deflate_f proc~h5ltmake_compressed_dataset_int_f_7->h5dcreate_f proc~h5ltmake_compressed_dataset_int_f_7->h5dwrite_f proc~h5ltmake_compressed_dataset_int_f_7->h5dclose_f proc~h5ltmake_compressed_dataset_int_f_7->h5pclose_f proc~h5ltmake_compressed_dataset_int_f_7->h5screate_simple_f proc~parallel_sum_d4->proc~parallel_sum_d1 proc~parallel_sum_d4->mpi_allreduce proc~parallel_sum_d3->proc~parallel_sum_d1 proc~parallel_sum_d3->mpi_allreduce 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->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_2->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_2->h5ltmake_dataset_double_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_int_f_5->h5ltmake_dataset_int_f proc~h5ltmake_compressed_dataset_int_f_5->proc~chunk_size proc~h5ltmake_compressed_dataset_int_f_5->h5sclose_f proc~h5ltmake_compressed_dataset_int_f_5->h5pset_chunk_f proc~h5ltmake_compressed_dataset_int_f_5->h5pcreate_f proc~h5ltmake_compressed_dataset_int_f_5->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_int_f_5->h5pset_deflate_f proc~h5ltmake_compressed_dataset_int_f_5->h5dcreate_f proc~h5ltmake_compressed_dataset_int_f_5->h5dwrite_f proc~h5ltmake_compressed_dataset_int_f_5->h5dclose_f proc~h5ltmake_compressed_dataset_int_f_5->h5pclose_f proc~h5ltmake_compressed_dataset_int_f_5->h5screate_simple_f proc~parallel_sum_d2->proc~parallel_sum_d1 proc~parallel_sum_d2->mpi_allreduce 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->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_6->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_6->h5ltmake_dataset_double_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

Called by

proc~~write_npa~~CalledByGraph proc~write_npa write_npa program~fidasim fidasim program~fidasim->proc~write_npa

Contents

Source Code


Source Code

subroutine write_npa
    !+ Writes [[libfida:npa]] to a HDF5 file
    integer(HID_T) :: fid, gid
    integer(HSIZE_T), dimension(3) :: dim3
    integer(HSIZE_T), dimension(2) :: dim2
    integer(HSIZE_T), dimension(1) :: d
    integer :: error

    integer, dimension(:), allocatable :: dcount
    real(Float64), dimension(:,:), allocatable :: ri, rf
    real(Float64), dimension(:), allocatable :: weight, energy, pitch
    integer, dimension(:), allocatable :: det, orbit_type
    integer :: i, npart, c, start_index, end_index
    character(charlim) :: filename = ''
    logical :: do_write = .True.

#ifdef _MPI
    integer :: rank
    integer, dimension(:), allocatable :: npart_image

    rank = my_rank()
    allocate(npart_image(0:num_ranks()-1))
#endif

#ifdef _MPI
    if(my_rank().ne.0) do_write = .False.
#endif

    filename=trim(adjustl(inputs%result_dir))//"/"//trim(adjustl(inputs%runid))//"_npa.h5"

    if(do_write) then
        !Open HDF5 interface
        call h5open_f(error)

        !Create file overwriting any existing file
        call h5fcreate_f(filename, H5F_ACC_TRUNC_F, fid, error)
    endif

    !! Active
    allocate(dcount(npa_chords%nchan))
    npart = npa%npart
    if(npart.gt.0) then
        do i=1,npa_chords%nchan
            dcount(i) = count(npa%part%detector.eq.i)
        enddo
    else
        dcount = 0
    endif

#ifdef _MPI
    call parallel_sum(dcount)
    call parallel_sum(npart)
#endif

    c = 0
#ifdef _MPI
    npart_image(:) = 0
    npart_image(rank) = npa%npart
    call parallel_sum(npart_image)

    do i=0,rank-1
        c = c +  npart_image(i)
    enddo
#endif
    start_index = 1 + c
    end_index = npa%npart + c

    if(npart.gt.0) then
        allocate(ri(3,npart),rf(3,npart))
        allocate(weight(npart),energy(npart),pitch(npart))
        allocate(det(npart),orbit_type(npart))
        ri = 0.d0     ; rf = 0.d0
        weight = 0.d0 ; energy = 0.d0
        pitch = 0.d0  ; det = 0
        orbit_type = 0
        c = 1
        do i=start_index, end_index
            ri(1,i) = npa%part(c)%xi
            ri(2,i) = npa%part(c)%yi
            ri(3,i) = npa%part(c)%zi
            rf(1,i) = npa%part(c)%xf
            rf(2,i) = npa%part(c)%yf
            rf(3,i) = npa%part(c)%zf
            weight(i) = npa%part(c)%weight
            energy(i) = npa%part(c)%energy
            pitch(i) = npa%part(c)%pitch
            det(i) = npa%part(c)%detector
            orbit_type(i) = npa%part(c)%class
            c = c + 1
        enddo
#ifdef _MPI
        call parallel_sum(ri)
        call parallel_sum(rf)
        call parallel_sum(weight)
        call parallel_sum(energy)
        call parallel_sum(pitch)
        call parallel_sum(det)
        call parallel_sum(orbit_type)
#endif
    endif

    if(do_write.and.(inputs%calc_npa.ge.1)) then
        !Write Active Flux
        d(1) = 1
        dim2 = [npa%nenergy, npa%nchan]
        dim3 = [npa%nenergy, npa%nchan, particles%nclass]
        if(particles%nclass.gt.1) then
            call h5ltmake_dataset_int_f(fid,"/nclass", 0, d, [particles%nclass], error)
            call h5ltmake_compressed_dataset_double_f(fid,"/flux",3,dim3,npa%flux, error)
            call h5ltset_attribute_string_f(fid,"/flux", "description", &
                 "Active Neutral flux: flux(energy,chan,class)", error)
        else
            call h5ltmake_compressed_dataset_double_f(fid,"/flux",2,dim3(1:2),npa%flux(:,:,1), error)
            call h5ltset_attribute_string_f(fid,"/flux", "description", &
                 "Active Neutral flux: flux(energy,chan)", error)
        endif
        call h5ltset_attribute_string_f(fid,"/flux", "units","neutrals/(s*dE)", error)

        call h5ltmake_dataset_int_f(fid,"/nenergy", 0, d, [npa%nenergy], error)
        call h5ltmake_dataset_int_f(fid,"/nchan", 0, d, [npa%nchan], error)
        call h5ltmake_compressed_dataset_double_f(fid,"/energy",1,dim2(1:1),&
             npa%energy, error)
        call h5ltmake_compressed_dataset_double_f(fid,"/radius",1,dim2(2:2),&
             npa_chords%radius, error)
        call h5ltmake_compressed_dataset_int_f(fid,"/count",1,dim2(2:2), dcount, error)

        !Add attributes
        call h5ltset_attribute_string_f(fid, "/", "version", version, error)
        call h5ltset_attribute_string_f(fid,"/","description", &
             "NPA flux calculated by FIDASIM",error)
        call h5ltset_attribute_string_f(fid,"/nenergy","description",&
             "Number of energy values",error)
        call h5ltset_attribute_string_f(fid,"/nchan","description",&
             "Number of channels",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,"/radius","description", &
             "Detector line of sight radius at midplane or tangency point", error)
        call h5ltset_attribute_string_f(fid,"/radius","units","cm",error)
        call h5ltset_attribute_string_f(fid,"/count","description", &
             "Number of particles that hit the detector: count(chan)", error)


        if((npart.gt.0).and.(inputs%calc_npa.ge.2)) then
            !Create Group
            call h5gcreate_f(fid,"/particles",gid, error)
            call h5ltmake_dataset_int_f(gid, "nparticle", 0, d, [npart], error)
            d(1) = npart
            dim2 = [3, npart]
            call h5ltmake_compressed_dataset_double_f(gid,"ri",2,dim2, ri, error)
            call h5ltmake_compressed_dataset_double_f(gid,"rf",2,dim2, rf, error)
            call h5ltmake_compressed_dataset_double_f(gid,"pitch",1,d, pitch, error)
            call h5ltmake_compressed_dataset_double_f(gid,"energy",1,d, energy, error)
            call h5ltmake_compressed_dataset_double_f(gid,"weight",1,d, weight, error)
            call h5ltmake_compressed_dataset_int_f(gid,"detector",1,d, det, error)
            call h5ltmake_compressed_dataset_int_f(gid,"class",1,d, orbit_type, error)

            !Add attributes
            call h5ltset_attribute_string_f(gid,"nparticle","description", &
                 "Number of particles that hit a detector", error)
            call h5ltset_attribute_string_f(gid,"ri","description", &
                 "Neutral particle's birth position in machine coordinates: ri([x,y,z],particle)", error)
            call h5ltset_attribute_string_f(gid,"ri","units", "cm", error)
            call h5ltset_attribute_string_f(gid,"rf","description", &
                 "Neutral particle's hit position in machine coordinates: rf([x,y,z],particle)", error)
            call h5ltset_attribute_string_f(gid,"rf","units", "cm", error)
            call h5ltset_attribute_string_f(gid,"pitch","description", &
                 "Pitch value of the neutral particle: p = v_parallel/v  w.r.t. the magnetic field", error)
            call h5ltset_attribute_string_f(gid,"energy","description", &
                 "Energy value of the neutral particle", error)
            call h5ltset_attribute_string_f(gid,"energy","units","keV",error)
            call h5ltset_attribute_string_f(gid,"weight","description", &
                 "Neutral particle's contribution to the flux", error)
            call h5ltset_attribute_string_f(gid,"weight","units","neutrals/s",error)
            call h5ltset_attribute_string_f(gid,"detector","description", &
                 "Detector that the neutral particle hit", error)
            call h5ltset_attribute_string_f(gid,"class","description", &
                 "Class of the neutral particle", error)

            call h5ltset_attribute_string_f(fid,"/particles","coordinate_system", &
                 "Right-handed cartesian",error)
            call h5ltset_attribute_string_f(fid,"/particles","description", &
                 "Active NPA Monte Carlo particles",error)

            !Close group
            call h5gclose_f(gid, error)
        endif
    endif

    deallocate(dcount)
    if(npart.gt.0) then
        deallocate(ri,rf)
        deallocate(energy,pitch,weight,det,orbit_type)
    endif

    !! Passive
    allocate(dcount(npa_chords%nchan))
    npart = pnpa%npart
    if(npart.gt.0) then
        do i=1,npa_chords%nchan
            dcount(i) = count(pnpa%part%detector.eq.i)
        enddo
    else
        dcount = 0
    endif

#ifdef _MPI
    call parallel_sum(dcount)
    call parallel_sum(npart)
#endif

    c = 0
#ifdef _MPI
    npart_image(:) = 0
    npart_image(rank) = pnpa%npart
    call parallel_sum(npart_image)

    do i=0,rank-1
        c = c +  npart_image(i)
    enddo
#endif
    start_index = 1 + c
    end_index = pnpa%npart + c

    if(npart.gt.0) then
        allocate(ri(3,npart),rf(3,npart))
        allocate(weight(npart),energy(npart),pitch(npart))
        allocate(det(npart),orbit_type(npart))
        ri = 0.d0     ; rf = 0.d0
        weight = 0.d0 ; energy = 0.d0
        pitch = 0.d0  ; det = 0
        orbit_type = 0
        c = 1
        do i=start_index, end_index
            ri(1,i) = pnpa%part(c)%xi
            ri(2,i) = pnpa%part(c)%yi
            ri(3,i) = pnpa%part(c)%zi
            rf(1,i) = pnpa%part(c)%xf
            rf(2,i) = pnpa%part(c)%yf
            rf(3,i) = pnpa%part(c)%zf
            weight(i) = pnpa%part(c)%weight
            energy(i) = pnpa%part(c)%energy
            pitch(i) = pnpa%part(c)%pitch
            det(i) = pnpa%part(c)%detector
            orbit_type(i) = pnpa%part(c)%class
            c = c + 1
        enddo
#ifdef _MPI
        call parallel_sum(ri)
        call parallel_sum(rf)
        call parallel_sum(weight)
        call parallel_sum(energy)
        call parallel_sum(pitch)
        call parallel_sum(det)
        call parallel_sum(orbit_type)
#endif
    endif

    if(do_write.and.(inputs%calc_pnpa.ge.1)) then
        !Write Passive Flux
        d(1) = 1
        dim2 = [pnpa%nenergy, pnpa%nchan]
        dim3 = [pnpa%nenergy, pnpa%nchan, particles%nclass]
        if(particles%nclass.gt.1) then
            call h5ltmake_compressed_dataset_double_f(fid,"/pflux",3,dim3,pnpa%flux, error)
            call h5ltset_attribute_string_f(fid,"/pflux", "description", &
                 "Passive Neutral flux: pflux(energy,chan,class)", error)
        else
            call h5ltmake_compressed_dataset_double_f(fid,"/pflux",2,dim3(1:2),pnpa%flux(:,:,1), error)
            call h5ltset_attribute_string_f(fid,"/pflux", "description", &
                 "Passive Neutral flux: pflux(energy,chan)", error)
        endif
        call h5ltset_attribute_string_f(fid,"/pflux", "units","neutrals/(s*dE)", error)
        call h5ltmake_compressed_dataset_int_f(fid,"/pcount",1,dim2(2:2), dcount, error)
        call h5ltset_attribute_string_f(fid,"/pcount","description", &
             "Number of passive particles that hit the detector: pcount(chan)", error)

        if(inputs%calc_npa.le.0) then
            call h5ltmake_dataset_int_f(fid,"/nenergy", 0, d, [pnpa%nenergy], error)
            call h5ltmake_dataset_int_f(fid,"/nchan", 0, d, [pnpa%nchan], error)
            call h5ltmake_compressed_dataset_double_f(fid,"/energy",1,dim2(1:1),&
                 pnpa%energy, error)
            call h5ltmake_compressed_dataset_double_f(fid,"/radius",1,dim2(2:2),&
                 npa_chords%radius, error)

            !Add attributes
            call h5ltset_attribute_string_f(fid, "/", "version", version, error)
            call h5ltset_attribute_string_f(fid,"/","description", &
                 "NPA flux calculated by FIDASIM",error)
            call h5ltset_attribute_string_f(fid,"/nenergy","description",&
                 "Number of energy values",error)
            call h5ltset_attribute_string_f(fid,"/nchan","description",&
                 "Number of channels",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,"/radius","description", &
                 "Detector line of sight radius at midplane or tangency point", error)
            call h5ltset_attribute_string_f(fid,"/radius","units","cm",error)
        endif

        if((npart.gt.0).and.(inputs%calc_pnpa.ge.2)) then
            !Create Group
            call h5gcreate_f(fid,"/passive_particles",gid, error)
            call h5ltmake_dataset_int_f(gid, "nparticle", 0, d, [npart], error)
            d(1) = npart
            dim2 = [3, npart]
            call h5ltmake_compressed_dataset_double_f(gid,"ri",2,dim2, ri, error)
            call h5ltmake_compressed_dataset_double_f(gid,"rf",2,dim2, rf, error)
            call h5ltmake_compressed_dataset_double_f(gid,"pitch",1,d, pitch, error)
            call h5ltmake_compressed_dataset_double_f(gid,"energy",1,d, energy, error)
            call h5ltmake_compressed_dataset_double_f(gid,"weight",1,d, weight, error)
            call h5ltmake_compressed_dataset_int_f(gid,"detector",1,d, det, error)
            call h5ltmake_compressed_dataset_int_f(gid,"class",1,d, orbit_type, error)

            !Add attributes
            call h5ltset_attribute_string_f(gid,"nparticle","description", &
                 "Number of particles that hit a detector", error)
            call h5ltset_attribute_string_f(gid,"ri","description", &
                 "Neutral particle's birth position in machine coordinates: ri([x,y,z],particle)", error)
            call h5ltset_attribute_string_f(gid,"ri","units", "cm", error)
            call h5ltset_attribute_string_f(gid,"rf","description", &
                 "Neutral particle's hit position in machine coordinates: rf([x,y,z],particle)", error)
            call h5ltset_attribute_string_f(gid,"rf","units", "cm", error)
            call h5ltset_attribute_string_f(gid,"pitch","description", &
                 "Pitch value of the neutral particle: p = v_parallel/v  w.r.t. the magnetic field", error)
            call h5ltset_attribute_string_f(gid,"energy","description", &
                 "Energy value of the neutral particle", error)
            call h5ltset_attribute_string_f(gid,"energy","units","keV",error)
            call h5ltset_attribute_string_f(gid,"weight","description", &
                 "Neutral particle's contribution to the flux", error)
            call h5ltset_attribute_string_f(gid,"weight","units","neutrals/s",error)
            call h5ltset_attribute_string_f(gid,"detector","description", &
                 "Detector that the neutral particle hit", error)
            call h5ltset_attribute_string_f(gid,"class","description", &
                 "Class of the neutral particle", error)

            call h5ltset_attribute_string_f(fid,"/passive_particles","coordinate_system", &
                 "Right-handed cartesian",error)
            call h5ltset_attribute_string_f(fid,"/passive_particles","description", &
                 "Passive NPA Monte Carlo particles",error)

            !Close group
            call h5gclose_f(gid, error)
        endif
    endif

    if(do_write) then
        !Close file
        call h5fclose_f(fid, error)

        !Close HDF5 interface
        call h5close_f(error)
    endif

    if(inputs%verbose.ge.1) then
        write(*,'(T4,a,a)') 'NPA data written to: ',trim(filename)
    endif

#ifdef _MPI
    deallocate(npart_image)
#endif

end subroutine write_npa