Skip to article frontmatterSkip to article content

Example: HDD Decoding

Based on Proakis (2007, Example 7.5-1)

Note: Compare to Proakis (2007, Table 7.5-1), vary the last two rows

import numpy as np
import matplotlib.pyplot as plt
from scipy.spatial.distance import hamming

# Generator matrix G
G = np.array([
    [1, 0, 1, 0, 1],
    [0, 1, 0, 1, 1]
])
# Generate all possible messages and corresponding codewords for the given generator matrix G
def generate_message_codeword_pairs(G):
    # Number of message bits
    k = G.shape[0]
    # Generate all possible messages
    messages = [list(format(i, '0{}b'.format(k))) for i in range(2**k)]
    # Generate codewords by multiplying messages with G
    codeword_pairs = [(m, list(np.dot(np.array(m, dtype=int), G) % 2)) for m in messages]
    return codeword_pairs

# Get all message and corresponding codewords
message_codeword_pairs = generate_message_codeword_pairs(G)
message_codeword_pairs
[(['0', '0'], [0, 0, 0, 0, 0]), (['0', '1'], [0, 1, 0, 1, 1]), (['1', '0'], [1, 0, 1, 0, 1]), (['1', '1'], [1, 1, 1, 1, 0])]
# import numpy as np

# # Generator matrix G
# G = np.array([
#     [1, 0, 1, 0, 1],
#     [0, 1, 0, 1, 1]
# ])

# Function to calculate the Hamming distance between two codewords
def hamming_distance(codeword1, codeword2):
    # Calculate the Hamming distance
    return np.sum(codeword1 != codeword2)

# Generate all possible codewords from the generator matrix
def generate_codewords(G):
    # Number of message bits
    k = G.shape[0]
    # Generate all possible messages
    messages = [np.array(list(format(i, '0{}b'.format(k))), dtype=int) for i in range(2**k)]
    # Generate codewords by multiplying messages with G
    codewords = [np.dot(m, G) % 2 for m in messages]
    return codewords

# Function to find the minimum distance of the block code
def minimum_distance(G):
    # Generate all codewords
    codewords = generate_codewords(G)
    # Initialize the minimum distance to the maximum possible, which is the length of the codeword
    min_dist = G.shape[1]
    # Compare each pair of codewords to find the minimum Hamming distance
    for i in range(len(codewords)):
        for j in range(i + 1, len(codewords)):
            dist = hamming_distance(codewords[i], codewords[j])
            # Update minimum distance if a smaller one is found
            if dist < min_dist:
                min_dist = dist
    return min_dist

# Find the minimum distance
d_min = minimum_distance(G)
print('Minimum distance = ', d_min)
Minimum distance =  3

Standard Array


# Step 1: Generate all codewords
k = G.shape[0]  # Number of information bits
n = G.shape[1]  # Total number of bits in a codeword
# Generating all possible information messages u_m
information_messages = [np.array(list(np.binary_repr(i, width=k)), dtype=int) for i in range(2**k)]
# Creating corresponding codewords from u_m
codewords = [(message @ G) % 2 for message in information_messages]

# Step 2: Find the coset leaders and construct the standard array
# Initialize the standard array with the codewords as the first row
standard_array = [codewords]

# Create a list of all possible n-bit binary sequences
all_sequences = [np.array(list(np.binary_repr(i, width=n)), dtype=int) for i in range(2**n)]

# Define a function to calculate the weight of a binary sequence
def weight(seq):
    return np.sum(seq)

# Sort all sequences by their weight (except for the all-zero sequence which is a codeword)
sorted_sequences = sorted(all_sequences[1:], key=weight)

# Initialize the list of coset leaders with the all-zero sequence
coset_leaders = [np.zeros(n, dtype=int)]

