CN103745185B - 一种识别探测器晶体单元位置的方法和装置 - Google Patents

Abstract

本申请公开了一种识别探测器晶体单元位置的方法和装置。该方法包括：获取基于辐射源辐射探测器而生成的位置散点图；计算位置散点图上的峰值点，并基于每两个相邻峰值点确定一个初始谷值点；基于峰值点和初始谷值点对位置散点图进行初次分割，以确定晶体单元的初始边界；分别对各行、列晶体单元初始边界内图像区域进行灰度投影，得到各行、列晶体单元的灰度投影曲线；依据各行、列晶体单元的灰度投影曲线上的极小值点确定目标谷值点；基于峰值点和目标谷值点对位置散点图进行再次分割，以便确定各晶体单元的目标边界。通过本申请的实施方式，可以避免各晶体单元物理特性不一致而导致的谷值点不准确，从而能够准确地识别各晶体单元的位置。

Description

一种识别探测器晶体单元位置的方法和装置

技术领域

本申请涉及核医学领域，特别是涉及一种识别核医学设备的探测器晶体位置的方法和装置。

背景技术

核医学设备是目前医学上常用的检测设备，例如，单电子发射计算机断层（SPECT）设备、正电子发射计算机断层（PET）设备等。核医学设备能够将含有放射性核素的药物在体内的分布形成图像，该图像可以反映人体代谢、组织功能和结构形态。

核医学设备中，最为核心的部件为探测器，该部件用于检测引入病患体内的放射性核素所发出的射线（例如γ射线）。如图1所示，常用的核医学设备探测器包括由多个晶体单元组成的晶体阵列和光电转换器。其中，晶体阵列用于检测病患体内释放出的射线光子（例如γ光子）并将其转换成可见光，光电转换器用于将可见光转换成电信号，依据该电信号可以确定射线光子在二维位置散点图上所形成的像素点的坐标，再基于该像素点的坐标就可以确定射线管子入射的晶体单元所在的位置，以便核医学设备进行整个图像的绘制。因此，为了实现基于位置散点图上像素点的坐标来确定射线光子入射的晶体单元位置，就需要在位置散点图上识别各个晶体单元所对应的坐标范围，即需要识别各晶体单元在位置散点图中所对应的位置。

由于受到光电转换器响应非线性、输出信号不确定性以及晶体条物理特性非均一不完全一致性、前端信号采集线路以及康普顿散射等因素的影响，由电信号确定出的像素点坐标与光子入射的晶体单元并非成线性的对应关系，这就导致位置散点图上晶体单元的坐标范围与晶体单元之间成非线性的对应关系。例如，假设位置散点图为m×m像素的图像，而晶体阵列为n×n的晶体单元组合，此时，按照线性的对应关系，图1中左下角第1个晶体单元在位置散点图上的位置坐标范围应为0<x<m/n、0<x<m/n，但实际上在位置散点图中该晶体单元的坐标范围却不是这样。

为了在位置散点图中识别这种非线性对应关系下各晶体单元对应的位置，现有技术中是先以射线源发出的大量均匀光子辐射探测器晶体单元阵列，然后通过探测器在每个光子入射所产生的电信号确定每个光子的像素点坐标，并基于大量光子的像素点坐标生成位置散点图，再基于位置散点图的灰度值分布对位置散点图进行分割来确定在位置散点图中各个晶体单元位置的边界，从而依据晶体单元位置边界识别出各个晶体单元位置的坐标范围。其中，位置散点图上每个像素点的灰度值表示该像素点坐标上入射光子的数量，即表示该像素点坐标上接收的辐射能量。如图2所示的是一幅二维晶体单元位置散点图的图像。

在基于灰度值对位置散点图进行分割时，需要先找出位置散点图上的灰度值的峰值点和谷值点，然后依据峰值点和谷值点来对图像进行分割，其中，峰值点为灰度极大值点，谷值点为灰度极小值点。现有技术中是先找出峰值点，再将相邻峰值点之间的距离中点作为谷值点。但是，由于各个晶体单元的物理特性不完全一致，相邻峰值点之间的距离中点通常都不是真实的谷值点。因此，由于现有技术并不能准确地确定谷值点，位置散点图的分割以及分割得到的晶体单元位置的边界就不准确，从而导致从位置散点图上识别出的各个晶体单元所对应的位置不够准确。

发明内容

本申请所要解决的技术问题是，提供一种识别探测器晶体单元位置的方法和装置，以解决按照现有技术中以相邻峰值点的距离中点作为谷值点对位置散点图进行分割所导致的晶体单元位置识别不准确的技术问题。

为解决上述技术问题，本申请提供了一种识别探测器晶体单元位置的方法。该方法包括：

获取基于辐射源辐射所述探测器所产生的电信号而生成的位置散点图；

计算所述位置散点图上的各个峰值点，并基于坐标轴方向上每两个相邻峰值确定一个初始谷值点；

基于各个所述峰值点和各个所述初始谷值点，对所述位置散点图进行初次分割，以便确定各个所述晶体单元在所述位置散点图上对应的初始边界；

