Code covered by the BSD License

# This is an implementation of (7,4) hamming code using belief propagation

### Shuang Wang (view profile)

An implementation of (7,4) hamming code using belief propagation

```clear;

pchkArray = [1, 1, 1, 0, 1, 0, 0;
0, 1, 1, 1, 0, 1, 0;
1, 0, 1, 1, 0, 0, 1];      % parity check matrix for (7,4) hamming code

noiseCodeword =  [1, 0, 1, 1, 0, 0, 0];      % input code word with 1 bit error
decodeword = zeros(1, numel(noiseCodeword));                  % array for storing decoded word
edgeArray = [eye(7); pchkArray];        % build the edge array which stores the connection between variable nodes and factor nodes
[i, j, s] = find(edgeArray);            % find the index of edgeArray;
edgeMap = sparse(i, j, [1:length(s)]);  % map edge id to node id;

max_iter = 10;                          % max number of iteration for bp;
p = 0.01;                               % error probability for BSC

[m,n] = size(pchkArray);                % get number of factor nodes and variable nodes, where the num of factor nodes is equal to n + m, and num of factor nodes is n;
numOfEdge = sum(edgeArray(:));          % get number of edges between variable nodes and factor nodes
numOfLabel = 2;                         % set the num of bp labels, where 2 means we only use 0 and 1 as the label in our bp algorithm.
LabelValue = [0, 1];                    % set the value
messageF2V = ones(numOfEdge, numOfLabel);  % message from variable nodes to factor nodes
messageV2F = ones(numOfEdge, numOfLabel);  % message from factor nodes to variable
belief     = ones(n, numOfLabel);          % belief of each variable node.
epsVal = 10E-5;                            % a small value to make sure there is no zero denominator;

%% update factor node 1 to 7. since they only have one variable node neighbor, they only need to be update once.
for i = 1:n
for j = i:i
edgeID = edgeMap(i,j); % find the edgeID of current neighbor j for factor node i.
for k = 1: numOfLabel
messageF2V(edgeID,k) = fun1(noiseCodeword(i), LabelValue(k), p);
end;
end;
end;

%% start bp iteration
for iter = 1:max_iter

%% update variable nodes 1 to 7
for j = 1:n
prodVal = messageF2V(edgeMap(j,j), :);  % get the message from Left side factor node;
edgeIDArray = full(edgeMap(find(edgeMap(j+1:n+m,j)~= 0) + j, j));  % find neighbor's edgeID;

% mutliply all the incoming message from its neighbor
for i = 1:length(edgeIDArray)
edgeID = edgeIDArray(i); % find the index of current neighbor i.
prodVal = prodVal .* messageF2V(edgeID,:);
end;

%send out the message, divide prodVal by the incoming message you
%get the sending out message;
for i = 1:length(edgeIDArray)
edgeID = edgeIDArray(i); % find the index of current neighbor i.
messageV2F(edgeID,:) = prodVal ./ messageF2V(edgeID,:);
% normalization
messageV2F(edgeID,:) = messageV2F(edgeID,:) / sum(messageV2F(edgeID,:));
end;

end;

%% update factor nodes  8 to 10
for i = n+1 : n+m
edgeIDArray = full(edgeMap(i,find(edgeMap(i, :)~= 0)));  % find neighbor's edgeID;
for j = 1: length(edgeIDArray)
edgeID = edgeIDArray(j); % find the index of current neighbor j.
for k = 1:numOfLabel
messageF2V(edgeID,k) = sumMsg(LabelValue(k), edgeIDArray, messageV2F, j);
end;
% normalization
messageF2V(edgeID,:) = messageF2V(edgeID,:) / sum(messageF2V(edgeID,:));
end;
end;

%% calculate belief
for j = 1:n
belief(j, :) = messageF2V(edgeMap(j,j), :);  % get the message from Left side factor node;
edgeIDArray = full(edgeMap(find(edgeMap(j+1:n+m,j)~= 0) + j, j));  % find neighbor's edgeID;
% mutliply all the incoming message
for i = 1 : length(edgeIDArray)
edgeID = edgeIDArray(i); % find the edgeID of current neighbor i for variable node j.
belief(j, :) = belief(j, :) .* messageF2V(edgeID,:);
end;
end;

%% get decoded codeword by belief
decodeword = fix(belief(:,2) ./ belief(:,1))';
decodeword(find(decodeword ~= 0)) = 1;
%% validate the decoded codeword
checkVal = mod(sum(pchkArray * decodeword'),2);
if(iter >= max_iter || checkVal == 0)
break;
end;

end;

if(checkVal ==0)
disp('decode successful');
else
disp('decode failed;')

end;
disp(['  input codeword: ', int2str(noiseCodeword)]);
disp(['decoded codeword: ', int2str(decodeword)]);

```