I. 서 론

악기음과 음성 중의 배음 신호에서는 각 배음의 주파수가 기본음의 주파수, 즉 기본주파수의 정수배가 된다. 따라서, 배음 신호에서 기본주파수를 정밀하게 분석하는 것은 배음 신호의 특성을 분석하기 위해 중요한 과정이다. 배음 신호에서 기본주파수를 정밀하게 분석하는 방법은 자기상관함수(Autocorrelation function)에 기반하거나,^[1] 혹은 Discrete Fourier Transform(DFT)에 기반하는 방법 등이 있다.

DFT에 기반해서 배음의 기본주파수를 분석하는 경우 일반적으로 2단계의 과정을 거친다. 1단계는 coarse frequency estimation으로서, 신호의 DFT 크기 데이터에 대해 배음 구조의 특징을 이용해서 기본주파수에 해당되는 DFT 인덱스 값을 찾는다.^[2] 이 값으로부터 산출한 주파수는 실제의 기본주파수와 같지 않기 때문에, 2단계인 fine frequency estimation에서 주파수 추정치의 정확도를 높이기 위한 추가적인 방법이 적용되어야 한다. 이의 하나의 방법은 Phase Difference Estimation(PDE)과 같이 DFT의 위상 정보를 이용해서 주파수를 정밀하게 추정하는 방법이다.^[3,4] 또 다른 방법은 DFT의 크기 정보를 이용하는 방법인데, 예를 들어 신호에 Zero-padding^[5]을 적용한 후 DFT를 계산하고 이의 피크 위치를 찾는 방법이 있다. 이 고전적인 방법은 Zero-padding을 하지 않은 경우보다 상대적으로 더 정확하게 주파수를 추정할 수 있다.

Zero-padding 기반의 방법은 신호의 끝에 추가되는 0 데이터의 양이 많을수록 주파수 추정 에러가 줄어드는 경향이 있는데, 이는 정밀한 주파수 추정을 위해서 0 데이터가 매우 많이 추가되어야 함을 의미한다. 이 경우 DFT 계산시에 연산량과 메모리 사용량이 과도하게 늘어나게 되는 문제가 있다. 이의 대안으로서 DFT 크기 데이터들을 보간하여 실제 피크의 위치를 추정하는 방식으로 주파수를 추정하는 DFT 보간 기법들이 있다.^[3,6,7,8] 이들은 연산량이 Zero-padding에 비해 상대적으로 적으면서도 성능은 우수하다는 장점이 있다.

주파수 분석 기법들은 일반적으로 신호의 특성이나 잡음 등에 따라서 성능이 달라질 수 있다. 본 논문에서는 4가지 DFT 보간 기법들에 대해 기본주파수, 신호대잡음비(Signal to Noise Ratio, SNR), 배음의 수, DFT 사이즈 등에 따른 배음 신호의 기본주파수 추정 성능을 분석 및 비교한다. DFT 보간 기법들에 대한 비교 연구는 잡음이 없고 배음의 수가 고정된 비브라토음 모델에 대해 수행된 바 있으나,^[9] 정현파 기반의 배음 신호에 대해서 기본주파수, 배음의 수, 신호대잡음비, DFT 사이즈 등 다양한 측면에서 DFT 보간 기법들의 성능을 분석 및 비교한 연구는 많지 않았다. 따라서, 본 논문의 결과는 배음 신호의 기본주파수를 분석함에 있어 신호의 특성 및 잡음의 정도 등에 따른 최적의 DFT 보간 기법을 선택하는 것에 도움이 될 수 있을 것이다.

II. DFT 보간 기법들

0≤n≤N-1에 대해 x[n]의 Discrete-Time Fourier Transform(DTFT)을 다음과 같이 정의하자.

여기에서 w[n]는 윈도우 함수이며, 윈도우의 길이는 N이고, 𝜔는 bin 단위로 정규화된 주파수이다. 이제 x[n]=a×cos2πf0fsn+θ라고 하면, 이의 DTFT는 다음과 같다. 단, fs는 샘플링 주파수이다.

여기에서 W(ω)는 w[n]의 DTFT이고, ω0=N×f0/fs이다. Eq. (2)에서 |X(ω)|의 피크는 ω=±ω0에서 발생하며, 따라서 Fig. 1에서처럼 ω0 근처에서 |X(ω)|의 피크 위치를 찾음으로써 주파수 f0를 구할 수 있다.

