[发明专利]一种Givens迭代的Prony低频振荡分析方法

摘要

本发明公开了一种Givens迭代的Prony低频振荡分析方法，采用Givens迭代法来减小信号模型的阶数，提高解算特征方程系数求解方程组的数值稳定性，包括以下步骤：读入滤波后的测量数据；用奇异值分解法计算信号模型实际阶数；用Givens迭代法计算特征方程系数；计算频率和衰减因子；计算幅值和初相。本发明通过Givens迭代的Prony法分析低频振荡，扩展Prony分析时的信号模型的阶数可取信号模型的实际阶数，有效避免了特征方程系数求解方程组的系数矩阵的奇异，提高了数值计算的稳定性。本发明在所取信号模型阶数为信号模型实际阶数的条件下，本发明比传统Prony法计算出的低频振荡信息的结果精度高。

主权项

1.一种Givens迭代的Prony低频振荡分析方法，其特征在于：包括以下步骤：A、读入滤波后的测量数据读入测量数据，并对测量数据滤波去掉高频信号及噪声；B、用奇异值分解法计算信号模型实际阶数根据测量数据形成扩展阶矩阵X_e，所取信号模型阶数p_e大于信号模型的实际阶数p，p_e取值为满足的整数，然后采用奇异值分解法计算扩展阶矩阵X_e的有效秩r，即得到信号模型的实际阶数p＝r；扩展阶矩阵Xe为：式中，pe为所取信号模型的阶数，x(i)表示第i个测量数据，i＝0、1、2、…、N‑2，N为测量数据总个数；C、用Givens迭代法计算特征方程系数根据计算出的信号模型的实际阶数p，列出特征方程系数求解方程组，由于特征方程系数求解方程组的方程数大于未知数个数，不能直接求解，用Givens变换对未知数进行最小二乘估计，计算出特征方程的系数a＝[a1 a2 … ap]T  (2)特征方程系数求解方程组为：式中：e(n)为估计值与实际测量数据x(n)的误差，定义如下：信号模型的阶数为p时的特征方程为：式中，a0＝1；为了提高计算精度，采用迭代法求解a；其步骤如下：C1、设迭代计数k＝0，收敛精度δ＝0.00001，最大迭代次数km＝30；C2、根据测量数据形成方程组的系数矩阵X和右端向量c；用最小二乘法对特征方程系数求解方程组(3)求解，就是把式(3)写成下式所示的不相容方程组，并令误差平方和最小，X和c分别为式(5)经过变化后得到的方程组的简写形式的系数矩阵和右端向量：Xa＝c  (6)式中，a＝[a1 a2 … ap]T  (8)c＝‑[x(p) x(p+1) … x(N‑1)]T  (9)C3、用Givens变换求解式(6)，得到a的初值a(0)；C4、把a(k)代入式(c′)(k)＝Xa(k)，得到(c′)(k)；C5、按照式Δc(k)＝c‑(c′)(k)，计算残差Δc(k)；C6、用Givens变换求解式XΔa(k)＝Δc(k)，得到Δa(k)；C7、求C8、判断Δamax是否小于δ，如果Δamax小于δ，转至步骤D；C9、令a(k+1)＝a(k)+Δa(k)；C10、令k＝k+1，如果k式中，arctan、ln、Im、Re分别为反正切函数、自然对数函数、取复数虚部函数、取复数实部函数，Δt为采样时间间隔；E、计算幅值和初相把求出的z代入以下公式(11)，并根据该式计算出b后，再计算信号分量的幅值和初相；式(11)的方程数也大于未知数个数，不能直接求解，采用Givens变换对未知数进行最小二乘估计；计算b的公式为：Zb＝x  (11)式中b＝[b1 b2 … bp]T  (13)x＝[x(0) x(1) … x(N‑1)]T  (14)幅值和初相计算公式为步骤C中所述的Givens变换方法如下：对式(6)两边进行Givens变换如下：GXa＝Gc  (16)式中，X为式(7)所示的(N‑p)×p阶矩阵，G为(N‑p)×(N‑p)阶矩阵，由一系列Givens变换矩阵相乘得到；G＝Gp,N‑p…Gp,p+1…Gi,N‑p…Gi,j…Gi,i+1…G1,N‑p…G1,2  (17)式中，Gij为对单位阵E的元素eii、eij、eji、ejj修改后得到的矩阵式中，cosθ和sinθ分别为 被变换矩阵X的元素的函数，为：式中，上标(n)表示前次Givens变换得到的新元素；由式(16)和式(17)可知，式(17)中各矩阵依次左乘X和c，对应的式(19)中cosθ和sinθ中X的元素都是进行前次Givens变换得到的新元素。