modwritemeta.f90

module modwritemeta

  use iso_fortran_env, only : OUTPUT_UNIT

  use modconstant, only : type_parameter_halofinder,&
       type_parameter_pfof,&
       type_parameter_pfof_snap,&
       type_parameter_pfof_cone,&
       type_parameter_conecreator,&
       type_parameter_conecreator_part,&
       type_parameter_conecreator_grav,&
       type_parameter_psod_snap,&
       type_info_ramses,&
       type_info_cone,&
       type_info_cone_part,&
       type_info_cone_grav,&
       type_common_metadata,&
       NAME_PFOF_SNAP,&
       NAME_PFOF_CONE,&
       NAME_PSOD_SNAP,&
       NAME_CONECREATOR_PART,&
       NAME_CONECREATOR_GRAV,&
       GIT_VERSION

  use modhdf5, only : hid_t,&
       H5_STR_LEN,&
       Hdf5_create_group,&
       Hdf5_open_group,&
       Hdf5_close_group,&
       Hdf5_read_attr,&
       Hdf5_read_data,&
       Hdf5_write_attr,&
       Hdf5_write_data

  use modmpicommons, only : procid

  implicit none

  private
  public :: Write_meta_common, &
       Write_meta_halofinder_parameter, &
       Write_meta_conecreator_parameter, &
       Write_meta_info_ramses, &
       Write_meta_info_cone

  ! Some metadata constants whose value may change for some software or options
  character(len=H5_STR_LEN) :: PARTICLE_TYPE = 'dark_matter'
  character(len=H5_STR_LEN) :: PFOF_CELL_ORDER = 'none'

contains
  !=======================================================================
  !> Write metadata common to every hdf5 ouput files of the pfof toolbox
  subroutine Write_meta_common(file_id, codename, npart, constant_mass, process_id)

    integer(kind=hid_t), intent(in) :: file_id
    character(len=H5_STR_LEN), intent(in) :: codename
    integer(kind=8), intent(in) :: npart
    integer(kind=4), intent(in), optional :: constant_mass
    integer(kind=4), intent(in), optional :: process_id

    character(len=H5_STR_LEN) :: groupname
    character(len=H5_STR_LEN) :: aname
    character(len=H5_STR_LEN) :: adata
    character(len=H5_STR_LEN) :: dsetname
    integer(kind=hid_t) :: gr_id
    integer(kind=4) :: tmpint4

#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Enter Write_meta_common on process ', procid
#endif

    groupname = 'metadata'
    call Hdf5_open_group(file_id,groupname, gr_id)
    ! Common metadata
    aname = 'created_by'
    call Hdf5_write_attr(gr_id, aname, codename)
    aname = 'version'
    call Hdf5_write_attr(gr_id, aname, GIT_VERSION)
    aname = 'simulation_code'
    adata = 'ramses'
    call Hdf5_write_attr(gr_id, aname,adata)
    aname = 'particle_type'
    adata = PARTICLE_TYPE
    call Hdf5_write_attr(gr_id, aname, adata)
    aname = 'process_id'
    if(present(process_id)) then
       tmpint4=process_id
    else
       tmpint4=0
    end if
    call Hdf5_write_attr(gr_id, aname, tmpint4)
    aname = 'constant_mass'
    if(present(constant_mass)) then
       tmpint4 = constant_mass
    else
       tmpint4 = -1
    end if
    call Hdf5_write_attr(gr_id, aname, tmpint4)
    aname = 'units'
    adata = 'ramses'
    call Hdf5_write_attr(gr_id, aname, adata)
    aname = 'pfof_cells_order'
    adata = PFOF_CELL_ORDER
    call Hdf5_write_attr(gr_id, aname, adata)
    ! Write the number of particles in this simulation as an integer(kind=8) dataset
    dsetname = 'npart_simulation'
    call Hdf5_write_data(gr_id, dsetname, npart)
    call Hdf5_close_group(gr_id)

