1. 全栈程序员必看首页

数字信号处理公式变程序（一）——DFT、FFT

全栈程序员-站长 • • 未分类 • 阅读 2

数字信号处理公式变程序（一）——DFT、FFT之前搞了一些数字信号处理算法编程 OC 一直没来得及整理 现在整理一下 陆续会更新 包括 FFT 巴特沃斯滤波器 高低带通 高低带阻 数据差值 线性 sinc 三次样条 数据压缩 等距 平均 峰值检测 和模仿 matlab 的 STFT 功能 spectrogram 函数三维绘图 注 可能会有不足或者理解偏差的地方 路过的高人请不吝赐教 好啦 进入正题

之前搞了一些数字信号处理算法编程（OC），一直没来得及整理，现在整理一下。陆续会更新，包括FFT、巴特沃斯滤波器（高低带通、高低带阻）、数据差值（线性、sinc、三次样条*）、数据压缩（等距、平均、峰值检测）和模仿matlab的STFT功能（spectrogram函数三维绘图）。

 

注：可能会有不足或者理解偏差的地方，路过的高人请不吝赐教。

 

好啦，进入正题。

—————————————————————————————

在数字世界中FFT的意义不言而喻（我曾转载一篇文章有提到：http://blog.csdn.net/shengzhadon/article/details/），这里就不再赘述了。FFT（快速傅里叶变换）是DFT的一种特殊情况，就是当运算点的个数是2的整数次幂的时候进行的运算（不够用0补齐），那就先从DFT开始吧。

 

一、DFT（本部分就是翻译公式）

定义（来自百科）：离散傅里叶变换（Discrete Fourier Transform，缩写为DFT），是傅里叶变换在时域和频域上都呈离散的形式，将信号的时域采样变换为其DTFT的频域采样。在形式上，变换两端（时域和频域上）的序列是有限长的，而实际上这两组序列都应当被认为是离散周期信号的主值序列。即使对有限长的离散信号作DFT，也应当将其看作其周期延拓的变换。在实际应用中通常采用快速傅里叶变换计算DFT。

 

DFT的公式是：设x(n)为M点的有限长序列，即在0≤n≤M-1内有值，则定义x(n)的N点（N≥M。当N＞M时，补N-M个零值点），离散傅里叶变化定义为

数字信号处理公式变程序（一）——DFT、FFT

其中，数字信号处理公式变程序（一）——DFT、FFT为旋转因子（为方便编辑后续记作W(N, nk)，特此说明），其计算公式为数字信号处理公式变程序（一）——DFT、FFT，具有以下性质：

①周期性：W(N, nk) = W(N, (n+rN)k) = W(N, n(k+rN))，其中r问整数。

②共轭对称性：(W(N, nk))* = W(N, -nk)。

③可约性：W(N, nk) = W(mN, mnk)，W(N, nk) = W(N/m, nk/m)，其中m为整数，N/m为整数。

④特殊值：

W(N, N/2) = -1；  W(N, (k+N/2)) = -W(N, k)；  W(N, (N-k)n) = W(N, (N-n)k) = W(N, -nk)。

所以，在计算旋转因子的过程中可以适当的使用特殊值来提高运算的效率。

 

在计算DFT时，如果数据的点数够计算的点数，则截取计算点数长的数据进行DFT运算，否则将数据点个数补0至计算点个数，然后进行计算，举例如下（举例只是为了说明问题，没有按照编程语言的书写格式）。

比如：数据点数组为 dataArray = {1， 2， 3， 4，5}，可以看出数据长度为5。

   如果要求做4点DFT运算，则只需截取前4个数作为运算数组进行运算即可，即为{1， 2， 3， 4}；

   如果要求做8点DFT运算，则需在原数组后补三个0，使长度为8后再进行计算，计算数组为{1， 2， 3， 4，5，0，0，0}。

 

数字信号处理公式变程序（一）——DFT、FFT

 

二、FFT

 

 

1.DFT的运算量

复数乘法次数为N*N，复数加法次数为N(N-1)，若N>>1，则这两者都近似为N*N，它随N增大为急速增大。

改进途径：①利用旋转因子性质减小计算量；②由于运算量和N*N成正比，因而可将N点DFT分解成小点数的DFT，以减少运算量（点数越小，计算量越小）；③改为用FFT计算，复数乘法次数为(N/2)*log(2,N)。

 

2.FFT可分为按照时间抽选的基-2算法（库利-图基算法DIT-FFT）和按频率抽选的基-2算法（桑德-图基算法DIF-FFT）。本文采用DIT-FFT算法。

 

3.FFT计算原理及流程图

FFT的计算要求点数必须为2的整数次幂，如果点数不够用0补齐。例如计算{2，3，5，8，4}的8点FFT，需要补3个0后进行计算，如果计算该数组的5点FFT，则先计算8点FFT后截取前5个值即可（不提倡）。

 

（1）原理

公式推导

设序列x(n)点数为N=2^L，L为整数。将N=2^L，先按照n的奇偶分为两组，其中r = 0, 1, …, N/2-1

x(2r) = x1(r)

x(2r-1) = x2(r)

则可将DFT化为（式1）

数字信号处理公式变程序（一）——DFT、FFT

由上式可以看出一个N点的DFT可以分为两个N/2点的DFT，按照上式右组合成N点DFT。但是这里的x1(r)、x2(r)以及X1(k)、X2(k)都是N/2点的序列，即r,k满足r,k=0, 1, …, N/2-1。而X(k)却有N点，上式计算的只是X(k)的前半项数的结果，因此还需要计算后半项的值。

数字信号处理公式变程序（一）——DFT、FFT

将①②③代入（式1），就可将X(k)表达为前后两部分：

前半部分，X(k)，当k=0, 1, …, N/2-1

数字信号处理公式变程序（一）——DFT、FFT

后半部分，X(k)，当k=N/2, …, N-1

数字信号处理公式变程序（一）——DFT、FFT

蝶形运算说明

蝶形运算，符号表示如下图所示：

数字信号处理公式变程序（一）——DFT、FFT

所以，FFT的蝶形运算图表示为（8点，2^3 = 8，运算级数最大为L=3）

数字信号处理公式变程序（一）——DFT、FFT

在蝶形运算中变化规律由W(N, p)推导，其中N为FFT计算点数，J为下角标的值

L = 1时，W(N, p) = W(N, J) = W(2^L, J)，其中J = 0;

L = 2时，W(N, p) = W(N, J) = W(2^L, J)，其中J = 0, 1;

L = 3时，W(N, p) = W(N, J) = W(2^L, J)，其中J = 0, 1, 2, 3;

所以，W(N, p) = W(2^L, J)，其中J = 0, 1, …, 2^(L-1)-1

又因为2^L = 2^M*2^(L-M) = N*2^(L-M)，这里N为2的整数次幂，即N=2^M，

W(N, p) = W(2^L, J) = W(N*2^(L-M), J) = W(N, J*2^(M-L))

所以，p = J*2^(M-L)，此处J = 0, 1, …, 2^(L-1)-1，当J遍历结束但计算点数不够N时，J=J+2^L，后继续遍历，直到计算点数为N时不再循环。

举例：N=8点的FFT计算

当L=2时，J = 0, 1两个值，因此p = J*2^(M-L) = 0, 2两个值，即旋转因子有两个值W(8, 0)和W(8, 2)，计算中两行之间的距离B = 2^(L-1)=2^(2-1)=2，

数字信号处理公式变程序（一）——DFT、FFT

代入J=0, 1可求得X(0)、X(0+2)和X(1)、X(1+2)，即可求出第2级蝶形运算的X(0), X(1), X(2), X(3)，也就是求出一半，此时J加步进2^L=4，即J=J+2^L=4, 5，

再代入J=4, 5可求出X(4), X(4+2)和X(5), X(5+2)，即可求出第2级蝶形运算的X(4), X(5), X(6), X(7)，已经全部求出，J循环结束。

 

二进制倒序说明

前面数的排列顺序是进行二进制倒序后的排序。二进制倒序是指将某数转化为二进制表示，将最高位看做最低位、次高位看做次低位…以此类推，计算后的纸进行排序。例如3的二进制表示为011b（3=0*2^2+1*2+1），二进制倒序值为6（011b，最高位看做最低位…即6=0+1*2+1*2^2）

即0到7的二进制排序是：

       0,                4,                 2,                6,                 1,                5,                3,                7

000→000    001→100    010→010    011→110    100→001    101→101    110→011    111→111

 

（2）流程图

①FFT运算总流程图，包括

数字信号处理公式变程序（一）——DFT、FFT

 

数字信号处理公式变程序（一）——DFT、FFT

说明：蝶形运算中又三层循环

第一层（最外层），完成M次迭代过程，即算出A0(k), A1(k), …, Am(k)，其中k=0, 1, …, N；A2(k)为蝶形运算第2级的结果，如A0(k)=x(k), Am(k)=X(k)；

第二层（中间层），完成旋转因子的变化，步进为2^L；

第三层（最里层），完成相同旋转因子的蝶形运算。

 

 

 

 

 

本文开头就已经说了，我个人比较喜欢使用matlab，所以我程序的对比对象当然是matlab。给定t1[8] = {
0,1,2,3,4,5,6,7};signal=sin(t1);对signal进行8、16、4、13、6点FFT/DFT计算输出的结果和给定128点数据进行128点FFT计算结果如下图所示：

数字信号处理公式变程序（一）——DFT、FFT

 

复制一遍：

FFT，8点计算结果

0. + 0.000000i   2. – 2.097012i   -1. + 0.i   -0. + 0.i   -0. + 0.000000i   -0. – 0.i   -1. – 0.i   2. + 2.097012i   

FFT，16点计算结果

0. + 0.000000i   1. + 0.i   2. – 2.097012i   -1. – 2.i   -1. + 0.i   0. – 0.i   -0. + 0.i   0. – 0.i   -0. + 0.000000i   0. + 0.i   -0. – 0.i   0. + 0.i   -1. – 0.i   -1. + 2.i   2. + 2.097012i   1. – 0.i   

FFT，4点计算结果

1. + 0.000000i   -0. – 0.i   -0.073294 + 0.000000i   -0. + 0.i   

FFT，13点计算结果

0. + 0.000000i   1. + 0.i   0. – 3.i   -2. + 0.i   0. – 0.i   -0. + 0.i   -0. – 0.i   -0. + 0.i   -0. – 0.i   0. + 0.i   -2. – 0.i   0. + 3.i   1. – 0.i   

FFT，6点计算结果

0. + 0.000000i   -0. – 3.002073i   0. – 0.i   0. + 0.000000i   0. + 0.i   -0. + 3.002073i   

FFT，128点计算结果

0. + 0.000000i   0. + 0.010352i   0. + 0.021364i   0. + 0.033697i   0. + 0.047746i   0. + 0.061398i   0. + 0.050647i   -6. – 0.i   54. + 5.062080i   -57. – 6.i   9. + 1.i   0. + 0.093026i   -0. – 0.020889i   -0. – 0.043726i   -0. – 0.047995i   -0. – 0.047167i   -0. – 0.044825i   -0. – 0.042109i   -0. – 0.039423i   -0. – 0.036897i   -0. – 0.034570i   -0. – 0.032440i   -0. – 0.030495i   -0. – 0.028715i   -0. – 0.027084i   -0.092007 – 0.025581i   -0.084229 – 0.024193i   -0.077390 – 0.022905i   -0.071342 – 0.021706i   -0.065971 – 0.020587i   -0.061182 – 0.019537i   -0.056895 – 0.018550i   -0.053044 – 0.017619i   -0.049575 – 0.016737i   -0.046441 – 0.015901i   -0.043601 – 0.015106i   -0.041024 – 0.014347i   -0.038679 – 0.013621i   -0.036542 – 0.012927i   -0.034593 – 0.012260i   -0.032808 – 0.011617i   -0.031181 – 0.010997i   -0.029686 – 0.010398i   -0.028318 – 0.009818i   -0.027066 – 0.009256i   -0.025917 – 0.008710i   -0.024865 – 0.008177i   -0.023903 – 0.007660i   -0.023024 – 0.007153i   -0.022220 – 0.006657i   -0.021489 – 0.006172i   -0.020825 – 0.005695i   -0.020224 – 0.005227i   -0.019682 – 0.004767i   -0.019197 – 0.004313i   -0.018766 – 0.003867i   -0.018383 – 0.003419i   -0.018055 – 0.002986i   -0.017772 – 0.002551i   -0.017534 – 0.002121i   -0.017342 – 0.001693i   -0.017193 – 0.001268i   -0.017087 – 0.000845i   -0.017025 – 0.000423i   -0.017003 + 0.000000i   -0.017025 + 0.000422i   -0.017087 + 0.000845i   -0.017193 + 0.001268i   -0.017342 + 0.001693i   -0.017534 + 0.002121i   -0.017772 + 0.002551i   -0.018055 + 0.002986i   -0.018383 + 0.003425i   -0.018766 + 0.003866i   -0.019197 + 0.004314i   -0.019682 + 0.004767i   -0.020224 + 0.005227i   -0.020825 + 0.005695i   -0.021489 + 0.006172i   -0.022220 + 0.006656i   -0.023024 + 0.007153i   -0.023902 + 0.007659i   -0.024865 + 0.008177i   -0.025917 + 0.008710i   -0.027066 + 0.009256i   -0.028318 + 0.009818i   -0.029686 + 0.010398i   -0.031179 + 0.010998i   -0.032808 + 0.011617i   -0.034592 + 0.012259i   -0.036542 + 0.012927i   -0.038679 + 0.013621i   -0.041024 + 0.014347i   -0.043601 + 0.015106i   -0.046441 + 0.015901i   -0.049574 + 0.016737i   -0.053044 + 0.017619i   -0.056894 + 0.018549i   -0.061182 + 0.019537i   -0.065971 + 0.020587i   -0.071342 + 0.021706i   -0.077390 + 0.022905i   -0.084229 + 0.024193i   -0.092005 + 0.025579i   -0. + 0.027084i   -0. + 0.028715i   -0. + 0.030496i   -0. + 0.032441i   -0. + 0.034570i   -0. + 0.036897i   -0. + 0.039423i   -0. + 0.042109i   -0. + 0.044825i   -0. + 0.047166i   -0. + 0.047996i   -0. + 0.043726i   -0. + 0.020889i   0. – 0.093026i   9. – 1.i   -57. + 6.i   54. – 5.062086i   -6. + 0.i   0. – 0.050648i   0. – 0.061398i   0. – 0.047746i   0. – 0.033698i   0. – 0.021364i   0. – 0.010354i 

 

测试结果显示：iOS程序与MATLAB程序运算结果一致，但是MATLAB运算效率高。以下是多次比较求平均后的时间值

数字信号处理公式变程序（一）——DFT、FFT

 

为了更直观的观察，将其他的信号（几个信号叠加后的信号）通过fft计算后并作频域处理导入matlab，绘制图如下：

数字信号处理公式变程序（一）——DFT、FFT

数字信号处理公式变程序（一）——DFT、FFT

 

===================================================================

OK，FFT随笔一份，方便自己以后查阅。其中一定会有很多不足，望路过的大神订正。谢过

 

 

版权声明：本文内容由互联网用户自发贡献，该文观点仅代表作者本人。本站仅提供信息存储空间服务，不拥有所有权，不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容， 请联系我们举报，一经查实，本站将立刻删除。