分别对位置散点图中各行晶体单元初始边界内的图像区域进行垂直灰度投影，得到各行晶体单元的灰度投影曲线，并分别对位置散点图中各列晶体单元初始边界内的图像区域进行水平灰度投影，得到各列晶体单元的灰度投影曲线；

依据各行及各列晶体单元的灰度投影曲线上计算出的各个极小值点，确定各个目标谷值点；

基于各个所述峰值点和各个所述目标谷值点，对所述位置散点图进行再次分割，以便确定各个所述晶体单元在所述位置散点图上对应的目标边界。

可选的，所述计算所述位置散点图上的各个峰值点，包括：

根据所述位置散点图的灰度最大值和灰度最小值，生成表示不同阈值的灰度等高线；

根据所述灰度等高线计算位置散点图上的各个峰值点。

可选的，采用基于欧式距离的等权重泰森多边形Voronio图方法执行所述对所述位置散点图进行初次分隔；采用基于椭圆距离的不等权重Voronio图方法执行所述对所述位置散点图进行再次分割。

可选的，所述依据各行及各列晶体单元的灰度投影曲线上计算出的各个极小值点，确定各个目标谷值点，包括：

对于每一行及每一列晶体单元，

从该行或该列晶体单元的灰度投影曲线上选取多个当前参考点；

将各个所述当前参考点拟合成该行或该列晶体单元的灰度拟合曲线；

以所述灰度拟合曲线上的极小值点确定该行或该列晶体单元彼此间的各个目标谷值点。

可选的，所述方法还包括：

对所述位置散点图进行平滑增强，并以平滑增强所得的图像重新作为位置散点图，执行所述计算所述位置散点图上的各个峰值点。

此外，本申请还提供了一种识别探测器晶体单元位置的装置。该装置包括：

图像获取模块，用于获取基于辐射源辐射所述探测器所产生的电信号而生成的位置散点图；

峰值点计算模块，用于计算所述位置散点图上的各个峰值点；

初始谷值点确定模块，用于基于坐标轴方向上每两个相邻峰值点确定一个初始谷值点；

初次分割模块，用于基于各个所述峰值点和各个所述初始谷值点，对所述位置散点图进行初次分割，以便确定各个所述晶体单元在所述位置散点图上对应的初始边界；

灰度投影模块，用于分别对位置散点图中各行晶体单元初始边界内的图像区域进行垂直灰度投影，得到各行晶体单元的灰度投影曲线，并分别对位置散点图中各列晶体单元初始边界内的图像区域进行水平灰度投影，得到各列晶体单元的灰度投影曲线；

目标谷值点确定模块，用于依据各行及各列晶体单元的灰度投影曲线上计算出的各个极小值点，确定各个目标谷值点；

再次分割模块，用于基于各个所述峰值点和各个所述目标谷值点，对所述位置散点图进行再次分割，以便确定各个所述晶体单元在所述位置散点图上对应的目标边界。

可选的，所述峰值点计算模块包括：

等高线生成子模块，用于根据所述位置散点图的灰度最大值和灰度最小值，生成表示不同阈值的灰度等高线；

峰值点计算子模块，用于根据所述灰度等高线计算位置散点图上的各个峰值点。

可选的，所述装置还包括：

欧式距离子模块，用于采用基于欧式距离的等权重泰森多边形Voronio图方法执行所述对所述位置散点图进行初次分隔；

椭圆距离子模块，用于采用基于椭圆距离的不等权重Voronio图方法执行所述对所述位置散点图进行再次分割。

可选的，所述目标谷值点确定模块包括：

各行列触发子模块，用于对于每一行及每一列晶体单元分别触发参考点选取子模块；

参考点选取子模块，用于从该行或该列晶体单元的灰度投影曲线上选取多个当前参考点；

拟合曲线子模块，用于将各个所述当前参考点拟合成该行或该列晶体单元的灰度拟合曲线；

目标谷值点确定子模块，用于以所述灰度拟合曲线上的极小值点确定该行或该列当前晶体单元彼此间的各个目标谷值点。

可选的，所述装置还包括：

去噪模块，用于对所述位置散点图进行平滑增强，并以平滑增强所得的图像重新作为位置散点图，执行所述计算所述位置散点图上的各个峰值点。

与现有技术相比，本发明具有以下优点：

本申请实施例的技术方案，先基于相邻峰值点确定出的初始谷值点对位置散点图进行初次分割得到各个晶体单元的初始边界，然后基于初始边界对各行晶体单元和各列晶体单元进行灰度投影得到各行及各列晶体单元的灰度投影曲线，再以各灰度投影曲线上极小值点作为目标谷值点，对位置散点图进行再次分割得到各个晶体单元的目标边界，以便以目标边界识别各个晶体单元在位置散点图上的位置。通过本申请实施例，由于目标谷值点是计算其所在行或列的晶体单元的灰度投影曲线上的极小值点而确定的，而并非是以相邻峰值点的距离中点而确定的，并且，最终识别晶体单元位置所依据的目标边界正是基于目标谷值点对位置散点图进行分割而确定的，因此，可以实现基于图像的灰度极小值点来确定用于最终分割出晶体单元边界的谷值点，避免各个晶体单元物理特性不一致对谷值点准确性的影响，从而使得位置散点图中各个晶体单元位置的识别更加准确。

