主要内容

基于粒子滤波块的Simulink参数和状态估计金宝app

这个示例使用:

打开脚本

此示例演示了系统识别工具箱™中的粒子过滤器块的使用。将离散时间传递函数参数估计问题进行重新格式化和递归地解决作为状态估计问题。

介绍

系统辨识工具箱有三个Simulink模块用于非线性状态估计:金宝app

  • 粒子滤波:实现离散时间粒子滤波算法。

  • 扩展卡尔曼滤波:实现一阶离散时间扩展卡尔曼滤波算法。

  • 无迹卡尔曼滤波:实现离散时间无迹卡尔曼滤波算法。

这些阻止使用以不同采样率运金宝app行的多个传感器支持状态估计。使用这些块的典型工作流程如下：

  1. 使用MATLAB或Simulink函数建模您的植物和传感器行为。金宝app

  2. 配置块的参数。

  3. 模拟过滤器并分析结果以获得对滤波器性能的信心。

  4. 在硬件上部署过滤器。您可以使用Simulink Coder™软件为这些过滤器生成代码。金宝app

此示例使用粒子过滤器块来演示此工作流程的前两个步骤。最后讨论了最后两个步骤下一个步骤部分。本例中的目标是递归地估计离散时间传递函数(输出误差模型)的参数，其中模型参数在每个时间步骤中随着新信息的到来而更新。

如果您对扩展卡尔曼滤波器感兴趣，请参阅示例“用多个、多速率传感器估计非线性系统的状态”。Unscented卡尔曼滤波器的使用遵循与扩展卡尔曼滤波器相似的步骤。

将示例文件夹添加到MATLAB路径中。

目录(fullfile (matlabroot,“例子”那“识别”那“主要”))；

植物建模

大多数状态估计算法依赖于植物模型（状态转换函数），该植物模型（状态转换函数）从一个时间步长描述植物状态的演变。此功能通常表示为$ x [k + 1] = f（x [k]，w [k]，u [k]）$其中x是状态，w是过程噪声，u是可选的附加输入，例如系统输入或参数。粒子过滤块需要你用稍微不同的语法提供这个函数，$ X [k + 1] = f {pf} u (X [k], [k])美元．的差异是:

  • 粒子过滤器的工作原理是跟踪许多状态假设(粒子)的轨迹，块将所有状态假设一次传递到你的函数。具体来说，如果你的状态向量x美元有N_s美元元素，你选择N_p美元使用粒子,X美元有尺寸美元[N_s \;N_p]美元每一列都是状态假设。

  • 计算过程噪声的影响w美元关于状态假设$ x [k + 1] $在你的函数$ f_ {pf}（...）$．该块不会对过程噪声的概率分布作出任何假设w美元，也不需要w美元作为输入。

这个函数$ f_ {pf}（...）$可以是MATLAB函数，符合MATLAB编码器™或SIMULINK功能块的限制。金宝app创建后$ f_ {pf}（...）$，则在Particle Filter块中指定函数名。

在这个例子中，您将离散时间传递函数参数估计问题重新表述为状态估计问题。这个传递函数可以表示离散时间过程的动力学，也可以表示与信号重构器(如零阶保持器)耦合的连续时间动力学。假设你有兴趣估计一阶离散时间传递函数的参数:

$ $ y [k] = \压裂{20问^ {1}}{1 - 0.7 - q ^ {1}} u e [k] + [k] $ $

这里美元$ y [k]为植物产量，u [k]美元为植物的输入，$ e [k] $是测量噪音，美元问^ {1}$时延算子是这样的吗美元问^ {1}u [k] = [t - 1]美元．参数化传递函数为$$ \ frac {n \; q ^ { - 1}} {1 + d \; q ^ { - 1}} $$,在那里n美元和$ d $是待估计的参数。通过选择状态向量，传递函数和参数可以多种方式表示为必要的状态空间形式。一个选择是$ x [k] = [\;Y [K];d [k];n [k] \;] $其中，第二和第三状态表示参数估计。然后可以等同地写入传输函数

$$ x[k+1] = \left[ c \开始{数组}{}& # xA;-x_2[k] [k] + x_3[k] [k] [k] [k];x_2 [k] \ \ & # xA;x_3 [k] \ \ & # xA;结束\{数组}& # xA;\] $ $

测量噪声项$ e [k] $是在传感器建模中处理的。在这个例子中，你在一个MATLAB函数中实现了上面的表达式，以矢量形式计算效率:

类型pfblockstateTransitionFcNexample.
函数xNext = pfBlockStateTransitionFcnExample (x, u) % pfBlockStateTransitionFcnExample粒子状态转换函数%过滤器,估算参数%的输出,一阶,离散时间传递% % %函数模型输入:% x -粒子,NumberOfStates-by-NumberOfParticles矩阵% u -系统输入,一个标量% %输出:实现状态转移函数% xNext = [x(1)*x(2) + x(3)*u;% x (2);% x (3)];%在矢量形式(为所有粒子)。xNext = x;xNext (1) = bsxfun (@times x (1:) - x (2:)) + x (3:) * u;添加一个小的过程噪声(相对于每个状态的预期大小)，以增加粒子多样性1飞行;1 e 1], randn(大小(xNext))); end

传感器建模

Particle Filter块要求您提供一个度量可能性函数，用于计算每个状态假设的可能性(概率)。这个函数有它的形式$ l [k] = h_ {pf}（x [k]，y [k]，u [k]）$．L [k]美元是一个N_p美元元素的向量,N_p美元是你选择的粒子数。$ m ^ {th} $元素L [k]美元可能性是$ m ^ {th} $粒子(列)美元$ X [k]．美元$ y [k]为传感器测量。u [k]美元是一个可选的输入参数，它可以与状态转换功能的输入不同。

传感器测量此示例中的第一个状态。该示例假设根据高斯分布分发实际和预测测量之间的错误，但可以使用任何任意概率分布或其他方法来计算似然性。你创造了美元$ h_ {pf}(…)，并在Particle Filter块中指定函数名。

类型pfBlockMeasurementLikelihoodFcnExample;
函数= pfBlockMeasurementLikelihoodFcnExample可能性(粒子,测量)% pfBlockMeasurementLikelihoodFcnExample测量粒子滤波% %的似然函数是第一个国家输入% %:% - % NumberOfStates-by-NumberOfParticles矩阵测量粒子系统输出,一个标量% %输出:% likelihood -一个具有NumberOfParticles元素的向量，其第n %元素是第n个粒子的可能性%#codegen % Predicted measurement yHat = particles(1，:);假设误差分布在每个多元正态分布上，%均值为0，方差为1。计算相应的概率%密度函数e = bsxfun(@minus, yHat, measurement(:)');% Error numberOfMeasurements = 1;μ= 0;%均值=眼(测量数);% Variance measurementErrorProd = dot((e-mu)， Sigma \ (e-mu)， 1);c = 1/√((2*pi)^numberOfMeasurements * det(Sigma));可能性= c * exp(-0.5 * measureenterrorprod); end

滤波器结构

配置粒子过滤器块以进行估计。您可以指定状态转换和测量似然函数名称，粒子数和这些颗粒的初始分布。

在系统模型选项卡，指定以下参数:

状态转换

  1. 指定状态转移函数，pfblockstateTransitionFcNexample.,在函数．当您输入函数名并单击时申请，块检测到您的功能具有额外的输入，你美元，并创建输入端口StateTransitionFcnInputs．您将系统输入连接到这个端口。

初始化

  1. 指定10000在粒子数量．粒子数越高，估计越准确，计算成本也就越高。

  2. 指定高斯在分布从多变量高斯分布获取初始粒子。然后指定[0;0;0]在的意思是因为你有三个状态，这是你的最佳猜测。指定诊断（[10 5 100]）在协方差为第三种状态指定较大的方差(不确定性)，为前两种状态指定较小的方差。至关重要的是，这一初始粒子集的分布范围足够广(大的方差)，以覆盖潜在的真实状态。

测量1

  1. 指定测量可能性函数的名称，pfBlockMeasurementLikelihoodFcnExample,在函数．

