离子速度切片成像图像的处理方法

摘要

本发明公开了一种离子速度切片成像图像的处理方法，包括以下步骤：将离子速度切片成像所得图像的强度矩阵进行优化处理，以得到优化的强度矩阵；对所述优化强度矩阵进行矩阵变换，从而获得图像的半径矩阵和角度矩阵；分别计算半径矩阵对应的速度积分和角度矩阵对应的角度积分，以得到离子解离过程中的速度分布和角度分布；以及利用所述角度矩阵对所述角度积分进行拟合，得到离子解离过程中的各向异性参数。本发明能够准确、高效得获得离子速度切片成像图像的预定物理量。

权利要求

1.一种离子速度切片成像图像的处理方法，其特征在于，所述处理方法包括以下步骤：a)将离子速度切片成像所得图像的强度矩阵进行优化处理，以得到优化的强度矩阵；

包括通过离子速度切片成像所得图像的强度矩阵中相邻像素点之间的强度差值确定并去除图像中的奇异点；或包括利用二维图形去卷积去除离子速度切片成像图像存在的速度展宽，以恢复图像原始宽度；

b)对所述优化强度矩阵进行矩阵变换，从而获得图像的半径矩阵和角度矩阵；包括：设置离子速度切片成像所得图像对应的中心坐标为原点；

根据所述优化强度矩阵对应的每个像素强度坐标，根据直角坐标到极坐标变换公式计算其对应的半径和角度；以及由每个像素点对应的半径和角度分别对应构成所述半径矩阵和所述角度矩阵；

c)分别计算所述半径矩阵对应的速度积分和所述角度矩阵对应的角度积分，以得到离子解离过程中的速度分布和角度分布；

以及

d)利用所述角度矩阵对所述角度积分进行拟合，得到离子解离过程中的各向异性参数。

2.如权利要求1所述的处理方法，其特征在于，计算所述半径矩阵对应的速度积分包括：将位于单位半径内且角度不同的所有强度进行求和。

3.如权利要求1所述的处理方法，其特征在于，计算所述角度矩阵对应的角度积分包括：将位于单位角度内且半径不同的所有强度进行求和。

4.如权利要求1所述的处理方法，其特征在于，d)采用最 小二乘算法进行拟合得到所述各向异性参数。

5.如权利要求1所述的处理方法，其特征在于，所述处理方法运行在matlab程序平台上。

说明书

离子速度切片成像图像的处理方法

技术领域

[0001] 本发明属于离子速度成像领域，尤其涉及一种离子速度切片成像数据的处理方法。

背景技术

[0002] 离子速度成像是在飞行时间质谱仪的装置中以位置敏感的微通道板作为探测器，其后连接荧光屏和CCD(电荷耦合器件)相机实现的。通过离子透镜技术，实现对光致电离解离过程中产生的位置不同，但速度、荷质比相同的离子的聚焦成像。离子速度切片成像技术是在离子速度成像的基础上发展起来的，它是通过调节离子透镜的电场，使离子云有一定的展宽，然后再利用超快质量门进行探测，进而得到反映每一时间光致电离解离产物的速度分布和角度分布的图像，即切片成像。

[0003] 现有软件中有很多都能对图像进行处理，但是没有专门针对离子速度切片成像图像的。并且实际中，CCD相机的像素不可能无穷大，所得到的强度图像是离散的矩阵化图形，因此只能通过矩阵化的数据处理方法对图形进行分析，得到其中包含的物理量信息。

[0004] 此外，离子信号在光致电离解离过程中会因速度分布造成图像加宽，这种加宽会对速度分布的计算带来误差。

发明内容

[0005] 本发明的目的旨在至少解决现有技术中的上述问题之一。

[0006] 为此，本发明的实施例提出一种能够高精度获得离子速度切片成像图像预定物理量的图像处理方法。

[0007] 根据本发明的一个方面，本发明实施例提出了一种离子速度切片成像图像的处理方法，所述处理方法包括以下步骤：a)将离子速度切片成像所得图像的强度矩阵进行优化处理，以得到优化的强度矩阵；b)对所述优化强度矩阵进行矩阵变换，从而获得图像的半径矩阵和角度矩阵；c)分别计算所述半径矩阵对应的速度积分和所述角度矩阵对应的角度积分，以得到离子解离过程中的速度分布和角度分布；以及d)利用所述角度矩阵对所述角度积分进行拟合，得到离子解离过程中的各向异性参数。

[0008] 根据本发明进一步的实施例，所述步骤a包括：通过离子速度切片成像所得图像的强度矩阵中相邻像素点之间的强度差值确定并去除图像中的奇异点。

[0009] 根据本发明进一步的实施例，所述步骤a包括：利用二维图形去卷积去除离子速度切片成像图像存在的速度展宽，以恢复图像原始宽度。

[0010] 根据本发明进一步的实施例，所述步骤b包括：设置离子速度切片成像所得图像对应的中心坐标为原点；根据所述优化强度矩阵对应的每个像素强度坐标，根据直角坐标到极坐标变换公式计算其对应的半径和角度；以及由每个像素点对应的半径和角度分别对应构成所述半径矩阵和所述角度矩阵。

[0011] 根据本发明进一步的实施例，计算所述半径矩阵对应的速度积分包括：将位于单位半径内且角度不同的所有强度进行求和。

[0012] 根据本发明进一步的实施例，计算所述角度矩阵对应的角度积分包括：将位于单位角度内且半径不同的所有强度进行求和。

[0013] 根据本发明进一步的实施例，所述步骤d采用最小二乘算法进行拟合得到所述各向异性参数。

[0014] 根据本发明进一步的实施例，所述处理方法运行在matlab程序平台上。

[0015] 本发明通过矩阵化的数据处理方法对图形进行分析，去除实验中的噪音，优化图像，积分并求解出光致电离解离过程中速度分布、角度分布以及各向异性参数β值。本发明不仅提供了一种对离散的矩阵化的离子速度切片成像图像进行处理的方法，还能够快速、准确地得到图像所反应出来的光致电离解离过程中的上述物理量。

[0016] 本发明附加的方面和优点将在下面的描述中部分给出，部分将从下面的描述中变得明显，或通过本发明的实践了解到。

附图说明

[0017] 本发明的上述和/或附加的方面和优点从下面结合附图对实施例的描述中将变得明显和容易理解，其中：

[0018] 图1为本发明离子速度切片成像图像的处理方法的步骤流程图；

[0019] 图2为本发明图像处理方法的优化预处理步骤的实施例流程图；

[0020] 图3为本发明图像处理方法的坐标变换步骤的实施例流程图；

[0021] 图4为本发明图像处理方法的速度积分步骤的实施例流程图；

[0022] 图5为本发明图像处理方法的角度积分步骤的实施例流程图；

[0023] 图6为本发明图像处理方法的各向异性参数β值拟合步骤的实施例流程图；

[0024] 图7为本发明实施例的离子速度切片成像的原始图像示意图；

[0025] 图8为本发明实施例的离子速度切片成像的优化处理图像示意图；

[0026] 图9为本发明实施例的图像处理方法输出的速度积分示意图；

[0027] 图10为本发明实施例的图像处理方法输出的角度积分示意图；以及[0028] 图11为本发明实施例的图像处理方法输出的拟合β值示意图。

具体实施方式

[0029] 下面详细描述本发明的实施例，所述实施例的示例在附图中示出，其中自始至终相同或类似的标号表示相同或类似的元件或具有相同或类似功能的元件。下面通过参考附图描述的实施例是示例性的，仅用于解释本发明，而不能解释为对本发明的限制。

[0030] 参考图1，该图为本发明离子速度切片成像图像的处理方法的步骤流程图。

