将时域信号转换成频域信号的方法与系统

摘要

本发明涉及数字信号处理领域，公开了一种将时域信号转换成频域信号的方法与系统。该方法将时域信号用曲线平滑后再生为近似模拟信号，传递给进行有阻尼的受迫简谐振动的谐振系统，并计算该谐振系统的谐振强度，以作为频域信号强度。因为可以任选频点计算，所以在只要计算特定的一个或几个频点的信号时，计算量远低于传统计算方法。进一步通过计算频域信号强度变化过程，得到时间-频率分布图，若将相应频点的所有采样点对应时刻的频域信号强度叠加即可得到该频点的频域信号强度，从而成功将时域信号转换成频域信号。由于是逐个采样点计算，时间响应最高可达到采样周期并且响应幅度可调。

说明书

技术领域

本发明涉及数字信号处理领域，特别涉及将时域信号转换成频域信号的方 法与系统。

背景技术

快速傅里叶变换是目前对数字信号数据最常用的处理方法之一，它可以将 信号从时域变换到频域，产生统计频谱，进行进一步研究。而为了对诸如语音 信息等频谱不断变化的信号进行研究，发展了时频变换方法。时频变换方法目 前一些主要的方法是短时傅里叶变换和小波变换，以及魏格纳-维尔变换等，但 最主要的方法是短时傅里叶变换。

短时傅里叶变换(Short Time Fourier Transform)具体做法是将数字信号划 分成固定长度的帧，相邻帧之间有50%～90%的重叠，对每一帧数据先加窗处理， 然后做快速傅里叶变换(Fast Fourier Transform)得到该帧数据的一个频率域分 布，将所有帧的频率域分布按时间顺序排列，就得到了三维的时间-频率分布 图。

短时傅里叶变换具有速度快，计算量小的特点，但是其产生的时间-频率 分布图，频率的分割点数永远等于数据帧长度的一半，是均匀分布的，并且由 于功率值反映的是每帧数据对应时间片内的平均值，对瞬时变化的时间响应精 度只能达到每帧时长的1/2-1/5。因此，给定一个采样频率后，变换所得到时频 分布图的时间精度和频率精度的乘积基本为定值，我们只能在时间精度和频率 精度之间被迫作折衷选择，难以得到同时具有高频率精度和高时间精度的时-频 分布图。

对人耳蜗的研究表明，耳蜗有大约24000个听觉传感器，每个传感器都 对特定频率的信号产生强烈反应并将信号幅度转换成电信号，众多听觉传感 器通过神经束以并行的方式将这些电信号传递到大脑，大脑处理这些信号产 生听觉。所以，从理论上我们也应该可以有一种类似听觉传感器的方法对信 号进行并行的时频变换，能够在时间响应和频率精度上均高于短时傅里叶变 换的转换方法。

发明内容

本发明的目的在于提供一种将时域信号转换成频域信号的方法与系统， 时间响应最高可达到采样周期并且响应幅度可调。

为解决上述技术问题，本发明的实施方式公开了一种将时域信号转换成 频域信号的方法，包括以下步骤：

根据时域信号中前一采样点和当前采样点确定平滑曲线f(t)的函数形式 和相应参数；

将f(t)的参数和前一采样点对应时刻的状态参数输入到一个偏微分方程 的解析解中，求得当前采样点对应时刻的状态参数，其中该偏微分方程代表 一个有阻尼的受迫谐振系统，形式为k是该谐振系 统无阻尼时的谐振角频率，n为该谐振系统的阻尼因子，k和n对应该谐振系 统的一个共振角频率ω_r，状态参数中包含该谐振系统的速度v和位移x；

根据当前采样点对应时刻的状态参数中的速度v和位移x计算谐振强度， 作为频域中与共振角频率ω_r对应的频点的在当前采样点对应时刻的频域信 号强度。

本发明的实施方式还公开了一种声纹识别方法，在声纹识别的过程中使 用上文的方法根据声音信号产生时间-频率分布图。

本发明的实施方式还公开了一种地震波分析方法，使用上文的方法根据 地震波产生时间-频率分布图，以对地震波进行时频分析。

本发明的实施方式还公开了一种电力线路谐波分析方法，使用上文的方 法根据电力线路谐波产生时间-频率分布图，以对电力线路谐波进行分析。

本发明的实施方式还公开了一种雷达的回波分析方法，使用上文的方法 根据雷达的回波信号产生时间-频率分布图，以对回波信号进行时频分析。

本发明的实施方式还公开了一种将时域信号转换成频域信号的系统，包 括：

平滑模块，用于根据时域信号中前一采样点和当前采样点确定平滑曲线 f(t)的函数形式和相应参数；

谐振模块，用于将f(t)的参数和前一采样点对应时刻的状态参数输入到 一个偏微分方程的解析解中，求得当前采样点对应时刻的状态参数，其中该 偏微分方程代表一个有阻尼的受迫谐振系统，形式为k是该谐振系统无阻尼时的谐振角频率，n为该谐振系统的阻尼因子， k和n对应该谐振系统的一个共振角频率ω_r，状态参数中包含该谐振系统的 速度v和位移x；

计算模块，用于根据当前采样点对应时刻的状态参数中的速度v和位移 x计算谐振强度，作为频域中与共振角频率ω_r对应的频点的在当前采样点对 应时刻的频域信号强度。

本发明实施方式与现有技术相比，主要区别及其效果在于：

将时域信号用曲线平滑后再生为近似模拟信号，传递给进行有阻尼的受 迫简谐振动的谐振系统，并计算该谐振系统的谐振强度，以作为频域信号强 度，因为是采用特定形式曲线来作为外强迫力，从而可以直接采用解析解进 行精确计算，避免了用诸如欧拉法、威尔逊-θ法或线性加速法等积分算法来 求近似解，在提高计算准确性的同时极大地降低了计算量，并且由于可以任 选频点计算，所以在只要计算特定的一个或几个频点的信号时，计算量远低 于像快速傅里叶变换这类传统计算方法。

该系统包括平滑模块、谐振模块和计算模块，其中该平滑模块将时域信 号用曲线平滑后再生为近似模拟信号，传递给该谐振模块进行有阻尼的受迫 简谐振动，计算模块计算该受迫简谐振动的谐振强度，以作为频域信号强度， 由于可以任选频点计算，所以在只要计算特定的一个或几个频点的信号时， 计算量远低于传统计算方法，并且由于是可以直接采用解析解进行计算，从 而不需要用诸如欧拉法、威尔逊-θ法或线性加速法等积分算法来求近似解， 在保持计算准确性的同时极大地降低了计算复杂度。

进一步地，通过计算频域信号强度变化过程，进而得到时间-频率分布 图，若将相应频点的所有采样点对应时刻的频域信号强度叠加即可得到该频 点的频域信号强度，从而成功将时域信号转换成频域信号，因为可以任选频 点计算，频率精度与采样频率无关，所以即使时域信号采样频率低，频率精 度也可以远高于快速傅里叶变换，而且还可以任意选择感兴趣的频率区间来 设置任意多个频率点进行转换，由于是逐个采样点计算，时间响应最高可接 近采样周期，并且响应幅度可调，因此可以得到同时具有高频率精度和高时 间精度的时间-频率分布图。

进一步地，对时域信号中各采样点或具有不同共振角频率的谐振系统采 用相同或不同的曲线函数形式进行平滑连接，可适应一些特殊需求的场合。

进一步地，通过对共振角频率ω_r对应的受迫谐振系统的阻尼因子n进行 调整，可以提高受迫谐振系统的频域分辨精度或时域分辨精度。

进一步地，通过不同方法得到与各共振角频率ω_r对应的谐振系统的阻尼 因子n和谐振系统无阻尼时的谐振角频率k，可根据实际需求选择相应的谐 振系统。

附图说明

图1是本发明第一实施方式中一种将时域信号转换成频域信号的方法的 流程示意图；

图2是本发明第一实施方式中有阻尼的谐振系统在受到外强迫力持续作 用后系统的能量变化的时间特性曲线；

图3是本发明第一实施方式中有阻尼的谐振系统在受到外强迫力持续作 用时系统的谐振强度的频率特性曲线；

图4是本发明第一实施方式中采样值点由平滑曲线分段连接的示意图；

