Help !!! I have been trying to do a complex FFT using the DSP_fft() function provided by the TMS320C64x DSP Library (see TMS320C64x DSP Library Programmers' Reference - SPRU565A). I have provided the "twiddle factor" array by running the tw_fft16x16.exe program to produce a *.c file, also provided in source-code by TI. I have left the "scale" parameter for this proram at the default value (32765.5). (Q1)Does this mean Q1.15 ?? My input data represents complex values (stored in the specified real/imag ordering), and has been quantized to Q1.15 values. Even when my input data represents a sine wave (created in MATLAB, quantized to 16-bits in Q1.15 format), the output of the FFT appears to be garbage. (Q2) Does anyone have any experience with using (any of) theFFT functions provided in this TMS320C64x DSP Library ?? Would really, appreciate any help ...:) Phil |
|
Fixed-point FFT on TMS320C6416
Started by ●November 6, 2002
Reply by ●November 6, 20022002-11-06
phil_levy1 schrieb: > Help !!! > I have been trying to do a complex FFT using the DSP_fft() function > provided by the TMS320C64x DSP Library (see TMS320C64x DSP Library > Programmers' Reference - SPRU565A). I have provided the "twiddle > factor" array by running the tw_fft16x16.exe program to produce a *.c > file, also provided in source-code by TI. I have left the "scale" > parameter for this proram at the default value (32765.5). > > (Q1)Does this mean Q1.15 ?? > > My input data represents complex values (stored in the specified > real/imag ordering), and has been quantized to Q1.15 values. Even > when my input data represents a sine wave (created in MATLAB, > quantized to 16-bits in Q1.15 format), the output of the FFT appears > to be garbage. > > (Q2) Does anyone have any experience with using (any of) theFFT > functions provided in this TMS320C64x DSP Library ?? > > Would really, appreciate any help ...:) > > Phil Hello Phil, assembly optimized FFT subroutines require some precautions: - alignment of input / output buffer(s) in memory - input and output data maybe presented in a special order (i.e. bit-reversed) - buffer / FFT size constraints - twiddle factor size and presentation thus, I would recommend to check the *.h file or the source code of the library for calling constraints. Maybe the output is just bit-reversed i.e. you need to run a bit-reverse function to get in correct order? Rgds. Phil -- Dipl.-Inf. Phil Alder - Field Application Engineer ------------------------------ DSPecialists GmbH Rotherstr. 22 10245 Berlin, Germany http://www.DSPecialists.de ------------------------------ DSPecialists - Making the impossible work! |
Reply by ●November 6, 20022002-11-06
Phil- > I have been trying to do a complex FFT using the DSP_fft() function > provided by the TMS320C64x DSP Library (see TMS320C64x DSP Library > Programmers' Reference - SPRU565A). I have provided the "twiddle > factor" array by running the tw_fft16x16.exe program to produce a *.c > file, also provided in source-code by TI. I have left the "scale" > parameter for this proram at the default value (32765.5). > > (Q1)Does this mean Q1.15 ?? Yes, seems like it; i.e. 1 sign bit and 15 mantissa bits. The scale value is a bit weird, where 32765.5 comes from I'm not sure, I would have expected 32767. I've attached some C code we've used before on C62xx to get very precise fixed-point FFT results. It's the "accurate" version, not the "fast" version, but it should help with understanding the optimized routine provided by TI. Jeff Brower DSP sw/hw engineer Signalogic void FFTRadix2(int xy[], int n, short dir) /******************************************************************************* * NAME * * FFTRadix2 * * USAGE * * This routine is C-callable and can be called as: * * * * void FFTRadix2(short xy[],short n,short w[]) * * * * xy[] --- input and output sequences (dim-n) (input/output) * * n --- FFT size (input) * * w[] --- FFT coefficients (dim-n/2) (input) * * * * ORIGINAL Copyright (C) Texas Instruments 1999 * * * * MODIFICATIONS Copyright (C) Signalogic 2001 * * * * DESCRIPTION * * This routine is used to compute FFT of a complex sequence * * of size n, a power of 2, with "decimation-in-frequency * * decomposition" method, ie, the output is in bit-reversed order. * * * * 4 <= n <= 2048 * * Each complex value has interleaved 16-bit real and imaginary parts. * * Input xy and Coefficients w are 16 bit data. * * sine and cosine data are stored in an integer lookup table as * * sin(n*pi/2048) and cos(n*pi/2048) coefficients multiplied by 32767. * * The sin/cos lookup table limits the maximum FFT size to 2048; this * * can be changed by increasing/decreasing the table size and adjusting * * the value of MAXFFT. Default value of MAXFFT is 11. * * Rounding is used in the butterfly calculation to improve output * * accuracy; can be removed for better performance. * * * TO USE TI'S EXACT METHOD, you need to add Bitrev_xtra.c to the * * project, and then call digitrev_index() and bitrev_cplx() functions * * prior to calling FFTRadix2(). However, please note that this method * * requires a *huge* Index array to hold the pre-calculated bit-reverse * * index values, either on the heap or on the stack (this is why we tend * * to avoid it.) * * * * MEMORY NOTE: * * For assembly language implementation, align xy and w on different word * * boundaries to minimize memory bank hits. * * Both input xy and coefficient w should be aligned on word boundary. * *******************************************************************************/ { int n1,n2,ie,ia,i,j,k,l; short xt,yt; short c,s; int* wcosine; int* wsine; int skip; int TR,TI; wcosine = (int*)&COSINE; wsine = (int*)&SINTAB; skip = 1 << (MAXFFT - ORDVAL); n2 = n; ie = 1; for (k=n; k>1; k>>=1) { n1 = n2; n2 >>= 1; ia = 0; for (j=0; j<n2; j++) { c = (short)-wcosine[skip*ia]; s = (short)wsine[skip*ia]; if (dir < 0) s = -s; // inverse FFT check ia = ia + ie; for (i=j; i < n; i += n1) { l = i + n2; xt = (short)xy[2*l] - (short)xy[2*i]; xy[2*i] = (short)xy[2*i] + (short)xy[2*l]; yt = (short)xy[2*l+1] - (short)xy[2*i+1]; xy[2*i+1] = (short)xy[2*i+1] + (short)xy[2*l+1]; TR = (c*xt - s*yt); if (TR < 0) xy[2*l] = (TR - 16384) / 32767; else xy[2*l] = (TR + 16384) / 32767; TI = (c*yt + s*xt); if (TI < 0) xy[2*l+1] = (TI - 16384) / 32767; else xy[2*l+1] = (TI + 16384) / 32767; } } ie <<= 1; } } |