dsp_fft16x16r.h

来自「dm642函数库」· C头文件 代码 · 共 702 行 · 第 1/4 页

H
702
字号
/* ======================================================================== */
/*  TEXAS INSTRUMENTS, INC.                                                 */
/*                                                                          */
/*  DSPLIB  DSP Signal Processing Library                                   */
/*                                                                          */
/*      Release:        Revision 1.04b                                      */
/*      CVS Revision:   1.9     Sun Sep 29 03:32:20 2002 (UTC)              */
/*      Snapshot date:  23-Oct-2003                                         */
/*                                                                          */
/*  This library contains proprietary intellectual property of Texas        */
/*  Instruments, Inc.  The library and its source code are protected by     */
/*  various copyrights, and portions may also be protected by patents or    */
/*  other legal protections.                                                */
/*                                                                          */
/*  This software is licensed for use with Texas Instruments TMS320         */
/*  family DSPs.  This license was provided to you prior to installing      */
/*  the software.  You may review this license by consulting the file       */
/*  TI_license.PDF which accompanies the files in this library.             */
/* ------------------------------------------------------------------------ */
/*          Copyright (C) 2003 Texas Instruments, Incorporated.             */
/*                          All Rights Reserved.                            */
/* ======================================================================== */
/* ======================================================================== */
/*  Assembler compatibility shim for assembling 4.30 and later code on      */
/*  tools prior to 4.30.                                                    */
/* ======================================================================== */
/* ======================================================================== */
/*  End of assembler compatibility shim.                                    */
/* ======================================================================== */
/* ========================================================================= */
/*                                                                           */
/*  TEXAS INSTRUMENTS, INC.                                                  */
/*                                                                           */
/*  NAME                                                                     */
/*        DSP_fft16x16r : 16x16 FFT with rounding and scaling                */
/*                                                                           */
/*                                                                           */
/*  REVISION DATE                                                            */
/*      23-May-2002                                                          */
/*                                                                           */
/*  USAGE                                                                    */
/*        This routine is C-callable and can be called as:                   */
/*                                                                           */
/*        void DSP_fft16x16r                                                 */
/*        (                                                                  */
/*            int                       n,                                   */
/*            short *restrict           ptr_x,                               */
/*            const short *restrict     ptr_w,                               */
/*            unsigned char *restrict   brev,                                */
/*            short *restrict           ptr_y,                               */
/*            int                       radix,                               */
/*            int                       offset,                              */
/*            int                       nmax                                 */
/*        )                                                                  */
/*                                                                           */
/*        N      = length of fft in complex samples, power of 2 <=16384      */
/*        ptr_x  = pointer to complex data input                             */
/*        ptr_w  = pointer to complex twiddle factor (see below)             */
/*        brev   = pointer to bit reverse table containing 64 entries        */
/*        n_min  = smallest fft butterfly used in computation                */
/*                 used for decomposing fft into subffts, see notes          */
/*        offset = index in complex samples of sub-fft from start of main ff */
/*        n_max  = size of main fft in complex samples                       */
/*                                                                           */
/*  DESCRIPTION                                                              */
/*       The benchmark performs a mixed radix forwards fft using             */
/*       a special sequence of coefficients generated in the following       */
/*       way:                                                                */
/*                                                                           */
/*         -* generate vector of twiddle factors for optimized algorithm *-  */
/*        void tw_gen(short * w, int N)                                      */
/*        {                                                                  */
/*          int j, k;                                                        */
/*          double x_t, y_t, theta1, theta2, theta3;                         */
/*          const double PI = 3.141592654, M = 32767.0;                      */
/*                                       -* M is 16383 for scale by 4 *-     */
/*                                                                           */
/*          for (j=1, k=0; j <= N>>2; j = j<<2)                              */
/*          {                                                                */
/*              for (i=0; i < N>>2; i+=j)                                    */
/*              {                                                            */
/*                  theta1 = 2*PI*i/N;                                       */
/*                  x_t = M*cos(theta1);                                     */
/*                  y_t = M*sin(theta1);                                     */
/*                  w[k]   =  (short)x_t;                                    */
/*                  w[k+1] =  (short)y_t;                                    */
/*                                                                           */
/*                  theta2 = 4*PI*i/N;                                       */
/*                  x_t = M*cos(theta2);                                     */
/*                  y_t = M*sin(theta2);                                     */
/*                  w[k+2] =  (short)x_t;                                    */
/*                  w[k+3] =  (short)y_t;                                    */
/*                                                                           */
/*                  theta3 = 6*PI*i/N;                                       */
/*                  x_t = M*cos(theta3);                                     */
/*                  y_t = M*sin(theta3);                                     */
/*                  w[k+4] =  (short)x_t;                                    */
/*                  w[k+5] =  (short)y_t;                                    */
/*                  k+=6;                                                    */
/*              }                                                            */
/*          }                                                                */
/*        }                                                                  */
/*       This redundent set of twiddle factors is size 2*N short samples.    */
/*       As pointed out later dividing these twiddle factors by 2 will give  */
/*       an effective divide by 4 at each stage to guarentee no overflow.    */
/*       The function is accurate to about 68dB of signal to noise ratio     */
/*       to the DFT function below:                                          */
/*                                                                           */
/*        void dft(int n, short x[], short y[])                              */
/*        {                                                                  */
/*           int k,i, index;                                                 */
/*           const double PI = 3  4159654;                                   */
/*           short * p_x;                                                    */
/*           double arg, fx_0, fx_1, fy_0, fy_1, co, si;                     */
/*                                                                           */
/*           for(k = 0; k<n; k++)                                            */
/*           {                                                               */
/*             p_x = x;                                                      */
/*             fy_0 = 0;                                                     */
/*             fy_1 = 0;                                                     */
/*             for(i=0; i<n; i++)                                            */
/*             {                                                             */
/*               fx_0 = (double)p_x[0];                                      */
/*               fx_1 = (double)p_x[1];                                      */
/*               p_x += 2;                                                   */
/*               index = (i*k) % n;                                          */
/*               arg = 2*PI*index/n;                                         */
/*               co = cos(arg);                                              */
/*               si = -sin(arg);                                             */
/*               fy_0 += ((fx_0 * co) - (fx_1 * si));                        */
/*               fy_1 += ((fx_1 * co) + (fx_0 * si));                        */
/*             }                                                             */
/*             y[2*k] = (short)2*fy_0/sqrt(N);                               */
/*             y[2*k+1] = (short)2*fy_1/sqrt(N);                             */
/*           }                                                               */
/*        }                                                                  */
/*       Scaling takes place at each stage except the last one.              */
/*       This is a divide by 2 to prevent overflow. All shifts are rounded t */
/*       reduce truncation noise power by 3dB.                               */
/*       The function takes the table and input data and calculates the fft  */
/*       producing the frequency domain data in the Y array.                 */
/*       As the fft allows every input point to effect every output point in */
/*       a cache based system such as the c6211, this causes cache thrashing */
/*       This is mitigated by allowing the main fft of size N to be divided  */
/*       into several steps, allowing as much data reuse as possible.        */
/*                                                                           */
/*       For example the following function:                                 */
/*                                                                           */
/*       DSP_fft16x16r  (1024, &x_asm[0],&w[0],y_asm,brev,4,  0,1024);       */
/*                                                                           */
/*       is equvalent to:                                                    */
/*                                                                           */
/*       DSP_fft16x16r  (1024,&x_asm[2*0],  &w[0]  ,y_asm,brev,256, 0,1024); */
/*       DSP_fft16x16r  (256, &x_asm[2*0],  &w[2*768],y_asm,brev,4, 0,1024); */
/*       DSP_fft16x16r (256, &x_asm[2*256],&w[2*768],y_asm,brev,4, 256,1024); */
/*       DSP_fft16x16r (256, &x_asm[2*512],&w[2*768],y_asm,brev,4, 512,1024); */
/*       DSP_fft16x16r (256, &x_asm[2*768],&w[2*768],y_asm,brev,4, 768,1024); */
/*                                                                           */
/*       Notice how the 1st fft function is called on the entire 1K data set */
/*       it covers the 1st pass of the fft until the butterfly size is 256.  */
/*       The following 4 ffts do 256 pt ffts 25% of the size. These continue */
/*       down to the end when the buttefly is of size 4. The use an index to */
/*       the main twiddle factor array of 0.75*2*N. This is because the      */
/*       twiddle factor array is composed of successively decimated versions */
/*       of the main array.                                                  */
/*                                                                           */
/*       N not equal to a power of 4 can be used, i.e. 512. In this case to  */
/*       decompose the fft the following would be needed :                   */
/*                                                                           */
/*       DSP_fft16x16r   (512, &x_asm[0],&w[0],y_asm,brev,2,  0,512);        */
/*                                                                           */
/*       is equvalent to:                                                    */
/*                                                                           */
/*       DSP_fft16x16r   (512, &x_asm[0],    &w[0],  y_asm,brev,128, 0,512); */
/*       DSP_fft16x16r   (128, &x_asm[2*0],  &w[2*384],y_asm,brev,2, 0,512); */
/*       DSP_fft16x16r (128, &x_asm[2*128],&w[2*384],y_asm,brev,2, 128,512); */

⌨️ 快捷键说明

复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?