test.cpp

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

#include <iostream>
#include <iomanip>
#include <cstdlib>
#include <time.h>
#include <stdlib.h>
#include <stdio.h>
#include <fstream>
#include <vector>
#include <bitset>
//#include <pthread.h>
#include <windows.h>

#include "REC_signal.h"
#include "REC_const.h"

using namespace std;

template<class T>
T abs(T _a)
{
	return (_a>=0)? _a: -_a;
}

struct settings{

	int Flag;                       // 0 --- simulation,  1 --- process real data
	int NumberOfChannels;
	int PRNNoForChannel[16];
    char FileName[1024];
    double IF;                          // Intermediate Frequency
    double Fs;                          // sampling frequency
    double CodeFreq;                    // code frequency
    int CodeLength;
    double acqSearchBand;               // Search band in frequency domain[kHz]
    int    acqSamples;                  // Samples over every milisecond for acquisition
    double acqThreshold;
    double dllDampingRatio;             // damping ratio for Code tracking loop
    double dllNoiseBandwidth;           // [Hz]
    double dllCorrelationSpacing;       // [Chips]

    double pllDampingRatio;             // damping ratio for phase track loop
    double pllNoiseBandwidth;           // Bandwidth of the noise in phase track loop

    double navSolPeriod;                // Period for calculating psuedo-ranges and position

    double elevationMask;               // Elevation mask to exclude signals from sat. at low elevation

    int   useTropCorr;                  // Enable/Disable tropospheric correction: 0 --- off, 1--- on

    double x, y, z;                     // initial coordinates of the receiver in WGS84

    // For simulation
    int NoiseVar;
};

CSignal CATable[32];
CFixedFFT recFFT(65536);
CSignal *SinCarr ;
CSignal *CosCarr ;
CSignal ** CodeFreqI;
CSignal ** CodeFreqQ;
CSignal *SigAcqNonCoherent;
void * pnt;
bool AcquiringFlag = false;   // == 1  means the acquisition is runing, ==0 means not
ofstream globalout;
int gAcquiring;                 // == 1 means there is a channel perform acquisition, 0 means otherwise


void GetSetup(ifstream& inp, settings& mySetting)
{
	cout<<" Reading setting file, i.e. rec.inp"<<endl;

	// Scenario, 0 -- simulation, 1 means real data
	 inp.ignore(25); inp >>mySetting.Flag;  inp.ignore(81,'\n');

	 // Number of Channels
	 inp.ignore(25); inp >>mySetting.NumberOfChannels;  inp.ignore(81,'\n');
	 inp.ignore(25);
	 for(int i=0; i<mySetting.NumberOfChannels; i++)
	 {
		 inp >>mySetting.PRNNoForChannel[i];
	 }
	 inp.ignore(81,'\n');

	 // File name of real data
	 inp.ignore(25); inp >>mySetting.FileName;       inp.ignore(81,'\n');
	 inp.ignore(25); inp >>mySetting.IF;             inp.ignore(81,'\n');
	 inp.ignore(25); inp >>mySetting.Fs;             inp.ignore(81,'\n');
	 inp.ignore(25); inp >>mySetting.CodeFreq;       inp.ignore(81,'\n');
	 inp.ignore(25); inp >>mySetting.CodeLength;     inp.ignore(81,'\n');
	 inp.ignore(25); inp >>mySetting.acqSearchBand;  inp.ignore(81,'\n');
	 inp.ignore(25); inp >>mySetting.acqSamples;     inp.ignore(81,'\n');
	 inp.ignore(25); inp >>mySetting.acqThreshold;   inp.ignore(81,'\n');
	 inp.ignore(25); inp >>mySetting.dllDampingRatio;  inp.ignore(81,'\n');

	 inp.ignore(25); inp >>mySetting.dllNoiseBandwidth;  inp.ignore(81,'\n');
	 inp.ignore(25); inp >>mySetting.dllCorrelationSpacing;  inp.ignore(81,'\n');

	 inp.ignore(25); inp >>mySetting.pllDampingRatio;  inp.ignore(81,'\n');
	 inp.ignore(25); inp >>mySetting.pllNoiseBandwidth;  inp.ignore(81,'\n');

	 inp.ignore(25); inp >>mySetting.navSolPeriod;  inp.ignore(81,'\n');
	 inp.ignore(25); inp >>mySetting.elevationMask;  inp.ignore(81,'\n');

	 inp.ignore(25); inp >>mySetting.useTropCorr;  inp.ignore(81,'\n');

	 inp.ignore(25); inp >>mySetting.x>>mySetting.y>>mySetting.z;  inp.ignore(81,'\n');

	 inp.ignore(256,'\n');
	 inp.ignore(25); inp >>mySetting.NoiseVar;  inp.ignore(81,'\n');

}


