/proc/self/cwd/libfaad/cfft.c
Line | Count | Source (jump to first uncovered line) |
1 | | /* |
2 | | ** FAAD2 - Freeware Advanced Audio (AAC) Decoder including SBR decoding |
3 | | ** Copyright (C) 2003-2005 M. Bakker, Nero AG, http://www.nero.com |
4 | | ** |
5 | | ** This program is free software; you can redistribute it and/or modify |
6 | | ** it under the terms of the GNU General Public License as published by |
7 | | ** the Free Software Foundation; either version 2 of the License, or |
8 | | ** (at your option) any later version. |
9 | | ** |
10 | | ** This program is distributed in the hope that it will be useful, |
11 | | ** but WITHOUT ANY WARRANTY; without even the implied warranty of |
12 | | ** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
13 | | ** GNU General Public License for more details. |
14 | | ** |
15 | | ** You should have received a copy of the GNU General Public License |
16 | | ** along with this program; if not, write to the Free Software |
17 | | ** Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. |
18 | | ** |
19 | | ** Any non-GPL usage of this software or parts of this software is strictly |
20 | | ** forbidden. |
21 | | ** |
22 | | ** The "appropriate copyright message" mentioned in section 2c of the GPLv2 |
23 | | ** must read: "Code from FAAD2 is copyright (c) Nero AG, www.nero.com" |
24 | | ** |
25 | | ** Commercial non-GPL licensing of this software is possible. |
26 | | ** For more info contact Nero AG through Mpeg4AAClicense@nero.com. |
27 | | ** |
28 | | ** $Id: cfft.c,v 1.35 2007/11/01 12:33:29 menno Exp $ |
29 | | **/ |
30 | | |
31 | | /* |
32 | | * Algorithmically based on Fortran-77 FFTPACK |
33 | | * by Paul N. Swarztrauber(Version 4, 1985). |
34 | | * |
35 | | * Does even sized fft only |
36 | | */ |
37 | | |
38 | | /* isign is +1 for backward and -1 for forward transforms */ |
39 | | |
40 | | #include "common.h" |
41 | | #include "structs.h" |
42 | | |
43 | | #include <stdlib.h> |
44 | | |
45 | | #include "cfft.h" |
46 | | #include "cfft_tab.h" |
47 | | |
48 | | |
49 | | /* static function declarations */ |
50 | | static void passf2pos(const uint16_t ido, const uint16_t l1, const complex_t *cc, |
51 | | complex_t *ch, const complex_t *wa); |
52 | | static void passf2neg(const uint16_t ido, const uint16_t l1, const complex_t *cc, |
53 | | complex_t *ch, const complex_t *wa); |
54 | | static void passf3(const uint16_t ido, const uint16_t l1, const complex_t *cc, |
55 | | complex_t *ch, const complex_t *wa1, const complex_t *wa2, const int8_t isign); |
56 | | static void passf4pos(const uint16_t ido, const uint16_t l1, const complex_t *cc, complex_t *ch, |
57 | | const complex_t *wa1, const complex_t *wa2, const complex_t *wa3); |
58 | | static void passf4neg(const uint16_t ido, const uint16_t l1, const complex_t *cc, complex_t *ch, |
59 | | const complex_t *wa1, const complex_t *wa2, const complex_t *wa3); |
60 | | static void passf5(const uint16_t ido, const uint16_t l1, const complex_t *cc, complex_t *ch, |
61 | | const complex_t *wa1, const complex_t *wa2, const complex_t *wa3, |
62 | | const complex_t *wa4, const int8_t isign); |
63 | | static void cffti1(uint16_t n, complex_t *wa, uint16_t *ifac); |
64 | | |
65 | | |
66 | | /*---------------------------------------------------------------------- |
67 | | passf2, passf3, passf4, passf5. Complex FFT passes fwd and bwd. |
68 | | ----------------------------------------------------------------------*/ |
69 | | |
70 | | static void passf2pos(const uint16_t ido, const uint16_t l1, const complex_t *cc, |
71 | | complex_t *ch, const complex_t *wa) |
72 | 559k | { |
73 | 559k | uint16_t i, k, ah, ac; |
74 | | |
75 | 559k | if (ido == 1) |
76 | 0 | { |
77 | 0 | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION |
78 | | // TODO: remove this code once fuzzer proves it is totally unreahable |
79 | | // For supported frame lengths that have odd number of factor 2 it is |
80 | | // never the last factor; consequently `ido` should never be 1. |
81 | 0 | __builtin_trap(); |
82 | | /* |
83 | | #endif |
84 | | for (k = 0; k < l1; k++) |
85 | | { |
86 | | ah = 2*k; |
87 | | ac = 4*k; |
88 | | |
89 | | RE(ch[ah]) = RE(cc[ac]) + RE(cc[ac+1]); |
90 | | RE(ch[ah+l1]) = RE(cc[ac]) - RE(cc[ac+1]); |
91 | | IM(ch[ah]) = IM(cc[ac]) + IM(cc[ac+1]); |
92 | | IM(ch[ah+l1]) = IM(cc[ac]) - IM(cc[ac+1]); |
93 | | } |
94 | | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION |
95 | | */ |
96 | 0 | #endif |
97 | 559k | } else { |
98 | 1.11M | for (k = 0; k < l1; k++) |
99 | 559k | { |
100 | 559k | ah = k*ido; |
101 | 559k | ac = 2*k*ido; |
102 | | |
103 | 142M | for (i = 0; i < ido; i++) |
104 | 141M | { |
105 | 141M | complex_t t2; |
106 | | |
107 | 141M | RE(ch[ah+i]) = RE(cc[ac+i]) + RE(cc[ac+i+ido]); |
108 | 141M | RE(t2) = RE(cc[ac+i]) - RE(cc[ac+i+ido]); |
109 | | |
110 | 141M | IM(ch[ah+i]) = IM(cc[ac+i]) + IM(cc[ac+i+ido]); |
111 | 141M | IM(t2) = IM(cc[ac+i]) - IM(cc[ac+i+ido]); |
112 | | |
113 | 141M | #if 1 |
114 | 141M | ComplexMult(&IM(ch[ah+i+l1*ido]), &RE(ch[ah+i+l1*ido]), |
115 | 141M | IM(t2), RE(t2), RE(wa[i]), IM(wa[i])); |
116 | | #else |
117 | | ComplexMult(&RE(ch[ah+i+l1*ido]), &IM(ch[ah+i+l1*ido]), |
118 | | RE(t2), IM(t2), RE(wa[i]), IM(wa[i])); |
119 | | #endif |
120 | 141M | } |
121 | 559k | } |
122 | 559k | } |
123 | 559k | } |
124 | | |
125 | | static void passf2neg(const uint16_t ido, const uint16_t l1, const complex_t *cc, |
126 | | complex_t *ch, const complex_t *wa) |
127 | 12.2k | { |
128 | 12.2k | uint16_t i, k, ah, ac; |
129 | | |
130 | 12.2k | if (ido == 1) |
131 | 0 | { |
132 | 0 | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION |
133 | | // TODO: remove this code once fuzzer proves it is totally unreahable |
134 | | // For supported frame lengths that have odd number of factor 2 it is |
135 | | // never the last factor; consequently `ido` should never be 1. |
136 | 0 | __builtin_trap(); |
137 | | /* |
138 | | #endif |
139 | | for (k = 0; k < l1; k++) |
140 | | { |
141 | | ah = 2*k; |
142 | | ac = 4*k; |
143 | | |
144 | | RE(ch[ah]) = RE(cc[ac]) + RE(cc[ac+1]); |
145 | | RE(ch[ah+l1]) = RE(cc[ac]) - RE(cc[ac+1]); |
146 | | IM(ch[ah]) = IM(cc[ac]) + IM(cc[ac+1]); |
147 | | IM(ch[ah+l1]) = IM(cc[ac]) - IM(cc[ac+1]); |
148 | | } |
149 | | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION |
150 | | */ |
151 | 0 | #endif |
152 | 12.2k | } else { |
153 | 24.4k | for (k = 0; k < l1; k++) |
154 | 12.2k | { |
155 | 12.2k | ah = k*ido; |
156 | 12.2k | ac = 2*k*ido; |
157 | | |
158 | 3.02M | for (i = 0; i < ido; i++) |
159 | 3.01M | { |
160 | 3.01M | complex_t t2; |
161 | | |
162 | 3.01M | RE(ch[ah+i]) = RE(cc[ac+i]) + RE(cc[ac+i+ido]); |
163 | 3.01M | RE(t2) = RE(cc[ac+i]) - RE(cc[ac+i+ido]); |
164 | | |
165 | 3.01M | IM(ch[ah+i]) = IM(cc[ac+i]) + IM(cc[ac+i+ido]); |
166 | 3.01M | IM(t2) = IM(cc[ac+i]) - IM(cc[ac+i+ido]); |
167 | | |
168 | 3.01M | #if 1 |
169 | 3.01M | ComplexMult(&RE(ch[ah+i+l1*ido]), &IM(ch[ah+i+l1*ido]), |
170 | 3.01M | RE(t2), IM(t2), RE(wa[i]), IM(wa[i])); |
171 | | #else |
172 | | ComplexMult(&IM(ch[ah+i+l1*ido]), &RE(ch[ah+i+l1*ido]), |
173 | | IM(t2), RE(t2), RE(wa[i]), IM(wa[i])); |
174 | | #endif |
175 | 3.01M | } |
176 | 12.2k | } |
177 | 12.2k | } |
178 | 12.2k | } |
179 | | |
180 | | |
181 | | static void passf3(const uint16_t ido, const uint16_t l1, const complex_t *cc, |
182 | | complex_t *ch, const complex_t *wa1, const complex_t *wa2, |
183 | | const int8_t isign) |
184 | 186k | { |
185 | 186k | static real_t taur = FRAC_CONST(-0.5); |
186 | 186k | static real_t taui = FRAC_CONST(0.866025403784439); |
187 | 186k | uint16_t i, k, ac, ah; |
188 | 186k | complex_t c2, c3, d2, d3, t2; |
189 | | |
190 | 186k | if (ido == 1) |
191 | 0 | { |
192 | 0 | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION |
193 | | // TODO: remove this code once fuzzer proves it is totally unreahable |
194 | | // 3 is never the the biggest factor for supported frame lengths; |
195 | | // consequently `ido` should never be 1. |
196 | 0 | __builtin_trap(); |
197 | | /* |
198 | | #endif |
199 | | if (isign == 1) |
200 | | { |
201 | | for (k = 0; k < l1; k++) |
202 | | { |
203 | | ac = 3*k+1; |
204 | | ah = k; |
205 | | |
206 | | RE(t2) = RE(cc[ac]) + RE(cc[ac+1]); |
207 | | IM(t2) = IM(cc[ac]) + IM(cc[ac+1]); |
208 | | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),taur); |
209 | | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),taur); |
210 | | |
211 | | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2); |
212 | | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2); |
213 | | |
214 | | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+1])), taui); |
215 | | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+1])), taui); |
216 | | |
217 | | RE(ch[ah+l1]) = RE(c2) - IM(c3); |
218 | | IM(ch[ah+l1]) = IM(c2) + RE(c3); |
219 | | RE(ch[ah+2*l1]) = RE(c2) + IM(c3); |
220 | | IM(ch[ah+2*l1]) = IM(c2) - RE(c3); |
221 | | } |
222 | | } else { |
223 | | for (k = 0; k < l1; k++) |
224 | | { |
225 | | ac = 3*k+1; |
226 | | ah = k; |
227 | | |
228 | | RE(t2) = RE(cc[ac]) + RE(cc[ac+1]); |
229 | | IM(t2) = IM(cc[ac]) + IM(cc[ac+1]); |
230 | | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),taur); |
231 | | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),taur); |
232 | | |
233 | | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2); |
234 | | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2); |
235 | | |
236 | | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+1])), taui); |
237 | | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+1])), taui); |
238 | | |
239 | | RE(ch[ah+l1]) = RE(c2) + IM(c3); |
240 | | IM(ch[ah+l1]) = IM(c2) - RE(c3); |
241 | | RE(ch[ah+2*l1]) = RE(c2) - IM(c3); |
242 | | IM(ch[ah+2*l1]) = IM(c2) + RE(c3); |
243 | | } |
244 | | } |
245 | | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION |
246 | | */ |
247 | 0 | #endif |
248 | 186k | } else { |
249 | 186k | if (isign == 1) |
250 | 179k | { |
251 | 443k | for (k = 0; k < l1; k++) |
252 | 263k | { |
253 | 15.7M | for (i = 0; i < ido; i++) |
254 | 15.4M | { |
255 | 15.4M | ac = i + (3*k+1)*ido; |
256 | 15.4M | ah = i + k * ido; |
257 | | |
258 | 15.4M | RE(t2) = RE(cc[ac]) + RE(cc[ac+ido]); |
259 | 15.4M | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),taur); |
260 | 15.4M | IM(t2) = IM(cc[ac]) + IM(cc[ac+ido]); |
261 | 15.4M | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),taur); |
262 | | |
263 | 15.4M | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2); |
264 | 15.4M | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2); |
265 | | |
266 | 15.4M | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+ido])), taui); |
267 | 15.4M | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+ido])), taui); |
268 | | |
269 | 15.4M | RE(d2) = RE(c2) - IM(c3); |
270 | 15.4M | IM(d3) = IM(c2) - RE(c3); |
271 | 15.4M | RE(d3) = RE(c2) + IM(c3); |
272 | 15.4M | IM(d2) = IM(c2) + RE(c3); |
273 | | |
274 | 15.4M | #if 1 |
275 | 15.4M | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), |
276 | 15.4M | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); |
277 | 15.4M | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), |
278 | 15.4M | IM(d3), RE(d3), RE(wa2[i]), IM(wa2[i])); |
279 | | #else |
280 | | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), |
281 | | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); |
282 | | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), |
283 | | RE(d3), IM(d3), RE(wa2[i]), IM(wa2[i])); |
284 | | #endif |
285 | 15.4M | } |
286 | 263k | } |
287 | 179k | } else { |
288 | 21.9k | for (k = 0; k < l1; k++) |
289 | 14.5k | { |
290 | 1.18M | for (i = 0; i < ido; i++) |
291 | 1.16M | { |
292 | 1.16M | ac = i + (3*k+1)*ido; |
293 | 1.16M | ah = i + k * ido; |
294 | | |
295 | 1.16M | RE(t2) = RE(cc[ac]) + RE(cc[ac+ido]); |
296 | 1.16M | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),taur); |
297 | 1.16M | IM(t2) = IM(cc[ac]) + IM(cc[ac+ido]); |
298 | 1.16M | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),taur); |
299 | | |
300 | 1.16M | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2); |
301 | 1.16M | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2); |
302 | | |
303 | 1.16M | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+ido])), taui); |
304 | 1.16M | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+ido])), taui); |
305 | | |
306 | 1.16M | RE(d2) = RE(c2) + IM(c3); |
307 | 1.16M | IM(d3) = IM(c2) + RE(c3); |
308 | 1.16M | RE(d3) = RE(c2) - IM(c3); |
309 | 1.16M | IM(d2) = IM(c2) - RE(c3); |
310 | | |
311 | 1.16M | #if 1 |
312 | 1.16M | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), |
313 | 1.16M | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); |
314 | 1.16M | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), |
315 | 1.16M | RE(d3), IM(d3), RE(wa2[i]), IM(wa2[i])); |
316 | | #else |
317 | | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), |
318 | | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); |
319 | | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), |
320 | | IM(d3), RE(d3), RE(wa2[i]), IM(wa2[i])); |
321 | | #endif |
322 | 1.16M | } |
323 | 14.5k | } |
324 | 7.34k | } |
325 | 186k | } |
326 | 186k | } Line | Count | Source | 184 | 61.8k | { | 185 | 61.8k | static real_t taur = FRAC_CONST(-0.5); | 186 | 61.8k | static real_t taui = FRAC_CONST(0.866025403784439); | 187 | 61.8k | uint16_t i, k, ac, ah; | 188 | 61.8k | complex_t c2, c3, d2, d3, t2; | 189 | | | 190 | 61.8k | if (ido == 1) | 191 | 0 | { | 192 | 0 | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION | 193 | | // TODO: remove this code once fuzzer proves it is totally unreahable | 194 | | // 3 is never the the biggest factor for supported frame lengths; | 195 | | // consequently `ido` should never be 1. | 196 | 0 | __builtin_trap(); | 197 | | /* | 198 | | #endif | 199 | | if (isign == 1) | 200 | | { | 201 | | for (k = 0; k < l1; k++) | 202 | | { | 203 | | ac = 3*k+1; | 204 | | ah = k; | 205 | | | 206 | | RE(t2) = RE(cc[ac]) + RE(cc[ac+1]); | 207 | | IM(t2) = IM(cc[ac]) + IM(cc[ac+1]); | 208 | | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),taur); | 209 | | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),taur); | 210 | | | 211 | | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2); | 212 | | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2); | 213 | | | 214 | | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+1])), taui); | 215 | | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+1])), taui); | 216 | | | 217 | | RE(ch[ah+l1]) = RE(c2) - IM(c3); | 218 | | IM(ch[ah+l1]) = IM(c2) + RE(c3); | 219 | | RE(ch[ah+2*l1]) = RE(c2) + IM(c3); | 220 | | IM(ch[ah+2*l1]) = IM(c2) - RE(c3); | 221 | | } | 222 | | } else { | 223 | | for (k = 0; k < l1; k++) | 224 | | { | 225 | | ac = 3*k+1; | 226 | | ah = k; | 227 | | | 228 | | RE(t2) = RE(cc[ac]) + RE(cc[ac+1]); | 229 | | IM(t2) = IM(cc[ac]) + IM(cc[ac+1]); | 230 | | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),taur); | 231 | | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),taur); | 232 | | | 233 | | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2); | 234 | | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2); | 235 | | | 236 | | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+1])), taui); | 237 | | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+1])), taui); | 238 | | | 239 | | RE(ch[ah+l1]) = RE(c2) + IM(c3); | 240 | | IM(ch[ah+l1]) = IM(c2) - RE(c3); | 241 | | RE(ch[ah+2*l1]) = RE(c2) - IM(c3); | 242 | | IM(ch[ah+2*l1]) = IM(c2) + RE(c3); | 243 | | } | 244 | | } | 245 | | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION | 246 | | */ | 247 | 0 | #endif | 248 | 61.8k | } else { | 249 | 61.8k | if (isign == 1) | 250 | 58.6k | { | 251 | 146k | for (k = 0; k < l1; k++) | 252 | 87.5k | { | 253 | 5.35M | for (i = 0; i < ido; i++) | 254 | 5.27M | { | 255 | 5.27M | ac = i + (3*k+1)*ido; | 256 | 5.27M | ah = i + k * ido; | 257 | | | 258 | 5.27M | RE(t2) = RE(cc[ac]) + RE(cc[ac+ido]); | 259 | 5.27M | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),taur); | 260 | 5.27M | IM(t2) = IM(cc[ac]) + IM(cc[ac+ido]); | 261 | 5.27M | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),taur); | 262 | | | 263 | 5.27M | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2); | 264 | 5.27M | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2); | 265 | | | 266 | 5.27M | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+ido])), taui); | 267 | 5.27M | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+ido])), taui); | 268 | | | 269 | 5.27M | RE(d2) = RE(c2) - IM(c3); | 270 | 5.27M | IM(d3) = IM(c2) - RE(c3); | 271 | 5.27M | RE(d3) = RE(c2) + IM(c3); | 272 | 5.27M | IM(d2) = IM(c2) + RE(c3); | 273 | | | 274 | 5.27M | #if 1 | 275 | 5.27M | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), | 276 | 5.27M | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); | 277 | 5.27M | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), | 278 | 5.27M | IM(d3), RE(d3), RE(wa2[i]), IM(wa2[i])); | 279 | | #else | 280 | | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), | 281 | | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); | 282 | | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), | 283 | | RE(d3), IM(d3), RE(wa2[i]), IM(wa2[i])); | 284 | | #endif | 285 | 5.27M | } | 286 | 87.5k | } | 287 | 58.6k | } else { | 288 | 9.51k | for (k = 0; k < l1; k++) | 289 | 6.33k | { | 290 | 512k | for (i = 0; i < ido; i++) | 291 | 506k | { | 292 | 506k | ac = i + (3*k+1)*ido; | 293 | 506k | ah = i + k * ido; | 294 | | | 295 | 506k | RE(t2) = RE(cc[ac]) + RE(cc[ac+ido]); | 296 | 506k | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),taur); | 297 | 506k | IM(t2) = IM(cc[ac]) + IM(cc[ac+ido]); | 298 | 506k | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),taur); | 299 | | | 300 | 506k | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2); | 301 | 506k | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2); | 302 | | | 303 | 506k | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+ido])), taui); | 304 | 506k | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+ido])), taui); | 305 | | | 306 | 506k | RE(d2) = RE(c2) + IM(c3); | 307 | 506k | IM(d3) = IM(c2) + RE(c3); | 308 | 506k | RE(d3) = RE(c2) - IM(c3); | 309 | 506k | IM(d2) = IM(c2) - RE(c3); | 310 | | | 311 | 506k | #if 1 | 312 | 506k | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), | 313 | 506k | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); | 314 | 506k | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), | 315 | 506k | RE(d3), IM(d3), RE(wa2[i]), IM(wa2[i])); | 316 | | #else | 317 | | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), | 318 | | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); | 319 | | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), | 320 | | IM(d3), RE(d3), RE(wa2[i]), IM(wa2[i])); | 321 | | #endif | 322 | 506k | } | 323 | 6.33k | } | 324 | 3.18k | } | 325 | 61.8k | } | 326 | 61.8k | } |
Line | Count | Source | 184 | 125k | { | 185 | 125k | static real_t taur = FRAC_CONST(-0.5); | 186 | 125k | static real_t taui = FRAC_CONST(0.866025403784439); | 187 | 125k | uint16_t i, k, ac, ah; | 188 | 125k | complex_t c2, c3, d2, d3, t2; | 189 | | | 190 | 125k | if (ido == 1) | 191 | 0 | { | 192 | 0 | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION | 193 | | // TODO: remove this code once fuzzer proves it is totally unreahable | 194 | | // 3 is never the the biggest factor for supported frame lengths; | 195 | | // consequently `ido` should never be 1. | 196 | 0 | __builtin_trap(); | 197 | | /* | 198 | | #endif | 199 | | if (isign == 1) | 200 | | { | 201 | | for (k = 0; k < l1; k++) | 202 | | { | 203 | | ac = 3*k+1; | 204 | | ah = k; | 205 | | | 206 | | RE(t2) = RE(cc[ac]) + RE(cc[ac+1]); | 207 | | IM(t2) = IM(cc[ac]) + IM(cc[ac+1]); | 208 | | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),taur); | 209 | | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),taur); | 210 | | | 211 | | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2); | 212 | | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2); | 213 | | | 214 | | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+1])), taui); | 215 | | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+1])), taui); | 216 | | | 217 | | RE(ch[ah+l1]) = RE(c2) - IM(c3); | 218 | | IM(ch[ah+l1]) = IM(c2) + RE(c3); | 219 | | RE(ch[ah+2*l1]) = RE(c2) + IM(c3); | 220 | | IM(ch[ah+2*l1]) = IM(c2) - RE(c3); | 221 | | } | 222 | | } else { | 223 | | for (k = 0; k < l1; k++) | 224 | | { | 225 | | ac = 3*k+1; | 226 | | ah = k; | 227 | | | 228 | | RE(t2) = RE(cc[ac]) + RE(cc[ac+1]); | 229 | | IM(t2) = IM(cc[ac]) + IM(cc[ac+1]); | 230 | | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),taur); | 231 | | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),taur); | 232 | | | 233 | | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2); | 234 | | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2); | 235 | | | 236 | | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+1])), taui); | 237 | | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+1])), taui); | 238 | | | 239 | | RE(ch[ah+l1]) = RE(c2) + IM(c3); | 240 | | IM(ch[ah+l1]) = IM(c2) - RE(c3); | 241 | | RE(ch[ah+2*l1]) = RE(c2) - IM(c3); | 242 | | IM(ch[ah+2*l1]) = IM(c2) + RE(c3); | 243 | | } | 244 | | } | 245 | | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION | 246 | | */ | 247 | 0 | #endif | 248 | 125k | } else { | 249 | 125k | if (isign == 1) | 250 | 120k | { | 251 | 296k | for (k = 0; k < l1; k++) | 252 | 176k | { | 253 | 10.3M | for (i = 0; i < ido; i++) | 254 | 10.1M | { | 255 | 10.1M | ac = i + (3*k+1)*ido; | 256 | 10.1M | ah = i + k * ido; | 257 | | | 258 | 10.1M | RE(t2) = RE(cc[ac]) + RE(cc[ac+ido]); | 259 | 10.1M | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),taur); | 260 | 10.1M | IM(t2) = IM(cc[ac]) + IM(cc[ac+ido]); | 261 | 10.1M | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),taur); | 262 | | | 263 | 10.1M | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2); | 264 | 10.1M | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2); | 265 | | | 266 | 10.1M | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+ido])), taui); | 267 | 10.1M | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+ido])), taui); | 268 | | | 269 | 10.1M | RE(d2) = RE(c2) - IM(c3); | 270 | 10.1M | IM(d3) = IM(c2) - RE(c3); | 271 | 10.1M | RE(d3) = RE(c2) + IM(c3); | 272 | 10.1M | IM(d2) = IM(c2) + RE(c3); | 273 | | | 274 | 10.1M | #if 1 | 275 | 10.1M | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), | 276 | 10.1M | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); | 277 | 10.1M | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), | 278 | 10.1M | IM(d3), RE(d3), RE(wa2[i]), IM(wa2[i])); | 279 | | #else | 280 | | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), | 281 | | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); | 282 | | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), | 283 | | RE(d3), IM(d3), RE(wa2[i]), IM(wa2[i])); | 284 | | #endif | 285 | 10.1M | } | 286 | 176k | } | 287 | 120k | } else { | 288 | 12.4k | for (k = 0; k < l1; k++) | 289 | 8.25k | { | 290 | 668k | for (i = 0; i < ido; i++) | 291 | 660k | { | 292 | 660k | ac = i + (3*k+1)*ido; | 293 | 660k | ah = i + k * ido; | 294 | | | 295 | 660k | RE(t2) = RE(cc[ac]) + RE(cc[ac+ido]); | 296 | 660k | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),taur); | 297 | 660k | IM(t2) = IM(cc[ac]) + IM(cc[ac+ido]); | 298 | 660k | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),taur); | 299 | | | 300 | 660k | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2); | 301 | 660k | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2); | 302 | | | 303 | 660k | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+ido])), taui); | 304 | 660k | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+ido])), taui); | 305 | | | 306 | 660k | RE(d2) = RE(c2) + IM(c3); | 307 | 660k | IM(d3) = IM(c2) + RE(c3); | 308 | 660k | RE(d3) = RE(c2) - IM(c3); | 309 | 660k | IM(d2) = IM(c2) - RE(c3); | 310 | | | 311 | 660k | #if 1 | 312 | 660k | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), | 313 | 660k | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); | 314 | 660k | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), | 315 | 660k | RE(d3), IM(d3), RE(wa2[i]), IM(wa2[i])); | 316 | | #else | 317 | | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), | 318 | | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); | 319 | | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), | 320 | | IM(d3), RE(d3), RE(wa2[i]), IM(wa2[i])); | 321 | | #endif | 322 | 660k | } | 323 | 8.25k | } | 324 | 4.16k | } | 325 | 125k | } | 326 | 125k | } |
|
327 | | |
328 | | |
329 | | static void passf4pos(const uint16_t ido, const uint16_t l1, const complex_t *cc, |
330 | | complex_t *ch, const complex_t *wa1, const complex_t *wa2, |
331 | | const complex_t *wa3) |
332 | 3.38M | { |
333 | 3.38M | uint16_t i, k, ac, ah; |
334 | | |
335 | 3.38M | if (ido == 1) |
336 | 879k | { |
337 | 68.2M | for (k = 0; k < l1; k++) |
338 | 67.4M | { |
339 | 67.4M | complex_t t1, t2, t3, t4; |
340 | | |
341 | 67.4M | ac = 4*k; |
342 | 67.4M | ah = k; |
343 | | |
344 | 67.4M | RE(t2) = RE(cc[ac]) + RE(cc[ac+2]); |
345 | 67.4M | RE(t1) = RE(cc[ac]) - RE(cc[ac+2]); |
346 | 67.4M | IM(t2) = IM(cc[ac]) + IM(cc[ac+2]); |
347 | 67.4M | IM(t1) = IM(cc[ac]) - IM(cc[ac+2]); |
348 | 67.4M | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+3]); |
349 | 67.4M | IM(t4) = RE(cc[ac+1]) - RE(cc[ac+3]); |
350 | 67.4M | IM(t3) = IM(cc[ac+3]) + IM(cc[ac+1]); |
351 | 67.4M | RE(t4) = IM(cc[ac+3]) - IM(cc[ac+1]); |
352 | | |
353 | 67.4M | RE(ch[ah]) = RE(t2) + RE(t3); |
354 | 67.4M | RE(ch[ah+2*l1]) = RE(t2) - RE(t3); |
355 | | |
356 | 67.4M | IM(ch[ah]) = IM(t2) + IM(t3); |
357 | 67.4M | IM(ch[ah+2*l1]) = IM(t2) - IM(t3); |
358 | | |
359 | 67.4M | RE(ch[ah+l1]) = RE(t1) + RE(t4); |
360 | 67.4M | RE(ch[ah+3*l1]) = RE(t1) - RE(t4); |
361 | | |
362 | 67.4M | IM(ch[ah+l1]) = IM(t1) + IM(t4); |
363 | 67.4M | IM(ch[ah+3*l1]) = IM(t1) - IM(t4); |
364 | 67.4M | } |
365 | 2.50M | } else { |
366 | 27.3M | for (k = 0; k < l1; k++) |
367 | 24.8M | { |
368 | 24.8M | ac = 4*k*ido; |
369 | 24.8M | ah = k*ido; |
370 | | |
371 | 242M | for (i = 0; i < ido; i++) |
372 | 217M | { |
373 | 217M | complex_t c2, c3, c4, t1, t2, t3, t4; |
374 | | |
375 | 217M | RE(t2) = RE(cc[ac+i]) + RE(cc[ac+i+2*ido]); |
376 | 217M | RE(t1) = RE(cc[ac+i]) - RE(cc[ac+i+2*ido]); |
377 | 217M | IM(t2) = IM(cc[ac+i]) + IM(cc[ac+i+2*ido]); |
378 | 217M | IM(t1) = IM(cc[ac+i]) - IM(cc[ac+i+2*ido]); |
379 | 217M | RE(t3) = RE(cc[ac+i+ido]) + RE(cc[ac+i+3*ido]); |
380 | 217M | IM(t4) = RE(cc[ac+i+ido]) - RE(cc[ac+i+3*ido]); |
381 | 217M | IM(t3) = IM(cc[ac+i+3*ido]) + IM(cc[ac+i+ido]); |
382 | 217M | RE(t4) = IM(cc[ac+i+3*ido]) - IM(cc[ac+i+ido]); |
383 | | |
384 | 217M | RE(c2) = RE(t1) + RE(t4); |
385 | 217M | RE(c4) = RE(t1) - RE(t4); |
386 | | |
387 | 217M | IM(c2) = IM(t1) + IM(t4); |
388 | 217M | IM(c4) = IM(t1) - IM(t4); |
389 | | |
390 | 217M | RE(ch[ah+i]) = RE(t2) + RE(t3); |
391 | 217M | RE(c3) = RE(t2) - RE(t3); |
392 | | |
393 | 217M | IM(ch[ah+i]) = IM(t2) + IM(t3); |
394 | 217M | IM(c3) = IM(t2) - IM(t3); |
395 | | |
396 | 217M | #if 1 |
397 | 217M | ComplexMult(&IM(ch[ah+i+l1*ido]), &RE(ch[ah+i+l1*ido]), |
398 | 217M | IM(c2), RE(c2), RE(wa1[i]), IM(wa1[i])); |
399 | 217M | ComplexMult(&IM(ch[ah+i+2*l1*ido]), &RE(ch[ah+i+2*l1*ido]), |
400 | 217M | IM(c3), RE(c3), RE(wa2[i]), IM(wa2[i])); |
401 | 217M | ComplexMult(&IM(ch[ah+i+3*l1*ido]), &RE(ch[ah+i+3*l1*ido]), |
402 | 217M | IM(c4), RE(c4), RE(wa3[i]), IM(wa3[i])); |
403 | | #else |
404 | | ComplexMult(&RE(ch[ah+i+l1*ido]), &IM(ch[ah+i+l1*ido]), |
405 | | RE(c2), IM(c2), RE(wa1[i]), IM(wa1[i])); |
406 | | ComplexMult(&RE(ch[ah+i+2*l1*ido]), &IM(ch[ah+i+2*l1*ido]), |
407 | | RE(c3), IM(c3), RE(wa2[i]), IM(wa2[i])); |
408 | | ComplexMult(&RE(ch[ah+i+3*l1*ido]), &IM(ch[ah+i+3*l1*ido]), |
409 | | RE(c4), IM(c4), RE(wa3[i]), IM(wa3[i])); |
410 | | #endif |
411 | 217M | } |
412 | 24.8M | } |
413 | 2.50M | } |
414 | 3.38M | } |
415 | | |
416 | | static void passf4neg(const uint16_t ido, const uint16_t l1, const complex_t *cc, |
417 | | complex_t *ch, const complex_t *wa1, const complex_t *wa2, |
418 | | const complex_t *wa3) |
419 | 35.7k | { |
420 | 35.7k | uint16_t i, k, ac, ah; |
421 | | |
422 | 35.7k | if (ido == 1) |
423 | 5.27k | { |
424 | 662k | for (k = 0; k < l1; k++) |
425 | 657k | { |
426 | 657k | complex_t t1, t2, t3, t4; |
427 | | |
428 | 657k | ac = 4*k; |
429 | 657k | ah = k; |
430 | | |
431 | 657k | RE(t2) = RE(cc[ac]) + RE(cc[ac+2]); |
432 | 657k | RE(t1) = RE(cc[ac]) - RE(cc[ac+2]); |
433 | 657k | IM(t2) = IM(cc[ac]) + IM(cc[ac+2]); |
434 | 657k | IM(t1) = IM(cc[ac]) - IM(cc[ac+2]); |
435 | 657k | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+3]); |
436 | 657k | IM(t4) = RE(cc[ac+1]) - RE(cc[ac+3]); |
437 | 657k | IM(t3) = IM(cc[ac+3]) + IM(cc[ac+1]); |
438 | 657k | RE(t4) = IM(cc[ac+3]) - IM(cc[ac+1]); |
439 | | |
440 | 657k | RE(ch[ah]) = RE(t2) + RE(t3); |
441 | 657k | RE(ch[ah+2*l1]) = RE(t2) - RE(t3); |
442 | | |
443 | 657k | IM(ch[ah]) = IM(t2) + IM(t3); |
444 | 657k | IM(ch[ah+2*l1]) = IM(t2) - IM(t3); |
445 | | |
446 | 657k | RE(ch[ah+l1]) = RE(t1) - RE(t4); |
447 | 657k | RE(ch[ah+3*l1]) = RE(t1) + RE(t4); |
448 | | |
449 | 657k | IM(ch[ah+l1]) = IM(t1) - IM(t4); |
450 | 657k | IM(ch[ah+3*l1]) = IM(t1) + IM(t4); |
451 | 657k | } |
452 | 30.4k | } else { |
453 | 464k | for (k = 0; k < l1; k++) |
454 | 434k | { |
455 | 434k | ac = 4*k*ido; |
456 | 434k | ah = k*ido; |
457 | | |
458 | 4.15M | for (i = 0; i < ido; i++) |
459 | 3.72M | { |
460 | 3.72M | complex_t c2, c3, c4, t1, t2, t3, t4; |
461 | | |
462 | 3.72M | RE(t2) = RE(cc[ac+i]) + RE(cc[ac+i+2*ido]); |
463 | 3.72M | RE(t1) = RE(cc[ac+i]) - RE(cc[ac+i+2*ido]); |
464 | 3.72M | IM(t2) = IM(cc[ac+i]) + IM(cc[ac+i+2*ido]); |
465 | 3.72M | IM(t1) = IM(cc[ac+i]) - IM(cc[ac+i+2*ido]); |
466 | 3.72M | RE(t3) = RE(cc[ac+i+ido]) + RE(cc[ac+i+3*ido]); |
467 | 3.72M | IM(t4) = RE(cc[ac+i+ido]) - RE(cc[ac+i+3*ido]); |
468 | 3.72M | IM(t3) = IM(cc[ac+i+3*ido]) + IM(cc[ac+i+ido]); |
469 | 3.72M | RE(t4) = IM(cc[ac+i+3*ido]) - IM(cc[ac+i+ido]); |
470 | | |
471 | 3.72M | RE(c2) = RE(t1) - RE(t4); |
472 | 3.72M | RE(c4) = RE(t1) + RE(t4); |
473 | | |
474 | 3.72M | IM(c2) = IM(t1) - IM(t4); |
475 | 3.72M | IM(c4) = IM(t1) + IM(t4); |
476 | | |
477 | 3.72M | RE(ch[ah+i]) = RE(t2) + RE(t3); |
478 | 3.72M | RE(c3) = RE(t2) - RE(t3); |
479 | | |
480 | 3.72M | IM(ch[ah+i]) = IM(t2) + IM(t3); |
481 | 3.72M | IM(c3) = IM(t2) - IM(t3); |
482 | | |
483 | 3.72M | #if 1 |
484 | 3.72M | ComplexMult(&RE(ch[ah+i+l1*ido]), &IM(ch[ah+i+l1*ido]), |
485 | 3.72M | RE(c2), IM(c2), RE(wa1[i]), IM(wa1[i])); |
486 | 3.72M | ComplexMult(&RE(ch[ah+i+2*l1*ido]), &IM(ch[ah+i+2*l1*ido]), |
487 | 3.72M | RE(c3), IM(c3), RE(wa2[i]), IM(wa2[i])); |
488 | 3.72M | ComplexMult(&RE(ch[ah+i+3*l1*ido]), &IM(ch[ah+i+3*l1*ido]), |
489 | 3.72M | RE(c4), IM(c4), RE(wa3[i]), IM(wa3[i])); |
490 | | #else |
491 | | ComplexMult(&IM(ch[ah+i+l1*ido]), &RE(ch[ah+i+l1*ido]), |
492 | | IM(c2), RE(c2), RE(wa1[i]), IM(wa1[i])); |
493 | | ComplexMult(&IM(ch[ah+i+2*l1*ido]), &RE(ch[ah+i+2*l1*ido]), |
494 | | IM(c3), RE(c3), RE(wa2[i]), IM(wa2[i])); |
495 | | ComplexMult(&IM(ch[ah+i+3*l1*ido]), &RE(ch[ah+i+3*l1*ido]), |
496 | | IM(c4), RE(c4), RE(wa3[i]), IM(wa3[i])); |
497 | | #endif |
498 | 3.72M | } |
499 | 434k | } |
500 | 30.4k | } |
501 | 35.7k | } |
502 | | |
503 | | static void passf5(const uint16_t ido, const uint16_t l1, const complex_t *cc, |
504 | | complex_t *ch, const complex_t *wa1, const complex_t *wa2, const complex_t *wa3, |
505 | | const complex_t *wa4, const int8_t isign) |
506 | 186k | { |
507 | 186k | static real_t tr11 = FRAC_CONST(0.309016994374947); |
508 | 186k | static real_t ti11 = FRAC_CONST(0.951056516295154); |
509 | 186k | static real_t tr12 = FRAC_CONST(-0.809016994374947); |
510 | 186k | static real_t ti12 = FRAC_CONST(0.587785252292473); |
511 | 186k | uint16_t k, ac, ah; |
512 | 186k | complex_t c2, c3, c4, c5, t2, t3, t4, t5; |
513 | 186k | (void)wa1; |
514 | 186k | (void)wa2; |
515 | 186k | (void)wa3; |
516 | 186k | (void)wa4; |
517 | | |
518 | 186k | if (ido == 1) |
519 | 186k | { |
520 | 186k | if (isign == 1) |
521 | 179k | { |
522 | 9.46M | for (k = 0; k < l1; k++) |
523 | 9.28M | { |
524 | 9.28M | ac = 5*k + 1; |
525 | 9.28M | ah = k; |
526 | | |
527 | 9.28M | RE(t2) = RE(cc[ac]) + RE(cc[ac+3]); |
528 | 9.28M | IM(t2) = IM(cc[ac]) + IM(cc[ac+3]); |
529 | 9.28M | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+2]); |
530 | 9.28M | IM(t3) = IM(cc[ac+1]) + IM(cc[ac+2]); |
531 | 9.28M | RE(t4) = RE(cc[ac+1]) - RE(cc[ac+2]); |
532 | 9.28M | IM(t4) = IM(cc[ac+1]) - IM(cc[ac+2]); |
533 | 9.28M | RE(t5) = RE(cc[ac]) - RE(cc[ac+3]); |
534 | 9.28M | IM(t5) = IM(cc[ac]) - IM(cc[ac+3]); |
535 | | |
536 | 9.28M | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2) + RE(t3); |
537 | 9.28M | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2) + IM(t3); |
538 | | |
539 | 9.28M | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); |
540 | 9.28M | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); |
541 | 9.28M | RE(c3) = RE(cc[ac-1]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); |
542 | 9.28M | IM(c3) = IM(cc[ac-1]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); |
543 | | |
544 | 9.28M | ComplexMult(&RE(c5), &RE(c4), |
545 | 9.28M | ti11, ti12, RE(t5), RE(t4)); |
546 | 9.28M | ComplexMult(&IM(c5), &IM(c4), |
547 | 9.28M | ti11, ti12, IM(t5), IM(t4)); |
548 | | |
549 | 9.28M | RE(ch[ah+l1]) = RE(c2) - IM(c5); |
550 | 9.28M | IM(ch[ah+l1]) = IM(c2) + RE(c5); |
551 | 9.28M | RE(ch[ah+2*l1]) = RE(c3) - IM(c4); |
552 | 9.28M | IM(ch[ah+2*l1]) = IM(c3) + RE(c4); |
553 | 9.28M | RE(ch[ah+3*l1]) = RE(c3) + IM(c4); |
554 | 9.28M | IM(ch[ah+3*l1]) = IM(c3) - RE(c4); |
555 | 9.28M | RE(ch[ah+4*l1]) = RE(c2) + IM(c5); |
556 | 9.28M | IM(ch[ah+4*l1]) = IM(c2) - RE(c5); |
557 | 9.28M | } |
558 | 179k | } else { |
559 | 707k | for (k = 0; k < l1; k++) |
560 | 699k | { |
561 | 699k | ac = 5*k + 1; |
562 | 699k | ah = k; |
563 | | |
564 | 699k | RE(t2) = RE(cc[ac]) + RE(cc[ac+3]); |
565 | 699k | IM(t2) = IM(cc[ac]) + IM(cc[ac+3]); |
566 | 699k | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+2]); |
567 | 699k | IM(t3) = IM(cc[ac+1]) + IM(cc[ac+2]); |
568 | 699k | RE(t4) = RE(cc[ac+1]) - RE(cc[ac+2]); |
569 | 699k | IM(t4) = IM(cc[ac+1]) - IM(cc[ac+2]); |
570 | 699k | RE(t5) = RE(cc[ac]) - RE(cc[ac+3]); |
571 | 699k | IM(t5) = IM(cc[ac]) - IM(cc[ac+3]); |
572 | | |
573 | 699k | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2) + RE(t3); |
574 | 699k | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2) + IM(t3); |
575 | | |
576 | 699k | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); |
577 | 699k | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); |
578 | 699k | RE(c3) = RE(cc[ac-1]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); |
579 | 699k | IM(c3) = IM(cc[ac-1]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); |
580 | | |
581 | 699k | ComplexMult(&RE(c4), &RE(c5), |
582 | 699k | ti12, ti11, RE(t5), RE(t4)); |
583 | 699k | ComplexMult(&IM(c4), &IM(c5), |
584 | 699k | ti12, ti11, IM(t5), IM(t4)); |
585 | | |
586 | 699k | RE(ch[ah+l1]) = RE(c2) + IM(c5); |
587 | 699k | IM(ch[ah+l1]) = IM(c2) - RE(c5); |
588 | 699k | RE(ch[ah+2*l1]) = RE(c3) + IM(c4); |
589 | 699k | IM(ch[ah+2*l1]) = IM(c3) - RE(c4); |
590 | 699k | RE(ch[ah+3*l1]) = RE(c3) - IM(c4); |
591 | 699k | IM(ch[ah+3*l1]) = IM(c3) + RE(c4); |
592 | 699k | RE(ch[ah+4*l1]) = RE(c2) - IM(c5); |
593 | 699k | IM(ch[ah+4*l1]) = IM(c2) + RE(c5); |
594 | 699k | } |
595 | 7.34k | } |
596 | 186k | } else { |
597 | 0 | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION |
598 | | // TODO: remove this code once fuzzer proves it is totally unreahable |
599 | | // 5 if the biggest factor and it never repeated for supported frame |
600 | | // lengths; consequently `ido` should always be 1. |
601 | 0 | __builtin_trap(); |
602 | | /* |
603 | | #else |
604 | | uint16_t i; |
605 | | complex_t d3, d4, d5, d2; |
606 | | #endif |
607 | | if (isign == 1) |
608 | | { |
609 | | for (k = 0; k < l1; k++) |
610 | | { |
611 | | for (i = 0; i < ido; i++) |
612 | | { |
613 | | ac = i + (k*5 + 1) * ido; |
614 | | ah = i + k * ido; |
615 | | |
616 | | RE(t2) = RE(cc[ac]) + RE(cc[ac+3*ido]); |
617 | | IM(t2) = IM(cc[ac]) + IM(cc[ac+3*ido]); |
618 | | RE(t3) = RE(cc[ac+ido]) + RE(cc[ac+2*ido]); |
619 | | IM(t3) = IM(cc[ac+ido]) + IM(cc[ac+2*ido]); |
620 | | RE(t4) = RE(cc[ac+ido]) - RE(cc[ac+2*ido]); |
621 | | IM(t4) = IM(cc[ac+ido]) - IM(cc[ac+2*ido]); |
622 | | RE(t5) = RE(cc[ac]) - RE(cc[ac+3*ido]); |
623 | | IM(t5) = IM(cc[ac]) - IM(cc[ac+3*ido]); |
624 | | |
625 | | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2) + RE(t3); |
626 | | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2) + IM(t3); |
627 | | |
628 | | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); |
629 | | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); |
630 | | RE(c3) = RE(cc[ac-ido]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); |
631 | | IM(c3) = IM(cc[ac-ido]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); |
632 | | |
633 | | ComplexMult(&RE(c5), &RE(c4), |
634 | | ti11, ti12, RE(t5), RE(t4)); |
635 | | ComplexMult(&IM(c5), &IM(c4), |
636 | | ti11, ti12, IM(t5), IM(t4)); |
637 | | |
638 | | IM(d2) = IM(c2) + RE(c5); |
639 | | IM(d3) = IM(c3) + RE(c4); |
640 | | RE(d4) = RE(c3) + IM(c4); |
641 | | RE(d5) = RE(c2) + IM(c5); |
642 | | RE(d2) = RE(c2) - IM(c5); |
643 | | IM(d5) = IM(c2) - RE(c5); |
644 | | RE(d3) = RE(c3) - IM(c4); |
645 | | IM(d4) = IM(c3) - RE(c4); |
646 | | |
647 | | #if 1 |
648 | | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), |
649 | | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); |
650 | | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), |
651 | | IM(d3), RE(d3), RE(wa2[i]), IM(wa2[i])); |
652 | | ComplexMult(&IM(ch[ah+3*l1*ido]), &RE(ch[ah+3*l1*ido]), |
653 | | IM(d4), RE(d4), RE(wa3[i]), IM(wa3[i])); |
654 | | ComplexMult(&IM(ch[ah+4*l1*ido]), &RE(ch[ah+4*l1*ido]), |
655 | | IM(d5), RE(d5), RE(wa4[i]), IM(wa4[i])); |
656 | | #else |
657 | | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), |
658 | | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); |
659 | | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), |
660 | | RE(d3), IM(d3), RE(wa2[i]), IM(wa2[i])); |
661 | | ComplexMult(&RE(ch[ah+3*l1*ido]), &IM(ch[ah+3*l1*ido]), |
662 | | RE(d4), IM(d4), RE(wa3[i]), IM(wa3[i])); |
663 | | ComplexMult(&RE(ch[ah+4*l1*ido]), &IM(ch[ah+4*l1*ido]), |
664 | | RE(d5), IM(d5), RE(wa4[i]), IM(wa4[i])); |
665 | | #endif |
666 | | } |
667 | | } |
668 | | } else { |
669 | | for (k = 0; k < l1; k++) |
670 | | { |
671 | | for (i = 0; i < ido; i++) |
672 | | { |
673 | | ac = i + (k*5 + 1) * ido; |
674 | | ah = i + k * ido; |
675 | | |
676 | | RE(t2) = RE(cc[ac]) + RE(cc[ac+3*ido]); |
677 | | IM(t2) = IM(cc[ac]) + IM(cc[ac+3*ido]); |
678 | | RE(t3) = RE(cc[ac+ido]) + RE(cc[ac+2*ido]); |
679 | | IM(t3) = IM(cc[ac+ido]) + IM(cc[ac+2*ido]); |
680 | | RE(t4) = RE(cc[ac+ido]) - RE(cc[ac+2*ido]); |
681 | | IM(t4) = IM(cc[ac+ido]) - IM(cc[ac+2*ido]); |
682 | | RE(t5) = RE(cc[ac]) - RE(cc[ac+3*ido]); |
683 | | IM(t5) = IM(cc[ac]) - IM(cc[ac+3*ido]); |
684 | | |
685 | | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2) + RE(t3); |
686 | | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2) + IM(t3); |
687 | | |
688 | | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); |
689 | | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); |
690 | | RE(c3) = RE(cc[ac-ido]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); |
691 | | IM(c3) = IM(cc[ac-ido]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); |
692 | | |
693 | | ComplexMult(&RE(c4), &RE(c5), |
694 | | ti12, ti11, RE(t5), RE(t4)); |
695 | | ComplexMult(&IM(c4), &IM(c5), |
696 | | ti12, ti11, IM(t5), IM(t4)); |
697 | | |
698 | | IM(d2) = IM(c2) - RE(c5); |
699 | | IM(d3) = IM(c3) - RE(c4); |
700 | | RE(d4) = RE(c3) - IM(c4); |
701 | | RE(d5) = RE(c2) - IM(c5); |
702 | | RE(d2) = RE(c2) + IM(c5); |
703 | | IM(d5) = IM(c2) + RE(c5); |
704 | | RE(d3) = RE(c3) + IM(c4); |
705 | | IM(d4) = IM(c3) + RE(c4); |
706 | | |
707 | | #if 1 |
708 | | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), |
709 | | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); |
710 | | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), |
711 | | RE(d3), IM(d3), RE(wa2[i]), IM(wa2[i])); |
712 | | ComplexMult(&RE(ch[ah+3*l1*ido]), &IM(ch[ah+3*l1*ido]), |
713 | | RE(d4), IM(d4), RE(wa3[i]), IM(wa3[i])); |
714 | | ComplexMult(&RE(ch[ah+4*l1*ido]), &IM(ch[ah+4*l1*ido]), |
715 | | RE(d5), IM(d5), RE(wa4[i]), IM(wa4[i])); |
716 | | #else |
717 | | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), |
718 | | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); |
719 | | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), |
720 | | IM(d3), RE(d3), RE(wa2[i]), IM(wa2[i])); |
721 | | ComplexMult(&IM(ch[ah+3*l1*ido]), &RE(ch[ah+3*l1*ido]), |
722 | | IM(d4), RE(d4), RE(wa3[i]), IM(wa3[i])); |
723 | | ComplexMult(&IM(ch[ah+4*l1*ido]), &RE(ch[ah+4*l1*ido]), |
724 | | IM(d5), RE(d5), RE(wa4[i]), IM(wa4[i])); |
725 | | #endif |
726 | | } |
727 | | } |
728 | | } |
729 | | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION |
730 | | */ |
731 | 0 | #endif |
732 | 0 | } |
733 | 186k | } Line | Count | Source | 506 | 61.8k | { | 507 | 61.8k | static real_t tr11 = FRAC_CONST(0.309016994374947); | 508 | 61.8k | static real_t ti11 = FRAC_CONST(0.951056516295154); | 509 | 61.8k | static real_t tr12 = FRAC_CONST(-0.809016994374947); | 510 | 61.8k | static real_t ti12 = FRAC_CONST(0.587785252292473); | 511 | 61.8k | uint16_t k, ac, ah; | 512 | 61.8k | complex_t c2, c3, c4, c5, t2, t3, t4, t5; | 513 | 61.8k | (void)wa1; | 514 | 61.8k | (void)wa2; | 515 | 61.8k | (void)wa3; | 516 | 61.8k | (void)wa4; | 517 | | | 518 | 61.8k | if (ido == 1) | 519 | 61.8k | { | 520 | 61.8k | if (isign == 1) | 521 | 58.6k | { | 522 | 3.22M | for (k = 0; k < l1; k++) | 523 | 3.16M | { | 524 | 3.16M | ac = 5*k + 1; | 525 | 3.16M | ah = k; | 526 | | | 527 | 3.16M | RE(t2) = RE(cc[ac]) + RE(cc[ac+3]); | 528 | 3.16M | IM(t2) = IM(cc[ac]) + IM(cc[ac+3]); | 529 | 3.16M | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+2]); | 530 | 3.16M | IM(t3) = IM(cc[ac+1]) + IM(cc[ac+2]); | 531 | 3.16M | RE(t4) = RE(cc[ac+1]) - RE(cc[ac+2]); | 532 | 3.16M | IM(t4) = IM(cc[ac+1]) - IM(cc[ac+2]); | 533 | 3.16M | RE(t5) = RE(cc[ac]) - RE(cc[ac+3]); | 534 | 3.16M | IM(t5) = IM(cc[ac]) - IM(cc[ac+3]); | 535 | | | 536 | 3.16M | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2) + RE(t3); | 537 | 3.16M | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2) + IM(t3); | 538 | | | 539 | 3.16M | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); | 540 | 3.16M | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); | 541 | 3.16M | RE(c3) = RE(cc[ac-1]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); | 542 | 3.16M | IM(c3) = IM(cc[ac-1]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); | 543 | | | 544 | 3.16M | ComplexMult(&RE(c5), &RE(c4), | 545 | 3.16M | ti11, ti12, RE(t5), RE(t4)); | 546 | 3.16M | ComplexMult(&IM(c5), &IM(c4), | 547 | 3.16M | ti11, ti12, IM(t5), IM(t4)); | 548 | | | 549 | 3.16M | RE(ch[ah+l1]) = RE(c2) - IM(c5); | 550 | 3.16M | IM(ch[ah+l1]) = IM(c2) + RE(c5); | 551 | 3.16M | RE(ch[ah+2*l1]) = RE(c3) - IM(c4); | 552 | 3.16M | IM(ch[ah+2*l1]) = IM(c3) + RE(c4); | 553 | 3.16M | RE(ch[ah+3*l1]) = RE(c3) + IM(c4); | 554 | 3.16M | IM(ch[ah+3*l1]) = IM(c3) - RE(c4); | 555 | 3.16M | RE(ch[ah+4*l1]) = RE(c2) + IM(c5); | 556 | 3.16M | IM(ch[ah+4*l1]) = IM(c2) - RE(c5); | 557 | 3.16M | } | 558 | 58.6k | } else { | 559 | 307k | for (k = 0; k < l1; k++) | 560 | 303k | { | 561 | 303k | ac = 5*k + 1; | 562 | 303k | ah = k; | 563 | | | 564 | 303k | RE(t2) = RE(cc[ac]) + RE(cc[ac+3]); | 565 | 303k | IM(t2) = IM(cc[ac]) + IM(cc[ac+3]); | 566 | 303k | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+2]); | 567 | 303k | IM(t3) = IM(cc[ac+1]) + IM(cc[ac+2]); | 568 | 303k | RE(t4) = RE(cc[ac+1]) - RE(cc[ac+2]); | 569 | 303k | IM(t4) = IM(cc[ac+1]) - IM(cc[ac+2]); | 570 | 303k | RE(t5) = RE(cc[ac]) - RE(cc[ac+3]); | 571 | 303k | IM(t5) = IM(cc[ac]) - IM(cc[ac+3]); | 572 | | | 573 | 303k | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2) + RE(t3); | 574 | 303k | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2) + IM(t3); | 575 | | | 576 | 303k | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); | 577 | 303k | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); | 578 | 303k | RE(c3) = RE(cc[ac-1]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); | 579 | 303k | IM(c3) = IM(cc[ac-1]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); | 580 | | | 581 | 303k | ComplexMult(&RE(c4), &RE(c5), | 582 | 303k | ti12, ti11, RE(t5), RE(t4)); | 583 | 303k | ComplexMult(&IM(c4), &IM(c5), | 584 | 303k | ti12, ti11, IM(t5), IM(t4)); | 585 | | | 586 | 303k | RE(ch[ah+l1]) = RE(c2) + IM(c5); | 587 | 303k | IM(ch[ah+l1]) = IM(c2) - RE(c5); | 588 | 303k | RE(ch[ah+2*l1]) = RE(c3) + IM(c4); | 589 | 303k | IM(ch[ah+2*l1]) = IM(c3) - RE(c4); | 590 | 303k | RE(ch[ah+3*l1]) = RE(c3) - IM(c4); | 591 | 303k | IM(ch[ah+3*l1]) = IM(c3) + RE(c4); | 592 | 303k | RE(ch[ah+4*l1]) = RE(c2) - IM(c5); | 593 | 303k | IM(ch[ah+4*l1]) = IM(c2) + RE(c5); | 594 | 303k | } | 595 | 3.18k | } | 596 | 61.8k | } else { | 597 | 0 | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION | 598 | | // TODO: remove this code once fuzzer proves it is totally unreahable | 599 | | // 5 if the biggest factor and it never repeated for supported frame | 600 | | // lengths; consequently `ido` should always be 1. | 601 | 0 | __builtin_trap(); | 602 | | /* | 603 | | #else | 604 | | uint16_t i; | 605 | | complex_t d3, d4, d5, d2; | 606 | | #endif | 607 | | if (isign == 1) | 608 | | { | 609 | | for (k = 0; k < l1; k++) | 610 | | { | 611 | | for (i = 0; i < ido; i++) | 612 | | { | 613 | | ac = i + (k*5 + 1) * ido; | 614 | | ah = i + k * ido; | 615 | | | 616 | | RE(t2) = RE(cc[ac]) + RE(cc[ac+3*ido]); | 617 | | IM(t2) = IM(cc[ac]) + IM(cc[ac+3*ido]); | 618 | | RE(t3) = RE(cc[ac+ido]) + RE(cc[ac+2*ido]); | 619 | | IM(t3) = IM(cc[ac+ido]) + IM(cc[ac+2*ido]); | 620 | | RE(t4) = RE(cc[ac+ido]) - RE(cc[ac+2*ido]); | 621 | | IM(t4) = IM(cc[ac+ido]) - IM(cc[ac+2*ido]); | 622 | | RE(t5) = RE(cc[ac]) - RE(cc[ac+3*ido]); | 623 | | IM(t5) = IM(cc[ac]) - IM(cc[ac+3*ido]); | 624 | | | 625 | | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2) + RE(t3); | 626 | | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2) + IM(t3); | 627 | | | 628 | | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); | 629 | | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); | 630 | | RE(c3) = RE(cc[ac-ido]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); | 631 | | IM(c3) = IM(cc[ac-ido]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); | 632 | | | 633 | | ComplexMult(&RE(c5), &RE(c4), | 634 | | ti11, ti12, RE(t5), RE(t4)); | 635 | | ComplexMult(&IM(c5), &IM(c4), | 636 | | ti11, ti12, IM(t5), IM(t4)); | 637 | | | 638 | | IM(d2) = IM(c2) + RE(c5); | 639 | | IM(d3) = IM(c3) + RE(c4); | 640 | | RE(d4) = RE(c3) + IM(c4); | 641 | | RE(d5) = RE(c2) + IM(c5); | 642 | | RE(d2) = RE(c2) - IM(c5); | 643 | | IM(d5) = IM(c2) - RE(c5); | 644 | | RE(d3) = RE(c3) - IM(c4); | 645 | | IM(d4) = IM(c3) - RE(c4); | 646 | | | 647 | | #if 1 | 648 | | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), | 649 | | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); | 650 | | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), | 651 | | IM(d3), RE(d3), RE(wa2[i]), IM(wa2[i])); | 652 | | ComplexMult(&IM(ch[ah+3*l1*ido]), &RE(ch[ah+3*l1*ido]), | 653 | | IM(d4), RE(d4), RE(wa3[i]), IM(wa3[i])); | 654 | | ComplexMult(&IM(ch[ah+4*l1*ido]), &RE(ch[ah+4*l1*ido]), | 655 | | IM(d5), RE(d5), RE(wa4[i]), IM(wa4[i])); | 656 | | #else | 657 | | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), | 658 | | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); | 659 | | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), | 660 | | RE(d3), IM(d3), RE(wa2[i]), IM(wa2[i])); | 661 | | ComplexMult(&RE(ch[ah+3*l1*ido]), &IM(ch[ah+3*l1*ido]), | 662 | | RE(d4), IM(d4), RE(wa3[i]), IM(wa3[i])); | 663 | | ComplexMult(&RE(ch[ah+4*l1*ido]), &IM(ch[ah+4*l1*ido]), | 664 | | RE(d5), IM(d5), RE(wa4[i]), IM(wa4[i])); | 665 | | #endif | 666 | | } | 667 | | } | 668 | | } else { | 669 | | for (k = 0; k < l1; k++) | 670 | | { | 671 | | for (i = 0; i < ido; i++) | 672 | | { | 673 | | ac = i + (k*5 + 1) * ido; | 674 | | ah = i + k * ido; | 675 | | | 676 | | RE(t2) = RE(cc[ac]) + RE(cc[ac+3*ido]); | 677 | | IM(t2) = IM(cc[ac]) + IM(cc[ac+3*ido]); | 678 | | RE(t3) = RE(cc[ac+ido]) + RE(cc[ac+2*ido]); | 679 | | IM(t3) = IM(cc[ac+ido]) + IM(cc[ac+2*ido]); | 680 | | RE(t4) = RE(cc[ac+ido]) - RE(cc[ac+2*ido]); | 681 | | IM(t4) = IM(cc[ac+ido]) - IM(cc[ac+2*ido]); | 682 | | RE(t5) = RE(cc[ac]) - RE(cc[ac+3*ido]); | 683 | | IM(t5) = IM(cc[ac]) - IM(cc[ac+3*ido]); | 684 | | | 685 | | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2) + RE(t3); | 686 | | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2) + IM(t3); | 687 | | | 688 | | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); | 689 | | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); | 690 | | RE(c3) = RE(cc[ac-ido]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); | 691 | | IM(c3) = IM(cc[ac-ido]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); | 692 | | | 693 | | ComplexMult(&RE(c4), &RE(c5), | 694 | | ti12, ti11, RE(t5), RE(t4)); | 695 | | ComplexMult(&IM(c4), &IM(c5), | 696 | | ti12, ti11, IM(t5), IM(t4)); | 697 | | | 698 | | IM(d2) = IM(c2) - RE(c5); | 699 | | IM(d3) = IM(c3) - RE(c4); | 700 | | RE(d4) = RE(c3) - IM(c4); | 701 | | RE(d5) = RE(c2) - IM(c5); | 702 | | RE(d2) = RE(c2) + IM(c5); | 703 | | IM(d5) = IM(c2) + RE(c5); | 704 | | RE(d3) = RE(c3) + IM(c4); | 705 | | IM(d4) = IM(c3) + RE(c4); | 706 | | | 707 | | #if 1 | 708 | | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), | 709 | | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); | 710 | | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), | 711 | | RE(d3), IM(d3), RE(wa2[i]), IM(wa2[i])); | 712 | | ComplexMult(&RE(ch[ah+3*l1*ido]), &IM(ch[ah+3*l1*ido]), | 713 | | RE(d4), IM(d4), RE(wa3[i]), IM(wa3[i])); | 714 | | ComplexMult(&RE(ch[ah+4*l1*ido]), &IM(ch[ah+4*l1*ido]), | 715 | | RE(d5), IM(d5), RE(wa4[i]), IM(wa4[i])); | 716 | | #else | 717 | | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), | 718 | | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); | 719 | | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), | 720 | | IM(d3), RE(d3), RE(wa2[i]), IM(wa2[i])); | 721 | | ComplexMult(&IM(ch[ah+3*l1*ido]), &RE(ch[ah+3*l1*ido]), | 722 | | IM(d4), RE(d4), RE(wa3[i]), IM(wa3[i])); | 723 | | ComplexMult(&IM(ch[ah+4*l1*ido]), &RE(ch[ah+4*l1*ido]), | 724 | | IM(d5), RE(d5), RE(wa4[i]), IM(wa4[i])); | 725 | | #endif | 726 | | } | 727 | | } | 728 | | } | 729 | | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION | 730 | | */ | 731 | 0 | #endif | 732 | 0 | } | 733 | 61.8k | } |
Line | Count | Source | 506 | 125k | { | 507 | 125k | static real_t tr11 = FRAC_CONST(0.309016994374947); | 508 | 125k | static real_t ti11 = FRAC_CONST(0.951056516295154); | 509 | 125k | static real_t tr12 = FRAC_CONST(-0.809016994374947); | 510 | 125k | static real_t ti12 = FRAC_CONST(0.587785252292473); | 511 | 125k | uint16_t k, ac, ah; | 512 | 125k | complex_t c2, c3, c4, c5, t2, t3, t4, t5; | 513 | 125k | (void)wa1; | 514 | 125k | (void)wa2; | 515 | 125k | (void)wa3; | 516 | 125k | (void)wa4; | 517 | | | 518 | 125k | if (ido == 1) | 519 | 125k | { | 520 | 125k | if (isign == 1) | 521 | 120k | { | 522 | 6.24M | for (k = 0; k < l1; k++) | 523 | 6.11M | { | 524 | 6.11M | ac = 5*k + 1; | 525 | 6.11M | ah = k; | 526 | | | 527 | 6.11M | RE(t2) = RE(cc[ac]) + RE(cc[ac+3]); | 528 | 6.11M | IM(t2) = IM(cc[ac]) + IM(cc[ac+3]); | 529 | 6.11M | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+2]); | 530 | 6.11M | IM(t3) = IM(cc[ac+1]) + IM(cc[ac+2]); | 531 | 6.11M | RE(t4) = RE(cc[ac+1]) - RE(cc[ac+2]); | 532 | 6.11M | IM(t4) = IM(cc[ac+1]) - IM(cc[ac+2]); | 533 | 6.11M | RE(t5) = RE(cc[ac]) - RE(cc[ac+3]); | 534 | 6.11M | IM(t5) = IM(cc[ac]) - IM(cc[ac+3]); | 535 | | | 536 | 6.11M | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2) + RE(t3); | 537 | 6.11M | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2) + IM(t3); | 538 | | | 539 | 6.11M | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); | 540 | 6.11M | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); | 541 | 6.11M | RE(c3) = RE(cc[ac-1]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); | 542 | 6.11M | IM(c3) = IM(cc[ac-1]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); | 543 | | | 544 | 6.11M | ComplexMult(&RE(c5), &RE(c4), | 545 | 6.11M | ti11, ti12, RE(t5), RE(t4)); | 546 | 6.11M | ComplexMult(&IM(c5), &IM(c4), | 547 | 6.11M | ti11, ti12, IM(t5), IM(t4)); | 548 | | | 549 | 6.11M | RE(ch[ah+l1]) = RE(c2) - IM(c5); | 550 | 6.11M | IM(ch[ah+l1]) = IM(c2) + RE(c5); | 551 | 6.11M | RE(ch[ah+2*l1]) = RE(c3) - IM(c4); | 552 | 6.11M | IM(ch[ah+2*l1]) = IM(c3) + RE(c4); | 553 | 6.11M | RE(ch[ah+3*l1]) = RE(c3) + IM(c4); | 554 | 6.11M | IM(ch[ah+3*l1]) = IM(c3) - RE(c4); | 555 | 6.11M | RE(ch[ah+4*l1]) = RE(c2) + IM(c5); | 556 | 6.11M | IM(ch[ah+4*l1]) = IM(c2) - RE(c5); | 557 | 6.11M | } | 558 | 120k | } else { | 559 | 400k | for (k = 0; k < l1; k++) | 560 | 396k | { | 561 | 396k | ac = 5*k + 1; | 562 | 396k | ah = k; | 563 | | | 564 | 396k | RE(t2) = RE(cc[ac]) + RE(cc[ac+3]); | 565 | 396k | IM(t2) = IM(cc[ac]) + IM(cc[ac+3]); | 566 | 396k | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+2]); | 567 | 396k | IM(t3) = IM(cc[ac+1]) + IM(cc[ac+2]); | 568 | 396k | RE(t4) = RE(cc[ac+1]) - RE(cc[ac+2]); | 569 | 396k | IM(t4) = IM(cc[ac+1]) - IM(cc[ac+2]); | 570 | 396k | RE(t5) = RE(cc[ac]) - RE(cc[ac+3]); | 571 | 396k | IM(t5) = IM(cc[ac]) - IM(cc[ac+3]); | 572 | | | 573 | 396k | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2) + RE(t3); | 574 | 396k | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2) + IM(t3); | 575 | | | 576 | 396k | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); | 577 | 396k | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); | 578 | 396k | RE(c3) = RE(cc[ac-1]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); | 579 | 396k | IM(c3) = IM(cc[ac-1]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); | 580 | | | 581 | 396k | ComplexMult(&RE(c4), &RE(c5), | 582 | 396k | ti12, ti11, RE(t5), RE(t4)); | 583 | 396k | ComplexMult(&IM(c4), &IM(c5), | 584 | 396k | ti12, ti11, IM(t5), IM(t4)); | 585 | | | 586 | 396k | RE(ch[ah+l1]) = RE(c2) + IM(c5); | 587 | 396k | IM(ch[ah+l1]) = IM(c2) - RE(c5); | 588 | 396k | RE(ch[ah+2*l1]) = RE(c3) + IM(c4); | 589 | 396k | IM(ch[ah+2*l1]) = IM(c3) - RE(c4); | 590 | 396k | RE(ch[ah+3*l1]) = RE(c3) - IM(c4); | 591 | 396k | IM(ch[ah+3*l1]) = IM(c3) + RE(c4); | 592 | 396k | RE(ch[ah+4*l1]) = RE(c2) - IM(c5); | 593 | 396k | IM(ch[ah+4*l1]) = IM(c2) + RE(c5); | 594 | 396k | } | 595 | 4.16k | } | 596 | 125k | } else { | 597 | 0 | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION | 598 | | // TODO: remove this code once fuzzer proves it is totally unreahable | 599 | | // 5 if the biggest factor and it never repeated for supported frame | 600 | | // lengths; consequently `ido` should always be 1. | 601 | 0 | __builtin_trap(); | 602 | | /* | 603 | | #else | 604 | | uint16_t i; | 605 | | complex_t d3, d4, d5, d2; | 606 | | #endif | 607 | | if (isign == 1) | 608 | | { | 609 | | for (k = 0; k < l1; k++) | 610 | | { | 611 | | for (i = 0; i < ido; i++) | 612 | | { | 613 | | ac = i + (k*5 + 1) * ido; | 614 | | ah = i + k * ido; | 615 | | | 616 | | RE(t2) = RE(cc[ac]) + RE(cc[ac+3*ido]); | 617 | | IM(t2) = IM(cc[ac]) + IM(cc[ac+3*ido]); | 618 | | RE(t3) = RE(cc[ac+ido]) + RE(cc[ac+2*ido]); | 619 | | IM(t3) = IM(cc[ac+ido]) + IM(cc[ac+2*ido]); | 620 | | RE(t4) = RE(cc[ac+ido]) - RE(cc[ac+2*ido]); | 621 | | IM(t4) = IM(cc[ac+ido]) - IM(cc[ac+2*ido]); | 622 | | RE(t5) = RE(cc[ac]) - RE(cc[ac+3*ido]); | 623 | | IM(t5) = IM(cc[ac]) - IM(cc[ac+3*ido]); | 624 | | | 625 | | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2) + RE(t3); | 626 | | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2) + IM(t3); | 627 | | | 628 | | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); | 629 | | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); | 630 | | RE(c3) = RE(cc[ac-ido]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); | 631 | | IM(c3) = IM(cc[ac-ido]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); | 632 | | | 633 | | ComplexMult(&RE(c5), &RE(c4), | 634 | | ti11, ti12, RE(t5), RE(t4)); | 635 | | ComplexMult(&IM(c5), &IM(c4), | 636 | | ti11, ti12, IM(t5), IM(t4)); | 637 | | | 638 | | IM(d2) = IM(c2) + RE(c5); | 639 | | IM(d3) = IM(c3) + RE(c4); | 640 | | RE(d4) = RE(c3) + IM(c4); | 641 | | RE(d5) = RE(c2) + IM(c5); | 642 | | RE(d2) = RE(c2) - IM(c5); | 643 | | IM(d5) = IM(c2) - RE(c5); | 644 | | RE(d3) = RE(c3) - IM(c4); | 645 | | IM(d4) = IM(c3) - RE(c4); | 646 | | | 647 | | #if 1 | 648 | | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), | 649 | | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); | 650 | | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), | 651 | | IM(d3), RE(d3), RE(wa2[i]), IM(wa2[i])); | 652 | | ComplexMult(&IM(ch[ah+3*l1*ido]), &RE(ch[ah+3*l1*ido]), | 653 | | IM(d4), RE(d4), RE(wa3[i]), IM(wa3[i])); | 654 | | ComplexMult(&IM(ch[ah+4*l1*ido]), &RE(ch[ah+4*l1*ido]), | 655 | | IM(d5), RE(d5), RE(wa4[i]), IM(wa4[i])); | 656 | | #else | 657 | | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), | 658 | | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); | 659 | | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), | 660 | | RE(d3), IM(d3), RE(wa2[i]), IM(wa2[i])); | 661 | | ComplexMult(&RE(ch[ah+3*l1*ido]), &IM(ch[ah+3*l1*ido]), | 662 | | RE(d4), IM(d4), RE(wa3[i]), IM(wa3[i])); | 663 | | ComplexMult(&RE(ch[ah+4*l1*ido]), &IM(ch[ah+4*l1*ido]), | 664 | | RE(d5), IM(d5), RE(wa4[i]), IM(wa4[i])); | 665 | | #endif | 666 | | } | 667 | | } | 668 | | } else { | 669 | | for (k = 0; k < l1; k++) | 670 | | { | 671 | | for (i = 0; i < ido; i++) | 672 | | { | 673 | | ac = i + (k*5 + 1) * ido; | 674 | | ah = i + k * ido; | 675 | | | 676 | | RE(t2) = RE(cc[ac]) + RE(cc[ac+3*ido]); | 677 | | IM(t2) = IM(cc[ac]) + IM(cc[ac+3*ido]); | 678 | | RE(t3) = RE(cc[ac+ido]) + RE(cc[ac+2*ido]); | 679 | | IM(t3) = IM(cc[ac+ido]) + IM(cc[ac+2*ido]); | 680 | | RE(t4) = RE(cc[ac+ido]) - RE(cc[ac+2*ido]); | 681 | | IM(t4) = IM(cc[ac+ido]) - IM(cc[ac+2*ido]); | 682 | | RE(t5) = RE(cc[ac]) - RE(cc[ac+3*ido]); | 683 | | IM(t5) = IM(cc[ac]) - IM(cc[ac+3*ido]); | 684 | | | 685 | | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2) + RE(t3); | 686 | | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2) + IM(t3); | 687 | | | 688 | | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); | 689 | | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); | 690 | | RE(c3) = RE(cc[ac-ido]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); | 691 | | IM(c3) = IM(cc[ac-ido]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); | 692 | | | 693 | | ComplexMult(&RE(c4), &RE(c5), | 694 | | ti12, ti11, RE(t5), RE(t4)); | 695 | | ComplexMult(&IM(c4), &IM(c5), | 696 | | ti12, ti11, IM(t5), IM(t4)); | 697 | | | 698 | | IM(d2) = IM(c2) - RE(c5); | 699 | | IM(d3) = IM(c3) - RE(c4); | 700 | | RE(d4) = RE(c3) - IM(c4); | 701 | | RE(d5) = RE(c2) - IM(c5); | 702 | | RE(d2) = RE(c2) + IM(c5); | 703 | | IM(d5) = IM(c2) + RE(c5); | 704 | | RE(d3) = RE(c3) + IM(c4); | 705 | | IM(d4) = IM(c3) + RE(c4); | 706 | | | 707 | | #if 1 | 708 | | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), | 709 | | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); | 710 | | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), | 711 | | RE(d3), IM(d3), RE(wa2[i]), IM(wa2[i])); | 712 | | ComplexMult(&RE(ch[ah+3*l1*ido]), &IM(ch[ah+3*l1*ido]), | 713 | | RE(d4), IM(d4), RE(wa3[i]), IM(wa3[i])); | 714 | | ComplexMult(&RE(ch[ah+4*l1*ido]), &IM(ch[ah+4*l1*ido]), | 715 | | RE(d5), IM(d5), RE(wa4[i]), IM(wa4[i])); | 716 | | #else | 717 | | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), | 718 | | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); | 719 | | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), | 720 | | IM(d3), RE(d3), RE(wa2[i]), IM(wa2[i])); | 721 | | ComplexMult(&IM(ch[ah+3*l1*ido]), &RE(ch[ah+3*l1*ido]), | 722 | | IM(d4), RE(d4), RE(wa3[i]), IM(wa3[i])); | 723 | | ComplexMult(&IM(ch[ah+4*l1*ido]), &RE(ch[ah+4*l1*ido]), | 724 | | IM(d5), RE(d5), RE(wa4[i]), IM(wa4[i])); | 725 | | #endif | 726 | | } | 727 | | } | 728 | | } | 729 | | #ifdef FUZZING_BUILD_MODE_UNSAFE_FOR_PRODUCTION | 730 | | */ | 731 | 0 | #endif | 732 | 0 | } | 733 | 125k | } |
|
734 | | |
735 | | |
736 | | /*---------------------------------------------------------------------- |
737 | | cfftf1, cfftf, cfftb, cffti1, cffti. Complex FFTs. |
738 | | ----------------------------------------------------------------------*/ |
739 | | |
740 | | static INLINE void cfftf1pos(uint16_t n, complex_t *c, complex_t *ch, |
741 | | const uint16_t *ifac, const complex_t *wa, |
742 | | const int8_t isign) |
743 | 1.05M | { |
744 | 1.05M | uint16_t i; |
745 | 1.05M | uint16_t k1, l1, l2; |
746 | 1.05M | uint16_t na, nf, ip, iw, ix2, ix3, ix4, ido; |
747 | | |
748 | 1.05M | nf = ifac[1]; |
749 | 1.05M | na = 0; |
750 | 1.05M | l1 = 1; |
751 | 1.05M | iw = 0; |
752 | | |
753 | 5.35M | for (k1 = 2; k1 <= nf+1; k1++) |
754 | 4.30M | { |
755 | 4.30M | ip = ifac[k1]; |
756 | 4.30M | l2 = ip*l1; |
757 | 4.30M | ido = n / l2; |
758 | | |
759 | 4.30M | switch (ip) |
760 | 4.30M | { |
761 | 3.38M | case 4: |
762 | 3.38M | ix2 = iw + ido; |
763 | 3.38M | ix3 = ix2 + ido; |
764 | | |
765 | 3.38M | if (na == 0) |
766 | 1.84M | passf4pos((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw], &wa[ix2], &wa[ix3]); |
767 | 1.53M | else |
768 | 1.53M | passf4pos((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)ch, c, &wa[iw], &wa[ix2], &wa[ix3]); |
769 | | |
770 | 3.38M | na = 1 - na; |
771 | 3.38M | break; |
772 | 559k | case 2: |
773 | 559k | if (na == 0) |
774 | 559k | passf2pos((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw]); |
775 | 0 | else |
776 | 0 | passf2pos((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)ch, c, &wa[iw]); |
777 | | |
778 | 559k | na = 1 - na; |
779 | 559k | break; |
780 | 179k | case 3: |
781 | 179k | ix2 = iw + ido; |
782 | | |
783 | 179k | if (na == 0) |
784 | 95.6k | passf3((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw], &wa[ix2], isign); |
785 | 83.9k | else |
786 | 83.9k | passf3((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)ch, c, &wa[iw], &wa[ix2], isign); |
787 | | |
788 | 179k | na = 1 - na; |
789 | 179k | break; |
790 | 179k | case 5: |
791 | 179k | ix2 = iw + ido; |
792 | 179k | ix3 = ix2 + ido; |
793 | 179k | ix4 = ix3 + ido; |
794 | | |
795 | 179k | if (na == 0) |
796 | 177k | passf5((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw], &wa[ix2], &wa[ix3], &wa[ix4], isign); |
797 | 2.11k | else |
798 | 2.11k | passf5((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)ch, c, &wa[iw], &wa[ix2], &wa[ix3], &wa[ix4], isign); |
799 | | |
800 | 179k | na = 1 - na; |
801 | 179k | break; |
802 | 4.30M | } |
803 | | |
804 | 4.30M | l1 = l2; |
805 | 4.30M | iw += (ip-1) * ido; |
806 | 4.30M | } |
807 | | |
808 | 1.05M | if (na == 0) |
809 | 4.05k | return; |
810 | | |
811 | 316M | for (i = 0; i < n; i++) |
812 | 315M | { |
813 | 315M | RE(c[i]) = RE(ch[i]); |
814 | 315M | IM(c[i]) = IM(ch[i]); |
815 | 315M | } |
816 | 1.05M | } |
817 | | |
818 | | static INLINE void cfftf1neg(uint16_t n, complex_t *c, complex_t *ch, |
819 | | const uint16_t *ifac, const complex_t *wa, |
820 | | const int8_t isign) |
821 | 12.6k | { |
822 | 12.6k | uint16_t i; |
823 | 12.6k | uint16_t k1, l1, l2; |
824 | 12.6k | uint16_t na, nf, ip, iw, ix2, ix3, ix4, ido; |
825 | | |
826 | 12.6k | nf = ifac[1]; |
827 | 12.6k | na = 0; |
828 | 12.6k | l1 = 1; |
829 | 12.6k | iw = 0; |
830 | | |
831 | 75.3k | for (k1 = 2; k1 <= nf+1; k1++) |
832 | 62.6k | { |
833 | 62.6k | ip = ifac[k1]; |
834 | 62.6k | l2 = ip*l1; |
835 | 62.6k | ido = n / l2; |
836 | | |
837 | 62.6k | switch (ip) |
838 | 62.6k | { |
839 | 35.7k | case 4: |
840 | 35.7k | ix2 = iw + ido; |
841 | 35.7k | ix3 = ix2 + ido; |
842 | | |
843 | 35.7k | if (na == 0) |
844 | 17.8k | passf4neg((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw], &wa[ix2], &wa[ix3]); |
845 | 17.8k | else |
846 | 17.8k | passf4neg((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)ch, c, &wa[iw], &wa[ix2], &wa[ix3]); |
847 | | |
848 | 35.7k | na = 1 - na; |
849 | 35.7k | break; |
850 | 12.2k | case 2: |
851 | 12.2k | if (na == 0) |
852 | 12.2k | passf2neg((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw]); |
853 | 0 | else |
854 | 0 | passf2neg((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)ch, c, &wa[iw]); |
855 | | |
856 | 12.2k | na = 1 - na; |
857 | 12.2k | break; |
858 | 7.34k | case 3: |
859 | 7.34k | ix2 = iw + ido; |
860 | | |
861 | 7.34k | if (na == 0) |
862 | 100 | passf3((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw], &wa[ix2], isign); |
863 | 7.24k | else |
864 | 7.24k | passf3((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)ch, c, &wa[iw], &wa[ix2], isign); |
865 | | |
866 | 7.34k | na = 1 - na; |
867 | 7.34k | break; |
868 | 7.34k | case 5: |
869 | 7.34k | ix2 = iw + ido; |
870 | 7.34k | ix3 = ix2 + ido; |
871 | 7.34k | ix4 = ix3 + ido; |
872 | | |
873 | 7.34k | if (na == 0) |
874 | 7.24k | passf5((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw], &wa[ix2], &wa[ix3], &wa[ix4], isign); |
875 | 100 | else |
876 | 100 | passf5((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)ch, c, &wa[iw], &wa[ix2], &wa[ix3], &wa[ix4], isign); |
877 | | |
878 | 7.34k | na = 1 - na; |
879 | 7.34k | break; |
880 | 62.6k | } |
881 | | |
882 | 62.6k | l1 = l2; |
883 | 62.6k | iw += (ip-1) * ido; |
884 | 62.6k | } |
885 | | |
886 | 12.6k | if (na == 0) |
887 | 377 | return; |
888 | | |
889 | 6.04M | for (i = 0; i < n; i++) |
890 | 6.03M | { |
891 | 6.03M | RE(c[i]) = RE(ch[i]); |
892 | 6.03M | IM(c[i]) = IM(ch[i]); |
893 | 6.03M | } |
894 | 12.2k | } |
895 | | |
896 | | void cfftf(cfft_info *cfft, complex_t *c) |
897 | 12.6k | { |
898 | 12.6k | cfftf1neg(cfft->n, c, cfft->work, (const uint16_t*)cfft->ifac, (const complex_t*)cfft->tab, -1); |
899 | 12.6k | } |
900 | | |
901 | | void cfftb(cfft_info *cfft, complex_t *c) |
902 | 1.05M | { |
903 | 1.05M | cfftf1pos(cfft->n, c, cfft->work, (const uint16_t*)cfft->ifac, (const complex_t*)cfft->tab, +1); |
904 | 1.05M | } |
905 | | |
906 | | static void cffti1(uint16_t n, complex_t *wa, uint16_t *ifac) |
907 | 64.7k | { |
908 | 64.7k | static uint16_t ntryh[4] = {3, 4, 2, 5}; |
909 | | #ifndef FIXED_POINT |
910 | | real_t arg, argh, argld, fi; |
911 | | uint16_t ido, ipm; |
912 | | uint16_t i1, k1, l1, l2; |
913 | | uint16_t ld, ii, ip; |
914 | | #endif |
915 | 64.7k | uint16_t ntry = 0, i, j; |
916 | 64.7k | uint16_t ib; |
917 | 64.7k | uint16_t nf, nl, nq, nr; |
918 | | |
919 | 64.7k | nl = n; |
920 | 64.7k | nf = 0; |
921 | 64.7k | j = 0; |
922 | | |
923 | 175k | startloop: |
924 | 175k | j++; |
925 | | |
926 | 175k | if (j <= 4) |
927 | 175k | ntry = ntryh[j-1]; |
928 | 0 | else |
929 | 0 | ntry += 2; |
930 | | |
931 | 175k | do |
932 | 369k | { |
933 | 369k | nq = nl / ntry; |
934 | 369k | nr = nl - ntry*nq; |
935 | | |
936 | 369k | if (nr != 0) |
937 | 110k | goto startloop; |
938 | | |
939 | 258k | nf++; |
940 | 258k | ifac[nf+1] = ntry; |
941 | 258k | nl = nq; |
942 | | |
943 | 258k | if (ntry == 2 && nf != 1) |
944 | 24.7k | { |
945 | 118k | for (i = 2; i <= nf; i++) |
946 | 93.6k | { |
947 | 93.6k | ib = nf - i + 2; |
948 | 93.6k | ifac[ib+1] = ifac[ib]; |
949 | 93.6k | } |
950 | 24.7k | ifac[2] = 2; |
951 | 24.7k | } |
952 | 258k | } while (nl != 1); |
953 | | |
954 | 64.7k | ifac[0] = n; |
955 | 64.7k | ifac[1] = nf; |
956 | | |
957 | | #ifndef FIXED_POINT |
958 | | argh = (real_t)2.0*(real_t)M_PI / (real_t)n; |
959 | | i = 0; |
960 | | l1 = 1; |
961 | | |
962 | 207k | for (k1 = 1; k1 <= nf; k1++) |
963 | 165k | { |
964 | 165k | ip = ifac[k1+1]; |
965 | 165k | ld = 0; |
966 | 165k | l2 = l1*ip; |
967 | 165k | ido = n / l2; |
968 | 165k | ipm = ip - 1; |
969 | | |
970 | 629k | for (j = 0; j < ipm; j++) |
971 | 463k | { |
972 | 463k | i1 = i; |
973 | 463k | RE(wa[i]) = 1.0; |
974 | 463k | IM(wa[i]) = 0.0; |
975 | | ld += l1; |
976 | | fi = 0; |
977 | | argld = ld*argh; |
978 | | |
979 | 11.9M | for (ii = 0; ii < ido; ii++) |
980 | 11.4M | { |
981 | 11.4M | i++; |
982 | 11.4M | fi++; |
983 | 11.4M | arg = fi * argld; |
984 | 11.4M | RE(wa[i]) = (real_t)cos(arg); |
985 | 11.4M | #if 1 |
986 | 11.4M | IM(wa[i]) = (real_t)sin(arg); |
987 | | #else |
988 | | IM(wa[i]) = (real_t)-sin(arg); |
989 | | #endif |
990 | 11.4M | } |
991 | | |
992 | 463k | if (ip > 5) |
993 | 0 | { |
994 | 0 | RE(wa[i1]) = RE(wa[i]); |
995 | 0 | IM(wa[i1]) = IM(wa[i]); |
996 | 0 | } |
997 | 463k | } |
998 | 165k | l1 = l2; |
999 | 165k | } |
1000 | | #else // FIXED_POINT |
1001 | | (void)wa; /* whoa! does it work for FIXED_POINT? */ |
1002 | | #endif |
1003 | 64.7k | } Line | Count | Source | 907 | 23.3k | { | 908 | 23.3k | static uint16_t ntryh[4] = {3, 4, 2, 5}; | 909 | | #ifndef FIXED_POINT | 910 | | real_t arg, argh, argld, fi; | 911 | | uint16_t ido, ipm; | 912 | | uint16_t i1, k1, l1, l2; | 913 | | uint16_t ld, ii, ip; | 914 | | #endif | 915 | 23.3k | uint16_t ntry = 0, i, j; | 916 | 23.3k | uint16_t ib; | 917 | 23.3k | uint16_t nf, nl, nq, nr; | 918 | | | 919 | 23.3k | nl = n; | 920 | 23.3k | nf = 0; | 921 | 23.3k | j = 0; | 922 | | | 923 | 62.6k | startloop: | 924 | 62.6k | j++; | 925 | | | 926 | 62.6k | if (j <= 4) | 927 | 62.6k | ntry = ntryh[j-1]; | 928 | 0 | else | 929 | 0 | ntry += 2; | 930 | | | 931 | 62.6k | do | 932 | 132k | { | 933 | 132k | nq = nl / ntry; | 934 | 132k | nr = nl - ntry*nq; | 935 | | | 936 | 132k | if (nr != 0) | 937 | 39.3k | goto startloop; | 938 | | | 939 | 93.2k | nf++; | 940 | 93.2k | ifac[nf+1] = ntry; | 941 | 93.2k | nl = nq; | 942 | | | 943 | 93.2k | if (ntry == 2 && nf != 1) | 944 | 8.25k | { | 945 | 39.5k | for (i = 2; i <= nf; i++) | 946 | 31.2k | { | 947 | 31.2k | ib = nf - i + 2; | 948 | 31.2k | ifac[ib+1] = ifac[ib]; | 949 | 31.2k | } | 950 | 8.25k | ifac[2] = 2; | 951 | 8.25k | } | 952 | 93.2k | } while (nl != 1); | 953 | | | 954 | 23.3k | ifac[0] = n; | 955 | 23.3k | ifac[1] = nf; | 956 | | | 957 | | #ifndef FIXED_POINT | 958 | | argh = (real_t)2.0*(real_t)M_PI / (real_t)n; | 959 | | i = 0; | 960 | | l1 = 1; | 961 | | | 962 | | for (k1 = 1; k1 <= nf; k1++) | 963 | | { | 964 | | ip = ifac[k1+1]; | 965 | | ld = 0; | 966 | | l2 = l1*ip; | 967 | | ido = n / l2; | 968 | | ipm = ip - 1; | 969 | | | 970 | | for (j = 0; j < ipm; j++) | 971 | | { | 972 | | i1 = i; | 973 | | RE(wa[i]) = 1.0; | 974 | | IM(wa[i]) = 0.0; | 975 | | ld += l1; | 976 | | fi = 0; | 977 | | argld = ld*argh; | 978 | | | 979 | | for (ii = 0; ii < ido; ii++) | 980 | | { | 981 | | i++; | 982 | | fi++; | 983 | | arg = fi * argld; | 984 | | RE(wa[i]) = (real_t)cos(arg); | 985 | | #if 1 | 986 | | IM(wa[i]) = (real_t)sin(arg); | 987 | | #else | 988 | | IM(wa[i]) = (real_t)-sin(arg); | 989 | | #endif | 990 | | } | 991 | | | 992 | | if (ip > 5) | 993 | | { | 994 | | RE(wa[i1]) = RE(wa[i]); | 995 | | IM(wa[i1]) = IM(wa[i]); | 996 | | } | 997 | | } | 998 | | l1 = l2; | 999 | | } | 1000 | | #else // FIXED_POINT | 1001 | 23.3k | (void)wa; /* whoa! does it work for FIXED_POINT? */ | 1002 | 23.3k | #endif | 1003 | 23.3k | } |
Line | Count | Source | 907 | 41.4k | { | 908 | 41.4k | static uint16_t ntryh[4] = {3, 4, 2, 5}; | 909 | 41.4k | #ifndef FIXED_POINT | 910 | 41.4k | real_t arg, argh, argld, fi; | 911 | 41.4k | uint16_t ido, ipm; | 912 | 41.4k | uint16_t i1, k1, l1, l2; | 913 | 41.4k | uint16_t ld, ii, ip; | 914 | 41.4k | #endif | 915 | 41.4k | uint16_t ntry = 0, i, j; | 916 | 41.4k | uint16_t ib; | 917 | 41.4k | uint16_t nf, nl, nq, nr; | 918 | | | 919 | 41.4k | nl = n; | 920 | 41.4k | nf = 0; | 921 | 41.4k | j = 0; | 922 | | | 923 | 112k | startloop: | 924 | 112k | j++; | 925 | | | 926 | 112k | if (j <= 4) | 927 | 112k | ntry = ntryh[j-1]; | 928 | 0 | else | 929 | 0 | ntry += 2; | 930 | | | 931 | 112k | do | 932 | 237k | { | 933 | 237k | nq = nl / ntry; | 934 | 237k | nr = nl - ntry*nq; | 935 | | | 936 | 237k | if (nr != 0) | 937 | 71.5k | goto startloop; | 938 | | | 939 | 165k | nf++; | 940 | 165k | ifac[nf+1] = ntry; | 941 | 165k | nl = nq; | 942 | | | 943 | 165k | if (ntry == 2 && nf != 1) | 944 | 16.4k | { | 945 | 78.8k | for (i = 2; i <= nf; i++) | 946 | 62.3k | { | 947 | 62.3k | ib = nf - i + 2; | 948 | 62.3k | ifac[ib+1] = ifac[ib]; | 949 | 62.3k | } | 950 | 16.4k | ifac[2] = 2; | 951 | 16.4k | } | 952 | 165k | } while (nl != 1); | 953 | | | 954 | 41.4k | ifac[0] = n; | 955 | 41.4k | ifac[1] = nf; | 956 | | | 957 | 41.4k | #ifndef FIXED_POINT | 958 | 41.4k | argh = (real_t)2.0*(real_t)M_PI / (real_t)n; | 959 | 41.4k | i = 0; | 960 | 41.4k | l1 = 1; | 961 | | | 962 | 207k | for (k1 = 1; k1 <= nf; k1++) | 963 | 165k | { | 964 | 165k | ip = ifac[k1+1]; | 965 | 165k | ld = 0; | 966 | 165k | l2 = l1*ip; | 967 | 165k | ido = n / l2; | 968 | 165k | ipm = ip - 1; | 969 | | | 970 | 629k | for (j = 0; j < ipm; j++) | 971 | 463k | { | 972 | 463k | i1 = i; | 973 | 463k | RE(wa[i]) = 1.0; | 974 | 463k | IM(wa[i]) = 0.0; | 975 | 463k | ld += l1; | 976 | 463k | fi = 0; | 977 | 463k | argld = ld*argh; | 978 | | | 979 | 11.9M | for (ii = 0; ii < ido; ii++) | 980 | 11.4M | { | 981 | 11.4M | i++; | 982 | 11.4M | fi++; | 983 | 11.4M | arg = fi * argld; | 984 | 11.4M | RE(wa[i]) = (real_t)cos(arg); | 985 | 11.4M | #if 1 | 986 | 11.4M | IM(wa[i]) = (real_t)sin(arg); | 987 | | #else | 988 | | IM(wa[i]) = (real_t)-sin(arg); | 989 | | #endif | 990 | 11.4M | } | 991 | | | 992 | 463k | if (ip > 5) | 993 | 0 | { | 994 | 0 | RE(wa[i1]) = RE(wa[i]); | 995 | 0 | IM(wa[i1]) = IM(wa[i]); | 996 | 0 | } | 997 | 463k | } | 998 | 165k | l1 = l2; | 999 | 165k | } | 1000 | | #else // FIXED_POINT | 1001 | | (void)wa; /* whoa! does it work for FIXED_POINT? */ | 1002 | | #endif | 1003 | 41.4k | } |
|
1004 | | |
1005 | | cfft_info *cffti(uint16_t n) |
1006 | 2.89k | { |
1007 | 2.89k | cfft_info *cfft = (cfft_info*)faad_malloc(sizeof(cfft_info)); |
1008 | | |
1009 | 2.89k | cfft->n = n; |
1010 | 2.89k | cfft->work = (complex_t*)faad_malloc(n*sizeof(complex_t)); |
1011 | | |
1012 | | #ifndef FIXED_POINT |
1013 | | cfft->tab = (complex_t*)faad_malloc(n*sizeof(complex_t)); |
1014 | | |
1015 | | cffti1(n, cfft->tab, cfft->ifac); |
1016 | | #else |
1017 | 2.89k | cffti1(n, NULL, cfft->ifac); |
1018 | | |
1019 | 2.89k | switch (n) |
1020 | 2.89k | { |
1021 | 1.00k | case 64: cfft->tab = (complex_t*)cfft_tab_64; break; |
1022 | 1.00k | case 512: cfft->tab = (complex_t*)cfft_tab_512; break; |
1023 | | #ifdef LD_DEC |
1024 | | case 256: cfft->tab = (complex_t*)cfft_tab_256; break; |
1025 | | #endif |
1026 | | |
1027 | 0 | #ifdef ALLOW_SMALL_FRAMELENGTH |
1028 | 439 | case 60: cfft->tab = (complex_t*)cfft_tab_60; break; |
1029 | 439 | case 480: cfft->tab = (complex_t*)cfft_tab_480; break; |
1030 | | #ifdef LD_DEC |
1031 | | case 240: cfft->tab = (complex_t*)cfft_tab_240; break; |
1032 | | #endif |
1033 | 0 | #endif |
1034 | 0 | case 128: cfft->tab = (complex_t*)cfft_tab_128; break; |
1035 | 2.89k | } |
1036 | 2.89k | #endif |
1037 | | |
1038 | 2.89k | return cfft; |
1039 | 2.89k | } |
1040 | | |
1041 | | void cfftu(cfft_info *cfft) |
1042 | 23.3k | { |
1043 | 23.3k | if (cfft->work) faad_free(cfft->work); |
1044 | | #ifndef FIXED_POINT |
1045 | | if (cfft->tab) faad_free(cfft->tab); |
1046 | | #endif |
1047 | | |
1048 | 23.3k | if (cfft) faad_free(cfft); |
1049 | 23.3k | } |
1050 | | |