利用GM-PHD跟踪器跟踪密集杂波中的点目标

打开实时脚本

这个例子向您展示了如何使用高斯混合概率假设密度(GM-PHD)跟踪器跟踪密集杂波中的点目标，该跟踪器使用恒定速度模型。

设置场景

本例中使用的场景是使用trackingScenario．该场景由五个以恒定速度移动的点目标组成。目标在静态2-D雷达传感器的视野内移动。你可以使用monostaticRadarSensor来模拟二维传感器，并将其安装在静态平台上。你可以使用FalseAlarmRate属性，以控制杂波的密度。的值FalseAlarmRate属性表示传感器的一个分辨率单元产生假报警的概率。基于的虚警率 $10^{- 3.}$ 根据本例中定义的传感器的分辨率，每一步大约会产生48个假警报。

%可重复的目标位置rng (2022);创建一个trackingScenario对象scene = trackingScenario(“UpdateRate”10“StopTime”10);numTgts = 5;初始化每个目标的位置和速度x = 100*(2*rand(numTgts,1) - 1);y = 100*(2*rand(numTgts,1) - 1);z = 0 (numTgts,1);vx = 5*randn(numTgts,1);vy = 5*randn(numTgts,1);vz = 0 (numTgts,1);在给定位置和速度的场景中添加平台。。为i = 1:numTgts thisTgt =平台(场景);thisTgt.Trajectory.Position = [x(i) y(i) z(i)];thisTgt.Trajectory.Velocity = [vx(i) vy(i) vz(i)];结束在场景中添加一个检测平台。detectingPlatform =平台(场景);detectingPlatform.Trajectory.Position = [-200 00];使用单站雷达传感器模拟二维雷达。雷达= monostaticRadarSensor(1，.．.“UpdateRate”场景。UpdateRate,.．.“DetectionProbability”, 0.9,.．.% Pd“FalseAlarmRate”1 e - 3,.．.% Pfa“FieldOfView”, 120年[1],.．.“ScanMode”，“没有扫描”，.．.“DetectionCoordinates”，“场景”，.．.%在场景框架中报告“之内”,真的,.．.使能INS在场景帧中报告“ReferenceRange”, 300,.．.%短程“ReferenceRCS”10.．.平台的默认RCS“MaxUnambiguousRange”, 400,.．.%短程“RangeResolution”, 1.．.“AzimuthResolution”, 1.．.“HasFalseAlarms”,真正的);%报告误报将雷达安装在探测平台上。。detectingPlatform。传感器=雷达;

设置跟踪器和指标

跟踪器

您使用GM-PHD点对象跟踪器来跟踪目标。配置PHD跟踪器的第一步是设置传感器的配置。定义配置trackingSensorConfiguration对象。

的SensorIndex的值设置为1以匹配模拟传感器的值。由于传感器是点对象传感器，每次扫描最多对每个对象输出一次检测，因此可以设置MaxNumDetsPerObject属性设置为1。

的SensorTransformFcn，SensorTransformParameters,SensorLimits一起允许您定义传感器可以检测到轨道的区域。的SensorTransformFcn定义履带状态的转换( $x_{跟踪}$ )进入传感器使用的中间空间( $x_{传感器}$ )定义航迹可探测性。计算检测概率的整体计算如下:

$x_{传感器} ＝ SensorTransformFcn （ x_{跟踪} ， SensorTransformParameters ）$

$P_{d} ＝$

$｛ \begin{array}{l} 配置． DetectionProbability & {SensorLimits （：， 1 ） \leq x}_{传感器} \leq SensorLimits （：， 2 ） \\ 配置． MinDetectionProbability & 否则 \end{array}$

为了计算具有给定状态协方差的不确定状态的检测概率，跟踪器使用类似于无味卡尔曼滤波器的sigma-point计算生成状态样本。

的签名SensorTransformFcn类似于典型的测量模型。因此，可以使用类似的函数cvmeas，cameas作为SensorTransformFcn．在本例中，您假设所有轨道都是可检测到的。因此,SensorTransformFcn定义为@ x (x, params)和SensorLimits定义为(负无穷到正无穷)对所有州都适用。

