%%%%%%5相位編組法實現紋理直線邊緣的檢測程序代碼(代碼以子函數形式給出)
function gradephase(x)
[m,n]=size(x);
bw=edge(x,'sobel');
gy=x(1:m-1,1:n-1)-x(2:m,1:n-1);
gx= x(1:m-1,1:n-1)-x(1:m-1,2:n);
g=gy./(gx+eps);
ph-atan(g)+(sign(gx)<0&sign(gy)>0)*pi+...
    (sign(gx)<0&sign(gy)>0)+ (sign(gx)>0&sign(gy)<0)*2*pi+(sign(gx)==0)*pi;
grdgp=floor(ph/pi*4);
cn=0;s=[ ];
In_spt=cell(1,1);
% cell數組類似c語言的指針，可以動態改變大小。S用來存放直線編組區域內的像％素坐標。
pline=[];
for  i=2:m-1
    for  j=2:n-1
        if  bw(i,j)==D&gradgp(i,j)==0 %搜索邊緣點。
            ph_cp=gradgp(i,j);gradgp(i,j)=0;bw(i,j)=0;
            cn=cn+1;p=[i,j];
            s=[s,p];in_spt(cn)-[p];
            while ~isempty(s) %為了在8個鄰域進行搜索，采用入棧和出棧操作。
                [cs,rs]=size(s);
                ps=s(:,rs);s=s(:,1:rs-1);
                col=ps(1,1);row=ps(2,1);
                if  legal(col+1,row,m,n)&grdgp(col+1,row)=ph_cp
                    s= [s,[col+1,row]’]; in_spt(cn)=[ in_spt(cn), [col+1,row]’];
                    bw(col+1,row)=0;gradgp(col+1,row)=0;
                end
                if  legal(col+1,row+1,m,n)&grdgp(col+1,row+1)=ph_cp
                    s= [s,[col+1,row+1]’]; in_spt(cn)=[ in_spt(cn), [col+1,row+1]’];
                    bw(col+1,row+1)=0;gradgp(col+1,row+1)=0;
                end
                if  legal(col,row+1,m,n)&grdgp(col,row+1)=ph_cp
                    s= [s,[col,row+1]’]; in_spt(cn)=[ in_spt(cn), [col,row+1]’];
                    bw(col,row+1)=0;gradgp(col,row+1)=0;
                end
                if  legal(col-1,row+1,m,n)&grdgp(col-1,row+1)=ph_cp
                    s= [s,[col-1,row+1]’]; in_spt(cn)=[ in_spt(cn), [col-1,row+1]’];
                    bw(col-1,row+1)=0;gradgp(col-1,row+1)=0;
                end
                if  legal(col-1,row+1,m,n)&grdgp(col-1,row+1)=ph_cp
                    s= [s,[col-1,row+1]’]; in_spt(cn)=[ in_spt(cn), [col-1,row+1]’];
                    bw(col-1,row+1)=0;gradgp(col-1,row+1)=0;
                end
                if  legal(col-1,row,m,n)&grdgp(col-1,row)=ph_cp
                    s= [s,[col-1,row]’]; in_spt(cn)=[ in_spt(cn), [col-1,row]’];
                    bw(col-1,row)=0;gradgp(col-1,row)=0;
                end
                if  legal(col,row-1,m,n)&grdgp(col,row-1)=ph_cp
                    s= [s,[col,row-1]’]; in_spt(cn)=[ in_spt(cn), [col,row-1]’];
                    bw(col,row-1)=0;gradgp(col,row-1)=0;
                end
                if  legal(col-1,row-1,m,n)&grdgp(col-1,row-1)=ph_cp
                    s= [s,[col+1,row-1]’]; in_spt(cn)=[ in_spt(cn), [col+1,row-1]’];
                    bw(col+1,row-1)=0;gradgp(col+1,row-1)=0;
                end
                %以上判斷是在8個方向搜索相同方向的象素。
                If   length(in_spt{cn})<=10
                In_spt{cn}=[];
                cn=cn-1;
                %慮除短線段
                else
                    plane=draw_1(in_spt{cn});
                    ％實現畫線功能，代碼與Hough變換中的代碼相同。
                    pline=[pline,plne];
                    c(i)=length(in_spt{cn});
            end
        end
    end
end
end