# Impulse response of Butterworth filter

**URL:** <https://discourse.julialang.org/t/impulse-response-of-butterworth-filter/51543>\
**Category:** Signal and Image Processing\
**Created:** [December 9, 2020, 7:10pm UTC](https://discourse.julialang.org/t/impulse-response-of-butterworth-filter/51543 "2020-12-09T19:10:01Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![rafael.guerra](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rafael.guerra/32/216610_2.png) [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)\
**Post date:** [December 9, 2020, 7:10pm UTC](https://discourse.julialang.org/t/impulse-response-of-butterworth-filter/51543/1 "2020-12-09T19:10:01Z")

</div>

In the following code example, a 6th-order Butterworth filter is defined over 4-50 Hz and the frequency response plotted. How to use [DSP.jl](https://docs.juliadsp.org/stable/contents/) to compute the corresponding zero-phase filter time response?  
Using `filt(fil,ht1)` in code example below, it seems to provide the minimum-phase response, while the workaround using `ifft` and `fftshift` outputs what I expect but it is rather convoluted.

```julia
using DSP, FFTW, Plots

fs = 1000; Δt = 1.0/fs; # sampling frequency [Hz] amd period [s]
fc1, fc2 = 4, 50; # low-cut / high-cut frequencies [Hz]
order = 6; # Butterworth filter order
fil = digitalfilter(Bandpass(fc1,fc2,fs=fs), Butterworth(order));
Hf = freqz(fil,0:1:200,fs); # 0-200 Hz response
plot(20*log10.(abs.(Hf)),title="Frequency Response", label="Butterworth 6th order (4-50 Hz)",
     xlabel="Frequency [Hz]", ylabel="Amplitude [dB]")

n = 500;
t = Δt * vec(-n÷2:n÷2-1);

# minimum phase impulse response?
ht1 = [1.0; zeros(n-1,1)]; # impulse at time 0
ht1 = filt(fil,ht1)
ht1 = [zeros(n÷2,1); ht1[1:n÷2]]
plot(t, ht1, title="Impulse Response", label="Min-phase (4-50 Hz)")

# 0-phase filter time response (using ifft)
ht0 = fftshift(ifft(fftshift(abs.(freqz(fil,-n÷2:n÷2-1,fs))))); # over -250 to +250 Hz
plot(t, real.(ht0), title="0-phase filter Time Response", label="Zero-phase (4-50 Hz)",
     xlabel="Time [s]", ylabel="Amplitude")

```

![Butterworth_6th_order_4_50Hz_frequency_response](https://global.discourse-cdn.com/julialang/original/3X/6/6/66497e61f93f66a71eb262a5354e10405328ee05.png)

![Butterworth_6th_order_4_50Hz_min_phase_impulse_response](https://global.discourse-cdn.com/julialang/original/3X/d/0/d0f541f571d7f9b3f0346b9f85f778dcf9a75fe9.png)

![Butterworth_6th_order_4_50Hz_zero-phase_time_response](https://global.discourse-cdn.com/julialang/original/3X/f/1/f15f977d9e99b976cb01d8e6fed9609fab92e283.png)

Thanks in advance.

---

<div class="post-metadata">

**Author:** ![mgkuhn](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mgkuhn/32/6276_2.png) [@mgkuhn](https://discourse.julialang.org/u/mgkuhn)\
**Post date:** [December 11, 2020, 10:58pm UTC](https://discourse.julialang.org/t/impulse-response-of-butterworth-filter/51543/2 "2020-12-11T22:58:57Z")

</div>

Try constructing an order-3 filter and then apply it with [filtfilt](https://docs.juliadsp.org/stable/filters/#DSP.Filters.filtfilt). This doubles the order while making the filter zero phase.

---

<div class="post-metadata">

**Author:** ![rafael.guerra](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rafael.guerra/32/216610_2.png) [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)\
**Post date:** [December 13, 2020, 3:56pm UTC](https://discourse.julialang.org/t/impulse-response-of-butterworth-filter/51543/3 "2020-12-13T15:56:56Z")

</div>

@mgkuhn, thank you for the feedback. I did not succeed in making DSP `filtfilt` to work in this case. However, the same principle of filtering in two-passes and flipping in time the first filter output, works if done manually:

```julia
t = Δt * vec(-(n-1):n-1);
ht = [zeros(n-1,1); 1.0; zeros(n-1,1)]; # impulse at time 0
fil3 = digitalfilter(Bandpass(fc1,fc2,fs=fs), Butterworth(3));
ht3 = filt(fil3,ht)
ht6 = filt(fil3, reverse(ht3, dims=1))

```

This approach seems to be more restrictive (even filter order) than the workaround using `ifft` and `fftshift`. But I am not a `DSP` expert, just an interested user. The documentation for [DSP.jl](https://docs.juliadsp.org/latest/contents/) would benefit from more examples.

---

<div class="post-metadata">

**Author:** ![mgkuhn](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mgkuhn/32/6276_2.png) [@mgkuhn](https://discourse.julialang.org/u/mgkuhn)\
**Post date:** [December 17, 2020, 5:21pm UTC](https://discourse.julialang.org/t/impulse-response-of-butterworth-filter/51543/4 "2020-12-17T17:21:52Z")

</div>

Works fine for me (Julia v1.5.3, DSP v0.6.9):

```julia
using Plots
using DSP

fs = 1000; Δt = 1.0/fs; # sampling frequency [Hz] amd period [s]
fc1, fc2 = 4, 50; # low-cut / high-cut frequencies [Hz]
fil3 = digitalfilter(Bandpass(fc1,fc2,fs=fs), Butterworth(3));
n = 500;
ht = [zeros(n-1,1); 1.0; zeros(n-1,1)]; # impulse at time 0
y = filtfilt(fil3, ht);
plot(y)

```

![imp](https://global.discourse-cdn.com/julialang/original/3X/f/1/f1bcf37a52a5f024acb55242e1f15c6bcef63566.png)

---

<div class="post-metadata">

**Author:** ![rafael.guerra](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rafael.guerra/32/216610_2.png) [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)\
**Post date:** [December 17, 2020, 5:31pm UTC](https://discourse.julialang.org/t/impulse-response-of-butterworth-filter/51543/5 "2020-12-17T17:31:54Z")

</div>

Looks good. I had only tested `filtfilt` with an impulse on first sample: `ht = [1.0; zeros(n-1,1)]`, which seems to chop the signal.