PWGLF/STRANGENESS/Hypernuclei/AliAnalysisTaskHelium3Pi.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 notice 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 //-----------------------------------------------------------------
  14 //                 AliAnalysisTaskHelium3Pi class
  15 //-----------------------------------------------------------------
  16
  17
  18 class TTree;
  19 class TParticle;
  20 class TVector3;
  21
  22 #include "AliAnalysisManager.h"
  23 #include <AliMCEventHandler.h>
  24 #include <AliMCEvent.h>
  25 #include <AliStack.h>
  26
  27 class AliESDVertex;
  28 class AliAODVertex;
  29 class AliESDv0;
  30 class AliAODv0;
  31 class AliCascadeVertexer;
  32
  33 #include <iostream>
  34 #include "AliAnalysisTaskSE.h"
  35 #include "TList.h"
  36 #include "TH1.h"
  37 #include "TH2.h"
  38 #include "TH3.h"
  39 #include "TNtuple.h"
  40 #include "TGraph.h"
  41 #include "TCutG.h"
  42 #include "TF1.h"
  43 #include "TCanvas.h"
  44 #include "TMath.h"
  45 #include "TChain.h"
  46 #include "Riostream.h"
  47 #include "AliLog.h"
  48 #include "AliCascadeVertexer.h"
  49 #include "AliESDEvent.h"
  50 #include "AliESDtrack.h"
  51 #include "AliExternalTrackParam.h"
  52 #include "AliAODEvent.h"
  53 #include "AliInputEventHandler.h"
  54 #include "AliESDcascade.h"
  55 #include "AliAODcascade.h"
  56 #include "AliAnalysisTaskHelium3Pi.h"
  57 #include "AliESDtrackCuts.h"
  58 #include "AliCentrality.h"
  59 #include "TString.h"
  60 #include <TDatime.h>
  61 #include <TRandom3.h>
  62 #include <TLorentzVector.h>
  63 //#include <AliVTrack.h>
  64
  65 const Int_t AliAnalysisTaskHelium3Pi::fgNrot = 15;
  66
  67 ClassImp(AliAnalysisTaskHelium3Pi)
  68
  69 //________________________________________________________________________
  70 AliAnalysisTaskHelium3Pi::AliAnalysisTaskHelium3Pi()
  71 : AliAnalysisTaskSE(),
  72   fAnalysisType("ESD"),
  73   fCollidingSystems(0),
  74   fESDtrackCuts(0),
  75   fDataType("REAL"),
  76   fListHist(0),
  77   fHistEventMultiplicity(0),
  78   fHistTrackMultiplicity(0),
  79   fHistTrackMultiplicityCent(0),
  80   fHistTrackMultiplicitySemiCent(0),
  81   fHistTrackMultiplicityMB(0),
  82   fHistTrackMultiplicityPVCent(0),
  83   fHistTrackMultiplicityPVSemiCent(0),
  84   fHistTrackMultiplicityPVMB(0),
  85   fhBB(0),
  86   fhTOF(0),
  87   fhMassTOF(0),
  88   fhBBPions(0),
  89   fhBBHe(0),
  90   fhNaPos(0),
  91   fhNaNeg(0),
  92   fBetavsTPCsignalPos(0),
  93   fBetavsTPCsignalNeg(0),
  94   fNtuple1(0),
  95   trunNumber(0),
  96   tbunchcross(0),
  97   torbit(0),
  98   tperiod(0),
  99   teventtype(0),
 100   tTrackNumber(0),
 101   tpercentile(0),
 102   txPrimaryVertex(0),
 103   tyPrimaryVertex(0),
 104   tzPrimaryVertex(0),
 105   txSecondaryVertex(0),
 106   tySecondaryVertex(0),
 107   tzSecondaryVertex(0),
 108   tdcaTracks(0),
 109   tCosPointingAngle(0),
 110   tDCAV0toPrimaryVertex(0),
 111   tHeSign(0),
 112   tHepInTPC(0),
 113   tHeTPCsignal(0),
 114   tDcaHeToPrimVertex(0),
 115   tHeEta(0),
 116   tmomHex(0),
 117   tmomHey(0),
 118   tmomHez(0),
 119   tmomHeAtSVx(0),
 120   tmomHeAtSVy(0),
 121   tmomHeAtSVz(0),
 122   tHeTPCNcls(0),
 123   tHeimpactXY(0),
 124   tHeimpactZ(0),
 125   tHeITSClusterMap(0),
 126   tIsHeITSRefit(0),
 127   tPionSign(0),
 128   tPionpInTPC(0),
 129   tPionTPCsignal(0),
 130   tDcaPionToPrimVertex(0),
 131   tPionEta(0),
 132   tmomPionx(0),
 133   tmomPiony(0),
 134   tmomPionz(0),
 135   tmomNegPionAtSVx(0),
 136   tmomNegPionAtSVy(0),
 137   tmomNegPionAtSVz(0),
 138   tPionTPCNcls(0),
 139   tPionimpactXY(0),
 140   tPionimpactZ(0),
 141   tPionITSClusterMap(0),
 142   tIsPiITSRefit(0),
 143   txn(0),
 144   txp(0),
 145   tchi2He(0),
 146   tchi2Pi(0),
 147   fNtuple4(0),
 148   tHelrunNumber(0),
 149   tHelBCNumber(0),
 150   tHelOrbitNumber(0),
 151   tHelPeriodNumber(0),
 152   tHeleventtype(0),
 153   tHelisHeITSrefit(0),
 154   tHelpercentile(0),
 155   tHelSign(0),
 156   tHelpinTPC(0),
 157   tHelGetTPCsignal(0),
 158   tHelPx(0),
 159   tHelPy(0),
 160   tHelPz(0),
 161   tHelEta(0),
 162   tHelisTOF(0),
 163   tHelpoutTPC(0),
 164   tHeltimeTOF(0),
 165   tHeltrackLenghtTOF(0),
 166   tHelimpactXY(0),
 167   tHelimpactZ(0),
 168   tHelmapITS(0),
 169   tHelTPCNcls(0),
 170   tHelTRDsignal(0),
 171   tHelxPrimaryVertex(0),
 172   tHelyPrimaryVertex(0),
 173   tHelzPrimaryVertex(0),
 174   tHelchi2PerClusterTPC(0),
 175   fPIDResponse(0)
 176
 177 {
 178   // Dummy Constructor
 179
 180   //  printf("Dummy Constructor");
 181
 182   fESDtrackCuts = new AliESDtrackCuts("fESDtrackCuts");
 183   fESDtrackCuts->SetRequireITSStandAlone(kFALSE);
 184   fESDtrackCuts->SetRequireITSPureStandAlone(kFALSE);
 185
 186   fESDtrackCuts->SetRequireTPCRefit(kTRUE);
 187   fESDtrackCuts->SetAcceptKinkDaughters(kFALSE);
 188   fESDtrackCuts->SetMinNClustersTPC(60);
 189   fESDtrackCuts->SetMaxChi2PerClusterTPC(5);
 190   fESDtrackCuts->SetEtaRange(-0.9,0.9);
 191
 192 }
 193
 194 //________________________________________________________________________
 195 AliAnalysisTaskHelium3Pi::AliAnalysisTaskHelium3Pi(const char *name)
 196 : AliAnalysisTaskSE(name),
 197   fAnalysisType("ESD"),
 198   fCollidingSystems(0),
 199   fESDtrackCuts(0),
 200   fDataType("REAL"),
 201   fListHist(0),
 202   fHistEventMultiplicity(0),
 203   fHistTrackMultiplicity(0),
 204   fHistTrackMultiplicityCent(0),
 205   fHistTrackMultiplicitySemiCent(0),
 206   fHistTrackMultiplicityMB(0),
 207   fHistTrackMultiplicityPVCent(0),
 208   fHistTrackMultiplicityPVSemiCent(0),
 209   fHistTrackMultiplicityPVMB(0),
 210   fhBB(0),
 211   fhTOF(0),
 212   fhMassTOF(0),
 213   fhBBPions(0),
 214   fhBBHe(0),
 215   fhNaPos(0),
 216   fhNaNeg(0),
 217   fBetavsTPCsignalPos(0),
 218   fBetavsTPCsignalNeg(0),
 219   fNtuple1(0),
 220   trunNumber(0),
 221   tbunchcross(0),
 222   torbit(0),
 223   tperiod(0),
 224   teventtype(0),
 225   tTrackNumber(0),
 226   tpercentile(0),
 227   txPrimaryVertex(0),
 228   tyPrimaryVertex(0),
 229   tzPrimaryVertex(0),
 230   txSecondaryVertex(0),
 231   tySecondaryVertex(0),
 232   tzSecondaryVertex(0),
 233   tdcaTracks(0),
 234   tCosPointingAngle(0),
 235   tDCAV0toPrimaryVertex(0),
 236   tHeSign(0),
 237   tHepInTPC(0),
 238   tHeTPCsignal(0),
 239   tDcaHeToPrimVertex(0),
 240   tHeEta(0),
 241   tmomHex(0),
 242   tmomHey(0),
 243   tmomHez(0),
 244   tmomHeAtSVx(0),
 245   tmomHeAtSVy(0),
 246   tmomHeAtSVz(0),
 247   tHeTPCNcls(0),
 248   tHeimpactXY(0),
 249   tHeimpactZ(0),
 250   tHeITSClusterMap(0),
 251   tIsHeITSRefit(0),
 252   tPionSign(0),
 253   tPionpInTPC(0),
 254   tPionTPCsignal(0),
 255   tDcaPionToPrimVertex(0),
 256   tPionEta(0),
 257   tmomPionx(0),
 258   tmomPiony(0),
 259   tmomPionz(0),
 260   tmomNegPionAtSVx(0),
 261   tmomNegPionAtSVy(0),
 262   tmomNegPionAtSVz(0),
 263   tPionTPCNcls(0),
 264   tPionimpactXY(0),
 265   tPionimpactZ(0),
 266   tPionITSClusterMap(0),
 267   tIsPiITSRefit(0),
 268   txn(0),
 269   txp(0),
 270   tchi2He(0),
 271   tchi2Pi(0),
 272   fNtuple4(0),
 273   tHelrunNumber(0),
 274   tHelBCNumber(0),
 275   tHelOrbitNumber(0),
 276   tHelPeriodNumber(0),
 277   tHeleventtype(0),
 278   tHelisHeITSrefit(0),
 279   tHelpercentile(0),
 280   tHelSign(0),
 281   tHelpinTPC(0),
 282   tHelGetTPCsignal(0),
 283   tHelPx(0),
 284   tHelPy(0),
 285   tHelPz(0),
 286   tHelEta(0),
 287   tHelisTOF(0),
 288   tHelpoutTPC(0),
 289   tHeltimeTOF(0),
 290   tHeltrackLenghtTOF(0),
 291   tHelimpactXY(0),
 292   tHelimpactZ(0),
 293   tHelmapITS(0),
 294   tHelTPCNcls(0),
 295   tHelTRDsignal(0),
 296   tHelxPrimaryVertex(0),
 297   tHelyPrimaryVertex(0),
 298   tHelzPrimaryVertex(0),
 299   tHelchi2PerClusterTPC(0),
 300   fPIDResponse(0)
 301 {
 302
 303   // Define input and output slots here
 304   // Input slot #0 works with a TChain
 305   //DefineInput(0, TChain::Class());
 306   // Output slot #0 writes into a TList container ()
 307
 308   DefineInput(0, TChain::Class());
 309
 310   DefineOutput(1, TList::Class());
 311   DefineOutput(2, TTree::Class());
 312   DefineOutput(3, TTree::Class());
 313
 314   fESDtrackCuts = new AliESDtrackCuts("fESDtrackCuts");
 315   fESDtrackCuts->SetRequireITSStandAlone(kFALSE);
 316   fESDtrackCuts->SetRequireITSPureStandAlone(kFALSE);
 317
 318   fESDtrackCuts->SetRequireTPCRefit(kTRUE);
 319   fESDtrackCuts->SetAcceptKinkDaughters(kFALSE);
 320   fESDtrackCuts->SetMinNClustersTPC(60);
 321   fESDtrackCuts->SetMaxChi2PerClusterTPC(5);
 322   fESDtrackCuts->SetEtaRange(-0.9,0.9);
 323
 324 }
 325 //_______________________________________________________
 326 AliAnalysisTaskHelium3Pi::~AliAnalysisTaskHelium3Pi()
 327 {
 328   // Destructor
 329   if (fListHist) {
 330     delete fListHist;
 331     fListHist = 0;
 332   }
 333
 334   if (fESDtrackCuts) delete fESDtrackCuts;
 335   if(fNtuple1) delete fNtuple1;
 336   if(fNtuple4) delete fNtuple4;
 337 }
 338 //=================DEFINITION BETHE BLOCH==============================
 339
 340 Double_t AliAnalysisTaskHelium3Pi::BetheBloch(Double_t betaGamma,Double_t charge,Bool_t isPbPb) {
 341
 342   Double_t kp1, kp2, kp3, kp4, kp5;
 343
 344   if(isPbPb){
 345
 346     //    pass1 2011
 347     kp1 = 4.7*charge*charge;
 348     kp2 = 8.98482806165147636e+00;
 349     kp3 = 1.54000000000000005e-05;
 350     kp4 = 2.30445734159456084e+00;
 351     kp5 = 2.25624744086878559e+00;
 352
 353   }
 354
 355   else{
 356
 357     // to be defined ...
 358     //pass1 2011
 359     kp1 = 4.7*charge*charge;
 360     kp2 = 8.98482806165147636e+00;
 361     kp3 = 1.54000000000000005e-05;
 362     kp4 = 2.30445734159456084e+00;
 363     kp5 = 2.25624744086878559e+00;
 364
 365   }
 366
 367   Double_t beta = betaGamma / TMath::Sqrt(1.0 + betaGamma * betaGamma);
 368
 369   Double_t aa = TMath::Power(beta, kp4);
 370   Double_t bb = TMath::Power(1.0 / betaGamma, kp5);
 371
 372   bb = TMath::Log(kp3 + bb);
 373
 374   Double_t out = (kp2 - aa - bb) * kp1 / aa;
 375
 376   return out;
 377
 378 }
 379
 380 //==================DEFINITION OF OUTPUT OBJECTS==============================
 381
 382 void AliAnalysisTaskHelium3Pi::UserCreateOutputObjects()
 383 {
 384
 385   fListHist = new TList();
 386   fListHist->SetOwner();  // IMPORTANT!
 387
 388   if(! fHistEventMultiplicity ){
 389     fHistEventMultiplicity   = new TH1F( "fHistEventMultiplicity" , "Nb of Events" , 12 , -0.5, 11.5);
 390     fHistEventMultiplicity->GetXaxis()->SetBinLabel(1 ,"All Events");
 391     fHistEventMultiplicity->GetXaxis()->SetBinLabel(2 ,"Events w/PV");
 392     fHistEventMultiplicity->GetXaxis()->SetBinLabel(3 ,"Events w/|Vz|<10cm");
 393     fHistEventMultiplicity->GetXaxis()->SetBinLabel(4 ,"Central Events");
 394     fHistEventMultiplicity->GetXaxis()->SetBinLabel(5 ,"SemiCentral Events");
 395     fHistEventMultiplicity->GetXaxis()->SetBinLabel(6 ,"MB Events");
 396     fHistEventMultiplicity->GetXaxis()->SetBinLabel(7 ,"Central Events  w/|Vz|<10cm");
 397     fHistEventMultiplicity->GetXaxis()->SetBinLabel(8 ,"SemiCentral Events  w/|Vz|<10cm");
 398     fHistEventMultiplicity->GetXaxis()->SetBinLabel(9 ,"MB Events w/|Vz|<10cm");
 399     fHistEventMultiplicity->GetXaxis()->SetBinLabel(10,"Any Events");
 400     fHistEventMultiplicity->GetXaxis()->SetBinLabel(11,"Any Events w/|Vz|<10cm");
 401
 402     fListHist->Add(fHistEventMultiplicity);
 403   }
 404
 405   if(! fHistTrackMultiplicity ){
 406     fHistTrackMultiplicity   = new TH2F( "fHistTrackMultiplicity" , "Nb of Tracks", 2500,0, 25000,210,-1,104);
 407     fHistTrackMultiplicity->GetXaxis()->SetTitle("Number of tracks");
 408     fHistTrackMultiplicity->GetYaxis()->SetTitle("Percentile");
 409     fListHist->Add(fHistTrackMultiplicity);
 410   }
 411
 412   if(! fHistTrackMultiplicityCent ){
 413     fHistTrackMultiplicityCent   = new TH2F( "fHistTrackMultiplicityCent", "Nb of Tracks Central Events", 2500,0, 25000,210,-1,104 );
 414     fHistTrackMultiplicityCent->GetXaxis()->SetTitle("Number of tracks");
 415     fHistTrackMultiplicityCent->GetYaxis()->SetTitle("Percentile");
 416     fListHist->Add(fHistTrackMultiplicityCent);
 417   }
 418
 419   if(! fHistTrackMultiplicitySemiCent ){
 420     fHistTrackMultiplicitySemiCent   = new TH2F( "fHistTrackMultiplicitySemiCent" , "Nb of Tracks SemiCentral Events", 2500,0, 25000 ,210,-1,104);
 421     fHistTrackMultiplicitySemiCent->GetXaxis()->SetTitle("Number of tracks");
 422     fHistTrackMultiplicitySemiCent->GetYaxis()->SetTitle("Percentile");
 423     fListHist->Add(fHistTrackMultiplicitySemiCent);
 424   }
 425
 426   if(! fHistTrackMultiplicityMB ){
 427     fHistTrackMultiplicityMB   = new TH2F( "fHistTrackMultiplicityMB" , "Nb of Tracks MBral Events", 2500,0, 25000,210,-1,104 );
 428     fHistTrackMultiplicityMB->GetXaxis()->SetTitle("Number of tracks");
 429     fHistTrackMultiplicityMB->GetYaxis()->SetTitle("Percentile");
 430     fListHist->Add(fHistTrackMultiplicityMB);
 431   }
 432
 433   if(! fHistTrackMultiplicityPVCent ){
 434     fHistTrackMultiplicityPVCent   = new TH2F( "fHistTrackMultiplicityPVCent" , "Nb of Tracks Central Events", 2500,0, 25000,210,-1,104 );
 435     fHistTrackMultiplicityPVCent->GetXaxis()->SetTitle("Number of tracks");
 436     fHistTrackMultiplicityPVCent->GetYaxis()->SetTitle("Percentile");
 437     fListHist->Add(fHistTrackMultiplicityPVCent);
 438   }
 439
 440   if(! fHistTrackMultiplicityPVSemiCent ){
 441     fHistTrackMultiplicityPVSemiCent   = new TH2F( "fHistTrackMultiplicityPVSemiCent" , "Nb of Tracks SemiCentral Events", 2500,0, 25000 ,210,-1,104);
 442     fHistTrackMultiplicityPVSemiCent->GetXaxis()->SetTitle("Number of tracks");
 443     fHistTrackMultiplicityPVSemiCent->GetYaxis()->SetTitle("Percentile");
 444     fListHist->Add(fHistTrackMultiplicityPVSemiCent);
 445   }
 446
 447   if(! fHistTrackMultiplicityPVMB ){
 448     fHistTrackMultiplicityPVMB   = new TH2F( "fHistTrackMultiplicityPVMB" , "Nb of Tracks MBral Events", 2500,0, 25000,210,-1,104 );
 449     fHistTrackMultiplicityPVMB->GetXaxis()->SetTitle("Number of tracks");
 450     fHistTrackMultiplicityPVMB->GetYaxis()->SetTitle("Percentile");
 451     fListHist->Add(fHistTrackMultiplicityPVMB);
 452   }
 453
 454   if(! fhBB ){
 455     fhBB = new TH2F( "fhBB" , "BetheBlochTPC" , 120,-6,6,150,0,1500);
 456     fhBB->GetXaxis()->SetTitle("p/z (GeV/#it{c})");
 457     fhBB->GetYaxis()->SetTitle("TPC Signal");
 458     fListHist->Add(fhBB);
 459   }
 460
 461   if(! fhTOF ){
 462     fhTOF = new TH2F( "fhTOF" , "Scatter Plot TOF" , 120,-6,6,100,0,1.2);
 463     fhTOF->GetXaxis()->SetTitle("p/z (GeV/#it{c})");
 464     fhTOF->GetYaxis()->SetTitle("#beta");
 465     fListHist->Add(fhTOF);
 466   }
 467
 468   if(! fhMassTOF){
 469     fhMassTOF=new TH1F ("fhMassTOF","Particle Mass - TOF", 300,0 ,5);
 470     fhMassTOF->GetXaxis()->SetTitle("Mass (GeV/#it{c}^{2})");
 471     fListHist->Add(fhMassTOF);
 472   }
 473
 474   if(! fhBBPions ){
 475     fhBBPions = new TH2F( "fhBBPions" , "Bethe-Bloch TPC Pions" , 120,-6,6,150,0,1500);
 476     fhBBPions->GetXaxis()->SetTitle("p/z (GeV/#it{c})");
 477     fhBBPions->GetYaxis()->SetTitle("TPC Signal");
 478     fListHist->Add(fhBBPions);
 479   }
 480
 481   if(! fhBBHe ){
 482     fhBBHe = new TH2F( "fhBBHe" , "Bethe-Bloch TPC He" , 120,-6,6,150,0,1500);
 483     fhBBHe->GetXaxis()->SetTitle("p/z (GeV/#it{c})");
 484     fhBBHe->GetYaxis()->SetTitle("TPC Signal");
 485     fListHist->Add(fhBBHe);
 486   }
 487
 488   if(! fhNaPos ){
 489     fhNaPos = new TH2F( "fhNaPos" , "Distribution Pos" , 500,0,5,40,-10,10);
 490     fhNaPos->GetXaxis()->SetTitle("p/z (GeV/#it{c})");
 491     fhNaPos->GetYaxis()->SetTitle("(TPCSignal-bbtheo)/bbtheo (He)");
 492     fListHist->Add(fhNaPos);
 493   }
 494
 495   if(! fhNaNeg ){
 496     fhNaNeg = new TH2F( "fhNaNeg" , "Distribution Neg" , 500,0,5,40,-10,10);
 497     fhNaNeg->GetXaxis()->SetTitle("p/z (GeV/#it{c})");
 498     fhNaNeg->GetYaxis()->SetTitle("(TPCSignal-bbtheo)/bbtheo (He)");
 499     fListHist->Add(fhNaNeg);
 500   }
 501
 502   if(! fBetavsTPCsignalPos ){
 503     fBetavsTPCsignalPos = new TH2F("fBetavsTPCsignalPos","fBetavsTPCsignalPos",100,0,1.2,150,0,1500);
 504     fBetavsTPCsignalPos->GetXaxis()->SetTitle("#beta");
 505     fBetavsTPCsignalPos->GetYaxis()->SetTitle("TPC Signal");
 506     fListHist->Add(fBetavsTPCsignalPos);
 507   }
 508
 509   if(! fBetavsTPCsignalNeg ){
 510     fBetavsTPCsignalNeg = new TH2F("fBetavsTPCsignalNeg","fBetavsTPCsignalNeg",100,0,1.2,150,0,1500);
 511     fBetavsTPCsignalNeg->GetXaxis()->SetTitle("#beta");
 512     fBetavsTPCsignalNeg->GetYaxis()->SetTitle("TPC Signal");
 513     fListHist->Add(fBetavsTPCsignalNeg);
 514   }
 515
 516
 517
 518   if(! fNtuple1 ) {
 519
 520     fNtuple1 = new TTree("fNtuple1","fNtuple1");
 521
 522     fNtuple1->Branch("trunNumber"           ,&trunNumber           ,"trunNumber/F");
 523     fNtuple1->Branch("tbunchcross"          ,&tbunchcross          ,"tbunchcross/F");
 524     fNtuple1->Branch("torbit"               ,&torbit               ,"torbit/F");
 525     fNtuple1->Branch("tperiod"              ,&tperiod              ,"tperiod/F");
 526     fNtuple1->Branch("teventtype"           ,&teventtype           ,"teventtype/F");
 527     fNtuple1->Branch("tTrackNumber"         ,&tTrackNumber         ,"tTrackNumber/F");
 528     fNtuple1->Branch("tpercentile"          ,&tpercentile          ,"tpercentile/F") ;
 529     fNtuple1->Branch("txPrimaryVertex"      ,&txPrimaryVertex      ,"txPrimaryVertex/F");
 530     fNtuple1->Branch("tyPrimaryVertex"      ,&tyPrimaryVertex      ,"tyPrimaryVertex/F");
 531     fNtuple1->Branch("tzPrimaryVertex"      ,&tzPrimaryVertex      ,"tzPrimaryVertex/F");
 532     fNtuple1->Branch("txSecondaryVertex"    ,&txSecondaryVertex    ,"txSecondaryVertex/F");
 533     fNtuple1->Branch("tySecondaryVertex"    ,&tySecondaryVertex    ,"tySecondaryVertex/F");
 534     fNtuple1->Branch("tzSecondaryVertex"    ,&tzSecondaryVertex    ,"tzSecondaryVertex/F");
 535     fNtuple1->Branch("tdcaTracks"           ,&tdcaTracks           ,"tdcaTracks/F");
 536     fNtuple1->Branch("tCosPointingAngle"    ,&tCosPointingAngle    ,"tCosPointingAngle/F");
 537     fNtuple1->Branch("tDCAV0toPrimaryVertex",&tDCAV0toPrimaryVertex,"tDCAV0toPrimaryVertex/F");
 538     fNtuple1->Branch("tHeSign"              ,&tHeSign              ,"tHeSign/F");
 539     fNtuple1->Branch("tHepInTPC"            ,&tHepInTPC            ,"tHepInTPC/F");
 540     fNtuple1->Branch("tHeTPCsignal"         ,&tHeTPCsignal         ,"tHeTPCsignal/F");
 541     fNtuple1->Branch("tDcaHeToPrimVertex"   ,&tDcaHeToPrimVertex   ,"tDcaHeToPrimVertex/F");
 542     fNtuple1->Branch("tHeEta"               ,&tHeEta               ,"tHeEta/F");
 543     fNtuple1->Branch("tmomHex"              ,&tmomHex              ,"tmomHex/F");
 544     fNtuple1->Branch("tmomHey"              ,&tmomHey              ,"tmomHey/F");
 545     fNtuple1->Branch("tmomHez"              ,&tmomHez              ,"tmomHez/F");
 546     fNtuple1->Branch("tmomHeAtSVx"          ,&tmomHeAtSVx          ,"tmomHeAtSVx/F");
 547     fNtuple1->Branch("tmomHeAtSVy"          ,&tmomHeAtSVy          ,"tmomHeAtSVy/F");
 548     fNtuple1->Branch("tmomHeAtSVz"          ,&tmomHeAtSVz          ,"tmomHeAtSVz/F");
 549     fNtuple1->Branch("tHeTPCNcls"           ,&tHeTPCNcls           ,"tHeTPCNcls/F");
 550     fNtuple1->Branch("tHeimpactXY"          ,&tHeimpactXY          ,"tHeimpactXY/F");
 551     fNtuple1->Branch("tHeimpactZ"           ,&tHeimpactZ           ,"tHeimpactZ/F");
 552     fNtuple1->Branch("tHeITSClusterMap"     ,&tHeITSClusterMap     ,"tHeITSClusterMap/F");
 553     fNtuple1->Branch("tIsHeITSRefit"        ,&tIsHeITSRefit        ,"tIsHeITSRefit/F");
 554     fNtuple1->Branch("tPionSign"            ,&tPionSign            ,"tPionSign/F");
 555     fNtuple1->Branch("tPionpInTPC"          ,&tPionpInTPC          ,"tPionpInTPC/F");
 556     fNtuple1->Branch("tPionTPCsignal"       ,&tPionTPCsignal       ,"tPionTPCsignal/F");
 557     fNtuple1->Branch("tDcaPionToPrimVertex" ,&tDcaPionToPrimVertex ,"tDcaPionToPrimVertex/F");
 558     fNtuple1->Branch("tPionEta"             ,&tPionEta             ,"tPionEta/F");
 559     fNtuple1->Branch("tmomPionx"            ,&tmomPionx            ,"tmomPionx/F");
 560     fNtuple1->Branch("tmomPiony"            ,&tmomPiony            ,"tmomPiony/F");
 561     fNtuple1->Branch("tmomPionz"            ,&tmomPionz            ,"tmomPionz/F");
 562     fNtuple1->Branch("tmomNegPionAtSVx"     ,&tmomNegPionAtSVx     ,"tmomNegPionAtSVx/F");
 563     fNtuple1->Branch("tmomNegPionAtSVy"     ,&tmomNegPionAtSVy     ,"tmomNegPionAtSVy/F");
 564     fNtuple1->Branch("tmomNegPionAtSVz"     ,&tmomNegPionAtSVz     ,"tmomNegPionAtSVz/F");
 565     fNtuple1->Branch("tPionTPCNcls"         ,&tPionTPCNcls         ,"tPionTPCNcls/F");
 566     fNtuple1->Branch("tPionimpactXY"        ,&tPionimpactXY        ,"tPionimpactXY/F");
 567     fNtuple1->Branch("tPionimpactZ"         ,&tPionimpactZ         ,"tPionimpactZ/F");
 568     fNtuple1->Branch("tPionITSClusterMap"   ,&tPionITSClusterMap   ,"tPionITSClusterMap/F");
 569     fNtuple1->Branch("tIsPiITSRefit"        ,&tIsPiITSRefit        ,"tIsPiITSRefit/F");
 570     fNtuple1->Branch("txn"                  ,&txn                  ,"txn/F");
 571     fNtuple1->Branch("txp"                  ,&txp                  ,"txp/F");
 572     fNtuple1->Branch("tchi2He"              ,&tchi2He              ,"tchi2He/F");
 573     fNtuple1->Branch("tchi2Pi"              ,&tchi2Pi              ,"tchi2Pi/F");
 574
 575   }
 576
 577   if(! fNtuple4 ) {
 578
 579     fNtuple4 = new TTree("fNtuple4","fNtuple4");
 580
 581     fNtuple4->Branch("tHelrunNumber"        ,&tHelrunNumber        ,"tHelrunNumber/F");
 582     fNtuple4->Branch("tHelBCNumber"         ,&tHelBCNumber         ,"tHelBCNumber/F");
 583     fNtuple4->Branch("tHelOrbitNumber"      ,&tHelOrbitNumber      ,"tHelOrbitNumber/F");
 584     fNtuple4->Branch("tHelPeriodNumber"     ,&tHelPeriodNumber     ,"tHelPeriodNumber/F");
 585     fNtuple4->Branch("tHeleventtype"        ,&tHeleventtype        ,"tHeleventtype/F");
 586     fNtuple4->Branch("tHelisHeITSrefit"     ,&tHelisHeITSrefit     ,"tHelisHeITSrefit/F");
 587     fNtuple4->Branch("tHelpercentile"       ,&tHelpercentile       ,"tHelpercentile/F");
 588     fNtuple4->Branch("tHelSign"             ,&tHelSign             ,"tHelSign/F");
 589     fNtuple4->Branch("tHelpinTPC"           ,&tHelpinTPC           ,"tHelpinTPC/F");
 590     fNtuple4->Branch("tHelGetTPCsignal"     ,&tHelGetTPCsignal     ,"tHelGetTPCsignal/F");
 591     fNtuple4->Branch("tHelPx"               ,&tHelPx               ,"tHelPx/F");
 592     fNtuple4->Branch("tHelPy"               ,&tHelPy               ,"tHelPy/F");
 593     fNtuple4->Branch("tHelPz"               ,&tHelPz               ,"tHelPz/F");
 594     fNtuple4->Branch("tHelEta"              ,&tHelEta              ,"tHelEta/F");
 595     fNtuple4->Branch("tHelisTOF"            ,&tHelisTOF            ,"tHelisTOF/F");
 596     fNtuple4->Branch("tHelpoutTPC"          ,&tHelpoutTPC          ,"tHelpoutTPC/F");
 597     fNtuple4->Branch("tHeltimeTOF"          ,&tHeltimeTOF          ,"tHeltimeTOF/F");
 598     fNtuple4->Branch("tHeltrackLenghtTOF"   ,&tHeltrackLenghtTOF   ,"tHeltrackLenghtTOF/F");
 599     fNtuple4->Branch("tHelimpactXY"         ,&tHelimpactXY         ,"tHelimpactXY/F");
 600     fNtuple4->Branch("tHelimpactZ"          ,&tHelimpactZ          ,"tHelimpactZ/F");
 601     fNtuple4->Branch("tHelmapITS"           ,&tHelmapITS           ,"tHelmapITS/F");
 602     fNtuple4->Branch("tHelTPCNcls"          ,&tHelTPCNcls          ,"tHelTPCNcls/F");
 603     fNtuple4->Branch("tHelTRDsignal"        ,&tHelTRDsignal        ,"tHelTRDsignal/F");
 604     fNtuple4->Branch("tHelxPrimaryVertex"   ,&tHelxPrimaryVertex   ,"tHelxPrimaryVertex/F");
 605     fNtuple4->Branch("tHelyPrimaryVertex"   ,&tHelyPrimaryVertex   ,"tHelyPrimaryVertex/F");
 606     fNtuple4->Branch("tHelzPrimaryVertex"   ,&tHelzPrimaryVertex   ,"tHelzPrimaryVertex/F");
 607     fNtuple4->Branch("tHelchi2PerClusterTPC",&tHelchi2PerClusterTPC,"tHelchi2PerClusterTPC/F");
 608
 609
 610   }
 611
 612   PostData(1,  fListHist);
 613   PostData(2,  fNtuple1);
 614   PostData(3,  fNtuple4);
 615 }// end UserCreateOutputObjects
 616
 617
 618 //====================== USER EXEC ========================
 619
 620 void AliAnalysisTaskHelium3Pi::UserExec(Option_t *)
 621 {
 622   //_______________________________________________________________________
 623
 624   //!*********************!//
 625   //!  Define variables   !//
 626   //!*********************!//
 627
 628   Double_t pinTPC=0.,poutTPC=0.,TPCSignal=0.;
 629   Double_t xPrimaryVertex=0.,yPrimaryVertex=0.,zPrimaryVertex=0.;
 630   Double_t massTOF=0.,timeTOF=0.,trackLenghtTOF=0.,betaTOF=0.;
 631
 632   ULong_t  status=0;
 633   //  ULong_t  statusT=0;
 634   ULong_t  statusPi=0;
 635
 636   Bool_t   isTPC=kFALSE,isTOF=kFALSE,IsHeITSRefit=kFALSE,IsPiITSRefit=kFALSE ;
 637
 638   Float_t nSigmaNegPion=0.;
 639
 640   Double_t cutNSigma = 3;
 641   Double_t bbtheoM=0.,bbtheo=0.;
 642   Double_t zNathashaNeg=0;
 643   Double_t zNathashaPos=0;
 644   Double_t fPos[3]={0.,0.,0.};
 645   Double_t runNumber=0.;
 646   //  Double_t evNumber=0.;
 647
 648   Double_t BCNumber=0.;
 649   Double_t OrbitNumber=0.;
 650   Double_t PeriodNumber=0.;
 651
 652   Double_t        Helium3Mass = 2.80839;
 653   Double_t        PionMass    = 0.13957;
 654   // TLORENTZ vectors
 655
 656   TLorentzVector  vPion,vHelium,vSum;
 657
 658   //!----------------------------------------------------------------
 659
 660   //! A set of very loose parameters for cuts
 661
 662   Double_t fgChi2max=33.;     //! max chi2
 663   Double_t fgDNmin=0.05;      //! min imp parameter for the 1st daughter = 500um
 664   Double_t fgDCAmax=1.0;      //! max DCA between the daughter tracks in cm
 665   Double_t fgCPAmin=0.99;     //! min cosine of V0's pointing angle
 666   //  Double_t fgRmin=0.2;    //! min radius of the fiducial volume //original
 667   Double_t fgRmin=0.1;        //! min radius of the fiducial volume = 1 mm
 668   Double_t fgRmax=200.;       //! max radius of the fiducial volume = 2 m
 669
 670   //------------------------------------------
 671
 672   // Main loop
 673   // Called for EACH event
 674
 675   AliVEvent *event = InputEvent();
 676   if (!event) { Printf("ERROR: Could not retrieve event"); return; }
 677
 678   Info("AliAnalysisTaskHelium3Pi","Starting UserExec");
 679
 680   SetDataType("REAL");
 681
 682   // create pointer to event
 683   AliESDEvent* lESDevent = dynamic_cast<AliESDEvent*>(event);
 684   if (!lESDevent) {
 685     AliError("Cannot get the ESD event");
 686     return;
 687   }
 688
 689   fHistEventMultiplicity->Fill(0);
 690
 691   Double_t lMagneticField=lESDevent->GetMagneticField();
 692   Int_t TrackNumber = -1;
 693
 694
 695   //*****************//
 696   //*   Centrality  *//
 697   //*****************//
 698
 699   AliCentrality *centrality = lESDevent->GetCentrality();
 700   Float_t percentile=centrality->GetCentralityPercentile("V0M");
 701
 702   TrackNumber = lESDevent->GetNumberOfTracks();
 703   if (TrackNumber<2) return;
 704
 705   //****************************************
 706
 707   // PID
 708
 709   AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
 710   AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
 711   fPIDResponse=inputHandler->GetPIDResponse(); // data member di tipo "const AliPIDResponse *fPIDResponse;"
 712   //  cout<<"fPIDResponse "<<fPIDResponse<<endl;
 713   //===========================================
 714
 715   Int_t eventtype=-99;
 716
 717   Bool_t isSelectedCentral     = (inputHandler->IsEventSelected() & AliVEvent::kCentral);
 718   Bool_t isSelectedSemiCentral = (inputHandler->IsEventSelected() & AliVEvent::kSemiCentral);
 719   Bool_t isSelectedMB          = (inputHandler->IsEventSelected() & AliVEvent::kMB);
 720   Bool_t isSelectedAny         = (inputHandler->IsEventSelected() & AliVEvent::kAny);
 721
 722   if(isSelectedCentral){
 723     fHistEventMultiplicity->Fill(3);
 724     fHistTrackMultiplicityCent->Fill(TrackNumber,percentile);
 725     eventtype=1;
 726   }
 727
 728   if(isSelectedSemiCentral){
 729     fHistEventMultiplicity->Fill(4);
 730     fHistTrackMultiplicitySemiCent->Fill(TrackNumber,percentile);
 731     eventtype=2;
 732   }
 733
 734   if(isSelectedMB){
 735     fHistEventMultiplicity->Fill(5);
 736     fHistTrackMultiplicityMB->Fill(TrackNumber,percentile);
 737     eventtype=3;
 738   }
 739
 740   if(!isSelectedCentral && !isSelectedSemiCentral && !isSelectedMB && isSelectedAny){
 741     fHistEventMultiplicity->Fill(9);
 742     fHistTrackMultiplicity->Fill(TrackNumber,percentile); //tracce per evento
 743     eventtype=4;
 744   }
 745
 746   //if(isSelectedCentral || isSelectedSemiCentral || isSelectedMB || isSelectedAny){
 747   if(eventtype ==1  || eventtype ==2  || eventtype==3  || eventtype==4){
 748
 749     // ANALISYS
 750
 751     // Primary vertex cut
 752
 753     const AliESDVertex *vtx = lESDevent->GetPrimaryVertexTracks();
 754
 755     if(vtx->GetNContributors()<1) {
 756
 757       // SPD vertex cut
 758       vtx = lESDevent->GetPrimaryVertexSPD();
 759
 760       if(vtx->GetNContributors()<1) {
 761         Info("AliAnalysisTaskHelium3Pi","No good vertex, skip event");
 762         return; // NO GOOD VERTEX, SKIP EVENT
 763       }
 764     }
 765
 766     fHistEventMultiplicity->Fill(1); // analyzed events with PV
 767
 768     xPrimaryVertex=vtx->GetX();
 769     yPrimaryVertex=vtx->GetY();
 770     zPrimaryVertex=vtx->GetZ();
 771
 772     if(TMath::Abs(zPrimaryVertex)>10) return;
 773
 774     if(eventtype==1){
 775       fHistTrackMultiplicityPVCent->Fill(TrackNumber,percentile);
 776       fHistEventMultiplicity->Fill(6);
 777     }
 778
 779     if(eventtype==2){
 780       fHistTrackMultiplicityPVSemiCent->Fill(TrackNumber,percentile);
 781       fHistEventMultiplicity->Fill(7);
 782     }
 783
 784     if(eventtype==3){
 785       fHistTrackMultiplicityPVMB->Fill(TrackNumber,percentile);
 786       fHistEventMultiplicity->Fill(8);
 787     }
 788
 789     if(eventtype==4){
 790       fHistEventMultiplicity->Fill(10);
 791     }
 792
 793
 794     fHistEventMultiplicity->Fill(2);
 795
 796     //Find Pair candidates
 797
 798     TArrayI PionsTPC(TrackNumber);        //Neg pions
 799     Int_t nPionsTPC=0;
 800
 801     TArrayI HeTPC(TrackNumber);        //helium3
 802     Int_t nHeTPC=0;
 803
 804     const Double_t speedOfLight =  TMath::C()*1E2*1E-12; // cm/ps
 805
 806     Float_t impactXY=-999, impactZ=-999;
 807     Float_t impactXYpi=-999, impactZpi=-999;
 808
 809
 810     //*************************************************************
 811
 812     runNumber = lESDevent->GetRunNumber();
 813     BCNumber    = lESDevent->GetBunchCrossNumber();
 814     OrbitNumber = lESDevent->GetOrbitNumber();
 815     PeriodNumber= lESDevent->GetPeriodNumber();
 816
 817     //*************************************************************
 818
 819     for (Int_t j=0; j<TrackNumber; j++) { //loop on tracks
 820
 821       AliESDtrack *esdtrack=lESDevent->GetTrack(j);
 822       //      AliVTrack*  esdtrack= (AliVTrack *) fEvent->GetTrack(iT);
 823
 824
 825       if(!esdtrack) {
 826         AliError(Form("ERROR: Could not retrieve esdtrack %d",j));
 827         continue;
 828       }
 829
 830       // ************** Track cuts ****************
 831
 832       if (!fESDtrackCuts->AcceptTrack(esdtrack)) continue;
 833
 834
 835       status  = (ULong_t)esdtrack->GetStatus();
 836       isTPC   = (((status) & AliESDtrack::kTPCin)  != 0);
 837       isTOF   = ((((status) & AliESDtrack::kTOFout) != 0) && (((status) & AliESDtrack::kTIME) != 0));
 838
 839       UInt_t mapITS=esdtrack->GetITSClusterMap();
 840
 841       //----------------------------------------------
 842
 843       //****** Cuts from  AliV0Vertex.cxx *************
 844
 845       Double_t d=esdtrack->GetD(xPrimaryVertex,yPrimaryVertex,lMagneticField);
 846       //    if (TMath::Abs(d)<fgDPmin) continue;
 847       if (TMath::Abs(d)>fgRmax) continue;
 848
 849       //---- (Usefull) Stuff
 850
 851       TPCSignal=esdtrack->GetTPCsignal();
 852
 853       if (TPCSignal<10)continue;
 854       if (TPCSignal>1000)continue;
 855
 856       if(!isTPC)continue;
 857       if(!esdtrack->GetTPCInnerParam())continue;
 858
 859       AliExternalTrackParam trackIn(*esdtrack->GetInnerParam());
 860       pinTPC= trackIn.GetP();
 861
 862       //pinTPC= esdtrack->GetTPCMomentum();
 863
 864       poutTPC=pinTPC;
 865
 866
 867       if((status) & (AliESDtrack::kITSrefit!=0)){
 868         fhBB->Fill(pinTPC*esdtrack->GetSign(),TPCSignal);
 869       }
 870
 871       timeTOF=esdtrack->GetTOFsignal();                 // ps
 872       trackLenghtTOF= esdtrack->GetIntegratedLength();  // cm
 873
 874       if(isTOF){
 875
 876         if(!esdtrack->GetOuterParam())continue;
 877
 878         AliExternalTrackParam trackOut(*esdtrack->GetOuterParam());
 879
 880         poutTPC = trackOut.GetP();
 881
 882         betaTOF= (trackLenghtTOF/timeTOF)/2.99792458e-2;
 883
 884         fhTOF->Fill(poutTPC*esdtrack->GetSign(),betaTOF);
 885
 886         Double_t mass2=(poutTPC*poutTPC)*((((speedOfLight*speedOfLight)*(timeTOF*timeTOF))-(trackLenghtTOF*trackLenghtTOF))/(trackLenghtTOF*trackLenghtTOF));
 887         if(mass2>0) massTOF=TMath::Sqrt(mass2);
 888         fhMassTOF->Fill(massTOF);
 889
 890         if(esdtrack->GetSign() < 0.)fBetavsTPCsignalNeg->Fill(betaTOF,TPCSignal);
 891         if(esdtrack->GetSign() > 0.)fBetavsTPCsignalPos->Fill(betaTOF,TPCSignal);
 892
 893       }
 894
 895       //pass2
 896
 897       // bbtheo =BetheBloch((2*pinTPC)/3.,2,kTRUE);    //! OK
 898       // bbtheoM=(1 - 0.08*5)*bbtheo;                  //! OK
 899       // bbtheoP=(1 + 0.08*5)*bbtheo;                  //! OK
 900
 901
 902       bbtheo = fPIDResponse->NumberOfSigmas((AliPIDResponse::EDetector)0,esdtrack,(AliPID::EParticleType) 7);
 903
 904       if(esdtrack->GetSign()<0){
 905         zNathashaNeg=bbtheo;//(TPCSignal-bbtheo)/bbtheo;
 906         //      cout<<"BBtheo 1 :"<<zNathashaNeg<<endl;
 907         fhNaNeg->Fill(pinTPC,zNathashaNeg);
 908       }
 909
 910       if(esdtrack->GetSign() > 0.){
 911         zNathashaPos=bbtheo;//(TPCSignal-bbtheo)/bbtheo;
 912         fhNaPos->Fill(pinTPC,zNathashaPos);
 913       }
 914
 915       nSigmaNegPion=TMath::Abs(fPIDResponse->NumberOfSigmasTPC(esdtrack,(AliPID::EParticleType) 2));
 916       //2 is pion
 917
 918       if ( (nSigmaNegPion < cutNSigma)){
 919
 920         //      cout<<"Nsigma pi: "<<nSigmaNegPion<<endl;
 921
 922         fhBBPions->Fill(pinTPC*esdtrack->GetSign(),TPCSignal);
 923
 924         if(pinTPC<3.){
 925           PionsTPC[nPionsTPC++]=j;
 926         }
 927       }
 928
 929       //      nSigmaNegPion=(fPIDResponse->NumberOfSigmasTPC(esdtrack,(AliPID::EParticleType) 2));
 930
 931       bbtheoM = TMath::Abs((fPIDResponse->NumberOfSigmasTPC(esdtrack,(AliPID::EParticleType) 7)));
 932
 933       //      if( TPCSignal > bbtheoM ) {
 934       //      if( bbtheoM > -3.) {
 935       if( bbtheoM < 3.) {
 936
 937         if(pinTPC>0.6){
 938
 939           fhBBHe->Fill(pinTPC*esdtrack->GetSign(),TPCSignal);
 940           HeTPC[nHeTPC++]=j;
 941
 942           Bool_t isHeITSrefit=((status) & (AliESDtrack::kITSrefit));
 943
 944           esdtrack->GetImpactParameters(impactXY, impactZ);
 945
 946           Int_t  fIdxInt[200]; //dummy array
 947           Int_t nClustersTPC = esdtrack->GetTPCclusters(fIdxInt);
 948
 949           Float_t chi2PerClusterTPC = esdtrack->GetTPCchi2()/(Float_t)(nClustersTPC);
 950
 951           tHelrunNumber         =(Float_t)runNumber;
 952           tHelBCNumber          =(Float_t)BCNumber;
 953           tHelOrbitNumber       =(Float_t)OrbitNumber;
 954           tHelPeriodNumber      =(Float_t)PeriodNumber;
 955           tHeleventtype         =(Float_t)eventtype;
 956           tHelisHeITSrefit      =(Float_t)isHeITSrefit;
 957           tHelpercentile        =(Float_t)percentile;
 958           tHelSign              =(Float_t)esdtrack->GetSign();
 959           tHelpinTPC            =(Float_t)pinTPC;
 960           tHelGetTPCsignal      =(Float_t)esdtrack->GetTPCsignal();
 961           tHelPx                =(Float_t)esdtrack->Px();
 962           tHelPy                =(Float_t)esdtrack->Py();
 963           tHelPz                =(Float_t)esdtrack->Pz();
 964           tHelEta               =(Float_t)esdtrack->Eta();
 965           tHelisTOF             =(Float_t)isTOF;
 966           tHelpoutTPC           =(Float_t)poutTPC;
 967           tHeltimeTOF           =(Float_t)timeTOF;
 968           tHeltrackLenghtTOF    =(Float_t)trackLenghtTOF;
 969           tHelimpactXY          =(Float_t)impactXY;
 970           tHelimpactZ           =(Float_t)impactZ;
 971           tHelmapITS            =(Float_t)mapITS;
 972           tHelTPCNcls           =(Float_t)esdtrack->GetTPCNcls();
 973           tHelTRDsignal         =(Float_t)esdtrack->GetTRDsignal();
 974           tHelxPrimaryVertex    =(Float_t)xPrimaryVertex;
 975           tHelyPrimaryVertex    =(Float_t)yPrimaryVertex;
 976           tHelzPrimaryVertex    =(Float_t)zPrimaryVertex;
 977           tHelchi2PerClusterTPC =(Float_t)chi2PerClusterTPC;
 978
 979           //      fNtuple4->Fill();
 980         }
 981       }
 982     }  //! track
 983
 984     PionsTPC.Set(nPionsTPC);
 985     HeTPC.Set(nHeTPC);
 986
 987     Double_t        DcaHeToPrimVertex=0;
 988     Double_t        DcaPionToPrimVertex=0;
 989
 990     impactXY=-999, impactZ=-999;
 991     impactXYpi=-999, impactZpi=-999;
 992
 993     // Track
 994
 995     // Vettors for il PxPyPz
 996
 997     Double_t momPionVett[3];
 998     for(Int_t i=0;i<3;i++)momPionVett[i]=0;
 999
1000     Double_t momHeVett[3];
1001     for(Int_t i=0;i<3;i++)momHeVett[i]=0;
1002
1003     //At SV
1004
1005     Double_t momPionVettAt[3];
1006     for(Int_t i=0;i<3;i++)momPionVettAt[i]=0;
1007
1008     Double_t momHeVettAt[3];
1009     for(Int_t i=0;i<3;i++)momHeVettAt[i]=0;
1010
1011     //---------------   LOOP PAIRS   ----------------
1012
1013     for (Int_t k=0; k < nPionsTPC; k++) {                           //! Pions Loop
1014
1015       DcaPionToPrimVertex=0.;
1016       DcaHeToPrimVertex=0;
1017
1018       Int_t PionIdx=PionsTPC[k];
1019
1020       AliESDtrack  *PionTrack=lESDevent->GetTrack(PionIdx);
1021
1022       statusPi = (ULong_t)PionTrack->GetStatus();
1023       //      isTOFPi  = ((((statusPi) & (AliESDtrack::kTOFout)) != 0) && (((statusPi) & (AliESDtrack::kTIME)) != 0));
1024       IsPiITSRefit = ((statusPi) & (AliESDtrack::kITSrefit));
1025
1026       if (PionTrack)
1027         DcaPionToPrimVertex = TMath::Abs(PionTrack->GetD(xPrimaryVertex, yPrimaryVertex,lMagneticField)); //OK
1028
1029       if(DcaPionToPrimVertex<0.2)continue;
1030
1031       AliExternalTrackParam trackInPion(*PionTrack);
1032
1033       for (Int_t i=0; i<nHeTPC; i++){                               //! Helium Loop
1034
1035         Int_t HeIdx=HeTPC[i];
1036
1037          AliESDtrack  *HeTrack=lESDevent->GetTrack(HeIdx);
1038
1039         //      statusT= (ULong_t)HeTrack->GetStatus();
1040         //      isTOFHe   = (((statusT & AliESDtrack::kTOFout) != 0) && ((statusT & AliESDtrack::kTIME) != 0));
1041         IsHeITSRefit = (status & AliESDtrack::kITSrefit);
1042
1043         if (HeTrack)
1044           DcaHeToPrimVertex = TMath::Abs(HeTrack->GetD(xPrimaryVertex, yPrimaryVertex,lMagneticField)); //OK
1045
1046         AliExternalTrackParam trackInHe(*HeTrack);
1047
1048         if ( DcaPionToPrimVertex < fgDNmin)                //OK
1049           if ( DcaHeToPrimVertex < fgDNmin) continue;    //OK
1050
1051         Double_t xn, xp;
1052         Double_t dca=0.;
1053
1054         dca= PionTrack->GetDCA(HeTrack,lMagneticField,xn,xp); //!dca (Neg to Pos)
1055
1056         if (dca > fgDCAmax) continue;
1057         if ((xn+xp) > 2*fgRmax) continue;
1058         if ((xn+xp) < 2*fgRmin) continue;
1059
1060         //CORRECTION from AliV0Vertex
1061
1062         Bool_t corrected=kFALSE;
1063         if ((trackInPion.GetX() > 3.) && (xn < 3.)) {
1064           //correct for the beam pipe material
1065           corrected=kTRUE;
1066         }
1067         if ((trackInHe.GetX() > 3.) && (xp < 3.)) {
1068           //correct for the beam pipe material
1069           corrected=kTRUE;
1070         }
1071         if (corrected) {
1072           dca=trackInPion.GetDCA(&trackInHe,lMagneticField,xn,xp);
1073           if (dca > fgDCAmax) continue;
1074           if ((xn+xp) > 2*fgRmax) continue;
1075           if ((xn+xp) < 2*fgRmin) continue;
1076         }
1077
1078         //=============================================//
1079         // Make "V0" with found tracks                 //
1080         //=============================================//
1081
1082         trackInPion.PropagateTo(xn,lMagneticField);
1083         trackInHe.PropagateTo(xp,lMagneticField);
1084
1085         AliESDv0 vertex(trackInPion,PionIdx,trackInHe,HeIdx);
1086         if (vertex.GetChi2V0() > fgChi2max) continue;
1087
1088         Float_t CosPointingAngle=vertex.GetV0CosineOfPointingAngle(xPrimaryVertex,yPrimaryVertex,zPrimaryVertex); //PointingAngle
1089         if (CosPointingAngle < fgCPAmin) continue;
1090
1091         vertex.SetDcaV0Daughters(dca);
1092         vertex.SetV0CosineOfPointingAngle(CosPointingAngle);
1093
1094         fPos[0]=vertex.Xv();
1095         fPos[1]=vertex.Yv();
1096         fPos[2]=vertex.Zv();
1097
1098         HeTrack->PxPyPz(momHeVett);
1099         PionTrack->PxPyPz(momPionVett);
1100
1101         Double_t raggio=TMath::Sqrt(fPos[0]*fPos[0]+fPos[1]*fPos[1]+fPos[2]*fPos[2]);
1102         HeTrack->GetPxPyPzAt(raggio,lMagneticField,momHeVettAt);
1103         PionTrack->GetPxPyPzAt(raggio,lMagneticField,momPionVettAt);
1104
1105         //------------------------------------------------------------------------//
1106
1107         HeTrack->GetImpactParameters(impactXY, impactZ);
1108
1109         PionTrack->GetImpactParameters(impactXYpi, impactZpi);
1110
1111         if(vertex.GetD(xPrimaryVertex,yPrimaryVertex,zPrimaryVertex)>3) continue;
1112
1113         //salvo solo fino a 3.1 GeV/c2
1114
1115         vHelium.SetXYZM(2*momHeVettAt[0],2*momHeVettAt[1],2*momHeVettAt[2],Helium3Mass);
1116         vPion.SetXYZM(momPionVettAt[0],momPionVettAt[1],momPionVettAt[2],PionMass);
1117         vSum=vHelium+vPion;
1118
1119         if(vSum.M()>3.2)
1120           continue;
1121
1122         Int_t  fIdxInt[200]; //dummy array
1123         Int_t nClustersTPCHe = HeTrack->GetTPCclusters(fIdxInt);
1124         Int_t nClustersTPCPi = PionTrack->GetTPCclusters(fIdxInt);
1125
1126         //----------------------------------------------------------------------//
1127
1128         trunNumber              =(Float_t)runNumber;
1129         tbunchcross             =(Float_t)BCNumber;
1130         torbit                  =(Float_t)OrbitNumber;
1131         tperiod                 =(Float_t)PeriodNumber;
1132         teventtype              =(Float_t)eventtype;
1133         tTrackNumber            =(Float_t)TrackNumber;
1134         tpercentile             =(Float_t)percentile;
1135         txPrimaryVertex         =(Float_t)xPrimaryVertex; //PRIMARY
1136         tyPrimaryVertex         =(Float_t)yPrimaryVertex;
1137         tzPrimaryVertex         =(Float_t)zPrimaryVertex;
1138         txSecondaryVertex       =(Float_t)fPos[0]; //SECONDARY
1139         tySecondaryVertex       =(Float_t)fPos[1];
1140         tzSecondaryVertex       =(Float_t)fPos[2];
1141         tdcaTracks              =(Float_t)dca;           //between 2 tracks
1142         tCosPointingAngle       =(Float_t)CosPointingAngle;          //cosPointingAngle da V0
1143         tDCAV0toPrimaryVertex   =(Float_t)vertex.GetD(xPrimaryVertex,yPrimaryVertex,zPrimaryVertex);
1144         tHeSign                 =(Float_t)HeTrack->GetSign(); //He
1145         tHepInTPC               =(Float_t)trackInHe.GetP();
1146         tHeTPCsignal            =(Float_t)HeTrack->GetTPCsignal();
1147         tDcaHeToPrimVertex      =(Float_t)DcaHeToPrimVertex;
1148         tHeEta                  =(Float_t)HeTrack->Eta();
1149         tmomHex                 =(Float_t)momHeVett[0];
1150         tmomHey                 =(Float_t)momHeVett[1];
1151         tmomHez                 =(Float_t)momHeVett[2];
1152         tmomHeAtSVx             =(Float_t)momHeVettAt[0];
1153         tmomHeAtSVy             =(Float_t)momHeVettAt[1];
1154         tmomHeAtSVz             =(Float_t)momHeVettAt[2];
1155         tHeTPCNcls              =(Float_t)HeTrack->GetTPCNcls();
1156         tHeimpactXY             =(Float_t)impactXY;
1157         tHeimpactZ              =(Float_t)impactZ;
1158         tHeITSClusterMap        =(Float_t)HeTrack->GetITSClusterMap();
1159         tIsHeITSRefit           =(Float_t)IsHeITSRefit;
1160         tPionSign               =(Float_t)PionTrack->GetSign(); //Pion
1161         tPionpInTPC             =(Float_t)trackInPion.GetP();
1162         tPionTPCsignal          =(Float_t)PionTrack->GetTPCsignal();
1163         tDcaPionToPrimVertex    =(Float_t)DcaPionToPrimVertex;
1164         tPionEta                =(Float_t)PionTrack->Eta();
1165         tmomPionx               =(Float_t)momPionVett[0];
1166         tmomPiony               =(Float_t)momPionVett[1];
1167         tmomPionz               =(Float_t)momPionVett[2];
1168         tmomNegPionAtSVx        =(Float_t)momPionVettAt[0];
1169         tmomNegPionAtSVy        =(Float_t)momPionVettAt[1];
1170         tmomNegPionAtSVz        =(Float_t)momPionVettAt[2];
1171         tPionTPCNcls            =(Float_t)PionTrack->GetTPCNcls();
1172         tPionimpactXY           =(Float_t)impactXYpi;
1173         tPionimpactZ            =(Float_t)impactZpi;
1174         tPionITSClusterMap      =(Float_t)PionTrack->GetITSClusterMap();
1175         tIsPiITSRefit           =(Float_t)IsPiITSRefit;
1176         txn                     =(Float_t)xn;
1177         txp                     =(Float_t)xp;
1178         tchi2He                 =(Float_t)HeTrack->GetTPCchi2()/(Float_t)(nClustersTPCHe);
1179         tchi2Pi                 =(Float_t)PionTrack->GetTPCchi2()/(Float_t)(nClustersTPCPi);
1180
1181
1182         fNtuple1->Fill();
1183         vertex.Delete();
1184       }// positive TPC
1185
1186     } //negative tpc
1187
1188   }
1189
1190   PostData(1,fListHist);
1191   PostData(2,fNtuple1);
1192   PostData(3,fNtuple4);
1193
1194 } //end userexec
1195
1196
1197 //________________________________________________________________________
1198
1199 void AliAnalysisTaskHelium3Pi::Terminate(Option_t *)
1200 {
1201   // Draw result to the screen
1202   // Called once at the end of the query
1203 }
1204
1205