fft of even number of elements of pulse yields one or multiple zeros in the array. How do I avoid this?

2 views (last 30 days)
I am trying to deconvolve two signals using fft(). The issue I'm having is with the pulse width. If the pulse width is an even number of ones with zero padding, the fft() yields a zero in the resulting array. Therefore, when I try to divide my Y_fft by H_fft, I get NaN as the result.
The code below zeros pads a pulse to some power of 2 because my signal I want deconvolve will be forced to this length. I have played around with finding the zeros and replacing them 1e-15 or replacing them with the average of the surrounding values. Both are quick fixes with no real math to prove why to do this. Therefore, my deconvolved signal is reconstructed but with some jitter. I would greatly appreciate feedback on the correct way to handle this. At the very least, some understand as to why even number of pulses yield zeros after the fft().
Note: The position of the zeros seem to also be at location of powers of 2 and shifted.
if true
clc
clear
close all
%original signal
x = [1 2 3 4 5 6 5 4 3 2 1];
lx = length(x);
%odd pulse case
pulses_odd = 3;
h1 = ones(1,pulses_odd);
%creation of convolved signal
y1 = conv(x,h1);
ly1 = length(y1);
y1 = [y1 zeros(1, pow2(nextpow2(ly1)) - ly1)];
ly1 = length(y1);
%zero padding impulse function
h1 = [ones(1,pulses_odd) zeros(1,ly1-pulses_odd)];
Y_fft1 = fft(y1);
H_fft1 = fft(h1);
x1_decon = Y_fft1./H_fft1;
x1_decon = ifft(x1_decon);
%even pulse case
pulses_even = 4;
h2 = ones(1,pulses_even);
%creation of convolved signal
y2 = conv(x,h2);
ly2 = length(y2);
y2 = [y2 zeros(1, pow2(nextpow2(ly2)) - ly2)];
ly2 = length(y2);
%zero padding impulse function
h2 = [ones(1,pulses_even) zeros(1,ly2-pulses_even)];
Y_fft2 = fft(y2);
H_fft2 = fft(h2);
x2_decon = Y_fft2./H_fft2;
x2_decon = ifft(x2_decon);
plot(x,'linewidth',2)
hold on
plot(x1_decon,'d')
plot(x2_decon,'s')
hold off
end

Accepted Answer

David Goodmanson
David Goodmanson on 1 Nov 2018
Hi Keith,
This is happening because the array length of 16 (after you zerofill) and a pulse length of 4 have a factor in common. For an array length of 16, a wave of period 4 fits (4 periods give 16) so that wave can contribute to the fft.
For the constant-height pulse of length 4, a wave of period 4 exactly averages out over the width of the pulse and the fft amplitude due to the nonzero part of the pulse is zero. It doesn't matter what is going on in the 12 points off the end of the pulse because the pulse height is zero there and can't contribute to the fft. So you get zero.
The same argument shows that you get zero for a wave of period 2. For a pulse length of 3 or 5 or anything that does not have a factor in common with 16, a similar argument shows that zero cannot happen.
One way to fix this, similar to what you are doing, is to introduce a small element into the pulse function that doesn't affect the convolution to any notable extent but does not allow a zero in the fft. The annotated line below does that. If the length of the pulse is an even number, including that line gives a good result and commenting it out gives NaNs.
%original signal
x = [1 2 3 4 5 6 5 4 3 2 1];
pulses = 4;
h = ones(1,pulses);
h(1) = h(1) + 1e-8; % assure nonzero fft
%creation of convolved signal
y = conv(x,h);
N = pow2(nextpow2(length(y)))
y(end+1:N) = 0 % zerofill
h(end+1:N) = 0; % zerofill
z = fft(y);
zh = fft(h);
x_decon = ifft(z./zh);
figure(1)
plot(x,'linewidth',2)
hold on
plot((x_decon),'d')
hold off
  1 Comment
Keith Hernandez
Keith Hernandez on 2 Nov 2018
Thank you. I appreciate the insight. I implemented the your suggestion and got acceptable results. Not sure why they didn't cover this in my undergrad class.

Sign in to comment.

More Answers (0)

Categories

Find more on Programming in Help Center and File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!