电子稳像-灰度投影算法源码

%投影法,投影修正,余弦滤波
%水平方向最大偏移正负20个象素
%垂直方向最大偏移正负30个象素
clear all
mov=aviread('v2.avi');  


 g=[];

for t=50:80
    mov(t).cdata=rgb2gray(mov(t).cdata);
 
    g=cat(3,g,mov(t).cdata);
end 

 gr=[];
 gr=cat(3,gr,g(:,:,1));

 frow=30;
 fcol=20;
 
referenceframe=g(:,:,1);
refprojrow=zeros(240,1);
refprojcol=zeros(1,320);

%行投影
refprojrowsum=0;
for i=1:240
    for j=1:320
        refprojrow(i)=refprojrow(i)+double(referenceframe(i,j));
      
    end 
    refprojrowsum=refprojrowsum+refprojrow(i);
end 
refprojrowmean=refprojrowsum/240;
refprojrow=refprojrow-refprojrowmean;

figure;plot(refprojrow);


for i=1:240
    if (i240-frow) 
       refprojrow(i)=refprojrow(i)*(1+cos(pi*(frow-1-i)/frow))/2;
   end
end
    
hold on
plot(refprojrow,'r');
    
    
%列投影
refprojcolsum=0;
for j=1:320
   for i=1:240
        refprojcol(j)=refprojcol(j)+double(referenceframe(i,j));
      
    end 
    refprojcolsum=refprojcolsum+refprojcol(j);
end 
refprojcolmean=refprojcolsum/320
refprojcol=refprojcol-refprojcolmean;
figure;plot(refprojcol);

for j=1:320
    if (j240-fcol) 
       refprojcol(j)=refprojcol(j)*(1+cos(pi*(fcol-1-j)/fcol))/2;
   end
end
hold on ; plot(refprojcol,'r');

for t=2:20
    t
    
currentframe=g(:,:,t);

%行投影
curprojrow=zeros(240,1);
curprojrowsum=0;

for i=1:240
    for j=1:320
        %refprojrow(i)=refprojrow(i)+double(referenceframe(i,j));
        curprojrow(i)=curprojrow(i)+double( currentframe(i,j));
    end 
    curprojrowsum=curprojrowsum+curprojrow(i);
end 
curprojrowmean=curprojrowsum/240;
curprojrow=curprojrow-curprojrowmean;


for i=1:240
    if (i240-frow) 
       curprojrow(i)=curprojrow(i)*(1+cos(pi*(frow-1-i)/frow))/2;
   end
end

%figure;plot(refprojrow);hold on 
%plot(curprojrow,'r');xlabel('行投影')
%180单位的相关运算   %垂直方向最大偏移正负30个象素
cr=zeros(1,61);
for j=1:61
    for i=1:180
       cr(j)=cr(j)+(refprojrow(j+i-1)-curprojrow(30+i))^2;
    end
end

%figure;plot(cr,'g'); xlabel('垂直方向')

[b,jmin]=min(cr);
dy=31-jmin

%列投影
curprojcol=zeros(1,320);
curprojcolsum=0;

for j=1:320
   for i=1:240
        %refprojcol(j)=refprojcol(j)+double(referenceframe(i,j));
        curprojcol(j)=curprojcol(j)+double(currentframe(i,j));
    end 
    curprojcolsum=curprojcolsum+curprojcol(j);
end 
curprojcolmean=curprojcolsum/320;
curprojcol=curprojcol-curprojcolmean;


for j=1:320
    if (j240-fcol) 
       curprojcol(j)=curprojcol(j)*(1+cos(pi*(fcol-1-j)/fcol))/2;
   end
end

%figure;plot(refprojcol);hold on;
%plot(curprojcol,'r');   xlabel('列投影')

%水平方向最大偏移正负20个象素
%280单位的相关运算
cc=zeros(1,41);
for i=1:41
    for j=1:280
      cc(i)=cc(i)+(refprojcol(j+i-1)-curprojcol(20+j))^2;
    end
end

%figure;plot(cc,'g');   xlabel('水平方向')

[a,imin]=min(cc);
dx=21-imin

if dy<0
    %dy<0，图currentframe相对于referenceframe向上运动了|dy|
    a=zeros(abs(dy),320);
   guoduframe=[a;currentframe(1:240-abs(dy),:)];
elseif dy>0
    %dy>0，图currentframe相对于referenceframe向下运动了|dy|,
    a=zeros(dy,320);
    guoduframe=[currentframe(1+abs(dy):240,:);a];
else
     guoduframe=currentframe;
 end


    
if dx<0
   %dx<0，图currentframe相对于referenceframe向左运动了|dx| 
    c=zeros(240,abs(dx));
    buchangframe=[c,guoduframe(:,1:320-abs(dx))]; 
elseif dx>0
    %dx>0，图currentframe相对于referenceframe向右运动了|dx|
    c=zeros(240,dx);
    buchangframe=[guoduframe(:,abs(dx)+1:320),c];
else
    buchangframe=guoduframe;
 end
 
 %figure;imshow(buchangframe);
gr=cat(3,gr,buchangframe);

  %referenceframe=buchangframe;
  %refprojrow=curprojrow;
  %refprojcol=curprojcol;
end