c	program main
c	call gset							!  enter graphic mode
c	call smoth
c	call aset							!  return to alpha mode
c	end	
c	-----------------------------------------------------------------------
	subroutine smoth
c
c	...demo of display of time series data
c
c	(C) Copyright 1988, 1989 by Jim Farrell      All Rights Reserved.
c
	parameter(n=256)
	real x(n),y(n)
c
	parameter (a2pi=6.283185307,rate=40.0)
c
	irk=1
	do 200 i=1,n
		wt=a2pi*real(i-1)/rate
		ty=10.0*sin(wt)+10.0*cos(2.0*wt)
		x(i)=ty
200	continue
	call putstr(10,20,' DEMONSTRATION OF TPLOT WINDOWS')
	call putstr(12,20,'  FOR A SIGNAL PROCESSING USAGE')
	call pause('WHEN READY ')		!  display message and pause
	call gcls
	xmn=0.0
	xmx=real(n)/rate
	call tplot(xmn,xmx,x,n,1,2)		!  plot clean signal
	call pause ('REVIEW PLOT OF PURE SIGNAL')
	do 300 i=1,n
		if(x(i).ge.0.0)then
			x(i)=x(i)+5.0*rand(irk)		!  add noise
		else
			x(i)=x(i)-5.0*rand(irk)		!  add noise
		endif
300	continue
	call tplot(xmn,xmx,x,n,2,2)		!  plot signal & noise
	call pause ('REVIEW PLOT OF SIGNAL AND NOISE')
	a=0.0
	b=0.0
	c=0.0
	call exsmoth(x,n,0.12,a,b,c,y)	!  smoothing function
	call tplot(xmn,xmx,y,n,3,2)		!  plot smoothed signal
	call pause ('REVIEW PLOT OF SMOOTHED SIGNAL')
	return
	end
c	--------------------------------------------------------------------
	subroutine exsmoth(x,nx,al,a,b,c,xs)
c
c	...smooth a time series
c
c		al - smoothing coefficient (0 < al < 1)
c
	real x(*),xs(*)
c
	if(a.eq.0.0.and.b.eq.0.0.and.c.eq.0.0)then	! compute coefficients
		c=x(1)-2.0*x(2)+x(3)
		b=x(2)-x(1)-1.5*c
		a=x(1)-b-0.5*c
	endif
	be=1.0-al								!  al - smoothing coefficient
	be3=be**3
	al2=al**2
	al3=al*al2
	do 300 i=1,nx
		xs(i)=a+b+0.5*c
		dif=xs(i)-x(i)
		a=x(i)+be3*dif
		b=b+c-1.5*al2*(2.0-al)*dif
		c=c-al3*dif
300	continue
	return
	end