完善资料让更多小伙伴认识你,还能领取20积分哦, 立即完善>
C语言版FFT简单测试
本次我们来自己封装一个FFT函数,进行简单的测试。 fft.c #include "math.h" #include "fft.h" //精度0.0001弧度 //复数的交换 void conjugate_complex(int n,complex in[],complex out[]) { int i = 0; for(i=0;i out.imag = -in.imag; out.real = in.real; } } //求所有复数的模 void c_abs(complex f[],float out[],int n) { int i = 0; float t; for(i=0;i t = f.real * f.real + f.imag * f.imag; out = sqrt(t); } } //求复数的和 void c_plus(complex a,complex b,complex *c) { c->real = a.real + b.real; c->imag = a.imag + b.imag; } //求复数的差 void c_sub(complex a,complex b,complex *c) { c->real = a.real - b.real; c->imag = a.imag - b.imag; } //求复数的积 void c_mul(complex a,complex b,complex *c) { c->real = a.real * b.real - a.imag * b.imag; c->imag = a.real * b.imag + a.imag * b.real; } //求复数的商 void c_div(complex a,complex b,complex *c) { c->real = (a.real * b.real + a.imag * b.imag)/(b.real * b.real +b.imag * b.imag); c->imag = (a.imag * b.real - a.real * b.imag)/(b.real * b.real +b.imag * b.imag); } #define SWAP(a,b) tempr=(a);(a)=(b);(b)=tempr void Wn_i(int n,int i,complex *Wn,char flag) { Wn->real = cos(2*PI*i/n); if(flag == 1) Wn->imag = -sin(2*PI*i/n); else if(flag == 0) Wn->imag = -sin(2*PI*i/n); } //傅里叶变化 void fft(int N,complex f[]) { complex t,wn;//中间变量 int i,j,k,m,n,l,r,M; int la,lb,lc; /*----计算分解的级数M=log2(N)----*/ for(i=N,M=1;(i=i/2)!=1;M++); /*----按照倒位序重新排列原信号----*/ for(i=1,j=N/2;i<=N-2;i++) { if(i t=f[j]; f[j]=f; f=t; } k=N/2; while(k<=j) { j=j-k; k=k/2; } j=j+k; } /*----FFT算法----*/ for(m=1;m<=M;m++) { la=pow(2,m); //la=2^m代表第m级每个分组所含节点数 lb=la/2; //lb代表第m级每个分组所含碟形单元数 //同时它也表示每个碟形单元上下节点之间的距离 /*----碟形运算----*/ for(l=1;l<=lb;l++) { r=(l-1)*pow(2,M-m); for(n=l-1;n lc=n+lb; //n,lc分别代表一个碟形单元的上、下节点编号 Wn_i(N,r,&wn,1);//wn=Wnr c_mul(f[lc],wn,&t);//t = f[lc] * wn复数运算 c_sub(f[n],t,&(f[lc]));//f[lc] = f[n] - f[lc] * Wnr c_plus(f[n],t,&(f[n]));//f[n] = f[n] + f[lc] * Wnr } } } } //傅里叶逆变换 void ifft(int N,complex f[]) { int i=0; conjugate_complex(N,f,f); fft(N,f); conjugate_complex(N,f,f); for(i=0;i f.imag = (f.imag)/N; f.real = (f.real)/N; } } fft.H #ifndef __FFT_H__ #define __FFT_H__ typedef struct complex //复数类型 { float real; //实部 float imag; //虚部 }complex; #define PI 3.1415926535897932384626433832795028841971 /// void conjugate_complex(int n,complex in[],complex out[]); void c_plus(complex a,complex b,complex *c);//复数加 void c_mul(complex a,complex b,complex *c) ;//复数乘 void c_sub(complex a,complex b,complex *c); //复数减法 void c_div(complex a,complex b,complex *c); //复数除法 void fft(int N,complex f[]);//傅立叶变换 输出也存在数组f中 void ifft(int N,complex f[]); // 傅里叶逆变换 void c_abs(complex f[],float out[],int n);//复数数组取模 #endif 测试代码部分 #include "sys.h" #include "delay.h" #include "usart.h" #include "led.h" #include "math.h" #include "fft.h" #define N 256 //采样点数 #define Fs 44800 //采样频率 #define F 175 //分辨率 #define P1 30 //测试相位 #define P2 60 #define P3 90 //FFT测试数据集 输入数组 complex FFT_256PointIn[N]; //FFT测试数据集 输出数组 float FFT_256PointOut[N/2]; //填入数组 void InitBufInArray() { unsigned short i; for(i=0; i FFT_256PointIn.real = 1500 * sin(2*PI * i * 350.0 / Fs+(PI*P1/180)) +2700 * sin(2*PI * i * 8400.0 / Fs+(PI*P2/180)) +4000 * sin(2*PI * i * 18725.0 / Fs+(PI*P3/180)); FFT_256PointIn.imag = 0; } } /****************************************************************** 函数名称:GetPowerMag() 函数功能:计算各次谐波幅值 参数说明: 备 注:先将FFT_256PointIn分解成实部(X)和虚部(Y), 然后计算幅值:(sqrt(X*X+Y*Y)*2/N 然后计算相位:atan2(Y/X) 作 者:土耳其冰激凌 *******************************************************************/ void GetPowerMag() { unsigned short i; float X,Y,P,Mag; c_abs(FFT_256PointIn,FFT_256PointOut,N/2); for(i=0; i X = FFT_256PointIn.real/N; //计算实部 Y = FFT_256PointIn.imag/N; //计算虚部 Mag = FFT_256PointOut*2/N; //计算幅值 P = atan2(Y,X)*180/PI; //计算相位 printf("%d ",i); printf("%d ",F*i); printf("%f ",Mag); printf("%f ",P); printf("%f ",X); printf("%f rn",Y); } } int main(void) { delay_init(); //延时函数初始化 LED_Init(); //初始化与LED连接的硬件接口 uart_init(9600); //初始化串口 9600波特率 printf("这是一个FFT 测试实验rn"); InitBufInArray(); fft(N,FFT_256PointIn); printf("点数 频率 幅值 实部 虚部n"); GetPowerMag(); while(1) { LED0=0; LED1=1; delay_ms(300); //延时300ms LED0=1; LED1=0; delay_ms(300); //延时300ms } } 输出打印结果 点数 频率 幅值 相位 实部 虚部 0 0 0.000031 0.000000 0.000015 0.000000 1 175 0.000025 36.658634 0.000010 0.000008 2 350 1500.000000 -60.000000 374.999969 -649.519043 3 525 0.000020 52.221466 0.000006 0.000008 4 700 0.000014 -28.152960 0.000006 -0.000003 5 875 0.000060 -17.167744 0.000029 -0.000009 6 1050 0.000012 -4.051142 0.000006 -0.000000 7 1225 0.000028 27.618505 0.000013 0.000007 8 1400 0.000036 56.088974 0.000010 0.000015 9 1575 0.000013 85.033867 0.000001 0.000007 10 1750 0.000019 -126.592506 -0.000006 -0.000008 11 1925 0.000099 155.458847 -0.000045 0.000021 12 2100 0.000021 -158.704437 -0.000010 -0.000004 13 2275 0.000018 -89.424667 0.000000 -0.000009 14 2450 0.000015 -87.367409 0.000000 -0.000007 15 2625 0.000019 -84.178368 0.000001 -0.000010 16 2800 0.000086 112.500000 -0.000017 0.000040 17 2975 0.000008 -30.908798 0.000003 -0.000002 18 3150 0.000028 -93.026222 -0.000001 -0.000014 19 3325 0.000032 -100.869904 -0.000003 -0.000016 20 3500 0.000028 -35.795570 0.000012 -0.000008 21 3675 0.000248 -9.752425 0.000122 -0.000021 22 3850 0.000036 -59.215248 0.000009 -0.000015 23 4025 0.000010 -76.748405 0.000001 -0.000005 24 4200 0.000050 154.164047 -0.000023 0.000011 25 4375 0.000015 -60.317295 0.000004 -0.000006 26 4550 0.000014 108.588715 -0.000002 0.000006 27 4725 0.000029 98.175682 -0.000002 0.000014 28 4900 0.000016 -124.238319 -0.000004 -0.000007 29 5075 0.000048 108.763763 -0.000008 0.000023 30 5250 0.000047 -21.094650 0.000022 -0.000008 31 5425 0.000009 -38.687473 0.000003 -0.000003 32 5600 0.000113 -87.260597 0.000003 -0.000056 33 5775 0.000016 -63.487057 0.000004 -0.000007 34 5950 0.000022 7.387201 0.000011 0.000001 35 6125 0.000027 63.663502 0.000006 0.000012 36 6300 0.000032 62.835876 0.000007 0.000014 37 6475 0.000076 -38.793385 0.000030 -0.000024 38 6650 0.000009 16.024696 0.000004 0.000001 39 6825 0.000012 -110.747261 -0.000002 -0.000006 40 7000 0.000028 82.453056 0.000002 0.000014 41 7175 0.000016 -27.415680 0.000007 -0.000004 42 7350 0.000027 160.582352 -0.000013 0.000004 43 7525 0.000087 -87.681641 0.000002 -0.000044 44 7700 0.000022 -135.357925 -0.000008 -0.000008 45 7875 0.000014 -169.563293 -0.000007 -0.000001 46 8050 0.000009 1.227116 0.000004 0.000000 47 8225 0.000016 42.631676 0.000006 0.000006 48 8400 2700.000000 -30.000000 1169.134277 -675.000000 49 8575 0.000008 175.232651 -0.000004 0.000000 50 8750 0.000017 147.675629 -0.000007 0.000005 51 8925 0.000056 179.124176 -0.000028 0.000000 52 9100 0.000021 126.577728 -0.000006 0.000008 53 9275 0.000040 28.367682 0.000018 0.000010 54 9450 0.000011 -37.718578 0.000004 -0.000003 55 9625 0.000010 131.671509 -0.000003 0.000004 56 9800 0.000026 172.892471 -0.000013 0.000002 57 9975 0.000018 -143.289139 -0.000007 -0.000005 58 10150 0.000012 -115.207047 -0.000003 -0.000006 59 10325 0.000008 -87.947166 0.000000 -0.000004 60 10500 0.000016 -127.038078 -0.000005 -0.000006 61 10675 0.000024 35.620342 0.000010 0.000007 62 10850 0.000000 0.000000 0.000000 0.000000 63 11025 0.000030 167.302933 -0.000015 0.000003 64 11200 0.000087 -37.874985 0.000034 -0.000027 65 11375 0.000009 -56.846428 0.000002 -0.000004 66 11550 0.000000 0.000000 0.000000 0.000000 67 11725 0.000016 65.517853 0.000003 0.000007 68 11900 0.000020 -62.873825 0.000005 -0.000009 69 12075 0.000049 146.105896 -0.000021 0.000014 70 12250 0.000025 144.522079 -0.000010 0.000007 71 12425 0.000022 83.036545 0.000001 0.000011 72 12600 0.000024 124.945473 -0.000007 0.000010 73 12775 0.000017 150.301559 -0.000008 0.000004 74 12950 0.000006 73.393593 0.000001 0.000003 75 13125 0.000064 101.329811 -0.000006 0.000031 76 13300 0.000022 21.395346 0.000010 0.000004 77 13475 0.000044 -22.719753 0.000020 -0.000008 78 13650 0.000002 67.932175 0.000000 0.000001 79 13825 0.000015 36.952290 0.000006 0.000004 80 14000 0.000061 -90.000000 0.000000 -0.000031 81 14175 0.000024 4.587837 0.000012 0.000001 82 14350 0.000027 -164.107590 -0.000013 -0.000004 83 14525 0.000026 94.922646 -0.000001 0.000013 84 14700 0.000017 113.708290 -0.000003 0.000008 85 14875 0.000049 -94.115646 -0.000002 -0.000025 86 15050 0.000020 3.759978 0.000010 0.000001 87 15225 0.000019 -4.117675 0.000010 -0.000001 88 15400 0.000017 -42.499950 0.000006 -0.000006 89 15575 0.000017 -79.784363 0.000002 -0.000009 90 15750 0.000016 -98.852882 -0.000001 -0.000008 91 15925 0.000054 139.931839 -0.000021 0.000017 92 16100 0.000026 165.618698 -0.000013 0.000003 93 16275 0.000015 -164.399490 -0.000007 -0.000002 94 16450 0.000016 -51.657200 0.000005 -0.000006 95 16625 0.000035 167.911438 -0.000017 0.000004 96 16800 0.000008 -132.260605 -0.000003 -0.000003 97 16975 0.000021 94.842651 -0.000001 0.000011 98 17150 0.000035 -114.642807 -0.000007 -0.000016 99 17325 0.000037 34.428101 0.000015 0.000011 100 17500 0.000016 -93.046417 -0.000000 -0.000008 101 17675 0.000022 32.025196 0.000009 0.000006 102 17850 0.000022 159.313171 -0.000010 0.000004 103 18025 0.000012 137.351868 -0.000005 0.000004 104 18200 0.000034 -31.864140 0.000014 -0.000009 105 18375 0.000028 -157.610428 -0.000013 -0.000005 106 18550 0.000020 56.170975 0.000006 0.000008 107 18725 3999.999756 -0.000001 1999.999878 -0.000040 108 18900 0.000024 70.995369 0.000004 0.000011 109 19075 0.000075 44.852394 0.000026 0.000026 110 19250 0.000026 -126.419708 -0.000008 -0.000010 111 19425 0.000031 -90.581345 -0.000000 -0.000016 112 19600 0.000086 67.500000 0.000017 0.000040 113 19775 0.000032 117.184265 -0.000007 0.000014 114 19950 0.000016 -152.757553 -0.000007 -0.000004 115 20125 0.000047 -12.603742 0.000023 -0.000005 116 20300 0.000009 -68.715248 0.000002 -0.000004 117 20475 0.000067 144.246552 -0.000027 0.000020 118 20650 0.000014 8.900421 0.000007 0.000001 119 20825 0.000013 150.707443 -0.000006 0.000003 120 21000 0.000051 67.017677 0.000010 0.000023 121 21175 0.000016 8.680163 0.000008 0.000001 122 21350 0.000013 -115.376968 -0.000003 -0.000006 123 21525 0.000073 7.219605 0.000036 0.000005 124 21700 0.000011 45.488987 0.000004 0.000004 125 21875 0.000051 79.440254 0.000005 0.000025 126 22050 0.000061 0.000000 0.000031 0.000000 127 22225 0.000012 -156.431137 -0.000006 -0.000002 |
|
|
|
只有小组成员才能发言,加入小组>>
调试STM32H750的FMC总线读写PSRAM遇到的问题求解?
1627 浏览 1 评论
X-NUCLEO-IHM08M1板文档中输出电流为15Arms,15Arms是怎么得出来的呢?
1550 浏览 1 评论
984 浏览 2 评论
STM32F030F4 HSI时钟温度测试过不去是怎么回事?
688 浏览 2 评论
ST25R3916能否对ISO15693的标签芯片进行分区域写密码?
1601 浏览 2 评论
1867浏览 9评论
STM32仿真器是选择ST-LINK还是选择J-LINK?各有什么优势啊?
650浏览 4评论
STM32F0_TIM2输出pwm2后OLED变暗或者系统重启是怎么回事?
518浏览 3评论
536浏览 3评论
stm32cubemx生成mdk-arm v4项目文件无法打开是什么原因导致的?
506浏览 3评论
小黑屋| 手机版| Archiver| 电子发烧友 ( 湘ICP备2023018690号 )
GMT+8, 2024-11-24 09:59 , Processed in 0.767259 second(s), Total 75, Slave 59 queries .
Powered by 电子发烧友网
© 2015 bbs.elecfans.com
关注我们的微信
下载发烧友APP
电子发烧友观察
版权所有 © 湖南华秋数字科技有限公司
电子发烧友 (电路图) 湘公网安备 43011202000918 号 电信与信息服务业务经营许可证:合字B2-20210191 工商网监 湘ICP备2023018690号