Inverse of filter function

Question

azzendix 2017-5-2

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/338325-inverse-of-filter-function

评论： Paul 2024-1-15

采纳的回答： Christoph F.

在 MATLAB Online 中打开

How to write the inverse of the filter function?

I use a filter function with input x to get output y. I want to get input x back.

y = filter(b,a,x,zi)
y = <inverse filter to get original input x>

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Christoph F. 2017-5-2

3
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/338325-inverse-of-filter-function#answer_265307

In the z domain, the transfer function of a filter H(z) is B(z)/A(z).

The inverse of the transfer function is A(z)/B(z). However, this only makes sense if the zeros of B(z) are inside the unit circle; if they are outside, the inverse filter will have poles outside the unit circle and hence be unstable.

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Paul 2024-1-15

在 MATLAB Online 中打开

In addition to the possiblity of the inverse filter having poles outside the unit circle, which is of practical concern, another consideration is the causality of the inverse filter.

Define the signal

rng(100);
N = 1000;
sig = sin(2*pi*40*(0:N-1)*0.001);
n = 0:N-1;

Define an IIR filter with one zero and four poles and filter the signal

[B,A] = zp2tf(0.6,[0.1 0.2 0.3 0.4],1)

B = 1×5

0 0 0 1.0000 -0.6000

A = 1×5

1.0000 -1.0000 0.3500 -0.0500 0.0024

sigf = filter(B,A,sig);

figure

plot(n,sig,n,sigf)

Attempting to filter with the inverse filter casues an error because the inverse filter is non-causal (more zeros than poles)

try
    sigr = filter(A,B,sigf);
catch ME
    ME.message
end
ans = 'First denominator filter coefficient must be non-zero.'

Instead, we mutiply the inverse filter by 1/z^3 (the power is the number of leading zeros in B), then filter

sigr = filter(A,B(end-1:end),sigf);

then shift the output back by three samples.

sigr = sigr(4:end);
nr   = numel(sigr);

This approach recovered the input signal except for the last three samples.

figure

plot(n,sig,n(1:nr),sigr)

figure

plot(n(1:nr),sig(1:nr)-sigr)

请先登录，再进行评论。

Answer 2

Jan 2017-5-2

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/338325-inverse-of-filter-function#answer_265309

This is not possible. Remember that filtering destroys the original information. If you e.g. remove the high frequency noise from a signal, it does not exist in the output anymore. Then there is no way to re-create it.

8 个评论
显示 6更早的评论隐藏 6更早的评论

John D'Errico 2017-5-4

编辑：John D'Errico 2017-5-4

在 MATLAB Online 中打开

Jan, consider this simple problem:

K = ones(1,3)/3;
S0 = rand(1,10);
S1 = conv(S0,K,'same');

Given only S1 and the kernel K, can you recover the original signal? (To within floating point precision.) I've written it in the form of a conv question, but since filter and conv are really the same thing under the hood, this still applies.

Lets try, and see what happens.

We can write the convolution in terms of a matrix multiply.

A = spdiags(repmat(1/3,10,3),-1:1,10,10);

A is full rank.

rank(full(A))
ans =
    10
S1m = S0*A;
S1 - S1m
ans =
   0     0     0     0     0     0     0     0     0     0

And as we see, we accomplish the same result by writing the convolution as a matrix multiplication. And as long as the convolution matrix is nonsingular, we can recover the input. Here, slash applies, since the problem has been posed in linear algebraic form.

S0m = S1/A;
S0 - S0m
ans =
  -2.2204e-16   2.2204e-16            0  -1.1102e-16   1.1102e-16  -2.2204e-16   1.1102e-16            0  -5.5511e-17   8.3267e-17

A is even pretty well conditioned.

cond(full(A))
ans =
     17.255

so we don't expect to lose much of the original signal. If it were true that A was not full rank, then the original signal is not recoverable. Or, if A has a large condition number, then there will be numerical issues. Both cases can be viewed from a signal processing point of view too.

John D'Errico 2017-5-5

I think the difference is, suppose you took the mean of your data. Then the original data is not recoverable. But here you have a linear transformation of the data. N points in, N points out. As long as the mapping is not singular, then there is an inverse transformation.

Paul 2024-1-15

在 MATLAB Online 中打开

@John D'Errico

Using the 'same' option with conv is not the same operation as using filter. In this example, the 'same' option results in an output that is shifte by 1 sample relative to the output of filter. Can the A matrix be changed appropriately?

Probably won't change the general conclusion ....

rng(100)
K = ones(1,3)/3;
S0 = rand(1,10);
S1 = conv(S0,K,'same');
S1a = conv(S0,K);
S1b = filter(K,1,S0);
[S1;S1a(1:numel(S0));S1b]
ans = 3×10
    0.2739    0.4154    0.5159    0.4247    0.3237    0.2657    0.5394    0.5444    0.5126    0.2373
    0.1811    0.2739    0.4154    0.5159    0.4247    0.3237    0.2657    0.5394    0.5444    0.5126
    0.1811    0.2739    0.4154    0.5159    0.4247    0.3237    0.2657    0.5394    0.5444    0.5126