2、。對于一個粒子經(jīng)歷以下過程：1．源抽樣：由于是點源，能量和位置的分布均為δ函數(shù)，抽樣得到。抽樣得：。2．到達閃爍體：如果1/2，無法到達閃爍體，丟棄，返回源抽樣重新產(chǎn)生粒子。粒子進入閃爍體的瞬間，狀態(tài)為：其他量不變。3．輸運過程：抽樣得到到下次碰撞的距離，根據(jù)當前粒子狀態(tài)中的算出下次碰撞的坐標，如果不在閃爍體區(qū)域（0碰撞有兩種可能：光電效應(yīng)和康普頓散射。根據(jù)粒子當前的能量，（由NaI(Tl)閃爍體宏觀界面數(shù)據(jù)）線性插值確定它的光電效應(yīng)截面和康普頓散射界面。抽樣得到本次反應(yīng)的類型。如果光電效應(yīng)，E=0，輸運過程結(jié)束。如果康普頓效應(yīng)，
3、抽樣獲取碰撞后的能量和運動方向（康普頓散射的能量分布密度函數(shù)知道，具體抽樣方法參考講義。）如果E<1KeV輸運過程結(jié)束，反之，重復(fù)本過程直到輸運過程結(jié)束。4．記錄與統(tǒng)計：記錄末態(tài)能量，計算沉積能量，考慮到測量系統(tǒng)分辨率，多道記錄能量為沉積能量的高斯展寬。記錄能量。其中，。由標準正態(tài)分布抽樣得到。五：程序如下：count=input('inputthecount:');%輸入模擬粒子數(shù)sigmaedata=[28370,13845,6908,2555,1223,2602,1925,905.3,479.7,164.5,74.24,23.86,66.60,36.62,22.29,9.978
4、,5.298,1.668,0.7378,0.2361,0.1099,0.06211,0.03939,0.02030];sigmacdata=[0.0220,0.0393,0.0568,0.0904,0.1209,0.1479,0.1722,0.2136,0.2480,0.3092,0.3486,0.3932,0.4153,0.4268,0.4319,0.4291,0.4215,0.3969,0.3691,0.3269,0.2944,0.2709,0.2512,0.2209];Edata=[1,1.5,2,3,4,5,6,8,10,15,20,30,40,50,60,80,100,15
5、0,200,300,400,500,600,800];%截面數(shù)據(jù)channel=zeros(1,ceil(662/5)+10);%多道數(shù)組nget=0;%探測到的總計數(shù)ntotal=0;%進入探測器的總計數(shù)forii=1:count%count個粒子循環(huán)collidetime=0;%當前粒子碰撞次數(shù)E0=input('inputtheenerg0:');%輸入模擬面源的一個能量E1=input('inputtheenergy1:');%輸入模擬面源的另一個能量E2=input('inputtheenergy2:')%輸入模擬面源的另另一個能量%粒子狀態(tài)初始化E0=622;E=E0;z
6、=-2;r=0;theta=2*pi*rand(1);%源抽樣，z,r,theta坐標miu=2*rand(1)-1;fai=2*pi*rand(1);%方向角抽樣ifmiu1%一個粒子在閃爍體中的輸運過程sigmae=interp1(Edata,sigmaedata,E,'linear');sigmac=interp1(Edata,sigmacdata,E,'linear');sigmat=sigmae+sigmac;%線性插值
7、得到截面數(shù)據(jù)L=-log(rand(1))/sigmat;%下次碰撞的距離%計算下次碰撞位置坐標rnew=sqrt(r^2+L^2*(1-miu^2)+2*r*L*sqrt(1-miu^2)*cos(fai-theta));z=z+L*miu;cdth=(rnew^2+r^2-L^2*(1-miu^2))/2/r/rnew;sdth=L*sqrt(1-miu^2)*sin(fai-theta)/rnew;dtheta=asin(sdth);ifcdth<0
