久久久久久久999_99精品久久精品一区二区爱城_成人欧美一区二区三区在线播放_国产精品日本一区二区不卡视频_国产午夜视频_欧美精品在线观看免费

標題: VC++6.0測試的FFT程序(快速傅里葉變換) [打印本頁]

作者: heise 時間: 2014-8-18 21:37
標題: VC++6.0測試的FFT程序(快速傅里葉變換)
                                                                                             #include"stdio.h"

#include "math.h"

typedef unsigned char u8;
typedef unsigned int u16;
typedef unsigned long u32;

#define  PI    3.141592653589793238462643       //定義圓周率
#define  FFT_N 64                            //采樣點數

typedef struct Compex                               //復數結構體
{
double real;
double image;
}COMPLEX;

COMPLEX FFT_result[FFT_N]={{1,0},{3,0},{2,0},{5,0},{8,0},{4,0},{1,0},{3,0},{2,0},{5,0},{8,0},{4,0},{1,0},{3,0},{2,0},{5,0},{8,0},{4,0},{1,0},{3,0},{2,0},{5,0},{8,0},{4,0},{1,0},{3,0},{2,0},{5,0},{8,0},{4,0},{1,0},{3,0},{2,0},{5,0},{8,0},{4,0},{1,0},{3,0},{2,0},{5,0},{8,0},{4,0},{1,0},{3,0},{2,0},{5,0},{8,0},{4,0},{1,0},{3,0},{2,0},{5,0},{8,0},{4,0},{1,0},{3,0},{2,0},{5,0},{8,0},{4,0},{1,0},{3,0},{2,0},{5,0}};                         //輸入輸出數據存儲數組

void Init_forward(void);                            //倒位序，采用所謂的雷德算法
COMPLEX MUL_EE(COMPLEX X1,COMPLEX X2);             //復數乘法公式
COMPLEX ADD_EE(COMPLEX X1,COMPLEX X2);             //復數加法公式
COMPLEX SUB_EE(COMPLEX X1,COMPLEX X2);             //復數減法公式

/**************************雷德算法思想*************************
* 自然順序       倒位序
* 000             000
* 001             100
* 010             010
* 011             110
* 100             001
* 101             101
* 110             011
* 111             000

由此可見到位序實際上就是鏡像運算，然而我們沒有采用鏡像算法，（據說可以用匯編來實現比較容易）

我們所要做的工作是：
1.如果我們已知自然順序的一個數想要知道下一個數只需要將當前數加1即可
2.再觀察倒位序之后數據的規律，，，
3.如果我們知道倒位序后的其中一個數，要想求出下一個數，同樣也可以采用加法，
  然而，此處的加法跟我們從小學習的不同，我們要實現的是向低位進位的加法（這里看不懂腦袋砸兩下）
   話不多說，此函數就是實現這個功能仔細分析，如果能看懂說明邏輯思維不錯
*/

/**********************************************************
*函數名稱：COMPLEX MUL_EE(COMPLEX X1,COMPLEX X2)
*函數功能：實現復數乘法
*輸入參數：COMPLEX X1,COMPLEX X2
*返回值：COMPLEX c
***********************************************************/
COMPLEX MUL_EE(COMPLEX X1,COMPLEX X2)
{
COMPLEX c;
c.real  = X1.real*X2.real - X1.image*X2.image;
c.image = X1.real*X2.image + X1.image*X2.real;
return  c;
}

/**********************************************************
*函數名稱：COMPLEX ADD_EE(COMPLEX X1,COMPLEX X2)
*函數功能：進行復數加法
*輸入參數：COMPLEX X1,COMPLEX X2
*返回值：
***********************************************************/
COMPLEX ADD_EE(COMPLEX X1,COMPLEX X2)
{
COMPLEX c;
c.real  = X1.real + X2.real;
c.image = X1.image + X2.image;
return  c;
}

/**********************************************************
*函數名稱：COMPLEX Dcc_EE(COMPLEX X1,COMPLEX X2)
*函數功能：進行復數減法
*輸入參數：COMPLEX X1,COMPLEX X2
*返回值：COMPLEX c
***********************************************************/
COMPLEX SUB_EE(COMPLEX X1,COMPLEX X2)
{
COMPLEX c;
c.real  = X1.real - X2.real;
c.image = X1.image - X2.image;
return  c;
}

/**********************************************************
*函數名稱：void Init_forward(void)
*函數功能：倒位序
*輸入函數：void
*返回值：void
***********************************************************/
void Init_forward(void)
{
u8 I,J,LH,N1,K ;
COMPLEX T;                            //替換結構體
LH = FFT_N/2;                         //N/2
J = LH;
N1 = FFT_N - 2;

for(I=1;I<N1;I++)                      //從1到N-2開始倒位序
{
      if(I<J)                         //此處的意思是當I不等于J時交換位置
      {                               //然而I>J時不交換因為之前已經交換了
         T = FFT_result[I];
         FFT_result[I] = FFT_result[J];
         FFT_result[J] = T;
         }
         K = LH;                         //將給K賦值N/2
      while(J>=K)                      //循環，，判斷所需判斷的位是否為1
      {
         J = J-K;
         K = K/2;
   }
         J = J+K;
   }
}

/**************************************************************
*函數名稱：void  FFT_Run(void)
*函數功能：進行快速傅里葉運算
*輸入參數：void
*返回值：void
***************************************************************/
void  FFT_Run(void)
{
u8 B,P,K;
u8 L = 0;                                           //蝶形變換級數
u8 M = 0;                                           //N = 2^M
u8 J;
u8 FFT_N1 = FFT_N;
COMPLEX  Result_Wn,Result_MUL,Result_ADD,Result_SUB;
Init_forward();                                     //進行倒位序運算

for(M=1; (FFT_N1=FFT_N1/2)!=1;M++);                   //計算蝶形級數
for(L=1;L<=M;L++)
{
   B = pow(2,L-1);                                  // 旋轉因子個數
   for(J=0;J<=B-1;J++)
   {
         P = pow(2,M-L)*J;                            //旋轉因子系數
         for(K=J;K<FFT_N;K=K+pow(2,L))
   {
   Result_Wn.real  =  cos((2*PI/FFT_N)*P);
            Result_Wn.image = -sin((2*PI/FFT_N)*P);

            Result_MUL    =  MUL_EE(Result_Wn,FFT_result[K+B]);  //復數乘法

            Result_ADD    =  ADD_EE(FFT_result[K],Result_MUL); //復數加法
            Result_SUB    =  SUB_EE(FFT_result[K],Result_MUL); //復數減法
            FFT_result[K] =  Result_ADD;                      //把加法后的結果放到 FFT_result[K]
            FFT_result[K+B] =  Result_SUB;                      //把減法之后的結果放到FFT_result[K+B]
   }
   }
}
}

void main(void)
{
u8 a;
u8 M;
u8 FFT_N1 = FFT_N;

FFT_Run();
for(a=0;a<FFT_N;a++)
{
   printf("%f",FFT_result[a].real/100);
   printf(" ");
   printf("%f",FFT_result[a].image/100);
   printf("\n");
}
  a= pow(2,3);
  printf("%d",a);
  printf("\n");
}

經過一星期已搞定，學弟學妹可以看看。。

作者: a7336615a1230 時間: 2015-5-16 23:18
謝謝學長。

作者: mokezhi 時間: 2018-7-22 15:49
謝謝學長咯

作者: jackyhuang0408 時間: 2020-6-7 00:58
學習，謝謝分享！

歡迎光臨 (http://www.zg4o1577.cn/bbs/)