adding cosmics selection (Marta)
[u/mrichter/AliRoot.git] / PWG4 / JetTasks / AliPWG4CosmicCandidates.cxx
CommitLineData
1f329128 1/**************************************************************************
2 * Copyright(c) 1998-1999, 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// Analysis task looking for cosmic candidates
18// Task checks if particles are back-to-back in eta and phi
19//
20//
21// Author : Marta Verweij - UU - marta.verweij@cern.ch
22//-----------------------------------------------------------------------
23
24#include "TVector3.h"
25#include <iostream>
26#include "TH1.h"
27#include "TH2.h"
28#include "TH3.h"
29#include "TList.h"
30#include "TChain.h"
31
32#include "AliAnalysisTask.h"
33#include "AliAnalysisManager.h"
34
35#include "AliESDEvent.h"
36#include "AliESDInputHandler.h"
37
38#include "AliESDtrack.h"
39#include "AliESDtrackCuts.h"
40#include "AliExternalTrackParam.h"
41
42#include "AliLog.h"
43
44#include "AliPWG4CosmicCandidates.h"
45
ffc1a68e 46//using namespace std; //required for resolving the 'cout' symbol
1f329128 47using namespace std;
1f329128 48
ffc1a68e 49ClassImp(AliPWG4CosmicCandidates)
1f329128 50
51//________________________________________________________________________
52AliPWG4CosmicCandidates::AliPWG4CosmicCandidates()
ffc1a68e 53: AliAnalysisTaskSE(),
54 fTrackCuts(0),
55 fPtMin(5.),
56 fMaxCosmicAngle(0.002),
57 fNEventAll(0),
58 fNEventSel(0),
59 fPtSignedCosmicCandidates(0),
60 fDeltaPtCosmicCandidates(0),
61 fDeltaPhiSumEta(0),
62 fDCAZCosmicCandidates(0),
63 fDCARCosmicCandidates(0),
64 fTheta(0),
65 fThetaZoom(0),
66 fThetaPt1Pt2(0),
67 fThetaPt1Pt2Signed(0),
68 fDeltaPhiSumEtaPt1(0),
69 fDeltaPhiSumEtaPt2(0),
70 fThetaDCAZ1DCAZ2(0),
71 fRisol(0),
72 fRisolTheta(0),
73 fHistListCosmics(0)
1f329128 74{
ffc1a68e 75 //
1f329128 76 // Default constructor
ffc1a68e 77 //
1f329128 78}
79
80//________________________________________________________________________
81AliPWG4CosmicCandidates::AliPWG4CosmicCandidates(const char *name)
82 : AliAnalysisTaskSE(name),
83 fTrackCuts(0),
84 fPtMin(5.),
85 fMaxCosmicAngle(0.002),
ffc1a68e 86 fNEventAll(0),
87 fNEventSel(0),
1f329128 88 fPtSignedCosmicCandidates(0),
89 fDeltaPtCosmicCandidates(0),
90 fDeltaPhiSumEta(0),
91 fDCAZCosmicCandidates(0),
92 fDCARCosmicCandidates(0),
93 fTheta(0),
94 fThetaZoom(0),
95 fThetaPt1Pt2(0),
96 fThetaPt1Pt2Signed(0),
97 fDeltaPhiSumEtaPt1(0),
98 fDeltaPhiSumEtaPt2(0),
99 fThetaDCAZ1DCAZ2(0),
100 fRisol(0),
101 fRisolTheta(0),
102 fHistListCosmics(0)
103{
104 // Constructor. Initialization of Inputs and Outputs
105
106 // Define input and output slots here
107 // Input slot #0 works with a TChain
108 // DefineInput(0, TChain::Class());
109 // Output slot #0 id reserved by the base class for AOD
110 // Output slot #1 writes into a TList
111 DefineOutput(1, TList::Class());
112 // Output slot #2 writes into a AliESDtrackCuts
113 DefineOutput(2, AliESDtrackCuts::Class());
114}
115//________________________________________________________________________
116AliPWG4CosmicCandidates::AliPWG4CosmicCandidates(const AliPWG4CosmicCandidates &res)
117 : AliAnalysisTaskSE(res),
118 fTrackCuts(0),
119 fPtMin(5.),
120 fMaxCosmicAngle(0.002),
ffc1a68e 121 fNEventAll(0),
122 fNEventSel(0),
1f329128 123 fPtSignedCosmicCandidates(0),
124 fDeltaPtCosmicCandidates(0),
125 fDeltaPhiSumEta(0),
126 fDCAZCosmicCandidates(0),
127 fDCARCosmicCandidates(0),
128 fTheta(0),
129 fThetaZoom(0),
130 fThetaPt1Pt2(0),
131 fThetaPt1Pt2Signed(0),
132 fDeltaPhiSumEtaPt1(0),
133 fDeltaPhiSumEtaPt2(0),
134 fThetaDCAZ1DCAZ2(0),
135 fRisol(0),
136 fRisolTheta(0),
137 fHistListCosmics(0)
138{
139 // Dummy copy constructor
140}
141
142//________________________________________________________________________
143AliPWG4CosmicCandidates& AliPWG4CosmicCandidates::operator=(const AliPWG4CosmicCandidates& /*trclass*/)
144{
145 // Dummy assignment operator
146 return *this;
147}
148
149//________________________________________________________________________
150void AliPWG4CosmicCandidates::LocalInit()
151{
152 //
153 // Only called once at beginning
154 //
155 PostData(2,fTrackCuts);
156}
157
158//________________________________________________________________________
159void AliPWG4CosmicCandidates::UserCreateOutputObjects()
160{
161 //Create output objects
162 // Called once
163 AliDebug(2,Form(">> AliPWG4CosmicCandidates::UserCreateOutputObjects \n"));
164
165 Bool_t oldStatus = TH1::AddDirectoryStatus();
166 TH1::AddDirectory(kFALSE);
167
168 OpenFile(1);
169 fHistListCosmics = new TList();
170
171 fNEventAll = new TH1F("fNEventAll","NEventAll",1,-0.5,0.5);
172 fHistListCosmics->Add(fNEventAll);
173 fNEventSel = new TH1F("fNEventSel","NEvent Selected for analysis",1,-0.5,0.5);
174 fHistListCosmics->Add(fNEventSel);
175
176 Float_t fgkPtMin=0.;
177 Float_t fgkPtMax=100.;
178 Int_t fgkNPtBins= (int)(fgkPtMax-fgkPtMin);
179
180 Int_t fgkNPhiBins=18;
181 Float_t kMinPhi = -0.5*TMath::Pi();
182 Float_t kMaxPhi = 3./2.*TMath::Pi();
183
184 Int_t fgkNThetaBins=fgkNPhiBins*8;
185 Float_t kMinTheta = -0.5*TMath::Pi();
186 Float_t kMaxTheta = 3./2.*TMath::Pi();
187
188 Int_t fgkNDCARBins=80;
189 Float_t fgkDCARMin = -0.2;
190 Float_t fgkDCARMax = 0.2;
60829498 191 Float_t *binsDCAR=new Float_t[fgkNDCARBins+1];
192 for(Int_t i=0; i<=fgkNDCARBins; i++) binsDCAR[i]=(Float_t)fgkDCARMin + (fgkDCARMax-fgkDCARMin)/fgkNDCARBins*(Float_t)i ;
1f329128 193
194 Int_t fgkNDCAZBins=80;
195 Float_t fgkDCAZMin = -2.;
196 Float_t fgkDCAZMax = 2.;
60829498 197 Float_t *binsDCAZ=new Float_t[fgkNDCAZBins+1];
198 for(Int_t i=0; i<=fgkNDCAZBins; i++) binsDCAZ[i]=(Float_t)fgkDCAZMin + (fgkDCAZMax-fgkDCAZMin)/fgkNDCAZBins*(Float_t)i ;
1f329128 199
200 fPtSignedCosmicCandidates = new TH1F("fPtSignedCosmicCandidates","fPtSignedCosmicCandidates",2*(int)(fgkPtMax-fgkPtMin), -1.*fgkPtMax, fgkPtMax);
201 fHistListCosmics->Add(fPtSignedCosmicCandidates);
202
203 fDeltaPtCosmicCandidates = new TH1F("fDeltaPtCosmicCandidates","fDeltaPtCosmicCandidates",fgkNPtBins, -50., 50.);
204 fHistListCosmics->Add(fDeltaPtCosmicCandidates);
205
206 fDeltaPhiSumEta = new TH2F("fDeltaPhiSumEta","fDeltaPhiSumEta",fgkNPhiBins*4,kMinPhi,kMaxPhi,80, -2.,2.);
207 fHistListCosmics->Add(fDeltaPhiSumEta);
208
209 fDCAZCosmicCandidates = new TH2F("fDCAZCosmicCandidates","fDCAZCosmicCandidates",fgkNDCAZBins,binsDCAZ,fgkNDCAZBins,binsDCAZ);
210 fHistListCosmics->Add(fDCAZCosmicCandidates);
211
212 fDCARCosmicCandidates = new TH2F("fDCARCosmicCandidates","fDCARCosmicCandidates",fgkNDCARBins,binsDCAR,fgkNDCARBins,binsDCAR);
213 fHistListCosmics->Add(fDCARCosmicCandidates);
214
215 fTheta = new TH1F("fTheta","fTheta",fgkNThetaBins,kMinTheta,kMaxTheta);
216 fHistListCosmics->Add(fTheta);
217
218 fThetaZoom = new TH1F("fThetaZoom","fThetaZoom",100,TMath::Pi()-1.,TMath::Pi()+1.);
219 fHistListCosmics->Add(fThetaZoom);
220
221 fThetaPt1Pt2 = new TH3F("fThetaPt1Pt2","fThetaPt1Pt2",fgkNThetaBins,kMinTheta,kMaxTheta,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax);
222 fHistListCosmics->Add(fThetaPt1Pt2);
223
224 fThetaPt1Pt2Signed = new TH3F("fThetaPt1Pt2Signed","fThetaPt1Pt2Signed",fgkNThetaBins,kMinTheta,kMaxTheta,4*(int)(fgkPtMax-fgkPtMin),-1.*fgkPtMax,fgkPtMax,4*(int)(fgkPtMax-fgkPtMin),-1.*fgkPtMax,fgkPtMax);
225 fHistListCosmics->Add(fThetaPt1Pt2Signed);
226
227 fDeltaPhiSumEtaPt1 = new TH3F("fDeltaPhiSumEtaPt1","fDeltaPhiSumEtaPt1",fgkNThetaBins,kMinTheta,kMaxTheta,80, -2.,2.,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax);
228 fHistListCosmics->Add(fDeltaPhiSumEtaPt1);
229
230 fDeltaPhiSumEtaPt2 = new TH3F("fDeltaPhiSumEtaPt2","fDeltaPhiSumEtaPt2",fgkNThetaBins,kMinTheta,kMaxTheta,80, -2.,2.,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax);
231 fHistListCosmics->Add(fDeltaPhiSumEtaPt2);
232
233 fThetaDCAZ1DCAZ2 = new TH3F("fThetaDCAZ1DCAZ2","fThetaDCAZ1DCAZ2",fgkNThetaBins,kMinTheta,kMaxTheta,fgkNDCAZBins,-2.,2.,fgkNDCAZBins,-2.,2.);
234 fHistListCosmics->Add(fThetaDCAZ1DCAZ2);
235
236 fRisol = new TH1F("fRisol","fRisol",100,0.,10.);
237 fHistListCosmics->Add(fRisol);
238
239 fRisolTheta = new TH2F("fRisolTheta","fRisolTheta",100,0.,10.,fgkNThetaBins,kMinTheta,kMaxTheta);
240 fHistListCosmics->Add(fRisolTheta);
241
242 TH1::AddDirectory(oldStatus);
243
244 if(binsDCAR) delete [] binsDCAR;
245 if(binsDCAZ) delete [] binsDCAZ;
246
247}
248
249//________________________________________________________________________
250void AliPWG4CosmicCandidates::UserExec(Option_t *)
251{
252// // Main loop
253// // Called for each event
254
255 // All events without selection
256 fNEventAll->Fill(0.);
257
258 if (!fInputEvent) {
259 AliDebug(2,Form("ERROR: fESD not available"));
260 cout << "ERROR: fESD not available" << endl;
261 PostData(1, fHistListCosmics);
262 return;
263 }
264
265 const AliVVertex *vtx = fInputEvent->GetPrimaryVertex();
266 // Need vertex cut
267 if (!vtx || vtx->GetNContributors() < 2) {
268 // Post output data
269 PostData(1, fHistListCosmics);
270 return;
271 }
272
273 // AliDebug(2,Form("Vertex title %s, status %d, nCont %d\n",vtx->GetTitle(), vtx->GetStatus(), vtx->GetNContributors()));
274 double primVtx[3];
275 vtx->GetXYZ(primVtx);
276 if(TMath::Sqrt(primVtx[0]*primVtx[0] + primVtx[1]*primVtx[1])>1. || TMath::Abs(primVtx[2]>10.)){
277 // Post output data
278 PostData(1, fHistListCosmics);
279 return;
280 }
281 if(!fInputEvent->GetNumberOfTracks() || fInputEvent->GetNumberOfTracks()<2){
282 // Post output data
283 PostData(1, fHistListCosmics);
284 return;
285 }
286 Int_t nTracks = fInputEvent->GetNumberOfTracks();
287
288 if(!fTrackCuts) {
289 // Post output data
290 PostData(1, fHistListCosmics);
291 return;
292 }
293
294 fNEventSel->Fill(0.);
295
296 Float_t dcaR[2] = {0.,0.};
297 Float_t dcaZ[2] = {0.,0.};
298
299 // Track loop to fill a pT spectrum
300 for (Int_t iTrack1 = 0; iTrack1 < nTracks; iTrack1++) {
301
302 AliESDtrack* track1 = (AliESDtrack*)fInputEvent->GetTrack(iTrack1);
303 if (!track1) continue;
304 if(!(fTrackCuts->AcceptTrack(track1))) { continue; }
305 if(track1->Pt()<fPtMin) continue;
306 //Start 2nd track loop to look for correlations
307 for (Int_t iTrack2 = iTrack1+1; iTrack2 < nTracks; iTrack2++) {
308 AliESDtrack *track2 = (AliESDtrack*)fInputEvent->GetTrack(iTrack2);
309 if(!track2) continue;
310 if(!(fTrackCuts->AcceptTrack(track2))) { continue; }
311
312 //Check if back-to-back
313 Double_t mom1[3],mom2[3];
314 track1->GetPxPyPz(mom1);
315 track2->GetPxPyPz(mom2);
316 // Double_t cosTheta = (mom1[0]*mom2[0]+mom1[1]*mom2[1]+mom1[2]*mom2[2])/( TMath::Sqrt(mom1[0]*mom1[0]+mom1[1]*mom1[1]+mom1[2]*mom1[2])*TMath::Sqrt(mom2[0]*mom2[0]+mom2[1]*mom2[1]+mom2[2]*mom2[2]) );
317 TVector3 momv1(mom1[0],mom1[1],mom1[2]);
318 TVector3 momv2(mom2[0],mom2[1],mom2[2]);
319 //Double_t theta = momv1.Angle(momv2);
320 Double_t theta = momv1.Phi()-momv2.Phi();
321 if(theta<-0.5*TMath::Pi()) theta+=2.*TMath::Pi();
322
323 fDeltaPtCosmicCandidates->Fill(track1->Pt()-track2->Pt());
324 Float_t deltaPhi = track1->Phi()-track2->Phi();
325 if(deltaPhi<-0.5*TMath::Pi()) deltaPhi+=2.*TMath::Pi();
326 fDeltaPhiSumEta->Fill(deltaPhi,track1->Eta()+track2->Eta());
327
328 track1->GetImpactParameters(dcaR[0],dcaZ[0]);
329 track2->GetImpactParameters(dcaR[1],dcaZ[1]);
330
331 if(track2->Pt()<0.5) continue;
332 Double_t rIsol = TMath::Sqrt( deltaPhi*deltaPhi+(track1->Eta()-track2->Eta())*(track1->Eta()-track2->Eta()) );
333 fRisol->Fill(rIsol); //Fill R histogram
334 if(track2->Pt()<fPtMin) continue;
335
336 fTheta->Fill(theta);
337 fThetaZoom->Fill(theta);
338 fThetaPt1Pt2->Fill(theta,track1->Pt(),track2->Pt());
339 fThetaPt1Pt2Signed->Fill(theta,track1->GetSign()*track1->Pt(),track2->GetSign()*track2->Pt());
340 fDeltaPhiSumEtaPt1->Fill(deltaPhi,track1->Eta()+track2->Eta(),track1->Pt());
341 fDeltaPhiSumEtaPt2->Fill(deltaPhi,track1->Eta()+track2->Eta(),track2->Pt());
342 fThetaDCAZ1DCAZ2->Fill(theta,dcaZ[0],dcaZ[1]);
343 fRisolTheta->Fill(rIsol,theta);
344 if(TMath::Abs(TMath::Pi()-theta)<fMaxCosmicAngle) {
345 fDCAZCosmicCandidates->Fill(dcaZ[0],dcaZ[1]);
346 fDCARCosmicCandidates->Fill(dcaR[0],dcaR[1]);
347 fPtSignedCosmicCandidates->Fill(track1->GetSign()*track1->Pt());
348 }
349
350 } // track2 loop
351
352 } // track1 loop
353
354 PostData(1,fHistListCosmics);
355}
356
357//________________________________________________________________________
358void AliPWG4CosmicCandidates::Terminate(Option_t *)
359{
360 // Called once at the end of the query
361
362
363}