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 
72 /**
73   @addtogroup DeprecatedRealFFT
74   @{
75  */
76 
77 /**
78   @brief         Processing function for the floating-point RFFT/RIFFT.
79                  Source buffer is modified by this function.
80 
81   @deprecated    Do not use this function.  It has been superceded by \ref arm_rfft_fast_f32 and will be removed in the future.
82   @param[in]     S    points to an instance of the floating-point RFFT/RIFFT structure
83   @param[in]     pSrc points to the input buffer
84   @param[out]    pDst points to the output buffer
85 
86  */
87 
arm_rfft_f32(const arm_rfft_instance_f32 * S,float32_t * pSrc,float32_t * pDst)88 void arm_rfft_f32(
89   const arm_rfft_instance_f32 * S,
90         float32_t * pSrc,
91         float32_t * pDst)
92 {
93   const arm_cfft_radix4_instance_f32 *S_CFFT = S->pCfft;
94 
95   /* Calculation of Real IFFT of input */
96   if (S->ifftFlagR == 1U)
97   {
98      /*  Real IFFT core process */
99      arm_split_rifft_f32 (pSrc, S->fftLenBy2, S->pTwiddleAReal, S->pTwiddleBReal, pDst, S->twidCoefRModifier);
100 
101 
102      /* Complex radix-4 IFFT process */
103      arm_radix4_butterfly_inverse_f32 (pDst, S_CFFT->fftLen, S_CFFT->pTwiddle, S_CFFT->twidCoefModifier, S_CFFT->onebyfftLen);
104 
105     /* Bit reversal process */
106     if (S->bitReverseFlagR == 1U)
107     {
108       arm_bitreversal_f32 (pDst, S_CFFT->fftLen, S_CFFT->bitRevFactor, S_CFFT->pBitRevTable);
109     }
110   }
111   else
112   {
113     /* Calculation of RFFT of input */
114 
115     /* Complex radix-4 FFT process */
116     arm_radix4_butterfly_f32 (pSrc, S_CFFT->fftLen, S_CFFT->pTwiddle, S_CFFT->twidCoefModifier);
117 
118     /* Bit reversal process */
119     if (S->bitReverseFlagR == 1U)
120     {
121       arm_bitreversal_f32 (pSrc, S_CFFT->fftLen, S_CFFT->bitRevFactor, S_CFFT->pBitRevTable);
122     }
123 
124     /*  Real FFT core process */
125     arm_split_rfft_f32 (pSrc, S->fftLenBy2, S->pTwiddleAReal, S->pTwiddleBReal, pDst, S->twidCoefRModifier);
126   }
127 
128 }
129 
130 /**
131   @} end of DeprecatedRealFFT group
132  */
133 
134 /**
135   @brief         Core Real FFT process
136   @param[in]     pSrc      points to input buffer
137   @param[in]     fftLen    length of FFT
138   @param[in]     pATable   points to twiddle Coef A buffer
139   @param[in]     pBTable   points to twiddle Coef B buffer
140   @param[out]    pDst      points to output buffer
141   @param[in]     modifier  twiddle coefficient modifier that supports different size FFTs with the same twiddle factor table
142  */
143 
arm_split_rfft_f32(float32_t * pSrc,uint32_t fftLen,const float32_t * pATable,const float32_t * pBTable,float32_t * pDst,uint32_t modifier)144 void arm_split_rfft_f32(
145         float32_t * pSrc,
146         uint32_t fftLen,
147   const float32_t * pATable,
148   const float32_t * pBTable,
149         float32_t * pDst,
150         uint32_t modifier)
151 {
152         uint32_t i;                                    /* Loop Counter */
153         float32_t outR, outI;                          /* Temporary variables for output */
154   const float32_t *pCoefA, *pCoefB;                    /* Temporary pointers for twiddle factors */
155         float32_t CoefA1, CoefA2, CoefB1;              /* Temporary variables for twiddle coefficients */
156         float32_t *pDst1 = &pDst[2], *pDst2 = &pDst[(4U * fftLen) - 1U];      /* temp pointers for output buffer */
157         float32_t *pSrc1 = &pSrc[2], *pSrc2 = &pSrc[(2U * fftLen) - 1U];      /* temp pointers for input buffer */
158 
159   /* Init coefficient pointers */
160   pCoefA = &pATable[modifier * 2];
161   pCoefB = &pBTable[modifier * 2];
162 
163   i = fftLen - 1U;
164 
165   while (i > 0U)
166   {
167      /*
168        outR = (  pSrc[2 * i]             * pATable[2 * i]
169                - pSrc[2 * i + 1]         * pATable[2 * i + 1]
170                + pSrc[2 * n - 2 * i]     * pBTable[2 * i]
171                + pSrc[2 * n - 2 * i + 1] * pBTable[2 * i + 1]);
172 
173        outI = (  pIn[2 * i + 1]         * pATable[2 * i]
174                + pIn[2 * i]             * pATable[2 * i + 1]
175                + pIn[2 * n - 2 * i]     * pBTable[2 * i + 1]
176                - pIn[2 * n - 2 * i + 1] * pBTable[2 * i]);
177       */
178 
179     /* read pATable[2 * i] */
180     CoefA1 = *pCoefA++;
181     /* pATable[2 * i + 1] */
182     CoefA2 = *pCoefA;
183 
184     /* pSrc[2 * i] * pATable[2 * i] */
185     outR = *pSrc1 * CoefA1;
186     /* pSrc[2 * i] * CoefA2 */
187     outI = *pSrc1++ * CoefA2;
188 
189     /* (pSrc[2 * i + 1] + pSrc[2 * fftLen - 2 * i + 1]) * CoefA2 */
190     outR -= (*pSrc1 + *pSrc2) * CoefA2;
191     /* pSrc[2 * i + 1] * CoefA1 */
192     outI += *pSrc1++ * CoefA1;
193 
194     CoefB1 = *pCoefB;
195 
196     /* pSrc[2 * fftLen - 2 * i + 1] * CoefB1 */
197     outI -= *pSrc2-- * CoefB1;
198     /* pSrc[2 * fftLen - 2 * i] * CoefA2 */
199     outI -= *pSrc2 * CoefA2;
200 
201     /* pSrc[2 * fftLen - 2 * i] * CoefB1 */
202     outR += *pSrc2-- * CoefB1;
203 
204     /* write output */
205     *pDst1++ = outR;
206     *pDst1++ = outI;
207 
208     /* write complex conjugate output */
209     *pDst2-- = -outI;
210     *pDst2-- = outR;
211 
212     /* update coefficient pointer */
213     pCoefB = pCoefB + (modifier * 2U);
214     pCoefA = pCoefA + ((modifier * 2U) - 1U);
215 
216     i--;
217 
218   }
219 
220   pDst[2U * fftLen] = pSrc[0] - pSrc[1];
221   pDst[(2U * fftLen) + 1U] = 0.0f;
222 
223   pDst[0] = pSrc[0] + pSrc[1];
224   pDst[1] = 0.0f;
225 
226 }
227 
228 
229 /**
230   @brief         Core Real IFFT process
231   @param[in]     pSrc      points to input buffer
232   @param[in]     fftLen    length of FFT
233   @param[in]     pATable   points to twiddle Coef A buffer
234   @param[in]     pBTable   points to twiddle Coef B buffer
235   @param[out]    pDst      points to output buffer
236   @param[in]     modifier  twiddle coefficient modifier that supports different size FFTs with the same twiddle factor table
237  */
238 
arm_split_rifft_f32(float32_t * pSrc,uint32_t fftLen,const float32_t * pATable,const float32_t * pBTable,float32_t * pDst,uint32_t modifier)239 void arm_split_rifft_f32(
240         float32_t * pSrc,
241         uint32_t fftLen,
242   const float32_t * pATable,
243   const float32_t * pBTable,
244         float32_t * pDst,
245         uint32_t modifier)
246 {
247         float32_t outR, outI;                          /* Temporary variables for output */
248   const float32_t *pCoefA, *pCoefB;                    /* Temporary pointers for twiddle factors */
249         float32_t CoefA1, CoefA2, CoefB1;              /* Temporary variables for twiddle coefficients */
250         float32_t *pSrc1 = &pSrc[0], *pSrc2 = &pSrc[(2U * fftLen) + 1U];
251 
252   pCoefA = &pATable[0];
253   pCoefB = &pBTable[0];
254 
255   while (fftLen > 0U)
256   {
257      /*
258        outR = (  pIn[2 * i]             * pATable[2 * i]
259                + pIn[2 * i + 1]         * pATable[2 * i + 1]
260                + pIn[2 * n - 2 * i]     * pBTable[2 * i]
261                - pIn[2 * n - 2 * i + 1] * pBTable[2 * i + 1]);
262 
263        outI = (  pIn[2 * i + 1]         * pATable[2 * i]
264                - pIn[2 * i]             * pATable[2 * i + 1]
265                - pIn[2 * n - 2 * i]     * pBTable[2 * i + 1]
266                - pIn[2 * n - 2 * i + 1] * pBTable[2 * i]);
267       */
268 
269      CoefA1 = *pCoefA++;
270      CoefA2 = *pCoefA;
271 
272      /* outR = (pSrc[2 * i] * CoefA1 */
273      outR = *pSrc1 * CoefA1;
274 
275      /* - pSrc[2 * i] * CoefA2 */
276      outI = -(*pSrc1++) * CoefA2;
277 
278      /* (pSrc[2 * i + 1] + pSrc[2 * fftLen - 2 * i + 1]) * CoefA2 */
279      outR += (*pSrc1 + *pSrc2) * CoefA2;
280 
281      /* pSrc[2 * i + 1] * CoefA1 */
282      outI += (*pSrc1++) * CoefA1;
283 
284      CoefB1 = *pCoefB;
285 
286      /* - pSrc[2 * fftLen - 2 * i + 1] * CoefB1 */
287      outI -= *pSrc2-- * CoefB1;
288 
289      /* pSrc[2 * fftLen - 2 * i] * CoefB1 */
290      outR += *pSrc2 * CoefB1;
291 
292      /* pSrc[2 * fftLen - 2 * i] * CoefA2 */
293      outI += *pSrc2-- * CoefA2;
294 
295      /* write output */
296      *pDst++ = outR;
297      *pDst++ = outI;
298 
299      /* update coefficient pointer */
300      pCoefB = pCoefB + (modifier * 2);
301      pCoefA = pCoefA + (modifier * 2 - 1);
302 
303      /* Decrement loop count */
304      fftLen--;
305   }
306 
307 }
308