ATLAS Offline Software
Loading...
Searching...
No Matches
SUSYToolsAlg.cxx
Go to the documentation of this file.
1/*
2 Copyright (C) 2002-2026 CERN for the benefit of the ATLAS collaboration
3*/
4
5// SUSYToolsAlg.cxx
6// Base class
7#include "SUSYTools/SUSYToolsAlg.h"
8
9// EDM includes
10#include "xAODEventInfo/EventInfo.h"
11#include "xAODMuon/MuonContainer.h"
12#include "xAODEgamma/ElectronContainer.h"
13#include "xAODEgamma/PhotonContainer.h"
14#include "xAODJet/JetContainer.h"
15#include "xAODTau/TauJetContainer.h"
16#include "xAODMissingET/MissingETContainer.h"
17#include "xAODMissingET/MissingETAuxContainer.h"
18#include "xAODTrigger/TrigCompositeContainer.h"
19
20// For shallow copy containers
21#include "xAODCore/ShallowCopy.h"
22
23// For the forcing of the tau truth container build
24#include "TauAnalysisTools/ITauTruthMatchingTool.h"
25
26// GRL
27#include "AsgAnalysisInterfaces/IGoodRunsListSelectionTool.h"
28
29// For handling systematics
30#include "AsgMessaging/StatusCode.h"
31#include "PATInterfaces/SystematicSet.h"
32
33// For finding calibration files
34#include "PathResolver/PathResolver.h"
35
36// For CBK access
37#include "xAODCutFlow/CutBookkeeper.h"
38#include "xAODCutFlow/CutBookkeeperContainer.h"
39#ifndef XAOD_STANDALONE // For now metadata is Athena-only
40#include "AthAnalysisBaseComps/AthAnalysisHelper.h"
41#include "EventInfo/EventStreamInfo.h"
42#endif
43
44// Common accessors
45#include "SUSYTools/SUSYObjDef_xAOD.h"
46
47// For output of histograms
48#include "TH1.h"
49#include "TH2.h"
50#include "TFile.h"
51#include "TEfficiency.h"
52//#include <AnaAlgorithm/IHistogramWorker.h>
53
54// For configuration
55#include "TString.h"
56#include "TEnv.h"
57#include "THashList.h"
58
59// For string manipulation
60#include <regex>
61#include <cmath>
62#include <atomic>
63
64//====================================================================================================
65namespace Cut {
66 const unsigned int NSel=5;
67 enum { all=0, baseline, passOR, signal, trigmatch };
68 static const TString SCut[] = {"All","Baseline","passOR","Signal","TrigM"};
69}
70
71//====================================================================================================
72std::string findInReg(SG::AuxTypeRegistry& reg, xAOD::IParticle* obj, const std::string& searchkey);
73int getSize(std::map<std::string,std::vector<std::string>> &collection, const std::string& object);
74
75//====================================================================================================
76// Constructor
77SUSYToolsAlg::SUSYToolsAlg(const std::string& name,
78 ISvcLocator* pSvcLocator )
79 : EL::AnaAlgorithm (name, pSvcLocator)
80 , m_SUSYTools("ST::SUSYObjDef_xAOD/SUSYTools",this)
81 , m_tauTruthMatchingTool("TauAnalysisTools::TauTruthMatchingTool/TauTruthMatchingTool")
82 , m_Nevts(0)
83 , m_kernel("StreamDAOD")
84{
85
86 declareProperty( "DoSyst", m_doSyst = false );
87 declareProperty( "GRLFiles", m_GRLFiles );
88 declareProperty( "maxEvts", m_maxEvts = -999 );
89 declareProperty( "LumiBlockFilter", m_lbfilter = 90 );
90 declareProperty( "ConfigFile", m_configFile = "SUSYTools/SUSYTools_Default.conf" );
91
92 // asg Tool Handles must be dealt with differently
93 m_tauTruthMatchingTool.declarePropertyFor( this, "TauTruthMatchingTool", "The TTMT" );
94 m_SUSYTools.declarePropertyFor( this, "SUSYTools", "The SUSYTools instance" );
95}
96
97//====================================================================================================
98// Destructor
99SUSYToolsAlg::~SUSYToolsAlg() { }
100
101//====================================================================================================
102// Initialize
103StatusCode SUSYToolsAlg::initialize() {
104 ATH_MSG_INFO("Initializing " << name() << "...");
105
106 ATH_CHECK(m_SUSYTools.retrieve());
107 ATH_MSG_INFO("Retrieved tool: " << m_SUSYTools->name() );
108
109 // read some of the property values (may have been autoconfigured at this point)
110 m_mcCampaign = *(m_SUSYTools->getProperty<std::string>("mcCampaign"));
111 m_isPHYSLITE = *(m_SUSYTools->getProperty<bool>("IsPHYSLITE"));
112
113
114 // Need truth matching for tau CP tools
115 if( !m_SUSYTools->isData() ){
116 m_tauTruthMatchingTool.setTypeAndName("TauAnalysisTools::TauTruthMatchingTool/TauTruthMatchingTool");
117 ATH_CHECK( m_tauTruthMatchingTool.setProperty("TruthJetContainerName", "AntiKt4TruthDressedWZJets") );
118 ATH_CHECK( m_tauTruthMatchingTool.retrieve() );
119 ATH_MSG_INFO("Retrieved tool: " << m_tauTruthMatchingTool->name() );
120 }
121
122 // GRL tool
123 if (m_SUSYTools->isData()) {
124 m_grl.setTypeAndName("GoodRunsListSelectionTool/grl");
125 m_grl.isUserConfigured();
126
127 std::vector<std::string> myGRLs;
128 for ( const auto& x : m_GRLFiles ) myGRLs.push_back(x);
129
130 ANA_CHECK( m_grl.setProperty("GoodRunsListVec", myGRLs) );
131 ANA_CHECK( m_grl.setProperty("PassThrough", (myGRLs.size()>0?false:true)) );
132 ANA_CHECK( m_grl.retrieve() );
133
134 ATH_MSG_INFO("GRL tool retrieved & initialized: " << m_grl->name() );
135 }
136
137 // Systematics
138 m_sysInfoList.clear();
139 // this is the nominal set
140 if (m_doSyst) {
141 m_sysInfoList = m_SUSYTools->getSystInfoList();
142 } else {
143 ST::SystInfo infodef;
144 infodef.affectsKinematics = false;
145 infodef.affectsWeights = false;
146 infodef.affectsType = ST::Unknown;
147 m_sysInfoList.push_back(infodef);
148 }
149 // order in object groups
150 groupSysts();
151
152 // setup triggers to be checked
153 m_triggers.clear();
154 if(m_mcCampaign.find("mc20") != std::string::npos){
155 // Trigger for Run2
156 m_triggers["el"] = {};
157 m_triggers["el"].push_back("HLT_e24_lhtight_nod0_ivarloose");
158 m_triggers["el"].push_back("HLT_e60_lhmedium_nod0");
159 m_triggers["el"].push_back("HLT_e26_lhtight_nod0_ivarloose");
160 m_triggers["el"].push_back("HLT_e140_lhloose_nod0");
161 m_triggers["el"].push_back("HLT_e24_lhtight_ivarloose_L1EM22VHI");
162 m_triggers["el"].push_back("HLT_e24_lhtight_ivarloose_L1eEM26M");
163 m_triggers["el"].push_back("HLT_e26_lhtight_ivarloose_L1EM22VHI");
164 m_triggers["el"].push_back("HLT_e26_lhtight_ivarloose_L1eEM26M");
165 m_triggers["el"].push_back("HLT_e60_lhvloose_L1eEM26M");
166 m_triggers["el"].push_back("HLT_e60_lhvloose_L1EM22VHI");
167
168 m_triggers["ph"] = {};
169 m_triggers["ph"].push_back("HLT_g120_loose");
170 m_triggers["ph"].push_back("HLT_g140_loose");
171 m_triggers["ph"].push_back("HLT_g120_loose_L1eEM26M");
172 m_triggers["ph"].push_back("HLT_g120_loose_L1EM22VHI");
173 m_triggers["ph"].push_back("HLT_g140_loose_L1eEM26M");
174 m_triggers["ph"].push_back("HLT_g140_loose_L1EM22VHI");
175
176 m_triggers["mu"] = {};
177 m_triggers["mu"].push_back("HLT_mu24_ivarmedium");
178 m_triggers["mu"].push_back("HLT_mu26_ivarmedium");
179 m_triggers["mu"].push_back("HLT_mu50");
180 m_triggers["mu"].push_back("HLT_mu24_ivarmedium_L1MU14FCH");
181 m_triggers["mu"].push_back("HLT_mu24_ivarmedium_L1MU18VFCH");
182 m_triggers["mu"].push_back("HLT_mu26_ivarmedium_L1MU14FCH");
183 m_triggers["mu"].push_back("HLT_mu26_ivarmedium_L1MU18VFCH");
184 }
185 else{
186 // Trigger for Run3
187 m_triggers["el"] = {};
188 m_triggers["el"].push_back("HLT_e26_lhtight_ivarloose_L1EM22VHI");
189 m_triggers["el"].push_back("HLT_e60_lhmedium_L1EM22VHI");
190 m_triggers["el"].push_back("HLT_e140_lhloose_L1EM22VHI");
191
192 m_triggers["ph"] = {};
193 m_triggers["ph"].push_back("HLT_g140_loose_L1EM22VHI");
194 m_triggers["ph"].push_back("HLT_g300_L1EM22VHI"); // HLT_g300_etcut_L1EM22VHI not working in mc21
195
196 m_triggers["mu"] = {};
197 m_triggers["mu"].push_back("HLT_mu24_ivarmedium_L1MU14FCH");
198 m_triggers["mu"].push_back("HLT_mu50_L1MU14FCH");
199 m_triggers["mu"].push_back("HLT_mu60_0eta105_msonly_L1MU14FCH");
200 m_triggers["mu"].push_back("HLT_mu60_L1MU14FCH");
201 m_triggers["mu"].push_back("HLT_mu80_msonly_3layersEC_L1MU14FCH");
202
203 }
204
205
206 // book histograms
207 ATH_CHECK( bookHistograms() );
208
209 // retrieve SUSYTools config file
210 TEnv rEnv;
211 int success = -1;
212 m_configFile = (PathResolverFindCalibFile(m_configFile)).c_str();
213 success = rEnv.ReadFile(m_configFile.c_str(), kEnvAll);
214 if(success != 0){
215 ATH_MSG_ERROR( "Cannot open config file!");
216 return StatusCode::FAILURE;
217 }
218 ATH_MSG_INFO( "Config file opened" );
219
220 m_configDict.clear();
221 m_configDict["Jet.LargeRcollection"] = rEnv.GetValue("Jet.LargeRcollection", "AntiKt10LCTopoTrimmedPtFrac5SmallR20Jets" );
222 m_configDict["TrackJet.Collection"] = rEnv.GetValue("TrackJet.Collection", "AntiKtVR30Rmax4Rmin02TrackJets" );
223 m_configDict["Jet.WtaggerConfig"] = rEnv.GetValue("Jet.WtaggerConfig", "None");
224 m_configDict["Jet.ZtaggerConfig"] = rEnv.GetValue("Jet.ZtaggerConfig", "None");
225 m_configDict["Jet.ToptaggerConfig"] = rEnv.GetValue("Jet.ToptaggerConfig", "None");
226
227 // Trim comments and extra spaces from config entries
228 std::regex comment("#.*$");
229 std::regex trimspaces("^ +| +$|( ) +");
230 for (const auto& keyval : m_configDict) {
231 m_configDict[keyval.first] = regex_replace(m_configDict[keyval.first], comment, "");
232 m_configDict[keyval.first] = regex_replace(m_configDict[keyval.first], trimspaces, "$1");
233 ATH_MSG_DEBUG("config " << keyval.first << " : " << m_configDict[keyval.first]);
234 }
235 m_FatJetCollection = m_configDict["Jet.LargeRcollection"];
236 m_TrkJetCollection = m_configDict["TrackJet.Collection"];
237 if (m_configDict["TrackJet.Collection"].find("None")!=std::string::npos) m_TrkJetCollection = "";
238
239 m_slices.clear();
240 m_slices["ele"] = bool(rEnv.GetValue("Slices.Ele", true));
241 m_slices["pho"] = bool(rEnv.GetValue("Slices.Pho", true));
242 m_slices["mu"] = bool(rEnv.GetValue("Slices.Mu", true));
243 m_slices["tau"] = bool(rEnv.GetValue("Slices.Tau", true));
244 m_slices["jet"] = bool(rEnv.GetValue("Slices.Jet", true));
245 m_slices["bjet"] = bool(rEnv.GetValue("Slices.BJet", true));
246 m_slices["fjet"] = bool(rEnv.GetValue("Slices.FJet", true));
247 m_slices["tjet"] = bool(rEnv.GetValue("Slices.TJet", true));
248 m_slices["met"] = bool(rEnv.GetValue("Slices.MET", true));
249 for (const auto& x : m_slices) { ATH_MSG_DEBUG( "Slice " << x.first << ": " << ((x.second)?"true":"false")); }
250
251 return StatusCode::SUCCESS;
252}
253
254//====================================================================================================
255// Finalize
256StatusCode SUSYToolsAlg::finalize() {
257 ATH_MSG_INFO ("Finalizing " << name() << "...");
258
259 ATH_MSG_INFO("Time per event: (CPU / WALL)");
260 ATH_MSG_INFO(" Entire loop - " << m_clock0.CpuTime() * 1e3 / m_Nevts << " / " << m_clock0.RealTime() * 1e3 / m_Nevts << " ms");
261 ATH_MSG_INFO(" First event - " << m_clock1.CpuTime() * 1e3 << " / " << m_clock1.RealTime() * 1e3 << " ms");
262 ATH_MSG_INFO(" Excluding first event - " << m_clock2.CpuTime() * 1e3 / (m_Nevts - 1) << " / " << m_clock2.RealTime() * 1e3 / (m_Nevts - 1) << " ms");
263
264
265 // label cutflow histograms
266 for (std::string hname : {"el_n_flow_nominal", "ph_n_flow_nominal", "mu_n_flow_nominal", "jet_n_flow_nominal", "bjet_n_flow_nominal", "fatjet_n_flow_nominal", "trkjet_n_flow_nominal", "tau_n_flow_nominal"} ) {
267 TH1* h = hist("Cutflow/"+hname);
268 for (unsigned int i=1; i < Cut::NSel+1; i++){ h->GetXaxis()->SetBinLabel(i, Cut::SCut[i-1].Data()); }
269 h->GetXaxis()->SetLabelSize(0.05);
270 }
271
272 // normalize and label efficiency histograms
273 for (std::string hname : {"el_trigmatch_eff_nominal", "ph_trigmatch_eff_nominal", "mu_trigmatch_eff_nominal"}) {
274 std::string obj = hname.substr(0,2);
275 TH1* h = hist("Trigger/"+hname);
276 for (unsigned int i=1; i < m_triggers[obj].size()+1; i++){ h->GetXaxis()->SetBinLabel(i, m_triggers[obj][i-1].substr(4,std::string::npos).c_str()); }
277 h->GetXaxis()->SetLabelSize(0.05);
278 if (m_obj_count[obj]["sig"]>0) h->Scale(1./static_cast<float>(m_obj_count[obj]["sig"]));
279 }
280
281 // label weights histograms
282 std::map<std::string,std::string> mapweights = {{"event","EventWeight"},{"electrons","Electron"},{"photons","Photon"},{"muons","Muon"},{"taus","Tau"},{"jets","Jet"},{"btags","BTag"}};
283 for (const auto& weight : mapweights) {
284 TH1* h = hist("Syst/weight_"+weight.first);
285 std::vector syslist = m_syst_weights[weight.second];
286 for (unsigned int i=1; i < syslist.size()+1; i++) { h->GetXaxis()->SetBinLabel(i, syslist.at(i-1).c_str()); }
287 h->GetXaxis()->SetLabelSize(0.05);
288 // average weight
289 h->Scale(1./static_cast<float>(h->GetBinContent(1)));
290
291 TH2* h2 = hist2d("Syst/syst_"+weight.first);
292 syslist = m_syst_all[weight.second];
293 for (unsigned int i=1; i < syslist.size()+1; i++) { h2->GetYaxis()->SetBinLabel(syslist.size()-i+1, syslist.at(i-1).c_str()); h2->Fill(0.5,syslist.size()-i,1.); }
294 h2->GetYaxis()->SetLabelSize(0.05);
295 }
296
297 return StatusCode::SUCCESS;
298}
299
300//====================================================================================================
301// Execute
302StatusCode SUSYToolsAlg::execute() {
303 // manual event limit
304 if (m_maxEvts>=0 && m_Nevts>=(unsigned int)m_maxEvts) return StatusCode::SUCCESS;
305
306 ATH_MSG_VERBOSE("Executing " << name() << "...");
307 m_clock0.Start(false);
308 if (m_Nevts == 0) m_clock1.Start(false);
309 else m_clock2.Start(false);
310
311 if (m_SUSYTools->resetSystematics() != StatusCode::SUCCESS) {
312 ATH_MSG_ERROR( "Cannot reset SUSYTools systematics" );
313 }
314
315 // METMaker crashes if there's no PV, check that before anything,
316 // and skip to the next event if there isn't one
317 if (m_SUSYTools->GetPrimVtx() == nullptr) {
318 ATH_MSG_WARNING("No PV found for this event! Skipping...");
319 return StatusCode::SUCCESS;
320 }
321
322 // call PRW tool to apply all relevant decorations
323 ATH_CHECK( m_SUSYTools->ApplyPRWTool());
324
325 // get EventInfo
326 const xAOD::EventInfo* evtInfo(0);
327 ATH_CHECK( evtStore()->retrieve(evtInfo, "EventInfo") );
328 bool isData = m_SUSYTools->isData();
329
330 // manual max events within lumiblock selection (on top of Nevts)
331 if (m_maxEvts>=-1 && (evtInfo->lumiBlock() != (unsigned int)m_lbfilter)) return StatusCode::SUCCESS;
332
333 // get file/kernel info on first event
334 if ( m_Nevts==0 ) {
335 std::string stream;
336#ifdef ROOTCORE
337 // Read the CutBookkeeper container on the first event
338 const xAOD::CutBookkeeperContainer* completeCBC = 0;
339 if (!evtStore()->event()->retrieveMetaInput(completeCBC, "CutBookkeepers").isSuccess()) {
340 ATH_MSG_ERROR("Failed to retrieve CutBookkeepers from MetaData! Exiting.");
341 }
342 std::string cbkname;
343 for ( const auto *cbk : *completeCBC ) {
344 cbkname = cbk->name();
345 stream = cbk->inputStream();
346 ATH_MSG_INFO("cbkname: " << cbkname << ", stream: " << stream);
347 if (cbkname.find("PHYSVAL")!=std::string::npos) m_kernel="PHYSVAL";
348 else if (cbkname.find("PHYSLITE")!=std::string::npos) m_kernel="PHYSLITE";
349 else if (cbkname.find("PHYS")!=std::string::npos) m_kernel="PHYS";
350 else if (stream.find("StreamDAOD")!=std::string::npos) m_kernel=stream.replace(stream.find("Stream"),6,"");
351 else m_kernel="other";
352 }
353#else
354 // check the stream info
355 const EventStreamInfo* esi = 0;
356 CHECK( inputMetaStore()->retrieve(esi) );
357 for ( auto stream : esi->getProcessingTags() ) {
358 if (stream.find("PHYSVAL")!=std::string::npos) m_kernel="PHYSVAL";
359 else if (stream.find("PHYSLITE")!=std::string::npos) m_kernel="PHYSLITE";
360 else if (stream.find("PHYS")!=std::string::npos) m_kernel="PHYS";
361 else if (stream.find("StreamDAOD")!=std::string::npos) m_kernel=stream.replace(stream.find("Stream"),6,"");
362 else m_kernel="other";
363 }
364#endif
365 ATH_MSG_INFO("Input file kernel type: " << m_kernel);
366 }
367
368 int period = 100;
369 if (m_Nevts==0 || (m_Nevts % period == 99)) ATH_MSG_INFO("===>>> start processing event #" << evtInfo->eventNumber() << ", run #" << evtInfo->runNumber() << " | " << m_Nevts << " events processed so far <<<===");
370
371 //==================================================
372 // General info
373 TH1* info_runNo = hist("EventInfo/info_runNo");
374 TH1* info_evtNo = hist("EventInfo/info_evtNo");
375 TH1* info_lumiBlock = hist("EventInfo/info_lumiBlock");
376 info_runNo->Fill(evtInfo->runNumber());
377 info_evtNo->Fill(evtInfo->eventNumber());
378 info_lumiBlock->Fill(evtInfo->lumiBlock());
379
380 // First, create the nominal containers
381 // Only create them if slice active & configured
382 // Only create them freshly if affected by systematics
383
384 xAOD::ElectronContainer* electrons_nominal(0);
385 xAOD::ShallowAuxContainer* electrons_nominal_aux(0);
386 if (m_slices["ele"]) {
387 ATH_CHECK( m_SUSYTools->GetElectrons(electrons_nominal, electrons_nominal_aux,true, m_isPHYSLITE?"AnalysisElectrons":"Electrons") );
388 ATH_MSG_DEBUG( "Number of electrons: " << electrons_nominal->size() );
389 }
390
391 xAOD::PhotonContainer* photons_nominal(0);
392 xAOD::ShallowAuxContainer* photons_nominal_aux(0);
393 if (m_slices["pho"]) {
394 ATH_CHECK( m_SUSYTools->GetPhotons(photons_nominal, photons_nominal_aux, true, m_isPHYSLITE?"AnalysisPhotons":"Photons") );
395 ATH_MSG_DEBUG( "Number of photons: " << photons_nominal->size() );
396 }
397
398 xAOD::MuonContainer* muons_nominal(0);
399 xAOD::ShallowAuxContainer* muons_nominal_aux(0);
400 if (m_slices["mu"]) {
401 ATH_CHECK( m_SUSYTools->GetMuons(muons_nominal, muons_nominal_aux, true, m_isPHYSLITE?"AnalysisMuons":"Muons") );
402 ATH_MSG_DEBUG( "Number of muons: " << muons_nominal->size() );
403 }
404
405 xAOD::JetContainer* jets_nominal(0);
406 xAOD::ShallowAuxContainer* jets_nominal_aux(0);
407 if (m_slices["jet"]) {
408 ATH_CHECK( m_SUSYTools->GetJets(jets_nominal, jets_nominal_aux, true, m_isPHYSLITE?"AnalysisJets":"") );
409 ATH_MSG_DEBUG( "Number of jets: " << jets_nominal->size() );
410 }
411
412 static std::atomic<bool> doFatJets = true; // only create if slice on && collection configured
413 if (m_FatJetCollection.empty()) doFatJets = false;
414 const xAOD::JetContainer* fatjets_test(0);
415 xAOD::JetContainer* fatjets_nominal(0);
416 xAOD::ShallowAuxContainer* fatjets_nominal_aux(0);
417 if( m_slices["fjet"] && doFatJets && evtStore()->contains<xAOD::JetContainer>( m_FatJetCollection ) ){
418 if( !evtStore()->retrieve( fatjets_test, m_FatJetCollection ).isSuccess() ){
419 ATH_MSG_ERROR("Failed to retrieve xAOD::JetContainer with name " << m_FatJetCollection );
420 return StatusCode::FAILURE;
421 } else {
422 ATH_CHECK( m_SUSYTools->GetFatJets(fatjets_nominal,fatjets_nominal_aux,true,m_FatJetCollection,true) );
423 ATH_MSG_DEBUG( "Number of Large Radius jets: " << fatjets_nominal->size() );
424 }
425 } else if( m_slices["jet"] && doFatJets) {
426 ATH_MSG_DEBUG("FatJets xAOD::JetContainer with name " << m_FatJetCollection <<" not found");
427 doFatJets=false;
428 }
429
430 static std::atomic<bool> doTrkJets = true; // only create if slice on && collection configured
431 if (m_TrkJetCollection.empty()) doTrkJets = false;
432 const xAOD::JetContainer* trkjets_test(0);
433 xAOD::JetContainer* trkjets_nominal(0);
434 xAOD::ShallowAuxContainer* trkjets_nominal_aux(0);
435 if( m_slices["tjet"] && doTrkJets && evtStore()->contains<xAOD::JetContainer>( m_TrkJetCollection ) ){
436 if( !evtStore()->retrieve( trkjets_test, m_TrkJetCollection ).isSuccess() ){
437 ATH_MSG_ERROR("Failed to retrieve xAOD::JetContainer with name " << m_TrkJetCollection );
438 return StatusCode::FAILURE;
439 } else {
440 ATH_CHECK( m_SUSYTools->GetTrackJets(trkjets_nominal,trkjets_nominal_aux,true,m_TrkJetCollection) );
441 ATH_MSG_DEBUG( "Number of track jets: " << trkjets_nominal->size() );
442 }
443 } else if(m_slices["jet"] && doTrkJets) {
444 ATH_MSG_DEBUG("TrkJets xAOD::JetContainer with name " << m_TrkJetCollection <<" not found");
445 doTrkJets=false;
446 }
447
449 //if( m_slices["tau"] && !isData) {
450 // const xAOD::TauJetContainer* p_TauJets = 0;
451 // ATH_CHECK( evtStore()->retrieve(p_TauJets, "TauJets") );
452 // //if( !evtStore()->contains<xAOD::TruthParticleContainer>("TruthTaus") && p_TauJets->size()==0 ){
453 // // // If there are no taus, then we need to force the building of the container
454 // // ATH_MSG_DEBUG("FIXME: TruthTaus hack for empty TauJets.");
455 // // xAOD::TruthParticleContainer* noTT = new xAOD::TruthParticleContainer;
456 // // xAOD::TruthParticleAuxContainer* noTTAux = new xAOD::TruthParticleAuxContainer;
457 // // noTT->setStore( noTTAux );
458 // // ATH_CHECK( evtStore()->record(noTT, "TruthTaus") );
459 // // ATH_CHECK( evtStore()->record(noTTAux, "TruthTausAux.") );
460 // // ATH_CHECK( m_tauTruthMatchingTool->retrieveTruthTaus() );
461 // //}
462 // for(const auto& tau : *p_TauJets) {
463 // const xAOD::TruthParticle* trueTau = 0;
464 // trueTau = m_tauTruthMatchingTool->getTruth(*tau);
465 // if( trueTau ){
466 // ATH_MSG_DEBUG("getTruth tau " <<tau->pt() <<" " <<tau->eta()
467 // <<" " <<tau->phi()
468 // <<" trueTau " <<trueTau->pt() <<" " <<trueTau->eta()
469 // <<" " <<trueTau->phi());
470 // } else {
471 // ATH_MSG_DEBUG("getTruth tau " <<tau->pt() <<" " <<tau->eta()
472 // <<" " <<tau->phi() << "trueTau not found");
473 // }
474 // }
475 //}
476
477 const xAOD::TauJetContainer* taus_gettruth(0);
478 xAOD::TauJetContainer* taus_nominal(0);
479 xAOD::ShallowAuxContainer* taus_nominal_aux(0);
480 if (m_slices["tau"]) {
481 if (!isData && !m_isPHYSLITE) {
482 ATH_CHECK( evtStore()->retrieve(taus_gettruth,"TauJets") );
483 for(const auto tau : *taus_gettruth) {
484 m_tauTruthMatchingTool->getTruth(*tau);
485 }
486 }
487 ATH_CHECK( m_SUSYTools->GetTaus(taus_nominal, taus_nominal_aux, true, m_isPHYSLITE?"AnalysisTauJets":"TauJets") );
488 ATH_MSG_DEBUG( "Number of taus: " << taus_nominal->size() );
489 }
490
491 m_metcst_nominal = new xAOD::MissingETContainer;
492 m_metcst_nominal_aux = new xAOD::MissingETAuxContainer;
493 m_metcst_nominal->setStore(m_metcst_nominal_aux);
494 m_metcst_nominal->reserve(10);
495 m_mettst_nominal = new xAOD::MissingETContainer;
496 m_mettst_nominal_aux = new xAOD::MissingETAuxContainer;
497 m_mettst_nominal->setStore(m_mettst_nominal_aux);
498 m_mettst_nominal->reserve(10);
499
500 double metsig_cst(0.);
501 double metsig_tst(0.);
502 if (m_slices["met"]) {
503 // (met, jets, electrons, muons, photons, taus, doTST, doJVT)
504 ATH_CHECK( m_SUSYTools->GetMET(*m_metcst_nominal, jets_nominal, electrons_nominal, muons_nominal, photons_nominal, 0, false, false) );
505 ATH_MSG_DEBUG("RefFinal CST etx=" << (*m_metcst_nominal)["Final"]->mpx()
506 << ", ety=" << (*m_metcst_nominal)["Final"]->mpy()
507 << ", et=" << (*m_metcst_nominal)["Final"]->met()
508 << ", sumet=" << (*m_metcst_nominal)["Final"]->sumet());
509 ATH_CHECK( m_SUSYTools->GetMETSig(*m_metcst_nominal, metsig_cst, false, false) );
510 ATH_MSG_DEBUG("METSignificance = " << metsig_cst);
511
512 ATH_CHECK( m_SUSYTools->GetMET(*m_mettst_nominal, jets_nominal, electrons_nominal, muons_nominal, photons_nominal, 0, true, true) ); // bugfix for PHYSLITE with p5511
513 ATH_MSG_DEBUG("RefFinal TST etx=" << (*m_mettst_nominal)["Final"]->mpx()
514 << ", ety=" << (*m_mettst_nominal)["Final"]->mpy()
515 << ", et=" << (*m_mettst_nominal)["Final"]->met()
516 << ", sumet=" << (*m_mettst_nominal)["Final"]->sumet());
517 ATH_CHECK( m_SUSYTools->GetMETSig(*m_mettst_nominal, metsig_tst, true, true) );
518 ATH_MSG_DEBUG("METSignificance = " << metsig_tst);
519
520 hist("MET/met_et") ->Fill( (*m_mettst_nominal)["Final"]->met()*0.001 );
521 hist("MET/met_sumet") ->Fill( (*m_mettst_nominal)["Final"]->sumet()*0.001 );
522 hist("MET/met_phi") ->Fill( (*m_mettst_nominal)["Final"]->phi() );
523 hist("MET/met_et_tst")->Fill( (*m_mettst_nominal)["PVSoftTrk"]->met() *0.001 );
524 hist("MET/met_et_el") ->Fill( (*m_mettst_nominal)["RefEle"]->met() *0.001 );
525 hist("MET/met_et_ph") ->Fill( (*m_mettst_nominal)["RefGamma"]->met() *0.001 );
526 hist("MET/met_et_mu") ->Fill( (*m_mettst_nominal)["Muons"]->met() *0.001 );
527 hist("MET/met_et_jet")->Fill( (*m_mettst_nominal)["RefJet"]->met() *0.001 );
528 //hist("MET/met_et_tau")->Fill( (*m_mettst_nominal)["RefTau"]->met() *0.001 );
529
530 hist("MET/met_significance")->Fill( metsig_tst );
531 }
532
533 //--- Overlap Removal
534 ATH_CHECK( m_SUSYTools->OverlapRemoval(electrons_nominal, muons_nominal, jets_nominal, photons_nominal, taus_nominal) );
535
536 //--- Trigger
537 bool isRun3Trig = false;
538 if (m_mcCampaign.find("mc21") != std::string::npos || m_mcCampaign.find("mc23") != std::string::npos) isRun3Trig = true;
539
540 hist("Trigger/isMETTrigPassed")->Fill( double(m_SUSYTools->IsMETTrigPassed()) );
541
542 //--- Monitoring
543 for (const auto& obj : m_objects) { for (const auto& lev : m_levels) { m_obj_count[obj][lev] = 0; } }
544
545 //----- Electrons
546
547 TH1* el_n_flow_nominal = hist("Cutflow/el_n_flow_nominal");
548 TH1* el_trigmatch_eff_nominal = hist("Trigger/el_trigmatch_eff_nominal");
549
550 if (m_slices["ele"]) {
551 ATH_MSG_DEBUG("Processing electrons");
552 m_obj_count["el"]["nom"] = electrons_nominal->size();
553 for(auto el : *electrons_nominal) {
554 el_n_flow_nominal->Fill(Cut::all);
555 if ( ST::acc_baseline(*el) == 1 ){
556 el_n_flow_nominal->Fill(Cut::baseline);
557 stdHistsForObj(el,"el","bsl");
558 m_obj_count["el"]["bsl"] += 1;
559
560 if ( ST::acc_passOR(*el) == 1 ){
561 el_n_flow_nominal->Fill(Cut::passOR);
562 if ( ST::acc_signal(*el) == 1 ){
563 el_n_flow_nominal->Fill(Cut::signal);
564 stdHistsForObj(el,"el","sig");
565 m_obj_count["el"]["sig"] += 1;
566
567 bool passTM=false;
568 unsigned int idx=1;
569 for(const auto& t : m_triggers["el"]){
570 ATH_MSG_DEBUG( " Processing trigger " << t );
571 bool passit = ((isRun3Trig||t.find("_L1")==std::string::npos) ? m_SUSYTools->IsTrigMatched(el, t) : false);
572 passTM |= passit;
573 if(passit) el_trigmatch_eff_nominal->SetBinContent(idx, el_trigmatch_eff_nominal->GetBinContent(idx)+1);
574 #ifdef XAOD_STANDALONE
575 m_heffs["Trigger/el_pt_"+t]->Fill(passit,el->pt()/1000.);
576 m_heffs["Trigger/el_eta_"+t]->Fill(passit,el->eta());
577 m_heffs["Trigger/el_phi_"+t]->Fill(passit,el->phi());
578 #endif
579 idx++;
580 }
581 if(passTM) el_n_flow_nominal->Fill(Cut::trigmatch);
582 }
583 }
584 }
585 stdHistsForObj(el,"el","nom");
586 }
587 hist("Electron/el_nom_N")->Fill(m_obj_count["el"]["nom"]);
588 hist("Electron/el_bsl_N")->Fill(m_obj_count["el"]["bsl"]);
589 hist("Electron/el_sig_N")->Fill(m_obj_count["el"]["sig"]);
590 }
591
592 //----- Photons
593
594 TH1* ph_n_flow_nominal = hist("Cutflow/ph_n_flow_nominal");
595 TH1* ph_trigmatch_eff_nominal = hist("Trigger/ph_trigmatch_eff_nominal");
596
597 if (m_slices["pho"]) {
598 ATH_MSG_DEBUG("Processing photons");
599 m_obj_count["ph"]["nom"] = photons_nominal->size();
600 for(auto ph : *photons_nominal) {
601 ph_n_flow_nominal->Fill(Cut::all);
602 if ( ST::acc_baseline(*ph) == 1 ){
603 ph_n_flow_nominal->Fill(Cut::baseline);
604 stdHistsForObj(ph,"ph","bsl");
605 m_obj_count["ph"]["bsl"] += 1;
606 if ( ST::acc_passOR(*ph) == 1 ){
607 ph_n_flow_nominal->Fill(Cut::passOR);
608 if ( ST::acc_signal(*ph) == 1 ){
609 ph_n_flow_nominal->Fill(Cut::signal);
610 stdHistsForObj(ph,"ph","sig");
611 m_obj_count["ph"]["sig"] += 1;
612
613 bool passTM=false;
614 unsigned int idx=1;
615 for(const auto& t : m_triggers["ph"]){
616 bool passit = ((isRun3Trig||t.find("_L1")==std::string::npos) ? m_SUSYTools->IsTrigMatched(ph, t) : false);
617 passTM |= passit;
618 if(passit) ph_trigmatch_eff_nominal->SetBinContent(idx, ph_trigmatch_eff_nominal->GetBinContent(idx)+1);
619 #ifdef XAOD_STANDALONE
620 m_heffs["Trigger/ph_pt_"+t]->Fill(passit,ph->pt()/1000.);
621 m_heffs["Trigger/ph_eta_"+t]->Fill(passit,ph->eta());
622 m_heffs["Trigger/ph_phi_"+t]->Fill(passit,ph->phi());
623 #endif
624 idx++;
625 }
626 if (passTM) ph_n_flow_nominal->Fill(Cut::trigmatch);
627 }
628 }
629 }
630 stdHistsForObj(ph,"ph","nom");
631 }
632 hist("Photon/ph_nom_N")->Fill(m_obj_count["ph"]["nom"]);
633 hist("Photon/ph_bsl_N")->Fill(m_obj_count["ph"]["bsl"]);
634 hist("Photon/ph_sig_N")->Fill(m_obj_count["ph"]["sig"]);
635 }
636
637 //----- Muons
638
639 TH1* mu_n_flow_nominal = hist("Cutflow/mu_n_flow_nominal");
640 TH1* mu_trigmatch_eff_nominal = hist("Trigger/mu_trigmatch_eff_nominal");
641
642 if (m_slices["mu"]) {
643 ATH_MSG_DEBUG("Processing muons");
644 m_obj_count["mu"]["nom"] = muons_nominal->size();
645 for(auto mu : *muons_nominal) {
646 mu_n_flow_nominal->Fill(Cut::all);
647 if ( ST::acc_baseline(*mu) == 1 ){
648 mu_n_flow_nominal->Fill(Cut::baseline);
649 stdHistsForObj(mu,"mu","bsl");
650 m_obj_count["mu"]["bsl"] += 1;
651 if ( ST::acc_passOR(*mu) == 1 ){
652 mu_n_flow_nominal->Fill(Cut::passOR);
653 if ( ST::acc_signal(*mu) == 1 ){
654 mu_n_flow_nominal->Fill(Cut::signal);
655 stdHistsForObj(mu,"mu","sig");
656 m_obj_count["mu"]["sig"] += 1;
657
658 bool passTM=false;
659 unsigned int idx=1;
660 for(const auto& t : m_triggers["mu"]){
661 bool passit = ((isRun3Trig||t.find("_L1")==std::string::npos) ? m_SUSYTools->IsTrigMatched(mu, t) : false);
662 passTM |= passit;
663 if(passit) mu_trigmatch_eff_nominal->SetBinContent(idx, mu_trigmatch_eff_nominal->GetBinContent(idx)+1);
664 #ifdef XAOD_STANDALONE
665 m_heffs["Trigger/mu_pt_"+t]->Fill(passit,mu->pt()/1000.);
666 m_heffs["Trigger/mu_eta_"+t]->Fill(passit,mu->eta());
667 m_heffs["Trigger/mu_phi_"+t]->Fill(passit,mu->phi());
668 #endif
669 idx++;
670 }
671 if(passTM) mu_n_flow_nominal->Fill(Cut::trigmatch);
672 }
673 }
674 }
675 stdHistsForObj(mu,"mu","nom");
676 }
677 hist("Muon/mu_nom_N")->Fill(m_obj_count["mu"]["nom"]);
678 hist("Muon/mu_bsl_N")->Fill(m_obj_count["mu"]["bsl"]);
679 hist("Muon/mu_sig_N")->Fill(m_obj_count["mu"]["sig"]);
680 }
681
682 //----- Jets
683
684 TH1* jet_n_flow_nominal = hist("Cutflow/jet_n_flow_nominal");
685 TH1* bjet_n_flow_nominal = hist("Cutflow/bjet_n_flow_nominal");
686 bool bjet = false;
687
688 if (m_slices["jet"]) {
689 m_obj_count["jet"]["nom"] = jets_nominal->size();
690 for(auto jet : *jets_nominal) {
691 ATH_MSG_DEBUG("Processing jets");
692 jet_n_flow_nominal->Fill(Cut::all);
693 bjet = false;
694 if ( ST::acc_bjet(*jet) >= 3 ) bjet = true;
695 if (bjet) bjet_n_flow_nominal->Fill(Cut::all);
696 //
697 if ( ST::acc_baseline(*jet) == 1 ){
698 jet_n_flow_nominal->Fill(Cut::baseline);
699 if (bjet) {
700 bjet_n_flow_nominal->Fill(Cut::baseline);
701 stdHistsForObj(jet,"bjet","bsl");
702 m_obj_count["bjet"]["bsl"] += 1;
703 }
704 stdHistsForObj(jet,"jet","bsl");
705 m_obj_count["jet"]["bsl"] += 1;
706
707 if ( ST::acc_passOR(*jet) == 1 ){
708 jet_n_flow_nominal->Fill(Cut::passOR);
709 if (bjet) bjet_n_flow_nominal->Fill(Cut::passOR);
710 if ( ST::acc_signal(*jet) == 1 ){
711 jet_n_flow_nominal->Fill(Cut::signal);
712 jet_n_flow_nominal->Fill(Cut::trigmatch); //no trig matching for jets
713 if (bjet) {
714 bjet_n_flow_nominal->Fill(Cut::signal);
715 bjet_n_flow_nominal->Fill(Cut::trigmatch); //no trig matching for (b)jets
716 stdHistsForObj(jet,"bjet","sig");
717 m_obj_count["bjet"]["sig"] += 1;
718 }
719 stdHistsForObj(jet,"jet","sig");
720 m_obj_count["jet"]["sig"] += 1;
721 }
722 }
723 }
724 if (bjet) {
725 stdHistsForObj(jet,"bjet","nom");
726 }
727 stdHistsForObj(jet,"jet","nom");
728 }
729 hist("Jet/jet_nom_N")->Fill(m_obj_count["jet"]["nom"]);
730 hist("Jet/jet_bsl_N")->Fill(m_obj_count["jet"]["bsl"]);
731 hist("Jet/jet_sig_N")->Fill(m_obj_count["jet"]["sig"]);
732 hist("bJet/bjet_nom_N")->Fill(m_obj_count["bjet"]["nom"]);
733 hist("bJet/bjet_bsl_N")->Fill(m_obj_count["bjet"]["bsl"]);
734 hist("bJet/bjet_sig_N")->Fill(m_obj_count["bjet"]["sig"]);
735 }
736
737 //----- Large R. Jets
738
739 TH1* fatjet_n_flow_nominal = hist("Cutflow/fatjet_n_flow_nominal");
740
741 if( m_slices["fjet"] && doFatJets ) {
742 ATH_MSG_DEBUG("Processing fatjets");
743 if (m_configDict.find("Jet.WtaggerConfig")!=m_configDict.end() && m_configDict.find("WtaggerKey")==m_configDict.end() && fatjets_nominal->size()>0) {
744 SG::AuxTypeRegistry& reg = SG::AuxTypeRegistry::instance();
745 std::string taggerKey;
746 taggerKey = findInReg(reg, (*fatjets_nominal)[0], "(SmoothW.*?)_.*");
747 if (!taggerKey.empty()) { m_configDict["WtaggerKey"] = taggerKey; }
748 }
749 if (m_configDict.find("Jet.ZtaggerConfig")!=m_configDict.end() && m_configDict.find("ZtaggerKey")==m_configDict.end() && fatjets_nominal->size()>0) {
750 SG::AuxTypeRegistry& reg = SG::AuxTypeRegistry::instance();
751 std::string taggerKey;
752 taggerKey = findInReg(reg, (*fatjets_nominal)[0], "(SmoothZ.*?)_.*");
753 if (!taggerKey.empty()) { m_configDict["ZtaggerKey"] = taggerKey; }
754 }
755 if (m_configDict.find("Jet.ToptaggerConfig")!=m_configDict.end() && m_configDict.find("ToptaggerKey")==m_configDict.end() && fatjets_nominal->size()>0) {
756 SG::AuxTypeRegistry& reg = SG::AuxTypeRegistry::instance();
757 std::string taggerKey;
758 taggerKey = findInReg(reg, (*fatjets_nominal)[0], "(.*Tag.*Top.*?)_.*");
759 if (!taggerKey.empty()) { m_configDict["ToptaggerKey"] = taggerKey; }
760 }
761 ATH_MSG_DEBUG("Tagger configuration retrieved.");
762 ATH_MSG_DEBUG("WTagger: " << m_configDict["WtaggerKey"]);
763 ATH_MSG_DEBUG("ZTagger: " << m_configDict["ZtaggerKey"]);
764 ATH_MSG_DEBUG("TopTagger: " << m_configDict["ToptaggerKey"]);
765
766 m_obj_count["fatjet"]["nom"] = fatjets_nominal->size();
767 for(auto fatjet : *fatjets_nominal) {
768 fatjet_n_flow_nominal->Fill(Cut::all);
769 if ( ST::acc_baseline(*fatjet) == 1 ){
770 fatjet_n_flow_nominal->Fill(Cut::baseline);
771 stdHistsForObj(fatjet,"fatjet","bsl",m_configDict);
772 m_obj_count["fatjet"]["bsl"] += 1;
773 if ( ST::acc_passOR(*fatjet) == 1 ){
774 fatjet_n_flow_nominal->Fill(Cut::passOR);
775 if ( ST::acc_signal(*fatjet) == 1 ){
776 fatjet_n_flow_nominal->Fill(Cut::signal);
777 fatjet_n_flow_nominal->Fill(Cut::trigmatch); //no trig matching for jets
778 stdHistsForObj(fatjet,"fatjet","sig",m_configDict);
779 m_obj_count["fatjet"]["sig"] += 1;
780 }
781 }
782 }
783 stdHistsForObj(fatjet,"fatjet","nom",m_configDict);
784 }
785 hist("LargeRJet/fatjet_nom_N")->Fill(m_obj_count["fatjet"]["nom"]);
786 hist("LargeRJet/fatjet_bsl_N")->Fill(m_obj_count["fatjet"]["bsl"]);
787 hist("LargeRJet/fatjet_sig_N")->Fill(m_obj_count["fatjet"]["sig"]);
788 }
789
790 //----- Track Jets
791
792 TH1* trkjet_n_flow_nominal = hist("Cutflow/trkjet_n_flow_nominal");
793
794 if( m_slices["tjet"] && doTrkJets ) {
795 m_obj_count["trkjet"]["nom"] = trkjets_nominal->size();
796 for(auto trkjet : *trkjets_nominal) {
797 trkjet_n_flow_nominal->Fill(Cut::all);
798 if ( ST::acc_baseline(*trkjet) == 1 ){
799 trkjet_n_flow_nominal->Fill(Cut::baseline);
800 stdHistsForObj(trkjet,"trkjet","bsl");
801 m_obj_count["trkjet"]["bsl"] += 1;
802 if ( ST::acc_passOR(*trkjet) == 1 ){
803 trkjet_n_flow_nominal->Fill(Cut::passOR);
804 if ( ST::acc_signal(*trkjet) == 1 ){
805 trkjet_n_flow_nominal->Fill(Cut::signal);
806 trkjet_n_flow_nominal->Fill(Cut::trigmatch); //no trig matching for jets
807 stdHistsForObj(trkjet,"trkjet","sig");
808 m_obj_count["trkjet"]["sig"] += 1;
809 }
810 }
811 }
812 stdHistsForObj(trkjet,"trkjet","nom");
813 }
814 hist("TrackJet/trkjet_nom_N")->Fill(m_obj_count["trkjet"]["nom"]);
815 hist("TrackJet/trkjet_bsl_N")->Fill(m_obj_count["trkjet"]["bsl"]);
816 hist("TrackJet/trkjet_sig_N")->Fill(m_obj_count["trkjet"]["sig"]);
817 }
818
819 //----- Taus
820
821 TH1* tau_n_flow_nominal = hist("Cutflow/tau_n_flow_nominal");
822
823 if (m_slices["tau"]) {
824 m_obj_count["tau"]["nom"] = taus_nominal->size();
825 for(auto tau : *taus_nominal) {
826 tau_n_flow_nominal->Fill(Cut::all);
827 if ( ST::acc_baseline(*tau) == 1 ){
828 tau_n_flow_nominal->Fill(Cut::baseline);
829 stdHistsForObj(tau,"tau","bsl");
830 m_obj_count["tau"]["bsl"] += 1;
831 if ( ST::acc_passOR(*tau) == 1 ){
832 tau_n_flow_nominal->Fill(Cut::passOR);
833
834 if ( ST::acc_signal(*tau) == 1 ){
835 tau_n_flow_nominal->Fill(Cut::signal);
836 tau_n_flow_nominal->Fill(Cut::trigmatch); //no trig matching for taus
837 stdHistsForObj(tau,"tau","sig");
838 m_obj_count["tau"]["sig"] += 1;
839 }
840 }
841 }
842 stdHistsForObj(tau,"tau","nom");
843 }
844 hist("Tau/tau_nom_N")->Fill(m_obj_count["tau"]["nom"]);
845 hist("Tau/tau_bsl_N")->Fill(m_obj_count["tau"]["bsl"]);
846 hist("Tau/tau_sig_N")->Fill(m_obj_count["tau"]["sig"]);
847 }
848
849 // Set up the event weights
850 // Base should include all weights that do not depend on individual objects
851 double base_event_weight(1.);
852 if (!isData) {
853 base_event_weight *= evtInfo->mcEventWeight();
854
855 if (m_slices["met"]) {
856 const xAOD::MissingETContainer* met_truth(0);
857 ATH_CHECK( evtStore()->retrieve(met_truth, "MET_Truth") );
858
859 ATH_MSG_DEBUG("Truth MET etx=" << (*met_truth)["NonInt"]->mpx()
860 << ", ety=" << (*met_truth)["NonInt"]->mpy()
861 << ", et=" << (*met_truth)["NonInt"]->met()
862 << ", sumet=" << (*met_truth)["NonInt"]->sumet());
863 ATH_MSG_DEBUG("CST residuals: detx=" << (*m_metcst_nominal)["Final"]->mpx() - (*met_truth)["NonInt"]->mpx()
864 << ", dety=" << (*m_metcst_nominal)["Final"]->mpy() - (*met_truth)["NonInt"]->mpy());
865 ATH_MSG_DEBUG("TST residuals: detx=" << (*m_mettst_nominal)["Final"]->mpx() - (*met_truth)["NonInt"]->mpx()
866 << ", dety=" << (*m_mettst_nominal)["Final"]->mpy() - (*met_truth)["NonInt"]->mpy());
867 }
868 }
869
870
871 // Additionally define a nominal weight for each object type
872 double electrons_weight_nominal(1.);
873 double photons_weight_nominal(1.);
874 double muons_weight_nominal(1.);
875 double jets_weight_nominal(1.);
876 //double fatjets_weight_nominal(1.);
877 //double trkjets_weight_nominal(1.);
878 double btag_weight_nominal(1.);
879 double taus_weight_nominal(1.);
880 double event_weight_nominal(1.);
881 double prw_weight_nominal(1.);
882
883 //--- Weights
884 TH1* weight_event = hist("Syst/weight_event");
885 TH1* weight_electrons = hist("Syst/weight_electrons");
886 TH1* weight_muons = hist("Syst/weight_muons");
887 TH1* weight_photons = hist("Syst/weight_photons");
888 TH1* weight_taus = hist("Syst/weight_taus");
889 TH1* weight_jets = hist("Syst/weight_jets");
890 TH1* weight_btags = hist("Syst/weight_btags");
891
892 bool isNominal(true);
893 for (const auto& sysInfo : m_sysInfoList) {
894 CP::SystematicSet sys = sysInfo.systset;
895 if (m_SUSYTools->applySystematicVariation(sys) != StatusCode::SUCCESS) {
896 ATH_MSG_ERROR( "Cannot configure SUSYTools for systematic var. %s" << sys.name() );
897 } else {
898 ATH_MSG_VERBOSE( "Variation \"" << sys.name() << "\" configured..." );
899 }
900 if (isNominal && (sysInfo.affectsKinematics || sysInfo.affectsWeights)) {
901 ATH_MSG_ERROR("Expected nominal variation!");
902 return StatusCode::FAILURE;
903 }
904 if (sysInfo.affectsKinematics || sysInfo.affectsWeights) isNominal = false;
905
906 double event_weight = base_event_weight;
907 bool syst_affectsEventWeight = (sysInfo.affectsType == ST::SystObjType::EventWeight);
908
909 // If nominal, compute the nominal weight, otherwise recompute the weight
910
911 double prw_weight=1.;
912 if(isNominal){
913 if(!isData) prw_weight_nominal = m_SUSYTools->GetPileupWeight();
914 event_weight *= prw_weight_nominal;
915 weight_event->SetBinContent(1, weight_event->GetBinContent(1)+event_weight);
916 }
917 else if (!syst_affectsEventWeight) {
918 event_weight *= prw_weight_nominal;
919 }
920 else {
921 if(!isData){
922 ATH_CHECK(m_SUSYTools->ApplyPRWTool());
923 prw_weight = m_SUSYTools->GetPileupWeight();
924 }
925 event_weight *= prw_weight;
926 size_t iwbin = find(m_syst_weights["EventWeight"].begin(), m_syst_weights["EventWeight"].end(), sys.name()) - m_syst_weights["EventWeight"].begin();
927 if(iwbin < m_syst_weights["EventWeight"].size()) { weight_event->SetBinContent(iwbin+1, weight_event->GetBinContent(iwbin+1)+event_weight); }
928 }
929
930
931 // Define the generic collection pointers
932 xAOD::ElectronContainer* electrons(electrons_nominal);
933 xAOD::PhotonContainer* photons(photons_nominal);
934 xAOD::MuonContainer* muons(muons_nominal);
935 xAOD::JetContainer* jets(jets_nominal);
936 xAOD::JetContainer* fatjets(fatjets_nominal);
937 xAOD::JetContainer* trkjets(trkjets_nominal);
938 xAOD::TauJetContainer* taus(taus_nominal);
939 xAOD::MissingETContainer *mettst(m_mettst_nominal), *metcst(m_metcst_nominal);
940 xAOD::MissingETAuxContainer *mettst_aux(m_mettst_nominal_aux), *metcst_aux(m_metcst_nominal_aux);
941
942 bool syst_affectsElectrons = ST::testAffectsObject(xAOD::Type::Electron, sysInfo.affectsType);
943 bool syst_affectsMuons = ST::testAffectsObject(xAOD::Type::Muon, sysInfo.affectsType);
944 bool syst_affectsPhotons = ST::testAffectsObject(xAOD::Type::Photon, sysInfo.affectsType);
945 bool syst_affectsTaus = ST::testAffectsObject(xAOD::Type::Tau, sysInfo.affectsType);
946 bool syst_affectsJets = ST::testAffectsObject(xAOD::Type::Jet, sysInfo.affectsType);
947 bool syst_affectsBTag = ST::testAffectsObject(xAOD::Type::BTag, sysInfo.affectsType);
948 bool syst_affectsLRT = ST::testAffectsObject(sysInfo.affectsType) == "LRT_Objects";
949
950 // Apply LRT uncertainty.
951 // Needs to be before retrieving leptons as they could be removed if the track is vetoed
952 if (syst_affectsLRT)
953 {
954 ATH_CHECK( m_SUSYTools->ApplyLRTUncertainty() );
955 }
956
957 // If necessary (kinematics affected), make a shallow copy with the variation applied
958 // otherwise set the collection pointers to the nominal
959 if (sysInfo.affectsKinematics) {
960 if (m_slices["el"] && syst_affectsElectrons) {
961 ATH_MSG_DEBUG("Get systematics-varied electrons");
962 xAOD::ElectronContainer* electrons_syst(0);
963 xAOD::ShallowAuxContainer* electrons_syst_aux(0);
964 ATH_CHECK( m_SUSYTools->GetElectrons(electrons_syst, electrons_syst_aux, true , m_isPHYSLITE?"AnalysisElectrons":"Electrons") );
965 electrons = electrons_syst;
966 }
967
968 if (m_slices["pho"] && syst_affectsPhotons) {
969 ATH_MSG_DEBUG("Get systematics-varied photons");
970 xAOD::PhotonContainer* photons_syst(0);
971 xAOD::ShallowAuxContainer* photons_syst_aux(0);
972 ATH_CHECK( m_SUSYTools->GetPhotons(photons_syst, photons_syst_aux, true, m_isPHYSLITE?"AnalysisPhotons":"Photons") );
973 photons = photons_syst;
974 }
975
976 if (m_slices["mu"] && syst_affectsMuons) {
977 ATH_MSG_DEBUG("Get systematics-varied muons");
978 xAOD::MuonContainer* muons_syst(0);
979 xAOD::ShallowAuxContainer* muons_syst_aux(0);
980 ATH_CHECK( m_SUSYTools->GetMuons(muons_syst, muons_syst_aux, true, m_isPHYSLITE?"AnalysisMuons":"Muons") );
981 muons = muons_syst;
982 }
983
984 if (syst_affectsJets) {
985 if (m_slices["jet"]) {
986 ATH_MSG_DEBUG("Get systematics-varied jets");
987 xAOD::JetContainer* jets_syst(0);
988 xAOD::ShallowAuxContainer* jets_syst_aux(0);
989 ATH_CHECK( m_SUSYTools->GetJetsSyst(*jets_nominal, jets_syst, jets_syst_aux, true, m_isPHYSLITE?"AnalysisJets":"") );
990 jets = jets_syst;
991 }
992 if (m_slices["fjet"] and doFatJets) {
993 ATH_MSG_DEBUG("Get systematics-varied fatjets");
994 xAOD::JetContainer* fatjets_syst(0);
995 xAOD::ShallowAuxContainer* fatjets_syst_aux(0);
996 ATH_CHECK( m_SUSYTools->GetFatJets(fatjets_syst, fatjets_syst_aux, true, m_FatJetCollection, true) );
997 fatjets = fatjets_syst;
998 }
999 if (m_slices["trkjet"]) {
1000 ATH_MSG_DEBUG("Get systematics-varied trkjets");
1001 xAOD::JetContainer* trkjets_syst(0);
1002 xAOD::ShallowAuxContainer* trkjets_syst_aux(0);
1003 ATH_CHECK( m_SUSYTools->GetJetsSyst(*trkjets_nominal, trkjets_syst, trkjets_syst_aux) );
1004 trkjets = trkjets_syst;
1005 }
1006 }
1007
1008 if (m_slices["tau"] && syst_affectsTaus) {
1009 ATH_MSG_DEBUG("Get systematics-varied taus");
1010 xAOD::TauJetContainer* taus_syst(0);
1011 xAOD::ShallowAuxContainer* taus_syst_aux(0);
1012 ATH_CHECK( m_SUSYTools->GetTaus(taus_syst, taus_syst_aux,true, m_isPHYSLITE?"AnalysisTauJets":"TauJets") );
1013 taus = taus_syst;
1014 }
1015
1016 ATH_MSG_DEBUG("Get systematics-varied MET");
1017 m_mettst_syst = new xAOD::MissingETContainer;
1018 m_mettst_syst_aux = new xAOD::MissingETAuxContainer;
1019 m_metcst_syst = new xAOD::MissingETContainer;
1020 m_metcst_syst_aux = new xAOD::MissingETAuxContainer;
1021 m_mettst_syst->setStore(m_mettst_syst_aux);
1022 m_metcst_syst->setStore(m_metcst_syst_aux);
1023 m_mettst_syst->reserve(10);
1024 m_metcst_syst->reserve(10);
1025 //
1026 if (m_slices["met"]) {
1027 ATH_CHECK( m_SUSYTools->GetMET(*m_metcst_syst, jets, electrons, muons, photons, 0) );
1028 ATH_CHECK( m_SUSYTools->GetMET(*m_mettst_syst, jets, electrons, muons, photons, 0, true, true) );
1029 }
1030 //
1031 mettst = m_mettst_syst;
1032 mettst_aux = m_mettst_syst_aux;
1033 metcst = m_metcst_syst;
1034 metcst_aux = m_metcst_syst_aux;
1035 }
1036
1037 //--- Overlap Removal
1038 ATH_CHECK( m_SUSYTools->OverlapRemoval(electrons, muons, jets, photons, taus) );
1039
1040 //--- Electrons
1041 if (m_slices["ele"]) {
1042 ATH_MSG_DEBUG("Working on electrons");
1043 float electrons_weight(1.);
1044 for ( const auto el : *electrons ) {
1045 if( !isData ){
1046 if (isNominal || syst_affectsElectrons) {
1047 if ((ST::acc_signal(*el) == 1) && (isNominal || sysInfo.affectsWeights)) {
1048 //electrons_weight *= m_SUSYTools->GetSignalElecSF( *el ); // (*el, true, true, false, true) to switch off trigger SF
1049 electrons_weight *= m_SUSYTools->GetSignalElecSF(*el, true, true, false, true);
1050 }
1051 }
1052 }
1053 ATH_MSG_VERBOSE( " Electron passing baseline selection? " << static_cast<int>( ST::acc_baseline(*el)));
1054 ATH_MSG_VERBOSE( " Electron passing signal selection? " << static_cast<int>( ST::acc_signal(*el)));
1055 if (ST::acc_signal(*el) == 1 && ST::acc_effscalefact.isAvailable(*el))
1056 ATH_MSG_VERBOSE( " Electron weight " << ST::acc_effscalefact(*el) );
1057
1058 }
1059 if (isNominal) {
1060 electrons_weight_nominal = electrons_weight;
1061 weight_electrons->SetBinContent(1, weight_electrons->GetBinContent(1)+electrons_weight);
1062 }
1063 else if (!syst_affectsElectrons) {
1064 electrons_weight = electrons_weight_nominal;
1065 }
1066 else if ( sysInfo.affectsWeights ){
1067 size_t iwbin = find(m_syst_weights["Electron"].begin(), m_syst_weights["Electron"].end(), sys.name()) - m_syst_weights["Electron"].begin();
1068 if(iwbin < m_syst_weights["Electron"].size()) { weight_electrons->SetBinContent(iwbin+1, weight_electrons->GetBinContent(iwbin+1)+electrons_weight); }
1069 }
1070
1071 event_weight *= electrons_weight;
1072 ATH_MSG_DEBUG("Combined electron scale factor: " << electrons_weight);
1073 }
1074
1075 //--- Photons
1076 if (m_slices["pho"]) {
1077 ATH_MSG_DEBUG("Working on photons");
1078 float photons_weight(1.);
1079 for ( const auto ph : *photons ) {
1080 if( !isData ){
1081 if (isNominal || syst_affectsPhotons) {
1082 if ((ST::acc_signal(*ph) == 1) && (isNominal || sysInfo.affectsWeights)) {
1083 photons_weight *= m_SUSYTools->GetSignalPhotonSF( *ph );
1084 }
1085 }
1086 }
1087 ATH_MSG_VERBOSE( " Photon passing baseline selection? " << static_cast<int>(ST::acc_baseline(*ph)));
1088 ATH_MSG_VERBOSE( " Photon passing signal selection? " << static_cast<int>(ST::acc_signal(*ph)));
1089 if (ST::acc_signal(*ph) == 1 && ST::acc_effscalefact.isAvailable(*ph))
1090 ATH_MSG_VERBOSE( " Photon weight " << ST::acc_effscalefact(*ph) );
1091 }
1092 if (isNominal) {
1093 photons_weight_nominal = photons_weight;
1094 weight_photons->SetBinContent(1, weight_photons->GetBinContent(1)+photons_weight);
1095 }
1096 else if (!syst_affectsPhotons) {
1097 photons_weight = photons_weight_nominal;
1098 }
1099 else if ( sysInfo.affectsWeights ){
1100 size_t iwbin = find(m_syst_weights["Photon"].begin(), m_syst_weights["Photon"].end(), sys.name()) - m_syst_weights["Photon"].begin();
1101 if(iwbin < m_syst_weights["Photon"].size()) { weight_photons->SetBinContent(iwbin+1, weight_photons->GetBinContent(iwbin+1)+photons_weight); }
1102 }
1103
1104 event_weight *= photons_weight;
1105 ATH_MSG_DEBUG("Combined photon scale factor: " << photons_weight);
1106 }
1107
1108 //--- Muons
1109 if (m_slices["mu"]) {
1110 ATH_MSG_DEBUG("Working on muons");
1111 float muons_weight(1.);
1112 for ( const auto mu : *muons ) {
1113 if( !isData ){
1114 if (isNominal || syst_affectsMuons) {
1115 if ((ST::acc_signal(*mu) == 1) && (isNominal || sysInfo.affectsWeights)) {
1116 muons_weight *= m_SUSYTools->GetSignalMuonSF(*mu);;
1117 }
1118 }
1119 }
1120 ATH_MSG_VERBOSE( " Muon passing baseline selection? " << static_cast<int>(ST::acc_baseline(*mu)));
1121 ATH_MSG_VERBOSE( " Muon passing signal selection? " << static_cast<int>(ST::acc_signal(*mu)));
1122 ATH_MSG_VERBOSE( " Muon is a cosmic ray? " << static_cast<int>(ST::acc_cosmic(*mu)));
1123 if (ST::acc_signal(*mu) == 1 && ST::acc_effscalefact.isAvailable(*mu))
1124 ATH_MSG_VERBOSE( " Muon weight " << ST::acc_effscalefact(*mu) );
1125 }
1126 if (isNominal) {
1127 muons_weight_nominal = muons_weight;
1128 weight_muons->SetBinContent(1, weight_muons->GetBinContent(1)+muons_weight);
1129 }
1130 else if (!syst_affectsMuons) {
1131 muons_weight = muons_weight_nominal;
1132 }
1133 else if ( sysInfo.affectsWeights ){
1134 size_t iwbin = find(m_syst_weights["Muon"].begin(), m_syst_weights["Muon"].end(), sys.name()) - m_syst_weights["Muon"].begin();
1135 if(iwbin < m_syst_weights["Muon"].size()) { weight_muons->SetBinContent(iwbin+1, weight_muons->GetBinContent(iwbin+1)+muons_weight); }
1136 }
1137
1138 event_weight *= muons_weight;
1139 ATH_MSG_DEBUG("Combined muon scale factor: " << muons_weight);
1140 }
1141
1142 //--- Jets
1143 if (m_slices["jet"]) {
1144 ATH_MSG_DEBUG("Working on jets");
1145 for ( const auto jet : *jets ) {
1146 ATH_MSG_VERBOSE( " Jet is bad? " << static_cast<int>(ST::acc_bad(*jet)));
1147 ATH_MSG_VERBOSE( " Jet is baseline ? " << static_cast<int>(ST::acc_baseline(*jet)));
1148 ATH_MSG_VERBOSE( " Jet passes OR ? " << static_cast<int>(ST::acc_passOR(*jet)));
1149 if (ST::acc_passOR(*jet)) {
1150 if (m_slices["bjet"] && (isNominal || syst_affectsJets || syst_affectsBTag)) {
1151 m_SUSYTools->IsBJetContinuous(*jet);
1152 }
1153 ATH_MSG_VERBOSE( " Jet is bjet ? " << static_cast<int>(ST::acc_bjet(*jet)));
1154 }
1155 }
1156
1157 float jet_weight(1.);
1158 float btag_weight(1.);
1159 if( !isData ) { //isMC
1160 if (isNominal) { //btagging
1161 btag_weight_nominal = btag_weight = m_SUSYTools->BtagSF(jets);
1162 weight_btags->SetBinContent(1, weight_btags->GetBinContent(1)+btag_weight);
1163 }
1164 else if (!syst_affectsBTag){
1165 btag_weight = btag_weight_nominal;
1166 }
1167 else{
1168 btag_weight = m_SUSYTools->BtagSF(jets);
1169 size_t iwbin = find(m_syst_weights["BTag"].begin(), m_syst_weights["BTag"].end(), sys.name()) - m_syst_weights["BTag"].begin();
1170 if(iwbin < m_syst_weights["BTag"].size()) { weight_btags->SetBinContent(iwbin+1, weight_btags->GetBinContent(iwbin+1)+btag_weight); }
1171 }
1172
1173 if(isNominal){ //JVT
1174 jet_weight = m_SUSYTools->JVT_SF(jets);
1175 jet_weight *= m_SUSYTools->FJVT_SF(jets);
1176 jets_weight_nominal = jet_weight;
1177 jets_weight_nominal *= m_SUSYTools->FJVT_SF(jets);
1178 weight_jets->SetBinContent(1, weight_jets->GetBinContent(1)+jet_weight);
1179 }
1180 else if (!syst_affectsJets || (syst_affectsJets && !sysInfo.affectsWeights)){
1181 jet_weight = jets_weight_nominal;
1182 }
1183 else if ( syst_affectsJets && sysInfo.affectsWeights ){
1184 jet_weight = m_SUSYTools->JVT_SF(jets);
1185 jet_weight *= m_SUSYTools->FJVT_SF(jets);
1186 size_t iwbin = find(m_syst_weights["Jet"].begin(), m_syst_weights["Jet"].end(), sys.name()) - m_syst_weights["Jet"].begin();
1187 if(iwbin < m_syst_weights["Jet"].size()) { weight_jets->SetBinContent(iwbin+1, weight_jets->GetBinContent(iwbin+1)+jet_weight); }
1188 }
1189 }
1190 else{ //data
1191 weight_jets->SetBinContent(1, 1.);
1192 weight_btags->SetBinContent(1, 1.);
1193 }
1194 event_weight *= jet_weight;
1195 event_weight *= btag_weight;
1196 ATH_MSG_DEBUG("Combined b-tagging scale factor: " << btag_weight);
1197 ATH_MSG_DEBUG("Combined jet scale factor: " << jet_weight);
1198 }
1199
1200 //--- Large R. Jets
1201
1202 if( m_slices["fjet"] && doFatJets ) {
1203 ATH_MSG_DEBUG("Working on fat jets");
1204 for ( const auto fatjet : *fatjets ) {
1205 ATH_MSG_VERBOSE( " Jet is bad? " << static_cast<int>(ST::acc_bad(*fatjet)));
1206 ATH_MSG_VERBOSE( " Jet is baseline ? " << static_cast<int>(ST::acc_baseline(*fatjet)));
1207 ATH_MSG_VERBOSE( " Jet passes OR ? " << static_cast<int>(ST::acc_passOR(*fatjet)));
1208 }
1209
1210 float fatjet_weight(1.);
1211 //disable - noJVTforfatjets if( !isData ) { //isMC
1212 //disable - noJVTforfatjets if(isNominal){ //JVT - no JVT cuts for fat jets!
1213 //disable - noJVTforfatjets fatjets_weight_nominal = fatjet_weight = m_SUSYTools->JVT_SF(fatjets);
1214 //disable - noJVTforfatjets weight_fatjets->SetBinContent(1, weight_fatjets->GetBinContent(1)+fatjet_weight);
1215 //disable - noJVTforfatjets }
1216 //disable - noJVTforfatjets if (!syst_affectsJets || (syst_affectsJets && !sysInfo.affectsWeights)){
1217 //disable - noJVTforfatjets fatjet_weight = fatjets_weight_nominal;
1218 //disable - noJVTforfatjets }
1219 //disable - noJVTforfatjets else if ( syst_affectsJets && sysInfo.affectsWeights ){
1220 //disable - noJVTforfatjets fatjet_weight = m_SUSYTools->JVT_SF(fatjets);
1221 //disable - noJVTforfatjets size_t iwbin = find(m_syst_weights["Jet"].begin(), m_syst_weights["Jet"].end(), sys.name()) - m_syst_weights["Jet"].begin();
1222 //disable - noJVTforfatjets if(iwbin < m_syst_weights["Jet"].size()) { weight_fatjets->SetBinContent(iwbin+1, weight_fatjets->GetBinContent(iwbin+1)+fatjet_weight); }
1223 //disable - noJVTforfatjets }
1224 //disable - noJVTforfatjets }
1225
1226 ATH_MSG_DEBUG("Combined large radius jet scale factor: " << fatjet_weight);
1227 }
1228
1229 //--- Track Jets
1230
1231 if( m_slices["tjet"] && doTrkJets ) {
1232 ATH_MSG_DEBUG("Working on trk jets");
1233 for ( const auto trkjet : *trkjets ) {
1234 ATH_MSG_VERBOSE( " Jet is bad? " << static_cast<int>(ST::acc_bad(*trkjet)));
1235 ATH_MSG_VERBOSE( " Jet is baseline ? " << static_cast<int>(ST::acc_baseline(*trkjet)));
1236 ATH_MSG_VERBOSE( " Jet passes OR ? " << static_cast<int>(ST::acc_passOR(*trkjet)));
1237 }
1238
1239 float trkjet_weight(1.);
1240 //disable - noJVTfortrkjets if( !isData ) { //isMC
1241 //disable - noJVTfortrkjets if(isNominal){ //JVT - no JVT cuts for trk jets!
1242 //disable - noJVTfortrkjets trkjets_weight_nominal = trkjet_weight = 1;//m_SUSYTools->JVT_SF(trkjets);
1243 //disable - noJVTfortrkjets weight_trkjets->SetBinContent(1, weight_trkjets->GetBinContent(1)+trkjet_weight);
1244 //disable - noJVTfortrkjets }
1245 //disable - noJVTfortrkjets else if (!syst_affectsJets || (syst_affectsJets && !sysInfo.affectsWeights)){
1246 //disable - noJVTfortrkjets trkjet_weight = trkjets_weight_nominal;
1247 //disable - noJVTfortrkjets }
1248 //disable - noJVTfortrkjets else if ( syst_affectsJets && sysInfo.affectsWeights ){
1249 //disable - noJVTfortrkjets trkjet_weight = 1;//m_SUSYTools->JVT_SF(trkjets);
1250 //disable - noJVTfortrkjets size_t iwbin = find(m_syst_weights["Jet"].begin(), m_syst_weights["Jet"].end(), sys.name()) - m_syst_weights["Jet"].begin();
1251 //disable - noJVTfortrkjets if(iwbin < m_syst_weights["Jet"].size()) { weight_trkjets->SetBinContent(iwbin+1, weight_trkjets->GetBinContent(iwbin+1)+trkjet_weight); }
1252 //disable - noJVTfortrkjets }
1253 //disable - noJVTfortrkjets }
1254
1255 ATH_MSG_DEBUG("Combined track jet scale factor: " << trkjet_weight);
1256 }
1257
1258 //--- Taus
1259 if (m_slices["tau"]) {
1260 ATH_MSG_DEBUG("Working on taus");
1261 float taus_weight(1.);
1262 for ( const auto ta : *taus ) {
1263 if( !isData ){
1264 if (isNominal || syst_affectsTaus) {
1265 if ((ST::acc_signal(*ta) == 1) && (isNominal || sysInfo.affectsWeights)) {
1266 taus_weight *= m_SUSYTools->GetSignalTauSF(*ta, true, false); //true, "HLT_tau25_medium1_tracktwo");;
1267 }
1268 }
1269 }
1270 ATH_MSG_VERBOSE( " Tau passing baseline selection? " << static_cast<int>(ST::acc_baseline(*ta)));
1271 ATH_MSG_VERBOSE( " Tau passing signal selection? " << static_cast<int>(ST::acc_signal(*ta)));
1272 if (ST::acc_signal(*ta) == 1 && ST::acc_effscalefact.isAvailable(*ta))
1273 ATH_MSG_VERBOSE( " Tau weight " << ST::acc_effscalefact(*ta) );
1274 }
1275
1276 if (isNominal) {
1277 taus_weight_nominal = taus_weight;
1278 weight_taus->SetBinContent(1, weight_taus->GetBinContent(1)+taus_weight);
1279 }
1280 else if (!syst_affectsTaus) {
1281 taus_weight = taus_weight_nominal;
1282 }
1283 else if ( sysInfo.affectsWeights ){
1284 size_t iwbin = find(m_syst_weights["Tau"].begin(), m_syst_weights["Tau"].end(), sys.name()) - m_syst_weights["Tau"].begin();
1285 if(iwbin < m_syst_weights["Tau"].size()) { weight_taus->SetBinContent(iwbin+1, weight_taus->GetBinContent(iwbin+1)+taus_weight); }
1286 }
1287
1288 event_weight *= taus_weight;
1289 }
1290
1291
1292 ATH_MSG_DEBUG("Full event weight: " << event_weight);
1293 if (isNominal) {event_weight_nominal = event_weight;}
1294 else if (sysInfo.affectsWeights) ATH_MSG_DEBUG("Difference with nominal weight: " << event_weight - event_weight_nominal);
1295
1296 // Clean up the systematics copies
1297 if (sysInfo.affectsKinematics) {
1298 delete metcst;
1299 delete mettst;
1300 delete metcst_aux;
1301 delete mettst_aux;
1302 }
1303 isNominal = false;
1304 }
1305
1306 m_clock0.Stop();
1307 if (m_Nevts == 0) {
1308 m_clock1.Stop();
1309 } else {
1310 m_clock2.Stop();
1311 }
1312 ++m_Nevts;
1313
1314 delete m_metcst_nominal;
1315 delete m_metcst_nominal_aux;
1316 delete m_mettst_nominal;
1317 delete m_mettst_nominal_aux;
1318
1319 return StatusCode::SUCCESS;
1320}
1321
1322//====================================================================================================
1323// Book histograms
1324StatusCode SUSYToolsAlg::bookHistograms(void) {
1325 // General info
1326 ATH_CHECK( book(TH1D("EventInfo/info_runNo", ";RunNumber;N", 1500, 250000, 400000)) );
1327 ATH_CHECK( book(TH1D("EventInfo/info_evtNo", ";EventNumber;N", 600, 100000000, 400000000)) );
1328 ATH_CHECK( book(TH1D("EventInfo/info_lumiBlock", ";LumiBlock;N", 1000, 0, 1000)) );
1329
1330 // Kinematic histograms
1331 m_objects = {"el", "ph", "mu", "tau", "jet", "bjet", "fatjet", "trkjet"};
1332 m_levels = {"nom", "bsl", "sig"};
1333 m_vars = {"pt", "eta", "phi", "e", "mass", "d0", "z0", "isolCalo", "isolTrack", "truthType", "truthOrigin", "N", "bweight", "bweightpb", "bweightpc", "bweightpu", "bweightptau", "pid", "parentpid", "wtagged", "ztagged", "toptagged", "toptaggedscore","nTracks","nTracksCharged","nTracksIsolation","TauIDScoreSigTrans"};
1334 std::map<std::string,std::vector<std::string>> cfg_hist_labels = {
1335 {"pt",{"p_{T} [GeV]","N"}}, {"eta",{"#eta","N"}}, {"phi",{"#phi","N"}}, {"e",{"energy","N"}}, {"mass",{"mass","N"}},
1336 {"d0",{"d0 significance","N"}}, {"z0",{"z0 * sin(#theta)","N"}}, {"isolCalo",{"isolation (calo) / pT","N"}}, {"isolTrack",{"isolation (track) / pT","N"}},
1337 {"truthType",{"truth type","N"}}, {"truthOrigin",{"truth origin","N"}}, {"N",{"count","N"}},
1338 {"bweight",{"b-weight","N"}}, {"bweightpb",{"b-weight (pb)","N"}}, {"bweightpc",{"b-weight (pc)","N"}}, {"bweightpu",{"b-weight (pu)","N"}}, {"bweightptau", {"b-weight (ptau)","N"}},
1339 {"pid", {"pdg ID", "N"}}, {"parentpid", {"parent pdg ID","N"}},
1340 {"wtagged", {"W tagged", "N"}}, {"ztagged", {"Z tagged", "N"}}, {"toptagged", {"Top tagged", "N"}}, {"toptaggedscore", {"Top tagger score", "N"}},
1341 {"nTracks", {"N tracks", "N"}}, {"nTracksCharged", {"N tracks charged", "N"}}, {"nTracksIsolation", {"N tracks isolation", "N"}}, {"TauIDScoreSigTrans", {"TauIDScoreSigTrans","N"}} };
1342 std::map<std::string,int> cfg_hist_nbins = {
1343 {"pt",100}, {"eta",40}, {"phi",64}, {"e",100}, {"mass",100},
1344 {"d0",50}, {"z0",60}, {"isolCalo",70}, {"isolTrack",70},
1345 {"truthType",50}, {"truthOrigin",50}, {"N",16},
1346 {"bweight",200}, {"bweightpb",200}, {"bweightpc",200}, {"bweightpu",200}, {"bweightptau",200},
1347 {"pid",61}, {"parentpid",61},
1348 {"wtagged",3}, {"ztagged",3}, {"toptagged",3}, {"toptaggedscore",50},
1349 {"nTracks",20}, {"nTracksCharged",20}, {"nTracksIsolation",5}, {"TauIDScoreSigTrans",50} };
1350 std::map<std::string,std::vector<float>> cfg_hist_minmax = {
1351 {"pt",{0,200}}, {"eta",{-4,4}}, {"phi",{-3.2,3.2}}, {"e",{0,200}}, {"mass",{0,200}},
1352 {"d0",{-10,10}}, {"z0",{-1.5,1.5}}, {"isolCalo",{-0.2,0.5}}, {"isolTrack",{-0.2,0.5}},
1353 {"truthType",{0,50}}, {"truthOrigin",{0,50}}, {"N",{0,16}},
1354 {"bweight",{-10,10}}, {"bweightpb",{-0.5,1.5}}, {"bweightpc",{-0.5,1.5}}, {"bweightpu",{-0.5,1.5}}, {"bweightptau",{-0.5,1.5}},
1355 {"pid",{-30.5,30.5}}, {"parentpid",{-30.5,30.5}},
1356 {"wtagged", {-1,2}}, {"ztagged", {-1,2}}, {"toptagged", {-1,2}}, {"toptaggedscore", {0,1}},
1357 {"nTracks", {0,20}}, {"nTracksCharged", {0,20}}, {"nTracksIsolation", {0,5}}, {"TauIDScoreSigTrans", {0,1}} };
1358 std::map<std::string,std::string> labels_objects = { {"el","Electron"}, {"ph","Photon"}, {"mu","Muon"}, {"jet","Jet"}, {"bjet","b-Jet"}, {"tau","Tau"}, {"fatjet","Large-R jet"}, {"trkjet","Track jet"} };
1359 std::map<std::string,std::string> labels_levels = { {"nom","Nominal"}, {"bsl","Baseline"}, {"sig","Signal"} };
1360 std::map<std::string,std::string> labels_dir = { {"el","Electron"}, {"ph","Photon"}, {"mu","Muon"}, {"jet","Jet"}, {"bjet","bJet"}, {"tau","Tau"}, {"fatjet","LargeRJet"}, {"trkjet","TrackJet"} };
1361
1362 for (const auto& obj : m_objects) {
1363 m_obj_count[obj] = std::map<std::string,int>();
1364 for (const auto& lev : m_levels) {
1365 m_obj_count[obj][lev] = 0;
1366 }
1367 }
1368
1369 for (const auto& obj : m_objects) {
1370 for (const auto& lev : m_levels) {
1371 for (const auto& var : m_vars) {
1372 if (var.find("bweight")!=std::string::npos && obj.compare("bjet")!=0) continue; // bweights only for bjets
1373 if (var.find("tagged")!=std::string::npos && obj.compare("fatjet")!=0) continue; // boson tagging only for fjets
1374 if ((var.find("nTracks")!=std::string::npos||var.find("ScoreSigTrans")!=std::string::npos) && obj.compare("tau")!=0) continue; // nTracks/RNN/GNTau score only for taus
1375 if (var.find("isol")!=std::string::npos && !(obj.compare("el")==0||obj.compare("mu")==0||obj.compare("ph")==0)) continue; // isol only for e/ph/mu
1376 std::string key = labels_dir[obj] + "/" + obj + "_" + lev + "_" + var;
1377 std::string labels = ";"+labels_levels[lev]+" "+labels_objects[obj]+" "+cfg_hist_labels[var][0]+";"+cfg_hist_labels[var][1];
1378 ATH_CHECK( book(TH1D(key.c_str(), labels.c_str(), cfg_hist_nbins[var], cfg_hist_minmax[var][0], cfg_hist_minmax[var][1])) );
1379 m_hists[key] = hist(key);
1380 ATH_MSG_INFO("Defined histogram: " << key.c_str() << ", " << m_hists[key]);
1381 }
1382 }
1383 }
1384
1386 #ifdef XAOD_STANDALONE
1387 for (std::string obj : {"el","mu","ph"}) {
1388 for (auto trg : m_triggers[obj]) {
1389 for (std::string var : {"pt","eta","phi"} ) {
1390 std::string key = "Trigger/"+obj+"_"+var+"_"+trg;
1391 std::string labels = ";"+labels_objects[obj]+" "+cfg_hist_labels[var][0]+";Efficiency "+trg;
1392 ATH_CHECK(book(TEfficiency(key.c_str(), labels.c_str(), cfg_hist_nbins[var], cfg_hist_minmax[var][0], cfg_hist_minmax[var][1])));
1393 m_heffs[key] = histeff(key);
1394 ATH_MSG_INFO("Defined histogram: " << key.c_str() << ", " << m_heffs[key]);
1395 }
1396 }
1397 }
1398 #endif
1399 ATH_CHECK( book(TH1D("Trigger/isMETTrigPassed", "isMETTrigPassed", 2, -0.5, 1.5)) );
1400
1401 ATH_CHECK( book(TH1D("Trigger/el_trigmatch_eff_nominal", "Electron Trigger Matching Efficiency (Nominal);Electron Trigger Matching Efficiency (Nominal);N", getSize(m_triggers,"el"), 0, getSize(m_triggers,"el")) ) );
1402 ATH_CHECK( book(TH1D("Trigger/ph_trigmatch_eff_nominal", "Photon Trigger Matching Efficiency (Nominal);Photon Trigger Matching Efficiency (Nominal);N", getSize(m_triggers,"ph"), 0, getSize(m_triggers,"ph")) ) );
1403 ATH_CHECK( book(TH1D("Trigger/mu_trigmatch_eff_nominal", "Muon Trigger Matching Efficiency (Nominal);Muon Trigger Matching Efficiency (Nominal);N", getSize(m_triggers,"mu"), 0, getSize(m_triggers,"mu")) ) );
1404
1405 // Extra histograms
1406 ATH_CHECK( book(TH1D("Cutflow/el_n_flow_nominal", "Electron Cutflow (Nominal);Electron Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1407 ATH_CHECK( book(TH1D("Cutflow/ph_n_flow_nominal", "Photon Cutflow (Nominal);Photon Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1408 ATH_CHECK( book(TH1D("Cutflow/mu_n_flow_nominal", "Muon Cutflow (Nominal);Muon Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1409 ATH_CHECK( book(TH1D("Cutflow/jet_n_flow_nominal", "Jet Cutflow (Nominal);Jet Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1410 ATH_CHECK( book(TH1D("Cutflow/bjet_n_flow_nominal", "b-jet Cutflow (Nominal);b-jet Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1411 ATH_CHECK( book(TH1D("Cutflow/fatjet_n_flow_nominal", "Large R. Jet Cutflow (Nominal);Large R. Jet Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1412 ATH_CHECK( book(TH1D("Cutflow/trkjet_n_flow_nominal", "Track Jet Cutflow (Nominal);Track Jet Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1413 ATH_CHECK( book(TH1D("Cutflow/tau_n_flow_nominal", "Tau Cutflow (Nominal);Tau Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1414
1415 ATH_CHECK( book(TH1D("MET/met_et", "MET (Nominal);MET (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) ); //MET (+Components)
1416 ATH_CHECK( book(TH1D("MET/met_phi", "MET_phi (Nominal);MET_phi (Nominal);N", 50, -5, 5) ) );
1417 ATH_CHECK( book(TH1D("MET/met_sumet", "MET_sumet (Nominal);MET_sumet (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1418 ATH_CHECK( book(TH1D("MET/met_et_tst", "MET [PVSoftTrk] (Nominal);MET [PVSoftTrk] (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1419 ATH_CHECK( book(TH1D("MET/met_et_el", "MET [RefEle] (Nominal);MET [RefEle] (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1420 ATH_CHECK( book(TH1D("MET/met_et_ph", "MET [RefGamma] (Nominal);MET [RefGamma] (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1421 ATH_CHECK( book(TH1D("MET/met_et_mu", "MET [Muons] (Nominal);MET [Muons] (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1422 ATH_CHECK( book(TH1D("MET/met_et_jet", "MET [RefJet] (Nominal);MET [RefJet] (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1423 //ATH_CHECK( book(TH1D("MET/met_et_tau", "MET [RefTau] (Nominal);MET [RefTau] (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1424 ATH_CHECK( book(TH1D("MET/met_significance", "MET Significance;MET Significance;N / [0.5]", 60, 0, 30) ) );
1425
1426 ATH_CHECK( book(TH1D("Syst/weight_event", "Event weights (Nom+Systematics) [MC*PRW];Event weights (Nom+Systematics) [MC*PRW];weight", getSize(m_syst_weights,"EventWeight"), 0, getSize(m_syst_weights,"EventWeight"))) ); //weights
1427 ATH_CHECK( book(TH1D("Syst/weight_electrons", "Electron total weights (Nom+Systematics);Electron total weights (Nom+Systematics);weight", getSize(m_syst_weights,"Electron") , 0, getSize(m_syst_weights,"Electron") )) ); //weights
1428 ATH_CHECK( book(TH1D("Syst/weight_muons", "Muon total weights (Nom+Systematics);Muon total weights (Nom+Systematics);weight", getSize(m_syst_weights,"Muon") , 0, getSize(m_syst_weights,"Muon") )) ); //weights
1429 ATH_CHECK( book(TH1D("Syst/weight_photons", "Photon total weights (Nom+Systematics);Photon total weights (Nom+Systematics);weight", getSize(m_syst_weights,"Photon") , 0, getSize(m_syst_weights,"Photon") )) ); //weights
1430 ATH_CHECK( book(TH1D("Syst/weight_taus", "Tau total weights (Nom+Systematics);Tau total weights (Nom+Systematics);weight", getSize(m_syst_weights,"Tau") , 0, getSize(m_syst_weights,"Tau") )) ); //weights
1431 ATH_CHECK( book(TH1D("Syst/weight_jets", "Jet total weights (Nom+Systematics);Jet total weights (Nom+Systematics);weight", getSize(m_syst_weights,"Jet") , 0, getSize(m_syst_weights,"Jet") )) ); //weights
1432 //ATH_CHECK( book(TH1D("Syst/weight_fatjets", "Large R. jet total weights (Nom+Systematics);Large R. jet total weights (Nom+Systematics);weight", getSize(m_syst_weights,"Jet") , 0, getSize(m_syst_weights,"Jet") )) ); //weights
1433 //ATH_CHECK( book(TH1D("Syst/weight_trkjets", "Track jet total weights (Nom+Systematics);Track jet total weights (Nom+Systematics);weight", getSize(m_syst_weights,"Jet") , 0, getSize(m_syst_weights,"Jet") )) ); //weights
1434 ATH_CHECK( book(TH1D("Syst/weight_btags", "Btagging total weights (Nom+Systematics);Btagging total weights (Nom+Systematics);weight", getSize(m_syst_weights,"BTag") , 0, getSize(m_syst_weights,"BTag") )) ); //weights
1435
1436 ATH_CHECK( book(TH2D("Syst/syst_event", "Systematics list [MC+PRW];;Syst name;", 1,0,1, static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"EventWeight"))),0,static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"EventWeight"))))) ); // syst names
1437 ATH_CHECK( book(TH2D("Syst/syst_electrons", "Systematics list [Electron];Syst name;", 1,0,1, static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"Electron"))), 0,static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"Electron"))))) ); // syst names
1438 ATH_CHECK( book(TH2D("Syst/syst_muons", "Systematics list [Muon];Syst name;", 1,0,1, static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"Muon"))), 0,static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"Muon"))))) ); // syst names
1439 ATH_CHECK( book(TH2D("Syst/syst_photons", "Systematics list [Photon];Syst name;", 1,0,1, static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"Photon"))), 0,static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"Photon"))))) ); // syst names
1440 ATH_CHECK( book(TH2D("Syst/syst_taus", "Systematics list [Tau];Syst name;", 1,0,1, static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"Tau"))), 0,static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"Tau"))))) ); // syst names
1441 ATH_CHECK( book(TH2D("Syst/syst_jets", "Systematics list [Jet];Syst name;", 1,0,1, static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"Jet"))), 0,static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"Jet"))))) ); // syst names
1442 ATH_CHECK( book(TH2D("Syst/syst_btags", "Systematics list [BTag];Syst name;", 1,0,1, static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"BTag"))), 0,static_cast<int>(std::ceil(1.2*getSize(m_syst_all,"BTag"))))) ); // syst names
1443
1444 return StatusCode::SUCCESS;
1445}
1446
1447//====================================================================================================
1448// collect systematics per object (all / weight-only)
1449void SUSYToolsAlg::groupSysts(void) {
1450 m_syst_all.clear();
1451 m_syst_weights.clear();
1452 const std::vector<std::string> Nominal{"Nominal"};
1453 for (const auto& sysInfo : m_sysInfoList) {
1454 std::string sys_name = sysInfo.systset.name();
1455 std::string sys_affects = ST::testAffectsObject(sysInfo.affectsType);
1456 ATH_MSG_DEBUG("Syst " << sys_name << " affects " << sys_affects);
1457
1458 bool sys_affects_weights = sysInfo.affectsWeights;
1459
1460 // collect all syst names
1461 // per affected object
1462 const auto & [pairPtr, inserted] = m_syst_all.try_emplace(sys_affects, Nominal);
1463 pairPtr->second.push_back(sys_name);
1464
1465 // weight related syst
1466 if (sys_affects_weights) {
1467 const auto & [pairPtr_w, inserted_w] = m_syst_weights.try_emplace(sys_affects, Nominal);
1468 pairPtr_w->second.push_back(sys_name);
1469 }
1470 }
1471}
1472
1473//====================================================================================================
1474// Fill object histograms consistently for different selection levels / object types
1475void SUSYToolsAlg::stdHistsForObj(xAOD::IParticle *obj, const std::string& objtype, const std::string& objlevel, std::map<std::string,std::string> config) {
1476 std::map<std::string,std::string> labels_dir = { {"el","Electron"}, {"ph","Photon"}, {"mu","Muon"}, {"jet","Jet"}, {"bjet","bJet"}, {"tau","Tau"}, {"fatjet","LargeRJet"}, {"trkjet","TrackJet"} };
1477 std::string dir = labels_dir[objtype]+"/";
1478
1479 ATH_MSG_VERBOSE("Filling histograms for " << obj->type() << ", " << objlevel);
1480 hist(dir+objtype+"_"+objlevel+"_pt")->Fill( obj->pt()/1000. );
1481 hist(dir+objtype+"_"+objlevel+"_eta")->Fill( obj->eta() );
1482 hist(dir+objtype+"_"+objlevel+"_phi")->Fill( obj->phi() );
1483 hist(dir+objtype+"_"+objlevel+"_e")->Fill( obj->e()/1000. );
1484 hist(dir+objtype+"_"+objlevel+"_mass")->Fill( obj->m()/1000. );
1485 //
1486 if (objtype=="el" || objtype=="mu" || objtype=="ph") {
1487 if (objtype=="el") {
1488 if (ST::acc_z0sinTheta.isAvailable(*obj))
1489 hist(dir+objtype+"_"+objlevel+"_z0")->Fill( ST::acc_z0sinTheta(*obj) );
1490 if (ST::acc_d0sig.isAvailable(*obj))
1491 hist(dir+objtype+"_"+objlevel+"_d0")->Fill( ST::acc_d0sig(*obj) );
1492 hist(dir+objtype+"_"+objlevel+"_isolCalo")->Fill( STAlg::acc_topoetcone20(*obj)/obj->pt() );
1493 hist(dir+objtype+"_"+objlevel+"_isolTrack")->Fill( STAlg::acc_ptvarcone30_TTVA_LooseCone(*obj)/obj->pt() );
1494 } else if (objtype=="mu") {
1495 if (ST::acc_z0sinTheta.isAvailable(*obj))
1496 hist(dir+objtype+"_"+objlevel+"_z0")->Fill( ST::acc_z0sinTheta(*obj) );
1497 if (ST::acc_d0sig.isAvailable(*obj))
1498 hist(dir+objtype+"_"+objlevel+"_d0")->Fill( ST::acc_d0sig(*obj) );
1499 hist(dir+objtype+"_"+objlevel+"_isolCalo")->Fill( STAlg::acc_topoetcone20(*obj)/obj->pt() );
1500 hist(dir+objtype+"_"+objlevel+"_isolTrack")->Fill( STAlg::acc_ptvarcone30_TTVA(*obj)/obj->pt() );
1501 } else if (objtype=="ph") {
1502 hist(dir+objtype+"_"+objlevel+"_isolCalo")->Fill( (STAlg::acc_topoetcone40(*obj)-2450)/obj->pt() );
1503 hist(dir+objtype+"_"+objlevel+"_isolTrack")->Fill( STAlg::acc_ptcone20(*obj)/obj->pt() );
1504 }
1505 }
1506 //
1507 if (objtype=="el" || objtype=="mu" || objtype=="ph" || objtype=="tau") {
1508 if ( objtype!="tau" ) {
1509 if (ST::acc_truthType.isAvailable(*obj) && ST::acc_truthOrigin.isAvailable(*obj)) {
1510 ATH_MSG_DEBUG(objtype << ", " << objlevel << ": " << ST::acc_truthType(*obj) << ", " << ST::acc_truthOrigin(*obj));
1511 hist(dir+objtype+"_"+objlevel+"_truthType")->Fill( ST::acc_truthType(*obj) );
1512 hist(dir+objtype+"_"+objlevel+"_truthOrigin")->Fill( ST::acc_truthOrigin(*obj) );
1513 }
1514 } else {
1515 if(!m_isPHYSLITE){
1516 bool istruthmatched = STAlg::acc_IsTruthMatched.isAvailable(*obj) ? (bool)STAlg::acc_IsTruthMatched(*obj) : false;
1517 int pid(0),ppid(0);
1518 if (istruthmatched && STAlg::acc_truthParticleLink.isAvailable(*obj)) {
1519 const auto *tp = *(STAlg::acc_truthParticleLink(*obj));
1520 if (tp) {
1521 pid = tp->pdgId();
1522 ppid = (tp->nParents()>0)?tp->parent(0)->pdgId():0;
1523 }
1524 }
1525 hist(dir+objtype+"_"+objlevel+"_pid")->Fill( pid );
1526 hist(dir+objtype+"_"+objlevel+"_parentpid")->Fill( ppid );
1527 }
1528 }
1529 }
1530
1531 if (objtype=="tau") {
1532 hist(dir+objtype+"_"+objlevel+"_nTracks")->Fill( dynamic_cast<xAOD::TauJet*>(obj)->nTracks() );
1533 hist(dir+objtype+"_"+objlevel+"_nTracksCharged")->Fill( dynamic_cast<xAOD::TauJet*>(obj)->nTracksCharged() );
1534 hist(dir+objtype+"_"+objlevel+"_nTracksIsolation")->Fill( dynamic_cast<xAOD::TauJet*>(obj)->nTracksIsolation() );
1535 if ( m_configDict.find("Tau.Tagger")!=m_configDict.end() && m_configDict["Tau.Tagger"].find("GNTau")!=std::string::npos) {
1536 hist(dir+objtype+"_"+objlevel+"_TauIDScoreSigTrans")->Fill( STAlg::acc_GNTauScoreSigTrans(*obj) );
1537 } else {
1538 hist(dir+objtype+"_"+objlevel+"_TauIDScoreSigTrans")->Fill( STAlg::acc_RNNJetScoreSigTrans(*obj) );
1539 }
1540 }
1541 //
1542
1543 if (objtype=="bjet") {
1544 hist(dir+objtype+"_"+objlevel+"_bweight")->Fill( ST::acc_btag_weight(*obj) );
1545 hist(dir+objtype+"_"+objlevel+"_bweightpb")->Fill( ST::acc_btag_pb(*obj) );
1546 hist(dir+objtype+"_"+objlevel+"_bweightpc")->Fill( ST::acc_btag_pc(*obj) );
1547 hist(dir+objtype+"_"+objlevel+"_bweightpu")->Fill( ST::acc_btag_pu(*obj) );
1548 hist(dir+objtype+"_"+objlevel+"_bweightptau")->Fill( ST::acc_btag_ptau(*obj) );
1549 }
1550
1551 //
1552 const static SG::ConstAccessor<char> acc_topValidPtRangeLow(config["ToptaggerKey"]+"_ValidPtRangeLow");
1553 const static SG::ConstAccessor<char> acc_topValidPtRangeHigh(config["ToptaggerKey"]+"_ValidPtRangeHigh");
1554 const static SG::ConstAccessor<char> acc_topValidEtaRange(config["ToptaggerKey"]+"_ValidEtaRange");
1555 const static SG::ConstAccessor<char> acc_wtagged(config["WtaggerKey"]+"_Tagged");
1556 const static SG::ConstAccessor<char> acc_ztagged(config["ZtaggerKey"]+"_Tagged");
1557 const static SG::ConstAccessor<char> acc_toptagged(config["ToptaggerKey"]+"_Tagged");
1558 const static SG::ConstAccessor<float> acc_topscore(config["ToptaggerKey"]+"_Score");
1559
1560 if (objtype=="fatjet" && acc_topValidPtRangeLow.isAvailable(*obj) && acc_topValidPtRangeHigh.isAvailable(*obj) && acc_topValidEtaRange.isAvailable(*obj)) {
1561 bool kin = acc_topValidPtRangeLow(*obj) && acc_topValidPtRangeHigh(*obj) && acc_topValidEtaRange(*obj);
1562 if (config.find("WtaggerKey")!=config.end()) hist(dir+objtype+"_"+objlevel+"_wtagged") ->Fill( acc_wtagged.isAvailable(*obj) ? kin && acc_wtagged(*obj) : -1);
1563 if (config.find("ZtaggerKey")!=config.end()) hist(dir+objtype+"_"+objlevel+"_ztagged") ->Fill( acc_ztagged.isAvailable(*obj) ? kin && acc_ztagged(*obj) : -1);
1564 if (config.find("ToptaggerKey")!=config.end()) hist(dir+objtype+"_"+objlevel+"_toptagged")->Fill( acc_toptagged.isAvailable(*obj) ? kin && acc_toptagged(*obj) : -1);
1565 if (config.find("ToptaggerKey")!=config.end()) hist(dir+objtype+"_"+objlevel+"_toptaggedscore")->Fill( acc_topscore.isAvailable(*obj) ? acc_topscore(*obj) : -1);
1566 }
1567}
1568
1569
1570//====================================================================================================
1571std::string findInReg(SG::AuxTypeRegistry& reg, xAOD::IParticle* obj, const std::string& searchkey) {
1572 std::regex re_tag(searchkey);
1573 std::smatch matches;
1574 for (auto x : obj->getAuxIDs()) {
1575 std::string xname = reg.getName(x);
1576 std::regex_search(xname, matches, re_tag);
1577 if (matches.size()>1 and !matches[1].str().empty()) {
1578 return matches[1].str();
1579 }
1580 }
1581 return "";
1582}
1583
1584//====================================================================================================
1585int getSize(std::map<std::string,std::vector<std::string>> &collection, const std::string& object) {
1586 if (collection.find(object) == collection.end()) return 1;
1587 return static_cast<int>(collection[object].size());
1588}
phi
Scalar phi() const
phi method
Definition AmgMatrixBasePlugin.h:67
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_VERBOSE
#define ATH_MSG_VERBOSE(x)
Definition AthMsgStreamMacros.h:28
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition AthMsgStreamMacros.h:32
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition AthMsgStreamMacros.h:29
Data
@ Data
Definition BaseObject.h:11
ANA_CHECK
#define ANA_CHECK(EXP)
check whether the given expression was successful
Definition Control/AthToolSupport/AsgMessaging/AsgMessaging/MessageCheck.h:324
CHECK
#define CHECK(...)
Evaluate an expression and check for errors.
Definition Control/AthenaKernel/AthenaKernel/errorcheck.h:422
EventStreamInfo.h
This file contains the class definition for the EventStreamInfo class.
taus
static Double_t taus
Definition LArPhysWaveHECTool.cxx:37
PathResolverFindCalibFile
std::string PathResolverFindCalibFile(const std::string &logical_file_name)
Definition PathResolver.cxx:325
SUSYObjDef_xAOD.h
getSize
int getSize(std::map< std::string, std::vector< std::string > > &collection, const std::string &object)
Definition SUSYToolsAlg.cxx:1585
findInReg
std::string findInReg(SG::AuxTypeRegistry &reg, xAOD::IParticle *obj, const std::string &searchkey)
Definition SUSYToolsAlg.cxx:1571
size
size_t size() const
Number of registered mappings.
x
#define x
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
h
Header file for AthHistogramAlgorithm.
AthHistogramming::book
StatusCode book(const TH1 &hist, const std::string &tDir="", const std::string &stream="")
Simplify the booking and registering (into THistSvc) of histograms.
Definition AthHistogramming.h:305
AthHistogramming::hist
TH1 * hist(const std::string &histName, const std::string &tDir="", const std::string &stream="")
Simplify the retrieval of registered histograms of any type.
Definition AthHistogramming.cxx:164
AthHistogramming::hist2d
TH2 * hist2d(const std::string &histName, const std::string &tDir="", const std::string &stream="")
Simplify the retrieval of registered 2-d histograms.
Definition AthHistogramming.h:371
CP::SystematicSet
Class to wrap a set of SystematicVariations.
Definition SystematicSet.h:31
DataVector::size
size_type size() const noexcept
Returns the number of elements in the collection.
EL::AnaAlgorithm::inputMetaStore
ConstMetaStorePtr_t inputMetaStore() const
Definition AnaAlgorithm.cxx:72
EL::AnaAlgorithm::AnaAlgorithm
AnaAlgorithm(const std::string &name, ISvcLocator *pSvcLocator)
constructor with parameters
Definition AnaAlgorithm.cxx:40
EventStreamInfo
This class provides the summary information stored for data written as a Event Stream.
Definition EventStreamInfo.h:28
EventStreamInfo::getProcessingTags
const std::set< std::string > & getProcessingTags() const
Definition EventStreamInfo.cxx:35
SG::AuxTypeRegistry
Handle mappings between names and auxid_t.
Definition AuxTypeRegistry.h:61
SG::AuxTypeRegistry::instance
static AuxTypeRegistry & instance()
Return the singleton registry instance.
Definition AuxTypeRegistry.cxx:640
SG::ConstAccessor
Helper class to provide constant type-safe access to aux data.
Definition ConstAccessor.h:55
SG::ConstAccessor::isAvailable
bool isAvailable(const ELT &e) const
Test to see if this variable exists in the store.
SUSYToolsAlg::m_vars
std::vector< std::string > m_vars
Definition SUSYToolsAlg.h:115
SUSYToolsAlg::m_heffs
std::map< std::string, TEfficiency * > m_heffs
Definition SUSYToolsAlg.h:82
SUSYToolsAlg::execute
StatusCode execute()
Definition SUSYToolsAlg.cxx:302
SUSYToolsAlg::m_grl
asg::AnaToolHandle< IGoodRunsListSelectionTool > m_grl
Definition SUSYToolsAlg.h:86
SUSYToolsAlg::m_mettst_nominal
xAOD::MissingETContainer * m_mettst_nominal
Definition SUSYToolsAlg.h:125
SUSYToolsAlg::m_clock2
TStopwatch m_clock2
Definition SUSYToolsAlg.h:107
SUSYToolsAlg::m_triggers
std::map< std::string, std::vector< std::string > > m_triggers
Definition SUSYToolsAlg.h:110
SUSYToolsAlg::m_mettst_syst
xAOD::MissingETContainer * m_mettst_syst
Definition SUSYToolsAlg.h:118
SUSYToolsAlg::m_syst_weights
std::map< std::string, std::vector< std::string > > m_syst_weights
Definition SUSYToolsAlg.h:102
SUSYToolsAlg::m_syst_all
std::map< std::string, std::vector< std::string > > m_syst_all
Definition SUSYToolsAlg.h:102
SUSYToolsAlg::m_metcst_nominal
xAOD::MissingETContainer * m_metcst_nominal
Definition SUSYToolsAlg.h:123
SUSYToolsAlg::m_slices
std::map< std::string, bool > m_slices
Definition SUSYToolsAlg.h:80
SUSYToolsAlg::m_metcst_syst
xAOD::MissingETContainer * m_metcst_syst
Definition SUSYToolsAlg.h:120
SUSYToolsAlg::m_obj_count
std::map< std::string, std::map< std::string, int > > m_obj_count
Definition SUSYToolsAlg.h:116
SUSYToolsAlg::m_kernel
std::string m_kernel
Definition SUSYToolsAlg.h:93
SUSYToolsAlg::m_TrkJetCollection
std::string m_TrkJetCollection
Definition SUSYToolsAlg.h:96
SUSYToolsAlg::m_objects
std::vector< std::string > m_objects
Definition SUSYToolsAlg.h:113
SUSYToolsAlg::m_mcCampaign
std::string m_mcCampaign
Definition SUSYToolsAlg.h:79
SUSYToolsAlg::m_mettst_nominal_aux
xAOD::MissingETAuxContainer * m_mettst_nominal_aux
Definition SUSYToolsAlg.h:126
SUSYToolsAlg::m_hists
std::map< std::string, TH1 * > m_hists
Definition SUSYToolsAlg.h:81
SUSYToolsAlg::m_SUSYTools
asg::AnaToolHandle< ST::ISUSYObjDef_xAODTool > m_SUSYTools
Definition SUSYToolsAlg.h:85
SUSYToolsAlg::m_metcst_nominal_aux
xAOD::MissingETAuxContainer * m_metcst_nominal_aux
Definition SUSYToolsAlg.h:124
SUSYToolsAlg::m_configFile
std::string m_configFile
Definition SUSYToolsAlg.h:94
SUSYToolsAlg::groupSysts
void groupSysts(void)
Definition SUSYToolsAlg.cxx:1449
SUSYToolsAlg::bookHistograms
StatusCode bookHistograms(void)
Definition SUSYToolsAlg.cxx:1324
SUSYToolsAlg::m_tauTruthMatchingTool
asg::AnaToolHandle< TauAnalysisTools::ITauTruthMatchingTool > m_tauTruthMatchingTool
Definition SUSYToolsAlg.h:87
SUSYToolsAlg::stdHistsForObj
void stdHistsForObj(xAOD::IParticle *obj, const std::string &objtype, const std::string &objlevel, std::map< std::string, std::string > config=std::map< std::string, std::string >())
Definition SUSYToolsAlg.cxx:1475
SUSYToolsAlg::m_mettst_syst_aux
xAOD::MissingETAuxContainer * m_mettst_syst_aux
Definition SUSYToolsAlg.h:119
SUSYToolsAlg::m_clock1
TStopwatch m_clock1
Definition SUSYToolsAlg.h:106
SUSYToolsAlg::m_levels
std::vector< std::string > m_levels
Definition SUSYToolsAlg.h:114
SUSYToolsAlg::m_metcst_syst_aux
xAOD::MissingETAuxContainer * m_metcst_syst_aux
Definition SUSYToolsAlg.h:121
SUSYToolsAlg::m_clock0
TStopwatch m_clock0
Definition SUSYToolsAlg.h:105
SUSYToolsAlg::initialize
StatusCode initialize()
Definition SUSYToolsAlg.cxx:103
SUSYToolsAlg::m_sysInfoList
std::vector< ST::SystInfo > m_sysInfoList
Definition SUSYToolsAlg.h:101
SUSYToolsAlg::m_Nevts
unsigned int m_Nevts
Definition SUSYToolsAlg.h:89
SUSYToolsAlg::m_configDict
std::map< std::string, std::string > m_configDict
Definition SUSYToolsAlg.h:83
SUSYToolsAlg::m_FatJetCollection
std::string m_FatJetCollection
Definition SUSYToolsAlg.h:95
SUSYToolsAlg::m_GRLFiles
std::vector< std::string > m_GRLFiles
Definition SUSYToolsAlg.h:78
SUSYToolsAlg::finalize
StatusCode finalize()
Definition SUSYToolsAlg.cxx:256
xAOD::EventInfo_v1::lumiBlock
uint32_t lumiBlock() const
The current event's luminosity block number.
xAOD::EventInfo_v1::runNumber
uint32_t runNumber() const
The current event's run number.
xAOD::EventInfo_v1::eventNumber
uint64_t eventNumber() const
The current event's event number.
xAOD::EventInfo_v1::mcEventWeight
float mcEventWeight(size_t i=0) const
The weight of one specific MC event used in the simulation.
Definition EventInfo_v1.cxx:203
xAOD::IParticle
Class providing the definition of the 4-vector interface.
Definition Event/xAOD/xAODBase/xAODBase/IParticle.h:41
xAOD::ShallowAuxContainer
Class creating a shallow copy of an existing auxiliary container.
Definition ShallowAuxContainer.h:50
contains
bool contains(const std::string &s, const std::string &regx)
does a string contain the substring
Definition hcg.cxx:116
find
std::string find(const std::string &s)
return a remapped string
Definition hcg.cxx:140
Cut::SCut
static const TString SCut[]
Definition SUSYToolsAlg.cxx:68
Cut::NSel
const unsigned int NSel
Definition SUSYToolsAlg.cxx:66
Cut::trigmatch
@ trigmatch
Definition SUSYToolsAlg.cxx:67
EL
This module defines the arguments passed from the BATCH driver to the BATCH worker.
Definition AsgComponentFactories.h:16
STAlg::acc_RNNJetScoreSigTrans
static const SG::ConstAccessor< float > acc_RNNJetScoreSigTrans("RNNJetScoreSigTrans")
STAlg::acc_ptvarcone30_TTVA
static const SG::ConstAccessor< float > acc_ptvarcone30_TTVA("ptvarcone30_Nonprompt_All_MaxWeightTTVA_pt1000")
STAlg::acc_IsTruthMatched
static const SG::ConstAccessor< char > acc_IsTruthMatched("IsTruthMatched")
STAlg::acc_GNTauScoreSigTrans
static const SG::ConstAccessor< float > acc_GNTauScoreSigTrans("GNTauScoreSigTrans_v0prune")
STAlg::acc_topoetcone20
static const SG::ConstAccessor< float > acc_topoetcone20("topoetcone20")
STAlg::acc_ptvarcone30_TTVA_LooseCone
static const SG::ConstAccessor< float > acc_ptvarcone30_TTVA_LooseCone("ptvarcone30_Nonprompt_All_MaxWeightTTVALooseCone_pt1000")
STAlg::acc_ptcone20
static const SG::ConstAccessor< float > acc_ptcone20("ptcone20")
STAlg::acc_topoetcone40
static const SG::ConstAccessor< float > acc_topoetcone40("topoetcone40")
STAlg::acc_truthParticleLink
static const SG::ConstAccessor< ElementLink< xAOD::TruthParticleContainer > > acc_truthParticleLink("truthParticleLink")
ST::acc_d0sig
static const SG::ConstAccessor< float > acc_d0sig("d0sig")
ST::acc_btag_weight
static const SG::ConstAccessor< double > acc_btag_weight("btag_weight")
ST::acc_effscalefact
static const SG::ConstAccessor< double > acc_effscalefact("effscalefact")
ST::acc_bad
static const SG::ConstAccessor< char > acc_bad("bad")
ST::acc_truthType
static const SG::ConstAccessor< int > acc_truthType("truthType")
ST::acc_truthOrigin
static const SG::ConstAccessor< int > acc_truthOrigin("truthOrigin")
ST::acc_z0sinTheta
static const SG::ConstAccessor< float > acc_z0sinTheta("z0sinTheta")
ST::acc_signal
static const SG::ConstAccessor< char > acc_signal("signal")
ST::acc_baseline
static const SG::ConstAccessor< char > acc_baseline("baseline")
ST::testAffectsObject
static bool testAffectsObject(xAOD::Type::ObjectType type, unsigned int test)
Definition ISUSYObjDef_xAODTool.h:134
ST::acc_bjet
static const SG::ConstAccessor< char > acc_bjet("bjet")
ST::acc_cosmic
static const SG::ConstAccessor< char > acc_cosmic("cosmic")
ST::acc_btag_pu
static const SG::ConstAccessor< float > acc_btag_pu("btag_pu")
ST::acc_passOR
static const SG::ConstAccessor< char > acc_passOR("passOR")
ST::acc_btag_pc
static const SG::ConstAccessor< float > acc_btag_pc("btag_pc")
ST::acc_btag_ptau
static const SG::ConstAccessor< float > acc_btag_ptau("btag_ptau")
ST::acc_btag_pb
static const SG::ConstAccessor< float > acc_btag_pb("btag_pb")
xAODType::Jet
@ Jet
The object is a jet.
Definition ObjectType.h:40
xAODType::Photon
@ Photon
The object is a photon.
Definition ObjectType.h:47
xAODType::Muon
@ Muon
The object is a muon.
Definition ObjectType.h:48
xAODType::BTag
@ BTag
The object is a b-tagging object.
Definition ObjectType.h:60
xAODType::Electron
@ Electron
The object is an electron.
Definition ObjectType.h:46
xAODType::Tau
@ Tau
The object is a tau (jet).
Definition ObjectType.h:49
xAOD::PhotonContainer
PhotonContainer_v1 PhotonContainer
Definition of the current "photon container version".
Definition Event/xAOD/xAODEgamma/xAODEgamma/PhotonContainer.h:17
xAOD::ElectronContainer
ElectronContainer_v1 ElectronContainer
Definition of the current "electron container version".
Definition Event/xAOD/xAODEgamma/xAODEgamma/ElectronContainer.h:17
xAOD::EventInfo
EventInfo_v1 EventInfo
Definition of the latest event info version.
Definition IEventInfoCnvTool.h:16
xAOD::Unknown
@ Unknown
Track fitter not defined.
Definition TrackingPrimitives.h:42
xAOD::TauJet
TauJet_v3 TauJet
Definition of the current "tau version".
Definition Event/xAOD/xAODTau/xAODTau/TauJet.h:17
xAOD::MissingETAuxContainer
MissingETAuxContainer_v1 MissingETAuxContainer
Definition MissingETAuxContainer.h:16
xAOD::JetContainer
JetContainer_v1 JetContainer
Definition of the current "jet container version".
Definition JetContainer.h:17
xAOD::TauJetContainer
TauJetContainer_v3 TauJetContainer
Definition of the current "taujet container version".
Definition Event/xAOD/xAODTau/xAODTau/TauJetContainer.h:17
xAOD::MuonContainer
MuonContainer_v1 MuonContainer
Definition of the current "Muon container version".
Definition Event/xAOD/xAODMuon/xAODMuon/MuonContainer.h:14
xAOD::CutBookkeeperContainer
CutBookkeeperContainer_v1 CutBookkeeperContainer
Define the latest version of the CutBookkeeperContainer class.
Definition CutBookkeeperContainer.h:18
ST::SystInfo::affectsType
unsigned int affectsType
Definition ISUSYObjDef_xAODTool.h:68
Generated on for ATLAS Offline Software by doxygen 1.16.1