纪振超, 李加福, 杜 华, 朱小平,张 辉, 禹 静

(1. 中国计量大学 计量测试工程学院,浙江 杭州 310018；2. 中国计量科学研究院,北京 100029)

1 引 言

在医用直线加速器周期检定和验收以及治疗计划系统数据采集工作中,三维模体是不可替代的质控工具,其测量结果对于放疗剂量的验证、治疗计划的制定和实施具有重要意义。在对加速器的检测过程中,三维模体放置于升降平台上,探测器初始位置对准射野中心点,由三维模体中的三维运动机构移动探测器测量射线束的水吸收剂量分布[1]。在此过程中,探测器在三维模体中位置的准确性会对测量结果产生影响[2,3],在YY/T 1538-2017标准中[4],规定了探测器沿各轴到达指定位置的误差不应超过0.2 mm。因此,根据三维模体的应用场景,如何实现对三维模体定位精度的便携式原位检测成为迫切需要解决的技术问题。

目前,对三维模体定位精度检测的研究较少,现有的检测方法主要是利用2个激光位移传感器和1个数显游标卡尺测量探测器运动轨迹延长装置的空间位移量[5],该方法可以在三维模体无水和有水状态下进行原位检测,但在检测时,需要搭建机械结构复杂的支架,对传感器的安装要求较高,检测过程较为繁琐,不适合在应用现场检测；因在有水和无水状态下三维模体的定位精度并无明显差异,故对三维模体定位精度的检测一般在无水状态下进行。

若不使用运动轨迹延长装置,则需要在三维模体的上方或玻璃壁外侧检测。若在其上方检测,则需要将其从升降平台上拆卸下来,将检测仪器悬挂于稳定性较高的支架上,该方法可以避免玻璃介质的遮挡,实现无水检测,但无法实现原位检测,且检测过程也较为繁琐；若在其玻璃壁外侧检测,可在有水或无水状态下实现便携式原位检测,但介质的遮挡会影响测量精度。对于多介质下的三维测量,常用的方法有摄影法[6,7]、激光扫描法[8,9]等,文献[6]提出了多介质下基于粒子群标定的双目视觉定位算法,考虑了成像光线在空气与水交界面处的折射,该算法得到的定位误差约1 mm。文献[9]建立了激光扫描仪在多介质下测量时的折射补偿模型,当只考虑玻璃介质时,平均测量误差为0.09 mm,同时考虑玻璃和水介质时,平均测量误差为0.14 mm。而三维模体的标称定位精度一般为0.1 mm,以上方法很难满足检测要求。

针对现有检测方法的不足,本文提出了基于激光跟踪仪折射补偿的原位检测方法,由于水对定位精度并无明显影响,而水对激光的衰减作用较大[10]且靶球在水中易受损,因此,为了提高检测精度和保护靶球,本文选择在无水状态下检测三维模体的定位精度。针对检测过程中激光的折射现象,建立了折射补偿模型,通过实验验证了补偿模型的可行性,并利用蒙特卡洛法评价了本文检测方法的测量不确定度,可实现激光跟踪仪在玻璃介质下的高精度测量与三维模体定位精度便携式原位检测。

2 测量过程及折射补偿模型

2.1 测量过程

激光跟踪仪是通过绝对测距仪(absolute distance meter,ADM)、干涉测距仪(interferometer,IFM)以及两个角度编码器进行跟踪测量。测量仪器坐标系原点到靶球中心的距离l,水平角α和垂直角β,构成球坐标[11]：

(1)

三维模体定位精度检测示意图如图1所示,在检测过程中,靶球在激光跟踪仪的鸟巢位置初始化后,首先在三维模体玻璃壁的外侧断光,然后将靶球固定在三维模体内的运动导轨上并续光,续光时由ADM测量绝对距离,以此作为基准距离,随后由IFM测量靶球移动时的相对位移量[12]。

由于ADM与IFM的激光波长不同,在由空气入射到玻璃后,会分散成折射角不同的2束激光,一般波长越短,相应的折射率越大,发生偏折的程度也越大。由光的折射规律可知,激光束的入射光路将与出射光路互相平行,所有的光路均在同一折射平面内,且该平面与玻璃平面互相垂直[13]。

图1 三维模体定位精度检测示意图Fig.1 3D phantom positioning accuracy detection diagram

图2 靶球位置动态测量过程Fig.2 Dynamic measurement process of target ball position

(2)

2.2 折射补偿模型

2.2.1 ADM折射补偿

ADM测距技术是利用相位法测量激光跟踪仪坐标原点到靶球中心的距离,相位法测距是通过测量调制光在往返过程中相位的变化量来计算位移量[14,15]。在多介质情况下,激光束在玻璃中的波长要小于在空气中的波长,其测量原理图见图3。

图3 相位法测距原理Fig.3 Principle of phase measurement

(3)

ADM激光束在真空中的传播速度为c,空气折射率为na,ADM,玻璃折射率为ng,ADM,则激光往返时间为

(4)

由光的传播规律可知,光在不同的介质中传播时,频率不变,速度发生变化,因此,在t时刻,相位的变化量为

2π N+Δφ1

(5)

式中:f为激光频率;N为相位变化量中包含2π 的整数倍;Δφ1为不足周期波的余相位。

由式(5)可求解待测距离：

(6)

为了解决相位测距的多值性问题,ADM采用零点法[15],即通过改变调制频率f,使Δφ1=0。

设第1个零点对应的调制频率为f1,则测得的距离为

(7)

设第k个零点对应的调制频率为fk,则测得的距离为

(8)

其中,Nk=N1+k-1,则由式(7)、(8)可求得ADM测量的理论值：

(9)

在实际测量过程中,ADM将使用空气折射率计算距离,这将对激光在玻璃中的传播距离产生测量误差,由式(5)和(9)可知,ADM实际测量得到的距离：

(10)

由式(9)和(10)可知,ADM的测量误差为

(11)

2.2.2 IFM折射补偿

激光跟踪仪在使用IFM测距时,由氦氖激光器发射的激光束经分光镜分成反射光束和透射光束,反射光束作为干涉参考光,透射光束作为测量光,测量光经玻璃到达靶球位置,反射回来时,与参考光汇合后产生干涉条纹,见图4。

图4 激光干涉测距原理Fig.4 Principle of laser interferometric ranging

当2束光的光程差为半波长的偶数倍时,合成光的振幅是两个分振幅之和,此时出现明条纹,反之,出现暗条纹。干涉条纹随着靶球的移动而变化,其变化情况通过干涉仪中的光电转换器被转换成电脉冲信号,计数器对处理后的电脉冲信号进行计算得到总脉冲数量,最后将脉冲数转换到靶球位移的变化量[15,16]。

根据激光干涉测距原理,并结合图2,靶球在P0位置时,测量光与参考光的光程差Δ0为

(12)

(13)

在该过程中,光程差的变化量为

(14)

激光的干涉条纹明暗交替变化1次,其光程差变化1个波长,设激光在真空中的波长为λ0,可得干涉条纹明暗交替的变化次数为

(15)

在实际测量过程中,IFM将使用空气中激光的波长计算靶球移动的距离,则IFM实际测得的相对位移量为

(16)

由式(16)可以看出,从P0到P1的IFM测量误差为

(17)

以此类推,当i>1时,从Pi-1到Pi的IFM测量误差为

(18)

2.2.3 靶球位置坐标的求解

设玻璃外表面的平面方程为

Ax+By+Cz+D=0

(19)

其单位法向量n=(A,B,C),A、B、C、D分别为拟合系数,玻璃内表面与外表面平行,玻璃厚度为h,则内表面的平面方程为

Ax+By+Cz+D+h=0

(20)

(21)

(22)

(23)

t=e×b

(24)

(25)

则激光在玻璃中的传播距离为

(26)

(27)

由式(27),可求解Ci坐标,即出射点坐标。

由ADM和IFM的测量误差可知,当i=0时,CiPi的长度为：

(28)

当i>0时,CiPi的长度为：

(29)

(30)

由式(30),可求解Pi坐标,即靶球位置坐标。

3 实 验

3.1 折射补偿模型验证

选用与三维模体玻璃壁材质相同的有机玻璃平板,利用API Radian激光跟踪仪分别测量加玻璃前后的靶球位置坐标,靶球以P0—P1—P2—P3—P0的路径移动,其中P0为初始位置,此位置的绝对距离由ADM测量,之后,由IFM测量靶球的相对位移量,折射补偿模型验证实验如图5所示。

图5 折射补偿模型验证实验Fig.5 Verification experiment of refraction compensation model

在加玻璃平板前,重复上述移动路径3次,取平均值可得到5组坐标测量值和4组相邻点距离测量值,即移动路径P0—P1—P2—P3—P0下每个位置的X、Y、Z坐标值和每两个相邻位置的距离,以此作为基准值。

靶球座的位置固定不变,玻璃平板以5种不同的姿态放置,每种姿态下重复上述移动路径3次,取平均值可得到25组坐标测量值和20组相相邻点距离测量值。

由折射补偿模型可知,在补偿有玻璃介质情况下的测量误差时,需要确定玻璃平板的外平面方程、厚度以及折射率。

为了减小玻璃表面不平整性引起的测量误差,利用激光跟踪仪扫描每种姿态下包含所有入射点的小范围玻璃外表面,平面方程拟合系数(A、B、C、D)见表1。

表1 不同姿态下玻璃平板外平面方程拟合系数Tab.1 Fitting coefficients of equation of the outer glass plane in different postures

利用激光跟踪仪扫描包含所有入射点的玻璃内表面与外表面,拟合平面方程后可求得玻璃的平均厚度为7.272 mm。

ADM激光束的空气折射率na,ADM和IFM激光束的空气折射率na,IFM均设为1,ADM激光束的玻璃折射率ng,ADM可利用式(11)估算,由ADM测量误差以及ADM激光束在玻璃中的传播距离可反解得到ng,ADM。

同样, IFM激光束的玻璃折射率ng,IFM可利用式(18)估算,由IFM的测量误差以及IFM激光束在玻璃中的传播距离变化量可反解得到ng,IFM。

对于ng,ADM的确定, 利用玻璃平板不同摆放姿态下在P0位置时的测距误差和玻璃中的传播距离估算； 对于ng,IFM的确定, 由式(18)可以看出, 当IFM光束在玻璃中的传播距离变化量较大时, 相应的测量误差也较大, 因此, 为了提高估算的准确性, 取5次玻璃中传播距离变化量较大的测量结果来估算ng,IFM, 玻璃折射率估算结果见表2。

表2 玻璃折射率Tab.2 Glass refractive index

基于玻璃平板外平面方程、厚度以及折射率,可求得补偿后的测量点坐标和相邻点距离,补偿前后X、Y、Z坐标以及相邻点距离与无玻璃介质情况下的测量偏差变化如图6所示。由图6可以看出: 补偿前,在玻璃平板同一摆放姿态下,X、Z坐标测量偏差的变化幅度较小,Y坐标和相邻点距离测量偏差的变化幅度则较大。在不同摆放姿态下,测量偏差的变化规律基本相同,其中,X、Z坐标和相邻点距离测量偏差的变化幅度较小,Y坐标测量偏差的变化幅度较大；补偿后,测量偏差在0附近变化,变化幅度都较小。在不同摆放姿态下,测量偏差的变化规律也基本相同。

图6 测量偏差对比Fig.6 Comparision of measurement deviations

补偿前后,X、Y、Z坐标以及相邻点距离的最大测量偏差Dmax和平均测量偏差Dave的变化情况见表3。

从表3可以看出,补偿前后,X、Y、Z坐标的最大测量偏差分别由-3.584 mm、1.366 mm、-1.845 mm减小到-0.058 mm、0.048 mm、-0.063 mm,平均测量偏差分别由3.410 mm、0.407 mm、1.732 mm减小到0.022 mm、0.015 mm、0.035 mm,相邻点距离的最大测量偏差由0.362 mm减小到-0.043 mm,平均测量偏差由0.266 mm减小到0.017 mm。

可见折射补偿模型具有较高的准确性,补偿后与空气中激光跟踪仪的测量精度相当。

3.2 三维模体定位精度检测

选用IBA公司的三维模体为检测对象,其标称定位精度为0.1 mm,X轴和Y轴方向的最大行程为478 mm,Z轴方向的最大行程为410 mm,利用2台API Radian激光跟踪仪分别在无玻璃介质和有玻璃介质情况下,同时测量三维模体运动机构的位移量,三维模体定位精度检测实验如图7所示。

图7 三维模体定位精度检测实验Fig.7 3D phantom positioning accuracy detection experiment

三维模体的定位精度主要分为单轴定位精度和空间定位精度,为了使检测范围能尽量覆盖三维模体的单轴和空间范围,移动路径和测量方案的规划如下：

(1)沿单轴移动,即分别沿X轴、Y轴和Z轴方向移动,移动间距为50 mm,单向移动次数为7次,重复去程与回程3次,取平均值可得到21组相邻点距离测量值；

(2)沿2条空间对角线移动,即同时沿X轴、Y轴和Z轴方向移动,移动间距为50 mm,单向移动次数为7次,重复去程与回程3次,取平均值可得到14组相邻点距离测量值。

在无玻璃介质的情况下,对于悬挂于支架上的激光跟踪仪,首先检测其测量的稳定性,在同一位置,每隔1 s测量1次靶球位置坐标,连续测量60次,将第1个点坐标置零,靶球位置坐标的变化情况如图8所示。

