东华理工核数据处理课程设计 - 图文

核数据处理课程设计

题 目：

关于γ能谱寻峰、边界道的确定以及峰面积计算

学生姓名：

班 级：

学 号：

指导教师：吴 和 喜

二零一二年六月

核数据处理课程设计报告

--γ能谱寻峰、边界道以及峰面积计算

摘要：本文简述两种不同的平滑过程与寻峰方法，简单比较法寻峰方法和导数法寻峰的对

比分析，试用采用八阶高斯函数拟合后公式求导后对寻峰以及峰面积结果分析。

关键词：能谱、峰面积、导数法、寻峰、边界道、拟合

Abstract:This method of derivative spectroscopy data after smoothing peak search, and

compare order derivative France to find the peak derivative method for the spectroscopy data the peak search and compare different fitting formula derivation on the peak search and peak areaanalysis of results.

Keywords: spectroscopy, peak area, the derivative method,peak search, boundary Road, fitting

一、实验原理： 1、采用高斯函数的平滑方法

如果把谱数据中的统计涨落看成是“白噪声”，当使用匹配滤波器进行滤波时，可以得到最佳的信号噪声比。所谓匹配滤波器，就是该滤波器在道域中冲击响应函数与峰函数互为镜象。一般情况下，谱中的峰函数可以近似为高斯函数。由于高斯函数是偶函数，所以匹配滤波器在道域中的响应函数也应该是高斯函数。

2、简单比较法寻峰

简单比较法是一种最直观又快速的寻峰方法。其基本思想是：在谱数据中，若某一道

的数据比其他相邻的几道大很多，则可认为该道存在一个峰。连续检索光滑后的谱数据，如果在第m道同时满足下列关系式

ym?i?ym?kym?ym?i

则第m道附近有一个峰。式中，k是寻峰阈值,通常取1-1.5：m=2i+1为寻峰窗口宽度，一般i去1或2。在第m道附近的谱数据中用二阶差值多项式计算出精确峰位。

简单比较法适于寻找强单峰，在高本底上导找弱峰和分辨重峰的能力都比较差。但算法简单，程序运行速度较快。

导数是一条光沿曲线上数值的微商或微商函数,在数学上称它为曲线上 各 点 的 斜

率。在导数光谱术中导数吸收光谱是指光强度或吸光强度对波长 的变化率曲线。当将原吸收曲线进 行一阶、二阶直至四 阶求导时,便可得到各阶导数光谱。导数光谱的基本特征(1)对原曲线的极值求一、二阶导数时通过原点,求二阶导数时为极小值,四阶导数为极大值,并恢复至原位置。据此特征可以准确的测定最大的吸收峰位。（2）极值数目随求导阶数的增加而增加。如当原曲线有一个极大值时 ,求n阶导数就有n+1个极大值和极小值,而原曲线有两个极大值时,则出现n+3个极值。(3)随求导阶数的增加,谱带变窄和复杂化。两组分组合的吸收曲线原光小的吸收谱带为大的吸收谱带所掩蔽,经二阶求导之后,原曲线中的小 “卫星”峰就变成极小值,四阶导数变为极大值。这时改善了谱线的分辨率,导数光谱更加明显。

在γ射线衍射分析中各衍射谱也有类似的特征,不过这时的衍射曲线是晶体衍射强度 4 相对于衍射角5的一条光滑曲线,通常它是高斯、柯西和修正的洛仑兹函数形式,或其复合形式的分布函数。其一阶导数曲线通过零点 ,对应于衍射线的峰位3而二阶导数分布曲线的极小值也对应于峰位。由于一阶导数的零点确定峰位比较困难,特别是在多相衍射线重叠的情况下,一阶导数定峰不如二阶导数明显、准确,因此一般采用经光滑处理求二阶、四阶以至偶数阶导数的定峰方法。这是一种简便而有效的精确测定衍射峰位的方法。导数方法可以准确测定衍射峰位和重叠衍射峰的数目。因此,我们也称它为导数定峰法。所谓导数定峰法就是利用一条光滑衍射曲线上各点的数值微商（或求导）方式来确定重叠衍射峰位的方法。 能谱分析对寻峰方法的基本要求如下： （1）识别弱峰 （2）剔除假峰

（3）分辨重峰

（4）不仅能计算出峰位的整数道址，还能计算出峰位的精确值，某些情况下要求峰位的误差小于0.1道。

一、导数法寻峰原理

一阶导数法：

经过平滑后的能谱我们可以很容易的计算其各阶导数。一阶导数计算过程如下： （1） 采用五点三次平滑公式的一阶导数公式： yi?，1（yi-2-8yi-1?8yi?1-yi?2） 12（2）沿道指i方向检索一阶导数从正变负， 其零点值对应的道址为峰位。 （3）再利用峰高判断条件： 0.8FWHM

这就是一阶导数寻峰程序中实际应用的峰高判定条件。

一阶导数法寻峰时最常用的一种方法，用这种方法能找出大部分峰，计算简单，运行书度也快，但主要缺点是不能分辨很近的重峰。 二阶导数法：

二阶导数和一阶导数法比较，二阶导数寻峰方法具有更高的重峰分辨能力。具体过程如下：

（1）先计算平滑数据二阶导数数值。采用五点三次平滑公式的二阶导数公式：

1,yi?（2yi-2-yi-1-2yi-yi?1?2yi?2）

7（2） 沿道指i增加方向检索二阶导数值， （3） 找出二阶导数局部负的极小值对应 的道指。如果这个负极大值的绝对值大于 改导数的标准差的若干倍，则认为该道指 是一个真峰。

（4） 为了进一步剔除假峰，采用判定条件 --\D（ymp-1）y”mp-1y\mp?1D(y\mp?1)?K2

二阶导数重叠峰导数图

对于二阶导数法，其峰边界为真正的极大值点对应的道址

三阶导数法

三阶导数计算方法其峰位为其三阶导数值从负变为正的零点所对应的道址，峰的边界为三阶导数从正变为负的零点对应的道址。

三阶导数可以确定重峰区的各个组分峰的峰位，而且能够狠灵敏的分辨重峰。 一阶导，峰位处由正到负过零点，边界处由负到正过零点 ；二阶导，峰位处负的局部最小值，边界道为正的极大值；三阶导，峰位处由负到正过零点，边界处由正到负过零点。

三、实验过程：

流程图：第一种方法，先高斯函数平滑，简单比较法寻峰

开始 获取能谱数据 高斯函数法对能谱进行平滑 得到平滑后的能谱数据 简单比较法寻峰 得到计算峰面积 结束

matlab编程程序如下：

clc; clear;

[Filename,Pathname]=uigetfile('*.txt','选择谱数据'); fid=fopen([Pathname Filename],'r'); [array,count]=fscanf(fid,'%d',inf); fclose(fid); x=1:count; plot(x,array)

peak_position=[]; %峰位置信息

%谱光滑 采用高斯函数光滑 y=0.2438 0.205 0.1218 0.05126 k1=4; for i=1:count

if(i<=k1)

for i=1:k1

array_smooth(i)=0.2438*array(i)+0.205*(array(i+1)+array(i+1))+0.1218*(array(i+2)+array(i+2))+0.05126*(array(i+3)+array(i+3)); %对称做镜像处理

end

elseif(i>k1&i<(count-k1))

array_smooth(i)=0.2438*array(i)+0.205*(array(i-1)+array(i+1))+0.1218*(array(i-2)+array(i+2))+0.05126*(array(i-3)+array(i+3)); else

array_smooth(i)=0.2438*array(i)+0.205*(array(i-1)+array(i-1))+0.1218*(array(i-2)+array(i-2))+0.05126*(array(i-3)+array(i-3)); end end

x2=1:count;

plot(x2,array_smooth)

p=ones(1,count);

m=1; %控制参数1 j=2; %控制参数2

for i=1+j:count-j-m

p(i)=(array_smooth(i)*array_smooth(i+m))/(array_smooth(i-j)*array_smooth(i+m+j)); end

k=1.5; %峰判断控制参数 for i=1:count

upper(i)=1+k/sqrt(array_smooth(i)); lower(i)=1-k/sqrt(array_smooth(i)); end%简单比较法寻峰 factor=1; for i=1+factor:count-factor

if(p(i)>upper(i)&p(i)>p(i+1)&p(i)>p(i-1)) peak_position=[peak_position i]; end

end

plot(array_smooth,'b*'); hold on;

peak=array(peak_position) plot(peak_position,peak,'r+'); 实验结果如下：

四、本文第二种思想：gauss8阶函数拟合法，此方法的主要思想是：不将数据平滑处理自然就没有了平滑公式，直接用八阶高斯函数拟合源数据，在对拟合后的函数进行求导，一阶导，峰位处由正到负过零点，边界处由负到正过零点 ；二阶导，峰位处负的局部最小值，边界道为正的极大值；三阶导，峰位处由负到正过零点，边界处由正到负过零点。其他过程都与前面的计算方法相同。

开始 获取能谱数据 高斯函数法对能谱进行拟合 对拟合后的函数求导 导数法的寻峰 得到计算峰面积 结束

