c ---------------------------
      subroutine demod_bloom1976(x,n,omega,d1,d2)
      implicit none
c
c Fourier Analysis of Time Series: An Introduction
c P. Bloomfield (1976: Wiley); p148
c Minor code changes to make it more 'modern'
c
      integer n                        ! INPUT
      real    x(n), omega, d1(n), d2(n)
      integer i                        ! LOCAL
      real    arg

      do i=1,n
         arg   = (i-1)*omega
         d1(i) =  x(i)*cos(arg)*2
         d2(i) = -x(i)*sin(arg)*2
      end do

      return
      end
c ---------------------------
      subroutine polar_bloom1976(x,y,n)
      implicit none
c
c Fourier Analysis of Time Series: An Introduction
c P. Bloomfield (1976: Wiley); p150
c Minor code changes to make it more 'modern'
c
      integer n                        ! INPUT
      real    x(n), y(n)
      integer i                        ! LOCAL
      real    amp, phi

      do i=1,n
         amp = sqrt(x(i)*x(i)+y(i)*y(i))
         phi = atan2(y(i),x(i))        ! radians
         x(i)= amp
         y(i)= phi
      end do

      return
      end
c ---------------------------
c CALL SEQUENCE
c
c omega  = (2*pi*nom)/np2
c cutoff = (pi*(npa+nst))/np2
c ns     = np2/(nst-npa)
c call demod_bloom1976 (x,n,omega,d1,d2) ! DEMODULATION
c call lopass_bloom1976(d1,cutoff,ns)    ! SMOOTH COSINE TERM
c call lopass_bloom1976(d2,cutoff,ns)    ! SMOOTH SINE TERM
c lim = n-2*ns
c call polar_bloom1976(d1,d2,lim)        ! AMPLITUDE & PHASES
c                                        ! USING SMOOTHED VALUES
c ---------------------------
c NCL uses a Lanczos Low Pass Filter
c ---------------------------