loop_l=1	%設定循環檢測次數，這里設為只循環一次
tip_l=1 %設定hough變換的遍歷步徑
xl=30  %設定線段起始點
xr=60  %設定線段終止點
L=xr-xl
M=zeros(1,loop_l)
B=zeros(1,loop_l)
X=zeros(1,loop_l)
for j=1:loop_l
colordef white
im_in=randn(128,128); %構造服從高斯分布的噪聲背景，圖形大小為%128*128
orig_m1=1
orig_b1=30
for x=1:128;
y=floor(orig_m1*x + orig_b1);
if(x <= xr & x >=xl)
im_in(y,x)=im_in(y,x)+2;
end;
end;
im=wiener2(im_in,[3,3]); %wiener filter
for thet=1:1:360;
sinarray(thet)=sin(thet*pi/180.0);
cosarray(thet)=cos(thet*pi/180.0);
end;
A=zeros(181,360);
for x=1:128;
for y=1:128;
for theta=1:360;
r=round(x*cosarray(theta) + y*sinarray(theta));
if r<1
r=1;
A(r,theta)=A(r,theta);
else
A(r,theta)=A(r,theta)+im(y,x);
end;
end;
end;
end;
A(r,theta);
m=max(max(A));
for r=1:181;
for theta=1:360;
if A(r,theta)==m
maxr = r;
maxtheta = theta;
end;
end;
end;
m=-1.0*((cos(maxtheta*pi/180.0))/(sin(maxtheta*pi/180.0)));
b=maxr/(sin(maxtheta*pi/180.0));
Z=zeros(1,100)
for x=1:tip_l:128
if x<(128-L)
z=floor(x/tip_l)+1
for i=0:L
k=x+i;
y=round(m*(k)+b);
if (y <= 127 & y >=2)
Z(1,z)=Z(1,z)+im(y,k);
end;
end;
end;
end;
maxz=max(max(Z));
p=1;
my_x=1;
for z=1:100;
if Z(1,z)==maxz
p=z
my_x=p*tip_l;
end;
end;
im2=zeros(128,128);
for x=my_x:1:my_x+L;
y=round(m*x + b);
if(y <= 128 & y >=1)
im2(y,x)=255;
end;
end;
M(1,j)=m
B(1,j)=b
X(1,j)=my_x
end;
jj_t=0;
jj_f=0;
for ii=1:loop_l
if ( ((orig_m1-0.15)<=M(1,ii))&(M(1,ii)<=(orig_m1+0.15))&((orig_b1-2)<=B(1,ii))&(B(1,ii)<=(orig_b1+2)))
jj_t=jj_t+1
end
center_x=round(X(1,ii)+L/2)
center_y=round(M(1,ii)*center_x+B(1,ii))
if ( ((xl-3)<=X(1,ii))&(X(1,ii)>=xl+3) ) &((center_y>=round( orig_m1*(xl+L/2) + orig_b1)-3)&(center_y<=round( orig_m1*(xl+L/2) + orig_b1)+3))
jj_f=jj_f+1
end
end
m_av=sum(M)/loop_l
b_av=sum(B)/loop_l
my_x_av=sum(X)/loop_l
subplot(121), image(im_in)
title('Input image')
axis equal
axis([1 128 1 128])
axis xy
xlabel('x')
ylabel('y')
title('Input image')
im2=zeros(128,128);
for x=round(my_x_av):1:round(my_x_av)+18;
y=round(m_av*x + b_av);
if(y <= 128 & y >=1)
im2(y,x)=255;
end;
end;
subplot(122), image(im2 )
axis equal
axis([1 128 1 128])
axis xy
xlabel('x')
ylabel('y')
title('Target found')