粒子物理蒙特卡罗模拟库Geant4之能谱制作
生活随笔
收集整理的這篇文章主要介紹了
粒子物理蒙特卡罗模拟库Geant4之能谱制作
小編覺得挺不錯的,現在分享給大家,幫大家做個參考.
γ衰變模型
抱怨!!
兩日光陰浪費
在GEANT4那些物理學家寫的奇怪邏輯文檔中我難以找到“粒子γ衰變模型”。所以這個物理問題變成了一個文檔搜尋問題。
找不到文檔只能自己查找源代碼,所以只能使用grep工具,結果由于我忘了使用shell通道
root@master:# grep -r gamma ./ | grep decay
以至于浪費兩日時間。所以歸根結底是一個Linux shell問題。:)
現實意義
能譜的制作在粒子物理中有著異常重要的意義。那么如何模擬粒子衰變以及收集能量沉積并進行數據分析變得非常必要。通過對被研究對象能譜的獲取和分析可以直接或間接地獲得物質的結構、組成元素的種類與含量等重要信息。
構建粒子類型、記錄能譜及定義衰變鏈
void TrackingAction::PreUserTrackingAction(const G4Track* track) {Run* run = static_cast<Run*>(G4RunManager::GetRunManager()->GetNonConstCurrentRun());G4ParticleDefinition* particle = track->GetDefinition();G4String name = particle->GetParticleName();fCharge = particle->GetPDGCharge();fMass = particle->GetPDGMass(); G4double Ekin = track->GetKineticEnergy();G4int ID = track->GetTrackID();G4bool condition = false;// check LifeTime//G4double meanLife = particle->GetPDGLifeTime();//count particles//if (ID>1) run->ParticleCount(name, Ekin, meanLife);//energy spectrum//G4int ih = 0;if (particle == G4Electron::Electron()||particle == G4Positron::Positron()) ih = 1;else if (particle == G4NeutrinoE::NeutrinoE()||particle == G4AntiNeutrinoE::AntiNeutrinoE()) ih = 2;else if (particle == G4Gamma::Gamma()) ih = 3;else if (particle == G4Alpha::Alpha()) ih = 4;else if (fCharge > 2.) ih = 5;if (ih) G4AnalysisManager::Instance()->FillH1(ih, Ekin);if (ih) G4CsvAnalysisManager::Instance()->FillH1(ih, Ekin);//Ion//if (fCharge > 2.) {//build decay chainif (ID == 1) fEvent->AddDecayChain(name);else fEvent->AddDecayChain(" ---> " + name);// //full chain: put at rest; if not: kill secondary G4Track* tr = (G4Track*) track;if (fFullChain) tr->SetTrackStatus(fStopButAlive);else if (ID>1) tr->SetTrackStatus(fStopAndKill);//fTime_birth = track->GetGlobalTime();}//example of saving random number seed of this fEvent, under condition//condition = (ih == 3);if (condition) G4RunManager::GetRunManager()->rndmSaveThisEvent(); }嚴格來說、應該叫獲取粒子類型,因為即將展示如何在交互命令中或者腳本中定義衰變粒子類型。
衰變腳本
錒系元素Am241的衰變腳本
/control/verbose 2 /run/verbose 1 # /gun/particle ion /gun/ion 95 241 #在這里定義了類型 # /tracking/verbose 2 /run/beamOn 1 # /control/cout/ignoreThreadsExcept 0 /tracking/verbose 0 # /analysis/setFileName Am241 #存儲數據文件為.root格式 /analysis/h1/set 1 150 0. 1500 keV #e+ e- /analysis/h1/set 2 150 0. 1500 keV #neutrino /analysis/h1/set 3 150 0. 1500 keV #gamma /analysis/h1/set 6 100 0. 2500 keV #EkinTot (Q) /analysis/h1/set 7 150 0. 15e3 keV #P balance /analysis/h1/set 8 100 0. 100. year #time of life /analysis/h1/set 9 100 1. 3. MeV #EvisTot # /run/printProgress 100000 # /run/beamOn 1000數據格式轉換
.root格式是適用于歐洲原子能機構CERN的ROOT分析工具的數據格式,不適于其他數據工具使用。將其轉換為.xml文件:
root@master:# sudo vi root2xml.C {TFile *f = TFile::Open("Example.xml","recreate");TFile *mf = TFile::Open("Co60.root");TH1F *h = (TH1F*)mf->Get("3");f->cd();h->Write(); }可以得到.xml格式能譜數據:
<?xml version="1.0"?> <root setup="2xoo" ref="null" created="2017-11-16 22:13:05" modified="2017-11-16 22:13:05" uuid="43680090-cad8-11e7-9717-0100007fbeef" version="3" file_version="53430"><XmlKey name="3" cycle="1" created="2017-11-16 22:13:45"><Object class="TH1D" ref="id0"><TH1D version="1"><TH1 version="7"><TNamed version="1"><TObject fUniqueID="0" fBits="3000008"/><fName str="3"/><fTitle str="energy spectrum (%): gamma"/></TNamed><TAttLine version="2"><fLineColor v="1"/><fLineStyle v="1"/><fLineWidth v="1"/></TAttLine><TAttFill version="2"><fFillColor v="0"/><fFillStyle v="101"/></TAttFill><TAttMarker version="2"><fMarkerColor v="1"/><fMarkerStyle v="1"/><fMarkerSize v="1.000000e+00"/></TAttMarker><fNcells v="152"/><fXaxis><TAxis version="9"><TNamed version="1"><TObject fUniqueID="0" fBits="3000000"/><fName str="xaxis"/><fTitle str=" [keV]"/> ... ... ...繪圖
繪圖腳本
{gROOT->Reset();// Draw histos filled by Geant4 simulation // TFile f("Am241.root");TCanvas* c1 = new TCanvas("c1", " ");c1->SetLogy(1);c1->cd();c1->Update();TH1D* hist = (TH1D*)f.Get("3");hist->Draw("HIST"); }繪圖命令
root@master:# root -l Am241.root總結
以上是生活随笔為你收集整理的粒子物理蒙特卡罗模拟库Geant4之能谱制作的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 零基础写java网络爬虫
- 下一篇: 基与坐标系(阅读《理解矩阵》笔记)