图5是本发明第二实施方式中一种将时域信号转换成频域信号的方法的 流程示意图；

图6是本发明第三实施方式中一种将时域信号转换成频域信号的方法的 流程示意图；

图7a至7g是本发明第四实施方式中对语音片段、粉红噪声、调频调幅 混合声，分别由短时傅里叶变换和本发明转换方法进行转换得到的频谱图结 果；

图8是本发明第八实施方式中一种将时域信号转换成频域信号的系统的 结构示意图；

图9是本发明第九实施方式中一种将时域信号转换成频域信号的系统的 结构示意图。

具体实施方式

在以下的叙述中，为了使读者更好地理解本申请而提出了许多技术细节。 但是，本领域的普通技术人员可以理解，即使没有这些技术细节和基于以下 各实施方式的种种变化和修改，也可以实现本申请各权利要求所要求保护的 技术方案。

为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发 明的实施方式作进一步地详细描述。

本发明第一实施方式涉及一种将时域信号转换成频域信号的方法。图1 是该将时域信号转换成频域信号的方法的流程示意图。

具体地说，如图1所示，该将时域信号转换成频域信号的方法包括以下 步骤：

在步骤101中，根据时域信号中前一采样点和当前采样点确定平滑曲线 f(t)的函数形式和相应参数。

可以理解，f(t)可以对于时域信号中各采样点和具有不同共振角频率ω_r的谐振系统采用相同平滑曲线函数形式，也可以对于时域信号中各采样点或 具有不同共振角频率ω_r的谐振系统采用不同平滑曲线函数形式。

例如，f(t)可以采用a+b cos(pt)的形式，其中p为时域信号的采样角频 率的一半，a为前一采样值与当前采样值的算数平均值，b为前一采样值与当 前采样值的差的一半，可以采用a′+b′p′t的形式，其中p′为时域信号的采样 频率，a′为前一采样值，b′为前一采样值与当前采样值的差，可以采用a₀+a₁p′t +a₂p′²t²+a₃p′³t³的形式，其中p′为时域信号的采样频率，也可以采用a+ bcos(pt)、a’+b’p’t和a₀+a₁p′t+a₂p′²t²+a₃p′³t³的形式组合。

此外，可以理解，在本发明的其他实施方式中，f(t)还可以采用其他的平 滑曲线函数形式。

对时域信号中各采样点或具有不同共振角频率的谐振系统采用相同或 不同的曲线函数形式进行平滑连接，可适应一些特殊需求的场合。

需注意的是，这些平滑曲线的参数可以通过直接输入前后两个或多个采 样点的值来确定，也可以只将相关采样点输入到一个算法，由该算法进行预 定曲线形式的曲线拟合或平滑，最终得到曲线参数，例如三次样条曲线、三 次拟合曲线等。只要这些平滑曲线可以通过或接近波形图中前一采样点和当 前采样点就可以了。

此后进入步骤102，将f(t)的参数和前一采样点对应时刻的状态参数输入 到一个偏微分方程的解析解中，求得当前采样点对应时刻的状态参数。其中 该偏微分方程代表一个有阻尼的受迫谐振系统，形式为k是该谐振系统无阻尼时的谐振角频率，n为该谐振系统的阻尼因子， k和n对应该谐振系统的一个共振角频率ω_r，状态参数中包含该谐振系统的 速度v和位移x。

可以理解，求一个偏微分方程的解析解是一种成熟的数字方法，在教科 书中可以找到，本申请中不进行详细说明了。

此外，在本发明的各个实施方式中，可以通过以下几个方式来得到谐振 系统的k和n：

1）谐振系统的阻尼因子n不变且为大于零的实数，由得 到与各共振角频率ω_r对应的谐振系统无阻尼时的谐振角频率k；

2）谐振系统的振幅衰减倍数C不变且为大于零的实数，由得 到与各共振角频率ω_r对应的谐振系统无阻尼时的谐振角频率k，再由 得到与各共振角频率ω_r对应的阻尼因子n，其中振幅衰减倍数C 为f(t)的振幅与谐振系统的谐振振幅A的比例，即定义稳态时表示外强迫力 的f(t)的振幅与谐振系统模型因为共振而产生的谐振振幅A的比例定义为振 幅衰减倍数C；

3）谐振系统的能量衰减倍数Q不变且为大于零的实数，由得到与各共振角频率ω_r对应的谐振系统无阻尼时的谐振角频率k， 再由得到与各共振角频率ω_r对应的阻尼因子n，其中能量衰减 倍数Q为f(t)的振幅的平方与谐振系统与质量无关的能量E’的比例，即定义 稳态时表示外强迫力的f(t)的振幅的平方与谐振系统模型因为共振而累积的 与质量无关的能量E’的比例定义为能量衰减倍数Q。

通过不同方法得到与各共振角频率ω_r对应的谐振系统的阻尼因子n和谐 振系统无阻尼时的谐振角频率k，可根据实际需求选择相应的谐振系统。

当然，可以理解，在本发明的其他实施方式中，也可以使用其他方式来 得到谐振系统的k和n。

此后进入步骤103，根据当前采样点对应时刻的状态参数中的速度v和 位移x计算谐振强度，作为频域中与共振角频率ω_r对应的频点的在当前采样 点对应时刻的频域信号强度。可以理解，在本发明的各个实施方式中，上述 谐振强度可以用与质量无关的能量E’=k²x²+v²来表示，也可以用谐振振幅 $A=\sqrt{x^2+\frac{v^2}{k^2}}$来表示。

将时域信号用曲线平滑后再生为近似模拟信号，传递给进行有阻尼的受 迫简谐振动的谐振系统，并计算该谐振系统的谐振强度，以作为频域信号强 度，因为是采用特定形式曲线来作为外强迫力，从而可以直接采用解析解进 行精确计算，避免了用诸如欧拉法、威尔逊-θ法或线性加速法等积分算法来 求近似解，在提高计算准确性的同时极大地降低了计算量，并且由于可以任 选频点计算，所以在只要计算特定的一个或几个频点的信号时，计算量远低 于像快速傅里叶变换这类传统计算方法。

更具体地说，该将时域信号转换成频域信号的过程主要如下：

众所周知，对于具有阻尼和外强迫力的弹簧简谐振动系统，定义m为谐 振子质量，x为谐振子的瞬时偏移平衡点的位置，v为瞬时速度，a为瞬时加 速度。

假设：

1)阻尼作用与瞬时速度v成正比，阻尼系数定为μ，阻力R=-μv；

2)弹簧的弹性系数为K，回复力F=-Kx；

3)外强迫力为Hsinβt，其中H为外力的振幅，β为外力的角频率；

则系统有以下瞬时力平衡关系：ma+μν+Kx＝Hsinβt

将上式以时间微分形式表示：

为阻尼率；

是谐振系统无阻尼时的固有谐振频率；

设：其中k就是ω₀，代入上式，则上式可转化 为：

这就是有阻尼强迫振动系统(即有阻尼的受迫谐振系统)的去量纲化微分 方程。

其中k为谐振系统无阻尼时的固有频率（即谐振角频率），n称为阻尼 因子，n/k称为阻尼比，hsinβt称为外强迫力。

可知这个有阻尼谐振系统的无外强迫力时的自由振动角频率ω为：

$\omega =\sqrt{k^2-n^2}---\left(2\right)$

而这个有阻尼谐振系统的有外强迫力时的共振角频率ω_r为：

${\omega }_r=\sqrt{k^2-{2n}^2}---\left(3\right)$

众所周知，有阻尼的谐振系统在受到外强迫力f(t)＝hsinβt持续作用后， 产生受迫振动的同时，也将外强迫力带来的能量不断的转化为系统自身的动 能和势能，同时系统的振幅不断增大。但因为阻尼作用的存在，随着系统振 幅的增大，阻尼消耗的能量也逐渐增大，最终系统的传入能量会与阻尼损耗 的能量达到平衡，使系统达到一个稳态；而当外强迫力消失，系统的能量会 因为阻尼衰减而快速降低，振幅也逐渐减弱到零。

