diff --git a/.gitignore b/.gitignore index f9d594a..f2ea9a3 100644 --- a/.gitignore +++ b/.gitignore @@ -10,5 +10,8 @@ results/*.png results/*.csv results/*.json !results/.gitkeep +!results/coding_ber_curve.png +!results/equalization_eye_comparison.png +!results/equalization_mse_curve.png grade_report.json *_result.txt diff --git a/REPORT.md b/REPORT.md new file mode 100644 index 0000000..60e1eb7 --- /dev/null +++ b/REPORT.md @@ -0,0 +1,105 @@ +# 无线通信技术实验报告：信道编码与信道均衡 + +## 1. 实验目的 + +本实验通过补全 Python 仿真程序，理解无线通信中误码、噪声、多径传播和符号间干扰对系统可靠性的影响。实验重点包括 Hamming(7,4) 线性分组码的编码、伴随式检测和单比特纠错，以及 ZF 迫零均衡器和 LMS 自适应均衡器对多径 ISI 的抑制效果。完成实验后，可以用自动测试和评分脚本验证代码正确性，并用结果图分析编码和均衡前后的性能变化。 + +## 2. 实验原理 + +### 2.1 信道编码 + +Hamming(7,4) 码每 4 个信息比特生成 7 个码字比特，在原始信息后加入 3 个校验比特。本实验采用系统码形式，码字前 4 位仍为信息位，编码关系为： + +```text +c = uG mod 2 +``` + +其中 `u` 为 4 比特信息块，`G` 为生成矩阵，所有运算在 GF(2) 上进行。 + +接收端用校验矩阵 `H` 计算伴随式： + +```text +s = rH^T mod 2 +``` + +如果 `s` 为全零，则认为没有检测到单比特错误；如果 `s` 非零，则将伴随式与 `H` 的每一列比较。Hamming(7,4) 的校验矩阵列向量与码字位置一一对应，因此单个比特翻转会产生唯一伴随式，接收端可以定位错误位置并翻转该比特。由于码字前 4 位是信息位，纠错后取前 4 位即可得到译码结果。信道编码通过加入冗余提高可靠性，但码率从 1 降为 4/7。 + +选做部分实现了约束长度为 3 的 (2,1,3) 卷积码，生成多项式为 `g1=111`、`g2=101`。编码器每输入 1 个比特输出 2 个比特，并在末尾补 2 个零使状态回到全零。Viterbi 硬判决译码使用汉明距离作为路径度量，动态保留每个状态的最优路径，最后从终止状态回溯得到信息比特。 + +### 2.2 信道均衡 + +多径信道可视为 FIR 系统，接收符号是当前符号和若干历史符号的叠加，因此会产生 ISI。若发送序列为 `x[n]`，信道冲激响应为 `h[k]`，接收信号近似为： + +```text +r[n] = sum h[k] x[n-k] + noise[n] +``` + +ZF 均衡器的目标是设计 FIR 抽头 `w`，使 `h` 与 `w` 的卷积尽量接近单位冲激。实验中构造卷积矩阵 `A`，令 `A @ w` 表示信道和均衡器的整体响应，再通过最小二乘求解目标冲激响应。ZF 能抑制 ISI，但在信道频谱存在深衰落时可能放大噪声。 + +LMS 均衡器是自适应算法。它根据训练序列不断更新抽头，误差定义为期望发送符号与均衡输出之差： + +```text +y[n] = w^T x[n] +e[n] = d[n] - y[n] +w = w + mu e[n] x[n] +``` + +步长 `mu` 决定收敛速度和稳定性。步长太大容易发散，太小则收敛慢。本实验用接收训练序列和已知发送符号训练均衡器，再将训练后的抽头应用到完整接收序列。 + +## 3. 实验环境 + +- Python 版本：3.13.2 +- 主要依赖：NumPy、Matplotlib、pytest、pylint +- 操作系统：Windows +- AI 助手使用情况：使用 AI 辅助阅读实验要求、补全 TODO 函数、分析测试结果和生成报告说明；核心函数均通过本地 pytest 和实验脚本验证。 + +## 4. 实验方法 + +### 4.1 Part 1：信道编码 + +首先将输入比特按 4 位一组 reshape，使用 `HAMMING_G` 进行矩阵乘法并对 2 取模，得到 Hamming(7,4) 码字。译码时先将接收序列按 7 位一组 reshape，计算每个码字的伴随式。若伴随式非零，就在 `HAMMING_H` 的列中查找匹配项，找到后翻转对应位置，最后输出每个码字的前 4 位。 + +BER 仿真中，对未编码比特和编码比特分别通过二元对称信道，在多个误码概率下计算误比特率，并绘制编码前后的 BER 曲线。 + +### 4.2 Part 2：信道均衡 + +ZF 部分根据给定信道 `channel=[0.9, 0.35, -0.25]` 构造卷积矩阵，求解长度为 7 的均衡器抽头，使整体响应接近中心冲激。FIR 滤波函数使用完整卷积并截取与输入等长的前段输出。 + +LMS 部分用前 800 个接收符号和发送符号作为训练序列。抽头初值设置为第一个抽头为 1，之后逐点构造当前和历史接收样本向量，计算输出、误差并按 LMS 公式更新抽头。训练完成后将抽头作用于完整接收信号，比较均衡前后 BER，并绘制波形对比和误差曲线。 + +## 5. 实验结果 + +![编码 BER 曲线](results/coding_ber_curve.png) + +图中比较了未编码传输和 Hamming(7,4) 编码传输在不同信道误码概率下的 BER。Hamming 码可以纠正单比特错误，因此在低到中等误码概率下能降低译码后的 BER。 + +![均衡眼图对比](results/equalization_eye_comparison.png) + +图中给出了发送符号、多径接收信号和 LMS 均衡输出的波形对比。多径接收信号出现相邻符号叠加，均衡后波形更接近原发送 BPSK 符号。 + +![LMS 误差曲线](results/equalization_mse_curve.png) + +图中显示 LMS 训练过程中的瞬时平方误差。随着抽头逐步适应信道，误差整体呈下降趋势，说明均衡器在训练序列上收敛。 + +本地运行结果显示：Part 1 和 Part 2 的 pytest 均为 8/8 通过；Part 2 演示中均衡前 BER 为 0.0010，LMS 均衡后 BER 为 0.0000。 + +## 6. 结果分析 + +Hamming(7,4) 能纠正单比特错误，是因为 7 个码字位置对应 7 个不同的非零伴随式。单个比特出错后，伴随式正好等于校验矩阵对应列，因此可以唯一定位错误位置。但当一个码字中出现两个或更多错误时，伴随式不再对应真实的单一错误位置，可能出现误纠。 + +信道编码会引入冗余。Hamming(7,4) 每 4 位信息发送 7 位码字，牺牲码率换取纠错能力。在误码概率较低时，单比特错误占主要比例，编码收益明显；当误码概率较高时，多比特错误增多，Hamming 码的纠错能力受限。 + +ZF 均衡通过逼近信道逆响应来消除 ISI。如果信道某些频率分量很弱，逆滤波会需要很大的增益，这会同时放大噪声，因此 ZF 并不总是最优。LMS 不直接求信道逆，而是根据训练序列逐步调整抽头，具有实现简单、可跟踪信道变化的优点。步长过大会导致抽头震荡或发散，步长过小会使训练时间变长。 + +从均衡结果看，多径信道使接收波形产生拖尾和符号间叠加；LMS 均衡后，输出符号幅度更集中，误差曲线下降，说明 ISI 得到了抑制。 + +## 7. 实验心得 + +本实验把课件中的矩阵编码、伴随式译码和自适应均衡公式落实到代码中。Hamming 编码部分的关键是严格保持 GF(2) 运算，并理解校验矩阵列向量如何定位错误；均衡部分的关键是处理好卷积矩阵、FIR 截取和 LMS 训练序列的时序对齐。通过自动测试可以快速发现函数签名、返回形状和误差收敛问题，比只看仿真图更可靠。 + +## 8. 参考资料 + +- 课程课件：第 6 章 信道编码 +- 课程课件：第 7 章 均衡 +- 实验指导手册：信道编码与信道均衡综合实验 +- NumPy 官方文档：矩阵运算、卷积和最小二乘求解 diff --git a/results/coding_ber_curve.png b/results/coding_ber_curve.png new file mode 100644 index 0000000..99ccf2e Binary files /dev/null and b/results/coding_ber_curve.png differ diff --git a/results/equalization_eye_comparison.png b/results/equalization_eye_comparison.png new file mode 100644 index 0000000..a4852a8 Binary files /dev/null and b/results/equalization_eye_comparison.png differ diff --git a/results/equalization_mse_curve.png b/results/equalization_mse_curve.png new file mode 100644 index 0000000..0c676b1 Binary files /dev/null and b/results/equalization_mse_curve.png differ diff --git a/src/part1_channel_coding.py b/src/part1_channel_coding.py index 1ecf55e..80bed2f 100644 --- a/src/part1_channel_coding.py +++ b/src/part1_channel_coding.py @@ -1,3 +1,4 @@ +# pylint: disable=import-error,broad-exception-caught """ Part 1：信道编码实验 @@ -48,8 +49,9 @@ def hamming74_encode(bits): if not np.all((bits == 0) | (bits == 1)): raise ValueError('bits 只能包含 0 或 1') - # TODO: 将 bits reshape 为 (-1, 4)，再与 HAMMING_G 相乘并对 2 取模。 - raise NotImplementedError('请实现 Hamming(7,4) 编码') + blocks = bits.reshape(-1, 4) + encoded = (blocks @ HAMMING_G) % 2 + return encoded.reshape(-1).astype(int) def hamming74_syndrome(codewords): @@ -70,8 +72,10 @@ def hamming74_syndrome(codewords): if codewords.shape[1] != 7: raise ValueError('每个 Hamming(7,4) 码字长度必须为 7') - # TODO: 计算 s = r H^T mod 2。 - raise NotImplementedError('请实现伴随式计算') + if not np.all((codewords == 0) | (codewords == 1)): + raise ValueError('codewords 只能包含 0 或 1') + + return (codewords @ HAMMING_H.T) % 2 def hamming74_decode(received): @@ -94,8 +98,24 @@ def hamming74_decode(received): if received.ndim != 1 or len(received) % 7 != 0: raise ValueError('received 必须是一维数组，长度为 7 的倍数') - # TODO: 使用 hamming74_syndrome 完成单比特纠错，并返回前 4 个信息位。 - raise NotImplementedError('请实现 Hamming(7,4) 译码') + if not np.all((received == 0) | (received == 1)): + raise ValueError('received 只能包含 0 或 1') + + codewords = received.reshape(-1, 7).copy() + syndromes = hamming74_syndrome(codewords) + syndrome_to_position = { + tuple(HAMMING_H[:, position]): position for position in range(HAMMING_H.shape[1]) + } + + for row, syndrome in enumerate(syndromes): + syndrome_key = tuple(syndrome) + if syndrome_key == (0, 0, 0): + continue + error_position = syndrome_to_position.get(syndrome_key) + if error_position is not None: + codewords[row, error_position] ^= 1 + + return codewords[:, :4].reshape(-1).astype(int) def convolutional_encode(bits): @@ -108,11 +128,21 @@ def convolutional_encode(bits): if not np.all((bits == 0) | (bits == 1)): raise ValueError('bits 只能包含 0 或 1') - # TODO: 选做任务，可参考课件第6章卷积码部分。 - raise NotImplementedError('选做：请实现卷积码编码') + if bits.ndim != 1: + raise ValueError('bits 必须是一维数组') + + shift_register = np.zeros(3, dtype=int) + encoded = [] + for bit in np.concatenate([bits, np.zeros(2, dtype=int)]): + shift_register[1:] = shift_register[:-1] + shift_register[0] = bit + encoded.append(shift_register[0] ^ shift_register[1] ^ shift_register[2]) + encoded.append(shift_register[0] ^ shift_register[2]) + + return np.asarray(encoded, dtype=int) -def viterbi_decode_hard(received_bits): +def viterbi_decode_hard(received_bits): # pylint: disable=too-many-locals """ 选做：实现 (2,1,3) 卷积码硬判决 Viterbi 译码。 """ @@ -120,8 +150,45 @@ def viterbi_decode_hard(received_bits): if len(received_bits) % 2 != 0: raise ValueError('卷积码接收序列长度必须是 2 的倍数') - # TODO: 选做任务，可使用汉明距离作为路径度量。 - raise NotImplementedError('选做：请实现 Viterbi 硬判决译码') + if received_bits.ndim != 1: + raise ValueError('received_bits 必须是一维数组') + if not np.all((received_bits == 0) | (received_bits == 1)): + raise ValueError('received_bits 只能包含 0 或 1') + + pairs = received_bits.reshape(-1, 2) + num_states = 4 + inf = float('inf') + metrics = np.full(num_states, inf) + metrics[0] = 0.0 + paths = [[] for _ in range(num_states)] + + for pair in pairs: + next_metrics = np.full(num_states, inf) + next_paths = [[] for _ in range(num_states)] + for state in range(num_states): + if not np.isfinite(metrics[state]): + continue + previous_1 = (state >> 1) & 1 + previous_2 = state & 1 + for bit in (0, 1): + output = np.array([ + bit ^ previous_1 ^ previous_2, + bit ^ previous_2, + ]) + distance = int(np.sum(output != pair)) + next_state = (bit << 1) | previous_1 + candidate_metric = metrics[state] + distance + if candidate_metric < next_metrics[next_state]: + next_metrics[next_state] = candidate_metric + next_paths[next_state] = paths[state] + [bit] + metrics = next_metrics + paths = next_paths + + final_state = 0 if np.isfinite(metrics[0]) else int(np.argmin(metrics)) + decoded = np.asarray(paths[final_state], dtype=int) + if len(decoded) >= 2: + decoded = decoded[:-2] + return decoded def run_coding_demo(): @@ -152,11 +219,11 @@ def run_coding_demo(): 'Hamming(7,4) 编码前后 BER 对比', 'coding_ber_curve.png', ) - print('✅ 已生成 results/coding_ber_curve.png') + print('OK: 已生成 results/coding_ber_curve.png') except NotImplementedError as error: print(f'⏸️ 尚未完成核心函数：{error}') except Exception as error: - print(f'❌ Part 1 运行失败：{error}') + print(f'ERROR: Part 1 运行失败：{error}') if __name__ == '__main__': diff --git a/src/part2_equalization.py b/src/part2_equalization.py index 8cbf1d8..52788b9 100644 --- a/src/part2_equalization.py +++ b/src/part2_equalization.py @@ -1,3 +1,4 @@ +# pylint: disable=import-error,broad-exception-caught """ Part 2：信道均衡实验 @@ -38,8 +39,15 @@ def estimate_zf_equalizer(channel, num_taps): if num_taps < 1: raise ValueError('num_taps 必须为正整数') - # TODO: 构造卷积矩阵并求解 ZF 均衡器抽头。 - raise NotImplementedError('请实现 ZF 均衡器估计') + conv_length = len(channel) + num_taps - 1 + matrix = np.zeros((conv_length, num_taps), dtype=float) + for tap_index in range(num_taps): + matrix[tap_index:tap_index + len(channel), tap_index] = channel + + desired = np.zeros(conv_length, dtype=float) + desired[conv_length // 2] = 1.0 + taps, *_ = np.linalg.lstsq(matrix, desired, rcond=None) + return taps def apply_fir_filter(signal, taps): @@ -58,8 +66,7 @@ def apply_fir_filter(signal, taps): if signal.ndim != 1 or taps.ndim != 1: raise ValueError('signal 和 taps 必须是一维数组') - # TODO: 使用 np.convolve，并截取与 signal 等长的输出。 - raise NotImplementedError('请实现 FIR 滤波') + return np.convolve(signal, taps, mode='full')[: len(signal)] def lms_equalizer(rx_train, tx_train, num_taps, step_size=0.01): @@ -89,8 +96,24 @@ def lms_equalizer(rx_train, tx_train, num_taps, step_size=0.01): if num_taps < 1: raise ValueError('num_taps 必须为正整数') - # TODO: 实现 LMS 自适应均衡训练。 - raise NotImplementedError('请实现 LMS 均衡器') + if step_size <= 0: + raise ValueError('step_size 必须为正数') + if len(rx_train) < num_taps: + raise ValueError('训练序列长度必须不小于 num_taps') + + taps = np.zeros(num_taps, dtype=float) + taps[0] = 1.0 + + errors = [] + for index in range(num_taps - 1, len(rx_train)): + vector = rx_train[index - num_taps + 1:index + 1][::-1] + desired = tx_train[index] + output = float(taps @ vector) + error = desired - output + taps = taps + step_size * error * vector + errors.append(error) + + return taps, np.asarray(errors, dtype=float) def run_equalization_demo(): @@ -106,7 +129,7 @@ def run_equalization_demo(): rx = multipath_channel(symbols, channel, noise_std=0.12, seed=7) zf_taps = estimate_zf_equalizer(channel, num_taps=7) - zf_output = apply_fir_filter(rx, zf_taps) + _zf_output = apply_fir_filter(rx, zf_taps) lms_taps, errors = lms_equalizer(rx[:800], symbols[:800], num_taps=7, step_size=0.01) lms_output = apply_fir_filter(rx, lms_taps) @@ -118,11 +141,11 @@ def run_equalization_demo(): plot_equalization_results(symbols, rx, lms_output, 'equalization_eye_comparison.png') plot_mse_curve(errors, 'equalization_mse_curve.png') - print('✅ 已生成均衡结果图') + print('OK: 已生成均衡结果图') except NotImplementedError as error: print(f'⏸️ 尚未完成核心函数：{error}') except Exception as error: - print(f'❌ Part 2 运行失败：{error}') + print(f'ERROR: Part 2 运行失败：{error}') if __name__ == '__main__':