Program g10abfe
! G10ABF Example Program Text
! Mark 26.2 Release. NAG Copyright 2017.
! .. Use Statements ..
Use nag_library, Only: g10abf, g10zaf, nag_wp
! .. Implicit None Statement ..
Implicit None
! .. Parameters ..
Integer, Parameter :: nin = 5, nout = 6
! .. Local Scalars ..
Integer :: i, ifail, j, lcomm, ldc, lwt, n, &
nord, nrho
Character (1) :: mode, weight
! .. Local Arrays ..
Real (Kind=nag_wp), Allocatable :: c(:,:), comm(:), df(:), h(:), &
res(:), rho(:), rss(:), wt(:), &
wwt(:), x(:), xord(:), y(:), &
yhat(:,:), yord(:)
Integer, Allocatable :: iwrk(:)
! .. Executable Statements ..
Write (nout,*) ' G10ABF Example Program Results'
Write (nout,*)
! Skip heading in data file
Read (nin,*)
! Read in problem size and control parameters
Read (nin,*) n, weight, nrho
If (weight=='W' .Or. weight=='w') Then
lwt = n
Else
lwt = 0
End If
lcomm = 9*n + 14
ldc = n - 1
Allocate (x(n),y(n),wt(lwt),xord(n),yord(n),wwt(n),yhat(n,nrho), &
c(ldc,3),res(n),h(n),comm(lcomm),iwrk(n),rho(nrho),rss(nrho),df(nrho))
! Read in the smoothing parameters
Read (nin,*) rho(1:nrho)
! Read in data
If (lwt>0) Then
Read (nin,*)(x(i),y(i),wt(i),i=1,n)
Else
Read (nin,*)(x(i),y(i),i=1,n)
End If
! Reorder data into increasing X and remove tied observations, weighting
! accordingly
ifail = 0
Call g10zaf(weight,n,x,y,wt,nord,xord,yord,wwt,rss(1),iwrk,ifail)
! Fit cubic spline the first time
! NB: These are weighted as G10ZAF creates weights
ifail = 0
mode = 'P'
Call g10abf(mode,'W',nord,xord,yord,wwt,rho(1),yhat(1,1),c,ldc,rss(1), &
df(1),res,h,comm,ifail)
! Fit cubic spline the remaining NRHO - 1 times
mode = 'Q'
Do i = 2, nrho
ifail = 0
Call g10abf(mode,'W',nord,xord,yord,wwt,rho(i),yhat(1,i),c,ldc,rss(i), &
df(i),res,h,comm,ifail)
End Do
! Display results
Write (nout,99999) 'Smoothing coefficient (rho) = ', rho(1:nrho)
Write (nout,99998) 'Residual sum of squares = ', rss(1:nrho)
Write (nout,99998) 'Degrees of freedom = ', df(1:nrho)
Write (nout,*)
Write (nout,*) ' X Y Fitted Values'
Do i = 1, nord
Write (nout,99997) xord(i), yord(i), (yhat(i,j),j=1,nrho)
End Do
99999 Format (1X,A,10(2X,F8.2))
99998 Format (1X,A,10(F10.3))
99997 Format (1X,2F8.4,14X,10(2X,F8.4))
End Program g10abfe