А.Б.Сергиенко. Пример использования пакета Communications - моделирование модема V.32bis

Одну из разновидностей систем связи, хорошо знакомых пользователям персональных компьютеров, представляют собой модемы, предназначенные для передачи данных по телефонным сетям. Рассмотрим моделирование одного из реальных модемных стандартов средствами пакета Communications. Речь пойдет о моделировании приемника и передатчика модема, работающего в соответствии с Рекомендацией ITU-T V.32bis [1]. Достаточно подробное описание процесса формирования сигнала имеется также в книге [2]. Конечно, данный протокол в настоящее время является устаревшим, но он обладает многими характерными чертами систем цифровой связи и благодаря своей простоте удобен для использования в учебных и демонстрационных целях. Рекомендация V.32bis предусматривает различные скорости передачи данных - от 4800 до 14 400 бит/с. В данном примере будет моделироваться обработка сигнала для максимально возможной скорости 14 400 бит/с.

Вы можете загрузить архив model32.zip, содержащий следующие файлы, относящиеся к данному примеру:

• mod32144.m - функция моделирования передатчика модема;
• demod32144.m - функция моделирования приемника модема;
• map32144.mat - карта сигнального созвездия.

Краткое описание стандарта V.32bis

Входным сигналом модема является последовательный битовый поток данных, для кодирования и модуляции которого в протоколе V.32bis используется следующая последовательность процедур.

1. Входной битовый поток данных проходит через скремблер, структурная схема которого показана на рис. 1. Скремблер представляет собой рекурсивный фильтр, операции сложения в котором производятся по модулю два. Порядок фильтра m = 23, положение отвода: для вызывающего модема n = 18, для отвечающего модема n = 5. Назначение скремблера - сделать сигнал псевдослучайным.
2. 
   

   

   Рис. 1. Структурная схема скремблера

3. Скремблированный битовый поток разбивается на символы длиной по 6 бит.
4. Два первых (во времени) бита сформированного символа подвергаются дифференциальному кодированию. Это кодирование математически описывается следующим образом: Y(n) = (Q(n) - Y(n - 1)) mod 4. Здесь Q(n) - целое число, представляемое двумя кодируемыми битами, Y(n) - результирующее двухбитовое целое число, Y(n - 1) - результат дифференциального кодирования предыдущего символа.
5. Дифференциально-кодированные биты Y₁(n) и Y₂(n) подаются на сверточный кодер, структурная схема которого приведена на рис. 2. Кодер формирует избыточный бит Y₀(n), пропуская на выход биты Y₁(n) и Y₂(n) без изменений. Следует обратить внимание на то, что избыточный бит Y₀(n) зависит только от предыдущих, но не от текущих значений входных битов Y₁(n) и Y₂(n).

   

   Рис. 2. Структурная схема сверточного кодера

6. Избыточный бит Y₀(n), дифференциально-кодированные биты Y₁(n) и Y₂(n) и четыре бита, не подвергавшиеся дифференциальному кодированию, образуют символ, передаваемый в течение одного символьного такта (символьная скорость - 2400 символов/с). Этот символ отображается на комплексную плоскость с использованием созвездия квадратурной манипуляции, показанного на рис. 3. В рассматриваемых далее MATLAB-примерах используется созвездие, сформированное заранее в виде вектора комплексных чисел и сохраненное в MAT-файле map32144.mat (идентификатор переменной - constellation). Для создания рис. 3 использована функция modmap, вызов которой осуществляется для этого следующим образом:

   load map32144
modmap('qask/arb', real(constellation), imag(constellation))

   

   Рис. 3. Созвездие квадратурной манипуляции, используемое модемами V.32bis при скорости 14 400 бит/с

7. Для сужения полосы модулирующего сигнала и ослабления межсимвольной интерференции полученный после отображения комплексный сигнал подвергается интерполяции с использованием формирующего фильтра нижних частот (ФНЧ). В качестве такого фильтра, как правило, используется ФНЧ с косинусоидальным сглаживанием АЧХ (raised cosine).
8. Интерполированный комплексный сигнал s_M(t) после фильтрации используется для модуляции несущего колебания с частотой f₀ = 1800 Гц, в результате чего получается сигнал передачи данных s(t) = Re[s_M(t) exp(j2pf₀t)].

Функция моделирования передатчика модема V.32bis

