Metrics#

compute_ser_awgn_psk(order: int, snr_per_bit: ndarray | float) → ndarray | float[source]#

Compute the theoretical Symbol Error Rate (SER) for PSK modulation over an AWGN channel.

Signal Model#

The receiver observes \(y[n] = x[n] + b[n]\) with \(b[n] \sim \mathcal{CN}(0, \sigma^2)\) and detects the closest \(M\)-PSK symbol. With \(k = \log_2 M\) bits per symbol, the signal-to-noise ratio per bit \(\gamma_b = E_b / N_0\) and the Gaussian tail function

\[Q(u) = \frac{1}{\sqrt{2\pi}} \int_u^{\infty} e^{-t^2/2} \, dt,\]

the closed-form symbol error probabilities are

  • \(M = 2\) (BPSK), exact:

\[P_s = Q\left(\sqrt{2 \gamma_b}\right)\]
  • \(M = 4\) (QPSK), exact:

\[P_s = 2 Q\left(\sqrt{2 \gamma_b}\right) \left[1 - \frac{1}{2} Q\left(\sqrt{2 \gamma_b}\right)\right]\]
  • \(M > 4\), the standard union-bound approximation, tight at high SNR:

\[P_s \simeq 2 Q\left(\sqrt{2 k \gamma_b} \, \sin\frac{\pi}{M}\right)\]

Note that the argument is the SNR per bit \(\gamma_b\), not the SNR per symbol \(\gamma_s = k \gamma_b\): pass snr_per_bit = 10**(snr_dB/10) / np.log2(order) when the sweep is parameterized by the symbol SNR (as AWGN(snr_dB=...) is).

param order:

Modulation order \(M\) (2, 4, 8, 16, …).

type order:

int

param snr_per_bit:

Signal-to-noise ratio per bit \(\gamma_b = E_b / N_0\) in linear scale.

type snr_per_bit:

float or np.ndarray

returns:

float or np.ndarray – Theoretical SER value(s) \(P_s\), with the same shape as snr_per_bit.

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Section 4.3 (error probability of M-ary PSK); the same result appears in Proakis, Digital Communications, McGraw-Hill, 2001, pp. 271-272.

Examples

>>> print(f"{compute_ser_awgn_psk(2, 10.0):.4e}")   # BPSK, Eb/N0 = 10 dB
3.8721e-06
>>> print(f"{compute_ser_awgn_psk(4, 10.0):.4e}")   # QPSK, Eb/N0 = 10 dB
7.7442e-06
>>> print(f"{compute_ser_awgn_psk(8, 10.0):.4e}")   # 8-PSK, Eb/N0 = 10 dB
3.0342e-03
compute_ser_awgn_qam(order: int, snr_per_bit: ndarray | float) → ndarray | float[source]#

Compute the theoretical Symbol Error Rate (SER) for square QAM modulation over an AWGN channel.

Signal Model#

A square \(M\)-QAM constellation (\(M = 4, 16, 64, 256, \dots\)) is the product of two independent \(\sqrt{M}\)-PAM constellations. With \(k = \log_2 M\) bits per symbol, the signal-to-noise ratio per bit \(\gamma_b = E_b / N_0\) and the Gaussian tail function \(Q(\cdot)\), the error probability of one PAM component is

\[P_{\sqrt{M}} = 2 \left(1 - \frac{1}{\sqrt{M}}\right) Q\left(\sqrt{\frac{3 k \gamma_b}{M-1}}\right)\]

and a symbol is correct only if both components are, hence

\[P_s = 1 - \left(1 - P_{\sqrt{M}}\right)^2\]

The argument is the SNR per bit \(\gamma_b\); the average SNR per symbol is \(\gamma_s = k \gamma_b\), so the argument of \(Q(\cdot)\) also reads \(\sqrt{3 \gamma_s / (M-1)}\).

param order:

Modulation order \(M\) (4, 16, 64, 256, …). Only square constellations, i.e. even \(k = \log_2 M\), are covered by this expression.

type order:

int

param snr_per_bit:

Signal-to-noise ratio per bit \(\gamma_b = E_b / N_0\) in linear scale.

type snr_per_bit:

float or np.ndarray

returns:

float or np.ndarray – Theoretical SER value(s) \(P_s\), with the same shape as snr_per_bit.

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Section 4.3 (error probability of QAM); the same result appears in Proakis, Digital Communications, McGraw-Hill, 2001, p. 280.

Examples

>>> print(f"{compute_ser_awgn_qam(4, 10.0):.4e}")    # QPSK, Eb/N0 = 10 dB
7.7442e-06
>>> print(f"{compute_ser_awgn_qam(16, 10.0):.4e}")   # 16-QAM, Eb/N0 = 10 dB
7.0043e-03
>>> gamma_b = 10 ** (np.array([10.0, 20.0]) / 10) / 4   # 16-QAM at 10 and 20 dB symbol SNR
>>> print(np.array2string(compute_ser_awgn_qam(16, gamma_b), precision=4))
[2.2203e-01 1.1616e-05]
compute_metric_awgn_theo(modulation: str, order: int, snr_per_bit: ndarray | float, *, metrics: Sequence[str] = ('ser', 'ber')) → dict[str, ndarray | float][source]#

Theoretical performance of a modulation over an AWGN channel.

Signal Model#

Front-end to the closed-form expressions of compute_ser_awgn_psk() and compute_ser_awgn_qam(), both parameterized by the signal-to-noise ratio per bit \(\gamma_b = E_b / N_0\). The bit error rate follows through the usual Gray-mapping approximation, in which a symbol error corrupts a single bit out of \(k = \log_2 M\):

\[P_b \simeq \frac{P_s}{k}\]

which is accurate when \(P_s \ll 1\) and the constellation is Gray-mapped, not an identity.

A figure usually wants several of these curves at once, so the return value is a dictionary: one call, one entry per metric, ready to hand to plot_error_rate().

Axes: element-wise.

param modulation:

Modulation family: "PSK" or "QAM".

type modulation:

str

param order:

Modulation order \(M\) (4, 16, 64, …).

type order:

int

param snr_per_bit:

Signal-to-noise ratio per bit \(\gamma_b = E_b / N_0\), linear.

type snr_per_bit:

float or np.ndarray

param metrics:

Which quantities to return, among "ser", "ber", "mi" and "gmi". Default ("ser", "ber") – the two closed forms. "mi" and "gmi" are the mutual information and the bit-interleaved (generalized) mutual information of the same constellation at the same SNR, from constellation_capacity() and bicm_capacity(); they are computed by quadrature rather than in closed form, and cost accordingly.

type metrics:

sequence of str, optional, keyword-only

returns:

dict of str to float or np.ndarray – One entry per requested metric, each with the shape of snr_per_bit.

raises ValueError:

If the modulation or one of the metric names is unknown.

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Section 4.3 (symbol and bit error probabilities with Gray mapping).

Examples

>>> theory = compute_metric_awgn_theo("QAM", 16, 10.0)
>>> print(f"{theory['ser']:.4e}, {theory['ber']:.4e}")
7.0043e-03, 1.7511e-03
>>> print(f"{compute_metric_awgn_theo('PSK', 8, 10.0)['ser']:.4e}")
3.0342e-03
>>> rates = compute_metric_awgn_theo("QAM", 4, 10.0, metrics=("mi", "gmi"))
>>> print(f"{rates['mi']:.3f}, {rates['gmi']:.3f}")
2.000, 2.000
compute_ser_rayleigh_psk(order: int, snr_per_bit: ndarray | float, *, diversity: int = 1) → ndarray | float[source]#

SER of M-PSK over Rayleigh fading with L-branch combining.

Signal Model#

Each of the \(L\) branches carries the symbol through an independent Rayleigh coefficient, and the receiver combines them coherently, so the post-combining SNR is a sum of \(L\) exponential variables. Averaging the AWGN error probability over that law is done with Craig’s finite-range form of the Gaussian tail function, which turns the average into

\[P_s = \frac{1}{\pi} \int_{0}^{\frac{M-1}{M}\pi} \left[ \frac{\sin^2\theta} {\sin^2\theta + g_{\mathrm{PSK}}\,\bar{\gamma}_s} \right]^{L} d\theta, \qquad g_{\mathrm{PSK}} = \sin^2\frac{\pi}{M}\]

with \(\bar{\gamma}_s = k \bar{\gamma}_b\) the average SNR per symbol per branch. The expression is exact, not a bound: the bracket is the moment generating function of one branch on Craig’s contour, and the exponent \(L\) is what independence buys.

At high SNR the error rate decays as \(\bar{\gamma}^{-L}\), i.e. \(L\) decades per decade – the diversity order is the exponent, and reading it off a simulated curve is how a scheme is checked.

Which scheme has which \(L\), at equal total transmit power:

  • one antenna each side: \(L = 1\), per-branch SNR \(\bar\gamma\);

  • maximum ratio combining over \(N_r\) receive antennas: \(L = N_r\), per-branch SNR \(\bar\gamma\);

  • an orthogonal design (Alamouti, OSTBC) on \(N_t \times N_r\): \(L = N_t N_r\), per-branch SNR \(\bar\gamma / N_t\);

  • zero forcing on \(N_r \times N_t\), per stream: \(L = N_r - N_t + 1\), per-branch SNR \(\bar\gamma\).

The transmit scheme divides the per-branch SNR by \(N_t\) because the power is split over the antennas: that division is the 3 dB an orthogonal design pays against receive diversity, and it is the whole difference between the two.

Axes: element-wise – broadcasts over snr_per_bit.

param order:

Modulation order \(M\) (2, 4, 8, …).

type order:

int

param snr_per_bit:

Average SNR per bit per branch, \(\bar{\gamma}_b\), linear.

type snr_per_bit:

float or np.ndarray

param diversity:

Number \(L\) of independent branches combined. Default 1, i.e. no diversity.

type diversity:

int, optional, keyword-only

returns:

float or np.ndarray – Symbol error rate, with the shape of snr_per_bit.

raises ValueError:

If order is below 2 or diversity below 1.

References

M. K. Simon, M.-S. Alouini, Digital Communication over Fading Channels, 2nd ed., Wiley, 2005, Sections 5.1 and 9.2 (the MGF method and Craig’s representation); J. W. Craig, “A new, simple and exact result for calculating the probability of error for two-dimensional signal constellations”, MILCOM 1991; J. G. Proakis, M. Salehi, Digital Communications, 5th ed., 2008, Section 13.4 (the closed form this quadrature reproduces for BPSK).

Examples

>>> snr = 10 ** (np.array([10.0, 20.0, 30.0]) / 10)
>>> print(" ".join(f"{v:.3e}" for v in compute_ser_rayleigh_psk(2, snr)))
2.327e-02 2.481e-03 2.498e-04
>>> two = compute_ser_rayleigh_psk(2, snr, diversity=2)
>>> print(" ".join(f"{v:.3e}" for v in two))
1.599e-03 1.844e-05 1.872e-07

One decade per decade against two: that ratio is the diversity order, and it is what a space-time code buys.

compute_ser_rayleigh_qam(order: int, snr_per_bit: ndarray | float, *, diversity: int = 1) → ndarray | float[source]#

SER of square M-QAM over Rayleigh fading with L-branch combining.

Signal Model#

Same average as compute_ser_rayleigh_psk(), over the exact AWGN error probability of a square constellation, which Craig’s form writes as a difference of two finite integrals:

\[P_s = \frac{4}{\pi}\left(1 - \frac{1}{\sqrt{M}}\right) \int_{0}^{\frac{\pi}{2}} \left[\frac{\sin^2\theta} {\sin^2\theta + g\,\bar{\gamma}_s}\right]^{L} d\theta - \frac{4}{\pi}\left(1 - \frac{1}{\sqrt{M}}\right)^{2} \int_{0}^{\frac{\pi}{4}} \left[\cdots\right]^{L} d\theta\]

with \(g = \frac{3}{2(M-1)}\) and \(\bar{\gamma}_s = k\bar{\gamma}_b\) per branch. The second term is the correction that stops the union of the two quadrature components from double-counting the corner events; dropping it gives the familiar upper bound instead of the exact value.

The scheme-to-\(L\) table of compute_ser_rayleigh_psk() applies unchanged.

Axes: element-wise – broadcasts over snr_per_bit.

param order:

Modulation order \(M\), a perfect square (4, 16, 64, …).

type order:

int

param snr_per_bit:

Average SNR per bit per branch, linear.

type snr_per_bit:

float or np.ndarray

param diversity:

Number \(L\) of combined branches. Default 1.

type diversity:

int, optional, keyword-only

returns:

float or np.ndarray – Symbol error rate.

raises ValueError:

If order is not a perfect square, or diversity below 1.

References

Simon & Alouini, Digital Communication over Fading Channels, 2nd ed., 2005, Section 9.2.3; Craig, MILCOM 1991.

Examples

>>> snr = 10 ** (np.array([10.0, 20.0, 30.0]) / 10)
>>> print(" ".join(f"{v:.3e}" for v in compute_ser_rayleigh_qam(16, snr)))
1.346e-01 1.587e-02 1.616e-03
>>> four = compute_ser_rayleigh_qam(16, snr, diversity=4)
>>> print(" ".join(f"{v:.3e}" for v in four))
7.043e-04 1.449e-07 1.570e-11
compute_metric_rayleigh_theo(modulation: str, order: int, snr_per_bit: ndarray | float, *, diversity: int = 1) → dict[str, ndarray | float][source]#

Theoretical error rate over Rayleigh fading, PSK or QAM.

Signal Model#

Front-end to compute_ser_rayleigh_psk() and compute_ser_rayleigh_qam(), the fading counterpart of compute_metric_awgn_theo() and parameterized the same way. The bit error rate follows through the Gray-mapping approximation \(P_b \simeq P_s / k\), which is an upper-SNR statement here exactly as it is over AWGN.

Axes: element-wise.

param modulation:

"PSK" or "QAM".

type modulation:

str

param order:

Modulation order \(M\).

type order:

int

param snr_per_bit:

Average SNR per bit per branch, linear.

type snr_per_bit:

float or np.ndarray

param diversity:

Number of combined branches \(L\). Default 1.

type diversity:

int, optional, keyword-only

returns:

dict of str to float or np.ndarray – {"ser": ..., "ber": ...}, each with the shape of snr_per_bit – the same shape of answer compute_metric_awgn_theo() returns, so the two can be drawn on one figure without unpacking rules that differ.

raises ValueError:

If the modulation is unknown.

References

Simon & Alouini, 2nd ed., 2005, Chapter 9.

Examples

>>> snr = 10 ** (np.array([10.0, 20.0]) / 10)
>>> theory = compute_metric_rayleigh_theo("PSK", 4, snr, diversity=2)
>>> bool(np.allclose(theory["ber"], theory["ser"] / 2))
True
compute_ser(X_target: ndarray, X_detected: ndarray, axis: int | None = None) → ndarray | float[source]#

Compute the Symbol Error Rate (SER) between the target and detected symbols.

Signal Model#

The transmitted symbol indices \(x[n]\) are compared with the detected ones \(\hat{x}[n]\) over \(N\) symbols; the SER is the normalized count of mismatches:

\[\mathrm{SER} = \frac{1}{N} \sum_{n=0}^{N-1} \mathbb{1}\left\{x[n] \neq \hat{x}[n]\right\}\]

Comparison is exact equality, so the two arrays must carry the same kind of quantity – symbol indices on both sides (the output of a SymbolGenerator and of a SymbolDemapper), not complex constellation points.

Axes: declared axis – errors are counted along axis, which must then hold the symbols; with axis=None (default) both arrays are truncated to their common length \(N = \min(N_{\text{target}}, N_{\text{detected}})\) along the last axis, then every element is pooled – so a batch of frames stays aligned row by row, and the pooled rate is the mean of the per-row rates.

param X_target:

Reference symbol indices \(x[n]\) (what was transmitted).

type X_target:

ndarray

param X_detected:

Detected symbol indices \(\hat{x}[n]\) (what the receiver decided).

type X_detected:

ndarray

param axis:

Axis along which to count the errors. If None, the arrays are raveled into 1-D arrays. Default is None.

type axis:

int or None, optional

returns:

float or ndarray – The SER. A float when axis is None, an array of SER values otherwise (one per position of the remaining axes).

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Section 4.2 (symbol error probability).

Examples

>>> x = np.array([0, 1, 2, 3, 0])
>>> x_hat = np.array([0, 1, 2, 0, 0])
>>> print(compute_ser(x, x_hat))
0.2
>>> X = np.array([[0, 1, 2, 3], [0, 1, 2, 3]])
>>> X_hat = np.array([[0, 1, 2, 0], [0, 0, 2, 3]])
>>> print(compute_ser(X, X_hat, axis=1))
[0.25 0.25]
compute_ber(X_target: ndarray, X_detected: ndarray, width: int, axis: int | None = None) → ndarray | float[source]#

Compute the Bit Error Rate (BER) between target and detected symbols.

Signal Model#

Both symbol streams are first expanded into bits by sym_2_bin(), which writes each symbol on \(w\) bits (natural binary, MSB first). The BER is the normalized count of differing bits:

\[\mathrm{BER} = \frac{N_e}{N_b}, \qquad N_e = \sum_{m=0}^{N_b-1} \mathbb{1}\left\{b[m] \neq \hat{b}[m]\right\}, \qquad N_b = N w\]

where \(b[m]\) and \(\hat{b}[m]\) are the transmitted and detected bits and \(N\) is the number of symbols. Because the bits come from the natural binary code, this BER equals the true one only if the constellation is labelled the same way; with a Gray-mapped constellation the bit labelling must match sym_2_bin for the count to be meaningful.

Axes: element-wise – with axis=None (default) both arrays are raveled and the errors are counted over the whole bit stream.

param X_target:

Reference symbol indices \(x[n]\). Should be convertible to a 1-D array.

type X_target:

ndarray

param X_detected:

Detected symbol indices \(\hat{x}[n]\) obtained from a transmission system or a decoder. Should be convertible to a 1-D array.

type X_detected:

ndarray

param width:

Number of bits per symbol \(w = \log_2 M\).

type width:

int

param axis:

Axis along which to count the errors. Default is None. Note that sym_2_bin() only accepts a 1-D input, so values other than None are not currently usable.

type axis:

int or None, optional

returns:

float or ndarray – The Bit Error Rate \(N_e / N_b\).

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Section 4.3 (bit error probability).

Examples

>>> x = np.array([0, 1, 2, 3])          # bits 00 01 10 11
>>> x_hat = np.array([0, 1, 2, 0])      # bits 00 01 10 00 -> 2 wrong bits out of 8
>>> print(compute_ber(x, x_hat, width=2))
0.25
compute_evm(X_target: ndarray, X_estimated: ndarray, axis: int | None = None) → ndarray | float[source]#

Compute the Error Vector Magnitude (EVM) between the target and estimated signals.

Signal Model#

The error vector is the difference between the estimated (measured) symbol \(\hat{x}[n]\) and the ideal reference symbol \(x[n]\). The EVM returned here is the RMS error normalized by the average power of the reference constellation:

\[\mathrm{EVM} = \sqrt{ \frac{\frac{1}{N} \sum_{n=0}^{N-1} \left|x[n] - \hat{x}[n]\right|^2} {\frac{1}{N} \sum_{n=0}^{N-1} \left|x[n]\right|^2}}\]

This is the average-power normalization used by the cellular and WLAN standards, not the peak-amplitude (outermost constellation point) normalization also found in the literature: a given impairment yields a larger figure here than under the peak convention. The result is a linear ratio (multiply by 100 for a percentage, or take \(20 \log_{10}(\mathrm{EVM})\) for dB).

Axes: declared axis – both mean powers are reduced along axis (default None: over the whole array).

param X_target:

Reference (ideal) constellation symbols \(x[n]\).

type X_target:

ndarray

param X_estimated:

Estimated (measured) symbols \(\hat{x}[n]\).

type X_estimated:

ndarray

param axis:

Axis along which the mean powers are computed. If None, the means are computed over all dimensions. Default is None.

type axis:

int or None, optional

returns:

float or ndarray – The computed EVM value. If axis is None, returns a float. Otherwise, returns an array.

References

3GPP TS 38.104, NR; Base Station (BS) radio transmission and reception – error vector magnitude, normalized by the average power of the reference signal.

Examples

>>> X_target = np.array([1, 2, 3, 4])
>>> X_estimated = np.array([1.1, 2.1, 3.1, 3.9])
>>> print(float(compute_evm(X_target, X_estimated)))
0.03651483716701111
>>> print(round(float(compute_evm(X_target, X_estimated, axis=0)), 6))
0.036515
compute_effective_snr(X_target: ndarray, X_estimated: ndarray, sigma2_s: float = 1.0, unit: str = 'natural') → ndarray | float[source]#

Compute the effective Signal-to-Noise Ratio (SNR) between a target and an estimated signal.

Signal Model#

All the impairments accumulated along the link (noise, residual distortion, nonlinear interference) are lumped into an equivalent additive term of variance

\[\sigma^2 = \frac{1}{N} \sum_{n=0}^{N-1} \left|x[n] - \hat{x}[n]\right|^2\]

where \(x[n]\) is the reference symbol and \(\hat{x}[n]\) the estimated one. The effective SNR is then

\[\mathrm{SNR}_{\mathrm{eff}} = \frac{\sigma_s^2}{\sigma^2}\]

with \(\sigma_s^2\) the nominal signal power, supplied by the caller (sigma2_s, default 1 for a unit-power constellation) and not measured on X_target. This is the reciprocal of the squared EVM when \(\sigma_s^2\) equals the reference power: \(\mathrm{SNR}_{\mathrm{eff}} = \mathrm{EVM}^{-2}\).

Axes: element-wise – both arrays are raveled and the error power is reduced over all elements.

param X_target:

Reference signal \(x[n]\).

type X_target:

ndarray

param X_estimated:

Estimated signal \(\hat{x}[n]\).

type X_estimated:

ndarray

param sigma2_s:

Nominal signal power \(\sigma_s^2\). The default is 1.

type sigma2_s:

float, optional

param unit:

Unit of the output SNR: linear scale, decibels, or decibel-milliwatts. The default is “natural”.

type unit:

{“natural”, “dB”, “dBm”}, optional

returns:

float – The effective SNR in the specified unit.

raises ValueError:

If the specified unit is not one of “natural”, “dB”, or “dBm”.

Notes

  • For "dB", the SNR is converted using \(\mathrm{SNR}_{dB} = 10 \log_{10}(\mathrm{SNR})\).

  • For "dBm", the SNR is converted using \(\mathrm{SNR}_{dBm} = 10 \log_{10}(\mathrm{SNR}) + 30\).

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Section 4.2 (signal-to-noise ratio of a memoryless additive-noise link).

Examples

>>> x = np.ones(4)
>>> x_hat = np.array([1.1, 0.9, 1.1, 0.9])   # error power 0.01
>>> print(round(float(compute_effective_snr(x, x_hat)), 3))
100.0
>>> print(round(float(compute_effective_snr(x, x_hat, unit="dB")), 3))
20.0
compute_power(x: ndarray, unit: str = 'natural') → float[source]#

Compute the mean power of an input array.

Signal Model#

The mean power of the signal \(x[n]\) of length \(N\) is the time average of its squared modulus:

\[P_x = \frac{1}{N} \sum_{n=0}^{N-1} \left|x[n]\right|^2\]

expressed in the requested unit:

\[P_{dB} = 10 \log_{10} \left(P_x\right), \qquad P_{dBm} = 10 \log_{10} \left(P_x\right) + 30\]

The dBm conversion assumes \(P_x\) is expressed in watts.

Axes: element-wise – the average is reduced over all elements of the array.

param x:

Input signal \(x[n]\).

type x:

ndarray

param unit:

Unit of the output power: watts, decibels, or decibel-milliwatts. The default is “natural”.

type unit:

{“natural”, “dB”, “dBm”}, optional

returns:

float – Mean power \(P_x\) in the specified unit.

raises ValueError:

If the specified unit is not one of “natural”, “dB”, or “dBm”.

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Section 2.2 (energy and power of a signal).

Examples

>>> x = np.array([1 + 1j, 1 - 1j])
>>> print(round(float(compute_power(x)), 6))
2.0
>>> print(round(float(compute_power(x, unit="dB")), 4))
3.0103
compute_ccdf(data: ndarray, axis: int = -1) → tuple[ndarray, ndarray][source]#

Compute the empirical Complementary Cumulative Distribution Function (CCDF) of a dataset.

Signal Model#

The CCDF of a random variable \(Z\) is the probability of exceeding a threshold:

\[\mathrm{CCDF}(z) = \Pr\left\{Z > z\right\} = 1 - F_Z(z)\]

Given \(N\) samples sorted in ascending order \(z_{(1)} \le \dots \le z_{(N)}\), the empirical estimate returned for the \(i\)-th sorted value is

\[\widehat{\mathrm{CCDF}}\left(z_{(i)}\right) = 1 - \frac{i}{N}, \qquad i = 1, \dots, N\]

i.e. the fraction of samples strictly above \(z_{(i)}\) (the last entry is therefore exactly 0). This is the standard way of plotting the PAPR distribution of an OFDM signal: feed comnumpy.ofdm.metrics.compute_papr() values in, then plot ccdf against sorted_data on a semilog-y axis.

Axes: declared axis – sorting and the sample count \(N\) use axis (default -1).

param data:

Input samples \(z\) whose CCDF is estimated.

type data:

array-like

param axis:

Axis along which to compute the CCDF. Default is the last axis.

type axis:

int, optional

returns:
  • sorted_data (np.ndarray) – The input data sorted in ascending order along the specified axis.

  • ccdf (np.ndarray) – The CCDF values \(1 - i/N\) matching the sorted data.

References

R. van Nee, R. Prasad, OFDM for Wireless Multimedia Communications, Artech House, 2000, Chapter 6 (PAPR distribution).

Examples

>>> data = np.array([3.0, 1.0, 2.0, 4.0])
>>> sorted_data, ccdf = compute_ccdf(data)
>>> print(sorted_data)
[1. 2. 3. 4.]
>>> print(ccdf)
[0.75 0.5  0.25 0.  ]
compute_acpr(signal: ndarray, bandwidth: float, sampling_rate: float) → tuple[float, float][source]#

Calculate the Adjacent Channel Power Ratio (ACPR) of a given signal.

Signal Model#

ACPR measures spectral regrowth: the power a nonlinear transmitter leaks into the neighbouring channels, relative to the power it emits in its own channel. With \(X(f)\) the DFT of the baseband signal \(x[n]\) sampled at \(f_s\), the power in a band \(\mathcal{B}\) is estimated by Parseval’s relation as

\[P(\mathcal{B}) = \frac{1}{N} \sum_{f_i \in \mathcal{B}} \left|X(f_i)\right|^2\]

The three bands are deduced from the channel bandwidth \(B\), the main channel being centred on DC:

\[\mathcal{B}_{\mathrm{main}} = \left[-\tfrac{B}{2}, \tfrac{B}{2}\right], \quad \mathcal{B}_{\mathrm{right}} = \left[\tfrac{B}{2}, \tfrac{3B}{2}\right], \quad \mathcal{B}_{\mathrm{left}} = \left[-\tfrac{3B}{2}, -\tfrac{B}{2}\right]\]

and each ratio is reported in decibels:

\[\mathrm{ACPR}_{dB} = 10 \log_{10} \frac{P(\mathcal{B}_{\mathrm{adj}})}{P(\mathcal{B}_{\mathrm{main}})}\]

A clean transmitter therefore gives a large negative value; the less negative the ACPR, the more spectral regrowth. The sampling rate must satisfy \(f_s \ge 3B\) for the two adjacent bands to lie inside the represented spectrum.

Axes: element-wise – the DFT and the band powers are computed over the whole 1-D input; no axis argument.

param signal:

Input time-domain baseband signal \(x[n]\).

type signal:

np.ndarray

param bandwidth:

Bandwidth \(B\) of the main channel, in Hz.

type bandwidth:

float

param sampling_rate:

Sampling rate \(f_s\) of the signal, in Hz.

type sampling_rate:

float

returns:
  • acpr_right (float) – ACPR (in dB) of the upper adjacent channel.

  • acpr_left (float) – ACPR (in dB) of the lower adjacent channel.

References

S. C. Cripps, RF Power Amplifiers for Wireless Communications, 2nd ed., Artech House, 2006, Chapter 9 (spectral regrowth and adjacent channel power).

Examples

>>> f_s, B = 100.0, 20.0
>>> n = np.arange(1000)
>>> x = np.cos(2 * np.pi * 4 * n / f_s)     # single tone inside the channel
>>> y = x - 0.2 * x ** 3                    # cubic nonlinearity -> IM3 at 12 Hz
>>> print(np.round(compute_acpr(y, B, f_s), 4))
[-27.6193 -27.6193]
signal_report(x: ndarray, compute_papr: bool = False, papr_unit: str = 'dB') → dict[str, float][source]#

Summary statistics of a signal, as plain data.

Replaces the former in-chain monitor blocks: extract the signal with Sequential(observations=...), compute the report, and let the caller decide how to present it (logging, table, assertion, …).

Signal Model#

All statistics are computed on the modulus \(|x[n]|\) over every element; the average power is

\[\widehat{P} = \frac{1}{N} \sum_{n=0}^{N-1} \left|x[n]\right|^2\]

Axes: element-wise – statistics are reduced over all elements.

param x:

Input signal.

type x:

np.ndarray

param compute_papr:

Also include the PAPR. Default is False.

type compute_papr:

bool, optional

param papr_unit:

Unit of the PAPR entry. Default is “dB”.

type papr_unit:

{“dB”, “natural”}, optional

returns:

dict of str to float – Keys: min, max, mean, std, rms, energy, avg_power and optionally papr.

References

J. G. Proakis, M. Salehi, Digital Communications, 5th ed., McGraw-Hill, 2008, Section 2.2 (signal energy and power).

Examples

>>> report = signal_report(np.array([1.0+0.0j, 0.0+1.0j]))
>>> print(report["avg_power"], report["rms"])
1.0 1.0
class ErrorCounter(width: int = 1)[source]#

Accumulate errors over the trials of one Monte-Carlo point.

Signal Model#

A Monte-Carlo point runs the chain several times and reports one error rate. The rate is a ratio of totals, not a mean of ratios:

\[\mathrm{rate} = \frac{\sum_t N_{e,t}}{\sum_t N_t}\]

Averaging the per-trial rates instead gives the same answer only when every trial has the same length, and it throws away the number the point’s credibility rests on – how many errors were actually seen. A rate of \(10^{-4}\) from three errors and one from three hundred are different claims, and only the counts tell them apart.

That is what this object is for. It does not own the loop: a Monte-Carlo over a fading channel, over several detectors or over a split-step link has a different shape each time, and that shape is what the reader of an example came to see. It owns the counting, which is the same everywhere and the part that hides mistakes.

param width:

Bits per symbol. The default 1 counts symbol errors, one per mismatched entry. Give the constellation’s bits_per_symbol to count bit errors instead, expanded by sym_2_bin() exactly as compute_ber() does.

type width:

int, optional

n_errors#

Errors accumulated so far – symbols or bits, per width.

Type:

int

n_symbols#

Symbols compared so far.

Type:

int

n_trials#

Number of update() calls.

Type:

int

raises ValueError:

If width is not positive.

Examples

Three trials of a hundred symbols, one error in the second:

>>> import numpy as np
>>> counter = ErrorCounter()
>>> sent = np.zeros(100, dtype=int)
>>> for trial in range(3):
...     received = sent.copy()
...     if trial == 1:
...         received[7] = 1
...     _ = counter.update(sent, received)
>>> print(counter)
ErrorCounter(1 error / 300 symbols, rate 3.33e-03)

The counts are what makes the rate readable, and they are what an example should print beside it:

>>> print(counter.n_errors, counter.n_symbols, counter.n_trials)
1 300 3

width switches to bits, on the same convention as compute_ber():

>>> bits = ErrorCounter(width=4)
>>> _ = bits.update(np.array([0, 0]), np.array([0, 15]))
>>> print(bits.n_errors, bits.n_symbols, f"{bits.rate:.3f}")
4 2 0.500
update(reference: ndarray, detected: ndarray) → ErrorCounter[source]#

Count one trial, and return self so calls can be chained.

Both arrays carry symbol indices, as compute_ser() requires: the output of a SymbolGenerator against that of a SymbolDemapper, not complex constellation points. The shorter of the two sets the length, which is what a chain with a filter delay produces.

property rate: float#

Errors over the units compared – symbols, or bits if width.

reset() → ErrorCounter[source]#

Forget every trial, keeping width. Returns self.