CN102564614B - 激光光斑动态测量方法及测量仪 - Google Patents

Abstract

本发明提供激光光斑动态测量方法及测量仪，其中，测量方法包括下列步骤：1)CCD各像素分别探测所接收到的光强；2)将各像素探测到的光强数据记为矩阵I；基于激光光斑光强分布公式，根据实测的光强数据矩阵I及相应各像素的横纵轴坐标值，利用线性最小二乘法得出σ_x，σ_y，x_c，y_c，A和I₀；其中，x_c与y_c分别为光斑中心的横纵轴坐标值，σ_x与σ_y分别为横纵轴1/e半径，A为幅度，I₀为暗噪声；4)通过牛顿迭代法得出参数向量P(k)的修正向量D，计算P(k+1)＝P(k)+D；5)重复执行步骤4)，直至迭代次数达到预设值K或第k次迭代得出的参数向量P(k)的误差小于预设的容忍值，将此时的参数向量P(k)的各元素作为所测得的激光光束参数。本发明还提供了相应的激光光斑动态测量仪。本发明测量速度快、精度高、动态范围大。

Description

激光光斑动态测量方法及测量仪

技术领域

本发明涉及激光光束诊断技术领域，具体地说，本发明涉及一种激光光斑动态测量方法及测量仪。

背景技术

激光光斑动态测量仪是诊断连续或者脉冲激光的激光光束的测量系统，它可以应用于激光光束的优化、激光参数控制、高斯拟合分析、光束准直等各个应用领域。现有技术中的激光光斑动态测量仪主要有两类，一类是用相机作为光束参数测量系统，即相机法测量仪；另一类是用一个可移动的机械狭缝或是用刀片法扫描入射的光束从而获得光束信息，即机械法测量仪。机械法的优势是动态范围很大，并且无论是高功率还是低功率都能准确测量光束参数(包括腰宽、形状、位置、功率、光强分布等)。但是，机械法测量仪的测量速度较慢，难以快速分析出光束的细节信息，难以实现对光束的动态测量。

因此，当前迫切需要一种测量速度快、精度高、动态范围大的激光光斑动态测量方法及测量仪。

发明内容

本发明的目的是提供一种测量速度快、精度高、动态范围大的激光光斑动态测量方法及测量仪。

为实现上述发明目的，本发明提供了一种激光光斑动态测量方法，包括下列步骤：

1)将激光光束照射到CCD上，CCD各像素分别探测所接收到的光强；

2)将CCD各像素探测到的光强数据记为矩阵I；矩阵I中的每个元素分别对应于CCD相应像素所探测的光强数据；

基于激光光斑光强分布公式根据实测的光强数据矩阵I及相应各像素的横纵轴坐标值，利用线性最小二乘法得出σ_x，σ_y，x_c，y_c，A和I₀；其中，I为激光光斑光强变量，x与y为光斑的横纵轴坐标变量，x_c与y_c分别为光斑中心的横纵轴坐标值，σ_x与σ_y分别为横纵轴1/e半径，A为幅度，I₀为暗噪声；

3)将步骤2)得出的σ_x，σ_y，x_c，y_c，A，I₀分别记为σ_x(0)，σ_y(0)，x_c(0)，y_c(0)，A(0)，1₀(0)；

令参数向量 $P\left(k\right)=\left(\begin{array}{c}{\mathrm{\σ}}_x\left(k\right)\\ {\mathrm{\σ}}_y\left(k\right)\\ x_c\left(k\right)\\ y_c\left(k\right)\\ A\left(k\right)\\ I_0\left(k\right)\end{array}\right)=\left(\begin{array}{c}{t}_1\left(k\right)\\ t_2\left(k\right)\\ t_3\left(k\right)\\ t_4\left(k\right)\\ t_5\left(k\right)\\ t_6\left(k\right)\end{array}\right),$

4)通过牛顿迭代法得出参数向量P(k)的修正向量D，

