modifications and additional classes for nonprimary particle flow
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowTasks / AliAnalysisTaskFlowK0Candidates.cxx
CommitLineData
6e214c87 1/*************************************************************************
2* Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
3* *
4* Author: The ALICE Off-line Project. *
5* Contributors are mentioned in the code where appropriate. *
6* *
7* Permission to use, copy, modify and distribute this software and its *
8* documentation strictly for non-commercial purposes is hereby granted *
9* without fee, provided that the above copyright notice appears in all *
10* copies and that both the copyright notice and this permission notice *
11* appear in the supporting documentation. The authors make no claims *
12* about the suitability of this software for any purpose. It is *
13* provided "as is" without express or implied warranty. *
14**************************************************************************/
15
16/////////////////////////////////////////////////////
17// AliAnalysisTaskFlowK0Candidates:
18// Analysis task to select K0 candidates for flow analysis.
19// Uses one AliESDtrackCuts object for both daughters and
20// QA histograms to monitor the reconstruction.
21// Author: Carlos Perez (cperez@cern.ch)
22//////////////////////////////////////////////////////
23
24#include "TChain.h"
25#include "TList.h"
26#include "TH1D.h"
27#include "TVector3.h"
28
29#include "AliAnalysisTaskSE.h"
30#include "AliAnalysisManager.h"
31
32#include "AliESDEvent.h"
33#include "AliESDInputHandler.h"
34#include "AliESDtrack.h"
35#include "AliESDtrackCuts.h"
36
37#include "AliAODEvent.h"
38#include "AliAODInputHandler.h"
39#include "AliAODTrack.h"
40
41#include "AliCFManager.h"
42
43#include "TObjArray.h"
44#include "AliFlowCandidateTrack.h"
45#include "AliFlowEventCuts.h"
46
47#include "AliAnalysisTaskFlowK0Candidates.h"
48
49
50ClassImp(AliAnalysisTaskFlowK0Candidates)
51
52//_____________________________________________________________________________
53AliAnalysisTaskFlowK0Candidates::AliAnalysisTaskFlowK0Candidates() :
54 AliAnalysisTaskSE(),
55 fCutsEvent(NULL),
56 fCuts(NULL),
57 fList(NULL),
58 fMass(NULL),
59 fDCA(NULL),
60 fDL(NULL),
61 fCTP(NULL),
62 fMassF(NULL),
63 fDCAF(NULL),
64 fDLF(NULL),
65 fCTPF(NULL),
66 fMassMin(0),
67 fMassMax(0)
68{
69}
70
71//_____________________________________________________________________________
72AliAnalysisTaskFlowK0Candidates::AliAnalysisTaskFlowK0Candidates(const char *name, AliFlowEventCuts *cutsEvent, AliESDtrackCuts *cuts, Double_t MassMin, Double_t MassMax) :
73 AliAnalysisTaskSE(name),
74 fCutsEvent(cutsEvent),
75 fCuts(cuts),
76 fList(NULL),
77 fMass(NULL),
78 fDCA(NULL),
79 fDL(NULL),
80 fCTP(NULL),
81 fMassF(NULL),
82 fDCAF(NULL),
83 fDLF(NULL),
84 fCTPF(NULL),
85 fMassMin(MassMin),
86 fMassMax(MassMax)
87{
88 DefineInput( 0, TChain::Class() );
89 DefineOutput( 1, TObjArray::Class() );
90 DefineOutput( 2, TList::Class() );
91}
92
93//_____________________________________________________________________________
94AliAnalysisTaskFlowK0Candidates::~AliAnalysisTaskFlowK0Candidates()
95{
96 if(fList) {
97 fList->Delete();
98 delete fList;
99 }
100}
101
102//_____________________________________________________________________________
103void AliAnalysisTaskFlowK0Candidates::UserCreateOutputObjects()
104{
105 fList = new TList();
106 fList->SetOwner();
107 fMass = new TH1D( "fMass","Mass;M_{#pi#pi} [GeV];Counts per MeV", 180, 0.41, 0.59); fList->Add( fMass );
108 fDCA = new TH1D( "fDCA" ,"DCA;[cm];Counts per 10 um", 180, 0.00, 0.18); fList->Add( fDCA );
109 fDL = new TH1D( "fDL" ,"Decay Length;[cm];Counts per 0.1 mm", 180, 0.00, 1.80); fList->Add( fDL );
110 fCTP = new TH1D( "fCTP" ,"Cos#theta_{p}", 180,-1.10, 1.10); fList->Add( fCTP );
111
112 fMassF= new TH1D( "fMassF","Mass after Cuts;M_{#pi#pi} [GeV];Counts per MeV", 180, 0.41, 0.59); fList->Add( fMassF );
113 fDCAF = new TH1D( "fDCAF" ,"DCA after Cuts;[cm];Counts per 10 um", 180, 0.00, 0.18); fList->Add( fDCAF );
114 fDLF = new TH1D( "fDLF" ,"Decay Length after Cuts;[cm];Counts per 0.1 mm", 180, 0.00, 1.80); fList->Add( fDLF );
115 fCTPF = new TH1D( "fCTPF" ,"Cos#theta_{p} after Cuts", 180,-1.10, 1.10); fList->Add( fCTPF );
116
117 PostData( 2, fList);
118}
119
120
121//_____________________________________________________________________________
122void AliAnalysisTaskFlowK0Candidates::NotifyRun()
123{
124 AliESDEvent *fESD = dynamic_cast<AliESDEvent*>(InputEvent());
125 if(!fESD) return;
126}
127
128//_____________________________________________________________________________
129void AliAnalysisTaskFlowK0Candidates::UserExec(Option_t *)
130{
131 AliESDEvent *fESD = dynamic_cast<AliESDEvent*>(fInputEvent);
132 if(!fESD) return;
133 if(fCutsEvent)
134 if( !(fCutsEvent->IsSelected(InputEvent())) ) return;
135
136 TObjArray *POIselection = new TObjArray();
137 POIselection->SetOwner();
138
139 int nTracks = fESD->GetNumberOfTracks();
140 for(int i=0; i!=nTracks; ++i) { // first particle
141 AliESDtrack *ioT = fESD->GetTrack(i);
142 if(fCuts)
143 if( !(fCuts->IsSelected(ioT)) )
144 continue;
145 for(int j=i+1; j!=nTracks; ++j) { // second particle
146 AliESDtrack *joT = fESD->GetTrack(j);
147 if( (ioT->Charge()*joT->Charge()) > 0 ) // only keeping opposite signs
148 continue;
149 if(fCuts)
150 if( !(fCuts->IsSelected(joT)) )
151 continue;
152 AliESDtrack *iT = new AliESDtrack(*ioT);
153 AliESDtrack *jT = new AliESDtrack(*joT);
154 // getting distance of closest approach
155 double DCA = iT->PropagateToDCA(jT,fESD->GetMagneticField());
156 fDCA->Fill( DCA );
157 // getting decay length
158 double DL;
159 double vx = fESD->GetPrimaryVertex()->GetX();
160 double vy = fESD->GetPrimaryVertex()->GetY();
161 double vz = fESD->GetPrimaryVertex()->GetZ();
162 double px = (iT->Xv()+jT->Xv())/2;
163 double py = (iT->Yv()+jT->Yv())/2;
164 double pz = (iT->Zv()+jT->Zv())/2;
165 DL = (vx-px)*(vx-px);
166 DL += (vy-py)*(vy-py);
167 DL += (vz-pz)*(vz-pz);
168 DL = sqrt(DL);
169 fDL->Fill( DL );
170 // cos pointing angle
171 double CTP;
172 TVector3 vi, vj, vs, vp;
173 vi = TVector3( iT->Px(), iT->Py(), iT->Pz() );
174 vj = TVector3( jT->Px(), jT->Py(), jT->Pz() );
175 vs = vi + vj;
176 vp = TVector3( px, py, pz );
177 CTP = TMath::Cos( vp.Angle(vs) );
178 fCTP->Fill(CTP);
179
180 // getting invariant mass
181 double sum12 = iT->P()*iT->P()+jT->P()*jT->P();
182 double pro12 = iT->P()*iT->P()*jT->P()*jT->P();
183 double InvMass = 0;
184 InvMass += TMath::Power(0.13957018,2);
185 InvMass += sqrt( TMath::Power(0.13957018,4) + TMath::Power(0.13957018,2)*(sum12) + pro12 );
186 InvMass -= ( iT->Px()*jT->Px()+iT->Py()*jT->Py()+iT->Pz()*jT->Pz() );
187 InvMass *= 2;
188 InvMass = sqrt(InvMass);
189 fMass->Fill( InvMass );
190
191 // evaluating cuts
192 bool bW = (InvMass>fMassMin) && (InvMass<fMassMax);
193 bool bCTP = (fabs(CTP)>0.95);
194 bool bDCA = (DCA<0.05);
195 bool bDL = (DL>2.5);
196 if( (!bW) || (!bDCA) || (!bDL) || (!bCTP) ) {
197 delete iT;
198 delete jT;
199 continue;
200 }
201 fDCAF->Fill( DCA );
202 fDLF->Fill( DL );
203 fCTPF->Fill(CTP);
204 fMassF->Fill( InvMass );
205
206 AliFlowCandidateTrack *sTrack = new AliFlowCandidateTrack();
207 //booking pt of selected candidates
208 sTrack->SetMass( InvMass );
209 sTrack->SetPt( vs.Pt() );
210 sTrack->SetPhi( vs.Phi() );
211 sTrack->SetEta( vs.Eta() );
212 sTrack->SetCharge( 0 );
213 sTrack->AddDaughter( i );
214 sTrack->AddDaughter( j );
215 //saving selected candidate
216 POIselection->AddLast( sTrack );
217 delete iT;
218 delete jT;
219 }
220 }
221 PostData( 1, POIselection );
222 PostData( 2, fList);
223}
224
225//_____________________________________________________________________________
226void AliAnalysisTaskFlowK0Candidates::Terminate(Option_t *)
227{
228}
229