Реализуем модель передатчика модема в виде функции mod32144, которая будет получать на входе вектор передаваемых битов x и возвращать вектор отсчетов модулированного сигнала y:

function y = mod32144(x)

На рис. 4 показана структурная схема передатчика модема с подписанным функциями пакета Communications, которые используются при реализации различных блоков.

Рис. 4. Блок-схема передатчика модема V.32bis

Скремблер, представляющий собой дискретный фильтр, математические операции в котором выполняются по модулю два, удобно реализуется с помощью функции gffilter:

xs = gffilter(1, [1 zeros(1,17) 1 zeros(1,4) 1], x(:)');

Разделение скремблированного битового потока xs на шестибитовые символы выполняется с помощью функции vec2mat:

sy = vec2mat(xs, 6);

Далее два первых бита символов sy подвергаются дифференциальному кодированию, при этом преобразование пар битов в целые числа и обратно осуществляется с помощью функций bi2de и de2bi соответственно:

bits12 = sy(:,1:2);
y12 = bi2de(bits12);
y12(2:2:end) = -y12(2:2:end);
y12m = cumsum(y12);
y12m(2:2:end) = -y12m(2:2:end);
y12m = mod(y12m, 4);
bits12m = de2bi(y12m, 2);

Сверточный код, соответствующий структурной схеме рис. 2, является нелинейным, поэтому он не может быть задан в виде коэффициентов полиномов. Заполняем поля структуры conv_code, описывающей код, вручную, предварительно приведя некоторые комментарии.

• Кодер содержит три элемента памяти (на рис. 2 они обозначены буквой T, поскольку осуществляют задержку сигнала на один такт) и, таким образом, имеет 2³ = 8 внутренних состояний. Это значение присваивается полю numStates.
• На вход кодера поступает два бита, для которых возможно 2² = 4 различных комбинации. Это значение присваивается полю numInputSymbols.
• Кодер добавляет к этим двум битам избыточный третий, что дает на выходе 2³ = 8 комбинаций. Это значение присваивается полю numOutputSymbols.
• Данный кодер является систематическим, поэтому выходная комбинация бит получается из входной путем добавления третьего (старшего) бита, равного младшему биту исходного внутреннего состояния. Поэтому поле outputs (оно должно представлять собой матрицу с числом строк, равным числу состояний, то есть 8, и с числом столбцов, равным числу возможных входных символов, то есть 4) должно содержать строки двух типов: [0 1 2 3], повторяющие значения входных символов, для четных значений внутреннего состояния, и [4 5 6 7] (входной символ с добавленным единичным старшим битом) для нечетных значений внутреннего состояния.
• Поле nextStates задает алгоритм переходов между внутренними состояниями кодера. Оно должно представлять собой матрицу с числом строк, равным числу состояний, то есть 8, и с числом столбцов, равным числу возможных входных символов, то есть 4. Это поле заполняется вручную в соответствии с логикой работы кодера.

Реализуем сказанное в виде кода MATLAB:

conv_code.numInputSymbols = 4;
conv_code.numOutputSymbols = 8;
conv_code.numStates = 8;
conv_code.nextStates = [0 2 3 1; 4 7 5 6; ...
                        1 3 2 0; 7 4 6 5; ...
                        2 0 1 3; 6 5 7 4; ...
                        3 1 0 2; 5 6 4 7];
conv_code.outputs = repmat([0 1 2 3; 4 5 6 7], 4, 1);

Теперь осуществляем собственно сверточное кодирование дифференциально-кодированных битов bits21m с помощью функции convenc:

bits12m = bits12m';
bits12m = bits12m(:);
bits012 = convenc(bits12m, conv_code);

Результат работы сверточного кодера bits012 представляет собой вектор, в котором последовательно расположены тройки интересующих нас битов. Преобразуем вектор в матрицу с помощью функции vec2mat и объединяем результат с теми четырьмя столбцами из матрицы sy, которые не подвергались дополнительным преобразованиям. В результате получаем семибитовые символы в битовом (sy7) и десятичном (sy_dec) представлении:

sy012 = vec2mat(bits012, 3);
sy7 = [(sy012) sy(:,3:6)];
sy_dec = bi2de(sy7, 'left-msb');

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

load map32144
sy_mapped = modmap(sy_dec, 1, 1, 'qask/arb', ...
                      real(constellation), imag(constellation));

Теперь производим формирование спектра, выполняя интерполяцию сигнала sy_mapped с использованием SQRT-варианта фильтра с косинусоидальным сглаживанием АЧХ. Частота дискретизации повышается при этом в восемь раз, до значения 2400 x 8 = 19200 Гц:

Fd = 2400;
Fs = Fd*8;
s_m = rcosflt(sy_mapped, Fd, Fs, 'sqrt');

Наконец, интерполированный сигнал s_m подвергается аналоговой квадратурной модуляции с помощью функции amod:

Fc = 1800;
y = amod(s_m, Fc, Fs, 'qam');

Модулированный сигнал y является результатом работы функции.

Тестирование функции mod32144

Протестируем функцию, подав на ее вход вектор из 10 000 бит, сформированный с помощью функции randint:

x = randint(1, 10000);
y = mod32144(x);

Поскольку при квадратурной модуляции амплитуда и фаза сигнала меняются одновременно, график полученного сигнала не слишком нагляден. Однако все же выведем его начальный фрагмент, а также покажем оценку спектральной плотности мощности полученного сигнала, произведя ее методом Уэлча с помощью функции pwelch пакета Signal Processing:

plot(y(1:200))
pwelch(y, [], [], [], 19200)

Полученные графики представлены на рис. 5.

Рис. 5. Сформированный сигнал модема V.32bis (слева) и его спектр мощности (справа)

При наличии звуковой карты можно прослушать полученный сигнал, например, с помощью функции soundsc:

soundsc(repmat(y, 10, 1), 19200)

Функция repmat использована здесь для того, чтобы повторить сформированный вектор 10 раз, иначе звук оказывается слишком коротким.

Если вы когда-нибудь слышали шуршащий звук модема, то, без сомнения, узнаете его.

Функция моделирования приемника модема V.32bis

Функция demod32144, реализующая приемник модема V.32bis, выполняет преобразования, обратные по отношению к процессам, происходящим в передатчике. Функция принимает на входе массив отсчетов модулированного сигнала y и возвращает вектор принятых битов z:

function z = demod32144(y)

На рис. 6 показана структурная схема приемника модема с подписанным функциями пакета Communications, которые используются при реализации различных блоков.

Рис. 6. Блок-схема приемника модема V.32bis

Прежде всего зададим необходимые константы: символьную скорость Fd, несущую частоту Fc и частоту дискретизации Fs:

Fd = 2400;
Fc = 1800;
Fs = Fd*8;

Как уже говорилось, функции цифровой модуляции пакета Communications не поддерживают формирование спектра при квадратурной манипуляции. Поэтому демодуляцию сигнала придется реализовать, вручную выполнив необходимые шаги. Сначала рассчитаем фильтр с косинусоидальным сглаживанием АЧХ, идентичный тому, который использовался в передатчике. Для этого служит функция rcosine:

b = rcosine(Fd, Fs, 'sqrt');

Теперь производим аналоговую демодуляцию с помощью функции ademod, указав рассчитанный фильтр в качестве используемого ФНЧ:

x_m = ademod(y, Fc, Fs, 'qam', b, 1);

Выводим диаграмму рассеяния с помощью функции scatterplot, чтобы посмотреть, как принятые точки ложатся на комплексную плоскость:

scatterplot(x_m, Fs/Fd)

Далее с помощью функции demodmap ищем точки используемого созвездия, ближайшие к принятым точкам (символьная скорость задана в виде двухэлементного вектора, второй элемент которого указывает положение оптимальных точек дискретизации сигнала):

load map32144
sy = demodmap(x_m, [Fd 1], Fs, 'qask/arb', ...
           real(constellation), imag(constellation));

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

for k = 1:8
    sy1(:,k) = demodmap(x_m, [Fd 1], Fs, 'qask/arb', ...
        real(constellation((k-1)*16+(1:16))), ...
        imag(constellation((k-1)*16+(1:16))));
    sy1(:,k) = sy1(:,k) + (k-1)*16;
end

При обработке фильтрами в приемнике и передатчике сигнал приобретает некоторую задержку, которая в нашем примере равна шести символам (по умолчанию функции rcosflt и rcosine рассчитывают фильтр, длина импульсной характеристики которого соответствует шести символьным тактам). Поэтому необходимо удалить первые шесть принятых символов, соответствующих начальному переходному процессу фильтрации. Далее вектор десятичных представлений символов sy преобразуется в матрицу битовых представлений с помощью функции de2bi:

sy = sy(7:end);
sy1 = sy1(7:end,:);
sy_bits = de2bi(sy, 7, 'left-msb');

Подготовка структуры conv_code, описывающей сверточный код, производится точно так же, как это делалось в передатчике:

conv_code.numInputSymbols = 4;
conv_code.numOutputSymbols = 8;
conv_code.numStates = 8;
conv_code.nextStates = [0 2 3 1; 4 7 5 6; ...
                        1 3 2 0; 7 4 6 5; ...
                        2 0 1 3; 6 5 7 4; ...
                        3 1 0 2; 5 6 4 7];
conv_code.outputs = repmat([0 1 2 3; 4 5 6 7], 4, 1);

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

y012 = sy_bits(:,1:3);
y012 = y012';
y012 = y012(:);
y12 = vitdec(y012, conv_code, 40, 'trunc', 'hard');

Алгоритм Витерби попытался исправить ошибки, имевшиеся в двухбитовых комбинациях, защищенных сверточным кодом. Однако в каждом символе имеются еще четыре бита, передававшиеся без такой защиты. Чтобы обеспечить коррекцию полных шестибитовых символов, мы подвергаем исправленную последовательность y12 сверточному кодированию, чтобы получить исправленные значения избыточного бита. Результирующие трехбитовые комбинации определяют номера частичных созвездий, в которых следует искать точки, соответствующие исправленным символам. Результаты отображения принятых точек с использованием частичных созвездий сохранены в матрице sy1, так что сейчас мы просто выбираем из нее нужные элементы, сформировав набор линейных индексов с помощью функции sub2ind. После этого исправленные символы преобразуются в битовые представления с помощью функции de2bi:

bits012_corr = convenc(y12, conv_code);
bits012_corr = vec2mat(bits012_corr, 3);
y012_corr = bi2de(bits012_corr, 'left-msb');
sy = sy1(sub2ind(size(sy1), (1:length(sy))', y012_corr+1));
sy_bits = de2bi(sy, 7, 'left-msb');

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

y12 = sy_bits(:,2:3);
z12 = bi2de(y12);
z12dec = bitand([0;z12(1:end-1,:)]+z12, 3);
sy_bits(:,2:3) = (de2bi(z12dec));

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

zs = sy_bits(:, 2:end);
zs = zs';
zs = zs(:);
z = gffilter([1 zeros(1,17) 1 zeros(1,4) 1], 1, zs');

Вектор дескремблированных битов z является результатом работы функции.

Тестирование функции demod32144

Протестируем функцию, подав на ее вход сформированный ранее вектор y отсчетов модулированного сигнала:

z = demod32144(y);

На рис. 7 показана диаграмма рассеяния, сформированная нашей функцией. Видно, что принятые точки соответствуют использованному созвездию (см. ранее рис. 3).

Рис. 7. Диаграмма рассеяния, выведенная функцией demod32144 при приеме сигнала

Проверим правильность приема, сравнив переданную и принятую последовательности бит с помощью функции symerr (поскольку мы сравниваем векторы, содержащие только нули и единицы, функции symerr и biterr в данном случае оказываются эквивалентными). Принятый вектор содержит на два элемента больше, чем переданный, поскольку при передаче входной сигнал был дополнен двумя нулевыми битами для формирования целого числа шестибитовых символов, поэтому при сравнении нужно выделить первые 10 000 элементов:

symerr(x, z(1:10000))

ans =
0

Как видите, сигнал принят правильно.

Заключение

Данный пример является демонстрационным и не реализует ряда функций, необходимых в реальной системе цифровой связи. К таким функциям относятся адаптивная фильтрация для компенсации искажений, вносимых каналом связи, временная синхронизация, передача служебных сигналов для установления соединения и т. п. Все это также может быть реализовано средствами MATLAB с использованием прикладных пакетов Communications, Signal Processing и Filter Design.

Литература

1. A Duplex Modem Operating at Data Signalling Rates of up to 14400 bit/s for Use on the General Switched Telephone Network and on Leased Point-to-Point 2-Wire Telephone-Type Circuits. ITU-T Recommendation V.32bis. ITU-T, Geneva, 1991.
2. Лагутенко О. И. Модемы. Справочник пользователя. СПб.: Лань, 1997. - 368 с.