; This batch file creates a plot a bandstop filter which suppresses
; frequencies between 7 cycles per second and 15 cycles per second for
; data sampled every 0.02 seconds, using the Hanning window. It is used
; in Chapter 13, "Signal Processing", of _Using IDL_.

@sigprc01   ; compute time data sequence u

; Compute the kaiser filter coefficients:

delt = 0.02 ; sampling period in seconds
f_low = 15. ; frequencies above f_low will be passed
f_high = 7. ; frequencies below f_high will be passed
a_ripple = 50. ; ripple amplitude will be less than -50 dB
nterms = 40 ; the order of the filter

; Compute the impulse response = the filter coefficients

bs_ir_k = DIGITAL_FILTER(f_low*2*delt, f_high*2*delt, a_ripple, nterms)

u_filt = BLK_CON(bs_ir_k, u) ; convolve the Kaiser filter with the signal

v = FFT(u) ; spectrum of original signal

v_filt = FFT(u_filt) ; spectrum of filtered signal

; log-log plot of power spectra

f = FINDGEN(N/2+1) / (N*delt) ; f = [0.0, 1.0/(N*delt), ... , 1.0/(2.0*delt)]

PLOT, f, ABS(v(0:N/2))^2, YTITLE='Power Spectrum', /YLOG, $
   XTITLE='Frequency in cycles / second', /XLOG, $
   XRANGE=[1.0,1.0/(2.0*delt)], XSTYLE=1, $
   TITLE='Spectrum of u(k) Before (solid) and!CAfter (dashed) Digital Filtering'
OPLOT, f, ABS(v_filt(0:N/2))^2, LINESTYLE=2