附图说明

为了更清楚地说明本申请实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本申请中记载的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为一种核医学设备探测器的组成结构示意图；

图2为一幅二维晶体单元位置散点图的图像示意图；

图3为本申请中识别探测器晶体单元位置的方法一实施例的流程图；

图4为本申请实施例中计算位置散点图上各个峰值点一实施方式的流程图；

图5为本申请实施例中在位置散点图上产生等高线的示意图；

图6为本申请实施例中所形成的灰度投影曲线的示意图；

图7为本申请实施例中确定目标谷值点一实施方式的流程图；

图8为本申请实施例中拟合曲线与灰度投影曲线的对比示意图；

图9为本申请实施例中标识出目标谷值点的位置散点图的示意图；

图10为本申请中识别探测器晶体单元位置的装置一实施例的结构图；

图11为本申请实施例中峰值点计算模块1002一实施方式的结构图；

图12为本申请中识别探测器晶体单元位置的装置又一实施例的结构图；

图13为本申请实施例中目标谷值点确定模块1006一实施方式的结构图；

图14为本申请中识别探测器晶体单元位置的装置再一实施例的结构图。

具体实施方式

为了使本技术领域的人员更好地理解本申请方案，下面将结合本申请实施例中的附图，对本申请实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅是本申请一部分实施例，而不是全部的实施例。基于本申请中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本申请保护的范围。

在本文中，所涉及的“晶体单元位置”，表示的是探测器上各个晶体单元在二维的位置散点图中对应的坐标范围，其中，位置散点图表示的是光子入射位置的图像，在本申请实施例中特指由辐射源辐射出的大量均匀光子入射到探测器上所形成的图像。此外，本文中的探测器，表示的是位于核医学设备中、用于形成医疗放射图像的探测器，如SPECT、PET等探测器。

发明人经过研究发现，现有技术中之所以不能准确地在位置散点图上识别出各个晶体单元对应的位置，是由于位置散点图分割时采用的谷值点不准确。而谷值点之所以不准确，是由于现有技术中并未从晶体单元整体的灰度考虑谷值点的坐标。例如，一种现有技术是将相邻峰值点的距离中点作为谷值点，这种方式并非基于位置散点图上的灰度分布，因此无法避免各个晶体单元不同的物理特性所导致的相邻峰值点之间灰度分布不均匀，显然所确定出的谷值点是不准确的。

基于上述发明的发现，本申请的主要思想是：首先，以基于相邻峰值点确定出的初始谷值点对位置散点图进行初次分割，粗略的确定各个晶体单元的边界位置，即初始边界；然后，以各晶体单元的初始边界所确定的各行晶体单元及各列晶体单元在位置散点图上的位置，对各行晶体单元和各列晶体单元进行灰度投影；再以各行及各列晶体单元的灰度投影曲线上的极小值点作为目标谷值点对位置散点图进行再次分割，以便更精确地确定各个晶体单元的边界位置，即目标边界。由此可见，由于最终得到的各个晶体单元目标边界是基于目标谷值点确定的，而目标谷值点是计算其所在行或列的晶体单元的灰度投影曲线上的极小值点而确定的，从而实现了基于晶体单元整体的灰度来确定用于最终分割出晶体单元边界的谷值点，避免了各个晶体单元物理特性不一致对谷值点准确性的影响，使得位置散点图中各个晶体单元位置的识别更加准确。

下面结合附图，通过实施例来详细说明本申请中识别探测器晶体单元位置的方法和装置的具体实现方式。

参见图3，示出了本申请中识别探测器晶体单元位置的方法一实施例的流程图。在本实施例中，例如可以包括以下步骤：

S301、获取基于辐射源辐射所述探测器所产生的电信号而生成的位置散点图。

其中，辐射源可以依据探测器的闪烁晶体来选择，其所辐射出的射线也需要是探测器实际所要探测的射线。例如，对于PET探测器来说，其闪烁晶体为BGO（分子式为Bi₄Ge₃O₁₂）晶体，该晶体一般地会被分割成10×10的晶体阵列，则辐射源可以选用γ射线泛源，射线源可以是Na-22、Ge-68或Cs-137等γ射线，能量范围可以为0.5Mev到10Mev之间。

当一个光子入射到某一晶体单元中时，晶体单元被激发出可见光，该可见光会按照一定的比例到达各个光电转换器产生电信号，依据各个光电转换器的电信号，就可以计算出该光子入射所形成的像素点在二维位置散点图上的坐标。例如，对于图1所示的具有4个光电倍增管（Photomultiplier Tube，简称PMT，一种光电转换器）的探测器来说，光电倍增管A、B、C、D分别产生的电信号电流为I_A、I_B、I_C、I_D，则此时表示光子入射的像素点坐标（x，y）可采用重心法计算，即：

$x=\frac{I_A+I_B-I_C-I_D}{I_A+I_B+I_C+I_D},y=\frac{I_A-I_B+I_C-I_C}{I_A+I_B+I_C+I_D}.$

