GETIRREPDECOMP Function

    LOGICAL FUNCTION GETIRREPDECOMP(CHARS, IRREPCHARS, NIRREPS, NROT, IDECOMP, CNORM, TAbelian)
        IMPLICIT NONE
        INTEGER NIRREPS, NROT
        TYPE(Symmetry) IDECOMP
        complex(dp) IRREPCHARS(NROT, NIRREPS), CHARS(NROT)
        real(dp) CNORM, NORM, DIFF
        complex(dp) TOT
        INTEGER I, J
        logical TAbelian
        character(*), parameter :: this_routine = 'GETIRREPDECOMP'
        if (TAbelian) then
            ! We shouldn't be here!  Using symmetry "quantum" numbers
            ! rather than irreps.
            call stop_all(this_routine, "Should not be decomposing irreps with Abelian sym")
        end if
        IDECOMP%s = 0
!,. First check norm of this state
        CNORM = 0
        DO J = 1, NROT
            CNORM = CNORM + real(CONJG(CHARS(J)) * CHARS(J), dp)
        end do
        DO I = 1, NIRREPS
            TOT = 0
            DO J = 1, NROT
                TOT = TOT + CONJG(IRREPCHARS(J, I)) * CHARS(J)
            end do
            IF (ABS(TOT) >= 1.0e-2_dp) THEN
!   Calculate the normalization of the state I which matches (if it's an irrep, this will be 1)
                NORM = 0
                DO J = 1, NROT
                    NORM = NORM + real(CONJG(IRREPCHARS(J, I)) * IRREPCHARS(J, I), dp)
                end do
!               write(stdout,*) "IRREP ",I,(TOT+0.0_dp)/NORM
!                CALL writechars(stdout,CHARS,NROT,"REP   ")
!                CALL writechars(stdout,IRREPCHARS(1,I),NROT,"IRREP ")
                DIFF = ABS(TOT - NINT(ABS(TOT / NORM)) * NORM)
                IF (DIFF >= 1.0e-2_dp .AND. abs(CNORM - NROT) < 1.0e-12_dp) THEN
!   The given representation CHARS has fewer irreps in it than the one in IRREPCHARS, and is an irrep
!   Hurrah!  Remove it from the one in IRREPCHARS, and keep on going)
!                  DO J=1,NROT
!                    IRREPCHARS(J,I)=IRREPCHARS(J,I)-CHARS(J)*TOT/CNORM
!                  end do
!                  CALL writechars(stdout,IRREPCHARS(1,I),NROT,"NOW   ")
                else if (DIFF < 1.0e-2_dp) THEN
!   We've found an (ir)rep which is wholly in CHARS
                    IDECOMP%s = IBSET(IDECOMP%s, I - 1)
                    CNORM = 0
                    DO J = 1, NROT
                        CHARS(J) = CHARS(J) - (IRREPCHARS(J, I) * TOT) / NORM
                        CNORM = CNORM + real(CONJG(CHARS(J)) * CHARS(J), dp)
                    end do
                end if
            end if
        end do
        GETIRREPDECOMP = .FALSE.
        DO J = 1, NROT
            IF (ABS(CHARS(J)) > 1.0e-2_dp) GETIRREPDECOMP = .TRUE.
        end do
    END FUNCTION GETIRREPDECOMP