[0031] 如图1所示，本发明的图像处理方法包括：在利用CCD相机得到的离子速度切片成像所得原始图像，即离散的矩阵化图形之后，首先对该矩阵化数据对应的强度矩阵进行优化处理(步骤102)，以得到优化的强度矩阵。

[0032] 接着，对所述优化强度矩阵对应的处理后数据进行矩阵变换(步骤104)，从而获得图像的半径矩阵和角度矩阵。关于矩阵变换的步骤可以参考下文的具体实施例。

[0033] 然后，对矩阵变换得到的半径矩阵和角度矩阵分别计算其对应的积分，即，对半径矩阵进行对应的速度积分(步骤106)，以及对角度矩阵进行对应的角度积分(步骤108)，以分别得到对应的速度积分矩阵和角度积分矩阵，进而得到速度切片成像图像中所反映出来的离子信号在解离过程中的速度分布和角度分布。

[0034] 在得到角度积分矩阵后，再将角度矩阵对所得的角度积分进行拟合，从而得到离子信号解离过程中的各向异性参数β。

[0035] 下面，将结合具体实施例，对本发明对离子速度切片成像图像的处理方法的原理作出详细说明。

[0036] 参考图2，该图为本发明图像处理方法的优化预处理步骤102的实施例流程图，图2的实施例是确定离子速度切片成像的原始图像I中的奇异点，并消除的步骤。

[0037] 由于大气中宇宙射线无处不在，导致在成像图像中存在十分突兀的点，称之为奇异点。这种奇异点对后续的图像处理带来严重误差，因此本发明需要从成像图像中找出这些奇异点并进行去除。

[0038] 具体来说，如图2所示，首先对原始图像I对应的数据矩阵的下角标i、j，即矩阵每点数据的横、纵坐标值进行初始化(步骤202)。然后，通过i、j加1的方式，从原始图像I矩阵的数据I1，1开始，逐一依次对原始图像I矩阵的各个数据进行处理(步骤204)，其中原图像矩阵可看作强度矩阵。

[0039] 对于矩阵的每点数据Ii，j，判断数据Ii，j与相邻点的差值是否大于预定阈值M(步骤206)，其中M是去宇宙射线的判断值，即差值大于M的认为是奇异点。这里，相邻点是以数据Ii，j为中心与其邻近的几个点，关于相邻点数量的确定，由使用者根据实际需要设定，设定值为相邻点横纵坐标与Ii，j坐标相差的数值，一般设为1。

[0040] 若在步骤206中判断差值大于M，则以该相邻点的平均值赋予数据Ii，j(步骤208)，即通过差值判断的方法，以相邻点之间的强度差值来确定奇异点，并用附近的平均值代替奇异点对应的数值，从而消除奇异点噪音。

[0041] 通过依次将强度矩阵每个点的下角标i，j与图像矩阵默认的大小N进行比较(步骤210)，从而重复步骤204到步骤210，直至完成整个矩阵数据的奇异点搜索和消除。

[0042] 此外，由于超声分子束的飞行方向与飞行时间质谱(TOF)的方向是垂直的，所以分子束在沿气束飞行方向的速度展宽会使成像图像沿气束方向变宽，使图像失真。所以除了图2实施例的奇异点消除步骤之外，本发明的矩阵优化预处理步骤还可以包括利用去卷积的方法来去除速度分布所导致的图像加宽，恢复图像的原始宽度，从而避免速度展宽对后续计算速度分布物理量造成的误差影响。

[0043] 具体地，用去卷积方法去除图像加宽，需要得到加宽矩阵h，这个矩阵h是通过高斯函数拟合某方向(这个方向上除了温度和仪器带来的速度加宽之外，不包括其他反映化学性质的加宽)上的图像强度随坐标的变化得到的，为简便起见，这里每一次取的加宽矩阵h都是原图像矩阵的一列，可以多次取平均，然后用高斯函数拟合，得到拟合函数存储到h.txt文件中，等待后续调用。

[0044] 本发明在现有matlab程序平台的基础上，利用该平台附带的二维图像去卷积工具‘deconvlucy’结合图像矩阵I和加宽矩阵h来去除速度展宽的影响。

[0045] 因此，通过上述去除图像中的噪音、去除卷积的原始图像优化处理，可获得更真实反映原始图像的矩阵数据，从而在后续速度分布、角度分布等物理量的求解处理，能够获得更精确的计算结果。

[0046] 现在参考图3，该图为本发明图像处理方法的坐标变换步骤的实施例流程图。

[0047] 如图所示，对于优化预处理后得到的图像I’，与图2的步骤202和204类似，仍结合步骤302和步骤304从图像I’矩阵的数据I’1，1开始，逐一依次对图像I’矩阵的各个数据进行后续步骤306的处理。即，在步骤306中，根据计算I’i，j这点的半径和角度分别赋予给半径矩阵R中的Ri，j和角度矩阵A中的Ai，j。其中，半径矩阵R和角度矩阵A是分别用于存放半径信息和角度信息的矩阵。

[0048] 具体来说，以图像I’的中心坐标为原点，根据每个像素I’i，j的坐标(x，y)，根据直角坐标到极坐标变换公式，例如由公式 和 计算它们的半径r和角度θ，并分别保存为半径矩阵和角度矩阵。原图像矩阵I’可看作强度矩阵，这三个矩阵的行列相互关联构成每个像素的强度、半径和角度。相当于用了三个二维矩阵，实现三维矩阵的作用

[0049] 类似于图2的步骤210，通过步骤308的判断，相应地重复步骤204到步骤210，从而依次将强度矩阵每个点的下角标i，j进行坐标变换处理，直至完成图像I’对应整个矩阵数据的坐标变换。

[0050] 图4和图5分别显示了本发明图像处理方法的速度积分和角度积分步骤的实施例流程图。

[0051] 如图4所示，对于矩阵坐标变换得到的半径矩阵R和处理后图像矩阵I’，与图2的步骤202和204类似，结合步骤402和步骤404从图像矩阵I’的数据I’1，1和半径矩阵R的数据R1，1开始，逐一依次对图像矩阵I’以及半径矩阵R的各个数据进行后续步骤的处理。即，在步骤406和408中，将0作为初始值赋值给n，或者在循环中对n进行自加1。

[0052] 然后，前进到步骤410，判断Ri，j是否位于循环次数n、n+1与积分单元的乘积之间。积分单元的单位是像素，可以根据使用者需求调节，一般设为1，即为单位半径。步骤410判断的目的是判断当前强度I’i，j所对应的半径Ri，j是在哪个积分单元，以便进行求和。

[0053] 若否，表示强度I’i，j所对应的半径Ri，j不在第n个积分单元内，则返回至步骤408，并重复后续步骤。

[0054] 若是，则将I’i，j累加到速度积分矩阵Vn中(步骤412)。即，把半径在同一个小的区域范围内的(即单位半径)、但角度不同的所有强度进行求和。

[0055] 接着，判断(n+1)与积分单元的乘积是否不小于图像最大半径，即N*N的图像中能形成的最大圆环半径(步骤414)。步骤414用于判断当前第n个积分单元所处的半径位置，以防止出现死循环。

[0056] 若否，表示当前第n个积分单元在图像中，则返回至步骤408。否则，表示当前第n个积分单元已经到达图像边缘，完成一次半径积分求和，并前进至步骤416，将矩阵下角标i，j与矩阵的大小N进行比较，以判断是否完成整个矩阵数据的速度积分处理(步骤416)。

[0057] 这样，通过上述步骤，将位于单位半径内且角度不同的所有强度进行求和，从而计算半径矩阵R对应的速度积分矩阵。

[0058] 在图5所示的角度积分步骤中，所需的是矩阵坐标变换得到的角度矩阵A和处理后图像矩阵I’。并且在图5中，步骤502到步骤508与图4的步骤402到步骤408类似，这里不再赘述。

