Detectors

The ComNumpy Library implements several classical MIMO detectors.

class MaximumLikelihoodDetector(alphabet: ndarray, *, H: ndarray | None = None, name: str = 'ML Detector')

Maximum-likelihood (ML) MIMO detector for white Gaussian noise.

Signal Model

The detector assumes the flat MIMO observation model

\[\mathbf{y}[n] = \mathbf{H} \, \mathbf{x}[n] + \mathbf{b}[n]\]

where \(\mathbf{H}\) is the channel matrix of size \(N_r \times N_t\), \(\mathbf{x}[n] \in \mathcal{M}^{N_t}\) the transmitted vector drawn from the constellation \(\mathcal{M}\) and \(\mathbf{b}[n]\) a circular white Gaussian noise. The ML decision under this model is the exhaustive search

\[\widehat{\mathbf{x}}[n] = \arg \min_{\mathbf{x} \in \mathcal{M}^{N_t}} \|\mathbf{y}[n] - \mathbf{H}\mathbf{x}\|^2\]

and the detector returns the integer indices of \(\widehat{\mathbf{x}}[n]\) in the alphabet (\(\widehat{\mathbf{x}}[n]\) = alphabet[S_[:, n]]).

Warning

The exhaustive search evaluates \(|\mathcal{M}|^{N_t}\) candidates and becomes expensive for a large number of transmit antennas or a high-order constellation.

Axes: declared axis – expects (ant, N) with antennas on axis -2. H may be a stack (K, N_r, N_t) against Y of shape (K, N_r, N): channel k decides frame k, through an internal per-draw loop – convenience, not a speedup (D51).

param alphabet:

Symbol constellation \(\mathcal{M}\) (1D array).

type alphabet:

np.ndarray

param H:

Channel matrix \(\mathbf{H}\) of size \(N_r \times N_t\). Must be set before calling the detector.

type H:

np.ndarray, optional, keyword-only

param name:

Name of the detector. Default is "ML Detector".

type name:

str, optional, keyword-only

S_

Decisions of the last call, as alphabet indices (data-dependent, hence the trailing underscore, D23).

Type:

np.ndarray

raises ValueError:

If H is not set when the detector is called.

References

• E. G. Larsson, P. Stoica and G. Ganesan, Space-Time Block Coding for Wireless Communications, Cambridge University Press, 2003.

• D. Tse and P. Viswanath, Fundamentals of Wireless Communication, Cambridge University Press, 2005, Chapter 8.

Examples

>>> alphabet = np.array([-1.0 + 0j, 1.0 + 0j])
>>> H = np.eye(2)
>>> Y = np.array([[0.9, -1.1], [-0.8, 1.2]])
>>> MaximumLikelihoodDetector(alphabet, H=H)(Y)
array([[1, 0],
       [0, 1]])

get_candidates(alphabet: ndarray, N_t: int) → tuple[ndarray, ndarray]

Every vector of \(\mathcal{M}^{N_t}\), as indices and symbols.

The enumeration is mixed-radix counting – the order itertools.product(range(M), repeat=N_t) produces – written as arithmetic on an index grid, because the loop that built it one candidate at a time was the dominant cost for a large alphabet.

class SphereDecoder(alphabet: ndarray, *, H: ndarray | None = None, name: str = 'Sphere Decoder')

Maximum likelihood by tree search instead of exhaustion.

Signal Model

The decision is the one of MaximumLikelihoodDetector – exactly, not approximately – reached without evaluating every candidate. Write the thin QR factorization of the channel, \(\mathbf{H} = \mathbf{Q}\mathbf{R}\) with \(\mathbf{Q}^H\mathbf{Q} = \mathbf{I}_{N_t}\) and \(\mathbf{R}\) upper triangular, and project the observation on \(\mathbf{z} = \mathbf{Q}^H \mathbf{y}\):

\[\left\|\mathbf{y} - \mathbf{H}\mathbf{x}\right\|^2 = \underbrace{\left\|\mathbf{z} - \mathbf{R}\mathbf{x}\right\|^2} _{\text{depends on } \mathbf{x}} + \left\|\left(\mathbf{I} - \mathbf{Q}\mathbf{Q}^H\right) \mathbf{y}\right\|^2\]

The second term is a constant, so minimizing the first is maximum likelihood. Its value is what makes the search possible, because \(\mathbf{R}\) is triangular:

