CN102105771A - 数据频率分析的方法和装置 - Google Patents

Abstract

数据频率分析方法包括：来自第一传感器的信号的输入步骤(310)；来自至少一个第二传感器的信号的输入步骤(315)，每个第二传感器被安置在第一传感器附近以便来自每个第二传感器的信号与来自第一传感器的信号高度相关；对于每个传感器估计基于第一传感器和每个第二传感器的所有信号建立的传递函数或模型的估计步骤；和利用每个经估计的模型抽取***结构特性的抽取步骤(320)。

Description

数据频率分析的方法和装置

技术领域

本发明涉及数据频率分析的方法和装置。特别地，本发明应用于飞机飞行包线打开试验(essais d′ouverture)的数据分析。

本发明尤其在航空领域中应用于飞行控制，例如结构的振动模态的飞行分析和控制；在汽车领域中，应用于运输工具的振动的研究和控制；在电动力学方面(发电机的控制)，尤其在核领域中，应用于反应堆芯的振动监测；应用于力学(活动部件的研究和控制)中、地震学(石油勘探中所使用的信号的研究)中和在动物学(动物发出的声音的研究)中。

本发明的目标是，在测试过程(在飞机的情况下，是在飞行过程中)中估计运输工具的特征，特别是共振频率和频谱特征。换言之，涉及从来自安装在运输工具上的传感器的非常大量的信息，非常快速地、甚至实时地抽取恰当的信号。

背景技术

如此确定的所有***特征允许设计者改善***的结构以便增加其舒适性、其飞行包线、其消耗等。

希望分析的信号由向其添加有一个或多个正弦信号的噪声组成，所述正弦信号的频率和振幅可随时间而变化。涉及实时地估计这些频率和这些振幅。

在作为研究空气动力、惯性力和弹性力之间相互作用的学科的空气弹性力学的领域中，颤振(flottement)现象是非常危险的飞机结构(机翼、机身、尾翼...)的震荡性不稳定，因为它可能通过损坏结构直至断裂而影响该结构的完整性。正是两种或多种具有不同性质的飞机结构运动的组合，以适当相位差允许空气动力学力向***提供能量。因而，将稳定现象转化为能量不再被耗散的不稳定现象：这是模态耦合(couplage de modes)。最经常被提及的例子是机翼的颤振，其产生自弯曲和扭曲的模态之间的耦合，这些模态在正交相位(quadrature de phase)导致具有与移动相同的方向的升力的空气动力学力，并因此导致发散性振荡。

在表征颤振中几个参数是有影响的：质量、刚度、结构的形状以及操作状态例如速度。为了预防该现象，飞机制造者应当研究并证明，如果存在的话，其出现阈值位于最大运行速度以外(增加15％)。首先进行风箱试验，并通过飞机结构的地面振动试验来补充。理论研究因而使得能够定义没有颤振的范围，由此将可以通过“激励”飞机结构而逐步实现飞行包线全面打开。

使用模态参数的识别方法以便准实时地抽取频率值和阻尼值，并研究其在飞行包线中的演变。来自颤振试验的时间数据的分析是复杂的：数据被噪声污染并需要通过信号处理(尤其是滤波和亚采样)的修整(mise en forme)。目前，使用数个传感器以“全面地”并自动地抽取空气弹性力学结构的模态参数。

已知文献WO03005229，其描述了对来自传感器的信号的频率分析***。然而，该分析的解析方案是有限的。

发明内容

本发明旨在克服这些缺点。

因此，根据第一方面，本发明旨在一种***频率分析方法，其特征在于，所述方法包括：

-来自第一传感器的信号的输入步骤，

-来自至少一个第二传感器的信号的输入步骤，将每个第二传感器安置在第一传感器附近以便来自每个第二传感器的信号与来自第一传感器的信号高度关联，

-对于每个传感器估计基于第一和每个第二传感器的所有信号建立的传递函数或模型的估计步骤，和

-利用每个经估计的模型抽取***的结构特性的抽取步骤。

模型将来自传感器的信号视作通过白噪声激发的滤波器的输出(sortie)。结构特性包括例如频谱特性、频率、振幅、初始相、阻尼、模态。

因此，将表示结构模态的模型视作是线性的。借助于本发明的实施，可以通过进行频率/阻尼耦合的在线监测实现实时处理。由于能够实时地掌握***的结构特性，因而本发明使得能够实时地确保***例如飞机的性能是令人满意的。因此，通过满足节省时间和因此减少成本的不断增长的约束，改善了所使用的分析方法。

