📄 16_qam.cpp
字号:
#include<stdio.h>
#include<math.h>
#include<ctime>
#include<cstdlib>
#define data_len 10000 //仿真所用数据总的比特数目
#define code_len data_len/4 //仿真所用数据总的码元数目
#define pi 3.1415926
//用结构体定义复数
struct complex
{
double Real; //实部
double Imag; //虚部
};
//A、B是复数，此函数用来计算|A-B|^2的值
double abs(struct complex a,struct complex b)
{
return((a.Real-b.Real)*(a.Real-b.Real)+(a.Imag-b.Imag)*(a.Imag-b.Imag));
};
void main()
{
float SNR=-18; //输入信噪比（dB）
double coef_AWGN=sqrt(5/pow(10,SNR/10)); //乘AWGN噪声（0，1）用于产生所需要的噪声功率
int data[data_len]; //存放仿真所用0 1 数据
int receive_data[data_len]; //存放经判决后的数据
srand(time(0)); //产生随机数的种子
//产生随机数存放在data数组里面
for(int i=0;i<data_len;i++)
{
data[i]=rand()%2;
}
//输出待测试的随机数
printf("The rest data:\n");
for(i=0;i<data_len;i++)
{
printf("%d\n",data[i]);
}
//printf("\n\n%d\n\n",code_len);
//定义结构体数组（由发送数据构成）
struct complex code[data_len/4];
/*-----constellation（星座图）------
Q
1011 1001 | 0001 0011
|
1010 1000 | 0000 0010
----------------|----------------I----
1110 1100 | 0100 0110
|
1111 1101 | 0101 0111 I1 Q1 I2 Q2 */
for(i=0;i<code_len;i++)//将数据进行QAM映射
{
if ((data[i*4]==1)&&(data[i*4+2]==1))
code[i].Real=-3;
else if((data[i*4]==1)&&(data[i*4+2]==0))
code[i].Real=-1;
else if((data[i*4]==0)&&(data[i*4+2]==0))
code[i].Real=1;
else
code[i].Real=3;
//printf("\n%f",code[i].Real);
if ((data[i*4+1]==1)&&(data[i*4+3]==1))
code[i].Imag=-3;
else if((data[i*4+1]==1)&&(data[i*4+3]==0))
code[i].Imag=-1;
else if((data[i*4+1]==0)&&(data[i*4+3]==0))
code[i].Imag=1;
else
code[i].Imag=3;
//printf("\n%f",code[i].Imag);
}
//输出映射后的数据
printf("the data after mapping\n");
for(i=0;i<code_len;i++)
{
printf("%f + j*%f;\n", code[i].Real,code[i].Imag);
}
//产生AWGN： Box-muller alog
double noise_I[2*code_len];
double noise_Q[2*code_len];
double Y[2*code_len];
double R[2*code_len];
noise_I[0]=1;
noise_Q[0]=2;
for(i=0;i<2*code_len-1;i++)
{
noise_I[i+1]=fmod(40014*(noise_I[i]),2147483563);
noise_Q[i+1]=fmod(noise_Q[i],2147483399);
Y[i]=fmod((noise_I[i]-noise_Q[i]),2147483562);
while(Y[i]<0)
{
Y[i]=Y[i]+2147483562;
}
if(Y[i]>0)
R[i]=Y[i]/2147483563;
else
R[i]=2147483562.0/2147483563.0;
//printf("\n%f",R[i]);
}
//i=2*code_len-1
Y[i]=fmod((noise_I[i]-noise_Q[i]),2147483562);
while(Y[i]<0)
{
Y[i]=Y[i]+2147483562;
}
if(Y[i]>0)
R[i]=Y[i]/2147483563;
else
R[i]=2147483562.0/2147483563.0;
for(i=0;i<code_len;i++)
{
noise_I[i]=sqrt(-2*log(R[i])/log(exp(1.0)))*cos(2*pi*R[i+code_len]);
noise_Q[i]=sqrt(-2*log(R[i])/log(exp(1.0)))*sin(2*pi*R[i+code_len]);
}
//printf("\n%f\n",coef_AWGN);
for(i=0;i<code_len;i++)
{
noise_I[i]=coef_AWGN*noise_I[i];
noise_Q[i]=coef_AWGN*noise_Q[i];
//printf("\n%f\n%f",noise_I[i],noise_Q[i]);
}
//定义结构体数组（由接收数据构成）
struct complex code_AWGN[code_len];
//在信号上叠加awgn噪声
for(i=0;i<code_len;i++)
{
code_AWGN[i].Real=code[i].Real+noise_I[i];
code_AWGN[i].Imag=code[i].Imag+noise_Q[i];
}
//将接收数据赋给结构体数组
struct complex receive[data_len/4];
for(i=0;i<code_len;i++)
{
receive[i].Real=code_AWGN[i].Real;
receive[i].Imag=code_AWGN[i].Imag;
}
//BI1_0:星座图中所有使得I1等于0的星座点
struct complex BI1_0[8]={{1,-3},{1,-1},{1,3},{1,1},{3,-3},{3,-1},{3,3},{3,1}};
//BI1_1:星座图中所有使得I1等于1的星座点
struct complex BI1_1[8]={{-3,3},{-3,1},{-3,-3},{-3,-1},{-1,3},{-1,1},{-1,-3},{-1,-1}};
//BI2_0:星座图中所有使得I1等于0的星座点
struct complex BI2_0[8]={{-1,-3},{-1,-1},{-1,1},{-1,3},{1,-3},{1,-1},{1,1},{1,3}};
//BI2_1:星座图中所有使得I1等于1的星座点
struct complex BI2_1[8]={{-3,-3},{-3,-1},{-3,1},{-3,3},{3,-3},{3,-1},{3,1},{3,3}};
//BQ1_0:星座图中所有使得I1等于0的星座点
struct complex BQ1_0[8]={{-3,3},{-3,1},{-1,3},{-1,1},{3,3},{3,1},{1,3},{1,1}};
//BQ1_1:星座图中所有使得I1等于1的星座点
struct complex BQ1_1[8]={{-3,-3},{-3,-1},{-1,-3},{-1,-1},{3,-3},{3,-1},{1,-3},{1,-1}};
//BQ2_0:星座图中所有使得I1等于0的星座点
struct complex BQ2_0[8]={{-3,1},{-3,-1},{-1,1},{-1,-1},{3,1},{3,-1},{1,-1},{1,1}};
//BQ2_1:星座图中所有使得I1等于1的星座点
struct complex BQ2_1[8]={{-3,-3},{-3,3},{-1,3},{-1,-3},{1,-3},{1,3},{3,-3},{3,3}};
//软输出
double LLR_BI1[code_len];
double LLR_BI2[code_len];
double LLR_BQ1[code_len];
double LLR_BQ2[code_len];
for(i=0;i<code_len;i++)
{
double temp;
double min_1;
double min_2;
min_1=abs(receive[i],BI1_1[0])*abs(receive[i],BI1_1[0]);
for(int j=1;j<8;j++)
{
temp=abs(receive[i],BI1_1[j])*abs(receive[i],BI1_1[j]);
if(temp<min_1)
min_1=temp;
}
min_2=abs(receive[i],BI1_0[0])*abs(receive[i],BI1_0[0]);
for(j=1;j<8;j++)
{
temp=abs(receive[i],BI1_0[j])*abs(receive[i],BI1_0[j]);
if(temp<min_2)
min_2=temp;
}
LLR_BI1[i]=(min_1-min_2)/4;
min_1=abs(receive[i],BI2_1[0])*abs(receive[i],BI2_1[0]);
for(j=1;j<8;j++)
{
temp=abs(receive[i],BI2_1[j])*abs(receive[i],BI2_1[j]);
if(temp<min_1)
min_1=temp;
}
min_2=abs(receive[i],BI2_0[0])*abs(receive[i],BI2_0[0]);
for(j=1;j<8;j++)
{
temp=abs(receive[i],BI2_0[j])*abs(receive[i],BI2_0[j]);
if(temp<min_2)
min_2=temp;
}
LLR_BI2[i]=(min_1-min_2)/4;
min_1=abs(receive[i],BQ1_1[0])*abs(receive[i],BQ1_1[0]);
for(j=1;j<8;j++)
{
temp=abs(receive[i],BQ1_1[j])*abs(receive[i],BQ1_1[j]);
if(temp<min_1)
min_1=temp;
}
min_2=abs(receive[i],BQ1_0[0])*abs(receive[i],BQ1_0[0]);
for(j=1;j<8;j++)
{
temp=abs(receive[i],BQ1_0[j])*abs(receive[i],BQ1_0[j]);
if(temp<min_2)
min_2=temp;
}
LLR_BQ1[i]=(min_1-min_2)/4;
min_1=abs(receive[i],BQ2_1[0])*abs(receive[i],BQ2_1[0]);
for(j=1;j<8;j++)
{
temp=abs(receive[i],BQ2_1[j])*abs(receive[i],BQ2_1[j]);
if(temp<min_1)
min_1=temp;
}
min_2=abs(receive[i],BQ2_0[0])*abs(receive[i],BQ2_0[0]);
for(j=1;j<8;j++)
{
temp=abs(receive[i],BQ2_0[j])*abs(receive[i],BQ2_0[j]);
if(temp<min_2)
min_2=temp;
}
LLR_BQ2[i]=(min_1-min_2)/4;
}
//判决，判决结构存入receive_data数组中
for(i=0;i<code_len;i++)
{
if(LLR_BI1[i]>0)receive_data[4*i]=0;
else receive_data[4*i]=1;
if(LLR_BQ1[i]>0)receive_data[4*i+1]=0;
else receive_data[4*i+1]=1;
if(LLR_BI2[i]>0)receive_data[4*i+2]=0;
else receive_data[4*i+2]=1;
if(LLR_BQ2[i]>0)receive_data[4*i+3]=0;
else receive_data[4*i+3]=1;
}
//将发送数据和接收数据同时输出
for(i=0;i<data_len;i++)
printf("\n%d%d",data[i],receive_data[i]);
//计算误码率
float error_rate; int j;
for(i=0,j=0;i<data_len;i++)
if(data[i]!=receive_data[i])j++;
error_rate=((float)j)/data_len;
printf("\n%f",error_rate);
}
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -