CbmRoot
CbmBeamGenerator.cxx
Go to the documentation of this file.
1 
7 #include "CbmBeamGenerator.h"
8 
9 #include "CbmEventGenerator.h"
10 #include <TDatabasePDG.h>
11 #include <TParticlePDG.h>
12 
13 #include "FairLogger.h"
14 #include "FairMCEventHeader.h"
15 #include "FairPrimaryGenerator.h"
16 #include "FairRunSim.h"
17 #include "TFile.h"
18 #include "TRandom.h"
19 #include "TTree.h"
20 #include "TVector3.h"
21 #include "UEvent.h"
22 #include "UParticle.h"
23 #include "URun.h"
24 #include <cassert>
25 #include <sstream>
26 
27 
28 // ----- Default constructor --------------------------------------------
30  : FairGenerator("BeamGenerator", "CBM generator")
31  , fP(0.)
32  , fStartZ(0.)
33  , fIon(nullptr) {}
34 // --------------------------------------------------------------------------
35 
36 
37 // ----- Standard constructor --------------------------------------------
39  UInt_t beamA,
40  UInt_t beamQ,
41  Double_t momentum,
42  Double_t startZ)
43  : FairGenerator("BeamGenerator", "CBM generator")
44  , fP(momentum * Double_t(beamA))
45  , fStartZ(startZ)
46  , fIon(nullptr) {
47 
48  // --- Create the ion species and add it to the particle list
49  char name[20];
50  sprintf(name, "Beam_%d_%d_%d", beamZ, beamA, beamQ);
51  fIon = new FairIon(name, beamZ, beamA, beamQ);
52  FairRunSim* run = FairRunSim::Instance();
53  assert(run);
54  run->AddNewIon(fIon);
55 
56  // --- Tell the event generator to force the vertex at the given z
57  FairPrimaryGenerator* primGen = run->GetPrimaryGenerator();
58  CbmEventGenerator* evGen = dynamic_cast<CbmEventGenerator*>(primGen);
59  assert(evGen);
60  evGen->ForceVertexAtZ(startZ);
61  evGen->SmearVertexZ(kFALSE);
62 }
63 // --------------------------------------------------------------------------
64 
65 
66 // ----- Destructor -----------------------------------------------------
68 // --------------------------------------------------------------------------
69 
70 
71 // ----- Print info -----------------------------------------------------
72 void CbmBeamGenerator::Print(Option_t*) const { LOG(info) << ToString(); }
73 // --------------------------------------------------------------------------
74 
75 
76 // ----- Generate input event -------------------------------------------
77 Bool_t CbmBeamGenerator::ReadEvent(FairPrimaryGenerator* primGen) {
78 
79  TParticlePDG* ion = TDatabasePDG::Instance()->GetParticle(fIon->GetName());
80  assert(ion);
81 
82  // Note that the beam position in x and y and the beam angle are generated
83  // by CbmEventGenerator. Here, we generate the ion thus in z direction
84  // and at zero coordinates. It will be properly placed and rotated
85  // according to the beam profile specified to CbmEventGenerator.
86  primGen->AddTrack(ion->PdgCode(), 0., 0., fP, 0., 0., 0.);
87 
88  return kTRUE;
89 }
90 // --------------------------------------------------------------------------
91 
92 
93 // ----- Info -----------------------------------------------------------
94 std::string CbmBeamGenerator::ToString() const {
95 
96  std::stringstream ss;
97  ss << GetName() << " ion: " << fIon->GetName() << " Z " << fIon->GetZ()
98  << " A " << fIon->GetA() << " Q " << fIon->GetQ() << " mass "
99  << fIon->GetMass() << ", momentum " << fP
100  << " GeV/u, start Z = " << fStartZ;
101 
102  return ss.str();
103 }
104 // --------------------------------------------------------------------------
105 
106 
CbmBeamGenerator::~CbmBeamGenerator
virtual ~CbmBeamGenerator()
Destructor.
Definition: CbmBeamGenerator.cxx:67
CbmBeamGenerator::fStartZ
Double_t fStartZ
z coordinate of start point
Definition: CbmBeamGenerator.h:79
CbmBeamGenerator
Definition: CbmBeamGenerator.h:38
ClassImp
ClassImp(CbmBeamGenerator)
CbmBeamGenerator.h
CbmEventGenerator
CbmEventGenerator.
Definition: CbmEventGenerator.h:34
UParticle.h
CbmBeamGenerator::fP
Double_t fP
Total momentum [GeV].
Definition: CbmBeamGenerator.h:78
CbmBeamGenerator::ReadEvent
virtual Bool_t ReadEvent(FairPrimaryGenerator *primGen)
Generate one event (abstract in base class)
Definition: CbmBeamGenerator.cxx:77
UEvent.h
CbmEventGenerator::ForceVertexAtZ
void ForceVertexAtZ(Double_t zVertex)
Force event vertex to be at a given z.
Definition: CbmEventGenerator.cxx:39
CbmBeamGenerator::fIon
FairIon * fIon
Ion type.
Definition: CbmBeamGenerator.h:80
CbmBeamGenerator::ToString
std::string ToString() const
Info to string.
Definition: CbmBeamGenerator.cxx:94
CbmBeamGenerator::Print
virtual void Print(Option_t *opt="") const
Print info to logger.
Definition: CbmBeamGenerator.cxx:72
URun.h
CbmEventGenerator.h
CbmBeamGenerator::CbmBeamGenerator
CbmBeamGenerator()
Default constructor (should not be used)
Definition: CbmBeamGenerator.cxx:29