//--------------------------------------------------------------------------------
// Channel class implementation


class Channel{
	public:

		enum ChannelState{
			INACTIVE = 1,
			ACTIVE   = 2,
			ACQUIRING= 3,
			FINISHACQUIRING = 4,
			ACQUIRED = 5,
			TRACKING = 6
		};

	      // Constructor and destructor
          Channel()
	      {
			  m_AcqCount = 0;
			  m_PRNNo = 1;
			  m_CACode = CATable[m_PRNNo-1];
			  m_Bins = 0;
			  m_State = ACTIVE;
			  m_SamplesPerChip = 1;
			  m_FreqIndxForAcquiring = 0;

	      };
	      Channel(const settings& _mySet):m_Setting(_mySet)
	      {
			  m_AcqCount = 0;m_PRNNo = 1; m_CACode = CATable[m_PRNNo-1];
	          m_Bins = (int)(_mySet.acqSearchBand*2+1);
	          m_isAcquired = 0;
	          m_State = ACTIVE;
	          m_SamplesPerChip = static_cast<int>(_mySet.Fs/_mySet.CodeFreq);
	          m_FreqIndxForAcquiring = 0;

	      };

	      ~Channel()
	      {

		  };

	      // Methods
          void SetChannelPara(const settings& _mySet){ m_Setting = _mySet;};
          void SetState(const ChannelState& _state){ m_State = _state; };
          ChannelState GetState(){ return m_State; };

	      void Acquisition();
	      void Tracking();
	      void Run(const CSignal& _sig);
	      void Activate(){ m_isActive = 1; m_AcqCount = 0; m_State = ACTIVE; };
	      void InActivate() { m_isActive = 0; m_isAcquired = 0;};

	      void SetPRNNo(unsigned int _PRN){ m_PRNNo = _PRN; m_CACode = CATable[_PRN-1];m_isActive = 1; m_AcqCount = 0;};
          unsigned int GetPRNNo(){ return m_PRNNo;};
          double GetCodePhase(){ return m_CodePhase; };
          ChannelState m_State;
          const static int MAXSEQ = 30;
          double m_AcqCodeSequence[MAXSEQ];
          double m_AcqFreqSequence[MAXSEQ];

	private:
	      CSignal m_sig;
	      CSignal m_sigForAcq;

	      CSignal m_CACode;
	      settings m_Setting;
	      unsigned int m_PRNNo;             // PRN No. of the satellite tracked by the channel
	      int    m_isActive;                // to tell whether it is active, 0 -- inactive, 1 -- active
          int    m_isAcquired;              // whether the signal is acquired, 0 -- no, 1 -- yes
          int    m_FreqIndxForAcquiring;

          int    m_CodeShift;
          double m_DopplerShift;
          int m_Bins;                             // frequency bins for doppler searching
          int m_SamplesPerChip;
          double m_CodePhase;

          // I and Q value of previous epoch

          int    m_PrevI;
          int    m_PrevQ;

          // I and Q value of punctual, late and early
          int m_I, m_Q, m_IL, m_QL, m_IE, m_QE;

          int m_AcqCount;                   // times of acquisition, if this number is too large,
                                            // set to zero and inactivate this channel

};

