bchpoly.m

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

function [pg, pm, cs, h, t] = bchpoly(n, k)
%BCHPOLY Produces GF(2) BCH code generator polynomial.
%       BCHPOLY lists the codeword length, message length, and
%       error-correction capability for BCH code in a figure. The codeword
%       lengths listed are 7, 15, 31, 63, 127, 255, and 511.
%
%       Pg = BCHPOLY outputs the list to a vector Pg. The first column is
%       the code word length, the second column is the message length, and the
%       third column is the error-correction capability.
%
%       Pg = BCHPOLY(N) outputs the code word length equal or larger than N,
%       its message lengths, and error-correction capability on the first,
%       second and third column of Pg respectively. When N > 511, the error
%       control capability is not listed.
%
%       Pg = BCHPOLY(N, K) produces BCH code generator polynomial Pg based on
%       two positive numbers: the codeword length N and the message length K.
%       The codeword length must be the solution of 2^M - 1, where M is a 
%       integer larger than 2. K should be a legible message length, which can
%       be found by using Pg = BCHPOLY(N) to get a list. 
%       When N is a binary vector instead of a positive number, N is a
%       primitive polynomial generating the field GF(2^DEG), where DEG is the
%       degree of polynomial N(X). The codeword length in this case should be 
%       2^DEG - 1.
%
%       [Pg, Pm] = BCHPOLY(N, K) outputs minimum polynomials Pm in the field
%       GF(2^M). Each row of Pm is a minimum polynomial.
%
%       [Pg, Pm, CS] = BCHPOLY(N, K) list the cosets that are the roots of
%       the polynomial Pg.
%
%       [Pg, Pm. CS, H] = BCHPOLY(N, K) outputs the parity-check matrix H.
%
%       [Pg, Pm. CS, H, T] = BCHPOLY(N, K) outputs the error-control
%       capability of the BCH code.
%
%       See also: CYCLGEN, GFMINPOL.

%       Wes Wang 7/1/94, 9/30/95
%       Copyright (c) 1995-96 by The MathWorks, Inc.
%       $Revision: 1.1 $  $Date: 1996/04/01 17:52:13 $

% parameter assignment
n_lst = [3    7   1;
         4   15   3;
         5   31   5;
         6   63  11;
         7  127  17;
         8  255  33;
         9  511  57;
        10 1023 105];

% known error-correction capability
t_lst = [1,...
         [1:3],...
         [1:3],5 7,...
         [1:7],10 11 13 15,...
         [1:7], 9,10 11 13 14 15 21 23 27 31,...
         [1:15], 18 19 [21:23],[25:27],[29:31],42 43 45 47 55 59 63,...
         [1:16],[18:23],[25:31],[36:39],[41:43],[45:47],51,[53:55],58 59,[61:63],...
            85 87 91 93 95 109 111 119 121,...
         [1:31],[34:39],[41:47],[49:55],[57:63],[73:75],[77:79],82,83,85+[0:2],...
         89+[0:2],93+[0:2],102,103,106,107,109,110,111,115,117+[0:2],122,123,...
        125+[0:2],170,171,173,175,181,183,187,189,191,219,223,239,247,255];

% give a list of all BCH code 
% code word length / message length / error-correction capability
% When there is no output parameter, draw a figure.
% When there is output parameter, output the list.

