filtimag1231.m

一个用matlab编写的 关于笔画提取的仿真试验程序

function [a]=filtimag1231(A)

[m,n]=size(A);
dd=5;
j1=round(m/dd);
a=A;
x=zeros(m,n);
 
  x6=zeros(1,6);
  
tic;
for j=3:1:(j1+2)
        for k=3:1:(n-2)
            
            t=0;an=zeros(1,5);
         
         for i=0:3
         if (0<A(j+i*j1,k)&&A(j+i*j1,k)<255)
            a(j+i*j1,k)=A(j+i*j1,k);
         else
            t=t+1;
            b=j+i*j1;
            an(1,t)=b;
         end
         end
           if (0<A(m-4-j1+j,k)&&A(m-4-j1+j,k)<255)
               a(m-4-j1+j,k)=A(m-4-j1+j,k);
          else
            t=t+1;
            b=m-4-j1+j;
            an(1,t)=b;
           end
               
        %%%
    
           if t==5
            x1=zeros(3,3);b1=an(1,1);
            x2=zeros(3,3);b2=an(1,2);
            x3=zeros(3,3);b3=an(1,3);
            x4=zeros(3,3);b4=an(1,4);
            x5=zeros(3,3);b5=an(1,5);
             for p=1:1:3
             for q=1:1:3
               
                x1(p,q)=A(b1-2+p,k-2+q);%2=(1+mod)/2
                x2(p,q)=A(b2-2+p,k-2+q);%2=(1+mod)/2
                x3(p,q)=A(b3-2+p,k-2+q);%2=(1+mod)/2       
                x4(p,q)=A(b4-2+p,k-2+q);%2=(1+mod)/2                           
                x5(p,q)=A(b5-2+p,k-2+q);%2=(1+mod)/2                           
             end
             end
        
           elseif t==4
              x1=zeros(3,3);b1=an(1,1);
            x2=zeros(3,3);b2=an(1,2);
            x3=zeros(3,3);b3=an(1,3);
            x4=zeros(3,3);b4=an(1,4);
               for p=1:1:3
             for q=1:1:3
               
                x1(p,q)=A(b1-2+p,k-2+q);%2=(1+mod)/2
                x2(p,q)=A(b2-2+p,k-2+q);%2=(1+mod)/2
                x3(p,q)=A(b3-2+p,k-2+q);%2=(1+mod)/2       
                x4(p,q)=A(b4-2+p,k-2+q);%2=(1+mod)/2                           
                               
             end
             end
              
             elseif t==3
              x1=zeros(3,3);b1=an(1,1);
            x2=zeros(3,3);b2=an(1,2);
            x3=zeros(3,3);b3=an(1,3);
               for p=1:1:3
             for q=1:1:3
               
                x1(p,q)=A(b1-2+p,k-2+q);%2=(1+mod)/2
                x2(p,q)=A(b2-2+p,k-2+q);%2=(1+mod)/2
                x3(p,q)=A(b3-2+p,k-2+q);%2=(1+mod)/2       
                                 
             end
             end
         elseif t==2
              x1=zeros(3,3);b1=an(1,1);
            x2=zeros(3,3);b2=an(1,2);
           
               for p=1:1:3
             for q=1:1:3
               
                x1(p,q)=A(b1-2+p,k-2+q);%2=(1+mod)/2
                x2(p,q)=A(b2-2+p,k-2+q);%2=(1+mod)/2
                              
             end
             end
                           
         elseif t==1 
             
          x1=zeros(3,3);b1=an(1,1);
               for p=1:1:3
             for q=1:1:3
               
                x1(p,q)=A(b1-2+p,k-2+q);%2=(1+mod)/2
                                
             end
             end
        
         end
         
          %%以下处理
   
         if t>=1
             x111=x1;
               mod=3;
           [av3]=elsm3gai1231(x111,mod);
              if av3==-1,
                 mod=5;
                x1=zeros(mod,mod);
                %5＊5大小的模板
                 for p=1:1:mod
                   for q=1:1:mod
                      x1(p,q)=A(b1-3+p,k-3+q);%3=(1+mod)/2!!!!!!
                  end
                end
                x111=x1;
                 [av5]=elsm3gai1231(x111,mod);
              if av5==-1,
                   l1=9;
                  sum=0;
             for p=1:1:3
             for q=1:1:3
                if (A(b1-2+p,k-2+q)==0||A(b1-2+p,k-2+q)==255)
                      l1=l1-1;
                  else
                  sum=sum+double(A(b1-2+p,k-2+q));
                end
             end
             end
               if l1==0,
                   a(b1,k)=50; 
               else 
                   a(b1,k)=round(sum/l1);  
               end
            else
                  a(b1,k)=av5;
              end 
            else 
            a(b1,k)=av3;
            end
            t=t-1;
        end
        
      
        if t>=1
             x111=x2;
               mod=3;
           [av3]=elsm3gai1231(x111,mod);
              if av3==-1,
                 mod=5;
                x1=zeros(mod,mod);
                %5＊5大小的模板
                 for p=1:1:mod
                   for q=1:1:mod
                      x1(p,q)=A(b2-3+p,k-3+q);%3=(1+mod)/2!!!!!!
                  end
                end
                x111=x1;
                 [av5]=elsm3gai1231(x111,mod);
              if av5==-1,
                   l1=9;
                  sum=0;
             for p=1:1:3
             for q=1:1:3
                if (A(b2-2+p,k-2+q)==0||A(b2-2+p,k-2+q)==255)
                      l1=l1-1;
                  else
                  sum=sum+double(A(b2-2+p,k-2+q));
                end
             end
             end
            if l1==0,
                   a(b2,k)=50; 
               else 
                   a(b2,k)=round(sum/l1);  
               end
            else
                  a(b2,k)=av5;
              end 
            else 
            a(b2,k)=av3;
            end
            t=t-1;
       
        end
        
        if t>=1
              mod=3;
             x111=x3;
           [av3]=elsm3gai1231(x111,mod);
              if av3==-1,
                 mod=5;
                x1=zeros(mod,mod);
                %5＊5大小的模板
                 for p=1:1:mod
                   for q=1:1:mod
                      x1(p,q)=A(b3-3+p,k-3+q);%3=(1+mod)/2!!!!!!
                  end
                end
                x111=x1;
                 [av5]=elsm3gai1231(x111,mod);
              if av5==-1,
                   l1=9;
                  sum=0;
             for p=1:1:3
             for q=1:1:3
                if (A(b3-2+p,k-2+q)==0||A(b3-2+p,k-2+q)==255)
                      l1=l1-1;
                  else
                  sum=sum+double(A(b3-2+p,k-2+q));
                end
             end
             end
            if l1==0,
                   a(b3,k)=50; 
               else 
                   a(b3,k)=round(sum/l1);  
               end
            else
                  a(b3,k)=av5;
              end 
            else 
            a(b3,k)=av3;
            end
            t=t-1;
        
        end
        
        if t>=1
             x111=x4;
               mod=3;
           [av3]=elsm3gai1231(x111,mod);
              if av3==-1,
                 mod=5;
                x1=zeros(mod,mod);
                %5＊5大小的模板
                 for p=1:1:mod
                   for q=1:1:mod
                      x1(p,q)=A(b4-3+p,k-3+q);%3=(1+mod)/2!!!!!!
                  end
                end
                x111=x1;
                 [av5]=elsm3gai1231(x111,mod);
              if av5==-1,
                   l1=9;
                  sum=0;
             for p=1:1:3
             for q=1:1:3
                if (A(b4-2+p,k-2+q)==0||A(b4-2+p,k-2+q)==255)
                      l1=l1-1;
                  else
                  sum=sum+double(A(b4-2+p,k-2+q));
                end
             end
             end
            if l1==0,
                   a(b4,k)=50; 
               else 
                   a(b4,k)=round(sum/l1);  
               end
            else
                  a(b4,k)=av5;
              end 
            else 
            a(b4,k)=av3;
            end
           
            t=t-1;
       
        end
           if t>=1
             x111=x5;
               mod=3;
           [av3]=elsm3gai1231(x111,mod);
              if av3==-1,
                 mod=5;
                x1=zeros(mod,mod);
                %5＊5大小的模板
                 for p=1:1:mod
                   for q=1:1:mod
                      x1(p,q)=A(b5-3+p,k-3+q);%3=(1+mod)/2!!!!!!
                  end
                end
                x111=x1;
                 [av5]=elsm3gai1231(x111,mod);
              if av5==-1,
                   l1=9;
                  sum=0;
             for p=1:1:3
             for q=1:1:3
                if (A(b5-2+p,k-2+q)==0||A(b5-2+p,k-2+q)==255)
                      l1=l1-1;
                  else
                  sum=sum+double(A(b5-2+p,k-2+q));
                end
             end
             end
            if l1==0,
                   a(b5,k)=50; 
               else 
                   a(b5,k)=round(sum/l1);  
               end
            else
                  a(b5,k)=av5;
              end 
            else 
            a(b5,k)=av3;
            end
           
    
        end
           end
      
       end

 
