1 简介
基于MATLAB工具对运动目标的检测进行了研究,利用改进的帧差算法对获取的图像做帧差并将帧差图像做自适应阈值分割,然后自适应滤波去燥,最后利用逻辑与运算将两幅二值化图像合并成一幅图像,利用形态学腐蚀和膨胀的检测方法,提取运动物体的轮廓达到检测运动目标的目的.
2 完整代码
function MovingTargetDetectionByMMI()
%Moving Target Detection
%Based on Maximun Mutual Information
%
%读文件
Im1=imread('001.jpg');
Im2=imread('002.jpg');
Im3=imread('003.jpg');
Im1=rgb2gray(Im1);
Im2=rgb2gray(Im2);
Im3=rgb2gray(Im3);
tic;
d12=GetDifferenceImg(Im2,Im1);
d23=GetDifferenceImg(Im2,Im3);
d=d12.*d23;
se =[0 0 0 1 0 0 0
0 0 1 2 1 0 0
0 1 2 4 2 1 0
1 2 4 8 4 2 1
0 1 2 4 2 1 0
0 0 1 2 1 0 0
0 0 0 1 0 0 0];
for i=1:4
d = imfilter(d,se);
end
for i=1:2
d = medfilt2(d,[4 4]);
end
%%d=abs((d12-d23).^0.7);
d=uint8(d/max(max(d))*255);
level = graythresh(d);
BW = im2bw(d,level);
s=regionprops(BW,'BoundingBox');
figure(1)
subplot(2,2,1);
imshow(uint8(d12/max(max(d12))*255));
title('参考帧与前一帧的差值')
subplot(2,2,2);
imshow(uint8(d23/max(max(d23))*255));
title('参考帧与后一帧的差值')
subplot(2,2,3);
imshow(BW);
title('由前后帧得出的差值')
subplot(2,2,4);
imshow(Im2);
%imshow(d);
rectangle('Position',s(1).BoundingBox,'Curvature',[0,0],'LineWidth',2,'LineStyle','--','EdgeColor','r')
title('参考帧与检测结果')
%求相邻两帧重合部分差值主函数
function outImg=GetDifferenceImg(R,F)
[CA1,~,~,CD1]=dwt2(R,'db1');
[CA2,~,~,CD2]=dwt2(F,'db1');
CA1=uint8(CA1);
CA2=uint8(CA2);
fprintf('\n------PSO start\n');
[pa,mi]=PSO(CA1,CA2);
while mi<1.2
[pa,mi]=PSO(CA1,CA2);
end
fprintf('tx:%f ty:%f ang:%f mi:%f\n',pa(1),pa(2),pa(3),mi);
fprintf('------PSO end\n\n');
%pa=[0,0,0];
fprintf('------Powell start\n');
mi_old=0;
while abs(mi-mi_old)>0.01
mi_old=mi;
[pa,mi]=powell(R,F,pa);
end
fprintf('------Powell end\n\n');
time=toc;
fprintf('tx:%.4f ty:%.4f ang:%.2f mi:%f\n',pa(1),pa(2),pa(3),mi);
fprintf('time:%f\n',time);
outImg=GetDifference(pa(1),pa(2),pa(3),R,F);
%figure(6);imshow(outImg);
%求相邻两帧重合部分差值
function outImg=GetDifference(tx,ty,ang,R,F)
[m,n]=size(R);
%
R=im2double(R);
F=im2double(F);
theta=ang*pi/180; %旋转角度转弧度
cx=floor(n/2); %旋转的中心点
cy=floor(m/2);
outImg=zeros(m,n);
for j=1:m
for i=1:n
%参考图像在浮动图像平移后的对应点
% x=i-tx; %列
% y=j-ty; %
x=(i-cx)*cos(theta)-(j-cy)*sin(theta)+cx-tx;
y=(i-cx)*sin(theta)+(j-cy)*cos(theta)+cy-ty;
x1=floor(x);
y1=floor(y);
rval=R(j,i);
%图像重合部分求差
if(x1>=1&&x1<n&&y1>=1&&y1<m)
dy=y1-y;dx=x1-x;
%双线性插值
fval=(F(y1+1,x1)-F(y1,x1))*dy+(F(y1,x1+1)-F(y1,x1))*dx+(F(y1+1,x1+1)+F(y1,x1)-F(y1,x1+1)-F(y1+1,x1))*dy*dx+F(y1,x1);
outImg(j,i)=abs((rval-fval).^0.7*exp(-min([rval fval])/20));
%outImg(j,i)=abs((rval-fval).^2.5/(min([rval fval])).^0.2);
end
end
end
%outImg=uint8(outImg/max(max(outImg))*255);
%双线性插值求互信息
function out=BI_mi(tx,ty,ang,R,F)
[m,n]=size(R);
hist=zeros(256,256);
ha = zeros(1,256);
hb = zeros(1,256);
%归一化到256级灰度
% if max(max(r))~=min(min(r)) %max(max(a))结果是A中最大的元素，max(A)结果一个行向量，元素分别是A的每个列向量的最大的元素
% r = (r-min(min(r)))/(max(max(r))-min(min(r)));
% else
% r = zeros(M,N);
% end
%
% if max(max(f))-min(min(f))
% f = (f-min(min(f)))/(max(max(f))-min(min(f)));
% else
% f = zeros(M,N);
% end
%
% r = double(int16(r*255))+1;
% f = double(int16(f*255))+1;
R=R+1;
F=F+1;
theta=ang*pi/180; %旋转角度转弧度
cx=floor(n/2); %旋转的中心点
cy=floor(m/2);
%求联合概率密度
for j=1:m
for i=1:n
%参考图像在浮动图像平移后的对应点
% x=i-tx; %列
% y=j-ty; %
x=(i-cx)*cos(theta)-(j-cy)*sin(theta)+cx-tx;
y=(i-cx)*sin(theta)+(j-cy)*cos(theta)+cy-ty;
x1=floor(x);
y1=floor(y);
rval=R(j,i);
%图像重合部分求差
if(x1>=1&&x1<n&&y1>=1&&y1<m)
dy=y1-y;dx=x1-x;
%双线性插值
fval=(F(y1+1,x1)-F(y1,x1))*dy+(F(y1,x1+1)-F(y1,x1))*dx+(F(y1+1,x1+1)+F(y1,x1)-F(y1,x1+1)-F(y1+1,x1))*dy*dx+F(y1,x1);
hist(fval,rval)=hist(fval,rval)+1;
end
end
end
%下面求边缘概率密度
for i=1:256
ha(i)=sum(hist(i,:));
hb(i)=sum(hist(:,i));
end
%调用互信息函数
out=MI(hist,ha,hb);
%下面是求互信息的函数
function mi=MI(hist,ha,hb)
%
hsum = sum(sum(hist));
index = find(hist~=0);
p=hist/hsum;
Hab=sum(sum(-p(index).*log(p(index)))); %联合熵
hsum = sum(sum(ha));
index = find(ha~=0);
p = ha/hsum;
Ha = sum(sum(-p(index).*log(p(index)))); %边缘熵
hsum = sum(sum(hb));
index = find(hb~=0);
p = hb/hsum;
Hb = sum(sum(-p(index).*log(p(index)))); %边缘熵
mi = Ha+Hb-Hab;
%粒子群算法
function [pa,mi]=PSO(R,F)
%粒子群算法
%初始化
D=3; %维数
ps=30; %种群规模
VRmin=ones(D,1)*-20; %最小速度
VRmax=ones(D,1)*20; %最大速度
VR=[VRmin,VRmax];
%minmax = 1;
pos=40*rand(ps,D)-20; %随机产生初始位置
vel=8*rand(ps,D)-4; %产生随机速度
% ps=15; %种群规模
% VRmin=ones(D,1)*-10; %最小速度
% VRmax=ones(D,1)*10; %最大速度
% VR=[VRmin,VRmax];
% %minmax = 1;
% pos=20*rand(ps,D)-10; %随机产生初始位置
% vel=4*rand(ps,D)-2; %产生随机速度
%一些参数
maxinterations=20; %最大迭代次数
iw=1; %固定权重
iw1 = 0.9; % 最大惯性权重
iw2 = 0.4;
iwe = 15;
ac1=2;
ac2=2;
flagg=1;
% ergrd=1e-5;
% ergrdep=5; %
% mv=4;%
ergrd=1e-4;
ergrdep=5; %
mv=4;%
%初始个体极值
pbest=pos;
%求初始全局极值
for i=1:ps
p=pos(i,:); %第i个粒子位置
out(i)=BI_mi(p(1),p(2),p(3),R,F); %求函数值
end
pbestval=out; %每个粒子当前函数值
[gbestval,idx]=max(pbestval); %全局最优函数值
gbest=pbest(idx,:); %全局极值
tr(1)=gbestval; %保存当前全局最优函数值
% start PSO iterative procedures
cnt=0; % counter used for updating display according to df in the options
cnt2=0; % counter used for the stopping subroutine based on error convergence
%开始迭代
for i=1:maxinterations
for j=1:ps
if flagg==1 % randomization control, one random set for each particle at each epoch
rannum1=rand(1);
rannum2=rand(1);
end
p=pos(j,:); %第i个粒子位置
out(j)=BI_mi(p(1),p(2),p(3),R,F); %求函数值
e(j)=out(j);
%更新pbest
if pbestval(j)<=e(j);%%%====
pbestval(j)=e(j);
pbest(j,:)=pos(j,:);
end
%更新gbest
[iterbestval,idx1]=max(pbestval);
if gbestval<=iterbestval %%%===
gbestval=iterbestval;
gbest=pbest(idx1,:);
end
tr(i+1)=gbestval;
te=i;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%更新速度，位置
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if i<=iwe %%%===
iwt(i)=((iw2-iw1)/(iwe-1))*(i-1)+iw1; %惯性权重
else
iwt(i)=iw2;
end
% iwt(i)=1;
%%%%%%%%%%速度%%%%%%%%%%%%%%%%%%%%
if flagg==2 % 粒子的每个参数的随机数不同
for dimcnt=1:D
rannum1=rand(1);
rannum2=rand(1);
vel(j,dimcnt)=iwt(i)*vel(j,dimcnt)...
+ac1*rannum1*(pbest(j,dimcnt)-pos(j,dimcnt))...
+ac2*rannum2*(gbest(1,dimcnt)-pos(j,dimcnt));
end
else % 粒子的每个参数的随机数相同
vel(j,:)=iwt(i)*vel(j,:)...
+ac1*rannum1*(pbest(j,:)-pos(j,:))...
+ac2*rannum2*(gbest(1,:)-pos(j,:));
end
% 固定权重
% vel(j,:)=iw*vel(j,:)...
% +ac1*rannum1*(pbest(j,:)-pos(j,:))...
% +ac2*rannum2*(gbest(1,:)-pos(j,:));
%%%%%%%%%%位置%%%%%%%%%%%%%%%%%%%%%%%
pos(j,:)=pos(j,:)+vel(j,:);
%%%%%%%%%速度和位置范围％％％％％％％％％
for dimcnt=1:D
if vel(j,dimcnt)>mv
vel(j,dimcnt)=mv;
end
if vel(j,dimcnt)<-mv
vel(j,dimcnt)=-mv;
end
if pos(j,dimcnt)>=VR(dimcnt,2)
pos(j,dimcnt)=VR(dimcnt,2);
end
if pos(j,dimcnt)<=VR(dimcnt,1) %%%===
pos(j,dimcnt)=VR(dimcnt,1);
end
end
end %结束一次粒子循环
temp=gbest';
fprintf('%f,%f,%f,%f\n',temp(1),temp(2),temp(3),gbestval);
Y(i)=gbestval;
X(i)=i;
% 收敛条件
tmp1=abs(tr(i)-gbestval);
if tmp1>ergrd
cnt2=0;
elseif tmp1<=ergrd %%%===
cnt2=cnt2+1;
if cnt2>=ergrdep
break
end
end
end %迭代结束
fprintf('total interations:%d\n',i);
%OUT=[gbest';gbestval];
pa=gbest;
mi=gbestval;
%POWELL优化算法
function [pa,mi]=powell(R,F,x)
len=5; %搜索区间
itmax=30; %最大循环次数
e=1e-3; %允许误差
%方向矢量矩阵存放d1,d2,d3三个方向矢量
D=[1 0 0;
0 1 0;
0 0 1];
%起始点
%x0=[-10,-12,0];
x0=x;
fx0=BI_mi(x0(1),x0(2),x0(3),R,F);
%fx0=pv_mi(x0(1),x0(2),-x0(3),R,F);
%循环搜索
for k=0:itmax
%从起始点出发，沿d1方向搜索，得到最大值fx1,对应点x1
d1=D(1,:);
[x1,fx1,step]=oneDimSearch(R,F,x0,d1,len);
fprintf('%f, %f, %f, %f\n',x1(1),x1(2),x1(3),fx1);
%从x1出发，沿d2方向搜索
d2=D(2,:);
[x2,fx2,step]=oneDimSearch(R,F,x1,d2,len);
fprintf('%f, %f, %f, %f\n',x2(1),x2(2),x2(3),fx2);
%从x2出发，沿d3方向搜索
d3=D(3,:);
[x3,fx3,step]=oneDimSearch(R,F,x2,d3,len);
%加速方向
dn=x3-x0;
y=sum(dn.*dn);
fprintf('%f\n',y);
fprintf('%f, %f, %f, %f\n',x3(1),x3(2),x3(3),fx3);
if y<=e %满足结束条件，退出循环 %%%===
pa=x3;
mi=fx3;
return;
end
%调整搜索方向
%计算目标值下降最多的那个方向
cha=[fx0-fx1,fx1-fx2,fx2-fx3];
[maxcha,j0]=max(cha);
%
xe=2*x3-x0;
fe=BI_mi(xe(1),xe(2),xe(3),R,F);
%fe=pv_mi(xe(1),xe(2),-xe(3),R,F);
%这里求极大值
if (fe<=fx0) || (2*(fx0-2*fx3+fe)*(fx0-fx3-maxcha)*(fx0-fx3-maxcha) >= (fx0-fe)*(fx0-fe)*maxcha)
%不引入新的方向 %%%===
x0=x3; %下次搜索的新起点
fx0=fx3;
else %引进新的方向
%以xn为起点沿dn进行搜索，得到下次搜索的新起点
[x0,fx0,step]=oneDimSearch(R,F,x3,dn,len);
%替换方向
D(4,:)=dn;
for i=j0:3
D(i,:)=D(i+1,:);
end
end
end %end for
pa=x3;
mi=fx3;
%一维搜索。从-len到len，找出最大值的点
function [Y,fY,step]=oneDimSearch(R,F,X,direction,len)
%一维brent搜索求函数最小值
CGOLD =0.3819660;
a=-len;
b=len;
%a=a0;
%b=b0;
tol=0.01;
e=0.0;
ITMAX=100;
v=a+CGOLD*(b-a);
Xv=X+direction*v;
%fv=pv_mi(Xv(1),Xv(2),-Xv(3),R,F);
fv=BI_mi(Xv(1),Xv(2),Xv(3),R,F);
w=v;
x=v;
fw=fv;
fx=fv;
for k=0:ITMAX
xm=0.5*(a+b); %中点
tol1=tol*abs(x);
tol2=2*tol1;
%if abs(x-xm)=tol;
if abs(x-xm)<=(tol2-0.5*(b-a)) %满足精度，退出 %%%===
step=x;
Y=X+direction*step;
fY=fx;
break;
end
if (abs(e) > tol1)
r=(x-w)*(fx-fv);
q=(x-v)*(fx-fw);
p=(x-v)*q-(x-w)*r; %分子
q=2*(q-r); %分母
if (q>0)
p=-p;
end
q=abs(q);
etemp=e;
e=d;
if (abs(p) >= abs(0.5*q*etemp) || p<=q*(a-x) || p >= q*(b-x)) %使用黄金分割 %%%===
if x >= xm
e=a-x;
else
e=b-x;
end
d=CGOLD*e;
else %抛物线内插
d=p/q;
u=x+d; %u
if(u-a<tol2 || b-u<tol2) %u靠近边界 %%%===
if x > xm;
d=-abs(tol1);
else
d=abs(tol1);
end
end
end
else %使用黄金分割
if x >= xm
e=a-x;
else
e=b-x;
end
d=CGOLD*e;
end %if (abs(e) > tol1)
if abs(d)>=tol1
u=x+d;
else
if d>=0
u=x+abs(tol1);
else
u=x-abs(tol1);
end
end
Xu=X+direction*u;
% fu=pv_mi(Xu(1),Xu(2),-Xu(3),R,F); %函数值
fu=BI_mi(Xu(1),Xu(2),Xu(3),R,F); %函数值
%更新参数
if fu>=fx
if u>=x
a=x;
else
b=x;
end
v=w;w=x;x=u;
fv=fw;fw=fx;fx=fu;
else
if u<x
a=u;
else
b=u;
end
if (fu>=fw || w==x)
v=w;w=u;
fv=fw;fw=fu;
elseif (fu>=fv || v==x || v==w)
v=u;
fv=fu;
end
end
end %end for
- 1.
- 2.
- 3.
- 4.
- 5.
- 6.
- 7.
- 8.
- 9.
- 10.
- 11.
- 12.
- 13.
- 14.
- 15.
- 16.
- 17.
- 18.
- 19.
- 20.
- 21.
- 22.
- 23.
- 24.
- 25.
- 26.
- 27.
- 28.
- 29.
- 30.
- 31.
- 32.
- 33.
- 34.
- 35.
- 36.
- 37.
- 38.
- 39.
- 40.
- 41.
- 42.
- 43.
- 44.
- 45.
- 46.
- 47.
- 48.
- 49.
- 50.
- 51.
- 52.
- 53.
- 54.
- 55.
- 56.
- 57.
- 58.
- 59.
- 60.
- 61.
- 62.
- 63.
- 64.
- 65.
- 66.
- 67.
- 68.
- 69.
- 70.
- 71.
- 72.
- 73.
- 74.
- 75.
- 76.
- 77.
- 78.
- 79.
- 80.
- 81.
- 82.
- 83.
- 84.
- 85.
- 86.
- 87.
- 88.
- 89.
- 90.
- 91.
- 92.
- 93.
- 94.
- 95.
- 96.
- 97.
- 98.
- 99.
- 100.
- 101.
- 102.
- 103.
- 104.
- 105.
- 106.
- 107.
- 108.
- 109.
- 110.
- 111.
- 112.
- 113.
- 114.
- 115.
- 116.
- 117.
- 118.
- 119.
- 120.
- 121.
- 122.
- 123.
- 124.
- 125.
- 126.
- 127.
- 128.
- 129.
- 130.
- 131.
- 132.
- 133.
- 134.
- 135.
- 136.
- 137.
- 138.
- 139.
- 140.
- 141.
- 142.
- 143.
- 144.
- 145.
- 146.
- 147.
- 148.
- 149.
- 150.
- 151.
- 152.
- 153.
- 154.
- 155.
- 156.
- 157.
- 158.
- 159.
- 160.
- 161.
- 162.
- 163.
- 164.
- 165.
- 166.
- 167.
- 168.
- 169.
- 170.
- 171.
- 172.
- 173.
- 174.
- 175.
- 176.
- 177.
- 178.
- 179.
- 180.
- 181.
- 182.
- 183.
- 184.
- 185.
- 186.
- 187.
- 188.
- 189.
- 190.
- 191.
- 192.
- 193.
- 194.
- 195.
- 196.
- 197.
- 198.
- 199.
- 200.
- 201.
- 202.
- 203.
- 204.
- 205.
- 206.
- 207.
- 208.
- 209.
- 210.
- 211.
- 212.
- 213.
- 214.
- 215.
- 216.
- 217.
- 218.
- 219.
- 220.
- 221.
- 222.
- 223.
- 224.
- 225.
- 226.
- 227.
- 228.
- 229.
- 230.
- 231.
- 232.
- 233.
- 234.
- 235.
- 236.
- 237.
- 238.
- 239.
- 240.
- 241.
- 242.
- 243.
- 244.
- 245.
- 246.
- 247.
- 248.
- 249.
- 250.
- 251.
- 252.
- 253.
- 254.
- 255.
- 256.
- 257.
- 258.
- 259.
- 260.
- 261.
- 262.
- 263.
- 264.
- 265.
- 266.
- 267.
- 268.
- 269.
- 270.
- 271.
- 272.
- 273.
- 274.
- 275.
- 276.
- 277.
- 278.
- 279.
- 280.
- 281.
- 282.
- 283.
- 284.
- 285.
- 286.
- 287.
- 288.
- 289.
- 290.
- 291.
- 292.
- 293.
- 294.
- 295.
- 296.
- 297.
- 298.
- 299.
- 300.
- 301.
- 302.
- 303.
- 304.
- 305.
- 306.
- 307.
- 308.
- 309.
- 310.
- 311.
- 312.
- 313.
- 314.
- 315.
- 316.
- 317.
- 318.
- 319.
- 320.
- 321.
- 322.
- 323.
- 324.
- 325.
- 326.
- 327.
- 328.
- 329.
- 330.
- 331.
- 332.
- 333.
- 334.
- 335.
- 336.
- 337.
- 338.
- 339.
- 340.
- 341.
- 342.
- 343.
- 344.
- 345.
- 346.
- 347.
- 348.
- 349.
- 350.
- 351.
- 352.
- 353.
- 354.
- 355.
- 356.
- 357.
- 358.
- 359.
- 360.
- 361.
- 362.
- 363.
- 364.
- 365.
- 366.
- 367.
- 368.
- 369.
- 370.
- 371.
- 372.
- 373.
- 374.
- 375.
- 376.
- 377.
- 378.
- 379.
- 380.
- 381.
- 382.
- 383.
- 384.
- 385.
- 386.
- 387.
- 388.
- 389.
- 390.
- 391.
- 392.
- 393.
- 394.
- 395.
- 396.
- 397.
- 398.
- 399.
- 400.
- 401.
- 402.
- 403.
- 404.
- 405.
- 406.
- 407.
- 408.
- 409.
- 410.
- 411.
- 412.
- 413.
- 414.
- 415.
- 416.
- 417.
- 418.
- 419.
- 420.
- 421.
- 422.
- 423.
- 424.
- 425.
- 426.
- 427.
- 428.
- 429.
- 430.
- 431.
- 432.
- 433.
- 434.
- 435.
- 436.
- 437.
- 438.
- 439.
- 440.
- 441.
- 442.
- 443.
- 444.
- 445.
- 446.
- 447.
- 448.
- 449.
- 450.
- 451.
- 452.
- 453.
- 454.
- 455.
- 456.
- 457.
- 458.
- 459.
- 460.
- 461.
- 462.
- 463.
- 464.
- 465.
- 466.
- 467.
- 468.
- 469.
- 470.
- 471.
- 472.
- 473.
- 474.
3 仿真结果
4 参考文献
[1]汪惠兰, 林航飞, 李赔龙. 基于Matlab实时运动目标跟踪检测系统[J]. 2022(3).
博主简介：擅长智能优化算法、神经网络预测、信号处理、元胞自动机、图像处理、路径规划、无人机等多种领域的Matlab仿真，相关matlab代码问题可私信交流。
部分理论引用网络文献，若有侵权联系博主删除。
1475
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