根据特别的特征，在估计步骤中，将来自传感器的信号视作多项式。由于这些规定，信号的表示是密集的，因为多项式的系数的数目比使用的信号采样的数目小得多。

根据特别的特征，在估计步骤的过程中，解析具有随时间和随传感器之间的空间缓慢变化的系数的线性递推方程。

根据特别的特征，所述估计步骤包括：

-关于时间、阶数和传感器的空间的递推适应性建模步骤，和

-根据适应性建模步骤的结果对于每个阶数的模态估计步骤。

根据特别的特征，来自传感器的信号的每个输入步骤包括先于适应性建模步骤的来自传感器的信号的噪声水平的实时降低步骤。

根据特别的特征，模态估计步骤包括根据适应性建模步骤的结果抽取模型参数的抽取步骤。

根据特别的特征，模型参数的抽取步骤包括：具有N阶和等于传感器数目的维数的多项式矩阵的求逆步骤。

根据特别的特征，模态估计步骤适合于提供构成一组冗余信息的每个模型的参数，该组冗余信息使得能够减小所估计的模态的方差。

根据特别的特征，适应性建模步骤进行ARMA“自回归移动平均”类型的建模。

根据特别的特征，所述ARMA类型的建模在每个时刻对于每个传感器和对于所考虑的所有阶数进行。

根据特别的特征，估计步骤包括多项式矩阵的求逆步骤，所述矩阵是对称互谱矩阵，其在其主对角线上表示每个传感器的功率谱密度，换言之，互谱。

根据特别的特征，适应性建模步骤包括下述步骤的关于时间、关于模型的阶数N＝[1，2，...，Nmax])和传感器的数目的递推：

-前和后线性预测误差向量的计算，

-后和前部分相关矩阵的计算，

-前和后线性预测误差的协方差矩阵的计算，

-后和前线性预测误差的幂矩阵的计算，

-增益向量的直接计算

${\underset{\‾}{\mathrm{\θ}}}_{a,n}^N=\left\{{\left(\underset{k=0}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\λ}}^{n-k}{\underset{\‾}{\mathrm{\ξ}}}_k^N\underset{\‾}{\mathrm{\ξ}}_k^N{}_t\right)}^{-1}{\underset{\‾}{\mathrm{\ξ}}}_n^N\right\}{\mathrm{\α}}_n^N$

α为标量，λ为遗忘因子，并且

-通过知道

递推计算向量

和

-矩阵A_k的计算，A_k表示模型，k＝1至N。

根据特别的特征，本发明的目标方法，例如上面所简要叙述的，包括通过应用下述两个约束中的一个进行来自模态估计步骤的模态的分类步骤：

-来自相同模型的的每类仅一个模态，和

-所有估计具有相同权重而不受估计来源的影响。

根据特别的特征，所述信号代表飞机结构的加速。

根据第二方面，本发明旨在可加载于信息***中的计算机程序，所述程序包含使得能够施行本发明的如上面所简要叙述的目标方法的指令。

该程序的特别的优点、目标和特征与本发明的如上面所简要叙述的目标方法的优点、目标和特征是相似的，在本文中不再提醒。

附图说明

本发明的其他优点、目标和特征将从下述的参照附图以说明而无限制性的目的进行的描述中体现出来，附图中：

-图1示意地表示具有适合于执行本发明的目标方法的装置的飞机，

-图2表示来自图1中所示的装置的两个传感器的信号，

-图3以逻辑图的形式表示在本发明的目标方法的第一实施方式中应用的一些步骤，

-图4以逻辑图的形式表示本发明的目标方法的第二实施方式中应用的步骤，

-图5表示在图4中所示的步骤之一过程中应用的滤波器的布置，

-图6示意地表示在噪声降低***的一实施方式中使用的相继函数，

-图7示意地表示在每一时刻来自传感器的采样，其构成图4所示实施方式的算法的输入值，

-图8示意地表示在图4中所示的第二实施方式中使用的递推，

-图9示意地表示“云动力学”类型的非监督分类方法的分类演变，和

-图10给出了在图4所示第二实施方式中使用的选通窗(

de validation)的说明。

具体实施方式

在图1中，观察到飞机105，其配有：在机翼120前部的两个彼此邻近的传感器110和115；和在机翼120后部的两个彼此邻近的传感器125和130。

以说明性的目的，在图1中仅表示了两对邻近的传感器。然而，注意到，在本发明的实际实施中，使用多于两对的传感器。

术语“邻近”这里指接收相互间高度关联的信号的传感器。

邻近的传感器大致接收随时间错移并有区别地减弱的但按照基本线性的传递函数的相同振动。所讨论的传感器为例如加速计。

在图2中观察到，来自一对传感器的第一传感器的信号205具有噪声210和两个峰215及220，并且来自相同一对传感器的第二传感器的信号255具有噪声260和两个峰265及270。峰265对应于随时间错移和减弱的峰215。峰220对应于随时间错移和减弱的峰270。

如容易理解的是，通过使用一些由至少两个相近传感器构成的组(在本文中，是一些对)，本发明使得能够分析飞机的结构特性。这些结构特性包括例如频谱特性、频率、振幅、初始相、阻尼、模态。

在图3中，观察到，在第一实施方式中，本发明的频率分析目标方法首先包括使一些由多个传感器构成的组在经受振动的机械***的结构上定位的步骤305。在每个传感器组中，将至少一个所谓“第二”传感器定位在所谓“第一”传感器的附近。

在机械***的运行中，进行来自一个传感器组的第一传感器的信号的输入步骤310和来自该相同传感器组的至少一个第二传感器的信号的输入步骤315。每个来自传感器的信号的输入步骤包括来自传感器的信号的噪声水平的降低步骤。该噪声降低可以以已知的方法或关于一个向量地逐个传感器地来进行，所述向量对于其每个坐标都具有来自传感器的信号。优选地，该去噪函数通过基于小波的分解(Stéphane Mallat算法)来确保。在步骤310和315中，例如，对于寻求其小于16Hz的固有频率的结构，取样来自具有十分高的频率、例如在图6中的256Hz的传感器的信号。注意到小波的使用允许简单和快速的处理。

在步骤320中，抽取通过来自第一传感器和每个第二传感器的信号的处理产生的传递函数。因此，将表示结构模态的模型视作是线性的，并且进行频率/阻尼耦合的在线监测而不考虑输入即注入的激励。

在步骤320中，在考虑来自不同传感器组的不同传感器的信号之间的所述信息的相和值的改变后，将来自传感器的信号视作多项式并且抽取来自传感器的信号的信息。

在步骤320中，进行具有随时间缓慢变化的系数的线性递推方程的解析步骤325，以便在足够稳定的时间间隔内能够估计模型。

步骤325包括：

-关于时间(在每个时刻)、阶数(对于每个所考虑的阶数)和传感器(对于每个传感器)的空间的ARMA类型的可递推适应性建模步骤330，和

-根据适应性建模步骤的结果对于每个阶数的模态估计步骤350。

步骤330实现ARMA类型的参数型的建模。将来自传感器的采样汇集为一个向量，该向量的分量数目为所考虑的传感器的数目(参见图7)。例如，如果有四个传感器，则所考虑的向量为四维的。更一般地，在说明书的下文中，称“p”为一组邻近传感器的传感器数目。

步骤330进行所谓“随时间递推的”建模，因为它使用所获得的最近估计值以更新其参数。因此，在该实施方式中，本发明使用两个接续时刻之间的关系，因为它们被视作相关的和相干的。在一些实施方式中，正是前面的这两个时刻用作时间递推。

建模的最佳阶数N不必是已知的。优选地，为了实时地获得结果，不确定该最佳阶数N。相反地，进行递推直至足够大的阶数Nmax。如此，获得足够大的一组模型，以便所寻求的信息(结构模态)就在其中被表示。

在图3所示的实施方式中，步骤330包括对称互谱矩阵的处理步骤335，所述对称互谱矩阵在其主对角线上表示每个传感器的功率谱密度并且，以额外对角线(extra diagonaux)的术语表示，即互谱。

步骤335包括以在时刻n对于阶数N＝0为起始条件的关于时间的、并具有关于模型阶数N＝[1，2，...Nmax]递推的递推步骤340，其包括下述步骤：

-前和后线性预测的误差向量的计算，

-后和前部分相关矩阵的计算，

-前和后线性预测误差的协方差矩阵的计算，

-前线性预测误差的幂矩阵的计算，

-增益向量的直接计算，其中维数是传感器的数目： ${\underset{\‾}{\mathrm{\θ}}}_{a,n}^N=\left\{{\left(\underset{k=0}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\λ}}^{n-k}{\underset{\‾}{\mathrm{\ξ}}}_k^N\underset{\‾}{\mathrm{\ξ}}_k^N{}_t\right)}^{-1}{\underset{\‾}{\mathrm{\ξ}}}_n^N\right\}{\mathrm{\α}}_n^N$

α为标量，λ为遗忘因子，并且

-通过利用已知

递推计算向量

和

-矩阵A_k的计算，A_k表示模型，k＝1至N。

基于步骤330结束时所获得的所有模型，进行分类以抽取结构模态。将模型视作多项式的比(所谓“ARMA”模型，就是“AutoRegressive Moving Average”(“自回归移动平均”)的首字母缩写，相应法文术语为moyenne mobile autorégressive或autorégressive àmoyenne ajustée)。为了抽取表示结构模态的参数值，使协方差矩阵Ea和Eb对于每个时刻和每个阶数最小化。为此，计算向量θ_a和θ_b，它们分别表示矩阵Ea的逆矩阵与向量ε_a之积和矩阵Eb的逆矩阵与向量ε_b之积。

为了实时地进行这些计算，对于所有1至Nmax的值，直接计算向量θ_a，n ^N。然后，通过知道θ_a，n ^N递推计算θ_b，n ^N。

矩阵A_k(维数p的平方)——对于范围1至Nmax的所有阶数在时刻n待估计的模型的系数矩阵——的数目等于Nmax x(Nmax-1)/2，即对于N＝15为105。

结构模态的分类原则为考虑如果在给定阶数的参数是恰当的，则在高一级的阶数找到该参数。参数的分类是不受监督的并且在于寻求来自不同模型的相似目标。如下面所述的，对于该分类，使用轨线构建(constructions de trajectoire)。

步骤350提供构成一组冗余信息的每个模型的参数，该组冗余信息使得能够减少所估计的模态的方差。步骤350包括步骤355，步骤355根据适应性建模步骤的结果抽取模型参数。步骤355包括阶数N且维数等于传感器数目的多项式矩阵的求逆步骤360。

步骤360包括Cholesky分解的步骤365。

在其第二特别的实施方式中，本发明的目标方法遵循按照参照图4所描述的实时过程，该过程专用于飞行包线打开试验的数据分析。该实施方式使得能够在下一个出现前处理每个信息，而不考虑注入结构的激发。

总之，该方法的特征满足了序言中所陈述的安全和减少成本的要求，并且使得能够通过提供性能更好的模型分析改善飞行包线打开过程。

待分析的信号是在飞机的初级结构上进行的加速测量值。

在下述段落中，详细描述参照图4中所述的加速类型的p个信号的实时分析过程的每个运算。

分析方法的第一步骤405在于：对来自同一组邻近传感器的p个传感器的p个信号进行“去噪”，例如借助于S.Mallat提出的金字塔算法，通过使用从小波规格化正交的基数。该算法追溯至1983年的Burt和Adelson的建立在图像的显像和压缩上的工作。该具有大的使用简单性的算法，具有与待处理的采样数目成比例的计算量。

处理的非线性特点使实施复杂化了，这更加是由于滤波器工作台不是随机(causal)的。将滤波器以在图5中所示的方式进行布置，同时使用小波。

Mallat算法是完全实时的，它普及于在每个时刻下同时处理p个采样。保留的原则建立在信号的超取样(sur-échantillonnage)上以便允许通过分析和合成进行“去噪”。

因此，在步骤405中，首先进行信号的子带分解的步骤406。在该分析或分解的过程中，将通过连续信号取样所形成的组列首先视作在与离散化相关的某个范围下该信号的近似(事实上采样对应于最精细的范围)。

通过约定，该近似范围相应于j＝0。因此，从属于子空间V₀的数据开始：取样组列{x₀，x₁，...，x_k，...}＝<f|φ_0，k>，k∈Z因而构成希望分析的数据全体。近似子空间之间的关系为

因此，只需分解离散信号即取样组列在两个子空间V₁和W₁上，以便有具有解析度(résolution)2^-l即<f|ψ_1，k>和<f|φ_1，k>的取样组列。下面的关系式给出了对于两个相继解析的递推：

$<f|{\mathrm{\&phi;}}_{j,k}>=\underset{n\&Element;Z}{\mathrm{\&Sigma;}}{h}_{n-2k}<f|{\mathrm{\&phi;}}_{j-1,n}>$ 和 $<f|{\mathrm{\&psi;}}_{j,k}>=\underset{n\&Element;Z}{\mathrm{\&Sigma;}}{g}_{n-2k}<f|{\mathrm{\&phi;}}_{j-1,n}>$

这些方程式分别表示卷积g_-n*<f|φ_j-1，n>和h_-n*<f|φ_j-1，n>，这些卷积随后每两个进行抽选(décimation par deux)。因此，通过随后实施抽选(décimation)的由滤波器

和的简单滤波运算基于在范围j-1所获得的那些来计算在范围j下的近似系数<f|φ_j，k>和细节系数<f|ψ_j，k>。

然后，进行分解的系数的阈值确定步骤407，在其过程中，只保留第一子带的系数。

随后，进行合成或重构的步骤408，在此过程中，基于已知在近似子空间内的投影的信号重构是以下述方式进行的：

即： $<f|{\mathrm{\&phi;}}_{j-1,n}>=\underset{k\&Element;Z}{\mathrm{\&Sigma;}}{h}_{n-2k}<f|{\mathrm{\&phi;}}_{j,k}>+\underset{k\&Element;Z}{\mathrm{\&Sigma;}}{g}_{n-2k}<f|{\mathrm{\&psi;}}_{j,k}>$

重构是前一个的对偶运算。它通过先于对近似系数和细节系数的插值(interpolation)的数值滤波来获得，所述近似系数和细节系数得自分解。

因此，“去噪器”的一般结构是在图6中所示的那种。注意到，所考虑的最大频率为取样频率f的一半。在未加工信号中，所寻求的信号位于最低的频率中。在从0至f/16、从f/16至f/8、从f/8至f/4和从f/4中f/2范围的频带内的信号分解后，进行通过将小波的系数置于零的阈值确定。然后，进行合成以提供包含所寻求信号的经去噪的信号。

然后，在步骤410中，进行ARMA(自回归移动平均)类型的适应性建模。该信号建模步骤是参数类型的，即使得能够获得所研究的信号的频谱分析的类型的，并且是在时间、阶数和关于传感器空间可递推的。

步骤410包括确定前线性预测向量和后线性预测向量的步骤411。扩展到向量情形，每次预测向量确定在于，以表达与采样的N个最近向量的线性组合，在图7中举例说明了表示法，其表示对p个传感器的空间的向量形式的前线性预测。

因此，涉及直接在时间的空间中使信号建模。如果x _n表示采样组列，则模型为：

${\underset{\&OverBar;}{\overset{^}{x}}}_n=\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{A}_k{\underset{\&OverBar;}{x}}_{n-k}$

注意，向量x _n对于分量而言具有p个传感器的当前的p个采样。在该表达式中，A_k为相应于传感器数目的维数p的矩阵，N为建模的阶数。

标量方法的推广给出了误差及其预测的新的写入：

${\underset{\&OverBar;}{e}}_n={\underset{\&OverBar;}{x}}_n-{\underset{\&OverBar;}{\overset{^}{x}}}_n\&DoubleRightArrow;{\underset{\&OverBar;}{e}}_n={\underset{\&OverBar;}{x}}_n-\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{A}_k{\underset{\&OverBar;}{x}}_{n-k}---\left(1\right)$

即： $\underset{\&OverBar;}{e}\left(z\right)=(1-\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{A}_k{z}^{-k})\underset{\&OverBar;}{x}\left(z\right)=\mathrm{\&Phi;}\left(z\right)\underset{\&OverBar;}{x}\left(z\right)---\left(2\right)$

从外部角度看，e _n是通过采样向量组列x _n所激发的滤波器RIF的输出。线性化特性允许逆过程：x _n因而作为由e _n激发的滤波器的输出而出现。该滤波器通过多项式矩阵Φ(z)求逆而获得，它是稳定的并且是RII(“àréponse impulsionnelle infinie(无限脉冲响应)”的首字母缩写)类型的。为了确定后线性预测，x _n-N的估值

表达为向量{x _n-N+k}k＝1，...，N的线性组合。特性与就前线性预测所获得的特性类似。

然后，在步骤412中，计算协方差矩阵、后和前部分相关矩阵以及后和前预测误差的幂矩阵。

然后，将先验和后验的后和前线性预测误差的协方差矩阵最小化，该矩阵基于以前定义的后和前线性预测概念构建。将后验的前线性预测误差的协方差矩阵定义为：

$E_{a,n}^N=\underset{k=1}{\overset{n}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}^{n-k}{{\underset{\&OverBar;}{e}}_k}^t{\underset{\&OverBar;}{e}}_k=\underset{k=1}{\overset{n}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}^{n-k}{({\underset{\&OverBar;}{x}}_k-A_n^N{}_t\left({\underset{\&OverBar;}{x}}_{k-1}^N\right))}^t({\underset{\&OverBar;}{x}}_k-A_n^N{}_t\left({\underset{\&OverBar;}{x}}_{k-1}^N\right))$

其中将0＜＜λ≤1，称为遗忘因子或适应因子。

在步骤415中，对于N的所有值，获得阶数N的模型的系数矩阵。因此，在步骤416中，最佳向量

在协方差矩阵

为最小时获得，即：

$\frac{\&PartialD;E_{a,n}^N}{\&PartialD;A_n^N}=0$

在步骤417中，基于估计先验的后线性预测误差的协方差矩阵所考虑的相似的推理给出：

$E_{b,n}^N=\underset{k=1}{\overset{n}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}^{n-k}{({\underset{\&OverBar;}{x}}_{k-N}-B_n^N{}_t\left({\underset{\&OverBar;}{x}}_k^N\right))}^t({\underset{\&OverBar;}{x}}_{k-N}-B_n^N{}_t\left({\underset{\&OverBar;}{x}}_k^N\right))$

在多传感器的范围下，名为“ESA”的算法进行下面定义的两个向量的计算：

${\underset{\&OverBar;}{\mathrm{\&theta;}}}_{b,n}^N={\left(E_{b,n}^N\right)}^{-1}{\underset{\&OverBar;}{\mathrm{\&epsiv;}}}_{a,n}^N$

${\underset{\&OverBar;}{\mathrm{\&theta;}}}_{a,n+1}^N={\left(E_{a,n+1}^N\right)}^{-1}{\underset{\&OverBar;}{\mathrm{\&epsiv;}}}_{a,n+1}^N$

通过“ESA”算法直接计算的向量

是向量和矩阵之积。矩阵是向量

的并向量积的直至时刻n+1的加权和，并且向量是在时刻n+1的该相同向量即

向量

通过利用

的递推来获得。

然后，在步骤420中，进行逆问题求解，即模型参数的抽取。利用以上描述的涵盖多传感器情况的模型方程(1)和(2)：

${\underset{\&OverBar;}{e}}_n={\underset{\&OverBar;}{x}}_n-{\underset{\&OverBar;}{\overset{^}{x}}}_n\&DoubleRightArrow;{\underset{\&OverBar;}{e}}_n={\underset{\&OverBar;}{x}}_n-\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{A}_k{\underset{\&OverBar;}{x}}_{n-k}---\left(1\right)$

即： $\underset{\&OverBar;}{e}\left(z\right)=(1-\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{A}_k{z}^{-k})\underset{\&OverBar;}{x}\left(z\right)=\mathrm{\&Phi;}\left(z\right)\underset{\&OverBar;}{x}\left(z\right)---\left(2\right)$

其中

其即为p×p维的阶数N＝[1，Nmax]的z的多项式矩阵，p为所分析的传感器的数目。

e(z)为未知的误差向量，相反地，其协方差矩阵已被估算出。因此为了得出它，计算：

e(z) ^t ′e(z)＝Φ(z)x(z) ^t x(z)^tΦ(z)

或互谱矩阵：

x(z) ^t x(z)^t＝Φ(z)^-1 e(z) ^t ′e(z)Φ(z)

该矩阵是对称的，每个传感器的功率谱密度在其主对角线上显现。额外对角线术语是互谱。模型参数的抽取在于，在每个时刻和对于每个阶数进行阶数N和维数p的多项式矩阵的求逆：Φ(z)。

为此，矩阵e(z) ^t ′e(z)＝L^tL的Cholesky分解给出：

x(z) ^t x(z)＝Φ(z)^-1L(z)^tL(z)^tΦ(z)^-1＝(Φ(z)^-1L(z))^t(Φ(z)^-1L(z))

＝(Λ(z)^-1U(z)L(z))^t(Λ(z)^-1U(z)L(z))

其中

$\underset{\&OverBar;}{x}\left(z\right)\underset{\&OverBar;}{x_t}\left(z\right)=\mathrm{\&Lambda;}{\left(z\right)}^{-1t}\mathrm{\&Lambda;}{\left(z\right)}^{-1}\left[\begin{array}{cc}{\left|\right|F_{\mathrm{1,1}}\left|\right|}^2& {\left|\right|F_{p,1}\left|\right|}^2\\ .& \\ .& \\ .& \\ .& \\ {\left|\right|F_{p,1}\left|\right|}^2& {\left|\right|F_{p,p}\left|\right|}^2\end{array}\right]$

多项式F表示不同传递函数的分子，分母(实际上唯一的分母，根据重合理论)为本征值(valeurs propres)。

假定 $F_{i,j}^N\left(z\right)=\underset{q=0}{\overset{N}{\mathrm{\&Sigma;}}}{a}_q^i{z}^{-q}\&DoubleRightArrow;F_{i,j}^N{\left(z\right)}^t{F}_{i,j}^N\left(z\right)=\underset{n=0}{\overset{N}{\mathrm{\&Sigma;}}}\underset{m=n}{\overset{n-N}{\mathrm{\&Sigma;}}}{a}_{n-m}{a}_n{z}^{-m}$