void Channel::Acquisition()
{
    gAcquiring = 1;
	ofstream out1("inout\\out3.dat");
	register int i, j, n;
	const double m = (double)m_Setting.acqSamples;
	int mm = (int)m;
	const int sigSize = m_sig.size();
	const double resample = m/sigSize;

	int maximum;
	int maxCodePhase;
	int maxFreq, secondMax, secondCodePhase, secondFreq;
 	int mean;
 	int acquired;
 	CSignal sig7[m_Bins];

	clock_t start, endt;

    CSignal sig1(m_sig.Resample(m/sigSize) );
    CSignal sig2(sig1.size() );
    const CSignal refSig= m_CACode.Resample(m/m_CACode.size());

    CSignal sig3, sig4;
    CSignal sig5, sig6;

    n = sig1.size();

    //recFFT.FFT(sig1, sig2);

     start = clock();
     fft(sig1, sig2);


    //i = m_FreqIndxForAcquiring;
    for( i=0; i<m_Bins; i++)
    {
        sig5 = sig1*CodeFreqI[m_PRNNo-1][i] + sig2*CodeFreqQ[m_PRNNo-1][i];
        sig6 = sig1*CodeFreqQ[m_PRNNo-1][i] - sig2*CodeFreqI[m_PRNNo-1][i];

        //recFFT.FFT(sig5, sig6, 1);
        fft(sig5, sig6, -1);

        SigAcqNonCoherent[i] = SigAcqNonCoherent[i]+Magnitude(sig5, sig6);

       // m_sigAcqSave[i] = sig7[i];
       // out1<<setw(15)<<SigAcqNonCoherent[i]<<endl;
       // sig7[i] = SigAcqNonCoherent[i];
	}

        //
   // m_FreqIndxForAcquiring ++;

    // find the maximum of the correlation

   // if(m_FreqIndxForAcquiring >= m_Bins){
	//	m_FreqIndxForAcquiring = 0;
		maximum = SigAcqNonCoherent[0](0);
		maxCodePhase = 0;
		maxFreq  = 0;
		mean = 0;

		secondMax = SigAcqNonCoherent[0](0);
		secondCodePhase = 0;
		secondFreq = 0;

		int const acqSamplesPerChip = static_cast<int>((double)m_Setting.acqSamples/m_Setting.CodeLength);

		for( i=0; i<m_Bins; i++)
		{
			mean += SigAcqNonCoherent[i].Mean();
			for( j=0; j<n; j++)
			{
				if(maximum < SigAcqNonCoherent[i](j) )
				{
					maximum = SigAcqNonCoherent[i](j);
					maxCodePhase = j;
					maxFreq = i;
				}
			}
		}

		for( i=0; i<m_Bins; i++)
		{
			for( j=0; j<n; j++)
			{
				if( ( secondMax < SigAcqNonCoherent[i](j) ) && abs(j-maxCodePhase) > acqSamplesPerChip )
				{
					secondCodePhase = j;
					secondFreq = i;
					secondMax = SigAcqNonCoherent[i](j);
				}
			}
			//m_sigAcqSave[i] = 0;
		}

		mean = mean/m_Bins;
		acquired = 0;
		if( maximum/secondMax > 4 ) acquired = 1;

		m_CodePhase = sigSize - maxCodePhase/resample;
		m_CodeShift = static_cast<int>(m_CodePhase*m_Setting.CodeLength/sigSize+0.5);
        endt = clock();

	//	cout<<m_AcqCount<<"Mean = "<<mean<<" ,  Max = "<<maximum<<", SecondMax = "<<secondMax<<" , maxCodePhase = "<<m_CodePhase <<" "<<secondCodePhase<<" , Max freqency = "<<m_Setting.IF - m_Setting.acqSearchBand*1e3/2.0 + 0.5e3*maxFreq<<endl;

        m_AcqCount ++;
		cout<<"Time duration Acquisition  "<<m_PRNNo<<" "<<(double)(endt-start)/CLOCKS_PER_SEC<<"  [s]"<<endl;
  //  }

    if(m_AcqCount >= MAXSEQ ){

		int *CodeHist, *FreqHist;
		const double spacing = 400.0;
		const int CodeSpan = static_cast<int>(sigSize/spacing + 0.5);    // dimension of histogram
		int indx, indx1;
		bool condition1=false, condition2 = false;

		CodeHist = new int[CodeSpan];
		FreqHist = new int[m_Bins];

		for( i=0; i<CodeSpan; i++)
		{
			CodeHist[i] = 0;
		}
		for( i=0; i<m_Bins; i++)
		{
 			FreqHist[i] = 0;
		}
		for( i=0; i<m_Bins; i++)
		{
 			SigAcqNonCoherent[i] = 0;
		}


		for( i=0; i<MAXSEQ; i++)
		{
			indx = static_cast<int>(m_AcqCodeSequence[i]/spacing);
			indx = indx>=CodeSpan? CodeSpan-1: indx;

            CodeHist[indx] ++;
            indx = static_cast<int>(m_AcqFreqSequence[i]);
            FreqHist[indx] ++;
		}
		for( i=0; i<MAXSEQ; i++)
		{
			indx = static_cast<int>(m_AcqCodeSequence[i]/spacing);
            indx = indx>=CodeSpan? CodeSpan-1: indx;
			indx1 = static_cast<int>(m_AcqFreqSequence[i]);

            if(CodeHist[indx] >=5)                            // if the histogram of the code phase is large enough
            {
				condition1 = true;
			}
			if(FreqHist[indx1] >=5)                          // if the histogram of the freqency is large enough
			{
				condition2 = true;
			}