Project
Loading...
Searching...
No Matches
MatchITSTPCQC.cxx
Go to the documentation of this file.
1// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
2// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
3// All rights not expressly granted are reserved.
4//
5// This software is distributed under the terms of the GNU General Public
6// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
7//
8// In applying this license CERN does not waive the privileges and immunities
9// granted to it by virtue of its status as an Intergovernmental Organization
10// or submit itself to any jurisdiction.
11
12#include <algorithm>
13#include <map>
14#include <unordered_map>
15
16#include "GLOQC/MatchITSTPCQC.h"
17#include "ReconstructionDataFormats/TrackTPCITS.h"
18#include "DataFormatsTPC/TrackTPC.h"
19#include "DetectorsBase/Propagator.h"
20#include "SimulationDataFormat/MCUtils.h"
21#include "GlobalTracking/TrackCuts.h"
22#include <DetectorsBase/GRPGeomHelper.h>
23#include <TEfficiency.h>
24#include "ReconstructionDataFormats/PrimaryVertex.h"
25#include "ReconstructionDataFormats/V0.h"
26#include "DetectorsVertexing/SVertexerParams.h"
27#include "Framework/InputRecord.h"
28#include "Framework/TimingInfo.h"
29#include "GPUO2InterfaceUtils.h"
30#include "CommonConstants/LHCConstants.h"
31#include "DetectorsCommonDataFormats/DetID.h"
32#include "GPUO2InterfaceRefit.h"
33
34using namespace o2::gloqc;
35using namespace o2::mcutils;
36using MCTrack = o2::MCTrackT<float>;
37using DetID = o2::detectors::DetID;
38
43
44//_______________________________________________________
45
46void MatchITSTPCQC::deleteHistograms()
47{
48
49 LOG(debug) << "Deleting histos...";
50 for (int i = 0; i < matchType::SIZE; ++i) {
51 // Pt
52 delete mPtNum[i];
53 delete mPtDen[i];
54 delete mFractionITSTPCmatch[i];
55 delete mPtNum_noEta0[i];
56 delete mPtDen_noEta0[i];
57 delete mFractionITSTPCmatch_noEta0[i];
58 delete mPtPhysPrimNum[i];
59 delete mPtPhysPrimDen[i];
60 delete mFractionITSTPCmatchPhysPrim[i];
61
62 // Phi
63 delete mPhiNum[i];
64 delete mPhiDen[i];
65 delete mFractionITSTPCmatchPhi[i];
66 delete mPhiPhysPrimNum[i];
67 delete mPhiPhysPrimDen[i];
68 delete mFractionITSTPCmatchPhiPhysPrim[i];
69 delete mPhiVsPtNum[i];
70 delete mPhiVsPtDen[i];
71 delete mFractionITSTPCmatchPhiVsPt[i];
72
73 // Eta
74 delete mEtaNum[i];
75 delete mEtaDen[i];
76 delete mFractionITSTPCmatchEta[i];
77 delete mEtaPhysPrimNum[i];
78 delete mEtaPhysPrimDen[i];
79 delete mFractionITSTPCmatchEtaPhysPrim[i];
80 delete mEtaVsPtNum[i];
81 delete mEtaVsPtDen[i];
82 delete mFractionITSTPCmatchEtaVsPt[i];
83
84 // Clusters
85 delete mClsVsPtNum[i];
86 delete mClsVsPtDen[i];
87 delete mFractionITSTPCmatchClsVsPt[i];
88
89 // Chi2
90 delete mChi2VsPtNum[i];
91 delete mChi2VsPtDen[i];
92 delete mFractionITSTPCmatchChi2VsPt[i];
93
94 for (int j = 0; j < o2::track::PID::NIDs; ++j) {
95 // Pt
96 delete mPtNumVsTrkPID[i][j];
97 delete mPtDenVsTrkPID[i][j];
98 delete mFractionITSTPCmatchPtVsTrkPID[i][j];
99 // Phi
100 delete mPhiNumVsTrkPID[i][j];
101 delete mPhiDenVsTrkPID[i][j];
102 delete mFractionITSTPCmatchPhiVsTrkPID[i][j];
103 // Eta
104 delete mEtaNumVsTrkPID[i][j];
105 delete mEtaDenVsTrkPID[i][j];
106 delete mFractionITSTPCmatchEtaVsTrkPID[i][j];
107 }
108
109 // 1/Pt
110 delete m1OverPtNum[i];
111 delete m1OverPtDen[i];
112 delete mFractionITSTPCmatch1OverPt[i];
113 delete m1OverPtPhysPrimNum[i];
114 delete m1OverPtPhysPrimDen[i];
115 delete mFractionITSTPCmatchPhysPrim1OverPt[i];
116
117 // 3D eta/phi/pt
118 delete mEtaPhiPtNum[i];
119 delete mEtaPhiPtDen[i];
120 }
121
122 // Residuals
123 delete mResidualPt;
124 delete mResidualPhi;
125 delete mResidualEta;
126 // Others
127 delete mChi2Matching;
128 delete mChi2Refit;
129 delete mTimeResVsPt;
130 delete mDCAr;
131 delete mDCArVsPtNum;
132 delete mDCArVsPtDen;
133 delete mFractionITSTPCmatchDCArVsPt;
134
135 // K0
136 delete mK0MassVsPtVsOccpp;
137 delete mK0MassVsPtVsOccPbPb;
138}
139
140//__________________________________________________________
141
142void MatchITSTPCQC::reset()
143{
144 for (int i = 0; i < matchType::SIZE; ++i) {
145 // Pt
146 mPtNum[i]->Reset();
147 mPtDen[i]->Reset();
148 mPtNum_noEta0[i]->Reset();
149 mPtDen_noEta0[i]->Reset();
150
151 // Phi
152 mPhiNum[i]->Reset();
153 mPhiDen[i]->Reset();
154 mPhiVsPtNum[i]->Reset();
155 mPhiVsPtDen[i]->Reset();
156
157 // Eta
158 mEtaNum[i]->Reset();
159 mEtaDen[i]->Reset();
160 mEtaVsPtNum[i]->Reset();
161 mEtaVsPtDen[i]->Reset();
162
163 // Clusters
164 mClsVsPtNum[i]->Reset();
165 mClsVsPtDen[i]->Reset();
166
167 // Chi2
168 mChi2VsPtNum[i]->Reset();
169 mChi2VsPtDen[i]->Reset();
170
171 // 1/Pt
172 m1OverPtNum[i]->Reset();
173 m1OverPtDen[i]->Reset();
174
175 // 3D eta/phi/pt
176 mEtaPhiPtNum[i]->Reset();
177 mEtaPhiPtDen[i]->Reset();
178
179 if (mUseTrkPID) { // Vs Tracking PID hypothesis
180 for (int j = 0; j < o2::track::PID::NIDs; ++j) {
181 // Pt
182 mPtNumVsTrkPID[i][j]->Reset();
183 mPtDenVsTrkPID[i][j]->Reset();
184 // Phi
185 mPhiNumVsTrkPID[i][j]->Reset();
186 mPhiDenVsTrkPID[i][j]->Reset();
187 // Eta
188 mEtaNumVsTrkPID[i][j]->Reset();
189 mEtaDenVsTrkPID[i][j]->Reset();
190 }
191 }
192
193 if (mUseMC) {
194 mPtPhysPrimNum[i]->Reset();
195 mPtPhysPrimDen[i]->Reset();
196
197 mPhiPhysPrimNum[i]->Reset();
198 mPhiPhysPrimDen[i]->Reset();
199
200 mEtaPhysPrimNum[i]->Reset();
201 mEtaPhysPrimDen[i]->Reset();
202
203 m1OverPtPhysPrimNum[i]->Reset();
204 m1OverPtPhysPrimDen[i]->Reset();
205 }
206 }
207
208 // Residuals
209 mResidualPt->Reset();
210 mResidualPhi->Reset();
211 mResidualEta->Reset();
212 // Others
213 mChi2Matching->Reset();
214 mChi2Refit->Reset();
215 mTimeResVsPt->Reset();
216 mDCAr->Reset();
217 mDCArVsPtNum->Reset();
218 mDCArVsPtDen->Reset();
219
220 // K0
221 if (mDoK0QC) {
222 mK0MassVsPtVsOccpp->Reset();
223 mK0MassVsPtVsOccPbPb->Reset();
224 }
225}
226
227//__________________________________________________________
228bool MatchITSTPCQC::init()
229{
230 LOGP(debug, "Creating Variable Binning");
231 std::array<std::string, 2> title{"TPC", "ITS"};
232 std::array<std::string, 2> etaSel{Form(", |#eta| < %.1f", mEtaTPCCut), Form(", |#eta| < %.1f", mEtaCut)};
233 std::array<int, 2> maxNCls{156, 7};
234 // log binning for pT
235 const Int_t nbinsPt = mPtBins;
236 const Double_t xminPt = (mPtCut > 0) ? mPtCut : 0.01;
237 const Double_t xmaxPt = mPtMaxCut;
238 Double_t* xbinsPt = new Double_t[nbinsPt + 1];
239 Double_t xlogminPt = TMath::Log10(xminPt);
240 Double_t xlogmaxPt = TMath::Log10(xmaxPt);
241 Double_t dlogxPt = (xlogmaxPt - xlogminPt) / nbinsPt;
242 for (int i = 0; i <= nbinsPt; i++) {
243 Double_t xlogPt = xlogminPt + i * dlogxPt;
244 xbinsPt[i] = TMath::Exp(TMath::Log(10) * xlogPt);
245 }
246
247 LOGP(debug, "Creating Histograms");
248 // Data and MC
249 for (int i = 0; i < matchType::SIZE; ++i) {
250 // Pt
251 mPtNum[i] = new TH1D(Form("mPtNum_%s", title[i].c_str()), Form("Pt distribution of ITSTPC matched tracks, wrt %s tracks %s; Pt [GeV/c]; dNdPt", title[i].c_str(), etaSel[i].c_str()), mPtBins, mPtCut, mPtMaxCut);
252 mPtNum[i]->Sumw2();
253 mPtNum[i]->SetOption("logy");
254 mPtNum[i]->GetYaxis()->SetTitleOffset(1.4);
255 mPtDen[i] = new TH1D(Form("mPtDen_%s", title[i].c_str()), Form("Pt distribution of %s tracks %s; Pt [GeV/c]; dNdPt", title[i].c_str(), etaSel[i].c_str()), mPtBins, mPtCut, mPtMaxCut);
256 mPtDen[i]->Sumw2();
257 mPtDen[i]->SetOption("logy");
258 mPtDen[i]->GetYaxis()->SetTitleOffset(1.4);
259 mFractionITSTPCmatch[i] = new TEfficiency(Form("mFractionITSTPCmatch_%s", title[i].c_str()), Form("Fraction of ITSTPC matched tracks wrt %s tracks vs Pt %s; Pt [GeV/c]; Eff", title[i].c_str(), etaSel[i].c_str()), mPtBins, mPtCut, mPtMaxCut);
260 mPtNum_noEta0[i] = new TH1D(Form("mPtNum_noEta0_%s", title[i].c_str()), Form("Pt distribution of ITSTPC matched tracks without |eta| < %.2f, wrt %s tracks %s; Pt [GeV/c]; dNdPt", mEtaNo0Cut, title[i].c_str(), etaSel[i].c_str()), mPtBins, mPtCut, mPtMaxCut);
261 mPtNum_noEta0[i]->Sumw2();
262 mPtNum_noEta0[i]->SetOption("logy");
263 mPtNum_noEta0[i]->GetYaxis()->SetTitleOffset(1.4);
264 mPtDen_noEta0[i] = new TH1D(Form("mPtDen_noEta0_%s", title[i].c_str()), Form("Pt distribution of %s tracks without |eta| < %.2f %s; Pt [GeV/c]; dNdPt", title[i].c_str(), mEtaNo0Cut, etaSel[i].c_str()), mPtBins, mPtCut, mPtMaxCut);
265 mPtDen_noEta0[i]->Sumw2();
266 mPtDen_noEta0[i]->SetOption("logy");
267 mPtDen_noEta0[i]->GetYaxis()->SetTitleOffset(1.4);
268 mFractionITSTPCmatch_noEta0[i] = new TEfficiency(Form("mFractionITSTPCmatch_noEta0_%s", title[i].c_str()), Form("Fraction of ITSTPC matched tracks wrt %s tracks vs Pt without |eta| < %.2f %s; Pt [GeV/c]; Eff", title[i].c_str(), mEtaNo0Cut, etaSel[i].c_str()), mPtBins, mPtCut, mPtMaxCut);
269
270 // Phi
271 mPhiNum[i] = new TH1F(Form("mPhiNum_%s", title[i].c_str()), Form("Phi distribution of ITSTPC matched tracks, wrt %s tracks %s; Phi [rad]; dNdPhi", title[i].c_str(), etaSel[i].c_str()), 100, 0.f, 2 * TMath::Pi());
272 mPhiNum[i]->Sumw2();
273 mPhiDen[i] = new TH1F(Form("mPhiDen_%s", title[i].c_str()), Form("Phi distribution of %s tracks %s; Phi [rad]; dNdPhi", title[i].c_str(), etaSel[i].c_str()), 100, 0.f, 2 * TMath::Pi());
274 mPhiDen[i]->Sumw2();
275 mFractionITSTPCmatchPhi[i] = new TEfficiency(Form("mFractionITSTPCmatchPhi_%s", title[i].c_str()), Form("Fraction of ITSTPC matched tracks vs Phi wrt %s tracks %s; Phi [rad]; Eff", title[i].c_str(), etaSel[i].c_str()), 100, 0.f, 2 * TMath::Pi());
276 mPhiVsPtNum[i] = new TH2F(Form("mPhiVsPtNum_%s", title[i].c_str()), Form("Phi vs Pt distribution of ITSTPC matched tracks wrt %s %s; #it{p}_{T} [GeV#it{c}]; Phi [rad]; dNdPhi", title[i].c_str(), etaSel[i].c_str()), mPtBins, mPtCut, mPtMaxCut, 100, 0.f, 2 * TMath::Pi());
277 mPhiVsPtNum[i]->Sumw2();
278 mPhiVsPtDen[i] = new TH2F(Form("mPhiVsPtDen_%s", title[i].c_str()), Form("Phi vs Pt distribution of %s tracks %s; #it{p}_{T} [GeV#it{c}]; Phi [rad]; dNdPhi", title[i].c_str(), etaSel[i].c_str()), mPtBins, mPtCut, mPtMaxCut, 100, 0.f, 2 * TMath::Pi());
279 mPhiVsPtDen[i]->Sumw2();
280 mFractionITSTPCmatchPhiVsPt[i] = new TEfficiency(Form("mFractionITSTPCmatchPhiVsPt_%s", title[i].c_str()), Form("Fraction of ITSTPC matched tracks wrt %s tracks %s, Phi vs Pt; #it{p}_{T} [GeV#it{c}]; Phi [rad]; Eff", title[i].c_str(), etaSel[i].c_str()), mPtBins, mPtCut, mPtMaxCut, 100, 0.f, 2 * TMath::Pi());
281
282 // Eta
283 mEtaNum[i] = new TH1F(Form("mEtaNum_%s", title[i].c_str()), Form("Eta distribution of ITSTPC matched tracks, wrt %s tracks; Eta; dNdEta", title[i].c_str()), 100, -2.f, 2.f);
284 mEtaNum[i]->Sumw2();
285 mEtaNum[i]->GetYaxis()->SetTitleOffset(1.4);
286 mEtaDen[i] = new TH1F(Form("mEtaDen_%s", title[i].c_str()), Form("Eta distribution of %s tracks; Eta; dNdEta", title[i].c_str()), 100, -2.f, 2.f);
287 mEtaDen[i]->Sumw2();
288 mEtaDen[i]->GetYaxis()->SetTitleOffset(1.4);
289 mFractionITSTPCmatchEta[i] = new TEfficiency(Form("mFractionITSTPCmatchEta_%s", title[i].c_str()), Form("Fraction of ITSTPC matched tracks , wrt %s tracks, vs Eta; Eta; Eff", title[i].c_str()), 100, -2.f, 2.f);
290 mEtaVsPtNum[i] = new TH2F(Form("mEtaVsPtNum_%s", title[i].c_str()), Form("Eta vs Pt distribution of ITSTPC matched tracks, wrt %s tracks; #it{p}_{T} [GeV#it{c}]; Eta", title[i].c_str()), mPtBins, mPtCut, mPtMaxCut, 100, -2.f, 2.f);
291 mEtaVsPtNum[i]->Sumw2();
292 mEtaVsPtDen[i] = new TH2F(Form("mEtaVsPtDen_%s", title[i].c_str()), Form("Eta vs Pt distribution of %s tracks; #it{p}_{T} [GeV#it{c}]; Eta", title[i].c_str()), mPtBins, mPtCut, mPtMaxCut, 100, -2.f, 2.f);
293 mEtaVsPtDen[i]->Sumw2();
294 mFractionITSTPCmatchEtaVsPt[i] = new TEfficiency(Form("mFractionITSTPCmatchEtaVsPt_%s", title[i].c_str()), Form("Fraction of ITSTPC matched tracks, wrt %s tracks, Eta vs Pt; #it{p}_{T} [GeV#it{c}]; Eta; Eff", title[i].c_str()), mPtBins, mPtCut, mPtMaxCut, 100, -2.f, 2.f);
295
296 // Clusters
297 mClsVsPtNum[i] = new TH2F(Form("mClsVsPtNum_%s", title[i].c_str()), Form("#Clusters vs Pt distribution of ITSTPC matched tracks, wrt %s tracks; #it{p}_{T} [GeV#it{c}]; #Clusters", title[i].c_str()), mPtBins, mPtCut, mPtMaxCut, maxNCls[i], 0, maxNCls[i]);
298 mClsVsPtNum[i]->Sumw2();
299 mClsVsPtDen[i] = new TH2F(Form("mClsVsPtDen_%s", title[i].c_str()), Form("#Clusters vs Pt distribution of %s tracks; #it{p}_{T} [GeV#it{c}]; #Clusters", title[i].c_str()), mPtBins, mPtCut, mPtMaxCut, maxNCls[i], 0, maxNCls[i]);
300 mClsVsPtDen[i]->Sumw2();
301 mFractionITSTPCmatchClsVsPt[i] = new TEfficiency(Form("mFractionITSTPCmatchClsVsPt_%s", title[i].c_str()), Form("Fraction of ITSTPC matched tracks, wrt %s tracks, #Clusters vs Pt; #it{p}_{T} [GeV#it{c}]; #Clusters; Eff", title[i].c_str()), mPtBins, mPtCut, mPtMaxCut, maxNCls[i], 0, maxNCls[i]);
302
303 // Chi2
304 mChi2VsPtNum[i] = new TH2F(Form("mChi2VsPtNum_%s", title[i].c_str()), Form("Chi2 vs Pt distribution of ITSTPC matched tracks, wrt %s tracks; #it{p}_{T} [GeV#it{c}]; Chi2", title[i].c_str()), mPtBins, mPtCut, mPtMaxCut, 200, 0, 300);
305 mChi2VsPtNum[i]->Sumw2();
306 mChi2VsPtDen[i] = new TH2F(Form("mChi2VsPtDen_%s", title[i].c_str()), Form("Chi2 vs Pt distribution of %s tracks; #it{p}_{T} [GeV#it{c}]; Chi2", title[i].c_str()), mPtBins, mPtCut, mPtMaxCut, 200, 0, 300);
307 mChi2VsPtDen[i]->Sumw2();
308 mFractionITSTPCmatchChi2VsPt[i] = new TEfficiency(Form("mFractionITSTPCmatchChi2VsPt_%s", title[i].c_str()), Form("Fraction of ITSTPC matched tracks, wrt %s tracks, Chi2 vs Pt; #it{p}_{T} [GeV#it{c}]; Chi2; Eff", title[i].c_str()), mPtBins, mPtCut, mPtMaxCut, 200, 0, 300);
309
310 // 1/pt
311 m1OverPtNum[i] = new TH1D(Form("m1OverPtNum_%s", title[i].c_str()), Form("1/Pt distribution of matched tracks, wrt %s tracks %s; 1/Pt [c/GeV]; dNdPt", title[i].c_str(), etaSel[i].c_str()), 100, -20.f, 20.f);
312 m1OverPtNum[i]->Sumw2();
313 m1OverPtDen[i] = new TH1D(Form("m1OverPtDen_%s", title[i].c_str()), Form("1/Pt distribution of %s tracks %s; 1/Pt [c/GeV]; dNdPt", title[i].c_str(), etaSel[i].c_str()), 100, -20.f, 20.f);
314 m1OverPtDen[i]->Sumw2();
315 mFractionITSTPCmatch1OverPt[i] = new TEfficiency(Form("mFractionITSTPCmatch1OverPt_%s", title[i].c_str()), Form("Fraction of ITSTPC matched tracks vs 1/Pt, wrt %s tracks %s; 1/Pt [c/GeV]; Eff", title[i].c_str(), etaSel[i].c_str()), 100, -20.f, 20.f);
316
317 // 3d eta/phi/pt
318 mEtaPhiPtNum[i] = new TH3F(Form("mEtaPhiPtNum_%s", title[i].c_str()), Form("Numerator #eta vs #varphi vs #it{p}_{T}, wrt %s;#eta %s;#varphi;#it{p}_{T} [GeV#it{c}];Entries", title[i].c_str(), etaSel[i].c_str()), 100, -2., 2., 100, 0., 2 * TMath::Pi(), 100, 0.01, 20.);
319 mEtaPhiPtNum[i]->Sumw2();
320 mEtaPhiPtDen[i] = new TH3F(Form("mEtaPhiPtDen_%s", title[i].c_str()), Form("Denominator #eta vs #varphi vs #it{p}_{T}, wrt %s;#eta %s;#varphi;#it{p}_{T} [GeV#it{c}];Entries", title[i].c_str(), etaSel[i].c_str()), 100, -2., 2., 100, 0., 2 * TMath::Pi(), 100, 0.01, 20.);
321 mEtaPhiPtDen[i]->Sumw2();
322
323 if (mUseTrkPID) { // Vs Tracking PID hypothesis
324 for (int j = 0; j < o2::track::PID::NIDs; ++j) {
325 // Pt
326 mPtNumVsTrkPID[i][j] = new TH1D(Form("mPtNumVsTrkPID_%s_PID%i", title[i].c_str(), j), Form("Pt distribution of ITSTPC matched tracks, wrt %s tracks %s, TrkPID %i; Pt [GeV/c]; dNdPt", title[i].c_str(), etaSel[i].c_str(), j), mPtBins, mPtCut, mPtMaxCut);
327 mPtNumVsTrkPID[i][j]->Sumw2();
328 mPtDenVsTrkPID[i][j] = new TH1D(Form("mPtDenVsTrkPID_%s_PID%i", title[i].c_str(), j), Form("Pt distribution of %s tracks %s, TrkPID %i; Pt [GeV/c]; dNdPt", title[i].c_str(), etaSel[i].c_str(), j), mPtBins, mPtCut, mPtMaxCut);
329 mPtDenVsTrkPID[i][j]->Sumw2();
330 mFractionITSTPCmatchPtVsTrkPID[i][j] = new TEfficiency(Form("mFractionITSTPCmatchPtVsTrkPID_%s_PID%i", title[i].c_str(), j), Form("Fraction of ITSTPC matched tracks wrt %s tracks vs Pt %s, TrkPID %i; Pt [GeV/c]; Eff", title[i].c_str(), etaSel[i].c_str(), j), mPtBins, mPtCut, mPtMaxCut);
331
332 // Phi
333 mPhiNumVsTrkPID[i][j] = new TH1D(Form("mPhiNumVsTrkPID_%s_PID%i", title[i].c_str(), j), Form("Phi distribution of ITSTPC matched tracks, wrt %s tracks %s, TrkPID %i; Phi [rad]; dNdPhi", title[i].c_str(), etaSel[i].c_str(), j), 100, 0.f, 2 * TMath::Pi());
334 mPhiNumVsTrkPID[i][j]->Sumw2();
335 mPhiDenVsTrkPID[i][j] = new TH1D(Form("mPhiDenVsTrkPID_%s_PID%i", title[i].c_str(), j), Form("Phi distribution of %s tracks %s, TrkPID %i; Phi [rad]; dNdPhi", title[i].c_str(), etaSel[i].c_str(), j), 100, 0.f, 2 * TMath::Pi());
336 mPhiDenVsTrkPID[i][j]->Sumw2();
337 mFractionITSTPCmatchPhiVsTrkPID[i][j] = new TEfficiency(Form("mFractionITSTPCmatchPhiVsTrkPID_%s_PID%i", title[i].c_str(), j), Form("Fraction of ITSTPC matched tracks wrt %s tracks vs Phi %s, TrkPID %i; Phi [rad]; Eff", title[i].c_str(), etaSel[i].c_str(), j), 100, 0.f, 2 * TMath::Pi());
338
339 // Eta
340 mEtaNumVsTrkPID[i][j] = new TH1D(Form("mEtaNumVsTrkPID_%s_PID%i", title[i].c_str(), j), Form("Eta distribution of ITSTPC matched tracks, wrt %s tracks %s, TrkPID %i; Eta; dNdEta", title[i].c_str(), etaSel[i].c_str(), j), 100, -2.f, 2.f);
341 mEtaNumVsTrkPID[i][j]->Sumw2();
342 mEtaDenVsTrkPID[i][j] = new TH1D(Form("mEtaDenVsTrkPID_%s_PID%i", title[i].c_str(), j), Form("Eta distribution of %s tracks %s, TrkPID %i; Eta; dNdEta", title[i].c_str(), etaSel[i].c_str(), j), 100, -2.f, 2.f);
343 mEtaDenVsTrkPID[i][j]->Sumw2();
344 mFractionITSTPCmatchEtaVsTrkPID[i][j] = new TEfficiency(Form("mFractionITSTPCmatchEtaVsTrkPID_%s_PID%i", title[i].c_str(), j), Form("Fraction of ITSTPC matched tracks wrt %s tracks vs Eta %s, TrkPID %i; Eta; Eff", title[i].c_str(), etaSel[i].c_str(), j), 100, -2.f, 2.f);
345 }
346 }
347 }
348
349 mResidualPt = new TH2F("mResidualPt", "Residuals of ITS-TPC matching in #it{p}_{T}; #it{p}_{T}^{ITS-TPC} [GeV/c]; #it{p}_{T}^{ITS-TPC} - #it{p}_{T}^{TPC} [GeV/c]", mPtBins, mPtCut, mPtMaxCut, 100, -1.f, 1.f);
350 mResidualPhi = new TH2F("mResidualPhi", "Residuals of ITS-TPC matching in #it{#phi}; #it{#phi}^{ITS-TPC} [rad]; #it{#phi}^{ITS-TPC} - #it{#phi}^{TPC} [rad]", 100, 0.f, 2 * TMath::Pi(), 100, -1.f, 1.f);
351 mResidualEta = new TH2F("mResidualEta", "Residuals of ITS-TPC matching in #it{#eta}; #it{#eta}^{ITS-TPC}; #it{#eta}^{ITS-TPC} - #it{#eta}^{TPC}", 100, -2.f, 2.f, 100, -1.f, 1.f);
352 mChi2Matching = new TH1F("mChi2Matching", "Chi2 of matching; chi2", 200, 0, 300);
353 mChi2Matching->SetOption("logy");
354 mChi2Matching->GetYaxis()->SetTitleOffset(1.4);
355 mChi2Refit = new TH1F("mChi2Refit", "Chi2 of refit; chi2", 200, 0, 300);
356 mChi2Refit->SetOption("logy");
357 mChi2Refit->GetYaxis()->SetTitleOffset(1.4);
358 mDCAr = new TH1F("mDCAr", "DCA of TPC tracks; DCAr", 100, -mDCATPCCutY, mDCATPCCutY);
359 mDCArVsPtNum = new TH2F("mDCArVsPtNum", "DCA of TPC tracks Vs Pt Num; #it{p}_{T} [GeV/c]; DCAr", 100, 0, 20., 100, -mDCATPCCutY, mDCATPCCutY);
360 mDCArVsPtNum->Sumw2();
361 mDCArVsPtDen = new TH2F("mDCArVsPtDen", "DCA of TPC tracks Vs Pt Den; #it{p}_{T} [GeV/c]; DCAr", 100, 0, 20., 100, -mDCATPCCutY, mDCATPCCutY);
362 mDCArVsPtDen->Sumw2();
363 mFractionITSTPCmatchDCArVsPt = new TEfficiency("mFractionITSTPCmatchDCArVsPt", "Fraction of ITSTPC matched tracks wrt TPC vs DCAr; #it{p}_{T} [GeV#it{c}]; DCAr; Eff", 100, 0, 20., 200, -30, 30);
364
365 mTimeResVsPt = new TH2F("mTimeResVsPt", "Time resolution vs Pt; Pt [GeV/c]; time res [us]", nbinsPt, xbinsPt, 100, 0.f, 2.f);
366 mTimeResVsPt->SetOption("colz logz logy logx");
367 mTimeResVsPt->GetYaxis()->SetTitleOffset(1.4);
368
369 if (mUseMC) {
370 mcReader.initFromDigitContext("collisioncontext.root");
371
372 for (int i = 0; i < matchType::SIZE; ++i) {
373 mPtPhysPrimNum[i] = new TH1F(Form("mPtPhysPrimNum_%s", title[i].c_str()), Form("Pt distribution of matched tracks (physical primary), wrt %s tracks %s; Pt [GeV/c]; dNdPt", title[i].c_str(), etaSel[i].c_str()), nbinsPt, xbinsPt);
374 mPtPhysPrimNum[i]->Sumw2();
375 mPtPhysPrimDen[i] = new TH1F(Form("mPtPhysPrimDen_%s", title[i].c_str()), Form("Pt distribution of %s tracks (physical primary) %s; Pt [GeV/c]; dNdPt", title[i].c_str(), etaSel[i].c_str()), nbinsPt, xbinsPt);
376 mPtPhysPrimDen[i]->Sumw2();
377 mFractionITSTPCmatchPhiPhysPrim[i] = new TEfficiency(Form("mFractionITSTPCmatchPhiPhysPrim_%s", title[i].c_str()), Form("Fraction of ITSTPC matched tracks vs Phi (physical primary), wrt %s tracks %s; Phi [rad]; Eff", title[i].c_str(), etaSel[i].c_str()), 100, 0.f, 2 * TMath::Pi());
378
379 mEtaPhysPrimNum[i] = new TH1F(Form("mEtaPhysPrimNum_%s", title[i].c_str()), Form("Eta distribution of matched tracks (physical primary), wrt %s tracks; Eta; dNdEta", title[i].c_str()), 100, -2.f, 2.f);
380 mEtaPhysPrimNum[i]->Sumw2();
381 mEtaPhysPrimDen[i] = new TH1F(Form("mEtaPhysPrimDen_%s", title[i].c_str()), Form("Eta distribution of %s tracks (physical primary); Eta; dNdEta", title[i].c_str()), 100, -2.f, 2.f);
382 mEtaPhysPrimDen[i]->Sumw2();
383 mFractionITSTPCmatchEtaPhysPrim[i] = new TEfficiency(Form("mFractionITSTPCmatchEtaPhysPrim_%s", title[i].c_str()), Form("Fraction of ITSTPC matched tracks vs Eta (physical primary), wrt %s tracks; Eta; Eff", title[i].c_str()), 100, -2.f, 2.f);
384
385 mPhiPhysPrimNum[i] = new TH1F(Form("mPhiPhysPrimNum_%s", title[i].c_str()), Form("Phi distribution of matched tracks (physical primary), wrt %s tracks %s; Phi [rad]; dNdPhi", title[i].c_str(), etaSel[i].c_str()), 100, 0.f, 2 * TMath::Pi());
386 mPhiPhysPrimNum[i]->Sumw2();
387 mPhiPhysPrimDen[i] = new TH1F(Form("mPhiPhysPrimDen_%s", title[i].c_str()), Form("Phi distribution of %s tracks (physical primary) %s; Phi [rad]; dNdPhi", title[i].c_str(), etaSel[i].c_str()), 100, 0.f, 2 * TMath::Pi());
388 mPhiPhysPrimDen[i]->Sumw2();
389 mFractionITSTPCmatchPhysPrim[i] = new TEfficiency(Form("mFractionITSTPCmatchPhysPrim_%s", title[i].c_str()), Form("Fraction of ITSTPC matched tracks vs Pt (physical primary), wrt %s tracks %s; Pt [GeV/c]; Eff", title[i].c_str(), etaSel[i].c_str()), nbinsPt, xbinsPt);
390
391 m1OverPtPhysPrimNum[i] = new TH1D(Form("m1OverPtPhysPrimNum_%s", title[i].c_str()), Form("1/Pt distribution of matched tracks (physical primary), wrt %s tracks %s; 1/Pt [c/GeV]; dNd1/Pt", title[i].c_str(), etaSel[i].c_str()), 2 * mPtBins, -20., 20.);
392 m1OverPtPhysPrimNum[i]->Sumw2();
393 m1OverPtPhysPrimDen[i] = new TH1D(Form("m1OverPtPhysPrimDen_%s", title[i].c_str()), Form("1/PtPt distribution of %s tracks (physical primary) %s; 1/Pt [c/GeV]; dNd1/Pt", title[i].c_str(), etaSel[i].c_str()), 2 * mPtBins, -20., 20.);
394 m1OverPtPhysPrimDen[i]->Sumw2();
395 mFractionITSTPCmatchPhysPrim1OverPt[i] = new TEfficiency(Form("mFractionITSTPCmatchPhysPrim1OverPt_%s", title[i].c_str()), Form("Fraction of ITSTPC matched tracks vs 1/Pt (physical primary), wrt %s tracks %s; 1/Pt [c/GeV]; Eff", title[i].c_str(), etaSel[i].c_str()), 2 * mPtBins, -20., 20.);
396 }
397 }
398
399 // log binning for pT for K0s
400 const Int_t nbinsPtK0 = 10;
401 const Double_t xminPtK0 = 0.01;
402 const Double_t xmaxPtK0 = 20;
403 Double_t* xbinsPtK0 = new Double_t[nbinsPtK0 + 1];
404 Double_t xlogminPtK0 = TMath::Log10(xminPtK0);
405 Double_t xlogmaxPtK0 = TMath::Log10(xmaxPtK0);
406 Double_t dlogxPtK0 = (xlogmaxPtK0 - xlogminPtK0) / nbinsPtK0;
407 for (int i = 0; i <= nbinsPtK0; i++) {
408 Double_t xlogPtK0 = xlogminPtK0 + i * dlogxPtK0;
409 xbinsPtK0[i] = TMath::Exp(TMath::Log(10) * xlogPtK0);
410 }
411 // the other bins
412 const Int_t nbinsMassK0 = 100;
413 Double_t* ybinsMassK0 = new Double_t[nbinsMassK0 + 1];
414 Double_t yminMassK0 = 0.4;
415 Double_t ymaxMassK0 = 0.6;
416 Double_t dyMassK0 = (ymaxMassK0 - yminMassK0) / nbinsMassK0;
417 for (int i = 0; i <= nbinsMassK0; i++) {
418 ybinsMassK0[i] = yminMassK0 + i * dyMassK0;
419 }
420 const Int_t nbinsMultK0pp = mNBinsTPCOccpp;
421 Double_t* zbinsMultK0pp = new Double_t[nbinsMultK0pp + 1];
422 Double_t zminMultK0pp = mMinTPCOccpp;
423 Double_t zmaxMultK0pp = mMaxTPCOccpp;
424 Double_t dzMultK0pp = (zmaxMultK0pp - zminMultK0pp) / nbinsMultK0pp;
425 for (int i = 0; i <= nbinsMultK0pp; i++) {
426 zbinsMultK0pp[i] = zminMultK0pp + i * dzMultK0pp;
427 }
428
429 const Int_t nbinsMultK0PbPb = mNBinsTPCOccPbPb;
430 Double_t* zbinsMultK0PbPb = new Double_t[nbinsMultK0PbPb + 1];
431 Double_t zminMultK0PbPb = mMinTPCOccPbPb;
432 Double_t zmaxMultK0PbPb = mMaxTPCOccPbPb;
433 Double_t dzMultK0PbPb = (zmaxMultK0PbPb - zminMultK0PbPb) / nbinsMultK0PbPb;
434 for (int i = 0; i <= nbinsMultK0PbPb; i++) {
435 zbinsMultK0PbPb[i] = zminMultK0PbPb + i * dzMultK0PbPb;
436 }
437
438 if (mDoK0QC) {
439 // V0s
440 mK0MassVsPtVsOccpp = new TH3F("mK0MassVsPtVsOccpp", "K0 invariant mass vs Pt vs TPC occupancy; Pt [GeV/c]; K0s mass [GeV/c^2]; TPC occ.", nbinsPtK0, xbinsPtK0, nbinsMassK0, ybinsMassK0, nbinsMultK0pp, zbinsMultK0pp);
441
442 mK0MassVsPtVsOccPbPb = new TH3F("mK0MassVsPtVsOccPbPb", "K0 invariant mass vs Pt vs TPC occupancy; Pt [GeV/c]; K0s mass [GeV/c^2]; TPC occ", nbinsPtK0, xbinsPtK0, nbinsMassK0, ybinsMassK0, nbinsMultK0PbPb, zbinsMultK0PbPb);
443 }
444
445 LOG(info) << "Printing configuration cuts";
446 printParams();
447
448 delete[] xbinsPt;
449 delete[] xbinsPtK0;
450 delete[] ybinsMassK0;
451 delete[] zbinsMultK0pp;
452 delete[] zbinsMultK0PbPb;
453
454 return true;
455}
456
457//__________________________________________________________
458
459void MatchITSTPCQC::initDataRequest()
460{
461
462 // initialize data request, if it was not already done
463
464 mSrc &= mAllowedSources;
465
466 if (!mSrc[GID::Source::ITSTPC] || !mSrc[GID::Source::TPC] || !mSrc[GID::Source::ITS]) {
467 LOG(fatal) << "We cannot do ITSTPC QC, some sources are missing, check sources in " << mSrc;
468 }
469
470 mDataRequest = std::make_shared<o2::globaltracking::DataRequest>();
471 mDataRequest->requestTracks(mSrc, mUseMC);
472 if (mDoK0QC) {
473 mDataRequest->requestPrimaryVertices(mUseMC);
474 mDataRequest->requestSecondaryVertices(mUseMC);
475 mDataRequest->requestTPCOccMap();
476 }
477}
478
479//__________________________________________________________
480
481void MatchITSTPCQC::run(o2::framework::ProcessingContext& ctx)
482{
483 // Getting the B field
484 mBz = o2::base::Propagator::Instance()->getNominalBz();
485
486 // Getting the SVertexer config params
487 if (mTimestamp == -1 && mDoK0QC) {
488 // we have not yet initialized the SVertexer params; let's do it
489 ctx.inputs().get<o2::vertexing::SVertexerParams*>("SVParam");
490 const auto& svparam = o2::vertexing::SVertexerParams::Instance();
491 mFitterV0.setUseAbsDCA(svparam.useAbsDCA);
492 mFitterV0.setMaxR(svparam.maxRIni);
493 mFitterV0.setMinParamChange(svparam.minParamChange);
494 mFitterV0.setMinRelChi2Change(svparam.minRelChi2Change);
495 mFitterV0.setMaxDZIni(svparam.maxDZIni);
496 mFitterV0.setMaxDXYIni(svparam.maxDXYIni);
497 mFitterV0.setMaxChi2(svparam.maxChi2);
498 mFitterV0.setMatCorrType(o2::base::Propagator::MatCorrType(svparam.matCorr));
499 mFitterV0.setUsePropagator(svparam.usePropagator);
500 mFitterV0.setRefitWithMatCorr(svparam.refitWithMatCorr);
501 mFitterV0.setMaxStep(svparam.maxStep);
502 mFitterV0.setMaxSnp(svparam.maxSnp);
503 mFitterV0.setMinXSeed(svparam.minXSeed);
504
505 mTimestamp = ctx.services().get<o2::framework::TimingInfo>().creation;
506 auto grplhcif = o2::base::GRPGeomHelper::instance().getGRPLHCIF();
507 if (grplhcif->getBeamZ(0) != 1 || grplhcif->getBeamZ(1) != 1) {
508 LOG(info) << "We are in Heavy Ion: Z for beam 0 = " << grplhcif->getBeamZ(0) << " ; Z for beam 1 = " << grplhcif->getBeamZ(1);
509 mIsHI = true;
510 } else {
511 LOG(info) << "We are not in Heavy Ion: Z for beam 0 = " << grplhcif->getBeamZ(0) << " ; Z for beam 1 = " << grplhcif->getBeamZ(1);
512 }
513 }
514
515 static int evCount = 0;
516 mRecoCont.collectData(ctx, *mDataRequest);
517 mTPCTracks = mRecoCont.getTPCTracks();
518 mITSTracks = mRecoCont.getITSTracks();
519 mITSTPCTracks = mRecoCont.getTPCITSTracks();
520
521 LOG(info) << "****** Number of found ITSTPC tracks = " << mITSTPCTracks.size();
522 LOG(info) << "****** Number of found TPC tracks = " << mTPCTracks.size();
523 LOG(info) << "****** Number of found ITS tracks = " << mITSTracks.size();
524
525 // cache selection for TPC and ITS tracks
526 std::vector<bool> isTPCTrackSelectedEntry(mTPCTracks.size(), false);
527 std::vector<bool> isITSTrackSelectedEntry(mITSTracks.size(), false);
528 TrackCuts cuts;
529 // ITS track
530 cuts.setMinPtITSCut(mPtITSCut);
531 cuts.setEtaITSCut(mEtaITSCut);
532 cuts.setMinNClustersITS(mMinNClustersITS);
533 cuts.setMaxChi2PerClusterITS(mMaxChi2PerClusterITS);
534 for (auto it = mRequiredITSHits.begin(); it != mRequiredITSHits.end(); it++) {
535 cuts.setRequireHitsInITSLayers((*it).first, (*it).second);
536 }
537 // TPC track
538 cuts.setMinPtTPCCut(mPtTPCCut);
539 cuts.setEtaTPCCut(mEtaTPCCut);
540 cuts.setMinNTPCClustersCut(mNTPCClustersCut);
541 cuts.setMaxDCATPCCut(mDCATPCCut);
542 cuts.setMaxDCATPCCutY(mDCATPCCutY);
543 // ITS-TPC track kinematics
544 cuts.setMinPtCut(mPtCut);
545 cuts.setMaxPtCut(mPtMaxCut);
546 cuts.setEtaCut(-mEtaCut, mEtaCut);
547
548 for (size_t itrk = 0; itrk < mTPCTracks.size(); ++itrk) {
549 auto const& trkTpc = mTPCTracks[itrk];
550 o2::dataformats::GlobalTrackID id(itrk, GID::TPC);
551 if (cuts.isSelected(id, mRecoCont)) {
552 // NB: same cuts for numerator and denominator tracks of ITS-TPC matching
553 // To change cuts only for numerator, something like o2::dataformats::GlobalTrackID id(itrk, GID::ITSTPC) is necessary
554 isTPCTrackSelectedEntry[itrk] = true;
555 }
556 }
557
558 for (size_t itrk = 0; itrk < mITSTracks.size(); ++itrk) {
559 auto const& trkIts = mITSTracks[itrk];
560 o2::dataformats::GlobalTrackID id(itrk, GID::ITS);
561 if (cuts.isSelected(id, mRecoCont)) {
562 // NB: same cuts for numerator and denominator tracks of ITS-TPC matching
563 // To change cuts only for numerator, something like o2::dataformats::GlobalTrackID id(itrk, GID::ITSTPC) is necessary
564 isITSTrackSelectedEntry[itrk] = true;
565 }
566 }
567
568 // physical-primary flag of all labelled tracks, per (source, event); the kinematics are read and released one event at a time
569 std::map<std::pair<int, int>, std::unordered_map<int, bool>> physPrimary;
570 auto isPhysPrimary = [&physPrimary](const o2::MCCompLabel& lbl) { return physPrimary.at({lbl.getSourceID(), lbl.getEventID()}).at(lbl.getTrackID()); };
571 if (mUseMC) {
572 for (auto const& labels : {mRecoCont.getTPCITSTracksMCLabels(), mRecoCont.getTPCTracksMCLabels(), mRecoCont.getITSTracksMCLabels()}) {
573 for (auto const& lbl : labels) {
574 if (lbl.isValid()) {
575 physPrimary[{lbl.getSourceID(), lbl.getEventID()}][lbl.getTrackID()] = false;
576 }
577 }
578 }
579 for (auto& [sourceAndEvent, flags] : physPrimary) {
580 const auto [source, event] = sourceAndEvent;
581 const std::vector<o2::MCTrack>& pcontainer = mcReader.getTracks(source, event);
582 for (auto& [trackID, flag] : flags) {
583 flag = static_cast<size_t>(trackID) < pcontainer.size() && MCTrackNavigator::isPhysicalPrimary(pcontainer[trackID], pcontainer);
584 }
585 mcReader.releaseTracksForSourceAndEvent(source, event);
586 }
587 }
588
589 // numerator + eta, chi2...
590 if (mUseMC) {
591 for (int i = 0; i < matchType::SIZE; ++i) {
592 mMapLabels[i].clear();
593 }
594 for (int itrk = 0; itrk < static_cast<int>(mITSTPCTracks.size()); ++itrk) {
595 auto const& trk = mITSTPCTracks[itrk];
596 auto idxTrkTpc = trk.getRefTPC().getIndex();
597 if (trk.getRefITS().getSource() != GID::ITS) {
598 continue;
599 }
600 if (isTPCTrackSelectedEntry[idxTrkTpc]) {
601 auto lbl = mRecoCont.getTrackMCLabel({(unsigned int)(itrk), GID::Source::ITSTPC});
602 if (!lbl.isValid()) {
603 continue;
604 }
605 if (mMapLabels[matchType::TPC].find(lbl) == mMapLabels[matchType::TPC].end()) {
606 if (isPhysPrimary(lbl)) {
607 mMapLabels[matchType::TPC].insert({lbl, {.mIdx = itrk, .mIsPhysicalPrimary = true}});
608 } else {
609 mMapLabels[matchType::TPC].insert({lbl, {.mIdx = itrk, .mIsPhysicalPrimary = false}});
610 }
611 } else {
612 // winner (if more tracks have the same label) has the highest pt
613 if (mITSTPCTracks[mMapLabels[matchType::TPC].at(lbl).mIdx].getPt() < trk.getPt()) {
614 mMapLabels[matchType::TPC].at(lbl).mIdx = itrk;
615 }
616 }
617 }
618 auto idxTrkIts = trk.getRefITS().getIndex();
619 if (isITSTrackSelectedEntry[idxTrkIts]) {
620 auto lbl = mRecoCont.getTrackMCLabel({(unsigned int)(itrk), GID::Source::ITSTPC});
621 if (!lbl.isValid()) {
622 continue;
623 }
624 if (mMapLabels[matchType::ITS].find(lbl) == mMapLabels[matchType::ITS].end()) {
625 if (isPhysPrimary(lbl)) {
626 mMapLabels[matchType::ITS].insert({lbl, {.mIdx = itrk, .mIsPhysicalPrimary = true}});
627 } else {
628 mMapLabels[matchType::ITS].insert({lbl, {.mIdx = itrk, .mIsPhysicalPrimary = false}});
629 }
630 } else {
631 // winner (if more tracks have the same label) has the highest pt
632 if (mITSTPCTracks[mMapLabels[matchType::ITS].at(lbl).mIdx].getPt() < trk.getPt()) {
633 mMapLabels[matchType::ITS].at(lbl).mIdx = itrk;
634 }
635 }
636 }
637 }
638 LOG(debug) << "number of entries in map for nominator (without duplicates) = " << mMapLabels.size();
639 // now we use only the tracks in the map to fill the histograms (--> tracks have passed the
640 // track selection and there are no duplicated tracks wrt the same MC label)
641 for (int i = 0; i < matchType::SIZE; ++i) {
642 for (auto const& el : mMapLabels[i]) {
643 auto const& trk = mITSTPCTracks[el.second.mIdx];
644 o2::track::TrackParCov trkDen;
645 bool isEtaITSOk = true;
646 if (i == matchType::TPC) {
647 trkDen = mTPCTracks[trk.getRefTPC()];
648 } else {
649 trkDen = mITSTracks[trk.getRefITS()];
650 if (std::abs(trkDen.getEta()) > mEtaITSCut) {
651 // ITS track outside |eta | < 0.9, we don't fill pt, nor phi , nor phi vs pt histos
652 isEtaITSOk = false;
653 }
654 }
655 if (isEtaITSOk) {
656 mPtNum[i]->Fill(trkDen.getPt());
657 if (std::abs(trkDen.getEta()) > mEtaNo0Cut) {
658 mPtNum_noEta0[i]->Fill(trkDen.getPt());
659 }
660 mPhiNum[i]->Fill(trkDen.getPhi());
661 mPhiVsPtNum[i]->Fill(trkDen.getPt(), trkDen.getPhi());
662 m1OverPtNum[i]->Fill(trkDen.getSign() * trkDen.getPtInv());
663 mEtaPhiPtNum[i]->Fill(trkDen.getEta(), trkDen.getPhi(), trkDen.getPt());
664 // we fill also the denominator
665 mPtDen[i]->Fill(trkDen.getPt());
666 if (std::abs(trkDen.getEta()) > mEtaNo0Cut) {
667 mPtDen_noEta0[i]->Fill(trkDen.getPt());
668 }
669 mPhiDen[i]->Fill(trkDen.getPhi());
670 mPhiVsPtDen[i]->Fill(trkDen.getPt(), trkDen.getPhi());
671 m1OverPtDen[i]->Fill(trkDen.getSign() * trkDen.getPtInv());
672 mEtaPhiPtDen[i]->Fill(trkDen.getEta(), trkDen.getPhi(), trkDen.getPt());
673 if (mUseTrkPID) { // Vs Tracking PID hypothesis
674 mPtNumVsTrkPID[i][trkDen.getPID()]->Fill(trkDen.getPt());
675 mPhiNumVsTrkPID[i][trkDen.getPID()]->Fill(trkDen.getPhi());
676 // we fill also the denominator
677 mPtDenVsTrkPID[i][trkDen.getPID()]->Fill(trkDen.getPt());
678 mPhiDenVsTrkPID[i][trkDen.getPID()]->Fill(trkDen.getPhi());
679 }
680 }
681 mEtaNum[i]->Fill(trkDen.getEta());
682 mEtaVsPtNum[i]->Fill(trkDen.getPt(), trkDen.getEta());
683 // we fill also the denominator
684 mEtaDen[i]->Fill(trkDen.getEta());
685 mEtaVsPtDen[i]->Fill(trkDen.getPt(), trkDen.getEta());
686 if (i == matchType::TPC) {
687 auto tpcTrk = mTPCTracks[trk.getRefTPC()];
688 mClsVsPtNum[i]->Fill(tpcTrk.getPt(), tpcTrk.getNClusters());
689 mChi2VsPtNum[i]->Fill(tpcTrk.getPt(), tpcTrk.getChi2());
690 mClsVsPtDen[i]->Fill(tpcTrk.getPt(), tpcTrk.getNClusters());
691 mChi2VsPtDen[i]->Fill(tpcTrk.getPt(), tpcTrk.getChi2());
692 math_utils::Point3D<float> v{};
693 std::array<float, 2> dca{};
694 if (tpcTrk.propagateParamToDCA(v, mBz, &dca)) {
695 mDCArVsPtNum->Fill(tpcTrk.getPt(), dca[0]);
696 mDCArVsPtDen->Fill(tpcTrk.getPt(), dca[0]);
697 }
698 } else {
699 const auto& itsTrk = mITSTracks[trk.getRefITS()];
700 mClsVsPtNum[i]->Fill(itsTrk.getPt(), itsTrk.getNClusters());
701 mChi2VsPtNum[i]->Fill(itsTrk.getPt(), itsTrk.getChi2());
702 mClsVsPtDen[i]->Fill(itsTrk.getPt(), itsTrk.getNClusters());
703 mChi2VsPtDen[i]->Fill(itsTrk.getPt(), itsTrk.getChi2());
704 }
705 if (mUseTrkPID) { // Vs Tracking PID hypothesis
706 mEtaNumVsTrkPID[i][trkDen.getPID()]->Fill(trkDen.getEta());
707 // we fill also the denominator
708 mEtaDenVsTrkPID[i][trkDen.getPID()]->Fill(trkDen.getEta());
709 }
710 if (el.second.mIsPhysicalPrimary) {
711 if (isEtaITSOk) {
712 mPtPhysPrimNum[i]->Fill(trkDen.getPt());
713 mPhiPhysPrimNum[i]->Fill(trkDen.getPhi());
714 m1OverPtPhysPrimNum[i]->Fill(trkDen.getSign() * trkDen.getPtInv());
715 // we fill also the denominator
716 mPtPhysPrimDen[i]->Fill(trkDen.getPt());
717 mPhiPhysPrimDen[i]->Fill(trkDen.getPhi());
718 m1OverPtPhysPrimDen[i]->Fill(trkDen.getSign() * trkDen.getPtInv());
719 }
720 mEtaPhysPrimNum[i]->Fill(trkDen.getEta());
721 // we fill also the denominator
722 mEtaPhysPrimDen[i]->Fill(trkDen.getEta());
723 }
724 ++mNITSTPCSelectedTracks[i];
725 }
726 }
727 }
728 int iITSTPC = 0;
729 for (auto const& trk : mITSTPCTracks) {
730 if (trk.getRefTPC().getIndex() >= mTPCTracks.size()) {
731 LOG(fatal) << "******************** ATTENTION! for TPC track associated to matched track: idx = " << trk.getRefTPC().getIndex() << ", size of container = " << mTPCTracks.size() << " in TF " << evCount;
732 }
733 std::array<std::string, 2> title{"TPC", "ITS"};
734 for (int i = 0; i < matchType::SIZE; ++i) {
735 o2::track::TrackParCov trkRef;
736 unsigned int idxTrkRef{0};
737 bool fillHisto = false;
738 bool isEtaITSOk = true;
739 if (i == matchType::TPC) {
740 trkRef = mTPCTracks[trk.getRefTPC()];
741 idxTrkRef = trk.getRefTPC().getIndex();
742 if (isTPCTrackSelectedEntry[idxTrkRef]) {
743 fillHisto = true;
744 ++mNITSTPCSelectedTracks[i];
745 }
746 } else {
747 idxTrkRef = trk.getRefITS().getIndex();
748 if (trk.getRefITS().getSource() == GID::ITSAB) {
749 // do not use afterburner tracks
750 LOG(debug) << "Track (ITS) with id " << idxTrkRef << " for ITSTPC track " << iITSTPC << " is from afterburner";
751 continue;
752 }
753 if (idxTrkRef >= mITSTracks.size()) {
754 LOG(fatal) << "******************** ATTENTION! for ITS track associated to matched track (NOT from AB): idx = " << trk.getRefITS().getIndex() << ", size of container = " << mITSTracks.size() << " in TF " << evCount;
755 }
756 trkRef = mITSTracks[trk.getRefITS()];
757 LOG(debug) << "Checking track (ITS) with id " << idxTrkRef << " for ITSTPC track " << iITSTPC << " and pt = " << trkRef.getPt();
758 if (isITSTrackSelectedEntry[idxTrkRef]) {
759 LOG(debug) << "Track was selected (ITS), with id " << idxTrkRef << " for ITSTPC track " << iITSTPC << " , we keep it in the numerator, pt = " << trkRef.getPt();
760 fillHisto = true;
761 ++mNITSTPCSelectedTracks[i];
762 } else {
763 LOG(debug) << "Track was not selected (ITS), with id " << idxTrkRef << " for ITSTPC track " << iITSTPC << " , we don't keep it in the numerator, pt = " << trkRef.getPt();
764 }
765 if (std::abs(trkRef.getEta()) > mEtaITSCut) {
766 // ITS track outside |eta | < 0.9, we don't fill pt, nor phi , nor phi vs pt histos
767 isEtaITSOk = false;
768 LOG(debug) << "Track (ITS), with id " << idxTrkRef << " for ITSTPC track " << iITSTPC << " will be discarded when filling pt of phi related histograms, since eta = " << trkRef.getEta() << " , we don't keep it in the numerator, pt = " << trkRef.getPt();
769 }
770 }
771 if (fillHisto) {
772 if (!mUseMC) {
773 LOG(debug) << "Filling num (" << title[i] << ") with track with id " << idxTrkRef << " for ITSTPC track " << iITSTPC << " with pt = " << trkRef.getPt();
774 if (isEtaITSOk) {
775 mPtNum[i]->Fill(trkRef.getPt());
776 if (std::abs(trkRef.getEta()) > mEtaNo0Cut) {
777 mPtNum_noEta0[i]->Fill(trkRef.getPt());
778 }
779 mPhiNum[i]->Fill(trkRef.getPhi());
780 if (mUseTrkPID) { // Vs Tracking PID hypothesis
781 mPtNumVsTrkPID[i][trkRef.getPID()]->Fill(trkRef.getPt());
782 mPhiNumVsTrkPID[i][trkRef.getPID()]->Fill(trkRef.getPhi());
783 }
784 mPhiVsPtNum[i]->Fill(trkRef.getPt(), trkRef.getPhi());
785 m1OverPtNum[i]->Fill(trkRef.getSign() * trkRef.getPtInv());
786 mEtaPhiPtNum[i]->Fill(trkRef.getEta(), trkRef.getPhi(), trkRef.getPt());
787 }
788 mEtaNum[i]->Fill(trkRef.getEta());
789 if (mUseTrkPID) { // Vs Tracking PID hypothesis
790 mEtaNumVsTrkPID[i][trkRef.getPID()]->Fill(trkRef.getEta());
791 }
792 mEtaVsPtNum[i]->Fill(trkRef.getPt(), trkRef.getEta());
793 if (i == matchType::TPC) {
794 const auto& tpcTrk = mTPCTracks[trk.getRefTPC()];
795 mClsVsPtNum[i]->Fill(tpcTrk.getPt(), tpcTrk.getNClusters());
796 mChi2VsPtNum[i]->Fill(tpcTrk.getPt(), tpcTrk.getChi2());
797 } else {
798 const auto& itsTrk = mITSTracks[trk.getRefITS()];
799 mClsVsPtNum[i]->Fill(itsTrk.getPt(), itsTrk.getNClusters());
800 mChi2VsPtNum[i]->Fill(itsTrk.getPt(), itsTrk.getChi2());
801 }
802 }
803 if (i == matchType::TPC) {
804 mResidualPt->Fill(trk.getPt(), trk.getPt() - trkRef.getPt());
805 mResidualPhi->Fill(trk.getPhi(), trk.getPhi() - trkRef.getPhi());
806 mResidualEta->Fill(trk.getEta(), trk.getEta() - trkRef.getEta());
807 mChi2Matching->Fill(trk.getChi2Match());
808 mChi2Refit->Fill(trk.getChi2Refit());
809 mTimeResVsPt->Fill(trkRef.getPt(), trk.getTimeMUS().getTimeStampError());
810 math_utils::Point3D<float> v{};
811 std::array<float, 2> dca{-999, -999};
812 if (trkRef.propagateParamToDCA(v, mBz, &dca)) {
813 mDCAr->Fill(dca[0]);
814 if (!mUseMC) {
815 mDCArVsPtNum->Fill(trkRef.getPt(), dca[0]);
816 }
817 }
818 LOG(debug) << "*** chi2Matching = " << trk.getChi2Match() << ", chi2refit = " << trk.getChi2Refit() << ", timeResolution = " << trk.getTimeMUS().getTimeStampError();
819 }
820 } else {
821 LOG(debug) << "Not filling num (" << title[i] << ") for ITSTPC track " << iITSTPC << " for track with pt " << trkRef.getPt();
822 }
823 }
824 ++iITSTPC;
825 }
826
827 // now filling the denominator for the efficiency calculation
828 if (mUseMC) {
829 for (int i = 0; i < matchType::SIZE; ++i) {
830 mMapRefLabels[i].clear();
831 }
832 // filling the map where we store for each MC label, the track id of the reconstructed
833 // track with the highest number of TPC clusters
834 for (int itrk = 0; itrk < static_cast<int>(mTPCTracks.size()); ++itrk) {
835 auto const& trk = mTPCTracks[itrk];
836 if (isTPCTrackSelectedEntry[itrk]) {
837 auto lbl = mRecoCont.getTrackMCLabel({(unsigned int)(itrk), GID::Source::TPC});
838 if (!lbl.isValid()) {
839 continue;
840 }
841 if (mMapLabels[matchType::TPC].find(lbl) != mMapLabels[matchType::TPC].end()) {
842 // the track was already added to the denominator
843 continue;
844 }
845 if (mMapRefLabels[matchType::TPC].find(lbl) == mMapRefLabels[matchType::TPC].end()) {
846 if (isPhysPrimary(lbl)) {
847 mMapRefLabels[matchType::TPC].insert({lbl, {itrk, true}});
848 } else {
849 mMapRefLabels[matchType::TPC].insert({lbl, {itrk, false}});
850 }
851 } else {
852 // winner (if more tracks have the same label) has the highest number of TPC clusters
853 if (mTPCTracks[mMapRefLabels[matchType::TPC].at(lbl).mIdx].getNClusters() < trk.getNClusters()) {
854 mMapRefLabels[matchType::TPC].at(lbl).mIdx = itrk;
855 }
856 }
857 }
858 }
859 // same for ITS
860 // filling the map where we store for each MC label, the track id of the reconstructed
861 // track with the highest number of ITS clusters
862 for (int itrk = 0; itrk < static_cast<int>(mITSTracks.size()); ++itrk) {
863 auto const& trk = mITSTracks[itrk];
864 if (isITSTrackSelectedEntry[itrk]) {
865 auto lbl = mRecoCont.getTrackMCLabel({(unsigned int)(itrk), GID::Source::ITS});
866 if (!lbl.isValid()) {
867 continue;
868 }
869 if (mMapLabels[matchType::ITS].find(lbl) != mMapLabels[matchType::ITS].end()) {
870 // the track was already added to the denominator
871 continue;
872 }
873 if (mMapRefLabels[matchType::ITS].find(lbl) == mMapRefLabels[matchType::ITS].end()) {
874 if (isPhysPrimary(lbl)) {
875 mMapRefLabels[matchType::ITS].insert({lbl, {itrk, true}});
876 } else {
877 mMapRefLabels[matchType::ITS].insert({lbl, {itrk, false}});
878 }
879 } else {
880 // winner (if more tracks have the same label) has the highest number of ITS clusters
881 if (mITSTracks[mMapRefLabels[matchType::ITS].at(lbl).mIdx].getNClusters() < trk.getNClusters()) {
882 mMapRefLabels[matchType::ITS].at(lbl).mIdx = itrk;
883 }
884 }
885 }
886 }
887 LOG(debug) << "number of entries in map for denominator of TPC tracks (without duplicates) = " << mMapRefLabels[matchType::TPC].size() + mMapLabels[matchType::TPC].size();
888 LOG(debug) << "number of entries in map for denominator of ITS tracks (without duplicates) = " << mMapRefLabels[matchType::ITS].size() + mMapLabels[matchType::ITS].size();
889 // now we use only the tracks in the map to fill the histograms (--> tracks have passed the
890 // track selection and there are no duplicated tracks wrt the same MC label)
891 for (auto const& el : mMapRefLabels[matchType::TPC]) {
892 auto const& trk = mTPCTracks[el.second.mIdx];
893 mPtDen[matchType::TPC]->Fill(trk.getPt());
894 if (std::abs(trk.getEta()) > mEtaNo0Cut) {
895 mPtDen_noEta0[matchType::TPC]->Fill(trk.getPt());
896 }
897 mPhiDen[matchType::TPC]->Fill(trk.getPhi());
898 mPhiVsPtDen[matchType::TPC]->Fill(trk.getPt(), trk.getPhi());
899 mEtaPhiPtDen[matchType::TPC]->Fill(trk.getEta(), trk.getPhi(), trk.getPt());
900 mEtaDen[matchType::TPC]->Fill(trk.getEta());
901 mEtaVsPtDen[matchType::TPC]->Fill(trk.getPt(), trk.getEta());
902 m1OverPtDen[matchType::TPC]->Fill(trk.getSign() * trk.getPtInv());
903 mClsVsPtDen[matchType::TPC]->Fill(trk.getPt(), trk.getNClusters());
904 mChi2VsPtDen[matchType::TPC]->Fill(trk.getPt(), trk.getChi2());
905 math_utils::Point3D<float> v{};
906 std::array<float, 2> dca{};
907 if (auto trc = trk; trc.propagateParamToDCA(v, mBz, &dca)) {
908 mDCArVsPtDen->Fill(trc.getPt(), dca[0]);
909 }
910 if (el.second.mIsPhysicalPrimary) {
911 mPtPhysPrimDen[matchType::TPC]->Fill(trk.getPt());
912 mPhiPhysPrimDen[matchType::TPC]->Fill(trk.getPhi());
913 mEtaPhysPrimDen[matchType::TPC]->Fill(trk.getEta());
914 m1OverPtPhysPrimDen[matchType::TPC]->Fill(trk.getSign() * trk.getPtInv());
915 }
916 ++mNTPCSelectedTracks;
917 }
918 for (auto const& el : mMapRefLabels[matchType::ITS]) {
919 auto const& trk = mITSTracks[el.second.mIdx];
920 if (std::abs(trk.getEta()) < mEtaITSCut) {
921 mPtDen[matchType::ITS]->Fill(trk.getPt());
922 if (std::abs(trk.getEta()) > mEtaNo0Cut) {
923 mPtDen_noEta0[matchType::ITS]->Fill(trk.getPt());
924 }
925 mPhiDen[matchType::ITS]->Fill(trk.getPhi());
926 mPhiVsPtDen[matchType::ITS]->Fill(trk.getPt(), trk.getPhi());
927 mEtaPhiPtDen[matchType::ITS]->Fill(trk.getEta(), trk.getPhi(), trk.getPt());
928 m1OverPtDen[matchType::ITS]->Fill(trk.getSign() * trk.getPtInv());
929 }
930 mEtaDen[matchType::ITS]->Fill(trk.getEta());
931 mEtaVsPtDen[matchType::ITS]->Fill(trk.getPt(), trk.getEta());
932 mClsVsPtDen[matchType::ITS]->Fill(trk.getPt(), trk.getNClusters());
933 mChi2VsPtDen[matchType::ITS]->Fill(trk.getPt(), trk.getChi2());
934 if (el.second.mIsPhysicalPrimary) {
935 if (std::abs(trk.getEta()) < mEtaITSCut) {
936 mPtPhysPrimDen[matchType::ITS]->Fill(trk.getPt());
937 mPhiPhysPrimDen[matchType::ITS]->Fill(trk.getPhi());
938 m1OverPtPhysPrimDen[matchType::ITS]->Fill(trk.getSign() * trk.getPtInv());
939 }
940 mEtaPhysPrimDen[matchType::ITS]->Fill(trk.getEta());
941 }
942 ++mNITSSelectedTracks;
943 }
944 } else {
945 // if we are in data, we loop over all tracks (no check on the label)
946 for (size_t itrk = 0; itrk < mTPCTracks.size(); ++itrk) {
947 auto const& trk = mTPCTracks[itrk];
948 if (isTPCTrackSelectedEntry[itrk]) {
949 LOG(debug) << "Filling den (TPC) with track with pt = " << trk.getPt();
950 mPtDen[matchType::TPC]->Fill(trk.getPt());
951 if (std::abs(trk.getEta()) > mEtaNo0Cut) {
952 mPtDen_noEta0[matchType::TPC]->Fill(trk.getPt());
953 } else {
954 LOG(debug) << "Track (ITS) " << itrk << " with pt = " << trk.getPt() << " and eta = " << trk.getEta() << " not used for den pt, phi, phi vs pt, 1.pt histos";
955 }
956 mPhiDen[matchType::TPC]->Fill(trk.getPhi());
957 mPhiVsPtDen[matchType::TPC]->Fill(trk.getPt(), trk.getPhi());
958 mEtaPhiPtDen[matchType::TPC]->Fill(trk.getEta(), trk.getPhi(), trk.getPt());
959 mEtaDen[matchType::TPC]->Fill(trk.getEta());
960 mEtaVsPtDen[matchType::TPC]->Fill(trk.getPt(), trk.getEta());
961 m1OverPtDen[matchType::TPC]->Fill(trk.getSign() * trk.getPtInv());
962 mClsVsPtDen[matchType::TPC]->Fill(trk.getPt(), trk.getNClusters());
963 mChi2VsPtDen[matchType::TPC]->Fill(trk.getPt(), trk.getChi2());
964 math_utils::Point3D<float> v{};
965 std::array<float, 2> dca{};
966 if (auto trc = trk; trc.propagateParamToDCA(v, mBz, &dca)) {
967 mDCArVsPtDen->Fill(trc.getPt(), dca[0]);
968 }
969 ++mNTPCSelectedTracks;
970 }
971 }
972 for (size_t itrk = 0; itrk < mITSTracks.size(); ++itrk) {
973 auto const& trk = mITSTracks[itrk];
974 LOG(debug) << "Checking den for track (ITS) " << itrk << " with pt " << trk.getPt() << " and eta = " << trk.getEta();
975 if (isITSTrackSelectedEntry[itrk]) {
976 if (std::abs(trk.getEta()) < mEtaITSCut) {
977 LOG(debug) << "Filling den for track (ITS) " << itrk << " with pt = " << trk.getPt() << " and eta = " << trk.getEta();
978 mPtDen[matchType::ITS]->Fill(trk.getPt());
979 if (std::abs(trk.getEta()) > mEtaNo0Cut) {
980 mPtDen_noEta0[matchType::ITS]->Fill(trk.getPt());
981 }
982 mPhiDen[matchType::ITS]->Fill(trk.getPhi());
983 mPhiVsPtDen[matchType::ITS]->Fill(trk.getPt(), trk.getPhi());
984 mEtaPhiPtDen[matchType::ITS]->Fill(trk.getEta(), trk.getPhi(), trk.getPt());
985 m1OverPtDen[matchType::ITS]->Fill(trk.getSign() * trk.getPtInv());
986 } else {
987 LOG(debug) << "Track (ITS) " << itrk << " with pt = " << trk.getPt() << " and eta = " << trk.getEta() << " not used for num pt, phi, phi vs pt, 1.pt histos";
988 }
989 mEtaDen[matchType::ITS]->Fill(trk.getEta());
990 mEtaVsPtDen[matchType::ITS]->Fill(trk.getPt(), trk.getEta());
991 mClsVsPtDen[matchType::ITS]->Fill(trk.getPt(), trk.getNClusters());
992 mChi2VsPtDen[matchType::ITS]->Fill(trk.getPt(), trk.getChi2());
993 ++mNITSSelectedTracks;
994 } else {
995 LOG(debug) << "Not filling for this track (ITS) " << itrk << " with pt = " << trk.getPt();
996 }
997 }
998 }
999
1000 if (mDoK0QC && mRecoCont.getPrimaryVertices().size() > 0) {
1001 // now doing K0S
1002 mFitterV0.setBz(mBz);
1003 const auto pvertices = mRecoCont.getPrimaryVertices();
1004 LOG(info) << "****** Number of PVs = " << pvertices.size();
1005
1006 // getting occupancy estimator
1007 mNHBPerTF = o2::base::GRPGeomHelper::instance().getGRPECS()->getNHBFPerTF();
1008 if (!mParam) {
1009 // for occupancy estimator
1010 mParam = o2::gpu::GPUO2InterfaceUtils::getFullParamShared(0.f, mNHBPerTF);
1011 }
1012 size_t occupancyMapSizeBytes = o2::gpu::GPUO2InterfaceRefit::fillOccupancyMapGetSize(mNHBPerTF, mParam.get());
1013 LOG(debug) << "occupancyMapSizeBytes = " << occupancyMapSizeBytes;
1014 mTPCRefitterOccMap = mRecoCont.occupancyMapTPC;
1015 o2::gpu::GPUO2InterfaceUtils::paramUseExternalOccupancyMap(mParam.get(), mNHBPerTF, mTPCRefitterOccMap.data(), occupancyMapSizeBytes);
1016
1017 std::vector<float> mTBinClOcc;
1018 mTBinClOcc.clear();
1019 int mNTPCOccBinLength = mParam->rec.tpc.occupancyMapTimeBins;
1020 LOG(debug) << "mNTPCOccBinLength = " << mNTPCOccBinLength;
1021 mNTPCOccBinLengthInv = 1. / mNTPCOccBinLength;
1022 if (mNTPCOccBinLength > 1 && mTPCRefitterOccMap.size()) {
1023 int nTPCBinsInTF = mNHBPerTF * o2::constants::lhc::LHCMaxBunches / 8; // number of TPC time bins in 1 TF, considering that 1 TPC time bin is 8 bunches
1024 int ninteg = 0;
1025 int nTPCOccBinsInTF = nTPCBinsInTF * mNTPCOccBinLengthInv; // how many occupancy bins in 1 TF; mNTPCOccBinLengthInv is the inverse of the length of an occupancy bin
1026 int sumBins = std::max(1, int(o2::constants::lhc::LHCMaxBunches / 8 * mNTPCOccBinLengthInv)); // we will integrate occupancy at max for this number of bins: the max between 1 and the number of occupancy bins in 1 orbit
1027 LOG(debug) << "number of TPC TB in 1 TF = nTPCBinsInTF = " << nTPCBinsInTF << " ; number of occupancy bins in 1 TF = nTPCOccBinsInTF = " << nTPCOccBinsInTF;
1028 LOG(debug) << "bins to integrate = sumBins = " << sumBins;
1029 mTBinClOcc.resize(nTPCOccBinsInTF);
1030 std::vector<float> mltHistTB(nTPCOccBinsInTF);
1031 float sm = 0., tb = 0.5 * mNTPCOccBinLength;
1032 bool foundNotZero = false;
1033 for (int i = 0; i < nTPCOccBinsInTF; i++) { // for every occupancy bin in the TF
1034 mltHistTB[i] = mParam->GetUnscaledMult(tb);
1035 if (mParam->GetUnscaledMult(tb) != 0) {
1036 LOG(debug) << "i = " << i << " tb = " << tb << " mltHistTB[" << i << "] = " << mltHistTB[i];
1037 foundNotZero = true;
1038 }
1039 tb += mNTPCOccBinLength;
1040 }
1041 if (!foundNotZero) {
1042 LOG(debug) << "No mult bin was found different from 0!";
1043 }
1044 foundNotZero = false;
1045 // now we fill the occupancy map; we integrate the sumBins after the current one, but when we are at the last 27 bins of the TF, where we integrate what we have left till the end of the TF; for practical reasons, we start from the end, adding all the time, and then also removing the last bin, when we have enough, so that we always add together sumBins bins (except, as said, for the last part of the TF)
1046 for (int i = nTPCOccBinsInTF; i--;) {
1047 if (mltHistTB[i] != 0) {
1048 foundNotZero = true;
1049 }
1050 LOG(debug) << "i = " << i << " sm before = " << sm;
1051 sm += mltHistTB[i];
1052 LOG(debug) << "i = " << i << " sm after = " << sm;
1053 if (i + sumBins < nTPCOccBinsInTF) {
1054 LOG(debug) << "i = " << i << " sumBins = " << sumBins << " nTPCOccBinsInTF = " << nTPCOccBinsInTF << " we have to decrease sm by = " << mltHistTB[i + sumBins];
1055 sm -= mltHistTB[i + sumBins];
1056 LOG(debug) << "i = " << i << " sm after 2 = " << sm;
1057 }
1058 mTBinClOcc[i] = sm;
1059 LOG(debug) << "i = " << i << " mTBinClOcc[" << i << "] = " << mTBinClOcc[i];
1060 }
1061 if (!foundNotZero) {
1062 LOG(debug) << "No mult bin was found different from 0! sm = " << sm;
1063 }
1064 } else {
1065 mTBinClOcc.resize(1);
1066 }
1067 auto v0IDs = mRecoCont.getV0sIdx();
1068 auto nv0 = v0IDs.size();
1069 if (nv0 > mRecoCont.getV0s().size()) {
1070 mRefit = true;
1071 }
1072 LOG(debug) << "Found " << mRecoCont.getV0s().size() << " V0s in reco container";
1073 LOG(debug) << "Found " << nv0 << " V0s ids";
1074 // associating sec vtxs to prim vtx
1075 std::map<int, std::vector<int>> pv2sv;
1076 static int tfID = 0;
1077 for (int iv = 0; iv < nv0; iv++) {
1078 const auto v0id = v0IDs[iv];
1079 o2::dataformats::GlobalTrackID id0 = v0id.getProngID(0), id1 = v0id.getProngID(1);
1080 if (id0.getSourceDetectorsMask()[DetID::ITS] && id1.getSourceDetectorsMask()[DetID::ITS] &&
1081 mRecoCont.isTrackSourceLoaded(id0.getSource()) && mRecoCont.isTrackSourceLoaded(id1.getSource())) { // taking only tracks for which there is the ITS in the track sources
1082 pv2sv[v0id.getVertexID()].push_back(iv);
1083 } else {
1084 LOG(debug) << "Skipping this V0: the track sources do not contain ITS";
1085 }
1086 }
1087 int nV0sOk = 0;
1088 // processing every sec vtx for each prim vtx
1089 int myCount = 0;
1090 for (auto it : pv2sv) {
1091 int pvID = it.first;
1092 auto& vv = it.second;
1093 if (pvID < 0 || vv.size() == 0) {
1094 continue;
1095 }
1096 const auto& pv = mRecoCont.getPrimaryVertex(pvID);
1097 float pvTime = pv.getTimeStamp().getTimeStamp(); // in \mus
1098 for (int iv0 : vv) {
1099 nV0sOk += processV0(iv0, mRecoCont, mTBinClOcc, pvTime) ? 1 : 0;
1100 }
1101 ++myCount;
1102 }
1103
1104 LOG(debug) << "Processed " << nV0sOk << " V0s";
1105 }
1106 evCount++;
1107}
1108
1109//__________________________________________________________
1110bool MatchITSTPCQC::processV0(int iv, o2::globaltracking::RecoContainer& recoData, std::vector<float>& mTBinClOcc, float pvTime)
1111{
1112 o2::dataformats::V0 v0;
1113 auto v0s = recoData.getV0s();
1114 auto v0IDs = recoData.getV0sIdx();
1115 static int tfID = 0;
1116
1117 const auto& v0id = v0IDs[iv];
1118 ++mNK0;
1119 if (mNK0 % int(1 / mK0Scaling) == 0) {
1120 LOG(debug) << "Checking " << mNK0 << "th V0: refitting it, since we keep " << mK0Scaling * 100 << "% of all V0s";
1121 } else {
1122 LOG(debug) << "Checking " << mNK0 << "th K0: NOT refitting it, but skipping it, since we keep " << mK0Scaling * 100 << "% of all V0s";
1123 return false;
1124 }
1125 if (mRefit && !refitV0(v0id, v0, recoData)) {
1126 return false;
1127 }
1128 const auto& v0sel = mRefit ? v0 : v0s[iv];
1129 if (mMaxEtaK0 < std::abs(v0sel.getEta())) {
1130 return false;
1131 }
1132
1133 if (mCutK0Mass > 0 && std::abs(std::sqrt(v0sel.calcMass2AsK0()) - 0.497) > mCutK0Mass || v0sel.getDCA() > mK0MaxDCA || v0sel.getCosPA() < mK0MinCosPA) {
1134 if (v0sel.getDCA() > mK0MaxDCA && v0sel.getCosPA() < mK0MinCosPA) {
1135 LOG(debug) << "v0sel.getDCA() = " << v0sel.getDCA() << " max is " << mK0MaxDCA << " returning ... ";
1136 }
1137 if (v0sel.getCosPA() > mK0MinCosPA) {
1138 LOG(debug) << "v0sel.getCosPA() = " << v0sel.getCosPA() << " min is " << mK0MinCosPA << " returning ... ";
1139 }
1140 return false;
1141 }
1142 // get the corresponding PV
1143 int tb = pvTime / (8 * o2::constants::lhc::LHCBunchSpacingMUS) * mNTPCOccBinLengthInv; // V0 time in TPC time bins
1144 LOG(debug) << "pvTime = " << pvTime << " tb = " << tb;
1145 float mltTPC = tb < 0 ? mTBinClOcc[0] : (tb >= mTBinClOcc.size() ? mTBinClOcc.back() : mTBinClOcc[tb]);
1146 ++mNK0;
1147 LOG(debug) << "Filling K0 histogram with pt = " << v0sel.getPt() << " mass = " << std::sqrt(v0sel.calcMass2AsK0()) << " mult TPC = " << mltTPC;
1148 if (!mIsHI) {
1149 mK0MassVsPtVsOccpp->Fill(v0sel.getPt(), std::sqrt(v0sel.calcMass2AsK0()), mltTPC);
1150 } else {
1151 mK0MassVsPtVsOccPbPb->Fill(v0sel.getPt(), std::sqrt(v0sel.calcMass2AsK0()), mltTPC);
1152 }
1153 return true;
1154}
1155
1156//__________________________________________________________
1157bool MatchITSTPCQC::refitV0(const o2::dataformats::V0Index& id, o2::dataformats::V0& v0, o2::globaltracking::RecoContainer& recoData)
1158{
1159 LOG(debug) << "Refitting V0";
1160 if (!recoData.isTrackSourceLoaded(id.getProngID(0).getSource()) || !recoData.isTrackSourceLoaded(id.getProngID(1).getSource())) {
1161 return false;
1162 }
1163 auto seedP = recoData.getTrackParam(id.getProngID(0));
1164 auto seedN = recoData.getTrackParam(id.getProngID(1));
1165 bool isTPConly = (id.getProngID(0).getSource() == o2::dataformats::GlobalTrackID::TPC) || (id.getProngID(1).getSource() == o2::dataformats::GlobalTrackID::TPC);
1166 const auto& svparam = o2::vertexing::SVertexerParams::Instance();
1167 if (svparam.mTPCTrackPhotonTune && isTPConly) {
1168 mFitterV0.setMaxDZIni(svparam.mTPCTrackMaxDZIni);
1169 mFitterV0.setMaxDXYIni(svparam.mTPCTrackMaxDXYIni);
1170 mFitterV0.setMaxChi2(svparam.mTPCTrackMaxChi2);
1171 mFitterV0.setCollinear(true);
1172 }
1173 int nCand = mFitterV0.process(seedP, seedN);
1174 if (svparam.mTPCTrackPhotonTune && isTPConly) { // restore
1175 // Reset immediately to the defaults
1176 mFitterV0.setMaxDZIni(svparam.maxDZIni);
1177 mFitterV0.setMaxDXYIni(svparam.maxDXYIni);
1178 mFitterV0.setMaxChi2(svparam.maxChi2);
1179 mFitterV0.setCollinear(false);
1180 }
1181 if (nCand == 0) { // discard this pair
1182 return false;
1183 }
1184 const int cand = 0;
1185 if (!mFitterV0.isPropagateTracksToVertexDone(cand) && !mFitterV0.propagateTracksToVertex(cand)) {
1186 return false;
1187 }
1188 const auto& trPProp = mFitterV0.getTrack(0, cand);
1189 const auto& trNProp = mFitterV0.getTrack(1, cand);
1190 std::array<float, 3> pP{}, pN{};
1191 trPProp.getPxPyPzGlo(pP);
1192 trNProp.getPxPyPzGlo(pN);
1193 std::array<float, 3> pV0 = {pP[0] + pN[0], pP[1] + pN[1], pP[2] + pN[2]};
1194 auto p2V0 = pV0[0] * pV0[0] + pV0[1] * pV0[1] + pV0[2] * pV0[2];
1195 const auto& pv = recoData.getPrimaryVertex(id.getVertexID());
1196 const auto v0XYZ = mFitterV0.getPCACandidatePos(cand);
1197 float dx = v0XYZ[0] - pv.getX(), dy = v0XYZ[1] - pv.getY(), dz = v0XYZ[2] - pv.getZ(), prodXYZv0 = dx * pV0[0] + dy * pV0[1] + dz * pV0[2];
1198 float cosPA = prodXYZv0 / std::sqrt((dx * dx + dy * dy + dz * dz) * p2V0);
1199 new (&v0) o2::dataformats::V0(v0XYZ, pV0, mFitterV0.calcPCACovMatrixFlat(cand), trPProp, trNProp);
1200 v0.setDCA(mFitterV0.getChi2AtPCACandidate(cand));
1201 v0.setCosPA(cosPA);
1202 return true;
1203}
1204
1205//__________________________________________________________
1206void MatchITSTPCQC::finalize()
1207{
1208
1209 std::array<std::string, 2> title{"TPC", "ITS"};
1210
1211 // first we use denominators and nominators to set the TEfficiency; later they are scaled
1212
1213 // filling the efficiency
1214 for (int ti = 0; ti < matchType::SIZE; ++ti) {
1215 setEfficiency(mFractionITSTPCmatch[ti], mPtNum[ti], mPtDen[ti]);
1216 setEfficiency(mFractionITSTPCmatch_noEta0[ti], mPtNum_noEta0[ti], mPtDen_noEta0[ti]);
1217 setEfficiency(mFractionITSTPCmatchPhi[ti], mPhiNum[ti], mPhiDen[ti]);
1218 setEfficiency(mFractionITSTPCmatchEta[ti], mEtaNum[ti], mEtaDen[ti]);
1219 setEfficiency<2>(mFractionITSTPCmatchPhiVsPt[ti], mPhiVsPtNum[ti], mPhiVsPtDen[ti]);
1220 setEfficiency<2>(mFractionITSTPCmatchEtaVsPt[ti], mEtaVsPtNum[ti], mEtaVsPtDen[ti]);
1221 setEfficiency(mFractionITSTPCmatch1OverPt[ti], m1OverPtNum[ti], m1OverPtDen[ti]);
1222 setEfficiency<2>(mFractionITSTPCmatchClsVsPt[ti], mClsVsPtNum[ti], mClsVsPtDen[ti]);
1223 setEfficiency<2>(mFractionITSTPCmatchChi2VsPt[ti], mChi2VsPtNum[ti], mChi2VsPtDen[ti]);
1224 if (mUseTrkPID) { // Vs Tracking PID hypothesis
1225 for (int j = 0; j < o2::track::PID::NIDs; ++j) {
1226 setEfficiency(mFractionITSTPCmatchPtVsTrkPID[ti][j], mPtNumVsTrkPID[ti][j], mPtDenVsTrkPID[ti][j]);
1227 setEfficiency(mFractionITSTPCmatchPhiVsTrkPID[ti][j], mPhiNumVsTrkPID[ti][j], mPhiDenVsTrkPID[ti][j]);
1228 setEfficiency(mFractionITSTPCmatchEtaVsTrkPID[ti][j], mEtaNumVsTrkPID[ti][j], mEtaDenVsTrkPID[ti][j]);
1229 }
1230 }
1231 if (mUseMC) {
1232 setEfficiency(mFractionITSTPCmatchPhysPrim[ti], mPtPhysPrimNum[ti], mPtPhysPrimDen[ti]);
1233 setEfficiency(mFractionITSTPCmatchPhiPhysPrim[ti], mPhiPhysPrimNum[ti], mPhiPhysPrimDen[ti]);
1234 setEfficiency(mFractionITSTPCmatchEtaPhysPrim[ti], mEtaPhysPrimNum[ti], mEtaPhysPrimDen[ti]);
1235 setEfficiency(mFractionITSTPCmatchPhysPrim1OverPt[ti], m1OverPtPhysPrimNum[ti], m1OverPtPhysPrimDen[ti]);
1236 }
1237 }
1238 setEfficiency<2>(mFractionITSTPCmatchDCArVsPt, mDCArVsPtNum, mDCArVsPtDen);
1239 /*
1240 mPtTPC->Scale(scaleFactTPC);
1241 mPt->Scale(scaleFactITSTPC);
1242 mPhiTPC->Scale(scaleFactTPC);
1243 mPhi->Scale(scaleFactITSTPC);
1244 if (mUseMC) {
1245 mPtTPCPhysPrim->Scale(scaleFactTPC);
1246 mPtPhysPrim->Scale(scaleFactITSTPC);
1247 mPhiTPCPhysPrim->Scale(scaleFactTPC);
1248 mPhiPhysPrim->Scale(scaleFactITSTPC);
1249 }
1250 mEta->Scale(scaleFactITSTPC);
1251 mChi2Matching->Scale(scaleFactITSTPC);
1252 mChi2Refit->Scale(scaleFactITSTPC);
1253 //mTimeResVsPt->Scale(scaleFactITSTPC); // if to few entries, one sees nothing after normalization --> let's not normalize
1254 */
1255}
1256
1257//__________________________________________________________
1258template <int DIM, bool DEBUG>
1259void MatchITSTPCQC::setEfficiency(TEfficiency* eff, TH1* hnum, TH1* hden)
1260{
1261 // Trivial check if we initalized
1262 if (eff == nullptr) {
1263 LOG(fatal) << "Cannot get TEfficiency object ";
1264 }
1265 if (hnum == nullptr) {
1266 LOG(fatal) << "Cannot get numerator histogram for TEfficiency object " << eff->GetName();
1267 }
1268 if (hden == nullptr) {
1269 LOG(fatal) << "Cannot get denominator histogram for TEfficiency object " << eff->GetName();
1270 }
1271
1272 // we need to force to replace the total histogram, otherwise it will compare it to the previous passed one, and it might get an error of inconsistency in the bin contents
1273 if constexpr (DEBUG) { // checking
1274 bool bad{false};
1275 LOG(info) << "Setting efficiency " << eff->GetName() << " from " << hnum->GetName() << " and " << hden->GetName();
1276 LOG(info) << "Num " << hnum->GetName() << " " << hnum->GetNbinsX() << " " << hnum->GetNbinsY() << " with " << hnum->GetEntries() << " entries";
1277 LOG(info) << "Den " << hden->GetName() << " " << hden->GetNbinsX() << " " << hden->GetNbinsY() << " with " << hden->GetEntries() << " entries";
1278 if (hnum->GetDimension() != hden->GetDimension()) {
1279 LOGP(warning, "Histograms have different dimensions (num={} to den={})", hnum->GetDimension(), hden->GetDimension());
1280 bad = true;
1281 }
1282 if (!TEfficiency::CheckBinning(*hnum, *hden)) {
1283 LOGP(warning, "Histograms do not have a compatible binning");
1284 bad = true;
1285 }
1286 if constexpr (DIM == 3) {
1287 for (int i = 1; i <= hden->GetNbinsX(); i++) {
1288 for (int j = 1; j <= hden->GetNbinsY(); j++) {
1289 for (int k = 1; k <= hden->GetNbinsZ(); k++) {
1290 if (hden->GetBinContent(i, j, k) < hnum->GetBinContent(i, j, k)) {
1291 LOGP(warning, "bin {}/{}/{} -> den: {} < num: {}", i, j, k, hden->GetBinContent(i, j, k), hnum->GetBinContent(i, j, k));
1292 bad = true;
1293 }
1294 }
1295 }
1296 }
1297 } else if constexpr (DIM == 2) {
1298 for (int i = 1; i <= hden->GetNbinsX(); i++) {
1299 for (int j = 1; j <= hden->GetNbinsY(); j++) {
1300 if (hden->GetBinContent(i, j) < hnum->GetBinContent(i, j)) {
1301 LOGP(warning, "bin {}/{} -> den: {} < num: {}", i, j, hden->GetBinContent(i, j), hnum->GetBinContent(i, j));
1302 bad = true;
1303 }
1304 }
1305 }
1306 } else {
1307 for (int i = 1; i <= hden->GetNbinsX(); i++) {
1308 if (hden->GetBinContent(i) < hnum->GetBinContent(i)) {
1309 LOG(warning) << "bin " << i << " den: " << hden->GetBinContent(i) << " < num: " << hnum->GetBinContent(i) << " should be the opposite";
1310 bad = true;
1311 }
1312 }
1313 }
1314 if (bad) {
1315 LOG(info) << " `--> Histogram is bad!";
1316 return;
1317 } else {
1318 LOG(info) << " `--> Histogram is good!";
1319 }
1320 }
1321 // we need to force to replace the total histogram, otherwise it will compare it to the previous passed one, and it might get an error of inconsistency in the bin contents
1322 if (!eff->SetTotalHistogram(*hden, "f")) {
1323 LOG(fatal) << "Something went wrong when defining the efficiency denominator " << eff->GetName() << " from " << hnum->GetName();
1324 }
1325 if (!eff->SetPassedHistogram(*hnum, "")) {
1326 LOG(fatal) << "Something went wrong when defining the efficiency numerator " << eff->GetName() << " from " << hnum->GetName();
1327 }
1328 if constexpr (DIM == 3) {
1329 eff->SetTitle(Form("%s;%s;%s;%s;%s", eff->GetTitle(), hnum->GetXaxis()->GetTitle(), hnum->GetYaxis()->GetTitle(), hnum->GetZaxis()->GetTitle(), "Efficiency"));
1330 } else if constexpr (DIM == 2) {
1331 eff->SetTitle(Form("%s;%s;%s;%s", eff->GetTitle(), hnum->GetXaxis()->GetTitle(), hnum->GetYaxis()->GetTitle(), "Efficiency"));
1332 } else {
1333 eff->SetTitle(Form("%s;%s;%s", eff->GetTitle(), hnum->GetXaxis()->GetTitle(), "Efficiency"));
1334 }
1335}
1336
1337//__________________________________________________________
1338
1339void MatchITSTPCQC::getHistos(TObjArray& objar)
1340{
1341
1342 for (int i = 0; i < matchType::SIZE; ++i) {
1343 objar.Add(mPtNum[i]);
1344 objar.Add(mPtDen[i]);
1345 objar.Add(mFractionITSTPCmatch[i]);
1346
1347 objar.Add(mPtNum_noEta0[i]);
1348 objar.Add(mPtDen_noEta0[i]);
1349 objar.Add(mFractionITSTPCmatch_noEta0[i]);
1350
1351 objar.Add(mPtPhysPrimNum[i]);
1352 objar.Add(mPtPhysPrimDen[i]);
1353 objar.Add(mFractionITSTPCmatchPhysPrim[i]);
1354
1355 objar.Add(mPhiNum[i]);
1356 objar.Add(mPhiDen[i]);
1357 objar.Add(mFractionITSTPCmatchPhi[i]);
1358
1359 if (mUseTrkPID) { // Vs Tracking PID hypothesis
1360 for (int j = 0; j < o2::track::PID::NIDs; ++j) {
1361 // Pt
1362 objar.Add(mPtNumVsTrkPID[i][j]);
1363 objar.Add(mPtDenVsTrkPID[i][j]);
1364 objar.Add(mFractionITSTPCmatchPtVsTrkPID[i][j]);
1365 // Phi
1366 objar.Add(mPhiNumVsTrkPID[i][j]);
1367 objar.Add(mPhiDenVsTrkPID[i][j]);
1368 objar.Add(mFractionITSTPCmatchPhiVsTrkPID[i][j]);
1369 // Eta
1370 objar.Add(mEtaNumVsTrkPID[i][j]);
1371 objar.Add(mEtaDenVsTrkPID[i][j]);
1372 objar.Add(mFractionITSTPCmatchEtaVsTrkPID[i][j]);
1373 }
1374 }
1375
1376 objar.Add(mPhiPhysPrimNum[i]);
1377 objar.Add(mPhiPhysPrimDen[i]);
1378 objar.Add(mFractionITSTPCmatchPhiPhysPrim[i]);
1379
1380 objar.Add(mPhiVsPtNum[i]);
1381 objar.Add(mPhiVsPtDen[i]);
1382 objar.Add(mFractionITSTPCmatchPhiVsPt[i]);
1383
1384 objar.Add(mEtaNum[i]);
1385 objar.Add(mEtaDen[i]);
1386 objar.Add(mFractionITSTPCmatchEta[i]);
1387
1388 objar.Add(mEtaPhysPrimNum[i]);
1389 objar.Add(mEtaPhysPrimDen[i]);
1390 objar.Add(mFractionITSTPCmatchEtaPhysPrim[i]);
1391
1392 objar.Add(mEtaVsPtNum[i]);
1393 objar.Add(mEtaVsPtDen[i]);
1394 objar.Add(mFractionITSTPCmatchEtaVsPt[i]);
1395
1396 objar.Add(mClsVsPtNum[i]);
1397 objar.Add(mClsVsPtDen[i]);
1398 objar.Add(mFractionITSTPCmatchClsVsPt[i]);
1399
1400 objar.Add(mChi2VsPtNum[i]);
1401 objar.Add(mChi2VsPtDen[i]);
1402 objar.Add(mFractionITSTPCmatchChi2VsPt[i]);
1403
1404 objar.Add(m1OverPtNum[i]);
1405 objar.Add(m1OverPtDen[i]);
1406 objar.Add(mFractionITSTPCmatch1OverPt[i]);
1407
1408 objar.Add(m1OverPtPhysPrimNum[i]);
1409 objar.Add(m1OverPtPhysPrimDen[i]);
1410 objar.Add(mFractionITSTPCmatchPhysPrim1OverPt[i]);
1411
1412 objar.Add(mEtaPhiPtNum[i]);
1413 objar.Add(mEtaPhiPtDen[i]);
1414 }
1415 objar.Add(mChi2Matching);
1416 objar.Add(mChi2Refit);
1417 objar.Add(mTimeResVsPt);
1418 objar.Add(mResidualPt);
1419 objar.Add(mResidualPhi);
1420 objar.Add(mResidualEta);
1421 objar.Add(mDCAr);
1422 objar.Add(mDCArVsPtNum);
1423 objar.Add(mDCArVsPtDen);
1424 objar.Add(mFractionITSTPCmatchDCArVsPt);
1425
1426 // V0
1427 objar.Add(mK0MassVsPtVsOccpp);
1428 objar.Add(mK0MassVsPtVsOccPbPb);
1429}
1430
1431void MatchITSTPCQC::printParams() const
1432{
1433 LOG(info) << "MatchITSTPCQC parameters:";
1434 LOG(info) << " - minPtBins = " << mPtBins;
1435 LOG(info) << " - minPtITSCut = " << mPtITSCut;
1436 LOG(info) << " - etaITSCut = " << mEtaITSCut;
1437 LOG(info) << " - minNITSClustersCut = " << mMinNClustersITS;
1438 LOG(info) << " - maxChi2PerClusterITS = " << mMaxChi2PerClusterITS;
1439 LOG(info) << " - minPtTPCCut = " << mPtTPCCut;
1440 LOG(info) << " - etaTPCCut = " << mEtaTPCCut;
1441 LOG(info) << " - minNTPCClustersCut = " << mNTPCClustersCut;
1442 LOG(info) << " - mEtaNo0Cut = " << mEtaNo0Cut;
1443 LOG(info) << " - minDCACut = " << mDCATPCCut;
1444 LOG(info) << " - minDCACutY = " << mDCATPCCutY;
1445 LOG(info) << " - minPtCut = " << mPtCut;
1446 LOG(info) << " - maxPtCut = " << mPtMaxCut;
1447 LOG(info) << " - etaCut = " << mEtaCut;
1448 LOG(info) << " - cutK0Mass = " << mCutK0Mass;
1449 LOG(info) << " - maxEtaK0 = " << mMaxEtaK0;
1450 LOG(info) << " - minTPCOccpp = " << mMinTPCOccpp;
1451 LOG(info) << " - maxTPCOccpp = " << mMaxTPCOccpp;
1452 LOG(info) << " - nBinsTPCOccpp = " << mNBinsTPCOccpp;
1453 LOG(info) << " - minTPCOccPbPb = " << mMinTPCOccPbPb;
1454 LOG(info) << " - maxTPCOccPbPb = " << mMaxTPCOccPbPb;
1455 LOG(info) << " - nBinsTPCOccPbPb = " << mNBinsTPCOccPbPb;
1456}
DEBUG
#define DEBUG
Definition AlignmentTrack.cxx:31
TrackCuts.h
Class to perform track cuts.
debug
std::ostringstream debug
Definition ExpressionJSONHelpers.cxx:83
MatchITSTPCQC.h
Class to perform QC for ITSTPC matching.
i
int32_t i
Definition GPUCommonAlgorithm.h:459
GRPGeomHelper.h
Helper for geometry and GRP related CCDB requests.
labels
std::vector< o2::MCCompLabel > labels
Definition CATrackerSpec.cxx:79
LHCConstants.h
Header to collect LHC related constants.
MCUtils.h
Utility functions for MC particles.
j
uint32_t j
Definition RawData.h:0
SVertexerParams.h
Configurable params for secondary vertexer.
TrackTPCITS.h
Result of refitting TPC-ITS matched track.
o2::TrackCuts::setMinPtCut
void setMinPtCut(float value)
ITS+TPC kinematics.
Definition TrackCuts.h:73
o2::TrackCuts::setMinNClustersITS
void setMinNClustersITS(int32_t value)
Definition TrackCuts.h:59
o2::TrackCuts::setMinNTPCClustersCut
void setMinNTPCClustersCut(int32_t value)
Definition TrackCuts.h:69
o2::TrackCuts::isSelected
bool isSelected(GID trackIndex, o2::globaltracking::RecoContainer &data)
Definition TrackCuts.h:82
o2::TrackCuts::setMaxPtCut
void setMaxPtCut(float value)
Definition TrackCuts.h:74
o2::TrackCuts::setEtaITSCut
void setEtaITSCut(float value)
Definition TrackCuts.h:58
o2::TrackCuts::setMaxChi2PerClusterITS
void setMaxChi2PerClusterITS(float value)
Definition TrackCuts.h:60
o2::TrackCuts::setEtaTPCCut
void setEtaTPCCut(float value)
Definition TrackCuts.h:68
o2::TrackCuts::setMaxDCATPCCut
void setMaxDCATPCCut(float value)
Definition TrackCuts.h:70
o2::TrackCuts::setMinPtITSCut
void setMinPtITSCut(float value)
Definition TrackCuts.h:57
o2::TrackCuts::setMinPtTPCCut
void setMinPtTPCCut(float value)
TPC.
Definition TrackCuts.h:67
o2::TrackCuts::setEtaCut
void setEtaCut(float valueMin, float valueMax)
Definition TrackCuts.h:75
o2::TrackCuts::setMaxDCATPCCutY
void setMaxDCATPCCutY(float value)
Definition TrackCuts.h:71
o2::TrackCuts::setRequireHitsInITSLayers
void setRequireHitsInITSLayers(int8_t minNRequiredHits, std::set< uint8_t > requiredLayers)
Definition TrackCuts.h:61
o2::base::GRPGeomHelper::instance
static GRPGeomHelper & instance()
Definition GRPGeomHelper.h:120
o2::base::PropagatorImpl< float >::Instance
GPUd() value_type estimateLTFast(o2 static GPUd() float estimateLTIncrement(const o2 PropagatorImpl * Instance(bool uninitialized=false)
Definition Propagator.h:180
o2::dataformats::V0::setDCA
void setDCA(float d)
Definition V0.h:47
o2::detectors::DetID
Static class with identifiers, bitmasks and names for ALICE detectors.
Definition DetID.h:60
o2::gloqc::MatchITSTPCQC::processV0
bool processV0(int iv, o2::globaltracking::RecoContainer &recoData, std::vector< float > &mTBinClOcc, float pvTime)
Definition MatchITSTPCQC.cxx:1110
o2::gloqc::MatchITSTPCQC::refitV0
bool refitV0(const o2::dataformats::V0Index &id, o2::dataformats::V0 &v0, o2::globaltracking::RecoContainer &recoData)
Definition MatchITSTPCQC.cxx:1157
o2::gloqc::MatchITSTPCQC::getHistos
void getHistos(TObjArray &objar)
Definition MatchITSTPCQC.cxx:1339
o2::gloqc::MatchITSTPCQC::run
void run(o2::framework::ProcessingContext &ctx)
Definition MatchITSTPCQC.cxx:481
o2::gpu::GPUO2InterfaceRefit::fillOccupancyMapGetSize
static size_t fillOccupancyMapGetSize(uint32_t nHbfPerTf, const GPUParam *param=nullptr)
Definition GPUO2InterfaceRefit.cxx:80
o2::gpu::GPUO2InterfaceUtils::paramUseExternalOccupancyMap
static void paramUseExternalOccupancyMap(GPUParam *param, uint32_t nHbfPerTf, const uint32_t *occupancymap, int32_t occupancyMapSize)
Definition GPUO2InterfaceUtils.cxx:123
o2::gpu::GPUO2InterfaceUtils::getFullParamShared
static std::shared_ptr< GPUParam > getFullParamShared(float solenoidBz, uint32_t nHbfPerTf=0, std::unique_ptr< GPUO2InterfaceConfiguration > *pConfiguration=nullptr, std::unique_ptr< GPUSettingsO2 > *pO2Settings=nullptr, bool *autoMaxTimeBin=nullptr)
Definition GPUO2InterfaceUtils.cxx:118
o2::mcutils::MCTrackNavigator::isPhysicalPrimary
static bool isPhysicalPrimary(o2::MCTrack const &p, std::vector< o2::MCTrack > const &pcontainer)
Definition MCUtils.cxx:73
o2::steer::MCKinematicsReader::initFromDigitContext
bool initFromDigitContext(std::string_view filename)
Definition MCKinematicsReader.cxx:212
o2::steer::MCKinematicsReader::releaseTracksForSourceAndEvent
void releaseTracksForSourceAndEvent(int source, int event)
API to ask releasing tracks and track references (freeing memory) for source + event.
Definition MCKinematicsReader.cxx:115
o2::steer::MCKinematicsReader::getTracks
std::vector< MCTrack > const & getTracks(int source, int event) const
variant returning all tracks for source and event at once
Definition MCKinematicsReader.h:175
o2::vertexing::DCAFitterN::getTrack
Track & getTrack(int i, int cand=0)
Definition DCAFitterN.h:221
o2::vertexing::DCAFitterN::calcPCACovMatrixFlat
std::array< float, 6 > calcPCACovMatrixFlat(int cand=0) const
Definition DCAFitterN.h:256
o2::vertexing::DCAFitterN::getChi2AtPCACandidate
float getChi2AtPCACandidate(int cand=0) const
Definition DCAFitterN.h:207
event
struct _cl_event * event
Definition glcorearb.h:2982
v
const GLdouble * v
Definition glcorearb.h:832
source
GLsizei GLsizei GLchar * source
Definition glcorearb.h:798
v0
GLfloat v0
Definition glcorearb.h:811
flags
GLbitfield flags
Definition glcorearb.h:1570
id
GLuint id
Definition glcorearb.h:650
o2::track::TrackParCov
TrackParCovF TrackParCov
Definition Track.h:33
o2::globaltracking::RecoContainer::getTrackMCLabel
o2::MCCompLabel getTrackMCLabel(GTrackID id) const
Definition RecoContainer.h:505
o2::globaltracking::RecoContainer::getTrackParam
const o2::track::TrackParCov & getTrackParam(GTrackID gidx) const
Definition RecoContainer.cxx:1402
o2::globaltracking::RecoContainer::collectData
void collectData(o2::framework::ProcessingContext &pc, const DataRequest &request)
Definition RecoContainer.cxx:639
o2::globaltracking::RecoContainer::getPrimaryVertex
const o2::dataformats::PrimaryVertex & getPrimaryVertex(int i) const
Definition RecoContainer.h:720
o2::globaltracking::RecoContainer::occupancyMapTPC
gsl::span< const unsigned int > occupancyMapTPC
externally set TPC clusters occupancy map
Definition RecoContainer.h:375
o2::vertexing::SVertexerParams::useAbsDCA
bool useAbsDCA
use abs dca minimization
Definition SVertexerParams.h:38
o2::vertexing::SVertexerParams::usePropagator
bool usePropagator
use external propagator
Definition SVertexerParams.h:53
o2::vertexing::SVertexerParams::maxRIni
float maxRIni
don't consider as a seed (circles intersection) if its R exceeds this
Definition SVertexerParams.h:45
o2::vertexing::SVertexerParams::minParamChange
float minParamChange
stop when tracks X-params being minimized change by less that this value
Definition SVertexerParams.h:41
o2::vertexing::SVertexerParams::maxSnp
float maxSnp
max snp when external propagator is used
Definition SVertexerParams.h:51
o2::vertexing::SVertexerParams::minXSeed
float minXSeed
minimal X of seed in prong frame (within the radial resolution track should not go to negative X)
Definition SVertexerParams.h:52
o2::vertexing::SVertexerParams::matCorr
int matCorr
material correction to use
Definition SVertexerParams.h:48
o2::vertexing::SVertexerParams::maxDZIni
float maxDZIni
don't consider as a seed (circles intersection) if Z distance exceeds this
Definition SVertexerParams.h:43
o2::vertexing::SVertexerParams::maxStep
float maxStep
max step size when external propagator is used
Definition SVertexerParams.h:50
o2::vertexing::SVertexerParams::minRelChi2Change
float minRelChi2Change
stop when chi2 changes by less than this value
Definition SVertexerParams.h:42
o2::vertexing::SVertexerParams::maxChi2
float maxChi2
max dca from prongs to vertex
Definition SVertexerParams.h:40
o2::vertexing::SVertexerParams::refitWithMatCorr
bool refitWithMatCorr
refit V0 applying material corrections
Definition SVertexerParams.h:54
o2::vertexing::SVertexerParams::maxDXYIni
float maxDXYIni
don't consider as a seed (circles intersection) if XY distance exceeds this
Definition SVertexerParams.h:44
LOG
LOG(info)<< "Compressed in "<< sw.CpuTime()<< " s"