求教，在matlab中利用hough变换检测出了符合条件的点之后怎样把直线提取出来？怎样继续进行三阶贝塞尔曲线的拟合？

HOUGH变换程序部分如下 万分感激
a=180; %角度的值为0到180度
d=round(sqrt(m^2+n^2)); %图像对角线长度为p的最大值
s=zeros(a,2*d); %存储每个(a,p)个数
z=cell(a,2*d);  %用元胞存储每个被检测的点的坐标
for i=1:m 
    for j=1:n%遍历图像每个点
        if(q(i,j)==255)%只检测图像边缘的白点，其余点不检测
            for k=1:a
                p = round(i*cos(pi*k/180)+j*sin(pi*k/180));%对每个点从1到180度遍历一遍，取得经过该点的所有直线的p值（取整）
                if(p > 0)%若p大于0，则将点存储在（d，2d）空间
                    s(k,d+p)=s(k,d+p)+1;%（a，p）相应的累加器单元加一
                    z{k,d+p}=[z{k,d+p},[i,j]'];%存储点坐标
                else
                    ap=abs(p)+1;%若p小于0，则将点存储在（0，d）空间
                    s(k,ap)=s(k,ap)+1;%(a，p)相应的累加器单元加一
                    z{k,ap}=[z{k,ap},[i,j]'];%存储点坐标
                end
            end
        end
    end
end
for i=1:a 
    for j=1:d*2 %检查每个累加器单元中存储数量
        if(s(i,j) >90) %将提取直线的阈值设为90
            lp=z{i,j};%提取对应点坐标
            for k=1:s(i,j)%对满足阈值条件的累加器单元中(a，p)对应的所有点进行操作
               o(lp(1,k),lp(2,k),1)=255; %每个点R分量=255，G分量=0，B分量=0
               o(lp(1,k),lp(2,k),2)=0;
               o(lp(1,k),lp(2,k),3)=0;  %结果为在原图上对满足阈值要求的直线上的点赋红色
            end  
        end 
    end
end