计算P(k+1)＝P(k)+D；

5)重复执行步骤4)，直至迭代次数达到预设值K或第k次迭代得出的参数向量P(k)的误差小于预设的容忍值，将此时的参数向量P(k)的各元素作为所测得的激光光束参数。

其中，所述步骤2)包括下列子步骤：

21)将光强矩阵I中的各个元素依次记为I_i，将I_i所对应的CCD像素的横坐标记为x_i，将I_i所对应的CCD像素的纵坐标记为y_i，其中i＝1，2，3，...，N-1，N；N为光强矩阵I的总元素个数；

22)

令 $S=\left(\begin{array}{ccccc}{x}_1^2& y_1^2& x_1& y_1& 1\\ x_2^2& y_2^2& x_2& y_2& 1\\ & & M& & \\ & & M& & \\ x_N^2& y_N^2& x_N& y_N& 1\end{array}\right),$ $B=\left(\begin{array}{c}{a}_1\\ a_2\\ a_3\\ a_4\\ a_5\end{array}\right),$ $C=\left(\begin{array}{c}\ln I_1\\ \ln I_2\\ M\\ M\\ \ln I_N\end{array}\right)$

根据公式B＝(S^TS)^-1S^TC，计算出B；

23)计算 ${\mathrm{\σ}}_x=\sqrt{-\frac{1}{a_1}},$ ${\mathrm{\σ}}_y=\sqrt{-\frac{1}{a_2}},$ $x_c=\frac{{\mathrm{\σ}}_x^2{a}_3}{2},$ $y_c=\frac{{\mathrm{\σ}}_y^2{a}_4}{2},$ $A=\exp (a_5+\frac{x_c^2}{{\mathrm{\σ}}_x^2}+\frac{y_c^2}{{\mathrm{\σ}}_y^2}),$

其中，x_c与y_c分别为光斑中心的横纵轴坐标值，σ_x与σ_y分别为横纵轴1/e半径；

24)根据激光光斑光强分布公式

进一步计算出I₀；I₀为暗噪声，x与y是像素的横、纵坐标。

其中，步骤4)包括下列子步骤：

41)根据P(k)计算P(k，j)，P(k，j)是在P(k)的基础上，用t_j(k)+δ_p替换t_j(k)后形成的新的参数向量，j为1到6的整数；其中，δ_p为预先设定用于进行参数修正的小量；

42)计算梯度矩阵Gd，Gd(i，j)为梯度矩阵Gd中第i行第j列的元素，其中i是对应于CCD各像素的编号；光强矩阵I中的各个元素依次记为I_i，将I_i所对应的CCD像素的横纵坐标分别记为x_i、y_i；I(x_i，y_i，P(k))为基于激光光斑光强分布公式，利用参数向量P(k)计算出的坐标(x_i、y_i)处的光强；I(x_i，y_i，P(k，j))为基于激光光斑光强分布公式，利用参数向量P(k，j)计算出的坐标(x_i、y_i)处的光强；计算梯度从而得出梯度矩阵Gd；

43)计算D＝(Gd^TGd)^-1Gd^T(I_拟-I_实)，其中I_实为CCD各像素实际测得的N个光强值所组成的N维列向量，I_拟为N个拟合值I(x_i，y_i，P(k))所组成的N维列向量；

44)P(k+1)＝P(k)+D。

其中，所述步骤1)中，所述激光光束经过多重衰减片衰减后再照射到CCD上。

本发明还提供了一种激光光斑动态测量仪，包括CCD和数据处理装置，所述CCD用于接收激光光束照射，CCD各像素分别探测所接收到的光强；所述数据处理装置用于根据CCD各像素所探测的光强，利用前述激光光斑动态测量方法进行数据处理，从而得出所测激光的光束参数。

其中，所述激光光斑动态测量仪还包括衰减片，所述衰减片设置在激光光源与CCD之间。

