% repetition code and convolutional code using Belief Propagation Inference% decoding on AWGN channel using matlab or octave if sparse matrix support% is added.% V 1.1 04/15/05% Mail comments to J.C. Olivier corne.olivier@eng.up.ac.za; and BP decoder % by I. Kozintsev (www.kozintsev.net) with permission.%% note that no special consideration was given for terminating effects% as is usually done by trellis termination in conventional VA trellis % based inference.  The decoded BER for the convolutional code can be improved% slightly if those effects are included in the BP inference.%% Type demo to execute, and results are on the % monitor, in the form [block number,raw BER,decoded BER]clear alln_bits = 100; % number of information bits per block (frame)SNR = 0; % in dB, SNR = Eb/No for BPSK n_bl = 1000; % number of blocks to simualtechoice = 2; % 1 for repetition code and 2 for GSM convolutional code% compute channel parameterslog_off = 10^(SNR/10);noise_std = sqrt(1/(2*log_off));% conv encoder (from GSM)gg = [1  0  1  1  0  1  1      1  1  0  0  1  0  1      1  1  1  1  1  0  1];% initialize variablesber = 0; FER = 0; BER_dec = 0;% make conv code generator matrix and parity check matrices G and H Q = n_bits;  M=2*Q; errors_acc = zeros(1,Q); P1 = zeros(Q);  for loop=1:7     P1 = P1 + diag(gg(1,loop)*ones(1,Q-loop+1),loop-1); end P2 = zeros(Q); for loop=1:7     P2 = P2 + diag(gg(2,loop)*ones(1,Q-loop+1),loop-1); end P3 = zeros(Q); for loop=1:7     P3 = P3 + diag(gg(3,loop)*ones(1,Q-loop+1),loop-1); end G = [P1 P2 P3]; % generator matrix rowpick = 2; % now use elementary operations to get into systematic form G = [P|I]; for loop = 2*Q+2:3*Q     for loop2=1:rowpick-1         if G(loop2,loop) > 0             G(loop2,:) = mod(G(loop2,:) + G(rowpick,:),2);         end     end     rowpick = rowpick + 1; end P = G(1:Q,1:2*Q); % make parity check matrix H = [eye(M)      P].'; verify = mod(G*H',2);  % make repeat code G and H P1 = eye(Q); P2 = eye(Q); P3 = eye(Q); G_r = [P1 P2 P3]; % generator matrix P_r = [P1 P2]; % make parity check matrix H_r = [eye(M)        P_r].'; if choice == 1 % repeat code, else convolutional code   G = G_r;   H = H_r;end    % initialize variables ber = 0; BER_dec = 0;% over all blocksfor blocks=1:n_bl    % make random bipolar binary info sequence vector s  s = [round(rand(n_bits,1)).'];   % encode  input_bits_coded = mod(s*G,2); % block code conv encoding    % make received vector in AWGN  awgn = noise_std*(randn(length(input_bits_coded),1)).'; % real noise  %awgn = noise_std*(randn(length(s),1) + i*randn(length(s),1)); % complex noise   y = (2*input_bits_coded-1) + awgn; % received in the presence of noise   Prob(:,2) = 1./(1+exp(-2*y/noise_std^2)).';  Prob(:,1) = 1-Prob(:,2);  llr = log(Prob(:,2)./Prob(:,1));  % log likelyhood    err = sum(xor(input_bits_coded,((sign(llr)+1)/2).'));  ber = ber + err; % Raw channel BER  BER = ber/Q/blocks/3;    % decode block (frame)  [z_hat, success, k] = belief_prop_I_Kozintsev(Prob(:,1).',Prob(:,2).',H,100); % BP decoder    s_hat = z_hat(2*Q+1:3*Q).';   errors = xor(s(1:Q),s_hat);    errors_acc = errors_acc + errors;  err_dec = sum(errors);  BER_dec = BER_dec + err_dec/n_bits;   [blocks,BER,BER_dec/blocks] end