傅里叶变换与FFT应用

一、傅里叶变换

1.1 变换

      我们先给例子，假设在直角坐标系上有A(2,1),B(1,2);数和图之间存在的关系，称作变换；在图上我们想找对角线C，通过计算我们就知道C(3,3)；我们知道，在坐标系上有单位向量，ex与ey，ex·ex=ey·ey=1，它们的内积是1，ex·ey=0，我们在空间中找到了一组内积是0，自身长度为1，这种向量是需要的标准正交积。在平面内的任何向量都可以用这种向量来组合

1.2 傅里叶级数

     任何一个周期性的函数f(t)，都可以变换成一系列正余弦函数的和

      从实际现象上看，我们只能看到 时间t和幅值y所显示出来的一个图，这就是时域上的观察；但是如果我们从纵向深度上看，就可以看到这个信号被拆成了很多个信号，这就是频域。

      我们把纵向轴重新画一遍，另新坐标x轴为频率w，纵轴为其傅里叶变换F(t)，我们就可以在轴上获得到其信号，每个对应的采样频率都相当在图上的一个冲击信号。

      如上图，相当于有1 3 5 7 9 11 13 15频率组合成的一个信号，而纵轴就是在组合信号中每个信号所占的强度，这样我们就可以表达出一个信号下每个信号的分量，这就是频域图，也就是一个傅里叶变换。但是我们漏了一个东西，也就是每个信号都有个起始点，也就是信号中的相位差；所以频域上也可以拓展成一个3维图，y轴设置为相位。

      所以整体来说，我们对一个信号做傅里叶变换，就可以得到三个内容，一个是一系列的频率，第二个是每个频率的振幅是多大，第三个是每个频率的相位，变化的过程是相互的，可以进行反变换。

      因此就可以推出一个公式

1.3 傅里叶变换

        在圆上任意一点，我们都可以利用正交基表示，根据欧拉公式

       假设角度=wt，随着时间逆时针变换。所以每时每刻都能代表一种正交基。所以对于任意一个信号，我们可以将其看作一种周期性无穷大的信号

      我们做一个摘，假定一种频率w，让他去乘以这个信号，如果包含这个频率，那他就会有对应的值，如果不包含，就都是0，这样就可以摘出对应的频率w。

二、 FFT在MATLAB下的应用

参考内容为MATLAB的HELP

2.1 FFT的传统使用

      上面的推论，我们明白傅里叶变换其实就是一个数学公式，用于将按时间或空间采样的信号变换为按时序或者空间频率采样信号。在信号处理中，傅里叶变换可以揭示信号的重要特征。对于包含n个均匀采样点的向量x，其傅里叶变换定义为

      w为n个复单位根之一，我们说的这个复单位根，其实就是前面说的基向量的组合，在圆上的表达。在MATLAB里面的FFT函数使用快速傅里叶变换来计算数据的傅里叶变换。我们举个例子，以正弦信号x为例子，信号为时间t的函数，假设频率分量为15hz和20hz，使用10s周期内以1/50秒为增量进行采样的时间向量。

我们得到初始信号图如下：

      根据计算信号的傅里叶变换，并在频率空间创建对应于信号采样的向量f，我们利用fft函数用于计算输入信号x的dft，x是时域，y是对应频域(复数数组)，每个元素代表一个频率分量幅度和相位。前面假设的Ts为采样周期，fs为采样频率(1/Ts);  

      f=(0:length(y)-1)*fs/length(y), length(y)返回信号长度，即频域信号点数，0:length(y)-1 表示创建一个从0 到 length(y)-1的向量，fs/length(y)为计算频率分辨率，也就是上面所说的每个点对应的频率间隔，所得到的f就是频域向量。我们给出代码与图：

得到的图如下：

      绘图部分 plot(f,abs(y))绘制频谱图，'f'是频率向量，'abs(y)'为频域信号'y'的幅度（绝对值），计算频谱图时通常是对称的，上图我们看到右半其实就是左半的镜像，看左半，我们发现和前面给出的公式其实是一样的，在组合信号里，FT把存在的频率w给提取出来了，这里我们看到的就是频率15 20和它们对应的幅值。为了直观看出，我给出图像与对应的代码，为了区分出幅值的变换，我们给20hz的信号的幅值为2

% 设置参数
Ts = 1/50;  % 采样周期 (每秒50个采样点)
t = 0:Ts:10-Ts;  % 采样时间点
x = sin(2*pi*15*t) + 2*sin(2*pi*20*t);  % 信号源 (15 Hz 和 20 Hz 的正弦波)

% 绘制时间域信号
figure;
plot(t, x);
xlabel('时间 (秒)');
ylabel('幅值');
title('时间域信号');

% 计算FFT
y = fft(x);  % 快速傅里叶变换
n = length(y);  % FFT结果的长度
fs = 1/Ts;  % 采样频率

% 只保留前半部分的频率分量
half_n = ceil(n/2);  % 保留前半部分的索引范围
f = (0:half_n-1)*(fs/n);  % 生成频率向量
y_half = y(1:half_n);  % 截取前半部分的FFT结果

% 绘制频谱图
figure;
plot(f, abs(y_half));
xlabel("频率 (Hz)");
ylabel("幅度");
title('Magnitude');

得到结果如下：

2.2 带噪声的信号

      正常情况下，获取的信号肯定是带有随机噪声破坏的，而FFT可以清除随机噪声并显示频率(当然你的电路上最好能自己滤除大部分噪声)。再给个例子，再原始信号里我们加入一个噪声。rng("default")重置随机生成器种子，其实就是生成一个随机数而已。 xnoise=x+2.5*rand(size(t)

      生成标准差为2.5，且与t具有相同数量的相同大小的随机数数组，这样生成了一个无任何规律的分量，我们将这个分量添加到x内，得到带噪声的信号xnoise。当然，我认为你可以把这个噪声理解成直流分量，这样我们可以得到以下图：

为了更直观我全部给出，可以看到噪声已经把原始信号几乎破坏掉了，但是FFT还是能提取出里面频率占比最多的。这就是我们需要的效果。

2.3 正弦波相位

      使用傅里叶变换，我们之前也说了还可以得到第三种数值-原始信号的相位频谱。再给个例子，一个信号由15hz与40hz两个正弦波组合而成的信号。第一个相位是 -pi/4的余弦波，第二个相位是pi/2的正弦波，假定以100hz的频率对信号采样。 关于采样定理，一个周期至少采两个点，所以采集频率要是原始信号的2倍以上，但是ft的逆变换可以把原始信号还原，所以如果你的采样周期刚好会是信号的整数倍，其实是可以避免混叠现象的。

% 设置参数
Fs = 100;              % 采样频率 (Hz)
T = 1/Fs;              % 采样周期 (秒)
L = 1000;              % 信号长度 (采样点数)
t = (0:L-1)*T;         % 时间向量 (秒)

% 创建信号
x = cos(2*pi*15*t - pi/4) + cos(2*pi*40*t+pi/2);  % 组合信号

% 计算傅里叶变换
Y = fft(x);

% 双边频谱到单边频谱的转换
P2 = abs(Y/L);         % 双边频谱
P1 = P2(1:L/2+1);      % 单边幅值谱
P1(2:end-1) = 2*P1(2:end-1);  % 将非直流分量和奈奎斯特频率的幅值乘以2

% 生成频率向量
f = Fs*(0:(L/2))/L;

% 计算相位
tol = 1e-6;                        % 相位计算时的容差
Y(abs(Y) < tol) = 0;               % 将幅值低于容差的部分设置为0
theta = angle(Y);                  % 计算相位
theta1 = theta(1:L/2+1);           % 单边相位谱

% 绘制幅值谱
figure;
subplot(2,1,1);
stem(f, P1, 'filled');
title('单边幅值谱');
xlabel('频率 (Hz)');
ylabel('幅值');

% 绘制相位谱
subplot(2,1,2);
stem(f, theta1/pi, 'filled');
title('单边相位谱');
xlabel('频率 (Hz)');
ylabel('相位 / \pi');

得到效果如图：

2.4 噪声效果下的相位频谱图

      利用带通滤波我们可以去掉很多噪声，但是相位无法像幅值那样突出，只要被检测出的相位存在，他就一定会存在分量，并且这种分量和幅值不同，他不存在占比的问题，只要有，那就会干扰，因此这里就对需要的相位进行标出，因为对于的频率肯定只有一个对应相位。

小白一枚，如有理解错误，请大家批评指正。