[0059] 与图4步骤不同的是，图5的步骤510中积分单元单位是角度，即单位角度，并且该步骤执行的是判断角度矩阵中像素数据Ai，j是否位于循环次数n、n+1与积分单元的乘积之间。

[0060] 若是，则将I’i，j累加到角度积分矩阵Tn中(步骤512)。即，把角度在同一个小的区域范围内的(即单位角度)、但半径不同的所有强度进行求和。

[0061] 并且，在步骤514中接着判断(n+1)与积分单元的乘积是否不小于360度，判断当前第n个积分单元所处的角度位置，以防止出现死循环。

[0062] 若否，表示当前第n个积分单元在一次圆周循环中，则返回至步骤508。否则，表示当前第n个积分单元已经完成一次圆周循环，完成一次角度积分求和，并前进至步骤516，将矩阵下角标i，j与矩阵的大小N进行比较，以判断是否完成整个矩阵数据的角度积分处理(步骤516)。

[0063] 这样，通过上述步骤，将位于单位角度内且半径不同的所有强度进行求和，从而计算角度矩阵A对应的角度积分矩阵T。

[0064] 下面，参考图6，该图为本发明图像处理方法的各向异性参数β值拟合步骤的实施例流程图。

[0065] 如图6所示，首先，拟合各向异性参数β值需要角度积分得到的矩阵T，即实验值。步骤602、604和606是对各向异性参数β、归一化系数P和角度θ赋予初始值，在循环中对这三个值进行加法运算。其中各向异性参数β的取值范围是-1到2之间，归一化系数P的取值范围是0.001到10之间，角度θ的取值范围是0到360度之间。

