1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 //-----------------------------------------------------------------------
17 // Analysis task looking for cosmic candidates embedded in a p-p event
18 // Task checks if particles are back-to-back in eta and phi
21 // Author : Marta Verweij - UU - marta.verweij@cern.ch
22 //-----------------------------------------------------------------------
32 #include "AliAnalysisTask.h"
33 #include "AliAnalysisManager.h"
35 #include "AliESDEvent.h"
36 #include "AliESDInputHandler.h"
38 #include "AliESDtrack.h"
39 #include "AliESDtrackCuts.h"
40 #include "AliExternalTrackParam.h"
44 #include "AliPWG4CosmicCandidates.h"
46 //using namespace std; //required for resolving the 'cout' symbol
49 ClassImp(AliPWG4CosmicCandidates)
51 //________________________________________________________________________
52 AliPWG4CosmicCandidates::AliPWG4CosmicCandidates()
53 : AliAnalysisTaskSE(),
56 fMaxCosmicAngle(0.002),
59 fPtSignedCosmicCandidates(0),
60 fDeltaPtCosmicCandidates(0),
62 fDCAZCosmicCandidates(0),
63 fDCARCosmicCandidates(0),
67 fThetaPt1Pt2Signed(0),
68 fDeltaPhiSumEtaPt1(0),
69 fDeltaPhiSumEtaPt2(0),
76 // Default constructor
80 //________________________________________________________________________
81 AliPWG4CosmicCandidates::AliPWG4CosmicCandidates(const char *name)
82 : AliAnalysisTaskSE(name),
85 fMaxCosmicAngle(0.002),
88 fPtSignedCosmicCandidates(0),
89 fDeltaPtCosmicCandidates(0),
91 fDCAZCosmicCandidates(0),
92 fDCARCosmicCandidates(0),
96 fThetaPt1Pt2Signed(0),
97 fDeltaPhiSumEtaPt1(0),
98 fDeltaPhiSumEtaPt2(0),
104 // Constructor. Initialization of Inputs and Outputs
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());
115 //________________________________________________________________________
116 AliPWG4CosmicCandidates::AliPWG4CosmicCandidates(const AliPWG4CosmicCandidates &res)
117 : AliAnalysisTaskSE(res),
120 fMaxCosmicAngle(0.002),
123 fPtSignedCosmicCandidates(0),
124 fDeltaPtCosmicCandidates(0),
126 fDCAZCosmicCandidates(0),
127 fDCARCosmicCandidates(0),
131 fThetaPt1Pt2Signed(0),
132 fDeltaPhiSumEtaPt1(0),
133 fDeltaPhiSumEtaPt2(0),
140 // Dummy copy constructor
144 //________________________________________________________________________
145 AliPWG4CosmicCandidates& AliPWG4CosmicCandidates::operator=(const AliPWG4CosmicCandidates& /*trclass*/)
148 // Dummy assignment operator
153 //________________________________________________________________________
154 void AliPWG4CosmicCandidates::LocalInit()
157 // Only called once at beginning
159 PostData(2,fTrackCuts);
162 //________________________________________________________________________
163 void AliPWG4CosmicCandidates::UserCreateOutputObjects()
166 // Create output objects
168 AliDebug(2,Form(">> AliPWG4CosmicCandidates::UserCreateOutputObjects \n"));
170 Bool_t oldStatus = TH1::AddDirectoryStatus();
171 TH1::AddDirectory(kFALSE);
174 fHistListCosmics = new TList();
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);
182 Float_t fgkPtMax=100.;
183 Int_t fgkNPtBins= (int)(fgkPtMax-fgkPtMin);
184 Int_t fgkNPtBins2D= (int)((fgkPtMax-fgkPtMin)/4.);
186 Int_t fgkNPhiBins=18;
187 Float_t kMinPhi = -0.5*TMath::Pi();
188 Float_t kMaxPhi = 3./2.*TMath::Pi();
190 Int_t fgkNThetaBins=fgkNPhiBins*8;
191 Float_t kMinTheta = -0.5*TMath::Pi();
192 Float_t kMaxTheta = 3./2.*TMath::Pi();
194 Int_t fgkNDCARBins=40;
195 Float_t fgkDCARMin = -0.2;
196 Float_t fgkDCARMax = 0.2;
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 ;
200 Int_t fgkNDCAZBins=40;
201 Float_t fgkDCAZMin = -2.;
202 Float_t fgkDCAZMax = 2.;
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 ;
206 fPtSignedCosmicCandidates = new TH1F("fPtSignedCosmicCandidates","fPtSignedCosmicCandidates",2*fgkNPtBins, -1.*fgkPtMax, fgkPtMax);
207 fHistListCosmics->Add(fPtSignedCosmicCandidates);
209 fDeltaPtCosmicCandidates = new TH1F("fDeltaPtCosmicCandidates","fDeltaPtCosmicCandidates",fgkNPtBins, -50., 50.);
210 fHistListCosmics->Add(fDeltaPtCosmicCandidates);
212 fDeltaPhiSumEta = new TH2F("fDeltaPhiSumEta","fDeltaPhiSumEta",fgkNPhiBins*4,kMinPhi,kMaxPhi,80, -2.,2.);
213 fHistListCosmics->Add(fDeltaPhiSumEta);
215 fDCAZCosmicCandidates = new TH2F("fDCAZCosmicCandidates","fDCAZCosmicCandidates",fgkNDCAZBins,binsDCAZ,fgkNDCAZBins,binsDCAZ);
216 fHistListCosmics->Add(fDCAZCosmicCandidates);
218 fDCARCosmicCandidates = new TH2F("fDCARCosmicCandidates","fDCARCosmicCandidates",fgkNDCARBins,binsDCAR,fgkNDCARBins,binsDCAR);
219 fHistListCosmics->Add(fDCARCosmicCandidates);
221 fTheta = new TH1F("fTheta","fTheta",fgkNThetaBins,kMinTheta,kMaxTheta);
222 fHistListCosmics->Add(fTheta);
224 fThetaZoom = new TH1F("fThetaZoom","fThetaZoom",100,TMath::Pi()-1.,TMath::Pi()+1.);
225 fHistListCosmics->Add(fThetaZoom);
227 fThetaPt1Pt2 = new TH3F("fThetaPt1Pt2","fThetaPt1Pt2",fgkNThetaBins,kMinTheta,kMaxTheta,fgkNPtBins2D,fgkPtMin,fgkPtMax,fgkNPtBins2D,fgkPtMin,fgkPtMax);
228 fHistListCosmics->Add(fThetaPt1Pt2);
230 fThetaPt1Pt2Signed = new TH3F("fThetaPt1Pt2Signed","fThetaPt1Pt2Signed",fgkNThetaBins,kMinTheta,kMaxTheta,fgkNPtBins2D,-1.*fgkPtMax,fgkPtMax,fgkNPtBins2D,-1.*fgkPtMax,fgkPtMax);
231 fHistListCosmics->Add(fThetaPt1Pt2Signed);
233 fDeltaPhiSumEtaPt1 = new TH3F("fDeltaPhiSumEtaPt1","fDeltaPhiSumEtaPt1",fgkNThetaBins,kMinTheta,kMaxTheta,80, -2.,2.,fgkNPtBins2D,fgkPtMin,fgkPtMax);
234 fHistListCosmics->Add(fDeltaPhiSumEtaPt1);
236 fDeltaPhiSumEtaPt2 = new TH3F("fDeltaPhiSumEtaPt2","fDeltaPhiSumEtaPt2",fgkNThetaBins,kMinTheta,kMaxTheta,80, -2.,2.,fgkNPtBins2D,fgkPtMin,fgkPtMax);
237 fHistListCosmics->Add(fDeltaPhiSumEtaPt2);
239 fThetaDCAZ1DCAZ2 = new TH3F("fThetaDCAZ1DCAZ2","fThetaDCAZ1DCAZ2",fgkNThetaBins,kMinTheta,kMaxTheta,fgkNDCAZBins,-2.,2.,fgkNDCAZBins,-2.,2.);
240 fHistListCosmics->Add(fThetaDCAZ1DCAZ2);
242 fRisol = new TH1F("fRisol","fRisol",100,0.,10.);
243 fHistListCosmics->Add(fRisol);
245 fRisolTheta = new TH2F("fRisolTheta","fRisolTheta",100,0.,10.,fgkNThetaBins,kMinTheta,kMaxTheta);
246 fHistListCosmics->Add(fRisolTheta);
248 TH1::AddDirectory(oldStatus);
250 PostData(1, fHistListCosmics);
252 if(binsDCAR) delete [] binsDCAR;
253 if(binsDCAZ) delete [] binsDCAZ;
257 //________________________________________________________________________
258 void AliPWG4CosmicCandidates::UserExec(Option_t *)
261 // Called for each event
264 // All events without selection
265 fNEventAll->Fill(0.);
268 AliDebug(2,Form("ERROR: fESD not available"));
269 cout << "ERROR: fESD not available" << endl;
270 PostData(1, fHistListCosmics);
274 const AliVVertex *vtx = fInputEvent->GetPrimaryVertex();
276 PostData(1, fHistListCosmics);
281 TString vtxName(vtx->GetName());
282 if( vtx->GetNContributors() < 2 || (vtxName.Contains("TPCVertex")) ) {
284 PostData(1, fHistListCosmics);
288 // AliDebug(2,Form("Vertex title %s, status %d, nCont %d\n",vtx->GetTitle(), vtx->GetStatus(), vtx->GetNContributors()));
290 vtx->GetXYZ(primVtx);
291 if(TMath::Sqrt(primVtx[0]*primVtx[0] + primVtx[1]*primVtx[1])>1. || TMath::Abs(primVtx[2]>10.)){
293 PostData(1, fHistListCosmics);
296 if(!fInputEvent->GetNumberOfTracks() || fInputEvent->GetNumberOfTracks()<2){
298 PostData(1, fHistListCosmics);
301 Int_t nTracks = fInputEvent->GetNumberOfTracks();
305 PostData(1, fHistListCosmics);
309 fNEventSel->Fill(0.);
311 Float_t dcaR[2] = {0.,0.};
312 Float_t dcaZ[2] = {0.,0.};
314 for (Int_t iTrack1 = 0; iTrack1 < nTracks; iTrack1++) {
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; }
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();
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());
342 track1->GetImpactParameters(dcaR[0],dcaZ[0]);
343 track2->GetImpactParameters(dcaR[1],dcaZ[1]);
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;
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());
368 PostData(1,fHistListCosmics);
371 //________________________________________________________________________
372 void AliPWG4CosmicCandidates::Terminate(Option_t *)
375 // Called once at the end of the query