极大似然技术编码译码
1 实验环境
语言:C++
2 实验内容
2.1 编程实现 $(m,n)$ 群码
2.1.1 输入输出
输入:$m,n$ 和一致性校验矩阵 $H$ 中的 $H_{m \times r}$ 部分,其中 $r=n-m$.
输出:所有的码字,按原码对应的十进制整数升序排列.
2.1.2 设计思路及核心算法
首先输入 $m,n$:
1
std::cin >> m >> n;
调用
inputH(m,n)
函数输入一致性校验矩阵 $H$:1
BitMatrix H = inputH(m, n);
inputH(m,n)
函数代码如下:1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17uint64_t inputBit() {
std::string line;
std::cin >> line;
uint64_t x = 0;
int k = 0;
for (size_t i = 0; i < line.size(); i++)
if (line[i] == '1' || line[i] == '0')
x |= (line[i] - '0') << k++;
return x;
}
BitMatrix inputH(int m, int n) {
BitMatrix H(n, n - m);
for (int i = 0; i < m; i++) H.d[i] = inputBit();
for (int i = m; i < n; i++) H.d[i] = 1 << (i - m);
return H;
}inputH
函数传入行列参数m,n,声明一个BitMatrix
类,用来存放一致性校验矩阵的 $m \times r$部分.调用
HtoH2(H)
函数生成 $m \times n$ 的编码矩阵:
$$
\left[
\begin{matrix}
I_m & H_{m\times r} \
\end{matrix}
\right]
$$1
2
3
4
5
6
7
8BitMatrix HtoH2(BitMatrix &H) {
BitMatrix H2(H.n - H.m, H.n);
for (int i = 0; i < H.n - H.m; i++) H2.d[i] = 1 << i;
for (int i = 0; i < H.n - H.m; i++)
for (int j = 0; j < H.m; j++)
H2.d[i] |= (H.d[i] >> j & 1) << (j + H.n - H.m);
return H2;
}定义一个
Bm
矩阵用于存放原码,矩阵每行从上到下按原码对应的十进制递增顺序排列,并进行运算:$$
e_H(B^m)=B^m\times
\left[
\begin{matrix}
I_m & H_{m\times r} \
\end{matrix}
\right]
$$求出群码对应的矩阵并输出.
1
2
3
4
5
6
7
8
9
10
11
12
13BitMatrix H2 = HtoH2(H), Bm(1 << m, m);
for (int i = 0; i < (1 << m); i++)
for (int j = 0; j < m; j++)
Bm.d[i] |= (i >> j & 1) << (m - j - 1);
BitMatrix N = Bm * H2;
std::cout << "1. e(B^m):" << std::endl;
for (int i = 0; i < (1 << m); i++) {
std::cout << " e(";
printBit(Bm.d[i], m);
std::cout << ") = ";
printBit(N.d[i], n);
std::cout << std::endl;
}
2.2 构造陪集表
2.2.1 设计思路及核心算法
生成 coset table
的算法原理如下:
首先,定义了一个空的哈希表
cosets
用于存储 coset 的信息,以及一个空的哈希表visited
用于标记已经访问过的向量;创建一个向量
numbers
,其中存储了所有可能的向量(共 2^n 个向量),并按照向量中 1 的个数进行排序,从少到多;遍历
numbers
中的每个向量i
,对于每个向量,执行以下步骤:a. 创建一个空的 coset 对象
coset
,用于存储属于同一个 coset 的向量;b. 对于矩阵
N
的每一行(表示每个 coset),将当前向量i
与该行的向量x
进行异或运算,得到新的向量x' = N.d[j] ^ uint64_t(i)
;c. 检查
visited
哈希表,如果在visited
中找到了向量x'
,说明该向量已经属于某个 coset,跳过当前循环;d. 如果在
visited
中没有找到向量x'
,则将x'
添加到当前 coset 中,并将x'
添加到visited
哈希表中,值为当前 coset 的索引;e. 将当前 coset 的陪集头
e
更新为 coset 中具有最少 1 的个数的向量;f. 如果当前 coset 的陪集头
e
在cosets
哈希表中不存在,则将当前 coset 添加到cosets
中,并将e
对应的特征值(通过一致性校验矩阵H
计算)添加到syncdrome
哈希表中,值为当前 coset 的索引;g. 检查
cosetLeaders
的大小,如果已经达到了 $2^{n-m}$(即余类的数量),则停止循环.最后,得到了所有的 coset 信息,并按照 coset 陪集头对应的十进制递增的顺序输出.
通过上述算法,可以生成所有的 coset,并找到每个 coset 的陪集头和特征值.
2.2.2 BitMatrix
类
该程序将译码时所需的比特矩阵进行了封装,能够实现乘法 (元素对应异或) 操作.
1 | class BitMatrix { |
2.2.3 Coset
类
Coset
类由陪集、陪集头及添加操作组成,陪集中的元素使用一维向量进行存储,陪集头用整型变量单独存储,方便调用. 添加操作在每一次执行时更新一次陪集头.
1 | class Coset { |
2.3 极大似然译码
基于 2.2 构造陪集表 求出的左陪集表,用和函数 $e_H$ 相关联的极大似然译码函数 $d$ ,对 $x_t\in B^n$ 进行译码.
2.3.1 输入输出格式
- 输入:首行是译码个数 $N$ ,整数格式;
- 然后是 $N$ 行比特串 $x_t \in B^n$ ,字符串格式;
- 输出:$N$行译码结果,格式为 $d(x_t) = b_1 b_2 … b_m$ . 如:
d(001011)=001
.
2.3.2 设计思路及核心算法
输入 $N$ ,然后输入待解码的 $N$ 个比特串. 输入后需要进行字符串转整数操作;
1
2
3
4
5
6
7std::cin >> k;
for (int i = 0; i < k; i++) {
std::cerr << "Enter the word #" << i + 1 << ": ";
uint64_t x = inputBit();
/*
......................
*/inputBit()
函数如下:1
2
3
4
5
6
7
8
9
10uint64_t inputBit() {
std::string line;
std::cin >> line;
uint64_t x = 0;
int k = 0;
for (size_t i = 0; i < line.size(); i++)
if (line[i] == '1' || line[i] == '0')
x |= (line[i] - '0') << k++;
return x;
}确定输入的 $x_t$ 位于哪个
coset
中,可通过visited[x]
得到所属陪集编号;取出所属陪集的陪集头 $\epsilon$ ,与 $x_t$ 进行异或运算,得到的结果必然在群码对应的矩阵 $N$ 中找到,由此得到解码后结果.
核心代码如下:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26std::vector<uint64_t > cosetLeaders;
for (auto i : numbers) {
Coset coset;
bool notFound = true;
for (int j = 0; j < N.n; j++) {
uint64_t x = N.d[j] ^ uint64_t(i);
if (visited.find(x) != visited.end()) {
notFound = false;
break;
}
visited[x] = cosetLeaders.size();
coset.add(x);
}
if (notFound) {
cosets[coset.e] = coset;
syndrome[getSyndrome(coset.e, H)] = cosetLeaders.size();
cosetLeaders.push_back(coset.e);
}
// if cosetLeaders.size() = 2^(n-m), then stop.
// Because there are only 2^(n-m) cosets.
if (cosetLeaders.size() == (1 << (n - m))) break;
}
/*
............
*/
printBit(Bm.d[NtoBm[cosetLeaders[visited[x]] ^ x]], m);
2.4 基于特征值的译码
基于 2.2 构造陪集表 求出的左陪集表,利用函数 $f_H:B^n \rightarrow B^r,f_H(x) = x\cdot H$ 求出每个陪集头对应的特征值 (Syndrome) ,并利用特征值,对$x_t \in B^n$ 进行译码.
2.4.1 输入输出格式
- 输入:首行是译码个数N,整数格式;
- 然后是N行比特串 $x_t \in B^n$ ,字符串格式;
- 输出:N行译码所用的特征值、陪集头和译码结果,格式为 $f(x_t)$ = *** ,$\epsilon$ = ***,$d(xt) = b_1b_2…b_m$.
注:本程序采用的输出方式,对于每一个待解码二进制串输入,输出结果分为两行,第一行为采用极大似然译码计算的结果,结果前用序号 (a) 进行区分,第二行为基于特征值译码的结果,结果前用序号 (b) 区分,例如:
1 | (a) d(10111) = 10 |
2.4.2 设计思路及核心算法
输入 N ,然后输入待解码的 N 个比特串. 输入后需要进行字符串转整数操作;
通过函数 $f_H(x_t) = x_t \cdot H$ 计算 $x_t$ 特征值;
1
2
3
4
5std::cout << " (b) f("; //特征值
printBit(x, n);
uint64_t s = getSyndrome(x, H);
std::cout << ") = ";
printBit(s, n - m);getSyndrome()
函数细节:1
2
3
4
5uint64_t getSyndrome(uint64_t x, BitMatrix &H) {
BitMatrix X(1, H.n);
X.d[0] = x;
return (X * H).d[0];
}根据所得特征值在哈希表中对应到陪集头;
1
2std::cout << ", ε = "; //陪集头
printBit(cosetLeaders[syndrome[s]], n);陪集头与待译码比特串进行异或操作,得到的结果必然在群码对应的矩阵 $N$ 中找到,由此得到解码后结果.
1
2
3
4
5std::cout << ", d("; //译码结果
printBit(x, n);
std::cout << ") = ";
printBit(Bm.d[NtoBm[cosetLeaders[syndrome[s]] ^ x]], m);
std::cout << std::endl;
3 代码优化
我们在四个地方优化了原本的代码,使得程序更加高效.
由于实验中满足 $1<m<n<7$,所以使用
uint64_t
替代数组,并使用位运算进行相应计算,同时降低了空间复杂度和时间复杂度;其中,矩阵乘法的核心代码如下:
1
2
3
4
5// main.cpp, Line 30-33
for (int i = 0; i < n; i++)
for (int k = 0; k < m; k++)
if (d[i] >> k & 1)
res.d[i] ^= rhs.d[k];普通的矩阵乘法的时间复杂度是 $O(n^3)$,而使用位运算后,时间复杂度降低到了 $O(n^2)$.
当已经找到了所有的陪集时,就可以停止循环,不需要再继续遍历;
1
2// main.cpp, Line 132
if (cosetLeaders.size() == (1 << (n - m))) break;使用哈希表记录每一个码字所属的陪集,在解码的时候可以直接查找,而不需要再次遍历;
1
2// main.cpp, Line 123, 记录每一个码字所属的陪集
visited[x] = cosetLeaders.size();1
2// main.cpp, Line 155, 直接解码
printBit(Bm.d[NtoBm[cosetLeaders[visited[x]] ^ x]], m);该哈希表还可以快速地判断当前陪集是否重复,减少了计算量;
1
2
3
4
5// main.cpp, Line 119-122, 判断当前陪集是否重复
if (visited.find(x) != visited.end()) {
notFound = false;
break;
}考虑到 $1<m<n<7$,所以码字的范围最大不超过 $2^8 = 256$,所以实际上可以使用数组代替哈希表,降低部分常数.
特征值解码方式也使用了类似的方式;
1
2// main.cpp, Line 128, 记录每一个特征值对应的陪集(头)
syndrome[getSyndrome(coset.e, H)] = cosetLeaders.size();1
2// main.cpp, Line 155, 直接获取特征值对应的陪集头并进行异或运算
printBit(Bm.d[NtoBm[cosetLeaders[syndrome[s]] ^ x]], m);使用了
lowbit
快速计算最低位的 1 的位置,加快了bitcount
的速度.1
2// main.cpp, Line 17
while (n) n -= n & -n, s++;
4 实验成果示例
1 | $ ./main # 运行 |
附录:源代码
1 |
|