**==aa0001.f    processed by SPAG 4.03F  at 13:27 on  2 Mar 1994
C --- QUARTERLY SEASONAL ADJUSTMENT PROGRAM
C     U.S. CENSUS BUREAU
C     ECONOMIC RESEARCH AND ANALYSIS DIVISION
C     WRITTEN BY NORMAN BAKKA
C
C --- FEBRUARY 1, 1968
C --- INPUT ROUTINE
C  (1)READS AND EDITS HEADER CARD
C  (2)READS AND EDITS INPUT SERIES
C  (3)CALLS QSAP, WHICH IS THE MAIN SEASONAL ADJUSTMENT PROGRAM
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,Ib,Ibreak,Ic,ic2,icer,Icq,Icyr,Ie,ifmt,Ifor,Igrz,ijk,
     &        indec,inm,inx,inx2,Io12,Iovr,Iq
      INTEGER Irang,irerr,irs,Irun,issap,Itbl,itblon,itc1,Itchar,Itfmt,
     &        Ixerr,Iyr,iz,j,kji,Kp,Kpage,Kpropt,Kswerr,Kswt
      INTEGER Ktcopt,kyr,L0,l2,Len,Len2,Length,lfq,llq,lm,Lmon,Lstob,
     &        Lstyr,Lswv,Lterm,ltss,Lyr,Mach,max,Meas
      INTEGER myr,Myr2,ncmax,Ncol,nd,Nerr,newpg,nl,Nlen,Np,Nr,Ns1,Nsea,
     &        Nterm,Numtbl,Nw,Nyear,Y2KLIM
      DOUBLE PRECISION ox,Rv,Sd5,Sigml,Sigmu,Stci,Vc,Vi,Vs,Xorg,Xtc
      PARAMETER(Y2KLIM=45)
C*** End of declarations inserted by SPAG
      REAL siglim(2)
      DIMENSION Xorg(125),ox(125)
      CHARACTER*1 reset,stop,sspan,seaadj,rewind,addmul,prtout,tsma,
     &            sma(4),asciix,chart,trend,Ch,prompt,y2kset
      CHARACTER*3 quar(4),Tbl
      CHARACTER*4 csma(4),vsk(5),blank
      CHARACTER*6 madj(3),aadj(3),summ(4),seas(4),prnt(4),adj(3),run(4),
     &            Ser,mser,tser
      CHARACTER*18 Chper
      CHARACTER*24 Cname
      CHARACTER*25 crerr(6),cre
      CHARACTER*40 nf4,fmt,Tf4
      CHARACTER*80 Title,ofile,otemp,dfile,dtemp,ifile,tfile,efile,f1,
     &             Tblfil
      INTEGER begob,Begyr,Outdec,ssyear,Add,ascibm
      LOGICAL lerout,lerdat
      DIMENSION nl(4)
      COMMON /ccci  / Stci(125)
      COMMON /cee   / Ic,Icyr,Icq
      COMMON /cfilt / Ktcopt,Lterm(4),Nterm
      COMMON /chame / Cname(4),Ch(4),Chper(2)
      COMMON /chrtm / Mach
      COMMON /com4b / Iyr,Iq,Nsea,Lyr,Lmon
      COMMON /inpt  / Xorg,Ib,Ie,Meas,Kp,Kpropt
      COMMON /iopt  / Ibreak,Irang
      COMMON /itcom / Itchar,Itfmt,Itbl(10),Ifor,Numtbl,Iovr,Io12
      COMMON /lzero / L0
      COMMON /ngrz  / Igrz,Myr2,Nw
      COMMON /ns4   / Nlen,Ncol,Ns1,Len,Len2
      COMMON /rite  / Np,Nyear,Title,Kpage,Outdec,Ser
      COMMON /sigm  / Sigml,Sigmu,Sd5(125)
      COMMON /tadj  / Add
      COMMON /tcom  / Tf4,Tblfil,Tbl(10)
      COMMON /xerr  / Nr,Nerr,Kswerr,Irun,Ixerr
      COMMON /xk    / Kswt,Lswv,Begyr
C  *****  CHANGE TO INCORPORATE M & Q DIAGNOSTICS
      COMMON /mq2   / Vi,Vc,Vs,Rv,Length,Lstob
      COMMON /mq4   / Xtc(125),Lstyr
C  *****  END OF CHANGE
      DATA newpg/'1'/
      DATA crerr/'NLX11Q                   ',
     &     'TITLE                    ','DATA FILE NAME           ',
     &     'OUTPUT FILE NAME         ','X-11Q.2 COMMAND FILE NAME',
     &     'X-11Q.2 ERROR FILE NAME  '/
      DATA(madj(i),i=1,3)/'MULTIP','LICATI','VE    '/
      DATA(aadj(i),i=1,3)/'      ','ADDITI','VE    '/
      DATA(summ(i),i=1,4)/'SUMMAR','Y MEAS','URES  ','      '/
      DATA(seas(i),i=1,4)/'SEASON','AL ADJ','USTMEN','T     '/
      DATA(prnt(i),i=1,4)/' SHORT','STNDRD','  LONG','  FULL'/
      DATA(quar(i),i=1,4)/'1ST','2ND','3RD','4TH'/
      DATA max,newpg/125,'1'/
      DATA nl/24,24,32,44/
      DATA blank/'    '/
      DATA vsk/'3X1 ','3X3 ','3X5 ','3X9 ','STA '/
      NAMELIST /nlx11q/ Ser,fmt,indec,begob,Begyr,Lstob,Lstyr,rewind,
     & addmul,seaadj,prtout,Outdec,chart,ascibm,siglim,sspan,ssyear,
     & trend,sma,reset,Numtbl,y2kset,prompt,stop
      Nr=5
      Np=6
      Nw=0
      nd=0
      Io12=0
      reset='Y'
      itblon=0
      Irun=1
      Ixerr=0
      irs=0
      lerout=.false.
      lerdat=.false.
      icer=5
   10 PRINT *,'ENTER CON (FOR INTERACTIVE PROCESSING) OR X-11Q.2 ',
     &      'COMMAND FILE NAME '
      PRINT *,' > '
      READ ,ifile
      IF(ifile(1:3).ne.'CON'.and.ifile(1:3).ne.'con')THEN
       Nr=10
       tfile=ifile
       OPEN(Nr,FILE=ifile,STATUS='OLD',ERR=180)
      END IF
      icer=6
      Nerr=25
   20 IF(Nr.eq.10)THEN
       PRINT *,'ENTER X-11Q.2 ERROR FILE NAME OR @ (ERRORS STORED IN',
     &       ' C:\X11Q.ERR) :'
       PRINT *,' > '
       READ ,efile
       IF(efile.eq.'@')efile='C:\X11Q.ERR'
       tfile=efile
       OPEN(25,FILE=efile,ERR=190)
      END IF
   30 IF(reset.eq.'Y')THEN
       stop='N'
       begob=1
       Begyr=0
       Lstob=4
       Lstyr=0
       Title=' '
       ofile=' '
       dfile=' '
       indec=0
       addmul='M'
       seaadj='Y'
       prtout='D'
       Ser='      '
       siglim(1)=1.5
       siglim(2)=2.5
       trend='S'
       chart='N'
       asciix='Y'
       y2kset='Y'
       rewind='N'
       fmt='1'
       Outdec=0
       sma(1)='D'
       DO i=2,4
        sma(i)=' '
       END DO
       sspan='N'
       ssyear=0
       prompt='N'
      END IF
C  --  READ MANDITORY NAMELIST
      IF(Nr.eq.10)THEN
       READ(Nr,nlx11q,ERR=130,END=140)
       IF(stop.ne.'N'.and.stop.ne.'n')GO TO 200
      ELSE
       PRINT *,' '
       PRINT *,'  CURRENT SETTINGS FOR NLX11Q NAMELIST ARE : '
       PRINT nlx11q
   40  PRINT *,' '
       PRINT *,'  ENTER VARIABLES AND THEIR SETTINGS FOR NLX11Q',
     &       ' NAMELIST : '
       READ(*,nlx11q,ERR=150)
       IF(stop.ne.'N'.and.stop.ne.'n')GO TO 200
       PRINT nlx11q
       IF(prompt.eq.'Y'.or.prompt.eq.'y')CALL prop(1,*40)
      END IF
   50 Kswerr=0
      Kswt=4
      IF(irs.ne.2)THEN
       irs=9
       IF(reset.eq.'Y'.or.reset.eq.'y')irs=0
       IF(reset.eq.'N'.or.reset.eq.'n')irs=1
       IF(irs.eq.9)THEN
        CALL inerr(reset,irs,*200)
        IF(Nr.eq.10)WRITE(Nerr,1070)'RESET '
        IF(Nr.eq.5)PRINT 1070,'RESET '
       END IF
      END IF
      IF(begob.lt.1.or.begob.gt.4)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1010)'BEGOB ',begob,1,4
 1010  FORMAT(10X,'THE ENTRY FOR ',A6,' IS ',I6,', WHICH IS NOT ',
     &        'DEFINED',/,10X,'THE ENTRY SHOULD BE BETWEEN ',I6,' AND ',
     &        I6)
       IF(Nr.eq.5)PRINT 1010,'BEGOB ',begob,1,4
      END IF
      IF(Lstob.lt.1.or.Lstob.gt.4)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1010)'LSTOB ',Lstob,1,4
       IF(Nr.eq.5)PRINT 1010,'LSTOB ',Lstob,1,4
      END IF
      IF(.not.(y2kset.eq.'Y'.or.y2kset.eq.'y'.or.y2kset.eq.'N'.or.
     &         y2kset.eq.'n'))THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1021)'Y2KSET',y2kset
 1021  FORMAT(10X,'THE ENTRY FOR ',A6,' IS ',I6,', WHICH IS NOT ',
     &        'DEFINED',/,10X,'THE ENTRY SHOULD BE EITHER Y OR N')
       IF(Nr.eq.5)PRINT 1021,'Y2KSET',y2kset
      END IF
      IF(Begyr.le.0)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1022)'BEGYR ',Begyr
 1022  FORMAT(10X,'THE ENTRY FOR ',A6,' IS ',I6,', WHICH IS NOT ',
     &        'DEFINED',/,10X,'THE ENTRY SHOULD BE GREATER THAN ZERO')
       IF(Nr.eq.5)PRINT 1022,'BEGYR ',Begyr
      ELSE IF(Begyr.lt.1900)THEN
       IF(y2kset.eq.'Y'.or.y2kset.eq.'y')THEN
        IF(Begyr.le.Y2KLIM)THEN
         Begyr=2000+Begyr
        ELSE
         Begyr=1900+Begyr
        END IF
       ELSE
        Begyr=1900+Begyr
       END IF 
      END IF
      IF(Lstyr.le.0)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1022)'LSTYR ',Lstyr
       IF(Nr.eq.5)PRINT 1022,'LSTYR ',Lstyr
      ELSE IF(Lstyr.lt.1900)THEN
       IF(y2kset.eq.'Y'.or.y2kset.eq.'y')THEN
        IF(Lstyr.le.Y2KLIM)THEN
         Lstyr=2000+Lstyr
        ELSE
         Lstyr=1900+Lstyr
        END IF
       ELSE
        Lstyr=1900+Lstyr
       END IF 
      END IF
C --- CHECK LENGTH OF SERIES
      Ie=(Lstyr-Begyr)*4+Lstob
      Length=Ie-begob+1
      IF(Length.lt.12)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1020)
 1020  FORMAT(10X,'THE SERIES IS LESS THAN 12 QUARTERS LONG')
       IF(Nr.eq.5)PRINT 1020
      ELSE IF(Ie.gt.max)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1030)
 1030  FORMAT(10X,'SERIES EXCEEDS 30 CALENDAR YEARS IN LENGTH')
       IF(Nr.eq.5)PRINT 1030
      END IF
      kyr=Begyr
C     CHECK NO. OF DECIMAL PLACES IN INPUT
      IF(indec.lt.0.or.indec.gt.5)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1010)'INDEC ',indec,0,5
       IF(Nr.eq.5)PRINT 1010,'INDEC ',indec,0,5
      END IF
      IF(fmt.eq.'2')THEN
       ifmt=0
       WRITE(nf4,1040)indec
 1040  FORMAT('(4F12.',I1,',24X,I2,A6)')
      ELSE IF(fmt.eq.'1')THEN
       ifmt=0
       WRITE(nf4,1050)indec
 1050  FORMAT('(4F18.',I1,',I2,A6)')
      ELSE IF(fmt.eq.' ')THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1060)
 1060  FORMAT(10X,'FORMAT FOR DATA INPUT NOT ENTERED')
       IF(Nr.eq.5)PRINT 1060
      ELSE
       ifmt=1
       nf4=fmt
      END IF
      IF(Outdec.lt.0.or.Outdec.gt.3)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1010)'OUTDEC',Outdec,0,3
       IF(Nr.eq.5)PRINT 1010,'OUTDEC',Outdec,0,3
      END IF
      Add=9
      IF(addmul.eq.'M'.or.addmul.eq.'m')Add=0
      IF(addmul.eq.'A'.or.addmul.eq.'a')Add=1
      IF(Add.eq.9)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1070)'ADDMUL'
 1070  FORMAT(10X,'THE ENTRY FOR ',A6,' IS NOT DEFINED')
       IF(Nr.eq.5)PRINT 1070,'ADDMUL'
      END IF
      Meas=9
      IF(seaadj.eq.'Y'.or.seaadj.eq.'y')Meas=1
      IF(seaadj.eq.'N'.or.seaadj.eq.'n')Meas=2
      IF(Meas.eq.9)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1070)'SEAADJ'
       IF(Nr.eq.5)PRINT 1070,'SEAADJ'
      END IF
      Kpropt=9
      IF(prtout.eq.'S'.or.prtout.eq.'s')Kpropt=0
      IF(prtout.eq.'D'.or.prtout.eq.'d')Kpropt=1
      IF(prtout.eq.'L'.or.prtout.eq.'l')Kpropt=2
      IF(prtout.eq.'F'.or.prtout.eq.'f')Kpropt=3
      IF(Kpropt.eq.9)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1070)'PRTOUT'
       IF(Nr.eq.5)PRINT 1070,'PRTOUT'
      END IF
      Kp=Kpropt+1
C --- TEST SIGMA LIMITS FOR GRADUATING EXTREMES
      IF(siglim(2).le.0.)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1080)'UPPER'
 1080  FORMAT(A5,' SIGMA LIMIT IS NEGATIVE OR ZERO')
       IF(Nr.eq.5)PRINT 1080,'UPPER'
      END IF
      IF(siglim(1).le.0)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1080)'LOWER'
       IF(Nr.eq.5)PRINT 1080,'LOWER'
      END IF
      IF(siglim(2).lt.siglim(1))THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1090)
 1090  FORMAT('UPPER SIGMA LIMIT IS LESS THAN LOWER SIGMA LIMIT')
       IF(Nr.eq.5)PRINT 1090
      END IF
      Sigml=dble(siglim(1))
      Sigmu=dble(siglim(2))
      Lswv=0
C --- TEST SEASONAL FACTOR CURVE MOVING AVERAGES
      DO i=1,4
       tsma=sma(i)
       IF(tsma.eq.' ')tsma=sma(1)
       Lterm(i)=9
       IF(tsma.eq.'D'.or.tsma.eq.'d')Lterm(i)=0
       IF(tsma.eq.'1')Lterm(i)=1
       IF(tsma.eq.'3')Lterm(i)=2
       IF(tsma.eq.'5')Lterm(i)=3
       IF(tsma.eq.'9')Lterm(i)=4
       IF(tsma.eq.'S'.or.tsma.eq.'s')Lterm(i)=5
       IF(Lterm(i).eq.9)THEN
        CALL inerr(reset,irs,*200)
        IF(Nr.eq.10)WRITE(Nerr,1100)i
 1100   FORMAT(10X,'SEASONAL MOVING AVERAGE SELECTED FOR MONTH ',I2,
     &         ' IS UNDEFINED')
        IF(Nr.eq.5)PRINT 1100,i
       ELSE IF(Lterm(i).eq.0)THEN
        csma(i)=blank
       ELSE
        Lswv=1
        j=Lterm(i)
        csma(i)=vsk(j)
       END IF
      END DO
      Ktcopt=9
      IF(trend.eq.'S'.or.trend.eq.'s')Ktcopt=0
      IF(trend.eq.'5')Ktcopt=1
      IF(trend.eq.'7')Ktcopt=2
      IF(Ktcopt.eq.9)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1070)'TREND '
       IF(Nr.eq.5)PRINT 1070,'TREND '
      END IF
      issap=9
      Ibreak=0
      IF(sspan.eq.'N'.or.sspan.eq.'n')issap=0
      IF(sspan.eq.'A'.or.sspan.eq.'a')issap=1
      IF(sspan.eq.'S'.or.sspan.eq.'s')THEN
       issap=1
       Ibreak=1
      END IF
      IF(issap.eq.9)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1070)'SSPAN '
       IF(Nr.eq.5)PRINT 1070,'SSPAN '
      END IF
C --- TEST GRAPHICS OPTIONS
      Igrz=9
      IF(chart.eq.'N'.or.chart.eq.'n')Igrz=0
      IF(chart.eq.'S'.or.chart.eq.'s')Igrz=1
      IF(chart.eq.'A'.or.chart.eq.'a')Igrz=2
      IF(Igrz.eq.9)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1070)'CHART '
       IF(Nr.eq.5)PRINT 1070,'CHART '
      END IF
      Mach=9
      IF(asciix.eq.'Y'.or.asciix.eq.'y')Mach=0
      IF(asciix.eq.'N'.or.asciix.eq.'n')Mach=1
      IF(Mach.eq.9)THEN
       CALL inerr(reset,irs,*200)
       IF(Nr.eq.10)WRITE(Nerr,1070)'ASCIIX'
       IF(Nr.eq.5)PRINT 1070,'ASCIIX'
      END IF
      IF(issap.ne.0)THEN
       ltss=0
       DO i=1,4
        IF(Lterm(i).eq.0.and.ltss.lt.3)THEN
         ltss=3
        ELSE IF(Lterm(i).gt.ltss)THEN
         ltss=Lterm(i)
        END IF
       END DO
       IF(ltss.eq.5)THEN
        Nlen=Length-12
        Ncol=4
        IF(Nlen.gt.48)THEN
         Nlen=48
         IF(Lstob.lt.4)Nlen=44
        ENDIF
        IF(Nlen.lt.12)THEN
         Ncol=((Length-12)/4)+1
         IF(Ncol.le.1)THEN
          Ncol=0
          issap=0
          IF(Nr.eq.10)THEN
           WRITE(Nerr,1110)
           IF(Ixerr.eq.0)Ixerr=1
          END IF
          IF(Nr.eq.5)PRINT 1110
          GO TO 60
         END IF
         Nlen=Length-4*(Ncol-1)
        END IF
       ELSE
        Nlen=nl(ltss)
        ncmax=((Length-Nlen)/4)+1
        IF(ncmax.gt.4)THEN
         Ncol=4
        ELSE IF(ncmax.lt.2)THEN
         Ncol=0
         issap=0
         IF(Nr.eq.10)THEN
          IF(Ixerr.eq.0)Ixerr=1
          WRITE(Nerr,1110)
         END IF
         IF(Nr.eq.5)PRINT 1110
 1110    FORMAT('    NOT ENOUGH DATA TO PERFORM SLIDING SPANS ',
     &          'ANALYSIS')
         GO TO 60
        ELSE
         Ncol=ncmax
        END IF
       END IF
       Len=Nlen+(Ncol-1)*4
       IF(Length.eq.Len)THEN
        Iq=begob
        Iyr=Begyr
       ELSE
        l2=Length-(Nlen+(Ncol-1)*4)
        Iyr=Begyr+l2/4
        Iq=begob+mod(l2,4)
        IF(Iq.gt.4)THEN
         Iq=Iq-4
         Iyr=Iyr+1
        END IF
        IF(ltss.lt.5)THEN
         IF((begob.eq.1.and.Iq.ne.1).or.(begob.lt.Iq.and.Iyr.eq.Begyr))
     &      THEN
          Nlen=(Iq-begob)+Nlen
          Iq=begob
         ELSE IF(Iq.gt.1.and.Iyr.gt.Begyr)THEN
          Nlen=(Iq-1)+Nlen
          Iq=1
         END IF
        END IF
        Len=Nlen+(Ncol-1)*4
       END IF
       lm=mod(Len,4)
       IF(lm.eq.0)lm=4
       IF(ssyear.le.Iyr)THEN
        Icq=Iq
        Icyr=Iyr
        Ic=Iq
       ELSE
        Icq=1
        Icyr=ssyear
        Ic=(Icyr-Iyr)*4+Icq
       END IF
      END IF
   60 IF(Kswt.le.0)THEN
       IF(Nr.eq.10)WRITE(Nerr,1120)ifile
 1120  FORMAT('    EXECUTION FOR THIS X-11Q.2 RUN CEASES.',/,4X,
     &        'WILL ATTEMPT TO EXECUTE REMAINDER OF COMMANDS IN ',A80)
       IF(Nr.eq.5)PRINT *,'ERRORS FOUND IN NAMELIST INPUT . . .',
     &                  ' TRY AGAIN'
       GO TO 30
      END IF
      irerr=0
      IF(Numtbl.gt.0)CALL tabset(nf4,ifmt,itblon,reset,prompt,*200)
      IF(Nr.eq.5)PRINT *,' '
      CALL inc(2,tfile,Title,ifile,crerr,*30,*200)
      CALL inc(3,dfile,dtemp,ifile,crerr,*30,*200)
      CALL inc(4,ofile,otemp,ifile,crerr,*30,*200)
      IF(tfile.ne.'@')Title=tfile
   70 IF(Nw.gt.0.and.ofile.eq.'@')THEN
       ofile=otemp
       IF(lerout)GO TO 190
      ELSE IF(otemp.ne.ofile.or.Nw.eq.0)THEN
       IF(Np.eq.16)CLOSE(16)
       otemp=ofile
       OPEN(16,FILE=ofile,ERR=190)
       IF(lerout)lerout=.false.
       IF(Nw.eq.0)THEN
        Np=16
        Nw=Np
       END IF
      END IF
   80 IF(Nr.eq.10.and.dfile.eq.ifile)THEN
       nd=Nr
       dtemp=dfile
      ELSE IF(nd.gt.0.and.dfile.eq.'@')THEN
       dfile=dtemp
       IF(lerdat)THEN
        dtemp=dfile
        GO TO 180
       ENDIF
      ELSE IF(dtemp.ne.dfile.or.nd.eq.0)THEN
       IF(nd.eq.15)CLOSE(15)
       tfile=dfile
       dtemp=dfile
       OPEN(15,FILE=dfile,STATUS='OLD',ERR=180)
       IF(lerdat)lerdat=.false.
       IF(nd.eq.0)nd=15
      END IF
      IF(rewind.eq.'Y'.or.rewind.eq.'y')REWIND nd
      tfile=dfile
      tser=Ser
      lfq=1
   90 llq=lfq+3
      DO WHILE (.true.)
C --- READ DATA CARDS. CHECK DATE AND IDENT
C --- REGARDLESS OF MONTH IN WHICH SERIES STARTS, XORG(1) IS EQUATED TO
C --- JANUARY, XORG(2) TO FEB , ETC.
       IF(ifmt.lt.1)THEN
        READ(nd,nf4,END=120,ERR=160)(Xorg(i),i=lfq,llq),myr,mser
        IF(y2kset.eq.'Y'.or.y2kset.eq.'y')THEN
         IF(Myr.le.Y2KLIM)THEN
          myr=2000+myr
         ELSE
          myr=1900+myr
         END IF
        ELSE
         myr=1900+myr
        END IF
        IF(Ser.ne.mser)THEN
         irerr=1
         GO TO 110
        ELSE IF(Ser.eq.mser.and.irerr.eq.1)THEN
         irerr=0
        END IF
        IF(kyr.ne.myr)THEN
         irerr=2
         GO TO 110
        ELSE IF(kyr.eq.myr.and.irerr.eq.2)THEN
         irerr=0
        END IF
       ELSE
        IF(llq.gt.Ie)llq=Ie
        IF(begob.gt.lfq)lfq=begob
        READ(nd,nf4,END=170,ERR=160)(Xorg(i),i=lfq,llq)
       END IF
       DO ijk=lfq,llq
        ox(ijk)=Xorg(ijk)
       END DO
       IF(kyr.lt.Lstyr)THEN
        kyr=kyr+1
        lfq=llq+1
        GO TO 90
       ELSE
        IF(Add.eq.0)THEN
C     TEST FOR ZERO, NEGATIVE DATA- MULTIPLICATIVE VERSION
         iz=0
         DO i=begob,Ie
          IF(Xorg(i).le.0.0)CALL zdat(i,Begyr,iz)
         END DO
         IF(iz.eq.1)THEN
          Add=1
          IF(Nr.eq.5)PRINT 1130
          IF(Nr.eq.10)THEN
           WRITE(Nerr,1130)
 1130      FORMAT(10X,'SWITCHING TO ADDITIVE ADJUSTMENT DUE TO ',
     &            'NEGATIVE DATA VALUES')
           IF(Ixerr.eq.0)Ixerr=1
          END IF
         END IF
        END IF
        IF(Kswt.le.0)THEN
C --- SERIES REJECTED
         IF(Nr.eq.5)PRINT 1140
 1140    FORMAT(/48X,'SERIES REJECTED')
         IF(Nr.eq.10)THEN
          WRITE(Nerr,1140)
          WRITE(Nerr,1120)ifile
         END IF
        ELSE
         IF(irs.eq.2)reset='Y'
C --- WRITE TITLE PAGE
         DO i=1,4
          IF(Add.eq.0.and.i.le.3)THEN
           adj(i)=madj(i)
          ELSE IF(Add.eq.1.and.i.le.3)THEN
           adj(i)=aadj(i)
          END IF
          IF(Meas.eq.1)THEN
           run(i)=seas(i)
          ELSE
           run(i)=summ(i)
          END IF
         END DO
         Ib=begob
         Nyear=Begyr
         Myr2=Begyr
         WRITE(Nw,1150)newpg,Title,Ser,quar(begob),Nyear,quar(Lstob),
     &                 Lstyr,(adj(i),i=1,3),(run(i),i=1,4),prnt(Kp),
     &                 (siglim(i),i=1,2)
 1150    FORMAT(A1,42X,'QUARTERLY SEASONAL ADJUSTMENT PROGRAM',/,49X,
     &          'U.S. BUREAU OF THE CENSUS',/,47X,
     &          'STATISTICAL RESEARCH DIVISION',/,53X,
     &          'SEPTEMBER 1, 1987 ',///,10X,'SERIES TITLE- ',A80,10X,
     &          'SERIES NO.  ',A6,/,37X,'PERIOD COVERED- ',A3,
     &          ' QUAR,',I4,' TO ',A3,' QUAR,',I4,/,36X,
     &          'TYPE OF RUN - ',3A6,1X,4A6,/,52X,A6,' PRINTOUT',/,31X,
     &          'SIGMA LIMITS FOR GRADUATING EXTREME VALUES ARE ',F3.1,
     &          ' AND ',F3.1)
         IF(Lswv.ne.0)THEN
          WRITE(Nw,1160)(csma(i),i=1,4)
 1160     FORMAT(15X,'MOVING AVERAGES FOR SEASONAL FACTOR CURVES- 1   2'
     &           ,'   3   4',/,58X,4A4)
         END IF
         PRINT *,' '
         PRINT *,'     X-11Q.2 IS NOW EXECUTING . . . ',Ser
         CALL qsap(issap)
         IF(issap.eq.1)THEN
          issap=2
          L0=Ie-(Nlen+Iq-2+(Ncol-1)*4)
          DO ijk=1,Ncol
           Ib=L0+(ijk-1)*4+Iq-1
           Ie=Ib+Nlen-1
           DO kji=Ib,Ie
            Xorg(kji)=ox(kji)
           END DO
           Kpropt=0
           Kp=1
           CALL qsap(issap)
          END DO
          CALL ssapq
          Ib=mod(Ib,4)
          IF(Ib.eq.0)Ib=4
         END IF
         IF(Numtbl.gt.0)THEN
          itc1=Itchar+1
          ic2=Itchar+2
          f1=' '
          Tblfil(Itchar:ic2)='TBL'
          DO inm=1,Numtbl
           IF(Iovr.lt.3)THEN
            Tblfil(Itchar:ic2)=Tbl(inm)
            IF(Tblfil(itc1:itc1).eq.' ')Tblfil(itc1:itc1)='_'
           END IF
           IF(Itbl(inm).ge.999)THEN
            IF(Iovr.lt.0)THEN
             WRITE(f1,1170)ic2
 1170        FORMAT('(5X,6HTABLE ,A3,20H WAS STORED IN FILE ,A',I2,
     &              ',25H BY THE OVERWRITE OPTION.)')
             IF(ic2.lt.10)f1(41:42)=' A'
            ELSE IF(Iovr.ge.2.and.Itbl(inm).eq.1999)THEN
             WRITE(f1,1180)ic2
 1180        FORMAT('(5X,6HTABLE ,A3,22H WAS APPENDED TO FILE ,A',I2,
     &              ')')
             IF(ic2.lt.10)f1(43:44)=' A'
            ELSE
             WRITE(f1,1190)ic2
 1190        FORMAT('(5X,6HTABLE ,A3,20H WAS STORED IN FILE ,A',I2,')')
             IF(ic2.lt.10)f1(41:42)=' A'
            END IF
            PRINT f1,Tbl(inm),Tblfil(1:ic2)
           END IF
           IF(Itbl(inm).lt.999.and.Itbl(inm).gt.0)
     &        CALL terr(Tbl(inm),reset,irs,0,Itbl(inm),*100,*200)
  100     END DO
         END IF
         PRINT *,'     X-11Q.2 RUN FOR ',Ser,' COMPLETED.'
         PRINT *,' '
         PRINT *,' '
         IF(Nr.eq.10)Irun=Irun+1
        END IF
        GO TO 30
       END IF
  110 END DO
