2023BUAA信息论期末大作业（Python实现信源-编码-信道-译码-信宿通信过程）

山舟

已于 2023-09-19 22:46:30 修改

阅读量1.1k

点赞数 5

文章标签： python 开发语言

于 2023-06-10 23:51:47 首次发布

本文链接：https://blog.csdn.net/weixin_51983604/article/details/130659193

版权

文章目录

前言
零、总述
一、信源
二、信源编码
三、信道编码
- 1.(n,k)线性码
- 2.简单重复编码
四、信道
五、信道译码
- 1.(n,k)线性码
- 2.简单重复编码
六、信源译码
七、整合
感谢阅读，如有错误请批评指正

前言

以下内容是我和 @Rg猿、 @murmurto 两位大佬共同完成的，所有内容仅是我们对作业要求和课程内容的理解，难免存在纰漏和错误，发布在这里仅供参考，切勿抄袭。

零、总述

代码主要由如下部分构成：

在这里插入图片描述

每部分由（简略的）原理说明、代码和注释以及各部分的测试构成。

正文开始。

一、信源

1.无记忆信源

无记忆信源要求序列中的每个符号在生成时是独立的，即序列中每个符号在生成时符号的概率分布互相没有影响。

程序中每次调用ZeroMemoryInformationSource函数生成一个符号，每次调用时符号的概率分布是通过随机数生成的，因此互相没有影响，是无记忆的。

def ZeroMemoryInformationSource(symbols):
    """
    :param
        symbols: 信源符号集(list)
    :return:
        生成的符号(np.array)、符号对应的概率分布(np.array)
    """
    sign_n = 1  # 生成符号的个数，每次调用生成1个符号
    seq_array = np.empty(sign_n, dtype=int)  # 生成的符号序列，长度=sign_n
    sign_and_probability = []  # 符号集和概率的对应矩阵

    # 生成当前符号的概率分布
    # 由于sign_n=1，因此下面的循环仅会执行一次
    for k in range(sign_n):
        sign_probability = []  # 产生序列中第k个符号时，信源生成每个符号对应的概率
        for i in range(len(symbols)):
            initial_probability = random.random()  # 随机生成的概率
            sign_probability.append(initial_probability)
        total_probability = sum(sign_probability)  # 所有概率之和
        for i in range(len(symbols)):
            sign_probability[i] /= total_probability  # 归一化

        new_sign = np.random.choice(symbols, size=1, p=sign_probability)  # 生成序列的第k个符号
        seq_array[k] = new_sign[0]  # 符号添加到序列

		# 为了后续编码方便，规定以这种格式输出
        sign_and_probability.append(symbols)
        sign_and_probability.append(sign_probability)

	# 						转置
    return seq_array, np.transpose(np.array(sign_and_probability))

上面的函数每调用一次生成1个符号，因此要得到一个N长的序列需要调用N次，下面以N=5，符号集为[1,2,3]为例进行说明。

print('--------------------Zero Memory Start----------------------------')
for i in range(5):
    seq, sp = ZeroMemoryInformationSource([1, 2, 3])
    print(seq)  # 信源输出的序列
    print(sp)  # 每个符号对应的概率
print('---------------------Zero Memory End-----------------------------')

调用3次函数，每次输出一个符号，最终生成的符号序列即为[1, 1, 1, 3, 2]。

关于输出概率分布的格式，使用的是numpy的array，由于内部数据的格式统一，为保证概率为浮点，符号也输出为浮点格式，即下图中的1. 2. 3.分别代表符号1、2、3.

在生成序列第一个符号1时符号集[1, 2, 3]的概率分布是[0.2386, 0.7420, 0.0193]；
在生成序列第二个符号1时符号集[1, 2, 3]的概率分布是[0.5420, 0.2584, 0.1994]；
以此类推，在生成序列第五个符号2时符号集[1, 2, 3]的概率分布是[0.4048, 0.4958, 0.0992]。

显然，在生成每个符号时符号集对应的概率分布是相互无关的，信源是无记忆信源。

在这里插入图片描述

2.平稳信源

平稳信源要求符号的概率分布在时间上是不变的。

在代码中设计一个StationaryInformationSource类，其成员为symbols和probability，分别代表符号集和对应的概率分布。实例化时类内的概率分布作为固有属性，不会发生变化。

同时类内包含一个Generate方法用于生成符号序列，每次按照固定的概率分布选出一个符号进行输出。

class StationaryInformationSource:
    def __init__(self, symbols, probability=None):
        """
        初始化
        :param symbols: 符号集
        :param probability: 每个符号的概率分布（稳定不变）
        """
        self.symbols = symbols
        # 如果没有传入概率分布，则随机生成
        if probability is None:
            sign_probability = []  # 符号的概率分布
            for i in range(len(symbols)):
                # 随机生成概率
                initial_probability = random.random()
                sign_probability.append(initial_probability)

            # 所有概率之和
            total_probability = sum(sign_probability)
            for i in range(len(symbols)):
                sign_probability[i] /= total_probability  # 归一化得到每个符号的概率

            self.probability = sign_probability
        else:
            self.probability = probability

    def Generate(self):
        """
        输出符号序列
        :param
        :return:
            信源序列(np.array)、每个符号对应的概率(np.array)
        """
        sign_n = 1  # 每次生成一个符号
        seq_array = np.empty(1, dtype=int)  # 序列
        sign_and_probability = []

        # 生成符号集
        for i in range(sign_n):
            new_sign = np.random.choice(self.symbols, size=sign_n, p=self.probability)
            seq_array[i] = new_sign[0]

		
        sign_and_probability.append(self.symbols)
        sign_and_probability.append(self.probability)

        return seq_array, np.transpose(np.array(sign_and_probability))

下面简单地测试一下这段代码，先实例化出一个对象s，然后循环调用Generate方法生成一段序列，循环的次数就是生成序列的长度。

print('--------------------Stationary Start----------------------------')
s = StationaryInformationSource([1, 2, 3])  # 实例化
for i in range(5):
    seq, sp = s.Generate()
    print(seq)
    print(sp)
print('---------------------Stationary End-----------------------------')

输入的符号集为[1, 2, 3]，输出的序列为[3, 2, 2, 3, 1]，每次输出符号时符号集的概率均是[0.2261, 0.2860, 0.4877]，不随时间变化，是平稳的。

在这里插入图片描述

3.马尔可夫信源

马尔可夫信源要求从一个状态（原概率分布）通过转移矩阵（条件概率）进入下一个状态（新的概率分布）。下面代码中实现的是一阶马尔可夫过程，即当前状态仅与前一个状态有关。

在代码中设计MarkovInformationSource类，类的成员为

states_list：符号集
probability_matrix：转移矩阵
initial_state：符号的概率分布

其中states_list为必须传入的参数，后两个参数为可选参数，如果没有传入就在类的初始化函数内自动生成，否则就使用传入的。

类内的Generate方法用于生成符号序列，同时通过当前的概率分布和概率转移矩阵生成下一时刻的概率分布，不断更新，最后整个符号集的分布会趋于稳定，信源进入稳态。

class MarkovInformationSource:
    def __init__(self, states_list, probability_matrix=None, initial_state=None):
        """
        初始化马尔可夫信源

        Args:
            states_list: 状态集(list)
            probability_matrix: 状态转移概率矩阵，如果为None会随机生成(list)
            initial_state: 每个状态的概率，如果为None会随机生成(list)
        """
        self.states_list = states_list
        n = len(states_list)  # 状态数量

        # 没有传入状态转移概率矩阵，在初始化函数内生成
        if probability_matrix is None:
            # n*n状态转移概率矩阵
            self.probability_matrix = np.random.random((n, n))

            for i in range(n):
                total = sum(self.probability_matrix[i])
                for j in range(n):
                    self.probability_matrix[i, j] /= total
        else:
            self.probability_matrix = probability_matrix

        # 没有传入初始状态，在初始化函数内生成
        if initial_state is None:
            temp = []
            # 每个状态生成随机概率
            for i in range(len(states_list)):
                temp.append(random.random())

            # 归一化
            total = sum(temp)
            for i in range(len(states_list)):
                temp[i] /= total
            self.initial_state = temp
        else:
            self.initial_state = initial_state

    def Generate(self):
        """
        每次生成一个符号
        :param
        :return:
            生成的符号(np.array)、符号集及其概率分布(np.array)
        """
        sign_n = 1
        seq_array = np.empty(sign_n, dtype=int)  # 序列
        probability_matrix = []
        sign_and_probability = []

        # 生成长度为sign_n的马尔可夫序列
        for k in range(sign_n):
            # 生成符号
            new_sign = np.random.choice(self.states_list, size=1, p=self.initial_state)
            # 状态数量
            n = len(self.states_list)
            new_state_probability = []  # 转移后的概率
            # 概率分布与转移状态矩阵相乘得到下一时刻的概率分布
            for i in range(n):
                temp = 0
                for j in range(n):
                    temp += self.initial_state[j] * self.probability_matrix[j, i]
                new_state_probability.append(temp)

            self.initial_state = new_state_probability

            # 生成符号和对应的概率
            seq_array[k] = new_sign[0]

			# 输出符号集及其概率分布
            sign_and_probability.append(self.states_list)
            sign_and_probability.append(self.initial_state)

        return seq_array, np.transpose(np.array(sign_and_probability))

仍旧用一段代码来测试这个类的功能，实例化一个符号集为[1, 2, 3, 4]的类对象m，循环调用其Generate方法生成20个符号的序列。

print('--------------------Markov Start----------------------------')
m = MarkovInformationSource([1, 2, 3, 4])
for i in range(20):
    seq, sp = m.Generate()
    print(seq)
    print(sp)
print('---------------------Markov End-----------------------------')

输出的符号及对应的概率分布如下图所示，由于输出过多，仅截取前几个输出的符号。

在这里插入图片描述

上图是前几个输出符号的概率分布，可以看出每一个概率值都在向一个值趋近，也就是说这个马尔可夫过程在逐渐进入稳态，上面代码中生成了长20的符号序列，下图是输出的最后几个符号的概率分布，可以看到已经进入稳态，完全平稳。

在这里插入图片描述

二、信源编码

信源编码将信源的各种符号编为不同的码符号序列，这里仅考虑二元码符号集0、1，即编码结果仅是0、1的序列。

编码的大致过程如下图所示：
在这里插入图片描述

下面的代码在整体框架上都是相同的，不同之处就在于各自编码的实现过程。

1.等长编码

（1）原理

假定符号集 $S$ 长度为 $n$ ，则编码输出集 $C$ 中各个码字长度 $l$ 满足

$\log_2 n \le l < 1 + log_2 n$
将信源发送的符号序列编码至整数集 $=1,2,3,\dots,n$ ，再将该整数集中各个元素值对应的二进制数作为编码输出。

对符号编至整数集，再将整数集中各个元素值对应的二进制码作为编码输出。

（2）等长编码节点类EqualNode

建立等长编码的节点类，节点中存放符号的编号Id、符号对应的等长编码结果Code，实例化时需要给定其编号。

class EqualNode(object):
    """
    等长编码节点
    """

    def __init__(self, id):
        """
        等长编码节点初始化函数__init__
        输入：符号id
        输出：无
        功能：初始化符号编码Code
        """
        self.Id = id
        self.Code = np.array([], dtype=int)

（3）等长编码接口函数EqualLength

下面根据输入的符号和概率分布矩阵（由信源输出）来生成EqualNode类型的节点。首先按照符号集的符号对其按升序排列（这里信源输入的符号全部为整型或浮点，可以直接排序；如果信源输入的符号为字符、字符串等不能直接排序的数，则可用其在符号集中的下标排序），然后遍历符号和概率分布矩阵，生成EqualNode节点放入数组中。之后即可调用EqualLengthCoding进行编码。

def EqualLength(symbol_matrix):
    """
    等长码接口EqualLength:
    输入：Nx2 信源概率矩阵symbol_matrix N为符号个数 第一列为符号名称 第二列为符号概率
    输出：EqualLengthNode类的数组arr
    功能：返回EqualLengthNode类数组arr，统一形式
    """

    # 按符号集的Id升序排列
    len_symbol_matrix = np.size(symbol_matrix, 0)

    # 存EqualNode类的数组
    arr = np.array([])
    for i in range(len_symbol_matrix):
        m = EqualNode(symbol_matrix[i][0])  # 生成节点
        arr = np.append(arr, m)  # 存入arr

	# 等长编码
    EqualLengthCoding(symbol_matrix, arr)

    return SortById(arr)

（4）等长编码函数EquallengthCoding

下面EqualLengthCoding函数根据符号集的个数进行编码，首先确定符号个数，以此来确定二进制编码的位数，将符号集的id（十进制）编码为二进制数，并将所有编码扩展到同样的二进制位数（等长编码）并存放到对应节点的Code中。

def EqualLengthCoding(symbol_matrix, arr):
    """
    等长码生成函数EqualLengthCoding:
    输入：Nx2 信源概率矩阵symbol_matrix N为符号个数 第一列为符号名称 第二列为符号概率
        EqualLengthNode类的数组arr
    功能：对输入概率矩阵中的符号进行等长编码
    """
    len_symbol_matrix = np.size(symbol_matrix, 0)  # 符号序列的行数
    symbol_max = len_symbol_matrix  # 取符号的最大值
    n = math.ceil(math.log2(symbol_max))  # 2进制数的位数
    for i in range(len_symbol_matrix):
        binary = bin(int(symbol_matrix[i][0])).replace("0b", "")  # 符号变为二进制数
        len_binary = len(binary)
        # 把二进制数补满n位
        for j in range(n):
            if j < n - len_binary:
                arr[i].Code = np.append(arr[i].Code, 0)
            else:
                arr[i].Code = np.append(arr[i].Code, int(binary[j - n + len_binary]))

（5）等长编码测试函数EqualLengthCodeDisplay

下面对等长编码进行测试，输出符号及其对应的编码序列。

def EqualLengthCodeDisplay(source_probability):
    """
    等长编码结果展示函数
    输入：Nx2 信源概率矩阵source_probability N为符号个数 第一列为符号名称 第二列为符号概率
    输出：无
    功能：展示信源等长编码结果
    """

    code_equal_length = EqualLength(source_probability)
    for k in range(0, np.size(source_probability, 0)):
        print(code_equal_length[k].Id, code_equal_length[k].Code)

# 测试用信源概率矩阵
p_source = np.array([[2, 0.10], [1, 0.08], [3, 0.18], [7, 0.01],
                     [5, 0.15], [6, 0.11], [8, 0.12], [4, 0.16],
                     [9, 0.09]])
print("===============EqualLengthCode===============")
print("---------------SourceAfterCode---------------")
EqualLengthCodeDisplay(p_source)

