home *** CD-ROM | disk | FTP | other *** search
/ Shareware Overload / ShartewareOverload.cdr / progm / grafx.zip / SMOTH.FOR < prev    next >
Text File  |  1989-03-14  |  2KB  |  76 lines
  1. c    program main
  2. c    call gset                            !  enter graphic mode
  3. c    call smoth
  4. c    call aset                            !  return to alpha mode
  5. c    end    
  6. c    -----------------------------------------------------------------------
  7.     subroutine smoth
  8. c
  9. c    ...demo of display of time series data
  10. c
  11. c    (C) Copyright 1988, 1989 by Jim Farrell      All Rights Reserved.
  12. c
  13.     parameter(n=256)
  14.     real x(n),y(n)
  15. c
  16.     parameter (a2pi=6.283185307,rate=40.0)
  17. c
  18.     irk=1
  19.     do 200 i=1,n
  20.         wt=a2pi*real(i-1)/rate
  21.         ty=10.0*sin(wt)+10.0*cos(2.0*wt)
  22.         x(i)=ty
  23. 200    continue
  24.     call putstr(10,20,' DEMONSTRATION OF TPLOT WINDOWS')
  25.     call putstr(12,20,'  FOR A SIGNAL PROCESSING USAGE')
  26.     call pause('WHEN READY ')        !  display message and pause
  27.     call gcls
  28.     xmn=0.0
  29.     xmx=real(n)/rate
  30.     call tplot(xmn,xmx,x,n,1,2)        !  plot clean signal
  31.     call pause ('REVIEW PLOT OF PURE SIGNAL')
  32.     do 300 i=1,n
  33.         if(x(i).ge.0.0)then
  34.             x(i)=x(i)+5.0*rand(irk)        !  add noise
  35.         else
  36.             x(i)=x(i)-5.0*rand(irk)        !  add noise
  37.         endif
  38. 300    continue
  39.     call tplot(xmn,xmx,x,n,2,2)        !  plot signal & noise
  40.     call pause ('REVIEW PLOT OF SIGNAL AND NOISE')
  41.     a=0.0
  42.     b=0.0
  43.     c=0.0
  44.     call exsmoth(x,n,0.12,a,b,c,y)    !  smoothing function
  45.     call tplot(xmn,xmx,y,n,3,2)        !  plot smoothed signal
  46.     call pause ('REVIEW PLOT OF SMOOTHED SIGNAL')
  47.     return
  48.     end
  49. c    --------------------------------------------------------------------
  50.     subroutine exsmoth(x,nx,al,a,b,c,xs)
  51. c
  52. c    ...smooth a time series
  53. c
  54. c        al - smoothing coefficient (0 < al < 1)
  55. c
  56.     real x(*),xs(*)
  57. c
  58.     if(a.eq.0.0.and.b.eq.0.0.and.c.eq.0.0)then    ! compute coefficients
  59.         c=x(1)-2.0*x(2)+x(3)
  60.         b=x(2)-x(1)-1.5*c
  61.         a=x(1)-b-0.5*c
  62.     endif
  63.     be=1.0-al                                !  al - smoothing coefficient
  64.     be3=be**3
  65.     al2=al**2
  66.     al3=al*al2
  67.     do 300 i=1,nx
  68.         xs(i)=a+b+0.5*c
  69.         dif=xs(i)-x(i)
  70.         a=x(i)+be3*dif
  71.         b=b+c-1.5*al2*(2.0-al)*dif
  72.         c=c-al3*dif
  73. 300    continue
  74.     return
  75.     end
  76.