利用Excel进行FFT和Fourier的分析的基本步骤.docVIP

利用Excel进行FFT和Fourier分析的基本步骤 实例：杭州市2000人口分布密度 [根据2000年人口普查的街道数据经环带(rings)平均计算得到的结果，数据由冯健博士处理]。下面的变换实质是一种空间自相关的分析过程。 第一步，录入数据 在Excel中录入数据不赘述（见表1）。 表1 原始数据序列 表2 补充后的数据序列 第二步，补充数据 由于Fourier变换（FT）一般是借助快速Fourier变换（Fast Fourier Transformation, FFT）算法，而这种算法的技术过程涉及到对称处理，故数据序列的长度必须是2N（N=1,2,3,…，）。如果数据序列长度不是2N，就必须对数据进行补充或者裁减。现在数据长度是26，介于24=16到25=32之间，而26到32更近一些，如果裁减数据，就会损失许多信息。因此，采用补充数据的方式。 补充的方法非常简单，在数据序列后面加0，直到序列长度为32=25为止（表2）。当然，延续到64=26也可以，总之必须是2的整数倍。不过，补充的“虚拟数据”越多，变换结果的误差也就越大。 第三步，Fourier变换的选项设置 沿着工具（Tools）→数据分析（Data Analysis）的路径打开数据分析复选框（图1）。 图1 数据分析（Data Analysis）的路径 在数据分析选项框中选择傅立叶分析（Fourier Analysis）（图2）。 图2 数据分析（Data Analysis） 在Fourier分析对话框中进行如下设置：在输入区域中输入数据序列的单元格范围“$B$1:$B$33”；选中“标志位于第一行（L）”；将输出区域设为“$C$2”或者“$C$2:$C$33”（图3a）。 a b 图3 傅立叶分析（Fourier Analysis） 注意：如果“输入区域”设为“$B$2:$B$33”，则不选“标志位于第一行（L）”（图3b）。 表3 FFT的结果 第四步，输出FFT结果 选项设置完毕以后，确定（OK），立即得到FFT结果（表3）。 显然，表3给出的都是复数（complex numbers）。假定一个数据序列表为f(t)，则理论上Fourier变换的结果为 =F[f(t)], () 表3中给出的正是相应于F(ω)的复数，这里ω为角频率。 第五步，计算功率谱 Excel好像不能自动计算功率谱，这需要我们利用有关函数进行计算。计算公式为 式中A为复数的实部（real number），B为虚部（imaginary number），T为假设的周期长度，实则补充后的数据序列长度。对于本例，T=32。注意复数的平方乃是一个复数与其共轭（conjugate）复数的乘积，若F(ω)=a+bj，则|F(ω)|2=(a+bj)*(a-bj)=a2+b2。这样，根据表3中的FFT结果，我们有 其余依次类推。 显然，这样计算非常繁琐。一个简单的办法是调用Excel的模数（modulus）计算函数ImAbs，方法是在函数类别中找“其他”，在其他类中找“工程”类，在工程类中容易找到ImAbs函数（图4）。 确定以后，弹出一个选项框，选中第一个FFT结果，确定，得到218701.857（图5）。我们知道，复数的模数计算公式为 图4 模数计算函数 对于第一个FFT结果，由于虚部为0，模数就是其自身，即 但对于后面真正的复数，就不一样了。抓住第一个模数所在的单元格的右下角往下一拉，或者用鼠标双击该单元格的右下角，立即得到全部模数。 图5 计算模数 最后，用模数的2次方除以数据长度32立即得到全部功率谱密度结果（表4）。 表4 功率谱密度 下表是利用Mathcad2000计算的功率谱密度（表5）。利用Mathcad进行FFT，过程要简单得多，只要调用FFT命令，可以直接给出各种结果（包括图表）。但Mathcad的计算不求精度，有一定误差。将Mathcad的变换结果copy到Excel中进行比较，可以看到，如果不计误差，二者是一致的（表4）。 表5 借助Mathcad2000进行FFT的结果 第六步，功率谱分析 功率谱分析目前主要用于两个方面，一是侦测系统变化的某种周期或者节律，据此寻找因果关系（解释）或者进行某种发展预测（应用）；二是寻找周期以外的某些规律，据此对系统的时空结构特征进行解释。 表6 以对称点（f=0.5）为界，从完整的数据序列中截取一半 上面基于杭州人口密度数据的FFT，实际上是一种空间自相关分析过程，属于FT的第二类应用。这种过程不以寻找周期为目标，实际上也不存在任何周期。 不论目标是什么，都必须借助频谱图（频率－功率谱密度图）进行分析和解释。下面第一步就是绘制频谱图。首先要计算频率，线频或角频都可以，因为二者相差常数倍（2π