Моделирование ЦОС цифровой обработки сигналов в MATLAB. Часть 6. Моделирование структур цифровых фильтров с фиксированной точкой программными средствами MATLAB: квантование воздействия и вычисление реакции

Все статьи цикла:

• Часть 1. Синтез оптимальных (по Чебышеву) КИХ-фильтров программными средствами MATLAB
• Часть 2. Синтез оптимальных цифровых БИХ-фильтров программными средствами MATLAB
• Часть 3. Описание структур КИХ- и БИХ-фильтров в MATLAB
• Часть 4. Моделирование структур цифровых фильтров с фиксированной точкой программными средствами MATLAB: анализ характеристик КИХ-фильтров
• Часть 5. Моделирование структур цифровых фильтров с фиксированной точкой программными средствами MATLAB: анализ характеристик БИХ-фильтров
• Часть 6. Моделирование структур цифровых фильтров с фиксированной точкой программными средствами MATLAB: квантование воздействия и вычисление реакции
• Часть 7. Моделирование цифровых фильтров средствами программ GUI MATLAB: GUI FDATool
• Часть 8. Моделирование цифровой фильтрации средствами программ GUI MATLAB: GUI SPTool
• Часть 9. Моделирование цифровых преобразователей Гильберта и дифференциаторов программными средствами и средствами GUI FDATool MATLAB

Моделирование квантования в АЦП

Квантование воздействия в АЦП (далее — квантование) в MATLAB может моделироваться одним из следующих способов:

• на основе объекта quantizer с помощью
  функции quantize;
• на основе объекта fi.

Моделирование квантования
на основе объекта quantizer

При моделировании квантования на основе объекта quantizer (Quantizer Object) сначала создается данный объект одним из следующих основных способов:

или

Здесь q — имя объекта quantizer;
′name1′,′name2′,… — свойства объекта quantizer;
value1,value2,… — значения свойств объекта quantizer; struct — имя массива записей
(struct array), поля которого отображают
свойства объекта quantizer.

Выходным параметром является объект q
(массив записей — struct array) со списком
свойств (полей). Их полный список выводится с помощью функции get(q), а основнойпо имени q.

Доступные пользователю свойства устанавливаются с помощью функции:

В таблице 1 представлены основные доступные пользователю свойства объекта quantizer,
а в таблице 2 — четыре важнейших из множества свойств, доступных только для чтения.

Таблица 1. Доступные пользователю свойства объекта quantizer

Таблица 2. Свойства объекта quantizer
только для чтения

Пример 1

Создать двумя способами объект quantizer с именем q и заданными свойствами и сохранить его на диске для дальнейшего использования.

1. Первый способ:

   >> q=quantizer(′Mode′,′fixed′,′RoundMode′,′convergent′,… ′OverflowMode′,′saturate′,′Format′,[16 15]) >> q =

   DataMode = fixed RoundMode = convergent OverflowMode = saturate Format = [16 15]

2. Второй способ (свойства объекта q точно
   такие же):

   >> struct.Mode=′fixed′; >> struct.RoundMode=′convergent′; >> struct.OverflowMode=′saturate′; >> struct.Format=[16 15]; >> q=quantizer(struct); >> save q

Объекты quantizer могут использоваться в качестве входных параметров многих функций MATLAB, основные из которых, связанные с копированием квантователя и преобразованием чисел, приведены в таблице 3 (использован сохраненный в примере 1 квантователь q).

Таблица 3. Функции MATLAB для объектов quantizer

На основе созданного объекта quantizer выполняется непосредственно квантование с помощью функции:

Здесь q — имя объекта quantizer; x — исходные неквантованные данные: вектор или матрица; xq — квантованные данные; размеры x и xq совпадают.

В результате квантования обновляются свойства объекта q, доступные только для чтения, которые рекомендуется проверять.

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

Здесь x— вектор с линейно меняющимися элементами; xb, xf — значения первого и последнего элементов вектора x; N— число элементов вектора (по умолчанию 100).

