Shaping#

Probabilistic shaping: choosing how often each symbol is sent.

A uniform constellation is not the best a channel can do. The capacity of the Gaussian channel is reached by a Gaussian input, and a uniform distribution over a square QAM grid loses up to 1.53 dB of it – the shaping gap, \(10 \log_{10}(\pi e / 6)\). Probabilistic shaping closes that gap by sending inner constellation points more often than outer ones, at the same average power.

Two things are needed for that, and this module provides both.

The target distribution. For an average-power constraint the entropy-maximizing distribution over a discrete constellation is the Maxwell-Boltzmann family maxwell_boltzmann(), and nothing else: the shape is forced by the constraint, only its parameter is free.

A way to produce it from data. Data is uniform, so a distribution matcher is an invertible map from uniform bit strings onto sequences with the wanted empirical distribution. Two classical constructions are implemented, and they differ in what they hold fixed:

ConstantCompositionMatcher (CCDM) fixes the composition – every output block contains exactly \(n_i\) copies of symbol \(i\) – and enumerates the permutations of that multiset;
SphereShaper (ESS) fixes an energy budget – every output block satisfies \(\sum_j e(a_j) \leq E_{\max}\) – and enumerates the sequences inside that sphere. It is the better of the two at short blocklengths, because it does not force each block to reproduce the target distribution exactly.

Both are exactly invertible: they are enumerative codes, ranking and unranking a finite set with integer arithmetic, so nothing is approximated and decode(encode(bits)) == bits holds by construction rather than by tolerance.

Where this sits in a transmitter (PAS). The standard architecture is probabilistic amplitude shaping: the matcher shapes the amplitudes, a systematic FEC encoder produces the parity bits, and those parity bits – which are uniform – become the signs. A sign is equiprobable, so it costs nothing in shaping and the composite constellation keeps the symmetric Maxwell-Boltzmann distribution. That is why the matchers here work on non-negative amplitudes and why their alphabet is half a PAM constellation.

References

G. Böcherer, F. Steiner, P. Schulte, “Bandwidth efficient and rate-matched low-density parity-check coded modulation”, IEEE Trans. Commun., vol. 63, no. 12, pp. 4651-4665, Dec. 2015 (PAS); P. Schulte, G. Böcherer, “Constant composition distribution matching”, IEEE Trans. Inf. Theory, vol. 62, no. 1, pp. 430-434, Jan. 2016 (CCDM); F. M. J. Willems, J. J. Wuijts, “A pragmatic approach to shaped coded modulation”, IEEE Symp. Commun. Veh. Technol., 1993 (enumerative shell mapping); Y. C. Gültekin, W. J. van Houtum, A. G. C. Koonen, F. M. J. Willems, “Enumerative sphere shaping for wireless communications with short packets”, IEEE Trans. Wireless Commun., vol. 19, no. 2, pp. 1098-1112, Feb. 2020; G. D. Forney, L.-F. Wei, “Multidimensional constellations – Part I”, IEEE J. Sel. Areas Commun., vol. 7, no. 6, pp. 877-892, 1989 (the 1.53 dB shaping gain).

Examples

>>> import numpy as np
>>> amplitudes = np.array([1.0, 3.0, 5.0, 7.0])       # half of 8-PAM
>>> target = maxwell_boltzmann(amplitudes, entropy=1.5)
>>> print(np.round(target, 4))
[0.5376 0.322  0.1156 0.0248]
>>> matcher = ConstantCompositionMatcher(amplitudes, distribution=target,
...                                      length=32)
>>> matcher.n_bits, round(matcher.rate, 4)
(42, 1.3125)
>>> bits = np.zeros(42, dtype=int); bits[::3] = 1
>>> block = matcher.encode(bits)
>>> bool(np.array_equal(matcher.decode(block), bits))    # exactly invertible
True

maxwell_boltzmann(alphabet: ndarray, *, lam: float | None = None, entropy: float | None = None) → ndarray[source]#

The distribution shaping uses, and the only one it can use.

Signal Model#

Among all distributions on a fixed constellation with a given average energy, the one of maximum entropy is

\[p_i = \frac{e^{-\lambda \left|a_i\right|^2}} {\sum_j e^{-\lambda \left|a_j\right|^2}}, \qquad \lambda \geq 0\]

