- fix names and units of axes
[u/mrichter/AliRoot.git] / HLT / QA / tasks / AliAnalysisTaskHLTCentralBarrel.cxx
CommitLineData
ff84363f 1// $Id$
2//**************************************************************************
3//* This file is property of and copyright by the ALICE HLT Project *
4//* ALICE Experiment at CERN, All rights reserved. *
5//* *
6//* Primary Authors: Kalliopi Kanaki <Kalliopi.Kanaki@ift.uib.no> *
7//* for The ALICE HLT Project. *
8//* *
9//* Permission to use, copy, modify and distribute this software and its *
10//* documentation strictly for non-commercial purposes is hereby granted *
11//* without fee, provided that the above copyright notice appears in all *
12//* copies and that both the copyright notice and this permission notice *
13//* appear in the supporting documentation. The authors make no claims *
14//* about the suitability of this software for any purpose. It is *
15//* provided "as is" without express or implied warranty. *
16//**************************************************************************
17
18
19/** @file AliAnalysisTaskHLTCentralBarrel.cxx
20 @author Per Ivar L√łnne, Hege Erdal, Kalliopi Kanaki
21 @date
22 @brief An analysis task containing
23 loops over HLT and offline ESD trees for comparing
24 event and track properties
25*/
26
27#include <iostream>
28
29#include "AliAnalysisTaskHLTCentralBarrel.h"
30#include "AliESDEvent.h"
31#include "AliESDtrack.h"
32#include "AliESDInputHandler.h"
33#include "AliTracker.h"
34#include "AliESDVZERO.h"
35#include "AliHLTGlobalTriggerDecision.h"
36
37#include "TAxis.h"
38#include "TString.h"
39#include "TText.h"
40#include "TTimeStamp.h"
2de8763b 41#include "TH1.h"
ff84363f 42
43ClassImp(AliAnalysisTaskHLTCentralBarrel)
44//_______________________________________________________________________________________________//
45
46AliAnalysisTaskHLTCentralBarrel::AliAnalysisTaskHLTCentralBarrel()
47:AliAnalysisTaskSE()
48 ,fUseHLTTrigger(kFALSE)
49 ,fOutputList(0)
50 ,fEventOFF(0)
51 ,fEventHLT(0)
52 ,fTrackOFF(0)
53 ,fTrackHLT(0)
54 ,fTextBox(0)
55{
56 // Constructor
57 // Define input and output slots here
58 // Input slot #0 works with a TChain
59 // DefineInput(0, TChain::Class());
60 // Output slot #0 writes into a TH1 container
61
62 //DefineOutput(1, TList::Class());
63}
64
65AliAnalysisTaskHLTCentralBarrel::AliAnalysisTaskHLTCentralBarrel(const char *name)
66:AliAnalysisTaskSE(name)
67 ,fUseHLTTrigger(kFALSE)
68 ,fOutputList(0)
69 ,fEventOFF(0)
70 ,fEventHLT(0)
71 ,fTrackOFF(0)
72 ,fTrackHLT(0)
73 ,fTextBox(0)
74{
75 // Constructor
76 // Define input and output slots here
77 // Input slot #0 works with a TChain
78 //DefineInput(0, TChain::Class());
79 // Output slot #0 writes into a TH1 container
80
81 DefineOutput(1, TList::Class());
82}
83
84AliAnalysisTaskHLTCentralBarrel::~AliAnalysisTaskHLTCentralBarrel(){
85 // destructor
86}
87
88void AliAnalysisTaskHLTCentralBarrel::UserCreateOutputObjects(){
89 // Create THnSpare objects
90
91 OpenFile(1);
92 fOutputList = new TList();
93 fOutputList->SetName(GetName());
94
95 static const int sizeEvent = 6;
96
97 int binsEvent[sizeEvent] = { 200, 200, 250, 2000, 2000, 2 };
98 double minEvent [sizeEvent] = { -2, -2, -30, 0, 0, 0 };
99 double maxEvent [sizeEvent] = { 2, 2, 30, 20000, 20000, 2 };
100
101 fEventHLT = CreateEventTHnSparse("fEventHLT",sizeEvent,binsEvent,minEvent,maxEvent);
102 fEventOFF = CreateEventTHnSparse("fEventOFF",sizeEvent,binsEvent,minEvent,maxEvent);
103
104 static const int sizeTrack = 13;
105
106 Int_t binsTrack[sizeTrack] = {1500, 200, 200, 200, 200, 400, 400, 3, 400, 400, 6, 2000, 2 };
107 Double_t minTrack [sizeTrack] = { 0, 0, -1, -3, -1, -100, -100, -1.5, -100, -100, 0, 0, 0 };
108 Double_t maxTrack [sizeTrack] = { 150, 200, 4, 3, 7, 100, 100, 1.5, 100, 100, 6, 20000, 2 };
109
110 fTrackHLT = CreateTrackTHnSparse("fTrackHLT",sizeTrack,binsTrack,minTrack,maxTrack);
111 fTrackOFF = CreateTrackTHnSparse("fTrackOFF",sizeTrack,binsTrack,minTrack,maxTrack);
112
113 fTextBox = new TText();
114 fOutputList->Add(fEventOFF);
115 fOutputList->Add(fEventHLT);
116 fOutputList->Add(fTrackOFF);
117 fOutputList->Add(fTrackHLT);
118 fOutputList->Add(fTextBox);
119}
120
121void AliAnalysisTaskHLTCentralBarrel::NotifyRun(){
122
123 AliESDEvent *esdOFF = dynamic_cast<AliESDEvent*>(InputEvent());
124 TTimeStamp* timestamp = new TTimeStamp(esdOFF->GetTimeStamp());
125 fTextBox->SetName("text");
126
127 TString s = Form("Run %d, Date %d", esdOFF->GetRunNumber(), timestamp->GetDate());
128 fTextBox->SetTitle(s);
129}
130
131void AliAnalysisTaskHLTCentralBarrel::UserExec(Option_t *){
132// see header for documentation
133
134 AliESDEvent *esdOFF = dynamic_cast<AliESDEvent*>(InputEvent());
135 if(!esdOFF){
136 printf("ERROR: offline ESD is not available.\n");
137 return;
138 }
139
140 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*>(fInputHandler);
141 AliESDEvent *esdHLT = NULL;
142
143 if(esdH) esdHLT = esdH->GetHLTEvent();
144 if(!esdHLT){
145 printf("ERROR: HLT ESD is not available.\n");
146 return;
147 }
148
149 if(fUseHLTTrigger && !((AliHLTGlobalTriggerDecision*)esdHLT->GetHLTTriggerDecision())->Result()) return;
150
151 if(fUseCentrality){
152 Int_t centbin = CalculateCentrality(esdHLT);
153 Printf("Centrality bin = %d", centbin);
154 }
155
156
157
158 //============================ OFFLINE =============================//
159
160 const AliESDVertex *vertOFF = esdOFF->GetPrimaryVertexTracks();
161
162 Double_t bfield = esdOFF->GetMagneticField();
163 Int_t nr_tracksOFF = 0;
164
165 if(esdOFF->GetEventSpecie()==16) return;
166
167 for(Int_t i=0; i<esdOFF->GetNumberOfTracks(); i++){
168 AliESDtrack *esdTrackOFF = esdOFF->GetTrack(i);
169 if (!esdTrackOFF) continue;
170 if(!(esdTrackOFF->GetStatus()&AliESDtrack::kTPCin)) continue;
171 nr_tracksOFF++;
172 }
173
174 for(Int_t i=0; i<esdOFF->GetNumberOfTracks(); i++){
175
176 AliESDtrack *esdTrackOFF = esdOFF->GetTrack(i);
177 if (!esdTrackOFF) continue;
178 if(!(esdTrackOFF->GetStatus()&AliESDtrack::kTPCin)) continue;
179
180 //DCA calculations(from offline)
181 Double_t x[3];
182 Double_t b[3];
183 esdTrackOFF->GetXYZ(x);
184
185 AliTracker::GetBxByBz(x,b);
186 Bool_t isOK = esdTrackOFF->RelateToVertexTPCBxByBz(vertOFF, b, kVeryBig);
187 if(!isOK) continue;
188
189 const AliExternalTrackParam *track = esdTrackOFF->GetTPCInnerParam();
190 if(!track) continue;
191
192 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
193 if(vertOFF->GetX()==0 && vertOFF->GetY()==0 && vertOFF->GetZ()==0 ){
194 dca[0]=-99;
195 dca[1]=-99;
196 }
197 else esdTrackOFF->GetImpactParametersTPC(dca,cov);
198
199 Float_t DCAr =-99, DCAz = -99.;
200 Double_t trackOFF[] = {
201 TMath::Abs(esdTrackOFF->Pt())
202 ,esdTrackOFF->GetTPCNcls()
203 ,esdTrackOFF->Theta()
204 ,esdTrackOFF->Eta()
205 ,esdTrackOFF->Phi()
206 ,dca[0]
207 ,dca[1]
208 ,esdTrackOFF->Charge()
209 ,DCAr
210 ,DCAz
211 ,esdTrackOFF->GetNcls(0)
212 ,nr_tracksOFF
213 ,vertOFF->GetStatus()
214 };
215 fTrackOFF->Fill(trackOFF);
216 }
217
218 Double_t eventOFF[] = { vertOFF->GetX(), vertOFF->GetY(), vertOFF->GetZ(), vertOFF->GetNContributors(), nr_tracksOFF, vertOFF->GetStatus()};
219 fEventOFF->Fill(eventOFF);
220
221
222
223 //======================================== HLT ==========================================//
224
225 Int_t nr_tracksHLT = 0;
226 if(esdHLT->GetEventSpecie()==16) return;
227 const AliESDVertex *vertHLT = esdHLT->GetPrimaryVertexTracks();
228
229 for(Int_t i=0; i<esdHLT->GetNumberOfTracks(); i++){
230 AliESDtrack *esdTrackHLT = esdHLT->GetTrack(i);
231 if (!esdTrackHLT) continue;
232 if(!(esdTrackHLT->GetStatus()&AliESDtrack::kTPCin)) continue;
233 nr_tracksHLT++;
234 }
235
236 Double_t eventHLT[] = { vertHLT->GetX(), vertHLT->GetY(), vertHLT->GetZ(), vertHLT->GetNContributors(), nr_tracksHLT, vertHLT->GetStatus()};
237 fEventHLT->Fill(eventHLT);
238
239 for(Int_t i=0; i<esdHLT->GetNumberOfTracks(); i++){
240
241 AliESDtrack *esdTrackHLT = esdHLT->GetTrack(i);
242 if(!esdTrackHLT) continue;
243 if(!(esdTrackHLT->GetStatus()&AliESDtrack::kTPCin)) continue;
244
245 //DCA calculations
246 Float_t DCAr=-99;
247 Float_t DCAz=-99;
248 Float_t dca[2];
249 if(vertHLT->GetX()==0 && vertHLT->GetY()==0 && vertHLT->GetZ() ==0 ){
250 dca[0]=-99;
251 dca[1]=-99;
252 }
253 else{
254 //Calculating DCA "old" fashion
255 esdTrackHLT->GetDZ(esdHLT->GetPrimaryVertex()->GetXv(), esdHLT->GetPrimaryVertex()->GetYv(), esdHLT->GetPrimaryVertex()->GetZv(), bfield, dca);
256 // plotting the DCA calculated by Sergey
257 esdTrackHLT->GetImpactParametersTPC(DCAr,DCAz);
258 }
259
260 Double_t trackHLT[] = {
261 TMath::Abs(esdTrackHLT->Pt())
262 ,esdTrackHLT->GetTPCNcls()
2de8763b 263 ,esdTrackHLT->Theta()
ff84363f 264 ,esdTrackHLT->Eta()
2de8763b 265 ,esdTrackHLT->Phi()
ff84363f 266 ,dca[0]
267 ,dca[1]
9d716a57 268 ,esdTrackHLT->Charge()
269 ,DCAr
270 ,DCAz
ff84363f 271 ,esdTrackHLT->GetNcls(0)
272 ,nr_tracksHLT
273 ,vertHLT->GetStatus()
274 };
275 fTrackHLT->Fill(trackHLT);
276 }
277 // Post output data.
278 PostData(1, fOutputList);
279}
280
281void AliAnalysisTaskHLTCentralBarrel::Terminate(Option_t *){
282 // see header file of AliAnalysisTask for documentation
283}
284
285Int_t AliAnalysisTaskHLTCentralBarrel::CalculateCentrality(AliESDEvent* esd){
286//see header for documentation
287
288 Int_t centrality =-1;
289
290 AliESDVZERO* esdV0 = esd->GetVZEROData();
291 //AliESDZDC* esdZDC = esd->GetZDCData();
292 //Int_t partZDC = esdZDC->GetZDCParticipants();
293
294 Float_t multV0 = esdV0->GetMTotV0A() + esdV0->GetMTotV0C();
295
296 if ( multV0 >= 0. && multV0 <= 124.5 ) centrality = 90;
297 else if ( multV0 > 124.5 && multV0 <= 274.5 ) centrality = 80;
298 else if ( multV0 > 274.5 && multV0 <= 574.5 ) centrality = 70;
299 else if ( multV0 > 574.5 && multV0 <= 1224.5 ) centrality = 60;
300 else if ( multV0 > 1224.5 && multV0 <= 2174.5 ) centrality = 50;
301 else if ( multV0 > 2174.5 && multV0 <= 3624.5 ) centrality = 40;
302 else if ( multV0 > 3624.5 && multV0 <= 5574.5 ) centrality = 30;
303 else if ( multV0 > 5574.5 && multV0 <= 8274.5 ) centrality = 20;
304 else if ( multV0 > 8274.5 && multV0 <= 12024.5 ) centrality = 10;
305 else if ( multV0 > 12024.5 && multV0 <= 14674.5 ) centrality = 5;
306 else if ( multV0 > 14674.5 && multV0 <= 19449.5 ) centrality = 0;
307
308 return centrality;
309}
310
311THnSparseF* AliAnalysisTaskHLTCentralBarrel::CreateEventTHnSparse(const char* name, Int_t size, const Int_t* bins, Double_t* min, Double_t* max){
312//see header for documentation
313
314 THnSparseF *thn = new THnSparseF(name,"",size,bins,min,max);
2de8763b 315 thn->GetAxis(0)->SetTitle("primary vertex x");
316 thn->GetAxis(1)->SetTitle("primary vertex y");
317 thn->GetAxis(2)->SetTitle("primary vertex z");
ff84363f 318 thn->GetAxis(3)->SetTitle("number of contributors");
319 thn->GetAxis(4)->SetTitle("track multiplicity");
320 thn->GetAxis(5)->SetTitle("vertex status");
321 return thn;
322}
323
324THnSparseF* AliAnalysisTaskHLTCentralBarrel::CreateTrackTHnSparse(const char* name, Int_t size, const Int_t* bins, Double_t* min, Double_t* max){
325//see header for documentation
326
327 THnSparseF *thn = new THnSparseF(name,"",size,bins,min,max);
2de8763b 328 thn->GetAxis(0)->SetTitle("p_{T}");
329 thn->GetAxis(1)->SetTitle("TPC clusters/track");
330 thn->GetAxis(2)->SetTitle("#theta");
331 thn->GetAxis(3)->SetTitle("#eta");
332 thn->GetAxis(4)->SetTitle("#phi");
ff84363f 333 thn->GetAxis(5)->SetTitle("DCAr");
334 thn->GetAxis(6)->SetTitle("DCAz");
2de8763b 335 thn->GetAxis(7)->SetTitle("polarity");
ff84363f 336 thn->GetAxis(8)->SetTitle("DCArSG");
337 thn->GetAxis(9)->SetTitle("DCAzSG");
2de8763b 338 thn->GetAxis(10)->SetTitle("ITS clusters/track");
ff84363f 339 return thn;
340}
341
342// void AliAnalysisTaskHLTCentralBarrel::Fill(AliESDevent *esd, THnSparseF* thn){
343// //see header for documentation
344//
345// int nTracks = 0;
346//
347// for(int i=0; i<esdHLT->GetNumberOfTracks(); i++){
348//
349// AliESDtrack *esdTrack = esd->GetTrack(i);
350// if(!esdTrack) continue;
351// if(!(esdTrack->GetStatus()&AliESDtrack::kTPCin)) continue;
352//
353// //DCA calculations
354// Float_t DCAr=-99;
355// Float_t DCAz=-99;
356// Float_t dca[2];
357// if(vertHLT->GetX()==0 && vertHLT->GetY()==0 && vertHLT->GetZ() ==0 ){
358// dca[0]=-99;
359// dca[1]=-99;
360// }
361// else{
362// //Calculating DCA "old" fashion
363// esdTrackHLT->GetDZ(esdHLT->GetPrimaryVertex()->GetXv(), esdHLT->GetPrimaryVertex()->GetYv(), esdHLT->GetPrimaryVertex()->GetZv(), bfield, dca);
364// // plotting the DCA calculated by Sergey
365// esdTrackHLT->GetImpactParametersTPC(DCAr,DCAz);
366// }
367//
368// Double_t trackHLT[] = {
369// TMath::Abs(esdTrackHLT->Pt())
370// ,esdTrackHLT->GetTPCNcls()
371// ,esdTrackHLT->Theta()
372// ,esdTrackHLT->Eta()
373// ,esdTrackHLT->Phi()
374// ,DCAr
375// ,DCAz
376// ,esdTrackHLT->Charge()
377// ,dca[0]
378// ,dca[1]
379// ,esdTrackHLT->GetStatus()
380// ,esdTrackHLT->GetNcls(0)
381// };
382// fTrackHLT->Fill(trackHLT);
383// nr_tracksHLT++;
384// }
385// Double_t eventHLT[] = { vertHLT->GetX(), vertHLT->GetY(), vertHLT->GetZ(), vertHLT->GetNContributors(), nr_tracksHLT, vertHLT->GetStatus()};
386// fEventHLT->Fill(eventHLT);
387// }