%%%%以下为边缘消噪！
    
   for j1=2:(-1):1
     for k1=2:(n-1)
           if (0<A(j1,k1)&&A(j1,k1)<255)
            a(j1,k1)=A(j1,k1);
        elseif (A(j1,k1)<=0||A(j1,k1)>=255)
          r=0;
        for p=j1:j1+1
            for q=k1-1:k1+1
              r=r+1;
              x6(1,r)=a(p,q);
            end
       end 
       [t13,e113]=edge3(x6);
     if e113==-1,
           sum=0;
           l2=0;
         for p=j1:j1+2
            for q=k1-1:k1+1
                if a(p,q)==0||a(p,q)==255
                    break;
                else
             sum=sum+double(a(p,q));l2=l2+1;
                end
            end
           end
           if l2==0,
               a(j2,k2)=50;
           else
             a(j1,k1)=round(sum/l2);
         end
      else
          a(j1,k1)=e113;
      end
  end
  end      
end

 
for j2=m-1:m
     for k2=2:(n-1)
           if (0<A(j2,k2)&&A(j2,k2)<255)
            a(j2,k2)=A(j2,k2);
        elseif (A(j2,k2)<=0||A(j2,k2)>=255)
          r=0;
        for p=j2-1:j2
            for q=k2-1:k2+1
                r=r+1;
              x6(1,r)=a(p,q);
            end
       end 
       [t13,e113]=edge3(x6);
     if e113==-1
           sum=0;
           l2=0;
         for p=j2:-1:j2-2
            for q=k2-1:k2+1
                if a(p,q)==0||a(p,q)==255
                   break; 
                else
             sum=sum+double(a(p,q));l2=l2+1;
                end
            end
           end
           if l2==0,
               a(j2,k2)=50;
           else
             a(j2,k2)=round(sum/l2);
         end
          %a(j,k)=(double(a(j,k-1))+double(a(j,k+1))+double(a(j+1,k+1)))/3;
      else
          a(j2,k2)=e113;
      end
  end
  end      
