Valid Shape FFT Deconvolution

Hi MathWorks community,

I am trying to wrap my head around FFT convolution and deconvolution, and the effect of different "shape" (full, valid, and same). I have sample code displaying my problem. I'd like to reconstruct the input, x, in both deconvolutions. I got it for full conv/deconv, but not for valid conv/deconv. I didn't even try it for same conv/deconv yet since I don't even know if it is possible for valid conv/deconv. I hope it is possible. I may be missing some theory!

Thank you for any help,

~BH

Answer 1

[Paul]

Hi Benjamin,

Here's an approach for the 'valid' convolution.

Problem parameters:

x = 1:7; nx = numel(x);
h = 1:3; nh = numel(h);

nyf = nx + nh - 1;
nyv = nx - nh + 1;

Define the window function such that conv(x,h,'valid) = w*conv(x,h,'full') with appropriate zero padding as will be shown

d = nh - 1;
w = [zeros(1,d) ones(1,nyv) zeros(1,nyf - nyv - d)];
nw = numel(w);

Verify

yf = conv(x,h)

yf = 1×9

1 4 10 16 22 28 34 32 21

w.*yf

ans = 1×9

0 0 10 16 22 28 34 0 0

conv(x,h,'valid')

ans = 1×5

10 16 22 28 34

Now, just like we could compute yf as

yf = ifft(fft(x,nw).*fft(h,nw))  % same results as above

yf = 1×9

1.0000 4.0000 10.0000 16.0000 22.0000 28.0000 34.0000 32.0000 21.0000

We can compute yv (the zero-padded form of the valid convolution) as

yv = ifft(cconv(fft(w),fft(x,nyf).*fft(h,nyf),nw)/nw)  % same as above  (*)

yv = 
   0.0000 + 0.0000i   0.0000 + 0.0000i  10.0000 + 0.0000i  16.0000 + 0.0000i  22.0000 - 0.0000i  28.0000 + 0.0000i  34.0000 - 0.0000i  -0.0000 - 0.0000i   0.0000 - 0.0000i

Verify the imaginary part is rounding error and get rid of it

max(abs(imag(yv)))

ans = 2.3614e-15

yv = real(yv)

yv = 1×9

0 0.0000 10.0000 16.0000 22.0000 28.0000 34.0000 -0.0000 0.0000

From which conv(x,h,'valid') can be extracted with appropriate indexing.

If we want to recover x starting from fft(yv), it may be possible to invert (*) to get to x. But I couldn't see a way to do that. It seems clear that

fft(yv) = fft(w.*conv(x,h)) = fft(w.*ifft(fft(x,nw).*fft(h,nw)))

but I don't see a way to isolate an expression for fft(x), which could then be ifft'd.

Will be interested to see if anyone comes up with a solution.

Answer 2

[Bruno Luong]

You just do something that is obviously impossible. By FFT with truncation

X = fft(x,L)

You don't use the n-1 last elements so The (n-1) last elements of X is not used to get y.

There is then no way to get back X from Y and h2/H2

%% Valid conv and deconv:
x = 1:7;
h = 1:3;
m = length(x);
n = length(h);
L = m-n+1 

L = 5

X = fft(x,L);
H = fft(h,L);
Y = X.*H;
y = ifft(Y,L)

y = 1×5

23 19 10 16 22

% This x3/h will give the same vector in Fourier domain as x/h
x3 = x; x3([6 7]) = 0;
X3 = fft(x3,L);
Y3 = X3.*H;
y3 = ifft(Y3,L)

y3 = 1×5

23 19 10 16 22

AFAIK the discrete FFT theorem applies for circular convolution. Only FULL convolution is can be convert to FT by Padding 0s and make pb circular without interfering the head and the tail of the original data.

One cannot do that with 'SAME' and 'VALID' convolution. This is impossible, period.

Another way to conceive that only FULL convolution is invertible (interm of x as unknown) is to write down the convolution matrix (in term of flip(h) as band diagonal). You'll see that the simple analysis is matrix size that only FULL convolution makes sense of possibme inversion.

But SAME and VALID convolution are simply truncations of the FULL convolution (that is how @Paul did in his reply).

[Benjamin Hezrony]

@Paul, @Bruno Luong I wish I could formaly accept both of your answers! With your input and some further digging, I've come to understand why what I had been trying was impossible. Thank you both!

~BH