1 
2 extern void cpass(complex *,const complex *,unsigned int);
3 extern void cpassbig(complex *,const complex *,unsigned int);
4 extern void upass(complex *,const complex *,unsigned int);
5 extern void upassbig(complex *,const complex *,unsigned int);
6 extern void rpass(real *,const complex *,unsigned int);
7 extern void rpassbig(real *,const complex *,unsigned int);
8 extern void vpass(real *,const complex *,unsigned int);
9 extern void vpassbig(real *,const complex *,unsigned int);
10 
11 extern const complex d16[];
12 extern const complex d32[];
13 extern const complex d64[];
14 extern const complex d128[];
15 extern const complex d256[];
16 extern const complex d512[];
17 extern const complex d1024[];
18 extern const complex d2048[];
19 extern const complex d4096[];
20 extern const complex d8192[];
21 extern const complex d16384[];
22 extern const complex d32768[];
23 extern const complex d65536[];
24 
25 #define sqrthalf (d16[1].re)
26 
27 #define VOL *(volatile real *)&
28 
29 /* assumes sizeof(complex) is a multiple of 4 */
30 #define A0 (a)
31 #define A1 ((complex *) (n * (sizeof(complex) / 4) + (char *) a))
32 #define A2 (b)
33 #define A3 ((complex *) (n * (sizeof(complex) / 4) + (char *) b))
34 #define X0 (x)
35 #define X1 ((complex *) (n * (sizeof(complex) / 4) + (char *) x))
36 #define X2 (y)
37 #define X3 ((complex *) (n * (sizeof(complex) / 4) + (char *) y))
38 
39 #define TRANSFORM(a0,a1,a2,a3,wre,wim) { \
40   t1 = a0.re + a2.re; \
41   t2 = a0.im + a2.im; \
42   t3 = a1.re + a3.re; \
43   t5 = VOL a0.re; \
44   a0.re = t1; \
45   t6 = VOL a0.im; \
46   t5 -= a2.re; \
47   t4 = VOL a1.im; \
48   t4 += a3.im; \
49   a0.im = t2; \
50   t8 = VOL a1.im; \
51   a1.im = t4; \
52   t7 = VOL a1.re; \
53   t8 -= a3.im; \
54   t6 -= a2.im; \
55   t1 = t5 - t8; \
56   t5 += t8; \
57   a1.re = t3; \
58   t3 = t1 * wim; \
59   t7 -= a3.re; \
60   t1 *= wre; \
61   t8 = t5; \
62   t5 *= wre; \
63   t2 = t6; \
64   t2 += t7; \
65   t8 *= wim; \
66   t6 -= t7; \
67   t4 = t2; \
68   t4 *= wim; \
69   t7 = t6; \
70   t2 *= wre; \
71   t1 -= t4; \
72   t6 *= wre; \
73   t2 += t3; \
74   t7 *= wim; \
75   t6 -= t8; \
76   a2.re = t1; \
77   t5 += t7; \
78   a3.im = t6; \
79   a2.im = t2; \
80   a3.re = t5; \
81   }
82 
83 #define TRANSFORMHALF(a0,a1,a2,a3) { \
84   t1 = a0.re + a2.re; \
85   t2 = a0.im + a2.im; \
86   t3 = a1.re + a3.re; \
87   t5 = VOL a0.re; \
88   a0.re = t1; \
89   t6 = VOL a0.im; \
90   t5 -= a2.re; \
91   t4 = VOL a1.im; \
92   t4 += a3.im; \
93   a0.im = t2; \
94   t8 = VOL a1.im; \
95   a1.im = t4; \
96   t7 = a1.re - a3.re; \
97   t6 -= a2.im; \
98   t1 = t5; \
99   t8 -= a3.im; \
100   a1.re = t3; \
101   t2 = t6; \
102   t6 -= t7; \
103   t5 += t8; \
104   t2 += t7; \
105   t4 = t6; \
106   t1 -= t8; \
107   t4 += t5; \
108   t6 -= t5; \
109   t7 = t1; \
110   t4 *= sqrthalf; \
111   t7 -= t2; \
112   t6 *= sqrthalf; \
113   t2 += t1; \
114   t7 *= sqrthalf; \
115   a3.re = t4; \
116   t2 *= sqrthalf; \
117   a3.im = t6; \
118   a2.re = t7; \
119   a2.im = t2; \
120   }
121 
122 #define TRANSFORMZERO(a0,a1,a2,a3) { \
123   t5 = a0.re + a2.re; \
124   t6 = a0.im + a2.im; \
125   t7 = a1.re + a3.re; \
126   t8 = a1.im + a3.im; \
127   t1 = a0.re - a2.re; \
128   a0.re = t5; \
129   t2 = a0.im - a2.im; \
130   a0.im = t6; \
131   t3 = a1.re - a3.re; \
132   a1.re = t7; \
133   t4 = a1.im - a3.im; \
134   a1.im = t8; \
135   t6 = t1 - t4; \
136   t7 = t2 + t3; \
137   t1 += t4; \
138   t2 -= t3; \
139   a2.re = t6; \
140   a3.im = t2; \
141   a3.re = t1; \
142   a2.im = t7; \
143   }
144 
145 #define UNTRANSFORM(a0,a1,a2,a3,wre,wim) { \
146   t1 = a2.re * wre; \
147   t3 = a2.im * wim; \
148   t5 = a3.re * wre; \
149   t7 = a3.im * wim; \
150   t2 = a2.im * wre; \
151   t4 = a2.re * wim; \
152   t6 = a3.im * wre; \
153   t8 = a3.re * wim; \
154   t1 += t3; \
155   t5 -= t7; \
156   t2 -= t4; \
157   t6 += t8; \
158   t3 = t5 + t1; \
159   t5 -= t1; \
160   t4 = t2 - t6; \
161   t6 += t2; \
162   t1 = a0.re - t3; \
163   t3 += a0.re; \
164   t2 = a0.im - t6; \
165   t6 += a0.im; \
166   t8 = a1.re - t4; \
167   t4 += a1.re; \
168   t7 = a1.im - t5; \
169   t5 += a1.im; \
170   a2.re = t1; \
171   a1.re = t4; \
172   a0.im = t6; \
173   a0.re = t3; \
174   a1.im = t5; \
175   a2.im = t2; \
176   a3.im = t7; \
177   a3.re = t8; \
178   }
179 
180 #define UNTRANSFORMHALF(a0,a1,a2,a3) { \
181   t1 = a2.re + a2.im; \
182   t2 = a2.im - a2.re; \
183   t1 *= sqrthalf; \
184   t3 = a3.re - a3.im; \
185   t2 *= sqrthalf; \
186   t4 = a3.im + a3.re; \
187   t3 *= sqrthalf; \
188   t7 = t2; \
189   t4 *= sqrthalf; \
190   t8 = t3; \
191   t8 -= t1; \
192   t7 -= t4; \
193   t1 += t3; \
194   t2 += t4; \
195   t4 = a1.im - t8; \
196   t8 += a1.im; \
197   t3 = a1.re - t7; \
198   t7 += a1.re; \
199   a3.im = t4; \
200   t5 = a0.re - t1; \
201   t1 += a0.re; \
202   a1.im = t8; \
203   t6 = a0.im - t2; \
204   t2 += a0.im; \
205   a2.re = t5; \
206   a0.re = t1; \
207   a2.im = t6; \
208   a1.re = t7; \
209   a3.re = t3; \
210   a0.im = t2; \
211   }
212 
213 #define UNTRANSFORMZERO(a0,a1,a2,a3) { \
214   t1 = a2.re + a3.re; \
215   t2 = a2.im + a3.im; \
216   t3 = a2.im - a3.im; \
217   t4 = a3.re - a2.re; \
218   t5 = a0.re - t1; \
219   t6 = a0.im - t2; \
220   t7 = a1.re - t3; \
221   t8 = a1.im - t4; \
222   t1 += a0.re; \
223   t2 += a0.im; \
224   t3 += a1.re; \
225   t4 += a1.im; \
226   a2.re = t5; \
227   a3.re = t7; \
228   a2.im = t6; \
229   a1.im = t4; \
230   a1.re = t3; \
231   a0.im = t2; \
232   a0.re = t1; \
233   a3.im = t8; \
234   }
235 
236 #define R(a0,a1,b0,b1,wre,wim) { \
237   t1 = a0 - a1; \
238   t2 = b0 - b1; \
239   t5 = t1 * wim; \
240   t6 = t2 * wim; \
241   t3 = VOL a0; \
242   t1 *= wre; \
243   t3 += a1; \
244   t2 *= wre; \
245   t1 -= t6; \
246   t4 = VOL b0; \
247   t2 += t5; \
248   t4 += b1; \
249   a0 = t3; \
250   b1 = t2; \
251   a1 = t4; \
252   b0 = t1; \
253   }
254 
255 #define RHALF(a0,a1,b0,b1) { \
256   t1 = a0 - a1; \
257   t2 = b0 - b1; \
258   t3 = a0 + a1; \
259   t5 = t1 - t2; \
260   t1 += t2; \
261   t4 = VOL b0; \
262   t5 *= sqrthalf; \
263   t4 += b1; \
264   t1 *= sqrthalf; \
265   a0 = t3; \
266   b1 = t1; \
267   a1 = t4; \
268   b0 = t5; \
269   }
270 
271 #define RZERO(a0,a1,b0,b1) { \
272   t1 = a0 - a1; \
273   t2 = b0 - b1; \
274   t3 = a0 + a1; \
275   t4 = b0 + b1; \
276   b0 = t1; \
277   a0 = t3; \
278   b1 = t2; \
279   a1 = t4; \
280   }
281 
282 #define V(a0,a1,b0,b1,wre,wim) { \
283   t5 = b0 * wre; \
284   t1 = b1 * wim; \
285   t6 = b1 * wre; \
286   t5 += t1; \
287   t3 = b0 * wim; \
288   t2 = a0 - t5; \
289   t6 -= t3; \
290   t5 += a0; \
291   t4 = a1 - t6; \
292   t6 += a1; \
293   a1 = t2; \
294   a0 = t5; \
295   b1 = t4; \
296   b0 = t6; \
297   }
298 
299 #define VHALF(a0,a1,b0,b1) { \
300   t5 = b0 + b1; \
301   t6 = b1 - b0; \
302   t5 *= sqrthalf; \
303   t2 = VOL a0; \
304   t6 *= sqrthalf; \
305   t2 -= t5; \
306   t5 += a0; \
307   t4 = a1 - t6; \
308   t6 += a1; \
309   a1 = t2; \
310   a0 = t5; \
311   b0 = t6; \
312   b1 = t4; \
313   }
314 
315 #define VZERO(a0,a1,b0,b1) { \
316   t1 = a0 + b0; \
317   t2 = a0 - b0; \
318   t3 = a1 + b1; \
319   t4 = a1 - b1; \
320   a0 = t1; \
321   b0 = t3; \
322   a1 = t2; \
323   b1 = t4; \
324   }
325