当某一光子入射到探测器上时形成了一个像素点，则将该像素点的计数值加1，而该像素点的灰度值即为该计数值，由此可见，位置散点图上像素点的灰度值表示的是入射形成该像素点的光子数量。需要说明的是，由于本实施例的目的是在位置散点图上识别各个晶体单元的位置，即在位置散点图上找到各个晶体单元的边界，因此，在形成位置散点图时，需要辐射源利用大量均匀的光子辐射探测器的晶体阵列，这样，位置散点图上像素点的灰度分布才能体现出晶体单元边界，其中，在晶体单元边界处的灰度相对于晶体单元内部的灰度要小。

在本实施例中，S301形成位置散点图之后，可以直接进入S302从该位置散点图上确定各个晶体单元的边界。但由于大量光子形成的位置散点图上像素点灰度分布是离散、毛糙而并非连续、平滑的，这会导致后续位置散点图上的图像处理不精确。为了后续对位置散点图的图像处理更加精确，本实施例在S301执行完成以后，还可以先对所述位置散点图进行平滑增强，然后再以平滑增强所得的图像重新作为位置散点图，进入执行S302。具体地，可以采用高斯函数差分（Difference of Gaussian，简称DoG）滤波算子对所述位置散点图进行平滑去噪。更具体地，可以采用二阶微分算法的DoG滤波算子来对位置散点图进行高斯平滑模糊处理、检测图像中的团块处理。可以理解的是，之所以选用DoG滤波算子，是因为本实施例后续需要对位置散点图进行边缘分割，而DoG滤波算子即是主要用于边缘特征提取及模式识别中的分割预处理。其中，DoG滤波算子的主要参数为两个高斯函数的方差。由于方差设计可以对不同的图像特征情况下有不同的表现，因此，可以将S301中得到的位置散点图在不同参数下的高斯滤波结果相减而得到DoG图，并以该DoG图作为去噪后的位置散点图。具体地，DoG图像的计算方式可以采用如下公式：

$\begin{array}{c}D(x,y,\mathrm{\&sigma;})=(G(x,y,\mathrm{k\&sigma;})-G(x,y,\mathrm{\&sigma;}))*I(x,y)\\ =G(x,y,\mathrm{k\&sigma;})*I(x,y)-G(x,y,\mathrm{\&sigma;})*I(x,y)\end{array};$

其中，D(x,y,σ)为去噪后的位置散点图的灰度分布函数，I(x,y)为去噪前的位置散点图（即S301中得到的位置散点图）的灰度分布函数，k为一常数，G(x,y,σ)为方差为σ的高斯分布函数，具体可以采用如下公式来定义：

$G(x,y,\mathrm{\&sigma;})=\frac{1}{{2\mathrm{\&pi;\&sigma;}}^2}{e}^{-(x+y)/2{\mathrm{\&sigma;}}^2}.$

S302、计算所述位置散点图上的各个峰值点，并基于坐标轴方向上每两个相邻峰值点确定一个初始谷值点。

其中，由于在位置散点图上各个晶体单元位置中灰度是由四周向中心逐渐增大的，因此，可以基于轮廓线对位置散点图进行阈值分割，以便在位置散点图上产生表示灰度值高低不同的轮廓线（即灰度等高线），再基于灰度等高线的排列即可计算到各个峰值点。具体地，例如可以采用图4所示的实施方式来计算位置散点图上的各个峰值点：

S401、根据所述位置散点图的灰度最大值和灰度最小值，生成表示不同阈值的灰度等高线。

其中，对于灰度等高线的生成，一种可能的实施方式可以是，先依据所述位置散点图的灰度最大值和灰度最小值，计算灰度等高线之间的阈值差值，再以所述阈值差值在所述位置散点图上生成表示不同灰度阈值的灰度等高线。

需要说明的是，在计算灰度等高线之前，需要先对位置散点图进行二值化处理，二值化所需的阈值即为灰度等高线处的灰度值。其中，对应于不同灰度值的灰度等高线，二值化所需的阈值也不同，例如，可以各个灰度等高线的阈值可以成一个等差数列，而相邻灰度等高线的阈值之间的阈值差值可以采用下式来计算：

$s_0=\frac{T_{\max }-T_{\min }}{n};$

其中，T_max表示位置散点图上的灰度最大值，T_min表示位置散点图上的灰度最小值，n表示位置散点图上的灰度等高线层数，s₀表示各相邻等高线之间的阈值差值。

在利用灰度最大值、灰度最小值以及阈值差值来计算灰度等高线时，首先确定灰度值最小的一层灰度等高线的阈值作为阈值初始值，再以阈值差值计算出其他各层灰度等高线的阈值，再利用各个灰度等高线的阈值对位置散点图做二值化处理，即可在位置散点图上产生灰度等高线。如图5，示意性地示出了一位置散点图上所产生的灰度等高线。可以理解的是，阈值初始值可以是采用一个预设值，或者，也可以是直接采用阈值差值来作为阈值初始值。

S402、根据所述灰度等高线计算位置散点图上的各个峰值点。

