cll_rgst.m

图像对齐（inage registration）工具箱

             SG = vrat('STDMenu');
             tp = 1; 
             
             [POINTS1,RESULT]=detectorHR(M1, SG,nCP,dCP,tp);

             zprava('Sensed image CP computing','b'); 
             [POINTS2,RESULT]=detectorHR(M2, SG,nCP,dCP,tp);

             
          case (3),
          % manual %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

             nCP = str2num(get(findobj(gcbf,'Tag','NumCP'),'String'));
                          
             POINTS1=zeros(size(M1));
             POINTS2=zeros(size(M2));
               
             set(gcbf,'CurrentAxes',h1);
             
             zbarvi('Ref_im','r');

             ramek('Axes1',[1 0 0]);
             axis on;
             hold on;
             
             ptr=1;
             gptr=1;CN = 0;
             clear hn;
             
             for i=1:nCP,
                btt=1;x=0;y=0;
                while btt~=3
                   if ptr>2+CN
                      delete(hn(1+CN));
                      hn(1+CN)=hn(2+CN);
                      ptr=2+CN;
                      refresh;
                   end;
                   
                   xx=x;yy=y;
                   [x,y,btt]=ginput(1);
                   x=floor(x);
                   y=floor(y);
                   x=min(max(x,1),size(M1,2));
                   y=min(max(y,1),size(M1,1));
                   hn(ptr)=text(x,y,'x','Color',[1 0 0]);
                   ptr=ptr+1;
                end;
                
                delete(hn(ptr-1));
                ptr=ptr-1;
                refresh;
                CN = CN +1;
                
                ghn(gptr)=text(xx,yy,'x','Color',[0 0 1]);
                ghn2(gptr)=text(xx+3,yy+3,int2str(gptr),'Color',[0 0 1]);
                POINTS1(yy,xx)=1;
                zprava(['Reference image CP: ' num2str(i)],'b'); 
                gptr=gptr+1;
             end;
             
             ramek('Axes1',[0 0 0]);
             axis off;
             zbarvi('Ref_im','k');
             zbarvi('Sens_im','r');

             zprava('Sensed image CP computing','b'); 
             
             set(gcbf,'CurrentAxes',h2);
             ramek('Axes2',[1 0 0]);
             axis on;
             
             hold on;
             bptr=1;
             bgptr=1;
             CN=0;
             clear bhn;
             
             for i=1:nCP,
                btt=1;x=0;y=0;
                while btt~=3
                   
                   if bptr>2+CN
                      delete(bhn(1+CN));
                      bhn(1+CN)=bhn(2+CN);
                      bptr=2+CN;
                      refresh;
                   end;
                   xx=x;yy=y;
                   [x,y,btt]=ginput(1);
                   x=floor(x);
                   y=floor(y);
                   x=min(max(x,1),size(M2,2));
                   y=min(max(y,1),size(M2,1));
                   bhn(bptr)=text(x,y,'x','Color',[1 0 0]);
                   bptr=bptr+1;
                end;
                
                delete(bhn(bptr-1));
                bptr=bptr-1;
                refresh;
                CN = CN+1;                
                bghn(bgptr)=text(xx,yy,'x','Color',[0 0 1]);
                bghn2(bgptr)=text(xx+3,yy+3,int2str(bgptr),'Color',[0 0 1]);
                POINTS2(yy,xx)=1;
                zprava(['Sensed image CP: ' num2str(i)],'b'); 
                bgptr=bgptr+1;
             end;
             

 
             while gptr>1
                delete(ghn(gptr-1));
                delete(ghn2(gptr-1));
                gptr=gptr-1;
             end;
                               
             while bgptr>1
                delete(bghn(bgptr-1));
                delete(bghn2(bgptr-1));
                bgptr=bgptr-1;
             end;
             
             while ptr>1
                delete(hn(ptr-1));
                ptr=ptr-1;
             end;
 
             while bptr>1
                delete(bhn(bptr-1));
                bptr=bptr-1;
             end;

             refresh;
          zbarvi('Sens_im','k');

          set(gcbf,'pointer','arrow');
          ramek('Axes2',[0 0 0]);
          axis off;
                 
     end;                     
          
       
          
          data.F1 = POINTS1;
          data.F2 = POINTS2;
          
         [x1,y1]=find(POINTS1>0);
   		 [x2,y2]=find(POINTS2>0);
   
          data.x1 = x1;
          data.x2 = x2;
          
          data.y1 = y1;
          data.y2 = y2;
          set(gcbf,'Userdata',data);  
  	      cll_rgst wind;          

       end;	
    end;
    
    
    % Window selection %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
case('wsel')
    zprava('','b');
   
   	zbarvi('Sens_im','r');

	   data 	= get(gcbf,'Userdata');
   

   	M2 = data.M2;
    szM2=size(M2);
    rect=[1 1 szM2(2) szM2(1)];
    
      val = get(findobj(gcbf,'Tag','wnsz'),'String');
 if isempty(val)
         zprava('SENSED IMAGE: Window selection, hold left mouse and plot square inside.','r'); 

         set(gcbf,'CurrentAxes',data.h2);
     while ((rect(1)+rect(3))>szM2(2))|((rect(2)+rect(4))>szM2(1))|(rect(1)<1)|(rect(2)<1)
         [vrz,rect]=imcrop;
     end;
     
         dx=floor(rect(2));
         dy=floor(rect(1));

 else
         zprava('SENSED IMAGE: Window selection, left - choice, right - OK','r'); 

         wsz=str2num(val);
    
			dx=1;
      	dy=1;
      	button=1;
      	clear hn;
			flg3=0;
      
      	set(gcbf,'CurrentAxes',data.h2);
	      ramek('Axes2',[1 0 0]);
   	          
      	hold on;
	      tmp=M2;
      
   	   vel=size(M2);
      
			
               vy=-2*vel(2); 
               vx=-2*vel(1); 
               while ((vy+wsz)>vel(2))|((vx+wsz)>vel(1))|(vy<1)|(vx<1)|(button == 1 )
  			     [vy2,vx2,button]=ginput(1);
			     vy=floor(vy2);
			     vx=floor(vx2);
               
    
			   if button == 1
			      pomidx1=max(1,vy):(min(vy+wsz-1,vel(2)));
                  pomidx2=max(vx,1):(min(vx+wsz-1,vel(1)));  
                  tmp(max(vx,1),pomidx1)=zeros(size(pomidx1));
                  tmp(min(vx+wsz-1,vel(1)),pomidx1)=zeros(size(pomidx1));
                  tmp(pomidx2,max(vy,1))=zeros(size(pomidx2))';
                  tmp(pomidx2,min(vy+wsz-1,vel(2)))=zeros(size(pomidx2))';



                  hold on;
			      if flg3 == 1
				        delete(hn);
			      end;
			      hn=imagesc(tmp);
			      flg3=1;
			      tmp=M2;
              end;
              
              
			   end;
               
               vrz=M2(vx:(vx+wsz-1),vy:(vy+wsz-1));
               rect=[vy,vx,wsz,wsz];
			   delete(hn);
		      
               dx=vx;
			   dy=vy;
			end;
       
 
      
  
      h=figure;
      imagesc(vrz);
      title('Space bar close the image.','Color','b','FontSize',16)
      
%      title(['Rectangle: ' int2str(floor(rect(1))) ', ' int2str(floor(rect(2))) ', 'int2str(floor(rect(3))) ', 'int2str(floor(rect(4))) ]); 

		zprava('Space bar close the image','r');

		axis image;
		axis off;
		colormap gray;
		pause;
		close(h);
		zprava(' ','r');

		data.wn = vrz;
		data.wnr = dx;
		data.wnc = dy;
       	zbarvi('Sens_im','k');
		set(gcbf,'Userdata',data);  

      
      
    % Normalized cross correlation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
 case('ncc')	
    	zprava('Normalized cross-correlation','b');
       
      data 	= get(gcbf,'Userdata');
   
    	M1 = data.M1;
        M2 = data.M2;
        wn = data.wn;
		wr = data.wnr;
		wcc = data.wnc;
      
        [mI,nI]=size(M1);
		[mW,nW]=size(wn);
      
        wn=wn-mean(wn(:));

		mm=0:mW-1;
		nn=0:nW-1;
		nrmW=sqrt(sum(sum(wn.*wn)));
		CS=zeros(mI-mW+1,nI-nW+1);
      
      for i=1:(mI-mW+1)
		   for j=1:(nI-nW+1)
            WI=M1(i+mm,j+nn);
            WI=WI-mean(WI(:));
		       CS(i,j)=sum(sum(WI.*wn))/(sqrt(sum(sum(WI.*WI)))*nrmW);
		   end;
		end;


		[a1,a2]=max(CS);
		[b1,b2]=max(a1);
        c=a2(b2);
       
        h=figure;
        mesh(CS);
        title('Space bar to continue.','Color','b');
      
        zprava(['Shift parameters: vertical ', int2str(c-wr),', horizontal ',int2str(b2-wcc),' .Space bar to continue.'],'b');
        pause;
      
        close(h);

		M3=transformuj(M1,M2,[1 0 c-wr b2-wcc ]);	
        zprava(['Space bar to continue. The recently opened window must be active.'],'b');
        pause;
        close;
        zprava('','r');
		data.param = [ c-wr b2-wcc 0 1];
		data.M3=M3;
      
      set(gcbf,'Userdata',data);  

          % Phase cross correlation %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
 case('phc')	
    	zprava('Phase cross-correlation','b');
  
      data 	= get(gcbf,'Userdata');
   
   	M1 = data.M1;
      M2 = data.M2;
      wn = data.wn;
		wr = data.wnr;
		wcc = data.wnc;
      
      [km,dx,dy]=phasekorel(M1,wn);
      
	h=figure;
	mesh(km);
	title('Phase correlation.    Space bar to continue.');

     
      zprava(['Shift parameters: vertical ', int2str(dx-wr),', horizontal ',int2str(dy-wcc)],'b');
      pause;
      
      close(h);

		M3=transformuj(M1,M2,[1 0 dx-wr dy-wcc ]);	
        
        zprava('','r');
		data.param = [ dx-wr dy-wcc  0 1];
		data.M3=M3;
      
      set(gcbf,'Userdata',data);  

      
      
   
	   
   
   % Registration %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   
case('registration')
   
   	data 	= get(gcbf,'Userdata');
      zprava('Registration','b');
      
      
      if data.flg == 'AUTO'
         cll_rgst CPdetection;
         cll_rgst inv_comp;
         cll_rgst twobest;
      end;
      
      mthr = vrat('ThrMenu');             
      rdm = str2num(get(findobj(gcbf,'Tag','RadiusM'),'String'));
      rddm = str2num(get(findobj(gcbf,'Tag','RRadM'),'String'));
      ord = str2num(get(findobj(gcbf,'Tag','Ord'),'String'));	

      invt= data.intp;      
      
      M1 = data.M1;
      M2 = data.M2;
      
   	x1		= data.x1;
   	x2		= data.x2;
   
   	y1		= data.y1;
   	y2		= data.y2;
      
      best = data.best;

      par1=ps(best([1 2],:));  
          
      bbx=best(:,3)*par1(1)-best(:,4)*par1(2)+par1(3);
      bby=best(:,4)*par1(1)+best(:,3)*par1(2)+par1(4);
      
      points(best(:,[1 2]),[bbx bby],'Control points: Ref | Sensed_transform');
      uiwait;
      refresh;
      tx2=x2*par1(1)-y2*par1(2)+par1(3);
      ty2=y2*par1(1)+x2*par1(2)+par1(4);
      
		pary=zparuj([x1, y1],[x2, y2],[tx2, ty2],M1,M2,rdm,rddm,mthr,invt,ord);      
      
      %points(pary(:,[1 2]),pary(:,[5 6]),'Control points: R S_transf');
      %uiwait;
      %points(pary(:,[1 2]),pary(:,[7 8]),'Control points: R S_refin');
      %uiwait;
      %points(pary(:,[3 4]),pary(:,[7 8]),'Control points: S S_refin');
      %uiwait;

      [par2,s,uhel]=mnc(pary(:,[7 8 1 2]));
      
      s=sqrt(par2(1)^2+par2(2)^2);
		u1=acos(par2(1)/s);
		u2=asin(par2(2)/s);
		znameni=sign(u2)<0;
		uhel=(~znameni)*u1+(znameni)*(2*pi-u1);
      uhel=uhel/(pi/180);
      
  %    zprava(['Registration parameters: Tx: 'int2str(par2(3)) '  Ty: 'int2str(par2(4)) ' alfa: 'int2str(uhel) ' scale: 'int2str(s)],'b');
      % pary contain: [sour bodu z ref obr, 1.radek, sour bodu z sens obr, 1.radek, sour bodu z transf sens obr, 1.radek, sour bodu z sen upravene obr, 1.radek, ]
      
      M3=transformuj(M1,M2,par2);
      
      data.param = [par2(3) par2(4) uhel s ];
  		data.M3 = M3;
       
      set(gcbf,'Userdata',data);  
     
      
      % Computation of invariants at the positions of points %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   case('inv_comp')

      data 	= get(gcbf,'Userdata');
      
      M1		= data.M1; %img
      M2		= data.M2;
      
      F1    = data.F1;  %points
  	  F2	= data.F2;
     
      flgW  = data.flgW; 
      
      zprava('Invariant computation','b');

      rd = str2num(get(findobj(gcbf,'Tag','Radius'),'String'));
      rdd = str2num(get(findobj(gcbf,'Tag','RRad'),'String'));
      tps = get(findobj(gcbf,'Tag','INVTypMenu'),'String');
      tpv = get(findobj(gcbf,'Tag','INVTypMenu'),'Value');
      tp = tps(tpv,:);
      ord = str2num(get(findobj(gcbf,'Tag','Ord'),'String'));	
      

      
      zprava([tp ' invariant computation with parameters:' num2str(rd) ' '   num2str(rdd)] ,'b');
      
      if flgW == 'OFF'
      
      [pom1,INV1,pts1 ] = match(M1, F1, rd,0,tp,ord);
      [pom2,INV2,pts2 ] = match(M2, F2, rd,rdd,tp,ord);
      
      krouzkuj(pts1,data.h1);
      krouzkuj(pts2,data.h2);

      else
      INV1=[];
      INV2=[];
      end;
      
      
      data.INV1 = INV1;
      data.INV2 = INV2;
      
      data.INVS1 = INV1;
      data.INVS2 = INV2;

      data.rd  = rd;
      data.rdd = rdd;
      data.intp = tp;
      data.pts1 = pts1;
      data.pts2 = pts2;

	  data.choice = 1:size(INV1,2);

      set(gcbf,'Userdata',data);  
      zprava('Invariant computation finished!' ,'b');