PWGLF/STRANGENESS/Hypernuclei/AliAnalysisTaskNucleiv2SP.cxx

   1 /**************************************************************************
   2  * Contributors are not mentioned at all.                                 *
   3  *                                                                        *
   4  * Permission to use, copy, modify and distribute this software and its   *
   5  * documentation strictly for non-commercial purposes is hereby granted   *
   6  * without fee, provided that the above copyright noticxse appears in all *
   7  * copies and that both the copyright notice and this permission notice   *
   8  * appear in the supporting documentation. The authors make no claims     *
   9  * about the suitability of this software for any purpose. It is          *
  10  * provided "as is" without express or implied warranty.                  *
  11  **************************************************************************/
  12 //-----------------------------------------------------------------
  13 //                 AliAnalysisTaskNucleiv2SP class
  14 //-----------------------------------------------------------------
  15
  16 class TTree;
  17 class TParticle;
  18 class TVector3;
  19
  20 #include "AliAnalysisManager.h"
  21 #include <AliMCEventHandler.h>
  22 #include <AliMCEvent.h>
  23 #include <AliStack.h>
  24
  25 class AliESDVertex;
  26 class AliAODVertex;
  27 class AliESDv0;
  28 class AliAODv0;
  29
  30 #include <iostream>
  31
  32 #include "TList.h"
  33 #include "TH1.h"
  34 #include "TH2.h"
  35 #include "TH3.h"
  36 #include "TNtuple.h"
  37 #include "TGraph.h"
  38 #include "TF1.h"
  39 #include "TCanvas.h"
  40 #include "TMath.h"
  41 #include "TChain.h"
  42 #include "Riostream.h"
  43 #include "AliLog.h"
  44 #include "AliCascadeVertexer.h"
  45 #include "AliESDEvent.h"
  46 #include "AliESDtrack.h"
  47 #include "AliExternalTrackParam.h"
  48 #include "AliAODEvent.h"
  49 #include "AliInputEventHandler.h"
  50 #include "AliESDcascade.h"
  51 #include "AliAODcascade.h"
  52 #include "AliAnalysisTaskNucleiv2SP.h"
  53 #include "AliESDtrackCuts.h"
  54 #include "AliCentrality.h"
  55
  56
  57 ClassImp(AliAnalysisTaskNucleiv2SP)
  58
  59 using std::cout;
  60 using std::endl;
  61
  62 //________________________________________________________________________
  63 AliAnalysisTaskNucleiv2SP::AliAnalysisTaskNucleiv2SP()
  64 : AliAnalysisTaskSE(),
  65   fisPrimCut(kFALSE),
  66   fptc(1),
  67   fListHist(0),
  68   fHistEventMultiplicity(0),
  69   fHistTrackMultiplicity(0),
  70   fHistTrackMultiplicityCentral(0),
  71   fHistTrackMultiplicitySemiCentral(0),
  72   fHistTrackMultiplicityMB(0),
  73   fhBB(0),
  74   fhBBDeu(0),
  75   fhTOF(0),
  76   fhMassTOF(0),
  77   EPVzAvsCentrality(0),
  78   EPVzCvsCentrality(0),
  79   EPTPCvsCentrality(0),
  80   EPVzvsCentrality(0),
  81   EPTPCpvsCentrality(0),
  82   EPTPCnvsCentrality(0),
  83   hEvPlaneTPCvsEvPVz05(0),
  84   hEvPlaneTPCvsEvPVz075(0),
  85   hEvPlaneTPCvsEvPVz1530(0),
  86   hEvPlaneTPCvsEvPVz3050(0),
  87   hEvPlaneTPCvsEvPVz2040(0),
  88   hEvPlaneTPCvsEvPVz4060(0),
  89   hCos2DeltaTPCVzAvsCentrality(0),
  90   hCos2DeltaTPCVzCvsCentrality(0),
  91   hCos2DeltaVzAVzCvsCentrality(0),
  92   hCos2DeltaVzMVzAvsCentrality(0),
  93   hCos2DeltaVzMVzCvsCentrality(0),
  94   hCos2DeltaVzATPCvsCentrality(0),
  95   hCos2DeltaVzCTPCvsCentrality(0),
  96   hCos2DeltaVzCVzAvsCentrality(0),
  97   hCos2DeltaVzMTPCpvsCentrality(0),
  98   hCos2DeltaVzMTPCnvsCentrality(0),
  99   hCos2DeltaTPCpTPCnvsCentrality(0),
 100   hQVzAQVzCvsCentrality(0),
 101   ftree(0),
 102   tCentrality(0),
 103   tType(0),
 104   tHasTOF(0),
 105   tpT(0),
 106   tMassTOF(0),
 107   tuqV0A(0),
 108   tuqV0C(0),
 109   tCharge(0),
 110   tCosdeltaphiTPC(0),
 111   tCosdeltaphiV0M(0),
 112   tCosdeltaphiV0A(0),
 113   tCosdeltaphiV0C(0),
 114   timpactXY(0),
 115   timpactZ(0),
 116   tpull(0),
 117   tphi(0),
 118   fESDtrackCuts(0),
 119   fESDtrackCutsEP(0),
 120   fPIDResponse(0)
 121 {
 122   // Dummy Constructor
 123   fESDtrackCuts   = new AliESDtrackCuts("AliESDtrackCuts","AliESDtrackCuts");
 124   fESDtrackCutsEP = new AliESDtrackCuts("AliESDtrackCutsEP","AliESDtrackCutsEP");
 125   Initialize();
 126   cout<<"Dummy constructor"<<endl;
 127 }
 128
 129 //________________________________________________________________________
 130 AliAnalysisTaskNucleiv2SP::AliAnalysisTaskNucleiv2SP(const char *name)
 131 : AliAnalysisTaskSE(name),
 132   fisPrimCut(kFALSE),
 133   fptc(1),
 134   fListHist(0),
 135   fHistEventMultiplicity(0),
 136   fHistTrackMultiplicity(0),
 137   fHistTrackMultiplicityCentral(0),
 138   fHistTrackMultiplicitySemiCentral(0),
 139   fHistTrackMultiplicityMB(0),
 140   fhBB(0),
 141   fhBBDeu(0),
 142   fhTOF(0),
 143   fhMassTOF(0),
 144   EPVzAvsCentrality(0),
 145   EPVzCvsCentrality(0),
 146   EPTPCvsCentrality(0),
 147   EPVzvsCentrality(0),
 148   EPTPCpvsCentrality(0),
 149   EPTPCnvsCentrality(0),
 150   hEvPlaneTPCvsEvPVz05(0),
 151   hEvPlaneTPCvsEvPVz075(0),
 152   hEvPlaneTPCvsEvPVz1530(0),
 153   hEvPlaneTPCvsEvPVz3050(0),
 154   hEvPlaneTPCvsEvPVz2040(0),
 155   hEvPlaneTPCvsEvPVz4060(0),
 156   hCos2DeltaTPCVzAvsCentrality(0),
 157   hCos2DeltaTPCVzCvsCentrality(0),
 158   hCos2DeltaVzAVzCvsCentrality(0),
 159   hCos2DeltaVzMVzAvsCentrality(0),
 160   hCos2DeltaVzMVzCvsCentrality(0),
 161   hCos2DeltaVzATPCvsCentrality(0),
 162   hCos2DeltaVzCTPCvsCentrality(0),
 163   hCos2DeltaVzCVzAvsCentrality(0),
 164   hCos2DeltaVzMTPCpvsCentrality(0),
 165   hCos2DeltaVzMTPCnvsCentrality(0),
 166   hCos2DeltaTPCpTPCnvsCentrality(0),
 167   hQVzAQVzCvsCentrality(0),
 168   ftree(0),
 169   tCentrality(0),
 170   tType(0),
 171   tHasTOF(0),
 172   tpT(0),
 173   tMassTOF(0),
 174   tuqV0A(0),
 175   tuqV0C(0),
 176   tCharge(0),
 177   tCosdeltaphiTPC(0),
 178   tCosdeltaphiV0M(0),
 179   tCosdeltaphiV0A(0),
 180   tCosdeltaphiV0C(0),
 181   timpactXY(0),
 182   timpactZ(0),
 183   tpull(0),
 184   tphi(0),
 185   fESDtrackCuts(0),
 186   fESDtrackCutsEP(0),
 187   fPIDResponse(0)
 188 {
 189   // Define input and output slots here
 190   // Input slot #0 works with a TChain
 191   //DefineInput(0, TChain::Class());
 192   // Output slot #0 writes into a TList container ()
 193
 194   //
 195   // create track cuts
 196   //
 197   fESDtrackCuts   = new AliESDtrackCuts("AliESDtrackCuts","AliESDtrackCuts");
 198   fESDtrackCutsEP = new AliESDtrackCuts("AliESDtrackCutsEP","AliESDtrackCutsEP");
 199   //
 200   cout<<"Real constructor"<<endl;
 201   Initialize();
 202
 203   DefineInput(0, TChain::Class());
 204   DefineOutput(1, TList::Class());
 205   DefineOutput(2, TTree::Class());
 206
 207 }
 208
 209 void AliAnalysisTaskNucleiv2SP::Initialize()
 210 {
 211   //
 212   // updating parameters in case of changes
 213   //
 214   fESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(fisPrimCut,kTRUE);
 215   fESDtrackCuts->SetMaxDCAToVertexXY(3);
 216   fESDtrackCuts->SetMaxDCAToVertexZ(2);
 217   fESDtrackCuts->SetEtaRange(-0.8,0.8);
 218
 219   fESDtrackCutsEP = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
 220
 221 }
 222
 223 //________________________________________________________________________
 224 Float_t AliAnalysisTaskNucleiv2SP::GetEventPlaneForCandidate(AliESDtrack* track0, const TVector2* q,AliEventplane *pl){
 225
 226   // remove autocorrelations
 227
 228   TArrayF* qx = 0x0;
 229   TArrayF* qy = 0x0;
 230   TVector2 qcopy;
 231   // if(!fEtaGap){
 232   qx = pl->GetQContributionXArray();
 233   qy = pl->GetQContributionYArray();
 234   qcopy = *q;
 235
 236   TVector2 q0;
 237   if((track0->GetID()) < qx->fN){
 238     q0.Set(qx->At(track0->GetID()),qy->At(track0->GetID()));
 239   }
 240
 241   qcopy = qcopy - q0;
 242
 243   return qcopy.Phi()/2.;
 244
 245 }
 246 //________________________________________________________________________
 247 Float_t AliAnalysisTaskNucleiv2SP::GetPhi0Pi(Float_t phi){
 248   // Sets the phi angle in the range 0-pi
 249   Float_t result=phi;
 250   while(result<0){
 251     result=result+TMath::Pi();
 252   }
 253   while(result>TMath::Pi()){
 254     result=result-TMath::Pi();
 255   }
 256    return result;
 257 }
 258
 259 //==================DEFINITION OF OUTPUT OBJECTS==============================
 260
 261 void AliAnalysisTaskNucleiv2SP::UserCreateOutputObjects()
 262 {
 263   //-------------------------------------------------------
 264   fListHist = new TList();
 265   fListHist->SetOwner();  // IMPORTANT!
 266
 267   if(! fHistEventMultiplicity ){
 268
 269     fHistEventMultiplicity   = new TH1F( "fHistEventMultiplicity" , "Nb of Events" , 12 , -0.5,11.5);
 270
 271     fHistEventMultiplicity->GetXaxis()->SetBinLabel(1,"All Events");
 272     fHistEventMultiplicity->GetXaxis()->SetBinLabel(2,"Events w/PV");
 273     fHistEventMultiplicity->GetXaxis()->SetBinLabel(3,"Events w/|Vz|<10cm");
 274     fHistEventMultiplicity->GetXaxis()->SetBinLabel(4,"Central Events");
 275     fHistEventMultiplicity->GetXaxis()->SetBinLabel(5,"Semi-Central Events");
 276     fHistEventMultiplicity->GetXaxis()->SetBinLabel(6,"MB Events");
 277     //from HF
 278     fHistEventMultiplicity->GetXaxis()->SetBinLabel(7,"nEventsAnal");
 279     fHistEventMultiplicity->GetXaxis()->SetBinLabel(8,"nEvSelected");
 280     fHistEventMultiplicity->GetXaxis()->SetBinLabel(9,"nCandidatesSelected");
 281     fHistEventMultiplicity->GetXaxis()->SetBinLabel(10,"out of pt bounds");
 282     fHistEventMultiplicity->GetXaxis()->SetBinLabel(11,"mismatch lab");
 283     fHistEventMultiplicity->GetXaxis()->SetBinLabel(12,"non valid TPC EP");
 284     fListHist->Add(fHistEventMultiplicity);
 285   }
 286
 287   if(! fHistTrackMultiplicity ){
 288     fHistTrackMultiplicity  = new TH2F( "fHistTrackMultiplicity", "Nb of Tracks MB Events |Vz| < 10", 250,0, 25000,105,-0.5,104.5);
 289     fHistTrackMultiplicity->GetXaxis()->SetTitle("Number of tracks");
 290     fHistTrackMultiplicity->GetYaxis()->SetTitle("Percentile");
 291     fListHist->Add(fHistTrackMultiplicity);
 292   }
 293
 294   if(! fHistTrackMultiplicityCentral ){
 295     fHistTrackMultiplicityCentral  = new TH2F( "fHistTrackMultiplicityCentral", "Nb of Tracks MB Events |Vz| < 10", 250,0, 25000,105,-0.5,104.5);
 296     fHistTrackMultiplicityCentral->GetXaxis()->SetTitle("Number of tracks");
 297     fHistTrackMultiplicityCentral->GetYaxis()->SetTitle("Percentile");
 298     fListHist->Add(fHistTrackMultiplicityCentral);
 299   }
 300   if(! fHistTrackMultiplicitySemiCentral ){
 301     fHistTrackMultiplicitySemiCentral  = new TH2F( "fHistTrackMultiplicitySemiCentral", "Nb of Tracks MB Events |Vz| < 10", 250,0, 25000,105,-0.5,104.5);
 302     fHistTrackMultiplicitySemiCentral->GetXaxis()->SetTitle("Number of tracks");
 303     fHistTrackMultiplicitySemiCentral->GetYaxis()->SetTitle("Percentile");
 304     fListHist->Add(fHistTrackMultiplicitySemiCentral);
 305   }
 306   if(! fHistTrackMultiplicityMB ){
 307     fHistTrackMultiplicityMB  = new TH2F( "fHistTrackMultiplicityMB", "Nb of Tracks MB Events |Vz| < 10", 250,0, 25000,105,-0.5,104.5);
 308     fHistTrackMultiplicityMB->GetXaxis()->SetTitle("Number of tracks");
 309     fHistTrackMultiplicityMB->GetYaxis()->SetTitle("Percentile");
 310     fListHist->Add(fHistTrackMultiplicityMB);
 311   }
 312
 313   if(! fhBB ){
 314     fhBB = new TH2F( "fhBB" , "BetheBlochTPC" , 240,-6,6,250,0,1000);
 315     fListHist->Add(fhBB);
 316   }
 317
 318   if(! fhBBDeu ){
 319     fhBBDeu = new TH2F( "fhBBDeu" , "BetheBlochTPC - Deuteron" , 240,-6,6,250,0,1000);
 320     fListHist->Add(fhBBDeu);
 321   }
 322
 323   if(! fhTOF ){
 324     fhTOF = new TH2F( "fhTOF" , "Scatter Plot TOF" , 240,-6,6,500,0,1.2);
 325     fListHist->Add(fhTOF);
 326   }
 327   if(! fhMassTOF){
 328     fhMassTOF=new TH1F ("fhMassTOF","Particle Mass - TOF", 100,0 ,10);
 329     fListHist->Add(fhMassTOF);
 330   }
 331
 332   EPVzAvsCentrality  = new TH2D("EPVzAvsCentrality" , "EPVzAvsCentrality" , 80, -2, 2,105,-0.5,105.5);
 333   EPVzCvsCentrality  = new TH2D("EPVzCvsCentrality" , "EPVzCvsCentrality" , 80, -2, 2,105,-0.5,105.5);
 334   EPTPCvsCentrality  = new TH2D("EPTPCvsCentrality" , "EPTPCvsCentrality" , 80, -2, 2,105,-0.5,105.5);
 335   EPVzvsCentrality   = new TH2D("EPVzvsCentrality"  , "EPVzvsCentrality"  , 80, -2, 2,105,-0.5,105.5);
 336   EPTPCpvsCentrality = new TH2D("EPTPCpvsCentrality", "EPTPCpvsCentrality", 80, -2, 2,105,-0.5,105.5);
 337   EPTPCnvsCentrality = new TH2D("EPTPCnvsCentrality", "EPTPCnvsCentrality", 80, -2, 2,105,-0.5,105.5);
 338
 339   fListHist->Add(EPVzAvsCentrality);
 340   fListHist->Add(EPVzCvsCentrality);
 341   fListHist->Add(EPTPCvsCentrality);
 342   fListHist->Add(EPVzvsCentrality);
 343   fListHist->Add(EPTPCpvsCentrality);
 344   fListHist->Add(EPTPCnvsCentrality);
 345
 346   hEvPlaneTPCvsEvPVz05   = new TH2F("hEvPlaneTPCvsEvPVz05"  ,"hEvPlaneTPCvsEvPVz05"  ,100,-TMath::Pi(),TMath::Pi(),100,-TMath::Pi(),TMath::Pi());
 347   hEvPlaneTPCvsEvPVz075  = new TH2F("hEvPlaneTPCvsEvPVz075" ,"hEvPlaneTPCvsEvPVz075" ,100,-TMath::Pi(),TMath::Pi(),100,-TMath::Pi(),TMath::Pi());
 348   hEvPlaneTPCvsEvPVz1530 = new TH2F("hEvPlaneTPCvsEvPVz1530","hEvPlaneTPCvsEvPVz1530",100,-TMath::Pi(),TMath::Pi(),100,-TMath::Pi(),TMath::Pi());
 349   hEvPlaneTPCvsEvPVz3050 = new TH2F("hEvPlaneTPCvsEvPVz3050","hEvPlaneTPCvsEvPVz3050",100,-TMath::Pi(),TMath::Pi(),100,-TMath::Pi(),TMath::Pi());
 350   hEvPlaneTPCvsEvPVz2040 = new TH2F("hEvPlaneTPCvsEvPVz2040","hEvPlaneTPCvsEvPVz2040",100,-TMath::Pi(),TMath::Pi(),100,-TMath::Pi(),TMath::Pi());
 351   hEvPlaneTPCvsEvPVz4060 = new TH2F("hEvPlaneTPCvsEvPVz4060","hEvPlaneTPCvsEvPVz4060",100,-TMath::Pi(),TMath::Pi(),100,-TMath::Pi(),TMath::Pi());
 352
 353   fListHist->Add(hEvPlaneTPCvsEvPVz05);
 354   fListHist->Add(hEvPlaneTPCvsEvPVz075);
 355   fListHist->Add(hEvPlaneTPCvsEvPVz1530);
 356   fListHist->Add(hEvPlaneTPCvsEvPVz3050);
 357   fListHist->Add(hEvPlaneTPCvsEvPVz2040);
 358   fListHist->Add(hEvPlaneTPCvsEvPVz4060);
 359
 360   hCos2DeltaTPCVzAvsCentrality   = new TH2F("hCos2DeltaTPCVzAvsCentrality"  ,"hCos2DeltaTPCVzAvsCentrality"  ,100,-1.1,1.1,105,-0.5,105.5);
 361   hCos2DeltaTPCVzCvsCentrality   = new TH2F("hCos2DeltaTPCVzCvsCentrality"  ,"hCos2DeltaTPCVzCvsCentrality"  ,100,-1.1,1.1,105,-0.5,105.5);
 362   hCos2DeltaVzAVzCvsCentrality   = new TH2F("hCos2DeltaVzAVzCvsCentrality"  ,"hCos2DeltaVzAVzCvsCentrality"  ,100,-1.1,1.1,105,-0.5,105.5);
 363   hCos2DeltaVzMVzAvsCentrality   = new TH2F("hCos2DeltaVzMVzAvsCentrality"  ,"hCos2DeltaVzMVzAvsCentrality"  ,100,-1.1,1.1,105,-0.5,105.5);
 364   hCos2DeltaVzMVzCvsCentrality   = new TH2F("hCos2DeltaVzMVzCvsCentrality"  ,"hCos2DeltaVzMVzCvsCentrality"  ,100,-1.1,1.1,105,-0.5,105.5);
 365   hCos2DeltaVzATPCvsCentrality   = new TH2F("hCos2DeltaVzATPCvsCentrality"  ,"hCos2DeltaVzATPCvsCentrality"  ,100,-1.1,1.1,105,-0.5,105.5);
 366   hCos2DeltaVzCTPCvsCentrality   = new TH2F("hCos2DeltaVzCTPCvsCentrality"  ,"hCos2DeltaVzCTPCvsCentrality"  ,100,-1.1,1.1,105,-0.5,105.5);
 367   hCos2DeltaVzCVzAvsCentrality   = new TH2F("hCos2DeltaVzCVzAvsCentrality"  ,"hCos2DeltaVzCVzAvsCentrality"  ,100,-1.1,1.1,105,-0.5,105.5);
 368   hCos2DeltaVzMTPCpvsCentrality  = new TH2F("hCos2DeltaVzMTPCpvsCentrality" ,"hCos2DeltaVzMTPCpvsCentrality" ,100,-1.1,1.1,105,-0.5,105.5);
 369   hCos2DeltaVzMTPCnvsCentrality  = new TH2F("hCos2DeltaVzMTPCnvsCentrality" ,"hCos2DeltaVzMTPCnvsCentrality" ,100,-1.1,1.1,105,-0.5,105.5);
 370   hCos2DeltaTPCpTPCnvsCentrality = new TH2F("hCos2DeltaTPCpTPCnvsCentrality","hCos2DeltaTPCpTPCnvsCentrality",100,-1.1,1.1,105,-0.5,105.5);
 371
 372   fListHist->Add(hCos2DeltaTPCVzAvsCentrality);
 373   fListHist->Add(hCos2DeltaTPCVzCvsCentrality);
 374   fListHist->Add(hCos2DeltaVzAVzCvsCentrality);
 375   fListHist->Add(hCos2DeltaVzMVzAvsCentrality);
 376   fListHist->Add(hCos2DeltaVzMVzCvsCentrality);
 377   fListHist->Add(hCos2DeltaVzATPCvsCentrality);
 378   fListHist->Add(hCos2DeltaVzCTPCvsCentrality);
 379   fListHist->Add(hCos2DeltaVzCVzAvsCentrality);
 380   fListHist->Add(hCos2DeltaVzMTPCpvsCentrality);
 381   fListHist->Add(hCos2DeltaVzMTPCnvsCentrality);
 382   fListHist->Add(hCos2DeltaTPCpTPCnvsCentrality);
 383
 384
 385   hQVzAQVzCvsCentrality = new TH2F("hQVzAQVzCvsCentrality","hQVzAQVzCvsCentrality",1000,-5,5,105,-0.5,105.5);
 386   fListHist->Add(hQVzAQVzCvsCentrality);
 387
 388
 389   if(!ftree){
 390
 391     ftree = new TTree("ftree","ftree");
 392
 393     ftree->Branch("tCentrality"      ,&tCentrality      ,"tCentrality/D"    );
 394     ftree->Branch("tType"            ,&tType            ,"tType/D"          );
 395     ftree->Branch("tHasTOF"          ,&tHasTOF          ,"tHasTOF/D"        );
 396     ftree->Branch("tpT"              ,&tpT              ,"tpT/D"            );
 397     ftree->Branch("tMassTOF"         ,&tMassTOF         ,"tMassTOF/D"       );
 398     ftree->Branch("tuqV0A"           ,&tuqV0A           ,"tuqV0A/D"         );
 399     ftree->Branch("tuqV0C"           ,&tuqV0C           ,"tuqV0C/D"         );
 400     ftree->Branch("tCharge"          ,&tCharge          ,"tCharge/D"        );
 401     ftree->Branch("tCosdeltaphiTPC"  ,&tCosdeltaphiTPC  ,"tCosdeltaphiTPC/D");
 402     ftree->Branch("tCosdeltaphiV0M"  ,&tCosdeltaphiV0M  ,"tCosdeltaphiV0M/D");
 403     ftree->Branch("tCosdeltaphiV0A"  ,&tCosdeltaphiV0A  ,"tCosdeltaphiV0A/D");
 404     ftree->Branch("tCosdeltaphiV0C"  ,&tCosdeltaphiV0C  ,"tCosdeltaphiV0C/D");
 405     ftree->Branch("timpactXY"        ,&timpactXY        ,"timpactXY/D"      );
 406     ftree->Branch("timpactZ"         ,&timpactZ         ,"timpactZ/D"       );
 407     ftree->Branch("tpull"            ,&tpull            ,"tpull/D"          );
 408     ftree->Branch("tphi"             ,&tphi             ,"tphi/D"          );
 409
 410   }
 411
 412   PostData(1,  fListHist);
 413   PostData(2,  ftree);
 414 }// end UserCreateOutputObjects
 415
 416
 417 //====================== USER EXEC ========================
 418
 419 void AliAnalysisTaskNucleiv2SP::UserExec(Option_t *)
 420 {
 421   // Main loop
 422   // Called for EACH event
 423   //  cout<<"AliAnalysisTaskNucleiv2SP Starting UserExec"<<endl;
 424
 425   Info("AliAnalysisTaskNucleiv2SP","Starting UserExec");
 426
 427   AliVEvent *event = InputEvent();
 428   if (!event) { Printf("ERROR: Could not retrieve event"); return; }
 429
 430
 431   AliESDEvent* lESDevent = dynamic_cast<AliESDEvent*>(event);
 432   if (!lESDevent) {
 433     AliError("Cannot get the ESD event");
 434     return;
 435   }
 436
 437   fHistEventMultiplicity->Fill(1);
 438   fHistEventMultiplicity->Fill(7);
 439
 440   //_____________________________________________________
 441   //   Centrality
 442
 443   AliCentrality *centrality = lESDevent->GetCentrality();
 444   Float_t percentile=centrality->GetCentralityPercentile("V0M");
 445
 446   Int_t TrackNumber = lESDevent->GetNumberOfTracks();
 447   fHistTrackMultiplicity->Fill(TrackNumber,percentile); //tracce per evento
 448
 449   //______________________________________________________
 450   // PID
 451
 452   AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
 453   AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
 454   fPIDResponse=inputHandler->GetPIDResponse();
 455
 456   //=================================================================
 457
 458   Float_t  impactXY=-999., impactZ=-999.;
 459   Double_t TPCSignal=0.;
 460
 461   ULong_t  status=0;
 462
 463   Double_t pmax  = 10;
 464   Double_t ptmax = 7;
 465   // Primary vertex cut
 466
 467   const AliESDVertex *vtx = lESDevent->GetPrimaryVertexTracks();
 468
 469   if(vtx->GetNContributors()<1) {
 470
 471     // SPD vertex cut
 472     vtx = lESDevent->GetPrimaryVertexSPD();
 473
 474     if(vtx->GetNContributors()<1) {
 475       Info("AliAnalysisTaskHelium3Pi","No good vertex, skip event");
 476       return; // NO GOOD VERTEX, SKIP EVENT
 477     }
 478   }
 479
 480   fHistEventMultiplicity->Fill(2); // analyzed events with PV
 481
 482   if(TMath::Abs(vtx->GetZ())>10) return;
 483   fHistEventMultiplicity->Fill(3);
 484
 485   Bool_t isSelectedCentral     = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kCentral);
 486   Bool_t isSelectedSemiCentral = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kSemiCentral);
 487   Bool_t isSelectedMB          = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB);
 488
 489   fHistTrackMultiplicity->Fill(TrackNumber,percentile);
 490
 491   Int_t eventtype = -999;
 492
 493   //  cout<<"ET 1: "<<eventtype<<endl;
 494
 495   if(isSelectedCentral){
 496     fHistEventMultiplicity->Fill(4);
 497     fHistTrackMultiplicityCentral->Fill(TrackNumber,percentile);
 498     eventtype =1;
 499   }
 500
 501   if(isSelectedSemiCentral){
 502     fHistEventMultiplicity->Fill(5);
 503     fHistTrackMultiplicitySemiCentral->Fill(TrackNumber,percentile);
 504     eventtype =2;
 505   }
 506
 507   if(isSelectedMB){
 508     if(percentile<0)return;
 509     if(percentile>=80)return;
 510     fHistEventMultiplicity->Fill(6);
 511     fHistTrackMultiplicityMB->Fill(TrackNumber,percentile);
 512     eventtype =3;
 513   }
 514
 515   //    cout<<"ET 2: "<<eventtype<<endl;
 516
 517   if(eventtype!=1 && eventtype!=2 && eventtype!=3 )return;
 518
 519   AliEventplane *pl=lESDevent->GetEventplane();
 520
 521
 522   if(!pl ){
 523     AliError("AliAnalysisTaskSENucleiv2SP::UserExec:no eventplane! v2 analysis without eventplane not possible!\n");
 524     fHistEventMultiplicity->Fill(12);
 525   }
 526
 527   //Event plane from FLOW
 528
 529   Double_t qxEPa = 0, qyEPa = 0;
 530   Double_t qxEPc = 0, qyEPc = 0;
 531   Double_t qxEP =  0 , qyEP = 0;
 532
 533   Double_t evPlAngV0A = pl->CalculateVZEROEventPlane(lESDevent, 8, 2, qxEPa, qyEPa);
 534   Double_t evPlAngV0C = pl->CalculateVZEROEventPlane(lESDevent, 9, 2, qxEPc, qyEPc);
 535   Double_t evPlAngV0  = pl->CalculateVZEROEventPlane(lESDevent,10, 2, qxEP,  qyEP);
 536
 537   Double_t Qx2  = 0, Qy2  = 0;
 538   Double_t Qx2p = 0, Qy2p = 0;
 539   Double_t Qx2n = 0, Qy2n = 0;
 540
 541   for (Int_t iT = 0; iT < TrackNumber; iT++){
 542
 543     AliESDtrack* track = lESDevent->GetTrack(iT);
 544
 545     if (!track)
 546       continue;
 547
 548     if ((TMath::Abs(track->Eta()) > 0.8) || (track->Pt() < 0.2) || (track->GetTPCNcls() < 70) || (track->Pt() >= 20.0))
 549       continue;
 550     if(!fESDtrackCutsEP->AcceptTrack(track))
 551       continue;
 552     if(track->Eta()>0 && track->Eta()<0.8){
 553
 554       Qx2p += TMath::Cos(2*track->Phi());
 555       Qy2p += TMath::Sin(2*track->Phi());
 556     }
 557     if(track->Eta()<0 && track->Eta()> -0.8){
 558
 559       Qx2n += TMath::Cos(2*track->Phi());
 560       Qy2n += TMath::Sin(2*track->Phi());
 561     }
 562
 563     if(track->Eta()>0 && track->Eta()<0.8){ //half TPC
 564       Qx2 += TMath::Cos(2*track->Phi());
 565       Qy2 += TMath::Sin(2*track->Phi());
 566     }
 567   }
 568
 569   Double_t evPlAngTPC  = TMath::ATan2(Qy2,  Qx2) /2.;
 570   Double_t evPlAngTPCn = TMath::ATan2(Qy2n, Qx2n)/2.;
 571   Double_t evPlAngTPCp = TMath::ATan2(Qy2p, Qx2p)/2.;
 572
 573   EPVzAvsCentrality  ->Fill(evPlAngV0A  , percentile);
 574   EPVzCvsCentrality  ->Fill(evPlAngV0C  , percentile);
 575   EPVzvsCentrality   ->Fill(evPlAngV0   , percentile);
 576   EPTPCvsCentrality  ->Fill(evPlAngTPC  , percentile);
 577   EPTPCpvsCentrality ->Fill(evPlAngTPCp , percentile);
 578   EPTPCnvsCentrality ->Fill(evPlAngTPCn , percentile);
 579
 580   if(percentile>=0 && percentile<=5)
 581     hEvPlaneTPCvsEvPVz05  ->Fill(evPlAngTPC,evPlAngV0);
 582   if(percentile>=0 && percentile<=7.5)
 583     hEvPlaneTPCvsEvPVz075 ->Fill(evPlAngTPC,evPlAngV0);
 584   if(percentile>=15 && percentile<=30)
 585     hEvPlaneTPCvsEvPVz1530->Fill(evPlAngTPC,evPlAngV0);
 586   if(percentile>=30 && percentile<50)
 587     hEvPlaneTPCvsEvPVz3050->Fill(evPlAngTPC,evPlAngV0);
 588   if(percentile>=20 && percentile<=40)
 589     hEvPlaneTPCvsEvPVz2040->Fill(evPlAngTPC,evPlAngV0);
 590   if(percentile>=40 && percentile<=60)
 591     hEvPlaneTPCvsEvPVz4060->Fill(evPlAngTPC,evPlAngV0);
 592
 593   // For TPC, V0M, V0c and V0A resolution
 594
 595   hCos2DeltaTPCVzAvsCentrality  ->Fill(TMath::Cos(2.*(evPlAngTPC - evPlAngV0A)) , percentile);
 596   hCos2DeltaTPCVzCvsCentrality  ->Fill(TMath::Cos(2.*(evPlAngTPC - evPlAngV0C)) , percentile);
 597   hCos2DeltaVzAVzCvsCentrality  ->Fill(TMath::Cos(2.*(evPlAngV0A - evPlAngV0C)) , percentile);
 598   hCos2DeltaVzMVzAvsCentrality  ->Fill(TMath::Cos(2.*(evPlAngV0  - evPlAngV0A)) , percentile);
 599   hCos2DeltaVzMVzCvsCentrality  ->Fill(TMath::Cos(2.*(evPlAngV0  - evPlAngV0C)) , percentile);
 600   hCos2DeltaVzATPCvsCentrality  ->Fill(TMath::Cos(2.*(evPlAngV0A - evPlAngTPC)) , percentile);
 601   hCos2DeltaVzCTPCvsCentrality  ->Fill(TMath::Cos(2.*(evPlAngV0C - evPlAngTPC)) , percentile);
 602   hCos2DeltaVzCVzAvsCentrality  ->Fill(TMath::Cos(2.*(evPlAngV0C - evPlAngV0A)) , percentile);
 603   hCos2DeltaVzMTPCpvsCentrality ->Fill(TMath::Cos(2.*(evPlAngV0  - evPlAngTPCp)), percentile);
 604   hCos2DeltaVzMTPCnvsCentrality ->Fill(TMath::Cos(2.*(evPlAngV0  - evPlAngTPCn)), percentile);
 605   hCos2DeltaTPCpTPCnvsCentrality->Fill(TMath::Cos(2.*(evPlAngTPCp- evPlAngTPCn)), percentile);
 606
 607   //Scalar Product
 608
 609   Double_t  QV0AQV0C = qxEPa * qxEPc + qyEPa*qyEPc;
 610   hQVzAQVzCvsCentrality->Fill(QV0AQV0C,percentile);
 611
 612   //====================================================================================================================
 613
 614   // To remove auto-correlation
 615   TVector2 *q = 0x0;
 616   q = pl->GetQVector();
 617
 618   Float_t ptcExp  = -999;
 619   Double_t pullTPC = -999;
 620   Float_t deltaphiTPC = -3;
 621   Float_t deltaphiV0  = -3;
 622   Float_t deltaphiV0A = -3;
 623   Float_t deltaphiV0C = -3;
 624
 625   Float_t  uqV0A = -999;
 626   Float_t  uqV0C = -999;
 627
 628   for (Int_t j=0; j<TrackNumber; j++) { //loop on tracks
 629
 630     AliESDtrack *esdtrack=lESDevent->GetTrack(j);
 631     if (!fESDtrackCuts->AcceptTrack(esdtrack)) continue;
 632
 633     status  = (ULong_t)esdtrack->GetStatus();
 634
 635     Bool_t hasTOFout  = status&AliESDtrack::kTOFout;
 636     Bool_t hasTOF     = kFALSE;
 637     if (hasTOFout) hasTOF = kTRUE;
 638     Float_t length = esdtrack->GetIntegratedLength();
 639     if (length < 350.) hasTOF = kFALSE;
 640
 641     TPCSignal=esdtrack->GetTPCsignal();
 642
 643     if(TPCSignal<10)continue;
 644     if(TPCSignal>1000)continue;
 645     if(!esdtrack->GetInnerParam()) continue;
 646
 647     Double_t ptot = esdtrack->GetInnerParam()->GetP(); // momentum for dEdx determination
 648
 649     if(ptot<0.2)continue;
 650     fhBB->Fill(ptot*esdtrack->GetSign(),TPCSignal);
 651     esdtrack->GetImpactParameters(impactXY, impactZ);
 652
 653     ptcExp = -999;
 654     if(fptc==1)
 655       ptcExp  = AliExternalTrackParam::BetheBlochAleph(ptot/(0.938*2),1.45802,27.4992,4.00313e-15,2.48485,8.31768);
 656     if(fptc==2)
 657       ptcExp  = AliExternalTrackParam::BetheBlochAleph(ptot/(0.938*3),1.45802,27.4992,4.00313e-15,2.48485,8.31768);
 658     if(fptc==3)
 659       ptcExp  = 4*AliExternalTrackParam::BetheBlochAleph(2*ptot/(0.938*3),1.74962,27.4992,4.00313e-15,2.42485,8.31768);
 660
 661     pullTPC  = (TPCSignal - ptcExp)/(0.07*ptcExp);
 662
 663     Double_t p    = esdtrack->P();
 664     Double_t tof  = esdtrack->GetTOFsignal()-fPIDResponse->GetTOFResponse().GetStartTime(p);
 665     Double_t tPhi = esdtrack->Phi();
 666
 667     Float_t  beta = 0;
 668     Float_t  gamma = 0;
 669     Float_t  mass  = -99;
 670     Double_t pt  = esdtrack->Pt();
 671
 672     if(fptc==3)
 673       pt = 2*pt;
 674
 675     if(TMath::Abs(ptot) < pmax && TMath::Abs(pullTPC) <= 3 && TMath::Abs(pt) < ptmax){
 676
 677       fhBBDeu->Fill(ptot*esdtrack->GetSign(),TPCSignal);
 678
 679       //
 680       // Process TOF information
 681       //
 682       if (hasTOF) {
 683         beta = length / (2.99792457999999984e-02 * tof);
 684         gamma = 1/TMath::Sqrt(1 - beta*beta);
 685         mass = ptot/TMath::Sqrt(gamma*gamma - 1); // using inner TPC mom. as approx.
 686
 687         fhTOF->Fill(ptot*esdtrack->GetSign(),beta);
 688         if(fptc==1){
 689           if(TMath::Abs(mass) > 2.7)continue;
 690           if(TMath::Abs(mass) < 1. )continue;
 691         }
 692         if(fptc==2){
 693           if(TMath::Abs(mass) > 5.0)continue;
 694           if(TMath::Abs(mass) < 1.8 )continue;
 695         }
 696         if(fptc==3){
 697           if(TMath::Abs(mass) > 5.0)continue;
 698           if(TMath::Abs(mass) < 1.8 )continue;
 699         }
 700         fhMassTOF->Fill(mass);
 701       }
 702
 703       // Event Plane
 704       // Remove AutoCorrelation
 705
 706       evPlAngTPC = GetEventPlaneForCandidate(esdtrack,q,pl);
 707
 708       deltaphiTPC=TMath::Cos(2*GetPhi0Pi(tPhi-evPlAngTPC));
 709       deltaphiV0 =TMath::Cos(2*GetPhi0Pi(tPhi-evPlAngV0 ));
 710       deltaphiV0A=TMath::Cos(2*GetPhi0Pi(tPhi-evPlAngV0A));
 711       deltaphiV0C=TMath::Cos(2*GetPhi0Pi(tPhi-evPlAngV0C));
 712
 713       // Scalar Product
 714
 715       uqV0A = TMath::Cos(2*tPhi)*qxEPa+TMath::Sin(2*tPhi)*qyEPa;
 716       uqV0C = TMath::Cos(2*tPhi)*qxEPc+TMath::Sin(2*tPhi)*qyEPc;
 717
 718       tCentrality      = percentile;
 719       tType            = eventtype;
 720       tHasTOF          = hasTOF;
 721       tpT              = pt;
 722       tMassTOF         = mass;
 723       tuqV0A           = uqV0A;
 724       tuqV0C           = uqV0C;
 725       tCharge          = esdtrack->GetSign();
 726       tCosdeltaphiTPC  = deltaphiTPC;
 727       tCosdeltaphiV0M  = deltaphiV0;
 728       tCosdeltaphiV0A  = deltaphiV0A;
 729       tCosdeltaphiV0C  = deltaphiV0C;
 730       timpactXY        = impactXY;
 731       timpactZ         = impactZ;
 732       tpull            = pullTPC;
 733       tphi             = tPhi;
 734
 735       ftree->Fill();
 736     }
 737   }  //track
 738
 739   PostData(1, fListHist);
 740   PostData(2, ftree);
 741 } //end userexec
 742
 743
 744 //________________________________________________________________________
 745
 746 void AliAnalysisTaskNucleiv2SP::Terminate(Option_t *)
 747 {
 748   // Draw result to the screen
 749   // Called once at the end of the query
 750 }
 751