图8 悬挂式激光跟踪仪的测量稳定性Fig.8 Measurement stability test of suspended laser tracker

由图8可以看出,在60 s内,X坐标的稳定性要优于Y、Z坐标,主要原因在于：根据球坐标到笛卡尔坐标的转换公式,在该测量范围内,水平角的波动对其正弦值影响较大,垂直角的波动对其余弦值影响较大。坐标值整体的波动范围在0.014 mm以内,稳定性较好,故以该检测方法为基准。

在有玻璃介质的情况下,同样地,需要确定三维模体玻璃壁的相关参数。首先,通过激光跟踪仪扫描玻璃壁外表面,拟合得到的平面方程系数A、B、C、D分别为0.998 9、0.045 3、0.007 2、854.888 1。其次,通过激光跟踪仪转站测量的方式扫描玻璃壁内表面与外表面,拟合平面方程后可计算得到玻璃壁的平均厚度为14.489 mm。最后,确定三维模体玻璃壁的折射率,三维模体的玻璃壁与上述玻璃平板材质相同,均为有机玻璃,利用转站测量的方式和相同的估算方法得到三维模体玻璃壁的折射率与上述玻璃平板的折射率一致,即ng,ADM取1.491,ng,IFM取1.531。

基于上述移动路径以及三维模体玻璃壁的相关参数,可得三维模体运动机构位移量,无玻璃介质测量值和补偿值,进而可得到两种检测方法测得的定位误差,如图9所示。

图9 定位误差对比Fig.9 Comparison of positioning error

由图9可以看出,本文检测方法与悬挂式检测方法测得的三维模体定位误差基本吻合。

两种检测方法测得的单轴和空间的最大定位误差ΔEmax与平均定位误差ΔEave见表4。

由表4可以看出,本文检测方法与悬挂式检测方法测得的最大定位误差与平均定位误差相当。可进一步说明本文检测方法的有效性。

表4 单轴与空间定位误差Tab.4 Single-axis and spatial positioning errors mm

根据文献[4]中三维模体定位精度的评价方法,以平均定位误差评价三维模体的定位精度,本文方法可得到三维模体的单轴定位精度为0.064 mm,空间定位精度为0.096 mm,可以看出,三维模体的单轴定位精度略高于空间定位精度,实际定位精度略高于标称定位精度。

3.3 测量不确定度分析

采用蒙特卡洛法[17]分析本文检测方法的测量不确定度,测量不确定度来源有激光跟踪仪测距误差εl、水平角测量误差εα、垂直角测量误差εβ、三维模体玻璃壁厚度测量误差εh、ADM激光束玻璃折射率测量误差εng,ADM、IFM激光束玻璃折射率测量误差εng,IFM和玻璃壁外平面方程测量误差εF,各项测量不确定度来源服从的概率分布见表5。

表5 测量不确定度来源的概率分布Tab.5 Probability distribution of sources of measurement uncertainty

输出量为相邻点距离dp,测量模型是由无玻璃介质测量点先求解加玻璃介质测量点再求解补偿点最后计算相邻补偿点距离的过程,测量模型设为fp,则输入量与输出量之间的关系可表示为：

dp=fp(εl,εα,εβ,εh,εng,ADM,εng,IFM,εF)

(31)

设定蒙特卡洛实验次数为10 000次,即各项输入量随机产生10 000个样本值,可得到10 000个输出值,输出量dp的测量误差概率分布如图10所示。

图10 输出量dp误差的概率分布图Fig.10 The probability distribution of dp error

输出量dp的误差估计值μ、标准差σ和95%概率包含区间如表6所示。

表6 蒙特卡洛试验结果Tab.6 Monte Carlo test results mm

4 结 论

本文提出了基于激光跟踪仪折射补偿的三维模体定位精度原位检测方法,建立了激光跟踪仪ADM和IFM测距误差补偿与靶球位置坐标求解模型,并利用蒙特卡洛法分析了本文检测方法的测量不确定度。

实验结果表明:折射补偿模型具有较高的准确性,可实现激光跟踪仪在玻璃介质下的高精度测量,补偿后距离偏差估计值μ=-0.020 mm,标准差σ=0.012 mm,95%概率包含区间为[-0.043 mm,0.003 mm],可满足三维模体定位精度的检测要求,为解决三维模体定位精度的计量需求和量值溯源问题提供了一种有效的方法。与无玻璃介质遮挡的悬挂式检测方法相比,本文的检测方法可以很好地解决三维模体定位精度难以原位检测的问题,与利用探测器运动轨迹延长装置的检测方法相比,本文的原位检测方法具有较高的便携性。