Added QA for digits during reconstruction (Yves)
[u/mrichter/AliRoot.git] / TPC / AliTPCQADataMakerRec.cxx
CommitLineData
44f32dd2 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>
be4f1702 52#include <TString.h>
53#include <TSystem.h>
44f32dd2 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"
adcdb395 63#include "AliTPCCalROC.h"
44f32dd2 64#include "AliTPCClustersRow.h"
65#include "AliTPCclusterMI.h"
44ed7a66 66#include "AliSimDigits.h"
44f32dd2 67
68ClassImp(AliTPCQADataMakerRec)
69
70//____________________________________________________________________________
71AliTPCQADataMakerRec::AliTPCQADataMakerRec() :
4e25ac79 72 AliQADataMakerRec(AliQAv1::GetDetName(AliQAv1::kTPC),
44f32dd2 73 "TPC Rec Quality Assurance Data Maker"),
57acd2d2 74 fTPCdataQA(NULL)
44f32dd2 75{
76 // ctor
57acd2d2 77 fTPCdataQA = new AliTPCdataQA*[AliRecoParam::kNSpecies] ;
78 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++)
79 fTPCdataQA[specie] = NULL ;
80
be4f1702 81 for(Int_t i = 0; i < 6; i++)
82 fMapping[i] = 0;
44f32dd2 83}
84
85//____________________________________________________________________________
86AliTPCQADataMakerRec::AliTPCQADataMakerRec(const AliTPCQADataMakerRec& qadm) :
336156cc 87 AliQADataMakerRec(),
57acd2d2 88 fTPCdataQA(NULL)
44f32dd2 89{
90 //copy ctor
91 // Does not copy the calibration object, instead InitRaws have to be
92 // called again
93 SetName((const char*)qadm.GetName()) ;
94 SetTitle((const char*)qadm.GetTitle());
95
57acd2d2 96 fTPCdataQA = new AliTPCdataQA*[AliRecoParam::kNSpecies] ;
97 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++)
98 fTPCdataQA[specie] = NULL ;
99
be4f1702 100 for(Int_t i = 0; i < 6; i++)
101 fMapping[i] = 0;
102
44f32dd2 103 //
104 // Associate class histogram objects to the copies in the list
105 // Could also be done with the indexes
106 //
57acd2d2 107
44f32dd2 108}
109
110//__________________________________________________________________
111AliTPCQADataMakerRec& AliTPCQADataMakerRec::operator = (const AliTPCQADataMakerRec& qadm )
112{
113 // Equal operator.
114 this->~AliTPCQADataMakerRec();
115 new(this) AliTPCQADataMakerRec(qadm);
116 return *this;
117}
be4f1702 118
119//__________________________________________________________________
120AliTPCQADataMakerRec::~AliTPCQADataMakerRec()
121{
122 // Destructor
57acd2d2 123 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++)
124 if ( fTPCdataQA[specie] != NULL )
125 delete fTPCdataQA[specie] ;
126 delete[] fTPCdataQA;
be4f1702 127
128 for(Int_t i = 0; i < 6; i++)
129 delete fMapping[i];
130}
44f32dd2 131
132//____________________________________________________________________________
4e25ac79 133void AliTPCQADataMakerRec::EndOfDetectorCycle(AliQAv1::TASKINDEX_t task, TObjArray ** list)
44f32dd2 134{
135 //Detector specific actions at end of cycle
57acd2d2 136
137 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
4e25ac79 138 if ( !AliQAv1::Instance()->IsEventSpecieSet(specie) )
57acd2d2 139 continue ;
140 if(fTPCdataQA[specie] != NULL) { // do the final step of the QA for Raw data
44f32dd2 141
57acd2d2 142 fTPCdataQA[specie]->Analyse(); // 31/1-08 Analyse is now protected against
c75bf2f1 143 // RAW data files with no TPC data
57acd2d2 144
5e232cd6 145 SetEventSpecie(AliRecoParam::ConvertIndex(specie)) ;
57acd2d2 146 TH1F * histRawsOccupancy = (TH1F*)GetRawsData(kOccupancy) ;
147 TH1F * histRawsOccupancyVsSector = (TH1F*)GetRawsData(kOccupancyVsSector) ;
148 TH1F * histRawsNClustersPerEventVsSector = (TH1F*)GetRawsData(kNClustersPerEventVsSector) ;
149 TH1F * histRawsQVsSector = (TH1F*)GetRawsData(kQVsSector) ;
150 TH1F * histRawsQmaxVsSector = (TH1F*)GetRawsData(kQmaxVsSector) ;
151 if ( !histRawsOccupancy ||
152 !histRawsOccupancyVsSector ||
153 !histRawsNClustersPerEventVsSector ||
154 !histRawsQVsSector ||
155 !histRawsQmaxVsSector) {
156 AliError("Something very wrong here, corrupted memory ?????. Please check\n") ;
157 continue ;
158 }
159
160 //Add2RawsList(fTPCdataQA, 0);
161 // get the histograms and add them to the output
162 // 31/8-08 Histogram is only added if the Calibration class
163 // receives TPC data
164 const Int_t eventCounter = fTPCdataQA[specie]->GetEventCounter();
165 if(eventCounter>0) { // some TPC data has been processed
166
167 // Reset histograms and refill them
168 histRawsOccupancy->Reset();
169 histRawsOccupancyVsSector->Reset();
170 histRawsNClustersPerEventVsSector->Reset();
171 histRawsQVsSector->Reset();
172 histRawsQmaxVsSector->Reset();
adcdb395 173
57acd2d2 174 TH1F* hNorm72 = new TH1F("hNorm72", "histogram to normalize 72 sectors",
175 72, 0, 72);
176 hNorm72->Sumw2();
177 TH1F* hNorm108 = new TH1F("hNorm108", "histogram to normalize 108 sectors (medium and long pads are split up)",
178 108, 0, 108);
179 hNorm108->Sumw2();
180
181 for (Int_t iSec = 0; iSec < 72; iSec++) {
adcdb395 182
57acd2d2 183 AliTPCCalROC* occupancyROC =
184 fTPCdataQA[specie]->GetNoThreshold()->GetCalROC(iSec);
185 AliTPCCalROC* nclusterROC =
186 fTPCdataQA[specie]->GetNLocalMaxima()->GetCalROC(iSec);
187 AliTPCCalROC* qROC =
188 fTPCdataQA[specie]->GetMeanCharge()->GetCalROC(iSec);
189 AliTPCCalROC* qmaxROC =
190 fTPCdataQA[specie]->GetMaxCharge()->GetCalROC(iSec);
191
192 const Int_t nRows = occupancyROC->GetNrows();
193 for (Int_t iRow = 0; iRow < nRows; iRow++) {
194
195 Int_t helpSector = iSec;
196 if(iRow>=64)
197 helpSector += 36; // OROC (long pads)
198
199 const Int_t nPads = occupancyROC->GetNPads(iRow);
200 for (Int_t iPad = 0; iPad < nPads; iPad++) {
adcdb395 201
57acd2d2 202 histRawsOccupancy->Fill(occupancyROC->GetValue(iRow, iPad));
203 hNorm72->Fill(iSec);
204 histRawsOccupancyVsSector
205 ->Fill(iSec, occupancyROC->GetValue(iRow, iPad));
adcdb395 206
57acd2d2 207 const Int_t nClusters = TMath::Nint(nclusterROC->GetValue(iRow, iPad));
adcdb395 208
57acd2d2 209 if(nClusters>0) {
adcdb395 210
57acd2d2 211 histRawsNClustersPerEventVsSector->Fill(iSec, nClusters);
212 hNorm108->Fill(helpSector, nClusters);
213 histRawsQVsSector->Fill(helpSector,
214 nClusters*qROC->GetValue(iRow, iPad));
215 histRawsQmaxVsSector->Fill(helpSector,
216 nClusters*qmaxROC->GetValue(iRow, iPad));
217 }
218 }
219 }
220 } // end loop over sectors
adcdb395 221
57acd2d2 222 // Normalize histograms
223 histRawsOccupancyVsSector->Divide(hNorm72);
224 histRawsNClustersPerEventVsSector->Scale(1.0/Float_t(eventCounter));
225 histRawsQVsSector->Divide(hNorm108);
226 histRawsQmaxVsSector->Divide(hNorm108);
227 delete hNorm72;
228 delete hNorm108;
229 }
c75bf2f1 230 }
44f32dd2 231 }
4e25ac79 232 AliQAChecker::Instance()->Run(AliQAv1::kTPC, task, list) ;
44f32dd2 233}
234
235//____________________________________________________________________________
236void AliTPCQADataMakerRec::InitESDs()
237{
238 //create ESDs histograms in ESDs subdir
7d297381 239 const Bool_t expert = kTRUE ;
240 const Bool_t image = kTRUE ;
241
57acd2d2 242 TH1F * histESDclusters =
44f32dd2 243 new TH1F("hESDclusters", "N TPC clusters per track; N clusters; Counts",
244 160, 0, 160);
57acd2d2 245 histESDclusters->Sumw2();
7d297381 246 Add2ESDsList(histESDclusters, KClusters, !expert, image);
44f32dd2 247
57acd2d2 248 TH1F * histESDratio =
44f32dd2 249 new TH1F("hESDratio", "Ratio: TPC clusters / findable; Ratio: cluster/findable; Counts",
250 100, 0, 1);
57acd2d2 251 histESDratio->Sumw2();
7d297381 252 Add2ESDsList(histESDratio, kRatio, !expert, image);
44f32dd2 253
57acd2d2 254 TH1F * histESDpt =
44f32dd2 255 new TH1F("hESDpt", "P_{T} distribution; p_{T} [GeV/c]; Counts",
256 50, 0, 5);
57acd2d2 257 histESDpt->Sumw2();
7d297381 258 Add2ESDsList(histESDpt, kPt, !expert, image);
44f32dd2 259}
260
261//____________________________________________________________________________
262void AliTPCQADataMakerRec::InitRaws()
263{
336156cc 264 //
265 // Adding the raw
be4f1702 266 //
267
268 // Modified: 7/7 - 2008
269 // Laurent Aphecetche pointed out that the mapping was read from file
270 // for each event, so now we read in the map here and set if for
271 // the raw data qa
7d297381 272 const Bool_t expert = kTRUE ;
273 const Bool_t saveCorr = kTRUE ;
274 const Bool_t image = kTRUE ;
275
57acd2d2 276 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
5e232cd6 277 fTPCdataQA[specie] = new AliTPCdataQA(AliRecoParam::ConvertIndex(specie));
57acd2d2 278 LoadMaps(); // Load Altro maps
279 fTPCdataQA[specie]->SetAltroMapping(fMapping); // set Altro mapping
280 fTPCdataQA[specie]->SetRangeTime(100, 920); // set time bin interval
7d297381 281// Add2RawsList(fTPCdataQA, kTPCdataQ, !expert, image, !saveCorrA); // This is used by the AMORE monitoring <------- THIS WILL FAIL (YS)
57acd2d2 282 }
be4f1702 283
57acd2d2 284 TH1F * histRawsOccupancy =
adcdb395 285 new TH1F("hRawsOccupancy", "Occupancy (all pads); Occupancy; Counts",
286 100, 0, 1);
57acd2d2 287 histRawsOccupancy->Sumw2();
7d297381 288 Add2RawsList(histRawsOccupancy, kOccupancy, !expert, image, !saveCorr);
57acd2d2 289
290 TH1F * histRawsOccupancyVsSector =
adcdb395 291 new TH1F("hRawsOccupancyVsSector", "Occupancy vs sector; Sector; Occupancy",
292 72, 0, 72);
57acd2d2 293 histRawsOccupancyVsSector->Sumw2();
7d297381 294 Add2RawsList(histRawsOccupancyVsSector, kOccupancyVsSector, !expert, image, !saveCorr);
adcdb395 295
57acd2d2 296 TH1F * histRawsNClustersPerEventVsSector =
adcdb395 297 new TH1F("hRawsNClustersPerEventVsSector", "Nclusters per event vs sector; Sector; Nclusters per event",
298 72, 0, 72);
57acd2d2 299 histRawsNClustersPerEventVsSector->Sumw2();
7d297381 300 Add2RawsList(histRawsNClustersPerEventVsSector, kNClustersPerEventVsSector, !expert, image, !saveCorr);
adcdb395 301
57acd2d2 302 TH1F * histRawsQVsSector =
adcdb395 303 new TH1F("hRawsQVsSector", "<Q> vs sector (OROC med: 36-71, long: 72-107); Sector; <Q>",
304 108, 0, 108);
57acd2d2 305 histRawsQVsSector->Sumw2();
7d297381 306 Add2RawsList(histRawsQVsSector, kQVsSector, !expert, image, !saveCorr);
adcdb395 307
57acd2d2 308 TH1F * histRawsQmaxVsSector =
adcdb395 309 new TH1F("hRawsQmaxVsSector", "<Qmax> vs sector (OROC med: 36-71, long: 72-107); Sector; <Qmax>",
310 108, 0, 108);
57acd2d2 311 histRawsQmaxVsSector->Sumw2();
7d297381 312 Add2RawsList(histRawsQmaxVsSector, kQmaxVsSector, !expert, image, !saveCorr);
44f32dd2 313}
314
315//____________________________________________________________________________
44ed7a66 316void AliTPCQADataMakerRec::InitDigits()
317{
318 const Bool_t expert = kTRUE ;
319 const Bool_t image = kTRUE ;
320 TH1F * histDigitsADC =
321 new TH1F("hDigitsADC", "Digit ADC distribution; ADC; Counts",
322 1000, 0, 1000);
323 histDigitsADC->Sumw2();
324 Add2DigitsList(histDigitsADC, kDigitsADC, !expert, image);
325}
326
327//____________________________________________________________________________
44f32dd2 328void AliTPCQADataMakerRec::InitRecPoints()
329{
7d297381 330 const Bool_t expert = kTRUE ;
331 const Bool_t image = kTRUE ;
332
57acd2d2 333 TH1F * histRecPointsQmaxShort =
44f32dd2 334 new TH1F("hRecPointsQmaxShort", "Qmax distrbution (short pads); Qmax; Counts",
adcdb395 335 100, 0, 300);
57acd2d2 336 histRecPointsQmaxShort->Sumw2();
7d297381 337 Add2RecPointsList(histRecPointsQmaxShort, kQmaxShort, !expert, image);
44f32dd2 338
57acd2d2 339 TH1F * histRecPointsQmaxMedium =
44f32dd2 340 new TH1F("hRecPointsQmaxMedium", "Qmax distrbution (medium pads); Qmax; Counts",
adcdb395 341 100, 0, 300);
57acd2d2 342 histRecPointsQmaxMedium->Sumw2();
7d297381 343 Add2RecPointsList(histRecPointsQmaxMedium, kQmaxMedium, !expert, image);
44f32dd2 344
57acd2d2 345 TH1F * histRecPointsQmaxLong =
44f32dd2 346 new TH1F("hRecPointsQmaxLong", "Qmax distrbution (long pads); Qmax; Counts",
adcdb395 347 100, 0, 300);
57acd2d2 348 histRecPointsQmaxLong->Sumw2();
7d297381 349 Add2RecPointsList(histRecPointsQmaxLong, kQmaxLong, !expert, image);
44f32dd2 350
57acd2d2 351 TH1F * histRecPointsQShort =
44f32dd2 352 new TH1F("hRecPointsQShort", "Q distrbution (short pads); Q; Counts",
adcdb395 353 100, 0, 2000);
57acd2d2 354 histRecPointsQShort->Sumw2();
7d297381 355 Add2RecPointsList(histRecPointsQShort, kQShort, !expert, image);
44f32dd2 356
57acd2d2 357 TH1F * histRecPointsQMedium =
44f32dd2 358 new TH1F("hRecPointsQMedium", "Q distrbution (medium pads); Q; Counts",
adcdb395 359 100, 0, 2000);
57acd2d2 360 histRecPointsQMedium->Sumw2();
7d297381 361 Add2RecPointsList(histRecPointsQMedium, kQMedium, !expert, image);
44f32dd2 362
57acd2d2 363 TH1F * histRecPointsQLong =
44f32dd2 364 new TH1F("hRecPointsQLong", "Q distrbution (long pads); Q; Counts",
adcdb395 365 100, 0, 2000);
57acd2d2 366 histRecPointsQLong->Sumw2();
7d297381 367 Add2RecPointsList(histRecPointsQLong, kQLong, !expert, image);
44f32dd2 368
57acd2d2 369 TH1F * histRecPointsRow =
44f32dd2 370 new TH1F("hRecPointsRow", "Clusters per row; Row; Counts",
371 159, 0, 159);
57acd2d2 372 histRecPointsRow->Sumw2();
7d297381 373 Add2RecPointsList(histRecPointsRow, kRow, !expert, image);
44f32dd2 374}
375
376//____________________________________________________________________________
377void AliTPCQADataMakerRec::MakeESDs(AliESDEvent * esd)
378{
379 // make QA data from ESDs
380
381 const Int_t nESDTracks = esd->GetNumberOfTracks();
382 Int_t nTPCtracks = 0;
383 for(Int_t i = 0; i < nESDTracks; i++) {
384
385 AliESDtrack * track = esd->GetTrack(i);
386
387 if ((track->GetStatus() & AliESDtrack::kTPCrefit)==0)
388 continue;
389
390 nTPCtracks++;
391
392 Int_t nTPCclusters = track->GetTPCNcls();
393 Int_t nTPCclustersFindable = track->GetTPCNclsF();
47a596d8 394 if ( nTPCclustersFindable<=0) continue;
57acd2d2 395 GetESDsData(KClusters)->Fill(nTPCclusters);
396 GetESDsData(kRatio)->Fill(Float_t(nTPCclusters)/Float_t(nTPCclustersFindable));
397 GetESDsData(kPt)->Fill(track->Pt());
44f32dd2 398 }
399}
400
401//____________________________________________________________________________
402void AliTPCQADataMakerRec::MakeRaws(AliRawReader* rawReader)
403{
404 //
405 // To make QA for the RAW data we use the TPC Calibration framework
406 // to handle the data and then in the end extract the data
407 //
78328afd 408 rawReader->Reset() ;
57acd2d2 409 fTPCdataQA[AliRecoParam::AConvert(fEventSpecie)]->ProcessEvent(rawReader);
410 }
44f32dd2 411
412//____________________________________________________________________________
44ed7a66 413void AliTPCQADataMakerRec::MakeDigits(TTree* digitTree)
414{
415 TBranch* branch = digitTree->GetBranch("Segment");
416 AliSimDigits* digArray = 0;
417 branch->SetAddress(&digArray);
418
419 Int_t nEntries = Int_t(digitTree->GetEntries());
420
421 for (Int_t n = 0; n < nEntries; n++) {
422
423 digitTree->GetEvent(n);
424
425 if (digArray->First())
426 do {
427 Float_t dig = digArray->CurrentDigit();
428
429 GetDigitsData(kDigitsADC)->Fill(dig);
430 } while (digArray->Next());
431 }
432}
433
434//____________________________________________________________________________
44f32dd2 435void AliTPCQADataMakerRec::MakeRecPoints(TTree* recTree)
436{
437 AliTPCClustersRow *clrow = new AliTPCClustersRow();
438 clrow->SetClass("AliTPCclusterMI");
439 clrow->SetArray(0);
440 clrow->GetArray()->ExpandCreateFast(10000);
441
442 TBranch* branch = recTree->GetBranch("Segment");
443 branch->SetAddress(&clrow);
444
445 const Int_t nEntries = Int_t(recTree->GetEntries());
446 for (Int_t i = 0; i < nEntries; i++) {
447
448 branch->GetEntry(i);
449
450 const Int_t nClusters = clrow->GetArray()->GetEntriesFast();
451 for (Int_t icl=0; icl < nClusters; icl++){
452
453 AliTPCclusterMI* cluster =
454 (AliTPCclusterMI*)clrow->GetArray()->At(icl);
455
456 Float_t Qmax = cluster->GetMax();
457 Float_t Q = cluster->GetQ();
458 Int_t row = cluster->GetRow();
459
460 if(cluster->GetDetector()<36) { // IROC (short pads)
461
57acd2d2 462 GetRecPointsData(kQmaxShort)->Fill(Qmax);
463 GetRecPointsData(kQShort)->Fill(Q);
44f32dd2 464 } else { // OROC (medium and long pads)
465 row += 63;
466 if(cluster->GetRow()<64) { // medium pads
467
57acd2d2 468 GetRecPointsData(kQmaxMedium)->Fill(Qmax);
469 GetRecPointsData(kQMedium)->Fill(Q);
44f32dd2 470 } else { // long pads
471
57acd2d2 472 GetRecPointsData(kQmaxLong)->Fill(Qmax);
473 GetRecPointsData(kQLong)->Fill(Q);
44f32dd2 474 }
475 }
476
57acd2d2 477 GetRecPointsData(kRow)->Fill(row);
44f32dd2 478 } // end loop over clusters
479 } // end loop over tree
480
481 delete clrow;
482}
be4f1702 483
484//____________________________________________________________________________
485void AliTPCQADataMakerRec::LoadMaps()
486{
487 TString path = gSystem->Getenv("ALICE_ROOT");
488 path += "/TPC/mapping/Patch";
489
490 for(Int_t i = 0; i < 6; i++) {
491 TString path2 = path;
492 path2 += i;
493 path2 += ".data";
494 fMapping[i] = new AliTPCAltroMapping(path2.Data());
495 }
496}
497