rscore.m

数字通信第四版原书的例程

function [msg, err, ccode] = rscore(code, k, tp, dim, pow_dim);
%RSCORE The core function in Reed-Solomon decode.
%       MSG = RSCORE(CODE, K, TP, M, POW_M, T2) decodes a single codeword
%       vector CODE using Reed-Solomon decoding technique. The message length
%       is K. The complete (and correct) list of all members in GF(2^M) is
%       in TP. The code word length is provided in POW_M, which equals to
%       2^M - 1. The decoding result is provided in the output variable MSG.
%
%       [MSG, ERR] = RSCORE(CODE, K, TP, M, POW_M, T2) outputs the error 
%       detected in the decoding.
%
%       [MSG, ERR, CCODE] = RSCORE(CODE, K, TP, M, POW_M, T2) outputs the
%       corrected codeword in CCODE.
%
%       NOTE. Different from all of the other encoding/decoding functions,
%       this function takes the exponential input instead of regular input for
%       processing. For example [-Inf, 0, 1, 2, ...] represents
%       [0 1 alpha, alpha^2, ...] in GF(2^m). There are 2^M elements in
%       GF(2^M). Hence, the input CODE represents 2^M * (2^M - 1) bits of
%       information. The decoded MSG represents 2^M * K bits of information.
%       For speeding computation, no error-check is placed in this function,
%       all input variables must be presented.
%       

%       Wes Wang 8/11/94, 10/11/95.
%       Copyright (c) 1995-96 by The MathWorks, Inc.
%       $Revision: 1.1 $  $Date: 1996/04/01 18:02:35 $

t2 = pow_dim - k;
t = floor(t2 / 2);
code=code(:)';

tp_num = tp * 2.^[0:dim-1]';
tp_inv(tp_num+1) = 0:pow_dim;

% (1) find syndrome
for sy_i = 1 : t2
    % The ith element of syndrome equals the result of dividing code(X) by alpha^i+X
    [tmp, syndrome(sy_i)] = gfdeconv(code, [sy_i, 0], tp);
end;

% (2) find error location polynomial
% For non-binary case, using Berlekamp's iterative method to do the computation.
err = 0;
[sigma, err] = errlocp(syndrome, t, tp, pow_dim, err, 1);

% (3) find solution for the polynomial
loc_err = zeros(1, pow_dim) - Inf;
num_err = length(sigma) - 1;
if num_err > t
    % in the case of more errors than possibly correction.
    err = 1;
end;
if (~err) & (num_err > 0)
    cnt_err = 0;
    pos_err = [];
    er_i = 0;
    while (cnt_err < num_err) & (er_i < pow_dim * dim)
        test_flag = sigma(1);
        for er_j = 1 : num_err
            if sigma(er_j + 1) >= 0
                % The following 6 lines is equivelent to
                % tmp = gfmul(er_i * er_j, sigma(er_j+1), tp);
                tmp = er_i * er_j;
                if (tmp < 0) | (sigma(er_j+1) < 0)
                    tmp = -1;
                else
                    tmp = rem(tmp + sigma(er_j + 1), pow_dim);
                end;

                test_flag = gfplus(test_flag, tmp, tp_num, tp_inv);
            end;
        end;
        if test_flag < 0
            cnt_err = cnt_err + 1;
            pos_err = [pos_err, rem(pow_dim-er_i, pow_dim)];
        end;
        er_i = er_i + 1;
    end;
    pos_err = rem(pow_dim+pos_err, pow_dim);
    pos_err = pos_err + 1; % shift one location because power zero is one.
%    loc_err(pos_err) = ones(1, cnt_err);
    err = num_err;
else
    if err
        err = -1;
    end;
end;

% (4) find the amplitude of the error
if (err > 0)
    % Construct Z(X) in (6.34)
    Z = 1;
    for am_i = 1 : num_err
        Z(am_i+1) = gfplus(syndrome(am_i), sigma(am_i + 1), tp_num, tp_inv);
        if am_i > 1
            for am_j = 1 : am_i - 1
                Z(am_i + 1) = gfplus(Z(am_i + 1), gfmul(sigma(am_j + 1), syndrome(am_i - am_j), tp), tp_num, tp_inv);
            end;
        end;
    end;

    % use pos_err here.
    pos_err_1 = pos_err - 1;    % considering position starting from zero.
    er_loc = [];
    for am_i = 1 : length(pos_err_1)
        num = 0;
        den = 0;
        pos_err_inv = rem(pow_dim - pos_err_1(am_i), pow_dim);
        for am_j = 1 : num_err
%            num = gfadd(num, gfmul(Z(am_j + 1), pos_err_inv * am_j, tp), tp);
            tmp = pos_err_inv * am_j;
            if (tmp < 0) | (Z(am_j + 1) < 0)
                tmp = -1;
            else
                tmp = rem(tmp + Z(am_j + 1), pow_dim);
            end;
            num = gfplus(num, tmp, tp_num, tp_inv);
            if am_i ~= am_j
%                den = gfmul(den, gfadd(0, gfmul(pos_err_1(am_j), pos_err_inv, tp), tp), tp);
                if (den < 0)
                    den = -1;
                else
                    if (pos_err_1(am_j) < 0) | (pos_err_inv < 0)
                        tmp = -1;
                    else
                        tmp = rem(pos_err_1(am_j) + pos_err_inv, pow_dim);
                    end;
                    tmp = gfplus(0, tmp, tp_num, tp_inv);
                    if (tmp < 0)
                        den = -1;
                    else
                        den = rem(den + tmp, pow_dim);
                    end;
                end;            
            end;
%
        end;
%        er_loc(am_i) = gfmul(num, pow_dim-den, tp);
        tmp = pow_dim - den;
        if (tmp < 0) | (num < 0)
            er_loc(am_i) = -1;
        else
            er_loc(am_i) = rem(tmp  + num, pow_dim);
        end;
%
    end;
    loc_err(pos_err) = er_loc;
end;

ccode = gfplus(loc_err, code, tp_num, tp_inv);
msg = ccode(pow_dim-k+1 : pow_dim);
%-- end of rscore--