system.cpp

Turbo码性能仿真源程序信道编码方面的

//**************************************************************//
//                                                              //
//                     August 2002                              //
//                         by                                   //
//                Leung, Raymond Wai-Kong                       // 
//                                                              //
//                                                              //
//**************************************************************//

#include <stdio.h>
#include <stdlib.h>
#include <fstream.h>
#include <math.h>
#include "Global.h"


extern Random	Rnd;

//Constant of Modular Arithmetic Generator
const int	Random::A	=	48271;
const int	Random::M	=	2147483647;
const int	Random::Q	=	44488;
const int	Random::R	=	3399;
//Constant of Wichmann/Hill Generator
const int	Random::A1	=	171;
const int	Random::A2	=	172;
const int	Random::A3	=	170;
const int	Random::M1	=	30269;
const int	Random::M2	=	30307;
const int	Random::M3	=	30323;
//Box-Muller Method
const double	Random::pi2 = 6.2831853071795864;




//------------------------------Interleaver------------------------------------
//function : permutation
//input    : *rule, length, scrnum
void Scramble( int *rule, int length, int scrnum )
{
	int i, j, n, r, q;
	int **buffer, *rbuffer;

    q = (int) sqrt( length );
	if( q*q < length ) q++;

	rbuffer = new int [q];
	buffer = new int* [q];
	for( i=0; i<q; i++ ) buffer[i] = new int [q];

	for( n=0, i=0; i<q; i++ ) for( j=0; j<q; j++ )
	{
		if( n < length ) buffer[i][j] = rule[n++];
		else buffer[i][j] = -1;
	}

	for( n=0; n<scrnum; n++ )
	{
		for( i=0; i<q; i++ )
		{			
			r = Rnd.RandomInt( ) % q;//r = rand( ) % q;//random variables
			for( j=0; j<q; j++ ) rbuffer[j] = buffer[i][j];
			for( j=0; j<q; j++ ) buffer[i][j] = rbuffer[(j+r)%q];
		}
		for( j=0; j<q; j++ )
		{			
			r = Rnd.RandomInt( ) % q;//r = rand( ) % q;//random variables
			for( i=0; i<q; i++ ) rbuffer[i] = buffer[i][j];
			for( i=0; i<q; i++ ) buffer[i][j] = rbuffer[(i+r)%q];
		}
	}

	for( n=0, i=0; i<q; i++ ) for( j=0; j<q; j++ )
		if( buffer[i][j] > -1 ) rule[n++] = buffer[i][j];

	for( i=0; i<q; i++ )
		delete [] buffer[i];

	delete [] buffer;
	delete [] rbuffer;
}

//function : used in the EF-interleaver
//input    : *x, length, shift
void CyclicShift( int *x, int length, int shift )
{
	int j, m;
	int *bf = new int [length];
	
	for( m=shift, j=0; j<length-shift; j++ ) bf[j] = x[m++];
	for( m=0, j=length-shift;j<length;j++) bf[j] = x[m++];
	for(j=0;j<length;j++) x[j]=bf[j];
	delete [] bf;
}



