MATLAB：如何正确应用 ifft 以将“过滤”信号带回时域？

Question

我正在尝试通过同轴电缆获得高斯脉冲的输出。我制作了一个代表同轴电缆的矢量；我在网上得到了衰减和相位延迟信息，并使用欧拉方程创建了一个复杂的数组。

我对我的高斯向量进行 FFT 运算，并用电缆对其进行卷积。问题是，我不知道如何正确地 iFFT 卷积。我在 MathWorks 中阅读了 iFFt 并查看了其他人的问题。有人遇到了类似的问题，在答案中，有人建议删除 n = 2^nextpow2(L) 和 FFT over length(t) 。我能够从中获得更合理的情节，并且为什么会这样是有道理的。我对是否应该在 iFFt 中使用对称选项感到困惑。这对我的情节产生了很大的影响。我添加对称性的主要原因是因为我在 iFFTed 卷积（timeHF）中得到了复数。非常感谢您的帮助，谢谢！

clc, clear
Fs = 14E12;          %1 sample per pico seconds
tlim = 4000E-12;
t = -tlim:1/Fs:tlim; %in pico seconds
ag = 0.5;            %peak of guassian 
bg = 0;              %peak location
wg = 50E-12;         %FWHM

x = ag.*exp(-4 .* log(2) .* (t-bg).^2 / (wg).^2); %Gauss. in terms of FWHM
Ly = x;

L = length(t);
%n = 2^nextpow2(L);  %test output in time domain with and without as suggested online
fNum = fft(Ly,L);

frange = Fs/L*(0:(L/2)); %half of the spectrum
fNumMag = abs(fNum/L);   %divide by n to normalize


% COAX modulation ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

%phase data
mu = 4*pi*1E-7;
sigma_a = 2.9*1E7;
sigma_b = 5.8*1E6;
a = 0.42E-3;
b = 1.75E-3;
er = 1.508;
vf = 0.66;
c = 3E8;
l = 1;
Lso = sqrt(mu) /(4*pi^3/2) * (1/(sqrt(sigma_a)*a) + 1/(b*sqrt(sigma_b)));
Lo = mu/(2*pi) * log(b/a);
%to = l/(vf*c);
to = 12E-9; %measured
phase = -pi*to*(frange + 1/2 * Lso/Lo * sqrt(frange));

%attenuation Data
k1 = 0.34190;
k2 = 0.00377;
len = 1;
mldb = (k1 .* sqrt(frange) + k2 .* frange) ./ 100 .* len ./1E6;
mldb1 = mldb ./ 0.3048; %original eqaution is in inch
tfMag = 10.^(mldb1./-10);
  
% combine to make in complex form
 tfC = [];
     for ii = 1: L/2 + 1
        tfC(ii) = tfMag(ii) * (cosd(phase(ii)) + 1j*sind(phase(ii)));
     end
     
%END ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

%convolute both h and signal
fNum = fNum(1:L/2+1);
convHF = tfC.*fNum;
convHFMag = abs(convHF/L);
timeHF = ifft(convHF, length(t), 'symmetric');  %this is the part im confused about

% Ignore, 
% tfC(numel(fNum)) = 0;
% convHF = tfC.*fNum;
% convHFMag = abs(convHF/n);
% timeHF = ifft(convHF);

%% plotting

% subplot(2, 2, 1);
% plot(t, Ly)
% title('Gaussian input');
% xlabel('time in seconds')
% ylabel('V')
% grid

subplot(2, 2, 1)
plot(frange, abs(tfC(1: L/2 + 1)));
set(gca, 'Xscale', 'log')
title('coax cable model')
xlabel('Hz')
ylabel('|H(s)|V/V')
grid
ylim([0 1.1])

subplot(2, 2, 2);
plot(frange, convHFMag(1:L/2+1), '.-', frange, fNumMag(1:L/2+1)) %make both range and function the same lenght
title('The input signal Vs its convolution with coax');
xlabel('Hz')
ylabel('V')
legend('Convolution','Lorentzian in frequecuency domain');
xlim([0, 5E10])
grid