if nargin < 1
    if nargout < 1
        % title of the figure.
        yy = 'BCH Code Generated by Default Primitive Polynomial';

        % In case the figure exist
        hand = get(0, 'Child');
        if ~isempty(hand)
            for i = 1:length(hand)
                if strcmp(yy, get(hand(i), 'Name'))
                    set(0,'Current',hand(i));
                    if length(get(get(hand(i),'Child'),'Child')) == 278
                        set(hand(i), 'Position', [10 10 700 700])
                        return;
                    else
                        close(hand(i));
                    end
                end;
            end;
        end;

        % figure frame construction
        h = figure('Position',    [10 10 700 700],...
                   'NumberTitle', 'off',....
                   'Clipping',    'off',...
                   'Name',        yy);
        hght = 43;
        x = axes('Xlim', [0 9], 'Ylim', [0 hght], 'Visible','off',...
            'Nextplot','add',...
            'Clipping','off');
        set(x, 'visible','off');
        text(0, -2, 'N: code word length; K: message length; T: error-correction capability');
        for i=1:3
            text(3*(i-1),  hght+1,'   N')
            text(3*(i-1)+1,hght+1,'   K')
            text(3*(i-1)+2,hght+1,'   T')
        end;
        plot([0, 9],   [hght, hght],  'b-');
        plot([0, 9],   [hght, hght]+2,'b-');
        plot([0, 9],   [-1 -1],       'b-');
        plot([2.5,2.5],[-1,hght+2],   'b-');
        plot([5.5,5.5],[-1,hght+2],   'b-');
    end;

    % list of all numbers in the figure or output to variable pg.
    next_coset = 1;
    n_m = 0;
    n_base = 0;

    for i = 0 : 126
        y_p = rem(i,hght);
        x_p = floor(i / hght) * 3;

        if ((next_coset) | (y_p == 0)) & (nargout < 1)
            drawnow;
            if next_coset
                xx = num2str(n_lst(n_m+1, 2));
            else
                xx = num2str(n_lst(n_m, 2));
            end;
            xx = [setstr(ones(1, 4-length(xx))*32) xx];            
            xx=text(x_p, hght-y_p-1, xx);
            set(xx,'Color',[1 1 0]);
        end;

        if next_coset
            n_m = n_m + 1;
            cs = gfcosets(n_lst(n_m, 1));
            next_coset = 0;
        end;
        if (nargout < 1)
              xx = num2str(n_lst(n_m, 2) - sum(sum(~isnan(cs(2:2+i-n_base,:)))));
            xx = [setstr(ones(1, 4-length(xx))*32) xx];            
            xx=text(x_p+1, hght-y_p-1, xx);
            set(xx,'Color',[1 0 1]);
            xx = num2str(t_lst(i+1));
            xx = [setstr(ones(1, 4-length(xx))*32) xx];
            xx=text(x_p+2, hght-y_p-1, xx);
            set(xx,'Color',[0 1 1]);
        else
            pg = [pg;
                  [n_lst(n_m, 2),...
                   n_lst(n_m, 2)-sum(sum(~isnan(cs(2:2+i-n_base,:)))),...
                   t_lst(i+1)]];
        end;
        if n_base + n_lst(n_m, 3) == i + 1
            next_coset = 1;
            n_base = n_base + n_lst(n_m, 3);
        end;
    end;
    set(h, 'Nextplot','new');
    return;
end;

% When specific code word length is specified, output a list of the
% code word length / message length / error-correction capability

if length(n) > 1
    m = gftrunc(n);
    n = 2^length(m) - 1;
    df_flag = 0;
else
    df_flag = 1;
end;
indx = find(n <= n_lst(:, 2));
if isempty(indx)
    % beyong the range
    dim = 0;
    k = 11;
    while (dim == 0)
        if (2^k-1 >= n)
            dim = k;
        else
            k = k + 1;
        end;
    end;
else
    dim = n_lst(indx(1), 1);
end;
n = 2^dim - 1;


if nargin < 2
    cs = gfcosets(dim);
    [n_cs, m_cs] = size(cs);
    if nargout < 1
        yy = ['BCH Code Generated by ', num2str(dim), 'th Order Primitive Polynomial'];

        % in case the figure exist
        hand = get(0, 'Child');
        if ~isempty(hand)
            for i = hand(:)'
                if strcmp(yy, get(i, 'Name'))
                    set(0, 'Current', i);
                    sz = get(i, 'UserData');
                    if length(get(get(i, 'Child'), 'Child')) == sz(1)
                        set(i, 'Position', sz(2:5));
                        return;
                    else
                        close(i);
                    end;
                end;
            end;
        end;
        % figure frame construction

        x_a = 1 + (dim>=8) + (dim>10);
        if x_a < 3
            hght = ceil((n_cs - 2) / x_a);
        else
            hght = ceil((n_cs - 2) / x_a / 2);
        end;
        h = figure('Position', [10, 10, 240*x_a, hght*12+100],...
                   'NumberTitle', 'off',...
                   'Clipping',    'off',...
                   'Name',        yy);
        x = axes('Xlim', [0 2*x_a], 'Ylim', [-1 hght+2], 'Visible','off',...
            'Nextplot','add',...
            'Clipping','off');
        set(x, 'visible','off');
        if n <= 1023
            text(0, -1, 'K: message length; T: error-correction capability.');
        else
            text(0, -1, 'K: message length.');
        end;

        for i=1 : x_a
            text(2*(i-1),hght+1,'   K')
            if (dim <= 10)
                text(2*(i-1)+1,hght+1,'   T')
            else
                text(2*(i-1)+1,hght+1,'   K')
            end;
        end;
        text(0, hght+2, ['Code word length', num2str(n)]);

        plot([0, 2*x_a],   [hght, hght]+.5,  'b-');
        plot([0, 2*x_a],   [hght, hght]+2.5,'b-');
        plot([0, 2*x_a],   [-.5 -.5],       'b-');

        if x_a > 1
            plot([1.5,1.5],[-1,hght+2],   'b-');