Fig. 1.

(Color available online) |X(ω)| and |X[k]| around ω=ω0.

DTFT는 실제로는 DFT로 구현되어야 한다. x[n]의 DFT를 다음과 같이 정의하자.

이제 |X[k]|의 피크가 Fig. 1에서처럼 k=m bin에서 발생한다고 하자. ω0와 m은 일반적으로 일치하지 않기 때문에, m으로부터 주파수를 산출하면 에러를 포함한다. DFT 보간 기법에서는 k=m bin 주위의 DFT 크기 샘플들을 보간해서 ω0와 m 사이의 차이 𝛿를 계산하고, 주파수 추정치를 다음처럼 산출한다.

본 논문에서는 𝛿를 계산하기 위해서 4가지의 DFT 보간 기법을 각각 적용한다. 우선, 2차 피팅(Quadratic Fit, QF)은 고전적인 방법으로서, 가장 큰 값인 |X[m]|과 그 좌우의 값들 |X[m-1]|, |X[m+1]|이 포물선을 이룬다고 가정해 다음과 같이 𝛿를 구한다.^[3] QF은 |W(ω)|에 대한 모델링에 기반하지 않으므로, 임의의 윈도우 함수를 사용했을 경우에도 적용 가능하다.

두 번째로, Interpolated Fast Fourier Transform(IFFT)에서는 Rectangular, Hanning, Rife-Vincent 윈도우 함수를 사용한 경우들에 대해서 𝛿를 구한다.^[6] Hanning 윈도우 함수를 사용한 경우, ω0 주위에서 가장 큰 2개의 피크값을 |X[m]|과 |X[m+1]|라고 할 때 다음과 같이 𝛿를 구한다.

세 번째로, Weighted Multipoint Interpolated Discrete Fourier Transform(WMIDFT)에서는 가장 큰 값 |X[m]| 및 그 좌우의 값들인 |X[m-1]|과 |X[m+1]|을 이용해 𝛿를 구한다.^[7] WMIDFT는 |W(ω)|의 모델링에 기반하며, Rectangular, Hanning 윈도우 함수를 사용하는 경우에 대해서 적용이 가능하다. Hanning 윈도우에 대한 𝛿 계산은 다음과 같다.

마지막으로 Corrected Quadratically Interpolated Fast Fourier Transform(CQIFFT)^[8]에서는 신호에 zero-padding을 적절히 적용하여 DFT를 수행하고, DFT 크기에 로그를 취한 값들 중 가장 큰 값과 그 좌우 값들을 이용해 QF 기반으로 δ^를 구한다. 또한, δ^를 다음과 같이 보정해서 최종적으로 𝛿를 계산한다.

단, ξZp는 Zero-padding factor Zp와 윈도우 함수의 종류에 의해 계산되는 상수이며,^[8]Zp는 zero-padding이 적용된 전체 신호의 길이를 N으로 나눈 값이다. CQIFFT는 Rectangular, Hamming, Hanning, Blackman, Kaiser-Bessel 윈도우에 대해 적용할 수 있다.

본 논문에서는 네 가지의 DFT 보간 기법 적용을 위해 DFT 계산 시 Hanning 윈도우 함수를 사용한다. Hanning 윈도우는 적절한 메인로브 대역폭 및 사이드로브 감쇠 특성으로 인해서 DFT 기반 주파수 추정에서 자주 사용되고 있다.^[10,11] 또한, 본 논문에서 비교하는 네 가지의 DFT 보간 기법들에 대해 모두 적용이 가능하다.

III. 주파수 추정 에러의 특성 분석

이제 배음 신호가 0≤n≤N-1에 대해 다음과 같이 L개의 배음으로 구성되어 있으며, 잡음도 존재한다고 가정하자.

여기에서 al,fl,θl는 l번째 배음의 크기, 주파수, 위상이며, z[n]은 평균이 0인 백색 잡음이다. f0를 기본주파수라고 하면, 일반적으로 fl=l×f0라고 할 수 있다. 이 신호의 DTFT를 Eq. (1)을 이용해 구하고 정리하면 다음과 같다.

Eq. (10)에서 첫째 항은 ω=ω0에서 피크가 발생하는 기본주파수 성분에 해당되며, 이를 타겟 성분인 Xt(ω)라고 하자. Xi(ω)는 배음들 중 타겟 성분을 제외한 일종의 간섭 성분을 의미하며 Eq. (11)과 같다.

Z(ω)는 z[n]×w[n]의 DTFT이며, 잡음에 의한 성분이다.

DFT 보간 기법을 이용해 f0를 구할 때, Xi(ω)와 Z(ω)는 추정 에러를 발생시킨다. 신호대잡음비를 다음과 같이 기본음과 잡음의 파워 비율로 정의한다.

여기에서 𝜎는 잡음의 표준편차이다. 이제 기본주파수 추정 시에 신호대잡음비에 따른 Xi(ω)와 Z(ω)의 영향을 분석하고, 이에 따른 주파수 추정치의 에러 경향을 정성적으로 분석한다.

우선 신호대잡음비가 낮은 경우, 잡음의 파워가 상대적으로 크므로 |Z(ω)|≫Xi(ω)라고 근사하자. Fig. 2(a)와 같이 X(ω)=Xt(ω)+Z(ω)로 가정할 경우, 이를 기반으로 주파수 추정 에러의 경향을 여러 측면에서 정성적으로 분석할 수 있다. 우선, DFT 보간 기법들이 이용하는 |X[m-1]|, |X[m]|, |X[m+1]|은 잡음 성분에 의한 왜곡을 포함하며, 이에 따른 주파수 추정 에러가 발생한다. 신호대잡음비와 잡음 성분의 파워는 반비례하므로, 신호대잡음비가 커질수록 주파수 추정 에러는 작아지는 형태가 될 것이다. 또한, 기본주파수 f0가 증가함에 따라서 -0.5bin<m-ω0≤0.5bin이면서 m-ω0의 값이 1 bin 단위로 반복되므로, 주파수 추정 에러도 𝜔축에서 1 bin 주기로 반복하게 될 것이다. 배음의 수 L은 Xt(ω) 및 Z(ω)와 상관이 없으므로, 에러에 큰 영향이 없을 것이다.

Fig. 2.

(Color available online) Approximation of |X(ω)| (a) at low SNR and (b) at high SNR.

신호대잡음비가 높은 경우, 잡음의 파워가 상대적으로 작으므로 Xi(ω)≫|Z(ω)|라고 근사하자. Fig. 2(b)와 같이 X(ω)=Xt(ω)+Xi(ω)로 가정할 경우, |X[m-1]|, |X[m]|, |X[m+1]|은 배음들의 간섭 성분에 의한 왜곡을 포함하며, 이에 따른 주파수 추정 에러가 발생한다. 윈도우 함수가 고정된 경우, 우선 DTFT에서 연속된 배음들의 피크는 서로 ω0(=N×f0/fs) bin만큼 떨어져 있으므로, f0가 커질수록 Xi(ω)에 의한 간섭은 약해지고 주파수 추정 에러도 작아지는 경향을 보일 것이다. 또한, 기본음의 크기 a1이 클수록, 배음들의 크기 al(l≠1)가 작을수록 Xi(ω)의 영향이 감소하여 에러가 작아질 것이다. 배음의 수 L이 클수록 Xi(ω)의 영향이 커지므로 에러가 커지는 경향이 있을 것이다. 한편, 신호대잡음비는 Xi(ω)와 상관이 없으므로, 에러의 크기에 큰 영향이 없을 것이다.

IV. 실험 결과

DFT 보간 기법들의 성능 비교를 위해서 모의 실험 및 실제 악기음의 분석을 시행하였다. 이 중 모의 실험에서는 기본주파수 f0, 신호대잡음비, 배음의 수 L을 변수로 하고, DFT 사이즈를 2가지 값으로 사용하여 각 경우에 대해 주파수 추정 에러를 산출하였다.

4.1 f0에 따른 기본주파수 추정 성능 분석

우선 첫 번째 실험에서는 f0의 값에 따른 기본주파수 추정 에러를 분석하였다. 실험을 위한 신호는 Eq. (9)의 x[n]을 이용했으며, L=5, 1≤l≤5에 대해 al=1, θl은 [0,2π) 구간에서 랜덤한 값을 가지도록 해서 동일한 크기의 5개 배음으로 구성된 배음 신호를 생성했다. fs=48 kHz, z[n]은 평균이 0이고 표준편차가 𝜎인 백색 가우시안 잡음을 이용했으며, 신호대잡음비에 따라서 𝜎를 변화시켰다. DFT 사이즈 N 은 1024와 2048을 사용했는데, 이는 각각 21.3 ms와 42.7 ms의 시간 길이에 해당된다. ω0는 4 bin부터 0.01 bin 간격으로 증가하여 전체 1,000개의 값을 가지도록 하였는데, 결과적으로 ω0는 (4 bin, 14 bin) 범위의 값을 가진다. 이는 N=1024인 경우 f0를 187.5 Hz부터 0.46875 Hz 간격으로 증가시켜서 (187.5 Hz, 656.25 Hz) 범위가 되도록 하고, N=2048인 경우 f0를 93.75 Hz부터 0.234375 Hz 간격으로 증가시켜서 (93.75 Hz, 328.125 Hz) 범위가 되도록 한 것을 의미한다. ω0의 최소 범위를 4 bin으로 한 이유는, Hanning 윈도우를 사용했을 때 X(ω)에서 기본음과 2배음의 메인로브들이 𝜔축에서 서로 겹치지 않도록 하기 위해서이다.

본 논문에서는 배음 신호의 정밀한 기본주파수 추정을 목적으로 한다. 따라서, Coarse frequency estimation은 기존의 방법^[2]을 이용해 성공적으로 수행할 수 있다고 가정하였고, 기본주파수의 대략적인 범위를 미리 정한 후 범위 내의 가장 큰 피크의 DFT 인덱스 값을 찾는 방식으로 구현하였다. Fine frequency estimation을 위해서 QF, IFFT, WMIDFT, CQIFFT를 적용하였고, DFT 계산 시 모두 Hanning 윈도우를 이용하였다. CQIFFT의 경우 Zp=2를 사용했다.

기본주파수의 추정 에러는 Root Mean Square Error (RMSE)로 Ef0-f^02와 같이 계산했다. 식에서 E{∙}는 기댓값이고, f0는 진짜 기본주파수, f^0는 추정 기본주파수이다. Figs. 3과 4에서 결과를 보이고 있는데, 그림에서 각각의 결과는 1,000번의 독립 시행의 평균 결과이다.

Fig. 3.

(Color available online) RMSEs of each DFT interpolation technique as a function of f0 at 0 dB SNR. (a) N = 1024. (b) N = 2048. Note that the y-axis is displayed on a log scale.

Fig. 4.

(Color available online) RMSEs of each DFT interpolation technique as a function of f0 at 60 dB SNR. (a) N = 1024. (b) N = 2048.

신호대잡음비가 낮을 때의 결과인 Fig. 3에서는 CQIFFT가 가장 작은 주파수 추정 에러를 보이고, IFFT, WMIDFT, QF 순서로 에러가 점점 증가함을 보인다. IFFT, WMIDFT, QF의 경우 3장에서 분석한 바와 같이 에러가 f0에 대해 1 bin마다 주기적인 형태를 보인다. 다만, CQIFFT는 f0에 상관없이 에러의 크기가 일정한 형태를 보인다. DFT 사이즈인 N이 커지면 에러도 작아지지만, 에러의 형태는 (a)와 (b)의 결과가 대체로 유사하다.

신호대잡음비가 높을 때의 결과인 Fig. 4에서는 f0에 따라서 각 기법의 에러 정도가 달라진다. f0가 낮을 때는 WMIDFT가 가장 에러가 작지만, f0가 높을 때는 CQIFFT가 가장 에러가 작은 경향을 보인다. 즉, f0의 값에 따라서 최적의 방법이 달라질 수 있음을 의미한다. 또한, WMIDFT, CQIFFT, IFFT는 f0가 증가함에 따라서 3장에서 분석한 바와 같이 에러가 줄어들거나, 혹은 줄어들다가 일정한 크기의 형태로 유지되는 것이 관측된다. 다만, QF의 경우는 에러가 가장 크며 f0에 따른 크기의 형태 변화는 별로 없다.

4.2 신호대잡음비에 따른 기본주파수 추정 성능 분석

두 번째 실험에서는 신호대잡음비의 크기에 따른 기본주파수 추정 에러를 분석하였다. 신호와 실험 조건은 첫 번째 실험과 동일하게 하였다. 특정 신호대잡음비의 값 별로 f0의 1,000개의 값에 대해 각각 1,000번의 독립시행 결과를 구해 전체를 평균하여 RMSE를 구했는데, 이는 RMSE를 1,000,000번의 평균에 의한 결과로 산출한 것에 해당된다.

Fig. 5에서 신호대잡음비가 낮을 때는 3장에서 분석한 바와 같이 RMSE의 크기가 신호대잡음비에 반비례한다. 에러는 CQIFFT가 가장 작으며, IFFT, WMIDFT, QF의 순서로 에러가 커지는데, 이는 Fig. 3의 결과와 같은 경향을 보인다. 신호대잡음비가 높을 때는 3장의 분석과 같이 신호대잡음비의 값과 상관없이 에러의 크기가 일정한 경향을 보인다. WMIDFT의 에러가 가장 작고, CQIFFT, IFFT, QF의 순서로 에러가 커진다. Fig. 4의 결과에서는 기본주파수가 낮을 때는 WMIDFT가 가장 에러가 작았고, 기본주파수가 높을 때는 CQIFFT가 가장 에러가 작았는데, 기본주파수에 대해 전체 평균을 구했을 때는 WMIDFT가 CQIFFT보다 더 에러가 작음을 의미한다. 이상의 결과에 의하면 신호대잡음비 및 기본주파수의 값에 따라서 최적의 DFT 보간 기법이 달라질 수 있음을 알 수 있다.

Fig. 5.

(Color available online) RMSEs of each DFT interpolation technique as a function of SNR. (a) N = 1024. (b) N = 2048.

한편, 신호대잡음비가 높을 때는 Fig. 4와 같이 f0의 값에 따라서 보간 기법들의 성능이 달라질 수 있으므로, 실험에서 f0의 범위를 어떻게 설정하느냐에 따라서 전체 평균 에러의 결과도 달라질 가능성이 있다. 추가적으로, Fig. 5에서 신호대잡음비가 높을 때는 WMIDFT, CQIFFT, IFFT의 에러가 모두 작다는 것은 언급할만하다. 예를 들어, 신호대잡음비 = 55 dB일 때 WMIDFT, CQIFFT, IFFT의 RMSE는 N = 1024일 때 0.040 Hz, 0.071 Hz, 0.084 Hz이고, N = 2048일 때 0.019 Hz, 0.036 Hz, 0.041 Hz이다.

4.3 배음의 수에 따른 기본주파수 추정 성능 분석

세 번째 실험에서는 배음의 수에 따른 기본주파수 추정 에러를 분석했다. Eq. (9)에서 L을 1부터 5까지 변화시켰으며, 다른 실험 조건은 두 번째 실험 조건과 같다. RMSE 결과는 1,000,000번의 평균 결과이다.

Fig. 6에서 신호대잡음비 = 0 dB인 경우 RMSE는 모든 보간 기법들에 대해 L과 상관없이 일정한 크기의 형태를 보이는데, 이는 3장에서 분석한 바와 같이 신호대잡음비가 낮을 때 배음의 영향은 크지 않음을 의미한다. 신호대잡음비 = 60 dB인 경우 L이 1에서 2로 증가할 때 3장에서 분석한 바와 같이 에러가 증가하고, 이후로는L이 증가해도 에러는 일정한 크기의 형태를 보인다. 이는 기본음만 있을 때는 에러가 상대적으로 작지만, 2배음이 추가되면 이의 사이드로브가 기본주파수 추정 시 영향을 주면서 에러가 커짐을 의미한다. 3배음 ~ 5배음이 추가되면 이들의 사이드로브 영향도 추가되지만, 이들 배음은 주파수 축에서 기본음으로부터 2배음에 비해 더 멀리 떨어져 있다. 따라서, 본 실험에서 사용한 신호 모델에서 3배음 ~ 5배음의 사이드로브에 의한 영향은 2배음의 사이드로브에 의한 영향보다 매우 작았음을 의미한다. 신호대잡음비가 높을 때 L=1이면 CQIFFT가 WMIDFT보다 더 에러가 작았는데, 이는 단일 정현파에 대해서는 전자가 후자보다 성능이 더 좋았음을 의미한다.

