IT++ Logo
Simulation of PCCCs in an AWGN channel

This program simulates Parallel Concatenated Convolutional Codes (PCCCs) of coding rate 1/3 using a turbo decoder with two SISO RSC modules.

Reference: S. Benedetto, D. Divsalar, G. Motorsi and F. Pollara, "A Soft-Input Soft-Output Maximum A posteriori (MAP) Module to Decode Parallel and Serial Concatenated Codes", TDA Progress Report, nov. 1996

#include "itpp/itcomm.h"
using namespace itpp;
using std::cout;
using std::endl;
using std::string;
int main(void)
{
//general parameters
double threshold_value = 10;
string map_metric="maxlogMAP";
ivec gen = "07 05";//octal form, feedback first
int constraint_length = 3;
int nb_errors_lim = 3000;
int nb_bits_lim = int(1e6);
int perm_len = (1<<14);//total number of bits in a block (with tail)
int nb_iter = 10;//number of iterations in the turbo decoder
vec EbN0_dB = "0:0.1:5";
double R = 1.0/3.0;//coding rate (non punctured PCCC)
double Ec = 1.0;//coded bit energy
//other parameters
int nb_bits = perm_len-(constraint_length-1);//number of bits in a block (without tail)
vec sigma2 = (0.5*Ec/R)*pow(inv_dB(EbN0_dB), -1.0);//N0/2
double Lc;//scaling factor
int nb_blocks;//number of blocks
int nb_errors;
ivec perm(perm_len);
ivec inv_perm(perm_len);
bvec bits(nb_bits);
int cod_bits_len = perm_len*gen.length();
bmat cod1_bits;//tail is added
bvec tail;
bvec cod2_input;
bmat cod2_bits;
int rec_len = int(1.0/R)*perm_len;
bvec coded_bits(rec_len);
vec rec(rec_len);
vec dec1_intrinsic_coded(cod_bits_len);
vec dec2_intrinsic_coded(cod_bits_len);
vec apriori_data(perm_len);//a priori LLR for information bits
vec extrinsic_coded(perm_len);
vec extrinsic_data(perm_len);
bvec rec_bits(perm_len);
int snr_len = EbN0_dB.length();
mat ber(nb_iter,snr_len);
ber.zeros();
register int en,n;
//Recursive Systematic Convolutional Code
cc.set_generator_polynomials(gen, constraint_length);//initial state should be the zero state
//BPSK modulator
BPSK bpsk;
//AWGN channel
AWGN_Channel channel;
//SISO modules
SISO siso;
siso.set_generators(gen, constraint_length);
siso.set_map_metric(map_metric);
//BER
BERC berc;
//Randomize generators
//main loop
for (en=0;en<snr_len;en++)
{
cout << "EbN0_dB = " << EbN0_dB(en) << endl;
channel.set_noise(sigma2(en));
Lc = -2/sigma2(en);//normalisation factor for intrinsic information (take into account the BPSK mapping)
nb_errors = 0;
nb_blocks = 0;
while ((nb_errors<nb_errors_lim) && (nb_blocks*nb_bits<nb_bits_lim))
{
//permutation
perm = sort_index(randu(perm_len));
//inverse permutation
inv_perm = sort_index(perm);
//bits generation
bits = randb(nb_bits);
//parallel concatenated convolutional code
cc.encode_tail(bits, tail, cod1_bits);//tail is added here to information bits to close the trellis
cod2_input = concat(bits, tail);
cc.encode(cod2_input(perm), cod2_bits);
for (n=0;n<perm_len;n++)//output with no puncturing
{
coded_bits(3*n) = cod2_input(n);//systematic output
coded_bits(3*n+1) = cod1_bits(n,0);//first parity output
coded_bits(3*n+2) = cod2_bits(n,0);//second parity output
}
//BPSK modulation (1->-1,0->+1) + AWGN channel
rec = channel(bpsk.modulate_bits(coded_bits));
//form input for SISO blocks
for (n=0;n<perm_len;n++)
{
dec1_intrinsic_coded(2*n) = Lc*rec(3*n);
dec1_intrinsic_coded(2*n+1) = Lc*rec(3*n+1);
dec2_intrinsic_coded(2*n) = 0.0;//systematic output of the CC is already used in decoder1
dec2_intrinsic_coded(2*n+1) = Lc*rec(3*n+2);
}
//turbo decoder
apriori_data.zeros();//a priori LLR for information bits
for (n=0;n<nb_iter;n++)
{
//first decoder (terminated trellis)
siso.rsc(extrinsic_coded, extrinsic_data, dec1_intrinsic_coded, apriori_data, true);
//interleave
apriori_data = extrinsic_data(perm);
//threshold
apriori_data = SISO::threshold(apriori_data, threshold_value);
//second decoder (unterminated trellis)
siso.rsc(extrinsic_coded, extrinsic_data, dec2_intrinsic_coded, apriori_data);
//decision
apriori_data += extrinsic_data;//a posteriori information
rec_bits = bpsk.demodulate_bits(-apriori_data(inv_perm));//take into account the BPSK mapping
//count errors
berc.clear();
berc.count(bits, rec_bits.left(nb_bits));
ber(n,en) += berc.get_errorrate();
//deinterleave for the next iteration
apriori_data = extrinsic_data(inv_perm);
}//end iterations
nb_errors += int(berc.get_errors());//get number of errors at the last iteration
nb_blocks++;
}//end blocks (while loop)
//compute BER over all tx blocks
ber.set_col(en, ber.get_col(en)/nb_blocks);
}
it_file ff("pccc_bersim_awgn.it");
ff << Name("EbN0_dB") << EbN0_dB;
ff << Name("BER") << ber;
ff.close();
return 0;
}
Ordinary AWGN Channel for cvec or vec inputs and outputs.
Definition: channel.h:1089
void set_noise(double noisevar)
Set noise variance (for complex-valued channels the sum of real and imaginary parts)
Definition: channel.h:1094
Bit Error Rate Counter (BERC) Class.
void clear()
Clears the bit error counter.
double get_errors() const
Returns the counted number of bit errors.
double get_errorrate() const
Returns the estimated bit error rate.
void count(const bvec &in1, const bvec &in2)
Cumulative error counter.
BPSK modulator with real symbols.
Definition: modulator.h:877
void modulate_bits(const bvec &bits, vec &output) const
Modulate bits into BPSK symbols in complex domain.
Definition: modulator.cpp:310
void demodulate_bits(const vec &signal, bvec &output) const
Demodulate noisy BPSK symbols in complex domain into bits.
Definition: modulator.cpp:326
Mat< double > mat
Default Matrix Type.
Definition: mat.h:482
Automatic naming when saving.
Definition: itfile.h:429
A Recursive Systematic Convolutional Encoder/Decoder class.
void set_generator_polynomials(const ivec &gen, int constraint_length)
Set generator polynomials.
void encode(const bvec &input, bmat &parity_bits)
Encode a binary vector of inputs starting from zero state without adding of a tail.
void encode_tail(const bvec &input, bvec &tail, bmat &parity_bits)
Encode a binary vector of inputs and also adds a tail of K-1 zeros to force the encoder into the zero...
Soft Input Soft Output (SISO) modules.
Definition: siso.h:72
void set_map_metric(const std::string &in_MAP_metric)
Sets the metric for MAP algorithm (convolutional codes and multipath channels)
Definition: siso.h:574
void set_generators(const itpp::bmat &in_gen)
Sets convolutional code generator polynomials.
Definition: siso.h:590
void rsc(itpp::vec &extrinsic_coded, itpp::vec &extrinsic_data, const itpp::vec &intrinsic_coded, const itpp::vec &apriori_data)
SISO decoder for RSC codes
Definition: siso.h:728
Vec< double > vec
Definition of double vector type.
Definition: vec.h:529
ivec sort_index(const Vec< T > &data, SORTING_METHOD method=INTROSORT)
Return an index vector corresponding to a sorted vector data in increasing order.
Definition: sort.h:203
Vec< int > ivec
Definition of integer vector type.
Definition: vec.h:541
Vec< bin > bvec
Definition of binary vector type.
Definition: vec.h:553
The IT++ file format reading and writing class.
Definition: itfile.h:246
vec pow(const double x, const vec &y)
Calculates x to the power of y (x^y)
Definition: log_exp.h:176
double inv_dB(double x)
Inverse of decibel of x.
Definition: log_exp.h:73
void RNG_randomize()
Set a random seed for all Random Number Generators in the current thread.
Definition: random.cpp:254
double randu(void)
Generates a random uniform (0,1) number.
Definition: random.h:804
bin randb(void)
Generates a random bit (equally likely 0s and 1s)
Definition: random.h:793
Include file for the IT++ communications module.
Mat< bin > bmat
bin matrix
Definition: mat.h:508
itpp namespace
Definition: itmex.h:37
const Array< T > concat(const Array< T > &a, const T &e)
Append element e to the end of the Array a.
Definition: array.h:486

When you run this program, the results (BER and EbN0_dB) are saved into pccc_bersim_awgn.it file. Using the following MATLAB script:

clear all
itload('pccc_bersim_awgn.it');
figure
semilogy(EbN0_dB, BER, 'o-')
grid on
xlabel('E_b/N_0 [dB]')
ylabel('BER')
double dB(double x)
Decibel of x (10*log10(x))
Definition: log_exp.h:71
ITPP_EXPORT bool all(const bvec &testvec)
Returns true if all elements are ones and false otherwise.

the results can be displayed.

Similarly, the results can be displayed using the following Python script (pyitpp, numpy and matplotlib modules are required):

#!/usr/bin/env python
from pyitpp import itload
from matplotlib.pyplot import *
out = itload('pccc_bersim_awgn.it')
semilogy(out['EbN0_dB'], out['BER'].T, 'o-')
grid()
xlabel('$E_b/N_0 [dB]$')
ylabel('$BER$')
show()

Generated on Tue Aug 17 2021 10:59:15 for IT++ by Doxygen 1.9.4