#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>
#include <string.h>
// BCH(7,4)编码参数
#define N_BCH 7
#define K_BCH 4
#define T_BCH 1 // 纠错能力
// 高斯随机数生成器（Box-Muller变换）
double gaussrand(double mean, double stddev) {
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.0 * U1 - 1.0;
V2 = 2.0 * U2 - 1.0;
S = V1 * V1 + V2 * V2;
} while (S >= 1.0 || S == 0.0);
X = V1 * sqrt(-2.0 * log(S) / S);
}
else {
X = V2 * sqrt(-2.0 * log(S) / S);
}
phase = 1 - phase;
return mean + X * stddev;
}
// BCH(7,4)编码函数
void bch_encode(const int* input, int* output, int num_blocks) {
// 生成矩阵 (系统形式)
int G[K_BCH][N_BCH] = {
{1, 0, 0, 0, 1, 1, 0},
{0, 1, 0, 0, 1, 0, 1},
{0, 0, 1, 0, 0, 1, 1},
{0, 0, 0, 1, 1, 1, 1}
};
for (int i = 0; i < num_blocks; i++) {
for (int j = 0; j < N_BCH; j++) {
output[i * N_BCH + j] = 0;
for (int k = 0; k < K_BCH; k++) {
output[i * N_BCH + j] ^= input[i * K_BCH + k] * G[k][j];
}
}
}
}
// BCH(7,4)解码函数（带纠错）
void bch_decode(const int* input, int* output, int num_blocks) {
// 校验矩阵
int H[3][7] = {
{1, 1, 0, 1, 1, 0, 0}, // 对应位置：0,1,3,4
{1, 0, 1, 1, 0, 1, 0}, // 对应位置：0,2,3,5
{0, 1, 1, 1, 0, 0, 1} // 对应位置：1,2,3,6
};
// 校正子到错误位置的映射（修正后的）
// 索引 = s0*4 + s1*2 + s2
int error_pattern[8] = { -1, 6, 5, 2, 4, 1, 0, 3 };
// -1: 无错误
// 0: 位置0错误
// 1: 位置1错误
// 2: 位置2错误
// 3: 位置3错误
// 4: 位置4错误
// 5: 位置5错误
// 6: 位置6错误
for (int i = 0; i < num_blocks; i++) {
int received[N_BCH];
memcpy(received, &input[i * N_BCH], N_BCH * sizeof(int));
// 计算校正子（异或运算）
int syndrome[3] = { 0 };
for (int j = 0; j < 3; j++) {
for (int k = 0; k < N_BCH; k++) {
if (H[j][k]) {
syndrome[j] ^= received[k];
}
}
}
// 计算错误位置索引 (0-7)
int error_index = syndrome[0] * 4 + syndrome[1] * 2 + syndrome[2];
int error_loc = error_pattern[error_index];
// 如果有错误且位置有效，则纠正
if (error_loc >= 0 && error_loc < N_BCH) {
received[error_loc] ^= 1; // 翻转错误位
}
// 提取信息位
for (int j = 0; j < K_BCH; j++) {
output[i * K_BCH + j] = received[j];
}
}
}
// 计算误码率
double calculate_ber(const int* original, const int* received, int length) {
if (length <= 0) return 0.0; // 避免除以零
int errors = 0;
for (int i = 0; i < length; i++) {
if (original[i] != received[i]) {
errors++;
}
}
return (double)errors / length;
}
int main() {
// 设置随机种子
srand((unsigned int)time(NULL));
// 参数初始化
const int n = 1000000; // 数据长度（比特数）
const double SNR_dB_start = 0.0;
const double SNR_dB_end = 10.0;
const double SNR_dB_step = 1.0; // 增加步长以加快测试
const int num_snr = (int)((SNR_dB_end - SNR_dB_start) / SNR_dB_step) + 1;
// 动态分配存储BER结果的数组
double* BER_encoded = (double*)malloc(num_snr * sizeof(double));
double* BER_theory = (double*)malloc(num_snr * sizeof(double));
double* BER_test = (double*)malloc(num_snr * sizeof(double));
if (!BER_encoded || !BER_theory || !BER_test) {
fprintf(stderr, "内存分配失败\n");
exit(1);
}
// 生成随机数据
int* data = (int*)malloc(n * sizeof(int));
for (int i = 0; i < n; i++) {
data[i] = rand() % 2;
}
// BCH编码参数
const int num_blocks = n / K_BCH; // 完整块数
// BCH编码
int* encoded_bits = (int*)malloc(num_blocks * N_BCH * sizeof(int));
if (!encoded_bits) {
fprintf(stderr, "编码内存分配失败\n");
exit(1);
}
bch_encode(data, encoded_bits, num_blocks);
const int encoded_length = num_blocks * N_BCH;
// BPSK调制（0->-1, 1->1）
double* BPSK_encoded = (double*)malloc(encoded_length * sizeof(double));
for (int i = 0; i < encoded_length; i++) {
BPSK_encoded[i] = encoded_bits[i] ? 1.0 : -1.0;
}
// 未编码数据的BPSK调制
double* BPSK_test = (double*)malloc(n * sizeof(double));
for (int i = 0; i < n; i++) {
BPSK_test[i] = data[i] ? 1.0 : -1.0;
}
printf("开始仿真...数据长度: %d bits, SNR范围: %.1f dB 到 %.1f dB\n",
n, SNR_dB_start, SNR_dB_end);
printf("编码方案: BCH(%d,%d), 纠错能力: %d bit\n", N_BCH, K_BCH, T_BCH);
// 主循环：针对不同的SNR计算BER
for (int z = 0; z < num_snr; z++) {
double snr_db = SNR_dB_start + z * SNR_dB_step;
double snr_lin = pow(10.0, snr_db / 10.0);
// 关键：确保相同Eb/N0下的公平比较
double noise_var = 1.0 / (2.0 * snr_lin); // 未编码的噪声方差
double noise_std = sqrt(noise_var);
// 编码系统的噪声方差（考虑码率损失）
double noise_var_encoded = 1.0 / (2.0 * snr_lin * (double)K_BCH / N_BCH);
double noise_std_encoded = sqrt(noise_var_encoded);
// ---------------------------------------------------
// 编码部分处理
// ---------------------------------------------------
double* noisy_encoded = (double*)malloc(encoded_length * sizeof(double));
int* demod_encoded = (int*)malloc(encoded_length * sizeof(int));
int* decoded_bits = (int*)malloc(num_blocks * K_BCH * sizeof(int));
// 添加高斯噪声
for (int i = 0; i < encoded_length; i++) {
noisy_encoded[i] = BPSK_encoded[i] + gaussrand(0.0, noise_std_encoded);
}
// BPSK解调（硬判决）
for (int i = 0; i < encoded_length; i++) {
demod_encoded[i] = (noisy_encoded[i] > 0.0) ? 1 : 0;
}
// BCH解码
bch_decode(demod_encoded, decoded_bits, num_blocks);
// 计算编码BER
BER_encoded[z] = calculate_ber(data, decoded_bits, num_blocks * K_BCH);
// 释放内存
free(noisy_encoded);
free(demod_encoded);
free(decoded_bits);
// ---------------------------------------------------
// 未编码部分处理
// ---------------------------------------------------
double* noisy_test = (double*)malloc(n * sizeof(double));
int* demod_test = (int*)malloc(n * sizeof(int));
// 添加高斯噪声
for (int i = 0; i < n; i++) {
noisy_test[i] = BPSK_test[i] + gaussrand(0.0, noise_std);
}
// BPSK解调（硬判决）
for (int i = 0; i < n; i++) {
demod_test[i] = (noisy_test[i] > 0.0) ? 1 : 0;
}
// 计算未编码BER
BER_test[z] = calculate_ber(data, demod_test, n);
// 计算理论BER (未编码)
BER_theory[z] = 0.5 * erfc(sqrt(snr_lin));
// 释放内存
free(noisy_test);
free(demod_test);
// 打印进度
printf("SNR = %4.1f dB: 编码BER = %8.3e, 未编码BER = %8.3e, 理论BER = %8.3e\n",
snr_db, BER_encoded[z], BER_test[z], BER_theory[z]);
}
// 输出结果文件（用于绘图）
FILE* fp = fopen("ber_results.txt", "w");
if (fp) {
fprintf(fp, "SNR_dB BER_theory BER_encoded BER_test\n");
for (int z = 0; z < num_snr; z++) {
double snr_db = SNR_dB_start + z * SNR_dB_step;
fprintf(fp, "%.1f %.6e %.6e %.6e\n",
snr_db, BER_theory[z], BER_encoded[z], BER_test[z]);
}
fclose(fp);
printf("\n结果已保存到 ber_results.txt\n");
// 分析编码增益
int coding_gain_point = -1;
for (int z = num_snr - 1; z >= 0; z--) {
if (BER_encoded[z] < BER_test[z]) {
coding_gain_point = z;
break;
}
}
if (coding_gain_point >= 0) {
double snr_db = SNR_dB_start + coding_gain_point * SNR_dB_step;
printf("-> 在 %.1f dB 处观察到编码增益: 编码BER=%.2e < 未编码BER=%.2e\n",
snr_db, BER_encoded[coding_gain_point], BER_test[coding_gain_point]);
}
else {
printf("-> 警告: 未检测到编码增益\n");
printf(" 可能原因:\n");
printf(" - 解码算法实现问题\n");
printf(" - SNR范围选择不当（需要在更高SNR区域仿真）\n");
printf(" - 数据长度不足（当前为%d比特）\n", n);
printf(" 建议: 尝试设置SNR_dB_end=15.0并增加数据长度\n");
}
}
else {
printf("保存结果到文件失败\n");
}
// 释放所有内存
free(BER_encoded);
free(BER_theory);
free(BER_test);
free(data);
free(encoded_bits);
free(BPSK_encoded);
free(BPSK_test);
return 0;
}
扫一扫
541
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
