Geant4学习一:写一个简单程序
1、 基本介绍
Geant4是由欧洲核子中心(CERN)主导开发的一款采用C++语言编写的蒙特卡洛通用粒子输运软件。Geant4能够模拟各种粒子与不同物质的相互作用以及输运过程。其开源性和C++优越的编程性能使得Geant4在核物理与辐射探测、放射性医学等领域都得到了广泛的应用,并且数据库非常全面,更新速度很快。Geant4拥有丰富、全面的物理模型以及适用于不同领域的物理模型列表可供用户选择,用户还可以根据需要自行定义。Geant4非常适合用来开发特定领域的专用模拟程序,如TOPAS、Gate都是基于Geant4开发的放射医学领域的专用程序。
在Geant4中有三个逻辑层次,分别是Step、Event和Run,Run包含整个模拟的总信息,Event包含每个源粒子的全部输运过程,Step则表示某一个粒子的某一次相互作用,可以从中获得该粒子的能量(Energy)、动量(Momentum)、位置(Position)、输运时间(GlobleTime)、物理过程的名字(ProcessName)等信息。Step中记录的信息可以累积到Event当中处理,Event中累积的信息则可以导入到Run当中。
Geant4可以通过代码编辑复杂的几何形状,可以使用布尔运算,支持从分子尺度到行星尺度的粒子输运计算。程序支持对各种同位素、混合物进行自定义,并提供了丰富的参考案例和模板,同时拥有强大的用户界面以及可视化引擎便于调试。
2、程序构成
当我们要构建一个简单的Geant4程序来进行计算时需要在程序中设定几何结构、材料、粒子源等信息,并对感兴趣的物理信息进行获取和处理。以14MeV中子入射聚乙烯,记录出射的中子的能谱为例。
2.1、构建材料
我们模拟的对象都由特定的材料构成,Geant4中提供了构建各种材料的多个类,同位素(G4Isotope)、元素(G4Element)、材料(G4Material)。其中元素可以由一个或者多个同位素构成,材料可以由一个或者多个元素构成,也可以是若干元素和若干材料合成。我们可以按照这一逻辑定义任意类型的材料。同时Geant4的源代码中已经提供了各种基本同位素、元素以及常用材料的定义可以通过指针调用。在中子入射聚乙烯这一模拟中我们需要构建聚乙烯这一材料,以及环境中的空气。聚乙烯由氢元素和氧元素按照2:1的数量比定义如下:
G4Element* elH = new G4Element("Hydrogen" ,"elH" , 1., 1.01 *g/mole);//氢元素
G4Element* elC = new G4Element("carbon" ,"elC" , 6., 12.01 *g/mole);//氧元素G4Material* mPE = new G4Material("PE",0.95* g/cm3,2);//核数比定义聚乙烯 mPE->AddElement(elH, 2); mPE->AddElement(elC, 1);
空气我们可以直接调用Geant4提供的材料库
G4NistManager* man = G4NistManager::Instance();G4Material *mAir = man->FindOrBuildMaterial("G4_AIR");//空气
G4NistManager中的常用材料可以在源代码中查看(geant4.10.05\source\materials\src\G4NistMaterialBuilder.cc),其中对空气的定义如下:
void G4NistMaterialBuilder::NistCompoundMaterials()
{AddMaterial("G4_AIR", 0.00120479, 0, 85.7, 4, kStateGas);AddElementByWeightFraction( 6, 0.000124);AddElementByWeightFraction( 7, 0.755267);AddElementByWeightFraction( 8, 0.231781);AddElementByWeightFraction(18, 0.012827);
]
c此时我们完成了两个材料的定义及空气(mAir)和聚乙烯(mPE)。随后我们将对定义好的几何结构添加这两种材料。
2.2、构建几何
Geant4的几何设定有一个母体与子体的概念(mother&daughter),每一个程序都要由一个基本的计算空间,也就是最大的母体,命名为World,一般我们把它作为环境空间,材料通常是空气或者真空。子体是设置在母体内的,World没有母体,其他 体都必须设置在一个唯一的母体内。每个体的内部都可以设置0或者若干个子体。在这个简单程序里我们需要设置World和聚乙烯板这两个体,World作为母体,聚乙烯板是World的一个子体。Geant4中有多种几何体的类,常用的有立方体(G4Box)、球体(G4Sphere)、圆柱体(G4Tubs)等。这里我们将World和聚乙烯板都定义为立方体,聚乙烯板的尺寸为长宽20cm,后附3cm,世界的尺寸为长宽高均为24cm,下面为具体定义:
// WorldG4double PE_sizeXY = 20*cm, PE_sizeZ = 3*cm;//定义尺寸参数G4double world_sizeXY = 1.2*PE_sizeXY;G4double world_sizeZ = 8*PE_sizeZ;G4Box* solidWorld = new G4Box("World", //its name0.5*world_sizeXY, 0.5*world_sizeXY, 0.5*world_sizeZ);//its sizeG4LogicalVolume* logicWorld = new G4LogicalVolume(solidWorld, //its solidmAir, //its material"World"); //its nameG4VPhysicalVolume* physWorld = new G4PVPlacement(0, //no rotationG4ThreeVector(), //Position at (0,0,0)logicWorld, //its logical volume"World", //its name0, //its mother volumefalse, //no boolean operation0, //copy numbercheckOverlaps); //overlaps checking//聚乙烯板G4Box* solidPE = new G4Box("PEBox", //its name0.5*PE_sizeXY, 0.5*PE_sizeXY, 0.5*PE_sizeZ); //its sizeG4LogicalVolume* logicPE = new G4LogicalVolume(solidPE, //its solidmPE, //its material"PEBox"); //its nameG4ThreeVector pos1 = G4ThreeVector(0*cm, 0*cm, -1*cm);//定义位置坐标new G4PVPlacement(0, //no rotationpos1, //PositionlogicPE, //its logical volume"PEBox", //its namelogicWorld, //its mother volumefalse, //no boolean operation0, //copy numbercheckOverlaps); //overlaps checking
2.3、构建粒子源
Geant4构建粒子源有两种方法,一个使用粒子枪(fParticleGun)在PrimaryGeneratorAction中直接定义,一个是使用宏命令(G4GeneralParticleSource,GPS),这里我们先使用fParticleGun,设置粒子源的形状、尺寸、位置、类型、能量以及发射方向,以下是一个各向同性的中子点源,能量为2.2MeV:
PrimaryGeneratorAction::PrimaryGeneratorAction()
: G4VUserPrimaryGeneratorAction(),fParticleGun(0)
{G4int n_particle = 1;fParticleGun = new G4ParticleGun(n_particle);
//&&&&&&&&&&&&&&&&&&&&&&&&&&&&&Gamma&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&G4ParticleDefinition* particle = G4ParticleTable::GetParticleTable()->FindParticle("neutron");fParticleGun->SetParticleDefinition(particle);G4double Energy=2.2 ;//gamma能量G4double Px=0,Py=0,Pz=-100;//位置参数fParticleGun->SetParticlePosition( G4ThreeVector(Px,Py,Pz) );//设定位置fParticleGun->SetParticleEnergy(Eg);//设定能量
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
PrimaryGeneratorAction::~PrimaryGeneratorAction()
{delete fParticleGun;
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
void PrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent)
{//设置各向同性源//生成随机数float Seed1=G4UniformRand();//生成随机数(0,1)float Seed2=G4UniformRand();G4double cosTheta = 2* Seed1 - 1., phi = twopi* Seed2;G4double sinTheta = std::sqrt(1. - cosTheta*cosTheta);G4double ux = sinTheta*std::cos(phi),uy = sinTheta*std::sin(phi),uz = cosTheta;fParticleGun->SetParticleMomentumDirection(G4ThreeVector(ux,uy,uz));fParticleGun->GeneratePrimaryVertex(anEvent);//每个Event重新设定
}
2.4、获取信息
这里需要说明的是Geant4中每个粒子从产生到消失(截断值控制),发生的每一次相互作用,以及每次相互作用的位置、时间、动量与能量变化均可以提取。从每个初始粒子的产生到该初始粒子以及其产生的次级粒子全部截止,为一个Event,各种粒子都有专属的代码来识别,比如22代表γ,11代表电子。初始粒子以及其产生的次级粒子所发生的每一个相互作用为一个Step,可以通过ProcessName来确定相互作用的类型。每一个粒子都有一条径迹,由一个专属的TrackID,初始粒子的TrackID为1。每个粒子的径迹由若干个Step的点串联而成,每一条径迹的Step都从1开始编号,程伟StepNumber。此外可以通过ParentID来确定产生该粒子的母粒子的TrackID。比如初始中子的TrackID为1,中子的代码为2112,该中子由粒子源产生之后第一次与氢原子核发生弹性碰撞,产生了一个反冲质子,质子的代码为2212,该质子的TrackID为2,这个质子的ParentID为1,说明这是由初始中子产生的。这一个Step的相互作用的ProcessName为“hadElastic”,StepNumber为1,这里这个StepNumber是属于初始中子的。此时还可以获取中子和质子的能量和动量。需要单独说明的的是,在Geant4中与一个强制的Step,其ProcessName为“Transportation”,在一个粒子的径迹穿过某个体的边界进入另一个体是就会产生这一过程,这便于我们分析粒子从实体中逃逸的情况。这里我们要记录穿过聚乙烯板的中子的能量,就需要用到这一过程。中子入射到聚乙烯板中会与氢核和碳核发生弹性散射、非弹性散射以及辐射俘获。一部分被吸收,一部分被散射后从其他方向出射,一部分不发生相互作用直接穿过。既有透射也会由背散射,我们想记录从聚乙烯板背面出射的中子的内能,除了要判断中子是从聚乙烯运动到外面的World,还要通过穿越点的坐标限定穿出的位置在背面:
G4double Energy01;using std::ofstream;ofstream outfile1("./data/中子E.txt");if((aStep->GetTrack()->GetDefinition()->GetPDGEncoding()==2112)//判断为中子&& aStep->GetPostStepPoint()->GetPhysicalVolume() != NULL)//Step未超出计算边界{if(aStep->GetPreStepPoint()->GetPhysicalVolume()->GetName() == "PEBox"//前一点在聚乙烯中&& aStep->GetPostStepPoint()->GetPhysicalVolume()->GetName()== "World"//后一点在世界&& aStep->GetTrack()->GetPosition0().z()>4.9)//Step的位置,边界z坐标为5{Energy01=(aStep->GetTrack()->GetKineticEnergy());outfile1 << Energy01 << G4endl;} }
2.5、程序运行
在程序目录新建data文件夹用于存储数据。在init_vis.mac中设定运行参数,可以设定运行的粒子数然后直接运行,也可以调用可视化界面进行查看和操作。
# Macro file for the initialization of example TestEm4
# in interactive session
#
# Set some default verbose
/control/verbose 2
/control/saveHistory
/run/verbose 2
#
# Change the default number of threads (in multi-threaded mode)
/run/numberOfThreads 1
#
# Initialize kernel
#设置截断值
/run/setCut 0.1 mm
/run/initialize
#
# Visualization setting#运行可视化界面
#/control/execute vis.mac
#直接运行run
/run/beamOn 20000000
在终端输入make -j进行编译,然后运行程序即可。
Geant4学习一:写一个简单程序相关推荐
- python123程序设计题说句心里话_用c++写一个简单的计算器程序
// 050305.cpp : 定义控制台应用程序的入口点. // // 050304.cpp : 定义控制台应用程序的入口点. // //四则运算 #include "stdafx.h&q ...
- 用java做一个简单记事本_用记事本写一个简单的java程序
用记事本写一个简单的java程序 第一步: 安装好jdk,并设置好环境变量. 桌面-计算机(右键)-属性-高级系统设置-环境变量-path-在变量值后加上:和jdk安装路径加上(路径即为C:\Prog ...
- 怎样用java写一个简单的文件复制程序
怎样用java写一个简单的文件复制程序 代码来源:https://jingyan.baidu.com/article/c35dbcb0d6f1398916fcbc07.html package Num ...
- 用c++写一个简单的钓鱼(集卡)程序
用c++写一个简单的钓鱼(集卡)程序 因为我们C++的老师要求我们写一个面向对象的C++程序,要求是扑克牌游戏: 1.分花色,牌值表示能显示一副完整的牌 2.洗牌,分牌 3.出牌(游戏规则自己定,可以 ...
- 写一个简单的Java界面程序
写一个简单的Java界面程序 有时候未免想写一些有界面的java小程序练练手,那么如何写一个比较好看的界面话程序呢?下面小编就带你一步一步来搭建这个小洋房. 实现界面化编程要用到的一个主要包impor ...
- Java 百度AI 写一个简单的手势识别程序
教程地址:请关注我的https://edu.csdn.net/course/detail/23001 Java写一个简单的手势识别程序,这里采用百度是AI,视觉技术中的人体分析中的手势识别,识别图片中 ...
- DuiVision开发教程(2)-如何写一个简单的界面程序
基于DuiVision界面库开发的界面程序主要包括如下几部分内容: 1.资源定义,包括图片资源.各个窗口界面的xml定义文件 2.事件处理类代码,用于处理界面响应消息 3.其他业务逻辑代码 下面举例说 ...
- java递归怎么写_什么是递归?用Java写一个简单的递归程序
什么是递归?用Java写一个简单的递归程序 递归的定义 递归(recursion):以此类推是递归的基本思想,将规模大的问题转化为规模小的问题来解决. 递归的要素 自定义递归函数,并确定函数的基本功能 ...
- 如何搭建python框架_从零开始:写一个简单的Python框架
原标题:从零开始:写一个简单的Python框架 Python部落(python.freelycode.com)组织翻译,禁止转载,欢迎转发. 你为什么想搭建一个Web框架?我想有下面几个原因: 有一个 ...
- 用python写一个简单的web服务器
人生苦短,我用python 简洁高效,这才是理想的语言啊 分享一点python的学习经验-----如何用python写一个简单的web服务器 首先,我们需要简单地了解一下网络通信协议,这里用白话介绍一 ...
最新文章
- Android7.0 Rild工作流程
- 算法--06年华为面试:求两个数组的最小差值(Java实现)
- 国人主导研发的 HAWQ® 成 Apache® 顶级项目
- java超出gc开销_通过这5个简单的技巧减少GC开销
- java计算加速减速_Javascript加速运动与减速运动
- c语言5版第10章答案,第10章 指 针 参考答案 c语言(1)
- 作者:陈兴鹏(1963-),男,兰州大学资源环境学院教授、博士生导师。
- 特斯拉电动皮卡量产时间还会推迟 内部人士称已被推迟到2023年
- 2021年中国乙醛市场趋势报告、技术动态创新及2027年市场预测
- SecureCRT学习之道:SecureCRT常用快捷键设置与字体设置方法
- linux文本编辑命令vim查找,Linux编辑器vi中文本搜索与替换操作
- 疯狂的程序员 21-30
- godaddy域名转入步骤
- 考研数学常用基础知识默写版
- Jetty9开发(1)
- 记录你生活的点滴,体会分享的快乐
- 阿里云客服机器人人工服务配置文档
- 如何在WPF中使用虚拟键盘
- Lotus Notes Send EMail from VB or VBA
- 验证码识别dll库,识别率95%