MODULE DIAGNOSTIC

      USE opt
      USE retvparam

      IMPLICIT NONE
      REAL(DOUBLE) :: DOF

      CONTAINS

!----------------------------------------------------------------------
      SUBROUTINE DOFS( M, N, NLD, NLEV )

      INTEGER, INTENT(IN)                :: N, M, NLD, NLEV
      INTEGER                            :: I, J
      REAL(DOUBLE), DIMENSION(NLEV)      :: WR, WI
      REAL(DOUBLE), DIMENSION(N,N)       :: A, SM, IDEN
      REAL(DOUBLE), DIMENSION(N,M)       :: NM, G
      REAL(DOUBLE), DIMENSION(NLEV,NLEV) :: H, Q

!      INTEGER                            :: INFO
!      REAL(DOUBLE), DIMENSION(4*NLEV)    :: WORK


!  --- CALCULATE AVERAGING KERNEL NxN MATRIX
!  --- KSK is an nmax*nmax array
      !print *, m, nld, nlev
      !print *, N, 'N'
      IF (LM) THEN
         CALL MULT( G_LM, KHAT, A, N, M, N)
      ELSE
         CALL MULT( SHAT(:N,:N), KSK, A, N, N, N)
      END IF

!  --- CALCULATE TRACE OF KERNEL
      DOF = 0.D0
      DO I = 1, NLEV
         DOF = DOF + A(I+NLD,I+NLD)
      END DO
      !PRINT *, 'DOF : ', DOF

!  --- CONVERT TO VMR
      !DO I=1, NLEV
      !  DO J=1, NLEV
      !     !A(I=NLD,J+NLD) = A(I+NLD,J+NLD) * X(1,I) / X(1,J) ! * 1.0D12
      !     DOF = DOF + A(I+NLD,I+NLD)
      !ENDDO
      !ENDDO

      CALL FILEOPEN( 81, 1 )
      WRITE(81,*) NLEV
      DO I=1,NLEV
          WRITE(81,'(255ES26.18)') ( A(I+NLD,J+NLD), J=1, NLEV )
      ENDDO
      CALL FILECLOSE( 81, 1 )

!  --- CALCULATE MEASUREMENT ERROR : Se IS DIAGONAL, BUT FIRST GAIN
      G(:N,:M) = MATMUL( SHAT(:N,:N), KS(:N,:M) )
      CALL MULTDIAG( G(:N,:M), SE, NM, N, M )
      SM(:N,:N) = MATMUL( NM(:N,:M), TRANSPOSE( G(:N,:M) ))

      CALL FILEOPEN( 82, 2 )
      WRITE(82,*) NLEV
      DO I=1,NLEV
          WRITE(82,'(255ES26.18)') ( SM(I+NLD,J+NLD), J=1, NLEV )
      ENDDO
      CALL FILECLOSE( 82, 1 )

!  --- CALCULATE SMOOTHING ERROR
!  --- FILL DIAGONAL ELEMENTS OF IDEN, CREATE IDENTITY MATRIX
      IDEN = 0.0D0
      FORALL( I=1 : N ) IDEN(I,I) = 1.D0

      SM = MATMUL( (A-IDEN), SA(:N,:N) )
      SM = MATMUL( SM, TRANSPOSE( A-IDEN ))

      CALL FILEOPEN( 83, 2 )
      WRITE(83,*) NLEV
      DO I=1,NLEV
          WRITE(83,'(255ES26.18)') ( SM(I+NLD,J+NLD), J=1, NLEV )
      ENDDO
      CALL FILECLOSE( 83, 1 )

!  --- AK - SMOOTH
      CALL FILEOPEN( 84, 2 )
      WRITE(84,*) NLEV
      DO I=1,NLEV
          WRITE(84,'(255ES26.18)') ( A(I+NLD,J+NLD)-SM(I+NLD,J+NLD), J=1, NLEV )
      ENDDO
      CALL FILECLOSE( 84, 1 )

! --- COMMENT OUT TO USE LAPACK AND MAKE EIGENVECTORS
      RETURN

!  --- CALCULATE THE EIGENVALUES/VECTORS OF THE AVERAGING KERNEL
!  --- REQUIRES COMPILING WITH LAPACK
!  --- COPY A SINCE DGEEV WILL DESTROY CONTENTS
!  --- EIGENVALUE DECOMPOSITION
      H = A(NLD+1:NLEV,NLD+1:NLEV)
!      CALL DGEEV('N','V',NLEV, H, NLEV, WR, WI, 0, 1, Q, NLEV, WORK, 4*NLEV, INFO)
      !PRINT *, 'INFO :', INFO

!  --- EIGENVALUES OF AK
      CALL FILEOPEN( 85, 2 )
      WRITE(85,*) NLEV, 'first 2 rows are real & imag eigenvalues'
      WRITE(85,'(255ES26.18)') ( WR(J), J=1, NLEV )      ! REAL E-VALUE
      WRITE(85,'(255ES26.18)') ( WI(J), J=1, NLEV )      ! IMAGINARY #-VALUE
      DO I=1,NLEV
          WRITE(85,'(255ES26.18)') ( Q(I,J), J=1, NLEV )
      ENDDO
      CALL FILECLOSE( 85, 1 )

      END SUBROUTINE DOFS

      END MODULE DIAGNOSTIC