slaon的个人博客分享 http://blog.sciencenet.cn/u/slaon

博文

使用FFTW做FFT变换

已有 11947 次阅读 2017-4-10 20:18 |系统分类:科研笔记| fftw, 频率分辨率, 幅值计算, 直流项

一般FFT分析,都是求一个时间历程的频率成分,

用FFTW中的

fftw_plan_dft_r2c_1d

输入是一个实数列向量in[N],

输出是一个复数空间二维矩阵out[N][2],第1列是实数、第2列是复数。


求幅值时,要对复数向量的每个,取向量模,乘以2,再除以点数N。

$\textup{Amp}=2\sqrt{(x^2_{R}+x^2_{I})}/N$

频率计算

$f=i\cdot \frac{1}{t}$

其中,

$i$      是采样点编号,从1到N;

$1/t$  是频率分辨率(参 陈鹏 博文 link1),也是采样时间倒数,这个很关键,它是FFT变换的最小频率值。

如果t =0.01 s,那么FFT变换后最小频率只能分辨到 100 Hz!


注意:

1. 如果采样时间不到 1 s, 就设法补成 1 s。

一种方法是补零,(参考 彭真明 老师博文 link2

还有一种方法,如果已知频率大概范围,折算出周期后,用正周期的时间历程估计。


2. 如果有直流项(平均值大于零),需要把直流项滤掉。





   fftw_complex *out = NULL;
   fftw_plan p;
   in  = (double *)fftw_malloc(sizeof(double) * len);
   out = (fftw_complex *)fftw_malloc(sizeof(fftw_complex) * len);

   p = fftw_plan_dft_r2c_1d(len, in, out, FFTW_ESTIMATE);
   fftw_execute(p);

   for (int i = 0; i < len; i++)
   {
         amp  = sqrt(out[i][0]*out[i][0] + out[i][1]*out[i][1]);
         amp  = 2.0* amp/len;
         freq  =  i*1.0/t1;
   }









https://blog.sciencenet.cn/blog-531760-1048005.html

上一篇:声功率及板的透射损失
下一篇:去除 cfd post 输出eps文件中的莫名其妙的点
收藏 IP: 114.213.245.*| 热度|

2 杨正瓴 彭真明

该博文允许注册用户评论 请点击登录 评论 (12 个评论)

数据加载中...

Archiver|手机版|科学网 ( 京ICP备07017567号-12 )

GMT+8, 2024-4-27 05:41

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部