テンソル伝達関数

アクロス解析入門
羽佐田葉子
2007年3月24日
アクロス研究会＠静岡大学
アクロスデータ解析必須項目
• 離散フーリエ変換
…デジタル信号処理
– 複素数の表現
• アクロス信号とは
• ノイズと誤差の統計学
• テンソル伝達関数
まだあるかも･･･
…確率・統計
…線形代数学
アクロスデータ解析の流れ（１）
時間領域データ
フーリエ変換
送信信号で除算
誤差つき伝達関数
周波数領域データ
正逆回転合成
信号成分
取り出し
ノイズレベル
推定
重みつきスタッキング
誤差つき周波数領域データ
座標変換
テンソル伝達関数
アクロスデータ解析の流れ（２）
テンソル伝達関数
周波数窓
フーリエ逆変換
時間領域波形
存否イベント解析
波の到着イベント
離散フーリエ変換
Discrete Fourier transform
時間領域
フーリエ変換
time domain
xn
周波数領域
frequency domain
フーリエ逆変換
Xk = Ck + i Sk
時間 t
周波数 f
xn = Sk { Ck cos (2p fk tn) + Sk sin (2p fk tn) }
時間領域データxnをcosとsinの重ね合わせで表現したときの
各周波数成分の係数を複素数の形で表したのがXk
複素数の表現
z = a + i b = r exp(i q) ←オイラーの公式
実数部分
虚数部分
絶対値
偏角
a = Re[z]
b = Im[z]
r=|z|
q = Arg[z]
↑Matlab表記
Im
b
real(z)
imag(z)
abs(z)
angle(z)
z
r
q
a
Re
離散フーリエ変換
xn = Sk { Ck cos (2p fk tn)
周波数 f
+ Sk sin (2p fk tn) }
実数部分がcosの係数
虚数部分がsinの係数
時間 t
時間 t
時間 t
時間 t
時間 t
データ長の中に整数周期入る
周波数の波で表現
…
それ以外の周波数は整数周期
入る波の合成で表される
離散フーリエ変換
時間領域 x ( t )
Dt
NDt
フーリエ変換
周波数領域 X ( f )
フーリエ変換
時間 t
フーリエ逆変換
周波数 f
Df = 1/(NDt)
Xk = S xn exp (-2pi fk tn) / N
フーリエ逆変換 xn = S Xk exp (2pi fk tn)
n = 1, …, N
k = 1, …, N
tn = (n-1)Dt
fk = (k-1)Df
1/Dt
X=fft(x)/N
x=ifft(X)*N
↑Matlab表記
Df = 1/(NDt) :基本周波数
fNyquist = 1/(2Dt) : ナイキスト周波数
ナイキスト周波数より高周波側は、低周波側の複素共役
X2 = X*N
X3 = X*N-1
X4 = X*N-2
…
いろいろなフーリエ変換対
cosの場合
2
1
0
0
-2
sinの場合
1
2
時間(s)
3
4
-1
2
1
0
0
-2
半端な周期
0
0
1
2
時間(s)
3
4
-1
2
1
0
0
-2
0
1
2
時間(s)
3
4
-1
0
2
4
6
周波数(Hz)
8
10
0
2
4
6
周波数(Hz)
8
10
0
2
4
6
周波数(Hz)
8
10
いろいろなフーリエ変換
インパルス
2
0.05
0
0
-2
0
1
2
時間(s)
3
4
-0.05
2
0.05
0
0
-2
0
1
2
時間(s)
3
4
-0.05
2
0.1
0
0
-2
0
1
2
時間(s)
3
4
-0.1
0
2
4
6
周波数(Hz)
8
10
0
2
4
6
周波数(Hz)
8
10
0
2
4
6
周波数(Hz)
8
10
アクロスデータ解析の流れ（１）
時間領域データ
フーリエ変換
送信信号で除算
今ココ
誤差つき伝達関数
周波数領域データ
正逆回転合成
信号成分
取り出し
ノイズレベル
推定
重みつきスタッキング
誤差つき周波数領域データ
座標変換
テンソル伝達関数
アクロス観測データの離散フーリエ
変換
例）100Hzサンプリング、200秒間のデータ
x1, … , x20000 Dt = 0.01sec
フーリエ変換すると、
X1, … , X20000 Df = 1/200 = 0.005Hz
ナイキスト周波数は 1/(2×0.01) = 50Hz
50Hzよりも高周波の成分は0～50Hzの
ところに現れる → aliasing （折り返し歪）
この場合、アクロス信号は0~50Hzの範囲で
0.005Hzの倍数の周波数でなければならな
い！
アクロス信号（弾性波FM回転の場
合）回転周波数
例）周波数を10~20Hzの間で
変調させる。変調周期は50秒
→アクロス信号は
搬送波周波数を基準として
1/50=0.02Hz間隔に現れる
振幅
0.02Hz
搬送波周波数
20Hz
10Hz
50sec
200sec
時間 (s)
振幅
フーリエ変換
ノイズ
チャンネル
0.005Hz
周波数 f
10 20
50
80 90
周波数 (Hz)
アクロスデータ解析の流れ（１）
時間領域データ
フーリエ変換
送信信号で除算
誤差つき伝達関数
周波数領域データ
正逆回転合成
信号成分
ノイズレベル
取り出し
推定
今ココ
次ココ
重みつきスタッキング
誤差つき周波数領域データ
座標変換
テンソル伝達関数
ノイズと誤差の統計学
• 正規白色雑音 xn~N(0,s2)
• 白色雑音のフーリエ変換
実数部分と虚数部分が
それぞれ独立に同じ分布に従う
Re[Xk]~N(0,s2/(2N)), Im[Xk]~N(0,s2/(2N))
標準偏差は s/(2N)1/2、データ長の-1/2乗に
比例 →スタッキングの効果
2
0.4
1
0.2
0
0
-1
-0.2
-2
0
1
2
時間(s)
3
4
-0.4
0
2
4
6
周波数(Hz)
8
10
ノイズと誤差の統計学
• アクロスデータに含まれる誤差の推定
– ノイズチャンネルからノイズレベルを推定
– ノイズを正規白色雑音と仮定すると、
ノイズチャンネルの振幅の２乗平均から
誤差の分散が推定できる
• 実際の自然のノイズレベルは周波数依存
– 適当な周波数範囲を区切って、その中では分
散が一定と仮定して推定
• 推定した誤差をスタッキングの重みに使用
アクロスデータ解析の流れ（１）
時間領域データ
送信信号で除算
フーリエ変換
誤差つき伝達関数
周波数領域データ
正逆回転合成
信号成分
取り出し
ノイズレベル
推定
座標変換
重みつきスタッキング
テンソル伝達関数
誤差つき周波数領域データ
今ココ
テンソル伝達関数
• 伝達関数 transfer function とは…
入力
X(w)
線形システム
出力
Y(w)
伝達関数 H(w) = Y(w)/X(w)
– 線形システムの特性を表す関数
– ある周波数の信号が、そのシステムによって
どのように変化するか
– 絶対値 | H(w) | が振幅の増幅率を、
偏角Arg[H(w)]が位相の進みを表す
テンソル伝達関数
• アクロスで求まる伝達関数
– 受信信号を送信信号で割り、地下の伝播特性
を求める
• 弾性波アクロスの場合
– 入力：震源の力ベクトル(N)
出力：地面の変位ベクトル(m) テンソル伝達関数
Ur
HrR HrT HrZ FR
Ut = HtR HtT HtZ FT
Uz
HzR HzT HzZ FZ
与えた力の向き、
大きさに対する
変位の向き、大きさが
分かる
正逆回転の合成、座標変換により求める
アクロスデータ解析の流れ（２）
テンソル伝達関数
今ココ
周波数窓
フーリエ逆変換
時間領域波形
存否イベント解析
波の到着イベント
伝達関数を時間領域で見る
伝達関数（周波数領域）
入力
X(w)
線形システム
出力
Y(w)
伝達関数 H(w) = Y(w)/X(w)
＝入力が全て１の場合の出力
フーリエ逆変換で時間領域に変換すると？
振幅が１の全ての周波数のcosを足し合わせる
と時刻0のインパルス
全ての周波数の出力を足し合わせると、時刻0
にインパルスを入力したときの出力
時刻0のインパルスがどんなふうに伝わるか？
伝達関数を時間領域で見る
• 時間領域で見ることで、
別の経路を通ってきた波を
分離できる
時間(s)
• 現実には全ての周波数を
観測できないので
インパルスにはならない
• なるべく広い周波数範囲を
反射波観測すればシャープなパルス
が見える
• 存否イベント解析を使うと
狭い周波数範囲でも
それなりにパルスを分離可能
震源
10km
20km
30km
P波
S波
アクロスデータ解析必須項目
• 離散フーリエ変換
…デジタル信号処理
– 複素数の表現
• アクロス信号とは
• ノイズと誤差の統計学
• テンソル伝達関数
まだあるかも･･･
…確率・統計
…線形代数学

Download Report