Adding C0SCO-D-NOPF-FAST.
[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//-----------------------------------------------------------------------
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);
184
185 Int_t fgkNPhiBins=18;
186 Float_t kMinPhi = -0.5*TMath::Pi();
187 Float_t kMaxPhi = 3./2.*TMath::Pi();
188
189 Int_t fgkNThetaBins=fgkNPhiBins*8;
190 Float_t kMinTheta = -0.5*TMath::Pi();
191 Float_t kMaxTheta = 3./2.*TMath::Pi();
192
193 Int_t fgkNDCARBins=80;
194 Float_t fgkDCARMin = -0.2;
195 Float_t fgkDCARMax = 0.2;
60829498 196 Float_t *binsDCAR=new Float_t[fgkNDCARBins+1];
197 for(Int_t i=0; i<=fgkNDCARBins; i++) binsDCAR[i]=(Float_t)fgkDCARMin + (fgkDCARMax-fgkDCARMin)/fgkNDCARBins*(Float_t)i ;
1f329128 198
199 Int_t fgkNDCAZBins=80;
200 Float_t fgkDCAZMin = -2.;
201 Float_t fgkDCAZMax = 2.;
60829498 202 Float_t *binsDCAZ=new Float_t[fgkNDCAZBins+1];
203 for(Int_t i=0; i<=fgkNDCAZBins; i++) binsDCAZ[i]=(Float_t)fgkDCAZMin + (fgkDCAZMax-fgkDCAZMin)/fgkNDCAZBins*(Float_t)i ;
1f329128 204
205 fPtSignedCosmicCandidates = new TH1F("fPtSignedCosmicCandidates","fPtSignedCosmicCandidates",2*(int)(fgkPtMax-fgkPtMin), -1.*fgkPtMax, fgkPtMax);
206 fHistListCosmics->Add(fPtSignedCosmicCandidates);
207
208 fDeltaPtCosmicCandidates = new TH1F("fDeltaPtCosmicCandidates","fDeltaPtCosmicCandidates",fgkNPtBins, -50., 50.);
209 fHistListCosmics->Add(fDeltaPtCosmicCandidates);
210
211 fDeltaPhiSumEta = new TH2F("fDeltaPhiSumEta","fDeltaPhiSumEta",fgkNPhiBins*4,kMinPhi,kMaxPhi,80, -2.,2.);
212 fHistListCosmics->Add(fDeltaPhiSumEta);
213
214 fDCAZCosmicCandidates = new TH2F("fDCAZCosmicCandidates","fDCAZCosmicCandidates",fgkNDCAZBins,binsDCAZ,fgkNDCAZBins,binsDCAZ);
215 fHistListCosmics->Add(fDCAZCosmicCandidates);
216
217 fDCARCosmicCandidates = new TH2F("fDCARCosmicCandidates","fDCARCosmicCandidates",fgkNDCARBins,binsDCAR,fgkNDCARBins,binsDCAR);
218 fHistListCosmics->Add(fDCARCosmicCandidates);
219
220 fTheta = new TH1F("fTheta","fTheta",fgkNThetaBins,kMinTheta,kMaxTheta);
221 fHistListCosmics->Add(fTheta);
222
223 fThetaZoom = new TH1F("fThetaZoom","fThetaZoom",100,TMath::Pi()-1.,TMath::Pi()+1.);
224 fHistListCosmics->Add(fThetaZoom);
225
226 fThetaPt1Pt2 = new TH3F("fThetaPt1Pt2","fThetaPt1Pt2",fgkNThetaBins,kMinTheta,kMaxTheta,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax);
227 fHistListCosmics->Add(fThetaPt1Pt2);
228
229 fThetaPt1Pt2Signed = new TH3F("fThetaPt1Pt2Signed","fThetaPt1Pt2Signed",fgkNThetaBins,kMinTheta,kMaxTheta,4*(int)(fgkPtMax-fgkPtMin),-1.*fgkPtMax,fgkPtMax,4*(int)(fgkPtMax-fgkPtMin),-1.*fgkPtMax,fgkPtMax);
230 fHistListCosmics->Add(fThetaPt1Pt2Signed);
231
232 fDeltaPhiSumEtaPt1 = new TH3F("fDeltaPhiSumEtaPt1","fDeltaPhiSumEtaPt1",fgkNThetaBins,kMinTheta,kMaxTheta,80, -2.,2.,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax);
233 fHistListCosmics->Add(fDeltaPhiSumEtaPt1);
234
235 fDeltaPhiSumEtaPt2 = new TH3F("fDeltaPhiSumEtaPt2","fDeltaPhiSumEtaPt2",fgkNThetaBins,kMinTheta,kMaxTheta,80, -2.,2.,(int)(fgkPtMax-fgkPtMin),fgkPtMin,fgkPtMax);
236 fHistListCosmics->Add(fDeltaPhiSumEtaPt2);
237
238 fThetaDCAZ1DCAZ2 = new TH3F("fThetaDCAZ1DCAZ2","fThetaDCAZ1DCAZ2",fgkNThetaBins,kMinTheta,kMaxTheta,fgkNDCAZBins,-2.,2.,fgkNDCAZBins,-2.,2.);
239 fHistListCosmics->Add(fThetaDCAZ1DCAZ2);
240
241 fRisol = new TH1F("fRisol","fRisol",100,0.,10.);
242 fHistListCosmics->Add(fRisol);
243
244 fRisolTheta = new TH2F("fRisolTheta","fRisolTheta",100,0.,10.,fgkNThetaBins,kMinTheta,kMaxTheta);
245 fHistListCosmics->Add(fRisolTheta);
246
247 TH1::AddDirectory(oldStatus);
248
249 if(binsDCAR) delete [] binsDCAR;
250 if(binsDCAZ) delete [] binsDCAZ;
251
252}
253
254//________________________________________________________________________
255void AliPWG4CosmicCandidates::UserExec(Option_t *)
256{
9b58297c 257 // Main loop
258 // Called for each event
259 //
1f329128 260
261 // All events without selection
262 fNEventAll->Fill(0.);
263
264 if (!fInputEvent) {
265 AliDebug(2,Form("ERROR: fESD not available"));
266 cout << "ERROR: fESD not available" << endl;
267 PostData(1, fHistListCosmics);
268 return;
269 }
270
271 const AliVVertex *vtx = fInputEvent->GetPrimaryVertex();
272 // Need vertex cut
273 if (!vtx || vtx->GetNContributors() < 2) {
274 // Post output data
275 PostData(1, fHistListCosmics);
276 return;
277 }
278
279 // AliDebug(2,Form("Vertex title %s, status %d, nCont %d\n",vtx->GetTitle(), vtx->GetStatus(), vtx->GetNContributors()));
280 double primVtx[3];
281 vtx->GetXYZ(primVtx);
282 if(TMath::Sqrt(primVtx[0]*primVtx[0] + primVtx[1]*primVtx[1])>1. || TMath::Abs(primVtx[2]>10.)){
283 // Post output data
284 PostData(1, fHistListCosmics);
285 return;
286 }
287 if(!fInputEvent->GetNumberOfTracks() || fInputEvent->GetNumberOfTracks()<2){
288 // Post output data
289 PostData(1, fHistListCosmics);
290 return;
291 }
292 Int_t nTracks = fInputEvent->GetNumberOfTracks();
293
294 if(!fTrackCuts) {
295 // Post output data
296 PostData(1, fHistListCosmics);
297 return;
298 }
299
9b58297c 300 fNEventSel->Fill(0.);
1f329128 301
302 Float_t dcaR[2] = {0.,0.};
303 Float_t dcaZ[2] = {0.,0.};
304
1f329128 305 for (Int_t iTrack1 = 0; iTrack1 < nTracks; iTrack1++) {
306
307 AliESDtrack* track1 = (AliESDtrack*)fInputEvent->GetTrack(iTrack1);
308 if (!track1) continue;
309 if(!(fTrackCuts->AcceptTrack(track1))) { continue; }
310 if(track1->Pt()<fPtMin) continue;
311 //Start 2nd track loop to look for correlations
312 for (Int_t iTrack2 = iTrack1+1; iTrack2 < nTracks; iTrack2++) {
313 AliESDtrack *track2 = (AliESDtrack*)fInputEvent->GetTrack(iTrack2);
314 if(!track2) continue;
315 if(!(fTrackCuts->AcceptTrack(track2))) { continue; }
316
317 //Check if back-to-back
318 Double_t mom1[3],mom2[3];
319 track1->GetPxPyPz(mom1);
320 track2->GetPxPyPz(mom2);
321 // 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]) );
322 TVector3 momv1(mom1[0],mom1[1],mom1[2]);
323 TVector3 momv2(mom2[0],mom2[1],mom2[2]);
324 //Double_t theta = momv1.Angle(momv2);
325 Double_t theta = momv1.Phi()-momv2.Phi();
326 if(theta<-0.5*TMath::Pi()) theta+=2.*TMath::Pi();
327
328 fDeltaPtCosmicCandidates->Fill(track1->Pt()-track2->Pt());
329 Float_t deltaPhi = track1->Phi()-track2->Phi();
330 if(deltaPhi<-0.5*TMath::Pi()) deltaPhi+=2.*TMath::Pi();
331 fDeltaPhiSumEta->Fill(deltaPhi,track1->Eta()+track2->Eta());
332
333 track1->GetImpactParameters(dcaR[0],dcaZ[0]);
334 track2->GetImpactParameters(dcaR[1],dcaZ[1]);
335
336 if(track2->Pt()<0.5) continue;
337 Double_t rIsol = TMath::Sqrt( deltaPhi*deltaPhi+(track1->Eta()-track2->Eta())*(track1->Eta()-track2->Eta()) );
338 fRisol->Fill(rIsol); //Fill R histogram
339 if(track2->Pt()<fPtMin) continue;
340
341 fTheta->Fill(theta);
342 fThetaZoom->Fill(theta);
343 fThetaPt1Pt2->Fill(theta,track1->Pt(),track2->Pt());
344 fThetaPt1Pt2Signed->Fill(theta,track1->GetSign()*track1->Pt(),track2->GetSign()*track2->Pt());
345 fDeltaPhiSumEtaPt1->Fill(deltaPhi,track1->Eta()+track2->Eta(),track1->Pt());
346 fDeltaPhiSumEtaPt2->Fill(deltaPhi,track1->Eta()+track2->Eta(),track2->Pt());
347 fThetaDCAZ1DCAZ2->Fill(theta,dcaZ[0],dcaZ[1]);
348 fRisolTheta->Fill(rIsol,theta);
349 if(TMath::Abs(TMath::Pi()-theta)<fMaxCosmicAngle) {
350 fDCAZCosmicCandidates->Fill(dcaZ[0],dcaZ[1]);
351 fDCARCosmicCandidates->Fill(dcaR[0],dcaR[1]);
352 fPtSignedCosmicCandidates->Fill(track1->GetSign()*track1->Pt());
353 }
354
355 } // track2 loop
356
357 } // track1 loop
358
359 PostData(1,fHistListCosmics);
360}
361
362//________________________________________________________________________
363void AliPWG4CosmicCandidates::Terminate(Option_t *)
364{
9b58297c 365 //
1f329128 366 // Called once at the end of the query
9b58297c 367 //
1f329128 368
369}