变换函数和极限sensorTransformFcn = @(x,params)x;sensorTransformParameters = struct;sensorLimits = [-inf inf].*ones(6,1);%传感器配置的检测概率Pd =雷达.探测概率;配置= trackingSensorConfiguration(“SensorIndex”, 1.．.“IsValidTime”,真的,.．.%更新每一步“MaxNumDetsPerObject”, 1.．.“SensorTransformFcn”sensorTransformFcn,.．.“SensorTransformParameters”sensorTransformParameters,.．.“SensorLimits”sensorLimits,.．.“DetectionProbability”, Pd);

的ClutterDensity配置属性是指测量空间单位体积的虚警率。在本例中，测量空间定义为场景框架中报告的检测的笛卡尔坐标。由于传感器在笛卡尔坐标下的分辨率体积随分辨率的方位角和范围的变化而变化，可以在传感器的平均范围内计算出近似的值。

%传感器参数，用于计算分辨率仓的体积Rm = radar.MaxUnambiguousRange/2;% mean -rangedTheta = radar.ElevationResolution;偏差分数降低了“有效”分辨率大小dPhi = radar.AzimuthBiasFraction*radar.AzimuthResolution;dR = radar.RangeBiasFraction*radar. rangerresolution;细胞体积%VCell = 2 * ((Rm +博士)^ 3 - Rm ^ 3) / 3 *(信德(dTheta / 2)) *函数(dPhi);误报率Pfa = radar.FalseAlarmRate;定义杂波密度配置。ClutterDensity = Pfa/VCell;

您还定义了FilterInitializationFcn指定由此传感器初始化的过滤器类型和过滤器中组件的分布。在本例中，设置FilterInitializationFcn来initcvgmphd，它创建了一个匀速GM-PHD滤波器，并为跟踪器的每个低似然检测添加一个组件。的initcvgmphd在没有检测的情况下调用时不添加任何组件。这意味着在此配置下，仅当检测落在AssignmentThreshold多目标文件。看到一个ssignmentThreshold的属性trackerPHD欲知详情。

配置。FilterInitializationFcn = @initcvgmphd;

对象使用此配置创建跟踪器trackerPHD系统对象™。配置跟踪器时，可以指定出生率属性定义单位时间内在视野中出现的目标数目。的FilterInitializationFcn与配置一起使用时，为每个未分配的检测添加一个组件。在每个时间步中，您可以期望组件的数量大约等于假警报和新目标的数量。跟踪器分发出生率对所有这些分量都一视同仁。

%用雷达单元数计算单位时间内虚警数。NCells = radar.MaxUnambiguousRange/radar. rangerresolution *radar.FieldOfView(1)/radar.AzimuthResolution;numFalse = Pfa*NCells;为每个新组件选择0.05的初始权重。作为新的数量%目标是不知道的，新数量的组件被简单地用作数字误报的百分比。0.1是更新速率。出生率= 0.05*(numFalse)/0.1;使用定义的传感器配置创建跟踪器。追踪器=追踪器“出生率”出生率,.．.“SensorConfigurations”、配置);

指标

要评估跟踪器的性能，还需要为性能评估设置一个度量。在本例中，您使用了广义最优子模式分配(GOSPA)度量。GOSPA指标旨在通过分配单个成本值来评估跟踪器的性能。跟踪性能越好，GOSPA成本越低。值为0表示完美跟踪。

创建gospa度量对象。gospa = trackGOSPAMetric;在每一步初始化存储度量的变量。= 0 (0,1);Loc = 0 (0,1);Mt = 0 (0,1);Ft = 0 (0,1);

运行仿真

接下来，您将推进场景，从场景中收集检测，并在模拟检测上运行PHD跟踪器。

%创建一个显示display = helperClutterTrackingDisplay(场景);Count = 1;用于存储度量数据的计数器rng (2018);%可重复运行而推进(场景)%当前时间time = scene.SimulationTime;电流检测百分比检测=检测(场景);%的踪迹Tracks =跟踪器(检测，时间);%更新显示显示(场景、探测、轨迹);计算GOSPA。getTruth函数定义如下。真相= getTruth(场景);[~，gospaMetric(count)，~，loc(count)，mt(count)，ft(count)] = gospa(跟踪，真相);Count = Count + 1;结束

