PWG3/vertexingHF/AliAnalysisTaskSEHFQA.cxx

   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
  16 /* $Id$ */
  17
  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>
  31 #include <TDatabasePDG.h>
  32
  33 #include <AliAnalysisDataSlot.h>
  34 #include <AliAnalysisDataContainer.h>
  35 #include "AliAnalysisManager.h"
  36 #include "AliESDtrack.h"
  37 #include "AliVertexerTracks.h"
  38 #include "AliPID.h"
  39 #include "AliTPCPIDResponse.h"
  40 #include "AliAODHandler.h"
  41 #include "AliAODEvent.h"
  42 #include "AliAODVertex.h"
  43 #include "AliAODTrack.h"
  44 #include "AliAODMCParticle.h"
  45 #include "AliAODMCHeader.h"
  46 #include "AliAODRecoDecayHF2Prong.h"
  47 #include "AliAODRecoCascadeHF.h"
  48 #include "AliAnalysisVertexingHF.h"
  49 #include "AliAnalysisTaskSE.h"
  50 #include "AliCounterCollection.h"
  51 #include "AliRDHFCuts.h"
  52 #include "AliRDHFCutsDplustoKpipi.h"
  53 #include "AliRDHFCutsD0toKpipipi.h"
  54 #include "AliRDHFCutsDstoKKpi.h"
  55 #include "AliRDHFCutsDStartoKpipi.h"
  56 #include "AliRDHFCutsD0toKpi.h"
  57 #include "AliRDHFCutsLctopKpi.h"
  58
  59 #include "AliAnalysisTaskSEHFQA.h"
  60
  61 ClassImp(AliAnalysisTaskSEHFQA)
  62
  63 //____________________________________________________________________________
  64
  65 AliAnalysisTaskSEHFQA::AliAnalysisTaskSEHFQA():AliAnalysisTaskSE(),
  66   fNEntries(0x0),
  67   fOutputPID(0x0),
  68   fOutputTrack(0x0),
  69   fOutputCounters(0x0),
  70   fOutputCheckCentrality(0x0),
  71   fOutputEvSelection(0x0),
  72   fDecayChannel(AliAnalysisTaskSEHFQA::kD0toKpi),
  73   fCuts(0x0),
  74   fEstimator(AliRDHFCuts::kCentTRK),
  75   fReadMC(kFALSE),
  76   fSimpleMode(kFALSE),
  77   fOnOff()
  78 {
  79   //default constructor
  80   fOnOff[0]=kTRUE;
  81   fOnOff[1]=kTRUE;
  82   fOnOff[2]=kTRUE;
  83   fOnOff[3]=kTRUE;
  84 }
  85
  86 //____________________________________________________________________________
  87 AliAnalysisTaskSEHFQA::AliAnalysisTaskSEHFQA(const char *name, AliAnalysisTaskSEHFQA::DecChannel ch,AliRDHFCuts* cuts):
  88   AliAnalysisTaskSE(name),
  89   fNEntries(0x0),
  90   fOutputPID(0x0),
  91   fOutputTrack(0x0),
  92   fOutputCounters(0x0),
  93   fOutputCheckCentrality(0x0),
  94   fOutputEvSelection(0x0),
  95   fDecayChannel(ch),
  96   fCuts(0x0),
  97   fEstimator(AliRDHFCuts::kCentTRK),
  98   fReadMC(kFALSE),
  99   fSimpleMode(kFALSE),
 100   fOnOff()
 101 {
 102   //constructor
 103
 104   //SetCutObject(cuts);
 105   fCuts=cuts;
 106
 107   fOnOff[0]=kTRUE;
 108   fOnOff[1]=kTRUE;
 109   fOnOff[2]=kTRUE;
 110   fOnOff[3]=kTRUE;
 111
 112   // Output slot #1 writes into a TH1F container (number of events)
 113   DefineOutput(1,TH1F::Class());  //My private output
 114   // Output slot #2 writes into a TList container (PID)
 115   if (fOnOff[1]) DefineOutput(2,TList::Class());  //My private output
 116   // Output slot #3 writes into a TList container (Tracks)
 117   if (fOnOff[0]) DefineOutput(3,TList::Class());  //My private output
 118   // Output slot #4 writes into a AliRDHFCuts container (cuts)
 119   switch(fDecayChannel){
 120   case 0:
 121     DefineOutput(4,AliRDHFCutsDplustoKpipi::Class());  //My private output
 122     break;
 123   case 1:
 124     DefineOutput(4,AliRDHFCutsD0toKpi::Class());  //My private output
 125     break;
 126   case 2:
 127     DefineOutput(4,AliRDHFCutsDStartoKpipi::Class());  //My private output
 128     break;
 129   case 3:
 130     DefineOutput(4,AliRDHFCutsDstoKKpi::Class());  //My private output
 131     break;
 132   case 4:
 133     DefineOutput(4,AliRDHFCutsD0toKpipipi::Class());  //My private output
 134     break;
 135   case 5:
 136     DefineOutput(4,AliRDHFCutsLctopKpi::Class());  //My private output
 137     break;
 138   }
 139   if (fOnOff[2]) {
 140     // Output slot #5 writes into a TList container (AliCounterCollection)
 141     DefineOutput(5,TList::Class());  //My private output
 142     // Output slot #6 writes into a TList container (TH1F)
 143     DefineOutput(6,TList::Class());  //My private output
 144   }
 145
 146   if(fOnOff[3]) DefineOutput(7,TList::Class());  //My private output
 147
 148 }
 149
 150 //___________________________________________________________________________
 151 AliAnalysisTaskSEHFQA::~AliAnalysisTaskSEHFQA()
 152 {
 153   //destructor
 154   delete fNEntries;
 155
 156   delete fOutputPID;
 157
 158   delete fOutputTrack;
 159
 160   delete fOutputCounters;
 161
 162   delete fOutputCheckCentrality;
 163
 164   delete fOutputEvSelection;
 165
 166 }
 167
 168 //___________________________________________________________________________
 169 void AliAnalysisTaskSEHFQA::Init(){
 170
 171   //initialization
 172   if(fDebug > 1) printf("AnalysisTaskSEHFQA::Init() \n");
 173
 174   switch(fDecayChannel){
 175   case 0:
 176     {
 177       AliRDHFCutsDplustoKpipi* copycut=new AliRDHFCutsDplustoKpipi(*(static_cast<AliRDHFCutsDplustoKpipi*>(fCuts)));
 178       // Post the data
 179       PostData(4,copycut);
 180     }
 181     break;
 182   case 1:
 183     {
 184       AliRDHFCutsD0toKpi* copycut=new AliRDHFCutsD0toKpi(*(static_cast<AliRDHFCutsD0toKpi*>(fCuts)));
 185       // Post the data
 186       PostData(4,copycut);
 187     }
 188     break;
 189   case 2:
 190     {
 191       AliRDHFCutsDStartoKpipi* copycut=new AliRDHFCutsDStartoKpipi(*(static_cast<AliRDHFCutsDStartoKpipi*>(fCuts)));
 192       // Post the data
 193       PostData(4,copycut);
 194     }
 195     break;
 196   case 3:
 197     {
 198       AliRDHFCutsDstoKKpi* copycut=new AliRDHFCutsDstoKKpi(*(static_cast<AliRDHFCutsDstoKKpi*>(fCuts)));
 199       // Post the data
 200       PostData(4,copycut);
 201     }
 202     break;
 203   case 4:
 204     {
 205       AliRDHFCutsD0toKpipipi* copycut=new AliRDHFCutsD0toKpipipi(*(static_cast<AliRDHFCutsD0toKpipipi*>(fCuts)));
 206       // Post the data
 207       PostData(4,copycut);
 208     }
 209     break;
 210   case 5:
 211     {
 212       AliRDHFCutsLctopKpi* copycut=new AliRDHFCutsLctopKpi(*(static_cast<AliRDHFCutsLctopKpi*>(fCuts)));
 213       // Post the data
 214       PostData(4,copycut);
 215     }
 216     break;
 217
 218   default:
 219     return;
 220   }
 221
 222
 223
 224 }
 225
 226 //___________________________________________________________________________
 227 void AliAnalysisTaskSEHFQA::UserCreateOutputObjects()
 228 {
 229
 230   //create the output container
 231   if(fDebug > 1) printf("AnalysisTaskSEHFQA::UserCreateOutputObjects() \n");
 232
 233   //count events
 234
 235   fNEntries=new TH1F(GetOutputSlot(1)->GetContainer()->GetName(), "Counts the number of events", 10,-0.5,9.5);
 236   fNEntries->GetXaxis()->SetBinLabel(1,"nEventsAnal");
 237   fNEntries->GetXaxis()->SetBinLabel(2,"Pile-up Rej");
 238   fNEntries->GetXaxis()->SetBinLabel(3,"No VertexingHF");
 239   fNEntries->GetXaxis()->SetBinLabel(4,"nCandidates(AnCuts)");
 240   fNEntries->GetXaxis()->SetBinLabel(5,"EventsWithGoodVtx");
 241   //fNEntries->GetXaxis()->SetBinLabel(6,"N. of 0SMH");
 242   fNEntries->GetXaxis()->SetBinLabel(6,"N. of CSH1");
 243   if(fReadMC){
 244     fNEntries->GetXaxis()->SetBinLabel(7,"MC Cand from c");
 245     fNEntries->GetXaxis()->SetBinLabel(8,"MC Cand from b");
 246     fNEntries->GetXaxis()->SetBinLabel(9,"N fake Trks");
 247     fNEntries->GetXaxis()->SetBinLabel(10,"N true Trks");
 248   } else{
 249     fNEntries->GetXaxis()->SetBinLabel(7,"N candidates");
 250   }
 251
 252   fNEntries->GetXaxis()->SetNdivisions(1,kFALSE);
 253
 254   //PID
 255   if(fOnOff[1]){
 256     fOutputPID=new TList();
 257     fOutputPID->SetOwner();
 258     fOutputPID->SetName(GetOutputSlot(2)->GetContainer()->GetName());
 259
 260     //TOF pid
 261     TString hname="hTOFsig";
 262     TH1F* hTOFsig=new TH1F(hname.Data(),"Distribution of TOF signal;TOF time [ps];Entries", 100, -2.e3,40.e3);
 263
 264     hname="hTOFtime";
 265     TH1F* hTOFtime=new TH1F(hname.Data(),"Distribution of TOF time Kaon;TOF time(Kaon) [ps];Entries", 1000, 0.,50000.);
 266
 267     hname="hTOFtimeKaonHyptime";
 268     TH2F* hTOFtimeKaonHyptime=new TH2F(hname.Data(),"TOFtime - timeHypothesisForKaon;p[GeV/c];TOFtime - timeHypothesisForKaon [ps]",200,0.,4.,1000,-20000.,20000.);
 269
 270     hname="hTOFtimeKaonHyptimeAC";
 271     TH2F* hTOFtimeKaonHyptimeAC=new TH2F(hname.Data(),"TOFtime - timeHypothesisForKaon;p[GeV/c];TOFtime - timeHypothesisForKaon [ps]",200,0.,4.,1000,-20000.,20000.);
 272
 273     hname="hTOFsigmaKSigPid";
 274     TH2F* hTOFsigmaKSigPid=new TH2F(hname.Data(),"(TOFsignal-timeK)/tofSigPid;p[GeV/c];(TOFsignal-timeK)/tofSigPid",200,0.,4.,100,-5,5);
 275
 276     hname="hTOFsigmaPionSigPid";
 277     TH2F* hTOFsigmaPionSigPid=new TH2F(hname.Data(),"(TOFsignal-time#pi)/tofSigPid;p[GeV/c];(TOFsignal-time#pi)/tofSigPid",200,0.,4.,100,-5,5);
 278
 279     hname="hTOFsigmaProtonSigPid";
 280     TH2F* hTOFsigmaProtonSigPid=new TH2F(hname.Data(),"(TOFsignal-timep)/tofSigPid;p[GeV/c];(TOFsignal-time p)/tofSigPid",200,0.,4.,100,-5,5);
 281
 282     hname="hTOFsigPid3sigPion";
 283     TH1F* hTOFsigPid3sigPion=new TH1F(hname.Data(),"TOF PID resolution (#pi) [ps]",500,0.,1000.);
 284
 285     hname="hTOFsigPid3sigKaon";
 286     TH1F* hTOFsigPid3sigKaon=new TH1F(hname.Data(),"TOF PID resolution (K) [ps]",500,0.,1000.);
 287
 288     hname="hTOFsigPid3sigProton";
 289     TH1F* hTOFsigPid3sigProton=new TH1F(hname.Data(),"TOF PID resolution (p) [ps]",500,0.,1000.);
 290
 291
 292     //TPC pid
 293     hname="hTPCsig";
 294     TH1F* hTPCsig=new TH1F(hname.Data(),"Distribution of TPC signal;TPC sig;Entries", 100, 35.,100.);
 295
 296     hname="hTPCsigvsp";
 297     TH2F* hTPCsigvsp=new TH2F(hname.Data(),"TPCsig vs p;TPC p[GeV/c];TPCsig",200,0.,4.,1000,35.,100.);
 298
 299     hname="hTPCsigvspAC";
 300     TH2F* hTPCsigvspAC=new TH2F(hname.Data(),"TPCsig vs p;TPCp[GeV/c];TPCsig",200,0.,4.,1000,35.,100.);
 301
 302     hname="hTPCsigmaK";
 303     TH2F* hTPCsigmaK=new TH2F(hname.Data(),"TPC Sigma for K as a function of momentum;p[GeV/c];Sigma Kaon",200,0.,4.,200,-5,5);
 304
 305     hname="hTPCsigmaPion";
 306     TH2F* hTPCsigmaPion=new TH2F(hname.Data(),"TPC Sigma for #pi as a function of momentum;p[GeV/c];Sigma #pi",200,0.,4.,200,-5,5);
 307
 308     hname="hTPCsigmaProton";
 309     TH2F* hTPCsigmaProton=new TH2F(hname.Data(),"TPC Sigma for proton as a function of momentum;p[GeV/c];Sigma Proton",200,0.,4.,200,-5,5);
 310
 311     fOutputPID->Add(hTOFsig);
 312     fOutputPID->Add(hTPCsig);
 313     fOutputPID->Add(hTOFtime);
 314     fOutputPID->Add(hTOFtimeKaonHyptime);
 315     fOutputPID->Add(hTOFtimeKaonHyptimeAC);
 316     fOutputPID->Add(hTOFsigmaKSigPid);
 317     fOutputPID->Add(hTOFsigmaPionSigPid);
 318     fOutputPID->Add(hTOFsigmaProtonSigPid);
 319     fOutputPID->Add(hTOFsigPid3sigPion);
 320     fOutputPID->Add(hTOFsigPid3sigKaon);
 321     fOutputPID->Add(hTOFsigPid3sigProton);
 322     fOutputPID->Add(hTPCsigvsp);
 323     fOutputPID->Add(hTPCsigvspAC);
 324     fOutputPID->Add(hTPCsigmaK);
 325     fOutputPID->Add(hTPCsigmaPion);
 326     fOutputPID->Add(hTPCsigmaProton);
 327   }
 328
 329   //quality of the tracks
 330   if(fOnOff[0]){
 331     fOutputTrack=new TList();
 332     fOutputTrack->SetOwner();
 333     fOutputTrack->SetName(GetOutputSlot(3)->GetContainer()->GetName());
 334
 335     TString hname="hnClsITS";
 336     TH1F* hnClsITS=new TH1F(hname.Data(),"Distribution of number of ITS clusters;nITScls;Entries",7,-0.5,6.5);
 337
 338     hname="hnClsITS-SA";
 339     TH1F* hnClsITSSA=new TH1F(hname.Data(),"Distribution of number of ITS clusters(ITS-SA);nITScls;Entries",7,-0.5,6.5);
 340
 341
 342     hname="hnLayerITS";
 343     TH1F* hnLayerITS=new TH1F(hname.Data(),"Number of tracks with point in layer;ITS layer;",7,-1.5,5.5);
 344     hnLayerITS->GetXaxis()->SetBinLabel(1,"n tracks");
 345
 346     hname="hnLayerITSsa";
 347     TH1F* hnLayerITSsa=new TH1F(hname.Data(),"Number of tracks with point in layer;ITS layer;",7,-1.5,5.5);
 348     hnLayerITSsa->GetXaxis()->SetBinLabel(1,"n tracks");
 349
 350     hname="hnClsSPD";
 351     TH1F* hnClsSPD=new TH1F(hname.Data(),"Distribution of number of SPD clusters;nSPDcls;Entries",3,-0.5,2.5);
 352
 353     hname="hptGoodTr";
 354     TH1F* hptGoodTr=new TH1F(hname.Data(),"Pt distribution of 'good' tracks;p_{t}[GeV];Entries/0.05 GeV/c",400,0.,20.);
 355     hptGoodTr->SetTitleOffset(1.3,"Y");
 356
 357     hname="hdistrGoodTr";
 358     TH1F* hdistrGoodTr=new TH1F(hname.Data(),"Distribution of number of 'good' tracks per event;no.good-tracks/ev;Entries",4000,-0.5,3999.5);
 359     hdistrGoodTr->SetTitleOffset(1.3,"Y");
 360
 361     hname="hd0";
 362     TH1F* hd0=new TH1F(hname.Data(),"Impact parameter distribution of 'good' tracks;d_{0}[cm];Entries/10^{3} cm",200,-0.1,0.1);
 363
 364     fOutputTrack->Add(hnClsITS);
 365     fOutputTrack->Add(hnClsITSSA);
 366     fOutputTrack->Add(hnLayerITS);
 367     fOutputTrack->Add(hnLayerITSsa);
 368     fOutputTrack->Add(hnClsSPD);
 369     fOutputTrack->Add(hptGoodTr);
 370     fOutputTrack->Add(hdistrGoodTr);
 371     fOutputTrack->Add(hd0);
 372
 373     if(fReadMC){
 374       hname="hdistrFakeTr";
 375       TH1F* hdistrFakeTr=new TH1F(hname.Data(),"Distribution of number of fake tracks per event;no.fake-tracks/ev;Entries",4000,-0.5,3999.5);
 376       hdistrGoodTr->SetTitleOffset(1.3,"Y");
 377
 378       hname="hd0f";
 379       TH1F* hd0f=new TH1F(hname.Data(),"Impact parameter distribution of fake tracks;d_{0}[cm];Entries/10^{3} cm",200,-0.1,0.1);
 380
 381       hname="hptFakeTr";
 382       TH1F* hptFakeTr=new TH1F(hname.Data(),"Pt distribution of fake tracks;p_{t}[GeV];Entries/0.05 GeV/c",400,0.,20.);
 383       hptFakeTr->SetTitleOffset(1.3,"Y");
 384
 385       fOutputTrack->Add(hptFakeTr);
 386       fOutputTrack->Add(hdistrFakeTr);
 387       fOutputTrack->Add(hd0f);
 388
 389     }
 390   }
 391
 392
 393   if(fOnOff[2] && fCuts->GetUseCentrality()){
 394
 395     //Centrality (Counters)
 396     fOutputCounters=new TList();
 397     fOutputCounters->SetOwner();
 398     fOutputCounters->SetName(GetOutputSlot(5)->GetContainer()->GetName());
 399
 400     AliCounterCollection *stdEstimator=new AliCounterCollection("stdEstimator");
 401     stdEstimator->AddRubric("run",500000);
 402     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");
 403     stdEstimator->Init();
 404     AliCounterCollection *secondEstimator=new AliCounterCollection("secondEstimator");
 405     secondEstimator->AddRubric("run",500000);
 406     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");
 407     secondEstimator->Init();
 408
 409     fOutputCounters->Add(stdEstimator);
 410     fOutputCounters->Add(secondEstimator);
 411
 412     //Centrality (Checks)
 413     fOutputCheckCentrality=new TList();
 414     fOutputCheckCentrality->SetOwner();
 415     fOutputCheckCentrality->SetName(GetOutputSlot(6)->GetContainer()->GetName());
 416
 417     TString hname="hNtrackletsIn";
 418     TH1F* hNtrackletsIn=new TH1F(hname.Data(),"Number of tracklets in Centrality range;ntracklets;Entries",5000,-0.5,4999.5);
 419
 420     hname="hMultIn";
 421     TH1F* hMultIn=new TH1F(hname.Data(),"Multiplicity;multiplicity in Centrality range;Entries",10000,-0.5,9999.5);
 422
 423     hname="hNtrackletsOut";
 424     TH1F* hNtrackletsOut=new TH1F(hname.Data(),"Number of tracklets out of Centrality range;ntracklets;Entries",5000,-0.5,4999.5);
 425
 426     hname="hMultOut";
 427     TH1F* hMultOut=new TH1F(hname.Data(),"Multiplicity out of Centrality range;multiplicity;Entries",10000,-0.5,9999.5);
 428
 429     hname="hMultvsPercentile";
 430     TH2F* hMultvsPercentile=new TH2F(hname.Data(),"Multiplicity vs Percentile;multiplicity;percentile",10000,-0.5,9999.5,12,-10.,110);
 431
 432     fOutputCheckCentrality->Add(hNtrackletsIn);
 433     fOutputCheckCentrality->Add(hNtrackletsOut);
 434     fOutputCheckCentrality->Add(hMultIn);
 435     fOutputCheckCentrality->Add(hMultOut);
 436     fOutputCheckCentrality->Add(hMultvsPercentile);
 437
 438     PostData(6,fOutputCheckCentrality);
 439
 440   } else{
 441     if(fOnOff[0]){
 442       TString hname="hNtracklets";
 443       TH1F* hNtracklets=new TH1F(hname.Data(),"Number of tracklets;ntracklets;Entries",5000,-0.5,4999.5);
 444
 445       hname="hMult";
 446       TH1F* hMult=new TH1F(hname.Data(),"Multiplicity;multiplicity;Entries",10000,-0.5,9999.5);
 447       fOutputTrack->Add(hNtracklets);
 448       fOutputTrack->Add(hMult);
 449     }
 450   }
 451
 452   //event selection (z vertex for the moment)
 453   if(fOnOff[3]){
 454     fOutputEvSelection=new TList();
 455     fOutputEvSelection->SetOwner();
 456     fOutputEvSelection->SetName(GetOutputSlot(7)->GetContainer()->GetName());
 457     AliCounterCollection *evselection=new AliCounterCollection("evselection");
 458     evselection->AddRubric("run",500000);
 459     evselection->AddRubric("evnonsel","zvtx");
 460     evselection->Init();
 461
 462     TH1F* hzvtx=new TH1F("hzvtx", "Distribution of z_{VTX};z_{VTX} [cm];Entries",100,-20,20);
 463
 464     fOutputEvSelection->Add(evselection);
 465     fOutputEvSelection->Add(hzvtx);
 466   }
 467
 468   // Post the data
 469   PostData(1,fNEntries);
 470   if(fOnOff[1]) PostData(2,fOutputPID);
 471   if(fOnOff[0]) PostData(3,fOutputTrack);
 472   PostData(4,fCuts);
 473   if(fOnOff[2]) PostData(5,fOutputCounters);
 474   if(fOnOff[3]) PostData(7,fOutputEvSelection);
 475
 476   if(!fOnOff[0] && !fOnOff[1] && !fOnOff[2]) AliError("Nothing will be filled!");
 477 }
 478
 479 //___________________________________________________________________________
 480 void AliAnalysisTaskSEHFQA::UserExec(Option_t */*option*/)
 481 {
 482   // Execute analysis for current event
 483
 484   AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
 485   if(fDebug>2) printf("Analysing decay %d\n",fDecayChannel);
 486   // Post the data already here
 487   PostData(1,fNEntries);
 488   if(fOnOff[1]) PostData(2,fOutputPID);
 489   if(fOnOff[0]) PostData(3,fOutputTrack);
 490   PostData(4,fCuts);
 491   if(fOnOff[2]) {
 492     PostData(5,fOutputCounters);
 493     if(fCuts->GetUseCentrality()) PostData(6,fOutputCheckCentrality);
 494   }
 495
 496   TClonesArray *arrayProng =0;
 497   Int_t pdg=0;
 498   Int_t *pdgdaughters=0x0;
 499
 500   if(!aod && AODEvent() && IsStandardAOD()) {
 501     // In case there is an AOD handler writing a standard AOD, use the AOD
 502     // event in memory rather than the input (ESD) event.
 503     aod = dynamic_cast<AliAODEvent*> (AODEvent());
 504     // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
 505     // have to taken from the AOD event hold by the AliAODExtension
 506     AliAODHandler* aodHandler = (AliAODHandler*)
 507       ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
 508     if(aodHandler->GetExtensions()) {
 509
 510       AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
 511       AliAODEvent *aodFromExt = ext->GetAOD();
 512
 513
 514
 515       switch(fDecayChannel){
 516       case 0:
 517         arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
 518         pdg=411;
 519         if(fReadMC){
 520           pdgdaughters =new Int_t[3];
 521           pdgdaughters[0]=211;//pi
 522           pdgdaughters[1]=321;//K
 523           pdgdaughters[2]=211;//pi
 524         }
 525         break;
 526       case 1:
 527         arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("D0toKpi");
 528         pdg=421;
 529         if(fReadMC){
 530           pdgdaughters =new Int_t[2];
 531           pdgdaughters[0]=211;//pi
 532           pdgdaughters[1]=321;//K
 533         }
 534         break;
 535       case 2:
 536         arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Dstar");
 537         pdg=413;
 538         if(fReadMC){
 539           pdgdaughters =new Int_t[3];
 540           pdgdaughters[1]=211;//pi
 541           pdgdaughters[0]=321;//K
 542           pdgdaughters[2]=211;//pi (soft?)
 543         }
 544         break;
 545       case 3:
 546         arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
 547         pdg=431;
 548         if(fReadMC){
 549           pdgdaughters =new Int_t[3];
 550           pdgdaughters[0]=321;//K
 551           pdgdaughters[1]=321;//K
 552           pdgdaughters[2]=211;//pi
 553         }
 554         break;
 555       case 4:
 556         arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm4Prong");
 557         pdg=421;
 558         if(fReadMC){
 559           pdgdaughters =new Int_t[4];
 560           pdgdaughters[0]=321;
 561           pdgdaughters[1]=211;
 562           pdgdaughters[2]=211;
 563           pdgdaughters[3]=211;
 564         }
 565         break;
 566       case 5:
 567         arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
 568         pdg=4122;
 569         if(fReadMC){
 570           pdgdaughters =new Int_t[3];
 571           pdgdaughters[0]=2212;//p
 572           pdgdaughters[1]=321;//K
 573           pdgdaughters[2]=211;//pi
 574         }
 575         break;
 576       }
 577     }
 578   } else if(aod) {
 579     switch(fDecayChannel){
 580     case 0:
 581       arrayProng=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
 582       pdg=411;
 583       if(fReadMC){
 584         pdgdaughters =new Int_t[3];
 585         pdgdaughters[0]=211;//pi
 586         pdgdaughters[1]=321;//K
 587         pdgdaughters[2]=211;//pi
 588       }
 589       break;
 590     case 1:
 591       arrayProng=(TClonesArray*)aod->GetList()->FindObject("D0toKpi");
 592       pdg=421;
 593       if(fReadMC){
 594         pdgdaughters =new Int_t[2];
 595         pdgdaughters[0]=211;//pi
 596         pdgdaughters[1]=321;//K
 597       }
 598       break;
 599     case 2:
 600       arrayProng=(TClonesArray*)aod->GetList()->FindObject("Dstar");
 601       pdg=413;
 602       if(fReadMC){
 603         pdgdaughters =new Int_t[3];
 604         pdgdaughters[1]=211;//pi
 605         pdgdaughters[0]=321;//K
 606         pdgdaughters[2]=211;//pi (soft?)
 607       }
 608       break;
 609     case 3:
 610       arrayProng=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
 611       pdg=431;
 612       if(fReadMC){
 613         pdgdaughters =new Int_t[3];
 614         pdgdaughters[0]=321;//K
 615         pdgdaughters[1]=321;//K
 616         pdgdaughters[2]=211;//pi
 617       }
 618       break;
 619     case 4:
 620       arrayProng=(TClonesArray*)aod->GetList()->FindObject("Charm4Prong");
 621       pdg=421;
 622       if(fReadMC){
 623         pdgdaughters =new Int_t[4];
 624         pdgdaughters[0]=321;
 625         pdgdaughters[1]=211;
 626         pdgdaughters[2]=211;
 627         pdgdaughters[3]=211;
 628       }
 629       break;
 630     case 5:
 631       arrayProng=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
 632       pdg=4122;
 633       if(fReadMC){
 634         pdgdaughters =new Int_t[3];
 635         pdgdaughters[0]=2212;//p
 636         pdgdaughters[1]=321;//K
 637         pdgdaughters[2]=211;//pi
 638       }
 639       break;
 640     }
 641   }
 642   Bool_t isSimpleMode=fSimpleMode;
 643   if(!arrayProng) {
 644     AliInfo("Branch not found! The output will contain only trak related histograms\n");
 645     isSimpleMode=kTRUE;
 646     fNEntries->Fill(2);
 647   }
 648
 649   TClonesArray *mcArray = 0;
 650   AliAODMCHeader *mcHeader = 0;
 651
 652   if(!aod) {
 653     delete [] pdgdaughters;
 654     return;
 655   }
 656
 657   //check if MC
 658   if(fReadMC) {
 659     // load MC particles
 660     mcArray = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
 661     if(!mcArray) {
 662       printf("AliAnalysisTaskSEHFQA::UserExec: MC particles branch not found!\n");
 663       delete [] pdgdaughters;
 664       return;
 665     }
 666
 667     // load MC header
 668     mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
 669     if(!mcHeader) {
 670       printf("AliAnalysisTaskSEHFQA::UserExec: MC header branch not found!\n");
 671       delete [] pdgdaughters;
 672       return;
 673     }
 674   }
 675   // fix for temporary bug in ESDfilter
 676   // the AODs with null vertex pointer didn't pass the PhysSel
 677   if(!aod->GetPrimaryVertex() || TMath::Abs(aod->GetMagneticField())<0.001) {
 678     delete [] pdgdaughters;
 679     return;
 680   }
 681
 682   // count event
 683   fNEntries->Fill(0);
 684
 685   //count events with good vertex
 686   // AOD primary vertex
 687   AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
 688   TString primTitle = vtx1->GetTitle();
 689   if(primTitle.Contains("VertexerTracks") && vtx1->GetNContributors()>0) fNEntries->Fill(4);
 690
 691   // trigger class for PbPb C0SMH-B-NOPF-ALLNOTRD, C0SMH-B-NOPF-ALL
 692   //TString trigclass=aod->GetFiredTriggerClasses();
 693   //if(trigclass.Contains("C0SMH-B-NOPF-ALLNOTRD") || trigclass.Contains("C0SMH-B-NOPF-ALL")) fNEntries->Fill(5); //tmp
 694
 695   Int_t runNumber = aod->GetRunNumber();
 696
 697   Bool_t evSelbyCentrality=kTRUE,evSelected=kTRUE,evSelByVertex=kTRUE,evselByPileup=kFALSE;
 698   //select event
 699   if(!fCuts->IsEventSelected(aod)) {
 700     evSelected=kFALSE;
 701     if(fCuts->GetWhyRejection()==1) {fNEntries->Fill(1); evselByPileup=kTRUE;}// rejected for pileup
 702     if(fCuts->GetWhyRejection()==2 || fCuts->GetWhyRejection()==3) evSelbyCentrality=kFALSE; //rejected by centrality
 703     if(fCuts->GetWhyRejection()==4) evSelByVertex=kFALSE; //rejected by vertex
 704     if(fCuts->GetWhyRejection()==5) fNEntries->Fill(5);//tmp
 705     if(fCuts->GetWhyRejection()==6 && fOnOff[3]) ((AliCounterCollection*)fOutputEvSelection->FindObject("evselection"))->Count("evnonsel:zvtx/Run:%d",runNumber);
 706   }
 707
 708   if(evSelected || (!evSelected && !evSelbyCentrality && evSelByVertex && !evselByPileup)){ //events selected or not selected because of vtx or pileup
 709     if(fOnOff[2] && fCuts->GetUseCentrality()){
 710
 711       Float_t stdCentf=fCuts->GetCentrality(aod);
 712       Int_t stdCent = (Int_t)(stdCentf+0.5);
 713       Float_t secondCentf =fCuts->GetCentrality(aod,fEstimator);
 714       Int_t secondCent = (Int_t)(secondCentf+0.5);
 715       Int_t mincent=stdCent-stdCent%10;
 716       if(stdCentf==-1) {
 717         mincent=-10;
 718         stdCent=-1;
 719       }
 720       ((AliCounterCollection*)fOutputCounters->FindObject("stdEstimator"))->Count(Form("centralityclass:%d_%d/Run:%d",mincent,mincent+10,runNumber));
 721
 722       mincent=secondCent-secondCent%10;
 723       if(secondCentf==-1) {
 724         mincent=-10;
 725         secondCent=-1;
 726       }
 727       ((AliCounterCollection*)fOutputCounters->FindObject("secondEstimator"))->Count(Form("centralityclass:%d_%d/Run:%d",mincent,mincent+10,runNumber));
 728
 729       if(stdCent<fCuts->GetMinCentrality() || stdCent>fCuts->GetMaxCentrality()){
 730         ((TH1F*)fOutputCheckCentrality->FindObject("hNtrackletsOut"))->Fill(aod->GetTracklets()->GetNumberOfTracklets());
 731         ((TH1F*)fOutputCheckCentrality->FindObject("hMultOut"))->Fill(aod->GetHeader()->GetRefMultiplicity());
 732       }else{
 733         ((TH1F*)fOutputCheckCentrality->FindObject("hNtrackletsIn"))->Fill(aod->GetTracklets()->GetNumberOfTracklets());
 734         ((TH1F*)fOutputCheckCentrality->FindObject("hMultIn"))->Fill(aod->GetHeader()->GetRefMultiplicity());
 735       }
 736       ((TH2F*)fOutputCheckCentrality->FindObject("hMultvsPercentile"))->Fill(aod->GetHeader()->GetRefMultiplicity(),stdCentf);
 737
 738       PostData(6,fOutputCheckCentrality);
 739
 740     } else{
 741       if(fOnOff[0]){
 742         ((TH1F*)fOutputTrack->FindObject("hNtracklets"))->Fill(aod->GetTracklets()->GetNumberOfTracklets());
 743         ((TH1F*)fOutputTrack->FindObject("hMult"))->Fill(aod->GetHeader()->GetRefMultiplicity());
 744       }
 745     }
 746   }
 747
 748   if(fOnOff[3]){
 749     const AliVVertex *vertex = aod->GetPrimaryVertex();
 750     Double_t zvtx=vertex->GetZ();
 751     if(evSelected || (!evSelected && zvtx > 10.))
 752     ((TH1F*)fOutputEvSelection->FindObject("hzvtx"))->Fill(zvtx);
 753   }
 754
 755   if(!evSelected) {
 756     delete [] pdgdaughters;
 757     return; //discard all events not selected (vtx and/or centrality)
 758   }
 759
 760
 761   AliAODPidHF* pidHF=fCuts->GetPidHF();
 762   AliTPCPIDResponse* tpcres=new AliTPCPIDResponse();
 763   if(pidHF) pidHF->SetBetheBloch(*tpcres);
 764
 765
 766   Int_t ntracks=0;
 767   Int_t isGoodTrack=0, isFakeTrack=0;
 768
 769   if(aod) ntracks=aod->GetNTracks();
 770
 771   if(fOnOff[0] || fOnOff[1]){
 772     //loop on tracks in the event
 773     for (Int_t k=0;k<ntracks;k++){
 774       AliAODTrack* track=aod->GetTrack(k);
 775       AliAODPid *pid = track->GetDetPid();
 776
 777
 778       if(fOnOff[1]){
 779         if(!pid)  {if (fDebug>1)cout<<"No AliAODPid found"<<endl; continue;}
 780         Double_t times[AliPID::kSPECIES];
 781         pid->GetIntegratedTimes(times);
 782
 783         Double_t tofRes[AliPID::kSPECIES];
 784         pid->GetTOFpidResolution(tofRes);
 785
 786         //check TOF
 787         if(pidHF && pidHF->CheckStatus(track,"TOF")){
 788           ((TH1F*)fOutputPID->FindObject("hTOFtime"))->Fill(times[AliPID::kProton]);
 789           ((TH2F*)fOutputPID->FindObject("hTOFtimeKaonHyptime"))->Fill(track->P(),pid->GetTOFsignal()-times[3]); //3 is kaon
 790           ((TH1F*)fOutputPID->FindObject("hTOFsig"))->Fill(pid->GetTOFsignal());
 791           if (pid->GetTOFsignal()< 0) ((TH1F*)fOutputPID->FindObject("hTOFsig"))->Fill(-1);
 792
 793
 794           // test TOF sigma PID
 795           if (tofRes[2] != 0.) {   // protection against 'old' AODs...
 796             ((TH2F*)fOutputPID->FindObject("hTOFsigmaKSigPid"))->Fill(track->P(),(pid->GetTOFsignal()-times[AliPID::kKaon])/tofRes[3]);
 797             ((TH2F*)fOutputPID->FindObject("hTOFsigmaPionSigPid"))->Fill(track->P(),(pid->GetTOFsignal()-times[AliPID::kPion])/tofRes[2]);
 798             ((TH2F*)fOutputPID->FindObject("hTOFsigmaProtonSigPid"))->Fill(track->P(),(pid->GetTOFsignal()-times[AliPID::kProton])/tofRes[4]);
 799             for (Int_t iS=2; iS<5; iS++){ //we plot TOF Pid resolution for 3-sigma identified particles
 800               if ( (TMath::Abs(times[iS]-pid->GetTOFsignal())/tofRes[iS])<3.){
 801                 switch (iS) {
 802                 case AliPID::kPion:
 803                   ((TH1F*)fOutputPID->FindObject("hTOFsigPid3sigPion"))->Fill(tofRes[iS]);
 804                   break;
 805                 case AliPID::kKaon:
 806                   ((TH1F*)fOutputPID->FindObject("hTOFsigPid3sigKaon"))->Fill(tofRes[iS]);
 807                   break;
 808                 case AliPID::kProton:
 809                   ((TH1F*)fOutputPID->FindObject("hTOFsigPid3sigProton"))->Fill(tofRes[iS]);
 810                   break;
 811                 default:
 812                   break;
 813                 }
 814               }
 815             }
 816           }
 817         }//if TOF status
 818
 819         if(pidHF && pidHF->CheckStatus(track,"TPC")){
 820
 821           Double_t TPCp=pid->GetTPCmomentum();
 822           Double_t TPCsignal=pid->GetTPCsignal();
 823           ((TH1F*)fOutputPID->FindObject("hTPCsig"))->Fill(TPCsignal);
 824           ((TH1F*)fOutputPID->FindObject("hTPCsigvsp"))->Fill(TPCp,TPCsignal);
 825           //if (pidHF->IsKaonRaw(track, "TOF"))
 826           ((TH2F*)fOutputPID->FindObject("hTPCsigmaK"))->Fill(TPCp,tpcres->GetNumberOfSigmas(TPCp,TPCsignal,track->GetTPCNcls(),AliPID::kKaon));
 827           //if (pidHF->IsPionRaw(track, "TOF"))
 828           ((TH2F*)fOutputPID->FindObject("hTPCsigmaPion"))->Fill(TPCp,tpcres->GetNumberOfSigmas(TPCp,TPCsignal,track->GetTPCNcls(),AliPID::kPion));
 829           //if (pidHF->IsProtonRaw(track,"TOF"))
 830           ((TH2F*)fOutputPID->FindObject("hTPCsigmaProton"))->Fill(TPCp,tpcres->GetNumberOfSigmas(TPCp,TPCsignal,track->GetTPCNcls(),AliPID::kProton));
 831         }//if TPC status
 832       } //end PID histograms
 833
 834       Int_t nclsTot=0,nclsSPD=0;
 835
 836       //check clusters of the tracks
 837       if(fOnOff[0]){
 838
 839         ((TH1F*)fOutputTrack->FindObject("hnLayerITS"))->Fill(-1);
 840         for(Int_t l=0;l<6;l++) {
 841           if(TESTBIT(track->GetITSClusterMap(),l)) {
 842             ((TH1F*)fOutputTrack->FindObject("hnLayerITS"))->Fill(l);
 843             nclsTot++; if(l<2) nclsSPD++;
 844           }
 845         }
 846         ((TH1F*)fOutputTrack->FindObject("hnClsITS"))->Fill(nclsTot);
 847         ((TH1F*)fOutputTrack->FindObject("hnClsSPD"))->Fill(nclsSPD);
 848
 849         if(!(track->GetStatus()&AliESDtrack::kTPCin) && track->GetStatus()&AliESDtrack::kITSrefit && !(track->GetStatus()&AliESDtrack::kITSpureSA)){//tracks retrieved in the ITS and not reconstructed in the TPC
 850           ((TH1F*)fOutputTrack->FindObject("hnClsITS-SA"))->Fill(nclsTot);
 851           ((TH1F*)fOutputTrack->FindObject("hnLayerITS"))->Fill(-1);
 852           for(Int_t l=0;l<6;l++) {
 853             if(TESTBIT(track->GetITSClusterMap(),l)) {
 854               ((TH1F*)fOutputTrack->FindObject("hnLayerITSsa"))->Fill(l);
 855             }
 856           }
 857         }
 858         Int_t label=0;
 859         if(fReadMC){
 860           label=track->GetLabel();
 861           if (label<0)fNEntries->Fill(8);
 862           else fNEntries->Fill(9);
 863         }
 864
 865         if(isSimpleMode){
 866
 867           if (track->Pt()>0.3 &&
 868               track->GetStatus()&AliESDtrack::kTPCrefit &&
 869               track->GetStatus()&AliESDtrack::kITSrefit &&
 870               /*nclsTot>3 &&*/
 871               nclsSPD>0) {//count good tracks
 872
 873
 874             if(fReadMC && label<0) {
 875               ((TH1F*)fOutputTrack->FindObject("hptFakeTr"))->Fill(track->Pt());
 876               isFakeTrack++;
 877             } else {
 878               ((TH1F*)fOutputTrack->FindObject("hptGoodTr"))->Fill(track->Pt());
 879               isGoodTrack++;
 880             }
 881           }
 882         } //simple mode: no IsSelected on tracks: use "manual" cuts
 883       } //fill track histos
 884     } //end loop on tracks
 885
 886       //fill once per event
 887     if(fOnOff[0]){
 888       if (fReadMC) ((TH1F*)fOutputTrack->FindObject("hdistrFakeTr"))->Fill(isFakeTrack);
 889       ((TH1F*)fOutputTrack->FindObject("hdistrGoodTr"))->Fill(isGoodTrack);
 890     }
 891
 892     if(!isSimpleMode){
 893       // loop over candidates
 894       Int_t nCand = arrayProng->GetEntriesFast();
 895       Int_t ndaugh=3;
 896       if(fDecayChannel==AliAnalysisTaskSEHFQA::kD0toKpi) ndaugh=2;
 897       if(fDecayChannel==AliAnalysisTaskSEHFQA::kD0toKpipipi) ndaugh=4;
 898
 899       for (Int_t iCand = 0; iCand < nCand; iCand++) {
 900         AliAODRecoDecayHF *d = (AliAODRecoDecayHF*)arrayProng->UncheckedAt(iCand);
 901         if(d->GetSelectionMap()) {
 902           if(fDecayChannel==AliAnalysisTaskSEHFQA::kD0toKpi && !d->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts)) continue; //skip the D0 from Dstar
 903           if(fDecayChannel==AliAnalysisTaskSEHFQA::kDplustoKpipi && !d->HasSelectionBit(AliRDHFCuts::kDplusCuts)) continue; //skip the 3 prong !D+
 904         }
 905
 906         if(fReadMC){
 907           Int_t labD = d->MatchToMC(pdg,mcArray,ndaugh,pdgdaughters);
 908           if(labD>=0){
 909             AliAODMCParticle *partD = (AliAODMCParticle*)mcArray->At(labD);
 910             Int_t label=partD->GetMother();
 911             AliAODMCParticle *mot = (AliAODMCParticle*)mcArray->At(label);
 912             while(label>=0){//get first mother
 913               mot = (AliAODMCParticle*)mcArray->At(label);
 914               label=mot->GetMother();
 915             }
 916             Int_t pdgMotCode = mot->GetPdgCode();
 917
 918             if(TMath::Abs(pdgMotCode)==4) fNEntries->Fill(6); //from primary charm
 919             if(TMath::Abs(pdgMotCode)==5) fNEntries->Fill(7); //from beauty
 920
 921           }
 922         }//end MC
 923         else fNEntries->Fill(6); //count the candidates (data)
 924
 925         for(Int_t id=0;id<ndaugh;id++){
 926
 927           //other histograms to be filled when the cut object is given
 928           AliAODTrack* track=(AliAODTrack*)d->GetDaughter(id);
 929
 930           //track quality
 931
 932           if (fCuts->IsInFiducialAcceptance(d->Pt(),d->Y(pdg)) && fCuts->IsSelected(d,AliRDHFCuts::kTracks,aod)) {
 933
 934             Int_t label=0;
 935             if(fReadMC)label=track->GetLabel();
 936             if(fOnOff[0]){
 937
 938               if(fReadMC && label<0) {
 939                 isFakeTrack++;
 940                 ((TH1F*)fOutputTrack->FindObject("hptFakeTr"))->Fill(track->Pt());
 941
 942                 ((TH1F*)fOutputTrack->FindObject("hd0f"))->Fill(d->Getd0Prong(id));
 943               } else {
 944                 isGoodTrack++;
 945                 ((TH1F*)fOutputTrack->FindObject("hptGoodTr"))->Fill(track->Pt());
 946
 947                 ((TH1F*)fOutputTrack->FindObject("hd0"))->Fill(d->Getd0Prong(id));
 948               }
 949             }
 950             if (fCuts->IsSelected(d,AliRDHFCuts::kAll,aod) && fOnOff[1]){
 951
 952               AliAODPid *pid = track->GetDetPid();
 953               Double_t times[5];
 954               pid->GetIntegratedTimes(times);
 955               if(pidHF && pidHF->CheckStatus(track,"TOF")) ((TH2F*)fOutputPID->FindObject("hTOFtimeKaonHyptimeAC"))->Fill(track->P(),pid->GetTOFsignal()-times[AliPID::kKaon]);
 956               if(pidHF && pidHF->CheckStatus(track,"TPC")) ((TH2F*)fOutputPID->FindObject("hTPCsigvspAC"))->Fill(pid->GetTPCmomentum(),pid->GetTPCsignal());
 957
 958               fNEntries->Fill(3);
 959             } //end analysis cuts
 960           } //end acceptance and track cuts
 961         } //end loop on tracks in the candidate
 962       } //end loop on candidates
 963       if(fOnOff[0]){
 964         if(fReadMC) ((TH1F*)fOutputTrack->FindObject("hdistrFakeTr"))->Fill(isFakeTrack);
 965         ((TH1F*)fOutputTrack->FindObject("hdistrGoodTr"))->Fill(isGoodTrack);
 966       }
 967     }
 968   } //end if on pid or track histograms
 969
 970   delete tpcres;
 971   delete [] pdgdaughters;
 972   PostData(1,fNEntries);
 973   if(fOnOff[1]) PostData(2,fOutputPID);
 974   if(fOnOff[0]) PostData(3,fOutputTrack);
 975   PostData(4,fCuts);
 976   if(fOnOff[2]) PostData(5,fOutputCounters);
 977   //Post data 6 done in case of centrality on
 978
 979 }
 980
 981 //____________________________________________________________________________
 982 void AliAnalysisTaskSEHFQA::Terminate(Option_t */*option*/){
 983   //terminate analysis
 984
 985   fNEntries = dynamic_cast<TH1F*>(GetOutputData(1));
 986   if(!fNEntries){
 987     printf("ERROR: %s not available\n",GetOutputSlot(1)->GetContainer()->GetName());
 988     return;
 989   }
 990
 991   fOutputPID = dynamic_cast<TList*> (GetOutputData(2));
 992   if (!fOutputPID && fOnOff[1]) {
 993     printf("ERROR: %s not available\n",GetOutputSlot(2)->GetContainer()->GetName());
 994     return;
 995   }
 996
 997   fOutputTrack = dynamic_cast<TList*> (GetOutputData(3));
 998   if (!fOutputTrack && fOnOff[0]) {
 999     printf("ERROR: %s not available\n",GetOutputSlot(3)->GetContainer()->GetName());
1000     return;
1001   }
1002
1003 }
1004