\[\left\|\mathbf{z} - \mathbf{R}\mathbf{x}\right\|^2 = \sum_{k=N_t-1}^{0} \left|R_{kk}\right|^2 \left|c_k(x_{k+1}, \ldots, x_{N_t-1}) - x_k\right|^2, \qquad c_k = \frac{z_k - \sum_{i>k} R_{ki} x_i}{R_{kk}}\]

Every term is non-negative and layer \(k\) depends only on the layers already decided, so a partial sum can only grow: the moment it exceeds the best full metric found so far, the whole subtree below it can be discarded. That is the sphere: the search visits only the lattice points inside a ball whose radius shrinks each time a better solution is found.

The enumeration order within a layer is the Schnorr-Euchner one – alphabet points by increasing \(|c_k - a|\) – for two reasons. The first leaf reached is then the successive-interference-cancelling (Babai) solution, which gives a finite radius immediately, so no initial radius has to be guessed; and since the terms are visited in increasing order, the first one that exceeds the bound ends the layer, no continue needed.

What it costs. The visited-node count is data dependent – that is the whole point – and is reported as nodes_, the average per detected vector. At high SNR it approaches \(N_t\) (the Babai point is already ML and nothing else survives the first bound); at low SNR, or on an ill-conditioned channel, it grows towards the exhaustive \(|\mathcal{M}|^{N_t}\). In Python the per-node overhead is large, so this block is slower than the vectorized exhaustive search on small problems and becomes the only option when \(|\mathcal{M}|^{N_t}\) stops fitting in memory – 16-QAM on four antennas is 65 536 candidates per symbol, 64-QAM on four is 16.7 million.

Axes: declared axis – expects (ant, N) with antennas on axis -2. H may be a stack (K, N_r, N_t) against Y of shape (K, N_r, N): channel k decides frame k, through an internal per-draw loop – convenience, not a speedup (D51) – and returns the alphabet indices, like every other detector here.

param alphabet:

Symbol constellation \(\mathcal{M}\) (1D array).

type alphabet:

np.ndarray

param H:

Channel matrix \(\mathbf{H}\) of size \(N_r \times N_t\), with \(N_r \geq N_t\). Must be set before the call.

type H:

np.ndarray, optional, keyword-only

param name:

Name of the detector. Default is "Sphere Decoder".

type name:

str, optional, keyword-only

nodes_

Average number of tree nodes visited per detected vector in the last call – data-dependent, hence the trailing underscore (D23). Compare it with \(|\mathcal{M}|^{N_t}\).

Type:

float

S_

Decisions of the last call, as alphabet indices (D23).

Type:

np.ndarray

raises ValueError:

If H is not set, if the channel has fewer receive than transmit antennas (the decomposition needs a full-column-rank \(\mathbf{R}\)), or if that rank is deficient.

References

• U. Fincke, M. Pohst, “Improved methods for calculating vectors of short length in a lattice”, Mathematics of Computation, vol. 44, no. 170, pp. 463-471, 1985.

• C. P. Schnorr, M. Euchner, “Lattice basis reduction: improved practical algorithms and solving subset sum problems”, Mathematical Programming, vol. 66, pp. 181-199, 1994.

• E. Agrell, T. Eriksson, A. Vardy, K. Zeger, “Closest point search in lattices”, IEEE Trans. Inf. Theory, vol. 48, no. 8, pp. 2201-2214, 2002.

• B. Hassibi, H. Vikalo, “On the sphere-decoding algorithm I. Expected complexity”, IEEE Trans. Signal Process., vol. 53, no. 8, pp. 2806-2818, 2005.

Examples

>>> alphabet = np.array([-1.0 + 0j, 1.0 + 0j])
>>> H = np.eye(2)
>>> Y = np.array([[0.9, -1.1], [-0.8, 1.2]])
>>> SphereDecoder(alphabet, H=H)(Y)
array([[1, 0],
       [0, 1]])

The decision is the exhaustive one, on far fewer candidates:

>>> from comnumpy.core.utils import get_alphabet
>>> from comnumpy.mimo.utils import rayleigh_channel
>>> alphabet = get_alphabet("QAM", 16)
>>> H = rayleigh_channel(4, 4, seed=0)
>>> rng = np.random.default_rng(1)
>>> sent = rng.integers(0, 16, size=(4, 200))
>>> Y = H @ alphabet[sent] + 0.05 * (rng.standard_normal((4, 200))
...                                  + 1j * rng.standard_normal((4, 200)))
>>> decoder = SphereDecoder(alphabet, H=H)
>>> exhaustive = MaximumLikelihoodDetector(alphabet, H=H)
>>> bool(np.array_equal(decoder(Y), exhaustive(Y)))
True
>>> print(f"{decoder.nodes_:.1f} nodes visited, against {16 ** 4}")
4.5 nodes visited, against 65536

class LinearDetector(alphabet: ndarray, *, H: ndarray | None = None, method: Literal['zf', 'mmse'] = 'zf', sigma2: float | None = None, name: str = 'ZF Detector')

Linear MIMO detector (equalization followed by hard decision).

Signal Model

The detector assumes the flat MIMO observation model

\[\mathbf{y}[n] = \mathbf{H} \, \mathbf{x}[n] + \mathbf{b}[n], \qquad \mathbf{b}[n] \sim \mathcal{CN}\left(\mathbf{0}, \sigma^2 \mathbf{I}_{N_r}\right)\]

where \(\mathbf{H}\) is the channel matrix of size \(N_r \times N_t\). Detection proceeds in two steps: linear equalization

\[\mathbf{z}[n] = \mathbf{W} \mathbf{y}[n], \qquad \mathbf{W}_{zf} = \mathbf{H}^{\dagger}, \qquad \mathbf{W}_{mmse} = \left(\mathbf{H}^H \mathbf{H} + \sigma^2 \mathbf{I}_{N_t}\right)^{-1} \mathbf{H}^H\]

then per-antenna projection of \(\mathbf{z}[n]\) onto the nearest point of the constellation \(\mathcal{M}\). The detector returns the integer indices of the decisions in the alphabet.

Axes: declared axis – expects (ant, N) with antennas on axis -2. H may be a stack (K, N_r, N_t) against Y of shape (K, N_r, N) – one draw per trial, decided through numpy’s stacked linear algebra in one call (D51).

param alphabet:

Symbol constellation \(\mathcal{M}\) (1D array).

type alphabet:

np.ndarray

param H:

Channel matrix \(\mathbf{H}\) of size \(N_r \times N_t\). Must be set before calling the detector.

type H:

np.ndarray, optional, keyword-only

param method:

Linear equalizer \(\mathbf{W}\). Default is "zf".

type method:

{“zf”, “mmse”}, keyword-only

param sigma2:

Noise variance \(\sigma^2\) (required by the MMSE equalizer only).

type sigma2:

float, optional, keyword-only

param name:

Name of the detector. Default is "ZF Detector".

type name:

str, optional, keyword-only

raises ValueError:

If H is not set when the detector is called.

References

• E. G. Larsson, P. Stoica and G. Ganesan, Space-Time Block Coding for Wireless Communications, Cambridge University Press, 2003.

• D. Tse and P. Viswanath, Fundamentals of Wireless Communication, Cambridge University Press, 2005, Chapter 8.

Examples

>>> alphabet = np.array([-1.0 + 0j, 1.0 + 0j])
>>> Y = np.array([[0.9, -1.1], [-0.8, 1.2]])
>>> LinearDetector(alphabet, H=np.eye(2))(Y)
array([[1, 0],
       [0, 1]])

linear_estimator(Y: ndarray) → ndarray

Perform Zero Forcing or MMSE linear equalization

class OrderedSuccessiveInterferenceCancellationDetector(alphabet: ndarray, *, osic_type: str = 'sinr', H: ndarray | None = None, method: Literal['zf', 'mmse'] = 'zf', sigma2: float | None = None, name: str = 'OSIC Detector')

Ordered Successive Interference Cancellation (OSIC) MIMO detector.

Signal Model

The detector assumes the flat MIMO observation model

\[\mathbf{y}[n] = \mathbf{H} \, \mathbf{x}[n] + \mathbf{b}[n], \qquad \mathbf{b}[n] \sim \mathcal{CN}\left(\mathbf{0}, \sigma^2 \mathbf{I}_{N_r}\right)\]

where \(\mathbf{H}\) is the channel matrix of size \(N_r \times N_t\) and \(\mathbf{x}[n] \in \mathcal{M}^{N_t}\). At each of the \(N_t\) stages, the stream selected by the ordering rule is equalized (ZF or MMSE), hard-detected on the constellation \(\mathcal{M}\), and its estimated contribution is cancelled:

\[\mathbf{y}[n] \leftarrow \mathbf{y}[n] - \mathbf{h}_{k} \, \widehat{x}_{k}[n]\]

where \(\mathbf{h}_{k}\) is the column of \(\mathbf{H}\) associated with the detected stream \(k\). The ordering rule is "sinr" (max post-equalization SINR, MMSE), "colnorm" (max column norm of \(\mathbf{H}\), ZF) or "snr" (max post-equalization SNR, ZF). The detector returns the integer indices of the decisions in the alphabet.

Axes: declared axis – expects (ant, N) with antennas on axis -2. H may be a stack (K, N_r, N_t) against Y of shape (K, N_r, N): channel k decides frame k, through an internal per-draw loop – convenience, not a speedup (D51).

param alphabet:

Symbol constellation \(\mathcal{M}\) (1D array).

type alphabet:

np.ndarray

param osic_type:

Ordering strategy. Default is "sinr".

type osic_type:

{“sinr”, “colnorm”, “snr”}, keyword-only

param H:

Channel matrix \(\mathbf{H}\) of size \(N_r \times N_t\). Must be set before calling the detector.

type H:

np.ndarray, optional, keyword-only

param method:

Per-stage equalizer; overridden in __post_init__ by the ordering rule ("sinr" forces MMSE, "colnorm" and "snr" force ZF).

type method:

{“zf”, “mmse”}, keyword-only

param sigma2:

Noise variance \(\sigma^2\) (required for the "sinr" ordering).

type sigma2:

float, optional, keyword-only

param name:

Name of the detector. Default is "OSIC Detector".

type name:

str, optional, keyword-only

raises ValueError:

If osic_type is unknown, or if H (or sigma2 for the "sinr" ordering) is not set when the detector is called.

References

• Y. S. Cho, J. Kim, W. Y. Yang and C. G. Kang, MIMO-OFDM Wireless Communications with MATLAB, John Wiley & Sons, 2010.

• D. Tse and P. Viswanath, Fundamentals of Wireless Communication, Cambridge University Press, 2005, Chapter 8.

Examples

>>> alphabet = np.array([-1.0 + 0j, 1.0 + 0j])
>>> H = np.array([[1.0, 0.5], [0.0, 1.0]])
>>> Y = H @ alphabet[[1, 0]][:, None]
>>> OrderedSuccessiveInterferenceCancellationDetector(alphabet, osic_type="snr", H=H)(Y)
array([[1],
       [0]])

class ApproximateMessagePassingDetector(alphabet: ndarray, *, H: ndarray | None = None, sigma2: float | None = None, alpha: float = 1, N_it: int = 100, name: str = 'AMP Detector')

Approximate Message Passing (AMP) MIMO detector.

Signal Model

The detector assumes the flat MIMO observation model

\[\mathbf{y}[n] = \mathbf{H} \, \mathbf{x}[n] + \mathbf{b}[n], \qquad \mathbf{b}[n] \sim \mathcal{CN}\left(\mathbf{0}, \sigma^2 \mathbf{I}_{N_r}\right)\]

where \(\mathbf{H}\) is the channel matrix of size \(N_r \times N_t\) and \(\mathbf{x}[n] \in \mathcal{M}^{N_t}\). For each received vector \(\mathbf{y}\), the AMP iteration (\(t = 1, \dots, N_{it}\)) is

\[\mathbf{z}^{(t)} = \mathbf{x}^{(t)} + \mathbf{H}^H \mathbf{r}^{(t)}, \qquad \mathbf{x}^{(t+1)} = F\left(\mathbf{z}^{(t)}, \sigma^2 (1 + \tau_t^2)\right)\]

\[\mathbf{r}^{(t+1)} = \mathbf{y} - \mathbf{H}\mathbf{x}^{(t+1)} + \frac{\tau_{t+1}^2}{1 + \tau_t^2} \, \mathbf{r}^{(t)}\]

where \(F\) is the posterior-mean (soft) denoiser over the constellation \(\mathcal{M}\) and the state \(\tau_t^2\) is tracked from the denoiser output variance with the system ratio \(\beta = N_t / N_r\). The final estimate is hard-projected on \(\mathcal{M}\) and the detector returns the integer indices of the decisions in the alphabet.

Axes: declared axis – expects (ant, N) with antennas on axis -2.

param alphabet:

Symbol constellation \(\mathcal{M}\) (1D array).

type alphabet:

np.ndarray

param H:

Channel matrix \(\mathbf{H}\) of size \(N_r \times N_t\). Must be set before calling the detector.

type H:

np.ndarray, optional, keyword-only

param sigma2:

Noise variance \(\sigma^2\). Must be set before calling the detector.

type sigma2:

float, optional, keyword-only

param alpha:

Damping factor (reserved; unused by the current iteration). Default is 1.

type alpha:

float, keyword-only

param N_it:

Number of iterations \(N_{it}\). Default is 100.

type N_it:

int, keyword-only

param name:

Name of the detector. Default is "AMP Detector".

type name:

str, optional, keyword-only