其中，所述衰减片所述衰减片至少为两片，所有衰减片均安装在套筒内，以在光路中形成多重衰减。

相对于现有技术，本发明具有下列技术效果：

1、本发明测量速度快；

2、本发明测量精度高；

3、本发明的动态范围大，可测量较大能量动态范围的激光光束。

附图说明

图1示出了本发明一个实施例的LBPA的工作示意图；

图2示出了本发明一个实施例中二维高斯拟合的流程图。

具体实施方式

下面结合附图和实施例对本发明作进一步地描述。

根据本发明的一个实施例，提供了一种激光光斑动态测量方法及测量仪。下文中将激光光斑动态测量仪简称为LBPA。如图1所示，LBPA主要包括衰减片、CCD、伺服电路和数据处理装置。本实施例中数据处理装置为PC机，其上安装应用程序，应用程序的核心是下文中将要描述的算法。

本实施例的LBPA使用一个视频相机来采集图像，用电脑来存储图像，分析光束二维强度分布，给出光束的信息。

CCD相机是由一系列二维马赛克状探测器(被称为像素)组合而成。二维马赛克状探测器能快速记录投射到其表面的能量强度，因而能记录光学图样。激光光束的强度分布能够被一个像素一个像素的记录下来，很容易恢复成原始的图像。

在获得CCD所采集激光光束的原始数据后，可以生成伪彩图以供使用者观察。同时，计算机还需要对原始数据进行数据处理，得到光束分布和光束模式。

下面详细介绍对CCD所采集的原始数据进行数据处理的方法。

一、二维高斯拟合

激光光斑形状一般都是高斯形，对采集到激光光斑图片作二维高斯拟合可以得到光斑的信息。下面给出具体的拟合步骤：

高斯形公式：

$I=\mathrm{Aexp}(-\frac{{(x-x_c)}^2}{{\mathrm{\σ}}_x^2}-\frac{{(y-y_c)}^2}{{\mathrm{\σ}}_y^2})+I_0......\left(1\right)$

其中I为CCD探测到的光强值，A为幅度，x_c与y_c分别为光斑中心的横纵轴坐标值，σ_x与σ_y分别为横纵轴1/e半径，I₀为暗噪声，x与y是像素坐标。假如用1024x1280像素的CCD作为探测器，曝一次光后，本实施例可以得到一个1024x1280矩阵形的数据，这就是公式1左边的矩阵I。

如图2所示，本实施例中，二维高斯拟合分两步完成，下面分别介绍。

拟合的第一步：取对数做多项式拟合。

先忽略I₀，对公式1两边求对数，得

$\ln I_0=-\frac{{(x-x_c)}^2}{{\mathrm{\σ}}_x^2}-\frac{{(y-y_c)}^2}{{\mathrm{\σ}}_y^2}+\ln A$

展开后

lnI₀＝a₁x²+a₂y²+a₃x+a₄y+a₅

其中：

$a_1=-1/{\mathrm{\σ}}_x^2$

$a_2=-1/{\mathrm{\σ}}_y^2$

$a_3=2x_c/{\mathrm{\σ}}_x^2$

$a_4=2y_c/{\mathrm{\σ}}_y^2$

$a_5=-\frac{x_c^2}{{\mathrm{\σ}}_x^2}-\frac{y_c^2}{{\mathrm{\σ}}_y^2}+\ln A$

进行多项式拟合：

SB＝C

其中

$S=\left(\begin{array}{ccccc}{x}_1^2& y_1^2& x_1& y_1& 1\\ x_2^2& y_2^2& x_2& y_2& 1\\ & & M& & \\ & & M& & \\ x_N^2& y_N^2& x_N& y_N& 1\end{array}\right),$ $B=\left(\begin{array}{c}{a}_1\\ a_2\\ a_3\\ a_4\\ a_5\end{array}\right),$ $C=\left(\begin{array}{c}\ln I_1\\ \ln I_2\\ M\\ M\\ \ln I_N\end{array}\right)$