图2展示了这种变化的时间特性。上部曲线为输入信号能量，下部曲线 为谐振系统振幅。输入信号的能量反映在谐振系统振幅上，有快速上升，平 稳维持和阻尼衰减三个阶段。快速上升阶段T_r，谐振系统接收到的信号能量 高于被阻尼衰减掉的能量，因此谐振幅度快速上升。随着谐振振幅逐渐增高， 阻尼作用也同时增强，当阻尼所衰减的能量与输入信号的能量达到平衡时， 谐振幅度不再增高而是达到一个稳定值。阻尼衰减阶段T_d，当输入信号消失 时，随着谐振系统积蓄的能量不断被衰减，谐振振幅也逐渐降低直至消失。

外强迫力的角频率β越接近于谐振系统的共振角频率ω_r，系统的能量积 蓄作用就越强，振幅越大。如果两个角频率相等，则谐振幅度将达到最大(即 共振状态)。图3展示了这种变化的频率特性，如图3所示，当外强迫力的角 频率β等于共振角频率ω_r时，该谐振系统产生最大谐振振幅A_m，而当外强 迫力的角频率β偏离共振角频率ω_r，如等于ω_a、ω_b时，该谐振系统的振幅 变小为

数字信号（即时域信号）在以时间为横坐标轴，信号幅度为纵坐标轴的 波形图上，可以表示为一系列水平间隔相等但垂直位置不同的采样值点，如 图4中实心三角形标识，相邻二个采样值点的水平间距就是采样周期；为了 从这些采样值点重建原始信号，可以将波形图上每一对相邻采样值点之间用 平滑曲线段相连，这些平滑曲线段首尾相接，就构成了近似的连续模拟信号， 如图4中连接各采样值点的虚线和实线，每一段平滑曲线可以用f(t)表示， 其参数由相邻两个或多个采样值以及必要的其它因素，通过一个与曲线形式 对应的算法计算得到。

将(1)式中的外强迫力hsinβt用由离散数据拟合所得的分段平滑曲线f(t) 代替，则新的微分方程就可以描述在受到f(t)影响下，谐振系统在一个采样周 期中的运动过程，即该系统就构成了一个谐振传感器。

用平滑曲线f(t)代替外强迫力hsinβt后（1）式变为：

$\frac{d^{2}x}{{\mathrm{dt}}^2}+2n\frac{\mathrm{dx}}{\mathrm{dt}}+k^{2}x=f\left(t\right)---\left(4\right)$

对于（4）式，它属于二阶常系数线性非齐次微分方程，它的特征方程 $\frac{d^{2}x}{{\mathrm{dt}}^2}+2n\frac{\mathrm{dx}}{\mathrm{dt}}+k^{2}x=0$的根为

因为有阻尼谐振系统的自由振动角频率：

$\omega ={\omega }_0\sqrt{1-{\zeta }^2}=\sqrt{\frac{K}{m}-\frac{{\mu }^2}{{4m}^2}}$

将代入，得到（2）式，

因为我们要求的是某个实频率ω下的解，这样根据(2)式就有k>n的限 制，即阻尼比因此特征方程的根为一对共轭复根：

$r_1=-n+i\sqrt{k^2-n^2},$$r_2=-n-i\sqrt{k^2-n^2},$或表示为：r_1,2＝-n±iω。

所以（4）式的暂态解为：

X＝e^-nt(C₁cosωt+C₂sinωt)   (5)

其中C₁,C₂为待定系数。根据微分方程的求解知识，只要(4)式等号右边 的曲线函数f(t)的形式能使微分方程可以得到稳态解y^*，则就可以最终得到x 的解析解:

x＝X+y^*   (6)

从而不需要用诸如欧拉法、威尔逊-θ法或线性加速法等积分算法来求近 似解，增大计算复杂度。

因为(4)式是二阶常系数非齐次线性方程，对应的可以使用待定系数法， 常数变易法以及拉普拉斯变换等方法求解该方程得到解析解。

以待定系数法为例：众所周知，当f(t)满足以下两种常见形式时，我们 可以直接得到稳态解y^*：

（1）f(x)＝P_m(x)e^λx，其中λ是常数，P_m(x)是x的一个m次多项式：

P_m(x)＝a₀x^m+a₁x^m-1+…+a_m-1x+a_m；

（2）f(x)＝e^λx[P_l(x)cosωx+P_n(x)sinωx]，其中λ、ω是常数，P_l(x)、P_n(x) 分别是x的l次、n次多项式，其中有一个可以为零；