the Maxwell-Boltzmann family. This is not a modelling choice: it is what maximizing \(H(P)\) subject to \(\sum_i p_i |a_i|^2 \leq E\) gives, \(\lambda\) being the Lagrange multiplier of the power constraint. \(\lambda = 0\) returns the uniform distribution; growing \(\lambda\) moves mass onto the inner points, lowering both the energy and the entropy.

Since the entropy is strictly decreasing in \(\lambda\), specifying a target entropy determines \(\lambda\) uniquely; that is the useful parameterization in practice, because the entropy is the rate the constellation carries.

Axes: element-wise – the alphabet is a 1-D array of constellation points, and the result is a distribution over it.

param alphabet:: Constellation points \(a_i\), real or complex.
type alphabet:: np.ndarray
param lam:: The parameter \(\lambda\) itself, non-negative.
type lam:: float, optional, keyword-only
param entropy:: Target entropy \(H(P)\) in bits per symbol, in \((0, \log_2 M]\). Exactly one of lam and entropy must be given (D41).
type entropy:: float, optional, keyword-only
returns:: np.ndarray – The distribution \(p_i\), same length as alphabet.
raises ValueError:: If both or neither parameterization is given, if lam is negative, or if entropy is outside \((0, \log_2 M]\).

References

F. R. Kschischang, S. Pasupathy, “Optimal nonuniform signaling for Gaussian channels”, IEEE Trans. Inf. Theory, vol. 39, no. 3, pp. 913-929, 1993; Böcherer, Steiner, Schulte, IEEE Trans. Commun. 63(12), 2015, Section III.

Examples

>>> amplitudes = np.array([1.0, 3.0, 5.0, 7.0])
>>> print(np.round(maxwell_boltzmann(amplitudes, lam=0.0), 4))
[0.25 0.25 0.25 0.25]
>>> shaped = maxwell_boltzmann(amplitudes, entropy=1.5)
>>> print(round(distribution_entropy(shaped), 6))
1.5
>>> float(np.sum(shaped * amplitudes ** 2)) < 21.0    # uniform: 21
True

blahut_arimoto(alphabet: ndarray, *, sigma2: float, lam: float | None = None, energy: float | None = None, n_nodes: int = 40, tol: float = 1e-09, max_iter: int = 5000) → ndarray[source]#

The distribution that maximizes the mutual information, computed.

Signal Model#

maxwell_boltzmann() answers “which distribution has the largest entropy at a given energy?”, and has a closed form. The question an engineer actually cares about is a different one – “which distribution carries the most bits over this channel at a given energy?” – and it has no closed form at all:

\[\max_{P_X} \; I(X;Y) \quad \text{subject to} \quad \sum_i p_i \left|a_i\right|^2 \leq E\]

over the real AWGN channel \(Y = X + Z\), \(Z \sim \mathcal{N}(0, \sigma^2)\), with \(X\) drawn from the fixed constellation. The problem is concave in \(P_X\), and the Blahut-Arimoto algorithm solves it by alternating maximization: with the Lagrange multiplier \(\lambda\) of the energy constraint,

\[q(i \mid y) = \frac{p_i \, f(y \mid a_i)} {\sum_j p_j \, f(y \mid a_j)}, \qquad p_i \;\leftarrow\; \frac{e^{D_i - \lambda |a_i|^2}} {\sum_j e^{D_j - \lambda |a_j|^2}},\]

where \(D_i = \int f(y \mid a_i) \log q(i \mid y)\,\mathrm{d}y\) is evaluated by Gauss-Hermite quadrature: the substitution \(y = a_i + \sqrt{2}\,\sigma t\) turns the Gaussian integral into the quadrature’s own weight function exactly, so \(D_i\) is a weighted sum over n_nodes points and nothing is sampled.

\(\lambda\) plays the same role as in the Maxwell-Boltzmann family and is measured in the same units, nats per unit energy, so the two can be compared directly – and should be, because the result is the reference this module’s closed form is judged against: Maxwell-Boltzmann is not the maximizer of the mutual information, only of the entropy, and how little that costs is a number worth measuring rather than assuming (Kschischang and Pasupathy, 1993).

Warning

The energy constraint is what makes the answer bell-shaped. Left unconstrained (lam=0) the maximizer may spend more energy than the uniform distribution, pushing mass onto the outer points to separate them further – at low SNR it does exactly that. A distribution that concentrates on the inner points is the answer to a constrained problem, never to an unconstrained one.

