/proc/self/cwd/libfaad/cfft.c
Line | Count | Source |
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 | 445k | { |
73 | 445k | uint16_t i, k, ah, ac; |
74 | | |
75 | 445k | 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 | 445k | } else { |
98 | 890k | for (k = 0; k < l1; k++) |
99 | 445k | { |
100 | 445k | ah = k*ido; |
101 | 445k | ac = 2*k*ido; |
102 | | |
103 | 113M | for (i = 0; i < ido; i++) |
104 | 112M | { |
105 | 112M | complex_t t2; |
106 | | |
107 | 112M | RE(ch[ah+i]) = RE(cc[ac+i]) + RE(cc[ac+i+ido]); |
108 | 112M | RE(t2) = RE(cc[ac+i]) - RE(cc[ac+i+ido]); |
109 | | |
110 | 112M | IM(ch[ah+i]) = IM(cc[ac+i]) + IM(cc[ac+i+ido]); |
111 | 112M | IM(t2) = IM(cc[ac+i]) - IM(cc[ac+i+ido]); |
112 | | |
113 | 112M | #if 1 |
114 | 112M | ComplexMult(&IM(ch[ah+i+l1*ido]), &RE(ch[ah+i+l1*ido]), |
115 | 112M | 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 | 112M | } |
121 | 445k | } |
122 | 445k | } |
123 | 445k | } |
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 | 11.5k | { |
128 | 11.5k | uint16_t i, k, ah, ac; |
129 | | |
130 | 11.5k | 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 | 11.5k | } else { |
153 | 23.0k | for (k = 0; k < l1; k++) |
154 | 11.5k | { |
155 | 11.5k | ah = k*ido; |
156 | 11.5k | ac = 2*k*ido; |
157 | | |
158 | 2.86M | for (i = 0; i < ido; i++) |
159 | 2.85M | { |
160 | 2.85M | complex_t t2; |
161 | | |
162 | 2.85M | RE(ch[ah+i]) = RE(cc[ac+i]) + RE(cc[ac+i+ido]); |
163 | 2.85M | RE(t2) = RE(cc[ac+i]) - RE(cc[ac+i+ido]); |
164 | | |
165 | 2.85M | IM(ch[ah+i]) = IM(cc[ac+i]) + IM(cc[ac+i+ido]); |
166 | 2.85M | IM(t2) = IM(cc[ac+i]) - IM(cc[ac+i+ido]); |
167 | | |
168 | 2.85M | #if 1 |
169 | 2.85M | ComplexMult(&RE(ch[ah+i+l1*ido]), &IM(ch[ah+i+l1*ido]), |
170 | 2.85M | 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 | 2.85M | } |
176 | 11.5k | } |
177 | 11.5k | } |
178 | 11.5k | } |
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 | 188k | { |
185 | 188k | static real_t taur = FRAC_CONST(-0.5); |
186 | 188k | static real_t taui = FRAC_CONST(0.866025403784439); |
187 | 188k | uint16_t i, k, ac, ah; |
188 | 188k | complex_t c2, c3, d2, d3, t2; |
189 | | |
190 | 188k | 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 | 188k | } else { |
249 | 188k | if (isign == 1) |
250 | 181k | { |
251 | 441k | for (k = 0; k < l1; k++) |
252 | 260k | { |
253 | 14.9M | for (i = 0; i < ido; i++) |
254 | 14.6M | { |
255 | 14.6M | ac = i + (3*k+1)*ido; |
256 | 14.6M | ah = i + k * ido; |
257 | | |
258 | 14.6M | RE(t2) = RE(cc[ac]) + RE(cc[ac+ido]); |
259 | 14.6M | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),taur); |
260 | 14.6M | IM(t2) = IM(cc[ac]) + IM(cc[ac+ido]); |
261 | 14.6M | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),taur); |
262 | | |
263 | 14.6M | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2); |
264 | 14.6M | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2); |
265 | | |
266 | 14.6M | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+ido])), taui); |
267 | 14.6M | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+ido])), taui); |
268 | | |
269 | 14.6M | RE(d2) = RE(c2) - IM(c3); |
270 | 14.6M | IM(d3) = IM(c2) - RE(c3); |
271 | 14.6M | RE(d3) = RE(c2) + IM(c3); |
272 | 14.6M | IM(d2) = IM(c2) + RE(c3); |
273 | | |
274 | 14.6M | #if 1 |
275 | 14.6M | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), |
276 | 14.6M | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); |
277 | 14.6M | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), |
278 | 14.6M | 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 | 14.6M | } |
286 | 260k | } |
287 | 181k | } else { |
288 | 19.8k | for (k = 0; k < l1; k++) |
289 | 13.2k | { |
290 | 1.07M | for (i = 0; i < ido; i++) |
291 | 1.05M | { |
292 | 1.05M | ac = i + (3*k+1)*ido; |
293 | 1.05M | ah = i + k * ido; |
294 | | |
295 | 1.05M | RE(t2) = RE(cc[ac]) + RE(cc[ac+ido]); |
296 | 1.05M | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),taur); |
297 | 1.05M | IM(t2) = IM(cc[ac]) + IM(cc[ac+ido]); |
298 | 1.05M | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),taur); |
299 | | |
300 | 1.05M | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2); |
301 | 1.05M | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2); |
302 | | |
303 | 1.05M | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+ido])), taui); |
304 | 1.05M | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+ido])), taui); |
305 | | |
306 | 1.05M | RE(d2) = RE(c2) + IM(c3); |
307 | 1.05M | IM(d3) = IM(c2) + RE(c3); |
308 | 1.05M | RE(d3) = RE(c2) - IM(c3); |
309 | 1.05M | IM(d2) = IM(c2) - RE(c3); |
310 | | |
311 | 1.05M | #if 1 |
312 | 1.05M | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), |
313 | 1.05M | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); |
314 | 1.05M | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), |
315 | 1.05M | 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.05M | } |
323 | 13.2k | } |
324 | 6.65k | } |
325 | 188k | } |
326 | 188k | } Line | Count | Source | 184 | 73.4k | { | 185 | 73.4k | static real_t taur = FRAC_CONST(-0.5); | 186 | 73.4k | static real_t taui = FRAC_CONST(0.866025403784439); | 187 | 73.4k | uint16_t i, k, ac, ah; | 188 | 73.4k | complex_t c2, c3, d2, d3, t2; | 189 | | | 190 | 73.4k | 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 | 73.4k | } else { | 249 | 73.4k | if (isign == 1) | 250 | 70.4k | { | 251 | 169k | for (k = 0; k < l1; k++) | 252 | 98.9k | { | 253 | 5.54M | for (i = 0; i < ido; i++) | 254 | 5.44M | { | 255 | 5.44M | ac = i + (3*k+1)*ido; | 256 | 5.44M | ah = i + k * ido; | 257 | | | 258 | 5.44M | RE(t2) = RE(cc[ac]) + RE(cc[ac+ido]); | 259 | 5.44M | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),taur); | 260 | 5.44M | IM(t2) = IM(cc[ac]) + IM(cc[ac+ido]); | 261 | 5.44M | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),taur); | 262 | | | 263 | 5.44M | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2); | 264 | 5.44M | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2); | 265 | | | 266 | 5.44M | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+ido])), taui); | 267 | 5.44M | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+ido])), taui); | 268 | | | 269 | 5.44M | RE(d2) = RE(c2) - IM(c3); | 270 | 5.44M | IM(d3) = IM(c2) - RE(c3); | 271 | 5.44M | RE(d3) = RE(c2) + IM(c3); | 272 | 5.44M | IM(d2) = IM(c2) + RE(c3); | 273 | | | 274 | 5.44M | #if 1 | 275 | 5.44M | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), | 276 | 5.44M | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); | 277 | 5.44M | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), | 278 | 5.44M | 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.44M | } | 286 | 98.9k | } | 287 | 70.4k | } else { | 288 | 9.02k | for (k = 0; k < l1; k++) | 289 | 6.01k | { | 290 | 486k | for (i = 0; i < ido; i++) | 291 | 480k | { | 292 | 480k | ac = i + (3*k+1)*ido; | 293 | 480k | ah = i + k * ido; | 294 | | | 295 | 480k | RE(t2) = RE(cc[ac]) + RE(cc[ac+ido]); | 296 | 480k | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),taur); | 297 | 480k | IM(t2) = IM(cc[ac]) + IM(cc[ac+ido]); | 298 | 480k | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),taur); | 299 | | | 300 | 480k | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2); | 301 | 480k | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2); | 302 | | | 303 | 480k | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+ido])), taui); | 304 | 480k | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+ido])), taui); | 305 | | | 306 | 480k | RE(d2) = RE(c2) + IM(c3); | 307 | 480k | IM(d3) = IM(c2) + RE(c3); | 308 | 480k | RE(d3) = RE(c2) - IM(c3); | 309 | 480k | IM(d2) = IM(c2) - RE(c3); | 310 | | | 311 | 480k | #if 1 | 312 | 480k | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), | 313 | 480k | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); | 314 | 480k | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), | 315 | 480k | 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 | 480k | } | 323 | 6.01k | } | 324 | 3.01k | } | 325 | 73.4k | } | 326 | 73.4k | } |
Line | Count | Source | 184 | 115k | { | 185 | 115k | static real_t taur = FRAC_CONST(-0.5); | 186 | 115k | static real_t taui = FRAC_CONST(0.866025403784439); | 187 | 115k | uint16_t i, k, ac, ah; | 188 | 115k | complex_t c2, c3, d2, d3, t2; | 189 | | | 190 | 115k | 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 | 115k | } else { | 249 | 115k | if (isign == 1) | 250 | 111k | { | 251 | 272k | for (k = 0; k < l1; k++) | 252 | 161k | { | 253 | 9.39M | for (i = 0; i < ido; i++) | 254 | 9.22M | { | 255 | 9.22M | ac = i + (3*k+1)*ido; | 256 | 9.22M | ah = i + k * ido; | 257 | | | 258 | 9.22M | RE(t2) = RE(cc[ac]) + RE(cc[ac+ido]); | 259 | 9.22M | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),taur); | 260 | 9.22M | IM(t2) = IM(cc[ac]) + IM(cc[ac+ido]); | 261 | 9.22M | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),taur); | 262 | | | 263 | 9.22M | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2); | 264 | 9.22M | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2); | 265 | | | 266 | 9.22M | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+ido])), taui); | 267 | 9.22M | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+ido])), taui); | 268 | | | 269 | 9.22M | RE(d2) = RE(c2) - IM(c3); | 270 | 9.22M | IM(d3) = IM(c2) - RE(c3); | 271 | 9.22M | RE(d3) = RE(c2) + IM(c3); | 272 | 9.22M | IM(d2) = IM(c2) + RE(c3); | 273 | | | 274 | 9.22M | #if 1 | 275 | 9.22M | ComplexMult(&IM(ch[ah+l1*ido]), &RE(ch[ah+l1*ido]), | 276 | 9.22M | IM(d2), RE(d2), RE(wa1[i]), IM(wa1[i])); | 277 | 9.22M | ComplexMult(&IM(ch[ah+2*l1*ido]), &RE(ch[ah+2*l1*ido]), | 278 | 9.22M | 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 | 9.22M | } | 286 | 161k | } | 287 | 111k | } else { | 288 | 10.8k | for (k = 0; k < l1; k++) | 289 | 7.21k | { | 290 | 584k | for (i = 0; i < ido; i++) | 291 | 577k | { | 292 | 577k | ac = i + (3*k+1)*ido; | 293 | 577k | ah = i + k * ido; | 294 | | | 295 | 577k | RE(t2) = RE(cc[ac]) + RE(cc[ac+ido]); | 296 | 577k | RE(c2) = RE(cc[ac-ido]) + MUL_F(RE(t2),taur); | 297 | 577k | IM(t2) = IM(cc[ac]) + IM(cc[ac+ido]); | 298 | 577k | IM(c2) = IM(cc[ac-ido]) + MUL_F(IM(t2),taur); | 299 | | | 300 | 577k | RE(ch[ah]) = RE(cc[ac-ido]) + RE(t2); | 301 | 577k | IM(ch[ah]) = IM(cc[ac-ido]) + IM(t2); | 302 | | | 303 | 577k | RE(c3) = MUL_F((RE(cc[ac]) - RE(cc[ac+ido])), taui); | 304 | 577k | IM(c3) = MUL_F((IM(cc[ac]) - IM(cc[ac+ido])), taui); | 305 | | | 306 | 577k | RE(d2) = RE(c2) + IM(c3); | 307 | 577k | IM(d3) = IM(c2) + RE(c3); | 308 | 577k | RE(d3) = RE(c2) - IM(c3); | 309 | 577k | IM(d2) = IM(c2) - RE(c3); | 310 | | | 311 | 577k | #if 1 | 312 | 577k | ComplexMult(&RE(ch[ah+l1*ido]), &IM(ch[ah+l1*ido]), | 313 | 577k | RE(d2), IM(d2), RE(wa1[i]), IM(wa1[i])); | 314 | 577k | ComplexMult(&RE(ch[ah+2*l1*ido]), &IM(ch[ah+2*l1*ido]), | 315 | 577k | 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 | 577k | } | 323 | 7.21k | } | 324 | 3.63k | } | 325 | 115k | } | 326 | 115k | } |
|
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 | 2.90M | { |
333 | 2.90M | uint16_t i, k, ac, ah; |
334 | | |
335 | 2.90M | if (ido == 1) |
336 | 757k | { |
337 | 54.0M | for (k = 0; k < l1; k++) |
338 | 53.2M | { |
339 | 53.2M | complex_t t1, t2, t3, t4; |
340 | | |
341 | 53.2M | ac = 4*k; |
342 | 53.2M | ah = k; |
343 | | |
344 | 53.2M | RE(t2) = RE(cc[ac]) + RE(cc[ac+2]); |
345 | 53.2M | RE(t1) = RE(cc[ac]) - RE(cc[ac+2]); |
346 | 53.2M | IM(t2) = IM(cc[ac]) + IM(cc[ac+2]); |
347 | 53.2M | IM(t1) = IM(cc[ac]) - IM(cc[ac+2]); |
348 | 53.2M | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+3]); |
349 | 53.2M | IM(t4) = RE(cc[ac+1]) - RE(cc[ac+3]); |
350 | 53.2M | IM(t3) = IM(cc[ac+3]) + IM(cc[ac+1]); |
351 | 53.2M | RE(t4) = IM(cc[ac+3]) - IM(cc[ac+1]); |
352 | | |
353 | 53.2M | RE(ch[ah]) = RE(t2) + RE(t3); |
354 | 53.2M | RE(ch[ah+2*l1]) = RE(t2) - RE(t3); |
355 | | |
356 | 53.2M | IM(ch[ah]) = IM(t2) + IM(t3); |
357 | 53.2M | IM(ch[ah+2*l1]) = IM(t2) - IM(t3); |
358 | | |
359 | 53.2M | RE(ch[ah+l1]) = RE(t1) + RE(t4); |
360 | 53.2M | RE(ch[ah+3*l1]) = RE(t1) - RE(t4); |
361 | | |
362 | 53.2M | IM(ch[ah+l1]) = IM(t1) + IM(t4); |
363 | 53.2M | IM(ch[ah+3*l1]) = IM(t1) - IM(t4); |
364 | 53.2M | } |
365 | 2.14M | } else { |
366 | 22.2M | for (k = 0; k < l1; k++) |
367 | 20.0M | { |
368 | 20.0M | ac = 4*k*ido; |
369 | 20.0M | ah = k*ido; |
370 | | |
371 | 194M | for (i = 0; i < ido; i++) |
372 | 174M | { |
373 | 174M | complex_t c2, c3, c4, t1, t2, t3, t4; |
374 | | |
375 | 174M | RE(t2) = RE(cc[ac+i]) + RE(cc[ac+i+2*ido]); |
376 | 174M | RE(t1) = RE(cc[ac+i]) - RE(cc[ac+i+2*ido]); |
377 | 174M | IM(t2) = IM(cc[ac+i]) + IM(cc[ac+i+2*ido]); |
378 | 174M | IM(t1) = IM(cc[ac+i]) - IM(cc[ac+i+2*ido]); |
379 | 174M | RE(t3) = RE(cc[ac+i+ido]) + RE(cc[ac+i+3*ido]); |
380 | 174M | IM(t4) = RE(cc[ac+i+ido]) - RE(cc[ac+i+3*ido]); |
381 | 174M | IM(t3) = IM(cc[ac+i+3*ido]) + IM(cc[ac+i+ido]); |
382 | 174M | RE(t4) = IM(cc[ac+i+3*ido]) - IM(cc[ac+i+ido]); |
383 | | |
384 | 174M | RE(c2) = RE(t1) + RE(t4); |
385 | 174M | RE(c4) = RE(t1) - RE(t4); |
386 | | |
387 | 174M | IM(c2) = IM(t1) + IM(t4); |
388 | 174M | IM(c4) = IM(t1) - IM(t4); |
389 | | |
390 | 174M | RE(ch[ah+i]) = RE(t2) + RE(t3); |
391 | 174M | RE(c3) = RE(t2) - RE(t3); |
392 | | |
393 | 174M | IM(ch[ah+i]) = IM(t2) + IM(t3); |
394 | 174M | IM(c3) = IM(t2) - IM(t3); |
395 | | |
396 | 174M | #if 1 |
397 | 174M | ComplexMult(&IM(ch[ah+i+l1*ido]), &RE(ch[ah+i+l1*ido]), |
398 | 174M | IM(c2), RE(c2), RE(wa1[i]), IM(wa1[i])); |
399 | 174M | ComplexMult(&IM(ch[ah+i+2*l1*ido]), &RE(ch[ah+i+2*l1*ido]), |
400 | 174M | IM(c3), RE(c3), RE(wa2[i]), IM(wa2[i])); |
401 | 174M | ComplexMult(&IM(ch[ah+i+3*l1*ido]), &RE(ch[ah+i+3*l1*ido]), |
402 | 174M | 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 | 174M | } |
412 | 20.0M | } |
413 | 2.14M | } |
414 | 2.90M | } |
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 | 34.0k | { |
420 | 34.0k | uint16_t i, k, ac, ah; |
421 | | |
422 | 34.0k | if (ido == 1) |
423 | 5.18k | { |
424 | 655k | for (k = 0; k < l1; k++) |
425 | 650k | { |
426 | 650k | complex_t t1, t2, t3, t4; |
427 | | |
428 | 650k | ac = 4*k; |
429 | 650k | ah = k; |
430 | | |
431 | 650k | RE(t2) = RE(cc[ac]) + RE(cc[ac+2]); |
432 | 650k | RE(t1) = RE(cc[ac]) - RE(cc[ac+2]); |
433 | 650k | IM(t2) = IM(cc[ac]) + IM(cc[ac+2]); |
434 | 650k | IM(t1) = IM(cc[ac]) - IM(cc[ac+2]); |
435 | 650k | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+3]); |
436 | 650k | IM(t4) = RE(cc[ac+1]) - RE(cc[ac+3]); |
437 | 650k | IM(t3) = IM(cc[ac+3]) + IM(cc[ac+1]); |
438 | 650k | RE(t4) = IM(cc[ac+3]) - IM(cc[ac+1]); |
439 | | |
440 | 650k | RE(ch[ah]) = RE(t2) + RE(t3); |
441 | 650k | RE(ch[ah+2*l1]) = RE(t2) - RE(t3); |
442 | | |
443 | 650k | IM(ch[ah]) = IM(t2) + IM(t3); |
444 | 650k | IM(ch[ah+2*l1]) = IM(t2) - IM(t3); |
445 | | |
446 | 650k | RE(ch[ah+l1]) = RE(t1) - RE(t4); |
447 | 650k | RE(ch[ah+3*l1]) = RE(t1) + RE(t4); |
448 | | |
449 | 650k | IM(ch[ah+l1]) = IM(t1) - IM(t4); |
450 | 650k | IM(ch[ah+3*l1]) = IM(t1) + IM(t4); |
451 | 650k | } |
452 | 28.8k | } else { |
453 | 440k | for (k = 0; k < l1; k++) |
454 | 411k | { |
455 | 411k | ac = 4*k*ido; |
456 | 411k | ah = k*ido; |
457 | | |
458 | 3.95M | for (i = 0; i < ido; i++) |
459 | 3.53M | { |
460 | 3.53M | complex_t c2, c3, c4, t1, t2, t3, t4; |
461 | | |
462 | 3.53M | RE(t2) = RE(cc[ac+i]) + RE(cc[ac+i+2*ido]); |
463 | 3.53M | RE(t1) = RE(cc[ac+i]) - RE(cc[ac+i+2*ido]); |
464 | 3.53M | IM(t2) = IM(cc[ac+i]) + IM(cc[ac+i+2*ido]); |
465 | 3.53M | IM(t1) = IM(cc[ac+i]) - IM(cc[ac+i+2*ido]); |
466 | 3.53M | RE(t3) = RE(cc[ac+i+ido]) + RE(cc[ac+i+3*ido]); |
467 | 3.53M | IM(t4) = RE(cc[ac+i+ido]) - RE(cc[ac+i+3*ido]); |
468 | 3.53M | IM(t3) = IM(cc[ac+i+3*ido]) + IM(cc[ac+i+ido]); |
469 | 3.53M | RE(t4) = IM(cc[ac+i+3*ido]) - IM(cc[ac+i+ido]); |
470 | | |
471 | 3.53M | RE(c2) = RE(t1) - RE(t4); |
472 | 3.53M | RE(c4) = RE(t1) + RE(t4); |
473 | | |
474 | 3.53M | IM(c2) = IM(t1) - IM(t4); |
475 | 3.53M | IM(c4) = IM(t1) + IM(t4); |
476 | | |
477 | 3.53M | RE(ch[ah+i]) = RE(t2) + RE(t3); |
478 | 3.53M | RE(c3) = RE(t2) - RE(t3); |
479 | | |
480 | 3.53M | IM(ch[ah+i]) = IM(t2) + IM(t3); |
481 | 3.53M | IM(c3) = IM(t2) - IM(t3); |
482 | | |
483 | 3.53M | #if 1 |
484 | 3.53M | ComplexMult(&RE(ch[ah+i+l1*ido]), &IM(ch[ah+i+l1*ido]), |
485 | 3.53M | RE(c2), IM(c2), RE(wa1[i]), IM(wa1[i])); |
486 | 3.53M | ComplexMult(&RE(ch[ah+i+2*l1*ido]), &IM(ch[ah+i+2*l1*ido]), |
487 | 3.53M | RE(c3), IM(c3), RE(wa2[i]), IM(wa2[i])); |
488 | 3.53M | ComplexMult(&RE(ch[ah+i+3*l1*ido]), &IM(ch[ah+i+3*l1*ido]), |
489 | 3.53M | 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.53M | } |
499 | 411k | } |
500 | 28.8k | } |
501 | 34.0k | } |
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 | 188k | { |
507 | 188k | static real_t tr11 = FRAC_CONST(0.309016994374947); |
508 | 188k | static real_t ti11 = FRAC_CONST(0.951056516295154); |
509 | 188k | static real_t tr12 = FRAC_CONST(-0.809016994374947); |
510 | 188k | static real_t ti12 = FRAC_CONST(0.587785252292473); |
511 | 188k | uint16_t k, ac, ah; |
512 | 188k | complex_t c2, c3, c4, c5, t2, t3, t4, t5; |
513 | 188k | (void)wa1; |
514 | 188k | (void)wa2; |
515 | 188k | (void)wa3; |
516 | 188k | (void)wa4; |
517 | | |
518 | 188k | if (ido == 1) |
519 | 188k | { |
520 | 188k | if (isign == 1) |
521 | 181k | { |
522 | 8.98M | for (k = 0; k < l1; k++) |
523 | 8.80M | { |
524 | 8.80M | ac = 5*k + 1; |
525 | 8.80M | ah = k; |
526 | | |
527 | 8.80M | RE(t2) = RE(cc[ac]) + RE(cc[ac+3]); |
528 | 8.80M | IM(t2) = IM(cc[ac]) + IM(cc[ac+3]); |
529 | 8.80M | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+2]); |
530 | 8.80M | IM(t3) = IM(cc[ac+1]) + IM(cc[ac+2]); |
531 | 8.80M | RE(t4) = RE(cc[ac+1]) - RE(cc[ac+2]); |
532 | 8.80M | IM(t4) = IM(cc[ac+1]) - IM(cc[ac+2]); |
533 | 8.80M | RE(t5) = RE(cc[ac]) - RE(cc[ac+3]); |
534 | 8.80M | IM(t5) = IM(cc[ac]) - IM(cc[ac+3]); |
535 | | |
536 | 8.80M | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2) + RE(t3); |
537 | 8.80M | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2) + IM(t3); |
538 | | |
539 | 8.80M | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); |
540 | 8.80M | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); |
541 | 8.80M | RE(c3) = RE(cc[ac-1]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); |
542 | 8.80M | IM(c3) = IM(cc[ac-1]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); |
543 | | |
544 | 8.80M | ComplexMult(&RE(c5), &RE(c4), |
545 | 8.80M | ti11, ti12, RE(t5), RE(t4)); |
546 | 8.80M | ComplexMult(&IM(c5), &IM(c4), |
547 | 8.80M | ti11, ti12, IM(t5), IM(t4)); |
548 | | |
549 | 8.80M | RE(ch[ah+l1]) = RE(c2) - IM(c5); |
550 | 8.80M | IM(ch[ah+l1]) = IM(c2) + RE(c5); |
551 | 8.80M | RE(ch[ah+2*l1]) = RE(c3) - IM(c4); |
552 | 8.80M | IM(ch[ah+2*l1]) = IM(c3) + RE(c4); |
553 | 8.80M | RE(ch[ah+3*l1]) = RE(c3) + IM(c4); |
554 | 8.80M | IM(ch[ah+3*l1]) = IM(c3) - RE(c4); |
555 | 8.80M | RE(ch[ah+4*l1]) = RE(c2) + IM(c5); |
556 | 8.80M | IM(ch[ah+4*l1]) = IM(c2) - RE(c5); |
557 | 8.80M | } |
558 | 181k | } else { |
559 | 641k | for (k = 0; k < l1; k++) |
560 | 634k | { |
561 | 634k | ac = 5*k + 1; |
562 | 634k | ah = k; |
563 | | |
564 | 634k | RE(t2) = RE(cc[ac]) + RE(cc[ac+3]); |
565 | 634k | IM(t2) = IM(cc[ac]) + IM(cc[ac+3]); |
566 | 634k | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+2]); |
567 | 634k | IM(t3) = IM(cc[ac+1]) + IM(cc[ac+2]); |
568 | 634k | RE(t4) = RE(cc[ac+1]) - RE(cc[ac+2]); |
569 | 634k | IM(t4) = IM(cc[ac+1]) - IM(cc[ac+2]); |
570 | 634k | RE(t5) = RE(cc[ac]) - RE(cc[ac+3]); |
571 | 634k | IM(t5) = IM(cc[ac]) - IM(cc[ac+3]); |
572 | | |
573 | 634k | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2) + RE(t3); |
574 | 634k | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2) + IM(t3); |
575 | | |
576 | 634k | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); |
577 | 634k | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); |
578 | 634k | RE(c3) = RE(cc[ac-1]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); |
579 | 634k | IM(c3) = IM(cc[ac-1]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); |
580 | | |
581 | 634k | ComplexMult(&RE(c4), &RE(c5), |
582 | 634k | ti12, ti11, RE(t5), RE(t4)); |
583 | 634k | ComplexMult(&IM(c4), &IM(c5), |
584 | 634k | ti12, ti11, IM(t5), IM(t4)); |
585 | | |
586 | 634k | RE(ch[ah+l1]) = RE(c2) + IM(c5); |
587 | 634k | IM(ch[ah+l1]) = IM(c2) - RE(c5); |
588 | 634k | RE(ch[ah+2*l1]) = RE(c3) + IM(c4); |
589 | 634k | IM(ch[ah+2*l1]) = IM(c3) - RE(c4); |
590 | 634k | RE(ch[ah+3*l1]) = RE(c3) - IM(c4); |
591 | 634k | IM(ch[ah+3*l1]) = IM(c3) + RE(c4); |
592 | 634k | RE(ch[ah+4*l1]) = RE(c2) - IM(c5); |
593 | 634k | IM(ch[ah+4*l1]) = IM(c2) + RE(c5); |
594 | 634k | } |
595 | 6.65k | } |
596 | 188k | } 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 | 188k | } Line | Count | Source | 506 | 73.4k | { | 507 | 73.4k | static real_t tr11 = FRAC_CONST(0.309016994374947); | 508 | 73.4k | static real_t ti11 = FRAC_CONST(0.951056516295154); | 509 | 73.4k | static real_t tr12 = FRAC_CONST(-0.809016994374947); | 510 | 73.4k | static real_t ti12 = FRAC_CONST(0.587785252292473); | 511 | 73.4k | uint16_t k, ac, ah; | 512 | 73.4k | complex_t c2, c3, c4, c5, t2, t3, t4, t5; | 513 | 73.4k | (void)wa1; | 514 | 73.4k | (void)wa2; | 515 | 73.4k | (void)wa3; | 516 | 73.4k | (void)wa4; | 517 | | | 518 | 73.4k | if (ido == 1) | 519 | 73.4k | { | 520 | 73.4k | if (isign == 1) | 521 | 70.4k | { | 522 | 3.33M | for (k = 0; k < l1; k++) | 523 | 3.26M | { | 524 | 3.26M | ac = 5*k + 1; | 525 | 3.26M | ah = k; | 526 | | | 527 | 3.26M | RE(t2) = RE(cc[ac]) + RE(cc[ac+3]); | 528 | 3.26M | IM(t2) = IM(cc[ac]) + IM(cc[ac+3]); | 529 | 3.26M | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+2]); | 530 | 3.26M | IM(t3) = IM(cc[ac+1]) + IM(cc[ac+2]); | 531 | 3.26M | RE(t4) = RE(cc[ac+1]) - RE(cc[ac+2]); | 532 | 3.26M | IM(t4) = IM(cc[ac+1]) - IM(cc[ac+2]); | 533 | 3.26M | RE(t5) = RE(cc[ac]) - RE(cc[ac+3]); | 534 | 3.26M | IM(t5) = IM(cc[ac]) - IM(cc[ac+3]); | 535 | | | 536 | 3.26M | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2) + RE(t3); | 537 | 3.26M | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2) + IM(t3); | 538 | | | 539 | 3.26M | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); | 540 | 3.26M | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); | 541 | 3.26M | RE(c3) = RE(cc[ac-1]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); | 542 | 3.26M | IM(c3) = IM(cc[ac-1]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); | 543 | | | 544 | 3.26M | ComplexMult(&RE(c5), &RE(c4), | 545 | 3.26M | ti11, ti12, RE(t5), RE(t4)); | 546 | 3.26M | ComplexMult(&IM(c5), &IM(c4), | 547 | 3.26M | ti11, ti12, IM(t5), IM(t4)); | 548 | | | 549 | 3.26M | RE(ch[ah+l1]) = RE(c2) - IM(c5); | 550 | 3.26M | IM(ch[ah+l1]) = IM(c2) + RE(c5); | 551 | 3.26M | RE(ch[ah+2*l1]) = RE(c3) - IM(c4); | 552 | 3.26M | IM(ch[ah+2*l1]) = IM(c3) + RE(c4); | 553 | 3.26M | RE(ch[ah+3*l1]) = RE(c3) + IM(c4); | 554 | 3.26M | IM(ch[ah+3*l1]) = IM(c3) - RE(c4); | 555 | 3.26M | RE(ch[ah+4*l1]) = RE(c2) + IM(c5); | 556 | 3.26M | IM(ch[ah+4*l1]) = IM(c2) - RE(c5); | 557 | 3.26M | } | 558 | 70.4k | } else { | 559 | 291k | for (k = 0; k < l1; k++) | 560 | 288k | { | 561 | 288k | ac = 5*k + 1; | 562 | 288k | ah = k; | 563 | | | 564 | 288k | RE(t2) = RE(cc[ac]) + RE(cc[ac+3]); | 565 | 288k | IM(t2) = IM(cc[ac]) + IM(cc[ac+3]); | 566 | 288k | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+2]); | 567 | 288k | IM(t3) = IM(cc[ac+1]) + IM(cc[ac+2]); | 568 | 288k | RE(t4) = RE(cc[ac+1]) - RE(cc[ac+2]); | 569 | 288k | IM(t4) = IM(cc[ac+1]) - IM(cc[ac+2]); | 570 | 288k | RE(t5) = RE(cc[ac]) - RE(cc[ac+3]); | 571 | 288k | IM(t5) = IM(cc[ac]) - IM(cc[ac+3]); | 572 | | | 573 | 288k | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2) + RE(t3); | 574 | 288k | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2) + IM(t3); | 575 | | | 576 | 288k | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); | 577 | 288k | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); | 578 | 288k | RE(c3) = RE(cc[ac-1]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); | 579 | 288k | IM(c3) = IM(cc[ac-1]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); | 580 | | | 581 | 288k | ComplexMult(&RE(c4), &RE(c5), | 582 | 288k | ti12, ti11, RE(t5), RE(t4)); | 583 | 288k | ComplexMult(&IM(c4), &IM(c5), | 584 | 288k | ti12, ti11, IM(t5), IM(t4)); | 585 | | | 586 | 288k | RE(ch[ah+l1]) = RE(c2) + IM(c5); | 587 | 288k | IM(ch[ah+l1]) = IM(c2) - RE(c5); | 588 | 288k | RE(ch[ah+2*l1]) = RE(c3) + IM(c4); | 589 | 288k | IM(ch[ah+2*l1]) = IM(c3) - RE(c4); | 590 | 288k | RE(ch[ah+3*l1]) = RE(c3) - IM(c4); | 591 | 288k | IM(ch[ah+3*l1]) = IM(c3) + RE(c4); | 592 | 288k | RE(ch[ah+4*l1]) = RE(c2) - IM(c5); | 593 | 288k | IM(ch[ah+4*l1]) = IM(c2) + RE(c5); | 594 | 288k | } | 595 | 3.01k | } | 596 | 73.4k | } 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 | 73.4k | } |
Line | Count | Source | 506 | 115k | { | 507 | 115k | static real_t tr11 = FRAC_CONST(0.309016994374947); | 508 | 115k | static real_t ti11 = FRAC_CONST(0.951056516295154); | 509 | 115k | static real_t tr12 = FRAC_CONST(-0.809016994374947); | 510 | 115k | static real_t ti12 = FRAC_CONST(0.587785252292473); | 511 | 115k | uint16_t k, ac, ah; | 512 | 115k | complex_t c2, c3, c4, c5, t2, t3, t4, t5; | 513 | 115k | (void)wa1; | 514 | 115k | (void)wa2; | 515 | 115k | (void)wa3; | 516 | 115k | (void)wa4; | 517 | | | 518 | 115k | if (ido == 1) | 519 | 115k | { | 520 | 115k | if (isign == 1) | 521 | 111k | { | 522 | 5.64M | for (k = 0; k < l1; k++) | 523 | 5.53M | { | 524 | 5.53M | ac = 5*k + 1; | 525 | 5.53M | ah = k; | 526 | | | 527 | 5.53M | RE(t2) = RE(cc[ac]) + RE(cc[ac+3]); | 528 | 5.53M | IM(t2) = IM(cc[ac]) + IM(cc[ac+3]); | 529 | 5.53M | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+2]); | 530 | 5.53M | IM(t3) = IM(cc[ac+1]) + IM(cc[ac+2]); | 531 | 5.53M | RE(t4) = RE(cc[ac+1]) - RE(cc[ac+2]); | 532 | 5.53M | IM(t4) = IM(cc[ac+1]) - IM(cc[ac+2]); | 533 | 5.53M | RE(t5) = RE(cc[ac]) - RE(cc[ac+3]); | 534 | 5.53M | IM(t5) = IM(cc[ac]) - IM(cc[ac+3]); | 535 | | | 536 | 5.53M | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2) + RE(t3); | 537 | 5.53M | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2) + IM(t3); | 538 | | | 539 | 5.53M | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); | 540 | 5.53M | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); | 541 | 5.53M | RE(c3) = RE(cc[ac-1]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); | 542 | 5.53M | IM(c3) = IM(cc[ac-1]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); | 543 | | | 544 | 5.53M | ComplexMult(&RE(c5), &RE(c4), | 545 | 5.53M | ti11, ti12, RE(t5), RE(t4)); | 546 | 5.53M | ComplexMult(&IM(c5), &IM(c4), | 547 | 5.53M | ti11, ti12, IM(t5), IM(t4)); | 548 | | | 549 | 5.53M | RE(ch[ah+l1]) = RE(c2) - IM(c5); | 550 | 5.53M | IM(ch[ah+l1]) = IM(c2) + RE(c5); | 551 | 5.53M | RE(ch[ah+2*l1]) = RE(c3) - IM(c4); | 552 | 5.53M | IM(ch[ah+2*l1]) = IM(c3) + RE(c4); | 553 | 5.53M | RE(ch[ah+3*l1]) = RE(c3) + IM(c4); | 554 | 5.53M | IM(ch[ah+3*l1]) = IM(c3) - RE(c4); | 555 | 5.53M | RE(ch[ah+4*l1]) = RE(c2) + IM(c5); | 556 | 5.53M | IM(ch[ah+4*l1]) = IM(c2) - RE(c5); | 557 | 5.53M | } | 558 | 111k | } else { | 559 | 349k | for (k = 0; k < l1; k++) | 560 | 346k | { | 561 | 346k | ac = 5*k + 1; | 562 | 346k | ah = k; | 563 | | | 564 | 346k | RE(t2) = RE(cc[ac]) + RE(cc[ac+3]); | 565 | 346k | IM(t2) = IM(cc[ac]) + IM(cc[ac+3]); | 566 | 346k | RE(t3) = RE(cc[ac+1]) + RE(cc[ac+2]); | 567 | 346k | IM(t3) = IM(cc[ac+1]) + IM(cc[ac+2]); | 568 | 346k | RE(t4) = RE(cc[ac+1]) - RE(cc[ac+2]); | 569 | 346k | IM(t4) = IM(cc[ac+1]) - IM(cc[ac+2]); | 570 | 346k | RE(t5) = RE(cc[ac]) - RE(cc[ac+3]); | 571 | 346k | IM(t5) = IM(cc[ac]) - IM(cc[ac+3]); | 572 | | | 573 | 346k | RE(ch[ah]) = RE(cc[ac-1]) + RE(t2) + RE(t3); | 574 | 346k | IM(ch[ah]) = IM(cc[ac-1]) + IM(t2) + IM(t3); | 575 | | | 576 | 346k | RE(c2) = RE(cc[ac-1]) + MUL_F(RE(t2),tr11) + MUL_F(RE(t3),tr12); | 577 | 346k | IM(c2) = IM(cc[ac-1]) + MUL_F(IM(t2),tr11) + MUL_F(IM(t3),tr12); | 578 | 346k | RE(c3) = RE(cc[ac-1]) + MUL_F(RE(t2),tr12) + MUL_F(RE(t3),tr11); | 579 | 346k | IM(c3) = IM(cc[ac-1]) + MUL_F(IM(t2),tr12) + MUL_F(IM(t3),tr11); | 580 | | | 581 | 346k | ComplexMult(&RE(c4), &RE(c5), | 582 | 346k | ti12, ti11, RE(t5), RE(t4)); | 583 | 346k | ComplexMult(&IM(c4), &IM(c5), | 584 | 346k | ti12, ti11, IM(t5), IM(t4)); | 585 | | | 586 | 346k | RE(ch[ah+l1]) = RE(c2) + IM(c5); | 587 | 346k | IM(ch[ah+l1]) = IM(c2) - RE(c5); | 588 | 346k | RE(ch[ah+2*l1]) = RE(c3) + IM(c4); | 589 | 346k | IM(ch[ah+2*l1]) = IM(c3) - RE(c4); | 590 | 346k | RE(ch[ah+3*l1]) = RE(c3) - IM(c4); | 591 | 346k | IM(ch[ah+3*l1]) = IM(c3) + RE(c4); | 592 | 346k | RE(ch[ah+4*l1]) = RE(c2) - IM(c5); | 593 | 346k | IM(ch[ah+4*l1]) = IM(c2) + RE(c5); | 594 | 346k | } | 595 | 3.63k | } | 596 | 115k | } 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 | 115k | } |
|
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 | 939k | { |
744 | 939k | uint16_t i; |
745 | 939k | uint16_t k1, l1, l2; |
746 | 939k | uint16_t na, nf, ip, iw, ix2, ix3, ix4, ido; |
747 | | |
748 | 939k | nf = ifac[1]; |
749 | 939k | na = 0; |
750 | 939k | l1 = 1; |
751 | 939k | iw = 0; |
752 | | |
753 | 4.64M | for (k1 = 2; k1 <= nf+1; k1++) |
754 | 3.71M | { |
755 | 3.71M | ip = ifac[k1]; |
756 | 3.71M | l2 = ip*l1; |
757 | 3.71M | ido = n / l2; |
758 | | |
759 | 3.71M | switch (ip) |
760 | 3.71M | { |
761 | 2.90M | case 4: |
762 | 2.90M | ix2 = iw + ido; |
763 | 2.90M | ix3 = ix2 + ido; |
764 | | |
765 | 2.90M | if (na == 0) |
766 | 1.59M | passf4pos((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw], &wa[ix2], &wa[ix3]); |
767 | 1.30M | else |
768 | 1.30M | passf4pos((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)ch, c, &wa[iw], &wa[ix2], &wa[ix3]); |
769 | | |
770 | 2.90M | na = 1 - na; |
771 | 2.90M | break; |
772 | 445k | case 2: |
773 | 445k | if (na == 0) |
774 | 445k | 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 | 445k | na = 1 - na; |
779 | 445k | break; |
780 | 181k | case 3: |
781 | 181k | ix2 = iw + ido; |
782 | | |
783 | 181k | if (na == 0) |
784 | 103k | passf3((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw], &wa[ix2], isign); |
785 | 78.1k | else |
786 | 78.1k | passf3((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)ch, c, &wa[iw], &wa[ix2], isign); |
787 | | |
788 | 181k | na = 1 - na; |
789 | 181k | break; |
790 | 181k | case 5: |
791 | 181k | ix2 = iw + ido; |
792 | 181k | ix3 = ix2 + ido; |
793 | 181k | ix4 = ix3 + ido; |
794 | | |
795 | 181k | if (na == 0) |
796 | 180k | passf5((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw], &wa[ix2], &wa[ix3], &wa[ix4], isign); |
797 | 1.65k | else |
798 | 1.65k | 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 | 181k | na = 1 - na; |
801 | 181k | break; |
802 | 3.71M | } |
803 | | |
804 | 3.71M | l1 = l2; |
805 | 3.71M | iw += (ip-1) * ido; |
806 | 3.71M | } |
807 | | |
808 | 939k | if (na == 0) |
809 | 3.22k | return; |
810 | | |
811 | 257M | for (i = 0; i < n; i++) |
812 | 256M | { |
813 | 256M | RE(c[i]) = RE(ch[i]); |
814 | 256M | IM(c[i]) = IM(ch[i]); |
815 | 256M | } |
816 | 935k | } |
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 | 11.8k | { |
822 | 11.8k | uint16_t i; |
823 | 11.8k | uint16_t k1, l1, l2; |
824 | 11.8k | uint16_t na, nf, ip, iw, ix2, ix3, ix4, ido; |
825 | | |
826 | 11.8k | nf = ifac[1]; |
827 | 11.8k | na = 0; |
828 | 11.8k | l1 = 1; |
829 | 11.8k | iw = 0; |
830 | | |
831 | 70.7k | for (k1 = 2; k1 <= nf+1; k1++) |
832 | 58.9k | { |
833 | 58.9k | ip = ifac[k1]; |
834 | 58.9k | l2 = ip*l1; |
835 | 58.9k | ido = n / l2; |
836 | | |
837 | 58.9k | switch (ip) |
838 | 58.9k | { |
839 | 34.0k | case 4: |
840 | 34.0k | ix2 = iw + ido; |
841 | 34.0k | ix3 = ix2 + ido; |
842 | | |
843 | 34.0k | if (na == 0) |
844 | 17.0k | passf4neg((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw], &wa[ix2], &wa[ix3]); |
845 | 17.0k | else |
846 | 17.0k | passf4neg((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)ch, c, &wa[iw], &wa[ix2], &wa[ix3]); |
847 | | |
848 | 34.0k | na = 1 - na; |
849 | 34.0k | break; |
850 | 11.5k | case 2: |
851 | 11.5k | if (na == 0) |
852 | 11.5k | 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 | 11.5k | na = 1 - na; |
857 | 11.5k | break; |
858 | 6.65k | case 3: |
859 | 6.65k | ix2 = iw + ido; |
860 | | |
861 | 6.65k | if (na == 0) |
862 | 79 | passf3((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw], &wa[ix2], isign); |
863 | 6.57k | else |
864 | 6.57k | passf3((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)ch, c, &wa[iw], &wa[ix2], isign); |
865 | | |
866 | 6.65k | na = 1 - na; |
867 | 6.65k | break; |
868 | 6.65k | case 5: |
869 | 6.65k | ix2 = iw + ido; |
870 | 6.65k | ix3 = ix2 + ido; |
871 | 6.65k | ix4 = ix3 + ido; |
872 | | |
873 | 6.65k | if (na == 0) |
874 | 6.57k | passf5((const uint16_t)ido, (const uint16_t)l1, (const complex_t*)c, ch, &wa[iw], &wa[ix2], &wa[ix3], &wa[ix4], isign); |
875 | 79 | else |
876 | 79 | 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 | 6.65k | na = 1 - na; |
879 | 6.65k | break; |
880 | 58.9k | } |
881 | | |
882 | 58.9k | l1 = l2; |
883 | 58.9k | iw += (ip-1) * ido; |
884 | 58.9k | } |
885 | | |
886 | 11.8k | if (na == 0) |
887 | 293 | return; |
888 | | |
889 | 5.71M | for (i = 0; i < n; i++) |
890 | 5.70M | { |
891 | 5.70M | RE(c[i]) = RE(ch[i]); |
892 | 5.70M | IM(c[i]) = IM(ch[i]); |
893 | 5.70M | } |
894 | 11.5k | } |
895 | | |
896 | | void cfftf(cfft_info *cfft, complex_t *c) |
897 | 11.8k | { |
898 | 11.8k | cfftf1neg(cfft->n, c, cfft->work, (const uint16_t*)cfft->ifac, (const complex_t*)cfft->tab, -1); |
899 | 11.8k | } |
900 | | |
901 | | void cfftb(cfft_info *cfft, complex_t *c) |
902 | 939k | { |
903 | 939k | cfftf1pos(cfft->n, c, cfft->work, (const uint16_t*)cfft->ifac, (const complex_t*)cfft->tab, +1); |
904 | 939k | } |
905 | | |
906 | | static void cffti1(uint16_t n, complex_t *wa, uint16_t *ifac) |
907 | 56.5k | { |
908 | 56.5k | 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 | 56.5k | uint16_t ntry = 0, i, j; |
916 | 56.5k | uint16_t ib; |
917 | 56.5k | uint16_t nf, nl, nq, nr; |
918 | | |
919 | 56.5k | nl = n; |
920 | 56.5k | nf = 0; |
921 | 56.5k | j = 0; |
922 | | |
923 | 155k | startloop: |
924 | 155k | j++; |
925 | | |
926 | 155k | if (j <= 4) |
927 | 155k | ntry = ntryh[j-1]; |
928 | 0 | else |
929 | 0 | ntry += 2; |
930 | | |
931 | 155k | do |
932 | 324k | { |
933 | 324k | nq = nl / ntry; |
934 | 324k | nr = nl - ntry*nq; |
935 | | |
936 | 324k | if (nr != 0) |
937 | 98.7k | goto startloop; |
938 | | |
939 | 226k | nf++; |
940 | 226k | ifac[nf+1] = ntry; |
941 | 226k | nl = nq; |
942 | | |
943 | 226k | if (ntry == 2 && nf != 1) |
944 | 22.2k | { |
945 | 106k | for (i = 2; i <= nf; i++) |
946 | 83.9k | { |
947 | 83.9k | ib = nf - i + 2; |
948 | 83.9k | ifac[ib+1] = ifac[ib]; |
949 | 83.9k | } |
950 | 22.2k | ifac[2] = 2; |
951 | 22.2k | } |
952 | 226k | } while (nl != 1); |
953 | | |
954 | 56.5k | ifac[0] = n; |
955 | 56.5k | 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 | 184k | for (k1 = 1; k1 <= nf; k1++) |
963 | 147k | { |
964 | 147k | ip = ifac[k1+1]; |
965 | 147k | ld = 0; |
966 | 147k | l2 = l1*ip; |
967 | 147k | ido = n / l2; |
968 | 147k | ipm = ip - 1; |
969 | | |
970 | 560k | for (j = 0; j < ipm; j++) |
971 | 413k | { |
972 | 413k | i1 = i; |
973 | 413k | RE(wa[i]) = 1.0; |
974 | 413k | IM(wa[i]) = 0.0; |
975 | | ld += l1; |
976 | | fi = 0; |
977 | | argld = ld*argh; |
978 | | |
979 | 10.6M | for (ii = 0; ii < ido; ii++) |
980 | 10.2M | { |
981 | 10.2M | i++; |
982 | 10.2M | fi++; |
983 | 10.2M | arg = fi * argld; |
984 | 10.2M | RE(wa[i]) = (real_t)cos(arg); |
985 | 10.2M | #if 1 |
986 | 10.2M | IM(wa[i]) = (real_t)sin(arg); |
987 | | #else |
988 | | IM(wa[i]) = (real_t)-sin(arg); |
989 | | #endif |
990 | 10.2M | } |
991 | | |
992 | 413k | if (ip > 5) |
993 | 0 | { |
994 | 0 | RE(wa[i1]) = RE(wa[i]); |
995 | 0 | IM(wa[i1]) = IM(wa[i]); |
996 | 0 | } |
997 | 413k | } |
998 | 147k | l1 = l2; |
999 | 147k | } |
1000 | | #else // FIXED_POINT |
1001 | | (void)wa; /* whoa! does it work for FIXED_POINT? */ |
1002 | | #endif |
1003 | 56.5k | } Line | Count | Source | 907 | 19.5k | { | 908 | 19.5k | 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 | 19.5k | uint16_t ntry = 0, i, j; | 916 | 19.5k | uint16_t ib; | 917 | 19.5k | uint16_t nf, nl, nq, nr; | 918 | | | 919 | 19.5k | nl = n; | 920 | 19.5k | nf = 0; | 921 | 19.5k | j = 0; | 922 | | | 923 | 53.8k | startloop: | 924 | 53.8k | j++; | 925 | | | 926 | 53.8k | if (j <= 4) | 927 | 53.8k | ntry = ntryh[j-1]; | 928 | 0 | else | 929 | 0 | ntry += 2; | 930 | | | 931 | 53.8k | do | 932 | 112k | { | 933 | 112k | nq = nl / ntry; | 934 | 112k | nr = nl - ntry*nq; | 935 | | | 936 | 112k | if (nr != 0) | 937 | 34.2k | goto startloop; | 938 | | | 939 | 78.2k | nf++; | 940 | 78.2k | ifac[nf+1] = ntry; | 941 | 78.2k | nl = nq; | 942 | | | 943 | 78.2k | if (ntry == 2 && nf != 1) | 944 | 7.12k | { | 945 | 33.8k | for (i = 2; i <= nf; i++) | 946 | 26.7k | { | 947 | 26.7k | ib = nf - i + 2; | 948 | 26.7k | ifac[ib+1] = ifac[ib]; | 949 | 26.7k | } | 950 | 7.12k | ifac[2] = 2; | 951 | 7.12k | } | 952 | 78.2k | } while (nl != 1); | 953 | | | 954 | 19.5k | ifac[0] = n; | 955 | 19.5k | 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 | 19.5k | (void)wa; /* whoa! does it work for FIXED_POINT? */ | 1002 | 19.5k | #endif | 1003 | 19.5k | } |
Line | Count | Source | 907 | 36.9k | { | 908 | 36.9k | static uint16_t ntryh[4] = {3, 4, 2, 5}; | 909 | 36.9k | #ifndef FIXED_POINT | 910 | 36.9k | real_t arg, argh, argld, fi; | 911 | 36.9k | uint16_t ido, ipm; | 912 | 36.9k | uint16_t i1, k1, l1, l2; | 913 | 36.9k | uint16_t ld, ii, ip; | 914 | 36.9k | #endif | 915 | 36.9k | uint16_t ntry = 0, i, j; | 916 | 36.9k | uint16_t ib; | 917 | 36.9k | uint16_t nf, nl, nq, nr; | 918 | | | 919 | 36.9k | nl = n; | 920 | 36.9k | nf = 0; | 921 | 36.9k | j = 0; | 922 | | | 923 | 101k | startloop: | 924 | 101k | j++; | 925 | | | 926 | 101k | if (j <= 4) | 927 | 101k | ntry = ntryh[j-1]; | 928 | 0 | else | 929 | 0 | ntry += 2; | 930 | | | 931 | 101k | do | 932 | 212k | { | 933 | 212k | nq = nl / ntry; | 934 | 212k | nr = nl - ntry*nq; | 935 | | | 936 | 212k | if (nr != 0) | 937 | 64.5k | goto startloop; | 938 | | | 939 | 147k | nf++; | 940 | 147k | ifac[nf+1] = ntry; | 941 | 147k | nl = nq; | 942 | | | 943 | 147k | if (ntry == 2 && nf != 1) | 944 | 15.1k | { | 945 | 72.3k | for (i = 2; i <= nf; i++) | 946 | 57.2k | { | 947 | 57.2k | ib = nf - i + 2; | 948 | 57.2k | ifac[ib+1] = ifac[ib]; | 949 | 57.2k | } | 950 | 15.1k | ifac[2] = 2; | 951 | 15.1k | } | 952 | 147k | } while (nl != 1); | 953 | | | 954 | 36.9k | ifac[0] = n; | 955 | 36.9k | ifac[1] = nf; | 956 | | | 957 | 36.9k | #ifndef FIXED_POINT | 958 | 36.9k | argh = (real_t)2.0*(real_t)M_PI / (real_t)n; | 959 | 36.9k | i = 0; | 960 | 36.9k | l1 = 1; | 961 | | | 962 | 184k | for (k1 = 1; k1 <= nf; k1++) | 963 | 147k | { | 964 | 147k | ip = ifac[k1+1]; | 965 | 147k | ld = 0; | 966 | 147k | l2 = l1*ip; | 967 | 147k | ido = n / l2; | 968 | 147k | ipm = ip - 1; | 969 | | | 970 | 560k | for (j = 0; j < ipm; j++) | 971 | 413k | { | 972 | 413k | i1 = i; | 973 | 413k | RE(wa[i]) = 1.0; | 974 | 413k | IM(wa[i]) = 0.0; | 975 | 413k | ld += l1; | 976 | 413k | fi = 0; | 977 | 413k | argld = ld*argh; | 978 | | | 979 | 10.6M | for (ii = 0; ii < ido; ii++) | 980 | 10.2M | { | 981 | 10.2M | i++; | 982 | 10.2M | fi++; | 983 | 10.2M | arg = fi * argld; | 984 | 10.2M | RE(wa[i]) = (real_t)cos(arg); | 985 | 10.2M | #if 1 | 986 | 10.2M | IM(wa[i]) = (real_t)sin(arg); | 987 | | #else | 988 | | IM(wa[i]) = (real_t)-sin(arg); | 989 | | #endif | 990 | 10.2M | } | 991 | | | 992 | 413k | if (ip > 5) | 993 | 0 | { | 994 | 0 | RE(wa[i1]) = RE(wa[i]); | 995 | 0 | IM(wa[i1]) = IM(wa[i]); | 996 | 0 | } | 997 | 413k | } | 998 | 147k | l1 = l2; | 999 | 147k | } | 1000 | | #else // FIXED_POINT | 1001 | | (void)wa; /* whoa! does it work for FIXED_POINT? */ | 1002 | | #endif | 1003 | 36.9k | } |
|
1004 | | |
1005 | | cfft_info *cffti(uint16_t n) |
1006 | 3.60k | { |
1007 | 3.60k | cfft_info *cfft = (cfft_info*)faad_malloc(sizeof(cfft_info)); |
1008 | | |
1009 | 3.60k | cfft->n = n; |
1010 | 3.60k | 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 | 3.60k | cffti1(n, NULL, cfft->ifac); |
1018 | | |
1019 | 3.60k | switch (n) |
1020 | 3.60k | { |
1021 | 1.16k | case 64: cfft->tab = (complex_t*)cfft_tab_64; break; |
1022 | 1.16k | 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 | 640 | case 60: cfft->tab = (complex_t*)cfft_tab_60; break; |
1029 | 640 | 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 | 3.60k | } |
1036 | 3.60k | #endif |
1037 | | |
1038 | 3.60k | return cfft; |
1039 | 3.60k | } |
1040 | | |
1041 | | void cfftu(cfft_info *cfft) |
1042 | 19.5k | { |
1043 | 19.5k | if (cfft->work) faad_free(cfft->work); |
1044 | | #ifndef FIXED_POINT |
1045 | | if (cfft->tab) faad_free(cfft->tab); |
1046 | | #endif |
1047 | | |
1048 | 19.5k | if (cfft) faad_free(cfft); |
1049 | 19.5k | } |
1050 | | |