rec_signal.cpp

这是本人编写的软件接收机

//------------------------------------------------------------------------------
//
// REC_signal.cpp
//
// Purpose:
//
//   signal implementation for DSP in Software receiver
//
// Reference:
//
//
// Notes:
//
//   This software is protected by national and international copyright.
//   Any unauthorized use, reproduction or modificaton is unlawful and
//   will be prosecuted. Commercial and non-private application of the
//   software in any form is strictly prohibited unless otherwise granted
//   by the authors.
//
//   The code is provided without any warranty; without even the implied
//   warranty of merchantibility or fitness for a particular purpose.
//
// Last modified:
//
//   2007/11/04  YI Weiyong (1st edition)
//
// (c) 2007-2009  YI Weiyong
//
//------------------------------------------------------------------------------

#include <cmath>
#include <iostream>
#include <iomanip>
#include <fstream>
#include <time.h>

#include "REC_signal.h"
#include "REC_Const.h"

using namespace std;

// Constructor

CSignal::CSignal(unsigned int _n):m_n(_n){
	data = new accuracy[m_n];

    memset(data, 0, m_n*(sizeof(accuracy)/sizeof(char)) );
};

 CSignal::CSignal(accuracy *p, unsigned int _n):m_n(_n){
	data = new accuracy[m_n];
    memcpy(data, p, m_n);
};

 CSignal::CSignal(const CSignal& _sig):m_n(_sig.m_n){
	data = new accuracy[_sig.m_n];
    memcpy(data, _sig.data, m_n*(sizeof(accuracy)/sizeof(char)) );
};

 CSignal::~CSignal(){
	if(m_n != 0){
		delete [] data;
	}
	data = 0;
};

 void CSignal::Set(int *p, unsigned int _n){
	m_n = _n;
	data = new accuracy[m_n];

    memcpy(data, p, m_n*(sizeof(accuracy)/sizeof(char)) );
};


//Assignment
 CSignal&  CSignal::operator=(const accuracy _val)
{
	for(int i=0; i<m_n; i++) data[i] = _val;

 	return(*this);
}

 CSignal&  CSignal::operator=(const CSignal& _sig)
{
	if(this == &_sig) return(*this);

	if(data == 0){
		m_n = _sig.m_n;
		data = new accuracy [m_n];
	}
	// check dimension
	if(m_n != _sig.m_n){
		Resize(_sig.m_n);
	}
	// Copy elements

	memcpy(data, _sig.data, m_n*(sizeof(accuracy)/sizeof(char)) );
	// for(int i=0; i<m_n; i++) data[i] = _sig.data[i];
	return(*this);
}


// Resize:
//         change the dimension of the signal

 void  CSignal::Resize(unsigned int _n){
    if(m_n == _n) return;

	if(_n == 0 ){
		delete [] data;

		m_n = _n;
		return;
	}

	m_n = _n;
	data = new accuracy [m_n];

    memset(data, 0, m_n*(sizeof(accuracy)/sizeof(char)) );
	// for(int i=0; i<m_n; i++) data[i] = 0;
};

