ATLAS Offline Software
Loading...
Searching...
No Matches
SecVertexTruthMatchAlg.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#include "TrackingAnalysisAlgorithms/SecVertexTruthMatchAlg.h"
6#include "InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool.h"
7#include "xAODTruth/TruthParticle.h"
8
9#include "TH1.h"
10#include "TEfficiency.h"
11#include "TLorentzVector.h"
12#include <algorithm>
13#include <cmath>
14#include <numbers>
15
16namespace {
17 constexpr float GeV = 1000.;
18}
19
20namespace CP {
21
22 SecVertexTruthMatchAlg::SecVertexTruthMatchAlg( const std::string& name, ISvcLocator* svcLoc )
23 : EL::AnaAlgorithm( name, svcLoc ) {}
24
25 StatusCode SecVertexTruthMatchAlg::initialize() {
26
27 // Initializing Keys
28 ATH_CHECK(m_secVtxContainerKey.initialize());
29 ATH_CHECK(m_truthVtxContainerKey.initialize());
30 ATH_CHECK(m_trackParticleContainerKey.initialize());
31
32 // Retrieving the tool
33 ATH_CHECK(m_matchTool.retrieve());
34
35 if(m_writeHistograms) {
36 std::vector<std::string> recoTypes{"All", "Matched", "Merged", "Fake", "Split", "Other"};
37 std::vector<std::string> truthTypes{"Inclusive", "Reconstructable", "Accepted", "Seeded", "Reconstructed", "ReconstructedSplit"};
38 // define SM origin categories if SM origin tracking is enabled
39 std::vector<std::string> smOriginTypes;
40 if(m_doSMOrigin) {
41 smOriginTypes = {"FakeOrigin", "Pileup", "KshortDecay", "StrangeMesonDecay", "LambdaDecay",
42 "StrangeBaryonDecay", "TauDecay", "GammaConversion", "OtherDecay",
43 "HadronicInteraction", "OtherSecondary", "BHadronDecay", "DHadronDecay",
44 "Fragmentation", "OtherOrigin", "Signal"};
45 }
46
47 //determine histogram ranges depending if we are in normal or MuSA mode
48 //total bin counts stay the same for simplicity -- MuSA has less precision in general
49 float maxX = m_doMuSA ? 8000 : 500;
50 float maxY = m_doMuSA ? 10000 : 500;
51 float maxZ = m_doMuSA ? 10000 : 1500;
52 float maxLxy = m_doMuSA ? 8000 : 500;
53 float maxR = m_doMuSA ? 8000 : 600;
54 float mind0 = m_doMuSA ? 2000 : 100;
55 float maxd0 = m_doMuSA ? 2000 : 100;
56 float maxTrackd0 = m_doMuSA ? 3000 : 300;
57 float maxTrackz0 = m_doMuSA ? 5000 : 500;
58 float maxErrd0 = m_doMuSA ? 300 : 30;
59 float maxErrz0 = m_doMuSA ? 500 : 50;
60
61 float maxResR = m_doMuSA ? 2000 : 20;
62 float maxResZ = m_doMuSA ? 2000 : 20;
63
64 ANA_CHECK (book(TH1F("RecoVertex/matchType", "Vertex Match Type", 65, -0.5, 64.5)));
65 if (m_doSMOrigin) {
66 ANA_CHECK (book(TH1F("RecoVertex/smOriginType", "Vertex SM Origin Type", 65537, -0.5, 65536.5)));
67 }
68
69 // book the reco vertex histograms of one category and cache their pointers
70 auto bookRecoVertexHistos = [&](const std::string& category, bool withMatchScore) -> StatusCode {
71 RecoVertexHists& h = m_recoHists[category];
72 auto bookHist = [&](TH1*& target, const std::string& suffix, const char* title,
73 int nbins, double low, double high) -> StatusCode {
74 const std::string name = "RecoVertex/" + category + suffix;
75 ANA_CHECK (book(TH1F(name.c_str(), title, nbins, low, high)));
76 target = hist(name);
77 return StatusCode::SUCCESS;
78 };
79
80 ANA_CHECK (bookHist(h.x, "_x", "Reco vertex x [mm]", 1000, -maxX, maxX));
81 ANA_CHECK (bookHist(h.y, "_y", "Reco vertex y [mm]", 1000, -maxY, maxY));
82 ANA_CHECK (bookHist(h.z, "_z", "Reco vertex z [mm]", 1000, -maxZ, maxZ));
83 ANA_CHECK (bookHist(h.Lxy, "_Lxy", "Reco vertex L_{xy} [mm]", 500, 0, maxLxy));
84 ANA_CHECK (bookHist(h.pT, "_pT", "Reco vertex p_{T} [GeV]", 100, 0, 100));
85 ANA_CHECK (bookHist(h.eta, "_eta", "Reco vertex #eta", 100, -5, 5));
86 ANA_CHECK (bookHist(h.phi, "_phi", "Reco vertex #phi", 100, -std::numbers::pi, std::numbers::pi));
87 ANA_CHECK (bookHist(h.mass, "_mass", "Reco vertex mass [GeV]", 500, 0, 100));
88 ANA_CHECK (bookHist(h.mu, "_mu", "Reco vertex Red. Mass [GeV]", 500, 0, 100));
89 ANA_CHECK (bookHist(h.chi2, "_chi2", "Reco vertex recoChi2", 100, 0, 10));
90 ANA_CHECK (bookHist(h.dir, "_dir", "Reco vertex recoDirection", 100, -1, 1));
91 ANA_CHECK (bookHist(h.charge, "_charge", "Reco vertex recoCharge", 20, -10, 10));
92 ANA_CHECK (bookHist(h.H, "_H", "Reco vertex H [GeV]", 100, 0, 100));
93 ANA_CHECK (bookHist(h.HT, "_HT", "Reco vertex Mass [GeV]", 100, 0, 100));
94 ANA_CHECK (bookHist(h.minOpAng, "_minOpAng", "Reco vertex minOpAng", 100, -1, 1));
95 ANA_CHECK (bookHist(h.maxOpAng, "_maxOpAng", "Reco vertex maxOpAng", 100, -1, 1));
96 ANA_CHECK (bookHist(h.maxdR, "_maxdR", "Reco vertex maxDR", 100, 0, 10));
97 ANA_CHECK (bookHist(h.mind0, "_mind0", "Reco vertex min d0 [mm]", 100, 0, mind0));
98 ANA_CHECK (bookHist(h.maxd0, "_maxd0", "Reco vertex max d0 [mm]", 100, 0, maxd0));
99 ANA_CHECK (bookHist(h.ntrk, "_ntrk", "Reco vertex n tracks", 30, 0, 30));
100
101 // tracks
102 ANA_CHECK (bookHist(h.Trk_qOverP, "_Trk_qOverP", "Reco track qOverP ", 100, 0, .01));
103 ANA_CHECK (bookHist(h.Trk_theta, "_Trk_theta", "Reco track theta ", 64, 0, 3.2));
104 ANA_CHECK (bookHist(h.Trk_E, "_Trk_E", "Reco track E ", 100, 0, 100));
105 ANA_CHECK (bookHist(h.Trk_M, "_Trk_M", "Reco track M ", 100, 0, 10));
106 ANA_CHECK (bookHist(h.Trk_Pt, "_Trk_Pt", "Reco track Pt ", 100, 0, 100));
107 ANA_CHECK (bookHist(h.Trk_Px, "_Trk_Px", "Reco track Px ", 100, 0, 100));
108 ANA_CHECK (bookHist(h.Trk_Py, "_Trk_Py", "Reco track Py ", 100, 0, 100));
109 ANA_CHECK (bookHist(h.Trk_Pz, "_Trk_Pz", "Reco track Pz ", 100, 0, 100));
110 ANA_CHECK (bookHist(h.Trk_Eta, "_Trk_Eta", "Reco track Eta ", 100, -5, 5));
111 ANA_CHECK (bookHist(h.Trk_Phi, "_Trk_Phi", "Reco track Phi ", 63, -3.2, 3.2));
112 ANA_CHECK (bookHist(h.Trk_D0, "_Trk_D0", "Reco track D0 ", 300, -maxTrackd0, maxTrackd0));
113 ANA_CHECK (bookHist(h.Trk_Z0, "_Trk_Z0", "Reco track Z0 ", 500, -maxTrackz0, maxTrackz0));
114 ANA_CHECK (bookHist(h.Trk_errD0, "_Trk_errD0", "Reco track errD0 ", 300, 0, maxErrd0));
115 ANA_CHECK (bookHist(h.Trk_errZ0, "_Trk_errZ0", "Reco track errZ0 ", 500, 0, maxErrz0));
116 ANA_CHECK (bookHist(h.Trk_Chi2, "_Trk_Chi2", "Reco track Chi2 ", 100, 0, 10));
117 ANA_CHECK (bookHist(h.Trk_nDoF, "_Trk_nDoF", "Reco track nDoF ", 100, 0, 100));
118 ANA_CHECK (bookHist(h.Trk_charge, "_Trk_charge", "Reco track charge ", 3, -1.5, 1.5));
119
120 if (withMatchScore) {
121 ANA_CHECK (bookHist(h.positionRes_R, "_positionRes_R", "Position resolution for vertices matched to truth decays", 400, -maxResR, maxResR));
122 ANA_CHECK (bookHist(h.positionRes_Z, "_positionRes_Z", "Position resolution for vertices matched to truth decays", 400, -maxResZ, maxResZ));
123 ANA_CHECK (bookHist(h.matchScore_weight, "_matchScore_weight", "Vertex Match Score (weight)", 101, 0, 1.01));
124 ANA_CHECK (bookHist(h.matchScore_pt, "_matchScore_pt", "Vertex Match Score (pT)", 101, 0, 1.01));
125 }
126 return StatusCode::SUCCESS;
127 };
128
129 for(const auto& recoType : recoTypes) {
130 // truth matching -- don't book for non-matched vertices
131 ANA_CHECK (bookRecoVertexHistos(recoType, recoType != "All" and recoType != "Fake"));
132 }
133
134 // do reco vertices by SM origin if enabled -- NOTE these types are not exclusive (d decays will also be b decays in a cascade etc)
135 for(const auto& smOriginType : smOriginTypes) {
136 ANA_CHECK (bookRecoVertexHistos(smOriginType, false));
137 }
138
139 for(const auto& truthType : truthTypes) {
140 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_x").c_str(), "Truth vertex x [mm]", 1000, -maxX, maxX)));
141 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_y").c_str(), "Truth vertex y [mm]", 500, -maxY, maxY)));
142 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_z").c_str(), "Truth vertex z [mm]", 500, -maxZ, maxZ)));
143 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_R").c_str(), "Truth vertex r [mm]", 6000, 0, maxR)));
144 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Eta").c_str(), "Truth vertex Eta", 100, -5, 5)));
145 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Phi").c_str(), "Truth vertex Phi", 64, -3.2, 3.2)));
146 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Ntrk_out").c_str(), "Truth vertex n outgoing tracks", 100, 0, 100)));
147 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Parent_E").c_str(), "Reco track E", 100, 0, 100)));
148 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Parent_M").c_str(), "Reco track M", 500, 0, 500)));
149 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Parent_Pt").c_str(), "Reco track Pt", 100, 0, 100)));
150 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Parent_Eta").c_str(), "Reco track Eta", 100, -5, 5)));
151 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Parent_Phi").c_str(), "Reco track Phi", 63, -3.2, 3.2)));
152 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Parent_charge").c_str(), "Reco track charge", 3, -1, 1)));
153 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_ParentProdX").c_str(), "truthParentProd vertex x [mm]", 500, -500, 500)));
154 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_ParentProdY").c_str(), "truthParentProd vertex y [mm]", 500, -500, 500)));
155 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_ParentProdZ").c_str(), "truthParentProd vertex z [mm]", 500, -500, 500)));
156 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_ParentProdR").c_str(), "truthParentProd vertex r [mm]", 6000, 0, 600)));
157 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_ParentProdEta").c_str(), "truthParentProd vertex Eta", 100, -5, 5)));
158 ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_ParentProdPhi").c_str(), "truthParentProd vertex Phi", 64, -3.2, 3.2)));
159 }
160
161 // now add the efficiencies
162 // Define two different bin arrays - one for standard mode, one for MuSA
163 const std::vector<double> bins = m_doMuSA
164 ? std::vector<double>{0.0, 1, 5, 10, 20, 50, 100, 200, 300, 500, 750, 1000, 1500, 2000, 2500, 3000, 4000, 5000, 6000, 7000, 8000}
165 : std::vector<double>{0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 15, 20, 25, 30, 35, 40, 50, 60, 70, 80, 90, 100, 125, 150, 200, 300, 500};
166 const int nbins = std::size(bins) - 1;
167
168 ANA_CHECK (book(TEfficiency("Acceptance", "Acceptance", nbins, bins.data())));
169 ANA_CHECK (book(TEfficiency("eff_seed", "Seed efficiency", nbins, bins.data())));
170 ANA_CHECK (book(TEfficiency("eff_core", "Core efficiency", nbins, bins.data())));
171 ANA_CHECK (book(TEfficiency("eff_total", "Total efficiency", nbins, bins.data())));
172
173 }
174
175
176 return StatusCode::SUCCESS;
177 }
178
179 StatusCode SecVertexTruthMatchAlg::execute(const EventContext& ctx) {
180
181 //Retrieve the vertices:
182 SG::ReadHandle<xAOD::VertexContainer> recoVertexContainer(m_secVtxContainerKey, ctx);
183 if (!recoVertexContainer.isValid()) {
184 ATH_MSG_ERROR("Failed to retrieve secondary vertex container " << m_secVtxContainerKey.key());
185 return StatusCode::FAILURE;
186 }
187 SG::ReadHandle<xAOD::TruthVertexContainer> truthVertexContainer(m_truthVtxContainerKey, ctx);
188 if (!truthVertexContainer.isValid()) {
189 ATH_MSG_ERROR("Failed to retrieve truth vertex container " << m_truthVtxContainerKey.key());
190 return StatusCode::FAILURE;
191 }
192 SG::ReadHandle<xAOD::TrackParticleContainer> trackParticleContainer(m_trackParticleContainerKey, ctx);
193 if (!trackParticleContainer.isValid()) {
194 ATH_MSG_ERROR("Failed to retrieve track particle container " << m_trackParticleContainerKey.key());
195 return StatusCode::FAILURE;
196 }
197
198 std::vector<const xAOD::Vertex*> recoVerticesToMatch;
199 std::vector<const xAOD::TruthVertex*> truthVerticesToMatch;
200
201 for(const auto *recoVertex : *recoVertexContainer) {
202 if(recoVertex->vertexType() != xAOD::VxType::SecVtx ){
203 ATH_MSG_DEBUG("Vertex not labeled as secondary");
204 continue;
205 }
206 recoVerticesToMatch.push_back(recoVertex);
207 }
208
209 for(const auto *truthVertex : *truthVertexContainer) {
210 if(truthVertex->nIncomingParticles() != 1) {
211 continue;
212 }
213 const xAOD::TruthParticle* truthPart = truthVertex->incomingParticle(0);
214 if(not truthPart) {
215 continue;
216 }
217 if(!std::ranges::contains(m_targetPDGIDs.value(), std::abs(truthPart->pdgId()))) {
218 continue;
219 }
220 if(truthVertex->nOutgoingParticles() < 2) {
221 continue;
222 }
223 truthVerticesToMatch.push_back(truthVertex);
224 }
225
226 //pass to the tool for decoration:
227 ATH_CHECK( m_matchTool->matchVertices( recoVerticesToMatch, truthVerticesToMatch, trackParticleContainer.cptr() ) );
228
229 if(m_writeHistograms) {
230 static const xAOD::Vertex::ConstAccessor<int> matchTypeAcc("vertexMatchType");
231 static const xAOD::Vertex::ConstAccessor<int> originTypeAcc("vertexMatchOriginType");
232
233 static const std::map<InDetSecVtxTruthMatchUtils::VertexMatchOriginType, std::string> originTypeMap = {
234 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::FakeOrigin, "FakeOrigin"},
235 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::Pileup, "Pileup"},
236 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::KshortDecay, "KshortDecay"},
237 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::StrangeMesonDecay, "StrangeMesonDecay"},
238 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::LambdaDecay, "LambdaDecay"},
239 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::StrangeBaryonDecay, "StrangeBaryonDecay"},
240 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::TauDecay, "TauDecay"},
241 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::GammaConversion, "GammaConversion"},
242 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::OtherDecay, "OtherDecay"},
243 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::HadronicInteraction, "HadronicInteraction"},
244 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::OtherSecondary, "OtherSecondary"},
245 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::BHadronDecay, "BHadronDecay"},
246 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::DHadronDecay, "DHadronDecay"},
247 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::Fragmentation, "Fragmentation"},
248 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::OtherOrigin, "OtherOrigin"},
249 {InDetSecVtxTruthMatchUtils::VertexMatchOriginType::Signal, "Signal"},
250 };
251 //initialise references before loop
252 const auto& matchedHists = m_recoHists.at("Matched");
253 const auto& mergedHists = m_recoHists.at("Merged");
254 const auto& fakeHists = m_recoHists.at("Fake");
255 const auto& splitHists = m_recoHists.at("Split");
256 const auto& otherHists = m_recoHists.at("Other");
257 const auto& allHists = m_recoHists.at("All");
258 //
259 std::vector<const RecoVertexHists*> categories;
260 //assuming the match types are mutually exclusive, but a reasonable guess in any case
261 categories.reserve(8);
262 //
263 for(const auto * secVtx : recoVerticesToMatch) {
264 categories.clear();
265 const int matchTypeBitset = matchTypeAcc(*secVtx);
266 hist("RecoVertex/matchType")->Fill(matchTypeBitset);
267
268 if(InDetSecVtxTruthMatchUtils::isMatched(matchTypeBitset)) {
269 categories.push_back(&matchedHists);
270 }
271 if(InDetSecVtxTruthMatchUtils::isMerged(matchTypeBitset)) {
272 categories.push_back(&mergedHists);
273 }
274 if(InDetSecVtxTruthMatchUtils::isFake(matchTypeBitset)) {
275 categories.push_back(&fakeHists);
276 }
277 if(InDetSecVtxTruthMatchUtils::isSplit(matchTypeBitset)) {
278 categories.push_back(&splitHists);
279 }
280 if(InDetSecVtxTruthMatchUtils::isOther(matchTypeBitset)) {
281 categories.push_back(&otherHists);
282 }
283 categories.push_back(&allHists);
284
285 if (m_doSMOrigin) {
286 const int smOriginTypeBitset = originTypeAcc(*secVtx);
287 hist("RecoVertex/smOriginType")->Fill(smOriginTypeBitset);
288
289 for(const auto& [originType, name] : originTypeMap) {
290 if(InDetSecVtxTruthMatchUtils::isOriginType(smOriginTypeBitset, originType)) {
291 categories.push_back(&m_recoHists.at(name));
292 }
293 }
294 }
295
296 fillRecoHistograms(secVtx, categories);
297 }
298
299 static const xAOD::TruthVertex::ConstAccessor<int> truthTypeAcc("truthVertexMatchType");
300 for(const auto * truthVtx : truthVerticesToMatch) {
301 const int truthTypeBitset = truthTypeAcc(*truthVtx);
302 if(InDetSecVtxTruthMatchUtils::isReconstructable(truthTypeBitset)) {
303 fillTruthHistograms(truthVtx, "Reconstructable");
304
305 // fill efficiencies
306 efficiency("Acceptance")->Fill(InDetSecVtxTruthMatchUtils::isAccepted(truthTypeBitset), truthVtx->perp());
307 efficiency("eff_total")->Fill(InDetSecVtxTruthMatchUtils::isReconstructed(truthTypeBitset), truthVtx->perp());
308 }
309 if(InDetSecVtxTruthMatchUtils::isAccepted(truthTypeBitset)) {
310 fillTruthHistograms(truthVtx, "Accepted");
311 efficiency("eff_seed")->Fill(InDetSecVtxTruthMatchUtils::isSeeded(truthTypeBitset), truthVtx->perp());
312 }
313 if(InDetSecVtxTruthMatchUtils::isSeeded(truthTypeBitset)) {
314 fillTruthHistograms(truthVtx, "Seeded");
315 efficiency("eff_core")->Fill(InDetSecVtxTruthMatchUtils::isReconstructed(truthTypeBitset), truthVtx->perp());
316 }
317 if(InDetSecVtxTruthMatchUtils::isReconstructed(truthTypeBitset)) {
318 fillTruthHistograms(truthVtx, "Reconstructed");
319 }
320 if(InDetSecVtxTruthMatchUtils::isReconstructedSplit(truthTypeBitset)) {
321 fillTruthHistograms(truthVtx, "ReconstructedSplit");
322 }
323 fillTruthHistograms(truthVtx, "Inclusive");
324
325 }
326
327 }
328
329 return StatusCode::SUCCESS;
330
331 }
332 void SecVertexTruthMatchAlg::fillRecoHistograms(const xAOD::Vertex* secVtx, const std::vector<const RecoVertexHists*>& categories) {
333
334 // set of accessors for tracks and truth matching info
335 xAOD::Vertex::ConstAccessor<xAOD::Vertex::TrackParticleLinks_t> trkAcc("trackParticleLinks");
336 const xAOD::Vertex::ConstAccessor<std::vector<InDetSecVtxTruthMatchUtils::VertexTruthMatchInfo> > matchInfoAcc("truthVertexMatchingInfos");
337
338 TVector3 reco_pos(secVtx->x(), secVtx->y(), secVtx->z());
339 float Lxy = reco_pos.Perp();
340
341 size_t ntracks;
342 const xAOD::Vertex::TrackParticleLinks_t & trkParts = trkAcc( *secVtx );
343 ntracks = trkParts.size();
344
345 TLorentzVector sumP4(0,0,0,0);
346 double H = 0.0;
347 double HT = 0.0;
348 int charge = 0;
349 // NOTE: minOpAng/maxOpAng hold the cosines of the minimum/maximum opening
350 // angle between two tracks, i.e. minOpAng is the largest cosine
351 double minOpAng = -1.0* 1.e10;
352 double maxOpAng = 1.0* 1.e10;
353 double minD0 = 1.0* 1.e10;
354 double maxD0 = 0.0;
355 double maxDR = 0.0;
356 size_t nValidTracks = 0;
357
358 ATH_MSG_DEBUG("Loop over tracks");
359 for(size_t t = 0; t < ntracks; t++){
360 if(!trkParts[t].isValid()){
361 ATH_MSG_DEBUG("Track " << t << " is bad!");
362 continue;
363 }
364 const xAOD::TrackParticle & trk = **trkParts[t];
365 ++nValidTracks;
366
367 double trk_d0 = std::abs(trk.definingParameters()[0]);
368 double trk_z0 = std::abs(trk.definingParameters()[1]);
369
370 if(trk_d0 < minD0){ minD0 = trk_d0; }
371 if(trk_d0 > maxD0){ maxD0 = trk_d0; }
372
373 TLorentzVector vv;
374 // TODO: use values computed w.r.t SV
375 vv.SetPtEtaPhiM(trk.pt(),trk.eta(), trk.phi0(), trk.m());
376 sumP4 += vv;
377 H += vv.Vect().Mag();
378 HT += vv.Pt();
379
380 TLorentzVector v_minus_iv(0,0,0,0);
381 for(size_t j = 0; j < ntracks; j++){
382 if (j == t){ continue; }
383 if(!trkParts[j].isValid()){
384 ATH_MSG_DEBUG("Track " << j << " is bad!");
385 continue;
386 }
387
388 const xAOD::TrackParticle & trk_2 = **trkParts[j];
389
390 TLorentzVector tmp;
391 // TODO: use values computed w.r.t. SV
392 tmp.SetPtEtaPhiM(trk_2.pt(),trk_2.eta(), trk_2.phi0(), trk_2.m());
393 v_minus_iv += tmp;
394
395 if( j > t ) {
396 double tm = vv * tmp / ( vv.Mag() * tmp.Mag() );
397 if( minOpAng < tm ) minOpAng = tm;
398 if( maxOpAng > tm ) maxOpAng = tm;
399 }
400 }
401 double DR = vv.DeltaR(v_minus_iv);
402 if( DR > maxDR ){ maxDR = DR;}
403
404 charge += trk.charge();
405 //coverity[UNNECESSARY_STRING_COPY]
406 static const SG::ConstAccessor<float> Trk_Chi2("chiSquared");
407 //coverity[UNNECESSARY_STRING_COPY]
408 static const SG::ConstAccessor<float> Trk_nDoF("numberDoF");
409
410 const bool hasChi2 = Trk_Chi2(trk) >0. && Trk_nDoF(trk) >0.;
411 const auto& covDiag = trk.definingParametersCovMatrixDiagVec();
412
413 for (const RecoVertexHists* h : categories) {
414 if ( hasChi2 ) {
415 h->Trk_Chi2->Fill(Trk_Chi2(trk) / Trk_nDoF(trk));
416 h->Trk_nDoF->Fill(Trk_nDoF(trk));
417 }
418 h->Trk_D0->Fill(trk_d0);
419 h->Trk_Z0->Fill(trk_z0);
420 h->Trk_theta->Fill(trk.definingParameters()[3]);
421 h->Trk_qOverP->Fill(trk.definingParameters()[4]);
422 h->Trk_Eta->Fill(trk.eta());
423 h->Trk_Phi->Fill(trk.phi0());
424 h->Trk_E->Fill(trk.e() / GeV);
425 h->Trk_M->Fill(trk.m() / GeV);
426 h->Trk_Pt->Fill(trk.pt() / GeV);
427 h->Trk_Px->Fill(trk.p4().Px() / GeV);
428 h->Trk_Py->Fill(trk.p4().Py() / GeV);
429 h->Trk_Pz->Fill(trk.p4().Pz() / GeV);
430 h->Trk_charge->Fill(trk.charge());
431 if (covDiag.size() > 1) {
432 h->Trk_errD0->Fill(std::sqrt(covDiag[0]));
433 h->Trk_errZ0->Fill(std::sqrt(covDiag[1]));
434 }
435 }
436
437 } // end loop over tracks
438
439 const double sumP3Mag = sumP4.Vect().Mag();
440 const double recoPosMag = reco_pos.Mag();
441
442 xAOD::Vertex::ConstAccessor<float> Chi2("chiSquared");
443 xAOD::Vertex::ConstAccessor<float> nDoF("numberDoF");
444
445 for (const RecoVertexHists* h : categories) {
446 h->x->Fill(secVtx->x());
447 h->y->Fill(secVtx->y());
448 h->z->Fill(secVtx->z());
449 h->Lxy->Fill(Lxy);
450 h->ntrk->Fill(ntracks);
451 h->pT->Fill(sumP4.Pt() / GeV);
452 h->eta->Fill(sumP4.Eta());
453 h->phi->Fill(sumP4.Phi());
454 h->mass->Fill(sumP4.M() / GeV);
455 if (maxDR > 0) {
456 h->mu->Fill(sumP4.M() / maxDR / GeV);
457 }
458 h->chi2->Fill(Chi2(*secVtx)/nDoF(*secVtx));
459 if (sumP3Mag > 0 && recoPosMag > 0) {
460 h->dir->Fill(sumP4.Vect().Dot( reco_pos ) / sumP3Mag / recoPosMag);
461 }
462 h->charge->Fill(charge);
463 h->H->Fill(H / GeV);
464 h->HT->Fill(HT / GeV);
465 if (nValidTracks > 1) {
466 h->minOpAng->Fill(minOpAng);
467 h->maxOpAng->Fill(maxOpAng);
468 }
469 if (nValidTracks > 0) {
470 h->mind0->Fill(minD0);
471 h->maxd0->Fill(maxD0);
472 }
473 h->maxdR->Fill(maxDR);
474
475 // This includes all matched vertices, including splits
476 if (h->matchScore_weight) {
477 const auto& truthmatchinfo = matchInfoAcc(*secVtx);
478 if(not truthmatchinfo.empty()){
479 float matchScore_weight = std::get<1>(truthmatchinfo.at(0));
480 float matchScore_pt = std::get<2>(truthmatchinfo.at(0));
481
482 ATH_MSG_DEBUG("Match Score and probability: " << matchScore_weight << " " << matchScore_pt/0.01);
483
484 const ElementLink<xAOD::TruthVertexContainer>& truthVertexLink = std::get<0>(truthmatchinfo.at(0));
485 const xAOD::TruthVertex& truthVtx = **truthVertexLink ;
486
487 h->positionRes_R->Fill(Lxy - truthVtx.perp());
488 h->positionRes_Z->Fill(secVtx->z() - truthVtx.z());
489 h->matchScore_weight->Fill(matchScore_weight);
490 h->matchScore_pt->Fill(matchScore_pt);
491 }
492 }
493 }
494 }
495
496 void SecVertexTruthMatchAlg::fillTruthHistograms(const xAOD::TruthVertex* truthVtx, const std::string& truthType) {
497
498 hist("TruthVertex/" + truthType + "_x")->Fill(truthVtx->x());
499 hist("TruthVertex/" + truthType + "_y")->Fill(truthVtx->y());
500 hist("TruthVertex/" + truthType + "_z")->Fill(truthVtx->z());
501 hist("TruthVertex/" + truthType + "_R")->Fill(truthVtx->perp());
502 hist("TruthVertex/" + truthType + "_Eta")->Fill(truthVtx->eta());
503 hist("TruthVertex/" + truthType + "_Phi")->Fill(truthVtx->phi());
504 hist("TruthVertex/" + truthType + "_Ntrk_out")->Fill(truthVtx->nOutgoingParticles());
505
506 ATH_MSG_DEBUG("Plotting truth parent");
507 const xAOD::TruthParticle& truthPart = *truthVtx->incomingParticle(0);
508
509 hist("TruthVertex/" + truthType + "_Parent_E")->Fill(truthPart.e() / GeV);
510 hist("TruthVertex/" + truthType + "_Parent_M")->Fill(truthPart.m() / GeV);
511 hist("TruthVertex/" + truthType + "_Parent_Pt")->Fill(truthPart.pt() / GeV);
512 hist("TruthVertex/" + truthType + "_Parent_Phi")->Fill(truthPart.phi());
513 hist("TruthVertex/" + truthType + "_Parent_Eta")->Fill(truthPart.eta());
514 hist("TruthVertex/" + truthType + "_Parent_charge")->Fill(truthPart.charge());
515
516 ATH_MSG_DEBUG("Plotting truth prod vtx");
517 if(truthPart.hasProdVtx()){
518 const xAOD::TruthVertex & vertex = *truthPart.prodVtx();
519
520 hist("TruthVertex/" + truthType + "_ParentProdX")->Fill(vertex.x());
521 hist("TruthVertex/" + truthType + "_ParentProdY")->Fill(vertex.y());
522 hist("TruthVertex/" + truthType + "_ParentProdZ")->Fill(vertex.z());
523 hist("TruthVertex/" + truthType + "_ParentProdR")->Fill(vertex.perp());
524 hist("TruthVertex/" + truthType + "_ParentProdEta")->Fill(vertex.eta());
525 hist("TruthVertex/" + truthType + "_ParentProdPhi")->Fill(vertex.phi());
526 }
527 }
528
529} // namespace CP
ATH_CHECK
#define ATH_CHECK
Evaluate an expression and check for errors.
Definition AthCheckMacros.h:40
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x,...)
Definition AthMsgStreamMacros.h:43
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x,...)
Definition AthMsgStreamMacros.h:47
ANA_CHECK
#define ANA_CHECK(EXP)
check whether the given expression was successful
Definition Control/AthToolSupport/AsgMessaging/AsgMessaging/MessageCheck.h:325
H
#define H(x, y, z)
Definition MD5.cxx:114
bins
static const std::vector< std::string > bins
Definition MuonTriggerSFFilesTest.cxx:22
h
Header file for AthHistogramAlgorithm.
AthHistogramming::hist
TH1 * hist(std::string_view histName, std::string_view tDir="", std::string_view stream="")
Simplify the retrieval of registered histograms of any type.
Definition AthHistogramming.cxx:182
AthHistogramming::book
StatusCode book(const TH1 &hist, std::string_view tDir="", std::string_view stream="")
Simplify the booking and registering (into THistSvc) of histograms.
Definition AthHistogramming.h:297
AthHistogramming::efficiency
TEfficiency * efficiency(std::string_view effName, std::string_view tDir="", std::string_view stream="")
Simplify the retrieval of registered TEfficiency.
Definition AthHistogramming.cxx:214
CP::SecVertexTruthMatchAlg::m_doMuSA
Gaudi::Property< bool > m_doMuSA
Definition SecVertexTruthMatchAlg.h:61
CP::SecVertexTruthMatchAlg::m_writeHistograms
Gaudi::Property< bool > m_writeHistograms
Definition SecVertexTruthMatchAlg.h:58
CP::SecVertexTruthMatchAlg::m_targetPDGIDs
Gaudi::Property< std::vector< int > > m_targetPDGIDs
Definition SecVertexTruthMatchAlg.h:56
CP::SecVertexTruthMatchAlg::m_doSMOrigin
Gaudi::Property< bool > m_doSMOrigin
Definition SecVertexTruthMatchAlg.h:62
CP::SecVertexTruthMatchAlg::initialize
virtual StatusCode initialize() override
Definition SecVertexTruthMatchAlg.cxx:25
CP::SecVertexTruthMatchAlg::m_trackParticleContainerKey
SG::ReadHandleKey< xAOD::TrackParticleContainer > m_trackParticleContainerKey
Definition SecVertexTruthMatchAlg.h:53
CP::SecVertexTruthMatchAlg::m_secVtxContainerKey
SG::ReadHandleKey< xAOD::VertexContainer > m_secVtxContainerKey
Definition SecVertexTruthMatchAlg.h:47
CP::SecVertexTruthMatchAlg::m_matchTool
ToolHandle< IInDetSecVtxTruthMatchTool > m_matchTool
Definition SecVertexTruthMatchAlg.h:60
CP::SecVertexTruthMatchAlg::SecVertexTruthMatchAlg
SecVertexTruthMatchAlg(const std::string &name, ISvcLocator *svcLoc)
Regular Algorithm constructor.
Definition SecVertexTruthMatchAlg.cxx:22
CP::SecVertexTruthMatchAlg::fillTruthHistograms
void fillTruthHistograms(const xAOD::TruthVertex *truthVtx, const std::string &truthType)
Definition SecVertexTruthMatchAlg.cxx:496
CP::SecVertexTruthMatchAlg::m_truthVtxContainerKey
SG::ReadHandleKey< xAOD::TruthVertexContainer > m_truthVtxContainerKey
Definition SecVertexTruthMatchAlg.h:50
CP::SecVertexTruthMatchAlg::m_recoHists
std::unordered_map< std::string, RecoVertexHists > m_recoHists
Definition SecVertexTruthMatchAlg.h:77
CP::SecVertexTruthMatchAlg::fillRecoHistograms
void fillRecoHistograms(const xAOD::Vertex *secVtx, const std::vector< const RecoVertexHists * > &categories)
Definition SecVertexTruthMatchAlg.cxx:332
EL::AnaAlgorithm::AnaAlgorithm
AnaAlgorithm(const std::string &name, ISvcLocator *pSvcLocator)
constructor with parameters
Definition AnaAlgorithm.cxx:40
EL::AnaAlgorithm::execute
virtual::StatusCode execute()
execute this algorithm
Definition AnaAlgorithm.cxx:315
SG::ConstAccessor
Helper class to provide constant type-safe access to aux data.
Definition ConstAccessor.h:55
SG::ReadHandle::isValid
virtual bool isValid() override final
Can the handle be successfully dereferenced?
SG::ReadHandle::cptr
const_pointer_type cptr()
Dereference the pointer.
xAOD::TrackParticle_v1::definingParametersCovMatrixDiagVec
const std::vector< float > & definingParametersCovMatrixDiagVec() const
Returns the diagonal elements of the defining parameters covariance matrix.
Definition TrackParticle_v1.cxx:326
xAOD::TrackParticle_v1::definingParameters
DefiningParameters_t definingParameters() const
Returns a SVector of the Perigee track parameters.
Definition TrackParticle_v1.cxx:175
xAOD::TrackParticle_v1::m
virtual double m() const override final
The invariant mass of the particle..
Definition TrackParticle_v1.cxx:85
xAOD::TrackParticle_v1::p4
virtual FourMom_t p4() const override final
The full 4-momentum of the particle.
Definition TrackParticle_v1.cxx:122
xAOD::TrackParticle_v1::pt
virtual double pt() const override final
The transverse momentum ( ) of the particle.
Definition TrackParticle_v1.cxx:75
xAOD::TrackParticle_v1::eta
virtual double eta() const override final
The pseudorapidity ( ) of the particle.
Definition TrackParticle_v1.cxx:79
xAOD::TrackParticle_v1::charge
float charge() const
Returns the charge.
Definition TrackParticle_v1.cxx:143
xAOD::TrackParticle_v1::phi0
float phi0() const
Returns the parameter, which has range to .
Definition TrackParticle_v1.cxx:150
xAOD::TrackParticle_v1::e
virtual double e() const override final
The total energy of the particle.
Definition TrackParticle_v1.cxx:111
xAOD::TruthParticle_v1::m
virtual double m() const override final
The mass of the particle.
xAOD::TruthParticle_v1::pdgId
int pdgId() const
PDG ID code.
xAOD::TruthParticle_v1::hasProdVtx
bool hasProdVtx() const
Check for a production vertex on this particle.
Definition TruthParticle_v1.cxx:69
xAOD::TruthParticle_v1::e
virtual double e() const override final
The total energy of the particle.
xAOD::TruthParticle_v1::pt
virtual double pt() const override final
The transverse momentum ( ) of the particle.
Definition TruthParticle_v1.cxx:161
xAOD::TruthParticle_v1::prodVtx
const TruthVertex_v1 * prodVtx() const
The production vertex of this particle.
Definition TruthParticle_v1.cxx:75
xAOD::TruthParticle_v1::charge
double charge() const
Physical charge.
xAOD::TruthParticle_v1::eta
virtual double eta() const override final
The pseudorapidity ( ) of the particle.
Definition TruthParticle_v1.cxx:169
xAOD::TruthParticle_v1::phi
virtual double phi() const override final
The azimuthal angle ( ) of the particle.
Definition TruthParticle_v1.cxx:176
xAOD::TruthVertex_v1::z
float z() const
Vertex longitudinal distance along the beam line form the origin.
xAOD::TruthVertex_v1::eta
float eta() const
Vertex pseudorapidity.
Definition TruthVertex_v1.cxx:171
xAOD::TruthVertex_v1::y
float y() const
Vertex y displacement.
xAOD::TruthVertex_v1::phi
float phi() const
Vertex azimuthal angle.
Definition TruthVertex_v1.cxx:177
xAOD::TruthVertex_v1::incomingParticle
const TruthParticle_v1 * incomingParticle(size_t index) const
Get one of the incoming particles.
Definition TruthVertex_v1.cxx:70
xAOD::TruthVertex_v1::perp
float perp() const
Vertex transverse distance from the beam line.
Definition TruthVertex_v1.cxx:164
xAOD::TruthVertex_v1::nOutgoingParticles
size_t nOutgoingParticles() const
Get the number of outgoing particles.
xAOD::TruthVertex_v1::x
float x() const
Vertex x displacement.
xAOD::Vertex_v1::z
float z() const
Returns the z position.
xAOD::Vertex_v1::TrackParticleLinks_t
std::vector< ElementLink< xAOD::TrackParticleContainer > > TrackParticleLinks_t
Type for the associated track particles.
Definition Vertex_v1.h:128
xAOD::Vertex_v1::y
float y() const
Returns the y position.
xAOD::Vertex_v1::x
float x() const
Returns the x position.
CP
Select isolated Photons, Electrons and Muons.
Definition Control/xAODRootAccess/xAODRootAccess/TEvent.h:27
EL
This module defines the arguments passed from the BATCH driver to the BATCH worker.
Definition AsgComponentFactories.h:16
InDetSecVtxTruthMatchUtils::isOriginType
bool isOriginType(int matchInfo, VertexMatchOriginType type)
Definition InDetSecVtxTruthMatchTool.h:118
xAOD::VxType::SecVtx
@ SecVtx
Secondary vertex.
Definition TrackingPrimitives.h:591
xAOD::TruthVertex
TruthVertex_v1 TruthVertex
Typedef to implementation.
Definition TruthVertex.h:15
xAOD::TrackParticle
TrackParticle_v1 TrackParticle
Reference the current persistent version:
Definition Event/xAOD/xAODTracking/xAODTracking/TrackParticle.h:13
xAOD::Vertex
Vertex_v1 Vertex
Define the latest version of the vertex class.
Definition Event/xAOD/xAODTracking/xAODTracking/Vertex.h:16
xAOD::TruthParticle
TruthParticle_v1 TruthParticle
Typedef to implementation.
Definition Event/xAOD/xAODTruth/xAODTruth/TruthParticle.h:15
CP::SecVertexTruthMatchAlg::RecoVertexHists
cached pointers to the histograms of one reco vertex category
Definition SecVertexTruthMatchAlg.h:65
Generated on for ATLAS Offline Software by doxygen 1.17.0