ATLAS Offline Software
LArG4ShowerLibSvc.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2025 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 #include "LArG4ShowerLibSvc.h"
6 
7 // shower lib classes
9 
11 
12 #include "AthenaKernel/Units.h"
13 
14 #include <sstream>
15 #include "TFile.h"
16 #include "TTree.h"
17 #include "Randomize.hh"
18 
19 namespace Units = Athena::Units;
20 
21 LArG4ShowerLibSvc::LArG4ShowerLibSvc(const std::string& name,ISvcLocator* svc)
22  : base_class(name,svc)
23 {
24  /* BE SURE THIS ONE IS THE SAME AS IN LArG4FastSimSvc!!! */
25  enum DETECTOR {EMB=100000,EMEC=200000,FCAL1=300000,FCAL2=400000,FCAL3=500000,HECLOC=600000,HEC=700000};
26 
27  m_detmap["EMB"]=EMB;
28  m_detmap["EMEC"]=EMEC;
29  m_detmap["FCAL1"]=FCAL1;
30  m_detmap["FCAL2"]=FCAL2;
31  m_detmap["FCAL3"]=FCAL3;
32  m_detmap["HECLOC"]=HECLOC;
33  m_detmap["HEC"]=HEC;
34 }
35 
37 {
38  ATH_MSG_INFO("Initializing");
39 
40  std::vector<std::string> ignoredLibraryFiles{};
41  // iterate through filenames in list
42  for (const std::string& fileName : m_fileNameList) {
43  std::string resolvedFilename = PathResolverFindCalibFile(fileName);
44  if (resolvedFilename.empty()) {
45  ATH_MSG_WARNING("Could not resolve input filename " << (fileName) << ". Ignoring!");
46  ignoredLibraryFiles.push_back(fileName);
47  continue;
48  } else {
49  ATH_MSG_INFO("Resolving input filename to " << resolvedFilename);
50  }
51 
52  TFile rootfile(resolvedFilename.c_str(),"READ");
53 
54  if (rootfile.IsZombie()) {
55  ATH_MSG_WARNING("File " << resolvedFilename << " is not a valid ROOT file");
56  ignoredLibraryFiles.push_back(fileName);
57  continue;
58  }
59 
60  const ShowerLib::IShowerLib* library = nullptr;
61 
62  // trying to create a library out of provided file
64 
65  // if no library can be created based on the file nullptr is returned
66  if (library == nullptr) {
67  ATH_MSG_WARNING("File " << resolvedFilename << " is not a valid library file");
68  ignoredLibraryFiles.push_back(fileName);
69  continue;
70  }
71 
72  if (m_detmap.find(library->detector()) == m_detmap.end()) {
73  ATH_MSG_WARNING("Library " << resolvedFilename << " is produced for unknown detector: " << library->detector());
74  ignoredLibraryFiles.push_back(fileName);
75  delete library;
76  continue;
77  }
78 
79  std::stringstream location;
80  location << library->detector() << "/" << library->particle_id();
81  int key = m_detmap.find(library->detector())->second + library->particle_id();
82 
83  // put the library into map
84  m_libraryMap[key] = library;
85  m_locations[key] = location.str();
86 
87  ATH_MSG_DEBUG("Filename: " << resolvedFilename.c_str());
88  ATH_MSG_DEBUG("Location: " << location.str());
89  ATH_MSG_DEBUG("Release: " << library->release().c_str());
90  ATH_MSG_DEBUG("Geant ver: " << library->geantVersion().c_str());
91  ATH_MSG_DEBUG("Phys list: " << library->physicsList().c_str());
92  }
93 
94  // no point in the service with no libraries
95  if (m_libraryMap.empty()) {
96  ATH_MSG_ERROR("No library files found. Please check the configuration of this job.");
97  return StatusCode::FAILURE;
98  }
99  // some of the configured library files did not work
100  if (!ignoredLibraryFiles.empty() ) {
101  for (const std::string& fileName : ignoredLibraryFiles) {
102  ATH_MSG_ERROR("Failed to create a library from filename: " << (fileName));
103  }
104  ATH_MSG_ERROR("Some library filenames were invalid. Please check the configuration of this job.");
105  return StatusCode::FAILURE;
106  }
107  ATH_MSG_INFO("List of loaded libraries:");
108  for (const auto& m : m_libraryMap) {
109  ATH_MSG_INFO(" " << m_locations[m.first] << ": " << m.second->comment());
110 #ifdef DEBUG_FrozenShowers
111  m_statisticsMap[m.second] = m.second->createStatistics();
112 #endif
113  }
114  ATH_MSG_INFO("Shower library successfully initialized.");
115 
116  return StatusCode::SUCCESS;
117 }
118 
120 {
121 
122  ATH_MSG_INFO("Finalizing shower library service.");
123 
124  libmap::const_iterator iter;
125 
126  for (iter = m_libraryMap.begin(); iter != m_libraryMap.end(); ++iter) {
127  ATH_MSG_INFO("Found ShowerLib at location " << m_locations[(*iter).first]);
128  ATH_MSG_INFO(std::endl << (*iter).second->statistics());
129  if (m_statisticsMap.find((*iter).second) != m_statisticsMap.end())
130  ATH_MSG_INFO(m_statisticsMap.find((*iter).second)->second->statistics());
131  else
132  ATH_MSG_INFO("No statistics available for this kind of library");
133  // delete the library:
134  delete (*iter).second ;
135  }
136 
137  return StatusCode::SUCCESS;
138 }
139 
140 /*
141  * Returns library from internal map based on the particle. Returns nullptr if there
142  * is no library for needed particle/detector.
143  */
144 const ShowerLib::IShowerLib* LArG4ShowerLibSvc::getShowerLib(G4int particleCode, int detectorTag) const
145 {
146  int location;
147 
148  location = detectorTag + particleCode;
149 
150  libmap::const_iterator iter = m_libraryMap.find(location);
151  if (iter != m_libraryMap.end()) {
152  return (*iter).second;
153  }
154 
155  return nullptr;
156 }
157 
158 bool
159 LArG4ShowerLibSvc::checkLibrary(G4int particleCode, int detectorTag)
160 {
161  const ShowerLib::IShowerLib* library = getShowerLib(particleCode, detectorTag);
162  if (library != nullptr) return true;
163 
164  return false;
165 }
166 
167 /*
168  * Returns a shower based on the particle. Return empty shower if there is no
169  * appropriate showers in the libraries.
170  */
171 
172 #ifdef DEBUG_FrozenShowers
173 std::vector<EnergySpot>
174 LArG4ShowerLibSvc::getShower(const G4FastTrack& track, int detectorTag)
175 #else
176 std::vector<EnergySpot>
177 LArG4ShowerLibSvc::getShower(const G4FastTrack& track, int detectorTag) const
178 #endif
179 {
180  // get shower lib from the map
181  const ShowerLib::IShowerLib* library = getShowerLib(track.GetPrimaryTrack()->GetDefinition()->GetPDGEncoding(), detectorTag);
182 
183  if (library == nullptr) {
184  // no library in map
185  ATH_MSG_ERROR("No library for location: " << detectorTag << "/" << track.GetPrimaryTrack()->GetDefinition()->GetPDGEncoding());
186  return std::vector<EnergySpot>();
187  }
188 
189  // starting point of the shower:
190  G4ThreeVector PositionShower = track.GetPrimaryTrack()->GetPosition();
191 
192  // get a shower from the library
193  int randomShift = 0;
194  randomShift = (int)(CLHEP::RandGaussZiggurat::shoot(G4Random::getTheEngine(), 0., 2.5)+0.5);
195 
196 #ifdef DEBUG_FrozenShowers
197  std::vector<EnergySpot>* shower = library->getShower(track.GetPrimaryTrack(), m_statisticsMap[library], randomShift);
198 #else
199  std::vector<EnergySpot>* shower = library->getShower(track.GetPrimaryTrack(), nullptr, randomShift);
200 #endif
201 
202 
203  if (shower == nullptr) {
204  return std::vector<EnergySpot>();
205  }
206 
207  // axis of the shower, in global reference frame:
208  G4ThreeVector DirectionShower = track.GetPrimaryTrack()->GetMomentumDirection();
209 
210  // time of the track (as in LArBarrelCalculator.cc)
211  G4double tof = track.GetPrimaryTrack()->GetGlobalTime() / Units::ns;
212  G4double time = tof - ( track.GetPrimaryTrack()->GetPosition().mag()/Units::c_light ) / Units::ns;
213 
215 
216  // Create energy spots
217  for (hit = shower->begin(); hit != shower->end(); ++hit) {
218  G4ThreeVector hitpos = (*hit).GetPosition();
219  hitpos.rotateUz(DirectionShower); // rotate the hit to the needed direction
220  (*hit).SetPosition(hitpos+PositionShower);
221  (*hit).SetTime((*hit).GetTime() + time);
222  }
223 
224  if (msgSvc()->outputLevel(name()) <= MSG::VERBOSE) {
225  ATH_MSG_VERBOSE("Prepared " << shower->size() << " EnergySpot");
226  for (std::vector<EnergySpot>::const_iterator iter = shower->begin(); iter != shower->end(); ++iter) {
227  ATH_MSG_VERBOSE("EnergySpot: " << iter->GetPosition().x() << " " << iter->GetPosition().y() << " " << iter->GetPosition().z()
228  << " " << iter->GetEnergy() << " " << iter->GetTime());
229  }
230  }
231 
232  return *std::unique_ptr< std::vector<EnergySpot> >(shower);
233 }
234 
235 double
236 LArG4ShowerLibSvc::getContainmentZ(const G4FastTrack& track, int detectorTag)
237 {
238  // get shower lib from the map
239  const ShowerLib::IShowerLib* library = getShowerLib(track.GetPrimaryTrack()->GetDefinition()->GetPDGEncoding(), detectorTag);
240 
241  if (library == nullptr) {
242  return 0.0;
243  }
244 
245  return library->getContainmentZ(track.GetPrimaryTrack());
246 }
247 
248 double
249 LArG4ShowerLibSvc::getContainmentR(const G4FastTrack& track, int detectorTag)
250 {
251  // get shower lib from the map
252  const ShowerLib::IShowerLib* library = getShowerLib(track.GetPrimaryTrack()->GetDefinition()->GetPDGEncoding(), detectorTag);
253 
254  if (library == nullptr) {
255  return 0.0;
256  }
257 
258  return library->getContainmentR(track.GetPrimaryTrack());
259 }
xAOD::iterator
JetConstituentVector::iterator iterator
Definition: JetConstituentVector.cxx:68
ShowerLib::IShowerLib::getContainmentZ
virtual double getContainmentZ(const G4Track *track) const =0
get average length of showers for the given energy
ShowerLib::IShowerLib::release
virtual const std::string release() const
get Release tag
Definition: IShowerLib.h:127
createLinkingScheme.iter
iter
Definition: createLinkingScheme.py:62
LArSamples::HEC
@ HEC
Definition: CaloId.h:26
LArG4ShowerLibSvc::getShower
virtual std::vector< EnergySpot > getShower(const G4FastTrack &track, int detectorTag) const override
return list of energy depositions for given track (interface implementation)
Definition: LArG4ShowerLibSvc.cxx:177
LArG4ShowerLibSvc::finalize
virtual StatusCode finalize() override
Definition: LArG4ShowerLibSvc.cxx:119
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition: AthMsgStreamMacros.h:31
ShowerLib::IShowerLib::getContainmentR
virtual double getContainmentR(const G4Track *track) const =0
get average lateral spread of the showers for the given energy
LArG4ShowerLibSvc.h
ShowerLib::IShowerLib::detector
virtual const std::string detector() const
get detector tag
Definition: IShowerLib.h:121
LArG4ShowerLibSvc::LArG4ShowerLibSvc
LArG4ShowerLibSvc(const std::string &name, ISvcLocator *svc)
Definition: LArG4ShowerLibSvc.cxx:21
LArG4ShowerLibSvc::getShowerLib
const ShowerLib::IShowerLib * getShowerLib(G4int particleCode, int detectorTag) const
get shower library from StoreGate by track (using current volume name)
Definition: LArG4ShowerLibSvc.cxx:144
ShowerLib::iterateTTree
IShowerLib * iterateTTree(TFile *fname)
Definition: ShowerLibList.cxx:55
CaloCell_ID_FCS::FCAL1
@ FCAL1
Definition: FastCaloSim_CaloCell_ID.h:41
LArG4ShowerLibSvc::m_libraryMap
libmap m_libraryMap
mapping StoreGate key to handle in StoreGate
Definition: LArG4ShowerLibSvc.h:65
ATH_MSG_VERBOSE
#define ATH_MSG_VERBOSE(x)
Definition: AthMsgStreamMacros.h:28
LArCalib_HVScale2NtupleConfig.rootfile
string rootfile
Definition: LArCalib_HVScale2NtupleConfig.py:83
ShowerLib::IShowerLib::particle_id
virtual int particle_id() const
get particle tag
Definition: IShowerLib.h:124
LArSamples::EMEC
@ EMEC
Definition: CaloId.h:25
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition: AthMsgStreamMacros.h:33
LArSamples::EMB
@ EMB
Definition: CaloId.h:25
StdJOSetup.msgSvc
msgSvc
Provide convenience handles for various services.
Definition: StdJOSetup.py:36
LArG4ShowerLibSvc::checkLibrary
virtual bool checkLibrary(G4int particleCode, int detectorTag) override
Definition: LArG4ShowerLibSvc.cxx:159
LArG4ShowerLibSvc::m_detmap
std::map< std::string, int > m_detmap
Definition: LArG4ShowerLibSvc.h:69
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
LArG4ShowerLibSvc::m_locations
std::map< int, std::string > m_locations
Definition: LArG4ShowerLibSvc.h:68
LArG4ShowerLibSvc::m_statisticsMap
statmap m_statisticsMap
Definition: LArG4ShowerLibSvc.h:67
ShowerLib::IShowerLib
Class for shower library shower lib interface.
Definition: IShowerLib.h:40
ShowerLib::IShowerLib::getShower
virtual std::vector< EnergySpot > * getShower(const G4Track *track, ShowerLibStatistics *stats, int randomShift) const =0
get shower for given G4 track
Handler::svc
AthROOTErrorHandlerSvc * svc
Definition: AthROOTErrorHandlerSvc.cxx:10
LArG4ShowerLibSvc::initialize
virtual StatusCode initialize() override
Definition: LArG4ShowerLibSvc.cxx:36
Athena::Units
Definition: Units.h:45
PathResolver.h
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:240
LArG4ShowerLibSvc::m_fileNameList
StringArrayProperty m_fileNameList
property, list of library files
Definition: LArG4ShowerLibSvc.h:71
ShowerLibList.h
LArG4ShowerLibSvc::getContainmentZ
virtual double getContainmentZ(const G4FastTrack &track, int detectorTag) override
Definition: LArG4ShowerLibSvc.cxx:236
python.PhysicalConstants.c_light
float c_light
Definition: PhysicalConstants.py:73
Units.h
Wrapper to avoid constant divisions when using units.
PathResolverFindCalibFile
std::string PathResolverFindCalibFile(const std::string &logical_file_name)
Definition: PathResolver.cxx:283
python.CaloAddPedShiftConfig.int
int
Definition: CaloAddPedShiftConfig.py:45
CaloSwCorrections.time
def time(flags, cells_name, *args, **kw)
Definition: CaloSwCorrections.py:242
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition: AthMsgStreamMacros.h:32
LArG4ShowerLibSvc::getContainmentR
virtual double getContainmentR(const G4FastTrack &track, int detectorTag) override
Definition: LArG4ShowerLibSvc.cxx:249
ShowerLib::IShowerLib::geantVersion
virtual const std::string geantVersion() const
get geant version tag
Definition: IShowerLib.h:133
CaloCell_ID_FCS::FCAL2
@ FCAL2
Definition: FastCaloSim_CaloCell_ID.h:42
xAOD::track
@ track
Definition: TrackingPrimitives.h:513
ShowerLib::IShowerLib::physicsList
virtual const std::string physicsList() const
get geant 4 physics list name
Definition: IShowerLib.h:136
python.Constants.VERBOSE
int VERBOSE
Definition: Control/AthenaCommon/python/Constants.py:13
jobOptions.fileName
fileName
Definition: jobOptions.SuperChic_ALP2.py:39
python.SystemOfUnits.ns
float ns
Definition: SystemOfUnits.py:146
python.SystemOfUnits.m
float m
Definition: SystemOfUnits.py:106
mapkey::key
key
Definition: TElectronEfficiencyCorrectionTool.cxx:37