arm_cfft_radix2_f32.c 11 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471
  1. /* ----------------------------------------------------------------------
  2. * Project: CMSIS DSP Library
  3. * Title: arm_cfft_radix2_f32.c
  4. * Description: Radix-2 Decimation in Frequency CFFT & CIFFT Floating point processing function
  5. *
  6. * $Date: 18. March 2019
  7. * $Revision: V1.6.0
  8. *
  9. * Target Processor: Cortex-M cores
  10. * -------------------------------------------------------------------- */
  11. /*
  12. * Copyright (C) 2010-2019 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. #include "arm_math.h"
  29. void arm_radix2_butterfly_f32(
  30. float32_t * pSrc,
  31. uint32_t fftLen,
  32. const float32_t * pCoef,
  33. uint16_t twidCoefModifier);
  34. void arm_radix2_butterfly_inverse_f32(
  35. float32_t * pSrc,
  36. uint32_t fftLen,
  37. const float32_t * pCoef,
  38. uint16_t twidCoefModifier,
  39. float32_t onebyfftLen);
  40. extern void arm_bitreversal_f32(
  41. float32_t * pSrc,
  42. uint16_t fftSize,
  43. uint16_t bitRevFactor,
  44. const uint16_t * pBitRevTab);
  45. /**
  46. @ingroup groupTransforms
  47. */
  48. /**
  49. @addtogroup ComplexFFT
  50. @{
  51. */
  52. /**
  53. @brief Radix-2 CFFT/CIFFT.
  54. @deprecated Do not use this function. It has been superseded by \ref arm_cfft_f32 and will be removed in the future
  55. @param[in] S points to an instance of the floating-point Radix-2 CFFT/CIFFT structure
  56. @param[in,out] pSrc points to the complex data buffer of size <code>2*fftLen</code>. Processing occurs in-place
  57. @return none
  58. */
  59. void arm_cfft_radix2_f32(
  60. const arm_cfft_radix2_instance_f32 * S,
  61. float32_t * pSrc)
  62. {
  63. if (S->ifftFlag == 1U)
  64. {
  65. /* Complex IFFT radix-2 */
  66. arm_radix2_butterfly_inverse_f32(pSrc, S->fftLen, S->pTwiddle,
  67. S->twidCoefModifier, S->onebyfftLen);
  68. }
  69. else
  70. {
  71. /* Complex FFT radix-2 */
  72. arm_radix2_butterfly_f32(pSrc, S->fftLen, S->pTwiddle,
  73. S->twidCoefModifier);
  74. }
  75. if (S->bitReverseFlag == 1U)
  76. {
  77. /* Bit Reversal */
  78. arm_bitreversal_f32(pSrc, S->fftLen, S->bitRevFactor, S->pBitRevTable);
  79. }
  80. }
  81. /**
  82. @} end of ComplexFFT group
  83. */
  84. /* ----------------------------------------------------------------------
  85. ** Internal helper function used by the FFTs
  86. ** ------------------------------------------------------------------- */
  87. /**
  88. brief Core function for the floating-point CFFT butterfly process.
  89. param[in,out] pSrc points to in-place buffer of floating-point data type
  90. param[in] fftLen length of the FFT
  91. param[in] pCoef points to twiddle coefficient buffer
  92. param[in] twidCoefModifier twiddle coefficient modifier that supports different size FFTs with the same twiddle factor table
  93. return none
  94. */
  95. void arm_radix2_butterfly_f32(
  96. float32_t * pSrc,
  97. uint32_t fftLen,
  98. const float32_t * pCoef,
  99. uint16_t twidCoefModifier)
  100. {
  101. uint32_t i, j, k, l;
  102. uint32_t n1, n2, ia;
  103. float32_t xt, yt, cosVal, sinVal;
  104. float32_t p0, p1, p2, p3;
  105. float32_t a0, a1;
  106. #if defined (ARM_MATH_DSP)
  107. /* Initializations for the first stage */
  108. n2 = fftLen >> 1;
  109. ia = 0;
  110. i = 0;
  111. // loop for groups
  112. for (k = n2; k > 0; k--)
  113. {
  114. cosVal = pCoef[ia * 2];
  115. sinVal = pCoef[(ia * 2) + 1];
  116. /* Twiddle coefficients index modifier */
  117. ia += twidCoefModifier;
  118. /* index calculation for the input as, */
  119. /* pSrc[i + 0], pSrc[i + fftLen/1] */
  120. l = i + n2;
  121. /* Butterfly implementation */
  122. a0 = pSrc[2 * i] + pSrc[2 * l];
  123. xt = pSrc[2 * i] - pSrc[2 * l];
  124. yt = pSrc[2 * i + 1] - pSrc[2 * l + 1];
  125. a1 = pSrc[2 * l + 1] + pSrc[2 * i + 1];
  126. p0 = xt * cosVal;
  127. p1 = yt * sinVal;
  128. p2 = yt * cosVal;
  129. p3 = xt * sinVal;
  130. pSrc[2 * i] = a0;
  131. pSrc[2 * i + 1] = a1;
  132. pSrc[2 * l] = p0 + p1;
  133. pSrc[2 * l + 1] = p2 - p3;
  134. i++;
  135. } // groups loop end
  136. twidCoefModifier <<= 1U;
  137. // loop for stage
  138. for (k = n2; k > 2; k = k >> 1)
  139. {
  140. n1 = n2;
  141. n2 = n2 >> 1;
  142. ia = 0;
  143. // loop for groups
  144. j = 0;
  145. do
  146. {
  147. cosVal = pCoef[ia * 2];
  148. sinVal = pCoef[(ia * 2) + 1];
  149. ia += twidCoefModifier;
  150. // loop for butterfly
  151. i = j;
  152. do
  153. {
  154. l = i + n2;
  155. a0 = pSrc[2 * i] + pSrc[2 * l];
  156. xt = pSrc[2 * i] - pSrc[2 * l];
  157. yt = pSrc[2 * i + 1] - pSrc[2 * l + 1];
  158. a1 = pSrc[2 * l + 1] + pSrc[2 * i + 1];
  159. p0 = xt * cosVal;
  160. p1 = yt * sinVal;
  161. p2 = yt * cosVal;
  162. p3 = xt * sinVal;
  163. pSrc[2 * i] = a0;
  164. pSrc[2 * i + 1] = a1;
  165. pSrc[2 * l] = p0 + p1;
  166. pSrc[2 * l + 1] = p2 - p3;
  167. i += n1;
  168. } while ( i < fftLen ); // butterfly loop end
  169. j++;
  170. } while ( j < n2); // groups loop end
  171. twidCoefModifier <<= 1U;
  172. } // stages loop end
  173. // loop for butterfly
  174. for (i = 0; i < fftLen; i += 2)
  175. {
  176. a0 = pSrc[2 * i] + pSrc[2 * i + 2];
  177. xt = pSrc[2 * i] - pSrc[2 * i + 2];
  178. yt = pSrc[2 * i + 1] - pSrc[2 * i + 3];
  179. a1 = pSrc[2 * i + 3] + pSrc[2 * i + 1];
  180. pSrc[2 * i] = a0;
  181. pSrc[2 * i + 1] = a1;
  182. pSrc[2 * i + 2] = xt;
  183. pSrc[2 * i + 3] = yt;
  184. } // groups loop end
  185. #else /* #if defined (ARM_MATH_DSP) */
  186. n2 = fftLen;
  187. // loop for stage
  188. for (k = fftLen; k > 1; k = k >> 1)
  189. {
  190. n1 = n2;
  191. n2 = n2 >> 1;
  192. ia = 0;
  193. // loop for groups
  194. j = 0;
  195. do
  196. {
  197. cosVal = pCoef[ia * 2];
  198. sinVal = pCoef[(ia * 2) + 1];
  199. ia += twidCoefModifier;
  200. // loop for butterfly
  201. i = j;
  202. do
  203. {
  204. l = i + n2;
  205. a0 = pSrc[2 * i] + pSrc[2 * l];
  206. xt = pSrc[2 * i] - pSrc[2 * l];
  207. yt = pSrc[2 * i + 1] - pSrc[2 * l + 1];
  208. a1 = pSrc[2 * l + 1] + pSrc[2 * i + 1];
  209. p0 = xt * cosVal;
  210. p1 = yt * sinVal;
  211. p2 = yt * cosVal;
  212. p3 = xt * sinVal;
  213. pSrc[2 * i] = a0;
  214. pSrc[2 * i + 1] = a1;
  215. pSrc[2 * l] = p0 + p1;
  216. pSrc[2 * l + 1] = p2 - p3;
  217. i += n1;
  218. } while (i < fftLen);
  219. j++;
  220. } while (j < n2);
  221. twidCoefModifier <<= 1U;
  222. }
  223. #endif /* #if defined (ARM_MATH_DSP) */
  224. }
  225. void arm_radix2_butterfly_inverse_f32(
  226. float32_t * pSrc,
  227. uint32_t fftLen,
  228. const float32_t * pCoef,
  229. uint16_t twidCoefModifier,
  230. float32_t onebyfftLen)
  231. {
  232. uint32_t i, j, k, l;
  233. uint32_t n1, n2, ia;
  234. float32_t xt, yt, cosVal, sinVal;
  235. float32_t p0, p1, p2, p3;
  236. float32_t a0, a1;
  237. #if defined (ARM_MATH_DSP)
  238. n2 = fftLen >> 1;
  239. ia = 0;
  240. // loop for groups
  241. for (i = 0; i < n2; i++)
  242. {
  243. cosVal = pCoef[ia * 2];
  244. sinVal = pCoef[(ia * 2) + 1];
  245. ia += twidCoefModifier;
  246. l = i + n2;
  247. a0 = pSrc[2 * i] + pSrc[2 * l];
  248. xt = pSrc[2 * i] - pSrc[2 * l];
  249. yt = pSrc[2 * i + 1] - pSrc[2 * l + 1];
  250. a1 = pSrc[2 * l + 1] + pSrc[2 * i + 1];
  251. p0 = xt * cosVal;
  252. p1 = yt * sinVal;
  253. p2 = yt * cosVal;
  254. p3 = xt * sinVal;
  255. pSrc[2 * i] = a0;
  256. pSrc[2 * i + 1] = a1;
  257. pSrc[2 * l] = p0 - p1;
  258. pSrc[2 * l + 1] = p2 + p3;
  259. } // groups loop end
  260. twidCoefModifier <<= 1U;
  261. // loop for stage
  262. for (k = fftLen / 2; k > 2; k = k >> 1)
  263. {
  264. n1 = n2;
  265. n2 = n2 >> 1;
  266. ia = 0;
  267. // loop for groups
  268. j = 0;
  269. do
  270. {
  271. cosVal = pCoef[ia * 2];
  272. sinVal = pCoef[(ia * 2) + 1];
  273. ia += twidCoefModifier;
  274. // loop for butterfly
  275. i = j;
  276. do
  277. {
  278. l = i + n2;
  279. a0 = pSrc[2 * i] + pSrc[2 * l];
  280. xt = pSrc[2 * i] - pSrc[2 * l];
  281. yt = pSrc[2 * i + 1] - pSrc[2 * l + 1];
  282. a1 = pSrc[2 * l + 1] + pSrc[2 * i + 1];
  283. p0 = xt * cosVal;
  284. p1 = yt * sinVal;
  285. p2 = yt * cosVal;
  286. p3 = xt * sinVal;
  287. pSrc[2 * i] = a0;
  288. pSrc[2 * i + 1] = a1;
  289. pSrc[2 * l] = p0 - p1;
  290. pSrc[2 * l + 1] = p2 + p3;
  291. i += n1;
  292. } while ( i < fftLen ); // butterfly loop end
  293. j++;
  294. } while (j < n2); // groups loop end
  295. twidCoefModifier <<= 1U;
  296. } // stages loop end
  297. // loop for butterfly
  298. for (i = 0; i < fftLen; i += 2)
  299. {
  300. a0 = pSrc[2 * i] + pSrc[2 * i + 2];
  301. xt = pSrc[2 * i] - pSrc[2 * i + 2];
  302. a1 = pSrc[2 * i + 3] + pSrc[2 * i + 1];
  303. yt = pSrc[2 * i + 1] - pSrc[2 * i + 3];
  304. p0 = a0 * onebyfftLen;
  305. p2 = xt * onebyfftLen;
  306. p1 = a1 * onebyfftLen;
  307. p3 = yt * onebyfftLen;
  308. pSrc[2 * i] = p0;
  309. pSrc[2 * i + 1] = p1;
  310. pSrc[2 * i + 2] = p2;
  311. pSrc[2 * i + 3] = p3;
  312. } // butterfly loop end
  313. #else /* #if defined (ARM_MATH_DSP) */
  314. n2 = fftLen;
  315. // loop for stage
  316. for (k = fftLen; k > 2; k = k >> 1)
  317. {
  318. n1 = n2;
  319. n2 = n2 >> 1;
  320. ia = 0;
  321. // loop for groups
  322. j = 0;
  323. do
  324. {
  325. cosVal = pCoef[ia * 2];
  326. sinVal = pCoef[(ia * 2) + 1];
  327. ia = ia + twidCoefModifier;
  328. // loop for butterfly
  329. i = j;
  330. do
  331. {
  332. l = i + n2;
  333. a0 = pSrc[2 * i] + pSrc[2 * l];
  334. xt = pSrc[2 * i] - pSrc[2 * l];
  335. yt = pSrc[2 * i + 1] - pSrc[2 * l + 1];
  336. a1 = pSrc[2 * l + 1] + pSrc[2 * i + 1];
  337. p0 = xt * cosVal;
  338. p1 = yt * sinVal;
  339. p2 = yt * cosVal;
  340. p3 = xt * sinVal;
  341. pSrc[2 * i] = a0;
  342. pSrc[2 * i + 1] = a1;
  343. pSrc[2 * l] = p0 - p1;
  344. pSrc[2 * l + 1] = p2 + p3;
  345. i += n1;
  346. } while ( i < fftLen ); // butterfly loop end
  347. j++;
  348. } while ( j < n2 ); // groups loop end
  349. twidCoefModifier = twidCoefModifier << 1U;
  350. } // stages loop end
  351. n1 = n2;
  352. n2 = n2 >> 1;
  353. // loop for butterfly
  354. for (i = 0; i < fftLen; i += n1)
  355. {
  356. l = i + n2;
  357. a0 = pSrc[2 * i] + pSrc[2 * l];
  358. xt = pSrc[2 * i] - pSrc[2 * l];
  359. a1 = pSrc[2 * l + 1] + pSrc[2 * i + 1];
  360. yt = pSrc[2 * i + 1] - pSrc[2 * l + 1];
  361. p0 = a0 * onebyfftLen;
  362. p2 = xt * onebyfftLen;
  363. p1 = a1 * onebyfftLen;
  364. p3 = yt * onebyfftLen;
  365. pSrc[2 * i] = p0;
  366. pSrc[2 * l] = p2;
  367. pSrc[2 * i + 1] = p1;
  368. pSrc[2 * l + 1] = p3;
  369. } // butterfly loop end
  370. #endif /* #if defined (ARM_MATH_DSP) */
  371. }