libdspl-2.0
Digital Signal Processing Algorithm Library
fft_subkernel.c
1/*
2* Copyright (c) 2015-2022 Sergey Bakhurin
3* Digital Signal Processing Library [http://dsplib.org]
4*
5* This file is part of libdspl-2.0.
6*
7* is free software: you can redistribute it and/or modify
8* it under the terms of the GNU Lesser General Public License as published by
9* the Free Software Foundation, either version 3 of the License, or
10* (at your option) any later version.
11*
12* DSPL is distributed in the hope that it will be useful,
13* but WITHOUT ANY WARRANTY; without even the implied warranty of
14* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15* GNU General Public License for more details.
16*
17* You should have received a copy of the GNU Lesser General Public License
18* along with Foobar. If not, see <http://www.gnu.org/licenses/>.
19*/
20
21#include <stdlib.h>
22#include <stdio.h>
23#include <string.h>
24
25#include "dspl.h"
26#include "dft.h"
27
28
29
30/*******************************************************************************
312 points DFT
32*******************************************************************************/
33void dft2 (complex_t *x, complex_t* y)
34{
35 RE(y[0]) = RE(x[0]) + RE(x[1]);
36 IM(y[0]) = IM(x[0]) + IM(x[1]);
37
38 RE(y[1]) = RE(x[0]) - RE(x[1]);
39 IM(y[1]) = IM(x[0]) - IM(x[1]);
40}
41
42
43
44/*******************************************************************************
453 points DFT (Winograd algorithm)
46*******************************************************************************/
47void dft3 (complex_t *x, complex_t* y)
48{
49 double a, b, c, d;
50
51 a = RE(x[0]) - 0.5 * (RE(x[1]) + RE(x[2]));
52 b = DFT3_W * (IM(x[1]) - IM(x[2]));
53 c = IM(x[0]) - 0.5 * (IM(x[1]) + IM(x[2]));
54 d = DFT3_W * (RE(x[2]) - RE(x[1]));
55
56 RE(y[0]) = RE(x[0]) + RE(x[1]) + RE(x[2]);
57 IM(y[0]) = IM(x[0]) + IM(x[1]) + IM(x[2]);
58
59 RE(y[1]) = a + b;
60 IM(y[1]) = c + d;
61
62 RE(y[2]) = a - b;
63 IM(y[2]) = c - d;
64}
65
66
67
68/*******************************************************************************
694 points DFT (Winograd algorithm)
70*******************************************************************************/
71void dft4 (complex_t *x, complex_t* y)
72{
73 double a0, b0, c0, d0, a1, b1, c1, d1;
74
75 a0 = RE(x[0]) + RE(x[2]);
76 b0 = RE(x[1]) + RE(x[3]);
77 c0 = IM(x[0]) + IM(x[2]);
78 d0 = IM(x[1]) + IM(x[3]);
79
80 a1 = RE(x[0]) - RE(x[2]);
81 b1 = RE(x[1]) - RE(x[3]);
82 c1 = IM(x[0]) - IM(x[2]);
83 d1 = IM(x[1]) - IM(x[3]);
84
85 RE(y[0]) = a0 + b0;
86 IM(y[0]) = c0 + d0;
87
88 RE(y[1]) = a1 + d1;
89 IM(y[1]) = c1 - b1;
90
91 RE(y[2]) = a0 - b0;
92 IM(y[2]) = c0 - d0;
93
94 RE(y[3]) = a1 - d1;
95 IM(y[3]) = c1 + b1;
96}
97
98
99
100
101/*******************************************************************************
1025 points DFT (Winograd algorithm)
103*******************************************************************************/
104void dft5 (complex_t *x, complex_t* y)
105{
106 complex_t sum[14];
107 complex_t mul[6];
108
109 RE(sum[1]) = RE(x[1]) + RE(x[4]);
110 IM(sum[1]) = IM(x[1]) + IM(x[4]);
111
112 RE(sum[2]) = RE(x[1]) - RE(x[4]);
113 IM(sum[2]) = IM(x[1]) - IM(x[4]);
114
115 RE(sum[3]) = RE(x[3]) + RE(x[2]);
116 IM(sum[3]) = IM(x[3]) + IM(x[2]);
117
118 RE(sum[4]) = RE(x[3]) - RE(x[2]);
119 IM(sum[4]) = IM(x[3]) - IM(x[2]);
120
121 RE(sum[5]) = RE(sum[1]) + RE(sum[3]);
122 IM(sum[5]) = IM(sum[1]) + IM(sum[3]);
123
124 RE(sum[6]) = RE(sum[1]) - RE(sum[3]);
125 IM(sum[6]) = IM(sum[1]) - IM(sum[3]);
126
127 RE(sum[7]) = RE(sum[2]) + RE(sum[4]);
128 IM(sum[7]) = IM(sum[2]) + IM(sum[4]);
129
130 RE(y[0]) = RE(sum[5]) + RE(x[0]);
131 IM(y[0]) = IM(sum[5]) + IM(x[0]);
132
133 RE(mul[1]) = RE(sum[5]) * DFT5_W1;
134 IM(mul[1]) = IM(sum[5]) * DFT5_W1;
135
136 RE(mul[2]) = RE(sum[6]) * DFT5_W2;
137 IM(mul[2]) = IM(sum[6]) * DFT5_W2;
138
139 RE(mul[3]) = -IM(sum[2]) * DFT5_W3;
140 IM(mul[3]) = RE(sum[2]) * DFT5_W3;
141
142 RE(mul[4]) = -IM(sum[7]) * DFT5_W4;
143 IM(mul[4]) = RE(sum[7]) * DFT5_W4;
144
145 RE(mul[5]) = -IM(sum[4]) * DFT5_W5;
146 IM(mul[5]) = RE(sum[4]) * DFT5_W5;
147
148 RE(sum[9]) = RE(y[0]) + RE(mul[1]);
149 IM(sum[9]) = IM(y[0]) + IM(mul[1]);
150
151 RE(sum[10]) = RE(sum[9]) + RE(mul[2]);
152 IM(sum[10]) = IM(sum[9]) + IM(mul[2]);
153
154 RE(sum[11]) = RE(sum[9]) - RE(mul[2]);
155 IM(sum[11]) = IM(sum[9]) - IM(mul[2]);
156
157 RE(sum[12]) = RE(mul[3]) - RE(mul[4]);
158 IM(sum[12]) = IM(mul[3]) - IM(mul[4]);
159
160 RE(sum[13]) = RE(mul[4]) + RE(mul[5]);
161 IM(sum[13]) = IM(mul[4]) + IM(mul[5]);
162
163 RE(y[4]) = RE(sum[10]) + RE(sum[12]);
164 IM(y[4]) = IM(sum[10]) + IM(sum[12]);
165
166 RE(y[3]) = RE(sum[11]) + RE(sum[13]);
167 IM(y[3]) = IM(sum[11]) + IM(sum[13]);
168
169 RE(y[2]) = RE(sum[11]) - RE(sum[13]);
170 IM(y[2]) = IM(sum[11]) - IM(sum[13]);
171
172 RE(y[1]) = RE(sum[10]) - RE(sum[12]);
173 IM(y[1]) = IM(sum[10]) - IM(sum[12]);
174}
175
176
177
178/*******************************************************************************
1797 points DFT (Winograd algorithm)
180*******************************************************************************/
181void dft7 (complex_t *x, complex_t* y)
182{
183 complex_t sum[31];
184 complex_t mul[9];
185
186 RE(sum[1]) = RE(x[1]) + RE(x[6]);
187 IM(sum[1]) = IM(x[1]) + IM(x[6]);
188
189 RE(sum[2]) = RE(x[1]) - RE(x[6]);
190 IM(sum[2]) = IM(x[1]) - IM(x[6]);
191
192 RE(sum[3]) = RE(x[4]) + RE(x[3]);
193 IM(sum[3]) = IM(x[4]) + IM(x[3]);
194
195 RE(sum[4]) = RE(x[4]) - RE(x[3]);
196 IM(sum[4]) = IM(x[4]) - IM(x[3]);
197
198 /* Winograd paper mistake?! */
199 RE(sum[5]) = RE(x[2]) + RE(x[5]);
200 IM(sum[5]) = IM(x[2]) + IM(x[5]);
201
202 /* Winograd paper mistake?! */
203 RE(sum[6]) = RE(x[2]) - RE(x[5]);
204 IM(sum[6]) = IM(x[2]) - IM(x[5]);
205
206 RE(sum[7]) = RE(sum[1]) + RE(sum[3]);
207 IM(sum[7]) = IM(sum[1]) + IM(sum[3]);
208
209 RE(sum[8]) = RE(sum[7]) + RE(sum[5]);
210 IM(sum[8]) = IM(sum[7]) + IM(sum[5]);
211
212 RE(y[0]) = RE(sum[9]) = RE(sum[8]) + RE(x[0]);
213 IM(y[0]) = IM(sum[9]) = IM(sum[8]) + IM(x[0]);
214
215 RE(sum[10]) = RE(sum[1]) - RE(sum[3]);
216 IM(sum[10]) = IM(sum[1]) - IM(sum[3]);
217
218 RE(sum[11]) = RE(sum[3]) - RE(sum[5]);
219 IM(sum[11]) = IM(sum[3]) - IM(sum[5]);
220
221 RE(sum[12]) = RE(sum[5]) - RE(sum[1]);
222 IM(sum[12]) = IM(sum[5]) - IM(sum[1]);
223
224 RE(sum[13]) = RE(sum[2]) + RE(sum[4]);
225 IM(sum[13]) = IM(sum[2]) + IM(sum[4]);
226
227 RE(sum[14]) = RE(sum[13]) + RE(sum[6]);
228 IM(sum[14]) = IM(sum[13]) + IM(sum[6]);
229
230 RE(sum[15]) = RE(sum[2]) - RE(sum[4]);
231 IM(sum[15]) = IM(sum[2]) - IM(sum[4]);
232
233 RE(sum[16]) = RE(sum[4]) - RE(sum[6]);
234 IM(sum[16]) = IM(sum[4]) - IM(sum[6]);
235
236 RE(sum[17]) = RE(sum[6]) - RE(sum[2]);
237 IM(sum[17]) = IM(sum[6]) - IM(sum[2]);
238
239 RE(mul[1]) = DFT7_W1 * RE(sum[8]);
240 IM(mul[1]) = DFT7_W1 * IM(sum[8]);
241
242 RE(mul[2]) = DFT7_W2 * RE(sum[10]);
243 IM(mul[2]) = DFT7_W2 * IM(sum[10]);
244
245 RE(mul[3]) = DFT7_W3 * RE(sum[11]);
246 IM(mul[3]) = DFT7_W3 * IM(sum[11]);
247
248 RE(mul[4]) = DFT7_W4 * RE(sum[12]);
249 IM(mul[4]) = DFT7_W4 * IM(sum[12]);
250
251 RE(mul[5]) = -DFT7_W5 * IM(sum[14]);
252 IM(mul[5]) = DFT7_W5 * RE(sum[14]);
253
254 RE(mul[6]) = -DFT7_W6 * IM(sum[15]);
255 IM(mul[6]) = DFT7_W6 * RE(sum[15]);
256
257 RE(mul[7]) = -DFT7_W7 * IM(sum[16]);
258 IM(mul[7]) = DFT7_W7 * RE(sum[16]);
259
260 RE(mul[8]) = -DFT7_W8 * IM(sum[17]);
261 IM(mul[8]) = DFT7_W8 * RE(sum[17]);
262
263 RE(sum[18]) = RE(y[0]) + RE(mul[1]);
264 IM(sum[18]) = IM(y[0]) + IM(mul[1]);
265
266 RE(sum[19]) = RE(sum[18]) + RE(mul[2]);
267 IM(sum[19]) = IM(sum[18]) + IM(mul[2]);
268
269 RE(sum[20]) = RE(sum[19]) + RE(mul[3]);
270 IM(sum[20]) = IM(sum[19]) + IM(mul[3]);
271
272 RE(sum[21]) = RE(sum[18]) - RE(mul[2]);
273 IM(sum[21]) = IM(sum[18]) - IM(mul[2]);
274
275 RE(sum[22]) = RE(sum[21]) - RE(mul[4]);
276 IM(sum[22]) = IM(sum[21]) - IM(mul[4]);
277
278 RE(sum[23]) = RE(sum[18]) - RE(mul[3]);
279 IM(sum[23]) = IM(sum[18]) - IM(mul[3]);
280
281 RE(sum[24]) = RE(sum[23]) + RE(mul[4]);
282 IM(sum[24]) = IM(sum[23]) + IM(mul[4]);
283
284 RE(sum[25]) = RE(mul[5]) + RE(mul[6]);
285 IM(sum[25]) = IM(mul[5]) + IM(mul[6]);
286
287 RE(sum[26]) = RE(sum[25]) + RE(mul[7]);
288 IM(sum[26]) = IM(sum[25]) + IM(mul[7]);
289
290 RE(sum[27]) = RE(mul[5]) - RE(mul[6]);
291 IM(sum[27]) = IM(mul[5]) - IM(mul[6]);
292
293 RE(sum[28]) = RE(sum[27]) - RE(mul[8]);
294 IM(sum[28]) = IM(sum[27]) - IM(mul[8]);
295
296 RE(sum[29]) = RE(mul[5]) - RE(mul[7]);
297 IM(sum[29]) = IM(mul[5]) - IM(mul[7]);
298
299 RE(sum[30]) = RE(sum[29]) + RE(mul[8]);
300 IM(sum[30]) = IM(sum[29]) + IM(mul[8]);
301
302 RE(y[1]) = RE(sum[20]) + RE(sum[26]);
303 IM(y[1]) = IM(sum[20]) + IM(sum[26]);
304
305 RE(y[2]) = RE(sum[22]) + RE(sum[28]);
306 IM(y[2]) = IM(sum[22]) + IM(sum[28]);
307
308 RE(y[3]) = RE(sum[24]) - RE(sum[30]);
309 IM(y[3]) = IM(sum[24]) - IM(sum[30]);
310
311 RE(y[4]) = RE(sum[24]) + RE(sum[30]);
312 IM(y[4]) = IM(sum[24]) + IM(sum[30]);
313
314 RE(y[5]) = RE(sum[22]) - RE(sum[28]);
315 IM(y[5]) = IM(sum[22]) - IM(sum[28]);
316
317 RE(y[6]) = RE(sum[20]) - RE(sum[26]);
318 IM(y[6]) = IM(sum[20]) - IM(sum[26]);
319}
320
321
322/*******************************************************************************
3238 points DFT
324*******************************************************************************/
325void dft8 (complex_t *x, complex_t* y)
326{
327 complex_t t0[8];
328 complex_t t1[8];
329 double tmp;
330
331 transpose2x4(x, t0);
332
333 dft4(t0, t1);
334 dft4(t0+4, t1+4);
335
336 /* 0.707106781186548 - 707106781186548i */
337 tmp = (RE(t1[5]) + IM(t1[5])) * DFT8_W;
338 IM(t1[5]) = (IM(t1[5]) - RE(t1[5])) * DFT8_W;
339 RE(t1[5]) = tmp;
340
341 /* 0.000000000000000 - 1.000000000000000i */
342 tmp = RE(t1[6]);
343 RE(t1[6]) = IM(t1[6]);
344 IM(t1[6]) = -tmp;
345
346 /* -0.707106781186548 - 707106781186548i */
347 tmp = (IM(t1[7]) - RE(t1[7])) * DFT8_W;
348 IM(t1[7]) = -(IM(t1[7]) + RE(t1[7])) * DFT8_W;
349 RE(t1[7]) = tmp;
350
351 transpose4x2(t1, t0);
352
353 RE(t1[0]) = RE(t0[0]) + RE(t0[1]);
354 IM(t1[0]) = IM(t0[0]) + IM(t0[1]);
355 RE(t1[1]) = RE(t0[0]) - RE(t0[1]);
356 IM(t1[1]) = IM(t0[0]) - IM(t0[1]);
357
358 RE(t1[2]) = RE(t0[2]) + RE(t0[3]);
359 IM(t1[2]) = IM(t0[2]) + IM(t0[3]);
360 RE(t1[3]) = RE(t0[2]) - RE(t0[3]);
361 IM(t1[3]) = IM(t0[2]) - IM(t0[3]);
362
363 RE(t1[4]) = RE(t0[4]) + RE(t0[5]);
364 IM(t1[4]) = IM(t0[4]) + IM(t0[5]);
365 RE(t1[5]) = RE(t0[4]) - RE(t0[5]);
366 IM(t1[5]) = IM(t0[4]) - IM(t0[5]);
367
368 RE(t1[6]) = RE(t0[6]) + RE(t0[7]);
369 IM(t1[6]) = IM(t0[6]) + IM(t0[7]);
370 RE(t1[7]) = RE(t0[6]) - RE(t0[7]);
371 IM(t1[7]) = IM(t0[6]) - IM(t0[7]);
372
373 transpose2x4(t1, y);
374}
375
376/*******************************************************************************
37716 points DFT (Winograd algorithm)
378*******************************************************************************/
379void dft16(complex_t *x, complex_t* y)
380{
381 complex_t t0[16];
382 complex_t t1[16];
383 double tmp;
384
385 transpose4x4(x, t0);
386
387 dft4(t0, t1);
388 dft4(t0+4, t1+4);
389 dft4(t0+8, t1+8);
390 dft4(t0+12, t1+12);
391
392 /* #define DFT16_W1 0.923879532511287 */
393 /* #define DFT16_W2 0.382683432365090 */
394 /* #define DFT16_W3 0.707106781186548 */
395
396 /* 0.923879532511287 - 0.382683432365090i */
397 tmp = RE(t1[5]) * DFT16_W1 + IM(t1[5]) * DFT16_W2;
398 IM(t1[5]) = -RE(t1[5]) * DFT16_W2 + IM(t1[5]) * DFT16_W1;
399 RE(t1[5]) = tmp;
400
401 /* 0.707106781186548 - 0.707106781186547i */
402 tmp = ( RE(t1[6]) + IM(t1[6])) * DFT16_W3;
403 IM(t1[6]) = (-RE(t1[6]) + IM(t1[6])) * DFT16_W3;
404 RE(t1[6]) = tmp;
405
406 /* 0.382683432365090 - 0.923879532511287i */
407 tmp = RE(t1[7]) * DFT16_W2 + IM(t1[7]) * DFT16_W1;
408 IM(t1[7]) = -RE(t1[7]) * DFT16_W1 + IM(t1[7]) * DFT16_W2;
409 RE(t1[7]) = tmp;
410
411 /* 0.707106781186548 - 0.707106781186547i */
412 tmp = ( RE(t1[9]) + IM(t1[9])) * DFT16_W3;
413 IM(t1[9]) = (-RE(t1[9]) + IM(t1[9])) * DFT16_W3;
414 RE(t1[9]) = tmp;
415
416 /* 0.000000000000000 - 1.000000000000000i */
417 tmp = RE(t1[10]);
418 RE(t1[10]) = IM(t1[10]);
419 IM(t1[10]) = -tmp;
420
421 /* -0.707106781186547 - 0.707106781186548i */
422 tmp = (-RE(t1[11]) + IM(t1[11])) * DFT16_W3;
423 IM(t1[11]) = (-RE(t1[11]) - IM(t1[11])) * DFT16_W3;
424 RE(t1[11]) = tmp;
425
426 /* 0.382683432365090 - 0.923879532511287i */
427 tmp = RE(t1[13]) * DFT16_W2 + IM(t1[13]) * DFT16_W1;
428 IM(t1[13]) = -RE(t1[13]) * DFT16_W1 + IM(t1[13]) * DFT16_W2;
429 RE(t1[13]) = tmp;
430
431 /* -0.707106781186547 - 0.707106781186548i */
432 tmp = (-RE(t1[14]) + IM(t1[14])) * DFT16_W3;
433 IM(t1[14]) = (-RE(t1[14]) - IM(t1[14])) * DFT16_W3;
434 RE(t1[14]) = tmp;
435
436
437 /* -0.923879532511287 + 0.382683432365090i */
438 tmp = -RE(t1[15]) * DFT16_W1 - IM(t1[15]) * DFT16_W2;
439 IM(t1[15]) = RE(t1[15]) * DFT16_W2 - IM(t1[15]) * DFT16_W1;
440 RE(t1[15]) = tmp;
441
442 transpose4x4(t1, t0);
443
444 dft4(t0, t1);
445 dft4(t0+4, t1+4);
446 dft4(t0+8, t1+8);
447 dft4(t0+12, t1+12);
448
449 transpose4x4(t1, y);
450
451}
452
453
454/*******************************************************************************
45532 points DFT (Winograd algorithm)
456*******************************************************************************/
457void dft32(complex_t *x, complex_t* y, complex_t* w)
458{
459 complex_t t0[32];
460 complex_t t1[32];
461
462 int i;
463
464 transpose4x8(x, t0);
465 dft8(t0, t1);
466 dft8(t0+8, t1+8);
467 dft8(t0+16, t1+16);
468 dft8(t0+24, t1+24);
469
470
471 for(i = 0; i < 32; i++)
472 {
473 RE(t0[i]) = CMRE(t1[i], w[i]);
474 IM(t0[i]) = CMIM(t1[i], w[i]);
475 }
476
477 transpose8x4(t0, t1);
478
479 for(i = 0; i < 8; i++)
480 dft4(t1 + i*4, t0 + i*4);
481
482 transpose4x8(t0, y);
483}
484
485
486/*******************************************************************************
48764 points DFT (Winograd algorithm)
488*******************************************************************************/
489void dft64(complex_t *x, complex_t* y, complex_t* w)
490{
491 complex_t t0[64];
492 complex_t t1[64];
493
494 int i;
495
496 transpose8x8(x, t0);
497
498 for(i = 0; i < 8; i++)
499 dft8(t0 + i*8, t1 + i*8);
500
501 for(i = 0; i < 64; i++)
502 {
503 RE(t0[i]) = CMRE(t1[i], w[i]);
504 IM(t0[i]) = CMIM(t1[i], w[i]);
505 }
506
507 transpose8x8(t0, t1);
508
509 for(i = 0; i < 8; i++)
510 dft8(t1 + i*8, t0 + i*8);
511 transpose8x8(t0, y);
512}
513
514
515
516/*******************************************************************************
517256 points DFT (Winograd algorithm)
518*******************************************************************************/
519void dft128(complex_t *x, complex_t* y, complex_t* w)
520{
521 complex_t t0[128];
522 complex_t t1[128];
523
524 int i;
525
526 matrix_transpose_cmplx(x,8,16,t0);
527
528 for(i = 0; i < 8; i++)
529 dft16(t0 + i*16, t1 + i*16);
530
531 for(i = 0; i < 128; i++)
532 {
533 RE(t0[i]) = CMRE(t1[i], w[i]);
534 IM(t0[i]) = CMIM(t1[i], w[i]);
535 }
536
537 matrix_transpose_cmplx(t0, 16, 8, t1);
538
539 for(i = 0; i < 16; i++)
540 dft8(t1 + i*8, t0 + i*8);
541
542 matrix_transpose_cmplx(t0, 8, 16, y);
543}
544
545
546
547/*******************************************************************************
548256 points DFT (Winograd algorithm)
549*******************************************************************************/
550void dft256(complex_t *x, complex_t* y, complex_t* w)
551{
552 complex_t t0[256];
553 complex_t t1[256];
554
555 int i;
556
557 transpose16x16(x, t0);
558
559 for(i = 0; i < 16; i++)
560 dft16(t0 + i*16, t1 + i*16);
561
562 for(i = 0; i < 256; i++)
563 {
564 RE(t0[i]) = CMRE(t1[i], w[i]);
565 IM(t0[i]) = CMIM(t1[i], w[i]);
566 }
567
568 transpose16x16(t0, t1);
569
570 for(i = 0; i < 16; i++)
571 dft16(t1 + i*16, t0 + i*16);
572
573 transpose16x16(t0, y);
574}
575
576
577/*******************************************************************************
578512 points DFT (Winograd algorithm)
579*******************************************************************************/
580void dft512(complex_t *x, complex_t* y, complex_t* w, complex_t* w32)
581{
582 complex_t t0[512];
583 complex_t t1[512];
584
585 int i;
586
587 matrix_transpose_cmplx(x,16,32,t0);
588
589 for(i = 0; i < 16; i++)
590 dft32(t0 + i*32, t1 + i*32, w32);
591
592 for(i = 0; i < 512; i++)
593 {
594 RE(t0[i]) = CMRE(t1[i], w[i]);
595 IM(t0[i]) = CMIM(t1[i], w[i]);
596 }
597
598 matrix_transpose_cmplx(t0, 32, 16, t1);
599
600 for(i = 0; i < 32; i++)
601 dft16(t1 + i*16, t0 + i*16);
602
603 matrix_transpose_cmplx(t0, 16, 32, y);
604}
605
606
607/*******************************************************************************
6081024 points DFT (Winograd algorithm)
609*******************************************************************************/
610void dft1024(complex_t *x, complex_t* y, complex_t* w, complex_t* w32)
611{
612 complex_t t0[1024];
613 complex_t t1[1024];
614
615 int i;
616
617 matrix_transpose_cmplx(x,32,32,t0);
618
619 for(i = 0; i < 32; i++)
620 dft32(t0 + i*32, t1 + i*32, w32);
621
622 for(i = 0; i < 1024; i++)
623 {
624 RE(t0[i]) = CMRE(t1[i], w[i]);
625 IM(t0[i]) = CMIM(t1[i], w[i]);
626 }
627
628 matrix_transpose_cmplx(t0, 32, 32, t1);
629
630 for(i = 0; i < 32; i++)
631 dft32(t1 + i*32, t0 + i*32, w32);
632
633 matrix_transpose_cmplx(t0, 32, 32, y);
634}
635
636/*******************************************************************************
6372048 points DFT (Winograd algorithm)
638*******************************************************************************/
639void dft2048(complex_t *x, complex_t* y, complex_t* w,
640 complex_t* w32, complex_t* w64)
641{
642 complex_t *t0 = NULL;
643 complex_t *t1 = NULL;
644
645 int i;
646
647 t0 = (complex_t*)malloc(2048*sizeof(complex_t));
648 t1 = (complex_t*)malloc(2048*sizeof(complex_t));
649
650 matrix_transpose_cmplx(x,32,64,t0);
651
652 for(i = 0; i < 32; i++)
653 dft64(t0 + i*64, t1 + i*64, w64);
654
655 for(i = 0; i < 2048; i++)
656 {
657 RE(t0[i]) = CMRE(t1[i], w[i]);
658 IM(t0[i]) = CMIM(t1[i], w[i]);
659 }
660
661 matrix_transpose_cmplx(t0, 64, 32, t1);
662
663 for(i = 0; i < 64; i++)
664 dft32(t1 + i*32, t0 + i*32, w32);
665
666 matrix_transpose_cmplx(t0, 32, 64, y);
667
668 free(t0);
669 free(t1);
670}
671
672
673
674/*******************************************************************************
6754096 points DFT (Winograd algorithm)
676*******************************************************************************/
677void dft4096(complex_t *x, complex_t* y, complex_t* w, complex_t* w256)
678{
679 complex_t *t0 = NULL;
680 complex_t *t1 = NULL;
681
682 int i;
683
684 t0 = (complex_t*)malloc(4096*sizeof(complex_t));
685 t1 = (complex_t*)malloc(4096*sizeof(complex_t));
686
687 matrix_transpose_cmplx(x,16,256,t0);
688
689 for(i = 0; i < 16; i++)
690 dft256(t0 + i*256, t1 + i*256, w256);
691
692 for(i = 0; i < 4096; i++)
693 {
694 RE(t0[i]) = CMRE(t1[i], w[i]);
695 IM(t0[i]) = CMIM(t1[i], w[i]);
696 }
697
698 matrix_transpose_cmplx(t0, 256, 16, t1);
699
700 for(i = 0; i < 256; i++)
701 dft16(t1 + i*16, t0 + i*16);
702
703 matrix_transpose_cmplx(t0, 16, 256, y);
704
705 free(t0);
706 free(t1);
707}
708
709
710
711/*******************************************************************************
7124 x 2 matrix transpose
713*******************************************************************************/
714void transpose4x2(complex_t *x, complex_t* y)
715{
716 RE(y[ 0]) = RE(x[ 0]); IM(y[ 0]) = IM(x[ 0]);
717 RE(y[ 1]) = RE(x[ 4]); IM(y[ 1]) = IM(x[ 4]);
718 RE(y[ 2]) = RE(x[ 1]); IM(y[ 2]) = IM(x[ 1]);
719 RE(y[ 3]) = RE(x[ 5]); IM(y[ 3]) = IM(x[ 5]);
720 RE(y[ 4]) = RE(x[ 2]); IM(y[ 4]) = IM(x[ 2]);
721 RE(y[ 5]) = RE(x[ 6]); IM(y[ 5]) = IM(x[ 6]);
722 RE(y[ 6]) = RE(x[ 3]); IM(y[ 6]) = IM(x[ 3]);
723 RE(y[ 7]) = RE(x[ 7]); IM(y[ 7]) = IM(x[ 7]);
724}
725
726
727/*******************************************************************************
7282 x 4 matrix transpose
729*******************************************************************************/
730void transpose2x4(complex_t *x, complex_t* y)
731{
732 RE(y[ 0]) = RE(x[ 0]); IM(y[ 0]) = IM(x[ 0]);
733 RE(y[ 1]) = RE(x[ 2]); IM(y[ 1]) = IM(x[ 2]);
734 RE(y[ 2]) = RE(x[ 4]); IM(y[ 2]) = IM(x[ 4]);
735 RE(y[ 3]) = RE(x[ 6]); IM(y[ 3]) = IM(x[ 6]);
736 RE(y[ 4]) = RE(x[ 1]); IM(y[ 4]) = IM(x[ 1]);
737 RE(y[ 5]) = RE(x[ 3]); IM(y[ 5]) = IM(x[ 3]);
738 RE(y[ 6]) = RE(x[ 5]); IM(y[ 6]) = IM(x[ 5]);
739 RE(y[ 7]) = RE(x[ 7]); IM(y[ 7]) = IM(x[ 7]);
740}
741
742
743
744/*******************************************************************************
7454 x 4 matrix transpose
746*******************************************************************************/
747void transpose4x4(complex_t *x, complex_t* y)
748{
749 RE(y[ 0]) = RE(x[ 0]); IM(y[ 0]) = IM(x[ 0]);
750 RE(y[ 1]) = RE(x[ 4]); IM(y[ 1]) = IM(x[ 4]);
751 RE(y[ 2]) = RE(x[ 8]); IM(y[ 2]) = IM(x[ 8]);
752 RE(y[ 3]) = RE(x[12]); IM(y[ 3]) = IM(x[12]);
753 RE(y[ 4]) = RE(x[ 1]); IM(y[ 4]) = IM(x[ 1]);
754 RE(y[ 5]) = RE(x[ 5]); IM(y[ 5]) = IM(x[ 5]);
755 RE(y[ 6]) = RE(x[ 9]); IM(y[ 6]) = IM(x[ 9]);
756 RE(y[ 7]) = RE(x[13]); IM(y[ 7]) = IM(x[13]);
757 RE(y[ 8]) = RE(x[ 2]); IM(y[ 8]) = IM(x[ 2]);
758 RE(y[ 9]) = RE(x[ 6]); IM(y[ 9]) = IM(x[ 6]);
759 RE(y[10]) = RE(x[10]); IM(y[10]) = IM(x[10]);
760 RE(y[11]) = RE(x[14]); IM(y[11]) = IM(x[14]);
761 RE(y[12]) = RE(x[ 3]); IM(y[12]) = IM(x[ 3]);
762 RE(y[13]) = RE(x[ 7]); IM(y[13]) = IM(x[ 7]);
763 RE(y[14]) = RE(x[11]); IM(y[14]) = IM(x[11]);
764 RE(y[15]) = RE(x[15]); IM(y[15]) = IM(x[15]);
765}
766
767
768/*******************************************************************************
7698 x 4 matrix transpose
770*******************************************************************************/
771void transpose8x4(complex_t *x, complex_t* y)
772{
773 RE(y[ 0]) = RE(x[ 0]); IM(y[ 0]) = IM(x[ 0]);
774 RE(y[ 1]) = RE(x[ 8]); IM(y[ 1]) = IM(x[ 8]);
775 RE(y[ 2]) = RE(x[ 16]); IM(y[ 2]) = IM(x[ 16]);
776 RE(y[ 3]) = RE(x[ 24]); IM(y[ 3]) = IM(x[ 24]);
777 RE(y[ 4]) = RE(x[ 1]); IM(y[ 4]) = IM(x[ 1]);
778 RE(y[ 5]) = RE(x[ 9]); IM(y[ 5]) = IM(x[ 9]);
779 RE(y[ 6]) = RE(x[ 17]); IM(y[ 6]) = IM(x[ 17]);
780 RE(y[ 7]) = RE(x[ 25]); IM(y[ 7]) = IM(x[ 25]);
781 RE(y[ 8]) = RE(x[ 2]); IM(y[ 8]) = IM(x[ 2]);
782 RE(y[ 9]) = RE(x[ 10]); IM(y[ 9]) = IM(x[ 10]);
783 RE(y[10]) = RE(x[ 18]); IM(y[10]) = IM(x[ 18]);
784 RE(y[11]) = RE(x[ 26]); IM(y[11]) = IM(x[ 26]);
785 RE(y[12]) = RE(x[ 3]); IM(y[12]) = IM(x[ 3]);
786 RE(y[13]) = RE(x[ 11]); IM(y[13]) = IM(x[ 11]);
787 RE(y[14]) = RE(x[ 19]); IM(y[14]) = IM(x[ 19]);
788 RE(y[15]) = RE(x[ 27]); IM(y[15]) = IM(x[ 27]);
789 RE(y[16]) = RE(x[ 4]); IM(y[16]) = IM(x[ 4]);
790 RE(y[17]) = RE(x[ 12]); IM(y[17]) = IM(x[ 12]);
791 RE(y[18]) = RE(x[ 20]); IM(y[18]) = IM(x[ 20]);
792 RE(y[19]) = RE(x[ 28]); IM(y[19]) = IM(x[ 28]);
793 RE(y[20]) = RE(x[ 5]); IM(y[20]) = IM(x[ 5]);
794 RE(y[21]) = RE(x[ 13]); IM(y[21]) = IM(x[ 13]);
795 RE(y[22]) = RE(x[ 21]); IM(y[22]) = IM(x[ 21]);
796 RE(y[23]) = RE(x[ 29]); IM(y[23]) = IM(x[ 29]);
797 RE(y[24]) = RE(x[ 6]); IM(y[24]) = IM(x[ 6]);
798 RE(y[25]) = RE(x[ 14]); IM(y[25]) = IM(x[ 14]);
799 RE(y[26]) = RE(x[ 22]); IM(y[26]) = IM(x[ 22]);
800 RE(y[27]) = RE(x[ 30]); IM(y[27]) = IM(x[ 30]);
801 RE(y[28]) = RE(x[ 7]); IM(y[28]) = IM(x[ 7]);
802 RE(y[29]) = RE(x[ 15]); IM(y[29]) = IM(x[ 15]);
803 RE(y[30]) = RE(x[ 23]); IM(y[30]) = IM(x[ 23]);
804 RE(y[31]) = RE(x[ 31]); IM(y[31]) = IM(x[ 31]);
805}
806
807
808/*******************************************************************************
8094 x 8 matrix transpose
810*******************************************************************************/
811void transpose4x8(complex_t *x, complex_t* y)
812{
813 RE(y[ 0]) = RE(x[ 0]); IM(y[ 0]) = IM(x[ 0]);
814 RE(y[ 1]) = RE(x[ 4]); IM(y[ 1]) = IM(x[ 4]);
815 RE(y[ 2]) = RE(x[ 8]); IM(y[ 2]) = IM(x[ 8]);
816 RE(y[ 3]) = RE(x[ 12]); IM(y[ 3]) = IM(x[ 12]);
817 RE(y[ 4]) = RE(x[ 16]); IM(y[ 4]) = IM(x[ 16]);
818 RE(y[ 5]) = RE(x[ 20]); IM(y[ 5]) = IM(x[ 20]);
819 RE(y[ 6]) = RE(x[ 24]); IM(y[ 6]) = IM(x[ 24]);
820 RE(y[ 7]) = RE(x[ 28]); IM(y[ 7]) = IM(x[ 28]);
821 RE(y[ 8]) = RE(x[ 1]); IM(y[ 8]) = IM(x[ 1]);
822 RE(y[ 9]) = RE(x[ 5]); IM(y[ 9]) = IM(x[ 5]);
823 RE(y[10]) = RE(x[ 9]); IM(y[10]) = IM(x[ 9]);
824 RE(y[11]) = RE(x[ 13]); IM(y[11]) = IM(x[ 13]);
825 RE(y[12]) = RE(x[ 17]); IM(y[12]) = IM(x[ 17]);
826 RE(y[13]) = RE(x[ 21]); IM(y[13]) = IM(x[ 21]);
827 RE(y[14]) = RE(x[ 25]); IM(y[14]) = IM(x[ 25]);
828 RE(y[15]) = RE(x[ 29]); IM(y[15]) = IM(x[ 29]);
829 RE(y[16]) = RE(x[ 2]); IM(y[16]) = IM(x[ 2]);
830 RE(y[17]) = RE(x[ 6]); IM(y[17]) = IM(x[ 6]);
831 RE(y[18]) = RE(x[ 10]); IM(y[18]) = IM(x[ 10]);
832 RE(y[19]) = RE(x[ 14]); IM(y[19]) = IM(x[ 14]);
833 RE(y[20]) = RE(x[ 18]); IM(y[20]) = IM(x[ 18]);
834 RE(y[21]) = RE(x[ 22]); IM(y[21]) = IM(x[ 22]);
835 RE(y[22]) = RE(x[ 26]); IM(y[22]) = IM(x[ 26]);
836 RE(y[23]) = RE(x[ 30]); IM(y[23]) = IM(x[ 30]);
837 RE(y[24]) = RE(x[ 3]); IM(y[24]) = IM(x[ 3]);
838 RE(y[25]) = RE(x[ 7]); IM(y[25]) = IM(x[ 7]);
839 RE(y[26]) = RE(x[ 11]); IM(y[26]) = IM(x[ 11]);
840 RE(y[27]) = RE(x[ 15]); IM(y[27]) = IM(x[ 15]);
841 RE(y[28]) = RE(x[ 19]); IM(y[28]) = IM(x[ 19]);
842 RE(y[29]) = RE(x[ 23]); IM(y[29]) = IM(x[ 23]);
843 RE(y[30]) = RE(x[ 27]); IM(y[30]) = IM(x[ 27]);
844 RE(y[31]) = RE(x[ 31]); IM(y[31]) = IM(x[ 31]);
845}
846
847
848/*******************************************************************************
8494 x 8 matrix transpose
850*******************************************************************************/
851void transpose8x8(complex_t *x, complex_t* y)
852{
853 RE(y[ 0]) = RE(x[ 0]); IM(y[ 0]) = IM(x[ 0]);
854 RE(y[ 1]) = RE(x[ 8]); IM(y[ 1]) = IM(x[ 8]);
855 RE(y[ 2]) = RE(x[16]); IM(y[ 2]) = IM(x[16]);
856 RE(y[ 3]) = RE(x[24]); IM(y[ 3]) = IM(x[24]);
857 RE(y[ 4]) = RE(x[32]); IM(y[ 4]) = IM(x[32]);
858 RE(y[ 5]) = RE(x[40]); IM(y[ 5]) = IM(x[40]);
859 RE(y[ 6]) = RE(x[48]); IM(y[ 6]) = IM(x[48]);
860 RE(y[ 7]) = RE(x[56]); IM(y[ 7]) = IM(x[56]);
861 RE(y[ 8]) = RE(x[ 1]); IM(y[ 8]) = IM(x[ 1]);
862 RE(y[ 9]) = RE(x[ 9]); IM(y[ 9]) = IM(x[ 9]);
863 RE(y[ 10]) = RE(x[17]); IM(y[ 10]) = IM(x[17]);
864 RE(y[ 11]) = RE(x[25]); IM(y[ 11]) = IM(x[25]);
865 RE(y[ 12]) = RE(x[33]); IM(y[ 12]) = IM(x[33]);
866 RE(y[ 13]) = RE(x[41]); IM(y[ 13]) = IM(x[41]);
867 RE(y[ 14]) = RE(x[49]); IM(y[ 14]) = IM(x[49]);
868 RE(y[ 15]) = RE(x[57]); IM(y[ 15]) = IM(x[57]);
869 RE(y[ 16]) = RE(x[ 2]); IM(y[ 16]) = IM(x[ 2]);
870 RE(y[ 17]) = RE(x[10]); IM(y[ 17]) = IM(x[10]);
871 RE(y[ 18]) = RE(x[18]); IM(y[ 18]) = IM(x[18]);
872 RE(y[ 19]) = RE(x[26]); IM(y[ 19]) = IM(x[26]);
873 RE(y[ 20]) = RE(x[34]); IM(y[ 20]) = IM(x[34]);
874 RE(y[ 21]) = RE(x[42]); IM(y[ 21]) = IM(x[42]);
875 RE(y[ 22]) = RE(x[50]); IM(y[ 22]) = IM(x[50]);
876 RE(y[ 23]) = RE(x[58]); IM(y[ 23]) = IM(x[58]);
877 RE(y[ 24]) = RE(x[ 3]); IM(y[ 24]) = IM(x[ 3]);
878 RE(y[ 25]) = RE(x[11]); IM(y[ 25]) = IM(x[11]);
879 RE(y[ 26]) = RE(x[19]); IM(y[ 26]) = IM(x[19]);
880 RE(y[ 27]) = RE(x[27]); IM(y[ 27]) = IM(x[27]);
881 RE(y[ 28]) = RE(x[35]); IM(y[ 28]) = IM(x[35]);
882 RE(y[ 29]) = RE(x[43]); IM(y[ 29]) = IM(x[43]);
883 RE(y[ 30]) = RE(x[51]); IM(y[ 30]) = IM(x[51]);
884 RE(y[ 31]) = RE(x[59]); IM(y[ 31]) = IM(x[59]);
885 RE(y[ 32]) = RE(x[ 4]); IM(y[ 32]) = IM(x[ 4]);
886 RE(y[ 33]) = RE(x[12]); IM(y[ 33]) = IM(x[12]);
887 RE(y[ 34]) = RE(x[20]); IM(y[ 34]) = IM(x[20]);
888 RE(y[ 35]) = RE(x[28]); IM(y[ 35]) = IM(x[28]);
889 RE(y[ 36]) = RE(x[36]); IM(y[ 36]) = IM(x[36]);
890 RE(y[ 37]) = RE(x[44]); IM(y[ 37]) = IM(x[44]);
891 RE(y[ 38]) = RE(x[52]); IM(y[ 38]) = IM(x[52]);
892 RE(y[ 39]) = RE(x[60]); IM(y[ 39]) = IM(x[60]);
893 RE(y[ 40]) = RE(x[ 5]); IM(y[ 40]) = IM(x[ 5]);
894 RE(y[ 41]) = RE(x[13]); IM(y[ 41]) = IM(x[13]);
895 RE(y[ 42]) = RE(x[21]); IM(y[ 42]) = IM(x[21]);
896 RE(y[ 43]) = RE(x[29]); IM(y[ 43]) = IM(x[29]);
897 RE(y[ 44]) = RE(x[37]); IM(y[ 44]) = IM(x[37]);
898 RE(y[ 45]) = RE(x[45]); IM(y[ 45]) = IM(x[45]);
899 RE(y[ 46]) = RE(x[53]); IM(y[ 46]) = IM(x[53]);
900 RE(y[ 47]) = RE(x[61]); IM(y[ 47]) = IM(x[61]);
901 RE(y[ 48]) = RE(x[ 6]); IM(y[ 48]) = IM(x[ 6]);
902 RE(y[ 49]) = RE(x[14]); IM(y[ 49]) = IM(x[14]);
903 RE(y[ 50]) = RE(x[22]); IM(y[ 50]) = IM(x[22]);
904 RE(y[ 51]) = RE(x[30]); IM(y[ 51]) = IM(x[30]);
905 RE(y[ 52]) = RE(x[38]); IM(y[ 52]) = IM(x[38]);
906 RE(y[ 53]) = RE(x[46]); IM(y[ 53]) = IM(x[46]);
907 RE(y[ 54]) = RE(x[54]); IM(y[ 54]) = IM(x[54]);
908 RE(y[ 55]) = RE(x[62]); IM(y[ 55]) = IM(x[62]);
909 RE(y[ 56]) = RE(x[ 7]); IM(y[ 56]) = IM(x[ 7]);
910 RE(y[ 57]) = RE(x[15]); IM(y[ 57]) = IM(x[15]);
911 RE(y[ 58]) = RE(x[23]); IM(y[ 58]) = IM(x[23]);
912 RE(y[ 59]) = RE(x[31]); IM(y[ 59]) = IM(x[31]);
913 RE(y[ 60]) = RE(x[39]); IM(y[ 60]) = IM(x[39]);
914 RE(y[ 61]) = RE(x[47]); IM(y[ 61]) = IM(x[47]);
915 RE(y[ 62]) = RE(x[55]); IM(y[ 62]) = IM(x[55]);
916 RE(y[ 63]) = RE(x[63]); IM(y[ 63]) = IM(x[63]);
917}
918
919
920
921/*******************************************************************************
92216 x 16 matrix transpose
923*******************************************************************************/
924void transpose16x16(complex_t* x, complex_t* y)
925{
926 RE(y[ 0]) = RE(x[ 0]); IM(y[ 0]) = IM(x[ 0]);
927 RE(y[ 1]) = RE(x[ 16]); IM(y[ 1]) = IM(x[ 16]);
928 RE(y[ 2]) = RE(x[ 32]); IM(y[ 2]) = IM(x[ 32]);
929 RE(y[ 3]) = RE(x[ 48]); IM(y[ 3]) = IM(x[ 48]);
930 RE(y[ 4]) = RE(x[ 64]); IM(y[ 4]) = IM(x[ 64]);
931 RE(y[ 5]) = RE(x[ 80]); IM(y[ 5]) = IM(x[ 80]);
932 RE(y[ 6]) = RE(x[ 96]); IM(y[ 6]) = IM(x[ 96]);
933 RE(y[ 7]) = RE(x[112]); IM(y[ 7]) = IM(x[112]);
934 RE(y[ 8]) = RE(x[128]); IM(y[ 8]) = IM(x[128]);
935 RE(y[ 9]) = RE(x[144]); IM(y[ 9]) = IM(x[144]);
936 RE(y[ 10]) = RE(x[160]); IM(y[ 10]) = IM(x[160]);
937 RE(y[ 11]) = RE(x[176]); IM(y[ 11]) = IM(x[176]);
938 RE(y[ 12]) = RE(x[192]); IM(y[ 12]) = IM(x[192]);
939 RE(y[ 13]) = RE(x[208]); IM(y[ 13]) = IM(x[208]);
940 RE(y[ 14]) = RE(x[224]); IM(y[ 14]) = IM(x[224]);
941 RE(y[ 15]) = RE(x[240]); IM(y[ 15]) = IM(x[240]);
942 RE(y[ 16]) = RE(x[ 1]); IM(y[ 16]) = IM(x[ 1]);
943 RE(y[ 17]) = RE(x[ 17]); IM(y[ 17]) = IM(x[ 17]);
944 RE(y[ 18]) = RE(x[ 33]); IM(y[ 18]) = IM(x[ 33]);
945 RE(y[ 19]) = RE(x[ 49]); IM(y[ 19]) = IM(x[ 49]);
946 RE(y[ 20]) = RE(x[ 65]); IM(y[ 20]) = IM(x[ 65]);
947 RE(y[ 21]) = RE(x[ 81]); IM(y[ 21]) = IM(x[ 81]);
948 RE(y[ 22]) = RE(x[ 97]); IM(y[ 22]) = IM(x[ 97]);
949 RE(y[ 23]) = RE(x[113]); IM(y[ 23]) = IM(x[113]);
950 RE(y[ 24]) = RE(x[129]); IM(y[ 24]) = IM(x[129]);
951 RE(y[ 25]) = RE(x[145]); IM(y[ 25]) = IM(x[145]);
952 RE(y[ 26]) = RE(x[161]); IM(y[ 26]) = IM(x[161]);
953 RE(y[ 27]) = RE(x[177]); IM(y[ 27]) = IM(x[177]);
954 RE(y[ 28]) = RE(x[193]); IM(y[ 28]) = IM(x[193]);
955 RE(y[ 29]) = RE(x[209]); IM(y[ 29]) = IM(x[209]);
956 RE(y[ 30]) = RE(x[225]); IM(y[ 30]) = IM(x[225]);
957 RE(y[ 31]) = RE(x[241]); IM(y[ 31]) = IM(x[241]);
958 RE(y[ 32]) = RE(x[ 2]); IM(y[ 32]) = IM(x[ 2]);
959 RE(y[ 33]) = RE(x[ 18]); IM(y[ 33]) = IM(x[ 18]);
960 RE(y[ 34]) = RE(x[ 34]); IM(y[ 34]) = IM(x[ 34]);
961 RE(y[ 35]) = RE(x[ 50]); IM(y[ 35]) = IM(x[ 50]);
962 RE(y[ 36]) = RE(x[ 66]); IM(y[ 36]) = IM(x[ 66]);
963 RE(y[ 37]) = RE(x[ 82]); IM(y[ 37]) = IM(x[ 82]);
964 RE(y[ 38]) = RE(x[ 98]); IM(y[ 38]) = IM(x[ 98]);
965 RE(y[ 39]) = RE(x[114]); IM(y[ 39]) = IM(x[114]);
966 RE(y[ 40]) = RE(x[130]); IM(y[ 40]) = IM(x[130]);
967 RE(y[ 41]) = RE(x[146]); IM(y[ 41]) = IM(x[146]);
968 RE(y[ 42]) = RE(x[162]); IM(y[ 42]) = IM(x[162]);
969 RE(y[ 43]) = RE(x[178]); IM(y[ 43]) = IM(x[178]);
970 RE(y[ 44]) = RE(x[194]); IM(y[ 44]) = IM(x[194]);
971 RE(y[ 45]) = RE(x[210]); IM(y[ 45]) = IM(x[210]);
972 RE(y[ 46]) = RE(x[226]); IM(y[ 46]) = IM(x[226]);
973 RE(y[ 47]) = RE(x[242]); IM(y[ 47]) = IM(x[242]);
974 RE(y[ 48]) = RE(x[ 3]); IM(y[ 48]) = IM(x[ 3]);
975 RE(y[ 49]) = RE(x[ 19]); IM(y[ 49]) = IM(x[ 19]);
976 RE(y[ 50]) = RE(x[ 35]); IM(y[ 50]) = IM(x[ 35]);
977 RE(y[ 51]) = RE(x[ 51]); IM(y[ 51]) = IM(x[ 51]);
978 RE(y[ 52]) = RE(x[ 67]); IM(y[ 52]) = IM(x[ 67]);
979 RE(y[ 53]) = RE(x[ 83]); IM(y[ 53]) = IM(x[ 83]);
980 RE(y[ 54]) = RE(x[ 99]); IM(y[ 54]) = IM(x[ 99]);
981 RE(y[ 55]) = RE(x[115]); IM(y[ 55]) = IM(x[115]);
982 RE(y[ 56]) = RE(x[131]); IM(y[ 56]) = IM(x[131]);
983 RE(y[ 57]) = RE(x[147]); IM(y[ 57]) = IM(x[147]);
984 RE(y[ 58]) = RE(x[163]); IM(y[ 58]) = IM(x[163]);
985 RE(y[ 59]) = RE(x[179]); IM(y[ 59]) = IM(x[179]);
986 RE(y[ 60]) = RE(x[195]); IM(y[ 60]) = IM(x[195]);
987 RE(y[ 61]) = RE(x[211]); IM(y[ 61]) = IM(x[211]);
988 RE(y[ 62]) = RE(x[227]); IM(y[ 62]) = IM(x[227]);
989 RE(y[ 63]) = RE(x[243]); IM(y[ 63]) = IM(x[243]);
990 RE(y[ 64]) = RE(x[ 4]); IM(y[ 64]) = IM(x[ 4]);
991 RE(y[ 65]) = RE(x[ 20]); IM(y[ 65]) = IM(x[ 20]);
992 RE(y[ 66]) = RE(x[ 36]); IM(y[ 66]) = IM(x[ 36]);
993 RE(y[ 67]) = RE(x[ 52]); IM(y[ 67]) = IM(x[ 52]);
994 RE(y[ 68]) = RE(x[ 68]); IM(y[ 68]) = IM(x[ 68]);
995 RE(y[ 69]) = RE(x[ 84]); IM(y[ 69]) = IM(x[ 84]);
996 RE(y[ 70]) = RE(x[100]); IM(y[ 70]) = IM(x[100]);
997 RE(y[ 71]) = RE(x[116]); IM(y[ 71]) = IM(x[116]);
998 RE(y[ 72]) = RE(x[132]); IM(y[ 72]) = IM(x[132]);
999 RE(y[ 73]) = RE(x[148]); IM(y[ 73]) = IM(x[148]);
1000 RE(y[ 74]) = RE(x[164]); IM(y[ 74]) = IM(x[164]);
1001 RE(y[ 75]) = RE(x[180]); IM(y[ 75]) = IM(x[180]);
1002 RE(y[ 76]) = RE(x[196]); IM(y[ 76]) = IM(x[196]);
1003 RE(y[ 77]) = RE(x[212]); IM(y[ 77]) = IM(x[212]);
1004 RE(y[ 78]) = RE(x[228]); IM(y[ 78]) = IM(x[228]);
1005 RE(y[ 79]) = RE(x[244]); IM(y[ 79]) = IM(x[244]);
1006 RE(y[ 80]) = RE(x[ 5]); IM(y[ 80]) = IM(x[ 5]);
1007 RE(y[ 81]) = RE(x[ 21]); IM(y[ 81]) = IM(x[ 21]);
1008 RE(y[ 82]) = RE(x[ 37]); IM(y[ 82]) = IM(x[ 37]);
1009 RE(y[ 83]) = RE(x[ 53]); IM(y[ 83]) = IM(x[ 53]);
1010 RE(y[ 84]) = RE(x[ 69]); IM(y[ 84]) = IM(x[ 69]);
1011 RE(y[ 85]) = RE(x[ 85]); IM(y[ 85]) = IM(x[ 85]);
1012 RE(y[ 86]) = RE(x[101]); IM(y[ 86]) = IM(x[101]);
1013 RE(y[ 87]) = RE(x[117]); IM(y[ 87]) = IM(x[117]);
1014 RE(y[ 88]) = RE(x[133]); IM(y[ 88]) = IM(x[133]);
1015 RE(y[ 89]) = RE(x[149]); IM(y[ 89]) = IM(x[149]);
1016 RE(y[ 90]) = RE(x[165]); IM(y[ 90]) = IM(x[165]);
1017 RE(y[ 91]) = RE(x[181]); IM(y[ 91]) = IM(x[181]);
1018 RE(y[ 92]) = RE(x[197]); IM(y[ 92]) = IM(x[197]);
1019 RE(y[ 93]) = RE(x[213]); IM(y[ 93]) = IM(x[213]);
1020 RE(y[ 94]) = RE(x[229]); IM(y[ 94]) = IM(x[229]);
1021 RE(y[ 95]) = RE(x[245]); IM(y[ 95]) = IM(x[245]);
1022 RE(y[ 96]) = RE(x[ 6]); IM(y[ 96]) = IM(x[ 6]);
1023 RE(y[ 97]) = RE(x[ 22]); IM(y[ 97]) = IM(x[ 22]);
1024 RE(y[ 98]) = RE(x[ 38]); IM(y[ 98]) = IM(x[ 38]);
1025 RE(y[ 99]) = RE(x[ 54]); IM(y[ 99]) = IM(x[ 54]);
1026 RE(y[100]) = RE(x[ 70]); IM(y[100]) = IM(x[ 70]);
1027 RE(y[101]) = RE(x[ 86]); IM(y[101]) = IM(x[ 86]);
1028 RE(y[102]) = RE(x[102]); IM(y[102]) = IM(x[102]);
1029 RE(y[103]) = RE(x[118]); IM(y[103]) = IM(x[118]);
1030 RE(y[104]) = RE(x[134]); IM(y[104]) = IM(x[134]);
1031 RE(y[105]) = RE(x[150]); IM(y[105]) = IM(x[150]);
1032 RE(y[106]) = RE(x[166]); IM(y[106]) = IM(x[166]);
1033 RE(y[107]) = RE(x[182]); IM(y[107]) = IM(x[182]);
1034 RE(y[108]) = RE(x[198]); IM(y[108]) = IM(x[198]);
1035 RE(y[109]) = RE(x[214]); IM(y[109]) = IM(x[214]);
1036 RE(y[110]) = RE(x[230]); IM(y[110]) = IM(x[230]);
1037 RE(y[111]) = RE(x[246]); IM(y[111]) = IM(x[246]);
1038 RE(y[112]) = RE(x[ 7]); IM(y[112]) = IM(x[ 7]);
1039 RE(y[113]) = RE(x[ 23]); IM(y[113]) = IM(x[ 23]);
1040 RE(y[114]) = RE(x[ 39]); IM(y[114]) = IM(x[ 39]);
1041 RE(y[115]) = RE(x[ 55]); IM(y[115]) = IM(x[ 55]);
1042 RE(y[116]) = RE(x[ 71]); IM(y[116]) = IM(x[ 71]);
1043 RE(y[117]) = RE(x[ 87]); IM(y[117]) = IM(x[ 87]);
1044 RE(y[118]) = RE(x[103]); IM(y[118]) = IM(x[103]);
1045 RE(y[119]) = RE(x[119]); IM(y[119]) = IM(x[119]);
1046 RE(y[120]) = RE(x[135]); IM(y[120]) = IM(x[135]);
1047 RE(y[121]) = RE(x[151]); IM(y[121]) = IM(x[151]);
1048 RE(y[122]) = RE(x[167]); IM(y[122]) = IM(x[167]);
1049 RE(y[123]) = RE(x[183]); IM(y[123]) = IM(x[183]);
1050 RE(y[124]) = RE(x[199]); IM(y[124]) = IM(x[199]);
1051 RE(y[125]) = RE(x[215]); IM(y[125]) = IM(x[215]);
1052 RE(y[126]) = RE(x[231]); IM(y[126]) = IM(x[231]);
1053 RE(y[127]) = RE(x[247]); IM(y[127]) = IM(x[247]);
1054 RE(y[128]) = RE(x[ 8]); IM(y[128]) = IM(x[ 8]);
1055 RE(y[129]) = RE(x[ 24]); IM(y[129]) = IM(x[ 24]);
1056 RE(y[130]) = RE(x[ 40]); IM(y[130]) = IM(x[ 40]);
1057 RE(y[131]) = RE(x[ 56]); IM(y[131]) = IM(x[ 56]);
1058 RE(y[132]) = RE(x[ 72]); IM(y[132]) = IM(x[ 72]);
1059 RE(y[133]) = RE(x[ 88]); IM(y[133]) = IM(x[ 88]);
1060 RE(y[134]) = RE(x[104]); IM(y[134]) = IM(x[104]);
1061 RE(y[135]) = RE(x[120]); IM(y[135]) = IM(x[120]);
1062 RE(y[136]) = RE(x[136]); IM(y[136]) = IM(x[136]);
1063 RE(y[137]) = RE(x[152]); IM(y[137]) = IM(x[152]);
1064 RE(y[138]) = RE(x[168]); IM(y[138]) = IM(x[168]);
1065 RE(y[139]) = RE(x[184]); IM(y[139]) = IM(x[184]);
1066 RE(y[140]) = RE(x[200]); IM(y[140]) = IM(x[200]);
1067 RE(y[141]) = RE(x[216]); IM(y[141]) = IM(x[216]);
1068 RE(y[142]) = RE(x[232]); IM(y[142]) = IM(x[232]);
1069 RE(y[143]) = RE(x[248]); IM(y[143]) = IM(x[248]);
1070 RE(y[144]) = RE(x[ 9]); IM(y[144]) = IM(x[ 9]);
1071 RE(y[145]) = RE(x[ 25]); IM(y[145]) = IM(x[ 25]);
1072 RE(y[146]) = RE(x[ 41]); IM(y[146]) = IM(x[ 41]);
1073 RE(y[147]) = RE(x[ 57]); IM(y[147]) = IM(x[ 57]);
1074 RE(y[148]) = RE(x[ 73]); IM(y[148]) = IM(x[ 73]);
1075 RE(y[149]) = RE(x[ 89]); IM(y[149]) = IM(x[ 89]);
1076 RE(y[150]) = RE(x[105]); IM(y[150]) = IM(x[105]);
1077 RE(y[151]) = RE(x[121]); IM(y[151]) = IM(x[121]);
1078 RE(y[152]) = RE(x[137]); IM(y[152]) = IM(x[137]);
1079 RE(y[153]) = RE(x[153]); IM(y[153]) = IM(x[153]);
1080 RE(y[154]) = RE(x[169]); IM(y[154]) = IM(x[169]);
1081 RE(y[155]) = RE(x[185]); IM(y[155]) = IM(x[185]);
1082 RE(y[156]) = RE(x[201]); IM(y[156]) = IM(x[201]);
1083 RE(y[157]) = RE(x[217]); IM(y[157]) = IM(x[217]);
1084 RE(y[158]) = RE(x[233]); IM(y[158]) = IM(x[233]);
1085 RE(y[159]) = RE(x[249]); IM(y[159]) = IM(x[249]);
1086 RE(y[160]) = RE(x[ 10]); IM(y[160]) = IM(x[ 10]);
1087 RE(y[161]) = RE(x[ 26]); IM(y[161]) = IM(x[ 26]);
1088 RE(y[162]) = RE(x[ 42]); IM(y[162]) = IM(x[ 42]);
1089 RE(y[163]) = RE(x[ 58]); IM(y[163]) = IM(x[ 58]);
1090 RE(y[164]) = RE(x[ 74]); IM(y[164]) = IM(x[ 74]);
1091 RE(y[165]) = RE(x[ 90]); IM(y[165]) = IM(x[ 90]);
1092 RE(y[166]) = RE(x[106]); IM(y[166]) = IM(x[106]);
1093 RE(y[167]) = RE(x[122]); IM(y[167]) = IM(x[122]);
1094 RE(y[168]) = RE(x[138]); IM(y[168]) = IM(x[138]);
1095 RE(y[169]) = RE(x[154]); IM(y[169]) = IM(x[154]);
1096 RE(y[170]) = RE(x[170]); IM(y[170]) = IM(x[170]);
1097 RE(y[171]) = RE(x[186]); IM(y[171]) = IM(x[186]);
1098 RE(y[172]) = RE(x[202]); IM(y[172]) = IM(x[202]);
1099 RE(y[173]) = RE(x[218]); IM(y[173]) = IM(x[218]);
1100 RE(y[174]) = RE(x[234]); IM(y[174]) = IM(x[234]);
1101 RE(y[175]) = RE(x[250]); IM(y[175]) = IM(x[250]);
1102 RE(y[176]) = RE(x[ 11]); IM(y[176]) = IM(x[ 11]);
1103 RE(y[177]) = RE(x[ 27]); IM(y[177]) = IM(x[ 27]);
1104 RE(y[178]) = RE(x[ 43]); IM(y[178]) = IM(x[ 43]);
1105 RE(y[179]) = RE(x[ 59]); IM(y[179]) = IM(x[ 59]);
1106 RE(y[180]) = RE(x[ 75]); IM(y[180]) = IM(x[ 75]);
1107 RE(y[181]) = RE(x[ 91]); IM(y[181]) = IM(x[ 91]);
1108 RE(y[182]) = RE(x[107]); IM(y[182]) = IM(x[107]);
1109 RE(y[183]) = RE(x[123]); IM(y[183]) = IM(x[123]);
1110 RE(y[184]) = RE(x[139]); IM(y[184]) = IM(x[139]);
1111 RE(y[185]) = RE(x[155]); IM(y[185]) = IM(x[155]);
1112 RE(y[186]) = RE(x[171]); IM(y[186]) = IM(x[171]);
1113 RE(y[187]) = RE(x[187]); IM(y[187]) = IM(x[187]);
1114 RE(y[188]) = RE(x[203]); IM(y[188]) = IM(x[203]);
1115 RE(y[189]) = RE(x[219]); IM(y[189]) = IM(x[219]);
1116 RE(y[190]) = RE(x[235]); IM(y[190]) = IM(x[235]);
1117 RE(y[191]) = RE(x[251]); IM(y[191]) = IM(x[251]);
1118 RE(y[192]) = RE(x[ 12]); IM(y[192]) = IM(x[ 12]);
1119 RE(y[193]) = RE(x[ 28]); IM(y[193]) = IM(x[ 28]);
1120 RE(y[194]) = RE(x[ 44]); IM(y[194]) = IM(x[ 44]);
1121 RE(y[195]) = RE(x[ 60]); IM(y[195]) = IM(x[ 60]);
1122 RE(y[196]) = RE(x[ 76]); IM(y[196]) = IM(x[ 76]);
1123 RE(y[197]) = RE(x[ 92]); IM(y[197]) = IM(x[ 92]);
1124 RE(y[198]) = RE(x[108]); IM(y[198]) = IM(x[108]);
1125 RE(y[199]) = RE(x[124]); IM(y[199]) = IM(x[124]);
1126 RE(y[200]) = RE(x[140]); IM(y[200]) = IM(x[140]);
1127 RE(y[201]) = RE(x[156]); IM(y[201]) = IM(x[156]);
1128 RE(y[202]) = RE(x[172]); IM(y[202]) = IM(x[172]);
1129 RE(y[203]) = RE(x[188]); IM(y[203]) = IM(x[188]);
1130 RE(y[204]) = RE(x[204]); IM(y[204]) = IM(x[204]);
1131 RE(y[205]) = RE(x[220]); IM(y[205]) = IM(x[220]);
1132 RE(y[206]) = RE(x[236]); IM(y[206]) = IM(x[236]);
1133 RE(y[207]) = RE(x[252]); IM(y[207]) = IM(x[252]);
1134 RE(y[208]) = RE(x[ 13]); IM(y[208]) = IM(x[ 13]);
1135 RE(y[209]) = RE(x[ 29]); IM(y[209]) = IM(x[ 29]);
1136 RE(y[210]) = RE(x[ 45]); IM(y[210]) = IM(x[ 45]);
1137 RE(y[211]) = RE(x[ 61]); IM(y[211]) = IM(x[ 61]);
1138 RE(y[212]) = RE(x[ 77]); IM(y[212]) = IM(x[ 77]);
1139 RE(y[213]) = RE(x[ 93]); IM(y[213]) = IM(x[ 93]);
1140 RE(y[214]) = RE(x[109]); IM(y[214]) = IM(x[109]);
1141 RE(y[215]) = RE(x[125]); IM(y[215]) = IM(x[125]);
1142 RE(y[216]) = RE(x[141]); IM(y[216]) = IM(x[141]);
1143 RE(y[217]) = RE(x[157]); IM(y[217]) = IM(x[157]);
1144 RE(y[218]) = RE(x[173]); IM(y[218]) = IM(x[173]);
1145 RE(y[219]) = RE(x[189]); IM(y[219]) = IM(x[189]);
1146 RE(y[220]) = RE(x[205]); IM(y[220]) = IM(x[205]);
1147 RE(y[221]) = RE(x[221]); IM(y[221]) = IM(x[221]);
1148 RE(y[222]) = RE(x[237]); IM(y[222]) = IM(x[237]);
1149 RE(y[223]) = RE(x[253]); IM(y[223]) = IM(x[253]);
1150 RE(y[224]) = RE(x[ 14]); IM(y[224]) = IM(x[ 14]);
1151 RE(y[225]) = RE(x[ 30]); IM(y[225]) = IM(x[ 30]);
1152 RE(y[226]) = RE(x[ 46]); IM(y[226]) = IM(x[ 46]);
1153 RE(y[227]) = RE(x[ 62]); IM(y[227]) = IM(x[ 62]);
1154 RE(y[228]) = RE(x[ 78]); IM(y[228]) = IM(x[ 78]);
1155 RE(y[229]) = RE(x[ 94]); IM(y[229]) = IM(x[ 94]);
1156 RE(y[230]) = RE(x[110]); IM(y[230]) = IM(x[110]);
1157 RE(y[231]) = RE(x[126]); IM(y[231]) = IM(x[126]);
1158 RE(y[232]) = RE(x[142]); IM(y[232]) = IM(x[142]);
1159 RE(y[233]) = RE(x[158]); IM(y[233]) = IM(x[158]);
1160 RE(y[234]) = RE(x[174]); IM(y[234]) = IM(x[174]);
1161 RE(y[235]) = RE(x[190]); IM(y[235]) = IM(x[190]);
1162 RE(y[236]) = RE(x[206]); IM(y[236]) = IM(x[206]);
1163 RE(y[237]) = RE(x[222]); IM(y[237]) = IM(x[222]);
1164 RE(y[238]) = RE(x[238]); IM(y[238]) = IM(x[238]);
1165 RE(y[239]) = RE(x[254]); IM(y[239]) = IM(x[254]);
1166 RE(y[240]) = RE(x[ 15]); IM(y[240]) = IM(x[ 15]);
1167 RE(y[241]) = RE(x[ 31]); IM(y[241]) = IM(x[ 31]);
1168 RE(y[242]) = RE(x[ 47]); IM(y[242]) = IM(x[ 47]);
1169 RE(y[243]) = RE(x[ 63]); IM(y[243]) = IM(x[ 63]);
1170 RE(y[244]) = RE(x[ 79]); IM(y[244]) = IM(x[ 79]);
1171 RE(y[245]) = RE(x[ 95]); IM(y[245]) = IM(x[ 95]);
1172 RE(y[246]) = RE(x[111]); IM(y[246]) = IM(x[111]);
1173 RE(y[247]) = RE(x[127]); IM(y[247]) = IM(x[127]);
1174 RE(y[248]) = RE(x[143]); IM(y[248]) = IM(x[143]);
1175 RE(y[249]) = RE(x[159]); IM(y[249]) = IM(x[159]);
1176 RE(y[250]) = RE(x[175]); IM(y[250]) = IM(x[175]);
1177 RE(y[251]) = RE(x[191]); IM(y[251]) = IM(x[191]);
1178 RE(y[252]) = RE(x[207]); IM(y[252]) = IM(x[207]);
1179 RE(y[253]) = RE(x[223]); IM(y[253]) = IM(x[223]);
1180 RE(y[254]) = RE(x[239]); IM(y[254]) = IM(x[239]);
1181 RE(y[255]) = RE(x[255]); IM(y[255]) = IM(x[255]);
1182}
1183
1184
sum
int DSPL_API sum(double *x, int n, double *s)
Definition: sum.c:43
RE
#define RE(x)
Macro sets real part of the complex number.
Definition: dspl.h:420
complex_t
double complex_t[2]
Complex data type.
Definition: dspl.h:86
IM
#define IM(x)
Macro sets imaginary part of the complex number.
Definition: dspl.h:478
Generated on Wed Jan 5 2022 12:44:36 for libdspl-2.0 by   doxygen 1.9.2
libdspl-2.0 documentation is licensed under a Creative Commons Attribution-ShareAlike 4.0 International License.