线性卷积计算,计算线性卷积的直接计算方法

#include<xxgc.h> void convol(double *x,int N,double *h,int M,double *y,int L) /*y[0]..y[L-1],x[0]..x[N-1],h[0]..h[M-1];general:L=M+N-1; as:0<=k<=M+N-2,0<=n-k<=M-1,0<=k<=N-1 so:K>=0,K>=n-M+1; k<=N-1,K<=n; */ { int n,k,kk; for(n=0;n<L;n++) for(k=(n>=M)?n-M+1:0,kk=(N>n)?n:N-1,y[n]=0;k<=kk;k++) y[n]+=x[k]*h[n-k]; } main() { int N=5,M=3,L=7,i,a;/*n,gd=VGA,gm=1,w,key,kv;*/ double x[8]={6,2,4,3,8},h[8]={1,2,3},y[20]; dcomplex f[8]={6,0,2,0,4,0,3,0,8,0}; dcomplex r[8]={1,0,2,0,3,0,0,0}; /***step 1 ***/ clrscr(); printf("(1) directly convolution... \n"); convol(h,M,x,N,y,L); for(i=0;i<L;i++){printf("y[%d]=%f\n",i,y[i]);} getch(); /***step 2 ***/ clrscr(); printf("(2) fft convolution\n"); for(i=M+N-1,L=1;i!=0;){i=i>>1;L=L<<1;} fft_ccc(f,L,1); fft_ccc(r,L,1); for(i=0;i<L;i++)f[i]=cmul(f[i],r[i]); fft_ccc(f,L,-1); for(i=0;i<L;i++)printf("f[%d]=%f+j%f,\n",i,f[i].r,f[i].i); getch(); /***step 3 ZhieJieFFT***/ clrscr(); printf("(3) 5 data fft convolution\n"); {dcomplex x[8]={6,0,2,0,4,0,3,0,8,0}; dcomplex h[8]={1,0,2,0,3,0,}; fft_ccc(x,8,1); fft_ccc(h,8,1); for(i=0;i<8;i++)x[i]=cmul(x[i],h[i]); fft_ccc(x,8,-1); for(i=0;i<7;i++)printf("x[%d]=%f+j%f,\n",i,x[i].r,x[i].i); getch(); } /***step 3-2-DieJieSheQuFa***/ clrscr(); printf("(4) 5 data fft convolution in DieJieSheQuFa\n"); {dcomplex xx[5+3+(8-(3+5)%5)]={0,0,0,0,0,0,6,0,2,0,4,0,3,0,8,0}; dcomplex h[8]={1,0,2,0,3,0},y[8]; int j; /*printf("%d===",18%5);getch();*/ fft_ccc(h,8,1); for(j=0;j<=13-8;j=j+5){ for(i=0;i<8;i++)y[i]=xx[i+j]; fft_ccc(y,8,1); for(i=0;i<8;i++)y[i]=cmul(y[i],h[i]); fft_ccc(y,8,-1); /***notie why i=3???****/ for(i=3;i<8;i++)printf("y[%d]=%f+%fj\n",i,y[i].r,y[i].i); } getch(); } /***step 3-2-DieJieXiangJiaFa ***/ clrscr(); printf("(4) 5 data fft convolution in DieJieXianJiaFa\n"); { dcomplex xx[5+3+(8-(3+5)%5)]={6,0,2,0,4,0,3,0,8,0}; dcomplex h[8]={1,0,2,0,3,0},y[8]; int j; fft_ccc(h,8,1); /***please you ansewer why below 13-8****/ for(j=0;j<=13-8;j=j+5){ for(i=0;i<8;i++){if(i>=5){y[i].r=0.0;y[i].i=0.0;} else y[i]=xx[i+j];} for(i=0;i<8;i++)printf("y#[%d]=%f+%fj,\n",i,y[i].r,y[i].i); fft_ccc(y,8,1); for(i=0;i<8;i++)y[i]=cmul(y[i],h[i]); fft_ccc(y,8,-1); for(i=0;i<8;i++)printf("y$[%d]=%f+%fj,\n",i,y[i].r,y[i].i); getch();} } }