ecal.for

德国人开发的地震处理分析软件

$DEBUG

      PROGRAM ECAL

      DIMENSION ACEL(40000)

      WRITE(*,*)
      WRITE(*,*)
      WRITE(*,*)
      WRITE(*,*) '     SE ESTAN CALCULANDO LOS ESPECTROS ELASTICOS...'
      WRITE(*,*)
      WRITE(*,*) '              CUANDO TERMINE EL CALCULO'
      WRITE(*,*) '            CIERRE LA VENTANA DOS Window'
      WRITE(*,*)
      WRITE(*,*)
      WRITE(*,*)

      OPEN(5,FILE='ACEL.DAT',STATUS='OLD')
      OPEN(6,FILE='SDAT.DAT',STATUS='OLD')
      OPEN(7,FILE='ASd.SAL',STATUS='UNKNOWN')
      OPEN(8,FILE='ASv.SAL',STATUS='UNKNOWN')
      OPEN(9,FILE='ASa.SAL',STATUS='UNKNOWN')

      READ(6,'(I8,F8.3,I8,3F8.3)')NACEL,DT,NNRES,DTSA,SAINI,XI

      DO I=1,NACEL
         READ(5,'(E12.4)') ACEL(I)
      END DO

      DO I=1,NNRES
         W=2.0*3.141592654/(SAINI+(I-1)*DTSA)
         WD=W*SQRT(1-XI*XI)
         XT0=0
         VT0=0
         SD=0
         DO K=1,NACEL
            A0=ACEL(K)
            PA=(ACEL(K+1)-ACEL(K))/DT
            D=-PA/W/W
            CC=(-A0+2*XI*PA/W)/W/W
            G1=XT0-CC
            G2=(VT0+XI*W*G1+PA/W/W)/WD
            XT0=EXP(-XI*W*DT)*(G1*COS(WD*DT)+G2*SIN(WD*DT))+CC+D*DT
            Vt0=EXP(-XI*W*DT)*(-G1*WD*SIN(WD*DT)+G2*WD*COS(WD*DT))
     1          -XI*W*EXP(-XI*W*DT)*(G1*COS(WD*DT)+G2*SIN(WD*DT))+D
            IF (SD .LT. ABS(XT0)) THEN
               SD=ABS(XT0)
            END IF
         END DO
         WRITE(7,'(E12.5)') SD
         WRITE(8,'(E12.5)') SD*W
         WRITE(9,'(E12.5)') SD*W*W
      END DO

      CLOSE (5)
      CLOSE (6)
      CLOSE (7)
      CLOSE (8)
      CLOSE (9)

      STOP
      END