Пример 2

Рисунок. Характеристики квантователей при использовании: а) арифметики насыщения; б) модульной арифметики

На основе объекта quantizer с именем q и заданными свойствами (пример 1) выполнить квантование элементов вектора x при xb = –2 и xf = 2. Повторить процедуру для объекта qq (копии объекта q) со свойством OverflowMode = wrap (модульной арифметикой). Соответствующим векторам с квантованными элементами присвоить имена xq1 и xq2. Построить графики зависимостей x от
xq1 и x от xq2, отображающие характеристики квантователей при разных типах арифметики и иллюстрирующие суть различий между арифметикой насыщения и модульной арифметикой (рисунок):

Моделирование квантования на основе объекта fi

При моделировании квантования на основе объекта fi (Fixed-point Numeric Object) достаточно создать данный объект:

Здесь xin — имя объекта fi; v, s, w, f— входные параметры объекта fi; v— исходные неквантованные данные: вектор или матрица;
w— длина слова для данных с ФТ (по умолчанию 16); f — длина дробной части слова
для данных с ФТ (по умолчанию 15); s— скаляр, принимающий значения 1 (true) или 0
(false); при s = 0 (по умолчанию) формируются числа с ФТ без знака (беззнаковые —
положительные), а s = 1— со знаком и появляется возможность управлять длиной f
дробной части в слове длиной w.

Пример 3

Создать объект fi с именем xin и входными параметрами w = 16, f = 15, s = 1 и выполнить квантование элементов вектора x:

Выходным параметром является объект xin — массив записей (struct array), поля которого отображают свойства объекта. Их удобно разделить на три группы:

• Data Properties (свойства данных) — это группа полей, отображающих различные виды представления данных с ФТ, основными из которых являются:
  • xin.data — квантованные согласно входным параметрам объекта fi десятичные
    числа;
  • xin.double — квантованные согласно
    входным параметрам объекта fi десятичные числа типа double (при выводе не
    отличаются от xin.data);
  • xin.bin — двоичные эквиваленты чисел
    xin.data;
  • xin.hex — шестнадцатеричные эквиваленты чисел xin.data;
  • xin.int — целочисленные эквиваленты чисел xin.data.
• fimath Properties, где fimath — имя поля, в свою очередь, объединяющего группу полей (свойств) объекта xin, связанных с собенностями арифметики при представлении данных с ФТ, а также при работе с двумя объектами fi.
• numerictype Properties, где numerictype — имя поля, в свою очередь, объединяющего группу полей (свойств) объекта xin, связанных с типом данных, форматом, масштабированием и т. п., предназначенных только для чтения.

Пример 4

Создать объект fi с именем xtest и входными параметрами w = 16, f = 15, s = 1 и выполнить квантование элементов вектора a:

Вывести квантованные значения, а также их двоичные, шестнадцатеричные и целочисленные эквиваленты, и сравнить результаты
с полученными в таблице 3:

Пример 5

Для созданного в примере 3 объекта xin вывести список основных свойств (полей) группы fimath Properties:

Для данных с ФТ в группе fimath Properties представляют интерес два поля (два свойства) объекта xin:

• xin.RoundMode — отображает метод округления при представлении данных с ФТ (свойство RoundMode в табл. 1);
• xin.OverflowMode — отображает режим переполнения при сохранении данных в формате слова (свойство OverflowMode в табл. 1).

Остальные поля отображают особенности арифметики при работе с двумя объектами fi.

Присвоим свойству RoundMode значение ′convergent′ и сохраним объект xin на диске:

Пример 6

Для сохраненного на диске объекта xin (пример 5) вывести список основных свойств (полей) группы numerictype Properties:

Свойство DataTypeMode отображает режим масштабирования данных с ФТ. Его значение Fixed-point: binary point scaling соответствует отсутствию масштабирования данных с ФТ (о нем см. далее).

Остальные значения полей (свойств) объекта xin согласованы с входными параметрами объекта fi и предназначены только для
чтения.

Полные списки свойств (полей) объекта xin (пример 5), относящихся к группам fimath Properties и numerictype Properties, выводятся с помощью функций get(xin.fimath) и get(xin.numerictype), а общий список свойств — с помощью функции get(xin). В целях экономии места эти списки не приводятся.

Среди дополнительных свойств в группе numerictype Properties представляют интерес дополнительные свойства, связанные с масштабированием данных с ФТ. Определяющим является свойство Scaling.

Пример 7

Для сохраненного на диске объекта xin (пример 5) вывести значение свойства Scaling в группе numerictype Properties:

При значении свойства Scaling: ′BinaryPoint′
(пример 6), когда масштабирование отсутствует, представляют интерес следующие дополнительные свойства:

• FixedExponent — отображает показатель степени при определении веса младшего бита 2^{FixedExponent} в заданном формате данных с ФТ;
• Slope — отображает точность представления дробной части данных с ФТ — вес
  младшего значащего бита 2^{FixedExponent} (шаг
  квантования).

Пример 8

Для сохраненного на диске объекта xin (пример 5) определить вес младшего бита и сравнить его с выведенным значением свойства Slope:

Свойства группы numerictype Properties невозможно изменить после создания объекта fi. Для того чтобы появилась возможность управлять ими при создании объекта, следует использовать другие форматы описания объекта fi.

Например, для масштабирования данных с ФТ, когда их значения можно будет изменять согласно формуле:

где a_data — десятичное число с ФТ (xin.data);
a_int — целочисленный эквивалент (xin.int) десятичного числа a_data; bias — смещение; slope—
коэффициент масштабирования, равный:

следует выбрать формат:

Здесь slope соответствует свойству Slope и управляет, согласно (1), масштабированием данных, а именно: распределением битов
в слове длиной WordLength (свойство объекта fi). При slope = 2^{FixedExponent} длина дробной части вWordLength автоматически становится равной FractionLength = –FixedExponent, а старшие биты, длина которых составляет (WordLength – FractionLength – 1), интерпретируются как целая часть двоичного числа; bias соответствует свойству Bias (выводится по команде get) и управляет, согласно (1), смещением bias.

При выборе данного формата описания объекта fi изменяются свойства создаваемого объекта xin, которые рекомендуется проверять. В частности, свойства DataTypeMode и Scaling принимают значения:

Пример 9

Создать объект fi с именем x1 и параметрами w = 8, f = 7, s = 1 и выполнить квантование элементов вектора a:

Как видим, число 2.5555555, большее единицы, в данном формате не может быть представлено, и для него использована арифметика насыщения: присвоено максимальное в данном формате значение 0,9922.

Рассмотрим, что изменится при значении slope = 2ˆ(–5), когда на дробную часть FractionLength будет отведено не 7, а 5 битов,
и, соответственно, на целую часть — 2 бита (значение bias остается равным 0):

Новое значение параметра slope изменило и сами десятичные числа с ФТ, и их двоичный и целочисленный эквиваленты. В двоичных эквивалентах биты, интерпретируемые как целая часть двоичного числа, выделены полужирным шрифтом. Используя правила перевода десятичных чисел (целых и дробных) в двоичные, нетрудно убедиться в соответствии чисел x2.data и их двоичных эквивалентов. Однако целочисленные эквиваленты вычисляются на основании новых двоичных эквивалентов по «старому» правилу, в чем также легко убедиться:

Покажем, что десятичные числа x2.data совпадают с рассчитанными по формуле (1):

Добавим параметр bias = 5 при том же значении slope = 2ˆ(–5):

Покажем, что десятичные числа x3.data совпадают с рассчитанными по формуле (1):

Изменение данных с ФТ с помощью параметров slope и bias называется масштабированием данных с ФТ. Его не следует путать ни с нормированием исходных неквантованных данных, ни масштабированием с целью минимизации переполнений.

Формирование воздействия с ФТ

Сформировать воздействие с квантованными значениями, моделирующее цифровой сигнал с ФТ на выходе АЦП, можно одним из следующих двух способов:

• на основе объекта quantizer с помощью
  функции quantize;
• на основе объекта fi.

Отметим, что если исходные неквантованные значения воздействия по модулю превосходят единицу, то предварительно их следует нормировать, а далее при вычислении реакции учитывать нормирующий множитель.

Пример 10

Сформировать воздействие xq1, моделирующее цифровой сигнал с ФТ на выходе 12-разрядного АЦП, на основе объекта quantizer
с помощью функции quantize. В качестве исходного, дискретного, сигнала x_in с неквантованными значениями (числами типа double)
использовать равномерный белый шум, генерируемый с помощью функции rand(100,1),
значения которой не превосходят единицу
и потому не требуют нормирования. Сохранить векторы x_in и xq1 на диске.

Для формирования воздействия xq1 выполним следующие действия:

1. Создадим и сохраним на диске объект
   quantizer с именем q1, моделирующий
   12-разрядное АЦП, со следующими свойствами (пример 1):

   >> q1=quantizer(′Mode′,′fixed′,′RoundMode′,′convergent′,… ′OverflowMode′,′saturate′,′Format′,[12 11]) q1 =

   DataMode = fixed RoundMode = convergent OverflowMode = saturate Format = [12 11]

   >> save q1

2. Создадим вектор x_in неквантованных значений — дискретный сигнал, который
   в данном случае представляет собой равномерный белый шум, и сохраним его на
   диске:

   >> x_in=rand(100,1); >> save x_in

3. С помощью функции quantize создадим
   вектор xq1 квантованных значений с ФТ—цифровой сигнал на выходе 12-разрядного АЦП — и сохраним его на диске:

   >> xq1=quantize(q1,x_in); >> save xq1

Пример 11

Сформировать воздействие xq2, моделирующее цифровой сигнал с ФТ на выходе
12-разрядного АЦП на основе объекта fi. В качестве исходного дискретного сигнала с неквантованными значениями использовать
сохраненный на диске равномерный белый
шум x_in (пример 10).

Для формирования воздействия xq2 выполним следующие действия:

1. Создадим вектор x_in неквантованных значений — дискретный сигнал, в данном случае представляющий собой равномерный
   белый шум x_in:

   >> load x_in

2. Создадим объект fi с именем xq2 и входными параметрами w = 12, f = 11, s = 1
   и выполним квантование в векторе x_in:

   >> xq2=fi(x_in,1,12,11);

   Квантованные значения объекта xq2 хранятся в поле data, поэтому далее используется вектор xq2.data.

3. Установим требуемые свойства объекта xq2.
   Вывести полный список свойств объекта
   xq2 по команде get(xq2) и установить их
   требуемые значения.
   Из соображений экономии места список
   свойств не приводится. В данном случае изменим только значение свойства RoundMode,
   а остальные установим по умолчанию (они
   аналогичны свойствам объекта в примере 10):

   >> xq2.RoundMode=′convergent′;

4. Сохраним объект xq2 с новыми свойствами
   на диске:

   >> save xq2

Критерии сравнения
квантованного сигнала
с неквантованным

Для сравнения векторов с квантованными
и неквантованными элементами (квантованных и дискретных¹ сигналов) в качестве объективных критериев обычно используют
нормы:

где a, b — сравниваемые векторы (квантованный и дискретный сигналы); i — номер элемента; m — количество элементов.

Пример 12

Сравнить векторы x_in и xq2.data на основе норм (3–5), используя функцию norm:

Полученные значения свидетельствуют о близости векторов.

Вычисление реакции ЦФ с ФТ

Для вычисления реакции ЦФ с ФТ, описанных в виде объектов dfilt, используется функция filter, простейший формат которой при
нулевых начальных условиях имеет вид:

Здесь Hd — объект dfilt, описывающий КИХ- или БИХ-фильтр с ФТ; x, y—воздействие и реакция.

Реакция КИХ-фильтров с ФТ

Вычисление реакции КИХ-фильтров с ФТ поясним на примере. Воздействие сформируем на основе объекта quantizer с помощью функции quantize (пример 10).

Пример 13

Необходимо:

• вычислить реакцию y_out исходного КИХ-фильтра (объекта Hf3, пример 1 в [9]) на воздействие, заданное вектором x_in неквантованных значений (пример 10);
• вычислить реакцию yq_out КИХ-фильтра с ФТ с 16-разрядными коэффициентами (объекта Hq3c1, пример 5 в [9]) на воздействие, заданное вектором xq1 квантованных 12-разрядных значений (пример 10);
• сравнить реакции y_out и yq_out на основе норм (3–5).

Для решения задачи выполним следующие действия:

1. Загрузим сохраненный КИХ-фильтр с ФТ (объект Hq3c1) и выведем список его свойств (для экономии места он не приводится, смысл свойств поясняется в [5]):

   >> load Hq3c1 >> get(Hq3c1)

   Изменим значения требуемых свойств объекта Hq3c1 и сохраним объект Hq3c1 с новыми свойствами на диске:

   >> set(Hq3c1,′OverflowMode′,′saturate′) >> set(Hq3c1,′ProductMode′,′KeepMSB′) >> set(Hq3c1,′OutputMode′,′SpecifyPrecision′) >> set(Hq3c1,′OutputFracLength′,15) >> save Hq3c1

   Если пользователь ориентируется на реализацию фильтра с ФТ на базе конкретного ЦПОС, то соответствующие значения
   свойств должны быть согласованы с характеристиками ЦПОС. В данном случае значения свойств выбраны из соображений
   разумной логики.

2. Вычислим реакции y_out и yq_out с помощью функции filter:

   >> load Hf3 >> load x_in >> y_out=filter(Hf3,x_in); >> load Hq3c1 >> load xq1 >> yq_out=filter(Hq3c1,xq1);

   После загрузки объектов dfilt в рабочее
   пространство Workspace рекомендуется
   просмотреть, какие переменные для объекта dfilt были использованы по умолчанию, во избежание их дублирования.
   Реакция y_out исходного КИХ-фильтра на
   воздействие x_in представляет собой вектор неквантованных значений. Длины векторов x_in и y_out совпадают.
   Реакция yq_out КИХ-фильтра с ФТ представляет собой объект fi. Квантованные значения объекта yq_out хранятся в поле data, поэтому далее используется вектор yq_out.data.
   Длины векторов xq1 и yq_out.data совпадают.

3. Сравним векторы y_out и yq_out.data на
   основе норм (3–5):

   >> norm((y_out–yq_out.data),1) ans = 0.0056 >> norm(y_out–yq_out.data) ans = 7.2107e–004 >> norm((y_out–yq_out.data),inf) ans = 2.1064e–004

Полученные значения свидетельствуют о близости векторов.

Реакция БИХ-фильтров с ФТ

Вычисление реакции БИХ-фильтров с ФТ поясним на примере. Воздействие сформируем на основе объекта fi (пример 10).

Пример 14

Необходимо:

• вычислить реакцию Y_OUT исходного
  БИХ-фильтра (объекта Hn, пример 4 в [10])
  на воздействие, заданное вектором x_in неквантованных значений (пример 11);
• вычислить реакцию YQ_OUT БИХ-фильтра с ФТ с 16-разрядными коэффициентами
  (объекта Hqn2, пример 5 в [10]) на воздействие, заданное вектором xq2 квантованных
  12-разрядных значений (пример 11);
• сравнить реакции Y_OUT и YQ_OUT на основе норм (3–5).

Для решения задачи выполним следующие
действия:

1. Загрузим сохраненный БИХ-фильтр с ФТ
   (объект Hqn2) и выведем список его свойств
   (для экономии места он не приводится,
   смысл свойств поясняется в [5]):

   >> load Hqn2 >> get(Hqn2)

   Изменим значения требуемых свойств объекта Hqn2 и сохраним созданный объект Hqn2 с новыми свойствами на диске:

   >> set(Hqn2,′OverflowMode′,′saturate′) >> set(Hqn2,′ProductMode′,′KeepMSB′) >> set(Hqn2,′OutputMode′,′SpecifyPrecision′) >> set(Hqn2,′OutputFracLength′,15) >> set(Hqn2,′StageInputAutoScale′,0) >> set(Hqn2,′StageInputFracLength′,15) >> set(Hqn2,′StageOutputAutoScale′,0) >> set(Hqn2,′StageOutputFracLength′,15) >> save Hqn2

   2

2. Вычислим реакции Y_OUT и YQ_OUT с помощью функции filter:

   >> load Hn >> load x_in >> Y_OUT=filter(Hn,x_in); >> load xq2 >> load Hqn2 >> YQ_OUT=filter(Hqn2,xq2);

   Реакция Y_OUT исходного БИХ-фильтра
   на воздействие x_in представляет собой
   вектор неквантованных значений. Длины
   векторов x_in и Y_OUT совпадают.
   Реакция YQ_OUT БИХ-фильтра с ФТ
   представляет собой объект fi. Квантованные значения объекта YQ_OUT хранятся
   в поле data, поэтому далее используется
   вектор YQ_OUT.data. Длины векторов xq2
   и YQ_OUT.data совпадают.

3. Сравним векторы Y и Y1.data на основе
   норм (3–5):

   >> norm((Y_OUT–YQ_OUT.data),1) ans = 0.0177 >> norm(Y_OUT–YQ_OUT.data) ans = 0.0023 >> norm((Y_OUT–YQ_OUT.data),inf) ans = 5.6329e–004

Полученные значения свидетельствуют о близости векторов.

Литература

1. Ingle V., Proakis J. Digital Signal Processing Using
   MATLAB. Second Edition. Thomson.
2. Оппенгейм А., Шафер Р. Цифровая обработка
   сигналов. М.: Техносфера, 2006.
3. Сергиенко А. Б. Цифровая обработка сигналов,
   2-е изд. СПб.: ПИТЕР, 2006.
4. Солонина А. И., Улахович Д. А., Арбузов С. М.,
   Соловьева Е. Б. Основы цифровой обработки
   сигналов. 2-е изд. СПб.: БХВ-Петербург, 2005.
5. Солонина А. И., Арбузов С. М. Цифровая обработка сигналов. Моделирование в MATLAB.
   СПб.: БХВ-Петербург, 2008.
6. Солонина А. Моделирование цифровой обработки сигналов в MATLAB. Часть 1. Синтез оптимальных (по Чебышеву) КИХ-фильтров программными средствами MATLAB // Компоненты и технологии. 2008. № 11.
7. Солонина А. Моделирование цифровой обработки сигналов в MATLAB. Часть 2. Синтез оптимальных БИХ-фильтров программными средствами MATLAB // Компоненты и технологии.
   2008. № 12.
8. Солонина А. Моделирование цифровой обработки сигналов в MATLAB. Часть 3. Описание
   структур КИХ- и БИХ-фильтров в MATLAB //
   Компоненты и технологии. 2009. № 1.
9. Солонина А. Моделирование цифровой обработки сигналов в MATLAB. Часть 4. Моделирование структур цифровых фильтров с фиксированной точкой программными средствами
   MATLAB: анализ характеристик КИХ-фильтров // Компоненты и технологии. 2009. № 2.
10. Солонина А. Моделирование цифровой обработки сигналов вMATLAB. Часть 5. Моделирование структур цифровых фильтров с фиксированной точкой программными средствами
    MATLAB: анализ характеристик БИХ-фильтров // Компоненты и технологии. 2009. № 3.