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