//----------------------------Input/Output function----------------------------------
//function : format an int data and place it in the string named string.
//			use .FnInit( ) to initialize the function  FnInit -- File name initialization
void IO::FnInit( )
{
//	FileN = new char[256];
	ptr = 0;
}
//input    : num -- Attach an int number to the file name
void IO::FnJoin( int num )
{
	ptr += sprintf( FileN + ptr, "%d", num );
}
//input    : num -- Attach a float number with only 2 decimal to the file name
void IO::FnJoin( float num )
{
	ptr += sprintf( FileN + ptr, "%1.3f", num );
}
//input    : *prefix -- Attach a string to the file name
void IO::FnJoin( char *prefix )
{
	ptr += sprintf( FileN + ptr, "%s", prefix );
}
//function : output the simulation parameters into the file
void IO::FileHeader( )
{
	fstream file;
	int i;

	file.open( FileN, ios::out );
	if( !file ) 
	{
		cout << "Error -- Cannot write the result into the file : Header" << endl;
		exit (-1);
	}
	file << _SimulationMode << endl;									cout << _SimulationMode << endl;
	file << _ChannelType << endl;										cout << _ChannelType << endl;
	file << "Nt					= " << (int)NT << endl;					cout << "Nt					= " << (int)NT << endl;
	file << "Data length		= " << (int)DATALEN << endl;			cout << "Data length		= " << (int)DATALEN << endl;
	file << "E-Fold Number		= " << (int)EFOLDNUM << endl;			cout << "E-Fold Number		= " << (int)EFOLDNUM << endl;
	file << "F-Fold Number		= " << (int)FFOLDNUM << endl;			cout << "F-Fold Number		= " << (int)FFOLDNUM << endl;
	file << "Fold Number		= " << (int)FOLDNUM << endl;			cout << "Fold Number		= " << (int)FOLDNUM << endl;
	file << "NestDim			= " << (int)NESTDIM << endl;			cout << "NestDim			= " << (int)NESTDIM << endl;
	file << "Punctured Bit		= " << (int)PUNCTURE << endl;			cout << "Punctured Bit		= " << (int)PUNCTURE << endl;
	file << "State Order		= " << (int)POLYSTATEORDER << endl;		cout << "State Order		= " << (int)POLYSTATEORDER << endl;
	file << "N-Polynomial		= " << (int)NPOLY << endl;				cout << "N-Polynomial		= " << (int)NPOLY << endl;
	file << "D-Polynomial		= " << (int)DPOLY << endl;				cout << "D-Polynomial		= " << (int)DPOLY << endl;
	file << "Spreading Length	= " << (int)SPREADLEN << endl;			cout << "Spreading Length	= " << (int)SPREADLEN << endl;
	file << "Stacking Length	= " << (int)STACKNUM << endl;			cout << "Stacking Length	= " << (int)STACKNUM << endl;
	file << "Frame Length		= " << _FrameLen << endl;				cout << "Frame Length		= " << _FrameLen << endl;
	file << "Chip Length		= " << _ChipNum << endl;				cout << "Chip Length		= " << _ChipNum << endl;
	file << "TranChip Length	= " << _TranChipNum << endl;		cout << "TranChip Length		= " << _TranChipNum << endl;
	file << "MinBlockNum		= " << _MinBlockNum << endl;			cout << "MinBlockNum		= " << _MinBlockNum << endl;
	file << "MaxBlockNum		= " << _MaxBlockNum << endl;			cout << "MaxBlockNum		= " << _MaxBlockNum << endl;

	file << "\t";
	for( i=1; i<=_ItNum; i++ )
		file << "It=" << i << "\t";
	file << "\t";
	for( i=1; i<=_ItNum; i++ )
		file << "It=" << i << "\t";
	file << endl;
	file.close( );
}

//function : output the intermediate result in the screen
//input    : Snr, BlockNum, *Error
void IO::ScreenOut( double Snr, int BlockNum, int *BitError, int *FrameError )
{
	cout << "e=" << BitError[_ItNum-1];
	cout << ", snr=" << _SnrStart << "dB, ";
	cout.setf(ios::scientific);

	if( BitError[_ItNum-1] > 0 )
		cout << "  BER=" << (float)BitError[_ItNum-1] / ( (float)BlockNum * (float)_DataLen );
	else
		cout << "  BER<" << 1 / ( (float)BlockNum * (float)_DataLen );

	if( FrameError[_ItNum-1] > 0 )
		cout << "  FER=" << (float)FrameError[_ItNum-1] / (float)BlockNum;
	else
		cout << "  FER<" << 1 / (float) BlockNum;

	cout.unsetf(ios::scientific);
	cout << "  #k=" << BlockNum << endl;
}
//function : output the result to the file of the simulation point
void IO::FileOut( float Snr, int BlockNum, int *BitError, int *FrameError )
{
	fstream file;
	int i;

	file.open( FileN, ios::app );
	if( !file )
	{
		cout << "Cannot write the result into the file : Content" <<endl;
		exit (-1);
	}

	file << Snr << "\t";
	file.setf(ios::scientific);
	for( i=0; i<_ItNum; i++ )
		file << ( BitError[i] / (float)BlockNum / (float)_DataLen ) << "\t";
	file << "\t";
	for( i=0; i<_ItNum; i++ )
		file << ( FrameError[i] / (float)BlockNum ) << "\t";
	file.unsetf(ios::scientific);
	file << "#k=" << BlockNum << endl;
	file.close( );
}

void IO::ErrorMessage( char *message )
{
	cout << "ERROR -- " << message << endl;
	exit(-1);
}

//-------------------------------Random number generator--------------------------------------
//Modular Arithmetic Generator
void Random::SetState( int S )
{
	stateS = S;
}

void Random::ReadState( int* S )
{
	*S = stateS;
}

void Random::IntNextSate( int *S )
{
	int tmpState = A * ( *S % Q ) - R * ( *S / Q );

	if( tmpState >= 0 )
		*S = tmpState;
	else
		*S = tmpState + M;
}

