最大熵谱法分析程序 总是数组越界

CUoldandnew

这是在站里下载的魏凤英老师的现代气候诊断与统计分析的Fortran源程序中关于最大熵谱分析的程序，但是我按照其中给的数据运行，总是出现数组越界的情况，求各位给看看，哪里出现的问题

C       THIS IS A PROGRAM OF MAXIMUM ENTROPY SPECTRUM
C       ANALYSIS (MESA) FOR A TIME SERIES
        PROGRAM MESA
        PARAMETER (N=44, M=22)
        DIMENSION X(N),T(N),SE(N,M)
        DIMENSION B1(N),B2(N),A(N),AA(N),FI(M),SIGMA(N)
        CHARACTER*9 FNAME1,FNAME2
C       ***************************************************
C       * N: SAMPLE SIZE OF THE TIME SERIES               *
C       * M: MAXIMUM LAG TIME LENGTH                      *
C       * X(N): TIME SERIES                               *
C       * T(N):LENGTH OF PERIOD                           *
C       * SE(N,M): MAXIMUM ENTROPY SPECTRUM VALUES        *
C       ***************************************************
        WRITE(*,*) 'PLEASE ENTER THE DATA FILE NAME:'
        READ(*,'(A)') FNAME1
        WRITE(*,*) 'PLEASE ENTER THE OUTPUT FILE NAME:'
        READ(*,'(A)') FNAME2
        NW2=1
        OPEN (7,FILE=FNAME1)
        READ (7,*,END=9)(X(I),I=1,N)
9       CONTINUE
        CLOSE (7)
        CALL SMESA(X,N,M,NW2,K0,T,SE,FI,SIGMA,A,AA,B1,B2)
        OPEN (8,FILE=FNAME2)
        WRITE(8,80)
  80    FORMAT(20X,'*******MAXIMUM ENTROPY SPECTRUM ANALYSIS*******'/)
        WRITE (8,90)
        WRITE (8,100) (L,T(L),SE(L,K0), L=1,M)
        CLOSE(8)
  90    FORMAT (1X,7X,1HL,4X,1HT,10X,2HSE/,
     &          8X,23H-----------------------)
100    FORMAT (1X,I8,F7.2,E15.7)
        STOP
        END

CCC------------------------------END------------------------------------CCC

        SUBROUTINE SMESA(X,N,K,NW,K0,T,SE,FI,SIGMA,A,AA,B1,B2)
        DIMENSION X(N),T(N),SE(N,K)
        DIMENSION B1(N),B2(N),A(N),AA(N),FI(K),SIGMA(N)
C
C----------------------------------------------------------------------C
C
        IF (NW.EQ.0) GOTO 10
        WRITE(*,'(/9X,1Hk,6X,5HFI(K))')
   10   S=0.
        DO 50 I=1,N
   50   S=S+X(I)
        S=S/N
        DO 60 I=1,N
        X(I)=X(I)-S
   60   CONTINUE
        S=0.
        DO 11 I=1,N
   11   S=S+X(I)*X(I)
        SIGMA(1)=S/N
        M=1
        B1(1)=X(1)
        B2(N-1)=X(N)
        DO 12 J=2,N-1
        B1(J)=X(J)
   12   B2(J-1)=X(J)
        IF (NW.EQ.0) GOTO 301
        WRITE(*,'(I10,F12.3)') M, SIGMA(1)
  301   GOTO 101
  102   M=M+1
        DO 13 J=1,M-1
   13   AA(J)=A(J)
        DO 14 J=1,N-M
        B1(J)=B1(J)-AA(M-1)*B2(J)
   14   B2(J)=B2(J+1)-AA(M-1)*B1(J+1)
  101   COM=0.
        DEN=0.
        DO 15 J=1,N-M
        COM=COM+B1(J)*B2(J)
   15   DEN=DEN+B1(J)**2+B2(J)**2
        A(M)=2.0*COM/DEN
        SIGMA(M+1)=SIGMA(M)*(1-A(M)**2)
        IF(M.EQ.1)GOTO 102
  103   DO 16 J=1,M-1
        A(J)=AA(J)-A(M)*AA(M-J)
   16   CONTINUE
        FI(M) = (FLOAT(N+M)/FLOAT(N-M))*SIGMA(M+1)
        IF (FI(M-1).lt.FI(M-2).and.FI(M-1).lt.FI(M)) then
        K0=M-1
        GOTO 555
          ELSE
            CONTINUE
          ENDIF
        IF (NW.EQ.0) GOTO 300
        WRITE(*,'(I10,F12.3)') M, FI(m)
  300   IF (M.EQ.K) GOTO 555
        DO 19 I=1,K
        F1=0.
        F2=0.
        DO 18 J=1,M
        E=2.0*3.1416*FLOAT(J)*FLOAT(I)/FLOAT(N)
        F1=F1+A(J)*COS(E)
   18   F2=F2+A(J)*SIN(E)
        SE(I,M) = SIGMA(M+1)/((1.0-F1)**2+F2**2)
   19   CONTINUE
        IF (M.LT.K) GOTO 102
  555   DO 150 I=1,K
        T(I) = 2.*K/float(I)
  150   CONTINUE
        IF (NW.EQ.0) GOTO 999
        WRITE (*,'(/,12x,4Hk0 =,i4/)') K0
  999   RETURN
        END
CCC-----------------------------END----------------------------------

言深深

楼主觉得这么长的程序会有人帮你看么？更不用说能够看出来的有多少了···

建议自己debug一下，然后定位到出错行。要是不知道debug在哪里的话，就按F5···

CUoldandnew

哦，好的，没太发过帖子，下次不会这样发了，谢谢您

CUoldandnew

哦，好的，没太发过帖子，下次不会这样发了，谢谢您

tszhaoweiwen

楼主 我最近也用魏老师的这个程序 出现一样问题 楼主解决了吗？

tszhaoweiwen

如果大家都遇到这个问题，可以这样解决，魏老师的程序比较旧了，CVF6可能出现向上不兼容问题，本人改用fortran4.0就可以解决了 希望对同样有问题的人有帮助

冰块馒头

楼主的问题解决了吗？最近也在学习，遇到了相同的问题