subplot(2, 2, [3, 4]);
plot(t, Ly, t, timeHF)
% plot(t, real(timeHF(1:length(t)))) %make both range and function the same lenght
legend('Input', 'Output')
title('Signal at the output')
xlabel('time in seconds')
ylabel('V')
grid

Answer 1

深入理解 FFT 的原理以正确使用它很重要。

当您将傅里叶变换应用于实信号时，负频率系数是正频率系数的共轭。当您将 FFT 应用于实数值信号时，您可以在数学上表明应该在负频率 (-f) 的系数的共轭现在将出现在 (Fsampling-f) 处，其中 Fsampling=1/dt 是采样频率，并且dt 采样周期。这种行为称为混叠，当您将 fft 应用于离散时间信号时会出现这种现象，并且应选择足够小的采样周期以使这两个光谱不与香农标准重叠。

当您想对信号应用频率滤波器时，我们说我们保留频谱的前半部分，因为高频 (>Fsampling/2) 是由混叠引起的，不是原始信号的特征。为此，在乘以滤波器之前，我们在光谱的后半部分置零。但是，这样做还会丢失原始信号幅度的一半，而使用 ifft 将无法恢复。选项“对称”可以通过在高频 (>Fsampling/2) 中添加较低的系数 (

我简化了代码以简要解释正在发生的事情，并在第 20 行为您实现了手工对称化。请注意，我将采样周期从 1 皮秒减少到 100 皮秒，以便正确显示频谱：

close all
clc, clear
Fs = 14E10;          %1 sample per pico seconds % CHANGED to 100ps
tlim = 4000E-12;
t = -tlim:1/Fs:tlim; %in pico seconds
ag = 0.5;            %peak of guassian 
bg = 0;              %peak location
wg = 50E-12;         %FWHM
NT = length(t);

x_i = ag.*exp(-4 .* log(2) .* (t-bg).^2 / (wg).^2); %Gauss. in terms of FWHM
fftx_i = fft(x_i);
f = 1/(2*tlim)*(0:NT-1);

fftx_r = fftx_i;
fftx_r(floor(NT/2):end) = 0; % The removal of high frequencies due to aliasing leads to losing half the amplitude
% HER YOU APPLY FILTER
x_r1 = ifft(fftx_r); % without symmetrisation (half the amplitude lost)
x_r2 = ifft(fftx_r, 'symmetric'); % with symmetrisation
x_r3 = ifft(fftx_r+[0, conj(fftx_r(end:-1:2))]); % hand-made symmetrisation

figure();
subplot(211)
hold on
plot(t, x_i, 'r')
plot(t, x_r2, 'r-+')
plot(t, x_r3, 'r-o')
plot(t, x_r1, 'k--')
hold off
legend('Initial', 'Matlab sym', 'Hand made sym', 'No sym')
title('Time signals')
xlabel('time in seconds')
ylabel('V')
grid
subplot(212)
hold on
plot(f, abs(fft(x_i)), 'r')
plot(f, abs(fft(x_r2)), 'r-+')
plot(f, abs(fft(x_r3)), 'r-o')
plot(f, abs(fft(x_r1)), 'k--')
hold off
legend('Initial', 'Matlab sym', 'Hand made sym', 'No sym')
title('Power spectra')
xlabel('frequency in hertz')
ylabel('V')
grid

绘制结果：

如果您还有其他问题，请不要犹豫。祝你好运！

--------- 编辑 ----------

离散傅里叶变换的幅度与连续傅里叶变换的幅度不同。如果您对在频域中显示信号感兴趣，则需要根据您选择的约定应用归一化。通常，您使用 Dirac delta 函数的傅里叶变换的幅度在任何地方都具有幅度 1 的约定。

数值 Dirac delta 函数在索引处的幅度为 1，在其他地方为零，并导致功率谱在任何地方都等于 1。但是，在您的情况下，时间轴具有采样周期 dt，在这种情况下，数字狄拉克随时间的积分不是 1，而是 dt。您必须通过将频域信号乘以因子 dt（在您的情况下为 = 1picoseceond）来规范化频域信号以遵守约定。您还可以注意到，这使得频域信号与 [原始单位乘以时间] 是均匀的，这是傅里叶变换的正确单位。