JavaScript で 2 元ゴレイ符号: 符号化と総当たりによる誤り訂正

本記事で説明する内容:

- いくつかの 2 元ゴレイ符号の例
    - (23, 12, 7) 完全 2 元ゴレイ符号
    - (24, 12, 8) 拡大 2 元ゴレイ符号
    - (18, 6, 8) 短縮 2 元ゴレイ符号
- 平方剰余符号としての 2 元ゴレイ符号の構成
- 符号化
- 総当たりによる誤り訂正

本記事で説明しない内容:

- 一般のゴレイ符号
- BCH 符号としての 2 元ゴレイ符号の構成
- その他の方法による誤り訂正

2 元ゴレイ符号自体の説明は別記事にしました。

参考「[2 元ゴレイ符号 入門: 符号化 - Qiita](https://qiita.com/kerupani129/items/58cf347ab5349d686bda)」

# 1. 有限体上の計算の準備

## 1.1. 指数表および対数表を作成する
有限体上の掛け算や割り算を効率良くするために、$\alpha^i$ の値を対応させる指数表および対数表を作成します。

ここでは、原始多項式を $x^{11} + x^2 + 1$ とする有限体 $GF(2^{11})$ を使用することにします。

```js
const GF_LOG_NEGATIVE_INFINITY = 2047;

// 指数表を作成
const generateExponents = () => {

	const exponents = new Uint16Array(2048);

	exponents[0] = 1;

	for (let log = 1; log < 2047; log++) {
		const antilog = exponents[log - 1] << 1;
		exponents[log] = (antilog >>> 11 !== 0 ? antilog ^ 0b100000000101 : antilog);
	}

	exponents[GF_LOG_NEGATIVE_INFINITY] = 0;

	return exponents;

};

const GF_EXP = generateExponents();

// 対数表を作成
const generateLogarithms = exponents => {

	const logarithms = new Uint16Array(2048);

	for (let i = 0; i <= 2047; i++) {
		logarithms[exponents[i]] = i;
	}

	return logarithms;

};

const GF_LOG = generateLogarithms(GF_EXP);

// TODO: 
console.log(Array.from(
	GF_EXP.filter((_, i) => 0 === i % 89),
	(antilog, i) => {
		const j = i * 89;
		const log = (j !== GF_LOG_NEGATIVE_INFINITY ? String(j).padEnd(4) : '-∞ ');
		const antilogBinary = antilog.toString(2).padStart(11, '0');
		return `α^${log} = 0b${antilogBinary} = ${antilog}`;
	},
).join('\n'));
```

```:実行結果
α^0    = 0b00000000001 = 1
α^89   = 0b00101000010 = 322
α^178  = 0b00010101110 = 174
α^267  = 0b10010001100 = 1164
α^356  = 0b10001111100 = 1148
α^445  = 0b00111100001 = 481
α^534  = 0b01001111101 = 637
α^623  = 0b11110010110 = 1942
α^712  = 0b11101011111 = 1887
α^801  = 0b10111101110 = 1518
α^890  = 0b10010000011 = 1155
α^979  = 0b00010100111 = 167
α^1068 = 0b11111011011 = 2011
α^1157 = 0b00110100010 = 418
α^1246 = 0b00100011001 = 281
α^1335 = 0b10111110101 = 1525
α^1424 = 0b00101111010 = 378
α^1513 = 0b11011000000 = 1728
α^1602 = 0b11011010011 = 1747
α^1691 = 0b00100111111 = 319
α^1780 = 0b01000101000 = 552
α^1869 = 0b11101010100 = 1876
α^1958 = 0b10000111101 = 1085
α^-∞  = 0b00000000000 = 0
```

JavaScript では負の無限を `-Infinity` で表すことができますが `Uint16Array` 等では扱えないため、ここでは他の数を負の無限として扱うように定数 (`GF_LOG_NEGATIVE_INFINITY = 2047`) を定義します。

ここで、`GF_EXP[i]` は $x = \alpha^i$ を表し、`GF_LOG[x]` は $i = \log_\alpha x$ を表します。ただし `GF_EXP[GF_LOG_NEGATIVE_INFINITY] === 0` および `GF_LOG[0] === GF_LOG_NEGATIVE_INFINITY` です。

ここでは

```math
\begin{align}
\alpha^0 &= 1 \\
\alpha^i &= (\alpha^{i-1} \cdot \alpha) \bmod (\alpha^{11} + \alpha^2 + 1)
\end{align}
```

のようにして計算を行います。

原始多項式 $x^{11} + x^2 + 1$ の係数を 2 進数で表すと `100000000101` になります。

## 1.2. 有限体上の四則演算

ここで使用する有限体の性質より、以下のようにして有限体上の四則演算を行えます:

- 足し算 $x + y$
    - `x ^ y`
- 引き算 $x - y$
    - `x ^ y`
- 掛け算 $x \times y$
    - $x \neq 0 \land y \neq 0$ のとき `GF_EXP[(GF_LOG[x] + GF_LOG[y]) % 2047]`
    - $x = 0 \lor y = 0$ のとき `0`
- 割り算 $x \div y$
    - $x \neq 0 \land y \neq 0$ のとき `GF_EXP[(15 + GF_LOG[x] - GF_LOG[y]) % 2047]`
    - $x = 0$ のとき `0`

# 2. 生成多項式を計算する

ゴレイ符号は BCH 符号としても構成できますが、ここではゴレイ符号を平方剰余符号として生成多項式を計算します。

## 2.1. (23, 12, 7) 完全 2 元ゴレイ符号

```js
// 
/**
 * 生成多項式を計算する
 */
const golayGenerator = () => {

	// 平方剰余を計算する
	const quadraticResidues = new Uint8Array(11);

	for (let i = 0; i < 11; i++) {
		const x = i + 1;
		quadraticResidues[i] = x * x % 23;
	}

	// TODO: 
	console.log(quadraticResidues.join(', '));

	// 
	const generatorTypedArray = new Uint16Array(quadraticResidues.length + 1);

	// 意味的には G_{-1}(x) = 1
	generatorTypedArray[0] = 1;

	for (let j = 0; j < quadraticResidues.length; j++) {

		// 意味的には
		// β = α^{89}
		// G_j(x) = G_{j-1}(x) * (x + β^{quadraticResidues[j]})
		//        = G_{j-1}(x) * x + G_{j-1}(x) * β^{quadraticResidues[j]}
		for (let k = j + 1; k > 0; k--) {

			if ( generatorTypedArray[k - 1] === 0 ) continue;

			generatorTypedArray[k] ^= GF_EXP[(
				GF_LOG[generatorTypedArray[k - 1]] + 89 * quadraticResidues[j]
			) % 2047];

		}

	}

	// 意味的には
	// G(x) = G_{quadraticResidues.length-1}(x)

	// 係数が 0 または 1 か確認
	if ( generatorTypedArray.some(coefficient => 0 !== coefficient && 1 !== coefficient) ) {
		throw new Error('Invalid generator');
	}

	// 配列を 2 進数として単一の数に変換
	const generator = generatorTypedArray.reduce((binary, coefficient) => binary << 1 ^ coefficient);

	return generator;

};
```

```js:例
// 
const generator = golayGenerator();

console.log(`0b${generator.toString(2)}`);

// 
const eccLength = 31 - Math.clz32(generator);

const n = 23;
const k = n - eccLength;

console.log(`(${n}, ${k})`);
```

```text:例の実行結果
1, 4, 9, 16, 2, 13, 3, 18, 12, 8, 6
0b101011100011
(23, 12)
```

生成多項式を計算するためにまず数 $x^2 \bmod 23 = q$ となる整数 $x$ が存在する平方剰余 $q$ を求めます。

ここでは

```math
\begin{align}
(23 \pm x)^2 \bmod 23 &= (23^2 \pm 2 \cdot 23 x + x^2) \bmod 23 \\
&= x^2 \bmod 23
\end{align}
```

より、最大 11 つの平方剰余を計算すればいいことになります (詳しくは前述の別記事参照) 。

実際に計算すると平方剰余は

```math
q \in Q = \{1, 2, 3, 4, 6, 8, 9, 12, 13, 16, 18\}
```

より、生成多項式は

```math
\begin{align}
G_Q(x) &= \prod_{i \in Q} (x - \beta^i) \\
&= x^{11} + x^9 + x^7 + x^6 + x^5 + x + 1
\end{align}
```

となります。

ここでは JavaScript で計算する都合上

```math
\begin{alignat}{2}
&G_{Q,-1} & (x) &= 1 \\
&G_{Q,0} & (x) &= (x + \beta^1) \\
&G_{Q,1} & (x) &= (x + \beta^1) (x + \beta^2) \\
&G_{Q,2} & (x) &= (x + \beta^1) (x + \beta^2) (x + \beta^3) \\
&G_{Q,3} & (x) &= (x + \beta^1) (x + \beta^2) (x + \beta^3) (x + \beta^4) \\
& & &\hspace{0.5em} \vdots \\
&G_{Q,10} & (x) &= (x + \beta^1) (x + \beta^2) (x + \beta^3) (x + \beta^4) (x + \beta^6) (x + \beta^8) \\
& & &\phantom{{} = {}} \qquad (x + \beta^9) (x + \beta^{12}) (x + \beta^{13}) (x + \beta^{16}) (x + \beta^{18}) \\
& & &= G_Q(x)
\end{alignat}
```

のように計算します。

ここでは生成多項式の係数を 2 進数で表し、生成多項式

```math
\begin{align}
G_Q(x) &= x^{11} + x^9 + x^7 + x^6 + x^5 + x + 1 \\
&= 1 x^{11} + 0 x^{10} + 1 x^9 + 0 x^8 + 1 x^7 + 1 x^6 + 1 x^5 + 0 x^4 + 0 x^3 + 0 x^2 + 1 x + 1
\end{align}
```

は `0b101011100011` になります。

0 または 1 のみを含む配列を 2 進数として単一の数に変換するのは以下のように行えます。

```js
(0b000000000001, 0) => 0b000000000001 << 1 ^ 0; // 0b000000000010
(0b000000000010, 1) => 0b000000000010 << 1 ^ 1; // 0b000000000101
(0b000000000101, 0) => 0b000000000101 << 1 ^ 0; // 0b000000001010
(0b000000001010, 1) => 0b000000001010 << 1 ^ 1; // 0b000000010101
(0b000000010101, 1) => 0b000000010101 << 1 ^ 1; // 0b000000101011
(0b000000101011, 1) => 0b000000101011 << 1 ^ 1; // 0b000001010111
(0b000001010111, 0) => 0b000001010111 << 1 ^ 0; // 0b000010101110
(0b000010101110, 0) => 0b000010101110 << 1 ^ 0; // 0b000101011100
(0b000101011100, 0) => 0b000101011100 << 1 ^ 0; // 0b001010111000
(0b001010111000, 1) => 0b001010111000 << 1 ^ 1; // 0b010101110001
(0b010101110001, 1) => 0b010101110001 << 1 ^ 1; // 0b101011100011
```

ここでは符号長は $n = 23$ であり、生成多項式 $G_Q(x)$ は 11 次多項式よりデータ長は

```math
\begin{align}
k &= 23 - 11 \\
&= 12
\end{align}
```

になります。

## 2.2. (24, 12, 8) 拡大 2 元ゴレイ符号

(23, 12, 7) 完全 2 元ゴレイ符号の生成多項式に $x + 1$ をかけることで (24, 12, 8) 拡大 2 元ゴレイ符号の生成多項式を計算できます。

```js:例
// 意味的には
// generatorExtended = (x + 1) generator
const generator = 0b101011100011;
const generatorExtended = generator << 1 ^ generator;

console.log(`0b${generatorExtended.toString(2)}`);

// 
const eccLength = 31 - Math.clz32(generatorExtended);

const n = 24;
const k = n - eccLength;

console.log(`(${n}, ${k})`);
```

```text:例の実行結果
0b1111100100101
(24, 12)
```

実際に計算すると生成多項式は

```math
\begin{align}
G'_Q(x) &= (x + 1) G_Q(x) \\
&= (x + 1) (x^{11} + x^9 + x^7 + x^6 + x^5 + x + 1) \\
&= x^{12} + x^{11} + x^{10} + x^9 + x^8 + x^5 + x^2 + 1
\end{align}
```

となります。

ここでは生成多項式の係数を 2 進数で表し、生成多項式

```math
\begin{align}
G'_Q(x) &= x^{12} + x^{11} + x^{10} + x^9 + x^8 + x^5 + x^2 + 1 \\
&= 1 x^{12} + 1 x^{11} + 1 x^{10} + 1 x^9 + 1 x^8 + 0 x^7 + 0 x^6 + 1 x^5 + 0 x^4 + 0 x^3 + 1 x^2 + 0 x + 1
\end{align}
```

は `0b1111100100101` になります。

ここでは符号長は $n = 24$ であり、生成多項式 $G'_Q(x)$ は 12 次多項式よりデータ長は

```math
\begin{align}
k &= 24 - 12 \\
&= 12
\end{align}
```

になります。

## 2.3. (18, 6, 8) 短縮 2 元ゴレイ符号

データ符号の先頭 6 ビットを除去することで、(24, 12, 8) 拡大 2 元ゴレイ符号を短縮して (18, 6, 8) 短縮 2 元ゴレイ符号を構成できます。

```js:例
// 
const generatorExtended = 0b1111100100101;

console.log(`0b${generatorExtended.toString(2)}`);

// 
const eccLength = 31 - Math.clz32(generatorExtended);

const n = 18;
const k = n - eccLength;

console.log(`(${n}, ${k})`);
```

```text:例の実行結果
0b1111100100101
(18, 6)
```

符号長は異なりますが、生成多項式は同じです。

## 2.4. まとめ

|$(n, k)$|生成多項式|最小ハミング距離|誤り訂正能力|誤り訂正符号の長さ|
|-|-|-|-|-|
|$(23, 12)$|`0b101011100011`|7|3|11|
|$(24, 12)$|`0b1111100100101`|8|3|12|
|$(18, 6)$|`0b1111100100101`|8|3|12|

# 3. 符号化

有限体上で情報多項式を生成多項式で割った余りを計算すると、その余りの多項式の係数が誤り訂正符号になります。

```js
// 
const encodeBinaryCyclic = (n, k, generator) => {

	const codeLastOffset = n - 1;
	const dataLastOffset = k - 1;
	const eccLength = n - k;

	return data => {

		const shifted = data << eccLength;

		// 有限体 GF(2^m) 上で多項式の割り算の余りを計算する
		let remainder = shifted;

		for (
			let i = codeLastOffset, subtrahend = generator << dataLastOffset;
			i >= eccLength;
			i--, subtrahend >>>= 1
		) {
			if ( remainder >>> i !== 0 ) {
				remainder ^= subtrahend;
			}
		}

		// data と remainder を合わせて符号とする
		// shifted ^ remainder === data << eccLength ^ remainder
		const code = shifted ^ remainder;

		return code;

	};

};
```

2 元ゴレイ符号の場合は、以下のようにして符号化できます。

```js:例
// 
const encodePerfectBinaryGolay = encodeBinaryCyclic(23, 12, 0b101011100011);
const encodeExtendedBinaryGolay = encodeBinaryCyclic(24, 12, 0b1111100100101);
const encodeShortenedBinaryGolay = encodeBinaryCyclic(18, 6, 0b1111100100101);

// 
const data = 0b000000000111;
const dataShortened = 0b000111;
const codePerfect = encodePerfectBinaryGolay(data);
const codeExtended = encodeExtendedBinaryGolay(data);
const codeShortened = encodeShortenedBinaryGolay(dataShortened);

console.log(`0b${codePerfect.toString(2).padStart(23, '0')}`);
console.log(`0b${codeExtended.toString(2).padStart(24, '0')}`);
console.log(`      0b${codeShortened.toString(2).padStart(18, '0')}`);
```

```text:例の実行結果
0b00000000011111001001010
0b000000000111110010010100
      0b000111110010010100
```

上記の例では、以下のようにして割り算の余りを計算しています。

```js
0b000000000111_00000000000  ^ 0b1_01011100011  << 2 === 0b000000000010_01110001100;
0b000000000010_01110001100  ^ 0b1_01011100011  << 1 === 0b000000000000_11001001010;

0b000000000111_000000000000 ^ 0b1_111100100101 << 2 === 0b000000000000_110010010100;

      0b000111_000000000000 ^ 0b1_111100100101 << 2 ===       0b000000_110010010100;
```

# 4. 総当たりによる誤り訂正

ゴレイ符号の誤り訂正はいくつかの方法がありますが、データ長が短い場合は総当たりで誤り訂正する方が計算量を減らせます。

(18, 6, 8) 短縮 2 元ゴレイ符号は以下のようにして誤り訂正できます。

```js
// 
const popCount = x => {
	const a = x - (x >>> 1 & 0x55555555);
	const b = (a & 0x33333333) + (a >>> 2 & 0x33333333);
	const c = b + (b >>> 4) & 0x0f0f0f0f;
	const y = c * 0x01010101 >>> 24;
	return y;
};

const correctBinaryCyclic = (t, codes) => code => (
	codes.find(codeCorrected => popCount(code ^ codeCorrected) <= t)
);
```

```js:例
// 
const generateCodes = () => {

	const codes = new Uint32Array(64);

	for (let data = 0; data < 64; data++) {
		codes[data] = encodeShortenedBinaryGolay(data);
	}

	return codes;

};

const correctBinaryGolay = correctBinaryCyclic(3, generateCodes());

// 
const dataShortened = 0b000111;
const codeShortened = encodeShortenedBinaryGolay(dataShortened);

console.log(`0b${codeShortened.toString(2).padStart(18, '0')}`);

// 
const codeCorrupted = codeShortened ^ 0b100000100000100000;
const codeCorrected = correctBinaryGolay(codeCorrupted);

console.log(`0b${codeCorrupted.toString(2).padStart(18, '0')}`);
console.log(`0b${codeCorrected.toString(2).padStart(18, '0')}`);
```

```text:例の実行結果
0b000111110010010100
0b100111010010110100
0b000111110010010100
```

正しい符号の一覧を作成し、受け取った符号とのハミング距離が誤り訂正能力 (ここでは 3) 以下の正しい符号を探します。

2 つの値のハミング距離は、2 つの値の排他的論理和のハミング重み (ビットカウント) と等しいです。

参考「[JavaScript でなるべく速くビットカウントする - Qiita](https://qiita.com/kerupani129/items/e14c2f3111ae09ecf7cf)」

# 5. 関連記事

参考「[JavaScript で 2 元 BCH 符号: 符号化と総当たりによる誤り訂正 - Qiita](https://qiita.com/kerupani129/items/3cb3237eccca654007d3)」
参考「[JavaScript でリード・ソロモン符号: 符号化と誤り検出 - Qiita](https://qiita.com/kerupani129/items/a9858a12b987e829f768)」