CN104656119A - 一种闪烁脉冲信息复原的方法及系统 - Google Patents
一种闪烁脉冲信息复原的方法及系统 Download PDFInfo
- Publication number
- CN104656119A CN104656119A CN201310585902.1A CN201310585902A CN104656119A CN 104656119 A CN104656119 A CN 104656119A CN 201310585902 A CN201310585902 A CN 201310585902A CN 104656119 A CN104656119 A CN 104656119A
- Authority
- CN
- China
- Prior art keywords
- pulse
- module
- energy value
- energy
- scintillation
- 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
Links
- 238000000034 method Methods 0.000 title claims abstract description 85
- 238000009826 distribution Methods 0.000 claims description 20
- 238000001228 spectrum Methods 0.000 claims description 18
- 238000012937 correction Methods 0.000 claims description 10
- 238000011084 recovery Methods 0.000 claims description 10
- 239000013078 crystal Substances 0.000 claims description 9
- 238000012935 Averaging Methods 0.000 claims description 6
- 239000012634 fragment Substances 0.000 claims description 6
- 230000000149 penetrating effect Effects 0.000 claims description 6
- 239000007787 solid Substances 0.000 claims description 6
- 238000004458 analytical method Methods 0.000 claims description 5
- 230000005855 radiation Effects 0.000 claims description 4
- 238000012549 training Methods 0.000 claims description 4
- 238000012850 discrimination method Methods 0.000 claims description 3
- 239000006185 dispersion Substances 0.000 claims description 3
- 238000010845 search algorithm Methods 0.000 claims description 3
- 238000004364 calculation method Methods 0.000 abstract description 10
- 238000013139 quantization Methods 0.000 abstract description 9
- 238000010586 diagram Methods 0.000 description 38
- 230000010354 integration Effects 0.000 description 21
- 230000008569 process Effects 0.000 description 8
- 238000005070 sampling Methods 0.000 description 8
- 239000000523 sample Substances 0.000 description 7
- 230000008878 coupling Effects 0.000 description 6
- 238000010168 coupling process Methods 0.000 description 6
- 238000005859 coupling reaction Methods 0.000 description 6
- 238000001514 detection method Methods 0.000 description 4
- 238000007493 shaping process Methods 0.000 description 4
- 238000013461 design Methods 0.000 description 3
- 230000005251 gamma ray Effects 0.000 description 3
- 239000002245 particle Substances 0.000 description 3
- 239000004809 Teflon Substances 0.000 description 2
- 229920006362 Teflon® Polymers 0.000 description 2
- 239000002390 adhesive tape Substances 0.000 description 2
- 238000005315 distribution function Methods 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 238000005259 measurement Methods 0.000 description 2
- 238000002600 positron emission tomography Methods 0.000 description 2
- 238000012545 processing Methods 0.000 description 2
- 230000000630 rising effect Effects 0.000 description 2
- 238000009825 accumulation Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 238000002059 diagnostic imaging Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000013213 extrapolation Methods 0.000 description 1
- JEIPFZHSYJVQDO-UHFFFAOYSA-N ferric oxide Chemical compound O=[Fe]O[Fe]=O JEIPFZHSYJVQDO-UHFFFAOYSA-N 0.000 description 1
- 238000001914 filtration Methods 0.000 description 1
- 238000012804 iterative process Methods 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 239000002184 metal Substances 0.000 description 1
- 238000013508 migration Methods 0.000 description 1
- 230000005012 migration Effects 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 238000010606 normalization Methods 0.000 description 1
- 238000004064 recycling Methods 0.000 description 1
- 230000001105 regulatory effect Effects 0.000 description 1
- 238000007789 sealing Methods 0.000 description 1
- 230000002123 temporal effect Effects 0.000 description 1
- 238000005303 weighing Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01T—MEASUREMENT OF NUCLEAR OR X-RADIATION
- G01T1/00—Measuring X-radiation, gamma radiation, corpuscular radiation, or cosmic radiation
- G01T1/16—Measuring radiation intensity
- G01T1/20—Measuring radiation intensity with scintillation detectors
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01T—MEASUREMENT OF NUCLEAR OR X-RADIATION
- G01T1/00—Measuring X-radiation, gamma radiation, corpuscular radiation, or cosmic radiation
- G01T1/16—Measuring radiation intensity
- G01T1/20—Measuring radiation intensity with scintillation detectors
- G01T1/208—Circuits specially adapted for scintillation detectors, e.g. for the photo-multiplier section
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01T—MEASUREMENT OF NUCLEAR OR X-RADIATION
- G01T1/00—Measuring X-radiation, gamma radiation, corpuscular radiation, or cosmic radiation
- G01T1/29—Measurement performed on radiation beams, e.g. position or section of the beam; Measurement of spatial distribution of radiation
- G01T1/2914—Measurement of spatial distribution of radiation
- G01T1/2985—In depth localisation, e.g. using positron emitters; Tomographic imaging (longitudinal and transverse section imaging; apparatus for radiation diagnosis sequentially in different planes, steroscopic radiation diagnosis)
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01T—MEASUREMENT OF NUCLEAR OR X-RADIATION
- G01T1/00—Measuring X-radiation, gamma radiation, corpuscular radiation, or cosmic radiation
- G01T1/36—Measuring spectral distribution of X-rays or of nuclear radiation spectrometry
- G01T1/362—Measuring spectral distribution of X-rays or of nuclear radiation spectrometry with scintillation detectors
Landscapes
- Physics & Mathematics (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- High Energy & Nuclear Physics (AREA)
- Molecular Biology (AREA)
- Measurement Of Radiation (AREA)
- Photometry And Measurement Of Optical Pulse Characteristics (AREA)
- Nuclear Medicine (AREA)
Abstract
一种闪烁脉冲信息复原的方法，其包括步骤：获取低计数下的非堆积符合单事件的闪烁脉冲数据库，再对各路单事件的闪烁脉冲数据库建立闪烁脉冲的噪声模型；根据闪烁脉冲的噪声模型计算给定能量值的后验概率对数值；通过计算反复地调用第二步骤，获得满足最大后验概率条件的能量值。一种闪烁脉冲信息复原的系统，其包括涨落模型模块、后验概率模块及能量值搜索模块。本发明闪烁脉冲信息复原的方法及系统，有效提高了系统能量计算精确度，特别适合于稀疏量化电平ADC数字式核仪器的能量计算。
Description
技术领域
本发明涉及数字信号处理、光电信号处理和辐射探测领域，尤其涉及一种闪烁脉冲信息复原的方法及系统。
背景技术
在正电子寿命谱仪、正电子角度-动量关联分析仪等核分析领域，双符合高能粒子甄别器等核探测领域，以及正电子发射断层成像（Positron EmissionTomography，以下简称PET）等医学影像领域，探测高能粒子的主要工作原理是把高能射线转化为电信号，再利用快电子学得到粒子事件的各种信息。在一些带有符合逻辑运算和时间标记的仪器应用中，对闪烁脉冲进行整形放慢会影响其时间性能，并增加事件脉冲发生堆积的概率，所以将闪烁脉冲直接进行数字化更加能够满足这种需求。
较为完整的对闪烁脉冲进行数字化是采用高带宽高采样率的模拟数字转换器（Analog to Digital Converter，以下简称ADC）直接对闪烁脉冲进行采样和量化。这种数字化方法的成本较高，不能满足多通道系统的实际需求。为了减少这种数字化方案的成本，有两种精简模式，一种是采用降低采样率的方式。由于采样频率较低，无整形的闪烁脉冲信号的高频部分不但不能被采样，反而会由于频域的混叠，影响到奈奎斯特频率附近的分量的精确度。为了缓解这种不足，会在闪烁脉冲输出和数字化之间添加一个整形环节。这个整形环节把闪烁脉冲高于奈奎斯特频率的分量进行衰减而保存了闪烁脉冲的低频分量，其作用是减少闪烁脉冲数字化在频率的混叠，但由于奈奎斯特频率本身的限制，超过该频率的分量始终是不能被采样到的。
除了在时间轴上对ADC进行精简以外，还有一种方式是对ADC在电压轴上进行精简。在电压轴上进行精简的ADC可以利用开环设计提高有效带宽，适合于高速信号处理。例如，若干个比较器/ADC单元的数字化方法，或者是开路时间交织的ADC都属于这个范畴。这种ADC的特点是高采样率高带宽，但量化精度有限。由于闪烁脉冲具有比较鲜明的先验知识，通过利用这些知识可以显著提高数字化后信号复原的精度和优化量化电平的设置。
因此，针对上述稀疏量化电平ADC，有必要提供一种新的烁脉冲信息复原的方法及系统，以解决现有技术中存在的问题。
发明内容
有鉴于此，本发明的目的在于提供一种闪烁脉冲信息复原的方法及系统，一方面，显著提高数字化后信号复原的精度；另一方面，优化量化电平的设置可以获取更多的事件信息。
为实现上述目的，本发明提供如下技术方案：
一种闪烁脉冲信息复原的方法，其具体包括步骤：
S1：获取低计数下的非堆积符合单事件的闪烁脉冲数据库，再对各路单事件的闪烁脉冲数据库建立闪烁脉冲的噪声模型，对各路非堆积符合单事件的闪烁脉冲数据库求平均脉冲，闪烁脉冲的形状信息由平均脉冲给定；
S2：根据闪烁脉冲的噪声模型计算给定能量值的后验概率对数值，具体包括步骤：
（2.1）载入一个闪烁脉冲片段S0，该闪烁脉冲片段的起点是向上过阈值Vl的时间点t0，终点是t0+Δt，其中Δt是一个大于闪烁晶体下降沿时间常数两倍的值，根据脉冲噪声模型，计算给定能量值的似然函数；
（2.2）将所有时间点取对数，再相加，得到一个单调性与后验概率相同的函数的具体值；
S3：通过计算反复地调用步骤S2，获得满足最大后验概率条件的能量值，具体包括步骤：
（3.1）将步骤S2作为一个模块反复调用，计算不同试探能量值的后验概率；
（3.2）线性搜索满足最大后验概率的能量值；
（3.3）计算得到能量值后，根据能量值校准2.1中的时间基点，再重复3.1和3.2。
优选的，在上述闪烁脉冲信息复原的方法中，所述步骤S1中闪烁脉冲数据库中的样本数在2000个以上。
优选的，在上述闪烁脉冲信息复原的方法中，所述步骤S1中闪烁脉冲的噪声模型包括特征：闪烁脉冲形状、闪烁脉冲高度、闪烁脉冲噪声水平。
优选的，在上述闪烁脉冲信息复原的方法中，所述步骤S1中求平均脉冲的具体步骤包括：
（1.1）降低射源的辐射剂量得到弱源，通过弱源或调整弱源对探测器的立体角，降低每个探测器捕获的高能光子数，每个探测器接收到的事件是一个泊松流，其平均计数率为
其中，mi和qi分别是弱源的剂量和弱源对探测器的立体角，i是弱源的序数，n是弱源的个数；
（1.2）按照恒定系数鉴别方法或前沿甄别方法对齐脉冲；
（1.3）对对齐后的脉冲取平均。
优选的，在上述闪烁脉冲信息复原的方法中，所述步骤S1中对齐脉冲后，同一个能量值的闪烁脉冲在一个时间点上能量值会有一个分散的分布，这个分布被定义为脉冲的噪声模型，噪声模型的参数留作计算给定能量值的后验概率值。
优选的，在上述闪烁脉冲信息复原的方法中，所述步骤S1还包括：将一能量区间选择的符合脉冲对认定为单事件脉冲对，将此数字化的电脉冲预先存储后，进行离线分析，通过对齐脉冲对，获得平均脉冲信号，将平均脉冲作为系统响应，将数字脉冲加和，得到单事件数据的能谱。
一种闪烁脉冲信息复原的系统，其包括涨落模型模块、后验概率模块及能量值搜索模块，其中，
所述涨落模型模块用于对预采集的低计数率数据进行训练，获得闪烁探测器系统的涨落模型，该涨落模型表征平均脉冲、噪声方差和噪声分布偏度，该噪声模型提供给后验概率模块；
所述后验概率模块用于后验概率的计算，获得给定能量值的后验概率对数值；
所述能量值搜索模块用于搜索满足最大后验概率条件的能量值。
优选的，在上述闪烁脉冲信息复原的系统中，所述涨落模型模块包括平均脉冲模块、脉冲噪声模块及概率密度分布模块，其中，
所述平均脉冲模块用于计算单事件的平均脉冲；
所述脉冲噪声模块用于计算单事件脉冲的噪声方差；
所述概率密度分布模块用于估计噪声的概率密度分布函数。
优选的，在上述闪烁脉冲信息复原的系统中，所述后验概率模块包括似然函数分布模块及时间点乘积模块，其中，
所述似然函数分布模块按照给定能量值和涨落模型模块输入的脉冲噪声模型计算每个点的似然函数值；
所述时间点乘积模块负责计算后验概率对数值。
优选的，在上述闪烁脉冲信息复原的系统中，所述能量值搜索模块包括试探能量计算模块、线性搜索模块及时间基点校正模块，其中，
所述试探能量计算模块通过调用后验概率模块计算出每个试探能量的后验概率对数值；
所述线性搜索模块通过线性搜索算法计算满足最大后验概率条件的能量值；
所述时间基点校正模块校正样本的时间基点，并重复调用试探能量计算模块和线性搜索模块。
从上述技术方案可以看出，本发明实施例通过在上述闪烁脉冲信息复原方法与系统中，在系统运行时可有效地计算出各个单事件的能量信息；增加了系统在能量窗内的计数率，提高了系统的能量分辨率，特别适合于离线环境下的闪烁事件处理。
从上述技术方案可以看出，本发明实施例通过计算给定系统的最大后验估计对应的闪烁事件参数，准确地估计了闪烁脉冲的脉冲高度或者能量信息。
与现有技术相比，本发明的有益效果是：
（1）较好的能量分辨率。
（2）给定能量窗内的计数率高。
（3）可避免量化电平较少时带来的量化误差。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的有关本发明的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1为本发明闪烁脉冲信息复原的系统结构图；
图2为本发明闪烁脉冲信息复原的方法流程图；
图3为一种典型的闪烁脉冲示意图；
图4为一个能量的闪烁脉冲集合示意图；
图5为图3的闪烁脉冲在四个阈值的ADC下的数字化结果；
图6为一个能量分辨率较差的能谱示意图，其能量计算方法为数字门控积分；
图7为一个能量分辨率较差的能谱示意图，其能量计算方法为本发明所提出的方法；
图8为一个能量分辨率较差的能谱示意图，其能量计算方法为最小二乘拟合算法；
图9为一个能量分辨率较好的能谱示意图，其能量计算方法为数字门控积分；
图10为一个能量分辨率较好的能谱示意图，其能量计算方法为本发明所提出的方法；
图11为一个能量分辨率较好的能谱示意图，其能量计算方法为最小二乘拟合算法；
图12为最小二乘拟合算法得到的能量值对数字门控积分得到的能量值的散点图，能量值为图4的能量值；
图13为本发明方法得到的能量值对数字门控积分得到的能量值的散点图，能量值为图4的能量值；
图14为最小二乘拟合算法得到的能量值对数字门控积分得到的能量值的散点图，能量值为图7的能量值；
图15为本发明方法得到的能量值对数字门控积分得到的能量值的散点图，能量值为图7的能量值；
图16为一个能量值一维解寻优搜索过程图；
图17为本发明一种典型系统的示意图；
图18为本发明另一种典型系统的示意图。
具体实施方式
本发明公开了一种闪烁脉冲信息复原方法与系统，该方法及系统能有效实现稀疏量化电平ADC采样下的脉冲能量信息复原，提高闪烁脉冲高度估计的准确度，增加系统在能量窗内的计数率。
本发明公开的闪烁脉冲信息复原的方法，其步骤包括：
（1）获取低计数下的非堆积符合单事件的脉冲数据库，再对各路单事件的脉冲数据库建立噪声模型。要求脉冲数据库中的样本数在2000个以上。脉冲的个数越多，噪声模型越精确；
（2）计算给定能量值的后验概率对数值；
（3）线性搜索满足最大后验概率条件的能量值。
本发明公开的闪烁脉冲信息复原的系统，包括噪声模型模块、后验概率模块、能量值搜索模块。其中，噪声模型模块用于在单事件脉冲集合中获得噪声模型，其他模块用于处理基于稀疏量化电平ADC的闪烁脉冲信号；后验概率模块用于对输入的一个能量值计算后验概率；能量值搜索模块用于搜索满足最大后验概率条件的能量值。
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行详细地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都属于本发明保护的范围。
本方法是根据贝叶斯原理，获得满足最大后验概率条件的能量值。
本发明提供的闪烁脉冲信息复原的方法首先利用预先获取好的数据库获得噪声模型。该噪声模型将作为对电信号的固有先验知识，计算稀疏量化的闪烁脉冲的后验概率。通过不断地调用后验概率模块，线性搜索满足最大后验概率条件的能量值。
如图2所示，本发明公开的闪烁脉冲信息复原的方法包括以下步骤：
S1：获取低计数下的非堆积符合单事件的闪烁脉冲数据库，再对各路单事件的闪烁脉冲数据库建立闪烁脉冲的噪声模型，对各路非堆积符合单事件的闪烁脉冲数据库求平均脉冲，闪烁脉冲的形状信息由平均脉冲给定。要求脉冲数据库中的样本数大于2000。脉冲数越多，其统计噪声越小。
所述步骤S1中闪烁脉冲的噪声模型包括特征：闪烁脉冲形状、闪烁脉冲高度、闪烁脉冲噪声水平。这几个特征是计算单个时间样本对能量的后验概率的参数，没有这几个值就不能计算对应于能量的后验概率参数。这几个参数决定于系统。建立噪声模型是所属领域内技术人员都能实现的，但各自的方法不同。一般地，噪声模型可以通过仿真或实验测量得到。仿真时，根据假设给出每个脉冲在不同时间点的电压值分布。这个分布就是噪声模型。实验测量时，对齐采样得到的闪烁脉冲，每个时间点的电压值分布（闪烁脉冲的噪声模型）即可得到。
所述步骤S1中求平均脉冲的具体步骤包括：
（1.1）降低射源的辐射剂量得到弱源，通过弱源或调整弱源对探测器的立体角，降低每个探测器捕获的高能光子数，在极低的计数率下，发生堆积的概率较小。每个探测器接收到的事件是一个泊松流，其平均计数率为
其中，mi和qi分别是弱源的剂量和弱源对探测器的立体角。i是弱源的序数，n是弱源的个数。在获得脉冲数据库时，可以令n=1。并令mi足够小。
（1.2）按照恒定系数鉴别方法或前沿甄别方法对齐脉冲；
（1.3）对对齐后的脉冲取平均。
所述步骤S1中对齐脉冲后，同一个能量值的闪烁脉冲在一个时间点上能量值会有一个分散的分布，这个分布被定义为脉冲的噪声模型，噪声模型的参数留作计算给定能量值的后验概率值。
如图4所示，图4为一个能量的闪烁脉冲集合示意图。所述步骤S1还包括：将一能量区间选择的符合脉冲对认定为单事件脉冲对，将此数字化的电脉冲预先存储后，进行离线分析，通过对齐脉冲对，获得平均脉冲信号和脉冲残差序列，将平均脉冲作为系统响应，将数字脉冲加和，得到单事件数据的能谱。
如未有特别注明，本发明所述的闪烁脉冲均为正值。将数字脉冲加和，可以得到单事件数据的能谱。如实际获取的脉冲是负脉冲，则默认地经过了反相操作。因此，上升沿即是电脉冲的前沿，下降沿即是电脉冲的尾部。处理原始脉冲数据库得到的所有中间信息可称作衍生脉冲数据库，该衍生脉冲数据库可以是在原始脉冲数据库基础上所做的滤波、插值、拟合、外推、组合等等。
S2：根据闪烁脉冲的噪声模型计算给定能量值的后验概率对数值。
最大后验概率条件的能量值反应了一个统计学上的基本原理：当前观测在一个能量出现的可能性比其他能量出现这种观测的可能性更大。这个可能性表达为能量的先验概率乘以能量出现观测的似然概率。
所述S2具体包括步骤：
（2.1）载入一个闪烁脉冲片段S0，该闪烁脉冲片段的起点是向上过阈值Vl的时间点t0，终点是t0+Δt，其中Δt是一个大于闪烁晶体下降沿时间常数2倍的值，根据上面得到的脉冲噪声噪声模型，计算给定能量值的似然函数；
（2.2）将所有时间点取对数，再相加，得到一个单调性与后验概率相同的函数的具体值。
这个似然函数在实际计算中是不需要直接计算出来的，而是通过计算似然函数的对数值来得到的。这个对数似然函数是等于脉冲每个时间点的对数似然函数的加和。由于这个似然函数可能具有别的形式，其形式也是不固定的。
这个时间点就是脉冲到达时间的意思。我们以这个脉冲到达时间作为脉冲不同位置的零点。由于取对数并不改变单调性，所以当一个函数的对数值达到极值时，他的对数也达到了极值。所以最大化某一个函数等于最大化某个函数的对数。所以与S2描述的内容是一样，最大化后验概率的对数就是最大化后验概率。
S3：通过计算反复地调用步骤S2，获得满足最大后验概率条件的能量值。
步骤S3中在每次计算能量之后需要对脉冲的时间基点进行能量校正，这个校正后的时间点除了使能量值更为准确，其本身也是一个准确的时间信息。在PET或者核测量设备中，时间信息和能量信息是表达一个单事件脉冲全部信息。形状信息已经由平均脉冲给定。这个形状信息是绝大多数核测量设备所不需要的。
所述步骤S3具体包括步骤：
（3.1）将步骤2作为一个模块反复调用，计算不同试探能量值的后验概率；
（3.2）线性搜索满足最大后验概率的能量值，现有的线性搜索算方法包括黄金分割法，二分法，斐波那契级数法等；如图16所示，图16为一个能量值一维解寻优搜索过程图。
最大化后验概率等于最大化似然函数乘以能量的概率分布函数。后验概率等于似然函数乘以能量的概率分布函数。
（3.3）计算得到能量值后，根据能量值校准2.1中的时间基点，再重复3.1和3.2。
这个能量值是基于数字化得样本来算的，即这个测得的能量值决定于时间。在没有能量值修正的时候，时间是不准的，会遇到“时间游走”问题。所以在得到了能量后需要对之前较为不准的时间进行修正。修正了时间后，再对能量进行计算。得到了能量后再次修正时间。重复时间-能量-时间-能量的迭代过程，就可以提高时间和能量测量的准确度。所以引入时间相关的迭代就是为了提高测量的准确度。
如图6至图11所示，图6至图11为在能量分辨率较差和能量分辨率较好两种情况下，采用不同的能量计算方法得到的能谱示意图。采用的能量计算方法包括数字门控积分、本发明所提出的方法、最小二乘拟合算法。
图6中为一个能量分辨率较差的能谱示意图，其能量计算方法为数字门控积分。能量分辨率的计算方法是用峰宽（半高宽）除以峰位。门控积分法是说利用上升沿的阈值翻转，得到一个门控信号。在这个门控区间内对脉冲进行积分，归一化后的积分值表示为能量值。归一化方法为将全能峰归一化到511keV（千电子伏）。
图7为一个能量分辨率较差的能谱示意图，其能量计算方法为本发明所提出的方法。能量分辨率和归一化的计算方法同图6。
图8为一个能量分辨率较差的能谱示意图，其能量计算方法为最小二乘拟合算法。能量分辨率和归一化的计算方法同图6。
图9为一个能量分辨率较好的能谱示意图，其能量计算方法为数字门控积分。能量分辨率和归一化的计算方法同图6。
图10为一个能量分辨率较好的能谱示意图，其能量计算方法为本发明所提出的方法。能量分辨率和归一化的计算方法同图6。
图11为一个能量分辨率较好的能谱示意图，其能量计算方法为最小二乘拟合算法。能量分辨率和归一化的计算方法同图6。
如图12至图15所示，图12至图15为不同能量值情况下，最小二乘拟合算法和本发明方法得到的能量值对数字门控积分得到的能量值的散点图。
图12为最小二乘拟合算法得到的能量值对数字门控积分得到的能量值的散点图，能量值为图4的能量值。散点图的横坐标为门控积分的能量值，纵坐标为脉冲拟合的能量值。以门控积分为金标准，其误差水平可以用来衡量能量计算方法的精确度。用误差均方根（Root Mean Square Error，RMSE）表示能量精确度，可以看出各种能量计算方法的性能。
图13为本发明方法得到的能量值对数字门控积分得到的能量值的散点图，能量值为图4的能量值。散点图的横坐标为门控积分的能量值，纵坐标为脉冲拟合的能量值。
图14为最小二乘拟合算法得到的能量值对数字门控积分得到的能量值的散点图，能量值为图7的能量值。散点图的横坐标为门控积分的能量值，纵坐标为本方法给出的能量值。
图15为本发明方法得到的能量值对数字门控积分得到的能量值的散点图，能量值为图7的能量值。散点图的横坐标为门控积分的能量值，纵坐标为本方法给出的能量值。
从能谱图6-11中可以看出，本方法相对于传统的脉冲拟合方法的能量分辨率更好，更有利于拒绝散射事件。而从散点图12-15中可以看出本方法的精确度高，误差较小，更适合于对脉冲高度进行估计，适用于一些位置灵敏型的探测器设计。图16的搜索过程是本方法定义能量区间的缩小过程，示意图刻画了其缩小过程，即每个区间缩小0.618倍。
本发明闪烁脉冲信息复原的方法及系统，有效提高了系统能量计算精确度，特别适合于稀疏量化电平ADC数字式核仪器的能量计算。
如图1所示，本发明公开的闪烁脉冲信息复原的系统包括噪声模型模块100、后验概率模块200、能量值搜索模块300。其中噪声模型模块100用于对预采集的低计数率数据进行训练，获得闪烁探测器系统的噪声模型；后验概率模块200用于后验概率的计算，获得给定能量值的后验概率对数值；能量值搜索模块300用于搜索满足最大后验概率条件的能量值。
继续如图1所示，噪声模型模块100用于对预先采集好的数据进行训练，得到平均脉冲、噪声方差和噪声分布偏度。该噪声模型提供给模块200。
噪声模型模块100可以分为三个子模块，分别为平均脉冲模块110、脉冲噪声模块120、概率密度分布模块130。平均脉冲模块110用于计算单事件的平均脉冲。脉冲噪声模块120用于计算单事件脉冲的噪声方差。概率密度分布模块130用于估计噪声的概率密度分布函数。
后验概率模块200用于计算给定能量值的后验概率对数值。
后验概率模块200可以分为两个子模块，分别为似然函数分布模块210、时间点乘积模块220。似然函数分布模块210按照给定能量值和模块100输入的脉冲噪声模型计算每个点的似然函数值。时间点乘积模块220负责计算后验概率对数值。
能量值搜索模块300用于计算满足最大后验概率条件的能量值。
能量值搜索模块300可以分为三个子模块，分别为试探能量计算模块310、线性搜索模块320、时间基点校正模块330。试探能量计算模块310通过调用模块200计算出每个试探能量的后验概率对数值。线性搜索模块320通过线性搜索算法计算满足最大后验概率条件的能量值。时间基点校正模块330校正样本的时间基点，并重复调用模块310和模块320。
如图3至图5所示，图3为本发明一种典型闪烁脉冲的示意图；图4为一个能量的闪烁脉冲集合示意图；图5为图3的闪烁脉冲在四个阈值的ADC下的数字化结果。图5中V1表示参考电压1，V2表示参考电压2，V3表示参考电压3，V4表示参考电压4。如图17及图18所示，图17为本发明一种符合工作模式典型系统的示意图；图18为本发明另一种单通道工作模式典型系统的示意图。其中，400代表闪烁脉冲，500代表射源，600代表光电倍增管，700代表数字示波器。结合图3至图5以及图17及图18，通过几个具体的实施例，对本发明闪烁脉冲信息复原的方法及系统做进一步描述。本发明提出的闪烁脉冲信息复原的方法及系统，其涉及到的参数、滤波器设计需要根据与获取数据的特点进行调节以达到良好的能量分辨性能和较短的脉冲持续时间。此处列出所涉及的应用实施例处理数据的参数。
实施例1：
此处列出本实施例处理数据的参数：
步骤（1）所用的实际系统为使用LYSO晶体和Hamamatsu R9800PMT。晶体尺寸为16.5mm×16.5mm×10.0mm。晶体和PMT耦合面为100面，除开耦合面外，其他面均用特氟龙胶带包裹。数据获得系统的采样率为50Ghz，带宽16Ghz。如图10所示，射源为511kev的正电子湮灭伽马光子。平均脉冲上升时间约为2ns，下降沿用指数拟合后，时间常数为42.5497ns。
步骤（2.1）采用的Δt=200ns.
实施例2：
此处列出本实施例处理数据的参数：
如图11所示，步骤（1）所用的实际系统为使用LYSO晶体和FM300035SIPM。晶体尺寸为2.0mm×2.0mm×10.0mm。晶体和PMT耦合面为100面，除开耦合面外，其他面均用特氟龙胶带包裹。数据获得系统的采样率为50Ghz，带宽16Ghz。如图10所示，射源为511kev的正电子湮灭伽马光子。平均脉冲上升时间约为5ns。
步骤（2.1）采用的Δt=300ns.
实施例3：
此处列出本应用实例处理数据的参数：
步骤（1）所用的实际系统为使用LaBr晶体和HamamatsuR9800PMT。晶体尺寸为3.5mm×3.5mm×5.0mm。晶体和PMT耦合面为100面，除开耦合面为光导外，其他面均用金属密封。数据获得系统的采样率为50Ghz，带宽16Ghz。如图10所示，射源为511kev的正电子湮灭伽马光子。平均脉冲上升时间约为2ns。
步骤（2.1）采用的Δt=100ns.
本发明闪烁脉冲信息复原的方法及系统可以用于要求对信号进行直接数字化的核探测、核分析、核医学仪器。
对于本领域技术人员而言，显然本发明不限于上述示范性实施例的细节，而且在不背离本发明的精神或基本特征的情况下，能够以其他的具体形式实现本发明。因此，无论从哪一点来看，均应将实施例看作是示范性的，而且是非限制性的，本发明的范围由所附权利要求而不是上述说明限定，因此旨在将落在权利要求的等同要件的含义和范围内的所有变化囊括在本发明内。不应将权利要求中的任何附图标记视为限制所涉及的权利要求。
此外，应当理解，虽然本说明书按照实施方式加以描述，但并非每个实施方式仅包含一个独立的技术方案，说明书的这种叙述方式仅仅是为清楚起见，本领域技术人员应当将说明书作为一个整体，各实施例中的技术方案也可以经适当组合，形成本领域技术人员可以理解的其他实施方式。
Claims (10)
1.一种闪烁脉冲信息复原的方法，其特征在于：具体包括步骤：
S1：获取低计数下的非堆积符合单事件的闪烁脉冲数据库，再对各路单事件的闪烁脉冲数据库建立闪烁脉冲的噪声模型，对各路非堆积符合单事件的闪烁脉冲数据库求平均脉冲，闪烁脉冲的形状信息由平均脉冲给定；
S2：根据闪烁脉冲的噪声模型计算给定能量值的后验概率对数值，具体包括步骤：
（2.1）载入一个闪烁脉冲片段S0，该闪烁脉冲片段的起点是向上过阈值Vl的时间点t0，终点是t0+Δt，其中Δt是一个大于闪烁晶体下降沿时间常数两倍的值，根据脉冲噪声模型，计算给定能量值的似然函数；
（2.2）将所有时间点取对数，再相加，得到一个单调性与后验概率相同的函数的具体值；
S3：通过计算反复地调用步骤S2，获得满足最大后验概率条件的能量值，具体包括步骤：
（3.1）将步骤S2作为一个模块反复调用，计算不同试探能量值的后验概率；
（3.2）线性搜索满足最大后验概率的能量值；
（3.3）计算得到能量值后，根据能量值校准2.1中的时间基点，再重复3.1和3.2。
2.根据权利要求1所述的闪烁脉冲信息复原的方法，其特征在于：所述步骤S1中闪烁脉冲数据库中的样本数在2000个以上。
3.根据权利要求1所述的闪烁脉冲信息复原的方法，其特征在于：所述步骤S1中闪烁脉冲的噪声模型包括特征：闪烁脉冲形状、闪烁脉冲高度、闪烁脉冲噪声水平。
4.根据权利要求1所述的闪烁脉冲信息复原的方法，其特征在于：所述步骤S1中求平均脉冲的具体步骤包括：
（1.1）降低射源的辐射剂量得到弱源，通过弱源或调整弱源对探测器的立体角，降低每个探测器捕获的高能光子数，每个探测器接收到的事件是一个泊松流，其平均计数率为
其中，mi和qi分别是弱源的剂量和弱源对探测器的立体角，i是弱源的序数，n是弱源的个数；
（1.2）按照恒定系数鉴别方法或前沿甄别方法对齐脉冲；
（1.3）对对齐后的脉冲取平均。
5.根据权利要求4所述的闪烁脉冲信息复原的方法，其特征在于：所述步骤S1中对齐脉冲后，同一个能量值的闪烁脉冲在一个时间点上能量值会有一个分散的分布，这个分布被定义为脉冲的噪声模型，噪声模型的参数留作计算给定能量值的后验概率值。
6.根据权利要求1所述的闪烁脉冲信息复原的方法，其特征在于：所述步骤S1还包括：将一能量区间选择的符合脉冲对认定为单事件脉冲对，将此数字化的电脉冲预先存储后，进行离线分析，通过对齐脉冲对，获得平均脉冲信号，将平均脉冲作为系统响应，将数字脉冲加和，得到单事件数据的能谱。
7.一种闪烁脉冲信息复原的系统，其特征在于：包括涨落模型模块、后验概率模块及能量值搜索模块，其中，
所述涨落模型模块用于对预采集的低计数率数据进行训练，获得闪烁探测器系统的涨落模型，该涨落模型表征平均脉冲、噪声方差和噪声分布偏度，该噪声模型提供给后验概率模块；
所述后验概率模块用于后验概率的计算，获得给定能量值的后验概率对数值；
所述能量值搜索模块用于搜索满足最大后验概率条件的能量值。
8.根据权利要求7所述的闪烁脉冲信息复原的系统，其特征在于：所述涨落模型模块包括平均脉冲模块、脉冲噪声模块及概率密度分布模块，其中，
所述平均脉冲模块用于计算单事件的平均脉冲；
所述脉冲噪声模块用于计算单事件脉冲的噪声方差；
所述概率密度分布模块用于估计噪声的概率密度分布函数。
9.根据权利要求1所述的闪烁脉冲信息复原的系统，其特征在于：所述后验概率模块包括似然函数分布模块及时间点乘积模块，其中，
所述似然函数分布模块按照给定能量值和涨落模型模块输入的脉冲噪声模型计算每个点的似然函数值；
所述时间点乘积模块负责计算后验概率对数值。
10.根据权利要求1所述的闪烁脉冲信息复原的系统，其特征在于：所述能量值搜索模块包括试探能量计算模块、线性搜索模块及时间基点校正模块，其中，
所述试探能量计算模块通过调用后验概率模块计算出每个试探能量的后验概率对数值；
所述线性搜索模块通过线性搜索算法计算满足最大后验概率条件的能量值；
所述时间基点校正模块校正样本的时间基点，并重复调用试探能量计算模块和线性搜索模块。
Priority Applications (5)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201310585902.1A CN104656119B (zh) | 2013-11-19 | 2013-11-19 | 一种闪烁脉冲信息复原的方法及系统 |
| US15/037,660 US10371827B2 (en) | 2013-11-19 | 2013-12-25 | Method and system for recovering scintillation pulse information |
| PCT/CN2013/090389 WO2015074312A1 (zh) | 2013-11-19 | 2013-12-25 | 一种闪烁脉冲信息复原的方法及系统 |
| JP2016532558A JP6336070B2 (ja) | 2013-11-19 | 2013-12-25 | シンチレーションパルス情報を回復するための方法およびシステム |
| EP13897953.9A EP3073292B1 (en) | 2013-11-19 | 2013-12-25 | Method and system for recovering scintillation pulse information |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201310585902.1A CN104656119B (zh) | 2013-11-19 | 2013-11-19 | 一种闪烁脉冲信息复原的方法及系统 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN104656119A true CN104656119A (zh) | 2015-05-27 |
| CN104656119B CN104656119B (zh) | 2018-06-05 |
Family
ID=53178852
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201310585902.1A Active CN104656119B (zh) | 2013-11-19 | 2013-11-19 | 一种闪烁脉冲信息复原的方法及系统 |
Country Status (5)
| Country | Link |
|---|---|
| US (1) | US10371827B2 (zh) |
| EP (1) | EP3073292B1 (zh) |
| JP (1) | JP6336070B2 (zh) |
| CN (1) | CN104656119B (zh) |
| WO (1) | WO2015074312A1 (zh) |
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105005068A (zh) * | 2015-06-25 | 2015-10-28 | 南京瑞派宁信息科技有限公司 | 一种脉冲分类的方法与系统 |
| CN108646285A (zh) * | 2018-06-12 | 2018-10-12 | 高超 | 伽马射线能谱测量方法、装置及系统 |
| CN110208846A (zh) * | 2019-04-30 | 2019-09-06 | 韶关学院 | 一种基于压缩感知的软x光能谱还原方法 |
| CN111025370A (zh) * | 2019-12-10 | 2020-04-17 | 南昌大学 | 一种双态功能材料辐射传感装置及其组装方法 |
Families Citing this family (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105223906B (zh) * | 2015-09-15 | 2017-10-03 | 华中科技大学 | 一种数控系统伺服驱动信号谐波频率的自动校正方法 |
| CN109171787B (zh) * | 2018-08-27 | 2021-02-26 | 苏州瑞派宁科技有限公司 | 脉冲信号的采样方法、装置和计算机程序介质 |
| CN111239797B (zh) * | 2020-02-10 | 2022-09-30 | 成都理工大学 | 一种基于辐射粒子事件的采集器及快速核素识别方法 |
| CN113095037B (zh) * | 2021-03-30 | 2023-11-03 | 上海华力微电子有限公司 | Mos器件闪烁噪声模型及提取方法 |
| CN116299629A (zh) * | 2022-12-30 | 2023-06-23 | 合肥锐世数字科技有限公司 | 闪烁脉冲的处理方法、装置、设备及存储介质 |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20110017918A1 (en) * | 2008-03-19 | 2011-01-27 | Koninklijke Philips Electronics N.V. | Single photon radiation detector |
| CN102073059A (zh) * | 2010-12-31 | 2011-05-25 | 华中科技大学 | 一种数字化pileup波形处理方法及系统 |
| CN102262238A (zh) * | 2011-04-19 | 2011-11-30 | 苏州瑞派宁科技有限公司 | 一种提取闪烁脉冲信息的方法及装置 |
Family Cites Families (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| EP1007988B1 (en) * | 1997-05-07 | 2003-07-23 | The Board Of Regents, The University Of Texas System | Method and apparatus to prevent pile-up when detecting the energy of incoming signals |
| JP4064009B2 (ja) * | 1999-07-30 | 2008-03-19 | 株式会社東芝 | 線種弁別型放射線検出装置 |
| JP4160275B2 (ja) * | 2001-05-28 | 2008-10-01 | 浜松ホトニクス株式会社 | エネルギー測定方法及び測定装置 |
| US8330115B2 (en) * | 2005-12-01 | 2012-12-11 | Innovative American Technology, Inc. | High performance neutron detector with near zero gamma cross talk |
| US8446308B2 (en) | 2011-04-21 | 2013-05-21 | Kabushiki Kaisha Toshiba | Apparatus for detection of a leading edge of a photo sensor output signal |
-
2013
- 2013-11-19 CN CN201310585902.1A patent/CN104656119B/zh active Active
- 2013-12-25 JP JP2016532558A patent/JP6336070B2/ja active Active
- 2013-12-25 EP EP13897953.9A patent/EP3073292B1/en active Active
- 2013-12-25 US US15/037,660 patent/US10371827B2/en active Active
- 2013-12-25 WO PCT/CN2013/090389 patent/WO2015074312A1/zh active Application Filing
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20110017918A1 (en) * | 2008-03-19 | 2011-01-27 | Koninklijke Philips Electronics N.V. | Single photon radiation detector |
| CN102073059A (zh) * | 2010-12-31 | 2011-05-25 | 华中科技大学 | 一种数字化pileup波形处理方法及系统 |
| CN102262238A (zh) * | 2011-04-19 | 2011-11-30 | 苏州瑞派宁科技有限公司 | 一种提取闪烁脉冲信息的方法及装置 |
Non-Patent Citations (2)
| Title |
|---|
| ZHENZHOU DENG ET AL: "Empirical Bayesian Energy Estimation for Multi-Voltage Threshold Digitizer in PET", 《2013 IEEE NUCLEAR SCIENCE SYMPOSIUM AND MEDICAL IMAGING CONFERENCE》 * |
| 龙岸文 等: "基于PET列表式单事件数据处理的能量校正", 《中国核科学技术进展报告》 * |
Cited By (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105005068A (zh) * | 2015-06-25 | 2015-10-28 | 南京瑞派宁信息科技有限公司 | 一种脉冲分类的方法与系统 |
| CN105005068B (zh) * | 2015-06-25 | 2018-04-06 | 南京瑞派宁信息科技有限公司 | 一种脉冲分类的方法与系统 |
| CN108646285A (zh) * | 2018-06-12 | 2018-10-12 | 高超 | 伽马射线能谱测量方法、装置及系统 |
| CN108646285B (zh) * | 2018-06-12 | 2020-05-22 | 高超 | 伽马射线能谱测量方法、装置及系统 |
| CN110208846A (zh) * | 2019-04-30 | 2019-09-06 | 韶关学院 | 一种基于压缩感知的软x光能谱还原方法 |
| CN111025370A (zh) * | 2019-12-10 | 2020-04-17 | 南昌大学 | 一种双态功能材料辐射传感装置及其组装方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| EP3073292B1 (en) | 2020-04-29 |
| JP2017503151A (ja) | 2017-01-26 |
| EP3073292A1 (en) | 2016-09-28 |
| WO2015074312A1 (zh) | 2015-05-28 |
| EP3073292A4 (en) | 2017-07-12 |
| US10371827B2 (en) | 2019-08-06 |
| JP6336070B2 (ja) | 2018-06-06 |
| US20160291167A1 (en) | 2016-10-06 |
| CN104656119B (zh) | 2018-06-05 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN104656119A (zh) | 一种闪烁脉冲信息复原的方法及系统 | |
| Imperiale et al. | On nuclear spectrometry pulses digital shaping and processing | |
| US9772408B2 (en) | Method and device for extracting scintillation pulse information | |
| CN105824817B (zh) | 一种闪烁脉冲的数字化方法 | |
| CN110082368B (zh) | 一种基于硅光电倍增器(SiPM)的正电子湮没寿命谱仪 | |
| EP3859403A1 (en) | Method and device for processing nuclear energy spectrum | |
| Warburton et al. | New algorithms for improved digital pulse arrival timing with sub-GSps ADCs | |
| CN109142403B (zh) | 一种用于地球元素测井仪解谱获取元素产额的方法 | |
| CN109828297A (zh) | 一种多重粒子事件的检测方法与装置 | |
| WO2024140347A1 (zh) | 闪烁脉冲的处理方法、装置、设备及存储介质 | |
| Xiao et al. | A new pulse model for NaI (Tl) detection systems | |
| CN105005068B (zh) | 一种脉冲分类的方法与系统 | |
| CN105181153B (zh) | 基于波形面积的雪崩光探测器增益测量方法 | |
| CN115420226B (zh) | 基于脉冲宽度的射线作用位置定位装置及定位方法 | |
| CN110716226B (zh) | 一种用于剔除无效事例的反符合电路及方法 | |
| CN104730566A (zh) | 一种标定光电转换模块平均探测器响应函数的方法和系统 | |
| CN104931997B (zh) | X射线能量的平场图像测量方法 | |
| Pollastrone et al. | A clustering algorithm for scintillator signals applied to neutron and gamma patterns identification | |
| CN114675320B (zh) | 一种混合β能谱的解谱方法、系统及存储介质 | |
| Ivanova et al. | An adaptive γ-ray spectrometer with a high event processing rate | |
| Ivanova et al. | Gamma-ray spectrometer with high event processing rate | |
| Nattress et al. | Detector Characterization for Accurate Monte Carlo Simulation and Neutron Energy-spectrum Unfolding | |
| Guan et al. | Study of Silicon Photomultiplier external cross-talk | |
| CN111789623A (zh) | 一种pet信号的全带压缩感知方法 | |
| Deng et al. | Pileup Scintillation Pulses Processing using Expectation Maximization Deconvolution |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |