yProcessingClub

すみません、許してください

2信号のウィグナービレ分布を求めよ

信号処理

ウィグナービレ分布

2信号のウィグナービレ分布は？

yuri-processing-club.hatenablog.com
この記事では
信号が1つ(s(t)=\cos(2 \pi f_c t))のウィグナービレ分布を求め，
W_z(t,f) = \delta(f-f_c)
となった．

それでは2信号ではどうなるだろうか．

問題

s(t)=\cos(2 \pi f_1 t) + \cos(2 \pi f_2 t)のウィグナービレ分布を求めよ．


解法

まずs(t)を解析信号z(t)に変換する．
z(t)=\mathrm{e}^{j 2 \pi f_1 t} + \mathrm{e}^{j 2 \pi f_2 t}

また，
z^*(t)=\mathrm{e}^{-j 2 \pi f_1 t} + \mathrm{e}^{-j 2 \pi f_2 t}

次にSignal Kernelを求める．

z\left(t+\frac{\tau}{2}\right)z^*\left(t-\frac{\tau}{2}\right)\\ =\left( \mathrm{e}^{j 2 \pi f_1 \left( t +\frac{\tau}{2} \right)} + \mathrm{e}^{j 2 \pi f_2 \left( t +\frac{\tau}{2} \right)} \right) \left(\mathrm{e}^{-j 2 \pi f_1 \left( t -\frac{\tau}{2} \right)} + \mathrm{e}^{-j 2 \pi f_2 \left( t -\frac{\tau}{2} \right)}\right)\\ = \mathrm{e}^{j 2 \pi f_1 \tau} + \mathrm{e}^{j 2 \pi f_2 \tau} + 2 \cos \left( 2 \pi \left( f_2 - f_1 \right) t \right) \mathrm{e}^{j 2 \pi \left( \frac{f_1+f_2}{2} \right) \tau}


よって，ウィグナービレ分布は，

\begin{eqnarray*} W_z(t,f) &=& \int z\left(t+\frac{\tau}{2}\right)z^*\left(t-\frac{\tau}{2}\right) \mathrm{e}^{-j 2 \pi f \tau} d\tau\\ &=& \int \left( \mathrm{e}^{j 2 \pi f_1 \tau} + \mathrm{e}^{j 2 \pi f_2 \tau} + 2 \cos \left( 2 \pi \left( f_2 - f_1 \right) t \right) \mathrm{e}^{j 2 \pi \left( \frac{f_1+f_2}{2} \right) \tau} \right) \mathrm{e}^{-j 2 \pi f \tau} d\tau\\ &=& \int \mathrm{e}^{j 2 \pi f_1 \tau}\mathrm{e}^{-j 2 \pi f \tau}d\tau + \int \mathrm{e}^{j 2 \pi f_2 \tau}\mathrm{e}^{-j 2 \pi f \tau}d\tau\\ && +\int \left( 2 \cos \left( 2 \pi \left( f_2 - f_1 \right) t \right) \mathrm{e}^{j 2 \pi \left( \frac{f_1+f_2}{2} \right) \tau} \right) \mathrm{e}^{-j 2 \pi f \tau}d\tau\\ &=& \int \mathrm{e}^{j 2 \pi f_1 \tau}\mathrm{e}^{-j 2 \pi f \tau}d\tau + \int \mathrm{e}^{j 2 \pi f_2 \tau}\mathrm{e}^{-j 2 \pi f \tau}d\tau \\ &&+ 2 \cos \left( 2 \pi \left( f_2 - f_1 \right) t \right) \int \mathrm{e}^{j 2 \pi \left( \frac{f_1+f_2}{2} \right) \tau} \mathrm{e}^{-j 2 \pi f \tau}d\tau\\ &=& \delta(f-f_1) + \delta(f-f_2) + 2 \cos \left( 2 \pi \left( f_2 - f_1 \right) t \right) \delta \left( f - \frac{f_1+f_2}{2} \right) \end{eqnarray*}
となり計算終了である．


答え

s(t)=\cos(2 \pi f_1 t) + \cos(2 \pi f_2 t)のウィグナービレ分布は，
W_z(t,f) = \delta(f-f_1) + \delta(f-f_2) + 2 \cos \left( 2 \pi \left( f_2 - f_1 \right) t \right) \delta \left( f - \frac{f_1+f_2}{2} \right)



