FairRoot
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
FairShieldGenerator.cxx
Go to the documentation of this file.
1 /********************************************************************************
2  * Copyright (C) 2014 GSI Helmholtzzentrum fuer Schwerionenforschung GmbH *
3  * *
4  * This software is distributed under the terms of the *
5  * GNU Lesser General Public Licence (LGPL) version 3, *
6  * copied verbatim in the file "LICENSE" *
7  ********************************************************************************/
8 // -------------------------------------------------------------------------
9 // ----- FairShieldGenerator source file -----
10 // ----- Created 15/09/06 by V. Friese -----
11 // -------------------------------------------------------------------------
12 #include "FairShieldGenerator.h"
13 
14 #include "FairIon.h" // for FairIon
15 #include "FairLogger.h" // for logging
16 #include "FairPrimaryGenerator.h" // for FairPrimaryGenerator
17 #include "FairRunSim.h" // for FairRunSim
18 
19 #include <TDatabasePDG.h> // for TDatabasePDG
20 #include <TParticlePDG.h> // for TParticlePDG
21 #include <climits> // for INT_MAX
22 #include <cstdio> // for sprintf
23 #include <fstream> // for ifstream
24 #include <utility> // for pair
25 
27  : FairGenerator()
28  , fInputFile(nullptr)
29  , fFileName(nullptr)
30  , fPDG(nullptr)
31  , fIonMap()
32 {}
33 
35  : FairGenerator()
36  , fInputFile(nullptr)
37  , fFileName(fileName)
38  , fPDG(TDatabasePDG::Instance())
39  , fIonMap()
40 {
41  LOG(info) << "FairShieldGenerator: Opening input file " << fileName;
42  fInputFile = new std::ifstream(fFileName);
43  if (!fInputFile->is_open()) {
44  LOG(fatal) << "Cannot open input file.";
45  }
46  LOG(info) << "FairShieldGenerator: Looking for ions...";
47  Int_t nIons = RegisterIons();
48  LOG(info) << "FairShieldGenerator: " << nIons << " ions registered.";
49  CloseInput();
50  LOG(info) << "FairShieldGenerator: Reopening input file " << fileName;
51  fInputFile = new std::ifstream(fFileName);
52 }
53 
55 
57 {
58  // Check for input file
59  if (!fInputFile->is_open()) {
60  LOG(error) << "FairShieldGenerator: Input file not open!";
61  return kFALSE;
62  }
63 
64  // Define event variable to be read from file
65  Int_t eventId = 0;
66  Int_t nTracks = 0;
67  Double_t pBeam = 0.;
68  Double_t b = 0.;
69 
70  // Define track variables to be read from file
71  Int_t iPid = 0;
72  Int_t iMass = 0;
73  Int_t iCharge = 0;
74  Double_t px = 0.;
75  Double_t py = 0.;
76  Double_t pz = 0.;
77  Int_t pdgType = 0;
78 
79  // Read event header line from input file
80  *fInputFile >> eventId;
81  *fInputFile >> nTracks;
82  if (nTracks < 0 || nTracks > (INT_MAX - 1))
83  LOG(fatal) << "Error reading the number of events from event header.";
84  *fInputFile >> pBeam >> b;
85 
86  // If end of input file is reached : close it and abort run
87  if (fInputFile->eof()) {
88  LOG(info) << "FairShieldGenerator: End of input file reached ";
89  CloseInput();
90  return kFALSE;
91  }
92 
93  LOG(info) << "FairShieldGenerator: Event " << eventId << ", pBeam = " << pBeam << "GeV, b = " << b
94  << " fm, multiplicity " << nTracks;
95 
96  // Loop over tracks in the current event
97  for (Int_t itrack = 0; itrack < nTracks; itrack++) {
98 
99  // Read PID and momentum from file
100  *fInputFile >> iPid >> iMass >> iCharge >> px >> py >> pz;
101 
102  // Case ion
103  if (iPid == 1000) {
104  char ionName[20];
105  sprintf(ionName, "Ion_%d_%d", iMass, iCharge);
106  TParticlePDG* part = fPDG->GetParticle(ionName);
107  if (!part) {
108  LOG(warn) << "FairShieldGenerator::ReadEvent: Cannot find " << ionName << " in database!";
109  continue;
110  }
111  pdgType = part->PdgCode();
112  } else {
113  pdgType = fPDG->ConvertGeant3ToPdg(iPid);
114  } // "normal" particle
115 
116  // Give track to PrimaryGenerator
117  primGen->AddTrack(pdgType, px, py, pz, 0., 0., 0.);
118  }
119 
120  return kTRUE;
121 }
122 
123 void FairShieldGenerator::CloseInput()
124 {
125  if (fInputFile) {
126  if (fInputFile->is_open()) {
127  LOG(info) << "FairShieldGenerator: Closing input file " << fFileName;
128  fInputFile->close();
129  }
130  delete fInputFile;
131  fInputFile = nullptr;
132  }
133 }
134 
135 Int_t FairShieldGenerator::RegisterIons()
136 {
137  Int_t nIons = 0;
138  Int_t eventId, nTracks, iPid, iMass, iCharge;
139  Double_t pBeam, b, px, py, pz;
140  fIonMap.clear();
141 
142  while (!fInputFile->eof()) {
143 
144  *fInputFile >> eventId;
145  *fInputFile >> nTracks;
146  if (nTracks < 0 || nTracks > (INT_MAX - 1))
147  LOG(fatal) << "Error reading the number of events from event header.";
148  *fInputFile >> pBeam >> b;
149  if (fInputFile->eof()) {
150  continue;
151  }
152  for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
153  *fInputFile >> iPid >> iMass >> iCharge >> px >> py >> pz;
154  if (iPid == 1000) { // ion
155  char buffer[20];
156  sprintf(buffer, "Ion_%d_%d", iMass, iCharge);
157  TString ionName(buffer);
158  if (fIonMap.find(ionName) == fIonMap.end()) { // new ion
159  FairIon* ion = new FairIon(ionName, iCharge, iMass, iCharge);
160  fIonMap[ionName] = ion;
161  nIons++;
162  } // new ion
163  } // ion
164  } // track loop
165  } // event loop
166 
168  for (const auto& mi : fIonMap) {
169  FairIon* ion = mi.second;
170  run->AddNewIon(ion);
171  }
172 
173  return nIons;
174 }
175 
static FairRunSim * Instance()
Definition: FairRunSim.cxx:116
ClassImp(FairEventBuilder)
virtual Bool_t ReadEvent(FairPrimaryGenerator *primGen)
virtual void AddTrack(Int_t pdgid, Double_t px, Double_t py, Double_t pz, Double_t vx, Double_t vy, Double_t vz, Int_t parent=-1, Bool_t wanttracking=true, Double_t e=-9e9, Double_t tof=0., Double_t weight=0., TMCProcess proc=kPPrimary)
void AddNewIon(FairIon *ion)
Definition: FairRunSim.h:51