Fixes for bug #49914: Compilation breaks in trunk, and bug #48629: Trunk cannot read...
[u/mrichter/AliRoot.git] / TPC / AliTPCQADataMakerRec.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-2007, 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 /* $Id: $ */
18
19 /*
20   Based on AliPHOSQADataMaker
21   Produces the data needed to calculate the quality assurance. 
22   All data must be mergeable objects.
23   P. Christiansen, Lund, January 2008
24 */
25
26 /*
27   Implementation:
28
29   We have chosen to have the histograms as non-persistent meber to
30   allow better debugging. In the copy constructor we then have to
31   assign the pointers to the existing histograms in the copied
32   list. This have been implemented but not tested.
33
34   For the QA of the RAW data we use the class, AliTPCdataQA, from the
35   existing TPC Calibration framework (which is more advanced than the
36   standard QA framework) and extract the histograms at the end. This
37   has been tested with zero-suppressed data. The Analyse method of the
38   AliTPCdataQA class is called in the method, EndOfDetectorCycle, and
39   there also: 1d histogram(s) are projected and added to the QA list.
40 */
41
42 /*
43   TODO:
44   Sumw2 for RAW histogram(s)?
45   RecPoints and ESD could have many more histograms
46 */
47
48 #include "AliTPCQADataMakerRec.h"
49
50 // --- ROOT system ---
51 #include <TClonesArray.h>
52 #include <TString.h>
53 #include <TSystem.h>
54
55 // --- Standard library ---
56
57 // --- AliRoot header files ---
58 #include "AliQAChecker.h"
59 #include "AliESDEvent.h"
60 #include "AliESDtrack.h"
61 #include "AliLog.h"
62 #include "AliTPCCalPad.h"
63 #include "AliTPCCalROC.h"
64 #include "AliTPCClustersRow.h"
65 #include "AliTPCclusterMI.h"
66
67 ClassImp(AliTPCQADataMakerRec)
68
69 //____________________________________________________________________________ 
70 AliTPCQADataMakerRec::AliTPCQADataMakerRec() : 
71   AliQADataMakerRec(AliQAv1::GetDetName(AliQAv1::kTPC), 
72                     "TPC Rec Quality Assurance Data Maker"),
73   fTPCdataQA(NULL)
74 {
75   // ctor
76   fTPCdataQA = new AliTPCdataQA*[AliRecoParam::kNSpecies] ;
77   for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) 
78     fTPCdataQA[specie] = NULL ; 
79   
80   for(Int_t i = 0; i < 6; i++)
81     fMapping[i] = 0;
82 }
83
84 //____________________________________________________________________________ 
85 AliTPCQADataMakerRec::AliTPCQADataMakerRec(const AliTPCQADataMakerRec& qadm) :
86   AliQADataMakerRec(),
87   fTPCdataQA(NULL)
88 {
89   //copy ctor 
90   // Does not copy the calibration object, instead InitRaws have to be
91   // called again
92   SetName((const char*)qadm.GetName()) ; 
93   SetTitle((const char*)qadm.GetTitle()); 
94
95   fTPCdataQA = new AliTPCdataQA*[AliRecoParam::kNSpecies] ;
96   for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) 
97     fTPCdataQA[specie] = NULL ; 
98   
99   for(Int_t i = 0; i < 6; i++)
100     fMapping[i] = 0;
101
102   //
103   // Associate class histogram objects to the copies in the list
104   // Could also be done with the indexes
105   //
106
107 }
108
109 //__________________________________________________________________
110 AliTPCQADataMakerRec& AliTPCQADataMakerRec::operator = (const AliTPCQADataMakerRec& qadm )
111 {
112   // Equal operator.
113   this->~AliTPCQADataMakerRec();
114   new(this) AliTPCQADataMakerRec(qadm);
115   return *this;
116 }
117
118 //__________________________________________________________________
119 AliTPCQADataMakerRec::~AliTPCQADataMakerRec()
120 {
121   // Destructor
122   for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) 
123     if ( fTPCdataQA[specie] != NULL )
124       delete fTPCdataQA[specie] ; 
125   delete[] fTPCdataQA; 
126
127   for(Int_t i = 0; i < 6; i++) 
128     delete fMapping[i];
129 }
130  
131 //____________________________________________________________________________ 
132 void AliTPCQADataMakerRec::EndOfDetectorCycle(AliQAv1::TASKINDEX_t task, TObjArray ** list)
133 {
134   //Detector specific actions at end of cycle
135   
136   for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
137     if ( !AliQAv1::Instance()->IsEventSpecieSet(specie) ) 
138       continue ; 
139     if(fTPCdataQA[specie] != NULL) { // do the final step of the QA for Raw data
140
141       fTPCdataQA[specie]->Analyse(); // 31/1-08 Analyse is now protected against
142                            //         RAW data files with no TPC data
143
144       SetEventSpecie(specie) ; 
145       TH1F * histRawsOccupancy                 = (TH1F*)GetRawsData(kOccupancy) ;
146       TH1F * histRawsOccupancyVsSector         = (TH1F*)GetRawsData(kOccupancyVsSector) ;
147       TH1F * histRawsNClustersPerEventVsSector = (TH1F*)GetRawsData(kNClustersPerEventVsSector) ;
148       TH1F * histRawsQVsSector                 = (TH1F*)GetRawsData(kQVsSector) ;
149       TH1F * histRawsQmaxVsSector              = (TH1F*)GetRawsData(kQmaxVsSector) ;
150       if ( !histRawsOccupancy ||
151           !histRawsOccupancyVsSector ||
152           !histRawsNClustersPerEventVsSector ||
153           !histRawsQVsSector ||
154           !histRawsQmaxVsSector) {
155         AliError("Something very wrong here, corrupted memory ?????. Please check\n") ; 
156         continue ; 
157       }
158         
159       //Add2RawsList(fTPCdataQA, 0);
160       // get the histograms and add them to the output
161       // 31/8-08 Histogram is only added if the Calibration class 
162       //         receives TPC data 
163       const Int_t eventCounter = fTPCdataQA[specie]->GetEventCounter();
164       if(eventCounter>0) { // some TPC data has been processed
165
166         // Reset histograms and refill them 
167         histRawsOccupancy->Reset();
168         histRawsOccupancyVsSector->Reset();
169         histRawsNClustersPerEventVsSector->Reset();
170         histRawsQVsSector->Reset();
171         histRawsQmaxVsSector->Reset();
172       
173         TH1F* hNorm72 = new TH1F("hNorm72", "histogram to normalize 72 sectors",
174                                  72, 0, 72);
175         hNorm72->Sumw2();
176         TH1F* hNorm108 = new TH1F("hNorm108", "histogram to normalize 108 sectors (medium and long pads are split up)",
177                                   108, 0, 108);
178         hNorm108->Sumw2();
179
180         for (Int_t iSec = 0; iSec < 72; iSec++) {
181         
182           AliTPCCalROC* occupancyROC = 
183           fTPCdataQA[specie]->GetNoThreshold()->GetCalROC(iSec); 
184           AliTPCCalROC* nclusterROC = 
185           fTPCdataQA[specie]->GetNLocalMaxima()->GetCalROC(iSec); 
186           AliTPCCalROC* qROC = 
187           fTPCdataQA[specie]->GetMeanCharge()->GetCalROC(iSec); 
188           AliTPCCalROC* qmaxROC = 
189           fTPCdataQA[specie]->GetMaxCharge()->GetCalROC(iSec); 
190
191           const Int_t nRows = occupancyROC->GetNrows(); 
192           for (Int_t iRow = 0; iRow < nRows; iRow++) {
193
194             Int_t helpSector = iSec;
195             if(iRow>=64)
196               helpSector += 36; // OROC (long pads)
197
198             const Int_t nPads = occupancyROC->GetNPads(iRow); 
199             for (Int_t iPad = 0; iPad < nPads; iPad++) {
200         
201               histRawsOccupancy->Fill(occupancyROC->GetValue(iRow, iPad));
202               hNorm72->Fill(iSec);
203               histRawsOccupancyVsSector
204               ->Fill(iSec, occupancyROC->GetValue(iRow, iPad));
205
206               const Int_t nClusters = TMath::Nint(nclusterROC->GetValue(iRow, iPad));
207             
208               if(nClusters>0) {
209               
210                 histRawsNClustersPerEventVsSector->Fill(iSec, nClusters);
211                 hNorm108->Fill(helpSector, nClusters);
212                 histRawsQVsSector->Fill(helpSector, 
213                                        nClusters*qROC->GetValue(iRow, iPad));
214                 histRawsQmaxVsSector->Fill(helpSector, 
215                                           nClusters*qmaxROC->GetValue(iRow, iPad));
216               }
217             }
218           }
219         } // end loop over sectors
220       
221         // Normalize histograms
222         histRawsOccupancyVsSector->Divide(hNorm72);
223         histRawsNClustersPerEventVsSector->Scale(1.0/Float_t(eventCounter));
224         histRawsQVsSector->Divide(hNorm108);
225         histRawsQmaxVsSector->Divide(hNorm108);
226         delete hNorm72;
227         delete hNorm108;
228       }
229     }
230   }
231   AliQAChecker::Instance()->Run(AliQAv1::kTPC, task, list) ;  
232 }
233
234 //____________________________________________________________________________ 
235 void AliTPCQADataMakerRec::InitESDs()
236 {
237   //create ESDs histograms in ESDs subdir  
238   TH1F * histESDclusters = 
239     new TH1F("hESDclusters", "N TPC clusters per track; N clusters; Counts",
240              160, 0, 160);
241   histESDclusters->Sumw2();
242   Add2ESDsList(histESDclusters, KClusters);
243
244   TH1F * histESDratio = 
245     new TH1F("hESDratio", "Ratio: TPC clusters / findable; Ratio: cluster/findable; Counts",
246              100, 0, 1);
247   histESDratio->Sumw2();
248   Add2ESDsList(histESDratio, kRatio);
249   
250   TH1F * histESDpt = 
251     new TH1F("hESDpt", "P_{T} distribution; p_{T} [GeV/c]; Counts",
252              50, 0, 5);
253   histESDpt->Sumw2();
254   Add2ESDsList(histESDpt, kPt);
255 }
256
257 //____________________________________________________________________________ 
258 void AliTPCQADataMakerRec::InitRaws()
259 {
260   //
261   // Adding the raw 
262   //  
263
264   // Modified: 7/7 - 2008
265   // Laurent Aphecetche pointed out that the mapping was read from file
266   // for each event, so now we read in the map here and set if for 
267   // the raw data qa
268   for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
269     fTPCdataQA[specie] = new AliTPCdataQA(AliRecoParam::Convert(specie));
270     LoadMaps(); // Load Altro maps
271     fTPCdataQA[specie]->SetAltroMapping(fMapping); // set Altro mapping
272     fTPCdataQA[specie]->SetRangeTime(100, 920); // set time bin interval 
273 //    Add2RawsList(fTPCdataQA, kTPCdataQA); // This is used by the AMORE monitoring <------- THIS WILL FAIL (YS)
274   }
275
276   TH1F * histRawsOccupancy = 
277     new TH1F("hRawsOccupancy", "Occupancy (all pads); Occupancy; Counts",
278              100, 0, 1);
279   histRawsOccupancy->Sumw2();
280   Add2RawsList(histRawsOccupancy, kOccupancy);
281   
282   TH1F * histRawsOccupancyVsSector = 
283     new TH1F("hRawsOccupancyVsSector", "Occupancy vs sector; Sector; Occupancy",
284              72, 0, 72);
285   histRawsOccupancyVsSector->Sumw2();
286   Add2RawsList(histRawsOccupancyVsSector, kOccupancyVsSector);
287
288   TH1F * histRawsNClustersPerEventVsSector = 
289     new TH1F("hRawsNClustersPerEventVsSector", "Nclusters per event vs sector; Sector; Nclusters per event",
290              72, 0, 72);
291   histRawsNClustersPerEventVsSector->Sumw2();
292   Add2RawsList(histRawsNClustersPerEventVsSector, kNClustersPerEventVsSector);
293   
294   TH1F * histRawsQVsSector = 
295     new TH1F("hRawsQVsSector", "<Q> vs sector (OROC med: 36-71, long: 72-107); Sector; <Q>",
296              108, 0, 108);
297   histRawsQVsSector->Sumw2();
298   Add2RawsList(histRawsQVsSector, kQVsSector);
299
300   TH1F * histRawsQmaxVsSector = 
301     new TH1F("hRawsQmaxVsSector", "<Qmax> vs sector (OROC med: 36-71, long: 72-107); Sector; <Qmax>",
302              108, 0, 108);
303   histRawsQmaxVsSector->Sumw2();
304   Add2RawsList(histRawsQmaxVsSector, kQmaxVsSector);
305 }
306
307 //____________________________________________________________________________ 
308 void AliTPCQADataMakerRec::InitRecPoints()
309 {
310   TH1F * histRecPointsQmaxShort = 
311     new TH1F("hRecPointsQmaxShort", "Qmax distrbution (short pads); Qmax; Counts",
312              100, 0, 300);
313   histRecPointsQmaxShort->Sumw2();
314   Add2RecPointsList(histRecPointsQmaxShort, kQmaxShort);
315
316   TH1F * histRecPointsQmaxMedium = 
317     new TH1F("hRecPointsQmaxMedium", "Qmax distrbution (medium pads); Qmax; Counts",
318              100, 0, 300);
319   histRecPointsQmaxMedium->Sumw2();
320   Add2RecPointsList(histRecPointsQmaxMedium, kQmaxMedium);
321
322   TH1F * histRecPointsQmaxLong = 
323     new TH1F("hRecPointsQmaxLong", "Qmax distrbution (long pads); Qmax; Counts",
324              100, 0, 300);
325   histRecPointsQmaxLong->Sumw2();
326   Add2RecPointsList(histRecPointsQmaxLong, kQmaxLong);
327
328   TH1F * histRecPointsQShort = 
329     new TH1F("hRecPointsQShort", "Q distrbution (short pads); Q; Counts",
330              100, 0, 2000);
331   histRecPointsQShort->Sumw2();
332   Add2RecPointsList(histRecPointsQShort, kQShort);
333
334   TH1F * histRecPointsQMedium = 
335     new TH1F("hRecPointsQMedium", "Q distrbution (medium pads); Q; Counts",
336              100, 0, 2000);
337   histRecPointsQMedium->Sumw2();
338   Add2RecPointsList(histRecPointsQMedium, kQMedium);
339
340   TH1F * histRecPointsQLong = 
341     new TH1F("hRecPointsQLong", "Q distrbution (long pads); Q; Counts",
342              100, 0, 2000);
343   histRecPointsQLong->Sumw2();
344   Add2RecPointsList(histRecPointsQLong, kQLong);
345
346   TH1F * histRecPointsRow = 
347     new TH1F("hRecPointsRow", "Clusters per row; Row; Counts",
348              159, 0, 159);
349   histRecPointsRow->Sumw2();
350   Add2RecPointsList(histRecPointsRow, kRow);
351 }
352
353 //____________________________________________________________________________
354 void AliTPCQADataMakerRec::MakeESDs(AliESDEvent * esd)
355 {
356   // make QA data from ESDs
357   
358   const Int_t nESDTracks = esd->GetNumberOfTracks();
359   Int_t nTPCtracks = 0; 
360   for(Int_t i = 0; i < nESDTracks; i++) {
361     
362     AliESDtrack * track = esd->GetTrack(i);
363     
364     if ((track->GetStatus() & AliESDtrack::kTPCrefit)==0)
365       continue;
366     
367     nTPCtracks++;
368     
369     Int_t nTPCclusters         = track->GetTPCNcls();
370     Int_t nTPCclustersFindable = track->GetTPCNclsF();
371     if ( nTPCclustersFindable<=0) continue;
372     GetESDsData(KClusters)->Fill(nTPCclusters);
373     GetESDsData(kRatio)->Fill(Float_t(nTPCclusters)/Float_t(nTPCclustersFindable));
374     GetESDsData(kPt)->Fill(track->Pt()); 
375   }
376 }
377
378 //____________________________________________________________________________
379 void AliTPCQADataMakerRec::MakeRaws(AliRawReader* rawReader)
380 {
381   //
382   // To make QA for the RAW data we use the TPC Calibration framework 
383   // to handle the data and then in the end extract the data
384   //
385         rawReader->Reset() ; 
386   fTPCdataQA[AliRecoParam::AConvert(fEventSpecie)]->ProcessEvent(rawReader);  
387  }
388
389 //____________________________________________________________________________
390 void AliTPCQADataMakerRec::MakeRecPoints(TTree* recTree)
391 {
392   AliTPCClustersRow *clrow = new AliTPCClustersRow();
393   clrow->SetClass("AliTPCclusterMI");
394   clrow->SetArray(0);
395   clrow->GetArray()->ExpandCreateFast(10000);
396
397   TBranch* branch = recTree->GetBranch("Segment");
398   branch->SetAddress(&clrow);
399
400   const Int_t nEntries = Int_t(recTree->GetEntries());
401   for (Int_t i = 0; i < nEntries; i++) {
402     
403     branch->GetEntry(i);
404     
405     const Int_t nClusters = clrow->GetArray()->GetEntriesFast();
406     for (Int_t icl=0; icl < nClusters; icl++){
407       
408       AliTPCclusterMI* cluster = 
409         (AliTPCclusterMI*)clrow->GetArray()->At(icl);
410       
411       Float_t Qmax = cluster->GetMax();
412       Float_t Q    = cluster->GetQ();
413       Int_t   row  = cluster->GetRow();
414
415       if(cluster->GetDetector()<36) { // IROC (short pads)
416
417         GetRecPointsData(kQmaxShort)->Fill(Qmax);
418         GetRecPointsData(kQShort)->Fill(Q);
419       } else { // OROC (medium and long pads)
420         row += 63;
421         if(cluster->GetRow()<64) { // medium pads
422
423           GetRecPointsData(kQmaxMedium)->Fill(Qmax);
424           GetRecPointsData(kQMedium)->Fill(Q);
425         } else { // long pads
426
427           GetRecPointsData(kQmaxLong)->Fill(Qmax);
428           GetRecPointsData(kQLong)->Fill(Q);
429         }
430       }
431       
432       GetRecPointsData(kRow)->Fill(row);
433     } // end loop over clusters
434   } // end loop over tree
435
436   delete clrow;
437 }
438
439 //____________________________________________________________________________
440 void AliTPCQADataMakerRec::LoadMaps()
441 {
442   TString path = gSystem->Getenv("ALICE_ROOT");
443   path += "/TPC/mapping/Patch";
444
445   for(Int_t i = 0; i < 6; i++) {
446     TString path2 = path;
447     path2 += i;
448     path2 += ".data";
449     fMapping[i] = new AliTPCAltroMapping(path2.Data());
450   }
451 }
452