trackerPHD
Multi-sensor, multi-object PHD tracker
Description
The trackerPHD
System object™ is a tracker capable of processing detections of multiple targets from multiple
sensors by using a multi-target probability hypothesis density (PHD) filter to estimate the
states of point targets and extended objects. PHD is a function defined over the state-space
of the tracking system, and its value at a state is defined as the expected number of targets
per unit state-space volume. The PHD is represented by a weighted summation (mixture) of
probability density functions, and peaks in the PHD correspond to possible targets. For an
overview of how the tracker functions, see Algorithms.
By default, the trackerPHD can track extended objects using the ggiwphd filter, which models
detections from an extended object as a parse points cloud. You can also use
trackerPHD with the gmphd filters, which tracks point targets and extended objects with designated
shapes. Inputs to the tracker are detection reports generated by objectDetection, fusionRadarSensor,
irSensor, or
sonarSensor
objects. The tracker outputs all maintained tracks and their analysis information.
To track targets using this object:
Create the
trackerPHDobject and set its properties.Call the object with arguments, as if it were a function.
To learn more about how System objects work, see What Are System Objects?
Creation
Description
tracker = trackerPHDtrackerPHD
System object with default property values.
tracker = trackerPHD(Name,Value)trackerPHD('MaxNumTracks',100) creates a PHD tracker that allows a
maximum of 100 tracks. Enclose each property name in quotes.
Properties
Usage
To process detections and update tracks, call the tracker with arguments, as if it were a function (described here).
Note
You must specify the SensorConfigurations property before using the
tracker.
Syntax
Description
confirmedTracks = tracker(detections,time)detections, at the update time, time.
Confirmed tracks are corrected and predicted to the update time.
confirmedTracks = tracker(detections,config,time)config. Use this syntax
when the configurations of sensors are changing with time. To enable this syntax, set the
HasSensorConfigurationsInput property to
true.
[
also returns a list of tentative tracks, confirmedTracks,tentativeTracks,allTracks] = tracker(___)tentativeTracks, and a list
of all tracks, allTracks. You can use this output syntax with any of
the previous input syntaxes.
[
also returns the analysis information, confirmedTracks,tentativeTracks,allTracks,analysisInformation] = tracker(___)analysisInformation, which can
be used for track analysis. You can use this output syntax with any of the previous input
syntaxes.
Input Arguments
Output Arguments
Object Functions
To use an object function, specify the
System object as the first input argument. For
example, to release system resources of a System object named obj, use
this syntax:
release(obj)
Examples
Set up the sensor configuration, create a PHD tracker, and feed the tracker with detections.
% Create sensor configuration. Specify clutter density of the sensor and % set the IsValidTime property to true. configuration = trackingSensorConfiguration(1); configuration.ClutterDensity = 1e-7; configuration.IsValidTime = true; % Create a PHD tracker. tracker = trackerPHD('SensorConfigurations',configuration); % Create detections near points [5;-5;0] and [-5;5;0] at t=0, and % update the tracker with these detections. detections = cell(20,1); for i = 1:10 detections{i} = objectDetection(0,[5;-5;0] + 0.2*randn(3,1)); end for j = 11:20 detections{j} = objectDetection(0,[-5;5;0] + 0.2*randn(3,1)); end tracker(detections,0);
Update the tracker again after 0.1 seconds by assuming that targets move at a constant velocity of [1;2;0] unit per second.
dT = 0.1;
for i = 1:20
detections{i}.Time = detections{i}.Time + dT;
detections{i}.Measurement = detections{i}.Measurement + [1;2;0]*dT;
end
[confTracks,tentTracks,allTracks] = tracker(detections,dT);Visualize detections and confirmed tracks.
% Obtain measurements from detections. d = [detections{:}]; measurements = [d.Measurement]; % Extract positions of confirmed tracking using getTrackPositions function. % Note that we used the default sensor configuration % FilterInitializationFcn, initcvggiwphd, which uses a constant velocity % model and defines the states as [x;vx;y;vy;z;vy]. positionSelector = [1 0 0 0 0 0;0 0 1 0 0 0;0 0 0 0 1 0]; positions = getTrackPositions(confTracks,positionSelector); figure() plot(measurements(1,:),measurements(2,:),'x','MarkerSize',5,'MarkerEdgeColor','b'); hold on; plot(positions(1,1),positions(1,2),'v','MarkerSize',5,'MarkerEdgeColor','r' ); hold on; plot(positions(2,1),positions(2,2),'^','MarkerSize',5,'MarkerEdgeColor','r' ); legend('Detections','Track 1','Track 2') xlabel('x') ylabel('y')
Create a tracking scenario and specify its StopTime and UpdateRate properties.
scenario = trackingScenario; scenario.StopTime = Inf; scenario.UpdateRate = 0;
Add a tower platform in the scenario and specify its dimensions.
Tower = platform(scenario,'ClassID',3); Tower.Dimensions = struct( ... 'Length',10, ... 'Width',10, ... 'Height',60, ... 'OriginOffset',[0 0 30]);
Add a car platform in the scenario. Specify its dimensions and trajectory.
Car = platform(scenario,'ClassID',2); Car.Dimensions = struct( ... 'Length',4.7, ... 'Width',1.8, ... 'Height',1.4, ... 'OriginOffset',[-0.6 0 0.7]); Car.Trajectory = waypointTrajectory( ... [0 -15 -0.23;0.3 -29.5 -0.23; 0.3 -42 -0.39;0.3 -56.5 -0.23; ... -0.3 -78.2 -0.23;4.4 -96.4 -0.23], ... [0; 1.4 ; 2.7; 4.1; 6.3; 8.2], ... 'Course',[-88; -89; -89; -92; -84; -71], ... 'GroundSpeed',[10; 10; 10; 10; 10; 10], ... 'ClimbRate',[0; 0; 0; 0; 0; 0], ... 'AutoPitch',true, ... 'AutoBank',true);
Create a non-scanning radar, specify its properties, and mount the sensor on the tower.
NoScanning = fusionRadarSensor('SensorIndex',1, ... 'UpdateRate',10, ... 'MountingAngles',[-90 0 0], ... 'FieldOfView',[20 10], ... 'ScanMode','No scanning', ... 'HasINS',true, ... 'DetectionCoordinates','Scenario', ... 'TargetReportFormat','Detections', ... 'HasElevation',true); Tower.Sensors = NoScanning;
Create a theater plot to visualize sensor coverage, tracks, and detections.
tp = theaterPlot('XLim',[-58 58],'YLim',[-104 12],'ZLim',[-109 8]); set(tp.Parent,'YDir','reverse','ZDir','reverse'); view(tp.Parent,-37.5,30); % Platform plotter for the car. platp = platformPlotter(tp,'DisplayName','Targets','MarkerFaceColor','k'); % Detection plotter for sensor detections. detp = detectionPlotter(tp,'DisplayName','Detections','MarkerSize',6, ... 'MarkerFaceColor',[0.85 0.325 0.098],'MarkerEdgeColor','k','History',10000); % Coverage plotter for sensor. covp = coveragePlotter(tp,'DisplayName','Sensor Coverage'); % Track plotter for tracks. tPlotter = trackPlotter(tp,'DisplayName','Tracks');
Extract the sensor configuration of the sensor and use it to specify a PHD tracker.
sensorConfig = trackingSensorConfiguration(scenario.Platforms{1}.Sensors{1}, ...
'SensorTransformFcn',@cvmeas,'FilterInitializationFcn',@initcvggiwphd);
tracker = trackerPHD('SensorConfigurations',sensorConfig, ...
'PartitioningFcn',@(x)partitionDetections(x,1,4.7),...
'AssignmentThreshold',20,'ExtractionThreshold',0.8,...
'ConfirmationThreshold',1.5,'MergingThreshold',20,...
'DeletionThreshold',2e-1,'BirthRate',1e-3,...
'HasSensorConfigurationsInput',true);Simulate the scenario, generate detections, and use the detections to track the car. Update the theater plot during simulation.
while advance(scenario) && ishghandle(tp.Parent) % Generate sensor data. [dets,configs,sensorConfigPIDs] = detect(scenario); % Read sensor data. allDets = [dets{:}]; if ~isempty(allDets) % Extract measurement positions. meas = cat(2,allDets.Measurement)'; % Extract measurement noise. measCov = cat(3,allDets.MeasurementNoise); else meas = zeros(0,3); measCov = zeros(3,3,0); end % Obtain true positions. truePoses = platformPoses(scenario); truePosition = vertcat(truePoses(:).Position); % Update tracker with the detections and sensor configuration. [cTracks,tTracks,allTracks] = tracker(dets,configs,scenario.SimulationTime); % Update the theater plot. plotPlatform(platp,truePosition); plotDetection(detp,meas,measCov); plotCoverage(covp,coverageConfig(scenario)); % Update the track plotter. Extract track positions. positionSelector = [1 0 0 0 0 0; 0 0 1 0 0 0; 0 0 0 0 1 0]; positions = getTrackPositions(cTracks,positionSelector); % Label and plot the tracks. if ~isempty(cTracks) labels = cell(numel(cTracks),1); for i =1:numel(cTracks) labels{i} = {['T',num2str(cTracks(i).TrackID)]}; end plotTrack(tPlotter,positions,labels); end drawnow end
More About
Algorithms
Probability hypothesis density (PHD) is a function defined over the state-space of the tracking system, and its value at a state is defined as the expected number of targets per unit state-space volume. The PHD is usually approximated by a mixture of components, and each component corresponds to an estimate of the state. The commonly used approximations of PHD are Gaussian mixture, SMC mixture, GGIW mixture, and GIW mixture.
To understand PHD, take the Gaussian mixture as an example. The Gaussian mixture can be represented by
where M is the total number of components, N(x|mi,Pi) is a normal distribution with mean mi and covariance Pi, and wi is the weight of the ith component. The weight wi denotes the number, which can be fractional, of targets represented by the ith component. Integration of D(x) over a state-space region results in the expected number of targets in that region. Integrating D(x) over the whole state space results in the total expected number of targets (∑ wi), since the integration of a normal distribution over the whole state space is 1. The x-coordinates of the peaks (local maximums) of D(x) represent the most likely states of targets.
For example, the following figure illustrates a PHD function given by D(x) = N(x|−4,2) + 0.5N(x|3,0.4) + 0.5N(x|4,0.4). The weight summation of these components is 2, which means that two targets probably exist. From the peaks of D(x), the possible positions of these targets are at x = −4, x = 3, and x = 4. Notice that the last two components are very close to each other, which means that these two components can possibly be attributed to one object.
References
[1] Granstorm, K., C. Lundquiest, and O. Orguner. " Extended target tracking using a Gaussian-mixture PHD filter." IEEE Transactions on Aerospace and Electronic Systems. Vol. 48, Number 4, 2012, pp. 3268-3286.
[2] Granstorm, K., and O. Orguner." A PHD filter for tracking multiple extended targets using random matrices." IEEE Transactions on Signal Processing. Vol. 60, Number 11, 2012, pp. 5657-5671.
[3] Granstorm, K., and A. Natale, P. Braca, G. Ludeno, and F. Serafino."Gamma Gaussian inverse Wishart probability hypothesis density for extended target tracking using X-band marine radar data." IEEE Transactions on Geoscience and Remote Sensing. Vol. 53, Number 12, 2015, pp. 6617-6631.
[4] Panta, Kusha, et al. “Data Association and Track Management for the Gaussian Mixture Probability Hypothesis Density Filter.” IEEE Transactions on Aerospace and Electronic Systems, vol. 45, no. 3, July 2009, pp. 1003–16.
[5] Ristic, B., et al. “Adaptive Target Birth Intensity for PHD and CPHD Filters.” IEEE Transactions on Aerospace and Electronic Systems, vol. 48, no. 2, 2012, pp. 1656–68.
