CN113826030A - 具有脉冲堆积非参数分解的辐射检测 - Google Patents

具有脉冲堆积非参数分解的辐射检测 Download PDF

Info

Publication number
CN113826030A
CN113826030A CN202080037871.4A CN202080037871A CN113826030A CN 113826030 A CN113826030 A CN 113826030A CN 202080037871 A CN202080037871 A CN 202080037871A CN 113826030 A CN113826030 A CN 113826030A
Authority
CN
China
Prior art keywords
pulse
statistics
estimate
noise
mapping
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.)
Pending
Application number
CN202080037871.4A
Other languages
English (en)
Inventor
C·麦克莱恩
M·波利
J·曼顿
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Southern Innovation International Pty Ltd
Original Assignee
Southern Innovation International Pty Ltd
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Priority claimed from AU2019900974A external-priority patent/AU2019900974A0/en
Application filed by Southern Innovation International Pty Ltd filed Critical Southern Innovation International Pty Ltd
Publication of CN113826030A publication Critical patent/CN113826030A/zh
Pending legal-status Critical Current

Links

Images

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/17—Circuit arrangements not adapted to a particular type of detector
    • G01T1/171—Compensation of dead-time counting losses
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R23/00—Arrangements for measuring frequencies; Arrangements for analysing frequency spectra
    • G01R23/02—Arrangements for measuring frequency, e.g. pulse repetition rate; Arrangements for measuring period of current or voltage
    • G01R23/15—Indicating that frequency of pulses is either above or below a predetermined value or within or outside a predetermined range of values, by making use of non-linear or digital elements (indicating that pulse width is above or below a certain limit)
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R23/00—Arrangements for measuring frequencies; Arrangements for analysing frequency spectra
    • G01R23/02—Arrangements for measuring frequency, e.g. pulse repetition rate; Arrangements for measuring period of current or voltage
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R23/00—Arrangements for measuring frequencies; Arrangements for analysing frequency spectra
    • G01R23/16—Spectrum analysis; Fourier analysis
    • 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
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R29/00—Arrangements for measuring or indicating electric quantities not covered by groups G01R19/00 - G01R27/00
    • G01R29/02—Measuring characteristics of individual pulses, e.g. deviation from pulse flatness, rise time or duration
    • G01R29/027—Indicating that a pulse characteristic is either above or below a predetermined value or within or beyond a predetermined range of values
    • G01R29/033—Indicating that a pulse characteristic is either above or below a predetermined value or within or beyond a predetermined range of values giving an indication of the number of times this occurs, i.e. multi-channel analysers (the characteristic being frequency)

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Health & Medical Sciences (AREA)
  • High Energy & Nuclear Physics (AREA)
  • Molecular Biology (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Physics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Nonlinear Science (AREA)
  • Probability & Statistics with Applications (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Operations Research (AREA)
  • Algebra (AREA)
  • Evolutionary Biology (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Measurement Of Radiation (AREA)

Abstract

公开了一种确定在辐射检测器中接收的单个辐射量子的能谱的方法。能谱敏感统计是根据来自辐射检测器的数字观察的时间序列计算的，定义从脉冲幅度密度到能谱敏感统计的映射。通过对能谱敏感统计应用映射的逆来估计脉冲幅度的密度，从而确定能谱。统计可以基于：在不考虑脉冲簇的整体的情况下恒定长度L的第一组不重叠的时间间隔，所述恒定长度L至少与脉冲的持续时间一样长；以及在也不考虑脉冲簇的整体的情况下小于L的恒定长度L1的第二组不重叠的时间间隔。还公开了一种估计计数率的方法。

Description

具有脉冲堆积非参数分解的辐射检测
技术领域
提供了一种用于估计辐射量子的能量分布的新奇方法，所述量子诸如是入射在诸如X射线或伽马射线光谱学的光谱系统中的检测器上的光子。该方法对于脉冲堆积是一个问题的计数率方法特别有用。推导本发明实施例的估计量的关键步骤是将问题新奇重新表述为复合泊松过程的分解问题。该方法可以应用于任何形式的辐射检测器，其检测辐射的量子或其他微粒，诸如x射线、伽马射线或其他光子、中子、原子、分子或地震脉冲。来自这种检测器的光谱学的应用是众所周知的。这样的申请在现有技术中被广泛描述，现有技术包括在国际专利申请PCT/AU2005/001423、PCT/AU2009/000393、PCT/AU2009/000394、PCT/AU2009/000395、PCT/AU2009/001648、PCT/AU2012/000678、PCT/AU2014/050420、PCT/AU2015/050752、PCT/AU2017/050514和PCT/AU2017/050512中，出于描述当前发明的潜在应用和理解当前发明所需要的任何其他背景材料的目的，将其中的每一项以其整体并入本文中。
背景技术
X射线和伽马射线光谱学支撑着宽范围的科学、工业和商业过程。光谱学的一个目标是估计入射在检测器上的光子的能量分布。从信号处理的角度来看，挑战在于将检测器输出的脉冲流转换成每个脉冲下面积的直方图。脉冲是根据泊松分布生成的，泊松分布的比率对应于用于照射样本的X射线或伽马射线的强度。增加强度导致平均每秒更多脉冲，并且因此导致在获得准确直方图之前的更短时间。在诸如机场行李扫描之类的应用中，这直接转化为更大的吞吐量。当两个或更多个脉冲在时域中重叠时，出现脉冲堆积。随着计数率（平均每秒脉冲数）增加，脉冲堆积的发生率也增加。这增加了在确定存在的脉冲数量和每个脉冲下的面积方面的难度。在极限下，该问题是病态的：如果两个脉冲实质上上同时开始，则它们的叠加就无法与单个脉冲区分。X射线或伽马射线检测器对入射光子的响应可以建模为脉冲形状函数
Figure 94662DEST_PATH_IMAGE001
与狄拉克-德尔塔（Dirac-delta）脉冲的卷积的叠加，
Figure 11802DEST_PATH_IMAGE002
。
到达时间
Figure 853856DEST_PATH_IMAGE003
未知，并且形成泊松过程。每个光子到达被建模为时间
Figure 249066DEST_PATH_IMAGE004
的狄拉克-德尔塔，并且具有与光子能量成比例并引起检测器脉冲形状响应
Figure 63438DEST_PATH_IMAGE005
的幅度
Figure 696020DEST_PATH_IMAGE006
。幅度
Figure 25370DEST_PATH_IMAGE006
是其公共概率密度函数
Figure 224270DEST_PATH_IMAGE007
未知的同分布随机变量
Figure 893149DEST_PATH_IMAGE008
的实现。脉冲形状函数
Figure 886513DEST_PATH_IMAGE009
由检测器的几何形状和光子相互作用确定。在一些系统中，脉冲形状的变化极小，并且可以忽略，而其他系统（例如HPGe）单个脉冲形状可能彼此显著不同[1]。假设所有脉冲形状都是有因果的，即对于t < 0，
Figure 453892DEST_PATH_IMAGE010
，单峰，有限能量，并且随着
Figure 190903DEST_PATH_IMAGE011
而指数地向零衰减。脉冲形状函数的积分被归一化为单位，即
Figure 714289DEST_PATH_IMAGE013
，使得脉冲下的面积由
Figure 206450DEST_PATH_IMAGE014
给出。观察的信号由被噪声破坏的检测器输出组成，即，
Figure 448075DEST_PATH_IMAGE015
。
脉冲堆积校正的数学目标是在给定
Figure 723199DEST_PATH_IMAGE016
的均匀采样、有限长度版本的情况下对
Figure 179719DEST_PATH_IMAGE017
进行估计。我们始终假设噪声分布
Figure 514886DEST_PATH_IMAGE018
是具有已知方差
Figure 243807DEST_PATH_IMAGE020
的零均值高斯分布。我们还假设光子到达时间形成具有已知比率的齐次泊松过程。由被噪声过程W破坏的检测器响应R组成，其中
Figure 384938DEST_PATH_IMAGE022
、
Figure 882916DEST_PATH_IMAGE024
和
Figure 654563DEST_PATH_IMAGE026
表示每个序列的第k个元素，并且其中
Figure 605201DEST_PATH_IMAGE028
。让S、R和W是对应于这些信号的均匀采样时间序列，
Figure 300756DEST_PATH_IMAGE029
其中
Figure 653240DEST_PATH_IMAGE030
。
脉冲处理方法概述：几十年来，已经提出了许多方法来解决脉冲堆积的问题。方法可以大体上分类为两个类型：基于时域的和基于能量域的。一种流行的策略是试图检测时域中何时发生堆积，并拒绝或补偿受影响的脉冲。早期的光谱系统采用了基于拒绝的方法连同匹配滤波。该方法的缺点是，随着堆积概率增长，越来越多比例的脉冲被拒绝。系统迅速陷入瘫痪，从而给计数率设置了上限[1]。随着廉价计算能力的增加，补偿或校正堆积的策略的数量已经增长。这些包括模板拟合[2]、基线减法[3]、自适应滤波[4，5]、稀疏回归[6，7]等。这些方法都试图标识和补偿时域中的堆积，并且通常最适合于具有低脉冲形状变化的系统。这些方法的复杂度随着脉冲形状
Figure 595788DEST_PATH_IMAGE031
之间可变性的增加而显著增加。可以示出，任何试图表征单个脉冲的方法都将受到堆积影响。这些方法能做的最好的事情是减少堆积的发生。基于能量的方法试图基于脉冲集而不是单个脉冲的统计来解决堆积。它们通常对估计能量（脉冲簇下的面积）的直方图进行操作。Wielopolski和Gardner的早期工作[8]和他们的想法的更近的扩展[9]主要在能量域中使用基于集合的策略来操作。Trigano等人[10，11]利用边缘密度估计入射光谱，该边缘密度来自可变长度脉冲簇的统计特性的联合分布，其中检测到每个簇的开始和结束。这避免了标识单个脉冲的需要，并且对脉冲形状变化是鲁棒的。Ilhe等人[12]考察了指数散粒噪声过程，将脉冲形状限制为简单的指数，以获得易于处理的结果。已经做了进一步的工作[13]以允许更宽范围的脉冲形状。在这两种情况下，都需要脉冲形状的知识，连同特征函数及其导数的估计。
发明内容
我们选择了基于能量的堆积校正方法，以便i）避免与单个脉冲检测相关联的限制[14]和ii）能够处理脉冲形状变化，而不过度增加计算复杂度。我们不是利用联合分布[10，11]或散粒过程[12]方法，而是将堆积问题重铸为复合泊松过程的“分解”问题。复合泊松过程是离散时间随机过程，其中每个分量由独立同分布随机变量的随机数的总和组成，其中每个总和中的随机变量数是泊松分布的[15]。复合泊松过程的“分解”是使用随机和来估计从其中已经提取随机变量的分布的任务。Buchmann和
Figure 96040DEST_PATH_IMAGE032
[16]在保险索赔和排队论的背景下阐述了复合泊松过程的分解。均匀采样复合泊松过程的分解在最近时间受到了一些关注[16，17，18，12，19]。这些推导的上下文经常（合理地）假设每个事件都是可检测的（即，关于事件的数量没有歧义），或者密度估计量以每个观察中至少发生一个事件为条件[20]。这些假设在解决光谱堆积问题时价值有限。
对事件检测的无条件非参数分解研究在文献中受到相对少的关注。Gugushvili[18]提出了一种非参数的、基于核的估计量，用于存在高斯噪声下的分解问题。在本发明的实施例中，发明人已经设想，一旦获得了用于选择核带宽的方法，连同用于将观察的检测器输出变换以拟合数学模型的方法，该估计量可以容易地扩展并应用于光谱堆积问题的重构。
根据本发明的第一广义方面，提供了一种确定在辐射检测器中接收的单个辐射量子的能谱的方法，该方法包括以下步骤：（1）从辐射检测器获得数字观察的时间序列，该时间序列包括与单个量子的检测相对应的脉冲；（2）根据检测器信号计算能谱敏感统计，该能谱敏感统计定义从脉冲幅度密度到能谱敏感统计的映射；（3）通过对能谱敏感统计应用映射的逆来估计脉冲幅度的密度，从而确定能谱。
在实施例中，能谱敏感统计可以基于在多个时间间隔上的数字观察的总和，并且可以使用近似复合泊松过程来定义映射，该近似复合泊松过程可以通过建模噪声来增强。该映射可以表述为幅度的特征函数、能谱敏感统计和建模噪声之间的关系。能谱敏感统计的特征函数可以通过使用数字观察总和的直方图来计算，对该直方图应用逆傅立叶变换。幅度的特征函数的计算可以包括使用低通滤波器。
在第一实施例中，多个时间间隔是不重叠的，并且具有恒定的长度L，并且每个间隔被选择为包含零个或更多个近似完整的脉冲簇。这可以通过在每个时间间隔的开始和结束时要求检测器信号的最大值来实现。在该实施例中，复合泊松过程可以定义为每个时间间隔内脉冲幅度的总和。该映射可以按照等式（40）和（41）中定义来表述，其可以通过加窗函数来增强。
在第二实施例中，多个间隔包括恒定长度L的第一组不重叠的时间间隔，其是在不考虑脉冲簇整体的情况下选择的，以及小于L的恒定长度L1的第二组不重叠的时间间隔，其也是在不考虑脉冲簇整体的情况下选择的。L至少与脉冲的持续时间一样长，并且优选地，L1小于脉冲的持续时间。在该实施例中，复合泊松过程可以如第6节中所定义。映射可以按照第6节中定义来表述。第二实施例可以利用针对如第一实施例中的时间间隔的组定义的每组时间间隔的处理和计算。
在实施例中，使用数据驱动策略，该策略导致对核参数的接近最优的选择，这最小化了入射光子能量的估计概率密度函数的积分平方误差（ISE）。
根据本发明的第二个广义方面，提供了一种估计在辐射检测器中接收的单个辐射量子计数率的方法，该方法包括以下步骤：（1）从辐射检测器获得数字观察的时间序列，该时间序列包括与单个量子的检测相对应的脉冲；（2）根据检测器信号计算能谱敏感统计，该能谱敏感统计使用如上文关于第一广义方面所述的恒定长度L和恒定长度L1的间隔；（3）通过公式（99）确定复合泊松过程的特征函数的估计；（4）根据特征函数的估计来估计计数率。上面的步骤（4）可以通过使用优化例程或一些其他手段来拟合曲线、估计特征函数的估计的对数的DC偏移、或者将曲线拟合到特征函数的估计的对数来实现。
该申请的其余部分组织如下。第3、4和5节涉及本发明第一方面的第一实施例。第3节提供了初步背景，定义了符号，概述了数学模型，并给出了第一实施例的估计量的推导，包括修改。第4节示出了第一实施例的修改估计量对于模拟和实验数据二者的性能，并讨论了结果。第5节提供了第一实施例的结论。第6节描述了第二实施例，在相关的地方参考第一实施例。第7节描述了本发明的第二方面，一种估计计数率的新奇方法。
3第一实施例的估计量的推导
我们采用来解决堆积的一般方法基于以下策略：i）从
Figure 516657DEST_PATH_IMAGE033
获得对入射光子能量分布敏感的统计，并使用
Figure 989226DEST_PATH_IMAGE033
的观察的有限长度采样版本来估计那些统计，ii）获得从入射光子能量密度到观察的统计的统计特性的映射，iii）通过对该映射取逆来估计入射光子能量的密度。第3.1节描述了我们对统计的选择。第3.2节认为，这些统计（近似）具有与复合泊松过程相同的分布。第3.3节介绍了用于从这些统计中恢复能谱的分解技术。它基于[18]中的分解算法，但进一步发展以获得在积分平方误差方面的接近最优性能。我们处理堆积问题的方法遵循的一般主题是找到一些对基础能谱敏感的
Figure 912795DEST_PATH_IMAGE033
统计，从
Figure 838026DEST_PATH_IMAGE033
的有限时间采样版本中估计这些统计，然后在给定基础能谱的情况下对描述这些统计的统计特性的映射取逆，从而产生能谱的估计。
3.1统计的选择
我们希望从（2）中给出的观察信号获得光子能量的估计。在典型的现代光谱系统中，检测器输出
Figure 62334DEST_PATH_IMAGE033
由ADC均匀采样。不失一般性，我们假设算法可用的原始观察是
Figure 186148DEST_PATH_IMAGE034
：由于标识单个脉冲可能困难，因此我们取而代之寻找包含零个或更多个脉冲簇的固定长度
Figure 470499DEST_PATH_IMAGE035
的间隔。精确地，我们将这些间隔定义为
Figure 883025DEST_PATH_IMAGE036
，其中
Figure 911024DEST_PATH_IMAGE037
。
这里，
Figure 905656DEST_PATH_IMAGE038
被选择作为能量估计中的误差和创建间隔的概率之间的折衷。
Figure 95329DEST_PATH_IMAGE038
的值应该足够小，以确保在每个间隔内到达的总光子能量的估计误差可接受地低，但是相对于噪声方差又足够大，以确保获得大量的间隔。尽管随着计数率趋于无穷大，将观察的数据分区为间隔的概率接近于零，但与基于单个脉冲的堆积拒绝策略相比，该方法在更高的计数率下陷入瘫痪，因为准许多个光子在每个间隔内堆积。第4.2节描述了真实数据的L和
Figure 995152DEST_PATH_IMAGE038
的选择。每个间隔包含未知的、随机数量的脉冲，并且可以包含零个脉冲。
我们使用采样的原始观察来估计间隔
Figure 623580DEST_PATH_IMAGE039
中的总光子能量
Figure 925248DEST_PATH_IMAGE040
。由于每个脉冲下的面积与（1）中定义的光子能量
Figure 285822DEST_PATH_IMAGE041
成比例，我们让
Figure 485990DEST_PATH_IMAGE042
。
假设每个间隔
Figure 855792DEST_PATH_IMAGE039
中的光子到达数量、每个到达光子的能量和检测器输出噪声是随机的，并且与其他间隔无关。对于具有指数衰减的脉冲形状，在间隔中到达的少量光子能量可以在下一个间隔中被记录。泄漏量与
Figure 746387DEST_PATH_IMAGE038
成比例并且对于足够小的
Figure 605759DEST_PATH_IMAGE038
可忽略。因此，估计
Figure 214595DEST_PATH_IMAGE043
可以被视为弱相关的平稳过程的实现，其中每个估计根据随机变量X相同地分布。该关系在图1中使用典型的脉冲形状针对无噪声情况进行了说明。
3.2利用复合泊松过程的近似
在本小节中，我们关于
Figure 388087DEST_PATH_IMAGE044
描述了X的分布。然后，我们将在第3.3节对此进行取逆，以获得密度
Figure 398769DEST_PATH_IMAGE044
的估计量。使用（9）、（2）、（1）和
Figure 179774DEST_PATH_IMAGE045
是因果的事实，我们有
Figure 275906DEST_PATH_IMAGE046
。
如下证明，这简化为
Figure 987510DEST_PATH_IMAGE047
其中
Figure 915015DEST_PATH_IMAGE048
并且
Figure 53872DEST_PATH_IMAGE049
。
Figure 637300DEST_PATH_IMAGE051
和
Figure 697135DEST_PATH_IMAGE052
二者都是随机变量的i.i.d.序列。我们用Y和Z表示它们的分布。Z的分布完全从
Figure 682409DEST_PATH_IMAGE053
的分布确定，假设
Figure 992168DEST_PATH_IMAGE053
的分布为零均值具有已知方差
Figure 125209DEST_PATH_IMAGE054
的高斯分布。此外，Y是复合泊松过程，因为总和中的项数（长度为L的间隔中的光子到达数）具有泊松统计。等式（11）-（13）证明如下。（10）的第一项表示来自早期间隔的泄漏，并且近似为零。对于高斯噪声，通过关于
Figure 178615DEST_PATH_IMAGE055
执行泰勒展开，这容易地示出，
Figure 18395DEST_PATH_IMAGE056
Figure 499055DEST_PATH_IMAGE057
。
因此，属于先前间隔的一些能量将被包括在当前估计中的概率是有限的但很小的。实际上，该贡献与足够小
Figure 604546DEST_PATH_IMAGE038
的噪声相当。（10）中的第三项是零，因为
Figure 461643DEST_PATH_IMAGE045
是因果的。（10）中的第二项可以写成
Figure 155930DEST_PATH_IMAGE058
其中，我们假设脉冲形状
Figure 869808DEST_PATH_IMAGE045
足够平滑，使得
Figure 915124DEST_PATH_IMAGE059
Figure 310333DEST_PATH_IMAGE060
。它近似于所有在间隔
Figure 672176DEST_PATH_IMAGE061
内到达的光子的总能量。让
Figure 760218DEST_PATH_IMAGE062
指定在间隔
Figure 27251DEST_PATH_IMAGE061
中光子到达的数量。我们假设
Figure 226151DEST_PATH_IMAGE062
是带有比率参数
Figure 691767DEST_PATH_IMAGE063
的齐次泊松过程的实现，其中
Figure 950710DEST_PATH_IMAGE063
用长度为L的每个间隔的光子预期量来表述。此后，我们将假设（11）确切地成立，并写道
Figure 705040DEST_PATH_IMAGE064
。
最后，我们将
Figure 255101DEST_PATH_IMAGE065
写成
Figure 778486DEST_PATH_IMAGE067
其中，我们假设Z具有已知方差
Figure 208331DEST_PATH_IMAGE054
。在本小节中，我们使用复合泊松过程对第3.1节的统计进行建模。这允许我们根据可观察的量推导密度
Figure 246694DEST_PATH_IMAGE068
的估计量。在间隔
Figure 787396DEST_PATH_IMAGE061
内到达的光子数是我们指定
Figure 430867DEST_PATH_IMAGE062
的泊松随机变量。间隔Y中的总能量可以建模为复合泊松过程，即，
Figure 576153DEST_PATH_IMAGE070
其中，
Figure 305075DEST_PATH_IMAGE071
是间隔中第一个光子到达时间的索引，到达时间假设是有序的，并且表示光子能量的
Figure 383889DEST_PATH_IMAGE072
是具有密度函数
Figure 944184DEST_PATH_IMAGE068
的随机变量A的独立实现。
Figure 450251DEST_PATH_IMAGE073
形成了带有比率参数
Figure 666469DEST_PATH_IMAGE063
的齐次泊松过程。泊松率
Figure 548975DEST_PATH_IMAGE063
用长度为L的每个间隔的光子的预期数量来表述
Figure 714508DEST_PATH_IMAGE074
。
Y的实现和采样检测器响应之间的关系在图1中图示。通过将（2）代入（9）中，可以用
Figure 657056DEST_PATH_IMAGE075
来近似观察的
Figure 94991DEST_PATH_IMAGE076
，
Figure 577924DEST_PATH_IMAGE077
其中
Figure 50494DEST_PATH_IMAGE078
是不可观察随机变量Y的实现，其表示离散时间检测器响应的间隔中的光子能量，
Figure 898364DEST_PATH_IMAGE080
其中
Figure 636644DEST_PATH_IMAGE081
是Z的实现，其是独立的随机变量，表示采样过程和估计
Figure 860952DEST_PATH_IMAGE083
中的误差。我们假设Z具有已知方差
Figure 188029DEST_PATH_IMAGE054
。在X和Y的这些定义下，在有限长度的检测器输出中可以找到的间隔数是随机变量N。在高计数率下，该方法陷入瘫痪，因为能够将观察数据分区为间隔的概率接近于零。与基于堆积拒绝的策略相比，瘫痪的发作出现在更高的计数率下，因为准许多个光子在每个间隔内堆积。假设（3）-（6）中定义的时间序列已被均匀采样。不失一般性，假设单位采样间隔从
Figure 472379DEST_PATH_IMAGE084
开始，即，
Figure 681644DEST_PATH_IMAGE085
，
Figure 444063DEST_PATH_IMAGE086
。让R为离散时间随机过程，表示（1）的采样检测器响应。让
Figure 625646DEST_PATH_IMAGE087
是离散时间随机过程，其分量
Figure 159527DEST_PATH_IMAGE088
表示在固定时间间隔期间到达的光子总能量。复合泊松过程可以用于对Y建模，即，
Figure 59350DEST_PATH_IMAGE089
其中，
Figure 625460DEST_PATH_IMAGE062
是独立的泊松随机变量，并且
Figure 723866DEST_PATH_IMAGE090
是具有密度函数
Figure 350020DEST_PATH_IMAGE068
的独立同分布随机变量。
Figure 737139DEST_PATH_IMAGE073
形成带有比率参数
Figure 651480DEST_PATH_IMAGE063
的齐次泊松过程。过程Y是不可直接观察的。假设脉冲形状
Figure 807655DEST_PATH_IMAGE091
具有有限的支持。让
Figure 604710DEST_PATH_IMAGE092
为集合A的指示符函数。让脉冲长度
Figure 275863DEST_PATH_IMAGE093
由
Figure 449355DEST_PATH_IMAGE095
给出。让
Figure 460036DEST_PATH_IMAGE096
是离散时间随机过程，其表示由（2）给出的观察的检测器输出。它由被噪声过程W破坏的检测器响应R组成。不失一般性，我们假设单位采样间隔。根据观察S，我们形成了过程X，其中
Figure 427992DEST_PATH_IMAGE097
并且其中，
Figure 71595DEST_PATH_IMAGE099
是来自已知方差
Figure 517619DEST_PATH_IMAGE054
的独立噪声过程的随机变量。当我们让（1）中的脉冲形状
Figure 648386DEST_PATH_IMAGE100
Figure 115140DEST_PATH_IMAGE101
时，获得了简单的测试理论模型，在这种情况下，我们让
Figure 698568DEST_PATH_IMAGE102
，并且N只是样本长度K。从S获得
Figure 948284DEST_PATH_IMAGE103
对于真实数据来说更复杂。在那种情况下，我们将过程S分区为长度为L的不重叠块，其中
Figure 746607DEST_PATH_IMAGE104
。泊松率L用每个块的光子来表述。选择每个块的开始
Figure 56365DEST_PATH_IMAGE105
，使得任何脉冲的总能量完全包含在它到达的块内
Figure 127089DEST_PATH_IMAGE106
。
图1示出
Figure 180496DEST_PATH_IMAGE107
。我们让
Figure 817014DEST_PATH_IMAGE108
Figure 297674DEST_PATH_IMAGE109
其中
Figure 590115DEST_PATH_IMAGE110
是
Figure 525841DEST_PATH_IMAGE112
的估计。第4.2节描述了真实数据的L和
Figure 220127DEST_PATH_IMAGE038
的选择。在
Figure 871688DEST_PATH_IMAGE114
的该定义下，对于给定的样本长度K，Y中的分量数变成了随机变量。在高计数率下，该方法陷入瘫痪，因为能够创建块的概率接近于零。与基于堆积拒绝的策略相比，瘫痪的发作出现在更高的计数率下，因为准许多个光子堆积在每个区块内。让
Figure 713743DEST_PATH_IMAGE116
是离散时间随机过程，其分量
Figure 108952DEST_PATH_IMAGE088
由下式给出
Figure 657745DEST_PATH_IMAGE117
其中
Figure 745787DEST_PATH_IMAGE118
是常数，其被选择使得h和d是接近于零的小阈值。因此，随机变量
Figure 88519DEST_PATH_IMAGE088
表示在长度为L的固定时间间隔期间到达的总光子能量。d的值确保与光子到达相关联的信号在每个间隔的开始和结束时非常小。这在图1中图示。复合泊松过程可以用于对Y建模，即，
Figure 21840DEST_PATH_IMAGE119
Figure 690718DEST_PATH_IMAGE121
其中
Figure 11978DEST_PATH_IMAGE122
是带有比率参数
Figure 500728DEST_PATH_IMAGE063
的齐次泊松过程，并且
Figure 503320DEST_PATH_IMAGE123
是具有密度函数
Figure 574175DEST_PATH_IMAGE068
的独立同分布随机变量。让
Figure 738440DEST_PATH_IMAGE124
是离散时间随机过程，其表示由（2）给出的采样检测器输出。它由被噪声过程W破坏的检测器响应R组成。过程Y是不可直接观察的。使用（2）、（25）和（32），我们通过过程
Figure 980065DEST_PATH_IMAGE125
对观察建模，即，
Figure 848664DEST_PATH_IMAGE126
其中Z是已知方差
Figure 226556DEST_PATH_IMAGE054
的噪声过程。在对给定观察
Figure 827302DEST_PATH_IMAGE127
的建模中涉及的所有随机变量
Figure 290644DEST_PATH_IMAGE128
都被假设为独立的。让
Figure 182508DEST_PATH_IMAGE129
是N个独立同分布的观察。让
Figure 680485DEST_PATH_IMAGE130
成为
Figure 186553DEST_PATH_IMAGE131
的收集。让对应的特征函数为
Figure 465087DEST_PATH_IMAGE132
。
3.3估计量的基本形式
我们寻求对从光子能量A的分布到X的分布的映射取逆。我们的策略是首先根据
Figure 347593DEST_PATH_IMAGE133
获得X的特征函数，然后假设计数率和噪声特征已知，对该映射取逆。让
Figure 700077DEST_PATH_IMAGE132
是
Figure 455674DEST_PATH_IMAGE130
的特征函数。众所周知[15]对于带有比率
Figure 893609DEST_PATH_IMAGE063
的复合泊松过程Y，
Figure 579805DEST_PATH_IMAGE135
并且由于X = Y + Z，所以
Figure 583533DEST_PATH_IMAGE137
。
给定观察
Figure 696983DEST_PATH_IMAGE138
，我们可以形成X的特征函数的经验估计
Figure 622213DEST_PATH_IMAGE140
。将此视为真实的特征函数，我们可以对（40）、（41）取逆以获得A的特征函数，然后进行傅立叶变换，以找到幅度谱
Figure 945658DEST_PATH_IMAGE133
。具体而言，使用（40）、（41）并利用Z是高斯分布的假设来确保
Figure 272734DEST_PATH_IMAGE141
将是非零，
Figure 557085DEST_PATH_IMAGE142
，我们让
Figure 704032DEST_PATH_IMAGE143
是由下式描述的曲线
Figure 528769DEST_PATH_IMAGE144
。
暂时假设
Figure 710351DEST_PATH_IMAGE145
，在取（43）的区别对数并重新排列后，我们有
Figure 165604DEST_PATH_IMAGE146
。
理想情况下，通过进行傅立叶变换恢复
Figure 878476DEST_PATH_IMAGE133
Figure 444586DEST_PATH_IMAGE147
。
我们提出的估计量的基本形式在（88）中给出，并经由一系列步骤从（45）导出。首先，根据数据对
Figure 480675DEST_PATH_IMAGE148
进行估计（步骤1）。在（42）中简单地用该估计量代替
Figure 434725DEST_PATH_IMAGE148
不产生
Figure 556265DEST_PATH_IMAGE149
的ISE最佳估计。近似ISE是从
Figure 926066DEST_PATH_IMAGE148
的误差分布的近似估计获得的（步骤2）。然后，我们确定合理的加窗函数
Figure 629711DEST_PATH_IMAGE150
（在步骤3中），并通过下式对
Figure 426766DEST_PATH_IMAGE149
进行估计
Figure 35602DEST_PATH_IMAGE151
。
加窗函数
Figure 943515DEST_PATH_IMAGE150
被设计为最小化
Figure 16513DEST_PATH_IMAGE133
与基于（44）、（45）和（46）我们对
Figure 984469DEST_PATH_IMAGE133
的估计之间的近似ISE，但是其中（44）中的
Figure 80601DEST_PATH_IMAGE149
被（46）代替。类似的思想被用于从（44）对
Figure 605254DEST_PATH_IMAGE152
进行估计：（在步骤4中）找到加权函数
Figure 736021DEST_PATH_IMAGE153
，使得用下式替换（45）中的
Figure 874879DEST_PATH_IMAGE152
Figure 255045DEST_PATH_IMAGE154
产生比使用未加权估计
Figure 770339DEST_PATH_IMAGE155
更好的f A 估计。最后，加权函数
Figure 755613DEST_PATH_IMAGE156
被修改（在步骤5中）以计及（45）中的积分实际上必须被有限总和替换。以下小节对这五个步骤进行展开。
3.4估计
Figure 65372DEST_PATH_IMAGE157
需要
Figure 946215DEST_PATH_IMAGE159
的估计以估计
Figure 999622DEST_PATH_IMAGE160
。在该小节中，我们定义了直方图模型，并描述了我们基于
Figure 839402DEST_PATH_IMAGE161
值的直方图对
Figure 116800DEST_PATH_IMAGE163
的估计。假设从有限长度的数据样本中获得了N个间隔（和对应的
Figure 674820DEST_PATH_IMAGE138
值）。尽管经验特征函数
Figure 266338DEST_PATH_IMAGE164
提供了特征函数的一致的、渐近正态的估计量[21]，但它具有的缺点是随着数据点数N和所需评估点数
Figure 39253DEST_PATH_IMAGE165
增加，计算负担快速增长。取而代之，我们使用基于直方图的估计量，它具有较低的计算负担。假设观察的X值的直方图由
Figure 690815DEST_PATH_IMAGE166
向量
Figure 470552DEST_PATH_IMAGE167
表示，其中第m个仓中的计数由下式给出
Figure 928078DEST_PATH_IMAGE168
。
直方图的所有仓具有相等的宽度。仓宽度是按
Figure 742450DEST_PATH_IMAGE138
值的量值来选择的。由于选择不同的仓宽度的效果简单地等同于缩放
Figure 564913DEST_PATH_IMAGE138
值，因此在不失一般性的情况下我们假设仓宽度是统一的。在非负数据值和负数据值之间平均分配仓。直方图仓的数量2M以各种方式影响估计量，如在后面的小节中讨论的。目前，假设2M足够大以确保直方图包括所有
Figure 644995DEST_PATH_IMAGE138
值就足够了。我们通过形成缩放
Figure 843895DEST_PATH_IMAGE169
值的直方图对
Figure 512774DEST_PATH_IMAGE171
进行估计，并且取得逆离散傅立叶变换，即，
Figure 506138DEST_PATH_IMAGE173
。
这是经验特征函数的接近的近似，但是其中
Figure 322784DEST_PATH_IMAGE138
项被舍入到最近的直方图仓中心（并且u缩小到1/
Figure 325375DEST_PATH_IMAGE175
倍）。项
Figure 848761DEST_PATH_IMAGE176
只是计数具有相同值的舍入项的数量。清楚地，可以使用快速傅立叶变换（FFT）在某些离散点
Figure 826075DEST_PATH_IMAGE178
Figure 67700DEST_PATH_IMAGE180
高效地评估该函数。
3.5
Figure 608403DEST_PATH_IMAGE182
的误差分布
在（46）和（47）中滤波器
Figure 48612DEST_PATH_IMAGE183
和
Figure 649357DEST_PATH_IMAGE184
的设计依赖于
Figure 112700DEST_PATH_IMAGE186
和真实特征函数之间的误差统计。在本小节中，我们定义并描述了这些误差的特征。我们假设密度函数
Figure 267213DEST_PATH_IMAGE187
足够平滑（即，
Figure 499611DEST_PATH_IMAGE188
），并且直方图仓的宽度足够小（相对于加性噪声Z的标准偏差），使得通过将
Figure 271258DEST_PATH_IMAGE138
值舍入到每个直方图仓的中心而引入的误差跨每个仓近似均匀分布，具有零均值，并且相对于Z引起的峰值扩散很小。换句话说，由
Figure 221897DEST_PATH_IMAGE138
值的装仓引起的误差源被认为是可忽略的。由于泊松计数的统计特性和每个仓中的预期计数二者都是非整数（
Figure 166719DEST_PATH_IMAGE189
），所以在任何给定直方图仓中观察的计数数量和该仓的预期计数数量之间存在差异。在我们的模型中，我们将这两个误差源组合起来，并且称之为“直方图噪声”。我们强调，该噪声不同于（11）中建模的加性噪声Z，它导致直方图中的峰值扩散。让X的实现落入第m个仓中的概率为
Figure 332252DEST_PATH_IMAGE191
让第m个仓中的归一化直方图误差
Figure 274800DEST_PATH_IMAGE192
是第m个仓中观察的计数
Figure 712735DEST_PATH_IMAGE193
和预期计数
Figure 461248DEST_PATH_IMAGE194
之间的差值相对于直方图中的总计数N，即，
Figure 668239DEST_PATH_IMAGE196
。
使用（50）、（51）和（52）我们有
Figure 781688DEST_PATH_IMAGE197
。
如果直方图被建模为泊松向量，则示出
Figure 706919DEST_PATH_IMAGE198
。
由于直方图噪声的特征可以用观察间隔的总数N来表述，因此使用有限长度的观察数据的影响可以通过将该信息结合到
Figure 744276DEST_PATH_IMAGE200
和
Figure 71352DEST_PATH_IMAGE201
的设计中来计及。
3.6估计
Figure 90124DEST_PATH_IMAGE202
获得
Figure 564967DEST_PATH_IMAGE204
后，下一个任务是估计
Figure 592966DEST_PATH_IMAGE202
。不是在（42）中用
Figure 774549DEST_PATH_IMAGE205
代替
Figure 42850DEST_PATH_IMAGE207
，我们取而代之用（46）作为估计量，这要求我们选择加窗函数
Figure 677094DEST_PATH_IMAGE208
。在这一小节中，我们试图找到接近最优的函数
Figure 243205DEST_PATH_IMAGE208
。当考虑
Figure 607190DEST_PATH_IMAGE209
中的误差分布时，导致（46）中给出的形式的最低ISE估计量的加窗函数
Figure 233343DEST_PATH_IMAGE210
为
Figure 165003DEST_PATH_IMAGE211
其中
Figure 534804DEST_PATH_IMAGE212
表示
Figure 690979DEST_PATH_IMAGE213
的实部。我们无法计算
Figure 488034DEST_PATH_IMAGE215
，因为
Figure 159186DEST_PATH_IMAGE217
未知，所以取而代之我们试图找到近似值。我们让
Figure 67100DEST_PATH_IMAGE218
。
这通过考虑函数
Figure 890830DEST_PATH_IMAGE220
和
Figure 921103DEST_PATH_IMAGE222
之间的相对误差的量值来证明，其中
Figure 17235DEST_PATH_IMAGE223
。
相对误差的量值由下式给出
Figure 541888DEST_PATH_IMAGE224
。
因为
Figure 407076DEST_PATH_IMAGE225
，我们看到（63）的右手边在
Figure DEST_PATH_IMAGE227AA
-1时最大。因此，相对误差由下式限界
Figure 421300DEST_PATH_IMAGE229
当
Figure 4728DEST_PATH_IMAGE063
很小时，或者当
Figure 520023DEST_PATH_IMAGE231
时，其证明了该近似值。此外，我们注意到以上界限相当保守。光谱系统中光子能量的分布通常可以建模为K个高斯峰值的总和，其中第k个峰值具有位置
Figure 239717DEST_PATH_IMAGE233
和尺度
Figure 611793DEST_PATH_IMAGE235
，即，
Figure 682517DEST_PATH_IMAGE237
其中
Figure 546043DEST_PATH_IMAGE239
。
因此，特征函数将具有以下形式
Figure 385823DEST_PATH_IMAGE240
即，在对于一些
Figure 928800DEST_PATH_IMAGE241
，如
Figure 486820DEST_PATH_IMAGE243
衰减的包络内的振荡。（64）给出的上界限相当保守，因为对于大多数
Figure 78339DEST_PATH_IMAGE244
值，
Figure 585674DEST_PATH_IMAGE245
。在跨能谱的大多数评估点处，近似误差将显著更小。选择
Figure 502815DEST_PATH_IMAGE246
之后，我们可以使用（46）形成对
Figure 344869DEST_PATH_IMAGE247
的估计。加窗函数减少了由有限数量的数据样本引起的直方图噪声的影响。对于大的
Figure 740078DEST_PATH_IMAGE248
值，加窗的影响可忽略，并且估计量与直接使用（42）实质上相同。然而，在下式成立的区域中
Figure 554450DEST_PATH_IMAGE250
加窗变得重要，并起作用以对我们的
Figure 189962DEST_PATH_IMAGE202
估计限界，即，使用噪声Z是高斯分布的事实（因此
Figure 456996DEST_PATH_IMAGE251
并且因此
Figure 718213DEST_PATH_IMAGE252
），并且因为
Figure 387091DEST_PATH_IMAGE253
，我们看到
Figure 193504DEST_PATH_IMAGE254
。
这确保了（47）中的区别对数的变元保持有限，即使
Figure 947834DEST_PATH_IMAGE255
。
3.7估计
Figure 684846DEST_PATH_IMAGE257
一旦获得
Figure 270548DEST_PATH_IMAGE258
，我们就继续使用（47）对
Figure 700392DEST_PATH_IMAGE260
进行估计。这需要另一个加窗函数
Figure 752137DEST_PATH_IMAGE201
。在该小节中，我们找到了接近ISE最优的用于估计
Figure 27261DEST_PATH_IMAGE262
的函数
Figure 670732DEST_PATH_IMAGE263
。我们从定义函数
Figure 68215DEST_PATH_IMAGE265
开始以便于标注
Figure 797137DEST_PATH_IMAGE266
。
当
Figure 689000DEST_PATH_IMAGE267
时，ISE最小化，其中最佳滤波器
Figure 186978DEST_PATH_IMAGE268
由下式给出
Figure 958625DEST_PATH_IMAGE270
。
再次，由于
Figure 971580DEST_PATH_IMAGE271
、
Figure 854085DEST_PATH_IMAGE272
和
Figure 206569DEST_PATH_IMAGE273
未知，我们无法通过使用（73）-（74）计算最佳滤波器。我们取而代之进行以下观察以获得ISE最佳滤波器的近似。
3.7.1初始观察
只要估计的
Figure 962167DEST_PATH_IMAGE274
保持接近
Figure 400101DEST_PATH_IMAGE275
的真实值，最佳滤波器就保持接近单位一。对于u的小值，情况将总是如此，因为
Figure 883035DEST_PATH_IMAGE277
对于小u，
Figure 355605DEST_PATH_IMAGE278
（76）。
此外，等式（73）示出，如果
Figure 469055DEST_PATH_IMAGE280
，那么
Figure 207335DEST_PATH_IMAGE281
Figure 431643DEST_PATH_IMAGE282
，因此
Figure 493139DEST_PATH_IMAGE283
。对于较大的u值，当
Figure 839807DEST_PATH_IMAGE284
的量值变得相当于或小于
Figure 252334DEST_PATH_IMAGE285
时，估计量
Figure 280333DEST_PATH_IMAGE286
由噪声主导，并且不再提供有用的
Figure 272035DEST_PATH_IMAGE287
估计。在极端情况下，
Figure 461708DEST_PATH_IMAGE289
，因此
Figure 361531DEST_PATH_IMAGE290
，并且因此
Figure 989958DEST_PATH_IMAGE292
。
窗
Figure 291627DEST_PATH_IMAGE293
应该从估计中排除这些区域，因为这样做引入的偏差将小于未滤波噪声的方差。不幸的是，在达到该边界条件之前，
Figure 652201DEST_PATH_IMAGE295
的估计可能会严重降级，因此（77）没有特别有帮助。用于检测噪声何时开始占主导地位的更有用的方法如下。
3.7.2滤波器设计函数
对（67）的进一步操纵示出，对于典型的光谱系统，
Figure 852369DEST_PATH_IMAGE296
的量值将具有以下形式
Figure 222171DEST_PATH_IMAGE297
即：根据峰值宽度
Figure 112766DEST_PATH_IMAGE298
衰减的平均分量，以及根据谱峰位置
Figure 972138DEST_PATH_IMAGE299
变化的更快速衰减的振荡分量。在设计窗
Figure 580974DEST_PATH_IMAGE301
时，我们感兴趣的是衰减
Figure 754466DEST_PATH_IMAGE303
中的区域，其中
Figure 578197DEST_PATH_IMAGE305
，即，其中信号功率小于直方图噪声，该噪声在
Figure 546153DEST_PATH_IMAGE306
的估计期间通过移除
Figure 642285DEST_PATH_IMAGE308
而得到增强。为了获得
Figure 353889DEST_PATH_IMAGE309
的估计，将低通高斯型滤波器
Figure 281393DEST_PATH_IMAGE310
与
Figure 420251DEST_PATH_IMAGE311
进行卷积，以衰减除
Figure 3679DEST_PATH_IMAGE312
的缓慢变化的大尺度特征之外的所有。我们表示该
Figure 66444DEST_PATH_IMAGE314
Figure 51717DEST_PATH_IMAGE316
。
我们看到
Figure 361476DEST_PATH_IMAGE318
具有带尺度参数
Figure 494517DEST_PATH_IMAGE319
的瑞利分布。因此
Figure 547924DEST_PATH_IMAGE321
。
众所周知，瑞利分布随机变量
Figure 387704DEST_PATH_IMAGE322
的累积分布函数由下式给出
Figure 678483DEST_PATH_IMAGE324
因此，为了帮助计算窗
Figure 970924DEST_PATH_IMAGE326
，我们将利用函数
Figure 828022DEST_PATH_IMAGE327
以控制
Figure 522309DEST_PATH_IMAGE328
的形状。该函数
Figure 236187DEST_PATH_IMAGE329
提供我们可以对估计
Figure 281503DEST_PATH_IMAGE330
包含的信号能量比噪声能量更多多么有信心的指示。（84）中的近似源于这样一个事实，即
Figure 676712DEST_PATH_IMAGE332
也是受噪声
Figure 38555DEST_PATH_IMAGE038
轻微影响的随机变量。有时——特别是对于较大的
Figure 126596DEST_PATH_IMAGE334
值——直方图噪声可能导致足够大的
Figure 393630DEST_PATH_IMAGE336
值，从而给出虚假的自信心，并潜在地允许噪声结果破坏
Figure 654847DEST_PATH_IMAGE338
的估计。为了克服该问题，该函数被修改为在u中为单峰的
Figure 58146DEST_PATH_IMAGE339
。
该修改在如下假设上被证明，即假设高斯噪声导致
Figure 317089DEST_PATH_IMAGE340
在
Figure 71419DEST_PATH_IMAGE341
中的减小。因此，我们预期
Figure 621480DEST_PATH_IMAGE342
将在
Figure 144865DEST_PATH_IMAGE343
中增加。如果我们忽略由于
Figure 637026DEST_PATH_IMAGE344
中的峰值位置而导致的
Figure 613072DEST_PATH_IMAGE345
中的局部振荡，则由平滑后
Figure 966825DEST_PATH_IMAGE346
所近似的包络将不在
Figure 610295DEST_PATH_IMAGE347
中增加。等式（74）指示最佳窗具有
Figure 945462DEST_PATH_IMAGE348
Figure 736700DEST_PATH_IMAGE350
的形式，因此总体窗形状将在
Figure 815515DEST_PATH_IMAGE352
中减小。因此，如果在某个
Figure 147049DEST_PATH_IMAGE354
区域（信噪比高的地方）中的估计特征函数已经确定窗值应该是
Figure 653117DEST_PATH_IMAGE355
，那么拒绝在区域
Figure 869335DEST_PATH_IMAGE356
（信噪比将更差的地方）中
Figure 814157DEST_PATH_IMAGE358
的建议是合理的。使用如下知识：
Figure 166641DEST_PATH_IMAGE359
对于小
Figure 922238DEST_PATH_IMAGE360
应该接近单位一，对于大
Figure 360173DEST_PATH_IMAGE360
接近0，并且应该随着信噪比降低而“滚降”——我们考虑两个潜在的加窗函数作为对
Figure 843107DEST_PATH_IMAGE361
的近似。
3.7.3矩形窗
指示符函数提供了非常简单的加窗函数
Figure 315677DEST_PATH_IMAGE362
。
阈值
Figure 163547DEST_PATH_IMAGE364
确定了截止出现的点，并且可以按照期望手动选择（例如
Figure 901827DEST_PATH_IMAGE366
）。一旦选择了阈值，无论入射光谱中的峰值位置如何，估计量都表现出相似的ISE性能。不需要用户取决于入射光谱选择窗宽度（Gugushvili [18]提出一种用于一般分解问题的非参数估计量，其中使用矩形加窗方案。它要求手动选择窗宽度，其随着
Figure 126135DEST_PATH_IMAGE367
变化），而是由数据通路（data via）
Figure 515528DEST_PATH_IMAGE369
自动选择窗宽度。虽然简单是矩形窗的主要优点，但突变区域为最佳滤波器的滚降区域提供了差的模型。第二个滤波器形状试图改进这一点。
3.7.4 logistic窗
基于logistic函数的窗试图对更平滑的滚降建模。它由下式给出
Figure 799879DEST_PATH_IMAGE371
其中
Figure 212406DEST_PATH_IMAGE372
再次充当接受信号能量大于估计
Figure 787874DEST_PATH_IMAGE373
中的噪声能量的假设的阈值。阈值区域附近的滤波器滚降比率由
Figure 969457DEST_PATH_IMAGE375
控制。这提供了比矩形窗更平滑的过渡区域，减少了
Figure 487026DEST_PATH_IMAGE367
的最终估计中的吉布斯振荡。一旦再次地，尽管参数
Figure 386849DEST_PATH_IMAGE377
是手动选择的，但它们对
Figure 513812DEST_PATH_IMAGE367
的依赖性小得多，并且可以用于提供对多种入射光谱的接近最佳滤波。使用的典型值是
Figure 612218DEST_PATH_IMAGE379
。矩形和logistic窗函数的性能在第4节中进行了比较。
3.8估计
Figure 238371DEST_PATH_IMAGE381
在设计了窗函数
Figure 438539DEST_PATH_IMAGE382
并且因此设计了估计量
Figure 542762DEST_PATH_IMAGE384
之后，最终的任务是通过逆傅立叶变换来估计
Figure 761253DEST_PATH_IMAGE386
。本小节描述了数字实现中出现的几个问题。首先，在整个实线上对
Figure 558308DEST_PATH_IMAGE387
和
Figure 977264DEST_PATH_IMAGE388
进行数值评估是不可行的。取而代之的是，我们在有限间隔上的离散点处估计它。有限间隔被选择得足够大，使得由于排除了间隔之外的信号值而导致相当小的误差。这对于
Figure 150756DEST_PATH_IMAGE390
作为高斯混合被证明，因为对于一些
Figure 161437DEST_PATH_IMAGE392
，
Figure 191710DEST_PATH_IMAGE393
和
Figure 100891DEST_PATH_IMAGE394
的量值将如
Figure 546916DEST_PATH_IMAGE395
衰减。快速傅立叶变换（FFT）用于在离散点处评估
Figure 740000DEST_PATH_IMAGE396
，并且因此也确定了评估
Figure 878857DEST_PATH_IMAGE397
和
Figure 275335DEST_PATH_IMAGE398
的点。同样，FFT用于在离散点处评估
Figure 525051DEST_PATH_IMAGE399
的最终估计。为了使用FFT，间隔之外的信号应该足够小，以减少混叠的影响。评估点也需要足够密集，以避免评估
Figure 510324DEST_PATH_IMAGE400
时的任何“相位卷绕”模糊。这两个目标都可以通过增加直方图中的仓数量2 M（零填充）直到获得足够多的仓数量来实现。随着M增加，
Figure 882400DEST_PATH_IMAGE401
的采样密度增加，这允许检测和管理相位卷绕。更大的M也允许混叠（由
Figure 953124DEST_PATH_IMAGE402
的高斯形尾部引起）可忽略。典型地，M的值被选择为二的最小幂，其足够大使得直方图的非零值被限制到“下半部分”索引，即，
Figure 819580DEST_PATH_IMAGE404
。其次，（47）中的区别对数是未定义的，如果
Figure 659360DEST_PATH_IMAGE405
的话。根据数据对
Figure 140020DEST_PATH_IMAGE406
进行估计时，估计将为零存在很小但不是零的概率。在这种情况下，（47）中的区别对数是未定义的，并且该技术失败。随着
Figure 494778DEST_PATH_IMAGE407
增加，
Figure 351875DEST_PATH_IMAGE408
减少并且可能接近
Figure 856281DEST_PATH_IMAGE409
。当
Figure 507843DEST_PATH_IMAGE410
和
Figure 615476DEST_PATH_IMAGE411
具有相似的量值时，
Figure 10685DEST_PATH_IMAGE412
（并且因此
Figure 559478DEST_PATH_IMAGE413
）接近于零的概率可能变得显著。滤波器
Figure 460569DEST_PATH_IMAGE414
应该滚降得快于
Figure 789919DEST_PATH_IMAGE415
接近
Figure 723240DEST_PATH_IMAGE416
，以减少这对估计可能具有的影响。理想情况下，在噪声可能导致
Figure 205168DEST_PATH_IMAGE417
接近零的区域中
Figure 464111DEST_PATH_IMAGE418
应该为零。Gugushvili已经示出[18]，对于矩形窗，随着数据集长度增加
Figure 15179DEST_PATH_IMAGE420
，取逆失败的概率接近零。
3.9离散符号
我们暂时离题以引入附加符号。在整篇文章的其余部分，粗体字体将用于指示对应于命名函数的离散采样版本的2M×1向量，例如，
Figure 17770DEST_PATH_IMAGE421
表示2m×1向量，其值由在点
Figure 541155DEST_PATH_IMAGE422
Figure 518469DEST_PATH_IMAGE423
处评估的特征函数
Figure 822412DEST_PATH_IMAGE424
给出。方括号符号[k]用于索引向量中的特定元素，例如，
Figure 173234DEST_PATH_IMAGE425
具有
Figure 551126DEST_PATH_IMAGE426
值。我们还使用负索引来访问向量的元素，所用的方式类似于python编程语言。负索引应该相对于向量的长度来解释，即
Figure 214188DEST_PATH_IMAGE427
指代向量中的最后一个元素（其等同于
Figure 677531DEST_PATH_IMAGE428
）。
3.10估计量概述
我们使用的估算过程可以在以下步骤中总结。
1.使用（8）将采样的时间序列分区为间隔。
2.根据（9）计算每个间隔的
Figure 834974DEST_PATH_IMAGE138
值。
3.根据
Figure 332951DEST_PATH_IMAGE138
值生成直方图
Figure 901336DEST_PATH_IMAGE429
。
4.使用逆FFT计算
Figure 117553DEST_PATH_IMAGE396
，以高效地在各个采样点处评估（50）。
5.在适当的点处计算
Figure 59DEST_PATH_IMAGE430
和
Figure 165592DEST_PATH_IMAGE431
。
6.经由（46）使用
Figure 842561DEST_PATH_IMAGE432
和
Figure 280495DEST_PATH_IMAGE433
计算
Figure 29009DEST_PATH_IMAGE434
。
7.计算
Figure 235999DEST_PATH_IMAGE435
，
Figure 349449DEST_PATH_IMAGE436
的低通滤波版本。
8.经由（83）和（85）计算
Figure 274679DEST_PATH_IMAGE437
。
9.使用
Figure 312037DEST_PATH_IMAGE438
和（86）或（87）计算
Figure 639113DEST_PATH_IMAGE439
。
10. 经由（47）使用
Figure 923463DEST_PATH_IMAGE440
和
Figure 132728DEST_PATH_IMAGE439
计算
Figure 160727DEST_PATH_IMAGE441
。如果
Figure 342309DEST_PATH_IMAGE442
的任何元素为零，并且
Figure 607681DEST_PATH_IMAGE439
的对应元素非零，则估计失败，因为区别对数未定义。
11.根据下式使用
Figure 507504DEST_PATH_IMAGE443
的FFT计算
Figure 73615DEST_PATH_IMAGE444
Figure 172021DEST_PATH_IMAGE445
。
3.11性能度量
使用误差的积分平方（ISE）来度量估计量的性能。ISE度量估计密度的全局拟合。
Figure 798174DEST_PATH_IMAGE446
离散ISE度量由下式给出
Figure 919714DEST_PATH_IMAGE447
其中
Figure 289515DEST_PATH_IMAGE448
是一个2M×1向量，其元素包含每个直方图仓区域中的概率质量，即
Figure 258739DEST_PATH_IMAGE449
。
该向量
Figure 55794DEST_PATH_IMAGE450
表示对应的估计概率质量向量。
4第一实施例的数值结果
使用模拟和真实数据执行实验。
4.1模拟
Trigano等人[11]使用的理想密度用于这些模拟。它由六个高斯分布和一个伽马分布的混合组成，以模拟康普顿背景。混合密度由下式给出
Figure 664630DEST_PATH_IMAGE451
其中
Figure 634860DEST_PATH_IMAGE452
是具有均值
Figure 645541DEST_PATH_IMAGE453
和方差
Figure 613497DEST_PATH_IMAGE454
的正态分布的密度。伽马分布的密度由
Figure 522679DEST_PATH_IMAGE456
给出。在8192个等间距的整数点对密度进行采样，以产生概率质量的离散向量
Figure 234283DEST_PATH_IMAGE457
。进行FFT以获得
Figure 365050DEST_PATH_IMAGE458
，
Figure 566224DEST_PATH_IMAGE459
值的采样向量。
为实验选择了特定的计数率
Figure 884073DEST_PATH_IMAGE063
，对应于每个观察间隔的预期事件数。预期堆积密度经由（40）获得，即，离散向量
Figure 399368DEST_PATH_IMAGE458
缩放了
Figure 384641DEST_PATH_IMAGE063
，取幂，然后缩放了
Figure 507449DEST_PATH_IMAGE461
，并且最后应用FFT
Figure 312594DEST_PATH_IMAGE462
。
等式（93）与高斯分布卷积，以模拟抹去观察的光谱的噪声Z的效果
Figure 366001DEST_PATH_IMAGE463
。
这表示观察光谱的预期密度，包括堆积和加性噪声。使用根据（94）分布的随机变量创建观察直方图。实验通过对
Figure 268098DEST_PATH_IMAGE464
而参数化，其中
Figure 748758DEST_PATH_IMAGE466
并且
Figure 306778DEST_PATH_IMAGE468
。对于每个参数对
Figure 708416DEST_PATH_IMAGE464
，制作一千个观察直方图。概率质量向量
Figure 668282DEST_PATH_IMAGE469
的估计是使用（88）进行的，其中将（86）和（87）二者用于
Figure 319843DEST_PATH_IMAGE471
。对于两种窗形状使用
Figure 161897DEST_PATH_IMAGE473
的阈值，并且对于logistic形状使用
Figure 557106DEST_PATH_IMAGE475
的阈值。记录每个估计
Figure 371478DEST_PATH_IMAGE477
和真实向量
Figure 6990DEST_PATH_IMAGE469
之间误差的离散ISE度量。为了与渐近带宽结果进行比较，使用矩形窗进行估计，该矩形窗的带宽是根据Gugushvili在[18]中规定的条件1.3——即
Figure 274024DEST_PATH_IMAGE479
，其中
Figure 472924DEST_PATH_IMAGE481
——选择的。我们强调Gugushvili的滤波器的
Figure 141802DEST_PATH_IMAGE483
不要与（87）的
Figure 197483DEST_PATH_IMAGE485
混淆。渐近带宽标准是通过使用下式来实现的
Figure 951813DEST_PATH_IMAGE487
其中
Figure 954404DEST_PATH_IMAGE488
试验了Gugushvilli的
Figure 290838DEST_PATH_IMAGE489
的三个值，即
Figure 455103DEST_PATH_IMAGE491
。
还使用矩形滤波器（95）进行估计，该滤波器具有各种值
Figure 696729DEST_PATH_IMAGE492
的固定带宽。最后，根据（1）创建时间序列数据，具有理想的矩形脉冲形状和107个脉冲，其能量根据（92）分布。选择脉冲长度和计数率以给出泊松率
Figure 299748DEST_PATH_IMAGE493
。Trigano等人[11]描述的算法用于从包含32×1024（持续时间×能量）个仓的二维直方图中估计基础幅度密度——据报道，该仓的选择给出了最佳的准确性和合理的执行时间。记录核心算法的性能和处理时间，以便与我们提出的算法进行比较。图2描绘了由数据驱动的logistic形滤波器对具有参数对（
Figure 677640DEST_PATH_IMAGE495
）的实验做出的典型估计
Figure 278386DEST_PATH_IMAGE497
。还绘制了真实向量
Figure 554777DEST_PATH_IMAGE499
（细实线）和观察的直方图
Figure 899171DEST_PATH_IMAGE501
（包含一些噪声的下曲线）。在观察的直方图中可以清楚地看到堆积峰。尽管估计的密度受到振铃（由于吉布斯现象）影响，但它以其他方式估计真实密度，并校正观察的直方图中存在的堆积。图3描绘了在与图2相同的操作点进行的典型估计，但是估计量具有矩形滤波器，其中使用（96）和
Figure 131569DEST_PATH_IMAGE503
选择带宽。这对应于图6中的操作区域，其中固定带宽滤波器（
Figure 903216DEST_PATH_IMAGE504
）的性能接近数据驱动滤波器的性能。明显，在也校正堆积的同时，结果所得的估计包含更多的噪声。图4示出了使用矩形滤波器和各种固定带宽时，ISE度量的分布密度作为样本计数的函数。在分布平均值（MISE）之间绘制了线，以帮助可视化。数据驱动的矩形滤波器（86）的结果也用较粗的曲线绘制、连接。这清楚地图示了固定带宽滤波的弱点。对于任何固定的带宽，ISE随着样本数增加而减小，最终随着偏差成为主要误差源而渐近。在这一点（其是噪声和带宽相关的）上，尽管样本计数增加，但ISE仍然保持大体不变。固定带宽排除了在最终计算中使用一些估计
Figure 916172DEST_PATH_IMAGE506
，即使在它们具有高信噪比（SNR）时。图4还示出了我们提出的数据驱动带宽选择的矩形滤波器给出的结果。该曲线位于每个固定带宽曲线的拐点附近。这指示跨样本计数范围，为数据驱动矩形滤波器选择的带宽接近最佳带宽值（对于矩形滤波器）。图5-图7示出了在三个计数率
Figure 798677DEST_PATH_IMAGE507
下，ISE度量的分布密度作为每个直方图中估计总数N的函数。logistic和矩形滤波器的MISE曲线低于使用（96）给出的带宽对于大部分应用感兴趣区域获得的曲线。存在各种区域，其中非数据驱动带宽（
Figure DEST_PATH_IMAGE508
）给出与数据驱动带宽的相似性能，然而跨整个采样计数范围，这种性能并未得到维持。logistic滤波器形状具有略好于矩形滤波器形状的性能，尽管两个滤波器之间的差异对于ISE度量来说相对较小。表1比较了提出的算法和最近在[11]中描述的算法之间的结果。两种方法的ISE在测试（
Figure DEST_PATH_IMAGE510
）下的操作点处是相似的，然而我们提出的算法需要相当少的计算。
表1：与[11]中描述的算法的比较
Figure 961281DEST_PATH_IMAGE511
4.2真实数据
估计量被应用于真实数据，以评估其在实际应用中的有用性。在（8）中发现的阈值
Figure 903829DEST_PATH_IMAGE038
被选择为加性噪声
Figure DEST_PATH_IMAGE512
的标准偏差的近似一半。这确保了创建间隔的合理高概率，又确保了间隔能量估计的误差低。间隔长度L的值被选择为典型脉冲“长度”的近似四倍，即间隔长度
Figure 404080DEST_PATH_IMAGE513
的四倍。从锰样本获得能量直方图，其中光子通量率标称为每秒105个事件。观察的直方图主峰的形状中存在轻微的负偏斜，表明复杂的噪声源影响了系统。这在图8中几乎不可见。噪声被建模为双峰高斯混合，而不是单一高斯峰。使用简单的最小二乘优化例程来拟合双峰高斯参数
Figure DEST_PATH_IMAGE514
到噪声峰位于仓索引零周围。手动选择了
Figure 903326DEST_PATH_IMAGE063
的合适值。使用具有数据驱动带宽的logistic滤波器来估计真实密度。图8示出了观察的和估计的概率质量向量的曲线图。主峰（仓450 ~ 600）增强，而堆积虽未完全消除，但已衰减。一阶堆积峰已经减少。峰堆积比（主峰高度与第一个堆积峰高度之比）从大约6增加到大约120。这些改进与其他最先进的系统（例如[11]）相当。估计量未能完全解决堆积存在几个可能的原因。估计量的准确度取决于正确建模高斯噪声峰。
双峰高斯混合建模了噪声峰，使得最大误差小于噪声密度峰的1%。给定估计谱中的残余堆积峰低于主峰的1%，估计量对噪声建模中误差的敏感性可能在某些部分对此有贡献。未解决的堆积的第二个原因可能是由于观察光谱估计中的不确定性。几个残余堆积峰相对靠近观察的直方图的底部。残余峰可能仅是估计量的噪声引起的假象。最后，数学模型可能是观察的光谱的过于简单的近似。检测过程包括模型中没有包括的许多二阶影响（例如，弹道亏损、供给电荷耗尽、相关噪声、非线性等...）。这些较小的影响可能限制堆积校正估计量的准确度。
5第一实施例的概述
我们采用了Gugushvili [18]提出的用于在高斯噪声下进行分解的估计量，并适配其以用于校正X射线光谱中的脉冲堆积。我们提出了易于实现的数据驱动带宽选择机制，并跨光谱应用的感兴趣的广泛样本计数（104~109计数）范围提供了ISE/MISE的显著降低。数据驱动的矩形带宽选择接近最优（对于矩形滤波器），并且在感兴趣的范围内优于基于渐近结果或固定带宽的带宽选择。
尽管最初的结果看起来很有希望，但是需要进一步的工作来改进实际实现的性能。该估计仍然包含与吉布斯现象相关联的“振铃”假象。附加的滤波器形状试图减少这一点，还存在更接近MSE最佳的其他形状。
6第二实施例
本节给出第二实施例的能谱估计量的概述。第二实施例解决了第一实施例的要求在每个间隔中近似包含整个簇的问题。在第二实施例中，如果期望，则可以使用整个数据序列，并且通过引入2个不同的间隔长度L和L1来补偿重叠。
我们需要包括一些在第一实施例中没有提到的附加术语。特别是
Figure DEST_PATH_IMAGE516
。能谱估计量基于
Figure 110316DEST_PATH_IMAGE517
。
滤波器
Figure DEST_PATH_IMAGE518
的引入允许我们解决出现的几个实现问题。我们使用的估计过程可以在以下步骤中总结。
1.将采样时间序列分区为固定长度的间隔
Figure DEST_PATH_IMAGE520
。
2.根据
Figure DEST_PATH_IMAGE522
计算每个间隔的
Figure 286083DEST_PATH_IMAGE138
值。
3.从
Figure 24363DEST_PATH_IMAGE138
值生成直方图
Figure 248671DEST_PATH_IMAGE523
。
4.使用
Figure 575747DEST_PATH_IMAGE523
的逆FFT计算
Figure 656835DEST_PATH_IMAGE396
。
5.将采样的时间序列分区为具有长度L1的不同间隔集，并跟着类似的计算以获得
Figure DEST_PATH_IMAGE524
。
6.计算
Figure 69362DEST_PATH_IMAGE525
和
Figure DEST_PATH_IMAGE526
。
7.使用
Figure 910410DEST_PATH_IMAGE528
和
Figure 91993DEST_PATH_IMAGE530
计算
Figure DEST_PATH_IMAGE531
。
Figure DEST_PATH_IMAGE533
8.计算
Figure 609562DEST_PATH_IMAGE534
，
Figure DEST_PATH_IMAGE536
的低通滤波版本。
9.计算
Figure 77363DEST_PATH_IMAGE537
。
10.使用
Figure DEST_PATH_IMAGE538
计算
Figure 643473DEST_PATH_IMAGE539
。
11.使用
Figure 7458DEST_PATH_IMAGE540
和
Figure 633612DEST_PATH_IMAGE539
计算
Figure DEST_PATH_IMAGE541
。如果
Figure 568201DEST_PATH_IMAGE542
的任何元素为零，并且
Figure 938002DEST_PATH_IMAGE539
的对应元素非零，则估计失败，因为区别对数未定义。
12.使用
Figure DEST_PATH_IMAGE543
的FFT根据下式计算
Figure 156494DEST_PATH_IMAGE544
Figure 953549DEST_PATH_IMAGE546
。
6.1算法细节
将检测器输出流分区为长度L的不重叠间隔集，即
Figure DEST_PATH_IMAGE547
。让
Figure 562385DEST_PATH_IMAGE138
是第j个间隔中检测器输出样本的总和，即，
Figure 283347DEST_PATH_IMAGE548
。
假设L大于脉冲长度，第j个间隔可能包含“完整”脉冲以及被间隔末端截断的脉冲。可以示出，
Figure 294028DEST_PATH_IMAGE138
由我们用
Figure DEST_PATH_IMAGE549
表示的“完整”脉冲的能量、我们用
Figure 324301DEST_PATH_IMAGE550
表示的截断脉冲和噪声
Figure DEST_PATH_IMAGE551
的能量的叠加组成。
让检测器输出流分区为第二不重叠间隔集
Figure 420433DEST_PATH_IMAGE552
，其中
Figure DEST_PATH_IMAGE553
。让
Figure 945087DEST_PATH_IMAGE554
由下式给出
Figure DEST_PATH_IMAGE555
。
如果L1被选择为略小于脉冲长度，则
Figure 872591DEST_PATH_IMAGE556
项将不包含“完整”的脉冲，而是仅由截断脉冲
Figure DEST_PATH_IMAGE557
和噪声
Figure 11449DEST_PATH_IMAGE551
的能量的叠加组成。任何间隔内的截断脉冲数都具有泊松分布。我们有
Figure 594877DEST_PATH_IMAGE558
。
我们可以将间隔
Figure 920291DEST_PATH_IMAGE560
中的总能量分解为来自被截断的脉冲的能量贡献
Figure 639986DEST_PATH_IMAGE562
和来自完全包含在间隔
Figure 949744DEST_PATH_IMAGE564
中的脉冲的能量贡献
Figure 82785DEST_PATH_IMAGE566
，即，
Figure DEST_PATH_IMAGE567
其中，
Figure DEST_PATH_IMAGE569
表示脉冲完全包含在间隔（长度为
Figure DEST_PATH_IMAGE571
）中的区域中的噪声，并且
Figure DEST_PATH_IMAGE573
表示脉冲被截断（长度为
Figure 949241DEST_PATH_IMAGE574
）的区域中的噪声。因此，
Figure 851338DEST_PATH_IMAGE576
。
通过将（103）和（105）组合，我们有
Figure DEST_PATH_IMAGE577
。
重新排列给出
Figure 331998DEST_PATH_IMAGE578
。
我们可以用类似于我们估计
Figure DEST_PATH_IMAGE579
的方式或一些其他方法来估计
Figure 703068DEST_PATH_IMAGE580
，例如，经由经验特征函数或通过对
Figure DEST_PATH_IMAGE581
值的归一化直方图执行FFT。
当执行分解操作时，减小的间隔长度
Figure DEST_PATH_IMAGE583
的泊松率
Figure DEST_PATH_IMAGE585
用于计及在其上发生复合泊松过程
Figure DEST_PATH_IMAGE587
的子间隔。
6.2内部数量的可视化
为了帮助读者理解，图9绘制了估计过程期间获得的各种数量。上方的蓝色曲线（在仓零处的值约为0.3）绘制了
Figure 356903DEST_PATH_IMAGE588
，观察光谱的估计特征函数。棕色曲线用于示出
Figure DEST_PATH_IMAGE589
的真实值，该值在区域[6000，10000]中作为具有周期性空值的下曲线清晰可见。数量
Figure DEST_PATH_IMAGE591
示出为透明的红色，并表现为“噪声”，其平均密度在
Figure DEST_PATH_IMAGE593
仓周围达到峰值。
Figure DEST_PATH_IMAGE595
的预期值用黑色虚线示出。这是如下获得的：使用（75）、已知的
Figure 926556DEST_PATH_IMAGE063
值，并假设具有已知的
Figure 843696DEST_PATH_IMAGE596
以获得
Figure DEST_PATH_IMAGE597
的高斯噪声。数量
Figure 420171DEST_PATH_IMAGE598
用透明的蓝色曲线示出。这几乎不可见，因为它在间隔[0，4000]，[12000，16000]中与
Figure DEST_PATH_IMAGE599
紧密重合，并且在间隔[5000，11000]中与
Figure 815380DEST_PATH_IMAGE600
紧密重合。注意，
Figure DEST_PATH_IMAGE601
的颜色看起来在间隔[5000，11000]内从红色变为紫色，因为两个透明绘图重叠。实心黑线示出
Figure 705451DEST_PATH_IMAGE602
，它是
Figure DEST_PATH_IMAGE603
的低通滤波版本。如本节开头关于平滑的段落中所述，低通滤波移除了由于峰位置引起的
Figure 340963DEST_PATH_IMAGE604
中的任何本地振荡。项
Figure 607997DEST_PATH_IMAGE606
用作
Figure DEST_PATH_IMAGE607
的估计。可以看出，在
Figure 869214DEST_PATH_IMAGE608
的区域中，
Figure DEST_PATH_IMAGE609
提供了
Figure 351142DEST_PATH_IMAGE610
合理良好估计。当这两个数量彼此接近时，估计的质量下降，直到最终由噪声主导。滤波器
Figure DEST_PATH_IMAGE611
应包括
Figure 344505DEST_PATH_IMAGE612
的良好估计，同时排除差的估计。为了找到获得
Figure 895572DEST_PATH_IMAGE612
的良好估计的区域，我们解决了如下问题：给定
Figure DEST_PATH_IMAGE613
，局部区域中
Figure 632584DEST_PATH_IMAGE614
的计算值主要从噪声出现的概率是多少。
7计数率估计
先前的估计量假设
Figure 234598DEST_PATH_IMAGE063
是已知的。在没有先验知识的情况下，可以如下获得对
Figure 664442DEST_PATH_IMAGE063
的估计。
1.使用来自前一节的
Figure DEST_PATH_IMAGE615
，计算
Figure 968385DEST_PATH_IMAGE616
。
2.使用
Figure DEST_PATH_IMAGE617
，估计计数率。这可以用多种方式完成。
3.一种方式是使用优化例程或一些其他手段将曲线拟合为
Figure 243508DEST_PATH_IMAGE618
。拟合的参数可以用于获得对计数率的估计。
4.另一种方式涉及估计
Figure DEST_PATH_IMAGE619
的DC偏移。这可以通过对合适数量的
Figure 683717DEST_PATH_IMAGE620
点进行平均来完成。通过前一节中的由
Figure DEST_PATH_IMAGE621
滤波获得的点通常是合适的，尽管较少的点也可以产生充足的估计。
5.另一种方式涉及使用优化引擎或一些其他手段来将曲线拟合到
Figure 18883DEST_PATH_IMAGE622
。拟合
Figure DEST_PATH_IMAGE623
的合适的参数化曲线由下式给出
Figure 823504DEST_PATH_IMAGE624
其中
Figure DEST_PATH_IMAGE625
并且其中
Figure 699056DEST_PATH_IMAGE626
被选择以允许曲线足够准确地拟合。参数
Figure 197033DEST_PATH_IMAGE063
提供了对计数率的估计。优化引擎不需要对
Figure 968680DEST_PATH_IMAGE620
中的每个点给予相等的加权。
8各图的描述
以下各图有助于理解该过程。图1示出了一种用于对检测器输出分区的可能方案。该图示描绘了采样检测器对三个入射光子的响应。为了使图更清晰，噪声的影响已经被移除。输出响应被分区成几个相等长度（L）的区域。在每个区域中到达的脉冲数量对于处理系统来说是未知的。在第一个间隔中一个脉冲到达。在第二个间隔中两个脉冲到达。在第三个间隔中没有脉冲到达。在每个间隔内到达的总光子能量被计算为感兴趣的统计，是每个间隔内所有样本值的总和。间隔与脉冲到达在时间上不对齐。图2图示了估计过程的输出。入射光子能量的真实概率密度被绘制成实心黑线。光子到达率是这样的，使得在任何给定的间隔
Figure DEST_PATH_IMAGE627
内平均有三个光子到达。检测器输出信号
Figure 997947DEST_PATH_IMAGE628
中加性噪声的标准偏差等于一个直方图仓宽度。收集了一百万个间隔。绘制了每个间隔中总能量的直方图。这用蓝线绘制出。堆积的影响清楚明显，特别是在仓75、150和225周围。红色轨迹绘制了在系统处理数据后真实入射能谱的估计量。尽管估算中出现了一些噪声，但堆积的影响已经被移除。预期该估计平均正确恢复真实入射光谱。该结果是使用内部滤波器获得的，其带宽是根据数据自动确定的。图3图示了在相同操作条件下与图2相同的数量，然而在这种情况下，内部滤波器的带宽使用来自文献的渐近结果来确定。尽管估计的入射能量概率密度已经恢复，但与图2相比，方差显著更大。图8图示了系统对真实数据的操作。蓝色轨迹绘制观察的能量值的概率密度，而红色轨迹绘制入射光子能量的估计真实概率密度。没有黑色轨迹，因为真实的概率密度未知。在该实验中，锰样本的X射线荧光被用作光子源。光子到达率约为每秒105个光子。间隔长度被选择为使得光子之间的平均时间对应于两个间隔的长度。收集并分区足够的数据以形成5.9xl06的间隔。加性噪声的标准偏差对应于4.7个直方图仓。估计过程清楚地减少了堆积峰，并且增强了真实峰。图9图示了在第二实施例中描述的系统的模拟期间获得的各种量。在第5.1节“内部量的可视化”中对它进行了描述。图9-12涉及第二实施例。图10图示了该实验的观察的和真实的输入光子能量的概率密度，从其导出图9-13。黑色轨迹绘制了真实的概率密度。红色轨迹绘制了在给定的间隔长度内平均有三个光子到达时预期观察的密度。蓝色轨迹绘制了实际观察的密度。在观察的密度中可以看到高达十数量级的堆积。图10包括几个产生于典型光谱系统的绘图。实际入射光子密度（“理想密度”）用实心黑线绘制。通过分区时间序列数据获得的观察直方图以深蓝色示出。脉冲堆积导致的能谱失真明显。图13使用对数垂直轴绘制了各种内部量。向绘图中心倾斜的深蓝色曲线是
Figure DEST_PATH_IMAGE629
。水平穿过绘图的绿色量是
Figure 880453DEST_PATH_IMAGE630
。向绘图中心倾斜的上青色曲线是
Figure DEST_PATH_IMAGE631
。图11图示了曲线
Figure 45986DEST_PATH_IMAGE202
在复平面中的轨迹。图12图示了类似于图9的内部量，然而还存在一些附加的信号。主要是噪声的水平红色轨迹和对应的黑色虚线表示
Figure 988534DEST_PATH_IMAGE632
，直方图噪声的特征函数的量值。在图中心形成“噪声峰”的透明绿色绘图是估计
Figure DEST_PATH_IMAGE633
。该量在图9中以蓝色绘制，并且几乎不可见，因为它被图12中未示出的
Figure 223206DEST_PATH_IMAGE634
遮挡住了。平均值为-3的水平轨迹是
Figure DEST_PATH_IMAGE635
的绘图。青色轨迹是
Figure 722452DEST_PATH_IMAGE636
，加性高斯噪声的特征函数的量值，它从零仓处的零值开始，并且在8000仓周围倾斜到最小值。图13涉及计算计数率的第二个广义方面。它图示了在计算
Figure DEST_PATH_IMAGE637
时使用的内部量。青色轨迹是
Figure 195022DEST_PATH_IMAGE638
，加性高斯噪声的特征函数的量值，其从仓零处的零值开始，并在8000仓周围倾斜到最小值。在图中心倾斜到最小值的深蓝色轨迹是
Figure DEST_PATH_IMAGE639
，它是观察数据的特征函数的估计。平均值为-3的黄色/绿色水平轨迹是
Figure 105209DEST_PATH_IMAGE640
的估计。
参考文献
Figure DEST_PATH_IMAGE641
Figure 92756DEST_PATH_IMAGE642
Figure DEST_PATH_IMAGE643
Figure 317064DEST_PATH_IMAGE644
Figure DEST_PATH_IMAGE645

Claims (19)

1.一种确定在辐射检测器中接收的单个辐射量子的能谱的方法，所述方法包括以下步骤：
（1）从辐射检测器获得数字观察的时间序列，所述时间序列包括对应于单个量子检测的脉冲；
（2）根据检测器信号计算能谱敏感统计，能谱敏感统计定义从脉冲幅度密度到能谱敏感统计的映射；和
（3）通过对能谱敏感统计应用映射的逆来估计脉冲幅度的密度，从而确定能谱。
2.根据权利要求1所述的方法，进一步包括使能谱敏感统计基于多个时间间隔上的数字观察的总和。
3.根据权利要求2所述的方法，进一步包括使用近似复合泊松过程来定义所述映射。
4.根据权利要求3所述的方法，进一步包括通过建模噪声来增强近似复合泊松过程。
5.根据权利要求4所述的方法，进一步包括将所述映射表述为幅度的特征函数、能谱敏感统计和建模噪声之间的关系。
6.根据权利要求5所述的方法，进一步包括通过对数字观察的总和的直方图应用逆傅立叶变换来计算所述能谱敏感统计的特征函数。
7.根据权利要求5或权利要求6所述的方法，进一步包括用低通滤波器计算幅度的特征函数。
8.根据权利要求2至7中任一项所述的方法，进一步包括选择所述多个时间间隔中的每一个以包含零个或更多个近似完整的脉冲簇，并且将所述多个时间间隔定义为不重叠并且具有恒定长度L。
9.根据权利要求8所述的方法，进一步包括在每个时间间隔的开始和结束时要求检测器信号的最大值。
10.根据权利要求8或9所述的方法，进一步包括将近似复合泊松过程定义为每个时间间隔中脉冲幅度的总和。
11.根据权利要求8至10中任一项所述的方法，进一步包括如由等式（40）和（41）定义的表述映射。
12.根据权利要求11所述的方法，进一步包括通过加窗函数来增强所述映射。
13.根据权利要求2至7中任一项所述的方法，进一步包括选择所述多个间隔以包括：在不考虑脉冲簇的整体的情况下恒定长度为L的第一组不重叠的时间间隔；以及在也不考虑脉冲簇的整体的情况下小于L的恒定长度L1的第二组不重叠的时间间隔；其中L至少与脉冲的持续时间一样长。
14.根据权利要求13所述的方法，进一步包括选择L1为小于脉冲的持续时间。
15.根据权利要求13或权利要求14所述的方法，进一步包括如第6节中阐述的定义近似复合泊松过程。
16.根据权利要求13至15中任一项所述的方法，进一步包括如第6节中阐述的表述映射。
17.根据权利要求1至16中任一项所述的方法，进一步包括使用数据驱动策略，所述数据驱动策略被选择为导致核参数的接近最优的选择，所述核参数最小化辐射的单个量子的能量的估计概率密度函数的误差的积分平方。
18.一种估计在辐射检测器中接收的单个辐射量子的计数率的方法，所述方法包括以下步骤：
（1）从辐射检测器获得数字观察的时间序列，所述时间序列包括对应于单个量子检测的脉冲；
（2）基于多个时间间隔上的数字观察的总和，从检测器信号计算能谱敏感统计，所述能谱敏感统计使用近似复合泊松过程定义从脉冲幅度密度到能谱敏感统计的映射，所述多个时间间隔包括：恒定长度为L的第一组不重叠的时间间隔，其是在不考虑脉冲簇的整体的情况下选择的；以及小于L的恒定长度L1的第二组不重叠的时间间隔，其也是在不考虑脉冲簇的整体的情况下选择的；其中L至少与脉冲间隔的持续时间一样长；
（3）通过等式（99）确定近似复合泊松过程的特征函数的估计；
（4）根据特征函数的估计来估计计数率。
19.根据权利要求18所述的方法，进一步包括通过以下来估计计数率：使用优化例程或其他手段拟合曲线，估计特征函数的估计的对数的DC偏移，或将曲线拟合到特征函数的估计的对数。
CN202080037871.4A 2019-03-22 2020-03-23 具有脉冲堆积非参数分解的辐射检测 Pending CN113826030A (zh)

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
AU2019900974 2019-03-22
AU2019900974A AU2019900974A0 (en) 2019-03-22 Radiation detection with non-parametric decompounding of pulse pile-up
PCT/AU2020/050275 WO2020191435A1 (en) 2019-03-22 2020-03-23 Radiation detection with non-parametric decompounding of pulse pile-up

Publications (1)

Publication Number Publication Date
CN113826030A true CN113826030A (zh) 2021-12-21

Family

ID=72608371

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202080037871.4A Pending CN113826030A (zh) 2019-03-22 2020-03-23 具有脉冲堆积非参数分解的辐射检测

Country Status (8)

Country Link
US (1) US20220137111A1 (zh)
EP (1) EP3942337A4 (zh)
JP (1) JP7291239B2 (zh)
CN (1) CN113826030A (zh)
AU (1) AU2020249184B2 (zh)
CA (1) CA3134143A1 (zh)
IL (1) IL286443A (zh)
WO (1) WO2020191435A1 (zh)

Family Cites Families (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
NL7110516A (zh) * 1971-07-30 1973-02-01
US6590957B1 (en) * 2002-03-13 2003-07-08 William K. Warburton Method and apparatus for producing spectra corrected for deadtime losses in spectroscopy systems operating under variable input rate conditions
FR2870603B1 (fr) * 2004-05-19 2006-07-28 Commissariat Energie Atomique Mesure et traitement d'un signal comprenant des empilements d'impulsions elementaires
WO2006029475A1 (en) * 2004-09-16 2006-03-23 Southern Innovation International Pty Ltd Method and apparatus for resolving individual signals in detector output data.
US10304217B2 (en) * 2012-07-30 2019-05-28 Toshiba Medical Systems Corporation Method and system for generating image using filtered backprojection with noise weighting and or prior in
DE102013213362A1 (de) * 2013-07-08 2015-01-08 Fraunhofer-Gesellschaft zur Förderung der angewandten Forschung e.V. Verfahren zur Identifizierung und Quantifizierung von emittierenden Teilchen in Systemen
US9801595B2 (en) * 2014-09-08 2017-10-31 Toshiba Medical Systems Corporation Count-weighted least squares parameter estimation for a photon-counting detector
US10621756B2 (en) 2017-01-13 2020-04-14 Canon Medical Systems Corporation Apparatus and method for correcting bias in low-count computed tomography projection data
JP6912304B2 (ja) * 2017-07-20 2021-08-04 株式会社日立製作所 波高頻度分布取得装置、波高頻度分布取得方法、波高頻度分布取得プログラム及び放射線撮像装置

Also Published As

Publication number Publication date
WO2020191435A1 (en) 2020-10-01
IL286443A (en) 2021-10-31
JP2022528512A (ja) 2022-06-14
EP3942337A1 (en) 2022-01-26
AU2020249184B2 (en) 2021-09-16
AU2020249184A1 (en) 2021-06-24
JP7291239B2 (ja) 2023-06-14
US20220137111A1 (en) 2022-05-05
EP3942337A4 (en) 2022-12-21
CA3134143A1 (en) 2020-10-01

Similar Documents

Publication Publication Date Title
Van Hateren et al. On the genetically meaningful decomposition of grain-size distributions: A comparison of different end-member modelling algorithms
JP5111105B2 (ja) メインインパルスの連鎖を含んだ信号処理を含む測定方法および分析装置
JP5022902B2 (ja) 検出器出力データ内の個別信号分離装置および方法
US20060157655A1 (en) System and method for detecting hazardous materials
Ansari et al. Correction of highly noisy strong motion records using a modified wavelet de-noising method
WO2015085372A1 (en) Method and apparatus for resolving signals in data
CN110542406B (zh) 基于emd-mpf改进的陀螺仪信号去噪方法
Sepulcre et al. Sparse Regression Algorithm for Activity Estimation in $\gamma $ Spectrometry
EP3637150A1 (en) Gamma-ray spectrum classification
EP2923220B1 (en) Method of spectral data detection and manipulation
EP2433155B1 (en) Radiation detection
CN113826030A (zh) 具有脉冲堆积非参数分解的辐射检测
Lavielle et al. A simulated annealing version of the EM algorithm for non-Gaussian deconvolution
Mclean et al. Non-parametric decompounding of pulse pile-up under gaussian noise with finite data sets
RU2802542C2 (ru) Детектирование излучения с помощью непараметрического разложения наложенных импульсов
US8676744B2 (en) Physics-based, Bayesian sequential detection method and system for radioactive contraband
Trigano et al. Statistical pileup correction method for HPGe detectors
Unser et al. Maximum likelihood estimation of liner signal parameters for poisson processes
Guo et al. Structure-oriented filtering for seismic images using nonlocal median filter
Lopatin et al. Pileup attenuation for spectroscopic signals using a sparse reconstruction
Candy et al. Threat detection of radioactive contraband incorporating compton scattering physics: A model-based processing approach
Fu et al. The application of wavelet denoising in material discrimination system
CN116541668B (zh) 一种游泳划水次数确定方法、装置、设备及存储介质
EP0990216A1 (en) Method of analyzing capabilities of multiple-suppression computer seismic data processing software
Neuer et al. A Cognitive Filter to Automatically Determine Scintillation Detector Materials and to Control Their Spectroscopic Resolution During Temperature Changes

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination