📄 ldpc_decode.cpp
字号:
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include "time.h"
#define N 9216
#define K 4608
#define M 4608
double gaussrand()
{
static double V1, V2, S;
static int phase = 0;
double X;
if(phase == 0) {
do {
double U1 = (double)rand() / RAND_MAX;
double U2 = (double)rand() / RAND_MAX;
V1 = 2 * U1 - 1;
V2 = 2 * U2 - 1;
S = V1 * V1 + V2 * V2;
} while(S >= 1 || S == 0);
X = V1 * sqrt(-2 * log(S) / S);
} else
X = V2 * sqrt(-2 * log(S) / S);
phase = 1 - phase;
return X;
}
void encode(char *rs_encoded,char **G,char *ldpc_encoded_data)
{
int j,t;
for(j=0;j<M;j++)
for(t=0;t<K;t++)
ldpc_encoded_data[j] ^= (rs_encoded[t]&&G[j][t]); //校验位先存 G里面存的是(K*M)的转置
for(j=0;j<K;j++)
ldpc_encoded_data[M+j] = rs_encoded[j]; //信息位后存
}
void ldpc_decode(double *llr,char *ldpc_decoded)
{
int i,j,k,n;
double ram_init[36][256]={0};
double ram_ext[3][36][256]={0};
char ram_hard[36][256]={0};
char decision;
int sum_decision;
int neg_num;
int max_iter = 100;
// int biterrnum;
double min[2],temp,alpha = 0.8;
double sum1,sum2,sum3,temp_ext0,temp_ext1,temp_ext2;
// double sum_hardware[9216];
// FILE* fp;
/***水平更新时确定从哪个ram组(cnu单元)的哪个ram中取哪一个数***/
char ext_i[18][6]={{0,0,0,0,0,0},
{1,0,0,0,0,2},
{2,0,0,0,0,1},
{0,1,0,0,2,2},
{1,0,0,0,1,2},
{2,0,0,0,1,1},
{0,2,0,0,2,2},
{1,0,0,2,2,2},
{2,0,0,2,1,1},
{0,1,2,2,2,1},
{1,1,2,1,1,1},
{2,0,2,1,1,1},
{0,2,0,2,2,2},
{1,2,2,2,1,1},
{2,2,1,2,1,1},
{0,1,2,1,1,1},
{1,2,2,2,1,1},
{2,0,2,1,1,1}};
char vnu_i[18][6]={{0,6,12,18,25,30},
{0,7,19,26,31,12},
{0,8,13,20,32,26},
{1,6,14,21,25,31},
{1,15,27,33,20,8},
{1,9,16,34,25,21},
{2,6,28,35,16,20},
{2,10,17,27,21,30},
{2,11,22,22,16,34},
{3,7,26,13,33,19},
{3,9,28,17,31,18},
{3,23,9,12,35,24},
{4,7,29,17,34,18},
{4,24,10,23,13,32},
{4,29,10,32,23,15},
{5,8,14,22,28,30},
{5,19,11,15,33,27},
{5,24,35,11,14,29}};
unsigned char one_i[18][6]={{0,0,0,0,0,0},
{0,0,0,0,0,157},
{0,0,0,0,0,229},
{0,0,0,0,85,248},
{0,0,0,0,253,255},
{0,0,0,0,235,252},
{0,0,0,0,115,215},
{0,0,0,203,209,253},
{0,0,0,146,242,248},
{0,0,69,132,239,246},
{0,129,131,209,254,255},
{0,0,65,250,252,254},
{0,0,0,164,215,248},
{0,200,224,231,242,255},
{0,115,240,243,250,254},
{0,0,184,249,251,255},
{0,58,216,240,253,254},
{0,0,164,236,247,254}};
/***********存储初始化似然比***************/
for(i=0;i<256;i++)
for(j=0;j<36;j++)
ram_init[j][i] = llr[i*36+j];
for(n=0;;n++)
{
// biterrnum = 0; //for err test
sum_decision = 0;
// printf("ram_init[14][154]=%f,ram_ext0[14][154]=%f,ram_ext1[14][154]=%f,ram_ext2[14][154]=%f\n",ram_init[14][154],ram_ext[0][14][154],ram_ext[1][14][154],ram_ext[2][14][154]);
/************* 垂直更新********************/
for(i=0;i<256;i++)
for(j=0;j<36;j++)
{
temp_ext0 = ram_ext[0][j][i];
temp_ext1 = ram_ext[1][j][i];
temp_ext2 = ram_ext[2][j][i];
sum1 = ram_init[j][i] + ram_ext[0][j][i];
sum2 = ram_ext[1][j][i] + ram_ext[2][j][i];
sum3 = sum1 +sum2;
ram_ext[0][j][i] = ram_init[j][i] + sum2;
ram_ext[1][j][i] = temp_ext2 + sum1;
ram_ext[2][j][i] = temp_ext1 + sum1;
/* sum3 = ram_init[j][i] + ram_ext[0][j][i] + ram_ext[1][j][i] + ram_ext[2][j][i];
ram_ext[0][j][i] = sum3 - ram_ext[0][j][i];
ram_ext[1][j][i] = sum3 - ram_ext[1][j][i];
ram_ext[2][j][i] = sum3 - ram_ext[2][j][i];*/
ram_hard[j][i] = sum3 >= 0;
// if(ram_hard[j][i] != ldpc_encoded_data1[i*36+j]) //for err test
// biterrnum++; //for err test
// sum_hardware[i*36+j]=sum3;
}
// printf("after NO.%d,errornum is %d.\n",n+1,biterrnum); //for err test
// printf("sum=%f\n",sum_hardware[5558]);
/* fp = fopen("sum_hardware.txt","wb");
if(fp == NULL)
{
printf("file open error!\n");
return;
}
fwrite(sum_hardware,sizeof(double),9216,fp);
fclose(fp);*/
/**************水平更新******************/
for(i=0;i<256;i++)
for(j=0;j<18;j++)
{
decision = 0;
/****求最小值和负数的个数**********/
neg_num = 0;
if(ram_ext[ext_i[j][0]][vnu_i[j][0]][(one_i[j][0]+i)%256] < 0)
neg_num++;
min[1] = fabs(ram_ext[ext_i[j][0]][vnu_i[j][0]][(one_i[j][0]+i)%256]);
if(ram_ext[ext_i[j][1]][vnu_i[j][1]][(one_i[j][1]+i)%256] < 0)
neg_num++;
min[0] = fabs(ram_ext[ext_i[j][1]][vnu_i[j][1]][(one_i[j][1]+i)%256]);
if(min[0]<min[1]) /*最小值存在min[1],次小值存在min[0]*/
{
temp=min[0];
min[0]=min[1];
min[1]=temp;
}
for(k=2;k<6;k++)
{
if(ram_ext[ext_i[j][k]][vnu_i[j][k]][(one_i[j][k]+i)%256] < 0)
neg_num++;
if(fabs(ram_ext[ext_i[j][k]][vnu_i[j][k]][(one_i[j][k]+i)%256])<min[0])
{
min[0]=fabs(ram_ext[ext_i[j][k]][vnu_i[j][k]][(one_i[j][k]+i)%256]);
if(min[0]<min[1]) /*最小值存在min[1],次小值存在min[0]*/
{
temp=min[0];
min[0]=min[1];
min[1]=temp;
}
}
}
/*******找到最小值后更新每行**********/
for(k=0;k<6;k++)
{
decision ^= ram_hard[vnu_i[j][k]][(one_i[j][k]+i)%256];//对校验式进行验算
if(neg_num%2==0)
{
if(fabs(ram_ext[ext_i[j][k]][vnu_i[j][k]][(one_i[j][k]+i)%256])!=min[1])
ram_ext[ext_i[j][k]][vnu_i[j][k]][(one_i[j][k]+i)%256] = ((ram_ext[ext_i[j][k]][vnu_i[j][k]][(one_i[j][k]+i)%256] > 0) ? min[1]*alpha:((-1)*min[1]*alpha));
else
ram_ext[ext_i[j][k]][vnu_i[j][k]][(one_i[j][k]+i)%256] = ((ram_ext[ext_i[j][k]][vnu_i[j][k]][(one_i[j][k]+i)%256] > 0) ? min[0]*alpha:((-1)*min[0]*alpha));
}
else
{
if(fabs(ram_ext[ext_i[j][k]][vnu_i[j][k]][(one_i[j][k]+i)%256])!=min[1])
ram_ext[ext_i[j][k]][vnu_i[j][k]][(one_i[j][k]+i)%256] = ((ram_ext[ext_i[j][k]][vnu_i[j][k]][(one_i[j][k]+i)%256] > 0) ? ((-1)*min[1]*alpha):min[1]*alpha);
else
ram_ext[ext_i[j][k]][vnu_i[j][k]][(one_i[j][k]+i)%256] = ((ram_ext[ext_i[j][k]][vnu_i[j][k]][(one_i[j][k]+i)%256] > 0) ? ((-1)*min[0]*alpha):min[0]*alpha);
}
}
sum_decision += decision;
}
if(n < max_iter && sum_decision == 0)
{
printf("success at NO.%d iteration!\n ",(n+1));
for(i=0;i<256;i++)
for(j=0;j<36;j++)
ldpc_decoded[i*36+j] = ram_hard[j][i];
break;
}
if(n == max_iter)
{
printf("fail to decode after No.%d iteration!\n",max_iter);
for(i=0;i<256;i++)
for(j=0;j<36;j++)
ldpc_decoded[i*36+j] = ram_hard[j][i];
break;
}
}
}
int main(int argc, char* argv[])
//int sc_main(int argc, char* argv[])
{
char info[K];
char rs_encoded[N];
char rs_encoded1[N];
char ldpc_encoded_data[N];
char ldpc_encoded_data1[N];
char **G;
int i,j;
int errframe,temperr,ErrorTimes = 10;
long framenum,errbit;
int alpha = 1,MQAM = 1;
FILE *fp;
// sc_fixed<total_bits,integral_bits,SC_RND,SC_SAT> ldpc_ch_info[N];
double ldpc_ch_info[N];
// double ldpc_ch_info1[N]; //for err test
clock_t start, finish; //计算时间用的
double duration;
double Per,Ber;
float LDPC_RATE = 0.5;
float snr,sigma,noise_var;
int order[N]; //码字重排
// int initerrnum; //for err test
printf("input the snr:");
scanf("%f",&snr);
noise_var = alpha*alpha/(LDPC_RATE*pow(10,snr/10)*MQAM);// %复信号不需要*2
sigma = sqrt(noise_var/2);
errframe = 0;
errbit = 0;
framenum=0;
G = (char **)malloc(M*sizeof(char *));
for(i=0 ; i<M ; i++)
{
G[i] = (char *)calloc(K,sizeof(char));
}
//不同码率的生成矩阵不同，要根据码率分别读入。这里只读入了码率为1/2时的矩阵
fp = fopen("GG12.txt","rb"); //GG12按列存，读的时候其实读的是(K*M)维矩阵的转置(M*K)
if(fp==NULL)
{
printf("can not open the file!\n");
return 0;
}
for(j=0;j<M;j++)
fread(G[j],sizeof(char),K,fp);
fclose(fp);
fp = fopen("rorder.txt","rb");
if(fp==NULL)
{
printf("can not open the file!\n");
return 0;
}
fread(order,sizeof(int),N,fp);
fclose(fp);
start = clock();
while (errframe<ErrorTimes)
{
temperr=0;
// initerrnum=0;
framenum = framenum+1;
srand( (unsigned)time( NULL ) );
for(i=0;i<K;i++)
{
info[i]=(char)(rand()%2);
}
for(i=0;i<N;i++) //每次编码前对上一次的码字清零，不然会累积
ldpc_encoded_data[i]=0;
encode(info,G,ldpc_encoded_data); //LDPC编码
/****编码后码字重排****/
for(i=0;i<N;i++)
{
ldpc_encoded_data1[order[i]] = ldpc_encoded_data[i];
}
/******加入高斯噪声*******/
for(i=0;i<N;i++)
{
// ldpc_ch_info[i]=(sc_fixed<total_bits,integral_bits,SC_RND,SC_SAT>)(2*ldpc_encoded_data.data[i]-1+sigma*gaussrand());
ldpc_ch_info[i]=(double)(2*ldpc_encoded_data1[i]-1+sigma*gaussrand());
// cout<<ldpc_ch_info[i].to_string()<<endl;
}
/****************************************************************/
/*****************************************************************/
/***********for test*************/
/* fp=fopen("coded_data.txt","rb");
if(fp==NULL)
{
printf("can not open the file!\n");
}
fread(ldpc_encoded_data1,sizeof(char),N,fp);
fclose(fp);
fp = fopen("llr.txt","rb");
if(fp==NULL)
{
printf("Can not open the file!\n");
return 0;
}
fread(ldpc_ch_info,sizeof(double),N,fp);
fclose(fp);
*/
/* for(i=0;i<N;i++)
{
ldpc_encoded_data1[order[i]] = ldpc_encoded_data[i];
}
for(i=0;i<N;i++)
{
ldpc_ch_info[order[i]] = ldpc_ch_info1[i];
}*/
/**************************************************************/
/**************************************************************/
/***************************************************************/
/********test over*********/
/* for(i=0;i<N;i++)
if((ldpc_ch_info[i]>=0)!=ldpc_encoded_data1[i])
initerrnum++;
printf("initerrnum=%d\n",initerrnum);
*/
ldpc_decode(ldpc_ch_info,rs_encoded); //ldpc译码
/******译码后码字重排*******/
for(i=0;i<N;i++)
rs_encoded1[i] = rs_encoded[order[i]];
/*******计算错误bit数，误码率等*********/
for(i=0;i<N;i++)
if(ldpc_encoded_data[i]!=rs_encoded1[i])
temperr++;
errframe = errframe + (temperr>0);
errbit = errbit + temperr;
if(temperr!=0)
{
printf("total errframe is %d,the current error framenum is %ld\n ",errframe,framenum) ;
Ber = errbit/(double)(framenum*9216);
Per = (double)errframe/framenum;
fp = fopen("datarecord.txt","w");
if(fp==NULL)
{
printf("Can not open the file!\n");
return 0;
}
fprintf(fp ,"snr=%f ,framenum=%ld ,errframe=%d ,errbit=%ld ,ber=%f ,per=%f \n",snr,framenum,errframe,errbit,Ber,Per);
fclose(fp);
}
// if(framenum == 1)
// break;
}
printf("snr=%f ,framenum=%ld ,errframe=%d ,errbit=%ld ,ber=%f ,per=%f\n",snr,framenum,errframe,errbit,Ber,Per);
for( j=0 ; j<M ; j++)
{
free(G[j]);
}
free(G);
finish = clock();
duration = (double)(finish - start) / CLOCKS_PER_SEC;
printf( "%f seconds\n", duration );
return 0;
}
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -