PDA

View Full Version : Forward-Backward Filtering



hcrisp
03-10-2009, 04:43 PM
FIRFILT() does not include the FORWARD_BACK keyword like the IIRFILT() function does. Would there be any advantage to filtering both forward and backward when using an FIR filter structure?



x = ((indgen(1024)+64) MOD 256) GT 128
h = firdesign(101, 0.1, /lowpass)
plot, firfilt(h, x), color=230
oplot, iirfilt(h, x, /forward_back), color=250


They look very similar. Is the raw convolution trimmed such that there is no lag? Why else would FIRFILT() not include the FORWARD_BACK keyword?

rwagner
03-11-2009, 11:37 AM
Hi Hcrisp,
firfilt() does not have a forward-backward keyword for removing phase lag
because it is removed by default; the lag can be retained with /full;

-Ryan

hcrisp
03-11-2009, 12:44 PM
Ok, I see. Actually it's the OFFSET keyword that will retain the full lag.

So the following shows how forward-backward filtering could be done for FIR filters. Notice that the last two curves do not equal each other (look at the far right end of the plot). Any idea why?


x = ((indgen(1024)+64) MOD 256) GT 128
h = firdesign(101, 0.1, /lowpass)
plot, firfilt(h, x, offset=0), color=210 ; with lag
oplot, x ; original
oplot, iirfilt(h, x, /forward_back), color=250 ; forward-backward filtering
; or alternatively
oplot, REVERSE(firfilt(h, REVERSE(firfilt(h, x, o=0)), o=0)), color=230

hcrisp
03-11-2009, 02:39 PM
Another comment I wanted to make was that the default behavior of FIRFILT does not match the default output of FREQRESP_Z. If you use FREQRESP_Z to get the frequency response of an FIR filter and then plot it magnitude-phase, you will see a phase lag -- the group delay. Unfortunately, the FIRFILT function when called without the OFFSET keyword, produces zero-lag by default. Perhaps this was by design, but it is not non-intuitive and inconsistent with many reference texts and other analysis tools in the industry.

To get FREQRESP_Z to match FIRFILT, you must do this:


!P.Multi = [0, 1, 2]

; use h from previous post

; This is the response of an FIR filter with group delay (NOT default FIRFILT behavior)
resp = freqresp_z(h, outfreq=f)
plot, f, arg(resp) / !DTOR, xtitle='Frequency (Hz)', ytitle='Phase (deg)'
plot, f, abs(resp), xtitle='Frequency (Hz)', ytitle='Magnitude'

; This is the response of an FIR filter with zero-phase (default FIRFILT behavior)
window, /free
resp = complex(abs(resp), 0.)
plot, f, arg(resp) / !DTOR, xtitle='Frequency (Hz)', ytitle='Phase (deg)', yrange=[-180,180]
plot, f, abs(resp), xtitle='Frequency (Hz)', ytitle='Magnitude'

rwagner
03-12-2009, 10:25 AM
Notice that the last two curves do not equal each other (look at the far right end of the plot). Any idea why?


The differences are likely due to numerical round-off
since the filter polynomial order is so high.