Geant4 note1 包含材料、几何体、粒子源设置
Geant4基础
路径
path/to/geant4.10.06/examples/basic/B1
程序架构
- main文件:exampleB1.cc,任何工程文件夹下包含唯一的main文件,在main文件中包含唯一的main函数。其既是C++工程的切入点妈也是结束点。在工程中用于串联其他程序文件。
一些C++的知识
- 头文件.h:变量、常量、函数、类的声明
- 源文件.c:变量的定义和函数的实现
具体解析
1.ActionInitialization.hh/cc:与Geant4的多线程运行有关
2.DetectorConstruction.hh/cc:定义类"DetectorConstruction",指定本次Run中的所有材料。
3.PrimaryGeneratorAction.hh/cc:定义类"PrimaryGeneratorAction",指定本次Run中的粒子源项。
4.RunAction.hh/cc:定义类"RunAction",指定在本次Run中,Run的初始化/开始/结束的所有操作。5.EventAction.hh/cc:定义类"EventAction",指定在本次Run中,每一个Event的初始化/开始/结束的所有操作。
6.SteppingAction.hh/cc:定义类"SteppingAction",指定在本次Run中,每一个Step的初始化/开始/结束的所有操作。
main文件解析
#include "B1DetectorConstruction.hh"
#include "B1ActionInitialization.hh"
#ifdef G4MULTITHREADED
#include "G4MTRunManager.hh"
#else
#include "G4RunManager.hh"
#endif
#include "G4UImanager.hh"
#include "QBBC.hh"
#include "G4VisExecutive.hh"
#include "G4UIExecutive.hh"
#include "Randomize.hh"///以上为include需要的头文件
int evtid=0;//定义外部变量,用于记录当前Event的序号
int main(int argc,char** argv){
G4UIExecutive* ui = 0;
if ( argc == 1 ) {
ui = new G4UIExecutive(argc, argv);
}
//判断为图形交互模式还是批处理模式
//图形交互模式:图形界面(在安装Geant4时安装了图形引擎),可在图形界面输入运行命令
//批处理模式:没有图形交互界面,输入指令后直接运行
//如果argc==1则为图形交互模式
G4Random::setTheEngine(new CLHEP::RanecuEngine);
//初始化Geant4产生随机数的引擎
#ifdef G4MULTITHREADED
G4MTRunManager* runManager = new G4MTRunManager;#else
G4RunManager* runManager = new G4RunManager;#endif
//初始化RunManager类的对象,此对象起到管理整个Run的作用,G4MTRunManager为多线程下的RunManager
runManager->SetUserInitialization(new B1DetectorConstruction());
//通过RunManager指定本次Run中的所有材料,即DetectorConstruction类的对象
G4VModularPhysicsList* physicsList = new QBBC;
physicsList->SetVerboseLevel(1);
runManager->SetUserInitialization(physicsList);
//通过RunManager指定本次Run中的相互作用列表
runManager->SetUserInitialization(new B1ActionInitialization());
//指定本次Run中粒子源(PrimaryGeneratorAction);
//指定在本次Run中,Run的初始化/开始/结束的所有操作。(RunAction);
//指定在本次Run中,每一个Event的初始化/开始/结束的所有操作。(EventAction);
//指定在本次Run中,每一个Step的初始化/开始/结束的所有操作。(SteppingAction);
//一般在多线程中使用,非多线程程序可直接在main函数中指定上述对象,方法与材料的指定方法相同
G4VisManager* visManager = new G4VisExecutive;
visManager->Initialize();
//初始化图形交互界面(Visualization UI)
G4UImanager* UImanager = G4UImanager::GetUIpointer();
//初始化负责执行用户输入命令的manager
if ( ! ui ) {
//批处理模式
G4String command = "/control/execute ";
G4String fileName = argv[1];
UImanager->ApplyCommand(command+fileName);
}
else {
//图形交互模式
UImanager->ApplyCommand("/control/execute init_vis.mac");
ui->SessionStart();
delete ui;
}
delete visManager;
delete runManager;
//释放内存,删除前面动态分配定义的runManager和visManager
}
- G4RunManager 是G4中唯一一个运行管理类,在单线程模式中,必须在main()函数中显示创建。在多线程模式中,G4MTRunManager 必须在main()函数中显示创建。
1.设置判断采用何种模式运行的条件
2.初始化RunManager
3.指定本次Run中的所有材料(DetectorConstruction)
4.指定本次Run中的相互作用列表(PhysicsList)
5.指定本次Run中的粒子源
6.指定本次Run中的Run/Event/Step的各项操作
7.执行用户输入命令8.释放内存并结束
Geant4材料
基础
用于确定材料的种类、物理状态、几何位置
材料的种类、物理状态
原子层面定义——G4Element,描述原子序数、核子数
材料定义——G4Material,描述温度、密度
(1)直接定义简单材料,用到原子序数z、摩尔质量a、材料密度density
G4double z , a , density ;
G4string name ;
density = 1.390 *g/cm3;
a = 39.95*g/mole ;
name = "test";
z = 18;
G4Material*Al = new G4Material ( name , z , a ,density );
(2) 通过组成分子来定义材料
- 分子的定义需要定义分子的名字、符号、原子序数、摩尔质量。
- 再用分子组成材料时,需要定义密度、成分数量、原子数
例:如何构建H2O
G4double z , a , density;
G4string name , symbol;
G4 int nchengfen , nyuanzs;
a = 1.01 *g/mole;
G4Element*elH = new G4Element(name = "Qing",symbol = "H" , z = 1. ,a);
a = 16 *g/mole;
G4Element*elO = new G4Element(name ="yang",symbol = "O" , z =8. , a );
density = 1 *g/cm3;
G4Material*H2O = new G4Material(name = "water" , density , nchengfen = 2);
H2O->AddElement(elH , nyuanzs = 2);
H2O->AddElement(elO , nyuanzs = 1);
(3)通过质量分数来定义混合物材料
- 分子的定义需要用到分子的名字、符号、原子序数、摩尔质量。
- 材料的定义需要用到成份数、百分含量、密度。
例:空气的定义,N2 70% ,O2 30%
G4double z , a , fractionmass , density;
G4double name , symbol ;
G4int nchengfen;
a = 14.01 *g/mole;
G4Element*elN = new G4Element(name = "dan" , symbol = "N" , z = 7. , a );
a = 16 *g/mole;
G4Element*elO = new G4Element(name = "yang",symbol = "O",z = 8. ,a);
density = 1.290 *g/cm3
G4Material*Air = new G4Material(name = "AIr" , density , nchengfen = 2);
Air->AddElement(elN , fractionmass = 70*perCent);
Air->AddElement(elO , fractionmass = 30*perCent);
(4)通过Geant4材料数据库中直接定义材料
- 首先要创建对象,以访问Geant4材料数据库
G4NistManager* man = G4NistManager::Instance(); //创建对象以访问Geant4数据库
G4Material*H2O = man->FindOrBuildMaterial("G4_WATER");
G4Material*AIR = man->FindOrBuildMaterial("G4_AIR");
(5)在已有材料的基础上创建新材料
- 需要定义密度
G4doubule density;
density = 1.05*mg/cm3;
G4Material*water1 = new G4Material("Water_1.05" , density , "G4_WATER");
//在数据库已有材料“G4_WATER”中创建密度不同的新材料
G4double density;
density = 1.03*mg/cm3;
G4NistManger*man = G4NistManager::Instance();
G4Material*water2 = man->BuildMaterialWithNewDensity("Water_1.03" ,"G4_WATER",density);
材料的几何位置
Geant4中,几何体分为以下三种:Solid(实体)、Logical Volume(逻辑几何体)、Physical Volume(物理几何体)
(1)Solid:直接创建的几何体,仅描述几何体的形状和尺寸。
(2)Logical Volume:在Solid基础上创建,指定其材料的种类等。
(3)Physical Volume:在逻辑几何体基础上创建,指定其位置,包含了真实世界几何体的所有属性。
创建实体
例:创建长方体Box,其中各个参数为真实长宽高的一半。在创建不同实体时调用不同的构造函数,具体参照https://link.zhihu.com/?target=http://geant4-userdoc.web.cern.ch/geant4-userdoc/UsersGuides/ForApplicationDeveloper/fo/BookForApplicationDevelopers.pdf
G4double world_hx = 3.0*m;
G4double world_hy = 1.0*m;
G4double world_hz = 1.0*m;
G4Box*worldBox = new G4Box("WorldBox" , world_hx , world_hy , world_hz);
G4Box*environmentBox
= new G4Box("environmentBox" , world_hx/2 , world_hy/2 , world_hz/2);
创建逻辑几何体
- 创建逻辑几何体不仅要有solid,还要定义材料。
G4LogicalVolume*worldLog = new G4LogicalVolume(worldBox , Air , "WorldLog");
G4LogicalVolume*envLog = new G4LogicalVolume(environmentBox , H2O , "environmentLog");
- G4LogicalVolume构造函数包含指定实体、填充材料、名称、磁场情况、是否为探测器的一部分等等。其实例如下:
G4LogicalVolume( G4VSolid* pSolid,//指定依赖的实体
G4Material* pMaterial,//指定材料的种类
const G4String& Name,//名称
G4FieldManager* pFieldMgr=0,//指定所在磁场的状况,默认为没有磁场
G4VSensitiveDetector* pSDetector=0,//指定是否为探测器的一部分,默认不是
//如果指定其为探测器的一部分,则可以在SteppingAction中获取粒子在该部位相互作用的信息
G4UserLimits* pULimits=0,//指定是否应用用户自定义限制
//用户自定义限制是指用户自行设置的关于粒子在逻辑几何体中的最大步长等逻辑几何体属性的设置
G4bool Optimise=true)//指定是否开启优化,默认开启
//此处的优化指的是Geant4在处理粒子在逻辑几何体中的track时将几何体体素化的方法
创建物理几何体
利用已创建的逻辑几何体创建物理几何体,使用G4PVPlacement类创建需要定义一个bool量检测有无重叠部分
bool checkOverlaps = true ;
G4PVPlacement*physWorld = new G4PVPlacement(0 ,
G4ThreeVector(),
worldLog,
"physWorld",
0,//world几何体不依附于任何母体
false,
0,//给这个物理体分配一个编号
checkOverlaps);
G4PVPlacement*physenv = new G4PVPlacement(0,
G4ThreeVector(),
envLog,
"phyenvironment",
worldLog, //依附于逻辑几何体“worldLog”
flase,
0,
checkOverlaps);
- G4PVPlacement类依次的参数为:旋转参数、位置、逻辑几何体、名称、依附母逻辑体、布尔运算、物理题编号、检查重叠。
- 创建LogVolume时可以指定该几何体是否为探测器的一部分,但实际指定哪一部分为探测器一般在StepppingAction中完成。
上述所有步骤创建几何体可视化如下所示
DetectorConstruction
与main文件的连接
将main文件与DetectorConstruction连接共需要两步
(1)引入相关头文件
#include "DetectorConstruction.hh"
(2)在main函数中利用runManager引入新DetectorConstruction()类,将其指定为程序中的材料
runManager->SetUserInitialization( DetectorConstruction() );
DetectorConstruction.hh
#ifndef B1DetectorConstruction_h
#define B1DetectorConstruction_h 1
//#ifndef与#endif防止头文件的重复包含和编译
#include "G4VUserDetectorConstruction.hh" //一般用户自定义的DetectorConstruction类继承于此
#include "globals.hh"//包含Geant4中的全局变量等,除Actioninitializtion外均应包含
class G4VPhysicalVolume;//DetectorConstruction类中成员函数返回值中包括此类,故应提前声明
class DetectorConstruction : public G4VUserDetectorConstruction
//一般用户自定义的DetectorConstruction继承于G4VUserDetectorConstruction
{
public:
DetectorConstruction();//构造函数声明,将在对应源文件中定义
virtual ~DetectorConstruction();//析构函数声明,将在对应源文件中定义
virtual G4VPhysicalVolume* Construct();//成员函数Construct()声明,将在对应源文件中定义
//一般在此函数中指定材料的相关信息
};
#endif//#ifndef与#endif防止头文件的重复包含和编译
- 在头文件中声明了G4PhysicalVolume类以及DetectorConstruction类
- 需要引入G4UserDetectorConstruction.hh(用户定义的DetectorConstruction继承于此)、globals.hh(包含Geant4中的全局变量)
- 可见DetectorConstruction中有三个成员函数,该三个成员函数将在.cc文件中一个个定义
DetectorConstruction.cc
#include "DetectorConstruction.hh"//包含该源文件对应的头文件
#include "G4NistManager.hh"//包含在材料种类定义时所使用的数据库
#include "G4Box.hh"//包含在实体定义时所使用的形状:长方体
#include "G4LogicalVolume.hh"
#include "G4PVPlacement.hh"
//包含逻辑几何体,物理几何体(G4PVPlacement类)
#include "G4SystemOfUnits.hh"//包含物理单位制
DetectorConstruction::DetectorConstruction()
: G4VUserDetectorConstruction()
{ }
//定义构造函数,一般为空
DetectorConstruction::~DetectorConstruction()
{ }
//定义析构函数,一般为空
G4VPhysicalVolume* DetectorConstruction::Construct()//定义成员函数Construct()
{
G4double z, a, density;
density = 1.390*g/cm3;
//指定材料的密度
G4String name;
a = 39.95*g/mole;
//指定材料的摩尔质量
G4Material* lAr = new G4Material(name="liquidArgon", z=18., a, density);
//创建名称为"liquidArgon"的材料,其材料的原子序数为18,摩尔质量为39.95*g/mole,密度为1.390g/cm3。
bool checkOverlaps=true;// 指定是否检查有无重叠部分,一般应设置为true
G4double world_hx = 3.0*m;
G4double world_hy = 1.0*m;
G4double world_hz = 1.0*m;
//指定几何体的长宽高
G4Box* worldBox
= new G4Box("WorldBox", world_hx, world_hy, world_hz);
//创建一个长方体(盒子),名称为"WorldBox",变量指针为"worldBox"
//在其自身的坐标系中,其X轴的范围是-3m到3m,总长6m;其Y轴的范围是-1m至1m,总长2m;
//其Z轴的范围是-1m至1m,总长2m
G4Box* envBox
= new G4Box("environmentBox",world_hx/2,world_hy/2,world_hz/2);
//创建一个长方体(盒子),名称为"environmentBox",变量指针为"env"
//在其自身的坐标系中,其X轴的范围是-1.5m到1.5m,总长3m;其Y轴的范围是-0.5m至0.5m,总长1m
G4LogicalVolume* worldLog
= new G4LogicalVolume(worldBox, lAr, "WorldLog");
//创建一个逻辑几何体,名称为"WorldLog",其实体为"worldBox",材料为5.1节中的"Air"
G4LogicalVolume* envLog
= new G4LogicalVolume(envBox, H2O, "environmentLog");
//创建一个逻辑几何体,名称为"environmentLog",其实体为"envBox",材料为5.1节中的"H2O"
G4VPhysicalVolume* physWorld
= new G4PVPlacement(0,G4ThreeVector(),worldLog,"physWorld",0,false,0,checkOverlaps);
//创建一个物理几何体,名称为"physWorld",其逻辑几何体为"worldLog";
//它不依附于任何母体,本身为Geant4中的World,这种不依附任何母体的World应该有且仅有一个。
G4VPhysicalVolume* physenv
= new G4PVPlacement(0,G4ThreeVector(),envLog,"physenvironment",worldLog,false,0,checkOverlaps);
//创建一个物理几何体,名称为"physenvironment",其逻辑几何体为"envLog",依附于逻辑几何体"worldLog",
//即依附于World
return physWorld;//返回World对应的物理几何体,实际模拟时即使用该物理几何体
}
- 可以看出,在DetectorConstruction中仅需对成员函数Construct()进行详细定义,其定义过程如下
(1)定义材料G4Mterial* name = new G4Material();从而定义一个新材料name
(2)定义几何体Solid->LogVolume->phyVolume - 需要引入的头文件
(1)DetectorConstruction.hh
(2)G4NistManager.hh(材料定义时所使用的数据库)
(3)定义相关形状引入相关的形状头文件
(4)G4LogicalVolume.hh以及G4PYPlacement.hh
(5)G4SystemOfUnits.hh 包含物理单位制 - construct()函数末尾return physWorld; 返回World对应的物理几何体,实际模拟时即使用该物理几何体
Geant4粒子源
简介
- 粒子源具有两个变量 (1)粒子源的物理信息 ;(2)几何信息
- 在“PrimaryGeneratorAction.hh/cc”文件中定义以上信息
粒子源的物理信息
- 即确定粒子的种类和能量,特别注意,对于gamma/e-/e+,为避免红外发散,需要设置第三个信息截断范围。在截断范围以下,三种粒子不再产生任何次级粒子
种类
例:以创建alpha粒子为例说明具体的程序
G4int n_particle = 1;
fParticleGun = new G4ParticleGun(n_particle);//首先设置粒子产生数
G4ParticleTable*particletable = G4ParticleTable::GetParticleTable();
G4String particleName;
G4ParticleDefinition*alp = ParticleTable->FindParticle(particleName = "alpha");
fParticleGun->SetParticleDefinition(alp);
- 需要一个int类型变量,设置每次event发射的粒子数
G4ParticleGun::G4ParticleGun(G4int numberofparticles)
{
SetInitialValues();
NumberOfParticlesToBeGenerated = numberofparticles;
}
- 类似材料的定义,对于粒子的种类的定义需要通过G4ParticleTable对象查找粒子。
能量
依旧通过G4ParticleGun成员函数进行定义。
fParticleGun->SetParticleEnergy(1000*keV);
- 过程总结
(1)定义粒子发射数量->定义粒子类型->定义粒子能量
(2)设置粒子相关物理信息参量通过G4ParticleGun的成员函数实现
(3)寻找粒子种类通过G4ParticleTable实现
粒子源的几何信息
即确定粒子的位置和发射动量方向
位置
真实世界中粒子源具有一定的尺寸,因此粒子源位置在一定的几何范围内随机分布。利用“* G4UniformRand()”来实现该功能。该随机数函数产生0-1均匀分布的随机数
例:产生一个再1m* 1m * 2m范围内随机分布的粒子源
G4double envSizeXY = 1*m;
G4double envSizeZ = 2*m;
G4double x0 = envSizeXY*(G4UniformRand()-0.5);
G4double y0 = envSizeXY*(G4UniformRand()-0.5);
G4double z0 = envSizeZ*(G4UniformRand()-0.5);
fParticleGun -> SetParticlePosition( G4ThreeVector(x0 , y0 , z0) );
发射动量方向
通过确定粒子动量再xyz三个方向的分配情况,就可以确定粒子的发射方向。
G4double x = 0;
G4double y = 0;
G4double z = 1;
fParticleGun -> SetParticleMomentumDirection(G4ThreeVector(x , y ,z));
//按照此方法确定的动量分配方式为:
//x方向:分配|a|/sqrt(a^2+b^2+c^2)的动量
//y方向:分配|b|/sqrt(a^2+b^2+c^2)的动量
//z方向:分配|c|/sqrt(a^2+b^2+c^2)的动量
//如果a/b/c为负数,则表示与对应坐标轴正方向相反
//也可以采用随机分布的方法确定动量方向,与确定位置类似
PrimaryGeneratorAction
与ActionIntialization文件的连接
不同于DetectorConstruction,PrimaryGeneratorAction需要再ActionInitialization.cc中进行指定。
(1)引入头文件
#include "PrimaryGeneratorAction.hh";
(2)创建该类,通过SetUserAction方法指定其为Geant4程序中的粒子源
SetUserAction(new PrimaryGeneratorAction);
PrimaryGeneratorAction.hh
#ifndef PrimaryGeneratorAction_h
#define PrimaryGeneratorAction_h 1
//#ifndef与#endif防止头文件的重复包含和编译
#include "G4VUserPrimaryGeneratorAction.hh" //一般用户自定义的PrimaryGeneratorAction类继承于此
#include "G4ParticleGun.hh"//使用G4ParticleGun作为粒子源
#include "globals.hh" //包含Geant4中的全局变量等,除Actioninitializtion外均应包含
class G4ParticleGun;//使用G4ParticleGun作为粒子源
class G4Event;//成员函数中的参数中含有G4Event类,故提前声明
class PrimaryGeneratorAction : public G4VUserPrimaryGeneratorAction//一般用户自定义的PrimaryGeneratorAction类需继承G4VUserPrimaryGeneratorAction
{
public:
PrimaryGeneratorAction();//构造函数声明,将在对应源文件中定义
virtual ~PrimaryGeneratorAction();//析构函数声明,将在对应源文件中定义
virtual void GeneratePrimaries(G4Event*);//成员函数声明,将在对应源文件中定义
//一般在此函数中指定粒子源的相关信息
const G4ParticleGun* GetParticleGun() const { return fParticleGun; }
//用于访问私有变量 fParticleGun
private:
G4ParticleGun* fParticleGun; // G4ParticleGun的指针fParticleGun
};
#endif //#ifndef与#endif防止头文件的重复包含和编译
PrimaryGeneratorAction.cc
#include "PrimaryGeneratorAction.hh"
#include "G4ParticleTable.hh"
#include "G4ParticleDefinition.hh"
//包含粒子数据库
#include "G4SystemOfUnits.hh"
#include "Randomize.hh"
PrimaryGeneratorAction::PrimaryGeneratorAction()
: G4VUserPrimaryGeneratorAction() , fParticleGun(0)
//在该函数中构造粒子类型、能量以及发射动量方向
{
G4int n_particle = 1;
fParticleGun = new G4ParticleGun(n_particle);
G4ParticleTable* particletable = G4ParticleTable::GetparticleTable();
G4String particlename;
G4ParticleDefinition* particle = particletable->FindParticle(particlename = "alpha");
fParticleGun -> SetParticleDefinition(particle);
//定义粒子类型
G4double a = 0;
G4double b = 0;
G4double c = 1;
fParticleGun->SetParticleMomentumDirtection(G4ThreeVector( a , b , c ));
//定义发射方向
fParticelGun->SetParticleEnergy(1000*keV)
}
PrimaryGeneratorAction::~PrimaryGeneratorAction()
{
delete fParticleGun;
}
void PrimaryGeneratorAction::GeneratorPrimaries(G4Event* anEvent)
{
G4double envSizeXY = 1*m;
G4double envSizeZ = 2*m;
G4double x0 = envSizeXY*(G4UniformRand()-0.5);
G4double y0 = envSizeXY*(G4UniformRand()-0.5);
G4double z0 = envSizeZ*(G4UniformRand()-0.5);
fParticleGun->SetParticlePosition(G4ThreeVector( x0 , y0 , z0 ));
//定义发射位置
fParticleGun->GeneratorPrimaryVertex(anEvent);
//指定本次Run中每个event开始时的粒子的产生信息
}