1 /* Compute the sum of the squares of a vector of signed shorts
2
3 * This is the Altivec SIMD version. It's a little hairy because Altivec
4 * does not do 64-bit operations directly, so we have to accumulate separate
5 * 32-bit sums and carries
6
7 * Copyright 2004 Phil Karn, KA9Q
8 * May be used under the terms of the GNU Lesser General Public License (LGPL)
9 */
10
11 #include <altivec.h>
12 #include "fec.h"
13
sumsq_av(signed short * in,int cnt)14 unsigned long long sumsq_av(signed short *in,int cnt){
15 long long sum;
16 vector signed short x;
17 vector unsigned int sums,carries,s1,s2;
18 int pad;
19 union { vector unsigned char cv; vector unsigned int iv; unsigned int w[4]; unsigned char c[16];} s;
20
21 carries = sums = (vector unsigned int){0};
22 if((pad = (int)in & 15)!=0){
23 /* Load unaligned leading word */
24 x = vec_perm(vec_ld(0,in),(vector signed short){0},vec_lvsl(0,in));
25 if(cnt < 8){ /* Shift right to chop stuff beyond end of short block */
26 s.c[15] = (8-cnt)<<4;
27 x = vec_sro(x,s.cv);
28 }
29 sums = (vector unsigned int)vec_msum(x,x,(vector signed int){0});
30 in += 8-pad/2;
31 cnt -= 8-pad/2;
32 }
33 /* Everything is now aligned, rip through most of the block */
34 while(cnt >= 8){
35 x = vec_ld(0,in);
36 /* A single vec_msum cannot overflow, but we have to sum it with
37 * the earlier terms separately to handle the carries
38 * The cast to unsigned is OK because squares are always positive
39 */
40 s1 = (vector unsigned int)vec_msum(x,x,(vector signed int){0});
41 carries = vec_add(carries,vec_addc(sums,s1));
42 sums = vec_add(sums,s1);
43 in += 8;
44 cnt -= 8;
45 }
46 /* Handle trailing fragment, if any */
47 if(cnt > 0){
48 x = vec_ld(0,in);
49 s.c[15] = (8-cnt)<<4;
50 x = vec_sro(x,s.cv);
51 s1 = (vector unsigned int)vec_msum(x,x,(vector signed int){0});
52 carries = vec_add(carries,vec_addc(sums,s1));
53 sums = vec_add(sums,s1);
54 }
55 /* Combine 4 sub-sums and carries */
56 s.c[15] = 64; /* Shift right two 32-bit words */
57 s1 = vec_sro(sums,s.cv);
58 s2 = vec_sro(carries,s.cv);
59 carries = vec_add(carries,vec_addc(sums,s1));
60 sums = vec_add(sums,s1);
61 carries = vec_add(carries,s2);
62
63 s.c[15] = 32; /* Shift right one 32-bit word */
64 s1 = vec_sro(sums,s.cv);
65 s2 = vec_sro(carries,s.cv);
66 carries = vec_add(carries,vec_addc(sums,s1));
67 sums = vec_add(sums,s1);
68 carries = vec_add(carries,s2);
69
70 /* Extract sum and carries from right-hand words and combine into result */
71 s.iv = sums;
72 sum = s.w[3];
73
74 s.iv = carries;
75 sum += (long long)s.w[3] << 32;
76
77 return sum;
78 }
79
80