可以看到符号1 ~ 9被对应编码成了对应二进制数0001 ~ 1001，注意一定要对不够位数的二进制数补齐位数，例如符号1的编码直接产生为1，一定要补齐为0001，这样才满足等长编码。

在这里插入图片描述

2.费诺编码

（1）原理

费诺编码先把当前符号集的概率降序排列，然后在某两个符号间划分成两部分，使两部分各自的概率和之差最小，此时对一部分编码为0，另一部分编码为1，如此循环不断分割，直到每一部分都只剩一个符号。具体过程可参考下图：

在这里插入图片描述

（2）费诺编码节点FanoNode

先建立FanoNode类，表示一个结点，其中包含符号的ID及其对应的编码，实例化时给定其编号。

class FanoNode(object):
    """
    费诺编码节点
    """

    def __init__(self, id):
        """
        费诺编码节点初始化函数__init__
        输入：符号id
        功能：初始化符号编码Code
        """
        self.Id = id
        self.Code = np.array([], dtype=int)

（3）费诺编码接口函数Fano

费诺编码方法是先把当前符号集的概率降序排列，然后在某两个符号间划分成两部分，使两部分各自的概率和之差最小，不断分割，直到每一部分都只剩一个符号。容易注意到这是很典型的递归思想，因此这部分做一些预处理，将排序好之后的编码过程分离出下一个函数FanoCoding来单独进行递归处理。

先对输入信源概率矩阵按照概率降序排列，然后根据排序后的信源概率矩阵创建以费诺编码节点为元素的数组arr，将信源概率矩阵及arr输入费诺编码函数FanoCoding得到编码结果，输出arr。

def Fano(symbol_matrix):
    """
    费诺编码函数Fano:
    输入：Nx2 信源概率矩阵symbol_matrix N为符号个数 第一列为符号名称 第二列为符号概率
    输出：FanoNode数组arr,可通过arr[i].Code取出符号对应的编码
    功能：对输入概率矩阵中的符号进行费诺编码
    """
    len_symbol_matrix = np.size(symbol_matrix, 0)

    # 按概率降序排列
    idx = np.lexsort([-1 * symbol_matrix[:, 1]])
    sorted_symbol_matrix = symbol_matrix[idx, :]

    # 创建FanoNode类型的结点
    arr = np.array([])
    for i in range(0, len_symbol_matrix):
        m = FanoNode(sorted_symbol_matrix[i][0])
        arr = np.append(arr, m)

    # 递归编码
    FanoCoding(symbol_matrix, arr, 0)
    
    return SortById(arr)  # 按照ID升序排列

（4）费诺编码生成函数FanoCoding

下面的代码思想就是上文提到的，依次在传入的符号集（多数情况下不是完整的符号集，而是经过分割后的符号集）每两个符号之间进行分割，计算两部分各自的概率之和，在差值最小的位置就是这次分割的结果。分割点两侧一部分编码0，另一部分编码1，这就完成了本次递归的编码。再用同样的逻辑递归处理分割后的两部分符号，直到最后每个分割的符号集都只剩一个符号，编码结束。

具体代码如下：

def FanoCoding(symbol_matrix, arr, start):
    """
    费诺编码生成函数FanoCoding:
    输入：Nx2 信源概率矩阵symbol_matrix N为符号个数 第一列为符号名称 第二列为符号概率
        费诺节点数组arr
        初始下标start
    输出：无 可从arr[i].Code中取出第i个符号的编码
    功能：对输入概率矩阵中的符号进行费诺编码
    """
    len_symbol_matrix = np.size(symbol_matrix, 0)  # 当前符号集的符号个数

    # 按概率降序排列
    idx = np.lexsort([-1 * symbol_matrix[:, 1]])
    sorted_symbol_matrix = symbol_matrix[idx, :]

    if len_symbol_matrix > 1:  # 当前符号集中只有一个符号时不再需要分割
        difference_prev = 1
        id_split = 0
        # 分块
        for i in range(0, len_symbol_matrix + 1):
            sum1 = 0  # 前i个符号的概率和
            sum2 = 0  # 剩余符号的概率和

            for j in range(0, i):
                sum1 += sorted_symbol_matrix[j][1]
            for j in range(i, len_symbol_matrix):
                sum2 += sorted_symbol_matrix[j][1]

            # 本次分割两侧的概率差
            difference_cur = abs(sum1 - sum2)
            if difference_cur > difference_prev:
                id_split = i - 1
                break
            difference_prev = difference_cur

        # 编码
        for i in range(0, id_split):
            arr[i + start].Code = np.append(arr[i + start].Code, 0)  # 一部分编码0
        for i in range(id_split, len_symbol_matrix):
            arr[i + start].Code = np.append(arr[i + start].Code, 1)  # 另一部分编码1

        # 在分割点处，将当前符号集分为前后两部分
        next_process_matrix = np.split(sorted_symbol_matrix, [id_split])

        # 对前后两部分各自进行递归处理
        FanoCoding(next_process_matrix[0], arr, start)
        FanoCoding(next_process_matrix[1], arr, start + id_split)
    else:
        pass

（5）费诺编码测试函数FanoCodeDisplay

下面同样用一段代码对费诺编码进行测试，测试数据集为上图中的符号集和对应概率。

def FanoCodeDisplay(source_probability):
    """
    费诺编码结果展示函数
    输入：Nx2 信源概率矩阵source_probability N为符号个数 第一列为符号名称 第二列为符号概率
    输出：无
    功能：展示费诺香农编码结果
    """

    code_fano = Fano(source_probability)
    for i in range(0, np.size(source_probability, 0)):
        print(code_fano[i].Id, code_fano[i].Code)


print("===============FanoCode===============")
print("---------------SourceAfterCode---------------")
symbol_and_probability = np.array([[2, 0.19], [1, 0.20], [3, 0.18], [7, 0.01], [5, 0.15], [6, 0.10], [4, 0.17]])

FanoCodeDisplay(symbol_and_probability)
print("---------------SymbolAfterCode---------------")

编码输出结果如下图，可以看到与上图的编码结果完全一致。

在这里插入图片描述

3.香农编码

（1）原理

对给定符号集 $S$ 及其概率集 $P$ ，按照概率对符号集进行降序排列，每个符号对应的香农码码长 $l$ 满足:

$log_2 \frac{1}{p} \le l < log_2 \frac{1}{p} + 1$

计算第一个符号到每个符号的累加概率，将累加概率转化成2进制小数，取小数点后的前 $l$ 位即可得到香农编码结果。

（2）香农编码节点ShannonNode

建立类来存放符号的id及其编码。

class ShannonNode(object):
    """
    香农编码节点类
    """

    def __init__(self, id):
        """
        香农编码节点初始化函数__init__
        输入：符号id
        功能：初始化符号编码Code
        """
        self.Id = id
        self.Code = np.array([], dtype=int)

（3）香农编码接口函数Shannon

下面的部分仍然与前面的编码方法中这一部分类似，不同编码用相同的框架结构来实现，会非常方便，也便于出bug时进行调试。

def Shannon(symbol_matrix):
    """
    香农编码函数Shannon:
    输入：Nx2 信源概率矩阵symbol_matrix
        N为符号个数
        第一列为符号名称
        第二列为符号概率
    输出：ShannonNode数组arr,可通过arr[i].Code取出符号对应的编码
    功能：对输入概率矩阵中的符号进行香农编码
    """
    # 符号个数
    len_symbol_matrix = np.size(symbol_matrix, 0)

    # 按概率对符号集降序排列
    idx = np.lexsort([-1 * symbol_matrix[:, 1]])
    sorted_symbol_matrix = symbol_matrix[idx, :]

    # 将所有结点添加到arr数组内
    arr = np.array([])
    for i in range(0, len_symbol_matrix):
        m = ShannonNode(sorted_symbol_matrix[i][0])
        arr = np.append(arr, m)

    # 编码
    ShannonCoding(symbol_matrix, arr)
    
    return SortById(arr)