//------------------------------------------------------
// ReSample
//        purpose: Resample the signal with a specific factor
// Return:
//        resampled signal

 CSignal  CSignal::Resample(double factor)
{
	if(factor <= 0){
		cerr<<"Error: resampling factor is less or equal to zero, exit"<<endl;
		exit(0);
	}

	int n = (int)(m_n*factor+0.5);
	CSignal Aux(n);
	int i, j;
	int mm = (int)(1/factor+0.5);
	int sum;

    if(factor<1)
    {
		for( i=0; i<n; i++)
		{
			sum = data[ (int)(i/factor  +0.5) ];
			//for( j=0; j<mm; j++){
			//	sum += data[ (int)(i/factor + j +0.5) ];
			//}
			Aux.data[i] = sum;
		}
	}
	else
	{
		for( i=0; i<n; i++)
		{
			Aux.data[i] = data[ (int)(i/factor+0.5) ];
		}
    }
	return (Aux);
}

 void CSignal::Resample(CSignal& ReSampledSig, double factor)
{
	if(factor <= 0){
		cerr<<"Error: resampling factor is less or equal to zero, exit"<<endl;
		exit(0);
	}

	int n = (int)(m_n*factor+0.5);
	CSignal Aux(n);
	int i, j;
	int mm = (int)(1/(factor+0.5));
	int sum;

    if(factor<1)
    {
		for( i=0; i<n; i++)
		{
			sum = 0;
			for( j=0; j<mm; j++){
				sum += data[ (int)(i/factor + j +0.5) ];
			}
			ReSampledSig.data[i] = sum;
		}
	}
	else
	{
		for( i=0; i<n; i++)
		{
			ReSampledSig.data[i] = data[ (int)(i/factor+0.5) ];
		}
    }
}

 void CSignal::ZeroPadding(const unsigned int& numberOfZeros)
{
	unsigned int n ;
	unsigned int i = 0;


    if(numberOfZeros == 0)
    {
		if(!(m_n&(m_n-1))) return;

		while( (1<<i)<m_n ) i++;
	    n = 1<<i;
	}
	else
	{
		n = m_n+numberOfZeros;
	}

	int *temp;
	int m;
	temp = new int [m_n];
	//memcpy(temp, data, m_n*(sizeof(int)/sizeof(char)) );
	for( i = 0; i<m_n; i++) temp[i] = data[i];
	m = m_n;
	Resize(n);
	for( i = 0; i<m; i++) data[i] = temp[i];
	//memcpy(data, temp, m_n*(sizeof(int)/sizeof(char)) );

	delete [] temp;
}

 CSignal  CSignal::slice(int first, int last) const
{
	CSignal Aux(last-first+1);
	for(int i=first; i<=last; i++) Aux.data[i-first] = data[i];

	return Aux;
}

 void  CSignal::shift(int n)
{
	n = n % (int)m_n;

	if(n<0) n = m_n + n;

	int i, tick;
	int *temp;

    tick = m_n - n;

	temp = new int [m_n];
	for( i=0; i<tick; i++) temp[i] = data[i+n];
	for( i=tick; i<m_n; i++) temp[i] = data[i-tick];

	for( i=0; i<m_n; i++) data[i] = temp[i];

	delete [] temp;
}


 CSignal  CSignal::PSD()
{
	CSignal PSD;
	CSignal Aux(*this);

    Aux.ZeroPadding();
    CSignal Aux_im(Aux.m_n);

    CFixedFFT myFFT(Aux.m_n);
    myFFT.FFT(Aux, Aux_im);

    Aux = Aux.slice(0, (Aux.m_n>>1) );
    Aux_im = Aux_im.slice(0, (Aux_im.m_n>>1) );

    PSD = Magnitude(Aux, Aux_im);

	return PSD;
}

 accuracy  CSignal::Max(int& _Index)
{
	_Index = 0;
	accuracy RetVal;
	RetVal = data[0];
	for(int i=1; i<m_n; i++)
	{
		if(data[i]>data[i-1]){
			RetVal = data[i];
			_Index = i;
		}
	}
	return RetVal;
}

 int  CSignal::MaxIndex()
{
	int Index = 0;
	int RetVal;
	RetVal = data[0];
	for(int i=1; i<m_n; i++)
	{
		if(data[i]>data[i-1]){
			RetVal = data[i];
			Index = i;
		}
	}
	return Index;
}

 accuracy  CSignal::Max()
{
	accuracy RetVal;
	RetVal = data[0];
	for(int i=1; i<m_n; i++)
	{
		if(data[i]>data[i-1]){
			RetVal = data[i];
		}
	}
	return RetVal;
}

 void CSignal::PassFIR(const CSignal& _coef_b)
{
	int i, j, half_m_n(_coef_b.m_n/2);
	int coefsum(0), sum;
	CSignal Aux(*this);

	if(m_n <= _coef_b.m_n){
		cerr<<"Warning: the size of the signal is less than filter in PassFIR, do nothing and return"<<endl;
		return;
	}

	for( j=0; j<_coef_b.m_n; j++ )
	{
		coefsum += _coef_b.data[j];
	}

	for( i=0; i<half_m_n; i++ )
	{
		sum = 0;
		for( j=0; j<_coef_b.m_n; j++ )
		{
			sum += Aux.data[i+j] * _coef_b.data[j];
		}
		data[i] = sum/coefsum;
	}

	for( i=half_m_n; i<m_n - half_m_n; i++ )
	{
		sum = 0;
		for( j=0; j<_coef_b.m_n; j++ )
		{
			sum += Aux.data[i+j-half_m_n] * _coef_b.data[j];
		}
		data[i] = sum/coefsum;
	}
	for( i=m_n - half_m_n; i<m_n; i++ )
	{
		sum = 0;
		for( j=0; j<_coef_b.m_n; j++ )
		{
			sum += Aux.data[i+j - _coef_b.m_n] * _coef_b.data[j];
		}
		data[i] = sum/coefsum;
	}
}

 void CSignal::ZeroDC()
{
	int sum(0);
	int mean;
	int i;

	for( i=0; i<m_n; i++) sum += data[i];

	mean = (int)((double)sum/m_n+0.5);

	for( i=0; i<m_n; i++) data[i] -= mean;
}

 accuracy CSignal::Mean()
{
	int mean(0);
	int sum(0);

	for( int i=0; i<m_n; i++)
	{
		sum += data[i];
	}
	mean = sum/m_n;

	return mean;
}

 ostream& operator << (ostream& os, const CSignal& Sig)
{
  int w = os.width();
  for (int i=0; i<Sig.size(); i++)
    os << setw(w) << Sig(i);
//os << endl;
  return os;
}

 CSignal ZeroDC(const CSignal& _sig)
{
	CSignal Aux(_sig.m_n);
	int sum(0);
	int mean;
	int i;
	int m_n(_sig.m_n);

	for( i=0; i<m_n; i++) sum += _sig.data[i];

	mean = (int)((double)sum/m_n +0.5);

	for( i=0; i<m_n; i++) Aux.data[i] =_sig.data[i] - mean;

	return Aux;
}

//-------------------------------------------------------------------
// Correlate
//          Purpose: Correlate two signals with delay

 long Correlate(const CSignal& _sig1, const CSignal& _sig2, int delay)
{
	if(_sig1.m_n != _sig2.m_n){
		cerr<<" Error: imcompatible size in CSignal::Correlate(), exit!! "<<endl;
	    exit(0);
	}

	delay = delay % (int)_sig1.m_n;

    if(delay < 0) delay = _sig1.m_n + delay;

	long sum = 0;
	unsigned int tick = _sig1.m_n-delay;
	int i;

	for(i = 0; i<tick; i++){
		sum += _sig1.data[i+delay] * _sig2.data[i];
	}

	for(i = tick; i<_sig1.m_n; i++){
		sum += _sig1.data[i-tick] * _sig2.data[i];
	}

	return sum;
}

 CSignal Stack(const CSignal& _sig1, const CSignal& _sig2)
{
	int i;
	CSignal Aux(_sig1.m_n+_sig2.m_n);
	for(i = 0; i<_sig1.m_n; i++) Aux.data[i] = _sig1.data[i];
	for(i = 0; i<_sig2.m_n; i++) Aux.data[i+_sig1.m_n] = _sig2.data[i];

	return Aux;
}

// product of two signals

 CSignal operator * (const CSignal& _sig1, const CSignal& _sig2)
{
	if(_sig1.m_n != _sig2.m_n){
		cerr<<" Error: imcompatible size in CSignal:: CSignal * CSignal , exit!! "<<endl;
		exit(1);
	}

	CSignal Aux(_sig1.m_n);

	for(int i=0; i<_sig1.m_n; i++) Aux.data[i] = _sig1.data[i]*_sig2.data[i];

	return Aux;
}

CSignal operator * (const CSignal& _sig1, const accuracy& _val)
{
	CSignal Aux(_sig1.m_n);

	for(int i=0; i<_sig1.m_n; i++) Aux.data[i] = _sig1.data[i]*_val;

	return Aux;
}


// sum of two signals

 CSignal operator + (const CSignal& _sig1, const CSignal& _sig2)
{
	if(_sig1.m_n != _sig2.m_n){
		cerr<<" Error: imcompatible size in CSignal:: CSignal + CSignal , exit!! "<<endl;
		exit(1);
	}

	CSignal Aux(_sig1.m_n);

	for(int i=0; i<_sig1.m_n; i++) Aux.data[i] = _sig1.data[i]+_sig2.data[i];

	return Aux;
}

// sum of two signals

 CSignal operator - (const CSignal& _sig1, const CSignal& _sig2)
{
	if(_sig1.m_n != _sig2.m_n){
		cerr<<" Error: imcompatible size in CSignal:: CSignal - CSignal , exit!! "<<endl;
		exit(1);
	}

	CSignal Aux(_sig1.m_n);

	for(int i=0; i<_sig1.m_n; i++) Aux.data[i] = _sig1.data[i]-_sig2.data[i];

	return Aux;
}

// sum of two signals

 CSignal operator - (const CSignal& _sig1, const accuracy& _val)
{

	CSignal Aux(_sig1.m_n);

	for(int i=0; i<_sig1.m_n; i++) Aux.data[i] = _sig1.data[i] - _val;

	return Aux;
}


// Negation of a CSignal (unary minus)

 CSignal operator - (const CSignal& _sig)
{
	CSignal Aux(_sig.m_n);
	for( int i=0; i<_sig.m_n; i++) Aux.data[i] = - _sig.data[i];
	return Aux;
}

// Return the magnitude of the two parts(real and image) of one signal
// Just return Abs(real)+Abs(image), not very strict, but just an approximation

 CSignal Magnitude(const CSignal& _sig1, const CSignal& _sig2)
{
	if(_sig1.m_n != _sig2.m_n )
	{
		cerr<<"Error: incompatible size in Magnitude(const CSignal& _sig1, const CSignal& _sig2), exit!"<<endl;
		exit(1);
	}

	CSignal Aux(_sig1.m_n);

	for(int i=0; i<Aux.m_n; i++)
	{
		Aux.data[i] = rss(_sig1.data[i], _sig2.data[i]);
	}
	return Aux;
}

// Generate a Gaussian distributed signal

 CSignal Randn(unsigned int n, double variance, double scale)
{
     double random1, random2;
     CSignal Aux(n);
     for(int i = 0; i<n; i++)
     {
          random1 = (rand() + 1e-15) * 1.0/(RAND_MAX+1e-15);
          random2 = (rand() + 1e-15) * 1.0/(RAND_MAX+1e-15);
          Aux.data[i] =  (int)(sqrt( - 2*log(random1) ) * cos(random2 * pi2) * variance * scale);
	  }
     return Aux;
 }


//------------------------------------------------------
// Sine
//        purpose: generate a sinusoidal signal with samples length n
// In:
//        freq ------- freqency of the signal
//        fs   ------- sampling frequency, must be large than two times freq according to Nyquist law
//        start ------- start phase of the signal(Normalized phase, i.e. 0 to 1
//        n     ------- number of samples
//        scale ------- scale factor for the signal
// Out:
//        None
// Return:
//        a sinusoidal signal with length n

 CSignal Sine(double freq, double fs, double start, int n, int scale)
{
	double omega = 2*pi*freq/fs*MAGIC;
	double phase = start*2*pi*freq*MAGIC;
	CSignal Aux(n);

	for(int i = 0; i<n; i++)
	{
		Aux.data[i] = static_cast<int>(SinTab[( (int)(omega*i+phase+0.5) ) & (MAX_TRIG)]>>8);
	}
	return Aux;
}

//-------------------------------------------------------------------
// ShiftProd
//          Purpose: Compute the product of two signals
//                   the second signal is shifted with value delay

 CSignal ShiftProd(const CSignal& _sig1, const CSignal& _sig2, int delay)
{
	if(_sig1.m_n != _sig2.m_n){
		cerr<<" Error: imcompatible size in CSignal::Correlate(), exit!! "<<endl;
	    exit(0);
	}
    CSignal Aux(_sig1.m_n);

	delay = delay % (int)_sig1.m_n;

    if(delay < 0) delay = _sig1.m_n + delay;

	unsigned int tick = _sig1.m_n-delay;
	int i;

	for(i = 0; i<tick; i++){
		Aux.data[i]= _sig1.data[i+delay] * _sig2.data[i];
	}

	for(i = tick; i<_sig1.m_n; i++){
		Aux.data[i]= _sig1.data[i-tick] * _sig2.data[i];
	}

	return Aux;
}


// Generate PRN code regarding PRNNo

 CSignal CACode(unsigned int PRNNo)
{
	int tap1, tap2;
	int reg1[10], reg2[10], check[10];