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