PWGCF/EBYE/NetParticle/AliAnalysisTaskNetParticle.cxx

   1 //-*- Mode: C++ -*-
   2
   3 #include "TFile.h"
   4 #include "TChain.h"
   5 #include "TTree.h"
   6 #include "TH1F.h"
   7 #include "TH2F.h"
   8 #include "TH3F.h"
   9 #include "TCanvas.h"
  10 #include "TStopwatch.h"
  11 #include "TMath.h"
  12 #include "THashList.h"
  13
  14 #include "AliAnalysisTask.h"
  15 #include "AliAnalysisManager.h"
  16 #include "AliTracker.h"
  17 #include "AliESDEvent.h"
  18 #include "AliESDInputHandler.h"
  19 #include "AliESDpid.h"
  20 #include "AliStack.h"
  21 #include "AliMCEvent.h"
  22 #include "AliMCEventHandler.h"
  23 #include "AliESDtrackCuts.h"
  24 #include "AliKineTrackCuts.h"
  25 #include "AliMCParticle.h"
  26 #include "AliESDVZERO.h"
  27 #include "AliAnalysisTaskNetParticle.h"
  28 #include "AliGenEventHeader.h"
  29 #include "AliCentrality.h"
  30
  31 using namespace std;
  32
  33 // Task for NetParticle checks
  34 // Author: Jochen Thaeder <jochen@thaeder.de>
  35
  36 ClassImp(AliAnalysisTaskNetParticle)
  37
  38 /*
  39  * ---------------------------------------------------------------------------------
  40  *                            Constructor / Destructor
  41  * ---------------------------------------------------------------------------------
  42  */
  43
  44 //________________________________________________________________________
  45 AliAnalysisTaskNetParticle::AliAnalysisTaskNetParticle(const char *name) :
  46   AliAnalysisTaskSE(name),
  47   fHelper(new AliAnalysisNetParticleHelper),
  48   fEffCont(NULL),
  49   fDCA(NULL),
  50   fDist(NULL),
  51
  52   fOutList(NULL),
  53   fOutListEff(NULL),
  54   fOutListCont(NULL),
  55   fOutListDCA(NULL),
  56   fOutListQA(NULL),
  57
  58   fESD(NULL),
  59   fESDHandler(NULL),
  60
  61   fESDTrackCutsBase(new AliESDtrackCuts),
  62   fESDTrackCuts(NULL),
  63   fESDTrackCutsBkg(NULL),
  64   fESDTrackCutsEff(NULL),
  65
  66   fIsMC(kFALSE),
  67   fESDTrackCutMode(0),
  68   fModeEffCreation(0),
  69   fModeDCACreation(0),
  70   fModeDistCreation(0),
  71
  72   fMCEvent(NULL),
  73   fMCStack(NULL),
  74
  75   fHnQA(NULL),
  76   fUseQATHnSparse(kFALSE),
  77
  78   fEtaMax(0.9),
  79   fPtRange(new Float_t[2]),
  80   fPtRangeEff(new Float_t[2]) {
  81   // Constructor
  82
  83   AliLog::SetClassDebugLevel("AliAnalysisTaskNetParticle",10);
  84
  85   fPtRange[0] = 0.4;
  86   fPtRange[1] = 0.8;
  87   fPtRangeEff[0] = 0.2;
  88   fPtRangeEff[1] = 1.6;
  89
  90   // -- Output slots
  91   // -------------------------------------------------
  92   DefineOutput(1, TList::Class());
  93   DefineOutput(2, TList::Class());
  94   DefineOutput(3, TList::Class());
  95   DefineOutput(4, TList::Class());
  96   DefineOutput(5, TList::Class());
  97 }
  98
  99 //________________________________________________________________________
 100 AliAnalysisTaskNetParticle::~AliAnalysisTaskNetParticle() {
 101   // Destructor
 102
 103   if (fPtRange)          delete[] fPtRange;
 104   if (fPtRangeEff)       delete[] fPtRangeEff;
 105
 106   if (fESDTrackCutsBase) delete fESDTrackCutsBase;
 107   if (fESDTrackCuts)     delete fESDTrackCuts;
 108   if (fESDTrackCutsBkg)  delete fESDTrackCutsBkg;
 109   if (fESDTrackCutsEff)  delete fESDTrackCutsEff;
 110
 111   if (fHelper)           delete fHelper;
 112   if (fEffCont)          delete fEffCont;
 113   if (fDCA)              delete fDCA;
 114   if (fDist)             delete fDist;
 115 }
 116
 117 /*
 118  * ---------------------------------------------------------------------------------
 119  *                                 Public Methods
 120  * ---------------------------------------------------------------------------------
 121  */
 122
 123 //________________________________________________________________________
 124 void AliAnalysisTaskNetParticle::UserCreateOutputObjects() {
 125   // Create histograms
 126
 127   Bool_t oldStatus = TH1::AddDirectoryStatus();
 128   TH1::AddDirectory(kFALSE);
 129
 130   fOutList = new TList;
 131   fOutList->SetName(GetName()) ;
 132   fOutList->SetOwner(kTRUE);
 133
 134   fOutListEff = new TList;
 135   fOutListEff->SetName(Form("%s_eff",GetName()));
 136   fOutListEff->SetOwner(kTRUE) ;
 137
 138   fOutListCont = new TList;
 139   fOutListCont->SetName(Form("%s_cont",GetName()));
 140   fOutListCont->SetOwner(kTRUE) ;
 141
 142   fOutListDCA = new TList;
 143   fOutListDCA->SetName(Form("%s_dca",GetName()));
 144   fOutListDCA->SetOwner(kTRUE) ;
 145
 146   fOutListQA = new TList;
 147   fOutListQA->SetName(Form("%s_qa",GetName()));
 148   fOutListQA->SetOwner(kTRUE) ;
 149
 150   // ------------------------------------------------------------------
 151   // -- Get event / trigger statistics histograms
 152   // ------------------------------------------------------------------
 153
 154   fOutList->Add(fHelper->GetHEventStat0());
 155   fOutList->Add(fHelper->GetHEventStat1());
 156   fOutList->Add(fHelper->GetHTriggerStat());
 157   fOutList->Add(fHelper->GetHCentralityStat());
 158
 159   // ------------------------------------------------------------------
 160   // -- Add histograms from distribution class
 161   // ------------------------------------------------------------------
 162   if (fModeDistCreation == 1)
 163     fDist->CreateHistograms(fOutList);
 164
 165   // ------------------------------------------------------------------
 166   // -- Add histograms from efficiency/contamination class
 167   // ------------------------------------------------------------------
 168   if (fIsMC && fModeEffCreation == 1) {
 169     fOutListEff->Add(fEffCont->GetHnEff());
 170     fOutListCont->Add(fEffCont->GetHnCont());
 171   }
 172
 173   // ------------------------------------------------------------------
 174   // -- Add histograms from DCA class
 175   // ------------------------------------------------------------------
 176   if (fModeDCACreation == 1)
 177     fOutListDCA->Add(fDCA->GetHnDCA());
 178
 179   // ------------------------------------------------------------------
 180   // -- Add Multiplicity Correlations
 181   // ------------------------------------------------------------------
 182   /*
 183   fOutList->Add(new TH2F("fHMultCorrTPCref0", "Centrality vs Multiplicity (TPCref);Centrality;Multiplicity (TPCref)",
 184                          24, -0.5, 11.49, 10001, 0.5, 10000.49));
 185   fOutList->Add(new TH2F("fHMultCorrTracks0", "Centrality vs Multiplicity (Tracks);Centrality;Multiplicity (Tracks)",
 186                          24, -0.5, 11.49, 10001, 0.5, 10000.49));
 187   fOutList->Add(new TH2F("fHMultCorrESDTracks0", "Centrality vs Multiplicity (ESDTracks);Centrality;Multiplicity (ESDTracks)",
 188                          24, -0.5, 11.49, 10001, 0.5, 10000.49));
 189
 190   fOutList->Add(new TH2F("fHMultCorrTPCref1", "Centrality vs Multiplicity (TPCref);Centrality;Multiplicity (TPCref)",
 191                          101, -0.5, 100.49, 10001, 0.5, 10000.49));
 192   fOutList->Add(new TH2F("fHMultCorrTracks1", "Centrality vs Multiplicity (Tracks);Centrality;Multiplicity (Tracks)",
 193                          101, -0.5, 100.49, 10001, 0.5, 10000.49));
 194   fOutList->Add(new TH2F("fHMultCorrESDTracks1", "Centrality vs Multiplicity (ESDTracks);Centrality;Multiplicity (ESDTracks)",
 195                          101, -0.5, 100.49, 10001, 0.5, 10000.49));
 196   */
 197   // ------------------------------------------------------------------
 198   // -- Create THnSparseF - QA
 199   // ------------------------------------------------------------------
 200   /*
 201   Double_t dCent[2] = {-0.5, 8.5};
 202   Int_t iCent       = 9;
 203
 204   Double_t dEta[2]  = {-0.9, 0.9}; // -> 37 bins
 205   Int_t iEta        = Int_t((dEta[1]-dEta[0]) / 0.075) +1 ;
 206
 207   Double_t dRap[2]  = {-0.5, 0.5};
 208   Int_t iRap        = Int_t((dRap[1]-dRap[0]) / 0.075) +1 ;
 209
 210   Double_t dPhi[2]  = {0.0, TMath::TwoPi()};
 211   Int_t iPhi        = 42;
 212
 213   Double_t dPt[2]   = {0.1, 3.0};
 214   Int_t iPt         = Int_t((dPt[1]-dPt[0]) / 0.075);
 215
 216   Double_t dSign[2] = {-1.5, 1.5};
 217   Int_t iSign       = 3;
 218
 219   //                      cent:isAccepted: pInner:     pt:     sign:tpcSignal:nSigmaTPC:nSigmaTOF:      eta:     phi:       y: dcar: dcaz: nSigmaCdd: nSigmaCzz
 220   Int_t    bin[15] = {   iCent,       2  ,    iPt,    iPt,    iSign,      500,     50  ,     50  ,     iEta,    iPhi,    iRap,  50 ,  50 ,       50 ,       50 };
 221   Double_t min[15] = {dCent[0],      -0.5, dPt[0], dPt[0], dSign[0],       30,     -5.0,     -5.0,  dEta[0], dPhi[0], dRap[0], -10., -10.,      -10.,      -10.};
 222   Double_t max[15] = {dCent[1],       1.5, dPt[1], dPt[1], dSign[1],      500,      5.0,      5.0,  dEta[1], dPhi[1], dRap[1],  10.,  10.,       10.,       10.};
 223
 224   if (fUseQATHnSparse) {
 225     fOutListQA->Add(new THnSparseF("fHnQA", "cent:isAccepted:pInner:pt:sign:tpcSignal:nSigmaTPC:nSigmaTOF:eta:phi:u:dcar:dcaz:nSigmaCdd:nSigmaCzz",
 226                                    15, bin, min, max));
 227
 228     fHnQA = static_cast<THnSparseF*>(fOutListQA->Last());
 229     fHnQA->Sumw2();
 230
 231     fHnQA->GetAxis(0)->SetTitle("centrality");
 232     fHnQA->GetAxis(1)->SetTitle("isAccepted");
 233     fHnQA->GetAxis(2)->SetTitle("#it{p}_{Inner} (GeV/#it{c})");
 234     fHnQA->GetAxis(3)->SetTitle("#it{p}_{T} (GeV/#it{c})");
 235     fHnQA->GetAxis(4)->SetTitle("sign");
 236
 237     fHnQA->GetAxis(5)->SetTitle("TPC signal");
 238     fHnQA->GetAxis(6)->SetTitle("n #sigma TPC");
 239     fHnQA->GetAxis(7)->SetTitle("n #sigma TOF");
 240
 241     fHnQA->GetAxis(8)->SetTitle("#eta");
 242     fHnQA->GetAxis(9)->SetTitle("#varphi");
 243     fHnQA->GetAxis(10)->SetTitle("#it{y}");
 244
 245     fHnQA->GetAxis(11)->SetTitle("DCAr");
 246     fHnQA->GetAxis(12)->SetTitle("DCAz");
 247
 248     fHnQA->GetAxis(13)->SetTitle("n #sigma #sqrt(Cdd)/DCAr");
 249     fHnQA->GetAxis(14)->SetTitle("n #sigma #sqrt(Czz)/DCAz");
 250
 251     fHelper->BinLogAxis(fHnQA, 2);
 252     fHelper->BinLogAxis(fHnQA, 3);
 253   }
 254   */
 255   // ------------------------------------------------------------------
 256
 257   TH1::AddDirectory(oldStatus);
 258
 259   return;
 260 }
 261
 262 //________________________________________________________________________
 263 void AliAnalysisTaskNetParticle::UserExec(Option_t *) {
 264   // Called for each event
 265
 266   // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 267   // -- Setup Event
 268   // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 269   if (SetupEvent() < 0) {
 270     PostData(1,fOutList);
 271     PostData(2,fOutListEff);
 272     PostData(3,fOutListCont);
 273     PostData(4,fOutListDCA);
 274     PostData(5,fOutListQA);
 275     return;
 276   }
 277
 278   // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 279   // -- Process Efficiency / Contamination Determination
 280   // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 281   if (fIsMC && fModeEffCreation == 1)
 282     fEffCont->Process();
 283
 284   // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 285   // -- Process DCA Determination
 286   // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 287   if (fModeDCACreation == 1)
 288     fDCA->Process();
 289
 290   // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 291   // -- Process Distributions
 292   // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 293   if (fModeDistCreation == 1)
 294     fDist->Process();
 295
 296   // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 297   // -- Post output data
 298   // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 299
 300   PostData(1,fOutList);
 301   PostData(2,fOutListEff);
 302   PostData(3,fOutListCont);
 303   PostData(4,fOutListDCA);
 304   PostData(5,fOutListQA);
 305
 306   return;
 307 }
 308
 309 //________________________________________________________________________
 310 void AliAnalysisTaskNetParticle::Terminate(Option_t *){
 311   // Terminate
 312 }
 313
 314 /*
 315  * ---------------------------------------------------------------------------------
 316  *                                 Public Methods
 317  * ---------------------------------------------------------------------------------
 318  */
 319
 320 //________________________________________________________________________
 321 Int_t AliAnalysisTaskNetParticle::Initialize() {
 322   // Initialize event
 323
 324   // ------------------------------------------------------------------
 325   // -- ESD TrackCuts
 326   // ------------------------------------------------------------------
 327   TString sModeName("");
 328
 329   // -- Create ESD track cuts
 330   // --------------------------
 331   fESDTrackCutsBase->SetMinNCrossedRowsTPC(70);                                             // TPC
 332   fESDTrackCutsBase->SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);                    // TPC
 333
 334   fESDTrackCutsBase->SetMaxChi2PerClusterTPC(4);                                            // TPC
 335   fESDTrackCutsBase->SetAcceptKinkDaughters(kFALSE);                                        // TPC
 336   fESDTrackCutsBase->SetRequireTPCRefit(kTRUE);                                             // TPC
 337
 338   fESDTrackCutsBase->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kOff); // ITS
 339   fESDTrackCutsBase->SetClusterRequirementITS(AliESDtrackCuts::kSDD,AliESDtrackCuts::kOff); // ITS
 340   fESDTrackCutsBase->SetClusterRequirementITS(AliESDtrackCuts::kSSD,AliESDtrackCuts::kOff); // ITS
 341
 342   fESDTrackCutsBase->SetDCAToVertex2D(kFALSE);                                              // VertexConstrained
 343   fESDTrackCutsBase->SetRequireSigmaToVertex(kFALSE);                                       // VertexConstrained
 344
 345   fESDTrackCutsBase->SetEtaRange(-1.*fEtaMax, fEtaMax);                                     // Acceptance
 346   fESDTrackCutsBase->SetPtRange(fPtRange[0],fPtRange[1]);                                   // Acceptance
 347
 348   // -- Mode : clean cuts -> small contamination
 349   if (fESDTrackCutMode == 0) {
 350     sModeName = "Clean";
 351     fESDTrackCutsBase->SetRequireITSRefit(kTRUE);                                           // ITS
 352     fESDTrackCutsBase->SetMaxChi2PerClusterITS(36);                                         // ITS
 353   }
 354   // -- Mode : dirty cuts -> high efficiency
 355   else if (fESDTrackCutMode == 1) {
 356     sModeName = "Dirty";
 357     fESDTrackCutsBase->SetRequireITSRefit(kFALSE);                                          // ITS
 358   }
 359   // -- Mode : Default
 360   else {
 361     sModeName = "Base";
 362   }
 363
 364   fESDTrackCutsBase->SetName(Form("NetParticleCuts2010_%s",sModeName.Data()));
 365
 366   // -- Create ESD BKG track cuts
 367   // ------------------------------
 368   fESDTrackCutsBkg = static_cast<AliESDtrackCuts*>(fESDTrackCutsBase->Clone());
 369   fESDTrackCutsBkg->SetName(Form("NetParticleCuts2010_%s_Bkg",sModeName.Data()));
 370   fESDTrackCutsBkg->SetMaxDCAToVertexZ(10.);                                                // VertexConstrained
 371
 372   // -- Create ESD track cuts
 373   // ------------------------------
 374   fESDTrackCuts = static_cast<AliESDtrackCuts*>(fESDTrackCutsBase->Clone());
 375   fESDTrackCuts->SetName(Form("NetParticleCuts2010_%s",sModeName.Data()));
 376   fESDTrackCuts->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");                         // VertexConstrained  ->  7*(0.0026+0.0050/pt^1.01)
 377   fESDTrackCuts->SetMaxChi2TPCConstrainedGlobal(36);                                        // VertexConstrained
 378   fESDTrackCuts->SetMaxDCAToVertexZ(2);                                                     // VertexConstrained
 379
 380   // -- Create ESD Eff track cuts
 381   // ------------------------------
 382   fESDTrackCutsEff = static_cast<AliESDtrackCuts*>(fESDTrackCuts->Clone());
 383   fESDTrackCutsEff->SetName(Form("NetParticleCuts2010_%s_Eff",sModeName.Data()));
 384   fESDTrackCutsEff->SetPtRange(fPtRangeEff[0],fPtRangeEff[1]);                              // Acceptance
 385
 386   // ------------------------------------------------------------------
 387   // -- Initialize Helper
 388   // ------------------------------------------------------------------
 389   if (fHelper->Initialize(fIsMC))
 390     return -1;
 391
 392   // ------------------------------------------------------------------
 393   // -- Create / Initialize Efficiency/Contamination
 394   // ------------------------------------------------------------------
 395   if (fIsMC && fModeEffCreation == 1) {
 396     fEffCont = new AliAnalysisNetParticleEffCont;
 397     fEffCont->Initialize(fESDTrackCutsEff, fHelper);
 398   }
 399
 400   // ------------------------------------------------------------------
 401   // -- Create / Initialize DCA Determination
 402   // ------------------------------------------------------------------
 403   if (fModeDCACreation == 1) {
 404     fDCA = new AliAnalysisNetParticleDCA;
 405     fDCA->Initialize(fESDTrackCuts, fESDTrackCutsBkg, fHelper);
 406   }
 407
 408   // ------------------------------------------------------------------
 409   // -- Create / Initialize DCA Determination
 410   // ------------------------------------------------------------------
 411   if (fModeDistCreation == 1) {
 412     fDist = new AliAnalysisNetParticleDistribution;
 413     fDist->Initialize(fHelper, fESDTrackCuts, fIsMC, fPtRange, fEtaMax);
 414   }
 415
 416   // ------------------------------------------------------------------
 417   // -- Reset Event
 418   // ------------------------------------------------------------------
 419   ResetEvent();
 420
 421   return 0;
 422 }
 423
 424 /*
 425  * ---------------------------------------------------------------------------------
 426  *                            Setup/Reset Methods - private
 427  * ---------------------------------------------------------------------------------
 428  */
 429
 430 //________________________________________________________________________
 431 Int_t AliAnalysisTaskNetParticle::SetupEvent() {
 432   // Setup Reading of event
 433   // > return 0 for success / accepted event
 434   // > return -1 for failed setup
 435   // > return -2 for rejected event
 436
 437   ResetEvent();
 438
 439   // -- ESD Event
 440   // ------------------------------------------------------------------
 441   if (SetupESDEvent() < 0) {
 442     AliError("Setup ESD Event failed");
 443     return -1;
 444   }
 445
 446   // -- Setup MC Event
 447   // ------------------------------------------------------------------
 448   if (fIsMC && SetupMCEvent() < 0) {
 449     AliError("Setup MC Event failed");
 450     return -1;
 451   }
 452
 453   // -- Setup Event for Helper / EffCont  / DCA / Dist classes
 454   // ------------------------------------------------------------------
 455   fHelper->SetupEvent(fESDHandler, fMCEvent);
 456
 457   if (fIsMC && fModeEffCreation)
 458     fEffCont->SetupEvent(fESDHandler, fMCEvent);
 459
 460   if (fModeDCACreation == 1)
 461     fDCA->SetupEvent(fESDHandler, fMCEvent);
 462
 463   if (fModeDistCreation == 1)
 464     fDist->SetupEvent(fESDHandler, fMCEvent);
 465
 466   // -- Evaluate Event cuts
 467   // ------------------------------------------------------------------
 468   return fHelper->IsEventRejected() ? -2 : 0;
 469 }
 470
 471 //________________________________________________________________________
 472 Int_t AliAnalysisTaskNetParticle::SetupESDEvent() {
 473   // -- Setup ESD Event
 474   // > return 0 for success
 475   // > return -1 for failed setup
 476
 477   fESDHandler= dynamic_cast<AliESDInputHandler*>
 478     (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
 479   if (!fESDHandler) {
 480     AliError("Could not get ESD input handler");
 481     return -1;
 482   }
 483
 484   fESD = fESDHandler->GetEvent();
 485   if (!fESD) {
 486     AliError("Could not get ESD event");
 487     return -1;
 488   }
 489
 490   // -- Check PID response
 491   // ------------------------------------------------------------------
 492   if (!fESDHandler->GetPIDResponse()) {
 493     AliError("Could not get PID response");
 494     return -1;
 495   }
 496
 497   // -- Check Vertex
 498   // ------------------------------------------------------------------
 499   if (!fESD->GetPrimaryVertexTracks()) {
 500     AliError("Could not get vertex from tracks");
 501     return -1;
 502   }
 503
 504   // -- Check Centrality
 505   // ------------------------------------------------------------------
 506   if (!fESD->GetCentrality()) {
 507     AliError("Could not get centrality");
 508     return -1;
 509   }
 510
 511   return 0;
 512 }
 513
 514 //________________________________________________________________________
 515 Int_t AliAnalysisTaskNetParticle::SetupMCEvent() {
 516   // -- Setup MC Event
 517   // > return 0 for success
 518   // > return -1 for failed setup
 519
 520   AliMCEventHandler *mcH = dynamic_cast<AliMCEventHandler*>
 521     (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
 522
 523   fMCEvent = mcH->MCEvent();
 524   if (!fMCEvent) {
 525     AliError("MC event not available");
 526     return -1;
 527   }
 528
 529   // -- Get MC header
 530   // ------------------------------------------------------------------
 531   AliHeader* header = fMCEvent->Header();
 532   if (!header) {
 533     AliError("MC header not available");
 534     return -1;
 535   }
 536
 537   // -- Check Stack
 538   // ------------------------------------------------------------------
 539   fMCStack = fMCEvent->Stack();
 540   if (!fMCStack) {
 541     AliError("MC stack not available");
 542     return -1;
 543   }
 544
 545   // -- Check GenHeader
 546   // ------------------------------------------------------------------
 547   if (!header->GenEventHeader()) {
 548     AliError("Could not retrieve genHeader from header");
 549     return -1;
 550   }
 551
 552   // -- Check primary vertex
 553   // ------------------------------------------------------------------
 554   if (!fMCEvent->GetPrimaryVertex()){
 555     AliError("Could not get MC vertex");
 556     return -1;
 557   }
 558
 559   return 0;
 560 }
 561
 562 //________________________________________________________________________
 563 void AliAnalysisTaskNetParticle::ResetEvent() {
 564   // -- Reset event
 565
 566   // -- Reset ESD Event
 567   fESD       = NULL;
 568
 569   // -- Reset MC Event
 570   if (fIsMC)
 571     fMCEvent = NULL;
 572
 573   // -- Reset Dist Creation
 574   if (fModeDistCreation == 1)
 575     fDist->ResetEvent();
 576
 577   return;
 578 }
 579
 580 /*
 581  * ---------------------------------------------------------------------------------
 582  *                           Process Methods - private
 583  * ---------------------------------------------------------------------------------
 584  */
 585
 586 #if 0
 587
 588 //________________________________________________________________________
 589 void AliAnalysisTaskNetParticle::ProcessESDTrackLoop() {
 590   // -- Process ESD tracks and fill histograms
 591
 592   Int_t nTracks = 0;
 593
 594   for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
 595     AliESDtrack *track = fESD->GetTrack(idxTrack);
 596
 597     // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 598     // -- Check track cuts
 599     // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 600
 601     // -- Check if charged track is accepted for basic parameters
 602     if (!fHelper->IsTrackAcceptedBasicCharged(track))
 603       continue;
 604
 605     // -- Check if accepted - BKG
 606     if (!fESDTrackCutsBkg->AcceptTrack(track))
 607       continue;
 608
 609     // -- Check if accepted in rapidity window
 610     Double_t yP;
 611     if (!fHelper->IsTrackAcceptedRapidity(track, yP))
 612       continue;
 613
 614     // -- Check if accepted bt PID from TPC or TPC+TOF
 615     Double_t pid[2];
 616     Bool_t isAcceptedPID = fHelper->IsTrackAcceptedPID(track, pid);
 617
 618     // -- Check if accepted with thighter DCA cuts
 619     Bool_t isAcceptedDCA = fHelper->IsTrackAcceptedDCA(track);
 620
 621     // -- Check track cuts
 622     Bool_t isAcceptedVertex = fESDTrackCuts->AcceptTrack(track);
 623
 624     // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 625     // -- Fill tracks in QA THnSparseF
 626     // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 627
 628     if (fUseQATHnSparse && fHnQA->GetEntries() < 5e5) {
 629
 630       // -- Get dca r/z
 631       Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
 632       track->GetImpactParameters(dca,cov);
 633
 634       Float_t dcaRoverCdd = ( TMath::Sqrt(cov[0]) != 0. )  ? dca[0]/TMath::Sqrt(cov[0]) : -9.99;
 635       Float_t dcaZoverCzz = ( TMath::Sqrt(cov[2]) != 0. )  ? dca[1]/TMath::Sqrt(cov[2]) : -9.99;
 636
 637       Double_t aQA[15] = {
 638         Double_t(fHelper->GetCentralityBin()),                     //  0 centrality
 639         Double_t(isAcceptedVertex || isAcceptedDCA),  //  1 isAccepted -> Vertex || isAcceptedDCA
 640         track->GetInnerParam()->GetP(),               //  2 p at InnerParam
 641         track->Pt(),                                  //  3 pt
 642         track->GetSign(),                             //  4 sign
 643         track->GetTPCsignal(),                        //  5 TPC dE/dx
 644         pid[0],                                       //  6 n Sigma TPC
 645         pid[1],                                       //  7 n Sigma TOF
 646         track->Eta(),                                 //  8 eta
 647         track->Phi(),                                 //  9 phi
 648         yP,                                           // 10 rapidity
 649         dca[0],                                       // 11 dca r
 650         dca[1],                                       // 12 dca z
 651         dcaRoverCdd,                                  // 13 sqrt(cov[dd])
 652         dcaZoverCzz,                                  // 14 sqrt(cov[zz])
 653       };
 654
 655       fHnQA->Fill(aQA);
 656     }
 657
 658     // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 659     // -- Apply track cuts
 660     // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 661
 662     if (!isAcceptedVertex || !isAcceptedDCA)
 663       continue;
 664
 665     ++nTracks;
 666
 667     if (!isAcceptedPID)
 668       continue;
 669
 670   } // for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
 671
 672   // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 673   // -- Fill Particle Fluctuation Histograms
 674   // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
 675
 676   /*
 677   Int_t tpcRef = fESDTrackCuts->GetReferenceMultiplicity(fESD,kTRUE);
 678
 679   (static_cast<TH2F*>(fOutList->FindObject("fHMultCorrTPCref0")))->Fill(fHelper->GetCentralityBin(),  tpcRef);
 680   (static_cast<TH2F*>(fOutList->FindObject("fHMultCorrTracks0")))->Fill(fHelper->GetCentralityBin(),  nTracks);
 681
 682   Float_t centPercentile = fHelper->GetCentralityPercentile();
 683   (static_cast<TH2F*>(fOutList->FindObject("fHMultCorrTPCref1")))->Fill(centPercentile,  tpcRef);
 684   (static_cast<TH2F*>(fOutList->FindObject("fHMultCorrTracks1")))->Fill(centPercentile,  nTracks);
 685 */
 686   return;
 687 }
 688
 689 #endif