Line data Source code
1 : /******************************************************************************
2 : * ETSI TS 103 634 V1.5.1 *
3 : * Low Complexity Communication Codec Plus (LC3plus) *
4 : * *
5 : * Copyright licence is solely granted through ETSI Intellectual Property *
6 : * Rights Policy, 3rd April 2019. No patent licence is granted by implication, *
7 : * estoppel or otherwise. *
8 : ******************************************************************************/
9 :
10 : /* guard against unindended includes */
11 : #ifndef INCLUDED_FROM_IISFFT_C
12 : #error "this file must not be included"
13 : #endif
14 :
15 0 : static __inline void butterfly(const LC3_FLOAT a, const LC3_FLOAT b, LC3_FLOAT* aPlusb, LC3_FLOAT* aMinusb)
16 : {
17 0 : *aPlusb = a + b;
18 0 : *aMinusb = a - b;
19 0 : }
20 :
21 0 : static void fft2(LC3_FLOAT* vec)
22 : {
23 : /*
24 : 1.0000 1.0000
25 : 1.0000 -1.0000
26 : */
27 0 : LC3_FLOAT re1 = vec[0];
28 0 : LC3_FLOAT im1 = vec[1];
29 0 : LC3_FLOAT re2 = vec[2];
30 0 : LC3_FLOAT im2 = vec[3];
31 :
32 0 : vec[0] = re1 + re2;
33 0 : vec[1] = im1 + im2;
34 0 : vec[2] = re1 - re2;
35 0 : vec[3] = im1 - im2;
36 0 : }
37 :
38 0 : static void fft3(LC3_FLOAT* vec)
39 : {
40 0 : const LC3_FLOAT C31 = 0.5; /* cos(PI/3); sin(2*PI/3) */
41 0 : const LC3_FLOAT C32 = 0.866025403784439; /* cos(PI/3); sin(2*PI/3) */
42 :
43 0 : LC3_FLOAT re1 = vec[0];
44 0 : LC3_FLOAT im1 = vec[1];
45 0 : LC3_FLOAT re2 = vec[2];
46 0 : LC3_FLOAT im2 = vec[3];
47 0 : LC3_FLOAT re3 = vec[4];
48 0 : LC3_FLOAT im3 = vec[5];
49 : /*
50 : 1.0000 1.0000 1.0000
51 : C31 C32
52 : 1.0000 -0.5000 - 0.8660i -0.5000 + 0.8660i
53 : 1.0000 -0.5000 + 0.8660i -0.5000 - 0.8660i
54 : */
55 0 : LC3_FLOAT tmp1 = re2 + re3;
56 0 : LC3_FLOAT tmp3 = im2 + im3;
57 0 : LC3_FLOAT tmp2 = re2 - re3;
58 0 : LC3_FLOAT tmp4 = im2 - im3;
59 :
60 0 : vec[0] = re1 + tmp1;
61 0 : vec[1] = im1 + tmp3;
62 0 : vec[2] = re1 - C31 * tmp1 + C32 * tmp4;
63 0 : vec[4] = re1 - C31 * tmp1 - C32 * tmp4;
64 0 : vec[3] = im1 - C32 * tmp2 - C31 * tmp3;
65 0 : vec[5] = im1 + C32 * tmp2 - C31 * tmp3;
66 0 : }
67 :
68 0 : static void fft4(LC3_FLOAT* vec)
69 : {
70 : LC3_FLOAT temp0, temp1, temp2, temp3, temp4, temp5, temp6, temp7;
71 :
72 : /* Pre-additions */
73 0 : temp0 = vec[0] + vec[4];
74 0 : temp2 = vec[0] - vec[4];
75 0 : temp1 = vec[1] + vec[5];
76 0 : temp3 = vec[1] - vec[5];
77 0 : temp4 = vec[2] + vec[6];
78 0 : temp7 = vec[2] - vec[6];
79 0 : temp5 = vec[7] + vec[3];
80 0 : temp6 = vec[7] - vec[3];
81 :
82 : /* Post-additions */
83 0 : vec[0] = temp0 + temp4;
84 0 : vec[1] = temp1 + temp5;
85 0 : vec[2] = temp2 - temp6;
86 0 : vec[3] = temp3 - temp7;
87 0 : vec[4] = temp0 - temp4;
88 0 : vec[5] = temp1 - temp5;
89 0 : vec[6] = temp2 + temp6;
90 0 : vec[7] = temp3 + temp7;
91 0 : }
92 :
93 0 : static void fft5(LC3_FLOAT* vec)
94 : {
95 0 : const LC3_FLOAT C51 = 0.309016994374947; /* cos(2*PI/5); */
96 0 : const LC3_FLOAT C52 = 0.951056516295154; /* sin(2*PI/5); */
97 0 : const LC3_FLOAT C53 = 0.809016994374947; /* cos( PI/5); */
98 0 : const LC3_FLOAT C54 = 0.587785252292473; /* sin( PI/5); */
99 :
100 : LC3_FLOAT re1, im1, re2, im2, re3, im3, re4, im4, re5, im5, tmp1, tmp2, tmp3, tmp4, tmp5, tmp6, tmp7, tmp8;
101 :
102 0 : re1 = vec[0];
103 0 : im1 = vec[1];
104 0 : re2 = vec[2];
105 0 : im2 = vec[3];
106 0 : re3 = vec[4];
107 0 : im3 = vec[5];
108 0 : re4 = vec[6];
109 0 : im4 = vec[7];
110 0 : re5 = vec[8];
111 0 : im5 = vec[9];
112 : /*
113 : 1.0000 1.0000 1.0000 1.0000 1.0000
114 : C51 C52 C53 C54
115 : 1.0000 0.3090 - 0.9511i -0.8090 - 0.5878i -0.8090 + 0.5878i 0.3090 + 0.9511i
116 : 1.0000 -0.8090 - 0.5878i 0.3090 + 0.9511i 0.3090 - 0.9511i -0.8090 + 0.5878i
117 : 1.0000 -0.8090 + 0.5878i 0.3090 - 0.9511i 0.3090 + 0.9511i -0.8090 - 0.5878i
118 : 1.0000 0.3090 + 0.9511i -0.8090 + 0.5878i -0.8090 - 0.5878i 0.3090 - 0.9511i
119 : */
120 0 : tmp1 = re2 + re5;
121 0 : tmp2 = re2 - re5;
122 0 : tmp3 = im2 + im5;
123 0 : tmp4 = im2 - im5;
124 0 : tmp5 = re3 + re4;
125 0 : tmp6 = re3 - re4;
126 0 : tmp7 = im3 + im4;
127 0 : tmp8 = im3 - im4;
128 :
129 0 : vec[0] = re1 + tmp1 + tmp5;
130 0 : vec[1] = im1 + tmp3 + tmp7;
131 0 : vec[2] = re1 + C51 * tmp1 - C53 * tmp5 + C52 * tmp4 + C54 * tmp8;
132 0 : vec[8] = re1 + C51 * tmp1 - C53 * tmp5 - C52 * tmp4 - C54 * tmp8;
133 0 : vec[3] = im1 - C52 * tmp2 - C54 * tmp6 + C51 * tmp3 - C53 * tmp7;
134 0 : vec[9] = im1 + C52 * tmp2 + C54 * tmp6 + C51 * tmp3 - C53 * tmp7;
135 0 : vec[4] = re1 - C53 * tmp1 + C51 * tmp5 + C54 * tmp4 - C52 * tmp8;
136 0 : vec[6] = re1 - C53 * tmp1 + C51 * tmp5 - C54 * tmp4 + C52 * tmp8;
137 0 : vec[5] = im1 - C54 * tmp2 + C52 * tmp6 - C53 * tmp3 + C51 * tmp7;
138 0 : vec[7] = im1 + C54 * tmp2 - C52 * tmp6 - C53 * tmp3 + C51 * tmp7;
139 0 : }
140 :
141 :
142 :
143 0 : static void fft8(LC3_FLOAT* vec)
144 : {
145 0 : const LC3_FLOAT INV_SQRT2 = 7.071067811865475e-1;
146 : LC3_FLOAT temp1[16], temp2[16];
147 :
148 : /* Pre-additions */
149 0 : temp1[0] = vec[0] + vec[8];
150 0 : temp1[2] = vec[0] - vec[8];
151 0 : temp1[1] = vec[1] + vec[9];
152 0 : temp1[3] = vec[1] - vec[9];
153 0 : temp1[4] = vec[2] + vec[10];
154 0 : temp1[6] = vec[2] - vec[10];
155 0 : temp1[5] = vec[3] + vec[11];
156 0 : temp1[7] = vec[3] - vec[11];
157 0 : temp1[8] = vec[4] + vec[12];
158 0 : temp1[10] = vec[4] - vec[12];
159 0 : temp1[9] = vec[5] + vec[13];
160 0 : temp1[11] = vec[5] - vec[13];
161 0 : temp1[12] = vec[6] + vec[14];
162 0 : temp1[14] = vec[6] - vec[14];
163 0 : temp1[13] = vec[7] + vec[15];
164 0 : temp1[15] = vec[7] - vec[15];
165 :
166 : /* Pre-additions and core multiplications */
167 0 : temp2[0] = temp1[0] + temp1[8];
168 0 : temp2[4] = temp1[0] - temp1[8];
169 0 : temp2[1] = temp1[1] + temp1[9];
170 0 : temp2[5] = temp1[1] - temp1[9];
171 0 : temp2[8] = temp1[2] - temp1[11];
172 0 : temp2[10] = temp1[2] + temp1[11];
173 0 : temp2[9] = temp1[3] + temp1[10];
174 0 : temp2[11] = temp1[3] - temp1[10];
175 0 : temp2[2] = temp1[4] + temp1[12];
176 0 : temp2[7] = temp1[4] - temp1[12];
177 0 : temp2[3] = temp1[5] + temp1[13];
178 0 : temp2[6] = temp1[13] - temp1[5];
179 :
180 0 : temp1[1] = temp1[6] + temp1[14];
181 0 : temp1[2] = temp1[6] - temp1[14];
182 0 : temp1[0] = temp1[7] + temp1[15];
183 0 : temp1[3] = temp1[7] - temp1[15];
184 0 : temp2[12] = (temp1[0] + temp1[2]) * INV_SQRT2;
185 0 : temp2[14] = (temp1[0] - temp1[2]) * INV_SQRT2;
186 0 : temp2[13] = (temp1[3] - temp1[1]) * INV_SQRT2;
187 0 : temp2[15] = (temp1[1] + temp1[3]) * -INV_SQRT2;
188 :
189 : /* Post-additions */
190 0 : vec[0] = temp2[0] + temp2[2];
191 0 : vec[8] = temp2[0] - temp2[2];
192 0 : vec[1] = temp2[1] + temp2[3];
193 0 : vec[9] = temp2[1] - temp2[3];
194 0 : vec[4] = temp2[4] - temp2[6];
195 0 : vec[12] = temp2[4] + temp2[6];
196 0 : vec[5] = temp2[5] - temp2[7];
197 0 : vec[13] = temp2[5] + temp2[7];
198 0 : vec[6] = temp2[8] + temp2[14];
199 0 : vec[14] = temp2[8] - temp2[14];
200 0 : vec[7] = temp2[9] + temp2[15];
201 0 : vec[15] = temp2[9] - temp2[15];
202 0 : vec[2] = temp2[10] + temp2[12];
203 0 : vec[10] = temp2[10] - temp2[12];
204 0 : vec[3] = temp2[11] + temp2[13];
205 0 : vec[11] = temp2[11] - temp2[13];
206 0 : }
207 :
208 :
209 0 : static void fft9(LC3_FLOAT* vec)
210 : {
211 0 : const LC3_FLOAT C91 = 0.766044443118978; /* cos(2*PI/5); */
212 0 : const LC3_FLOAT C92 = 0.642787609686539;
213 0 : const LC3_FLOAT C93 = 0.17364817766693;
214 0 : const LC3_FLOAT C94 = 0.984807753012208;
215 0 : const LC3_FLOAT C95 = 0.5;
216 0 : const LC3_FLOAT C96 = 0.866025403784439;
217 0 : const LC3_FLOAT C97 = 0.939692620785908;
218 0 : const LC3_FLOAT C98 = 0.342020143325669;
219 :
220 : LC3_FLOAT re1, im1, re2_9p, re2_9m, im2_9p, im2_9m, re3_8p, re3_8m, im3_8p, im3_8m, re4_7p, re4_7m, im4_7p, im4_7m,
221 : re5_6p, re5_6m, im5_6p, im5_6m;
222 :
223 0 : re1 = vec[0];
224 0 : im1 = vec[1];
225 :
226 0 : butterfly(vec[1 * 2], vec[8 * 2], &re2_9p, &re2_9m);
227 0 : butterfly(vec[1 * 2 + 1], vec[8 * 2 + 1], &im2_9p, &im2_9m);
228 0 : butterfly(vec[2 * 2], vec[7 * 2], &re3_8p, &re3_8m);
229 0 : butterfly(vec[2 * 2 + 1], vec[7 * 2 + 1], &im3_8p, &im3_8m);
230 0 : butterfly(vec[3 * 2], vec[6 * 2], &re4_7p, &re4_7m);
231 0 : butterfly(vec[3 * 2 + 1], vec[6 * 2 + 1], &im4_7p, &im4_7m);
232 0 : butterfly(vec[4 * 2], vec[5 * 2], &re5_6p, &re5_6m);
233 0 : butterfly(vec[4 * 2 + 1], vec[5 * 2 + 1], &im5_6p, &im5_6m);
234 :
235 : /*
236 : 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000
237 : 1.0000 C91 - C92i C93 - C94i -C95 - C96i -C97 - C98i -C97 + C98i -C95 + C96i C93 + C94i C91 + C92i
238 : 1.0000 C93 - C94i -C97 - C98i -C95 + C96i C91 + C92i C91 - C92i -C95 - C96i -C97 + C98i C93 + C94i
239 : 1.0000 -C95 - C96i -C95 + C96i 1.0000 -C95 - C96i -C95 + C96i 1.0000 -C95 - C96i -C95 + C96i
240 : 1.0000 -C97 - C98i C91 + C92i -C95 - C96i C93 + C94i C93 - C94i -C95 + C96i C91 - C92i -C97 + C98i
241 : 1.0000 -C97 + C98i C91 - C92i -C95 + C96i C93 - C94i C93 + C94i -C95 - C96i C91 + C92i -C97 - C98i
242 : 1.0000 -C95 + C96i -C95 - C96i 1.0000 -C95 + C96i -C95 - C96i 1.0000 -C95 + C96i -C95 - C96i
243 : 1.0000 C93 + C94i -C97 + C98i -C95 - C96i C91 - C92i C91 + C92i -C95 + C96i -C97 - C98i C93 - C94i
244 : 1.0000 C91 + C92i C93 + C94i -C95 + C96i -C97 + C98i -C97 - C98i -C95 - C96i C93 - C94i C91 - C92i
245 : */
246 0 : vec[0] = re1 + re2_9p + re3_8p + re4_7p + re5_6p;
247 0 : vec[1] = im1 + im2_9p + im3_8p + im4_7p + im5_6p;
248 0 : vec[2] = re1 + C91 * re2_9p + C93 * re3_8p - C95 * re4_7p - C97 * re5_6p + C92 * im2_9m + C94 * im3_8m +
249 0 : C96 * im4_7m + C98 * im5_6m;
250 0 : vec[16] = re1 + C91 * re2_9p + C93 * re3_8p - C95 * re4_7p - C97 * re5_6p - C92 * im2_9m - C94 * im3_8m -
251 0 : C96 * im4_7m - C98 * im5_6m;
252 0 : vec[3] = im1 - C92 * re2_9m - C94 * re3_8m - C96 * re4_7m - C98 * re5_6m + C91 * im2_9p + C93 * im3_8p -
253 0 : C95 * im4_7p - C97 * im5_6p;
254 0 : vec[17] = im1 + C92 * re2_9m + C94 * re3_8m + C96 * re4_7m + C98 * re5_6m + C91 * im2_9p + C93 * im3_8p -
255 0 : C95 * im4_7p - C97 * im5_6p;
256 0 : vec[4] = re1 + C93 * re2_9p - C97 * re3_8p - C95 * re4_7p + C91 * re5_6p + C94 * im2_9m + C98 * im3_8m -
257 0 : C96 * im4_7m - C92 * im5_6m;
258 0 : vec[14] = re1 + C93 * re2_9p - C97 * re3_8p - C95 * re4_7p + C91 * re5_6p - C94 * im2_9m - C98 * im3_8m +
259 0 : C96 * im4_7m + C92 * im5_6m;
260 0 : vec[5] = im1 - C94 * re2_9m - C98 * re3_8m + C96 * re4_7m + C92 * re5_6m + C93 * im2_9p - C97 * im3_8p -
261 0 : C95 * im4_7p + C91 * im5_6p;
262 0 : vec[15] = im1 + C94 * re2_9m + C98 * re3_8m - C96 * re4_7m - C92 * re5_6m + C93 * im2_9p - C97 * im3_8p -
263 0 : C95 * im4_7p + C91 * im5_6p;
264 0 : vec[6] = re1 - C95 * (re2_9p + re3_8p + re5_6p) + re4_7p + C96 * (im2_9m - im3_8m + im5_6m);
265 0 : vec[12] = re1 - C95 * (re2_9p + re3_8p + re5_6p) + re4_7p - C96 * (im2_9m - im3_8m + im5_6m);
266 0 : vec[7] = im1 - C96 * (re2_9m - re3_8m + re5_6m) - C95 * (im2_9p + im3_8p + im5_6p) + im4_7p;
267 0 : vec[13] = im1 + C96 * (re2_9m - re3_8m + re5_6m) - C95 * (im2_9p + im3_8p + im5_6p) + im4_7p;
268 0 : vec[8] = re1 - C97 * re2_9p + C91 * re3_8p - C95 * re4_7p + C93 * re5_6p + C98 * im2_9m - C92 * im3_8m +
269 0 : C96 * im4_7m - C94 * im5_6m;
270 0 : vec[10] = re1 - C97 * re2_9p + C91 * re3_8p - C95 * re4_7p + C93 * re5_6p - C98 * im2_9m + C92 * im3_8m -
271 0 : C96 * im4_7m + C94 * im5_6m;
272 0 : vec[9] = im1 - C98 * re2_9m + C92 * re3_8m - C96 * re4_7m + C94 * re5_6m - C97 * im2_9p + C91 * im3_8p -
273 0 : C95 * im4_7p + C93 * im5_6p;
274 0 : vec[11] = im1 + C98 * re2_9m - C92 * re3_8m + C96 * re4_7m - C94 * re5_6m - C97 * im2_9p + C91 * im3_8p -
275 0 : C95 * im4_7p + C93 * im5_6p;
276 0 : }
277 :
|