Страница 1 из 1

АЧХ полученного КИХ фильтра не соответсвует действительности.

Добавлено: 17 мар 2023, 12:52
harkut
Параметры КИХ фильтра:
1) Частота дискретизации - 10 МГц;
2) Частота среза - 1.1 МГц;
3) Симметричный.

Получил коэффициенты при помощи питоновского SciPy.

Код на питоне:

Код: Выделить всё

#-----------------------------------------------------------------------------
# Библиотеки.
from scipy import signal
import matplotlib.pyplot as plt
import matplotlib
import numpy as np
import math
#-----------------------------------------------------------------------------
# Функции.

# FirFiniteInt(koeffs, precision).
# Получает koeffs и precision. Возвращает koeffs_int с разрядностью precision.
# Аргументы:
# 1)koeffs    - исходные коэффициенты фильтра;
# 2)precision - разрядность итоговых коэффициентов.
def FirFiniteInt(koeffs, precision):
 min_koeff = abs(min(koeffs))
 max_koeff = abs(max(koeffs))
 scale = 2**(precision - 1) / min_koeff if min_koeff > max_koeff\
  else (2**(precision - 1) - 1) / max_koeff
 print(scale)
 koeffs_int = []
 for item in koeffs:
  koeffs_int.append(round(item * scale))
 return koeffs_int

# FirFiniteInt(koeffs, precision).
# Получает koeffs и precision. Возвращает koeffs_int с разрядностью precision.
# Аргументы:
# 1)koeffs    - исходные коэффициенты фильтра;
# 2)precision - разрядность итоговых коэффициентов.
def FirFiniteInt2(koeffs, precision):
 sum_koeffs = sum(np.abs(koeffs)) 
 scale = (2**(precision - 1) - 1) / sum_koeffs
 print(scale)
 koeffs_int = []
 for item in koeffs:
  koeffs_int.append(round(item * scale))
 return koeffs_int