# Find the minimum weight sequence that is not a codeword and add it to the coset leaders
for seq in sorted_sequences:
    if not any((seq == cw).all() for cw in codewords):
        coset_leaders.append(seq)
        # Construct the new coset and add it to the standard array
        new_coset = [(seq + cw) % 2 for cw in codewords]
        standard_array.append(new_coset)
        # Stop when we have enough coset leaders
        if len(coset_leaders) == 2**(n-k):
            break

# Convert the standard array to a NumPy array for better display
standard_array_np = np.array(standard_array, dtype=int)

# Print out the coset leaders
coset_leaders_np = np.array(coset_leaders, dtype=int)

standard_array_np, coset_leaders_np
(array([[[0, 0, 0, 0, 0], [0, 1, 0, 1, 1], [1, 0, 1, 0, 1], [1, 1, 1, 1, 0]], [[0, 0, 0, 0, 1], [0, 1, 0, 1, 0], [1, 0, 1, 0, 0], [1, 1, 1, 1, 1]], [[0, 0, 0, 1, 0], [0, 1, 0, 0, 1], [1, 0, 1, 1, 1], [1, 1, 1, 0, 0]], [[0, 0, 1, 0, 0], [0, 1, 1, 1, 1], [1, 0, 0, 0, 1], [1, 1, 0, 1, 0]], [[0, 1, 0, 0, 0], [0, 0, 0, 1, 1], [1, 1, 1, 0, 1], [1, 0, 1, 1, 0]], [[1, 0, 0, 0, 0], [1, 1, 0, 1, 1], [0, 0, 1, 0, 1], [0, 1, 1, 1, 0]], [[0, 0, 0, 1, 1], [0, 1, 0, 0, 0], [1, 0, 1, 1, 0], [1, 1, 1, 0, 1]], [[0, 0, 1, 0, 1], [0, 1, 1, 1, 0], [1, 0, 0, 0, 0], [1, 1, 0, 1, 1]]]), array([[0, 0, 0, 0, 0], [0, 0, 0, 0, 1], [0, 0, 0, 1, 0], [0, 0, 1, 0, 0], [0, 1, 0, 0, 0], [1, 0, 0, 0, 0], [0, 0, 0, 1, 1], [0, 0, 1, 0, 1]]))

Example 7.5-2


H = np.array([
    [1, 0, 1, 0, 0],
    [0, 1, 0, 1, 0],
    [1, 1, 0, 0, 1]
])


# Compute the syndrome for the coset leaders using the parity-check matrix H
syndrome_coset_leaders = (coset_leaders @ H.T) % 2

# syndrome_coset_leaders, coset_leaders

syndrome_coset_leader_pairs = [(syndrome, leader) for syndrome, leader in zip(syndrome_coset_leaders, coset_leaders)]

syndrome_coset_leader_pairs
[(array([0, 0, 0]), array([0, 0, 0, 0, 0])), (array([0, 0, 1]), array([0, 0, 0, 0, 1])), (array([0, 1, 0]), array([0, 0, 0, 1, 0])), (array([1, 0, 0]), array([0, 0, 1, 0, 0])), (array([0, 1, 1]), array([0, 1, 0, 0, 0])), (array([1, 0, 1]), array([1, 0, 0, 0, 0])), (array([0, 1, 1]), array([0, 0, 0, 1, 1])), (array([1, 0, 1]), array([0, 0, 1, 0, 1]))]
e_7_5_2 = np.array([1, 0, 1, 0, 0])
syndrome_7_5_2 = (e_7_5_2 @ H.T) % 2
syndrome_7_5_2
array([0, 0, 1])
# Find the corresponding syndrome in syndrome_coset_leaders and point out the corresponding coset leader
index = np.where(np.all(syndrome_coset_leaders == syndrome_7_5_2, axis=1))[0]
corresponding_coset_leader = coset_leaders[index[0]] if index.size > 0 else None

syndrome_7_5_2, corresponding_coset_leader
(array([0, 0, 1]), array([0, 0, 0, 0, 1]))
References
  1. Proakis, J. (2007). Digital Communications (5th ed.). McGraw-Hill Professional.