Tuesday, July 10, 2018

Find reverse one pole lowpass filter


I need to find a filter that revert the one pole filter of the current signal, a function using Python (or MATLAB) scipy.signal.filtfilt or scipy.signal.lfilter



The "original" signal that i want to recover is :


[0, 0, 0, 0, 0, 0, 0, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2] 

The "filtered" signal that i have is:


[ 0., 0., 0., 0., 0., 0., 0., 0.25956786, 0.48537242, 0.68182373, 0.85275161, 1.00148237, 1.1309067, 1.24353671, 1.34155607, 1.42686391, 1.50111103, 1.62198067, 1.67093897, 1.71355379, 1.81104326, 1.83550966, 1.85680735, 1.87534702, 1.89148581, 1.90553486, 1.9177649, 1.92841125, 1.93767929, 1.95277071]

I want find the (b, a) values to execute the "original=signal.filtfilt(b,a,filtered)"


Original = Red Filtered = Blue


filtered=blue, original=red



Answer




The blue curve does not seem accurate and is a little non-smooth. However, I will give a single-pole filter to approximate it.


It is not difficult to show that the step response of the filter $$H(z)=\frac{1-p}{1-pz^{-1}}$$ is in the form of


$$s[n]=(1-p^{n+1})u[n]$$


To find the step response, it suffices to shift the signal given in the question to the left and divide it by $2$. Hence, the output at zero is $0.25956786/2= 0.1298$. To find the parameter $p$, substitute it in the given equation for $s[n]$ at $n=0$:


$$0.1298=1-p\Rightarrow p=0.8702$$ So the filter we were looking for is $$H(z)=\frac{0.1298}{1-0.8702z^{-1}}$$ More clearly, $b=0.1298$ and $a=[1,-0.8702]$.


Here is a comparison in Matlab: enter image description here


It can be seen that the first $10$ outputs are very close to the actual output. However, after the output becomes non-smooth, there is a small error. Alternatively, we could use the term at $n=22$ for our reference. In such case, $$1.95277071/2=1-p^{23}\Rightarrow p=0.8497.$$ This will lead to an accurate approximation at large $n$, but an error at small $n$.


Using the code


y = [0.25956786, 0.48537242, 0.68182373, 0.85275161, 1.00148237, 1.1309067, 1.24353671, 1.34155607, 1.42686391, 1.50111103, 1.62198067, 1.67093897, 1.71355379, 1.81104326, 1.83550966, 1.85680735, 1.87534702, 1.89148581, 1.90553486, 1.9177649, 1.92841125, 1.93767929, 1.95277071]
fun = @(P)2*(1-P.^(1:23))-y;

p = lsqnonlin(fun,.5);

$p=0.8644$ gives the least squares (LS) solution to this problem which is compared to the original output in the following plot ($\text{SSE}=0.02796$):


enter image description here


Note that, if the forward filter is $(b,a)$, then $(a,b)$ is the inverse filter which gives you the original input. If the error in the step response cannot be fixed, you can apply the LS approach on the overall output (taking into account the inverse filter) to get a better result.


No comments:

Post a Comment

periodic trends - Comparing radii in lithium, beryllium, magnesium, aluminium and sodium ions

Apparently the of last four, $\ce{Mg^2+}$ is closest in radius to $\ce{Li+}$. Is this true, and if so, why would a whole larger shell ($\ce{...