lebronze 发表于 2016-7-5 21:56

长度随机的信号FFT如何实现

想用C语言实现FFT、IFFT,但是遇到的问题是:如果信号的长度不是2的整数幂,那应该如何处理?
尝试了下面两种方法:
(1)截断原信号到2的整数幂。
(2)在原信号的后面补零到2的整数幂。
但发现结果都不正确。
这个应该有成熟的处理方法吧,希望懂的大神指导指导,多谢

lebronze 发表于 2016-7-13 10:21

多谢各位回复,目前的解决方法是:对原始信号末尾补零到2的整数幂,在进行FFT。
时域补零,相当于频域插值,不会改变频谱的包络。

impulse 发表于 2016-7-5 22:34

本帖最后由 impulse 于 2016-7-5 22:38 编辑

直接调用FFTW库-几个公认的比较好的、效率比较高的FFT算法之一
http://www.fftw.org/

lebronze 发表于 2016-7-5 22:58

impulse 发表于 2016-7-5 22:34
直接调用FFTW库-几个公认的比较好的、效率比较高的FFT算法之一
http://www.fftw.org/

多谢回复,不过我这边是要搞嵌入式开发,FFT,IFFT只能自己写,不能调用库

hcharlie 发表于 2016-7-6 07:56

从编程角度看,FFT,IFFT 的精髓就是必须是2的整数幂,加零等方法都是凑数的近似方法。

eastar 发表于 2016-7-6 08:18

参考fft matlab代码改成c行吗

impulse 发表于 2016-7-6 12:36

eastar 发表于 2016-7-6 08:18
参考fft matlab代码改成c行吗

matlab fft用的是FFTW

TestGuru 发表于 2016-7-6 12:56

C语言算法可参考: http://forum.vibunion.com/thread-136567-1-5.html

ZH----过客 发表于 2016-7-6 13:11

这是我找到的C++的代码,不知道有没有用你看一下时域抽取法FFT的C++实现。
   下面是算法的流程图

   倒序的流程图

