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

1 背景介绍

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

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

使用的方法：fft分析。

2 代码分析

案例代码

% 加载温度数据：tempratureData.mat中包含一个变量名为temp的5584*1数组（约16.5周，30min一个点的温度数据，16.5*7*24*2≈5584）

load tempratureData.mat 

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天。说明办公楼温度在周和日尺度下的变化具有明显的周期性。（后续可进一步绘制一下“周期-幅值”的图像，参考下面官方给出的另一个例子。）