calc_pgen_guga_pchb Function

    function calc_pgen_guga_pchb(this, ilutI, csf_i, ilutJ, excitInfo_in) result(pgen)
        class(GugaAliasSampler_t), intent(in) :: this
        integer(n_int), intent(in) :: ilutI(0:GugaBits%len_tot), ilutJ(GugaBits%len_tot)
        type(CSF_Info_t), intent(in) :: csf_i
        type(ExcitationInformation_t), intent(in), optional :: excitInfo_in
        type(ExcitationInformation_t) :: excitInfo
        real(dp) :: pgen
        integer :: ic
        integer :: nI(nel), nJ(nel)

        if (present(excitInfo_in)) then
            excitInfo = excitInfo_in
        else
            excitInfo = identify_excitation(ilutI, ilutJ)
        end if

        ic = get_excit_level_from_excitInfo(excitInfo)

        if (ic == 1) then
            if (t_guga_pchb_weighted_singles) then
                call decode_bit_det(nI, ilutI)
                call decode_bit_det(nJ, ilutJ)
                pgen = calc_pgen_mol_guga_single( &
                            ilutI, nI, csf_i, ilutJ, nJ, excitInfo)
            else
                pgen = calc_orb_pgen_uniform_singles(csf_i, excitInfo)
            end if

            pgen = pgen * pSingles

        else if (ic == 2) then

            pgen = pDoubles * this%calc_orb_pgen_guga_pchb_double(excitInfo)

        else
            pgen = 0.0_dp
        end if

    end function calc_pgen_guga_pchb