Home | History | Annotate | Download | only in arm
      1 /*
      2  *  Copyright (c) 2013 The WebRTC project authors. All Rights Reserved.
      3  *
      4  *  Use of this source code is governed by a BSD-style license
      5  *  that can be found in the LICENSE file in the root of the source
      6  *  tree. An additional intellectual property rights grant can be found
      7  *  in the file PATENTS.  All contributing project authors may
      8  *  be found in the AUTHORS file in the root of the source tree.
      9  *
     10  * Some code in this file was originally from file omxSP_FFTInit_R_S16S32.c
     11  * which was licensed as follows.
     12  * It has been relicensed with permission from the copyright holders.
     13  */
     14 
     15 /*
     16  * OpenMAX DL: v1.0.2
     17  * Last Modified Revision:
     18  * Last Modified Date:
     19  *
     20  * (c) Copyright 2007-2008 ARM Limited. All Rights Reserved.
     21  */
     22 
     23 #include <stdint.h>
     24 
     25 #include "dl/api/arm/armOMX.h"
     26 #include "dl/api/omxtypes.h"
     27 #include "dl/sp/api/armSP.h"
     28 #include "dl/sp/api/omxSP.h"
     29 
     30 /**
     31  * Function: omxSP_FFTInit_R_S16
     32  *
     33  * Description:
     34  * Initialize the real forward-FFT specification information struct.
     35  *
     36  * Remarks:
     37  * This function is used to initialize the specification structures
     38  * for functions <ippsFFTFwd_RToCCS_S16_Sfs> and
     39  * <ippsFFTInv_CCSToR_S16_Sfs>. Memory for *pFFTSpec must be
     40  * allocated prior to calling this function. The number of bytes
     41  * required for *pFFTSpec can be determined using
     42  * <FFTGetBufSize_R_S16>.
     43  *
     44  * Parameters:
     45  * [in]  order       base-2 logarithm of the desired block length;
     46  *			   valid in the range [1,12].
     47  * [out] pFFTFwdSpec pointer to the initialized specification structure.
     48  *
     49  * Return Value:
     50  * Standard omxError result. See enumeration for possible result codes.
     51  *
     52  */
     53 
     54 OMXResult omxSP_FFTInit_R_S16(OMXFFTSpec_R_S16* pFFTSpec, OMX_INT order) {
     55   OMX_INT i = 0, j = 0;
     56   OMX_SC16 *pTwiddle = NULL, *pTwiddle1 = NULL, *pTwiddle2 = NULL;
     57   OMX_SC16 *pTwiddle3 = NULL, *pTwiddle4 = NULL;
     58   OMX_S16 *pBuf = NULL;
     59   OMX_U16 *pBitRev = NULL;
     60   OMX_U32 pTmp = 0;
     61   OMX_INT Nby2 = 0, N = 0, M = 0, diff = 0, step = 0;
     62   OMX_S16 x = 0, y = 0, xNeg = 0;
     63   OMX_S32 xS32 = 0, yS32 = 0;
     64   ARMsFFTSpec_R_SC16 *pFFTStruct = NULL;
     65 
     66   /* Order zero not allowed */
     67   if (order == 0) {
     68     return OMX_Sts_BadArgErr;
     69   }
     70 
     71   /* Do the initializations */
     72   pFFTStruct = (ARMsFFTSpec_R_SC16*) pFFTSpec;
     73   Nby2 = 1 << (order - 1);
     74   N = Nby2 << 1;
     75   pBitRev = NULL ;  /* optimized implementations don't use bitreversal */
     76   pTwiddle = (OMX_SC16*) (sizeof(ARMsFFTSpec_R_SC16) + (OMX_S8*)pFFTSpec);
     77 
     78   /* Align to 32 byte boundary */
     79   pTmp = ((uintptr_t)pTwiddle)&31;  /* (uintptr_t)pTwiddle % 32 */
     80   if(pTmp != 0) {
     81     pTwiddle = (OMX_SC16*) ((OMX_S8*)pTwiddle + (32 - pTmp));
     82   }
     83 
     84   pBuf = (OMX_S16*) (sizeof(OMX_SC16) * (5 * N / 8) + (OMX_S8*)pTwiddle);
     85 
     86   /* Align to 32 byte boundary */
     87   pTmp = ((OMX_U32)pBuf)&31;                 /* (OMX_U32)pBuf % 32 */
     88   if(pTmp != 0) {
     89     pBuf = (OMX_S16*)((OMX_S8*)pBuf + (32 - pTmp));
     90   }
     91 
     92   /*
     93    * Filling Twiddle factors : exp^(-j*2*PI*k/ (N/2) ) ; k=0,1,2,...,3/4(N/2).
     94    * N/2 point complex FFT is used to compute N point real FFT.
     95    * The original twiddle table "armSP_FFT_S32TwiddleTable" is of size
     96    * (MaxSize/8 + 1). Rest of the values i.e., up to MaxSize are calculated
     97    * using the symmetries of sin and cos.
     98    * The max size of the twiddle table needed is 3/4(N/2) for a radix-4 stage.
     99    *
    100    * W = (-2 * PI) / N
    101    * N = 1 << order
    102    * W = -PI >> (order - 1)
    103    *
    104    * Note we use S32 twiddle factor table and round the values to 16 bits.
    105    */
    106 
    107   M = Nby2 >> 3;
    108   diff = 12 - (order - 1);
    109   step = 1 << diff;  /* Step into the twiddle table for the current order */
    110 
    111   xS32 = armSP_FFT_S32TwiddleTable[0];
    112   yS32 = armSP_FFT_S32TwiddleTable[1];
    113   x = (xS32 + 0x8000) >> 16;
    114   y = (yS32 + 0x8000) >> 16;
    115   xNeg = 0x7FFF;
    116 
    117   if((order-1) >= 3) {
    118     /* i = 0 case */
    119     pTwiddle[0].Re = x;
    120     pTwiddle[0].Im = y;
    121     pTwiddle[2 * M].Re = -y;
    122     pTwiddle[2 * M].Im = xNeg;
    123     pTwiddle[4 * M].Re = xNeg;
    124     pTwiddle[4 * M].Im = y;
    125 
    126     for (i=1; i<=M; i++){
    127       OMX_S16 x_neg = 0, y_neg = 0;
    128       j = i * step;
    129 
    130       xS32 = armSP_FFT_S32TwiddleTable[2 * j];
    131       yS32 = armSP_FFT_S32TwiddleTable[2 * j + 1];
    132       x = (xS32 + 0x8000) >> 16;
    133       y = (yS32 + 0x8000) >> 16;
    134       /* |x_neg = -x| doesn't work when x is 0x8000. */
    135       x_neg = (-(xS32 + 0x8000)) >> 16;
    136       y_neg = (-(yS32 + 0x8000)) >> 16;
    137 
    138       pTwiddle[i].Re = x;
    139       pTwiddle[i].Im = y;
    140       pTwiddle[2 * M - i].Re = y_neg;
    141       pTwiddle[2 * M - i].Im = x_neg;
    142       pTwiddle[2 * M + i].Re = y;
    143       pTwiddle[2 * M + i].Im = x_neg;
    144       pTwiddle[4 * M - i].Re = x_neg;
    145       pTwiddle[4 * M - i].Im = y;
    146       pTwiddle[4 * M + i].Re = x_neg;
    147       pTwiddle[4 * M + i].Im = y_neg;
    148       pTwiddle[6 * M - i].Re = y;
    149       pTwiddle[6 * M - i].Im = x;
    150     }
    151   }
    152   else {
    153     if ((order - 1) == 2) {
    154       pTwiddle[0].Re = x;
    155       pTwiddle[0].Im = y;
    156       pTwiddle[1].Re = -y;
    157       pTwiddle[1].Im = xNeg;
    158       pTwiddle[2].Re = xNeg;
    159       pTwiddle[2].Im = y;
    160     }
    161     if ((order-1) == 1) {
    162       pTwiddle[0].Re = x;
    163       pTwiddle[0].Im = y;
    164     }
    165   }
    166 
    167   /*
    168    * Now fill the last N/4 values : exp^(-j*2*PI*k/N);  k=1,3,5,...,N/2-1.
    169    * These are used for the final twiddle fix-up for converting complex to
    170    * real FFT.
    171    */
    172 
    173   M = N >> 3;
    174   diff = 12 - order;
    175   step = 1 << diff;
    176 
    177   pTwiddle1 = pTwiddle + 3 * N / 8;
    178   pTwiddle4 = pTwiddle1 + (N / 4 - 1);
    179   pTwiddle3 = pTwiddle1 + N / 8;
    180   pTwiddle2 = pTwiddle1 + (N / 8 - 1);
    181 
    182   xS32 = armSP_FFT_S32TwiddleTable[0];
    183   yS32 = armSP_FFT_S32TwiddleTable[1];
    184   x = (xS32 + 0x8000) >> 16;
    185   y = (yS32 + 0x8000) >> 16;
    186   xNeg = 0x7FFF;
    187 
    188   if((order) >= 3) {
    189     for (i = 1; i <= M; i += 2 ) {
    190       OMX_S16 x_neg = 0, y_neg = 0;
    191 
    192       j = i*step;
    193 
    194       xS32 = armSP_FFT_S32TwiddleTable[2 * j];
    195       yS32 = armSP_FFT_S32TwiddleTable[2 * j + 1];
    196       x = (xS32 + 0x8000) >> 16;
    197       y = (yS32 + 0x8000) >> 16;
    198       /* |x_neg = -x| doesn't work when x is 0x8000. */
    199       x_neg = (-(xS32 + 0x8000)) >> 16;
    200       y_neg = (-(yS32 + 0x8000)) >> 16;
    201 
    202       pTwiddle1[0].Re = x;
    203       pTwiddle1[0].Im = y;
    204       pTwiddle1 += 1;
    205       pTwiddle2[0].Re = y_neg;
    206       pTwiddle2[0].Im = x_neg;
    207       pTwiddle2 -= 1;
    208       pTwiddle3[0].Re = y;
    209       pTwiddle3[0].Im = x_neg;
    210       pTwiddle3 += 1;
    211       pTwiddle4[0].Re = x_neg;
    212       pTwiddle4[0].Im = y;
    213       pTwiddle4 -= 1;
    214     }
    215   }
    216   else {
    217     if (order == 2) {
    218       pTwiddle1[0].Re = -y;
    219       pTwiddle1[0].Im = xNeg;
    220     }
    221   }
    222 
    223   /* Update the structure */
    224   pFFTStruct->N = N;
    225   pFFTStruct->pTwiddle = pTwiddle;
    226   pFFTStruct->pBitRev = pBitRev;
    227   pFFTStruct->pBuf = pBuf;
    228 
    229   return OMX_Sts_NoErr;
    230 }
    231 /*****************************************************************************
    232  *                              END OF FILE
    233  *****************************************************************************/
    234 
    235