Reed–Solomon error correction
Reed–So 4}} (t 3 x^{2} + 2 x + 1}}, then the codeword is calculated as follows.
Errors in transmission might cause this to be received instead.
The syndromes are calculated by evaluating r at powers of α.
To correct the errors, first use the Berlekamp–Massey algorithm to calculate the error locator polynomial.
n | S_{n+1} | d | C | B | b | m |
---|---|---|---|---|---|---|
0 | 732 | 732 | 197 x + 1 | 1 | 732 | 1 |
1 | 637 | 846 | 173 x + 1 | 1 | 732 | 2 |
2 | 762 | 412 | 634 x^{2} + 173 x + 1 | 173 x + 1 | 412 | 1 |
3 | 925 | 576 | 329 x^{2} + 821 x + 1 | 173 x + 1 | 412 | 2 |
The final value of C is the error locator polynomial, Λ(x). The zeros can be found by trial substitution. They are x_{1} = 757 = 3^{−3} and x_{2} = 562 = 3^{−4}, corresponding to the error locations. To calculate the error values, apply the Forney algorithm.
Subtracting e_{1}x^{3} and e_{2}x^{4} from the received polynomial r reproduces the original codeword s.
Contents
Euclidean decoder[edit]
Another iterative method for calculating both the error locator polynomial and the error value polynomial is based on Sugiyama's adaptation of the Extended Euclidean algorithm .
Define S(x), Λ(x), and Ω(x) for t syndromes and e errors:
The key equation is:
For t = 6 and e = 3:
The middle terms are zero due to the relationship between Λ and syndromes.
The extended Euclidean algorithm can find a series of polynomials of the form
A_{i}(x) S(x) + B_{i}(x) x^{t} = R_{i}(x)
where the degree of R decreases as i increases. Once the degree of R_{i}(x) < t/2, then
A_{i}(x) = Λ(x)
B_{i}(x) = -Q(x)
R_{i}(x) = Ω(x).
B(x) and Q(x) don't need to be saved, so the algorithm becomes:
- R_{−1} = x^{t}
- R_{0} = S(x)
- A_{−1} = 0
- A_{0} = 1
- i = 0
- while degree of R_{i} >= t/2
- i = i + 1
- Q = R_{i-2} / R_{i-1}
- R_{i} = R_{i-2} - Q R_{i-1}
- A_{i} = A_{i-2} - Q A_{i-1}
to set low order term of Λ(x) to 1, divide Λ(x) and Ω(x) by A_{i}(0):
- Λ(x) = A_{i} / A_{i}(0)
- Ω(x) = R_{i} / A_{i}(0)
A_{i}(0) is the constant (low order) term of A_{i}.
Example[edit]
Using the same data as the Berlekamp Massey example above:
i | R_{i} | A_{i} |
---|---|---|
-1 | 001 x^{4} + 000 x^{3} + 000 x^{2} + 000 x + 000 | 000 |
0 | 925 x^{3} + 762 x^{2} + 637 x + 732 | 001 |
1 | 683 x^{2} + 676 x + 024 | 697 x + 396 |
2 | 673 x + 596 | 608 x^{2} + 704 x + 544 |
- Λ(x) = A_{2} / 544 = 329 x^{2} + 821 x + 001
- Ω(x) = R_{2} / 544 = 546 x + 732
Decoder using discrete Fourier transform[edit]
A discrete Fourier transform can be used for decoding. To avoid conflict with syndrome names, let c(x) = s(x) the encoded codeword. r(x) and e(x) are the same as above. Define C(x), E(x), and R(x) as the discrete Fourier transforms of c(x), e(x), and r(x). Since r(x) = c(x) + e(x), and since a discrete Fourier transform is a linear operator, R(x) = C(x) + E(x).
Transform r(x) to R(x) using discrete Fourier transform. Since the calculation for a discrete Fourier transform is the same as the calculation for syndromes, t coefficients of R(x) and E(x) are the same as the syndromes:
Use through as syndromes (they're the same) and generate the error locator polynomial using the methods from any of the above decoders.
Let v = number of errors. Generate E(x) using the known coefficients to , the error locator polynomial, and these formulas
Then calculate C(x) = R(x) - E(x) and take the inverse transform of C(x) to produce c(x).
Decoding beyond the error-correction bound[edit]
The Singleton bound states that the minimum distance d of a linear block code of size (n,k) is upper-bounded by n − k + 1. The distance d was usually understood to limit the error-correction capability to ⌊d/2⌋. The Reed–Solomon code achieves this bound with equality, and can thus correct up to ⌊(n − k + 1)/2⌋ errors. However, this error-correction bound is not exact.
In 1999, Madhu Sudan and Venkatesan Guruswami at MIT published "Improved Decoding of Reed–Solomon and Algebraic-Geometry Codes" introducing an algorithm that allowed for the correction of errors beyond half the minimum distance of the code. It applies to Reed–Solomon codes and more generally to algebraic geometric codes. This algorithm produces a list of codewords (it is a list-decoding algorithm) and is based on interpolation and factorization of polynomials over and its extensions.
Soft-decoding[edit]
The algebraic decoding methods described above are hard-decision methods, which means that for every symbol a hard decision is made about its value. For example, a decoder could associate with each symbol an additional value corresponding to the channel demodulator's confidence in the correctness of the symbol. The advent of LDPC and turbo codes, which employ iterated soft-decision belief propagation decoding methods to achieve error-correction performance close to the theoretical limit, has spurred interest in applying soft-decision decoding to conventional algebraic codes. In 2003, Ralf Koetter and Alexander Vardy presented a polynomial-time soft-decision algebraic list-decoding algorithm for Reed–Solomon codes, which was based upon the work by Sudan and Guruswami. In 2016, Steven J. Franke and Joseph H. Taylor published a novel soft-decision decoder.
Matlab Example[edit]
Encoder[edit]
Here we present a simple Matlab implementation for an encoder.
function [ encoded ] = rsEncoder( msg, m, prim_poly, n, k )
%RSENCODER Encode message with the Reed-Solomon algorithm
% m is the number of bits per symbol
% prim_poly: Primitive polynomial p(x). Ie for DM is 301
% k is the size of the message
% n is the total size (k+redundant)
% Example: msg = uint8('Test')
% enc_msg = rsEncoder(msg, 8, 301, 12, numel(msg));
% Get the alpha
alpha = gf(2, m, prim_poly);
% Get the Reed-Solomon generating polynomial g(x)
g_x = genpoly(k, n, alpha);
% Multiply the information by X^(n-k), or just pad with zeros at the end to
% get space to add the redundant information
msg_padded = gf([msg zeros(1, n-k)], m, prim_poly);
% Get the remainder of the division of the extended message by the
% Reed-Solomon generating polynomial g(x)
[~, reminder] = deconv(msg_padded, g_x);
% Now return the message with the redundant information
encoded = msg_padded - reminder;
end
% Find the Reed-Solomon generating polynomial g(x), by the way this is the
% same as the rsgenpoly function on matlab
function g = genpoly(k, n, alpha)
g = 1;
% A multiplication on the galois field is just a convolution
for k = mod(1 : n-k, n)
g = conv(g, [1 alpha .^ (k)]);
end
end
Decoder[edit]
Now the decoding part:
function [ decoded, error_pos, error_mag, g, S ] = rsDecoder( encoded, m, prim_poly, n, k )
%RSDECODER Decode a Reed-Solomon encoded message
% Example:
% [dec, ~, ~, ~, ~] = rsDecoder(enc_msg, 8, 301, 12, numel(msg))
max_errors = floor((n-k)/2);
orig_vals = encoded.x;
% Initialize the error vector
errors = zeros(1, n);
g = [];
S = [];
% Get the alpha
alpha = gf(2, m, prim_poly);
% Find the syndromes (Check if dividing the message by the generator
% polynomial the result is zero)
Synd = polyval(encoded, alpha .^ (1:n-k));
Syndromes = trim(Synd);
% If all syndromes are zeros (perfectly divisible) there are no errors
if isempty(Syndromes.x)
decoded = orig_vals(1:k);
error_pos = [];
error_mag = [];
g = [];
S = Synd;
return;
end
% Prepare for the euclidean algorithm (Used to find the error locating
% polynomials)
r0 = [1, zeros(1, 2*max_errors)]; r0 = gf(r0, m, prim_poly); r0 = trim(r0);
size_r0 = length(r0);
r1 = Syndromes;
f0 = gf([zeros(1, size_r0-1) 1], m, prim_poly);
f1 = gf(zeros(1, size_r0), m, prim_poly);
g0 = f1; g1 = f0;
% Do the euclidean algorithm on the polynomials r0(x) and Syndromes(x) in
% order to find the error locating polynomial
while true
% Do a long division
[quotient, remainder] = deconv(r0, r1);
% Add some zeros
quotient = pad(quotient, length(g1));
% Find quotient*g1 and pad
c = conv(quotient, g1);
c = trim(c);
c = pad(c, length(g0));
% Update g as g0-quotient*g1
g = g0 - c;
% Check if the degree of remainder(x) is less than max_errors
if all(remainder(1:end - max_errors) == 0)
break;
end
% Update r0, r1, g0, g1 and remove leading zeros
r0 = trim(r1); r1 = trim(remainder);
g0 = g1; g1 = g;
end
% Remove leading zeros
g = trim(g);
% Find the zeros of the error polynomial on this galois field
evalPoly = polyval(g, alpha .^ (n-1 : -1 : 0));
error_pos = gf(find(evalPoly == 0), m);
% If no error position is found we return the received work, because
% basically is nothing that we could do and we return the received message
if isempty(error_pos)
decoded = orig_vals(1:k);
error_mag = [];
return;
end
% Prepare a linear system to solve the error polynomial and find the error
% magnitudes
size_error = length(error_pos);
Syndrome_Vals = Syndromes.x;
b(:, 1) = Syndrome_Vals(1:size_error);
for idx = 1 : size_error
e = alpha .^ (idx*(n-error_pos.x));
err = e.x;
er(idx, :) = err;
end
% Solve the linear system
error_mag = (gf(er, m, prim_poly) \ gf(b, m, prim_poly))';
% Put the error magnitude on the error vector
errors(error_pos.x) = error_mag.x;
% Bring this vector to the galois field
errors_gf = gf(errors, m, prim_poly);
% Now to fix the errors just add with the encoded code
decoded_gf = encoded(1:k) + errors_gf(1:k);
decoded = decoded_gf.x;
end
% Remove leading zeros from galois array
function gt = trim(g)
gx = g.x;
gt = gf(gx(find(gx, 1) : end), g.m, g.prim_poly);
end
% Add leading zeros
function xpad = pad(x,k)
len = length(x);
if (len<k)
xpad = [zeros(1, k-len) x];
end
end
See Also on BitcoinWiki[edit]
- BCH code
- Chien search
- Berlekamp–Massey algorithm
- Forward error correction
- Berlekamp–Welch algorithm
- Folded Reed–Solomon code
Source[edit]