DFT与IDFT

DFT与IDFT

一.方法简介

序列x(n)(n=0,1,…N-1)的DFT定义为
$$X（k）=\sum _{n=0}^{N-1}x(n)e^{-j\frac{2\pi nk}{N}}$$

$$设x(n)=a(n)+jb(n),X(k)=A(k)+jB(K),Q=2\pi /N$$

则上式变为：
$$A（k）+jB(k)=\sum _{n=0}^{N-1}[a(n)+jb(n)][\cos (Qnk)-j\sin (Qnk)]$$
即
$$A（k）=\sum _{n=0}^{N-1}[a(n)\cos (Qnk)+b(n)\sin (Qnk)]$$

$$B(k)=\sum _{n=0}^{N-1}[b(n)\cos (Qnk)-a(n)\sin (Qnk)]$$

序列X（k）的IDFT定义为：
$$x(n)=\frac{1}{N}\sum _{k=0}^{N-1}X(k)W_N^{-nk},n=0,1,...N-1$$

$$它与DFT的区别在于将W_N改变为改变为W_N^{-1},并多了一个除以N的运算，公式如下$$

$$a(n)=\frac{1}{N}\sum _{k=0}^{N-1}[A(k)\cos (Qnk)-B(k)\sin (Qnk)]$$

$$b(n)=\frac{1}{N}\sum _{k=0}^{N-1}[B(k)\cos (Qnk)+A(k)\sin (Qnk]$$

二.编程实现

考滤到DFT和IDFT算法过程中有部分相似，可以把它们合成到一个算法。

DFT.c

/*
	x-存放要变换数据的实部
	y-存放要变换数据的虚部
	a-存放变换结果的实部
	b-存放变换结果的虚部
	n-数据长度
	sign-为1时执行DFT，为-1时执行IDFT
*/
#include "math.h"
void dft(x,y,a,b,n,sign)
int n, sign;
double x[],y[],a[],b[];
{
	int i,k;
	double c,d,q,w,s;
	q = 6.28318530718/n;
	for (k=0;k<n;k++)
	{
		w=k*q;
		a[k]=b[k]=0.0;
		for(i=0;i<n;i++)
		{
			d=i*w;
			c=cos(d);
			s=sin(d)*sign;
			a[k]+=c*x[i] + s*y[i];
			b[k]+=c*y[i] - s*x[i];
		}
	}
	if(sign == -1)
	{
		c=1.0/n;
		for (k=0;k<n;k++)
		{
			a[k]=c*a[k];
			b[k]=c*b[k];
		}
	}
}

下面验证此算法，对X(n)=(0,1,2,3,4,5,6,7)，做DFT和IDFT算法

dft_d.c

#include "stdio.h"
#include "math.h"
#include "dft.c"
#define N 4
static double  x[N],y[N],a[N],b[N],c[N];
main(){
	int k;
	int i=0;
	for(i=0; i<N; i++)
	{
		x[i]=i;
		y[i]=0;
		
		
	}
	dft(x,y,a,b,N,1);	//DFT变换
	for(i=0; i<N; i++)
	{
		c[i]=sqrt(a[i]*a[i]+b[i]*b[i]);	//算出模
		printf("%lf + j  %lf \n",a[i],b[i]);//输出变换后结果				
		printf("%lf \n",c[i]); //输出模值
		printf("\n");		
	}
	dft(a,b,x,y,N,-1); //IDFT变换
	for(i=0; i<N; i++)
	{
		printf("%lf \n",x[i]); //输出x(n)的实部
	}
	
}

运行结果：