CN116187064A - 一种连续信号时间序列二阶导数的数值仿真方法 - Google Patents

Abstract

本发明涉及一个连续信号的观测时间序列二阶导数仿真领域，特别涉及一种连续信号时间序列二阶导数的数值仿真方法。本发明方法包括：对连续信号进行采样，获得连续信号的观测时间序列；选取观测时间序列中的多个采样点作为观测数据，并对观测数据的二阶导数进行离散近似；对二阶导数进行Z变换，并进一步化为实部和虚部的组合形式后再进行化简；在某一频段，对化简后的二阶导数进行泰勒展开；比较泰勒展开后的式子和二阶求导算子的傅里叶变换，得到二阶导数中的各系数与自由参数的关系式；通过对自由参数取值，获得不同的二阶导数中的各系数，以相应频段内各采样点的估计值与真实值之间差值的最小二乘平方和最小为约束条件，实现高精度数值仿真。

Description

一种连续信号时间序列二阶导数的数值仿真方法

技术领域

本发明涉及一个连续信号的观测时间序列二阶导数高精度仿真领域，特别涉及一种连续信号时间序列二阶导数的高精度数值仿真方法，具体涉及到给定一个连续信号的观测时间序列，通过插值算法获得连续信号的二阶导数离散近似，保证在低频区间，与原连续信号的二阶导数的线性功率谱密度差异很小。

背景技术

针对以空间引力波探测为代表的许多空间科学卫星高精度探测任务中，由于卫星研制成本极高，风险巨大，需在地面通过数值仿真对整个系统进行全链条的高精度建模与仿真，仿真得到科学家所关注的观测信号，以及工程信号。为了仿真获得与真实系统高精度匹配的信号，需根据科学任务特点，对仿真方法进行研究，使得仿真方法引入的仿真误差对科学信号的影响可以忽略。在这些仿真方法中，二阶导数的数值仿真算法是一个重要的待解决的问题。例如，在空间引力波探测任务全链条仿真过程中，传感器仿真得到的数据是距离数据，但是科学目标数据是加速度数据，需要对距离数据进行二次求导获得。考虑到在空间引力波探测任务中科学家比较关注的是低频段科学数据高精度仿真，那么怎样进行导数的离散数值仿真使得仿真方法在低频段科学数据的线性功率谱密度中的影响可以忽略就是一个具有重要的研究意义的课题。

因此，我们的问题是，若给定一个连续信号的观测时间序列，那么怎么从给定连续信号的观测时间序列出发，通过差值得到连续信号的二阶导数的离散近似，使得在低频区间，得到连续信号的二阶导数的离散近似值的线性功率谱密度与相应的真实值的线性功率谱密度差异很小。故本发明将聚焦于连续信号二阶导数的高精度数值仿真方法。

发明内容

本发明的目的在于克服目前对于连续信号二阶导数算子离散仿真的实现过程中在低频段仿真精度不高的问题，提出了一种基于多个采样点的包含一个自由参数的二阶导数算子的仿真方法。

为达到上述目的，本发明通过下述技术方案实现。

本发明提出了一种连续信号时间序列二阶导数的数值仿真方法，所述方法通过差值算法获得连续信号二阶导数的离散近似值，包括：

对连续信号进行采样，获得连续信号的观测时间序列；

选取观测时间序列中的多个采样点作为观测数据，并对观测数据的二阶导数进行离散近似；

对二阶导数进行Z变换，并进一步化为实部和虚部的组合形式后再进行化简；

在某一频段，对化简后的二阶导数进行泰勒展开；

比较泰勒展开后的式子和二阶求导算子的傅里叶变换，得到二阶导数中的各系数与自由参数的关系式；

通过对自由参数取值，获得不同的二阶导数中的各系数，实现对连续信号时间序列二阶导数的数值仿真。

作为上述技术方案的改进之一，所述方法还包括：以相应频段内各采样点的估计值与真实值之间差值的最小二乘平方和最小为约束条件，获取最优的二阶导数的一组系数，实现对连续信号时间序列二阶导数的高精度数值仿真。

作为上述技术方案的改进之一，所述对连续信号进行采样，获得连续信号的观测时间序列，具体包括：

设定采样间隔为ΔT，对连续信号进行采样，并记所有观测时刻t的集合为{t_k}，观测数据y的集合为{y_k}，采样次数k＝1,2,3...，获得连续信号的观测时间序列。