具体地，可以根据所述灰度等高线的排列，采用重心法计算所述位置散点图上的各个峰值点。如图5所示的“+”，即为各个峰值点。

接着返回图3。

需要说明的是，确定初始谷值点并非是为了最终确定各晶体单元的目标边界，而是为了更准确地计算目标谷值点。其中，初始谷值点为基于相邻峰值点而确定出的，在坐标轴方向上的每两个相邻峰值点可以确定出一个初始谷值点，即同一行上两个相邻的峰值点之间产生一个初始谷值点，同一列上两个相邻的峰值点之间产生一个初始谷值点，而不同行且不同列的相邻峰值点之间不产生初始谷值点；如图5所示的位置散点图上，图像对角线方向上相邻的峰值点之间是不确定初始谷值点的。可以理解的是，初始谷值点可以采用现有技术中任意一种谷值点的确定方式，本申请实施例对此不做限定。例如，初始谷值点可以是相邻峰值点之间的距离中点或者某个等分点，如三等分点、四等分点；又如，初始谷值点也可以是在相邻峰值点之间连线上的搜索出的最小值点。

S303、基于各个所述峰值点和各个所述初始谷值点，对所述位置散点图进行初次分割，以便确定各个所述晶体单元在所述位置散点图上对应的初始边界。

其中，基于峰值点和初始谷值点对位置散点图的初次分割，是粗略地对位置散点图进行边缘分割而得到各晶体单元粗略的初始边界。例如，可以采用泰森多边形Voronio图方法，执行所述对所述位置散点图进行初次分割。Voronio图是一种点的分布图案，通常是由一组由连续两邻点直线的垂直平分线组成的连续多边形。

可以理解的是，由于仅需要粗略地确定各晶体单元的边界来进一步更准确地计算目标谷值点，因此，在分割时可以为节省图像处理所需的时间而采用相对不太准确的分割方法。例如，在采用Voronio图方法时，可以选用基于欧式距离的各点等权重Voronio图。其中，每个点的欧式距离可以采用下式计算：

${\left|\right|d\left|\right|}_i=\sqrt{x^2+y^2},$ 其中，d≡(x,y)。

需要说明的是，对位置散点图进行分割，实际上即是使用图像分割的方法对位置散点图的灰度分布函数进行处理，从而得到表示晶体单元边界的方程。由于此处属于图像处理中常用的现有技术，本申请在此不再赘述。

S304、分别对位置散点图中各行晶体单元初始边界内的图像区域进行垂直灰度投影，得到各行晶体单元的灰度投影曲线，并分别对位置散点图中各列晶体单元初始边界内的图像区域进行水平灰度投影，得到各列晶体单元的灰度投影曲线。

其中，一行晶体单元即为位置散点图中纵坐标范围大体相当的晶体单元，一列晶体单元即为位置散点图上横坐标范围大体相当的晶体单元，其中，这里对于晶体单元对应的坐标范围可以是以初次分割所得到的初始边界来确定的。

在本实施例中，垂直灰度投影是将每一行晶体单元范围内的像素点灰度向横坐标上投影，从而形成一条表示灰度与横坐标间对应关系的曲线，作为每一行晶体单元各自的灰度投影曲线。其中，各行晶体单元的灰度投影曲线上，每一个横坐标上的灰度值，是由该行晶体单元内该横坐标上的所有像素点的灰度值进行叠加而得到的。相似地，水平灰度投影是将每一列晶体单元范围内的像素点灰度向纵坐标上投影，从而形成一条表示灰度与纵坐标间对应关系的曲线，作为每一列晶体单元各自的灰度投影曲线。其中，各列晶体单元的灰度投影曲线上，每一个纵坐标上的灰度值，是由该列晶体单元内该纵坐标上的所有像素点的灰度值进行叠加而得到的。如图6，示意性地示出了一条灰度投影曲线的形状。

可以理解的是，由于灰度投影是基于各行及各列晶体单元的坐标范围来进行的，所以，需要基于各个晶体单元的初始边界来确定各行及各列晶体单元的坐标范围。此时，一种可能的确定方式可以是，对于各行或各列的晶体单元，先确定该行或该列所包括的晶体单元，然后将这些晶体单元的初始边界内的坐标范围合并成该行或该列晶体单元的坐标范围。但由于同一行或同一列中各个晶体单元的初始边界一般并不均匀，一行晶体单元的每个横坐标上对应的纵坐标范围都不同，一列晶体单元的每个纵坐标上对应的横坐标范围也都不同，这显然会增加灰度投影时的系统处理时间。为了简化灰度投影过程，缩短系统处理时间，另一种可能的确定方式可以是，对于各行或各列的晶体单元，先确定该行或该列所包括的晶体单元，然后以该行晶体单元的初始边界为该行粗略地确定一个固定不变的纵坐标范围，或者以该列晶体单元的初始边界为该列粗略地确定一个固定不变的横坐标范围。由于灰度投影之后还需要再次对位置散点图进行更加精确地再次分割，所以在灰度投影时各行及各列晶体单元的坐标范围可以不必十分精确。

S305、依据各行及各列晶体单元的灰度投影曲线上计算出的各个极小值点，确定各个目标谷值点。

其中，对于目标谷值点的确定，是在各行及各列晶体单元中分别进行的；对于每一行和每一列晶体单元来说，目标谷值点即位于其中相邻两个晶体单元之间。

需要说明的是，在基于灰度投影曲线来确定目标谷值点时，一种可能的确定方式是采用最小值搜索的方式来将搜索到的极小值直接作为目标谷值点。但是，由于灰度投影曲线上往往会在真正的谷值点上叠加噪声，并且前续步骤中也可能存在误差而导致灰度投影曲线上具有毛刺，这样直接将极小值作为目标谷值点就会导致目标谷值点并非真正的谷值点。为了使目标谷值点更加准确，另一种可能的确定方式是先对灰度投影曲线采样并拟合，以拟合的方式来消除灰度投影曲线上可能带有的毛刺和噪声，然后再以拟合曲线上的极小值点作为目标谷值点，这样就可以避免毛刺和噪声对目标谷值点准确性的影响。具体地，对于每一行及每一列晶体单元，例如可以采用图7所示的方式来确定该行或该列上的目标谷值点：

S701、从该行或该列晶体单元的灰度投影曲线上选取多个当前参考点。

一般地，各行的灰度投影曲线所选取的当前参考点可以都是相同横坐标的参考点。其中，当前参考点的个数可以是以位置散点图上像素点的分布来确定的，例如，假设位置散点图上每行具有m个像素点，则可以在各行灰度投影曲线上均匀地选取m个当前参考点，假设位置散点图上每列也具有m个像素点，则可以在各列灰度投影曲线上均匀地选取m个当前参考点。此外，在选取当前参考点时，还可以将位置散点图上的坏点删除，只保留没有问题的像素点，这样可以避免坏点给目标谷值点确定带来偏差。可以理解的是，当前参考点的数量越多，所拟合出曲线就越接近于真实的灰度投影曲线。

S702、将各个所述当前参考点拟合成该行或该列晶体单元的灰度拟合曲线。

其中，灰度拟合曲线的拟合手段，可以采用任意一种用于拟合峰谷相间特性曲线的拟合手段。例如，灰度拟合曲线的拟合手段，可以为傅里叶级数、多项式函数、抛物线函数等等。进一步而言，由于灰度投影曲线峰谷相间的弦波特性，本实施方式可以优选傅里叶级数来对各行、各列的当前参考点进行拟合。其中，所选用的傅里叶级数的级数越高，则灰度拟合曲线越精确，但花费在拟合上的系统处理时间也会相应越长。综合考虑拟合精度和系统处理时间，可以选用三级傅里叶级数来进行拟合。例如，对于各行的灰度投影曲线进行拟合时，三级傅里叶级数拟合的曲线可以表示如下：

$\begin{array}{c}h=a_0+a_1*\cos (x*w)+b_1*\sin (x*w)+a_2*\cos (x*2w)+b_2*\sin (x*2w)\\ +a_3*\cos (x*3w)+b_3*\sin (x*3w)\end{array};$

其中，h表示拟合曲线的灰度值，x表示横坐标值，a0～a3、b0～b3以及w均是需要拟合来确定的未知参量。可以理解的是，对于各列的灰度投影曲线进行拟合时，可以将上述傅里叶级数中表示横坐标值的x换成表示纵坐标值的y。

S703、以所述灰度拟合曲线上的极小值点确定该行或该列晶体单元彼此间的各个目标谷值点。

如图8，示意性地示出了拟合曲线与灰度投影曲线的对比，其中，拟合曲线的极小值点即为一个目标谷值点。如图9，示意性地示出了一幅位置散点图上各行晶体单元之间的目标谷值点（图9中的三角块）。

接着返回图3。S305执行完成以后，进入执行S306。

S306、基于各个所述峰值点和各个所述目标谷值点，对所述位置散点图进行再次分割，以便确定各个所述晶体单元在所述位置散点图上对应的目标边界。

其中，基于峰值点和目标谷值点对位置散点图的再次分割，是相对于初次分割更加精确地对位置散点图进行边缘分割而得到各晶体单元精确的目标边界。例如，与初次分割相似的是，也可以采用Voronio图方法来执行所述对所述位置散点图进行再次分割。但为了更准确地确定各晶体单元的目标边界，与初次分割不同的是，再次分割时可以选用基于椭圆距离的不等权重Voronio图方法执行所述对所述位置散点图进行再次分割。其中，每个点的团员距离可以采用下式计算：

${\left|\right|d\left|\right|}_i=\sqrt{\frac{x^2}{a_i^2}+\frac{y^2}{b_i^2}};$

其中，a_i、b_i分别为团员的长轴和短轴，其数值可由相邻谷值点来确定。

在S306执行完成以后，由于已经确定了各个晶体单元的目标边界，进而就可以在位置散点图上确定出各个晶体单元对应的坐标范围，也即可以实现对探测器晶体单元位置的识别。此时，可以建立一个表示位置散点图上的坐标范围与晶体单元之间对应关系的探测器晶体编码表，这样，当某个光子入射到探测器的晶体上而产生电信号时，可以依据由电信号计算得到的坐标，对照探测器晶体编码表查找该坐标位于哪一晶体单元的坐标范围，从而就可以确定出该光子是从哪一个晶体单元入射的。

