art.m

来自「几个图象处理程序」· M 代码 · 共 40 行

     sl=phantom;
     theta = 0:179;
     [R,xp] = radon(sl,theta);
      figure(1);
      imshow(sl);      
      title('原图');
      [width,length]=size(R); 
       width1=256;
       length1=256;    
rec=zeros(width1,length1);
for l=1:1%叠代次数
   for k=1:360%从180个方向扫描
    num=zeros(width,1);%初始一列量，用来存每条线扫过的格数
      y=zeros(width,1); %初始一列量，用来存每条线扫过的格数的值    
           for i=1:width1               
                 for j=1:length1            
                         n= round((i-128)*cos((k-1)*pi/180)+(-j+128)*sin((k-1)*pi/180)+184);               
                         num(n)=num(n)+1;
                         y(n)=y(n)+rec(i,j);
                 end
           end
    y=R(:,k)-y;%用已知的投影值减去计算的投影值之和
    for n=1:width
      if num(n)~=0
    y(n)=y(n)/num(n);
     else
          y(n)=y(n);
      end
   end
        for i=1:width1
              for j=1:length1            
                  n=round((i-128)*cos((k-1)*pi/180)+(-j+128)*sin((k-1)*pi/180))+184;     
                  rec(i,j)= rec(i,j)+y(n);       
              end
       end
 end
end
figure(2);
imshow(rec');
title('art');