C++实现代码:
1. FFT.h
#pragma once
#ifndef FFT_H
#define FFT_H
#include <cstring>
#include <cmath>
using namespace std;
#define pi 3.141592
class FFT
{
private:
void changeOrder(double* xr,double* xi,int n);
public:
FFT(void);
void FFT_1D(double* ctxr,double* ctxi,double* cfxr,double* cfxi,int len);
void rFFT_1D(double* ctxr,double* cfxr,double* cfxi,int len);
void IFFT_1D(double* cfxr,double* cfxi,double* ctxr,double* ctxi,int len);
void rIFFT_1D(double* cfxr,double* cfxi,double* ctxr,int len);
~FFT(void);
};
#endif
2.FFT.cpp
#include ".fft.h"
FFT::FFT()
{
}
///////////////////////////
//倒序实现
//xr实部,xi虚部,n为2的幂
///////////////////////////
void FFT::changeOrder(double *xr,double *xi,int n)
{
double temp;
int k;
int lh=n/2;
int j=lh;
int n1=n-2;
for(int i=1;i<=n1;i++)
{
if(i<j)
{
   temp=xr;
   xr=xr;
   xr=temp;
   temp=xi;
   xi=xi;
   xi=temp;
}
k=lh;
while(j>=k)
{
   j=j-k;
   k=(int)(k/2+0.5);
}
j=j+k;
}
}
//////////////////////////
//复数FFT
//ctxr和ctxi的长度为len
//cfxr和cfxi的长度为2的幂
//////////////////////////
void FFT::FFT_1D(double *ctxr,double *ctxi,double *cfxr,double *cfxi,int len)
{
int m=ceil(log((double)len)/log(2.0));
int l,b,j,p,k;
double rkb,ikb;
int n=1<<m;
double* rcos=new double;
double* isin=new double;
for(l=0;l<n/2;l++)         
{
rcos=cos(l*pi*2/n);
isin=sin(l*pi*2/n);
}                           
memcpy(cfxr,ctxr,sizeof(double)*len);
memcpy(cfxi,ctxi,sizeof(double)*len);
for(l=len;l<n;l++)
{
cfxr=0;
cfxi=0;
}
changeOrder(cfxr,cfxi,n);//倒序
for(l=1;l<=m;l++)
{
b=(int)(pow(2,l-1)+0.5);
for(j=0;j<b;j++)
{
   p=j*(int)(pow(2,m-l)+0.5);
   for(k=j;k<n;k+=(int)(pow(2,l)+0.5))
   {
    rkb=cfxr*rcos+cfxi*isin;
    ikb=cfxi*rcos-cfxr*isin;
    cfxr=cfxr-rkb;
    cfxi=cfxi-ikb;
    cfxr=cfxr+rkb;
    cfxi=cfxi+ikb;
   }
}
}
delete []rcos;
delete []isin;
}
/////////////////////////////
//实数FFT
//ctxr的长度为len
//cfxr和cfxi的长度为2的幂
////////////////////////////
void FFT::rFFT_1D(double *ctxr,double *cfxr,double *cfxi,int len)
{
int m=ceil(log((double)len)/log(2.0));
int n=1<<m;
double* rcos=new double;
double* isin=new double;
for(int l=0;l<n/2;l++)         
{
rcos=cos(l*pi*2/n);
isin=sin(l*pi*2/n);
}   
double* txr=new double[(len+1)/2];
double* txi=new double[(len+1)/2];
for(int i=0;i<len/2;i++)
{
txr=ctxr;
txi=ctxr;
}
if(len%2==1)
{
txr[(len-1)/2]=ctxr;
txi[(len-1)/2]=0;
}
FFT_1D(txr,txi,cfxr,cfxi,(len+1)/2);
double* x1r=new double;
double* x1i=new double;
double* x2r=new double;
double* x2i=new double;
x1r=cfxr;
x1i=0;
x2r=cfxi;
x2i=0;
for(int k=1;k<n/2;k++)
{
x1r=(cfxr+cfxr)/2.0;
x1i=(cfxi-cfxi)/2.0;
x2r=(cfxi+cfxi)/2.0;
x2i=(-cfxr+cfxr)/2.0;
}
double rkb,ikb;
for(i=0;i<n/2;i++)
{
rkb=x2r*rcos+x2i*isin;
ikb=x2i*rcos-x2r*isin;
cfxr=x1r-rkb;
cfxi=x1i-ikb;
cfxr=x1r+rkb;
cfxi=x1i+ikb;
}
delete []txr;
delete []txi;
delete []x1r;
delete []x1i;
delete []x2r;
delete []x2i;
delete []rcos;
delete []isin;
}
///////////////////////////////
//复数IFFT
//cfxr和cfxi的长度为n(2的幂)
//ctxr和ctxi的长度为len
///////////////////////////////
void FFT::IFFT_1D(double *cfxr,double *cfxi,double *ctxr,double *ctxi,int len)
{
int m=ceil(log((double)len)/log(2.0));
int n=1<<m;
double* txr=new double;
double* txi=new double;
for(int i=0;i<n;i++)
cfxi=-cfxi;
FFT_1D(cfxr,cfxi,txr,txi,n);
for(i=0;i<len;i++)
{
ctxr=txr/n;
ctxi=-txi/n;
}
delete []txr;
delete []txi;
}
//////////////////////////////
//实数IFFT
//cfxr和cfxi的长度为n(2的幂)
//ctxr的长度为len
//////////////////////////////
void FFT::rIFFT_1D(double *cfxr,double *cfxi,double *ctxr,int len)
{
int m=ceil(log((double)len)/log(2.0));
int n=1<<m;
double* txr=new double;
double* txi=new double;
for(int i=0;i<n;i++)
cfxi=-cfxi;
FFT_1D(cfxr,cfxi,txr,txi,n);
for(int i=0;i<len;i++)
{
ctxr=txr/n;
}
delete []txr;
delete []txi;
}
FFT::~FFT(void)
{
}
3.test.cpp
#include <iostream>
#include "FFT.h"
using namespace std;
void main()
{
double xr={1,2,3,4,5,6,7,8,9,10};   //实部
double xi={0,0,0,0,0,0,0,0,0,0};    //虚部
double *cxr,*cxi;
cxr=new double;
cxi=new double;
FFT f;
f.rFFT_1D(xr,cxr,cxi,10);
for(int i=0;i<16;i++)
{
cout<<cxr<<"+j"<<cxi;
cout<<endl;
}
cout<<endl;
double *fxr,*fxi;
fxr=new double;
fxi=new double;
f.rIFFT_1D(cxr,cxi,fxr,10);
for(i=0;i<10;i++)
{
cout<<fxr;
cout<<endl;
}
delete []fxr;
delete []fxi;
delete []cxr;
delete []cxi;

}

lebronze 发表于 2016-7-6 15:17

多谢各位的热心回复,很感动。
不过2的整数幂的fft我也会做,现在的问题是如果信号长度是任意的,可能不是2的整数幂,那要
怎么做fft??
matlab的fft函数就可以实现这个功能,但它内部是封装好的看不到代码。我查网上说是通过多种
混合基(不同长度的基混合)来做的fft。
希望有经验的人能指点迷津

lebronze 发表于 2016-7-6 15:18

eastar 发表于 2016-7-6 08:18
参考fft matlab代码改成c行吗

matlab的fft是封装好的,看不到啊

lebronze 发表于 2016-7-6 15:19

impulse 发表于 2016-7-6 12:36
matlab fft用的是FFTW

好的,那我还是下个FFTW看看他处理非2整数幂fft的思路吧

impulse 发表于 2016-7-6 20:03

lebronze 发表于 2016-7-6 15:19
好的,那我还是下个FFTW看看他处理非2整数幂fft的思路吧

FFTW好像没有嵌入式版本

eastar 发表于 2016-7-7 08:38

lebronze 发表于 2016-7-6 15:18
matlab的fft是封装好的,看不到啊

我看论坛里之前有人问fft源代码的

Generation 发表于 2016-7-7 08:46

你可以参考一下张明的《任意长度FFT算法的C++实现》

lebronze 发表于 2016-7-8 22:36

Generation 发表于 2016-7-7 08:46
你可以参考一下张明的《任意长度FFT算法的C++实现》

多谢回复,我这两天看了一下张明的代码。
他相当于使用多组基的fft来做,但是问题是它的代码是C++的,我本身就看不太懂,在自己用C写真的很困难。。。
我在想可不可以直接用DFT来做
页: [1] 2 3
查看完整版本: 长度随机的信号FFT如何实现