クロス項

\delta(f-f_1)及び\delta(f-f_2)は期待通りの項である．2 \cos \left( 2 \pi \left( f_2 - f_1 \right) t \right) \delta \left( f - \frac{f_1+f_2}{2} \right)は本来存在しないスペクトルであり，これをクロス項と言う．クロス項は2信号の中間の位置に出ると言われているが，計算結果を見ると，確かにf_1及びf_2の中間の周波数\frac{f_1+f_2}{2}の部分に発生していることが分かる．また，\cos \left( 2 \pi \left( f_2 - f_1 \right) t \right)を見ると分かるが，クロス項は時間方向に周波数f_2-f_1で振動している．


図解

f:id:Yuri-Processing-Club:20181115161108p:plain
ウィグナービレ分布を図にするとこんな感じだろう（多分）．

AM変調信号のウィグナービレ分布

信号処理

今回も引き続き情報の精度については自信が無い．
メモとしてブログに残しておく．

yuri-processing-club.hatenablog.com

この記事で，z(t)=A(t)\mathrm{e}^{j \phi (t)}のウィグナービレ分布は

W_z(t,f) = \mathcal{W}_A(t, f-f_i(t))
になることを計算で確かめた．

本記事ではこれの意味を考えていく．


f:id:Yuri-Processing-Club:20181114193536p:plain
周波数一定（f_0）の信号をf_cで変調しよう．



f:id:Yuri-Processing-Club:20181114193459p:plain
スペクトルで見ると，このようにf_cだけ平行移動する形となる．

この例では元の信号の周波数f_0も搬送波f_cも時間変化しない．




f:id:Yuri-Processing-Club:20181114193558p:plain
時間周波数平面で書き直してみるとこんな感じ．周波数の時間変化は無いのでグラフの傾きは0である．



次にA(t)の周波数が時間変化する信号を考える．
搬送波周波数は一定なので，
f:id:Yuri-Processing-Club:20181114193824p:plain
同じようにf_cだけ平行移動した形となる．



次にA(t)もf_cも変調する場合を考える．
f:id:Yuri-Processing-Club:20181114194006p:plain
ある時点t_1におけるもとの信号の周波数\pm f(t_1)はその時点での搬送波の周波数f_c(t_1)だけずれるので，変調後の周波数f'=f_c(t_1) \pm f(t_1)



f:id:Yuri-Processing-Club:20181114200201p:plain
時間-周波数平面W(t,f)で考えると，変調後のある点W_z(t_1,f_1)は，変調前はf_c(t_1)だけ下にあった．つまり変調前は\mathcal{W}_A(t_1,f_1-f_c(t_1))の位置にあった．

変調する搬送波周波数f_c(t)を瞬時周波数f_i(t)に書き直して，W_z(t,f) = \mathcal{W}_A(t, f-f_i(t))である．



振幅A(t)を瞬時周波数f_i(t)を使ってAM変調した信号のウィグナービレ分布W_zは，f_i(t)とA(t)のウィグナー分布 \mathcal{W}_Aを使って表せることが分かる．
若干煙に巻いた部分はあるが，何となくの解釈はこんな感じでよいかもしれないし，よくないかもしれない．

メモ：周波数シフト計算

信号処理


72ページより計算メモ．


ぜんぜんわからん・・・
とりあえず計算メモだけ残しておく．

\begin{eqnarray*} k_z(\nu,f)&=&\int W_z(t,f) \mathrm{e}^{-j 2 \pi \nu t} dt\\ &=&\int \left( \int z(t+\tau/2)z^*(t-\tau/2) \mathrm{e}^{-j 2 \pi f \tau} d\tau \right) \mathrm{e}^{-j 2 \pi \nu t} dt\\ &=& \int \int z(t+\tau/2)z^*(t-\tau/2)\mathrm{e}^{-j 2 \pi (f \tau + \nu t)} dt d\tau \end{eqnarray*}


u = t+\frac{\tau}{2}=U(t, \tau)，v=t-\frac{\tau}{2}=V(t,\tau)と置くと，

t=\frac{1}{2}(u+v)=T(u,v)，\tau=u-v=\mathcal{T}(u,v)

