《粒子物理与核物理试验中的数据分析报告》由会员分享，可在线阅读，更多相关《粒子物理与核物理试验中的数据分析报告（40页珍藏版）》请在金锄头文库上搜索。
1、2019-11-08,1,粒子物理与核物理实验中的数据分析,杨振伟 清华大学 第四讲:蒙特卡罗方法,2019-11-08,2,上一讲回顾,概率的基本概念,随机变量与概率密度函数,随机变量的平均值与方差,能不通过实验对随机变量进行研究吗？,2019-11-08,3,本讲要点,蒙特卡罗方法 随机数产生子 任意分布抽样之函数变换法与舍选法 蒙特卡罗方法中的精度问题 在粒子物理与核物理中的应用,2019-11-08,4,蒙特卡罗方法简介,蒙特卡罗方法就是利用一系列随机数来计算各种概率大小和随机变量均值等等的数值分析技术。通常的步骤为：,产生一系列在0，1之间均匀分布的随机数 。 利用这些随机数按某些概
2、率密度函数 抽样生成我们感兴趣的另一随机序列 。 利用这些 值来估计 的一些特性，例如：通过找到在区间 的 比例，给出积分值 。,第一层面上的应用：,蒙特卡罗计算 = 积分,第二层面上的应用：,蒙特卡罗变量 = “模拟的数据”,2019-11-08,5,随机数的产生,用物理方法产生 真正的随机数,不可重复 产生速度慢,用数学方法产生 伪随机数,可以重复 产生的速度快,2019-11-08,6,真随机数与伪随机数,美国兰德（RAND）公司在1950年代，利用真空管中产生的噪音制作了一个含十万个真正的随机数表，并运用于其开展的所有模拟研究中。,真正的随机数与伪随机数之间的区别在于：数据串是否具有可
3、压缩性，即能否用更短的形式来表示。,真正的随机数是不可压缩的，非常不规则，以至于无法用更短的形式来表示它。,在粒子物理与核物理研究中，随机数的可重复性经常也是非常有用的，尤其是程序的调试(debugging)。,2019-11-08,7,随机数产生子,目的是使在 0，1 范围内产生的伪随机数满足：,均匀性；相互独立性；长周期性,乘同余法,友情推荐,2019-11-08,8,CERN库的随机数产生子,PAW用户, gRandom-SetSeed(); Float_t random = gRandom-Rndm(1); , Real random(1) Call Rmarin(ISEED,0,0)
4、 Call Ranmar(random,1) ,注意: 用于产生子的 随机数种子还可以用 来保证后续进程的随 机数不重复。,Root 用户,粒子物理与核物理研究中，大都采用CERN程序库提供的随机数产生子。,2019-11-08,9,随机数均匀性与相关性检验,subroutine mc double precision lamda,M,x,x0,y call hbook1(10,r,100,0.,1.,0.) call hbook2(20,r(i+1) vs. r(i), &100,0.,1.,100,0.,1.,0.) x0=1. lamda=1220703125 ! 5*13 M=4294
5、967296. ! 2*32 do i=1,10000 x=Mod(lamda*x0,M) y=x/M call hfill(10,real(y),0.,1.0) if(i.gt.1)call & hfill(20,real(y_old),real(y),1.0) x0=x y_old=y end do return end,随机变量,第 I 个随机变量,第 I+1 个随机变量,频数,均匀性,相关性,2019-11-08,10,随机数均匀性与相关性检验,随机变量,第 I 个随机变量,第 I+1 个随机变量,频数,均匀性,相关性,void random() UInt_t lambda,M,x0;
6、 TH1F *h1 = new TH1F(“h1“,“,100,0,1); TH2F *h2 = new TH2F(“h2“,“,100,0,1,100,0,1); lambda=1220703125; /513 M =4294967296; /232 x0 =1; double y, y_old; for (int i=0; iFill(y); if (i1) h2-Fill(y_old,y); x0=x; y_old=y; ,2019-11-08,11,用蒙特卡罗法计算积分,对于计算积分值,解析解：,数值解：,函数必须解析可积,自变量不能太多,对函数是否解析可积和 是否太多自变量无要求,在
7、AB区间均匀投总数为N个点。,2019-11-08,12,蒙特卡罗方法中的精度问题,采用蒙特卡罗方法(MC)计算积分 与传统的梯形法相比有如下特点,一维积分:,多维积分:,对于维数大于4的积分, 用蒙特卡罗方计算积分总是最好。,2019-11-08,13,从均匀分布到任意分布的随机数,函数变换法,舍选法,寻找某个函数，当函数的自变量取均匀分布值时，对应的函数值自动满足给定分布。,从一个随机变量与对应概率密度函数最大值构成的二维均 匀分布中，按概率密度函数与自变量关系曲线切割得到。,2019-11-08,14,函数变换法,均匀分布,任意分布,2019-11-08,15,例子:指数分布抽样,抽样效
8、率为100%。,可采用函数变换法抽样的分布,指数分布 三维各向同性分布 二维随机角度的正、余弦分布 高斯分布 n 个自由度的 2 分布 伽马分布 二项式分布 泊松分布 Student 分布 （http:/pdg.lbl.gov/2008/reviews/monterpp.pdf）,2019-11-08,16,2019-11-08,17,舍选法,问题: 如何找到函数的最大值?,2019-11-08,18,舍选法举例,subroutine acc_rej real rvec(1) call hbook1(10,x(r),100,0.,10.,0.) call hbook1(20,x(r),100,
9、0.,10.,0.) call hbook2(30,f(x) vs. x(r),100,0.,10.,100,0.,1.1,0.) fmax=-999. do i=1,100 call ranmar(rvec,1) r=0+rvec(1)*(10.-0.) f=0.5*exp(-r/2.) if(fmax.lt.f)fmax=f end do fmax=1.2*fmax ntot=0 do i=1,10000 call ranmar(rvec,1) r=0+rvec(1)*(10.-0.) z=0.5*exp(-r/2.),if(z.gt.fmax)then fmax=z*1.2 write(
10、6,*)z greater than fmax end if call hfill(10,r,0.,1.0) call ranmar(rvec,1) u=rvec(1)*fmax if(u.lt.z)then call hfill(20,r,0.,1.0) call hfill(30,r,u,1.0) ntot=ntot+1 end if end do write(6,*)ntot=,ntot return end,2019-11-08,19,舍选法举例,void acc_rej() TH1F *h11 = new TH1F(“h11“,“,100,0,10); TH1F *h12 = new
11、 TH1F(“h12“,“,100,0,10); TH2F *h2 = new TH2F(“h2“,“,100,0,10,100,0,1); double fmax=-999.; for (int i=0;iUniform(0,10); double f = 0.5*exp(-r/2.); if (fmax f) fmax=f; fmax *= 1.2; cout “fmax=“ fmax endl;,int ntot=0; for (int i=0;iUniform(0,10); double z=0.5*exp(-r/2.); if (zfmax) fmax=1.2*z; h11-Fill
12、(r); double u=gRandom-Uniform(0,fmax); if (uFill(r); h2-Fill(r,u); ntot += 1; cout “ntot=“ ntot endl; ,ROOT脚本,2019-11-08,20,舍选法举例(续),频数,频数,舍选法存在效率问题。,二维均匀分布,2019-11-08,21,函数变换法与舍选法,函数变换法,优点：100%的抽样效率,缺点：函数须解析可积,舍选法,优点：方法简单，可用于非常复杂的函数,缺点：需要估计函数最大值，而且抽样效率低,粒子物理与核物理中，对常用的概率密度函数有各种建议采用的方法(见http:/pdg.lbl
13、.gov/2008/reviews/monterpp.pdf)。除此之外，舍选法最为常用。,初学者常犯的错误：对同一个过程做计算机模拟批处理， 没有考虑在批处理结果中存在随机数的重复性。,2019-11-08,22,常用概率密度分布函数的抽样,高斯(正态)分布, gROOT-Reset(); hx = new TH1F(“hx“,“x dis.“, 100,-10,10); gRandom-SetSeed(); Double_t x; const Double_t sigma=2.0; const Double_t mean=1.0; const Int_t kUPDATE = 1000; f
14、or ( Int_t i=0; iGaus(mean,sigma); hx-Fill(x); ,产生平均值为mean 标准偏差为sigma的 高斯分布。,在ROOT环境下采用已有的分布，可以容易完成布置的练习。,2019-11-08,23,蒙特卡罗统计检验,例如常用来检验理论与实验符合好坏的2 分布。,四个服从 N(0,1) 正态 分布的且相互独立的 随机变量平方和,一定符合自由度 为 4 的 2 分布,思考：如果出现不符合的情况，该如何解释？,2019-11-08,24,Toy 蒙特卡罗方法,粒子物理与核物理在实验的早期设计阶段，通常利用Toy 蒙特卡罗来估计可达到的测量精度（也称黑盒子方法
15、）。,在不做探测器模拟的情况下，可以对稳定的末态粒子动量 各分量进行含高斯分辨率的抽样，能损大小进行朗道分布 抽样，寿命进行指数分布抽样，等等，然后在所有末态中 寻找中间不稳定态E，根据能动量关系计算其对应的质量， 得到的质量分布称为Toy 蒙特卡罗结果。,2019-11-08,25,蒙特卡罗物理产生子,目的： 将理论用于某种物理过程的事例产生,输出量： 为对应某一物理过程的事例。对于每个事例，给出过程 产生的末态粒子和对应的动量,在粒子物理与核物理实验数据分析中，为了验证某一理论或模型，常常需要理论家提供蒙特卡罗物理产生子。,2019-11-08,26,蒙特卡罗物理产生子(续),简单情形,产
16、生 与 ,粒子物理与核物理中常用的产生子程序包,JETSET(PYTHIA) HERWIG ARIADNE,ISAJET PYTHIA HERWIG,KORALW EXCALIBUR ERATO,2019-11-08,27,蒙特卡罗探测器模拟,从产生子中输入粒子种类与动量，然后模拟粒子的输运过程,模拟探测器响应,多重散射(产生散射角) 粒子衰变(产生寿命) 电离能损(产生能损) 电磁与强子簇射 产生信号，电子学响应 ,输出量 = 模拟的数据,输入重建分析软件,用途: 预测“物理产生子层面上的”给定假设在“探测器层面上”应该观测到的响应。,通用软件包：GEANT3(FORTRAN)，GEANT4(C+),粒子与核物理中模拟的应用,用于实验初期的设计阶段建模分析 用于了解实验可能遇到物理过程的基本特征 用于了解实验仪器自身所受到的各种影
