3 //**************************************************************************
4 //* This file is property of and copyright by the ALICE HLT Project *
5 //* ALICE Experiment at CERN, All rights reserved. *
7 //* Primary Authors: Zhong-Bao Yin <Zhong-Bao.Yin@cern.ch> *
8 //* Matthias Richter <Matthias.Richter@ift.uib.no> *
9 //* for The ALICE HLT Project. *
11 //* Permission to use, copy, modify and distribute this software and its *
12 //* documentation strictly for non-commercial purposes is hereby granted *
13 //* without fee, provided that the above copyright notice appears in all *
14 //* copies and that both the copyright notice and this permission notice *
15 //* appear in the supporting documentation. The authors make no claims *
16 //* about the suitability of this software for any purpose. It is *
17 //* provided "as is" without express or implied warranty. *
18 //**************************************************************************
20 /** @file AliHLTTPCQADataMaker.cxx
21 @author Zhongbao Yin, Matthias Richter
23 @brief Container for the HLT offline QA
25 #include "AliHLTTPCQADataMaker.h"
26 #include "AliESDEvent.h"
32 #include "AliESDtrack.h"
36 /** ROOT macro for the implementation of ROOT specific class methods */
37 ClassImp(AliHLTTPCQADataMaker)
39 AliHLTTPCQADataMaker::AliHLTTPCQADataMaker()
40 : AliHLTQADataMakerBase()
42 // see header file for class documentation
44 // refer to README to build package
46 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
49 AliHLTTPCQADataMaker::~AliHLTTPCQADataMaker()
51 // see header file for class documentation
54 void AliHLTTPCQADataMaker::StartOfDetectorCycle()
56 // see header file for class documentation
59 void AliHLTTPCQADataMaker::EndOfDetectorCycle(AliQAv1::TASKINDEX_t, TObjArray** /*list*/)
61 // see header file for class documentation
63 if(GetESDsData(kPHLT)->GetEntries()){
64 GetESDsData(kPRatio)->Sumw2();
65 GetESDsData(kPRatio)->Add(GetESDsData(kPOffline));
66 GetESDsData(kPRatio)->Divide(GetESDsData(kPHLT));
69 if(GetESDsData(kPHLTFired)->GetEntries()){
70 GetESDsData(kPRatioFired)->Sumw2();
71 GetESDsData(kPRatioFired)->Add(GetESDsData(kPOfflineFired));
72 GetESDsData(kPRatioFired)->Divide(GetESDsData(kPHLTFired));
78 void AliHLTTPCQADataMaker::MakeRaws(AliRawReader * rawReader)
80 // see header file for class documentation
81 if (!rawReader) return;
84 void AliHLTTPCQADataMaker::InitESDs(){
86 //create ESDs histograms in ESDs subdir
87 const Bool_t expert = kTRUE ;
88 const Bool_t image = kTRUE ;
90 TH2F * histESDMultiplicity =
91 new TH2F("hESDMultiplicity",
92 "Number of tracks from all events; Number of offline tracks; Number of HLT tracks",
93 300, 0, 300, 300, 0, 300);
94 histESDMultiplicity->Sumw2();
95 Add2ESDsList(histESDMultiplicity, kMultiplicity, !expert, image);
97 TH2F * histESDMultiplicityFired =
98 new TH2F("hESDMultiplicityFired",
99 "Number of tracks from HLT triggered events; Number of offline tracks; Number of HLT tracks",
100 300, 0, 300, 300, 0, 300);
101 histESDMultiplicityFired->Sumw2();
102 Add2ESDsList(histESDMultiplicityFired, kMultiplicityFired, !expert, image);
105 new TH2F("hESDNCls", "Mean number of TPC clusters from all events; Offline; HLT",
106 200, 0, 200, 200, 0, 200);
107 histESDNCls->Sumw2();
108 Add2ESDsList(histESDNCls, kNCls, !expert, image);
110 TH2F * histESDNClsFired =
111 new TH2F("hESDNClsFired", "Mean number of TPC clusters from triggered events; Offline; HLT",
112 200, 0, 200, 200, 0, 200);
113 histESDNClsFired->Sumw2();
114 Add2ESDsList(histESDNClsFired, kNClsFired, !expert, image);
117 = new TH1F("hPHLT", "P distribution for all events from HLT; P [GeV/c]",
120 Add2ESDsList(histPHLT, kPHLT, !expert, image);
123 = new TH1F("hPOffline",
124 "P distribution for all events from offline; P [GeV/c]",
126 histPOffline->Sumw2();
127 Add2ESDsList(histPOffline, kPOffline, !expert, image);
130 = new TH1F("hPHLTFired",
131 "P distribution for fired events from HLT; P [GeV/c]",
133 histPHLTFired->Sumw2();
134 Add2ESDsList(histPHLTFired, kPHLTFired, !expert, image);
136 TH1F * histPOfflineFired
137 = new TH1F("hPOfflineFired",
138 "P distribution for fired events from offline; P [GeV/c]",
140 histPOfflineFired->Sumw2();
141 Add2ESDsList(histPOfflineFired, kPOfflineFired, !expert, image);
144 = new TH1F("hPRatio",
145 "Ratio of P distribution for all events; P [GeV/c]",
148 Add2ESDsList(histPRatio, kPRatio, !expert, image);
150 TH1F * histPRatioFired
151 = new TH1F("hPRatioFired",
152 "Ratio of P distribution for fired events; P [GeV/c]",
154 histPRatioFired->Sumw2();
155 Add2ESDsList(histPRatioFired, kPRatioFired, !expert, image);
159 "P_{T} distribution for all events from HLT; P_{T} [GeV/c]",
162 Add2ESDsList(histPtHLT, kPtHLT, !expert, image);
165 = new TH1F("hPtOffline",
166 "P_{T} distribution for all events from offline; P_{T} [GeV/c]",
168 histPtOffline->Sumw2();
169 Add2ESDsList(histPtOffline, kPtOffline, !expert, image);
171 TH1F * histPtHLTFired
172 = new TH1F("hPtHLTFired",
173 "P_{T} distribution for fired events from HLT; P_{T} [GeV/c]",
175 histPtHLTFired->Sumw2();
176 Add2ESDsList(histPtHLTFired, kPtHLTFired, !expert, image);
178 TH1F * histPtOfflineFired
179 = new TH1F("hPtOfflineFired",
180 "P_{T} distribution for fired events from offline; P_{T} [GeV/c]",
182 histPtOfflineFired->Sumw2();
183 Add2ESDsList(histPtOfflineFired, kPtOfflineFired, !expert, image);
185 TH1F * histNClsPerTrkHLT
186 = new TH1F("hNClsPerTrkHLT", "Clusters per HLT track; N cluster; Counts",
188 histNClsPerTrkHLT->Sumw2();
189 Add2ESDsList(histNClsPerTrkHLT, kNClsPerTrkHLT, !expert, image);
191 TH1F * histNClsPerTrkOffline
192 = new TH1F("hNClsPerTrkOffline",
193 "Clusters per offline track; N cluster; Counts",
195 histNClsPerTrkOffline->Sumw2();
196 Add2ESDsList(histNClsPerTrkOffline, kNClsPerTrkOffline, !expert, image);
198 TH1F * histNClsPerTrkHLTFired
199 = new TH1F("hNClsPerTrkHLTFired",
200 "Clusters per HLT track from HLT fired events; N cluster; Counts",
202 histNClsPerTrkHLTFired->Sumw2();
203 Add2ESDsList(histNClsPerTrkHLTFired, kNClsPerTrkHLTFired, !expert, image);
205 TH1F * histNClsPerTrkOfflineFired
206 = new TH1F("hNClsPerTrkOfflineFired",
207 "Clusters per offline track from HLT fired events; N cluster; Counts",
209 histNClsPerTrkOfflineFired->Sumw2();
210 Add2ESDsList(histNClsPerTrkOfflineFired, kNClsPerTrkOfflineFired, !expert, image);
213 new TH1F("hPhiHLT", "Phi distribution of HLT tracks; Phi; Counts",
216 Add2ESDsList(histPhiHLT, kPhiHLT, !expert, image);
218 TH1F * histPhiOffline =
219 new TH1F("hPhiOffline",
220 "Phi distribution of offline tracks; Phi; Counts",
222 histPhiOffline->Sumw2();
223 Add2ESDsList(histPhiOffline, kPhiOffline, !expert, image);
225 TH1F * histPhiHLTFired =
226 new TH1F("hPhiHLTFired", "Phi distribution of HLT tracks from HLT fired event ; Phi; Counts",
228 histPhiHLTFired->Sumw2();
229 Add2ESDsList(histPhiHLTFired, kPhiHLTFired, !expert, image);
231 TH1F * histPhiOfflineFired =
232 new TH1F("hPhiOfflineFired",
233 "Phi distribution of offline tracks from HLT fired events; Phi; Counts",
235 histPhiOfflineFired->Sumw2();
236 Add2ESDsList(histPhiOfflineFired, kPhiOfflineFired, !expert, image);
239 new TH1F("hEtaHLT", "Eta distribution of HLT tracks; Eta; Counts",
242 Add2ESDsList(histEtaHLT, kEtaHLT, !expert, image);
244 TH1F * histEtaOffline =
245 new TH1F("hEtaOffline", "Eta distribution of offline tracks; Eta; Counts",
248 Add2ESDsList(histEtaOffline, kEtaOffline, !expert, image);
250 TH1F * histEtaHLTFired =
251 new TH1F("hEtaHLTFired",
252 "Eta distribution of HLT tracks from HLT fired events; Eta; Counts",
254 histEtaHLTFired->Sumw2();
255 Add2ESDsList(histEtaHLTFired, kEtaHLTFired, !expert, image);
257 TH1F * histEtaOfflineFired =
258 new TH1F("hEtaOfflineFired",
259 "Eta distribution of offline tracks from HLT fired events; Eta; Counts",
261 histEtaHLTFired->Sumw2();
262 Add2ESDsList(histEtaOfflineFired, kEtaOfflineFired, !expert, image);
266 void AliHLTTPCQADataMaker::MakeESDs(AliESDEvent * esd, AliESDEvent* hltesd)
269 if (!esd || !hltesd) {
270 AliError("invalid parameter: missing ESDs");
274 // make QA data from ESDs
276 const Int_t nESDTracks = esd->GetNumberOfTracks();
277 const Int_t nHLTesdTracks = hltesd->GetNumberOfTracks();
278 GetESDsData(kMultiplicity)->Fill(nESDTracks, nHLTesdTracks);
281 Int_t nClsOffline = 0;
283 for(Int_t i = 0; i < nESDTracks; i++){
284 AliESDtrack * esdTrk = esd->GetTrack(i);
285 GetESDsData(kPOffline)->Fill(esdTrk->P());
286 GetESDsData(kPtOffline)->Fill(esdTrk->GetSignedPt());
287 GetESDsData(kNClsPerTrkOffline)->Fill(esdTrk->GetTPCNcls());
288 GetESDsData(kPhiOffline)->Fill(esdTrk->Phi()*TMath::RadToDeg());
289 GetESDsData(kEtaOffline)->Fill(esdTrk->Eta());
290 nClsOffline+=esdTrk->GetTPCNcls();
293 for(Int_t i = 0; i < nHLTesdTracks; i++){
294 AliESDtrack * hltEsdTrk = hltesd->GetTrack(i);
295 GetESDsData(kPHLT)->Fill(hltEsdTrk->P());
296 GetESDsData(kPtHLT)->Fill(hltEsdTrk->GetSignedPt());
297 GetESDsData(kNClsPerTrkHLT)->Fill(hltEsdTrk->GetTPCNcls());
298 GetESDsData(kPhiHLT)->Fill(hltEsdTrk->Phi()*TMath::RadToDeg());
299 GetESDsData(kEtaHLT)->Fill(hltEsdTrk->Eta());
300 nClsHLT += hltEsdTrk->GetTPCNcls();
304 nClsOffline /= Float_t(nESDTracks);
306 nClsHLT /= Float_t(nHLTesdTracks);
308 GetESDsData(kNCls)->Fill(nClsOffline, nClsHLT);
311 if(hltesd->IsHLTTriggerFired()){
312 GetESDsData(kMultiplicityFired)->Fill(nESDTracks, nHLTesdTracks);
313 GetESDsData(kNClsFired)->Fill(nClsOffline, nClsHLT);
315 for(Int_t i = 0; i < nESDTracks; i++){
316 AliESDtrack * esdTrk = esd->GetTrack(i);
317 GetESDsData(kPOfflineFired)->Fill(esdTrk->P());
318 GetESDsData(kPtOfflineFired)->Fill(esdTrk->GetSignedPt());
319 GetESDsData(kNClsPerTrkOfflineFired)->Fill(esdTrk->GetTPCNcls());
320 GetESDsData(kPhiOfflineFired)->Fill(esdTrk->Phi()*TMath::RadToDeg());
321 GetESDsData(kEtaOfflineFired)->Fill(esdTrk->Eta());
324 for(Int_t i = 0; i < nHLTesdTracks; i++){
325 AliESDtrack * hltEsdTrk = hltesd->GetTrack(i);
326 GetESDsData(kPHLTFired)->Fill(hltEsdTrk->P());
327 GetESDsData(kPtHLTFired)->Fill(hltEsdTrk->GetSignedPt());
328 GetESDsData(kNClsPerTrkHLTFired)->Fill(hltEsdTrk->GetTPCNcls());
329 GetESDsData(kPhiHLTFired)->Fill(hltEsdTrk->Phi()*TMath::RadToDeg());
330 GetESDsData(kEtaHLTFired)->Fill(hltEsdTrk->Eta());