Проектирование КИХ фильтра методом взвешивания

Методическое указание к лабораторной работе на учебном стенде LESO7. Шауэрман Александр А. shamrel@yandex.ru

Цель работы
Теоретические сведения
Ход работы
Вопросы для самопроверки
Содержание отчета

1 Цель работы

Научиться проектировать цифровой нерекурсивный фильтр (фильтр с конечной импульсной характеристикой) методом взвешивания, исследовать аппаратную реализацию фильтра на базе ПЛИС.

2 Теоретические сведения

2.1 Аппаратная реализация фильтра

Цифровые фильтры могут быть созданы двумя методами: БИХ и КИХ. Фильтры БИХ (Бесконечная Импульсная Характеристика) – это фильтры в которых значение на выходе зависит от значения на входе и от предыдущих значений на выходе. Это фильтры с обратной связью. В противоположность к ним КИХ (Конечная Импульсная Характеристика) фильтры не используют обратную связь и значение на их выходе связано только с текущими и предыдущими значениями на входе. Что касается стабильности, то КИХ фильтры всегда стабильны, с другой стороны, они нуждаются в более высоком порядке, чтобы соответствовать тем же параметрам, что и БИХ.

Разностное уравнение, описывающее связь между входным и выходным сигналами фильтра для БИХ фильтра, имеет вид:

$y(k) = \sum_{n=0}^{P-1}b(n)\cdot x(k-n)-\sum_{n=0}^{Q-1}a(n)\cdot x(k-n),$

(1)

где P – порядок входного сигнала, Q – порядок обратной связи, b(n) –коэффициенты входного сигнала, a(n) – коэффициенты обратной связи, x(n) – входной сигнал, y(n) – выходной сигнал.

На рисунке 1 показана структурная схема БИХ фильтра, реализующая уравнение (4.1).

Рисунок 1 – Схема БИХ фильтра

Разностное уравнение для КИХ фильтра можно получить, приравняв все коэффициенты обратной связи к нулю:

$y(k)=\sum_{n=0}^{P-1}b(n)\cdot x(k-n)$

(2)

Структурная схема КИХ фильтра показана на рисунке ниже:

Рисунок 2 – Схема КИХ фильтра

Для того чтобы разработать описание схемы фильтра на языке Verilog, присвоим имена всем промежуточным цепям фильтра. В качестве линии задержки будем использовать цепочку параллельных регистров r. Рассмотрим КИХ фильтр с числом коэффициентов N = 5 (Рисунок 3).

Рисунок 3 – Схема КИХ фильтра с подписанными цепями

Отсчет входного сигнала data_i поступает на вход линии задержки r, значение на выходе каждого элемента умножается с соответствующим коэффициентом h. Результаты умножения m с помощью многовходового сумматора складываются, полученное значение data_o – выход фильтра.

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

  always @(posedge clk_i) begin
    r[0] <= data_i;
    r[1] <= r[0];
    r[2] <= r[1];
    r[3] <= r[2];
  end

Передача сигнала по цепочке происходит в поведенческом блоке always. Обратите внимание, что запись значений должна происходить одновременно, поэтому в коде использовано неблокирующее присваивание.

Умножаем отсчеты сигнала на коэффициенты фильтра:

  assign    m[0] = data_i * h[0];
  assign    m[1] = r[0] * h[1];
  assign    m[2] = r[1] * h[2];
  assign    m[3] = r[2] * h[3];
  assign    m[4] = r[3] * h[4];

Реализуем многовходовый сумматор:

  always @(posedge clk_i) 
    data_o <= m[0] + m[1] + m[2] + m[3] + m[4];

Результат суммирования поступает на выход data_o.

Значения коэффициентов фильтра проще всего задать в блоке инициализации (initial):

  initial begin
    h[0] = -53;	// пример ФНЧ
    h[1] = 138;
    h[2] = 255;
    h[3] = 138; 
    h[4] = -53;
  end

Осталось объявить все цепи и регистры, и оформит фильтр в качестве параметризированного модуля:

Модуль фильтра КИХ

(* multstyle = "dsp" *) module  simple_fir
#(
  parameter     N = 5,      // количество коэффициентов 
  parameter     H_W = 9,    // разрядность h 
  parameter     DI_W = 8,   // разрядность входных данных
  parameter     DO_W = 18   // разрядность выходных данных
)
(
  input wire                    clk_i,
  input wire  signed [DI_W-1:0] data_i,
  output  reg        [DO_W-1:0] data_o 
);
// Линию задержки реализуем в виде массива регистров,
// разрядность слова DI_W, число элементов N-1 (элементов задержки
// должно быть на единицу меньше, чем число отсчетов в h).   
  reg  signed [DI_W-1:0] r [N-2:0];   // N-1 регистров                     
 
  // Коэффициенты фильтра
  reg  signed [H_W-1:0] h [N-1:0];   // N регистров  
 
  // Результат умножения
  wire  signed [H_W + DI_W -1:0] m [N-1:0];
 
  // Прописываем значения отсчетов h
  initial begin
    h[0] = -53;
    h[1] = 138;
    h[2] = 255;
    h[3] = 138; 
    h[4] = -53;
  end
 
  // Перемещаем входной сигнал по линии задержки 
  always @(posedge clk_i) begin
    r[0] <= data_i;
    r[1] <= r[0];
    r[2] <= r[1];
    r[3] <= r[2];
  end
 
  // Умножаем на отсчеты h
  assign    m[0] = data_i * h[0];
  assign    m[1] = r[0] * h[1];
  assign    m[2] = r[1] * h[2];
  assign    m[3] = r[2] * h[3];
  assign    m[4] = r[3] * h[4];
 
  always @(posedge clk_i) 
    data_o <= m[0] + m[1] + m[2] + m[3] + m[4];  
 
endmodule

Следует обратить внимание на атрибут (* multstyle = "dsp" *) перед заглавием модуля. Этот атрибут указывает компилятору Quartus, что для синтеза данного модуля следует использовать аппаратные множители (dsp-блоки).

Ниже на листинге экземпляр модуля фильтра в модуле верхнего уровня:

Модуль верхнего уровня

  wire   signed  [7:0] signal1;   // Сигнал с АЦП  
  wire   signed  [17:0] signal2;  // Сигнал после фильтра
 
  wire clk_pll;    // объявляем линию для выхода pll
 
  pll  pll_inst (           // экземпляр модуля pll
    .inclk0 ( clk_50MHz ),  // На вход модуля pll подаем сигнал с тактового генератора
    .c0 ( clk_pll )         // Подключаем линию clk_pll к выходу pll
  );
 
  assign signal1 = $signed (adc_d -128);
 
  // Экземпляр модуля фильтра
  simple_fir 
  #(
    .N(5),    // количество коэффициентов фильтра
    .H_W(9),  // разрядность коэффициентов
    .DI_W(8), // разрядность входных данных
    .DO_W(18) // разрядность выходных данных
    ) fir_inst
  (
    .clk_i (clk_pll),
    .data_i (signal1),
    .data_o (signal2)
   ); 
 
  always @(posedge clk_pll)
  begin
    dac1_d <= (signal1<<6) + $signed({1'b1, {(13){1'b0}}});
    dac2_d <= (signal2[17-:14]) + $signed({1'b1, {(13){1'b0}}});
  end
 
  assign dac1_clk = ~clk_pll;  // Подаем тактовые импульсы на вход CLK ЦАП
  assign dac1_wrt = ~clk_pll;  // Подаем тактовые импульсы на вход WRT ЦАП
 
  assign dac2_clk = ~clk_pll;  // Подаем тактовые импульсы на вход CLK ЦАП
  assign dac2_wrt = ~clk_pll;  // Подаем тактовые импульсы на вход WRT ЦАП
 
  assign adc_pwrdwn = 0;        // Переводим АЦП в рабочий режим
  assign adc_encode = clk_pll;  // Подаем тактовые импульсы на АЦП
 
endmodule

В Verilog все арифметические операции контекстно-зависимые, это значит, что разрядность синтезируемого элемента (сумматор, умножитель) зависит от разрядности операндов. Если требуется перемножить 8-ми битное число a и 6-ти битное число b, то число b будет дополнено до 8-ми битного и будет использована 8-ми битная арифметика. Однако наибольше число, полученное в результате такого умножения, будет занимать не более 14-ти разрядов (8 + 6 = 14). Поэтому для хранения результата перемножения 8-ми битного и 6-ти битных чисел достаточно 14-ти битной переменной (и это только в том случае, если перемножены максимальные по модулю отрицательные числа, в остальных во всех случаях достаточно 13 бит). Разрядность результата перемножения n-битного и k-битного чисел будет равна n + k.

При сложении n-битного и k-битного чисел, компилятор выберет максимальную разрядность операнда и увеличит ее на единицу, таким образом разрядность результата сложения будет равна:

(3)

где max – операция взятия наибольшего.

В реализации КИХ фильтра требуется сложить N чисел, как узнать разрядность результата? Предположим, что все N слагаемых одинаковы, тогда их сложение эквивалентно умножению одного из них на число N. Следовательно разрядность р результата сложения N n-разрядных чисел можно вычислить по формуле:

(4)

где ceil – округление до большего целого.

В примере модуля фильтра разрядность входных данных равна восьми (задается параметром DI_W = 8), разрядность отсчетов фильтра равна девяти (H_W = 9). Разрядность цепей m равна H_W + DI_W:

wire  signed [H_W + D_W -1:0] m [N-1:0];

Разрядность выходных данных определяет параметр DO_W и на основании (4) может быть найден как:

$DO\_W=H\_W+DI\_W+ceil(〖log〗_2 N)$

(5)

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

$DO\_W=H\_W+DI\_W+ceil(log_2 \alpha)-1$

(6)

где

$\alpha=\sum_{i=0}^{N-1}\frac{|b(n)|}{max(|b(n)|)}$

(7)

Результат синтеза модуля можно изучить в RTL Viewer (Рисунок 4).

Рисунок 4 – Модуль КИХ фильтра. Quartus RTL Viewer

2.2 Масштабируемый модуль фильтра

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

Рассмотрим листинг линии задержки:

  always @(posedge clk_i) begin
    r[0] <= data_i;
    r[1] <= r[0];
    r[2] <= r[1];
    r[3] <= r[2];
  end

Линия содержит 4 ячейки и предназначена для фильтра с 5-ю коэффициентами. Очевидно, что если понадобится сделать фильтр с N = 6, то придется дописать строку:

  r[4] <= r[3];

Это не является проблемой, пока число N не велико, но с увеличением N до нескольких десятков это не только увеличит код и уменьшит его читаемость, но и может стать источником серьезных ошибок. Для генерации линии задержки разумно использовать оператор for. Оператор for известен многим из классических языков программирования и предназначен для организации циклов. Синтез использования практически аналогичен с синтаксисом C/C++:

for ( initial_assignment ; expression ; step_assignment )
     function_statement

где initial_assignment – оператор назначения, который выполняется один раз при старте for; expression – выражение, значение, которое проверяется перед выполнением каждой итерации; step_assignment – оператор назначения, который выполняется в конце каждой итерации.

При старте оператора for первым выполняется оператор initial_assignment, который присваивает начальное значение переменной цикла. Далее проверяется условие (выражение expression) и если его значение истина, то выполняется function_statement, иначе цикл завершается. Проверка условия происходит перед началом каждой итерации. После выполнения function_statement выполняется step_assignment. Как правило step_assignment присваивает новое значение переменной цикла. Далее цикл повторяется.

Введем переменную цикла i типа integer, и выполним генерацию цепочки с помощью цикла:

  integer i;
  always @(posedge clk_i) begin
    for(i = 0; i < N; i = i + 1)
      r[i] <= r[i-1];
  end

В результате компилятор выдаст ошибку, так как элемент r[i-1] не определен для i = 0, действительно, в нулевой регистр должен записываться входной отсчет (data_i). Изменим код следующим образом:

  always @(posedge clk_i) begin
    for(i = 0; i < N-1; i = i + 1)
      if (i == 0)
        r[i] <= data_i;
      else
        r[i] <= r[i-1];
  end

Теперь при изменении параметра N компилятор автоматически сформирует код с нужным количеством присваивания.

Аналогичным образом с помощь цикла for оформим умножители:

Было:

  assign    m[0] = data_i * h[0];
  assign    m[1] = r[0] * h[1];
  assign    m[2] = r[1] * h[2];
  assign    m[3] = r[2] * h[3];
  assign    m[4] = r[3] * h[4];

Стало:

  always @(*) begin
    for(i = 0; i < N; i = i + 1)  
      if (i == 0)
        m[i] <= data_i * h[i];
      else
        m[i] <= r[i-1] * h[i];
  end

Заметим, в блоке always звездочка означает, что список чувствительности включает все входные сигналы данного оператора. Такой формат чувствительности используется для описания комбинационных схем. В этом случае m должен быть объявлен как регистр.

Преобразуем участок кода с помощь цикла for:

  always @(posedge clk_i) 
    data_o <= m[0] + m[1] + m[2] + m[3] + m[4];

Разобьем суммирование на шаги и введем вспомогательную переменную acc:

  reg  signed    [DO_W-1:0]    acc;
  always @(*) begin
    for(i = 0; i < N; i = i + 1)
      if (i == 0)
        acc = m[i];
      else
        acc = acc + m[i];
  end

Ниже на листинге представлен модуль КИХ фильтра с параметрической генерацией кода:

Модуль фильтра

(* multstyle = "dsp" *) module  fir
#(
  parameter     N = 5,      // количество коэффициентов 
  parameter     H_W = 9,    // разрядность h 
  parameter     DI_W = 8,   // разрядность входных данных
  parameter     DO_W = 18   // разрядность выходных данных
)
(  
  input wire                    clk_i,
  input wire  signed [DI_W-1:0] data_i,
  output  reg        [DO_W-1:0] data_o 
);
// Линию задержки реализуем в виде массива регистров,
// разрядность слова DI_W, число элементов N-1 (элементов задержки
// должно быть на единицу меньше, чем число отсчетов в h).   
  reg  signed [DI_W-1:0] r [N-2:0];   // N-1 регистров                     
 
  // Коэффициенты фильтра
  reg  signed [H_W-1:0] h [N-1:0];   // N регистров  
 
  // Результат умножения
  reg  signed [H_W + DI_W -1:0] m [N-1:0];
 
  initial begin
    h[0] = -53;
    h[1] = 138;
    h[2] = 255;
    h[3] = 138; 
    h[4] = -53;
  end
 
  integer i;
  reg  signed    [DO_W-1:0]    acc;
 
  always @(posedge clk_i) begin
    for(i = 0; i < N-1; i = i + 1)
      if (i == 0)
        r[i] <= data_i;
      else
        r[i] <= r[i-1];
  end
 
  always @(*) begin
    for(i = 0; i < N; i = i + 1)  
      if (i == 0)
        m[i] <= data_i * h[i];
      else
        m[i] <= r[i-1] * h[i];
  end
 
  always @(*) begin
    for(i = 0; i < N; i = i + 1)
      if (i == 0)
        acc = m[i];
      else
        acc = acc + m[i];
  end
 
 
always @(posedge clk_i) 
    data_o = acc; 
 
endmodule

Стоит отметит, что коэффициенты фильтра по-прежнему должны быть заданы вручную.

2.3 Расчет коэффициентов фильтра методом взвешивания

КИХ фильтр может быть однозначно задан коэффициентами импульсной характеристики {h(n)}. Дискретное преобразование импульсной характеристики (ДПФ) представляет собой выборку частотной характеристики фильтра {H(k)}. Связь между двумя последовательностями однозначная и может быть записана формулами:

$H(k)=\sum_{n=0}^{N-1}h(n)e^\frac{-j\cdot 2\pi\cdot n\cdot k}{N}$ ДПФ;

(8)

$h(n)=\frac{1}{N}\cdot\sum_{k=0}^{N-1}H(k)e^\frac{j\cdot 2\pi\cdot n\cdot k}{N},$ ОДПФ;

(9)

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

$H(\Omega)=\sum_{n=0}^{N-1}h(n)e^{-j\cdot 2\pi\cdot n\cdot \Omega}$

(10)

При конечном значении N аппроксимирующая функция H(Ω) не обеспечивает хорошей аппроксимации заданной АЧХ, возникает хорошо известное явление Гиббса, которое проявляется в виде выбросов и пульсаций определенного уровня около точек разрыва аппроксимирующей частотной функции. Один из методов борьбы с этим явлением – это изменение значения коэффициентов импульсной характеристики с помощью весовой последовательности w(k). Последовательность должна быть выбрана таким образом, чтобы полученная импульсная характеристики hw(k), будучи подставленной в выражение (10), обеспечивала бы требования к АЧХ. Модификация импульсной характеристики может быть произведена по формуле:

$hw_k=h_k\cdot w_k$

(11)

где w_k – последовательность конечной длины (окно).

Метод проектирования КИХ фильтра с применением окон называется методом взвешивания.

Желательно, чтобы окно обладало следующими свойствами:

Главный лепесток частотной характеристики окна должен содержать по возможности большую часть энергии спектра, а его ширина должна быть по возможности минимальной.
Энергия в боковых лепестках должна быстро уменьшаться при приближении Ω к 0.5.

Рассмотрим окна, применяемые при синтезе фильтров.

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

$w_k = \begin{cases} 1, & \text{if } 0\leq k\leq N-1,\\ 0, & \text{if } k < 0, k > N-1,\\ \end{cases}$

(12)

Рисунок 5 – Прямоугольное окно для случая N=13. а) Отсчеты, б) Модуль частотной характеристики

2. Окно Бартлетта или треугольное окно. Описывается выражением:

w_{k} = {\begin{cases} \frac{2 \cdot k}{N - 1}, & 0 \leq k \leq \frac{N - 1}{2}, \\ 2 - \frac{2 \cdot k}{N - 1}, & \frac{N - 1}{2} < k \leq N - 1, \\ 0 & k < 0, k > N - 1, \end{cases}

(13)

Рисунок 6 – Треугольное окно для случая N=13. а) Отсчеты, б) Модуль частотной характеристики

