我们今天继续学习短时傅里叶变换,看清它的特性和不足。今天要引入新的工具来弥补短时傅里叶变换的不足,尽管它自身也有不足之处(没有缺点的工具是不存在的)。
短时傅里叶变换(STFT / Windowed FT)
引入窗函数
g
(
t
)
g(t)
g(t) 做时间上的局部化:
g
(
t
)
,
Localization (Time)
V
g
f
(
x
,
w
)
=
∫
−
∞
+
∞
f
(
t
)
g
(
t
−
x
)
‾
exp
(
−
j
w
t
)
d
t
g(t) ,\ \text{Localization (Time)}\\ V_g f(x,w)=\int_{-\infty}^{+\infty}f(t)\overline{g(t-x)}\exp(-jwt)dt\\
g(t), Localization (Time)Vgf(x,w)=∫−∞+∞f(t)g(t−x)exp(−jwt)dt
- 我们非常关心短时傅里叶变换是否能够恰如其分地反映出信号在时域上的平移和频域上的调制:
f ( t ) → M w ′ T x ′ f ( t ) = f ( t − x ′ ) exp ( j w ′ t ) f(t)\to M_{w'}T_{x'}f(t)=f(t-x')\exp(jw't) f(t)→Mw′Tx′f(t)=f(t−x′)exp(jw′t)
根据短时傅里叶变换的定义:
V
g
(
M
w
′
T
x
′
f
)
(
x
,
w
)
=
∫
−
∞
+
∞
exp
(
j
w
′
t
)
f
(
t
−
x
′
)
g
(
t
−
x
)
‾
exp
(
−
j
w
t
)
d
t
=
t
′
=
t
−
x
′
∫
−
∞
+
∞
f
(
t
′
)
g
(
t
′
+
x
′
−
x
)
‾
exp
(
−
j
(
w
−
w
′
)
(
t
′
+
x
′
)
)
d
t
′
=
exp
(
−
j
(
w
−
w
′
)
x
′
)
∫
−
∞
+
∞
f
(
t
′
)
g
(
t
′
−
(
x
−
x
′
)
)
‾
exp
(
−
j
(
w
−
w
′
)
t
′
)
d
t
′
=
exp
(
−
j
(
w
−
w
′
)
x
′
)
V
g
f
(
x
−
x
′
,
w
−
w
′
)
V_g (M_{w'}T_{x'}f)(x,w)=\int_{-\infty}^{+\infty}\exp(jw't)f(t-x')\overline{g(t-x)}\exp(-jwt)dt\\ \xlongequal{t'=t-x'}\int_{-\infty}^{+\infty}f(t')\overline{g(t'+x'-x)}\exp(-j(w-w')(t'+x'))dt'\\ =\exp(-j(w-w')x')\int_{-\infty}^{+\infty}f(t')\overline{g(t'-(x-x'))}\exp(-j(w-w')t')dt'\\ =\exp(-j(w-w')x')V_g f(x-x',w-w')
Vg(Mw′Tx′f)(x,w)=∫−∞+∞exp(jw′t)f(t−x′)g(t−x)exp(−jwt)dtt′=t−x′∫−∞+∞f(t′)g(t′+x′−x)exp(−j(w−w′)(t′+x′))dt′=exp(−j(w−w′)x′)∫−∞+∞f(t′)g(t′−(x−x′))exp(−j(w−w′)t′)dt′=exp(−j(w−w′)x′)Vgf(x−x′,w−w′)
可见,短时傅里叶变换确实能够做到这一点。只不过相差了一个指数因子,这是一个小缺点。(后面的分析工具能弥补这个缺点)
如果一个信号有自己的傅里叶变换,我们通过Parseval关系知道傅里叶变换是“保能量”的:(注意, 这里及后面都抛开
π
2
\frac{\pi}{2}
2π不谈)
f
→
F
f
^
,
∣
∣
f
∣
∣
2
=
∣
∣
f
^
∣
∣
2
f\stackrel{\mathcal F}{\to}\widehat f, \quad ||f||^2=||\widehat f||^2
f→Ff
,∣∣f∣∣2=∣∣f
∣∣2
这说明傅里叶变换本身只是一个正交变换,是一个旋转。这也是对的,举个例子,时域频域就像同一块宝石的不同面,从不同的角度看,它的重量并没有变。
- 现在我们比较关心短时傅里叶变换是否也能够保能量?虽然我们很容易注意到,“加窗”这一操作对能量的损耗应该是很大的。
Parseval关系有一个简单的推广:
<
f
,
g
>
=
<
f
^
,
g
^
>
<f,g>=<\hat f, \widehat g>
<f,g>=<f^,g
>
我们关心Parseval关系是否能够推广到短时傅里叶变换,并且我们猜测有如下的形式:
<
V
g
1
f
1
,
V
g
2
f
2
>
=
<
f
1
,
f
2
>
<
g
1
,
g
2
>
‾
<V_{g_1}f_1,V_{g_2}f_2>=<f_1, f_2>\overline{<g_1,g_2>}
<Vg1f1,Vg2f2>=<f1,f2><g1,g2>
下面我们来推导这件事,为了减轻积分的难度,我们可以把短时傅里叶变换看做某种傅里叶变换,然后利用Parseval关系:
<
V
g
1
f
1
,
V
g
2
f
2
>
=
<
(
f
1
(
t
)
g
1
(
t
−
x
)
‾
)
∧
,
(
f
2
(
t
)
g
2
(
t
−
x
)
‾
)
∧
>
(
x
,
w
)
双重积分
=
<
(
f
1
(
t
)
g
1
(
t
−
x
)
‾
)
,
(
f
2
(
t
)
g
2
(
t
−
x
)
‾
)
>
(
x
,
t
)
=
∫
−
∞
+
∞
∫
−
∞
+
∞
f
1
(
t
)
g
1
(
t
−
x
)
‾
)
f
2
(
t
)
g
2
(
t
−
x
)
‾
‾
d
t
d
x
=
∫
−
∞
+
∞
f
1
(
t
)
f
2
(
t
)
‾
(
∫
−
∞
+
∞
g
1
(
t
−
x
)
‾
)
g
2
(
t
−
x
)
d
x
)
d
t
=
x
′
=
t
−
x
∫
−
∞
+
∞
f
1
(
t
)
f
2
(
t
)
‾
d
t
∫
−
∞
+
∞
g
1
(
x
)
‾
)
g
2
(
x
)
d
x
=
<
f
1
,
f
2
>
<
g
1
,
g
2
>
‾
<V_{g_1}f_1,V_{g_2}f_2>=<(f_1(t)\overline{g_1(t-x)})^{\wedge},(f_2(t)\overline{g_2(t-x)})^{\wedge}>_{(x,w)\text{双重积分}}\\ =<(f_1(t)\overline{g_1(t-x)}),(f_2(t)\overline{g_2(t-x)})>_{(x,t)}\\ =\int_{-\infty}^{+\infty}\int_{-\infty}^{+\infty}f_1(t)\overline{g_1(t-x)})\overline{f_2(t)\overline{g_2(t-x)}} dtdx\\ =\int_{-\infty}^{+\infty}f_1(t)\overline{f_2(t)}\left( \int_{-\infty}^{+\infty}\overline{g_1(t-x)})g_2(t-x) dx\right)dt\\ \xlongequal{x'=t-x}\int_{-\infty}^{+\infty}f_1(t)\overline{f_2(t)}dt\int_{-\infty}^{+\infty}\overline{g_1(x)})g_2(x) dx\\ =<f_1, f_2>\overline{<g_1,g_2>}
<Vg1f1,Vg2f2>=<(f1(t)g1(t−x))∧,(f2(t)g2(t−x))∧>(x,w)双重积分=<(f1(t)g1(t−x)),(f2(t)g2(t−x))>(x,t)=∫−∞+∞∫−∞+∞f1(t)g1(t−x))f2(t)g2(t−x)dtdx=∫−∞+∞f1(t)f2(t)(∫−∞+∞g1(t−x))g2(t−x)dx)dtx′=t−x∫−∞+∞f1(t)f2(t)dt∫−∞+∞g1(x))g2(x)dx=<f1,f2><g1,g2>
这个结果叫做Moyal公式。是Parseval的一个推广。它从能量的角度对短时傅里叶变换进行了说明:
∣
∣
V
g
f
∣
∣
2
=
<
V
g
f
,
V
g
f
>
=
<
f
,
f
>
<
g
,
g
>
‾
=
∣
∣
f
∣
∣
2
∣
∣
g
∣
∣
2
||V_gf||^2=<V_{g}f_,\ V_{g}f>=<f,f>\overline{<g,g>}=||f||^2||g||^2
∣∣Vgf∣∣2=<Vgf, Vgf>=<f,f><g,g>=∣∣f∣∣2∣∣g∣∣2
因此,当人们在选择这个窗
g
g
g 的时候,有意地将这个能量归一化:
∣
∣
g
∣
∣
2
=
1
||g||^2=1
∣∣g∣∣2=1,从而和傅里叶变换本身能够比较。
- 另外Moyal公式还有一个用途,就是做短时傅里叶变化的逆变换(Inversion)
f ( t ) = C ∫ R 2 V g f ( x , w ) ( M w T x h ( t ) ) d w d x f(t)=C\int_{\R^2} V_gf(x,w)(M_wT_xh(t))dwdx f(t)=C∫R2Vgf(x,w)(MwTxh(t))dwdx
这个
h
h
h 的加入就像一个药引子。
z
(
t
)
z(t)
z(t) 为任意关于t的函数
KaTeX parse error: Invalid delimiter type 'ordgroup' at position 140: …V_gf(x,w)\biggr{̲(̲}̲\int_{-\infty}^…
于是我们就得到了短时傅里叶变换的逆转公式:
∀
h
.
V
g
f
(
x
,
w
)
↔
1
<
g
,
h
>
‾
∫
R
2
V
g
f
(
x
,
w
)
(
M
w
T
x
h
(
t
)
)
d
w
d
x
=
f
(
t
)
\forall h.\quad V_gf(x,w)\leftrightarrow \frac{1}{\overline{<g,h>}}\int_{\R^2} V_gf(x,w)(M_wT_xh(t))dwdx=f(t)
∀h.Vgf(x,w)↔<g,h>1∫R2Vgf(x,w)(MwTxh(t))dwdx=f(t)
信号处理三部曲:
- Analysis(Transform、Decomposition, …)把复杂问题变简单。
- Process(Filtering,Operating, …)优化
- Synthesizing(Reconstruct,Anti-Transform, …)Restoration
有一个非常高的观点是Gabor提出的。Gabor认为通信就是把信号的能量,在不同的时间,不同的频率,换句话讲,在时频平面上做适当的安排。
Gabor发明了全息技术,因此而得诺贝尔物理学奖。全息照片是立体的。
我们人眼之所以能感知到立体,是因为不同物体散射的光进入我们人眼,眼睛所接受到的相位信息是不一样的。并不是因为两只眼睛我们才能感受到立体,你闭上一只眼睛同样能感受到立体。而照片我们之所以感受不到立体,是因为它把相位信息全部抹除了,只留下了强度信息。灰度不就是强度么。
全息技术,就是不但保留强度信息,更是通过技术手段保留相位信息,并通过特殊的播放装置,把相位信息播放出来。
我们继续,即将介绍下一个时频分析工具。我们来看一个信号:
k
(
x
,
w
)
=
δ
(
w
−
w
0
)
(Stationary)
⇒
k
(
x
,
w
)
=
δ
(
w
−
w
0
(
x
)
)
(Non-Stationary)
k(x,w)=\delta(w-w_0)\text{ (Stationary)}\\ \Rightarrow k(x,w)=\delta(w-w_0(x))\text{ (Non-Stationary)}
k(x,w)=δ(w−w0) (Stationary)⇒k(x,w)=δ(w−w0(x)) (Non-Stationary)
前者是平稳的情况,不随时间x发生变化,后者是非平稳的情况。我们现在主要分析非平稳的情况。
k
(
x
,
w
)
k(x,w)
k(x,w)是信号真实的表现(ground truth)。我们现在对
k
(
x
,
w
)
k(x,w)
k(x,w) 做一个反傅里叶变换。(注意,这里不在乎
1
2
π
\frac{1}{2\pi}
2π1)
∫
−
∞
+
∞
k
(
x
,
w
)
exp
(
j
w
t
)
d
w
=
∫
−
∞
+
∞
δ
(
w
−
w
0
(
x
)
)
exp
(
j
w
t
)
d
w
=
exp
(
j
w
0
(
x
)
t
)
\int_{-\infty}^{+\infty}k(x,w)\exp(jwt)dw=\int_{-\infty}^{+\infty}\delta(w-w_0(x))\exp(jwt)dw =\exp(jw_0(x)t)
∫−∞+∞k(x,w)exp(jwt)dw=∫−∞+∞δ(w−w0(x))exp(jwt)dw=exp(jw0(x)t)
注意,这里得到的
exp
(
j
w
0
(
x
)
t
)
\exp(jw_0(x)t)
exp(jw0(x)t) 是我的 ground-truth 经过傅里叶反变换得到的。这意味着我要分析的工具应当也长成这么个模样。这样就可以直接傅里叶变换得到我的 ground-truth 然而这样的分析工具只存在于理想中没法用,因为你并不知道
w
0
(
x
)
w_0(x)
w0(x)。因此我们需要跳出来看这个信号到底是什么?
- 如果 w 0 ( x ) w_0(x) w0(x) 是作为频域随时间x的变化,那么它指的就应该是**频率(Frequency)**了,只不过这个频率是时变的,此其一。
- 其二,对于信号而言,频率是什么意思?——相位求导就是频率: exp ( j w 0 x ) ⟹ w 0 x = ϕ ( x ) ϕ ′ ( x ) = w 0 \exp(jw_0x)\stackrel{w_0x=\phi(x)}{\Longrightarrow}\phi'(x)=w_0 exp(jw0x)⟹w0x=ϕ(x)ϕ′(x)=w0
所以,我们就可以把信号
f
(
x
)
=
exp
(
j
w
0
(
x
)
t
)
f(x)=\exp(jw_0(x)t)
f(x)=exp(jw0(x)t) 变成如下模样:
f
(
x
)
=
exp
(
j
ϕ
(
x
)
)
=
exp
(
j
ϕ
′
(
x
)
t
)
f(x)=\exp(j\phi(x))=\exp(j\phi'(x)t)
f(x)=exp(jϕ(x))=exp(jϕ′(x)t)
再使用一个差分近似来得到我想要的目标:
exp
(
j
ϕ
′
(
x
)
t
)
≈
差分近似
exp
(
j
ϕ
(
x
+
t
2
)
−
ϕ
(
x
−
t
2
)
t
t
)
=
exp
(
j
ϕ
(
x
+
t
2
)
)
exp
(
j
ϕ
(
x
−
t
2
)
)
‾
=
f
(
x
+
t
2
)
f
(
x
−
t
2
)
‾
\exp(j\phi'(x)t)\stackrel{\text{差分近似}}\approx\exp(j\frac{\phi(x+\frac{t}{2})-\phi(x-\frac{t}{2})}{t}t)\\ =\exp(j\phi(x+\frac{t}{2}))\overline{\exp(j\phi(x-\frac{t}{2}))}=f(x+\frac{t}{2})\overline{f(x-\frac{t}{2})}
exp(jϕ′(x)t)≈差分近似exp(jtϕ(x+2t)−ϕ(x−2t)t)=exp(jϕ(x+2t))exp(jϕ(x−2t))=f(x+2t)f(x−2t)
如此一来,我只要把基于单信号的某种操作
f
(
x
+
t
2
)
f
(
x
−
t
2
)
‾
f(x+\frac{t}{2})\overline{f(x-\frac{t}{2})}
f(x+2t)f(x−2t) 傅里叶变换过去,就能得到我ground truth的某种近似。
∫
−
∞
+
∞
f
(
x
+
t
2
)
f
(
x
−
t
2
)
‾
exp
(
−
j
w
t
)
d
t
\int_{-\infty}^{+\infty}f(x+\frac{t}{2})\overline{f(x-\frac{t}{2})}\exp(-jwt)dt
∫−∞+∞f(x+2t)f(x−2t)exp(−jwt)dt
这个新的时频分析工具记作
W
f
(
x
,
w
)
W_f(x,w)
Wf(x,w)(纪念 Wingner,这个变换统称Wingner-Ville Transform)和我们短时傅里叶变换的依赖性一样,他给出的是时频平面上的结果。但是它的做法跟短时傅里叶变换完全不同。是因为它是**二次的(Quadratic)**所以他没有线性特性,即
W
f
+
g
(
x
,
w
)
=
W
f
(
x
,
w
)
+
W
g
(
x
,
w
)
+
{Cross term}
W_{f+g}(x,w)= W_f(x,w)+W_g(x,w) +\text{\{Cross term\}}
Wf+g(x,w)=Wf(x,w)+Wg(x,w)+{Cross term},这也是它为人诟病的一点。但它也有很多优点:
-
Real 是实的。(如果把信号看做能量,那就应该是实的)
W f ( x , w ) ‾ = ∫ − ∞ + ∞ f ( x + t 2 ) ‾ f ( x − t 2 ) exp ( − j w t ) d t = t ′ = − t W f ( x , w ) \overline{W_f(x,w)}=\int_{-\infty}^{+\infty}\overline{f(x+\frac{t}{2})}f(x-\frac{t}{2})\exp(-jwt)dt\xlongequal{t'=-t}W_f(x,w) Wf(x,w)=∫−∞+∞f(x+2t)f(x−2t)exp(−jwt)dtt′=−tWf(x,w) -
W M w ′ T x ′ f ( x , w ) W_{M_{w'}T_{x'}f}(x,w) WMw′Tx′f(x,w)能够处理的相较于短时傅里叶变换更加干净:
W M w ′ T x ′ f ( x , w ) = ∫ − ∞ + ∞ ( M w ′ T x ′ f ) ( x + t 2 ) ( M w ′ T x ′ f ) ( x − t 2 ) ‾ exp ( − j w t ) d t = ∫ − ∞ + ∞ exp ( j w ′ ( x + t 2 ) ) f ( x − x ′ + t 2 ) exp ( − j w ′ ( x − t 2 ) ) f ( x − x ′ − t 2 ) ‾ exp ( − j w t ) d t = ∫ − ∞ + ∞ exp ( j w ′ t ) ) f ( x − x ′ + t 2 ) f ( x − x ′ − t 2 ) ‾ exp ( − j w t ) d t = ∫ − ∞ + ∞ f ( x − x ′ + t 2 ) f ( x − x ′ − t 2 ) ‾ exp ( − j ( w − w ′ ) t ) d t = W f ( x − x ′ , w − w ′ ) W_{M_{w'}T_{x'}f}(x,w)=\int_{-\infty}^{+\infty}(M_{w'}T_{x'}f)(x+\frac{t}{2})\overline{(M_{w'}T_{x'}f)(x-\frac{t}{2})}\exp(-jwt)dt\\ =\int_{-\infty}^{+\infty}\exp(jw'(x+\frac{t}{2}))f(x-x'+\frac{t}{2})\exp(-jw'(x-\frac{t}{2}))\overline{f(x-x'-\frac{t}{2})}\exp(-jwt)dt\\ =\int_{-\infty}^{+\infty}\exp(jw't))f(x-x'+\frac{t}{2})\overline{f(x-x'-\frac{t}{2})}\exp(-jwt)dt\\ =\int_{-\infty}^{+\infty}f(x-x'+\frac{t}{2})\overline{f(x-x'-\frac{t}{2})}\exp(-j(w-w')t)dt\\ =W_f(x-x',w-w') WMw′Tx′f(x,w)=∫−∞+∞(Mw′Tx′f)(x+2t)(Mw′Tx′f)(x−2t)exp(−jwt)dt=∫−∞+∞exp(jw′(x+2t))f(x−x′+2t)exp(−jw′(x−2t))f(x−x′−2t)exp(−jwt)dt=∫−∞+∞exp(jw′t))f(x−x′+2t)f(x−x′−2t)exp(−jwt)dt=∫−∞+∞f(x−x′+2t)f(x−x′−2t)exp(−j(w−w′)t)dt=Wf(x−x′,w−w′) -
V g ^ f ^ ( x , w ) = V g f ( − w , x ) exp ( . . . ) V_{\hat g}\hat f(x,w)=V_g f(-w,x)\exp(...) Vg^f^(x,w)=Vgf(−w,x)exp(...) 同样如此,能够去除相位因子。我们先对Wagner-Ville做一个推广:
W f , g ( x , w ) = ∫ − ∞ + ∞ f ( x + t 2 ) g ( x − t 2 ) ‾ exp ( − j w t ) d t W_{f,g}(x,w)=\int_{-\infty}^{+\infty}f(x+\frac{t}{2})\overline{g(x-\frac{t}{2})}\exp(-jwt)dt Wf,g(x,w)=∫−∞+∞f(x+2t)g(x−2t)exp(−jwt)dt
能够看到,他跟短时傅里叶变换非常相像,这一点点区别却造成了本质不同。
W f ^ , g ^ ( x , w ) = ∫ − ∞ + ∞ f ^ ( x + t 2 ) g ^ ( x − t 2 ) ‾ exp ( − j w t ) d t = h ^ 2 ( t ) = g ^ ( x − t 2 ) exp ( j w t ) h ^ 1 ( t ) = f ^ ( x + t 2 ) ∫ − ∞ + ∞ h ^ 1 ( t ) h ^ 2 ( t ) ‾ d t = Parseval关系 ∫ − ∞ + ∞ h 1 ( s ) h 2 ( s ) ‾ d s = 证明略 W f , g ( − w , x ) W_{\widehat f,\widehat g}(x,w)=\int_{-\infty}^{+\infty}\widehat f(x+\frac{t}{2})\overline{\widehat g(x-\frac{t}{2})}\exp(-jwt)dt\\ \xlongequal[\widehat h_2(t)=\widehat g(x-\frac{t}{2})\exp(jwt)]{\widehat h_1(t)=\widehat f(x+\frac{t}{2})}\int_{-\infty}^{+\infty}\widehat h_1(t)\overline{\widehat h_2(t)}dt\\ \xlongequal{\text{Parseval关系}}\int_{-\infty}^{+\infty}h_1(s)\overline {h_2(s)}ds\xlongequal{\text{证明略}}W_{f,g}(-w,x) Wf ,g (x,w)=∫−∞+∞f (x+2t)g (x−2t)exp(−jwt)dth 1(t)=f (x+2t)h 2(t)=g (x−2t)exp(jwt)∫−∞+∞h 1(t)h 2(t)dtParseval关系∫−∞+∞h1(s)h2(s)ds证明略Wf,g(−w,x)- Inversion Formula
f ( x + t 2 ) f ( x − t 2 ) ‾ = ∫ − ∞ + ∞ W f ( x , w ) exp ( j w t ) d t = x = t 2 f ( t ) f ( 0 ) ‾ ⟹ f ( t ) = 1 f ( 0 ) ‾ ∫ − ∞ + ∞ W f ( t 2 , w ) exp ( j w t ) d t f(x+\frac{t}{2})\overline{f(x-\frac{t}{2})}=\int_{-\infty}^{+\infty}W_f(x,w)\exp(jwt)dt\xlongequal{x=\frac{t}{2}}f(t)\overline{f(0)}\\ \Longrightarrow f(t)=\frac{1}{\overline{f(0)}}\int_{-\infty}^{+\infty}W_f(\frac{t}{2},w)\exp(jwt)dt f(x+2t)f(x−2t)=∫−∞+∞Wf(x,w)exp(jwt)dtx=2tf(t)f(0)⟹f(t)=f(0)1∫−∞+∞Wf(2t,w)exp(jwt)dt
- Inversion Formula