1 /*
2  * Copyright (c) 2007 - 2015 Joseph Gaeddert
3  *
4  * Permission is hereby granted, free of charge, to any person obtaining a copy
5  * of this software and associated documentation files (the "Software"), to deal
6  * in the Software without restriction, including without limitation the rights
7  * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
8  * copies of the Software, and to permit persons to whom the Software is
9  * furnished to do so, subject to the following conditions:
10  *
11  * The above copyright notice and this permission notice shall be included in
12  * all copies or substantial portions of the Software.
13  *
14  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
15  * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
16  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
17  * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
18  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
19  * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
20  * THE SOFTWARE.
21  */
22 
23 //
24 // Golay(24,12) half-rate forward error-correction code
25 //
26 // References:
27 //  [Lin:2004] Lin, Shu and Costello, Daniel L. Jr., "Error Control
28 //      Coding," Prentice Hall, New Jersey, 2nd edition, 2004.
29 //
30 
31 #include <stdio.h>
32 #include <stdlib.h>
33 #include <assert.h>
34 
35 #include "liquid.internal.h"
36 
37 #define DEBUG_FEC_GOLAY2412 0
38 
39 // P matrix [12 x 12]
40 unsigned int golay2412_P[12] = {
41     0x08ed, 0x01db, 0x03b5, 0x0769,
42     0x0ed1, 0x0da3, 0x0b47, 0x068f,
43     0x0d1d, 0x0a3b, 0x0477, 0x0ffe};
44 
45 #if 0
46 // generator matrix [12 x 24]
47 unsigned int golay2412_G[12] = {
48     0x008ed800, 0x001db400, 0x003b5200, 0x00769100,
49     0x00ed1080, 0x00da3040, 0x00b47020, 0x0068f010,
50     0x00d1d008, 0x00a3b004, 0x00477002, 0x00ffe001};
51 #endif
52 
53 // generator matrix transposed [24 x 12]
54 unsigned int golay2412_Gt[24] = {
55     0x08ed, 0x01db, 0x03b5, 0x0769, 0x0ed1, 0x0da3, 0x0b47, 0x068f,
56     0x0d1d, 0x0a3b, 0x0477, 0x0ffe, 0x0800, 0x0400, 0x0200, 0x0100,
57     0x0080, 0x0040, 0x0020, 0x0010, 0x0008, 0x0004, 0x0002, 0x0001};
58 
59 // parity check matrix [12 x 24]
60 unsigned int golay2412_H[12] = {
61     0x008008ed, 0x004001db, 0x002003b5, 0x00100769,
62     0x00080ed1, 0x00040da3, 0x00020b47, 0x0001068f,
63     0x00008d1d, 0x00004a3b, 0x00002477, 0x00001ffe};
64 
65 // multiply input vector with parity check matrix, H
golay2412_matrix_mul(unsigned int _v,unsigned int * _A,unsigned int _n)66 unsigned int golay2412_matrix_mul(unsigned int   _v,
67                                   unsigned int * _A,
68                                   unsigned int   _n)
69 {
70     unsigned int x = 0;
71     unsigned int i;
72     for (i=0; i<_n; i++) {
73         x <<= 1;
74 #if 0
75         // compute dot product mod 2
76         x |= liquid_count_ones_mod2( _A[i] & _v );
77 #else
78         // same as above, but exploit the fact that vectors are at
79         // most 24 bits long; liquid_count_ones_mod2() assumes full
80         // 32- or 64-bit integer
81         unsigned int c = 0;
82         unsigned int p = _A[i] & _v;
83         c += liquid_c_ones[ p & 0xff ]; p >>= 8;
84         c += liquid_c_ones[ p & 0xff ]; p >>= 8;
85         c += liquid_c_ones[ p & 0xff ];
86 
87         x |= c & 0x0001;    // mod 2
88 #endif
89     }
90     return x;
91 }
92 
fec_golay2412_encode_symbol(unsigned int _sym_dec)93 unsigned int fec_golay2412_encode_symbol(unsigned int _sym_dec)
94 {
95     // validate input
96     if (_sym_dec >= (1<<12)) {
97         fprintf(stderr,"error, fec_golay2412_encode_symbol(), input symbol too large\n");
98         exit(1);
99     }
100 
101     // compute encoded/transmitted message: v = m*G
102     return golay2412_matrix_mul(_sym_dec, golay2412_Gt, 24);
103 }
104 
105 // search for p[i] such that w(v+p[i]) <= 2, return -1 on fail
golay2412_parity_search(unsigned int _v)106 int golay2412_parity_search(unsigned int _v)
107 {
108     //assert( _v < (1<<12) );
109 
110     unsigned int i;
111     for (i=0; i<12; i++) {
112 #if 0
113         unsigned int wj = liquid_count_ones(_v ^ golay2412_P[i]);
114 #else
115         // same as above but faster, exploiting fact that P has
116         // only 12 bits of resolution
117         unsigned int wj = 0;
118         unsigned int p  = _v ^ golay2412_P[i];
119         wj += liquid_c_ones[ (p     ) & 0xff ];
120         wj += liquid_c_ones[ (p >> 8) & 0xff ];
121 #endif
122         if (wj <= 2)
123             return i;
124     }
125 
126     // could not find p[i] to satisfy criteria
127     return -1;
128 }
129 
fec_golay2412_decode_symbol(unsigned int _sym_enc)130 unsigned int fec_golay2412_decode_symbol(unsigned int _sym_enc)
131 {
132     // validate input
133     if (_sym_enc >= (1<<24)) {
134         fprintf(stderr,"error, fec_golay2412_decode_symbol(), input symbol too large\n");
135         exit(1);
136     }
137 
138     // state variables
139     unsigned int s=0;       // syndrome vector
140     unsigned int e_hat=0;   // estimated error vector
141     unsigned int v_hat=0;   // estimated transmitted message
142     unsigned int m_hat=0;   // estimated original message
143 
144     // compute syndrome vector, s = r*H^T = ( H*r^T )^T
145     s = golay2412_matrix_mul(_sym_enc, golay2412_H, 12);
146 #if DEBUG_FEC_GOLAY2412
147     printf("s (syndrome vector): "); liquid_print_bitstring(s,12); printf("\n");
148 #endif
149 
150     // compute weight of s (12 bits)
151     unsigned int ws = liquid_count_ones_uint16(s);
152 #if DEBUG_FEC_GOLAY2412
153     printf("w(s) = %u\n", ws);
154 #endif
155 
156     // step 2:
157     e_hat = 0;
158     if (ws <= 3) {
159 #if DEBUG_FEC_GOLAY2412
160         printf("    w(s) <= 3: estimating error vector as [s, 0(12)]\n");
161 #endif
162         // set e_hat = [s 0(12)]
163         e_hat = (s << 12) & 0xfff000;
164     } else {
165         // step 3: search for p[i] s.t. w(s+p[i]) <= 2
166 #if DEBUG_FEC_GOLAY2412
167         printf("    searching for w(s + p_i) <= 2...\n");
168 #endif
169         int s_index = golay2412_parity_search(s);
170 
171         if (s_index >= 0) {
172             // vector found!
173 #if DEBUG_FEC_GOLAY2412
174             printf("    w(s + p[%2u]) <= 2: estimating error vector as [s+p[%2u],u[%2u]]\n", s_index, s_index, s_index);
175 #endif
176             // NOTE : uj = 1 << (12-j-1)
177             e_hat = ((s ^ golay2412_P[s_index]) << 12) | (1 << (11-s_index));
178         } else {
179             // step 4: compute s*P
180             unsigned int sP = golay2412_matrix_mul(s, golay2412_P, 12);
181 #if DEBUG_FEC_GOLAY2412
182             printf("s*P: "); liquid_print_bitstring(sP,12); printf("\n");
183 #endif
184 
185             // compute weight of sP (12 bits)
186             unsigned int wsP = liquid_count_ones_uint16(sP);
187 #if DEBUG_FEC_GOLAY2412
188             printf("w(s*P) = %u\n", wsP);
189 #endif
190 
191             if (wsP == 2 || wsP == 3) {
192                 // step 5: set e = [0, s*P]
193 #if DEBUG_FEC_GOLAY2412
194                 printf("    w(s*P) in [2,3]: estimating error vector as [0(12), s*P]\n");
195 #endif
196                 e_hat = sP;
197             } else {
198                 // step 6: search for p[i] s.t. w(s*P + p[i]) == 2...
199 
200 #if DEBUG_FEC_GOLAY2412
201                 printf("    searching for w(s*P + p_i) == 2...\n");
202 #endif
203                 int sP_index = golay2412_parity_search(sP);
204 
205                 if (sP_index >= 0) {
206                     // vector found!
207 #if DEBUG_FEC_GOLAY2412
208                     printf("    w(s*P + p[%2u]) == 2: estimating error vector as [u[%2u],s*P+p[%2u]]\n", sP_index, sP_index, sP_index);
209 #endif
210                     // NOTE : uj = 1 << (12-j-1)
211                     //      [      uj << 1 2    ] [    sP + p[j]    ]
212                     e_hat = (1 << (23-sP_index)) | (sP ^ golay2412_P[sP_index]);
213                 } else {
214                     // step 7: decoding error
215 #if DEBUG_FEC_GOLAY2412
216                     printf("  **** decoding error\n");
217 #endif
218                 }
219             }
220         }
221     }
222 
223     // step 8: compute estimated transmitted message: v_hat = r + e_hat
224     v_hat = _sym_enc ^ e_hat;
225 #if DEBUG_FEC_GOLAY2412
226     printf("r (recevied vector):            "); liquid_print_bitstring(_sym_enc,24); printf("\n");
227     printf("e-hat (estimated error vector): "); liquid_print_bitstring(e_hat,24);    printf("\n");
228     printf("v-hat (estimated tx vector):    "); liquid_print_bitstring(v_hat,24);    printf("\n");
229 #endif
230 
231     // compute estimated original message: (last 12 bits of encoded message)
232     m_hat = v_hat & 0x0fff;
233 
234     return m_hat;
235 }
236 
237 // create Golay(24,12) codec object
fec_golay2412_create(void * _opts)238 fec fec_golay2412_create(void * _opts)
239 {
240     fec q = (fec) malloc(sizeof(struct fec_s));
241 
242     // set scheme
243     q->scheme = LIQUID_FEC_GOLAY2412;
244     q->rate = fec_get_rate(q->scheme);
245 
246     // set internal function pointers
247     q->encode_func      = &fec_golay2412_encode;
248     q->decode_func      = &fec_golay2412_decode;
249     q->decode_soft_func = NULL;
250 
251     return q;
252 }
253 
254 // destroy Golay(24,12) object
fec_golay2412_destroy(fec _q)255 void fec_golay2412_destroy(fec _q)
256 {
257     free(_q);
258 }
259 
260 // encode block of data using Golay(24,12) encoder
261 //
262 //  _q              :   encoder/decoder object
263 //  _dec_msg_len    :   decoded message length (number of bytes)
264 //  _msg_dec        :   decoded message [size: 1 x _dec_msg_len]
265 //  _msg_enc        :   encoded message [size: 1 x 2*_dec_msg_len]
fec_golay2412_encode(fec _q,unsigned int _dec_msg_len,unsigned char * _msg_dec,unsigned char * _msg_enc)266 void fec_golay2412_encode(fec _q,
267                           unsigned int _dec_msg_len,
268                           unsigned char *_msg_dec,
269                           unsigned char *_msg_enc)
270 {
271     unsigned int i=0;           // decoded byte counter
272     unsigned int j=0;           // encoded byte counter
273     unsigned int s0, s1, s2;    // three 8-bit symbols
274     unsigned int m0, m1;        // two 12-bit symbols (uncoded)
275     unsigned int v0, v1;        // two 24-bit symbols (encoded)
276 
277     // determine remainder of input length / 3
278     unsigned int r = _dec_msg_len % 3;
279 
280     for (i=0; i<_dec_msg_len-r; i+=3) {
281         // strip three input bytes (two uncoded symbols)
282         s0 = _msg_dec[i+0];
283         s1 = _msg_dec[i+1];
284         s2 = _msg_dec[i+2];
285 
286         // pack into two 12-bit symbols
287         m0 = ((s0 << 4) & 0x0ff0) | ((s1 >> 4) & 0x000f);
288         m1 = ((s1 << 8) & 0x0f00) | ((s2     ) & 0x00ff);
289 
290         // encode each 12-bit symbol into a 24-bit symbol
291         v0 = fec_golay2412_encode_symbol(m0);
292         v1 = fec_golay2412_encode_symbol(m1);
293 
294         // unpack two 24-bit symbols into six 8-bit bytes
295         // retaining order of bits in output
296         _msg_enc[j+0] = (v0 >> 16) & 0xff;
297         _msg_enc[j+1] = (v0 >>  8) & 0xff;
298         _msg_enc[j+2] = (v0      ) & 0xff;
299         _msg_enc[j+3] = (v1 >> 16) & 0xff;
300         _msg_enc[j+4] = (v1 >>  8) & 0xff;
301         _msg_enc[j+5] = (v1      ) & 0xff;
302 
303         j += 6;
304     }
305 
306     // if input length isn't divisible by 3, encode last 1 or two bytes
307     for (i=_dec_msg_len-r; i<_dec_msg_len; i++) {
308         // strip last input symbol
309         s0 = _msg_dec[i];
310 
311         // extend as 12-bit symbol
312         m0 = s0;
313 
314         // encode into 24-bit symbol
315         v0 = fec_golay2412_encode_symbol(m0);
316 
317         // unpack one 24-bit symbol into three 8-bit bytes, and
318         // append to output array
319         _msg_enc[j+0] = ( v0 >> 16 ) & 0xff;
320         _msg_enc[j+1] = ( v0 >>  8 ) & 0xff;
321         _msg_enc[j+2] = ( v0       ) & 0xff;
322 
323         j += 3;
324     }
325 
326     assert( j == fec_get_enc_msg_length(LIQUID_FEC_GOLAY2412,_dec_msg_len) );
327     assert( i == _dec_msg_len);
328 }
329 
330 // decode block of data using Golay(24,12) decoder
331 //
332 //  _q              :   encoder/decoder object
333 //  _dec_msg_len    :   decoded message length (number of bytes)
334 //  _msg_enc        :   encoded message [size: 1 x 2*_dec_msg_len]
335 //  _msg_dec        :   decoded message [size: 1 x _dec_msg_len]
336 //
337 //unsigned int
fec_golay2412_decode(fec _q,unsigned int _dec_msg_len,unsigned char * _msg_enc,unsigned char * _msg_dec)338 void fec_golay2412_decode(fec _q,
339                           unsigned int _dec_msg_len,
340                           unsigned char *_msg_enc,
341                           unsigned char *_msg_dec)
342 {
343     unsigned int i=0;                       // decoded byte counter
344     unsigned int j=0;                       // encoded byte counter
345     unsigned int r0, r1, r2, r3, r4, r5;    // six 8-bit bytes
346     unsigned int v0, v1;                    // two 24-bit encoded symbols
347     unsigned int m0_hat, m1_hat;            // two 12-bit decoded symbols
348 
349     // determine remainder of input length / 3
350     unsigned int r = _dec_msg_len % 3;
351 
352     for (i=0; i<_dec_msg_len-r; i+=3) {
353         // strip six input bytes (two encoded symbols)
354         r0 = _msg_enc[j+0];
355         r1 = _msg_enc[j+1];
356         r2 = _msg_enc[j+2];
357         r3 = _msg_enc[j+3];
358         r4 = _msg_enc[j+4];
359         r5 = _msg_enc[j+5];
360 
361         // pack six 8-bit symbols into two 24-bit symbols
362         v0 = ((r0 << 16) & 0xff0000) | ((r1 <<  8) & 0x00ff00) | ((r2     ) & 0x0000ff);
363         v1 = ((r3 << 16) & 0xff0000) | ((r4 <<  8) & 0x00ff00) | ((r5 << 0) & 0x0000ff);
364 
365         // decode each symbol into a 12-bit symbol
366         m0_hat = fec_golay2412_decode_symbol(v0);
367         m1_hat = fec_golay2412_decode_symbol(v1);
368 
369         // unpack two 12-bit symbols into three 8-bit bytes
370         _msg_dec[i+0] = ((m0_hat >> 4) & 0xff);
371         _msg_dec[i+1] = ((m0_hat << 4) & 0xf0) | ((m1_hat >> 8) & 0x0f);
372         _msg_dec[i+2] = ((m1_hat     ) & 0xff);
373 
374         j += 6;
375     }
376 
377     // if input length isn't divisible by 3, decode last 1 or two bytes
378     for (i=_dec_msg_len-r; i<_dec_msg_len; i++) {
379         // strip last input symbol (three bytes)
380         r0 = _msg_enc[j+0];
381         r1 = _msg_enc[j+1];
382         r2 = _msg_enc[j+2];
383 
384         // pack three 8-bit symbols into one 24-bit symbol
385         v0 = ((r0 << 16) & 0xff0000) | ((r1 <<  8) & 0x00ff00) | ((r2     ) & 0x0000ff);
386 
387         // decode into a 12-bit symbol
388         m0_hat = fec_golay2412_decode_symbol(v0);
389 
390         // retain last 8 bits of 12-bit symbol
391         _msg_dec[i] = m0_hat & 0xff;
392 
393         j += 3;
394     }
395 
396     assert( j== fec_get_enc_msg_length(LIQUID_FEC_GOLAY2412,_dec_msg_len) );
397     assert( i == _dec_msg_len);
398 
399     //return num_errors;
400 }
401 
402