modulatormex.c

本程序利用改进的CHOW算法

#include "mex.h"
#include "mat.h"
#include <math.h>
#define  bits_out plhs[0]
#define	bits_in prhs[0]
#define baud_in prhs[1]
#define	tag_in prhs[2]

//function bits_out=transmittermex(bits_in,baudrate,tag);

//Pass in a message to be digitally modulated
//For FSK,BPSK {0,1} for QPSK {0,1,2,3} for 16PSK {0,1,...,15}. (bits_in)
//Also, specify the baudrate. Finally, specify the type of modulation to use.
//tag = 0 (FSK), 1 (BPSK), 2 (QPSK), 3 (16PSK), 4 (16QAM) - unimplemented, 5 (DPSK) unimplemented
//returns modulated waveform (output_bits)

// This function was written by me, Joe Williams, a Stanford EE Grad Student. It is used for my 
// EE359 project (Wireless Comm). You may use this code in whatever manner you want as long as
// you give me credit.

void mexFunction(int nlhs, 
                     mxArray *plhs[],
                     int nrhs, 
                     const mxArray *prhs[])
    {
	
	  double *output_bits,*output_temp;
	  double *input_bits; 
      MATFile *ph;
	  double baudrate,fa,fb,Fs;
      int tag,samplesperbit,i,j,inputlength;
	  float f0,f1,wave0[100000],wave1[100000],wave2[100000],wave3[100000],wave[16][100000];
      double *var;
    
	  baudrate=mxGetScalar(baud_in);
      fa=900000000;
	  fb=899900000;
	  tag=mxGetScalar(tag_in);
	  Fs=9000000*20; // * 100
      input_bits = mxGetPr(bits_in);
	  inputlength = mxGetN(bits_in);
      //samplesperbit=9000; 
      samplesperbit=Fs/(baudrate/100);
	  mexPrintf("%d",samplesperbit);
      bits_out = mxCreateDoubleMatrix(1, samplesperbit*inputlength, mxREAL);
      output_bits = mxGetPr(bits_out);
	   
	  f0=(float)fa/(float)Fs*.02*3.1415926539;
	  f1=(float)fb/(float)Fs*.02*3.1415926539;
	 
	  switch(tag){//construct 1 period of signal
	  case 0://FSK
	    for(i=0;i<samplesperbit;i++)
	      wave0[i]=sin(f0*(float)i);
	    for(i=0;i<samplesperbit;i++)
	      wave1[i]=sin(f1*(float)i);
	    break;
	  case 1://BPSK
		for(i=0;i<samplesperbit;i++)
	      wave0[i]=sin(f0*(float)i)*sin(2*3.1415926539*.5/samplesperbit*i);
	    for(i=0;i<samplesperbit;i++)
	      wave1[i]=sin(f0*(float)i+3.1415926539)*sin(2*3.1415926539*.5/samplesperbit*i);
	    break;
	  case 2://QPSK
		for(i=0;i<samplesperbit;i++)
		  wave0[i]=sin(f0*(float)i)+cos(f0*(float)i);
		for(i=0;i<samplesperbit;i++)
		  wave1[i]=-sin(f0*(float)i)+cos(f0*(float)i);
		for(i=0;i<samplesperbit;i++)
		  wave2[i]=-sin(f0*(float)i)-cos(f0*(float)i);
		for(i=0;i<samplesperbit;i++)
		  wave3[i]=sin(f0*(float)i)-cos(f0*(float)i);
        break;
	  case 3://16PSK
		for(j=0;j<16;j++){
		  for(i=0;i<samplesperbit;i++)
	        //wave[j][i]=sin(f0*(float)i+(j)*3.1415926539/8);
		    wave[j][i]=sin(f0*(float)i+(j)*3.1415926539/8)*(sin(2*3.1415926539*.5/samplesperbit*i));
		}
      }		

	  //construct signal
	  if((tag==0)||(tag==1)){//FSK and BPSK
        for(i=0;i<inputlength;i++){
		  for(j=0;j<samplesperbit;j++){
			if(input_bits[i]==0)
              output_bits[i*samplesperbit+j]=wave0[j];
			else if(input_bits[i]==1)
			  output_bits[i*samplesperbit+j]=wave1[j];
          }
		}
	  }
	  else if(tag==2){//QPSK
        for(i=0;i<inputlength;i++){
		  for(j=0;j<samplesperbit;j++){
            if(input_bits[i]==0)
			  output_bits[i*samplesperbit+j]=wave0[j];
			else if(input_bits[i]==1)
			  output_bits[i*samplesperbit+j]=wave1[j];
			else if(input_bits[i]==2)
			  output_bits[i*samplesperbit+j]=wave2[j];
			else if(input_bits[i]==3)
			  output_bits[i*samplesperbit+j]=wave3[j];
		  }
		}
	  }
	  else if(tag==3){//16PSK
        for(i=0;i<inputlength;i++){
		  for(j=0;j<samplesperbit;j++)
		    output_bits[i*samplesperbit+j]=wave[(int)(input_bits[i])][j];
		}
	  }
}