coding conventions fix
[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
41dc4195 43#include "TMath.h"
6e214c87 44#include "TObjArray.h"
45#include "AliFlowCandidateTrack.h"
46#include "AliFlowEventCuts.h"
47
48#include "AliAnalysisTaskFlowK0Candidates.h"
49
50
51ClassImp(AliAnalysisTaskFlowK0Candidates)
52
53//_____________________________________________________________________________
54AliAnalysisTaskFlowK0Candidates::AliAnalysisTaskFlowK0Candidates() :
55 AliAnalysisTaskSE(),
56 fCutsEvent(NULL),
57 fCuts(NULL),
41dc4195 58 fQAList(NULL),
6e214c87 59 fMassMin(0),
41dc4195 60 fMassMax(0),
61 fDLcut(0),
26baa5fb 62 fEvent(NULL),
63 fMulti(NULL)
6e214c87 64{
26baa5fb 65 //
41dc4195 66 for(int i=0; i!=4; ++i) {
26baa5fb 67 fMass[i] = NULL;
68 fDCA[i] = NULL;
69 fDL[i] = NULL;
70 fCTP[i] = NULL;
71 fd0d0[i] = NULL;
72 fPhi[i] = NULL;
73 fEta[i] = NULL;
74 fPt[i] = NULL;
75 fAPhi[i] = NULL;
76 fAEta[i] = NULL;
77 fAPt[i] = NULL;
78 fBPhi[i] = NULL;
79 fBEta[i] = NULL;
80 fBPt[i] = NULL;
41dc4195 81 }
6e214c87 82}
83
84//_____________________________________________________________________________
85AliAnalysisTaskFlowK0Candidates::AliAnalysisTaskFlowK0Candidates(const char *name, AliFlowEventCuts *cutsEvent, AliESDtrackCuts *cuts, Double_t MassMin, Double_t MassMax) :
86 AliAnalysisTaskSE(name),
87 fCutsEvent(cutsEvent),
88 fCuts(cuts),
41dc4195 89 fQAList(NULL),
6e214c87 90 fMassMin(MassMin),
41dc4195 91 fMassMax(MassMax),
92 fDLcut(0),
26baa5fb 93 fEvent(NULL),
94 fMulti(NULL)
6e214c87 95{
26baa5fb 96 //
41dc4195 97 for(int i=0; i!=4; ++i) {
26baa5fb 98 fMass[i] = NULL;
99 fDCA[i] = NULL;
100 fDL[i] = NULL;
101 fCTP[i] = NULL;
102 fd0d0[i] = NULL;
103 fPhi[i] = NULL;
104 fEta[i] = NULL;
105 fPt[i] = NULL;
106 fAPhi[i] = NULL;
107 fAEta[i] = NULL;
108 fAPt[i] = NULL;
109 fBPhi[i] = NULL;
110 fBEta[i] = NULL;
111 fBPt[i] = NULL;
41dc4195 112 }
6e214c87 113 DefineInput( 0, TChain::Class() );
114 DefineOutput( 1, TObjArray::Class() );
115 DefineOutput( 2, TList::Class() );
116}
117
118//_____________________________________________________________________________
119AliAnalysisTaskFlowK0Candidates::~AliAnalysisTaskFlowK0Candidates()
120{
26baa5fb 121 //
41dc4195 122 if(fQAList) delete fQAList;
6e214c87 123}
124
125//_____________________________________________________________________________
126void AliAnalysisTaskFlowK0Candidates::UserCreateOutputObjects()
127{
26baa5fb 128 //
41dc4195 129 fQAList = new TList();
130 fQAList->SetOwner();
131 AddQAEvents();
132 AddQACandidates();
133 PostData( 2, fQAList);
134}
135//_____________________________________________________________________________
136void AliAnalysisTaskFlowK0Candidates::AddQAEvents()
137{
26baa5fb 138 //
41dc4195 139 TList *tQAEvents = new TList();
140 tQAEvents->SetName("Events");
141 tQAEvents->SetOwner();
26baa5fb 142 fEvent = new TH1D("Event", "Number of Events", 2, 0, 2);
143 tQAEvents->Add(fEvent);
144 fMulti = new TH1D("Multiplicity", "Multiplicity", 180, 0, 10000);
145 tQAEvents->Add(fMulti);
41dc4195 146 fQAList->Add(tQAEvents);
147}
148//_____________________________________________________________________________
149void AliAnalysisTaskFlowK0Candidates::AddQACandidates()
150{
26baa5fb 151 //
41dc4195 152 TList *tQACandidates[4];
153 TList *tQADaughters[4];
154 for(int i=0; i!=4; ++i) {
155 tQACandidates[i] = new TList();
156 tQACandidates[i]->SetOwner();
157 tQACandidates[i]->SetName(Form("Candidates%d",i));
26baa5fb 158 fMass[i] = new TH1D( Form("Mass%i",i), "Mass;M_{#pi#pi} [GeV];Counts per MeV", 180, 0.41, 0.59); tQACandidates[i]->Add( fMass[i] );
159 fDCA[i] = new TH1D( Form("DCA%i" ,i), "DCA;[cm];Counts per 10 um", 180, 0.00, 0.18); tQACandidates[i]->Add( fDCA[i] );
160 fDL[i] = new TH1D( Form("DL%i" ,i), "Decay Length;[cm];Counts per 0.1 mm", 180, 0.00, 1.80); tQACandidates[i]->Add( fDL[i] );
161 fCTP[i] = new TH1D( Form("CTP%i" ,i), "Cos#theta_{p}", 180,-1.10, 1.10); tQACandidates[i]->Add( fCTP[i] );
162 fd0d0[i] = new TH1D( Form("d0d0%i",i), "d_{0}xd_{0};[cm^{2}];Cnts 0.01 mm^{2}",180,-0.009,0.009);tQACandidates[i]->Add( fd0d0[i] );
163 fPhi[i] = new TH1D( Form("Phi%i" ,i), "Phi;[rad];Counts per degree", 180,0,TMath::TwoPi()); tQACandidates[i]->Add( fPhi[i] );
164 fEta[i] = new TH1D( Form("Eta%i" ,i), "Eta;;Counts per 0.04", 180,-3.60, 3.60); tQACandidates[i]->Add( fEta[i] );
165 fPt[i] = new TH1D( Form("Pt%i" ,i), "Pt;[GeV];Counts per 0.1 GeV", 180, 0.00,18.00); tQACandidates[i]->Add( fPt[i] );
41dc4195 166 tQADaughters[i] = new TList();
167 tQADaughters[i]->SetOwner();
168 tQADaughters[i]->SetName(Form("Daughters%d",i));
26baa5fb 169 fAPhi[i] = new TH1D( Form("PhiBef%i",i), "Phi prePropagation;[rad];Cnts per degree",180,0,TMath::TwoPi()); tQADaughters[i]->Add( fAPhi[i] );
170 fAEta[i] = new TH1D( Form("EtaBef%i",i), "Eta prePropagation;;Counts per 0.04" ,180,-3.6,3.6); tQADaughters[i]->Add( fAEta[i] );
171 fAPt[i] = new TH1D( Form("PtBef%i" ,i), "Pt prePropagation;[GeV];Counts per 0.1 GeV",180, 0, 18); tQADaughters[i]->Add( fAPt[i] );
172 fBPhi[i] = new TH1D( Form("PhiAft%i",i), "Phi posPropagation;[rad];Cnts per degree",180,0,TMath::TwoPi()); tQADaughters[i]->Add( fBPhi[i] );
173 fBPhi[i] = new TH1D( Form("EtaAft%i",i), "Eta posPropagation;;Counts per 0.04" ,180,-3.6,3.6); tQADaughters[i]->Add( fBPhi[i] );
174 fBPt[i] = new TH1D( Form("PtAft%i" ,i), "Pt posPropagation;[GeV];Counts per 0.1 GeV",180, 0, 18); tQADaughters[i]->Add( fBPt[i] );
41dc4195 175 tQACandidates[i]->Add(tQADaughters[i]);
176 fQAList->Add(tQACandidates[i]);
177 }
6e214c87 178}
6e214c87 179//_____________________________________________________________________________
180void AliAnalysisTaskFlowK0Candidates::NotifyRun()
181{
26baa5fb 182 //
41dc4195 183// AliESDEvent *fESD = dynamic_cast<AliESDEvent*>(InputEvent());
184// if(!fESD) return;
6e214c87 185}
186
187//_____________________________________________________________________________
188void AliAnalysisTaskFlowK0Candidates::UserExec(Option_t *)
189{
26baa5fb 190 //
41dc4195 191 AliESDEvent *fESD = dynamic_cast<AliESDEvent*>(InputEvent());
192 AliAODEvent *fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
193 if( (!fESD)&&(!fAOD) )
194 return;
195// printf("\nEvent found");
26baa5fb 196 fEvent->Fill( 0 );
6e214c87 197 if(fCutsEvent)
41dc4195 198 if( !(fCutsEvent->IsSelected(InputEvent())) )
199 return;
26baa5fb 200 fEvent->Fill( 1 );
41dc4195 201 if(fESD)
202 ReadFromESD(fESD);
203 else if(fAOD)
204 ReadFromAOD(fAOD);
205}
206//_____________________________________________________________________________
26baa5fb 207void AliAnalysisTaskFlowK0Candidates::ReadFromESD(const AliESDEvent *fESD)
41dc4195 208{
26baa5fb 209 //
210 TObjArray *gPOIselection = new TObjArray();
211 gPOIselection->SetOwner();
6e214c87 212 int nTracks = fESD->GetNumberOfTracks();
26baa5fb 213 fMulti->Fill( nTracks );
41dc4195 214 for(int i=0; i!=nTracks; ++i) {
6e214c87 215 AliESDtrack *ioT = fESD->GetTrack(i);
216 if(fCuts)
217 if( !(fCuts->IsSelected(ioT)) )
218 continue;
41dc4195 219 // printf("\n first particle OK...");
220 for(int j=i+1; j!=nTracks; ++j) {
6e214c87 221 AliESDtrack *joT = fESD->GetTrack(j);
41dc4195 222 if( (ioT->Charge()*joT->Charge()) > 0 )
6e214c87 223 continue;
224 if(fCuts)
225 if( !(fCuts->IsSelected(joT)) )
226 continue;
41dc4195 227 // printf("\n second also...");
6e214c87 228 AliESDtrack *iT = new AliESDtrack(*ioT);
229 AliESDtrack *jT = new AliESDtrack(*joT);
230 // getting distance of closest approach
231 double DCA = iT->PropagateToDCA(jT,fESD->GetMagneticField());
41dc4195 232 // printf(" propagated...");
6e214c87 233 // getting decay length
41dc4195 234 TVector3 vp, vv, vl;
235 vp = TVector3( (iT->Xv()+jT->Xv())/2, (iT->Yv()+jT->Yv())/2, (iT->Zv()+jT->Zv())/2 ); // candidate position
236 vv = TVector3( fESD->GetPrimaryVertex()->GetX(), fESD->GetPrimaryVertex()->GetY(), fESD->GetPrimaryVertex()->GetZ() ); // vertex position
237 vl = vp - vv; // flight line
26baa5fb 238 double gDL;
239 gDL = vl.Mag();
41dc4195 240 // getting cos pointing angle
26baa5fb 241 double gCTP;
41dc4195 242 TVector3 vi, vj, vs;
6e214c87 243 vi = TVector3( iT->Px(), iT->Py(), iT->Pz() );
244 vj = TVector3( jT->Px(), jT->Py(), jT->Pz() );
245 vs = vi + vj;
26baa5fb 246 gCTP = TMath::Cos( vl.Angle(vs) ); // Marta says: "Carlos, it is faster if you calculate this by hand!"
247 double gD0D0;
248 gD0D0 = iT->GetD(vv.X(),vv.Y(),fESD->GetMagneticField())*jT->GetD(vv.X(),vv.Y(),fESD->GetMagneticField());
6e214c87 249 // getting invariant mass
250 double sum12 = iT->P()*iT->P()+jT->P()*jT->P();
251 double pro12 = iT->P()*iT->P()*jT->P()*jT->P();
26baa5fb 252 double gInvMass = 0;
253 gInvMass += TMath::Power(0.13957018,2);
254 gInvMass += sqrt( TMath::Power(0.13957018,4) + TMath::Power(0.13957018,2)*(sum12) + pro12 );
255 gInvMass -= ( iT->Px()*jT->Px()+iT->Py()*jT->Py()+iT->Pz()*jT->Pz() );
256 gInvMass *= 2;
257 gInvMass = sqrt(gInvMass);
41dc4195 258 // filtering
259 int iLevel = 3;
26baa5fb 260 if( gCTP<0.8 )
41dc4195 261 iLevel = 2;
26baa5fb 262 if( gDL<fDLcut ) // 0.5
41dc4195 263 iLevel = 1;
264 if( DCA>0.05 )
265 iLevel = 0;
266 // printf(" candidate at level %d...",iLevel);
267 for(int h=0; h!=iLevel+1; ++h) {
268 // printf(" %d",h);
269 // candidate
26baa5fb 270 fDCA[h]->Fill( DCA );
271 fDL[h]->Fill( gDL );
272 fCTP[h]->Fill( gCTP );
273 fd0d0[h]->Fill( gD0D0 );
274 fMass[h]->Fill( gInvMass );
275 fPhi[h]->Fill( vs.Phi()+TMath::Pi() );
276 fEta[h]->Fill( vs.Eta() );
277 fPt[h]->Fill( vs.Pt() );
41dc4195 278 // daughters
26baa5fb 279 fAPhi[h]->Fill( iT->Phi() );
280 fAEta[h]->Fill( iT->Eta() );
281 fAPt[h]->Fill( iT->Pt() );
282 fAPhi[h]->Fill( jT->Phi() );
283 fAEta[h]->Fill( jT->Eta() );
284 fAPt[h]->Fill( jT->Pt() );
285 fBPhi[h]->Fill( ioT->Phi());
286 fBPhi[h]->Fill( ioT->Eta());
287 fBPt[h]->Fill( ioT->Pt() );
288 fBPhi[h]->Fill( joT->Phi());
289 fBPhi[h]->Fill( joT->Eta());
290 fBPt[h]->Fill( joT->Pt() );
41dc4195 291 }
292 // building candidate
26baa5fb 293 if( (gInvMass>fMassMin) && (gInvMass<fMassMax) && (iLevel==3) ) {
41dc4195 294 AliFlowCandidateTrack *sTrack = new AliFlowCandidateTrack();
26baa5fb 295 sTrack->SetMass( gInvMass );
41dc4195 296 sTrack->SetPt( vs.Pt() );
297 sTrack->SetPhi( vs.Phi()+TMath::Pi() );
298 sTrack->SetEta( vs.Eta() );
299 sTrack->SetCharge( 0 );
300 sTrack->AddDaughter( iT->GetID() );
301 sTrack->AddDaughter( jT->GetID() );
26baa5fb 302 gPOIselection->AddLast( sTrack );
6e214c87 303 }
6e214c87 304 delete iT;
305 delete jT;
41dc4195 306 } // endfor j
307 } // endfor i
26baa5fb 308 PostData( 1, gPOIselection );
41dc4195 309 PostData( 2, fQAList);
310}
311
312//_____________________________________________________________________________
26baa5fb 313void AliAnalysisTaskFlowK0Candidates::ReadFromAOD(const AliAODEvent *fAOD)
41dc4195 314{
26baa5fb 315 //
316 TObjArray *gPOIselection = new TObjArray();
317 gPOIselection->SetOwner();
41dc4195 318 int nTracks = fAOD->GetNumberOfTracks();
26baa5fb 319 fMulti->Fill( nTracks );
41dc4195 320 for(int i=0; i!=fAOD->GetNumberOfV0s(); ++i) {
321 AliAODv0 *iV0 = (AliAODv0*) fAOD->GetV0(i);
322 if ( iV0->GetOnFlyStatus() ) continue;
323 if ( (iV0->ChargeProng(0)*iV0->ChargeProng(1))>0 ) continue;
324 // getting distance of closest approach
325 double DCA = iV0->DcaV0Daughters();
326 // getting decay length
26baa5fb 327 double gDL;
41dc4195 328 double vv[3];
329 vv[0] = fAOD->GetPrimaryVertex()->GetX();
330 vv[1] = fAOD->GetPrimaryVertex()->GetY();
331 vv[2] = fAOD->GetPrimaryVertex()->GetZ();
26baa5fb 332 gDL = iV0->DecayLengthV0(vv);
41dc4195 333 // cos pointing angle
26baa5fb 334 double gCTP;
335 gCTP = iV0->CosPointingAngle(vv);
336 double gD0D0;
337 gD0D0 = iV0->Prodd0d0();
41dc4195 338 // getting invariant mass
26baa5fb 339 double gInvMass = iV0->MassK0Short();
41dc4195 340 // filtering
341 int iLevel = 3;
26baa5fb 342 if( gCTP<0.8 )
41dc4195 343 iLevel = 2;
26baa5fb 344 if( gDL<fDLcut ) // 0.5
41dc4195 345 iLevel = 1;
346 if( DCA>0.05 )
347 iLevel = 0;
348 for(int h=0; h!=iLevel+1; ++h) {
349 // candidate
26baa5fb 350 fDCA[h]->Fill( DCA );
351 fDL[h]->Fill( gDL );
352 fCTP[h]->Fill( gCTP );
353 fd0d0[h]->Fill( gD0D0 );
354 fMass[h]->Fill( gInvMass );
355 fPhi[h]->Fill( iV0->Phi() );
356 fEta[h]->Fill( iV0->Eta() );
357 fPt[h]->Fill( iV0->Pt() );
41dc4195 358 // daughters
26baa5fb 359 fAPhi[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(0) )->Phi() );
360 fAEta[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(0) )->Eta() );
361 fAPt[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(0) )->Pt() );
362 fAPhi[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(1) )->Phi() );
363 fAEta[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(1) )->Eta() );
364 fAPt[h]->Fill( ( (AliAODTrack*) iV0->GetDaughter(1) )->Pt() );
365 fBPhi[h]->Fill( iV0->PhiProng(0) );
366 fBPhi[h]->Fill( iV0->EtaProng(0) );
367 fBPt[h]->Fill( iV0->PtProng(0) );
368 fBPhi[h]->Fill( iV0->PhiProng(1) );
369 fBPhi[h]->Fill( iV0->EtaProng(1) );
370 fBPt[h]->Fill( iV0->PtProng(1) );
6e214c87 371 }
41dc4195 372 AliFlowCandidateTrack *sTrack = new AliFlowCandidateTrack();
26baa5fb 373 sTrack->SetMass( gInvMass );
41dc4195 374 sTrack->SetPt( iV0->Pt() );
375 sTrack->SetPhi( iV0->Phi() );
376 sTrack->SetEta( iV0->Eta() );
377 sTrack->SetCharge( 0 );
378 sTrack->AddDaughter( iV0->GetPosID() );
379 sTrack->AddDaughter( iV0->GetNegID() );
26baa5fb 380 gPOIselection->AddLast( sTrack );
6e214c87 381 }
26baa5fb 382 PostData( 1, gPOIselection );
41dc4195 383 PostData( 2, fQAList);
6e214c87 384}
385
386//_____________________________________________________________________________
387void AliAnalysisTaskFlowK0Candidates::Terminate(Option_t *)
388{
26baa5fb 389 //
6e214c87 390}
391