MATLAB静止背景下的多目标追踪

随着计算机技术以及智能汽车行业的发展，多目标的检测与追踪的实用性与研究价值逐渐提高。在计算机视觉的三层结构中，目标跟踪属于中间层，是其他高层任务，例如动作识别以及行为分析等的基础。其主要应用可包括视频监控，检测异常行为人机交互，对复杂场景中目标交互的识别与处理，以及虚拟现实及医学图像。

目标跟踪又包括单目标跟踪和多目标跟踪。单目标跟踪可以通过目标的表观建模或者运动建模，以处理光照、形变、遮挡等问题，而多目标跟踪问题则更加复杂，除了单目标跟踪回遇到的问题外，还需要目标间的关联匹配。另外在多目标跟踪任务中经常会碰到 目标的频繁遮挡、轨迹开始终止时刻未知、目标太小、表观相似、目标间交互、低帧率等等问题。

静止背景下的多目标追踪可分为两步来实现，第一步是在视频文件的每帧中检测出移动的目标，第二步是将检测到的目标与跟踪轨迹实时匹配。在本次实验中，利用混合高斯模型进行背景减除，使用形态学操作消除噪声，通过卡尔曼滤波预测目标位置，最后利用匈牙利算法进行匹配，实现静止背景下的多目标追踪。

1 实验原理

1.1 混合高斯模型

单高斯模型是利用高维高斯分布概率来进行模式分类：

1.2 卡尔曼滤波

卡尔曼滤波是一种利用线性系统状态方程，通过系统输入输出观测数据，对系统状态进行最优估计的算法。其核心思想是，根据当前的测量值、上一时间的预测值以及误差，计算得到当前值，并可以持续预测下一时间的值。

试验中利用卡尔曼滤波计算并预测每个轨迹在下一帧中的位置，建立帧间轨迹的关系。卡尔曼滤波将跟踪分为5种状态：新目标出现、目标匹配、目标遮挡、目标分离以及目标消失。其状态方程如下所示:

定义好了观测方程与状态方程之后就可以用卡尔曼滤波器实现运动目标的跟踪，步骤如下：

1）计算运动目标的特征信息（运动质心，以及外接矩形）。

2）用得到的特征信息初始化卡尔曼滤波器。

3）用卡尔曼滤波器对下一帧中对应的目标区域进行预测，当下一帧到来时，在预测区域内进行目标匹配。

4）如果匹配成功，则更新卡尔曼滤波器。

1.3 匈牙利匹配算法

匈牙利匹配算法是一种利用增广路径求取二分图最大匹配的算法。匈牙利树一般由宽度优先搜索（BFS）构成。从一个未匹配点出发运行 BFS，且必须走交替路，直至不能再扩展为止。交替路指的是，从一个未匹配点出发，依次经过非匹配边、匹配边...循环往复，形成的路径称为交替路。

示意图如图1所示：

图1 匈牙利树匹配

匈牙利算法的要点如下：

1. 从左边第一个顶点开始，挑选未匹配点进行搜索，寻找增广路。
2. 如果经过一个未匹配点，则寻找成功。更新路径信息，匹配变数+1，停止搜索。
3. 如果一直没有找到增广路，则不再从这个点开始搜索。
4. 找到增广路之后需要沿着路径更新匹配，通过prev 数组来记录路径上的点。

在实验中，用匈牙利算法将新一帧图片中检测到的运动物体匹配到对应的轨迹。匹配的过程是通过最小化卡尔曼预测得到的质心与检测到的质心之间的欧氏距离之和实现的。

通过卡尔曼滤波计算并预测每个轨迹在下一帧中的位置，然后计算预测的轨迹位置和每个新检测到的目标之间的欧几里得距离，将度量结果作为损失函数矩阵。损失矩阵的大小为（M,N），其中M是轨迹数目，N是检测到的运动物体数目。

2 实验内容

