The parallel rays of a direction¶

Each cell's ray and its step along it, a shear of its coordinates.

Requires

the lattice's extents (cells, its axes in order), the Direction v the rays run (1, a nonzero vector of the lattice's dimension, the way the light travels; any length).

Provides

each cell's ray (an integer, shared by the cells of one ray and no other) and its step along it (an integer, growing the way v runs): the keys of ordered(it.ray, it.step).

Law

m the major axis, the largest |v_m| (the first on a tie); the step is c_m sign(v_m); on every other axis a the cell's offset from the line through the origin, s_a = c_a - floor(c_m v_a / v_m + 1/2), the ray the offsets numbered (s_a + n_m in base n_a + 2 n_m + 1, the axes in order). The rays are digital straight lines of slope v_a / v_m, each crossing the major axis one cell a step, so they partition the lattice: every cell on exactly one ray, at exactly one step.

Source

Bresenham, Algorithm for computer control of a digital plotter, IBM Systems Journal 4(1), 1965; Reveilles, Geometrie discrete, calcul en nombres entiers et algorithmique (thesis, 1991), the arithmetic line.

Proof

invariant, limit.

No Trait form: waits on the middle: a key kernel, each cell's ray and step the keys a directional fold walks, which draft 5 does not write yet; until then a reference.

rays.parallel(shape, direction)

Each cell's ray and its step along it, a shear of its coordinates: the keys a directional law folds over, ordered(it.ray, it.step) with @ before.

In [1]:

import os, sys
root = os.getcwd()
while not os.path.isfile(os.path.join(root, 'docsite.json')):
    root = os.path.dirname(root)
sys.path.insert(0, os.path.join(root, 'tools'))
import numpy as np
import matplotlib.pyplot as plt
from lusine import library
from lusine.kernels._rays import ordered
from lusine.kernels.rays.parallel import parallel

import os, sys root = os.getcwd() while not os.path.isfile(os.path.join(root, 'docsite.json')): root = os.path.dirname(root) sys.path.insert(0, os.path.join(root, 'tools')) import numpy as np import matplotlib.pyplot as plt from lusine import library from lusine.kernels._rays import ordered from lusine.kernels.rays.parallel import parallel

The reference¶

The kernel's reference, in numpy; its Trait form is not written yet.

In [2]:

library.reference(parallel)

library.reference(parallel)

Out[2]:

the reference, parallel.py

def parallel(shape, direction):
    """THE PARALLEL RAYS OF A DIRECTION: each cell's ray and its step along
    it, a shear of its coordinates.

    Requires: the lattice's extents (cells, its axes in order), the
    `Direction` v the rays run (1, a nonzero vector of the lattice's
    dimension, the way the light travels; any length).
    Provides: each cell's `ray` (an integer, shared by the cells of one ray
    and no other) and its `step` along it (an integer, growing the way v
    runs): the keys of `ordered(it.ray, it.step)`.
    Law: m the major axis, the largest |v_m| (the first on a tie); the step
    is c_m sign(v_m); on every other axis a the cell's offset from the line
    through the origin, s_a = c_a - floor(c_m v_a / v_m + 1/2), the ray the
    offsets numbered (s_a + n_m in base n_a + 2 n_m + 1, the axes in order).
    The rays are digital straight lines of slope v_a / v_m, each crossing
    the major axis one cell a step, so they partition the lattice: every
    cell on exactly one ray, at exactly one step.
    Source: Bresenham, Algorithm for computer control of a digital plotter,
    IBM Systems Journal 4(1), 1965; Reveilles, Geometrie discrete, calcul en
    nombres entiers et algorithmique (thesis, 1991), the arithmetic line.
    Proof: invariant, limit."""
    shape = tuple(int(n) for n in shape)
    v, m = _major(direction, len(shape))
    idx = np.indices(shape, dtype=np.int64)
    cm = idx[m]
    step = cm * (1 if v[m] > 0 else -1)
    ray = np.zeros(shape, dtype=np.int64)
    nm = shape[m]
    for a in range(len(shape)):
        if a == m:
            continue
        s = idx[a] - np.floor(cm * (v[a] / v[m]) + 0.5).astype(np.int64)
        ray = ray * (shape[a] + 2 * nm + 1) + (s + nm)
    return ray, step

