通过“中”频率，我假设您想要滤除信号，以便仅存在 159 Hz 的分量。这意味着输出结果应该只包含一个 159 Hz 的正弦波。首先，我想提一下傅里叶变换显示了信号的频率分解。您可以将信号视为不同频率的许多正弦曲线的总和，它们可能因不同的相移而异相。对于一个频率的每个正弦波，都有一个幅度分量和一个相位分量。幅度告诉您正弦曲线在该频率下有多大，相位告诉您正弦曲线正在经历多少延迟。
现在，在计算 FFT 时，它执行Cooley-Tukey 算法，这是一种计算傅里叶变换的非常有效的方法。它的计算方式是，信号的前半部分显示来自0 <= f < fn 的频谱，而信号的后半部分显示来自-fn <= f < 0 的频谱。 fn 就是所谓的Nyquist frequency，它是采样频率的一半。请注意与后半部分相比，前半部分光谱范围内的排他性。我们不在前半部分包含fn，但我们在后半部分包含-fn。
这仅适用于 真实 信号，即您将三个正弦信号加在一起的情况。为了使这更有意义，您可以使用fftshift，以便重新组织频谱，使 0 Hz / DC 频率出现在中间，这将允许您绘制频谱，使其介于 -fn <= f < fn .
此外，信号的绘制方式是标准化，这意味着您看到的频率实际上在-1 <= f <= 1 之间。要将其转换为实际频率，请使用以下关系：
freq = i * fs / N
i 是您想要的 FFT 上的 bin 编号或点，在您的情况下是从 0 到 500。请记住，信号的前半部分表示从 0 到 fn 的频率分布和信号中从 0 到 500 是您所需要的。只需将i 替换为-i 即可找到负频率。 fs 是采样频率，N 是 FFT 的大小，在您的情况下是信号的总长度，即 1001。要生成对应于正确点的正确频率，您可以使用 linspace并在-fn 和fn 之间生成N+1 点以确保间距正确，但由于我们不在正端包含fn，因此我们将其从末尾删除范围。
另外，要绘制信号的幅度和相位，请分别使用abs 和angle。
因此，请尝试以这种方式绘制它，并同时关注幅度和相位。仅仅绘制光谱是不确定的，因为一般情况下，有实部和虚部。
%// Your code
fs = 1000;
t = (0 : 1/fs : 1)';
f = [ 86, 159, 392 ];
x = sum( cos(2*pi * t * f), 2 );
y = fft(x);
%// New code
%// Shift spectrum
ys = fftshift(y);
N = numel(x); %// Determine number of points
mag = abs(ys); %// Magnitude
pha = angle(ys); %// Phase
%// Generate frequencies
freq = linspace(-fs/2, fs/2, N+1);
freq(end) = [];
%// Draw stuff
figure;
subplot(2,1,1);
plot(freq, mag);
xlabel('Frequency (Hz)');
ylabel('Magnitude');
subplot(2,1,2);
plot(freq, pha);
xlabel('Frequency (Hz)');
ylabel('Phase (radians)');
这是我得到的：
如您所见，三个尖峰对应于三个正弦分量。你想要中间那个，那是159赫兹。要创建带通滤波器，您需要滤除除 +/- 159 Hz 的分量之外的所有分量。如果您想自动执行此操作，您将找到最接近 +/- 159 Hz 的 bin 位置，然后扩展围绕这两个点的邻域，并确保它们不受影响，同时将其余组件归零。
因为您有精确的正弦曲线，所以在这方面使用带通滤波器是完全可以接受的。通常，由于振铃和混叠效应，您不会这样做，因为带通滤波器的截止锐度会以这种方式在时域中引入不需要的混叠效应。见Wikipedia article on aliasing for more details。
因此，要找出我们需要过滤掉的位置，请尝试使用 min 并找出生成的频率与 159 Hz 之间的绝对差值 - 特别是找到位置。一旦你找到了 +159 Hz 和 -159 Hz 的这些点，在这些点周围扩展一个邻域，确保它们没有被触及，而其余点在频谱中设置为 0：
[~,min_pt_pos] = min(abs(freq - f(2))); %// Find location where +159 Hz is located
[~,min_pt_neg] = min(abs(freq + f(2))); %// Find location of where -159 Hz is
%// Neighbourhood size
ns = 100; %// Play with this parameter
%// Filtered signal
yfilt = zeros(1,numel(y));
%// Extract out the positive and negative frequencies centered at
%// 159 Hz
yfilt(min_pt_pos-ns/2 : min_pt_pos+ns/2) = ys(min_pt_pos-ns/2 : min_pt_pos+ns/2);
yfilt(min_pt_neg-ns/2 : min_pt_neg+ns/2) = ys(min_pt_neg-ns/2 : min_pt_neg+ns/2);
yfilt 现在包含过滤后的信号，去除了除 159 Hz 之外的所有分量。如果你想显示这个滤波信号的幅度和相位，我们可以这样做：
mag2 = abs(yfilt);
pha2 = phase(yfilt);
figure;
subplot(2,1,1);
plot(freq, mag2);
xlabel('Frequency (Hz)');
ylabel('Magnitude');
subplot(2,1,2);
plot(freq, pha2);
xlabel('Frequency (Hz)');
ylabel('Phase (radians)');
这是我们得到的：
正如我们所料，只有一个强频率会分解此信号，即 159 Hz。现在要重建此信号，您必须撤消我们所做的居中，并且您必须在此过滤结果上使用ifftshift，然后通过ifft 进行逆运算。您可能还会得到小的剩余虚部，因此最好在输出结果上使用real。
out = real(ifft(ifftshift(yfilt)));
如果我们绘制这个，我们得到：
plot(t, out);
xlabel('Time (seconds)');
ylabel('Height');
如您所见，有一个频率为 159 Hz 的正弦曲线。不过不要介意幅度。这仅仅是由于您选择绘制信号的点数不精确，因此某些时间点可能与信号的真实峰值不完全一致。请记住，如果存在多个正弦波，那么所有正弦波的峰值都会在某个点相遇，您将获得更高的峰值，而不是单个正弦波提供的峰值 1。因为幅度在 -1 到 1 之间徘徊，您可以确定只有一个正弦波存在。如果您要选择更细粒度的步长，从而在 FFT 中选择更多点，则可以避免这种情况。
希望这足以让您入门。祝你好运！