发布者：全栈程序员-站长，转载请注明出处：https://javaforall.net/223468.html原文链接：https://javaforall.net

赞 (0)
0 0

关于作者

全栈程序员-站长的头像

全栈程序员-站长

133.5K 文章
3 粉丝
本网站汇聚当前互联网主流语音，持续更新，欢迎关注公众号“全栈程序员社区”
上一篇 2026年3月17日 下午2:08
下一篇 2026年3月17日 下午2:08


相关推荐

  • 高效的NoSQL数据库服务Amozon DynamoDB体验分享

    高效的NoSQL数据库服务Amozon DynamoDB体验分享

    高效的NoSQL数据库服务Amozon DynamoDB体验分享AmazonDynamo 是一种全托管 NoSQL 数据库服务 相当于在 NoSQL 的基础上做了很多扩展 附加了很多降本增效的功能 让开发人员不需要额外去维护什么 开箱即用 它主要提供的服务有 无缝扩展 快速可预测的性能 减少繁琐的管理分布式数据库的工作负担 提供静态加密 按需备份 自动维护过期项目等 极大的减少开发人员的时间成本和提高建设效率

    全栈程序员-站长的头像 全栈程序员-站长
    2026年3月19日
    3
  • 时间戳格式化「建议收藏」

    时间戳格式化「建议收藏」

    时间戳格式化「建议收藏」须知：1. 时间戳分2种，一种是10位的，只包含年月日时分秒，也就是说，只精确到秒。一种是13位的，包含毫秒。这2种都叫时间戳，并不是只有精确到毫秒的才叫时间戳。10位时间戳就是从1970-01-01到当前的秒数，注意，不是毫秒数，所以需要按毫秒解析时，要*100013位时间戳就是从1970-01-01到当前的毫秒数，在java中用Instant对象对应。2. timestamp的格式化串用大写的S来表示毫秒数。S的个数和毫秒的位数严格对应，否则报错。如果规范中要求精确到毫秒，那么给的时间字符串

    全栈程序员-站长的头像 全栈程序员-站长
    2022年4月19日
    474
  • 共享内存函数（shmget、shmat、shmdt、shmctl）及其范例

    共享内存函数（shmget、shmat、shmdt、shmctl）及其范例

    共享内存函数（shmget、shmat、shmdt、shmctl）及其范例共享内存函数由 shmget shmat shmdt shmctl 四个函数组成 下面的表格列出了这四个函数的函数原型及其具体说明 1 nbsp shmget 函数原型 shmget 得到一个共享内存标识符或创建一个共享内存对象 所需头文件 include includ

    全栈程序员-站长的头像 全栈程序员-站长
    2026年3月17日
    2
  • pycharm进入和退出控制台console

    pycharm进入和退出控制台console

    pycharm进入和退出控制台consolePyCharm 是有交互式界面的进入 PyCharmconso 在代码页面右击 选择如下然后就可以再 python 控制台进行编程了 和 matlab 很像哈哈 退出 PyCharmconso 适合控制台编写代码 但是不适合运行脚本代码 因为每次运行他都会跳出来一个界面 有时候不小心进入了 PyCharmconso 还很尴尬 退出步骤如下取消掉 runwithPycho 并勾选 Emulatatermi

    全栈程序员-站长的头像 全栈程序员-站长
    2026年3月27日
    2
  • 函数声明[通俗易懂]

    函数声明[通俗易懂]

    函数声明[通俗易懂]语法描述通过函数声明构造的函数是Function对象，所以拥有一切Function对象所有的属性，方法和行为。函数默认返回undefined，如果想返回其他值，函数必须使用return语句来返回

    全栈程序员-站长的头像 全栈程序员-站长
    2022年8月4日
    9
  • 半年学习安排

    半年学习安排

    半年学习安排

    全栈程序员-站长的头像 全栈程序员-站长
    2021年5月7日
    127

发表回复

您的邮箱地址不会被公开。 必填项已用 * 标注

关注全栈程序员社区公众号