C --  ERRORS IN DATA INPUT AND OUTPUT
  120 inx=index(tfile,' ')-1
      f1=' '
      IF(irerr.eq.1)THEN
       CALL inerr(reset,irs,*200)
       WRITE(f1,1200)inx
 1200  FORMAT('(10X,A',I2,',14H SEARCHED FOR ,A6,18H, SERIES NOT ',
     &        'FOUND)')
       IF(inx.lt.10)f1(6:7)=' A'
       IF(Nr.eq.10)WRITE(Nerr,f1)tfile(1:inx),Ser
       IF(Nr.eq.5)PRINT f1,tfile(1:inx),Ser
      ELSE IF(irerr.eq.2)THEN
       CALL inerr(reset,irs,*200)
       WRITE(f1,1210)inx
 1210  FORMAT('(10X,A',I2,'," SEARCHED FOR ",A6,", DATA FOR "',
     &        ',I4," NOT FOUND")')
       IF(inx.lt.10)f1(6:7)=' A'
       IF(Nr.eq.10)WRITE(Nerr,f1)tfile(1:inx),tser,kyr
       IF(Nr.eq.5)PRINT f1,tfile,tser,kyr
      ELSE IF((kyr-1).lt.Lstyr)THEN
       CALL inerr(reset,irs,*200)
       WRITE(f1,1220)inx
 1220  FORMAT('(10X,A',I2,'," SEARCHED FOR ",A6,", REQUESTED ',
     &        'DATA TO ",I4," NOT FOUND")')
       IF(inx.lt.10)f1(6:7)=' A'
       IF(Nr.eq.10)WRITE(Nerr,f1)tfile(1:inx),tser,kyr
       IF(Nr.eq.5)PRINT f1,tfile,tser,kyr
      END IF
      IF(Nr.eq.10)THEN
       Irun=Irun+1
       WRITE(Nerr,1120)ifile
      END IF
      GO TO 30
  130 CALL inerr(reset,irs,*200)
      WRITE(Nerr,1230)ifile
 1230 FORMAT(10X,'ERROR IN READING NAMELIST NLX11Q FROM ',A80)
      WRITE(Nerr,1120)ifile
      Irun=Irun+1
      GO TO 30
  140 CALL inerr(reset,irs,*200)
      WRITE(Nerr,1240)ifile
 1240 FORMAT(10X,'END OF FILE ENCONTERED IN READING NAMELIST NLX11Q',
     &       ' FROM ',A80)
      IF(Ixerr.eq.0)Ixerr=1
      GO TO 200
  150 PRINT *,' '
      PRINT *,'  ERROR IN NAMELIST NLX11Q ... TRY AGAIN '
      PRINT *,'  CURRENT SETTINGS FOR NLX11Q NAMELIST ARE : '
      PRINT nlx11q
      READ(*,nlx11q,ERR=150)
      IF(stop.ne.'N'.and.stop.ne.'n')GO TO 200
      PRINT nlx11q
      GO TO 50
  160 CALL inerr(reset,irs,*200)
      IF(Nr.eq.10)WRITE(Nerr,1250)tfile
 1250 FORMAT(10X,'ERROR WHILE READING ',A80)
      IF(Nr.eq.5)PRINT 1250,tfile
      IF(Nr.eq.10)THEN
       WRITE(Nerr,1120)ifile
       Irun=Irun+1
      END IF
      GO TO 30
  170 CALL inerr(reset,irs,*200)
      IF(Nr.eq.10)WRITE(Nerr,1260)tfile
 1260 FORMAT(10X,'END OF FILE ENCOUNTERED WHILE READING ',A80)
      IF(Nr.eq.5)PRINT 1260,tfile
      IF(Nr.eq.10)THEN
       WRITE(Nerr,1120)ifile
       Irun=Irun+1
      END IF
      GO TO 30
  180 cre=crerr(icer)
      inx=index(cre,' N')
      inx2=index(tfile,' ')-1
      f1=' '
      WRITE(f1,1270)inx,inx2
 1270 FORMAT('(10X,A',I2,',1X,A',I2,',15H DOES NOT EXIST)')
      IF(inx.lt.10)f1(6:7)=' A'
      IF(inx2.lt.10)f1(13:14)=' A'
      IF(Nr.eq.5.or.icer.eq.5)THEN
       PRINT f1,cre(1:inx),tfile(1:inx2)
       IF(icer.eq.5)THEN
        IF(Nr.ne.5)Nr=5
        GO TO 10
       END IF
       CALL inc(13,dfile,dtemp,ifile,crerr,*80,*200)
      ELSE
       CALL inerr(reset,irs,*200)
       WRITE(Nerr,f1)cre(1:inx),tfile(1:inx2)
       WRITE(Nerr,1120)ifile
       Irun=Irun+1
       GO TO 30
      END IF
  190 inx=index(tfile,' ')-1
      lerout=.true.
      WRITE(f1,1280)inx
 1280 FORMAT('(10X,5HFILE ,A',I2,',28H CANNOT BE OPENED FOR OUTPUT)')
      IF(inx.lt.10)f1(14:15)=' A'
      IF(Nr.eq.5.or.icer.eq.6)THEN
       PRINT f1,tfile(1:inx)
       IF(icer.eq.6)GO TO 20
       CALL inc(14,ofile,otemp,ifile,crerr,*70,*200)
      ELSE
       CALL inerr(reset,irs,*200)
       WRITE(Nerr,f1)tfile(1:inx)
       WRITE(Nerr,1120)ifile
       Irun=Irun+1
       GO TO 30
      END IF
C --- FINISH
  200 WRITE(6,1290)
 1290 FORMAT(1X,'   END OF X-11Q.2 RUN')
      IF(Np.ne.6)THEN
       REWIND Np
       CLOSE(Np)
      END IF
      IF(Nr.eq.10)THEN
       CLOSE(Nr)
       CLOSE(25)
       IF(Ixerr.eq.1)THEN
        inx=index(ifile,' ')-1
        PRINT *,'  ERRORS FOUND IN X-11Q.2 RUN(S) STORED IN ',
     &        ifile(1:inx)
        inx=index(efile,' ')-1
        PRINT *,'  CHECK ',efile(1:inx),' FOR ERROR MESSAGES'
       END IF
      END IF
      IF(Nr.ne.nd.and.nd.ne.0)CLOSE(nd)
      STOP
      END
**==inerr.f    processed by SPAG 4.03F  at 13:27 on  2 Mar 1994
      SUBROUTINE inerr(Reset,Irs,*)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER Irs,Irun,Ixerr,Kpage,Kswerr,Kswt,Lswv,Lyr,Nerr,Np,Nr,Nyear
C*** End of declarations inserted by SPAG
      INTEGER Outdec
      CHARACTER*1 Reset
      CHARACTER*6 Lser
      CHARACTER*80 Title
      COMMON /rite  / Np,Nyear,Title,Kpage,Outdec,Lser
      COMMON /xerr  / Nr,Nerr,Kswerr,Irun,Ixerr
      COMMON /xk    / Kswt,Lswv,Lyr
      IF(Kswt.gt.0)THEN
       Kswt=0
       IF(Nr.eq.10)THEN
        WRITE(Nerr,1010)Irun
 1010   FORMAT(5X,'ERROR MESSAGES FOR X-11Q.2 RUN ',I3,' : ')
        WRITE(Nerr,1020)Lser
        IF(Ixerr.eq.0)Ixerr=1
       END IF
       IF(Nr.eq.5)PRINT 1020,Lser
 1020  FORMAT(5X,' REJECT ',A6)
       IF(Reset.eq.'Y'.or.Irs.eq.9)THEN
        Reset='N'
        Irs=2
       END IF
      END IF
      Kswerr=Kswerr+1
      IF(Nr.eq.10.and.Kswerr.eq.31)THEN
       WRITE(Nerr,1030)
 1030  FORMAT(10X,'ERROR MESSAGES FOR THIS SERIES DISCONTINUED')
       RETURN 1
      END IF
      RETURN
      END
**==zdat.f    processed by SPAG 4.03F  at 13:27 on  2 Mar 1994
      SUBROUTINE zdat(Idat,Ibyr,Iz)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER Ibyr,Idat,Irun,Ixerr,Iz,Kswerr,Nerr,nq,Nr,nyr
C*** End of declarations inserted by SPAG
      COMMON /xerr  / Nr,Nerr,Kswerr,Irun,Ixerr
      Iz=1
      nyr=(Idat-1)/4
      nq=Idat-nyr*4
      nyr=nyr+Ibyr
      IF(Nr.eq.10)WRITE(Nerr,1010)nq,nyr
      IF(Nr.eq.5)PRINT 1010,nq,nyr
 1010 FORMAT(10X,'THE DATA FOR ',I2,'/',I2,' IS ZERO OR NEGATIVE')
      RETURN
      END
**==inc.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE inc(Icer,File,File2,Ifile,Crerr,*,*)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER Icer,Irun,Ixerr,Kswerr,Nerr,Nr
C*** End of declarations inserted by SPAG
      CHARACTER*8 ce
      CHARACTER*25 Crerr(6)
      CHARACTER*80 File,File2,Ifile
      COMMON /xerr  / Nr,Nerr,Kswerr,Irun,Ixerr
      IF(Nr.eq.10)THEN
       READ(Nr,1010,ERR=30,END=40)File
 1010  FORMAT(A80)
      ELSE
       ce='   ENTER'
       IF(Icer.gt.10)THEN
        Icer=Icer-10
        ce='RE-ENTER'
       END IF
       DO WHILE (.true.)
        IF(Icer.eq.2)PRINT *,ce,' TITLE FOR X-11Q.2 RUN :'
        IF(Icer.eq.3)PRINT *,ce,' DATA FILE NAME :'
        IF(Icer.eq.4)PRINT *,ce,' OUTPUT FILE NAME :'
        PRINT *,' > '
        READ 1010,File
        IF(File.eq.'*')THEN
         PRINT *,File2
         GO TO 10
        END IF
        IF(ce.eq.'RE-ENTER')RETURN 1
        GO TO 20
   10  END DO
      END IF
   20 RETURN
   30 WRITE(25,1020)Crerr(Icer),Ifile
 1020 FORMAT(10X,'ERROR IN READING ',A25,' FOUND IN ',A80)
      WRITE(25,1030)
 1030 FORMAT('    EXECUTION FOR THIS X-11Q.2 RUN CEASES.',/,
     &       '    WILL ATTEMPT TO EXECUTE REMAINDER OF INSTRUCTIONS ',
     &       'IN THIS COMMAND FILE.')
      IF(Ixerr.eq.0)Ixerr=1
      RETURN 1
   40 WRITE(25,1040)Crerr(Icer),Ifile
 1040 FORMAT(10X,'END OF FILE ENCONTERED IN READING ',A25,' FOUND IN ',
     &       A80)
      IF(Ixerr.eq.0)Ixerr=1
      RETURN 2
      END
**==prop.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE prop(I,*)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER I
C*** End of declarations inserted by SPAG
      CHARACTER*1 yn
      CHARACTER*6 nl(2)
      DATA nl/'NLX11Q','NLTBLQ'/
      DO WHILE (.true.)
       PRINT *,' '
       PRINT *,' DO YOU WISH TO CHANGE ANY OF THE CURRENT SETTINGS IN ',
     &       nl(I),'? (Y=YES,N=NO)'
       PRINT *,' > '
       READ *,yn
       IF(yn.eq.'Y'.or.yn.eq.'y')RETURN 1
       IF(yn.eq.'N'.or.yn.eq.'n')RETURN
       PRINT *,'IMPROPER INPUT . . . TRY AGAIN'
      END DO
      END
**==writ.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
 
C --- SUBROUTINE TO WRITE OUTPUT
      SUBROUTINE writ(X,Kb,Ke,Ktab,Mtab,Kdec,Lop)
C --- X = ARRAY TO BE PRINTED FROM KB TO KE
C --- KTAB = LOCATION OF WRITE STATEMENT FOR CORRECT TITLE FORMAT
C --- MTAB = TABLE NUMBER
C --- KDEC = NUMBER OF DECIMAL PLACES, 1 = NONE, 2 = ONE PLACE
C --- NOP   = 1,TOTAL. = 2,AVG. = 3,S.D., = 4,MOVING 5-YEAR S.D.
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION c,count,part,parte,partf,Sd5,sdev,Sigml,Sigmu,
     &                 space,subtot,sum,table,total,X,xbar,Xorg,xx,Y,
     &                 yravg
      DOUBLE PRECISION yrsd,Z
      INTEGER i,Ib,ibeg2,Ie,Ifor,igo,Igrz,Io12,Iovr,ist,Itbl,Itchar,
     &        Itfmt,j,j1,j2,jb,je,k,Kb
      INTEGER Kdec,Ke,Kp,Kpage,Kpart,Kpropt,Kswt,Ktab,kyear,l,Lop,Lswv,
     &        Lyr,m,Mcd,Meas,mmod,Mtab,Myr,n
      INTEGER Ndec,newpg,no,nop,Nr,nt,Numtbl,Nw,Nyear
C*** End of declarations inserted by SPAG
      INTEGER Add
      INTEGER Ktcopt,Lterm,Nterm
      DIMENSION X(125),part(3),mmod(9),Xorg(125),xx(3)
      CHARACTER*6 tavg(4),siname(4),Sfc,Bpf,Serno
      CHARACTER*49 nform
      CHARACTER*32 mform
      CHARACTER*80 title2,Title
      COMMON /cfilt / Ktcopt,Lterm(4),Nterm
      COMMON /inpt  / Xorg,Ib,Ie,Meas,Kp,Kpropt
      COMMON /itcom / Itchar,Itfmt,Itbl(10),Ifor,Numtbl,Iovr,Io12
      COMMON /more  / Sfc,Bpf
      COMMON /mor2  / Kpart,Mcd
      COMMON /ngrz  / Igrz,Myr,Nw
      COMMON /rite  / Nr,Nyear,Title,Kpage,Ndec,Serno
      COMMON /sigm  / Sigml,Sigmu,Sd5(125)
      COMMON /tadj  / Add
      COMMON /xk    / Kswt,Lswv,Lyr
      COMMON /xybk  / Y(30),Z(30)
      DATA newpg,space/'1','0'/
      DATA(tavg(i),i=1,4)/'TOTAL ','   AVG','  S.D.','  S.D.'/
      DATA(siname(i),i=1,4)/'RATIOS','DIFFER','      ','ENCES '/
      DATA mmod(3),mmod(4),mmod(8),mmod(9)/'UN',' ','UN',' '/
      DATA(part(i),i=1,3),parte,partf/'B','C','D','E','F'/
      nop=1
      IF(Lop.eq.1)nop=2
      IF(Lop.eq.2)nop=4
      IF(Lop.eq.3)nop=3
      xbar=100.0
      IF(Add.gt.0)xbar=0.0
      k=Kpart
      table=Mtab
      IF(Ktab.ne.11)THEN
C --- WRITE SERIES TITLE
       Kpage=Kpage+1
       WRITE(Nr,1010)newpg,Title,Kpage,Serno
 1010  FORMAT(A1,15X,A80,2X,'PAGE',I3,', SERIES ',A6)
C --- WRITE TABLE TITLE
       GO TO(10,10,10,120,190),Kpart
      ELSE
       WRITE(Nr,1020)
 1020  FORMAT(' D10A. SEASONAL FACTORS ONE YEAR AHEAD')
       GO TO 220
      END IF
C --- PARTS B, C, AND D.
   10 GO TO(20,30,40,50,60,70,80,90,100,110),Ktab
   20 WRITE(Nr,1030)part(k),table
 1030 FORMAT(1X,A1,F3.0,' ORIGINAL SERIES')
      GO TO 220
   30 WRITE(Nr,1040)part(k),table
 1040 FORMAT(1X,A1,F3.0,' TREND-CYCLE, CENTERED 4-TERM MOVING AVERAGE')
      GO TO 220
   40 WRITE(Nr,1050)part(k),table,Bpf,mmod(Mtab),siname(Add+1),
     &              siname(Add+3)
 1050 FORMAT(1X,A1,F3.0,A6,1X,A2,'MODIFIED SI ',2A6)
      GO TO 220
   50 WRITE(Nr,1060)part(k),table,Bpf,siname(Add+1),siname(Add+3)
 1060 FORMAT(1X,A1,F3.0,A6,' REPLACEMENT VALUES FOR EXTREME SI ',2A6)
      GO TO 220
   60 WRITE(Nr,1070)part(k),table,Bpf,Sfc
 1070 FORMAT(1X,A1,F3.0,A6,' SEASONAL FACTORS, ',A6,' MOVING AVERAGE')
      GO TO 220
   70 WRITE(Nr,1080)part(k),table,Bpf
 1080 FORMAT(1X,A1,F3.0,A6,' SEASONALLY ADJUSTED SERIES')
      GO TO 220
   80 WRITE(Nr,1090)part(k),table,Bpf,Nterm
 1090 FORMAT(1X,A1,F3.0,A6,' TREND-CYCLE, ',i1,'-TERM HENDERSON CURVE')
      GO TO 220
   90 WRITE(Nr,1100)part(k),table,Bpf
 1100 FORMAT(1X,A1,F3.0,A6,' IRREGULAR SERIES')
      GO TO 220
  100 WRITE(Nr,1110)part(k),table,Bpf,Sigml,Sigmu
 1110 FORMAT(1X,A1,F3.0,A6,' WEIGHTS FOR IRREGULAR SERIES',/,16X,
     &       'GRADUATION RANGE FROM ',F4.1,' TO ',F4.1,' SIGMA')
      GO TO 220
  110 WRITE(Nr,1120)part(k),table,Bpf
 1120 FORMAT(1X,A1,F3.0,' ORIGINAL SERIES MODIFIED BY ',A6,' WEIGHTS')
      GO TO 220
C --- PART E
  120 GO TO(130,140,150,160,170,180),Ktab
  130 WRITE(Nr,1130)parte,table
 1130 FORMAT(1X,A1,F3.0,' ORIGINAL SERIES MODIFIED FOR EXTREMES')
      GO TO 220
  140 WRITE(Nr,1140)parte,table
 1140 FORMAT(1X,A1,F3.0,' MODIFIED SEASONALLY ADJUSTED SERIES')
      GO TO 220
  150 WRITE(Nr,1150)parte,table
 1150 FORMAT(1X,A1,F3.0,' MODIFIED IRREGULAR SERIES')
      GO TO 220
  160 WRITE(Nr,1160)parte,table,siname(Add+1),siname(Add+3)
 1160 FORMAT(1X,A1,F3.0,2A6,
     &       'OF ANNUAL TOTALS, ORIGINAL AND SEASONALLY ADJUSTED SERIES'
     &       /23X,'UNMODIFIED',26X,'MODIFIED')
      RETURN
  170 WRITE(Nr,1170)parte,table
 1170 FORMAT(1X,A1,F3.0,' QUARTER-TO-QUARTER CHANGES IN ORIGINAL SERIES'
     &       )
      GO TO 220
  180 WRITE(Nr,1180)parte,table
 1180 FORMAT(1X,A1,F3.0,
     &' QUARTER-TO-QUARTER CHANGES IN FINAL SEASONALLY ADJUSTED SERIES (
     &D11.)')
      GO TO 220
C --- PART F
  190 GO TO(200,210),Ktab
  200 WRITE(Nr,1190)partf,table,Mcd
 1190 FORMAT(1X,A1,F3.0,' QCD MOVING AVERAGE'/11X,'QCD IS',I2)
      IF(Mcd.gt.1)GO TO 220
      WRITE(Nr,1200)
 1200 FORMAT(6X,' SEE TABLE D11.')
      RETURN
  210 WRITE(Nr,1210)partf,table
 1210 FORMAT(1X,A1,F3.0,' SUMMARY MEASURES')
      RETURN
  220 WRITE(Nr,1220)space,tavg(nop)
 1220 FORMAT(A1,3X,'YEAR',10X,'1ST QUAR',16X,'2ND QUAR',16X,'3RD QUAR',
     &       16X,'4TH QUAR',16X,A6)
      jb=mod(Kb,4)
      IF(jb.eq.0)jb=4
      jb=jb-1
      je=mod(Ke,4)
      IF(je.gt.0)je=4-je
      IF(jb.gt.0.or.je.gt.0)THEN
       DO i=1,3
        xx(i)=9999999999.0
       END DO
      END IF
      j1=(Kb+3)/4*4-3
      j2=j1+3
      no=(Ke-1)/4*4
      total=0.0
      nt=0
      DO j=j1,Ke,4
       IF(j.le.j1)THEN
        m=Kb
        n=j2
       ELSE IF(j.lt.no)THEN
        m=j
        n=j+3
       ELSE
        m=j
        n=Ke
       END IF
       nt=nt+1
       subtot=0.0
       DO i=m,n
        subtot=subtot+X(i)
       END DO
       total=total+subtot
       count=n-m+1
       yravg=subtot/count
       GO TO(230,240,250,260),nop
  230  Y(nt)=subtot
       GO TO 270
  240  Y(nt)=yravg
       GO TO 270
  250  yrsd=0.0
       DO i=m,n
        yrsd=yrsd+(X(i)-xbar)**2
       END DO
       Y(nt)=sqrt(yrsd/count)
       GO TO 270
  260  Y(nt)=Sd5(n)
  270 END DO
      kyear=Nyear+j1/4
      n=1
C --- WRITE LINE(YEAR) OF DATA
      j=j1
      IF(Mtab.ne.17)THEN
       WRITE(nform,1230)Kdec,Kdec,Kdec
 1230  FORMAT('(A1,3X,I4, 8X,F10.',I1,',3(14X,F10.',I1,'), 12X',
     &        ',F10.',I1,')')
      ELSE
C --- LINE FOR TABLE OF WEIGHTS FOR IRREGULAR
       WRITE(nform,1240)Kdec,Kdec,Kdec
 1240  FORMAT('(A1,3X,I4,10X,2PF8.',I1,',3(16X,F8.',I1,'),14X',
     &        ',0PF8.',I1,')')
      END IF
      WRITE(mform,1250)Kdec,Kdec
 1250 FORMAT('(2X,A6,   8X,F10.',I1,',3(14X,F10.',I1,'))')
      DO WHILE (.true.)
       k=j+3
       IF(j.eq.j1.and.jb.gt.0)THEN
        WRITE(Nr,nform)space,kyear,(xx(l),l=1,jb),(X(i),i=j+jb,k),Y(n)
       ELSE IF(k.ge.Ke.and.je.gt.0)THEN
        WRITE(Nr,nform)space,kyear,(X(i),i=j,k-je),(xx(l),l=1,je),Y(n)
       ELSE
        WRITE(Nr,nform)space,kyear,(X(i),i=j,k),Y(n)
       END IF
       j=j+4
       kyear=kyear+1
       n=n+1
       IF(k.ge.Ke)THEN
        WRITE(Nr,1260)
 1260   FORMAT(' ')
C --- COMPUTE AVERAGE OR STND DEV FOR EACH QUARTER
C
C     TABLES FOR WHICH ABOVE IS NOT PRINTED
        IF(Kpart.lt.3)THEN
        ELSE IF(Kpart.eq.3)THEN
         IF(Mtab.eq.9)GO TO 420
         IF(Ktab.eq.11)GO TO 420
        ELSE
         GO TO 280
        END IF
        IF(Ktab.eq.5)GO TO 360
        GO TO 280
       END IF
      END DO
  280 GO TO(290,290,290,420),nop
  290 DO n=1,4
       ist=n+4*((Kb+3-n)/4)
       sum=0.0
       c=0.0
       DO i=ist,Ke,4
        c=c+1.0
        GO TO(300,300,310),nop
  300   sum=sum+X(i)
        GO TO 320
  310   sum=sum+(X(i)-xbar)**2
  320  END DO
       GO TO(330,330,340),nop
  330  Y(n)=sum/c
       GO TO 350
  340  Y(n)=sqrt(sum/c)
  350 END DO
      m=(nop+1)/2+1
      WRITE(Nr,mform)tavg(m),(Y(i),i=1,4)
C
C     DETERMINE WHICH TABLES ALSO HAVE TABLE AVG AND STND DEV
  360 IF(Mtab.eq.13)GO TO 410
      IF(Mtab.eq.12)GO TO 410
      GO TO(370,400,380,390,400),Kpart
  370 IF(Mtab.eq.1)GO TO 410
      GO TO 400
  380 IF(Mtab.eq.11)GO TO 410
      GO TO 400
  390 IF(Mtab.le.3)GO TO 410
C --- TABLE TOTAL ONLY
  400 WRITE(Nr,1270)total
 1270 FORMAT(20X,'TABLE TOTAL ',F12.3)
      GO TO 420
C --- TOTAL,AVERAGE, AND STANDARD DEVIATION OF TABLE
  410 c=Ke-Kb+1
      xbar=total/c
      sdev=0.0
      DO i=Kb,Ke
       sdev=sdev+(X(i)-xbar)**2
      END DO
      sdev=sqrt(sdev/c)
      WRITE(Nr,1280)total,xbar,sdev
 1280 FORMAT(20X,'TABLE TOTAL ',F12.3,3X,'MEAN ',F12.3,3X,
     &       'STD. DEVIATION ',F12.3)
  420 igo=0
      IF(Mtab.eq.10.and.Kpart.eq.3.and.Ifor.eq.1)igo=1
      IF(Numtbl.gt.0.and.igo.eq.0)THEN
       IF(Ktab.eq.11)Kb=Ib
       CALL trite(X,Kb,Ke,Lyr,Serno,Kpart,Mtab)
      END IF
      IF(Igrz.eq.0)RETURN
      IF(Kpart.eq.1.and.Mtab.eq.1)THEN
       title2='B1. ORIGINAL SERIES'
       CALL grzlst(Kb,ibeg2,Ke,Myr,X,X,108)
       CALL chrt(title2,6,1,ibeg2,Ke)
      ELSE IF(Kpart.eq.3.and.Mtab.eq.11)THEN
       title2='D11. FINAL SEASONALLY ADJUSTED SERIES'
       CALL grzlst(Kb,ibeg2,Ke,Myr,X,X,108)
       CALL chrt(title2,6,1,ibeg2,Ke)
      END IF
      IF(Igrz.eq.1)RETURN
      IF(Kpart.eq.3.and.Mtab.eq.13)THEN
       title2='D13. FINAL IRREGULAR SERIES'
       CALL grzlst(Kb,ibeg2,Ke,Myr,X,X,108)
       CALL chrt(title2,29,1,ibeg2,Ke)
      END IF
      RETURN
      END
**==xtrmes.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C --- SUBROUTINE TO FIND EXTEME VALUES IN ARRAY XI
      SUBROUTINE xtrmes(Xi,Kb,Ke,Xwt)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION c,Sd5,sdev,Sigml,Sigmu,temp,xbar,Xi,Xwt
      INTEGER i,iter,j,k,Kb,Ke,m,mb,me
C*** End of declarations inserted by SPAG
      INTEGER Add
      DIMENSION Xwt(125),Xi(125)
      COMMON /sigm  / Sigml,Sigmu,Sd5(125)
      COMMON /tadj  / Add
      iter=1
      xbar=100.0
      IF(Add.gt.0)xbar=0.0
   10 mb=(Kb+2)/4*4+1
      me=Ke/4*4-19
      IF(me.le.mb)me=mb
c      IF(me.le.mb)me=1
      DO i=mb,me,4
       sdev=0.0
       IF(me.le.mb)THEN
        j=Kb
        k=Ke
       ELSE
        j=i
        k=i+19
       END IF
       c=0.0
       DO m=j,k
        GO TO(30,20),iter
   20   IF(Xwt(m).le.0)GO TO 40
   30   c=c+1.0
        sdev=sdev+(Xi(m)-xbar)**2
   40  END DO
       sdev=sqrt(sdev/c)
       IF(me.gt.mb)THEN
        IF(i.eq.mb)THEN
         j=Kb
         k=mb+11
        ELSE IF(i.ne.me)THEN
         j=i+8
         k=i+11
        ELSE
         j=me+8
         k=Ke
        END IF
       END IF
       DO m=j,k
        Sd5(m)=sdev
        IF(sdev.gt.0)THEN
         temp=(abs(Xi(m)-xbar))/sdev
         IF(temp.gt.Sigmu)THEN
          Xwt(m)=0.0
          GO TO 50
         ELSE IF(temp.gt.Sigml)THEN
          Xwt(m)=(Sigmu-temp)/(Sigmu-Sigml)
          GO TO 50
         END IF
        END IF
        Xwt(m)=1.0
   50  END DO
      END DO
      GO TO(60,70),iter
   60 iter=2
      GO TO 10
   70 RETURN
      END
**==ftestq.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C --- F-TEST SUBROUTINE
      SUBROUTINE ftestq(X,Kb,Ke,Ind,Issap,Icol)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION cf,F,fratio,fvalue,P,prob,qms,rms,ssbq,Ssf,Ssmf,
     &                 ssq,ssr,sst,t,t1,temp,Test1,Test2,X
      INTEGER i,Icol,Ind,Iq,Iqfail,Issap,Issqf,ist,itype,Iyr,Kb,kdfbq,
     &        kdfr,kdft,Ke,Kpage,L0,lb,Lmon,Lyr
      INTEGER n,Np,Nsea,Nyear
