1 /**************************************************************************
2 * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
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
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.
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.
44 Sumw2 for RAW histogram(s)?
45 RecPoints and ESD could have many more histograms
48 #include "AliTPCQADataMakerRec.h"
50 // --- ROOT system ---
51 #include <TClonesArray.h>
55 // --- Standard library ---
57 // --- AliRoot header files ---
58 #include "AliQAChecker.h"
59 #include "AliESDEvent.h"
60 #include "AliESDtrack.h"
62 #include "AliTPCCalPad.h"
63 #include "AliTPCCalROC.h"
64 #include "AliTPCClustersRow.h"
65 #include "AliTPCclusterMI.h"
66 #include "AliSimDigits.h"
68 ClassImp(AliTPCQADataMakerRec)
70 //____________________________________________________________________________
71 AliTPCQADataMakerRec::AliTPCQADataMakerRec() :
72 AliQADataMakerRec(AliQAv1::GetDetName(AliQAv1::kTPC),
73 "TPC Rec Quality Assurance Data Maker"),
77 fTPCdataQA = new AliTPCdataQA*[AliRecoParam::kNSpecies] ;
78 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++)
79 fTPCdataQA[specie] = NULL ;
81 for(Int_t i = 0; i < 6; i++)
85 //____________________________________________________________________________
86 AliTPCQADataMakerRec::AliTPCQADataMakerRec(const AliTPCQADataMakerRec& qadm) :
91 // Does not copy the calibration object, instead InitRaws have to be
93 SetName((const char*)qadm.GetName()) ;
94 SetTitle((const char*)qadm.GetTitle());
96 fTPCdataQA = new AliTPCdataQA*[AliRecoParam::kNSpecies] ;
97 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++)
98 fTPCdataQA[specie] = NULL ;
100 for(Int_t i = 0; i < 6; i++)
104 // Associate class histogram objects to the copies in the list
105 // Could also be done with the indexes
110 //__________________________________________________________________
111 AliTPCQADataMakerRec& AliTPCQADataMakerRec::operator = (const AliTPCQADataMakerRec& qadm )
114 this->~AliTPCQADataMakerRec();
115 new(this) AliTPCQADataMakerRec(qadm);
119 //__________________________________________________________________
120 AliTPCQADataMakerRec::~AliTPCQADataMakerRec()
123 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++)
124 if ( fTPCdataQA[specie] != NULL )
125 delete fTPCdataQA[specie] ;
128 for(Int_t i = 0; i < 6; i++)
132 //____________________________________________________________________________
133 void AliTPCQADataMakerRec::EndOfDetectorCycle(AliQAv1::TASKINDEX_t task, TObjArray ** list)
135 //Detector specific actions at end of cycle
137 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
138 if ( !AliQAv1::Instance()->IsEventSpecieSet(specie) )
140 if(fTPCdataQA[specie] != NULL) { // do the final step of the QA for Raw data
142 fTPCdataQA[specie]->Analyse(); // 31/1-08 Analyse is now protected against
143 // RAW data files with no TPC data
145 SetEventSpecie(AliRecoParam::ConvertIndex(specie)) ;
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") ;
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
164 const Int_t eventCounter = fTPCdataQA[specie]->GetEventCounter();
165 if(eventCounter>0) { // some TPC data has been processed
167 // Reset histograms and refill them
168 histRawsOccupancy->Reset();
169 histRawsOccupancyVsSector->Reset();
170 histRawsNClustersPerEventVsSector->Reset();
171 histRawsQVsSector->Reset();
172 histRawsQmaxVsSector->Reset();
174 TH1F* hNorm72 = new TH1F("hNorm72", "histogram to normalize 72 sectors",
177 TH1F* hNorm108 = new TH1F("hNorm108", "histogram to normalize 108 sectors (medium and long pads are split up)",
181 for (Int_t iSec = 0; iSec < 72; iSec++) {
183 AliTPCCalROC* occupancyROC =
184 fTPCdataQA[specie]->GetNoThreshold()->GetCalROC(iSec);
185 AliTPCCalROC* nclusterROC =
186 fTPCdataQA[specie]->GetNLocalMaxima()->GetCalROC(iSec);
188 fTPCdataQA[specie]->GetMeanCharge()->GetCalROC(iSec);
189 AliTPCCalROC* qmaxROC =
190 fTPCdataQA[specie]->GetMaxCharge()->GetCalROC(iSec);
192 const Int_t nRows = occupancyROC->GetNrows();
193 for (Int_t iRow = 0; iRow < nRows; iRow++) {
195 Int_t helpSector = iSec;
197 helpSector += 36; // OROC (long pads)
199 const Int_t nPads = occupancyROC->GetNPads(iRow);
200 for (Int_t iPad = 0; iPad < nPads; iPad++) {
202 histRawsOccupancy->Fill(occupancyROC->GetValue(iRow, iPad));
204 histRawsOccupancyVsSector
205 ->Fill(iSec, occupancyROC->GetValue(iRow, iPad));
207 const Int_t nClusters = TMath::Nint(nclusterROC->GetValue(iRow, iPad));
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));
220 } // end loop over sectors
222 // Normalize histograms
223 histRawsOccupancyVsSector->Divide(hNorm72);
224 histRawsNClustersPerEventVsSector->Scale(1.0/Float_t(eventCounter));
225 histRawsQVsSector->Divide(hNorm108);
226 histRawsQmaxVsSector->Divide(hNorm108);
232 AliQAChecker::Instance()->Run(AliQAv1::kTPC, task, list) ;
235 //____________________________________________________________________________
236 void AliTPCQADataMakerRec::InitESDs()
238 //create ESDs histograms in ESDs subdir
239 const Bool_t expert = kTRUE ;
240 const Bool_t image = kTRUE ;
242 TH1F * histESDclusters =
243 new TH1F("hESDclusters", "N TPC clusters per track; N clusters; Counts",
245 histESDclusters->Sumw2();
246 Add2ESDsList(histESDclusters, KClusters, !expert, image);
248 TH1F * histESDratio =
249 new TH1F("hESDratio", "Ratio: TPC clusters / findable; Ratio: cluster/findable; Counts",
251 histESDratio->Sumw2();
252 Add2ESDsList(histESDratio, kRatio, !expert, image);
255 new TH1F("hESDpt", "P_{T} distribution; p_{T} [GeV/c]; Counts",
258 Add2ESDsList(histESDpt, kPt, !expert, image);
261 //____________________________________________________________________________
262 void AliTPCQADataMakerRec::InitRaws()
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
272 const Bool_t expert = kTRUE ;
273 const Bool_t saveCorr = kTRUE ;
274 const Bool_t image = kTRUE ;
276 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
277 fTPCdataQA[specie] = new AliTPCdataQA(AliRecoParam::ConvertIndex(specie));
278 LoadMaps(); // Load Altro maps
279 fTPCdataQA[specie]->SetAltroMapping(fMapping); // set Altro mapping
280 fTPCdataQA[specie]->SetRangeTime(100, 920); // set time bin interval
281 // Add2RawsList(fTPCdataQA, kTPCdataQ, !expert, image, !saveCorrA); // This is used by the AMORE monitoring <------- THIS WILL FAIL (YS)
284 TH1F * histRawsOccupancy =
285 new TH1F("hRawsOccupancy", "Occupancy (all pads); Occupancy; Counts",
287 histRawsOccupancy->Sumw2();
288 Add2RawsList(histRawsOccupancy, kOccupancy, !expert, image, !saveCorr);
290 TH1F * histRawsOccupancyVsSector =
291 new TH1F("hRawsOccupancyVsSector", "Occupancy vs sector; Sector; Occupancy",
293 histRawsOccupancyVsSector->Sumw2();
294 Add2RawsList(histRawsOccupancyVsSector, kOccupancyVsSector, !expert, image, !saveCorr);
296 TH1F * histRawsNClustersPerEventVsSector =
297 new TH1F("hRawsNClustersPerEventVsSector", "Nclusters per event vs sector; Sector; Nclusters per event",
299 histRawsNClustersPerEventVsSector->Sumw2();
300 Add2RawsList(histRawsNClustersPerEventVsSector, kNClustersPerEventVsSector, !expert, image, !saveCorr);
302 TH1F * histRawsQVsSector =
303 new TH1F("hRawsQVsSector", "<Q> vs sector (OROC med: 36-71, long: 72-107); Sector; <Q>",
305 histRawsQVsSector->Sumw2();
306 Add2RawsList(histRawsQVsSector, kQVsSector, !expert, image, !saveCorr);
308 TH1F * histRawsQmaxVsSector =
309 new TH1F("hRawsQmaxVsSector", "<Qmax> vs sector (OROC med: 36-71, long: 72-107); Sector; <Qmax>",
311 histRawsQmaxVsSector->Sumw2();
312 Add2RawsList(histRawsQmaxVsSector, kQmaxVsSector, !expert, image, !saveCorr);
315 //____________________________________________________________________________
316 void AliTPCQADataMakerRec::InitDigits()
318 const Bool_t expert = kTRUE ;
319 const Bool_t image = kTRUE ;
320 TH1F * histDigitsADC =
321 new TH1F("hDigitsADC", "Digit ADC distribution; ADC; Counts",
323 histDigitsADC->Sumw2();
324 Add2DigitsList(histDigitsADC, kDigitsADC, !expert, image);
327 //____________________________________________________________________________
328 void AliTPCQADataMakerRec::InitRecPoints()
330 const Bool_t expert = kTRUE ;
331 const Bool_t image = kTRUE ;
333 TH1F * histRecPointsQmaxShort =
334 new TH1F("hRecPointsQmaxShort", "Qmax distrbution (short pads); Qmax; Counts",
336 histRecPointsQmaxShort->Sumw2();
337 Add2RecPointsList(histRecPointsQmaxShort, kQmaxShort, !expert, image);
339 TH1F * histRecPointsQmaxMedium =
340 new TH1F("hRecPointsQmaxMedium", "Qmax distrbution (medium pads); Qmax; Counts",
342 histRecPointsQmaxMedium->Sumw2();
343 Add2RecPointsList(histRecPointsQmaxMedium, kQmaxMedium, !expert, image);
345 TH1F * histRecPointsQmaxLong =
346 new TH1F("hRecPointsQmaxLong", "Qmax distrbution (long pads); Qmax; Counts",
348 histRecPointsQmaxLong->Sumw2();
349 Add2RecPointsList(histRecPointsQmaxLong, kQmaxLong, !expert, image);
351 TH1F * histRecPointsQShort =
352 new TH1F("hRecPointsQShort", "Q distrbution (short pads); Q; Counts",
354 histRecPointsQShort->Sumw2();
355 Add2RecPointsList(histRecPointsQShort, kQShort, !expert, image);
357 TH1F * histRecPointsQMedium =
358 new TH1F("hRecPointsQMedium", "Q distrbution (medium pads); Q; Counts",
360 histRecPointsQMedium->Sumw2();
361 Add2RecPointsList(histRecPointsQMedium, kQMedium, !expert, image);
363 TH1F * histRecPointsQLong =
364 new TH1F("hRecPointsQLong", "Q distrbution (long pads); Q; Counts",
366 histRecPointsQLong->Sumw2();
367 Add2RecPointsList(histRecPointsQLong, kQLong, !expert, image);
369 TH1F * histRecPointsRow =
370 new TH1F("hRecPointsRow", "Clusters per row; Row; Counts",
372 histRecPointsRow->Sumw2();
373 Add2RecPointsList(histRecPointsRow, kRow, !expert, image);
376 //____________________________________________________________________________
377 void AliTPCQADataMakerRec::MakeESDs(AliESDEvent * esd)
379 // make QA data from ESDs
381 // Check id histograms already created for this Event Specie
382 if ( ! GetESDsData(KClusters) )
385 const Int_t nESDTracks = esd->GetNumberOfTracks();
386 Int_t nTPCtracks = 0;
387 for(Int_t i = 0; i < nESDTracks; i++) {
389 AliESDtrack * track = esd->GetTrack(i);
391 if ((track->GetStatus() & AliESDtrack::kTPCrefit)==0)
396 Int_t nTPCclusters = track->GetTPCNcls();
397 Int_t nTPCclustersFindable = track->GetTPCNclsF();
398 if ( nTPCclustersFindable<=0) continue;
399 GetESDsData(KClusters)->Fill(nTPCclusters);
400 GetESDsData(kRatio)->Fill(Float_t(nTPCclusters)/Float_t(nTPCclustersFindable));
401 GetESDsData(kPt)->Fill(track->Pt());
405 //____________________________________________________________________________
406 void AliTPCQADataMakerRec::MakeRaws(AliRawReader* rawReader)
409 // To make QA for the RAW data we use the TPC Calibration framework
410 // to handle the data and then in the end extract the data
412 if ( ! GetRawsData(kOccupancy) )
416 fTPCdataQA[AliRecoParam::AConvert(fEventSpecie)]->ProcessEvent(rawReader);
419 //____________________________________________________________________________
420 void AliTPCQADataMakerRec::MakeDigits(TTree* digitTree)
423 // Check id histograms already created for this Event Specie
424 if ( ! GetDigitsData(kDigitsADC) )
427 TBranch* branch = digitTree->GetBranch("Segment");
428 AliSimDigits* digArray = 0;
429 branch->SetAddress(&digArray);
431 Int_t nEntries = Int_t(digitTree->GetEntries());
433 for (Int_t n = 0; n < nEntries; n++) {
435 digitTree->GetEvent(n);
437 if (digArray->First())
439 Float_t dig = digArray->CurrentDigit();
441 GetDigitsData(kDigitsADC)->Fill(dig);
442 } while (digArray->Next());
446 //____________________________________________________________________________
447 void AliTPCQADataMakerRec::MakeRecPoints(TTree* recTree)
449 // Check id histograms already created for this Event Specie
450 if ( ! GetRecPointsData(kQmaxShort) )
453 AliTPCClustersRow *clrow = new AliTPCClustersRow();
454 clrow->SetClass("AliTPCclusterMI");
456 clrow->GetArray()->ExpandCreateFast(10000);
458 TBranch* branch = recTree->GetBranch("Segment");
459 branch->SetAddress(&clrow);
461 const Int_t nEntries = Int_t(recTree->GetEntries());
462 for (Int_t i = 0; i < nEntries; i++) {
466 const Int_t nClusters = clrow->GetArray()->GetEntriesFast();
467 for (Int_t icl=0; icl < nClusters; icl++){
469 AliTPCclusterMI* cluster =
470 (AliTPCclusterMI*)clrow->GetArray()->At(icl);
472 Float_t Qmax = cluster->GetMax();
473 Float_t Q = cluster->GetQ();
474 Int_t row = cluster->GetRow();
476 if(cluster->GetDetector()<36) { // IROC (short pads)
478 GetRecPointsData(kQmaxShort)->Fill(Qmax);
479 GetRecPointsData(kQShort)->Fill(Q);
480 } else { // OROC (medium and long pads)
482 if(cluster->GetRow()<64) { // medium pads
484 GetRecPointsData(kQmaxMedium)->Fill(Qmax);
485 GetRecPointsData(kQMedium)->Fill(Q);
486 } else { // long pads
488 GetRecPointsData(kQmaxLong)->Fill(Qmax);
489 GetRecPointsData(kQLong)->Fill(Q);
493 GetRecPointsData(kRow)->Fill(row);
494 } // end loop over clusters
495 } // end loop over tree
500 //____________________________________________________________________________
501 void AliTPCQADataMakerRec::LoadMaps()
503 TString path = gSystem->Getenv("ALICE_ROOT");
504 path += "/TPC/mapping/Patch";
506 for(Int_t i = 0; i < 6; i++) {
507 TString path2 = path;
510 fMapping[i] = new AliTPCAltroMapping(path2.Data());