Geant4中B1实例学习

注意:

  • 该笔记为学习B站UID298509409视频时的个人记录
  • 该笔记仅为个人学习、理解的备忘录,权威解释请看官网资料

编译运行Geant4中B1实例

  • 拷贝B1到目标文件夹(根据自己的geant4安装路径修改)
    • cp /opt/Geant4/geant41004/share/Geant4-10.4.0/examples/basic /home/username/
  • 进入复制的B1文件夹
    • cd /home/username/B1
  • 创建build文件夹
    • mkdir build
  • 进入build文件夹,并编译
    • cd build
    • cmake ..
  • make&&make install
    • make -j8
    • make install
  • 运行./exampleB1
    • ./exampleB1

建议学习一下cmake语法及其相关功能

运行结果如图: avatar

B1主函数以及运行脚本文件(.mac)

  • 主函数(exampleB1.cc)
    • 给不同的输入参数,执行不同功能
    • Detect interactive mode
    • batch mode
    • interactive mode
  • run1.mac

    • 初始化
    • 发射粒子、能量
  • init_vis.mac

    • 初始化交互界面
  • vis.mac

    • 设置了图形界面的信息:大小、颜色等等

Geant4是啥以及基本的概念

  • Geant4是函数库
    • 我们要做的就是继承以及调用其写好的代码,不要造轮子
    • 先研习Geant4自带的事例,查询是否有自己需要的
    • 问师兄,看论坛问询有无需要的功能
    • 看Geant4的说明书,查找需要的功能
  • run、Event、Step、Track
    • run:可以理解成一次测试或者一次实验
    • Event:run中的一次反应,包括入射粒子、所有产物,物质响应等等
    • Step:粒子与物质相互作用过程,如alpha被Si探测,其可能有多个step,类似粒子在Si中碰撞输运
    • Track:记录step中的信息,如位置、能量、时间、动量
  • 自己编写代码的思路
    • 构建几何模型:探测器材料、大小、位置等等
    • 设置粒子种类,反应类型、能量、位置等等
    • 确定目标,如能量沉积、能谱、分布等等
    • 物理过程,你能想到的Geant4基本搞定了,如有特别需要,去修改对应的代码

核心:我们是使用Geant4代码,不要造轮子,关键在于找到Geant4中需要的代码!!!

B1实例代码内容与作用

avatar

B1-构建几何模型-B1DetectorConstruction

  • 注意
    • 继承G4VUserDetectorConstruction
    • G4VPhysicalVolume* B1DetectorConstruction::Construct(){}中自己写实现
    • 有返回值,返回G4VPhysicalVolume

B1如何构建模型

  • World
    • 类似实验中的实验室坐标系,所有物体的位置均可设置成相对World的位置
    • vis.mac中World不可见:/vis/geometry/set/visibility World 0 false
    • 相关代码可以查看vis.mac,根据其代码查看Geant4说明书
  • Envelop
  • Shape1 and Shape2

  • 即World中有Envelop,Envelop中有两个Shape

每个模型都需要确定三个部分

  • solid volume

    • 作用:确定几何形状、尺寸
    • 调用类:G4Box
    • 语法:G4Box * [实例化name] = new G4Box(参数)
  • logcial volume

    • 作用:设置材料
    • 调用类:G4LogicalVolume
    • 语法:G4LogicalVolume * [实例化name] = new G4LogicalVolume(参数)
  • physical volume

    • 作用:物体的相对位置、旋转、平移
    • 调用类:G4VPhysicalVolume
    • 语法:G4VPhysicalVolume * [实例化name] = new G4PVPlacement(参数)
    • 注意:
      • 是否相对mother volume
      • G4VPhysicalVolume构造函数很多,查看说明书找到自己需要的
  • 相关具体代码以及使用说明请查看Geant4说明书

  • 知识点-如何判断粒子已经死亡

    • 粒子速度或能量为0
    • 手动设置条件,让粒子死亡
    • 粒子逃出介质,跑到了World中

B1-设置材料

  • 调用类:G4Material

    • Define a Material from the Geant4 Material Database-利用G4的材料库(B1用的此方法)

      • 具体内容查看说明书中Material中的Table
      • 示例:
      G4NistManager* man = G4NistManager::Instance();
      
      G4Material* H2O  = man->FindOrBuildMaterial("G4_WATER");
      G4Material* Air  = man->FindOrBuildMaterial("G4_AIR");
      
    • Define a Material from the Base Material-利用Material类自己定义

      • 示例:
      G4double density;
      
      density = 1.05*mg/cm3;
      G4Material* water1 = new G4Material("Water_1.05",density,"G4_WATER");
      
      density = 1.03*mg/cm3;
      G4NistManager* man = G4NistManager::Instance();
      G4Material* water2 = man->BuildMaterialWithNewDensity("Water_1.03","G4_WATER",density);
      
  • 调用类:G4Isotope(同位素)和G4Element(元素),根据需求自己定义

    • Define by number of atoms per each element in the molecule

      • 示例:
      a = 1.01*g/mole;
      G4Element* ele_H = new G4Element("Hydrogen",symbol="H",z=1.,a);
      a = 16.00*g/mole;
      G4Element* ele_O = new G4Element("Oxygen",symbol="O",z=8.,a);
      
      density = 1.000*g/cm3;
      G4Material* H2O = new G4Material("Water",density,ncomp=2);
      G4int natoms;
      H2O->AddElement(ele_H, natoms = 2);
      H2O->AddElement(ele_O, natoms = 1);
      
    • Define by fraction of mass of each element in the material

      • 示例:
      a = 14.01*g/mole;
      G4Element* ele_N = new G4Element(name="Nitrogen",symbol="N",z=7.,a);
      a = 16.00*g/mole;
      G4Element* ele_O = new G4Element(name="Oxygen",symbol="O",z=8.,a);
      
      density = 1.290*g/cm3;
      G4Material* Air = new G4Material(name="Air",density,ncomponents=2);
      G4double fracMass;
      Air->AddElement(ele_N, fracMass = 70.0*perCent);
      Air->AddElement(ele_O, fracMass = 30.0*perCent);
      
    • other
    //G4Isotope* [Isoname] = new G4Isotope(...);
    //Air->AddMaterial(...);
    
  • 可查看G4中自带的示例:examples/extended/electromagnetic/TestEm0/

  • 自己定义的材料G4会自动存入Material table中,并在模拟时调用

B1-设置粒子源(ParticleGun)

  • B1通过G4VUserPrimaryGeneratorAction实现粒子发射,代码为B1PrimaryGeneratorAction
    • 该代码继承:G4VUserPrimaryGeneratorAction
      • 动作类(action class)-如何理解?
      • 强制类,程序必须有
      • G4VUserPrimaryGeneratorAction注释提示:该类的作用并不是发射粒子的类
      • G4VUserPrimaryGeneratorAction注释提示:G4VPrimaryGenerator才是实现发射粒子的那个类,并且在G4VUserPrimaryGeneratorAction中GeneratePrimaries()函数实现
  • G4VPrimaryGenerator产生初级粒子的发射器有三个子类(子类的理解):
    • G4ParticleGun
    • G4GeneralParticleSource
    • G4HEPEvtInterface
  • B1PrimaryGeneratorAction.cc实现粒子发射

    • 设置粒子种类、动量方向、粒子能量
    • 设置粒子发射位置-加入了空间均匀分布
    • 调用fParticleGun->GeneratePrimaryVertex发射该粒子
  • B1中发射粒子代码

    fParticleGun->SetParticleDefinition(particle);//粒子种类
    fParticleGun->SetParticleMomentumDirection(G4ThreeVector(0.,0.,1.));//粒子方向
    fParticleGun->SetParticleEnergy(6.*MeV);//粒子能量
    fParticleGun->SetParticlePosition(G4ThreeVector(x0,y0,z0));//粒子发射位置
    
  • 在run1.mac中也可以修改发射粒子种类、能量,前提需要在代码中留出对应的代码功能!

B1-能量沉积

回顾

  • run, event, step, track
    • 一个run会有多个event
    • 每个event会有多个step,直至粒子消亡
    • 一个event的能量沉积相当于多个step的能损之和
    • 每个step的信息通过track获取
    • 进而可以通过循环求和得到粒子在介质中的能损

B1中通过Action class实现获取能量沉积

  • B1RunAction
    • 继承了G4UserRunAction
    • BeginOfRunAction()
    • EndOfRunAction()
    • AddEdep()
  • B1EventAction
    • 继承了G4UserEventAction
    • BeginOfEvnetAction()
    • EndOfEventAction()
    • AddEdep()
  • B1SteppingAction
    • 继承了G4UserSteppingAction
    • UserSteppingAction()
      • 判断粒子是否在logicvolume中:step->GerPreSetPoint()->GetTouchableHandle()->GetVolume()->GetLogicalVolume();
  • 注意
    • BeginOfRunAction() vs BeginOfEvnetAction()
      • BeginOfRunAction()是一个run调用一次
      • BeginOfEvnetAction()是一个Event调用一次
      • 如果一个run里有多个粒子,即多个evnet,那么一个BeginOfRunAction()会对应多个BeginOfEvnetAction()
    • 一个Event中每个step结束就会调用UserSteppingAction()
    • ybsim的能损实现是将step存入map中,然后在对map求和
    • B1RunAction->AddEdep() vs B1EventAction->AddEdep()
      • B1EventAction->AddEdep():多个step的Edep的叠加,单个粒子能损
      • B1RunAction->AddEdep():多个Event的Edep的叠加,多个粒子能损之和

B1-多线程

  • 线程 vs 进程
  • 例如:run 10 gammas
    • 单线程
      • 判断是否为有效step
      • 得到有效step的能量沉积,并且加和给EventAction.fEdep
      • EventAction.fEdep加和给RunAction.fEdep
      • 10个gamma逐次模拟
      • 得到最终的RunAction.fEdep
    • 多线程(难在同一过程,多线程分步后如何合并)
      • 同样上述步骤,分为两条线,一条线运行5个gamma,另一条线运行5个gamma
      • 对于EventAction.fEdep没有问题,自己和自己加和即可
      • 但两条线会有两个RunAction.fEdep,二者如何加和呢?
      • GetRunAction中对fEdep做了处理,即定义了两个量,一个为fEdep,一个为fEdep2
      • B1已经写好了合并代码(G4AccumulableManager* accumulableManager)
      • G4有个Merge功能,需要自己写代码将多个线程的变量合并