Invariant: the rays partition the lattice¶

Every cell lies on exactly one ray at exactly one step, whatever the direction, on a plane or in space; a ray's cells are a digital straight line, one cell a step along the major axis.

In [3]:

for shape, d in (((24, 31), (1, 0.37)), ((20, 20), (1, -1)), ((15, 40), (0.2, 1)), ((6, 7, 8), (1, 0.5, -0.25))):
    ray, step = parallel(shape, d)
    assert len(set(zip(ray.ravel().tolist(), step.ravel().tolist()))) == ray.size
ray, step = parallel((24, 31), (1, 0.37))
fig, ax = plt.subplots(1, 2, figsize=(8, 3.6))
ax[0].imshow(ray % 7, cmap='tab10'); ax[0].set_title('ray (coloured mod 7)')
ax[1].imshow(step, cmap='viridis'); ax[1].set_title('step along it')
plt.show()
print('one ray and one step a cell')

for shape, d in (((24, 31), (1, 0.37)), ((20, 20), (1, -1)), ((15, 40), (0.2, 1)), ((6, 7, 8), (1, 0.5, -0.25))): ray, step = parallel(shape, d) assert len(set(zip(ray.ravel().tolist(), step.ravel().tolist()))) == ray.size ray, step = parallel((24, 31), (1, 0.37)) fig, ax = plt.subplots(1, 2, figsize=(8, 3.6)) ax[0].imshow(ray % 7, cmap='tab10'); ax[0].set_title('ray (coloured mod 7)') ax[1].imshow(step, cmap='viridis'); ax[1].set_title('step along it') plt.show() print('one ray and one step a cell')

one ray and one step a cell

Invariant: the walk is the keys', not the storage's¶

The order ordered(it.ray, it.step) walks is fixed by the keys alone: the cells shuffled in storage are walked in the same order, so a float fold over it has the same bits however the lattice is laid out.

In [4]:

ray, step = parallel((30, 30), (1, 0.6))
order = ordered(ray, step)
for seed in range(3):
    perm = np.random.default_rng(seed).permutation(ray.size)
    assert np.array_equal(perm[ordered(ray.ravel()[perm], step.ravel()[perm])], order)
print('the same walk whatever the storage order')

ray, step = parallel((30, 30), (1, 0.6)) order = ordered(ray, step) for seed in range(3): perm = np.random.default_rng(seed).permutation(ray.size) assert np.array_equal(perm[ordered(ray.ravel()[perm], step.ravel()[perm])], order) print('the same walk whatever the storage order')

the same walk whatever the storage order

Limit: along an axis, along a diagonal¶

A beam down an axis has the columns for its rays and the rows for its steps; a diagonal one the diagonals.

In [5]:

ray, step = parallel((5, 4), (1, 0))
assert np.array_equal(step, np.indices((5, 4))[0]) and all(len(set(ray[:, x])) == 1 for x in range(4))
ray, _ = parallel((12, 12), (1, 1))
assert all(len(set(y - x for y, x in zip(*np.nonzero(ray == q)))) == 1 for q in np.unique(ray))
print('columns, and diagonals')

ray, step = parallel((5, 4), (1, 0)) assert np.array_equal(step, np.indices((5, 4))[0]) and all(len(set(ray[:, x])) == 1 for x in range(4)) ray, _ = parallel((12, 12), (1, 1)) assert all(len(set(y - x for y, x in zip(*np.nonzero(ray == q)))) == 1 for q in np.unique(ray)) print('columns, and diagonals')

columns, and diagonals