$F_{i,j}^N{\left(z\right)}^t{F}_{i,j}^N\left(z\right)=\left[\begin{array}{cccccccccc}{a}_0{a}_N& a_1{a}_N& a_2{a}_N& .& .& a_N{a}_N& 0& .& .& 0\\ & a_0{a}_{N-1}& & & & .& & & & \\ & & a_0{a}_{N-2}& & & .& & & & \\ & & & & & a_0{a}_0& & & & \\ 0& .& .& .& .& 0& a_0{a}_0& a_0{a}_0& & a_0{a}_N\end{array}\right]\left[\begin{array}{c}{z}^N\\ z^{N-1}\\ .\\ .\\ .\\ .\\ .\\ z^{-N}\end{array}\right]$

分子的系数为互谱矩阵(没有符号)的第N+1列的N+1个元素的平方根。

适应性建模的处理可以如在图8中那样表示。在该多传感器建模处理的示意图中，可见三个递推805、810和815的交错排列。

然后，在步骤425中，进行模态的分类。在每个时刻，估计{1:p:N^＊p}维的N个模型，其中p为传感器数。每个模型的参数(频率值和阻尼值)的抽取是使得能够减小估计的方差的一组冗余信息。目前力求在下述约束下按照某一标准对这些来自建模的数据进行分类：

-每类仅一个来自相同模型的模态，

-所有估计独立于估计来源都具有相同权重。

先验认识的缺乏(特别是密度)致使优先采用“云动力学”类型的非监督分类方法，该方法在于找到自然的(隐含的)类别以便集中未经划分的数据。

在后面图9中所示意的该方法根据下述算法令人满意地满足了所有要求：

-在步骤426中，选择按K个类别的初始划分，

-在步骤427中，对于每个数据***地寻求最好的类别；计算数据至重心的距离并且将元素分配给中心离其最近的类别(通过使用例如欧几里得距离或Kullback-Leibler距离)，和

-在步骤428中，进行“发送”类别和“接受”类别的重心的更新。

然后，返回步骤427和428直至收敛。

在步骤430中，进行模态的轨线构建。

问题在于，实时监测相应于其数目随时间演变的模态的频率的所有靶标的轨线。算法结构是围绕每个监测轨迹一个卡尔曼(Kalman)滤波器构建的。当通过模型提供一组“测量值”时，力求将其与现有轨迹相关联。在这里目的是：在所接受的测量值中选择可出自基于此预测测量值的靶标的那些测量值。经常使用的原则在于围绕所做的预测来定义窗，通常称为“选通(gating)”。

处理的一般结构基于下述原则：

-利用已知的轨线预测在时刻n+1(知道n)的状态，

-测量-靶标的联合在于：将测量值与利用已知轨线预测的那些测量值进行比较。该处理应当使得能够，不仅维持已经存在的轨迹而且还初设新的轨迹，并且必要时消除对应于已经偏离观察空间的靶标的那些轨迹。轨迹***性能取决的正是这些函数的性质。

-卡尔曼增益以及协方差矩阵的具有状态更新的滤波。

在该过程中，卡尔曼滤波器使得能够追踪多个靶标，其中预测占据基本的作用。对于每个靶标，它按照协方差最小提供状态的经过滤的估计，预测状态并允许“选通”的计算。

解由预测方程和滤波方程的两个方程***的全体构成，即：

对于每个靶标，选通窗使得能够选择可属于靶标的测量值。原则是，围绕预测的测量值确定一个区、观察空间中的一个体积。该区的大小通过预测的测量值(在这种情况下，高斯法)的统计学特性来确定。一般地，该体积的“维数”应当合理地进行选择。这是因为，测量值分类和概率事实上取决于它，以便来自靶标的测量值在界定该体积的表面之内。

图10给出了选通窗的举例说明。

测量-标靶联合技术是靶标监测过程的中心部分。存在许多技术，但其中的一些技术并不处理轨迹的出现和消失。因此，应该考虑用于实现该管理的补充机制。一个简单的方式在于采用下述法则：

-法则1：将任何不与现有轨迹相关的测量值视作新轨迹的起始。

-法则2：如果至少相继的Nd个测量值已与其相关，则确认(检测到)轨迹。

-法则3：如果至少连续的Nl个测量值本身与其不相关，则将轨迹视作消失的。

匈牙利方法使得能够借助于下述特别的解析方法来解决估计-测量的分配问题以寻求最小成本：

假定需分配m项任务给m个人，而C即是分配成本矩阵。任一个分配由记为(1，x)，(2，y)，..(k，t)，...(m，u)的m对定义，并且(x，y，...u)：{1，2，...，m}的排列。花费或总成本对应于特别的分配A：

