NAG Library Manual, Mark 30.1
Interfaces:  FL   CL   CPP   AD 

NAG FL Interface Introduction
Example description
    Program g02bzfe
!     G02BZF Example Program Text

!     Mark 30.1 Release. NAG Copyright 2024.

!     .. Use Statements ..
      Use nag_library, Only: dscal, g02buf, g02bzf, nag_wp, x04ccf
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Integer, Parameter               :: nin = 5, nout = 6
!     .. Local Scalars ..
      Real (Kind=nag_wp)               :: alpha, xsw, ysw
      Integer                          :: b, i, ierr, ifail, lc, ldx, m, n
      Character (1)                    :: mean, weight
!     .. Local Arrays ..
      Real (Kind=nag_wp), Allocatable  :: wt(:), x(:,:), xc(:), xmean(:),      &
                                          yc(:), ymean(:)
!     .. Executable Statements ..
      Write (nout,*) 'G02BZF Example Program Results'
      Write (nout,*)
      Flush (nout)

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

!     Read in the problem defining variables
      Read (nin,*) mean, m

!     Allocate memory for output arrays
      lc = (m*m+m)/2
      Allocate (xmean(m),ymean(m),xc(lc),yc(lc))

!     Loop over each block of data
      b = 0
d_lp: Do
!       Read in the number of observations in this block and the weight flag
        Read (nin,*,Iostat=ierr) n, weight
        If (ierr/=0) Then
          Exit d_lp
        End If

!       Keep a running total of the number of blocks of data
        b = b + 1

!       Allocate arrays to hold data and read the current block of data in
        ldx = n
        Allocate (x(ldx,m))
        If (weight=='W' .Or. weight=='w') Then
!         Weighted
          Allocate (wt(n))
          Do i = 1, n
            Read (nin,*) x(i,1:m), wt(i)
          End Do
        Else
!         Unweighted
          Allocate (wt(0))
          Do i = 1, n
            Read (nin,*) x(i,1:m)
          End Do
        End If

!       Summarise this block of data
        If (b==1) Then
!         This is the first block of data, so summarise the data into XMEAN
!         and XC
          ifail = 0
          Call g02buf(mean,weight,n,m,x,ldx,wt,xsw,xmean,xc,ifail)

        Else
!         This is not the first block of data, so summarise the data into
!         YMEAN and YC
          ifail = 0
          Call g02buf(mean,weight,n,m,x,ldx,wt,ysw,ymean,yc,ifail)

!         Update the running summaries
          ifail = 0
          Call g02bzf(mean,m,xsw,xmean,xc,ysw,ymean,yc,ifail)
        End If

        Deallocate (x,wt)
      End Do d_lp

!     Display results
      Write (nout,*) 'Means'
      Write (nout,99999) xmean(1:m)

      Write (nout,*)
      Flush (nout)
      ifail = 0
      Call x04ccf('Upper','Non-unit',m,xc,'Sums of squares and cross-products' &
        ,ifail)

      If (xsw>1.0_nag_wp .And. (mean=='M' .Or. mean=='m')) Then
!       Use DSCAL (F06EDF) to scale the sums of squares and cross-products
!       matrix XC, and so convert it to a covariance matrix
        alpha = 1.0_nag_wp/(xsw-1.0_nag_wp)
        Call dscal(lc,alpha,xc,1)

        Write (nout,*)
        Flush (nout)
        ifail = 0
        Call x04ccf('Upper','Non-unit',m,xc,'Covariance matrix',ifail)
      End If

99999 Format (1X,6F14.4)
    End Program g02bzfe