Ответ IIR в Python

Вы не увидите ничего симпатичного с линейной шкалой x. Я не знаю numpy, но я знаком с Matlab, и есть некоторые функции, которые можно делать в журналах. Попробуйте использовать масштаб x-log с:

 import matplotlib.pyplot as pyplot fig = pyplot.figure() ax = fig.add_subplot(2,1,1) line, = ax.plot(w/np.max(w), h_dB, color='blue', lw=2) ax.set_xscale('log') show() 

Я не тестировал его, я не установил python 🙁

Редактировать:

Я попытался смоделировать фильтр butterworth в matlab для порядка 4 фильтра IIR и одного порядка 20 фильтра IIR.

 %!/usr/local/bin/matlab %% Inputs fs = 40e6; fc = 1e6; BW = 20e3; fl = (fc - BW/2); fh = (fc + BW/2); %% Build bandpass filter IIR Butterworth order 4 N = 4; % Filter Order h = fdesign.bandpass('N,F3dB1,F3dB2', N, fl, fh, fs); Hd1 = design(h, 'butter'); %% Build bandpass filter IIR Butterworth order 50 N = 20; % Filter Order h = fdesign.bandpass('N,F3dB1,F3dB2', N, fl, fh, fs); Hd2 = design(h, 'butter'); %% Compare fvtool(Hd1,Hd2); 

И здесь коэффициенты A и B для первого фильтра:

 FilterStructure: 'Direct-Form II Transposed' A: [2.46193004641106e-06 0 -4.92386009282212e-06 0 2.46193004641106e-06] B: [1 -3.94637005453608 5.88902106889851 -3.93761314372475 0.995566972065978] 

Если я получу какое-то время, я попытаюсь сделать то же самое с numpy!

Я не думаю, что проблема заключается в том, как вы анализируете ответ – это в вашем выборе фильтра. Вы пытаетесь создать очень узкий ответ фильтра, используя только фильтр низкого порядка IIR. Я думаю, вам нужен либо фильтр более высокого порядка, либо ослаблять ваши ограничения.

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

  b, a = signal.butter(4, [1.0/4-1.0/2e2,1.0/4+1.0/2e2], 'bandpass', analog=False) w, h = signal.freqs(b, a) import matplotlib.pyplot as plt fig = plt.figure() plt.title('Digital filter frequency response') ax1 = fig.add_subplot(111) plt.semilogy(w, np.abs(h), 'b') plt.ylabel('Amplitude (dB)', color='b') plt.xlabel('Frequency (rad/sample)') ax2 = ax1.twinx() angles = np.unwrap(np.angle(h)) plt.plot(w, angles, 'g') plt.ylabel('Angle (radians)', color='g') plt.grid() plt.axis('tight') plt.show() 

который дает:

Проблема в том, что signal.freqz возвращает точки на полукруге … поэтому вы не можете расширять до большего диапазона x , если только вы не делаете это через signal.freqz . Я попытался немного подтолкнуть его, и я увидел, что вы можете использовать pass whole=True для signal.freqz , и вы получите то, что у вас есть выше, но зеркально отразились на отрицательном x . Вот и все. Тем не менее, есть еще один аргумент ключевого слова, который позволяет вам передать массив x точек, которые вы хотите signal.freqz на … поэтому я попытался использовать np.arange(-5., 5., 0.1) … и он не выглядел на всех, как сюжет, который вы ожидаете справа, это выглядело как куча отражений исходного сюжета. Так что это заставляет меня думать … Может быть, у сюжета, который у вас есть справа, а у левой – разные топоры? В частности, это одна угловая частота, а другая простая старая частота?

При дальнейшем signal.freqz возвращает w,h , где w – нормализованная угловая частота в радианах / образце. Таким образом, вам не нужно выполнять нормализацию по np.max(w) в вашем коде, чтобы сделать сюжет. Тем не менее это все еще не решает проблему. Ваш график справа находится в единицах fc , а fc – в МГц (например, 1 / образец).

Поэтому, чтобы сделать график слева согласным справа, я предполагаю, что это означает, что вам нужно «ненормализовать» вашу ось x , тогда вам нужно будет преобразовать данные из единиц угловой частоты, используемой в МГц.

Или, может быть, более вероятно, используйте другую функцию, чем signal.freqz .