sergeyhaki
commented
1 year ago Для определения сдвига фазы all-pass фильтра на частотах в диапазоне 200-3500 Гц можно вычислить фазовую характеристику фильтра и обработать ее взвешенным интегрированием. Фазовая характеристика фильтра представляет собой зависимость фазового угла от частоты. На основе этой характеристики можно вычислить задержку фильтра и пересчитать ее в сдвиг фазы на заданном диапазоне частот.
Для начала необходимо вычислить фазовую характеристику H(ω) фильтра на всем диапазоне частот. Для этого можно использовать функцию atan2, которая вычисляет арктангенс по двум аргументам, чтобы получить фазовый угол на каждой частоте. Затем нужно вычислить задержку фильтра на каждой частоте. Задержка фильтра может быть вычислена как производная фазовой характеристики по частоте. Для реализации этого вычисления может использоваться численное дифференцирование. Наконец, полученные задержки фильтра можно пересчитать в сдвиг фазы на диапазоне частот от 200 до 3500 Гц.
Пример кода на языке Python для определения сдвига фазы all-pass фильтра на диапазоне частот от 200 до 3500 Гц:
import numpy as np
def allpass_filter(x, b, a, xn=None, yn=None):
if xn is None:
xn = np.zeros(len(b)-1)
if yn is None:
yn = np.zeros(len(a)-1)
y = np.zeros(len(x))
for n in range(len(x)):
y_ = b[0]*x[n] + b[1]*xn[0] + b[2]*xn[1] - a[1]*yn[0] - a[2]*yn[1]
xn[1], xn[0] = xn[0], x[n]
yn[1], yn[0] = yn[0], y_
y[n] = y_
return y, xn, yn
def phase_delay(w, h):
delta_w = np.diff(w)
delta_h = np.diff(h)
delta_phi = np.arctan2(delta_h, delta_w)
delta_phi = np.concatenate(([delta_phi[0]], delta_phi))
tau = -np.cumsum(delta_phi)/w
return tau
# Коэффициенты фильтра
b = np.array([1.0, -1.9315, 0.9315])
a = np.array([0.9647, -1.9315, 0.9647])
# Диапазон частот
w = np.logspace(np.log10(200), np.log10(3500), 1000)
# Вычисляем фазовую характеристику фильтра
h = np.angle(np.polyval(b, np.exp(1j*w)) / np.polyval(a, np.exp(1j*w)))
# Вычисляем задержку фильтра
tau = phase_delay(w, h)
# Вычисляем сдвиг фазы на диапазоне частот от 200 до 3500 Гц
f_range = np.where((w>=200) & (w<=3500))[0]
phi = -np.angle(np.exp(-1j*np.cumsum(tau[f_range])*w[f_range]))
phi_deg = np.rad2deg(phi)
# Выводим результат
print(f"Сдвиг фазы на диапазоне частот от 200 до 3500 Гц:")
print(f"Минимальный сдвиг фазы: {np.min(phi_deg):.2f} градусов")
print(f"Максимальный сдвиг фазы: {np.max(phi_deg):.2f} градусов")В результате выполнения этого кода будет определен сдвиг фазы all-pass фильтра на диапазоне частот от 200 до 3500 Гц.
void freq(int32_t fout, uint16_t i, uint16_t q){ // Set a CLK0,1,2 to fout Hz with phase i, q (on PLLA) uint8_t rdiv = 0; // CLK pin sees fout/(2^rdiv) if(fout > 300000000){ i/=3; q/=3; fout/=3; } // for higher freqs, use 3rd harmonic if(fout < 500000){ rdiv = 7; fout = 128; } // Divide by 128 for fout 4..500kHz uint16_t d; if(fout < 30000000) d = (16 fxtal) / fout; else d = (32 fxtal) / fout; // Integer part .. maybe 44? if(fout < 3500000) d = (7 fxtal) / fout; // PLL at 189MHz to cover 160m (freq>1.48MHz) when using 27MHz crystal if(fout > 140000000) d = 4; // for f=140..300MHz; AN619; 4.1.3, this implies integer mode if(d % 2) d++; // even numbers preferred for divider (AN619 p.4 and p.6) if( (d (fout - 5000) / fxtal) != (d (fout + 5000) / fxtal) ) d += 2; // Test if multiplier remains same for freq deviation +/- 5kHz, if not use different divider to make same uint32_t fvcoa = d fout; // Variable PLLA VCO frequency at integer multiple of fout at around 27MHz16 = 432MHz // si5351 spectral purity considerations: https://groups.io/g/QRPLabs/message/42662
ifdef F_CLK2
else
endif
ifdef x
endif
}