2.1 目标检测

要实现目标检测，首先利用混合高斯模型区分前景背景。通过调用函数vision.ForegroundDetector设置检测子为混合高斯模型，其中参数分别为高斯核数目、训练背景帧数以及背景阈值。函数的返回值为一个二进制掩码，其中1的像素值对应前景，0的像素值对应背景。在实验中，选取前150帧图像作为背景帧，并设置阈值为0.6，高斯核数目为3。

完成背景减除后，通过设置blob分析子来寻找连通域，函数设置的参数为最小区域面积，返回值为目标面积、质心和边界框。实验中，设置最小区域面积为400，代码如下所示：

obj.detector = vision.ForegroundDetector('NumGaussians', 3, ...

'NumTrainingFrames', 150, 'MinimumBackgroundRatio', 0.6);

obj.blobAnalyser = vision.BlobAnalysis('BoundingBoxOutputPort', true,'AreaOutputPort', true, 'CentroidOutputPort', true, 'MinimumBlobArea', 400);

完成混合高斯混合模型以及blob分析子后，逐帧读取视频，为后续目标检测以及追踪的过程做准备：

frame = obj.reader.step();

在检测目标的过程中，利用形态学运算中的开运算以及闭运算可以消除噪声，使目标检测更为准确。开运算是通过先腐蚀再膨胀，去除孤立的像素点、总的位置和结构不变。而闭运算是先膨胀再腐蚀，弥合小裂缝，而总的位置和形状不变，通过填充图像的凹角来滤波图像。其效果如图2所展示：

图2.1 原始像素点 图2.2 开运算效果图

图2.3 闭运算效果图

图2 形态学运算效果图

实验中，首先使用检测子，即混合高斯模型得到前景图，对前景图使用8*8矩形进行开运算，切断临近物体间的联系。再使用15*15矩形进行闭运算，消除细小物体，最后填补物体中间的空洞。滤除噪声后，使用blob分析得到所有连通域的中心以及边界框的大小。代码如下所示：

function [centroids, bboxes, mask] = detectObjects(frame)

% Detect foreground.

mask = obj.detector.step(frame);

% Apply morphological operations to remove noise and fill in holes.

mask = imopen(mask, strel('rectangle', [8,8]));

mask = imclose(mask, strel('rectangle', [15,15]));

mask = imfill(mask, 'holes');

% Perform blob analysis to find connected components.

[~, centroids, bboxes] = obj.blobAnalyser.step(mask);

end

视频中对移动目标的检测结果如图3.1,图3.2所示：

图3.1 移动目标检测结果1

图3.2 移动目标检测结果2

2.2 目标跟踪

要进行目标跟踪，首先需要进行轨迹初始化，通过函数initializeTracks()来进行初始化，每一个轨迹代表视频中一个移动的目标。轨迹的结构包含如下信息：

1. ID，轨迹编号；
2. Bbox，目标的边界框；
3. kalmanFilter，用于预测目标位置的卡尔曼滤波器；
4. Age，目标被检测到的总帧数；
5. totalVisibleCount，目标可被检测到的全部帧数；
6. consecutiveInvisibleCount：连续未检测到目标的帧数。

代码如下所示：

function tracks = initializeTracks()

% create an empty array of tracks

tracks = struct(...

'id', {}, ...

'bbox', {}, ...

'kalmanFilter', {}, ...

'age', {}, ...

'totalVisibleCount', {}, ...

'consecutiveInvisibleCount', {});

end

为消除噪声对目标追踪的影响，仅在totalVisibleCount超过阈值时才显示目标的轨迹。当连续几帧没有检测到与跟踪相关的信息时，则假设该对象已经离开了可视图画面。通过参数consecutiveinvisiblecount可判断这种情况，当其超过阈值时，删除跟踪轨迹。如果跟踪时间较短，并且在大多数帧中标记为不可见，那么轨迹也可能作为噪声被删除。

轨迹初始化完成后，通过卡尔曼滤波计算并预测每个轨迹在下一帧的位置，函数输出为边框预测中心。然后调整目标边界框的位置，使其中心到达预测位置，并将结果作为轨迹的跟踪矩形框，代码如下所示：

function predictNewLocationsOfTracks()

for i = 1:length(tracks)

bbox = tracks(i).bbox;

% Predict the current location of the track.

predictedCentroid = predict(tracks(i).kalmanFilter);

%根据以前的轨迹，预测当前位置

% Shift the bounding box so that its center is at

% the predicted location.

predictedCentroid = int32(predictedCentroid) - bbox(3:4) / 2;

tracks(i).bbox = [predictedCentroid, bbox(3:4)];

end

end

完成位置预测后，创建损失函数矩阵，航代表轨迹，列代表检测到的目标。损失矩阵的大小为（M,N），其中M是轨迹数目，N是检测到的运动物体数目。对每个轨迹计算其卡尔曼滤波预测的轨迹位置和每个新检测到的目标之间的欧几里得距离，将度量结果作为损失函数矩阵。

function [assignments, unassignedTracks, unassignedDetections] = ...

detectionToTrackAssignment()

nTracks = length(tracks);

nDetections = size(centroids, 1);

% Compute the cost of assigning each detection to each track.

cost = zeros(nTracks, nDetections);

for i = 1:nTracks

cost(i, :) = distance(tracks(i).kalmanFilter, centroids);

end

通过函数assignDetectionsToTracks(cost, costOfNonAssignment)利用匈牙利匹配算法将新一帧图片中检测到的运动物体匹配到对应的轨迹。其中输入的参数为损失矩阵以及阈值，低于阈值时，取消匹配。返回值为匹配的结果以及未匹配成功的轨迹以及目标。

完成匹配后，对已分配的轨迹，将其更新至当前帧目标所在位置，对未分配的轨迹，增加其连续不可见帧数。设置两个阈值，invisibleForLong代表当连续不可见帧数大于它时，删除轨迹；ageThreshold代表当总出现帧数小于它时，当该参数与总可见帧数的比值小于0.6时，删除轨迹。

目标检测与轨迹跟踪的操作循环进行，直至视频结束，显示最终跟踪结果。

2.3 总结

实验程序可分解为11个部分，创建系统对象、初始化轨迹、读取视频帧、检测目标、预测已跟踪轨迹的新位置、分配新检测目标给轨迹、更新已分配的轨迹、更新未分配的轨迹、删除丢失的轨迹、创建新轨迹以及显示跟踪结果。

其中创建系统对象包含创建视频对象，设置检测子为高斯混合模型以及设置blob分析子。在预测已跟踪轨迹的新位置时，使用卡尔曼滤波器进行预测，并调整位置，显示预测结果矩形框。在更新已分配的轨迹部分中，根据轨迹对应的检测目标位置中心修正其卡尔曼滤波器，并修正轨迹存在帧数、目标检测到的总帧数以及连续未检测到目标的帧数。在删除丢失轨迹中，删除连续不可见帧数大于阈值或当轨迹存在帧数小于10时，根据总可见帧数与轨迹存在帧数的比值丢弃轨迹。

主函数代码如下：

obj = setupSystemObjects(); %创建系统对象

tracks = initializeTracks(); % 初始化轨迹

nextId = 1; % ID of the next track

% Detect moving objects, and track them across video frames.

while ~isDone(obj.reader)

frame = readFrame(); %读取一帧

[centroids, bboxes, mask] = detectObjects(frame);

predictNewLocationsOfTracks();

[assignments, unassignedTracks, unassignedDetections] = ...

detectionToTrackAssignment();

updateAssignedTracks();

updateUnassignedTracks();

deleteLostTracks();

createNewTracks();

displayTrackingResults();

end