#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Exit Write_meta_common on process ', procid
#endif

  end subroutine Write_meta_common

  !=======================================================================
  !> Write halofinder input parameters as metadata
  subroutine Write_meta_halofinder_parameter(file_id, param_halofinder)

    integer(kind=hid_t), intent(in) :: file_id
    class(type_parameter_halofinder), intent(in) :: param_halofinder

    character(len=H5_STR_LEN) :: groupname
    character(len=H5_STR_LEN) :: aname
    integer(kind=hid_t) :: gr_id
    integer(kind=hid_t) :: gr_pfof_id
    integer(kind=hid_t) :: gr_input_id
    integer(kind=hid_t) :: gr_fof_id
    integer(kind=hid_t) :: gr_output_id
    integer(kind=4) :: tmpint4
    character(len=16) :: c_b_name

#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Enter Write_meta_halofinder_parameter on process ', procid
#endif

    select type (param_halofinder)
    type is (type_parameter_pfof_snap)
       c_b_name = NAME_PFOF_SNAP
    type is (type_parameter_pfof_cone)
       c_b_name = NAME_PFOF_CONE
    type is (type_parameter_psod_snap)
       c_b_name = NAME_PSOD_SNAP
    end select
    groupname = 'metadata'
    call Hdf5_open_group(file_id,groupname, gr_id)
    ! pfof parameters:
    groupname = trim(c_b_name)//'_parameters'
    call Hdf5_create_group(gr_id, groupname, gr_pfof_id)
    groupname = 'input_parameters'
    call Hdf5_create_group(gr_pfof_id, groupname, gr_input_id)
    aname = 'input_path'
    call Hdf5_write_attr(gr_input_id, aname, param_halofinder%input_path)
    aname = 'part_input_file'
    call Hdf5_write_attr(gr_input_id, aname, param_halofinder%part_input_file)
    ! Write potential logical as an integer attribute (1=true, 0=false)
    aname = 'do_read_potential'
    if(param_halofinder%do_read_potential) then
       tmpint4=1
    else
       tmpint4=0
    end if
    call Hdf5_write_attr(gr_input_id, aname, tmpint4)
    ! Write force logical as an integer attribute (1=true, 0=false)
    aname = 'do_read_gravitational_field'
    if(param_halofinder%do_read_gravitational_field) then
       tmpint4=1
    else
       tmpint4=0
    end if
    call Hdf5_write_attr(gr_input_id, aname, tmpint4)
    select type (param_halofinder)
    type is (type_parameter_pfof_snap)
       aname = 'info_input_file'
       call Hdf5_write_attr(gr_input_id, aname, param_halofinder%info_input_file)
       aname = 'grpsize'
       call Hdf5_write_attr(gr_input_id, aname, param_halofinder%grpsize)
       aname = 'code_index'
       call Hdf5_write_attr(gr_input_id, aname, param_halofinder%code_index)
       aname = 'do_skip_star'
       if(param_halofinder%do_skip_star) then
          tmpint4=1
       else
          tmpint4=0
       end if
       call Hdf5_write_attr(gr_input_id, aname, tmpint4)
       aname = 'do_skip_metal'
       if(param_halofinder%do_skip_metal) then
          tmpint4=1
       else
          tmpint4=0
       end if
       call Hdf5_write_attr(gr_input_id, aname, tmpint4)
       aname='do_read_from_cube'
       tmpint4 = 0
       if(param_halofinder%do_read_from_cube) tmpint4 = 1
       call Hdf5_write_attr(gr_input_id, aname, tmpint4)
       aname = 'gatherread_factor'
       call Hdf5_write_attr(gr_input_id, aname, param_halofinder%gatherread_factor)
    type is (type_parameter_pfof_cone)
       aname = 'shell_first_id'
       call Hdf5_write_attr(gr_input_id, aname, param_halofinder%shell_first_id)
       aname = 'shell_last_id'
       call Hdf5_write_attr(gr_input_id, aname, param_halofinder%shell_last_id)
    end select
    call Hdf5_close_group(gr_input_id)
    groupname = 'fof_parameters'
    call Hdf5_create_group(gr_pfof_id, groupname, gr_fof_id)
    ! Write minimum halo mass Mmin as attribute
    aname = 'npart_halo_min'
    call Hdf5_write_attr(gr_fof_id, aname, param_halofinder%mmin)
    aname = 'npart_halo_max'
    call Hdf5_write_attr(gr_fof_id, aname, param_halofinder%mmax)
    ! Write doUnbinding as attribute (not implemented yet)
    aname = 'do_unbinding'
    if(param_halofinder%do_unbinding) then
       tmpint4=1
    else
       tmpint4=0
    end if
    call Hdf5_write_attr(gr_fof_id, aname, tmpint4)
    ! Write doSubHalo as attribute (not implemented yet)
    aname = 'do_subhalo'
    if(param_halofinder%do_subhalo) then
       tmpint4=1
    else
       tmpint4=0
    end if
    call Hdf5_write_attr(gr_fof_id, aname, tmpint4)
    select type (param_halofinder)
    type is (type_parameter_pfof_snap)
       ! Write percolation parameter b as attribute
       aname = 'percolation_length'
       call Hdf5_write_attr(gr_fof_id, aname, param_halofinder%percolation_length)
       aname = 'do_fof'
       if(param_halofinder%do_fof) then
          tmpint4=1
       else
          tmpint4=0
       end if
       call Hdf5_write_attr(gr_fof_id, aname, tmpint4)
    type is (type_parameter_pfof_cone)
       ! Write percolation parameter b as attribute
       aname = 'percolation_length'
       call Hdf5_write_attr(gr_fof_id, aname, param_halofinder%percolation_length)
    type is (type_parameter_psod_snap)
       aname = 'do_sod'
       if(param_halofinder%do_sod) then
          tmpint4=1
       else
          tmpint4=0
       end if
       call Hdf5_write_attr(gr_fof_id, aname, tmpint4)
    end select
    call Hdf5_close_group(gr_fof_id)
    groupname = 'output_parameters'
    call Hdf5_create_group(gr_pfof_id, groupname, gr_output_id)
    ! Write the simulation name as attribute
    aname = 'simulation_name'
    call Hdf5_write_attr(gr_output_id, aname, param_halofinder%simulation_name)
    ! Write do_timings as attribute
    aname = 'do_timings'
    if(param_halofinder%do_timings) then
       tmpint4=1
    else
       tmpint4=0
    end if
    call Hdf5_write_attr(gr_output_id, aname, tmpint4)
    select type (param_halofinder)
    type is (type_parameter_pfof_snap)
       ! write the snapshot number: if halo computed from a lightcone then outputNB = -1
       aname = 'snapshot'
       call Hdf5_write_attr(gr_output_id, aname, param_halofinder%snapshot)
       tmpint4 = 0
       if(param_halofinder%do_write_cube) tmpint4 = 1
       aname = 'do_write_cube'
       call Hdf5_write_attr(gr_output_id, aname, tmpint4)
       tmpint4 = 0
       if(param_halofinder%do_sort_cube) tmpint4 = 1
       aname = 'do_sort_cube'
       call Hdf5_write_attr(gr_output_id, aname, tmpint4)
       aname = 'gatherwrite_factor'
       call Hdf5_write_attr(gr_output_id, aname, param_halofinder%gatherwrite_factor)
    end select
    call Hdf5_close_group(gr_output_id)
    call Hdf5_close_group(gr_pfof_id)
    call Hdf5_close_group(gr_id)

