Imported Upstream version 3.0
[debian/gnuradio] / gnuradio-core / src / lib / filter / ccomplex_dotprod_3dnow.S
1 #
2 # Copyright 2002 Free Software Foundation, Inc.
3
4 # This file is part of GNU Radio
5
6 # GNU Radio is free software; you can redistribute it and/or modify
7 # it under the terms of the GNU General Public License as published by
8 # the Free Software Foundation; either version 2, or (at your option)
9 # any later version.
10
11 # GNU Radio is distributed in the hope that it will be useful,
12 # but WITHOUT ANY WARRANTY; without even the implied warranty of
13 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
14 # GNU General Public License for more details.
15
16 # You should have received a copy of the GNU General Public License
17 # along with GNU Radio; see the file COPYING.  If not, write to
18 # the Free Software Foundation, Inc., 51 Franklin Street,
19 # Boston, MA 02110-1301, USA.
20
21
22
23 # input and taps are guarenteed to be 16 byte aligned.
24 # n_2_ccomplex_blocks is != 0
25 #       
26 #
27 #  ccomplex_dotprod_generic (const float *input,
28 #                         const float *taps, unsigned n_2_ccomplex_blocks, float *result)
29 #  {
30 #    float sum0 = 0;
31 #    float sum1 = 0;
32 #    float sum2 = 0;
33 #    float sum3 = 0;
34 #  
35 #    do {
36 #  
37 #      sum0 += input[0] * taps[0] - input[1] * taps[1];
38 #      sum1 += input[0] * taps[1] + input[1] * taps[0];
39 #      sum2 += input[2] * taps[2] - input[3] * taps[3];
40 #      sum3 += input[2] * taps[3] + input[3] * taps[2];
41 #  
42 #      input += 4;
43 #      taps += 4;
44 #  
45 #    } while (--n_2_ccomplex_blocks != 0);
46 #  
47 #  
48 #    result[0] = sum0 + sum2;
49 #    result[1] = sum1 + sum3;
50 #  }
51 #               
52
53 # TODO: prefetch and better scheduling
54
55 #include "assembly.h"
56
57         .file   "ccomplex_dotprod_3dnow.S"
58         .version        "01.01"
59 .text
60         .p2align 4
61 .globl GLOB_SYMB(ccomplex_dotprod_3dnow)
62         DEF_FUNC_HEAD(ccomplex_dotprod_3dnow)
63 GLOB_SYMB(ccomplex_dotprod_3dnow):
64         pushl   %ebp
65         movl    %esp, %ebp
66         movl    8(%ebp), %eax           # input
67         movl    12(%ebp), %edx          # taps
68         movl    16(%ebp), %ecx          # n_2_ccomplex_blocks
69
70         # zero accumulators
71         
72         pxor    %mm6, %mm6              # mm6 = 0 0 
73
74         movq    0(%eax), %mm0
75
76         pxor    %mm7, %mm7              # mm7 = 0 0 
77
78         movq    0(%edx), %mm2
79
80         movq    8(%eax), %mm1
81
82         shrl    $1, %ecx                # ecx = n_2_ccomplex_blocks / 2
83
84         movq    8(%edx), %mm3
85
86         jmp     .L1_test
87
88         #
89         # 4 taps / loop
90         # something like ?? cycles / loop
91         #
92         
93         .p2align 4
94 .loop1: 
95
96 # complex prod: C += A * B,  w/ temp Z, mmPN=$80000000
97 #
98 #       movq    (%eax), %mmA
99 #       movq    (%edx), %mmB
100 #
101 #       # 3DNow! replacement for: pswapd  %mmA, %mmZ
102 #       # TODO: optimize the punpckhdq
103 #       movq    %mmA, %mmZ
104 #       punpckhdq       %mmZ, %mmZ
105 #       punpckldq       %mmA, %mmZ
106 #
107 #       pfmul   %mmB, %mmA
108 #       pfmul   %mmZ, %mmB
109 #
110 #       # 3DNow! replacement for: pfpnacc %mmB, %mmA
111 #       pxor    %mmPN, %mmA
112 #       pfacc   %mmB, %mmA
113 #
114 #       pfadd   %mmA, %mmC
115
116
117 # A=mm0, B=mm2, Z=mm4
118 # A'=mm1, B'=mm3, Z'=mm5
119
120         movq    %mm0, %mm4
121         movq    %mm1, %mm5
122         punpckhdq       %mm4, %mm4
123         punpckhdq       %mm5, %mm5
124         punpckldq       %mm0, %mm4
125         pfmul   %mm2, %mm0
126         punpckldq       %mm1, %mm5
127         pfmul   %mm4, %mm2
128         pfadd   %mm0, %mm6
129         movq    16(%edx), %mm0
130         pfmul   %mm3, %mm1
131         pfadd   %mm2, %mm7
132         movq    16(%eax), %mm2
133         pfadd   %mm1, %mm6
134         pfmul   %mm5, %mm3
135         movq    24(%edx), %mm1
136
137         movq    %mm0, %mm4
138         movq    %mm1, %mm5
139
140         pfadd   %mm3, %mm7
141         movq    24(%eax), %mm3
142
143 # unroll
144
145         punpckhdq       %mm4, %mm4
146         punpckhdq       %mm5, %mm5
147         punpckldq       %mm0, %mm4
148         pfmul   %mm2, %mm0
149         punpckldq       %mm1, %mm5
150         pfmul   %mm4, %mm2
151         pfadd   %mm0, %mm6
152         movq    32(%edx), %mm0
153         pfmul   %mm3, %mm1
154         pfadd   %mm2, %mm7
155         movq    32(%eax), %mm2
156         pfadd   %mm1, %mm6
157         pfmul   %mm5, %mm3
158         movq    40(%edx), %mm1
159
160         addl    $32, %eax
161         addl    $32, %edx
162
163         pfadd   %mm3, %mm7
164         movq    8(%eax), %mm3
165
166 .L1_test:
167         decl    %ecx
168         jge     .loop1
169
170         # We've handled the bulk of multiplies up to here.
171         # Let's see if original n_2_ccomplex_blocks was odd.
172         # If so, we've got 2 more taps to do.
173         
174         movl    16(%ebp), %ecx          # n_2_ccomplex_blocks
175         andl    $1, %ecx
176         je      .Leven
177         
178         # The count was odd, do 2 more taps.
179         # Note that we've already got mm0/mm2 & mm1/mm3 preloaded
180         # from the main loop.
181
182         movq    %mm0, %mm4
183         movq    %mm1, %mm5
184         punpckhdq       %mm4, %mm4
185         punpckhdq       %mm5, %mm5
186         punpckldq       %mm0, %mm4
187         pfmul   %mm2, %mm0
188         punpckldq       %mm1, %mm5
189         pfmul   %mm4, %mm2
190         pfadd   %mm0, %mm6
191         pfmul   %mm3, %mm1
192         pfadd   %mm2, %mm7
193         pfmul   %mm5, %mm3
194         pfadd   %mm1, %mm6
195         pfadd   %mm3, %mm7
196
197 .Leven:
198         # mmNP: negative inversor
199
200         pcmpeqd %mm0, %mm0              # set all bits to 1
201         psllq   $63, %mm0               # keep only hsb
202
203         pxor    %mm0, %mm6
204         pfacc   %mm7, %mm6
205
206         movl    20(%ebp), %eax          # result
207         movq    %mm6, (%eax)
208
209         femms
210
211         popl    %ebp
212         ret
213
214 FUNC_TAIL(ccomplex_dotprod_3dnow)
215         .ident  "Hand coded x86 3DNow! assembly"
216