3. Окно Хемминга. Описывается выражением:

$w_k = \begin{cases} 0.54-0.46\cdot \cos\left(\cfrac{2\pi \cdot k}{N-1}\right), & \text{} 0\leq k\leq N-1,\\ 0, & \text{} k < 0, k > N-1,\\ \end{cases}$

(14)

Рисунок 7 – Окно Хемминга для случая N=13. а) Отсчеты, б) Модуль частотной характеристики

4. Окно Хеннинга (Хана). Описывается выражением:

> $w_k = \begin{cases} 0.5-0.5\cdot \cos\left(\cfrac{2\pi \cdot k}{N-1}\right), & \text{} 0\leq k\leq N-1,\\ 0, & \text{} k < 0, k > N-1,\\ \end{cases}$

(15)

Рисунок 8 – Окно Хеннинга для случая N=13. а) Отсчеты, б) Модуль частотной характеристики

5. Окно Блэкмана. Описывается выражением:

$w_k = \begin{cases} 0.45-0.5\cdot \cos\left(\cfrac{2\pi \cdot k}{N-1}\right)+ 0.08\cdot \cos\left(\cfrac{4\pi \cdot k}{N-1}\right), & \text{} 0\leq k\leq N-1,\\ 0, & \text{} k < 0, k > N-1.\\ \end{cases}$

(16)

Рисунок 9 – Окно Блэкмана для случая N=13. а) Отсчеты, б) Модуль частотной характеристики

2.4 Алгоритм расчета КИХ фильтра

Рассмотрим алгоритм расчета на примере: пусть требуется рассчитать нерекурсивный фильтр нижних частот с линейной фазой методом взвешивания если полоса пропускания (0 ‑ 16) кГц, полосы непропускания (задержки): (23 ‑ ∞) кГц. Число отсчетов импульсной характеристики N = 13. Вид окна для взвешивания – окно Хемминга. Частота дискретизации – 100 кГц.

Вычисления рекомендуется проводить в прикладных программных пакетах для математических расчетов типа Mathcad, MatLab, SMath Studio, Mathematica. В примере расчеты произведены в среде Mathcad.

1.Аппроксимация частотной характеристики. Определяем аппроксимирующую частотную характеристику с учетом частоты дискретизации:

где f_pass – граничная частота полосы пропускания, f_stop – граничная частота полосы непропускания, f_s – частота дискретизации.

Построим график аппроксимирующей частотной характеристики фильтра:

Рисунок 10 – График аппроксимирующей функции. Mathcad

Следует отметить, что при проектировании цифровых фильтров удобно пользоваться нормированными к частоте дискретизации частотами, поэтому все частотные характеристики фильтров задаются на интервале от 0 до 1. Зададим аппроксимирующую функцию для нормированной частоты Ω, учтем что функция должна быть симметрична относительно f_s/2 (или Ω = 0.5):

Рисунок 11 – Желаемая амплитудно-частотная характеристика фильтра. Mathcad

2. Дискретизация аппроксимирующей функции. Для заданного числа N (количество отсчетов импульсной характеристики) определяем дискретную выборку (вектор) отсчетов модуля частотной характеристики (АЧХ) Hm_k и вектор нормированных частот Ω_k.

$\Omega_k=\cfrac{k}{N},$

$Hm_k=D_{norm}(\Omega_k)$

где k = 0, 1, 2 .. N

Рисунок 12 – График дискретной АЧХ. Mathcad

3. Дискретизация ФЧХ. На данном этапе у нас задана только АЧХ, то есть модуль комплексного коэффициента передачи , чтобы выполнить ОДПФ мы должны добавить ФЧХ. Для фильтра с четным количеством коэффициентов N ФЧХ рассчитывается следующим образом:

$\phi_k = \begin{cases} -\dfrac{N-1}{N}\cdot \pi\cdot k, & \text{} k=0..\dfrac{N}{2}-1,\\ -\dfrac{(N-1)\cdot(k-N)}{N}\cdot \pi, & \text{} k=\dfrac{N}{2}..N-1.\\ \end{cases}$

(17)

Для фильтра с нечетным N выражение примет вид:

$\phi_k = \begin{cases} -\dfrac{N-1}{N}\cdot \pi\cdot k, & \text{} k=0..\dfrac{N-1}{2},\\ -\dfrac{(N-1)\cdot(k-N)}{N}\cdot \pi, & \text{} k=\dfrac{N+1}{2}..N-1.\\ \end{cases}$

(18)

В рассматриваемом примере число отсчетов нечетное и равно 13. В среде Mathcad выражение для расчета ФЧХ может быть записано следующим образом:

Для визуального контроля построим на одном графике полученные отсчеты АЧХ и ФЧХ:

Рисунок 13 – Отсчеты АЧХ (синий) и ФЧХ (красный). Mathcad

Рассчитанные значение сведем в таблицу:

Таблица 1 – Значения отсчетов частотной выборки
Ω	0	0.077	0.154	0.231	0.308	0.385	0.462	0.538	0.615	0.692	0.769	0.846	0.923
Hm	1	1	1	0	0	0	0	0	0	0	0	1	1
φ, рад	0	-2.9	-5.8	-8.7	-11.6	-14.5	-17.4	17.4	14.5	11.6	8.7	5.8	2.9

4. Расчет импульсной характеристики. Используя обратное дискретное преобразование Фурье (ОДПФ) вычисляем коэффициенты импульсной характеристики.

$h_k=\frac{1}{N}\sum_{n=0}^{N-1}Hm_n\cdot e^{j\cdot \phi_n}\cdot e^{j \frac{2\pi\cdot n\cdot k}{N}}$

(19)

В общем случае результат обратного дискретного преобразования Фурье является комплексным. Однако, если дискретизация АЧХ и ФЧХ выполнена правильно, и в результате получена симметричная частотная функция, то мнимая часть в отсчетах импульсной характеристики будет равна нулю. На практике она отличается от нуля ввиду ошибок округления при расчете, но имеет очень маленькие значения, поэтому h_k можно считать вещественной.

Результат расчета сведем в таблицу и визуализируем на графике.

Таблица 2 – Значения отсчетов импульсной характеристики
k	0	1	2	3	4	5	6	7	8	9	10	11	12
h	0.064	-0.02	-0.093	-0.054	0.11	0.301	0.385	0.301	0.11	-0.054	-0.093	-0.02	0.064

Рисунок 14 – Отсчеты импульсной характеристики фильтра. Mathcad

5. Оконное сглаживание. Определяем отсчеты окна Хемминга для N = 13:

$w_k=0.54-0.46\cdot \cos\left(\cfrac{2\pi \cdot k}{N-1}\right)$

Рисунок 15 – Отсчеты оконной функции. Mathcad

Отсчеты взвешенной импульсной функции hm можно рассчитать, умножив каждый отсчет исходной импульсной функции на отсчет оконной функции с тем же порядковым номером:

$hw_k=h_k\cdot w_k$

(20)

Результат расчета сведем в таблицу и визуализируем с помощью графика.

Таблица 3 – Значения отсчетов импульсной характеристики после взвешивания
k	0	1	2	3	4	5	6	7	8	9	10	11	12
w	0.08	0.142	0.31	0.54	0.77	0.938	1	0.938	0.77	0.54	0.31	0.142	0.08
h	0.064	-0.02	-0.093	-0.054	0.11	0.301	0.385	0.301	0.11	-0.054	-0.093	-0.02	0.064
hw	0.005	-0.003	-0.029	-0.029	0.085	0.282	0.385	0.282	0.085	-0.029	-0.029	-0.003	0.005

Рисунок 16 – Отсчеты импульсной характеристики фильтра до (красный) и после (синий) взвешивания. Mathcad

6. По импульсной характеристики восстанавливаем частотную характеристику. Используя связь z-преобразования с частотной функцией, можно восстановить частотную характеристику по отсчетам импульсной характеристики:

$H(\Omega)=\sum_{n=0}^{N-1}h_n\cdot e^{-j\cdot 2 \pi \cdot\Omega\cdot n}}.$

(21)

Построим АЧХ КИХ фильтра с взвешенными и исходными коэффициентами на одном графике:

Рисунок 17 – АЧХ фильтров с взвешиванием (зеленый) и без (синий). Mathcad

7. Определяем длину кодового слова для отсчетов импульсной характеристики.

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

$D=2\cdot \cfrac{\max|h|}{\min|h|}$

(22)

где и max|h| и min|h|– максимальное и минимальное значение модуля отсчета импульсной характеристики. С учетом того, что отсчеты импульсной характеристики принимает отрицательные значения, то отношение максимального значения к минимальному должно быть умножено на два.

$D=2\cdot \cfrac{|0.385|}{|-0.003|}\approx 257.$

Знак приближенного равенства в выражении означает, что значения отсчетов взяты с округлением и результат вычисления не может быть точным, а лишь ориентировочным. Таким образом, для того чтобы обеспечить динамический диапазон значений для хранения отсчетов импульсной характеристики нужно использовать как минимум 9 разрядов (динамический диапазон 512).

Чтобы получить целочисленные коэффициенты фильтра, нормируем импульсную характеристику и полученные значения умножим на 255 (максимальное 9-ти разрядное положительное число):

$b=round\left(\cfrac{h\cdot 255}{\max|h|}\right)$

(23)

где round – операция округления до целого.

Выполняем расчет для коэффициентов до и после взвешивания. Результат сведем в таблицу:

Таблица 4 – Целочисленные значения отсчетов импульсной характеристики до и после взвешивания
k	0	1	2	3	4	5	6	7	8	9	10	11	12
h	0.064	-0.02	-0.093	-0.054	0.11	0.301	0.385	0.301	0.11	-0.054	-0.093	-0.02	0.064
b	42	-13	-62	-36	73	199	255	199	73	-36	-62	-13	42
hw	0.005	-0.003	-0.029	-0.029	0.085	0.282	0.385	0.282	0.085	-0.029	-0.029	-0.003	0.005
bw	3	-2	-19	-19	56	187	255	187	56	-19	-19	-2	3

3 Ход работы

Этап 1. Предварительный расчет

Для заданного варианта рассчитать коэффициенты КИХ фильтра. Результатом расчета должны быть целочисленные коэффициенты фильтра.

Таблица 5 – Задания на проектирование КИХ фильтра
Вар.	Число отсчетов	Полоса пропускания, кГц	Полоса непропускания, кГц	Вид окна для взвешивания
1	13	0 – 8.0	15.0 – ∞	Блэкмана
2	13	0 – 23.1	30.7 – ∞	Хемминга
3	15	0 – 7.0	13.0 – ∞	Хана
4	15	0 – 13.5	20.0 – ∞	Блэкмана
5	15	0 – 20.0	26.6 – ∞	Хемминга
6	15	0 – 26.7	33.3 – ∞	Хана
7	17	0 – 6.0	11.7 – ∞	Блэкмана
8	17	0 – 12.0	17.0 – ∞	Хемминга
9	17	0 – 18.0	23.0 – ∞	Хана
10	17	0 – 23.6	29.4 – ∞	Блэкмана

Расчет произвести для варианта со взвешиванием и без. Расчеты свести в таблицу:

Таблица 6 – Результат расчета фильтра
k	0	1	2	…	N-1
h (импульсная х-ка)
b (отсчеты фильтра)
w (оконная функция)
hw (взвешенная х-ка)
bw (отсчеты фильтра)

Построить график АЧХ и импульсной характеристики фильтра с невзвешенными и с взвешенными коэффициентами.

3.2 Этап 2. Создание проекта, настройка стенда

Создать проект в Quartus Prime.

В модуле верхнего уровня подключить выход АЦП на вход ЦАП. С помощью IP-функции ALTPLL установить частоту дискретизации АЦП и преобразования ЦАП 100 кГц.

Скомпилировать проект и загрузить файл конфигурации в учебный стенд LESO7.