#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Exit Write_meta_halofinder_parameter on process ', procid
#endif

  end subroutine Write_meta_halofinder_parameter


  !=======================================================================
  !> Write conecreator input parameters as metadata
  subroutine Write_meta_conecreator_parameter(file_id, param_cone)

    integer(kind=hid_t), intent(in) :: file_id
    class(type_parameter_conecreator), intent(in) :: param_cone

    character(len=H5_STR_LEN) :: groupname
    character(len=H5_STR_LEN) :: aname
    integer(kind=hid_t) :: gr_id
    integer(kind=hid_t) :: gr_param_id
    integer(kind=hid_t) :: gr_input_id
    integer(kind=hid_t) :: gr_output_id
    integer(kind=4) :: tmpint4
    character(len=16) :: c_b_name

#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Enter Write_meta_conecreator_parameter on process ', procid
#endif

    select type(param_cone)
    type is(type_parameter_conecreator_part)
       c_b_name = NAME_CONECREATOR_PART
    type is(type_parameter_conecreator_grav)
       c_b_name = NAME_CONECREATOR_GRAV
    end select
    groupname = 'metadata'
    call Hdf5_open_group(file_id, groupname, gr_id)
    groupname = c_b_name//'_parameters'
    call Hdf5_create_group(gr_id, groupname, gr_param_id)
    groupname = 'input_parameters'
    call Hdf5_create_group(gr_param_id, groupname, gr_input_id)
    aname = 'input_path'
    call Hdf5_write_attr(gr_input_id, aname, param_cone%input_path)
    aname = 'cone_input_file'
    call Hdf5_write_attr(gr_input_id, aname, param_cone%cone_input_file)
    aname = 'info_cone_input_file'
    call Hdf5_write_attr(gr_input_id, aname, param_cone%info_cone_input_file)
    aname = 'info_ramses_input_file'
    call Hdf5_write_attr(gr_input_id, aname, param_cone%info_ramses_input_file)
    aname = 'nfile'
    call Hdf5_write_attr(gr_input_id, aname, param_cone%nfile)
    aname = 'first_file'
    call Hdf5_write_attr(gr_input_id, aname, param_cone%first_file)
    aname = 'last_file'
    call Hdf5_write_attr(gr_input_id, aname, param_cone%last_file)
    aname = 'cone_max_radius'
    call Hdf5_write_attr(gr_input_id, aname, param_cone%cone_max_radius)
    select type(param_cone)
    type is (type_parameter_conecreator_part)
       aname = 'do_read_ramses_part_id'
       tmpint4 = 0
       if(param_cone%do_read_ramses_part_id) tmpint4 = 1
       call Hdf5_write_attr(gr_input_id, aname, tmpint4)
       aname = 'do_read_potential'
       tmpint4 = 0
       if(param_cone%do_read_potential) tmpint4 = 1
       call Hdf5_write_attr(gr_input_id, aname, tmpint4)
       aname = 'do_read_gravitational_field'
       tmpint4 = 0
       if(param_cone%do_read_gravitational_field) tmpint4 = 1
       call Hdf5_write_attr(gr_input_id, aname, tmpint4)
    type is (type_parameter_conecreator_grav)
       aname = 'nlevel'
       call Hdf5_write_attr(gr_input_id, aname, param_cone%nlevel)
       aname = 'levelmin'
       call Hdf5_write_attr(gr_input_id, aname, param_cone%levelmin)
    end select
    call Hdf5_close_group(gr_input_id)
    groupname = 'output_parameters'
    call Hdf5_create_group(gr_param_id, groupname, gr_output_id)
    aname = 'simulation_name'
    call Hdf5_write_attr(gr_output_id, aname, param_cone%simulation_name)
    aname = 'cube_size'
    call Hdf5_write_attr(gr_output_id, aname, param_cone%cube_size)
    call Hdf5_close_group(gr_output_id)
    call Hdf5_close_group(gr_param_id)
    call Hdf5_close_group(gr_id)

