]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/vertexingHF/AliAnalysisTaskSEHFQA.cxx
o updates for AOD MC (Julian)
[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliAnalysisTaskSEHFQA.cxx
CommitLineData
0dbb51a0 1/**************************************************************************
2 * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
27de2dfb 16/* $Id$ */
17
0dbb51a0 18/////////////////////////////////////////////////////////////
19//
20// AliAnalysisTaskSE for HF quality assurance
21//
22// Author: Chiara Bianchin, chiara.bianchin@pd.infn.it
23/////////////////////////////////////////////////////////////
24
25#include <Riostream.h>
26#include <TClonesArray.h>
27#include <TCanvas.h>
28#include <TList.h>
29#include <TH1F.h>
30#include <TH2F.h>
d48d81b5 31#include <TH3F.h>
32#include <TProfile2D.h>
0dbb51a0 33#include <TDatabasePDG.h>
34
35#include <AliAnalysisDataSlot.h>
36#include <AliAnalysisDataContainer.h>
37#include "AliAnalysisManager.h"
38#include "AliESDtrack.h"
824b448a 39#include "AliESDVertex.h"
0dbb51a0 40#include "AliVertexerTracks.h"
41#include "AliPID.h"
42#include "AliTPCPIDResponse.h"
43#include "AliAODHandler.h"
44#include "AliAODEvent.h"
45#include "AliAODVertex.h"
46#include "AliAODTrack.h"
47#include "AliAODMCParticle.h"
a4ef4383 48#include "AliAODMCHeader.h"
0dbb51a0 49#include "AliAODRecoDecayHF2Prong.h"
50#include "AliAODRecoCascadeHF.h"
51#include "AliAnalysisVertexingHF.h"
52#include "AliAnalysisTaskSE.h"
a4ef4383 53#include "AliCounterCollection.h"
0dbb51a0 54#include "AliRDHFCuts.h"
55#include "AliRDHFCutsDplustoKpipi.h"
56#include "AliRDHFCutsD0toKpipipi.h"
57#include "AliRDHFCutsDstoKKpi.h"
58#include "AliRDHFCutsDStartoKpipi.h"
59#include "AliRDHFCutsD0toKpi.h"
60#include "AliRDHFCutsLctopKpi.h"
a4b0ffc7 61#include "AliInputEventHandler.h"
0dbb51a0 62
d48d81b5 63#include "AliFlowEvent.h"
64#include "AliFlowTrackCuts.h"
65#include "AliFlowTrackSimple.h"
66#include "AliFlowVector.h"
67
0dbb51a0 68#include "AliAnalysisTaskSEHFQA.h"
69
0dbb51a0 70ClassImp(AliAnalysisTaskSEHFQA)
71
72//____________________________________________________________________________
73
74AliAnalysisTaskSEHFQA::AliAnalysisTaskSEHFQA():AliAnalysisTaskSE(),
75 fNEntries(0x0),
76 fOutputPID(0x0),
77 fOutputTrack(0x0),
a4ef4383 78 fOutputCounters(0x0),
79 fOutputCheckCentrality(0x0),
818c1271 80 fOutputEvSelection(0x0),
d48d81b5 81 fOutputFlowObs(0x0),
0dbb51a0 82 fDecayChannel(AliAnalysisTaskSEHFQA::kD0toKpi),
a4ef4383 83 fCuts(0x0),
d48d81b5 84 fFlowEvent(0x0),
85 fRFPcuts(0x0),
a4ef4383 86 fEstimator(AliRDHFCuts::kCentTRK),
19222b31 87 fReadMC(kFALSE),
5fc4893f 88 fSimpleMode(kFALSE),
89 fOnOff()
0dbb51a0 90{
91 //default constructor
5fc4893f 92 fOnOff[0]=kTRUE;
93 fOnOff[1]=kTRUE;
94 fOnOff[2]=kTRUE;
818c1271 95 fOnOff[3]=kTRUE;
d48d81b5 96 fOnOff[4]=kTRUE;
0dbb51a0 97}
98
99//____________________________________________________________________________
100AliAnalysisTaskSEHFQA::AliAnalysisTaskSEHFQA(const char *name, AliAnalysisTaskSEHFQA::DecChannel ch,AliRDHFCuts* cuts):
a4ef4383 101 AliAnalysisTaskSE(name),
102 fNEntries(0x0),
103 fOutputPID(0x0),
104 fOutputTrack(0x0),
105 fOutputCounters(0x0),
106 fOutputCheckCentrality(0x0),
818c1271 107 fOutputEvSelection(0x0),
d48d81b5 108 fOutputFlowObs(0x0),
a4ef4383 109 fDecayChannel(ch),
110 fCuts(0x0),
d48d81b5 111 fFlowEvent(0x0),
112 fRFPcuts(0x0),
a4ef4383 113 fEstimator(AliRDHFCuts::kCentTRK),
19222b31 114 fReadMC(kFALSE),
5fc4893f 115 fSimpleMode(kFALSE),
116 fOnOff()
0dbb51a0 117{
118 //constructor
119
120 //SetCutObject(cuts);
121 fCuts=cuts;
122
5fc4893f 123 fOnOff[0]=kTRUE;
124 fOnOff[1]=kTRUE;
125 fOnOff[2]=kTRUE;
818c1271 126 fOnOff[3]=kTRUE;
d48d81b5 127 fOnOff[4]=kTRUE;
5fc4893f 128
0dbb51a0 129 // Output slot #1 writes into a TH1F container (number of events)
130 DefineOutput(1,TH1F::Class()); //My private output
131 // Output slot #2 writes into a TList container (PID)
5fc4893f 132 if (fOnOff[1]) DefineOutput(2,TList::Class()); //My private output
0dbb51a0 133 // Output slot #3 writes into a TList container (Tracks)
5fc4893f 134 if (fOnOff[0]) DefineOutput(3,TList::Class()); //My private output
0dbb51a0 135 // Output slot #4 writes into a AliRDHFCuts container (cuts)
136 switch(fDecayChannel){
137 case 0:
138 DefineOutput(4,AliRDHFCutsDplustoKpipi::Class()); //My private output
139 break;
140 case 1:
141 DefineOutput(4,AliRDHFCutsD0toKpi::Class()); //My private output
142 break;
143 case 2:
144 DefineOutput(4,AliRDHFCutsDStartoKpipi::Class()); //My private output
145 break;
146 case 3:
147 DefineOutput(4,AliRDHFCutsDstoKKpi::Class()); //My private output
148 break;
149 case 4:
150 DefineOutput(4,AliRDHFCutsD0toKpipipi::Class()); //My private output
151 break;
a4ef4383 152 case 5:
0dbb51a0 153 DefineOutput(4,AliRDHFCutsLctopKpi::Class()); //My private output
154 break;
a4ef4383 155 }
5fc4893f 156 if (fOnOff[2]) {
157 // Output slot #5 writes into a TList container (AliCounterCollection)
158 DefineOutput(5,TList::Class()); //My private output
159 // Output slot #6 writes into a TList container (TH1F)
160 DefineOutput(6,TList::Class()); //My private output
161 }
818c1271 162
163 if(fOnOff[3]) DefineOutput(7,TList::Class()); //My private output
d48d81b5 164 if(fOnOff[4]) DefineOutput(8,TList::Class()); //My private output
818c1271 165
0dbb51a0 166}
167
168//___________________________________________________________________________
169AliAnalysisTaskSEHFQA::~AliAnalysisTaskSEHFQA()
170{
171 //destructor
a4ef4383 172 delete fNEntries;
173
174 delete fOutputPID;
175
176 delete fOutputTrack;
177
178 delete fOutputCounters;
179
180 delete fOutputCheckCentrality;
181
818c1271 182 delete fOutputEvSelection;
183
d48d81b5 184 if(fOnOff[4]) {
185 delete fOutputFlowObs;
186 delete fFlowEvent;
187 }
0dbb51a0 188}
189
190//___________________________________________________________________________
191void AliAnalysisTaskSEHFQA::Init(){
192
193 //initialization
194 if(fDebug > 1) printf("AnalysisTaskSEHFQA::Init() \n");
195
196 switch(fDecayChannel){
197 case 0:
198 {
199 AliRDHFCutsDplustoKpipi* copycut=new AliRDHFCutsDplustoKpipi(*(static_cast<AliRDHFCutsDplustoKpipi*>(fCuts)));
200 // Post the data
201 PostData(4,copycut);
202 }
203 break;
204 case 1:
205 {
206 AliRDHFCutsD0toKpi* copycut=new AliRDHFCutsD0toKpi(*(static_cast<AliRDHFCutsD0toKpi*>(fCuts)));
207 // Post the data
208 PostData(4,copycut);
209 }
210 break;
211 case 2:
212 {
213 AliRDHFCutsDStartoKpipi* copycut=new AliRDHFCutsDStartoKpipi(*(static_cast<AliRDHFCutsDStartoKpipi*>(fCuts)));
214 // Post the data
215 PostData(4,copycut);
216 }
217 break;
218 case 3:
219 {
220 AliRDHFCutsDstoKKpi* copycut=new AliRDHFCutsDstoKKpi(*(static_cast<AliRDHFCutsDstoKKpi*>(fCuts)));
221 // Post the data
222 PostData(4,copycut);
223 }
224 break;
225 case 4:
226 {
227 AliRDHFCutsD0toKpipipi* copycut=new AliRDHFCutsD0toKpipipi(*(static_cast<AliRDHFCutsD0toKpipipi*>(fCuts)));
228 // Post the data
229 PostData(4,copycut);
230 }
231 break;
232 case 5:
233 {
234 AliRDHFCutsLctopKpi* copycut=new AliRDHFCutsLctopKpi(*(static_cast<AliRDHFCutsLctopKpi*>(fCuts)));
235 // Post the data
236 PostData(4,copycut);
237 }
238 break;
239
240 default:
241 return;
242 }
243
244
245
246}
247
248//___________________________________________________________________________
249void AliAnalysisTaskSEHFQA::UserCreateOutputObjects()
250{
251
252 //create the output container
253 if(fDebug > 1) printf("AnalysisTaskSEHFQA::UserCreateOutputObjects() \n");
254
255 //count events
256
acc33385 257 fNEntries=new TH1F(GetOutputSlot(1)->GetContainer()->GetName(), "Counts the number of events", 10,-0.5,9.5);
0dbb51a0 258 fNEntries->GetXaxis()->SetBinLabel(1,"nEventsAnal");
259 fNEntries->GetXaxis()->SetBinLabel(2,"Pile-up Rej");
260 fNEntries->GetXaxis()->SetBinLabel(3,"No VertexingHF");
9af24f46 261 fNEntries->GetXaxis()->SetBinLabel(4,"nCandidates(AnCuts)");
262 fNEntries->GetXaxis()->SetBinLabel(5,"EventsWithGoodVtx");
5fc4893f 263 //fNEntries->GetXaxis()->SetBinLabel(6,"N. of 0SMH");
264 fNEntries->GetXaxis()->SetBinLabel(6,"N. of CSH1");
a4ef4383 265 if(fReadMC){
266 fNEntries->GetXaxis()->SetBinLabel(7,"MC Cand from c");
267 fNEntries->GetXaxis()->SetBinLabel(8,"MC Cand from b");
acc33385 268 fNEntries->GetXaxis()->SetBinLabel(9,"N fake Trks");
269 fNEntries->GetXaxis()->SetBinLabel(10,"N true Trks");
a4ef4383 270 } else{
271 fNEntries->GetXaxis()->SetBinLabel(7,"N candidates");
272 }
273
0dbb51a0 274 fNEntries->GetXaxis()->SetNdivisions(1,kFALSE);
275
276 //PID
5fc4893f 277 if(fOnOff[1]){
278 fOutputPID=new TList();
279 fOutputPID->SetOwner();
280 fOutputPID->SetName(GetOutputSlot(2)->GetContainer()->GetName());
0dbb51a0 281
5fc4893f 282 //TOF pid
23dfe729 283 TH1F* hTOFflags=new TH1F("hTOFflags","TOF flags",6,-0.5,5.5);
b906dda8 284 hTOFflags->SetMinimum(0.);
285 hTOFflags->GetXaxis()->SetBinLabel(1,"All Tracks");
23dfe729 286 hTOFflags->GetXaxis()->SetBinLabel(2,"kTPCout");
287 hTOFflags->GetXaxis()->SetBinLabel(3,"kTOFout");
288 hTOFflags->GetXaxis()->SetBinLabel(4,"kTIME");
289 hTOFflags->GetXaxis()->SetBinLabel(5,"kTOFpid");
290 hTOFflags->GetXaxis()->SetBinLabel(6,"kTOFmismatch");
b906dda8 291
5fc4893f 292 TString hname="hTOFsig";
293 TH1F* hTOFsig=new TH1F(hname.Data(),"Distribution of TOF signal;TOF time [ps];Entries", 100, -2.e3,40.e3);
0dbb51a0 294
5fc4893f 295 hname="hTOFtime";
296 TH1F* hTOFtime=new TH1F(hname.Data(),"Distribution of TOF time Kaon;TOF time(Kaon) [ps];Entries", 1000, 0.,50000.);
0dbb51a0 297
5fc4893f 298 hname="hTOFtimeKaonHyptime";
b906dda8 299 TH2F* hTOFtimeKaonHyptime=new TH2F(hname.Data(),"TOFtime - timeHypothesisForKaon;p[GeV/c];TOFtime - timeHypothesisForKaon [ps]",500,0.,10.,1000,-20000.,20000.);
0dbb51a0 300
5fc4893f 301 hname="hTOFtimeKaonHyptimeAC";
b906dda8 302 TH2F* hTOFtimeKaonHyptimeAC=new TH2F(hname.Data(),"TOFtime - timeHypothesisForKaon;p[GeV/c];TOFtime - timeHypothesisForKaon [ps]",500,0.,10.,1000,-20000.,20000.);
0dbb51a0 303
5fc4893f 304 hname="hTOFsigmaKSigPid";
23dfe729 305 TH2F* hTOFsigmaKSigPid=new TH2F(hname.Data(),"(TOFsignal-timeK)/tofSigPid;p[GeV/c];(TOFsignal-timeK)/tofSigPid",500,0.,10.,400,-20,20);
0dbb51a0 306
5fc4893f 307 hname="hTOFsigmaPionSigPid";
23dfe729 308 TH2F* hTOFsigmaPionSigPid=new TH2F(hname.Data(),"(TOFsignal-time#pi)/tofSigPid;p[GeV/c];(TOFsignal-time#pi)/tofSigPid",500,0.,10.,400,-20,20);
a4ef4383 309
5fc4893f 310 hname="hTOFsigmaProtonSigPid";
23dfe729 311 TH2F* hTOFsigmaProtonSigPid=new TH2F(hname.Data(),"(TOFsignal-timep)/tofSigPid;p[GeV/c];(TOFsignal-time p)/tofSigPid",500,0.,10.,400,-20,20);
a4ef4383 312
5fc4893f 313 hname="hTOFsigPid3sigPion";
314 TH1F* hTOFsigPid3sigPion=new TH1F(hname.Data(),"TOF PID resolution (#pi) [ps]",500,0.,1000.);
a4ef4383 315
5fc4893f 316 hname="hTOFsigPid3sigKaon";
317 TH1F* hTOFsigPid3sigKaon=new TH1F(hname.Data(),"TOF PID resolution (K) [ps]",500,0.,1000.);
0dbb51a0 318
5fc4893f 319 hname="hTOFsigPid3sigProton";
320 TH1F* hTOFsigPid3sigProton=new TH1F(hname.Data(),"TOF PID resolution (p) [ps]",500,0.,1000.);
a4ef4383 321
0dbb51a0 322
5fc4893f 323 //TPC pid
324 hname="hTPCsig";
325 TH1F* hTPCsig=new TH1F(hname.Data(),"Distribution of TPC signal;TPC sig;Entries", 100, 35.,100.);
0dbb51a0 326
5fc4893f 327 hname="hTPCsigvsp";
b906dda8 328 TH2F* hTPCsigvsp=new TH2F(hname.Data(),"TPCsig vs p;TPC p[GeV/c];TPCsig",500,0.,10.,1000,35.,100.);
0dbb51a0 329
5fc4893f 330 hname="hTPCsigvspAC";
b906dda8 331 TH2F* hTPCsigvspAC=new TH2F(hname.Data(),"TPCsig vs p;TPCp[GeV/c];TPCsig",500,0.,10.,1000,35.,100.);
5fc4893f 332
333 hname="hTPCsigmaK";
b906dda8 334 TH2F* hTPCsigmaK=new TH2F(hname.Data(),"TPC Sigma for K as a function of momentum;p[GeV/c];Sigma Kaon",500,0.,10.,400,-20,20);
5fc4893f 335
336 hname="hTPCsigmaPion";
b906dda8 337 TH2F* hTPCsigmaPion=new TH2F(hname.Data(),"TPC Sigma for #pi as a function of momentum;p[GeV/c];Sigma #pi",500,0.,10.,400,-20,20);
5fc4893f 338
339 hname="hTPCsigmaProton";
b906dda8 340 TH2F* hTPCsigmaProton=new TH2F(hname.Data(),"TPC Sigma for proton as a function of momentum;p[GeV/c];Sigma Proton",500,0.,10.,400,-20,20);
5fc4893f 341
b906dda8 342 fOutputPID->Add(hTOFflags);
5fc4893f 343 fOutputPID->Add(hTOFsig);
344 fOutputPID->Add(hTPCsig);
345 fOutputPID->Add(hTOFtime);
346 fOutputPID->Add(hTOFtimeKaonHyptime);
347 fOutputPID->Add(hTOFtimeKaonHyptimeAC);
5fc4893f 348 fOutputPID->Add(hTOFsigmaKSigPid);
349 fOutputPID->Add(hTOFsigmaPionSigPid);
350 fOutputPID->Add(hTOFsigmaProtonSigPid);
351 fOutputPID->Add(hTOFsigPid3sigPion);
352 fOutputPID->Add(hTOFsigPid3sigKaon);
353 fOutputPID->Add(hTOFsigPid3sigProton);
354 fOutputPID->Add(hTPCsigvsp);
355 fOutputPID->Add(hTPCsigvspAC);
356 fOutputPID->Add(hTPCsigmaK);
357 fOutputPID->Add(hTPCsigmaPion);
358 fOutputPID->Add(hTPCsigmaProton);
24748c57 359
360 if(fReadMC){
361 //TOF
362 hname="hTOFsigmaMCKSigPid";
363 TH2F* hTOFsigmaMCKSigPid=new TH2F(hname.Data(),"(TOFsignal-timeK)/tofSigPid;p[GeV/c];(TOFsignal-timeK)/tofSigPid",500,0.,10.,400,-20,20);
364
365 hname="hTOFsigmaMCPionSigPid";
366 TH2F* hTOFsigmaMCPionSigPid=new TH2F(hname.Data(),"(TOFsignal-time#pi)/tofSigPid;p[GeV/c];(TOFsignal-time#pi)/tofSigPid",500,0.,10.,400,-20,20);
367
368 hname="hTOFsigmaMCProtonSigPid";
369 TH2F* hTOFsigmaMCProtonSigPid=new TH2F(hname.Data(),"(TOFsignal-timep)/tofSigPid;p[GeV/c];(TOFsignal-time p)/tofSigPid",500,0.,10.,400,-20,20);
370
371 //TPC
372 hname="hTPCsigmaMCK";
373 TH2F* hTPCsigmaMCK=new TH2F(hname.Data(),"TPC Sigma for K as a function of momentum;p[GeV/c];Sigma Kaon",500,0.,10.,400,-20,20);
374
375 hname="hTPCsigmaMCPion";
376 TH2F* hTPCsigmaMCPion=new TH2F(hname.Data(),"TPC Sigma for #pi as a function of momentum;p[GeV/c];Sigma #pi",500,0.,10.,400,-20,20);
377
378 hname="hTPCsigmaMCProton";
379 TH2F* hTPCsigmaMCProton=new TH2F(hname.Data(),"TPC Sigma for proton as a function of momentum;p[GeV/c];Sigma Proton",500,0.,10.,400,-20,20);
380
381 fOutputPID->Add(hTOFsigmaMCKSigPid);
382 fOutputPID->Add(hTOFsigmaMCPionSigPid);
383 fOutputPID->Add(hTOFsigmaMCProtonSigPid);
384 fOutputPID->Add(hTPCsigmaMCK);
385 fOutputPID->Add(hTPCsigmaMCPion);
386 fOutputPID->Add(hTPCsigmaMCProton);
387
388 }
5fc4893f 389 }
0dbb51a0 390
391 //quality of the tracks
5fc4893f 392 if(fOnOff[0]){
393 fOutputTrack=new TList();
394 fOutputTrack->SetOwner();
395 fOutputTrack->SetName(GetOutputSlot(3)->GetContainer()->GetName());
0dbb51a0 396
5fc4893f 397 TString hname="hnClsITS";
398 TH1F* hnClsITS=new TH1F(hname.Data(),"Distribution of number of ITS clusters;nITScls;Entries",7,-0.5,6.5);
0dbb51a0 399
39dd297f 400 hname="hnClsITSselTr";
401 TH1F* hnClsITSselTr=new TH1F(hname.Data(),"Distribution of number of ITS clusters selected tracks;nITScls;Entries",7,-0.5,6.5);
402
5fc4893f 403 hname="hnClsITS-SA";
404 TH1F* hnClsITSSA=new TH1F(hname.Data(),"Distribution of number of ITS clusters(ITS-SA);nITScls;Entries",7,-0.5,6.5);
0dbb51a0 405
bc6846db 406
407 hname="hnLayerITS";
f3e1ad2f 408 TH1F* hnLayerITS=new TH1F(hname.Data(),"Number of tracks with point in layer;ITS layer;",7,-1.5,5.5);
409 hnLayerITS->GetXaxis()->SetBinLabel(1,"n tracks");
bc6846db 410
411 hname="hnLayerITSsa";
f3e1ad2f 412 TH1F* hnLayerITSsa=new TH1F(hname.Data(),"Number of tracks with point in layer;ITS layer;",7,-1.5,5.5);
413 hnLayerITSsa->GetXaxis()->SetBinLabel(1,"n tracks");
414
5fc4893f 415 hname="hnClsSPD";
416 TH1F* hnClsSPD=new TH1F(hname.Data(),"Distribution of number of SPD clusters;nSPDcls;Entries",3,-0.5,2.5);
0dbb51a0 417
5fc4893f 418 hname="hptGoodTr";
419 TH1F* hptGoodTr=new TH1F(hname.Data(),"Pt distribution of 'good' tracks;p_{t}[GeV];Entries/0.05 GeV/c",400,0.,20.);
420 hptGoodTr->SetTitleOffset(1.3,"Y");
0dbb51a0 421
824b448a 422 if(!fSimpleMode){
423 hname="hptGoodTrFromDaugh";
424 TH1F* hptGoodTrFromDaugh=new TH1F(hname.Data(),"Pt distribution of 'good' candidate's daughters;p_{t}[GeV];Entries/0.05 GeV/c",400,0.,20.);
425 hptGoodTrFromDaugh->SetTitleOffset(1.3,"Y");
426 fOutputTrack->Add(hptGoodTrFromDaugh);
427 }
428
5fc4893f 429 hname="hdistrGoodTr";
824b448a 430 TH1F* hdistrGoodTr=new TH1F(hname.Data(),"Distribution of number of 'good' candidate's daughters per event;no.good-tracks/ev;Entries",4000,-0.5,3999.5);
5fc4893f 431 hdistrGoodTr->SetTitleOffset(1.3,"Y");
0dbb51a0 432
824b448a 433 hname="hdistrSelTr";
434 TH1F* hdistrSelTr=new TH1F(hname.Data(),"Distribution of number of Selected tracks per event;no.good-tracks/ev;Entries",4000,-0.5,3999.5);
435 hdistrSelTr->SetTitleOffset(1.3,"Y");
436
5fc4893f 437 hname="hd0";
4563d678 438 TH1F* hd0=new TH1F(hname.Data(),"Impact parameter (rphi) distribution of 'good' tracks;d_{0rphi}[cm];Entries/10^{3} cm",200,-0.1,0.1);
439
440 hname="hd0z";
441 TH1F* hd0z=new TH1F(hname.Data(),"Impact parameter (z) distribution of 'good' tracks;d_{0z}[cm];Entries/10^{3} cm",200,-0.1,0.1);
0dbb51a0 442
5fc4893f 443 fOutputTrack->Add(hnClsITS);
39dd297f 444 fOutputTrack->Add(hnClsITSselTr);
5fc4893f 445 fOutputTrack->Add(hnClsITSSA);
bc6846db 446 fOutputTrack->Add(hnLayerITS);
447 fOutputTrack->Add(hnLayerITSsa);
5fc4893f 448 fOutputTrack->Add(hnClsSPD);
449 fOutputTrack->Add(hptGoodTr);
450 fOutputTrack->Add(hdistrGoodTr);
824b448a 451 fOutputTrack->Add(hdistrSelTr);
5fc4893f 452 fOutputTrack->Add(hd0);
4563d678 453 fOutputTrack->Add(hd0z);
acc33385 454
5fc4893f 455 if(fReadMC){
456 hname="hdistrFakeTr";
457 TH1F* hdistrFakeTr=new TH1F(hname.Data(),"Distribution of number of fake tracks per event;no.fake-tracks/ev;Entries",4000,-0.5,3999.5);
824b448a 458 hdistrFakeTr->SetTitleOffset(1.3,"Y");
acc33385 459
5fc4893f 460 hname="hd0f";
461 TH1F* hd0f=new TH1F(hname.Data(),"Impact parameter distribution of fake tracks;d_{0}[cm];Entries/10^{3} cm",200,-0.1,0.1);
acc33385 462
5fc4893f 463 hname="hptFakeTr";
464 TH1F* hptFakeTr=new TH1F(hname.Data(),"Pt distribution of fake tracks;p_{t}[GeV];Entries/0.05 GeV/c",400,0.,20.);
ad42e35b 465 hptFakeTr->SetTitleOffset(1.3,"Y");
824b448a 466 if(!fSimpleMode){
467 hname="hptFakeTrFromDaugh";
468 TH1F* hptFakeTrFromDaugh=new TH1F(hname.Data(),"Pt distribution of fake tracks from daughters;p_{t}[GeV];Entries/0.05 GeV/c",400,0.,20.);
469 hptFakeTrFromDaugh->SetTitleOffset(1.3,"Y");
470 fOutputTrack->Add(hptFakeTrFromDaugh);
471 }
acc33385 472
5fc4893f 473 fOutputTrack->Add(hptFakeTr);
474 fOutputTrack->Add(hdistrFakeTr);
475 fOutputTrack->Add(hd0f);
476
477 }
acc33385 478 }
479
0dbb51a0 480
5fc4893f 481 if(fOnOff[2] && fCuts->GetUseCentrality()){
a4ef4383 482
483 //Centrality (Counters)
484 fOutputCounters=new TList();
485 fOutputCounters->SetOwner();
486 fOutputCounters->SetName(GetOutputSlot(5)->GetContainer()->GetName());
487
488 AliCounterCollection *stdEstimator=new AliCounterCollection("stdEstimator");
489 stdEstimator->AddRubric("run",500000);
5fc4893f 490 stdEstimator->AddRubric("centralityclass","-10_0/0_10/10_20/20_30/30_40/40_50/50_60/60_70/70_80/80_90/90_100/-990_-980");
a4ef4383 491 stdEstimator->Init();
492 AliCounterCollection *secondEstimator=new AliCounterCollection("secondEstimator");
493 secondEstimator->AddRubric("run",500000);
5fc4893f 494 secondEstimator->AddRubric("centralityclass","-10_0/0_10/10_20/20_30/30_40/40_50/50_60/60_70/70_80/80_90/90_100/-990_-980");
a4ef4383 495 secondEstimator->Init();
1b765719 496
a4ef4383 497 fOutputCounters->Add(stdEstimator);
498 fOutputCounters->Add(secondEstimator);
1b765719 499
a4ef4383 500 //Centrality (Checks)
501 fOutputCheckCentrality=new TList();
502 fOutputCheckCentrality->SetOwner();
503 fOutputCheckCentrality->SetName(GetOutputSlot(6)->GetContainer()->GetName());
504
5fc4893f 505 TString hname="hNtrackletsIn";
a4ef4383 506 TH1F* hNtrackletsIn=new TH1F(hname.Data(),"Number of tracklets in Centrality range;ntracklets;Entries",5000,-0.5,4999.5);
507
508 hname="hMultIn";
509 TH1F* hMultIn=new TH1F(hname.Data(),"Multiplicity;multiplicity in Centrality range;Entries",10000,-0.5,9999.5);
510
511 hname="hNtrackletsOut";
512 TH1F* hNtrackletsOut=new TH1F(hname.Data(),"Number of tracklets out of Centrality range;ntracklets;Entries",5000,-0.5,4999.5);
513
514 hname="hMultOut";
515 TH1F* hMultOut=new TH1F(hname.Data(),"Multiplicity out of Centrality range;multiplicity;Entries",10000,-0.5,9999.5);
516
ac0c2841 517 hname="hMultvsPercentile";
7a5f2502 518 TH2F* hMultvsPercentile=new TH2F(hname.Data(),"Multiplicity vs Percentile;multiplicity;percentile",10000,-0.5,9999.5,240,-10.,110);
ac0c2841 519
7a5f2502 520 hname="hntrklvsPercentile";
521 TH2F* hntrklvsPercentile=new TH2F(hname.Data(),"N tracklets vs Percentile;ntracklets;percentile",5000,-0.5,4999.5,240,-10.,110);
23dfe729 522
138e0d89 523 hname="hnTPCTracksvsPercentile";
524 TH2F* hnTPCTracksvsPercentile=new TH2F(hname.Data(),"N TPC tracks vs Percentile;nTPCTracks;percentile",5000,-0.5,9999.5,240,-10.,110);
525
526 hname="hnTPCITSTracksvsPercentile";
527 TH2F* hnTPCITSTracksvsPercentile=new TH2F(hname.Data(),"N TPC+ITS tracks vs Percentile;nTPCITSTracks;percentile",5000,-0.5,9999.5,240,-10.,110);
528
529 hname="hnTPCITS1SPDTracksvsPercentile";
530 TH2F* hnTPCITS1SPDTracksvsPercentile=new TH2F(hname.Data(),"N TPC+ITS+1SPD tracks vs Percentile;nTPCITS1SPDTracks;percentile",5000,-0.5,9999.5,240,-10.,110);
531
b68a0f6a 532 hname="hV0MultiplicityPercentile";
533 TH2F*hV0MultiplicityPercentile = new TH2F(hname.Data(),"V0 Multiplicity vs Percentile;V0 multiplicity;percentile",1000,-0.5,9999.5,120,-10.,110);
534
535 hname="hV0MultiplicityNtrackletsIn";
536 TH2F*hV0MultiplicityNtrackletsIn = new TH2F(hname.Data(),"V0 Multiplicity vs Number of tracklets in the CC;V0 multiplicity;percentile",1000,-0.5,9999.5,5000,-0.5,4999.5);
537
538 hname="hStdPercentileSPDPercentile";
539 TH2F* hStdPercentileSPDPercentile = new TH2F(hname.Data(),"Std estimator Percentile Vs SPD Percentile;Std estimator percentile;SPD percentile",120,-10.,110,120,-10.,110);
540
a4ef4383 541 fOutputCheckCentrality->Add(hNtrackletsIn);
542 fOutputCheckCentrality->Add(hNtrackletsOut);
543 fOutputCheckCentrality->Add(hMultIn);
544 fOutputCheckCentrality->Add(hMultOut);
ac0c2841 545 fOutputCheckCentrality->Add(hMultvsPercentile);
7a5f2502 546 fOutputCheckCentrality->Add(hntrklvsPercentile);
138e0d89 547 fOutputCheckCentrality->Add(hnTPCTracksvsPercentile);
548 fOutputCheckCentrality->Add(hnTPCITSTracksvsPercentile);
549 fOutputCheckCentrality->Add(hnTPCITS1SPDTracksvsPercentile);
b68a0f6a 550 fOutputCheckCentrality->Add(hV0MultiplicityPercentile);
551 fOutputCheckCentrality->Add(hV0MultiplicityNtrackletsIn);
552 fOutputCheckCentrality->Add(hStdPercentileSPDPercentile);
82f89d0d 553
554 PostData(6,fOutputCheckCentrality);
555
a4ef4383 556 } else{
5fc4893f 557 if(fOnOff[0]){
558 TString hname="hNtracklets";
559 TH1F* hNtracklets=new TH1F(hname.Data(),"Number of tracklets;ntracklets;Entries",5000,-0.5,4999.5);
560
561 hname="hMult";
562 TH1F* hMult=new TH1F(hname.Data(),"Multiplicity;multiplicity;Entries",10000,-0.5,9999.5);
563 fOutputTrack->Add(hNtracklets);
564 fOutputTrack->Add(hMult);
565 }
a4ef4383 566 }
567
818c1271 568 //event selection (z vertex for the moment)
569 if(fOnOff[3]){
570 fOutputEvSelection=new TList();
571 fOutputEvSelection->SetOwner();
572 fOutputEvSelection->SetName(GetOutputSlot(7)->GetContainer()->GetName());
573 AliCounterCollection *evselection=new AliCounterCollection("evselection");
574 evselection->AddRubric("run",500000);
575 evselection->AddRubric("evnonsel","zvtx");
576 evselection->Init();
577
578 TH1F* hzvtx=new TH1F("hzvtx", "Distribution of z_{VTX};z_{VTX} [cm];Entries",100,-20,20);
579
23dfe729 580
2a844546 581 TH2F* hTrigCent=new TH2F("hTrigCent","Centrality vs. Trigger types",12,-1.5,10.5,12,-10,110);
23dfe729 582 hTrigCent->GetXaxis()->SetBinLabel(1,"All");
2a844546 583 hTrigCent->GetXaxis()->SetBinLabel(2,"kAny");
584 hTrigCent->GetXaxis()->SetBinLabel(3,"kMB");
585 hTrigCent->GetXaxis()->SetBinLabel(4,"kINT7");
586 hTrigCent->GetXaxis()->SetBinLabel(5,"kCINT5");
587 hTrigCent->GetXaxis()->SetBinLabel(6,"kCent");
588 hTrigCent->GetXaxis()->SetBinLabel(7,"kSemiCent");
589 hTrigCent->GetXaxis()->SetBinLabel(8,"kEMC1+7");
590 hTrigCent->GetXaxis()->SetBinLabel(9,"kEMCJET+GAMMA");
591 hTrigCent->GetXaxis()->SetBinLabel(10,"Muons");
592 hTrigCent->GetXaxis()->SetBinLabel(11,"PHOS");
593 hTrigCent->GetXaxis()->SetBinLabel(12,"Others");
594
595 TH2F* hTrigMul=new TH2F("hTrigMul","Multiplicity vs. Trigger types",12,-1.5,10.5,100,0.,10000.);
23dfe729 596 hTrigMul->GetXaxis()->SetBinLabel(1,"All");
2a844546 597 hTrigMul->GetXaxis()->SetBinLabel(2,"kAny");
598 hTrigMul->GetXaxis()->SetBinLabel(3,"kMB");
599 hTrigMul->GetXaxis()->SetBinLabel(4,"kINT7");
600 hTrigMul->GetXaxis()->SetBinLabel(5,"kCINT5");
601 hTrigMul->GetXaxis()->SetBinLabel(6,"kCent");
602 hTrigMul->GetXaxis()->SetBinLabel(7,"kSemiCent");
603 hTrigMul->GetXaxis()->SetBinLabel(8,"kEMC1+7");
604 hTrigMul->GetXaxis()->SetBinLabel(9,"kEMCJET+GAMMA");
605 hTrigMul->GetXaxis()->SetBinLabel(10,"Muons");
606 hTrigMul->GetXaxis()->SetBinLabel(11,"PHOS");
607 hTrigMul->GetXaxis()->SetBinLabel(12,"Others");
608
609 TH2F* hTrigCentSel=new TH2F("hTrigCentSel","Trigger types",12,-1.5,10.5,12,-10,110);
23dfe729 610 hTrigCentSel->GetXaxis()->SetBinLabel(1,"All");
2a844546 611 hTrigCentSel->GetXaxis()->SetBinLabel(2,"kAny");
612 hTrigCentSel->GetXaxis()->SetBinLabel(3,"kMB");
613 hTrigCentSel->GetXaxis()->SetBinLabel(4,"kINT7");
614 hTrigCentSel->GetXaxis()->SetBinLabel(5,"kCINT5");
615 hTrigCentSel->GetXaxis()->SetBinLabel(6,"kCent");
616 hTrigCentSel->GetXaxis()->SetBinLabel(7,"kSemiCent");
617 hTrigCentSel->GetXaxis()->SetBinLabel(8,"kEMC1+7");
618 hTrigCentSel->GetXaxis()->SetBinLabel(9,"kEMCJET+GAMMA");
619 hTrigCentSel->GetXaxis()->SetBinLabel(10,"Muons");
620 hTrigCentSel->GetXaxis()->SetBinLabel(11,"PHOS");
621 hTrigCentSel->GetXaxis()->SetBinLabel(12,"Others");
23dfe729 622
1b765719 623 AliCounterCollection *trigCounter=new AliCounterCollection("trigCounter");
624 trigCounter->AddRubric("run",500000);
5e67f11a 625 trigCounter->AddRubric("triggerType","All/Any/MB/Cent/SemiCent/EMCAL/MUON/NoPhysSelMUON/NoPhysSelEvNot7/NoPhysSelCMUP1/NoPhysSelMB/NoPhysSelCent/NoPhysSelSemiCent");
1b765719 626 trigCounter->Init();
627
818c1271 628 fOutputEvSelection->Add(evselection);
629 fOutputEvSelection->Add(hzvtx);
23dfe729 630 fOutputEvSelection->Add(hTrigCent);
631 fOutputEvSelection->Add(hTrigMul);
632 fOutputEvSelection->Add(hTrigCentSel);
1b765719 633 fOutputEvSelection->Add(trigCounter);
818c1271 634 }
d48d81b5 635 if(fOnOff[4]){ // FLOW OBSERVABLES
636 fOutputFlowObs=new TList();
637 fOutputFlowObs->SetOwner();
638 fOutputFlowObs->SetName(GetOutputSlot(8)->GetContainer()->GetName());
639
640 fFlowEvent = new AliFlowEvent(3000);
641 fRFPcuts = new AliFlowTrackCuts("rfpCuts");
642
65274ad6 643 TH2F *hFEvents = new TH2F("hFlowEvents","FlowEvent Selection",7,0,7,7,-10,60);
d48d81b5 644 hFEvents->GetXaxis()->SetBinLabel(1,"REACHED");
645 hFEvents->GetXaxis()->SetBinLabel(2,"TRIGGERED");
65274ad6 646 hFEvents->GetXaxis()->SetBinLabel(3,"kMB");
647 hFEvents->GetXaxis()->SetBinLabel(4,"kCent");
648 hFEvents->GetXaxis()->SetBinLabel(5,"kSemiC");
649 hFEvents->GetXaxis()->SetBinLabel(6,"Triggered + vtx cut");
650 hFEvents->GetXaxis()->SetBinLabel(7,"UnexpectedBehaviour");
d48d81b5 651 fOutputFlowObs->Add(hFEvents);
652
653 TProfile2D *hQ[3];
654 TH2F *hAngleQ[3];
655 TH3F *hPhiEta[3];
656 TString ref[3] = {"FB1","FB128","VZE"};
657 Int_t etabin[3] = {40,40,20};
658 Int_t etamax[3] = { 1, 1, 5};
659 for(Int_t i=0; i<3; ++i) {
660 hQ[i]= new TProfile2D( Form("h%s_Q",ref[i].Data()),
661 Form("Q_{2} components for %s",ref[i].Data()),
662 4,0,4,12,0,60,"s");
663 hQ[i]->GetXaxis()->SetBinLabel(1,"Qx^{-}");
664 hQ[i]->GetXaxis()->SetBinLabel(2,"Qy^{-}");
665 hQ[i]->GetXaxis()->SetBinLabel(3,"Qx^{+}");
666 hQ[i]->GetXaxis()->SetBinLabel(4,"Qy^{+}");
667 hQ[i]->GetYaxis()->SetTitle("Centrality");
668 fOutputFlowObs->Add(hQ[i]);
669
670 hAngleQ[i] = new TH2F( Form("h%s_AngleQ",ref[i].Data()),
671 Form("#Psi_{2} for %s",ref[i].Data()),
672 72,0,TMath::Pi(),12,0,60);
673 hAngleQ[i]->GetXaxis()->SetTitle( Form("#Psi_{2}^{%s}",ref[i].Data()) );
674 hAngleQ[i]->GetYaxis()->SetTitle("Centrality");
675 fOutputFlowObs->Add(hAngleQ[i]);
676
677 hPhiEta[i] = new TH3F( Form("h%s_PhiEta",ref[i].Data()),
678 Form("Eta vs Phi for %s",ref[i].Data()),
65274ad6 679 144,0,TMath::TwoPi(),etabin[i],-1.0*etamax[i],+1.0*etamax[i],12,0,60);
d48d81b5 680 hPhiEta[i]->GetXaxis()->SetTitle("Phi");
681 hPhiEta[i]->GetYaxis()->SetTitle("Eta");
682 hPhiEta[i]->GetZaxis()->SetTitle("Centrality");
683 fOutputFlowObs->Add(hPhiEta[i]);
684
685 }
686 TH3F *hTPCVZE_AngleQ = new TH3F("hTPCVZE_AngleQ","#Psi_{2}^{VZERO} vs #Psi_{2}^{TPC}", 72,0,TMath::Pi(),72,0,TMath::Pi(),12,0,60);
687 hTPCVZE_AngleQ->GetXaxis()->SetTitle("#Psi_{2}^{TPC}");
688 hTPCVZE_AngleQ->GetYaxis()->SetTitle("#Psi_{2}^{VZE}");
689 hTPCVZE_AngleQ->GetZaxis()->SetTitle("Centrality");
690 fOutputFlowObs->Add(hTPCVZE_AngleQ);
65274ad6 691
692 TH2F *hCentVsMultRPS = new TH2F("hCentVsMultRPS", " Centrality Vs. Multiplicity RPs",5000, 0, 5000.,12,0,60 );
693 hCentVsMultRPS->GetXaxis()->SetTitle("Multiplicity RPs");
694 hCentVsMultRPS->GetYaxis()->SetTitle("Centrality");
695 fOutputFlowObs->Add(hCentVsMultRPS);
d48d81b5 696 }
a4b0ffc7 697// AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
698// AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
699// AliPIDResponse *pidResp=inputHandler->GetPIDResponse();
700// fCuts->GetPidHF()->SetPidResponse(pidResp);
0dbb51a0 701 // Post the data
702 PostData(1,fNEntries);
d48d81b5 703
5fc4893f 704 if(fOnOff[1]) PostData(2,fOutputPID);
705 if(fOnOff[0]) PostData(3,fOutputTrack);
0dbb51a0 706 PostData(4,fCuts);
5fc4893f 707 if(fOnOff[2]) PostData(5,fOutputCounters);
818c1271 708 if(fOnOff[3]) PostData(7,fOutputEvSelection);
d48d81b5 709 if(fOnOff[4]) PostData(8,fOutputFlowObs);
818c1271 710
5fc4893f 711 if(!fOnOff[0] && !fOnOff[1] && !fOnOff[2]) AliError("Nothing will be filled!");
0dbb51a0 712}
713
714//___________________________________________________________________________
715void AliAnalysisTaskSEHFQA::UserExec(Option_t */*option*/)
716{
717 // Execute analysis for current event
718
719 AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
720 if(fDebug>2) printf("Analysing decay %d\n",fDecayChannel);
721 // Post the data already here
722 PostData(1,fNEntries);
5fc4893f 723 if(fOnOff[1]) PostData(2,fOutputPID);
724 if(fOnOff[0]) PostData(3,fOutputTrack);
0dbb51a0 725 PostData(4,fCuts);
5fc4893f 726 if(fOnOff[2]) {
727 PostData(5,fOutputCounters);
728 if(fCuts->GetUseCentrality()) PostData(6,fOutputCheckCentrality);
729 }
0dbb51a0 730
731 TClonesArray *arrayProng =0;
732 Int_t pdg=0;
a4ef4383 733 Int_t *pdgdaughters=0x0;
734
0dbb51a0 735 if(!aod && AODEvent() && IsStandardAOD()) {
736 // In case there is an AOD handler writing a standard AOD, use the AOD
737 // event in memory rather than the input (ESD) event.
738 aod = dynamic_cast<AliAODEvent*> (AODEvent());
739 // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
740 // have to taken from the AOD event hold by the AliAODExtension
741 AliAODHandler* aodHandler = (AliAODHandler*)
742 ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
743 if(aodHandler->GetExtensions()) {
744
745 AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
746 AliAODEvent *aodFromExt = ext->GetAOD();
747
748
749
750 switch(fDecayChannel){
751 case 0:
752 arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
753 pdg=411;
a4ef4383 754 if(fReadMC){
755 pdgdaughters =new Int_t[3];
756 pdgdaughters[0]=211;//pi
757 pdgdaughters[1]=321;//K
758 pdgdaughters[2]=211;//pi
759 }
0dbb51a0 760 break;
761 case 1:
762 arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("D0toKpi");
763 pdg=421;
a4ef4383 764 if(fReadMC){
765 pdgdaughters =new Int_t[2];
766 pdgdaughters[0]=211;//pi
767 pdgdaughters[1]=321;//K
768 }
0dbb51a0 769 break;
770 case 2:
771 arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Dstar");
772 pdg=413;
a4ef4383 773 if(fReadMC){
fb9311e8 774 pdgdaughters =new Int_t[3];
a4ef4383 775 pdgdaughters[1]=211;//pi
776 pdgdaughters[0]=321;//K
777 pdgdaughters[2]=211;//pi (soft?)
778 }
0dbb51a0 779 break;
780 case 3:
781 arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
782 pdg=431;
a4ef4383 783 if(fReadMC){
784 pdgdaughters =new Int_t[3];
785 pdgdaughters[0]=321;//K
786 pdgdaughters[1]=321;//K
787 pdgdaughters[2]=211;//pi
788 }
0dbb51a0 789 break;
790 case 4:
791 arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm4Prong");
792 pdg=421;
a4ef4383 793 if(fReadMC){
794 pdgdaughters =new Int_t[4];
795 pdgdaughters[0]=321;
796 pdgdaughters[1]=211;
797 pdgdaughters[2]=211;
798 pdgdaughters[3]=211;
799 }
0dbb51a0 800 break;
801 case 5:
802 arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
803 pdg=4122;
a4ef4383 804 if(fReadMC){
805 pdgdaughters =new Int_t[3];
806 pdgdaughters[0]=2212;//p
807 pdgdaughters[1]=321;//K
808 pdgdaughters[2]=211;//pi
809 }
0dbb51a0 810 break;
811 }
812 }
dc222f77 813 } else if(aod) {
0dbb51a0 814 switch(fDecayChannel){
815 case 0:
816 arrayProng=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
817 pdg=411;
a4ef4383 818 if(fReadMC){
819 pdgdaughters =new Int_t[3];
820 pdgdaughters[0]=211;//pi
821 pdgdaughters[1]=321;//K
822 pdgdaughters[2]=211;//pi
823 }
0dbb51a0 824 break;
825 case 1:
826 arrayProng=(TClonesArray*)aod->GetList()->FindObject("D0toKpi");
827 pdg=421;
a4ef4383 828 if(fReadMC){
829 pdgdaughters =new Int_t[2];
830 pdgdaughters[0]=211;//pi
831 pdgdaughters[1]=321;//K
832 }
0dbb51a0 833 break;
834 case 2:
835 arrayProng=(TClonesArray*)aod->GetList()->FindObject("Dstar");
836 pdg=413;
a4ef4383 837 if(fReadMC){
24bcb2f3 838 pdgdaughters =new Int_t[3];
a4ef4383 839 pdgdaughters[1]=211;//pi
840 pdgdaughters[0]=321;//K
841 pdgdaughters[2]=211;//pi (soft?)
842 }
0dbb51a0 843 break;
844 case 3:
845 arrayProng=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
846 pdg=431;
a4ef4383 847 if(fReadMC){
848 pdgdaughters =new Int_t[3];
849 pdgdaughters[0]=321;//K
850 pdgdaughters[1]=321;//K
851 pdgdaughters[2]=211;//pi
852 }
0dbb51a0 853 break;
854 case 4:
855 arrayProng=(TClonesArray*)aod->GetList()->FindObject("Charm4Prong");
856 pdg=421;
a4ef4383 857 if(fReadMC){
858 pdgdaughters =new Int_t[4];
859 pdgdaughters[0]=321;
860 pdgdaughters[1]=211;
861 pdgdaughters[2]=211;
862 pdgdaughters[3]=211;
863 }
0dbb51a0 864 break;
865 case 5:
866 arrayProng=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
867 pdg=4122;
a4ef4383 868 if(fReadMC){
869 pdgdaughters =new Int_t[3];
870 pdgdaughters[0]=2212;//p
871 pdgdaughters[1]=321;//K
872 pdgdaughters[2]=211;//pi
873 }
0dbb51a0 874 break;
875 }
876 }
19222b31 877 Bool_t isSimpleMode=fSimpleMode;
0dbb51a0 878 if(!arrayProng) {
879 AliInfo("Branch not found! The output will contain only trak related histograms\n");
880 isSimpleMode=kTRUE;
9af24f46 881 fNEntries->Fill(2);
0dbb51a0 882 }
883
a4ef4383 884 TClonesArray *mcArray = 0;
885 AliAODMCHeader *mcHeader = 0;
886
3655ecf2 887 if(!aod) {
888 delete [] pdgdaughters;
889 return;
890 }
891
a4ef4383 892 //check if MC
893 if(fReadMC) {
894 // load MC particles
895 mcArray = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
896 if(!mcArray) {
897 printf("AliAnalysisTaskSEHFQA::UserExec: MC particles branch not found!\n");
fb9311e8 898 delete [] pdgdaughters;
a4ef4383 899 return;
900 }
901
902 // load MC header
903 mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
904 if(!mcHeader) {
905 printf("AliAnalysisTaskSEHFQA::UserExec: MC header branch not found!\n");
fb9311e8 906 delete [] pdgdaughters;
a4ef4383 907 return;
908 }
909 }
23dfe729 910
911
912 UInt_t evSelMask=((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
913 Double_t centrality=fCuts->GetCentrality(aod);
914 Double_t multiplicity=aod->GetHeader()->GetRefMultiplicity();
17531bbd 915 Int_t runNumber = aod->GetRunNumber();
5e67f11a 916 TString trigClass=aod->GetFiredTriggerClasses();
138e0d89 917 Int_t nAODtracks=aod->GetNTracks();
918 Int_t nSelTracksTPCOnly=0;
919 Int_t nSelTracksTPCITS=0;
920 Int_t nSelTracksTPCITS1SPD=0;
921
922 for (Int_t k=0;k<nAODtracks;k++){
923 AliAODTrack* track=aod->GetTrack(k);
924 if(track->GetID()<0) continue;
925 Int_t nclsTot=0,nclsSPD=0;
926 for(Int_t l=0;l<6;l++) {
927 if(TESTBIT(track->GetITSClusterMap(),l)) {
928 nclsTot++; if(l<2) nclsSPD++;
929 }
930 }
931 UShort_t nTPCClus=track->GetTPCClusterMap().CountBits();
932 if(TMath::Abs(track->Eta())<0.8 && nTPCClus>=70 && track->GetStatus()&AliESDtrack::kTPCrefit){
933 if(track->TestFilterBit(1)) nSelTracksTPCOnly++;
934 if(track->GetStatus()&AliESDtrack::kITSrefit){
935 nSelTracksTPCITS++;
936 if(nclsSPD>0) nSelTracksTPCITS1SPD++;
937 }
938 }
939 }
940
d48d81b5 941 if(fOnOff[4]) {
942 FillFlowObs(aod);
943 PostData(8,fOutputFlowObs);
944 }
1b765719 945 if(fOnOff[3]){
946 TH2F* hTrigC=(TH2F*)fOutputEvSelection->FindObject("hTrigCent");
947 TH2F* hTrigM=(TH2F*)fOutputEvSelection->FindObject("hTrigMul");
948 AliCounterCollection* trigCount=(AliCounterCollection*)fOutputEvSelection->FindObject("trigCounter");
17531bbd 949
2a844546 950 hTrigC->Fill(-1.,centrality);
951 hTrigM->Fill(-1.,multiplicity);
da3b34af 952 trigCount->Count(Form("triggerType:All/Run:%d",runNumber));
5e67f11a 953 if(evSelMask==0){
954 if(aod->GetEventType()!=7){
955 trigCount->Count(Form("triggerType:NoPhysSelEvNot7/Run:%d",runNumber));
956 }else if(trigClass.Contains("CMUP1")){
957 trigCount->Count(Form("triggerType:NoPhysSelCMUP1/Run:%d",runNumber));
958 }else if(trigClass.Contains("MUON")){
959 trigCount->Count(Form("triggerType:NoPhysSelMUON/Run:%d",runNumber));
960 }else if(trigClass.Contains("CPBI2_B1-B") || trigClass.Contains(" CPBI2WU_B1-B")){
961 trigCount->Count(Form("triggerType:NoPhysSelMB/Run:%d",runNumber));
962 }else if(trigClass.Contains("CCENT") || trigClass.Contains("CVHN")){
963 trigCount->Count(Form("triggerType:NoPhysSelCent/Run:%d",runNumber));
964 }else if(trigClass.Contains("CSEMI") || trigClass.Contains("CVLN")){
965 trigCount->Count(Form("triggerType:NoPhysSelSemiCent/Run:%d",runNumber));
966 }
967 }
2a844546 968 if(evSelMask & AliVEvent::kAny){
969 hTrigC->Fill(0.,centrality);
970 hTrigM->Fill(0.,multiplicity);
971 trigCount->Count(Form("triggerType:Any/Run:%d",runNumber));
972 }
1b765719 973 if(evSelMask & AliVEvent::kMB){
974 hTrigC->Fill(1.,centrality);
975 hTrigM->Fill(1.,multiplicity);
976 trigCount->Count(Form("triggerType:MB/Run:%d",runNumber));
977 }
978 if(evSelMask & AliVEvent::kINT7){
979 hTrigC->Fill(2.,centrality);
980 hTrigM->Fill(2.,multiplicity);
981 }
982 if(evSelMask & AliVEvent::kCINT5){
983 hTrigC->Fill(3.,centrality);
984 hTrigM->Fill(3.,multiplicity);
985 }
986 if(evSelMask & AliVEvent::kCentral){
987 hTrigC->Fill(4.,centrality);
988 hTrigM->Fill(4.,multiplicity);
989 trigCount->Count(Form("triggerType:Cent/Run:%d",runNumber));
990 }
991 if(evSelMask & AliVEvent::kSemiCentral){
992 hTrigC->Fill(5.,centrality);
993 hTrigM->Fill(5.,multiplicity);
994 trigCount->Count(Form("triggerType:SemiCent/Run:%d",runNumber));
995 }
996 if(evSelMask & (AliVEvent::kEMC1 | AliVEvent::kEMC7)){
997 hTrigC->Fill(6.,centrality);
998 hTrigM->Fill(6.,multiplicity);
999 }
1000 if(evSelMask & (AliVEvent::kEMCEJE | AliVEvent::kEMCEGA)){
1001 hTrigC->Fill(7.,centrality);
1002 hTrigM->Fill(7.,multiplicity);
1003 trigCount->Count(Form("triggerType:EMCAL/Run:%d",runNumber));
1004 }
1005 if(evSelMask & (((AliVEvent::kCMUS5 | AliVEvent::kMUSH7) | (AliVEvent::kMUL7 | AliVEvent::kMUU7)) | (AliVEvent::kMUS7 | AliVEvent::kMUON))){
1006 hTrigC->Fill(8.,centrality);
1007 hTrigM->Fill(8.,multiplicity);
5e67f11a 1008 trigCount->Count(Form("triggerType:MUON/Run:%d",runNumber));
1b765719 1009 }
1010 if(evSelMask & (AliVEvent::kPHI1 | AliVEvent::kPHI7)){
1011 hTrigC->Fill(9.,centrality);
1012 hTrigM->Fill(9.,multiplicity);
1013 }
1014 if(evSelMask & (AliVEvent::kDG5 | AliVEvent::kZED)){
1015 hTrigC->Fill(10.,centrality);
1016 hTrigM->Fill(10.,multiplicity);
1017 }
23dfe729 1018 }
1019
1020
0dbb51a0 1021 // fix for temporary bug in ESDfilter
1022 // the AODs with null vertex pointer didn't pass the PhysSel
3655ecf2 1023 if(!aod->GetPrimaryVertex() || TMath::Abs(aod->GetMagneticField())<0.001) {
1024 delete [] pdgdaughters;
1025 return;
1026 }
0dbb51a0 1027
1028 // count event
1029 fNEntries->Fill(0);
9af24f46 1030
1031 //count events with good vertex
1032 // AOD primary vertex
1033 AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
824b448a 1034
1035 Double_t pos[3],cov[6];
1036 vtx1->GetXYZ(pos);
1037 vtx1->GetCovarianceMatrix(cov);
1038 const AliESDVertex vESD(pos,cov,100.,100);
1039
9af24f46 1040 TString primTitle = vtx1->GetTitle();
1041 if(primTitle.Contains("VertexerTracks") && vtx1->GetNContributors()>0) fNEntries->Fill(4);
1042
a3aa1279 1043 // trigger class for PbPb C0SMH-B-NOPF-ALLNOTRD, C0SMH-B-NOPF-ALL
5fc4893f 1044 //TString trigclass=aod->GetFiredTriggerClasses();
1045 //if(trigclass.Contains("C0SMH-B-NOPF-ALLNOTRD") || trigclass.Contains("C0SMH-B-NOPF-ALL")) fNEntries->Fill(5); //tmp
a3aa1279 1046
818c1271 1047
23dfe729 1048
1049
2a844546 1050 Bool_t evSelbyCentrality=kTRUE,evSelected=kTRUE,evSelByVertex=kTRUE,evselByPileup=kTRUE,evSelByPS=kTRUE;
0dbb51a0 1051 //select event
1052 if(!fCuts->IsEventSelected(aod)) {
82f89d0d 1053 evSelected=kFALSE;
2a844546 1054 if(fCuts->GetWhyRejection()==1) {fNEntries->Fill(1); evselByPileup=kFALSE;}// rejected for pileup
ac0c2841 1055 if(fCuts->GetWhyRejection()==2 || fCuts->GetWhyRejection()==3) evSelbyCentrality=kFALSE; //rejected by centrality
1056 if(fCuts->GetWhyRejection()==4) evSelByVertex=kFALSE; //rejected by vertex
5fc4893f 1057 if(fCuts->GetWhyRejection()==5) fNEntries->Fill(5);//tmp
00c00d8e 1058 if(fCuts->GetWhyRejection()==6 && fOnOff[3]) ((AliCounterCollection*)fOutputEvSelection->FindObject("evselection"))->Count(Form("evnonsel:zvtx/Run:%d",runNumber));
2a844546 1059 if(fCuts->GetWhyRejection()==7) { evSelByPS=kFALSE; }
0dbb51a0 1060 }
23dfe729 1061 if(evSelected){
1062 TH2F* hTrigS=(TH2F*)fOutputEvSelection->FindObject("hTrigCentSel");
2a844546 1063 hTrigS->Fill(-1.,centrality);
23dfe729 1064
2a844546 1065 if(evSelMask & AliVEvent::kAny) hTrigS->Fill(0.,centrality);
23dfe729 1066 if(evSelMask & AliVEvent::kMB) hTrigS->Fill(1.,centrality);
1067 if(evSelMask & AliVEvent::kINT7) hTrigS->Fill(2.,centrality);
1068 if(evSelMask & AliVEvent::kCINT5) hTrigS->Fill(3.,centrality);
1069 if(evSelMask & AliVEvent::kCentral) hTrigS->Fill(4.,centrality);
1070 if(evSelMask & AliVEvent::kSemiCentral) hTrigS->Fill(5.,centrality);
1071 if(evSelMask & (AliVEvent::kEMC1 | AliVEvent::kEMC7)) hTrigS->Fill(6.,centrality);
1072 if(evSelMask & (AliVEvent::kEMCEJE | AliVEvent::kEMCEGA)) hTrigS->Fill(7.,centrality);
1073 if(evSelMask & (((AliVEvent::kCMUS5 | AliVEvent::kMUSH7) | (AliVEvent::kMUL7 | AliVEvent::kMUU7)) | (AliVEvent::kMUS7 | AliVEvent::kMUON))) hTrigS->Fill(8.,centrality);
1074 if(evSelMask & (AliVEvent::kPHI1 | AliVEvent::kPHI7)) hTrigS->Fill(9.,centrality);
1075 if(evSelMask & (AliVEvent::kDG5 | AliVEvent::kZED)) hTrigS->Fill(10.,centrality);
1076 }
1077
2a844546 1078 if(evSelected || (!evSelbyCentrality && evSelByVertex && evselByPileup && evSelByPS)){ //events selected or not selected because of centrality
5fc4893f 1079 if(fOnOff[2] && fCuts->GetUseCentrality()){
818c1271 1080
ac0c2841 1081 Float_t stdCentf=fCuts->GetCentrality(aod);
1082 Int_t stdCent = (Int_t)(stdCentf+0.5);
5fc4893f 1083 Float_t secondCentf =fCuts->GetCentrality(aod,fEstimator);
1084 Int_t secondCent = (Int_t)(secondCentf+0.5);
82f89d0d 1085 Int_t mincent=stdCent-stdCent%10;
b68a0f6a 1086 AliAODVZERO *vzeroAOD = (AliAODVZERO*)aod->GetVZEROData();
1087 Float_t vzeroMult = vzeroAOD->GetMTotV0A() + vzeroAOD->GetMTotV0C();
1088 AliCentrality *aodcent = aod->GetCentrality();
1089 Float_t spdCentf = aodcent->GetCentralityPercentile("CL1");
5fc4893f 1090 if(stdCentf==-1) {
5587ae1a 1091 mincent=-10;
5fc4893f 1092 stdCent=-1;
1093 }
c40c91c6 1094 if(mincent==100)mincent--;
82f89d0d 1095 ((AliCounterCollection*)fOutputCounters->FindObject("stdEstimator"))->Count(Form("centralityclass:%d_%d/Run:%d",mincent,mincent+10,runNumber));
5fc4893f 1096
82f89d0d 1097 mincent=secondCent-secondCent%10;
5fc4893f 1098 if(secondCentf==-1) {
5587ae1a 1099 mincent=-10;
5fc4893f 1100 secondCent=-1;
1101 }
c40c91c6 1102 if(mincent==100)mincent--;
82f89d0d 1103 ((AliCounterCollection*)fOutputCounters->FindObject("secondEstimator"))->Count(Form("centralityclass:%d_%d/Run:%d",mincent,mincent+10,runNumber));
1104
1105 if(stdCent<fCuts->GetMinCentrality() || stdCent>fCuts->GetMaxCentrality()){
1106 ((TH1F*)fOutputCheckCentrality->FindObject("hNtrackletsOut"))->Fill(aod->GetTracklets()->GetNumberOfTracklets());
1107 ((TH1F*)fOutputCheckCentrality->FindObject("hMultOut"))->Fill(aod->GetHeader()->GetRefMultiplicity());
1108 }else{
1109 ((TH1F*)fOutputCheckCentrality->FindObject("hNtrackletsIn"))->Fill(aod->GetTracklets()->GetNumberOfTracklets());
1110 ((TH1F*)fOutputCheckCentrality->FindObject("hMultIn"))->Fill(aod->GetHeader()->GetRefMultiplicity());
1111 }
ac0c2841 1112 ((TH2F*)fOutputCheckCentrality->FindObject("hMultvsPercentile"))->Fill(aod->GetHeader()->GetRefMultiplicity(),stdCentf);
7a5f2502 1113 ((TH2F*)fOutputCheckCentrality->FindObject("hntrklvsPercentile"))->Fill(aod->GetTracklets()->GetNumberOfTracklets(),stdCentf);
138e0d89 1114 ((TH2F*)fOutputCheckCentrality->FindObject("hnTPCTracksvsPercentile"))->Fill(nSelTracksTPCOnly,stdCentf);
1115 ((TH2F*)fOutputCheckCentrality->FindObject("hnTPCITSTracksvsPercentile"))->Fill(nSelTracksTPCITS,stdCentf);
1116 ((TH2F*)fOutputCheckCentrality->FindObject("hnTPCITS1SPDTracksvsPercentile"))->Fill(nSelTracksTPCITS1SPD,stdCentf);
b68a0f6a 1117 ((TH2F*)fOutputCheckCentrality->FindObject("hV0MultiplicityPercentile"))->Fill(vzeroMult,stdCentf);
1118 ((TH2F*)fOutputCheckCentrality->FindObject("hV0MultiplicityNtrackletsIn"))->Fill(vzeroMult,aod->GetTracklets()->GetNumberOfTracklets());
1119 ((TH2F*)fOutputCheckCentrality->FindObject("hStdPercentileSPDPercentile"))->Fill(stdCentf,spdCentf);
0dbb51a0 1120
82f89d0d 1121 PostData(6,fOutputCheckCentrality);
1122
1123 } else{
5fc4893f 1124 if(fOnOff[0]){
1125 ((TH1F*)fOutputTrack->FindObject("hNtracklets"))->Fill(aod->GetTracklets()->GetNumberOfTracklets());
1126 ((TH1F*)fOutputTrack->FindObject("hMult"))->Fill(aod->GetHeader()->GetRefMultiplicity());
1127 }
a4ef4383 1128 }
82f89d0d 1129 }
1130
818c1271 1131 if(fOnOff[3]){
1132 const AliVVertex *vertex = aod->GetPrimaryVertex();
1133 Double_t zvtx=vertex->GetZ();
a4b0ffc7 1134 if(evSelected || (!evSelected && TMath::Abs(zvtx) > 10.))
818c1271 1135 ((TH1F*)fOutputEvSelection->FindObject("hzvtx"))->Fill(zvtx);
1136 }
1137
82f89d0d 1138 if(!evSelected) {
1139 delete [] pdgdaughters;
1140 return; //discard all events not selected (vtx and/or centrality)
a4ef4383 1141 }
1142
818c1271 1143
476db98d 1144 AliAODPidHF* pidHF=fCuts->GetPidHF();
362242bf 1145 if(!pidHF) {
1146 delete [] pdgdaughters;
1147 return;
1148 }
4940d5bf 1149 //AliPIDResponse* respF=pidHF->GetPidResponse();
476db98d 1150 AliTPCPIDResponse* tpcres=new AliTPCPIDResponse();
4940d5bf 1151 Bool_t oldPID=pidHF->GetOldPid();
1152 if(oldPID){
70e398d3 1153 Double_t alephParameters[5];
1154 pidHF->GetTPCBetheBlochParams(alephParameters);
1155 tpcres->SetBetheBlochParameters(alephParameters[0],alephParameters[1],alephParameters[2],alephParameters[3],alephParameters[4]);
1156 }
476db98d 1157
1158
0dbb51a0 1159 Int_t ntracks=0;
824b448a 1160 Int_t isGoodTrack=0, isFakeTrack=0, isSelTrack=0;
0dbb51a0 1161
1162 if(aod) ntracks=aod->GetNTracks();
1163
5fc4893f 1164 if(fOnOff[0] || fOnOff[1]){
1165 //loop on tracks in the event
1166 for (Int_t k=0;k<ntracks;k++){
1167 AliAODTrack* track=aod->GetTrack(k);
138e0d89 1168 if(track->GetID()<0) continue;
5fc4893f 1169 AliAODPid *pid = track->GetDetPid();
138e0d89 1170 if(!pid && fDebug>1) cout<<"No AliAODPid found"<<endl;
5fc4893f 1171
138e0d89 1172 if(pid && fOnOff[1]){
5fc4893f 1173 Double_t times[AliPID::kSPECIES];
1174 pid->GetIntegratedTimes(times);
0dbb51a0 1175
5fc4893f 1176 Double_t tofRes[AliPID::kSPECIES];
1177 pid->GetTOFpidResolution(tofRes);
1178
1179 //check TOF
b906dda8 1180 TH1F* htmpfl=((TH1F*)fOutputPID->FindObject("hTOFflags"));
1181 htmpfl->Fill(0.);
23dfe729 1182 if (track->GetStatus()&AliESDtrack::kTPCout) htmpfl->Fill(1.);
1183 if (track->GetStatus()&AliESDtrack::kTOFout) htmpfl->Fill(2.);
1184 if (track->GetStatus()&AliESDtrack::kTIME) htmpfl->Fill(3.);
1185 if (track->GetStatus()&AliESDtrack::kTOFpid) htmpfl->Fill(4.);
1186 if (track->GetStatus()&AliESDtrack::kTOFmismatch) htmpfl->Fill(5.);
b906dda8 1187
5fc4893f 1188 if(pidHF && pidHF->CheckStatus(track,"TOF")){
1189 ((TH1F*)fOutputPID->FindObject("hTOFtime"))->Fill(times[AliPID::kProton]);
1190 ((TH2F*)fOutputPID->FindObject("hTOFtimeKaonHyptime"))->Fill(track->P(),pid->GetTOFsignal()-times[3]); //3 is kaon
1191 ((TH1F*)fOutputPID->FindObject("hTOFsig"))->Fill(pid->GetTOFsignal());
1192 if (pid->GetTOFsignal()< 0) ((TH1F*)fOutputPID->FindObject("hTOFsig"))->Fill(-1);
1193
4940d5bf 1194 Double_t nsigma[3]={-10,-10,-10};
1195 pidHF->GetnSigmaTOF(track,(Int_t)AliPID::kPion,nsigma[0]);
1196 pidHF->GetnSigmaTOF(track,(Int_t)AliPID::kKaon,nsigma[1]);
1197 pidHF->GetnSigmaTOF(track,(Int_t)AliPID::kProton,nsigma[2]);
1198
1199 ((TH2F*)fOutputPID->FindObject("hTOFsigmaKSigPid"))->Fill(track->P(),nsigma[1]);
1200 ((TH2F*)fOutputPID->FindObject("hTOFsigmaPionSigPid"))->Fill(track->P(),nsigma[0]);
1201 ((TH2F*)fOutputPID->FindObject("hTOFsigmaProtonSigPid"))->Fill(track->P(),nsigma[2]);
1202 if(fReadMC){
1203 Int_t label=track->GetLabel();
1204 if(label<=0) continue;
1205 AliMCParticle* mcpart=(AliMCParticle*)mcArray->At(label);
1206 if(mcpart){
1207 Int_t abspdgcode=TMath::Abs(mcpart->PdgCode());
1208 if(abspdgcode==211) ((TH2F*)fOutputPID->FindObject("hTOFsigmaMCPionSigPid"))->Fill(track->P(),nsigma[0]);
1209 if(abspdgcode==321) ((TH2F*)fOutputPID->FindObject("hTOFsigmaMCKSigPid"))->Fill(track->P(),nsigma[1]);
1210 if(abspdgcode==2212) ((TH2F*)fOutputPID->FindObject("hTOFsigmaMCProtonSigPid"))->Fill(track->P(),nsigma[2]);
24748c57 1211
24748c57 1212 }
4940d5bf 1213 }
24748c57 1214
4940d5bf 1215 for (Int_t iS=2; iS<5; iS++){ //we plot TOF Pid resolution for 3-sigma identified particles
1216 if ( TMath::Abs(nsigma[iS-2])<3.){
1217 switch (iS) {
1218 case AliPID::kPion:
1219 ((TH1F*)fOutputPID->FindObject("hTOFsigPid3sigPion"))->Fill(tofRes[iS]);
1220 break;
1221 case AliPID::kKaon:
1222 ((TH1F*)fOutputPID->FindObject("hTOFsigPid3sigKaon"))->Fill(tofRes[iS]);
1223 break;
1224 case AliPID::kProton:
1225 ((TH1F*)fOutputPID->FindObject("hTOFsigPid3sigProton"))->Fill(tofRes[iS]);
1226 break;
1227 default:
1228 break;
5fc4893f 1229 }
a4ef4383 1230 }
1231 }
5fc4893f 1232 }//if TOF status
4940d5bf 1233 //}
5fc4893f 1234
1235 if(pidHF && pidHF->CheckStatus(track,"TPC")){
476db98d 1236
5fc4893f 1237 Double_t TPCp=pid->GetTPCmomentum();
1238 Double_t TPCsignal=pid->GetTPCsignal();
1239 ((TH1F*)fOutputPID->FindObject("hTPCsig"))->Fill(TPCsignal);
1240 ((TH1F*)fOutputPID->FindObject("hTPCsigvsp"))->Fill(TPCp,TPCsignal);
1241 //if (pidHF->IsKaonRaw(track, "TOF"))
4940d5bf 1242 Double_t nsigma[3]={-10,-10,-10};
1243 pidHF->GetnSigmaTPC(track,(Int_t)AliPID::kPion,nsigma[0]);
1244 pidHF->GetnSigmaTPC(track,(Int_t)AliPID::kKaon,nsigma[1]);
1245 pidHF->GetnSigmaTPC(track,(Int_t)AliPID::kProton,nsigma[2]);
24748c57 1246
4940d5bf 1247 ((TH2F*)fOutputPID->FindObject("hTPCsigmaK"))->Fill(TPCp,nsigma[1]);
1248
1249 ((TH2F*)fOutputPID->FindObject("hTPCsigmaPion"))->Fill(TPCp,nsigma[0]);
1250 ((TH2F*)fOutputPID->FindObject("hTPCsigmaProton"))->Fill(TPCp,nsigma[2]);
1251
24748c57 1252 if(fReadMC){
1253 Int_t label=track->GetLabel();
1254 if(label<=0) continue;
1255 AliMCParticle* mcpart=(AliMCParticle*)mcArray->At(label);
1256 if(mcpart){
1257 Int_t abspdgcode=TMath::Abs(mcpart->PdgCode());
1258 if(abspdgcode==211) ((TH2F*)fOutputPID->FindObject("hTPCsigmaMCPion"))->Fill(track->P(),nsigma[0]);
1259 if(abspdgcode==321) ((TH2F*)fOutputPID->FindObject("hTPCsigmaMCK"))->Fill(track->P(),nsigma[1]);
1260 if(abspdgcode==2212) ((TH2F*)fOutputPID->FindObject("hTPCsigmaMCProton"))->Fill(track->P(),nsigma[2]);
1261
1262 }
1263
1264 }
1265
1266
5fc4893f 1267 }//if TPC status
1268 } //end PID histograms
1269
1270 Int_t nclsTot=0,nclsSPD=0;
1271
1272 //check clusters of the tracks
1273 if(fOnOff[0]){
1274
bc6846db 1275 ((TH1F*)fOutputTrack->FindObject("hnLayerITS"))->Fill(-1);
5fc4893f 1276 for(Int_t l=0;l<6;l++) {
1277 if(TESTBIT(track->GetITSClusterMap(),l)) {
bc6846db 1278 ((TH1F*)fOutputTrack->FindObject("hnLayerITS"))->Fill(l);
5fc4893f 1279 nclsTot++; if(l<2) nclsSPD++;
1280 }
1281 }
1282 ((TH1F*)fOutputTrack->FindObject("hnClsITS"))->Fill(nclsTot);
1283 ((TH1F*)fOutputTrack->FindObject("hnClsSPD"))->Fill(nclsSPD);
39dd297f 1284 if(track->Pt()>0.3 &&
1285 TMath::Abs(track->Eta())<0.8 &&
1286 track->GetStatus()&AliESDtrack::kITSrefit &&
1287 track->GetStatus()&AliESDtrack::kTPCrefit &&
1288 nclsSPD>0){
1289 ((TH1F*)fOutputTrack->FindObject("hnClsITSselTr"))->Fill(nclsTot);
1290 }
5fc4893f 1291 if(!(track->GetStatus()&AliESDtrack::kTPCin) && track->GetStatus()&AliESDtrack::kITSrefit && !(track->GetStatus()&AliESDtrack::kITSpureSA)){//tracks retrieved in the ITS and not reconstructed in the TPC
1292 ((TH1F*)fOutputTrack->FindObject("hnClsITS-SA"))->Fill(nclsTot);
bc6846db 1293 ((TH1F*)fOutputTrack->FindObject("hnLayerITS"))->Fill(-1);
1294 for(Int_t l=0;l<6;l++) {
1295 if(TESTBIT(track->GetITSClusterMap(),l)) {
1296 ((TH1F*)fOutputTrack->FindObject("hnLayerITSsa"))->Fill(l);
1297 }
1298 }
acc33385 1299 }
5fc4893f 1300 Int_t label=0;
1301 if(fReadMC){
1302 label=track->GetLabel();
1303 if (label<0)fNEntries->Fill(8);
1304 else fNEntries->Fill(9);
acc33385 1305 }
0dbb51a0 1306
5fc4893f 1307
824b448a 1308 if (track->Pt()>0.3 &&
1309 track->GetStatus()&AliESDtrack::kTPCrefit &&
1310 track->GetStatus()&AliESDtrack::kITSrefit &&
1311 /*nclsTot>3 &&*/
1312 nclsSPD>0) {//count good tracks
5fc4893f 1313
1314
824b448a 1315 if(fReadMC && label<0) {
1316 ((TH1F*)fOutputTrack->FindObject("hptFakeTr"))->Fill(track->Pt());
1317 isFakeTrack++;
1318 } else {
1319 ((TH1F*)fOutputTrack->FindObject("hptGoodTr"))->Fill(track->Pt());
1320 isGoodTrack++;
ad42e35b 1321 }
824b448a 1322 }
1323 if(fCuts->IsDaughterSelected(track,&vESD,fCuts->GetTrackCuts())){
1324 isSelTrack++;
1325 }//select tracks for our analyses
818c1271 1326 } //fill track histos
5fc4893f 1327 } //end loop on tracks
0dbb51a0 1328
4940d5bf 1329 //fill once per event
5fc4893f 1330 if(fOnOff[0]){
1331 if (fReadMC) ((TH1F*)fOutputTrack->FindObject("hdistrFakeTr"))->Fill(isFakeTrack);
1332 ((TH1F*)fOutputTrack->FindObject("hdistrGoodTr"))->Fill(isGoodTrack);
824b448a 1333 ((TH1F*)fOutputTrack->FindObject("hdistrSelTr"))->Fill(isSelTrack);
5fc4893f 1334 }
0dbb51a0 1335
5fc4893f 1336 if(!isSimpleMode){
1337 // loop over candidates
1338 Int_t nCand = arrayProng->GetEntriesFast();
1339 Int_t ndaugh=3;
1340 if(fDecayChannel==AliAnalysisTaskSEHFQA::kD0toKpi) ndaugh=2;
1341 if(fDecayChannel==AliAnalysisTaskSEHFQA::kD0toKpipipi) ndaugh=4;
1342
1343 for (Int_t iCand = 0; iCand < nCand; iCand++) {
1344 AliAODRecoDecayHF *d = (AliAODRecoDecayHF*)arrayProng->UncheckedAt(iCand);
1345 if(d->GetSelectionMap()) {
1346 if(fDecayChannel==AliAnalysisTaskSEHFQA::kD0toKpi && !d->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts)) continue; //skip the D0 from Dstar
1347 if(fDecayChannel==AliAnalysisTaskSEHFQA::kDplustoKpipi && !d->HasSelectionBit(AliRDHFCuts::kDplusCuts)) continue; //skip the 3 prong !D+
1348 }
acc33385 1349
5fc4893f 1350 if(fReadMC){
1351 Int_t labD = d->MatchToMC(pdg,mcArray,ndaugh,pdgdaughters);
1352 if(labD>=0){
1353 AliAODMCParticle *partD = (AliAODMCParticle*)mcArray->At(labD);
1354 Int_t label=partD->GetMother();
1355 AliAODMCParticle *mot = (AliAODMCParticle*)mcArray->At(label);
1356 while(label>=0){//get first mother
1357 mot = (AliAODMCParticle*)mcArray->At(label);
1358 label=mot->GetMother();
1359 }
24748c57 1360 if(mot){
1361 Int_t pdgMotCode = mot->GetPdgCode();
5fc4893f 1362
24748c57 1363 if(TMath::Abs(pdgMotCode)==4) fNEntries->Fill(6); //from primary charm
1364 if(TMath::Abs(pdgMotCode)==5) fNEntries->Fill(7); //from beauty
1365 }
acc33385 1366 }
5fc4893f 1367 }//end MC
1368 else fNEntries->Fill(6); //count the candidates (data)
1369
1370 for(Int_t id=0;id<ndaugh;id++){
1371
1372 //other histograms to be filled when the cut object is given
1373 AliAODTrack* track=(AliAODTrack*)d->GetDaughter(id);
1374
1375 //track quality
1376
1377 if (fCuts->IsInFiducialAcceptance(d->Pt(),d->Y(pdg)) && fCuts->IsSelected(d,AliRDHFCuts::kTracks,aod)) {
1378
1379 Int_t label=0;
1380 if(fReadMC)label=track->GetLabel();
ad42e35b 1381 if(fOnOff[0]){
818c1271 1382
ad42e35b 1383 if(fReadMC && label<0) {
1384 isFakeTrack++;
824b448a 1385 ((TH1F*)fOutputTrack->FindObject("hptFakeTrFromDaugh"))->Fill(track->Pt());
5fc4893f 1386
ad42e35b 1387 ((TH1F*)fOutputTrack->FindObject("hd0f"))->Fill(d->Getd0Prong(id));
1388 } else {
824b448a 1389 ((TH1F*)fOutputTrack->FindObject("hptGoodTrFromDaugh"))->Fill(track->Pt());
ad42e35b 1390 ((TH1F*)fOutputTrack->FindObject("hd0"))->Fill(d->Getd0Prong(id));
4563d678 1391 Double_t d0rphiz[2],covd0[3];
1392 Bool_t isDCA=track->PropagateToDCA(aod->GetPrimaryVertex(),aod->GetMagneticField(),9999.,d0rphiz,covd0);
1393 if(isDCA) ((TH1F*)fOutputTrack->FindObject("hd0z"))->Fill(d0rphiz[1]);
ad42e35b 1394 }
818c1271 1395 }
5fc4893f 1396 if (fCuts->IsSelected(d,AliRDHFCuts::kAll,aod) && fOnOff[1]){
0dbb51a0 1397
5fc4893f 1398 AliAODPid *pid = track->GetDetPid();
138e0d89 1399 if(pid){
1400 Double_t times[5];
1401 pid->GetIntegratedTimes(times);
1402 if(pidHF && pidHF->CheckStatus(track,"TOF")) ((TH2F*)fOutputPID->FindObject("hTOFtimeKaonHyptimeAC"))->Fill(track->P(),pid->GetTOFsignal()-times[AliPID::kKaon]);
1403 if(pidHF && pidHF->CheckStatus(track,"TPC")) ((TH2F*)fOutputPID->FindObject("hTPCsigvspAC"))->Fill(pid->GetTPCmomentum(),pid->GetTPCsignal());
1404 }
5fc4893f 1405 fNEntries->Fill(3);
1406 } //end analysis cuts
1407 } //end acceptance and track cuts
1408 } //end loop on tracks in the candidate
1409 } //end loop on candidates
824b448a 1410
5fc4893f 1411 }
1412 } //end if on pid or track histograms
fb9311e8 1413
476db98d 1414 delete tpcres;
fb9311e8 1415 delete [] pdgdaughters;
a4ef4383 1416 PostData(1,fNEntries);
5fc4893f 1417 if(fOnOff[1]) PostData(2,fOutputPID);
1418 if(fOnOff[0]) PostData(3,fOutputTrack);
a4ef4383 1419 PostData(4,fCuts);
5fc4893f 1420 if(fOnOff[2]) PostData(5,fOutputCounters);
1421 //Post data 6 done in case of centrality on
82f89d0d 1422
0dbb51a0 1423}
1424
d48d81b5 1425//____________________________________________________________________________
1426void AliAnalysisTaskSEHFQA::FillFlowObs(AliAODEvent *aod){
1427 //fills the flow observables
65274ad6 1428 Double_t cc;
1429 cc = fCuts->GetCentrality(aod);
1430 ((TH2F*) fOutputFlowObs->FindObject("hFlowEvents"))->Fill(0., cc);
d48d81b5 1431
1432 UInt_t mask=((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
1433 UInt_t trigger=AliVEvent::kMB | AliVEvent::kCentral | AliVEvent::kSemiCentral;
d48d81b5 1434 if(mask & trigger) {
65274ad6 1435 ((TH2F*) fOutputFlowObs->FindObject("hFlowEvents"))->Fill(1.,cc); // fired
1436 if (mask & AliVEvent::kMB) ((TH2F*) fOutputFlowObs->FindObject("hFlowEvents"))->Fill(2.,cc);
1437 if (mask & AliVEvent::kCentral) ((TH2F*) fOutputFlowObs->FindObject("hFlowEvents"))->Fill(3.,cc);
1438 if (mask & AliVEvent::kSemiCentral) ((TH2F*) fOutputFlowObs->FindObject("hFlowEvents"))->Fill(4.,cc);
d48d81b5 1439 Bool_t rejected=false;
d48d81b5 1440 if(cc<0 || cc>60) rejected=true;
1441 const AliVVertex *vertex = aod->GetPrimaryVertex();
1442 Double_t zvtx=vertex->GetZ();
65274ad6 1443 if(TMath::Abs(zvtx)>fCuts->GetMaxVtxZ()) rejected=true;
d48d81b5 1444 if(rejected) return; //not interesting for flow QA
1445 } else {
1446 return;
1447 }
1448
1449 // event accepted
65274ad6 1450 ((TH2F*) fOutputFlowObs->FindObject("hFlowEvents"))->Fill(5.,cc);
d48d81b5 1451 fRFPcuts->SetParamType(AliFlowTrackCuts::kGlobal);
1452 fRFPcuts->SetPtRange(0.2,5.);
1453 fRFPcuts->SetEtaRange(-0.8,0.8);
1454 fRFPcuts->SetMinNClustersTPC(70);
1455 fRFPcuts->SetMinChi2PerClusterTPC(0.2);
1456 fRFPcuts->SetMaxChi2PerClusterTPC(4.0);
1457 fRFPcuts->SetAcceptKinkDaughters(kFALSE);
1458 fRFPcuts->SetEvent(aod);
1459
1460 TString ref[3] = {"FB1","FB128","VZE"};
1461 Double_t psi[3];
1462 for(Int_t i=0; i!=3; ++i) {
1463 if(i==0) { // switching to bit 1
1464 fRFPcuts->SetMinimalTPCdedx(10.);
1465 fRFPcuts->SetAODfilterBit(1);
1466 } else { // switching to bit 128
1467 fRFPcuts->SetMinimalTPCdedx(-1);
1468 fRFPcuts->SetAODfilterBit(128);
1469 }
1470 if(i>1) {
1471 fRFPcuts->SetParamType(AliFlowTrackCuts::kV0);
1472 fRFPcuts->SetEtaRange(-5,+5);
1473 fRFPcuts->SetPhiMin(0);
1474 fRFPcuts->SetPhiMax(TMath::TwoPi());
1475 }
1476 fFlowEvent->Fill(fRFPcuts,fRFPcuts);
1477 fFlowEvent->TagSubeventsInEta(-5,0,0,+5);
65274ad6 1478 // getting informationt
d48d81b5 1479 AliFlowVector vQ, vQaQb[2];
1480 fFlowEvent->Get2Qsub(vQaQb,2);
1481 vQ = vQaQb[0]+vQaQb[1];
1482 Double_t dMa=vQaQb[0].GetMult();
1483 Double_t dMb=vQaQb[1].GetMult();
1484 if( dMa<2 || dMb<2 ) {
65274ad6 1485 ((TH2F*) fOutputFlowObs->FindObject("hFlowEvents"))->Fill(6.,cc); //???
d48d81b5 1486 continue;
1487 }
1488 psi[i] = vQ.Phi()/2;
1489 // publishing
1490 ((TProfile2D*) fOutputFlowObs->FindObject( Form("h%s_Q",ref[i].Data())))->Fill(0,cc,vQaQb[0].X()/dMa,dMa); // Qx-
1491 ((TProfile2D*) fOutputFlowObs->FindObject( Form("h%s_Q",ref[i].Data())))->Fill(1,cc,vQaQb[0].Y()/dMa,dMa); // Qy-
1492 ((TProfile2D*) fOutputFlowObs->FindObject( Form("h%s_Q",ref[i].Data())))->Fill(2,cc,vQaQb[1].X()/dMb,dMb); // Qx+
1493 ((TProfile2D*) fOutputFlowObs->FindObject( Form("h%s_Q",ref[i].Data())))->Fill(3,cc,vQaQb[1].Y()/dMb,dMb); // Qy+
1494 ((TH2F*) fOutputFlowObs->FindObject( Form("h%s_AngleQ",ref[i].Data()) ))->Fill(psi[i],cc); // Psi
1495 AliFlowTrackSimple *track;
1496 for(Int_t t=0; t!=fFlowEvent->NumberOfTracks(); ++t) {
1497 track = (AliFlowTrackSimple*) fFlowEvent->GetTrack(t);
1498 if(!track) continue;
1499 if(!track->InRPSelection()) continue;
1500 ((TH3F*) fOutputFlowObs->FindObject( Form("h%s_PhiEta",ref[i].Data()) ))->Fill(track->Phi(),track->Eta(),cc,track->Weight()); //PhiEta
1501 }
65274ad6 1502
1503 //histo filled only for TPCFB1
1504 if (i==0) {
1505 ((TH2F*) fOutputFlowObs->FindObject("hCentVsMultRPS"))->Fill(fFlowEvent->GetNumberOfRPs(),cc);
1506 }
d48d81b5 1507 }
1508 // TPC vs VZERO
1509 ((TH3F*) fOutputFlowObs->FindObject( "hTPCVZE_AngleQ" ))->Fill(psi[0],psi[2],cc);
1510}
1511
a4ef4383 1512//____________________________________________________________________________
0dbb51a0 1513void AliAnalysisTaskSEHFQA::Terminate(Option_t */*option*/){
1514 //terminate analysis
1515
1516 fNEntries = dynamic_cast<TH1F*>(GetOutputData(1));
1517 if(!fNEntries){
1518 printf("ERROR: %s not available\n",GetOutputSlot(1)->GetContainer()->GetName());
1519 return;
1520 }
1521
1522 fOutputPID = dynamic_cast<TList*> (GetOutputData(2));
5fc4893f 1523 if (!fOutputPID && fOnOff[1]) {
0dbb51a0 1524 printf("ERROR: %s not available\n",GetOutputSlot(2)->GetContainer()->GetName());
1525 return;
1526 }
1527
1528 fOutputTrack = dynamic_cast<TList*> (GetOutputData(3));
5fc4893f 1529 if (!fOutputTrack && fOnOff[0]) {
0dbb51a0 1530 printf("ERROR: %s not available\n",GetOutputSlot(3)->GetContainer()->GetName());
1531 return;
1532 }
1533
1534}
1535