（4）香农编码生成函数ShannonCoding

下面的代码可以参考书上135页表5.11来看，matrix矩阵的各列分别代表符号、符号概率、累积概率和编码长度。依次计算累积概率，计算 $log_2p(x_i)$ 后向上取整即得到编码长度 $l$ ，对累积概率求二进制数，取前 $l$ 位即为符号的香农编码。

def ShannonCoding(symbol_matrix, arr):
    """
    香农码生成函数ShannonCoding:
    输入：Nx2 信源概率矩阵symbol_matrix
        N为符号个数
        第一列为符号名称
        第二列为符号概率
    输出：无 可从arr[i].Code中取出第i个符号的编码
    功能：对输入概率矩阵中的符号进行香农编码
    """
    # 符号个数
    len_symbol_matrix = np.size(symbol_matrix, 0)

    # 生成香农编码矩阵
    # 每一行代表一个符号的4项内容
    # 第一列是符号，第二列是该符号概率，第三列是累计概率，第四列是编码长度
    matrix = np.zeros((len_symbol_matrix, 4))
    for i in range(0, len_symbol_matrix):
        matrix[i][0] = symbol_matrix[i][0]
        matrix[i][1] = symbol_matrix[i][1]
        if i == 0:
            matrix[i][2] = 0
        elif i == 1:
            matrix[i][2] = matrix[i - 1][1]
        else:
            matrix[i][2] = matrix[i - 1][1] + matrix[i - 1][2]
        matrix[i][3] = math.ceil(-math.log2(matrix[i][1]))
    # 长度最长的编码
    lmax = int(max(matrix[:, 3]))

    # 把概率变为二进制码，取码长位即生成对应的香农编码
    for i in range(0, len_symbol_matrix):
        tmp = matrix[i][2]
        for j in range(0, lmax):
            if j < matrix[i][3]:
                tmp = tmp * 2
                arr[i].Code = np.append(arr[i].Code, math.floor(tmp))
                if tmp >= 1:
                    tmp -= 1

（5）香农编码测试函数ShannonCodeDisplay

用书上135页表5.11的数据来测试香农编码。

def ShannonCodeDisplay(source_probability):
    """
    香农编码结果展示函数
    输入：Nx2 信源概率矩阵source_probability N为符号个数 第一列为符号名称 第二列为符号概率
    输出：无
    功能：展示信源香农编码结果
    """

    code_shannon = Shannon(source_probability)
    for i in range(0, np.size(source_probability, 0)):
        print(code_shannon[i].Id, code_shannon[i].Code)

print("===============FanoCode===============")
print("---------------SourceAfterCode---------------")
symbol_and_probability = np.array([[2, 0.19], [1, 0.20], [3, 0.18], [7, 0.01], [5, 0.15], [6, 0.10], [4, 0.17]])

FanoCodeDisplay(symbol_and_probability)

编码结果如下，可以看到与书上的结果完全一致。

在这里插入图片描述

4.哈夫曼编码

哈夫曼编码需要用到二叉树相关的知识，如果对二叉树不熟悉的话，可以看这里：数据结构（四）：二叉树

下面的原理和代码可参考书136页例5.4.4的信源和表5.12、表5.13的哈夫曼编码理解，下面的代码和书上的可能有细微区别，我会在最后测试编码时说明。

（1）原理

对给定符号集 $S$ 及其概率集 $P$ ，按照概率对符号集进行降序排列，并创建二叉树节点数组，创建一个新的节点，将末尾两个符号节点分别挂在新节点的左右子节点，小的在左，大的在右，新节点的概率为左右子节点的概率之和。删除节点数组中的末尾两节点，将新节点加入节点数组，对节点数组按照节点的Value属性（节点的概率）降序排列，重复上述操作直至节点数组长度为1。以最后一个节点为根节点创建二叉树，即可得到哈夫曼树。对哈夫曼树遍历，设置临时数组存放哈夫曼编码，遍历左子节点就在临时数组末尾接0，遍历右子节点就在临时数组末尾接1，若当前节点无子节点，则将临时数组值赋给当前节点的Code属性，即可得到各个符号对应的哈夫曼编码。

（2）创建二叉树

设置每个二叉树节点类的成员变量分别为符号ID、节点概率、左子节点、右子节点以及编码。

class BinaryTreeNode(object):
    """
    二叉树节点类
    """

    def __init__(self, id, value, left, right):
        """
        初始化函数__init__：
           	输入：节点名称Id、节点概率值Value、左子节点Left、右子节点Right、编码结果Code
            功能：初始化二叉树节点
        """
        self.Id = id
        self.Value = value
        self.Left = left
        self.Right = right
        self.Code = np.array([], dtype=int)

建立二叉树类，成员变量为二叉树的根（类型为BinaryTreeNode），成员函数为哈夫曼编码输出，从根节点沿着树枝一路生成叶子节点的编码。

class BinaryTree(object):
    """
    二叉树类
    """

    def __init__(self, root):
        """
        初始化函数__init__：
            输入：BinaryTreeNode类的根节点root
            功能：初始化二叉树
        """
        self.Root = root
        
    def HuffmanOutput(self, root, arr):
        """
        哈夫曼编码函数HuffmanOutput：
            输入：BinaryTreeNode类的根节点root、父节点哈夫曼码数组arr
            功能：从根节点沿着二叉树生成各个叶子节点的哈夫曼码
        """
        if root:  # 当前节点不为空
            self.HuffmanOutput(root.Left, arr=np.append(arr, 0))
            self.HuffmanOutput(root.Right, arr=np.append(arr, 1, ))

            # 没有子节点
            if root.Left is None and root.Right is None:
                root.Code = arr

（3）哈夫曼编码接口函数Huffman

输入信源概率矩阵，使用哈夫曼树生成函数HuffmanCoding得到节点数组并返回，各节点变量的Code成员即为生成的编码（此函数从功能上可以省略，仅是为与其他编码方式格式保持统一）。

def Huffman(symbol_matrix):
    """
    哈夫曼编码接口函数Huffman:
    输入：Nx2 信源概率矩阵symbol_matrix
        N为符号个数
        第一行为符号名称
        第二行为符号概率
    输出：BinaryTreeNode数组arr,可通过arr[i].Code取出符号对应的编码
    功能：对输入概率矩阵中的符号进行香农编码
    """
    arr = HuffmanCoding(symbol_matrix)

    return SortById(arr)

（4）哈夫曼树生成函数HuffmanCoding

1、对输入信源概率矩阵按照概率降序排列；
2、创建节点数组；
3、创建新节点，取出节点数组末尾两节点分别接在新节点左右子节点上，将两节点概率和作为新节点的概率，删除节点数组末尾两节点，将新节点放入数组末尾，并重新对数组按照节点Value属性降序排列，重复上述步骤直至节点数组长度为1；
4、将最后一个节点作为根节点建立哈夫曼树；
5、根据哈夫曼树编码并返回上述节点数组；

