본문으로 바로가기

전 아주 예전부터 제 블로그에 필터에 관한 글을 올리고 있었습니다. 물론 단편적인 글들이었지만요. 그 중에는 실제 자이로센서와 가속도센서를 융합하는 상보필터[바로가기]를 다루기 시작했고, 그 후 MATLAB이나 Python으로 구현하는 1차 혹은 2차 필터들에 대한 글들을 올렸습니다. 물론 이전에 올린 글들이 많지만, 문득 필터들에 대한 이야기를 1차 저역/고역 통과필터, 2차 필터, Band Pass, Band Stop 필터등등에 대해 한 번 쫘~악 정리하고 싶어지더라구요. 그래서 이런 것들을 정리할 겸 필터를 대상으로 연재를 시작할까 합니다.^^

오늘 다룰 이야기는 1차 저역 통과필터와 1차 고역 통과필터 흔히 LPF(Low Pass Filter), HPF(High Pass Filter)에 대한 이야기 입니다.

흔히 Hz 단위의 차단주파수(f_cut)가 결정되면 위의 간단한 변환으로 차단 각주파수가 계산됩니다.

이를 1차식으로 s-domain에서 표현하면 LPF는 위와 같이...

HPF는 위와 같이 표현됩니다.

만약 흔히 알려진 시정수(tau)로 LPF와 HPF를 표현하고 싶다면... 위 관계식을 이용해서

이렇게 표현할 수도 있습니다. 뭐 여하튼.. 1차의 저역통과필터식과 고역통과필터식은 합하면 '1'이 됩니다.^^.

1차 LPF와 HPF의 Bode선도를 그려보면 위와 같습니다.

위 보드선도를 그린 Python 코드 보기

f_cut 지점에서 -3dB가 떨어지는 것을 확인할 수 있습니다.

이제 이를 디지털 필터, Digital Filter로 변환하는 과정을 한 번 보도록 하겠습니다. 위 수식처럼 애초 s-domain에서의 LPF를 이산영역에서 표현할 수 있습니다. 위 과정의 상세한 전개는 [바로가기]에서 다루었습니다.

위 식은 HPF의 이산영역 수식입니다.

식 전개의 편의성을 위해 tau_ts를 위와 같이 두고, LPF와 HPF의 이산영역에서의 수식을 다시 편집하고 이 수식을 z-domain에서 표현하면...

와 같습니다. 즉, 설정하고 싶은 주파수 f_cut을 이용해서 알 수 있는 tau와 이산영역이 되면서 당연히 생겨난 ts(sample time)을 공식으로 위 수식을 완성할 수 있을 겁니다.

일반적으로 biquad라고 알려진 Direct2Form에 일반화해서 적용하기 위해 위의 수식을 사용하겠습니다.

필터는 위와 같은 모양이 될 거구요^^

그러면 Low Pass Filter의 계수는 위의 공식으로~~~

High Pass Filter의 계수는 위의 공식으로 결정됩니다. 이제 digital 영역에서 손쉽게 1차 저역/고역 통과필터의 계수를 찾을 수 있겠죠^^ 여기서... 

tau_ts는 위의 식입니다.

위 시험신호를 생성한 Python 코드보기

이제 시험신호를 만들어서 한 번 적용해 봐야죠^^& 위 시험신호는 10kHz 샘플 주파수에서 생성된 것으로 가정했구요. 기본적으로 2Hz 주 성분에 10Hz부터 500Hz까지 50Hz 간격으로 신호를 생성해서 합했습니다. 이 신호의 FFT결과를 보면

위와 같이 의도한 대로 깔끔하게 나타나네요^^

위 FFT 결과를 만들어낸 Python 코드 보기

앞부분에서 사용한 대로 공식화하고 이 필터의 계수를 찾고 적용하는 부분에 대해 Python으로 코드를 확인해 보겠습니다.

