S变换 (s-transform)是一种时频分析 的工具。
和其他时频分析工具一样,透过S变换,我们可以同时从时域 以及频域 观察一个信号的能量分布。S变换的特别之处在它既保持与傅里叶变换 的直接关系,又可在不同频率有不同的分辨率。此外,S变换与小波变换 (wavelet transform)有密切的关系,或可视为连续小波变换 (continuous wavelet transform)的变形。S变换的清晰度略优于加伯变换 (Gabor transform),而不如韦格纳分布 (Wigner distribution function)、科恩克莱斯分布 、改良式韦格纳分布 (Modified Wigner distribution function)。
定义
一个信号x(t)的S变换为
S
x
(
t
,
f
)
=
∫
−
∞
∞
x
(
τ
)
|
f
|
e
−
π
(
t
−
τ
)
2
f
2
e
−
j
2
π
f
τ
d
τ
{\displaystyle S_{x}(t,\,f)=\int _{-\infty }^{\infty }x(\tau )\,|f|\,e^{-\pi (t-\tau )^{2}f^{2}}e^{-j2\pi f\tau }\,d\tau }
其中窗函数 为高斯窗函数
w
(
t
,
f
)
=
|
f
|
e
−
π
t
2
f
2
{\displaystyle w(t,f)=|f|e^{-\pi t^{2}f^{2}}}
另种表示-频谱表示式
借着折积定理
x
(
t
)
∗
h
(
t
)
=
F
−
1
(
X
(
f
)
⋅
H
(
f
)
)
{\displaystyle x(t)\ast h(t)={\mathcal {F}}^{-1}(X(f)\cdot H(f))}
S变换能以频域
X
(
f
)
{\displaystyle X(f)}
表示,
S
x
(
t
,
f
)
=
∫
−
∞
∞
(
x
(
τ
)
e
−
j
2
π
f
τ
)
(
|
f
|
e
−
π
(
t
−
τ
)
2
f
2
)
d
τ
{\displaystyle S_{x}(t,f)=\int _{-\infty }^{\infty }(x(\tau )e^{-j2\pi f\tau })(|f|e^{-\pi (t-\tau )^{2}f^{2}})\,d\tau }
这里可将
S
x
(
t
,
f
)
{\displaystyle S_{x}(t,f)}
看成
x
(
t
)
e
−
j
2
π
f
t
{\displaystyle x(t)e^{-j2\pi ft}}
与
|
f
|
e
−
π
t
2
f
2
{\displaystyle |f|e^{-\pi t^{2}f^{2}}}
的卷积,
将
x
(
τ
)
e
−
j
2
π
f
τ
{\displaystyle x(\tau )e^{-j2\pi f\tau }}
以及
|
f
|
e
−
π
(
t
−
τ
)
2
f
2
{\displaystyle |f|e^{-\pi (t-\tau )^{2}f^{2}}}
分别取傅里叶变换 可得
S
x
(
t
,
f
)
=
∫
−
∞
∞
X
(
f
+
α
)
e
−
π
α
2
/
f
2
e
j
2
π
α
t
d
α
{\displaystyle S_{x}(t,f)=\int _{-\infty }^{\infty }X(f+\alpha )\,e^{-\pi \alpha ^{2}/f^{2}}\,e^{j2\pi \alpha t}\,d\alpha }
窗函数(window function)
对一讯号
x
(
t
)
{\displaystyle x(t)}
进行S变换,写为
S
x
(
t
,
f
)
=
∫
−
∞
∞
x
(
τ
)
w
(
t
−
τ
,
f
)
e
−
j
2
π
f
τ
d
τ
{\displaystyle S_{x}(t,f)=\int _{-\infty }^{\infty }x(\tau )\,w(t-\tau ,f)e^{-j2\pi f\tau }\,d\tau }
其中
w
(
t
,
f
)
{\displaystyle w(t,f)}
为S变换之窗函数 ,常见之窗函数为高斯函数,即
w
(
t
,
f
)
=
|
f
|
e
−
π
t
2
f
2
{\displaystyle w(t,f)=|f|e^{-\pi t^{2}f^{2}}}
事实上,S变换之窗函数并不局限于高斯函数,考虑以下一般化型态的S变换(generalized S transform):
S
x
(
t
,
f
,
P
)
=
∫
−
∞
∞
x
(
τ
)
w
(
t
−
τ
,
f
,
P
)
e
−
j
2
π
f
τ
d
τ
{\displaystyle S_{x}(t,f,P)=\int _{-\infty }^{\infty }x(\tau )\,w(t-\tau ,f,P)e^{-j2\pi f\tau }\,d\tau }
其中
w
(
t
,
f
,
P
)
{\displaystyle w(t,f,P)}
为 一般化之窗函数(generalized window function), 参数 P 为一系列数值而组成之矩阵,可控制window function,举例来说
w
(
t
,
f
,
γ
)
=
|
f
|
e
−
π
t
2
f
2
/
γ
2
{\displaystyle w(t,f,\gamma )=|f|e^{-\pi t^{2}f^{2}/\gamma ^{2}}}
,此时
γ
{\displaystyle \gamma }
为矩阵 P 内之唯一数值
通过选取较小的
γ
{\displaystyle \gamma }
值,可得到时域上较窄的窗函数,提高时间轴的分辨率,然而,频率轴上的分辨率将因此变差。
改善时频两域分辨率的问题,可选择非对称的窗函数(asymmetric window),在正方向(forward direction)上有较陡峭之斜率,反方向(backward direction)上有较缓之斜率,借由牺牲终端时间的分辨率(一般而言较不重要)换取初始时间上更好的分辨率。[ 1]
S变换可以沿着时间轴方向积分,将可以得到
x
(
t
)
{\displaystyle x(t)}
的频谱
X
(
f
)
{\displaystyle X(f)}
。推导如下,
利用Gaussian window所包含面积等于1的特性,
∫
−
∞
∞
|
f
|
e
−
π
(
t
−
τ
)
2
f
2
d
t
=
|
f
|
∫
−
∞
∞
e
−
π
(
t
−
τ
)
2
f
2
d
t
=
1
{\displaystyle \int _{-\infty }^{\infty }|f|e^{-\pi (t-\tau )^{2}f^{2}}\,dt=|f|\int _{-\infty }^{\infty }e^{-\pi (t-\tau )^{2}f^{2}}\,dt=1}
因此,沿着时间轴t积分,
∫
−
∞
∞
S
x
(
t
,
f
)
d
t
=
∫
−
∞
∞
x
(
τ
)
[
∫
−
∞
∞
|
f
|
e
−
π
(
t
−
τ
)
2
f
2
d
t
]
e
−
j
2
π
f
τ
d
τ
=
X
(
f
)
{\displaystyle \int _{-\infty }^{\infty }S_{x}(t,f)\,dt=\int _{-\infty }^{\infty }x(\tau )\left[\int _{-\infty }^{\infty }|f|\,e^{-\pi (t-\tau )^{2}f^{2}}\,dt\right]\,e^{-j2\pi f\tau }\,d\tau =X(f)}
这表示S频谱是可逆的,同时也提供一个简单的逆变换。
x
(
τ
)
=
∫
−
∞
∞
[
∫
−
∞
∞
S
x
(
t
,
f
)
d
t
]
e
j
2
π
f
τ
d
f
{\displaystyle x(\tau )=\int _{-\infty }^{\infty }\left[\int _{-\infty }^{\infty }S_{x}(t,f)\,dt\right]\,e^{j2\pi f\tau }\,df}
=
∫
−
∞
∞
X
(
f
)
e
j
2
π
f
τ
d
f
{\displaystyle =\int _{-\infty }^{\infty }X(f)\,e^{j2\pi f\tau }\,df}
滤波应用(Filtering)
S变换如同其他时频分析变换,皆可以设计波器来达到消除噪声留下讯号的功用,
利用逆S变换,我们可以设计一个S域的滤波器U(t,f) ,对x(t)进行讯号处理
x
f
i
l
t
e
r
(
τ
)
=
∫
−
∞
∞
[
∫
−
∞
∞
S
x
(
t
,
f
)
⋅
U
(
t
,
f
)
d
t
]
e
j
2
π
f
τ
d
f
{\displaystyle x_{filter}(\tau )=\int _{-\infty }^{\infty }\left[\int _{-\infty }^{\infty }S_{x}(t,f)\cdot U(t,f)\,dt\right]\,e^{j2\pi f\tau }\,df}
离散时间S变换
S变换 相较于加伯变换 ,虽在清晰度有较好的改善,但也有其缺点,就是运算复杂度变高,积分的范围会随着
f
{\displaystyle f\,}
的增加而增加。
因此,这里利用上面推导的频谱 表示式来推导离散时间S变换
频谱表示式
S
x
(
t
,
f
)
=
∫
−
∞
∞
X
(
f
+
α
)
e
−
π
α
2
/
f
2
e
j
2
π
α
t
d
α
{\displaystyle S_{x}(t,f)=\int _{-\infty }^{\infty }X(f+\alpha )\,e^{-\pi \alpha ^{2}/f^{2}}\,e^{j2\pi \alpha t}\,d\alpha }
令
t
=
n
Δ
T
f
=
m
Δ
F
α
=
p
Δ
F
{\displaystyle t=n\Delta _{T}\,\,f=m\Delta _{F}\,\,\alpha =p\Delta _{F}}
Δ
T
{\displaystyle \Delta _{T}\,}
表示取样时间间隔
,
Δ
F
{\displaystyle ,\Delta _{F}\,}
表示取样频率
如果要使用FFT的方式来实作,必须另加条件
Δ
T
⋅
Δ
F
=
1
/
N
{\displaystyle \Delta _{T}\cdot \Delta _{F}=1/N}
首先先对
x
(
t
)
{\displaystyle x(t)}
做傅里叶变换得到
X
(
f
)
{\displaystyle X(f)}
X
[
m
Δ
F
]
=
1
N
∑
k
=
0
N
−
1
x
[
k
Δ
T
]
e
−
j
2
π
m
k
N
{\displaystyle X[m\Delta _{F}]={\frac {1}{N}}\,\sum _{k=0}^{N-1}x[k\Delta _{T}]\,e^{\frac {-j2\pi mk}{N}}}
接着带入频谱表示式中,
S
x
(
n
Δ
T
,
m
Δ
F
)
=
∑
p
=
0
N
−
1
X
[
(
p
+
m
)
Δ
F
]
e
−
π
p
2
m
2
e
j
2
p
n
N
{\displaystyle S_{x}(n\Delta _{T}\,,m\Delta _{F})=\sum _{p=0}^{N-1}X[(p+m)\,\Delta _{F}]\,e^{-\pi {\frac {p^{2}}{m^{2}}}}\,e^{\frac {j2pn}{N}}}
当 m=0 时,S变换就定义成
S
x
(
n
Δ
T
,
0
)
=
1
N
∑
k
=
0
N
−
1
x
[
k
Δ
F
{\displaystyle S_{x}(n\Delta _{T}\,,0)={\frac {1}{N}}\,\sum _{k=0}^{N-1}x[k\Delta _{F}}
]
流程
Step1 : 计算
X
[
p
Δ
F
]
{\displaystyle X[p\Delta _{F}]\,}
,这个步骤只需要计算一次。
Step2 : 计算
e
−
π
p
2
m
2
{\displaystyle e^{-\pi {\frac {p^{2}}{m^{2}}}}}
for
f
=
m
Δ
F
{\displaystyle f=m\Delta _{F}}
Step3 : 将
X
[
p
Δ
F
]
{\displaystyle X[p\Delta _{F}]}
移动至
X
[
(
p
+
m
)
Δ
F
]
{\displaystyle X[(p+m)\Delta _{F}]}
Step4 : 将Step2,Step3的结果相乘得到
B
[
m
,
p
]
=
X
[
(
p
+
m
)
Δ
F
]
⋅
e
−
π
p
2
m
2
{\displaystyle B[m,p]=X[(p+m)\Delta _{F}]\cdot e^{-\pi {\frac {p^{2}}{m^{2}}}}}
Step5 : 对B[m,p]取逆离散傅里叶变换 (IDFT)可得到,
S
x
(
n
Δ
T
,
m
Δ
F
)
{\displaystyle S_{x}(n\Delta _{T}\,,m\Delta _{F})}
在
f
=
m
Δ
F
{\displaystyle f=m\Delta _{F}}
的行向量
Step6 : 重复Step2~5直到
S
x
(
n
Δ
T
,
m
Δ
F
)
{\displaystyle S_{x}(n\Delta _{T}\,,m\Delta _{F})}
全部定义完成。
S变换特性
S变换与加伯变换(Gabor Transform)很相似,
G
x
(
t
,
f
)
=
∫
−
∞
∞
x
(
τ
)
e
−
π
(
t
−
τ
)
2
e
−
j
2
π
f
τ
d
τ
{\displaystyle G_{x}(t,f)=\int _{-\infty }^{\infty }x(\tau )\,e^{-\pi (t-\tau )^{2}}e^{-j2\pi f\tau }\,d\tau }
S
x
(
t
,
f
)
=
∫
−
∞
∞
x
(
τ
)
|
f
|
e
−
π
(
t
−
τ
)
2
f
2
e
−
j
2
π
f
τ
d
τ
{\displaystyle S_{x}(t,f)=\int _{-\infty }^{\infty }x(\tau )\,|f|e^{-\pi (t-\tau )^{2}f^{2}}e^{-j2\pi f\tau }\,d\tau }
唯一的不同就在于S变换的Gaussian Window的宽度会随着
f
{\displaystyle f}
改变。
低频
时域分辨率差
频域分辨率佳
高频
频域分辨率差
时域分辨率佳
原因就是
f
{\displaystyle f\,}
在高频时,Gaussian Window宽度变小,时域分辨率好;反之,
f
{\displaystyle f\,}
在低频时,Gaussian Window宽度变宽,频域分辨率好。
但是,当
f
→
0
{\displaystyle f\rightarrow 0}
时,Gaussian Window会无穷无尽的变宽,就丧失时频分析只做局部分析的精神。
一种解决的方式是:使Gaussian Window宽度不再因
f
2
{\displaystyle f^{2}\,}
改变\,产生带宽剧烈的变化,
S变换一般式
S
x
(
t
,
f
)
=
|
s
(
f
)
|
∫
−
∞
∞
x
(
τ
)
e
−
π
(
t
−
τ
)
2
s
2
(
f
)
e
−
j
2
π
f
τ
d
τ
{\displaystyle S_{x}(t,f)=|s(f)|\int _{-\infty }^{\infty }x(\tau )\,e^{-\pi (t-\tau )^{2}s^{2}(f)}e^{-j2\pi f\tau }\,d\tau }
s
(
f
)
{\displaystyle s(f)}
是一个相对平缓的曲线(见底下示意图),当
f
→
0
{\displaystyle f\rightarrow 0}
时,
s
(
f
)
≠
0
{\displaystyle s(f)\neq 0}
S变换是一种运算量高的时频分析工具,尤其在低频部分,Gaussian Window宽度变宽,频域分辨率比加伯变换来的好,所以S变换对于低频讯号分析比较有优势
例如:声音讯号,人耳对高频的部分没有太特别的感觉,但在低频部分却比较敏感,如:中央Do = 262Hz,高八度Do = 512Hz 可以很清楚的听出两个不同的音,
但10000Hz 和 10170Hz对人来说差别不大,再说人耳对3KHz以内的声音最敏感,所以能分析低频讯号就显得重要。
此时,就可以使用S变换,来强调低频讯号,而牺牲高频讯号。
与韦格纳分布的比较
韦格纳分布是时频分析工具中,具有高清晰度的一个,但最大的缺点是有交叉项(cross-term)的问题。若一个信号是由数个信号成分组合而成,那么使用韦格纳分布来分析时就会受到两两信号成分之间的交叉项干扰,这将会产生一些不必要的噪声。一个信号x的韦格纳分布为
W
x
(
t
,
f
)
=
∫
−
∞
∞
x
(
t
+
τ
2
)
x
∗
(
t
−
τ
2
)
e
−
j
2
π
τ
f
d
τ
{\displaystyle W_{x}(t,f)=\int _{-\infty }^{\infty }x(t+{\frac {\tau }{2}})x^{*}(t-{\frac {\tau }{2}})e^{-j2\pi \tau f}\,d\tau }
交叉项是在积分中两个x项相乘时产生的。S变换的计算原理与韦格纳分布不同,是直接对
x
(
τ
)
{\displaystyle x(\tau )}
进行变换,不会有交叉项的问题。
与加伯变换的比较
加伯变换的定义为
G
x
(
t
,
f
)
=
∫
−
∞
∞
x
(
τ
)
e
[
−
π
(
t
−
τ
)
2
f
2
]
e
(
−
j
2
π
f
τ
)
d
τ
{\displaystyle G_{x}(t,f)=\int _{-\infty }^{\infty }x(\tau )e^{\left[-\pi (t-\tau )^{2}f^{2}\right]}e^{\left(-j2\pi f\tau \right)}d\tau }
我们知道加伯变换是短时距傅里叶变换的一种特殊形式,其中只要把短时距傅里叶变换的窗函数用高斯函数来替代就成了加伯变换;S变换则可视为一种窗函数会随f变化的加伯变换;随着频率的升高,高斯函数在时域上的宽度会越来越窄,使得时域上的分辨率会增加,反之牺牲频域上的分辨率。
加伯变换和S变换原理相同,两者唯一不同的地方就是窗函数的
f
2
{\displaystyle f^{2}}
和强度
f
{\displaystyle f}
,基本上都是由短时距傅里叶变换延伸而来;两者共有的好处是不会像韦格纳分布一样会有交叉项;又S变换在低频时的频率分辨率会优于加伯变换。
与小波变换的关系
连续小波变换可以视为将一个信号对小波做相关 (correlation):
W
(
τ
,
d
)
=
∫
−
∞
∞
x
(
t
)
W
(
t
−
τ
,
d
)
d
t
{\displaystyle W(\tau ,\,d)=\int _{-\infty }^{\infty }x(t)\,W(t-\tau ,\,d)\,dt}
而S变换可以视为连续小波变换乘上一个相位项:
S
(
τ
,
f
)
=
e
−
j
2
π
τ
f
W
(
τ
,
d
)
{\displaystyle S(\tau ,\,f)=e^{-j2\pi \tau f}W(\tau ,\,d)}
而其用的母小波为:
w
(
t
,
f
)
=
|
f
|
exp
[
−
π
t
2
f
2
]
exp
[
−
j
2
π
f
t
]
{\displaystyle w(t,\,f)=|f|\,\exp[-\pi t^{2}f^{2}]\exp[-j2\pi ft]}
S变换,加伯变换和短时距傅里叶变换的比较
假设我们想要对一讯号
x
(
t
)
{\displaystyle x(t)}
做时频分析
x
(
t
)
=
c
o
s
(
2
π
t
)
,
t
<
10
,
{\displaystyle x(t)=cos(2\pi t),t<10,}
x
(
t
)
=
c
o
s
(
2
π
t
)
,
10
≤
t
<
20
,
{\displaystyle x(t)=cos(2\pi t),10\leq t<20,}
x
(
t
)
=
c
o
s
(
2
π
t
)
,
20
≤
t
.
{\displaystyle x(t)=cos(2\pi t),20\leq t.}
如果使用一个宽度为1秒的矩形函数来做短时距傅里叶变换会得到
如果使用一个加伯变换会得到时频分析图
如果用S变换会得到时频分析图
由上三图可知,S变换是一个在低频时频域分辨率高,高频时时域分辨率高的时频分析;举例来说,对于100Hz与300Hz和1000Hz与1200Hz这两组声音,哪一组对于人耳会有明显的差异呢? 答案是100Hz和300Hz,低频时的频率差异对于人耳较明显,当频率越高时,人耳就越难分辨出频率的差异;同样的道理,S变换便符合我们的需求,低频时讯号变化慢拥有低的时域分辨率和拥有高的频域分辨率,高频时因为讯号变化很快则拥有高的时域分辨率和低的频域分辨率。
参考文献
R. G. Stockwell, L. Mansinha, and R. P. Lowe, "Localization of the complex spectrum: the S transform," IEEE Trans. Signal Processing, vol. 44, no. 4, pp. 998–1001, Apr. 1996.
Jian-Jiun Ding, Time frequency analysis and wavelet transform class note, the Department of Electrical Engineering, National Taiwan University (NTU), Taipei, Taiwan, 2010.
Sitanshu Sekhar, Ganapati Panda and Nithin V George, "An Improved S-Transform for Time-Frequency Analysis," "IACC2009", pp. 315-319, March 2009.
^ The S-transform with windows of arbitrary and varying shape, C. Robert Pinnegar, Lalu Mansinha, GEOPHYSICS, VOL.68, NO.1 (JANUARY-FEBRUARY 2003); P.381–385