PROGRAM BT5_14
	PARAMETER (NMAX=100,MMAX=10)
	DIMENSION X(NMAX*MMAX),R(MMAX*MMAX),RT(MMAX*MMAX)
	CHARACTER VARNAME(4)*2
	INTEGER N,M
	DATA VARNAME /'X1','X2','X3','X4'/

	OPEN(1,FILE='C:\MASTER\G_TR\Bai Tap TK trong KH\New\BANG17.TXT')
	READ(1,*)
	READ(1,*)
	READ(1,*) N,M
	DO I=1,N
	   L=(M-1)*N+I
	   READ(1,*) TMP, (X(J),J=I,L,N)
	ENDDO
	
	CALL CORRE(X,N,M,R)
	WRITE(*,*)'   MA TRAN TUONG QUAN '
	WRITE(*,'("       ",4A7)') (VARNAME(J),J=1,M)
	JK=0
	DO J=1,M
	   JK=(J-1)*M+1
	   WRITE(*,'(A7,4F7.4)') VARNAME(J),(R(K),K=JK,JK+M-1)
	ENDDO
	DO K=2,4
	   J=1
	   RT=R
	   R1K=HSTQR(RT,M,J,K)
	   WRITE(*,'("HE SO TQ RIENG GIUA ",A2,"-",A2,"=",F7.4)')
     &         VARNAME(J),VARNAME(K),R1K
	ENDDO
	J=1
	RT=R
	R1K=HSTQB(RT,M,J)
	WRITE(*,'("HE SO TQ BOI GIUA ",A2," VA ",A2,"...",A2,"=",
     &      F7.4)') VARNAME(J),VARNAME(J+1),VARNAME(M),R1K

	END

!+++++++++++++++++++++++++++++++++

	SUBROUTINE CORRE(X,N,M,R)
! CHUONG TRINH NAY TINH MA TRAN TUONG QUAN CHUAN HOA
! CUA TAP M BIEN TU SO LIEU BAN DAU CO DUNG LUONG N
! INPUT: + X MANG MOT CHIEU KICH THUOC N*M LUU TRU SO
!          LIEU BAN DAU DANG MA TRAN N HANG M COT
!          (X(1,1), X(2,1),...,X(N,1),X(1,2),...)
!        + N DUNG LUONG MAU
!        + M SO BIEN
! OUTPUT: + R MANG MOT CHIEU KICH THUOC M*M LUU TRU MA TRAN
!            TUONG QUAN CHUAN HOA CUA M BIEN
!
	DIMENSION X(N*M),XX(N,M),R(M*M),TB(M),SX(M)
	K=0
	DO J=1,M
	   TB(J)=0.0
	   SX(J)=0.0
	   DO I=1,N
	      K=K+1
	      TB(J)=TB(J)+X(K)
	      SX(J)=SX(J)+X(K)*X(K)
	      XX(I,J)=X(K)
	   ENDDO
	   TB(J)=TB(J)/REAL(N)
	   SX(J)=SQRT(SX(J)/REAL(N)-TB(J)*TB(J))
	ENDDO
	JK=0
	DO J=1,M
	   DO K=1,M
	      JK=JK+1
	      R(JK)=0.0
	      DO I=1,N
	         R(JK)=R(JK)+XX(I,J)*XX(I,K)
	      ENDDO
	      R(JK)=R(JK)/REAL(N)-TB(J)*TB(K)
	      R(JK)=R(JK)/(SX(J)*SX(K))
	   ENDDO
	ENDDO
	RETURN
	END

	FUNCTION HSTQR(R,M,J,K)