def HuffmanCoding(symbol_matrix):
    """
    哈夫曼树生成函数HuffmanCoding:
    输入：Nx2 信源概率矩阵
        N为符号个数
        第一列为符号名称
        第二列为符号概率
    输出：编码后二叉树节点
        可通过node.Code取出对应节点的哈夫曼编码
    功能：对输入概率矩阵中的符号进行哈夫曼编码
    """
    huffman_output = np.array([], dtype=int)
    len_row = np.size(symbol_matrix, 0)
    id = 100
    i = 0
    node = np.array([])  # 包含全部节点的数组

    # 对符号概率矩阵按概率值降序排列
    idx = np.lexsort([-1 * symbol_matrix[:, 1]])
    sorted_symbol_matrix = symbol_matrix[idx, :]

    # 创建BinaryTreeNode数组
    while i < len_row:
        n = BinaryTreeNode(sorted_symbol_matrix[i][0], sorted_symbol_matrix[i][1], None, None)
        node = np.append(node, n)
        i += 1
    process_array = node.copy()

    # 创建哈夫曼树
    while np.size(process_array, 0) > 1:
        idx_last = np.size(process_array, 0)
        # 创建新节点m
        # m.Id = id
        # m.Value = process_array[idx_last - 1].Value + process_array[idx_last - 2].Value
        # m.Left = process_array[idx_last - 1]
        # m.Right = process_array[idx_last - 2]
        m = BinaryTreeNode(id, process_array[idx_last - 1].Value + process_array[idx_last - 2].Value,
                           process_array[idx_last - 1], process_array[idx_last - 2])

        id += 1
        # 删除process_array尾部两个元素
        process_array = np.delete(process_array, -1)
        process_array = np.delete(process_array, -1)

        # process_array尾部添加新生成的节点m
        process_array = np.append(process_array, m)

        # 对process_array按Value进行重新排序
        process_array = SortByValue(process_array)

    # 建立哈夫曼树
    tree = BinaryTree(process_array[0])
    # 根据哈夫曼树编码
    tree.HuffmanOutput(process_array[0], huffman_output)
    
    return node

（5）哈夫曼编码测试函数HuffmanCodeDisplay

下面测试哈夫曼编码部分，使用书上138页表5.13给出的样例。

def HuffmanCodeDisplay(source_probability):
    """
    哈夫曼编码结果展示函数
    输入：Nx2 信源概率矩阵source_probability N为符号个数 第一列为符号名称 第二列为符号概率
    输出：无
    功能：展示信源哈夫曼编码结果
    """

    code_huffman = Huffman(source_probability)
    for i in range(0, np.size(source_probability, 0)):
        print(code_huffman[i].Id, code_huffman[i].Code)

print("===============HuffmanCode===============")
print("---------------SourceAfterCode---------------")
symbol_and_probability = np.array([[2, 0.2], [1, 0.4], [5, 0.1], [3, 0.2], [4, 0.1]])

HuffmanCodeDisplay(symbol_and_probability)

编码结果如下，但结果与书上并不一致。

在这里插入图片描述

但用这里的编码结果和书上表5.12、5.13的编码结果具体计算一下各自的平均码长如下：

在这里插入图片描述
可以看到三种哈夫曼编码的实现方法最后的平均码长是相同的，下面说明这里实现的代码和书上两种编码方式的不同。

表5.12编码时概率较大的编码0，概率较小的编码1，具体如s2->s3时0.4概率编码为0，而0.2概率编码为1。上面的代码在编码时与这里刚好相反，概率较大的编码1，概率较小的编码0。
表5.13编码时新生成的节点概率如果和已有的相等，则放在最前面。而上面的代码对于新生成的概率结点，如果已有相同的概率，在对概率降序排列时并不能保证新生成的节点在最上面。

由上，尽管这三种编码结果并不相同，但都是合理的，最后的平均码长也是相同的。

5.排序函数

下面是两个排序函数，在上面的代码中出现，第一个是按照符号的id排升序，第二个是按照二叉树节点的概率值排降序。

def SortById(arr):
    """
    排序辅助函数SortById：
        输入：待排序数组
        输出：排序后数组
        功能：对输入类数组按Id属性降序排列
    """
    cmp_key = lambda p: p.Id
    return sorted(arr, key=cmp_key)


def SortByValue(arr):
    """
    排序辅助函数SortByValue：
        输入：待排序数组
        输出：排序后数组
        功能：对输入类数组按Value属性降序排列
    """
    cmp_key = lambda p: p.Value
    return sorted(arr, key=cmp_key, reverse=True)

6.信道编码接口函数

该部分代码是将信道编码的各个函数统一调用，SourceCodeMethod类内是各种编码的枚举，

class SourceCodeMethod(object):
    EQUAL_LENGTH = 0
    FANO         = 1
    SHANNON      = 2
    HUFFMAN      = 3
    
def Code(source_probability, symbol_sequence, method):
    """
    编码统一接口函数Code
    输入：Nx2 信源概率矩阵source_probability |N为符号个数|第一行为符号名称|第二行为符号概率
        符号序列symbol_sequence
        编码方式选择参数method |0:
    输出：符号与编码的字典、编码平均长度、0的概率、1的概率、编码输出Code
    功能：对输入符号序列按照信源概率进行编码
    """
    if method == 0:
        arr = EqualLength(source_probability)
    elif method == 1:
        arr = Fano(source_probability)
    elif method == 2:
        arr = Shannon(source_probability)
    else:
        arr = Huffman(source_probability)
    code_dict = {}
    sum0 = 0
    sum1 = 0
    len_average = 0
    for i in range(0, np.size(source_probability, 0)):
        code_dict.update({np.array2string(arr[i].Code).replace('[', '').replace(' ', '').replace(']', ''): str(arr[i].Id)})
        len_average += np.size(arr[i].Code, 0) * source_probability[i][1]
        sum0 += source_probability[i][1] * str(arr[i].Code).count("0")
        sum1 += source_probability[i][1] * str(arr[i].Code).count("1")
    p0 = sum0 / (sum1 + sum0)
    p1 = sum1 / (sum1 + sum0)
    code = np.array([], dtype=int)
    len_symbol = np.size(symbol_sequence, 0)
    for i in range(0, len_symbol):
        code = np.append(code, arr[int(symbol_sequence[i]) - 1].Code)
    return code_dict, len_average, p0, p1, code

三、信道编码

1.(n,k)线性码

输入序列为信源编码得到的0、1序列，设置一个生成矩阵 $G$ ， $G=[I_k|Q]$ ，将输入序列与 $Q$ 矩阵相乘记得校正子，将校正子拼接在输入序列右边即可得到长度为n的线性码。

def ChaEncodeNK(str_in):  # nk码编码程序
    array_out = str_in
    length = np.size(str_in, 0)  # 查看输入序列长度
    steps = 10 - length  # 计算需要补充的0个数
    # str用于与(n, k)矩阵相乘求校验码，为str_in补0至长度为10的序列
    str = np.zeros(10, dtype=int)
    for i in range(0, 10 - steps):
        str[i] = str_in[i]
    # 生成矩阵G=Ik|Q
    # Q矩阵
    matrix = np.array([
        [1, 0, 1, 1],
        [1, 1, 0, 1],
        [1, 1, 1, 0],
        [0, 1, 1, 1],
        [1, 0, 0, 1],
        [1, 0, 1, 0],
        [0, 1, 0, 1],
        [1, 1, 0, 0],
        [0, 1, 1, 0],
        [0, 0, 1, 1],
    ])
    str_check = np.dot(str, matrix)  # 点乘，计算校验位
    for i in range(0, np.size(str_check, 0)):
        str_check[i] = str_check[i] % 2  # 除2取余
    array_out = np.append(array_out, str_check)  # 在输出序列后追加校验码
    return array_out

