多道γ能谱分析软件中寻峰算法比较总结.docx
《多道γ能谱分析软件中寻峰算法比较总结.docx》由会员分享，可在线阅读，更多相关《多道γ能谱分析软件中寻峰算法比较总结.docx（15页珍藏版）》请在冰豆网上搜索。
多道γ能谱分析软件中寻峰算法比较总结
自动寻峰
由于谱结构的复杂和统计涨落的影响，从谱中正确地找到全部存在的峰是比较困难的。
尤其是找到位于很高本底上的弱峰，分辨出相互靠得很近的重峰更为困难。
谱分析对寻峰方法的基本要求如下：
（1）比较高的重峰分辨能力。
能确定相互距离很近的峰的峰位。
（2）能识别弱峰，特别是位于高本底上的弱峰。
（3）假峰出现的几率要小。
（4）不仅能计算出峰位的整数道址，还能计算出峰位的精确值，某些情况下要求峰位的误差小于0.2道。
很多作者对寻峰方法进行了研究，提出了很多有效的寻峰方法。
目的：
判断有没有峰存在
确定峰位（高斯分布的数学期望），以便把峰位对应的道址，转换成能量
确定峰边界——为计算峰面积服务（峰边界道的确定，直接影响峰面积的计算）分为两个步骤：
谱变换和峰判定
要求：
支持手动/自动寻峰，参数输入，同时计算并显示峰半高宽、精确峰位、峰宽等信息，能够区分康普顿边沿和假峰
感兴区内寻峰
人工设置感兴趣大小，然后在感兴区内采用简单方法寻峰重点研究：
对感兴区内的弱峰寻峰、重峰的分解
对于一个单峰区，当峰形在峰位两侧比较对称时，可以由峰的FWHM计算峰区的左、右
边界道址。
峰区的宽度取为3FWHM，FWHM的值可以根据峰位mp由测量系统的FWHM
刻度公式计算。
由于峰形对称，左、右边界道和峰位的距离都是1.5FWHNM。
mL
INT（mp
1.5FWHM
0.5）
mR
INT（mp
1.5FWHM
0.5）
式中mp是峰位，INT的含义是取整数。
对于存在有低能尾部的峰，其峰形函数描述（参见图）。
ym
HEXP[（m
mp）2/2
2]，m≥mp－J
ym
HEXP[J（2m
2mp
J）/22]，m≤mp－J
式中H为峰高，mp为峰位，是高斯函数的标准偏差，J为接点的道址和峰位之间的距
离。
在峰位的左侧，有一个接点，其道址为mp－J。
在接点的右侧，峰函数是高斯函数。
在接点的左侧，峰函数用指数曲线来描述。
这时峰区的左、右边界道址为
m
INT（m
p
1.12FWHM2
/J0.5J0.5）
L
mR
INT（mp
1.5FWHM
0.5）
带有低能尾部的峰函数的图形
全谱自动寻峰
基于核素库法：
能量刻度完成后，根据核素库中的能量计算对应的道址，在各个道址附
近（左右10道附近）采用简单的寻峰方法（导数法）
方法：
根据仪器选择开发
IF函数法/简单比较法（适于寻找强单峰，速度快）
满足条件：
dataimdataikdataidataim可认为有峰存在
然后在datai-m至datai+m中找最大值，对应的道值即为峰位
k：
找峰阈值，根据高斯分布，一般k取值1—1.5常用5点、7点极大值法（m取2，3）判定峰是否有意义
一般，用R=N0/Nb≥R0确定峰是否有意义
R为峰谷比，R0为设定值（经验值）
N0为净峰幅度与基底之和
Nb为基底计数
intCMmcaView:
:
SearPeakCompare（intBeginch,intEndch,intm,floatk）
高斯乘积函数找峰法（可靠性差，不建议采用）
描述谱峰形状的函数主要是高斯函数
G（i）
A
exp
（i
i0）2/22则由相邻的数据点定义
2
一个新的函数（第一高斯乘积函数，只与
FWHM
2.3556
有关）：
Pm（i）
G（i）G（im1）
exp（11.092m）m
2
G（i2）G（im）
H2
m是步长（用道表示），是高斯乘积函数的阶数，则
Pm（i）称为第m阶高斯乘积函数。
找
峰的灵敏度与m有关，随m的增加灵敏度提高。
为避免基线参数的影响，最好扣除本底后，再应用高斯乘积函数找峰。
考虑统计涨落的影响，把判断无峰存在的
1变为一个“单位带”。
即峰的判断为：
1
k/
yi
无峰
Pm（i）
k/
yi
k3
（1
）有峰
峰位的确定：
由
Pm（i）过1的两点求平均来确定；峰边界的确定：
“单位带”下限的两个最
端点；半高宽的确定：
函数Pm（i）在“1上”的截距；组合峰的确定：
在乘积函数的两个峰之
间没有处于“带内”的乘积函数值
导数法（一阶、二阶、三阶）
yi'1
m
Cjyij
Nmj
m
Nm为规范化常数，Cj平滑的变换系数。
3次多项式5点光滑一阶导数公式：
（可以采用）
yi'1（yi28yi18yi1yi2）峰位确定：
一阶导数值由正变负=0处；峰边界确定：
一阶导数
12
由负变正=0处
CalculateDifferential（0,size,m,differ）;
for（intj=m;j<=size-m;j++）
{
for（inti=1;i<=m;i++）
{
if（differ[j-i]）>0&&differ[j-i]>maxtemp）{maxtemp=differ[j-i];nmax=j-i;}
if（differ[j+i]）<0&&differ[j+i]}
if（（nmin-nmax）>0.8*fwhm&&（nmin-nmax）<3*fwhm）
//FWHM参数根据仪器能量分辨率可人工确定，fwhm~20
peakposition[p++]=j+0.5;//保持峰位对应的道址
}
5点光滑二阶导数公式（软件中推荐采用）
''1
yi（2yi2yi12yiyi12yi2）
//7点二阶导数
（5*（countsdata[j-3]+countsdata[j+3]）-3*（countsdata[j-1]+countsdata[j+1]）-4*countsdata[j]）/42;
yi'
1
（22.0yi367.0yi258.0yi158.0yi167.0yi222.0yi3）
252.0
软件中推荐采用11点以上的公式
峰位确定：
二阶导数最小值对应的道址；峰边界确定：
二阶导数正极大值点
for（intj=m;j<=size-m;j++）//m~30
{
intmaxtemp=-0.5,mintemp=-0.5;
If（differ[j]<-0.05）
for（inti=1;i<=m;i++）
{
if（differ[j-i]>maxtemp）{maxtemp=differ[j-i];nmax=j-i;}
if（differ[j+i]>mintemp）{mintemp=differ[j+i];nmin=j+i;}
}
if（（nmin-nmax）>0.8*fwhm&&（nmin-nmax）<3*fwhm）
//FWHM参数根据仪器能量分辨率可人工确定，fwhm~20
peakposition[p++]=j+0.5;//保持峰位对应的道址
}
1
c0jni
Di
ni1
2ni
ni
j
j
1
j
k
k
1/2
Di
cjnij;
c2jni
j
j
k
j
k
Significanceof2ndDerivative:
S
Di
/
cj
100
p2
j2
exp
21j2
p2
p2
wherepistheassumedpeakwidth.
k:
Gouptocj
10-6
Peak
‘found’
when
S>Threshold
试验：
系列
1为处理后的原始能谱，系列
2为5点一阶导数，系列
3为5点二阶导数，系列
4为对称零面积法寻峰
只要选择好合适的寻峰阈值，足以满足准确寻找到全能峰，并剔除假峰（如康普顿边沿，反散射峰）
5点光滑三阶导数公式判定各感兴区是单峰还是重峰
yi'''1（yi22yi12yi1yi2）
2
峰位确定：
三阶导数由负变正=0处；峰边界确定：
三阶导数由正变负
=0处
判定峰是否有意义0.8FWHM≤N≤3FWHM
峰高判定条件
|ym|maxTRHympe0.5/
这个公式就是在一阶导数法寻峰程序中实际应用的峰高判定条件。
CalculateDifferential（Beginch,Endch,m,differ）;
intCMmcaView:
:
SearPeakDifferential（intBeginch,intEndch,intfwhm,intdiffer[],intm）
{
intn1=0,differ[Endch-Beginch+1],nmax=0,nmin=0,maxtemp,mintemp,temp;maxtemp=differ[0];mintemp=differ[0];
for（intj=1;j<=Endch-Beginch;j++）
{
temp=differ[j-1];
if（_copysign（temp,differ[j]）!
=differ[j-1]&&differ[j]<0）n1=j+Beginch;if（differ[j]maxtemp）{maxtemp=differ[j];nmax=j+Beginch;}
}
if（（nmin-nmax）>0.8*fwhm&&（nmin-nmax）<3*fwhm）
elsereturn（0）;
returnn1;
}
对称零面积法（推荐自动寻峰中采用，可探测弱峰和重峰）
面积为零的“窗”函数与实验谱数据进行褶积变换，且要求“窗”函数为对称函数。
对
线性基底的褶积变换将为零，只有存在峰的地方不为零。
yi'
m
m
Cjyij
Cj
0Cj
Cj
j
m
j
m
匹配滤波器法（类峰形函数）Cj
j
2
1
m
k
2
2]
exp[
2]
exp[
2m
2
2
1km
m
~
Cjdataij
yi
jm
1f
峰判定准则Ri~
yi
m
2
C2jdataij
j
m
2m+1为变换宽度，
FWHM2.3556为峰宽参数，若变换后的y'和其均方根误差的比值超
过预先给定的寻峰阈值（f），则认为找到了一个峰。
峰位的确定：
Ri的正极值对应的道址；峰边界的确定：
Ri的正峰两边相邻的两个极小值
之间的距离可以作为峰的宽度信息；半宽度：
两过零截距。
CalculateArea（0,size,m,fwhm,area,R）;
for（intj=m;j<=size-m;j++）
{if（area[j]>0&&R[j]>fh）
for（inti=1;i<=m;i++）
{
if（area[j-i]）>0&&area[j-i]0&&area[j+i]}
if（（nmin2-nmin1）>0.6*fwhm&&（nmin2-nmin1）<=2*fwhm）
peakposition[p++]=j+0.5;//保持峰位对应的道址
}
协方差法（曲线拟合寻峰，计算机寻峰中采用，可分辨重峰，比较好的寻峰方法，但计
算较为复杂，运算速度较慢）
1975年H.P.BLOK等提出了一种新的寻峰方法，称为协方差法。
用一个峰形函数与实验
谱数据逐段拟合（一个高斯形函数与实验谱yi的协方差）
yijyi'Cjbi，Cj为峰形/高斯函数CjEXP[2.773（j/H2）]H为峰FWHM，y'i为拟合峰高，bi
为本底常数（在峰区内假定不变）mjm
m
m
m
m
gj
gjCjyi
j
gjCj
gjyij
yi'
jm
jm
jm
jm
m
m
m
gj
gjCj2
（
gjCj）2
j
mj
m
j
m
m
m
m
m
yi'
gj
gjCjyij
gjCj
gjyij
f（f判峰阈值）判定是否存在峰
用Ri
j
m
jm
j
m
jm
yi'
m
m
m
m
12
gj
gj
gjCj2
（
gjCj）2
j
m
j
m
jm
j
m
Cj通常为纯峰形函数高斯函数：
Cj
exp
4ln2（
j）2
，H为峰的FWHM
H
gj为各道计数的权重因子gj
1
或
gj
exp[
2（j
H）4]
yij
yij
参数选择：
H的取值最好与实验谱峰的半宽度接近，
2m+1一般取2H左右最好，f一般取
2~5
峰位确定：
当Ri为极大值对应的道址；峰边界确定：
Ri为负极大值处对应的道址
为了更好地分辨出落在一个强峰‘肩部’上的弱峰，可以在一个峰的左半部分和右半部
分别计算Ri值，寻找相互靠得很近的组分峰。
线性拟合寻峰方法（适合于在峰区内分辨重峰）
吸取匹配滤波器方法的优点，同时用一阶导数法和线性拟合双重峰的技术来提高分辨重
峰的能力，形成了一种新的寻峰方法，称为线性拟合寻峰方法。
Deconvolutionmethod
Firstthebackgroundisremoved（ifdesired）,thenMarkovspectrumiscalculated（ifdesired）,thentheresponsefunctionisgeneratedaccordingtogivensigmaanddeconvolutioniscarriedout.
可以提供多种算法，方便自行选择
总结
1．对于弱峰，数据光滑前，高斯乘积函数法和协方差法不能使用，若先光滑再找峰，又
容易影响重叠峰的分辨；而导数法和对称零面积变换法，无论峰的统计质量如何，均可
使用。
2．从统计假峰及高基底的抑制能力及重峰的分辨能力来看，一、三阶导数法和对称零面
积变换法是较好的。
对于一、三阶导数法，可先用适当多数据点的一阶导数法找峰，选
取适当的灵敏度常数，以抑制假峰；然后用少点的三阶导数法（或用一阶导数法重复三
次）检查是否有漏峰和重峰。
对称零面积变换法同理。
3．从高基底的抑制能力和弱峰识别的准确度来看，对称零面积变换法最好。
（在计算机
自动找峰程序中，最好采用对称零面积变换法。
）
参考资料
4.对找到的峰进行净面积判定是降低假峰出现几率的有效方法。
当峰的净面积比峰的总面
积（峰的净面积和本底面积之和）的标准偏差大若干倍时，才确认该峰是一个真峰，否
则认为它是假峰，予以剔除。
峰的判弃主要是利用峰面积来进行判定真假峰。
对于给定的灵敏因子S，若峰的净面积为NetAREA，峰的宽度为Width。
这些参数满足下式认为峰有意义，应保留，否则将找到的此峰丢弃。
此式为：
NetAREAWidth
S（AREA/Width）12
S越大灵敏度越高，一般情况下S=3。
参考文献
SPECTRAN-FVersion2ListingsVolume.4CommonSubroutinesandDatastructures,CANBERRAIndustries,Inc.1980能谱的数据处理（原文）.doc
M.A.Mariscotti,Nucl.Instr.Meth.50（1967）309MariscottiAlgorithmmodifiedbyRouttiandPrussin:
J.T.RouttiandS.G.Prussin,Nucl.Instr.Meth.72（1969）125
升级会员 