!
! An example of cubic spline interpolation using spine.f and seval.f
! from netlib.org
!
	parameter(n=11,ns=200)
	dimension xdat(n), ydat(n), qa(n), qb(n), b(n), c(n), d(n)

	print *,'#m=0, S=4'        !this is a command for graph
	do i=0,n-1
		xdat(i+1) = i+sin(i*1.0)/2
		ydat(i+1) = i+cos(i*i*1.0)
		print *, xdat(i+1), ydat(i+1)
	end do
	call qspline (xdat, ydat, n, qa, qb) 
	call spline (n, xdat, ydat, b, c, d) 
	step = ( xdat(n)-xdat(1) )/(ns-1)
	print *,'#m=4, S=0'        !this is a command for graph
	do i=1,ns
		x = xdat(1) + (i-1)*step
		y = qseval(xdat,ydat,qa,qb,n,x) 
		print *, x, y
	end do
	print *,'#m=5, S=0'        !this is a command for graph
	do i=1,ns
		x = xdat(1) + (i-1)*step
		y = seval(n,x,xdat,ydat,b,c,d) 
		print *, x, y
	end do
	print *,'#m=6, S=0'        !this is a command for graph
	do i=1,ns
		x = xdat(1) + (i-1)*step
		y = pinterp(xdat,ydat,n,x) 
		print *, x, y
	end do
	end