2.简单重复编码

采用3次简单重复进行编码。

def RepeatEncode(str_in):  # 重复码编码程序
    length = np.size(str_in, 0)
    array_out = np.zeros(0, dtype=int)
    for i in range(0, length):
        if str_in[i] == 1:
            array_out = np.append(array_out, [1, 1, 1])  # 在输出序列后追加校验码
        else:
            array_out = np.append(array_out, [0, 0, 0])  # 在输出序列后追加校验码
    return array_out[:length * 3]

四、信道

1.单一信道

在这里插入图片描述

输入符号根据信道的输出矩阵直接生成输出符号即可，由于是单一信道，直接在一系列信道（channel）中任意选择一个当做这里的信道。

def SingleChannel(symbol_probability, channel, symbol_sequence):
    """
    单一信道编码
    :param symbol_sequence: 1xN 符号序列
    :param symbol_probability: 1x2 符号概率矩阵 分别表示符号0 1的概率
    :param channel: 1x2x2 信道矩阵
    :return: 信道输出符号序列、信道输出概率分布、等效信道概率矩阵
    """

    channel_matrix = channel[0]  # 由于是单一信道，取第一个信道
    channel_output = np.array([], dtype=int)  # 初始化符号输出序列

    # 确定符号输出序列
    symbol = np.array([0, 1])
    for i in range(np.size(symbol_sequence, 0)):
        if symbol_sequence[i] == 0:
            # 按照信道概率分布输出符号
            channel_output = np.append(channel_output, np.random.choice(symbol, size=1, p=channel_matrix[0]))
        else:
            channel_output = np.append(channel_output, np.random.choice(symbol, size=1, p=channel_matrix[1]))

    return channel_output, np.dot(symbol_probability, channel_matrix), channel_matrix

测试如下，输入序列为[1 0 1 1 0 0 1]，信道传输矩阵为

在这里插入图片描述

这里将噪声的影响体现为输入0时会有1%的概率输出1，而输入1时不受影响。误码率很小，输出与输入相同，0、1概率分布也不变。

在这里插入图片描述

2.级联信道

在这里插入图片描述

级联信道是将多个信道串联使用，总的传递矩阵为各信道传递矩阵相乘。

def CascadedChannel(symbol_probability, channel, symbol_sequence):
    """
    级联信道编码
    :param symbol_sequence: 1xN 符号序列
    :param symbol_probability: 1x2 符号概率矩阵 分别表示符号01的概率
    :param channel: 1x2x2 信道矩阵
    :return: 信道输出符号序列、信道输出概率分布、等效信道概率矩阵
    """
    channel_matrix = channel[0]  # 初始化级联信道矩阵
    n = np.size(channel, 0)  # 信道个数
    channel_output = np.array([], dtype=int)  # 初始化符号输出序列

    # 确定级联后信道矩阵
    for i in range(1, n):
        channel_matrix = np.dot(channel_matrix, channel[i])

    # 确定输出符号序列
    symbol = np.array([0, 1])
    for i in range(np.size(symbol_sequence, 0)):
        if symbol_sequence[i] == 0:
            channel_output = np.append(channel_output, np.random.choice(symbol, size=1, p=channel_matrix[0]))
        else:
            channel_output = np.append(channel_output, np.random.choice(symbol, size=1, p=channel_matrix[1]))
    return channel_output, np.dot(symbol_probability, channel_matrix), channel_matrix

测试如下，输入序列同上，信道要经过的信道可以为任意个，这里给出级联的三个信道依次为：

在这里插入图片描述

级联后经计算等效的信道如下图，第一个符号和第三个符号经过信道后从1变为0，出现了误码。

在这里插入图片描述

3.和信道

在这里插入图片描述

和信道是将多个信道组合后，在每一个时刻，以一定概率选择其中一个信道来传输。

def SumChannel(symbol_probability, channel, channel_probability, channel_id, symbol_sequence):
    """
    和信道
    :param symbol_sequence: 1xN 符号序列
    :param symbol_probability: 1x2 符号概率矩阵 分别表示符号01的概率
    :param channel: 1x2x2 信道矩阵
    :param channel_probability: 信道编号概率集合
    :param channel_id: 信道编号集合
    :return: 信道输出符号序列、信道输出概率分布、等效信道概率矩阵
    """

    channel_output = np.array([], dtype=int)  # 初始化符号输出序列
    id = np.random.choice(channel_id, size=1, p=channel_probability)  # 确定用哪个信道
    channel_matrix = channel[id].reshape((2, 2))  # 信道矩阵

    # 确定输出符号序列
    symbol = np.array([0, 1])
    for i in range(np.size(symbol_sequence, 0)):
        if symbol_sequence[i] == 0:
            channel_output = np.append(channel_output, np.random.choice(symbol, size=1, p=channel_matrix[0]))
        else:
            channel_output = np.append(channel_output, np.random.choice(symbol, size=1, p=channel_matrix[1]))

    return channel_output, np.dot(symbol_probability, channel_matrix), channel_matrix

信道输入同前，提供的信道同前，下面是两次的运行结果，可以看到两次选择的信道并不一样（也可能一样），但经过信道后的输出序列相同，输出的概率分布改变。

在这里插入图片描述

4.并联信道

在这里插入图片描述

多个信道并联，每个信道的输出仅与该信道的输入相关。

def ParallelChannel(symbol_probability, channel, symbol_sequence):
    """
    并联信道编码
    :param symbol_sequence: 1xN 符号序列
    :param symbol_probability: 1x2 符号概率矩阵 分别表示符号01的概率
    :param channel: Nx2x2 信道矩阵
    :return: 平均互信息、输入平均符号熵、输出平均符号熵、信道输出符号序列
    """
    n = np.size(channel, 0)  # 信道个数
    l_channel_matrix = int(math.pow(2, n))  # 并联等效信道矩阵长度
    channel_matrix = np.ones((l_channel_matrix, l_channel_matrix))  # 初始化并联等效信道信源矩阵
    input_probability_parallel = np.zeros(pow(2, n))
    tmp_id = np.array([])
    H_YX = 0
    H_total = 0
    H_input_avr = 0
    H_output_avr = 0

    # 得到并联输入概率矩阵
    for i in range(pow(2, n)):
        i_binary = bin(i).replace("0b", "").zfill(n)
        input_probability_parallel[i] = pow(symbol_probability[0], i_binary.count('0')) * pow(symbol_probability[1], i_binary.count('1'))
    # 得到信道概率矩阵
    for i in range(l_channel_matrix):
        i_binary = bin(i).replace("0b", "").zfill(n)
        for j in range(l_channel_matrix):
            j_binary = bin(j).replace("0b", "").zfill(n)
            for k in range(n):
                idx = int(i_binary[k])
                idy = int(j_binary[k])
                channel_matrix[i][j] *= channel[k][idx][idy]

    # 得到并联符号集和概率
    symbol = np.array([])
    symbol_probility_n = np.ones(l_channel_matrix)  # 并联后符号概率
    for i in range(l_channel_matrix):
        i_binary = bin(i).replace("0b", "").zfill(n)
        symbol = np.append(symbol, i_binary)
        for j in range(n):
            idx = int(i_binary[j])
            symbol_probility_n[i] *= symbol_probability[idx]

    # 根据信道矩阵得到输出
    channel_output = np.array([], dtype=int)
    for i in range(int(np.size(symbol_sequence) / n)):
        id = 0
        for j in range(n):
            id = id * 2 + symbol_sequence[i * n + j]
        tmp_id = np.append(tmp_id, id)

        y = np.random.choice(symbol, size=1, p=channel_matrix[id])  # 输出 为n长01字符串
        # 拼接所有输出
        for j in range(n):
            channel_output = np.append(channel_output, int(y[0][j]))
    output_probability_parallel = np.dot(symbol_probility_n, channel_matrix)
    for j in range((i + 1) * n, np.size(symbol_sequence)):
        channel_output = np.append(channel_output, symbol_sequence[j])
    # 并联熵计算
    # 损失熵
    for i in range(pow(2, n)):
        for j in range(pow(2, n)):
            H_YX -= input_probability_parallel[i] * channel_matrix[i][j] * np.log2(channel_matrix[i][j] + 1e-8)
    # 输入符号熵
    H_input_avr = -np.sum(input_probability_parallel * np.log2(input_probability_parallel + 1e-8)) / pow(2, n)
    # 输出符号熵
    for i in range(int(np.size(symbol_sequence) / n)):
        H_total += -output_probability_parallel[int(tmp_id[i])] * np.log2(output_probability_parallel[int(tmp_id[i])] + 1e-8)
    H_output_avr = H_total / int(np.size(symbol_sequence) / n)
    # HY
    H_Y = -np.sum(output_probability_parallel * np.log2(output_probability_parallel + 1e-8))
    # 平均互信息
    I = H_Y - H_YX
    return I, H_input_avr, H_output_avr, channel_output

