一、什么是DFT？
離散傅里葉變換（Discrete Fourier Transform，縮寫為DFT），是傅里葉變換在時域和頻域上都呈離散的形式，將信號的時域采樣變換為其DTFT的頻域采樣。在形式上，變換兩端（時域和頻域上）的序列是有限長的，而實際上這兩組序列都應當被認為是離散周期信號的主值序列。即使對有限長的離散信號作DFT，也應當將其看作其周期延拓的變換。在實際應用中通常采用快速傅里葉變換計算DFT

二、DFT的公式
離散傅里葉變化定義為
$$X(k)=\sum _{n=0}^{N?1}x(n)W_N^{nk},k=0,1,\dots ,N?1$$

其中，$$W_N^{nk}$$為旋轉因子，特此說明，其計算公式為$$W_N^{nk}=e^{?j\frac{2\pi }{N}nk},\mathrm{N}=0,1_{\mu },\dots ,\mathrm{N}?1$$
1.利用歐拉公式$$e^{j\theta }=\cos \theta +j\sin \theta$$得：$$X(k)=\sum _{n=0}^{N?1}x(n)\cos \left(\frac{2\pi }{N}kn\right)?jx(n)\sin \left(\frac{2\pi }{N}kn\right)$$
2.幅值得計算
$$\mathrm{asin}(x)+b\cos (x)=\sqrt{a^2+b^2}\sin (x+\mathrm{arctanb}/a)$$

三，DFT源程序
抽樣點為：64
原始序列為：
0.6*sin(2*pi*500*i)+0.6*sin(2*pi*50*i)

#include <stdio.h> #include <stdlib.h> #include <math.h> #define pi 3.14159 #define N 64typedef struct {double real;double imag; }complex;void Wn(int ,int ,int ,complex *);//定義旋轉因子 void c_jiafa(complex ,complex ,complex *); //復數加法運算 void c_chengfa(complex ,complex ,complex *); //復數乘法運算int main() {double y[N];complex A[N],a[N],f[N]; //a[N]為輸入的信號，A[N]是DFT后的序列 for(int i=0;i<N;i++) //提供初始信號 {a[i].real=0.6*sin(2*pi*500*i)+0.6*sin(2*pi*50*i);//0.6*sin(2*pi*500*i)+0.6*sin(2*pi*50*i)a[i].imag=0;}for(int k=0;k<N;k++){//int n;complex sum={0,0};for(int n=0;n<N;n++){complex wn,t;Wn(N,n,k,&wn);c_chengfa(a[n],wn,&t);c_jiafa(sum,t,&sum);}A[k].real=sum.real;A[k].imag=sum.imag;} for(int i=0;i<N;i++)//幅值計算{y[i]=sqrt(A[i].real*A[i].real+A[i].imag*A[i].imag);printf("%d\t%lf\n",i,y[i]);}system("pause");return 0; }void Wn(int n1,int i,int k,complex *Wn) //定義旋轉因子 {//用歐拉公式分成實部和虛部 Wn->real=cos(2*pi*i*k/n1);Wn->imag=-sin(2*pi*i*k/n1); }void c_jiafa(complex a,complex b,complex *c) //復數加法運算 {c->real=a.real+b.real; c->imag=a.imag+b.imag; }void c_chengfa(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; }

四，程序編譯，產生dat文件

五，gnuplot繪圖

結果為：

總結

以上是生活随笔為你收集整理的C语言编写DFT计算程序， 并绘制幅度谱的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。