在这个场景中，可以看到PHD跟踪器跟踪了所有目标。这也可以使用GOSPA度量和相关组件进行定量评估。在下面的图中，请注意GOSPA指标在几个步骤后下降。由于每条航迹的建立延迟，GOSPA度量的初始值较高。

图(“单位”，“归一化”，“位置”，[0.1 0.1 0.8 0.8]) subplot(2,2，[1 2]);情节(gospaMetric“线宽”2);标题(“总GOSPA”）;ylabel (“GOSPA指标”）;包含(“时间步”）;网格在；次要情节(2、2、3);情节(太“线宽”2);标题(“错失目标GOSPA”）;ylabel (“未达到的目标指标”）;包含(“时间步”）;网格在；次要情节(2、2、4);图(见英国《金融时报》“线宽”2);标题(“假轨GOSPA”）;ylabel (“虚假目标度量”）;包含(“时间步”）;网格在；

分析性能

确定跟踪器性能的典型方法是在不同的场景实现上运行多个模拟。蒙特卡洛模拟有助于消除随机事件的影响，如假警报的位置，目标偏离事件和测量中的噪声。

在本节中，您将运行该场景的不同实现和具有不同误报率的跟踪器，并为每个实现计算系统的平均GOSPA。运行场景和计算系统的平均GOSPA的过程封装在helper函数中helperRunMonteCarloAnalysis．为加速蒙特卡罗模拟，可以为monostaticRadarSensor模型以及trackerPHD使用MATLAB®Coder™工具箱。生成代码的过程被包装在辅助函数中helperGenerateCode．要为算法生成代码，需要将代码组装成一个独立的函数。这个函数被命名为clutterSimTracker_kernel在这个例子中。的clutterSimTracker_kernel函数被编写为支持四种误报率。金宝app要重新生成具有不同误报率的代码，可以使用以下命令。

Pfas = 10.^(-4，-3,4);选择您的假警报设置helperGenerateCode(场景,pfa);生成代码

要重新生成用于更多假警报设置的代码，可以修改函数以包含更多持久变量。有关如何为系统对象™生成代码的详细信息，请参阅如何为跟踪器生成C代码的例子。

的helperRunMonteCarloAnalysis函数使用parfor来执行每次蒙特卡洛运行。如果您拥有Parallel Computing Toolbox™许可证，那么蒙特卡洛运行可以分布在多个工作程序上，以进一步加速模拟。

打开标志以运行蒙特卡罗模拟Pfas = 10.^linspace(-4，-3,4);runMonteCarlo = false;如果runMonteCarlo numRunsPerPfa = 50;% #好< UNRCH >= 0 (4,numRunsPerPfa,4);为i = 1:4 gospaMCs(i，:，:) = helperRunMonteCarloAnalysis(scene, Pfas, Pfas(i)， numRunsPerPfa);结束保存clutterTrackingMCRuns.matgospaMCs；其他的%重载结果来自上次运行负载(“clutterTrackingMCRuns.mat”，“gospaMCs”）;结束%绘图结果图(“单位”，“归一化”，“位置”，[0.1 0.1 0.8 0.8]) subplot(2,2，[1 2]);情节(gospaMCs(:: 1)”,“线宽”2);标题(“GOSPA”）;传奇(strcat ('Pfa = '字符串(pfa)),“定位”，“水平”，“位置”，“NorthOutside”）;ylabel (“平均GOSPA指标”）;包含(“蒙特卡洛跑”）;网格在；次要情节(2、2、3);情节(gospaMCs(:: 3)”,“线宽”2);标题(“错失目标GOSPA”）;ylabel (“平均未达到目标指标”）;包含(“蒙特卡洛跑”）;网格在；次要情节(2、2、4);情节(gospaMCs(:: 4)”,“线宽”2);标题(“假轨GOSPA”）;ylabel (“平均虚假跟踪指标”）;包含(“蒙特卡洛跑”）;网格在；

上面的图显示了跟踪器在这种情况下的性能，通过对每个虚警率运行50个蒙特卡洛实现。随着虚警率的增加，产生假航迹的概率也会增加。在传感器附近，这种可能性甚至更高，因为那里分辨率单元的密度要高得多。由于假警报在这一区域间隔较近且出现频率较高，因此可能被错误地归类为低速假轨。跟踪器的这种行为可以在每个场景运行的GOSPA度量的平均“假跟踪组件”中观察到。注意，随着误报率的增加，图中的峰值数量也会增加。这也会导致总GOSPA指标的增加。“遗漏的目标组件”在除一次运行外的所有运行中都为零。这种类型的事件是由传感器多次未能击中目标引起的。

总结

在本例中，您学习了如何配置和初始化GM-PHD跟踪器，以跟踪给定虚警率的点目标。您还学习了如何使用GOSPA度量及其相关组件评估跟踪器的性能。此外，您还学习了如何在不同的假警报设置下运行该场景的多个实现，以限定跟踪器的性能特征。

效用函数

`helperRunMonteCarloAnalysis`

函数明确的gospaMC = helperRunMonteCarloAnalysis(场景，Pfas, Pfai, numRuns)clutterSimTracker_kernel；根据提供的Pfa计算要运行的跟踪器settingToRun = find(pas == Pfai,1，“第一”）;为每次Monte Carlo运行初始化ospa高斯amc = 0 (numRuns,4);使用parfor进行并行模拟parfori = 1:numRuns rng(i，“旋风”）;%在每次运行前重新启动场景重启(现场);在此场景中，%度量vs时间= 0 (0,4);%计数器计数= 0;%前进场景并运行而提前(场景)计数=计数+ 1;%使用雷达传感器使用targetpose功能own = scene.Platforms{end};tgtpositions = targetpositions(自己的，“rotmat”）;insPose =姿势(自己的);姿势=平台姿势(场景，“rotmat”）;真相=摆(1:end-1);time = scene.SimulationTime;在步骤结束时重置跟踪器，以便在下一次调用跟踪器时%重置到时间0。systemToReset = settingToRun*(abs(time - scene.StopTime) < 1e-3);存储ospa度量。跟踪和检测可以输出运行%场景的代码生成，而不需要蒙特卡罗运行。[~，~，gospaMetric(count，:)] = clutterSimTracker_kernel(Pfas, tgtpositions, insPose, truth, time, settingToRun, systemToReset);结束计算本次运行中GOSPA的平均值%从时间步骤20开始，允许轨迹建立。gospaMC(i，:) = mean(gospaMetric(20:end，:));结束结束

`helperGenerateCode`

函数helperGenerateCode(场景,pfa)% #好< DEFNU >使用platformpose生成样本姿势姿势=平台姿势(场景，“rotmat”）;为内核函数生成样例输入％% Pfas必须是编译时常数，不能随时间变化Pfas = code . constant (Pfas);%目标姿态为一个最大20个元素的可变大小数组tgtpositions = code .typeof(pose (1)，[20 1]，[10 0]);% Scalar ins poseinsPose =姿势(场景。平台{1},“真正的”）;%真值信息，最多20个目标真相= code .typeof(pose (1)，[20 1]，[10 0]);%的时间时间= 0;%要运行哪个假警报设置和要重置哪个。重置是必要的%以重新初始化跟踪器systemToRun = 1;systemToReset = 1;input = {Pfas, tgtpositions, insPose, truth, time, systemToRun, systemToReset};% #好< NASGU >%与MATLAB文件相同的名称，允许阴影MATLAB函数时% MEX文件可用，并自动在MEX中执行代码。codegenclutterSimTracker_kernelarg游戏输入- oclutterSimTracker_kernel；结束

`getTruth`

函数真相= getTruth(场景)platpose = platformpose(场景);%真实信息= platpose (1:end-1);最后一个对象是所有权结束

`clutterSimTracker_kernel`

此函数定义在一个名为clutterSimTracker_kernel，与此脚本在同一工作文件夹中可用。

函数[detects, tracks, ospaMetric] = clutterSimTracker_kernel(Pfas, tgtpositions, insPose, truth, time, systemToRun, systemToReset) assert(nummel (Pfas) == 4，只支持4个假警报设置。金宝app重写更多持久变量以添加更多设置）;持续的Tracker1 tracker2 tracker3 tracker4.．.雷达1雷达2雷达3雷达4.．.．Reset1 reset2 reset3 reset4.．.．gospa如果|| isempty(reset1) || isempty(reset2) || isempty(reset3) || isempty(reset4) gospa = trackGOSPAMetric(“CutoffDistance”, 50);结束如果isempty(tracker1) || isempty(radar1) || isempty(reset1) tracker1 = setupTracker(Pfas(1));radar1 = setupRadar(Pfas(1));Reset1 = 0 (1,1);结束如果isempty(tracker2) || isempty(radar2) || isempty(reset2) tracker2 = setupTracker(Pfas(2));radar2 = setupRadar(Pfas(2));Reset2 = 0 (1,1);结束如果isempty(tracker3) || isempty(radar3) || isempty(reset3) tracker3 = setupTracker(Pfas(3));radar3 = setupRadar(Pfas(3));Reset3 = 0 (1,1);结束如果isempty(tracker4) || isempty(radar4) || isempty(reset4) tracker4 = setupTracker(Pfas(4));radar4 = setupRadar(Pfas(4));Reset4 = 0 (1,1);结束开关systemToRun情况下1 detection = radar1(tgtpose, insPose, time);Tracks = tracker1(检测，时间);情况下2 detection = radar2(tgtpose, insPose, time);Tracks = tracker2(检测，时间);情况下3 detection = radar3(tgtpose, insPose, time);Tracks = tracker3(检测，时间);情况下4 detection = radar4(tgtpose, insPose, time);Tracks = tracker4(检测，时间);否则错误(“Idx出界”）;结束[~，gp，~，loc,mT,fT] = gospa(轨道，真理);ospaMetric = [gp loc mT fT];开关systemToReset情况下1 reset1 = 0 (0,1);情况下2 reset2 = 0 (0,1);情况下3 reset3 = 0 (0,1);情况下4 reset4 = 0 (0,1);结束结束函数跟踪器= setupTracker(Pfa) Pd = 0.9;VCell = 0.0609;单元格体积内联值%;Kc = Pfa/VCell;杂波密度%出生率出生率= 0.05*(5 + Pfa*48000)/0.1;变换函数和极限sensorTransformFcn = @(x,params)x;sensorTransformParameters = struct;sensorLimits = bsxfun(@times，[-inf inf]，ones(6,1));将信息组装到trackingSensorConfiguration配置= trackingSensorConfiguration(“SensorIndex”, 1.．.“IsValidTime”,真的,.．.%更新每一步“DetectionProbability”帕金森病,.．.% Detection可检测状态的概率“ClutterDensity”Kc,.．.杂波密度“SensorTransformFcn”sensorTransformFcn,.．.“SensorTransformParameters”sensorTransformParameters,.．.“SensorLimits”sensorLimits,.．.“MaxNumDetsPerObject”, 1.．.“FilterInitializationFcn”, @initcvgmphd);追踪器=追踪器“SensorConfigurations”配置,.．.“出生率”,出生率);结束函数雷达= setupRadar(Pfa) Pd = 0.9;雷达= monostaticRadarSensor(1，.．.“UpdateRate”10.．.“DetectionProbability”帕金森病,.．.“FalseAlarmRate”Pfa,.．.“FieldOfView”, 120年[1],.．.“ScanMode”，“没有扫描”，.．.“DetectionCoordinates”，“场景”，.．.%在场景框架中报告“之内”,真的,.．.使能INS在场景帧中报告“ReferenceRange”, 300,.．.%短程“ReferenceRCS”10.．.“MaxUnambiguousRange”, 400,.．.%短程“RangeResolution”, 1.．.“AzimuthResolution”, 1.．.“HasFalseAlarms”,真正的);结束