\begin{eqnarray*} dtd\tau &=& \mathrm{ABS} \left( \mathrm{det} \left( \begin{array}{cc} \partial T / \partial u & \partial T / \partial v \\ \partial \mathcal{T} / \partial u & \partial \mathcal{T} / \partial v \\ \end{array} \right) \right) dudv\\ &=& \mathrm{ABS} \left( \mathrm{det} \left( \begin{array}{cc} 1/2 & 1/2 \\ 1 & -1 \\ \end{array} \right) \right) dudv\\ &=& \mathrm{ABS} \left( -1/2 - 1/2 \right)dudv\\ &=& \mathrm{ABS} \left( -1 \right)dudv\\ &=& dudv\ \end{eqnarray*}
↑この辺かなり微妙である．ヤコビアン行列が云々らしいが・・・？

よって，

\begin{eqnarray*} k_z(\nu,f)&=& \int \int z(t+\tau/2)z^*(t-\tau/2)\mathrm{e}^{-j 2 \pi (f \tau + \nu t)} dt d\tau \\ &=& \int \int z(u)z^*(v)\mathrm{e}^{-j 2 \pi (f (u-v)+ \nu \frac{1}{2}(u+v))} dudv \\ &=& \int \int z(u)z^*(v)\mathrm{e}^{-j 2 \pi (fu-fv+ \frac{1}{2}\nu u + \frac{1}{2} \nu v)} dudv \\ &=& \int \int z(u)z^*(v)\mathrm{e}^{-j 2 \pi (f + \frac{1}{2}\nu)u } \mathrm{e}^{j 2 \pi (f - \frac{1}{2}\nu)v } dudv \\ &=& \int z(u)\mathrm{e}^{-j 2 \pi (f + \frac{1}{2}\nu)u } du \int z^*(v)\mathrm{e}^{j 2 \pi (f - \frac{1}{2}\nu)v } dv \\ &=& Z \left( f+\frac{\nu}{2} \right) Z^* \left( f-\frac{\nu}{2} \right) \end{eqnarray*}

よって，
W_z(t,f)=\int Z \left( f+\frac{\nu}{2} \right) Z^* \left( f-\frac{\nu}{2} \right) \mathrm{e}^{j 2 \pi \nu t} d\nu


\nuを周波数シフト（またはドップラー）と言うそうだが，何がどうなってるのかぜんぜんわからん

ウィグナー分布の求め方は2通りあって，
1.時間ドメインで相関関数的なやつ（instantaneous autocorrelation function, IAF）K_z(t,f)を求めてフーリエ変換
2.周波数ドメインで相関関数的なやつ（spectral correlation function, SCF）k_z(\nu,f)を求めて逆フーリエ変換（本記事の話）
みたいな感じかな・・・・・・？

ふわふわしたことを言うが，おそらく時間-周波数平面として(t - f)平面と(\tau - \nu)平面があるっぽくて，後者は曖昧度関数が持つ平面っぽい．

LFM信号のウィグナービレ分布

信号処理


71ページから解説．

yuri-processing-club.hatenablog.com
ここと併せて読まれたし．


上の記事では信号としてz(t)=\mathrm{e}^{j \phi (t)}を考えた．

今回は振幅が時間変化する信号z(t)=A(t)\mathrm{e}^{j \phi (t)}を考える．

Signal Kernelは
\begin{eqnarray*} K_z(t,\tau) &=& z(t+\tau/2)z^*(t-\tau/2)\\ &=& A(t+\tau/2)\mathrm{e}^{j \phi (t+\tau/2)} A(t-\tau/2)\mathrm{e}^{-j \phi (t-\tau/2)}\\ &=& A(t+\tau/2) A(t-\tau/2) \mathrm{e}^{j \phi (t+\tau/2)}\mathrm{e}^{-j \phi (t-\tau/2)}\\ &=& A(t+\tau/2) A(t-\tau/2) \mathrm{e}^{j (\phi (t+\tau/2) - \phi (t-\tau/2))}\\ &=& K_A(t,\tau)\mathrm{e}^{j\psi(t,\tau)} \end{eqnarray*}

ここで\phi(t)が二次式ならば，limの近似が行え，
\psi (t, \tau) = \phi ' (t) \tau = 2 \pi f_i (t) \tauである．

