ATLAS Offline Software
PixelRadSimFluenceMapAlg.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
3 */
4 
6 #include "GaudiKernel/EventIDRange.h"
7 #include "TFile.h"
8 #include <memory>
9 PixelRadSimFluenceMapAlg::PixelRadSimFluenceMapAlg(const std::string& name, ISvcLocator* pSvcLocator):
10  ::AthReentrantAlgorithm(name, pSvcLocator)
11 {
12 }
13 
15  ATH_MSG_DEBUG("PixelRadSimFluenceMapAlg::initialize()");
18  return StatusCode::SUCCESS;
19 }
20 
21 StatusCode PixelRadSimFluenceMapAlg::execute(const EventContext& ctx) const {
22  ATH_MSG_DEBUG("PixelRadSimFluenceMapAlg::execute()");
23 
25  if (writeFluenceMapHandle.isValid()) {
26  ATH_MSG_DEBUG("CondHandle " << writeFluenceMapHandle.fullKey() << " is already valid.. In theory this should not be called, but may happen if multiple concurrent events are being processed out of order.");
27  return StatusCode::SUCCESS;
28  }
29 
30  // Construct the output Cond Object and fill it in
31  std::unique_ptr<PixelRadiationDamageFluenceMapData> writeFluenceCdo(std::make_unique<PixelRadiationDamageFluenceMapData>());
32 
33  const EventIDBase start{EventIDBase::UNDEFNUM, EventIDBase::UNDEFEVT, 0,
34  0, EventIDBase::UNDEFNUM, EventIDBase::UNDEFNUM};
35  const EventIDBase stop {EventIDBase::UNDEFNUM, EventIDBase::UNDEFEVT, EventIDBase::UNDEFNUM-1,
36  EventIDBase::UNDEFNUM-1, EventIDBase::UNDEFNUM, EventIDBase::UNDEFNUM};
37 
38  EventIDRange rangeW{start, stop};
39 
41  const PixelModuleData *moduleData = *moduleDataHandle;
42 
43  // mapping files for radiation damage simulation
44  writeFluenceCdo -> setFluenceLayer(moduleData->getFluenceLayer());
45  writeFluenceCdo -> setFluenceLayer3D(moduleData->getFluenceLayer3D());
46 
47  // Create mapping file for radiation damage simulation
48  std::vector<PixelHistoConverter> ramoPotentialMap;
49  std::vector<PixelHistoConverter> lorentzMap_e;
50  std::vector<PixelHistoConverter> lorentzMap_h;
51  std::vector<PixelHistoConverter> distanceMap_e;
52  std::vector<PixelHistoConverter> distanceMap_h;
53  for (unsigned int i=0; i<moduleData->getRadSimFluenceMapList().size(); i++) {
54  ATH_MSG_INFO("Using maps located in: "<< moduleData->getRadSimFluenceMapList().at(i) << " for layer No." << i);
55  std::unique_ptr<TFile> mapsFile(TFile::Open((moduleData->getRadSimFluenceMapList().at(i)).c_str(), "READ")); //this is the ramo potential
56 
57  if (!mapsFile) {
58  ATH_MSG_FATAL("Cannot open file: " << moduleData->getRadSimFluenceMapList().at(i));
59  return StatusCode::FAILURE;
60  }
61 
62  //Setup ramo weighting field map
63  std::unique_ptr<TH3F> ramoPotentialMap_hold(mapsFile->Get<TH3F>("hramomap1"));
64  if (!ramoPotentialMap_hold) {
65  ramoPotentialMap_hold.reset(mapsFile->Get<TH3F>("ramo3d"));
66  }
67  if (!ramoPotentialMap_hold) {
68  ATH_MSG_FATAL("Did not find a Ramo potential map and an approximate form is available yet. Exit...");
69  return StatusCode::FAILURE;
70  }
71 
72  ramoPotentialMap_hold->SetDirectory(nullptr);
73  std::unique_ptr<TH2F> lorentzMap_e_hold(mapsFile->Get<TH2F>("lorentz_map_e"));
74  std::unique_ptr<TH2F> lorentzMap_h_hold(mapsFile->Get<TH2F>("lorentz_map_h"));
75  std::unique_ptr<TH2F> distanceMap_e_hold(mapsFile->Get<TH2F>("edistance"));
76  std::unique_ptr<TH2F> distanceMap_h_hold(mapsFile->Get<TH2F>("hdistance"));
77 
78  if (!lorentzMap_e_hold || !lorentzMap_h_hold || !distanceMap_e_hold || !distanceMap_h_hold) {
79  ATH_MSG_FATAL("Cannot read one of the histograms needed");
80  return StatusCode::FAILURE;
81  }
82 
83  lorentzMap_e_hold->SetDirectory(nullptr);
84  lorentzMap_h_hold->SetDirectory(nullptr);
85  distanceMap_e_hold->SetDirectory(nullptr);
86  distanceMap_h_hold->SetDirectory(nullptr);
87 
88  ramoPotentialMap.emplace_back();
89  ATH_CHECK(ramoPotentialMap.back().setHisto3D(ramoPotentialMap_hold.get()));
90  lorentzMap_e.emplace_back();
91  lorentzMap_h.emplace_back();
92  distanceMap_e.emplace_back();
93  distanceMap_h.emplace_back();
94  ATH_CHECK(lorentzMap_e.back().setHisto2D(lorentzMap_e_hold.get()));
95  ATH_CHECK(lorentzMap_h.back().setHisto2D(lorentzMap_h_hold.get()));
96  ATH_CHECK(distanceMap_e.back().setHisto2D(distanceMap_e_hold.get()));
97  ATH_CHECK(distanceMap_h.back().setHisto2D(distanceMap_h_hold.get()));
98 
99  mapsFile->Close();
100  }
101  writeFluenceCdo -> setLorentzMap_e(std::move(lorentzMap_e));
102  writeFluenceCdo -> setLorentzMap_h(std::move(lorentzMap_h));
103  writeFluenceCdo -> setDistanceMap_e(std::move(distanceMap_e));
104  writeFluenceCdo -> setDistanceMap_h(std::move(distanceMap_h));
105  writeFluenceCdo -> setRamoPotentialMap(std::move(ramoPotentialMap));
106 
107  // Create mapping file for radiation damage simulation for 3D sensor
108  std::vector<PixelHistoConverter> ramoPotentialMap3D;
109  std::vector<PixelHistoConverter> eFieldMap3D;
110  std::vector<PixelHistoConverter> xPositionMap3D_e;
111  std::vector<PixelHistoConverter> xPositionMap3D_h;
112  std::vector<PixelHistoConverter> yPositionMap3D_e;
113  std::vector<PixelHistoConverter> yPositionMap3D_h;
114  std::vector<PixelHistoConverter> timeMap3D_e;
115  std::vector<PixelHistoConverter> timeMap3D_h;
116  PixelHistoConverter avgChargeMap3D_e;
117  PixelHistoConverter avgChargeMap3D_h;
118 
119  for (unsigned int i=0; i<moduleData->getRadSimFluenceMapList3D().size(); i++) {
120  ATH_MSG_INFO("Using maps located in: "<< moduleData->getRadSimFluenceMapList3D().at(i) << " for 3D sensor layer No." << i);
121  std::unique_ptr<TFile> mapsFile3D(TFile::Open((moduleData->getRadSimFluenceMapList3D().at(i)).c_str(), "READ")); //this is the ramo potential
122 
123  if (!mapsFile3D) {
124  ATH_MSG_FATAL("Cannot open file: " << moduleData->getRadSimFluenceMapList3D().at(i));
125  return StatusCode::FAILURE;
126  }
127 
128  //Setup ramo weighting field map
129  std::unique_ptr<TH2F> ramoPotentialMap3D_hold(mapsFile3D->Get<TH2F>("ramo"));
130  std::unique_ptr<TH2F> eFieldMap3D_hold(mapsFile3D->Get<TH2F>("efield"));
131  if (!ramoPotentialMap3D_hold || !eFieldMap3D_hold) {
132  ATH_MSG_FATAL("Did not find a Ramo potential or e-field map for 3D and an approximate form is available yet. Exit...");
133  return StatusCode::FAILURE;
134  }
135  ramoPotentialMap3D_hold->SetDirectory(nullptr);
136  eFieldMap3D_hold->SetDirectory(nullptr);
137  ramoPotentialMap3D.emplace_back();
138  eFieldMap3D.emplace_back();
139  ATH_CHECK(ramoPotentialMap3D.back().setHisto2D(ramoPotentialMap3D_hold.get()));
140  ATH_CHECK(eFieldMap3D.back().setHisto2D(eFieldMap3D_hold.get()));
141 
142  //Now setup the E-field.
143  std::unique_ptr<TH3F> xPositionMap3D_e_hold(mapsFile3D->Get<TH3F>("xPosition_e"));
144  std::unique_ptr<TH3F> xPositionMap3D_h_hold(mapsFile3D->Get<TH3F>("xPosition_h"));
145  std::unique_ptr<TH3F> yPositionMap3D_e_hold(mapsFile3D->Get<TH3F>("yPosition_e"));
146  std::unique_ptr<TH3F> yPositionMap3D_h_hold(mapsFile3D->Get<TH3F>("yPosition_h"));
147  std::unique_ptr<TH2F> timeMap3D_e_hold(mapsFile3D->Get<TH2F>("etimes"));
148  std::unique_ptr<TH2F> timeMap3D_h_hold(mapsFile3D->Get<TH2F>("htimes"));
149 
150  if (!xPositionMap3D_e_hold || !xPositionMap3D_h_hold || !yPositionMap3D_e_hold || !yPositionMap3D_h_hold || !timeMap3D_e_hold || !timeMap3D_h_hold) {
151  ATH_MSG_FATAL("Cannot find one of the maps.");
152  return StatusCode::FAILURE;
153  }
154 
155  xPositionMap3D_e_hold->SetDirectory(nullptr);
156  xPositionMap3D_h_hold->SetDirectory(nullptr);
157  yPositionMap3D_e_hold->SetDirectory(nullptr);
158  yPositionMap3D_h_hold->SetDirectory(nullptr);
159  timeMap3D_e_hold->SetDirectory(nullptr);
160  timeMap3D_h_hold->SetDirectory(nullptr);
161 
162  //Now, determine the time to reach the electrode and the trapping position.
163  xPositionMap3D_e.emplace_back();
164  xPositionMap3D_h.emplace_back();
165  yPositionMap3D_e.emplace_back();
166  yPositionMap3D_h.emplace_back();
167  timeMap3D_e.emplace_back();
168  timeMap3D_h.emplace_back();
169  ATH_CHECK(xPositionMap3D_e.back().setHisto3D(xPositionMap3D_e_hold.get()));
170  ATH_CHECK(xPositionMap3D_h.back().setHisto3D(xPositionMap3D_h_hold.get()));
171  ATH_CHECK(yPositionMap3D_e.back().setHisto3D(yPositionMap3D_e_hold.get()));
172  ATH_CHECK(yPositionMap3D_h.back().setHisto3D(yPositionMap3D_h_hold.get()));
173  ATH_CHECK(timeMap3D_e.back().setHisto2D(timeMap3D_e_hold.get()));
174  ATH_CHECK(timeMap3D_h.back().setHisto2D(timeMap3D_h_hold.get()));
175 
176  std::unique_ptr<TH2F> avgCharge3D_e_hold(mapsFile3D->Get<TH2F>("avgCharge_e"));
177  std::unique_ptr<TH2F> avgCharge3D_h_hold(mapsFile3D->Get<TH2F>("avgCharge_h"));
178 
179  if (!avgCharge3D_e_hold || !avgCharge3D_h_hold) {
180  ATH_MSG_ERROR("Cannot find one of the charge maps.");
181  return StatusCode::FAILURE;
182  }
183 
184  avgCharge3D_e_hold->SetDirectory(nullptr);
185  avgCharge3D_h_hold->SetDirectory(nullptr);
186 
187  // Get average charge data (for charge chunk effect correction)
188  ATH_CHECK(avgChargeMap3D_e.setHisto2D(avgCharge3D_e_hold.get()));
189  ATH_CHECK(avgChargeMap3D_h.setHisto2D(avgCharge3D_h_hold.get()));
190 
191  mapsFile3D->Close();
192  }
193  writeFluenceCdo -> setRamoPotentialMap3D(std::move(ramoPotentialMap3D));
194  writeFluenceCdo -> setEFieldMap3D(std::move(eFieldMap3D));
195  writeFluenceCdo -> setXPositionMap3D_e(std::move(xPositionMap3D_e));
196  writeFluenceCdo -> setXPositionMap3D_h(std::move(xPositionMap3D_h));
197  writeFluenceCdo -> setYPositionMap3D_e(std::move(yPositionMap3D_e));
198  writeFluenceCdo -> setYPositionMap3D_h(std::move(yPositionMap3D_h));
199  writeFluenceCdo -> setTimeMap3D_e(std::move(timeMap3D_e));
200  writeFluenceCdo -> setTimeMap3D_h(std::move(timeMap3D_h));
201  writeFluenceCdo -> setAvgChargeMap3D_e(avgChargeMap3D_e);
202  writeFluenceCdo -> setAvgChargeMap3D_h(avgChargeMap3D_h);
203 
204  if (rangeW.stop().isValid() && rangeW.start()>rangeW.stop()) {
205  ATH_MSG_FATAL("Invalid intersection rangeW: " << rangeW);
206  return StatusCode::FAILURE;
207  }
208  if (writeFluenceMapHandle.record(rangeW, std::move(writeFluenceCdo)).isFailure()) {
209  ATH_MSG_FATAL("Could not record PixelRadiationDamageFluenceMapData " << writeFluenceMapHandle.key() << " with EventRange " << rangeW << " into Conditions Store");
210  return StatusCode::FAILURE;
211  }
212  ATH_MSG_INFO("recorded new CDO " << writeFluenceMapHandle.key() << " with range " << rangeW << " into Conditions Store");
213 
214  return StatusCode::SUCCESS;
215 }
216 
217 
PixelRadSimFluenceMapAlg::m_writeFluenceMapKey
SG::WriteCondHandleKey< PixelRadiationDamageFluenceMapData > m_writeFluenceMapKey
Definition: PixelRadSimFluenceMapAlg.h:33
ATH_MSG_FATAL
#define ATH_MSG_FATAL(x)
Definition: AthMsgStreamMacros.h:34
SG::ReadCondHandle
Definition: ReadCondHandle.h:44
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition: AthMsgStreamMacros.h:31
PixelModuleData::getRadSimFluenceMapList
const std::vector< std::string > & getRadSimFluenceMapList() const
Definition: PixelModuleData.cxx:343
PixelRadSimFluenceMapAlg.h
PixelModuleData
Definition: PixelModuleData.h:22
mergePhysValFiles.start
start
Definition: DataQuality/DataQualityUtils/scripts/mergePhysValFiles.py:14
PixelModuleData::getRadSimFluenceMapList3D
const std::vector< std::string > & getRadSimFluenceMapList3D() const
Definition: PixelModuleData.cxx:349
PixelModuleFeMask_create_db.stop
int stop
Definition: PixelModuleFeMask_create_db.py:76
SG::WriteCondHandle::record
StatusCode record(const EventIDRange &range, T *t)
record handle, with explicit range DEPRECATED
Definition: WriteCondHandle.h:157
python.TrigEgammaMonitorHelper.TH2F
def TH2F(name, title, nxbins, bins_par2, bins_par3, bins_par4, bins_par5=None, bins_par6=None, path='', **kwargs)
Definition: TrigEgammaMonitorHelper.py:45
AthReentrantAlgorithm
An algorithm that can be simultaneously executed in multiple threads.
Definition: AthReentrantAlgorithm.h:83
python.setupRTTAlg.size
int size
Definition: setupRTTAlg.py:39
SG::WriteCondHandle::fullKey
const DataObjID & fullKey() const
Definition: WriteCondHandle.h:41
PixelRadSimFluenceMapAlg::m_moduleDataKey
SG::ReadCondHandleKey< PixelModuleData > m_moduleDataKey
Definition: PixelRadSimFluenceMapAlg.h:30
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition: AthMsgStreamMacros.h:33
PixelHistoConverter
Definition: PixelHistoConverter.h:25
lumiFormat.i
int i
Definition: lumiFormat.py:85
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
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
PixelModuleData::getFluenceLayer3D
const std::vector< double > & getFluenceLayer3D() const
Definition: PixelModuleData.cxx:346
PixelRadSimFluenceMapAlg::initialize
virtual StatusCode initialize() override
Definition: PixelRadSimFluenceMapAlg.cxx:14
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:228
SG::WriteCondHandle::key
const std::string & key() const
Definition: WriteCondHandle.h:40
SG::CondHandleKey::initialize
StatusCode initialize(bool used=true)
PixelHistoConverter::setHisto2D
StatusCode setHisto2D(const TH2 *histo)
Definition: PixelHistoConverter.cxx:41
PixelRadSimFluenceMapAlg::execute
virtual StatusCode execute(const EventContext &ctx) const override
Definition: PixelRadSimFluenceMapAlg.cxx:21
SG::WriteCondHandle::isValid
bool isValid() const
Definition: WriteCondHandle.h:248
PixelRadSimFluenceMapAlg::PixelRadSimFluenceMapAlg
PixelRadSimFluenceMapAlg(const std::string &name, ISvcLocator *pSvcLocator)
Definition: PixelRadSimFluenceMapAlg.cxx:9
SG::WriteCondHandle
Definition: WriteCondHandle.h:26
PixelModuleData::getFluenceLayer
const std::vector< double > & getFluenceLayer() const
Definition: PixelModuleData.cxx:340