int Random::RandomInt( )
{
	Random::IntNextSate( &stateS );
	return stateS;
}
//Wichmann/Hill Generator
void Random::SetState( int X, int Y, int Z )
{
	stateX = X;
	stateY = Y;
	stateZ = Z;
}
void Random::ReadState( int* X, int* Y, int* Z )
{
	*X = stateX;
	*Y = stateY;
	*Z = stateZ;
}
void Random::DoubleNextState( int* X, int* Y, int* Z )
{
	*X = ( *X * A1 ) % M1;
	*Y = ( *Y * A2 ) % M2;
	*Z = ( *X * A3 ) % M3;
}
double Random::RandomDouble( )
{
	double tmp;

	Random::DoubleNextState( &stateX, &stateY, &stateZ );
	tmp = (double)stateX / (double)M1 + (double)stateY / (double)M2 + (double)stateZ / (double)M3;
	return tmp - (double) ( (int)tmp );//mod 1 operatoion
}
//Both Wichmann/Hill and Modular Arithmetic Generator
void Random::SetState( int S, int X, int Y, int Z )
{
	stateS = S;
	stateX = X;
	stateY = Y;
	stateZ = Z;
}
void Random::ReadState( int* S, int* X, int* Y, int* Z )
{
	*S = stateS;
	*X = stateX;
	*Y = stateY;
	*Z = stateZ;
}
//Box-Muller Method
float Random::RealNoise( double sigma )
{
	double Amp = sigma * sqrt( -2.0 * log( Random::RandomDouble( ) ) );
	double Angle = pi2 * Random::RandomDouble( );

	return (float) ( Amp * cos( Angle ) );
}

Complex Random::ComplexNoise( double sigma )
{
//	struct Complex symbol;
	Complex symbol;
	double Amp = sigma * sqrt( -2.0 * log( Random::RandomDouble( ) ) );
	double Angle = pi2 * Random::RandomDouble( );

	symbol.Re = (float) ( Amp * cos( Angle ) );
	symbol.Im = (float) ( Amp * sin( Angle ) );

	return symbol;
}

Complex Random::RayleighGenerator( )
{
	Complex symbol;
	double Amp = sqrt( -log( Random::RandomDouble( ) ) );
	double Angle = pi2 * Random::RandomDouble( );

	symbol.Re = (float) ( Amp * cos( Angle ) );
	symbol.Im = (float) ( Amp * sin( Angle ) );

	return symbol;
}





Complex operator + ( const Complex& a, const Complex& b )
{
	Complex result;

	result.Re = a.Re + b.Re;
	result.Im = a.Im + b.Im;
	return result;
}

Complex operator - ( const Complex& a, const Complex& b )
{
	Complex result;

	result.Re = a.Re - b.Re;
	result.Im = a.Im - b.Im;
	return result;
}

Complex operator * ( const Complex& a, const Complex& b )
{
	Complex result;

	result.Re = a.Re * b.Re - a.Im * b.Im;
	result.Im = a.Re * b.Im + a.Im * b.Re;
	return result;
}





float LClip1( float x )
{
	if( x > (float) LTOOLARGE1 ) return (float) LTOOLARGE1;
	else if( x < (float) LTOOSMALL1 ) return (float) LTOOSMALL1;
	else return x;
}

float LClip2( float x )
{
	if( x > (float) LTOOLARGE2 ) return (float) LTOOLARGE2;
	else if( x < (float) LTOOSMALL2 ) return (float) LTOOSMALL2;
	else return x;
}

float LClip3( float x )
{
	if( x > (float) LTOOLARGE3 ) return (float) LTOOLARGE3;
	else if( x < (float) LTOOSMALL3 ) return (float) LTOOSMALL3;
	else return x;
}

float RClip1( double x )
{
	if( x > (float) RTOOLARGE1 ) return (float) RTOOLARGE1;
	else if( x < (float) RTOOSMALL1 ) return (float) RTOOSMALL1;
	else return (float) x;
}

float RClip2( double x )
{
	if( x > (float) RTOOLARGE2 ) return (float) RTOOLARGE2;
	else if( x < (float) RTOOSMALL2 ) return (float) RTOOSMALL2;
	else return (float) x;
}

float RClip3( double x )
{
	if( x > (float) RTOOLARGE3 ) return (float) RTOOLARGE3;
	else if( x < (float) RTOOSMALL3 ) return (float) RTOOSMALL3;
	else return (float) x;
}

float RClipT( double x )
{
	if( x < (float) TINY ) return (float) TINY;
	else return (float) x;
}