24小时热门版块排行榜

返回列表

【奖励】本帖被评价1次，作者baobiao007增加金币 0.8 个

baobiao007

木虫 (职业作家)

应助: 201 (大学生)
金币: 6482.7
帖子: 3050
在线: 1009.9小时
虫号: 505962

[资源] 【分享】时域褶积带通滤波c程序【原创，无重复】

#include
#include
#include
#include"FFT.h"

/*带通滤波因子
  num--滤波因子长度（奇数）
  f0--通频带中心频率
  hf--通频带半宽度
  dt--时间采样间隔(与输入信号相同)
*/
void FilterH(float h[],int num,float f0,float hf,float dt)
{
int i,j;
double pi=3.1415926;
float t;
float *w;
            //产生汉明时窗w[]
            w=(float *)calloc(num,sizeof(float));
for(i=0; i<=num/2; i++)
      w[ i ]=0.54+0.46*cos(2*pi*i/num);
for(i=num/2+1; i       w[ i ]=w[num-i];
//产生滤波因子
h[0]=4.0*hf;
for(i=1; i {
      h[ i ]=0.0;
                  if(i<=num/2)
      { t=i*dt;j=i;}
      else
      { j=num-i; t= -j*dt;}
      h[ i ]=2.0*sin(2*pi*hf*t)*cos(2*pi*f0*t)*w[ j ]/(pi*t);
   }
}
/*输入信号与滤波因子褶积
  h[]--滤波因子,x[]--输入信号
  保证x[]更长,结果存放在y[],长度也为m
*/
void ConvFilter(float h[],int n,float x[],int m,float y[],float dt)
{
int i,j,k;
            int half;
half=(n-1)/2;
for(i=0; i {
                  y[ i ]=0.0;
      for(j= -half; j<=half; j++)
      {
k=i-j;
if(j<0 && k          y[ i ]+=h[n+j]*x[k]*dt;
if(j>=0 && k>=0)
         y[ i ]+=h[j]*x[k]*dt;
      }
   }
}
//产生合成信号
//30 60 90hz正弦信号相加
void Generate(float x[],int num)
{
int i;
double pi=3.1415926;
float dt=0.002;
float f1,f2,f3;
f1=30.0;
f2=60.0;
f3=90.0;
for(i=0; i       x[ i ]=sin(2*pi*f1*i*dt)+sin(2*pi*f2*dt*i)+sin(2*pi*f3*i*dt);
}
void main()
{
const int N=101;//滤波因子长度
const int M=128;//信号长度
float Inisig[M]={0};//要滤波的信号
float Afsig2[M]={0};//滤波后的信号
float Afsig[M]={0};//理想情况下滤波后的信号
float h[N];//滤波因子
int i;
FILE *fp1,*fp2,*fp3;
fp1=fopen("inisig.txt","w" );//合成信号
      fp2=fopen("afsig.txt","w" );//理想50hz信号
fp3=fopen("afsig2.txt","w" );//滤波得到的50hz信号

//产生信号,20hz,50hz,80hz正弦信号叠加
      double pi=3.1415926;
float dt=0.002;
float f1,f2,f3;
f1=20.0;
f2=50.0;
f3=80.0;
for(i=0; i { Afsig[ i ]=sin(2*pi*f2*dt*i);
Inisig[ i ]=Afsig+sin(2*pi*f1*i*dt)+sin(2*pi*i*f3*dt);
}

FilterH(h,N,50.0,10.0,dt);
            //滤去低于70hz和高于80hz信号
ConvFilter(h,N,Inisig,M,Afsig2,dt);

for(i=0; i {
         fprintf(fp1,"%f\n",Inisig[ i ]);
                     fprintf(fp2,"%f\n",Afsig[ i ]);
                     fprintf(fp3,"%f\n",Afsig2[ i ]);
}
            fclose(fp1);
fclose(fp2);
fclose(fp3);
}
滤波前：
滤波后：

[ Last edited by baobiao007 on 2011-2-14 at 20:42 ]