这里N为像素矩阵中像素的个数。N＝1024x1280。

求矩阵B有

S^TAS＝S^TC

B＝(S^TS)^-1S^TC

通过开源库函数OPENCV里面的矩阵运算，很容易求得矩阵B，这样就可求得：

${\mathrm{\σ}}_x=\sqrt{-\frac{1}{a_1}},$ ${\mathrm{\σ}}_y=\sqrt{-\frac{1}{a_2}},$ $x_c=\frac{{\mathrm{\σ}}_x^2{a}_3}{2},$ $y_c=\frac{{\mathrm{\σ}}_y^2{a}_4}{2},$ $A=\exp (a_5+\frac{x_c^2}{{\mathrm{\σ}}_x^2}+\frac{y_c^2}{{\mathrm{\σ}}_y^2})$

这些参数只是近似准确，只是作为下一步拟合的初值。

第二步拟合：非线性最小二乘法拟合。

第(1)步：取i为0到N的数，I写为I(x_i，y_i)

设小量δ_p＝1×10^-6

误差值为ee；

拟合误差的容忍值为tol；

设 $P\left(k\right)=\left(\begin{array}{c}{\mathrm{\σ}}_x\left(k\right)\\ {\mathrm{\σ}}_y\left(k\right)\\ x_c\left(k\right)\\ y_c\left(k\right)\\ A\left(k\right)\\ I_0\left(k\right)\end{array}\right)=\left(\begin{array}{c}{t}_1\left(k\right)\\ t_2\left(k\right)\\ t_3\left(k\right)\\ t_4\left(k\right)\\ t_5\left(k\right)\\ t_6\left(k\right)\end{array}\right)$

第(2)步：

用P(k，j)表示在P(k)中用t_j(k)+δ_p替代t_j(k)后的新参数向量，j为1到6的整数。

求 $\mathrm{Gd}(i,j)=\frac{I(x_i,y_i,P(k,j))-I(x_i,y_i,P\left(k\right))}{{\mathrm{\δ}}_p}$

Gd(i，j)为梯度矩阵Gd中第i行第j列的元素，其中i是对应于CCD各像素的编号；光强矩阵I中的各个元素依次记为I_i，将I_i所对应的CCD像素的横纵坐标分别记为x_i、y_i；I(x_i，y_i，P(k))为基于激光光斑光强分布公式，利用参数向量P(k)计算出的坐标(x_i、y_i)处的光强；I(x_i，y_i，P(k，j))为基于激光光斑光强分布公式，利用参数向量P(k，j)计算出的坐标(x_i、y_i)处的光强；

计算每个梯度Gd(i，j)，即可得出梯度矩阵Gd；当CCD像素个数为1024x1280时，梯度矩阵Gd如下：

$\mathrm{Gd}=\left(\begin{array}{ccccc}\mathrm{Gd}\left(\mathrm{1,1}\right)& \mathrm{Gd}\left(\mathrm{1,2}\right)& ...& ...& \mathrm{Gd}\left(\mathrm{1,6}\right)\\ \mathrm{Gd}\left(\mathrm{2,1}\right)& \mathrm{Gd}\left(\mathrm{2,2}\right)& ...& ...& \mathrm{Gd}\left(\mathrm{2,6}\right)\\ & & M& & \\ & & M& & \\ \mathrm{Gd}\left(1024x\mathrm{1280,1}\right)& \mathrm{Gd}\left(1024x\mathrm{1280,2}\right)& ...& ...& \mathrm{Gd}\left(1024x\mathrm{1280,6}\right)\end{array}\right)$

第(3)步：通过梯度矩阵Gd进一步求出指向参量准确值的方向向量：

D＝(Gd^TGd)^-1Gd^T(I_拟-I_实)，

其中I_实为CCD各像素实际测得的N个光强值所组成的N维列向量，I_拟为N个拟合值I(x_i，y_i，P(k))所组成的N维列向量；