Axes: element-wise – the alphabet is a 1-D array of constellation points, and the result is a distribution over it.

param alphabet:: Constellation points \(a_i\), real (a complex alphabet with a zero imaginary part is accepted). A square QAM is the product of two independent PAM axes; pass one of them.
type alphabet:: np.ndarray
param sigma2:: Noise variance \(\sigma^2\) of the real AWGN channel, positive. Only the ratio of the alphabet’s energy to \(\sigma^2\) matters, so this is where the SNR enters.
type sigma2:: float
param lam:: The Lagrange multiplier \(\lambda\) itself, non-negative. Useful when the multiplier is the object of interest; note that the energy it lands on is an output, not a choice.
type lam:: float, optional, keyword-only
param energy:: Average energy budget \(\sum_i p_i a_i^2\) the answer must spend, reached by a root-find on \(\lambda\). This is the parameterization with an operational meaning – a transmitter has a power budget, not a Lagrange multiplier – and the one to use when comparing against maxwell_boltzmann(), which only makes sense at equal energy. Exactly one of lam and energy must be given (D41).
type energy:: float, optional, keyword-only
param n_nodes:: Gauss-Hermite nodes for \(D_i\). Default 40. Measured on 16-PAM at 10 dB against a 200-node reference, the largest probability moves by 3.3e-7 at 20 nodes, 7.7e-10 at 40 and 3.2e-13 at 80; a cleaner channel converges faster still, so 40 is well inside what the result is read to.
type n_nodes:: int, optional, keyword-only
param tol:: Stop when Blahut’s bound puts the objective within this many nats of its maximum. Default 1e-9.
type tol:: float, optional, keyword-only
param max_iter:: Iteration ceiling. Default 5000. Reaching it is not silent: the remaining bound is logged, so the answer always comes with the distance it may still be from the maximum.
type max_iter:: int, optional, keyword-only
returns:: np.ndarray – The maximizing distribution \(p_i\), same length as alphabet.
raises ValueError:: If both or neither parameterization is given, if sigma2 is not positive, if lam is negative, if the alphabet is complex, or if energy is outside the range the constraint can reach.

References

R. E. Blahut, “Computation of channel capacity and rate-distortion functions”, IEEE Trans. Inf. Theory, vol. 18, no. 4, pp. 460-473, July 1972; S. Arimoto, “An algorithm for computing the capacity of arbitrary discrete memoryless channels”, IEEE Trans. Inf. Theory, vol. 18, no. 1, pp. 14-20, Jan. 1972; F. R. Kschischang, S. Pasupathy, “Optimal nonuniform signaling for Gaussian channels”, IEEE Trans. Inf. Theory, vol. 39, no. 3, pp. 913-929, May 1993 (that Maxwell-Boltzmann is within a hundredth of a bit of this).

Examples

The uniform distribution is optimal when the constraint is loose and the channel is clean – a 4-PAM at 20 dB has nothing to gain:

>>> pam4 = np.array([-3.0, -1.0, 1.0, 3.0])
>>> clean = blahut_arimoto(pam4, sigma2=5.0 / 100, lam=0.0)
>>> print(np.round(clean, 4))
[0.25 0.25 0.25 0.25]

Impose a budget and the answer becomes the bell shape shaping is named after. The uniform law spends 5 on this alphabet, so asking for 3 is asking for 2.2 dB less power:

>>> best = blahut_arimoto(pam4, sigma2=5.0 / 100, energy=3.0)
>>> print(np.round(best, 4))
[0.125 0.375 0.375 0.125]
>>> print(round(float(np.sum(best * pam4 ** 2)), 6))
3.0

At that energy the closed form is the same law, which is the whole reason the closed form is used everywhere else:

>>> print(np.round(maxwell_boltzmann(pam4, lam=0.1373), 4))
[0.125 0.375 0.375 0.125]

distribution_entropy(distribution: ndarray) → float[source]#

Entropy of a symbol distribution, in bits.

Signal Model#

\[H(P) = -\sum_{i} p_i \log_2 p_i\]

It is the number of bits a symbol of that distribution carries, hence the ceiling any distribution matcher can approach and the quantity its rate loss is measured against.

Axes: element-wise – the distribution is a 1-D array of probabilities.