符合上述要求的平滑曲线有多种形式，比如最简单的a+bx直线型，三 阶多项式曲线a₀+a₁x+a₂x²+a₃x³型，以及a+be^λx型，b cos px型等。针对上述 每一种类型，都是可以得到一个解析解的。

得到了解析解x，也就得到了瞬时位置x和瞬时速度v(即x的导数)的计 算公式，就可以由任意一个采样值点所对应时刻的谐振系统的初始位置x₀和 初始速度v₀，以及连接下一采样值点的平滑曲线的参数，计算出下一采样值 点所对应时刻的瞬时位置x和瞬时速度v。对数字信号连续做上述计算，就 可以模拟出谐振系统的振动全过程。

谐振系统的谐振强度可以从每一步计算得到的瞬时位置x和瞬时速度v 计算出来。

因为谐振系统积蓄的总能量E等于弹性势能与瞬时动能的和：

$E=\frac{1}{2}{\mathrm{Kx}}^2+\frac{1}{2}{\mathrm{mv}}^2---\left(7\right)$

将代入(7)式，并令则上式转化为与质量m无关的能量 值E′＝k²x²+ν²   (8)

因此可以用去质量化的与质量无关的能量E’表示在某一时刻谐振系统 的瞬时能量，称为与质量无关的能量。

因为谐振系统的最大位移（即谐振振幅A）出现在瞬时速度v=0的时刻， 此时的系统能量和一起代入(7)式后得到谐振振幅A的算法：

$A=\sqrt{x^2+\frac{v^2}{k^2}}---\left(9\right)$

与质量无关的能量E’和谐振振幅A之间的关系为：E’＝k²A²

为了让谐振模型满足不同频率的转换目的，需要根据转换目标来设计谐 振模型每一个共振角频率ω_r的参数n和k，这里以hsinβt作为外强迫力来说 明参数设计方案。

根据谐振系统的特性，当外强迫力hsinβt的角频率β与谐振系统的共振 角频率ω_r相等时，即：${\omega }_r=\beta =\sqrt{k^2-{2n}^2}---\left(10\right)$

此时系统有最大谐振振幅A：

第一种谐振模型设计方案：所有共振角频率下的谐振振幅衰减速率相同， 即阻尼因子n固定不变；

该设计的特性是谐振振幅的衰减倍数随着固有频率的升高而升高。

设给定的共振角频率为ω_r则按如下步骤确定n和k：

a)给n指定>0的任意实数；

b)由$k^2={\omega }_r^2+{2n}^2$计算出k。

第二种谐振模型设计方案：要求外强迫力振幅的平方与谐振系统模型因 为共振而累积的与质量无关的能量E’的比值对所有频率保持不变；

定义外强迫力h sin βt的振幅的平方h²与与质量无关的能量E’的比值为能 量衰减倍数Q，即与(8)（9）式一同代入(11)式,则有如下关系：

4n²(k²-n²)＝k²Q   (12)

因此由(10)式和(12)式可以得到以下设定方法：

设给定的共振角频率为ω_r按如下步骤确定n和k：

a)给Q指定>0的任意实数；

b)由$k^2=\frac{Q+\sqrt{{4\omega }_r^4+Q^2}}{2}$计算出k²；

c)由$n=\sqrt{\frac{k^2-{\omega }_r^2}{2}}$计算出n。

当共振角频率ω_r逐渐增大时，k²趋近于n趋近于

该设计的特性是谐振振幅的衰减速度随着频率升高而稳定直至趋近于

第三种谐振模型设计方案：外强迫力振幅与谐振振幅的比例对所有频率 保持不变；该设计的特性是谐振振幅的衰减速度随着频率升高而减慢，直至 趋近于零。

定义外强迫力h sin βt的振幅与谐振系统模型因为共振而产生的谐振振幅 A的比例定义为振幅衰减倍数C，即代入(11)式，则有如下关系：

4n²(k²-n²)＝C²   (13)

因此由(10)式和(13)式可以得到以下设定方法：

设给定的共振角频率为ω_r则按如下步骤确定n和k：

a)给C指定>0的任意实数；

b)由$k^2=\sqrt{{\omega }_r^2+C^2}$计算出k；

