Introduction to Digital Signal Processing.md

Introduction to Digital Signal Processing

1. 디지털 신호처리 기초

신호의 종류

신호는 시공간적으로 변하는 물리량으로, 형태에 따라 아날로그 신호와 디지털 신호로 나눌 수 있다.

신호를 목적에 맞게 가공하는 과정을 신호처리라고 한다.

하나의 소스에서 발생되는 신호를 스칼라 신호라고 하고, 시간에 따라 두 개 이상의 소스에서 발생하는 신호를 벡터 신호라고 한다. 스칼라 신호는 하나의 독립 변수로 표현할 수 있으므로 1차원 신호라고 하고, 벡터 신호는 두 개 이상의 독립 변수로 표현되므로 다차원 신호라고 한다. 예를 들어 영상(image)신호는 2차원 신호이며, 동영상(video) 신호는 3차원 신호이다.

디지털화

아날록 신호를 디지털 시스템에서 처리하기 위한 과정을 디지털화 과정이라고 한다. 디지털화는 기본적으로 샘플링과 양자화를 거쳐서 수행된다. 샘플링은 연속 시간 신호를 이산 시간 신호로 변환하는 과정을 말하고, 양자화는 신호의 크기를 연속된 크기에서 몇 개의 크기로 매핑하는 과정을 말한다.

일반적으로 연속 신호를 샘플링하여 얻은 신호를 이산 신호라고 하고, 이를 양자화해서 얻은 신호를 디지털 신호라고 한다.

나이퀴스트 샘플링 정리

샘플링을 촘촘하게 수행할수록 이산 신호는 연속 신호로 손실 없이 복원될 수 있다. 나이퀴스트가 증명한 샘플링 정리에 따르면 대역 제한된 신호(band-limited signal)에 대해, 몇 개의 이산 표본으로 원래의 아날로그 신호를 완벽하게 표현할 수 있다.

2. 시간 영역에서의 분석

신호

단위 임펄스 신호는 하나의 샘플 신호를 가지는 신호로, 다음과 같이 정의된다. $$ \delta[n] = \begin{cases} 1, &(n=0) \ 0, &(n \ne 0) \end{cases} $$ 단위 계단 신호는 신호가 계단 형태여서 붙여진 이름으로 다음과 같이 정의된다. $$ u[n] = \begin{cases} 1, &(n \ge 0) \ 0, &(n < 0) \end{cases} $$ 단위 계단 신호를 1 지연시킨 후 빼서 단위 임펄스 신호를 얻을 수 있다.

시스템

이산 신호 $x$는 샘플 인덱스 $n$에 대하여 $x[n]$와 같이 표현한다.

$x[n]$은 일종의 함수이기 때문에 상수로 착각하지 않도록 하자.

이산 신호를 입력으로 받을 때 이산 출력 신호를 생성하는 함수를 시스템이라고 한다.

선형 $\cdot$ 비선형 시스템

중첩 합 정리를 만족시키는 시스템을 선형 시스템, 그렇지 않은 시스템을 비선형 시스템이라고 한다.

중첩 합 정리는 합 성질과 스칼라 곱 성질을 모두 만족시키는 것으로, 시스템 $T$에 대하여 간단히 다음과 같이 나타낼 수 있다. $$ T{\alpha x[n] + \beta y[n]} = \alpha T{x[n]} + \beta T{ y[n] } $$

시불변 $\cdot$ 시가변 시스템

어떤 입력 신호 $x[n]$에 대하여 출력 신호 $y[n]$이 생성되는 시스템 $T$과 모든 정수 $k$에 대하여 $x[n-k]$를 입력 신호로 했을 때 출력 신호가 $y[n-k]$가 생성되는 시스템을 시불변 시스템이라고 하며 그렇지 않은 시스템을 시가변 시스템이라고 한다.

예제

다음 시스템이 시불변 시스템인지 시가변 시스템인지 밝히시오. $$ y[n]= \sum _{k= -\infty} ^{n} x[k] $$ 신호 $x_1[n]$을 $n_0$만큼 지연시켜서 얻은 데이터를 시스템에 통과시키면 $$ \sum _{k= - \infty } ^{n} x[k-n_0 ] = \sum _{ k = -\infty} ^{ n - n_0} x[k] = y[n-n_0] $$ 을 얻을 수 있다. 따라서 이 시스템은 시불변 시스템이다.

위에 주어진 시스템은 누산기 시스템이라 불린다.

임펄스 응답 시스템

임의의 시스템이 갖는 고유한 성질을 알고자 할 때 주로 사용하는 방법은 시스템에 임퍼스 입력 $\delta[n]$을 넣어 나오는 출력 $h[n]$의 형태로 파악하는 것이다. 이를 임펄스 응답이라고 한다. 임펄스 응답은 이산 시스템의 고유한 성질을 반영하게 되어 임펄스 응답 시스템이라고 부른다.

LTI 시스템

선형성과 시불변성을 동시에 만족하는 시스템을 LTI(Linear time-invariant) 시스템이라고 한다. LTI 시스템은 입력과 출력의 관계를 쉽게 해석할 수 있어서 매우 중요한 시스템이다. LTI 시스템에서 임펄스 응답은 아주 큰 역할을 한다. $$ x[n] = \sum _{k=-\infty} ^{\infty} x[k] \delta [n-k] $$

$$ y[n] = T{x[n]} = \sum _{k = -\infty} ^{\infty} x[n] T { \delta [n-k] };;\cdots;\textrm{linearity} $$

$$ y[n] = \sum _{ k = -\infty } ^{\infty} x[n]h[n-k] = x[n]*h[n] ;;\cdots;\textrm{time invariant} $$

위의 과정을 통해 LTI시스템의 임펄스 응답 $h[n]$을 가지고 있으면, 그 출력 $y[n]$을 곧바로 알 수 있다는 사실을 알 수 있다. 또한 LTI의 시스템의 이러한 특성 때문에, 여러 LTI 시스템을 직렬로 통과시키는 것은 여러 개의 임펄스 응답 $h_2[x]$, $h_2[n]$, $\cdots$을 연쇄적으로 컨벌루션을 수행하는 것으로 계산할 수 있다. 이는 결국 여러 개의 LTI 시스템의 임펄스 응답을 가지고 있으면 이를 단 한 번의 컨벌루션으로 줄일 수 있다는 것을 함의한다. 따라서 시스템을 설계할 때, 여러 LTI 시스템을 직렬로 연결해서 마음에 드는 시스템을 얻어내고 나서 그 최종 임펄스 응답만을 기록하면 된다. 이와 같은 방법으로 뮤직 플레이어에서 가끔 볼 수 있는 이퀄라이저가 만들어진다.

인과성(causality)과 안정성(Stability)

LTI 시스템의 임펄스 응답이 절대 수렴할 때, LTI 시스템은 안정성을 갖는다고 말한다. $$ S = \sum _{ k = -\infty} ^{\infty} |h[k]| < \infty $$ LTI 시스템이 과거와 현재에만 영향을 받고, 미래 신호에 대해서는 영향을 받지 않을 때, LTI 시스템은 인과성을 갖는다고 말한다. $$ y(t) = x(t) * h(t) = \int _{-\infty} ^{\infty} x(\tau)h(t-\tau)d\tau $$

$$ =\int _{-\infty} ^t x(\tau)h(t-\tau)d\tau + \int _t ^{\infty} x(\tau)h(t-\tau) d\tau $$

$$ 위의;적분식에서;오른쪽;항은;미래;시간에;대한;식이므로;인과성을;갖기;위해서는 $$

$$ h(t-\tau) = 0 ;for;t \le \tau < \infty $$

따라서, 인과성의 필요충분조건은 다음과 같다. $$ \therefore h(t)=0 ;;\textrm{if} ; t<0 $$

3. 신호의 변환(transform)

신호의 확장

신호를 확장한다는 것은 같은 신호를 다른 기저로 표현한다는 것이다. 어떤 공간의 모든 원소를 선형 독립인 몇 개의 원소의 선형 조합으로 표현 가능할 때, 그 몇개의 원소를 그 공간의 기저라고 한다.

정현$\cdot$여현파의 선형 독립

정현파 함수와 여현파 함수에는 다음과 같은 성질이 있다. $$ <\sin(mx), \sin(nx)> = \int {-\pi} ^{\pi} \sin(mx)\sin(nx)dx = \pi \delta{mn} $$

$$ <\cos(mx), \cos(nx)> = \int {-\pi} ^{\pi} \cos(mx)\cos(nx)dx = \pi \delta{mn} $$

$$ <\sin(mx), \cos(nx)> = \int _{-\pi} ^{\pi} \sin(mx)\cos(nx)dx = 0 $$

$$ <\sin(mx), 1> = \int _{-\pi} ^{\pi} \sin(mx)dx = 0 $$

$$ <\cos(nx), 1> = \int _{-\pi} ^{\pi} \cos(nx)dx = 0 $$

푸리에 계수

따라서 폐구간 $[-\pi, \pi]$ 위에서 직교(orthogonal)하는 정현파 및 여현파 함수들을 이용해 다음과 같이 원본 연속 신호 $f$를 확장할 수 있다. $$ f(x) = \frac{1}{2} a_0 + \sum _{n=1} ^{\infty} a_n \cos(nx) + b_n \sin(nx) $$ 각 계수는 다음과 같이 구할 수 있다. $$ a_0 = \frac{1}{\pi} \int _{-\pi} ^{\pi} f(x)dx $$

$$ a_n = \frac{1}{\pi} \int _{-\pi} ^{\pi} f(x)\cos(nx)dx $$

$$ b_n = \frac{1}{\pi} \int _{-\pi} ^{\pi} f(x)\sin(nx)dx $$

여기서 오일러 전개를 이용하면 $$ \cos(nx) = \frac{1}{2} (e^{jnx} + e^{-jnx} ) $$

$$ \sin(nx) = \frac{1}{2j} (e^{jnx} - e^{-jnx} ) $$

이므로 $c_n$을 다음과 같이 정의하면 $$ c_n = \frac{1}{2}a_n - \frac{1}{2}jb_n $$

$$ c_{-n} = \frac{1}{2}a_n + \frac{1}{2}jb_n $$

다음 식에 의해서 $$ a_n\cos(nx) + b_n \sin(nx) = \left(\frac{a_n}{2}-\frac{b_n}{2}j \right) e^{jnx} + \left(\frac{a_n}{2}+\frac{b_n}{2}j \right) e^{-jnx} $$ $f$를 다음과 같이 다시 쓸 수 있다. $$ f(x) = \sum _{n= -\infty} ^{\infty} c_n e^{jnx} $$

푸리에 계수 $c_n$은 $a_n$과 $b_n$ 성분을 모두 가지고 있다. 따라서 실수가 아닌 복소수로 존재할 수 있다. $\vert c_n \vert$은 각 주파수 성분의 진폭이 되고, $\angle c_n$은 각 주파수 성분의 위상이 된다.

푸리에 변환

일반적으로 연속 시간 위에서 정의되는 비주기 신호의 푸리에 분석을 연속 시간 푸리에 변환(CTFT: Continuous Time Fourier Transform)이라고 한다.

비주기 연속 신호의 푸리에 변환은 다음과 같이 계산한다. $$ X(j\Omega) = \int _{-\infty} ^{\infty} x(t)e^{-j\Omega t} dt $$

이산 푸리에 변환

이산 푸리에 변환은 연속 푸리에 변환을 영차 홀드 근사를 이용해서 리만 근사한 것으로 다음과 같이 계산한다. $$ X[k] = \sum _{n=0} ^{N-1} x[n] e ^ { - j \cdot 2\pi\cdot \frac{2\pi}{N} n k} $$

주파수 응답 함수

임펄스 응답 $h[n]$을 갖는 선형 시불변 이산 시스템은 다음과 같이 나타낼 수 있었다. $$ y[n] = x[n]*h[n] = \sum _{k=-\infty} ^{\infty} h[k]x[n-k] $$ 여기에 $x[n]$을 복소 지수 함수라고 생각하면 선형 시불변 시스템의 출력은 다음과 같이 표현된다. $$ x[n] = e^{j\omega n} $$

$$ y[n] = \sum _{k=-\infty} ^{\infty} e^{j\omega (n-k)} = \left ( \sum _{k=-\infty} ^{\infty} h[k]e^{-j\omega k} \right ) e^{j\omega n} $$

여기서 괄호로 묶은 부분은 임펄스 응답 $h[k]$에 대한 푸리에 변환이다. 따라서 이 계수를 $H(e^{j\omega})$로 대체하면 다음과 같이 다시 쓸 수 있다. $$ y[n] = H(e^{j\omega}) e^{j\omega n} $$ 이 계수는 $n$에 대하여 독립이므로 이는 복소 상수로 볼 수 있고, 이렇게 입력의 상수배가 되면서 출력을 만들어낼 수 있는 함수를 그 함수의 고유 함수(eigen function)라고 한다. 따라서 복소 지수 함수 $e^{j\omega n}$은 LTI 시스템의 고유 함수가 된다.

이 때, $H(e^{j\omega})$ 를 주파수 응답 함수라고 한다.

무한 임펄스 응답 저주파 패스 필터(IIR LPF: Infinite Impulse Response Low Pass Filter)

주파수 대역에서 이상적인 저주파 패스 필터를 설계하면 다음과 같은 주파수 응답 함수를 생각할 수 있다. $$ H(e^{j\omega}) = \begin{cases} 1 & (|\omega| < \omega_c)\ 0 & otherwise \end{cases} $$ 이 함수의 푸리에 역변환인 임펄스 응답 $h[n]$은 다음과 같다. (이는 구형파 신호의 푸리에 변환이 $sinc$함수인 것과 푸리에 변환 결과에 한 번 더 푸리에 변환을 수행한 결과와 원본 신호 간의 관계를 따져 보면 알 수 있다.) $$ h[n] =\frac{\sin(\omega_c n)}{\pi n} $$ 이 때, $w_c$는 $0&lt;w_c &lt; \pi$를 만족하도록 정규화되어 있기 때문에 $w_c =\pi f_c$에서 , 컷오프 주파수 $f_c$는 나이퀴스트 주파수 $f_s = \frac{N}{2}$에 대하여 컷오프 될 주파수 샘플인 $N_c$를 정하고 $f_c = \frac{N_c}{f_s}$ 으로 잡아야 한다. 따라서 실습에서는 임펄스 응답 $h_n$의 $w_c$를 다음과 같이 잡았다. $$ w_c = \frac{2\pi N_c}{N} $$

# IIR 임펄스 응답
def h(n: int):

    f_c = 600
    w_c = 2*pi*f_c / SAMPLE_RATE

    if n==0:
        return w_c / pi
    else:
        return sin(w_c * n) / pi / n