基于帧图像序列的目标检测与跟踪MATLAB实现

基于帧图像序列的目标检测与跟踪MATLAB程序。该系统实现了从视频读取、背景建模、目标检测、目标跟踪到结果可视化的全流程处理。

%% 基于帧图像序列的目标检测与跟踪系统
% 描述: 实现目标检测、跟踪与轨迹分析的完整流程

%% 主程序入口
function objectDetectionAndTracking()
    % 清空环境
    clear; close all; clc;
    
    % 参数设置
    params = struct();
    params.videoSource = 'pedestrians.mp4'; % 视频文件路径
    params.useWebcam = false;                % 是否使用摄像头
    params.samplingRate = 5;                 % 采样率（每N帧处理一次）
    params.minBlobArea = 500;                % 最小目标面积
    params.maxBlobArea = 20000;              % 最大目标面积
    params.minTrackLength = 10;              % 最小轨迹长度
    params.maxTrackGap = 5;                  % 最大轨迹中断帧数
    params.kalmanProcessNoise = 1e-2;        % 卡尔曼过程噪声
    params.kalmanMeasurementNoise = 1e-1;    % 卡尔曼测量噪声
    params.iouThreshold = 0.3;               % IOU匹配阈值
    params.displayResults = true;            % 是否显示结果
    params.saveResults = true;               % 是否保存结果
    
    % 初始化系统
    trackerSystem = initializeTrackingSystem(params);
    
    % 处理视频序列
    processVideoSequence(trackerSystem);
    
    % 分析结果
    analyzeTrackingResults(trackerSystem);
    
    % 保存结果
    if params.saveResults
        saveTrackingResults(trackerSystem);
    end
end

%% 初始化跟踪系统
function system = initializeTrackingSystem(params)
    % 创建系统结构体
    system.params = params;
    system.tracks = {};              % 存储所有轨迹
    system.nextTrackId = 1;          % 下一个轨迹ID
    system.frameCount = 0;           % 已处理帧数
    system.detectionHistory = {};    % 检测结果历史
    system.trackingHistory = {};     % 跟踪结果历史
    
    % 初始化背景建模器
    system.backgroundModel = vision.ForegroundDetector(...
        'NumGaussians', 5, ...
        'NumTrainingFrames', 100, ...
        'MinimumBackgroundRatio', 0.7);
    
    % 初始化卡尔曼滤波器参数
    system.kalmanParams = struct(...
        'A', [1 0 0 0 1 0; 0 1 0 0 0 1; 0 0 1 0 0 0; 0 0 0 1 0 0; 0 0 0 0 1 0; 0 0 0 0 0 1], ... % 状态转移矩阵
        'H', [1 0 0 0 0 0; 0 1 0 0 0 0], ... % 观测矩阵
        'Q', diag([params.kalmanProcessNoise*ones(1,4), params.kalmanProcessNoise*ones(1,2)]), ... % 过程噪声协方差
        'R', eye(2)*params.kalmanMeasurementNoise, ... % 测量噪声协方差
        'P', eye(6)*10, ... % 初始估计误差协方差
        'B', []); % 控制输入矩阵（无控制输入）
    
    % 创建显示窗口
    if params.displayResults
        system.fig = figure('Name', '目标检测与跟踪系统', 'NumberTitle', 'off', ...
            'Position', [100, 100, 1200, 600]);
    end
    
    % 视频源初始化
    if params.useWebcam
        system.videoSource = webcam;
        system.videoInfo = imaqhwinfo(system.videoSource);
        system.isLiveVideo = true;
    else
        if exist(params.videoSource, 'file')
            system.videoSource = VideoReader(params.videoSource);
            system.isLiveVideo = false;
        else
            error('视频文件不存在: %s', params.videoSource);
        end
    end
    
    fprintf('目标检测与跟踪系统已初始化\n');
    fprintf('视频源: %s\n', params.videoSource);
    fprintf('参数设置:\n');
    disp(params);
end

%% 处理视频序列
function processVideoSequence(system)
    % 主处理循环
    while true
        % 获取下一帧
        if system.isLiveVideo
            frame = snapshot(system.videoSource);
        else
            if hasFrame(system.videoSource)
                frame = readFrame(system.videoSource);
            else
                break; % 视频结束
            end
        end
        
        % 更新帧计数
        system.frameCount = system.frameCount + 1;
        
        % 按采样率处理帧
        if mod(system.frameCount, system.params.samplingRate) == 0
            fprintf('处理帧: %d\n', system.frameCount);
            
            % 1. 背景建模与前景检测
            foregroundMask = detectForeground(system, frame);
            
            % 2. 目标检测（连通区域分析）
            detections = detectObjects(system, foregroundMask, frame);
            
            % 3. 目标跟踪
            tracks = trackObjects(system, detections, frame);
            
            % 4. 更新系统状态
            system.detectionHistory{system.frameCount} = detections;
            system.trackingHistory{system.frameCount} = tracks;
            system.tracks = tracks; % 更新当前轨迹
            
            % 5. 可视化结果
            if system.params.displayResults
                visualizeResults(system, frame, detections, tracks);
            end
        end
        
        % 退出条件（实时视频限制帧数）
        if system.isLiveVideo && system.frameCount > 300
            break;
        end
    end
    
    fprintf('视频处理完成，共处理 %d 帧\n', system.frameCount);
end

%% 前景检测
function foregroundMask = detectForeground(system, frame)
    % 转换为灰度图像
    if size(frame, 3) == 3
        grayFrame = rgb2gray(frame);
    else
        grayFrame = frame;
    end
    
    % 使用背景建模器检测前景
    foregroundMask = step(system.backgroundModel, grayFrame);
    
    % 形态学操作去除噪声
    seOpen = strel('disk', 2);
    seClose = strel('rectangle', [5, 5]);
    foregroundMask = imopen(foregroundMask, seOpen);
    foregroundMask = imclose(foregroundMask, seClose);
    
    % 填充小孔
    foregroundMask = imfill(foregroundMask, 'holes');
end

%% 目标检测（连通区域分析）
function detections = detectObjects(system, foregroundMask, frame)
    % 连通区域分析
    cc = bwconncomp(foregroundMask);
    stats = regionprops(cc, 'BoundingBox', 'Area', 'Centroid', 'PixelIdxList');
    
    % 筛选目标
    detections = [];
    for i = 1:length(stats)
        area = stats(i).Area;
        if area >= system.params.minBlobArea && area <= system.params.maxBlobArea
            bbox = stats(i).BoundingBox;
            centroid = stats(i).Centroid;
            
            % 计算边界框宽高比
            aspectRatio = bbox(3)/bbox(4);
            if aspectRatio > 0.2 && aspectRatio < 5.0 % 合理的宽高比范围
                detection = struct(...
                    'id', i, ...
                    'bbox', bbox, ...
                    'centroid', centroid, ...
                    'area', area, ...
                    'frame', system.frameCount, ...
                    'timestamp', tic);
                detections = [detections; detection];
            end
        end
    end
end

%% 目标跟踪
function tracks = trackObjects(system, detections, frame)
    % 预测现有轨迹的新位置
    predictedTracks = predictTracks(system);
    
    % 关联检测与轨迹
    [matchedPairs, unmatchedDetections, unmatchedTracks] = associateDetectionsToTracks(...
        system, detections, predictedTracks);
    
    % 更新匹配的轨迹
    updatedTracks = updateMatchedTracks(system, matchedPairs, detections, predictedTracks);
    
    % 处理未匹配的轨迹（标记为丢失）
    lostTracks = handleUnmatchedTracks(system, unmatchedTracks, predictedTracks);
    
    % 创建新轨迹
    newTracks = createNewTracks(system, unmatchedDetections, detections);
    
    % 合并所有轨迹
    tracks = [updatedTracks; lostTracks; newTracks];
    
    % 移除旧轨迹
    tracks = removeOldTracks(system, tracks);
end

%% 预测轨迹位置（卡尔曼滤波）
function predictedTracks = predictTracks(system)
    predictedTracks = {};
    for i = 1:length(system.tracks)
        track = system.tracks{i};
        
        % 初始化卡尔曼滤波器（如果需要）
        if isempty(track.kalmanFilter)
            track.kalmanFilter = configureKalmanFilter(...
                'ConstantVelocity', ...
                track.centroid, ...
                [100, 30], ... % 初始位置不确定性
                system.params.kalmanMeasurementNoise, ...
                system.params.kalmanProcessNoise);
        end
        
        % 预测下一步位置
        predict(track.kalmanFilter);
        predictedPos = track.kalmanFilter.State(1:2);
        
        % 更新轨迹预测位置
        track.predictedCentroid = predictedPos';
        track.age = track.age + 1;
        track.timeSinceUpdate = track.timeSinceUpdate + 1;
        track.history = [track.history; predictedPos'];
        
        predictedTracks{end+1} = track;
    end
end

%% 配置卡尔曼滤波器
function kalmanFilter = configureKalmanFilter(model, initialLocation, initialEstimateError, measurementNoise, processNoise)
    % 创建卡尔曼滤波器
    kalmanFilter = vision.KalmanFilter(...
        'State', [initialLocation, zeros(2,1)]', ...
        'StateTransitionModel', [1 0 1 0; 0 1 0 1; 0 0 1 0; 0 0 0 1], ...
        'MeasurementModel', [1 0 0 0; 0 1 0 0], ...
        'ProcessNoise', eye(4)*processNoise, ...
        'MeasurementNoise', eye(2)*measurementNoise, ...
        'InitialStateEstimateError', eye(4)*initialEstimateError(1), ...
        'InitialErrorCovariance', eye(4)*initialEstimateError(2));
end

%% 关联检测与轨迹（IOU匹配）
function [matchedPairs, unmatchedDetections, unmatchedTracks] = associateDetectionsToTracks(system, detections, tracks)
    % 初始化输出
    matchedPairs = [];
    unmatchedDetections = 1:length(detections);
    unmatchedTracks = 1:length(tracks);
    
    % 如果没有检测或轨迹，直接返回
    if isempty(detections) || isempty(tracks)
        return;
    end
    
    % 计算代价矩阵（IOU距离）
    costMatrix = zeros(length(tracks), length(detections));
    for i = 1:length(tracks)
        track = tracks{i};
        for j = 1:length(detections)
            detection = detections(j);
            costMatrix(i, j) = 1 - computeIOU(track.predictedBbox, detection.bbox);
        end
    end
    
    % 使用匈牙利算法求解最优匹配
    [assignment, cost] = assignDetectionsToTracks(costMatrix, system.params.iouThreshold);
    
    % 解析匹配结果
    matchedPairs = [];
    for i = 1:size(assignment, 1)
        trackIdx = assignment(i, 1);
        detIdx = assignment(i, 2);
        
        if trackIdx > 0 && detIdx > 0
            matchedPairs = [matchedPairs; trackIdx, detIdx];
            
            % 从未匹配列表中移除
            unmatchedTracks(unmatchedTracks == trackIdx) = [];
            unmatchedDetections(unmatchedDetections == detIdx) = [];
        end
    end
end

%% 计算IOU（交并比）
function iou = computeIOU(bbox1, bbox2)
    % 计算交集区域
    x1_min = bbox1(1); y1_min = bbox1(2);
    x1_max = bbox1(1) + bbox1(3); y1_max = bbox1(2) + bbox1(4);
    
    x2_min = bbox2(1); y2_min = bbox2(2);
    x2_max = bbox2(1) + bbox2(3); y2_max = bbox2(2) + bbox2(4);
    
    inter_x_min = max(x1_min, x2_min);
    inter_y_min = max(y1_min, y2_min);
    inter_x_max = min(x1_max, x2_max);
    inter_y_max = min(y1_max, y2_max);
    
    inter_width = max(0, inter_x_max - inter_x_min);
    inter_height = max(0, inter_y_max - inter_y_min);
    intersection = inter_width * inter_height;
    
    % 计算并集区域
    area1 = bbox1(3) * bbox1(4);
    area2 = bbox2(3) * bbox2(4);
    union = area1 + area2 - intersection;
    
    % 计算IOU
    iou = intersection / union;
end

%% 更新匹配的轨迹
function updatedTracks = updateMatchedTracks(system, matchedPairs, detections, tracks)
    updatedTracks = {};
    for i = 1:size(matchedPairs, 1)
        trackIdx = matchedPairs(i, 1);
        detIdx = matchedPairs(i, 2);
        
        track = tracks{trackIdx};
        detection = detections(detIdx);
        
        % 更新卡尔曼滤波器
        correct(track.kalmanFilter, detection.centroid');
        
        % 更新轨迹状态
        track.centroid = detection.centroid;
        track.bbox = detection.bbox;
        track.predictedBbox = detection.bbox; % 简化处理
        track.age = 1;
        track.timeSinceUpdate = 0;
        track.hits = track.hits + 1;
        track.history = [track.history; detection.centroid];
        
        updatedTracks{end+1} = track;
    end
end

%% 处理未匹配的轨迹
function lostTracks = handleUnmatchedTracks(system, unmatchedTracks, tracks)
    lostTracks = {};
    for i = 1:length(unmatchedTracks)
        trackIdx = unmatchedTracks(i);
        track = tracks{trackIdx};
        
        % 标记为丢失
        track.timeSinceUpdate = track.timeSinceUpdate + 1;
        
        % 如果丢失时间过长，标记为不可见
        if track.timeSinceUpdate <= system.params.maxTrackGap
            lostTracks{end+1} = track;
        end
    end
end

%% 创建新轨迹
function newTracks = createNewTracks(system, unmatchedDetections, detections)
    newTracks = {};
    for i = 1:length(unmatchedDetections)
        detIdx = unmatchedDetections(i);
        detection = detections(detIdx);
        
        % 创建新轨迹
        track = struct(...
            'id', system.nextTrackId, ...
            'centroid', detection.centroid, ...
            'bbox', detection.bbox, ...
            'predictedBbox', detection.bbox, ...
            'kalmanFilter', [], ...
            'age', 1, ...
            'timeSinceUpdate', 0, ...
            'hits', 1, ...
            'history', detection.centroid, ...
            'color', rand(1, 3)); % 随机颜色用于可视化
        
        newTracks{end+1} = track;
        system.nextTrackId = system.nextTrackId + 1;
    end
end

%% 移除旧轨迹
function tracks = removeOldTracks(system, tracks)
    validTracks = {};
    for i = 1:length(tracks)
        track = tracks{i};
        if track.timeSinceUpdate <= system.params.maxTrackGap || track.age < system.params.minTrackLength
            validTracks{end+1} = track;
        end
    end
    tracks = validTracks;
end

%% 可视化结果
function visualizeResults(system, frame, detections, tracks)
    % 复制帧用于绘制
    displayFrame = frame;
    
    % 绘制检测结果
    for i = 1:length(detections)
        bbox = detections(i).bbox;
        centroid = detections(i).centroid;
        
        % 绘制边界框
        displayFrame = insertShape(displayFrame, 'Rectangle', bbox, ...
            'Color', 'yellow', 'LineWidth', 2);
        
        % 绘制质心
        displayFrame = insertMarker(displayFrame, centroid, ...
            'o', 'Color', 'cyan', 'Size', 8);
    end
    
    % 绘制轨迹
    for i = 1:length(tracks)
        track = tracks{i};
        
        % 绘制边界框
        displayFrame = insertShape(displayFrame, 'Rectangle', track.bbox, ...
            'Color', track.color, 'LineWidth', 2);
        
        % 绘制轨迹历史
        if size(track.history, 1) > 1
            for j = 2:size(track.history, 1)
                pt1 = track.history(j-1, :);
                pt2 = track.history(j, :);
                displayFrame = insertShape(displayFrame, 'Line', ...
                    [pt1(1), pt1(2), pt2(1), pt2(2)], ...
                    'Color', track.color, 'LineWidth', 2);
            end
        end
        
        % 绘制轨迹ID
        displayFrame = insertText(displayFrame, track.bbox(1:2), ...
            num2str(track.id), 'FontSize', 14, ...
            'BoxColor', track.color, 'TextColor', 'white');
    end
    
    % 显示帧信息
    infoStr = sprintf('帧: %d\n检测目标: %d\n跟踪目标: %d', ...
        system.frameCount, length(detections), length(tracks));
    displayFrame = insertText(displayFrame, [10, 10], infoStr, ...
        'FontSize', 14, 'BoxColor', 'black', 'TextColor', 'white');
    
    % 显示图像
    if ishandle(system.fig)
        figure(system.fig);
        imshow(displayFrame);
        title(sprintf('目标检测与跟踪 (帧 %d)', system.frameCount));
        drawnow;
    end
end

%% 分析跟踪结果
function analyzeTrackingResults(system)
    % 计算统计信息
    totalTracks = system.nextTrackId - 1;
    activeTracks = length(system.tracks);
    avgTrackLength = 0;
    maxTrackLength = 0;
    
    % 计算轨迹长度
    trackLengths = [];
    for i = 1:length(system.trackingHistory)
        if ~isempty(system.trackingHistory{i})
            for j = 1:length(system.trackingHistory{i})
                track = system.trackingHistory{i}{j};
                if isfield(track, 'hits')
                    trackLengths(end+1) = track.hits;
                end
            end
        end
    end
    
    if ~isempty(trackLengths)
        avgTrackLength = mean(trackLengths);
        maxTrackLength = max(trackLengths);
    end
    
    % 显示统计信息
    fprintf('\n===== 跟踪结果分析 =====\n');
    fprintf('总轨迹数: %d\n', totalTracks);
    fprintf('活动轨迹数: %d\n', activeTracks);
    fprintf('平均轨迹长度: %.2f 帧\n', avgTrackLength);
    fprintf('最长轨迹长度: %d 帧\n', maxTrackLength);
    
    % 绘制轨迹长度分布
    if ~isempty(trackLengths)
        figure;
        histogram(trackLengths, 20);
        xlabel('轨迹长度 (帧)');
        ylabel('数量');
        title('轨迹长度分布');
        grid on;
    end
    
    % 绘制目标运动轨迹
    plotTrajectories(system);
end

%% 绘制目标运动轨迹
function plotTrajectories(system)
    figure('Name', '目标运动轨迹', 'NumberTitle', 'off', 'Position', [100, 100, 1000, 800]);
    
    % 收集所有轨迹点
    allTrajectories = {};
    for i = 1:length(system.trackingHistory)
        if ~isempty(system.trackingHistory{i})
            for j = 1:length(system.trackingHistory{i})
                track = system.trackingHistory{i}{j};
                if isfield(track, 'history') && size(track.history, 1) > 1
                    allTrajectories{end+1} = struct(...
                        'id', track.id, ...
                        'color', track.color, ...
                        'points', track.history);
                end
            end
        end
    end
    
    % 绘制轨迹
    if ~isempty(allTrajectories)
        for i = 1:length(allTrajectories)
            traj = allTrajectories{i};
            plot(traj.points(:,1), traj.points(:,2), ...
                'Color', traj.color, 'LineWidth', 1.5);
            hold on;
            plot(traj.points(1,1), traj.points(1,2), 'o', ...
                'Color', traj.color, 'MarkerSize', 8, 'MarkerFaceColor', traj.color);
            plot(traj.points(end,1), traj.points(end,2), 's', ...
                'Color', traj.color, 'MarkerSize', 8, 'MarkerFaceColor', traj.color);
        end
        
        xlabel('X 坐标 (像素)');
        ylabel('Y 坐标 (像素)');
        title('目标运动轨迹');
        grid on;
        axis equal;
        legend('轨迹', '起点', '终点', 'Location', 'bestoutside');
    else
        text(0.5, 0.5, '无有效轨迹数据', 'HorizontalAlignment', 'center');
    end
end

%% 保存跟踪结果
function saveTrackingResults(system)
    % 准备数据
    results = struct();
    results.params = system.params;
    results.frameCount = system.frameCount;
    results.totalTracks = system.nextTrackId - 1;
    results.tracks = system.tracks;
    results.detectionHistory = system.detectionHistory;
    results.trackingHistory = system.trackingHistory;
    
    % 保存为MAT文件
    save('tracking_results.mat', 'results');
    
    % 导出轨迹数据为CSV
    exportTrajectoriesToCSV(system);
    
    fprintf('跟踪结果已保存:\n');
    fprintf('  - tracking_results.mat (完整数据)\n');
    fprintf('  - trajectories.csv (轨迹数据)\n');
end

%% 导出轨迹数据到CSV
function exportTrajectoriesToCSV(system)
    % 创建表格数据
    data = [];
    
    for frameIdx = 1:length(system.trackingHistory)
        if ~isempty(system.trackingHistory{frameIdx})
            for trackIdx = 1:length(system.trackingHistory{frameIdx})
                track = system.trackingHistory{frameIdx}{trackIdx};
                if isfield(track, 'id') && isfield(track, 'centroid')
                    data = [data; struct(...
                        'Frame', frameIdx, ...
                        'TrackID', track.id, ...
                        'X', track.centroid(1), ...
                        'Y', track.centroid(2), ...
                        'Width', track.bbox(3), ...
                        'Height', track.bbox(4), ...
                        'Age', track.age, ...
                        'Hits', track.hits, ...
                        'TimeSinceUpdate', track.timeSinceUpdate)];
                end
            end
        end
    end
    
    % 写入CSV文件
    if ~isempty(data)
        % 转换为表格
        T = struct2table(data);
        
        % 写入CSV
        writetable(T, 'trajectories.csv');
    end
end

%% 辅助函数：使用匈牙利算法分配检测与轨迹
function [assignments, unassignedTracks, unassignedDetections] = assignDetectionsToTracks(costMatrix, costOfNonAssignment)
    % 简化版的匈牙利算法实现
    numTracks = size(costMatrix, 1);
    numDets = size(costMatrix, 2);
    
    % 初始化输出
    assignments = zeros(0, 2);
    unassignedTracks = 1:numTracks;
    unassignedDets = 1:numDets;
    
    % 如果没有轨迹或检测，直接返回
    if numTracks == 0 || numDets == 0
        return;
    end
    
    % 构建扩展代价矩阵
    costMatrixExtended = [costMatrix, ones(numTracks, 1)*costOfNonAssignment];
    costMatrixExtended = [costMatrixExtended; ones(1, numDets+1)*costOfNonAssignment];
    
    % 使用匈牙利算法求解
    [assignment, cost] = hungarian(costMatrixExtended);
    
    % 解析结果
    for i = 1:size(assignment, 1)
        if assignment(i, 1) <= numTracks && assignment(i, 2) <= numDets
            trackIdx = assignment(i, 1);
            detIdx = assignment(i, 2);
            
            if costMatrix(trackIdx, detIdx) <= costOfNonAssignment
                assignments = [assignments; trackIdx, detIdx];
                
                % 从未分配列表中移除
                unassignedTracks(unassignedTracks == trackIdx) = [];
                unassignedDets(unassignedDets == detIdx) = [];
            end
        end
    end
    
    % 剩余的是未分配的轨迹和检测
    unassignedTracks = [unassignedTracks; (numTracks+1):size(costMatrixExtended, 1)];
    unassignedDets = [unassignedDets; (numDets+1):size(costMatrixExtended, 2)];
    unassignedTracks(unassignedTracks > numTracks) = [];
    unassignedDets(unassignedDets > numDets) = [];
end

%% 匈牙利算法实现（简化版）
function [assignment, cost] = hungarian(costMatrix)
    % 简化实现 - 实际应用中应使用更高效版本
    [numRows, numCols] = size(costMatrix);
    assignment = zeros(min(numRows, numCols), 2);
    cost = 0;
    
    % 贪心匹配（实际应用中应使用完整匈牙利算法）
    usedRows = false(numRows, 1);
    usedCols = false(numCols, 1);
    
    for k = 1:min(numRows, numCols)
        minCost = inf;
        minRow = 0;
        minCol = 0;
        
        for i = 1:numRows
            if ~usedRows(i)
                for j = 1:numCols
                    if ~usedCols(j) && costMatrix(i, j) < minCost
                        minCost = costMatrix(i, j);
                        minRow = i;
                        minCol = j;
                    end
                end
            end
        end
        
        if minCost < inf
            assignment(k, :) = [minRow, minCol];
            cost = cost + minCost;
            usedRows(minRow) = true;
            usedCols(minCol) = true;
        end
    end
end

%% 演示函数
function demoObjectTracking()
    % 创建示例视频
    createSampleVideo();
    
    % 运行目标检测与跟踪
    objectDetectionAndTracking();
end

%% 创建示例视频
function createSampleVideo()
    % 创建视频写入器
    videoFile = 'pedestrians.mp4';
    writer = VideoWriter(videoFile, 'MPEG-4');
    writer.FrameRate = 10;
    open(writer);
    
    % 创建背景
    bg = zeros(480, 640, 3, 'uint8');
    bg(:,:,1) = 100; % R
    bg(:,:,2) = 120; % G
    bg(:,:,3) = 150; % B
    
    % 添加一些静态对象
    bg(200:300, 100:150, :) = 200; % 矩形
    bg(150:200, 400:500, 1) = 200; % 红色区域
    
    % 创建移动目标
    for frameIdx = 1:100
        frame = bg;
        
        % 目标1：从左向右移动
        x1 = 50 + frameIdx * 5;
        y1 = 200;
        if x1 < 600
            frame(y1:y1+50, x1:x1+30, 1) = 255; % 红色目标
        end
        
        % 目标2：从上向下移动
        x2 = 300;
        y2 = 50 + frameIdx * 3;
        if y2 < 400
            frame(y2:y2+40, x2:x2+40, 2) = 255; % 绿色目标
        end
        
        % 目标3：对角线移动
        x3 = 100 + frameIdx * 4;
        y3 = 100 + frameIdx * 2;
        if x3 < 550 && y3 < 380
            frame(y3:y3+60, x3:x3+30, 3) = 255; % 蓝色目标
        end
        
        % 添加噪声
        noise = randn(480, 640, 3) * 10;
        frame = uint8(double(frame) + noise);
        
        % 写入帧
        writeVideo(writer, frame);
    end
    
    % 关闭写入器
    close(writer);
    disp(['示例视频已创建: ', videoFile]);
end

程序功能说明

1. 系统初始化模块

• 参数配置（视频源、检测阈值、跟踪参数等）
• 背景建模器初始化（高斯混合模型）
• 卡尔曼滤波器参数配置
• 视频源初始化（文件或摄像头）

2. 视频处理流水线

1. 前景检测：使用高斯混合模型分离前景和背景

2. 目标检测：连通区域分析提取目标边界框和质心

3. 目标跟踪：

   • 卡尔曼滤波预测目标位置
   • IOU匹配关联检测与轨迹
   • 更新轨迹状态
   • 创建新轨迹
   • 移除旧轨迹

4. 结果可视化：实时显示检测结果和跟踪轨迹

3. 轨迹分析模块

• 轨迹统计（总数、长度分布）
• 运动轨迹可视化
• 轨迹数据导出（CSV格式）

4. 辅助功能

• 示例视频生成
• 匈牙利算法实现（数据关联）
• IOU计算（交并比）
• 结果保存（MAT文件和CSV）

关键技术实现

1. 背景建模与目标检测

function foregroundMask = detectForeground(system, frame)
    % 转换为灰度图像
    grayFrame = rgb2gray(frame);
    
    % 使用背景建模器检测前景
    foregroundMask = step(system.backgroundModel, grayFrame);
    
    % 形态学操作去除噪声
    seOpen = strel('disk', 2);
    seClose = strel('rectangle', [5, 5]);
    foregroundMask = imopen(foregroundMask, seOpen);
    foregroundMask = imclose(foregroundMask, seClose);
    foregroundMask = imfill(foregroundMask, 'holes');
end

2. 卡尔曼滤波跟踪

function predictedTracks = predictTracks(system)
    for i = 1:length(system.tracks)
        track = system.tracks{i};
        
        % 初始化卡尔曼滤波器
        if isempty(track.kalmanFilter)
            track.kalmanFilter = configureKalmanFilter(...
                'ConstantVelocity', track.centroid, [100, 30], ...
                system.params.kalmanMeasurementNoise, ...
                system.params.kalmanProcessNoise);
        end
        
        % 预测下一步位置
        predict(track.kalmanFilter);
        predictedPos = track.kalmanFilter.State(1:2);
        
        % 更新轨迹预测位置
        track.predictedCentroid = predictedPos';
        track.age = track.age + 1;
        track.timeSinceUpdate = track.timeSinceUpdate + 1;
        track.history = [track.history; predictedPos'];
        
        predictedTracks{end+1} = track;
    end
end

3. 数据关联（IOU匹配）

function [matchedPairs, unmatchedDetections, unmatchedTracks] = associateDetectionsToTracks(system, detections, tracks)
    % 计算代价矩阵（IOU距离）
    costMatrix = zeros(length(tracks), length(detections));
    for i = 1:length(tracks)
        track = tracks{i};
        for j = 1:length(detections)
            detection = detections(j);
            costMatrix(i, j) = 1 - computeIOU(track.predictedBbox, detection.bbox);
        end
    end
    
    % 使用匈牙利算法求解最优匹配
    [assignment, cost] = assignDetectionsToTracks(costMatrix, system.params.iouThreshold);
    
    % 解析匹配结果
    matchedPairs = [];
    for i = 1:size(assignment, 1)
        trackIdx = assignment(i, 1);
        detIdx = assignment(i, 2);
        if trackIdx > 0 && detIdx > 0 && costMatrix(trackIdx, detIdx) <= system.params.iouThreshold
            matchedPairs = [matchedPairs; trackIdx, detIdx];
        end
    end
end

4. 轨迹可视化

function visualizeResults(system, frame, detections, tracks)
    displayFrame = frame;
    
    % 绘制检测结果
    for i = 1:length(detections)
        bbox = detections(i).bbox;
        displayFrame = insertShape(displayFrame, 'Rectangle', bbox, ...
            'Color', 'yellow', 'LineWidth', 2);
    end
    
    % 绘制轨迹
    for i = 1:length(tracks)
        track = tracks{i};
        displayFrame = insertShape(displayFrame, 'Rectangle', track.bbox, ...
            'Color', track.color, 'LineWidth', 2);
        
        % 绘制轨迹历史
        if size(track.history, 1) > 1
            for j = 2:size(track.history, 1)
                pt1 = track.history(j-1, :);
                pt2 = track.history(j, :);
                displayFrame = insertShape(displayFrame, 'Line', ...
                    [pt1(1), pt1(2), pt2(1), pt2(2)], ...
                    'Color', track.color, 'LineWidth', 2);
            end
        end
    end
    
    imshow(displayFrame);
    title(sprintf('目标检测与跟踪 (帧 %d)', system.frameCount));
end

算法原理与数学基础

1. 背景建模（高斯混合模型）

高斯混合模型（GMM）将每个像素表示为多个高斯分布的混合：

$P(xt)=i=1∑Kωi,t⋅N(xt∣μi,t,Σi,t)$

其中K是高斯分量数量，ω是权重，μ是均值，Σ是协方差矩阵。

2. 卡尔曼滤波

卡尔曼滤波是一种递归估计算法，包含预测和更新两个步骤：

预测步骤：

$x^k∣k−1=Fkx^k−1∣k−1$
$Pk∣k−1=FkPk−1∣k−1FkT+Qk$

更新步骤：

$Kk=Pk∣k−1HkT(HkPk∣k−1HkT+Rk)−1$
$x^k∣k=x^k∣k−1+Kk(zk−Hkx^k∣k−1)$
$Pk∣k=(I−KkHk)Pk∣k−1$

3. 数据关联（匈牙利算法）

匈牙利算法解决指派问题，最小化总匹配成本：

$mini=1∑nj=1∑nCijxij$

约束条件：$∑j=1nxij=1, ∑i=1nxij=1$

4. IOU（交并比）

IOU衡量两个边界框的重叠程度：

$IoU=∣B1∪B2∣∣B1∩B2∣$

其中B1和B2是两个边界框。

使用说明

1. 基本使用

% 运行主程序
objectDetectionAndTracking();

% 使用摄像头实时处理
system.params.useWebcam = true;
system.params.videoSource = '';
objectDetectionAndTracking();

2. 参数调整

% 调整检测参数
system.params.minBlobArea = 300;    % 减小最小目标面积
system.params.maxBlobArea = 30000;  % 增大最大目标面积

% 调整跟踪参数
system.params.maxTrackGap = 10;     % 允许更长的轨迹中断
system.params.iouThreshold = 0.5;   % 提高匹配阈值

% 调整卡尔曼滤波参数
system.params.kalmanProcessNoise = 0.01;
system.params.kalmanMeasurementNoise = 0.1;

3. 处理自定义视频

% 指定视频文件
system.params.videoSource = 'my_video.mp4';
system.params.useWebcam = false;

% 调整采样率（处理每N帧）
system.params.samplingRate = 2; % 处理每2帧

参考代码 基于帧图像序列的目标检测与跟踪 www.youwenfan.com/contentalh/96326.html

扩展功能

1. 深度学习目标检测

function detections = detectObjectsWithYOLO(system, frame)
    % 使用预训练的YOLOv4模型
    net = yolov4ObjectDetector('csp-darknet53-coco');
    
    % 检测目标
    [bboxes, scores, labels] = detect(net, frame);
    
    % 转换为统一格式
    detections = [];
    for i = 1:size(bboxes, 1)
        if scores(i) > 0.5 % 置信度阈值
            detection = struct(...
                'id', i, ...
                'bbox', [bboxes(i,1), bboxes(i,2), bboxes(i,3)-bboxes(i,1), bboxes(i,4)-bboxes(i,2)], ...
                'centroid', [(bboxes(i,1)+bboxes(i,3))/2, (bboxes(i,2)+bboxes(i,4))/2], ...
                'area', prod(bboxes(i,3:4)-bboxes(i,1:2)), ...
                'frame', system.frameCount, ...
                'timestamp', tic, ...
                'label', labels(i));
            detections = [detections; detection];
        end
    end
end

2. 多目标跟踪算法（DeepSORT）

function tracks = deepSortTracking(system, detections, tracks)
    % 提取外观特征（使用ReID模型）
    features = extractReIDFeatures(detections);
    
    % 预测轨迹位置（卡尔曼滤波）
    predictedTracks = predictTracks(system);
    
    % 级联匹配（优先考虑近期匹配的轨迹）
    [matches, unmatchedTracks, unmatchedDetections] = cascadeMatching(...
        predictedTracks, detections, features);
    
    % 更新匹配的轨迹
    updatedTracks = updateMatchedTracks(system, matches, detections, predictedTracks);
    
    % 创建新轨迹
    newTracks = createNewTracks(system, unmatchedDetections, detections);
    
    % 合并轨迹
    tracks = [updatedTracks; newTracks];
end

3. 轨迹预测与行为分析

function analyzeBehavior(tracks)
    % 分析轨迹运动模式
    for i = 1:length(tracks)
        track = tracks{i};
        positions = track.history;
        
        % 计算速度和加速度
        velocities = diff(positions);
        accelerations = diff(velocities);
        
        % 识别运动模式
        if std(velocities(:,1)) < 0.5 && std(velocities(:,2)) < 0.5
            behavior = '静止';
        elseif std(velocities(:,1)) > 2*std(velocities(:,2))
            behavior = '水平运动';
        elseif std(velocities(:,2)) > 2*std(velocities(:,1))
            behavior = '垂直运动';
        else
            behavior = '任意方向运动';
        end
        
        % 识别路径类型
        if isStraightLine(positions)
            pathType = '直线';
        elseif isCircular(positions)
            pathType = '圆形';
        else
            pathType = '复杂路径';
        end
        
        fprintf('轨迹 %d: 行为=%s, 路径=%s\n', track.id, behavior, pathType);
    end
end

4. 实时性能优化

function optimizedProcessing(system)
    % 使用GPU加速
    if gpuDeviceCount > 0
        system.backgroundModel.UseGPU = true;
    end
    
    % 降低分辨率
    system.params.processingResolution = [320, 240];
    
    % 跳帧处理
    system.params.skipFrames = 2; % 每3帧处理1帧
    
    % 并行处理
    if isempty(gcp('nocreate'))
        parpool('local', 4);
    end
    parfor frameIdx = 1:totalFrames
        % 处理帧
    end
end

常见问题

1. 目标检测不准确

   • 调整形态学操作参数（结构元素大小）
   • 修改目标面积阈值
   • 使用更高级的检测算法（如YOLO、SSD）

2. 轨迹ID切换频繁

   • 调整IOU匹配阈值
   • 增加卡尔曼滤波器的过程噪声
   • 实现轨迹平滑算法

3. 处理速度慢

   • 降低处理分辨率
   • 跳帧处理
   • 使用GPU加速
   • 优化算法实现

4. 漏检和误检

   • 调整背景建模参数
   • 添加阴影去除步骤
   • 使用多特征融合（颜色、纹理等）

实际应用建议

1. 智能监控系统

   • 异常行为检测（徘徊、奔跑、打架）
   • 人数统计和密度分析
   • 入侵检测

2. 自动驾驶

   • 行人检测和跟踪
   • 车辆检测和跟踪
   • 交通流量分析

3. 体育分析

   • 球员跟踪和热力图
   • 战术分析
   • 动作识别

4. 零售分析

   • 顾客行为分析
   • 货架关注度分析
   • 排队检测