霍夫变换与最小二乘法相结合的直线拟合

2003年12月第17卷　第4期南昌航空工业学院学报(自然科学版)

Journal of Nanchang Institute of Aeronautical T echnology (Natural Science ) Dec. 2003　V ol. 17　N o. 4

霍夫变换与最小二乘法相结合的直线拟合

曾接贤, 张桂梅, 储　, 鲁宇明

(南昌航空工业学院土木建筑系, 江西南昌　330034)

[关键词]　霍夫变换; 最小二乘法; 直线拟合

[摘　要]　将霍夫变换与最小二乘法相结合, 研究对实验数据和图像处理中的二值边缘图进行直线拟合的方法。

首先, 用霍夫变换剔除数据点集中的干扰点或噪声, 并将分布在不同直线附近的点分离出来; 然后, 用最小二乘法拟合各直线。该方法既解决了直接使用最小二乘法拟合时, 拟合直线易受干扰点或噪声的影响和数据点分布在多条直线附近而无法拟合的两个问题; 同时也解决了直接使用霍夫变换时, 拟合直线精度不高和直线段有效区间不容易控制的问题。

[中图分类号]　TP391. 41　　[文献标识码]　A 　　[文章编号]　1001-4926(2003) 04-0009-05

Fit Line Using A Method Combined H ough Transform Square

ZE NGJie 2xian , ZHANG G ui 2mei CHU , (Department o f Architecture and Civil Engineering Nanchang f Nanchang 330034)

K ey w ords :H ough trans form ; Least square ; Fit line

Abstract :A new approach to fit line is , H and least square have been combined to process experiment data and the contour of noise in the set of data points have been deleted by using H ough trans 2form , meanwhile lines are separated ; Secondly , lines have been fit by using least square. When fitting line using least square , it s ome problems. On the one hand , it can ’t overcome the interference of interferential points and noise , on the other hand , when data points distributed in the vicinity of a few lines , it is difficult to fit line based on these points. When fitting lines using H ough trans form , it tends to have low precision and difficult to control line in the valid trivial. The presented method can overcome these problems. The effectiveness of the alg orithm is dem onstrated by the experiment.

1　引言

工程上经常要进行直线拟合, 如实验数据分析、统计数据分析、图像处理中将离散形式的二值边缘图拟

合成连续的直线、C AD 系统中获取三维物体图像的线图等。直线拟合在计算数学中有多种方法, 最常用的是平方逼近, 即最小二乘法。最小二乘法的定义[1]是

i =1

∑[y i -φ(x i ) ]2=min

N

(1)

其中(x i , y i ) 为已知数据点; φ(x ) 为拟合函数, 直线拟合时为一次多项式。

显然, 最小二乘法考虑的是已知数据点到拟合函数(直线拟合时为直线) 的距离平方和为最小。这里存在两个问题, 一是已知数据点集中若存在干扰点或噪声时, 拟合函数并不通过最多的数据点, 拟合误差较大; 二是已知数据点集趋向于分布在多条直线附近时, 需事先对已知数据进行分离预处理, 否则拟合结果没有意

[收稿日期]　2003-09-20

[基金项目]　国家自然科学基金(N o :60275037) ; 江西省自然科学基金(N o :0311019) ; 南昌航空工业学院测试技术与控制工程研究中心开

放基金(N o :2002.007)

[作者简介]　曾接贤(1958-) , 男, 硕士学位。现为南昌航空工业学院副教授。主要从事工程图学、计算机图形学和计算机视觉等方面的

研究。

　　　　　　　　　　　　　　　　南昌航空工业学院学报(自然科学版) 　　　　　　　　　　　　　　　　2003年10

义。对数据点集进行分离预处理是一项既难又费时的工作, 有时甚至不可能。

霍夫变换(H ough Trans form ) 是模式识别领域中对二值图像进行直线检测的有效方法[2], 它检测已知点的共线性, 是一种全局性的检测方法。当已知数据点集中存在干扰点和噪声时, 它可以很好地抑制干扰和噪声, 同时它还可以将已知数据点集拟合成多条直线。但是, 霍夫变换的精度却不容易控制, 当实际问题对拟合直线的精度要求较高时, 则不能用霍夫变换。霍夫变换的输出是共线点直线方程的参数, 当需要得到线段时, 还需进一步处理。

本文将霍夫变换与最小二乘法相结合, 研究对实验数据和图像处理中的边缘图进行直线拟合的方法。首先, 用霍夫变换剔除数据点集中的干扰点和噪声, 并将分布在不同直线附近的点分离出来; 然后, 用最小二乘法拟合各直线。该方法既解决了直接使用最小二乘法拟合直线时存在的两个问题; 同时也解决了直接使用霍夫变换时, 拟合直线精度不高和直线段有效区间不容易控制的问题。

2　霍夫变换原理与实现方法

2. 1　霍夫变换原理

1962年,P. V. C. H ough 根据数学对偶性原理提出了检测图像直线的方法[3], 此后, 该方法被不断地研究

和发展, 主要应用于模式识别领域中对二值图像进行直线检测。其原理如图1所示, , 平面直角坐标系中的直线l 表达为

ρ=x cos ; +y sin ; , ρ02π

其中, ρ为l 相对于原点的距离, ; 为ρ与x (2)

根据式

(2) , 直线l p 点的正弦曲线。显然, 若能确定参数空间中的p 点(, 由上述霍夫变换原理可知, 霍夫变换具有如下性质[4]:

(1) 直角坐标系中的一个点映射到参数空间中为一条曲线; (2) 参数空间中的一个点对应直角坐标系中的一条直线; (3) 直角坐标系中的共点线映射到参数空间中为一条曲线;

(4) 直角坐标系中的共线点映射到参数空间中后为一个交于同一点的曲线族。

2. 2　霍夫变换实现方法

工程中的实验数据和图像处理中的二值边缘图, 通常都是离散数据, 因此, 根据霍夫变换性质, 可按下列步骤实现霍夫变换。

(1) 将参数空间量化成m ×n (m 为; 的等份数, n 为ρ的等份数) 个单元, 并设置累加器矩阵Q m ×n ; (2) 给参数空间中的每个单元分配一个累加器Q (i , j ) , 并把累加器的初始值置为零;

(3) 取出直角坐标系中的点(x i , y i ) (i =1, 2, 3…s , s 为直角坐标系中的已知点数) 代入式(2) , 并以量化

后的; 值计算出ρ;

(4) 在参数空间中, 找到; 和ρ所对应的单元, 并将该单元的累加器加1, 即Q (i , j ) =Q (i , j ) +1;

第4期　　　　　　　　曾接贤、张桂梅、储、鲁宇明:霍夫变换与最小二乘法相结合的直线拟合　　　　　　　　　11

(5) 当直角坐标系中的点都经过(3) 、(4) 两步遍历后, 检验参数空间中每个累加器的值, 累加器值最大的

单元所对应的; 和ρ即为直角坐标系中直线方程式(2) 的参数。

由上述霍夫变换过程可知, 如果参数空间中的; 和ρ量化过粗, 则参数空间中的凝聚效果较差, 找不出直线的准确参数; 和ρ; 反之, ; 和ρ量化过细, 那么计算量将增大。

另外, 当直角坐标系中的点分布在R 条直线附近时, 可在第5步检测累加器时, 取出累加器中前R 个值最大的单元所对应的; k 和ρk (k =1, 2, …, R ) , 以; k 和ρk (k =1, 2, …, R ) 为直角坐标系中直线方程式(2) 的参数, 即可同时实现多条直线的拟合。

3　霍夫变换与最小二乘法相结合的直线拟合方法

假设采集到的数据集为M =(x i , y i ) T , (i =1, 2, …, s , s 为数据集中的数据点数) , M 中的数据点分布在

R 条直线附近, 根据实验目的要求给定误差阈值为d k 。

(1) 霍夫变换。将式(2) 改写成

ρk =x i cos ; k +y i sin ; k (k =1, 2, …, s ; k =1, 2, …, R )

根据式(3) , 对M 作霍夫变换, 可得拟合直线的参数(; k , ρk ) 。

3

(2) 找出拟合直线(; k , ρ将式(3) k ) 附近的点集M k 。

(3)

y i =a k x i +(4)

ρ, b k =-sin ; k sin ; k

计算M 中的点到由式(4) 其中a k =-d +a k

2

, (x i , y i ) ∈M , (i , =1, 2, …, s , k =1, 2, …, R )

d ki

(5) (6) (7)

如果则

3

3

(x i , y i ) ]M k (x kj , y kj ) , (j =1, 2, …s 3, s 3≤s )

M k 为符合误差阈值要求的第k 条霍夫变换直线附近的点集。

333

(3) 以点集M k 为拟合数据, 分别拟合各直线, 可得直线方程式(4) 的参数(a k , b k 。以((x kj ) min , y kj ) ∈

M k 和((x kj ) max , y kj ) ∈M k 为端点, 可确定各直线段的区间, 即

33

y kj =a k x kj +b k , (x kj )

min

33

≤x kj ≤(x kj ) max (8)

4　实验

4. 1　实验数据分布于一条直线附近时的直线拟合

已知某实验中得到的实验数据如表1所示。

表1　实验数据

序号

x y

10. 000015. 000013

22. 50005. 000014

34. 000015

45. 000016

59. 000017

6789101112

10. 000013. 000014. 000015. 000020. 000024. 000025. 000018

19

20

21

22

23

24

17. 886817. 886819. 618820. 773524. 010023. 660323. 660326. 547028. 279129. 4338

序号

x y

30. 000034. 000035. 000037. 000040. 000045. 000044. 000017. 500027. 500032. 500042. 500050. 000032. 320535. 207335. 207339. 000038. 094040. 980839. 826127. 990427. 990436. 650636. 650643. 8675

　　(1) 首先对数据进行预处理,

检查它们是否分布于一条直线附近。如图2(a ) 所示, 实验数据分布在一条直线附近。

(2) 将实验数据直接作最小二乘法拟合, 得直线的斜率a 直=0. 6142, 截距b 直=13. 8112, 拟合结果如图2

　　　　　　　　　　　　　　　　南昌航空工业学院学报(自然科学版) 　　　　　　　　　　　　　　　　2003年12

(b ) 所示。

(3) 将实验数据作霍夫变换。霍夫变换时, 参数空间的量化值为:0≤;

1000个单元, 累加器最大值为11, 输出参数为; =120°, ρ=13。将直线方程式(2) 写成斜截式, 代入; =120°, ρ=13, 得直线的斜率a 霍=0. 5774, 截距b 霍=15. 0111。拟合结果如图2(c ) 所示。

(4) 利用第3步霍夫变换结果, 设误差阈值d =5, 剔除干扰点2后, 用最小二乘法拟合, 得直线的斜率a

=0. 5650, 截距b =15. 4886。拟合结果如图2(d ) 所示

。

霍夫变换拟合的直线与最小二乘法拟合的直线斜率和截距不同, 这是因为, 霍夫变换考虑的是拟合直线通过最多的数据点, 而最小二乘法考虑的是所有数据点距离拟合直线的平方和最小, 当数据集中存在干扰点或噪声时, 直接采用最小二乘法拟合, 噪声, 然后采用最小二乘法拟合。, 拟合。

3. 2　2表2　实验数据

序号

x y

110. 000010. 000012110. 00080. 02082340. 562725. 29643426. 723758. 07964541. 961152. 16905673. 923224. 4781

17. 686815. 38241313. 196215. 29002473. 241448. 17833536. 382252. 50334664. 647033. 29785793. 753525. 7308

324. 539820. 18031429. 586620. 66282591. 532373. 19343644. 384047. 88344782. 431228. 80365859. 638219. 3361

432. 462325. 72831535. 297130. 76512652. 802843. 02283762. 278737. 551948

540. 747931. 52991646. 743932. 67642757. 580640. 26443869. 405333. 437349

649. 269937. 49711762. 559649. 85462845. 259540. 79283976. 563229. 304750

759. 603444. 73271876. 859953. 76392955. 385047. 88284085. 698124. 03075121. 926954. 49826289. 163841. 2621

870. 015552. 02331985. 278465. 76243054. 998435. 40434198. 818216. 45585233. 646160. 4338

982. 160560. 527320

10

11

91. 7063101. 803467. 211421

74. 281422

序号

x y

97. 4590104. 162322. 467268. 15733165. 243642. 190942110. 00010. 00005339. 938944. 0990

78. 98513210. 000067. 73504326. 698963. 29015451. 052150. 3844

24. 83353317. 995263. 11904430. 053253. 27055575. 335236. 3646

序号

x y

序号

x y

序号

x y

90. 0434104. 803517. 118918. 63525921. 715440. 7206

15. 88706053. 408068. 6076

69. 97586175. 620075. 1056

序号

x y

　　当对实验数据进行预处理后, 发现其分布在两条直线附近, 如图3(a ) 所示。根据实验目的的要求, 采用两种直线拟合方法。

(1) 直接采用霍夫变换进行直线拟合。霍夫变换时, 参数空间的量化值为:0≤;

第4期　　　　　　　　曾接贤、张桂梅、

储、鲁宇明:霍夫变换与最小二乘法相结合的直线拟合　　　　　　　　　13

360×1000个单元, 累加器最大值分别为13和12。L 1的拟合参数; 1=60°, ρ1=63. 7, L 2的拟合参数; 2=125°, ρ2=2. 5拟合结果如图3(b ) 所示。

(2) 首先采用霍夫变换找出分布在拟合直线附近的点, 然后采用最小二乘法拟合。在寻找拟合直线附近

的点时, 取误差阈值d =6。图3(c ) 为L2的拟合结果, 图3(d ) 为L1的拟合结果

。

4. 3　图像数据的直线拟合

在模式识别、计算机视觉和C AD 系统中, 经常需要获得三维物体图像的线图, 以便后续处理。图4所示是L 型多面体的综合图像及其各棱线的编号, 图像大小为360×400像素

。

(1) 图像预处理

加入1%的模拟噪声后, 采用图5所示的3×3模板算子, 得到如图6(a ) 所示的二值边缘图。

(2) 霍夫变换

ρ

所示, 其中; 和ρ为棱线的法线式参数, a 和b 为棱线的斜截式参数。

(3) 利用霍夫变换结果检测棱线周围点

根据式(5) , 取误差阈值d k =4(像素) , 检测到的棱线周围点如图6(c ) 所示。(4) 最小二乘法拟合

根据第3步检测到的棱线周围点, 采用最小二乘法拟合。拟合结果如图6(d ) 所示, 各棱线的参数如表3所示

。

(下转第40页)

　　　　　　　　　　　　　　　　南昌航空工业学院学报(自然科学版) 　　　　　　　　　　　　　　　　2003年40

[33]S prangle P , Penao J R , and Hafizi B. Phys. Rev. E ,2001,64:026504.[34]Huang C G, and Zhang Y Z , Phys. Rev. A , 2001, 65:015802.[35]Ribeiro F J , and C ohen M L. Phys. Rev. E , 2001,64:046602.[36]Marang os J. Nature ,2000, 406:243.

[37]D ogariu A , K uzmich A , and Wang L J. Phys. Rev. A ,2001,63:053806.[38]D ogariu A , K uzmich A , Cao H , et al. Optics Express ,2001,8(6) :344.[39]K uzmich A , D ogariu A , Wang L J , et al. Phys. Rev. Lett. ,2000,86:3925.[40]Blaauboer M , K ofman A G, K ozhekin A E , et al. Phys. Rev. A ,1998, 57:4905.

[41]Bortman -Arbiv D , Wils on -G ordon A D , and Friedmann H. Phys. Rev. A ,2001,63:043818. [42]Nian -Hua Liu , Shi -Y ao Zhu , H ong Chen , and X iang Wu , Phys. Rev. E , 2002, 65:046607.

[43]Li -G ang Wang , Nian -Hua Liu , Qiang Lin , and Shi -Y ao Zhu , Europhys. Lett. ,2002, 60(6) :834-840. [44]Mandel L , and Walf E. Optical C oherence and Quantum Optics ,Cambridge University Press 1995.

(上接第13页)

3

表3　L 型多面体各棱线拟合数据(; 的单位为度, ρ、b 、b 的单位为像素)

棱线号

;

霍夫变换

ρ

a b a

3

11503771. 7321-404. 0001. 75849

21503111. 7321-272. 0001. 10102670. 2126281. 50320. 2076284. 3248

31502571. 7321-164. 1. 11-28208-1. 8807793. 0513-1. 8952788. 0255

41501901. 30. . -35. 370112-28136-1. 8807639. 6874-1. 8930635. 2823

51. . 1. 772373. 105913-28109-1. 8807582. 1759-1. 9006587. 3746

6. 8. 53820. 211310. 473114-2838-1. 8807430. 9421-1. 8939433. 4725

70. 212686. 23610. 202288. 5856152833-1. 8807279. 7082-1. 8974281. 3009

81021750. 2126171. 09040. 2067172. 7944

最小二乘法

棱线号

b 3

-. -1840. 2126161. 88930. 2103162. 7951

;

霍夫变换

ρ

a b a 3b

3

最小二乘法

5　结论

本文所提出的霍夫变换与最小二乘法相结合的直线拟合, 充分利用了霍夫变换抗噪声能力强和最小二

乘法拟合精度高的特性。该方法既解决了直接使用最小二乘法拟合时, 拟合直线易受干扰点或噪声的影响和数据点分布在多条直线附近而无法拟合的两个问题; 同时也解决了直接使用霍夫变换时, 拟合直线精度不高和直线段有效区间不容易控制的问题。

[参考文献]

[1]聂铁军等. 数值计算方法[M].西安:西北工业大学出版社,1990.

[2]孙丰荣, 刘积仁. 快速霍夫变换算法[J].计算机学报,V ol. 24,N o. 10(2001) :1102～1109. [3]H ough P V C. Method and means of recognizing complex patterns. U. S. Patent 3069654,1962.

[4]王四龙, 宁书年等. 利用霍夫变换机助提取重磁数据构造信息[J].煤炭学报,V ol. 23,N o. 4(1998) :342～346.