よって，
K_z(t,\tau) =K_A(t,\tau)\mathrm{e}^{j2 \pi f_i (t) \tau}

これをフーリエ変換したものがウィグナービレ分布であるから，
\begin{eqnarray*} W_z(t,f) &=& \int K_z(t,\tau) \mathrm{e}^{-j 2 \pi f \tau} d\tau\\ &=& \int K_A(t,\tau)\mathrm{e}^{j2 \pi f_i (t) \tau} \mathrm{e}^{-j 2 \pi f \tau} d\tau\\ &=& \int K_A(t,\tau)\mathrm{e}^{-j2 \pi (f - f_i (t)) \tau} d\tau\\ &=& \mathcal{W}_A(t, f-f_i(t)) \end{eqnarray*}

（区別のため，ウィグナービレ分布をW，ウィグナー分布を\mathcal{W}と書く）

Signal Kernelの計算方法

信号処理


66-67ページから解説．


ある時間信号z(t)の時間-周波数表現\rho_z(t,f)があるとして，\rho_z(t,f)を逆フーリエ変換したものをK_z(t,\tau)と書く．ここでK_z(t,\tau)をSignal Kernelと呼ぶことにする．K_z(t,\tau)をフーリエ変換すると\rho_z(t,f)になるわけだから，数式で書くと，
\rho_z(t,f) = \int K_z(t,\tau) \mathrm{e}^{-j 2 \pi f \tau} d\tau


ここでK_z(t,\tau)はどのような形になるか．

ある信号z(t)の例として，解析信号z(t)=\mathrm{e}^{j \phi (t)}を考える．