输入和并联的三个信道矩阵同前，可以看到输出序列与输入相同，没有误码，输出为000（即三个信道均输出0）、001、……、111的概率如下图d数组所示。

在这里插入图片描述

5.信道接口函数

与信道编码接口函数相同，该部分是为了将函数的调用统一，并计算单一信道、和信道以及串联信道的熵（并联信道的熵在单独的代码内）。

class ChannelMethod(object):
    SINGLE = 0
    CASCADE = 1
    SUM = 2
    PARALLEL = 3


def CalculateChannelH(input_probability, output_probability, channel, symbol_sequence):
    """
    计算单一信道、和信道、串联信道的熵
    :param input_probability: 输入符号概率矩阵
    :param output_probability: 输出符号概率
    :param channel: 信道
    :param symbol_sequence: 符号序列
    :return: 平均互信息、输入符号熵、输出符号熵
    """

    H_YX = 0
    H_total = 0
    l_channel = np.size(channel, 0)
    l_symbol = np.size(symbol_sequence, 0)
    # 损失熵
    for i in range(l_channel):
        for j in range(l_channel):
            H_YX -= input_probability[i] * channel[i][j] * np.log2(channel[i][j] + 1e-8)
    # 输入符号熵
    H_input_avr = np.sum(-input_probability * np.log2(input_probability + 1e-8)) / l_symbol
    # HY
    H_Y = np.sum(-output_probability * np.log2(output_probability + 1e-8))
    # 平均互信息
    I = H_Y - H_YX
    # 输出符号熵
    for i in range(l_symbol):
        H_total += -output_probability[int(symbol_sequence[i])] * np.log2(output_probability[int(symbol_sequence[i])] + 1e-8)
    H_output_avr = H_total / l_symbol
    return I, H_input_avr, H_output_avr

6.信道部分功能和熵的测试

使用的三个信道为前文提到的三个信道，输入序列为[1 0 1 1 0 0]，运行查看各个熵的值以及输出序列。

a = np.array([[[0.99, 0.01], [0, 1]], [[0.95, 0.05], [0.1, 0.9]], [[1, 0], [0.04, 0.96]]])
p_channel = np.array([0.1, 0.2, 0.7])  # 选择信道的概率
id_channel = np.array([0, 1, 2])
symbol_prob = np.array([0.6, 0.4])  # 符号概率
symbol_seq = np.array([1, 0, 1, 1, 0, 0])  # 输入序列

print("==============SingleChannel==============")
I_single, H_input_avr_single, H_output_avr_single, symbol_out_single = Channel(symbol_prob, a, p_channel, id_channel, symbol_seq, ChannelMethod.SINGLE)
print("I_single = ", I_single, "\nH_input_avr_single = ", H_input_avr_single, "\nH_output_avr_single = ", H_output_avr_single)
print("symbol_out_single = ", symbol_out_single)

print("==============CascadedChannel==============")
I_cascaded, H_input_avr_cascaded, H_output_avr_cascaded, symbol_out_cascaded = Channel(symbol_prob, a, p_channel, id_channel, symbol_seq, ChannelMethod.CASCADE)
print("I_cascaded = ", I_cascaded, "\nH_input_avr_cascaded = ", H_input_avr_cascaded, "\nH_output_avr_cascaded = ", H_output_avr_cascaded)
print("symbol_out_cascaded = ", symbol_out_cascaded)

print("==============SumChannel==============")
I_sum, H_input_avr_sum, H_output_avr_sum, symbol_out_sum = Channel(symbol_prob, a, p_channel, id_channel, symbol_seq, ChannelMethod.SUM)
print("I_sum = ", I_sum, "\nH_input_avr_sum = ", H_input_avr_sum, "\nH_output_avr_sum = ", H_output_avr_sum)
print("symbol_out_sum = ", symbol_out_sum)

print("==============ParallelChannel==============")
I_parallel, H_input_avr_parallel, H_output_avr_parallel, symbol_out_parallel = Channel(symbol_prob, a, p_channel, id_channel, symbol_seq, ChannelMethod.PARALLEL)
print("I_parallel = ", I_parallel, "\nH_input_avr_parallel = ", H_input_avr_parallel, "\nH_output_avr_parallel = ", H_output_avr_parallel)
print("symbol_out_parallel = ", symbol_out_parallel)

运行结果如下：结果依次表示输入和输出的互信息量、输入符号的平均信息熵、输出符号的平均信息熵以及输出序列。

在这里插入图片描述

五、信道译码

1.(n,k)线性码

根据信道编码和信道传递矩阵计算出通过信道后的序列，与信道编码比较可获得错误图样，将错误图样矩阵与校验矩阵相乘即可得到校正子，用校正子可将通过信道引起的误码纠正。

# 校验子和错误图样的字典
NKDict = {
    "0000": 100,  # 全部正确
    "0001": 13,  # 第13位错误
    "0010": 12,
    "0100": 11,
    "1000": 10,
    "0011": 9,
    "0110": 8,
    "1100": 7,
    "0101": 6,
    "1010": 5,
    "1001": 4,
    "0111": 3,
    "1110": 2,
    "1101": 1,
    "1011": 0,
    "1111": 50,  # 错误过多
}

def ChaDecodeNK(str_in):  # nk码译码程序
    print('信道输出', str_in)
    matrix = np.array([
        [1, 0, 1, 1],
        [1, 1, 0, 1],
        [1, 1, 1, 0],
        [0, 1, 1, 1],
        [1, 0, 0, 1],
        [1, 0, 1, 0],
        [0, 1, 0, 1],
        [1, 1, 0, 0],
        [0, 1, 1, 0],
        [0, 0, 1, 1],
    ])
    matrix_a = np.array([
        [1, 0, 0, 0],
        [0, 1, 0, 0],
        [0, 0, 1, 0],
        [0, 0, 0, 1]
    ])
    matrix = np.transpose(matrix)
    matrix = np.concatenate((matrix, matrix_a), axis=1)  # 组合成H矩阵
    length = np.size(str_in, 0)  # 查看输入序列长度
    steps = 14 - length  # 计算需要补充的0个数
    # str为str_in在中间补0得到的长度为14的序列，用于与H矩阵相乘
    str = np.zeros(14, dtype=int)
    for i in range(0, 10 - steps):
        str[i] = str_in[i]
    for i in range(10 - steps, 10):
        str[i] = 0
    for i in range(10, 14):
        str[i] = str_in[i - steps]
    str_out = np.dot(matrix, np.transpose(str))  # 计算校验子
    str_check = np.transpose(str_out)  # 转置
    for i in range(0, np.size(str_check, 0)):
        str_check[i] = str_check[i] % 2
    str_check = np.array_str(str_check)  # 得到校验子
    str_check = str_check.replace(' ', '').replace('[', '').replace(']', '')

    change_pos = NKDict[str_check]  # 查字典中校验子位置，找到错误类型
    if change_pos != 100:
        print('error has been corrected!, 反转位置：', change_pos)
        str[change_pos] = 1 - str[change_pos]  # 将此位反转
    elif change_pos == 50:
        print('信道误码率过高，误码类型未查询')
    else:
        print('no error!')  # 表示译码正确
    print('str_out', str[:length - 4])  # 输出译码序列
    return str[:length - 4]

2.简单重复编码

按照最佳译码规则选出每个编码对应的译码符号输出。

RepeatDict = {
    "000": 0,
    "001": 0,
    "010": 0,
    "011": 1,
    "100": 0,
    "101": 1,
    "110": 1,
    "111": 1
}

def RepeatDecode(str_in):  # 重复码译码程序
    str_check = np.array_str(str_in)  # 将输入序列强制类型转换为str
    str_handle = str_check.replace(' ', '').replace('[', '').replace(']', '')  # 去符号
    array_out = np.zeros(0, dtype=int)
    for i in range(0, len(str_handle), 3):
        stra = str_handle[i] + str_handle[i + 1] + str_handle[i + 2]
        array_out = np.append(array_out, RepDict[stra])
    return array_out

六、信源译码

根据信源编码的对应关系直接进行译码。

def SourceDecode(str_in, _dict):  # 信源译码
    str_check = np.array_str(str_in)  # 将输入序列强制类型转换为str
    str_handle = str_check.replace(' ', '').replace('[', '').replace(']', '')  # 去符号
    print(str_handle)  # 得到字符串
    return _dict[str_handle]  # 根据字典译码

七、整合

将以上所有文件整合，在main内分别调用实现一次完整的信息传输过程，同时各部分的接口函数均有符号集概率分布，因此也可对应计算出对应的熵，如下所示。由于传输的每一个环节我们都设置了多个方法，下面只展示从信源开始到信宿结束众多组合中的一种。

import numpy as np

import Source # 信源
import Coding  # 信源编码
import Channel  # 信道
import Decoding  # 信道编码、信道解码以及信源译码

# 信源
n = 3
symbols = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20]
source1 = Source.MarkovInformationSource(symbols)
source2 = Source.StationaryInformationSource(symbols)
source3 = Source.ZeroMemoryInformationSource(symbols)

# 信源编码
source_code_method = Coding.SourceCodeMethod()  # 编码方式
source_code_dict = dict()  # 符号与编码字典
len_average = 0  # 平均编码长度
p0 = 0  # 0的概率
p1 = 0  # 1的概率

# 信道矩阵
ChannelMethod = channel.ChannelMethod()  # 编码方式

# 信道编码
ChannelCodeMethod = Decoding.ChannelCodeMethod()

channel_matrix = np.array([[[0.95, 0.05], [0.05, 0.95]], [[0.95, 0.05], [0.05, 0.95]], [[0.95, 0.05], [0.05, 0.95]]])
channel_probability = np.array([0.1, 0.2, 0.7])
channel_id = np.array([0, 1, 2])

channel_method = ChannelMethod.PARALLEL
channel_code_method = ChannelCodeMethod.REP
source_sequence = np.array([], dtype=int)  # 信源序列
source_code_sequence = np.array([], dtype=int)  # 信源编码结果
source_decode_sequence = np.array([], dtype=int)  # 信源译码结果
channel_decode_sequence = np.array([], dtype=int)  # 信道译码结果

for i in range(n):
    # 信源生成
    sign, probability_matrix = source1.Generate()
    source_sequence = np.append(source_sequence, sign)
    H_source_avr = Coding.CalculateSourceH(probability_matrix, sign)

    # 信源编码
    source_code_dict, len_average, p0, p1, source_code = Coding.Code(probability_matrix, sign,
                                                                     method=source_code_method.HUFFMAN)
    H_source_code_avr = Coding.CalculateH(p0, p1, source_code)

    # 信道编码
    if channel_method == ChannelCodeMethod.NK:
        channel_in = Decoding.ChaEncodeNK(source_code)
    else:
        channel_in = Decoding.RepeatEncode(source_code)

    # 信道
    p01 = np.array([p0, p1])
    if channel_method != ChannelMethod.PARALLEL:
        if channel_code_method != ChannelCodeMethod.REP:
            I, H_channel_input_avr, H_channel_output_avr, channel_out = channel.ChannelForNK(p01, channel_matrix,
                                                                                             channel_probability,
                                                                                             channel_id, channel_in,
                                                                                             channel_method)
        else:
            I, H_channel_input_avr, H_channel_output_avr, channel_out = channel.ChannelForRep(p01, channel_matrix,
                                                                                              channel_probability,
                                                                                              channel_id, channel_in,
                                                                                              channel_method)
    else:
        I, H_channel_input_avr, H_channel_output_avr, channel_out = channel.ParallelChannel(p01, channel_matrix,
                                                                                            channel_in)

    # 信道解码
    if channel_method == ChannelCodeMethod.NK:
        channel_decode = Decoding.ChaDecodeNK(channel_out)
    else:
        channel_decode = Decoding.RepeatDecode(channel_out)

    # 信源解码
    source_decode = Decoding.SourceDecode(channel_decode, source_code_dict)
    source_decode_sequence = np.append(source_decode_sequence, source_decode)

将测试代码放入UI展示，由于整个过程中可以观测的信息过多，下面只展示一小部分。

符号集为[1,20]的20个整数，用马尔可夫信源、霍夫曼编码、（n，k）码编码后进入信道，再经过解码得到最终结果。

在这里插入图片描述
由于马尔可夫信源在最开始时不平稳，经过一定时刻后平稳，所以UI展示了各时刻符号集每个符号的概率（数据过多，仅以时刻1为例）。

然后以最后一个时刻（时刻10）生成的符号为例，对其经过通信系统的整个流程进行分析。
在这里插入图片描述

通信系统经过10个时刻信源发出的序列、信源编码后的序列、信道译码后的序列以及最终收到的序列如下图。

在这里插入图片描述

信宿输出的符号-1代表经过误码后编码结果不在编码字典中（找不到对应的信源符号），因误码也可能导致某一符号被译为其它已存在的符号，但这里并没有发生。

以上即为我们做的所有内容。

感谢阅读，如有错误请批评指正

山舟

关注

5
点赞
踩
12

收藏

觉得还不错? 一键收藏
打赏
0
评论
2023BUAA信息论期末大作业（Python实现信源-编码-信道-译码-信宿通信过程）

以下内容是我和@Rg猿、@murmurto两位大佬共同完成的，所有内容仅是我们对作业要求和课程内容的理解，难免存在纰漏和错误，发布在这里仅供参考，切勿抄袭。
复制链接

扫一扫