Подключить к задействованному в проекте входу АЦП выход генератора LESO5. К выходу ЦАП подключить вход анализатора сигналов LESO4 (Рисунок 18).

Рисунок 18 – Схема стенда для экспериментального исследования КИХ фильтра

Подать гармонический сигнал (размах: 1.8 В, частота: 1000 Гц, смещение: 0 В), убедиться, что на выходе ЦАП сигнал без искажений. Наличие искажений контролировать визуально по осциллограмме и по спектру (паразитные гармоники не должны быть выше -40 дБ.).

Подать импульсный сигнал.

Установить шаг частоты импульсов равный Δf = 1000 Гц , соответствующий период следования рассчитать по формуле:

$T=\cfrac{1}{\Delta f}$

Длительность импульса выбрать равным периоду дискретизации АЦП. В этом случае спектр сигнала в диапазоне частот от 0 до 0.5f_д, будет наиболее равномерным.

На Рисунке 19 показано окно генератора с настройками импульсов для исследования цепи с частотой дискретизации 100 кГц. На Рисунке 20 и Рисунке 21 показаны осциллограмма и спектр полученного сигнала.

Рисунок 19 – Настройки генератора LESO5. Импульсный сигнал

Рисунок 20 – Осциллограмма импульсного сигнала

Рисунок 21 – Спектр импульсного сигнала

Полученные настройки генератора следует использовать для исследования АЧХ цепей в последующих этапах работы.

Листинг кода для первого этапа

module leso7_fir(
 
  (* chip_pin = "23" *) 
  input          clk_50MHz,  // Вход с тактового генератора
 
  (* chip_pin = "136, 135, 133, 132, 129, 128, 127, 126" *) 
  input   [7:0]    adc_d,    // Шина данных АЦП
 
  (* chip_pin = "137" *)
  output         adc_encode, // Сигнал ENCODE АЦП
 
  (* chip_pin = "125" *)
  output        adc_pwrdwn, // Сигнал PWRDWN АЦП
 
  (* chip_pin = "106, 105, 104, 103, 100, 99, 98, 87, 86, 85, 84, 83, 80, 77" *)
  output  reg [13:0]  dac1_d,    // Шина данных ЦАП1
 
  (* chip_pin = "76" *)
  output        dac1_wrt,    // Сигнал WRT ЦАП1
 
  (* chip_pin = "75" *)
  output        dac1_clk,    // Сигнал CLK ЦАП1
 
  (* chip_pin = "72, 71, 70, 69, 68, 67, 66, 65, 60, 59, 58, 55, 54, 53" *)
  output reg  [13:0]  dac2_d,    // Шина данных ЦАП2
 
  (* chip_pin = "73" *)
  output        dac2_wrt,    // Сигнал WRT ЦАП2
 
  (* chip_pin = "74" *)
  output        dac2_clk    // Сигнал CLK ЦАП2
  // после объявления последнего порта запятая не ставится
 
);
 
  wire   signed  [7:0] signal1;   // Сигнал с АЦП  
  wire clk_pll;    // объявляем линию для выхода pll
 
  pll  pll_inst (           // экземпляр модуля pll
    .inclk0 ( clk_50MHz ),  // На вход модуля pll подаем сигнал с тактового генератора
    .c0 ( clk_pll )         // Подключаем линию clk_pll к выходу pll
  );
 
  assign signal1 = $signed (adc_d -128);
 
 
  always @(posedge clk_pll)
  begin
    dac1_d <= (signal1<<6) + $signed({1'b1, {(13){1'b0}}});
  end
 
  assign dac1_clk = ~clk_pll;  // Подаем тактовые импульсы на вход CLK ЦАП
  assign dac1_wrt = ~clk_pll;  // Подаем тактовые импульсы на вход WRT ЦАП
 
  assign dac2_clk = ~clk_pll;  // Подаем тактовые импульсы на вход CLK ЦАП
  assign dac2_wrt = ~clk_pll;  // Подаем тактовые импульсы на вход WRT ЦАП
 
  assign adc_pwrdwn = 0;        // Переводим АЦП в рабочий режим
  assign adc_encode = clk_pll;  // Подаем тактовые импульсы на АЦП
 
endmodule

3.3 Этап 3. Исследование КИХ фильтра без взвешивания

В модуль верхнего уровня вставить экземпляр модуля фильтра с учетом разрядности полученных коэффициентов и выходного сигнала.

В модуле фильтра изменить коэффициенты на невзвешенные коэффициенты фильтра согласно предварительному расчету.

Сигнал с АЦП должен быть переведен в знаковую форму.

На первый выход ЦАП подать сигнал с АЦП, на второй сигнал подать сигнал с выхода фильтра. Сигналы подавать с учетом масштабирования.

Листинг модуля фильтра:

(* multstyle = "dsp" *) module  fir
#(
  parameter     N = 5,      // количество коэффициентов 
  parameter     H_W = 9,    // разрядность h 
  parameter     DI_W = 8,   // разрядность входных данных
  parameter     DO_W = 18   // разрядность выходных данных
)
(  
  input wire                    clk_i,
  input wire  signed [DI_W-1:0] data_i,
  output  reg        [DO_W-1:0] data_o 
);
// Линию задержки реализуем в виде массива регистров,
// разрядность слова DI_W, число элементов N-1 (элементов задержки
// должно быть на единицу меньше, чем число отсчетов в h).   
  reg  signed [DI_W-1:0] r [N-2:0];   // N-1 регистров                     
 
  // Коэффициенты фильтра
  reg  signed [H_W-1:0] h [N-1:0];   // N регистров  
 
  // Результат умножения
  reg  signed [H_W + DI_W -1:0] m [N-1:0];
 
  initial begin
    h[0] = 42;
    h[1] = -13;
    h[2] = -62;
    h[3] = -36; 
    h[4] = 73;
    h[5] = 199;
    h[6] = 255;
    h[7] = 199;
    h[8] = 73; 
    h[9] = -36;
    h[10] = -62;
    h[11] = -13;
    h[12] = 42;
 
  end
 
  integer i;
  reg  signed    [DO_W-1:0]    acc;
 
  always @(posedge clk_i) begin
    for(i = 0; i < N-1; i = i + 1)
      if (i == 0)
        r[i] <= data_i;
      else
        r[i] <= r[i-1];
  end
 
  always @(*) begin
    for(i = 0; i < N; i = i + 1)  
      if (i == 0)
        m[i] <= data_i * h[i];
      else
        m[i] <= r[i-1] * h[i];
  end
 
  always @(*) begin
    for(i = 0; i < N; i = i + 1)
      if (i == 0)
        acc = m[i];
      else
        acc = acc + m[i];
  end
 
always @(posedge clk_i) 
    data_o = acc; 
 
endmodule

Листинг модуля верхнего уровня:

module leso7_fir(
 
  (* chip_pin = "23" *) 
  input          clk_50MHz,  // Вход с тактового генератора
 
  (* chip_pin = "136, 135, 133, 132, 129, 128, 127, 126" *) 
  input   [7:0]    adc_d,    // Шина данных АЦП
 
  (* chip_pin = "137" *)
  output         adc_encode, // Сигнал ENCODE АЦП
 
  (* chip_pin = "125" *)
  output        adc_pwrdwn, // Сигнал PWRDWN АЦП
 
  (* chip_pin = "106, 105, 104, 103, 100, 99, 98, 87, 86, 85, 84, 83, 80, 77" *)
  output  reg [13:0]  dac1_d,    // Шина данных ЦАП1
 
  (* chip_pin = "76" *)
  output        dac1_wrt,    // Сигнал WRT ЦАП1
 
  (* chip_pin = "75" *)
  output        dac1_clk,    // Сигнал CLK ЦАП1
 
  (* chip_pin = "72, 71, 70, 69, 68, 67, 66, 65, 60, 59, 58, 55, 54, 53" *)
  output reg  [13:0]  dac2_d,    // Шина данных ЦАП2
 
  (* chip_pin = "73" *)
  output        dac2_wrt,    // Сигнал WRT ЦАП2
 
  (* chip_pin = "74" *)
  output        dac2_clk    // Сигнал CLK ЦАП2
  // после объявления последнего порта запятая не ставится
 
);
 
  wire   signed  [7:0] signal1;   // Сигнал с АЦП  
  wire   signed  [18:0] signal2;  // Сигнал после фильтра
 
  wire clk_pll;    // объявляем линию для выхода pll
 
  pll  pll_inst (           // экземпляр модуля pll
    .inclk0 ( clk_50MHz ),  // На вход модуля pll подаем сигнал с тактового генератора
    .c0 ( clk_pll )         // Подключаем линию clk_pll к выходу pll
  );
 
  assign signal1 = $signed (adc_d -128);
 
  // Экземпляр модуля фильтра
  fir 
  #(
    .N(13),   // количество коэффициентов фильтра
    .H_W(9),  // разрядность коэффициентов
    .DI_W(8), // разрядность входных данных
    .DO_W(19) // разрядность выходных данных
    ) fir_inst
  (
    .clk_i (clk_pll),
    .data_i (signal1),
    .data_o (signal2)
   ); 
 
 
  always @(posedge clk_pll)
  begin
    dac1_d <= (signal1 << 6) + $signed({1'b1, {(13){1'b0}}});
    dac2_d <= (signal2[18-:14] << 2) + $signed({1'b1, {(13){1'b0}}});
  end
 
  assign dac1_clk = ~clk_pll;  // Подаем тактовые импульсы на вход CLK ЦАП
  assign dac1_wrt = ~clk_pll;  // Подаем тактовые импульсы на вход WRT ЦАП
 
  assign dac2_clk = ~clk_pll;  // Подаем тактовые импульсы на вход CLK ЦАП
  assign dac2_wrt = ~clk_pll;  // Подаем тактовые импульсы на вход WRT ЦАП
 
  assign adc_pwrdwn = 0;        // Переводим АЦП в рабочий режим
  assign adc_encode = clk_pll;  // Подаем тактовые импульсы на АЦП
 
endmodule

Скомпилировать проект и загрузить файл конфигурации в учебный стенд LESO7.

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

Сохранить осциллограмму импульсной характеристики (Рисунок 22) и спектр импульсной характеристики в линейном (Рисунок 23) и логарифмическом масштабах (Рисунок 24).

Рисунок 22 – Импульсная характеристика фильтра. Зеленый трек – вход фильтра, желтый трек – выход фильтра

Рисунок 23 – Спектр импульсной характеристики фильтра. Линейный масштаб

Рисунок 24 – Спектр импульсной характеристики фильтра. Логарифмический масштаб

По спектру в логарифмическом масштабе определить затухание в полосе пропускания на границах диапазона и в полосе задержки. В качестве нулевого уровня выбрать максимальный уровень спектра в полосе пропускания (на Рисункне 24 – это -15.71 дБ). Для измерения уровня удобно использовать инструмент «курсор».

С помощью RTL-Viewer изучить схему реализации фильтра. Для этого открыть RTL-Viewer, выбрать блок fir:fir_inst, раскрыть его (значок «+»).

Для следующих пунктов работы следует сохранить бинарный файл конфигурации.

Сделать выводы.

3.4 Этап 4. Исследование КИХ фильтра с взвешиванием

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

Сохранить осциллограмму импульсной характеристики (Рисунок 25) и спектр импульсной характеристики в линейном (Рисунок 26) и логарифмическом масштабах (Рисунок 27).

Рисунок 25 – Импульсная характеристика фильтра. Зеленый трек – вход фильтра, желтый трек – выход фильтра

Рисунок 26 – Спектр импульсной характеристики фильтра. Линейный масштаб

Рисунок 27 – Спектр импульсной характеристики фильтра. Логарифмический масштаб

Для следующих пунктов работы следует сохранить бинарный файл конфигурации.

Сделать выводы.

3.5 Этап 5. Исследование прохождения сигнала через цифровой фильтр

Загрузить в ПЛИС файл конфигурации фильтра с коэффициентами без взвешивания.

Поочередно подавать с генератора LESO5 на вход стенда сигналы разных форм, наблюдать формы и спектр сигналов на входе и выходе фильтра. Исследуемые формы: синусоидальная, прямоугольная. Размах сигнала для каждой формы установить 1 В.

Для каждой формы произвести исследование на частотах f_s/100 и f_s/20, где f_s– частота дискретизации.

При исследовании спектра время развертки должно быть выбрано таким образом, чтобы по диапазон частот был от 0 Гц до f_s/2.

Сохранить скриншоты осциллограмм (Рисунок 28, Рисунок 29, Рисунок 30, Рисунок 31) и спектров сигналов.

Рисунок 28 – Сигнал синусоидальной (гармонической) формы (1 кГц) на входе и выходе фильтра. Зеленый трек – вход фильтра, желтый трек – выход фильтра

Рисунок 29 – Сигнал синусоидальной (гармонической) формы (5 кГц) на входе и выходе фильтра. Зеленый трек – вход фильтра, желтый трек – выход фильтра

Рисунок 30 – Сигнал прямоугольной формы (1 кГц) на входе и выходе фильтра. Зеленый трек – вход фильтра, желтый трек – выход фильтра

Рисунок 31 – Сигнал прямоугольной формы (5 кГц) на входе и выходе фильтра. Зеленый трек – вход фильтра, желтый трек – выход фильтра

Повторить исследование для фильтра с взвешенными коэффициентами.

Рисунок 32 – Сигнал прямоугольной формы (1 кГц) на входе и выходе фильтра с взвешенными коэффициентами. Зеленый трек – вход фильтра, желтый трек – выход фильтра

Рисунок 33 – Сигнал прямоугольной формы (5 кГц) на входе и выходе фильтра с взвешенными коэффициентами. Зеленый трек – вход фильтра, желтый трек – выход фильтра

Сделать выводы.

4 Вопросы для самопроверки:

Дать определение КИХ фильтра.
Назвать причины синтеза КИХ фильтров с линейной фазой.
Дать определение метода окон (взвешивания). На примере фильтров (ФНЧ, ФВЧ, ПФ, РФ) перечислить этапы синтеза КИХ фильтров с использованием метода взвешивания.
Что такое «окно»? Привести примеры «окон».
ДПФ и его применение при проектировании КИХ фильтра.
В чем суть явление Гиббса?

5 Содержание отчета:

Название и цель работы.
Структурная схема рассчитанного КИХ фильтра.
Таблица с расчетами коэффициентов фильтра.
Расчет разрядности коэффициентов фильтра и выходного сигнала.
Рассчитанные и измеренные графики АЧХ фильтра с взвешенными и с невзвешенными коэффициентами.
Исходные коды модуля фильтра и модуля верхнего уровня.
Осциллограммы и скриншоты измерений (для пунктов, где это указано).
Для каждого этапа сделать выводы.

Литература

16 июня 2019

Добавить комментарий

При компиляции проекта жалуется синтаксисе в pll на верхнем уровне. Можно подробнее узнать о этом модуле и как исправить ошибку. Прикрепляю скрин с ошибкой компилятора

Прикрепленные файлы:

scrin.png

Заговори, чтобы я тебя увидел.

ответить

shamrel
сб, 11/30/2019 - 18:16
Постоянная ссылка (Permalink)

Приветствую! Первое, что бросается в глаза, это то, что файл имеет расширение .sv, это значит, что вы используете SystemVerilog, в то время как все примеры были в Verilog. Но это не должно было вызвать проблем с пропущенной ";". Следующая странность, судя по логам в скриншоте, происходит компиляция файла verilog5.sv, а ошибка в файле Kinf.sv. Не очень понятно, какая архитектура вашего проекта.

Вам нужно разобраться с алгоритмом создания проекта. Рекомендую почитать статю Пишем "демку" для LESO2 на Verilog

ответить

dimatsss
сб, 11/30/2019 - 18:30
Постоянная ссылка (Permalink)

Файл переименовал в sv по причине жалоб Quartus

Прикрепленные файлы:

scrin_2.png

Заговори, чтобы я тебя увидел.

ответить

shamrel
сб, 11/30/2019 - 19:38
Постоянная ссылка (Permalink)

В приложении проект фильтра. Там много лишнего. При желании можно разобраться. Открывать Quartus Prime 18

Прикрепленные файлы:

Проект фильтра

ответить

Сообщение об ошибке

Новости ЛЭСО

Контактные лица

Проектирование КИХ фильтра методом взвешивания

1 Цель работы

2 Теоретические сведения

2.1 Аппаратная реализация фильтра

2.2 Масштабируемый модуль фильтра

2.3 Расчет коэффициентов фильтра методом взвешивания

2.4 Алгоритм расчета КИХ фильтра

3 Ход работы

Этап 1. Предварительный расчет

3.2 Этап 2. Создание проекта, настройка стенда

3.3 Этап 3. Исследование КИХ фильтра без взвешивания

3.4 Этап 4. Исследование КИХ фильтра с взвешиванием

3.5 Этап 5. Исследование прохождения сигнала через цифровой фильтр

4 Вопросы для самопроверки:

5 Содержание отчета:

Литература

Комментарии:

Похожие статьи

Сообщение об ошибке

Вы здесь

Форма поиска

Новости ЛЭСО

Контактные лица

Проектирование КИХ фильтра методом взвешивания

1 Цель работы

2 Теоретические сведения

2.1 Аппаратная реализация фильтра

2.2 Масштабируемый модуль фильтра

2.3 Расчет коэффициентов фильтра методом взвешивания

2.4 Алгоритм расчета КИХ фильтра

3 Ход работы

Этап 1. Предварительный расчет

3.2 Этап 2. Создание проекта, настройка стенда

3.3 Этап 3. Исследование КИХ фильтра без взвешивания

3.4 Этап 4. Исследование КИХ фильтра с взвешиванием

3.5 Этап 5. Исследование прохождения сигнала через цифровой фильтр

4 Вопросы для самопроверки:

5 Содержание отчета:

Литература

Комментарии:

Похожие статьи