Program c09ecfe

!     C09ECF Example Program Text
!     Mark 26.1 Release. NAG Copyright 2016.

!     .. Use Statements ..
      Use nag_library, Only: c09abf, c09ecf, c09edf, c09eyf, nag_wp
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Integer, Parameter               :: nin = 5, nout = 6
!     .. Local Scalars ..
      Integer                          :: i, i1, ifail, ilevel, itype_coeffs,  &
                                          j1, lda, ldb, ldcoefs, lenc, m, n,   &
                                          nf, nwcn, nwct, nwl, nwlinv, nwlmax
      Character (10)                   :: mode, wavnam, wtrans
!     .. Local Arrays ..
      Real (Kind=nag_wp), Allocatable  :: a(:,:), b(:,:), c(:), coefs(:,:)
      Integer, Allocatable             :: dwtlvm(:), dwtlvn(:)
      Integer                          :: icomm(180)
!     .. Executable Statements ..
      Write (nout,*) 'C09ECF Example Program Results'
      Write (nout,*)
!     Skip heading in data file
      Read (nin,*)
!     Read problem parameters
      Read (nin,*) m, n
      Read (nin,*) wavnam, mode
      lda = m
      ldb = m
      Allocate (a(lda,n),b(ldb,n))

      Write (nout,99999) wavnam, mode, m, n

!     Read data array and write it out

      Do i = 1, m
        Read (nin,*) a(i,1:n)
      End Do

      Write (nout,*) ' Input Data     A :'
      Do i = 1, m
        Write (nout,99998) a(i,1:n)
      End Do

!     Query wavelet filter dimensions
!     For Multi-Resolution Analysis, decomposition, wtrans = 'M'
      wtrans = 'Multilevel'
      ifail = 0
      Call c09abf(wavnam,wtrans,mode,m,n,nwlmax,nf,nwct,nwcn,icomm,ifail)

      lenc = nwct
      Allocate (c(lenc),dwtlvm(nwlmax),dwtlvn(nwlmax))

      nwl = nwlmax

!     Perform Discrete Wavelet transform
      ifail = 0
      Call c09ecf(m,n,a,lda,lenc,c,nwl,dwtlvm,dwtlvn,icomm,ifail)

      Write (nout,99997) nwl
      Write (nout,99996)
      Write (nout,99995) dwtlvm(1:nwl)
      Write (nout,99994)
      Write (nout,99995) dwtlvn(1:nwl)

!     Allocate an array in which to extract coefficients
      ldcoefs = dwtlvm(nwl)
      Allocate (coefs(ldcoefs,dwtlvn(nwl)))

!     Extract each set of coefficients, working from the deepest level
      Write (nout,99993)
      Do ilevel = nwl, 1, -1
        Write (nout,99992) ilevel, dwtlvm(nwl-ilevel+1), dwtlvn(nwl-ilevel+1)

        Do itype_coeffs = 0, 3
          Select Case (itype_coeffs)
          Case (0)
            If (ilevel==nwl) Then
              Write (nout,99991) 'Approximation coefficients '
            End If
          Case (1)
            Write (nout,99991) 'Vertical coefficients      '
          Case (2)
            Write (nout,99991) 'Horizontal coefficients    '
          Case (3)
            Write (nout,99991) 'Diagonal coefficients      '
          End Select
          If (itype_coeffs>0 .Or. ilevel==nwl) Then
!           Call the 2D extraction routine c09eaf
            Call c09eyf(ilevel,itype_coeffs,lenc,c,coefs,ldcoefs,icomm,ifail)
            Do i1 = 1, dwtlvm(nwl-ilevel+1)
              Write (nout,99989)(coefs(i1,j1),j1=1,dwtlvn(nwl-ilevel+1))
            End Do
          End If
        End Do
      End Do

      nwlinv = nwl

!     Reconstruct original data
      ifail = 0
      Call c09edf(nwlinv,lenc,c,m,n,b,ldb,icomm,ifail)

      Write (nout,99990)
      Do i = 1, m
        Write (nout,99998) b(i,1:n)
      End Do

99999 Format (1X,' MLDWT :: Wavelet  : ',A,/,1X,'          End mode : ',A,/,   &
        1X,'          M        : ',I4,/,1X,'          N        : ',I4,/)
99998 Format (8(F8.4,1X),:)
99997 Format (/,1X,' Number of Levels : ',I10)
99996 Format (1X,' Number of coefficients in 1st dimension for each level :')
99995 Format (8(I8,1X),:)
99994 Format (1X,' Number of coefficients in 2nd dimension for each level :')
99993 Format (/,1X,' Wavelet coefficients C : ')
99992 Format (1X,55('-'),/,1X,' Level : ',I10,'; output is ',I10,' by ',I10,/, &
        1X,55('-'))
99991 Format (1X,A28,': ')
99990 Format (/,1X,' Reconstruction           B : ')
99989 Format (4X,5(F8.4,1X),:)
    End Program c09ecfe