Sunday 6 September 2015

filters - Exponential moving average cut-off frequency


I am trying to implement low pass filter from this example.


What is the cut-off frequency for this type of filter? Is it $F_s \frac{1-\alpha}{2\pi\alpha}$, where $F_s$ is sampling frequency?



Answer




The discrete recurrence relation given on the linked page is


$y[n] = (1-\alpha)y[n-1] + \alpha x[n]$


with $x[n]$ being the input samples and $y[n]$ being the output samples. So taking the Z-transform


$Y(z) = (1-\alpha)z^{-1}Y(z) + \alpha X(z)$


$H(z) = \dfrac{Y(z)}{X(z)}= \dfrac{\alpha}{1-(1-\alpha)z^{-1}}$


To find the -3 dB corner of the filter, in terms of the normalized digital radian frequency, $\Omega$ in the range $[0, \pi]$, one solves


$\left|H(z=e^{j\Omega_{3dB}})\right|^2 = \dfrac{1}{2} = \left|\dfrac{\alpha}{1-(1-\alpha)e^{-j\Omega_{3dB}}}\right|^2$


$ = \dfrac{\alpha^2}{\left|1-(1-\alpha)\cos(-\Omega_{3dB}) -j (1-\alpha)\sin(-\Omega_{3dB})\right|^2}$


which results in this equation


$\left[1-(1-\alpha)\cos(\Omega_{3dB})\right]^2+\left[(1-\alpha)\sin(\Omega_{3dB})\right]^2 = 2\alpha^2$



After some more algebra and trigonometry, I come up with


$\Omega_{3dB} = \cos^{-1}\left[1-\dfrac{\alpha^2}{2(1-\alpha)}\right]$


and since


$f_{3dB} = \dfrac{\Omega_{3dB}}{\pi} \cdot \dfrac{F_s}{2}$


I come up with


$f_{3dB} = \dfrac{F_s}{2\pi}\cos^{-1}\left[1-\dfrac{\alpha^2}{2(1-\alpha)}\right]$


Here's some Octave/MatLab code so you can verify the result on a plot; just change your Fs and alpha as required:


Fs = 10000.0
alpha = 0.01


b = [alpha];
a = [1 -(1-alpha)];

[H, W] = freqz(b, a, 1024, Fs);

f3db = Fs/(2*pi)*acos(1-(alpha^2)/(2*(1-alpha)))

plot(W, 20*log10(abs(H)), [f3db, f3db], [-40, 0]);
grid on


Update in response to comment:


The graph is a plot of the filter's magnitude response in dB (on the y-axis) vs. input frequency in Hz (on the x-axis). The vertical line designates the -3 dB corner frequency.


Since you want to specify your -3 dB corner frequency to find the value for $\alpha$, let's start with this equation from above:


$\left[1-(1-\alpha)\cos(\Omega_{3dB})\right]^2+\left[(1-\alpha)\sin(\Omega_{3dB})\right]^2 = 2\alpha^2$


and after some algebra and trigonometry, one can get a quadratic equation in $\alpha$


$\alpha^2 +2(1- \cos(\Omega_{3dB}))\alpha - 2(1- \cos(\Omega_{3dB})) = 0$


which has solutions


$\alpha = \cos(\Omega_{3db}) - 1 \pm \sqrt{\cos^2(\Omega_{3dB}) -4\cos(\Omega_{3dB}) +3}$


of which only the $+$ answer from the $\pm$ can yield a positive answer for $\alpha$.


Using the above solution for $\alpha$ and the relation



$\Omega_{3dB} = \dfrac{\pi}{F_s/2}f_{3dB}$


one can find $\alpha$ given $f_{3dB}$ and $F_s$.


Here is an updated Octave script:


Fs = 40000
f3db = 1

format long;
omega3db = f3db * pi/(Fs/2)

alpha = cos(omega3db) - 1 + sqrt(cos(omega3db).^2 - 4*cos(omega3db) + 3)


b = [alpha];
a = [1 -(1-alpha)];

[H, W] = freqz(b, a, 32768, Fs);

figure(1);
plot(W, 20*log10(abs(H)), [f3db, f3db], [-40, 0]);
xlabel('Frequency (Hz)');
ylabel('Magnitude Response (dB)');

title('EWMA Filter Frequency Response');
grid on;


W2 = [0:75] * pi/(Fs/2); % 0 to 75 Hz
H2 = freqz(b, a, W2);
W2 = W2 / (pi/(Fs/2));

figure(2);
plot(W2, 20*log10(abs(H2)), [f3db, f3db], [-20, 0]);

xlabel('Frequency (Hz)');
ylabel('Magnitude Response (dB)');
title('EWMA Filter Frequency Response near DC');
grid on;

No comments:

Post a Comment

readings - Appending 内 to a company name is read ない or うち?

For example, if I say マイクロソフト内のパートナーシップは強いです, is the 内 here read as うち or ない? Answer 「内」 in the form: 「Proper Noun + 内」 is always read 「ない...