Fixes for coverity.
[u/mrichter/AliRoot.git] / PWGJE / 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//-----------------------------------------------------------------------
9b58297c 17// Analysis task looking for cosmic candidates embedded in a p-p event
1f329128 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{
9b58297c 139 //
140 // Dummy copy constructor
141 //
1f329128 142}
143
144//________________________________________________________________________
145AliPWG4CosmicCandidates& AliPWG4CosmicCandidates::operator=(const AliPWG4CosmicCandidates& /*trclass*/)
146{
9b58297c 147 //
148 // Dummy assignment operator
149 //
150 return *this;
1f329128 151}
152
153//________________________________________________________________________
154void AliPWG4CosmicCandidates::LocalInit()
155{
156 //
157 // Only called once at beginning
158 //
159 PostData(2,fTrackCuts);
160}
161
162//________________________________________________________________________
163void AliPWG4CosmicCandidates::UserCreateOutputObjects()
164{
9b58297c 165 //
166 // Create output objects
1f329128 167 // Called once
168 AliDebug(2,Form(">> AliPWG4CosmicCandidates::UserCreateOutputObjects \n"));
169
170 Bool_t oldStatus = TH1::AddDirectoryStatus();
171 TH1::AddDirectory(kFALSE);
172
173 OpenFile(1);
174 fHistListCosmics = new TList();
175
176 fNEventAll = new TH1F("fNEventAll","NEventAll",1,-0.5,0.5);
177 fHistListCosmics->Add(fNEventAll);
178 fNEventSel = new TH1F("fNEventSel","NEvent Selected for analysis",1,-0.5,0.5);
179 fHistListCosmics->Add(fNEventSel);
180
181 Float_t fgkPtMin=0.;
182 Float_t fgkPtMax=100.;
183 Int_t fgkNPtBins= (int)(fgkPtMax-fgkPtMin);
21a0adae 184 Int_t fgkNPtBins2D= (int)((fgkPtMax-fgkPtMin)/4.);
1f329128 185
186 Int_t fgkNPhiBins=18;
187 Float_t kMinPhi = -0.5*TMath::Pi();
188 Float_t kMaxPhi = 3./2.*TMath::Pi();
189
190 Int_t fgkNThetaBins=fgkNPhiBins*8;
191 Float_t kMinTheta = -0.5*TMath::Pi();
192 Float_t kMaxTheta = 3./2.*TMath::Pi();
193
21a0adae 194 Int_t fgkNDCARBins=40;
1f329128 195 Float_t fgkDCARMin = -0.2;
196 Float_t fgkDCARMax = 0.2;
60829498 197 Float_t *binsDCAR=new Float_t[fgkNDCARBins+1];
198 for(Int_t i=0; i<=fgkNDCARBins; i++) binsDCAR[i]=(Float_t)fgkDCARMin + (fgkDCARMax-fgkDCARMin)/fgkNDCARBins*(Float_t)i ;
1f329128 199
21a0adae 200 Int_t fgkNDCAZBins=40;
1f329128 201 Float_t fgkDCAZMin = -2.;
202 Float_t fgkDCAZMax = 2.;
60829498 203 Float_t *binsDCAZ=new Float_t[fgkNDCAZBins+1];
204 for(Int_t i=0; i<=fgkNDCAZBins; i++) binsDCAZ[i]=(Float_t)fgkDCAZMin + (fgkDCAZMax-fgkDCAZMin)/fgkNDCAZBins*(Float_t)i ;
1f329128 205
21a0adae 206 fPtSignedCosmicCandidates = new TH1F("fPtSignedCosmicCandidates","fPtSignedCosmicCandidates",2*fgkNPtBins, -1.*fgkPtMax, fgkPtMax);
1f329128 207 fHistListCosmics->Add(fPtSignedCosmicCandidates);
208
209 fDeltaPtCosmicCandidates = new TH1F("fDeltaPtCosmicCandidates","fDeltaPtCosmicCandidates",fgkNPtBins, -50., 50.);
210 fHistListCosmics->Add(fDeltaPtCosmicCandidates);
211
212 fDeltaPhiSumEta = new TH2F("fDeltaPhiSumEta","fDeltaPhiSumEta",fgkNPhiBins*4,kMinPhi,kMaxPhi,80, -2.,2.);
213 fHistListCosmics->Add(fDeltaPhiSumEta);
214
215 fDCAZCosmicCandidates = new TH2F("fDCAZCosmicCandidates","fDCAZCosmicCandidates",fgkNDCAZBins,binsDCAZ,fgkNDCAZBins,binsDCAZ);
216 fHistListCosmics->Add(fDCAZCosmicCandidates);
217
218 fDCARCosmicCandidates = new TH2F("fDCARCosmicCandidates","fDCARCosmicCandidates",fgkNDCARBins,binsDCAR,fgkNDCARBins,binsDCAR);
219 fHistListCosmics->Add(fDCARCosmicCandidates);
220
221 fTheta = new TH1F("fTheta","fTheta",fgkNThetaBins,kMinTheta,kMaxTheta);
222 fHistListCosmics->Add(fTheta);
223
224 fThetaZoom = new TH1F("fThetaZoom","fThetaZoom",100,TMath::Pi()-1.,TMath::Pi()+1.);
225 fHistListCosmics->Add(fThetaZoom);
226
21a0adae 227 fThetaPt1Pt2 = new TH3F("fThetaPt1Pt2","fThetaPt1Pt2",fgkNThetaBins,kMinTheta,kMaxTheta,fgkNPtBins2D,fgkPtMin,fgkPtMax,fgkNPtBins2D,fgkPtMin,fgkPtMax);
1f329128 228 fHistListCosmics->Add(fThetaPt1Pt2);
229
21a0adae 230 fThetaPt1Pt2Signed = new TH3F("fThetaPt1Pt2Signed","fThetaPt1Pt2Signed",fgkNThetaBins,kMinTheta,kMaxTheta,fgkNPtBins2D,-1.*fgkPtMax,fgkPtMax,fgkNPtBins2D,-1.*fgkPtMax,fgkPtMax);
1f329128 231 fHistListCosmics->Add(fThetaPt1Pt2Signed);
232
21a0adae 233 fDeltaPhiSumEtaPt1 = new TH3F("fDeltaPhiSumEtaPt1","fDeltaPhiSumEtaPt1",fgkNThetaBins,kMinTheta,kMaxTheta,80, -2.,2.,fgkNPtBins2D,fgkPtMin,fgkPtMax);
1f329128 234 fHistListCosmics->Add(fDeltaPhiSumEtaPt1);
235
21a0adae 236 fDeltaPhiSumEtaPt2 = new TH3F("fDeltaPhiSumEtaPt2","fDeltaPhiSumEtaPt2",fgkNThetaBins,kMinTheta,kMaxTheta,80, -2.,2.,fgkNPtBins2D,fgkPtMin,fgkPtMax);
1f329128 237 fHistListCosmics->Add(fDeltaPhiSumEtaPt2);
238
239 fThetaDCAZ1DCAZ2 = new TH3F("fThetaDCAZ1DCAZ2","fThetaDCAZ1DCAZ2",fgkNThetaBins,kMinTheta,kMaxTheta,fgkNDCAZBins,-2.,2.,fgkNDCAZBins,-2.,2.);
240 fHistListCosmics->Add(fThetaDCAZ1DCAZ2);
241
242 fRisol = new TH1F("fRisol","fRisol",100,0.,10.);
243 fHistListCosmics->Add(fRisol);
244
245 fRisolTheta = new TH2F("fRisolTheta","fRisolTheta",100,0.,10.,fgkNThetaBins,kMinTheta,kMaxTheta);
246 fHistListCosmics->Add(fRisolTheta);
247
11245558 248 TH1::AddDirectory(oldStatus);
249
250 PostData(1, fHistListCosmics);
1f329128 251
252 if(binsDCAR) delete [] binsDCAR;
253 if(binsDCAZ) delete [] binsDCAZ;
254
255}
256
257//________________________________________________________________________
258void AliPWG4CosmicCandidates::UserExec(Option_t *)
259{
9b58297c 260 // Main loop
261 // Called for each event
262 //
1f329128 263
264 // All events without selection
265 fNEventAll->Fill(0.);
266
267 if (!fInputEvent) {
268 AliDebug(2,Form("ERROR: fESD not available"));
269 cout << "ERROR: fESD not available" << endl;
270 PostData(1, fHistListCosmics);
271 return;
272 }
273
274 const AliVVertex *vtx = fInputEvent->GetPrimaryVertex();
b85b355c 275 if(!vtx){
276 PostData(1, fHistListCosmics);
277 return;
278 }
279
280// Need vertex cut
75946d5d 281 TString vtxName(vtx->GetName());
b85b355c 282 if( vtx->GetNContributors() < 2 || (vtxName.Contains("TPCVertex")) ) {
1f329128 283 // Post output data
284 PostData(1, fHistListCosmics);
285 return;
286 }
287
288 // AliDebug(2,Form("Vertex title %s, status %d, nCont %d\n",vtx->GetTitle(), vtx->GetStatus(), vtx->GetNContributors()));
289 double primVtx[3];
290 vtx->GetXYZ(primVtx);
291 if(TMath::Sqrt(primVtx[0]*primVtx[0] + primVtx[1]*primVtx[1])>1. || TMath::Abs(primVtx[2]>10.)){
292 // Post output data
293 PostData(1, fHistListCosmics);
294 return;
295 }
296 if(!fInputEvent->GetNumberOfTracks() || fInputEvent->GetNumberOfTracks()<2){
297 // Post output data
298 PostData(1, fHistListCosmics);
299 return;
300 }
301 Int_t nTracks = fInputEvent->GetNumberOfTracks();
302
303 if(!fTrackCuts) {
304 // Post output data
305 PostData(1, fHistListCosmics);
306 return;
307 }
308
9b58297c 309 fNEventSel->Fill(0.);
1f329128 310
311 Float_t dcaR[2] = {0.,0.};
312 Float_t dcaZ[2] = {0.,0.};
313
1f329128 314 for (Int_t iTrack1 = 0; iTrack1 < nTracks; iTrack1++) {
315
316 AliESDtrack* track1 = (AliESDtrack*)fInputEvent->GetTrack(iTrack1);
317 if (!track1) continue;
318 if(!(fTrackCuts->AcceptTrack(track1))) { continue; }
319 if(track1->Pt()<fPtMin) continue;
320 //Start 2nd track loop to look for correlations
321 for (Int_t iTrack2 = iTrack1+1; iTrack2 < nTracks; iTrack2++) {
322 AliESDtrack *track2 = (AliESDtrack*)fInputEvent->GetTrack(iTrack2);
323 if(!track2) continue;
324 if(!(fTrackCuts->AcceptTrack(track2))) { continue; }
325
326 //Check if back-to-back
327 Double_t mom1[3],mom2[3];
328 track1->GetPxPyPz(mom1);
329 track2->GetPxPyPz(mom2);
330 // 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]) );
331 TVector3 momv1(mom1[0],mom1[1],mom1[2]);
332 TVector3 momv2(mom2[0],mom2[1],mom2[2]);
333 //Double_t theta = momv1.Angle(momv2);
334 Double_t theta = momv1.Phi()-momv2.Phi();
335 if(theta<-0.5*TMath::Pi()) theta+=2.*TMath::Pi();
336
337 fDeltaPtCosmicCandidates->Fill(track1->Pt()-track2->Pt());
338 Float_t deltaPhi = track1->Phi()-track2->Phi();
339 if(deltaPhi<-0.5*TMath::Pi()) deltaPhi+=2.*TMath::Pi();
340 fDeltaPhiSumEta->Fill(deltaPhi,track1->Eta()+track2->Eta());
341
342 track1->GetImpactParameters(dcaR[0],dcaZ[0]);
343 track2->GetImpactParameters(dcaR[1],dcaZ[1]);
344
345 if(track2->Pt()<0.5) continue;
346 Double_t rIsol = TMath::Sqrt( deltaPhi*deltaPhi+(track1->Eta()-track2->Eta())*(track1->Eta()-track2->Eta()) );
347 fRisol->Fill(rIsol); //Fill R histogram
348 if(track2->Pt()<fPtMin) continue;
349
350 fTheta->Fill(theta);
351 fThetaZoom->Fill(theta);
352 fThetaPt1Pt2->Fill(theta,track1->Pt(),track2->Pt());
353 fThetaPt1Pt2Signed->Fill(theta,track1->GetSign()*track1->Pt(),track2->GetSign()*track2->Pt());
354 fDeltaPhiSumEtaPt1->Fill(deltaPhi,track1->Eta()+track2->Eta(),track1->Pt());
355 fDeltaPhiSumEtaPt2->Fill(deltaPhi,track1->Eta()+track2->Eta(),track2->Pt());
356 fThetaDCAZ1DCAZ2->Fill(theta,dcaZ[0],dcaZ[1]);
357 fRisolTheta->Fill(rIsol,theta);
358 if(TMath::Abs(TMath::Pi()-theta)<fMaxCosmicAngle) {
359 fDCAZCosmicCandidates->Fill(dcaZ[0],dcaZ[1]);
360 fDCARCosmicCandidates->Fill(dcaR[0],dcaR[1]);
361 fPtSignedCosmicCandidates->Fill(track1->GetSign()*track1->Pt());
362 }
363
364 } // track2 loop
365
366 } // track1 loop
367
368 PostData(1,fHistListCosmics);
369}
370
371//________________________________________________________________________
372void AliPWG4CosmicCandidates::Terminate(Option_t *)
373{
9b58297c 374 //
1f329128 375 // Called once at the end of the query
9b58297c 376 //
1f329128 377
378}