c)由$n=\sqrt{\frac{k^2-{\omega }_r^2}{2}}$计算出n。

由于瞬时位置x和瞬时速度v是以相邻两个采样值点为区间段进行计算 得到的，影响其每一段区间计算结果的因素包括区间的初始状态值x₀,v₀和 连接相邻采样值点的表示外强迫力的f(t)的曲线参数，以及与f(t)形式对应的 x和v的解析解。所以，在某些特殊需求的场合下，对不同共振频率的谐振 模型可以应用不同形式的曲线进行平滑连接；或者更进一步，对同一个共振 频率的谐振模型的不同采样值点之间应用不同的表示外强迫力的f(t)的曲线 形式，即使用多类型平滑曲线的组合来平滑连接所有采样值点，只要确保始 终应用与当前连接相邻采样值点的平滑曲线f(t)形式相对应的x和v的解析解， 就可以正确计算x和v，从而实现振动能量累积。

本发明第二实施方式涉及一种将时域信号转换成频域信号的方法。图5 是该将时域信号转换成频域信号的方法的流程示意图。

具体地说，如图5所示，第二实施方式在第一实施方式的基础上进行了 改进，主要改进之处在于，上述方法还包括以下步骤：

在步骤304中，按时序获取与共振角频率ω_r对应的频点的在时域信号 中所有采样点对应时刻的频域信号强度，以形成与共振角频率ω_r对应的频点 的时间-频域信号强度曲线。

此后进入步骤305，将预定的频率序列ω₁，…，ω_m乘以2π并依次作为 谐振系统的共振角频率ω_r，以得到与各共振角频率ω_r对应的频点的时间- 频域信号强度曲线。

此后进入步骤306，将与各共振角频率ω_r对应的频点的时间-频域信号 强度曲线按预定的频率序列排列，即得到时间-频率分布图。

通过计算频域信号强度变化过程，进而得到时间-频率分布图，若将相 应频点的所有采样点对应时刻的频域信号强度叠加即可得到该频点的频域信 号强度，从而成功将时域信号转换成频域信号，因为可以任选频点计算，频 率精度与采样频率无关，所以即使时域信号采样频率低，频率精度也可以远 高于快速傅里叶变换，而且还可以任意选择感兴趣的频率区间来设置任意多 个频率点进行转换，由于是逐个采样点计算，时间响应最高可达到采样周期 并且响应幅度可调，因此可以得到同时具有高频率精度和高时间精度的时间- 频率分布图。

经过上述改进，转换方法的基本步骤为：

选取一种平滑曲线的函数形式，作为所述谐振系统模型的表示外强迫力 的f(t)部分，求解所述谐振系统模型得出所述的瞬时位移x和所述的瞬时速 度v的解析解；构造一组能让所述谐振系统模型在指定频率（相当于对应共 振角频率ω_r的频点）产生共振的k和n，然后计算出连接数字信号每一对相 邻采样值点的平滑曲线段的参数，代入所述的解析解，连续计算所述的指定 频率下的每一个所述采样值点对应时刻的所述谐振强度，得到所述指定频率 下的时间-谐振强度曲线（即时间-频域信号强度曲线）；计算同一组数字信 号在不同所述指定频率下的所述时间-谐振强度曲线，按频率高低顺序排列这 些曲线，就得到以谐振强度表示的时间-频率分布图。

作为本发明的一个优选例，以余弦函数作为平滑曲线来进行说明，但不 以此为限。

为了使曲线能够在每一个采样值点平滑，并且减少平滑后的波形产生的 额外高次谐波，选取余弦函数作为平滑曲线的基本形式，具体形式定为f(t)= a+bcos(pt)，p为数字信号的采样角频率的一半，a为两相邻采样值的算数平 均值，b为前一采样值与当前采样值的差的一半。这样相邻的二采样值点就 可以用1/2周余弦曲线相连，并且在连接点处的一阶导数连续。如图4中实 线所示，三角符号代表采样值，实线代表平滑曲线。

具体计算步骤如下：

i.以单一的a+bcos(pt)形式作为外强迫力，构造一个具有阻尼和强迫力的 谐振系统模型，该谐振系统的微分方程为：