通过本实施例的技术方案，由于目标谷值点是依据其所在行或列的晶体单元的灰度投影曲线上的极小值点而确定的，而并非是以相邻峰值点的距离中点而确定的，并且，最终识别晶体单元位置所依据的目标边界正是基于目标谷值点对位置散点图进行分割而确定的，因此，可以实现基于图像的灰度极小值点来确定用于最终分割出晶体单元边界的谷值点，避免各个晶体单元物理特性不一致对谷值点准确性的影响，从而使得位置散点图中各个晶体单元位置的识别更加准确。

对应于方法实施例，本申请实施例还提供了一种识别探测器晶体单元位置的装置的实现方式。

参见图10，示出了本申请中识别探测器晶体单元位置的装置一实施例的结构图。在本实施例中，所述装置例如可以包括：

图像获取模块1001，用于获取基于辐射源辐射所述探测器所产生的电信号而生成的位置散点图；

峰值点计算模块1002，用于计算所述位置散点图上的各个峰值点；

初始谷值点确定模块1003，用于基于坐标轴方向上每两个相邻峰值点确定一个初始谷值点；

初次分割模块1004，用于基于各个所述峰值点和各个所述初始谷值点，对所述位置散点图进行初次分割，以便确定各个所述晶体单元在所述位置散点图上对应的初始边界；

灰度投影模块1005，用于分别对位置散点图中各行晶体单元初始边界内的图像区域进行垂直灰度投影，得到各行晶体单元的灰度投影曲线，并分别对位置散点图中各列晶体单元初始边界内的图像区域进行水平灰度投影，得到各列晶体单元的灰度投影曲线；

目标谷值点确定模块1006，用于依据各行及各列晶体单元的灰度投影曲线上计算出的各个极小值点，确定各个目标谷值点；

再次分割模块1007，用于基于各个所述峰值点和各个所述目标谷值点，对所述位置散点图进行再次分割，以便确定各个所述晶体单元在所述位置散点图上对应的目标边界。

参见图11，示出了本申请实施例中峰值点计算模块1002一实施方式的结构图。其中，所述峰值点计算模块1002例如可以包括：

等高线生成子模块1101，用于根据所述位置散点图的灰度最大值和灰度最小值，生成表示不同阈值的灰度等高线；

峰值点计算子模块1102，用于根据所述灰度等高线计算位置散点图上的各个峰值点。

参见图12，示出了本申请中识别探测器晶体单元位置的装置又一实施例的结构图。在本实施例中，除了图10所示的所有结构外，所述装置例如还可以包括：

欧式距离分割模块1201，用于采用基于欧式距离的等权重泰森多边形Voronio图方法执行所述对所述位置散点图进行初次分隔；

椭圆距离分割模块1202，用于采用基于椭圆距离的不等权重Voronio图方法执行所述对所述位置散点图进行再次分割。

参见图13，示出了本申请实施例中目标谷值点确定模块1006一实施方式的结构图。其中，所述目标谷值点确定模块1006例如可以包括：

各行列触发子模块1301，用于对于每一行及每一列晶体单元分别触发参考点选取子模块1302；

参考点选取子模块1302，用于从该行或该列晶体单元的灰度投影曲线上选取多个当前参考点；

拟合曲线子模块1303，用于将各个所述当前参考点拟合成该行或该列晶体单元的灰度拟合曲线；

目标谷值点确定子模块1304，用于以所述灰度拟合曲线上的极小值点确定该行或该列当前晶体单元彼此间的各个目标谷值点。

参见图14，示出了本申请中识别探测器晶体单元位置的装置再一实施例的结构图。在本实施例中，除了图10所示的所有结构外，所述装置例如还可以包括：

去噪模块1401，用于对所述位置散点图进行平滑增强，并以平滑增强所得的图像重新作为位置散点图，执行所述计算所述位置散点图上的各个峰值点。

通过本申请中装置实施例的技术方案，由于目标谷值点是依据其所在行或列的晶体单元的灰度投影曲线上的极小值点而确定的，而并非是以相邻峰值点的距离中点而确定的，并且，最终识别晶体单元位置所依据的目标边界正是基于目标谷值点对位置散点图进行分割而确定的，因此，可以实现基于图像的灰度极小值点来确定用于最终分割出晶体单元边界的谷值点，避免各个晶体单元物理特性不一致对谷值点准确性的影响，从而使得位置散点图中各个晶体单元位置的识别更加准确。

需要说明的是，在本文中，诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括所述要素的过程、方法、物品或者设备中还存在另外的相同要素。

对于装置实施例而言，由于其基本对应于方法实施例，所以相关之处参见方法实施例的部分说明即可。以上所描述的装置实施例仅仅是示意性的，其中所述作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部模块来实现本实施例方案的目的。本领域普通技术人员在不付出创造性劳动的情况下，即可以理解并实施。

以上所述仅是本申请的具体实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本申请原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本申请的保护范围。

