ATLAS Offline Software
Loading...
Searching...
No Matches
SUSYToolsAlg.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// 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 }
733
734 //----- Large R. Jets
735
736 TH1* fatjet_n_flow_nominal = hist("Cutflow/fatjet_n_flow_nominal");
737
738 if( m_slices["fjet"] && doFatJets ) {
739 ATH_MSG_DEBUG("Processing fatjets");
740 if (m_configDict.find("Jet.WtaggerConfig")!=m_configDict.end() && m_configDict.find("WtaggerKey")==m_configDict.end() && fatjets_nominal->size()>0) {
741 SG::AuxTypeRegistry& reg = SG::AuxTypeRegistry::instance();
742 std::string taggerKey;
743 taggerKey = findInReg(reg, (*fatjets_nominal)[0], "(SmoothW.*?)_.*");
744 if (!taggerKey.empty()) { m_configDict["WtaggerKey"] = taggerKey; }
745 }
746 if (m_configDict.find("Jet.ZtaggerConfig")!=m_configDict.end() && m_configDict.find("ZtaggerKey")==m_configDict.end() && fatjets_nominal->size()>0) {
747 SG::AuxTypeRegistry& reg = SG::AuxTypeRegistry::instance();
748 std::string taggerKey;
749 taggerKey = findInReg(reg, (*fatjets_nominal)[0], "(SmoothZ.*?)_.*");
750 if (!taggerKey.empty()) { m_configDict["ZtaggerKey"] = taggerKey; }
751 }
752 if (m_configDict.find("Jet.ToptaggerConfig")!=m_configDict.end() && m_configDict.find("ToptaggerKey")==m_configDict.end() && fatjets_nominal->size()>0) {
753 SG::AuxTypeRegistry& reg = SG::AuxTypeRegistry::instance();
754 std::string taggerKey;
755 taggerKey = findInReg(reg, (*fatjets_nominal)[0], "(.*Tag.*Top.*?)_.*");
756 if (!taggerKey.empty()) { m_configDict["ToptaggerKey"] = taggerKey; }
757 }
758 ATH_MSG_DEBUG("Tagger configuration retrieved.");
759 ATH_MSG_DEBUG("WTagger: " << m_configDict["WtaggerKey"]);
760 ATH_MSG_DEBUG("ZTagger: " << m_configDict["ZtaggerKey"]);
761 ATH_MSG_DEBUG("TopTagger: " << m_configDict["ToptaggerKey"]);
762
763 m_obj_count["fatjet"]["nom"] = fatjets_nominal->size();
764 for(auto fatjet : *fatjets_nominal) {
765 fatjet_n_flow_nominal->Fill(Cut::all);
766 if ( ST::acc_baseline(*fatjet) == 1 ){
767 fatjet_n_flow_nominal->Fill(Cut::baseline);
768 stdHistsForObj(fatjet,"fatjet","bsl",m_configDict);
769 m_obj_count["fatjet"]["bsl"] += 1;
770 if ( ST::acc_passOR(*fatjet) == 1 ){
771 fatjet_n_flow_nominal->Fill(Cut::passOR);
772 if ( ST::acc_signal(*fatjet) == 1 ){
773 fatjet_n_flow_nominal->Fill(Cut::signal);
774 fatjet_n_flow_nominal->Fill(Cut::trigmatch); //no trig matching for jets
775 stdHistsForObj(fatjet,"fatjet","sig",m_configDict);
776 m_obj_count["fatjet"]["sig"] += 1;
777 }
778 }
779 }
780 stdHistsForObj(fatjet,"fatjet","nom",m_configDict);
781 }
782 hist("LargeRJet/fatjet_nom_N")->Fill(m_obj_count["fatjet"]["nom"]);
783 hist("LargeRJet/fatjet_bsl_N")->Fill(m_obj_count["fatjet"]["bsl"]);
784 hist("LargeRJet/fatjet_sig_N")->Fill(m_obj_count["fatjet"]["sig"]);
785 }
786
787 //----- Track Jets
788
789 TH1* trkjet_n_flow_nominal = hist("Cutflow/trkjet_n_flow_nominal");
790
791 if( m_slices["tjet"] && doTrkJets ) {
792 m_obj_count["trkjet"]["nom"] = trkjets_nominal->size();
793 for(auto trkjet : *trkjets_nominal) {
794 trkjet_n_flow_nominal->Fill(Cut::all);
795 if ( ST::acc_baseline(*trkjet) == 1 ){
796 trkjet_n_flow_nominal->Fill(Cut::baseline);
797 stdHistsForObj(trkjet,"trkjet","bsl");
798 m_obj_count["trkjet"]["bsl"] += 1;
799 if ( ST::acc_passOR(*trkjet) == 1 ){
800 trkjet_n_flow_nominal->Fill(Cut::passOR);
801 if ( ST::acc_signal(*trkjet) == 1 ){
802 trkjet_n_flow_nominal->Fill(Cut::signal);
803 trkjet_n_flow_nominal->Fill(Cut::trigmatch); //no trig matching for jets
804 stdHistsForObj(trkjet,"trkjet","sig");
805 m_obj_count["trkjet"]["sig"] += 1;
806 }
807 }
808 }
809 stdHistsForObj(trkjet,"trkjet","nom");
810 }
811 hist("TrackJet/trkjet_nom_N")->Fill(m_obj_count["trkjet"]["nom"]);
812 hist("TrackJet/trkjet_bsl_N")->Fill(m_obj_count["trkjet"]["bsl"]);
813 hist("TrackJet/trkjet_sig_N")->Fill(m_obj_count["trkjet"]["sig"]);
814 }
815
816 //----- Taus
817
818 TH1* tau_n_flow_nominal = hist("Cutflow/tau_n_flow_nominal");
819
820 if (m_slices["tau"]) {
821 m_obj_count["tau"]["nom"] = taus_nominal->size();
822 for(auto tau : *taus_nominal) {
823 tau_n_flow_nominal->Fill(Cut::all);
824 if ( ST::acc_baseline(*tau) == 1 ){
825 tau_n_flow_nominal->Fill(Cut::baseline);
826 stdHistsForObj(tau,"tau","bsl");
827 m_obj_count["tau"]["bsl"] += 1;
828 if ( ST::acc_passOR(*tau) == 1 ){
829 tau_n_flow_nominal->Fill(Cut::passOR);
830
831 if ( ST::acc_signal(*tau) == 1 ){
832 tau_n_flow_nominal->Fill(Cut::signal);
833 tau_n_flow_nominal->Fill(Cut::trigmatch); //no trig matching for taus
834 stdHistsForObj(tau,"tau","sig");
835 m_obj_count["tau"]["sig"] += 1;
836 }
837 }
838 }
839 stdHistsForObj(tau,"tau","nom");
840 }
841 hist("Tau/tau_nom_N")->Fill(m_obj_count["tau"]["nom"]);
842 hist("Tau/tau_bsl_N")->Fill(m_obj_count["tau"]["bsl"]);
843 hist("Tau/tau_sig_N")->Fill(m_obj_count["tau"]["sig"]);
844 }
845
846 // Set up the event weights
847 // Base should include all weights that do not depend on individual objects
848 double base_event_weight(1.);
849 if (!isData) {
850 base_event_weight *= evtInfo->mcEventWeight();
851
852 if (m_slices["met"]) {
853 const xAOD::MissingETContainer* met_truth(0);
854 ATH_CHECK( evtStore()->retrieve(met_truth, "MET_Truth") );
855
856 ATH_MSG_DEBUG("Truth MET etx=" << (*met_truth)["NonInt"]->mpx()
857 << ", ety=" << (*met_truth)["NonInt"]->mpy()
858 << ", et=" << (*met_truth)["NonInt"]->met()
859 << ", sumet=" << (*met_truth)["NonInt"]->sumet());
860 ATH_MSG_DEBUG("CST residuals: detx=" << (*m_metcst_nominal)["Final"]->mpx() - (*met_truth)["NonInt"]->mpx()
861 << ", dety=" << (*m_metcst_nominal)["Final"]->mpy() - (*met_truth)["NonInt"]->mpy());
862 ATH_MSG_DEBUG("TST residuals: detx=" << (*m_mettst_nominal)["Final"]->mpx() - (*met_truth)["NonInt"]->mpx()
863 << ", dety=" << (*m_mettst_nominal)["Final"]->mpy() - (*met_truth)["NonInt"]->mpy());
864 }
865 }
866
867
868 // Additionally define a nominal weight for each object type
869 double electrons_weight_nominal(1.);
870 double photons_weight_nominal(1.);
871 double muons_weight_nominal(1.);
872 double jets_weight_nominal(1.);
873 //double fatjets_weight_nominal(1.);
874 //double trkjets_weight_nominal(1.);
875 double btag_weight_nominal(1.);
876 double taus_weight_nominal(1.);
877 double event_weight_nominal(1.);
878 double prw_weight_nominal(1.);
879
880 //--- Weights
881 TH1* weight_event = hist("Syst/weight_event");
882 TH1* weight_electrons = hist("Syst/weight_electrons");
883 TH1* weight_muons = hist("Syst/weight_muons");
884 TH1* weight_photons = hist("Syst/weight_photons");
885 TH1* weight_taus = hist("Syst/weight_taus");
886 TH1* weight_jets = hist("Syst/weight_jets");
887 TH1* weight_btags = hist("Syst/weight_btags");
888
889 bool isNominal(true);
890 for (const auto& sysInfo : m_sysInfoList) {
891 CP::SystematicSet sys = sysInfo.systset;
892 if (m_SUSYTools->applySystematicVariation(sys) != StatusCode::SUCCESS) {
893 ATH_MSG_ERROR( "Cannot configure SUSYTools for systematic var. %s" << sys.name() );
894 } else {
895 ATH_MSG_VERBOSE( "Variation \"" << sys.name() << "\" configured..." );
896 }
897 if (isNominal && (sysInfo.affectsKinematics || sysInfo.affectsWeights)) {
898 ATH_MSG_ERROR("Expected nominal variation!");
899 return StatusCode::FAILURE;
900 }
901 if (sysInfo.affectsKinematics || sysInfo.affectsWeights) isNominal = false;
902
903 double event_weight = base_event_weight;
904 bool syst_affectsEventWeight = (sysInfo.affectsType == ST::SystObjType::EventWeight);
905
906 // If nominal, compute the nominal weight, otherwise recompute the weight
907
908 double prw_weight=1.;
909 if(isNominal){
910 if(!isData) prw_weight_nominal = m_SUSYTools->GetPileupWeight();
911 event_weight *= prw_weight_nominal;
912 weight_event->SetBinContent(1, weight_event->GetBinContent(1)+event_weight);
913 }
914 else if (!syst_affectsEventWeight) {
915 event_weight *= prw_weight_nominal;
916 }
917 else {
918 if(!isData){
919 ATH_CHECK(m_SUSYTools->ApplyPRWTool());
920 prw_weight = m_SUSYTools->GetPileupWeight();
921 }
922 event_weight *= prw_weight;
923 size_t iwbin = find(m_syst_weights["EventWeight"].begin(), m_syst_weights["EventWeight"].end(), sys.name()) - m_syst_weights["EventWeight"].begin();
924 if(iwbin < m_syst_weights["EventWeight"].size()) { weight_event->SetBinContent(iwbin+1, weight_event->GetBinContent(iwbin+1)+event_weight); }
925 }
926
927
928 // Define the generic collection pointers
929 xAOD::ElectronContainer* electrons(electrons_nominal);
930 xAOD::PhotonContainer* photons(photons_nominal);
931 xAOD::MuonContainer* muons(muons_nominal);
932 xAOD::JetContainer* jets(jets_nominal);
933 xAOD::JetContainer* fatjets(fatjets_nominal);
934 xAOD::JetContainer* trkjets(trkjets_nominal);
935 xAOD::TauJetContainer* taus(taus_nominal);
936 xAOD::MissingETContainer *mettst(m_mettst_nominal), *metcst(m_metcst_nominal);
937 xAOD::MissingETAuxContainer *mettst_aux(m_mettst_nominal_aux), *metcst_aux(m_metcst_nominal_aux);
938
939 bool syst_affectsElectrons = ST::testAffectsObject(xAOD::Type::Electron, sysInfo.affectsType);
940 bool syst_affectsMuons = ST::testAffectsObject(xAOD::Type::Muon, sysInfo.affectsType);
941 bool syst_affectsPhotons = ST::testAffectsObject(xAOD::Type::Photon, sysInfo.affectsType);
942 bool syst_affectsTaus = ST::testAffectsObject(xAOD::Type::Tau, sysInfo.affectsType);
943 bool syst_affectsJets = ST::testAffectsObject(xAOD::Type::Jet, sysInfo.affectsType);
944 bool syst_affectsBTag = ST::testAffectsObject(xAOD::Type::BTag, sysInfo.affectsType);
945 bool syst_affectsLRT = ST::testAffectsObject(sysInfo.affectsType) == "LRT_Objects";
946
947 // Apply LRT uncertainty.
948 // Needs to be before retrieving leptons as they could be removed if the track is vetoed
949 if (syst_affectsLRT)
950 {
951 ATH_CHECK( m_SUSYTools->ApplyLRTUncertainty() );
952 }
953
954 // If necessary (kinematics affected), make a shallow copy with the variation applied
955 // otherwise set the collection pointers to the nominal
956 if (sysInfo.affectsKinematics) {
957 if (m_slices["el"] && syst_affectsElectrons) {
958 ATH_MSG_DEBUG("Get systematics-varied electrons");
959 xAOD::ElectronContainer* electrons_syst(0);
960 xAOD::ShallowAuxContainer* electrons_syst_aux(0);
961 ATH_CHECK( m_SUSYTools->GetElectrons(electrons_syst, electrons_syst_aux, true , m_isPHYSLITE?"AnalysisElectrons":"Electrons") );
962 electrons = electrons_syst;
963 }
964
965 if (m_slices["pho"] && syst_affectsPhotons) {
966 ATH_MSG_DEBUG("Get systematics-varied photons");
967 xAOD::PhotonContainer* photons_syst(0);
968 xAOD::ShallowAuxContainer* photons_syst_aux(0);
969 ATH_CHECK( m_SUSYTools->GetPhotons(photons_syst, photons_syst_aux, true, m_isPHYSLITE?"AnalysisPhotons":"Photons") );
970 photons = photons_syst;
971 }
972
973 if (m_slices["mu"] && syst_affectsMuons) {
974 ATH_MSG_DEBUG("Get systematics-varied muons");
975 xAOD::MuonContainer* muons_syst(0);
976 xAOD::ShallowAuxContainer* muons_syst_aux(0);
977 ATH_CHECK( m_SUSYTools->GetMuons(muons_syst, muons_syst_aux, true, m_isPHYSLITE?"AnalysisMuons":"Muons") );
978 muons = muons_syst;
979 }
980
981 if (syst_affectsJets) {
982 if (m_slices["jet"]) {
983 ATH_MSG_DEBUG("Get systematics-varied jets");
984 xAOD::JetContainer* jets_syst(0);
985 xAOD::ShallowAuxContainer* jets_syst_aux(0);
986 ATH_CHECK( m_SUSYTools->GetJetsSyst(*jets_nominal, jets_syst, jets_syst_aux, true, m_isPHYSLITE?"AnalysisJets":"") );
987 jets = jets_syst;
988 }
989 if (m_slices["fjet"] and doFatJets) {
990 ATH_MSG_DEBUG("Get systematics-varied fatjets");
991 xAOD::JetContainer* fatjets_syst(0);
992 xAOD::ShallowAuxContainer* fatjets_syst_aux(0);
993 ATH_CHECK( m_SUSYTools->GetFatJets(fatjets_syst, fatjets_syst_aux, true, m_FatJetCollection, true) );
994 fatjets = fatjets_syst;
995 }
996 if (m_slices["trkjet"]) {
997 ATH_MSG_DEBUG("Get systematics-varied trkjets");
998 xAOD::JetContainer* trkjets_syst(0);
999 xAOD::ShallowAuxContainer* trkjets_syst_aux(0);
1000 ATH_CHECK( m_SUSYTools->GetJetsSyst(*trkjets_nominal, trkjets_syst, trkjets_syst_aux) );
1001 trkjets = trkjets_syst;
1002 }
1003 }
1004
1005 if (m_slices["tau"] && syst_affectsTaus) {
1006 ATH_MSG_DEBUG("Get systematics-varied taus");
1007 xAOD::TauJetContainer* taus_syst(0);
1008 xAOD::ShallowAuxContainer* taus_syst_aux(0);
1009 ATH_CHECK( m_SUSYTools->GetTaus(taus_syst, taus_syst_aux,true, m_isPHYSLITE?"AnalysisTauJets":"TauJets") );
1010 taus = taus_syst;
1011 }
1012
1013 ATH_MSG_DEBUG("Get systematics-varied MET");
1014 m_mettst_syst = new xAOD::MissingETContainer;
1015 m_mettst_syst_aux = new xAOD::MissingETAuxContainer;
1016 m_metcst_syst = new xAOD::MissingETContainer;
1017 m_metcst_syst_aux = new xAOD::MissingETAuxContainer;
1018 m_mettst_syst->setStore(m_mettst_syst_aux);
1019 m_metcst_syst->setStore(m_metcst_syst_aux);
1020 m_mettst_syst->reserve(10);
1021 m_metcst_syst->reserve(10);
1022 //
1023 if (m_slices["met"]) {
1024 ATH_CHECK( m_SUSYTools->GetMET(*m_metcst_syst, jets, electrons, muons, photons, 0) );
1025 ATH_CHECK( m_SUSYTools->GetMET(*m_mettst_syst, jets, electrons, muons, photons, 0, true, true) );
1026 }
1027 //
1028 mettst = m_mettst_syst;
1029 mettst_aux = m_mettst_syst_aux;
1030 metcst = m_metcst_syst;
1031 metcst_aux = m_metcst_syst_aux;
1032 }
1033
1034 //--- Overlap Removal
1035 ATH_CHECK( m_SUSYTools->OverlapRemoval(electrons, muons, jets, photons, taus) );
1036
1037 //--- Electrons
1038 if (m_slices["ele"]) {
1039 ATH_MSG_DEBUG("Working on electrons");
1040 float electrons_weight(1.);
1041 for ( const auto el : *electrons ) {
1042 if( !isData ){
1043 if (isNominal || syst_affectsElectrons) {
1044 if ((ST::acc_signal(*el) == 1) && (isNominal || sysInfo.affectsWeights)) {
1045 //electrons_weight *= m_SUSYTools->GetSignalElecSF( *el ); // (*el, true, true, false, true) to switch off trigger SF
1046 electrons_weight *= m_SUSYTools->GetSignalElecSF(*el, true, true, false, true);
1047 }
1048 }
1049 }
1050 ATH_MSG_VERBOSE( " Electron passing baseline selection? " << static_cast<int>( ST::acc_baseline(*el)));
1051 ATH_MSG_VERBOSE( " Electron passing signal selection? " << static_cast<int>( ST::acc_signal(*el)));
1052 if (ST::acc_signal(*el) == 1 && ST::acc_effscalefact.isAvailable(*el))
1053 ATH_MSG_VERBOSE( " Electron weight " << ST::acc_effscalefact(*el) );
1054
1055 }
1056 if (isNominal) {
1057 electrons_weight_nominal = electrons_weight;
1058 weight_electrons->SetBinContent(1, weight_electrons->GetBinContent(1)+electrons_weight);
1059 }
1060 else if (!syst_affectsElectrons) {
1061 electrons_weight = electrons_weight_nominal;
1062 }
1063 else if ( sysInfo.affectsWeights ){
1064 size_t iwbin = find(m_syst_weights["Electron"].begin(), m_syst_weights["Electron"].end(), sys.name()) - m_syst_weights["Electron"].begin();
1065 if(iwbin < m_syst_weights["Electron"].size()) { weight_electrons->SetBinContent(iwbin+1, weight_electrons->GetBinContent(iwbin+1)+electrons_weight); }
1066 }
1067
1068 event_weight *= electrons_weight;
1069 ATH_MSG_DEBUG("Combined electron scale factor: " << electrons_weight);
1070 }
1071
1072 //--- Photons
1073 if (m_slices["pho"]) {
1074 ATH_MSG_DEBUG("Working on photons");
1075 float photons_weight(1.);
1076 for ( const auto ph : *photons ) {
1077 if( !isData ){
1078 if (isNominal || syst_affectsPhotons) {
1079 if ((ST::acc_signal(*ph) == 1) && (isNominal || sysInfo.affectsWeights)) {
1080 photons_weight *= m_SUSYTools->GetSignalPhotonSF( *ph );
1081 }
1082 }
1083 }
1084 ATH_MSG_VERBOSE( " Photon passing baseline selection? " << static_cast<int>(ST::acc_baseline(*ph)));
1085 ATH_MSG_VERBOSE( " Photon passing signal selection? " << static_cast<int>(ST::acc_signal(*ph)));
1086 if (ST::acc_signal(*ph) == 1 && ST::acc_effscalefact.isAvailable(*ph))
1087 ATH_MSG_VERBOSE( " Photon weight " << ST::acc_effscalefact(*ph) );
1088 }
1089 if (isNominal) {
1090 photons_weight_nominal = photons_weight;
1091 weight_photons->SetBinContent(1, weight_photons->GetBinContent(1)+photons_weight);
1092 }
1093 else if (!syst_affectsPhotons) {
1094 photons_weight = photons_weight_nominal;
1095 }
1096 else if ( sysInfo.affectsWeights ){
1097 size_t iwbin = find(m_syst_weights["Photon"].begin(), m_syst_weights["Photon"].end(), sys.name()) - m_syst_weights["Photon"].begin();
1098 if(iwbin < m_syst_weights["Photon"].size()) { weight_photons->SetBinContent(iwbin+1, weight_photons->GetBinContent(iwbin+1)+photons_weight); }
1099 }
1100
1101 event_weight *= photons_weight;
1102 ATH_MSG_DEBUG("Combined photon scale factor: " << photons_weight);
1103 }
1104
1105 //--- Muons
1106 if (m_slices["mu"]) {
1107 ATH_MSG_DEBUG("Working on muons");
1108 float muons_weight(1.);
1109 for ( const auto mu : *muons ) {
1110 if( !isData ){
1111 if (isNominal || syst_affectsMuons) {
1112 if ((ST::acc_signal(*mu) == 1) && (isNominal || sysInfo.affectsWeights)) {
1113 muons_weight *= m_SUSYTools->GetSignalMuonSF(*mu);;
1114 }
1115 }
1116 }
1117 ATH_MSG_VERBOSE( " Muon passing baseline selection? " << static_cast<int>(ST::acc_baseline(*mu)));
1118 ATH_MSG_VERBOSE( " Muon passing signal selection? " << static_cast<int>(ST::acc_signal(*mu)));
1119 ATH_MSG_VERBOSE( " Muon is a cosmic ray? " << static_cast<int>(ST::acc_cosmic(*mu)));
1120 if (ST::acc_signal(*mu) == 1 && ST::acc_effscalefact.isAvailable(*mu))
1121 ATH_MSG_VERBOSE( " Muon weight " << ST::acc_effscalefact(*mu) );
1122 }
1123 if (isNominal) {
1124 muons_weight_nominal = muons_weight;
1125 weight_muons->SetBinContent(1, weight_muons->GetBinContent(1)+muons_weight);
1126 }
1127 else if (!syst_affectsMuons) {
1128 muons_weight = muons_weight_nominal;
1129 }
1130 else if ( sysInfo.affectsWeights ){
1131 size_t iwbin = find(m_syst_weights["Muon"].begin(), m_syst_weights["Muon"].end(), sys.name()) - m_syst_weights["Muon"].begin();
1132 if(iwbin < m_syst_weights["Muon"].size()) { weight_muons->SetBinContent(iwbin+1, weight_muons->GetBinContent(iwbin+1)+muons_weight); }
1133 }
1134
1135 event_weight *= muons_weight;
1136 ATH_MSG_DEBUG("Combined muon scale factor: " << muons_weight);
1137 }
1138
1139 //--- Jets
1140 if (m_slices["jet"]) {
1141 ATH_MSG_DEBUG("Working on jets");
1142 for ( const auto jet : *jets ) {
1143 ATH_MSG_VERBOSE( " Jet is bad? " << static_cast<int>(ST::acc_bad(*jet)));
1144 ATH_MSG_VERBOSE( " Jet is baseline ? " << static_cast<int>(ST::acc_baseline(*jet)));
1145 ATH_MSG_VERBOSE( " Jet passes OR ? " << static_cast<int>(ST::acc_passOR(*jet)));
1146 if (ST::acc_passOR(*jet)) {
1147 if (m_slices["bjet"] && (isNominal || syst_affectsJets || syst_affectsBTag)) {
1148 m_SUSYTools->IsBJetContinuous(*jet);
1149 }
1150 ATH_MSG_VERBOSE( " Jet is bjet ? " << static_cast<int>(ST::acc_bjet(*jet)));
1151 }
1152 }
1153
1154 float jet_weight(1.);
1155 float btag_weight(1.);
1156 if( !isData ) { //isMC
1157 if (isNominal) { //btagging
1158 btag_weight_nominal = btag_weight = m_SUSYTools->BtagSF(jets);
1159 weight_btags->SetBinContent(1, weight_btags->GetBinContent(1)+btag_weight);
1160 }
1161 else if (!syst_affectsBTag){
1162 btag_weight = btag_weight_nominal;
1163 }
1164 else{
1165 btag_weight = m_SUSYTools->BtagSF(jets);
1166 size_t iwbin = find(m_syst_weights["BTag"].begin(), m_syst_weights["BTag"].end(), sys.name()) - m_syst_weights["BTag"].begin();
1167 if(iwbin < m_syst_weights["BTag"].size()) { weight_btags->SetBinContent(iwbin+1, weight_btags->GetBinContent(iwbin+1)+btag_weight); }
1168 }
1169
1170 if(isNominal){ //JVT
1171 jet_weight = m_SUSYTools->JVT_SF(jets);
1172 jet_weight *= m_SUSYTools->FJVT_SF(jets);
1173 jets_weight_nominal = jet_weight;
1174 jets_weight_nominal *= m_SUSYTools->FJVT_SF(jets);
1175 weight_jets->SetBinContent(1, weight_jets->GetBinContent(1)+jet_weight);
1176 }
1177 else if (!syst_affectsJets || (syst_affectsJets && !sysInfo.affectsWeights)){
1178 jet_weight = jets_weight_nominal;
1179 }
1180 else if ( syst_affectsJets && sysInfo.affectsWeights ){
1181 jet_weight = m_SUSYTools->JVT_SF(jets);
1182 jet_weight *= m_SUSYTools->FJVT_SF(jets);
1183 size_t iwbin = find(m_syst_weights["Jet"].begin(), m_syst_weights["Jet"].end(), sys.name()) - m_syst_weights["Jet"].begin();
1184 if(iwbin < m_syst_weights["Jet"].size()) { weight_jets->SetBinContent(iwbin+1, weight_jets->GetBinContent(iwbin+1)+jet_weight); }
1185 }
1186 }
1187 else{ //data
1188 weight_jets->SetBinContent(1, 1.);
1189 weight_btags->SetBinContent(1, 1.);
1190 }
1191 event_weight *= jet_weight;
1192 event_weight *= btag_weight;
1193 ATH_MSG_DEBUG("Combined b-tagging scale factor: " << btag_weight);
1194 ATH_MSG_DEBUG("Combined jet scale factor: " << jet_weight);
1195 }
1196
1197 //--- Large R. Jets
1198
1199 if( m_slices["fjet"] && doFatJets ) {
1200 ATH_MSG_DEBUG("Working on fat jets");
1201 for ( const auto fatjet : *fatjets ) {
1202 ATH_MSG_VERBOSE( " Jet is bad? " << static_cast<int>(ST::acc_bad(*fatjet)));
1203 ATH_MSG_VERBOSE( " Jet is baseline ? " << static_cast<int>(ST::acc_baseline(*fatjet)));
1204 ATH_MSG_VERBOSE( " Jet passes OR ? " << static_cast<int>(ST::acc_passOR(*fatjet)));
1205 }
1206
1207 float fatjet_weight(1.);
1208 //disable - noJVTforfatjets if( !isData ) { //isMC
1209 //disable - noJVTforfatjets if(isNominal){ //JVT - no JVT cuts for fat jets!
1210 //disable - noJVTforfatjets fatjets_weight_nominal = fatjet_weight = m_SUSYTools->JVT_SF(fatjets);
1211 //disable - noJVTforfatjets weight_fatjets->SetBinContent(1, weight_fatjets->GetBinContent(1)+fatjet_weight);
1212 //disable - noJVTforfatjets }
1213 //disable - noJVTforfatjets if (!syst_affectsJets || (syst_affectsJets && !sysInfo.affectsWeights)){
1214 //disable - noJVTforfatjets fatjet_weight = fatjets_weight_nominal;
1215 //disable - noJVTforfatjets }
1216 //disable - noJVTforfatjets else if ( syst_affectsJets && sysInfo.affectsWeights ){
1217 //disable - noJVTforfatjets fatjet_weight = m_SUSYTools->JVT_SF(fatjets);
1218 //disable - noJVTforfatjets size_t iwbin = find(m_syst_weights["Jet"].begin(), m_syst_weights["Jet"].end(), sys.name()) - m_syst_weights["Jet"].begin();
1219 //disable - noJVTforfatjets if(iwbin < m_syst_weights["Jet"].size()) { weight_fatjets->SetBinContent(iwbin+1, weight_fatjets->GetBinContent(iwbin+1)+fatjet_weight); }
1220 //disable - noJVTforfatjets }
1221 //disable - noJVTforfatjets }
1222
1223 ATH_MSG_DEBUG("Combined large radius jet scale factor: " << fatjet_weight);
1224 }
1225
1226 //--- Track Jets
1227
1228 if( m_slices["tjet"] && doTrkJets ) {
1229 ATH_MSG_DEBUG("Working on trk jets");
1230 for ( const auto trkjet : *trkjets ) {
1231 ATH_MSG_VERBOSE( " Jet is bad? " << static_cast<int>(ST::acc_bad(*trkjet)));
1232 ATH_MSG_VERBOSE( " Jet is baseline ? " << static_cast<int>(ST::acc_baseline(*trkjet)));
1233 ATH_MSG_VERBOSE( " Jet passes OR ? " << static_cast<int>(ST::acc_passOR(*trkjet)));
1234 }
1235
1236 float trkjet_weight(1.);
1237 //disable - noJVTfortrkjets if( !isData ) { //isMC
1238 //disable - noJVTfortrkjets if(isNominal){ //JVT - no JVT cuts for trk jets!
1239 //disable - noJVTfortrkjets trkjets_weight_nominal = trkjet_weight = 1;//m_SUSYTools->JVT_SF(trkjets);
1240 //disable - noJVTfortrkjets weight_trkjets->SetBinContent(1, weight_trkjets->GetBinContent(1)+trkjet_weight);
1241 //disable - noJVTfortrkjets }
1242 //disable - noJVTfortrkjets else if (!syst_affectsJets || (syst_affectsJets && !sysInfo.affectsWeights)){
1243 //disable - noJVTfortrkjets trkjet_weight = trkjets_weight_nominal;
1244 //disable - noJVTfortrkjets }
1245 //disable - noJVTfortrkjets else if ( syst_affectsJets && sysInfo.affectsWeights ){
1246 //disable - noJVTfortrkjets trkjet_weight = 1;//m_SUSYTools->JVT_SF(trkjets);
1247 //disable - noJVTfortrkjets size_t iwbin = find(m_syst_weights["Jet"].begin(), m_syst_weights["Jet"].end(), sys.name()) - m_syst_weights["Jet"].begin();
1248 //disable - noJVTfortrkjets if(iwbin < m_syst_weights["Jet"].size()) { weight_trkjets->SetBinContent(iwbin+1, weight_trkjets->GetBinContent(iwbin+1)+trkjet_weight); }
1249 //disable - noJVTfortrkjets }
1250 //disable - noJVTfortrkjets }
1251
1252 ATH_MSG_DEBUG("Combined track jet scale factor: " << trkjet_weight);
1253 }
1254
1255 //--- Taus
1256 if (m_slices["tau"]) {
1257 ATH_MSG_DEBUG("Working on taus");
1258 float taus_weight(1.);
1259 for ( const auto ta : *taus ) {
1260 if( !isData ){
1261 if (isNominal || syst_affectsTaus) {
1262 if ((ST::acc_signal(*ta) == 1) && (isNominal || sysInfo.affectsWeights)) {
1263 taus_weight *= m_SUSYTools->GetSignalTauSF(*ta, true, false); //true, "HLT_tau25_medium1_tracktwo");;
1264 }
1265 }
1266 }
1267 ATH_MSG_VERBOSE( " Tau passing baseline selection? " << static_cast<int>(ST::acc_baseline(*ta)));
1268 ATH_MSG_VERBOSE( " Tau passing signal selection? " << static_cast<int>(ST::acc_signal(*ta)));
1269 if (ST::acc_signal(*ta) == 1 && ST::acc_effscalefact.isAvailable(*ta))
1270 ATH_MSG_VERBOSE( " Tau weight " << ST::acc_effscalefact(*ta) );
1271 }
1272
1273 if (isNominal) {
1274 taus_weight_nominal = taus_weight;
1275 weight_taus->SetBinContent(1, weight_taus->GetBinContent(1)+taus_weight);
1276 }
1277 else if (!syst_affectsTaus) {
1278 taus_weight = taus_weight_nominal;
1279 }
1280 else if ( sysInfo.affectsWeights ){
1281 size_t iwbin = find(m_syst_weights["Tau"].begin(), m_syst_weights["Tau"].end(), sys.name()) - m_syst_weights["Tau"].begin();
1282 if(iwbin < m_syst_weights["Tau"].size()) { weight_taus->SetBinContent(iwbin+1, weight_taus->GetBinContent(iwbin+1)+taus_weight); }
1283 }
1284
1285 event_weight *= taus_weight;
1286 }
1287
1288
1289 ATH_MSG_DEBUG("Full event weight: " << event_weight);
1290 if (isNominal) {event_weight_nominal = event_weight;}
1291 else if (sysInfo.affectsWeights) ATH_MSG_DEBUG("Difference with nominal weight: " << event_weight - event_weight_nominal);
1292
1293 // Clean up the systematics copies
1294 if (sysInfo.affectsKinematics) {
1295 delete metcst;
1296 delete mettst;
1297 delete metcst_aux;
1298 delete mettst_aux;
1299 }
1300 isNominal = false;
1301 }
1302
1303 m_clock0.Stop();
1304 if (m_Nevts == 0) {
1305 m_clock1.Stop();
1306 } else {
1307 m_clock2.Stop();
1308 }
1309 ++m_Nevts;
1310
1311 delete m_metcst_nominal;
1312 delete m_metcst_nominal_aux;
1313 delete m_mettst_nominal;
1314 delete m_mettst_nominal_aux;
1315
1316 return StatusCode::SUCCESS;
1317}
1318
1319//====================================================================================================
1320// Book histograms
1321StatusCode SUSYToolsAlg::bookHistograms(void) {
1322 // General info
1323 ATH_CHECK( book(TH1D("EventInfo/info_runNo", ";RunNumber;N", 1500, 250000, 400000)) );
1324 ATH_CHECK( book(TH1D("EventInfo/info_evtNo", ";EventNumber;N", 600, 100000000, 400000000)) );
1325 ATH_CHECK( book(TH1D("EventInfo/info_lumiBlock", ";LumiBlock;N", 1000, 0, 1000)) );
1326
1327 // Kinematic histograms
1328 m_objects = {"el", "ph", "mu", "tau", "jet", "bjet", "fatjet", "trkjet"};
1329 m_levels = {"nom", "bsl", "sig"};
1330 m_vars = {"pt", "eta", "phi", "e", "mass", "d0", "z0", "isolCalo", "isolTrack", "truthType", "truthOrigin", "N", "bweight", "bweightpb", "bweightpc", "bweightpu", "pid", "parentpid", "wtagged", "ztagged", "toptagged", "toptaggedscore","nTracks","nTracksCharged","nTracksIsolation","RNNJetScoreSigTrans"};
1331 std::map<std::string,std::vector<std::string>> cfg_hist_labels = {
1332 {"pt",{"p_{T} [GeV]","N"}}, {"eta",{"#eta","N"}}, {"phi",{"#phi","N"}}, {"e",{"energy","N"}}, {"mass",{"mass","N"}},
1333 {"d0",{"d0 significance","N"}}, {"z0",{"z0 * sin(#theta)","N"}}, {"isolCalo",{"isolation (calo) / pT","N"}}, {"isolTrack",{"isolation (track) / pT","N"}},
1334 {"truthType",{"truth type","N"}}, {"truthOrigin",{"truth origin","N"}}, {"N",{"count","N"}},
1335 {"bweight",{"b-weight","N"}}, {"bweightpb",{"b-weight (pb)","N"}}, {"bweightpc",{"b-weight (pc)","N"}}, {"bweightpu",{"b-weight (pu)","N"}},
1336 {"pid", {"pdg ID", "N"}}, {"parentpid", {"parent pdg ID","N"}},
1337 {"wtagged", {"W tagged", "N"}}, {"ztagged", {"Z tagged", "N"}}, {"toptagged", {"Top tagged", "N"}}, {"toptaggedscore", {"Top tagger score", "N"}},
1338 {"nTracks", {"N tracks", "N"}}, {"nTracksCharged", {"N tracks charged", "N"}}, {"nTracksIsolation", {"N tracks isolation", "N"}}, {"RNNJetScoreSigTrans", {"RNNJetScoreSigTrans","N"}} };
1339 std::map<std::string,int> cfg_hist_nbins = {
1340 {"pt",100}, {"eta",40}, {"phi",64}, {"e",100}, {"mass",100},
1341 {"d0",50}, {"z0",60}, {"isolCalo",70}, {"isolTrack",70},
1342 {"truthType",50}, {"truthOrigin",50}, {"N",16},
1343 {"bweight",200}, {"bweightpb",200}, {"bweightpc",200}, {"bweightpu",200},
1344 {"pid",61}, {"parentpid",61},
1345 {"wtagged",3}, {"ztagged",3}, {"toptagged",3}, {"toptaggedscore",50},
1346 {"nTracks",20}, {"nTracksCharged",20}, {"nTracksIsolation",5}, {"RNNJetScoreSigTrans",50} };
1347 std::map<std::string,std::vector<float>> cfg_hist_minmax = {
1348 {"pt",{0,200}}, {"eta",{-4,4}}, {"phi",{-3.2,3.2}}, {"e",{0,200}}, {"mass",{0,200}},
1349 {"d0",{-10,10}}, {"z0",{-1.5,1.5}}, {"isolCalo",{-0.2,0.5}}, {"isolTrack",{-0.2,0.5}},
1350 {"truthType",{0,50}}, {"truthOrigin",{0,50}}, {"N",{0,16}},
1351 {"bweight",{-10,10}}, {"bweightpb",{-0.5,1.5}}, {"bweightpc",{-0.5,1.5}}, {"bweightpu",{-0.5,1.5}},
1352 {"pid",{-30.5,30.5}}, {"parentpid",{-30.5,30.5}},
1353 {"wtagged", {-1,2}}, {"ztagged", {-1,2}}, {"toptagged", {-1,2}}, {"toptaggedscore", {0,1}},
1354 {"nTracks", {0,20}}, {"nTracksCharged", {0,20}}, {"nTracksIsolation", {0,5}}, {"RNNJetScoreSigTrans", {0,1}} };
1355 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"} };
1356 std::map<std::string,std::string> labels_levels = { {"nom","Nominal"}, {"bsl","Baseline"}, {"sig","Signal"} };
1357 std::map<std::string,std::string> labels_dir = { {"el","Electron"}, {"ph","Photon"}, {"mu","Muon"}, {"jet","Jet"}, {"bjet","bJet"}, {"tau","Tau"}, {"fatjet","LargeRJet"}, {"trkjet","TrackJet"} };
1358
1359 for (const auto& obj : m_objects) {
1360 m_obj_count[obj] = std::map<std::string,int>();
1361 for (const auto& lev : m_levels) {
1362 m_obj_count[obj][lev] = 0;
1363 }
1364 }
1365
1366 for (const auto& obj : m_objects) {
1367 for (const auto& lev : m_levels) {
1368 for (const auto& var : m_vars) {
1369 if (var.find("bweight")!=std::string::npos && obj.compare("bjet")!=0) continue; // bweights only for bjets
1370 if (var.find("tagged")!=std::string::npos && obj.compare("fatjet")!=0) continue; // boson tagging only for fjets
1371 if ((var.find("nTracks")!=std::string::npos||var.find("RNNJetScore")!=std::string::npos) && obj.compare("tau")!=0) continue; // nTracks/RNN score only for taus
1372 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
1373 std::string key = labels_dir[obj] + "/" + obj + "_" + lev + "_" + var;
1374 std::string labels = ";"+labels_levels[lev]+" "+labels_objects[obj]+" "+cfg_hist_labels[var][0]+";"+cfg_hist_labels[var][1];
1375 ATH_CHECK( book(TH1D(key.c_str(), labels.c_str(), cfg_hist_nbins[var], cfg_hist_minmax[var][0], cfg_hist_minmax[var][1])) );
1376 m_hists[key] = hist(key);
1377 ATH_MSG_INFO("Defined histogram: " << key.c_str() << ", " << m_hists[key]);
1378 }
1379 }
1380 }
1381
1383 #ifdef XAOD_STANDALONE
1384 for (std::string obj : {"el","mu","ph"}) {
1385 for (auto trg : m_triggers[obj]) {
1386 for (std::string var : {"pt","eta","phi"} ) {
1387 std::string key = "Trigger/"+obj+"_"+var+"_"+trg;
1388 std::string labels = ";"+labels_objects[obj]+" "+cfg_hist_labels[var][0]+";Efficiency "+trg;
1389 ATH_CHECK(book(TEfficiency(key.c_str(), labels.c_str(), cfg_hist_nbins[var], cfg_hist_minmax[var][0], cfg_hist_minmax[var][1])));
1390 m_heffs[key] = histeff(key);
1391 ATH_MSG_INFO("Defined histogram: " << key.c_str() << ", " << m_heffs[key]);
1392 }
1393 }
1394 }
1395 #endif
1396 ATH_CHECK( book(TH1D("Trigger/isMETTrigPassed", "isMETTrigPassed", 2, -0.5, 1.5)) );
1397
1398 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")) ) );
1399 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")) ) );
1400 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")) ) );
1401
1402 // Extra histograms
1403 ATH_CHECK( book(TH1D("Cutflow/el_n_flow_nominal", "Electron Cutflow (Nominal);Electron Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1404 ATH_CHECK( book(TH1D("Cutflow/ph_n_flow_nominal", "Photon Cutflow (Nominal);Photon Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1405 ATH_CHECK( book(TH1D("Cutflow/mu_n_flow_nominal", "Muon Cutflow (Nominal);Muon Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1406 ATH_CHECK( book(TH1D("Cutflow/jet_n_flow_nominal", "Jet Cutflow (Nominal);Jet Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1407 ATH_CHECK( book(TH1D("Cutflow/bjet_n_flow_nominal", "b-jet Cutflow (Nominal);b-jet Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1408 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) ) );
1409 ATH_CHECK( book(TH1D("Cutflow/trkjet_n_flow_nominal", "Track Jet Cutflow (Nominal);Track Jet Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1410 ATH_CHECK( book(TH1D("Cutflow/tau_n_flow_nominal", "Tau Cutflow (Nominal);Tau Cutflow (Nominal);N", Cut::NSel, 0, Cut::NSel) ) );
1411
1412 ATH_CHECK( book(TH1D("MET/met_et", "MET (Nominal);MET (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) ); //MET (+Components)
1413 ATH_CHECK( book(TH1D("MET/met_phi", "MET_phi (Nominal);MET_phi (Nominal);N", 50, -5, 5) ) );
1414 ATH_CHECK( book(TH1D("MET/met_sumet", "MET_sumet (Nominal);MET_sumet (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1415 ATH_CHECK( book(TH1D("MET/met_et_tst", "MET [PVSoftTrk] (Nominal);MET [PVSoftTrk] (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1416 ATH_CHECK( book(TH1D("MET/met_et_el", "MET [RefEle] (Nominal);MET [RefEle] (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1417 ATH_CHECK( book(TH1D("MET/met_et_ph", "MET [RefGamma] (Nominal);MET [RefGamma] (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1418 ATH_CHECK( book(TH1D("MET/met_et_mu", "MET [Muons] (Nominal);MET [Muons] (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1419 ATH_CHECK( book(TH1D("MET/met_et_jet", "MET [RefJet] (Nominal);MET [RefJet] (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1420 //ATH_CHECK( book(TH1D("MET/met_et_tau", "MET [RefTau] (Nominal);MET [RefTau] (Nominal) [GeV];N / [10 GeV]", 50, 0, 500) ) );
1421 ATH_CHECK( book(TH1D("MET/met_significance", "MET Significance;MET Significance;N / [0.5]", 60, 0, 30) ) );
1422
1423 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
1424 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
1425 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
1426 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
1427 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
1428 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
1429 //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
1430 //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
1431 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
1432
1433 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
1434 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
1435 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
1436 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
1437 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
1438 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
1439 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
1440
1441 return StatusCode::SUCCESS;
1442}
1443
1444//====================================================================================================
1445// collect systematics per object (all / weight-only)
1446void SUSYToolsAlg::groupSysts(void) {
1447 m_syst_all.clear();
1448 m_syst_weights.clear();
1449 const std::vector<std::string> Nominal{"Nominal"};
1450 for (const auto& sysInfo : m_sysInfoList) {
1451 std::string sys_name = sysInfo.systset.name();
1452 std::string sys_affects = ST::testAffectsObject(sysInfo.affectsType);
1453 ATH_MSG_DEBUG("Syst " << sys_name << " affects " << sys_affects);
1454
1455 bool sys_affects_weights = sysInfo.affectsWeights;
1456
1457 // collect all syst names
1458 // per affected object
1459 const auto & [pairPtr, inserted] = m_syst_all.try_emplace(sys_affects, Nominal);
1460 pairPtr->second.push_back(sys_name);
1461
1462 // weight related syst
1463 if (sys_affects_weights) {
1464 const auto & [pairPtr_w, inserted_w] = m_syst_weights.try_emplace(sys_affects, Nominal);
1465 pairPtr_w->second.push_back(sys_name);
1466 }
1467 }
1468}
1469
1470//====================================================================================================
1471// Fill object histograms consistently for different selection levels / object types
1472void SUSYToolsAlg::stdHistsForObj(xAOD::IParticle *obj, const std::string& objtype, const std::string& objlevel, std::map<std::string,std::string> config) {
1473 std::map<std::string,std::string> labels_dir = { {"el","Electron"}, {"ph","Photon"}, {"mu","Muon"}, {"jet","Jet"}, {"bjet","bJet"}, {"tau","Tau"}, {"fatjet","LargeRJet"}, {"trkjet","TrackJet"} };
1474 std::string dir = labels_dir[objtype]+"/";
1475
1476 ATH_MSG_VERBOSE("Filling histograms for " << obj->type() << ", " << objlevel);
1477 hist(dir+objtype+"_"+objlevel+"_pt")->Fill( obj->pt()/1000. );
1478 hist(dir+objtype+"_"+objlevel+"_eta")->Fill( obj->eta() );
1479 hist(dir+objtype+"_"+objlevel+"_phi")->Fill( obj->phi() );
1480 hist(dir+objtype+"_"+objlevel+"_e")->Fill( obj->e()/1000. );
1481 hist(dir+objtype+"_"+objlevel+"_mass")->Fill( obj->m()/1000. );
1482 //
1483 if (objtype=="el" || objtype=="mu" || objtype=="ph") {
1484 if (objtype=="el") {
1485 if (ST::acc_z0sinTheta.isAvailable(*obj))
1486 hist(dir+objtype+"_"+objlevel+"_z0")->Fill( ST::acc_z0sinTheta(*obj) );
1487 if (ST::acc_d0sig.isAvailable(*obj))
1488 hist(dir+objtype+"_"+objlevel+"_d0")->Fill( ST::acc_d0sig(*obj) );
1489 hist(dir+objtype+"_"+objlevel+"_isolCalo")->Fill( STAlg::acc_topoetcone20(*obj)/obj->pt() );
1490 hist(dir+objtype+"_"+objlevel+"_isolTrack")->Fill( STAlg::acc_ptvarcone30_TTVA_LooseCone(*obj)/obj->pt() );
1491 } else if (objtype=="mu") {
1492 if (ST::acc_z0sinTheta.isAvailable(*obj))
1493 hist(dir+objtype+"_"+objlevel+"_z0")->Fill( ST::acc_z0sinTheta(*obj) );
1494 if (ST::acc_d0sig.isAvailable(*obj))
1495 hist(dir+objtype+"_"+objlevel+"_d0")->Fill( ST::acc_d0sig(*obj) );
1496 hist(dir+objtype+"_"+objlevel+"_isolCalo")->Fill( STAlg::acc_topoetcone20(*obj)/obj->pt() );
1497 hist(dir+objtype+"_"+objlevel+"_isolTrack")->Fill( STAlg::acc_ptvarcone30_TTVA(*obj)/obj->pt() );
1498 } else if (objtype=="ph") {
1499 hist(dir+objtype+"_"+objlevel+"_isolCalo")->Fill( (STAlg::acc_topoetcone40(*obj)-2450)/obj->pt() );
1500 hist(dir+objtype+"_"+objlevel+"_isolTrack")->Fill( STAlg::acc_ptcone20(*obj)/obj->pt() );
1501 }
1502 }
1503 //
1504 if (objtype=="el" || objtype=="mu" || objtype=="ph" || objtype=="tau") {
1505 if ( objtype!="tau" ) {
1506 if (ST::acc_truthType.isAvailable(*obj) && ST::acc_truthOrigin.isAvailable(*obj)) {
1507 ATH_MSG_DEBUG(objtype << ", " << objlevel << ": " << ST::acc_truthType(*obj) << ", " << ST::acc_truthOrigin(*obj));
1508 hist(dir+objtype+"_"+objlevel+"_truthType")->Fill( ST::acc_truthType(*obj) );
1509 hist(dir+objtype+"_"+objlevel+"_truthOrigin")->Fill( ST::acc_truthOrigin(*obj) );
1510 }
1511 } else {
1512 if(!m_isPHYSLITE){
1513 bool istruthmatched = STAlg::acc_IsTruthMatched.isAvailable(*obj) ? (bool)STAlg::acc_IsTruthMatched(*obj) : false;
1514 int pid(0),ppid(0);
1515 if (istruthmatched && STAlg::acc_truthParticleLink.isAvailable(*obj)) {
1516 const auto *tp = *(STAlg::acc_truthParticleLink(*obj));
1517 if (tp) {
1518 pid = tp->pdgId();
1519 ppid = (tp->nParents()>0)?tp->parent(0)->pdgId():0;
1520 }
1521 }
1522 hist(dir+objtype+"_"+objlevel+"_pid")->Fill( pid );
1523 hist(dir+objtype+"_"+objlevel+"_parentpid")->Fill( ppid );
1524 }
1525 }
1526 }
1527
1528 if (objtype=="tau") {
1529 hist(dir+objtype+"_"+objlevel+"_nTracks")->Fill( dynamic_cast<xAOD::TauJet*>(obj)->nTracks() );
1530 hist(dir+objtype+"_"+objlevel+"_nTracksCharged")->Fill( dynamic_cast<xAOD::TauJet*>(obj)->nTracksCharged() );
1531 hist(dir+objtype+"_"+objlevel+"_nTracksIsolation")->Fill( dynamic_cast<xAOD::TauJet*>(obj)->nTracksIsolation() );
1532 hist(dir+objtype+"_"+objlevel+"_RNNJetScoreSigTrans")->Fill( STAlg::acc_RNNJetScoreSigTrans(*obj) );
1533 }
1534 //
1535
1536 if (objtype=="bjet") {
1537 hist(dir+objtype+"_"+objlevel+"_bweight")->Fill( ST::acc_btag_weight(*obj) );
1538 hist(dir+objtype+"_"+objlevel+"_bweightpb")->Fill( ST::acc_btag_dl1pb(*obj) );
1539 hist(dir+objtype+"_"+objlevel+"_bweightpc")->Fill( ST::acc_btag_dl1pc(*obj) );
1540 hist(dir+objtype+"_"+objlevel+"_bweightpu")->Fill( ST::acc_btag_dl1pu(*obj) );
1541 }
1542
1543 //
1544 const static SG::ConstAccessor<bool> acc_topValidPtRangeLow(config["ToptaggerKey"]+"_ValidPtRangeLow");
1545 const static SG::ConstAccessor<bool> acc_topValidPtRangeHigh(config["ToptaggerKey"]+"_ValidPtRangeHigh");
1546 const static SG::ConstAccessor<bool> acc_topValidEtaRange(config["ToptaggerKey"]+"_ValidEtaRange");
1547 const static SG::ConstAccessor<bool> acc_wtagged(config["WtaggerKey"]+"_Tagged");
1548 const static SG::ConstAccessor<bool> acc_ztagged(config["ZtaggerKey"]+"_Tagged");
1549 const static SG::ConstAccessor<bool> acc_toptagged(config["ToptaggerKey"]+"_Tagged");
1550 const static SG::ConstAccessor<float> acc_topscore(config["ToptaggerKey"]+"_Score");
1551
1552 if (objtype=="fatjet" && acc_topValidPtRangeLow.isAvailable(*obj) && acc_topValidPtRangeHigh.isAvailable(*obj) && acc_topValidEtaRange.isAvailable(*obj)) {
1553 bool kin = acc_topValidPtRangeLow(*obj) && acc_topValidPtRangeHigh(*obj) && acc_topValidEtaRange(*obj);
1554 if (config.find("WtaggerKey")!=config.end()) hist(dir+objtype+"_"+objlevel+"_wtagged") ->Fill( acc_wtagged.isAvailable(*obj) ? kin && acc_wtagged(*obj) : -1);
1555 if (config.find("ZtaggerKey")!=config.end()) hist(dir+objtype+"_"+objlevel+"_ztagged") ->Fill( acc_ztagged.isAvailable(*obj) ? kin && acc_ztagged(*obj) : -1);
1556 if (config.find("ToptaggerKey")!=config.end()) hist(dir+objtype+"_"+objlevel+"_toptagged")->Fill( acc_toptagged.isAvailable(*obj) ? kin && acc_toptagged(*obj) : -1);
1557 if (config.find("ToptaggerKey")!=config.end()) hist(dir+objtype+"_"+objlevel+"_toptaggedscore")->Fill( acc_topscore.isAvailable(*obj) ? acc_topscore(*obj) : -1);
1558 }
1559}
1560
1561
1562//====================================================================================================
1563std::string findInReg(SG::AuxTypeRegistry& reg, xAOD::IParticle* obj, const std::string& searchkey) {
1564 std::regex re_tag(searchkey);
1565 std::smatch matches;
1566 for (auto x : obj->getAuxIDs()) {
1567 std::string xname = reg.getName(x);
1568 std::regex_search(xname, matches, re_tag);
1569 if (matches.size()>1 and !matches[1].str().empty()) {
1570 return matches[1].str();
1571 }
1572 }
1573 return "";
1574}
1575
1576//====================================================================================================
1577int getSize(std::map<std::string,std::vector<std::string>> &collection, const std::string& object) {
1578 if (collection.find(object) == collection.end()) return 1;
1579 return static_cast<int>(collection[object].size());
1580}
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:1577
findInReg
std::string findInReg(SG::AuxTypeRegistry &reg, xAOD::IParticle *obj, const std::string &searchkey)
Definition SUSYToolsAlg.cxx:1563
x
#define x
AthCommonDataStore< AthCommonMsg< Algorithm > >::declareProperty
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)
Definition AthCommonDataStore.h:145
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:303
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:198
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:369
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:114
SUSYToolsAlg::m_heffs
std::map< std::string, TEfficiency * > m_heffs
Definition SUSYToolsAlg.h:81
SUSYToolsAlg::execute
StatusCode execute()
Definition SUSYToolsAlg.cxx:302
SUSYToolsAlg::m_grl
asg::AnaToolHandle< IGoodRunsListSelectionTool > m_grl
Definition SUSYToolsAlg.h:85
SUSYToolsAlg::m_mettst_nominal
xAOD::MissingETContainer * m_mettst_nominal
Definition SUSYToolsAlg.h:124
SUSYToolsAlg::m_clock2
TStopwatch m_clock2
Definition SUSYToolsAlg.h:106
SUSYToolsAlg::m_triggers
std::map< std::string, std::vector< std::string > > m_triggers
Definition SUSYToolsAlg.h:109
SUSYToolsAlg::m_mettst_syst
xAOD::MissingETContainer * m_mettst_syst
Definition SUSYToolsAlg.h:117
SUSYToolsAlg::m_syst_weights
std::map< std::string, std::vector< std::string > > m_syst_weights
Definition SUSYToolsAlg.h:101
SUSYToolsAlg::m_syst_all
std::map< std::string, std::vector< std::string > > m_syst_all
Definition SUSYToolsAlg.h:101
SUSYToolsAlg::m_metcst_nominal
xAOD::MissingETContainer * m_metcst_nominal
Definition SUSYToolsAlg.h:122
SUSYToolsAlg::m_slices
std::map< std::string, bool > m_slices
Definition SUSYToolsAlg.h:79
SUSYToolsAlg::m_metcst_syst
xAOD::MissingETContainer * m_metcst_syst
Definition SUSYToolsAlg.h:119
SUSYToolsAlg::m_obj_count
std::map< std::string, std::map< std::string, int > > m_obj_count
Definition SUSYToolsAlg.h:115
SUSYToolsAlg::m_kernel
std::string m_kernel
Definition SUSYToolsAlg.h:92
SUSYToolsAlg::m_TrkJetCollection
std::string m_TrkJetCollection
Definition SUSYToolsAlg.h:95
SUSYToolsAlg::m_objects
std::vector< std::string > m_objects
Definition SUSYToolsAlg.h:112
SUSYToolsAlg::m_mcCampaign
std::string m_mcCampaign
Definition SUSYToolsAlg.h:78
SUSYToolsAlg::m_mettst_nominal_aux
xAOD::MissingETAuxContainer * m_mettst_nominal_aux
Definition SUSYToolsAlg.h:125
SUSYToolsAlg::m_hists
std::map< std::string, TH1 * > m_hists
Definition SUSYToolsAlg.h:80
SUSYToolsAlg::m_SUSYTools
asg::AnaToolHandle< ST::ISUSYObjDef_xAODTool > m_SUSYTools
Definition SUSYToolsAlg.h:84
SUSYToolsAlg::m_metcst_nominal_aux
xAOD::MissingETAuxContainer * m_metcst_nominal_aux
Definition SUSYToolsAlg.h:123
SUSYToolsAlg::m_configFile
std::string m_configFile
Definition SUSYToolsAlg.h:93
SUSYToolsAlg::groupSysts
void groupSysts(void)
Definition SUSYToolsAlg.cxx:1446
SUSYToolsAlg::bookHistograms
StatusCode bookHistograms(void)
Definition SUSYToolsAlg.cxx:1321
SUSYToolsAlg::m_tauTruthMatchingTool
asg::AnaToolHandle< TauAnalysisTools::ITauTruthMatchingTool > m_tauTruthMatchingTool
Definition SUSYToolsAlg.h:86
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:1472
SUSYToolsAlg::m_mettst_syst_aux
xAOD::MissingETAuxContainer * m_mettst_syst_aux
Definition SUSYToolsAlg.h:118
SUSYToolsAlg::m_clock1
TStopwatch m_clock1
Definition SUSYToolsAlg.h:105
SUSYToolsAlg::m_levels
std::vector< std::string > m_levels
Definition SUSYToolsAlg.h:113
SUSYToolsAlg::m_metcst_syst_aux
xAOD::MissingETAuxContainer * m_metcst_syst_aux
Definition SUSYToolsAlg.h:120
SUSYToolsAlg::m_clock0
TStopwatch m_clock0
Definition SUSYToolsAlg.h:104
SUSYToolsAlg::initialize
StatusCode initialize()
Definition SUSYToolsAlg.cxx:103
SUSYToolsAlg::m_sysInfoList
std::vector< ST::SystInfo > m_sysInfoList
Definition SUSYToolsAlg.h:100
SUSYToolsAlg::m_Nevts
unsigned int m_Nevts
Definition SUSYToolsAlg.h:88
SUSYToolsAlg::m_configDict
std::map< std::string, std::string > m_configDict
Definition SUSYToolsAlg.h:82
SUSYToolsAlg::m_FatJetCollection
std::string m_FatJetCollection
Definition SUSYToolsAlg.h:94
SUSYToolsAlg::m_GRLFiles
std::vector< std::string > m_GRLFiles
Definition SUSYToolsAlg.h:77
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:114
find
std::string find(const std::string &s)
return a remapped string
Definition hcg.cxx:138
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_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:133
ST::acc_btag_dl1pu
static const SG::ConstAccessor< float > acc_btag_dl1pu("btag_dl1pu")
ST::acc_btag_dl1pb
static const SG::ConstAccessor< float > acc_btag_dl1pb("btag_dl1pb")
ST::acc_bjet
static const SG::ConstAccessor< char > acc_bjet("bjet")
ST::acc_cosmic
static const SG::ConstAccessor< char > acc_cosmic("cosmic")
ST::acc_passOR
static const SG::ConstAccessor< char > acc_passOR("passOR")
ST::acc_btag_dl1pc
static const SG::ConstAccessor< float > acc_btag_dl1pc("btag_dl1pc")
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:18
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:67
object
Generated on for ATLAS Offline Software by doxygen 1.14.0