CECCO: Linear Codes and Decoding
CECCO supports linear codes over finite fields and (partially) over the rationals. For details about finite field construction consult the demo on prime fields, the demo on finite extension fields, or the documentation
The main example here is based on the binary field F4 (as extension of F2 with irreducible polynomial 1 + x + x^2, runtime-generated LUTs). It can easily be modified to any other finite field.
At the end of the demo we show an example of pseudo-soft GMD decoding, which (in the form shown here) requires binary extension fields. CECCO also supports soft-input decoding for linear codes over arbitrary finite fields, cf. the demo on soft decoding.
using F = Ext<Fp<2>, MOD{1, 1, 1}>;
After selecting some code length n and code dimension k we calculate a random k x n matrix with full rank:
const size_t n = 13;
const size_t k = 6;
auto G = ZeroMatrix<F>(k, n);
do {
G.randomize();
} while (G.rank() < k);
Based on the matrix we can set up a (random) linear code C'[F_2; n, k]:
auto Cp = LinearCode(n, k, G);
CECCO can determine a variety of properties of a linear code. For some properties, this requires extensive
calculations. For that reason, the user can choose the level of detail for outputting a code to the console.
With the IO manipulator CECCO::showbasic only field, length, and dimension are shown. The
manipulator
CECCO::showmost adds generator matrix G, parity check matrix H as well as some other properties.
The
manipulator CECCO::showall shows all properties, even the ones with expensive calculations such
as the
weight enumerator A(x):
// std::cout << showbasic << Cp << std::endl;
// std::cout << showmost << Cp << std::endl;
std::cout << showall << Cp << std::endl;
A possible output (with the CECCO::showall manipulator) is:
[F_2; 17, 9], dmin = 4
Linear code with properties: {
G =
⌈1 1 0 1 0 0 0 1 1 0 1 1 0 1 1 1 0⌉
|0 0 1 1 0 0 1 0 1 1 0 0 1 0 1 0 1|
|0 1 1 1 0 0 1 1 1 1 0 1 1 1 0 0 0|
|0 0 0 0 0 1 1 0 1 1 1 0 0 0 1 1 1|
|1 0 0 1 0 0 0 0 0 0 1 0 1 0 0 0 0|
|1 1 1 1 1 1 0 0 1 0 1 1 1 1 1 1 1|
|0 0 0 0 1 0 0 1 1 0 0 0 1 0 1 1 0|
|1 1 1 1 1 0 0 0 1 1 1 1 0 0 1 1 1|
⌊1 1 1 1 1 1 1 1 0 0 0 0 0 0 1 1 0⌋
H =
⌈1 0 0 0 0 0 0 0 1 1 1 1 0 1 0 1 0⌉
|0 1 0 0 0 0 0 0 0 1 1 1 1 0 1 0 1|
|0 0 1 0 0 0 0 0 0 0 1 1 1 1 0 1 0|
|0 0 0 1 0 0 0 0 0 0 0 1 1 1 1 0 1|
|0 0 0 0 1 0 0 0 1 1 1 1 1 0 1 0 0|
|0 0 0 0 0 1 0 0 0 1 1 1 1 1 0 1 0|
|0 0 0 0 0 0 1 0 0 0 1 1 1 1 1 0 1|
⌊0 0 0 0 0 0 0 1 1 1 1 0 1 0 1 0 0⌋
A(x) = 1 + 24x^4 + 87x^6 + 185x^8 + 166x^10 + 46x^12 + 3x^14 tmax = 1 polynomial(gamma = 1 + x^4 + x^7 + x^8)
}
We see the code parameters, the minimum distance, G and H, the weight enumerator A(x). The number of correctable errors with a BD decoder is shown as tmax=1 and CECCO also found that the random code happens to be polynomial with generator polynomial 1 + x^4 + x^7 + x^8 (cf. the demo on polynomial and cyclic codes).
We typically prefer systematic encoding with the message in the leftmost k components of each codeword (this requires a generator matrix in standard form, with a k x k identity matrix in the leftmost k columns). We can realize this by moving from code C' to code C:
auto C = Cp.get_equivalent_code_in_standard_form();
assert(C.is_equivalent(Cp));
std::cout << showall << C << std::endl;
[F_2; 17, 9], dmin = 4
Linear code with properties: {
G =
⌈1 0 0 0 0 0 0 0 0 0 0 0 1 0 1 1 0⌉
|0 1 0 0 0 0 0 0 0 0 0 0 0 1 0 1 1|
|0 0 1 0 0 0 0 0 0 1 0 0 0 1 1 0 0|
|0 0 0 1 0 0 0 0 0 0 1 0 0 0 1 1 0|
|0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 1 1|
|0 0 0 0 0 1 0 0 0 1 0 0 1 1 0 0 0|
|0 0 0 0 0 0 1 0 0 0 1 0 0 1 1 0 0|
|0 0 0 0 0 0 0 1 0 0 0 1 0 0 1 1 0|
⌊0 0 0 0 0 0 0 0 1 0 0 0 1 0 0 1 1⌋
H =
⌈1 0 0 0 0 0 0 0 1 1 1 1 0 1 0 1 0⌉
|0 1 0 0 0 0 0 0 0 1 1 1 1 0 1 0 1|
|0 0 1 0 0 0 0 0 0 0 1 1 1 1 0 1 0|
|0 0 0 1 0 0 0 0 0 0 0 1 1 1 1 0 1|
|0 0 0 0 1 0 0 0 1 1 1 1 1 0 1 0 0|
|0 0 0 0 0 1 0 0 0 1 1 1 1 1 0 1 0|
|0 0 0 0 0 0 1 0 0 0 1 1 1 1 1 0 1|
⌊0 0 0 0 0 0 0 1 1 1 1 0 1 0 1 0 0⌋
A(x) = 1 + 24x^4 + 87x^6 + 185x^8 + 166x^10 + 46x^12 + 3x^14 tmax = 1 polynomial(gamma = 1 + x^4 + x^7 + x^8)
}
CECCO can deliver the dual code and applies MacWilliams' identity in order to calculate the weight enumerator. We see that the dual code is not polynomial (other than cyclicity, polynomiality is not preserved by the dual code).
auto Cd = C.get_dual();
std::cout << showall << Cd << std::endl;
[F_2; 17, 8], dmin = 3
Linear code with properties: {
G =
⌈1 0 0 0 0 0 0 0 1 1 1 1 0 1 0 1 0⌉
|0 1 0 0 0 0 0 0 0 1 1 1 1 0 1 0 1|
|0 0 1 0 0 0 0 0 0 0 1 1 1 1 0 1 0|
|0 0 0 1 0 0 0 0 0 0 0 1 1 1 1 0 1|
|0 0 0 0 1 0 0 0 1 1 1 1 1 0 1 0 0|
|0 0 0 0 0 1 0 0 0 1 1 1 1 1 0 1 0|
|0 0 0 0 0 0 1 0 0 0 1 1 1 1 1 0 1|
⌊0 0 0 0 0 0 0 1 1 1 1 0 1 0 1 0 0⌋
H =
⌈1 0 0 0 0 0 0 0 0 0 0 0 1 0 1 1 0⌉
|0 1 0 0 0 0 0 0 0 0 0 0 0 1 0 1 1|
|0 0 1 0 0 0 0 0 0 1 0 0 0 1 1 0 0|
|0 0 0 1 0 0 0 0 0 0 1 0 0 0 1 1 0|
|0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 1 1|
|0 0 0 0 0 1 0 0 0 1 0 0 1 1 0 0 0|
|0 0 0 0 0 0 1 0 0 0 1 0 0 1 1 0 0|
|0 0 0 0 0 0 0 1 0 0 0 1 0 0 1 1 0|
⌊0 0 0 0 0 0 0 0 1 0 0 0 1 0 0 1 1⌋
A(x) = 1 + 3x^3 + 3x^4 + 7x^5 + 21x^6 + 40x^7 + 53x^8 + 53x^9 + 40x^10 + 21x^11 + 7x^12 + 3x^13 + 3x^14 + x^17 tmax = 1
}
Besides outputting code properties we can also query them one by one. Some examples are shown below:
const auto A = C.get_weight_enumerator();
const size_t dmin = C.get_dmin();
const bool is_perfect = C.is_perfect();
if (C.is_polynomial()) {
std::cout << "C is polynomial ";
if (C.is_cyclic()) std::cout << " and cyclic ";
std::cout << " with generator polynomial " << C.get_gamma() << std::endl;
} else {
std::cout << "C is not polynomial" << std::endl;
}
C is polynomial with generator polynomial 1 + x^4 + x^7 + x^8
Some properties are only defined in special cases. One example is the Bhattacharryya bound for binary codes. It requires the Bhattacharyya parameter of a binary channel:
if constexpr (std::is_same_v<F, Fp<2>>) {
const double pe = 0.1;
BSC channel(pe);
std::cout << "Word error probability for BSC with pe=" << pe << " is at most "
<< C.Bhattacharyya_bound(channel.get_Bhattacharyya_param()) << std::endl;
}
Word error probability of ML decoding for BSC with pe=0.01 is at most 0.0435097
CECCO can calculate minimal trellises for codes over arbitrary finite fields. These trellises are used by the Viterbi and BCJR decoding algorithms (cf. the demo on soft decoding). They can be exported in TiKZ format for use in LaTeX documents:
auto T = C.get_minimal_trellis();
T.export_as_tikz("trellis.tex");
For general linear codes, CECCO supports a variety of hard-input decoders. Standard array-based
bounded-distance BD and ML decoding with and without erasures (dec_BD, dec_BD_EE, dec_ML, dec_ML_EE),
hard-input Viterbi ML decoding with and without erasures (dec_Viterbi, dec_Viterbi_EE), and pseudo-soft-input GMD decoding
(dec_GMD). Note that the EE decoders as well as GMD decoding require
#define CECCO_ERASURE_SUPPORT.
First, we generate a random codeword of C':
auto u = Vector<F>(k).randomize();
std::cout << "Random message: " << u << std::endl;
auto c = C.enc(u);
Random message: (0, 1, 0, 0, 1, 0, 1, 0, 1)
Then we transmit over a discrete memoryless channel (a BSC in case F = Fp<2>) for a
hard-input decoding scenario:
const double pe = 0.05;
SDMC<F> channel(pe);
// RX
auto r = channel(c);
BD decoding can fail with a CECCO::decoding_failure so we have to wrap the decoder in a
try/catch block:
Vector<F> c_est, u_est;
try {
c_est = C.dec_BD(r);
u_est = C.encinv(c_est);
std::cout << "BD decoding message estimate: " << u_est << std::endl;
if (u_est == u) {
std::cout << "... this is correct decoding." << std::endl;
} else {
std::cout << "... this is wrong decoding/a word error with " << dH(u, u_est)
<< " (message) symbol errors" << std::endl;
}
} catch (const decoding_failure& e) {
std::cout << "BD decoding failure!" << std::endl;
}
ML decoding can never fail:
c_est = C.dec_ML(r);
u_est = C.encinv(c_est);
std::cout << "ML decoding message estimate: " << u_est << std::endl;
if (u_est == u) {
std::cout << "... this is correct decoding." << std::endl;
} else {
std::cout << "... this is wrong decoding/a word error with " << dH(u, u_est) << " (message) symbol errors"
<< std::endl;
}
For and exemplary transmission, both decoders may provide the correct message estimate:
BD decoding message estimate: ( 0, 1, 0, 0, 1, 0, 1, 0, 1 )
... this is correct decoding.
ML decoding message estimate: ( 0, 1, 0, 0, 1, 0, 1, 0, 1 )
... this is correct decoding.
In case the received vector r contains more than tmax errors, different ML decoders may come up with different codeword estimates (due to random tie breaking or implementational details). However, it is guaranteed that the Hamming distance of both codeword estimates to r is identical. Thus, we can assert the following (since Viterbi is an alternative ML decoder):
assert(dH(r, c_est) == dH(r, C.dec_Viterbi(r)));
Transmitting over a discrete memoryless channel with erasures looks almost identical (using the corresponding EE decoder):
const double pe = 0.025;
const double px = 0.05;
SDMEC<F> channel(pe, px);
// RX
auto r = channel(c);
Vector<F> c_est, u_est;
try {
c_est = C.dec_BD_EE(r);
u_est = C.encinv(c_est);
std::cout << "BD EE decoding message estimate: " << u_est << std::endl;
if (u_est == u) {
std::cout << "... this is correct decoding." << std::endl;
} else {
std::cout << "... this is wrong decoding/a word error with " << dH(u, u_est)
<< " (message) symbol errors" << std::endl;
}
} catch (const decoding_failure& e) {
std::cout << "BD EE decoding failure!" << std::endl;
}
c_est = C.dec_ML_EE(r);
u_est = C.encinv(c_est);
std::cout << "ML EE decoding message estimate: " << u_est << std::endl;
if (u_est == u) {
std::cout << "... this is correct decoding." << std::endl;
} else {
std::cout << "... this is wrong decoding/a word error with " << dH(u, u_est) << " (message) symbol errors"
<< std::endl;
}
assert(dH(r, c_est) == dH(r, C.dec_Viterbi_EE(r)));
Besides the hard received vector r, pseudo-soft GMD decoding also requires a vector of reliabilities. Reliability 1 means full reliability, 0 means completely unreliable. One way to obtain r and the reliability vector is to calculate it based on the per-symbol LLR matrix obtained from transmitting the binary image of the codeword over a BI-AWGN channel. In this approach, we find the two F-symbols with the smallest cost (smallest entry in the per-symbol LLR matrix column) and then map their difference to the interval [0, 1] using the tanh function.
// note: this only works if F has characteristic 2
DEMUX<F, Fp<2>> demux;
const double EbNodB = 0;
BI_AWGN channel(EbNodB);
LLRCalculator<F> llrcalculator(channel);
// RX
auto y = channel(demux(c));
auto llrs = llrcalculator(y);
std::cout << "Per-symbol LLR matrix: " << std::endl << llrs << std::endl;
auto r = Vector<F>(n);
auto reliabilities = Vector<double>(n);
for (size_t i = 0; i < n; ++i) {
size_t best_label = 0; // initialization, F(0) is best
double best = 0.0;
double second = std::numeric_limits<double>::infinity();
for (size_t a = 1; a < F::get_q(); ++a) {
const double cost = llrs(a - 1, i);
if (cost < best) {
second = best;
best = cost;
best_label = a;
} else if (cost < second) {
second = cost;
}
}
// hard decision: most probable symbol
r.set_component(i, F(best_label));
// reliability in [0, 1]: best-vs-second through tanh
const double lambda = second - best;
reliabilities.set_component(i, std::tanh(lambda / 2.0));
}
std::cout << "Hard decision r: " << std::endl << r << std::endl;
std::cout << "Reliabilities: " << std::endl << reliabilities << std::endl;
try {
auto c_est = C.dec_GMD(r, reliabilities);
auto u_est = C.encinv(c_est);
std::cout << "GMD decoding message estimate: " << u_est << std::endl;
if (u_est == u) {
std::cout << "... this is correct decoding." << std::endl;
} else {
std::cout << "... this is wrong decoding/a word error with " << dH(u, u_est) << " (message) symbol errors"
<< std::endl;
}
} catch (const decoding_failure& e) {
std::cout << "GMD decoding failure: " << e.what() << std::endl;
}
We stress that GMD decoding can fail in case all of its internal BD EE decoders fail.
A complete, compilable demo is shown below:
#include <iostream>
#define CECCO_ERASURE_SUPPORT
#include "cecco.hpp"
using namespace CECCO;
int main(void) {
using F = Ext<Fp<2>, MOD{1, 1, 1}>;
const size_t n = 13;
const size_t k = 6;
auto G = ZeroMatrix<F>(k, n);
do {
G.randomize();
} while (G.rank() < k);
auto Cp = LinearCode(n, k, G);
// std::cout << showbasic << Cp << std::endl;
// std::cout << showmost << Cp << std::endl;
std::cout << showall << Cp << std::endl;
auto C = Cp.get_equivalent_code_in_standard_form();
assert(C.is_equivalent(Cp));
std::cout << showall << C << std::endl;
auto Cd = C.get_dual();
std::cout << showall << Cd << std::endl;
const auto A = C.get_weight_enumerator();
const size_t dmin = C.get_dmin();
const bool is_perfect = C.is_perfect();
if (C.is_polynomial()) {
std::cout << "C is polynomial ";
if (C.is_cyclic()) std::cout << " and cyclic ";
std::cout << " with generator polynomial " << C.get_gamma() << std::endl;
} else {
std::cout << "C is not polynomial" << std::endl;
}
if constexpr (std::is_same_v<F, Fp<2>>) {
const double pe = 0.1;
BSC channel(pe);
std::cout << "Word error probability for BSC with pe=" << pe << " is at most "
<< C.Bhattacharyya_bound(channel.get_Bhattacharyya_param()) << std::endl;
}
auto T = C.get_minimal_trellis();
T.export_as_tikz("trellis.tex");
// TX
auto u = Vector<F>(k).randomize();
std::cout << "Random message: " << u << std::endl;
auto c = C.enc(u);
{
const double pe = 0.05;
SDMC<F> channel(pe);
// RX
auto r = channel(c);
Vector<F> c_est, u_est;
try {
c_est = C.dec_BD(r);
u_est = C.encinv(c_est);
std::cout << "BD decoding message estimate: " << u_est << std::endl;
if (u_est == u) {
std::cout << "... this is correct decoding." << std::endl;
} else {
std::cout << "... this is wrong decoding/a word error with " << dH(u, u_est)
<< " (message) symbol errors" << std::endl;
}
} catch (const decoding_failure& e) {
std::cout << "BD decoding failure!" << std::endl;
}
c_est = C.dec_ML(r);
u_est = C.encinv(c_est);
std::cout << "ML decoding message estimate: " << u_est << std::endl;
if (u_est == u) {
std::cout << "... this is correct decoding." << std::endl;
} else {
std::cout << "... this is wrong decoding/a word error with " << dH(u, u_est) << " (message) symbol errors"
<< std::endl;
}
assert(dH(r, c_est) == dH(r, C.dec_Viterbi(r)));
}
{
const double pe = 0.025;
const double px = 0.05;
SDMEC<F> channel(pe, px);
// RX
auto r = channel(c);
Vector<F> c_est, u_est;
try {
c_est = C.dec_BD_EE(r);
u_est = C.encinv(c_est);
std::cout << "BD EE decoding message estimate: " << u_est << std::endl;
if (u_est == u) {
std::cout << "... this is correct decoding." << std::endl;
} else {
std::cout << "... this is wrong decoding/a word error with " << dH(u, u_est)
<< " (message) symbol errors" << std::endl;
}
} catch (const decoding_failure& e) {
std::cout << "BD EE decoding failure!" << std::endl;
}
c_est = C.dec_ML_EE(r);
u_est = C.encinv(c_est);
std::cout << "ML EE decoding message estimate: " << u_est << std::endl;
if (u_est == u) {
std::cout << "... this is correct decoding." << std::endl;
} else {
std::cout << "... this is wrong decoding/a word error with " << dH(u, u_est) << " (message) symbol errors"
<< std::endl;
}
assert(dH(r, c_est) == dH(r, C.dec_Viterbi_EE(r)));
}
// note: this only works if F has characteristic 2
DEMUX<F, Fp<2>> demux;
const double EbNodB = 0;
BI_AWGN channel(EbNodB);
LLRCalculator<F> llrcalculator(channel);
// RX
auto y = channel(demux(c));
auto llrs = llrcalculator(y);
std::cout << "Per-symbol LLR matrix: " << std::endl << llrs << std::endl;
auto r = Vector<F>(n);
auto reliabilities = Vector<double>(n);
for (size_t i = 0; i < n; ++i) {
size_t best_label = 0; // initialization, F(0) is best
double best = 0.0;
double second = std::numeric_limits<double>::infinity();
for (size_t a = 1; a < F::get_q(); ++a) {
const double cost = llrs(a - 1, i);
if (cost < best) {
second = best;
best = cost;
best_label = a;
} else if (cost < second) {
second = cost;
}
}
// hard decision: most probable symbol
r.set_component(i, F(best_label));
// reliability in [0, 1]: best-vs-second through tanh
const double lambda = second - best;
reliabilities.set_component(i, std::tanh(lambda / 2.0));
}
std::cout << "Hard decision r: " << std::endl << r << std::endl;
std::cout << "Reliabilities: " << std::endl << reliabilities << std::endl;
try {
auto c_est = C.dec_GMD(r, reliabilities);
auto u_est = C.encinv(c_est);
std::cout << "GMD decoding message estimate: " << u_est << std::endl;
if (u_est == u) {
std::cout << "... this is correct decoding." << std::endl;
} else {
std::cout << "... this is wrong decoding/a word error with " << dH(u, u_est) << " (message) symbol errors"
<< std::endl;
}
} catch (const decoding_failure& e) {
std::cout << "GMD decoding failure: " << e.what() << std::endl;
}
return 0;
}