def get1stFilterCoeffi(f_cut, ts, isLPF):
    from numpy import pi
    
    w_cut = 2*pi*f_cut
    tau = 1/w_cut
    
    tau_ts = 1/(tau+ts)
    if isLPF=='LPF':
        return -tau*tau_ts, ts*tau_ts, 0
    elif isLPF=='HPF':
        return -tau*tau_ts, tau*tau_ts, -tau*tau_ts

먼저 위 코드는 1차 LPFHPF의 계수를 direct2form에 맞춰 계산하는 코드입니다. 뭐 크게 어렵지 않습니다.^^.

def direct2FormModel(data, a1, a2, b0, b1, b2):
    from numpy import zeros, arange
    
    result = zeros((len(data),))
    timeZone = zeros((len(data),))
    
    for n in arange(2, len(data)):
        sum0 = -a1*timeZone[n-1] - a2*timeZone[n-2]
        timeZone[n] = data[n] + sum0
        result[n] = b0*timeZone[n] + b1*timeZone[n-1] + b2*timeZone[n-2]
        
    return result


def d2f_1st(data, a1, b0, b1):
    return direct2FormModel(data, a1, 0, b0, b1, 0)

먼저 direct2FormModel 함수는 2차 필터까지 고려한 것으로 d2f_1st라는 함수로 2차에서만 사용하는 계수는 0을 넣은 것입니다. 아무튼... 위 함수는 direct2form에 맞춰 필털를 동작시키는 함수입니다. 계수는 아까 get1stFilterCoeffi함수로 찾았으니까 적용할 수 있는거죠^^...

이제 이 함수들을 이용해서 LPF를 적용해보면

이런 결과가 나타납니다. 

a1, b0, b1 = get1stFilterCoeffi(200, Ts, 'LPF')
dataLPF = d2f_1st(inputSig, a1, b0, b1)

위 저역통과필터의 결과는 이미 만들어둔 함수로 인해 위 코드 딱 두줄이네요^^ 200Hz를 차단주파수로 가지게 해서 필터 계수를 얻고 이를 적용해서 데이터를 얻었네요.

LPF가 적용된 결과를 그린 Python 코드 보기

실제 FFT 결과도 LPF에 어울리게 나타납니다.

어때요? 그렇죠?^^ 실제 200Hz 다음에 나오는 값이 대략 0.7쯤인데 이를 20log10을 적용시켜보면 -3dB가 되거든요^^. 

위 FFT 결과를 그린 Python 코드 확인하기

이제 HPF를 적용한 경우를 볼까요?

a1, b0, b1 = get1stFilterCoeffi(200, Ts, 'HPF')
dataHPF = d2f_1st(inputSig, a1, b0, b1)

위와 같이 HPF라고 알려주고 계수를 받아 적용해보면~

요렇게 나타납니다. 저주파 영역이 사라졌기 때문에 필터의 결과가 센터(0)에 딱 붙어 있네요^^

위 HPF의 결과를 그린 Python 코드 보기

그리고... FFT의 결과도 예상대로 입니다.^^

FFT 결과 그래프 그린 Python 코드 보기

이렇게 해서 시간영역에서 1차 저역통과필터와 1차 고역통과필터의 개념을 이야기했고, 이를 digital로 구현하기 위해 각 계수를 찾는 공식을 도출하고 이를 code로 표현하는 것 까지 다루어 보았습니다. 실제 예시로 Python을 사용했구요^^ 오늘 글 전체에서 다룬 Python의 전체 예시는 Github[바로가기]에 두었습니다....


댓글을 달아 주세요

  1. BlogIcon IT넘버원 2017.01.06 04:28 신고

    아무리봐도 천재이신듯 합니다.^^

  2. BlogIcon 핑구야 날자 2017.01.07 07:33 신고

    필터를 연재 하시는군요 어렵지만 잘 보고 갑니다

  3. BlogIcon 지후대디 2017.01.08 22:35 신고

    수식보다 예시로 남긴 코드가 더 보기가 편한게 참 저도 직업병인가 봅니다. ^^

  4. BlogIcon PEACE- 2017.02.02 20:45 신고

    이 글을 읽으며 또 좋은 지식얻어갑니다~^^
    차단 주파수를 LPF와 HPF를 통해 찾아내셨는데
    혹시 f_cut(차단주파수)를 찾는데 사용된 데이터는 어떻게 되나요?..
    단순히 LPF와 HPF로 저렇게 만들어 진다는게 이해가 안갑니다.!

    • BlogIcon PinkWink 2017.02.03 11:01 신고

      혼돈이 있을 수 있습니다.^^
      이 부분은 (적어도 이 문서는) 판단에 대한 부분입니다. 차단주파수를 사용자가 결정한 이후,
      결정된 차단주파수를 이용해서 어떻게 필터를 설계할 것인가를 다루거든요.^^

  5. 핑크님 2017.12.17 21:40 신고

    핑크님 오랜만입니다.
    너무 해결안되는 궁금중이 있어 찾아왔습니다 .총 두가지의 질문이 있습니다.
    1.상보필터든 머든 cutoff freq 는 어떻게 구해내는 건가요 ?
    아무리 고민해봐도 잘모르겠네요 ..ㅠㅠ
    그리고, 제어공학 책을 몇 번을 반복하고 있는데
    Closed loop 를 이용할 때의 장점이나 Open loop sys이용할 때 장점이라고있는데
    이용한다라는 게 뭘 이용한다는지 잘 와닿진 않네요.

    제가 실제로 핑크님 강의를 평소 많이 찾아보면서 LQR 제어로 외륜로봇 구현에도 성공하였습니다.
    그런데도 그 LQR제어 마저 저런 Closed loop 활용 그런 개념이 잘 이해가 가지않네요.
    2. Closed loop TF 과 open loop TF을 이용한다라 함은 .. 멀 어떻게 활용할 수 있는걸가요 ?


    읽어주셔서 감사합니다.
    아차. 핑크님 핑크님은 혹시 메일 주소도 비공개이신지요 ?
    한번씩 관련 자료로 보여드리고 여쭤보고 싶은 적이 너무 많던지라 ..
    아무튼 항상감사드립니다

    • BlogIcon PinkWink 2017.12.18 14:03 신고

      cutoff freq.는 사용자가 정하는 겁니다.^^. 내가 원하는 주파수를 차단하는 거죠^^ 뭐 통과시키든 말든...

      cutoff freq를 어떻게 결정하느냐... 시스템의 노이즈 제거를 예로 든다면... 어디까지 노이즈를 잘라버릴까..를 디자이너가 결정해야죠.. 원신호를 좀 손해를 보더라도 여기까지 잘라야지.. 하고 결정하는 겁니다.^^

      그러나 이어진 질문은 이해가 안됩니다. cutoff 이야기하다가 갑자기 왜 closed/open loop으로 가는 지를 모르겠네요.ㅠㅠ.

      LQR도 closed loop controller입니다. 현재와 과거의 상태를 받아와서 현재의 제어기어 넣어서 현재의 제어입력을 결정하는 모든 제어기는 closed-loop입니다. 뭐 PID도 그렇구요...

      open-loop도 있습니다. 그냥 현재값이 얼마면 이정도 값으로 출력해라``고 해주면 되죠... open-loop도 정말 많이 사용하는 방식입니다.

      두 방식 모두 장단점이 있지만, 그건 대상 시스템의 영향을 많이 받구요...

      그리고 이메일...은... 음.. 한번 공개한 적이 있는데.. 제가 일을 못할 정도였습니다. 너무 많이 와서ㅠㅠ. 답변을 못하면 심지어 항의 메일까지ㅠㅠ. 그래서 댓글로만 질문을 받는답니다. 또 저는 댓글도 공개(open)의 일부로 봐서 댓글과 제 답변을 또 공개하고 싶거든요...^^ 좀 더 많은 자료는 별도로 블로그에 만드셔서 저에게 질문하셔도 됩니다.^^