[0066] 在步骤608中，由公式I(θ)＝P(4π)-1[1+βP2(cosθ)]求出对应不同各向异性参数β和归一化系数P的角度积分的理论值(其中P2是二阶Legendre多项式，计算理论值I(θ)与实验值T(θ)的差值平方。并把这些差值平方放入到关于β、P参数的二维向量矩阵C(β，P)中求和，即通过预定间隔依次循环β、P和θ参数，求出每一个β和P对应的所有角度θ的求和。

[0067] 然后，判断角度θ是否小于360度(步骤610)。

[0068] 若否，表示在当前各向异性参数β和归一化系数P下的理论值和实验值的差值平方和还未完成，则返回至步骤606。

[0069] 否则，表示在当前各向异性参数β和归一化系数P下的理论值和实验值的差值平方和已经完成，前进至步骤612，比较各向异性参数β和归一化系数P是否已经超出其取值范围，以判断是否完成对于所有不同各向异性参数β和归一化系数P的角度积分理论值I(θ)与实验值T(θ)差值平方和的计算。

[0070] 这样，通过上述步骤得到了对于所有不同各向异性参数β和归一化系数P的角度积分理论值I(θ)与实验值T(θ)差值平方和，并在步骤614中，选取C矩阵中最小值对应的β值即为拟合值。

[0071] 下面，结合具体的示例以及matlab平台对本发明的离子速度切片成像图像的处理方法进行说明。

[0072] 以氧气离子的某个电子态的振动基态B2∑g-，v＝0为例，在此电子态、振动态上+ +O2 离子解离生成O，即：

[0073]

[0074] 1、导入数据，进行初步处理。

[0075] 将该图像对应的数据文件输入到matlab程序中并运行，从而数据导入界面如图7所示，其中原始图像如图7左图框所示。

[0076] 通过改变“图像最小值增加值”和“图像最大值增加值”，调节图像显示的对比度；“图像旋转角度”是以角度为单位旋转图像，输入完这三个参数后，点击“绘制图像”可在左图框中观察到调节后的结果。为达到比较好的图像效果，以上三步可重复多次。接着，点击“结果输出”按钮，这步是将调整后的数据导入到临时文件，以备后续数据处理。

[0077] 然后，从图像中读出必要的各个参数，包括：

[0078] a、环圆心坐标，用以确定环的中心；

[0079] b、误符合范围坐标，做速度积分和角度积分时需要去掉的范围；

[0080] c、母体离子误符合范围坐标，去宇宙射线是避免去掉的范围；

[0081] d、背景值。

[0082] 至此，数据初步处理已完成。

[0083] 2、图像处理

[0084] 点击上图7中的“图像处理”按钮，进入图像处理界面，如下图8所示：

[0085] 通过适当设置“去宇宙射线的步长”、“去宇宙射线的幅度”、“图像最小值增加值”和“图像最大值减小值”来去除原始实验数据中的宇宙射线，将环移至图像中心，并利用matlab二维图像去卷积。点击“图像处理”，从左图框中可以观察到设置不同参数所得到的结果。达到满意图像后，点击“结果输出”，将数据导入到临时文件，以备后续数据处理，同时跳转界面到图7。至此，图像处理已完成。

[0086] 3、速度积分

[0087] 点击图7中的“速度积分”，进入速度积分界面，如下图9所示：

[0088] 通过设置“速度积分单位长度”、“速度积分角度最小值”“速度积分角度最大值”和“横坐标转化因子”来得到所要求精度和范围的速度积分。“速度积分单位长度”控制精度，“速度积分角度最小值”和“速度积分角度最大值”共同控制速度积分的角度范围。“横坐标转化因子”是用来选择横纵坐标单位，为pixels、m/s或者eV，对结果无本质影响，此因子是由图像校准得到的。以上设置完成后，点击“速度积分”，从左图框中直接可看到积分结果。然后点击“结果输出”，将速度积分结果输出到文件中，同时跳转界面到图7。此时，速度积分已完成。

[0089] 4、角度分布

[0090] 点击图7中的“角度积分”将得到角度分布界面，如下图10所示：

[0091] 与速度积分类似，通过设置“角度积分单位长度”、“角度积分半径最小值”和“角度积分半径最大值”来得到所要求精度和范围的角度积分。以上设置完成后，点击“角度积分”，从左图框中直接可看到积分结果。然后点击“结果输出”，将角度积分结果输出到文件中，同时跳转界面到图7。此时，角度积分已完成。

[0092] 5、β值的拟合

[0093] 点击图7中的“β值拟合”，进入β值拟合界面，如下图11所示：

[0094] 根据拟合β值的需要，选择不同的角度范围。包括：0°到360°、0°到180°和180°到360°。只需点击图框上面的按钮即可完成选择。接着点击“β拟合”，将在图框中出现拟合结果。最后点击“结果输出”，将数值数据输出到文件中。对于不同的速度，可以重复4、5两步，得到不同速度的角度积分及β值。至此，β值的拟合就完成了。

[0095] 本发明能够准确、高效地处理离子速度切片成像图像，并获得图像中反应的物理量。

[0096] 另外，将本发明的处理方法在matlab程序上运行，采用窗口可视化界面，使得各参数的输入、调节简单易行，并且可以时刻看到调节结果。程序中注意保存了每一步的结果输出文件，使得输入调节过程中不需要额外进行记录。程序运行完，可以从输出文件中轻松找到所要的结果。集成了多种处理程序，使对离子速度切片成像图像的处理能够一次完成，大大减少了实验数据处理的工作量，使实验者能轻松、快速、准确得到图像中反应的物理量。

[0097] 当然，显然可知本发明的处理方法不局限于运行在matlab程序平台上，本领域普通技术人员可以依照本发明的处理方法的技术方案相应编程算法来进行图像处理。

[0098] 尽管已经示出和描述了本发明的实施例，对于本领域的普通技术人员而言，可以理解在不脱离本发明的原理和精神的情况下可以对这些实施例进行多种变化、修改、替换和变型，本发明的范围由所附权利要求及其等同限定。