周波数領域での積は、時間領域では畳み込みに対応する。この関係性を利用すると、周波数を好みにカットする関数を作って、それを逆フーリエ変換するという方法で、フィルタを作ることができる。

逆フーリエ変換_フィルタ

このようなフィルタの設計法を窓関数設計という。

理想的なフィルタは無限に続きまた因果性も満たさない

ローパスフィルタを考えてみる。

理想的なローパスフィルタは、ある周波数以上の信号をバッサリ遮断するものだから、矩形をかけることにあたる。

矩形を逆フーリエ変換することでフィルタが得られるわけだが、矩形の逆フーリエ変換はsinc関数になる。

cut_sinc_sinc

sinc関数は無限に続いているし、時刻が負の時も値をもっているので因果性をみたさない。したがって理想的に周波数を遮断するフィルタは物理的に不可能といえる。

sinc関数を途中でカットする

妥協案として考えられるのが、無限に続くsinc関数をあるところでカットして有限にすること。あと、時刻が負の時に値も持つ件を解決するために、時間が正の方向にシフトする。そうすることで、実現可能なフィルタとなる。ただ、当然ながら周波数を遮断する性能は落ちる。では、上記のようなカットしたりシフトしたりする影響はどこにでてくるのか。

その影響は2つある。ひとつは遮断特性、もうひとつは振幅が一定でなくなる(リプル)。

cut_sinc_rec

理想的なフィルタでは通過する周波数と遮断される周波数が境界で完全にわかれていたが、途中でカットしたフィルタではなだらかに遮断されていくことになる。

また、通過する周波数において振幅が一定ではなく波を打つようになる。

遮断特性はフィルタの次数をあげることで急峻になっていくが、リプルは次数をあげて消える気配はない。

cut_sinc_cut_sinc

この点を考えてみる。

そもそもsinc関数をある範囲でバッサリカットするということは、矩形窓関数をかけることにあたる。したがって、周波数領域で矩形窓関数のフーリエ変換したものを理想的フィルタの形状に畳み込んでいることになる。

しだって、窓関数の次数をあげるとメインローブの幅が狭くなるので、通過と遮断の立ち上がりが急峻になる。

一方、リプルはサイドローブの影響。したがって、次数をあげてもサイドローブはなくならないので、リプルもなくならない。リプルを小さくするには、サイドローブが小さい窓関数でsinc関数をカットしてやればよいということになる。

例えば、ハミング窓を作用させることでリプルを小さくすることができる。ただし、遮断特性はなだらかになる。

cut_sinc_hamming
"""
Design filter by window
"""
# -*- coding: utf-8 -*-
import numpy as np
import scipy.signal as signal
import matplotlib.pyplot as plt
import sys, os
sys.path.append(os.pardir)  

fftpt = 1024
tap_num = 20
#sinc function
x = np.linspace(-tap_num,tap_num,fftpt)
s = np.sinc(x)
bc.showPlot(s, title='sinc func')
bc.showSpectrum(s)
#fft sinc in varing tap
tap_list = [3,5,10,15,20]
bin_per_tap = fftpt/(2*tap_num)
cut_sinc = np.zeros(fftpt, 'float')
f_list = np.linspace(-np.pi, np.pi, fftpt)
w_hamming = signal.hamming(fftpt)
for tap in tap_list:
    cut_sinc = s.copy()
    cent_bin = int(fftpt/2)
    cut_bin = int(tap * bin_per_tap)
    cut_sinc[:cent_bin - cut_bin] = 0
    cut_sinc[cent_bin + cut_bin:] = 0
    #bc.showPlot(cut_sinc, title='sinc func(cut)')
    cut_sinc = w_hamming * cut_sinc
    fftdt = np.fft.fft(cut_sinc, fftpt)/np.sqrt(fftpt)
    fftdt = np.fft.fftshift(fftdt)
    plt.plot(f_list, 20*np.log10(np.abs(fftdt)))
plt.legend(tap_list)
plt.title('fft sinc func (cut)')
plt.show()