New support for QA histos
[u/mrichter/AliRoot.git] / CORRFW / AliCFEventRecCuts.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 // Cut on the Event at reconstructed level: for the moment 
16 // the requirements on the number of charged tracks and on 
17 // the vertex position and resolution are implemented
18 // The argument of IsSelected member function (passed object) is cast into 
19 // an AliESDEvent. In the future may be modified to use AliVEvent interface
20 // and include more cut variables.
21 // The class derives from AliCFCutBase
22 // Author:S.Arcelli Silvia.Arcelli@cern.ch
23 //
24 //
25 #include "TH1F.h"
26 #include "TList.h"
27 #include "TBits.h"
28 #include "AliLog.h"
29 #include "AliESDEvent.h"
30 #include "AliESDVertex.h"
31 #include "AliCFEventRecCuts.h"
32 ClassImp(AliCFEventRecCuts) 
33 //____________________________________________________________________
34 AliCFEventRecCuts::AliCFEventRecCuts() : 
35   AliCFCutBase(),
36   fNTracksMin(-1),
37   fNTracksMax(1000000),
38   fRequireVtxCuts(kFALSE),
39   fVtxXMax(1.e99),
40   fVtxYMax(1.e99),
41   fVtxZMax(1.e99),
42   fVtxXMin(-1.e99),
43   fVtxYMin(-1.e99),
44   fVtxZMin(-1.e99),
45   fVtxXResMax(1.e99),
46   fVtxYResMax(1.e99),
47   fVtxZResMax(1.e99),
48   fBitMap(0x0)
49 {
50   //
51   //ctor
52   //
53   fBitMap=new TBits(0);
54   Initialise();
55 }
56
57 //____________________________________________________________________
58 AliCFEventRecCuts::AliCFEventRecCuts(Char_t* name, Char_t* title) : 
59   AliCFCutBase(name,title),
60   fNTracksMin(-1),
61   fNTracksMax(1000000),
62   fRequireVtxCuts(kFALSE),
63   fVtxXMax(1.e99),
64   fVtxYMax(1.e99),
65   fVtxZMax(1.e99),
66   fVtxXMin(-1.e99),
67   fVtxYMin(-1.e99),
68   fVtxZMin(-1.e99),
69   fVtxXResMax(1.e99),
70   fVtxYResMax(1.e99),
71   fVtxZResMax(1.e99),
72   fBitMap(0x0)
73  {
74   //
75   //ctor
76   //
77   fBitMap=new TBits(0);
78   Initialise();
79  }
80
81 //____________________________________________________________________
82 AliCFEventRecCuts::AliCFEventRecCuts(const AliCFEventRecCuts& c) : 
83   AliCFCutBase(c),
84   fNTracksMin(c.fNTracksMin),
85   fNTracksMax(c.fNTracksMax),
86   fRequireVtxCuts(c.fRequireVtxCuts),
87   fVtxXMax(c.fVtxXMax),
88   fVtxYMax(c.fVtxYMax),
89   fVtxZMax(c.fVtxZMax),
90   fVtxXMin(c.fVtxXMin),
91   fVtxYMin(c.fVtxYMin),
92   fVtxZMin(c.fVtxZMin),
93   fVtxXResMax(c.fVtxXResMax),
94   fVtxYResMax(c.fVtxYResMax),
95   fVtxZResMax(c.fVtxZResMax),
96   fBitMap(c.fBitMap)
97 {
98   //
99   //copy constructor
100   //
101   for (Int_t i=0; i<c.kNCuts; i++){
102     for (Int_t j=0; j<c.kNStepQA; j++){
103       if(c.fhQA[i][j]) fhQA[i][j] = (TH1F*)c.fhQA[i][j]->Clone();
104     }
105   }
106
107 }
108
109 //____________________________________________________________________
110 AliCFEventRecCuts::~AliCFEventRecCuts() {
111   //
112   //dtor
113   //
114
115   for (Int_t i=0; i<kNCuts; i++){
116     for (Int_t j=0; j<kNStepQA; j++){
117       if(fhQA[i][j]) delete fhQA[i][j];
118     }
119   }
120
121   if(fBitMap)delete fBitMap;
122
123 }
124 //_____________________________________________________________________________
125 void AliCFEventRecCuts::Initialise()
126 {
127
128   //
129   //initialization
130   //
131
132   //
133   // sets pointers to histos to zero
134   //
135
136   for(Int_t i=0; i<kNCuts; i++){
137     for(Int_t j =0; j<kNStepQA; j++){
138       fhQA[i][j]=0x0;
139     }
140   }
141 }
142
143 //____________________________________________________________________
144 AliCFEventRecCuts& AliCFEventRecCuts::operator=(const AliCFEventRecCuts& c)
145 {
146   //
147   // Assignment operator
148   //
149   if (this != &c) {
150     AliCFCutBase::operator=(c) ;
151     fNTracksMin=c.fNTracksMin;
152     fNTracksMax=c.fNTracksMax;
153     fRequireVtxCuts=c.fRequireVtxCuts;
154     fVtxXMax=c.fVtxXMax;
155     fVtxYMax=c.fVtxYMax;
156     fVtxZMax=c.fVtxZMax;
157     fVtxXMin=c.fVtxXMin;
158     fVtxYMin=c.fVtxYMin;
159     fVtxZMin=c.fVtxZMin;
160     fVtxXResMax=c.fVtxXResMax;
161     fVtxYResMax=c.fVtxYResMax;
162     fVtxZResMax=c.fVtxZResMax;
163     fBitMap=c.fBitMap;
164   }
165
166   for (Int_t i=0; i<c.kNCuts; i++){
167     for (Int_t j=0; j<c.kNStepQA; j++){
168       if(c.fhQA[i][j]) fhQA[i][j] = (TH1F*)c.fhQA[i][j]->Clone();
169     }
170   }
171
172
173   return *this ;
174 }
175
176 //____________________________________________________________________
177 Bool_t AliCFEventRecCuts::IsSelected(TObject* obj) {
178   //
179   //Check if the requested cuts are passed
180   //
181
182
183   SelectionBitMap(obj);
184
185   if (fIsQAOn) FillHistograms(obj,0);
186   Bool_t isSelected = kTRUE;
187
188   for (UInt_t icut=0; icut<fBitMap->GetNbits();icut++)
189         if(!fBitMap->TestBitNumber(icut)) isSelected = kFALSE;
190
191   if (!isSelected) return kFALSE ;
192   if (fIsQAOn) FillHistograms(obj,1);
193   return kTRUE;
194
195 }
196 //____________________________________________________________________
197 void AliCFEventRecCuts::SelectionBitMap(TObject* obj) {
198   //
199   //cut on the number of charged tracks and on the event vertex.
200   //so far specific to AliESDEvents
201   //
202
203   //Check if the requested cuts are passed and return a bitmap
204   for(Int_t j=0;j<kNCuts;j++)fBitMap->SetBitNumber(j,kFALSE);
205   AliESDEvent* esd = dynamic_cast<AliESDEvent *>(obj);
206   if ( !esd ) return;
207
208   //now start checking the cuts,
209   //first assume the event will be accepted: 
210   for(Int_t j=0;j<kNCuts;j++)fBitMap->SetBitNumber(j,kTRUE);
211
212   //Number of charged tracks:
213   Int_t nTracks = esd->GetNumberOfTracks();
214   if(nTracks<fNTracksMin || nTracks>fNTracksMax)
215     fBitMap->SetBitNumber(0,kFALSE); 
216   
217   if(fRequireVtxCuts){
218     const AliESDVertex* vtxESD = esd->GetVertex();
219     if(!vtxESD){
220       for(Int_t j=1;j<kNCuts;j++)fBitMap->SetBitNumber(j,kFALSE); 
221       return;
222     }
223     // Require the vertex to have been reconstructed successfully
224     if (strcmp(vtxESD->GetName(), "default")==0){
225       AliWarning(Form(" No reconstructed vertex found, skip event"));    
226       for(Int_t j=1;j<kNCuts;j++)fBitMap->SetBitNumber(j,kFALSE); 
227       return;
228     }    
229     // Pick up the position and uncertainties
230     
231     Double_t vtxPos[3];
232     vtxPos[0] = vtxESD->GetXv();
233     vtxPos[1] = vtxESD->GetYv();
234     vtxPos[2] = vtxESD->GetZv();
235     
236     Double_t vtxRes[3];
237     vtxRes[0] = vtxESD->GetXRes();
238     vtxRes[1] = vtxESD->GetYRes();
239     vtxRes[2] = vtxESD->GetZRes();
240  
241     // Apply the cut
242     
243     if (vtxPos[0]>fVtxXMax || vtxPos[0]<fVtxXMin)
244       fBitMap->SetBitNumber(1,kFALSE); 
245     if (vtxPos[1]>fVtxYMax || vtxPos[1]<fVtxYMin)
246       fBitMap->SetBitNumber(2,kFALSE); 
247     if (vtxPos[2]>fVtxZMax || vtxPos[2]<fVtxZMin)
248       fBitMap->SetBitNumber(3,kFALSE); 
249     if (vtxRes[0]==0 || vtxRes[0]>fVtxXResMax)
250       fBitMap->SetBitNumber(4,kFALSE); 
251     if (vtxRes[1]==0 || vtxRes[1]>fVtxYResMax)
252       fBitMap->SetBitNumber(5,kFALSE); 
253     if (vtxRes[2]==0 || vtxRes[2]>fVtxZResMax)
254       fBitMap->SetBitNumber(6,kFALSE); 
255   }  
256   return;
257 }
258
259 //_____________________________________________________________________________
260 void AliCFEventRecCuts::FillHistograms(TObject* obj, Bool_t b)
261 {
262   //
263   // fill the QA histograms
264   //
265
266   if(!fIsQAOn) return;
267   // cast TObject into VParticle
268   AliESDEvent* esd = dynamic_cast<AliESDEvent *>(obj);
269   if (!esd ) return  ;
270
271   // index = 0: fill histograms before cuts
272   // index = 1: fill histograms after cuts
273   Int_t index = -1;
274   index = ((b) ? 1 : 0);
275
276
277   //number of charged tracks:
278   Int_t nTracks = esd->GetNumberOfTracks();
279   fhQA[kNTracks][index]->Fill(nTracks);
280
281   //look at vertex parameters:
282   const AliESDVertex* vtxESD = esd->GetVertex();
283   if(!vtxESD)return;
284   // Require the vertex to have been reconstructed successfully
285   if (strcmp(vtxESD->GetName(), "default")==0)return;
286   // vertex position and uncertainties
287   fhQA[kVtxPosX][index]->Fill(vtxESD->GetXv());
288   fhQA[kVtxPosY][index]->Fill(vtxESD->GetYv());
289   fhQA[kVtxPosZ][index]->Fill(vtxESD->GetZv());
290   fhQA[kVtxResX][index]->Fill(vtxESD->GetXRes());
291   fhQA[kVtxResY][index]->Fill(vtxESD->GetYRes());
292   fhQA[kVtxResZ][index]->Fill(vtxESD->GetZRes());
293   
294 }
295
296 //____________________________________________________________________
297 void AliCFEventRecCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t *bins)
298 {
299   //
300   //setting x-axis bin limits of QA histogram fhQA[index] 
301   // 
302
303   for(Int_t i=0;i<kNStepQA;i++){
304     if(!fhQA[index][i]){AliWarning("non-existing histogram!");
305     return;
306     }
307     fhQA[index][i]->GetXaxis()->Set(nbins,bins);
308   }
309 }
310 //____________________________________________________________________
311 void AliCFEventRecCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t xmin, Double_t xmax)
312 {
313   //
314   //setting x-axis bins and range of QA histogram fhQA[index] 
315   // 
316
317   for(Int_t i=0;i<kNStepQA;i++){
318     if(!fhQA[index][i]){AliWarning("non-existing histogram!");
319     return;
320     }
321     fhQA[index][i]->GetXaxis()->Set(nbins,xmin,xmax);
322   }
323 }
324
325 //_____________________________________________________________________________
326  void AliCFEventRecCuts::DefineHistograms() {
327   //
328   // histograms for cut variables
329   //
330   Int_t color = 2;
331
332   if(!fIsQAOn) {
333     AliInfo(Form("No QA histos requested, Please first set the QA flag on!"));
334     return;
335   }  
336   
337   // book QA histograms
338
339   Char_t str[256];
340   for (Int_t i=0; i<kNStepQA; i++) {
341     if (i==0) sprintf(str," ");
342     else sprintf(str,"_cut");
343
344     fhQA[kNTracks][i]   = new  TH1F(Form("%s_NTracks%s",GetName(),str),                 "",501,-0.5,500.5);
345     fhQA[kVtxPosX][i]   = new  TH1F(Form("%s_Vtx_Pos_X%s",GetName(),str),               "",100,-5.,5.);
346     fhQA[kVtxPosY][i]   = new  TH1F(Form("%s_Vtx_Pos_Y%s",GetName(),str),               "",100,-5.,5.);
347     fhQA[kVtxPosZ][i]   = new  TH1F(Form("%s_Vtx_Pos_Z%s",GetName(),str),               "",200,-50.,50.);
348
349     fhQA[kVtxResX][i]   = new  TH1F(Form("%s_Vtx_Res_X%s",GetName(),str),               "",100,-1.,1.);
350     fhQA[kVtxResY][i]   = new  TH1F(Form("%s_Vtx_Res_Y%s",GetName(),str),               "",100,-1.,1.);
351     fhQA[kVtxResZ][i]   = new  TH1F(Form("%s_Vtx_Res_Z%s",GetName(),str),               "",100,-1.,1.);
352  
353     fhQA[kNTracks][i]   ->SetXTitle("Number of ESD tracks");
354     fhQA[kVtxPosX][i]   ->SetXTitle("Vertex Position X (cm)");
355     fhQA[kVtxPosY][i]   ->SetXTitle("Vertex Position Y (cm)");
356     fhQA[kVtxPosZ][i]   ->SetXTitle("Vertex Position Z (cm)");
357     fhQA[kVtxResX][i]   ->SetXTitle("Vertex Resolution X (cm)");
358     fhQA[kVtxResY][i]   ->SetXTitle("Vertex Resolution Y (cm)");
359     fhQA[kVtxResZ][i]   ->SetXTitle("Vertex Resolution Z (cm)");
360
361   }
362
363   for(Int_t i=0; i<kNCuts; i++) fhQA[i][1]->SetLineColor(color);
364
365 }
366
367 //_____________________________________________________________________________
368 void AliCFEventRecCuts::AddQAHistograms(TList *qaList) {
369   //
370   // saves the histograms in a TList
371   //
372
373   DefineHistograms();
374
375   for (Int_t j=0; j<kNStepQA; j++) {
376     for(Int_t i=0; i<kNCuts; i++)
377         qaList->Add(fhQA[i][j]);
378   }
379 }