    PROGRAM g13bafe

!      G13BAF Example Program Text

!      Mark 23 Release. NAG Copyright 2011.

!      .. Use Statements ..
       USE nag_library, ONLY : g13ajf, g13baf, nag_wp
!      .. Implicit None Statement ..
       IMPLICIT NONE
!      .. Parameters ..
       INTEGER, PARAMETER              :: nin = 5, nout = 6
!      .. Local Scalars ..
       REAL (KIND=nag_wp)              :: cx, cy, rms
       INTEGER                         :: i, idd, ifail, ifv, ii, ij, ipar,    &
                                          iqxd, ist, iw, nb, nmr, npar, nparx, &
                                          nst, nwa, nx, ny, pp, qp, sy
!      .. Local Arrays ..
       REAL (KIND=nag_wp), ALLOCATABLE :: b(:), fsd(:), fva(:), par(:),        &
                                          parx(:), st(:), w(:), wa(:), x(:),   &
                                          y(:)
       INTEGER                         :: isf(4), mrx(7)
       INTEGER, ALLOCATABLE            :: mr(:)
!      .. Intrinsic Functions ..
       INTRINSIC                          max, min, mod
!      .. Executable Statements ..
       WRITE (nout,*) 'G13BAF Example Program Results'
       WRITE (nout,*)

!      Skip heading in data file
       READ (nin,*)

!      Read in the problem size
       READ (nin,*) nx

!      Read univariate ARIMA for series
       READ (nin,*) mrx(1:7)
       READ (nin,*) cx

!      Calculate number of backforecasts required
       iqxd = mrx(3) + mrx(6)*mrx(7)
       IF (iqxd/=0) THEN
          nmr = 14
       ELSE
          nmr = 7
       END IF

!      Back forecasts will be stored in first IQXD elements
!      of Y, the series will be stored in last NX elements of
!      Y, so calculate start point for the series
       sy = iqxd + 1

!      Calculate length of series with back forecasts
       ny = nx + iqxd

       ALLOCATE (y(ny),mr(nmr))

!      Read in the series into the end of Y
       READ (nin,*) y(sy:ny)

!      Get back forecasts if required
       IF (iqxd/=0) THEN

!         Calculate number of parameters in ARIMA model
          nparx = mrx(1) + mrx(3) + mrx(4) + mrx(6)

          ist = mrx(4) + mrx(7) + mrx(2) + mrx(5) + mrx(3) + &
             max(mrx(1),mrx(6)*mrx(7))
          ifv = max(1,iqxd)
          qp = mrx(6)*mrx(7) + mrx(3)
          pp = mrx(4)*mrx(7) + mrx(1)
          iw = 6*nx + 5*nparx + qp*(qp+11) + 3*pp + 7
          ALLOCATE (parx(nparx),x(nx),st(ist),fva(ifv),fsd(ifv),w(iw))

!         Read in initial values
          READ (nin,*) parx(1:nparx)

!         Reverse series
          x(nx:1:-1) = y(sy:ny)

!         Possible sign reversal for ARIMA constant
          idd = mrx(2) + mrx(5)
          IF (mod(idd,2)/=0) THEN
             cx = -cx
          END IF

!         Calculate back forecasts
          ifail = 0
          CALL g13ajf(mrx,parx,nparx,cx,1,x,nx,rms,st,ist,nst,iqxd,fva,fsd, &
             ifv,isf,w,iw,ifail)

!         Move back forecasts into Y, in reverse order
          y(1:iqxd) = fva(iqxd:1:-1)

!         Reverse sign for ARIMA constant back again
          IF (mod(idd,2)/=0) THEN
             cx = -cx
          END IF
       END IF

!      Read model by which to filter series
       READ (nin,*) mr(1:7)

!      Calculate NPAR
       ipar = mr(1) + mr(3) + mr(4) + mr(6)
       npar = ipar + nparx

       ALLOCATE (par(npar))

!      Read in initial parameter values
       READ (nin,*) par(1:ipar)

       IF (iqxd/=0) THEN
!         Move ARIMA series into MR
          mr(8:14) = mrx(1:7)

!         Move parameters of ARIMA for Y into PAR
          par((ipar+1):(ipar+nparx)) = parx(1:nparx)
       END IF

!      Move constant
       cy = cx

!      Set parameters for call to filter routine G13BAF
       IF (nmr==14) THEN
          nwa = mr(3) + mr(6)*mr(7) + mr(8) + mr(9) + (mr(11)+mr(12))*mr(14)
          nwa = nwa*(nwa+2)
          nb = ny + max(mr(3)+mr(6)*mr(7),mr(1)+mr(2)+(mr(4)+mr(5))*mr(7))
       ELSE
          nwa = 1
          nb = ny
       END IF
       ALLOCATE (wa(nwa),b(nb))

!      Filter series by call to G13BAF
       ifail = 0
       CALL g13baf(y,ny,mr,nmr,par,npar,cy,wa,nwa,b,nb,ifail)

!      Display results
       IF (iqxd/=0) THEN
          WRITE (nout,*) '                 Original        Filtered'
          WRITE (nout,*) 'Backforecasts    y-series         series'
          ij = -iqxd
          DO i = 1, iqxd
             WRITE (nout,99999) ij, y(i), b(i)
             ij = ij + 1
          END DO
          WRITE (nout,*)
       END IF
       WRITE (nout,*) &
          '       Filtered        Filtered        Filtered        Filtered'
       WRITE (nout,*) &
          '        series          series          series          series'
       DO i = iqxd + 1, ny, 4
          WRITE (nout,99998) (ii-iqxd,b(ii),ii=i,min(ny,i+3))
       END DO

99999  FORMAT (1X,I8,F17.4,F15.4)
99998  FORMAT (1X,I5,F9.4,I7,F9.4,I7,F9.4,I7,F9.4)
    END PROGRAM g13bafe
