📄 channelfading.c
字号:
/*********************************************************************
* 文件名称： ChannelFading.c
* 主要内容： GSM系统六径衰落信道仿真。采用Jakes衰落模型。
* 本程序Model A选用第二组数据。
The reduced setting (6 taps) is defined thereunder.
model A (Typical case for hilly terrain)
number Relative Average relative doppler
time (μs) power (dB) spectrum
(1) (2) (1) (2)
1 0.0 0.0 0.0 0.0 CLASS
2 0.1 0.2 -1.5 -2.0 CLASS
3 0.3 0.4 -4.5 -4.0 CLASS
4 0.5 0.6 -7.5 -7.0 CLASS
5 15.0 15.0 -8.0 -6.0 CLASS
6 17.2 17.2 -17.7 -12.0 CLASS
model B (Profile for equalization test)
Tap Relative Average relative doppler
number time (μs) power (dB) spectrum
1 0.0 0.0 CLASS
2 3.2 0.0 CLASS
3 6.4 0.0 CLASS
4 9.6 0.0 CLASS
5 12.8 0.0 CLASS
6 16.0 0.0 CLASS
*=====================================================================
* 版权：作者个人所有
*
* 当前版本： 0.1.3
* 作 者： 金利忠，蒋登峰, 高霁
* 完成日期： 2004年5月25日
*=====================================================================
**********************************************************************/
#include "ChannelFading.h"
/*void main()
{
double in_Xdata[FRAME_LEN];
double in_Ydata[FRAME_LEN];
double out_Xdata[FRAME_LEN+maxDELAY];
double out_Ydata[FRAME_LEN+maxDELAY];
int i;
in_Xdata[0]=1;
in_Ydata[0]=1;
for (i = 1; i<FRAME_LEN; i++)
{
in_Xdata[i] =0; // generate the random data
in_Ydata[i] =0;
//printf("%f,%f\n", in_Xdata[i],in_Ydata[i]);
}
JakesModel(in_Xdata,in_Ydata,out_Xdata,out_Ydata);
for (i = 0; i<FRAME_LEN+maxDELAY; i++)
{
printf("%f,,,%f---Time:%d\n", out_Xdata[i],out_Ydata[i],i);
}
}*/
/*****************************************************
函数名：void JakesModel(const double *Xcomplex, const double *Ycomplex,
double *oXcomplex, double *oYcomplex)
函数功能：每一径平坦衰落和延时后合成函数。
输入参数说明：
Xcomplex: 输入数据实部分量，double型指针
Ycomplex: 输入数据虚部分量，double型指针
输出参数说明：
oXcomplex:输出数据实部分量，double型指针
返回值：oYcomplex:输出数据虚部分量，double型指针
无
其它说明：这个程序中，直接对输出数据进行处理，没有在数据处理是额外的寄存器的需求。
具体可参考程序文档。
修改日期 版本号 修改人 修改内容
------------------------------------------------------
2004/05/25 V1.0 金利忠 创建
2004/05/26 V1.1 金利忠
2004/05/29 V1.2 金利忠
2004/05/30 V1.3 金利忠
******************************************************/
void JakesModel(const double *Xcomplex, const double *Ycomplex,
double *oXcomplex, double *oYcomplex, int lasttime[1],
unsigned int gseed[2],double SNR,const int chchoice,
const int vchoice,const double chfac)
{
/*****************************************
先对模型中用到的参数初始化。
*****************************************/
int jj,n,i;
int j;
int delay[]={0,0,0,1,16,19};//延时初始化
double power_rela[]={0,-2.0,-4.0,-7.0,-6.0,-12.0};//相对功率的初始化
//For debug
//int delay[TAP_NUMBER]={0,0,0,1,16,19};//延时初始化
//double power_rela[6]={0,-200,-200,-200,-200,-1200};//相对功率的初始化
double power[6];
int p;
double sum_power=0.0;
double gM[FRAME_LEN];
double gN[FRAME_LEN];
double Xcc[1];
double Xss[1];
double outX;
double outY;
int tap_num; //TAP计数器
double sigma;
int maxDELAY = 19;
double jakesWM;
jakesWM = (2*PI*vchoice*1000.0*9)/(3600.0*3.0);
if (chchoice == 2)
{
delay[0] = 0;
delay[1] = 3;
delay[2] = 7;
delay[3] = 10;
delay[4] = 14;
delay[5] = 17;
power_rela[0] = 0;
power_rela[1] = 0;
power_rela[2] = 0;
power_rela[3] = 0;
power_rela[4] = 0;
power_rela[5] = 0;
maxDELAY = 17;
}
/* for(jj=1;jj<JAKES_N0+1;jj++) //取得gama和sita
{
for(i=1;i<TAP_NUMBER+1;i++)
{
beta[jj-1][i-1]=PI*jj*i/(JAKES_N0+1);
gama[i-1]=(2*PI*(i-1))/(JAKES_N0+1);
sita[jj-1][i-1]=beta[jj-1][i-1]+gama[i-1];
}
}
for(n=1;n<JAKES_N+1;n++) //取得alfa和Wn
{
alfa[n-1]=(2*PI*n)/JAKES_N;
JAKES_Wn[n-1]=jakesWM*cos(alfa[n-1]);
}
*/
for (p=0; p<TAP_NUMBER; p++)
{
power[p]=pow(10.0,power_rela[p]/10);
sum_power=sum_power+power[p]; //取相对功率的和。
}
/*****************************************
产生高斯噪声。
*****************************************/
G_rand(gM,gN,gseed);
sigma = pow(10.0,(SNR/10));
/*****************************************
开始进行多径衰落合成。
*****************************************/
//初始化输出数据前maxDELAY位
for(j=0;j<maxDELAY;j++)
{
oXcomplex[j]=0;
oYcomplex[j]=0;
}
Xcc[0]=0;
Xss[0]=0;
if (lasttime[0]>=1073741824) //prevent overflow
{
lasttime[0] = 0;
}
for(j=0;j<FRAME_LEN;j++) //j时刻
{
//kkkk=Xcomplex[j];
oXcomplex[j+maxDELAY]=Xcomplex[j];
oYcomplex[j+maxDELAY]=Ycomplex[j];
//最大延时的一径先算
h_gen(TAP_NUMBER, j+lasttime[0], Xcc, Xss, jakesWM);
outX=sqrt(power[TAP_NUMBER-1]/sum_power)*oXcomplex[j];
outY=sqrt(power[TAP_NUMBER-1]/sum_power)*oYcomplex[j];
oXcomplex[j]=outX*(*Xcc)-outY*(*Xss);
oYcomplex[j]=outY*(*Xcc)+outX*(*Xss);
//再从第一径开始算，一共剩下5径
for(tap_num=1;tap_num<TAP_NUMBER;tap_num++)
{
h_gen(tap_num, j+lasttime[0], Xcc, Xss, jakesWM);
outX=sqrt(power[tap_num-1]/sum_power)*oXcomplex[j+maxDELAY-delay[tap_num-1]];
outY=sqrt(power[tap_num-1]/sum_power)*oYcomplex[j+maxDELAY-delay[tap_num-1]];
oXcomplex[j]=(outX*(Xcc[0])-outY*(Xss[0]))+oXcomplex[j];
oYcomplex[j]=(outY*(Xcc[0])+outX*(Xss[0]))+oYcomplex[j];
}
oXcomplex[j]=gM[j]/(sqrt(2*sigma))+oXcomplex[j]/(sqrt(chfac));//加入噪声M[j]
oYcomplex[j]=gN[j]/(sqrt(2*sigma))+oYcomplex[j]/(sqrt(chfac));//加入噪声N[j]
}
lasttime[0]=j+lasttime[0];
}
/*****************************************************
函数名：void h_gen(const int tap, const int tim, double *Xc, double *Xs)
函数功能：产生每一径的衰落系数h(j,tap)，使用的是Jakes模型。
输入参数说明：
tap: 输入多径的径数，int型
tim: 输入时间，int型
输出参数说明：
oXcomplex:输出数据实部分量，double型指针
oYcomplex:输出数据虚部分量，double型指针
返回值：
无
其它说明
修改日期 版本号 修改人 修改内容
------------------------------------------------------
2004/05/25 V1.0 金利忠 创建
2004/05/26 V1.1 金利忠
2004/05/29 V1.2 金利忠
2004/05/30 V1.3 金利忠
******************************************************/
void h_gen(const int tap, const int tim, double Xc[1], double Xs[1], double jakesWM)
{
double XsS=0;
double XcS=0;
double cos_w;
double cos_alfa;
double factor_n;
double Nalpha=PI; //待定
//double alpha_hat=0.5*PI; //待定
double betan[JAKES_N0];
int i,jj,n;
//Classical Jakes
/* for(jj=1;jj<JAKES_N0+1;jj++) //取得gama和sita
{
betan[jj-1] = (PI*jj)/(JAKES_N0);
}
for(jj=1;jj<JAKES_N0+1;jj++) //取得gama和sita
{
for(i=1;i<TAP_NUMBER+1;i++)
{
beta[jj-1][i-1]=(PI*jj)/(JAKES_N0+1);
gama[i-1]=(2*PI*(i-1))/(JAKES_N0+1);
sita[jj-1][i-1]=beta[jj-1][i-1]+gama[i-1];
}
}
for(n=1;n<JAKES_N+1;n++) //取得alfa和Wn
{
alfa[n-1]=(2*PI*n)/JAKES_N;
JAKES_Wn[n-1]=(jakesWM)*cos(alfa[n-1]);
}
for(i=1;i<JAKES_N0+1;i++)
{
cos_w=cos(JAKES_Wn[i-1]*tim*(SYMBLE_T1)+sita[i-1][tap-1]);
//XcS=XcS+2*cos(beta[i-1][tap-1])*cos_w;
//XsS=XsS+2*sin(beta[i-1][tap-1])*cos_w;
XcS=XcS+2*cos(betan[i-1])*cos_w;
XsS=XsS+2*sin(betan[i-1])*cos_w;
}
cos_alfa=cos(jakesWM*tim*(SYMBLE_T1)+(tap)*(PI/2));
factor_n=sqrt(2.0/JAKES_N);
//XcS=2*XcS+sqrt(2)*fabs(cos(Nalpha*tap))*cos_alfa;
//XsS=2*XsS+sqrt(2)*sin(Nalpha*tap)*cos_alfa;
XcS = XcS+cos_alfa;
XsS = XsS+cos_alfa;
Xc[0]=factor_n*XcS;
Xs[0]=factor_n*XsS;
*/
//Modified (Ref to IEEE)
for(jj=1;jj<JAKES_N0+1;jj++) //取得gama和sita
{
for(i=1;i<TAP_NUMBER+1;i++)
{
beta[jj-1][i-1]=(PI*jj*i)/(JAKES_N0+1);
gama[i-1]=(2*PI*(i-1))/(JAKES_N0+1);
sita[jj-1][i-1]=beta[jj-1][i-1]+gama[i-1];
}
}
for(n=1;n<JAKES_N+1;n++) //取得alfa和Wn
{
alfa[n-1]=(2*PI*n)/JAKES_N;
JAKES_Wn[n-1]=(jakesWM)*cos(alfa[n-1]);
}
for(i=1;i<JAKES_N0+1;i++)
{
cos_w=cos(JAKES_Wn[i-1]*tim*(SYMBLE_T1)+(tap*PI/2));
XcS=XcS+2*cos(beta[i-1][tap-1])*cos_w;
XsS=XsS+2*sin(beta[i-1][tap-1])*cos_w;
//XcS=XcS+2*cos(betan[i-1])*cos_w;
//XsS=XsS+2*sin(betan[i-1])*cos_w;
}
cos_alfa=cos(jakesWM*tim*(SYMBLE_T1)+(tap)*(PI/2));
factor_n=sqrt(2.0/JAKES_N);
XcS=XcS+2*cos_alfa*fabs(cos(PI*tap));
XsS=XsS+2*cos_alfa*sin(PI*tap);
//XcS = XcS+cos_alfa;
//XsS = XsS+cos_alfa;
Xc[0]=factor_n*XcS;
Xs[0]=factor_n*XsS;
}
/*void h_gen(const int tap, const int tim, double Xc[1], double Xs[1])
{
double XsS=0;
double XcS=0;
double cos_w;
double cos_alfa;
double factor_n;
double Nalpha=PI; //待定
double alpha_hat=0.5*PI; //待定
int i,jj,n;
double betan[JAKES_N0];
//printf("%d\n",JAKES_N0);
for(jj=1;jj<JAKES_N0+1;jj++) //取得gama和sita
{
betan[jj-1] = (PI*jj)/(JAKES_N0);
}
for(n=1;n<JAKES_N+1;n++) //取得alfa和Wn
{
alfa[n-1]=(2*PI*n)/JAKES_N;
JAKES_Wn[n-1]=jakesWM*cos(alfa[n-1]);
}
for(i=1;i<JAKES_N0+1;i++)
{
cos_w=cos(JAKES_Wn[i-1]*tim*(SYMBLE_T1));//+alpha_hat*tap);//sita[i-1][tap-1]);
XcS=XcS+2*cos(betan[i-1])*cos_w;
XsS=XsS+2*sin(betan[i-1])*cos_w;
}
cos_alfa=cos(jakesWM*tim*SYMBLE_T1);//+alpha_hat*tap);//gama[tap-1]
factor_n=sqrt(2.0/JAKES_N);
//XcS=2*XcS+2*fabs(cos(Nalpha*tap))*cos_alfa;//sqrt(2)
//XcS=2*XcS+2*cos(Nalpha*tap)*cos_alfa;//sqrt(2)
//XsS=2*XsS+2*sin(Nalpha*tap)*cos_alfa;//sqrt(2)
XcS = XcS+cos_alfa;
XsS = XsS+cos_alfa;
Xc[0]=factor_n*XcS;
Xs[0]=factor_n*XsS;
}*/
/*****************************************************
函数名：void G_rand(double *M, double *N)
函数功能：产生高斯随机变量，采用Box-Muller方法。
无输出参数
输出参数说明：
M: 输出复数X变量，double型
N: 输出复数Y变量，double型
返回值：
无
其它说明：均值为0，方差为1
修改日期 版本号 修改人 修改内容
------------------------------------------------------
2004/05/27 V1.0 高霁 创建
2004/05/30 V1.1 金利忠 原型：if (ii == 1)
{
X[1] = 999987256;
Y[1] = 198102025;
}
******************************************************/
void G_rand(double *M, double *N, unsigned int gseed[2])
{
unsigned long int X[FRAME_LEN];
unsigned long int Y[FRAME_LEN];
double XX[FRAME_LEN];
double YY[FRAME_LEN];
unsigned long int a = 16807;
unsigned long int b = 0;
unsigned long int c = 2147483647;
int ii;
double d,e,f,g;
for(ii = 0;ii < FRAME_LEN;ii++)
{
if (ii == 0)
{
X[0] = gseed[0];
Y[0] = gseed[1];
}
else
{
X[ii] = (a * X[ii-1] + b)%c;
Y[ii] = (a * Y[ii-1] + b)%c;
}
XX[ii] = X[ii] / (double)c;
YY[ii] = Y[ii] / (double)c;
d=sqrt(-2 * log(XX[ii]));
e=cos(2 * PI * YY[ii]);
f=sqrt(-2 * log(XX[ii]));
g=sin(2 * PI * YY[ii]);
M[ii] = d * e;
N[ii] = f * g;
}
gseed[0] = X[FRAME_LEN - 1];
gseed[1] = Y[FRAME_LEN - 1];
}
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -