Кодирование Рида-Соломона для чайников

Есть способ передавать данные, теряя часть по пути, но так, чтобы потерянное можно было вернуть по прибытии. Это третья, завершающая часть моего простого изложения алгоритма избыточного кодирования по Риду-Соломону. Реализовать это в коде не прочитав первую, или хотя бы вторую часть на эту тему будет проблематично, но чтобы понять для себя что можно сделать с использованием кодировки Рида-Соломона, можно ограничиться прочтением этой статьи.

Что может этот код?

Итак, что из себя представляет избыточный код Рида-Соломона с практической точки зрения? Допустим, есть у нас сообщение – "DON'T PANIC". Если добавить к нему несколько избыточных байт, допустим 6 штук: "rrrrrrDON'T PANIC" (каждый r – это рассчитанный по алгоритму байт), а затем передать через какую-нибудь среду с помехами, или сохранить там, где данные могут понемногу портиться, то по окончании передачи или хранения у нас может остаться такое, например: "rrrrrrDON'AAAAAAA" (6 байт оказались с ошибкой). Если мы знаем номера байтов, где вместо букв, которые были при создании кода, вдруг оказались какие-нибудь "A", то мы можем полностью восстановить сообщение в исходное "rrrrrrDON'T PANIC". После этого можно для красоты убрать избыточные символы. Теперь текст можно печатать на обложку.

Вообще, избыточных символов к сообщению мы можем добавить сколько угодно. Количество избыточных символов равно количеству исправляемых ошибок (это верно лишь в том случае, когда нам известны номера позиций ошибок). Как правило, ошибки, положение которых известно, называют erasures. Благозвучного перевода найти не могу ("стирание" мне не кажется благозвучным), так что в дальнейшем я буду применять термин "опечатки" и ставить его в кавычки (прекрасно понимаю, что этот термин обычно несёт похожий, но другой смысл). Исправление "опечаток" полезно, например, при восстановлении блоков QR кода, которые по какой-либо причине не удалось прочитать.

Также код Рида-Соломона позволяет исправлять ошибки, положение которых неизвестно, но тогда на каждую одну исправляемую ошибку должно приходиться 2 избыточных символа. "rrrrrrDON'T PANIC", принятые как "rrrrrrDO___ PANIC" легко будут исправлены без дополнительной информации. Неправильно принятый байт, положение которого неизвестно, в дальнейшем я буду называть "ошибкой" и тоже брать в кавычки.

Можно комбинировать исправление "ошибок" и "опечаток". Если, например, есть 3 избыточных символа, то можно исправить одну "ошибку" и одну "опечатку". Ещё раз обращу внимание на то, что чтобы исправить "опечатку", нужно каким-то образом (не связанным с алгоритмом Рида-Соломона) узнать номер байта "опечатки". Что важно, и "ошибки" и "опечатки" могут быть исправлены алгоритмом и в избыточных байтах тоже.

Стоит отметить, что если количество переданных и принятых байт отличается, то здесь код Рида-Соломона практически бессилен. То есть, если на расшифровку попадёт такое: "rrrrrrDO'AIC", то ничего сделать не получится, если, конечно, неизвестно какие позиции у пропавших букв.

Как закодировать сообщение?

Здесь уже не обойтись без понимания арифметики с полиномами в полях Галуа. Ранее мы научились представлять сообщения в виде полиномов и проводить операции сложения, умножения и деления над ними. Уже этого почти достаточно, чтобы создать код Рида-Соломона из сообщения. Единственно, для того, чтобы это сделать понадобится ещё полином-генератор. Это результат такого произведения:

Где – это примитивный член поля (как правило, выбирают 2), а – это количество избыточных символов. То есть, прежде чем создавать код Рида-Соломона из сообщения, нужно определиться с количеством избыточных символов, которое мы считаем достаточным, затем перемножить биномы вида в количестве штук по правилам перемножения полиномов. Для любого сообщения можно использовать один и тот же полином-генератор, и любое сообщение в таком случае будет закодировано с одним и тем же количеством избыточных символов.

Пример: Мы решили использовать 4 избыточных символа, тогда нужно составить такое выражение:

Так как мы работаем с полем Галуа с характеристикой 2, то вместо минуса можно смело писать плюс, не боясь никаких последствий. Жаль, что это не работает с количеством денег после похода в магазин. Итак, возводим в степень, и перемножаем (по правилам поля Галуа GF[256], порождающий полином 285):

Выражение выше есть полином-генератор, который необходим для того, чтобы закодировать сообщение любой длины, добавив к нему 4 избыточных символа, которые позволят скорректировать 2 "ошибки" или 4 "опечатки".

Прежде чем приводить пример кодирования, нужно договориться об обозначениях. Полиномы, записанные "по-математически" с иксами и степенями выглядят довольно-таки громоздко. На самом деле, при написании программы достаточно знать коэффициенты полинома, а степени можно узнать из положения этих коэффициентов. Таким образом полученный в примере выше полином-генератор можно записать так: {116, 231, 216, 30, 1}. Также, для ещё большей компактности, можно опустить скобки и запятые и записать всё в шестнадцатеричном представлении: 74 E7 D8 1E 01. Выходит в 2 раза короче. Надо отметить, что если в "математической" записи мы не пишем члены, коэффициенты которых равны нулю, то при принятой здесь шестнадцатеричной записи они обязательны, и, например, нужно записывать так: или 00 00 00 00 0A. Там, где "математическая" запись позволит более понятно объяснить суть, я буду прибегать к ней.

Итак, чтобы представить сообщение «DON'T PANIC» в полиномиальной форме, с учётом соглашения выше достаточно просто записать его байты:

44 4F 4E 27 54 20 50 41 4E 49 43.

Чтобы создать код Рида-Соломона с 4 избыточными символами, сдвигаем полином вправо на 4 позиции (что эквивалентно умножению его на ):

00 00 00 00 44 4F 4E 27 54 20 50 41 4E 49 43

Теперь делим полученный полином на полином-генератор (74 E7 D8 1E 01), берём остаток от деления (DB 22 58 5C) и записываем вместо нулей к полиному, который мы делили. (это эквивалентно операции сложения):

DB 22 58 5C 44 4F 4E 27 54 20 50 41 4E 49 43

Вот эта строка как раз и будет кодом Рида-Соломона для сообщения «DON'T PANIC» с 4 избыточными символами.

Как раскодировать сообщение?

Здесь всё посложнее будет. Ненамного, но всё же. Вопрос про раскодировать, собственно "не вопрос!" – убираем избыточные символы и остаётся исходное сообщение. Вопрос в том, как узнать, были ли ошибки при передаче, и если были, то как их исправить.

В первую очередь нужно отметить, что при проверке на наличие ошибок нужно знать количество избыточных символов. А во-вторую – надо научиться считать значение полинома при определённом . Про количество избыточных символов нам должен заранее сообщить тот, кто кодировал сообщение, а вот чтобы вычислить значение полинома нужно написать ещё одну функцию для работы с полиномами. Это элементарщина – просто вместо подставляется нужное значение. Но пример, всё же, никогда не помешает.

Пример: Нужно вычислить полиномпри . Подставляем, возводим в степень: , перемножаем, , складываем и получаем число . Сложение, умножение и возведение в степень здесь по правилам поля Галуа GF[256] (порождающий полином 285)

Код приводить не буду, оставлю ссылку на гитхаб. Там всё что я описывал в этой и предыдущих статьях реализовано на C#, в виде демо-приложения (собирается под win в VS2019, бинарник тоже выложен). Можно посмотреть как работает арифметика в поле Галуа, а также посмотреть, как работает кодирование Рида-Соломона.

Итак, прежде чем исправлять "ошибки" или "опечатки" нужно узнать есть ли они. Элементарно. Нужно вычислить полином принятого сообщения с избыточными символами при равном степеням примитивного члена. Это те же числа, которые мы использовали при составлении полинома-генератора: , – количество избыточных символов, – примитивный член. Если ошибок нет, то все вычисленные значения будут равны нулю. Закодированное ранее сообщение «DON'T PANIC» с 4 избыточными символами, в виде полинома в шестнадцатеричном представлении:

DB 22 58 5C 44 4F 4E 27 54 20 50 41 4E 49 43,

если вычислить этот полином при равном 2, 4, 8, 16, то получатся значения: 0, 0, 0, 0, ведь здесь сообщение точно в таком же виде, в котором оно и было закодировано. Если изменить хотя бы один байт, например, последний символ сделаем более правильным: 42 вместо 43:

DB 22 58 5C 44 4F 4E 27 54 20 50 41 4E 49 42,

то результат такого же вычисления станет равным 13, 18, B5, 5D. Эти значения называются синдромами. Их тоже можно принять за полином. Тогда это будет полином синдромов.

Итак, чтобы узнать есть ли ошибки в принятом сообщении, нужно посчитать полином синдромов. Если он состоит из одних нулей (также можно говорить, что он равен нулю), то ошибок нет.

Важное, но совсем занудное дополнение

Может случиться так, что сообщение с ошибками будет иметь синдром равным нулю. Это случится в том случае, когда полином амплитуд ошибок (о нём будет ниже) кратен полиному-генератору. Так что проверку ошибок по полиному синдромов кода Рида-Соломона нельзя считать 100% гарантией отсутствия ошибок. Можно даже посчитать вероятность такого случая.

Допустим мы кодируем сообщение из 4 символов четырьмя же избыточными символами, то есть передаём 8 байт. Также возьмём для примера вероятность ошибки при передаче одного символа в 10%. То есть, в среднем на каждые 10 символов приходится один, который передался как случайное число от 00 до FF. Это, конечно же совсем синтетическая ситуация, которая вряд ли будет в реальности, но здесь можно точно вычислить вероятности.

Для рассчёта я рассуждаю так: Полиномы, кратные полиному-генератору получаются умножением генератора на другие полиномы. Пятизначный кратный полином - получается умножением на константу от 1 до 255. Шестизначный - умножением на бином первой степени а их, без нулей ровно Те же рассуждения для 7 и 8 -значных полиномов, кратных генератору. Затем надо найти вероятности выпадения 5, 6, 7 и 8 ошибок подряд, и для каждой из них вычислить вероятность, что такая случайная последовательность ошибок окажется кратной полиному-генератору. Сложить их, и тогда мы получим вероятность того, что при передаче 4 байт с 4 избыточными символами, при вероятности ошибки при передаче одного символа 10% получится не обнаруживаемая кодом Рида-Соломона ошибочная передача. Рассчёт в маткаде:

Итого, на каждые ~500 Тб при такой передаче окажется один блок из 4 ошибочных символов, которые алгоритм посчитает корректными. Цифры большие, но вероятность не 0. При вероятности ошибки в 1% речь идёт об эксабайтах. Рассчёт, конечно не эталон, может быть даже с ошибками, но даёт понять об порядках чисел.

Что же делать, если синдром не равен нулю? Конечно же исправлять ошибки! Для начала рассмотрим случай с "опечатками", когда мы точно знаем номера позиций некорректно принятых байт. Ошибёмся намеренно в нашем закодированном сообщении 4 раза, столько же, сколько у нас избыточных символов:

DB 22 58 5C 44 4F 4E 27 54 20 41 41 41 41 41

41 – это буква A, поэтому их 5 подряд получилось. Позиции ошибок считаются слева направо, начиная с 0. Для удобства используем шестнадцатеричную систему при нумерации:

00 01 02 03 04 05 06 07 08 09 0A 0B 0C 0D 0E
DB 22 58 5C 44 4F 4E 27 54 20 50 41 4E 49 43
DB 22 58 5C 44 4F 4E 27 54 20 41 41 41 41 41

Позиции ошибок: 0A 0C 0D 0E.

Итак, если мы находимся на стороне приёмника, то у нас есть следующая информация:

• Сообщение с 4 избыточными символами;

• само сообщение: DB 22 58 5C 44 4F 4E 27 54 20 41 41 41 41 41;

• В сообщении есть ошибки в позициях 0A 0C 0D 0E.

Этого достаточно, чтобы восстановить сообщение в исходное состояние. Но обо всём по порядку.

Для продолжения необходимо разучить ещё одну операцию с полиномами в полях Галуа - взятие формальной производной от полинома. Формальная производная полинома в поле Галуа похожа на обычную производную. Формальной она называется потому, что в полях вроде GF[256] нет дробных чисел, и соответственно нельзя определить производную, как отношение бесконечно малых величин. Вычисляется похоже на обычную производную, но с особенностями. Если при обычном дифференцировании , то для формальной производной в поле Галуа с основанием 2, формула для дифференцирования члена такая: . Это значит, что достаточно просто переписать полином, начиная с первой степени (нулевая выкидывается) и у оставшегося убрать (обнулить, извиняюсь) члены с нечётными степенями. Пример:

Необходимо найти производную

(Это рандомный полином, не связан с примером). Производная суммы равна сумме производных, соответственно применяем формулу для производной члена и получаем:

Или, если записывать в шестнадцатеричном виде, то это же самое выглядит так:

(01 2D A5 C6 8C DF )' = 2D 00 C6 00 DF .

Думаю, что из примера в шестнадцатеричном виде проще всего составить алгоритм нахождения формальной производной.

Теперь можно уже исправить "опечатки"? Как бы не так! Нужно ещё два полинома. Полином-локатор и полином ошибок.

Полином-локатор – это полином, корнями которого являются числа обратные примитивному члену в степени позиции ошибки. Сложно? Можно проще. Полином-локатор это произведение вида

где – это примитивный член, и так далее – это позиции ошибок.

Пример: у нас есть позиции ошибок 10, 12, 13, 14; примитивный член тогда полином локатор будет таким:

Перемножаем и получаем полином-локатор для позиций ошибок 10, 12, 13, 14:

Или в шестнадцатеричной записи: 01 2D A5 C6 8C.

Про полином-локатор нужно понять следующее: из него можно получить позиции ошибок, и наоборот – из позиций ошибок можно получить полином-локатор. По сути, это две разные записи одного и того же – позиций ошибок.

Полином ошибок – его по-разному называют в разных статьях, он не так уж и сложен. Представляет из себя произведение полинома синдромов и полином-локатора, с отброшенными старшими степенями. Продолжая пример, найдём полином ошибок для искажённого сообщения:

DB 22 58 5C 44 4F 4E 27 54 20 41 41 41 41 41

Полином синдромов: 72 BD 22 5B

Произведение полинома синдромов и полинома-локатора не буду расписывать в "математическом" виде, напишу так:

(72 BD 22 5B)(01 2D A5 C6 8C) = 72 4B 10 22 D9 C0 57 15

У результата оставляем количество младших членов, равное количеству избыточных символов, в нашем случае их 4, старшие степени просто выбрасываем, они не нужны. Остаётся

72 4B 10 22

Это и есть полином ошибок.

Осталось посчитать амплитуды ошибок. Звучит угрожающе, но на деле это просто значения, которые нужно прибавить к искажённым символам сообщения чтобы получились неискажённые символы. Для этого воспользуемся алгоритмом Форни. Здесь придётся привести фрагмент кода, словами расписать так, чтобы было понятно, очень сложно.

Функция принимает на входе

• полином синдромов (Syndromes),

• полином, в котором члены – позиции ошибок (ErrPos),

• количество избыточных символов (NumOfErCorrSymbs).

Класс GF_Byte - это просто байт, для которого переопределены арифметические операции так, чтобы они выполнялись по правилам поля Галуа GF[256], класс GF_Poly – Это полином в поле Галуа. По сути, массив GF_Byte. Для него также переопределны арифметические операции так, чтобы они выполнялись по правилам арифметики с полиномами в полях Галуа.

public static GF_Poly FindMagnitudesFromErrPos(
   GF_Poly Syndromes,
   GF_Poly ErrPos,
   uint NumOfErCorrSymbs)
{
 	//Вычисление локатора из позиций ошибок
	GF_Poly Locator = CalcLocatorPoly(ErrPos);
	//Произведение для вычисления полинома ошибок
	GF_Poly Product = Syndromes * Locator;
	//Полином ошибок. DiscardHiDeg оставляет указаное количество младших степеней
	GF_Poly ErrPoly = Product.DiscardHiDeg(NumOfErCorrSymbs);
	//Производная локатора
	GF_Poly LocatorDer = Locator.FormalDerivative();
	//Здесь будут амплитуды ошибок. Количество членов - это самая большая позиция ошибки
	GF_Poly Magnitudes = new GF_Poly(ErrPos.GetMaxCoef());

	//Перебор каждой заданной позиции ошибки
	for (uint i = 0; i < ErrPos.Len; i++) {
		//число обратное примитивному члену в степени позиции ошибки
		GF_Byte Xi = 1 / GF_Byte.Pow_a(ErrPos[i]);
		//значение полинома ошибок при x = Xi
		GF_Byte W = ErrPoly.Eval(Xi);
		//значение производной локатора при x = Xi
		GF_Byte L = LocatorDer.Eval(Xi);
		//Это как раз и будет найденное значение ошибки,
		//которое надо вычесть из ошибочного символа, чтобы он стал не ошибочным
		GF_Byte Magnitude = W / L;
		//запоминаем найденную амплитуду в текущей позиции ошибки
		Magnitudes[ErrPos[i]] = Magnitude;
	}            
	return Magnitudes;
}

Если скормить функции следующие параметры:

• полином синдромов 72 BD 22 5B

• полином, в котором члены - позиции ошибок 0A 0C 0D 0E

• количество символов коррекции ошибок 4,

то на выходе она даст полином амплитуд ошибок:

00 00 00 00 00 00 00 00 00 00 11 00 0F 08 02.

Теперь можно прибавить полученное к искажённому сообщению

DB 22 58 5C 44 4F 4E 27 54 20 41 41 41 41 41

(по правилам сложения полиномов, конечно же), и на выходе получится исходное сообщение:

DB 22 58 5C 44 4F 4E 27 54 20 50 41 4E 49 43.

Первые 4 байта - это избыточные символы. Если бы в них оказались "опечатки", то разницы никакой для алгоритма нет, разве что они нам не нужны после исправления. Можно их просто отбросить:

44 4F 4E 27 54 20 50 41 4E 49 43 Это исходное сообщение «DON'T PANIC».

Здесь должно быть понятно, как исправлять ошибки, положение которых известно. Само по себе уже это может нести практическую пользу. В QR кодах на обшарпанных стенах могут стереться некоторые квадратики, и программа, которая их расшифровывает сможет определить в каких именно местах находятся байты, которые не удалось прочитать, которые "стёрлись" – erasures, или как мы договорились писать по-русски "опечатки". Но нам этого, конечно же недостаточно. Мы хотим уметь выявлять испорченные байты без дополнительной информации, чтобы передавать их по радио, или по лазерному лучу, или записывать на диски (кого я обманываю? CD давно мертвы), может быть, захотим реализовать передачу через ультразвук под водой, чтобы управлять моделью подводной лодки, а какие-нибудь неблагодарные дельфины будут портить случайные данные своими песнями. Для всего этого нам понадобится уметь выявлять, в каких именно байтах при передаче попортились биты.

Как найти позиции ошибок?

Вспомним про полином-локатор. Его можно составить из заранее известных позиций ошибок, а ещё его можно вычислить из полинома-синдромов и количества избыточных символов. Есть не один алгоритм, который позволяет это сделать. Здесь будет алгоритм алгоритм Берлекэмпа-Мэсси. Если хочется много математики, то гугл с википедией на неё не скупятся. Я, если честно, не вник до конца в циклические полиномы и прочее-прочее-прочее. Стыдно, немножко, конечно, но я взял реализацию этого алгоритма с сайта Wikiversity переписал его на C#, и постарался сделать его более доходчивым и читаемым:

public static GF_Poly CalcLocatorPoly(GF_Poly Syndromes, uint NumOfErCorrSymbs) {
	//Алгоритм Берлекэмпа-Мэсси
	GF_Poly Locator;
	GF_Poly Locator_old;
	
	//Присваиваем локатору инициализирующее значение (1*X^0)
	Locator = new GF_Poly(new byte[] { 1 });
	Locator_old = new GF_Poly(Locator);

	uint Synd_Shift = 0;

	for (uint i = 0; i < NumOfErCorrSymbs; i++) {
		uint K = i + Synd_Shift;
		GF_Byte Delta = Syndromes[K];

		for (uint j = 1; j < Locator.Len; j++) {
			Delta += Locator[j] * Syndromes[K - j];
		}
		//Умножение полинома на икс (эквивалентно сдвигу вправо на 1 байт)
		Locator_old = Locator_old.MultiplyByXPower(1);
		if (Delta.val != 0) {
			if (Locator_old.Len > Locator.Len) {
				GF_Poly Locator_new = Locator_old.Scale(Delta);
				Locator_old = Locator.Scale(Delta.Inverse());
				Locator = Locator_new;
			}
			//Scale – умножение на константу. Можно было бы
			//вместо использования Scale 
			//умножить на полином нулевой степени. Разницы нет, но так короче:
			Locator += Locator_old.Scale(Delta);
		}
	}
	return Locator;
}

Приведённый алгоритм считает локатор. Если количество "ошибок" больше, чем количество избыточных символов, поделённое на 2, то алгоритм не сработает правильно.

Если в сообщении, которое мы используем для примера –

DB 22 58 5C 44 4F 4E 27 54 20 50 41 4E 49 43,

ошибиться в нулевом и последнем символе (2 "ошибки", мы притворяемся, что не знаем в каких позициях ошиблись), получится такой полином:

02 22 58 5C 44 4F 4E 27 54 20 50 41 4E 49 01,

Полином синдромов для него 4B A7 E8 BD. Если выполнить функцию, приведённую выше с параметрами 4B A7 E8 BD, и 4 (количество избыточных символов), то она вернёт нам такой полином: 01 12 13. Это не похоже на позиции ошибок, которые мы ожидаем, но полином-локатор содержит в себе информацию о позициях ошибок, ведь это "полином, корнями которого являются числа обратные примитивному члену в степени позиции ошибки". Из этого, если немного поскрипеть мозгами или ручкой по бумаге следует, что позиция ошибки – это логарифм числа по основанию примитивного члена, обратного корню полинома.

E – позиция ошибки, a – примитивный член (2, как правило), R – корень полинома.

Что-ж, будем искать корни в поле. Поиск корней полинома в поле Галуа занятие лёгкое и непыльное. В GF[256] может быть 256 чисел всего, так что иксу негде разгуляться. Просто считаем полином 256 раз, подставляя вместо x число, и если полином посчитался как нуль, то записываем к массиву с корнями текущее значение x. Дальше считаем по формуле и получаем позиции ошибок 00 и 0E, именно там где они и были допущены. Теперь эти значения вместе с синдромами и цифрой 4 можно скармливать алгоритму Форни, чтобы он исправил "ошибки" также, как он исправлял "опечатки".

Если у нас есть 4 избыточных символа, при этом мы знаем что есть 2 "опечатки" в известных позициях, то алгоритм Берлекэмпа-Мэсси сможет найти ещё одну "ошибку". Но для этого его нужно будет совсем немного модифицировать. Всего то надо там где мы писали

	//Присваиваем локатору инициализирующее значение (1*X^0)
	Locator = new GF_Poly(new byte[] { 1 });

нужно локатор инициализировать не единичным полиномом, а полиномом-локатором, рассчитанным из известных позиций ошибок. И ещё изменить пару строчек. Весь код, напомню, есть на гитхабе.

Надеюсь материал в этой статье поможет тем, кто захочет в каком-нибудь своём проекте реализовать избыточное кодирование без сторонних библиотек. Просьба: Если что-то не понятно, не стесняйтесь комментировать. Постараюсь ответить на вопросы, или внести правки в статью.