]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Adding Tasks to Read AOD for (3He,pi) analysis
authorlramona <ramona.lea@cern.ch>
Wed, 19 Mar 2014 14:30:03 +0000 (15:30 +0100)
committerhristov <Peter.Hristov@cern.ch>
Thu, 27 Mar 2014 15:25:15 +0000 (16:25 +0100)
PWGLF/CMakelibPWGLFSTRANGENESS.pkg
PWGLF/PWGLFSTRANGENESSLinkDef.h
PWGLF/STRANGENESS/Hypernuclei/AddTask_Helium3PiAOD.C [new file with mode: 0644]
PWGLF/STRANGENESS/Hypernuclei/AliAnalysisTaskHelium3PiAOD.cxx [new file with mode: 0644]
PWGLF/STRANGENESS/Hypernuclei/AliAnalysisTaskHelium3PiAOD.h [new file with mode: 0644]

index 4aa817cf2c7a1d9b7182a06ea1f5a8f52497a549..e9437e47b1ac940bc3eedbb48dc422ad25087436 100644 (file)
@@ -53,6 +53,7 @@ set ( SRCS
   STRANGENESS/Hypernuclei/AliAnalysisTaskLambdaNAOD.cxx
   STRANGENESS/Hypernuclei/AliAnalysisTaskHdibaryonLPpi.cxx
   STRANGENESS/Hypernuclei/AliAnalysisTaskHelium3Pi.cxx
+  STRANGENESS/Hypernuclei/AliAnalysisTaskHelium3PiAOD.cxx
   STRANGENESS/Hypernuclei/AliAnalysisTaskHelium3PiMC.cxx
   STRANGENESS/Hypernuclei/AliAODRecoDecayLF.cxx
   STRANGENESS/Hypernuclei/AliAODRecoDecayLF2Prong.cxx
index 8af5fb4aef341bd9e1175bc97f89701ea1f3b717..8b07f33dc9e427abc3a3d5bd53276541048c46a7 100644 (file)
@@ -32,6 +32,7 @@
 #pragma link C++ class AliAnalysisTaskLambdaNAOD+;
 #pragma link C++ class AliAnalysisTaskHdibaryonLPpi+;
 #pragma link C++ class AliAnalysisTaskHelium3Pi+;
+#pragma link C++ class AliAnalysisTaskHelium3PiAOD+;
 #pragma link C++ class AliAnalysisTaskHelium3PiMC+;
 #pragma link C++ class AliAODRecoDecayLF+;
 #pragma link C++ class AliAODRecoDecayLF2Prong+;
diff --git a/PWGLF/STRANGENESS/Hypernuclei/AddTask_Helium3PiAOD.C b/PWGLF/STRANGENESS/Hypernuclei/AddTask_Helium3PiAOD.C
new file mode 100644 (file)
index 0000000..4c85e8a
--- /dev/null
@@ -0,0 +1,38 @@
+AliAnalysisTask *AddTask_Helium3Pi(TString name="name"){
+
+  //get the current analysis manager
+  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+  if (!mgr) {
+    Error("AddTask_Helium3Pi", "No analysis manager found.");
+    return 0;
+  }
+  
+  //========= Add task to the ANALYSIS manager =====
+
+  AliAnalysisTaskHelium3PiAOD *taskHelium3Pi = new AliAnalysisTaskHelium3PiAOD(name);
+   
+  mgr->AddTask(taskHelium3Pi);
+  
+  //================================================
+  //              data containers
+  //================================================
+  //            find input container
+
+  AliAnalysisDataContainer *cinput   = mgr->GetCommonInputContainer();
+  
+  TString outputFileName = AliAnalysisManager::GetCommonFileName();
+  //AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("Helium3Pi_tree", TTree::Class(), AliAnalysisManager::kOutputContainer, "AnalysisResults.root");  
+  //AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("Helium3Pi_tree", TTree::Class(), AliAnalysisManager::kOutputContainer, "AnalysisResults.root");  
+  AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("clisthistHyper", TList::Class(), AliAnalysisManager::kOutputContainer, outputFileName);
+   
+  AliAnalysisDataContainer *coutput2 = mgr->CreateContainer("treeHyper", TTree::Class(),AliAnalysisManager::kOutputContainer,outputFileName);
+  AliAnalysisDataContainer *coutput3 = mgr->CreateContainer("treeHelium" , TTree::Class(),AliAnalysisManager::kOutputContainer,outputFileName);
+  //           connect containers
+  mgr->ConnectInput  (taskHelium3Pi,  0, cinput );
+  mgr->ConnectOutput (taskHelium3Pi,  1, coutput1);
+  mgr->ConnectOutput (taskHelium3Pi,  2, coutput2);
+  mgr->ConnectOutput (taskHelium3Pi,  3, coutput3);
+
+  return taskHelium3Pi;
+}
diff --git a/PWGLF/STRANGENESS/Hypernuclei/AliAnalysisTaskHelium3PiAOD.cxx b/PWGLF/STRANGENESS/Hypernuclei/AliAnalysisTaskHelium3PiAOD.cxx
new file mode 100644 (file)
index 0000000..7c4973c
--- /dev/null
@@ -0,0 +1,1196 @@
+/**************************************************************************
+ * Contributors are not mentioned at all.                                 *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+//
+//-----------------------------------------------------------------
+//                 AliAnalysisTaskHelium3PiAOD class
+//-----------------------------------------------------------------
+
+
+class TTree;
+class TParticle;
+class TVector3;
+
+#include "AliAnalysisManager.h"
+#include <AliMCEventHandler.h>
+#include <AliMCEvent.h>
+#include <AliStack.h>
+
+class AliESDVertex;
+class AliAODVertex;
+class AliESDv0;
+class AliAODv0; 
+class AliCascadeVertexer;
+
+#include <iostream>
+#include "AliAnalysisTaskSE.h"
+#include "TList.h"
+#include "TH1.h"
+#include "TH2.h"
+#include "TH3.h"
+#include "TNtuple.h"
+#include "TGraph.h"
+#include "TCutG.h"
+#include "TF1.h"
+#include "TCanvas.h"
+#include "TMath.h"
+#include "TChain.h"
+#include "Riostream.h"
+#include "AliLog.h"
+#include "AliCascadeVertexer.h"
+#include "AliESDEvent.h"
+#include "AliESDtrack.h"
+#include "AliExternalTrackParam.h"
+#include "AliAODEvent.h"
+#include "AliInputEventHandler.h"
+#include "AliESDcascade.h"
+#include "AliAODcascade.h"
+#include "AliAnalysisTaskHelium3PiAOD.h"
+#include "AliESDtrackCuts.h"
+#include "AliCentrality.h"
+#include "TString.h"
+#include <TDatime.h>
+#include <TRandom3.h>
+#include <TLorentzVector.h>
+#include <AliVTrack.h>
+
+
+ClassImp(AliAnalysisTaskHelium3PiAOD)
+
+//________________________________________________________________________
+AliAnalysisTaskHelium3PiAOD::AliAnalysisTaskHelium3PiAOD() 
+: AliAnalysisTaskSE(),
+  fAnalysisType("AOD"), 
+  fCollidingSystems(0), 
+  fESDtrackCuts(0),
+  fDataType("REAL"),
+  fListHist(0), 
+  fHistEventMultiplicity(0),         
+  fHistTrackMultiplicity(0),      
+  fHistTrackMultiplicityCent(0),      
+  fHistTrackMultiplicitySemiCent(0),  
+  fHistTrackMultiplicityMB(0),        
+  fHistTrackMultiplicityPVCent(0),      
+  fHistTrackMultiplicityPVSemiCent(0),  
+  fHistTrackMultiplicityPVMB(0),        
+  fhBB(0),    
+  fhTOF(0),   
+  fhMassTOF(0),
+  fhBBPions(0),
+  fhBBHe(0),   
+  fhNaPos(0),  
+  fhNaNeg(0),  
+  fBetavsTPCsignalPos(0),  
+  fBetavsTPCsignalNeg(0),  
+  fNtuple1(0),
+  trunNumber(0),
+  tbunchcross(0),
+  torbit(0),
+  tperiod(0),
+  teventtype(0),
+  tTrackNumber(0),
+  tpercentile(0),
+  txPrimaryVertex(0),
+  tyPrimaryVertex(0),
+  tzPrimaryVertex(0),
+  txSecondaryVertex(0),
+  tySecondaryVertex(0),
+  tzSecondaryVertex(0),
+  tdcaTracks(0),
+  tCosPointingAngle(0),
+  tDCAV0toPrimaryVertex(0),
+  tHeSign(0),
+  tHepInTPC(0),
+  tHeTPCsignal(0),
+  tDcaHeToPrimVertex(0),
+  tHeEta(0),
+  tmomHex(0),
+  tmomHey(0),
+  tmomHez(0),
+  tmomHeAtSVx(0),
+  tmomHeAtSVy(0),
+  tmomHeAtSVz(0),
+  tHeTPCNcls(0),
+  tHeimpactXY(0),
+  tHeimpactZ(0),
+  tHeITSClusterMap(0),
+  tIsHeITSRefit(0),
+  tPionSign(0),
+  tPionpInTPC(0),
+  tPionTPCsignal(0),
+  tDcaPionToPrimVertex(0),
+  tPionEta(0),
+  tmomPionx(0),
+  tmomPiony(0),
+  tmomPionz(0),
+  tmomNegPionAtSVx(0),
+  tmomNegPionAtSVy(0),
+  tmomNegPionAtSVz(0),
+  tPionTPCNcls(0),
+  tPionimpactXY(0),
+  tPionimpactZ(0),
+  tPionITSClusterMap(0),
+  tIsPiITSRefit(0),
+  txn(0),
+  txp(0),
+  tchi2He(0),
+  tchi2Pi(0),
+  fNtuple4(0),
+  tHelrunNumber(0),
+  tHelBCNumber(0),
+  tHelOrbitNumber(0),
+  tHelPeriodNumber(0),
+  tHeleventtype(0),
+  tHelisHeITSrefit(0),
+  tHelpercentile(0),
+  tHelSign(0),
+  tHelpinTPC(0),
+  tHelGetTPCsignal(0),
+  tHelPx(0),
+  tHelPy(0),
+  tHelPz(0),
+  tHelEta(0),
+  tHelisTOF(0),
+  tHelpoutTPC(0),
+  tHeltimeTOF(0),
+  tHeltrackLenghtTOF(0),
+  tHelimpactXY(0),
+  tHelimpactZ(0),
+  tHelmapITS(0),
+  tHelTPCNcls(0),
+  tHelTRDsignal(0),
+  tHelxPrimaryVertex(0),
+  tHelyPrimaryVertex(0),
+  tHelzPrimaryVertex(0),
+  tHelchi2PerClusterTPC(0),
+  fPIDResponse(0)
+  
+{
+  // Dummy Constructor
+  
+  //  printf("Dummy Constructor");
+  
+  fESDtrackCuts = new AliESDtrackCuts("fESDtrackCuts");
+  fESDtrackCuts->SetRequireITSStandAlone(kFALSE);
+  fESDtrackCuts->SetRequireITSPureStandAlone(kFALSE);
+      
+  fESDtrackCuts->SetRequireTPCRefit(kTRUE);
+  fESDtrackCuts->SetAcceptKinkDaughters(kFALSE);
+  fESDtrackCuts->SetMinNClustersTPC(60);
+  fESDtrackCuts->SetMaxChi2PerClusterTPC(5);
+  fESDtrackCuts->SetEtaRange(-0.9,0.9);
+
+}
+
+//________________________________________________________________________
+AliAnalysisTaskHelium3PiAOD::AliAnalysisTaskHelium3PiAOD(const char *name) 
+: AliAnalysisTaskSE(name), 
+  fAnalysisType("AOD"), 
+  fCollidingSystems(0), 
+  fESDtrackCuts(0),
+  fDataType("REAL"),
+  fListHist(0), 
+  fHistEventMultiplicity(0),         
+  fHistTrackMultiplicity(0),      
+  fHistTrackMultiplicityCent(0),      
+  fHistTrackMultiplicitySemiCent(0),  
+  fHistTrackMultiplicityMB(0),        
+  fHistTrackMultiplicityPVCent(0),      
+  fHistTrackMultiplicityPVSemiCent(0),  
+  fHistTrackMultiplicityPVMB(0),        
+  fhBB(0),    
+  fhTOF(0),   
+  fhMassTOF(0),
+  fhBBPions(0),
+  fhBBHe(0),   
+  fhNaPos(0),  
+  fhNaNeg(0),  
+  fBetavsTPCsignalPos(0),  
+  fBetavsTPCsignalNeg(0),  
+  fNtuple1(0),
+  trunNumber(0),
+  tbunchcross(0),
+  torbit(0),
+  tperiod(0),
+  teventtype(0),
+  tTrackNumber(0),
+  tpercentile(0),
+  txPrimaryVertex(0),
+  tyPrimaryVertex(0),
+  tzPrimaryVertex(0),
+  txSecondaryVertex(0),
+  tySecondaryVertex(0),
+  tzSecondaryVertex(0),
+  tdcaTracks(0),
+  tCosPointingAngle(0),
+  tDCAV0toPrimaryVertex(0),
+  tHeSign(0),
+  tHepInTPC(0),
+  tHeTPCsignal(0),
+  tDcaHeToPrimVertex(0),
+  tHeEta(0),
+  tmomHex(0),
+  tmomHey(0),
+  tmomHez(0),
+  tmomHeAtSVx(0),
+  tmomHeAtSVy(0),
+  tmomHeAtSVz(0),
+  tHeTPCNcls(0),
+  tHeimpactXY(0),
+  tHeimpactZ(0),
+  tHeITSClusterMap(0),
+  tIsHeITSRefit(0),
+  tPionSign(0),
+  tPionpInTPC(0),
+  tPionTPCsignal(0),
+  tDcaPionToPrimVertex(0),
+  tPionEta(0),
+  tmomPionx(0),
+  tmomPiony(0),
+  tmomPionz(0),
+  tmomNegPionAtSVx(0),
+  tmomNegPionAtSVy(0),
+  tmomNegPionAtSVz(0),
+  tPionTPCNcls(0),
+  tPionimpactXY(0),
+  tPionimpactZ(0),
+  tPionITSClusterMap(0),
+  tIsPiITSRefit(0),
+  txn(0),
+  txp(0),
+  tchi2He(0),
+  tchi2Pi(0),
+  fNtuple4(0),
+  tHelrunNumber(0),
+  tHelBCNumber(0),
+  tHelOrbitNumber(0),
+  tHelPeriodNumber(0),
+  tHeleventtype(0),
+  tHelisHeITSrefit(0),
+  tHelpercentile(0),
+  tHelSign(0),
+  tHelpinTPC(0),
+  tHelGetTPCsignal(0),
+  tHelPx(0),
+  tHelPy(0),
+  tHelPz(0),
+  tHelEta(0),
+  tHelisTOF(0),
+  tHelpoutTPC(0),
+  tHeltimeTOF(0),
+  tHeltrackLenghtTOF(0),
+  tHelimpactXY(0),
+  tHelimpactZ(0),
+  tHelmapITS(0),
+  tHelTPCNcls(0),
+  tHelTRDsignal(0),
+  tHelxPrimaryVertex(0),
+  tHelyPrimaryVertex(0),
+  tHelzPrimaryVertex(0),
+  tHelchi2PerClusterTPC(0),
+  fPIDResponse(0)
+{                                        
+
+  // Define input and output slots here
+  // Input slot #0 works with a TChain
+  //DefineInput(0, TChain::Class());
+  // Output slot #0 writes into a TList container ()
+
+  DefineInput(0, TChain::Class());
+
+  DefineOutput(1, TList::Class());
+  DefineOutput(2, TTree::Class());
+  DefineOutput(3, TTree::Class());
+
+  fESDtrackCuts = new AliESDtrackCuts("fESDtrackCuts");
+  fESDtrackCuts->SetRequireITSStandAlone(kFALSE);
+  fESDtrackCuts->SetRequireITSPureStandAlone(kFALSE);
+      
+  fESDtrackCuts->SetRequireTPCRefit(kTRUE);
+  fESDtrackCuts->SetAcceptKinkDaughters(kFALSE);
+  fESDtrackCuts->SetMinNClustersTPC(60);
+  fESDtrackCuts->SetMaxChi2PerClusterTPC(5);
+  fESDtrackCuts->SetEtaRange(-0.9,0.9);
+
+}
+//_______________________________________________________
+AliAnalysisTaskHelium3PiAOD::~AliAnalysisTaskHelium3PiAOD() 
+{ 
+  // Destructor
+  if (fListHist) {
+    delete fListHist;
+    fListHist = 0;
+  }
+  
+  if (fESDtrackCuts) delete fESDtrackCuts;
+  if(fNtuple1) delete fNtuple1;
+  if(fNtuple4) delete fNtuple4;
+}
+//=================DEFINITION BETHE BLOCH==============================
+
+Double_t AliAnalysisTaskHelium3PiAOD::BetheBloch(Double_t betaGamma,Double_t charge,Bool_t isPbPb) {
+
+  Double_t kp1, kp2, kp3, kp4, kp5;
+  
+  if(isPbPb){
+
+    //    pass1 2011
+    kp1 = 4.7*charge*charge;
+    kp2 = 8.98482806165147636e+00;
+    kp3 = 1.54000000000000005e-05;
+    kp4 = 2.30445734159456084e+00;
+    kp5 = 2.25624744086878559e+00;
+
+  }
+  
+  else{
+
+    // to be defined ...
+    //pass1 2011
+    kp1 = 4.7*charge*charge;
+    kp2 = 8.98482806165147636e+00;
+    kp3 = 1.54000000000000005e-05;
+    kp4 = 2.30445734159456084e+00;
+    kp5 = 2.25624744086878559e+00;
+
+  }
+
+  Double_t beta = betaGamma / TMath::Sqrt(1.0 + betaGamma * betaGamma);
+  
+  Double_t aa = TMath::Power(beta, kp4);
+  Double_t bb = TMath::Power(1.0 / betaGamma, kp5);
+  
+  bb = TMath::Log(kp3 + bb);
+  
+  Double_t out = (kp2 - aa - bb) * kp1 / aa;
+
+  return out;
+}
+
+//==================DEFINITION OF OUTPUT OBJECTS==============================
+
+void AliAnalysisTaskHelium3PiAOD::UserCreateOutputObjects()
+{
+
+  fListHist = new TList();
+  fListHist->SetOwner();  // IMPORTANT!
+
+  if(! fHistEventMultiplicity ){
+    fHistEventMultiplicity   = new TH1F( "fHistEventMultiplicity" , "Nb of Events" , 10 , -0.5, 9.5);
+    fHistEventMultiplicity->GetXaxis()->SetBinLabel(1,"All Events");
+    fHistEventMultiplicity->GetXaxis()->SetBinLabel(2,"Events w/PV");
+    fHistEventMultiplicity->GetXaxis()->SetBinLabel(3,"Events w/|Vz|<10cm");
+    fHistEventMultiplicity->GetXaxis()->SetBinLabel(4,"Central Events");
+    fHistEventMultiplicity->GetXaxis()->SetBinLabel(5,"SemiCentral Events");
+    fHistEventMultiplicity->GetXaxis()->SetBinLabel(6,"MB Events");
+    fHistEventMultiplicity->GetXaxis()->SetBinLabel(7,"Central Events  w/|Vz|<10cm");
+    fHistEventMultiplicity->GetXaxis()->SetBinLabel(8,"SemiCentral Events  w/|Vz|<10cm");
+    fHistEventMultiplicity->GetXaxis()->SetBinLabel(9,"MB Events   w/|Vz|<10cm");
+
+    fListHist->Add(fHistEventMultiplicity);
+  }
+
+  if(! fHistTrackMultiplicity ){
+    fHistTrackMultiplicity   = new TH2F( "fHistTrackMultiplicity" , "Nb of Tracks", 2500,0, 25000,210,-1,104);
+    fHistTrackMultiplicity->GetXaxis()->SetTitle("Number of tracks");
+    fHistTrackMultiplicity->GetYaxis()->SetTitle("Percentile");
+    fListHist->Add(fHistTrackMultiplicity);
+  } 
+
+  if(! fHistTrackMultiplicityCent ){
+    fHistTrackMultiplicityCent   = new TH2F( "fHistTrackMultiplicityCent", "Nb of Tracks Central Events", 2500,0, 25000,210,-1,104 );
+    fHistTrackMultiplicityCent->GetXaxis()->SetTitle("Number of tracks");
+    fHistTrackMultiplicityCent->GetYaxis()->SetTitle("Percentile");
+    fListHist->Add(fHistTrackMultiplicityCent);
+  } 
+
+  if(! fHistTrackMultiplicitySemiCent ){
+    fHistTrackMultiplicitySemiCent   = new TH2F( "fHistTrackMultiplicitySemiCent" , "Nb of Tracks SemiCentral Events", 2500,0, 25000 ,210,-1,104);
+    fHistTrackMultiplicitySemiCent->GetXaxis()->SetTitle("Number of tracks");
+    fHistTrackMultiplicitySemiCent->GetYaxis()->SetTitle("Percentile");
+    fListHist->Add(fHistTrackMultiplicitySemiCent);
+  } 
+  if(! fHistTrackMultiplicityMB ){
+    fHistTrackMultiplicityMB   = new TH2F( "fHistTrackMultiplicityMB" , "Nb of Tracks MBral Events", 2500,0, 25000,210,-1,104 );
+    fHistTrackMultiplicityMB->GetXaxis()->SetTitle("Number of tracks");
+    fHistTrackMultiplicityMB->GetYaxis()->SetTitle("Percentile");
+    fListHist->Add(fHistTrackMultiplicityMB);
+  } 
+
+  if(! fHistTrackMultiplicityPVCent ){
+    fHistTrackMultiplicityPVCent   = new TH2F( "fHistTrackMultiplicityPVCent" , "Nb of Tracks Central Events", 2500,0, 25000,210,-1,104 );
+    fHistTrackMultiplicityPVCent->GetXaxis()->SetTitle("Number of tracks");
+    fHistTrackMultiplicityPVCent->GetYaxis()->SetTitle("Percentile");
+    fListHist->Add(fHistTrackMultiplicityPVCent);
+  } 
+
+  if(! fHistTrackMultiplicityPVSemiCent ){
+    fHistTrackMultiplicityPVSemiCent   = new TH2F( "fHistTrackMultiplicityPVSemiCent" , "Nb of Tracks SemiCentral Events", 2500,0, 25000 ,210,-1,104);
+    fHistTrackMultiplicityPVSemiCent->GetXaxis()->SetTitle("Number of tracks");
+    fHistTrackMultiplicityPVSemiCent->GetYaxis()->SetTitle("Percentile");
+    fListHist->Add(fHistTrackMultiplicityPVSemiCent);
+  } 
+  if(! fHistTrackMultiplicityPVMB ){
+    fHistTrackMultiplicityPVMB   = new TH2F( "fHistTrackMultiplicityPVMB" , "Nb of Tracks MBral Events", 2500,0, 25000,210,-1,104 );
+    fHistTrackMultiplicityPVMB->GetXaxis()->SetTitle("Number of tracks");
+    fHistTrackMultiplicityPVMB->GetYaxis()->SetTitle("Percentile");
+    fListHist->Add(fHistTrackMultiplicityPVMB);
+  } 
+
+  if(! fhBB ){
+    fhBB = new TH2F( "fhBB" , "BetheBlochTPC" , 120,-6,6,150,0,1500);
+    fhBB->GetXaxis()->SetTitle("p/z (GeV/#it{c})");
+    fhBB->GetYaxis()->SetTitle("TPC Signal");
+    fListHist->Add(fhBB);
+  }
+
+  if(! fhTOF ){
+    fhTOF = new TH2F( "fhTOF" , "Scatter Plot TOF" , 120,-6,6,100,0,1.2);
+    fhTOF->GetXaxis()->SetTitle("p/z (GeV/#it{c})");
+    fhTOF->GetYaxis()->SetTitle("#beta");
+    fListHist->Add(fhTOF);
+  }
+
+  if(! fhMassTOF){
+    fhMassTOF=new TH1F ("fhMassTOF","Particle Mass - TOF", 300,0 ,5);
+    fhMassTOF->GetXaxis()->SetTitle("Mass (GeV/#it{c}^{2})");
+    fListHist->Add(fhMassTOF);
+  }
+
+  if(! fhBBPions ){
+    fhBBPions = new TH2F( "fhBBPions" , "Bethe-Bloch TPC Pions" , 120,-6,6,150,0,1500);
+    fhBBPions->GetXaxis()->SetTitle("p/z (GeV/#it{c})");
+    fhBBPions->GetYaxis()->SetTitle("TPC Signal");
+    fListHist->Add(fhBBPions);
+  }
+  
+  if(! fhBBHe ){
+    fhBBHe = new TH2F( "fhBBHe" , "Bethe-Bloch TPC He" , 120,-6,6,150,0,1500);
+    fhBBHe->GetXaxis()->SetTitle("p/z (GeV/#it{c})");
+    fhBBHe->GetYaxis()->SetTitle("TPC Signal");
+    fListHist->Add(fhBBHe);
+  }
+  
+  if(! fhNaPos ){
+    fhNaPos = new TH2F( "fhNaPos" , "Distribution Pos" , 500,0,5,40,-10,10);
+    fhNaPos->GetXaxis()->SetTitle("p/z (GeV/#it{c})");
+    fhNaPos->GetYaxis()->SetTitle("(TPCSignal-bbtheo)/bbtheo (He)");
+    fListHist->Add(fhNaPos);
+  }
+  
+  if(! fhNaNeg ){
+    fhNaNeg = new TH2F( "fhNaNeg" , "Distribution Neg" , 500,0,5,40,-10,10);
+    fhNaNeg->GetXaxis()->SetTitle("p/z (GeV/#it{c})");
+    fhNaNeg->GetYaxis()->SetTitle("(TPCSignal-bbtheo)/bbtheo (He)");
+    fListHist->Add(fhNaNeg);
+  }
+
+  if(! fBetavsTPCsignalPos ){
+    fBetavsTPCsignalPos = new TH2F("fBetavsTPCsignalPos","fBetavsTPCsignalPos",100,0,1.2,150,0,1500);
+    fBetavsTPCsignalPos->GetXaxis()->SetTitle("#beta");
+    fBetavsTPCsignalPos->GetYaxis()->SetTitle("TPC Signal");
+    fListHist->Add(fBetavsTPCsignalPos);
+  }
+  
+  if(! fBetavsTPCsignalNeg ){
+    fBetavsTPCsignalNeg = new TH2F("fBetavsTPCsignalNeg","fBetavsTPCsignalNeg",100,0,1.2,150,0,1500);
+    fBetavsTPCsignalNeg->GetXaxis()->SetTitle("#beta");
+    fBetavsTPCsignalNeg->GetYaxis()->SetTitle("TPC Signal");
+    fListHist->Add(fBetavsTPCsignalNeg);
+  }
+  
+
+  
+  if(! fNtuple1 ) {
+    
+    fNtuple1 = new TTree("fNtuple1","fNtuple1");
+    
+    fNtuple1->Branch("trunNumber"           ,&trunNumber           ,"trunNumber/F");
+    fNtuple1->Branch("tbunchcross"          ,&tbunchcross          ,"tbunchcross/F");
+    fNtuple1->Branch("torbit"               ,&torbit               ,"torbit/F");
+    fNtuple1->Branch("tperiod"              ,&tperiod              ,"tperiod/F");
+    fNtuple1->Branch("teventtype"           ,&teventtype           ,"teventtype/F");
+    fNtuple1->Branch("tTrackNumber"         ,&tTrackNumber         ,"tTrackNumber/F");
+    fNtuple1->Branch("tpercentile"          ,&tpercentile          ,"tpercentile/F") ;
+    fNtuple1->Branch("txPrimaryVertex"      ,&txPrimaryVertex      ,"txPrimaryVertex/F");
+    fNtuple1->Branch("tyPrimaryVertex"      ,&tyPrimaryVertex      ,"tyPrimaryVertex/F");
+    fNtuple1->Branch("tzPrimaryVertex"      ,&tzPrimaryVertex      ,"tzPrimaryVertex/F");
+    fNtuple1->Branch("txSecondaryVertex"    ,&txSecondaryVertex    ,"txSecondaryVertex/F");
+    fNtuple1->Branch("tySecondaryVertex"    ,&tySecondaryVertex    ,"tySecondaryVertex/F");
+    fNtuple1->Branch("tzSecondaryVertex"    ,&tzSecondaryVertex    ,"tzSecondaryVertex/F");
+    fNtuple1->Branch("tdcaTracks"           ,&tdcaTracks           ,"tdcaTracks/F");
+    fNtuple1->Branch("tCosPointingAngle"    ,&tCosPointingAngle    ,"tCosPointingAngle/F");
+    fNtuple1->Branch("tDCAV0toPrimaryVertex",&tDCAV0toPrimaryVertex,"tDCAV0toPrimaryVertex/F");
+    fNtuple1->Branch("tHeSign"              ,&tHeSign              ,"tHeSign/F");
+    fNtuple1->Branch("tHepInTPC"            ,&tHepInTPC            ,"tHepInTPC/F");
+    fNtuple1->Branch("tHeTPCsignal"         ,&tHeTPCsignal         ,"tHeTPCsignal/F");
+    fNtuple1->Branch("tDcaHeToPrimVertex"   ,&tDcaHeToPrimVertex   ,"tDcaHeToPrimVertex/F");
+    fNtuple1->Branch("tHeEta"               ,&tHeEta               ,"tHeEta/F");
+    fNtuple1->Branch("tmomHex"              ,&tmomHex              ,"tmomHex/F");
+    fNtuple1->Branch("tmomHey"              ,&tmomHey              ,"tmomHey/F");
+    fNtuple1->Branch("tmomHez"              ,&tmomHez              ,"tmomHez/F");
+    fNtuple1->Branch("tmomHeAtSVx"          ,&tmomHeAtSVx          ,"tmomHeAtSVx/F");
+    fNtuple1->Branch("tmomHeAtSVy"          ,&tmomHeAtSVy          ,"tmomHeAtSVy/F");
+    fNtuple1->Branch("tmomHeAtSVz"          ,&tmomHeAtSVz          ,"tmomHeAtSVz/F");
+    fNtuple1->Branch("tHeTPCNcls"           ,&tHeTPCNcls           ,"tHeTPCNcls/F");
+    fNtuple1->Branch("tHeimpactXY"          ,&tHeimpactXY          ,"tHeimpactXY/F");
+    fNtuple1->Branch("tHeimpactZ"           ,&tHeimpactZ           ,"tHeimpactZ/F");
+    fNtuple1->Branch("tHeITSClusterMap"     ,&tHeITSClusterMap     ,"tHeITSClusterMap/F");
+    fNtuple1->Branch("tIsHeITSRefit"        ,&tIsHeITSRefit        ,"tIsHeITSRefit/F");
+    fNtuple1->Branch("tPionSign"            ,&tPionSign            ,"tPionSign/F");
+    fNtuple1->Branch("tPionpInTPC"          ,&tPionpInTPC          ,"tPionpInTPC/F");
+    fNtuple1->Branch("tPionTPCsignal"       ,&tPionTPCsignal       ,"tPionTPCsignal/F");
+    fNtuple1->Branch("tDcaPionToPrimVertex" ,&tDcaPionToPrimVertex ,"tDcaPionToPrimVertex/F");
+    fNtuple1->Branch("tPionEta"             ,&tPionEta             ,"tPionEta/F");
+    fNtuple1->Branch("tmomPionx"            ,&tmomPionx            ,"tmomPionx/F");
+    fNtuple1->Branch("tmomPiony"            ,&tmomPiony            ,"tmomPiony/F");
+    fNtuple1->Branch("tmomPionz"            ,&tmomPionz            ,"tmomPionz/F");
+    fNtuple1->Branch("tmomNegPionAtSVx"     ,&tmomNegPionAtSVx     ,"tmomNegPionAtSVx/F");
+    fNtuple1->Branch("tmomNegPionAtSVy"     ,&tmomNegPionAtSVy     ,"tmomNegPionAtSVy/F");
+    fNtuple1->Branch("tmomNegPionAtSVz"     ,&tmomNegPionAtSVz     ,"tmomNegPionAtSVz/F");
+    fNtuple1->Branch("tPionTPCNcls"         ,&tPionTPCNcls         ,"tPionTPCNcls/F");
+    fNtuple1->Branch("tPionimpactXY"        ,&tPionimpactXY        ,"tPionimpactXY/F");
+    fNtuple1->Branch("tPionimpactZ"         ,&tPionimpactZ         ,"tPionimpactZ/F");
+    fNtuple1->Branch("tPionITSClusterMap"   ,&tPionITSClusterMap   ,"tPionITSClusterMap/F");
+    fNtuple1->Branch("tIsPiITSRefit"        ,&tIsPiITSRefit        ,"tIsPiITSRefit/F");
+    fNtuple1->Branch("txn"                  ,&txn                  ,"txn/F");
+    fNtuple1->Branch("txp"                  ,&txp                  ,"txp/F");
+    fNtuple1->Branch("tchi2He"              ,&tchi2He              ,"tchi2He/F");
+    fNtuple1->Branch("tchi2Pi"              ,&tchi2Pi              ,"tchi2Pi/F");
+    
+  }
+  
+  if(! fNtuple4 ) {
+    fNtuple4 = new TTree("fNtuple4","fNtuple4");
+    
+    fNtuple4->Branch("tHelrunNumber"        ,&tHelrunNumber        ,"tHelrunNumber/F");
+    fNtuple4->Branch("tHelBCNumber"         ,&tHelBCNumber         ,"tHelBCNumber/F");
+    fNtuple4->Branch("tHelOrbitNumber"      ,&tHelOrbitNumber      ,"tHelOrbitNumber/F");
+    fNtuple4->Branch("tHelPeriodNumber"     ,&tHelPeriodNumber     ,"tHelPeriodNumber/F");
+    fNtuple4->Branch("tHeleventtype"        ,&tHeleventtype        ,"tHeleventtype/F");
+    fNtuple4->Branch("tHelisHeITSrefit"     ,&tHelisHeITSrefit     ,"tHelisHeITSrefit/F");
+    fNtuple4->Branch("tHelpercentile"       ,&tHelpercentile       ,"tHelpercentile/F");
+    fNtuple4->Branch("tHelSign"             ,&tHelSign             ,"tHelSign/F");
+    fNtuple4->Branch("tHelpinTPC"           ,&tHelpinTPC           ,"tHelpinTPC/F");
+    fNtuple4->Branch("tHelGetTPCsignal"     ,&tHelGetTPCsignal     ,"tHelGetTPCsignal/F");
+    fNtuple4->Branch("tHelPx"               ,&tHelPx               ,"tHelPx/F");
+    fNtuple4->Branch("tHelPy"               ,&tHelPy               ,"tHelPy/F");
+    fNtuple4->Branch("tHelPz"               ,&tHelPz               ,"tHelPz/F");
+    fNtuple4->Branch("tHelEta"              ,&tHelEta              ,"tHelEta/F");
+    fNtuple4->Branch("tHelisTOF"            ,&tHelisTOF            ,"tHelisTOF/F");
+    fNtuple4->Branch("tHelpoutTPC"          ,&tHelpoutTPC          ,"tHelpoutTPC/F");
+    fNtuple4->Branch("tHeltimeTOF"          ,&tHeltimeTOF          ,"tHeltimeTOF/F");
+    fNtuple4->Branch("tHeltrackLenghtTOF"   ,&tHeltrackLenghtTOF   ,"tHeltrackLenghtTOF/F");
+    fNtuple4->Branch("tHelimpactXY"         ,&tHelimpactXY         ,"tHelimpactXY/F");
+    fNtuple4->Branch("tHelimpactZ"          ,&tHelimpactZ          ,"tHelimpactZ/F");
+    fNtuple4->Branch("tHelmapITS"           ,&tHelmapITS           ,"tHelmapITS/F");
+    fNtuple4->Branch("tHelTPCNcls"          ,&tHelTPCNcls          ,"tHelTPCNcls/F");
+    fNtuple4->Branch("tHelTRDsignal"        ,&tHelTRDsignal        ,"tHelTRDsignal/F");
+    fNtuple4->Branch("tHelxPrimaryVertex"   ,&tHelxPrimaryVertex   ,"tHelxPrimaryVertex/F");
+    fNtuple4->Branch("tHelyPrimaryVertex"   ,&tHelyPrimaryVertex   ,"tHelyPrimaryVertex/F");
+    fNtuple4->Branch("tHelzPrimaryVertex"   ,&tHelzPrimaryVertex   ,"tHelzPrimaryVertex/F");
+    fNtuple4->Branch("tHelchi2PerClusterTPC",&tHelchi2PerClusterTPC,"tHelchi2PerClusterTPC/F");
+    
+
+  } 
+
+  PostData(1,  fListHist);
+  PostData(2,  fNtuple1);
+  PostData(3,  fNtuple4);
+}// end UserCreateOutputObjects
+
+
+//====================== USER EXEC ========================
+
+void AliAnalysisTaskHelium3PiAOD::UserExec(Option_t *) 
+{
+  //_______________________________________________________________________
+  
+  //!*********************!//
+  //!  Define variables   !//
+  //!*********************!//
+
+  Double_t pinTPC=0.,poutTPC=0.,TPCSignal=0.;
+  Double_t xPrimaryVertex=0.,yPrimaryVertex=0.,zPrimaryVertex=0.;
+  Double_t massTOF=0.,timeTOF=0.,trackLenghtTOF=0.,betaTOF=0.;
+
+  ULong_t  status=0;
+  //  ULong_t  statusT=0;
+  ULong_t  statusPi=0;
+
+  Bool_t   isTPC=kFALSE,isTOF=kFALSE,IsHeITSRefit=kFALSE,IsPiITSRefit=kFALSE ;
+
+  Float_t nSigmaNegPion=0.;
+
+  Double_t cutNSigma = 3;
+  Double_t bbtheoM=0.,bbtheo=0.;
+  Double_t zNathashaNeg=0;
+  Double_t zNathashaPos=0;
+  Double_t fPos[3]={0.,0.,0.};
+  Double_t runNumber=0.;
+  //  Double_t evNumber=0.;
+
+  Double_t BCNumber=0.;
+  Double_t OrbitNumber=0.;
+  Double_t PeriodNumber=0.;
+
+  Double_t        Helium3Mass = 2.80839; 
+  Double_t        PionMass    = 0.13957; 
+  // TLORENTZ vectors
+  
+  TLorentzVector  vPion,vHelium,vSum;
+
+  //!----------------------------------------------------------------
+
+  //! A set of very loose parameters for cuts 
+  
+  Double_t fgChi2max=33.;     //! max chi2
+  Double_t fgDNmin=0.05;      //! min imp parameter for the 1st daughter = 500um
+  Double_t fgDCAmax=1.0;      //! max DCA between the daughter tracks in cm
+  Double_t fgCPAmin=0.99;     //! min cosine of V0's pointing angle  
+  //  Double_t fgRmin=0.2;    //! min radius of the fiducial volume //original
+  Double_t fgRmin=0.1;        //! min radius of the fiducial volume = 1 mm 
+  Double_t fgRmax=200.;       //! max radius of the fiducial volume = 2 m
+
+  //------------------------------------------
+
+  // Main loop
+  // Called for EACH event
+
+  AliVEvent *event = InputEvent();
+  if (!event) { Printf("ERROR: Could not retrieve event"); return; }
+    
+  Info("AliAnalysisTaskHelium3PiAOD","Starting UserExec");  
+
+  SetDataType("REAL");
+  
+  // create pointer to event
+
+  AliAODEvent *lAODevent=(AliAODEvent *)InputEvent();  
+  if (!lAODevent) {
+    Printf("ERROR: aod not available");
+    //    fHistLog->Fill(98);
+    return;
+  }
+
+  fHistEventMultiplicity->Fill(0);
+
+  Double_t lMagneticField=lAODevent->GetMagneticField();
+  Int_t TrackNumber = -1;
+
+   
+  //*****************//  
+  //*   Centrality  *//
+  //*****************//
+  AliCentrality *centrality = lAODevent->GetCentrality();
+  Float_t percentile=centrality->GetCentralityPercentile("V0M");
+
+  TrackNumber = lAODevent->GetNumberOfTracks();
+  if (TrackNumber<2) return;  
+
+  fHistTrackMultiplicity->Fill(TrackNumber,percentile); //tracce per evento
+
+  //****************************************
+  
+  // PID
+  
+  AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
+  AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
+  fPIDResponse=inputHandler->GetPIDResponse(); // data member di tipo "const AliPIDResponse *fPIDResponse;"
+  //  cout<<"fPIDResponse "<<fPIDResponse<<endl;
+  //===========================================
+
+  Int_t eventtype=-99;
+  
+  Bool_t isSelectedCentral     = (inputHandler->IsEventSelected() & AliVEvent::kCentral);
+  Bool_t isSelectedSemiCentral = (inputHandler->IsEventSelected() & AliVEvent::kSemiCentral);
+  Bool_t isSelectedMB          = (inputHandler->IsEventSelected() & AliVEvent::kMB);
+  if(isSelectedCentral){
+    fHistEventMultiplicity->Fill(3);
+    fHistTrackMultiplicityCent->Fill(TrackNumber,percentile); 
+    eventtype=1;
+  }
+
+  if(isSelectedSemiCentral){
+    fHistEventMultiplicity->Fill(4);
+    fHistTrackMultiplicitySemiCent->Fill(TrackNumber,percentile); 
+    eventtype=2;
+  }
+
+  if(isSelectedMB){
+    fHistEventMultiplicity->Fill(5);
+    fHistTrackMultiplicityMB->Fill(TrackNumber,percentile); 
+    eventtype=3;
+  }
+  
+  if(isSelectedCentral || isSelectedSemiCentral || isSelectedMB){
+    
+    // ANALISYS
+    
+    // Primary vertex cut
+
+    Double_t lBestPrimaryVtxPos[3]          = {-100.0, -100.0, -100.0};
+
+    AliAODVertex*  vtx= lAODevent->GetPrimaryVertex();
+    if(!vtx){
+      cout<<"No PV"<<endl;
+      return;
+    }
+
+    vtx->GetXYZ( lBestPrimaryVtxPos );
+    
+    fHistEventMultiplicity->Fill(1); // analyzed events with PV
+    xPrimaryVertex=lBestPrimaryVtxPos[0];
+    yPrimaryVertex=lBestPrimaryVtxPos[1];
+    zPrimaryVertex=lBestPrimaryVtxPos[2];  
+    
+    if(TMath::Abs(zPrimaryVertex)>10) return;
+    
+    if(eventtype==1){
+      fHistTrackMultiplicityPVCent->Fill(TrackNumber,percentile); 
+      fHistEventMultiplicity->Fill(6); 
+    }
+    
+    if(eventtype==2){
+      fHistTrackMultiplicityPVSemiCent->Fill(TrackNumber,percentile); 
+      fHistEventMultiplicity->Fill(7); 
+    }
+    
+    if(eventtype==3){
+      fHistTrackMultiplicityPVMB->Fill(TrackNumber,percentile); 
+      fHistEventMultiplicity->Fill(8); 
+    }
+    
+    
+    fHistEventMultiplicity->Fill(2);
+    
+    //Find Pair candidates
+    
+    TArrayI PionsTPC(TrackNumber);        //Neg pions
+    Int_t nPionsTPC=0;
+    
+    TArrayI HeTPC(TrackNumber);        //helium3
+    Int_t nHeTPC=0;
+    
+    const Double_t speedOfLight =  TMath::C()*1E2*1E-12; // cm/ps
+    
+    Float_t impactXY=-999, impactZ=-999;
+    Float_t impactXYpi=-999, impactZpi=-999;
+
+    
+    //*************************************************************
+    
+    runNumber = lAODevent->GetRunNumber();
+    BCNumber    = lAODevent->GetBunchCrossNumber();
+    OrbitNumber = lAODevent->GetOrbitNumber();
+    PeriodNumber= lAODevent->GetPeriodNumber();
+    
+    //*************************************************************
+
+    for (Int_t j=0; j<TrackNumber; j++) { //loop on tracks
+
+      AliVTrack *track = lAODevent->GetTrack(j);
+      AliESDtrack *esdtrack=new AliESDtrack(track);
+      
+      //      AliVTrack*  esdtrack= (AliVTrack *) fEvent->GetTrack(iT);
+
+     
+      if(!esdtrack) { 
+       AliError(Form("ERROR: Could not retrieve esdtrack %d",j)); 
+       continue; 
+      }
+
+      // ************** Track cuts ****************
+      
+      if (!fESDtrackCuts->AcceptTrack(esdtrack)) continue;
+
+      
+      status  = (ULong_t)esdtrack->GetStatus();
+      isTPC   = (((status) & AliESDtrack::kTPCin)  != 0);
+      isTOF   = ((((status) & AliESDtrack::kTOFout) != 0) && (((status) & AliESDtrack::kTIME) != 0));
+
+      
+      UInt_t mapITS=esdtrack->GetITSClusterMap();
+            
+      //----------------------------------------------
+      
+      //****** Cuts from  AliV0Vertex.cxx *************
+      
+      Double_t d=esdtrack->GetD(xPrimaryVertex,yPrimaryVertex,lMagneticField);
+      //    if (TMath::Abs(d)<fgDPmin) continue;
+      if (TMath::Abs(d)>fgRmax) continue;
+      
+      //---- (Usefull) Stuff
+      
+      TPCSignal=esdtrack->GetTPCsignal(); 
+      
+      if (TPCSignal<10)continue;
+      if (TPCSignal>1000)continue;
+      
+      if(!isTPC)continue;
+      if(!esdtrack->GetTPCInnerParam())continue;
+      
+      AliExternalTrackParam trackIn(*esdtrack->GetInnerParam()); 
+      pinTPC= trackIn.GetP(); 
+      
+      //pinTPC= esdtrack->GetTPCMomentum();
+
+      poutTPC=pinTPC;
+      
+      
+      if((status) & (AliESDtrack::kITSrefit!=0)){
+       fhBB->Fill(pinTPC*esdtrack->GetSign(),TPCSignal);
+      }
+      
+      timeTOF=esdtrack->GetTOFsignal();                 // ps
+      trackLenghtTOF= esdtrack->GetIntegratedLength();  // cm
+      
+      if(isTOF){
+       
+       if(!esdtrack->GetOuterParam())continue;    
+       
+       AliExternalTrackParam trackOut(*esdtrack->GetOuterParam()); 
+       
+       poutTPC = trackOut.GetP(); 
+       
+       betaTOF= (trackLenghtTOF/timeTOF)/2.99792458e-2;
+       
+       fhTOF->Fill(poutTPC*esdtrack->GetSign(),betaTOF);
+       
+       Double_t mass2=(poutTPC*poutTPC)*((((speedOfLight*speedOfLight)*(timeTOF*timeTOF))-(trackLenghtTOF*trackLenghtTOF))/(trackLenghtTOF*trackLenghtTOF));
+       if(mass2>0) massTOF=TMath::Sqrt(mass2);
+       fhMassTOF->Fill(massTOF);
+       
+       if(esdtrack->GetSign() < 0.)fBetavsTPCsignalNeg->Fill(betaTOF,TPCSignal);
+       if(esdtrack->GetSign() > 0.)fBetavsTPCsignalPos->Fill(betaTOF,TPCSignal);
+       
+      }
+      
+      //pass2
+      
+      // bbtheo =BetheBloch((2*pinTPC)/3.,2,kTRUE);    //! OK
+      // bbtheoM=(1 - 0.08*5)*bbtheo;                  //! OK 
+      // bbtheoP=(1 + 0.08*5)*bbtheo;                  //! OK
+      
+
+      bbtheo = fPIDResponse->NumberOfSigmas((AliPIDResponse::EDetector)0,esdtrack,(AliPID::EParticleType) 7);
+      
+      if(esdtrack->GetSign()<0){
+       zNathashaNeg=bbtheo;//(TPCSignal-bbtheo)/bbtheo;
+       //      cout<<"BBtheo 1 :"<<zNathashaNeg<<endl;
+       fhNaNeg->Fill(pinTPC,zNathashaNeg); 
+      }
+      
+      if(esdtrack->GetSign() > 0.){
+               zNathashaPos=bbtheo;//(TPCSignal-bbtheo)/bbtheo;
+       fhNaPos->Fill(pinTPC,zNathashaPos); 
+      }
+      
+      nSigmaNegPion=TMath::Abs(fPIDResponse->NumberOfSigmasTPC(esdtrack,(AliPID::EParticleType) 2));
+      //2 is pion
+      
+      if ( (nSigmaNegPion < cutNSigma)){ 
+       
+       //      cout<<"Nsigma pi: "<<nSigmaNegPion<<endl;
+       
+       fhBBPions->Fill(pinTPC*esdtrack->GetSign(),TPCSignal);
+       
+       if(pinTPC<3.){
+         PionsTPC[nPionsTPC++]=j;
+       }
+      }
+    
+      //      nSigmaNegPion=(fPIDResponse->NumberOfSigmasTPC(esdtrack,(AliPID::EParticleType) 2));
+      
+      bbtheoM = TMath::Abs((fPIDResponse->NumberOfSigmasTPC(esdtrack,(AliPID::EParticleType) 7)));
+      
+      //      if( TPCSignal > bbtheoM ) {
+      //      if( bbtheoM > -3.) {
+      if( bbtheoM < 3.) {
+       
+       if(pinTPC>0.6){
+         
+         fhBBHe->Fill(pinTPC*esdtrack->GetSign(),TPCSignal);
+         HeTPC[nHeTPC++]=j;
+         
+         Bool_t isHeITSrefit=((status) & (AliESDtrack::kITSrefit));
+         
+         esdtrack->GetImpactParameters(impactXY, impactZ);
+         
+         Int_t  fIdxInt[200]; //dummy array
+         Int_t nClustersTPC = esdtrack->GetTPCclusters(fIdxInt);
+         
+         Float_t chi2PerClusterTPC = esdtrack->GetTPCchi2()/(Float_t)(nClustersTPC);
+         
+         tHelrunNumber         =(Float_t)runNumber;
+         tHelBCNumber          =(Float_t)BCNumber;
+         tHelOrbitNumber       =(Float_t)OrbitNumber;
+         tHelPeriodNumber      =(Float_t)PeriodNumber;
+         tHeleventtype         =(Float_t)eventtype;
+         tHelisHeITSrefit      =(Float_t)isHeITSrefit;
+         tHelpercentile        =(Float_t)percentile;
+         tHelSign              =(Float_t)esdtrack->GetSign();
+         tHelpinTPC            =(Float_t)pinTPC;
+         tHelGetTPCsignal      =(Float_t)esdtrack->GetTPCsignal();
+         tHelPx                =(Float_t)esdtrack->Px();
+         tHelPy                =(Float_t)esdtrack->Py();
+         tHelPz                =(Float_t)esdtrack->Pz();
+         tHelEta               =(Float_t)esdtrack->Eta();
+         tHelisTOF             =(Float_t)isTOF;
+         tHelpoutTPC           =(Float_t)poutTPC;
+         tHeltimeTOF           =(Float_t)timeTOF;
+         tHeltrackLenghtTOF    =(Float_t)trackLenghtTOF;
+         tHelimpactXY          =(Float_t)impactXY;
+         tHelimpactZ           =(Float_t)impactZ;
+         tHelmapITS            =(Float_t)mapITS;
+         tHelTPCNcls           =(Float_t)esdtrack->GetTPCNcls();
+         tHelTRDsignal         =(Float_t)esdtrack->GetTRDsignal();
+         tHelxPrimaryVertex    =(Float_t)xPrimaryVertex;
+         tHelyPrimaryVertex    =(Float_t)yPrimaryVertex;
+         tHelzPrimaryVertex    =(Float_t)zPrimaryVertex;
+         tHelchi2PerClusterTPC =(Float_t)chi2PerClusterTPC;            
+                 
+         fNtuple4->Fill();
+       }
+      }
+    }  //! track
+         
+    PionsTPC.Set(nPionsTPC);
+    HeTPC.Set(nHeTPC);
+    
+    Double_t        DcaHeToPrimVertex=0;
+    Double_t        DcaPionToPrimVertex=0;
+    
+    impactXY=-999, impactZ=-999;
+    impactXYpi=-999, impactZpi=-999;
+    
+    // Track 
+    
+    // Vettors for il PxPyPz
+    
+    Double_t momPionVett[3];
+    for(Int_t i=0;i<3;i++)momPionVett[i]=0;
+    
+    Double_t momHeVett[3];
+    for(Int_t i=0;i<3;i++)momHeVett[i]=0;
+    
+    //At SV
+    
+    Double_t momPionVettAt[3];
+    for(Int_t i=0;i<3;i++)momPionVettAt[i]=0;
+    
+    Double_t momHeVettAt[3];
+    for(Int_t i=0;i<3;i++)momHeVettAt[i]=0;
+    
+    //---------------   LOOP PAIRS   ----------------
+    
+    for (Int_t k=0; k < nPionsTPC; k++) {                           //! Pions Loop
+      
+      DcaPionToPrimVertex=0.;
+      DcaHeToPrimVertex=0;
+      
+      Int_t PionIdx=PionsTPC[k];
+     
+      AliVTrack *trackpi = (AliVTrack*)lAODevent->GetTrack(PionIdx);
+      AliESDtrack *PionTrack = new AliESDtrack(trackpi); 
+      
+      statusPi = (ULong_t)PionTrack->GetStatus();
+      //      isTOFPi  = ((((statusPi) & (AliESDtrack::kTOFout)) != 0) && (((statusPi) & (AliESDtrack::kTIME)) != 0));
+      IsPiITSRefit = ((statusPi) & (AliESDtrack::kITSrefit)); 
+      
+      if (PionTrack) 
+       DcaPionToPrimVertex = TMath::Abs(PionTrack->GetD(xPrimaryVertex, yPrimaryVertex,lMagneticField)); //OK
+      
+      if(DcaPionToPrimVertex<0.2)continue; 
+      
+      AliExternalTrackParam trackInPion(*PionTrack);  
+      
+      for (Int_t i=0; i<nHeTPC; i++){                               //! Helium Loop
+       
+       Int_t HeIdx=HeTPC[i];
+        
+       AliVTrack *trackhe = (AliVTrack*)lAODevent->GetTrack(HeIdx);
+       AliESDtrack *HeTrack = new AliESDtrack(trackhe);
+       
+       //      statusT= (ULong_t)HeTrack->GetStatus();
+       //      isTOFHe   = (((statusT & AliESDtrack::kTOFout) != 0) && ((statusT & AliESDtrack::kTIME) != 0));
+       IsHeITSRefit = (status & AliESDtrack::kITSrefit); 
+       
+       if (HeTrack) 
+         DcaHeToPrimVertex = TMath::Abs(HeTrack->GetD(xPrimaryVertex, yPrimaryVertex,lMagneticField)); //OK
+       
+       AliExternalTrackParam trackInHe(*HeTrack); 
+    
+       if ( DcaPionToPrimVertex < fgDNmin)                //OK
+         if ( DcaHeToPrimVertex < fgDNmin) continue;    //OK
+       
+       Double_t xn, xp;
+       Double_t dca=0.;
+       
+       dca= PionTrack->GetDCA(HeTrack,lMagneticField,xn,xp); //!dca (Neg to Pos)
+       
+       if (dca > fgDCAmax) continue;
+       if ((xn+xp) > 2*fgRmax) continue;
+       if ((xn+xp) < 2*fgRmin) continue;
+       
+       //CORRECTION from AliV0Vertex
+       
+       Bool_t corrected=kFALSE;
+       if ((trackInPion.GetX() > 3.) && (xn < 3.)) {
+         //correct for the beam pipe material
+         corrected=kTRUE;
+       }
+       if ((trackInHe.GetX() > 3.) && (xp < 3.)) {
+         //correct for the beam pipe material
+         corrected=kTRUE;
+       }
+       if (corrected) {
+         dca=trackInPion.GetDCA(&trackInHe,lMagneticField,xn,xp);
+         if (dca > fgDCAmax) continue;
+         if ((xn+xp) > 2*fgRmax) continue;
+         if ((xn+xp) < 2*fgRmin) continue;
+       }
+       
+       //=============================================//
+       // Make "V0" with found tracks                 //
+       //=============================================//
+       
+       trackInPion.PropagateTo(xn,lMagneticField); 
+       trackInHe.PropagateTo(xp,lMagneticField);
+       
+       AliESDv0 vertex(trackInPion,PionIdx,trackInHe,HeIdx);
+       if (vertex.GetChi2V0() > fgChi2max) continue;
+       
+       Float_t CosPointingAngle=vertex.GetV0CosineOfPointingAngle(xPrimaryVertex,yPrimaryVertex,zPrimaryVertex); //PointingAngle
+       if (CosPointingAngle < fgCPAmin) continue;
+       
+       vertex.SetDcaV0Daughters(dca);
+       vertex.SetV0CosineOfPointingAngle(CosPointingAngle);
+       
+       fPos[0]=vertex.Xv();
+       fPos[1]=vertex.Yv(); 
+       fPos[2]=vertex.Zv(); 
+       
+       HeTrack->PxPyPz(momHeVett);
+       PionTrack->PxPyPz(momPionVett); 
+       
+       Double_t raggio=TMath::Sqrt(fPos[0]*fPos[0]+fPos[1]*fPos[1]+fPos[2]*fPos[2]);
+       HeTrack->GetPxPyPzAt(raggio,lMagneticField,momHeVettAt);
+       PionTrack->GetPxPyPzAt(raggio,lMagneticField,momPionVettAt); 
+       
+       //------------------------------------------------------------------------//
+       
+       HeTrack->GetImpactParameters(impactXY, impactZ);
+       
+       PionTrack->GetImpactParameters(impactXYpi, impactZpi);
+       
+       if(vertex.GetD(xPrimaryVertex,yPrimaryVertex,zPrimaryVertex)>3) continue;
+       
+       //salvo solo fino a 3.1 GeV/c2
+       
+       vHelium.SetXYZM(2*momHeVettAt[0],2*momHeVettAt[1],2*momHeVettAt[2],Helium3Mass); 
+       vPion.SetXYZM(momPionVettAt[0],momPionVettAt[1],momPionVettAt[2],PionMass);       
+       vSum=vHelium+vPion;
+       
+       if(vSum.M()>3.2)
+         continue;
+
+       Int_t  fIdxInt[200]; //dummy array
+       Int_t nClustersTPCHe = HeTrack->GetTPCclusters(fIdxInt);
+       Int_t nClustersTPCPi = PionTrack->GetTPCclusters(fIdxInt);
+       
+       //----------------------------------------------------------------------//
+
+       trunNumber              =(Float_t)runNumber;
+       tbunchcross             =(Float_t)BCNumber;
+       torbit                  =(Float_t)OrbitNumber;
+       tperiod                 =(Float_t)PeriodNumber;
+       teventtype              =(Float_t)eventtype;
+       tTrackNumber            =(Float_t)TrackNumber;
+       tpercentile             =(Float_t)percentile;
+       txPrimaryVertex         =(Float_t)xPrimaryVertex; //PRIMARY
+       tyPrimaryVertex         =(Float_t)yPrimaryVertex;
+       tzPrimaryVertex         =(Float_t)zPrimaryVertex;
+       txSecondaryVertex       =(Float_t)fPos[0]; //SECONDARY
+       tySecondaryVertex       =(Float_t)fPos[1];
+       tzSecondaryVertex       =(Float_t)fPos[2];
+       tdcaTracks              =(Float_t)dca;           //between 2 tracks
+       tCosPointingAngle       =(Float_t)CosPointingAngle;          //cosPointingAngle da V0
+       tDCAV0toPrimaryVertex   =(Float_t)vertex.GetD(xPrimaryVertex,yPrimaryVertex,zPrimaryVertex);
+       tHeSign                 =(Float_t)HeTrack->GetSign(); //He
+       tHepInTPC               =(Float_t)trackInHe.GetP();
+       tHeTPCsignal            =(Float_t)HeTrack->GetTPCsignal();
+       tDcaHeToPrimVertex      =(Float_t)DcaHeToPrimVertex;
+       tHeEta                  =(Float_t)HeTrack->Eta();
+       tmomHex                 =(Float_t)momHeVett[0];
+       tmomHey                 =(Float_t)momHeVett[1];
+       tmomHez                 =(Float_t)momHeVett[2];
+       tmomHeAtSVx             =(Float_t)momHeVettAt[0];
+       tmomHeAtSVy             =(Float_t)momHeVettAt[1];
+       tmomHeAtSVz             =(Float_t)momHeVettAt[2];
+       tHeTPCNcls              =(Float_t)HeTrack->GetTPCNcls();
+       tHeimpactXY             =(Float_t)impactXY;
+       tHeimpactZ              =(Float_t)impactZ;
+       tHeITSClusterMap        =(Float_t)HeTrack->GetITSClusterMap();
+       tIsHeITSRefit           =(Float_t)IsHeITSRefit;
+       tPionSign               =(Float_t)PionTrack->GetSign(); //Pion
+       tPionpInTPC             =(Float_t)trackInPion.GetP();
+       tPionTPCsignal          =(Float_t)PionTrack->GetTPCsignal();
+       tDcaPionToPrimVertex    =(Float_t)DcaPionToPrimVertex;
+       tPionEta                =(Float_t)PionTrack->Eta();
+       tmomPionx               =(Float_t)momPionVett[0];
+       tmomPiony               =(Float_t)momPionVett[1];
+       tmomPionz               =(Float_t)momPionVett[2];
+       tmomNegPionAtSVx        =(Float_t)momPionVettAt[0];
+       tmomNegPionAtSVy        =(Float_t)momPionVettAt[1];
+       tmomNegPionAtSVz        =(Float_t)momPionVettAt[2];
+       tPionTPCNcls            =(Float_t)PionTrack->GetTPCNcls();
+       tPionimpactXY           =(Float_t)impactXYpi;
+       tPionimpactZ            =(Float_t)impactZpi;
+       tPionITSClusterMap      =(Float_t)PionTrack->GetITSClusterMap();
+       tIsPiITSRefit           =(Float_t)IsPiITSRefit;
+       txn                     =(Float_t)xn;
+       txp                     =(Float_t)xp;
+       tchi2He                 =(Float_t)HeTrack->GetTPCchi2()/(Float_t)(nClustersTPCHe);
+       tchi2Pi                 =(Float_t)PionTrack->GetTPCchi2()/(Float_t)(nClustersTPCPi);
+               
+
+       fNtuple1->Fill();  
+       vertex.Delete();
+      }// positive TPC
+      
+    } //negative tpc
+    
+  }
+  
+  PostData(1,fListHist);
+  PostData(2,fNtuple1);
+  PostData(3,fNtuple4);
+  
+} //end userexec
+
+
+//________________________________________________________________________
+
+void AliAnalysisTaskHelium3PiAOD::Terminate(Option_t *) 
+{
+  // Draw result to the screen
+  // Called once at the end of the query
+}
+
+
diff --git a/PWGLF/STRANGENESS/Hypernuclei/AliAnalysisTaskHelium3PiAOD.h b/PWGLF/STRANGENESS/Hypernuclei/AliAnalysisTaskHelium3PiAOD.h
new file mode 100644 (file)
index 0000000..660c305
--- /dev/null
@@ -0,0 +1,165 @@
+#ifndef ALIANALYSISTASKHELIUM3PIAOD_H
+#define ALIANALYSISTASKHELIUM3PIAOD_H
+
+/*  See cxx source for full Copyright notice */
+
+//-----------------------------------------------------------------
+//                 AliAnalysisTaskHelium3Pion class
+//-----------------------------------------------------------------
+
+class TList;
+class TH1F;
+class TH2F;
+class TH3F;
+class TTree;
+class AliESDtrackCuts;
+
+#include <AliPIDResponse.h>
+
+#include "TString.h"
+
+#include "AliAnalysisTaskSE.h"
+
+class AliAnalysisTaskHelium3PiAOD : public AliAnalysisTaskSE {
+ public:
+  AliAnalysisTaskHelium3PiAOD();
+  AliAnalysisTaskHelium3PiAOD(const char *name);
+  virtual ~AliAnalysisTaskHelium3PiAOD();
+  
+  virtual void  UserCreateOutputObjects();
+  virtual void  UserExec(Option_t *option);
+  virtual void  Terminate(Option_t *);
+  
+  void SetCollidingSystems(Short_t collidingSystems = 0)     {fCollidingSystems = collidingSystems;}
+  void SetAnalysisType    (const char* analysisType = "ESD") {fAnalysisType = analysisType;}
+  void SetDataType    (const char* dataType = "REAL") {fDataType = dataType;}
+  
+  Double_t BetheBloch(Double_t bg,Double_t Charge,Bool_t isPbPb);
+
+  //  Bool_t IsTrackAccepted(AliVTrack *track);
+  
+  private:
+  
+  TString fAnalysisType;            //! "ESD" or "AOD" analysis type   
+  
+
+  Short_t fCollidingSystems;        //! 0 = pp collisions or 1 = AA collisions
+  
+  AliESDtrackCuts *fESDtrackCuts; 
+
+  TString fDataType;                //! "REAL" or "SIM" data type      
+
+  TList        *fListHist;                  //! List of  histograms
+
+  TH1F *fHistEventMultiplicity;
+  TH2F *fHistTrackMultiplicity;
+  TH2F *fHistTrackMultiplicityCent;
+  TH2F *fHistTrackMultiplicitySemiCent;
+  TH2F *fHistTrackMultiplicityMB;
+  TH2F *fHistTrackMultiplicityPVCent;
+  TH2F *fHistTrackMultiplicityPVSemiCent;
+  TH2F *fHistTrackMultiplicityPVMB;
+  TH2F *fhBB;
+  TH2F *fhTOF;
+  TH1F *fhMassTOF;
+  TH2F *fhBBPions;
+  TH2F *fhBBHe;
+  TH2F *fhNaPos;
+  TH2F *fhNaNeg;
+  TH2F *fBetavsTPCsignalPos;
+  TH2F *fBetavsTPCsignalNeg;
+   
+  TTree *fNtuple1;                  //! Tree Pairs Pi/Proton "standard"
+  
+  Float_t trunNumber;
+  Float_t tbunchcross;
+  Float_t torbit;
+  Float_t tperiod;
+  Float_t teventtype;
+  Float_t tTrackNumber;
+  Float_t tpercentile;
+  Float_t txPrimaryVertex;
+  Float_t tyPrimaryVertex;
+  Float_t tzPrimaryVertex;
+  Float_t txSecondaryVertex;
+  Float_t tySecondaryVertex;
+  Float_t tzSecondaryVertex;
+  Float_t tdcaTracks;
+  Float_t tCosPointingAngle;
+  Float_t tDCAV0toPrimaryVertex;
+  Float_t tHeSign;
+  Float_t tHepInTPC;
+  Float_t tHeTPCsignal;
+  Float_t tDcaHeToPrimVertex;
+  Float_t tHeEta;
+  Float_t tmomHex;
+  Float_t tmomHey;
+  Float_t tmomHez;
+  Float_t tmomHeAtSVx;
+  Float_t tmomHeAtSVy;
+  Float_t tmomHeAtSVz;
+  Float_t tHeTPCNcls;
+  Float_t tHeimpactXY;
+  Float_t tHeimpactZ;
+  Float_t tHeITSClusterMap;
+  Float_t tIsHeITSRefit;
+  Float_t tPionSign;
+  Float_t tPionpInTPC;
+  Float_t tPionTPCsignal;
+  Float_t tDcaPionToPrimVertex;
+  Float_t tPionEta;
+  Float_t tmomPionx;
+  Float_t tmomPiony;
+  Float_t tmomPionz;
+  Float_t tmomNegPionAtSVx;
+  Float_t tmomNegPionAtSVy;
+  Float_t tmomNegPionAtSVz;
+  Float_t tPionTPCNcls;
+  Float_t tPionimpactXY;
+  Float_t tPionimpactZ;
+  Float_t tPionITSClusterMap;
+  Float_t tIsPiITSRefit;
+  Float_t txn;
+  Float_t txp;
+  Float_t tchi2He;
+  Float_t tchi2Pi;
+  
+  TTree *fNtuple4;                  //! Tree He 
+
+  Float_t tHelrunNumber;
+  Float_t tHelBCNumber;
+  Float_t tHelOrbitNumber;
+  Float_t tHelPeriodNumber;
+  Float_t tHeleventtype;
+  Float_t tHelisHeITSrefit;
+  Float_t tHelpercentile;
+  Float_t tHelSign;
+  Float_t tHelpinTPC;
+  Float_t tHelGetTPCsignal;
+  Float_t tHelPx;
+  Float_t tHelPy;
+  Float_t tHelPz;
+  Float_t tHelEta;
+  Float_t tHelisTOF;
+  Float_t tHelpoutTPC;
+  Float_t tHeltimeTOF;
+  Float_t tHeltrackLenghtTOF;
+  Float_t tHelimpactXY;
+  Float_t tHelimpactZ;
+  Float_t tHelmapITS;
+  Float_t tHelTPCNcls;
+  Float_t tHelTRDsignal;
+  Float_t tHelxPrimaryVertex;
+  Float_t tHelyPrimaryVertex;
+  Float_t tHelzPrimaryVertex;
+  Float_t tHelchi2PerClusterTPC;
+  
+  AliPIDResponse *fPIDResponse;     //! pointer to PID response
+
+  AliAnalysisTaskHelium3PiAOD(const AliAnalysisTaskHelium3PiAOD&);            // not implemented
+  AliAnalysisTaskHelium3PiAOD& operator=(const AliAnalysisTaskHelium3PiAOD&); // not implemented
+  
+  ClassDef(AliAnalysisTaskHelium3PiAOD, 0);
+};
+
+#endif