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

1.1. 指数表および対数表を作成する

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

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

const GF_LOG_NEGATIVE_INFINITY = 255;

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

	const exponents = new Uint8Array(256);

	exponents[0] = 1;

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

	exponents[GF_LOG_NEGATIVE_INFINITY] = 0;

	return exponents;

};

const GF_EXP = generateExponents();

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

	const logarithms = new Uint8Array(256);

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

	return logarithms;

};

const GF_LOG = generateLogarithms(GF_EXP);

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

実行結果

α^0   = 0b00000001 = 1
α^1   = 0b00000010 = 2
α^2   = 0b00000100 = 4
α^3   = 0b00001000 = 8
α^4   = 0b00010000 = 16
α^5   = 0b00100000 = 32
α^6   = 0b01000000 = 64
α^7   = 0b10000000 = 128
α^8   = 0b00011101 = 29
α^9   = 0b00111010 = 58

(中略)

α^253 = 0b01000111 = 71
α^254 = 0b10001110 = 142
α^-∞ = 0b00000000 = 0

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

ここで、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 です。

ここでは

\begin{align}
\alpha^0 &= 1 \\
\alpha^i &= (\alpha^{i-1} \cdot \alpha) \bmod (\alpha^8 + \alpha^4 + \alpha^3 + \alpha^2 + 1)
\end{align}

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

原始多項式 $x^8 + x^4 + x^3 + x^2 + 1$ の係数を 2 進数で表すと 100011101 になります。

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]) % 255]
- $x = 0 \lor y = 0$ のとき 0
割り算 $x \div y$
- $x \neq 0 \land y \neq 0$ のとき GF_EXP[(255 + GF_LOG[x] - GF_LOG[y]) % 255]
- $x = 0$ のとき 0

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

ここでは様々な誤り訂正符号の長さで符号化する想定で、各誤り訂正符号の長さ (eccLength) ごとに生成多項式を計算して生成多項式の一覧を作成します。

// 生成多項式の一覧を作成
// 
// 意味的には
// G_{eccLength}(x) = (x - α^0) (x - α) (x - α^2) ... (x - α^{eccLength-1})
const generateGenerators = () => {

	const generators = new Map()

	const generator = new Uint8Array(255).fill(GF_LOG_NEGATIVE_INFINITY);

	// eccLength === 0 のとき
	// generator === [0]
	// 意味的には G_0(x) = α^0 = 1
	generator[0] = 0;

	for (let eccLength = 1; eccLength < 255; eccLength++) { 

		// 意味的には
		// G_{eccLength}(x) = G_{eccLength-1}(x) * (x + α^{eccLength-1})
		//                  = G_{eccLength-1}(x) * x + G_{eccLength-1}(x) * α^{eccLength-1}
		for (let i = eccLength; i > 0; i--) {

			// GF_EXP[generator[i - 1]] === 0 のとき
			if ( generator[i - 1] === GF_LOG_NEGATIVE_INFINITY ) continue;

			const a = GF_EXP[generator[i]];
			const b = GF_EXP[(generator[i - 1] + eccLength - 1) % 255];
			generator[i] = GF_LOG[a ^ b];

		}

		generators.set(
			eccLength,
			generator.slice(0, eccLength + 1),
		);

	}

	return generators;

};

const RS_GENERATORS = generateGenerators();

// TODO: 
console.log(RS_GENERATORS);

ここでは生成多項式は

G_{d'}(x) = \prod_{i=0}^{d'-1} (x - \alpha^i)

より

\begin{align}
G_0(x) &= 1 \\
G_{d'}(x) &= G_{d'-1}(x) \cdot (x + \alpha^{d'-1}) \\
&= G_{d'-1}(x) \cdot x + G_{d'-1}(x) \cdot \alpha^{d'-1}
\end{align}

のようにして計算を行います ($G_0(x)$ は符号化には使いません) 。

生成多項式は符号化時に掛け算の計算から始めるため、あらかじめべき乗 $\alpha^i$ の指数 $i$ の状態で扱うことにします。

例えば、生成多項式 $G_4(x) = \alpha^0 x^4 + \alpha^{75} x^3 + \alpha^{249} x^2 + \alpha^{78} x + \alpha^6$ は 0, 75, 249, 78, 6 で表します。

3. 符号化

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

/**
 * 誤り訂正符号を計算する
 */
const rsECC = (data, eccLength) => {

	const generator = RS_GENERATORS.get(eccLength);

	const remainder = new Uint8Array(data.length + generator.length - 1);
	remainder.set(data);

	for (let offset = 0; offset < data.length; offset++) {

		// GF_EXP[generator[0]] === 1 より
		// 意味的には remainder[offset] / GF_EXP[generator[0]] === remainder[offset]
		// すなわち
		// generator[0] === 0 より
		// 意味的には (255 + GF_LOG[remainder[offset]] - generator[0]) mod 255 === GF_LOG[remainder[offset]]
		const quotient = GF_LOG[remainder[offset]];

		// remainder[offset] === 0 のとき
		if ( quotient === GF_LOG_NEGATIVE_INFINITY ) continue;

		for (let i = 0; i < generator.length; i++) {

			// GF_EXP[generator[0]] === 0 のとき
			if ( generator[i] === GF_LOG_NEGATIVE_INFINITY ) continue;

			// 意味的には remainder[offset + i] -= GF_EXP[quotient] * GF_EXP[generator[i]]
			remainder[offset + i] ^= GF_EXP[(quotient + generator[i]) % 255];

		}

	}

	return remainder.slice(data.length);

};

例

// 
const data = new Uint8Array([16, 240, 80]);
const ecc = rsECC(data, 4);

const code = new Uint8Array(data.length + ecc.length);
code.set(data);
code.set(ecc, data.length);

// 
console.log(data.join(', '));
console.log(ecc.join(', '));

console.log(code.join(', '));

例の実行結果

16, 240, 80
14, 177, 166, 169
16, 240, 80, 14, 177, 166, 169

生成多項式 $G_{d'}(x)$ の $x^{d'}$ の係数は常に $\alpha^0 = 1$ なため、上記の例では

\begin{alignat}{4}
(16 x^6 + 240 x^5 &+ & 80 x^4) &- 16 x^2 & G_4(x) &= 23 x^4 + 211 x^3 + 116 x^2 & \\
(23 x^4 + 211 x^3 &+ & 116 x^2) &- 23 & G_4(x) &= 14 x^3 + 177 x^2 + 166 x & + 169
\end{alignat}

のようにして割り算の余り

(16 x^6 + 240 x^5 + 80 x^4) \bmod G_4(x) = 14 x^3 + 177 x^2 + 166 x + 169

を計算しています。

4. 誤り検出

符号多項式 $C(x)$ に $G_{d'}(x) = 0$ となる $x = \alpha^i$ ($i = 0, 1, \dots , d' - 2, d' - 1$) を代入することでシンドロームを計算します。

シンドロームが全て 0 ならば誤りを含まず、シンドロームが 1 つでも 0 でないなら誤りを含むと判定します。

/**
 * シンドロームを計算する
 */
const rsSyndromes = (code, eccLength) => {

	const syndromes = new Uint8Array(eccLength);

	// 生成多項式 G_{eccLength}(x) の根 α^i
	for (let i = 0; i < eccLength; i++) {
		// 意味的にはシンドローム S_i = C(α^i)
		syndromes[i] = code
			.map(antilog => GF_EXP[(GF_LOG[antilog] + i) % 255])
			.reduce((a, b) => a ^ b);
	}

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

	return syndromes;

};

/**
 * 誤りを検出する
 */
// 符号多項式 C(x) が正しいとき、C(x) は生成多項式 G(x) で割り切れる
// 符号多項式 C(x) が正しいとき、G(x) = 0 ならば C(x) = 0
// G_{eccLength}(x) = (x - α^0) (x - α) (x - α^2) ... (x - α^{eccLength-1}) より
// x = α^0, α, α^2 ... α^{eccLength-1} ならば G_{eccLength}(x) = 0
const rsIsCorrupted = (code, eccLength) => rsSyndromes(code, eccLength).some(syndrome => syndrome !== 0);

例

// 
const data = new Uint8Array([16, 240, 80]);
const ecc = rsECC(data, 4);

const code = new Uint8Array(data.length + ecc.length);
code.set(data);
code.set(ecc, data.length);

// 
console.log(code.join(', '));
console.log(rsIsCorrupted(code, 4));

code[1] ^= 1;

console.log(code.join(', '));
console.log(rsIsCorrupted(code, 4));

例の実行結果

16, 240, 80, 14, 177, 166, 169
0, 0, 0, 0
false
16, 241, 80, 14, 177, 166, 169
1, 2, 4, 8
true

5. 関連記事

参考「JavaScript で 2 元 BCH 符号: 符号化と総当たりによる誤り訂正 - Qiita」
参考「JavaScript で 2 元ゴレイ符号: 符号化と総当たりによる誤り訂正 - Qiita」