当前位置:文档之家› DSP实验(fft算法的C语言实现)——njust

DSP实验(fft算法的C语言实现)——njust

/*main()函数*/
#include
#include
#include "fft.h"

// float result;
int N=64;
float s[128]={0.494875,0.038333,0.227436,0.327883,0.899469,0.313730,0.251676,0.432989,
0.842382,0.184489,0.508179,0.452240,0.325584,0.380076,0.886480,0.761261,
0.883766,0.457406,0.799202,0.134077,0.065314,0.375145,0.373523,0.484022,
0.969459,0.342061,0.252689,0.584887,0.523704,0.163419,0.486398,0.496061,
0.843194,0.806198,0.857786,0.609754,0.565730,0.611899,0.102977,0.158316,
0.413650,0.560410,0.268677,0.784254,0.387871,0.030984,0.585502,0.558559,
0.200696,0.087422,0.933230,0.259380,0.204171,0.049208,0.606161,0.546349,
0.095837,0.636996,0.442948,0.066382,0.374293,0.249103,0.924875,0.629499,
0.878309,0.641674,0.798391,0.435026,0.981140,0.095958,0.527482,0.545646,
0.284343,0.370803,0.064693,0.544809,0.836376,0.145322,0.171520,0.068047,
0.824012,0.133971,0.884786,0.514737,0.963636,0.120495,0.048290,0.380152,
0.412791,0.401391,0.420997,0.376954,0.907337,0.670162,0.961839,0.162979,
0.748649,0.374066,0.454237,0.038561,0.562432,0.372312,0.792784,0.795231,
0.382914,0.252790,0.342928,0.967804,0.479808,0.368328,0.764567,0.377149,
0.900306,0.183432,0.368317,0.917457,0.515916,0.090307,0.735311,0.004712,
0.603123,0.956867,0.397432,0.731551,0.684639,0.978503,0.203785,0.593303};

main()
{
int i;
void FFT(struct compx * ,int );
struct compx xin[64];
for(i=0;i<64;i++)
{
xin[i].real=s[2*i];
xin[i].imag=s[2*i+1];
}

FFT(xin,64);

for(i=0;i<64;i++)
{
printf("%.6f",xin[i].real);
printf("+%.6fj\n",xin[i].imag);
// result[i]=sqrt(pow(s[i].real,2)+pow(s[i].imag,2));
}
return 0;
}

void FFT( struct compx * xin,int N) /*FFT变换函数*/
{
int f,m,LH,nm,i,k,j,L;
double p,ps;
int le,B,ip;
float pi;
struct compx w,t;
LH=N/2;
f=N;
for(m=1;(f=f/2)!=1;m++)
{
;
}
nm=N-2;
j=N/2;
/*以下是位倒序*/
for(i=1;i<=nm;i++)
{
if(i{
t.real=xin[j].real;
t.imag=xin[j].imag;
xin[j].real=xin[i].real;
xin[j].imag=xin[i].imag;
xin[i].real=t.real;
xin[i].imag=t.imag;
}
k=LH;
while(j>=k)
{
j=j-k;
k=k/2;
}
j=j+k;
}
/*以下是fft运算*/
{
for(L=1;L<=m;L++)
{
le=pow(2,L);
B=le/2;
pi=3.1415926;
for(j=0;j<=B-1;j++)
{
p=pow(2,m-L)*j;
ps=2*pi/N*p;
w.real=cos(ps);
w.imag=-sin(ps);

for(i=j;i<=N-1;i=i+le)
{
ip=i+B;
t=EE(xin[ip],w);
xin[ip].real=xin[i].real-t.real;
xin[ip].imag=xin[i].imag-t.imag;
xin[i].real=xin[i].real+t.real;
xin[i].imag=xin[i].imag+t.imag;
}
}
}
}
return;
}



vector.asm 文件

.global _c_int00
.sect ".vec"

RESET:
MVKL _c_int00, B0
MVKH _c_int00, B0

B B0
nop
nop
nop
nop
nop

NMI:
B IRP
NOP
NOP
NOP
NOP
NOP
NOP
NOP

.space 8*4*2

IE4
B IRP
NOP
NOP
NOP
NOP
NOP
NOP
NOP

IE5
B IRP
NOP
NOP
NOP
NOP
NOP
NOP
NOP

IE6
; B _RecFromPC
B IRP
NOP
NOP
NOP
NOP
NOP
NOP
NOP

IE7
B IRP
NOP
NOP
NOP
NOP
NOP
NOP
NOP

;
IE8
B IRP
NOP
NOP
NOP
NOP
NOP
NOP
NOP

IE9
; MV B0,B20
; MV B1,B21
; MV B2,B22
; MV B3,B23
B IRP
NOP
NOP
NOP
NOP
NOP
NOP
NOP

IE10
B IRP
NOP
NOP
NOP
NOP
NOP
NOP
NOP

IE11
B IRP
NOP
NOP
NOP
NOP
NOP
NOP
NOP

IE12
B IRP
NOP
NOP
NOP
NOP
NOP
NOP
NOP

IE13
B IRP
NOP
NOP
NOP
NOP
NOP
NOP
NOP

IE14
B IRP
NOP
NOP
NOP
NOP
NOP
NOP
NOP

IE15
B IRP
NOP
NOP
NOP
NOP
NOP
NOP
NOP

fft.h文件

#ifndef _fft_h_
#define _fft_h_
struct compx /*定义一个复数结构*/
{
float real;
float imag;
}compx;
struct compx b1,b2,b3;
struct compx EE(struct compx b1,struct compx b2) /*定义复数相乘结构*/
{
struct compx b3;
b3.real=b1.real*b2.real-b1.imag*b2.imag;
b3.imag=b1.real*b2.imag+b1.imag*b2.real;
return(b3);
}
#endif

相关主题
文本预览
相关文档 最新文档