spyrit.misc.walsh_hadamard.fwht

spyrit.misc.walsh_hadamard.fwht(x, order=True)

Fast Walsh-Hadamard transform of x.

Due to the inherent numerical instability of the Hadamard transform (lots of additions and subtractions), it is recommended to use float64 whenever possible.

This is computed using Amit Portnoy’s algorithm available in the package hadamard-transform at https://github.com/amitport/hadamard-transform.

Args:

x (np.ndarray): batched input signal of shape \((* , n)\), where \(n\) is a power of two. The transform applies to the last dimension.

order (bool | list, optional): True for sequency (default), False for natural. It is also possible to provide a list of permutation indices.

Returns:

np.ndarray: transformed signal of shape \((* , n)\)

Example:

Example 1: Fast sequency-ordered (i.e., Walsh) Hadamard transform

>>> import numpy as np
>>> import spyrit.misc.walsh_hadamard as wh
>>> x = np.array([1, 3, 0, -1, 7, 5, 1, -2])
>>> y = wh.fwht(x)
>>> print(y)
[14 -8 -8 18 -4 -2 -6  4]

Example 2: Fast Hadamard transform

>>> import numpy as np
>>> import spyrit.misc.walsh_hadamard as wh
>>> x = np.array([1, 3, 0, -1, 7, 5, 1, -2])
>>> y = wh.fwht(x, False)
>>> print(y)
[14  4 18 -4 -8 -6 -8 -2]

Example 3: Permuted fast Hadamard transform

>>> import numpy as np
>>> import spyrit.misc.walsh_hadamard as wh
>>> x = np.array([1, 3, 0, -1, 7, 5, 1, -2])
>>> ind = [1, 0, 3, 2, 7, 4, 5, 6]
>>> y = wh.fwht(x, ind)
>>> print(y)
[ 4 14 -4 18 -2 -8 -6 -8]

Example 4: Comparison with Walsh-Hadamard transform via matrix-vector product

>>> from spyrit.misc.walsh_hadamard import fwht, walsh_matrix
>>> x = np.array([1, 3, 0, -1, 7, 5, 1, -2])
>>> y1 = fwht(x)
>>> H = walsh_matrix(8)
>>> y2 = H @ x
>>> print(f"Diff = {np.linalg.norm(y1-y2)}")
Diff = 0.0

Example 5: Comparison with the fast Walsh-Hadamard transform from sympy

>>> import spyrit.misc.walsh_hadamard as wh
>>> import sympy as sp
>>> x = np.array([1, 3, 0, -1, 7, 5, 1, -2])
>>> y1 = wh.fwht(x)
>>> y2 = sp.fwht(x)
>>> y3 = wh.sequency_perm(np.array(y2,dtype=x.dtype))
>>> print(y1)
[14 -8 -8 18 -4 -2 -6  4]
>>> print(f"Diff = {np.linalg.norm(y1-y3)}")
Diff = 0.0

Example 6: Computation times for 100 signals of length 2**12

>>> import timeit
>>> import numpy as np
>>> import spyrit.misc.walsh_hadamard as wh
>>> x = np.random.rand(100, 2**12)
>>> t = timeit.timeit(lambda: wh.fwht(x), number=100)
>>> print(f"Fast Walsh transform, no ind (100x): {t:.3f} seconds")
Fast Walsh transform, no ind (100x): ... seconds
>>> t = timeit.timeit(lambda: wh.fwht(x,False), number=100)
>>> print(f"Fast Hadamard transform (100x): {t:.3f} seconds")
Fast Hadamard transform (100x): ... seconds
>>> ind = wh.sequency_perm_ind(x.shape[-1])
>>> t = timeit.timeit(lambda: wh.fwht(x,ind), number=100)
>>> print(f"Fast Walsh transform, with ind (100x): {t:.3f} seconds")
Fast Walsh transform, with ind (100x): ... seconds

>>> import sympy as sp
>>> t = timeit.timeit(lambda: sp.fwht(x[0,:]), number=10)
>>> print(f"Fast Hadamard transform from sympy (only 1 signal, 10x): {t:.3f} seconds")
Fast Hadamard transform from sympy (only 1 signal, 10x): ... seconds