第(4)步：计算P(k+1)＝P(k)+D

上述(2)至(4)步即为一次迭代过程。

通过设定迭代步数及容忍的误差值，就可以控制最后参数准确度。其中，每次叠代后的误差值为

$\mathrm{ee}=\sqrt{\begin{array}{c}{({\mathrm{\σ}}_x(k+1)-{\mathrm{\σ}}_x\left(k\right))}^2+{({\mathrm{\σ}}_y(k+1)-{\mathrm{\σ}}_x\left(k\right))}^2+{(x_c(k+1)-x_c\left(k\right))}^2\\ +{(y_c(k+1)-y_c\left(k\right))}^2+{(A(k+1)-A\left(k\right))}^2+{(I_0(k+1)-I_0\left(k\right))}^2\end{array}}.$

二、测定光强

首先从CCD的手册中找到不同波长下，CCD的量子效率。用一束已知功率的激光照射在CCD上，读出光子数的总和，用功率除以光子数再除以该激光波长所对应的量子效率，最后除以曝光时间，即可得出用于光斑测量测量时的单位值。实际应用中，首先在应用程序的界面输入激光波长，然后输入所用的衰减片的衰减比，在软件的背后程序里就会查表找到相应的量子效率，并且将衰减片在该波长下的衰减比找到，这样用单位值乘以量子效率乘以曝光时间再除以衰减比就得到要测的光强值。

根据本发明的另一个优选实施例，还提供了一种基于上述激光光斑动态测量方法的测量仪(即LBPA)。该LBPA中，为提高动态范围，设计了一组螺纹套筒，这些套筒既能够安放衰减片，又能够将互相拧在一起，使多个衰减片形成级联，从而方便地实现各种衰减组合，如此一来，对于更强的光强，只要选择的衰减片足够好，就都可以用该LBPA进行探测。

本发明能够高速高精度地测出激光的输出模式，可以应用于激光光束的优化、激光参数控制、光束准直等各个领域。

最后，上述的实施例仅用来说明本说明，它不应该理解为是对本发明的保护范围进行任何限制。而且，本领域的技术人员可以明白，在不脱离上述实施例精神和原理下，对上述实施例所进行的各种等效变化、变型以及在文中没有描述的各种改进均在本专利的保护范围之内。

Claims (1)

1.一种激光光斑动态测量方法，包括下列步骤：

1）将激光光束照射到CCD上，CCD各像素分别探测所接收到的光强；其中所述激光光束经过多重衰减片衰减后再照射到CCD上；

2）将CCD各像素探测到的光强数据记为矩阵I；

基于激光光斑光强分布公式 $I=\mathrm{Aexp}(-\frac{{(x-x_c)}^2}{{\mathrm{\σ}}_x^2}-\frac{{(y-y_c)}^2}{{\mathrm{\σ}}_y^2})+I_0,$ 根据实测的光强数据矩阵I及相应各像素的横纵轴坐标值，利用线性最小二乘法得出σ_x，σ_y，x_c，y_c，A和I₀；激光光斑光强分布公式中，I为激光光斑光强变量，x与y为光斑的横纵轴坐标变量，x_c与y_c分别为光斑中心的横纵轴坐标值，σ_x与σ_y分别为横纵轴1/e半径，A为幅度，I₀为暗噪声；

3）将步骤2）得出的σ_x，σ_y，x_c，y_c，A，I₀分别记为σ_x(0)，σ_y(0)，x_c(0)，y_c(0)，A(0)，I₀(0)；

令参数向量 $P\left(k\right)=\left(\begin{array}{c}{\mathrm{\σ}}_x\left(k\right)\\ {\mathrm{\σ}}_x\left(k\right)\\ x_c\left(k\right)\\ y_c\left(k\right)\\ A\left(k\right)\\ I_0\left(k\right)\end{array}\right)=\left(\begin{array}{c}{t}_1\left(k\right)\\ t_2\left(k\right)\\ t_3\left(k\right)\\ t_4\left(k\right)\\ t_5\left(k\right)\\ t_6\left(k\right)\end{array}\right),$

