FFT修改

FFT-dit的算法讨论

级的概念

M=log2NM=log_2NM=log2N,所以N点DFT可分成M级。

下图为8点FFT时间抽取算法信号流程图，方框分别为m=0,m=1,m=2;

码位倒置

交换后输出序列x(k)依照正序排列，但输入序列x(n)的次序不再是原来的自然次序，这是由于将x(n)按奇,偶分开产生的。

雷德算法：对于自然顺序（二进制）我们是在低位加 1 得到下一位数，对于倒位序我们是在高位加 1 向低位进位。比如已知一个倒位序数是J求其下一个倒位序数，N位总数，把J与N/2比较若J<N/2则J的最高位为 0 ，把最高位置 1 ，就得到了J的下一个倒位序数；若J>=N/2则说明J的最高位为1 ，把最高位置0 ，比较次高位，若次高位为0 ，则把次高位置1，得到J的下一个倒位序，若次高位为1 ，则把次高位置0,以此类推…

参考：https://blog.csdn.net/corcplusplusorjava/article/details/17119567

Rader参考程序：

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
int main(void)
{int array[8]={0,1,2,3,4,5,6,7};int i,j,k;int N = 8;int temp;j = 0;for(i = 0; i < N -1; i ++){if(i < j){temp = array[i];array[i] = array[j];array[j] = temp;}k = N >> 1;while( k <= j){j = j - k;k >>= 1;}j = j + k;}for( i = 0; i < N; i ++)printf("%d %d\n",i,array[i]);printf("\n");return 0;
}

结果：

WrW^rWr因子的分布

规律：第m级，W2m+1rW_{2^{m+1}}^rW2m+1r

在fft中，乘法主要来自旋转因子。因为Wr=cos(2πrN)−jsin(2πrN),所以在对WrW^r=cos(\frac{2{\pi}r}{N})-jsin(\frac{2{\pi}r}{N}),所以在对W^rWr=cos(N2πr)−jsin(N2πr),所以在对Wr相乘时，必须产生相应的正，余弦函数。编程时，产生正，余弦函数两个方法，一个是在每一步直接产生，二是在程序开始前预先计算W^r存于一个数组中，等效一个正，余弦函数的“表”。

蝶形单元

第m级，有

Xm+1(p)=xm(p)+WNrXm(q)X_{m+1}(p)=x_{m}(p)+W_N^{r}X_m(q)Xm+1(p)=xm(p)+WNrXm(q)

Xm+1(q)=xm(p)−WNrXm(q)X_{m+1}(q)=x_{m}(p)-W_N^{r}X_m(q)Xm+1(q)=xm(p)−WNrXm(q)

"组"的概念

每一级的N/2个蝶形单元可以分成若干组，每一组有着相同的结构及W^r因子分布。

第m级的组数是N/2m+1N/2^{m+1}N/2m+1.

程序

A.1<<i是把1左移i位，每次左移以为就是乘以2，所以1<<i的结果是1乘以2的i次方 i<<1就是把i左移一位，即i乘以2，假如i=5,最后结果就是5*2=10

B.FFT蝶形算法使用三重循环，下一层的数据都是由上一层计算得到的。

