Program g13bjfe

!     G13BJF Example Program Text

!     Mark 26.1 Release. NAG Copyright 2016.

!     .. Use Statements ..
      Use nag_library, Only: g13bjf, nag_wp, x04caf
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Integer, Parameter               :: nin = 5, nout = 6
!     .. Local Scalars ..
      Integer                          :: dp, i, ifail, imwa, isttf, iwa, kfc, &
                                          kzef, ldparx, ldxxy, mx, n, ncf,     &
                                          ncg, nch, nci, nev, nfv, nis, npara, &
                                          nparx, nser, nsttf, qp, qx, smx
!     .. Local Arrays ..
      Real (Kind=nag_wp), Allocatable  :: fsd(:), fva(:), para(:), parx(:,:),  &
                                          rmsxy(:), sttf(:), wa(:), xxy(:,:)
      Integer                          :: mr(7)
      Integer, Allocatable             :: mrx(:,:), mt(:,:), mwa(:)
!     .. Intrinsic Procedures ..
      Intrinsic                        :: max, sum
!     .. Executable Statements ..
      Write (nout,*) 'G13BJF Example Program Results'
      Write (nout,*)

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

!     Read in problem size
      Read (nin,*) kfc, nev, nfv, nser, kzef

!     Number of input series
      nis = nser - 1

      Allocate (mt(4,nser))

!     Read in the orders for the output noise
      Read (nin,*) mr(1:7)

!     Read in transfer function
      Do i = 1, nis
        Read (nin,*) mt(1:4,i)
      End Do

!     Calculate NPARA
      npara = 0
      Do i = 1, nis
        npara = npara + mt(2,i) + mt(3,i)
      End Do
      npara = npara + mr(1) + mr(3) + mr(4) + mr(6) + nser

!     Calculate array sizes
      n = nev + nfv
      ldxxy = n
      ncf = 0
      Do i = 1, nis
        If (mt(4,i)>1) Then
          ncf = sum(mt(1:3,i))
        End If
      End Do
      isttf = mr(4)*mr(7) + mr(2) + mr(5)*mr(7) + mr(3) +                      &
        max(mr(1),mr(6)*mr(7)) + ncf
      qp = mr(3) + mr(6)*mr(7)
      dp = mr(2) + mr(5)*mr(7)
      smx = 0
      qx = qp
      nci = nser
      Do i = 1, nis
        If (mt(4,i)==3) Then
          mx = max(mt(1,i)+mt(2,i),mt(3,i))
          nci = nci + 1
        Else
          mx = 0
        End If
        If (mt(3,i)>0) Then
          nci = nci + 1
        End If
        smx = smx + mx
        qx = max(qx,mx)
      End Do
      ncg = npara + qx + smx
      nch = dp + 6*qx + nev
      If (qx>0) Then
        nci = nci + 1
      End If
      If (mr(1)>0) Then
        nci = nci + 1
      End If
      If (mr(3)>0) Then
        nci = nci + 1
      End If
      If (mr(4)>0) Then
        nci = nci + 1
      End If
      If (mr(6)>0) Then
        nci = nci + 1
      End If
      iwa = 2*(ncg**2) + nch*(nci+4)
      imwa = 16*nser + 7*ncg + 3*npara + 27
      Allocate (para(npara),xxy(ldxxy,nser),rmsxy(nser),mrx(7,nser),fva(nfv),  &
        fsd(nfv),sttf(isttf),wa(iwa),mwa(imwa))

!     Read in multi-input model parameters
      Read (nin,*) para(1:npara)

!     Read in the observed values for the input and output series
      Read (nin,*)(xxy(i,1:nser),i=1,nev)

!     Read in the future values for the input series
      Read (nin,*)(xxy(nev+i,1:nis),i=1,nfv)

      If (nis>=1) Then
!       Read in residual variance of input series
        Read (nin,*) rmsxy(1:nis)

!       Read in orders for input series ARIMA where available
!       (i.e. where residual variance is not zero)
        ldparx = 0
        Do i = 1, nis
          If (rmsxy(i)/=0.0E0_nag_wp) Then
            Read (nin,*) mrx(1:7,i)
            nparx = mrx(1,i) + mrx(3,i) + mrx(4,i) + mrx(6,i)
            ldparx = max(ldparx,nparx)
          Else
            mrx(1:7,i) = 0
          End If
        End Do
      Else
!       No input series
        ldparx = 1
      End If

      Allocate (parx(ldparx,nser))

!     Read in parameters for each input series ARIMA
      If (nis>0) Then
        Do i = 1, nis
          If (rmsxy(i)/=0.0E0_nag_wp) Then
            nparx = mrx(1,i) + mrx(3,i) + mrx(4,i) + mrx(6,i)
            If (nparx>0) Then
              Read (nin,*) parx(1:nparx,i)
            End If
          End If
        End Do
      End If

      ifail = 0
      Call g13bjf(mr,nser,mt,para,npara,kfc,nev,nfv,xxy,ldxxy,kzef,rmsxy,mrx,  &
        parx,ldparx,fva,fsd,sttf,isttf,nsttf,wa,iwa,mwa,imwa,ifail)

!     Display results
      Write (nout,99999) 'After processing', nev, ' sets of observations'
      Write (nout,99998) nsttf, ' values of the state set are derived'
      Write (nout,*)
      Write (nout,99997) sttf(1:nsttf)
      Write (nout,*)
      Write (nout,*) 'The residual mean square for the output'
      Write (nout,99996) 'series is also derived and its value is',            &
        rmsxy(nser)
      Write (nout,*)
      Write (nout,*) 'The forecast values and their standard errors are'
      Write (nout,*)
      Write (nout,*) '   I       FVA       FSD'
      Write (nout,*)
      Write (nout,99995)(i,fva(i),fsd(i),i=1,nfv)
      Write (nout,*)
      Flush (nout)
      ifail = 0
      Call x04caf('General',' ',n,nser,xxy,ldxxy,                              &
        'The values of z(t) and n(t) are',ifail)
      Write (nout,99994) 'The first ', nis,                                    &
        ' columns hold the z(t) and the last column the n(t)'

99999 Format (1X,A,I3,A)
99998 Format (1X,I3,A)
99997 Format (1X,6F10.4)
99996 Format (1X,A,F10.4)
99995 Format (1X,I4,F10.3,F10.4)
99994 Format (1X,A,I0,A)
    End Program g13bjfe