1 /**************************************************************************
2 * Copyright(c) 1998-2007, 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 **************************************************************************/
18 // Analysis Task for the Quality Assurance of Cosmic Data
19 // Two track segments in the are matched in angle and charged. T
20 // The quality of the matching in is checked by comparing
21 // the tarnsverse momenta and starting points of the track segments
25 // andreas.morsch@cern.ch
35 #include "AliESDEvent.h"
36 #include "AliESDtrack.h"
37 #include "AliAnalysisTaskCosmic.h"
39 ClassImp(AliAnalysisTaskCosmic)
41 //________________________________________________________________________
42 AliAnalysisTaskCosmic::AliAnalysisTaskCosmic(const char *name)
43 : AliAnalysisTaskSE(name),
56 for (Int_t i = 0; i < 6; i++) {
72 DefineOutput(1, TList::Class());
75 AliAnalysisTaskCosmic::AliAnalysisTaskCosmic(const AliAnalysisTaskCosmic& task)
76 : AliAnalysisTaskSE(task),
87 fhDZvsZ = (TH2F*) ((task.fhDZvsZ) ->Clone());
88 fhDZvsPhi = (TH2F*) ((task.fhDZvsPhi) ->Clone());
89 fhCh1Ch2 = (TH2F*) ((task.fhCh1Ch2) ->Clone());
90 fhPh1Ph2 = (TH2F*) ((task.fhPh1Ph2) ->Clone());
91 fhCl1Cl2G = (TH2F*) ((task.fhCl1Cl2G) ->Clone());
92 fhCl1Cl2B = (TH2F*) ((task.fhCl1Cl2B) ->Clone());
93 for (Int_t i = 0; i < 6; i++) {
94 fhPt[i] = (TH1F*) ((task.fhPt[i]) ->Clone());
95 fhTheta[i] = (TH1F*) ((task.fhTheta[i]) ->Clone());
96 fhPhi[i] = (TH1F*) ((task.fhPhi[i]) ->Clone());
97 fhDPhi[i] = (TH1F*) ((task.fhDPhi[i]) ->Clone());
98 fhDTheta[i] = (TH1F*) ((task.fhDTheta[i]) ->Clone());
99 fhDZ[i] = (TH1F*) ((task.fhDZ[i]) ->Clone());
100 fhDX[i] = (TH1F*) ((task.fhDX[i]) ->Clone());
101 fhDY[i] = (TH1F*) ((task.fhDY[i]) ->Clone());
102 fhDPt[i] = (TH1F*) ((task.fhDPt[i]) ->Clone());
103 fhD1ovPt[i] = (TH1F*) ((task.fhD1ovPt[i]) ->Clone());
104 fhDPtovPt[i] = (TH1F*) ((task.fhDPtovPt[i]) ->Clone());
106 fpDPt[i] = (TProfile*) ((task.fpDPt[i]) ->Clone());
107 fpDPtS[i] = (TProfile*) ((task.fpDPtS[i]) ->Clone());
111 //______________________________________________________________________________
112 AliAnalysisTaskCosmic& AliAnalysisTaskCosmic::operator=(const AliAnalysisTaskCosmic& task)
114 // Assignment operator
117 AliAnalysisTaskSE::operator=(task);
119 fHists = new TList();
120 fhDZvsZ = (TH2F*) ((task.fhDZvsZ) ->Clone());
121 fhDZvsPhi = (TH2F*) ((task.fhDZvsPhi) ->Clone());
122 fhCh1Ch2 = (TH2F*) ((task.fhCh1Ch2) ->Clone());
123 fhPh1Ph2 = (TH2F*) ((task.fhPh1Ph2) ->Clone());
124 fhCl1Cl2G = (TH2F*) ((task.fhCl1Cl2G) ->Clone());
125 fhCl1Cl2B = (TH2F*) ((task.fhCl1Cl2B) ->Clone());
126 for (Int_t i = 0; i < 6; i++) {
127 fhPt[i] = (TH1F*) ((task.fhPt[i]) ->Clone());
128 fhTheta[i] = (TH1F*) ((task.fhTheta[i]) ->Clone());
129 fhPhi[i] = (TH1F*) ((task.fhPhi[i]) ->Clone());
130 fhDPhi[i] = (TH1F*) ((task.fhDPhi[i]) ->Clone());
131 fhDTheta[i] = (TH1F*) ((task.fhDTheta[i]) ->Clone());
132 fhDZ[i] = (TH1F*) ((task.fhDZ[i]) ->Clone());
133 fhDX[i] = (TH1F*) ((task.fhDX[i]) ->Clone());
134 fhDY[i] = (TH1F*) ((task.fhDY[i]) ->Clone());
135 fhDPt[i] = (TH1F*) ((task.fhDPt[i]) ->Clone());
136 fhD1ovPt[i] = (TH1F*) ((task.fhD1ovPt[i]) ->Clone());
137 fhDPtovPt[i] = (TH1F*) ((task.fhDPtovPt[i]) ->Clone());
139 fpDPt[i] = (TProfile*) ((task.fpDPt[i]) ->Clone());
140 fpDPtS[i] = (TProfile*) ((task.fpDPtS[i]) ->Clone());
147 //________________________________________________________________________
148 void AliAnalysisTaskCosmic::UserCreateOutputObjects()
152 TString ext[6] = {"PC", "NC", "PZ", "NZ", "_Good", "_Bad"};
156 fHists = new TList();
158 for (Int_t i = 0; i < 6; i++) {
160 sprintf(name, "fhPt%2s", ext[i].Data());
161 fhPt[i] = new TH1F(name, " pT distribution", 800, 0., 200.);
162 fhPt[i]->SetXTitle("p_{T} [Gev]");
164 sprintf(name, "fhPhi%2s", ext[i].Data());
165 fhPhi[i] = new TH1F(name, "Phi distribution", 62, 0., 2. * TMath::Pi());
166 fhPhi[i]->SetXTitle("#phi [rad]");
168 sprintf(name, "fhTheta%2s", ext[i].Data());
169 fhTheta[i] = new TH1F(name, "Theta distribution", 62, 0., TMath::Pi());
170 fhTheta[i]->SetXTitle("#theta [rad]");
172 sprintf(name, "fhDPhi%2s", ext[i].Data());
173 fhDPhi[i] = new TH1F(name, "DeltaPhi distribution", 320, -0.4, 0.4);
174 fhDPhi[i]->SetXTitle("#Delta#phi [rad]");
176 sprintf(name, "fhDTheta%2s", ext[i].Data());
177 fhDTheta[i] = new TH1F(name, "DeltaTheta distribution", 320, -0.4, 0.4);
178 fhDTheta[i]->SetXTitle("#Delta#theta [rad]");
180 sprintf(name, "fhDZ%2s", ext[i].Data());
181 fhDZ[i] = new TH1F(name, "DeltaZ distribution", 200, -10., 10.);
182 fhDZ[i]->SetXTitle("#DeltaZ [cm]");
184 sprintf(name, "fhDX%2s", ext[i].Data());
185 fhDX[i] = new TH1F(name, "DeltaX distribution", 200, -10., 10.);
186 fhDX[i]->SetXTitle("#DeltaX [cm]");
188 sprintf(name, "fhDY%2s", ext[i].Data());
189 fhDY[i] = new TH1F(name, "DeltaY distribution", 200, -10, 10.);
190 fhDY[i]->SetXTitle("#DeltaY [cm]");
192 sprintf(name, "fhDPt%2s", ext[i].Data());
193 fhDPt[i] = new TH1F(name, "DeltaPt distribution", 200, -20., 20.);
194 fhDPt[i]->SetXTitle("#Delta p_{T} [GeV]");
197 sprintf(name, "fhD1ovPt%2s", ext[i].Data());
198 fhD1ovPt[i] = new TH1F(name, "Delta 1/Pt distribution", 200, -1., 1.);
199 fhD1ovPt[i]->SetXTitle("#Delta 1/Pt");
202 sprintf(name, "fhDPtovPt%2s", ext[i].Data());
203 fhDPtovPt[i] = new TH1F(name, "DeltaPt/Pt distribution", 200, -2., 2.);
204 fhDPtovPt[i]->SetXTitle("#DeltaPt/Pt");
208 // Delta Pt/ Pt vs Pt
209 sprintf(name, "fpDPt%2s", ext[i].Data());
210 fpDPt[i] = new TProfile(name, "#Delta Pt / Pt", 20, 0., 20., -1, 1., "S");
211 fpDPt[i]->SetXTitle("p_{T} [GeV]");
212 fpDPt[i]->SetYTitle("#Delta 1/p_{T} [GeV^{-1}]");
214 sprintf(name, "fpDPtS%2s", ext[i].Data());
215 fpDPtS[i] = new TProfile(name, "#Delta Pt / <sigma>", 20, 0., 20., 0., 10.);
216 fpDPtS[i]->SetXTitle("p_{T}");
217 fpDPtS[i]->SetYTitle("#Delta p_{T} / <#sigma_{p_{T}}>");
220 fHists->Add(fhPt[i]);
221 fHists->Add(fhPhi[i]);
222 fHists->Add(fhTheta[i]);
223 fHists->Add(fhDPhi[i]);
224 fHists->Add(fhDPt[i]);
225 fHists->Add(fhD1ovPt[i]);
226 fHists->Add(fhDPtovPt[i]);
227 fHists->Add(fhDTheta[i]);
228 fHists->Add(fhDZ[i]);
229 fHists->Add(fhDX[i]);
230 fHists->Add(fhDY[i]);
231 fHists->Add(fpDPt[i]);
232 fHists->Add(fpDPtS[i]);
235 fhDZvsZ = new TH2F("fhDZvsZ", "dz vs z", 500, -250., 250., 100, -10., 10.);
236 fhDZvsZ->SetXTitle("z_{in} * sign(z_{in}) * sign(z_{out}) [cm]");
237 fhDZvsZ->SetYTitle("#Delta z [cm]");
240 fhDZvsPhi = new TH2F("fhDZvsPhi", "dz vs phi", 36, 0., TMath::Pi(), 50, -2., 2.);
242 fhCh1Ch2 = new TH2F("fCh1Ch2", "ch1 vs ch2", 8, -2., 2., 8, -2., 2.);
243 fhPh1Ph2 = new TH2F("fPh1Ph2", "ph1 vs ph2", 128, 0., 2. * TMath::Pi(), 128, 0., 2. * TMath::Pi());
244 fhCl1Cl2G = new TH2F("fCl1Cl2G", "#cl vs #cl", 200, 0., 200., 200, 0., 200);
245 fhCl1Cl2B = new TH2F("fCl1Cl2B", "#cl vs #cl", 200, 0., 200., 200, 0., 200);
247 fHists->Add(fhDZvsZ);
248 fHists->Add(fhDZvsPhi);
249 fHists->Add(fhCh1Ch2);
250 fHists->Add(fhPh1Ph2);
251 fHists->Add(fhCl1Cl2G);
252 fHists->Add(fhCl1Cl2B);
257 //________________________________________________________________________
258 void AliAnalysisTaskCosmic::UserExec(Option_t *)
261 // Called for each event
264 Printf("ERROR: fESD not available");
268 AliESDEvent* esdE = (AliESDEvent*)fInputEvent;
270 if (esdE->GetNumberOfTracks() != 2) {
276 for (Int_t iTracks = 0; iTracks < esdE->GetNumberOfTracks(); iTracks++) {
277 AliESDtrack* track = esdE->GetTrack(iTracks);
281 // const AliExternalTrackParam * track = trackG->GetTPCInnerParam();
282 // if (!track) continue;
284 Float_t pt = track->Pt();
285 Short_t charge = track->Charge();
286 Float_t phi = track->Phi();
287 if (phi > 0. && phi < TMath::Pi()) charge*=(-1);
292 UInt_t status = track->GetStatus();
293 if ((status&AliESDtrack::kTPCrefit) ==0) continue;
295 Int_t nClustersTPC = track->GetTPCclusters(0);
296 if (nClustersTPC < 50) continue;
300 Float_t z = track->GetZ();
301 Float_t x = track->Xv();
302 Float_t y = track->Yv();
303 Float_t theta = track->Theta();
306 const AliExternalTrackParam * trackOut = track->GetOuterParam();
308 if (trackOut)zOut = trackOut->Zv();
312 fhPt[kPosC] ->Fill(pt);
314 fhPt[kNegC] ->Fill(pt);
317 // if ((TMath::Abs(esdE->GetCurrentL3()) > 1.e-3) && (pt < 1. || pt > 50.)) continue;
322 fhPhi[kPosC] ->Fill(phi);
323 fhTheta[kPosC]->Fill(theta);
325 fhPhi[kNegC] ->Fill(phi);
326 fhTheta[kNegC]->Fill(theta);
331 fhPt[kPosZ] ->Fill(pt);
332 fhPhi[kPosZ] ->Fill(phi);
333 fhTheta[kPosZ]->Fill(theta);
335 fhPt[kNegZ] ->Fill(pt);
336 fhPhi[kNegZ] ->Fill(phi);
337 fhTheta[kNegZ]->Fill(theta);
342 // Tracks coming from above
343 if (phi > 0. && phi < TMath::Pi()) {
345 // printf("Track#1 %5d %5d %13.3f %13.3f \n", Int_t(Entry()), iTracks, phi, zOut);
347 Float_t dphiMin = 999.;
350 // Search for a matching track (in dphi-dtheta)
351 for (Int_t jTracks = 0; jTracks < esdE->GetNumberOfTracks(); jTracks++) {
352 if (jTracks == iTracks) continue;
354 AliESDtrack* track2 = esdE->GetTrack(jTracks);
356 UInt_t status2 = track2->GetStatus();
357 if ((status2&AliESDtrack::kTPCrefit) ==0) continue;
359 Int_t nClustersTPC2 = track2->GetTPCclusters(0);
360 if (nClustersTPC2 < 50) continue;
361 // if ((TMath::Abs(esdE->GetCurrentL3()) > 1.e-3) && (track2->Pt() < 1. || track2->Pt() > 50.)) continue;
362 // const AliExternalTrackParam * track2 = trackG2->GetTPCInnerParam();
363 // if (!track2) continue;
365 Float_t phi2 = track2->Phi() - TMath::Pi();
366 Float_t theta2 = TMath::Pi() - track2->Theta();
367 Float_t dphi = phi2 - phi;
368 Float_t dtheta = theta2 - theta;
370 if (dphi > TMath::Pi()) dphi -= 2. * TMath::Pi();
371 if (dphi < -TMath::Pi()) dphi += 2. * TMath::Pi();
373 Float_t dR = TMath::Sqrt(dphi * dphi + dtheta * dtheta);
383 // we found a matching track candidate ...
384 AliESDtrack* track2 = esdE->GetTrack(jMin);
385 const AliExternalTrackParam * trackOut2 = track2->GetOuterParam();
387 if (trackOut2) zOut2 = trackOut2->Zv();
390 Float_t theta2 = - track2->Theta() + TMath::Pi();
391 Float_t z2 = track2->GetZ();
392 Float_t x2 = track2->Xv();
393 Float_t y2 = track2->Yv();
394 Short_t charge2 = track2->Charge();
398 Float_t pt1 = track->Pt();
399 Float_t pt2 = track2->Pt();
400 Float_t dpt = pt2 - pt1;
401 Float_t d1pt = 1./pt2 - 1./pt1;
403 Float_t ptm = 0.5 * (pt1 + pt2);
405 Int_t nClustersTPC2 = track2->GetTPCclusters(0);
408 fhDPhi[kPosC] ->Fill(dphiMin);
409 fhDTheta[kPosC]->Fill(theta2 - theta);
411 fhDPhi[kNegC] ->Fill(dphiMin);
412 fhDTheta[kNegC]->Fill(theta2 - theta);
416 fhDPhi[kPosZ] ->Fill(dphiMin);
417 fhDTheta[kPosZ]->Fill(theta2 - theta);
419 fhDPhi[kNegZ] ->Fill(dphiMin);
420 fhDTheta[kNegZ]->Fill(theta2 - theta);
423 if (TMath::Abs(dpt)/ptm > 0.5) {
424 fhDPhi[kBad] ->Fill(dphiMin);
425 fhDTheta[kBad]->Fill(theta2 - theta);
427 fhDPhi[kGood] ->Fill(dphiMin);
428 fhDTheta[kGood]->Fill(theta2 - theta);
431 // if (TMath::Abs(rMin < 0.04) && (charge == charge2) && TMath::Abs(dz) < 0.5 && TMath::Abs(dy) && TMath::Abs(dx) < 0.5 )
432 // if (TMath::Abs(rMin < 0.04) && (charge == charge2) && (zOut * zOut2) < 0.)
433 if (TMath::Abs(rMin < 0.04) && (charge == charge2))
436 if (TMath::Abs(dpt)/ptm < .1) {
437 fhCl1Cl2G->Fill(nClustersTPC, nClustersTPC2);
440 if (TMath::Abs(dpt)/ptm > .5) fhCl1Cl2B->Fill(nClustersTPC, nClustersTPC2);
441 fhPh1Ph2->Fill(track->Phi(), track2->Phi());
444 Double_t sigmaPt1 = TMath::Sqrt(track->GetSigma1Pt2());
445 Double_t sigmaPt2 = TMath::Sqrt(track->GetSigma1Pt2());
446 Double_t sigmaPt = 0.5 * TMath::Sqrt(sigmaPt1 * sigmaPt1 + sigmaPt2 * sigmaPt2);
447 if (TMath::Abs(dpt)/ptm > 0.2 && pt > 10. && charge > 0.)
448 // printf("Track#2 %5d %5d %13.3f %13.3f %13.3f %13.3f\n", Entry(), jMin, track2->Phi(), dz, dpt, zOut2);
449 if (TMath::Abs(zOut) > 1.e-10 && TMath::Abs(zOut2) > 1.e-10) fhDZvsZ->Fill(TMath::Abs(zOut) * zOut/zOut2, dz);
450 if (zOut * zOut2 > 0. && zOut < 0) fhDZvsPhi->Fill(phi, dz);
452 fhCh1Ch2->Fill(charge, track2->Charge());
456 fhDPt[kPosC]->Fill(dpt);
457 fhD1ovPt[kPosC]->Fill(d1pt);
458 fpDPt[kPosC]->Fill(ptm, d1pt);
459 if (ptm > 5. && ptm < 15.) fhDPtovPt[kPosC]->Fill(dpt/ptm);
460 fhDZ[kPosC]->Fill(dz);
461 fhDX[kPosC]->Fill(dx);
462 fhDY[kPosC]->Fill(dy);
463 fpDPtS[kPosC]->Fill(ptm, 2. * TMath::Abs(1./pt1 - 1./pt2)/sigmaPt);
465 fhDPt[kNegC]->Fill(dpt);
466 fhD1ovPt[kNegC]->Fill(d1pt);
467 fpDPt[kNegC]->Fill(ptm, d1pt);
468 if (ptm > 5. && ptm < 15.) fhDPtovPt[kNegC]->Fill(dpt/ptm);
469 fhDZ[kNegC]->Fill(dz);
470 fhDX[kNegC]->Fill(dx);
471 fhDY[kNegC]->Fill(dy);
472 fpDPtS[kNegC]->Fill(ptm, 2. * TMath::Abs(1./pt1 - 1./pt2)/sigmaPt);
476 fhDPt[kPosZ]->Fill(dpt);
477 fhD1ovPt[kPosZ]->Fill(d1pt);
478 fpDPt[kPosZ]->Fill(ptm, d1pt);
479 fhDZ[kPosZ]->Fill(dz);
480 fhDX[kPosZ]->Fill(dx);
481 fhDY[kPosZ]->Fill(dy);
482 fpDPtS[kPosZ]->Fill(ptm, 2. * TMath::Abs(1./pt1 - 1./pt2)/sigmaPt);
484 fhDPt[kNegZ]->Fill(dpt);
485 fhD1ovPt[kNegZ]->Fill(d1pt);
486 fpDPt[kNegZ]->Fill(ptm, d1pt);
487 fhDZ[kNegZ]->Fill(dz);
488 fhDX[kNegZ]->Fill(dx);
489 fhDY[kNegZ]->Fill(dy);
490 fpDPtS[kNegZ]->Fill(ptm, 2. * TMath::Abs(1./pt1 - 1./pt2)/sigmaPt);
494 if (TMath::Abs(dpt)/ptm > 0.5) {
495 fhDPt[kBad]->Fill(dpt);
496 fpDPt[kBad]->Fill(ptm, TMath::Abs(dpt)/ptm);
497 fhDZ[kBad]->Fill(dz);
498 fhDX[kBad]->Fill(dx);
499 fhDY[kBad]->Fill(dy);
500 fpDPtS[kBad]->Fill(ptm, 2. * TMath::Abs(1./pt1 - 1./pt2)/sigmaPt);
502 fhDPt[kGood]->Fill(dpt);
503 fpDPt[kGood]->Fill(ptm, TMath::Abs(dpt)/ptm);
504 fhDZ[kGood]->Fill(dz);
505 fhDX[kGood]->Fill(dx);
506 fhDY[kGood]->Fill(dy);
507 fpDPtS[kGood]->Fill(ptm, 2. * TMath::Abs(1./pt1 - 1./pt2)/sigmaPt);
511 } // found possible match
517 //________________________________________________________________________
518 void AliAnalysisTaskCosmic::Terminate(Option_t *)