#include <stdio.h>
#include<stdlib.h>
#define pi 3.1415
#define N  64
typedef struct
{double real,imag;
} complex; //复数结构
complex fft_outa[100];//a的全部数据
complex fft_outao[100];//a的单个数据
complex fft_outb[100];//b的全部数据
complex fft_outbo[100];//b的单个数据
double x_r[N],x_i[N];
double y_r[N],y_[N];void Rader(double x_r[N],double x_i[N])//雷德算法排序
{int i,j,k;j=0;double t_r,t_i;for(i=0; i<N-1; i++){if(i<j){t_r=x_r[i];t_i=x_i[i];x_r[i]=x_r[j];x_i[i]=x_i[j];x_r[j]=t_r;x_i[j]=t_i;}k=N>>1;while(k<=j){j=j-k;k>>=1;}j=j+k;}
}
void dit(double x_r[N],double x_i[N])//三层循环
{int m,r,p,q;int step,factor_step;int a,b,i;double t_real,t_imag;double w_re,w_im;for(m=0; m<log2(N); m++) //第一层循环{a=1;for(i=0;i<m+1;i++)a=a*2;b=a/2;w_re=cos(pi/b);w_im=-sin(pi/b);step=1<<(m+1);// factor_step=N>>(m+1);//旋转因素变化速度//初始化旋转因子double factor_real=1.0;double factor_imag=0.0;for(r=0; r<step/2; r++){//   w_re=cos(2*pi*(r+1)*factor_step/N);//   w_im=-sin(2*pi*(r+1)*factor_step/N);for(p=r; p<N; p+=step){q=p+step/2;//蝶形运算的两个输入t_real=x_r[q]*factor_real-x_i[q]*factor_imag;t_imag=x_r[q]*factor_imag+x_i[q]*factor_real;x_r[q]=x_r[p]-t_real;x_i[q]=x_i[p]-t_imag;x_r[p]=x_r[p]+t_real;x_i[p]=x_i[p]+t_imag;}t_real=factor_real*w_re-factor_imag*w_im;t_imag=factor_real*w_im+factor_imag*w_re;factor_real=t_real;factor_imag=t_imag;}}}int main()
{int n,i;double y[N];for(n=0; n<N; n++)//输入信号{x_r[n]=cos(n*pi/6);x_i[n]=0;}Rader(x_r,x_i);dit(x_r,x_i);for(i=0; i<N; i++)//输出fft{y[i]=sqrt(x_r[i]*x_r[i]+x_i[i]*x_i[i]);printf("%d %f\n",i,y[i]);}
}

c结果：

matlab的程序结果：

N=64;
n=0:1:N;
x=cos(n*pi/6);
X=fft(x);
figure
stem(n,abs(X))

对比结果还是有差异。

FFT过程中，对C理解还不够，对蝶形运算算法还不理解，对三重循环还不熟悉。

修改了，matlab的程序

>> n=0:1:63;
>> x=cos(n*pi/6);
>> y=fft(x);
>> figure
>> stem(n,abs(y))

得出和C一样的结果

FFT C语言修改了matlab相关推荐

快速傅里叶变换(FFT)c语言实现
快速傅里叶变换(FFT)c语言实现:(参考:FFT多种编程语言实现).注意:输入数据个数必须为2的n次方,数据不够可以用0补齐. #include <stdio.h> #include ...
R语言修改dataframe的列名（column name）实战
R语言修改dataframe的列名(column name)实战目录 R语言修改dataframe的列名(column name)实战 #使用colnames函数修改dataframe的列名
c语言一行代码太长,C语言修改一行代码，运行效率居然提升数倍，这个技巧你知道吗...
对编译.链接.OS内核.系统调优等技术感兴趣的童鞋,不妨右上角关注一下吧,近期会持续更新相关方面的专题文章!引言近日,网上看到一篇文章,分析数组访问的性能问题.文章经过一系列"有理有据&q ...
Win11系统语言修改不了中文怎么办
一些升级了Win11系统的朋友发现升级后发现是英文版的,怎么把英文版的换成中文版的呢?下面为大家带来如何把Win11系统语言从英文变成中文,方法非常简单. Win11系统语言修改不了中文怎么办 1.首 ...
c语言time.h时区不对,用C语言修改系统时区，发现一堆问题，请各位大侠不吝赐教。...
用C语言修改系统时区,发现一堆问题,请各位大侠不吝赐教. (2012-06-13 03:14:10) 标签: 系统 c语言杂谈用C语言修改系统时区,发现一堆问题,请各位大侠不吝赐教.已经实现,用s ...
易语言修改IP和DNS
易语言修改IP和DNS 修改IP和DNS的时候,有一个参数叫"连接名称", 不清楚这个参数导致走了不少弯路 "连接名称"对应的是"网络和共享中心&qu ...
闽江学院c语言期末试卷,Matlab期末复习08_闽江学院：matlab6.5（周赢武）_ppt_大学课件预览_高等教育资讯网...
Matlab期末复习 2008.06.04 第 1章 MATLAB语言概述第 2章基本语法第 4章 Matlab的其它函数库第 6章 Matlab在信号与系统中的应用第 9章 Matlab工 ...
c语言实现1024点fft程序,数字信号处理的步骤与注意事项,并编写1024个采样点的FFT C语言程序...
数字信号处理的步骤与注意事项,并编写1024个采样点的FFT C语言程序 1. 数字信号处理 1.1 数字信号处理概述数字信号处理是研究如何用数字或符号序列来表示信号以及如何对这些序列进行处理的一门 ...
linux c语言修改文件的时间属性,请教一个关于用标准C语言修改文件创建时间、修改时间和访问时间的问题。...
请教一个关于用标准C语言修改文件创建时间.修改时间和访问时间的问题. 标准C里面有没有这种方法呢?我需要在Unix下运行,具体的是Solaris 9下,把某个目录下的所有文件的这三个时间属性都改成某年 ...

FFT C语言修改了matlab

目录

FFT修改

FFT-dit的算法讨论

程序

c结果：

matlab的程序结果：

修改了，matlab的程序

FFT C语言修改了matlab相关推荐

最新文章

热门文章

FFT C语言 修改了matlab

目录

FFT修改

FFT-dit的算法讨论

程序

c结果：

matlab的程序结果：

修改了，matlab的程序

FFT C语言 修改了matlab相关推荐

最新文章

热门文章

FFT C语言修改了matlab

FFT C语言修改了matlab相关推荐