ONOSOKKI -- info channel 5月号 第68号「計測コラム」

━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━
計測コラム
emm68 号用
━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━
波形と FFT−6
6.フーリエ級数・フーリエ変換
(1) フーリエ級数
フーリエ係数
a0
+ a1 cos x + b1 sin x + a 2 cos 2x + b2 sin 2x +・・・
2
f (x ) =
M
=
∑ (a
m
cos mx + bm sin mx )
m =0
am =
bm =
1
π
1
π
2π
∫ f (x ) cos mxdx
0
2π
∫ f (x ) sin mxdx
0
(2)1 周期を T としたフーリエ級数
フーリエ係数
上記(1)の x を x = 2πft = ωt、1 周期 2π = T と置き換えて;
f (t ) =
a0
+ a1 cos ω t + b1 sin ω t + a2 cos 2 ωt + b2 sin 2ω t +・・・
2
M
=
∑ (a
m
cos mω t + bm sin mω t )
m =0
2
T
2
bm =
T
am =
T
1、・・・
∫ f (t ) cos mω tdt m = 0、 0
T
1、 2、・・・
∫ f (t ) sin mω tdt m = 0、 0
(3)複素フーリエ級数・複素フーリエ係数
オイラーの公式
e jθ + e − jθ
2
jθ
e − e − jθ
sin θ =
2j
cos θ =
を考慮して(2)は;
-1-
f (t ) =
a0 M ⎧ 1
1
⎫
jmω t
+
+ (am + j bm ) e − jmω t ⎬
⎨ (am − j bm ) e
2
2 m =1 ⎩ 2
⎭
∑
M
∑
= C m e jmω t M → ∞
m=−M
C 0 =
Cm =
1
T
T
a0
a − j bm
a + j bm
、 C m = m
、 C - m = m
2
2
2
∫ f (t) e
0
− jmω t
dt
<補足>
1 T
1
f (t ) e − j ( − m ) ω t dt =
T 0
T
一方、f (t )は実関数なので
C−m =
* ⎛ 1
C m =⎜
⎝T
よって
∫
∫
T
0
*
1
⎞
f (t ) e − jmωt dt ⎟ =
T
⎠
∫
T
∫
T
0
0
f (t ) e jm ω t dt
(
)
f*(t ) e − jmωt dt =
*
1
T
∫
T
0
f (t ) e jmωt dt
*
C−m = Cm
実関数の複素フーリエ係数は Cm、m = 1 ∼ M を計算すれば、あとは複素共役 Cm*
をとればよいことになります。
なお、 Cmの実数部と虚数部は、複素フーリエ級数展開することにより、正負の周波数
の 1/2となります。
成分にわかれ、フーリエ係数 am、bm (4)フーリエ変換、フーリエ逆変換
T を無限大に拡大することで、周期性のない波形【非周期関数】でもフーリエ級数
に展開できるようにしたものです。
フーリエ係数 Cm を求める式も 1 つの関数と考え、これをフーリエ変換といいます。
f (t)はフーリエ変換の逆の計算をする関数と考えて、フーリエ逆変換といい、フー
リエ変換、フーリエ逆変換を一対の式として考えます。
-2-
フーリエ変換
F (ω ) =
∫
∞
−∞
f (t ) e − jω t dt =
∫
∞
−∞
f (t ) cos ωtdt − j
∫
∞
−∞
f (t )sin ωtdt
= Re[F (ω )] + j Im[F (ω )]
フーリエ逆変換
f (t ) =
1
T
∫
∞
−∞
F (ω ) e jω t dω
実部 Re 、虚部 Im は
Re[F(ω )]=
∞
1
T
∫
Im[F (ω )] = −
1
T
f (t ) cos ω tdt
−∞
∫
∞
f (t )sin ω tdt
−∞
(5)離散的フーリエ変換
フーリエ逆変換
実用性を考えると T は有限をとることになります。また、FFT アナライザでは AD
変換してサンプリングされたデータ列の有限個 N をとり、フーリエ変換が実行さ
れます。サンプリングすると【離散的】
(連続していない飛び飛びの値)になります
から、T = Nh(h はサンプリング間隔の時間)より T の代わりに N 個とすると、離散
的フーリエ変換、フーリエ逆変換の第k項を考える(周波数は kωになるので)
。
∫
∞
−∞
N −1
は
∑
、t は nh に置き換えて;
n =0
F (kω ) =
1
N
N −1
∑ f (nh ) e
− jkω nh
n =0
=
N −1
1 N −1
f (nh ) cos kωnh − j
f (nh )sin kωnh
N n =0
n =0
fk ( nh ) =
1
N
∑
∑
N −1
∑ F (kω ) e
jkω nh
n =0
実部 Re 、虚部 Imは
Re[F (kω )] =
1
N
Im[F (kω )] = −
N −1
∑ f (nh )cos kωnh
n =0
1
N
N −1
∑ f (nh )sin kωnh
n =0
-3-
(6)離散的フーリエ級数、フーリエ係数
サンプリングされた離散値でフーリエ係数を求めるには、上記(5)を展開して次
式で計算します。フーリエ級数、フーリエ係数と実部、虚部の関係を対比して記
します。
<フーリエ級数、フーリエ係数>
f (t ) =
a0
+ a1 cos ω t + b1 sin ω t + a2 cos 2 ωt + b2 sin 2ω t +・・・
2
=
M −1
∑ (a
m
cos mω t + bm sin mω t )
m =0
2
T
2
bm =
T
am =
T
1、・・・
∫ f (t ) cos mω tdt m = 0、 0
T
1、 2、・・・
∫ f (t ) sin mω tdt m = 0、 0
第k項のフーリエ係数は;
ak =
bk =
2
N
2
N
N −1
∑ f (nh )cos kωnh
n =0
N −1
∑ f (nh )sin kωnh
n =0
第 k 項のフーリエ変換をフーリエ係数 an、bn に一致するように 2/N で考えると;
F (kω ) = Re[F (kω )] + j Im[F (kω )]
Re[F (kω )] =
2
N
Im[F (kω )] = −
N −1
∑ f (nh ) cos kωnh = a
k
n =0
2
N
N −1
∑ f (nh ) sin kωnh = −b
k
n =0
第 k 項のフーリエ係数の式または Re [ F (kω) ]、Im [F (kω)]の式を使うと表計算ソ
フトで計算することができます。
ω = 2πfo の fo は1周期 T の逆数です。T = Nh、N はサンプル数、h はサンプル時
間間隔ですから、h が一定なら N は周波数分解能 fo に関係します。
-4-
(7)フーリエ係数の計算例
前号で 10Hz と 20Hz の cos を合成した;
f (t ) = f ( nh ) = cos (20π n h + θ ) + cos (40π n h + φ )
30 × 2π
( rad)
θ = 30 (deg) = 360
45 × 2π
φ = 45 (deg) =
( rad)
360
のサンプルデータを作ってみました。
f (t)の式から時系列データを作り、f (t)の 10Hz のフーリエ係数の計算をするにあた
り、次のような数値を代入し計算してみましょう。
・ サンプル周波数 1000Hz、
・ サンプル時間h = 1/1000 (s)、
・ サンプル数 N = 2048
・ t = nh
・ n = 0、1、2、・・・2047
として
n、nh、f (nh)、cos (20πnh)、sin (20πnh)、f (nh) cos (20πnh)、f (nh) sin (20πnh)を
順に A ∼ G 列で計算します。
そして;
a (10Hz) =
b (10Hz) =
2
N
2
N
N −1
∑ f (nh )cos(20π nh ) ( F列の合計から計算する )
n =0
N −1
∑ f (nh )sin(20π nh ) ( G列の合計から計算する )
n =0
にて 10Hz のスペクトル a、b を N = 500(n = 0 ∼ 499)と N = 2048(n = 0 ∼ 2047)
の場合を求めると;
N = 500(Nh = 0.5s)の場合
a = 0.866
b = -0.500
N = 2048(Nh = 2.048s)の場合
a = 0.853
b = -0.508
-5-
Microsoft-EXCEL での計算の様子
f (t)の式から 10Hz のスペクトルは;
a = cos 30 deg = 0.866
b = sin 30 deg = 0.5
であることがわかっていますから、N = 500 のときは一致しますが、N = 2048 では
少し差があります。これはフーリエ級数の仮定である「周期性」によります。
N = 500(n = 0 ∼ 499 間で 0.5 s 間)の場合は 10Hz のちょうど 5 周期分になりま
すが、N = 2048 では 10Hz の整数倍になっていません。
-6-
次図 A のようにちょうど 1 周期(またはその整数倍)になるように N をとると、
同じ波形が連続で続いて描くことができますが、1周期(又はその整数倍)とな
らない N では不連続になります(図 B)。
不連続点があると、その中間点をとることで連続となりますから、中間点で滑ら
かに繋がった波形として近似計算がおこなわれます。
MUXO
<H5A + QN 0* DXODKL3J
0+ CY7@[VCBIF;2
A8?TIP=>G[VCBIF;2
R.
`]
`\
R/
`]
+ QN?EB4 0, DXODKL3J
0- CY7@Z[V@BI1Z[VW
6S:>9F4F;2
Z[V
Z[V
`_
Z[V
`^
Z[V
`_
-7-
さて、10Hz のフーリエ係数;
a = 0.866
b = -0.5
より、
f10 Hz (t ) = 0.866 cos (2π × 10t ) − 0.5 sin (2π × 10t )
となります。さらに、この式は又、Re = 0.866、Im = -b = 0.5 より;
C = a 2 + b 2 = 0.866 2 + 0.52 = 1
Im
0.5
π
(30 deg )
= tan −1
= Re
0.866 6
π⎞
⎛
f10 Hz (t ) = cos⎜ 2π × 10t + ⎟
6⎠
⎝
θ = tan −1
であることがわかります。
同様に、他の周波数成分を求めると、観測波形 f (t)は;
f (t ) =
∑f
n
となります。
これが、FFT アナライザの概念です。
-8-
(t )