#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Exit Write_meta_conecreator_parameter on process ', procid
#endif

  end subroutine Write_meta_conecreator_parameter

  !=======================================================================
  !> Write ramses info as metadata
  subroutine Write_meta_info_ramses(file_id, inforamses, islast)

    integer(kind=hid_t), intent(in) :: file_id
    type(type_info_ramses), intent(in) :: inforamses
    logical(kind=4), intent(in) :: islast

    character(len=H5_STR_LEN) :: groupname
    character(len=H5_STR_LEN) :: aname
    integer(kind=hid_t) :: gr_id
    integer(kind=hid_t) :: gr_ramses_id

#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Enter Write_meta_info_ramses on process ', procid
#endif

    groupname = 'metadata'
    call Hdf5_open_group(file_id, groupname, gr_id)
    ! Ramses Info Metadata
    groupname = 'ramses_info'
    if(islast) then
       groupname = 'ramses_info_last'
    end if
    call Hdf5_create_group(gr_id,groupname,gr_ramses_id)
    aname = 'ncpu'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%ncpu)
    aname = 'ndim'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%ndim)
    aname = 'levelmin'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%levelmin)
    aname = 'levelmax'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%levelmax)
    aname = 'ngridmax'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%ngridmax)
    aname = 'nstep_coarse'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%nstep_coarse)
    aname = 'boxlen'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%boxlen)
    aname = 'time'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%time)
    aname = 'aexp'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%aexp)
    aname = 'h0'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%h0)
    aname = 'omega_m'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%omega_m)
    aname = 'omega_l'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%omega_l)
    aname = 'omega_k'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%omega_k)
    aname = 'omega_b'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%omega_b)
    aname = 'unit_l'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%unit_l)
    aname = 'unit_d'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%unit_d)
    aname = 'unit_t'
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%unit_t)
    call Hdf5_close_group(gr_ramses_id)
    call Hdf5_close_group(gr_id)

