デジタル信号処理④ - AIST Digital Human Research

デジタル信号処理④
2002.6.4
前回までの内容
1. サンプリング定理
2. フーリエ変換
1.
2.
3.
4.
5.
実フーリエ級数
複素フーリエ級数
フーリエ変換
離散フーリエ変換 (DFT(FFT))
DFTを使った実際のフーリエ変換
①②
3. 線形システム
1. ラプラス変換
2. Z変換
•
インパルス応答
③
今回の内容
1. Z変換の復習
•
•
•
Z変換
Z変換を行う上での重要な公式
Z変換の具体例
2. デジタルフィルタ
1. フィルタの種類
2. FIRフィルタ
Z変換の復習①
Z変換とは、パルス状の信号系列を一種の周波数のパラメータ
と考えられる変数z(複素数領域)の関数へ変換すること。
これにより、デジタル信号処理を行うシステムの設計と解析が
容易になる。
デジタル信号系列をx(n)とするとき、z変換は、以下。

X ( z )   x ( n) z  n

(Tは、サンプリング周期)
nを0,1,2,3からなる自然数とすると、

X ( z )   x ( n) z
n 0
n
Z変換の復習②
• ラプラス変換の式で、Z  esT なる置き換えを
行い、有理多項式に変換したものであるが、
ラプラス変換を使うことなく、デジタル信号
系列x(n)から容易に求めることが可能。
Z変換の復習③
Z変換を行う上で重要な公式
Z [ x(n)]  X ( z ) のとき、
• σ
Z[(n)]  1
• 線形性
Z[ax1 (n)  bx2 (n)]  aX1[ z]  bX2 [ z]
• 推移定理
Z[ x(n  m)]  X ( z) z
m
• たたみ込み

Z [  x ( k ) h ( n  k )]  X ( z )  H ( z )
k 0
 H ( z)  X ( z)
Z変換の復習④
x
Z変換の具体例
2
1
n
0
X ( z )  1.0  1.0 z 1  1.0 z 2  0.0 z 3  0.0 z 4
5
6
 2.0 z  2.0 z  2.0 z
1
 1.0  1.0 z  1.0 z
6
 2.0 z  2.0 z
7
2
7
 2.0 z
5
たたみ込みの復習 x(3)  (3  3)
インパルス応答
R1
時不変性
線形性
x
x(t)
C
R2
インパルス応答は、以下。
 R1  R2 
1
h(t ) 
exp 
t 
CR1
 CR1R2 
y(t)
 0.7  (0)
1.0
0.5
0.7
0.8
0.3
x(3)  h(3  3)
 0.7  h(0)
y
01 2 3 4
重ね合わせ
y(3)  x(0)h(3)  x(1)h(2)  x(2)h(1)  x(3)h(0)
 0.5h(3  0)  1.0h(3  1)  0.3h(3  2)  0.7h(3  3)
たたみ込みの具体例①
たたみ込みから出力パルスを求める
線形時不変システムにおいて、
インパルス応答がh(n)で、入力がx(n)であるとき
出力y(n)は以下のようになる。

y ( n)   h( n  k ) x ( k )
k 0
例えば、インパルス応答と入力データが以下のとき
{h(n)}  {0.2,0.1,0.5}
{x(n)}  {3,2,1,5}
たたみ込みの具体例②
たたみ込みから出力パルスを求める

y (0)   h(0  k ) x(k )
k 0
 x(0)h(0)  x(1)h(1)  
 3  0.2  0  0.6
{h(n)}  {0.2,0.1,0.5}
{x(n)}  {3,2,1,5}


y ( n)   h( n  k ) x ( k )
y (1)   h(1  k ) x(k )
k 0
k 0
 x(0)h(1)  x(1)h(0)  x(2)h(1) |  
 3  0.1  2  0.2  0.7
{ y(n)}  {0.6,0.7,1.9,2.1,1.0,2.5}

y (2)   h(2  k ) x(k )
k 0
と求めることが可能
 x(0)h(2)  x(1)h(1)  x(2)h(0)  x(3)h(1) 
 3  0.5  2  0.1  1 0.2  1.9
たたみ込みの具体例③
Z変換を利用して出力パルスを求める


{h(n)}  {0.2,0.1,0.5}
Z [ y(n)]  Z  h(n  k ) x(k )
 k 0
 H ( z)  X ( z)
{x(n)}  {3,2,1,5}
x
h
n
X ( z)  3  2z 1  z 2  5z 3
Y ( z)  X ( z)  H ( z)
n
H ( z)  0.2  0.1z 1  0.5z 2
 0.6  0.7 z 1  1.9 z 2  2.1z 3  1.0 z 4  2.5 z 5
{ y(n)}  {0.6,0.7,1.9,2.1,1.0,2.5}

デジタルフィルタ①
• フィルタ
様々な信号から望みの信号だけを取り出
す技術 (濾紙やコーヒーのフィルタ)
• デジタルフィルタ
信号をパルス列として表現し、数値演算処
理に基づいて行うフィルタ
デジタルフィルタ② フィルタの種類
通過
フィルタの呼称
除去
ローパスフィルタ
(低域通過フィルタ)
ハイカットフィルタ
(高域除去フィルタ)
f
バンドパスフィルタ
(帯域通過フィルタ)
色がついた部分を除去
するフィルタを
バンドエリミネーション
フィルタ
(帯域除去フィルタ)
f
振幅スペクトル
振幅スペクトル
振幅スペクトル
特性
ハイパスフィルタ
(高域通過フィルタ)
f
ローカットフィルタ
(低域除去フィルタ)
デジタルフィルタ③
フィルタの周波数特性の表し方
振幅スペクトル
通過域
減衰量
カットオフ
周波数
カットオフ
周波数
f
FIRフィルタ①
• FIRフィルタ(Finite Impulse Response)
そのインパルス応答が有限時間長で表されるも
の
• フィルタの設計とは、フィルタ係数h(n)を決めるこ
と
N 1
y (n)   h(i ) x(n  i )
i 0
N 1
H ( z )   h(i ) z
i 0
i
y(n),x(n)とフィルタ係数h(n)
の関係
FIRフィルタの伝達関数
FIRフィルタ②
4次(フィルタ係数5個)のFIRフィルタのブロック線図
1サンプリング時間の遅延
1個前のデータ
x(n)
Z 1
Z 1
h0
Z 1
h1
+
Z 1
h2
+
h3
+
h4
+
y(n)
4
y ( n)   h( k ) x ( n  k )
k 0
 h(0) x(n)  h(1) x(n  1)  h(2) x(n  2)  h(3) x(n  3)  h(4) x(n  4)
FIRフィルタ③
FIRフィルタの具体例(フィルタ長13)
h(n)
•
デジタルフィルタの
インパルス応答の例
{h(n)}={0.0000 0.0637 -0.0000 -0.1061 0.0000 0.3183 0.5000 0.3183 0.0000
-0.1061 0.0000 0.0637 0.0000 -0.1061 -0.0000 0.0637 0.0000 -0.0455 -0.0000
0.0354 0.0000}
FIRフィルタ④
振幅スペクトル
FIRフィルタの具体例
フィルタの周波数特性
ローパスフィルター

2f

fs:サンプリング周波数
横軸は、
fs なる正規化がされている。
FIRフィルタ⑤
FIRフィルタの具体例
実際に使ってみる
フィルタのインパルス応答h(n)と入力x(n)から、
N 1
y (n)   h(i ) x(n  i )
i 0
の計算を行う。
FIRフィルタ⑥
FIRフィルタの具体例(チャープ信号)
• 入力x(n)として、「チャープ信号」と呼ばれる信号
を与える。チャープ信号とは、時間の経過ととも
に周波数が高くなる特殊な信号。
だんだん高い音になる
「全ての周波数を
ほぼ均等に含んだ
信号」であることがわかる
FFTで
振幅スペクトル
を出す
波形
FIRフィルタ⑦
FIRフィルタの具体例
チャープ信号を利用してフィルタの周波数特性を調べてみる
• チャープ信号を入力することで、どの周波数がど
のように影響を受けるかが直感的に理解できる。
FIRフィルタ
0[Hz]
4[kHz]
2[kHz]ぐらいから
急に小さくなっている。
FIRフィルタの具体例
MATLABを使ったプログラム例
%チャープ信号の作成
• Ts=1;
%計測時間1[s]
• Fs=8192;
%サンプリング周波数8192[Hz]
• delta_t=Ts/Fs;
%サンプリング周期
• t=linspace(0,Ts-delta_t,Ts*Fs); %時間
• y=sin(2*pi*Fs/4*t.*t);
%チャープ信号の作成
• sound(y,Fs);
%音を聞く
FIRフィルタの具体例
MATLABを使ったプログラム例
%チャープ信号のスペクトルを見る
• ffty=fft(y,Fs);
%FFT
• delta_f=1/Ts;
%周波数分解能
• f=linspace(0,Fs-delta_f,Fs);
%グラフの横軸データ
• plot(f,abs(ffty),‘LineWidth’,2.0);axis([0 4096 0 70]) %描画
FIRフィルタの具体例
MATLABを使ったプログラム例
•
•
•
for n=-6:6,
h(n+7)=1/2*sinc(n/2)
end
%FIRフィルタの係数の作成
•
•
•
•
•
•
•
y(1:13)=0;
for i=14:8192,
y(i)=0;
for j=1:13,
y(i)=y(i)+h(j)*x(i-j);
end
end
%出力y(n)の作成
•
plot(t,y)
%y(n)を描画する
•
sound(y,Fs)
%y(n)を音として聞く
フーリエ級数法①
希望のデジタルフィルタの周波数特性H(ω)を
フーリエ逆変換することでFIRフィルタのフィルタ係数
を求める方法
FIRフィルタが無限のフィルタ係数を持つと考えると、

