import numpy as np
from PIL import Image

def get_matrix():
    img = Image.open("output2g.png").convert("L")
    matrix = (np.array(img) < 128).astype(np.uint8)
    return matrix

def check_rank(mat):
    # Gaussian elimination over GF(2)
    m, n = mat.shape
    A = mat.copy()
    rank = 0
    for j in range(n):
        if rank >= m: break
        # Find pivot
        pivot = -1
        for i in range(rank, m):
            if A[i, j] == 1:
                pivot = i
                break
        if pivot != -1:
            # Swap rows
            A[[rank, pivot]] = A[[pivot, rank]]
            # Eliminate
            for i in range(rank + 1, m):
                if A[i, j] == 1:
                    A[i] ^= A[rank]
            rank += 1
    return rank

def main():
    P = get_matrix() # 16384 x 2048
    print(f"P shape: {P.shape}")
    
    # Check rank of P
    # P has 2048 columns. Rank is at most 2048.
    # Since H = [P^T | I], the rank of H is exactly 2048 (because of I).
    # The correction capability depends on d_min.
    
    # Let's check some properties of the columns of P.
    # Each column of P is a syndrome for a data bit error.
    # Column j of P corresponds to data bit j.
    # The syndrome for data bit j error is P[:, j].
    # Wait, if G = [I | P], then codeword c = [d | dP] = [d | p].
    # Syndrome s = H c^T = [P^T | I] [d | p]^T = P^T d^T + p^T.
    # For a codeword, s = 0.
    # If error e = [e_d | e_p], syndrome s = P^T e_d^T + e_p^T.
    # If a single data bit j is in error, e_d = bit j set, e_p = 0.
    # s = P^T [0...1...0]^T = j-th column of P^T = j-th row of P.
    # Wait. If s = H e^T, then s is a 2048-bit vector.
    # P^T is 2048 x 16384. Its columns (length 2048) are the syndromes of data bits.
    # So the ROWS of P are the syndromes of the 16384 data bits.
    # The columns of I (length 2048) are the syndromes of the 2048 parity bits.
    
    # Check if all rows of P are unique and non-zero
    rows_p = ["".join(map(str, row)) for row in P]
    unique_rows = set(rows_p)
    print(f"Unique rows in P: {len(unique_rows)} out of {len(rows_p)}")
    
    # If all rows are unique and non-zero, it can correct at least 1 error.
    # What about d_min?
    # If any two rows are the same, d_min=2 (two data bit errors give 0 syndrome).
    # If any row is 0, d_min=2 (one data bit error and one parity bit error give 0 syndrome).
    # Check for zero rows
    zero_rows = [i for i, row in enumerate(rows_p) if "1" not in row]
    print(f"Zero rows in P: {len(zero_rows)}")

    # Check weights
    weights = [row.count("1") for row in rows_p]
    print(f"Row weights: min={min(weights)}, max={max(weights)}, avg={sum(weights)/len(weights)}")

    # Check XORs of 2 rows (t=2 capability)
    # This is hard to do for all 16384^2.
    # But we can check a few.
    
if __name__ == "__main__":
    main()