采样时间

  1. 在块对话框的底部，输入1采样时间．如果您在状态转换和测量似然函数之间具有不同的采样时间，或者如果您有多个具有不同样本时间的传感器，则可以在其中配置这些传感器块输出,多重速率的选项卡。

颗粒过滤器涉及使用具有更高可能性的诸如较高的似然和播种新颗粒的颗粒。这由选项控制重采样组。本示例使用默认设置。

默认情况下，区块只输出状态假设的平均值，并根据它们的可能性加权。要查看所有的粒子、权重，或选择提取状态估计的不同方法，请查看块输出,多重速率的选项卡。

仿真和结果

对于简单的测试，真正的工厂模型是用白噪声输入模拟的。从工厂的输入和噪声测量被送入粒子过滤器块。以下Simulink模型表示金宝app此设置。

open_system (“pfBlockExample”）

模拟系统，比较估计的和真实的参数:

SIM（“pfBlockExample”) open_system ('pfblockexample /参数估计'）

图中显示了真实的分子和分母参数，以及它们的粒子滤波估计。经过10个时间步长，估计值近似收敛于真值。即使初始状态猜想与真实值相差甚远，也能获得收敛性。

故障排除

这里列出了一些潜在的实现问题和故障排除思路，以防粒子过滤器没有按照您的应用程序的预期执行。

粒子滤波器的故障排除通常是通过查看粒子集和它们的权重来执行的，这些权重可以通过选择得到输出所有粒子和输出所有重量在块输出,多重速率的选项卡的块对话框。

第一个完整性检查是确保状态转换和度量可能性函数能够很好地捕获系统的行为。如果您的系统有一个仿真模型(因此可以访问仿真中的真实状态)，那么可以尝试使用真实状态初始化过滤器。然后，您可以验证状态转移函数是否准确地计算了真实状态的时间传播，以及测量似然函数是否为这些粒子计算了很高的可能性。

初始粒子集很重要。在你的模拟开始时，确保至少有一些粒子有很高的可能性。如果真实状态在状态假设的初始传播范围之外，估计可能是不准确的，甚至会出现分歧。

如果状态估计精度最初很好，但随着时间的推移而恶化，问题可能是粒子退化或粒子枯竭[2]。当粒子分布太广时，粒子会发生简并;当粒子重新采样后，粒子会聚集在一起，粒子会发生贫化。由于直接重采样，粒子退化导致粒子贫化。在本例中使用的状态转移函数中加入人工过程噪声是一种实用的方法。有大量关于解决这些问题的文献，并且基于您的应用程序可以使用更系统的方法。[1]和[2]是两个有用的参考文献。

下一个步骤

  1. 验证状态估计:一旦滤波器在模拟中按预期执行，通常使用广泛的蒙特卡罗模拟进一步验证性能。有关更多信息，请参见在Simulink中验证在线状态估计金宝app．你可以使用下面的选项随机性组的粒子过滤器对话框，以促进这些模拟。

  2. 生成代码：粒子过滤器块使用Simulink Coder™软件支持C和C ++代码。金宝app金宝app您提供给此块的功能必须符合MATLAB Coder™软件的限制（如果您使用的是使用MATLAB函数来建模您的系统）和Simulink编码器软件（如果您使用Simulink功能块来建模您的系统）。金宝app

总结

这个例子展示了如何在系统识别工具箱中使用粒子过滤器块。你递归地估计一个离散时间传递函数的参数，其中参数在每个时间步随着新信息的到来而更新。

close_system (“pfBlockExample”, 0)

从MATLAB路径中移除示例文件文件夹。

rmpath (fullfile (matlabroot,“例子”那“识别”那“主要”))；

参考文献

[1]西蒙，丹。最优状态估计：卡尔曼，H无限远和非线性方法。John Wiley＆Sons，2006。

[2] Doucet，Arnaud和Adam M. Johansen。“粒子过滤和平滑的教程：十五年后。”非线性滤波手册12.656-704（2009）：3。

另请参阅

||

相关的话题

系统识别工具箱文档

金宝app

尝试matlab，sim金宝appulink等产品下载188bet金宝搏

得到审判现在