Soft Demapping of Spherical Codes from
Cartesian Powers of PAM Constellations

We denote the set of the first $n$ positive integers by $[n]$. A permutation of length $n$ is a bijection from $[n]$ to itself. The set of all permutations of length $n$ forms a group under function composition. This group is called the symmetric group on $[n]$ and is denoted by $\mathcal{S}_{n}$.

There are different ways to represent a permutation $\alpha\in\mathcal{S}_{n}$. In the two-line form, the permutation $\alpha$ is given as

The matrix form of $\alpha$ is given as

where $\bm{e}_{j}$ denotes a column vector of length $n$ with $1$ in the $j^{\text{th}}$ position and $0$ in every other position. The set of permutation matrices $\{\mathsf{P}_{\alpha}\mid\alpha\in\mathcal{S}_{n}\}$ forms a group under matrix multiplication and is usually denoted as $\mathcal{A}_{n-1}$. This group is, of course, isomorphic to $\mathcal{S}_{n}$. A signed permutation matrix is a generalized permutation matrix whose nonzero entries are $\pm 1$. The set of all $n\times n$ signed permutation matrices forms a group under matrix multiplication and is usually denoted as $\mathcal{B}_{n}$.

A Variant I permutation code $\mathcal{C}$ of length $n$ [2] is characterized as the set of points in real Euclidean $n$-space $\mathbb{R}^{n}$ that are generated by the group $\mathcal{A}_{n-1}$ acting on a vector $\bm{x}$:

This means that the code $\mathcal{C}$ is the set of all distinct vectors that can be formed by permuting the order of the $n$ elements of a vector $\bm{x}$. Each vector in $\mathcal{C}$ is called a codeword. Likewise, a Variant II permutation code $\mathcal{C}$ of length $n$ is characterized as the set of points in $\mathbb{R}^{n}$ that are generated by the group $\mathcal{B}_{n}$ acting on a vector $\bm{x}$:

Here, the code $\mathcal{C}$ is the set of all distinct vectors that can be formed by permuting the order or changing the sign of the $n$ elements of a vector $\bm{x}$.

The initial vector of a Variant I or II permutation code is the codeword

that satisfies $0<x_{1}\leq x_{2}\leq x_{3}\leq\dots\leq x_{n}$. The left most inequality is strict to be consistent with the use of a PAM constituent constellation.

We consider the communication system of PAS as described in [3]. As shown in Fig. 1, at the first stage in the transmitter, a sequence of binary information inputs is mapped to a Variant I permutation code whose initial vector has an empirical distribution that resembles a desired input distribution. This first step is called permutation encoding and is described further in Section IV. Each symbol in the constituent constellation is labeled according to some Gray labeling of the used PAM. Therefore, corresponding to each Variant I codeword, there is a binary label obtained by concatenating the Gray label of each of the elements in the codeword. If the constituent constellation²²2We assume that the constituent constellation has a power-of-two size. is $\mathcal{M}$, the length of the label of each Variant I codeword of blocklength $n$ is exactly $n\log_{2}\lvert\mathcal{M}\rvert$. The binary label of a codeword should not be confused with the information bits at the input of the permutation encoder that produced the codeword. Each Variant I permutation codeword is the output of the permutation encoder applied to a sequence of typically less than $n\log_{2}\lvert\mathcal{M}\rvert$ input information bits.

To protect the labels, a soft-decodable channel code is used. We emphasize that instead of encoding the input information bits that was used by the permutation encoder, the codeword labels are encoded by the channel code. This is fundamental in the development of PAS that uses bit-interleaved coded modulation as it allows the permutation codewords to have a semi-Gray labeling: if codewords are close in Euclidean distance, their labels are also close in the Hamming distance. If the input information bits are used to label the codewords, a semi-Gray labeling, even for small blocklengths, appears to be challenging to find. The channel encoder produces a sequence of parity bits which is used to assign signs to the elements of the Variant I permutation codeword³³3Sometimes, some binary inputs are directly used as sign bits of some of the elements of the permutation codeword. This is done when the ratio of the blocklength of the permutation code to the length of its binary label, $1/\log_{2}\lvert\mathcal{M}\rvert$, is larger than the overhead of the channel code. See [3] for details.. This will produce a Variant II permutation codeword which will be transmitted over an AWGN channel with noise variance $\sigma^{2}$.

A noisy version of each transmitted Variant II permutation codeword is received at the receiver. By demapping, we mean the detection of the binary label and sign bits of a permutation codeword from its received noisy version. The first block at the receiver estimates the log-likelihood ratio (LLR) for each bit in the label of the received word—a process we refer to as soft demapping. The estimated LLRs are then passed to the decoder of the channel code to approximate the label of the transmitted Variant I permutation codeword at the receiver. If an acceptable Variant I codeword is obtained, the inverse of permutation encoder produces the corresponding binary input sequence. Conversely, if the reproduced symbols do not form a valid Variant I permutation codeword, a block error is declared.

A bipartite graph $G=(\mathcal{U},\mathcal{V},\mathcal{E})$ is a triple consisting of two disjoint set of vertices $\mathcal{U}$ and $\mathcal{V}$ and a set of edges $\mathcal{E}\subseteq\mathcal{U}\times\mathcal{V}$. Each edge in $\mathcal{E}$ connects a vertex in $\mathcal{U}$ to a vertex in $\mathcal{V}$. The graph $G$ is called balanced if $\lvert\mathcal{U}\rvert=\lvert\mathcal{V}\rvert=n$. The bipartite graph $G$ is complete if $\mathcal{E}=\mathcal{U}\times\mathcal{V}$. All bipartite graphs that we consider are assumed complete and balanced. The graph $G$ is weighted if there exist a weight function $w:\mathcal{E}\to\mathbb{R}\cup\{-\infty\}$. A matching is a set of edges without common vertices. The weight of a matching is the sum of the weight of the edges in that matching. A perfect matching is a matching that matches all vertices of the graph. The maximum weight perfect matching problem is the problem of finding a perfect matching with maximum weight. This problem is usually called the assignment problem and is stated as follows.

The assignment problem is typically solved using the Hungarian method [13, 14]. The time complexity of the Hungarian algorithm is $O\left(n^{3}\right)$. However, the weight function may have additional structure that can be exploited to reduce the complexity of solving the assignment problem.

A weight function $w$ is said to be multiplicative if there exist two functions

such that

If the weight function is multiplicative the assignment problem can be solved using sorting as follows: Assume

and that there are two permutations $\sigma_{\mathcal{U}}$ and $\sigma_{\mathcal{V}}$ on $[n]$ such that

In this case, the assignment problem can be immediately solved using the upper bound in the rearrangement inequality [15]. The bijection $f$ that solves the assignment problem is defined by

The solution can be found by sorting in $O(n\log n)$ time.

We consider a vector Gaussian noise channel

where the input $\bm{X}$ is a random vector uniformly distributed over a permutation code $\mathcal{C}$, the noise term $\bm{Z}$ is a standard normal random vector and $\sigma$ is the standard deviation of the additive white Gaussian noise. The output $\bm{Y}$ is a noisy version of $\bm{X}$. The joint distribution⁴⁴4Instead of distinguishing between discrete, continuous or mixed random vectors, we always refer to $p$ as a distribution for convenience. of the input and the output is denoted by $p(\bm{x},\bm{y})$. Similarly, $p(\bm{y}\mid\bm{x})$ is the conditional distribution of the output conditional on the input and $p(\bm{y})$ is the output distribution. List decoding of a Variant I permutation code is to find a set $\mathcal{L}=\{\bm{c}_{0},\bm{c}_{1},\dots,\bm{c}_{L-1}\}$ of codewords with the highest likelihoods so that

where $L=\lvert\mathcal{L}\rvert$ is the list size.

To be able to find such a list of candidate codewords, we need to define a weight function for the maximum likelihood decoding of permutation codes. Because the noise is Gaussian, a codeword $\bm{x}$ has a higher likelihood than a codeword $\bm{x}^{\prime}$ if $\bm{x}$ is closer to the received word $\bm{y}$ than $\bm{x}^{\prime}$. That is,

Using the constant-energy property of permutation codes, we can simplify further to get

Therefore, list decoding with maximum likelihood reward is equivalent to list decoding with maximum correlation reward. For Variant I permutation codes, we can define the following assignment problem:

This is an assignment problem with

and weight function

The weight function is clearly multiplicative with $w_{\mathcal{U}}(x_{i})=x_{i}$ and $w_{\mathcal{V}}(y_{j})=y_{j}$. The list decoding problem, then, can be solved by a direct application of Murty’s algorithm as described in the Appendix.

The assignment problem defined for Variant I permutation codes is not suitable for Variant II permutation codes as in a Variant II code, each symbol can have a plus or minus sign. It is not possible to accommodate such arbitrary sign flips in the assignment problem defined for Variant I permutation codes. In what follows, we introduce maximum orbit-likelihood decoding of Variant II permutation codes that allows us to reformulate the decoding problem in the framework of an assignment problem.

A trivial subcode of a Variant II permutation code $\mathcal{C}$ corresponding to a codeword $\bm{x}\in\mathcal{C}$ is the subcode containing $\bm{x}$ and all codewords obtained by changing the sign of elements in $\bm{x}$ in all possible ways. We denote the group of $n$-dimensional diagonal unitary matrices (i.e., $n\times n$ diagonal integer matrices with $\pm 1$ diagonal entries) as $\mathcal{D}_{n}$. A trivial subcode corresponding to $\bm{x}$ is the orbit of $\bm{x}$ under the group action $\phi$ defined by

Similar to finding the codeword in a Variant I permutation code with the maximum likelihood, we can think about finding the orbit in a Variant II permutation code with the maximum orbit likelihood—a process that we call maximum orbit-likelihood decoding. Also, similar to a list of $L$ most likely codewords for a Variant I permutation code, we can think about finding a list of $L$ most likely orbits for a Variant II permutation code—a process that we refer to as orbit list-decoding. In order to use Murty’s algorithm, a weight function for the orbits needs to be defined. We denote the likelihood of the orbit of $\bm{x}\in\mathcal{C}$ by

The key result of this section is the following theorem.

Theorem 1 suggests the weight function

for Murty’s algorithm in orbit list-decoding of Variant II permutation codes. This weight function, however, is not multiplicative. Using a generalized version of the rearrangement inequality [16], it turns out that this weight function has all the required properties of a multiplicative weight function that we used to simplify Murty’s algorithm in the Appendix.

Mutual information of Variant II permutation codes can be computed using the orbit likelihoods by the method discussed in the Appendix. Note that the contribution of the likelihood of one orbit in $p(\bm{y})$ in (7) is equivalent to the contribution of all $2^{n}$ codewords in the orbit.

A Variant I permutation code $\mathcal{C}$ is uniquely identified by its initial vector $\bm{x}$. Note that when the length of the code is $n$, the number of distinct elements in $\bm{x}$ is at most $n$. That is,

We denote the number of unique elements of $\bm{x}$ by $u$. We define $\bm{\mu}$ as the vector whose elements are the distinct elements of $\bm{x}$ in ascending order. That is, if

then $\mu_{1}<\mu_{2}<\mu_{3}<\dots<\mu_{u}$. The number of times each $\mu_{i}$ occurs in $\bm{x}$ is denoted as $m_{i}$. This means that

The size of this code is

and its rate is $r=\frac{\log_{2}\lvert\mathcal{C}\rvert}{n}$.

The problem of encoding permutation codes has been considered by many [19, 4, 20, 18]. Here, we adapt the encoder provided in [18]. A pseudocode for this algorithm is provided in Algorithm 1. A pseudocode for the inverse encoder is given in Algorithm 2. The computational complexity of both algorithms is $O(nu)$ arithmetic operations.

A Variant II permutation code $\mathcal{C}$, similar to a Variant I permutation code, is uniquely identified by its initial vector $\bm{x}$. Note that all elements of the initial vector $\bm{x}$ are positive. We use the same notation $u$ to denote the number of distinct elements of $\bm{x}$. Similarly, the vector $\bm{\mu}$ is defined as the vector whose elements are the distinct elements of $\bm{x}$ in ascending order. The number of times each $\mu_{i}$ occurs in $\bm{x}$ is denoted as $m_{i}$. The rate of this code is

Let $k=\lfloor nr\rfloor$. Encoding an integer in the range $0$ to $2^{k}-1$ can be done using an encoder for its constituent Variant I permutation code. The $k-n$ least significant bits are mapped to a codeword $\bm{c}$ in the constituent Variant I permutation code using an encoder as in Section IV-A. The signs of the symbols in the resulting codeword $\bm{c}$ are chosen using the $n$ most significant bits. If the encoder for the constituent Variant I permutation code is $\operatorname{enc}_{\mathrm{I}}$, the encoder for the Variant II permutation code will be the function

such that

where, base-$2$ representation of $q$ is $b_{n+1}b_{n+2}\dots b_{k}$. Similarly, an inverse encoder can be defined using the inverse encoder for the constituent Variant I permutation code.

We are now ready to explicitly describe the first block in the PAS system of Fig. 1. Assume we want to use a Variant I permutation code with initial vector $\bm{x}$ and blocklength $n$ whose rate is $r$. The first step for us is to expurgate the corresponding Variant I code so that its size is a power of $2$. Let $k_{\mathrm{a}}=\lfloor nr\rfloor$. To identify a sequence of amplitudes of length $n$, i.e., a Variant I codeword of blocklength $n$, one can use the encoder defined in Section IV-A. The input is an integer $i$ in the range $[0,2^{k_{\mathrm{a}}})$. Therefore, each codeword in this Variant I permutation code can be indexed by a binary vector of length $k_{\mathrm{a}}$. This way, if $nr$ is not an integer, the Variant I permutation code is expurgated so that only the first $2^{k_{\mathrm{a}}}$ codewords, based on the lexicographical ordering, are utilized. This can potentially introduce a considerable correlation between the elements of codewords and breaks the geometric uniformity of the permutation code [21].

In order to properly randomize so that each codeword is surrounded by approximately equal number of codewords at any given distance, we first try to randomly spread the set of integers in the range of $[0,2^{k_{\mathrm{a}}})$ on the integers in the range of $[0,2^{nr})$. One way to achieve this is to take the sampled integer $i$ from $[0,2^{k_{\mathrm{a}}})$ and multiply it by a randomly chosen integer $e$ coprime with $2^{nr}$ and add another randomly chosen integer $d$ in $[0,2^{k_{\mathrm{a}}})$ modulo $2^{nr}$. The coprime factor $e$ and the dither $d$ can be generated based on some shared randomness so that the receiver can also reproduce them synchronously. Consequently, the first block in Fig. 1 takes in an integer $i$ from $[0,2^{k_{\mathrm{a}}})$ at the input and produces a Variant I codeword $\bm{c}$ by

The codeword $\bm{c}$ serves as the set of amplitudes in PAS and is labeled based on the Gray labeling of the constituent PAM constellation. The labels are then encoded using a systematic soft-decodable channel code. The resulting parity bits are used to set the signs of each of the amplitudes and as a result, a Variant II permutation code is produced.

The input to the first block at the receiver in Fig. 1 contains noisy version of Variant II permutation codewords. Assume that the received word is $\bm{y}$ and that we want to estimate the LLR for the $i^{\mathrm{th}}$ bit, $\mathrm{LLR}_{i}$, which is a part of the label in the $j^{\mathrm{th}}$ symbol, $c_{j}$, of the Variant II codeword $\bm{c}$. Four different estimators are discussed in what follows.

If the $i^{\mathrm{th}}$ bit is a sign bit, we only take into account the contribution of the half of the codewords in each orbit that have positive sign in the numerator and the other half in the denominator. Using a technique similar to the proof of Theorem 1, one can show that for a sign bit

To find the corresponding orbit of each codeword on the list, the assignment problem needs to be solved with the frozen symbol as an inclusion constraint (See the Appendix). Once the received word is sorted by the magnitude of its elements, the orbit likelihood of each item on the list can be found using $O(n)$ operations (see Theorem 1). There are a total of $O(n\lvert\mathcal{M}\rvert)$ orbits to consider. Therefore, the orbit likelihoods can be computed with $O(n^{2}\lvert\mathcal{M}\rvert)$ complexity. Once the orbit likelihoods are computed, each of the $n\log_{2}\lvert\mathcal{M}\rvert$ LLRs can be found with $O(\lvert\mathcal{M}\rvert)$ complexity. Therefore, the computational complexity of this method is $O(n^{2}\lvert\mathcal{M}\rvert+n\lvert\mathcal{M}\rvert\log\lvert\mathcal{M}\rvert)$ operations per permutation codeword.

Consider the abstract channel for which the channel input is the label of the transmitted permutation codeword $\bm{B}=(B_{1},B_{2},\dots,B_{k_{\mathrm{a}}+n})$ and its output is the noisy permutation codeword $\bm{Y}=(Y_{1},Y_{2},\dots,Y_{n})$. The bit-metric decoding (BMD) rate⁵⁵5Here, $\left[x\right]^{+}$ is the same as $\max\{x,0\}$. [7]

is an achievable rate [24] and can be used to compare different methods of soft demapping. This is particularly useful as it provides an achievable rate when a binary soft-decodable channel code is intended to be used, without specifically restricting the system to a particular implementation of the channel code. The BMD rate of permutations codes, as well as generalizations discussed in the next section, for different methods of estimating LLRs are shown in Fig. 4.

To better compare the symbol-by-symbol method with our proposed method using orbit decoding as well as decoding using the BCJR algorithm, we use a permutation code in a PAS scheme. We compare the permutation block error rate (BLER) for Code 2 from Example 2 for these three estimation methods. By each block, we mean a permutation codeword (or, in general, a spherical codeword) and BLER is the ratio of blocks received in error after the channel decoder to the total number of blocks sent. For the channel code, we use the $(10\,860,8\,448)$ low density parity check (LDPC) code from the 5G NR standard [25]. Each LDPC code frame contains $67$ permutation codewords. Among the sign bits, $14$ are used for systematic information bits [3]. In Fig. 5, we show the BLER for these three LLR estimation methods. Note that we think of the case of using a permutation code with symbol-by-symbol decoding [3] as the baseline for comparison with other spherical codes and decoding methods that we study. For instance, it is clear that at BLER $=10^{-3}$, about $0.3$ dB improvement is achieved by using the orbit decoding with frozen symbols compared with the symbol-by-symbol method. Also, about $0.2$ dB improvement is achieved compared with the BCJR method. In Fig. 5, similar results are shown for more general spherical codes as well. The structure of these codes will be explained in Section V.

We design codes with blocklength $n$ that are shells of the product of $n$ $2p$-PAM constellations, constrained to have total energy $\mathtt{E}$ in each block. This is formalized as follows.

We refer to the complete shell code of length $n$ and squared radius $\mathtt{E}$ over the $2p$-PAM constellation as the $(n,\mathtt{E},p)$ code. Complete shell code automatically deliver a shaping gain, as the symmetry of the construction induces an approximate Boltzmann distribution over the symbols. That is, if codewords are selected uniformly, the probability distribution over the symbols $x_{i}$ is $p(x_{i})\propto\exp[\lambda x_{i}^{2}]$, where $\lambda$ is a function of the code parameters $n$ and $\mathtt{E}$. Note that this is a discretized version of a Gaussian distribution.

Complete shell codes can be decomposed into a union of permutation subcodes as $\mathcal{C}=\cup_{i=1}^{t(\mathcal{C})}\mathcal{C}^{i}$, where each $\mathcal{C}^{i}$ is a permutation code, and the number of distinct permutation subcodes $t(\mathcal{C})$ is a function of the code parameters. We will also write $t(\mathcal{C})=t(n,\mathtt{E},p)$ when $\mathcal{C}$ is the $(n,\mathtt{E},p)$ code. We refer to each subcode as a type class following the terminology of, e.g., [26, Ch. 11.1]. Describing the code in this union form is useful when describing and analyzing the code performance.

We call the initial vector in a type class the type class representative. For instance, the $(8,32,4)$ code has two type classes, with representatives $(1,1,1,1,1,3,3,3)$ and $(1,1,1,1,1,1,1,5)$. The first class contains $14\,336\approx 2^{13.8}$ codewords, and the second $2048=2^{11}$. In general, the number of type classes $t(n,\mathtt{E},p)$ of a given energy $\mathtt{E}$ is a complicated function, which we usually compute numerically. A loose bound for the maximum value of $t(n,\mathtt{E},p)$ when varying $n$ and $p$ is $O(n^{p})$. This becomes intractable quickly for larger values of $p$.

For large $n$, orbit decoding is impractical. Even BCJR decoding becomes unwieldy as $\mathtt{E}$ grows. Here we define a new family of codes that has competitive rate and shaping properties, with vastly improved encoding and decoding complexities. In particular, keeping only a small number of the largest type classes in a complete shell code is often sufficient to create a high-performing code. We formalize this as follows.

The maximal $k$-class partial shell code consists of the $k$ largest type classes. In many cases, relatively small values of $k$ result in codes whose rates are only marginally smaller than that of the complete shell code, with the benefit that the encoding and decoding algorithms now have a constant number of type classes, i.e., we have replaced the function $t(\mathcal{C})$ with the constant $k$. Note that $k=1$ corresponds to a permutation code.

The acceptable reduction in rate will depend on the use-case for the code, but we already see the complexity benefits of using selective codes instead of complete codes. Conversely, if we want to create a selective code from a union of permutation codes, we see that the rate increases for adding each subsequent code to the union decrease quickly.

As previously noted, one of the appealing features of complete codes over PAM constellations is that they induce an approximate marginal Boltzmann distribution. Partial shell codes also have this property. This is a result of the fact that they are constructed from the largest type classes available, whose symbols $m_{i}$ are distributed proportionally to $\exp[\lambda(2i-1)^{2}]n$, a Boltzmann distribution, where parameter $\lambda$ depends on $\mathtt{E}$. This result (more specifically its dual) is shown in [29].

The exact Boltzmann distribution does not give an integer solution in general, but the non-integer solution is often a near approximation for solutions. In fact, since we are often interested not only in the largest type class, but in the $k$ largest type classes, an alternative method for constructing $k$-class partial shell codes is to start with the non-integer Boltzmann solution and search for several nearby solutions, e.g., by iterating over all combinations of integer values within some neighborhood of each non-integer entry, and keeping the resulting vectors that are the $k$ best solutions.

The primary benefit of partial codes over complete codes is their improved complexity of decoding. Encoding can be done using Algorithm 3 with using $k$ type-classes.

Fig. 5 also shows the SNR values for the 4-maximal partial shell code. We see that they are intermediate between the results for the permutation and complete codes. It is important to note that, compared with orbit decoding, the complete shell code with symbol-by-symbol decoding performs very well and is less than 0.15 dB worse than decoding with the BCJR algorithm. In terms of computational complexity, this is much cheaper than the BCJR algorithm and, therefore, the complete shell code with symbol-by-symbol decoding is still practically competitive.