Something is wrong with my FFT results
31 views (last 30 days)
Show older comments
Nathan Lively
on 9 Mar 2022
Answered: Nathan Lively
on 10 Mar 2022
I've been over this for a couple of hours and I can't find where the error is. I just want to confirm that at each step of the process that I'm doing it correctly so all of the outputs should match.
- Import loudspeaker response data in magnitude and phase.
- Convert to complex number.
- Interpolate.
- IFFT
The first three steps return matching plots, but the plot of the fourth step has something wrong with it. What am I doing wrong??
And here's the magPhase2complex function I made:
function Z = magPhase2complex(mag_dB,phase_deg)
Z = [db2mag(mag_dB) .* exp(1j*(deg2rad(phase_deg)))];
end
some n
% Get some Data
CSVmainFile = '/Users/nathanlively/Downloads/PRX615M_Main.csv';
M_TF = readtable(CSVmainFile,'ReadVariableNames',false);
M_TF.Properties.VariableNames = {'frequency','magnitude','phase','coherence'};
% Convert to complex
M_TF.Z = magPhase2complex(M_TF.magnitude,M_TF.phase);
% Interpolate and fill missing
newFs = 192000;
main = timetable('SampleRate',newFs);
main.Hz = linspace(0,newFs,newFs)';
main.Z = interp1(M_TF.frequency,M_TF.Z,main.Hz,'pchip',NaN);
main.Z = fillmissing(main.Z,'nearest','EndValues','extrap');
main.magnitude = mag2db(abs(main.Z));
% IFFT
main.IR = real(ifft(main.Z)) * height(main);
main.Zrecon = fft(main.IR) / height(main);
main.magnitudeRecon = mag2db(2*abs(main.Zrecon));
semilogx(M_TF.frequency,M_TF.magnitude, main.Hz,main.magnitude, main.Hz,main.magnitudeRecon)
xlim([20 20000]);ylim([max(M_TF.magnitude)-10 max(M_TF.magnitude)]);
legend('original','interp1','reconstructed magnitude','FontSize',22)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/919469/image.png)
16 Comments
Bjorn Gustavsson
on 10 Mar 2022
@Chris Turnes - fair points, the question just struck me as I read your description, and thought that "the optimal procedure surely ought to be the averages of the real and the conjugated complex parts". But I see the use-case and that my flicker of an idea doesn't really apply...
Accepted Answer
More Answers (0)
See Also
Categories
Find more on Logical 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!