param distribution:: Probabilities \(p_i\), non-negative and summing to one.
type distribution:: np.ndarray
returns:: float – Entropy in bits per symbol, in \([0, \log_2 M]\).
raises ValueError:: If the probabilities are negative or do not sum to one.

References

T. M. Cover, J. A. Thomas, Elements of Information Theory, 2nd ed., Wiley, 2006, Section 2.1.

Examples

>>> print(round(distribution_entropy(np.full(8, 1 / 8)), 4))
3.0
>>> print(round(distribution_entropy(np.array([0.5, 0.25, 0.25])), 4))
1.5

composition_from_distribution(distribution: ndarray, length: int) → tuple[int, ...][source]#

Round a distribution to an integer composition of a block.

Signal Model#

A constant-composition matcher emits blocks holding exactly \(n_i\) copies of symbol \(i\), so it needs integers, not probabilities. The composition closest to \(n p_i\) under the constraint \(\sum_i n_i = n\) is obtained by the largest-remainder rule: take \(\lfloor n p_i \rfloor\) for every symbol, then give the \(n - \sum_i \lfloor n p_i \rfloor\) remaining slots to the symbols with the largest fractional parts.

The composition is not the distribution: it is a quantization of it, and the difference is one reason a finite block loses rate. It vanishes as \(n\) grows.

Axes: element-wise – a 1-D distribution in, a tuple of counts out.

param distribution:: Target probabilities, summing to one.
type distribution:: np.ndarray
param length:: Block length \(n\), the total number of symbols.
type length:: int
returns:: tuple of int – The counts \(n_i\), summing to length.
raises ValueError:: If length is not positive, or the distribution is invalid.

References

P. Schulte, G. Böcherer, IEEE Trans. Inf. Theory 62(1), 2016, Section IV (choice of the composition).

Examples

>>> composition_from_distribution(np.array([0.5, 0.25, 0.25]), 8)
(4, 2, 2)
>>> composition_from_distribution(np.array([0.4, 0.35, 0.25]), 10)
(4, 4, 2)

shaping_gain_dB(alphabet: ndarray, distribution: ndarray) → float[source]#

Energy saved by a shaped distribution, at equal rate.

Signal Model#

Shaping is only worth something if it lowers the energy without lowering the rate, so the comparison has to be made at equal entropy. On a one-dimensional grid of spacing \(\Delta\), a source of entropy \(H\) spread uniformly would occupy an interval of width \(2^{H}\Delta\) and cost

\[E_{\mathrm{unif}} = \frac{\left(2^{H} \Delta\right)^2}{12}\]

against the \(E_P = \sum_i p_i a_i^2\) the shaped distribution actually costs. The gain is their ratio,

\[G_s = 10 \log_{10} \frac{E_{\mathrm{unif}}}{E_P}\]

A uniform distribution over \(M\) points gives \(10\log_{10}(M^2/(M^2-1))\), i.e. essentially 0 dB, which is the calibration this definition needs to be worth anything. The supremum over all one-dimensional distributions is \(10\log_{10}(\pi e/6) = 1.53\) dB, reached in the limit of a Gaussian over an infinitely fine, infinitely wide grid.

Axes: element-wise – a 1-D real alphabet and its distribution.

param alphabet:: Constellation points, real and on a regular grid. A complex constellation is two independent dimensions; pass one of them.
type alphabet:: np.ndarray
param distribution:: Probability of each point.
type distribution:: np.ndarray
returns:: float – Shaping gain in dB, at most 1.53.
raises ValueError:: If the alphabet is complex, or has fewer than two points.

References

G. D. Forney, L.-F. Wei, “Multidimensional constellations – Part I: Introduction, figures of merit, and generalized cross constellations”, IEEE J. Sel. Areas Commun., vol. 7, no. 6, pp. 877-892, 1989, Section III (shaping gain and the 1.53 dB ultimate gain).

Examples

>>> pam8 = np.arange(-7, 8, 2).astype(float)
>>> uniform = np.full(8, 1 / 8)
>>> print(round(shaping_gain_dB(pam8, uniform), 4))     # the calibration
0.0684
>>> shaped = maxwell_boltzmann(pam8, entropy=2.5)
>>> print(round(shaping_gain_dB(pam8, shaped), 4))
1.5055

class ConstantCompositionMatcher(alphabet: ndarray, *, composition: tuple[int, ...] = (), distribution: ndarray | None = None, length: int = 0)[source]#

Constant composition distribution matching (CCDM).

Signal Model#

The matcher is an enumerative code on the permutations of a fixed multiset. Given a composition \((n_1, \ldots, n_M)\) with \(\sum_i n_i = n\), the set of output blocks is every sequence containing exactly \(n_i\) copies of symbol \(i\), and there are

\[N = \binom{n}{n_1, \ldots, n_M} = \frac{n!}{n_1! \, n_2! \cdots n_M!}\]

of them. The matcher maps \(k = \lfloor \log_2 N \rfloor\) input bits onto them by unranking: the integer the bits spell out is the lexicographic index of the block. Reading the index back off the block – ranking – inverts it exactly, which is why the pair is a bijection and not an approximation.

Every block has the same empirical distribution \(n_i / n\) by construction, so the shaped statistics are exact at every block length. What a finite block costs is rate: the matcher carries \(k/n\) bits per symbol where the composition itself is worth \(H(n_i/n)\), and

\[R_{\mathrm{loss}} = H\left(\frac{n_i}{n}\right) - \frac{k}{n} \;\xrightarrow[n \to \infty]{}\; 0\]

Axes: element-wise – a block is a 1-D sequence of length symbol indices; batches are handled by DistributionMatcher.

param alphabet:: The symbols the matcher emits, usually the non-negative amplitudes of a PAM constellation (see the module docstring on PAS). Only its length is used by the enumeration.
type alphabet:: np.ndarray
param composition:: The counts \(n_i\) directly.
type composition:: tuple of int, optional, keyword-only
param distribution:: A target distribution, quantized to a composition by composition_from_distribution(). Requires length. Exactly one of composition and distribution is given (D41).
type distribution:: np.ndarray, optional, keyword-only
param length:: Block length \(n\), required with distribution and forbidden with composition, which already carries it.
type length:: int, optional, keyword-only

n_bits#

\(k\), the number of input bits per block.

Type:: int

length#

\(n\), the number of output symbols per block.

Type:: int

raises ValueError:: If the parameterization is ambiguous, or if the composition supports fewer than one bit.

References

P. Schulte, G. Böcherer, “Constant composition distribution matching”, IEEE Trans. Inf. Theory, vol. 62, no. 1, pp. 430-434, Jan. 2016.

Examples

>>> matcher = ConstantCompositionMatcher(np.array([1.0, 3.0]),
...                                      composition=(6, 2))
>>> matcher.n_bits, matcher.length
(4, 8)
>>> block = matcher.encode(np.array([1, 0, 1, 1]))
>>> print(block)
[0 0 1 0 0 0 1 0]
>>> print(np.bincount(block))            # the composition, exactly
[6 2]
>>> print(matcher.decode(block))
[1 0 1 1]

property rate: float#: Bits per output symbol, \(k/n\).

property rate_loss: float#: What the finite block costs, in bits per symbol.

encode(bits: ndarray) → ndarray[source]#: Map n_bits bits onto one block of symbol indices.

decode(block: ndarray) → ndarray[source]#: Read the n_bits bits back off one block of indices.

class SphereShaper(alphabet: ndarray, *, length: int, max_energy: int = -1, n_bits: int = -1)[source]#

Enumerative sphere shaping (ESS).

Signal Model#

Instead of forcing every block to have the same composition, ESS keeps every block inside an energy sphere: with an integer energy \(e_i\) attached to each amplitude, the code is

\[\mathcal{C} = \left\{ (s_1, \ldots, s_n) : \sum_{j=1}^{n} e_{s_j} \leq E_{\max} \right\}\]

and the matcher enumerates it. Counting is a one-dimensional recursion over the remaining budget,

\[N_t(E) = \sum_{i} N_{t-1}\left(E - e_i\right), \qquad N_0(E) = 1 \;\; \text{for } E \geq 0\]

from which the lexicographic unranking follows exactly as for CCDM.

Compared with a constant composition, the sphere is a larger set at the same average energy, because it lets a block spend more on one symbol and less on another. That is worth up to a few tenths of a dB at the block lengths where shaping is actually deployed: CCDM needs thousands of symbols to be efficient, ESS is already good at a few dozen. What it gives up is the exact per-block distribution, which is only reproduced on average.

The energies are derived from the alphabet: \(|a_i|^2\) shifted so the cheapest amplitude costs zero, then divided by the greatest common divisor. Both steps are exact on the odd-integer grids used in practice (amplitudes 1, 3, 5, 7 give energies 0, 1, 3, 6) and only change the labelling of the budget, not the code.

Axes: element-wise – a block is a 1-D sequence of length symbol indices.

param alphabet:: Amplitudes, on a grid whose squared values are integers after shifting (odd integers, as produced by get_alphabet("PAM", M, norm=False) restricted to its positive half).
type alphabet:: np.ndarray
param length:: Block length \(n\).
type length:: int, keyword-only
param max_energy:: The budget \(E_{\max}\), in the reduced integer unit.
type max_energy:: int, optional, keyword-only
param n_bits:: Target number of input bits; the smallest budget carrying that many bits is used. Exactly one of the two is given (D41).
type n_bits:: int, optional, keyword-only

n_bits#

\(k = \lfloor \log_2 |\mathcal{C}| \rfloor\).

Type:: int

max_energy#

The budget actually used.

Type:: int

energies#

The reduced integer energy of each amplitude. The rate property gives \(k/n\).

Type:: tuple of int

raises ValueError:: If the parameterization is ambiguous, if the squared amplitudes are not integers after shifting, or if the budget carries no bit.

References

F. M. J. Willems, J. J. Wuijts, “A pragmatic approach to shaped coded modulation”, IEEE Symp. Commun. Veh. Technol., 1993; Y. C. Gültekin et al., “Enumerative sphere shaping for wireless communications with short packets”, IEEE Trans. Wireless Commun., vol. 19, no. 2, pp. 1098-1112, Feb. 2020.

Examples

>>> shaper = SphereShaper(np.array([1.0, 3.0]), length=8, max_energy=3)
>>> shaper.energies, shaper.n_bits
((0, 1), 6)
>>> block = shaper.encode(np.array([1, 0, 1, 1, 0, 1]))
>>> print(block)
[0 1 0 0 0 0 1 1]
>>> int(np.sum(block))                   # at most max_energy = 3
3
>>> print(shaper.decode(block))
[1 0 1 1 0 1]

encode(bits: ndarray) → ndarray[source]#: Map n_bits bits onto one block of symbol indices.

decode(block: ndarray) → ndarray[source]#: Read the n_bits bits back off one block of indices.

property rate: float#: Bits per output symbol, \(k/n\).

class DistributionMatcher(shaper: ConstantCompositionMatcher | SphereShaper, *, name: str = 'distribution matcher')[source]#

Shape a bit stream into symbol indices, block by block.

Signal Model#

The chain-level wrapper of ConstantCompositionMatcher and SphereShaper: it cuts the input into blocks of \(k\) bits, maps each one onto \(n\) symbol indices, and concatenates. Its output is what SymbolMapper consumes, so a shaped transmitter is a chain like any other.

Axes: axis -1 – the bit stream is read in blocks of n_bits along the last axis, whose length changes (\(k\) bits in, \(n\) amplitudes out); leading axes are carried through unchanged.

param shaper:: The enumerative code doing the work.
type shaper:: ConstantCompositionMatcher or SphereShaper
param name:: Block name.
type name:: str, optional, keyword-only
raises ShapeError:: If the input length is not a multiple of shaper.n_bits.

References

G. Böcherer, F. Steiner, P. Schulte, IEEE Trans. Commun. 63(12), 2015 (the PAS transmitter this block sits in).

Examples

>>> shaper = ConstantCompositionMatcher(np.array([1.0, 3.0]),
...                                     composition=(6, 2))
>>> bits = np.array([1, 0, 1, 1, 0, 0, 0, 1])
>>> indices = DistributionMatcher(shaper)(bits)
>>> indices.size, int(np.sum(indices == 1))
(16, 4)
>>> print(DistributionDematcher(shaper)(indices))
[1 0 1 1 0 0 0 1]

class DistributionDematcher(shaper: ConstantCompositionMatcher | SphereShaper, *, name: str = 'distribution dematcher')[source]#

Read the bit stream back off shaped symbol indices.

Signal Model#

The exact inverse of DistributionMatcher: blocks of \(n\) indices become blocks of \(k\) bits. Since the enumerative code is a bijection onto its block set, this is an identity on any sequence the matcher produced – and an error on any other, which is the honest answer after a detection mistake.

