6
$\begingroup$

I have 2 vectors, $a$ & $c$, both of length M. I know they are related by $a*b=c$. My goal is to recover $b$. Obviously $b=$deconv$(c,a)$. I am only interested in the first M elements of the convolution (ie discarding the last M-1 elements). I have two methods of doing this deconvolution in MATLAB and I get different answers. I want to know why:

Method1: using MATLABs default algorithm: b= deconv([c;zeros(M-1,1)],a);

Method2: Using toeplitz matrices:

A=toeplitz(a,[a(1);zeros(M-1,1)]);
b=A^-1*c

where A is the toeplitz matrix form of the vector a with 0s on upper triangle.

I believe in theory these methods should give me identical answers, but the answers are somewhat different (but still close). Why is this?

Thanks!

$\endgroup$
8
$\begingroup$

I think you are constructing the convolution matrix incorrectly. Let's do an easy example with some numbers to see how it is constructed

Let's imagine that youhave the following vectors of dimensions $M$ and $N$, And let's for the sake of it, say that $M = 5$ and $N = 4$

$\mathbf{a} = [a_{1},a_{2},a_{3},a_{4}]^{T}$

$\mathbf{b} = [b_{1},b_{2},b_{3},b_{4},b_{5}]^{T}$

$ \mathbf{c} = \mathbf{a}*\mathbf{b} $

In matrix form, if we assume $\mathbf{a}$ is the system response, the convolution is given by

\begin{equation} c = \mathbf{A}\mathbf{b} \end{equation}

with

\begin{equation} \mathbf{A} = \left[ \begin{array}{ccccc} a_{1} & 0 & 0 & 0 & 0 \\ a_{2} & a_{1} & 0 & 0 & 0 \\ a_{3} & a_{2} & a_{1} & 0 & 0 \\ a_{4} & a_{3} & a_{2} & a_{1} & 0 \\ 0 & a_{4} & a_{3} & a_{2} & a_{1} \\ 0 & 0 & a_{4} & a_{3} & a_{2} \\ 0 & 0 & 0 & a_{4} & a_{3} \\ 0 & 0 & 0 & 0 & a_{4} \\ \end{array}\right] \end{equation}

By doing the convolution with this matrix, you get the full operation including the transients. So the correct matlab code would be

A = toeplitz([a,zeros(1,M-1)],[a(1),zeros(1,M-1)])

(Check help toepliz to verify how it constructs the matrix)

Use this matrix and the result of the convolution is exactly as with the matrix multiplication.

Now, a word of wisdom, as you can clearly see, the full convolution matrix with the transients is not square, so in order to "invert" it and do deconvolution you will need to use the Moore-Penrose Pseudoinverse. And you will get some differences in the beginning and end of the resulting vector, but it is the best you can do, as using the pseudoinverse would be the same as performing least squares, since the deconvolution is given by

$\mathbf{b} = \left(\mathbf{A}^{T}\mathbf{A}\right)^{-1}\mathbf{A}^{T}\mathbf{c}$

This result will always hold, because the dimensions of $A$ are always $(M+N-1)\times M$, so the matrix is always "skinny" and the equation holds.

Try this, it should give you the result you were looking for.

$\endgroup$
  • 1
    $\begingroup$ Perfect answer. Thnx! My mistake was that in the toeplitz matrix I left out the last rows as I didnt think I neeed the last transient part of the convolution, but these rows make a difference when it comes time to take the inverse! $\endgroup$ – DankMasterDan Mar 19 '14 at 18:53

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.