4）通过牛顿迭代法得出参数向量P(k)的修正向量D，

计算P(k+1)=P(k)+D；

5）重复执行步骤4），直至迭代次数达到预设值K或第k次迭代得出的参数向量P(k)的误差小于预设的容忍值，将此时的参数向量P(k)的各元素作为所测得的激光光束参数，从而得出激光光束的模式；

其中，所述步骤2）包括下列子步骤：

21）将光强矩阵I中的各个元素依次记为I_i，将I_i所对应的CCD像素的横坐标记为x_i，将I_i所对应的CCD像素的纵坐标记为y_i，其中i=1，2，3，…，N-1，N；N为光强矩阵I的总元素个数；

22）

令 $S=\left(\begin{array}{ccccc}{x}_1^2& y_1^2& x_1& y_1& 1\\ x_2^2& y_2^2& x_2& y_2& 1\\ & & \·& & \\ & & \·& & \\ & & \·& & \\ & & \·& & \\ & & \·& & \\ & & \·& & \\ x_N^2& y_N^2& x_N& y_N& 1\end{array}\right),$ $\left(\begin{array}{c}{a}_1\\ a_2\\ a_3\\ a_4\\ a_5\end{array}\right),$ $C=\left(\begin{array}{c}\ln I_1\\ \ln I_2\\ \·\\ \·\\ \·\\ \·\\ \·\\ \·\\ \ln I_N\end{array}\right)$

根据公式B=(S^TS)^-1S^TC，计算出B；

23）计算 ${\mathrm{\σ}}_x=\sqrt{-\frac{1}{a_1},}$ ${\mathrm{\σ}}_y=\sqrt{-\frac{1}{a_2},}$ $x_c=\frac{{\mathrm{\σ}}_x^2{a}_3}{2},$ $y_c=\frac{{\mathrm{\σ}}_y^2{a}_4}{2},$ $A=\exp (a_5+\frac{x_c^2}{{\mathrm{\σ}}_x^2}+\frac{y_c^2}{{\mathrm{\σ}}_y^2}),$

其中，x_c与y_c分别为光斑中心的横纵轴坐标值，σ_x与σ_y分别为横纵轴1/e半径；

24）根据激光光斑光强分布公式，进一步计算出I₀；I₀为暗噪声，x与y是像素的横、纵坐标；

其中，所述步骤4）包括下列子步骤：

41）根据P(k)计算P(k,j)，P(k,j)是在P(k)的基础上，用t_j(k)+δ_p替换t_j(k)后形成的新的参数向量，j为1到6的整数；其中，δ_p为预先设定用于进行参数修正的小量；

42）计算梯度矩阵Gd，Gd(i,j)为梯度矩阵Gd中第i行第j列的元素，其中i是对应于CCD各像素的编号；光强矩阵I中的各个元素依次记为I_i，将I_i所对应的CCD像素的横纵坐标分别记为x_i、y_i；I(x_i,y_i,P(k))为基于激光光斑光强分布公式，利用参数向量P(k)计算出的坐标(x_i、y_i)处的光强；I(x_i,y_i,P(k,j))为基于激光光斑光强分布公式，利用参数向量P(k,j)计算出的坐标(x_i、y_i)处的光强；计算梯度 $\mathrm{Gd}(i,j)=\frac{I(x_i,y_i,P(k,j))-I(x_i,y_i,P\left(k\right))}{{\mathrm{\δ}}_p},$ 从而得出梯度矩阵绷；

43)计算D=(Gd^TGd)^-1Gd^T(I_拟-I_实)，其中I_实为CCD各像素实际测得的N个光强值所组成的N维列向量，I_拟为N个拟合值I(x_i,y_i,P(k))所组成的N维列向量；

44)P(k+1)=P(k)+D。

Images