五、MATLAB实验程序如下： % y 原始数据，1025*1 % x = (1:1025)';

% y2为前600个数据 也可用全部数据 y2 = y(1:600); % 这四个峰出现在y(1:600)中

% x2 = (1:600)'; x2 = (1:length(y2))';

% Set up fittype and options. ft = fittype( 'gauss8' ); opts = fitoptions( ft ); opts.Display = 'Off';

opts.Lower = [-Inf -Inf 0 -Inf -Inf 0 -Inf -Inf 0 -Inf -Inf 0 -Inf -Inf 0 -Inf -Inf 0 -Inf -Inf 0 -Inf -Inf 0];

opts.StartPoint = [1236 282 3.92319043722374 907 336 4.22548549053934 338.232193942756 342 5.30830398585762 308 230 6.26590931040574 254.442173171231 277 6.49267774242084 135.84970058518 329 7.80555875890908 112 1 21.9443863743451 102 398 10.6865668950446];

opts.Upper = [Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf Inf];

% Fit model to data.

[fitresult, gof] = fit( x2, y2, ft, opts );

% Plot fit with data. figure( 'Name', 'fit gauss8' ); h = plot( fitresult, x2, y2 );

legend( h, 'y2 vs. x2', 'fit gauss8', 'Location', 'NorthEast' ); % Label axes xlabel( 'x2' );

ylabel( 'y2' );

title('gauss8拟合的效果图'); grid on

dd111 = diff(fitresult(x2));

x20 = x2(dd111(1:end-1).*dd111(2:end)<0) %导数值跨过0的值

% x20 = [228; 238; 259; 267; 281; 308; 336; 376;

figure

plot(x2,fitresult(x2),x20,fitresult(x20),'o') title('边界道址、峰位');

xz11 = 228 - (238-228); % 第一个峰左侧的边界道址（对称点） Z1 = trapz(x2(xz11:238),y2(xz11:238)) % 第一个峰的面积

Z2 = trapz(x2(267:308),y2(267:308)) % 第二个峰的面积

Z3 = trapz(x2(308:376),y2(308:376)) % 第三个峰的面积

xz42 = 397 + (397-376); % 第四个峰右侧的边界道址（对称点） Z4 = trapz(x2(376:xz42),y2(376:xz42)) % 第四个峰的面积

六、实验结果输出：

397];

七、各阶导数寻峰效果比较：

(1)一阶导数法寻峰时最常用的一种方法，用这种方法能找出大部分峰，计算简单，运行书度也快，对孤立峰有良好分辨能力，但主要缺点是不能分辨很近的重峰

(2) 二阶导数寻峰方法具有更高的重峰分辨能力，但灵敏度差。

(3) 三阶导数可以确定重峰区的各个组分峰的峰位，有很灵敏的重峰分辨能力，但是寻峰速度慢。

八、不同拟合函数的效果比较：

采用平滑后的多项式曲线平滑公式的寻峰效果更好，由于数据跳跃性太大，采用八阶高斯函数拟合法可能剔除弱峰，无法分辨重峰，误差比较大，但是平滑本来就会丢失弱峰和重峰，所以相比平滑公式求导寻峰高斯函数拟合法也是可以的。虽然老师上课没有讲这种拟合函数的方法，但是我认为在一般的数据分析中（单个或者几个强峰的分析）效果还是可行的。

九、实验心得体会： 在本次《核数据处理课程设计》中我深刻体会到了这门课的难度，由于对MATLAB软件不熟悉，我几乎是从零学起，每一个语句都来之不易，上课的理论知识还算比较扎实，编程的思想我很清楚，我知道自己要干什么，但却无法用MATLAB语句写出来。最终，在这些天不断的摸索中我终于初步掌握了MATLAB的编程，学会了用MATLAB对核数据的一般处理，平滑（拟合曲线）—求导寻峰—确定边界道—计算峰面积等等。总之，在这次课程设计中我学到了很多，学会了将理论知识应用到实践当中去，理论知识掌握得比较扎实。同时，我也深刻认识到自己的不足,我的MATLAB编程我能力确实比较差，在以后的学习中一定会加强这方面的训练，努力掌握软件编程方法这个处理核数据最关键的步骤。 参考文献：

1、《核数据处理》 东华理工大学，吴和喜 2、BASIC语言程序库，张巨洪，1984

3、北京地区X射线分析专业会议论文集，陶琨，1983 4、V. J Hammond el al.，J. Opt. Soc. Amer.43(1953) 5、J.Steinier et al，Aual，Chem，44（1972）

本文来源：https://www.bwwdw.com/article/w6m3.html