# PrintKoeffsVerilog(koeffs, precision, filename).
# Создает filename. Форматирует koeffs и записывает в filename.
# Аргументы:
# 1)koeffs    - коэффициенты фильтра с разрядностью precision;
# 2)precision - разрядность koeffs;
# 3)filename  - название файла в который будут записываться coeffs.
def PrintKoeffsVerilog(koeffs, precision, filename):
 with open(filename, 'w') as f:  
  for i in range(len(koeffs)//2, -1, -1):
   sign = '-' if koeffs[i] < 0 else ''
   comma = ',' if i!=0 else ''
   print(f'{sign}{precision}\'d{abs(koeffs[i])}{comma}', file = f)


def FirOutLength(in_length, koeffs_precision, koeffs):
 k = 0
 max_koeff = abs(max(koeffs))
 for i in koeffs:
  k += i/max_koeff
 alpha = math.ceil(k)
 out_length = in_length + koeffs_precision + alpha
 return out_length
 
#-----------------------------------------------------------------------------
# Формат графиков.

# Шрифты на графиках.
font = {'family' : 'GOST Type BU',
        'weight' : 'bold',
        'size'   : 30}
matplotlib.rc('font', **font)
plt.rc('legend', fontsize=25)

#matplotlib.rcParams["axes.formatter.limits"] = (0, 2)
#-----------------------------------------------------------------------------
# Параметры фильтра.

# Порядок фильтра.
order = 40
# Частота дисктретизации.
fd = 10e6
# Частота начала среза.
fb = 1.1e6
# Частота конца среза.
fe = 1.3e6
# Частотные границы фильтра.
freq_borders = [0,
                fb/(fd/2),
                fe/(fd/2),
                1]
# Уровень сигнала до начала среза.
B1 = 1
# Уровень сигнала после среза.
B2 = 0
# Уровни ослабления сигнала.
level_borders = [B1,
                 B1,
                 B2,
                 B2]
# Разрядность входа фильтра.
in_length = 2
# Разрядность коэффициентов фильтра.
precision = 16

#-----------------------------------------------------------------------------
#Находим коэффициенты фильтра.
y1 = signal.firwin2(order + 1,
                    freq_borders,
                    level_borders,
                    window=('kaiser', 5))
#Целочисленные коэффициенты фильтра.
y1_i = FirFiniteInt2(y1, precision)
sum_y1_i = sum(np.abs(y1_i))
print(f'sum_y1_i = {sum_y1_i}')
out_length = FirOutLength(in_length, precision, y1_i)
#-----------------------------------------------------------------------------
#Формируем файл с коэффициентами для Verilog.
PrintKoeffsVerilog(y1_i, precision, 'data_out/fir_filter_coef_int.vh')

print('------------------' + '\nКоэффициенты фильтра:')
print(y1)
print('------------------' + '\nКоэффициенты с ограниченной точностью:')
print(y1_i)
print('------------------' + '\nПорядок фильтра:')
print(order)
print('------------------' + '\nРазрядность входа фильтра:')
print(in_length)
print('------------------' + '\nРазрядность коэффициентов фильтра:')
print(precision)
print('------------------' + '\nРазрядность выхода фильтра:')
print(out_length)
#-----------------------------------------------------------------------------
#Параметры графика.
fig, axes = plt.subplots(2)

# # Создание форматера
# formatter = matplotlib.ticker.ScalarFormatter()
# formatter.set_powerlimits((-3, 2))

# # Установка форматера для оси X
# ax.yaxis.set_major_formatter(formatter)

az = np.array([1.0,0])
f1,h1 = signal.freqz(y1, az, fs = fd)
f1_i,h1_i = signal.freqz(y1_i, fs = fd)

#for i in range(len(f1)):
# print(f1[i], np.abs(h1[i]))

axes[0].plot(f1,
             (np.abs(h1)),
             '--g',
             linewidth=2,
             label = 'label_1')
axes[1].plot(f1_i,
            (np.abs(h1_i)),
            '-.y',
            linewidth=4,
            label = 'label_2')

for ax in axes.flat:
    #ax.set_xlim(left = 0, right = 200e3)
    ax.set(xlabel='x-label', ylabel='y-label')
    ax.label_outer()
    ax.minorticks_on()
    ax.grid(which='major',
            color = 'b', 
            linewidth = 1)
    plt.grid(which='minor', 
             color = 'g', 
             linestyle = ':')



    
# plt.legend()
# plt.minorticks_on()
# plt.grid(which='major',
#          color = 'b', 
#          linewidth = 1)
# plt.grid(which='minor', 
#         color = 'b', 
#         linestyle = ':')
# axes[0].set_xlabel('Частота, Гц')
# axes[0].set_ylabel('Амплитуда')

plt.show()

Полученное АЧХ и параметры:
Снимок экрана от 2023-03-17 16-51-16.png
Снимок экрана от 2023-03-17 16-50-57.png
Итоговые полученные целочисленные коэффициенты для Verilog:

Код: Выделить всё

16'd5299,
16'd4780,
16'd3415,
16'd1697,
16'd190,
-16'd712,
-16'd916,
-16'd615,
-16'd138,
16'd227,
16'd348,
16'd257,
16'd78,
-16'd65,
-16'd116,
-16'd88,
-16'd31,
16'd12,
16'd26,
16'd18,
16'd6
Реализация фильтра на Verilog:

Код: Выделить всё

//----------------------------------------------------------------------------
module FIR_FILTER (C, D, CLR, Y);

  parameter kOrder                                  = 40;
  parameter kKoeffWidth                             = 16;
  localparam kKoeffAmount = kOrder / 2 + 1;
  localparam [(kKoeffAmount * kKoeffWidth) - 1:0] kKoeffs = {
`include "src/data_out/fir_filter_coef_int.vh" 
                                                           };

  parameter kInWidth                                = 2;
  parameter kOutWidth                               = 23;

  localparam kShiftRegWidth = (kOrder + 1) * kInWidth;

  
  input                     C;
  input [kInWidth-1:0]      D;
  input                     CLR;
  
  output bit signed [kOutWidth-1:0] Y;
  
  reg [kShiftRegWidth-1:0]   st;

  initial begin
    integer i;
    bit signed [kKoeffWidth-1:0] koeff;
    for (i = 0; kKoeffAmount > i; ++i) begin
      koeff = kKoeffs[i * kKoeffWidth+:kKoeffWidth];
      $display("koeff = %d", koeff);
    end
    Y = 0;
  end
  
  // always @(posedge C) st <= {st[kShiftRegWidth-1-kInWidth:0],D};
  always @(posedge C) st <= {D,st[kShiftRegWidth-1:kInWidth]};
  
  always @(posedge C) Y <= fir_summing();

  function bit signed [kOutWidth-1:0] fir_summing;
    integer i;
    bit signed [kInWidth-1:0] st_item;
    bit signed [kKoeffWidth-1:0] koeff;
    // integer st_item;
    fir_summing = 0;
    st_item = st[kOrder / 2 * kInWidth+:kInWidth];
    koeff = kKoeffs[(kKoeffAmount-1)*kKoeffWidth+:kKoeffWidth];
    fir_summing += koeff * st_item;
    //$display("i = %d st_item = %d fir_summing = %d koeff = %d\n",
    //         i, st_item, fir_summing, koeff);
    for (i = 0; kKoeffAmount - 1 > i; ++i) begin : linear_combination
      st_item = st[i * kInWidth+:kInWidth];
      koeff = kKoeffs[i * kKoeffWidth+:kKoeffWidth];
      fir_summing += koeff * st_item;
      //$display("i = %d st_item = %d fir_summing = %d koeff = %d\n",
      //         i, st_item, fir_summing, koeff);

      st_item = st[(kOrder - i) * kInWidth+:kInWidth];
      fir_summing += koeff * st_item;
      //$display("i = %d st_item = %d fir_summing = %d koeff = %d\n",
      //         i, st_item, fir_summing, koeff);

    end // block: linear_combination
  endfunction // fir_summing
endmodule // FIR_FILTER
//---------------------------------------------------------------------------
Полученный результат:
1) Сигнал 0.5 МГц.
Снимок экрана от 2023-03-17 16-49-04.png
2) Сигнал 2.1 МГц.
Снимок экрана от 2023-03-17 16-48-29.png
Как видно ослабление сигнала не соответствует полученному АЧХ.
Вопрос. Почему? Алиасинг? Или фильтр сделан/рассчитан не так?




И ещё вопрос. Возможно ли получить такую АЧХ используя одноразрядные коэффициенты.