]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TRD/TRDqaAnalysis/AliTRDqaJPsi.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
TRD module
[u/mrichter/AliRoot.git] / TRD / TRDqaAnalysis / AliTRDqaJPsi.cxx
CommitLineData
8e2f611a 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 iLs" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id$ */
17
18//
19// This class is a part of a package of high level QA monitoring for TRD.
20//
21// S. Radomski
22// radomski@physi.uni-heidelberg.de
23// March 2008
24//
25
26#include "AliTRDqaJPsi.h"
27#include "AliTRDqaAT.h"
28
29#include "TFile.h"
30#include "TH1D.h"
31#include "TH2D.h"
32#include "TChain.h"
33
34#include "AliESDEvent.h"
35#include "AliESDtrack.h"
36#include "AliKFParticle.h"
37
38#include "TLorentzVector.h"
39
8e2f611a 40//______________________________________________________________________________
41
42AliTRDqaJPsi::AliTRDqaJPsi()
43 : AliAnalysisTask("",""),
44 fChain(0),
45 fESD(0),
46 fOutputContainer(0),
47 fnKFtracks(0)
48{
49 //
50 // default dummy constructor
51 //
02f3bfcc 52
ce42378b 53 for (Int_t k = 0; k < fgknSteps; k++) {
54 fStatus[k] = 0x0;
55 fAngleSM[k] = 0x0;
56 fInvMass[k] = 0x0;
57 fInvMassVec[k] = 0x0;
58 fInvMassDiff[k] = 0x0;
59 fPtAngle[k] = 0x0;
60 }
61
62 for (Int_t l = 0; l < 2*fgknSteps; l++) {
63 fnTracks[l] = 0x0;
64 fPt[l] = 0x0;
65 fPID[l] = 0x0;
66 }
67
02f3bfcc 68 for (Int_t i = 0; i < 1000; i++) {
ce42378b 69 fSM[i] = 0;
70 fTracks[i] = 0x0;
71 fVec[i] = 0x0;
72 for (Int_t j = 0; j < fgknSteps; j++) {
73 fInSample[i][j] = 0;
74 }
02f3bfcc 75 }
ce42378b 76
8e2f611a 77}
8e2f611a 78
c336eab0 79//______________________________________________________________________________
80AliTRDqaJPsi:: AliTRDqaJPsi(const AliTRDqaJPsi & /*trd*/)
8e2f611a 81 : AliAnalysisTask("",""),
82 fChain(0),
83 fESD(0),
84 fOutputContainer(0),
85 fnKFtracks(0)
86{
87 //
c336eab0 88 // Copy constructor
8e2f611a 89 //
90
ce42378b 91 for (Int_t k = 0; k < fgknSteps; k++) {
92 fStatus[k] = 0x0;
93 fAngleSM[k] = 0x0;
94 fInvMass[k] = 0x0;
95 fInvMassVec[k] = 0x0;
96 fInvMassDiff[k] = 0x0;
97 fPtAngle[k] = 0x0;
98 }
99
100 for (Int_t l = 0; l < 2*fgknSteps; l++) {
101 fnTracks[l] = 0x0;
102 fPt[l] = 0x0;
103 fPID[l] = 0x0;
104 }
105
a987273c 106 for (Int_t i = 0; i < 1000; i++) {
ce42378b 107 fSM[i] = 0;
108 fTracks[i] = 0x0;
109 fVec[i] = 0x0;
110 for (Int_t j = 0; j < fgknSteps; j++) {
111 fInSample[i][j] = 0;
112 }
a987273c 113 }
114
8e2f611a 115}
116
8e2f611a 117//______________________________________________________________________________
118AliTRDqaJPsi::AliTRDqaJPsi(const char *name)
119 : AliAnalysisTask(name,""),
120 fChain(0),
121 fESD(0),
122 fOutputContainer(0),
123 fnKFtracks(0)
124
125{
126 // Constructor.
127 // Input slot #0 works with an Ntuple
128 DefineInput(0, TChain::Class());
129 // Output slot #0 writes into a TH1 container
130 DefineOutput(0, TObjArray::Class()) ;
ce42378b 131
132 for (Int_t k = 0; k < fgknSteps; k++) {
133 fStatus[k] = 0x0;
134 fAngleSM[k] = 0x0;
135 fInvMass[k] = 0x0;
136 fInvMassVec[k] = 0x0;
137 fInvMassDiff[k] = 0x0;
138 fPtAngle[k] = 0x0;
139 }
140
141 for (Int_t l = 0; l < 2*fgknSteps; l++) {
142 fnTracks[l] = 0x0;
143 fPt[l] = 0x0;
144 fPID[l] = 0x0;
145 }
146
02f3bfcc 147 for (Int_t i = 0; i < 1000; i++) {
ce42378b 148 fSM[i] = 0;
149 fTracks[i] = 0x0;
150 fVec[i] = 0x0;
151 for (Int_t j = 0; j < fgknSteps; j++) {
152 fInSample[i][j] = 0;
153 }
02f3bfcc 154 }
155
8e2f611a 156}
157
158//______________________________________________________________________________
159void AliTRDqaJPsi::ConnectInputData(const Option_t *)
160{
161 // Initialisation of branch container and histograms
162
163 //AliInfo(Form("*** Initialization of %s", GetName())) ;
164
165 fChain = (TChain*)GetInputData(0);
166 fESD = new AliESDEvent();
167 fESD->ReadFromTree(fChain);
168}
169
170//________________________________________________________________________
171void AliTRDqaJPsi::CreateOutputObjects()
172{
c336eab0 173 //
174 // Create the objects that contain the analysis output
175 //
8e2f611a 176
177 Int_t c = 0;
178 fOutputContainer = new TObjArray(100);
179
180 const char *charge[2] = {"Neg", "Pos"};
181
182 // build histograms
183
c336eab0 184 for(Int_t i=0; i<fgknSteps; i++) {
8e2f611a 185
186 fStatus[i] = new TH1D(Form("status_%d", i), "status", 32, -0.5, 31.5);
187 fOutputContainer->AddAt(fStatus[i], c++);
188
189 fInvMass[i] = new TH1D(Form("mass_%d", i), ";m_{inv} (GeV);", 100, 0, 5);
190 fOutputContainer->AddAt(fInvMass[i], c++);
191
192 fInvMassVec[i] = new TH1D(Form("massVec_%d", i), ";m_{inv} (GeV);", 100, 0, 5);
193 fOutputContainer->AddAt(fInvMassVec[i], c++);
194
195 fInvMassDiff[i] = new TH1D(Form("massDiff_%d", i), ";m_{inv} (GeV);", 100, -1, 1);
196 fOutputContainer->AddAt(fInvMassDiff[i], c++);
197
198 fAngleSM[i] = new TH1D(Form("angleSM_%d", i), ";#delta SM", 19, -0.5, 18.5);
199 fOutputContainer->AddAt(fAngleSM[i], c++);
200
201 fPtAngle[i] = new TH2D(Form("ptAngle_%d", i), ";p_{T} (GeV/c);#delta SM",
202 20, 0, 5, 10, -0.5, 9.5);
203 fOutputContainer->AddAt(fPtAngle[i], c++);
204
205
206 for(Int_t j=0; j<2; j++) {
c336eab0 207 fnTracks[j*fgknSteps+i] =
8e2f611a 208 new TH1D(Form("nTracks%s_%d", charge[j],i), Form("%s;number of tracks",charge[j]), 100, -0.5, 99.5);
c336eab0 209 fPt[j*fgknSteps+i] = new TH1D(Form("pt%s_%d", charge[j], i), Form("%s;p_{T} (GeV/c)", charge[j]), 100, 0, 5);
210 fPID[j*fgknSteps+i] = new TH1D(Form("pid%s_%d", charge[j], i), ";electron LQ", 100, 0, 1);
8e2f611a 211
c336eab0 212 fOutputContainer->AddAt(fnTracks[j*fgknSteps+i], c++);
213 fOutputContainer->AddAt(fPt[j*fgknSteps+i], c++);
214 fOutputContainer->AddAt(fPID[j*fgknSteps+i], c++);
8e2f611a 215 }
216 }
217
218 //TH2D *fnGoodTracks;
219
220 printf("n hist = %d\n", c);
221}
222//______________________________________________________________________________
223void AliTRDqaJPsi::Exec(Option_t *)
224{
225 /*
226 Selection steps:
227 - Parameters In and Out
228 - TRDrefit bit
229 - TRDpid and quality
230 - pt > 0.8 GeV/c
231 - PID > 0.9
232 */
233
234
235 // Process one event
236 Long64_t entry = fChain->GetReadEntry() ;
20a6e9c2 237 if (!(entry%100)) Info("Exec", "Entry = %lld", entry);
8e2f611a 238
239 // Processing of one event
240
241 if (!fESD) {
242 //AliError("fESD is not connected to the input!") ;
243 return ;
244 }
245
246
247 Int_t nTracks = fESD->GetNumberOfTracks();
c336eab0 248 Int_t cTracks[2*fgknSteps] = {0,0,0,0,0,0,0,0,0,0};
8e2f611a 249 fnKFtracks = 0;
250
251 // track loop
252 for(Int_t i=0; i<nTracks; i++) {
253
254 //
255 // track selection
256 //
257 // param in and Out
258 // TRDrefit and TRDPid bit
259 //
260
261 AliESDtrack *track = fESD->GetTrack(i);
262 const AliExternalTrackParam *paramOut = track->GetOuterParam();
263 const AliExternalTrackParam *paramIn = track->GetInnerParam();
264
265 // long track ..
266 if (!paramIn) continue;
267 if (!paramOut) continue;
268
269 Int_t step = 0;
270 Int_t charge = (track->Charge() > 0) ? 1 : 0;
271 UInt_t status = track->GetStatus();
272 Double_t pt = track->Pt();
273 Double_t pid = track->GetTRDpid(AliPID::kElectron);
274
275 Double_t esdPid[5];
276 track->GetESDpid(esdPid);
277
278 // create a kalman particle
279 Int_t pdg = (charge == 0)? -11 : 11;
c336eab0 280 for(Int_t k=0; k<fgknSteps; k++) fInSample[fnKFtracks][k] = 0;
8e2f611a 281
7bce990c 282 fVec[fnKFtracks] = CreateVector(track);
283 fTracks[fnKFtracks] = new AliKFParticle(*track, pdg);
8e2f611a 284 fSM[fnKFtracks] = AliTRDqaAT::GetSector(paramOut->GetAlpha());
285 fnKFtracks++;
286
287 //AliTRDqaAT::PrintPID(track);
288
289 // apply the cuts
290
c336eab0 291 cTracks[fgknSteps *charge + step]++;
8e2f611a 292 FillHist(track, step++);
293
294 if (!(status & AliESDtrack::kTRDrefit)) continue;
295
c336eab0 296 cTracks[fgknSteps *charge + step]++;
8e2f611a 297 FillHist(track, step++);
298
299 if (!(status & AliESDtrack::kTRDpid)) continue;
ed15ef4f 300 if (track->GetTRDntracklets() < 6) continue;
8e2f611a 301
c336eab0 302 cTracks[fgknSteps *charge + step]++;
8e2f611a 303 FillHist(track, step++);
304
305 if (pt < 0.8) continue;
306
c336eab0 307 cTracks[fgknSteps *charge + step]++;
8e2f611a 308 FillHist(track, step++);
309
310 if (pid < 0.3) continue; //
311 //if (esdPid[AliPID::kElectron] < 0.5) continue;
312
c336eab0 313 cTracks[fgknSteps *charge + step]++;
8e2f611a 314 FillHist(track, step);
315
c336eab0 316 for(Int_t k=0; k<2*fgknSteps; k++) fnTracks[k]->Fill(cTracks[k]);
8e2f611a 317 }
318
319 // calculate invariant mass
320
c336eab0 321 for(Int_t k=0; k<fgknSteps; k++) {
8e2f611a 322 for(Int_t i=0; i<fnKFtracks; i++) {
323 if (!fInSample[i][k]) continue;
324 for(Int_t j=i+1; j<fnKFtracks; j++) {
325 if (!fInSample[j][k]) continue;
7bce990c 326 AliKFParticle jpsi(*(fTracks[i]), *(fTracks[j]));
327 TLorentzVector jpsiVec = (*(fVec[i])) + (*fVec[j]);
8e2f611a 328 fInvMass[k]->Fill(jpsi.GetMass());
329 fInvMassVec[k]->Fill(jpsiVec.M());
330 fInvMassDiff[k]->Fill(jpsiVec.M() - jpsi.GetMass());
331
332 if (jpsi.GetMass() > 2.5 && jpsi.GetMass() < 3.5) {
333 Int_t dSM = TMath::Abs(fSM[i] - fSM[j]);
334 if (dSM > 9) dSM = (18-dSM);
335 fAngleSM[k]->Fill(dSM);
336 fPtAngle[k]->Fill(TMath::Hypot(jpsi.GetPx(), jpsi.GetPy()), dSM);
337 }
338
339 }
340 }
341 }
342
343 PostData(0, fOutputContainer);
344}
345
346//______________________________________________________________________________
347void AliTRDqaJPsi::Terminate(Option_t *)
348{
349 // save histograms
350 fOutputContainer = (TObjArray*)GetOutputData(0);
351
352 TFile *file = new TFile("outJPsi.root", "RECREATE");
353 fOutputContainer->Write();
354
355 file->Flush();
356 file->Close();
357 delete file;
358
359 //for(Int_t i=0; i<fOutputContainer->GetEntries(); i++) {
360 // TObject *obj = fOu
361 // }
362}
363
364//______________________________________________________________________________
365
366void AliTRDqaJPsi::FillHist(AliESDtrack *track, Int_t step) {
c336eab0 367 //
368 // Fill the histograms
369 //
8e2f611a 370
371 Int_t charge = (track->Charge() > 0) ? 1 : 0;
372 UInt_t status = track->GetStatus();
373 Double_t pt = track->Pt();
374 Double_t pid = track->GetTRDpid(AliPID::kElectron);
375
376 Double_t esdPid[5];
377 track->GetESDpid(esdPid);
378
c336eab0 379 Int_t id = charge * fgknSteps + step;
8e2f611a 380 AliTRDqaAT::FillStatus(fStatus[step], status);
381 fPt[id]->Fill(pt);
382 fPID[id]->Fill(pid);
383 //fPID[id]->Fill(esdPid[AliPID::kElectron]);
384
385 fInSample[fnKFtracks-1][step] = 1;
386}
387
388//______________________________________________________________________________
389
390TLorentzVector *AliTRDqaJPsi::CreateVector(AliESDtrack *track) {
391
392 TLorentzVector *vec = new TLorentzVector();
393 vec->SetXYZM(track->Px(), track->Py(), track->Pz(), 0);
394 return vec;
395}
396
397//______________________________________________________________________________