! HAM NAY TINH HE SO TUONG QUAN RIENG GIUA BIEN XJ VA XK
! KHI XET DONG THOI M BIEN X1,X2,...,XM
! INPUT: + R MANG MOT CHIEU CHUA MA TRAN TUONG QUAN CHUAN HOA
!          (KICH THUOC M*M ==> R(M*M)
!        + M KICH THUOC MA TRAN TUONG QUAN, DONG THOI CUNG LA
!          SO BIEN DUOC XET
!        + J,K CHI SO BIEN DUOC TINH TUONG QUAN RIENG
! OUTPUT: HE SO TUONG QUAN RIENG GIUA BIEN THU J VA THU K
!
	DIMENSION R(M*M),RT(M*M)
	RT=R
	RJK=PPDS(RT,M,J,K)
	RT=R
	RJJ=PPDS(RT,M,J,J)
	RT=R
	RKK=PPDS(RT,M,K,K)
	HSTQR=-RJK/(SQRT(RJJ*RKK))
	RETURN
	END

	FUNCTION HSTQB(R,M,J)
! HAM NAY TINH HE SO TUONG QUAN BOI GIUA BIEN XJ VA TAP HOP
!        M-1 BIEN X1,X2,X(J-1),X(J+1),..,XM 
! INPUT: + R MANG MOT CHIEU CHUA MA TRAN TUONG QUAN CHUAN HOA
!          (KICH THUOC M*M ==> R(M*M)
!        + M KICH THUOC MA TRAN TUONG QUAN, DONG THOI CUNG LA
!          SO BIEN DUOC XET
!        + J CHI SO BIEN DUOC TINH TUONG QUAN BOI VOI TAP M-1 
!            BIEN CON LAI
! OUTPUT: HE SO TUONG QUAN BOI GIUA BIEN THU J VA CAC BIEN KHAC
!
	DIMENSION R(M*M),RT(M*M)
	RT=R
	RJJ=PPDS(RT,M,J,J)
	RT=R
	RR=DET(RT,M)
	HSTQB=SQRT(1-RR/RJJ)
	RETURN
	END

!
	FUNCTION PPDS(A,N,IR,JC)
!
! HAM NAY TINH PHAN PHU DAI SO CUA PHAN TU B(IR,JC)
! (HANG IR, COT JC) CUA MA TRAN B(NxN) MA CAC PHAN TU
! CUA NO DUOC LUU TRONG MANG A(N*N) THEO QUI CACH COT
! TRUOC DONG SAU
! INPUT: + MANG A(N*N) CHUA MA TRAN DAU VAO CUA B
!        + N KICH THUOC CUA B
!        + IR CHI SO HANG
!        + JC CHI SO COT
! OUTPUT: PHAN PHU DAI SO CUA PHAN TU HANG IR COT JC
!
	DIMENSION A(N*N), B(N,N)
	K=0
	DO J=1,N
	   DO I=1,N
	      K=K+1
	      B(I,J)=A(K)
	   ENDDO
	ENDDO
	
	K=0
	DO J=1,N
	   IF (J.NE.JC) THEN
	   DO I=1,N
	      IF (I.NE.IR) THEN
	         K=K+1
	         A(K)=B(I,J)
	      ENDIF
	   ENDDO
	   ENDIF
	ENDDO
	
	PPDS=DET(A,N-1)*(-1)**(IR+JC)
	RETURN
	END


	FUNCTION DET(X,N)
! HAM NAY TINH DINH THUC CUA MA TRAN A(N,N)
! INPUT: + MANG X DO DAI N*N CHUA CAC PHAN TU CUA MA TRAN A
!        + N KICH THUOC MA TRAN A
!   CHU Y: X(1)=A(1,1), X(2)=A(2,1),..., X(N)=A(N,1),X(N+1)=A(1,2),
!          X(N+2)=A(2,2),...
! OUTPUT: DINH THUC CUA A
!
	PARAMETER (EP=1.0E-6)
	DIMENSION X(N*N),A(N,N)
	K=0
	DO J=1,N
	   DO I=1,N
	      K=K+1
	      A(I,J)=X(K)
	   ENDDO
	ENDDO
	D=1.0
	N1=N-1
	DO 10 K=1,N1
	   AM=0.0
	   DO 11 I=K,N
	      T=A(I,K)
	      IF (ABS(T).GE.ABS(AM)) THEN
	         AM=T
	         J=I
	      ENDIF
  11     CONTINUE
         IF (ABS(AM).LE.EP) THEN
	      DT=0.0
	      DET=DT
	      RETURN
	   ELSE
	      IF (J.NE.K) THEN
	         D=-D
	         DO I=K,N
	            T=A(J,I)
	            A(J,I)=A(K,I)
	            A(K,I)=T
	         ENDDO
	      ENDIF
	   ENDIF
	   M=K+1
	   DO I=M,N
	      T=A(I,K)/AM
	      DO J=M,N
	         A(I,J)=A(I,J)-T*A(K,J)
		  ENDDO
	   ENDDO
	   D=D*A(K,K)
  10  CONTINUE
      DT=D*A(N,N)
	DET=DT
	RETURN
	END