y ( n) 
 h(i) x(n  i)
i  
H ( z) 

 h(i) z
i  
i
ze
j
なる置き換えを行い
H(ω)を求めると、
H () 

 jk
h
(
i
)
e

k  
複素フーリエ級数展開の式
に似ている
フーリエ級数法②
希望のデジタルフィルタの周波数特性H(ω)を
フーリエ逆変換することでFIRフィルタのフィルタ係数
を求めることが可能
1 
jn
h( n) 
H ()e d

2  
フーリエ級数法③
問題例
• 次のようなローパスフィルタを設計する。
遮断周波数100[Hz], サンプリング周波数
400Hzとする。
H(f)
1
0
100
f[Hz]
フーリエ級数法④
問題例
1. 周波数をサンプリング周波数で正規化し
ておく fc=100/400
2. 各周波数ωc =2πfc=0.5π
H(f)
1
H ()  1, |  | 0.5
 0,0.5 |  | 
フーリエ変換を
行うためy軸対称
にしておく。
フーリエ級数法⑤
問題例
• フーリエ逆変換する
k
sin
1 0.5  jk
1
2
h( k ) 
e
d


2 0.5
2 k
2
k=0,±1,±2, ±3, ±4,…を代入してh(k)を求める。
フーリエ級数法⑥
問題例
実際に、h(k)を使うときは、例えば、k=0, ±1 ±2, ±3,…,±Nのとき、
h(-N) →h(0)
h(-N+1)→h(1)
…
h(N-1) →h(2N-1)
h(N) →h(2N)
と係数をずらし、以下のFIRの式
で使える形にする。
N 1
y (n)   h(i ) x(n  i )
i 0
FIRフィルタの特徴
利点
1. 設計が容易である
•
フーリエ級数を利用した設計法がある
2. 直線位相特性をもつ
•
フィルタ処理を行った後の波形は、位相に関しては、
元波形を平行移動させた結果が得られる
3. 常に安定なフィルタを実現できる
•
ある周波数について発散することがない
欠点
• 鋭いカットオフ特性を得るには、フィルタの係数
を多くしなければならない。ソフトウェアの処理
の負荷を増大させる原因となる。
これまでの内容
1. サンプリング定理
2. フーリエ変換
1.
2.
3.
4.
5.
実フーリエ級数
複素フーリエ級数
フーリエ変換
離散フーリエ変換 (DFT(FFT))
DFTを使った実際のフーリエ変換
①②
3. 線形システム
1. ラプラス変換
2. Z変換(インパルス応答・たたみ込み)
4. デジタルフィルタ
1.
FIRフィルタ
③④
これまでの内容
試験日程のお知らせ
• 2002.7.9に実施
• 持ち込みの可否に関しては、後日、連絡がある
• デジタル信号処理については、
– FFTの結果得られるフーリエ係数の解釈の方法
– デジタル信号系列のz変換
– z変換を利用した周波数応答の解析
z exp( jT )
– FIRフィルタ
の範囲から出題予定。
資料:http://www.dh.aist.go.jp/~ynishida/lecture
質問:[email protected]