1  /* ----------------------------------------------------------------------
2   * Project:      CMSIS DSP Library
3   * Title:        arm_rfft_f32.c
4   * Description:  RFFT & RIFFT Floating point process function
5   *
6   * $Date:        23 April 2021
7   * $Revision:    V1.9.0
8   *
9   * Target Processor: Cortex-M and Cortex-A cores
10   * -------------------------------------------------------------------- */
11  /*
12   * Copyright (C) 2010-2021 ARM Limited or its affiliates. All rights reserved.
13   *
14   * SPDX-License-Identifier: Apache-2.0
15   *
16   * Licensed under the Apache License, Version 2.0 (the License); you may
17   * not use this file except in compliance with the License.
18   * You may obtain a copy of the License at
19   *
20   * www.apache.org/licenses/LICENSE-2.0
21   *
22   * Unless required by applicable law or agreed to in writing, software
23   * distributed under the License is distributed on an AS IS BASIS, WITHOUT
24   * WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
25   * See the License for the specific language governing permissions and
26   * limitations under the License.
27   */
28  
29  #include "dsp/transform_functions.h"
30  
31  /* ----------------------------------------------------------------------
32   * Internal functions prototypes
33   * -------------------------------------------------------------------- */
34  
35  extern void arm_radix4_butterfly_f32(
36          float32_t * pSrc,
37          uint16_t fftLen,
38    const float32_t * pCoef,
39          uint16_t twidCoefModifier);
40  
41  extern void arm_radix4_butterfly_inverse_f32(
42          float32_t * pSrc,
43          uint16_t fftLen,
44    const float32_t * pCoef,
45          uint16_t twidCoefModifier,
46          float32_t onebyfftLen);
47  
48  extern void arm_bitreversal_f32(
49          float32_t * pSrc,
50          uint16_t fftSize,
51          uint16_t bitRevFactor,
52    const uint16_t * pBitRevTab);
53  
54  void arm_split_rfft_f32(
55          float32_t * pSrc,
56          uint32_t fftLen,
57    const float32_t * pATable,
58    const float32_t * pBTable,
59          float32_t * pDst,
60          uint32_t modifier);
61  
62  void arm_split_rifft_f32(
63          float32_t * pSrc,
64          uint32_t fftLen,
65    const float32_t * pATable,
66    const float32_t * pBTable,
67          float32_t * pDst,
68          uint32_t modifier);
69  
70  /**
71    @ingroup groupTransforms
72   */
73  
74  /**
75    @addtogroup RealFFT
76    @{
77   */
78  
79  /**
80    @brief         Processing function for the floating-point RFFT/RIFFT.
81                   Source buffer is modified by this function.
82  
83    @deprecated    Do not use this function.  It has been superceded by \ref arm_rfft_fast_f32 and will be removed in the future.
84    @param[in]     S    points to an instance of the floating-point RFFT/RIFFT structure
85    @param[in]     pSrc points to the input buffer
86    @param[out]    pDst points to the output buffer
87    @return        none
88  
89    @par
90                     For the RIFFT, the source buffer must at least have length
91                     fftLenReal + 2.
92                     The last two elements must be equal to what would be generated
93                     by the RFFT:
94                       (pSrc[0] - pSrc[1]) and 0.0f
95   */
96  
arm_rfft_f32(const arm_rfft_instance_f32 * S,float32_t * pSrc,float32_t * pDst)97  void arm_rfft_f32(
98    const arm_rfft_instance_f32 * S,
99          float32_t * pSrc,
100          float32_t * pDst)
101  {
102    const arm_cfft_radix4_instance_f32 *S_CFFT = S->pCfft;
103  
104    /* Calculation of Real IFFT of input */
105    if (S->ifftFlagR == 1U)
106    {
107       /*  Real IFFT core process */
108       arm_split_rifft_f32 (pSrc, S->fftLenBy2, S->pTwiddleAReal, S->pTwiddleBReal, pDst, S->twidCoefRModifier);
109  
110  
111       /* Complex radix-4 IFFT process */
112       arm_radix4_butterfly_inverse_f32 (pDst, S_CFFT->fftLen, S_CFFT->pTwiddle, S_CFFT->twidCoefModifier, S_CFFT->onebyfftLen);
113  
114      /* Bit reversal process */
115      if (S->bitReverseFlagR == 1U)
116      {
117        arm_bitreversal_f32 (pDst, S_CFFT->fftLen, S_CFFT->bitRevFactor, S_CFFT->pBitRevTable);
118      }
119    }
120    else
121    {
122      /* Calculation of RFFT of input */
123  
124      /* Complex radix-4 FFT process */
125      arm_radix4_butterfly_f32 (pSrc, S_CFFT->fftLen, S_CFFT->pTwiddle, S_CFFT->twidCoefModifier);
126  
127      /* Bit reversal process */
128      if (S->bitReverseFlagR == 1U)
129      {
130        arm_bitreversal_f32 (pSrc, S_CFFT->fftLen, S_CFFT->bitRevFactor, S_CFFT->pBitRevTable);
131      }
132  
133      /*  Real FFT core process */
134      arm_split_rfft_f32 (pSrc, S->fftLenBy2, S->pTwiddleAReal, S->pTwiddleBReal, pDst, S->twidCoefRModifier);
135    }
136  
137  }
138  
139  /**
140    @} end of RealFFT group
141   */
142  
143  /**
144    @brief         Core Real FFT process
145    @param[in]     pSrc      points to input buffer
146    @param[in]     fftLen    length of FFT
147    @param[in]     pATable   points to twiddle Coef A buffer
148    @param[in]     pBTable   points to twiddle Coef B buffer
149    @param[out]    pDst      points to output buffer
150    @param[in]     modifier  twiddle coefficient modifier that supports different size FFTs with the same twiddle factor table
151    @return        none
152   */
153  
arm_split_rfft_f32(float32_t * pSrc,uint32_t fftLen,const float32_t * pATable,const float32_t * pBTable,float32_t * pDst,uint32_t modifier)154  void arm_split_rfft_f32(
155          float32_t * pSrc,
156          uint32_t fftLen,
157    const float32_t * pATable,
158    const float32_t * pBTable,
159          float32_t * pDst,
160          uint32_t modifier)
161  {
162          uint32_t i;                                    /* Loop Counter */
163          float32_t outR, outI;                          /* Temporary variables for output */
164    const float32_t *pCoefA, *pCoefB;                    /* Temporary pointers for twiddle factors */
165          float32_t CoefA1, CoefA2, CoefB1;              /* Temporary variables for twiddle coefficients */
166          float32_t *pDst1 = &pDst[2], *pDst2 = &pDst[(4U * fftLen) - 1U];      /* temp pointers for output buffer */
167          float32_t *pSrc1 = &pSrc[2], *pSrc2 = &pSrc[(2U * fftLen) - 1U];      /* temp pointers for input buffer */
168  
169    /* Init coefficient pointers */
170    pCoefA = &pATable[modifier * 2];
171    pCoefB = &pBTable[modifier * 2];
172  
173    i = fftLen - 1U;
174  
175    while (i > 0U)
176    {
177       /*
178         outR = (  pSrc[2 * i]             * pATable[2 * i]
179                 - pSrc[2 * i + 1]         * pATable[2 * i + 1]
180                 + pSrc[2 * n - 2 * i]     * pBTable[2 * i]
181                 + pSrc[2 * n - 2 * i + 1] * pBTable[2 * i + 1]);
182  
183         outI = (  pIn[2 * i + 1]         * pATable[2 * i]
184                 + pIn[2 * i]             * pATable[2 * i + 1]
185                 + pIn[2 * n - 2 * i]     * pBTable[2 * i + 1]
186                 - pIn[2 * n - 2 * i + 1] * pBTable[2 * i]);
187        */
188  
189      /* read pATable[2 * i] */
190      CoefA1 = *pCoefA++;
191      /* pATable[2 * i + 1] */
192      CoefA2 = *pCoefA;
193  
194      /* pSrc[2 * i] * pATable[2 * i] */
195      outR = *pSrc1 * CoefA1;
196      /* pSrc[2 * i] * CoefA2 */
197      outI = *pSrc1++ * CoefA2;
198  
199      /* (pSrc[2 * i + 1] + pSrc[2 * fftLen - 2 * i + 1]) * CoefA2 */
200      outR -= (*pSrc1 + *pSrc2) * CoefA2;
201      /* pSrc[2 * i + 1] * CoefA1 */
202      outI += *pSrc1++ * CoefA1;
203  
204      CoefB1 = *pCoefB;
205  
206      /* pSrc[2 * fftLen - 2 * i + 1] * CoefB1 */
207      outI -= *pSrc2-- * CoefB1;
208      /* pSrc[2 * fftLen - 2 * i] * CoefA2 */
209      outI -= *pSrc2 * CoefA2;
210  
211      /* pSrc[2 * fftLen - 2 * i] * CoefB1 */
212      outR += *pSrc2-- * CoefB1;
213  
214      /* write output */
215      *pDst1++ = outR;
216      *pDst1++ = outI;
217  
218      /* write complex conjugate output */
219      *pDst2-- = -outI;
220      *pDst2-- = outR;
221  
222      /* update coefficient pointer */
223      pCoefB = pCoefB + (modifier * 2U);
224      pCoefA = pCoefA + ((modifier * 2U) - 1U);
225  
226      i--;
227  
228    }
229  
230    pDst[2U * fftLen] = pSrc[0] - pSrc[1];
231    pDst[(2U * fftLen) + 1U] = 0.0f;
232  
233    pDst[0] = pSrc[0] + pSrc[1];
234    pDst[1] = 0.0f;
235  
236  }
237  
238  
239  /**
240    @brief         Core Real IFFT process
241    @param[in]     pSrc      points to input buffer
242    @param[in]     fftLen    length of FFT
243    @param[in]     pATable   points to twiddle Coef A buffer
244    @param[in]     pBTable   points to twiddle Coef B buffer
245    @param[out]    pDst      points to output buffer
246    @param[in]     modifier  twiddle coefficient modifier that supports different size FFTs with the same twiddle factor table
247    @return        none
248   */
249  
arm_split_rifft_f32(float32_t * pSrc,uint32_t fftLen,const float32_t * pATable,const float32_t * pBTable,float32_t * pDst,uint32_t modifier)250  void arm_split_rifft_f32(
251          float32_t * pSrc,
252          uint32_t fftLen,
253    const float32_t * pATable,
254    const float32_t * pBTable,
255          float32_t * pDst,
256          uint32_t modifier)
257  {
258          float32_t outR, outI;                          /* Temporary variables for output */
259    const float32_t *pCoefA, *pCoefB;                    /* Temporary pointers for twiddle factors */
260          float32_t CoefA1, CoefA2, CoefB1;              /* Temporary variables for twiddle coefficients */
261          float32_t *pSrc1 = &pSrc[0], *pSrc2 = &pSrc[(2U * fftLen) + 1U];
262  
263    pCoefA = &pATable[0];
264    pCoefB = &pBTable[0];
265  
266    while (fftLen > 0U)
267    {
268       /*
269         outR = (  pIn[2 * i]             * pATable[2 * i]
270                 + pIn[2 * i + 1]         * pATable[2 * i + 1]
271                 + pIn[2 * n - 2 * i]     * pBTable[2 * i]
272                 - pIn[2 * n - 2 * i + 1] * pBTable[2 * i + 1]);
273  
274         outI = (  pIn[2 * i + 1]         * pATable[2 * i]
275                 - pIn[2 * i]             * pATable[2 * i + 1]
276                 - pIn[2 * n - 2 * i]     * pBTable[2 * i + 1]
277                 - pIn[2 * n - 2 * i + 1] * pBTable[2 * i]);
278        */
279  
280       CoefA1 = *pCoefA++;
281       CoefA2 = *pCoefA;
282  
283       /* outR = (pSrc[2 * i] * CoefA1 */
284       outR = *pSrc1 * CoefA1;
285  
286       /* - pSrc[2 * i] * CoefA2 */
287       outI = -(*pSrc1++) * CoefA2;
288  
289       /* (pSrc[2 * i + 1] + pSrc[2 * fftLen - 2 * i + 1]) * CoefA2 */
290       outR += (*pSrc1 + *pSrc2) * CoefA2;
291  
292       /* pSrc[2 * i + 1] * CoefA1 */
293       outI += (*pSrc1++) * CoefA1;
294  
295       CoefB1 = *pCoefB;
296  
297       /* - pSrc[2 * fftLen - 2 * i + 1] * CoefB1 */
298       outI -= *pSrc2-- * CoefB1;
299  
300       /* pSrc[2 * fftLen - 2 * i] * CoefB1 */
301       outR += *pSrc2 * CoefB1;
302  
303       /* pSrc[2 * fftLen - 2 * i] * CoefA2 */
304       outI += *pSrc2-- * CoefA2;
305  
306       /* write output */
307       *pDst++ = outR;
308       *pDst++ = outI;
309  
310       /* update coefficient pointer */
311       pCoefB = pCoefB + (modifier * 2);
312       pCoefA = pCoefA + (modifier * 2 - 1);
313  
314       /* Decrement loop count */
315       fftLen--;
316    }
317  
318  }
319