C*** End of declarations inserted by SPAG
      INTEGER Outdec
      DIMENSION X(125),temp(125)
      CHARACTER*2 test,blank,star1,star2
      CHARACTER*6 Serno
      CHARACTER*80 Title
      COMMON /com4b / Iyr,Iq,Nsea,Lyr,Lmon
      COMMON /lzero / L0
      COMMON /rite  / Np,Nyear,Title,Kpage,Outdec,Serno
      COMMON /ssft  / Ssf(4),Ssmf(4),Issqf(4)
      COMMON /tests / F(2),P(2),Test1,Test2,Iqfail
      DATA blank,star1,star2/'  ','* ','**'/
      IF(Ind.gt.0.and.Issap.eq.2)RETURN
      itype=0
      Icol=0
      IF(Issap.eq.2)Icol=(Kb-L0-Iq+5)/4
      IF(Ind.eq.1)THEN
       lb=Kb+1
       DO i=lb,Ke
        temp(i)=X(i)-X(i-1)
       END DO
      ELSE
       lb=Kb
       DO i=Kb,Ke
        temp(i)=X(i)
       END DO
      END IF
      DO WHILE (.true.)
       t=0.0
       cf=0.0
       ssbq=0.0
       sst=0.0
       DO n=1,4
        t1=0.0
        ssq=0.0
        ist=n
        IF(lb.gt.n)THEN
         DO WHILE (.true.)
          ist=ist+4
          IF(lb.le.ist)GO TO 10
         END DO
        END IF
   10   DO i=ist,Ke,4
         t1=t1+1.0
         ssq=ssq+temp(i)
         sst=sst+temp(i)**2
        END DO
        ssbq=ssbq+ssq*ssq/t1
        cf=cf+ssq
        t=t+t1
       END DO
       cf=cf*cf/t
       kdfbq=3
       kdfr=t-4.0
       kdft=t-1.0
       ssbq=ssbq-cf
       sst=sst-cf
       ssr=sst-ssbq
       qms=ssbq/3.0
       cf=dble(kdfr)
       rms=ssr/cf
       fratio=qms/rms
       prob=fvalue(fratio,3,kdfr)*100.0
       IF(Ind.eq.0)THEN
        F(1)=fratio
        P(1)=prob
        IF(Issap.eq.2)THEN
         Ssf(Icol)=fratio
        ELSE
         Kpage=Kpage+1
         WRITE(Np,1010)Title,Kpage,Serno
         WRITE(Np,1020)ssbq,kdfbq,qms,fratio,ssr,kdfr,rms,sst,kdft
        END IF
        RETURN
       END IF
       test=blank
       IF(prob.lt.1.0)test=star1
       IF(prob.lt.0.1)test=star2
       IF(itype.eq.0)THEN
        WRITE(Np,1030)
        IF(test.eq.blank)WRITE(Np,1050)test,fratio
        IF(test.ne.blank)WRITE(Np,1040)test,fratio
        lb=Ke-11
        IF(lb.le.Kb)RETURN
        itype=1
        GO TO 20
       ELSE IF(test.ne.blank)THEN
        WRITE(Np,1090)test,fratio
        WRITE(Np,1100)
       ELSE
        WRITE(Np,1060)test,fratio
        IF(prob.le.5.0)WRITE(Np,1070)test
        IF(prob.gt.5.0)WRITE(Np,1080)test
        WRITE(Np,1100)
       END IF
       RETURN
 1010  FORMAT('1',15X,A80,2X,'PAGE',I3,', SERIES ',A6,//)
 1020  FORMAT(35X,'STABLE SEASONALITY TEST',/,43X,
     &        'SUM OF   DEGREES OF    MEAN',/,42X,
     &        'SQUARES    FREEDOM     SQUARE        F',/,20X,
     &        'BETWEEN QUARTERS',F13.3,8X,I2,F12.3,F12.3,/,20X,
     &        'RESIDUAL',6X,F15.3,I10,F12.3,/,20X,'TOTAL',9X,F15.3,I10)
 1030  FORMAT(/,42X,'TEST FOR THE PRESENCE OF RESIDUAL SEASONALITY')
 1040  FORMAT(/,18X,A2,'RESIDUAL SEASONALITY PRESENT IN THE ENTIRE ',
     &        'SERIES AT THE 1 PER CENT LEVEL.        F =',F10.2)
 1050  FORMAT(/,18X,A2,'NO EVIDENCE OF RESIDUAL SEASONALITY IN THE ',
     &        'ENTIRE SERIES AT THE 1 PER CENT LEVEL. F =',F10.2)
 1060  FORMAT(/,18X,A2,'NO EVIDENCE OF RESIDUAL SEASONALITY IN THE ',
     &        'LAST 3 YEARS AT THE 1 PER CENT LEVEL.  F =',F10.2)
 1070  FORMAT(/,18X,A2,'RESIDUAL SEASONALITY PRESENT IN THE LAST 3 ',
     &        'YEARS AT THE 5 PER CENT LEVEL.')
 1080  FORMAT(/,18X,A2,'NO EVIDENCE OF RESIDUAL SEASONALITY IN THE ',
     &        'LAST 3 YEARS AT THE 5 PER CENT LEVEL.')
 1090  FORMAT(/,18X,A2,'RESIDUAL SEASONALITY PRESENT IN THE LAST 3 ',
     &        'YEARS AT THE 1 PER CENT LEVEL.         F =',F10.2)
 1100  FORMAT(/,1X,'NOTE: SUDDEN LARGE CHANGES IN THE LEVEL OF THE ',
     &      'SEASONALLY ADJUSTED SERIES WILL INVALIDATE THE RESULTS OF '
     &      ,'THIS TEST FOR THE',/,50X,'LAST THREE YEAR PERIOD.')
   20 END DO
      END
**==mstest.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE mstest(Array,Jfda,Jlda,Nw,Nyr,Issap,Icol)
C --- AN F TEST FOR MOVING SEASONALITY
C
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION Array,array1,F,fvalue,P,Ssf,Ssmf,Test1,Test2
      INTEGER i,i1,Icol,ifmo,Iqfail,Issap,Issqf,j,j1,Jfda,Jlda,k,k1,l,
     &        l1,lmo,m,Muladd,n,n1
      INTEGER ndgfre,nmin1,noyrs,Nw,Nyr
C*** End of declarations inserted by SPAG
      DOUBLE PRECISION suma1,xbar,colss,colmn,rowss,rowmn,fnyr,totss,
     &                 errss,degfre,fnoyrs,rowssn,errssn,c
      COMMON /ssft  / Ssf(4),Ssmf(4),Issqf(4)
      COMMON /tadj  / Muladd
      COMMON /tests / F(2),P(2),Test1,Test2,Iqfail
      DIMENSION Array(125),array1(125)
      c=1.0D0
      ifmo=(Jfda+Nyr-2)/Nyr*Nyr+1
      lmo=Jlda/Nyr*Nyr
      noyrs=(lmo-ifmo)/Nyr+1
      fnoyrs=noyrs
      IF(Muladd.eq.0)THEN
C -- MULTIPLICATIVE MODEL
       c=10000.0D0
       DO j=ifmo,lmo
        array1(j)=abs(Array(j)-100.0)
       END DO
      ELSE
C --- ADDITIVE MODEL
       DO i=ifmo,lmo
        array1(i)=abs(Array(i))
       END DO
      END IF
C --- ANALYSIS OF VARIANCE TEST
      suma1=0.0D0
      DO k=ifmo,lmo
       suma1=suma1+array1(k)
      END DO
      fnyr=Nyr
      xbar=suma1/(fnyr*fnoyrs)
      colss=0.0D0
      DO l=1,Nyr
       colmn=0.0D0
       k1=ifmo+l-1
       DO m=k1,lmo,Nyr
        colmn=colmn+array1(m)
       END DO
       colmn=colmn/fnoyrs
       colss=colss+(colmn-xbar)*(colmn-xbar)
      END DO
      colss=colss*fnoyrs*c
      rowss=0.0D0
      DO n=ifmo,lmo,Nyr
       rowmn=0.0D0
       l1=n+Nyr-1
       DO i1=n,l1
        rowmn=rowmn+array1(i1)
       END DO
       rowmn=rowmn/fnyr
       rowss=rowss+(rowmn-xbar)*(rowmn-xbar)
      END DO
      rowss=rowss*fnyr*c
      totss=0.0D0
      DO j1=ifmo,lmo
       totss=totss+(array1(j1)-xbar)*(array1(j1)-xbar)
      END DO
      errss=totss*c-colss-rowss
      degfre=fnoyrs-1.0D0
      rowssn=rowss/degfre
      errssn=errss/((fnyr-1.0D0)*degfre)
      ndgfre=(Nyr-1)*(noyrs-1)
      F(2)=rowssn/errssn
      n1=noyrs-1
      P(2)=fvalue(F(2),n1,ndgfre)*100.0
      IF(Issap.eq.2)THEN
       Ssmf(Icol)=F(2)
       RETURN
      END IF
      nmin1=noyrs-1
      WRITE(Nw,1010)rowss,nmin1,rowssn,F(2),errss,ndgfre,errssn
 1010 FORMAT(//,20X,'MOVING SEASONALITY TEST',/,43X,'SUM OF',5X,
     &       'DGRS.OF',9X,'MEAN',/,42X,'SQUARES',5X,'FREEDOM',8X,
     &       'SQUARE',7X,'F-VALUE',/,21X,'BETWEEN YEARS',F17.4,I9,F17.6,
     &       F12.3,/,29X,'ERROR',F17.4,I9,F17.6,/)
      RETURN
      END
**==combft.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE combft(Issap,Icol,Nw)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION F,P,Ssf,Ssmf,test,Test1,Test2
      INTEGER Icol,Iqfail,Issap,Issqf,Nw
C*** End of declarations inserted by SPAG
C --- THIS SUBROUTINUE PRODUCES THE COMBINED TEST FOR THE PRESENCE OF
C --- IDENTIFIABLE SEASONALITY.
C
      COMMON /ssft  / Ssf(4),Ssmf(4),Issqf(4)
      COMMON /tests / F(2),P(2),Test1,Test2,Iqfail
      Iqfail=0
      Test1=9.0
      IF(F(1)*9.0.gt.7.0)Test1=7.0/F(1)
      Test2=(3.0*F(2))/F(1)
      IF(Test2.gt.9.0)Test2=9.0
      test=(Test1+Test2)/2.0
      IF(Test1.ge.1.0.or.Test2.ge.1.0)Iqfail=1
      IF(P(1).ge.0.1)Iqfail=2
      IF(P(2).le.5.0.and.test.ge.1.0)Iqfail=2
      IF(Issap.eq.2)THEN
       Issqf(Icol)=Iqfail
       RETURN
      END IF
      WRITE(Nw,1010)
      IF(Iqfail.eq.0)WRITE(Nw,1030)
      IF(Iqfail.eq.1)WRITE(Nw,1020)
      IF(Iqfail.eq.2)WRITE(Nw,1040)
 1010 FORMAT(//,20X,
     &    'COMBINED TEST FOR THE PRESENCE OF IDENTIFIABLE SEASONALITY :'
     &    )
 1020 FORMAT(/,30X,'IDENTIFIABLE SEASONALITY PROBABLY NOT PRESENT')
 1030 FORMAT(/,30X,'IDENTIFIABLE SEASONALITY PRESENT')
 1040 FORMAT(/,30X,'IDENTIFIABLE SEASONALITY NOT PRESENT')
      RETURN
      END
**==fvalue.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      DOUBLE PRECISION FUNCTION fvalue(X,M,N)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,i1,j,j1,k,l,M,N
C*** End of declarations inserted by SPAG
C --- THIS FUNCTION CALCULATES F-DISTRIBUTION PROBABILITY LEVELS FOR
C ---   PR(R.V. WITH F-DISTRIBUTION WITH M AND N DEGREES OF FREEDOM > X)
      DOUBLE PRECISION w,z,p,y,d,zk,X
      IF(X.gt.0.0)THEN
       IF(X.le.90.0)THEN
        IF(X.le.40.0.or.N.le.150)THEN
         l=(M/2)*2-M+2
         k=(N/2)*2-N+2
         w=X*dble(float(M))/dble(float(N))
         z=1.0D0/(1.0D0+w)
         IF(l.ne.1)THEN
          IF(k.ne.1)THEN
           d=z*z
           p=w*z
          ELSE
           p=dsqrt(z)
           d=0.5D0*z*p
           p=1.0D0-p
          END IF
         ELSE IF(k.ne.1)THEN
          p=dsqrt(w*z)
          d=0.5D0*p*z/w
         ELSE
          p=dsqrt(w)
          y=0.31830988618379D0
          d=y*z/p
          p=2.0D0*y*datan(p)
         END IF
         y=2.0D0*w/z
         j1=k+2
         IF(N.ge.j1)THEN
          IF(l.ne.1)THEN
           zk=z**((N-1)/2)
           d=d*zk*dble(float(N))/dble(float(k))
           p=p*zk+w*z*(zk-1.0D0)/(z-1.0D0)
          ELSE
           DO j=j1,N,2
            d=(1.0D0+dble(float(l))/dble(float(j-2)))*d*z
            p=p+d*y/dble(float(j-1))
           END DO
          END IF
         END IF
         y=w*z
         i1=l+2
         IF(M.ge.i1)THEN
          z=2.0D0/z
          k=N-2
          DO i=i1,M,2
           zk=dble(float(i+k))
           d=y*d*zk/dble(float(i-2))
           p=p-z*d/zk
          END DO
         END IF
         IF(p.lt.1.0D0)THEN
          IF(p.le.0.0D0)GO TO 10
          GO TO 20
         END IF
        END IF
       END IF
       fvalue=0.0
       RETURN
      END IF
   10 fvalue=1.0
      X=0.0
      RETURN
   20 fvalue=sngl(1.0D0-p)
      RETURN
      END
**==replce.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE replce(X,Y,J,K,Xwt)
C --- THIS SUBROUTINE REPLACES VALUES IN ARRAY -X- WHICH HAVE A WEIGHT
C     LESS THAN 1.0
C --- THE REPLACEMENT VALUES ARE STORED IN ARRAY -Y-
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION ct,r,sum,X,Xwt,Y
      INTEGER i,inc,J,K,k2,kno,l,m,m2,mm,n,no
C*** End of declarations inserted by SPAG
      DIMENSION X(125),Y(125),Xwt(125)
      DO i=J,K
       r=X(i)
       IF(Xwt(i).ge.1.0)GO TO 90
       sum=X(i)*Xwt(i)
       kno=0
       IF(i-1.le.J)THEN
C --- EXTREME VALUES IN FIRST 2 LOCATIONS AT EITHER END OF ARRAY ARE
C      REPLACED BY THE THREE NEAREST NONEXTREME VALUES
        m=J
        n=K
        inc=1
       ELSE IF(K-i.le.1)THEN
        m=K
        n=J
        inc=-1
       ELSE
C --- EXTREME CENTER VALUES ARE REPLACED BY TWO NEAREST NONEXTREME
C      VALUES ON EACH SIDE
        k2=2
        mm=i
        GO TO 20
       END IF
   10  IF(Xwt(m).ge.1.0)THEN
        sum=sum+X(m)
        kno=kno+1
        IF(kno.ge.3)GO TO 50
       END IF
       IF(m.eq.n)GO TO 50
       m=m+inc
       GO TO 10
C
   20  inc=1
       n=K
C
   30  no=0
       m=mm
C
       DO WHILE (m.ne.n)
        m=m+inc
        IF(Xwt(m).ge.1.0)THEN
C
         sum=sum+X(m)
         no=no+1
         IF(no.ge.k2)GO TO 40
        END IF
       END DO
C
   40  kno=kno+no
       IF(kno.lt.4)THEN
        IF(mm.eq.i)THEN
C
         IF(inc.le.0)THEN
C
          k2=4-kno
          mm=m2
          GO TO 20
         ELSE
          m2=m
          n=J
          k2=4-no
          inc=-1
          GO TO 30
         END IF
        END IF
       END IF
C
       IF(kno.ge.4)GO TO 70
       GO TO 60
   50  IF(kno.ge.3)GO TO 70
   60  sum=0.0
       DO l=J,K
        sum=sum+X(l)
       END DO
       ct=K-J+1
       GO TO 80
   70  ct=kno
       ct=ct+Xwt(i)
   80  r=sum/ct
   90  Y(i)=r
      END DO
      RETURN
      END
**==avec.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C --- SUBROUTINE TO CALCULATE AVERAGE QUARTER-TO-QUARTER CHANGE
      SUBROUTINE avec(X,J,K,Xbar,Xbar2,Sdc)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION count,Sdc,X,Xbar,Xbar2
      INTEGER i,J,K,km,m,n
C*** End of declarations inserted by SPAG
C --- INPUT IS ARRAY X - LOCATION OF 1ST VALUE IS  J, OF LAST IS K.
C --- OUTPUT
C           SDC  = STD. DEV. OF CHANGES WITH REGARD TO SIGN
C           XBAR = AVERAGE CHANGE WITHOUT REGARD TO SIGN
C           XBAR2 = AVG CHANGE WITH REGARD TO SIGN
C --- IF ADDITIVE ADJUSTMENT, CHANGES ARE QUARTER-TO-QUARTER DIFFERENCE
C     IF MULTIPLICATIVE, QUARTER-TO-QUARTER PERCENT CHANGE
      INTEGER Add
      DIMENSION X(125),Xbar(4),Xbar2(4),Sdc(4)
      COMMON /tadj  / Add
      DO m=1,4
       Xbar(m)=0.0
       Xbar2(m)=0.0
       Sdc(m)=0.0
       km=K-m
       count=km-J+1
       DO i=J,km
        n=i+m
        IF(Add.le.0)THEN
         Xbar(m)=Xbar(m)+abs(X(n)-X(i))/X(i)*100.0
         Xbar2(m)=Xbar2(m)+(X(n)-X(i))/X(i)*100.0
        ELSE
         Xbar(m)=Xbar(m)+abs(X(n)-X(i))
         Xbar2(m)=Xbar2(m)+(X(n)-X(i))
        END IF
       END DO
       Xbar(m)=Xbar(m)/count
       Xbar2(m)=Xbar2(m)/count
       IF(Add.le.0)THEN
        DO i=J,km
         n=i+m
         Sdc(m)=Sdc(m)+((X(n)-X(i))/X(i)*100.0-Xbar2(m))**2
        END DO
       ELSE
        DO i=J,km
         n=i+m
         Sdc(m)=Sdc(m)+(X(n)-X(i)-Xbar2(m))**2
        END DO
       END IF
       IF(Sdc(m).gt.0)Sdc(m)=sqrt(Sdc(m)/count)
      END DO
      RETURN
      END
**==avedur.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C --- AVERAGE DURATION OF RUN SUBROUTINE
      SUBROUTINE avedur(Y,L,M,Adr)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION Adr,count,runs,Y
      INTEGER i,L,M
C*** End of declarations inserted by SPAG
      DIMENSION Y(125)
      i=L
      runs=1.
      IF(Y(i).lt.Y(i+1))THEN
      ELSE IF(Y(i).eq.Y(i+1))THEN
       DO WHILE (.true.)
        i=i+1
        IF(i.ge.M)GO TO 30
        IF(Y(i).lt.Y(i+1))GO TO 10
        IF(Y(i).ne.Y(i+1))THEN
C***** SPAG has made duplicate copies of the following statement 
         runs=runs+1.0
         GO TO 20
        END IF
       END DO
      ELSE
       GO TO 20
      END IF
   10 DO WHILE (.true.)
       i=i+1
       IF(i.ge.M)GO TO 30
       IF(Y(i).gt.Y(i+1))THEN
C***** The following statement is a duplicate copy made by SPAG
        runs=runs+1.0
        GO TO 20
       END IF
      END DO
   20 DO WHILE (.true.)
       i=i+1
       IF(i.ge.M)GO TO 30
       IF(Y(i).lt.Y(i+1))THEN
        runs=runs+1.0
        GO TO 10
       END IF
      END DO
   30 count=M-L
      Adr=count/runs
      RETURN
      END
**==si.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C --- SUBROUTINE TO COMPUTE SI RATIOS FOR PART B (TABLES 4 AND 9)
      SUBROUTINE si(Xsi,J,K,Xsit,Iter)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,Iter,J,K,l,m,n
      DOUBLE PRECISION wt,X,Xi,Xs,Xsi,Xsit,Xwt,Y
C*** End of declarations inserted by SPAG
      DIMENSION Xsi(125),Xsit(125),wt(30)
      COMMON /extr  / Xs(129),Xi(125),Xwt(125)
      COMMON /xybk  / X(30),Y(30)
      CALL qsf(Xsi,J,K,Iter)
      CALL divsub(Xsi,Xs,Xi,J,K)
      CALL xtrmes(Xi,J,K,Xwt)
      l=J-1
      DO n=1,4
       l=l+1
       m=0
       DO i=l,K,4
        m=m+1
        wt(m)=Xwt(i)
        X(m)=Xsi(i)
       END DO
       CALL replce(X,Y,1,m,wt)
       m=0
       DO i=l,K,4
        m=m+1
        Xsit(i)=Y(m)
       END DO
      END DO
      DO i=J,K
       IF(Xsit(i).ne.Xsi(i))THEN
        Xsi(i)=Xsit(i)
       ELSE
        Xsit(i)=9999999999.
       END IF
      END DO
      RETURN
      END
**==divsub.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C --- SUBROUTINE TO DIVIDE/SUBTRACT TWO SERIES
      SUBROUTINE divsub(X,Y,Z,J,K)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,J,K
      DOUBLE PRECISION X,Y,Z
C*** End of declarations inserted by SPAG
      INTEGER Add
      DIMENSION X(125),Y(125),Z(125)
      COMMON /tadj  / Add
      DO i=J,K
       IF(Add.le.0)THEN
        Z(i)=(X(i)/Y(i))*100.
       ELSE
        Z(i)=X(i)-Y(i)
       END IF
      END DO
      RETURN
      END
**==chrt.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE chrt(Title,Icode,Noser,Ibeg2,Ienda)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Add,Xdata3,Xy,Y1,Y2,Ydiff,Ymax,Ymin
      INTEGER ialphc,Ibeg2,icod2,Icode,Ienda,Ifrst,Ixy,jx,jxo,jyo,l,
     &        Last,Llyr,Logsw,Lstyr,Mach,N1,Noser,Npts,Nseas
      INTEGER Nyr
C*** End of declarations inserted by SPAG
      CHARACTER*1 Ia(110,55),Iq(4),I1,I2,I3,I4,I5,I6,I7,I8,I9,Ip(10),
     &            I10,I11,I12
      CHARACTER*80 Title
      COMMON /arim1 / Nseas,Xy,Ixy
      COMMON /chrtm / Mach
      COMMON /chrtp / Ia,Iq,Ip
      COMMON /chrtq / I1,I2,I3,I4,I5,I6,I7,I8,I9
      COMMON /chrtq3/ I10,I11,I12
      COMMON /grzar1/ Y1(125),Y2(125),Xdata3(31),Ymin,Ymax,Ydiff
      COMMON /grzar2/ Nyr,Npts,N1,Ifrst,Last,Llyr,Lstyr
      COMMON /lg1   / Logsw
      COMMON /tadj  / Add
      IF(Mach.eq.0)THEN
       I3=char(205)
       I4=char(186)
       I8=char(187)
       I9=char(188)
       I10=char(202)
       I11=char(203)
       I12=char(223)
      ELSE
       I3='+'
       I4='I'
       I8=I3
       I9=I3
       I10=I4
       I11=I4
       I12='@'
      END IF
      DO l=1,10
       Ip(l)=I3
      END DO
      IF(Icode.eq.0)RETURN
      Logsw=0
      Nseas=4
      IF(Ymax.le.Ymin)RETURN
      icod2=0
      IF(Icode.eq.15)THEN
       icod2=Icode
       Icode=5
      END IF
      IF(Icode.eq.29)THEN
       IF(Add.eq.0)THEN
        icod2=Icode
       ELSE
        icod2=19
       END IF
       Icode=9
      END IF
      CALL setup(Icode,icod2,jxo,jyo,ialphc)
      jx=2
      IF(Icode.eq.9)CALL aver(Y1,N1,'I',Icode,icod2,jx)
      IF(Icode.eq.7)CALL month(Icode,Ibeg2,Ienda,jx)
      IF(Icode.le.6)CALL yrly(Icode,icod2,jx,Noser)
      CALL out(Title,Icode)
      RETURN
      END
**==setup.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C**************
      SUBROUTINE setup(Icode,Icod2,Jxo,Jyo,Ialphc)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Fact1,Xdata3,Xy,Y1,Y2,Ydiff,Ydifl,Ymax,Ymaxl,Ymid,Ymin,Yminl
      INTEGER i,Ialphc,Ibottm,Icmax,Icod2,Icode,Ifact2,Ifrst,Imid,Inyr,
     &        Ixy,Iyr,j,j2,Jxo,Jyo,k,k2,Last,Llyr
      INTEGER Logsw,Lstyr,N1,Npts,Nseas,Nyr
C*** End of declarations inserted by SPAG
      CHARACTER*1 Ia(110,55),I1,I2,I3,I4,I5,I6,I7,I8,I9,Iq(4)
      CHARACTER*1 Ip(10),I10,I11,I12
      COMMON /arim1 / Nseas,Xy,Ixy
      COMMON /chrtp / Ia,Iq,Ip
      COMMON /chrtp2/ Iyr(30),Inyr
      COMMON /chrtq / I1,I2,I3,I4,I5,I6,I7,I8,I9
      COMMON /chrtq2/ Icmax
      COMMON /chrtq3/ I10,I11,I12
      COMMON /chrtr / Yminl,Ymaxl,Ydifl,Ymid(14)
      COMMON /chrts / Ibottm,Imid(14),Fact1,Ifact2
      COMMON /grzar1/ Y1(125),Y2(125),Xdata3(31),Ymin,Ymax,Ydiff
      COMMON /grzar2/ Nyr,Npts,N1,Ifrst,Last,Llyr,Lstyr
      COMMON /lg1   / Logsw
C  SET UP CHART BOUNDARIES
C***************
      DO i=1,110
       DO j=1,55
        Ia(i,j)=' '
       END DO
      END DO
      Inyr=0
      DO i=1,30,2
       IF(Ifrst.eq.1)THEN
        j=((i-1)/2)+1
        Iyr(j)=int(Xdata3(i))
       ELSE
        j=(i+1)/2
        Iyr(j)=int(Xdata3(i+1))
       END IF
       IF(j.le.Nyr)Inyr=Inyr+1
      END DO
      Icmax=N1+1
      IF(Icode.eq.19)Icmax=Icmax+1
      IF(Icode.eq.11)Icmax=1+(Nseas+1)*4
      Jxo=1
      Jyo=1
      Ialphc=Ifrst-1
      Ydiff=Ymax-Ymin
      IF(Logsw.eq.1)THEN
       Ymaxl=alog10(Ymax)
       Yminl=alog10(Ymin)
       Ydifl=Ymaxl-Yminl
      END IF
      DO j=1,53,4
       j2=((53-j)/4)+1
       IF(Logsw.eq.0)THEN
        Ymid(j2)=(Ydiff*(j-1))/52+Ymin
       ELSE
        Ymid(j2)=10**((Ydifl*(j-1))/52+Yminl)
       END IF
      END DO
      Ia(1,1)=I11
      Ia(1,Ibottm)=I10
      Ia(Icmax+1,1)=I8
      Ia(Icmax+1,Ibottm)=I9
      DO i=2,Ibottm-1
       Ia(1,i)=I4
       Ia(Icmax+1,i)=I4
      END DO
      IF((Icode-7)*(Icode-9)*(Icod2-15).ne.0)THEN
       k=Nseas-Ifrst+2
       DO i=1,Ibottm
        DO j=k,Icmax,Nseas
         Ia(j,i)=I7
        END DO
       END DO
      END IF
      DO i=2,Icmax
       Ia(i,1)=I3
       Ia(i,Ibottm)=I3
       IF((Icod2-15)*(Icode-7)*(Icode-9).ne.0)THEN
        DO k2=3,12,3
         Ia(i,Imid(k2))=I7
        END DO
       END IF
      END DO
      RETURN
      END
**==yrly.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C*****************
      SUBROUTINE yrly(Icode,Icod2,Jx,Noser)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Fact1,Xdata3,Xy,Y1,Y2,Ydiff,Ymax,Ymin
      INTEGER i,Ibottm,Icod2,Icode,Ifact2,Ifrst,Imid,Ixy,j,Jx,jy1,jy2,l,
     &        Last,Llyr,Lstyr,N1,Noser,Npts,Nseas
      INTEGER Nyr
C*** End of declarations inserted by SPAG
      CHARACTER*1 Ia(110,55),I1,I2,I3,I4,I5,I6,I7,I8,I9,Iq(4)
      CHARACTER*1 Ip(10),I10,I11,I12
      COMMON /arim1 / Nseas,Xy,Ixy
      COMMON /chrtp / Ia,Iq,Ip
      COMMON /chrtq / I1,I2,I3,I4,I5,I6,I7,I8,I9
      COMMON /chrtq3/ I10,I11,I12
      COMMON /chrts / Ibottm,Imid(14),Fact1,Ifact2
      COMMON /grzar1/ Y1(125),Y2(125),Xdata3(31),Ymin,Ymax,Ydiff
      COMMON /grzar2/ Nyr,Npts,N1,Ifrst,Last,Llyr,Lstyr
C  PLOT CHART TYPES 1-6
C******************
      CHARACTER*1 itype
      DO i=1,N1
       Xy=Y1(i)
       CALL value(*10)
       itype=I1
       l=mod(Ifrst+i-2,Nseas)+1
       IF((Icode-5)*(Icode-6).eq.0)itype=Iq(l)
       Ia(Jx,Ixy)=itype
       IF(Icode.eq.6)THEN
        DO j=Ibottm-1,Ixy,-1
         Ia(Jx,j)=itype
        END DO
       END IF
       IF(Noser.eq.2)THEN
        jy1=Ixy
        Xy=Y2(i)
        CALL value(*10)
        jy2=Ixy
        Ia(Jx,jy2)=I12
        IF(Icod2.eq.15)CALL dot(jy1,jy2,Jx)
       END IF
   10  Jx=Jx+1
      END DO
      IF(Icod2.ge.15)Icode=Icod2
      RETURN
      END
**==month.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C***************
      SUBROUTINE month(Icode,Ibeg2,Ienda,Jx)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Aa,Bb,Xy
      INTEGER Ibeg2,Icode,Ienda,Ifrst,Ixy,Jx,l,Last,Llyr,Lstyr,N1,Npts,
     &        Nseas,Nyr,nyr2
C*** End of declarations inserted by SPAG
      CHARACTER*1 Ia(110,55),Iq(4)
      CHARACTER*1 Ip(10)
      COMMON /arim1 / Nseas,Xy,Ixy
      COMMON /chrtp / Ia,Iq,Ip
      COMMON /grzar2/ Nyr,Npts,N1,Ifrst,Last,Llyr,Lstyr
      COMMON /grzg2 / Aa(31),Bb(31)
C  PLOT CHART TYPES 7 AND 8 (MONTHLY CHARTS)
C***************
      CHARACTER*1 itype
      CALL grzmth(Ibeg2,Ienda)
      DO l=1,Nseas
       nyr2=Nyr
       IF((l.lt.Ifrst).or.(l.gt.Last))nyr2=Nyr-1
       CALL grzmyr(l)
       itype=Iq(l)
       CALL aver(Aa,nyr2,itype,Icode,0,Jx)
      END DO
      RETURN
      END
**==out.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C***************
      SUBROUTINE out(Title,Icode)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Fact1,Xdata3,Y1,Y2,Ydiff,Ydifl,Ymax,Ymaxl,Ymid,Ymin,Yminl
      INTEGER i,Ibottm,Icmax,Icode,Ifact2,Ifrst,Igrz,iline,Imid,Inyr,
     &        iyf,iyl,Iyr,j,j2,Last,Llyr,Lstyr,n,N1
      INTEGER Npts,Nw,Nyr,Nyr2
C*** End of declarations inserted by SPAG
      CHARACTER*1 Ia(110,55),Iq(4)
      CHARACTER*18 for
      CHARACTER*1 Ip(10)
      COMMON /chrtp / Ia,Iq,Ip
      COMMON /chrtp2/ Iyr(30),Inyr
      COMMON /chrtr / Yminl,Ymaxl,Ydifl,Ymid(14)
      COMMON /chrtq2/ Icmax
      COMMON /chrts / Ibottm,Imid(14),Fact1,Ifact2
      COMMON /grzar1/ Y1(125),Y2(125),Xdata3(31),Ymin,Ymax,Ydiff
      COMMON /grzar2/ Nyr,Npts,N1,Ifrst,Last,Llyr,Lstyr
      COMMON /ngrz  / Igrz,Nyr2,Nw
C  WRITE OUT THE CHART
C***************
      CHARACTER*80 Title
      WRITE(Nw,1010)Title
 1010 FORMAT(' ',A80)
 1020 FORMAT(1X,120A1)
 1030 FORMAT(1X,F10.2,110A1)
 1040 FORMAT(1X,10X,110A1)
      iline=Icmax+1
      j2=1
      DO j=1,55
       IF(j.eq.1.or.j.eq.55)THEN
        WRITE(Nw,1020)(Ip(i),i=1,10),(Ia(i,j),i=1,iline)
       ELSE IF(j.eq.Imid(j2))THEN
        WRITE(Nw,1030)Ymid(j2),(Ia(i,j),i=1,iline)
        j2=j2+1
        IF(j2.gt.14)j2=1
       ELSE
        WRITE(Nw,1040)(Ia(i,j),i=1,iline)
       END IF
      END DO
      IF(Icode.eq.7)THEN
       iyf=int(Xdata3(1))
       iyl=int(Xdata3(Nyr))
       WRITE(Nw,1050)iyf,iyl
      ELSE
       n=min0(14,Inyr)
       IF(Icmax.le.61)n=min0(8,Inyr)
       for='(12X,I4,12(4X,I4))'
       IF(Ifrst.eq.1)THEN
        for(10:10)='3'
       ELSE IF(Ifrst.eq.2)THEN
        for(3:3)='5'
       ELSE IF(Ifrst.eq.3)THEN
        for(3:3)='4'
       ELSE IF(Ifrst.eq.4)THEN
        for(3:3)='3'
       END IF
       WRITE(Nw,for)(Iyr(i),i=1,n)
      END IF
 1050 FORMAT(T20,I4,'  TO  ',I4)
      RETURN
      END
**==aver.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C************************
      SUBROUTINE aver(Val,N,I,Icode,Icod2,Jx)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER Icod2,Icode,Ifrst,imnth,istep,Ixy,j,Jx,jyave,jyval,k,l,
     &        Last,Llyr,Lstyr,N,N1,Npts,Nseas,Nyr
      REAL sum,Val,Xy
C*** End of declarations inserted by SPAG
      CHARACTER*1 Ia(110,55),I1,I2,I3,I4,I5,I6,I7,I8,I9,Iq(4)
      CHARACTER*1 Ip(10)
      COMMON /arim1 / Nseas,Xy,Ixy
      COMMON /chrtp / Ia,Iq,Ip
      COMMON /chrtq / I1,I2,I3,I4,I5,I6,I7,I8,I9
      COMMON /grzar2/ Nyr,Npts,N1,Ifrst,Last,Llyr,Lstyr
C PLOT BY MONTH AROUND MEAN OF EACH MONTH
C************************
      DIMENSION Val(*)
      CHARACTER*1 itype
      CHARACTER*1 I
      imnth=Ifrst
      itype=Iq(imnth)
      IF(Icode.eq.9)I=itype
      IF(Icod2.eq.19)THEN
       Xy=0.
      ELSE IF(Icod2.eq.29)THEN
       Xy=100.
      ELSE
       sum=0
       DO l=1,N
        sum=sum+Val(l)
       END DO
       Xy=sum/float(N)
      END IF
      CALL value(*20)
      jyave=Ixy
      DO k=1,N
       Xy=Val(k)
       CALL value(*10)
       jyval=Ixy
       istep=1
       IF(jyval.lt.jyave)istep=-1
       DO j=jyave,jyval,istep
        Ia(Jx,j)=I
       END DO
       IF((Icode-7)*(Icode-9).eq.0)Ia(Jx,jyave)=I1
       Jx=Jx+1
       IF(Icode.eq.9)THEN
        imnth=imnth+1
        IF(imnth.eq.Nseas+1)imnth=1
        I=Iq(imnth)
       END IF
   10 END DO
   20 RETURN
      END
**==value.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE value(*)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Fact1,Xdata3,Xy,Y1,Y2,Ydiff,Ydifl,yl,Ymax,Ymaxl,Ymid,Ymin,
     &     Yminl,yvalue
      INTEGER Ibottm,Ifact2,Imid,Jy,Logsw,Nseas
C*** End of declarations inserted by SPAG
      COMMON /arim1 / Nseas,Xy,Jy
      COMMON /chrtr / Yminl,Ymaxl,Ydifl,Ymid(14)
      COMMON /chrts / Ibottm,Imid(14),Fact1,Ifact2
      COMMON /grzar1/ Y1(125),Y2(125),Xdata3(31),Ymin,Ymax,Ydiff
      COMMON /lg1   / Logsw
C*************
      IF(Logsw.eq.0)THEN
       IF(Xy.gt.Ymax)RETURN 1
       yvalue=Fact1*(Xy-Ymin)/Ydiff
       Jy=Ifact2-int(yvalue)
      ELSE
       yl=alog10(Xy)
       IF(yl.gt.Ymaxl)RETURN 1
       yvalue=Fact1*(yl-Yminl)/Ydifl
       Jy=Ifact2-int(yvalue)
      END IF
      RETURN
      END
**==dot.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C*************
      SUBROUTINE dot(Jy1,Jy2,Jx)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER j,j1,j2,Jx,Jy1,Jy2
C*** End of declarations inserted by SPAG
      CHARACTER*1 Ia(110,55),I1,I2,I3,I4,I5,I6,I7,I8,I9,Iq(4)
      CHARACTER*1 Ip(10)
      COMMON /chrtp / Ia,Iq,Ip
      COMMON /chrtq / I1,I2,I3,I4,I5,I6,I7,I8,I9
C*************
      j1=min(Jy1,Jy2)+1
      j2=max(Jy1,Jy2)-1
      IF(j2.ge.j1)THEN
       DO j=j1,j2
        Ia(Jx,j)=I7
       END DO
      END IF
      RETURN
      END
**==grzlst.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE grzlst(Ibeg,Ibeg2,Iend,Lyr,A,B,Nmon)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL A1,B1,Grzmax,Grzmin,Xdata3,Ydiff
      INTEGER i,Ibeg,Ibeg2,Iend,Ifrst,k,Kpage,l,Last,Llyr,Lstyr,Lyr,N1,
     &        Ndec,nfreq,Nmon,Npts,Nw,Nyear,Nyr
C*** End of declarations inserted by SPAG
      CHARACTER*80 Title
      CHARACTER*6 Serno
      COMMON /grzar1/ A1(125),B1(125),Xdata3(31),Grzmin,Grzmax,Ydiff
      COMMON /grzar2/ Nyr,Npts,N1,Ifrst,Last,Llyr,Lstyr
      COMMON /rite  / Nw,Nyear,Title,Kpage,Ndec,Serno
      DOUBLE PRECISION A(125),B(125)
      Kpage=Kpage+1
      WRITE(Nw,1010)Title,Kpage,Serno
 1010 FORMAT('1',17X,A80,1X,'PAGE ',I2,', SERIES ',A6)
      Grzmin=999999999.
      Grzmax=-999999999.
      Ibeg2=Iend-Nmon+1
      IF(Ibeg2.lt.Ibeg.or.Nmon.eq.0)Ibeg2=Ibeg
      l=0
      nfreq=4
      Ifrst=mod(Ibeg2-1,nfreq)+1
      Last=mod(Iend-1,nfreq)+1
      DO k=Ibeg2,Iend
       Npts=Npts+1
       l=l+1
       A1(l)=dble(A(k))
       B1(l)=dble(B(k))
       Grzmin=amin1(Grzmin,A1(l),B1(l))
       Grzmax=amax1(Grzmax,A1(l),B1(l))
      END DO
      Llyr=Lyr+Ibeg2/nfreq
      Lstyr=Lyr+(Iend+nfreq-1)/nfreq-1
      Nyr=Lstyr-Llyr+1
      N1=Iend-Ibeg2+1
      DO i=1,Nyr
       Xdata3(i)=Llyr+i-1
      END DO
      RETURN
      END
**==grzmth.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE grzmth(Ibeg,Iend)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL A1,B1,Ser1,Ser2,Xdata3,Ydiff,Ymax,Ymin
      INTEGER i,Ibeg,Iend,iyear,j,month,nfreq
C*** End of declarations inserted by SPAG
      COMMON /grzmon/ Ser1(31,12),Ser2(31,12)
      COMMON /grzar1/ A1(125),B1(125),Xdata3(31),Ymin,Ymax,Ydiff
C---- INITIALIZE THE ARRAYS
      DO i=1,31
       DO j=1,12
        Ser1(i,j)=9999999999.
        Ser2(i,j)=9999999999.
       END DO
      END DO
C---PUT THE PROPER DATA VALUES TO BE PLOTTED INTO THE ARRAYS
      i=0
      j=0
      nfreq=4
      DO iyear=1,31
       DO month=1,nfreq
        i=i+1
        IF(i.ge.Ibeg)THEN
         IF(i.gt.Iend)GO TO 10
         j=j+1
         Ser1(iyear,month)=A1(j)
         Ser2(iyear,month)=B1(j)
        END IF
       END DO
      END DO
   10 RETURN
      END
**==grzmyr.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE grzmyr(L)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Aa,Bb,Ser1,Ser2
      INTEGER i,iyear,L
C*** End of declarations inserted by SPAG
      COMMON /grzg2 / Aa(31),Bb(31)
      COMMON /grzmon/ Ser1(31,12),Ser2(31,12)
      i=0
      DO iyear=1,31
       IF(Ser1(iyear,L).le.999999998.)THEN
        i=i+1
        Aa(i)=Ser1(iyear,L)
        Bb(i)=Ser2(iyear,L)
       END IF
      END DO
      RETURN
      END
**==aa0002.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      BLOCK DATA 
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Cut,Fact1
      INTEGER Ibottm,Ifact2,ii,Imid,j,Kount
C*** End of declarations inserted by SPAG
      CHARACTER*1 Ia(110,55),Ialpha(4),I1,I2,I3,I4,I5,I6,I7,I8,I9,Ch,
     &            Ip(10)
      CHARACTER*18 Chper
      CHARACTER*24 Cname
      CHARACTER*50 F1
      CHARACTER*35 F2
      COMMON /chame / Cname(4),Ch(4),Chper(2)
      COMMON /chrtp / Ia,Ialpha,Ip
      COMMON /chrtq / I1,I2,I3,I4,I5,I6,I7,I8,I9
      COMMON /chrts / Ibottm,Imid(14),Fact1,Ifact2
      COMMON /com5  / Kount(4),Cut(4,4)
      COMMON /mform / F1,F2
      DATA Ialpha/'1','2','3','4'/
      DATA I1,I2,I4,I5,I6,I7/'*','-','I','\','/','.'/
      DATA Ibottm,Fact1,Ifact2/55,52.,54/
      DATA Imid/2,6,10,14,18,22,26,30,34,38,42,46,50,54/
      DATA F1/'(1X,I2,A1,I4,2X, (F10.3,5X),3X,2A1,3X,F10.3,5X,A4)'/
      DATA F2/'(10X, (1X,I2,A1,I4,1X,A1,5X),8X,A7)'/
      DATA(Cut(1,ii),ii=1,4)/3.,4.,5.,6./
      DATA(Cut(2,ii),ii=1,4)/3.,4.,5.,6./
      DATA(Cut(3,ii),ii=1,4)/3.,5.,7.,10./
      DATA(Cut(4,ii),ii=1,4)/3.,4.,5.,6./
      DATA(Chper(j),j=1,2)/'QUARTER-TO-QUARTER','   YEAR-TO-YEAR   '/
      DATA(Cname(j),j=1,4)/'SEASONAL FACTORS        ',
     &     'SEASONALLY ADJUSTED DATA','CHANGES IN S. A. DATA   ',
     &     'CHANGES IN S. A. DATA   '/
      DATA(Ch(j),j=1,4)/'%','#','$','@'/
      END
**==qsap.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE qsap(Issap)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION adrc,adrci,Adri,adrqcd,cbar,cbar2,cibar,cibar2,
     &                 et1,et3,F,obar,obar2,P,qbar,qbar2,R,r1,r2,Rati
      DOUBLE PRECISION Ratic,ratio,rato,ric,rii,ris,rmsr,rtot,Rv,sbar,
     &                 sbar2,sdcc,sdcci,sdci,sdco,sdcq,sdcs,Stci,temps,
     &                 tmp,Test1
      DOUBLE PRECISION Test2,ti,torg,tsa,Vc,Vi,Vs,Xi,xmi,xmorg,
     &                 xmtci,Xorg,Xs,xsi,xsit,Xtc,xtci,xtcsi,xtem
      DOUBLE PRECISION xtem2,xtemp,Xwt,xmcd
      INTEGER i,i2,Ib,ib1,ibeg2,icol,Ie,ie1,Igrz,il,il2,Iqfail,ir,Issap,
     &        iyear,j,j3,jtble,kde,Kp,Kpage
      INTEGER Kpart,Kswt,Ktcopt,kyear,l,Length,llda1,Lstob,Lstyr,Lswv,
     &        Lterm,ltest,Mach,Mcd,Meas,modq,Myr,nd2,Ndec,Np,Nr
      INTEGER nt1,nt2,Nterm,Nw,Nyear
C*** End of declarations inserted by SPAG
      INTEGER Add,Begyr
      DOUBLE PRECISION ibar(4),ibar2(4),mobar(4),mcibar(4),mibar(4)
      DIMENSION Xorg(125),xtcsi(125),Xtc(125),xsi(125),xtci(125),
     &          xtemp(125),xmorg(125),xmtci(125),xmi(125),xsit(125),
     &          ratio(4),xtem(4),xtem2(4),obar(4),cibar(4),
     &          cbar(4),sbar(4),qbar(4),obar2(4),cibar2(4),cbar2(4),
     &          sbar2(4),qbar2(4),sdco(4),sdcci(4),sdcc(4),sdcs(4),
     &          sdcq(4),sdci(4),ric(4),rii(4),ris(4),rtot(4),rato(4),
     &          rmsr(8)
      CHARACTER*6 word(6),pre,fin,Bpf,bla,Sfc,Serno,cmsr(6)
      CHARACTER*80 title2,Title,kse
      COMMON /ccci  / Stci(125)
      COMMON /cfilt / Ktcopt,Lterm(4),Nterm
      COMMON /chrtm / Mach
      COMMON /extr  / Xs(129),Xi(125),Xwt(125)
      COMMON /inpt  / Xorg,Ib,Ie,Meas,Kp,Np
      COMMON /more  / Sfc,Bpf
      COMMON /mor2  / Kpart,Mcd
      COMMON /ngrz  / Igrz,Myr,Nw
      COMMON /qrat  / Rati(12)
      COMMON /rite  / Nr,Nyear,Title,Kpage,Ndec,Serno
      COMMON /tadj  / Add
      COMMON /tests / F(2),P(2),Test1,Test2,Iqfail
      COMMON /xk    / Kswt,Lswv,Begyr
C  *****  CHANGE TO INCORPORATE M & Q DIAGNOSTICS
      COMMON /mq1   / Adri,R,Ratic
      COMMON /mq2   / Vi,Vc,Vs,Rv,Length,Lstob
      COMMON /mq4   / Xtc,Lstyr
      CHARACTER*7 mq
      DATA mq/'QUARTER'/
C  *****  END OF CHANGE
      DATA pre,fin,bla/'PRELIM','FINAL ','      '/
      DATA(word(i),i=1,6)/' PERCE',' DIFFE','NT CHA','RENCE ','NGE   ',
     &     '      '/
      DATA rmsr/2.3,4.1,5.2,6.5,2.1,3.8,5.0,6.9/
      DATA cmsr/'3-TERM',' 3 X 3',' 3 X 5',' 3 X 9','STABLE',' DIFF.'/
      Bpf=bla
      Kpart=1
      nd2=Ndec
      IF(Add.eq.0.and.Ndec.eq.0)nd2=2
      DO i=Ib,Ie
       xtcsi(i)=Xorg(i)
      END DO
C --- WRITE ORIGINAL SERIES - TABLE B1
      IF(Issap.lt.2)CALL writ(Xorg,Ib,Ie,1,1,Ndec,4)
C --- TREND-CYCLE, CENTERED 4-TERM MOVING AVERAGE  (TABLE 2)
   10 ib1=Ib+2
      ie1=Ie-2
      DO i=ib1,ie1
       Xtc(i)=(xtcsi(i-2)+xtcsi(i+2)+2.0*(xtcsi(i-1)+xtcsi(i)+xtcsi(i+1)
     &        ))/8.0
      END DO
      GO TO(30,30,20,20,30,30,30,20,30,30,30,20),Kp
   20 IF(Issap.lt.2)CALL writ(Xtc,ib1,ie1,2,2,Ndec,4)
C --- SI RATIOS -- PRELIMINARY  (TABLE 3 OR 4)
   30 ltest=Lterm(1)
      DO il=2,4
       IF(ltest.ne.Lterm(il))ltest=6
      END DO
      Sfc=cmsr(2)
      IF(ltest.gt.0)Sfc=cmsr(ltest)
      jtble=Kpart/2+3
      CALL divsub(xtcsi,Xtc,xsi,ib1,ie1)
      GO TO(50,50,50,40,50,50,50,40,50,50,50,40),Kp
   40 IF(Issap.lt.2)CALL writ(xsi,ib1,ie1,3,jtble,nd2,1)
   50 GO TO(60,80,80),Kpart
C --- PART B - REPLACE  EXTREME SI RATIOS - TABLE B4
   60 CALL si(xsi,ib1,ie1,xtemp,1)
      GO TO(80,80,80,70),Kp
   70 IF(Issap.lt.2)CALL writ(xtemp,ib1,ie1,4,4,nd2,2)
C --- SEASONAL FACTORS -- PRELIMINARY (TABLE 5)
   80 CALL qsf(xsi,ib1,ie1,1)
C        FILL IN MISSING TERMS
      Xs(Ib)=Xs(Ib+4)
      Xs(Ib+1)=Xs(Ib+5)
      Xs(Ie)=Xs(Ie-4)
      Xs(Ie-1)=Xs(Ie-5)
      GO TO(100,100,100,90,100,100,100,90,100,100,100,90),Kp
   90 IF(Issap.lt.2)CALL writ(Xs,Ib,Ie,5,5,nd2,1)
C --- SEASONALLY ADJUST DATA -- PRELIMINARY  (TABLE 6)
  100 CALL divsub(xtcsi,Xs,xtci,Ib,Ie)
      GO TO(120,120,120,110,120,120,120,110,120,120,120,110),Kp
  110 IF(Issap.lt.2)CALL writ(xtci,Ib,Ie,6,6,Ndec,4)
C --- TREND-CYCLE CURVE - TABLE 7
  120 CALL vtcq(Xtc,xtci,Ib,Ie)
      GO TO(140,140,140,130,140,140,130,130,140,140,130,130),Kp
  130 IF(Issap.lt.2)CALL writ(Xtc,Ib,Ie,7,7,Ndec,4)
C --- SI RATIOS  (TABLE 8 OR 9)
  140 Sfc=cmsr(3)
      IF(ltest.gt.0)Sfc=cmsr(ltest)
      GO TO(150,190,210),Kpart
C        PART B
  150 CALL divsub(Xorg,Xtc,xsi,Ib,Ie)
      GO TO(170,170,170,160),Kp
  160 IF(Issap.lt.2)CALL writ(xsi,Ib,Ie,3,8,nd2,1)
C --- PART B - REPLACE  EXTREME SI RATIOS - TABLE B9
  170 CALL si(xsi,Ib,Ie,xtemp,2)
      GO TO(220,220,220,180),Kp
  180 IF(Issap.lt.2)CALL writ(xtemp,Ib,Ie,4,9,nd2,2)
      GO TO 220
C        PART C
  190 CALL divsub(xtcsi,Xtc,xsi,Ib,Ie)
      GO TO(220,220,220,220,220,220,220,200),Kp
  200 IF(Issap.lt.2)CALL writ(xsi,Ib,Ie,3,9,nd2,1)
      GO TO 220
C        PART D
  210 Bpf=fin
C --- FINAL UNMODIFIED SI RATIOS - TABLE D8
      CALL divsub(Xorg,Xtc,xsit,Ib,Ie)
      IF(Issap.lt.2)CALL writ(xsit,Ib,Ie,3,8,nd2,1)
      CALL ftestq(xsit,Ib,Ie,0,Issap,icol)
      CALL mstest(xsit,Ib,Ie,Nr,4,Issap,icol)
      CALL combft(Issap,icol,Nr)
      IF(Igrz.eq.2.and.Issap.lt.2)THEN
       title2='D8. FINAL UNMODIFIED SI-RATIOS'
       CALL grzlst(Ib,ibeg2,Ie,Myr,xsit,xsit,108)
       CALL chrt(title2,7,1,ibeg2,Ie)
      END IF
C --- FINAL REPLACEMENT VALUES FOR EXTREME SI RATIOS - TABLE D9
      CALL divsub(xtcsi,Xtc,xsi,Ib,Ie)
      DO i=Ib,Ie
       IF(xsit(i).ne.xsi(i))THEN
        xtemp(i)=xsi(i)
       ELSE
        xtemp(i)=9999999999.
       END IF
      END DO
      IF(Issap.lt.2)CALL writ(xtemp,Ib,Ie,4,9,nd2,1)
C --- SEASONAL FACTORS  (TABLE 10)
  220 CALL qsf(xsi,Ib,Ie,2)
      IF(Kpart.eq.3.and.Issap.lt.2)THEN
       WRITE(Nw,1010)
 1010  FORMAT(' D 9A. YEAR TO YEAR CHANGE IN IRREGULAR AND SEASONAL',
     &        ' COMPONENTS AND MOVING SEASONALITY RATIO',/,23X,'1ST',
     &        21X,'2ND',21X,'3RD',21X,'4TH')
       kde=3
       IF(Add.ne.0)kde=2
       WRITE(kse,1020)kde,kde
 1020  FORMAT('(4X,6H  I   ,6X,F10.',I1,',3(14X,F10.',I1,'))')
       WRITE(Nw,kse)(Rati(i),i=1,4)
       kse(9:9)='S'
       WRITE(Nw,kse)(Rati(i),i=5,8)
       kse(7:12)='RATIO '
       WRITE(Nw,kse)(Rati(i),i=9,12)
       r1=0
       r2=0
       DO ir=1,4
        r1=r1+Rati(ir)
        r2=r2+Rati(ir+4)
       END DO
       R=r1/r2
       WRITE(Nw,1030)R
 1030  FORMAT(/,4X,'OVERALL MOVING SEASONALITY RATIO : ',12F9.3)
       IF(Lswv.eq.0)THEN
        il=0
        llda1=Ie-Ib+1
        IF(llda1.ge.180)il=4
        IF(R.lt.rmsr(il+1))il2=1
        IF(R.lt.rmsr(il+2).and.R.ge.rmsr(il+1))il2=2
        IF(R.lt.rmsr(il+3).and.R.ge.rmsr(il+2))il2=3
        IF(R.lt.rmsr(il+4).and.R.ge.rmsr(il+3))il2=4
        IF(R.ge.rmsr(il+4))il2=5
        WRITE(Nw,1040)cmsr(il2)
 1040   FORMAT(4X,'TRY A ',A6,' MOVING AVERAGE FOR ALL QUARTERS')
       END IF
      END IF
      GO TO(240,240,230,230,240,240,230,230,230,230,230,230),Kp
  230 IF(Issap.lt.2)CALL writ(Xs,Ib,Ie,5,10,nd2,1)
      IF(Kpart.eq.3)THEN
C --- SEASONAL FACTORS ONE YEAR AHEAD
       j=Ie+1
       Xs(j)=(3.0*Xs(Ie-3)-Xs(Ie-7))/2.0
       Xs(j+1)=(3.0*Xs(Ie-2)-Xs(Ie-6))/2.0
       Xs(j+2)=(3.0*Xs(Ie-1)-Xs(Ie-5))/2.0
       Xs(j+3)=(3.0*Xs(Ie)-Xs(Ie-4))/2.0
       iyear=Nyear
       modq=mod(j,4)+1
       IF(modq.gt.1)modq=0
       Nyear=Nyear+modq
       j3=j+3
       IF(Issap.lt.2)CALL writ(Xs,j,j3,11,0,nd2,1)
       Nyear=iyear
       IF(Igrz.eq.2.and.Issap.lt.2)THEN
        title2='D10. FINAL SEASONAL FACTORS WITH ONE YEAR FORECAST'
        CALL grzlst(Ib,ibeg2,Ie+4,Myr,Xs,Xs,108)
        CALL chrt(title2,7,1,ibeg2,Ie+4)
       END IF
       IF(Issap.eq.2)CALL ssrit(Xs,Ib,Ie,1,129,Xorg,Lstob)
      END IF
C --- SEASONALLY ADJUST  (TABLE 11)
  240 CALL divsub(Xorg,Xs,xtci,Ib,Ie)
      GO TO(270,270,250),Kpart
  250 GO TO(270,260),Meas
C --- IF PART D AND SUMMARY MEASURES RUN, REPLACE SEASONALLY ADJUSTED
C     SERIES BY ORIGINAL SERIES
  260 DO i=Ib,Ie
       xtci(i)=Xorg(i)
      END DO
  270 GO TO(320,320,320,280,320,320,320,280,280,280,280,280),Kp
  280 IF(Issap.lt.2)CALL writ(xtci,Ib,Ie,6,11,Ndec,4)
      CALL ftestq(xtci,Ib,Ie,1,Issap,icol)
      IF(Issap.eq.2)CALL ssrit(xtci,Ib,Ie,2,125,Xorg,Lstob)
      GO TO(320,320,290),Kpart
C --- FINAL TREND-CYCLE, TABLE D12
C --- DIVIDE MODIFIED ORIGINAL SERIES BY FINAL SEASONAL FACTORS
  290 CALL divsub(xtcsi,Xs,xtemp,Ib,Ie)
      GO TO(310,300),Meas
C --- IF SUMMARY MEASURES RUN, COMPUTE MOVING AVERAGE OF
C     OF MODIFIED ORIGINAL
  300 DO i=Ib,Ie
       xtemp(i)=xtcsi(i)
      END DO
  310 CALL vtcq(Xtc,xtemp,Ib,Ie)
      IF(Issap.lt.2.and.Kp.ne.9)CALL writ(Xtc,Ib,Ie,7,12,Ndec,4)
C --- IRREGULAR  (TABLE 13)
  320 CALL divsub(xtci,Xtc,Xi,Ib,Ie)
      GO TO(350,350,330,330,350,350,330,330,340,330,330,330),Kp
  330 IF(Issap.lt.2)CALL writ(Xi,Ib,Ie,8,13,nd2,3)
  340 GO TO(350,350,410),Kpart
C --- WEIGHTS FOR IRREGULAR - TABLES B17 AND C17
  350 CALL xtrmes(Xi,Ib,Ie,Xwt)
      Bpf=pre
      GO TO(380,380,380,370,360,360,360,360),Kp
  360 Bpf=fin
  370 IF(Issap.lt.2)CALL writ(Xwt,Ib,Ie,9,17,nd2,2)
C --- MODIFY ORIGINAL SERIES BY WEIGHTS - TABLES C1 AND D1
  380 DO i=Ib,Ie
       xtcsi(i)=Xorg(i)
       IF(Xwt(i).lt.1.0)THEN
        IF(Add.le.0)THEN
         xtcsi(i)=(Xorg(i)*(100.+Xwt(i)*(Xi(i)-100.)))/Xi(i)
        ELSE
         xtcsi(i)=Xorg(i)-(Xi(i)*(1.0-Xwt(i)))
        END IF
       END IF
      END DO
      Kpart=Kpart+1
      GO TO(400,400,390,390,400,400,390,390),Kp
  390 IF(Issap.lt.2)CALL writ(xtcsi,Ib,Ie,10,1,Ndec,4)
  400 Bpf=bla
      Kp=Kp+4
      GO TO 10
C
C --- PART E
C
C --- MODIFY ORIGINAL,SEASONALLY ADJUSTED,AND IRREGULAR FOR EXTREMES
C --- TABLES E1, E2, AND E3
  410 IF(Kp.eq.9)RETURN
      IF(Igrz.ne.0.or.Issap.eq.2)THEN
       title2=
     &'E0.ORIGINAL SERIES,@-SEASONALLY ADJUSTED SERIES,COINCIDENCE OF PO
     &INTS'
       IF(Mach.eq.0)title2(20:20)=char(223)
       CALL grzlst(Ib,ibeg2,Ie,Myr,Xorg,xtci,108)
       CALL chrt(title2,15,2,ibeg2,Ie)
      END IF
      DO i=Ib,Ie
       torg=Xorg(i)
       tsa=xtci(i)
       ti=Xi(i)
       IF(Xwt(i).ne.0)GO TO 430
       IF(Add.le.0)THEN
        torg=Xtc(i)*Xs(i)/100.
        ti=100.
       ELSE
        torg=Xtc(i)+Xs(i)
        ti=0.0
       END IF
       tsa=Xtc(i)
       GO TO(430,420),Meas
  420  torg=Xtc(i)
  430  xmtci(i)=tsa
       xmi(i)=ti
       xmorg(i)=torg
      END DO
      Kpart=4
      IF(Kp.eq.10)GO TO 470
      IF(Issap.lt.2)CALL writ(xmorg,Ib,Ie,1,1,Ndec,4)
      GO TO(440,450),Meas
  440 IF(Issap.lt.2)CALL writ(xmtci,Ib,Ie,2,2,Ndec,4)
      CALL ftestq(xmtci,Ib,Ie,1,Issap,icol)
      IF(Issap.eq.2)RETURN
  450 CALL writ(xmi,Ib,Ie,3,3,nd2,3)
C --- RATIOS OF ANNUAL TOTALS - TABLE E4
      GO TO(460,470),Meas
  460 IF(Kp.ne.11)THEN
       CALL writ(Xorg,Ib,Ie,4,4,Ndec,4)
       ib1=(Ib+2)/4*4+1
       ie1=Ie/4*4
       kyear=Nyear+ib1/4
       DO i=ib1,ie1,4
        IF(Add.le.0)THEN
         et1=(Xorg(i)+Xorg(i+1)+Xorg(i+2)+Xorg(i+3))
     &       /(xtci(i)+xtci(i+1)+xtci(i+2)+xtci(i+3))
         et3=(xmorg(i)+xmorg(i+1)+xmorg(i+2)+xmorg(i+3))
     &       /(xmtci(i)+xmtci(i+1)+xmtci(i+2)+xmtci(i+3))
        ELSE
         et1=(Xorg(i)+Xorg(i+1)+Xorg(i+2)+Xorg(i+3))
     &       -(xtci(i)+xtci(i+1)+xtci(i+2)+xtci(i+3))
         et3=(xmorg(i)+xmorg(i+1)+xmorg(i+2)+xmorg(i+3))
     &       -(xmtci(i)+xmtci(i+1)+xmtci(i+2)+xmtci(i+3))
        END IF
        WRITE(Nr,1050)kyear,et1,et3
 1050   FORMAT(4X,I4,15X,F7.2,29X,F7.2)
        kyear=kyear+1
       END DO
      END IF
C --- CHANGES IN ORGINAL SERIES  (TABLE E5)
  470 ib1=Ib+1
      ie1=Ie-1
      DO i=Ib,ie1
       xtemp(i+1)=Xorg(i+1)-Xorg(i)
       IF(Add.le.0)xtemp(i+1)=xtemp(i+1)/Xorg(i)*100.
      END DO
      CALL writ(xtemp,ib1,Ie,5,5,nd2,2)
      GO TO(480,490),Meas
C --- CHANGES IN FINAL SEASONALLY ADJUSTED SERIES  (TABLE E6)
  480 DO i=Ib,ie1
       xtemp(i+1)=xtci(i+1)-xtci(i)
       IF(Add.le.0)xtemp(i+1)=xtemp(i+1)/xtci(i)*100.0
      END DO
      CALL writ(xtemp,ib1,Ie,6,6,nd2,2)
C
C --- PART F
C
  490 Kpart=5
      DO i=1,125
       Stci(i)=xtci(i)
      END DO
      CALL avec(xtci,Ib,Ie,cibar,cibar2,sdcci)
      CALL avedur(xtci,Ib,Ie,adrci)
      CALL avedur(Xtc,Ib,Ie,adrc)
      CALL avedur(Xi,Ib,Ie,Adri)
C --- FIND QCD AND COMPUTE QCD MOVING AVERAGE
      CALL avec(Xtc,Ib,Ie,cbar,cbar2,sdcc)
      CALL avec(Xi,Ib,Ie,ibar,ibar2,sdci)
      DO i=1,4
       ratio(i)=ibar(i)/cbar(i)
      END DO
      i=3
      DO WHILE (ratio(i).le..995)
       i=i-1
       IF(i.le.0)THEN
        Mcd=1
        DO i=Ib,Ie
         xtemp(i)=xtci(i)
        END DO
        nt1=Ib
        nt2=Ie
        GO TO 500
       END IF
      END DO
      Mcd=i+1
      nt1=Ib+Mcd/2
      nt2=Ie-(Mcd-1)/2
      j=(Mcd+2)/2
      xmcd=Mcd
      DO i=nt1,nt2
       tmp=0.
       DO i2=1,Mcd
        l=i+i2-j
        tmp=tmp+Stci(l)
       END DO
       xtemp(i)=tmp/xmcd
      END DO
C --- COMPUTE AVGE MO-TO-MO CHANGES IN MCD FOR SUMMARY MEASURES TABLE
  500 CALL avedur(xtemp,nt1,nt2,adrqcd)
      CALL avec(xtemp,nt1,nt2,qbar,qbar2,sdcq)
      IF(Kp.eq.12)CALL writ(xtemp,nt1,nt2,1,1,Ndec,4)
C --- SUMMARY MEASURES
      CALL writ(Xorg,0,0,2,2,Ndec,4)
      CALL avec(Xorg,Ib,Ie,obar,obar2,sdco)
      CALL avec(Xs,Ib,Ie,sbar,sbar2,sdcs)
      CALL avec(xmorg,Ib,Ie,mobar,xtem,xtem2)
      CALL avec(xmtci,Ib,Ie,mcibar,xtem,xtem2)
      CALL avec(xmi,Ib,Ie,mibar,xtem,xtem2)
      WRITE(Nr,1060)word(Add+1),word(Add+3),word(Add+5)
 1060 FORMAT(/,' F 2.A:',2X,'AVERAGE',2A6,A3,' WITHOUT REGARD TO SIGN ',
     &       'OVER INDICATED SPAN',/,14X,'SPAN',22X,
     &       'TABLE NO. AND SYMBOL',37X,'MODIFIED SERIES',/,15X,'IN',9X,
     &       'B1',8X,'D11',7X,'D13',7X,'D12',7X,'D10',7X,'F1',18X,'E1',
     &       8X,'E2',8X,'E3',/,12X,'QUARTERS',7X,'O',8X,'CI',9X,'I',9X,
     &       'C',9X,'S',8X,'QCD',18X,'O',8X,'CI',9X,'I')
      DO i=1,4
       WRITE(Nr,1070)i,obar(i),cibar(i),ibar(i),cbar(i),sbar(i),qbar(i),
     &               mobar(i),mcibar(i),mibar(i)
      END DO
 1070 FORMAT(16X,I1,3X,6F10.2,10X,3F10.2)
      GO TO(510,520),Meas
  510 DO i=1,4
       temps=cbar(i)**2+ibar(i)**2+sbar(i)**2
       ric(i)=cbar(i)**2/temps*100.
       rii(i)=ibar(i)**2/temps*100.
       ris(i)=sbar(i)**2/temps*100.
       rtot(i)=ric(i)+rii(i)+ris(i)
       rato(i)=(temps/obar(i)**2)*100.
      END DO
      WRITE(Nr,1080)word(Add+1),word(Add+3),word(Add+5)
 1080 FORMAT(/,' F 2.B:',2X,'RELATIVE CONTRIBUTIONS OF COMPONENTS TO',
     &       2A6,A4,'IN ORIGINAL SERIES',/,14X,'SPAN',7X,
     &       'TABLE NO. AND ','SYMBOL',/,15X,'IN',9X,'D13',7X,'D12',7X,
     &       'D10',/,12X,'QUARTERS',7X,'I',9X,'C',9X,'S',7X,'TOTAL',5X,
     &       'RATIO')
      DO i=1,4
       WRITE(Nr,1090)i,rii(i),ric(i),ris(i),rtot(i),rato(i)
      END DO
 1090 FORMAT(16X,I1,3X,5F10.2)
  520 WRITE(Nr,1100)adrci,Adri,adrc,adrqcd
 1100 FORMAT(/,' F 2.C:',2X,'AVERAGE DURATION OF RUN',6X,'CI',9X,'I',9X,
     &       'C',7X,'QCD',/,30X,4F10.2)
      WRITE(Nr,1110)(ratio(i),i=1,4),Mcd
 1110 FORMAT(/,' F 2.D:',2X,'I/C RATIO FOR QUARTERS SPAN',/,19X,'1',9X,
     &       '2',9X,'3',9X,'4',/,10X,4F10.2,//,7X,
     &       'QUARTERS FOR CYCLICAL DOMINANCE',11X,I1)
      WRITE(Nr,1120)word(Add+1),word(Add+3),word(Add+5)
 1120 FORMAT(/,' F 2.E:',2X,'AVERAGE',2A6,A4,'WITH REGARD TO SIGN AND ',
     &       'STANDARD DEVIATION OVER INDICATED SPAN',/,14X,'SPAN',11X,
     &       ' B1',13X,'D13',13X,'D12',13X,'D10',13X,'D11',13X,' F1',/,
     &       15X,'IN',13X,'O',15X,'I',15X,'C',15X,'S',15X,'CI',13X,
     &       'QCD',/,12X,'QUARTERS     AVG    S.D.     AVG    S.D.    ',
     &    ' AVG    S.D.     AVG    S.D.     AVG    S.D.     AVG    S.D.'
     &    )
      DO i=1,4
       WRITE(Nr,1130)i,obar2(i),sdco(i),ibar2(i),sdci(i),cbar2(i),
     &               sdcc(i),sbar2(i),sdcs(i),cibar2(i),sdcci(i),
     &               qbar2(i),sdcq(i)
      END DO
 1130 FORMAT(16X,I1,3X,12F8.2)
C  *****  CHANGES TO INCORPORATE M & Q STATISTICS
      CALL vcomp(4,xmi,xmorg)
      Kpage=Kpage+1
      WRITE(Nw,1140)Title,Kpage,Serno
 1140 FORMAT('1',17X,A80,1X,'PAGE ',I2,', SERIES ',A6)
      WRITE(Nw,1150)Vi,Vc,Vs,Rv
 1150 FORMAT(//,' F 2.F:',2X,'RELATIVE CONTRIBUTION OF THE ',
     &       'COMPONENTS TO THE STATIONARY PORTION OF THE VARIANCE IN ',
     &       'THE ORIGINAL SERIES',/,24X,'I',7X,'C',7X,'S',4X,'TOTAL',/,
     &       19X,4F8.2,/)
      WRITE(Nw,1160)Ratic,R
 1160 FORMAT(/,' F 2.G:',2X,'THE FINAL I/C RATIO FROM TABLE D12:',F12.2,
     &       /,9X,'THE FINAL I/S RATIO FROM TABLE D10:',F12.2)
      WRITE(Nw,1170)
 1170 FORMAT(/,' F 2.H:',77X,'STATISTIC PROBABILITY',/,100X,'LEVEL')
      WRITE(Nw,1180)F(1),P(1),F(2),P(2)
 1180 FORMAT(10X,'F-TEST FOR STABLE SEASONALITY FROM TABLE D 8.',26X,
     &       ':',F11.3,F8.2,'%',/,10X,
     &       'F-TEST FOR MOVING SEASONALITY FROM TABLE D 8.',26X,':',
     &       F11.3,F8.2,'%')
      CALL f3cal(rii,ratio,4)
      CALL f3gen(4,mq)
C  *****  END OF CHANGES
      RETURN
      END
**==ssrit.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE ssrit(X,L1,L2,Isec,Nd,Stor,Lstob)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,i0,icol,Isec,L0,L1,l10,L2,l20,la,Len,Len2,ll0,Lstob,
     &        Muladd,Ncol,Nd,Nlen,Ns1
      REAL S,Sa,Ten,xss
C*** End of declarations inserted by SPAG
C  *****
C  *****  Stores results from X-11Q.2 runs into sliding spans variables.
C  *****
      DOUBLE PRECISION X(Nd),Stor(125)
      COMMON /kten  / Ten
      COMMON /lzero / L0
      COMMON /ns4   / Nlen,Ncol,Ns1,Len,Len2
      COMMON /sser  / S(60,4),Sa(60,4)
      COMMON /tadj  / Muladd
      Ten=10.0**10
      l10=L1-L0
      icol=(l10/4)+1
      IF(l10.ne.0)THEN
       DO i=1,l10
        IF(Isec.eq.1)S(i,icol)=Ten
        IF(Isec.eq.2)Sa(i,icol)=Ten
       END DO
      END IF
      DO i0=L1,L2
       i=i0-L0+1
       xss=sngl(X(i0))
       IF(Isec.eq.1)THEN
        S(i,icol)=xss
       ELSE
        Sa(i,icol)=xss
        IF(Muladd.eq.1)S(i,icol)=(sngl(Stor(i0))/xss)*100.0
       END IF
      END DO
      la=mod(L1,4)
      IF(la.eq.0)la=4
      ll0=Len+L0+la-2
      IF(mod(ll0,4).ne.0)THEN
       l20=((ll0/4)+1)*4
      ELSE
       l20=ll0
      END IF
      IF(l20.eq.L2.and.Lstob.lt.4)l20=(4-Lstob)+l20
      IF(l20.ne.L2)THEN
       DO i=L2+1,l20
        i0=i-L0+1
        IF(Isec.eq.1)S(i0,icol)=Ten
        IF(Isec.eq.2)Sa(i0,icol)=Ten
       END DO
      END IF
      RETURN
      END
**==mvavg.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C --- LIP FOR MVAVG,MVAVG,MVAVGS
C --- THIS ROUTINE COMPUTES WEIGHTED MOVING AVERAGES FOR VSFC AND VTC.
C --- WEIGHTS ARE ORDERED SO THAT FIRST N WEIGHTS ARE APPLIED TO VALUES
C --- FOR FIRST AND LAST MOVING AVERAGE TERMS. LAST N VALUES ARE FOR
C --- MIDDLE TERMS. WITHIN EACH SET OF N VALUES FOR END TERMS, FIRST
C --- WEIGHT IS APPLIED TO FIRST AND LAST INPUT VALUE.
      SUBROUTINE mvavg(Stb,Sta,Lf,Ll,N,Wt)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,ii,j,jj,kk,Lf,Ll,m,mm,N,nn
      DOUBLE PRECISION Sta,Stb,swt,Wt
C*** End of declarations inserted by SPAG
      DIMENSION Sta(*),Stb(*),Wt(*)
      m=N/2
      DO i=1,m
       ii=Lf+i-1
       jj=Ll-i+1
       Stb(ii)=0.
       Stb(jj)=0.
       kk=Lf
       mm=Ll
       nn=N*(i-1)
       swt=0.
       DO j=1,N
        nn=nn+1
        Stb(ii)=Stb(ii)+Sta(kk)*Wt(nn)
        Stb(jj)=Stb(jj)+Sta(mm)*Wt(nn)
        swt=swt+Wt(nn)
        kk=kk+1
        mm=mm-1
       END DO
       Stb(ii)=Stb(ii)/swt
       Stb(jj)=Stb(jj)/swt
      END DO
      IF(N-m.ne.m)THEN
       ii=ii+1
       jj=jj-1
       DO i=ii,jj
        kk=i-m-1
        Stb(i)=0.
        swt=0.
        DO j=1,N
         kk=kk+1
         mm=nn+j
         swt=swt+Wt(mm)
         Stb(i)=Stb(i)+Sta(kk)*Wt(mm)
        END DO
        Stb(i)=Stb(i)/swt
       END DO
      END IF
      RETURN
      END
**==qsf.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C --- LIP FOR SUMRY,SUMRY,SUMRYS
      SUBROUTINE qsf(Stsi,Lfda,Llda,Ksect)
C --- SEASONAL FACTOR CURVE ROUTINE.
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,j,k,kfda,ki,Kpart,Ksect,Ktcopt,l,Lfda,Llda,Lterm,m,Mcd,
     &        mfda,mlda,Muladd,Nterm
      DOUBLE PRECISION p,Rati,stavg,stmp,Sts,Stsi,stsp,sum,tmp,w3a3,
     &                 w3a5,w3a9,w3b3,w3b5,w3b9,w3c5,w3c9,w3d9,w3e9,w3x1
      DOUBLE PRECISION w3x3,w3x5,w3x9,Xi,Xwt
C*** End of declarations inserted by SPAG
      COMMON /cfilt / Ktcopt,Lterm(4),Nterm
      COMMON /extr  / Sts(129),Xi(125),Xwt(125)
      COMMON /mor2  / Kpart,Mcd
      COMMON /tadj  / Muladd
      COMMON /qrat  / Rati(12)
      DIMENSION w3x9(66),w3a9(50),w3b9(45),w3c9(32),w3d9(28),w3e9(18),
     &          w3x5(28),w3a5(18),w3b5(15),w3c5(8),w3x3(15),w3a3(8),
     &          w3b3(6),w3x1(6),stsp(57),stavg(57),stmp(300),Stsi(125)
C --- WEIGHTS FOR 3X9 MOVING AVGE. WEIGHTS IN W3A9-W3E9 ARE FOR SERIES
C --- SHORTER THAN 11 YEARS.
      DATA w3x9/.246,.221,.197,.173,.112,.051,.0,.0,.0,.0,.0,.208,.192,
     &     .176,.160,.144,.092,.028,.0,.0,.0,.0,.173,.163,.154,.143,
     &     .133,.123,.079,.032,.0,.0,.0,.141,.137,.132,.128,.123,.117,
     &     .113,.075,.034,.0,.0,.084,.120,.118,.117,.116,.114,.113,.111,
     &     .073,.034,.0,.037,.074,.111,.111,.111,.111,.111,.111,.111,
     &     .074,.037/
      DATA w3a9/.246,.221,.197,.173,.112,.051,.0,.0,.0,.0,.208,.192,
     &     .176,.160,.144,.092,.028,.0,.0,.0,.173,.163,.154,.143,.133,
     &     .123,.079,.032,.0,.0,.141,.137,.132,.128,.123,.117,.113,.075,
     &     .034,.0,.084,.120,.118,.117,.116,.114,.113,.111,.073,.034/
      DATA w3b9/.246,.221,.197,.173,.112,.051,.0,.0,.0,.208,.192,.176,
     &     .160,.144,.092,.028,.0,.0,.173,.163,.154,.143,.133,.123,.079,
     &     .032,.0,.141,.137,.132,.128,.123,.117,.113,.075,.034,.111,
     &     .111,.111,.111,.111,.111,.111,.111,.111/
      DATA w3c9/.246,.221,.197,.173,.112,.051,.0,.0,.208,.192,.176,.160,
     &     .144,.092,.028,.0,.173,.163,.154,.143,.133,.123,.079,.032,
     &     .125,.125,.125,.125,.125,.125,.125,.125/
      DATA w3d9/.246,.221,.197,.173,.112,.051,.0,.208,.192,.176,.160,
     &     .144,.092,.028,.143,.143,.143,.143,.143,.143,.143,.143,.143,
     &     .143,.143,.143,.143,.143/
      DATA w3e9/.246,.221,.197,.173,.112,.051,.167,.167,.167,.167,.167,
     &     .167,.167,.167,.167,.167,.167,.167/
C --- WEIGHTS FOR 3X5 MOVING AVERAGE.
      DATA w3x5/.283,.283,.283,.150,.0,.0,.0,.250,.250,.250,.183,.067,
     &     .0,.0,.150,.217,.217,.217,.133,.067,.0,.067,.133,.200,.200,
     &     .200,.133,.067/
      DATA w3a5/.283,.283,.283,.150,.0,.0,.250,.250,.250,.183,.067,.0,
     &     .150,.217,.217,.217,.133,.067/
      DATA w3b5/.283,.283,.283,.150,.0,.250,.250,.250,.183,.067,.200,
     &     .200,.200,.200,.200/
      DATA w3c5/.283,.283,.283,.150,.250,.250,.250,.250/
C --- WEIGHTS FOR 3X3 MOVING AVERAGE.
      DATA w3x3/.407,.407,.185,.0,.0,.259,.370,.259,.111,.0,.111,.222,
     &     .333,.222,.111/
      DATA w3a3/.407,.407,.185,.0,.259,.370,.259,.111/
      DATA w3b3/.407,.407,.185,.333,.333,.333/
C --- WEIGHTS FOR 3-TERM MOVING AVERAGE.
      DATA w3x1/.610,.390,.0,.333,.333,.333/
      mfda=Lfda
      mlda=Llda
      kfda=mfda+3
C --- FOR EACH MONTH SEPARATELY.
      DO j=mfda,kfda
       m=j-(j-1)/4*4
       k=3
C --- MOVE MONTH TO *STSP*.
       DO i=j,mlda,4
        k=k+1
        stsp(k)=Stsi(i)
       END DO
C --- COMPUTE IBAR,SBAR,AND RATIO IN FINAL ROUND ONLY.
C --- COMPUTE 7-TERM MOVING AVERAGE FOR ESTIMATE OF S.
       IF(Kpart.eq.3)THEN
        IF(Ksect.eq.2)THEN
         DO i=1,3
          stsp(i)=(stsp(4)+stsp(5)+stsp(6))/3.
          ki=k+i
          stsp(ki)=(stsp(k)+stsp(k-1)+stsp(k-2))/3.
         END DO
         DO i=4,k
          stavg(i)=(stsp(i-3)+stsp(i-2)+stsp(i-1)+stsp(i)+stsp(i+1)
     &             +stsp(i+2)+stsp(i+3))/7.
         END DO
C --- COMPUTE SBAR.
         Rati(m+4)=0.
         DO i=5,k
          IF(Muladd.ne.0)THEN
           Rati(m+4)=Rati(m+4)+abs(stavg(i-1)-stavg(i))
          ELSE
           Rati(m+4)=Rati(m+4)+abs(stavg(i-1)-stavg(i))/stavg(i)*100.
          END IF
         END DO
         p=k-4
         Rati(m+4)=Rati(m+4)/p
C --- DIVIDE SI/S FOR ESTIMATE OF I
         CALL addmul(stavg,stsp,stavg,4,k)
C --- COMPTUE IBAR
         Rati(m)=0.
         DO i=5,k
          IF(Muladd.ne.0)THEN
           Rati(m)=Rati(m)+abs(stavg(i-1)-stavg(i))
          ELSE
           Rati(m)=Rati(m)+abs(stavg(i-1)-stavg(i))/stavg(i)*100.
          END IF
         END DO
         Rati(m)=Rati(m)/p
C --- COMPUTE IBAR/SBAR RATIO. LARGEST RATIO PRINTED IS 999.99.
         Rati(m+8)=999.99
         IF(Rati(m).lt.Rati(m+4)*999.99)Rati(m+8)=Rati(m)/Rati(m+4)
        END IF
       END IF
       Nterm=Lterm(m)
C --- CHECK IF MOVING AVERAGE FOR THIS MONTH PRESELECTED.
       IF(Nterm.le.0)THEN
        IF(Ksect.eq.1)GO TO 40
        GO TO 30
       ELSE
        GO TO(50,40,30,20,10),Nterm
       END IF
C --- STABLE SEASONAL. AVERAGE OF ALL SI RATIOS FOR THIS MONTH.
   10  tmp=0.
       DO i=4,k
        tmp=tmp+stsp(i)
       END DO
       sum=k-3
       tmp=tmp/sum
       DO i=4,k
        stavg(i)=tmp
       END DO
       GO TO 60
C --- COMPUTE 3X9 MOVING AVERAGE.
   20  IF(k.lt.13)THEN
        IF(k.lt.11)THEN
         IF(k.lt.9)GO TO 10
         IF(k.eq.9)THEN
          CALL mvavg(stavg,stsp,4,k,6,w3e9)
         ELSE
          CALL mvavg(stavg,stsp,4,k,7,w3d9)
         END IF
        ELSE IF(k.eq.11)THEN
         CALL mvavg(stavg,stsp,4,k,8,w3c9)
        ELSE
         CALL mvavg(stavg,stsp,4,k,9,w3b9)
        END IF
       ELSE IF(k.eq.13)THEN
        CALL mvavg(stavg,stsp,4,k,10,w3a9)
       ELSE
        CALL mvavg(stavg,stsp,4,k,11,w3x9)
       END IF
       GO TO 60
C --- COMPUTE 3X5 MOVING AVERAGE.
   30  IF(k.lt.9)THEN
        IF(k.lt.7)GO TO 10
        IF(k.eq.7)THEN
         CALL mvavg(stavg,stsp,4,k,4,w3c5)
        ELSE
         CALL mvavg(stavg,stsp,4,k,5,w3b5)
        END IF
       ELSE IF(k.eq.9)THEN
        CALL mvavg(stavg,stsp,4,k,6,w3a5)
       ELSE
        CALL mvavg(stavg,stsp,4,k,7,w3x5)
       END IF
       GO TO 60
C --- COMPUTE 3X3 MOVING AVERAGE.
   40  IF(k.lt.7)THEN
        IF(k.le.5)GO TO 10
        CALL mvavg(stavg,stsp,4,k,3,w3b3)
        GO TO 60
       ELSE IF(k.eq.7)THEN
        CALL mvavg(stavg,stsp,4,k,4,w3a3)
        GO TO 60
       ELSE
        CALL mvavg(stavg,stsp,4,k,5,w3x3)
        GO TO 60
       END IF
C --- COMPUTE 3-TERM MOVING AVERAGE
   50  CALL mvavg(stavg,stsp,4,k,3,w3x1)
   60  l=j
       DO i=4,k
        stmp(l)=stavg(i)
        l=l+4
       END DO
      END DO
C --- COMPUTE CENTERED 4-TERM MOVING AVERAGE.
      j=mfda+2
      k=mlda-2
      DO i=j,k
       Sts(i)=(stmp(i-2)+stmp(i+2)+2.*(stmp(i-1)+stmp(i)+stmp(i+1)))/8.
      END DO
C --- FILL IN MISSING END TERMS BY REPEATING FIRST AND LAST AVAILABLE
C --- MOVING AVERAGE VALUES 2 TIMES.
      Sts(mfda)=Sts(mfda+2)
      Sts(mfda+1)=Sts(mfda+2)
      Sts(mlda)=Sts(mlda-2)
      Sts(mlda-1)=Sts(mlda-2)
C --- DIVIDE SEASONALS BY 4-TERM MOVING AVERAGE.
      CALL addmul(Sts,stmp,Sts,mfda,mlda)
      RETURN
      END
**==vtcq.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
C --- LIP FOR VTC,VTC,VTCS
C --- VARIABLE TREND-CYCLE ROUTINE.
      SUBROUTINE vtcq(Stc,Stci,Lfda,Llda)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION Adri,apcc,apci,r,Ratic,Ratis,Stc,Stci,Sti,w5h,
     &                 w7h,Xs,Xwt
      INTEGER i,k,Ktcopt,Lfda,Llda,Lterm,mfda,mlda,Muladd,Nterm
C*** End of declarations inserted by SPAG
      DIMENSION w5h(15),w7h(28),Stc(125),Stci(125)
      COMMON /cfilt / Ktcopt,Lterm(4),Nterm
      COMMON /extr  / Xs(129),Sti(125),Xwt(125)
      COMMON /tadj  / Muladd
C  *****  CHANGE TO INCORPORATE M & Q DIAGNOSTICS
      COMMON /mq1   / Adri,Ratis,Ratic
C  *****  END OF CHANGE
C --- 55-TERM HENDERSON WEIGHTS.
      DATA(w5h(i),i=1,15)/.670,.403,-.073,.0,.0,.257,.522,.294,-.073,.0,
     &     -.073,.294,.558,.294,-.073/
C --- 7-TERM HENDERSON WEIGHTS.
      DATA(w7h(i),i=1,28)/.670,.403,-.073,.0,.0,.0,.0,.257,.522,.294,
     &     -.073,.0,.0,.0,-.073,.294,.558,.294,-.073,.0,.0,-.059,.059,
     &     .294,.413,.294,.059,-.059/
C --- CALCULATE 5-TERM HENDERSON. DROP END TERMS. CALC  IRREGULAR SER.
C --- CHECK IF TREND-CYCLE MOVING AVERAGE PRESELECTED.
      IF(Ktcopt.ne.1)THEN
       IF(Ktcopt.eq.2)GO TO 10
       mfda=Lfda+2
       mlda=Llda-2
       DO i=mfda,mlda
        Stc(i)=.559441*Stci(i)+.293706*(Stci(i+1)+Stci(i-1))
     &         -.073427*(Stci(i+2)+Stci(i-2))
       END DO
       CALL addmul(Sti,Stci,Stc,mfda,mlda)
C --- CALCULATE IBAR/CBAR RATIO
       apcc=0.
       apci=0.
       k=mlda-1
       DO i=mfda,k
        IF(Muladd.ne.0)THEN
         apcc=apcc+abs(Stc(i+1)-Stc(i))
         apci=apci+abs(Sti(i+1)-Sti(i))
        ELSE
         apcc=apcc+abs(Stc(i+1)-Stc(i))/Stc(i)
         apci=apci+abs(Sti(i+1)-Sti(i))/Sti(i)
        END IF
       END DO
       Ratic=apci/apcc
C --- CHECK IF TREND-CYCLE MOVING AVERAGE PRESELECTED.
       IF(Ktcopt.ne.1)THEN
        IF(Ktcopt.eq.2)GO TO 10
        r=Ratic*3
        IF(r.gt.3.495)GO TO 10
       END IF
      END IF
      Nterm=5
      CALL mvavg(Stc,Stci,Lfda,Llda,5,w5h)
      RETURN
   10 Nterm=7
      CALL mvavg(Stc,Stci,Lfda,Llda,7,w7h)
      RETURN
      END
**==addmul.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE addmul(Result,Array1,Array2,Jfda,Jlda)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION Array1,Array2,Result
      INTEGER i,Jfda,Jlda,Muladd
C*** End of declarations inserted by SPAG
      DIMENSION Result(130),Array1(130),Array2(130)
      COMMON /tadj  / Muladd
      IF(Muladd.eq.0)THEN
       DO i=Jfda,Jlda
        Result(i)=Array1(i)/Array2(i)*100.
       END DO
       RETURN
      END IF
      DO i=Jfda,Jlda
       Result(i)=Array1(i)-Array2(i)
      END DO
      RETURN
      END
**==tabset.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE tabset(Nf4,Ifmt,Irun,Reset,Prompt,*)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,i0,i1,ic,ic2,idx,Ifmt,Ifor,io,Io12,Iovr,irs,Irun,Irunt,
     &        it,itb,Itbl,itc1,Itchar,Itfmt
      INTEGER Ixerr,j,k2,Kpage,Kswerr,Kswt,ktbl,Lswv,Lyr,Nerr,Np,Nr,
     &        Numtbl,Nyear
C*** End of declarations inserted by SPAG
      INTEGER tbldec,Outdec
      LOGICAL fxst
      CHARACTER*1 forcst,cpart,stop,treset,overit,Reset,Prompt
      CHARACTER*3 Tbl
      CHARACTER*6 Ser
      CHARACTER*40 tblfmt,drive,Tf4,Nf4
      CHARACTER*80 Tblfil,tbf,f1,Title
      COMMON /itcom / Itchar,Itfmt,Itbl(10),Ifor,Numtbl,Iovr,Io12
      COMMON /rite  / Np,Nyear,Title,Kpage,Outdec,Ser
      COMMON /tcom  / Tf4,Tblfil,Tbl(10)
      COMMON /xerr  / Nr,Nerr,Kswerr,Irunt,Ixerr
      COMMON /xk    / Kswt,Lswv,Lyr
      NAMELIST /nltblq/ Numtbl,Tbl,tblfmt,drive,tbldec,forcst,treset,
     & overit,stop
      irs=0
      k2=Kswt
      IF(Irun.eq.0)treset=' '
      IF(Irun.eq.0.or.treset.eq.'Y'.or.treset.eq.'y')THEN
       treset=' '
       stop='N'
       tblfmt=' '
       tbldec=0
       forcst='N'
       drive=' '
       overit='Y'
       DO i=1,10
        Tbl(i)='   '
       END DO
      END IF
   10 Irun=1
      IF(Nr.eq.10)THEN
       READ(Nr,nltblq,ERR=80,END=90)
       IF(stop.ne.'N'.and.stop.ne.'n')RETURN 1
      ELSE
       IF(Kswt.eq.0)Kswt=k2
       PRINT *,' '
       PRINT *,'  CURRENT SETTINGS FOR NLTBLQ NAMELIST ARE : '
       PRINT nltblq
   20  PRINT *,' '
       PRINT *,'  ENTER VARIABLES AND THEIR SETTINGS FOR NLTBLQ',
     &       ' NAMELIST : '
       READ(*,nltblq,ERR=100)
       IF(stop.ne.'N'.and.stop.ne.'n')RETURN 1
       PRINT nltblq
       IF(Prompt.eq.'Y'.or.Prompt.eq.'y')CALL prop(2,*20)
      END IF
      IF(irs.ne.2)THEN
       irs=9
       IF(treset.eq.' ')treset=Reset
       IF(treset.eq.'Y'.or.treset.eq.'y')irs=0
       IF(treset.eq.'N'.or.treset.eq.'n')irs=1
       IF(irs.eq.9)THEN
        CALL inerr(treset,irs,*40)
        IF(Nr.eq.10)WRITE(Nerr,1010)'TRESET'
 1010   FORMAT(10X,'THE ENTRY FOR ',A6,' IS NOT DEFINED')
        IF(Nr.eq.5)PRINT 1010,'TRESET '
       END IF
      END IF
      IF(tbldec.lt.0.or.tbldec.gt.5)THEN
       CALL inerr(treset,irs,*40)
       IF(Nr.eq.10)WRITE(Nerr,1020)'TBLDEC',tbldec,0,5
 1020  FORMAT(10X,'THE ENTRY FOR ',A6,' IS ',I6,', WHICH IS NOT ',
     &        'DEFINED',/,10X,'THE ENTRY SHOULD BE BETWEEN ',I6,' AND ',
     &        I6)
       IF(Nr.eq.5)PRINT 1020,'TBLDEC',tbldec,0,5
      END IF
      Ifor=9
      IF(forcst.eq.'Y'.or.forcst.eq.'y')Ifor=1
      IF(forcst.eq.'N'.or.forcst.eq.'n')Ifor=0
      IF(Ifor.eq.9)THEN
       CALL inerr(treset,irs,*40)
       IF(Nr.eq.10)WRITE(Nerr,1010)'FORCST'
       IF(Nr.eq.5)PRINT 1010,'FORCST'
      END IF
      Iovr=9
      IF(overit.eq.'N'.or.overit.eq.'n')Iovr=0
      IF(overit.eq.'Y'.or.overit.eq.'y')Iovr=1
      IF(overit.eq.'A'.or.overit.eq.'a')Iovr=2
      IF(overit.eq.'O'.or.overit.eq.'o')Iovr=3
      IF(Iovr.eq.9)THEN
       CALL inerr(treset,irs,*40)
       IF(Nr.eq.10)WRITE(Nerr,1010)'OVERIT'
       IF(Nr.eq.5)PRINT 1010,'OVERIT'
      END IF
      DO j=1,Numtbl
       IF(Tbl(j).eq.' ')CALL terr(Tbl(j),treset,irs,3,Itbl(j),*30,*40)
       READ(Tbl(j),1030)cpart,ktbl
 1030  FORMAT(A1,I2)
       IF(ktbl.gt.17.or.ktbl.le.0)CALL terr(Tbl(j),treset,irs,1,Itbl(j),
     &    *30,*40)
       ic=ichar(cpart)-ichar('a')
       IF(ic.lt.0)ic=ichar(cpart)-ichar('A')
       IF(ic.le.0.and.ic.gt.5)CALL terr(Tbl(j),treset,irs,1,Itbl(j),*30,
     &                                  *40)
       Itbl(j)=ic*100+ktbl
       IF((Itbl(j).ge.314).or.(Itbl(j).ge.407.and.Itbl(j).le.417).or.
     &    (Itbl(j).ge.503.and.Itbl(j).le.517).or.
     &    (Itbl(j).ge.114..and.Itbl(j).le.116).or.
     &    (Itbl(j).ge.214.and.Itbl(j).le.216).or.
     &    (Itbl(j).eq.112.or.Itbl(j).eq.203.or.Itbl(j).eq.208.or.Itbl(j)
     &    .eq.212.or.Itbl(j).eq.303))
     &    CALL terr(Tbl(j),treset,irs,1,Itbl(j),*30,*40)
       IF(Itbl(j).eq.309.or.Itbl(j).eq.404.or.Itbl(j).eq.502)
     &    CALL terr(Tbl(j),treset,irs,2,Itbl(j),*30,*40)
   30 END DO
   40 IF(Kswt.eq.0)THEN
       IF(Nr.eq.10)THEN
        PRINT *,'  ERROR(S) FOUND IN NLTBLQ NAMELIST ENTRIES'
        PRINT *,'  NO TABLES WILL BE STORED.  X-11.2 RUN CONTINUES'
        Numtbl=0
        RETURN
       END IF
       GO TO 10
      END IF
      IF(tblfmt.eq.'2')THEN
       Itfmt=0
       WRITE(Tf4,1040)tbldec
 1040  FORMAT('(4F12.',I1,',24X,I2.2,A6)')
      ELSE IF(tblfmt.eq.'1')THEN
       Itfmt=0
       WRITE(Tf4,1050)tbldec
 1050  FORMAT('(4F18.',I1,',I2.2,A6)')
      ELSE IF(tblfmt.eq.' ')THEN
       Itfmt=Ifmt
       Tf4=Nf4
      ELSE
       Itfmt=1
       Tf4=tblfmt
      END IF
      itb=1
      Itchar=0
      IF(drive.ne.' ')THEN
       IF(drive(1:1).eq.' ')THEN
        DO WHILE (.true.)
         itb=itb+1
         IF(drive(itb:itb).ne.' ')GO TO 50
        END DO
       END IF
   50  DO WHILE (.true.)
        Itchar=Itchar+1
        it=Itchar+itb-1
        IF(drive(it:it).eq.' ')THEN
         Tblfil(1:(Itchar-1))=drive(itb:(it-1))
         GO TO 60
        END IF
       END DO
      END IF
   60 IF(Itchar.eq.0)Itchar=1
      IF(Iovr.eq.3)THEN
       Tblfil(Itchar:Itchar+5)='X11Q2.'
       Itchar=Itchar+6
       RETURN
      END IF
      Tblfil(Itchar:(Itchar+5))=Ser
      idx=index(Ser,' ')
      IF(idx.eq.0)idx=7
      Itchar=Itchar+idx-1
      Tblfil(Itchar:Itchar)='.'
      Itchar=Itchar+1
      IF(Iovr.eq.0)THEN
       tbf=Tblfil
       i1=1
       DO i=1,Numtbl
        tbf(Itchar:Itchar+2)=Tbl(i)
        itc1=Itchar+1
        IF(tbf(itc1:itc1).eq.' ')tbf(itc1:itc1)='_'
        io=0
        i0=0
   70   INQUIRE(FILE=tbf,EXIST=fxst,ERR=110)
        i0=i0+1
        IF(fxst)CALL novrit(i,io,tbf,*70,*40)
        IF(.not.fxst.and.io.gt.0.and.io.gt.i1)THEN
         i1=i0
         Tblfil=tbf
        END IF
        tbf(Itchar-1:Itchar+4)='.     '
       END DO
       IF(i1.gt.1)THEN
        Itchar=Itchar+2
        Iovr=-1
       END IF
      END IF
      RETURN
   80 CALL inerr(treset,irs,*40)
      WRITE(Nerr,1060)
 1060 FORMAT(10X,'ERROR IN READING NLTBLQ')
      Numtbl=0
      RETURN 1
   90 CALL inerr(treset,irs,*40)
      WRITE(Nerr,1070)
 1070 FORMAT(10X,'END OF FILE ENCOUNTERED IN READING NLTBLQ')
      RETURN 1
  100 PRINT *,' ERROR IN READING NLTBLQ.  RE-ENTER NLTBLQ NAMELIST.'
      GO TO 10
  110 ic2=Itchar+2
      WRITE(f1,1080)(Itchar+2)
 1080 FORMAT('(10X,17HUNABLE TO ACCESS ,A',I2,',18H FOR TABLE STORAGE)')
      IF(ic2.lt.10)f1(27:28)=' A'
      WRITE(Nerr,f1)tbf(1:Itchar+2)
      IF(Nr.eq.10)THEN
       CALL inerr(treset,irs,*40)
       WRITE(Nerr,1090)
 1090  FORMAT(10X,'X-11.2 RUN WILL CONTINUE WITHOUT TABLE STORAGE')
       IF(Ixerr.eq.0)Ixerr=1
       Numtbl=0
       RETURN
      ELSE IF(Nr.eq.5)THEN
       PRINT f1,tbf(1:Itchar+2)
       GO TO 10
      END IF
      END
**==terr.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE terr(T,Treset,Irs,Ix,Ij,*,*)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER Ij,Irs,Irun,Ix,Ixerr,Kswerr,Nerr,Nr
C*** End of declarations inserted by SPAG
      CHARACTER*1 Treset
      CHARACTER*3 T
      COMMON /xerr  / Nr,Nerr,Kswerr,Irun,Ixerr
      IF(Ix.eq.1)THEN
       CALL inerr(Treset,Irs,*10)
       IF(Nr.eq.10)WRITE(Nerr,1010)T
 1010  FORMAT(10X,'TABLE ',A3,' DOES NOT EXIST')
       IF(Nr.eq.5)PRINT 1010,T
       Ij=0
      ELSE IF(Ix.eq.0)THEN
       IF(Nr.eq.10)WRITE(Nerr,1020)T
 1020  FORMAT(10X,'TABLE ',A3,' NOT AVAILABLE ON THIS LEVEL OF ',
     &        'PRINTOUT')
       IF(Nr.eq.5)PRINT 1020,T
      ELSE IF(Ix.eq.2)THEN
       CALL inerr(Treset,Irs,*10)
       IF(Nr.eq.10)WRITE(Nerr,1030)T
 1030  FORMAT(10X,'TABLE ',A3,' NOT AVAILABLE FOR PRINTING')
       IF(Nr.eq.5)PRINT 1030,T
       Ij=0
      ELSE
       CALL inerr(Treset,Irs,*10)
       IF(Nr.eq.10)WRITE(Nerr,1040)
 1040  FORMAT(10X,'TABLE NUMBER NOT ENTERED.  CHECK TBL, NUMTBL ',
     &        'VARIABLES IN NLTBLQ NAMELIST.')
       IF(Nr.eq.5)PRINT 1040
       Ij=0
      END IF
      RETURN 1
   10 RETURN 2
      END
**==trite.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE trite(X,L1,L2,Lyr,Ser,Kpart,Ktbl)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,ic2,Ifor,Io12,Iovr,iret,Irun,Itbl,itc1,Itchar,Itfmt,
     &        Ixerr,iyr,j,j1,j2,k,Kpart,Kswerr,Ktbl
      INTEGER l,L1,L2,l3,l4,Lyr,Nerr,NL2,Nr,Numtbl
      DOUBLE PRECISION X
C*** End of declarations inserted by SPAG
      PARAMETER(NL2=125)
      LOGICAL fxst
      CHARACTER*3 Tbl
      CHARACTER*6 Ser
      CHARACTER*40 Tf4
      CHARACTER*80 Tblfil,f1
      DIMENSION X(NL2)
      COMMON /itcom / Itchar,Itfmt,Itbl(10),Ifor,Numtbl,Iovr,Io12
      COMMON /tcom  / Tf4,Tblfil,Tbl(10)
      COMMON /xerr  / Nr,Nerr,Kswerr,Irun,Ixerr
      k=Kpart*100+Ktbl
      IF(k.eq.300.and.Ifor.eq.0)RETURN
      IF(k.eq.300.and.Ifor.eq.1)k=310
      iret=0
      DO j=1,Numtbl
       IF(k.eq.Itbl(j))THEN
        Itbl(j)=999
        iret=j
       END IF
      END DO
      IF(iret.eq.0)RETURN
      itc1=Itchar+1
      Tblfil(Itchar:Itchar+2)='TBL'
      IF(Iovr.lt.3)THEN
       Tblfil(Itchar:Itchar+2)=Tbl(iret)
       IF(Tblfil(itc1:itc1).eq.' ')Tblfil(itc1:itc1)='_'
      END IF
      fxst=.false.
      IF(Iovr.eq.2.or.Iovr.eq.3)THEN
       INQUIRE(FILE=Tblfil,EXIST=fxst,ERR=10)
       IF(fxst)THEN
        IF(Iovr.eq.2.or.(Iovr.eq.3.and.Io12.eq.1))
     &     OPEN(12,FILE=Tblfil,ACCESS='APPEND',STATUS='OLD',ERR=10)
        IF(Iovr.eq.3.and.Io12.eq.0)fxst=.false.
       END IF
      END IF
      IF(.not.fxst)OPEN(12,FILE=Tblfil,ERR=10)
      IF(Iovr.eq.3.and.Io12.eq.0)Io12=1
      IF(fxst)Itbl(iret)=1999
      IF(L1.ne.1)THEN
       DO j=1,L1-1
        X(j)=0.
       END DO
      END IF
      l3=L2
      l=l3/4
      IF(mod(l3,4).ne.0)l=l+1
      l4=l*4
      IF(l3.ne.l4)THEN
       DO j=l3+1,l4
        X(j)=0.
       END DO
      END IF
      IF(Itfmt.eq.0)THEN
       DO i=1,l
        j1=(i-1)*4+1
        j2=i*4
        iyr=Lyr+i-1
        IF(iyr.ge.2000)THEN
         iyr=iyr-2000
        ELSE
         iyr=iyr-1900
        END IF
        WRITE(12,Tf4)(X(j),j=j1,j2),iyr,Ser
       END DO
      ELSE
       WRITE(12,Tf4)(X(j),j=1,l4)
      END IF
      CLOSE(12)
      RETURN
   10 Itbl(iret)=0
      ic2=Itchar+2
      WRITE(f1,1010)ic2
 1010 FORMAT('(10X,15HUNABLE TO OPEN ,A',I2,',21H, CANNOT STORE TABLE ',
     &       ',A3)')
      IF(ic2.lt.10)f1(25:26)=' A'
      IF(Nr.eq.10)THEN
       WRITE(Nerr,f1)Tblfil(1:ic2),Tbl(iret)
       IF(Ixerr.eq.0)Ixerr=1
      END IF
      IF(Nr.eq.5)PRINT f1,Tblfil(1:ic2),Tbl(iret)
      RETURN
      END
**==novrit.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE novrit(J,Iasc,Tbf,*,*)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,Iasc,ic,Ifor,Io12,Iovr,Irun,Itbl,Itchar,Itfmt,Ixerr,J,
     &        Kswerr,Nerr,Nr,Numtbl
C*** End of declarations inserted by SPAG
      CHARACTER*3 Tbl
      CHARACTER*40 Tf4
      CHARACTER*80 Tblfil,Tbf
      COMMON /itcom / Itchar,Itfmt,Itbl(10),Ifor,Numtbl,Iovr,Io12
      COMMON /tcom  / Tf4,Tblfil,Tbl(10)
      COMMON /xerr  / Nr,Nerr,Kswerr,Irun,Ixerr
      IF(Iasc.eq.0)THEN
       Iasc=47
       DO i=2,-1,-1
        ic=Itchar+i
        Tbf(ic+2:ic+2)=Tbf(ic:ic)
       END DO
       Tbf(ic:ic)='$'
      END IF
      Iasc=Iasc+1
      IF(Iasc.eq.58)Iasc=65
      IF(Iasc.eq.91)Iasc=35
      IF(Iasc.eq.39.and.Nr.eq.10)THEN
       IF(Ixerr.eq.0)Ixerr=1
       WRITE(Nerr,1010)Tbl(J)
 1010  FORMAT(10X,'OVERWRITE OPTION ACTIVATED; CANNOT ASSIGN NEW ',
     &        'FILENAME.',/,10X,'TABLE ',A3,' WILL NOT BE STORED')
      ELSE IF(Iasc.eq.39.and.Nr.eq.5)THEN
       PRINT 1010,Tbl(J)
      END IF
      IF(Iasc.eq.39)RETURN 2
      Tbf(Itchar:Itchar)=char(Iasc)
      RETURN 1
      END
**==vcomp.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE vcomp(Ny,Xmi,Xmorg)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION b,bx,flen,fny,fsum1,fsum2,ftemp,Rv,sum1,sum2,
     &                 temp,trend,vars,Vc,Vi,vo,Vs,w,Xi,xmat
      DOUBLE PRECISION Xmi,Xmorg,Xorg,Xs,Xtc,Xwt
      INTEGER i,Ib,Ie,iflag,j,j1,j2,k,Kp,Kswt,Length,lf,lfda,llda,Lstob,
     &        Lstyr,Lswv,Lyr,Meas,Np
      INTEGER nvar,Ny
C*** End of declarations inserted by SPAG
      DIMENSION xmat(31,31),b(31),bx(31),trend(125),temp(125)
      DIMENSION Xmi(125),Xmorg(125)
      INTEGER Add
      COMMON /extr  / Xs(129),Xi(125),Xwt(125)
      COMMON /inpt  / Xorg(125),Ib,Ie,Meas,Kp,Np
      COMMON /mq2   / Vi,Vc,Vs,Rv,Length,Lstob
      COMMON /mq4   / Xtc(125),Lstyr
      COMMON /tadj  / Add
      COMMON /xk    / Kswt,Lswv,Lyr
C *** INITIALIZE VARIABLES, TAKE LOG OF TREND IF NOT ADDITIVE ADJUSTMENT
      lfda=Ib
      llda=Ie
      IF(Add.ne.1)CALL logar(Xtc,lfda,llda)
      sum1=fsum1(Ny)
      sum2=fsum2(Ny)
      nvar=Lstyr-Lyr+2
      fny=float(Ny)
      flen=float(Length)
C *** SET UP X-PRIME-X MATRIX FOR SPLINE REGRESSION
      xmat(nvar,nvar)=fsum2(Lstob)
      DO i=1,nvar-1
       IF(i.eq.1)THEN
        xmat(i,i)=flen
       ELSE
        ftemp=flen-float(Ny*(i-2)+(Ny-lfda+1))
        xmat(i,i)=sum2+fny*fny*ftemp
        IF(i.eq.2.and.lfda.gt.1)xmat(i,i)=xmat(i,i)-fsum2(lfda-1)
       END IF
       DO j=i+1,nvar
        ftemp=flen-float(Ny*(j-2)+(Ny-lfda+1))
        IF(j.eq.nvar)ftemp=0.0
        IF(i.eq.1)THEN
         xmat(i,j)=sum1+fny*ftemp
         IF(j.eq.nvar)xmat(i,j)=fsum1(Lstob)
         IF(j.eq.2.and.lfda.gt.1)xmat(i,j)=xmat(i,j)-fsum1(lfda-1)
        ELSE
         xmat(i,j)=fny*xmat(1,j)
         IF(i.eq.2.and.lfda.gt.1)xmat(i,j)=xmat(i,j)-fsum1(lfda-1)
        END IF
        xmat(j,i)=xmat(i,j)
       END DO
      END DO
C *** SET UP X-PRIME-Y MATRIX FOR SPLINE REGRESSION
      j2=llda
      lf=mod(lfda,Ny)
      IF(lf.eq.0)lf=Ny
      DO j=1,nvar
       b(j)=0.0
       j1=lfda
       IF(j.gt.2)j1=((lfda/Ny)+j-2)*Ny+1
       DO k=j1,j2
        w=1.0
        IF(j.eq.1)trend(k)=0.0
        IF(j.gt.1)w=float(k-j1+1)
        IF(j.eq.2)w=w+(lf-1)
        IF(w.gt.Ny)w=Ny
        b(j)=b(j)+w*Xtc(k)
       END DO
      END DO
C *** USE GAUSSIAN ELIMINATION AND PIVOTING SUBROUTINE TO ESTIMATE
C     COEFFICIENTS TO BE USED IN GENERATING THE SPLINE TREND
      CALL lassol(nvar,xmat,b,31,bx,iflag)
      IF(iflag.eq.2)THEN
       PRINT *,'MATRIX APPEARS SINGULAR, WILL NOT PROCEED FURTHER'
       STOP 10
      END IF
C *** USING SPLINE REGRESSION COEFFICIENTS, COMPUTE SPLINE TREND
C     AND STORE IN VARIABLE TREND
      DO j=1,nvar
       j1=lfda
       IF(j.gt.2)j1=((lfda/Ny)+j-2)*Ny+1
       j2=llda
       DO k=j1,j2
        w=1.0
        IF(j.gt.1)w=float(k-j1+1)
        IF(w.gt.Ny)w=Ny
        trend(k)=trend(k)+w*bx(j)
       END DO
      END DO
C *** CALCULATE RELATIVE CONTRIBUTION OF EACH COMPONENT
C     TO STATIONARY PORTION OF THE VARIANCE, USING SPLINE
C     TREND AS DETRENDER
      IF(Add.ne.1)THEN
       CALL antilg(Xtc,lfda,llda)
       CALL antilg(trend,lfda,llda)
      END IF
      Vs=vars(Xs,lfda,llda,1)
      Vi=vars(Xmi,lfda,llda,1)
      CALL divsub(Xtc,trend,temp,lfda,llda)
      Vc=vars(temp,lfda,llda,0)
      CALL divsub(Xmorg,trend,temp,lfda,llda)
      vo=vars(temp,lfda,llda,0)/100.0
      Vs=Vs/vo
      Vi=Vi/vo
      Vc=Vc/vo
      Rv=Vc+Vs+Vi
      RETURN
      END
**==fsum1.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      DOUBLE PRECISION FUNCTION fsum1(N)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER N
C*** End of declarations inserted by SPAG
      fsum1=dble(N*(N+1))/2.0
      RETURN
      END
**==fsum2.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      DOUBLE PRECISION FUNCTION fsum2(N)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER N
C*** End of declarations inserted by SPAG
      fsum2=dble(N*(N+1)*(2*N+1))/6.0
      RETURN
      END
**==lassol.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE lassol(N,A,B,M,X,Iflag)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,ib,idxpiv,Iflag,ip1,j,k,kp1,M,N,nm1,np1
C*** End of declarations inserted by SPAG
C
C  LASSOL SOLVES A SYSTEM OF N LINEAR EQUATIONS IN N UNKOWNS, AX = B,
C  USING GAUSSIAN ELIMINATION WITH PARITAL PIVOTING AND ROW
C  EQUILIBRATION.
C
C  TAKEN FROM NUMERICAL COMPUTING: AN INTRODUCTION BY SHAMPINE AND
C  ALLEN.  ADAPTED BY BRIAN MONSELL, 12-6-88.
C
      DOUBLE PRECISION A(M,M),B(N),X(N),ab(31,32),rowmax,big,quot,sum,
     &                 tempb,tempi,scale,one
      EQUIVALENCE(big,quot,sum,scale),(rowmax,tempb,tempi)
      DATA one/1.0/
      np1=N+1
      nm1=N-1
C
C  FORM THE N BY (N+1) MATRIX AB, THE FIRST N COLUMNS OF WHICH ARE
C  A AND THE REMAINING COLUMN B.  CALCULATE SCALE FACTORS AND SCALE
C  AB.
C
      DO i=1,N
       rowmax=0.
       DO j=1,N
        rowmax=dmax1(rowmax,dabs(A(i,j)))
       END DO
       scale=one/rowmax
       DO j=1,N
        ab(i,j)=A(i,j)*scale
       END DO
       ab(i,np1)=B(i)*scale
      END DO
C
C  BEGIN BASIC ELIMINATION LOOP.  ROWS OF AB ARE PHYSICALLY
C  INTERCHANGED IN ORDER TO BRING ELEMENT OF LARGEST MAGNITUDE
C  INTO PIVOTAL POSITION
C
      DO k=1,nm1
       big=0.0
       DO i=k,N
        tempb=dabs(ab(i,k))
        IF(big.lt.tempb)THEN
         big=tempb
         idxpiv=i
        END IF
       END DO
       IF(big.eq.0)GO TO 10
       IF(idxpiv.ne.k)THEN
C
C  PIVOT IS IN ROW IDXPIV.  INTERCHANGE ROUW IDXPIV WITH ROW K.
C
        DO i=k,np1
         tempi=ab(k,i)
         ab(k,i)=ab(idxpiv,i)
         ab(idxpiv,i)=tempi
        END DO
       END IF
       kp1=k+1
C
C  ELIMINATE X(K) FROM EQUATIONS K+1,K+2,...,K+N.
C
       DO i=kp1,N
        quot=ab(i,k)/ab(k,k)
        DO j=kp1,np1
         ab(i,j)=ab(i,j)-(quot*ab(k,j))
        END DO
       END DO
      END DO
C
C  BEGIN CALCULATION OF SOLUTION X USING BACK SUBSTITUTION.
C
      IF(ab(N,N).ne.0.0)THEN
       X(N)=ab(N,np1)/ab(N,N)
       DO ib=2,N
        i=np1-ib
        ip1=i+1
        sum=0.0
        DO j=ip1,N
         sum=sum+ab(i,j)*X(j)
        END DO
        X(i)=(ab(i,np1)-sum)/ab(i,i)
       END DO
C
C  SET IFLAG = 1 FOR NORMAL RETURN
C            = 2 IF MATRIX APPEARS SINGULAR TO THE CODE.
C
       Iflag=1
       RETURN
      END IF
   10 Iflag=2
      RETURN
      END
**==vars.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      DOUBLE PRECISION FUNCTION vars(X,I,J,Iopt)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER I,Iopt,J
      DOUBLE PRECISION varian,varlog,X
C*** End of declarations inserted by SPAG
      INTEGER Add
      COMMON /tadj  / Add
      DIMENSION X(*)
      IF(Add.eq.1)vars=varian(X,I,J,Iopt)
      IF(Add.eq.0)vars=varlog(X,I,J,Iopt)
      RETURN
      END
**==varlog.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      DOUBLE PRECISION FUNCTION varlog(X,I,J,Iopt)
C --- THIS FUNCTION COMPUTES THE LOG VARIANCE OF X. IF IOPT IS EQUAL
C --- TO 1 THE LOG MEAN IS ASSUMED TO BE 1.
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER I,Iopt,J,k
      DOUBLE PRECISION tmp,tmp2,X
C*** End of declarations inserted by SPAG
      DIMENSION X(*)
      tmp=0.0
      IF(Iopt.eq.0)THEN
       DO k=I,J
        tmp=tmp+dlog(X(k))
       END DO
       tmp=tmp/float(J-I+1)
      END IF
      varlog=0.0
      DO k=I,J
       IF(Iopt.eq.0)tmp2=dlog(X(k))-tmp
       IF(Iopt.eq.1)tmp2=dlog(X(k)/100.0)
       varlog=varlog+tmp2*tmp2
      END DO
      RETURN
      END
**==varian.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      DOUBLE PRECISION FUNCTION varian(X,I,J,Iopt)
C --- THIS FUNCTION COMPUTES THE VARIANCE OF X. IF IOPT = 0 COMPUTE THE
C --- MEAN , IF IOPT = 1 THE MEAN IS ASSUMED TO BE ZERO, AND IF IOPT = 2
C --- THE MEAN IS ASSUMED TO BE ONE HUNDRED.
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION ave,X
      INTEGER I,Iopt,J,k
C*** End of declarations inserted by SPAG
      DIMENSION X(*)
      IF(Iopt.eq.2)ave=100.0
      IF(Iopt.eq.1)ave=0.0
      IF(Iopt.eq.0)THEN
       DO k=I,J
        ave=ave+X(k)
       END DO
       ave=ave/(J-I+1)
      END IF
      varian=0.0
      DO k=I,J
       varian=varian+(X(k)-ave)*(X(k)-ave)
      END DO
      RETURN
      END
**==logar.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE logar(X,I,J)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER I,J,k
      DOUBLE PRECISION tmp,X
C*** End of declarations inserted by SPAG
      DIMENSION X(*)
      DO k=I,J
       tmp=X(k)
       X(k)=dlog(tmp)
      END DO
      RETURN
      ENTRY antilg(X,I,J)
      DO k=I,J
       tmp=X(k)
       X(k)=dexp(tmp)
      END DO
      RETURN
      END
**==v15.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE v15(X,Trend,Iopt,Jopt,Lfd1,Length,Vx,Iy15)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER Iopt,Iy15,Jopt,l,l2,l3,Length,Lfd1
      DOUBLE PRECISION temp,Trend,vars,Vx,X
C*** End of declarations inserted by SPAG
      DIMENSION X(*),Trend(*),temp(64)
      DO l=1,Iy15
       l3=Length+l+Lfd1-Iy15-1
       IF(Jopt.eq.1)THEN
        l2=l+Lfd1-1
        temp(l)=X(l3)/Trend(l2)
       ELSE
        temp(l)=X(l3)
       END IF
      END DO
      Vx=vars(temp,1,Iy15,Iopt)
      RETURN
      END
**==f3cal.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE f3cal(Isq,Smic,Ny)
C --- THIS SUBROUTINE CALCULATES THE QUALITY CONTROL STATISTICS IN TABLE
C --- F3.
C
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION Adri,ave,ave1,ave2,ave3,ave4,count,ct,ct1,ct2,
     &                 diff,F,fn,P,Qual,Ratic,Ratis,rmcd,Rv,sd
      DOUBLE PRECISION sdev,Smic,Sti,t1,t2,temp,Test1,Test2,twt,Vc,Vi,
     &                 Vs,wt,Xorg,Xs,Xwt
      INTEGER i,i1,Ib,Ie,ifail,Iqfail,j,k,Kfail,klda,kny,Kp,Kpart,
     &        Ktcopt,l,Length,Lstob,Lterm,Mcd,Meas
      INTEGER mtype,n,Nn,Np,Nterm,Ny,Nyrs
C*** End of declarations inserted by SPAG
      DOUBLE PRECISION Isq(4),Nm
      INTEGER Add
      DIMENSION wt(11),Smic(4),temp(372)
      COMMON /cfilt / Ktcopt,Lterm(4),Nterm
      COMMON /extr  / Xs(129),Sti(125),Xwt(125)
      COMMON /inpt  / Xorg(125),Ib,Ie,Meas,Kp,Np
      COMMON /mor2  / Kpart,Mcd
      COMMON /mq1   / Adri,Ratis,Ratic
      COMMON /mq2   / Vi,Vc,Vs,Rv,Length,Lstob
      COMMON /mq3   / Qual,Nm(11),Nn,Kfail,Nyrs
      COMMON /tadj  / Add
      COMMON /tests / F(2),P(2),Test1,Test2,Iqfail
      DATA wt/10.,11.,10.,8.,11.,10.,18.,7.,7.,4.,4./
      mtype=Lterm(1)
      DO i=2,Ny
       IF(mtype.ne.Lterm(i))mtype=6
      END DO
      ifail=0
      kny=Ny/4
      Nm(1)=Isq(kny)/10
      kny=12/Ny
      Nm(2)=Vi/10
      Nm(3)=(Ratic*kny-1.0)/2.0
      fn=Ie-Ib+1
      t1=dabs(3.0*(fn-1.0)/Adri-2.0*(fn-1.0))
      t2=dsqrt(1.6*fn-2.9)*2.577
      Nm(4)=t1/t2
      IF(Mcd.eq.1)THEN
       rmcd=1+(Smic(1)-1.0)/(Smic(1)-Smic(2))
       IF(rmcd.lt.0.5)rmcd=0.5
       IF(rmcd.gt.1.0)rmcd=1.0
      ELSE
       rmcd=Mcd+(Smic(Mcd)-1.0)/(Smic(Mcd-1)-Smic(Mcd))
       IF(Mcd.eq.Ny.and.Smic(Ny-1).le.Smic(Ny))rmcd=kny*15.5
      END IF
      Nm(5)=(rmcd*kny-0.5)/5.0
      Nm(6)=dabs(Ratis-4.0)/2.5
      Nm(7)=dsqrt((Test1+Test2)/2.0)
      Nyrs=(Ie-Ib+1)/Ny
      Nn=7
      IF(mtype.ne.5.and.Nyrs.ge.6)THEN
       Nn=11
       n=2-Add
       DO i=Ib,Ie
        temp(i)=Xs(i)/100
       END DO
       sd=sdev(temp,Ib,Ie,1,n)
       ave=1-dble(Add)
       DO i=Ib,Ie
        temp(i)=(temp(i)-ave)/sd
       END DO
       ct1=0.0
       ct2=0.0
       count=0.0
       ave1=0.0
       ave2=0.0
       ave3=0.0
       ave4=0.0
       klda=Ib+Ny-1
       DO i=Ib,klda
        ct=0.0
        i1=i+Ny
        DO j=i1,Ie,Ny
         ct=ct+1.0
         count=count+1.0
         diff=dabs(temp(j)-temp(j-Ny))
         k=j
         ave1=ave1+diff
        END DO
        ave2=ave2+dabs(temp(k)-temp(i))/ct
        k=k-2*Ny
        j=k-3*Ny
        IF(j.ge.Ib)THEN
         ave3=ave3+dabs(temp(k)-temp(j))/3.0
         ct1=ct1+1.0
         j=j+Ny
         DO l=j,k,Ny
          i1=l-Ny
          IF(i1.ge.Ib)THEN
           ct2=ct2+1.0
           diff=dabs(temp(l)-temp(i1))
           ave4=ave4+diff
          END IF
         END DO
        END IF
       END DO
       ave1=ave1/count
       ave2=ave2/Ny
       IF(ct1.gt.0.0)THEN
        ave3=ave3/ct1
        ave4=ave4/ct2
       END IF
       Nm(8)=10.0*ave1
       Nm(9)=10.0*ave2
       Nm(10)=10.0*ave4
       Nm(11)=10.0*ave3
      END IF
      Qual=0.0
      twt=0.0
      DO i=1,11
       IF(i.le.Nn)THEN
        IF(Nm(i).lt.0.0)Nm(i)=0.0
        IF(Nm(i).gt.3.0)Nm(i)=3.0
        IF(Nm(i).ge.1.0)ifail=ifail+1
        IF(i.ne.6.or.(mtype.eq.3.and.i.eq.6))Qual=Qual+Nm(i)*wt(i)
       ELSE IF(i.eq.10)THEN
        Qual=Qual+Nm(1)*wt(i)
       ELSE IF(i.eq.11)THEN
        Qual=Qual+Nm(2)*wt(i)
       ELSE
        Qual=Qual+Nm(7)*wt(i)
       END IF
       IF(i.ne.6.or.(mtype.eq.3.and.i.eq.6))twt=twt+wt(i)
      END DO
      Qual=Qual/twt
      Kfail=ifail
      RETURN
      END
**==sdev.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      DOUBLE PRECISION FUNCTION sdev(X,I,J,K,Iopt)
C --- THIS FUNCTION CALCULATES THE STANDARD DEVIATION OF X. IF IOPT = 0
C --- THE MEAN OF X IS COMPUTED, IF IOPT = 1 THE MEAN IS ASSUMED TO BE
C --- ZERO, AND IF IOPT = 2 THE MEAN IS ASSUMED TO BE ONE.
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION ave,fn,X
      INTEGER I,Iopt,J,K,l
C*** End of declarations inserted by SPAG
      DOUBLE PRECISION totals
      DIMENSION X(*)
      fn=dble((J-I)/K+1)
      ave=0.0
      IF(Iopt.eq.2)ave=1.0
      IF(Iopt.eq.0)ave=totals(X,I,J,K,1)
      sdev=0.0
      DO l=I,J,K
       sdev=sdev+(X(l)-ave)*(X(l)-ave)
      END DO
      sdev=sqrt(sdev/fn)
      RETURN
      END
**==f3gen.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE f3gen(Ny,Mq)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER Igrz,k,Kfail,Kpage,l,Myr,Ndec,Nn,Nr,Nw,Ny,Nyear,Nyrs
      DOUBLE PRECISION Qu,Qual
C*** End of declarations inserted by SPAG
      CHARACTER*6 span(4),Serno
      CHARACTER*7 Mq
      CHARACTER*80 Title
C --- THIS SUBROUTINE PRINTS THE F3 TABLE.
      COMMON /mq3   / Qual,Qu(11),Nn,Kfail,Nyrs
      COMMON /ngrz  / Igrz,Myr,Nw
      COMMON /rite  / Nr,Nyear,Title,Kpage,Ndec,Serno
      DATA span/'THREE ','MONTHS','ONE QU','ARTER '/
      k=1
      IF(Ny.eq.4)k=3
      l=k+1
      Kpage=Kpage+1
      WRITE(Nw,1010)Title,Kpage,Serno
 1010 FORMAT('1',17X,A80,1X,'PAGE ',I2,', SERIES ',A6)
      WRITE(Nw,1020)
 1020 FORMAT(//,' F 3. MONITORING AND QUALITY ASSESSMENT STATISTICS')
      WRITE(Nw,1030)
 1030 FORMAT('0',6X,
     &'ALL THE MEASURES BELOW ARE IN THE RANGE FROM 0 TO 3 WITH AN ACCEP
     &TANCE REGION FROM 0 TO 1.')
      WRITE(Nw,1040)span(k),span(l),Qu(1)
 1040 FORMAT('0',3X,
     &       '1. THE RELATIVE CONTRIBUTION OF THE IRREGULAR OVER ',2A6,
     &       ' SPAN (FROM TABLE F 2.B).',T99,'M1  = ',F5.3,//)
      WRITE(Nw,1050)Qu(2)
 1050 FORMAT('0',3X,
     &'2. THE RELATIVE CONTRIBUTION OF THE IRREGULAR COMPONENT TO THE ST
     &ATIONARY PORTION OF',T99,'M2  = ',F5.3,/,7X,
     &'THE VARIANCE (FROM TABLE F 2.F).',/)
      WRITE(Nw,1060)Mq,Mq,Qu(3),Mq,Mq
 1060 FORMAT('0',3X,'3. THE AMOUNT OF ',A7,' TO ',A7,
     &       ' CHANGE IN THE IRREGULAR COMPONENT AS COMPARED TO THE',
     &       T99,'M3  = ',F5.3,/,7X,'AMOUNT OF ',A7,' TO ',A7,
     &       ' CHANGE IN THE TREND-CYCLE (FROM TABLE F2.G).',/)
      WRITE(Nw,1070)Qu(4)
 1070 FORMAT('0',3X,
     &'4. THE AMOUNT OF AUTOCORRELATION IN THE IRREGULAR AS DESCRIBED BY
     & THE AVERAGE DURATION',T99,'M4  = ',F5.3,/,7X,
     &'OF RUN (TABLE F 2.C).',/)
      WRITE(Nw,1080)Mq,Qu(5)
 1080 FORMAT('0',3X,'5. THE NUMBER OF ',A7,
     &   'S IT TAKES THE CHANGE IN THETREND-CYCLE TO SURPASS THE AMOUNT'
     &   ,T99,'M5  = ',F5.3,/,6X,
     &   ' OF CHANGE IN THE IRREGULAR (FROM TABLE F 2.D).',/)
      WRITE(Nw,1090)Qu(6)
 1090 FORMAT('0',3X,
     &'6. THE AMOUNT OF YEAR TO YEAR CHANGE IN THE IRREGULAR AS COMPARED
     & TO THE AMOUNT OF YEAR',T99,'M6  = ',F5.3,/,7X,
     &'TO YEAR CHANGE IN THE SEASONAL (FROM TABLE F 2.G).',/)
      WRITE(Nw,1100)Qu(7)
 1100 FORMAT('0',3X,
     &'7. THE AMOUNT OF MOVING SEASONALITY PRESENT RELATIVE TO THE AMOUN
     &T OF STABLE',T99,'M7  = ',F5.3,/,7X,
     &'SEASONALITY (FROM TABLE F 2.H).',/)
      IF(Nn.gt.7)THEN
       WRITE(Nw,1110)Qu(8)
 1110  FORMAT('0',3X,
     &'8. THE SIZE OF THE FLUCTUATIONS IN THE SEASONAL COMPONENT THROUGH
     &OUT THE WHOLE SERIES.',T99,'M8  = ',F5.3,//)
       WRITE(Nw,1120)Qu(9)
 1120  FORMAT('0',3X,
     &'9. THE AVERAGE LINEAR MOVEMENT IN THE SEASONAL COMPONENT THROUGHO
     &UT THE WHOLE SERIES.',T99,'M9  = ',F5.3,//)
       WRITE(Nw,1130)Qu(10)
 1130  FORMAT('0',2X,'10. SAME AS 8, CALCULATED FOR RECENT YEARS ',
     &        'ONLY.',T99,'M10 = ',F5.3,//)
       WRITE(Nw,1140)Qu(11)
 1140  FORMAT('0',2X,'11. SAME AS 9, CALCULATED FOR RECENT YEARS ONLY .'
     &        ,T99,'M11 = ',F5.3)
      END IF
      IF(Qual.lt.0.8)WRITE(Nw,1150)Qual
      IF(Qual.lt.1.0.and.Qual.ge.0.8)WRITE(Nw,1160)Qual
      IF(Qual.lt.1.2.and.Qual.ge.1.0)WRITE(Nw,1170)Qual
      IF(Qual.gt.1.2)WRITE(Nw,1180)Qual
 1150 FORMAT('0',/,42X,'*** ACCEPTED *** AT THE LEVEL  ',F5.2)
 1160 FORMAT('0',/,28X,'*** CONDITIONALLY ACCEPTED *** AT THE LEVEL  ',
     &       F5.2)
 1170 FORMAT('0',/,28X,'*** CONDITIONALLY REJECTED *** AT THE LEVEL ',
     &       F5.2)
 1180 FORMAT('0',/,42X,'*** REJECTED *** AT THE LEVEL ',F5.2)
      IF(Kfail.gt.0)THEN
       WRITE(Nw,1190)Kfail
 1190  FORMAT(/,42X,'*** CHECK THE ',I2,' ABOVE MEASURES WHICH ',
     &        'FAILED.')
      END IF
      RETURN
      END
**==totals.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      DOUBLE PRECISION FUNCTION totals(X,I,J,K,Iopt)
C --- THIS FUNCTION CALCULATES TOTALS AND AVERAGES. A TOTAL IS RETURNED
C --- IF IOPT IS ZERO AND AN AVERAGE IF IOPT IS ONE.
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      DOUBLE PRECISION fn,X
      INTEGER I,Iopt,J,K,l
C*** End of declarations inserted by SPAG
      DIMENSION X(*)
      totals=0.0
      DO l=I,J,K
       totals=totals+X(l)
      END DO
      IF(Iopt.ne.0)THEN
       fn=(J-I)/K+1
       totals=totals/fn
      END IF
      RETURN
      END
**==ssapq.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
 
      SUBROUTINE ssapq
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL c,dmax,S,Sa,yy
      INTEGER i1,i2,Ibreak,Iq,iquar,Irang,Iyr,Kpage,l,l0,Len,Len2,
     &        length,Lquar,Lyr,Muladd,n2,n48,Ncol,Nlen
      INTEGER Np,Ns1,Nsea,Nyear,nyears,nyearz
C*** End of declarations inserted by SPAG
C  *****
C  *****  Main Subroutine for the Sliding Spans Analysis.  This subroutine
C  *****  initializes variables and calls subroutines which generate
C  *****  breakdown tables, a listing of each observation in each span,
C  *****  and a range analysis of the seasonal factors
C  *****
      DIMENSION yy(60,4),c(60,4),dmax(60,4)
      INTEGER Y,Qter,Outdec
      CHARACTER*1 fc(3)
      CHARACTER*6 Lser
      CHARACTER*35 F2
      CHARACTER*50 F1
      CHARACTER*80 Title
      COMMON /com4  / Y(8),Qter(8)
      COMMON /com4b / Iyr,Iq,Nsea,Lyr,Lquar
      COMMON /iopt  / Ibreak,Irang
      COMMON /mform / F1,F2
      COMMON /ns4   / Nlen,Ncol,Ns1,Len,Len2
      COMMON /rite  / Np,Nyear,Title,Kpage,Outdec,Lser
      COMMON /sser  / S(60,4),Sa(60,4)
      COMMON /tadj  / Muladd
      DATA fc/'2','3','4'/
      Nsea=4
      CALL header
      CALL ssftst(Ncol,Np)
      l=Iq+(Iyr*Nsea)+Nlen-1
      Lquar=mod(l,Nsea)
      IF(Lquar.eq.0)THEN
       Lquar=Nsea
       Lyr=(l/Nsea)-1
      ELSE
       Lyr=l/Nsea
      END IF
      nyears=Lyr-Iyr+1
      Len2=Len-Nsea+1
      Ns1=Ncol+1
      nyearz=nyears+Ncol-3
      iquar=Iq+Nsea-1
      Kpage=Kpage+1
      WRITE(Np,1010)Title,Kpage,Lser
 1010 FORMAT('1',17X,A80,1X,'PAGE ',I2,', SERIES ',A6)
      IF(Muladd.eq.0)WRITE(Np,1020)Lser
      IF(Muladd.eq.1)WRITE(Np,1030)Lser
 1020 FORMAT(//,' S 1.  RANGE ANALYSIS OF SEASONAL FACTORS FOR ',A6,//,
     &       9X,'MEANS OF SEASONAL FACTORS FOR EACH QUARTER',/,9X,
     &       '(MOVEMENTS WITHIN A QUARTER SHOULD BE SMALL)')
 1030 FORMAT(//,' S 1.  RANGE ANALYSIS OF IMPLIED ADJUSTMENT FACTORS ',
     &       'FOR ',A6,//,9X,
     &       'MEANS OF IMPLIED ADJUSTMENT FACTORS FOR EACH',' QUARTER',
     &       /,9X,'(MOVEMENTS WITHIN A QUARTER SHOULD BE SMALL)')
      CALL range(Np,S)
      IF(Muladd.eq.0)CALL btab(S,1,0,iquar,nyearz,dmax)
      IF(Muladd.eq.1)CALL btab(Sa,2,0,iquar,nyearz,dmax)
      CALL xchng(Sa,c,yy,nyears)
      iquar=iquar+1
      CALL btab(c,3,1,iquar,nyearz,dmax)
      iquar=iquar+Nsea-1
      CALL btab(yy,4,2,iquar,nyearz,dmax)
      IF(Ibreak.eq.0)THEN
       DO l=1,Ncol
        i1=l
        i2=l+Ncol
        Y(i1)=Iyr+l-1
        Y(i2)=Lyr+l-1
        Qter(i1)=Iq
        Qter(i2)=Lquar
       END DO
       n2=Ncol*2
       length=Len-Iq+1
       n48=(length/48)+1
       IF(mod(length,48).eq.0)n48=n48-1
       l0=Ncol-1
       F1(17:17)=fc(l0)
       F2(6:6)=fc(l0)
       IF(Muladd.eq.0)CALL qlist(S,1,0,dmax,n2,n48)
       IF(Muladd.eq.1)CALL qlist(Sa,2,0,dmax,n2,n48)
       CALL qlist(c,3,1,dmax,n2,n48)
       CALL qlist(yy,4,2,dmax,n2,n48)
      END IF
      RETURN
      END
**==btab.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE btab(X,Icol,Ic2,Iqu,Nyearz,Dmax)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Dmax,X
      INTEGER Ibreak,Ic2,Icol,Iqu,Irang,Kpage,Len,Len2,Ncol,Nlen,Np,Ns1,
     &        Nyear,Nyearz
C*** End of declarations inserted by SPAG
C  *****
C  *****  This subroutine is called for each seasonally adjusted estimate
C  *****  (seasonal factors, quarter-to-quarter changes, etc.).  Calls other
C  *****  subroutines which determine which observations are flagged as
C  *****  extreme, prints out breakdown table, calculate and print out
C  *****  sliding spans histogram, and calculate and print out hinge
C  *****  values.  Prints heading on top of page.
C  *****
      DIMENSION Dmax(60,4),X(60,4)
      INTEGER Outdec
      CHARACTER*6 Lser
      CHARACTER*80 Title
      COMMON /iopt  / Ibreak,Irang
      COMMON /ns4   / Nlen,Ncol,Ns1,Len,Len2
      COMMON /rite  / Np,Nyear,Title,Kpage,Outdec,Lser
      Kpage=Kpage+1
      IF(Irang.eq.0)WRITE(Np,1010)Title,Kpage,Lser
 1010 FORMAT('1',17X,A80,1X,'PAGE ',I2,', SERIES ',A6)
      CALL mflag(X,Icol,Ic2,Iqu,Dmax)
      CALL calca(icol,ic2)
      IF(Irang.eq.1)RETURN
      CALL btrit(Nyearz,Icol,Ic2)
      CALL hinge(Np,Dmax,Len2,Icol,Ic2)
      RETURN
      END
**==mflag.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE mflag(X,Iopt,Iop2,Kq,Dmax)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Aquar,Ayr,Cut,Dmax,Ten,X
      INTEGER i,Idquar,Idy,Iop2,Iopt,Iq,Iscf,Isign,Itf,Itot,Iturn,Iyr,j,
     &        k,Kount,Kq,Len,Len2
      INTEGER Lquar,Lyr,Ncol,Nlen,Nquar,Ns1,Nsea,Ntot,Nyr
C*** End of declarations inserted by SPAG
C  *****
C  *****  This subroutine controls the flagging of individual quarters
C  *****  for a given set of seasonally adjusted estimates.  In addition,
C  *****  the values for the sliding spans histogram (KOUNT) are generated
C  *****  here.
C  *****
      DIMENSION X(60,4),Dmax(60,4)
      CHARACTER*1 Csign,Cturn
      CHARACTER*4 Per
      COMMON /ave   / Aquar(4),Ayr(15)
      COMMON /com1  / Isign,Iscf,Iturn,Itf,Itot
      COMMON /com1c / Csign(60,4),Per(60,4),Cturn(60,4)
      COMMON /com3  / Nquar(4),Nyr(15),Ntot
      COMMON /com4b / Iyr,Iq,Nsea,Lyr,Lquar
      COMMON /com5  / Kount(4),Cut(4,4)
      COMMON /ctot  / Idquar(4),Idy(15)
      COMMON /kten  / Ten
      COMMON /ns4   / Nlen,Ncol,Ns1,Len,Len2
      Iscf=0
      Isign=0
      Itf=0
      Itot=0
      Iturn=0
      Ntot=0
      DO j=1,15
       Ayr(j)=0.
       Idy(j)=0
       Nyr(j)=0
      END DO
      DO j=1,Nsea
       Aquar(j)=0.
       Idquar(j)=0
       Nquar(j)=0
      END DO
      DO k=1,4
       Kount(k)=0
      END DO
      DO i=1,Len
       Csign(i,Iopt)=' '
       Cturn(i,Iopt)=' '
       Per(i,Iopt)='    '
       IF(i.le.Kq.or.i.ge.Len2)THEN
        Dmax(i,Iopt)=Ten
       ELSE
        CALL rplus(X,i,Iopt,Iop2,Dmax(i,Iopt))
       END IF
      END DO
      DO k=1,3
       Kount(k)=Kount(k)-Kount(k+1)
      END DO
      RETURN
      END
**==rplus.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE rplus(X,I,Nopt,Nop2,R)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Aquar,Ayr,Cut,R,Ten,X,xbase,xj,xmn,xmx
      INTEGER I,ibase,Ic,Icq,Icyr,Idquar,Idy,ij,Iq,Iscf,Isign,Itf,Itot,
     &        ittot,Iturn,iy,Iyr,j,jay,jay2,jbase1
      INTEGER jay3,jbase,k,Kount,Len,Len2,Lquar,Lyr,Ncol,Nlen,Nop2,Nopt,
     &        Nquar,Ns1,Nsea,Ntot,Nyr
C*** End of declarations inserted by SPAG
C  *****
C  *****  This subroutine computes the maximum percentage difference R
C  *****  for observation I over NCOL spans, and determines if this quarter
C  *****  should be flagged.  Quarterly (AQUAR) and yearly (YMON) averages of
C  *****  the maximum percentage difference are incremented, as well as
C  *****  the sliding spans histogram values (KOUNT) and the number of
C  *****  quarters flagged for each quarter (NQUAR) and year (NYR).
C  *****  Determines the total number of observations tested (ITOT), the
C  *****  number of obserations tested for each calendar quarter (IDQUAR)
C  *****  and year (IDY).  Determines if the estimates for observation I
C  *****  have undergone a change of direction.
C  *****
      DIMENSION X(60,4)
      INTEGER qmod
      LOGICAL chnge
      CHARACTER*1 Ch,Csign,Cturn
      CHARACTER*4 Per,temp
      CHARACTER*18 Chper
      CHARACTER*24 Cname
      COMMON /ave   / Aquar(4),Ayr(15)
      COMMON /cee   / Ic,Icyr,Icq
      COMMON /chame / Cname(4),Ch(4),Chper(2)
      COMMON /com1  / Isign,Iscf,Iturn,Itf,Itot
      COMMON /com1c / Csign(60,4),Per(60,4),Cturn(60,4)
      COMMON /com3  / Nquar(4),Nyr(15),Ntot
      COMMON /com5  / Kount(4),Cut(4,4)
      COMMON /com4b / Iyr,Iq,Nsea,Lyr,Lquar
      COMMON /ctot  / Idquar(4),Idy(15)
      COMMON /kten  / Ten
      COMMON /ns4   / Nlen,Ncol,Ns1,Len,Len2
      qmod=mod(I,Nsea)
      IF(qmod.eq.0)qmod=Nsea
      iy=(I-1)/Nsea
      j=1
      DO WHILE (X(I,j).eq.Ten)
       j=j+1
      END DO
      jay=j
      xmx=X(I,j)
      xmn=X(I,j)
      DO WHILE (.true.)
       j=j+1
       IF(j.gt.Ncol.or.X(I,j).eq.Ten)THEN
        temp='    '
        jay2=j-1
        R=xmx-xmn
        IF(Nop2.eq.0)R=(R/abs(xmn))*100.0
        IF(I.lt.Ic)THEN
         Per(I,Nopt)=' -- '
         RETURN
        END IF
        IF(R.ge.Cut(Nopt,1))THEN
         DO k=1,4
          IF(R.ge.Cut(Nopt,k))THEN
           temp(k:k)=Ch(Nopt)
           Kount(k)=Kount(k)+1
          END IF
         END DO
         Per(I,Nopt)=temp
         Nquar(qmod)=Nquar(qmod)+1
         Nyr(iy)=Nyr(iy)+1
         Ntot=Ntot+1
        END IF
        Aquar(qmod)=Aquar(qmod)+R
        Ayr(iy)=Ayr(iy)+R
        Idquar(qmod)=Idquar(qmod)+1
        Idy(iy)=Idy(iy)+1
        Itot=Itot+1
        Csign(I,Nopt)=' '
        IF(Nopt.eq.2)GO TO 10
        j=jay
        jbase=0
        xbase=100.
        IF(Nop2.gt.0)xbase=0.
        DO WHILE (.true.)
         IF(X(I,j).eq.xbase)THEN
          ibase=0
         ELSE IF(X(I,j).lt.xbase)THEN
          ibase=-1
         ELSE
          ibase=1
         END IF
         IF(jbase.eq.0)THEN
          jbase=ibase
         ELSE
          ij=ibase+jbase
          IF(ij.eq.0)THEN
           Isign=Isign+1
           Csign(I,Nopt)='*'
           RETURN
          END IF
         END IF
         j=j+1
         IF(j.gt.Ncol.or.X(I,j).eq.Ten)GO TO 10
        END DO
       ELSE
        IF(xmx.lt.X(I,j))xmx=X(I,j)
        IF(xmn.gt.X(I,j))xmn=X(I,j)
       END IF
      END DO
   10 IF(Ncol.eq.2)RETURN
      jay3=jay2-jay+1
      Cturn(I,Nopt)=' '
      IF(jay3.lt.3)RETURN
      ibase=0
      jbase=0
      chnge=.false.
      ittot=ittot+1
      DO j=jay+1,jay2
c-----------------------------------------------------------------------
c     Check to see if observation I has already been flagged as a 
c     turning point.
c-----------------------------------------------------------------------
       IF(Cturn(I,Nopt).eq.' ')THEN
c-----------------------------------------------------------------------
c     Compute change between adjoining spans.  Check to see if the
c     change is greater than one percent
c-----------------------------------------------------------------------
        xj=X(I,j)-X(I,j-1)
        IF(Nop2.eq.0)xj=(xj/X(I,j-1))*100.0
        IF(xj.gt.1.0)chnge=.true.
c-----------------------------------------------------------------------
c     set ibase to be 1 if the change of direction is up, -1 if it
c     is down
c-----------------------------------------------------------------------
        IF(xj.lt.0)ibase=-1
        IF(xj.gt.0)ibase=1
c-----------------------------------------------------------------------
c     If 4 sliding spans, save direction of first change
c-----------------------------------------------------------------------
        IF(jay3.eq.4.and.j.eq.2)jbase1=ibase
        IF(jbase.eq.0)THEN
c-----------------------------------------------------------------------
c     If jbase hasn't been set yet, set jbase to be the direction of the
c     current change
c-----------------------------------------------------------------------
         jbase=ibase
        ELSE IF(ibase.ne.jbase)THEN
c-----------------------------------------------------------------------
c     If direction has changed, test to see if the current change is 
c     greater than one percent.
c-----------------------------------------------------------------------
         IF(abs(xj).gt.1)THEN
c-----------------------------------------------------------------------
c     If so, declare a turning point
c-----------------------------------------------------------------------
          Iturn=Iturn+1
          Cturn(I,Nopt)='^'
         ELSE
c-----------------------------------------------------------------------
c     If not, update the direction indicator variable jbase
c-----------------------------------------------------------------------
          jbase=ibase
         END IF
        ELSE IF(jay3.eq.4.and.j.eq.4)THEN
c-----------------------------------------------------------------------
c     If there is no change in direction for the last span, determine
c     the direction of the change from the second to the fourth span.
c-----------------------------------------------------------------------
         xj=X(I,j)-X(I,j-2)
         IF(xj.lt.0)ibase=-1
         IF(xj.gt.0)ibase=1
c-----------------------------------------------------------------------
c     If the direction of this change is different than that of the 
c     first two spans, test to see if this change is greater than one
c     percent.
c-----------------------------------------------------------------------
         IF(ibase.ne.jbase1)THEN
          IF(Nop2.eq.0)xj=(xj/X(I,j-2))*100.0
          IF(abs(xj).gt.1.and.chnge)THEN
c-----------------------------------------------------------------------
c     If so, declare a turning point, if at least one of the individal
c     changes is greater than one percent.
c-----------------------------------------------------------------------
           Iturn=Iturn+1
           Cturn(I,Nopt)='^'
          END IF
         END IF
        END IF
       END IF
      END DO
      RETURN
      END
**==calca.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE calca(I,I2)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Aquar,Ayr
      INTEGER I,I2,Idquar,Idy,Iq,Iscf,Isign,Itf,Itot,Iturn,Iyr,j,j1,jm,
     &        jy,Len,Len2,Lquar,Lyr,Ncol
      INTEGER Nlen,Ns1,Nsea
C*** End of declarations inserted by SPAG
C  *****
C  *****  Computes quarterly (AQUAR) and yearly (AYR) averages of the
C  *****  maximum percentage difference.  Determines number of observations
C  *****  with changes in direction which were flagged as extreme (ISCF).
C  *****
      INTEGER t2
      CHARACTER*4 Per
      CHARACTER*1 Csign,Cturn
      COMMON /ave   / Aquar(4),Ayr(15)
      COMMON /com1  / Isign,Iscf,Iturn,Itf,Itot
      COMMON /com1c / Csign(60,4),Per(60,4),Cturn(60,4)
      COMMON /com4b / Iyr,Iq,Nsea,Lyr,Lquar
      COMMON /ctot  / Idquar(4),Idy(15)
      COMMON /ns4   / Nlen,Ncol,Ns1,Len,Len2
      t2=Iq+((Iyr+Ncol-2)*Nsea)+Nlen-1
      jm=mod(t2,Nsea)
      jy=t2/Nsea
      IF(jm.eq.0)jy=jy-1
      jy=jy-Iyr
      DO j=1,jy
       IF(Idy(j).gt.0)Ayr(j)=Ayr(j)/Idy(j)
      END DO
      DO j=1,Nsea
       Aquar(j)=Aquar(j)/Idquar(j)
      END DO
      j1=Nsea+1
      IF(I2.eq.2)j1=j1+Nsea
      DO j=j1,Len2
       IF(Per(j,I).ne.'    '.and.Csign(j,I).eq.'*'.and.I.ne.3)
     &    Iscf=Iscf+1
       IF(Per(j,I).ne.'    '.and.Cturn(j,I).eq.'^')Itf=Itf+1
      END DO
      RETURN
      END
**==btrit.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE btrit(Nyearz,Nopt,No2)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Aquar,Ayr,Cut,fper
      INTEGER i,Ic,Icq,Icyr,ij,Iq,Iscf,Isign,Itf,Itot,Iturn,iy,Iyr,j,
     &        Kount,Kpage,Lquar,Lyr,n2,No2
      INTEGER Nopt,Np,Nquar,Nsea,Ntot,Nyear,Nyearz,Nyr
C*** End of declarations inserted by SPAG
C  *****
C  *****  Prints out summary breakdown tables for total, quarters, years.
C  *****  Calculates percentage of quarters flagged (FPER), percentage of
C  *****  quarters with change in direction (SPER), percentage of quarters
C  *****  with change of direction that were flagged (SFPER).
C  *****
      INTEGER Outdec
      CHARACTER*1 Ch
      CHARACTER*4 c
      CHARACTER*6 Lser
      CHARACTER*18 Chper
      CHARACTER*24 Cname
      CHARACTER*80 Title
      COMMON /ave   / Aquar(4),Ayr(15)
      COMMON /cee   / Ic,Icyr,Icq
      COMMON /chame / Cname(4),Ch(4),Chper(2)
      COMMON /com1  / Isign,Iscf,Iturn,Itf,Itot
      COMMON /com3  / Nquar(4),Nyr(15),Ntot
      COMMON /com4b / Iyr,Iq,Nsea,Lyr,Lquar
      COMMON /com5  / Kount(4),Cut(4,4)
      COMMON /rite  / Np,Nyear,Title,Kpage,Outdec,Lser
      fper=(real(Ntot)/real(Itot))*100.
      n2=Nopt+1
      IF(No2.eq.0)WRITE(Np,1010)n2,Cut(Nopt,1),Cname(Nopt),Lser
      IF(No2.gt.0)WRITE(Np,1020)n2,Cut(Nopt,1),Chper(No2),Cname(Nopt),
     &                          Lser
      WRITE(Np,1030)Ntot,Itot,fper
      IF(Nopt.eq.1)WRITE(Np,1040)
      IF(Nopt.eq.3)WRITE(Np,1050)
      IF(Nopt.eq.4)WRITE(Np,1060)
      WRITE(Np,1070)
      DO j=1,Nsea
       WRITE(Np,1080)j,Nquar(j),Aquar(j)
      END DO
      WRITE(Np,1070)
      ij=Icyr-Iyr+1
      IF((No2.eq.2.or.(No2.eq.1.and.Iq.eq.Nsea)).and.Icyr.eq.Iyr)ij=ij+1
      DO j=ij,Nyearz
       iy=Iyr+j
       WRITE(Np,1090)iy,Nyr(j),Ayr(j)
      END DO
      WRITE(Np,1100)
      CALL prst(Nopt,0)
C  *****
C  *****  Prints out sliding spans histogram.
C  *****
      c='    '
      IF(No2.eq.0)WRITE(Np,1110)n2,Cname(Nopt)
      IF(No2.gt.0)WRITE(Np,1120)n2,Chper(No2),Cname(Nopt)
      DO i=1,3
       c(i:i)=Ch(Nopt)
       WRITE(Np,1130)c,Cut(Nopt,i),Cut(Nopt,i+1),Kount(i)
      END DO
      c(4:4)=Ch(Nopt)
      WRITE(Np,1140)c,Cut(Nopt,4),Kount(4)
 1010 FORMAT(' S ',I1,'.A  BREAKDOWN OF ',F3.1,' PERCENT OR MORE ',
     &       'DIFFERENCES IN ',A24,' FOR ',A6,/)
 1020 FORMAT(' S ',I1,'.A  BREAKDOWN OF ',F3.1,' PERCENT OR MORE ',
     &       'DIFFERENCES IN ',A18,' ',A24,' FOR ',A6,/)
 1030 FORMAT(2X,'TOTAL : ',I3,' OUT OF ',I3,' (',F5.1,' %)')
 1040 FORMAT(2X,'(USUALLY, 15% IS TOO HIGH, AND 25% IS MUCH TOO HIGH)')
 1050 FORMAT(2X,'(USUALLY, 35% IS TOO HIGH, AND 40% IS MUCH TOO HIGH)')
 1060 FORMAT(2X,'(USUALLY, 10% IS TOO HIGH)')
 1070 FORMAT(/)
 1080 FORMAT(10X,'QUARTER ',I1,' : ',I3,5X,'(AMPD = ',F5.1,')')
 1090 FORMAT(10X,I4,5X,' : ',I3,5X,'(AMPD = ',F5.1,')')
 1100 FORMAT(/,'  AMPD = AVERAGE MAXIMUM PERCENTAGE DIFFERENCE')
 1110 FORMAT(/,' S ',I1,'.B  HISTOGRAM OF MAXIMUM PERCENTAGE ',
     &       'DIFFERENCE IN ',A24,' FOR FLAGGED QUARTERS',/)
 1120 FORMAT(/,' S ',I1,'.B  HISTOGRAM OF MAXIMUM PERCENTAGE ',
     &       'DIFFERENCE IN ',A18,' ',A24,' FOR FLAGGED QUARTERS',/)
 1130 FORMAT('  ',A4,' : GREATER THAN OR EQUAL TO ',F4.1,'% BUT LESS ',
     &       'THAN ',F4.1,'% :',1X,I3)
 1140 FORMAT('  ',A4,' : GREATER THAN OR EQUAL TO ',F4.1,'%',T61,':',1X,
     &       I3)
      RETURN
      END
**==hinge.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE hinge(Np,Dmax,Len2,Nopt,No2)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Dmax,temp,Ten,ts,ts65,ts85,xo
      INTEGER i,Ibreak,Ic,Icq,Icyr,Irang,it,l2,l65,l85,Len2,n,n1,n2,n3,
     &        No2,Nopt,Np
C*** End of declarations inserted by SPAG
C  *****
C  *****  Calculates and prints out hinge values for the maximum
C  *****  percentage difference (TS).  Includes the 65th (TS65) and
C  *****  85th (TS85) percentiles.
C  *****
      DIMENSION Dmax(60,4),ts(5),xo(60),temp(60)
      CHARACTER*1 Ch
      CHARACTER*18 Chper
      CHARACTER*24 Lab
      COMMON /cee   / Ic,Icyr,Icq
      COMMON /chame / Lab(4),Ch(4),Chper(2)
      COMMON /iopt  / Ibreak,Irang
      COMMON /kten  / Ten
      l2=Len2-1
      it=0
      DO i=Ic,l2
       IF(Dmax(i,Nopt).eq.Ten)THEN
        it=it+1
       ELSE
        n=i-it-Ic+1
        xo(n)=Dmax(i,Nopt)
       END IF
      END DO
      CALL ssort(xo,temp,n,1)
      ts(1)=xo(1)
      ts(5)=xo(n)
      IF(mod(n,2).eq.1)THEN
       ts(3)=xo((n+1)/2)
      ELSE
       ts(3)=(xo(n/2)+xo((n/2)+1))/2
      END IF
      n2=(n+1)/2
      IF(mod(n2,2).eq.1)THEN
       n1=(n2+1)/2
       n3=n-n1+1
       ts(2)=xo(n1)
       ts(4)=xo(n3)
      ELSE
       n1=n2/2
       n3=n1+1
       ts(2)=(xo(n1)+xo(n3))/2
       n1=n-n1+1
       n3=n-n3+1
       ts(4)=(xo(n1)+xo(n3))/2
      END IF
      l85=n-int(real(n)*0.15+0.5)
      ts85=xo(l85)
      l65=n-int(real(n)*0.35+0.5)
      ts65=xo(l65)
      n2=Nopt+1
      IF(Irang.eq.1)THEN
       IF(No2.eq.0)WRITE(Np,1010)Lab(Nopt)
       IF(No2.gt.0)WRITE(Np,1020)Chper(No2),Lab(Nopt)
      ELSE
       IF(No2.eq.0)WRITE(Np,1030)n2,Lab(Nopt)
       IF(No2.gt.0)WRITE(Np,1040)n2,Chper(No2),Lab(Nopt)
      END IF
 1010 FORMAT(/,5X,'STATISTICS FOR MAXIMUM PERCENTAGE DIFFERENCE OF ',
     &       'THE ',A24,/)
 1020 FORMAT(/,5X,'STATISTICS FOR MAXIMUM PERCENTAGE DIFFERENCE OF ',
     &       'THE ',A14,1X,A24,/)
 1030 FORMAT(/,' S ',I1,'.C  STATISTICS FOR MAXIMUM PERCENTAGE ',
     &       'DIFFERENCE OF THE ',A24,/)
 1040 FORMAT(/,' S ',I1,'.C  STATISTICS FOR MAXIMUM PERCENTAGE ',
     &       'DIFFERENCE OF THE ',A14,1X,A24,/)
      WRITE(Np,1050)(ts(i),i=1,3),ts65,ts(4),ts85,ts(5)
 1050 FORMAT(5X,'MINIMUM',T22,': ',F10.3,/,5X,'25TH PERCENTILE',T22,
     &       ': ',F10.3,/,5X,'MEDIAN',T22,': ',F10.3,/,3X,
     &       '->65TH PERCENTILE',T22,': ',F10.3,'<-',/,5X,
     &       '75TH PERCENTILE',T22,': ',F10.3,/,3X,'->85TH PERCENTILE',
     &       T22,': ',F10.3,'<-',/,5X,'MAXIMUM',T22,': ',F10.3)
      RETURN
      END
**==ssort.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE ssort(X,Y,N,Kflag)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,ij,il,iu,j,k,Kflag,kk,l,m,N,nn
      REAL r,t,tt,tty,ty,X,Y
C*** End of declarations inserted by SPAG
C***BEGIN PROLOGUE  SSORT
C***DATE WRITTEN   761101   (YYMMDD)
C***REVISION DATE  820801   (YYMMDD)
C***CATEGORY NO.  N6A2B1
C***KEYWORDS  QUICKSORT,SINGLETON QUICKSORT,SORT,SORTING
C***AUTHOR  JONES, R. E., (SNLA)
C           WISNIEWSKI, J. A., (SNLA)
C***PURPOSE  SSORT SORTS ARRAY X AND OPTIONALLY MAKES THE SAME
C            INTERCHANGES IN ARRAY Y.  THE ARRAY X MAY BE SORTED IN
C            INCREASING ORDER OR DECREASING ORDER.  A SLIGHTLY MODIFIED
C            QUICKSORT ALGORITHM IS USED.
C***END PROLOGUE  SSORT
      DIMENSION X(N),Y(N),il(21),iu(21)
C***FIRST EXECUTABLE STATEMENT  SSORT
      nn=N
      kk=iabs(Kflag)
C
C ALTER ARRAY X TO GET DECREASING ORDER IF NEEDED
C
      IF(Kflag.lt.1)THEN
       DO i=1,nn
        X(i)=-X(i)
       END DO
      END IF
      GO TO(10,80),kk
C
C SORT X ONLY
C
   10 m=1
      i=1
      j=nn
      r=.375
   20 IF(i.eq.j)GO TO 50
      IF(r.gt..5898437)THEN
       r=r-.21875
      ELSE
       r=r+3.90625E-2
      END IF
   30 k=i
C                                  SELECT A CENTRAL ELEMENT OF THE
C                                  ARRAY AND SAVE IT IN LOCATION T
      ij=i+ifix(float(j-i)*r)
      t=X(ij)
C                                  IF FIRST ELEMENT OF ARRAY IS GREATER
C                                  THAN T, INTERCHANGE WITH T
      IF(X(i).gt.t)THEN
       X(ij)=X(i)
       X(i)=t
       t=X(ij)
      END IF
      l=j
C                                  IF LAST ELEMENT OF ARRAY IS LESS THAN
C                                  T, INTERCHANGE WITH T
      IF(X(j).lt.t)THEN
       X(ij)=X(j)
       X(j)=t
       t=X(ij)
C                                  IF FIRST ELEMENT OF ARRAY IS GREATER
C                                  THAN T, INTERCHANGE WITH T
       IF(X(i).gt.t)THEN
        X(ij)=X(i)
        X(i)=t
        t=X(ij)
       END IF
      END IF
      DO WHILE (.true.)
C                                  FIND AN ELEMENT IN THE SECOND HALF OF
C                                  THE ARRAY WHICH IS SMALLER THAN T
       l=l-1
       IF(X(l).le.t)THEN
        DO WHILE (.true.)
C                                  FIND AN ELEMENT IN THE FIRST HALF OF
C                                  THE ARRAY WHICH IS GREATER THAN T
         k=k+1
         IF(X(k).ge.t)THEN
C                                  INTERCHANGE THESE ELEMENTS
          IF(k.le.l)THEN
           tt=X(l)
           X(l)=X(k)
           X(k)=tt
           GO TO 40
          ELSE
C                                  SAVE UPPER AND LOWER SUBSCRIPTS OF
C                                  THE ARRAY YET TO BE SORTED
           IF(l-i.le.j-k)THEN
            il(m)=k
            iu(m)=j
            j=l
            m=m+1
           ELSE
            il(m)=i
            iu(m)=l
            i=k
            m=m+1
           END IF
           GO TO 60
          END IF
         END IF
        END DO
       END IF
   40 END DO
C                                  BEGIN AGAIN ON ANOTHER PORTION OF
C                                  THE UNSORTED ARRAY
   50 m=m-1
      IF(m.eq.0)GO TO 150
      i=il(m)
      j=iu(m)
   60 IF(j-i.ge.1)GO TO 30
      IF(i.eq.1)GO TO 20
      i=i-1
      DO WHILE (.true.)
       i=i+1
       IF(i.eq.j)GO TO 50
       t=X(i+1)
       IF(X(i).gt.t)THEN
        k=i
        DO WHILE (.true.)
         X(k+1)=X(k)
         k=k-1
         IF(t.ge.X(k))THEN
          X(k+1)=t
          GO TO 70
         END IF
        END DO
       END IF
   70 END DO
C
C SORT X AND CARRY Y ALONG
C
   80 m=1
      i=1
      j=nn
      r=.375
   90 IF(i.eq.j)GO TO 120
      IF(r.gt..5898437)THEN
       r=r-.21875
      ELSE
       r=r+3.90625E-2
      END IF
  100 k=i
C                                  SELECT A CENTRAL ELEMENT OF THE
C                                  ARRAY AND SAVE IT IN LOCATION T
      ij=i+ifix(float(j-i)*r)
      t=X(ij)
      ty=Y(ij)
C                                  IF FIRST ELEMENT OF ARRAY IS GREATER
C                                  THAN T, INTERCHANGE WITH T
      IF(X(i).gt.t)THEN
       X(ij)=X(i)
       X(i)=t
       t=X(ij)
       Y(ij)=Y(i)
       Y(i)=ty
       ty=Y(ij)
      END IF
      l=j
C                                  IF LAST ELEMENT OF ARRAY IS LESS THAN
C                                  T, INTERCHANGE WITH T
      IF(X(j).lt.t)THEN
       X(ij)=X(j)
       X(j)=t
       t=X(ij)
       Y(ij)=Y(j)
       Y(j)=ty
       ty=Y(ij)
C                                  IF FIRST ELEMENT OF ARRAY IS GREATER
C                                  THAN T, INTERCHANGE WITH T
       IF(X(i).gt.t)THEN
        X(ij)=X(i)
        X(i)=t
        t=X(ij)
        Y(ij)=Y(i)
        Y(i)=ty
        ty=Y(ij)
       END IF
      END IF
      DO WHILE (.true.)
C                                  FIND AN ELEMENT IN THE SECOND HALF OF
C                                  THE ARRAY WHICH IS SMALLER THAN T
       l=l-1
       IF(X(l).le.t)THEN
        DO WHILE (.true.)
C                                  FIND AN ELEMENT IN THE FIRST HALF OF
C                                  THE ARRAY WHICH IS GREATER THAN T
         k=k+1
         IF(X(k).ge.t)THEN
C                                  INTERCHANGE THESE ELEMENTS
          IF(k.le.l)THEN
           tt=X(l)
           X(l)=X(k)
           X(k)=tt
           tty=Y(l)
           Y(l)=Y(k)
           Y(k)=tty
           GO TO 110
          ELSE
C                                  SAVE UPPER AND LOWER SUBSCRIPTS OF
C                                  THE ARRAY YET TO BE SORTED
           IF(l-i.le.j-k)THEN
            il(m)=k
            iu(m)=j
            j=l
            m=m+1
           ELSE
            il(m)=i
            iu(m)=l
            i=k
            m=m+1
           END IF
           GO TO 130
          END IF
         END IF
        END DO
       END IF
  110 END DO
C                                  BEGIN AGAIN ON ANOTHER PORTION OF
C                                  THE UNSORTED ARRAY
  120 m=m-1
      IF(m.eq.0)GO TO 150
      i=il(m)
      j=iu(m)
  130 IF(j-i.ge.1)GO TO 100
      IF(i.eq.1)GO TO 90
      i=i-1
      DO WHILE (.true.)
       i=i+1
       IF(i.eq.j)GO TO 120
       t=X(i+1)
       ty=Y(i+1)
       IF(X(i).gt.t)THEN
        k=i
        DO WHILE (.true.)
         X(k+1)=X(k)
         Y(k+1)=Y(k)
         k=k-1
         IF(t.ge.X(k))THEN
          X(k+1)=t
          Y(k+1)=ty
          GO TO 140
         END IF
        END DO
       END IF
  140 END DO
C
C CLEAN UP
C
  150 IF(Kflag.ge.1)RETURN
      DO i=1,nn
       X(i)=-X(i)
      END DO
      RETURN
      END
**==qlist.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE qlist(X,Nopt,Nop2,Dmax,N2,N48)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL Dmax,X
      INTEGER i,Ibreak,Ic,Icq,Icyr,Iq,Irang,iy,Iyr,Kpage,l,l0,l1,l2,Len,
     &        Len2,Lquar,Lyr,m,N2
      INTEGER N48,Ncol,Nlen,Nop2,Nopt,Np,Ns1,Nsea,Nyear
C*** End of declarations inserted by SPAG
C  *****
C  *****  Prints out each observation in sliding spans analysis, with
C  *****  date, estimates (example, seasonal factors) for each span,
C  *****  maximum percentage difference (DMAX), and an indication of
C  *****  whether the observation was flagged as an extreme (PER)
C  *****  or changed direction (CSIGN).
C  *****
      DIMENSION X(60,4),Dmax(60,4)
      INTEGER Y,Qter,Outdec
      CHARACTER*1 dash(3),Csign,Ch,cpart(5),Cturn
      CHARACTER*4 Per
      CHARACTER*7 f(2)
      CHARACTER*6 Lser
      CHARACTER*18 Chper
      CHARACTER*24 Cname
      CHARACTER*50 F1
      CHARACTER*35 F2
      CHARACTER*80 Title
      COMMON /cee   / Ic,Icyr,Icq
      COMMON /chame / Cname(4),Ch(4),Chper(2)
      COMMON /com1c / Csign(60,4),Per(60,4),Cturn(60,4)
      COMMON /com4  / Y(8),Qter(8)
      COMMON /com4b / Iyr,Iq,Nsea,Lyr,Lquar
      COMMON /iopt  / Ibreak,Irang
      COMMON /mform / F1,F2
      COMMON /ns4   / Nlen,Ncol,Ns1,Len,Len2
      COMMON /rite  / Np,Nyear,Title,Kpage,Outdec,Lser
      DATA dash/'-','/',' '/
      DATA f/'MAXIMUM','% DIFF.'/
      DATA cpart/'A','B','C','D','E'/
      iy=Iyr
      m=Iq-1
      DO l0=1,N48
       l1=(l0-1)*48+Iq
       l2=l0*48+Iq-1
       IF(l0.eq.N48)l2=Len
       Kpage=Kpage+1
       WRITE(Np,1010)Title,Kpage,Lser
 1010  FORMAT('1',17X,A80,1X,'PAGE ',I2,', SERIES ',A6)
       IF(Nop2.eq.0)WRITE(Np,1030)cpart(Nopt),Cname(Nopt),Lser
       IF(Nop2.gt.0)WRITE(Np,1040)cpart(Nopt),Chper(Nop2),Cname(Nopt),
     &                            Lser
       WRITE(Np,1020)
       WRITE(Np,F2)(Qter(i),dash(2),Y(i),dash(1),i=1,Ncol),f(1)
       WRITE(Np,F2)(Qter(i),dash(2),Y(i),dash(3),i=Ncol+1,N2),f(2)
       DO l=l1,l2
        m=m+1
        IF(m.gt.Nsea)THEN
         m=1
         iy=iy+1
        END IF
        WRITE(Np,F1)m,dash(1),iy,(X(l,i),i=1,Ncol),Csign(l,Nopt),
     &              Cturn(l,Nopt),Dmax(l,Nopt),Per(l,Nopt)
       END DO
      END DO
      IF(Irang.eq.1)THEN
       CALL hinge(Np,Dmax,Len2,Nopt,Nop2)
       CALL prst(Nopt,Irang)
      END IF
      WRITE(Np,1020)
      IF(Nopt.ne.3)WRITE(Np,1050)
      IF(Ncol.gt.2)WRITE(Np,1060)
      IF(Iq.ne.Ic)WRITE(Np,1070)
 1020 FORMAT(' ')
 1030 FORMAT(' S 7.',A1,'  SLIDING SPANS ANALYSIS OF ',A26,' FOR ',A6)
 1040 FORMAT(' S 7.',A1,'  SLIDING SPANS ANALYSIS OF ',A18,' ',A26,
     &       ' FOR ',A6)
 1050 FORMAT('  * - SIGN CHANGE FOR OBSERVATION OVER ALL SPANS')
 1060 FORMAT('  ^ - TURNING POINT FOR OBSERVATION OVER ALL SPANS')
 1070 FORMAT(' -- - OBSERVATION NOT INCLUDED IN SLIDING SPANS ',
     &       'COMPARISONS')
      RETURN
      END
**==range.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE range(Np,X)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,Ibreak,Ic,Icq,Icyr,Iq,Irang,Iyr,j,j4,k,Kountr,l,Len,
     &        Len2,Lquar,Lyr,m1,Muladd,Ncol
      INTEGER Nlen,Np,Ns1,Nsea
      REAL Ten,X,xavg,xmn,xmnx,xmpd,xmx,xran,xranr
C*** End of declarations inserted by SPAG
C  *****
C  *****  Performs range analysis of the seasonal factors.  Calculates
C  *****  quarterly means of the seasonal factors for each span and for
C  *****  all spans (XAVG), as well as the maximum percentage difference
C  *****  for the quarterly means (XMPD).  Computes values of the range
C  *****  (XRAN) and the range ratio (XRANR) for quarterly means and
C  *****  seasonal factors for each span and for all spans.  Prints out
C  *****  results.
C  *****
      DIMENSION X(60,4),xavg(5,4),xmnx(5,2),xran(5,2),xranr(5,2)
      DIMENSION xmpd(4)
      CHARACTER*1 f(3)
      CHARACTER*3 Xcm,xc(2)
      CHARACTER*6 clab(6)
      CHARACTER*57 f1
      CHARACTER*27 f2
      COMMON /cee   / Ic,Icyr,Icq
      COMMON /charx / Xcm(5,4)
      COMMON /com4b / Iyr,Iq,Nsea,Lyr,Lquar
      COMMON /iopt  / Ibreak,Irang
      COMMON /kcom  / Kountr(4,4)
      COMMON /kten  / Ten
      COMMON /ns4   / Nlen,Ncol,Ns1,Len,Len2
      COMMON /tadj  / Muladd
      DATA xc/' SF','IAF'/
      DATA(f(i),i=1,3)/'2','3','4'/
      DATA(clab(i),i=1,6)/'SPAN 1','SPAN 2','SPAN 3','SPAN 4','  MPD ',
     &     ' TOTAL'/
      DO j=1,4
       DO i=1,5
        IF(i.le.4)Kountr(i,j)=0
        IF(j.le.2)THEN
         xran(i,j)=0.
         xranr(i,j)=0.
         xmnx(i,j)=100.
        END IF
        xavg(i,j)=0.
       END DO
      END DO
      DO i=1,Ncol
       DO j=1,Len
        IF(X(j,i).ne.Ten.and.j.ge.Ic)THEN
         j4=mod(j,Nsea)
         IF(j4.eq.0)j4=Nsea
         xavg(i,j4)=xavg(i,j4)+X(j,i)
         xavg(Ns1,j4)=xavg(Ns1,j4)+X(j,i)
         IF(xmnx(i,1).gt.X(j,i))xmnx(i,1)=X(j,i)
         IF(xmnx(i,2).lt.X(j,i))xmnx(i,2)=X(j,i)
         Kountr(i,j4)=Kountr(i,j4)+1
        END IF
       END DO
       CALL compb(xavg,xmnx,xran,xranr,i)
      END DO
      CALL compb(xavg,xmnx,xran,xranr,Ns1)
      DO i=1,Nsea
       xmx=xavg(1,i)
       xmn=xavg(1,i)
       DO j=2,Ncol
        IF(xmx.lt.xavg(j,i))xmx=xavg(j,i)
        IF(xmn.gt.xavg(j,i))xmn=xavg(j,i)
       END DO
       xmpd(i)=((xmx-xmn)/xmn)*100.
      END DO
      f1='(2X,8HQUARTER ,I1,3X, (F10.2,2X,A3),F10.2,5X,F10.2,2X,A3)'
      f2='(//,8X, (9X,A6),2(9X,A6),/)'
      f1(22:22)=f(Ncol-1)
      f2(8:8)=f(Ncol-1)
      WRITE(Np,f2)(clab(j),j=1,Ncol),clab(5),clab(6)
      DO k=1,Nsea
       WRITE(Np,f1)k,(xavg(j,k),Xcm(j,k),j=1,Ncol),xmpd(k),xavg(Ns1,k),
     &             Xcm(Ns1,k)
      END DO
      WRITE(Np,1010)Ncol
 1010 FORMAT(/,10X,'MPD = MAXIMUM PERCENT DIFFERENCE = (MAX - MIN) ',
     &       '/ MIN',/,10X,'TOTAL = AVERAGE TAKEN OVER ALL ',I1,
     &       ' SPANS',//,10X,'SUMMARY OF RANGE MEASURES')
      m1=Muladd+1
      WRITE(Np,1020)xc(m1),xc(m1),xc(m1),xc(m1)
 1020 FORMAT(/,20X,'RANGE',9X,'R-R',14X,'MIN',12X,'MAX',11X,'RANGE',9X,
     &       'R-R',/,20X,'MEANS',8X,'MEANS',13X,A3,12X,A3,12X,A3,10X,A3,
     &       /)
      DO k=1,Ncol
       WRITE(Np,1030)k,xran(k,1),xranr(k,1),(xmnx(k,l),l=1,2),xran(k,2),
     &               xranr(k,2)
 1030  FORMAT(4X,'SPAN ',I1,5X,F10.2,5X,F10.4,3(5X,F10.2),5X,F10.4,/)
      END DO
      WRITE(Np,1040)xran(Ns1,1),xranr(Ns1,1),(xmnx(Ns1,l),l=1,2),
     &              xran(Ns1,2),xranr(Ns1,2)
 1040 FORMAT(2X,'ALL SPANS',4X,F10.2,5X,F10.4,3(5X,F10.2),5X,F10.4,/)
      WRITE(Np,1050)
      IF(m1.eq.1)WRITE(Np,1060)
      IF(m1.eq.2)WRITE(Np,1070)
 1050 FORMAT(10X,'R-R = RANGE RATIO = MAX / MIN')
 1060 FORMAT(10X,'SF = SEASONAL FACTORS')
 1070 FORMAT(10X,'IAF = IMPLIED ADJUSTMENT FACTORS')
      Irang=0
      IF(xran(Ns1,2).lt.10)THEN
       Irang=1
       Ibreak=0
       WRITE(Np,1080)
 1080  FORMAT(/,5X,'WARNING : RANGE OF SEASONAL FACTORS IS TOO LOW ',
     &        'FOR SUMMARY SLIDING SPANS MEASURES TO BE RELIABLE.',/,
     &        15X,'SUMMARY SLIDING SPANS STATISTICS NOT PRINTED OUT')
      END IF
      RETURN
      END
**==compb.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE compb(Xavg,Xmnx,Xran,Xranr,I)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER I,i2,ij,Iq,Iyr,j,k,Kountr,Len,Len2,Lquar,Lyr,Ncol,Nlen,
     &        Ns1,Nsea
      REAL Ten,Xavg,xmax,xmin,Xmnx,Xran,Xranr
C*** End of declarations inserted by SPAG
C  *****
C  *****  Computes quarterly means of seasonal factors for each span
C  *****  and all spans.  Determines which quarters have the maximum
C  *****  and minimum averages for each span and all spans (XMNX).
C  *****
      CHARACTER*3 Xcm
      DIMENSION Xavg(5,4),Xmnx(5,2),Xran(5,2),Xranr(5,2)
      COMMON /charx / Xcm(5,4)
      COMMON /com4b / Iyr,Iq,Nsea,Lyr,Lquar
      COMMON /kcom  / Kountr(4,4)
      COMMON /ns4   / Nlen,Ncol,Ns1,Len,Len2
      COMMON /kten  / Ten
      Xran(I,2)=Xmnx(I,2)-Xmnx(I,1)
      IF(I.ne.Ns1)THEN
       IF(Xmnx(Ns1,1).gt.Xmnx(I,1))Xmnx(Ns1,1)=Xmnx(I,1)
       IF(Xmnx(Ns1,2).lt.Xmnx(I,2))Xmnx(Ns1,2)=Xmnx(I,2)
       DO j=1,Nsea
        Xavg(I,j)=Xavg(I,j)/Kountr(I,j)
       END DO
      ELSE
       DO j=1,Nsea
        k=0
        DO i2=1,Ncol
         k=Kountr(i2,j)+k
        END DO
        Xavg(I,j)=Xavg(I,j)/k
       END DO
      END IF
      Xranr(I,2)=Ten
      IF(Xmnx(I,1).ne.0.0)Xranr(I,2)=Xmnx(I,2)/Xmnx(I,1)
      xmax=100.
      xmin=100.
      ij=1
      DO j=ij,Nsea
       IF(xmin.gt.Xavg(I,j))xmin=Xavg(I,j)
       IF(xmax.lt.Xavg(I,j))xmax=Xavg(I,j)
      END DO
      Xran(I,1)=xmax-xmin
      Xranr(I,1)=Ten
      IF(xmin.ne.0.0)Xranr(I,1)=xmax/xmin
      DO j=1,Nsea
       Xcm(I,j)='   '
       IF(xmax.eq.Xavg(I,j))THEN
        Xcm(I,j)='MAX'
       ELSE IF(xmin.eq.Xavg(I,j))THEN
        Xcm(I,j)='MIN'
       END IF
      END DO
      RETURN
      END
**==xchng.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE xchng(X,C,Yy,Nyears)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      REAL C,Ten,X,Yy
      INTEGER i,i1,i2,Iq,iten,Iyr,iyy,iyy2,Len,Len2,Lquar,Lyr,Ncol,Nlen,
     &        Ns1,Nsea,Nyears
C*** End of declarations inserted by SPAG
C  *****
C  *****  Calculate quarter-to-quarter and year-to-year changes of the
C  *****  seasonally adjusted data
C  *****
      DIMENSION X(60,4),Yy(60,4),C(60,4)
      COMMON /com4b / Iyr,Iq,Nsea,Lyr,Lquar
      COMMON /kten  / Ten
      COMMON /ns4   / Nlen,Ncol,Ns1,Len,Len2
      DO i=1,Ncol
       i1=((i-1)*Nsea)+1
       i2=(Nyears+i-1)*Nsea
       IF(i1.ne.1)THEN
        DO iten=1,i1-1
         C(iten,i)=Ten
         Yy(iten,i)=Ten
        END DO
       END IF
       IF(i2.ne.Len)THEN
        DO iten=i2+1,Len
         C(iten,i)=Ten
         Yy(iten,i)=Ten
        END DO
       END IF
       DO iyy=i1,i2
        C(iyy,i)=Ten
        Yy(iyy,i)=Ten
        IF(X(iyy,i).ne.Ten.and.iyy.ge.(i1+1))THEN
         iyy2=iyy-1
         IF(X(iyy2,i).ne.Ten)C(iyy,i)=((X(iyy,i)-X(iyy2,i))/abs(X(iyy2,i
     &                                )))*100.
        END IF
        IF(X(iyy,i).ne.Ten.and.iyy.ge.(i1+Nsea))THEN
         iyy2=iyy-Nsea
         IF(X(iyy2,i).ne.Ten)Yy(iyy,i)=((X(iyy,i)-X(iyy2,i))/abs(X(iyy2,
     &                                 i)))*100.
        END IF
       END DO
      END DO
      RETURN
      END
**==header.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE header
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER Ibreak,Ic,Icq,Icyr,Iq,Irang,Iyr,Kpage,Len,Len2,Lquar,Lyr,
     &        Ncol,Nlen,Np,Ns1,Nsea,Nyear
C*** End of declarations inserted by SPAG
C  *****
C  *****  Prints heading which identifies which options are being used
C  *****  in this sliding spans analysis.
C  *****
      INTEGER Outdec
      CHARACTER*6 Lser
      CHARACTER*80 Title
      COMMON /cee   / Ic,Icyr,Icq
      COMMON /com4b / Iyr,Iq,Nsea,Lyr,Lquar
      COMMON /iopt  / Ibreak,Irang
      COMMON /ns4   / Nlen,Ncol,Ns1,Len,Len2
      COMMON /rite  / Np,Nyear,Title,Kpage,Outdec,Lser
      Kpage=Kpage+1
      WRITE(Np,1010)Title,Kpage,Lser
 1010 FORMAT('1',17X,A80,1X,'PAGE ',I2,', SERIES ',A6)
      WRITE(Np,1020)
 1020 FORMAT(//,' SLIDING SPANS ANALYSIS',//,
     &       ' S 0.   SUMMARY OF OPTIONS SELECTED FOR THIS RUN',//)
      IF(Ibreak.eq.1)WRITE(Np,1030)
      IF(Ibreak.eq.0)WRITE(Np,1040)
 1030 FORMAT('  PRINT ONLY BREAKDOWN TABLES')
 1040 FORMAT('  PRINT ALL TABLES')
      WRITE(Np,1050)Ncol,Nlen,Iq,Iyr
 1050 FORMAT('  YEAR-TO-YEAR CHANGES ANALYZED',/,'  NUMBER OF SPANS : ',
     &       I5,/,'  LENGTH OF SPANS : ',I5,/,'  QUARTER OF FIRST ',
     &       'OBSERVATION IN FIRST SPAN : ',I5,/,'  YEAR    OF FIRST ',
     &       'OBSERVATION IN FIRST SPAN : ',I5)
      IF(Ic.ne.Iq)WRITE(Np,1060)Icyr
 1060 FORMAT('  YEAR    OF FIRST OBSERVATION USED IN SLIDING SPANS ',
     &       'COMPARISON : ',I5)
      WRITE(Np,1070)Lser
 1070 FORMAT('  NAME OF SERIES BEING ADJUSTED : ',A6)
      IF(Ncol.lt.4)WRITE(Np,1080)
 1080 FORMAT(/,'   WARNING: DUE TO THE LENGTH OF THE SERIES, FEWER ',
     &       'THAN FOUR SPANS HAVE BEEN USED TO COMPILE THE MEASURES ',
     &       'GENERATED BELOW.',/,12X,
     &       'IN THIS SITUATION, THE THRESHOLD ',
     &       'VALUES USED TO DETERMINE ADJUSTABILITY (15%, 25%, 40%) ',
     &       'WHICH APPEAR WITH THE',/,12X,'SUMMARY TABLES SHOULD BE ',
     &       'LOWERED.',/)
      RETURN
      END
**==ssftst.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE ssftst(Ncol,Kout)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER i,Issqf,Kout,Ncol
      DOUBLE PRECISION Ssf,Ssmf
C*** End of declarations inserted by SPAG
      CHARACTER*3 qf(4)
      CHARACTER*5 table
      CHARACTER*6 span(4)
      CHARACTER*18 label(2)
      CHARACTER*46 ftfmt
      COMMON /ssft  / Ssf(4),Ssmf(4),Issqf(4)
      DATA table/'(D 8)'/
      DATA span/'SPAN 1','SPAN 2','SPAN 3','SPAN 4'/
      DATA label/'STABLE SEASONALITY','MOVING SEASONALITY'/
      WRITE(Kout,1010)
      WRITE(Kout,1010)
 1010 FORMAT(' ')
      WRITE(Kout,1020)
 1020 FORMAT(5X,'SUMMARY OF F-TESTS FOR STABLE AND MOVING SEASONALITY ',
     &       'FOR EACH SPAN',/)
      ftfmt=' '
      WRITE(ftfmt,1030)Ncol
 1030 FORMAT('(27X,',I1,'(9X,A6))')
      WRITE(Kout,ftfmt)(span(i),i=1,Ncol)
      ftfmt=' '
      WRITE(ftfmt,1040)Ncol
 1040 FORMAT('(5X,A18,3X,A5,',I1,'(3X,F8.2,4X))')
      WRITE(Kout,ftfmt)label(1),table,(Ssf(i),i=1,Ncol)
      WRITE(Kout,ftfmt)label(2),table,(Ssmf(i),i=1,Ncol)
      ftfmt=' '
      WRITE(ftfmt,1050)Ncol
 1050 FORMAT('(5X,25HIDENTIFIABLE SEASONALITY?,9X,',I1,'(A3,12X))')
      DO i=1,Ncol
       qf(i)='   '
       IF(Issqf(i).eq.0)qf(i)='YES'
       IF(Issqf(i).eq.1)qf(i)='???'
       IF(Issqf(i).eq.2)qf(i)=' NO'
      END DO
      WRITE(Kout,ftfmt)(qf(i),i=1,Ncol)
      WRITE(Kout,1010)
      WRITE(Kout,1060)
 1060 FORMAT(/,10X,'YES = IDENTIFIABLE SEASONALITY PRESENT',/,10X,
     &       '??? = IDENTIFIABLE SEASONALITY PROBABLY NOT PRESENT',/,
     &       10X,' NO = IDENTIFIABLE SEASONALITY NOT PRESENT',//)
      RETURN
      END
**==prst.f    processed by SPAG 4.03F  at 13:28 on  2 Mar 1994
      SUBROUTINE prst(Nopt,Irang)
      IMPLICIT NONE
C*** Start of declarations inserted by SPAG
      INTEGER Im,Irang,Iscf,Isign,Itf,Itot,itp,Iturn,Iyr,Kpage,Len,Len2,
     &        Lobs,Lyr,Ncol,Nlen,Nopt,Np,Ns1,Nsea
      INTEGER Nyear
      REAL Outdec,sfper,sper,tfper,tper
C*** End of declarations inserted by SPAG
      CHARACTER*6 Lser
      CHARACTER*80 Title
      COMMON /com1  / Isign,Iscf,Iturn,Itf,Itot
      COMMON /com4b / Iyr,Im,Nsea,Lyr,Lobs
      COMMON /nsl   / Nlen,Ncol,Ns1,Len,Len2
      COMMON /rite  / Np,Nyear,Title,Kpage,Outdec,Lser
      IF(Nopt.ne.2)THEN
       sper=(real(Isign)/real(Itot))*100.
       WRITE(Np,1010)Isign,Itot,sper
       IF(Irang.eq.0)THEN
        sfper=0.0
        IF(Isign.gt.0)sfper=(real(Iscf)/real(Isign))*100.
        WRITE(Np,1020)Iscf,Isign,sfper
       END IF
      END IF
      IF(Ncol.eq.2)RETURN
      itp=Itot-2*Nsea
      IF(itp.gt.0)THEN
       tper=(real(Iturn)/real(itp))*100.
       WRITE(Np,1030)Iturn,itp,tper,Ncol
       IF(Irang.eq.1)RETURN
       tfper=0.0
       IF(Iturn.gt.0)tfper=(real(Itf)/real(Iturn))*100.
       WRITE(Np,1040)Itf,Iturn,tfper
      END IF
      RETURN
 1010 FORMAT(/,2X,I3,' OUT OF ',I3,' (',F5.1,' %)',' OF THE ',
     &       'QUARTERS TESTED HAD A CHANGE OF DIRECTION.')
 1020 FORMAT(2X,I3,' OUT OF ',I3,' (',F5.1,' %)',' OF THE ',
     &       'QUARTER(S) WITH CHANGES OF DIRECTION WERE FLAGGED.',/)
 1030 FORMAT(/,2X,I3,' OUT OF ',I3,' (',F5.1,' %)',' OF THE ',
     &       'QUARTERS TESTED HAD A TURNING POINT OVER ALL ',I2,
     &       ' SPANS.')
 1040 FORMAT(2X,I3,' OUT OF ',I3,' (',F5.1,' %)',' OF THE ',
     &       'QUARTER(S) WITH TURNING POINTS WERE FLAGGED.',//)
      END

