ATLAS Offline Software
Loading...
Searching...
No Matches
GetLCSinglePionsPerf.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//-----------------------------------------------------------------------
6// File and Version Information:
7// $Id: CheckSinglePionsReco.cxx,v 1.1.1.1 2008/11/04 08:56:11 menke Exp $
8//
9// Description: see CheckSinglePionsReco.h
10//
11// Environment:
12// Software developed for the ATLAS Detector at CERN LHC
13//
14// Author List:
15// Gennady Pospelov
16//
17//-----------------------------------------------------------------------
18
19//-----------------------
20// This Class's Header --
21//-----------------------
22#include "CaloLocalHadCalib/GetLCSinglePionsPerf.h"
23
24#include <iterator>
25#include <cmath>
26
27#include "AthenaKernel/errorcheck.h"
28#include "GaudiKernel/ISvcLocator.h"
29#include "GaudiKernel/StatusCode.h"
30
31#include "CaloEvent/CaloCell.h"
32#include "CaloSimEvent/CaloCalibrationHit.h"
33#include "CaloSimEvent/CaloCalibrationHitContainer.h"
34#include "CaloDmDetDescr/CaloDmDescrManager.h"
35#include "CaloIdentifier/CaloCell_ID.h"
36#include "GeneratorObjects/McEventCollection.h"
37#include <CLHEP/Vector/LorentzVector.h>
38#include "CaloEvent/CaloCompositeKineBase.h"
39#include "CaloEvent/CaloRecoStatus.h"
40
41#include "xAODCaloEvent/CaloCluster.h"
42#include "TruthUtils/HepMCHelpers.h"
43
44#include "TProfile.h"
45#include "TProfile2D.h"
46#include "TFile.h"
47#include "TH1F.h"
48#include "TVector3.h"
49#include "TTree.h"
50
51
52using CLHEP::HepLorentzVector;
53using CLHEP::MeV;
54
55
56/* ****************************************************************************
57
58**************************************************************************** */
59GetLCSinglePionsPerf::GetLCSinglePionsPerf(const std::string& name, ISvcLocator* pSvcLocator)
60 : AthAlgorithm(name, pSvcLocator),
61 m_id_helper(nullptr),
62 m_calo_id(nullptr),
63 m_caloDmDescrManager(nullptr),
64 m_outputFileName("CheckSinglePionsReco.root"),
65 m_outputFile(nullptr),
66 m_distance_cut(1.5),
67 m_netabin(25), m_etamin(0.0), m_etamax(5.0),
68 m_deta(0),
69 m_nphibin(1), m_phimin(-M_PI), m_phimax(M_PI),
70 m_dphi(0),
71 m_nlogenerbin(22), m_logenermin(2.0), m_logenermax(6.4),
72 m_nlogenerbin2(3),
73 m_dlogener(0),
74 m_nnormtype(3),
75 m_ncluscoll(4),
76 m_ntagcases(3),
77 m_nmoments(0),
78 m_nmomsums(6),
79 m_truthPiEngFraction(0.1),
80 m_doEngRecOverTruth(true),
81 m_doEngTag(true),
82 m_doEngRecSpect(true),
83 m_doEngNoiseClus(true),
84 m_doClusMoments(true),
85 m_doRecoEfficiency(true),
86 m_useRecoEfficiency(false),
87 m_useGoodClus(false),
88 m_doCalibHitsValidation(true),
89 m_mc_ener(0),
90 m_mc_eta(0),
91 m_mc_phi(0),
92 m_mc_etabin(0),
93 m_mc_enerbin(0),
94 m_mc_phibin(0)
95{
96 // collection name where cluster always have moments
97 declareProperty("ClusterBasicCollName",m_clusterBasicCollName);
98
99 // collections names after different correction steps
100 declareProperty("ClusterCollectionNames",m_clusterCollNames);
101
102 // distance cut for clusters to calcuate engReco
103 declareProperty("DistanceCut",m_distance_cut);
104
105 // Name of output file to save histograms in
106 declareProperty("OutputFileName",m_outputFileName);
107
108 // Names of CalibrationHitContainers to use
109 declareProperty("CalibrationHitContainerNames",m_CalibrationHitContainerNames);
110 declareProperty("DMCalibrationHitContainerNames",m_DMCalibrationHitContainerNames);
111
112 // to create different types of validation histograms
113 declareProperty("doEngRecOverTruth",m_doEngRecOverTruth);
114 declareProperty("doEngTag",m_doEngTag);
115 declareProperty("doEngRecSpect",m_doEngRecSpect);
116 declareProperty("doEngNoiseClus",m_doEngNoiseClus);
117 declareProperty("doClusMoments",m_doClusMoments);
118 declareProperty("doRecoEfficiency",m_doRecoEfficiency);
119 declareProperty("useRecoEfficiency",m_useRecoEfficiency);
120 declareProperty("useGoodClus",m_useGoodClus);
121
122 declareProperty("etamin",m_etamin);
123 declareProperty("etamax",m_etamax);
124 declareProperty("netabin",m_netabin);
125
126 declareProperty("phimin",m_phimin);
127 declareProperty("phimax",m_phimax);
128 declareProperty("nphibin",m_nphibin);
129
130 declareProperty("logenermin",m_logenermin);
131 declareProperty("logenermax",m_logenermax);
132 declareProperty("nlogenerbin",m_nlogenerbin);
133
134 declareProperty("doCalibHitsValidation",m_doCalibHitsValidation);
135// Only use clusters that contain at least m_truthPiEngFraction % of the true pion energy
136 declareProperty("usePionClustersOnly",m_usePionClusters);
137
138 // dead material identifier description manager
139 m_caloDmDescrManager = CaloDmDescrManager::instance();
140
141 m_validMoments.resize(0);
142 m_validMoments.push_back(moment_name_pair(std::string("ENG_CALIB_TOT"),xAOD::CaloCluster::ENG_CALIB_TOT));
143 m_validMoments.push_back(moment_name_pair(std::string("ENG_CALIB_OUT_L"),xAOD::CaloCluster::ENG_CALIB_OUT_L));
144 m_validMoments.push_back(moment_name_pair(std::string("ENG_CALIB_DEAD_TOT"),xAOD::CaloCluster::ENG_CALIB_DEAD_TOT));
145 m_validMoments.push_back(moment_name_pair(std::string("ENG_CALIB_DEAD_EMB0"),xAOD::CaloCluster::ENG_CALIB_DEAD_EMB0));
146 m_validMoments.push_back(moment_name_pair(std::string("ENG_CALIB_DEAD_TILE0"),xAOD::CaloCluster::ENG_CALIB_DEAD_TILE0));
147 m_validMoments.push_back(moment_name_pair(std::string("ENG_CALIB_DEAD_TILEG3"),xAOD::CaloCluster::ENG_CALIB_DEAD_TILEG3));
148 m_validMoments.push_back(moment_name_pair(std::string("ENG_CALIB_DEAD_EME0"),xAOD::CaloCluster::ENG_CALIB_DEAD_EME0));
149 m_validMoments.push_back(moment_name_pair(std::string("ENG_CALIB_DEAD_HEC0"),xAOD::CaloCluster::ENG_CALIB_DEAD_HEC0));
150 m_validMoments.push_back(moment_name_pair(std::string("ENG_CALIB_DEAD_FCAL"),xAOD::CaloCluster::ENG_CALIB_DEAD_FCAL));
151 m_validMoments.push_back(moment_name_pair(std::string("ENG_CALIB_DEAD_LEAKAGE"),xAOD::CaloCluster::ENG_CALIB_DEAD_LEAKAGE));
152 m_validMoments.push_back(moment_name_pair(std::string("ENG_CALIB_DEAD_UNCLASS"),xAOD::CaloCluster::ENG_CALIB_DEAD_UNCLASS));
153
154 m_nmoments = m_validMoments.size();
155
156}
157
158
159
160/* ****************************************************************************
161
162**************************************************************************** */
167
168
169
170/* ****************************************************************************
171
172**************************************************************************** */
173StatusCode GetLCSinglePionsPerf::initialize()
174{
175 ATH_CHECK( detStore()->retrieve(m_calo_id,"CaloCell_ID") );
176 ATH_CHECK( detStore()->retrieve(m_id_helper) );
177
178
179 m_deta = (m_etamax - m_etamin)/m_netabin;
180 m_dphi = (m_phimax - m_phimin)/m_nphibin;
181 m_dlogener = (m_logenermax - m_logenermin)/m_nlogenerbin;
182
183 m_outputFile = new TFile(m_outputFileName.c_str(),"RECREATE");
184 m_outputFile->cd();
185
186 // defining histograms
187 char hname[256];
188 float h_eta_min = -5.0;
189 float h_eta_max=5.0;
190 int h_nch_eta = 50;
191 if(m_etamin != 0.0) {
192 h_eta_min = m_etamin;
193 h_eta_max = m_etamax;
194 }
195 const int n_logener_bins = m_nlogenerbin*2;
196 double *xbins = new double[n_logener_bins+1];
197 const double inv_n_logener_bins = 1. / static_cast<double> (n_logener_bins);
198 for(int i_bin=0; i_bin<=n_logener_bins; i_bin++) {
199 xbins[i_bin] = pow(10,m_logenermin+i_bin*(m_logenermax-m_logenermin)*inv_n_logener_bins);
200 }
201 m_ncluscoll = m_clusterCollNames.size();
202
203 /* **************************************************************************
204 histograms to check classification performance
205 *************************************************************************** */
206 if(m_doEngTag) {
207 // as a function of pion eta for several energy bins
208 m_engTag_vs_eta.resize(m_ntagcases);
209 for(int i_tag=0; i_tag<m_ntagcases; i_tag++) {
210 m_engTag_vs_eta[i_tag].resize(m_nlogenerbin, nullptr);
211 for(int i_ener=0; i_ener<m_nlogenerbin; i_ener++){
212 sprintf(hname, "hp_engTag_vs_eta_tag%d_ener%d",i_tag, i_ener);
213 TProfile *hp = new TProfile(hname, hname, h_nch_eta, h_eta_min, h_eta_max);
214 hp->GetXaxis()->SetTitle("#eta");
215 hp->GetYaxis()->SetTitle("E_{class}/E_{tot}");
216 m_engTag_vs_eta[i_tag][i_ener] = hp;
217 } // i_ener
218 } // i_tag
219 // as a function of pion energy for several eta bins
220 m_engTag_vs_ebeam.resize(m_ntagcases);
221 for(int i_tag=0; i_tag<m_ntagcases; i_tag++){
222 m_engTag_vs_ebeam[i_tag].resize(m_netabin, nullptr);
223 for(int i_eta=0; i_eta<m_netabin; i_eta++) {
224 sprintf(hname, "hp_engTag_vs_ebeam_tag%d_eta%d",i_tag, i_eta);
225 TProfile *hp = new TProfile(hname, hname, n_logener_bins, xbins);
226 hp->GetXaxis()->SetTitle("E_{#pi}, MeV");
227 hp->GetYaxis()->SetTitle("E_{class}/E_{tot}");
228 m_engTag_vs_ebeam[i_tag][i_eta] = hp;
229 }
230 }
231 } // m_doEngTag
232
233
234 /* **************************************************************************
235 histograms to check reconstructed energy after each step of local hadronic
236 calibration, reco energy is defined as a sum of clusters energies within certain
237 angle around particle direction
238 *************************************************************************** */
239 if(m_doEngRecOverTruth) {
240 // as a function of pion eta for several energy bins
241 m_engRecOverTruth_vs_eta.resize(m_ncluscoll);
242 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++){
243 m_engRecOverTruth_vs_eta[i_coll].resize(m_nnormtype);
244 for(int i_norm=0; i_norm<m_nnormtype; i_norm++){
245 m_engRecOverTruth_vs_eta[i_coll][i_norm].resize(m_nlogenerbin, nullptr);
246 for(int i_ener=0; i_ener<m_nlogenerbin; i_ener++){
247 double ylow(0), yup(0);
248 if(i_norm != kLEVEL_PARTICLE) {
249 ylow = 0.0;
250 yup = 100.0;
251 }
252 sprintf(hname, "hp_engRecOverTruth_vs_eta_coll%d_norm%d_ener%d",i_coll, i_norm, i_ener);
253 TProfile *hp = new TProfile(hname, hname, h_nch_eta, h_eta_min, h_eta_max, ylow, yup);
254 hp->GetXaxis()->SetTitle("#eta");
255 if(i_norm==kLEVEL_PARTICLE) {
256 hp->GetYaxis()->SetTitle("E_{reco}/E_{#pi}");
257 }else{
258 hp->GetYaxis()->SetTitle("E_{reco}/E_{calib}");
259 }
260 m_engRecOverTruth_vs_eta[i_coll][i_norm][i_ener] = hp;
261 } // i_ener
262 } // i_norm
263 } // i_coll
264 // as a function of pion energy for several eta bins
265 m_engRecOverTruth_vs_ebeam.resize(m_ncluscoll);
266 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++){
267 m_engRecOverTruth_vs_ebeam[i_coll].resize(m_nnormtype);
268 for(int i_norm=0; i_norm<m_nnormtype; i_norm++){
269 m_engRecOverTruth_vs_ebeam[i_coll][i_norm].resize(m_netabin, nullptr);
270 for(int i_eta=0; i_eta<m_netabin; i_eta++) {
271 double ylow(0), yup(0);
272 if(i_norm != kLEVEL_PARTICLE) {
273 ylow = 0.0;
274 yup = 100.0;
275 }
276 sprintf(hname, "hp_m_engRecOverTruth_vs_ebeam_coll%d_norm%d_eta%d",i_coll, i_norm, i_eta);
277 TProfile *hp = new TProfile(hname, hname, n_logener_bins, xbins, ylow, yup);
278 hp->GetXaxis()->SetTitle("E_{#pi}, MeV");
279 if(i_norm==kLEVEL_PARTICLE) {
280 hp->GetYaxis()->SetTitle("E_{reco}/E_{#pi}");
281 }else{
282 hp->GetYaxis()->SetTitle("E_{reco}/E_{calib}");
283 }
284 m_engRecOverTruth_vs_ebeam[i_coll][i_norm][i_eta] = hp;
285 } // i_eta
286 } // i_norm
287 } // i_coll
288 } // m_doEngRecOverTruth
289
290
291 /* **************************************************************************
292 histograms to save reco and noise clusters energy spectras
293 *************************************************************************** */
294 if(m_doEngRecSpect) {
295 m_engRecSpect.resize(m_ncluscoll);
296 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++){
297 m_engRecSpect[i_coll].resize(m_nlogenerbin);
298 for(int i_ener=0; i_ener<m_nlogenerbin; i_ener++){
299 m_engRecSpect[i_coll][i_ener].resize(m_netabin, nullptr);
300 for(int i_eta=0; i_eta<m_netabin; i_eta++){
301 sprintf(hname, "hp_engRecSpect_coll%d_ener%d_eta%d",i_coll, i_ener, i_eta);
302 TH1F *h1 = new TH1F(hname, hname, 110, -0.2, 2.0);
303 h1->GetXaxis()->SetTitle("E_{reco}/E_{#pi}");
304 m_engRecSpect[i_coll][i_ener][i_eta] = h1;
305 } // i_eta
306 } // i_ener
307 } // i_coll
308 // energy spectra of incident pions
309 m_engPionSpect.resize(m_nlogenerbin);
310 for(int i_ener=0; i_ener<m_nlogenerbin; i_ener++){
311 m_engPionSpect[i_ener].resize(m_netabin, nullptr);
312 for(int i_eta=0; i_eta<m_netabin; i_eta++){
313 sprintf(hname, "hp_engPionSpect_ener%d_eta%d",i_ener, i_eta);
314 double x1 = pow(10,m_logenermin+i_ener*m_dlogener);
315 double x2 = pow(10,m_logenermin+(i_ener+1)*m_dlogener);
316 TH1F *h1 = new TH1F(hname, hname, 20, x1, x2);
317 h1->GetXaxis()->SetTitle("E_{#pi}, MeV");
318 m_engPionSpect[i_ener][i_eta] = h1;
319 } // i_eta
320 } // i_ener
321 } // m_doEngRecSpect
322
323 /* **************************************************************************
324 histograms to save noise clusters spectras
325 *************************************************************************** */
326 if(m_doEngNoiseClus) {
327 // noise clusters spectras
328 m_engNoiseClusSpect.resize(m_ncluscoll);
329 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++){
330 m_engNoiseClusSpect[i_coll].resize(m_netabin, nullptr);
331 for(int i_eta=0; i_eta<m_netabin; i_eta++){
332 sprintf(hname, "hp_engNoiseClusSpect_coll%d_eta%d",i_coll, i_eta);
333 TH1F *h1 = new TH1F(hname, hname, 100, -2000., 2000.);
334 h1->GetXaxis()->SetTitle("E_{clus}, MeV");
335 m_engNoiseClusSpect[i_coll][i_eta] = h1;
336 } // i_eta
337 } // i_coll
338 // average
339 m_engNoiseClus_vs_eta.resize(m_ncluscoll, nullptr);
340 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++){
341 sprintf(hname, "hp_engNoiseClus_vs_eta_coll%d",i_coll);
342 TProfile *hp = new TProfile(hname, hname, h_nch_eta, h_eta_min, h_eta_max);
343 hp->GetXaxis()->SetTitle("#eta");
344 hp->GetYaxis()->SetTitle("E_{clus}, MeV");
345 m_engNoiseClus_vs_eta[i_coll] = hp;
346 //hp->SetErrorOption("s");
347 }
348 } // m_doEngNoiseClus
349
350 /* **************************************************************************
351 histograms to save cluster moments
352 *************************************************************************** */
353 if(m_doClusMoments) {
354 m_clusMoment_vs_eta.resize(m_nmoments);
355 for(int i_mom=0; i_mom<m_nmoments; i_mom++){
356 m_clusMoment_vs_eta[i_mom].resize(m_nmomsums);
357 for(int i_sum=0; i_sum<m_nmomsums; i_sum++){
358 m_clusMoment_vs_eta[i_mom][i_sum].resize(m_nlogenerbin, nullptr);
359 for(int i_ener=0; i_ener<m_nlogenerbin; i_ener++){
360 sprintf(hname, "hp_clusMoment_vs_eta_mom%d_sum%d_ener%d",i_mom, i_sum, i_ener);
361 TProfile *hp = new TProfile(hname, hname, h_nch_eta, h_eta_min, h_eta_max);
362 hp->GetXaxis()->SetTitle("#eta");
363 hp->GetYaxis()->SetTitle("Mom_{cls}/Mom_{sum}");
364 hp->SetTitle(m_validMoments[i_mom].first.c_str());
365 m_clusMoment_vs_eta[i_mom][i_sum][i_ener] = hp;
366 } // i_ener
367 } // i_sum
368 } // i_mom
369
370 m_clusMoment_vs_ebeam.resize(m_nmoments);
371 for(int i_mom=0; i_mom<m_nmoments; i_mom++){
372 m_clusMoment_vs_ebeam[i_mom].resize(m_nmomsums);
373 for(int i_sum=0; i_sum<m_nmomsums; i_sum++){
374 m_clusMoment_vs_ebeam[i_mom][i_sum].resize(m_netabin, nullptr);
375 for(int i_eta=0; i_eta<m_netabin; i_eta++){
376 sprintf(hname, "hp_clusMoment_vs_ebeam_mom%d_sum%d_eta%d",i_mom, i_sum, i_eta);
377 TProfile *hp = new TProfile(hname, hname, n_logener_bins, xbins);
378 hp->GetXaxis()->SetTitle("E_{#pi}");
379 hp->GetYaxis()->SetTitle("Mom_{cls}/Mom_{sum}");
380 hp->SetTitle(m_validMoments[i_mom].first.c_str());
381 m_clusMoment_vs_ebeam[i_mom][i_sum][i_eta] = hp;
382 } // i_eta
383 } // i_sum
384 } // i_mom
385
386 } // m_doClusMoments
387
388
389 /* **************************************************************************
390 histograms to save reconstruction efficiency
391 *************************************************************************** */
392 if(m_doRecoEfficiency) {
393 // as a function of eta for several energy bins
394 for(int i_eff=0; i_eff<2; i_eff++){ // 0-reconstructed, 1-all events
395 m_RecoEfficiency_vs_eta[i_eff].resize(m_nlogenerbin, nullptr);
396 for(int i_ener=0; i_ener<m_nlogenerbin; i_ener++){
397 sprintf(hname, "hp_m_RecoEfficiency_eff%d_ener%d", i_eff, i_ener);
398 TH1F *h1 = new TH1F(hname, hname, h_nch_eta, h_eta_min, h_eta_max);
399 h1->GetXaxis()->SetTitle("#eta");
400 h1->GetYaxis()->SetTitle("Efficiency");
401 m_RecoEfficiency_vs_eta[i_eff][i_ener] = h1;
402 } // i_eta
403 } // i_eff
404 // as a function of energy for several eta bins
405 for(int i_eff=0; i_eff<2; i_eff++){ // 0-reconstructed, 1-all events
406 m_RecoEfficiency_vs_ebeam[i_eff].resize(m_netabin, nullptr);
407 for(int i_eta=0; i_eta<m_netabin; i_eta++){
408 sprintf(hname, "hp_m_RecoEfficiency_eff%d_eta%d", i_eff, i_eta);
409 TH1F *h1 = new TH1F(hname, hname, n_logener_bins, xbins);
410 h1->GetXaxis()->SetTitle("E_{#pi}, MeV");
411 h1->GetYaxis()->SetTitle("Efficiency");
412 m_RecoEfficiency_vs_ebeam[i_eff][i_eta] = h1;
413 } // i_eta
414 } // i_eff
415 }
416
417 /* **************************************************************************
418 histograms to check calibration hits machinery
419 *************************************************************************** */
420 if(m_doCalibHitsValidation) {
421 m_SumCalibHitOverEbeam_vs_eta.resize(m_nlogenerbin2, nullptr);
422 m_DefaultCalculatorOverEbeam_vs_eta.resize(m_nlogenerbin2, nullptr);
423 m_SumCalibHitAssignedOverEbeam_vs_eta.resize(m_nlogenerbin2, nullptr);
424 m_SumCalibHitAssignedOverEbeam_vs_etaphi.resize(m_nlogenerbin2, nullptr);
425 for(int i_ener=0; i_ener<m_nlogenerbin2; i_ener++){
426 sprintf(hname, "hp_SumCalibHitOverEbeam_vs_eta_ener%d", i_ener);
427 TProfile *hp = new TProfile(hname, hname, h_nch_eta, h_eta_min, h_eta_max);
428 hp->GetXaxis()->SetTitle("#eta");
429 hp->GetYaxis()->SetTitle("#sum E_{calib}/E_{#pi}");
430 m_SumCalibHitOverEbeam_vs_eta[i_ener] = hp;
431 sprintf(hname, "hp_DefaultCalculatorOverEbeam_vs_eta_ener%d", i_ener);
432 hp = new TProfile(hname, hname, h_nch_eta, h_eta_min, h_eta_max);
433 hp->GetXaxis()->SetTitle("#eta");
434 hp->GetYaxis()->SetTitle("#sum E_{DefCalc}/E_{#pi}");
435 m_DefaultCalculatorOverEbeam_vs_eta[i_ener] = hp;
436 sprintf(hname, "hp_SumCalibHitAssignedOverEbeam_vs_eta_ener%d", i_ener);
437 hp = new TProfile(hname, hname, h_nch_eta, h_eta_min, h_eta_max);
438 hp->GetXaxis()->SetTitle("#eta");
439 hp->GetYaxis()->SetTitle("#sum E_{assigned}/E_{#pi}");
440 m_SumCalibHitAssignedOverEbeam_vs_eta[i_ener] = hp;
441 sprintf(hname, "hp_SumCalibHitAssignedOverEbeam_vs_etaphi_ener%d", i_ener);
442 TProfile2D *hp2 = new TProfile2D(hname, hname, 25, h_eta_min, h_eta_max, 25, -M_PI, M_PI);
443 hp2->GetXaxis()->SetTitle("#eta");
444 hp2->GetYaxis()->SetTitle("#phi");
445 hp2->GetZaxis()->SetTitle("#sum E_{assigned}/E_{#pi}");
446 m_SumCalibHitAssignedOverEbeam_vs_etaphi[i_ener] = hp2;
447 } // i_ener
448
449// m_SumCalibHitOverEbeam_vs_ebeam.resize(m_netabin, 0);
450// // m_DefaultCalculatorOverEbeam_vs_ebeam.resize(m_netabin, 0);
451// for(int i_eta=0; i_eta<m_netabin; i_eta++){
452// sprintf(hname, "hp_SumCalibHitOverEbeam_vs_ebeam_eta%d", i_eta);
453// TProfile *hp = new TProfile(hname, hname, n_logener_bins, xbins);
454// hp->GetXaxis()->SetTitle("#eta");
455// hp->GetYaxis()->SetTitle("#sum E_{calib}/E_{#pi}");
456// m_SumCalibHitOverEbeam_vs_ebeam[i_eta] = hp;
457// // sprintf(hname, "hp_DefaultCalculatorOverEbeam_vs_ebeam_eta%d", i_eta);
458// // hp = new TProfile(hname, hname, n_logener_bins, xbins);
459// // hp->GetXaxis()->SetTitle("#eta");
460// // hp->GetYaxis()->SetTitle("#sum E_{DefCalc}/E_{#pi}");
461// // m_DefaultCalculatorOverEbeam_vs_ebeam[i_eta] = hp;
462// } // i_eta
463 } // m_doCalibHitsValidation
464
465
466 delete [] xbins;
467
468 ATH_CHECK( m_clusterBasicCollName.initialize() );
469 ATH_CHECK( m_clusterCollNames.initialize() );
470 ATH_CHECK( m_CalibrationHitContainerNames.initialize() );
471 ATH_CHECK( m_DMCalibrationHitContainerNames.initialize() );
472
473 return StatusCode::SUCCESS;
474}
475
476
477
478/* ****************************************************************************
479writing histograms into the file
480**************************************************************************** */
481StatusCode GetLCSinglePionsPerf::finalize()
482{
483 ATH_MSG_INFO( "Writing out histograms" );
484
485 // ereco over truth after each step of calibration
486 if(m_doEngRecOverTruth) {
487 m_outputFile->cd();
488 m_outputFile->mkdir("engRecOverTruth");
489 m_outputFile->cd("engRecOverTruth");
490 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++){
491 for(int i_norm=0; i_norm<m_nnormtype; i_norm++){
492 for(int i_ener=0; i_ener<m_nlogenerbin; i_ener++){
493 if(m_engRecOverTruth_vs_eta[i_coll][i_norm][i_ener]) m_engRecOverTruth_vs_eta[i_coll][i_norm][i_ener]->Write();
494 } // i_ener
495 } // i_norm
496 } // i_coll
497 //
498 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++){
499 for(int i_norm=0; i_norm<m_nnormtype; i_norm++){
500 for(int i_eta=0; i_eta<m_netabin; i_eta++) {
501 if(m_engRecOverTruth_vs_ebeam[i_coll][i_norm][i_eta]) m_engRecOverTruth_vs_ebeam[i_coll][i_norm][i_eta]->Write();
502 } // i_eta
503 } // i_norm
504 } // i_coll
505 } // m_doEngRecOverTruth
506
507 // classification performance
508 if(m_doEngTag) {
509 m_outputFile->cd();
510 m_outputFile->mkdir("engTag");
511 m_outputFile->cd("engTag");
512 for(int i_tag=0; i_tag<m_ntagcases; i_tag++) {
513 for(int i_ener=0; i_ener<m_nlogenerbin; i_ener++){
514 if(m_engTag_vs_eta[i_tag][i_ener]) m_engTag_vs_eta[i_tag][i_ener]->Write();
515 } // i_ener
516 } // i_tag
517 // as a function of pion energy for several eta bins
518 for(int i_tag=0; i_tag<m_ntagcases; i_tag++){
519 for(int i_eta=0; i_eta<m_netabin; i_eta++) {
520 if(m_engTag_vs_ebeam[i_tag][i_eta]) m_engTag_vs_ebeam[i_tag][i_eta]->Write();
521 }
522 }
523 } // m_doEngTag
524
525 // reco energy spactras
526 if(m_doEngRecSpect) {
527 m_outputFile->cd();
528 m_outputFile->mkdir("engRecSpect");
529 m_outputFile->cd("engRecSpect");
530 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++){
531 for(int i_ener=0; i_ener<m_nlogenerbin; i_ener++){
532 for(int i_eta=0; i_eta<m_netabin; i_eta++){
533 if(m_engRecSpect[i_coll][i_ener][i_eta]) m_engRecSpect[i_coll][i_ener][i_eta]->Write();
534 } // i_eta
535 } // i_ener
536 } // i_coll
537 // pione energy spectra
538 for(int i_ener=0; i_ener<m_nlogenerbin; i_ener++){
539 for(int i_eta=0; i_eta<m_netabin; i_eta++){
540 if(m_engPionSpect[i_ener][i_eta]) m_engPionSpect[i_ener][i_eta]->Write();
541 } // i_eta
542 } // i_ener
543 } // m_doEngRecSpect
544
545 // noise clusters spectra
546 if(m_doEngNoiseClus) {
547 m_outputFile->cd();
548 m_outputFile->mkdir("engNoiseClus");
549 m_outputFile->cd("engNoiseClus");
550 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++){
551 for(int i_eta=0; i_eta<m_netabin; i_eta++){
552 if(m_engNoiseClusSpect[i_coll][i_eta]) m_engNoiseClusSpect[i_coll][i_eta]->Write();
553 } // i_eta
554 } // i_coll
555 // average noise
556 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++){
557 if(m_engNoiseClus_vs_eta[i_coll]) m_engNoiseClus_vs_eta[i_coll] -> Write();
558 }
559 } // m_doEngRecSpect
560
561 if(m_doClusMoments) {
562 m_outputFile->cd();
563 m_outputFile->mkdir("clusMoment");
564 m_outputFile->cd("clusMoment");
565 for(int i_mom=0; i_mom<m_nmoments; i_mom++){
566 for(int i_sum=0; i_sum<m_nmomsums; i_sum++){
567 for(int i_ener=0; i_ener<m_nlogenerbin; i_ener++){
568 if(m_clusMoment_vs_eta[i_mom][i_sum][i_ener]) m_clusMoment_vs_eta[i_mom][i_sum][i_ener]->Write();
569 } // i_ener
570 } // i_sum
571 } // i_mom
572 for(int i_mom=0; i_mom<m_nmoments; i_mom++){
573 for(int i_sum=0; i_sum<m_nmomsums; i_sum++){
574 for(int i_eta=0; i_eta<m_netabin; i_eta++){
575 if(m_clusMoment_vs_ebeam[i_mom][i_sum][i_eta]) m_clusMoment_vs_ebeam[i_mom][i_sum][i_eta]->Write();
576 } // i_eta
577 } // i_sum
578 } // i_mom
579
580 } // m_doClusMoments
581
582 if(m_doRecoEfficiency) {
583 m_outputFile->cd();
584 m_outputFile->mkdir("RecoEff");
585 m_outputFile->cd("RecoEff");
586 for(int i_ener=0; i_ener<m_nlogenerbin; i_ener++){
587 if(m_RecoEfficiency_vs_eta[0][i_ener] && m_RecoEfficiency_vs_eta[1][i_ener]) {
588 m_RecoEfficiency_vs_eta[0][i_ener]->Write();
589 m_RecoEfficiency_vs_eta[1][i_ener]->Write();
590 }
591 } // i_eta
592 for(int i_eta=0; i_eta<m_netabin; i_eta++){
593 if(m_RecoEfficiency_vs_ebeam[0][i_eta] && m_RecoEfficiency_vs_ebeam[1][i_eta]) {
594 m_RecoEfficiency_vs_ebeam[0][i_eta]->Write();
595 m_RecoEfficiency_vs_ebeam[1][i_eta]->Write();
596 }
597 } // i_eta
598 } // m_doRecoEfficiency
599
600 /* **************************************************************************
601 histograms to check calibration hits machinery
602 *************************************************************************** */
603 if(m_doCalibHitsValidation) {
604 m_outputFile->cd();
605 m_outputFile->mkdir("CalibHitsValidation");
606 m_outputFile->cd("CalibHitsValidation");
607 for(int i_ener=0; i_ener<m_nlogenerbin2; i_ener++){
608 if(m_SumCalibHitOverEbeam_vs_eta[i_ener]) m_SumCalibHitOverEbeam_vs_eta[i_ener]->Write();
609 if(m_DefaultCalculatorOverEbeam_vs_eta[i_ener]) m_DefaultCalculatorOverEbeam_vs_eta[i_ener]->Write();
610 if(m_SumCalibHitAssignedOverEbeam_vs_eta[i_ener]) m_SumCalibHitAssignedOverEbeam_vs_eta[i_ener]->Write();
611 if(m_SumCalibHitAssignedOverEbeam_vs_etaphi[i_ener]) m_SumCalibHitAssignedOverEbeam_vs_etaphi[i_ener]->Write();
612 } // i_ener
613// for(int i_eta=0; i_eta<m_netabin; i_eta++){
614// if(m_SumCalibHitOverEbeam_vs_ebeam[i_eta]) m_SumCalibHitOverEbeam_vs_ebeam[i_eta]->Write();
615// // if(m_DefaultCalculatorOverEbeam_vs_ebeam[i_eta]) m_DefaultCalculatorOverEbeam_vs_ebeam[i_eta]->Write();
616// } // i_eta
617
618 } // m_doCalibHitsValidation
619
620 /* ********************************************
621 writing tree with histogram parameters
622 ******************************************** */
623 m_outputFile->cd();
624 m_outputFile->mkdir("Parameters");
625 m_outputFile->cd("Parameters");
626 TTree *tree = new TTree("ParamTree","ParamTree");
627 tree->Branch("m_netabin",&m_netabin,"m_netabin/I");
628 tree->Branch("m_etamin",&m_etamin,"m_etamin/F");
629 tree->Branch("m_etamax",&m_etamax,"m_etamax/F");
630 tree->Branch("m_nphibin",&m_nphibin,"m_nphibin/I");
631 tree->Branch("m_phimin",&m_phimin,"m_phimin/F");
632 tree->Branch("m_phimax",&m_phimax,"m_phimax/F");
633 tree->Branch("m_nlogenerbin",&m_nlogenerbin,"m_nlogenerbin/I");
634 tree->Branch("m_logenermin",&m_logenermin,"m_logenermin/F");
635 tree->Branch("m_logenermax",&m_logenermax,"m_logenermax/F");
636 tree->Branch("m_nnormtype",&m_nnormtype,"m_nnormtype/I");
637 tree->Branch("m_ncluscoll",&m_ncluscoll,"m_ncluscoll/I");
638 tree->Branch("m_ntagcases",&m_ntagcases,"m_ntagcases/I");
639 tree->Branch("m_nmoments",&m_nmoments,"m_nmoments/I");
640 tree->Branch("m_nmomsums",&m_nmomsums,"m_nmomsums/I");
641 tree->Fill();
642 tree->Write();
643
644 m_outputFile->Close();
645
646 return StatusCode::SUCCESS;
647}
648
649
650
651/* ****************************************************************************
652
653**************************************************************************** */
654StatusCode GetLCSinglePionsPerf::execute(const EventContext& ctx)
655{
656
657 /* ********************************************
658 reading particles
659 ******************************************** */
660 const McEventCollection* truthEvent=nullptr;
661 ATH_CHECK( evtStore()->retrieve(truthEvent, "TruthEvent") );
662
663 if( truthEvent->at(0)->particles_empty() ){
664 ATH_MSG_ERROR( "No particles in McEventCollection" );
665 return StatusCode::FAILURE;
666 }
667 HepMC::ConstGenParticlePtr gen{nullptr};
668 for (auto p: *truthEvent->at(0)) {
669 if (MC::isGenStable(p)) {
670 if (std::abs(p->pdg_id()) == MC::PIPLUS ||
671 std::abs(p->pdg_id()) == MC::PI0) {
672 gen = std::move(p);
673 break;
674 }
675 }
676 }
677 if(!gen){
678 ATH_MSG_ERROR( "No final stable pion in McEventCollection" );
679 return StatusCode::FAILURE;
680 }
681 m_mc_eta = gen->momentum().pseudoRapidity();
682 m_mc_phi = gen->momentum().phi();
683 m_mc_ener = gen->momentum().e();
684
685 // energy and eta bins for particle
686 m_mc_enerbin = int( (log10(m_mc_ener) - m_logenermin)/m_dlogener);
687 m_mc_etabin = int((fabs(m_mc_eta)-m_etamin)/m_deta);
688 m_mc_phibin = int((m_mc_phi-m_phimin)/m_dphi);
689
690 if(m_mc_etabin <0 || m_mc_etabin >= m_netabin || m_mc_enerbin <0 || m_mc_enerbin >= m_nlogenerbin || m_mc_phibin !=0) return StatusCode::SUCCESS;
691
692 SG::ReadHandle<xAOD::CaloClusterContainer> clusColl (m_clusterBasicCollName, ctx);
693 if(clusColl->empty()) {
694 ATH_MSG_WARNING( "Empty ClusterContainer " << m_clusterBasicCollName );
695 return StatusCode::SUCCESS;
696 }
697
698 // check that we have at least one meaningfull cluster
699 if(m_doRecoEfficiency) {
700 bool PionWasReconstructed = false;
701 float engClusCalibThs = 1.0*MeV;
702 HepLorentzVector hlv_pion(1,0,0,1);
703 hlv_pion.setREtaPhi(1./cosh(m_mc_eta),m_mc_eta,m_mc_phi);
704 for (const xAOD::CaloCluster* theCluster : *clusColl) {
705 double mx_calib_tot;
706 if( !theCluster->retrieveMoment( xAOD::CaloCluster::ENG_CALIB_TOT, mx_calib_tot) ) {
707 ATH_MSG_ERROR( "Moment ENG_CALIB_TOT is absent" );
708 }
709 double mx_dens, mx_center_lambda;
710 if (!theCluster->retrieveMoment(xAOD::CaloCluster::CENTER_LAMBDA, mx_center_lambda))
711 mx_center_lambda = 0;
712 if (!theCluster->retrieveMoment(xAOD::CaloCluster::FIRST_ENG_DENS, mx_dens))
713 mx_dens = 0;
714
715 HepLorentzVector hlv_cls(1,0,0,1);
716 hlv_cls.setREtaPhi(1./cosh(theCluster->eta()), theCluster->eta(), theCluster->phi());
717 double r = hlv_pion.angle(hlv_cls.vect());
718 if(r < m_distance_cut*angle_mollier_factor(m_mc_eta)
719 && mx_calib_tot > engClusCalibThs
720 && theCluster->size() > 1
721 && mx_dens > 0
722 && mx_center_lambda != 0.0 ) {
723 PionWasReconstructed = true;
724 break;
725 }
726 }
727 if(PionWasReconstructed){
728 m_RecoEfficiency_vs_ebeam[0][m_mc_etabin]->Fill(m_mc_ener);
729 m_RecoEfficiency_vs_eta[0][m_mc_enerbin]->Fill(m_mc_eta);
730 }
731 m_RecoEfficiency_vs_ebeam[1][m_mc_etabin]->Fill(m_mc_ener);
732 m_RecoEfficiency_vs_eta[1][m_mc_enerbin]->Fill(m_mc_eta);
733 if(m_useRecoEfficiency && !PionWasReconstructed) return StatusCode::SUCCESS;
734 }
735
736 if(m_doEngRecSpect) {
737 m_engPionSpect[m_mc_enerbin][m_mc_etabin]->Fill(m_mc_ener);
738 }
739
740 if(m_doEngRecOverTruth || m_doEngTag || m_doEngRecSpect || m_doEngNoiseClus){
741 fill_reco (*clusColl, ctx);
742 }
743
744 if(m_doClusMoments) {
745 fill_moments (*clusColl, ctx);
746 }
747
748 if(m_doCalibHitsValidation) {
749 fill_calibhits (*clusColl, ctx);
750 }
751
752 return StatusCode::SUCCESS;
753}
754
755
756
757/* ****************************************************************************
758+ average reconstructed energy after different steps of calibration, energy spectras
759+ classification performance
760+ noise clusters
761**************************************************************************** */
762int GetLCSinglePionsPerf::fill_reco (const xAOD::CaloClusterContainer& clusColl,
763 const EventContext& ctx)
764{
765 /* ********************************************
766 reading cluster collections for each calibration level
767 ******************************************** */
768 std::vector<const xAOD::CaloClusterContainer *> clusCollVector;
769 for (const SG::ReadHandleKey<xAOD::CaloClusterContainer> & k : m_clusterCollNames) {
770 SG::ReadHandle<xAOD::CaloClusterContainer> pColl (k, ctx);
771 clusCollVector.push_back(pColl.cptr());
772 }
773
774 HepLorentzVector hlv_pion(1,0,0,1);
775 hlv_pion.setREtaPhi(1./cosh(m_mc_eta),m_mc_eta,m_mc_phi);
776
777 std::vector<double > engClusSum; // for 4 collections
778 engClusSum.resize(m_ncluscoll, 0.0);
779 std::vector<double > engClusSumPresOnly; // for 4 collections
780 engClusSumPresOnly.resize(m_ncluscoll, 0.0);
781
782 double engClusSumCalib(0);
783 double engClusSumCalibPresOnly(0);
784 double engClusSumTrueOOC(0);
785 double engClusSumTrueDM(0);
786 std::vector<double> engClusSumCalibTagged;
787 engClusSumCalibTagged.resize(m_ntagcases, 0.0);
788
789 int nGoodClus = 0;
790 unsigned int iClus = -1;
791 for (const xAOD::CaloCluster* theCluster : clusColl) {
792 ++iClus;
793 double mx_calib_tot;
794 if( !theCluster->retrieveMoment( xAOD::CaloCluster::ENG_CALIB_TOT, mx_calib_tot) ) {
795 ATH_MSG_ERROR( "Moment ENG_CALIB_TOT is absent" );
796 }
797 double mx_calib_emb0 = 0;
798 double mx_calib_eme0 = 0;
799 double mx_calib_tileg3 = 0;
800 if( !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_EMB0, mx_calib_emb0)
801 || !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_EME0, mx_calib_eme0)
802 || !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_TILEG3, mx_calib_tileg3)){
803 ATH_MSG_ERROR( "One of the moment ENG_CALIB_EMB0, ENG_CALIB_EME0, ENG_CALIB_TILEG3 is absent" );
804 }
805 double mx_calib_out = 0;
806 double mx_calib_dead_tot = 0;
807 if( !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_OUT_L, mx_calib_out)
808 || !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_DEAD_TOT, mx_calib_dead_tot) ){
809 ATH_MSG_ERROR( "One of the moment ENG_CALIB_OUT_L, ENG_CALIB_DEAD_TOT is absent" );
810 }
811
812 double mx_dens, mx_center_lambda;
813 if (!theCluster->retrieveMoment(xAOD::CaloCluster::CENTER_LAMBDA, mx_center_lambda))
814 mx_center_lambda = 0;
815 if (!theCluster->retrieveMoment(xAOD::CaloCluster::FIRST_ENG_DENS, mx_dens))
816 mx_dens = 0;
817
818 bool clusIsGood(true);//Fabian
819 if(m_useGoodClus){
820 if(theCluster->size() <= 1 || mx_dens <=0 || mx_center_lambda ==0 || ( m_usePionClusters && mx_calib_tot <= m_truthPiEngFraction * m_mc_ener )) {
821 clusIsGood = false;
822 }
823 }
824
825 HepLorentzVector hlv_cls(1,0,0,1);
826 hlv_cls.setREtaPhi(1./cosh(theCluster->eta()), theCluster->eta(), theCluster->phi());
827 double r = hlv_pion.angle(hlv_cls.vect());
828
829 if((!m_usePionClusters && r < m_distance_cut*angle_mollier_factor(m_mc_eta)) ||
830 ( ( m_usePionClusters && r < m_distance_cut*angle_mollier_factor(m_mc_eta)) && clusIsGood ) ) {
831
832 theCluster = clusCollVector[kTOPO]->at(iClus); // this collection has RecoStatus
833 CaloRecoStatus recoStatus = theCluster->recoStatus();
834
835 if(recoStatus.checkStatus(CaloRecoStatus::TAGGEDEM) ) {
836 engClusSumCalibTagged[kTAGEM] += mx_calib_tot;
837 }else if(recoStatus.checkStatus(CaloRecoStatus::TAGGEDHAD) ) {
838 engClusSumCalibTagged[kTAGHAD] += mx_calib_tot;
839 }else if(recoStatus.checkStatus(CaloRecoStatus::TAGGEDUNKNOWN) ) {
840 engClusSumCalibTagged[kTAGUNK] += mx_calib_tot;
841 }else{
842 ATH_MSG_ERROR("CheckSinglePionsReco::execute() -> Error! Unkown classification status " << recoStatus.getStatusWord());
843 }
844
845 if(clusIsGood){
846 nGoodClus++;
847 engClusSumCalib += mx_calib_tot;
848 engClusSumCalibPresOnly += (mx_calib_emb0+mx_calib_eme0+mx_calib_tileg3);
849 engClusSumTrueOOC += mx_calib_out;
850 engClusSumTrueDM += mx_calib_dead_tot;
851 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++){
852 theCluster = clusCollVector[i_coll]->at(iClus);
853 engClusSum[i_coll] += theCluster->e();
854 engClusSumPresOnly[i_coll] += (theCluster->eSample(CaloSampling::PreSamplerB)+theCluster->eSample(CaloSampling::PreSamplerE)+theCluster->eSample(CaloSampling::TileGap3));
855 }
856 }
857
858 }else{
859 if(m_doEngNoiseClus){
860 // to save noise clusters
861 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++){
862 theCluster = clusCollVector[i_coll]->at(iClus);
863 m_engNoiseClusSpect[i_coll][m_mc_etabin]->Fill(theCluster->e());
864 m_engNoiseClus_vs_eta[i_coll]->Fill(theCluster->eta(), theCluster->e());
865 }
866 } // m_doEngNoiseClus
867 } // r
868 // xxx
869 } // iclus
870
871 if(nGoodClus == 0) return 0;
872
873 /* ********************************************
874 filling histograms for ereco over truth after each step of calibration
875 ******************************************** */
876 if(m_doEngRecOverTruth){
877 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++) {
878 // with respect to the initial pion energy
879 int i_norm = 0;
880 double enom = engClusSum[i_coll];
881 if(i_coll != kTOPO) {
882 enom = enom - engClusSumPresOnly[i_coll];
883 }
884 if(m_mc_ener != 0.0) {
885 m_engRecOverTruth_vs_eta[i_coll][i_norm][m_mc_enerbin]->Fill(m_mc_eta, enom/m_mc_ener);
886 m_engRecOverTruth_vs_ebeam[i_coll][i_norm][m_mc_etabin]->Fill(m_mc_ener, enom/m_mc_ener);
887 }
888 // with respect to the calibration level
889 i_norm = 1;
890 double edenom = engClusSumCalib - engClusSumCalibPresOnly; // for em or had
891 if(i_coll == kTOPO_OOC){ // ooc
892 edenom += engClusSumTrueOOC;
893 }else if(i_coll == kTOPO){ // dm
894 edenom += (engClusSumTrueOOC + engClusSumTrueDM);
895 }
896 if(edenom != 0.0) {
897 m_engRecOverTruth_vs_eta[i_coll][i_norm][m_mc_enerbin]->Fill(m_mc_eta, enom/edenom);
898 m_engRecOverTruth_vs_ebeam[i_coll][i_norm][m_mc_etabin]->Fill(m_mc_ener, enom/edenom);
899 }
900 // with respect to the calibration level with "ideal" previous step
901 i_norm = 2;
902 enom = engClusSum[i_coll] - engClusSumPresOnly[i_coll];
903 if(i_coll == kTOPO_OOC){ // ooc
904 double engClusSumRecoOOC = engClusSum[kTOPO_OOC] - engClusSum[kTOPO_W];
905 enom = engClusSumCalib - engClusSumCalibPresOnly + engClusSumRecoOOC;
906 }else if(i_coll == kTOPO){ // dm
907 double engClusSumRecoDM = engClusSum[kTOPO] - (engClusSum[kTOPO_OOC]-engClusSumPresOnly[kTOPO_OOC]);
908 enom = (engClusSumCalib - engClusSumCalibPresOnly) + engClusSumTrueOOC + engClusSumRecoDM;
909 }
910 if(edenom != 0.0) {
911 m_engRecOverTruth_vs_eta[i_coll][i_norm][m_mc_enerbin]->Fill(m_mc_eta, enom/edenom);
912 m_engRecOverTruth_vs_ebeam[i_coll][i_norm][m_mc_etabin]->Fill(m_mc_ener, enom/edenom);
913 }
914 } // i_coll
915 } // m_doEngRecOverTruth
916
917 /* ********************************************
918 filling histograms with classification results
919 ******************************************** */
920 if(m_doEngTag) {
921 if(engClusSumCalib!=0.0) {
922 const double inv_engClusSumCalib = 1. / engClusSumCalib;
923 for(int i_tag=0; i_tag<m_ntagcases; i_tag++) {
924 m_engTag_vs_eta[i_tag][m_mc_enerbin]->Fill(m_mc_eta, engClusSumCalibTagged[i_tag]*inv_engClusSumCalib);
925 m_engTag_vs_ebeam[i_tag][m_mc_etabin]->Fill(m_mc_ener, engClusSumCalibTagged[i_tag]*inv_engClusSumCalib);
926 }
927 }
928 } // m_doEngTag
929
930 /* ********************************************
931 filling histograms for reco energy spectras
932 ******************************************** */
933 if(m_doEngRecSpect && m_mc_ener > 0){
934 for(int i_coll=0; i_coll<m_ncluscoll; i_coll++){
935 double enom = engClusSum[i_coll];
936 if(i_coll != kTOPO) {
937 enom = enom - engClusSumPresOnly[i_coll];
938 }
939 m_engRecSpect[i_coll][m_mc_enerbin][m_mc_etabin]->Fill(enom/m_mc_ener);
940 }
941 } // m_doEngRecSpect
942
943 return 0;
944}
945
946
947
948/* ****************************************************************************
949to check moments assignment
950**************************************************************************** */
951int GetLCSinglePionsPerf::fill_moments (const xAOD::CaloClusterContainer& clusColl,
952 const EventContext& ctx)
953{
954 ATH_MSG_DEBUG("GetLCSinglePionsPerf::fill_moments got cont. size: "<<clusColl.size());
955 /* ********************************************
956 reading calibration containers
957 ******************************************** */
958 std::vector<const CaloCalibrationHitContainer *> v_cchc;
959 for (const SG::ReadHandleKey<CaloCalibrationHitContainer>& k : m_CalibrationHitContainerNames) {
960 SG::ReadHandle<CaloCalibrationHitContainer> cchc (k, ctx);
961 v_cchc.push_back(cchc.cptr());
962 }
963
964 std::vector<const CaloCalibrationHitContainer *> v_dmcchc;
965 for (const SG::ReadHandleKey<CaloCalibrationHitContainer>& k : m_DMCalibrationHitContainerNames) {
966 SG::ReadHandle<CaloCalibrationHitContainer> cchc (k, ctx);
967 v_dmcchc.push_back(cchc.cptr());
968 }
969
970 // calculate total calibration energy int active+inactive hits
971 double engCalibTot = 0.0;
972 std::vector<double > engCalibTotSmp;
973 engCalibTotSmp.resize(CaloSampling::Unknown, 0.0);
974 for (const CaloCalibrationHitContainer* cchc : v_cchc) {
975 //loop over cells in calibration container
976 for (const CaloCalibrationHit* hit : *cchc) {
977 Identifier myId = hit->cellID();
978 if(!myId.is_valid()) {
979 ATH_MSG_WARNING("Error! Bad identifier " << myId << " in container '" << cchc->Name() << "',"
980 << " AtlasDetectorID says '" << m_id_helper->show_to_string(myId) << "'");
981 continue;
982 }
983 if(!m_id_helper->is_lar(myId) && !m_id_helper->is_tile(myId)) {
984 ATH_MSG_WARNING("Error! Bad identifier (nor lar or tile) " << myId << " in container '" << cchc->Name() << "',"
985 << " AtlasDetectorID says '" << m_id_helper->show_to_string(myId) << "'");
986 continue;
987 }
988
989 CaloSampling::CaloSample nsmp = CaloSampling::CaloSample(m_calo_id->calo_sample(myId));
990 engCalibTot += hit->energyTotal();
991 engCalibTotSmp[nsmp] += hit->energyTotal();
992 }
993 }
994
995
996 double engCalibDeadTot = 0.0;
997 std::vector<double > engCalibDeadTotInArea;
998 engCalibDeadTotInArea.resize(CaloDmDescrArea::DMA_MAX,0.0);
999 for (const CaloCalibrationHitContainer* dmcchc : v_dmcchc) {
1000 //loop over cells in calibration container
1001 for (const CaloCalibrationHit* hit : *dmcchc) {
1002 Identifier myId = hit->cellID();
1003 if(!myId.is_valid()) {
1004 ATH_MSG_ERROR("GetLCSinglePionsPerf::fill_moments() -> Error! Bad dead material identifier " << myId << " in container '" << dmcchc->Name() << "',"
1005 << " AtlasDetectorID says '" << m_id_helper->show_to_string(myId) << "'");
1006 continue;
1007 }
1008 if(!m_id_helper->is_lar_dm(myId) && !m_id_helper->is_tile_dm(myId)) {
1009 ATH_MSG_ERROR("GetLCSinglePionsPerf::fill_moments() -> Error! Bad dead material identifier (nor lar_dm or tile_dm) " << myId << " in container '" << dmcchc->Name() << "',"
1010 << " AtlasDetectorID says '" << m_id_helper->show_to_string(myId) << "'");
1011 continue;
1012 }
1013 engCalibDeadTot += hit->energyTotal();
1014 //
1015 CaloDmDescrElement* myCDDE(nullptr);
1016 myCDDE = m_caloDmDescrManager->get_element(myId);
1017 if( !myCDDE ) {
1018 ATH_MSG_ERROR("GetLCSinglePionsPerf::fill_moments() -> Error! Absent DM description element!");
1019 continue;
1020 }
1021 int nDmArea = m_caloDmDescrManager->get_dm_area(myId);
1022 engCalibDeadTotInArea[nDmArea] += hit->energyTotal();
1023 engCalibDeadTotInArea[CaloDmDescrArea::DMA_ALL] += hit->energyTotal();
1024 }
1025 }
1026// double engCalibDeadTotWithPres = engCalibDeadTot + engCalibTotSmp[CaloSampling::PreSamplerB]
1027// + engCalibTotSmp[CaloSampling::PreSamplerE]
1028// + engCalibTotSmp[CaloSampling::TileGap3];
1029
1030 /* ********************************************
1031 reading clusters moments
1032 ******************************************** */
1033 double engClusSumCalib(0);
1034 double engClusSumCalibPresOnly(0);
1035 std::vector<std::vector<double > > clsMoments; // [iClus][m_nmoments]
1036 clsMoments.resize(clusColl.size());
1037 std::vector<double > clsMomentsSum; // [m_nmoments]
1038 clsMomentsSum.resize(m_nmoments,0.0);
1039
1040 // retrieving cluster moments
1041 unsigned int iClus = -1;
1042 for (const xAOD::CaloCluster* theCluster : clusColl) {
1043 ++iClus;
1044 clsMoments[iClus].resize(m_validMoments.size(),0.0);
1045 int iMoment=0;
1046 for(moment_name_vector::const_iterator im=m_validMoments.begin(); im!=m_validMoments.end(); ++im, ++iMoment){
1047 double mx;
1048 if( theCluster->retrieveMoment((*im).second, mx) ) {
1049 clsMoments[iClus][iMoment] = mx;
1050 }else{
1051 ATH_MSG_ERROR("GetLCSinglePionsPerf::fill_moments() -> Error! Can't retrieve moment " << (*im).first << " " << (*im).second);
1052 }
1053 clsMomentsSum[iMoment] += clsMoments[iClus][iMoment];
1054 }
1055 double mx;
1056 if(theCluster->retrieveMoment( xAOD::CaloCluster::ENG_CALIB_TOT, mx ) ){
1057 engClusSumCalib += mx;
1058 }else{
1059 ATH_MSG_ERROR("GetLCSinglePionsPerf::fill_moments() -> Error! Can't retrieve moment xAOD::CaloCluster::ENG_CALIB_TOT ");
1060 }
1061 double mx_calib_emb0, mx_calib_eme0, mx_calib_tileg3;
1062 if( !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_EMB0, mx_calib_emb0)
1063 || !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_EME0, mx_calib_eme0)
1064 || !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_TILEG3, mx_calib_tileg3)){
1065 ATH_MSG_WARNING("One of the moment ENG_CALIB_EMB0, ENG_CALIB_EME0, ENG_CALIB_TILEG3 is absent");
1066 }else{
1067 engClusSumCalibPresOnly += (mx_calib_emb0+mx_calib_eme0+mx_calib_tileg3);
1068 }
1069 } // iClus
1070
1071 double engCalibDeadTotWithPres = engCalibDeadTot + engClusSumCalibPresOnly;
1072
1073 // DM energy before barrel presampler, inside it, and between presampler and strips
1074 double eng_calib_dead_emb0 = engCalibDeadTotInArea[CaloDmDescrArea::DMA_EMB0]
1075 + engCalibDeadTotInArea[CaloDmDescrArea::DMA_EMB1]
1076 + engCalibTotSmp[CaloSampling::PreSamplerB];
1077 // DM energy between barrel and tile
1078 double eng_calib_dead_tile0 = engCalibDeadTotInArea[CaloDmDescrArea::DMA_EMB3_TILE0];
1079 // DM energy before scintillator and inside scintillator
1080 double eng_calib_dead_tileg3 = engCalibDeadTotInArea[CaloDmDescrArea::DMA_SCN]
1081 + engCalibTotSmp[CaloSampling::TileGap3];
1082 // DM energy beforee endcap presampler, inside it and between presampler and strips
1083 double eng_calib_dead_eme0 = engCalibDeadTotInArea[CaloDmDescrArea::DMA_EME0]
1084 + engCalibDeadTotInArea[CaloDmDescrArea::DMA_EME12]
1085 + engCalibTotSmp[CaloSampling::PreSamplerE];
1086 // DM energy between emec and hec
1087 double eng_calib_dead_hec0 = engCalibDeadTotInArea[CaloDmDescrArea::DMA_EME3_HEC0];
1088 // DM energy before FCAL and between HEC and FCAL
1089 double eng_calib_dead_fcal = engCalibDeadTotInArea[CaloDmDescrArea::DMA_FCAL0]
1090 + engCalibDeadTotInArea[CaloDmDescrArea::DMA_HEC_FCAL];
1091 // DM leakage behind the calorimeter
1092 double eng_calib_dead_leakage = engCalibDeadTotInArea[CaloDmDescrArea::DMA_LEAK];
1093 // the rest of DM energy which remains unclassified
1094 double eng_calib_dead_unclass = engCalibDeadTotWithPres - eng_calib_dead_emb0 - eng_calib_dead_tile0
1095 - eng_calib_dead_tileg3 - eng_calib_dead_eme0 - eng_calib_dead_hec0 - eng_calib_dead_fcal
1096 - eng_calib_dead_leakage;
1097
1098 // filling histograms
1099 int iMoment=0;
1100 for(moment_name_vector::const_iterator im=m_validMoments.begin(); im!=m_validMoments.end(); ++im, ++iMoment){
1101 double xnorm = 0.0;
1102 switch ( (*im).second ) {
1103 case xAOD::CaloCluster::ENG_CALIB_TOT:
1104 xnorm = clsMomentsSum[iMoment];
1105 break;
1106 case xAOD::CaloCluster::ENG_CALIB_OUT_T:
1107 xnorm = engCalibTot - engClusSumCalib;
1108 break;
1109 case xAOD::CaloCluster::ENG_CALIB_OUT_L:
1110 xnorm = engCalibTot - engClusSumCalib;
1111 break;
1112 case xAOD::CaloCluster::ENG_CALIB_OUT_M:
1113 xnorm = engCalibTot - engClusSumCalib;
1114 break;
1115 case xAOD::CaloCluster::ENG_CALIB_DEAD_TOT:
1116 xnorm = engCalibDeadTotWithPres;
1117 break;
1118 case xAOD::CaloCluster::ENG_CALIB_DEAD_EMB0:
1119 xnorm = eng_calib_dead_emb0;
1120 break;
1121 case xAOD::CaloCluster::ENG_CALIB_DEAD_TILE0:
1122 xnorm = eng_calib_dead_tile0;
1123 break;
1124 case xAOD::CaloCluster::ENG_CALIB_DEAD_TILEG3:
1125 xnorm = eng_calib_dead_tileg3;
1126 break;
1127 case xAOD::CaloCluster::ENG_CALIB_DEAD_EME0:
1128 xnorm = eng_calib_dead_eme0;
1129 break;
1130 case xAOD::CaloCluster::ENG_CALIB_DEAD_HEC0:
1131 xnorm = eng_calib_dead_hec0;
1132 break;
1133 case xAOD::CaloCluster::ENG_CALIB_DEAD_FCAL:
1134 xnorm = eng_calib_dead_fcal;
1135 break;
1136 case xAOD::CaloCluster::ENG_CALIB_DEAD_LEAKAGE:
1137 xnorm = eng_calib_dead_leakage;
1138 break;
1139 case xAOD::CaloCluster::ENG_CALIB_DEAD_UNCLASS:
1140 xnorm = eng_calib_dead_unclass;
1141 break;
1142 default:
1143 ATH_MSG_ERROR("GetLCSinglePionsPerf::fill_moments() -> Error! Not implemented for " << (*im).first << " " << (*im).second);
1144 break;
1145 }
1146 if(m_doClusMoments && xnorm > m_mc_ener*0.0001 && xnorm > 0) {
1147 // moments assigned to first 3 maximum clusters
1148 const double inv_xnorm = 1. / xnorm;
1149 for(unsigned int i_cls=0; i_cls<clusColl.size(); i_cls++){
1150 m_clusMoment_vs_eta[iMoment][i_cls][m_mc_enerbin]->Fill(m_mc_eta, clsMoments[i_cls][iMoment]*inv_xnorm);
1151 m_clusMoment_vs_ebeam[iMoment][i_cls][m_mc_etabin]->Fill(m_mc_ener, clsMoments[i_cls][iMoment]*inv_xnorm);
1152 if(i_cls>=2) break;
1153 }
1154 // sum of moments assigned to clusters wrt total sum of moment
1155 m_clusMoment_vs_eta[iMoment][kMOM_SUMCLS][m_mc_enerbin]->Fill(m_mc_eta, clsMomentsSum[iMoment]/xnorm);
1156 m_clusMoment_vs_ebeam[iMoment][kMOM_SUMCLS][m_mc_etabin]->Fill(m_mc_ener, clsMomentsSum[iMoment]/xnorm);
1157 // unassigned energy wrt pion energy
1158 m_clusMoment_vs_eta[iMoment][kMOM_MISS][m_mc_enerbin]->Fill(m_mc_eta, (xnorm-clsMomentsSum[iMoment])/m_mc_ener);
1159 m_clusMoment_vs_ebeam[iMoment][kMOM_MISS][m_mc_etabin]->Fill(m_mc_ener, (xnorm-clsMomentsSum[iMoment])/m_mc_ener);
1160 // unassigned energy wrt pion energy
1161 m_clusMoment_vs_eta[iMoment][kMOM_MCENER][m_mc_enerbin]->Fill(m_mc_eta, xnorm/m_mc_ener);
1162 m_clusMoment_vs_ebeam[iMoment][kMOM_MCENER][m_mc_etabin]->Fill(m_mc_ener, xnorm/m_mc_ener);
1163 }
1164 } // i_mom
1165
1166 return 0;
1167}
1168
1169
1170
1171/* ****************************************************************************
1172calibration hits validation
1173**************************************************************************** */
1174int GetLCSinglePionsPerf::fill_calibhits (const xAOD::CaloClusterContainer& clusColl,
1175 const EventContext& ctx)
1176{
1177 /* ********************************************
1178 reading calibration containers
1179 ******************************************** */
1180 std::vector<const CaloCalibrationHitContainer *> v_cchc;
1181 for (const SG::ReadHandleKey<CaloCalibrationHitContainer>& k : m_CalibrationHitContainerNames) {
1182 SG::ReadHandle<CaloCalibrationHitContainer> cchc (k, ctx);
1183 v_cchc.push_back(cchc.cptr());
1184 }
1185
1186 std::vector<const CaloCalibrationHitContainer *> v_dmcchc;
1187 for (const SG::ReadHandleKey<CaloCalibrationHitContainer>& k : m_DMCalibrationHitContainerNames) {
1188 SG::ReadHandle<CaloCalibrationHitContainer> cchc (k, ctx);
1189 v_dmcchc.push_back(cchc.cptr());
1190 }
1191
1192 // calculate total calibration energy int active+inactive hits
1193 double engCalibTot = 0.0;
1194 for (const CaloCalibrationHitContainer* cchc : v_cchc) {
1195 //loop over cells in calibration container
1196 for (const CaloCalibrationHit* hit : *cchc) {
1197 Identifier myId = hit->cellID();
1198 if(!myId.is_valid()) {
1199 ATH_MSG_ERROR("Error! Bad identifier " << myId << " in container '" << cchc->Name() << "',"
1200 << " AtlasDetectorID says '" << m_id_helper->show_to_string(myId) << "'");
1201 continue;
1202 }
1203 if(!m_id_helper->is_lar(myId) && !m_id_helper->is_tile(myId)) {
1204 ATH_MSG_ERROR("Error! Bad identifier (nor lar or tile) " << myId << " in container '" << cchc->Name() << "',"
1205 << " AtlasDetectorID says '" << m_id_helper->show_to_string(myId) << "'");
1206 continue;
1207 }
1208 engCalibTot += hit->energyTotal();
1209 }
1210 }
1211
1212 // calibration energy in dead material
1213 double engCalibDeadTot = 0.0;
1214 double engDefaultCalculator = 0.0;
1215 for (const CaloCalibrationHitContainer* dmcchc : v_dmcchc) {
1216 //loop over cells in calibration container
1217 for (const CaloCalibrationHit* hit : *dmcchc) {
1218 Identifier myId = hit->cellID();
1219 if(!myId.is_valid()) {
1220 ATH_MSG_ERROR("GetLCSinglePionsPerf::fill_moments() -> Error! Bad dead material identifier " << myId << " in container '" << dmcchc->Name() << "',"
1221 << " AtlasDetectorID says '" << m_id_helper->show_to_string(myId) << "'");
1222 continue;
1223 }
1224 if(!m_id_helper->is_lar_dm(myId) && !m_id_helper->is_tile_dm(myId)) {
1225 ATH_MSG_ERROR("GetLCSinglePionsPerf::fill_moments() -> Error! Bad dead material identifier (nor lar_dm or tile_dm) " << myId << " in container '" << dmcchc->Name() << "',"
1226 << " AtlasDetectorID says '" << m_id_helper->show_to_string(myId) << "'");
1227 continue;
1228 }
1229 engCalibDeadTot += hit->energyTotal();
1230 //
1231 CaloDmDescrElement* myCDDE(nullptr);
1232 myCDDE = m_caloDmDescrManager->get_element(myId);
1233 if( !myCDDE ) {
1234 ATH_MSG_ERROR("GetLCSinglePionsPerf::fill_moments() -> Error! Absent DM description element!");
1235 continue;
1236 }
1237 int nDmArea = m_caloDmDescrManager->get_dm_area(myId);
1238 if(nDmArea == CaloDmDescrArea::DMA_DEFCALC) engDefaultCalculator += hit->energyTotal();
1239 }
1240 }
1241
1242 // calibration energy assigned to one or another cluster
1243 // retrieving cluster moments
1244 double engCalibAssigned = 0.0;
1245 double engClusSumCalibPresOnly = 0.0;
1246 for (const xAOD::CaloCluster* theCluster : clusColl) {
1247 double mx_calib_tot, mx_calib_ooc, mx_calib_dm;
1248// if( !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_TOT, mx_calib_tot)
1249// || !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_OUT_L, mx_calib_ooc)
1250// || !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_DEAD_TOT mx_calib_dm)){
1251// ATH_MSG_ERROR("One of the moment ENG_CALIB_TOT, ENG_CALIB_OUT_L, ENG_CALIB_DEAD_TOT is absent");
1252// }
1253 if( !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_TOT, mx_calib_tot) ){
1254 ATH_MSG_WARNING("Moment ENG_CALIB_TOT is absent");
1255 }
1256 if( !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_OUT_L, mx_calib_ooc) ){
1257 ATH_MSG_WARNING("Moment ENG_CALIB_OUT_L is absent" );
1258 }
1259 if( !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_DEAD_TOT, mx_calib_dm) ){
1260 ATH_MSG_WARNING("Moment ENG_CALIB_DEAD_TOT is absent");
1261 }
1262 engCalibAssigned += (mx_calib_tot + mx_calib_ooc + mx_calib_dm);
1263
1264 double mx_calib_emb0, mx_calib_eme0, mx_calib_tileg3;
1265 if( !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_EMB0, mx_calib_emb0)
1266 || !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_EME0, mx_calib_eme0)
1267 || !theCluster->retrieveMoment(xAOD::CaloCluster::ENG_CALIB_TILEG3, mx_calib_tileg3)){
1268 ATH_MSG_WARNING("One of the moment ENG_CALIB_EMB0, ENG_CALIB_EME0, ENG_CALIB_TILEG3 is absent");
1269 }else{
1270 engClusSumCalibPresOnly += (mx_calib_emb0+mx_calib_eme0+mx_calib_tileg3);
1271 }
1272
1273 } // iClus
1274 // energy in presamplers is included already in the dead material
1275 engCalibAssigned -= engClusSumCalibPresOnly;
1276
1277 int n_enerbin = int( (log10(m_mc_ener) - m_logenermin)/ ((m_logenermax - m_logenermin)/m_nlogenerbin2));
1278 // sum of moments assigned to clusters wrt total sum of moment
1279 m_SumCalibHitOverEbeam_vs_eta[n_enerbin]->Fill(m_mc_eta, (engCalibTot+engCalibDeadTot)/m_mc_ener);
1280// m_SumCalibHitOverEbeam_vs_ebeam[m_mc_etabin]->Fill(m_mc_ener, (engCalibTot+engCalibDeadTot)/m_mc_ener);
1281 m_SumCalibHitAssignedOverEbeam_vs_eta[n_enerbin]->Fill(m_mc_eta, engCalibAssigned/m_mc_ener);
1282 m_SumCalibHitAssignedOverEbeam_vs_etaphi[n_enerbin]->Fill(m_mc_eta, m_mc_phi, engCalibAssigned/m_mc_ener);
1283
1284 m_DefaultCalculatorOverEbeam_vs_eta[n_enerbin]->Fill(m_mc_eta, engDefaultCalculator/m_mc_ener);
1285// m_DefaultCalculatorOverEbeam_vs_ebeam[m_mc_etabin]->Fill(m_mc_ener, engDefaultCalculator/m_mc_ener);
1286
1287 return 0;
1288}
1289
1290
1291
1292/* ****************************************************************************
1293
1294**************************************************************************** */
1295double GetLCSinglePionsPerf::angle_mollier_factor(double x)
1296{
1297 double eta = fabs(x);
1298 double ff;
1299 if(eta<1.6){
1300 ff = atan(5.0*1.7/(200.0*cosh(eta)));
1301 }else if(eta<3.2){
1302 ff = atan(5.0*1.6/(420./tanh(eta)));
1303 }else{
1304 ff = atan(5.0*0.95/(505./tanh(eta)));
1305 }
1306 return ff*(1./atan(5.0*1.7/200.0));
1307}
1308
1309
1310
M_PI
#define M_PI
Definition ActiveFraction.h:11
eta
Scalar eta() const
pseudorapidity method
Definition AmgMatrixBasePlugin.h:83
ATH_CHECK
#define ATH_CHECK
Evaluate an expression and check for errors.
Definition AthCheckMacros.h:40
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition AthMsgStreamMacros.h:33
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition AthMsgStreamMacros.h:31
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition AthMsgStreamMacros.h:32
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition AthMsgStreamMacros.h:29
CaloCell_ID.h
errorcheck.h
Helpers for checking error return status codes and reporting errors.
MeV
static const double MeV
Definition HIEventSelectionToolRun3.cxx:72
HepMCHelpers.h
ATLAS-specific HepMC functions.
x
#define x
AthAlgorithm::AthAlgorithm
AthAlgorithm(const std::string &name, ISvcLocator *pSvcLocator)
Constructor.
Definition AthAlgorithm.cxx:10
AthCommonAlgorithm< Gaudi::Algorithm >::declareProperty
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)
Definition AthCommonDataStore.h:145
AthCommonAlgorithm< Gaudi::Algorithm >::evtStore
ServiceHandle< StoreGateSvc > & evtStore()
Definition AthCommonDataStore.h:85
AthCommonAlgorithm< Gaudi::Algorithm >::detStore
const ServiceHandle< StoreGateSvc > & detStore() const
Definition AthCommonDataStore.h:95
CaloCalibrationHit
Class to store calorimeter calibration hit.
Definition CaloCalibrationHit.h:24
CaloDmDescrManager::instance
static const CaloDmDescrManager * instance()
Definition CaloDmDescrManager.cxx:67
CaloRecoStatus
reconstruction status indicator
Definition CaloRecoStatus.h:12
CaloRecoStatus::getStatusWord
virtual const store_type & getStatusWord() const
retrieve the entire status word
Definition CaloRecoStatus.h:118
CaloRecoStatus::checkStatus
virtual bool checkStatus(const StatusIndicator &statusIndicator) const
Check status.
Definition CaloRecoStatus.h:123
DataVector::at
const T * at(size_type n) const
Access an element, as an rvalue.
DataVector::size
size_type size() const noexcept
Returns the number of elements in the collection.
GetLCSinglePionsPerf::m_SumCalibHitAssignedOverEbeam_vs_etaphi
std::vector< TProfile2D * > m_SumCalibHitAssignedOverEbeam_vs_etaphi
Definition GetLCSinglePionsPerf.h:129
GetLCSinglePionsPerf::m_DMCalibrationHitContainerNames
SG::ReadHandleKeyArray< CaloCalibrationHitContainer > m_DMCalibrationHitContainerNames
Definition GetLCSinglePionsPerf.h:67
GetLCSinglePionsPerf::m_clusterCollNames
SG::ReadHandleKeyArray< xAOD::CaloClusterContainer > m_clusterCollNames
Definition GetLCSinglePionsPerf.h:63
GetLCSinglePionsPerf::finalize
virtual StatusCode finalize()
Definition GetLCSinglePionsPerf.cxx:481
GetLCSinglePionsPerf::fill_moments
int fill_moments(const xAOD::CaloClusterContainer &clusColl, const EventContext &ctx)
Definition GetLCSinglePionsPerf.cxx:951
GetLCSinglePionsPerf::m_DefaultCalculatorOverEbeam_vs_eta
std::vector< TProfile * > m_DefaultCalculatorOverEbeam_vs_eta
Definition GetLCSinglePionsPerf.h:126
GetLCSinglePionsPerf::initialize
virtual StatusCode initialize()
Definition GetLCSinglePionsPerf.cxx:173
GetLCSinglePionsPerf::m_calo_id
const CaloCell_ID * m_calo_id
Definition GetLCSinglePionsPerf.h:59
GetLCSinglePionsPerf::m_validMoments
moment_name_vector m_validMoments
Definition GetLCSinglePionsPerf.h:68
GetLCSinglePionsPerf::m_RecoEfficiency_vs_eta
std::vector< TH1F * > m_RecoEfficiency_vs_eta[2]
Definition GetLCSinglePionsPerf.h:121
GetLCSinglePionsPerf::m_RecoEfficiency_vs_ebeam
std::vector< TH1F * > m_RecoEfficiency_vs_ebeam[2]
Definition GetLCSinglePionsPerf.h:122
GetLCSinglePionsPerf::m_id_helper
const AtlasDetectorID * m_id_helper
Definition GetLCSinglePionsPerf.h:58
GetLCSinglePionsPerf::m_engRecOverTruth_vs_ebeam
std::vector< std::vector< std::vector< TProfile * > > > m_engRecOverTruth_vs_ebeam
Definition GetLCSinglePionsPerf.h:106
GetLCSinglePionsPerf::m_engPionSpect
std::vector< std::vector< TH1F * > > m_engPionSpect
Definition GetLCSinglePionsPerf.h:113
GetLCSinglePionsPerf::m_engRecOverTruth_vs_eta
std::vector< std::vector< std::vector< TProfile * > > > m_engRecOverTruth_vs_eta
Definition GetLCSinglePionsPerf.h:105
GetLCSinglePionsPerf::angle_mollier_factor
static double angle_mollier_factor(double x)
Definition GetLCSinglePionsPerf.cxx:1295
GetLCSinglePionsPerf::m_SumCalibHitAssignedOverEbeam_vs_eta
std::vector< TProfile * > m_SumCalibHitAssignedOverEbeam_vs_eta
Definition GetLCSinglePionsPerf.h:128
GetLCSinglePionsPerf::m_clusterBasicCollName
SG::ReadHandleKey< xAOD::CaloClusterContainer > m_clusterBasicCollName
Definition GetLCSinglePionsPerf.h:62
GetLCSinglePionsPerf::m_engNoiseClusSpect
std::vector< std::vector< TH1F * > > m_engNoiseClusSpect
Definition GetLCSinglePionsPerf.h:115
GetLCSinglePionsPerf::moment_name_pair
std::pair< std::string, xAOD::CaloCluster::MomentType > moment_name_pair
Definition GetLCSinglePionsPerf.h:50
GetLCSinglePionsPerf::execute
virtual StatusCode execute(const EventContext &ctx)
Execute method.
Definition GetLCSinglePionsPerf.cxx:654
GetLCSinglePionsPerf::m_engTag_vs_eta
std::vector< std::vector< TProfile * > > m_engTag_vs_eta
Definition GetLCSinglePionsPerf.h:108
GetLCSinglePionsPerf::m_clusMoment_vs_eta
std::vector< std::vector< std::vector< TProfile * > > > m_clusMoment_vs_eta
Definition GetLCSinglePionsPerf.h:118
GetLCSinglePionsPerf::m_CalibrationHitContainerNames
SG::ReadHandleKeyArray< CaloCalibrationHitContainer > m_CalibrationHitContainerNames
Definition GetLCSinglePionsPerf.h:66
GetLCSinglePionsPerf::m_SumCalibHitOverEbeam_vs_eta
std::vector< TProfile * > m_SumCalibHitOverEbeam_vs_eta
Definition GetLCSinglePionsPerf.h:124
GetLCSinglePionsPerf::m_engRecSpect
std::vector< std::vector< std::vector< TH1F * > > > m_engRecSpect
Definition GetLCSinglePionsPerf.h:111
GetLCSinglePionsPerf::GetLCSinglePionsPerf
GetLCSinglePionsPerf(const std::string &name, ISvcLocator *pSvcLocator)
Definition GetLCSinglePionsPerf.cxx:59
GetLCSinglePionsPerf::fill_reco
int fill_reco(const xAOD::CaloClusterContainer &clusColl, const EventContext &ctx)
Definition GetLCSinglePionsPerf.cxx:762
GetLCSinglePionsPerf::m_caloDmDescrManager
const CaloDmDescrManager * m_caloDmDescrManager
Definition GetLCSinglePionsPerf.h:60
GetLCSinglePionsPerf::m_engNoiseClus_vs_eta
std::vector< TProfile * > m_engNoiseClus_vs_eta
Definition GetLCSinglePionsPerf.h:116
GetLCSinglePionsPerf::m_clusMoment_vs_ebeam
std::vector< std::vector< std::vector< TProfile * > > > m_clusMoment_vs_ebeam
Definition GetLCSinglePionsPerf.h:119
GetLCSinglePionsPerf::fill_calibhits
int fill_calibhits(const xAOD::CaloClusterContainer &clusColl, const EventContext &ctx)
Definition GetLCSinglePionsPerf.cxx:1174
GetLCSinglePionsPerf::m_engTag_vs_ebeam
std::vector< std::vector< TProfile * > > m_engTag_vs_ebeam
Definition GetLCSinglePionsPerf.h:109
Identifier::is_valid
bool is_valid() const
Check if id is in a valid state.
McEventCollection
This defines the McEventCollection, which is really just an ObjectVector of McEvent objectsFile: Gene...
Definition McEventCollection.h:31
SG::ReadHandleKey
Property holding a SG store/key/clid from which a ReadHandle is made.
Definition StoreGate/StoreGate/ReadHandleKey.h:40
SG::ReadHandle::cptr
const_pointer_type cptr()
Dereference the pointer.
xAOD::CaloCluster_v1::ENG_CALIB_OUT_M
@ ENG_CALIB_OUT_M
Attached Calibration Hit energy outside clusters but inside the calorimeter with medium matching (Ang...
Definition CaloCluster_v1.h:203
xAOD::CaloCluster_v1::ENG_CALIB_OUT_L
@ ENG_CALIB_OUT_L
Attached Calibration Hit energy outside clusters but inside the calorimeter with loose matching (Angl...
Definition CaloCluster_v1.h:199
xAOD::CaloCluster_v1::ENG_CALIB_DEAD_UNCLASS
@ ENG_CALIB_DEAD_UNCLASS
Attached Calibration Hit energy in dead material in unclassified areas of the detector.
Definition CaloCluster_v1.h:249
xAOD::CaloCluster_v1::ENG_CALIB_DEAD_HEC0
@ ENG_CALIB_DEAD_HEC0
Attached Calibration Hit energy in dead material between EME3 and HEC0.
Definition CaloCluster_v1.h:241
xAOD::CaloCluster_v1::FIRST_ENG_DENS
@ FIRST_ENG_DENS
First Moment in E/V.
Definition CaloCluster_v1.h:146
xAOD::CaloCluster_v1::ENG_CALIB_DEAD_TILEG3
@ ENG_CALIB_DEAD_TILEG3
Attached Calibration Hit energy in dead material before scintillator.
Definition CaloCluster_v1.h:235
xAOD::CaloCluster_v1::CENTER_LAMBDA
@ CENTER_LAMBDA
Shower depth at Cluster Centroid.
Definition CaloCluster_v1.h:139
xAOD::CaloCluster_v1::ENG_CALIB_DEAD_EME0
@ ENG_CALIB_DEAD_EME0
Attached Calibration Hit energy in dead material before EME0, between EME0 and EME1.
Definition CaloCluster_v1.h:238
xAOD::CaloCluster_v1::ENG_CALIB_DEAD_TILE0
@ ENG_CALIB_DEAD_TILE0
Attached Calibration Hit energy in dead material between EMB3 and TILE0.
Definition CaloCluster_v1.h:233
xAOD::CaloCluster_v1::ENG_CALIB_DEAD_FCAL
@ ENG_CALIB_DEAD_FCAL
Attached Calibration Hit energy in dead material before FCAL, between FCAL and HEC.
Definition CaloCluster_v1.h:244
xAOD::CaloCluster_v1::ENG_CALIB_TOT
@ ENG_CALIB_TOT
Calibration Hit energy inside the cluster.
Definition CaloCluster_v1.h:198
xAOD::CaloCluster_v1::ENG_CALIB_OUT_T
@ ENG_CALIB_OUT_T
Attached Calibration Hit energy outside clusters but inside the calorimeter with tight matching (Angl...
Definition CaloCluster_v1.h:207
xAOD::CaloCluster_v1::ENG_CALIB_DEAD_LEAKAGE
@ ENG_CALIB_DEAD_LEAKAGE
Attached Calibration Hit energy in dead material behind calorimeters.
Definition CaloCluster_v1.h:246
xAOD::CaloCluster_v1::ENG_CALIB_EMB0
@ ENG_CALIB_EMB0
Calibration Hit energy inside the cluster barrel presampler.
Definition CaloCluster_v1.h:221
xAOD::CaloCluster_v1::ENG_CALIB_EME0
@ ENG_CALIB_EME0
Calibration Hit energy inside the cluster endcap presampler.
Definition CaloCluster_v1.h:223
xAOD::CaloCluster_v1::ENG_CALIB_DEAD_EMB0
@ ENG_CALIB_DEAD_EMB0
Attached Calibration Hit energy in dead material before EMB0, between EMB0 and EMB1.
Definition CaloCluster_v1.h:230
xAOD::CaloCluster_v1::ENG_CALIB_DEAD_TOT
@ ENG_CALIB_DEAD_TOT
Attached Calibration Hit energy in dead material.
Definition CaloCluster_v1.h:227
xAOD::CaloCluster_v1::ENG_CALIB_TILEG3
@ ENG_CALIB_TILEG3
Calibration Hit energy inside the cluster scintillator.
Definition CaloCluster_v1.h:225
r
int r
Definition globals.cxx:22
HepMC::ConstGenParticlePtr
HepMC3::ConstGenParticlePtr ConstGenParticlePtr
Definition GenParticle.h:20
MC::PI0
static const int PI0
Definition HepMCHelpers.h:112
MC::PIPLUS
static const int PIPLUS
Definition HepMCHelpers.h:110
MC::isGenStable
bool isGenStable(const T &p)
Determine if the particle is stable at the generator (not det-sim) level,.
Definition HepMCHelpers.h:60
xAOD::CaloCluster
CaloCluster_v1 CaloCluster
Define the latest version of the calorimeter cluster class.
Definition Event/xAOD/xAODCaloEvent/xAODCaloEvent/CaloCluster.h:19
xAOD::CaloClusterContainer
CaloClusterContainer_v1 CaloClusterContainer
Define the latest version of the calorimeter cluster container.
Definition Event/xAOD/xAODCaloEvent/xAODCaloEvent/CaloClusterContainer.h:17
tree
TChain * tree
Definition tile_monitor.h:30
Generated on for ATLAS Offline Software by doxygen 1.16.1