1 #include <assert.h>
2 #include <stdio.h>
3 #include <stdint.h>
4 #include <stdlib.h>
5 #include <string.h>
6
7 #define BUFSZ (64 * 64 * sizeof(uint16_t))
8 #define NUM_COEFF_BUCKETS (4)
9 #define NUM_OTHER_BUCKETS (0)
10 #define NUM_TOTAL_BUCKETS ((NUM_COEFF_BUCKETS) + (NUM_OTHER_BUCKETS))
11 #define MAX_COEFF_BUCKET ((NUM_COEFF_BUCKETS) - 1)
12
13 #define clz(x) __builtin_clz(x)
14 #define ilog2(x) (sizeof(x) * 8 - clz(x) - 1)
15
print_coeffs(const int16_t * buf,uint32_t size,uint32_t ccc)16 void print_coeffs(const int16_t *buf, uint32_t size, uint32_t ccc)
17 {
18 uint32_t i;
19 printf("Buf size %u, ccc %u\n", size, ccc);
20 for (i = 0; i < size; i++)
21 printf("%i ", buf[i]);
22 printf("\n");
23 }
24
count_coeffs(const int16_t * buf,uint32_t size,uint64_t * buckets,uint64_t * num_signs,uint16_t * excess)25 void count_coeffs(const int16_t *buf, uint32_t size, uint64_t *buckets, uint64_t *num_signs, uint16_t *excess)
26 {
27 *excess = 0;
28 uint32_t i;
29
30 for (i = 0; i < size; i++) {
31 int16_t curr = buf[i];
32 int16_t is_signed = curr >> 15;
33 *num_signs += (is_signed & 1);
34
35 uint16_t abs = (curr ^ is_signed) - is_signed;
36 if (abs > MAX_COEFF_BUCKET) {
37 *excess += abs - MAX_COEFF_BUCKET;
38 abs = MAX_COEFF_BUCKET;
39 }
40
41 buckets[abs]++;
42 }
43 }
44
print_buckets(const uint64_t * buckets,uint64_t num_signs)45 void print_buckets(const uint64_t *buckets, uint64_t num_signs)
46 {
47 uint32_t i;
48 for (i = 0; i < NUM_COEFF_BUCKETS; i++)
49 printf("%3u: %lu\n", i, buckets[i]);
50 printf("Signs: %lu\n", num_signs);
51 }
52
update_matrix(const uint64_t * buckets,uint64_t * mat)53 void update_matrix(const uint64_t *buckets, uint64_t *mat)
54 {
55 for (int y = 0; y < NUM_TOTAL_BUCKETS; y++) {
56 for (int x = 0; x < NUM_TOTAL_BUCKETS; x++) {
57 int curr_pos = y * NUM_TOTAL_BUCKETS + x;
58 mat[curr_pos] += buckets[x] * buckets[y];
59 }
60 }
61 }
62
is_power_of_two(uint32_t u)63 static inline int is_power_of_two(uint32_t u)
64 {
65 return (u & (u - 1)) == 0;
66 }
67
process_rdcosts(FILE * in,FILE * out)68 int process_rdcosts(FILE *in, FILE *out)
69 {
70 void *buf = malloc(BUFSZ);
71 uint32_t *u32buf = (uint32_t *)buf;
72 int16_t *i16buf = (int16_t *)buf;
73 int rv = 0;
74
75 float weights[NUM_TOTAL_BUCKETS] = {0.0f};
76
77 uint64_t mat[NUM_TOTAL_BUCKETS * NUM_TOTAL_BUCKETS] = {0};
78
79 while (!feof(in)) {
80 uint32_t size, ccc, size_sqrt;
81 uint64_t cg_buckets[NUM_TOTAL_BUCKETS] = {0};
82 uint64_t cg_num_signs = 0;
83 uint16_t excess = 0;
84 size_t n_read;
85
86 n_read = fread(buf, sizeof(uint32_t), 2, in);
87 size = u32buf[0];
88 ccc = u32buf[1];
89
90 // Can't rely on feof() alone when reading from a pipe that might only get
91 // closed long after the last data has been poured in
92 if (n_read == 0) {
93 break;
94 }
95 if (feof(in) || n_read < 2) {
96 fprintf(stderr, "Unexpected EOF when reading header, managed still to read %u u32's\n", n_read);
97 rv = 1;
98 goto out;
99 }
100 if (!is_power_of_two(size)) {
101 fprintf(stderr, "Errorneous block size %u\n", size);
102 rv = 1;
103 goto out;
104 }
105
106 size_sqrt = 1 << (ilog2(size) >> 1);
107 n_read = fread(buf, sizeof(int16_t), size, in);
108 if (n_read != size) {
109 fprintf(stderr, "Unexpected EOF when reading block, managed still to read %u i16's\n", n_read);
110 rv = 1;
111 goto out;
112 }
113
114 count_coeffs(i16buf, size, cg_buckets, &cg_num_signs, &excess);
115 update_matrix(cg_buckets, mat);
116 }
117 for (int y = 0; y < NUM_TOTAL_BUCKETS; y++) {
118 for (int x = 0; x < NUM_TOTAL_BUCKETS; x++) {
119 int curr_pos = y * NUM_TOTAL_BUCKETS + x;
120 printf("%lu ", mat[curr_pos]);
121 }
122 printf("\n");
123 }
124 fflush(stdout);
125
126 out:
127 free(buf);
128 return rv;
129 }
130
main(int ar,char ** av)131 int main(int ar, char **av)
132 {
133 return process_rdcosts(stdin, stdout);
134 }
135