0
$\begingroup$

I have been trying to calculate the inverse of a sweep (by inverse I mean something that when convolved with my original signal will yield a unit impulse) and although I managed to find a time-domain solution, in the process of trying to do it in the frequency domain I encountered an issue which puzzles me a bit.

So, in order to calculate the inverse of the signal I tried the same thing in two slightly different ways. The first was to synthesize the inverse filter by setting its magnitude and phase responses manually and then taking it back to the time domain. I calculated the magnitude as the reciprocal of the original signals magnitude and the phase as the negative of the original signal's phase.

The next try was to calculate the complex spectrum of the original signal, calculate the reciprocal of that and then take this back to the time-domain.

The plots below show the results. I have also plotted the differences between the two approaches (the time domain difference doesn't provide much info, it is just there for visualization purpose).

enter image description here

So, now my question is, why do I get such a different time-domain signals? I know that the correct is the second on but, assuming that most of the residuals originate from round-off errors of the implementation (I used MATLAB's build-in functions for that), I can't think of a good reason.

Any insights would be most welcome here.

UPDATE: For completeness I also provide the MATLAB script I used to calculate the signals.

%% Create signal and calculate the inverses
time = 0:1/1000:10; % Create a time vector
sweep = chirp(time, 1, 10, 500, 'logarithmic'); % Create the sweep

% Artificial inverse filter
artMag = 1./(abs(fft(sweep))); % Magnitude of filter
artPhase = -angle(fft(sweep)); % Phase of filter
artInv = pol2cart(artPhase, artMag); % Convert to Cartesian to IFFT it
artInv = ifft(artInv); % Time-domain signal
artInv = ifftshift(artInv); % Shift the signal to centre it

% Original inverse filter
invSig = 1./fft(sweep); % Calculate the inverse filter
invSigMag = abs(invSig); % Calculate its magnitude
invSigPhase = angle(invSig); % Calculate its phase
invTime = ifft(invSig); % Time-domain signal

%% Plot results
figure(1)
subplot(3, 3, 1)
plot(artMag); grid on
title('Artificial inverse magnitude', 'FontSize', 16)
grid on

subplot(3, 3, 2)
plot(artPhase); grid on
title('Artificial negated phase', 'FontSize', 16)

subplot(3, 3, 3)
plot(artInv); grid on
title('Artificial inverse IFFT', 'FontSize', 16)

subplot(3, 3, 4)
plot(invSigMag); grid on
title('Reciprocal signal magnitude', 'FontSize', 16)

subplot(3, 3, 5)
plot(invSigPhase); grid on
title('Reciprocal signal phase', 'FontSize', 16)

subplot(3, 3, 6)
plot(invTime); grid on
title('Reciprocal signal IFFT', 'FontSize', 16)

subplot(3, 3, 7); grid on
plot(artMag - invSigMag)
title('Magnitude residuals', 'FontSize', 16)

subplot(3, 3, 8); grid on
plot(artPhase - invSigPhase)
title('Phase residuals', 'FontSize', 16)

subplot(3, 3, 9); grid on
plot(artInv - invTime)
title('Residuals in time domain', 'FontSize', 16)
$\endgroup$
  • 1
    $\begingroup$ If magnitude and phase are identical (or almost identical), so should be the time domain responses. It's probably a bug in your code, but we can't tell without seeing it $\endgroup$ – Hilmar Oct 11 at 13:43
  • $\begingroup$ Yep, I added the code. I added the plotting lines just to make sure you know which plot show what. Thanks for the comment, I should had done that from the beginning. $\endgroup$ – ZaellixA Oct 11 at 13:58
1
$\begingroup$

Two bugs in your code.

  1. pol2cart() returns the real and imaginary part, but you only use the real part and discard the imaginary part
  2. You apply ifftshift() on one method but not the other

Try instead

[a,b] = pol2cart(artPhase, artMag);
artInv = ifft(a+1i*b); % Time-domain signal
plot(artInv-invSig);
$\endgroup$
  • $\begingroup$ Yep, this is it. I really thought that pol2cart returns an imaginary number and was treating it as such. Thanks a lot for your time here. As for the ifftshift() I wasn't using it originally but added it at some later stage to visualize the differences. Thanks again... $\endgroup$ – ZaellixA Oct 11 at 19:16

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.