#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Exit Write_meta_info_ramses on process ', procid
#endif

  end subroutine Write_meta_info_ramses


  !=======================================================================
  !> Write light cone info as metadata
  subroutine Write_meta_info_cone(file_id, infocone, islast)

    integer(kind=hid_t), intent(in) :: file_id
    class(type_info_cone), intent(in) :: infocone
    logical(kind=4), intent(in) :: islast

    character(len=H5_STR_LEN) :: groupname
    character(len=H5_STR_LEN) :: aname

    integer(kind=hid_t) :: gr_id
    integer(kind=hid_t) :: gr_cone_id

#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Enter Write_meta_info_cone on process ', procid
#endif

    groupname = 'metadata'
    call Hdf5_open_group(file_id, groupname, gr_id)

    groupname = 'cone_info'
    if(islast) then
       groupname = 'cone_info_last'
    end if
    call Hdf5_create_group(gr_id, groupname, gr_cone_id)
    aname = 'ncpu'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%ncpu)
    aname = 'nstride'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%nstride)
    aname = 'nstep_coarse'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%nstep_coarse)
    aname = 'cone_id'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%cone_id)
    aname = 'nglobalfile'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%nglobalfile)
    aname = 'isfullsky'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%isfullsky)
    aname = 'aexp'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%aexp)
    aname = 'observer_x'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%observer_x)
    aname = 'observer_y'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%observer_y)
    aname = 'observer_z'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%observer_z)
    aname = 'observer_rds'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%observer_rds)
    aname = 'cone_zlim'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%cone_zlim)
    aname = 'amax'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%amax)
    aname = 'amin'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%amin)
    aname = 'zmax'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%zmax)
    aname = 'zmin'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%zmin)
    aname = 'dmax'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%dmax)
    aname = 'dmin'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%dmin)
    aname = 'dtol'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%dtol)
    aname = 'thetay'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%thetay)
    aname = 'thetaz'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%thetaz)
    aname = 'theta'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%theta)
    aname = 'phi'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%phi)
    aname = 'aendconem2'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%aendconem2)
    aname = 'aendconem1'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%aendconem1)
    aname = 'aendcone'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%aendcone)
    aname = 'zendconem2'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%zendconem2)
    aname = 'zendconem1'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%zendconem1)
    aname = 'zendcone'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%zendcone)
    aname = 'dendconem2'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%dendconem2)
    aname = 'dendconem1'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%dendconem1)
    aname = 'dendcone'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%dendcone)
    aname = 'future'
    call Hdf5_write_attr(gr_cone_id, aname, infocone%future)
    select type(infocone)
    type is(type_info_cone_part)
       aname = 'npart'
       call Hdf5_write_data(gr_cone_id, aname, infocone%npart)
       aname = 'aexpold'
       call Hdf5_write_attr(gr_cone_id, aname, infocone%aexpold)
       aname = 'zexpold'
       call Hdf5_write_attr(gr_cone_id, aname, infocone%zexpold)
       aname = 'zexp'
       call Hdf5_write_attr(gr_cone_id, aname, infocone%zexp)
       aname = 'dexpold'
       call Hdf5_write_attr(gr_cone_id, aname, infocone%dexpold)
       aname = 'dexp'
       call Hdf5_write_attr(gr_cone_id, aname, infocone%dexp)
    type is(type_info_cone_grav)
       aname = 'nglobalcell'
       call Hdf5_write_data(gr_cone_id, aname, infocone%nglobalcell)
       aname = 'nlevel'
       call Hdf5_write_attr(gr_cone_id, aname, infocone%nlevel)
       aname = 'levelmin'
       call Hdf5_write_attr(gr_cone_id, aname, infocone%levelmin)
       aname = 'levelmax'
       call Hdf5_write_attr(gr_cone_id, aname, infocone%levelmax)
    end select
    call Hdf5_close_group(gr_cone_id)
    call Hdf5_close_group(gr_id)

#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Exit Write_meta_info_cone on process ', procid
#endif

  end subroutine Write_meta_info_cone

end module modwritemeta