NAG Library Manual, Mark 28.5
```    Program g13bbfe

!     G13BBF Example Program Text

!     Mark 28.5 Release. NAG Copyright 2022.

!     .. Use Statements ..
Use nag_library, Only: g13ajf, g13bbf, 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, iwa, nb, nmr, npar,   &
nparx, nst, 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 Procedures ..
Intrinsic                        :: max, min, mod
!     .. Executable Statements ..
Write (nout,*) 'G13BBF Example Program Results'
Write (nout,*)

!     Skip heading in data file

!     Read univariate ARIMA for series

!     Calculate number of backforecasts required
iqxd = mrx(3) + mrx(6)*mrx(7)
If (iqxd/=0) Then
nmr = 10
Else
nmr = 3
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))

!     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))

!       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,0,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

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

Allocate (par(npar))

!     Read in initial parameter values

If (iqxd/=0) Then
!       Move ARIMA for series into MR
mr(4:10) = 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 G13BBF
If (nmr==10) Then
iwa = mr(3) + mr(4) + mr(5) + (mr(7)+mr(8))*mr(10)
iwa = npar + iwa*(iwa+2)
nb = ny + max(mr(1)+mr(2),mr(3))
Else
iwa = mr(1) + npar
nb = ny
End If

Allocate (wa(iwa),b(nb))

!     Filter series by call to G13BBF
ifail = 0
Call g13bbf(y,ny,mr,nmr,par,npar,cy,wa,iwa,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.1,F16.1)
99998 Format (1X,I5,F10.1,I6,F10.1,I6,F10.1,I6,F10.1)
End Program g13bbfe
```