CN102661742A - 基于曲率特征加权质心点约束的自适应标志点布局方法 - Google Patents

基于曲率特征加权质心点约束的自适应标志点布局方法 Download PDF

Info

Publication number
CN102661742A
CN102661742A CN2012101744310A CN201210174431A CN102661742A CN 102661742 A CN102661742 A CN 102661742A CN 2012101744310 A CN2012101744310 A CN 2012101744310A CN 201210174431 A CN201210174431 A CN 201210174431A CN 102661742 A CN102661742 A CN 102661742A
Authority
CN
China
Prior art keywords
coordinate
point
coordinate system
error
under
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN2012101744310A
Other languages
English (en)
Other versions
CN102661742B (zh
Inventor
吕乃光
谭启蒙
娄小平
燕必希
孙鹏
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beijing Information Science and Technology University
Original Assignee
Beijing Information Science and Technology University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beijing Information Science and Technology University filed Critical Beijing Information Science and Technology University
Priority to CN201210174431.0A priority Critical patent/CN102661742B/zh
Publication of CN102661742A publication Critical patent/CN102661742A/zh
Application granted granted Critical
Publication of CN102661742B publication Critical patent/CN102661742B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明提供一种基于曲率特征加权质心点约束的自适应标志点布局方法,所述方法根据特征标志点规划布局的原则,在规定的测量空间中合理规划布局尽可能少的转换点的空间位置信息,使得这些转换点的分布情况既能满足预先给定的转换精度要求,又可以显著地提高特征标志点规划布局效率。

Description

基于曲率特征加权质心点约束的自适应标志点布局方法
技术领域
本发明属于光电测量领域,涉及一种基于曲率特征加权质心点约束的自适应特征标志点布局方法,尤其适用于坐标系转换中的特征标志点规划布局设计法。
背景技术
多年来,坐标转换问题一直是地理测绘、遥感、视觉测量、生物医学、机器人导航等众多研究领域中普遍存在的问题之一。尤其是在飞机、轮船、汽车、涡轮机、发电机等大型设备的制造装配领域,为了获得大尺寸机械部件的完整3D面形,就需要经纬仪、摄影测量***、关节式测量臂、激光跟踪仪等多种测量设备及仪器来组成大尺寸测量***。这是一个十分复杂而繁琐的过程,其中涉及了诸多研究方面,例如:摄像机标定、不同3D点云数据的拼接、测量点云与CAD模型配准、多传感器融合等,鉴于此处仅考虑刚体运动,因此,上述所有问题就是不同坐标系的转换问题,即求解坐标系之间的旋转矩阵和平移矢量。
坐标转换问题的关键之处在于如何解算出准确的转换参数,而转换参数的求解通常需要人为设置若干在各个不同坐标系下具有已知坐标值的特征标志点,通过联立线性方程组来求得转换参数。根据不同的计算功能,特征标志点主要可以分为两部分:转换点集和测试点集。其中,转换点集主要用于计算不同坐标系之间的坐标转换参数;然后,将所得的转换参数值代入测试点集,根据误差评价模型,利用不同的坐标转换误差评价方法,即可得到相应的转换误差评价指标及相关参数。
通常,在工业现场布设大量的标志点参与测量环节,这就需要在全部标志点已知的基础上,利用最小二乘法获取坐标转换最优值。其优点在于,基于大量数据获取全局最优解,可以有效消除随机误差干扰;而缺点是需要大量的标志点,不仅提高了测量成本,同时增加了计算量。然而,在实际装配现场中,由于装配仪器、测量设备众多,极易造成对标志点的遮挡,从而限制了测量标志点的数量;同时,相同数量的标志点布设在操作现场中的分布不同,对于测量精度也会产生不同的影响。
因此,特征标志点的规划布局研究具有十分重要的研究意义和应用价值,作为坐标转换问题的首要前提,它将直接决定着坐标转换误差水平,进而影响最终的测量精度。
目前,许多文献都致力于研究如何分布标志点以提高坐标转换精度的问题。耿娜、邾继贵等在《基于刚体运动学的坐标系配准理论及算法》中指出,标志点对坐标转换精度的影响主要在于标志点的测量误差以及标志点的空间布局方式。该文献在转换参数求解中,推导出标志点测量误差与转换误差之间的关系,并针对于现场测量中标志点分布十分复杂的情况,给出了标志点布设原则:将标志点选定在转换点附近,且均匀分布在较广的范围内,构成不规则的形状,适当地增加标志点数量,可有效提高转换精度。该文献剖析了标志点测量误差与转换误差之间的关系,但对于标志点布局问题,并未进行定量的分析与推导,仅以结论的方式给出。
针对旋转微小角度情况下的坐标系转换问题,Hakan S.Kutoglu等指出了标志点分布对转换精度的影响明显重要于标志点数量对转换精度的影响。此外,转换参数之间的相关性与转换精度无关。然而,此项研究仅适用于旋转微小角度,推导中的一些假设条件无法满足于一般情况,因此,具有一定局限性。
王玉成等提出了坐标转换精度受到标志点测量误差及其分布的影响,从测量误差传播角度,推导了标志点测量误差对于转换精度的影响,而对于标志点的分布情况,则建议标志点应该均匀分布于测量范围内,但这种建议只是概念性的提示,并没有严密的理论推导作为依据。
Tevfik Ayan等指出:在坐标转换参数的计算中,其中一部分误差被估计的转换参数值所吸收,误差的吸收量依赖于标志点的几何分布。利用外部可靠性调整识别并剔除粗大误差对转换参数的影响。Tevfik Ayan等同时指出,误差吸收量的多少主要依赖于标志点个数的冗余量。这里,将标志点的几何分布与粗大误差联系在一起,主要解决了如何识别并剔除粗大误差的影响,但并未深入探讨标志点的几何分布对转换误差的影响。
从上述文献可知,多数文献更多地强调了标志点的测量精度对坐标转换精度的影响;而对于标志点自身的空间分布情况,仅凭借实际经验给出一些定性的结论分析,缺乏***的、完整的理论推导过程。存在这一现象的主要原因在于:标志点的测量误差对转换精度的影响十分直观,易于进行公式推导;而标志点集分布本身存在多种复杂因素,其对转换精度的影响往往是多种因素交织在一起,这极大地增加了推导的难度。然而,标志点的测量精度与标志点的空间分布是一种从属关系,它无法完全表征标志点空间分布的其它特性。
发明内容
本发明要解决的技术问题在于:确保在规定的测量空间范围内,合理规划布局尽可能少的转换点就可满足预先给定的测量精度要求。
根据本发明的一方面,提供一种基于曲率特征加权质心点约束的自适应标志点布局方法,所述方法包括:(1)将实际测量现场作为立方体空间,以所述立方体的几何中心作为坐标原点建立测量空间坐标系,并且选择两个不同的第一观测站位和第二观测站位分别独立获取立方体空间的三维形貌特征,第一观测站位的第一坐标系和第二观测站位的第二坐标系与测量空间坐标系互不相同;(2)选取立方体中的两条体对角线方向处的四个顶点构成初始转换点集,使得初始转换点集的质心点坐标位于测量空间坐标系的坐标原点处;(3)选择不同于转换点且包含在初始转换点集内部的多个标志点构成测试点集,并且使得测试点集的质心点坐标与初始转换点集的质心点坐标彼此重合;(4)利用初始转换点集在第一坐标系和第二坐标系下的坐标值求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量的初值,将测试点集在第一坐标系下的坐标值代入旋转矩阵和平移矢量的初值来计算误差评价指标初值;(5)在立方体测量空间中布设等间距的密集转换点集,计算每个转换点处的曲率特征函数,利用密集转换点集在第一坐标系和第二坐标系下的坐标值求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量,将测试点集在第一坐标系下的坐标值代入旋转矩阵和平移矢量来计算误差评价指标目标值;(6)保持立方体几何中心位置不变,以密集转换点集中各个标志点之间的间距作为预定步长,依次缩小立方体测量空间的三个坐标分量取值,获得多个由大至小的子立方体空间;在每个子立方体中,将所述子立方体中所包含的各个密集标志点的曲率特征函数作为权重因子来求解所述子立方体的加权质心点坐标值,并将加权质心点坐标值与测量空间坐标系原点进行比较;如果加权质心点坐标值与测量空间坐标系原点之间存在偏移量,则将求得的加权质心点作为中心,继续根据所述步长逐步压缩所述子立方体空间,直到所述子立方体空间缩小至预定大小的空间区域为止,将由此获得的子立方体空间作为局部立方体区域;如果加权质心点坐标值与测量空间坐标系原点之间不存在偏移量,则搜索下一个子立方体空间中的局部立方体区域;(7)对于在步骤(6)中所获得的一个局部立方体区域,以所述局部立方体区域的加权质心点为中心,以所述步长的整数倍为球半径,构造一组由小至大的搜索球直至所述局部立方体区域边缘;按照球半径由小至大的顺序,依次对每个球表面及内部且不属于前一个球所包含的密集标志点进行遍历搜索;在每个球的搜索过程中,将每个密集标志点及其关于坐标原点呈中心对称的对称标志点作为一对标志点,每次只将一对标志点引入初始转换点集构成当前转换点集,并按照步骤(4)的方式针对当前转换点集计算误差评价指标当前值,如果计算的误差评价指标当前值小于误差评价指标初值,则将该对标志点引入转换点集而形成更新的转换点集;利用更新的转换点集,继续针对下一对标志点按照相同的方式计算误差评价指标当前值,直至该球体所包含的每个标志点结束,将该次搜索中满足要求的标志点引入到初始转换点集中而构成新转换点集作为下一个球体搜索的初始转换点集,并计算此时的误差评价指标作为下一个球体搜索的误差评价指标初值;重复上述搜索过程,直至该局部立方体区域的每个球体全部搜索完毕为止;(8)对所有局部立方体区域重复执行步骤(7)的遍历搜索操作,寻找满足步骤(7)的条件的转换标志点对,直至计算出的误差评价指标趋近或低于误差评价指标目标值,此时布设的标志点的数目以及相应的坐标位置构成最优转换点集。
在步骤(4)中,利用初始转换点集在第一坐标系和第二坐标系下的坐标值求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量的初值的步骤包括:计算初始转换点集的质心点在第一坐标系和第二坐标系下的坐标值,作为初始转换点集在第一坐标系和第二坐标系下的基准点坐标值;计算初始转换点集中各个标志点与基准点在第一坐标系下的坐标差值,计算初始转换点集中各个标志点与基准点在第二坐标系下的坐标差值;利用所述第一坐标系下的坐标差值以及第二坐标系下的坐标差值,通过最小二乘法计算旋转矩阵的初值;根据计算的旋转矩阵的初值以及初始转换点集在第一坐标系和第二坐标系下的基准点坐标,计算平移矢量的初值。
在步骤(4)中,将测试点集在第一坐标系下的坐标值代入旋转矩阵和平移矢量的初值来计算误差评价指标初值的步骤包括:将测试点集在第一坐标系下的坐标值代入旋转矩阵和平移矢量的初值来计算测试点集中各个测试点在第二坐标系下的计算坐标值;计算测试点集中每一个测试点在第二坐标系下的计算坐标值与所述测试点在第二坐标系下的实际坐标值之差;根据下列方式中的一种来计算误差评价指标初值:(a1)计算所述坐标值之差的标准差和最大绝对误差,作为误差评价指标初值;(a2)计算所述坐标值之差的均方根误差,作为误差评价指标初值;(a3)计算所述坐标值之差的标准差和最大绝对误差,并且计算所述坐标值之差的均方根误差,将所述坐标值之差的标准差和最大绝对误差以及均方根误差作为误差评价指标初值。
在步骤(5)中,利用密集转换点集在第一坐标系和第二坐标系下的坐标值求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量的步骤包括:计算密集转换点集的质心点在第一坐标系和第二坐标系下的坐标值,作为密集转换点集在第一坐标系和第二坐标系下的基准点坐标值;计算密集转换点集中各个标志点与基准点在第一坐标系下的坐标差值,计算密集转换点集中各个标志点与基准点在第二坐标系下的坐标差值;利用所述第一坐标系下的坐标差值以及第二坐标系下的坐标差值,通过最小二乘法计算旋转矩阵;根据计算的旋转矩阵的初值以及密集转换点集在第一坐标系和第二坐标系下的基准点坐标,计算平移矢量。
在步骤(5)中,将测试点集在第一坐标系下的坐标值代入旋转矩阵和平移矢量来计算误差评价指标目标值的步骤包括:将测试点集在第一坐标系下的坐标值代入旋转矩阵和平移矢量来计算测试点集中各个测试点在第二坐标系下的计算坐标值;计算测试点集中每一个测试点在第二坐标系下的计算坐标值与所述测试点在第二坐标系下的实际坐标值之差;根据下列方式中的一种来计算误差评价指标目标值:(b1)计算所述坐标值之差的标准差和最大绝对误差,作为误差评价指标目标值;(b2)计算所述坐标值之差的均方根误差,作为误差评价指标目标值;(b3)计算所述坐标值之差的标准差和最大绝对误差,并且计算所述坐标值之差的均方根误差,将所述坐标值之差的标准差和最大绝对误差以及均方根误差作为误差评价指标目标值。
在步骤(6)中,所述预定大小的空间区域的边长是立方体测量空间的边长的1/5。
在步骤(6)中,如果两个或多个子立方体空间同时包含等于或大于预定百分比的转换标志点,则将所述两个或多个子立方体空间作为同一子立方体空间。
所述预定百分比是85%。
所述自适应标志点布局方法还包括:(9)利用最优转换点集在第一坐标系和第二坐标系下的坐标值求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量的最优值。
在步骤(9)中,利用最优转换点集在第一坐标系和第二坐标系下的坐标值求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量的最优值的步骤包括:计算最优转换点集的质心点在第一坐标系和第二坐标系下的坐标值,作为最优转换点集在第一坐标系和第二坐标系下的基准点坐标值;计算最优转换点集中各个标志点与基准点在第一坐标系下的坐标差值,计算最优转换点集中各个标志点与基准点在第二坐标系下的坐标差值;利用所述第一坐标系下的坐标差值以及第二坐标系下的坐标差值,通过最小二乘法计算旋转矩阵的最优值;根据计算的旋转矩阵的最优值以及最优转换点集在第一坐标系和第二坐标系下的基准点坐标,计算平移矢量的最优值。
附图说明
通过结合附图,从下面的实施例的描述中,本发明这些和/或其它方面及优点将会变得清楚,并且更易于理解,其中:
图1是根据本发明的坐标转换示意图;
图2是根据本发明的基于曲率特征加权质心点约束的自适应特征标志点布局方法的流程图。
具体实施方式
以下,参照附图来详细说明本发明的实施例。
图1是根据本发明的坐标转换示意图。图1中的黑点表示特征标志点。
如图1所示,假设存在两个空间直角坐标系分别为第一坐标系CA和第二坐标系CB,从CA至CB的坐标转换参数可以归结为旋转矩阵R和平移矢量T;鉴于特征标志点主要包括转换点集和测试点集,假定转换点集包含N个转换标志点,它们在上述两个坐标系中具有各自不同的坐标值可以分别表示为:XCAi和XCBi(i=1,2,…,N);测试点集包含M个检验标志点,它们在这两个坐标系中的坐标值均可表示为:XTAj和XTBj(i=1,2,…,M)。因此,在上述两个坐标系下,转换点集和测试点集中的全部标志点均满足等式(1)中的刚体变换坐标转换模型,其具体形式可以表示为:
XCBi=T+R·XCAi,XTBj=T+R·XTAj                     (1)
1.布局参数
转换点集和测试点集的布局参数,其主要可以分为点集内参数和点集间参数:前者包括点个数、基准点坐标和坐标差值;后者则定义初始转换点集和测试点集之间的包络融合度。
1.1点集内参数
1.1.1点个数
根据前面的假设可知,转换点集中包含的点个数为N,测试点集包含的点个数为M。
1.1.2基准点坐标
由于利用转换点集来估算坐标转换参数,这就需要将转换点集的质心点作为基准点。该质心点在第一坐标系CA和第二坐标系CB下的坐标值分别表示为:
X AO = Σ i = 1 N X CAi N , X BO = Σ i = 1 N X CBi N - - - ( 2 )
且上述坐标值可满足如下关系:
XBO=T+R·XAO                               (3)
1.1.3坐标差值
坐标差值是指同一坐标系下各标志点与基准点之间的坐标值之差,其符号即代表标志点与基准点之间的相对位置关系。
在第一坐标系CA下,转换点集和测试点集中各标志点与基准点之间的坐标差值可以表示为:
Xcai=XCAi-XAO,Xtaj=XTAj-XAO              (4)
同理,在第二坐标系CB下,转换点集和测试点集中各标志点与基准点之间的坐标差值可以表示为:
Xcbi=XCBi-XBO,Xtbj=XTBj-XBO              (5)
将等式(4)和(5)代入等式(3)可以推导出:
Xcbi=R·Xcai                      (6)
等式(6)反映了旋转矩阵R与转换点集分别在第一坐标系CA和第二坐标系CB下的坐标差值Xcai和Xcbi之间的函数关系。
1.2点集间布局参数
点集间布局参数则主要表现为转换点集与测试点集之间的包络融合度,且该参数仅与这两个点集各自的分布情况及其相对的位置关系密切相关,而与坐标系无关。现仅以第一坐标系CA为例进行说明。
1.2.1转换点集包络球
在第一坐标系CA下,转换点集包含N个转换点,其中每个转换点的坐标值为XCAi(i=1,2,…,N)。现以转换点集的质心点OCA为中心,以各个转换点到质心点的最大欧氏距离为半径RC,将其所生成的一个能够囊括所有转换点的最小包络球定义为转换点集包络球,它在第一坐标系CA下的球心坐标和球半径RC分别可以写成下列形式:
XA0=(xA0,yA0,zA0)T
R C = max i ( ( x CAi - x A 0 ) 2 + ( y CAi - y A 0 ) 2 + ( z CAi - z A 0 ) 2 ) , ( i = 1,2 , . . . , N ) - - - ( 7 )
1.2.2测试点集包络球
在第一坐标系CA下,测试点集包含M个转换点,其中每个测试点的坐标值为XTAj(i=1,2,…,M)。这里,以测试点集的质心OTA为中心,以各个测试点到点OTA的最大欧氏距离为半径RT,生成一个能够包含所有测试点的最小包络球可定义为测试点集包络球,其在第一坐标系CA下的球心坐标值和球半径RT分别表示为:
X TA 0 = ( x TA 0 , y TA 0 , z TA 0 ) T = ( Σ j = 1 M x TAj M , Σ j = 1 M y TAj M , Σ j = 1 M z TAj M ) T ;
R T = max j ( ( x TAj - x TA 0 ) 2 + ( y TAj - y TA 0 ) 2 + ( z TAj - z TA 0 ) 2 ) , ( j = 1,2 , . . . , M ) - - - ( 8 )
1.2.3包络融合度
转换点集和测试点集之间的包络融合度问题,可以转化为两个包络球之间的相交、重叠程度,现定义两个点集之间的包络交融度为φ,数学表示如下:
φ = ln ( R C + R T O CA O TA ‾ ) = ln ( R C + R T ( x AO - x TAO ) 2 + ( y AO - y TAO ) 2 + ( z AO - z TAO ) 2 ) - - - ( 9 )
等式(9)中,
Figure BDA00001702748500092
表示2个包络球心点之间的欧氏距离。
包络融合度参数Φ可以合理地解释转换点集和测试点集相互交融的程度。由等式(9)可知,包络融合度参数Φ的取值范围是(-∞,+∞)。包络融合度参数Φ具有以下性质:
(1)当RC=RT=0时,φ→-∞;
(2)当RC+RT>0, O CA O TA ‾ = 0 时, φ = ln ( R C + R T O CA O TA ‾ ) = lg ( R C + R T 0 ) → + ∞ ;
(3)当 R C + R T = O CA O TA ‾ 时,φ=ln(1)=0;
(4)当 R C + R T > O CA O TA ‾ 时, R C + R T O CA O TA ‾ > 1 , φ = ln ( R C + R T O CA O TA ‾ ) > 0 ;
(5)当 R C + R T < O CA O TA &OverBar; 时, 0 < R C + R T O CA O TA &OverBar; < 1 , &phi; = ln ( R C + R T O CA O TA &OverBar; ) < 0 .
根据包络融合度的定义及相关性质可知,如果转换点集包络球与测试点集包络球互不相交时,则包络融合度为负数,且随着两个包络球的球心之间的欧氏距离不断增大而逐步减小,直至φ→-∞;如果转换点集包络球与测试点集包络球呈现相切状态,则包络融合度为零;如果转换点集包络球与测试点集包络球出现相交时,则包络融合度为正数,且随着两个包络球的球心之间欧氏距离的不断减小而逐渐增大。当且仅当两个包络球心完全重合时,即转换点集与测试点集此时能够完全相互融合,则包络融合度φ→+∞。
2.坐标转换参数求解
基于刚体变换的坐标转换参数求解主要涉及对旋转矩阵R和平移矢量T的数值估计。通常,将转换点集在第一坐标系CA和第二坐标系CB下的坐标值作为已知条件,来解算坐标转换参数R和T。在计算过程中,需要先解算旋转矩阵R,然后再求取平移矢量T。
2.1旋转矩阵
主要利用最小二乘法求解旋转矩阵R,该方法直接利用最小二乘法原理,求解超定线性方程组的全局最优解,其具体的求解过程现描述如下:
将等式(6)进行向量展开可得:
x cbi y cbi z cbi = r 1 r 2 r 3 r 4 r 5 r 6 r 7 r 8 r 9 &CenterDot; x cai y cai z cai - - - ( 10 )
改写成线性方程组可以表示为:
x cbi = r 1 &CenterDot; x cai + r 2 &CenterDot; y cai + r 3 &CenterDot; z cai y cbi = r 4 &CenterDot; x cai + r 5 &CenterDot; y cai + r 6 &CenterDot; z cai z cbi = r 7 &CenterDot; x cai + r 8 &CenterDot; y cai + r 9 &CenterDot; z cai - - - ( 11 )
现假设待求未知列向量为P,可令P=(r1,r2,r3,r4,r5,r6,r7,r8,r9)T,则等式(11)随即可以整理为以下形式:
x ca 1 y ca 1 z ca 1 0 0 0 0 0 0 0 0 0 x ca 1 y ca 1 z ca 1 0 0 0 0 0 0 0 0 0 x ca 1 y ca 1 z ca 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . x caN y caN z caN 0 0 0 0 0 0 0 0 0 x caN y caN z caN 0 0 0 0 0 0 0 0 0 x caN y caN z caN &CenterDot; r 1 r 2 r 3 r 4 r 5 r 6 r 7 r 8 r 9 = x cb 1 y cb 1 z cb 1 . . . x cbN y cbN z cbN - - - ( 12 )
等式(12)的矩阵形式可以表示为:
U·P=F       (13)
其中,U代表常系数矩阵,其矩阵维数为(3N×9),该矩阵仅与各转换点在第一坐标系CA下的坐标差值有关;而列向量F的维数为(3N×1),它仅与转换点在第二坐标系CB下的坐标差值有关;而列向量P即为待求的未知向量。只有当3N>9,也就是说N>3时,根据最小二乘法,方可求出等式(11)中表示的超定线性方程组最优解为:
P=(UT·U)-1·UT·F                           (14)
随后,应按照旋转矩阵的表示形式对最优解P重新整理可得:
R = P ( 1 ) P ( 2 ) P ( 3 ) P ( 4 ) P ( 5 ) P ( 6 ) P ( 7 ) P ( 8 ) P ( 9 ) - - - ( 15 )
其中,P(i)表示列向量的第i个元素。
注意,根据最小二乘法的前提要求,N应大于3,考虑到N所表示的转换点个数应为整数,这就说明至少需要4个转换点,才能保证等式(14)中的(UT·U)存在逆矩阵,从而解得最优解。
2.2平移矢量
根据等式(3)可知,平移矢量T主要受基准点坐标和旋转矩阵R这两个方面的制约。参照等式(3),平移矢量T可以根据以下等式进行计算:
T=XBO-R·XAO                                (16)
3.坐标转换误差评价体系的建立
为了能够有效衡量前面所解算出的坐标转换参数的精度水平,需要建立一套坐标转换误差评价体系,该坐标转换误差评价体系主要包括:坐标转换误差评价模型和误差评价方法。
3.1坐标转换误差评价模型
3.1.1转换点集-转换点集评价模型
转换点集-转换点集评价模型是一种普遍应用的坐标转换精度评价方式。转换点集作为其中唯一的核心因素,主要涉及两个方面:即旋转矩阵、平移矢量的求解以及坐标转换精度评价。前者采用最小二乘法以得到转换参数最优解;后者则是在前者的基础上,利用随后描述的三种不同的误差评价方法,进行坐标转换精度评价。
实际上,转换点集-转换点集评价模型通过对下列目标函数进行优化,以解决一个非线性的最小二乘问题,则该目标函数可以表示为:
min H ( R , T ) = min &Sigma; i = 1 N | | R ( X CAi ) + T - X CBi | | 2 - - - ( 17 )
其中,H即为需要优化的目标函数,||·||表示该矩阵的欧式范数。
该模型简单明了,操作方便,该模型能够将坐标转换参数的求解与坐标转换精度评价有效地统一在等式(17)中显示的目标函数之中,但坐标转换精度的评价结果容易受到转换点集自身因素(例如,转换点个数、转换点的空间位置分布等因素)的影响,具有一定的局限性。
3.1.2转换点集-测试点集评价模型
转换点集-测试点集评价模型同样反映了一种常用的坐标转换评价方法,它在前一种评价模型的基础之上,额外引入一组固定的且区别于转换点集的测试点集着重致力于坐标转换精度的评价工作。该模型主要涉及两个部分:一方面,基于转换点集利用最小二乘法、奇异值分解法和单位四元数法等方法解算旋转矩阵和平移矢量;另一方面,将解算出的转换参数代入测试点集,利用随后描述的三种误差评价方法,计算相应的误差评价指标,即评价误差、标准差和最大绝对误差等。
转换点集-测试点集评价模型的数学表达式如下:
min H ( R , T ) = &alpha; &CenterDot; min &Sigma; i = 1 N | | R &CenterDot; ( X CAi ) + T - X CBi | | 2 + ( 1 - &alpha; ) &CenterDot; min &Sigma; j = 1 M [ R ( X TAj ) + T - X TBj ] - - - ( 18 )
等式(18)集中体现了上述两部分内容,其中,α表示上述两项在评价模型中所占的权重系数,0<α<1,通常取α=0.5;H即为需要优化的目标函数;||·||表示该矩阵的欧式范数;等号右边的第一项主要表现为利用转换点集求解旋转平移矩阵的最小二乘原理,而第二项则体现了测试点集评价坐标转换精度。该模型基于最小二乘原则,能够在所有转换标志点基础上获得转换参数最优解,并且该组最优解能够使得测试点集中的坐标重投影误差呈现出最小化。
该模型的优点主要体现在:固定选取一组区别于转换点集的测试点集作为评价基准,可以对不同转换点集所求得的相应转换参数,实现统一的精度评价与比较,从而避免了转换点集自身对转换误差评价的影响。
3.2坐标转换误差评价方法
坐标转换精度的高低主要依赖于对坐标转换误差的大小进行衡量和评价。这里,对坐标转换误差的评价主要采用三种评价方法:坐标值误差法,均方根误差法和相对欧氏距离误差法。
为方便推导,现做出以下假设:经计算解得的从第一坐标系CA到第二坐标系CB的转换参数分别为:旋转矩阵Rcal和平移矢量Tcal;将各测试点在第一坐标系CA下的坐标值XTAj代入解算出的转换参数可得出其在第二坐标系CB下的计算坐标值,即
X TBj cal = T cal + R cal &CenterDot; X TAj , ( j = 1,2 , . . . , M ) - - - ( 19 )
3.2.1坐标值误差法
当求解出从第一坐标系CA到第二坐标系CB的转换参数Rcal和Tcal后,可将测试点集在第一坐标系CA下的坐标值代入上述转换参数值,从而求得每个测试点在第二坐标系CB下的计算坐标值再取测试点的计算坐标值与实际坐标值XTBj之差即为测试点集的坐标值误差矩阵,可以表示为:
error co = ( error co 1 , error co 2 , . . . , error coj , . . . , error coM ) T
= ( X TB 1 cal - X TB 1 , X TB 2 cal - X TB 2 , . . . , X TBj cal - X TBj , . . . , X TBM cal - X TBM ) T - - - ( 20 )
等式(20)中,坐标值误差矩阵errorco是一个M×3的矩阵,其中每一列分别代表了测试点关于X、Y、Z轴的坐标值误差,现对每个列使用平均误差、标准差、最大误差、最小误差以及最大绝对误差作为评价坐标转换精度的参数指标。现以X轴为例,列出上述参数指标的计算公式,其余各坐标轴均与此类似。
xmean error co = &Sigma; j = 1 M error co ( j , 1 ) M ;
xstd error co = &Sigma; j = 1 M [ error co ( j , 1 ) - xmean error co ] 2 M ;
x max error co = max j [ error co ( j , 1 ) ] ;
x min error co = min j [ error co ( j , 1 ) ] ;
xabsm error co = max j | error co ( j , 1 ) | ; - - - ( 21 )
3.2.2均方根误差法
均方根误差(RMS Error)作为一种衡量精度的重要参数,它集中反映了误差的离散程度。其计算等式表示如下:
error rms = ( X TB 1 cal - X TB 1 ) 2 + ( X TB 2 cal - X TB 2 ) 2 + . . . + ( X TBj cal - X TBj ) 2 + . . . + ( X TBM cal - X TBM ) 2 M - - - ( 22 )
3.2.3相对欧氏距离误差法
在第二坐标系CB下,在测试点集中任意选择某一个测试点,分别利用各测试点的实际坐标值XTBj及各自的计算坐标值
Figure BDA00001702748500137
求取其余测试点相对于所选择的测试点的欧氏距离向量,可以表示为:
dis t = ( x TBj - x TB 1 ) 2 + ( y TBj - y TB 1 ) 2 + ( z TBj - z TB 1 ) 2 , j=2,3,…,M,t=j-1;
dis=(dis1,dis2,…,dist,…,disM-1)T
dis t cal = ( x TBj cal - x TB 1 cal ) 2 + ( y TBj cal - y TB 1 cal ) 2 + ( z TBj cal - z TB 1 ca 1 ) 2 , j=2,3,…,M,t=j-1;
dis cal = ( dis 1 cal , dis 2 cal , . . . , d t cal , . . . , dis M - 1 cal ) T ; - - - ( 23 )
则相对欧氏距离误差向量即可表示为:
error dis = ( dis 1 cal - dis 1 , dis 2 cal - dis 2 , . . . , dis t cal - dis t , . . . , dis M - 1 cal - dis M - 1 ) - - - ( 24 )
该误差向量的维数为(M-1)×1,对该向量取平均误差、标准差、最大误差、最小误差以及最大绝对误差作为评价坐标转换精度的参数指标,分别表示如下:
mean error dis = &Sigma; t = 1 M - 1 error dis ( t , 1 ) M - 1 ;
std error dis = &Sigma; t = 1 M - 1 [ error dis ( t , 1 ) - mean error dis ] 2 M - 1 ;
max error dis = max t [ error dis ( t , 1 ) ] ;
min error dis = min t [ error dis ( t , 1 ) ] ;
dbsm error dis = max t | error dis ( t , 1 ) | ; - - - ( 25 )
然而,在实际操作中发现,相对欧氏距离误差法普遍适用于转换点集-转换点集误差评价模型,但其对于转换点集-测试点集误差评价模型具有一定的适用范围:当测试点集的相对空间位置固定时,无论转换点集发生怎样的变化,相对欧氏距离误差法无效,其中的误差评价参数指标均为恒定值。下面,将对上述结论做出以下证明:
一方面,在第一坐标系CA下,选择其中某一个测试点作为参考点,其坐标值为XTA0,其余M-1个测试点分别与该参考点之间的相对欧氏距离
Figure BDA00001702748500146
即可表示为:
( d jj A ) 2 = ( X TAjj - X TA 0 ) T &CenterDot; ( X TAjj - X TA 0 ) , ( jj = 1,2 , . . . , M - 1 ) - - - ( 26 )
另一方面,在第二坐标系CB下,所选择的参考点的估计坐标值为
Figure BDA00001702748500148
同理,其余M-1个测试点分别与该参考点之间的相对欧氏距离
Figure BDA00001702748500149
则可表示为:
( d jj Bcal ) 2 = ( X TBjj cal - X TB 0 cal ) T &CenterDot; ( X TBjj cal - X TB 0 cal ) , ( jj = 1,2 , . . . , M - 1 ) - - - ( 27 )
现对两个相对欧氏距离之间的关系作出如下推导:
( d jj Bcal ) 2 = ( X TBjj cal - X TB 0 cal ) T &CenterDot; ( X TBjj cal - X TB 0 cal )
= [ ( T cal + R cal &CenterDot; X TAjj ) - ( T cal + R cal &CenterDot; X TA 0 ) ] T &CenterDot; [ ( T cal + R cal &CenterDot; X TAjj ) - ( T cal + R cal &CenterDot; X TA 0 ) ]
= [ R cal &CenterDot; ( X TAjj - X TA 0 ) ] T &CenterDot; [ R cal &CenterDot; ( X TAjj - X TA 0 ) ]
= ( X TAjj - X TA 0 ) T &CenterDot; ( R cal ) T &CenterDot; R cal &CenterDot; ( X TAjj - X TA 0 )
= ( X TAjj - X TA 0 ) T &CenterDot; ( R cal &CenterDot; ( R cal ) T ) T &CenterDot; ( X TAjj - X TA 0 )
= ( X TAjj - X TA 0 ) T &CenterDot; I &CenterDot; ( X TAjj - X TA 0 ) = ( d jj A ) 2 - - - ( 28 )
在上述推导中,涉及旋转矩阵其中一个性质,即Rcal·(Rcal)T=I。I为单位矩阵。
另外,在第二坐标系CB下,所选择的参考测试点的实际坐标值为XTB0,则其余M-1个测试点分别与该参考点之间的相对欧氏距离
Figure BDA00001702748500151
即可表示为:
( d jj B ) 2 = ( X TBjj - X TB 0 ) T &CenterDot; ( X TBjj - X TB 0 ) , ( jj = 1,2 , . . . , M - 1 ) - - - ( 29 )
根据相对欧氏距离的误差计算公式,可得到相对欧氏距离误差参量为:
error jj = d jj B - d jj Bcal = d jj B - d jj A - - - ( 30 )
从等式(30)可以看出,测试点集中的M-1个测试点的相对欧氏距离误差仅与测试点分别在两个坐标系下的坐标值密切相关,而与坐标转换参数无关。
基于上述推导可以看出,无论转换点集解算出的旋转矩阵、平移矢量具有何种解,则只要测试点集确定,则其相对的欧氏距离就可以获得,也就是说,无论转换点集是否发生改变,都不会影响最终计算出的相对欧氏距离误差评价参数值,即,相对欧氏距离误差评价方法无效。
4.特征标志点规划布局设计的指导原则
特征标志点规划布局的指导原则主要包含初始转换点集选取原则、测试点集选取原则、新转换点集规划布局原则和标志点规划布局约束条件四个方面。
4.1初始转换点集的选取原则
(1)初始转换点集中包含的点个数为4个
(2)将转换点集的质心点作为转换点集的基准点,从而转换点集的基准点坐标位于坐标原点处。也就是说,各初始转换点关于坐标原点呈中心对称分布。
4.2测试点集的选取原则
(1)测试点坐标与转换点坐标必须互不相同;
(2)测试点必须包含在初始点集坐标值范围之内;
(3)测试点集质心坐标与转换点集质心坐标相互重合;
(4)测试点集与转换点集的包络融合度趋于无穷大,即φ→+∞。
4.3新转换点集的规划布局原则
(1)新引入的转换点个数ΔN应为大于或等于1的正整数;
(2)新引入的ΔN个转换点的坐标均值应与初始转换点集的基准点坐标尽可能地相互重合,使得新转换点集仍关于坐标原点呈中心对称分布;其中,当ΔN为奇数时,至少有一个新引入的转换点位于坐标原点处;当ΔN为偶数时,对于每一个新引入的转换点总能找到其关于坐标原点呈中心对称的对应的新引入的转换点。
4.4标志点规划布局的约束条件
4.4.1曲率特征加权质心点约束
鉴于曲率特征是3D测量空间自身所固有的属性,与外界因素无关,因此,曲率特征分布将直接决定着标志点规划布局的疏密程度。为了直观体现曲率特征对标志点规划布局的影响,现引入曲率特征加权质心点概念,可以作为标志点规划布局的其中一个重要的约束条件。
曲率特征加权质心点的定义过程如下:假设存在一组标志点集可表示为X={Xi|i=1,2,…,n},Xi=(xi,yi,zi)T,其中每个标志点处的主曲率设为k1,k2,且每个标志点处的高斯曲率和平均曲率分别可设定为kgas,kavg,则该组标志点集的曲率特征加权质心点坐标XO可以计算得到为如下:
X O = &Sigma; i = 1 n Q ( X i ) &CenterDot; X i &Sigma; i = 1 n Q ( X i ) - - - ( 31 )
其中,Q(Xi)表示每个标志点处的曲率特征函数,其数学定义公式为:
Q ( X i ) = k 1 2 + k 2 2 = ( 2 &CenterDot; k avg ) 2 - 2 k gas
k avg = k 1 + k 2 2
kgas=k1·k2                           (32)
从等式(31)和式(32)可以看出,将各标志点处的曲率特征作为求解质心点坐标的权重系数,加权质心点位置将会明显靠近曲率变化明显的区域。基于刚体运动学原理,需要将标志点集视为刚体,而具有剧烈曲率变化的标志点所在位置对于坐标转换精度的影响更为显著,这就需要在曲率变化明显的区域布置相对较多的标志点,以保证坐标转换精度。因此,曲率特征加权质心点约束可用于确定具有显著曲率特征的空间区域,从而为标志点布局奠定有利基础。
4.4.2坐标转换误差的平方和最小约束
根据最小二乘法基本原理可知,在求解坐标转换参数的过程中,要求获得的最优解是基于所有参与计算的转换点,使得全部转换点的坐标转换误差的平方和最小。
4.4.3测试点的绝对误差指标最小约束
这里主要采用坐标值误差法、均方根误差法和相对欧氏距离误差法等三种误差评价方法对坐标转换精度进行评定。坐标值误差法和相对欧氏距离误差法主要涉及两个参数指标(即,标准差和最大绝对误差),而均方根误差法主要涉及均方根误差。这里研究特征标志点规划布局的最终目的就在于使得标准差和最大绝对误差以及/或者均方根误差降至最小值。然而,特征标志点规划布局设计是基于同一组测试点集进行的,因此,相对欧氏距离误差法处于无效状态,将不能够参与坐标转换误差评价,因此,本发明仅采用坐标值误差法和均方根误差法。
下面将详细描述根据本发明的基于曲率特征加权质心点约束的自适应特征标志点布局方法。基于曲率特征加权质心点约束的自适应特征标志点布局方法的目的在于,根据前面归纳出的特征标志点规划布局的原则,在规定的测量空间中合理规划布局尽可能少的转换点的空间位置信息,使得这些转换点的分布情况既能满足预先给定的转换精度要求(即,位于允许的转换误差范围之内),又可以显著地提高特征标志点规划布局效率。在此基础上,利用最小二乘法确定坐标系之间的转换参数即可视为最优解。
图2是根据本发明的基于曲率特征加权质心点约束的自适应特征标志点布局方法的流程图。
参照图2,在操作201,将实际的测量现场作为立方体空间,以该立方体的几何中心作为坐标原点建立测量空间坐标系,并确定该立方体空间的各个坐标分量的取值范围。接着,选择两个不同的第一观测站位和第二观测站位分别独立获取立方体空间的3D形貌特征,第一观测站位的坐标系为第一坐标系CA,第二观测站位的坐标系为第二坐标系CB,且第一坐标系CA和第二坐标系CB与测量空间坐标系互不相同,从第一坐标系CA到第二坐标系CB的坐标转换参数即为所求,坐标转换参数主要包括旋转矩阵R和平移矢量T。
在操作202,选取立方体中两条体对角线方向处的四个顶点构成初始转换点集,从而使得初始转换点集的质心点坐标位于测量空间坐标系的坐标原点处。因此,初始转换点集中的转换点在第一坐标系CA和第二坐标系CB下的坐标值分别表示为XCAi和XCBi(i=1,2,3,4),且满足等式(1)中的数学关系。
在操作203,选择不同于转换点且包含于初始转换点集内部的M个标志点构成测试点集,并且使得测试点集的质心点坐标与初始转换点集的质心点坐标彼此重合。
各个测试点在第一坐标系CA和第二坐标系CB下的坐标值分别表示为XTAj和XTBj(i=1,2,…,M),这两组坐标值具有一一对应关系,它们同样满足等式(1)中的关系式。
在操作204,利用初始转换点集在第一坐标系CA和第二坐标系CB下的坐标值求解出第一坐标系CA和第二坐标系CB之间的旋转矩阵R和平移矢量T的初值,将测试点集在第一坐标系CA下的坐标值代入旋转矩阵R和平移矢量T的初值来计算误差评价指标初值。
可分别计算初始转换点集和测试点集的布局参数,布局参数主要可以分为点集内参数和点集间参数:前者包括标志点个数、基准点坐标和坐标差值;后者定义初始转换点集和测试点集之间的包络融合度。其中,初始转换点集包含的点个数为4个,测试点集包含的点个数为M个;基准点坐标在第一坐标系CA和第二坐标系CB下的坐标值分别如等式(2)所示,且二者之间同样满足等式(3)中的关系;初始转换点集、测试点集各自在第一坐标系CA和第二坐标系CB下的坐标差值应分别按照等式(4)和等式(5)进行计算;初始转换点集与测试点集之间的包络融合度则应根据等式(9)进行解算。然后,根据计算出的初始转换点集的布局参数,利用最小二乘法计算旋转矩阵R(参见等式(14)),并且根据旋转矩阵R计算平移矢量T(参见等式(16))。
利用初始转换点集在第一坐标系CA和第二坐标系CB下的坐标值求解出从第一坐标系CA到第二坐标系CB的旋转矩阵R和平移矢量T的初值的具体操作如下:计算初始转换点集的质心点在第一坐标系CA和第二坐标系CB下的坐标值,作为初始转换点集在第一坐标系CA和第二坐标系CB下的基准点坐标值;计算初始转换点集中各个标志点与基准点在第一坐标系CA下的坐标差值,计算初始转换点集中各个标志点与基准点在第二坐标系CB下的坐标差值;利用所述第一坐标系CA下的坐标差值以及第二坐标系CB下的坐标差值,通过最小二乘法计算旋转矩阵R的初值(参见等式(14));根据计算的旋转矩阵R的初值以及初始转换点集在第一坐标系CA和第二坐标系CB下的基准点坐标,计算平移矢量T的初值(参见等式(16))。
在计算误差评价指标初值的过程中,可根据转换点集测试点集误差评价模型(参见等式(18)),采用坐标值误差法和/或均方根误差法,计算相应的标准差和最大绝对误差以及/或者均方根误差,作为误差评价指标初值。
具体地,将测试点集在第一坐标系CA下的坐标值代入旋转矩阵R和平移矢量T的初值来计算测试点集中各个测试点在第二坐标系CB下的计算坐标值;计算测试点集中每一个测试点在第二坐标系CB下的计算坐标值与所述测试点在第二坐标系CB下的实际坐标值之差。然后,根据下列方式中的一种来计算误差评价指标初值:(a1)采用坐标值误差法计算所述坐标值之差的标准差和最大绝对误差(参见等式(21)),作为误差评价指标初值;(a2)采用均方根误差法计算所述坐标值之差的均方根误差(参见等式(22)),作为误差评价指标初值;(a3)采用坐标值误差法计算所述坐标值之差的标准差和最大绝对误差,采用均方根误差法计算所述坐标值之差的均方根误差,所述坐标值之差的标准差和最大绝对误差以及均方根误差作为误差评价指标初值。
如果对误差的精度要求较高,可采用上述方式(a1);如果对误差的离散程度要求较小,可采用上述方式(a2);如果对误差的精度要求较高,并且对误差的离散程度要求较小,可采用上述方式(a3)。
在操作205,在立方体测量空间中布设等间距的密集转换点集,计算每个转换点处的曲率特征函数(参见等式(32)),利用密集转换点集在第一坐标系CA和第二坐标系CB下的坐标值求解从第一坐标系CA到第二坐标系CB的旋转矩阵R和平移矢量T,将测试点集在第一坐标系CA下的坐标值代入旋转矩阵R和平移矢量T来计算误差评价指标目标值。
同样地,可分别计算密集转换点集和测试点集的布局参数。
同样地,利用密集转换点集在第一坐标系CA和第二坐标系CB下的坐标值求解从第一坐标系CA到第二坐标系CB的旋转矩阵R和平移矢量T的具体操作如下:计算密集转换点集的质心点在第一坐标系CA和第二坐标系CB下的坐标值,作为密集转换点集在第一坐标系CA和第二坐标系CB下的基准点坐标值;计算密集转换点集中各个标志点与基准点在第一坐标系CA下的坐标差值,计算密集转换点集中各个标志点与基准点在第二坐标系CB下的坐标差值;利用所述第一坐标系CA下的坐标差值以及第二坐标系CB下的坐标差值,通过最小二乘法计算旋转矩阵R(参见等式(14));根据计算的旋转矩阵R以及密集转换点集在第一坐标系CA和第二坐标系CB下的基准点坐标,计算平移矢量T(参见等式(16))。
同样地,在计算误差评价指标目标值的过程中,可根据转换点集-测试点集误差评价模型(参见等式(18)),采用坐标值误差法和/或均方根误差法,计算相应的标准差和最大绝对误差以及/或者均方根误差,作为误差评价指标目标值。
具体地,将测试点集在第一坐标系CA下的坐标值代入旋转矩阵R和平移矢量T来计算测试点集中各个测试点在第二坐标系CB下的计算坐标值;计算测试点集中每一个测试点在第二坐标系CB下的计算坐标值与所述测试点在第二坐标系CB下的实际坐标值之差。然后,根据下列方式中的一种来计算误差评价指标目标值:(b1)采用坐标值误差法计算所述坐标值之差的标准差和最大绝对误差(参见等式(21)),作为误差评价指标目标值;(b2)采用均方根误差法计算所述坐标值之差的均方根误差(参见等式(22)),作为误差评价指标目标值;(b3)采用坐标值误差法计算所述坐标值之差的标准差和最大绝对误差,采用均方根误差法计算所述坐标值之差的均方根误差,所述坐标值之差的标准差和最大绝对误差以及均方根误差作为误差评价指标目标值。
在操作206,保持立方体几何中心位置不变,以密集转换点集中各个标志点之间的间距作为预定步长,依次缩小立方体测量空间的三个坐标分量取值,获得多个由大至小的子立方体空间;在每个子立方体中,将所述子立方体中所包含的各个密集标志点的曲率特征函数作为权重因子来求解所述子立方体的加权质心点坐标值(参见等式(31)),并将加权质心点坐标值与测量空间坐标系原点进行比较,如果加权质心点坐标值与测量空间坐标系原点之间存在偏移量,则表示该子立方体空间中存在显著的曲率特征区域,随即将求得的加权质心点作为中心,继续根据所述步长逐步压缩所述子立方体空间,直到所述子立方体空间缩小至预定大小的空间区域为止,此时获得的空间区域可称之为局部立方体区域;如果加权质心点坐标值与测量空间坐标系原点之间不存在偏移量,则搜索下一个子立方体空间中可能存在的局部立方体区域。至此,根据给定立方体空间中的曲率特征变化,可以通过上述操作自适应确定具有显著曲率特征的不同局部立方体区域,并计算出每个局部立方体区域的加权质心点坐标。
所述预定大小的空间区域也是一个子立方体空间,所述预定大小的空间区域的边长可以是立方体测量空间的边长的1/5。因此,所述预定大小的空间区域的三个坐标分量的阈值范围均为a/5≤ξ≤a,其中,a为立方体测量空间的边长,ξ是所述预定大小的空间区域的边长。
优选的是,如果两个或多个子立方体空间同时包含等于或大于预定百分比的转换标志点,即可视为同一子立方体空间,因此可加快压缩子立方体空间的速度。优选的是,所述预定百分比可以是85%。
在步骤207,对于在步骤206中所获得的一个局部立方体区域,以所述局部立方体区域的加权质心点为中心,以所述步长的整数倍为球半径,构造一组由小至大的搜索球直至所述局部立方体区域边缘;按照球半径由小至大的顺序,依次对每个球表面及内部且不属于前一个球所包含的密集标志点进行遍历搜索;在每个球的搜索过程中,将每个密集标志点及其关于坐标原点呈中心对称的对称标志点作为一对标志点,每次只将一对标志点引入初始转换点集构成当前转换点集,并按照步骤204的方式针对当前转换点集计算误差评价指标当前值,如果计算的误差评价指标当前值小于误差评价指标初值,则将该对标志点引入转换点集而形成更新的转换点集,如果计算的误差评价指标当前值不小于误差评价指标初值,则不将该对标志点引入转换点集;利用更新的转换点集,继续针对下一对标志点按照相同的方式计算误差评价指标当前值,直至该球体所包含的每个标志点结束,将该次搜索中满足要求的标志点引入到初始转换点集中而构成新转换点集作为下一个球体搜索的初始转换点集,并计算此时的误差评价指标作为下一个球体搜索的误差评价指标初值;重复上述搜索过程,直至该局部立方体区域的每个球体全部搜索完毕为止。因此,如果引入了满足条件的一对标志点,则转换点集被更新,然后利用更新的转换点集,针对下一对标志点按照相同的方式计算误差评价指标当前值。
在步骤208,对所有局部立方体区域重复执行步骤207的遍历搜索操作,寻找满足步骤207的条件的转换标志点对,直至计算出的误差评价指标趋近或低于误差评价指标目标值,并且不再发生显著变化,此时布设的标志点的数目以及相应的坐标位置构成最优转换点集。在此基础上,可利用最优转换点集在第一坐标系和第二坐标系下的坐标值,通过等式(14)和(16)求解从第一坐标系CA到第二坐标系CB的旋转矩阵和平移矢量的最优值。
同样地,利用最优转换点集在第一坐标系CA和第二坐标系CB下的坐标值求解从第一坐标系CA到第二坐标系CB的旋转矩阵R和平移矢量T的最优值的具体操作如下:计算最优转换点集的质心点在第一坐标系CA和第二坐标系CB下的坐标值,作为最优转换点集在第一坐标系CA和第二坐标系CB下的基准点坐标值;计算最优转换点集中各个标志点与基准点在第一坐标系CA下的坐标差值,计算最优转换点集中各个标志点与基准点在第二坐标系CB下的坐标差值;利用所述第一坐标系CA下的坐标差值以及第二坐标系CB下的坐标差值,通过最小二乘法计算旋转矩阵R的最优值(参见等式(14));根据计算的旋转矩阵R的最优值以及最优转换点集在第一坐标系CA和第二坐标系CB下的基准点坐标,计算平移矢量T的最优值(参见等式(16))。
与现有技术相比,根据本发明的基于曲率特征加权质心点约束的自适应标志点布局方法操作简单、快速,无需人为布设大量的特征标志点,即可获得较高的坐标转换精度,可用于坐标转换的特征标志点规划布局设计,进而可普遍应用于多种实际工程领域,具有显著的研究意义和应用价值。
虽然本发明是参照其示例性的实施例被具体描述和显示的,但是本领域的普通技术人员应该理解,在不脱离由权利要求限定的本发明的精神和范围的情况下,可以对其进行形式和细节的各种改变。

Claims (10)

1.一种基于曲率特征加权质心点约束的自适应标志点布局方法,所述方法包括:
(1)将实际测量现场作为立方体空间,以所述立方体的几何中心作为坐标原点建立测量空间坐标系,并且选择两个不同的第一观测站位和第二观测站位分别独立获取立方体空间的三维形貌特征,第一观测站位的第一坐标系和第二观测站位的第二坐标系与测量空间坐标系互不相同;
(2)选取立方体中的两条体对角线方向处的四个顶点构成初始转换点集,使得初始转换点集的质心点坐标位于测量空间坐标系的坐标原点处;
(3)选择不同于转换点且包含在初始转换点集内部的多个标志点构成测试点集,并且使得测试点集的质心点坐标与初始转换点集的质心点坐标彼此重合;
(4)利用初始转换点集在第一坐标系和第二坐标系下的坐标值求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量的初值,将测试点集在第一坐标系下的坐标值代入旋转矩阵和平移矢量的初值来计算误差评价指标初值;
(5)在立方体测量空间中布设等间距的密集转换点集,计算每个转换点处的曲率特征函数,利用密集转换点集在第一坐标系和第二坐标系下的坐标值求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量,将测试点集在第一坐标系下的坐标值代入旋转矩阵和平移矢量来计算误差评价指标目标值;
(6)保持立方体几何中心位置不变,以密集转换点集中各个标志点之间的间距作为预定步长,依次缩小立方体测量空间的三个坐标分量取值,获得多个由大至小的子立方体空间;在每个子立方体中,将所述子立方体中所包含的各个密集标志点的曲率特征函数作为权重因子来求解所述子立方体的加权质心点坐标值,并将加权质心点坐标值与测量空间坐标系原点进行比较;如果加权质心点坐标值与测量空间坐标系原点之间存在偏移量,则将求得的加权质心点作为中心,继续根据所述步长逐步压缩所述子立方体空间,直到所述子立方体空间缩小至预定大小的空间区域为止,将由此获得的子立方体空间作为局部立方体区域;如果加权质心点坐标值与测量空间坐标系原点之间不存在偏移量,则搜索下一个子立方体空间中的局部立方体区域;
(7)对于在步骤(6)中所获得的一个局部立方体区域,以所述局部立方体区域的加权质心点为中心,以所述步长的整数倍为球半径,构造一组由小至大的搜索球直至所述局部立方体区域边缘;按照球半径由小至大的顺序,依次对每个球表面及内部且不属于前一个球所包含的密集标志点进行遍历搜索;在每个球的搜索过程中,将每个密集标志点及其关于坐标原点呈中心对称的对称标志点作为一对标志点,每次只将一对标志点引入初始转换点集构成当前转换点集,并按照步骤(4)的方式针对当前转换点集计算误差评价指标当前值,如果计算的误差评价指标当前值小于误差评价指标初值,则将该对标志点引入转换点集而形成更新的转换点集;利用更新的转换点集,继续针对下一对标志点按照相同的方式计算误差评价指标当前值,直至该球体所包含的每个标志点结束,将该次搜索中满足要求的标志点引入到初始转换点集中而构成新转换点集作为下一个球体搜索的初始转换点集,并计算此时的误差评价指标作为下一个球体搜索的误差评价指标初值;重复上述搜索过程,直至该局部立方体区域的每个球体全部搜索完毕为止;
(8)对所有局部立方体区域重复执行步骤(7)的遍历搜索操作,寻找满足步骤(7)的条件的转换标志点对,直至计算出的误差评价指标趋近或低于误差评价指标目标值,此时布设的标志点的数目以及相应的坐标位置构成最优转换点集。
2.根据权利要求1所述的自适应标志点布局方法,其中,在步骤(4)中,利用初始转换点集在第一坐标系和第二坐标系下的坐标值求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量的初值的步骤包括:
计算初始转换点集的质心点在第一坐标系和第二坐标系下的坐标值,作为初始转换点集在第一坐标系和第二坐标系下的基准点坐标值;
计算初始转换点集中各个标志点与基准点在第一坐标系下的坐标差值,计算初始转换点集中各个标志点与基准点在第二坐标系下的坐标差值;
利用所述第一坐标系下的坐标差值以及第二坐标系下的坐标差值,通过最小二乘法计算旋转矩阵的初值;
根据计算的旋转矩阵的初值以及初始转换点集在第一坐标系和第二坐标系下的基准点坐标,计算平移矢量的初值。
3.根据权利要求2所述的自适应标志点布局方法,其中,在步骤(4)中,将测试点集在第一坐标系下的坐标值代入旋转矩阵和平移矢量的初值来计算误差评价指标初值的步骤包括:
将测试点集在第一坐标系下的坐标值代入旋转矩阵和平移矢量的初值来计算测试点集中各个测试点在第二坐标系下的计算坐标值;
计算测试点集中每一个测试点在第二坐标系下的计算坐标值与所述测试点在第二坐标系下的实际坐标值之差;
根据下列方式中的一种来计算误差评价指标初值:(a1)计算所述坐标值之差的标准差和最大绝对误差,作为误差评价指标初值;(a2)计算所述坐标值之差的均方根误差,作为误差评价指标初值;(a3)计算所述坐标值之差的标准差和最大绝对误差,并且计算所述坐标值之差的均方根误差,将所述坐标值之差的标准差和最大绝对误差以及均方根误差作为误差评价指标初值。
4.根据权利要求1所述的自适应标志点布局方法,其中,在步骤(5)中,利用密集转换点集在第一坐标系和第二坐标系下的坐标值求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量的步骤包括:
计算密集转换点集的质心点在第一坐标系和第二坐标系下的坐标值,作为密集转换点集在第一坐标系和第二坐标系下的基准点坐标值;
计算密集转换点集中各个标志点与基准点在第一坐标系下的坐标差值,计算密集转换点集中各个标志点与基准点在第二坐标系下的坐标差值;
利用所述第一坐标系下的坐标差值以及第二坐标系下的坐标差值,通过最小二乘法计算旋转矩阵;
根据计算的旋转矩阵的初值以及密集转换点集在第一坐标系和第二坐标系下的基准点坐标,计算平移矢量。
5.根据权利要求4所述的自适应标志点布局方法,其中,在步骤(5)中,将测试点集在第一坐标系下的坐标值代入旋转矩阵和平移矢量来计算误差评价指标目标值的步骤包括:
将测试点集在第一坐标系下的坐标值代入旋转矩阵和平移矢量来计算测试点集中各个测试点在第二坐标系下的计算坐标值;
计算测试点集中每一个测试点在第二坐标系下的计算坐标值与所述测试点在第二坐标系下的实际坐标值之差;
根据下列方式中的一种来计算误差评价指标目标值:(b1)计算所述坐标值之差的标准差和最大绝对误差,作为误差评价指标目标值;(b2)计算所述坐标值之差的均方根误差,作为误差评价指标目标值;(b3)计算所述坐标值之差的标准差和最大绝对误差,并且计算所述坐标值之差的均方根误差,将所述坐标值之差的标准差和最大绝对误差以及均方根误差作为误差评价指标目标值。
6.根据权利要求1所述的自适应标志点布局方法,其中,在步骤(6)中,所述预定大小的空间区域的边长是立方体测量空间的边长的1/5。
7.根据权利要求6所述的自适应标志点布局方法,其中,在步骤(6)中,如果两个或多个子立方体空间同时包含等于或大于预定百分比的转换标志点,则将所述两个或多个子立方体空间作为同一子立方体空间。
8.根据权利要求7所述的自适应标志点布局方法,其中,所述预定百分比是85%。
9.根据权利要求1所述的自适应标志点布局方法,还包括:
(9)利用最优转换点集在第一坐标系和第二坐标系下的坐标值求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量的最优值。
10.根据权利要求9所述的自适应标志点布局方法,其中,在步骤(9)中,利用最优转换点集在第一坐标系和第二坐标系下的坐标值求解从第一坐标系到第二坐标系的旋转矩阵和平移矢量的最优值的步骤包括:
计算最优转换点集的质心点在第一坐标系和第二坐标系下的坐标值,作为最优转换点集在第一坐标系和第二坐标系下的基准点坐标值;
计算最优转换点集中各个标志点与基准点在第一坐标系下的坐标差值,计算最优转换点集中各个标志点与基准点在第二坐标系下的坐标差值;
利用所述第一坐标系下的坐标差值以及第二坐标系下的坐标差值,通过最小二乘法计算旋转矩阵的最优值;
根据计算的旋转矩阵的最优值以及最优转换点集在第一坐标系和第二坐标系下的基准点坐标,计算平移矢量的最优值。
CN201210174431.0A 2012-05-30 2012-05-30 基于曲率特征加权质心点约束的自适应标志点布局方法 Expired - Fee Related CN102661742B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201210174431.0A CN102661742B (zh) 2012-05-30 2012-05-30 基于曲率特征加权质心点约束的自适应标志点布局方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201210174431.0A CN102661742B (zh) 2012-05-30 2012-05-30 基于曲率特征加权质心点约束的自适应标志点布局方法

Publications (2)

Publication Number Publication Date
CN102661742A true CN102661742A (zh) 2012-09-12
CN102661742B CN102661742B (zh) 2014-04-02

Family

ID=46771277

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201210174431.0A Expired - Fee Related CN102661742B (zh) 2012-05-30 2012-05-30 基于曲率特征加权质心点约束的自适应标志点布局方法

Country Status (1)

Country Link
CN (1) CN102661742B (zh)

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107860346A (zh) * 2017-09-30 2018-03-30 北京卫星制造厂 一种测量坐标系配准方法
CN109489553A (zh) * 2018-12-27 2019-03-19 中国科学院长春光学精密机械与物理研究所 一种空间标志点库的生成方法、装置、设备及存储介质
CN109859274A (zh) * 2018-12-24 2019-06-07 深圳市银星智能科技股份有限公司 机器人、其物体标定方法及视教交互方法
CN111819838A (zh) * 2018-03-06 2020-10-23 富士胶片株式会社 摄影评价图、摄影评价图生成装置、摄影评价图生成方法及摄影评价图生成程序
CN111819839A (zh) * 2018-03-06 2020-10-23 富士胶片株式会社 摄影装置、摄影方法、摄影程序以及摄影***
CN117492409A (zh) * 2024-01-03 2024-02-02 深圳市钧诚精密制造有限公司 一种五轴倾斜单角度取零点坐标的方法、***及介质

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2007170821A (ja) * 2005-12-19 2007-07-05 Enzan Kobo:Kk 三次元変位計測方法
CN101331379A (zh) * 2005-12-16 2008-12-24 株式会社Ihi 自身位置辨认方法和装置以及三维形状的计测方法和装置
CN101694370A (zh) * 2009-09-15 2010-04-14 北京信息科技大学 大尺寸工业摄影测量***的精度评价方法及基准装置

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101331379A (zh) * 2005-12-16 2008-12-24 株式会社Ihi 自身位置辨认方法和装置以及三维形状的计测方法和装置
JP2007170821A (ja) * 2005-12-19 2007-07-05 Enzan Kobo:Kk 三次元変位計測方法
CN101694370A (zh) * 2009-09-15 2010-04-14 北京信息科技大学 大尺寸工业摄影测量***的精度评价方法及基准装置

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
王玉成等: "坐标转换中公共点选取对于转换精度的影响", 《现代测绘》 *
耿娜等: "基于刚体运动学的坐标系配准理论及算法", 《传感技术学报》 *
赵宝锋等: "坐标转换模型及公共点选取对转换成果精度的影响", 《淮海工学院学报(自然科学版)》 *

Cited By (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107860346A (zh) * 2017-09-30 2018-03-30 北京卫星制造厂 一种测量坐标系配准方法
CN107860346B (zh) * 2017-09-30 2019-12-20 北京卫星制造厂 一种测量坐标系配准方法
CN111819838A (zh) * 2018-03-06 2020-10-23 富士胶片株式会社 摄影评价图、摄影评价图生成装置、摄影评价图生成方法及摄影评价图生成程序
CN111819839A (zh) * 2018-03-06 2020-10-23 富士胶片株式会社 摄影装置、摄影方法、摄影程序以及摄影***
CN111819838B (zh) * 2018-03-06 2022-03-01 富士胶片株式会社 摄影评价图、摄影评价图生成装置、摄影评价图生成方法及摄影评价图生成程序
CN111819839B (zh) * 2018-03-06 2022-06-14 富士胶片株式会社 摄影装置、摄影方法、摄影程序以及摄影***
US11494889B2 (en) 2018-03-06 2022-11-08 Fujifilm Corporation Imaging evaluation map, imaging evaluation map generating device, imaging evaluation map generating method, and imaging evaluation map generating program
US11546506B2 (en) 2018-03-06 2023-01-03 Fujifilm Corporation Imaging apparatus, imaging method, imaging program, and imaging system
CN109859274A (zh) * 2018-12-24 2019-06-07 深圳市银星智能科技股份有限公司 机器人、其物体标定方法及视教交互方法
CN109489553A (zh) * 2018-12-27 2019-03-19 中国科学院长春光学精密机械与物理研究所 一种空间标志点库的生成方法、装置、设备及存储介质
CN117492409A (zh) * 2024-01-03 2024-02-02 深圳市钧诚精密制造有限公司 一种五轴倾斜单角度取零点坐标的方法、***及介质
CN117492409B (zh) * 2024-01-03 2024-03-08 深圳市钧诚精密制造有限公司 一种五轴倾斜单角度取零点坐标的方法、***及介质

Also Published As

Publication number Publication date
CN102661742B (zh) 2014-04-02

Similar Documents

Publication Publication Date Title
CN102661742B (zh) 基于曲率特征加权质心点约束的自适应标志点布局方法
CN106767780B (zh) 基于Chebyshev多项式插值逼近的扩展椭球集员滤波方法
CN104061932B (zh) 一种利用引力矢量和梯度张量进行导航定位的方法
CN104007700B (zh) 一种基于全局敏感度分析的三轴数控机床的关键性几何误差辨识方法
CN102305608B (zh) 多目标二维交叉运动模拟***误差测量补偿方法
CN104050316B (zh) 一种数控机床空间加工误差分布特征分析方法
CN106646645B (zh) 一种重力正演加速方法
CN104156519A (zh) 一种面向加工精度可靠度提升的多轴数控机床几何精度设计方法
CN104504240B (zh) 航天器总装精度测量计算方法
CN104729481B (zh) 一种基于pnp透视模型的合作目标位姿精度测量方法
CN103411589B (zh) 一种基于四维实数矩阵的三维图像匹配导航方法
CN101707026A (zh) 数字地图线状要素化简的组合优化方法
CN101413793A (zh) 复合空间型面几何误差评定方法
Gao et al. Development and calibration of an accurate 6-degree-of-freedom measurement system with total station
CN102735204A (zh) 一种基于弦线的航空薄壁叶片加工扭曲度误差测量方法
CN105953795A (zh) 一种用于航天器表面巡视的导航装置及方法
Krivec et al. Local properties of three-body atomic wave functions
CN105738915A (zh) 三维雷达测量方法及装置
CN101847262A (zh) 一种快速三维点云搜索匹配方法
CN110516350A (zh) 一种基于各向异性加权的ers点误差修正方法
CN104950805B (zh) 一种基于Floyd算法的空间误差补偿方法
CN103808286A (zh) 一种基于全站仪的钢结构三维精度检测分析方法及其应用
CN102087117A (zh) 飞船交会对接用测距敏感器精度的地面测量方法
CN108508463A (zh) 基于Fourier-Hermite正交多项式扩展椭球集员滤波方法
CN104200063A (zh) 机床空间加工误差的非确定性描述及预测方法

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20140402

Termination date: 20160530