Claims (10)

1.一种识别探测器晶体单元位置的方法，其特征在于，包括：

获取基于辐射源辐射所述探测器所产生的电信号而生成的位置散点图；位置散点图表示光子入射位置的图像；

计算所述位置散点图上的各个峰值点，并基于坐标轴方向上每两个相邻峰值确定一个初始谷值点；

基于各个所述峰值点和各个所述初始谷值点，对所述位置散点图进行初次分割，以便确定各个所述晶体单元在所述位置散点图上对应的初始边界；

分别对位置散点图中各行晶体单元初始边界内的图像区域进行垂直灰度投影，得到各行晶体单元的灰度投影曲线，并分别对位置散点图中各列晶体单元初始边界内的图像区域进行水平灰度投影，得到各列晶体单元的灰度投影曲线；

依据各行及各列晶体单元的灰度投影曲线上计算出的各个极小值点，确定各个目标谷值点；

基于各个所述峰值点和各个所述目标谷值点，对所述位置散点图进行再次分割，以便确定各个所述晶体单元在所述位置散点图上对应的目标边界。

2.根据权利要求1所述的方法，其特征在于，所述计算所述位置散点图上的各个峰值点，包括：

根据所述位置散点图的灰度最大值和灰度最小值，生成表示不同阈值的灰度等高线；

根据所述灰度等高线计算位置散点图上的各个峰值点。

3.根据权利要求1所述的方法，其特征在于，采用基于欧式距离的等权重泰森多边形Voronio图方法执行所述对所述位置散点图进行初次分隔；

采用基于椭圆距离的不等权重Voronio图方法执行所述对所述位置散点图进行再次分割。

4.根据权利要求1所述的方法，其特征在于，所述依据各行及各列晶体单元的灰度投影曲线上计算出的各个极小值点，确定各个目标谷值点，包括：

对于每一行及每一列晶体单元，

从该行或该列晶体单元的灰度投影曲线上选取多个当前参考点；

将各个所述当前参考点拟合成该行或该列晶体单元的灰度拟合曲线；

以所述灰度拟合曲线上的极小值点确定该行或该列晶体单元彼此间的各个目标谷值点。

5.根据权利要求1所述的方法，其特征在于，还包括：

对所述位置散点图进行平滑增强，并以平滑增强所得的图像重新作为位置散点图，执行所述计算所述位置散点图上的各个峰值点。

6.一种识别探测器晶体单元位置的装置，其特征在于，包括：

图像获取模块，用于获取基于辐射源辐射所述探测器所产生的电信号而生成的位置散点图；位置散点图表示光子入射位置的图像；

峰值点计算模块，用于计算所述位置散点图上的各个峰值点；

初始谷值点确定模块，用于基于坐标轴方向上每两个相邻峰值点确定一个初始谷值点；

初次分割模块，用于基于各个所述峰值点和各个所述初始谷值点，对所述位置散点图进行初次分割，以便确定各个所述晶体单元在所述位置散点图上对应的初始边界；

灰度投影模块，用于分别对位置散点图中各行晶体单元初始边界内的图像区域进行垂直灰度投影，得到各行晶体单元的灰度投影曲线，并分别对位置散点图中各列晶体单元初始边界内的图像区域进行水平灰度投影，得到各列晶体单元的灰度投影曲线；

目标谷值点确定模块，用于依据各行及各列晶体单元的灰度投影曲线上计算出的各个极小值点，确定各个目标谷值点；

再次分割模块，用于基于各个所述峰值点和各个所述目标谷值点，对所述位置散点图进行再次分割，以便确定各个所述晶体单元在所述位置散点图上对应的目标边界。

7.根据权利要求6所述的装置，其特征在于，所述峰值点计算模块包括：

等高线生成子模块，用于根据所述位置散点图的灰度最大值和灰度最小值，生成表示不同阈值的灰度等高线；

峰值点计算子模块，用于根据所述灰度等高线计算位置散点图上的各个峰值点。

8.根据权利要求6所述的装置，其特征在于，还包括：

欧式距离子模块，用于采用基于欧式距离的等权重泰森多边形Voronio图方法执行所述对所述位置散点图进行初次分隔；

椭圆距离子模块，用于采用基于椭圆距离的不等权重Voronio图方法执行所述对所述位置散点图进行再次分割。

9.根据权利要求6所述的装置，其特征在于，所述目标谷值点确定模块包括：

各行列触发子模块，用于对于每一行及每一列晶体单元分别触发参考点选取子模块；

参考点选取子模块，用于从该行或该列晶体单元的灰度投影曲线上选取多个当前参考点；

拟合曲线子模块，用于将各个所述当前参考点拟合成该行或该列晶体单元的灰度拟合曲线；

目标谷值点确定子模块，用于以所述灰度拟合曲线上的极小值点确定该行或该列当前晶体单元彼此间的各个目标谷值点。

10.根据权利要求6所述的装置，其特征在于，还包括：

去噪模块，用于对所述位置散点图进行平滑增强，并以平滑增强所得的图像重新作为位置散点图，执行所述计算所述位置散点图上的各个峰值点。