Axes: axis -1 – indices are read in blocks of length along the last axis, whose length changes (\(n\) amplitudes in, \(k\) bits out); leading axes are carried through unchanged.

param shaper:: The same code used to shape.
type shaper:: ConstantCompositionMatcher or SphereShaper
param name:: Block name.
type name:: str, optional, keyword-only
raises ShapeError:: If the input length is not a multiple of shaper.length.

References

Schulte, G. Böcherer, IEEE Trans. Inf. Theory 62(1), 2016.

Examples

>>> shaper = SphereShaper(np.array([1.0, 3.0]), length=8, max_energy=3)
>>> bits = np.array([1, 0, 1, 1, 0, 1])
>>> indices = DistributionMatcher(shaper)(bits)
>>> print(DistributionDematcher(shaper)(indices))
[1 0 1 1 0 1]

class AmplitudeMapper(alphabet: ndarray, *, seed: int | None = None, name: str = 'amplitude mapper')[source]#

Amplitude indices to signed PAM symbols: the sign of PAS.

Signal Model#

A distribution matcher emits amplitudes; a PAM constellation is signed. The missing half is one bit per symbol,

\[y[n] = \varepsilon[n] \, a_{x[n]}, \qquad \Pr(\varepsilon[n] = \pm 1) = \tfrac{1}{2}\]

and it is free: an equiprobable sign leaves the amplitude distribution untouched, so the composite constellation carries the symmetric Maxwell-Boltzmann law

\[P_Y(\pm a_i) = \tfrac{1}{2} P_A(a_i)\]

at the same energy, and adds exactly one bit per symbol to the rate. That is the whole reason probabilistic amplitude shaping works: the systematic FEC encoder’s parity bits, which are uniform, are spent on the signs, and only the amplitudes go through the matcher.

This block draws the signs itself, from its own generator, which is the same distribution a parity stream has. It is a stochastic block, so seed() reaches it (D6).

Axes: element-wise – applied pointwise, shape-agnostic.

param alphabet:: Non-negative amplitudes \(a_i\), indexed by the input.
type alphabet:: np.ndarray
param seed:: Local RNG seed.
type seed:: int, optional, keyword-only
param name:: Name of the block. Default "amplitude mapper".
type name:: str, optional, keyword-only
raises ValueError:: If an amplitude is negative – the sign is this block’s business, and a signed input would make the constellation asymmetric.

References

Examples

>>> mapper = AmplitudeMapper(np.array([1.0, 3.0, 5.0, 7.0]), seed=0)
>>> print(mapper(np.array([0, 1, 0, 3, 2])))
[-1. -3. -1.  7.  5.]
>>> print(AmplitudeDemapper(mapper.alphabet)(mapper(np.array([0, 3]))))
[0 3]

class AmplitudeDemapper(alphabet: ndarray, *, name: str = 'amplitude demapper')[source]#

Signed PAM samples back to amplitude indices.

Signal Model#

The inverse of AmplitudeMapper on the half that carries the shaped data:

\[\hat{x}[n] = \arg\min_i \left| \left|\Re e\, y[n]\right| - a_i \right|\]

The real part is taken first because the amplitudes live on one quadrature lane of the constellation – each lane of a QAM symbol is demapped on its own. On that lane, deciding on the magnitude is not a shortcut but the maximum likelihood decision: the ASK points are symmetric about the origin, so the nearest point always has the amplitude nearest to \(|\Re e\, y|\). The sign is dropped because in PAS it carries parity, not data – a real receiver feeds it to the FEC decoder instead.

Axes: element-wise – applied pointwise, shape-agnostic.

param alphabet:: The same non-negative amplitudes the mapper was given.
type alphabet:: np.ndarray
param name:: Name of the block. Default "amplitude demapper".
type name:: str, optional, keyword-only
raises ValueError:: If an amplitude is negative.

References

G. Böcherer, F. Steiner, P. Schulte, IEEE Trans. Commun. 63(12), 2015, Section IV.

Examples

>>> demapper = AmplitudeDemapper(np.array([1.0, 3.0, 5.0, 7.0]))
>>> print(demapper(np.array([-0.8, 3.4, -5.2, 6.2])))
[0 1 2 3]