raises ValueError:

If H or sigma2 is not set when the detector is called.

References

• C. Jeon, R. Ghods, A. Maleki and C. Studer, “Optimality of large MIMO detection via approximate message passing,” Proc. IEEE Int. Symp. Information Theory (ISIT), 2015, pp. 1227-1231.

• D. Tse and P. Viswanath, Fundamentals of Wireless Communication, Cambridge University Press, 2005, Chapter 8.

Examples

>>> alphabet = np.array([-1.0 + 0j, 1.0 + 0j])
>>> Y = np.array([[0.9], [-1.1]])
>>> ApproximateMessagePassingDetector(alphabet, H=np.eye(2), sigma2=0.1)(Y)
array([[1],
       [0]])

class OrthogonalApproximateMessagePassingDetector(alphabet: ndarray, *, H: ndarray | None = None, sigma2: float | None = None, alpha: float = 1, N_it: int = 100, type: Literal['H', 'pinv', 'MMSE'] = 'MMSE', name: str = 'OAMP Detector')

Orthogonal Approximate Message Passing (OAMP) MIMO detector.

Signal Model

The detector assumes the flat MIMO observation model

\[\mathbf{y}[n] = \mathbf{H} \, \mathbf{x}[n] + \mathbf{b}[n], \qquad \mathbf{b}[n] \sim \mathcal{CN}\left(\mathbf{0}, \sigma^2 \mathbf{I}_{N_r}\right)\]

where \(\mathbf{H}\) is the channel matrix of size \(N_r \times N_t\) and \(\mathbf{x}[n] \in \mathcal{M}^{N_t}\). For each received vector \(\mathbf{y}\), the OAMP iteration (\(t = 1, \dots, N_{it}\)) alternates a linear estimator and a posterior-mean denoiser \(F\) over the constellation \(\mathcal{M}\):

\[\mathbf{z}^{(t)} = \mathbf{x}^{(t)} + \mathbf{W}_t \left(\mathbf{y} - \mathbf{H}\mathbf{x}^{(t)}\right), \qquad \mathbf{x}^{(t+1)} = F\left(\mathbf{z}^{(t)}, \tau_t^2\right)\]

where \(\mathbf{W}_t\) is the matched filter \(\mathbf{H}^H\) (type="H"), the pseudo-inverse \(\mathbf{H}^{\dagger}\) (type="pinv") or the MMSE matrix \(v_t^2 \mathbf{H}^H (v_t^2 \mathbf{H}\mathbf{H}^H + \sigma^2 \mathbf{I}_{N_r})^{-1}\) (type="MMSE"). The error variances \(v_t^2\) and \(\tau_t^2\) are tracked from the residual and from \(\mathbf{B}_t = \mathbf{I}_{N_t} - \mathbf{W}_t \mathbf{H}\). The final estimate is hard-projected on \(\mathcal{M}\) and the detector returns the integer indices of the decisions in the alphabet.

Axes: declared axis – expects (ant, N) with antennas on axis -2.

param alphabet:

Symbol constellation \(\mathcal{M}\) (1D array).

type alphabet:

np.ndarray

param H:

Channel matrix \(\mathbf{H}\) of size \(N_r \times N_t\). Must be set before calling the detector.

type H:

np.ndarray, optional, keyword-only

param sigma2:

Noise variance \(\sigma^2\). Must be set before calling the detector.

type sigma2:

float, optional, keyword-only

param alpha:

Damping factor (reserved; unused by the current iteration). Default is 1.

type alpha:

float, keyword-only

param N_it:

Number of iterations \(N_{it}\). Default is 100.

type N_it:

int, keyword-only

param type:

Linear estimator \(\mathbf{W}_t\). Default is "MMSE".

type type:

{“H”, “pinv”, “MMSE”}, keyword-only

param name:

Name of the detector. Default is "OAMP Detector".

type name:

str, optional, keyword-only

raises ValueError:

If H or sigma2 is not set when the detector is called, or if type is unknown.

References

• J. Ma and L. Ping, “Orthogonal AMP,” IEEE Access, vol. 5, pp. 2020-2033, 2017.

• D. Tse and P. Viswanath, Fundamentals of Wireless Communication, Cambridge University Press, 2005, Chapter 8.

Examples

>>> alphabet = np.array([-1.0 + 0j, 1.0 + 0j])
>>> Y = np.array([[0.9], [-1.1]])
>>> OrthogonalApproximateMessagePassingDetector(alphabet, H=np.eye(2), sigma2=0.1)(Y)
array([[1],
       [0]])