/* ---------------------------------------------------------------------- * Copyright (C) 2010-2013 ARM Limited. All rights reserved. * * $Date: 17. January 2013 * $Revision: V1.4.1 * * Project: CMSIS DSP Library * Title: arm_cfft_radix4_f32.c * * Description: Radix-4 Decimation in Frequency CFFT & CIFFT Floating point processing function * * * Target Processor: Cortex-M4/Cortex-M3/Cortex-M0 * * Redistribution and use in source and binary forms, with or without * modification, are permitted provided that the following conditions * are met: * - Redistributions of source code must retain the above copyright * notice, this list of conditions and the following disclaimer. * - Redistributions in binary form must reproduce the above copyright * notice, this list of conditions and the following disclaimer in * the documentation and/or other materials provided with the * distribution. * - Neither the name of ARM LIMITED nor the names of its contributors * may be used to endorse or promote products derived from this * software without specific prior written permission. * * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS * FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE * COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, * BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER * CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN * ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE * POSSIBILITY OF SUCH DAMAGE. * -------------------------------------------------------------------- */ #include "arm_math.h" extern void arm_bitreversal_f32( float32_t * pSrc, uint16_t fftSize, uint16_t bitRevFactor, uint16_t * pBitRevTab); /** * @ingroup groupTransforms */ /* ---------------------------------------------------------------------- ** Internal helper function used by the FFTs ** ------------------------------------------------------------------- */ /* * @brief Core function for the floating-point CFFT butterfly process. * @param[in, out] *pSrc points to the in-place buffer of floating-point data type. * @param[in] fftLen length of the FFT. * @param[in] *pCoef points to the twiddle coefficient buffer. * @param[in] twidCoefModifier twiddle coefficient modifier that supports different size FFTs with the same twiddle factor table. * @return none. */ void arm_radix4_butterfly_f32( float32_t * pSrc, uint16_t fftLen, float32_t * pCoef, uint16_t twidCoefModifier) { float32_t co1, co2, co3, si1, si2, si3; uint32_t ia1, ia2, ia3; uint32_t i0, i1, i2, i3; uint32_t n1, n2, j, k; #ifndef ARM_MATH_CM0_FAMILY_FAMILY /* Run the below code for Cortex-M4 and Cortex-M3 */ float32_t xaIn, yaIn, xbIn, ybIn, xcIn, ycIn, xdIn, ydIn; float32_t Xaplusc, Xbplusd, Yaplusc, Ybplusd, Xaminusc, Xbminusd, Yaminusc, Ybminusd; float32_t Xb12C_out, Yb12C_out, Xc12C_out, Yc12C_out, Xd12C_out, Yd12C_out; float32_t Xb12_out, Yb12_out, Xc12_out, Yc12_out, Xd12_out, Yd12_out; float32_t *ptr1; float32_t p0,p1,p2,p3,p4,p5; float32_t a0,a1,a2,a3,a4,a5,a6,a7; /* Initializations for the first stage */ n2 = fftLen; n1 = n2; /* n2 = fftLen/4 */ n2 >>= 2u; i0 = 0u; ia1 = 0u; j = n2; /* Calculation of first stage */ do { /* index calculation for the input as, */ /* pSrc[i0 + 0], pSrc[i0 + fftLen/4], pSrc[i0 + fftLen/2], pSrc[i0 + 3fftLen/4] */ i1 = i0 + n2; i2 = i1 + n2; i3 = i2 + n2; xaIn = pSrc[(2u * i0)]; yaIn = pSrc[(2u * i0) + 1u]; xbIn = pSrc[(2u * i1)]; ybIn = pSrc[(2u * i1) + 1u]; xcIn = pSrc[(2u * i2)]; ycIn = pSrc[(2u * i2) + 1u]; xdIn = pSrc[(2u * i3)]; ydIn = pSrc[(2u * i3) + 1u]; /* xa + xc */ Xaplusc = xaIn + xcIn; /* xb + xd */ Xbplusd = xbIn + xdIn; /* ya + yc */ Yaplusc = yaIn + ycIn; /* yb + yd */ Ybplusd = ybIn + ydIn; /* index calculation for the coefficients */ ia2 = ia1 + ia1; co2 = pCoef[ia2 * 2u]; si2 = pCoef[(ia2 * 2u) + 1u]; /* xa - xc */ Xaminusc = xaIn - xcIn; /* xb - xd */ Xbminusd = xbIn - xdIn; /* ya - yc */ Yaminusc = yaIn - ycIn; /* yb - yd */ Ybminusd = ybIn - ydIn; /* xa' = xa + xb + xc + xd */ pSrc[(2u * i0)] = Xaplusc + Xbplusd; /* ya' = ya + yb + yc + yd */ pSrc[(2u * i0) + 1u] = Yaplusc + Ybplusd; /* (xa - xc) + (yb - yd) */ Xb12C_out = (Xaminusc + Ybminusd); /* (ya - yc) + (xb - xd) */ Yb12C_out = (Yaminusc - Xbminusd); /* (xa + xc) - (xb + xd) */ Xc12C_out = (Xaplusc - Xbplusd); /* (ya + yc) - (yb + yd) */ Yc12C_out = (Yaplusc - Ybplusd); /* (xa - xc) - (yb - yd) */ Xd12C_out = (Xaminusc - Ybminusd); /* (ya - yc) + (xb - xd) */ Yd12C_out = (Xbminusd + Yaminusc); co1 = pCoef[ia1 * 2u]; si1 = pCoef[(ia1 * 2u) + 1u]; /* index calculation for the coefficients */ ia3 = ia2 + ia1; co3 = pCoef[ia3 * 2u]; si3 = pCoef[(ia3 * 2u) + 1u]; Xb12_out = Xb12C_out * co1; Yb12_out = Yb12C_out * co1; Xc12_out = Xc12C_out * co2; Yc12_out = Yc12C_out * co2; Xd12_out = Xd12C_out * co3; Yd12_out = Yd12C_out * co3; /* xb' = (xa+yb-xc-yd)co1 - (ya-xb-yc+xd)(si1) */ //Xb12_out -= Yb12C_out * si1; p0 = Yb12C_out * si1; /* yb' = (ya-xb-yc+xd)co1 + (xa+yb-xc-yd)(si1) */ //Yb12_out += Xb12C_out * si1; p1 = Xb12C_out * si1; /* xc' = (xa-xb+xc-xd)co2 - (ya-yb+yc-yd)(si2) */ //Xc12_out -= Yc12C_out * si2; p2 = Yc12C_out * si2; /* yc' = (ya-yb+yc-yd)co2 + (xa-xb+xc-xd)(si2) */ //Yc12_out += Xc12C_out * si2; p3 = Xc12C_out * si2; /* xd' = (xa-yb-xc+yd)co3 - (ya+xb-yc-xd)(si3) */ //Xd12_out -= Yd12C_out * si3; p4 = Yd12C_out * si3; /* yd' = (ya+xb-yc-xd)co3 + (xa-yb-xc+yd)(si3) */ //Yd12_out += Xd12C_out * si3; p5 = Xd12C_out * si3; Xb12_out += p0; Yb12_out -= p1; Xc12_out += p2; Yc12_out -= p3; Xd12_out += p4; Yd12_out -= p5; /* xc' = (xa-xb+xc-xd)co2 + (ya-yb+yc-yd)(si2) */ pSrc[2u * i1] = Xc12_out; /* yc' = (ya-yb+yc-yd)co2 - (xa-xb+xc-xd)(si2) */ pSrc[(2u * i1) + 1u] = Yc12_out; /* xb' = (xa+yb-xc-yd)co1 + (ya-xb-yc+xd)(si1) */ pSrc[2u * i2] = Xb12_out; /* yb' = (ya-xb-yc+xd)co1 - (xa+yb-xc-yd)(si1) */ pSrc[(2u * i2) + 1u] = Yb12_out; /* xd' = (xa-yb-xc+yd)co3 + (ya+xb-yc-xd)(si3) */ pSrc[2u * i3] = Xd12_out; /* yd' = (ya+xb-yc-xd)co3 - (xa-yb-xc+yd)(si3) */ pSrc[(2u * i3) + 1u] = Yd12_out; /* Twiddle coefficients index modifier */ ia1 += twidCoefModifier; /* Updating input index */ i0++; } while(--j); twidCoefModifier <<= 2u; /* Calculation of second stage to excluding last stage */ for (k = fftLen >> 2u; k > 4u; k >>= 2u) { /* Initializations for the first stage */ n1 = n2; n2 >>= 2u; ia1 = 0u; /* Calculation of first stage */ j = 0; do { /* index calculation for the coefficients */ ia2 = ia1 + ia1; ia3 = ia2 + ia1; co1 = pCoef[ia1 * 2u]; si1 = pCoef[(ia1 * 2u) + 1u]; co2 = pCoef[ia2 * 2u]; si2 = pCoef[(ia2 * 2u) + 1u]; co3 = pCoef[ia3 * 2u]; si3 = pCoef[(ia3 * 2u) + 1u]; /* Twiddle coefficients index modifier */ ia1 += twidCoefModifier; i0 = j; do { /* index calculation for the input as, */ /* pSrc[i0 + 0], pSrc[i0 + fftLen/4], pSrc[i0 + fftLen/2], pSrc[i0 + 3fftLen/4] */ i1 = i0 + n2; i2 = i1 + n2; i3 = i2 + n2; xaIn = pSrc[(2u * i0)]; yaIn = pSrc[(2u * i0) + 1u]; xbIn = pSrc[(2u * i1)]; ybIn = pSrc[(2u * i1) + 1u]; xcIn = pSrc[(2u * i2)]; ycIn = pSrc[(2u * i2) + 1u]; xdIn = pSrc[(2u * i3)]; ydIn = pSrc[(2u * i3) + 1u]; /* xa - xc */ Xaminusc = xaIn - xcIn; /* (xb - xd) */ Xbminusd = xbIn - xdIn; /* ya - yc */ Yaminusc = yaIn - ycIn; /* (yb - yd) */ Ybminusd = ybIn - ydIn; /* xa + xc */ Xaplusc = xaIn + xcIn; /* xb + xd */ Xbplusd = xbIn + xdIn; /* ya + yc */ Yaplusc = yaIn + ycIn; /* yb + yd */ Ybplusd = ybIn + ydIn; /* (xa - xc) + (yb - yd) */ Xb12C_out = (Xaminusc + Ybminusd); /* (ya - yc) - (xb - xd) */ Yb12C_out = (Yaminusc - Xbminusd); /* xa + xc -(xb + xd) */ Xc12C_out = (Xaplusc - Xbplusd); /* (ya + yc) - (yb + yd) */ Yc12C_out = (Yaplusc - Ybplusd); /* (xa - xc) - (yb - yd) */ Xd12C_out = (Xaminusc - Ybminusd); /* (ya - yc) + (xb - xd) */ Yd12C_out = (Xbminusd + Yaminusc); pSrc[(2u * i0)] = Xaplusc + Xbplusd; pSrc[(2u * i0) + 1u] = Yaplusc + Ybplusd; Xb12_out = Xb12C_out * co1; Yb12_out = Yb12C_out * co1; Xc12_out = Xc12C_out * co2; Yc12_out = Yc12C_out * co2; Xd12_out = Xd12C_out * co3; Yd12_out = Yd12C_out * co3; /* xb' = (xa+yb-xc-yd)co1 - (ya-xb-yc+xd)(si1) */ //Xb12_out -= Yb12C_out * si1; p0 = Yb12C_out * si1; /* yb' = (ya-xb-yc+xd)co1 + (xa+yb-xc-yd)(si1) */ //Yb12_out += Xb12C_out * si1; p1 = Xb12C_out * si1; /* xc' = (xa-xb+xc-xd)co2 - (ya-yb+yc-yd)(si2) */ //Xc12_out -= Yc12C_out * si2; p2 = Yc12C_out * si2; /* yc' = (ya-yb+yc-yd)co2 + (xa-xb+xc-xd)(si2) */ //Yc12_out += Xc12C_out * si2; p3 = Xc12C_out * si2; /* xd' = (xa-yb-xc+yd)co3 - (ya+xb-yc-xd)(si3) */ //Xd12_out -= Yd12C_out * si3; p4 = Yd12C_out * si3; /* yd' = (ya+xb-yc-xd)co3 + (xa-yb-xc+yd)(si3) */ //Yd12_out += Xd12C_out * si3; p5 = Xd12C_out * si3; Xb12_out += p0; Yb12_out -= p1; Xc12_out += p2; Yc12_out -= p3; Xd12_out += p4; Yd12_out -= p5; /* xc' = (xa-xb+xc-xd)co2 + (ya-yb+yc-yd)(si2) */ pSrc[2u * i1] = Xc12_out; /* yc' = (ya-yb+yc-yd)co2 - (xa-xb+xc-xd)(si2) */ pSrc[(2u * i1) + 1u] = Yc12_out; /* xb' = (xa+yb-xc-yd)co1 + (ya-xb-yc+xd)(si1) */ pSrc[2u * i2] = Xb12_out; /* yb' = (ya-xb-yc+xd)co1 - (xa+yb-xc-yd)(si1) */ pSrc[(2u * i2) + 1u] = Yb12_out; /* xd' = (xa-yb-xc+yd)co3 + (ya+xb-yc-xd)(si3) */ pSrc[2u * i3] = Xd12_out; /* yd' = (ya+xb-yc-xd)co3 - (xa-yb-xc+yd)(si3) */ pSrc[(2u * i3) + 1u] = Yd12_out; i0 += n1; } while(i0 < fftLen); j++; } while(j <= (n2 - 1u)); twidCoefModifier <<= 2u; } j = fftLen >> 2; ptr1 = &pSrc[0]; /* Calculations of last stage */ do { xaIn = ptr1[0]; yaIn = ptr1[1]; xbIn = ptr1[2]; ybIn = ptr1[3]; xcIn = ptr1[4]; ycIn = ptr1[5]; xdIn = ptr1[6]; ydIn = ptr1[7]; /* xa + xc */ Xaplusc = xaIn + xcIn; /* xa - xc */ Xaminusc = xaIn - xcIn; /* ya + yc */ Yaplusc = yaIn + ycIn; /* ya - yc */ Yaminusc = yaIn - ycIn; /* xb + xd */ Xbplusd = xbIn + xdIn; /* yb + yd */ Ybplusd = ybIn + ydIn; /* (xb-xd) */ Xbminusd = xbIn - xdIn; /* (yb-yd) */ Ybminusd = ybIn - ydIn; /* xa' = xa + xb + xc + xd */ a0 = (Xaplusc + Xbplusd); /* ya' = ya + yb + yc + yd */ a1 = (Yaplusc + Ybplusd); /* xc' = (xa-xb+xc-xd) */ a2 = (Xaplusc - Xbplusd); /* yc' = (ya-yb+yc-yd) */ a3 = (Yaplusc - Ybplusd); /* xb' = (xa+yb-xc-yd) */ a4 = (Xaminusc + Ybminusd); /* yb' = (ya-xb-yc+xd) */ a5 = (Yaminusc - Xbminusd); /* xd' = (xa-yb-xc+yd)) */ a6 = (Xaminusc - Ybminusd); /* yd' = (ya+xb-yc-xd) */ a7 = (Xbminusd + Yaminusc); ptr1[0] = a0; ptr1[1] = a1; ptr1[2] = a2; ptr1[3] = a3; ptr1[4] = a4; ptr1[5] = a5; ptr1[6] = a6; ptr1[7] = a7; /* increment pointer by 8 */ ptr1 += 8u; } while(--j); #else float32_t t1, t2, r1, r2, s1, s2; /* Run the below code for Cortex-M0 */ /* Initializations for the fft calculation */ n2 = fftLen; n1 = n2; for (k = fftLen; k > 1u; k >>= 2u) { /* Initializations for the fft calculation */ n1 = n2; n2 >>= 2u; ia1 = 0u; /* FFT Calculation */ j = 0; do { /* index calculation for the coefficients */ ia2 = ia1 + ia1; ia3 = ia2 + ia1; co1 = pCoef[ia1 * 2u]; si1 = pCoef[(ia1 * 2u) + 1u]; co2 = pCoef[ia2 * 2u]; si2 = pCoef[(ia2 * 2u) + 1u]; co3 = pCoef[ia3 * 2u]; si3 = pCoef[(ia3 * 2u) + 1u]; /* Twiddle coefficients index modifier */ ia1 = ia1 + twidCoefModifier; i0 = j; do { /* index calculation for the input as, */ /* pSrc[i0 + 0], pSrc[i0 + fftLen/4], pSrc[i0 + fftLen/2], pSrc[i0 + 3fftLen/4] */ i1 = i0 + n2; i2 = i1 + n2; i3 = i2 + n2; /* xa + xc */ r1 = pSrc[(2u * i0)] + pSrc[(2u * i2)]; /* xa - xc */ r2 = pSrc[(2u * i0)] - pSrc[(2u * i2)]; /* ya + yc */ s1 = pSrc[(2u * i0) + 1u] + pSrc[(2u * i2) + 1u]; /* ya - yc */ s2 = pSrc[(2u * i0) + 1u] - pSrc[(2u * i2) + 1u]; /* xb + xd */ t1 = pSrc[2u * i1] + pSrc[2u * i3]; /* xa' = xa + xb + xc + xd */ pSrc[2u * i0] = r1 + t1; /* xa + xc -(xb + xd) */ r1 = r1 - t1; /* yb + yd */ t2 = pSrc[(2u * i1) + 1u] + pSrc[(2u * i3) + 1u]; /* ya' = ya + yb + yc + yd */ pSrc[(2u * i0) + 1u] = s1 + t2; /* (ya + yc) - (yb + yd) */ s1 = s1 - t2; /* (yb - yd) */ t1 = pSrc[(2u * i1) + 1u] - pSrc[(2u * i3) + 1u]; /* (xb - xd) */ t2 = pSrc[2u * i1] - pSrc[2u * i3]; /* xc' = (xa-xb+xc-xd)co2 + (ya-yb+yc-yd)(si2) */ pSrc[2u * i1] = (r1 * co2) + (s1 * si2); /* yc' = (ya-yb+yc-yd)co2 - (xa-xb+xc-xd)(si2) */ pSrc[(2u * i1) + 1u] = (s1 * co2) - (r1 * si2); /* (xa - xc) + (yb - yd) */ r1 = r2 + t1; /* (xa - xc) - (yb - yd) */ r2 = r2 - t1; /* (ya - yc) - (xb - xd) */ s1 = s2 - t2; /* (ya - yc) + (xb - xd) */ s2 = s2 + t2; /* xb' = (xa+yb-xc-yd)co1 + (ya-xb-yc+xd)(si1) */ pSrc[2u * i2] = (r1 * co1) + (s1 * si1); /* yb' = (ya-xb-yc+xd)co1 - (xa+yb-xc-yd)(si1) */ pSrc[(2u * i2) + 1u] = (s1 * co1) - (r1 * si1); /* xd' = (xa-yb-xc+yd)co3 + (ya+xb-yc-xd)(si3) */ pSrc[2u * i3] = (r2 * co3) + (s2 * si3); /* yd' = (ya+xb-yc-xd)co3 - (xa-yb-xc+yd)(si3) */ pSrc[(2u * i3) + 1u] = (s2 * co3) - (r2 * si3); i0 += n1; } while( i0 < fftLen); j++; } while(j <= (n2 - 1u)); twidCoefModifier <<= 2u; } #endif /* #ifndef ARM_MATH_CM0_FAMILY_FAMILY */ } /* * @brief Core function for the floating-point CIFFT butterfly process. * @param[in, out] *pSrc points to the in-place buffer of floating-point data type. * @param[in] fftLen length of the FFT. * @param[in] *pCoef points to twiddle coefficient buffer. * @param[in] twidCoefModifier twiddle coefficient modifier that supports different size FFTs with the same twiddle factor table. * @param[in] onebyfftLen value of 1/fftLen. * @return none. */ void arm_radix4_butterfly_inverse_f32( float32_t * pSrc, uint16_t fftLen, float32_t * pCoef, uint16_t twidCoefModifier, float32_t onebyfftLen) { float32_t co1, co2, co3, si1, si2, si3; uint32_t ia1, ia2, ia3; uint32_t i0, i1, i2, i3; uint32_t n1, n2, j, k; #ifndef ARM_MATH_CM0_FAMILY_FAMILY float32_t xaIn, yaIn, xbIn, ybIn, xcIn, ycIn, xdIn, ydIn; float32_t Xaplusc, Xbplusd, Yaplusc, Ybplusd, Xaminusc, Xbminusd, Yaminusc, Ybminusd; float32_t Xb12C_out, Yb12C_out, Xc12C_out, Yc12C_out, Xd12C_out, Yd12C_out; float32_t Xb12_out, Yb12_out, Xc12_out, Yc12_out, Xd12_out, Yd12_out; float32_t *ptr1; float32_t p0,p1,p2,p3,p4,p5,p6,p7; float32_t a0,a1,a2,a3,a4,a5,a6,a7; /* Initializations for the first stage */ n2 = fftLen; n1 = n2; /* n2 = fftLen/4 */ n2 >>= 2u; i0 = 0u; ia1 = 0u; j = n2; /* Calculation of first stage */ do { /* index calculation for the input as, */ /* pSrc[i0 + 0], pSrc[i0 + fftLen/4], pSrc[i0 + fftLen/2], pSrc[i0 + 3fftLen/4] */ i1 = i0 + n2; i2 = i1 + n2; i3 = i2 + n2; /* Butterfly implementation */ xaIn = pSrc[(2u * i0)]; yaIn = pSrc[(2u * i0) + 1u]; xcIn = pSrc[(2u * i2)]; ycIn = pSrc[(2u * i2) + 1u]; xbIn = pSrc[(2u * i1)]; ybIn = pSrc[(2u * i1) + 1u]; xdIn = pSrc[(2u * i3)]; ydIn = pSrc[(2u * i3) + 1u]; /* xa + xc */ Xaplusc = xaIn + xcIn; /* xb + xd */ Xbplusd = xbIn + xdIn; /* ya + yc */ Yaplusc = yaIn + ycIn; /* yb + yd */ Ybplusd = ybIn + ydIn; /* index calculation for the coefficients */ ia2 = ia1 + ia1; co2 = pCoef[ia2 * 2u]; si2 = pCoef[(ia2 * 2u) + 1u]; /* xa - xc */ Xaminusc = xaIn - xcIn; /* xb - xd */ Xbminusd = xbIn - xdIn; /* ya - yc */ Yaminusc = yaIn - ycIn; /* yb - yd */ Ybminusd = ybIn - ydIn; /* xa' = xa + xb + xc + xd */ pSrc[(2u * i0)] = Xaplusc + Xbplusd; /* ya' = ya + yb + yc + yd */ pSrc[(2u * i0) + 1u] = Yaplusc + Ybplusd; /* (xa - xc) - (yb - yd) */ Xb12C_out = (Xaminusc - Ybminusd); /* (ya - yc) + (xb - xd) */ Yb12C_out = (Yaminusc + Xbminusd); /* (xa + xc) - (xb + xd) */ Xc12C_out = (Xaplusc - Xbplusd); /* (ya + yc) - (yb + yd) */ Yc12C_out = (Yaplusc - Ybplusd); /* (xa - xc) + (yb - yd) */ Xd12C_out = (Xaminusc + Ybminusd); /* (ya - yc) - (xb - xd) */ Yd12C_out = (Yaminusc - Xbminusd); co1 = pCoef[ia1 * 2u]; si1 = pCoef[(ia1 * 2u) + 1u]; /* index calculation for the coefficients */ ia3 = ia2 + ia1; co3 = pCoef[ia3 * 2u]; si3 = pCoef[(ia3 * 2u) + 1u]; Xb12_out = Xb12C_out * co1; Yb12_out = Yb12C_out * co1; Xc12_out = Xc12C_out * co2; Yc12_out = Yc12C_out * co2; Xd12_out = Xd12C_out * co3; Yd12_out = Yd12C_out * co3; /* xb' = (xa+yb-xc-yd)co1 - (ya-xb-yc+xd)(si1) */ //Xb12_out -= Yb12C_out * si1; p0 = Yb12C_out * si1; /* yb' = (ya-xb-yc+xd)co1 + (xa+yb-xc-yd)(si1) */ //Yb12_out += Xb12C_out * si1; p1 = Xb12C_out * si1; /* xc' = (xa-xb+xc-xd)co2 - (ya-yb+yc-yd)(si2) */ //Xc12_out -= Yc12C_out * si2; p2 = Yc12C_out * si2; /* yc' = (ya-yb+yc-yd)co2 + (xa-xb+xc-xd)(si2) */ //Yc12_out += Xc12C_out * si2; p3 = Xc12C_out * si2; /* xd' = (xa-yb-xc+yd)co3 - (ya+xb-yc-xd)(si3) */ //Xd12_out -= Yd12C_out * si3; p4 = Yd12C_out * si3; /* yd' = (ya+xb-yc-xd)co3 + (xa-yb-xc+yd)(si3) */ //Yd12_out += Xd12C_out * si3; p5 = Xd12C_out * si3; Xb12_out -= p0; Yb12_out += p1; Xc12_out -= p2; Yc12_out += p3; Xd12_out -= p4; Yd12_out += p5; /* xc' = (xa-xb+xc-xd)co2 - (ya-yb+yc-yd)(si2) */ pSrc[2u * i1] = Xc12_out; /* yc' = (ya-yb+yc-yd)co2 + (xa-xb+xc-xd)(si2) */ pSrc[(2u * i1) + 1u] = Yc12_out; /* xb' = (xa+yb-xc-yd)co1 - (ya-xb-yc+xd)(si1) */ pSrc[2u * i2] = Xb12_out; /* yb' = (ya-xb-yc+xd)co1 + (xa+yb-xc-yd)(si1) */ pSrc[(2u * i2) + 1u] = Yb12_out; /* xd' = (xa-yb-xc+yd)co3 - (ya+xb-yc-xd)(si3) */ pSrc[2u * i3] = Xd12_out; /* yd' = (ya+xb-yc-xd)co3 + (xa-yb-xc+yd)(si3) */ pSrc[(2u * i3) + 1u] = Yd12_out; /* Twiddle coefficients index modifier */ ia1 = ia1 + twidCoefModifier; /* Updating input index */ i0 = i0 + 1u; } while(--j); twidCoefModifier <<= 2u; /* Calculation of second stage to excluding last stage */ for (k = fftLen >> 2u; k > 4u; k >>= 2u) { /* Initializations for the first stage */ n1 = n2; n2 >>= 2u; ia1 = 0u; /* Calculation of first stage */ j = 0; do { /* index calculation for the coefficients */ ia2 = ia1 + ia1; ia3 = ia2 + ia1; co1 = pCoef[ia1 * 2u]; si1 = pCoef[(ia1 * 2u) + 1u]; co2 = pCoef[ia2 * 2u]; si2 = pCoef[(ia2 * 2u) + 1u]; co3 = pCoef[ia3 * 2u]; si3 = pCoef[(ia3 * 2u) + 1u]; /* Twiddle coefficients index modifier */ ia1 = ia1 + twidCoefModifier; i0 = j; do { /* index calculation for the input as, */ /* pSrc[i0 + 0], pSrc[i0 + fftLen/4], pSrc[i0 + fftLen/2], pSrc[i0 + 3fftLen/4] */ i1 = i0 + n2; i2 = i1 + n2; i3 = i2 + n2; xaIn = pSrc[(2u * i0)]; yaIn = pSrc[(2u * i0) + 1u]; xbIn = pSrc[(2u * i1)]; ybIn = pSrc[(2u * i1) + 1u]; xcIn = pSrc[(2u * i2)]; ycIn = pSrc[(2u * i2) + 1u]; xdIn = pSrc[(2u * i3)]; ydIn = pSrc[(2u * i3) + 1u]; /* xa - xc */ Xaminusc = xaIn - xcIn; /* (xb - xd) */ Xbminusd = xbIn - xdIn; /* ya - yc */ Yaminusc = yaIn - ycIn; /* (yb - yd) */ Ybminusd = ybIn - ydIn; /* xa + xc */ Xaplusc = xaIn + xcIn; /* xb + xd */ Xbplusd = xbIn + xdIn; /* ya + yc */ Yaplusc = yaIn + ycIn; /* yb + yd */ Ybplusd = ybIn + ydIn; /* (xa - xc) - (yb - yd) */ Xb12C_out = (Xaminusc - Ybminusd); /* (ya - yc) + (xb - xd) */ Yb12C_out = (Yaminusc + Xbminusd); /* xa + xc -(xb + xd) */ Xc12C_out = (Xaplusc - Xbplusd); /* (ya + yc) - (yb + yd) */ Yc12C_out = (Yaplusc - Ybplusd); /* (xa - xc) + (yb - yd) */ Xd12C_out = (Xaminusc + Ybminusd); /* (ya - yc) - (xb - xd) */ Yd12C_out = (Yaminusc - Xbminusd); pSrc[(2u * i0)] = Xaplusc + Xbplusd; pSrc[(2u * i0) + 1u] = Yaplusc + Ybplusd; Xb12_out = Xb12C_out * co1; Yb12_out = Yb12C_out * co1; Xc12_out = Xc12C_out * co2; Yc12_out = Yc12C_out * co2; Xd12_out = Xd12C_out * co3; Yd12_out = Yd12C_out * co3; /* xb' = (xa+yb-xc-yd)co1 - (ya-xb-yc+xd)(si1) */ //Xb12_out -= Yb12C_out * si1; p0 = Yb12C_out * si1; /* yb' = (ya-xb-yc+xd)co1 + (xa+yb-xc-yd)(si1) */ //Yb12_out += Xb12C_out * si1; p1 = Xb12C_out * si1; /* xc' = (xa-xb+xc-xd)co2 - (ya-yb+yc-yd)(si2) */ //Xc12_out -= Yc12C_out * si2; p2 = Yc12C_out * si2; /* yc' = (ya-yb+yc-yd)co2 + (xa-xb+xc-xd)(si2) */ //Yc12_out += Xc12C_out * si2; p3 = Xc12C_out * si2; /* xd' = (xa-yb-xc+yd)co3 - (ya+xb-yc-xd)(si3) */ //Xd12_out -= Yd12C_out * si3; p4 = Yd12C_out * si3; /* yd' = (ya+xb-yc-xd)co3 + (xa-yb-xc+yd)(si3) */ //Yd12_out += Xd12C_out * si3; p5 = Xd12C_out * si3; Xb12_out -= p0; Yb12_out += p1; Xc12_out -= p2; Yc12_out += p3; Xd12_out -= p4; Yd12_out += p5; /* xc' = (xa-xb+xc-xd)co2 - (ya-yb+yc-yd)(si2) */ pSrc[2u * i1] = Xc12_out; /* yc' = (ya-yb+yc-yd)co2 + (xa-xb+xc-xd)(si2) */ pSrc[(2u * i1) + 1u] = Yc12_out; /* xb' = (xa+yb-xc-yd)co1 - (ya-xb-yc+xd)(si1) */ pSrc[2u * i2] = Xb12_out; /* yb' = (ya-xb-yc+xd)co1 + (xa+yb-xc-yd)(si1) */ pSrc[(2u * i2) + 1u] = Yb12_out; /* xd' = (xa-yb-xc+yd)co3 - (ya+xb-yc-xd)(si3) */ pSrc[2u * i3] = Xd12_out; /* yd' = (ya+xb-yc-xd)co3 + (xa-yb-xc+yd)(si3) */ pSrc[(2u * i3) + 1u] = Yd12_out; i0 += n1; } while(i0 < fftLen); j++; } while(j <= (n2 - 1u)); twidCoefModifier <<= 2u; } /* Initializations of last stage */ j = fftLen >> 2; ptr1 = &pSrc[0]; /* Calculations of last stage */ do { xaIn = ptr1[0]; yaIn = ptr1[1]; xbIn = ptr1[2]; ybIn = ptr1[3]; xcIn = ptr1[4]; ycIn = ptr1[5]; xdIn = ptr1[6]; ydIn = ptr1[7]; /* Butterfly implementation */ /* xa + xc */ Xaplusc = xaIn + xcIn; /* xa - xc */ Xaminusc = xaIn - xcIn; /* ya + yc */ Yaplusc = yaIn + ycIn; /* ya - yc */ Yaminusc = yaIn - ycIn; /* xb + xd */ Xbplusd = xbIn + xdIn; /* yb + yd */ Ybplusd = ybIn + ydIn; /* (xb-xd) */ Xbminusd = xbIn - xdIn; /* (yb-yd) */ Ybminusd = ybIn - ydIn; /* xa' = (xa+xb+xc+xd) * onebyfftLen */ a0 = (Xaplusc + Xbplusd); /* ya' = (ya+yb+yc+yd) * onebyfftLen */ a1 = (Yaplusc + Ybplusd); /* xc' = (xa-xb+xc-xd) * onebyfftLen */ a2 = (Xaplusc - Xbplusd); /* yc' = (ya-yb+yc-yd) * onebyfftLen */ a3 = (Yaplusc - Ybplusd); /* xb' = (xa-yb-xc+yd) * onebyfftLen */ a4 = (Xaminusc - Ybminusd); /* yb' = (ya+xb-yc-xd) * onebyfftLen */ a5 = (Yaminusc + Xbminusd); /* xd' = (xa-yb-xc+yd) * onebyfftLen */ a6 = (Xaminusc + Ybminusd); /* yd' = (ya-xb-yc+xd) * onebyfftLen */ a7 = (Yaminusc - Xbminusd); p0 = a0 * onebyfftLen; p1 = a1 * onebyfftLen; p2 = a2 * onebyfftLen; p3 = a3 * onebyfftLen; p4 = a4 * onebyfftLen; p5 = a5 * onebyfftLen; p6 = a6 * onebyfftLen; p7 = a7 * onebyfftLen; /* xa' = (xa+xb+xc+xd) * onebyfftLen */ ptr1[0] = p0; /* ya' = (ya+yb+yc+yd) * onebyfftLen */ ptr1[1] = p1; /* xc' = (xa-xb+xc-xd) * onebyfftLen */ ptr1[2] = p2; /* yc' = (ya-yb+yc-yd) * onebyfftLen */ ptr1[3] = p3; /* xb' = (xa-yb-xc+yd) * onebyfftLen */ ptr1[4] = p4; /* yb' = (ya+xb-yc-xd) * onebyfftLen */ ptr1[5] = p5; /* xd' = (xa-yb-xc+yd) * onebyfftLen */ ptr1[6] = p6; /* yd' = (ya-xb-yc+xd) * onebyfftLen */ ptr1[7] = p7; /* increment source pointer by 8 for next calculations */ ptr1 = ptr1 + 8u; } while(--j); #else float32_t t1, t2, r1, r2, s1, s2; /* Run the below code for Cortex-M0 */ /* Initializations for the first stage */ n2 = fftLen; n1 = n2; /* Calculation of first stage */ for (k = fftLen; k > 4u; k >>= 2u) { /* Initializations for the first stage */ n1 = n2; n2 >>= 2u; ia1 = 0u; /* Calculation of first stage */ j = 0; do { /* index calculation for the coefficients */ ia2 = ia1 + ia1; ia3 = ia2 + ia1; co1 = pCoef[ia1 * 2u]; si1 = pCoef[(ia1 * 2u) + 1u]; co2 = pCoef[ia2 * 2u]; si2 = pCoef[(ia2 * 2u) + 1u]; co3 = pCoef[ia3 * 2u]; si3 = pCoef[(ia3 * 2u) + 1u]; /* Twiddle coefficients index modifier */ ia1 = ia1 + twidCoefModifier; i0 = j; do { /* index calculation for the input as, */ /* pSrc[i0 + 0], pSrc[i0 + fftLen/4], pSrc[i0 + fftLen/2], pSrc[i0 + 3fftLen/4] */ i1 = i0 + n2; i2 = i1 + n2; i3 = i2 + n2; /* xa + xc */ r1 = pSrc[(2u * i0)] + pSrc[(2u * i2)]; /* xa - xc */ r2 = pSrc[(2u * i0)] - pSrc[(2u * i2)]; /* ya + yc */ s1 = pSrc[(2u * i0) + 1u] + pSrc[(2u * i2) + 1u]; /* ya - yc */ s2 = pSrc[(2u * i0) + 1u] - pSrc[(2u * i2) + 1u]; /* xb + xd */ t1 = pSrc[2u * i1] + pSrc[2u * i3]; /* xa' = xa + xb + xc + xd */ pSrc[2u * i0] = r1 + t1; /* xa + xc -(xb + xd) */ r1 = r1 - t1; /* yb + yd */ t2 = pSrc[(2u * i1) + 1u] + pSrc[(2u * i3) + 1u]; /* ya' = ya + yb + yc + yd */ pSrc[(2u * i0) + 1u] = s1 + t2; /* (ya + yc) - (yb + yd) */ s1 = s1 - t2; /* (yb - yd) */ t1 = pSrc[(2u * i1) + 1u] - pSrc[(2u * i3) + 1u]; /* (xb - xd) */ t2 = pSrc[2u * i1] - pSrc[2u * i3]; /* xc' = (xa-xb+xc-xd)co2 - (ya-yb+yc-yd)(si2) */ pSrc[2u * i1] = (r1 * co2) - (s1 * si2); /* yc' = (ya-yb+yc-yd)co2 + (xa-xb+xc-xd)(si2) */ pSrc[(2u * i1) + 1u] = (s1 * co2) + (r1 * si2); /* (xa - xc) - (yb - yd) */ r1 = r2 - t1; /* (xa - xc) + (yb - yd) */ r2 = r2 + t1; /* (ya - yc) + (xb - xd) */ s1 = s2 + t2; /* (ya - yc) - (xb - xd) */ s2 = s2 - t2; /* xb' = (xa+yb-xc-yd)co1 - (ya-xb-yc+xd)(si1) */ pSrc[2u * i2] = (r1 * co1) - (s1 * si1); /* yb' = (ya-xb-yc+xd)co1 + (xa+yb-xc-yd)(si1) */ pSrc[(2u * i2) + 1u] = (s1 * co1) + (r1 * si1); /* xd' = (xa-yb-xc+yd)co3 - (ya+xb-yc-xd)(si3) */ pSrc[2u * i3] = (r2 * co3) - (s2 * si3); /* yd' = (ya+xb-yc-xd)co3 + (xa-yb-xc+yd)(si3) */ pSrc[(2u * i3) + 1u] = (s2 * co3) + (r2 * si3); i0 += n1; } while( i0 < fftLen); j++; } while(j <= (n2 - 1u)); twidCoefModifier <<= 2u; } /* Initializations of last stage */ n1 = n2; n2 >>= 2u; /* Calculations of last stage */ for (i0 = 0u; i0 <= (fftLen - n1); i0 += n1) { /* index calculation for the input as, */ /* pSrc[i0 + 0], pSrc[i0 + fftLen/4], pSrc[i0 + fftLen/2], pSrc[i0 + 3fftLen/4] */ i1 = i0 + n2; i2 = i1 + n2; i3 = i2 + n2; /* Butterfly implementation */ /* xa + xc */ r1 = pSrc[2u * i0] + pSrc[2u * i2]; /* xa - xc */ r2 = pSrc[2u * i0] - pSrc[2u * i2]; /* ya + yc */ s1 = pSrc[(2u * i0) + 1u] + pSrc[(2u * i2) + 1u]; /* ya - yc */ s2 = pSrc[(2u * i0) + 1u] - pSrc[(2u * i2) + 1u]; /* xc + xd */ t1 = pSrc[2u * i1] + pSrc[2u * i3]; /* xa' = xa + xb + xc + xd */ pSrc[2u * i0] = (r1 + t1) * onebyfftLen; /* (xa + xb) - (xc + xd) */ r1 = r1 - t1; /* yb + yd */ t2 = pSrc[(2u * i1) + 1u] + pSrc[(2u * i3) + 1u]; /* ya' = ya + yb + yc + yd */ pSrc[(2u * i0) + 1u] = (s1 + t2) * onebyfftLen; /* (ya + yc) - (yb + yd) */ s1 = s1 - t2; /* (yb-yd) */ t1 = pSrc[(2u * i1) + 1u] - pSrc[(2u * i3) + 1u]; /* (xb-xd) */ t2 = pSrc[2u * i1] - pSrc[2u * i3]; /* xc' = (xa-xb+xc-xd)co2 - (ya-yb+yc-yd)(si2) */ pSrc[2u * i1] = r1 * onebyfftLen; /* yc' = (ya-yb+yc-yd)co2 + (xa-xb+xc-xd)(si2) */ pSrc[(2u * i1) + 1u] = s1 * onebyfftLen; /* (xa - xc) - (yb-yd) */ r1 = r2 - t1; /* (xa - xc) + (yb-yd) */ r2 = r2 + t1; /* (ya - yc) + (xb-xd) */ s1 = s2 + t2; /* (ya - yc) - (xb-xd) */ s2 = s2 - t2; /* xb' = (xa+yb-xc-yd)co1 - (ya-xb-yc+xd)(si1) */ pSrc[2u * i2] = r1 * onebyfftLen; /* yb' = (ya-xb-yc+xd)co1 + (xa+yb-xc-yd)(si1) */ pSrc[(2u * i2) + 1u] = s1 * onebyfftLen; /* xd' = (xa-yb-xc+yd)co3 - (ya+xb-yc-xd)(si3) */ pSrc[2u * i3] = r2 * onebyfftLen; /* yd' = (ya+xb-yc-xd)co3 + (xa-yb-xc+yd)(si3) */ pSrc[(2u * i3) + 1u] = s2 * onebyfftLen; } #endif /* #ifndef ARM_MATH_CM0_FAMILY_FAMILY */ } /** * @addtogroup ComplexFFT * @{ */ /** * @details * @brief Processing function for the floating-point Radix-4 CFFT/CIFFT. * @deprecated Do not use this function. It has been superceded by \ref arm_cfft_f32 and will be removed * in the future. * @param[in] *S points to an instance of the floating-point Radix-4 CFFT/CIFFT structure. * @param[in, out] *pSrc points to the complex data buffer of size 2*fftLen. Processing occurs in-place. * @return none. */ void arm_cfft_radix4_f32( const arm_cfft_radix4_instance_f32 * S, float32_t * pSrc) { if(S->ifftFlag == 1u) { /* Complex IFFT radix-4 */ arm_radix4_butterfly_inverse_f32(pSrc, S->fftLen, S->pTwiddle, S->twidCoefModifier, S->onebyfftLen); } else { /* Complex FFT radix-4 */ arm_radix4_butterfly_f32(pSrc, S->fftLen, S->pTwiddle, S->twidCoefModifier); } if(S->bitReverseFlag == 1u) { /* Bit Reversal */ arm_bitreversal_f32(pSrc, S->fftLen, S->bitRevFactor, S->pBitRevTable); } } /** * @} end of ComplexFFT group */