Fig. 6.

(Color available online) RMSEs of each DFT interpolation technique as a function of the number of harmonics. (a) N = 1024. (b) N = 2048.

4.4 악기음에 적용한 결과 예

DFT 보간법을 University of Iowa musical instrument samples^[12]에서 제공하는 ‘Flute.nonvib.ff.A4.stereo.aif’ 파일에 적용해서 기본주파수의 값들을 분석하였다. 이는 무향실에서 녹음한 A4음의 플루트 연주음으로서, 24 bit로 되어 있고 44.1 kHz의 샘플링주파수를 가진다. 스테레오 신호로 되어있지만 좌우의 데이터가 같아서 실제로는 모노 신호라고 할 수 있다. 이의 파형은 Fig. 7(a)에 보이는데, 파형이 깨끗하고 소리를 들었을 때 잡음이 거의 없어서 신호대잡음비가 높은 것으로 판단된다.

Fig. 7.

(Color available online) Results of each DFT interpolation technique for the flute sound. (a) Waveform. (b) Fundamental frequency trajectories.

이 신호에 대해 프레임의 길이를 2048 샘플로 하고, 연속되는 두 프레임 사이의 시간 간격을 1024 샘플로 해서 전체 프레임들을 구성하였다. 각 프레임별로 Hanning 윈도우 함수를 적용해서 DFT를 계산한 후에, DFT 보간 기법들을 적용해서 구한 기본주파수 결과들을 Fig. 7(b)에서 보이고 있다. 그림에서 IFFT, WMIDFT, CQIFFT의 결과는 대체로 유사하며, QF의 결과는 다소 다름을 관측할 수 있다. DFT 보간 기법을 적용하지 않았을 경우, 즉 Eq. (4)에서 𝛿=0인 경우, 기본주파수는 전반적으로 452.1973 Hz로 측정되었으며, 이는 Fig. 7(b)의 결과들과 차이가 있다.

4.5 논의 사항

이상의 실험 결과에 근거해 기본주파수 분석을 위한 DFT 보간 기법에 대해 논의한다. 우선 QF, IFFT, WMIDFT는 N-point DFT 수행 후 Eqs. (5), (6), (7)을 적용했고, CQIFFT는 Zp=2로 Zero-padding을 해서 2N-point DFT를 수행한 후 로그값을 이용해 QF 및 Eq. (8)을 적용했다. 각 방법의 전체 연산량에서 DFT 연산량이 가장 주요하며, 따라서 CQIFFT가 다른 기법들에 비해서 연산량이 더 많은 편이다. 추가적인 연산량을 감수할 수 있다면, 신호대잡음비가 낮은 경우 혹은 신호대잡음비가 높고 기본주파수가 상대적으로 높거나 단일 정현파로 구성되는 등 배음의 영향이 적은 경우 CQIFFT가 좋은 성능을 보였다. 신호대잡음비가 높고 기본주파수가 상대적으로 낮을 때는 WMIDFT가 좋은 성능을 보였다. 추가적인 연산량을 감수할 수 없고 신호대잡음비가 낮으면, IFFT가 WMIDFT보다 성능이 더 좋았다.

이상의 논의 사항은 본 논문의 신호 및 실험 조건에 대해서 유효하다는 것을 밝힌다. 이와 관련하여, 추후 다양한 신호의 특성 및 실험 조건에 대한 추가적인 연구가 필요하다.

V. 결 론

배음 신호의 정밀한 기본주파수 분석을 위해서 4가지의 DFT 보간 기법들에 대해 주파수 추정 에러의 특성을 정성적으로 분석하고, 모의 실험 및 녹음된 플루트음을 이용해 각 기법의 성능을 비교하였다. 또한, 실험 결과에 기반해서 기본주파수, 신호대잡음비, 배음의 수 등에 따른 최적의 기법을 논의하였다. 본 논문의 결과가 배음 신호의 기본주파수 분석을 위해 DFT 보간 기법을 선택하는 데 도움이 되기를 기대한다.