Adding the TPC QA task
[u/mrichter/AliRoot.git] / PWG1 / AliTPCtaskQA.cxx
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 //
17 // 
18 //Task for analysis if TPC base QA 
19
20 // QA histogram
21 // Dimensions
22 // 0 - chi2
23 // 1 - number of clusters
24 // 2 - number of findable clusters
25 // 3 - number of clusters/ findable clusters  
26 // 4 - pt          - at the entrance of the TPC
27 // 5 - eta         - at the entrance of the TPC
28 // 6 - phi         - at the entrance of the TPC
29
30 //-----------------------------------------------------------------------
31 // Author : M.Ivanov  marian.ivanov@cern.ch - 
32 //-----------------------------------------------------------------------
33
34
35
36 //
37 // ROOT includes
38 #include <TChain.h>
39 #include <TMath.h>
40 #include <TVectorD.h>
41 #include <TSystem.h>
42 #include <TFile.h>
43 #include <TParticle.h>
44
45 // ALIROOT includes
46 #include <AliAnalysisManager.h>
47 #include <AliESDInputHandler.h>
48 #include "AliStack.h"
49 #include "AliMCEvent.h"
50 #include "AliMCEventHandler.h"
51 #include "AliMathBase.h"
52
53 #include <AliESD.h>
54 #include "AliExternalTrackParam.h"
55 #include "AliTracker.h"
56 #include "AliTPCseed.h"
57 //
58 #include "AliTPCtaskQA.h"
59 //
60 #include <THnSparse.h>
61
62 //
63
64 // STL includes
65 #include <iostream>
66
67 using namespace std;
68
69
70 ClassImp(AliTPCtaskQA)
71
72 //________________________________________________________________________
73 AliTPCtaskQA::AliTPCtaskQA() : 
74   AliAnalysisTask(), 
75   fMCinfo(0),     //! MC event handler
76   fESD(0),
77   fList(0),
78   fTPCqa(0)
79 {
80   //
81   // Default constructor (should not be used)
82   //
83 }
84
85 AliTPCtaskQA::AliTPCtaskQA(const AliTPCtaskQA& info) : 
86   AliAnalysisTask(info), 
87   fMCinfo(info.fMCinfo),     //! MC event handler
88   fESD(info.fESD),        //!
89   fList(0),
90   fTPCqa(0)
91 {
92   //
93   // Dummy Copy  constructor - no copy constructor for THnSparse 
94   //
95   fList = (TObjArray*)(info.fList->Clone());
96 }
97
98
99
100 //________________________________________________________________________
101 AliTPCtaskQA::AliTPCtaskQA(const char *name) : 
102   AliAnalysisTask(name, "AliTPCtaskQA"), 
103   fMCinfo(0),     //! MC event handler
104   fESD(0),
105   fList(0),
106   fTPCqa(0)
107 {
108   //
109   // Normal constructor
110   //
111   // Input slot #0 works with a TChain
112   DefineInput(0, TChain::Class());
113   // Output slot #0 writes into a TList
114   DefineOutput(0, TObjArray::Class());
115   //
116   //make histos
117   Init(); 
118 }
119
120 void AliTPCtaskQA::Init(){
121   //
122   // Init qa histogram
123   // Dimensions
124   //
125   // 0 - chi2
126   // 1 - number of clusters
127   // 2 - number of findable clusters
128   // 3 - number of clusters/ findable clusters  
129   // 4 - pt          - at the entrance of the TPC
130   // 5 - eta         - at the entrance of the TPC
131   // 6 - phi         - at the entrance of the TPC
132   
133
134
135   Double_t xmin[7],  xmax[7];
136   Int_t    nbins[7];
137   // 
138   nbins[0]=100;                // chi2
139   xmin[0]=0; xmax[0]=10;
140   //
141   nbins[1]=80;                 // ncls
142   xmin[1]=0; xmax[1]=160;
143   //
144   nbins[2]=80;                 // nclsf
145   xmin[2]=0; xmax[2]=160;
146   //
147   nbins[3]=40;                 // ncls/nclsf
148   xmin[3] =-0.1; xmax[3]=1.1;
149   //
150   nbins[4]=50;                 // pt
151   xmin[4] =0.1; xmax[4]=100;
152
153   nbins[5]=40;                 // eta
154   xmin[5] =-2; xmax[5]=2;
155
156   nbins[6]= 360;                 // phi - 10 bins per sector
157   xmin[6] = -TMath::Pi(); xmax[6]=TMath::Pi();
158
159
160
161   fTPCqa = new THnSparseS("TPC qa","TPC qa",7,nbins,xmin,xmax);
162   //
163   //
164   BinLogX(fTPCqa->GetAxis(4));
165   
166   char *hisAxisName[7] ={"chi2/N_{cl}","N_{cl}","N_{clF}","N_{clR}","p_{t}","#eta","#phi"};
167   //  
168   for (Int_t i=0;i<7;i++) {
169     fTPCqa->GetAxis(i)->SetTitle(hisAxisName[i]);
170     fTPCqa->GetAxis(i)->SetName(hisAxisName[i]);
171   }
172   fList = new TObjArray(3);
173   fList->AddAt(fTPCqa,0);
174 }
175
176
177
178
179 AliTPCtaskQA::~AliTPCtaskQA(){
180   //
181   //
182   //
183   delete fTPCqa;
184 }
185
186
187 //________________________________________________________________________
188 void AliTPCtaskQA::ConnectInputData(Option_t *) 
189 {
190   //
191   // Connect the input data
192   //
193
194   TTree* tree=dynamic_cast<TTree*>(GetInputData(0));
195   if (!tree) {
196     //Printf("ERROR: Could not read chain from input slot 0");
197   }
198   else {
199     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
200     if (!esdH) {
201       //Printf("ERROR: Could not get ESDInputHandler");
202     }
203     else {
204       fESD = esdH->GetEvent();
205       //Printf("*** CONNECTED NEW EVENT ****");
206     }  
207   }
208   AliMCEventHandler* mcinfo = (AliMCEventHandler*) (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());  
209   mcinfo->SetReadTR(kTRUE);
210   
211   fMCinfo = mcinfo->MCEvent();
212
213
214 }
215
216
217
218
219
220
221 //________________________________________________________________________
222 void AliTPCtaskQA::CreateOutputObjects() 
223 {
224   //
225   // Connect the output objects
226   //
227
228 }
229
230
231 //________________________________________________________________________
232 void AliTPCtaskQA::Exec(Option_t *) {
233   //
234   // Execute analysis for current event 
235   //
236
237
238   // If MC has been connected   
239
240   if (!fMCinfo){
241     cout << "Not MC info\n" << endl;
242   }else{
243     ProcessMCInfo();
244   }
245   //
246   PostData(0, fList);
247 }      
248
249
250
251
252
253
254
255
256
257
258 void  AliTPCtaskQA::ProcessMCInfo(){
259   //
260   //
261   //
262   //
263
264   if (!fTPCqa) Init();
265   Int_t npart   = fMCinfo->GetNumberOfTracks();
266   Int_t ntracks = fESD->GetNumberOfTracks(); 
267   if (npart<=0) return;
268   if (ntracks<=0) return;
269   //
270   //
271   TParticle * particle= new TParticle;
272   TClonesArray * trefs = new TClonesArray("AliTrackReference");
273   
274   for (Int_t itrack=0;itrack<ntracks;itrack++){
275     AliESDtrack *track = fESD->GetTrack(itrack);
276     //
277     if ((track->GetStatus()&AliESDtrack::kTPCrefit)==0) continue;   // only refited tracks
278     Int_t ipart = TMath::Abs(track->GetLabel());
279     //
280     Int_t status = fMCinfo->GetParticleAndTR(ipart, particle, trefs);
281     if (status<0) continue;
282     if (!particle) continue;
283     if (!trefs) continue;
284     //
285     //
286     AliTrackReference *tpcRef=0;
287     for (Int_t iref=0; iref<trefs->GetEntries(); iref++){
288       AliTrackReference *ref = (AliTrackReference*)trefs->At(iref);
289       if (ref->DetectorId()== AliTrackReference::kTPC){
290         tpcRef=ref;
291         break;
292       }
293     }
294     if (!tpcRef) continue;
295     
296     //
297     // Fill histos
298     //
299     Double_t x[7];
300     x[0]= track->GetTPCchi2()/track->GetTPCNcls();
301     x[1]= track->GetTPCNcls();
302     x[2]= track->GetTPCNclsF();
303     x[3]= Float_t(track->GetTPCNcls())/Float_t(track->GetTPCNclsF());
304     x[4]= tpcRef->Pt();
305     x[5]= -0.5*TMath::Log((tpcRef->P()+tpcRef->Pz())/(tpcRef->P()-tpcRef->Pz()));
306     x[6]= TMath::ATan2(tpcRef->Y(),tpcRef->X());
307     //
308     fTPCqa->Fill(x);
309   }
310 }
311
312
313
314
315
316
317 void AliTPCtaskQA::BinLogX(TAxis *axis) {
318   //
319   //
320   //
321   Int_t bins = axis->GetNbins();
322   
323   Double_t from = axis->GetXmin();
324   Double_t to   = axis->GetXmax();
325   Double_t *new_bins = new Double_t[bins + 1];
326   
327   new_bins[0] = from;
328   Double_t factor = pow(to/from, 1./bins);
329   
330   for (Int_t i = 1; i <= bins; i++) {
331     new_bins[i] = factor * new_bins[i-1];
332   }
333   axis->Set(bins, new_bins);
334   delete new_bins;
335 }
336
337
338 AliTPCtaskQA* AliTPCtaskQA::ReadFromFile(const char *fname){
339   //
340   //
341   //
342   AliTPCtaskQA *qa = new AliTPCtaskQA;
343   TFile *f = new TFile(fname);
344   TObjArray *array= (TObjArray*)f->Get("tpcTaskQA");
345   qa->fTPCqa =  (THnSparse*)array->At(0);  
346   delete f;
347   return qa;
348 }
349