可见完成创建系统对象以及初始化轨迹后，循环剩余9个步骤，直至处理完整个视频。

3 实验结果分析

在本次试验中，通过混合高斯模型背景减除，形态学操作消除噪声，利用卡尔曼滤波预测每个轨迹在下一帧中的位置，最后通过匈牙利匹配算法完成目标与轨迹之间的匹配。实现了在静止背景下的多目标检测与跟踪。

如图4所示可见，对于匀速移动的目标，检测效果较好，可以实现较好的目标检测与跟踪。

图4 匀速移动目标检测结果

但静止背景下，基于动态的多目标追踪，很容易受到环境的影响。当检测目标为行人时，风吹动树叶或是有车辆经过，甚至是光照导致的影子变化都会很大程度地影响跟踪效果。图5展示了明显的失败样例：

图5 目标追踪失败样例

在实验中使用的卡尔曼滤波预测目标下一帧所在位置的模型，只适用于匀速变化，而汽车存在加速运动。从右边的二值图像可以看出，汽车的车灯造成像素点大量变化，只使用帧间相减以及形态学操作来得到检测目标的本实验并不适用于此。造成了较大的误差，导致检测以及跟踪目标失败。要规避由前景中不同物体造成的实验误差，可使用其他方法在前景检测中对检测到的物体进行分类。

通过多次试验发现，程序参数的鲁棒性较差，在进行形态学计算时，不同的视频调整开运算以及闭运算的结构元素，目标检测的效果差异很大。

在后续的学习过程中可以继续对本次实验进行更新修正，得到更好的实验效果。

附录

function objectracking()

obj = setupSystemObjects(); %初始化

tracks = initializeTracks(); % 初始化轨迹

nextId = 1; % ID of the next track

% Detect moving objects, and track them across video frames.

while ~isDone(obj.reader)

frame = readFrame(); %读取一帧

[centroids, bboxes, mask] = detectObjects(frame);

predictNewLocationsOfTracks();

[assignments, unassignedTracks, unassignedDetections] = ...

detectionToTrackAssignment();

updateAssignedTracks();

updateUnassignedTracks();

deleteLostTracks();

createNewTracks();

displayTrackingResults();

end

function obj = setupSystemObjects()

obj.reader = vision.VideoFileReader('test.mp4');

obj.maskPlayer = vision.VideoPlayer('Position', [740, 400, 700, 400]);

obj.videoPlayer = vision.VideoPlayer('Position', [20, 400, 700, 400]);%创建视频播放对象

obj.detector = vision.ForegroundDetector('NumGaussians', 3, ...

'NumTrainingFrames', 150, 'MinimumBackgroundRatio', 0.6);

%GMM前景检测，高斯核数目3，前150帧为背景帧，阈值为0.6，返回值为背景

obj.blobAnalyser = vision.BlobAnalysis('BoundingBoxOutputPort', true, ...

'AreaOutputPort', true, 'CentroidOutputPort', true, ...

'MinimumBlobArea', 400);%返回值为目标面积、质心和边界框

end

function tracks = initializeTracks()

% create an empty array of tracks

tracks = struct(...

'id', {}, ...

'bbox', {}, ...

'kalmanFilter', {}, ...

'age', {}, ...

'totalVisibleCount', {}, ...

'consecutiveInvisibleCount', {});

end

function frame = readFrame()

frame = obj.reader.step();

end

function [centroids, bboxes, mask] = detectObjects(frame)

% Detect foreground.

mask = obj.detector.step(frame);

% Apply morphological operations to remove noise and fill in holes.

mask = imopen(mask, strel('rectangle', [8,8]));

mask = imclose(mask, strel('rectangle', [15,15]));

mask = imfill(mask, 'holes');

% Perform blob analysis to find connected components.

[~, centroids, bboxes] = obj.blobAnalyser.step(mask);

end

function predictNewLocationsOfTracks()

for i = 1:length(tracks)

bbox = tracks(i).bbox;

% Predict the current location of the track.

predictedCentroid = predict(tracks(i).kalmanFilter);

%根据以前的轨迹，预测当前位置

% Shift the bounding box so that its center is at

% the predicted location.

predictedCentroid = int32(predictedCentroid) - bbox(3:4) / 2;

tracks(i).bbox = [predictedCentroid, bbox(3:4)];

end

end

function [assignments, unassignedTracks, unassignedDetections] = ...

detectionToTrackAssignment()

nTracks = length(tracks);

nDetections = size(centroids, 1);

% Compute the cost of assigning each detection to each track.

cost = zeros(nTracks, nDetections);

for i = 1:nTracks

cost(i, :) = distance(tracks(i).kalmanFilter, centroids);

end

% Solve the assignment problem.

costOfNonAssignment = 20;

[assignments, unassignedTracks, unassignedDetections] = ...

assignDetectionsToTracks(cost, costOfNonAssignment);

end

function updateAssignedTracks()

numAssignedTracks = size(assignments, 1);

for i = 1:numAssignedTracks

trackIdx = assignments(i, 1);

detectionIdx = assignments(i, 2);

centroid = centroids(detectionIdx, :);

bbox = bboxes(detectionIdx, :);

correct(tracks(trackIdx).kalmanFilter, centroid);

tracks(trackIdx).bbox = bbox;

tracks(trackIdx).age = tracks(trackIdx).age + 1;

tracks(trackIdx).totalVisibleCount = ...

tracks(trackIdx).totalVisibleCount + 1;

tracks(trackIdx).consecutiveInvisibleCount = 0;

end

end

function updateUnassignedTracks()

for i = 1:length(unassignedTracks)

ind = unassignedTracks(i);

tracks(ind).age = tracks(ind).age + 1;

tracks(ind).consecutiveInvisibleCount = ...

tracks(ind).consecutiveInvisibleCount + 1;

end

end

%% Delete Lost Tracks

function deleteLostTracks()

if isempty(tracks)

return;

end

invisibleForTooLong = 20;

ageThreshold = 10;

ages = [tracks(:).age];

totalVisibleCounts = [tracks(:).totalVisibleCount];

visibility = totalVisibleCounts ./ ages;

lostInds = (ages < ageThreshold & visibility < 0.6) | ...

[tracks(:).consecutiveInvisibleCount] >= invisibleForTooLong;

tracks = tracks(~lostInds);

end

function createNewTracks()

centroids = centroids(unassignedDetections, :);

bboxes = bboxes(unassignedDetections, :);

for i = 1:size(centroids, 1)

centroid = centroids(i,:);

bbox = bboxes(i, :);

kalmanFilter = configureKalmanFilter('ConstantVelocity', ...

centroid, [200, 50], [100, 25], 100);

end

end

function displayTrackingResults()

frame = im2uint8(frame);

mask = uint8(repmat(mask, [1, 1, 3])) .* 255;

minVisibleCount = 20;

if ~isempty(tracks)

reliableTrackInds = ...

[tracks(:).totalVisibleCount] > minVisibleCount;

reliableTracks = tracks(reliableTrackInds);

if ~isempty(reliableTracks)

bboxes = cat(1, reliableTracks.bbox);

ids = int32([reliableTracks(:).id]);

labels = cellstr(int2str(ids'));

predictedTrackInds = ...

[reliableTracks(:).consecutiveInvisibleCount] > 0;

isPredicted = cell(size(labels));

isPredicted(predictedTrackInds) = {' predicted'};

labels = strcat(labels, isPredicted);

frame = insertObjectAnnotation(frame, 'rectangle', ...

bboxes, labels);

mask = insertObjectAnnotation(mask, 'rectangle', ...

bboxes, labels);

end

end

obj.maskPlayer.step(mask);

obj.videoPlayer.step(frame);

end

end