最近畢設(shè)需要用到希爾伯特變換的知識(shí)，今天做完之后決定還是記錄下思路。

當(dāng)然是數(shù)字信號(hào)的希爾伯特變換

上面是連續(xù)信號(hào)的希爾伯特變換，離散的應(yīng)該也能根據(jù)上面寫(xiě)（沒(méi)現(xiàn)成的圖片，懶得編輯公式了）。

這里打算采用使用卷積的方法來(lái)計(jì)算。由于希爾伯特變換的傳輸函數(shù)的傅里葉變換是  H(w)=        -j   w>=0

                                                                                                                                                                         j     w<0

所以我們可以先求原始信號(hào)的離散傅里葉變換E（K），然后按照下面的公式就可以求出希爾伯特變換之后的信號(hào)的離散傅里葉變換Eh（K）。

 然后對(duì)Eh（K）求反傅里葉變換就可得到我們需要的信號(hào)的希爾伯特變換信號(hào)。

下面貼代碼思路

先建立一個(gè)復(fù)信號(hào)的結(jié)構(gòu)體：

typedef struct {
  Float64 r;  /* 實(shí)部 */
  Float64 i;  /* 虛部 */
} CPX;

接著是離散傅里葉變換的函數(shù)  第一個(gè)參數(shù)dir代表正變換和反變換。

void DFT(int dir,int framelen,CPX *signal,CPX *dft_s)
{
int i,k;
double arg;
double cosarg,sinarg;

for(i=0;i<framelen;i++)
{
arg=-dir*2.0*3.141592654*(double)i/(double)framelen;

     for(k=0;k<framelen;k++)
{
        cosarg=cos(k*arg);
        sinarg=sin(k*arg);
        dft_s[i].r+=(signal[k].r*cosarg-signal[k].i*sinarg);
        dft_s[i].i+=(signal[k].r*sinarg+signal[k].i*cosarg);
}
}
/*返回?cái)?shù)據(jù)*/
if(dir==-1)
{
      for(i=0;i<framelen;i++)
 {
           dft_s[i].r=dft_s[i].r/(double)framelen;
           dft_s[i].i=dft_s[i].i/(double)framelen;
 }

}
}

上主函數(shù)

int main（void）

{

CPX *signal;
CPX *dft_s；

CPX *hdft_s;

CPX *hsignal;

signal=calloc(framelen,sizeof(CPX));  //   原始信號(hào)
   dft_s=calloc(framelen,sizeof(CPX));    //  原始信號(hào)的傅里葉變換
   hdft_s=calloc(framelen,sizeof(CPX));  //  希爾伯特變換的離散 傅里葉變換
   hsignal=calloc(framelen,sizeof(CPX)); //  希爾伯特變換后信號(hào)

DFT(1,framelen,signal, dft_s);  //求原始信號(hào) 傅里葉變換

for(i=0;i<framelen;i++)              //求出希爾伯特變換信號(hào)的傅里葉變換
   {
  if(i<=framelen/2)
  {
 hdft_s[i].r=dft_s[i].i;
 hdft_s[i].i=-dft_s[i].r;
  }
  else
  {
  hdft_s[i].r=-dft_s[i].i;
  hdft_s[i].i=dft_s[i].r;
  }
   }                                  

  DFT(-1,framelen, hdft_s,hsignal);    //利用反傅里葉變換求出希爾伯特變換信號(hào)

}

最后這個(gè)hsignal就是我們要的希爾伯特變換信號(hào)了。