数字信号处理的几个前沿课题

第10章 数字信号处理的几个前沿课题

前面介绍了数字信号处理的基本知识，本章我们将介绍时谱分析、小波变换、地震观测系统仿真与地面运动恢复等几个数字信号前沿课题，以便大家在实际工作中参考。

10.1 时谱（倒谱）分析

时谱分析(Cepstrum analysis)是一种非线性信号处理技术，它在语言、图像、和噪声处理领域中都有广泛的应用。

时谱可分为两类：复时谱和功率时谱。MATLAB 信号处理工具箱提供复时谱分析的工具函数。

复时谱（Complex cepstrum ）的定义为：

[]{}ωπ

ωπ

π

ω

d e

e X n x

n

j j ⎰-=)(ln 21

)(ˆ （10-1）

由上式可见，复时谱实际上是序列x(n)的Fourier 变换取自然对数，再取Fourier 逆

变换，得到的复时谱仍然是一个序列。也就是说，复时谱是x(n)从时间域至频率域、频率域至频率域、频率域至时间域的三次变换。

MATLAB 信号处理工具箱函数cceps 用于估计一个序列x 的复时谱，调用格式为：

xhat=cceps(x)

式中，x 为输入序列（实序列）；xhat 为复时谱（复序列）。

MATLAB 信号处理工具箱还提供了序列实时（倒）谱的计算程序rceps ，调用格式为Y=rceps(x)，其中x 为实序列；y 为实时谱，执行的操作为：

ωπ

π

πω

d e

X C j x ⎰-

=

)(ln 21 （10-2）

由此可知，我们不能从序列x 的实时谱重构原始序列，因为实时谱是根据序列Fourier

变换的幅值计算的，丢失了相位方面的信息。但如果需要，可采用最小相位模式估计原始序列。

由于复时谱从复频谱计算得到，不损失相位信息，因此复时谱是可逆的，实时谱过程是不可逆的。

时谱分析技术广泛地应用于语言信号分析、同态滤波技术中。 这里举一个说明复时谱在具有回声信号测量中的应用。

【例10-1】设原信号是一个45Hz 的正弦波，在传播过程中遇到障碍产生回声，回声振幅衰减为原信号的0.5，并与原信号有0.2s 的延迟。在某测点测到的信号是原信号和回声信号的叠加。试使用复时谱分析该测点的信号。

%Samp10_1

t=0:0.01:1.49; %时间信号，采样间隔为0.01s sig=sin(2*pi*45*t); %原始信号

echo=0.5*[zeros(1,20),sig(1:length(t)-20)];

%在信号前面补加20个零，并使其振幅衰减一半作为回声信号 sigecho=sig+echo;

c=cceps(sigecho); %求原始信号与回声信号叠加的复时谱

subplot(2,1,1),plot(t,sigecho); %绘制原始信号与回声的叠加信号

xlabel('时间/s');axis([0 1.5 -1 1]);grid on

title('时间域信号')

subplot(2,1,2),plot(t,c); %绘制其复时谱

xlabel('时间/s');axis([0 1.5 -1 1]);grid on

title('时谱域信号')

程序运行结果为图10-1。可以看到，信号的复时谱在t=0.2秒处有一个峰值，这就是回声信号。表示运用复时谱可以明确检测到回声信号。在数字地震仪监测核爆破的工作中，如果采用两次核爆在几乎同一时间、同一地点发生，爆破产生的地震波经过地球内部传到同一地震台上的传播路径几乎一致，并且爆破源产生的波形相差不大，因此在地震台上表现为顺序到达的两个波形类似的波，若后面的波形被淹没，可以采用复时谱的方法进行检测，明确给出是否是两次爆破事件，如果是两次爆破事件，就可以明确给出两次事件的时间间隔。

图10-1 正弦波与回声信号叠加的波形和时谱形状

10.2 地震观测系统的仿真和地面运动的恢复

由前面的线性系统的讲解可知，一个系统可以用它的系统函数或脉冲响应来表示：

)()()(t x t h t y ⊗= （10-3）

式中，x(t)为输入信号，对于地震观测系统来讲为地面运动；y(t)为系统的输出，对于地震观测系统来讲为地震记录；h(t)为系统的脉冲响应。在频率域内，根据卷积定理，该式可以表示为：

)()()(ωωωX H Y = （10-4）

式中，)(ωH 为系统的传递函数，()ωX 和()ωY 分别为x(t)和y(t)的Fourier 变换。

现在设想一个频带范围很宽的线性系统，比如宽频带地震仪，其系统函数为()ωH ;另一个频带较窄的系统，比如短周期地震仪，其系统函数为()ω1H ，对于同样的输入()ωX 有：

)()()(ωωωX H Y =，)()()(11ωωωX H Y = （10-5）

式中，()ω1Y 为频带较窄的系统记录的频谱；()ω1H 为频带较窄系统的传递函数。比较前面两式可得：

()()()()

ωωωωH X H Y 11=

（10-6）

将上式变换到时间域就得到频带较窄系统的输出y 1(t)。也就是说，如果知道了宽频带和窄频带系统的传递函数()ωH 和()ω1H ，原则上可以从宽频带系统的输出推测出频带较窄系统的输出。但如果我们知道窄频带系统的输出及其两种系统的传递函数，无法得到宽频带系统的输出。这样就使得我们在记录某种信号时采用宽频带记录，然后仿真到其他各种窄频带的记录仪器上对信号进行分析。

如果已知地震仪的输出和地震仪的传递函数，我们可以求出地面运动为：

)()()(ωωωH Y X = （10-7）

下面我们就举例给出将宽频带系统的输出仿真到窄频带输出及地面运动恢复的方法。 【例10-2】 设计一个Butterworth 模拟宽频带带通滤波器，设计指标为：通带频率：1~40Hz ，低端阻带边界0.2Hz ，高端阻带边界42.5Hz ，通带波纹1dB ，阻带衰减20dB 。再设计一个窄频带带通滤波器，设计指标为通带频率：5~7Hz ，低频段过渡带宽 1.5Hz ，高频段过渡带宽0.5Hz ，通带波纹1dB ，阻带衰减20dB 。假设一个信号

t f t f t f t x 3212sin 5.02cos 5.02sin )(πππ++=，其中f1=6Hz,f2=10Hz,f3=19Hz 。信号的采

样频率为50Hz 。（1）模拟输入信号通过宽频带滤波器的输出信号并与原来信号进行比较。（2）运用宽带滤波器的传递函数和输出得到输入信号，并与原输入信号进行比较（恢复地面运动）（3）模拟输入信号通过窄频带滤波器的输出信号并与原来信号进行比较。（4）将在宽频带