运动目标检测MATLAB代码实现——差分法、GMM、ViBe
此文章会使用简练而清晰的语言描述下列三种算法，分别为帧间差分法(Temporal Difference)、混合高斯法(GMM)、ViBe算法，并使用MATLAB进行实现。
帧间差分法
帧间差分法分为两种：二帧差分法和三帧差分法，原理类似，这里一并介绍。
摄像机采集的视频序列具有连续性的特点，如果场景内没有运动目标，则连续帧的变化很微弱，如果存在运动目标，则连续的帧和帧之间会有明显地变化。帧间差分法就是借鉴了这个思想。该类算法对时间上连续的两帧或三帧图像进行差分运算，不同帧对应的像素点相减，判断灰度差的绝对值，当绝对值超过一定阈值时，即可判断为运动目标，从而实现运动目标的检测功能。
差分示意图
算法整体的流程图如下所示
三帧差分法
下面是二帧差分法的MATLAB代码实现
clear;clc;%读取视频name = 'fanfan.mp4';%文件名字可自行修改v = VideoReader(name);videoSource = vision.VideoFileReader(name,'ImageColorSpace', 'RGB', 'VideoOutputDataType', 'uint8');n = v.NumberOfFrames;frame_last = rgb2gray(step(videoSource));frame = step(videoSource);figure('units','normalized','outerposition',[0 0 1 1])for i=1:n
frame = step(videoSource);
frame_now = rgb2gray(frame);
frame_now = medfilt2(frame_now);
frame_diff = abs(frame_now - frame_last);
fgMask = imbinarize(frame_diff);
fgMask = imdilate(fgMask,strel('rectangle',[3 3]));
fgMask = imerode(fgMask,strel('arbitrary',eye(3)));
frame_last = frame_now;
figure(1)
subplot(122);imshow(fgMask);title('Original Image')
subplot(121);imshow(frame);title('Segmentation')end
下面是三帧差分法的MATLAB代码实现
clear;clc%读取视频name = 'fanfan.mp4';%文件名字可自行修改v = VideoReader(name);videoSource = vision.VideoFileReader(name,'ImageColorSpace', 'RGB', 'VideoOutputDataType', 'uint8');n = v.NumberOfFrames;frame_last = rgb2gray(step(videoSource));frame_first = rgb2gray(step(videoSource));frame = step(videoSource);figure('units','normalized','outerposition',[0 0 1 1])for i=1:n
figure(1)
frame = step(videoSource);
subplot(121);imshow(frame);title('Original Image')
frame_second = rgb2gray(frame);
frame = step(videoSource);
frame_third = rgb2gray(frame);
frame_diff1 = abs(frame_second - frame_first);
frame_diff2 = abs(frame_third - frame_second);
fgMask = imbinarize(min(frame_diff1,frame_diff2));
fgMask = imdilate(fgMask,strel('rectangle',[3 3]));
fgMask = imerode(fgMask,strel('arbitrary',eye(3)));
fgMask = imfill(fgMask, 'hole');
frame_first = frame_second;
subplot(122);imshow(fgMask);title('Segmentation')end
效果如下
原图像
二帧差分法
三帧差分法
GMM——高斯混合建模
假设每个像素在时域上符合正态分布，在一定阈值范围内的像素判定为背景，并用来更新模型，不符合该分布的像素即为前景。
算法流程如下：
- 初始化背景模型，令 Ix,y∼N(μ0,σ0)
- 按照如下公式检测像素 Ix,y 属于前景还是背景，其中 λ 为阈值参数 {if|Ix,y−μ(x,y)|<λ×σ(x,y)Ix,y∈bgelseIx,y∈fg
- 按照如下公式更新参数，其中 α 为更新速率 μ(x,y)=(1−α)×μ(x,y)+α×Ix,yσ(x,y)=(1−α)×σ2(x,y)+α×(Ix,y−μ(x,y))2
- 重复步骤2、3直到结束。
下面是GMM算法的MATLAB代码实现
clc;clear;filename = 'car.avi'; %文件名字可自行修改vidObj = VideoReader(filename);I=readFrame(vidObj); % 读入第一帧作为背景帧I = rgb2gray(I);fr_bw = I; [height,width] = size(fr_bw); %求每帧图像大小fg = zeros(height, width); %定义前景和背景矩阵bg_bw = zeros(height, width);
C = 3; % 单高斯模型的个数(通常为3-5)M = 3; % 代表背景的模型个数D = 2.5; % 偏差阈值alpha = 0.01; % 学习率thresh = 0.25; % 前景阈值sd_init = 15; % 初始化标准差w = zeros(height,width,C); % 初始化权重矩阵mean = zeros(height,width,C); % 像素均值sd = zeros(height,width,C); % 像素标准差u_diff = zeros(height,width,C); % 像素与某个高斯模型均值的绝对距离p = alpha/(1/C); % 初始化p变量，用来更新均值和标准差rank = zeros(1,C); %各个高斯分布的优先级（w/sd)
pixel_depth = 8; % 每个像素8bit分辨率pixel_range = 2^pixel_depth -1; % 像素值范围[0,255]figure('units','normalized','outerposition',[0 0 1 1])for i=1:height
for j=1:width
for k=1:C
mean(i,j,k) = rand*pixel_range; %初始化第k个高斯分布的均值
w(i,j,k) = 1/C; % 初始化第k个高斯分布的权重
sd(i,j,k) = sd_init; % 初始化第k个高斯分布的标准差
end
endend
while hasFrame(vidObj)
fr=readFrame(vidObj); % 依次读入各帧图像
fr_bw =rgb2gray(fr);
% 计算新像素与第m个高斯模型均值的绝对距离
for m=1:C
u_diff(:,:,m) = abs(double(fr_bw) - double(mean(:,:,m)));
end
% 更新高斯模型的参数
for i=1:height
for j=1:width
match = 0; %匹配标记;
for k=1:C
if (abs(u_diff(i,j,k)) <= D*sd(i,j,k)) % 像素与第k个高斯模型匹配
match = 1; %将匹配标记置为1
% 更新权重、均值、标准差、p
w(i,j,k) = (1-alpha)*w(i,j,k) + alpha;
p = alpha/w(i,j,k);
mean(i,j,k) = (1-p)*mean(i,j,k) + p*double(fr_bw(i,j));
sd(i,j,k) = sqrt((1-p)*(sd(i,j,k)^2) + p*((double(fr_bw(i,j)) - mean(i,j,k)))^2);
else % 像素与第k个高斯模型不匹配
w(i,j,k) = (1-alpha)*w(i,j,k); %略微减少权重
end
end
bg_bw(i,j)=0;
for k=1:C
bg_bw(i,j) = bg_bw(i,j)+ mean(i,j,k)*w(i,j,k);
end
% 像素值与任一高斯模型都不匹配，则创建新的模型
if (match == 0)
[min_w, min_w_index] = min(w(i,j,:)); %寻找最小权重
mean(i,j,min_w_index) = double(fr_bw(i,j));%初始化均值为当前观测像素的均值
sd(i,j,min_w_index) = sd_init; %初始化标准差为6
end
rank = w(i,j,:)./sd(i,j,:); % 计算模型优先级
rank_ind = [1:1:C];%优先级索引
% 计算前景
fg(i,j) = 0;
while ((match == 0)&&(k<=M))
if (abs(u_diff(i,j,rank_ind(k))) <= D*sd(i,j,rank_ind(k)))% 像素与第k个高斯模型匹配
fg(i,j) = 0; %该像素为背景，置为黑色
else
fg(i,j) = 255; %否则为前景，置为白色
end
k = k+1;
end
end
end
figure(1)
subplot(1,2,1),imshow(fr);title('Original Image') %显示最后一帧图像
fg = imdilate(fg,strel('rectangle',[6 6]));
fg = imerode(fg,strel('arbitrary',eye(5)));
subplot(1,2,2),imshow(fg);title('Segmentation') %显示前景end
效果如下
GMM
ViBe算法——随机背景更新策略
像素的变化存在不确定性，很难用一个固定的模型来表征。ViBe算法基于这样一个假设：在无法确定像素变化的模型时，随机模型在一定程度上更适合模拟像素变化的不确定性。
记图片中的每一个像素为 v(x) ，样本集中的值为 vi ，样本的个数为 N ，那么样本可以表示为： M(x)={v1,v2,…,vN}
为了比较像素点与本集合中的值，定义 SR(v(x)) 是以为中心，半径为 R 的区域，若样本集中的值与当前像素点的距离在该范围内的个数大于预先规定的阈值 #min ，则将该像素点分类为背景点，即如下公式
{if#{SR(v(x))}>#minv(x)∈bgelsev(x)∈fg
算法流程如下
1.模型初始化
对于一个像素点，结合周围像素的分布特性，随机选取 N 次作为它的模型样本值。
8邻域分布抽取
假设视频每一帧的分辨率为 X×Y ，那么一共需要选取 X×Y×N 模型样本值。
2.背景判定
对于一个像素点 (x,y) ，计算其和背景样本值的距离 d ，判断它是否在 SR(v(x)) 圆内
颜色相近即距离d比较小
计算得出一共有 #SR(v(x)) 个点在阈值内，若 #SR(v(x)) 大于 #min ，则判定这个像素点是背景，否则判定为前景。
3.模型更新
对背景模型的更新就是使得背景模型能够适应背景的不断变化，比如光照的变化，背景物体的变更。
每一个背景点 (x,y) 有 1φ 的概率去更新自己的模型样本值，同时也有 1φ 的概率去更新它的邻居点 (x′,y′) 的模型样本值。
更新邻居的样本值利用了像素值的空间传播特性，背景模型逐渐向外扩散，同时当前景点计数达到临界值时也将其变为背景，并有 1φ 的概率去更新自己的模型样本值。
代码如下，一共有4个程序
%主程序main.mclear; clc; %参数param.numberOfSamples = 20;param.matchingThreshold = 20;param.matchingNumber = 2;param.updateFactor = 5;param.numberOfHistoryImages = 20;param.lastHistoryImageSwapped = 0;
filename = 'car.avi';vidObj = VideoReader(filename);
firstFrame = true;height = vidObj.Height;width = vidObj.Width;
param.height = height;param.width = width;
figure('units','normalized','outerposition',[0 0 1 1])while hasFrame(vidObj)
vidFrame = readFrame(vidObj);
figure(1)
subplot(121),imshow(vidFrame), title('Original Image');
vidFrame = rgb2gray(vidFrame);
vidFrame = double(vidFrame);
if firstFrame
firstFrame = false;
initViBe;
end
segmentationMap = vibeSegmentation(vidFrame, historyImages, historyBuffer, param);
[historyImages, historyBuffer] = vibeUpdate(vidFrame, segmentationMap, historyImages, historyBuffer, param, ...
jump, neighborX, neighborY, position);
segmentationMap = medfilt2(segmentationMap);
subplot(122),imshow(segmentationMap), title('Segmentation');end
%initViBe.m%% ParametersnumberOfSamples = param.numberOfSamples;matchingThreshold = param.matchingThreshold;matchingNumber = param.matchingNumber;updateFactor = param.updateFactor;numberOfHistoryImages = param.numberOfHistoryImages;
%% Initialize ViBehistoryImages = cell(1, numberOfHistoryImages);for ii = 1:length(historyImages)
historyImages{ii} = vidFrame;end
historyBuffer = cell(1, numberOfSamples - numberOfHistoryImages);for ii = 1:length(historyBuffer)
historyBuffer{ii} = vidFrame + double(floor(rand(height, width))*20 - 10);end
%% Random Partsize = 2*max(height, width) + 1;jump = floor(rand(1, size)*2*updateFactor) + 1;neighborX = floor(rand(1, size)*3) - 1;neighborY = floor(rand(1, size)*3) - 1;position = floor(rand(1, size)*numberOfSamples) + 1;
%vibeSegmentation.mfunction segmentationMap = vibeSegmentation(buffer, historyImages, historyBuffer, param) %% Parameters
height = param.height;
width = param.width;
numberOfSamples = param.numberOfSamples;
matchingThreshold = param.matchingThreshold;
matchingNumber = param.matchingNumber;
numberOfHistoryImages = param.numberOfHistoryImages;
%% Segmentation
segmentationMap = uint8(ones(height, width)*(matchingNumber - 1));
% First and Second history Image structure
distance1 = abs(buffer - historyImages{1}) <= matchingThreshold;
distance2 = abs(buffer - historyImages{2}) <= matchingThreshold;
for ii = 1:height
for jj = 1:width
if ~distance1(ii, jj)
segmentationMap(ii, jj) = matchingNumber;
end
if distance2(ii, jj)
segmentationMap(ii, jj) = segmentationMap(ii, jj) - 1;
end
end
end
% match the image and samples
numberOfTests = numberOfSamples - numberOfHistoryImages;
for kk = 1:numberOfTests
distance3 = uint8(abs(buffer - historyBuffer{kk}) <= matchingThreshold);
segmentationMap = segmentationMap - distance3;
end
segmentationMap = uint8(segmentationMap*255);end
%vibeUpdate.mfunction [historyImages, historyBuffer] = vibeUpdate(buffer, updatingMask, historyImages, historyBuffer, param, ...
jump, neighborX, neighborY, position)
%% Parameters
height = param.height;
width = param.width;
numberOfHistoryImages = param.numberOfHistoryImages;
%% Update Model
for indY = 2:height - 1
shift = floor(rand()*width) + 1;
indX = jump(shift) + 1;
while indX < width
if updatingMask(indY, indX) == 0
value = buffer(indY, indX);
if position(shift) <= numberOfHistoryImages
historyImages{position(shift)}(indY, indX) = value;
historyImages{position(shift)}...
(indY + neighborY(shift), indX + neighborX(shift)) = value;
else
pos = position(shift) - numberOfHistoryImages;
historyBuffer{pos}(indY, indX) = value;
historyBuffer{pos}...
(indY + neighborY(shift), indX + neighborX(shift)) = value;
end
end
shift = shift + 1;
indX = indX + jump(shift);
end
end
end
效果如下
1620
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
