Line data Source code
1 : /*
2 : * libpoporon - bch.c
3 : *
4 : * This file is part of libpoporon.
5 : *
6 : * SPDX-License-Identifier: BSD-3-Clause
7 : */
8 :
9 : #include "internal/common.h"
10 :
11 : #define BCH_MAX_POLY 64
12 : #define BCH_MAX_T 16
13 :
14 : struct _poporon_bch_t {
15 : poporon_gf_t *gf;
16 : uint8_t correction_capability;
17 : uint16_t codeword_length;
18 : uint16_t data_length;
19 : uint16_t parity_bits;
20 : uint32_t gen_poly;
21 : uint16_t gen_poly_deg;
22 : };
23 :
24 321 : static inline int32_t bch_compute_syndromes(poporon_bch_t *bch, uint32_t codeword, uint16_t *syndromes)
25 : {
26 : poporon_gf_t *gf;
27 : uint16_t exp_val;
28 : int32_t i, j, has_nonzero, syndrome_count;
29 :
30 321 : gf = bch->gf;
31 321 : has_nonzero = 0;
32 321 : syndrome_count = 2 * bch->correction_capability;
33 :
34 2239 : for (i = 0; i < syndrome_count; i++) {
35 1918 : syndromes[i] = 0;
36 :
37 30944 : for (j = 0; j < bch->codeword_length; j++) {
38 29026 : if (codeword & (1U << j)) {
39 15376 : exp_val = (uint16_t)(((i + 1) * j) % gf->field_size);
40 15376 : syndromes[i] ^= gf->exp[exp_val];
41 : }
42 : }
43 :
44 1918 : if (syndromes[i] != 0) {
45 702 : has_nonzero = 1;
46 : }
47 : }
48 :
49 321 : return has_nonzero;
50 : }
51 :
52 1348 : static inline uint16_t bch_poly_eval(poporon_bch_t *bch, const uint16_t *poly, int32_t degree, uint16_t x)
53 : {
54 : poporon_gf_t *gf;
55 : uint16_t sum, log_x, exp_val;
56 : int32_t i;
57 :
58 1348 : gf = bch->gf;
59 :
60 1348 : if (x == 0) {
61 0 : return poly[0];
62 : }
63 :
64 1348 : sum = 0;
65 1348 : log_x = gf->log[x];
66 :
67 5327 : for (i = 0; i <= degree; i++) {
68 3979 : if (poly[i] != 0) {
69 3979 : exp_val = (gf->log[poly[i]] + (log_x * i) % gf->field_size) % gf->field_size;
70 3979 : sum ^= gf->exp[exp_val];
71 : }
72 : }
73 :
74 1348 : return sum;
75 : }
76 :
77 128 : static inline int32_t bch_berlekamp_massey(poporon_bch_t *bch, const uint16_t *syndromes, uint16_t *error_locator)
78 : {
79 : poporon_gf_t *gf;
80 : uint16_t current[BCH_MAX_POLY], prev[BCH_MAX_POLY], temp[BCH_MAX_POLY], prev_discrepancy, discrepancy, log_mult,
81 : multiplier, log_sum, log_product;
82 : int32_t error_count, shift, syndrome_count, iteration, i;
83 :
84 128 : gf = bch->gf;
85 128 : error_count = 0;
86 128 : shift = 1;
87 128 : prev_discrepancy = 1;
88 128 : syndrome_count = 2 * bch->correction_capability;
89 :
90 128 : pmemset(current, 0, sizeof(current));
91 128 : pmemset(prev, 0, sizeof(prev));
92 128 : pmemset(temp, 0, sizeof(temp));
93 128 : current[0] = 1;
94 128 : prev[0] = 1;
95 :
96 892 : for (iteration = 0; iteration < syndrome_count; iteration++) {
97 764 : discrepancy = syndromes[iteration];
98 :
99 1740 : for (i = 1; i <= error_count; i++) {
100 976 : if (current[i] != 0 && syndromes[iteration - i] != 0) {
101 914 : log_sum = (gf->log[current[i]] + gf->log[syndromes[iteration - i]]) % gf->field_size;
102 914 : discrepancy ^= gf->exp[log_sum];
103 : }
104 : }
105 :
106 764 : if (discrepancy == 0) {
107 518 : shift++;
108 : } else {
109 246 : log_mult = (gf->field_size - gf->log[prev_discrepancy] + gf->log[discrepancy]) % gf->field_size;
110 246 : multiplier = gf->exp[log_mult];
111 :
112 246 : if (2 * error_count <= iteration) {
113 246 : pmemcpy(temp, current, sizeof(temp));
114 :
115 15626 : for (i = 0; i < BCH_MAX_POLY - shift; i++) {
116 15380 : if (prev[i] != 0) {
117 252 : log_product = (gf->log[prev[i]] + gf->log[multiplier]) % gf->field_size;
118 252 : current[i + shift] ^= gf->exp[log_product];
119 : }
120 : }
121 :
122 246 : pmemcpy(prev, temp, sizeof(prev));
123 246 : error_count = iteration + 1 - error_count;
124 246 : prev_discrepancy = discrepancy;
125 246 : shift = 1;
126 : } else {
127 0 : for (i = 0; i < BCH_MAX_POLY - shift; i++) {
128 0 : if (prev[i] != 0) {
129 0 : log_product = (gf->log[prev[i]] + gf->log[multiplier]) % gf->field_size;
130 0 : current[i + shift] ^= gf->exp[log_product];
131 : }
132 : }
133 0 : shift++;
134 : }
135 : }
136 : }
137 :
138 128 : pmemcpy(error_locator, current, BCH_MAX_POLY * sizeof(uint16_t));
139 :
140 128 : return error_count;
141 : }
142 :
143 128 : static inline int32_t bch_chien_search(poporon_bch_t *bch, const uint16_t *error_locator, int32_t error_count,
144 : uint16_t *error_pos)
145 : {
146 : poporon_gf_t *gf;
147 : uint16_t alpha_inv;
148 128 : int32_t found = 0, i;
149 :
150 128 : gf = bch->gf;
151 :
152 1349 : for (i = 0; i < bch->codeword_length; i++) {
153 1348 : alpha_inv = gf->exp[(gf->field_size - i) % gf->field_size];
154 :
155 1348 : if (bch_poly_eval(bch, error_locator, error_count, alpha_inv) == 0) {
156 243 : error_pos[found++] = (uint16_t)i;
157 :
158 243 : if (found >= error_count) {
159 127 : break;
160 : }
161 : }
162 : }
163 :
164 128 : return found;
165 : }
166 :
167 38 : static inline uint32_t bch_get_minimal_polynomial(poporon_gf_t *gf, int32_t exp)
168 : {
169 : uint32_t binary_poly;
170 : uint16_t poly[BCH_MAX_POLY], root, log_prod;
171 : int32_t poly_deg, conjugate, i, j;
172 :
173 38 : pmemset(poly, 0, sizeof(poly));
174 38 : poly[0] = 1;
175 38 : poly_deg = 0;
176 :
177 38 : conjugate = exp;
178 :
179 : do {
180 135 : root = gf->exp[conjugate];
181 :
182 463 : for (j = poly_deg; j >= 0; j--) {
183 328 : if (j + 1 < BCH_MAX_POLY) {
184 328 : poly[j + 1] ^= poly[j];
185 : }
186 328 : if (poly[j] != 0 && root != 0) {
187 328 : log_prod = (gf->log[poly[j]] + gf->log[root]) % gf->field_size;
188 328 : poly[j] = gf->exp[log_prod];
189 : } else {
190 0 : poly[j] = 0;
191 : }
192 : }
193 135 : poly_deg++;
194 :
195 135 : conjugate = (conjugate * 2) % gf->field_size;
196 135 : } while (conjugate != exp);
197 :
198 38 : binary_poly = 0;
199 211 : for (i = 0; i <= poly_deg; i++) {
200 173 : if (poly[i] == 1) {
201 142 : binary_poly |= (1U << i);
202 : }
203 : }
204 :
205 38 : return binary_poly;
206 : }
207 :
208 38 : static inline uint32_t bch_poly_multiply_binary(uint32_t a, int32_t deg_a, uint32_t b)
209 : {
210 : int32_t i;
211 38 : uint32_t result = 0;
212 :
213 228 : for (i = 0; i <= deg_a; i++) {
214 190 : if (a & (1U << i)) {
215 114 : result ^= (b << i);
216 : }
217 : }
218 :
219 38 : return result;
220 : }
221 :
222 76 : static inline int32_t bch_poly_degree_binary(uint32_t poly)
223 : {
224 76 : int32_t deg = 0, i;
225 :
226 76 : if (poly == 0) {
227 0 : return -1;
228 : }
229 :
230 2010 : for (i = 31; i >= 0; i--) {
231 2010 : if (poly & (1U << i)) {
232 76 : deg = i;
233 76 : break;
234 : }
235 : }
236 :
237 76 : return deg;
238 : }
239 :
240 13 : static inline bool bch_build_generator(poporon_bch_t *bch)
241 : {
242 : poporon_gf_t *gf;
243 : uint32_t gen, min_poly;
244 : int32_t gen_deg, i, root_exp, conj, min_deg;
245 : bool *used;
246 :
247 13 : gf = bch->gf;
248 13 : used = (bool *)pcalloc(gf->field_size + 1, sizeof(bool));
249 :
250 13 : if (!used) {
251 0 : return false;
252 : }
253 :
254 13 : gen = 1;
255 13 : gen_deg = 0;
256 :
257 89 : for (i = 1; i <= 2 * bch->correction_capability; i++) {
258 76 : root_exp = i % gf->field_size;
259 :
260 76 : if (used[root_exp]) {
261 38 : continue;
262 : }
263 :
264 38 : conj = root_exp;
265 : do {
266 135 : used[conj] = true;
267 135 : conj = (conj * 2) % gf->field_size;
268 135 : } while (conj != root_exp);
269 :
270 38 : min_poly = bch_get_minimal_polynomial(gf, root_exp);
271 38 : min_deg = bch_poly_degree_binary(min_poly);
272 :
273 38 : gen = bch_poly_multiply_binary(gen, gen_deg, min_poly);
274 38 : gen_deg = bch_poly_degree_binary(gen);
275 : }
276 :
277 13 : pfree(used);
278 :
279 13 : bch->gen_poly = gen;
280 13 : bch->gen_poly_deg = (uint16_t)gen_deg;
281 13 : bch->parity_bits = (uint16_t)gen_deg;
282 13 : bch->data_length = bch->codeword_length - bch->parity_bits;
283 :
284 13 : return true;
285 : }
286 :
287 17 : extern poporon_bch_t *poporon_bch_create(uint8_t symbol_size, uint16_t generator_polynomial,
288 : uint8_t correction_capability)
289 : {
290 : poporon_bch_t *bch;
291 :
292 17 : if (symbol_size < 3 || symbol_size > 5) {
293 3 : return NULL;
294 : }
295 :
296 14 : if (correction_capability < 1 || correction_capability > BCH_MAX_T) {
297 1 : return NULL;
298 : }
299 :
300 13 : bch = (poporon_bch_t *)pcalloc(1, sizeof(poporon_bch_t));
301 13 : if (!bch) {
302 0 : return NULL;
303 : }
304 :
305 13 : bch->gf = poporon_gf_create(symbol_size, generator_polynomial);
306 13 : if (!bch->gf) {
307 0 : pfree(bch);
308 0 : return NULL;
309 : }
310 :
311 13 : bch->correction_capability = correction_capability;
312 13 : bch->codeword_length = (uint16_t)((1 << symbol_size) - 1);
313 :
314 13 : if (!bch_build_generator(bch)) {
315 0 : poporon_gf_destroy(bch->gf);
316 0 : pfree(bch);
317 0 : return NULL;
318 : }
319 :
320 13 : return bch;
321 : }
322 :
323 14 : extern void poporon_bch_destroy(poporon_bch_t *bch)
324 : {
325 14 : if (!bch) {
326 1 : return;
327 : }
328 :
329 13 : if (bch->gf) {
330 13 : poporon_gf_destroy(bch->gf);
331 : }
332 :
333 13 : pfree(bch);
334 : }
335 :
336 7 : extern uint16_t poporon_bch_get_codeword_length(const poporon_bch_t *bch)
337 : {
338 7 : return bch ? bch->codeword_length : 0;
339 : }
340 :
341 6 : extern uint16_t poporon_bch_get_data_length(const poporon_bch_t *bch)
342 : {
343 6 : return bch ? bch->data_length : 0;
344 : }
345 :
346 1 : extern uint8_t poporon_bch_get_correction_capability(const poporon_bch_t *bch)
347 : {
348 1 : return bch ? bch->correction_capability : 0;
349 : }
350 :
351 74 : extern bool poporon_bch_encode(poporon_bch_t *bch, uint32_t data, uint32_t *codeword)
352 : {
353 : uint32_t shifted, remainder, gen;
354 : int32_t gen_deg, i;
355 :
356 74 : if (!bch || !codeword) {
357 2 : return false;
358 : }
359 :
360 72 : if (data >= (1U << bch->data_length)) {
361 1 : return false;
362 : }
363 :
364 71 : shifted = data << bch->parity_bits;
365 :
366 71 : remainder = shifted;
367 71 : gen = bch->gen_poly;
368 71 : gen_deg = bch->gen_poly_deg;
369 :
370 442 : for (i = bch->codeword_length - 1; i >= gen_deg; i--) {
371 371 : if (remainder & (1U << i)) {
372 181 : remainder ^= (gen << (i - gen_deg));
373 : }
374 : }
375 :
376 71 : *codeword = shifted ^ remainder;
377 :
378 71 : return true;
379 : }
380 :
381 196 : extern bool poporon_bch_decode(poporon_bch_t *bch, uint32_t received, uint32_t *corrected, int32_t *num_errors)
382 : {
383 : uint32_t corrected_word;
384 : uint16_t syndromes[BCH_MAX_POLY], error_locator[BCH_MAX_POLY], error_positions[BCH_MAX_T];
385 : int32_t error_count, found, i;
386 :
387 196 : if (!bch || !corrected) {
388 2 : return false;
389 : }
390 :
391 194 : pmemset(syndromes, 0, sizeof(syndromes));
392 194 : pmemset(error_locator, 0, sizeof(error_locator));
393 194 : pmemset(error_positions, 0, sizeof(error_positions));
394 :
395 194 : received &= ((1U << bch->codeword_length) - 1);
396 194 : *corrected = received;
397 :
398 194 : if (num_errors) {
399 194 : *num_errors = 0;
400 : }
401 :
402 194 : if (!bch_compute_syndromes(bch, received, syndromes)) {
403 66 : return true;
404 : }
405 :
406 128 : error_count = bch_berlekamp_massey(bch, syndromes, error_locator);
407 :
408 128 : if (error_count > bch->correction_capability) {
409 0 : return false;
410 : }
411 :
412 128 : found = bch_chien_search(bch, error_locator, error_count, error_positions);
413 :
414 128 : if (found != error_count) {
415 1 : return false;
416 : }
417 :
418 127 : corrected_word = received;
419 370 : for (i = 0; i < found; i++) {
420 243 : corrected_word ^= (1U << error_positions[i]);
421 : }
422 :
423 127 : if (bch_compute_syndromes(bch, corrected_word, syndromes)) {
424 0 : return false;
425 : }
426 :
427 127 : *corrected = corrected_word;
428 :
429 127 : if (num_errors) {
430 127 : *num_errors = found;
431 : }
432 :
433 127 : return true;
434 : }
435 :
436 158 : extern uint32_t poporon_bch_extract_data(const poporon_bch_t *bch, uint32_t codeword)
437 : {
438 158 : if (!bch) {
439 0 : return 0;
440 : }
441 :
442 158 : return (codeword >> bch->parity_bits) & ((1U << bch->data_length) - 1);
443 : }
|