PWGHF/hfe/AliAnalysisTaskHFEemcQA.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, 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
  17 //////////////////////////////////////////////
  18 //    QA task for EMCAL electron analysis   //
  19 //  Author: Deepa Thomas, Shingo Sakai      //
  20 //////////////////////////////////////////////
  21
  22 #include "TChain.h"
  23 #include "TTree.h"
  24 #include "TH1F.h"
  25 #include "TCanvas.h"
  26 #include "THnSparse.h"
  27
  28 #include "AliAnalysisTask.h"
  29 #include "AliAnalysisManager.h"
  30
  31 #include "AliESDEvent.h"
  32 #include "AliESDInputHandler.h"
  33 #include "AliESDtrackCuts.h"
  34 #include "AliAODEvent.h"
  35 #include "AliAODHandler.h"
  36
  37 #include "AliPID.h"
  38 #include "AliESDpid.h"
  39 #include "AliAODPid.h"
  40 #include "AliPIDResponse.h"
  41 #include "AliHFEcontainer.h"
  42 #include "AliHFEcuts.h"
  43 #include "AliHFEpid.h"
  44 #include "AliHFEpidBase.h"
  45 #include "AliHFEpidQAmanager.h"
  46 #include "AliHFEtools.h"
  47 #include "AliCFContainer.h"
  48 #include "AliCFManager.h"
  49
  50 #include "AliAnalysisTaskHFEemcQA.h"
  51
  52 using std::cout;
  53 using std::endl;
  54
  55 ClassImp(AliAnalysisTaskHFEemcQA)
  56   //________________________________________________________________________
  57   AliAnalysisTaskHFEemcQA::AliAnalysisTaskHFEemcQA(const char *name)
  58 : AliAnalysisTaskSE(name),
  59   fVevent(0),
  60   fESD(0),
  61   fAOD(0),
  62   fpidResponse(0),
  63   fFlagSparse(kFALSE),
  64   fUseTender(kTRUE),
  65   fTracks_tender(0),
  66   fCaloClusters_tender(0),
  67   fOutputList(0),
  68   fNevents(0),
  69   fVtxZ(0),
  70   fVtxX(0),
  71   fVtxY(0),
  72   fTrigMulti(0),
  73   fHistClustE(0),
  74   fEMCClsEtaPhi(0),
  75   fHistoNCls(0),
  76   fHistoNClsE1(0),
  77   fHistoNClsE2(0),
  78   fHistoNClsE3(0),
  79   fHistoNCells(0),
  80   fHistoCalCell(0),
  81   fNegTrkIDPt(0),
  82   fTrkPt(0),
  83   fTrketa(0),
  84   fTrkphi(0),
  85   fdEdx(0),
  86   fTPCNpts(0),
  87   fTPCnsig(0),
  88   fHistPtMatch(0),
  89   fEMCTrkMatch(0),
  90   fEMCTrkPt(0),
  91   fEMCTrketa(0),
  92   fEMCTrkphi(0),
  93   fEMCdEdx(0),
  94   fEMCTPCnsig(0),
  95   fEMCTPCNpts(0),
  96   fClsEAftMatch(0),
  97   fHistdEdxEop(0),
  98   fHistNsigEop(0),
  99   fHistEop(0),
 100   fM20(0),
 101   fM02(0),
 102   fM20EovP(0),
 103   fM02EovP(0),
 104   fEleCanTPCNpts(0),
 105   fEleCanTPCNCls(0),
 106   fEleCanITSNCls(0),
 107   fEleCanITShit(0),
 108   fEleCanSPD1(0),
 109   fEleCanSPD2(0),
 110   fEleCanSPDBoth(0),
 111   fEleCanSPDOr(0),
 112   fSparseElectron(0),
 113   fvalueElectron(0)
 114 {
 115   // Constructor
 116
 117   fvalueElectron = new Double_t[6];
 118   // Define input and output slots here
 119   // Input slot #0 works with a TChain
 120   DefineInput(0, TChain::Class());
 121   // Output slot #0 id reserved by the base class for AOD
 122   // Output slot #1 writes into a TH1 container
 123   DefineOutput(1, TList::Class());
 124 }
 125 //________________________________________________________________________
 126 AliAnalysisTaskHFEemcQA::AliAnalysisTaskHFEemcQA()
 127   : AliAnalysisTaskSE("DefaultTask_HfeEMCQA"),
 128   fVevent(0),
 129   fESD(0),
 130   fAOD(0),
 131   fpidResponse(0),
 132   fFlagSparse(kFALSE),
 133   fUseTender(kTRUE),
 134   fTracks_tender(0),
 135   fCaloClusters_tender(0),
 136   fOutputList(0),
 137   fNevents(0),
 138   fVtxZ(0),
 139   fVtxX(0),
 140   fVtxY(0),
 141   fTrigMulti(0),
 142   fHistClustE(0),
 143   fEMCClsEtaPhi(0),
 144   fHistoNCls(0),
 145   fHistoNClsE1(0),
 146   fHistoNClsE2(0),
 147   fHistoNClsE3(0),
 148   fHistoNCells(0),
 149   fHistoCalCell(0),
 150   fNegTrkIDPt(0),
 151   fTrkPt(0),
 152   fTrketa(0),
 153   fTrkphi(0),
 154   fdEdx(0),
 155   fTPCNpts(0),
 156   fTPCnsig(0),
 157   fHistPtMatch(0),
 158   fEMCTrkMatch(0),
 159   fEMCTrkPt(0),
 160   fEMCTrketa(0),
 161   fEMCTrkphi(0),
 162   fEMCdEdx(0),
 163   fEMCTPCnsig(0),
 164   fEMCTPCNpts(0),
 165   fClsEAftMatch(0),
 166   fHistdEdxEop(0),
 167   fHistNsigEop(0),
 168   fHistEop(0),
 169   fM20(0),
 170   fM02(0),
 171   fM20EovP(0),
 172   fM02EovP(0),
 173   fEleCanTPCNpts(0),
 174   fEleCanTPCNCls(0),
 175   fEleCanITSNCls(0),
 176   fEleCanITShit(0),
 177   fEleCanSPD1(0),
 178   fEleCanSPD2(0),
 179   fEleCanSPDBoth(0),
 180   fEleCanSPDOr(0),
 181   fSparseElectron(0),
 182   fvalueElectron(0)
 183 {
 184   //Default constructor
 185
 186   fvalueElectron = new Double_t[6];
 187   // Define input and output slots here
 188   // Input slot #0 works with a TChain
 189   DefineInput(0, TChain::Class());
 190   // Output slot #0 id reserved by the base class for AOD
 191   // Output slot #1 writes into a TH1 container
 192   // DefineOutput(1, TH1I::Class());
 193   DefineOutput(1, TList::Class());
 194   //DefineOutput(3, TTree::Class());
 195 }
 196 //________________________________________________________________________
 197 AliAnalysisTaskHFEemcQA::~AliAnalysisTaskHFEemcQA()
 198 {
 199   //Destructor
 200   delete fOutputList;
 201   delete fTracks_tender;
 202   delete fCaloClusters_tender;
 203   delete fSparseElectron;
 204   delete []fvalueElectron;
 205 }
 206 //________________________________________________________________________
 207 void AliAnalysisTaskHFEemcQA::UserCreateOutputObjects()
 208 {
 209   // Create histograms
 210   // Called once
 211   AliDebug(3, "Creating Output Objects");
 212
 213   /////////////////////////////////////////////////
 214   //Automatic determination of the analysis mode//
 215   ////////////////////////////////////////////////
 216   AliVEventHandler *inputHandler = dynamic_cast<AliVEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
 217   if(!TString(inputHandler->IsA()->GetName()).CompareTo("AliAODInputHandler")){
 218     SetAODAnalysis();
 219   } else {
 220     SetESDAnalysis();
 221   }
 222   printf("Analysis Mode: %s Analysis\n", IsAODanalysis() ? "AOD" : "ESD");
 223
 224   ////////////////
 225   //Output list//
 226   ///////////////
 227   fOutputList = new TList();
 228   fOutputList->SetOwner();
 229
 230   fNevents = new TH1F("fNevents","No of events",3,-0.5,2.5);
 231   fOutputList->Add(fNevents);
 232   fNevents->GetYaxis()->SetTitle("counts");
 233   fNevents->GetXaxis()->SetBinLabel(1,"All");
 234   fNevents->GetXaxis()->SetBinLabel(2,"With >2 Trks");
 235   fNevents->GetXaxis()->SetBinLabel(3,"Vtx_{z}<10cm");
 236
 237   fVtxZ = new TH1F("fVtxZ","Z vertex position;Vtx_{z};counts",1000,-50,50);
 238   fOutputList->Add(fVtxZ);
 239
 240   fVtxY = new TH1F("fVtxY","Y vertex position;Vtx_{y};counts",1000,-50,50);
 241   fOutputList->Add(fVtxY);
 242
 243   fVtxX = new TH1F("fVtxX","X vertex position;Vtx_{x};counts",1000,-50,50);
 244   fOutputList->Add(fVtxX);
 245
 246   fTrigMulti = new TH2F("fTrigMulti","Multiplicity distribution for different triggers; Trigger type; multiplicity",11,-1,10,2000,0,2000);
 247   fOutputList->Add(fTrigMulti);
 248
 249   fHistClustE = new TH1F("fHistClustE", "EMCAL cluster energy distribution; Cluster E;counts", 5000, 0.0, 50.0);
 250   fOutputList->Add(fHistClustE);
 251
 252   fEMCClsEtaPhi = new TH2F("fEMCClsEtaPhi","EMCAL cluster #eta and #phi distribution;#eta;#phi",100,-0.9,0.9,200,0,6.3);
 253   fOutputList->Add(fEMCClsEtaPhi);
 254
 255   fHistoNCls = new TH1F("fHistoNCls","No of EMCAL cluster in the event;N^{EMC}_{cls};counts",150,0,150);
 256   fOutputList->Add(fHistoNCls);
 257
 258   fHistoNClsE1 = new TH1F("fHistoNClsE1","No of EMCAL cluster in the event (E>0.1 GeV);N^{EMC}_{cls};counts",150,0,150);
 259   fOutputList->Add(fHistoNClsE1);
 260
 261   fHistoNClsE2 = new TH1F("fHistoNClsE2","No of EMCAL cluster in the event (E>0.2 GeV);N^{EMC}_{cls};counts",150,0,150);
 262   fOutputList->Add(fHistoNClsE2);
 263
 264   fHistoNClsE3 = new TH1F("fHistoNClsE3","No of EMCAL cluster in the event (E>0.5 GeV);N^{EMC}_{cls};counts",150,0,150);
 265   fOutputList->Add(fHistoNClsE3);
 266
 267   //fHistoNCells = new TH1F("fHistoNCells","No of EMCAL cells in a cluster;N^{EMC}_{cells};counts",30,0,30);
 268   fHistoNCells = new TH2F("fHistoNCells","No of EMCAL cells in a cluster;Cluster E;N^{EMC}_{cells}",300,0,30,30,0,30);
 269   fOutputList->Add(fHistoNCells);
 270
 271   fHistoCalCell = new TH2F("fHistoCalCell","EMCAL cells in a cluster;cell ID;E (GeV)",15000,-0.5,14999.5,600,0,30);
 272   fOutputList->Add(fHistoCalCell);
 273
 274   fNegTrkIDPt = new TH1F("fNegTrkIDPt", "p_{T} distribution of tracks with negative track id;p_{T} (GeV/c);counts", 500, 0.0, 50.0);
 275   fOutputList->Add(fNegTrkIDPt);
 276
 277   fTrkPt = new TH1F("fTrkPt","p_{T} distribution of all tracks;p_{T} (GeV/c);counts",1000,0,100);
 278   fOutputList->Add(fTrkPt);
 279
 280   fTrketa = new TH1F("fTrketa","All Track #eta distribution;#eta;counts",100,-1.5,1.5);
 281   fOutputList->Add(fTrketa);
 282
 283   fTrkphi = new TH1F("fTrkphi","All Track #phi distribution;#phi;counts",100,0,6.3);
 284   fOutputList->Add(fTrkphi);
 285
 286   fdEdx = new TH2F("fdEdx","All Track dE/dx distribution;p (GeV/c);dE/dx",200,0,20,500,0,160);
 287   fOutputList->Add(fdEdx);
 288
 289   fTPCNpts = new TH2F("fTPCNpts","All track TPC Npoints used for dE/dx calculation;p (GeV/c);N points",200,0,20,200,0.,200.);
 290   fOutputList->Add(fTPCNpts);
 291
 292   fTPCnsig = new TH2F("fTPCnsig","All Track TPC Nsigma distribution;p (GeV/c);#sigma_{TPC-dE/dx}",1000,0,50,200,-10,10);
 293   fOutputList->Add(fTPCnsig);
 294
 295   fHistPtMatch = new TH1F("fHistPtMatch", "p_{T} distribution of tracks matched to EMCAL;p_{T} (GeV/c);counts",1000, 0.0, 100.0);
 296   fOutputList->Add(fHistPtMatch);
 297
 298   fEMCTrkMatch = new TH2F("fEMCTrkMatch","Distance of EMCAL cluster to its closest track;#phi;z",100,-0.3,0.3,100,-0.3,0.3);
 299   fOutputList->Add(fEMCTrkMatch);
 300
 301   fEMCTrkPt = new TH1F("fEMCTrkPt","p_{T} distribution of tracks with EMCAL cluster;p_{T} (GeV/c);counts",1000,0,100);
 302   fOutputList->Add(fEMCTrkPt);
 303
 304   fEMCTrketa = new TH1F("fEMCTrketa","#eta distribution of tracks matched to EMCAL;#eta;counts",100,-1.5,1.5);
 305   fOutputList->Add(fEMCTrketa);
 306
 307   fEMCTrkphi = new TH1F("fEMCTrkphi","#phi distribution of tracks matched to EMCAL;#phi;counts",100,0,6.3);
 308   fOutputList->Add(fEMCTrkphi);
 309
 310   fEMCdEdx = new TH2F("fEMCdEdx","dE/dx distribution of tracks matched to EMCAL;p (GeV/c);dE/dx",200,0,20,500,0,160);
 311   fOutputList->Add(fEMCdEdx);
 312
 313   fEMCTPCnsig = new TH2F("fEMCTPCnsig","TPC Nsigma distribution of tracks matched to EMCAL;p (GeV/c);#sigma_{TPC-dE/dx}",1000,0,50,200,-10,10);
 314   fOutputList->Add(fEMCTPCnsig);
 315
 316   fEMCTPCNpts = new TH2F("fEMCTPCNpts","TPC Npoints used for dE/dx for tracks matched to EMCAL;p (GeV/c);N points",200,0,20,200,0.,200.);
 317   fOutputList->Add(fEMCTPCNpts);
 318
 319   fClsEAftMatch = new TH1F("fClsEAftMatch", "EMCAL cluster energy distribution after track matching; Cluster E;counts", 500, 0.0, 50.0);
 320   fOutputList->Add(fClsEAftMatch);
 321
 322   fHistEop = new TH2F("fHistEop", "E/p distribution;p_{T} (GeV/c);E/p", 200,0,20,60, 0.0, 3.0);
 323   fOutputList->Add(fHistEop);
 324
 325   fHistdEdxEop = new TH2F("fHistdEdxEop", "E/p vs dE/dx;E/p;dE/dx", 60, 0.0, 3.0, 500,0,160);
 326   fOutputList->Add(fHistdEdxEop);
 327
 328   fHistNsigEop = new TH2F ("fHistNsigEop", "E/p vs TPC nsig",60, 0.0, 3.0, 200, -10,10);
 329   fOutputList->Add(fHistNsigEop);
 330
 331   fM20 = new TH2F ("fM20","M20 vs pt distribution",200,0,20,400,0,2);
 332   fOutputList->Add(fM20);
 333
 334   fM02 = new TH2F ("fM02","M02 vs pt distribution",200,0,20,400,0,2);
 335   fOutputList->Add(fM02);
 336
 337   fM20EovP = new TH2F ("fM20EovP","M20 vs E/p distribution",400,0,3,400,0,2);
 338   fOutputList->Add(fM20EovP);
 339
 340   fM02EovP = new TH2F ("fM02EovP","M02 vs E/p distribution",400,0,3,400,0,2);
 341   fOutputList->Add(fM02EovP);
 342
 343   fEleCanTPCNpts = new TH2F("fEleCanTPCNpts","TPC Npoints used for dE/dx for electron candidates;p_{T} (GeV/c);N points",200,0,20,200,0,200);
 344   fOutputList->Add(fEleCanTPCNpts);
 345
 346   fEleCanTPCNCls = new TH2F("fEleCanTPCNCls","TPC N clusters for electron candidates;p_{T} (GeV/c);N TPC clusters",200,0,20,171,-0.5,170.5);
 347   fOutputList->Add(fEleCanTPCNCls);
 348
 349   fEleCanITSNCls = new TH2F("fEleCanITSNCls","ITS N clusters for electron candidates;p_{T} (GeV/c);N ITS clusters",200,0,20,8,-0.5,7.5);
 350   fOutputList->Add(fEleCanITSNCls);
 351
 352   fEleCanITShit = new TH1F("fEleCanITShit","ITS hit map;ITS layer;counts",7,-0.5,6.5);
 353   fOutputList->Add(fEleCanITShit);
 354
 355   fEleCanSPD1 = new TH2F("fEleCanSPD1","Hit on SPD layer 1;p_{T} (GeV/c);Hit",200,0,20,1,0,1);
 356   fOutputList->Add(fEleCanSPD1);
 357
 358   fEleCanSPD2 = new TH2F("fEleCanSPD2","Hit on SPD layer 2;p_{T} (GeV/c);Hit",200,0,20,1,0,1);
 359   fOutputList->Add(fEleCanSPD2);
 360
 361   fEleCanSPDBoth = new TH2F("fEleCanSPDBoth","Tracks with hits on both SPD layer;p_{T} (GeV/c);Hit",200,0,20,1,0,1);
 362   fOutputList->Add(fEleCanSPDBoth);
 363
 364   fEleCanSPDOr = new TH2F("fEleCanSPDOr","Tracks with hits on both SPD layer;p_{T} (GeV/c);Hit",200,0,20,1,0,1);
 365   fOutputList->Add(fEleCanSPDOr);
 366
 367   Int_t bins[6]={8,500,200,400,400,400}; //trigger, pt, TPCnsig, E/p, M20, M02
 368   Double_t xmin[6]={-0.5,0,-10,0,0,0};
 369   Double_t xmax[6]={7.5,25,10,2,2,2};
 370   fSparseElectron = new THnSparseD ("Electron","Electron",6,bins,xmin,xmax);
 371   fOutputList->Add(fSparseElectron);
 372
 373   PostData(1,fOutputList);
 374 }
 375
 376 //________________________________________________________________________
 377 void AliAnalysisTaskHFEemcQA::UserExec(Option_t *)
 378 {
 379   // Main loop
 380   // Called for each event
 381   // Post output data.
 382
 383   UInt_t evSelMask=((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
 384
 385   fVevent = dynamic_cast<AliVEvent*>(InputEvent());
 386   if (!fVevent) {
 387     printf("ERROR: fVEvent not available\n");
 388     return;
 389   }
 390
 391   fESD = dynamic_cast<AliESDEvent*>(InputEvent());
 392   if (fESD) {
 393     //   printf("fESD available\n");
 394     //return;
 395   }
 396
 397   //////////////
 398   //if Tender //
 399   //////////////
 400   if(fUseTender){
 401     //new branches with calibrated tracks and clusters
 402     if(IsAODanalysis()) fTracks_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("AODFilterTracks"));
 403     if(!IsAODanalysis()) fTracks_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("ESDFilterTracks"));
 404
 405     fCaloClusters_tender = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("EmcCaloClusters"));
 406   }
 407
 408   ////////////////////
 409   //cuts initialised//
 410   ///////////////////
 411   AliESDtrackCuts* esdTrackCutsH = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(kFALSE);
 412   esdTrackCutsH->SetMaxDCAToVertexXY(2.4);
 413   esdTrackCutsH->SetMaxDCAToVertexZ(3.2);
 414   esdTrackCutsH->SetDCAToVertex2D(kTRUE);
 415
 416   fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
 417   if (fAOD) {
 418     // printf("fAOD available\n");
 419     //return;
 420   }
 421
 422   ///////////////////
 423   //PID initialised//
 424   //////////////////
 425   fpidResponse = fInputHandler->GetPIDResponse();
 426
 427   ////////////////
 428   //Event vertex//
 429   ///////////////
 430   Int_t ntracks = -999;
 431   if(!fUseTender)ntracks = fVevent->GetNumberOfTracks();
 432   if(fUseTender) ntracks = fTracks_tender->GetEntries();
 433   if(ntracks < 1) printf("There are %d tracks in this event\n",ntracks);
 434
 435   fNevents->Fill(0); //all events
 436   Double_t Zvertex = -100, Xvertex = -100, Yvertex = -100;
 437   const AliVVertex *pVtx = fVevent->GetPrimaryVertex();
 438   Double_t NcontV = pVtx->GetNContributors();
 439   if(NcontV<2)return;
 440   fNevents->Fill(1); //events with 2 tracks
 441
 442   Zvertex = pVtx->GetZ();
 443   Yvertex = pVtx->GetY();
 444   Xvertex = pVtx->GetX();
 445   fVtxZ->Fill(Zvertex);
 446   fVtxX->Fill(Xvertex);
 447   fVtxY->Fill(Yvertex);
 448
 449   /////////////////
 450   //trigger check//
 451   /////////////////
 452   TString firedTrigger;
 453   TString TriggerEG1("EG1");
 454   TString TriggerEG2("EG2");
 455   fVevent->GetFiredTriggerClasses();
 456
 457   Bool_t EG1tr = kFALSE;
 458   Bool_t EG2tr = kFALSE;
 459
 460   if(firedTrigger.Contains(TriggerEG1))EG1tr = kTRUE;
 461   if(firedTrigger.Contains(TriggerEG2))EG2tr = kTRUE;
 462
 463   Int_t trigger = -1;
 464   if (fAOD){
 465     //Double_t multiplicity=fAOD->GetHeader()->GetRefMultiplicity();
 466     AliAODHeader *header = dynamic_cast<AliAODHeader*>(fAOD->GetHeader());
 467     if(!header) AliFatal("Not a standard AOD");
 468     Double_t multiplicity = header->GetRefMultiplicity();
 469
 470     fTrigMulti->Fill(-0.5, multiplicity);
 471     if(evSelMask & AliVEvent::kAny) fTrigMulti->Fill(0.5, multiplicity);
 472     if(evSelMask & AliVEvent::kMB) fTrigMulti->Fill(1.5, multiplicity);
 473     if(evSelMask & AliVEvent::kINT7) fTrigMulti->Fill(2.5, multiplicity);
 474     if(evSelMask & AliVEvent::kINT8) fTrigMulti->Fill(3.5, multiplicity);
 475     if(evSelMask & AliVEvent::kEMC1) fTrigMulti->Fill(4.5, multiplicity);
 476     if(evSelMask & AliVEvent::kEMC7) fTrigMulti->Fill(5.5, multiplicity);
 477     if(evSelMask & AliVEvent::kEMC8) fTrigMulti->Fill(6.5, multiplicity);
 478     if(evSelMask & AliVEvent::kEMCEJE) fTrigMulti->Fill(7.5, multiplicity);
 479     if(evSelMask & AliVEvent::kEMCEGA) fTrigMulti->Fill(8.5, multiplicity);
 480     if(evSelMask & AliVEvent::kEMCEGA & EG2tr) fTrigMulti->Fill(9.5, multiplicity);
 481
 482     if(evSelMask & AliVEvent::kMB) trigger =0;
 483     if(evSelMask & AliVEvent::kINT7) trigger =1;
 484     if(evSelMask & AliVEvent::kINT8) trigger =2;
 485     if(evSelMask & AliVEvent::kEMC1) trigger =3;
 486     if(evSelMask & AliVEvent::kEMC7) trigger =4;
 487     if(evSelMask & AliVEvent::kEMC8) trigger =5;
 488     if(evSelMask & AliVEvent::kEMCEJE) trigger =6;
 489     if(evSelMask & AliVEvent::kEMCEGA) trigger =7;
 490   }
 491
 492   ////////////////////
 493   //event selection//
 494   ///////////////////
 495   if(fabs(Zvertex>10.0))return;
 496   fNevents->Fill(2); //events after z vtx cut
 497
 498   /////////////////////////////
 499   //EMCAL cluster information//
 500   ////////////////////////////
 501   Int_t Nclust = -999;
 502   if(!fUseTender) Nclust = fVevent->GetNumberOfCaloClusters();
 503   if(fUseTender) Nclust = fCaloClusters_tender->GetEntries();
 504
 505   int NclustAll= 0;
 506   int NclustE1 = 0; //# of clust E>0.1
 507   int NclustE2 = 0; //# of clust E>0.2
 508   int NclustE3 = 0; //# of clust E>0.5
 509
 510   for(Int_t icl=0; icl<Nclust; icl++)
 511   {
 512     AliVCluster *clust = 0x0;
 513     if(!fUseTender) clust = fVevent->GetCaloCluster(icl);
 514     if(fUseTender) clust = dynamic_cast<AliVCluster*>(fCaloClusters_tender->At(icl));
 515     if(!clust)  printf("ERROR: Could not receive cluster matched calibrated from track %d\n", icl);
 516
 517     if(clust && clust->IsEMCAL())
 518     {
 519       Double_t clustE = clust->E();
 520       Float_t  emcx[3]; // cluster pos
 521       clust->GetPosition(emcx);
 522       TVector3 clustpos(emcx[0],emcx[1],emcx[2]);
 523       Double_t emcphi = clustpos.Phi();
 524       Double_t emceta = clustpos.Eta();
 525       fHistClustE->Fill(clustE);
 526       fEMCClsEtaPhi->Fill(emceta,emcphi);
 527       fHistoNCells->Fill(clustE,clust->GetNCells());
 528       //fHistoNCells->Fill(clust->GetNCells());
 529
 530       NclustAll++;
 531       if(clustE>0.1)NclustE1++;
 532       if(clustE>0.2)NclustE2++;
 533       if(clustE>0.5)NclustE3++;
 534     }
 535   }
 536
 537   fHistoNCls->Fill(NclustAll);
 538   fHistoNClsE1->Fill(NclustE1);
 539   fHistoNClsE2->Fill(NclustE2);
 540   fHistoNClsE3->Fill(NclustE3);
 541
 542   // cell information
 543   AliVCaloCells *fCaloCells = fVevent->GetEMCALCells();
 544
 545   //Int_t nSACell, iSACell, mclabel;
 546   Short_t cellAddr, nSACell;
 547   Int_t  mclabel;
 548   Short_t iSACell;
 549   Double_t cellAmp=-1., cellTimeT=-1., clusterTime=-1., efrac=-1.;
 550
 551   nSACell = fCaloCells->GetNumberOfCells();
 552   for(iSACell = 0; iSACell < nSACell; iSACell++ ){
 553     Bool_t haveCell = fCaloCells->GetCell(iSACell, cellAddr, cellAmp, cellTimeT , mclabel, efrac);
 554     //virtual Bool_t   GetCell(Short_t pos, Short_t &cellNumber, Double_t &amplitude, Double_t &time, Int_t &mclabel,    Double_t  &efrac)
 555     if(haveCell)fHistoCalCell->Fill(cellAddr,cellAmp);
 556
 557   }
 558
 559   /////////////////////////////////
 560   //Look for kink mother for AOD//
 561   /////////////////////////////////
 562   Int_t numberofvertices = 100;
 563   if(fAOD) numberofvertices = fAOD->GetNumberOfVertices();
 564   Double_t listofmotherkink[numberofvertices];
 565   Int_t numberofmotherkink = 0;
 566   if(IsAODanalysis())
 567   {
 568     for(Int_t ivertex=0; ivertex < numberofvertices; ivertex++) {
 569       AliAODVertex *aodvertex = fAOD->GetVertex(ivertex);
 570       if(!aodvertex) continue;
 571       if(aodvertex->GetType()==AliAODVertex::kKink) {
 572         AliAODTrack *mother = (AliAODTrack *) aodvertex->GetParent();
 573         if(!mother) continue;
 574         Int_t idmother = mother->GetID();
 575         listofmotherkink[numberofmotherkink] = idmother;
 576         numberofmotherkink++;
 577       }
 578     }
 579   } //+++
 580
 581
 582   ///////////////
 583   //Track loop///
 584   ///////////////
 585   for (Int_t iTracks = 0; iTracks < ntracks; iTracks++) {
 586
 587     AliVParticle* Vtrack = 0x0;
 588     if(!fUseTender) Vtrack  = fVevent->GetTrack(iTracks);
 589     if(fUseTender) Vtrack = dynamic_cast<AliVTrack*>(fTracks_tender->At(iTracks));
 590
 591     if (!Vtrack) {
 592       printf("ERROR: Could not receive track %d\n", iTracks);
 593       continue;
 594     }
 595     AliVTrack *track = dynamic_cast<AliVTrack*>(Vtrack);
 596     AliESDtrack *etrack = dynamic_cast<AliESDtrack*>(Vtrack);
 597     AliAODTrack *atrack = dynamic_cast<AliAODTrack*>(Vtrack);
 598
 599     ////////////////////
 600     //Apply track cuts//
 601     ////////////////////
 602     if(fAOD)
 603       if(!atrack->TestFilterMask(AliAODTrack::kTrkGlobalNoDCA)) continue; //mimimum cuts
 604
 605     if(fESD)
 606       if(!esdTrackCutsH->AcceptTrack(etrack))continue;
 607
 608     //reject kink
 609     if(IsAODanalysis()){
 610       Bool_t kinkmotherpass = kTRUE;
 611       for(Int_t kinkmother = 0; kinkmother < numberofmotherkink; kinkmother++) {
 612         if(track->GetID() == listofmotherkink[kinkmother]) {
 613           kinkmotherpass = kFALSE;
 614           continue;
 615         }
 616       }
 617       if(!kinkmotherpass) continue;
 618     }
 619     else{
 620       if(etrack->GetKinkIndex(0) != 0) continue;
 621     }
 622
 623     ////////////////////
 624     //Track properties//
 625     ///////////////////
 626     Double_t dEdx =-999, fTPCnSigma=-999;
 627     dEdx = track->GetTPCsignal();
 628     fTPCnSigma = fpidResponse->NumberOfSigmasTPC(track, AliPID::kElectron);
 629
 630     if(track->GetID()<0) fNegTrkIDPt->Fill(track->Pt());
 631     fTrkPt->Fill(track->Pt());
 632     fTrketa->Fill(track->Eta());
 633     fTrkphi->Fill(track->Phi());
 634     fdEdx->Fill(track->P(),dEdx);
 635     fTPCNpts->Fill(track->P(),track->GetTPCsignalN());
 636     fTPCnsig->Fill(track->P(),fTPCnSigma);
 637
 638     ///////////////////////////
 639     //Track matching to EMCAL//
 640     //////////////////////////
 641     Int_t EMCalIndex = -1;
 642     EMCalIndex = track->GetEMCALcluster();
 643     if(EMCalIndex < 0) continue;
 644     fHistPtMatch->Fill(track->Pt());
 645
 646     AliVCluster *clustMatch=0x0;
 647     if(!fUseTender) clustMatch = (AliVCluster*)fVevent->GetCaloCluster(EMCalIndex);
 648     if(fUseTender) clustMatch = dynamic_cast<AliVCluster*>(fCaloClusters_tender->At(EMCalIndex));
 649
 650     if(clustMatch && clustMatch->IsEMCAL())
 651     {
 652       /////////////////////////////////////////////
 653       //Properties of tracks matched to the EMCAL//
 654       /////////////////////////////////////////////
 655       fEMCTrkMatch->Fill(clustMatch->GetTrackDx(),clustMatch->GetTrackDz());
 656       if(TMath::Abs(clustMatch->GetTrackDx())>0.05 || TMath::Abs(clustMatch->GetTrackDz())>0.05) continue;
 657
 658       fEMCTrkPt->Fill(track->Pt());
 659       fEMCTrketa->Fill(track->Eta());
 660       fEMCTrkphi->Fill(track->Phi());
 661       fEMCdEdx->Fill(track->P(),dEdx);
 662       fEMCTPCnsig->Fill(track->P(),fTPCnSigma);
 663       fEMCTPCNpts->Fill(track->P(),track->GetTPCsignalN());
 664
 665       Double_t clustMatchE = clustMatch->E();
 666       fClsEAftMatch->Fill(clustMatchE);
 667
 668       //EMCAL EID info
 669       Double_t eop = -1.0;
 670       Double_t m02 = -99999;
 671       if(track->P()>0)eop = clustMatchE/track->P();
 672       m02 =clustMatch->GetM02();
 673
 674       if(track->Pt()>1.0){
 675         fHistdEdxEop->Fill(eop,dEdx);
 676         fHistNsigEop->Fill(eop,fTPCnSigma);
 677         fM20EovP->Fill(eop,clustMatch->GetM20());
 678         fM02EovP->Fill(eop,clustMatch->GetM02());
 679       }
 680
 681       fHistEop->Fill(track->Pt(),eop);
 682       fM20->Fill(track->Pt(),clustMatch->GetM20());
 683       fM02->Fill(track->Pt(),clustMatch->GetM02());
 684
 685       //EID THnsparse
 686       fvalueElectron[0] = trigger;
 687       fvalueElectron[1] = track->Pt();
 688       fvalueElectron[2] = fTPCnSigma;
 689       fvalueElectron[3] = eop;
 690       fvalueElectron[4] = clustMatch->GetM20();
 691       fvalueElectron[5] = clustMatch->GetM02();
 692
 693       if(fFlagSparse){
 694         //cout << "filling sparse"<<endl;
 695         fSparseElectron->Fill(fvalueElectron);
 696       }
 697
 698       ////////////////////////////////////////////////
 699       //Track properties of EMCAL electron cadidates//
 700       ///////////////////////////////////////////////
 701       if(fTPCnSigma > -1 && fTPCnSigma < 3 && eop>0.9 && eop<1.2 && m02 > 0.006 && m02 < 0.35){
 702         fEleCanTPCNpts->Fill(track->Pt(),track->GetTPCsignalN());
 703         fEleCanTPCNCls->Fill(track->Pt(),track->GetTPCNcls());
 704
 705         Int_t fITSncls=0;
 706         for(Int_t l=0;l<6;l++) {
 707           if(TESTBIT(track->GetITSClusterMap(),l)) {
 708             fEleCanITShit->Fill(l);
 709             if(l==0) fEleCanSPD1->Fill(track->Pt(),0.5);
 710             if(l==1) fEleCanSPD2->Fill(track->Pt(),0.5);
 711             if(l==0 && l==1) fEleCanSPDBoth->Fill(track->Pt(),0.5);
 712             if(l==0 || l==1) fEleCanSPDOr->Fill(track->Pt(),0.5);
 713             fITSncls++;
 714           }
 715         }
 716         fEleCanITSNCls->Fill(track->Pt(),fITSncls++);
 717       }
 718     }
 719   } //track loop
 720
 721   PostData(1, fOutputList);
 722 }
 723 //________________________________________________________________________
 724 void AliAnalysisTaskHFEemcQA::Terminate(Option_t *)
 725 {
 726   // Draw result to the screen
 727   // Called once at the end of the query
 728
 729   fOutputList = dynamic_cast<TList*> (GetOutputData(1));
 730   if (!fOutputList) {
 731     printf("ERROR: Output list not available\n");
 732     return;
 733   }
 734
 735 }