고역 통과 필터 설계
다른 글에서 사용했던 유사한 방법으로 Butterworth HPF를 설계하겠습니다. 이번엔 필터의 Q를 이용하지 않고, 필터의 차단 주파수와 감쇄량을 이용하는 방법을 알아보겠습니다.
다음과 같은 조건으로 고역 통과 필터를 만들겠습니다.
- 통과대역의 특성으로 차단 주파수와 통과대역 주파수$f_{c}=F_{p}=1000Hz$로 하겠습니다.
- 그러면 $F_{p}$ 에서 감쇄량은 -3dB가 됩니다.
- 그리고 저지 대역의 주파수 $F_{s}=80Hz$이며 감쇄량이 -80dB입니다.
- 이때 샘플링률 $f_{s}=44100Hz$, 샘플 구간 $T_{s}=1/f_{s}$입니다.
필터의 전달 함수 구하기
먼저 필터의 조건으로 다음과 같이 필터를 계산하기 위한 사전 조치를 취합니다.
$$A_{p}=10^{\left | \delta_{p}/10 \right |}$$
$$A_{s}=10^{\left | \delta_{s}/10 \right |}$$
$$\omega_{p}=2\pi F_{p}$$
$$\omega_{s}=2\pi F_{s}$$
$$\Omega_{p}= 2/T_{s} tan\left ( \omega_{p}T_{s}/2 \right )$$
$$\Omega_{s}= 2/T_{s} tan\left ( \omega_{s}T_{s}/2 \right )$$
그런 다음 조건에 부합하는 필터 차수 N을 다음식으로 구합니다. 여기서 ceil()은 올림 정수를 구하는 함수입니다.
$$N=ceil\left (\frac{1}{2} \frac{log((A_{s}-1)/(A_{p}-1))}{log(\Omega_{p}/\Omega_{s}))} \right )$$
그리고 N에 의해 $\omega_{c}=1$인 조건에 정규화된 버터워스 필터의 전달 함수를 선택합니다.
$$N==1\rightarrow H(s)=\frac{1}{s+1}$$
$$N==2\rightarrow H(s)=\frac{1}{s^{2}+\sqrt{2}s+1}$$
$$N==3\rightarrow H(s)=\frac{1}{(s+1)(s^{2}+s+1)}$$
$$N==4\rightarrow H(s)=\frac{1}{(s^{2}+0.76537s+1)(s^{2}+1.8477s+1)}$$
$H(s)$를 $s\rightarrow \Omega_{c}/s$로 치환해서 정규화된 필터를 HPF로 변경합니다.
이를 디지털 필터로 변환하기 위해서 또다시
$$s\rightarrow \frac{2}{T_{s}}\frac{z-1}{z+1}$$로 치환해서 $H(s)$를 $H(z)$로 변경합니다.
이 방법은 주워진 조건을 최대한 만족하는 필터의 전달 함수를 구해줍니다.
필터 계수 구하기
다음과 같이 필터 계수를 구하고 특성을 확인하는 GNU Octave의 스크립 파일을 만듭니다.
##
# butterworth HPF filter using BLT
##
clc;
clear all;
## FFT
NFFT = 2^12;
## load package
pkg load signal;
pkg load symbolic;
## symbol
z = sym('z','complex');
s = sym('s','complex');
## sampling rate
fs = 44100.0; #Hz
# sample time
Ts = 1.0/fs; #second
#frequency
fp = 1000.0; #Hz
fs = 80.0; #Hz
wp = 2.0*pi*(fp) #rad/s
ws = 2.0*pi*(fs) #rad/s
wc = wp
#prewaping frequency
Wp = 2.0/Ts * tan(wp*Ts/2.0);
Ws = 2.0/Ts * tan(ws*Ts/2.0);
Wc = Wp;
#attenuation
delta_p = -3.0; #dB
delta_s = -80.0; #dB
Ap = 10.0^abs(delta_p/10.0);
As = 10.0^abs(delta_s/10.0);
N_num = log10((As-1)/(Ap-1));
N_den = log10(Wp/Ws);
N = 0.5*N_num/N_den;
N=ceil(N)
#Hs_den = 1;
if(N==1)
Hs_den = s+1;
elseif(N==2)
Hs_den = s^2+ sqrt(2)*s + 1;
elseif(N==3)
Hs_den = (s+1)*(s^2+s+1);
elseif(N==4)
Hs_den = (s^2+0.76536*s+1)*(s^2+1.8477*s+1);
else
return
endif
Hs=1/Hs_den;
#subsitution
Hs = subs(Hs,s,(Wc/s));
Hs = simplify(expand(Hs));
#BLT
zb = 2.0/Ts * (z-1) /(z+1);
Hz = subs(Hs,s,(zb));
Hz = simplify(expand(Hz));
#filter coeffieient
[num_z,den_z]= numden(Hz);
coeff_b = coeffs(num_z,z,'all');
coeff_a = coeffs(den_z,z,'all');
coeff_num = double(coeff_b);
coeff_dev = double(coeff_a);
##Coefficient normalization
a0 = coeff_dev(1);
coeff_num /= a0
coeff_dev /= a0
#display
figure(1);
clf;
zplane(coeff_num,coeff_dev);
figure(2);
clf;
subplot(2,1,1);
[h,w] = freqz(coeff_num,coeff_dev,NFFT,fs);
m = 20.0*log10(abs(h));
semilogx(w,m);
#xlim([20 fs/2]);
ylim([-60 10])
xlabel("frequency(Hz)")
ylabel("magitude(dB)")
grid();
subplot(2,1,2);
semilogx(w,arg(h)*180/pi);
#xlim([20 fs/2]);
ylim([-200 200])
xlabel("frequency(Hz)");
ylabel("phase(deg)");
grid();
필터 특성 검토
z-영역에서 필터의 특성은 다음과 같습니다. 모든 극점이 단위원 안에 있어 안정적입니다.
주파수 특성은 다음과 같습니다. 고역 통과 필터가 형성되어 있습니다.
광고좀 꾹 눌러주시면 고맙겠습니다.
위의 내용을 참조용으로만 사용해주세요. 무단 도용이나 무단 복제는 불허합니다.
기타 문의 사항은 gigasound@naver.com에 남겨 주시면 고맙겠습니다.
'Audio Processing' 카테고리의 다른 글
오디오 필터(16)-IIR 필터 실행, 오디오 볼륨 처리 (0) | 2021.10.29 |
---|---|
오디오 필터 실험(15)-원시 데이터 얻기 (0) | 2021.10.29 |
오디오 필터(13)-IIR 버터워스 저역 통과 필터(Butterworth Low Pass Filter) 설계 (0) | 2021.10.29 |
오디오 필터(12)-IIR 노치 필터(Notch Filter, NF) (0) | 2021.10.29 |
오디오 필터(11)-IIR 전대역 통과필터(All Pass Filter, APF) (0) | 2021.10.29 |