作为上述技术方案的改进之一，所述选取观测时间序列中的多个采样点作为观测数据，并对观测数据的二阶导数进行离散近似，具体包括：

在时刻t_k左、右两侧分别选择m个采样点作为观测数据；

对观测数据的二阶导数

进行离散近似:

其中，T为采样周期，a_1,j为二阶导数的各系数，j表示t_k时刻左右两侧第j个采样点，j＜0为左侧，j＞0为右侧，j为整数且-m≤j≤m，y_k+j表示t_k时刻左、右两侧的第j个采样点的观测值。

作为上述技术方案的改进之一，所述对二阶导数进行Z变换，并进一步化为实部和虚部的组合形式后再进行化简，其中，

对二阶导数

进行Z变换，得到：

其中，

表示对二阶导数/>

进行Z变换，Z[y_k]表示观测值y_k进行Z变换，z表示Z变换的算子。

进一步化为实部和虚部的组合形式，具体包括：

取z＝e^iΩ，代入

得到实部和虚部的组合形式：

其中，i为虚数符号，变量Ω＝ωT，ω＝2πf为角频率，f为观测频率，a_2,0表示二阶导数中的一个系数；

再进行化简，具体包括：

将

的实部取0，化简得到：

作为上述技术方案的改进之一，所述在低频段，对化简后的二阶导数进行泰勒展开，得到：

其中，l为泰勒展开的级数。

作为上述技术方案的改进之一，所述二阶求导算子的傅里叶变换，为：

作为上述技术方案的改进之一，当m＝2时，二阶导数中的各系数与自由参数的关系式，包括：

其中，α为自由参数，a_2,-2，a_2,-1，a₂₀，a₂₁和a₂₂均为二阶导数中的系数。

本发明与现有技术相比优点在于：

1、本发明的二阶导数算子估计方法，可以设计5项系数及以上二阶导数估计器。

2、本发明设计得到的二阶导数算子估计方法具有广泛适用性，适用于其他数学物理计算中的信号导数算子计算。

3、本发明的导数算子估计方法，可以通过改变α的值，改变初始条件，根据不同的仿真需求设计不同的二阶导数算子估计方案。

4、对导数算子的估计，可以看出对于关注的低频段内的二阶导数算子，本发明可以获得很好的估计。

附图说明

图1是二阶导数算子估计流程图；

图2是单参数二阶离散导数算子仿真算法示例图；

图3是二阶离散导数算子仿真算法1mHz～0,1Hz关注频段仿真最优结果图；

图4是二阶离散导数算子仿真算法1mHz～0,5Hz关注频段仿真最优结果图。

具体实施方式

下面结合附图和实施例对本发明的技术方案进行详细的说明。

对于连续信号二阶导数的数值仿真，仿真精度评判标准有很多。根据不同的仿真精度判断标准，可以给出不同的仿真算法。在本发明中，针对一个连续信号，已知已获得其观测时间序列，基于此观测时间序列，怎么通过差值算法构建连续信号本身的二阶导数的离散近似为仿真目标，采用在低频区间线性功率谱一致作为评判标准，即仿真算法获得的连续信号二阶导数的离散近似值的线性功率谱密度在低频区间与连续信号二阶导数的真实值的线性功率谱密度一致作为仿真精度标准。在本发明中，将基于原始连续信号的5个采样点的值，通过差值，获得其二阶导数的离散近似值。目前，针对连续信号二阶导数算法，基于5个采样点的仿真方法包含泰勒级数法和PI方法。但是这些仿真方法的求解精度已经固定，不能根据所研究的问题进行仿真精度的改善。在本发明中，将针对连续信号二阶导数仿真算法，给出包含一个自由参数的一类高精度二阶导数仿真算法。通过选取不同的参数，获得不同的仿真精度。同时，针对连续信号二阶导数算法，泰勒级数法和PI方法也是本发明二阶导数仿真算法的两个特例。因此，本发明获得的二阶导数仿真算法使用范围更广。

实施例

设在观测时刻t_k，获得一个连续信号的观测值记为y_k，这里下标k为自然数。并记所有观测时刻的集合为{t_k}，观测数据的集合为{y_k}。假设采样时间均匀，考虑采样间隔为ΔT，则观测时刻t_k和t_k-1之间关系为t_k＝t_k-1+(k-1)ΔT,k＝1,2,3...。在计算二阶导数时，可以在时刻t_k左右两侧分别选择m个的观测数据对二阶导数进行离散近似:

其中，a_2,j为二阶导数估计过程中的待求系数，本发明基于5个采样点法进行仿真方法研究，m取值为2。

将(1)进行Z变换，得到

为采样周期，f_s＝1/T为采样频率。

取z＝e^iΩ，代入(2)，化为实部和虚部的组合形式：

由于在频域上，二阶求导算子的傅里叶变换对应于

考虑到实际数据为实数，因此(3)中虚部为0，则有

a_2,-j＝a_2j,j＝1,...,m (5)

因此(3)可以化简为：

在低频处，对(6)进行泰勒展开，得到：

根据(4)和(7)，比较Ω的零次方和二次方系数，分别得到：

a₂₀+2(a₂₁+a₂₂)＝0, (8)

a₂₁+4a₂₂＝1. (9)

对于二阶导数，方程(8)和(9)给出三个系数满足两个方程，因此包含一个自由参数。若假设当Ω＝π时，有

其中α为一个实数，则有：

a₂₀-2a₂₁+2a₂₂＝α (11)

联立求解方程(8)，(9)和(11)，可以得到

以上结果给出了一个单参数二阶导数离散求解方法，为本发明的主要结果，我们称之为PI-α方法。

从(12)可知，当α＝-16/3时，得到泰勒级数方法：

当α＝0时，得到PI方法：

上述可以看出，泰勒级数法和PI方法只是本研究获得的二阶单参导数在参数取特定值时的特例。

以关注频段1mHz～0.1Hz内的离散点的估计值与真实信号值差值的最小二乘平方和最小为优化目标，得到的最优的一组系数为：

a_2,-2＝-0.083，a_2,-1＝1.333，a₂₀＝-2.50，a₂₁＝1.333，a₂₂＝-0.083。 (15)

该组系数与泰勒级数法参数相同，因此可看到泰勒级数法既可以是当参数α＝-16/3时得到的导数算子，也可以是从最小二乘法求得。对空间引力波探测任务，其探测频段为0.1mHz到1Hz，所以当我们关注频段放宽至1mHz～0.5Hz后，得到的最优的一组系数为：

a_2,-2＝-0.081，a_2,-1＝1.324，a₂₀＝-2.485，a₂₁＝1.324，a₂₂＝-0.081。 (16)

(15)和(16)说明本发明提出的方法可适用于不同的关注频段需求，应用更加广泛且更具针对性。

本发明的二阶导数算子估计流程如图1所示，下面结合图1说明本发明的流程。(图1所示的方案是两条路线，左侧路线是对传统仿真方法和本发明提出的仿真方法进行的一个仿真结果比较。右侧路线是对本发明仿真方法寻找最优解，包含但不仅仅包含左侧路线的几组仿真结果。)

1、确定想要使用的采样点数(2m+1)，并选取m的值。为了获得较好的导数算子计算结果，通常需要选择5个以上的采样点，才可以应用更多的系数约束条件，而上限没有限制，为了兼顾优化速度，选择m＝2。

2、确定泰勒展开的阶数l,通常需要选择3阶以上的系数，才可以应用更多的系数约束条件，而上限没有限制，为了优化计算速度，选择l＝10。

3、确定导数算子估计的采样周期T，为简化计算选择T＝1。

4、设定离散的频率矩阵freq，该矩阵为两个矩阵，在我们关注的低频段区间1mHz和0.1Hz之间采用均匀对数采样的方法取50个点，在0.1Hz到0.5Hz，间隔0.01Hz均匀采样。

5、确定离散频率对应的离散点弧度Omega。

6、确定方法想要使用的单参数α值，将其存入矩阵alpha。

7、根据参数α的值确定系数矩阵a2，并计算离散点的导数算子，并画图得到图2。

8、根据4中频域段，确定我们关注的低频段离散点。

9、确定二阶导数的最小二乘误差，确定约束条件，使用最小二乘法进行优化选择。

10、对最优化选择得到的系数代入导数估计矩阵，对离散二阶导数进行仿真，得到图3和图4。

对于传统的二阶导数仿真方法而言，他们的精度是固定的，而我们采取的方法相对于传统方法而言，可以实现传统方法，且根据不同的关注频段选择最优的，精度最高的系数。例如：在1mHz～0.1Hz频段，得到的最优系数为泰勒展开的系数，说明泰勒展开法在该频段精度更好。在1mHz～0.5Hz频段，得到该仿真方法的系数与其他方法都不一致，但是我们可以寻找到最优系数。可通过最小二乘的平方误差比较优劣性，见表1。

表1 1mHz～0.1Hz频段二阶离散化导数仿真精度

最后所应说明的是，以上实施例仅用以说明本发明的技术方案而非限制。尽管参照实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，对本发明的技术方案进行修改或者等同替换，都不脱离本发明技术方案的精神和范围，其均应涵盖在本发明的权利要求范围当中。

Claims (8)

1.一种连续信号时间序列二阶导数的数值仿真方法，所述方法通过差值算法获得连续信号二阶导数的离散近似值，包括：

对连续信号进行采样，获得连续信号的观测时间序列；

选取观测时间序列中的多个采样点作为观测数据，并对观测数据的二阶导数进行离散近似；

对二阶导数进行Z变换，并进一步化为实部和虚部的组合形式后再进行化简；

在某一频段，对化简后的二阶导数进行泰勒展开；

比较泰勒展开后的式子和二阶求导算子的傅里叶变换，得到二阶导数中的各系数与自由参数的关系式；

通过对自由参数取值，获得不同的二阶导数中的各系数，实现对连续信号时间序列二阶导数的数值仿真。

2.根据权利要求1所述的连续信号时间序列二阶导数的数值仿真方法，其特征在于，所述方法还包括：以相应频段内各采样点的估计值与真实值之间差值的最小二乘平方和最小为约束条件，获取最优的二阶导数的一组系数，实现对连续信号时间序列二阶导数的高精度数值仿真。

3.根据权利要求1所述的连续信号时间序列二阶导数的数值仿真方法，其特征在于，所述对连续信号进行采样，获得连续信号的观测时间序列，具体包括：

设定采样间隔为ΔT，对连续信号进行采样，并记所有观测时刻t的集合为{t_k}，观测数据y的集合为{y_k}，采样次数k＝1,2,3...，获得连续信号的观测时间序列。

4.根据权利要求3所述的连续信号时间序列二阶导数的数值仿真方法，其特征在于，所述选取观测时间序列中的多个采样点作为观测数据，并对观测数据的二阶导数进行离散近似，具体包括：

在时刻t_k左、右两侧分别选择m个采样点作为观测数据；

对观测数据的二阶导数

进行离散近似:

其中，T为采样周期，a_1,j为二阶导数的各系数，j表示t_k时刻左右两侧第j个采样点，j＜0为左侧，j＞0为右侧，j为整数且-m≤j≤m，y_k+j表示t_k时刻左、右两侧的第j个采样点的观测值。

5.根据权利要求4所述的连续信号时间序列二阶导数的数值仿真方法，其特征在于，所述对二阶导数进行Z变换，并进一步化为实部和虚部的组合形式后再进行化简，其中，

对二阶导数

进行Z变换，得到：

其中，

表示对二阶导数/>

进行Z变换，Z[y_k]表示观测值y_k进行Z变换，z表示Z变换的算子；

进一步化为实部和虚部的组合形式，具体包括：

取z＝e^iΩ，代入

得到实部和虚部的组合形式：

其中，i为虚数符号，变量Ω＝ωT，ω＝2πf为角频率，f为观测频率，a_2,0表示二阶导数中的一个系数；

再进行化简，具体包括：

将

的实部取0，化简得到：

6.根据权利要求5所述的连续信号时间序列二阶导数的数值仿真方法，其特征在于，所述在低频段，对化简后的二阶导数进行泰勒展开，得到：

其中，l为泰勒展开的级数。

7.根据权利要求6所述的连续信号时间序列二阶导数的数值仿真方法，其特征在于，所述二阶求导算子的傅里叶变换，为：

8.根据权利要求7所述的连续信号时间序列二阶导数的数值仿真方法，其特征在于，当m＝2时，二阶导数中的各系数与自由参数的关系式，包括：

其中，α为自由参数，a_2,-2，a_2,-1，a₂₀，a₂₁和a₂₂均为二阶导数中的系数。

Images