Sample code for the MaxwellRules post "An unexpected matrix trick" showing how to reduce the dimension of a homography given a fixed coordinate along a dimension.

homography_reduction.py

import numpy as np

DIM = 4  # Dimension of matrix which represents the homography
IGNORE_DIM = 2  # This would be t on the blog post example

# We can't remove the last dim of the homography because it's
# special.
assert 0 <= IGNORE_DIM <= (DIM - 2)

assert IGNORE_DIM

rng = np.random.default_rng(42)


def normalize_homography(h: np.ndarray) -> np.ndarray:
    return h / h[-1, -1]


def normalize_vec(vec: np.ndarray) -> np.ndarray:
    return vec / vec[-1]


def reduce_h_dim(
    h: np.ndarray, ignored_component_value: float, ignored_dim: int
) -> np.ndarray:
    # ASSUME: mat is a homography acting on a homogenous coordinates vector

    new_col = ignored_component_value * np.delete(h[:, ignored_dim], ignored_dim)

    reduced_h = np.delete(np.delete(h, ignored_dim, axis=1), ignored_dim, axis=0)
    reduced_h[:, -1] = reduced_h[:, -1] + new_col

    return reduced_h


def main(rng):
    H = normalize_homography(rng.normal(3, 2, (DIM, DIM)))
    vec = normalize_vec(rng.normal(2, 1, (DIM,)))

    # Eliminate the dimension we didn't care about and
    # normalize the homogenous vector
    res = np.delete(normalize_vec(H @ vec), IGNORE_DIM)

    reduced_H = reduce_h_dim(H, vec[IGNORE_DIM], IGNORE_DIM)
    reduced_vec = np.delete(vec, IGNORE_DIM)

    # No need to eliminate the dimension we didn't care about in the
    # result as it was taken care of by reduced_H and reduced_ved
    # then normalize the homogenous vector
    reduced_res = normalize_vec(reduced_H @ reduced_vec)

    assert np.allclose(reduced_res, res)


if __name__ == "__main__":
    # Run a few times to check we weren't just lucky
    for _ in range(2**8):
        main(rng)