return_perturbed_ground_spec Subroutine

    subroutine return_perturbed_ground_spec(perturbs, ilut_list, pert_ground_local, all_norm_pert)

        ! This routine takes a perturbation operator, the list of all
        ! possible determinants (ilut_list). It reads in a popsfile, applies
        ! the perturbation to the read in state, and then copies the amplitudes
        ! of this perturbed popsfile state to pert_ground_local, using an
        ! ordering and addressing defined by ilut_list. It all divides these
        ! amplitudes by the norm of the read in wave function *BEFORE* the
        ! perturbation is applied (as needed for spectral calculations) and
        ! returns the norm of the perturbed popsfile as an output.

        use bit_rep_data, only: nifd, extract_sign
        use bit_reps, only: decode_bit_det
        use DetBitOps, only: DetBitEq, EncodeBitDet, ilut_lt, ilut_gt
        use FciMCData, only: TotWalkers, CurrentDets
        use PopsfileMod, only: read_popsfile_wrapper
        use sort_mod, only: sort
        use SystemData, only: nel

        type(perturbation), intent(in), allocatable :: perturbs(:)
        integer(n_int), intent(in) :: ilut_list(:, :)
        real(dp), intent(inout) :: pert_ground_local(:)
        real(dp), intent(out) :: all_norm_pert

        integer :: i, j
        integer :: nI(nel)
        real(dp) :: real_sign(lenof_sign)
        real(dp) :: norm_pert
        character(len=*), parameter :: t_r = 'return_perturbed_ground_spec'

        pert_ground_local = 0.0_dp

        ! Read in the POPSFILE and apply the perturbation operator in perturbs,
        ! as we do so. Afterwards, the perturbed wave function will be stored in
        ! CurrentDets.
        call read_popsfile_wrapper(perturbs)

        j = 0
        norm_pert = 0.0_dp
        call sort(CurrentDets(:, 1:TotWalkers), ilut_lt, ilut_gt)
        do i = 1, int(TotWalkers)
            call extract_sign(CurrentDets(:, i), real_sign)
            ! If the POPSFILE was generated from a calculation with the linear-scaling
            ! algorithm then it is possible to have the same determinant twice in a
            ! row, once with zero weight. If this is the case, skip the second case.
            if (i > 1) then
                if (DetBitEq(CurrentDets(:, i - 1), CurrentDets(:, i), nifd)) cycle
            end if
            norm_pert = norm_pert + real_sign(1) * real_sign(1)
            do
                j = j + 1
                if (DetBitEq(CurrentDets(:, i), ilut_list(:, j), nifd)) then
                    pert_ground_local(j) = real_sign(1)
                    exit
                end if
            end do
        end do

        ! The contributions to pops_norm are added in as the determinants are
        ! read in from the popsfile, which is only done on the root process, so
        ! brooadcast the value to other processors.
        call MPIBCast(pops_norm)
        pert_ground_local = pert_ground_local / sqrt(pops_norm)

        call MPISumAll(norm_pert, all_norm_pert)
        all_norm_pert = sqrt(all_norm_pert)

    end subroutine return_perturbed_ground_spec