使用fft分析采样数据的周期性（matlab案例解析）

1 背景介绍

这是一个对matlab官方案例的分析笔记。

案例要解决的问题是：需要根据测得的办公楼温度数据，发现办公楼温度变化的周期性。

使用的方法：fft分析。

2 代码分析

案例代码

?

1 2 3 4 5 6 7 8 9 10 11 12 13 14

% 加载温度数据：tempratureData.mat中包含一个变量名为temp的5584*1数组（约16.5周，30min一个点的温度数据，16.5*7*24*2≈5584）load tempratureData.mat <br><br>Fs = 1/(60*30);   %  采样频率为1/(60*30) hz (即每30min采样一个点) t = (0:length(temp)-1)/Fs; % 时间   drawPlotFunction(t/(60*60*24*7),temp,...   'Time in weeks','Temperature (Fahrenheit)') %内部绘图函数   NFFT = length(temp);              % fft点数，即5584≈16.5*7*24*2 F = (0 : 1/NFFT : 1/2-1/NFFT)*Fs; % 频率向量，即最小分辨率为1/NFFT*Fs（即1/(16.5*7*24*60*60) hz, 以16.5周作为一个周期），最大分辨率约为1/2*Fs（即2*0.5h=1h作为一个周期）   TEMP = fft(temp,NFFT); % FFT处理后得到的数据，为5584*1复数数组 TEMP(1) = 0; % 清除直流分量方便可视化   drawPlotFunction(F*60*60*24*7,abs(TEMP(1:NFFT/2)),...   'Frequency (cycles/week)','Magnitude',[],[],[0 10])

最后一行绘图代码解析：

1. “F*60*60*24*7”解释：

　　对于  “1/60*60*24*7 hz”(原来"每1周一次"对应的频率)   乘系数(*60*60*24*7)转换后变为  “1 hz”(新的"每1周一次"对应的频率)  。F中的其它频率 乘系数后 也进行了线性缩放，可以把“每1周一次对应的频率（1hz）”看作是一个对应实际周期值的参照点，按比例推算其它频率值对应的实际周期值是多少。

2. “abs(TEMP(1:NFFT/2)”解释：
　　1）fft变换后，TEMP为复数数组，abs(TEMP)：将 复数数据 转换为 幅值数据 ，可用于进一步幅值绘图。
　　2）TEMP(1:NFFT/2)：取TEMP前一半的数据（5584/2），这是因为FFT变换的结果前一半数据与后一半数据幅值相同（共轭），且上面的频率向量F也是只取了5584/2个数据。（观察数据可以发现TEMP的第2792与2794个数据、第2791与2795个数据...分别互为共轭（实部相同，虚部相差一负号），绘制abs(TEMP)的图像可以发现它是关于中心对称的。）

绘图结果：

 绘图结果分析：

可以看到频谱图在Frequency f1=1和f2=7处具有明显峰值：f1=1对应每1周一次，即周期T1为7天；f2=7对应f1的7倍，则周期T2=1/7*T1，也即周期T2为1天。说明办公楼温度在周和日尺度下的变化具有明显的周期性。 （后续可进一步绘制一下“周期-幅值”的图像，参考下面官方给出的另一个例子。）

3 延伸阅读

另一个官方例子：

太阳黑子周期预测。

如何理解FT、DFT、FFT：

FT：  纯干货数学推导_傅里叶级数与傅里叶变换

DFT：离散傅里叶变换零基础入门-中文1（针对工科生，无需连续傅立叶变换知识）

FFT：快速傅里叶变换(FFT)——有史以来最巧妙的算法？