瞬時周波数f_i (t) = \frac{\phi ' (t)}{2 \pi}より，

\rho_z(t,f) = \delta (f - f_i(t))である．

例えば5秒時点において100Hzの部分にスペクトルがあるとすれば\rho_z(5,100) = 1となるわけである．

K_z(t,\tau)は\rho_z(t,f) = \delta (f - f_i(t))の逆フーリエ変換なので，
\begin{eqnarray*} K_z(t,\tau) &=& \int \rho_z(t,f) \mathrm{e}^{j 2 \pi f \tau} df\\ &=& \int \delta (f - f_i(t)) \mathrm{e}^{j 2 \pi f \tau} df\\ &=& \mathrm{e}^{j 2 \pi f_i(t) \tau}\\ &=& \mathrm{e}^{j \phi ' (t) \tau}\\ \end{eqnarray*}

ここの計算は \mathcal{F}^{-1} [ \delta (f - f_0) ] = \mathrm{e}^{j 2 \pi f_0 \tau}を使えば行ける．

さて，瞬時周波数\phi ' (t)だが，
\displaystyle \phi ' (t) = \lim _{\tau \to 0} \frac{ \phi(t + \tau / 2) - \phi(t - \tau / 2)}{\tau}
と書ける．これは中心有限差分などで検索すると色々出てくるが，ある点 \phi(t)における傾き（ =\phi ' (t)）は， \phi(t)の前後の2点 \phi(t + \tau / 2)及び\phi(t - \tau / 2)を通る直線の傾きに等しい（\tau \to 0において）と言える．ここで\phi ' (t) が直線ならば、言い換えると\phi ' (t) が一次式（または定数）ならば，さらに言い換えると \bf {\phi (t)}が二次式（または一次式）ならば，\tau \to 0でなくても傾きは等しくなる．つまり，

\phi ' (t) \approx \frac{ \phi(t + \tau / 2) - \phi(t - \tau / 2)}{\tau} \ \ (\mathrm{if} \ \ \phi (t) \mathrm{\ \ is \ \ quadratic \ \ or \ \ linear})

これを代入すると，
\begin{eqnarray*} K_z(t,\tau) &=& \mathrm{e}^{j \phi ' (t) \tau}\\ &=& \mathrm{e}^{j(\phi(t + \tau / 2) - \phi(t - \tau / 2))}\\ &=& \mathrm{e}^{j\phi(t + \tau / 2)}\mathrm{e}^{-j\phi(t + \tau / 2)}\\ &=& z(t+\tau / 2)z^*(t-\tau / 2) \end{eqnarray*}


よって，
\begin{eqnarray*} \rho_z(t,f) &=& \int K_z(t,\tau) \mathrm{e}^{-j 2 \pi f \tau} d\tau\\ &=& \int z(t+\tau / 2)z^*(t-\tau / 2) \mathrm{e}^{-j 2 \pi f \tau} d\tau\\ \end{eqnarray*}

この形はなんじゃらほいと考えると，まさしくウィグナービレ分布なわけである．

くどいようだが，この話が正しくなる条件として\phi ' (t) が直線である必要がある．そうでなければ極限の近似が上手くいかないからである．

s(t)=cos(2πft)のウィグナービレ分布を求めよ

信号処理

ウィグナービレ分布とは？

ウィグナービレ分布とは時間周波数解析法の一種である．「解析信号z(t)=\mathcal{A}(s(t))のウィグナー分布」をウィグナービレ分布といい，ウィグナー分布と比較してクロス項の抑制が期待できる．

問題

s(t)=\cos(2 \pi f_c t)のウィグナービレ分布を求めよ．


解法

まず解析信号z(t)を求める．
オイラーの公式より，
s(t)=\cos(2 \pi f_c t) = \frac{1}{2} \left( \mathrm{e}^{j 2 \pi f t} + \mathrm{e}^{-j 2 \pi f t} \right)

また，\mathrm{e}^{j 2 \pi ft}の解析信号は
\mathcal{A}(\mathrm{e}^{j 2 \pi ft})= \begin{cases} 0 \ \ (\mathrm{if} \ \ f < 0)\\ 2 \mathrm{e}^{j 2 \pi ft} \ \ (\mathrm{if} \ \ f > 0)\\ \end{cases}
であるので，

\begin{eqnarray*} z(t)&=&\mathcal{A}(s)=\frac{1}{2} \mathcal{A}(\mathrm{e}^{j 2 \pi f_c t})+\frac{1}{2} \mathcal{A}(\mathrm{e}^{-j 2 \pi f_c t})=\frac{1}{2} 2 \mathrm{e}^{j 2 \pi f_c t}+\frac{1}{2}0\\ &=&\mathrm{e}^{j 2 \pi f_c t} \end{eqnarray*}


また，z^*(t)=\mathrm{e}^{-j 2 \pi f_c t}．
よって，
\begin{eqnarray*} z\left(t+\frac{\tau}{2}\right)z^*\left(t-\frac{\tau}{2}\right) &=& \mathrm{e}^{j 2 \pi f_c \left( t + \frac{\tau}{2} \right)} \mathrm{e}^{-j 2 \pi f_c \left( t - \frac{\tau}{2} \right)}\\ &=& \mathrm{e}^{j 2 \pi f_c t} \mathrm{e}^{j 2 \pi f_c \frac{\tau}{2}} \mathrm{e}^{-j 2 \pi f_c t} \mathrm{e}^{j 2 \pi f_c \frac{\tau}{2}}\\ &=& \mathrm{e}^{j 2 \pi f_c \tau} \end{eqnarray*}



次にz(t)のウィグナー分布（＝ウィグナービレ分布）W_z(t, f)を計算する．
\begin{eqnarray*} W_z(t, f)&=&\int z\left(t+\frac{\tau}{2}\right)z^*\left(t-\frac{\tau}{2}\right) \mathrm{e}^{-j 2 \pi f \tau}d\tau\\ &=& \int \mathrm{e}^{j 2 \pi f_c \tau} \mathrm{e}^{-j 2 \pi f \tau}d\tau\\ &=& \delta(f-f_c) \end{eqnarray*}
となり，計算終了である．


答え

s(t)=\cos(2 \pi f_c t)のウィグナービレ分布W_z(t,f)は，
W_z(t,f) = \delta(f-f_c)


クロス項抑制

解析信号は周波数f_cの信号で負の周波数を持たないものであり，そのウィグナービレ分布はW_z(t,f) = \delta(f-f_c)と期待通りのスペクトルを持つ．ウィグナー分布ではf=\pm f_cの中間のf=0の部分にクロス項が出ていたが，ウィグナー分布は信号を解析信号に変換することで負のスペクトルを消したため，2信号の中間に出るクロス項の発生を抑制することができた．



yuri-processing-club.hatenablog.com

s(t)=cos(2πft)のウィグナー分布を求めよ

信号処理

ウィグナー分布とは？

ウィグナー分布とは時間周波数解析法の一種である．

問題

s(t)=\cos(2 \pi f_c t)のウィグナー分布を求めよ．

解法

ウィグナー分布は以下で定義される．
W_s(t,f)=\int s\left(t+\frac{\tau}{2}\right)s^*\left(t-\frac{\tau}{2}\right) \mathrm{e}^{-j 2 \pi f \tau}d\tau


s(t)は実信号なので，s^*(t)=s(t)．
よって，
s\left(t+\frac{\tau}{2}\right)s^*\left(t-\frac{\tau}{2}\right) = s\left(t+\frac{\tau}{2}\right)s\left(t-\frac{\tau}{2}\right)

s(t)=\cos(2 \pi f_c t)を代入して，
s\left(t+\frac{\tau}{2}\right)s\left(t-\frac{\tau}{2}\right) = \cos 2 \pi f_c \left(t+\frac{\tau}{2}\right) \cos 2 \pi f_c \left(t-\frac{\tau}{2}\right)


ここで\alpha = 2 \pi f_c \left( t + \frac{\tau}{2} \right)，\beta = 2 \pi f_c \left( t - \frac{\tau}{2} \right)とおくと，
\alpha - \beta = 2 \pi f_c \tau，\alpha + \beta = 2 \pi f_c 2t
加法定理： \cos \alpha \cos \beta = \frac{\cos(\alpha - \beta) + \cos(\alpha + \beta)}{2}を用いて，

\cos 2 \pi f_c \left(t+\frac{\tau}{2}\right) \cos 2 \pi f_c \left(t-\frac{\tau}{2}\right) = \frac{1}{2}cos ( 2 \pi f_c \tau ) + \frac{1}{2}cos ( 2 \pi f_c 2t )


これをW_s(t,f)=\int s\left(t+\frac{\tau}{2}\right)s^*\left(t-\frac{\tau}{2}\right)\mathrm{e}^{-j 2 \pi f \tau}d\tauに代入して，

\begin{eqnarray*} W_s(t,f)&=&\int {\left( \frac{1}{2}\cos ( 2 \pi f_c \tau ) + \frac{1}{2}\cos ( 2 \pi f_c 2t ) \right) \mathrm{e}^{-j 2 \pi f \tau}d\tau}\\ &=& \frac{1}{2} \int \cos ( 2 \pi f_c \tau )\mathrm{e}^{-j 2 \pi f \tau}d\tau + \frac{1}{2}\cos ( 2 \pi f_c 2t ) \int \mathrm{e}^{-j 2 \pi f \tau}d\tau \end{eqnarray*}

ここで\cos(2 \pi f_c t)のフーリエ変換は \frac{\delta(f-f_c) + \delta(f+f_c)}{2}
1のフーリエ変換は\delta(f)より，

\begin{eqnarray*} W_s(t,f)&=& \frac{1}{2}\frac{\delta(f-f_c) + \delta(f+f_c)}{2} + \frac{1}{2}\cos ( 2 \pi f_c 2t )\delta(f)\\ &=& \frac{1}{4}\delta(f-f_c) + \frac{1}{4}\delta(f+f_c) + \frac{1}{2}\cos ( 2 \pi f_c 2t )\delta(f) \end{eqnarray*} 　
となり，計算終了である．


答え

s(t)=\cos(2 \pi f_c t)のウィグナー分布W_s(t,f)は，
W_s(t,f) = \frac{1}{4}\delta(f-f_c) + \frac{1}{4}\delta(f+f_c) + \frac{1}{2}\cos ( 2 \pi f_c 2t )\delta(f)


クロス項

周波数f_cの信号であるので，両側スペクトルを考えると，f=\pm f_cにスペクトルがあるべきで，\frac{1}{4}\delta(f-f_c)及び\frac{1}{4}\delta(f+f_c)の項は期待通りである．しかし，\frac{1}{2}\cos ( 2 \pi f_c 2t )\delta(f)という謎の項が出ている．これをクロス項といい，アーティファクトである．クロス項の発生個所を見ると，f=0の位置，つまりf=f_c及びf=-f_cの中間の部分にスペクトルが出ていることが分かる．このように，クロス項は信号と信号の中間の位置に出る．



yuri-processing-club.hatenablog.com