end
  


for j3=2:m-1
    for k3=2:-1:1
      if (0<A(j3,k3)&&A(j3,k3)<255)
            a(j3,k3)=A(j3,k3);
        elseif (A(j3,k3)<=0||A(j3,k3)>=255)
          r=0;
        for p=j3-1:j3+1
            for q=k3:k3+1
              r=r+1;
              x6(1,r)=a(p,q);
            end
       end 
       [t13,e113]=edge3(x6);
     if e113==-1
             sum=0;
           l2=0;
         for p=j3-1:j3+1
            for q=k3:k3+2
                if (a(p,q)==0||a(p,q)==255)
                    break;
                else
             sum=sum+double(a(p,q));l2=l2+1;
                end
            end
           end
           if l2==0,
               a(j3,k3)=50;
           else
             a(j3,k3)=round(sum/l2);
         end
          %a(j,k)=(double(a(j,k-1))+double(a(j,k+1))+double(a(j+1,k+1)))/3;
      else
        
          a(j3,k3)=e113;
      end
  end
  end      
end
   

for j=2:m-1
    for k=n-1:n
      if (0<A(j,k)&&A(j,k)<255)
            a(j,k)=A(j,k);
        elseif (A(j,k)<=0||A(j,k)>=255)
          r=0;
        for p=j-1:j+1
            for q=k-1:k
       
              r=r+1;
              x6(1,r)=a(p,q);
            end
       end 
       [t13,e113]=edge3(x6);
     if e113==-1
            sum=0;
           l2=0;
         for p=j:-1:j-1
            for q=k-2:k
                if (a(p,q)==0||a(p,q)==255)
                    break;
                else
             sum=sum+double(a(p,q));l2=l2+1;
                end
            end
           end
           if l2==0,
               a(j,k)=50;
           else
             a(j,k)=round(sum/l2);
         end
          %a(j,k)=(double(a(j,k-1))+double(a(j,k+1))+double(a(j+1,k+1)))/3;
      else
          a(j,k)=e113;
      end
  end
  end      
end

 toc;   
 %bw1=edge(a,'canny');
 %   imwrite(a,'as2.jpg');   
 %   imwrite(bw1,'aaaaabw1122.jpg');   
   
     %   figure(3),imshow(a);