zoukankan      html  css  js  c++  java
  • Geant4 入射粒子设置

    Geant4 入射粒子设置


    在 PrimaryGeneratorAction 中设置 G4ParticleGun.

    注意:关于粒子种类的定义,有两个位置可以定义,第一是在 PrimaryGeneratorAction 类的构造函数中定义,第二是在 GeneratePrimaries() 函数中定义。

    • 在 PrimaryGeneratorAction 类的构造函数中定义,粒子种类可以在 .mac 文件中修改,举例代码取自 example B2:
     1 B2PrimaryGeneratorAction::B2PrimaryGeneratorAction()
     2  : G4VUserPrimaryGeneratorAction()
     3 {
     4   G4int nofParticles = 1;
     5   fParticleGun = new G4ParticleGun(nofParticles);
     6 
     7   // default particle kinematic
     8 
     9   G4ParticleDefinition* particleDefinition 
    10     = G4ParticleTable::GetParticleTable()->FindParticle("proton");
    11 
    12   fParticleGun->SetParticleDefinition(particleDefinition);
    13   fParticleGun->SetParticleMomentumDirection(G4ThreeVector(0.,0.,1.));
    14   fParticleGun->SetParticleEnergy(3.0*GeV);
    15 }
    • 在 GeneratePrimaries() 函数中定义,粒子种类可以在 .mac 文件中修改,举例代码取自 example B5:
      1 #include "B5PrimaryGeneratorAction.hh"
      2 
      3 #include "G4Event.hh"
      4 #include "G4ParticleGun.hh"
      5 #include "G4ParticleTable.hh"
      6 #include "G4ParticleDefinition.hh"
      7 #include "G4GenericMessenger.hh"
      8 #include "G4SystemOfUnits.hh"
      9 #include "Randomize.hh"
     10 
     11 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
     12 
     13 B5PrimaryGeneratorAction::B5PrimaryGeneratorAction()
     14 : G4VUserPrimaryGeneratorAction(),     
     15   fParticleGun(nullptr), fMessenger(nullptr), 
     16   fPositron(nullptr), fMuon(nullptr), fPion(nullptr), 
     17   fKaon(nullptr), fProton(nullptr),
     18   fMomentum(1000.*MeV),
     19   fSigmaMomentum(50.*MeV),
     20   fSigmaAngle(2.*deg),
     21   fRandomizePrimary(true)
     22 {
     23   G4int nofParticles = 1;
     24   fParticleGun  = new G4ParticleGun(nofParticles);
     25   
     26   auto particleTable = G4ParticleTable::GetParticleTable();
     27   fPositron = particleTable->FindParticle("e+");
     28   fMuon = particleTable->FindParticle("mu+");
     29   fPion = particleTable->FindParticle("pi+");
     30   fKaon = particleTable->FindParticle("kaon+");
     31   fProton = particleTable->FindParticle("proton");
     32   
     33   // default particle kinematics
     34   fParticleGun->SetParticlePosition(G4ThreeVector(0.,0.,-8.*m));
     35   fParticleGun->SetParticleDefinition(fPositron);
     36   
     37   // define commands for this class
     38   DefineCommands();
     39 }
     40 
     41 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
     42 
     43 B5PrimaryGeneratorAction::~B5PrimaryGeneratorAction()
     44 {
     45   delete fParticleGun;
     46   delete fMessenger;
     47 }
     48 
     49 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
     50 
     51 void B5PrimaryGeneratorAction::GeneratePrimaries(G4Event* event)
     52 {
     53   G4ParticleDefinition* particle;  
     54   if (fRandomizePrimary) {
     55     G4int i = (int)(5.*G4UniformRand());
     56     switch(i) {
     57       case 0:
     58           particle = fPositron;
     59           break;
     60       case 1:
     61           particle = fMuon;
     62           break;
     63       case 2:
     64           particle = fPion;
     65           break;
     66       case 3:
     67           particle = fKaon;
     68           break;
     69       default:
     70           particle = fProton;
     71           break;
     72     }
     73     fParticleGun->SetParticleDefinition(particle);
     74   }
     75   else {
     76     particle = fParticleGun->GetParticleDefinition();
     77   }
     78   
     79   auto pp = fMomentum + (G4UniformRand()-0.5)*fSigmaMomentum;
     80   auto mass = particle->GetPDGMass();
     81   auto ekin = std::sqrt(pp*pp+mass*mass)-mass;
     82   fParticleGun->SetParticleEnergy(ekin);
     83   
     84   auto angle = (G4UniformRand()-0.5)*fSigmaAngle;
     85   fParticleGun->SetParticleMomentumDirection(
     86                   G4ThreeVector(std::sin(angle),0.,std::cos(angle)));
     87   
     88   fParticleGun->GeneratePrimaryVertex(event);
     89 }
     90 
     91 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
     92 
     93 void B5PrimaryGeneratorAction::DefineCommands()
     94 {
     95   // Define /B5/generator command directory using generic messenger class
     96   fMessenger 
     97     = new G4GenericMessenger(this, 
     98                              "/B5/generator/", 
     99                              "Primary generator control");
    100             
    101   // momentum command
    102   auto& momentumCmd
    103     = fMessenger->DeclarePropertyWithUnit("momentum", "GeV", fMomentum, 
    104         "Mean momentum of primaries.");
    105   momentumCmd.SetParameterName("p", true);
    106   momentumCmd.SetRange("p>=0.");                                
    107   momentumCmd.SetDefaultValue("1.");
    108   // ok
    109   //momentumCmd.SetParameterName("p", true);
    110   //momentumCmd.SetRange("p>=0.");                                
    111   
    112   // sigmaMomentum command
    113   auto& sigmaMomentumCmd
    114     = fMessenger->DeclarePropertyWithUnit("sigmaMomentum",
    115         "MeV", fSigmaMomentum, "Sigma momentum of primaries.");
    116   sigmaMomentumCmd.SetParameterName("sp", true);
    117   sigmaMomentumCmd.SetRange("sp>=0.");                                
    118   sigmaMomentumCmd.SetDefaultValue("50.");
    119 
    120   // sigmaAngle command
    121   auto& sigmaAngleCmd
    122     = fMessenger->DeclarePropertyWithUnit("sigmaAngle", "deg", fSigmaAngle, 
    123         "Sigma angle divergence of primaries.");
    124   sigmaAngleCmd.SetParameterName("t", true);
    125   sigmaAngleCmd.SetRange("t>=0.");                                
    126   sigmaAngleCmd.SetDefaultValue("2.");
    127 
    128   // randomizePrimary command
    129   auto& randomCmd
    130     = fMessenger->DeclareProperty("randomizePrimary", fRandomizePrimary);
    131   G4String guidance
    132     = "Boolean flag for randomizing primary particle types.
    ";   
    133   guidance
    134     += "In case this flag is false, you can select the primary particle
    ";
    135   guidance += "  with /gun/particle command.";                               
    136   randomCmd.SetGuidance(guidance);
    137   randomCmd.SetParameterName("flg", true);
    138   randomCmd.SetDefaultValue("true");
    139 }
    140 
    141 //..oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
    • 第 54 - 72 行,粒子源随机出射一种粒子。
    • DefineCommands() 函数中定义了 macro 文件中使用的一些命令,有待研究

    可以使用 G4GeneralParticleSource 替代 G4ParticleGun

    • “For many applications G4ParticleGun is a suitable particle generator. However if you want to generate primary particles in more sophisticated manner, you can utilize G4GeneralParticleSource, the GEANT4 General Particle Source module (GPS), discussed in the next section ( General Particle Source).” Geant4 手册原话。

    代码举例,取自 HPGe_60Co

    PrimaryGeneratorAction.hh 文件

     1 #ifndef PrimaryGeneratorAction_h
     2 #define PrimaryGeneratorAction_h 1
     3 
     4 #include "G4VUserPrimaryGeneratorAction.hh"
     5 #include "globals.hh"
     6 
     7 class G4GeneralParticleSource;
     8 class G4Event;
     9 
    10 class PrimaryGeneratorAction : public G4VUserPrimaryGeneratorAction
    11 {
    12 public:
    13   PrimaryGeneratorAction();    
    14   virtual ~PrimaryGeneratorAction();
    15 
    16   virtual void GeneratePrimaries(G4Event* );
    17 
    18   const G4GeneralParticleSource* GetParticleGun() const {return fParticleGun;}
    19 
    20   // Set methods
    21   void SetRandomFlag(G4bool );
    22 
    23 private:
    24   G4GeneralParticleSource*          fParticleGun; // G4 particle gun
    25 };
    26 
    27 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
    28 
    29 #endif
     1 #include "PrimaryGeneratorAction.hh"
     2 
     3 #include "G4LogicalVolumeStore.hh"
     4 #include "G4LogicalVolume.hh"
     5 #include "G4Box.hh"
     6 #include "G4Event.hh"
     7 #include "G4ParticleGun.hh"
     8 #include "G4GeneralParticleSource.hh"
     9 #include "G4ParticleTable.hh"
    10 #include "G4ParticleDefinition.hh"
    11 #include "G4SystemOfUnits.hh"
    12 #include "G4RandomDirection.hh"
    13 #include "G4IonTable.hh"
    14 #include "G4Geantino.hh"
    15 
    16 #include "Randomize.hh"
    17 
    18 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
    19 
    20 PrimaryGeneratorAction::PrimaryGeneratorAction()
    21  : G4VUserPrimaryGeneratorAction()
    22 {
    23   fParticleGun = new G4GeneralParticleSource();
    24 }
    25 
    26 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
    27 
    28 PrimaryGeneratorAction::~PrimaryGeneratorAction()
    29 {
    30   delete fParticleGun;
    31 }
    32 
    33 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
    34 
    35 void PrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent)
    36 {
    37   fParticleGun->GeneratePrimaryVertex(anEvent);
    38 }
    39 
    40 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
    • 使用 G4GeneralParticleSource 的 PrimaryGeneratorAction 很简单。
    1 /run/initialize
    2 /tracking/verbose 0
    3 
    4 
    5 /gps/particle ion
    6 /gps/ion 27 60 0 0
    7 /gps/pos/centre 0 0 0 mm
    8 
    9 /run/beamOn 10000
    • 有关 gps的 macro 文件

    得到 60Co 的相关能谱,这个能谱比较好理解。

    得到 241Am 能谱,这个就不好理解,241Am alpha 能量约为 5.486 MeV,但得到的能谱却很奇怪。

    我认为这是我对使用 gps 产生 alpha 源不熟悉导致的,在 G4 的论坛中有类似的问题,回答中提供的图也有类似的问题,我不理解。https://geant4-forum.web.cern.ch/top/monthly

  • 相关阅读:
    北风设计模式课程---单一职责原则
    北风设计模式课程---访问者模式(Visitor)
    北风设计模式课程---访问者(Visitor)模式
    北风设计模式课程---命令模式
    社交专题---3、装逼的那些事儿
    社交专题---2、男生的套路有多深
    MHA 日常管理
    【屌丝程序的口才逆袭演讲稿50篇】第十篇:程序猿们请看看外面的世界吧【张振华.Jack】
    Win32 Windows编程 九
    java实现大数相加问题
  • 原文地址:https://www.cnblogs.com/kurrrr/p/13673142.html
Copyright © 2011-2022 走看看