Dear all,
I am running a program for DSK6711 that calculates the FFT-radix 4. As is shown
in the attached code, data length is defined as a constant value of 128 samples
(NUMDATA)and twiddle arrays (TR,TI), which values are defined in twiddleI.h and
twiddleR.h, are 64 samples each one (all files are avalaible in the CD "DSP
System Design Using the TMS320C6000"). At this point the program works ok!
My next step is to increase the number of points of FFT but before do it, I want
to know how was the twiddle factor generated. To understand it I removed
twiddleI.h and twiddleR.h files from my project and replaced code in the main
function with the following loop:
for (i=0;i
{
W4[i].real = cos(2*PI*i/NUMDATA);
W4[i].imag = sin(2*PI*i/NUMDATA);
}
but now the output values are wrong for the same 128 data record. So I have
compared the real and imaginary values generated with loop with values defined
in twiddleI/R.h files and they are different and I don't understand why!
So my questions are:
(1) How were generated TR and TI values in header files to make correct the FFT
calculation?
(2) How can be possible that TR,TI were 64 samples length when in Chassaing
(page 204 of TMS320C6x book) the length of the set up for a FFT radix 4 is
N/RADIX length (whichs means in my case 128/42!)?
The source code is listed below.
Any help would be really appreciated! =)
Thanks in advance,
Benjamin Sanchez
Electronic Department
Technical University of Catalonia UPC
Barcelona (SPAIN)
/******************/
fftradix4.c
/******************/
#include
#include "twiddleR.h"
#include "twiddleI.h"
#define NUMDATA 128 // number of real data samples 128
#define NUMPOINTS 64 // number of point in the DFT,
// NUMDATA/2
#define PI 3.141592653589793 // definition of pi //
typedef struct { //data type for the radix-4 twiddle factors
short imag;
short real;
} COEFF;
typedef struct {
short real;
short imag;
} COMPLEX;
#pragma DATA_ALIGN(x,NUMPOINTS);
COMPLEX x[NUMPOINTS+1];// array of complex DFT data
COEFF W4[NUMPOINTS];
short g[NUMDATA];
COMPLEX A[NUMPOINTS]; // array of complex A coefficients //
COMPLEX B[NUMPOINTS]; // array of complex B coefficients //
COMPLEX G[2*NUMPOINTS]; // array of complex DFT result //
unsigned short IIndex[NUMPOINTS], JIndex[NUMPOINTS];
int count;
int magR[NUMDATA];
int magI[NUMDATA];
int mag[NUMDATA];
...
int data[NUMDATA]={...}; //signal waveform
int output[NUMDATA/2];
void main(void)
{
int i,j,k;
short tr[NUMPOINTS], ti[NUMPOINTS];
//Read Twiddle factors to COMPLEX array and make Q-15;
make_q15(tr, TR, NUMPOINTS); //Data in Header files from Matlab
make_q15(ti, TI, NUMPOINTS);
for(i=0;i
{
W4[i].real = tr[i];
W4[i].imag = ti[i];
}
// Initialize A,B arrays //
for(k=0; k
{
A[k].imag = (short)(16383.0*(-cos(2*PI/(double)(2*NUMPOINTS)*(double)k)));
A[k].real = (short)(16383.0*(1.0 -
sin(2*PI/(double)(2*NUMPOINTS)*(double)k)));
B[k].imag = (short)(16383.0*(cos(2*PI/(double)(2*NUMPOINTS)*(double)k)));
B[k].real = (short)(16383.0*(1.0 +
sin(2*PI/(double)(2*NUMPOINTS)*(double)k)));
}
// Initialize tables for FFT digit reversal function //
R4DigitRevIndexTableGen(NUMPOINTS, &count, IIndex, JIndex);
for(i=0;i
g[i]a[i];
// Call FFT algorithm //
fft();
for (k=0; k
{
magR[k] = (G[k].real*G[k].real) << 1;
magI[k] = (G[k].imag*G[k].imag) << 1;
mag[k]=magR[k]+magI[k];
}
}
void fft()
{
int n;
for (n=0; n
{
x[n].imag = g[2*n + 1];
x[n].real = g[2*n];
}
radix4(NUMPOINTS, (short *)x, (short *)W4);
digit_reverse((int *)x, IIndex, JIndex, count);
x[NUMPOINTS].real = x[0].real;
x[NUMPOINTS].imag = x[0].imag;
split1(NUMPOINTS, x, A, B, G);
G[NUMPOINTS].real = x[0].real - x[0].imag;
G[NUMPOINTS].imag = 0;
for (n=1; n
{
G[2*NUMPOINTS-n].real = G[n].real;
G[2*NUMPOINTS-n].imag = -G[n].imag;
}
}
void radix4(int n, short x[], short w[])
{
int n1, n2, ie, ia1, ia2, ia3, i0, i1, i2, i3, j, k;
short t, r1, r2, s1, s2, co1, co2, co3, si1, si2, si3;
n2 = n;
ie = 1;
for (k = n; k > 1; k >>= 2)
{
n1 = n2;
n2 >>= 2;
ia1 = 0;
for (j = 0; j < n2; j++)
{
ia2 = ia1 + ia1;
ia3 = ia2 + ia1;
co1 = w[ia1 * 2 + 1];
si1 = w[ia1 * 2];
co2 = w[ia2 * 2 + 1];
si2 = w[ia2 * 2];
co3 = w[ia3 * 2 + 1];
si3 = w[ia3 * 2];
ia1 = ia1 + ie;
for (i0 = j; i0 < n; i0 += n1)
{
i1 = i0 + n2;
i2 = i1 + n2;
i3 = i2 + n2;
r1 = x[2 * i0] + x[2 * i2];
r2 = x[2 * i0] - x[2 * i2];
t = x[2 * i1] + x[2 * i3];
x[2 * i0] = r1 + t;
r1 = r1 - t;
s1 = x[2 * i0 + 1] + x[2 * i2 + 1];
s2 = x[2 * i0 + 1] - x[2 * i2 + 1];
t = x[2 * i1 + 1] + x[2 * i3 + 1];
x[2 * i0 + 1] = s1 + t;
s1 = s1 - t;
x[2 * i2] = (r1 * co2 + s1 * si2) >> 15;
x[2 * i2 + 1] = (s1 * co2-r1 * si2)>>15;
t = x[2 * i1 + 1] - x[2 * i3 + 1];
r1 = r2 + t;
r2 = r2 - t;
t = x[2 * i1] - x[2 * i3];
s1 = s2 - t;
s2 = s2 + t;
x[2 * i1] = (r1 * co1 + s1 * si1) >>15;
x[2 * i1 + 1] = (s1 * co1-r1 * si1)>>15;
x[2 * i3] = (r2 * co3 + s2 * si3) >>15;
x[2 * i3 + 1] = (s2 * co3-r2 * si3)>>15;
}
}
ie <<= 2;
}
}
void digit_reverse(int *yx, unsigned short *JIndex, unsigned short *IIndex, int
count)
{
int i;
unsigned short I, J;
int YXI, YXJ;
for (i = 0; i
{
I = IIndex[i];
J = JIndex[i];
YXI = yx[I];
YXJ = yx[J];
yx[J] = YXI;
yx[I] = YXJ;
}
}
void R4DigitRevIndexTableGen(int n, int *count, unsigned short *IIndex, unsigned
short
*JIndex)
{
int j, n1, k, i;
j = 1;
n1 = n - 1;
*count = 0;
for(i=1; i<=n1; i++)
{
if(i < j)
{
IIndex[*count] = (unsigned short)(i-1);
JIndex[*count] = (unsigned short)(j-1);
*count = *count + 1;
}
k = n >> 2;
while(k*3 < j)
{
j = j - k*3;
k = k >> 2;
}
j = j + k;
}
}
void split1(int N, COMPLEX *X, COMPLEX *A, COMPLEX *B, COMPLEX *G)
{
int k;
int Tr, Ti;
for (k=0; k
{
Tr = (int)X[k].real * (int)A[k].real - (int)X[k].imag * (int)A[k].imag +
(int)X[N-k].real * (int)B[k].real + (int)X[N-k].imag * (int)B[k].imag;
G[k].real = (short)(Tr>>15);
Ti = (int)X[k].imag * (int)A[k].real + (int)X[k].real * (int)A[k].imag +
(int)X[N-k].real * (int)B[k].imag - (int)X[N-k].imag * (int)B[k].real;
G[k].imag = (short)(Ti>>15);
}
}
void make_q15(short out[], float in[], int N)
{
int i;
for(i=0;i
{
out[i]=0x7fff*in[i]; //Convert to Q-15, good approximate
}
}
/******************/
twiddleR.h
/******************/
float TR[64] ={0.9999,
0.9952,
0.9808,
0.9569,
0.9239,
0.8819,
0.8315,
0.7730,
0.7071,
0.6344,
0.5556,
0.4714,
0.3827,
0.2903,
0.1951,
0.0980,
0.0000,
-0.0980,
-0.1951,
-0.2903,
-0.3827,
-0.4714,
-0.5556,
-0.6344,
-0.7071,
-0.7730,
-0.8315,
-0.8819,
-0.9239,
-0.9569,
-0.9808,
-0.9952,
-1.0000,
-0.9952,
-0.9808,
-0.9569,
-0.9239,
-0.8819,
-0.8315,
-0.7730,
-0.7071,
-0.6344,
-0.5556,
-0.4714,
-0.3827,
-0.2903,
-0.1951,
-0.0980,
0.0000,
0.0980,
0.1951,
0.2903,
0.3827,
0.4714,
0.5556,
0.6344,
0.7071,
0.7730,
0.8315,
0.8819,
0.9239,
0.9569,
0.9808,
0.9952};
/******************/
twiddleI.h
/******************/
float TI[64] = {0.0000,
0.0980,
0.1951,
0.2903,
0.3827,
0.4714,
0.5556,
0.6344,
0.7071,
0.7730,
0.8315,
0.8819,
0.9239,
0.9569,
0.9808,
0.9952,
0.9999,
0.9952,
0.9808,
0.9569,
0.9239,
0.8819,
0.8315,
0.7730,
0.7071,
0.6344,
0.5556,
0.4714,
0.3827,
0.2903,
0.1951,
0.0980,
0.0000,
-0.0980,
-0.1951,
-0.2903,
-0.3827,
-0.4714,
-0.5556,
-0.6344,
-0.7071,
-0.7730,
-0.8315,
-0.8819,
-0.9239,
-0.9569,
-0.9808,
-0.9952,
-1.0000,
-0.9952,
-0.9808,
-0.9569,
-0.9239,
-0.8819,
-0.8315,
-0.7730,
-0.7071,
-0.6344,
-0.5556,
-0.4714,
-0.3827,
-0.2903,
-0.1951,
-0.0980};