F(x，y，..t，...，u)＝C1，x+C2，y+...+Ck，t+...+Cm，u

因而问题在于确定(x，y，..t，..，u)，以使F最小。

为了应用本发明的目标方法，在本发明的一个优选的实施方式中，提出一般用途计算机，在该计算机中配有可加载的计算机程序，所述程序包含执行上面详细描述的步骤和算法的指令。

Claims (15)

1.数据频率分析方法，其特征在于，所述方法包括：

-来自第一传感器的信号的输入步骤(310)，

-来自至少一个第二传感器的信号的输入步骤(315)，每个第二传感器定位在所述第一传感器附近以便来自每个第二传感器的信号与来自所述第一传感器的信号高度相关，

-对于每个传感器估计基于所述第一传感器和每个第二传感器的所有信号建立的传递函数或模型的估计步骤，和

-利用每个估计的模型抽取***结构特性的抽取步骤(320)。

2.根据权利要求1所述的方法，其特征在于，在所述抽取步骤(320)的过程中，将来自所述传感器的信号视为多项式。

3.根据权利要求1或2所述的方法，其特征在于，在所述抽取步骤(320)的过程中，解析具有随时间和随传感器之间的空间缓慢变化的系数的线性递推方程。

4.根据权利要求1至3中任一项所述的方法，其特征在于，所述抽取步骤(320)包括：

-关于时间、阶数和传感器的空间的递推适应性建模步骤，和

-根据所述适应性建模步骤对于每个阶数的模态估计步骤。

5.根据权利要求4所述的方法，其特征在于，每个来自传感器的信号的输入步骤包括在所述适应性建模步骤之前的来自传感器的信号的噪声水平的实时降低步骤。

6.根据权利要求4或5所述的方法，其特征在于，所述模态估计步骤包括根据所述适应性建模步骤的结果抽取模型参数的抽取步骤。

7.根据权利要求6所述的方法，其特征在于，所述的模型参数的抽取步骤包括阶数为N和维数等于传感器数目的多项式矩阵的求逆步骤。

8.根据权利要求4至7中任一项所述的方法，其特征在于，所述模态估计步骤适合于提供构成一组冗余信息的每个模型的参数，该组冗余信息允许减小被估计模态的方差。

9.根据权利要求4至8中任一项所述的方法，其特征在于，所述适应性建模步骤进行参数类型的建模。

10.根据权利要求4至9中任一项所述的方法，其特征在于，所述适应性建模步骤进行ARMA也就是“自回归移动平均”类型的建模。

11.根据权利要求10的方法，其特征在于，所述ARMA类型的建模在每个时刻对于每个传感器和对于所有所考虑的阶数来进行。

12.根据权利要求4至11中任一项所述的方法，其特征在于，所述方法包括多项式矩阵的求逆的步骤，所述多项式矩阵是对称互谱矩阵，该对称互谱矩阵在其主对角线上表示每个传感器的功率谱密度，换言之，互谱。

13.根据权利要求4至12中任一项所述的方法，其特征在于，所述适应性建模步骤包括下述步骤的以在时刻n对于阶数N＝0为起始条件的关于时间、包括关于模型阶数N＝[1，2，...，Nmax]递推的递推：

-前和后线性预测误差向量的计算，

-后和前部分相关矩阵的计算，

-前和后线性预测误差的协方差矩阵的计算，

-前线性预测误差的幂矩阵的计算，

-增益向量的直接计算，

${\underset{\&OverBar;}{\mathrm{\&theta;}}}_{a,n}^N=\left\{{\left(\underset{k=0}{\overset{n}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}^{n-k}{\underset{\&OverBar;}{\mathrm{\&xi;}}}_k^N\underset{\&OverBar;}{\mathrm{\&xi;}}_k^N{}_t\right)}^{-1}{\underset{\&OverBar;}{\mathrm{\&xi;}}}_n^N\right\}{\mathrm{\&alpha;}}_n^N$

其中α为标量，λ为遗忘因子，并且

-通过知道

递推计算向量

和

-矩阵A_k的计算，A_k表示模型，k＝1至N。

14.根据权利要求1至13中任一项所述的方法，其特征在于，所述方法包括通过应用以下两个约束之一对来自所述模态估计步骤的模态进行分类的步骤：

-每类仅一个来自相同模型的模态，和

-所有估计不受估计来源影响地具有相同权重。

15.可加载于信息***中的计算机程序，所述程序包含允许实施根据权利要求1至14中任一项所述的方法的指令。

Images