B1-动作类的初始化类

B1ActionInitialization

  • 继承了G4VUserActionInitialization
    • 其中的Build()用于实例化所有的动作类,会被G4RunManager调用
      • invoked by G4RunManager for sequential execution
      • invoked by G4WorkerRunManager for multi-threaded execution
    • 其中的BuildForMaster(),实例化UserRunAction,z只在多线程模式被调用(对于多线程,这个是必须的类)
      • invoked from G4MTRunManager for multi-threaded execution
    • 其中初始化函数有
      • SetUserAction(G4VUserPrimaryGeneratorAction*)
      • SetUserAction(G4UserRunAction*)
      • SetUserAction(G4UserEventAction*)
      • SetUserAction(G4UserStackingAction*)
      • SetUserAction(G4UserTrackingAction*)
      • SetUserAction(G4UserSteppingAction*)
    • 上述初始化函数,查看对应的.cc文件可以看到,就是调用的G4RunManager
  • 对比B1ActionInitialization::BuildForMaster()和B1ActionInitialization::Build()

exampleB1.cc

  • 调用runManager初始化B1DetectorConstruction()
    • G4RunManager是唯一的管理类
    • 在main函数中调用

B1-Physics list

建模过程中我们用户需要

  • 设置各种粒子属性
  • 设置探测器,介质等等
  • 确定所需的物理量
  • 设置Physicslist

Physicslist作用

  • 反应过程中每一步具体的动作在其中定义
  • 即每个step考虑什么物理过程,如何输运均在其中定义
  • B1中的Physicslist是调用的G4库中写好的
    • G4VModularPhysicsList* physicsList = new QBBC
    • 可在源码(.tar.gz解压)中查看G4中的PhysicsList或这官网查看Physics List Guide
  • 需要在main函数中定义
  • physicslist、physics process、process
    • process:一种具体反应过程,如瑞利散射
    • physics constructor:反应过程的集合
      • electromagnetic
      • hadronic
      • decay
      • photolepton-hadron
      • optical
      • parameterization
      • transportation
    • physicslist
      • Physics Model = final state generator
      • Physics Process = cross section + model
      • Physics List = list of processes for each particle

特殊需求需要自己去写physicslist

北京大学杨彪同学编写的SIM代码

注意

  • 该代码适用于北京大学实验核物理团队有关cluster实验中物理过程以及探测器阵列响应的模拟
  • 编写了课题组实验中常用的探测器类
    • 硅条探测器:W1,BB7,方硅,扇形硅
    • 闪烁体探测器:NaI, CsI,BGO,BC408
  • 物理过程包括:
    • 初级反应
    • 初级产物衰变或破裂
    • 不变质量分析
    • 等等
  • 根据保密、核不扩散条约、知识产权保护等原则,源代码不予公开,如有需要烦请联系北京大学实验核物理团队以及北京大学博士毕业生杨彪同学

代码结构(个人学习总结

  • 出于保密(调侃)、核不扩散条约(玩笑)、知识产权以及成果保护(真正的原因)等原因,不予公开
  • 如有需要烦请联系北京大学实验核物理团队以及北京大学博士毕业生杨彪同学

编译、运行

  • 复制完整的sim包(包括相关参数文件、数据文件等等)
  • 修改/sim/sim/CMakeLists.txt文件夹相关路径
    • 特别是yblib:set(YB_LIB_DIR /[filepath]/sim/yblib)
    • root库:set(ROOT_INCLUDE_DIR /opt/root61600/include)
  • 进入/[filepath]/sim/sim/build文件夹,执行cmake ..
    • 如若报错需要确认上步路径是否修改正确,或删除该文件夹下非G4的cmake生成的文件夹(可与B1对比)
  • 运行`make -j8',生成可执行文件sim

  • 运行./sim

    • ./sim phy:物理模拟, 生成文件存入/simfile/
    • ./sim vis:可视化
    • ./sim kin:运动学计算,生成文件存入/kinfile/

运行 ./sim vis结果如图: avatar

BUG

Error1:/usr/bin/ld: 找不到 -lMathMore(参考网站

  • 意思是找不到boost_serialization共享库,这个库的文件名应该为“libboost_serializatio.so”,其命名规则是:lib+库名(即xxx)+.so。
  • 用locate命令定位XXX库文件
    • locate libMathMore.so
    • 输出结果:/opt/root61600/lib/libMathMore.so
  • 再用软链接将两者链接起来
    • sudo ln -s /opt/root61600/lib/libMathMore.so /usr/lib/libMathMore.so

to be continued

In [2]:
!jupyter nbconvert geant4note_xi.ipynb --to html
[NbConvertApp] Converting notebook geant4note_xi.ipynb to html

[NbConvertApp] Writing 312784 bytes to geant4note_xi.html

In [ ]: