TOF + macros to each detector folder
[u/mrichter/AliRoot.git] / TRD / qaAnalysis / AliTRDqaEnergyDeposit.cxx
1
2 /**************************************************************************
3  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4  *                                                                        *
5  * Author: The ALICE Off-line Project.                                    *
6  * Contributors are mentioned in the code where appropriate.              *
7  *                                                                        *
8  * Permission to use, copy, modify and distribute this software and its   *
9  * documentation strictly for non-commercial purposes is hereby granted   *
10  * without fee, provided that the above copyright notice appears in all   *
11  * copies and that both the copyright notice and this permission notice   *
12  * appear in the supporting documentation. The authors make no claims     *
13  * about the suitability of this software for any purpose. It is          *
14  * provided "as is" without express or implied warranty.                  *
15  **************************************************************************/
16
17 /* $Id: AliTRDqaEnergyDeposit.cxx $ */
18
19 //
20 // This class is a part of a package of high level QA monitoring for TRD.
21 // In this class the dEdX is analyzed as a function of the particle type,
22 // momentum and chamber.
23 //
24 // S. Radomski
25 // radomski@physi.uni-heidelberg.de
26 // March 2008
27 //
28
29 #include "AliTRDqaEnergyDeposit.h"
30
31 #include "TMath.h"
32 #include "TH1D.h"
33 #include "TH2D.h"
34 #include "TFile.h"
35 #include "TChain.h"
36
37 #include "AliESDEvent.h"
38 #include "AliESDtrack.h"
39
40 //______________________________________________________________________________
41
42 AliTRDqaEnergyDeposit::AliTRDqaEnergyDeposit() 
43   : AliAnalysisTask("",""),  
44     fChain(0),
45     fESD(0),
46     fOutputContainer(0)
47 {
48   //
49   // Default constructor
50   //
51
52   for (Int_t i = 0; i < 2; i++) {
53     fSignalPtSum[i] = 0x0;
54   }
55
56   for (Int_t j = 0; j < 2*5; j++) {
57     fSignalPtType[j] = 0x0;
58     fSignalPtPure[j] = 0x0;
59     fProb[j]         = 0x0;
60   }
61
62 }
63
64 //______________________________________________________________________________
65
66 AliTRDqaEnergyDeposit:: AliTRDqaEnergyDeposit(const AliTRDqaEnergyDeposit & /*trd*/)
67   : AliAnalysisTask("",""),  
68     fChain(0),
69     fESD(0),
70     fOutputContainer(0)
71 {
72   //
73   // Copy constructor
74   //
75
76   for (Int_t i = 0; i < 2; i++) {
77     fSignalPtSum[i] = 0x0;
78   }
79
80   for (Int_t j = 0; j < 2*5; j++) {
81     fSignalPtType[j] = 0x0;
82     fSignalPtPure[j] = 0x0;
83     fProb[j]         = 0x0;
84   }
85
86 }
87
88 //______________________________________________________________________________
89 AliTRDqaEnergyDeposit::AliTRDqaEnergyDeposit(const char *name) 
90   : AliAnalysisTask(name,""),  
91     fChain(0),
92     fESD(0),
93     fOutputContainer(0)
94 {
95   //
96   // Constructor.
97   //
98
99   // Input slot #0 works with an Ntuple
100   DefineInput(0, TChain::Class());
101   // Output slot #0 writes into a TH1 container
102   DefineOutput(0,  TObjArray::Class()) ; 
103
104   for (Int_t i = 0; i < 2; i++) {
105     fSignalPtSum[i] = 0x0;
106   }
107
108   for (Int_t j = 0; j < 2*5; j++) {
109     fSignalPtType[j] = 0x0;
110     fSignalPtPure[j] = 0x0;
111     fProb[j]         = 0x0;
112   }
113
114 }
115
116 //______________________________________________________________________________
117 void AliTRDqaEnergyDeposit::ConnectInputData(const Option_t *)
118 {
119   // Initialisation of branch container and histograms 
120
121   //AliInfo(Form("*** Initialization of %s", GetName())) ; 
122
123   fChain = (TChain*)GetInputData(0);
124   fESD = new AliESDEvent();
125   fESD->ReadFromTree(fChain);
126 }
127
128 //________________________________________________________________________
129 void AliTRDqaEnergyDeposit::CreateOutputObjects()
130 {
131   // build histograms
132   
133   // prepare the scale from 0.5 to 10 GeV
134   const Int_t knbinsx = 50;
135   Double_t scalex[knbinsx+1];
136   Double_t dd = (TMath::Log(10) - TMath::Log(0.5)) / knbinsx;
137   for(Int_t ix=0; ix<knbinsx+1; ix++) {
138     scalex[ix] = 0.5 * TMath::Exp(dd * ix);
139   }
140
141   const Int_t knbinsy = 50;
142   Double_t scaley[knbinsy+1];
143   for(Int_t iy=0; iy<knbinsy+1; iy++) {
144     scaley[iy] = iy * (3e3/100.);
145   }
146   
147   const char *title = ";p_{T};dEdX (a. u.)";
148   const char *charge[2] = {"Pos", "Neg"};
149
150   // build histograms
151   fOutputContainer = new TObjArray(50);
152   Int_t c=0;
153
154   for(Int_t i=0; i<2; i++) {
155
156     fSignalPtSum[i] = new TH2D(Form("ptSig%s", charge[i]), title, knbinsx, scalex, knbinsy, scaley);
157     fOutputContainer->AddAt(fSignalPtSum[i], c++);
158
159     for(Int_t j=0; j<AliPID::kSPECIES; j++) {
160       Int_t idx = AliPID::kSPECIES*i+j;
161
162       //
163       fSignalPtType[idx] = 
164         new TH2D(Form("ptSig%s%d", charge[i], j), title, knbinsx, scalex, knbinsy, scaley);
165       fOutputContainer->AddAt(fSignalPtType[idx], c++);
166
167       //
168       fSignalPtPure[idx] = 
169         new TH2D(Form("ptSigPure%s%d", charge[i], j), title, knbinsx, scalex, knbinsy, scaley);
170       fOutputContainer->AddAt(fSignalPtPure[idx], c++);
171       
172       fProb[idx] = new TH1D(Form("prob%s%d", charge[i], j), ";LQ", 100, 0, 1);
173       fOutputContainer->AddAt(fProb[idx], c++);
174     }  
175   }
176
177   printf("n hist = %d\n", c);
178 }
179
180 //______________________________________________________________________________
181 void AliTRDqaEnergyDeposit::Exec(Option_t *) 
182 {
183   //
184   // Process one event
185   //
186
187   //  Long64_t entry = fChain->GetReadEntry() ;
188   
189   // Processing of one event 
190    
191   if (!fESD) {
192     //AliError("fESD is not connected to the input!") ; 
193     return ; 
194   }
195   
196   Int_t nTracks = fESD->GetNumberOfTracks();
197   //fNTracks->Fill(nTracks); 
198
199   // track loop
200   for(Int_t i=0; i<nTracks; i++) {
201     
202     //
203     // track selection 
204     //
205     // param in and Out
206     // TRDrefit and TRDPid bit
207     //
208  
209     AliESDtrack *track = fESD->GetTrack(i);
210     const AliExternalTrackParam *paramOut = track->GetOuterParam();
211     const AliExternalTrackParam *paramIn = track->GetInnerParam();
212
213     // long track ..
214     if (!paramIn) continue;
215     if (!paramOut) continue;
216     
217     UInt_t status = track->GetStatus();
218     if (!(status & AliESDtrack::kTRDrefit)) continue;
219     if (!(status & AliESDtrack::kTRDpid)) continue;
220     if (track->GetTRDntrackletsPID() < 6) continue;
221
222     // standard selection
223     
224     Int_t idx = (track->GetSign() > 0) ? 0 : 1;
225     Double_t pt = paramOut->Pt();
226
227     Double_t signal = 0;
228     for(Int_t k=0; k<6; ++k)
229       signal += track->GetTRDslice(k, -1);
230     signal /= 6;
231
232     fSignalPtSum[idx]->Fill(pt, signal);
233     
234     for(Int_t k=0; k<AliPID::kSPECIES; ++k) {
235       
236       Double_t lq = track->GetTRDpid(k);
237       fProb[AliPID::kSPECIES*idx+k]->Fill(lq);
238       if (lq > 0.8) fSignalPtType[AliPID::kSPECIES*idx+k]->Fill(pt, signal);
239     }
240   }
241
242   PostData(0, fOutputContainer);
243 }
244
245 //______________________________________________________________________________
246 void AliTRDqaEnergyDeposit::FillElectrons() const 
247 {
248   //
249   // Fill the electrons
250   //  
251
252 }
253
254 //______________________________________________________________________________
255 void AliTRDqaEnergyDeposit::Terminate(Option_t *)
256 {
257   //
258   // Retrieve histograms
259   //
260
261   fOutputContainer = (TObjArray*)GetOutputData(0);
262   
263   // post processing
264
265
266   // save the results
267   TFile *file = new TFile("outEnergyDeposit.root", "RECREATE");
268   fOutputContainer->Write();
269  
270   file->Flush();
271   file->Close();
272   delete file;
273
274   //for(Int_t i=0; i<fOutputContainer->GetEntries(); i++) {
275   //  TObject *obj = fOu
276   // }
277 }
278
279 //______________________________________________________________________________