1 /**************************************************************************
2 * Copyright(c) 1998-1999, 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 **************************************************************************/
17 Produces the data needed to calculate the quality assurance.
18 All data must be mergeable objects.
20 Based on PHOS code written by
21 Y. Schutz CERN July 2007
24 // --- ROOT system ---
25 #include <TClonesArray.h>
31 // --- Standard library ---
33 // --- AliRoot header files ---
34 #include "AliESDCaloCluster.h"
35 #include "AliESDEvent.h"
37 #include "AliEMCALDigit.h"
38 #include "AliEMCALHit.h"
39 #include "AliEMCALQADataMaker.h"
40 #include "AliQAChecker.h"
41 #include "AliEMCALRecPoint.h"
42 #include "AliEMCALRawUtils.h"
43 #include "AliEMCALReconstructor.h"
44 #include "AliEMCALRecParam.h"
46 ClassImp(AliEMCALQADataMaker)
48 //____________________________________________________________________________
49 AliEMCALQADataMaker::AliEMCALQADataMaker() :
50 AliQADataMaker(AliQA::GetDetName(AliQA::kEMCAL), "EMCAL Quality Assurance Data Maker")
55 //____________________________________________________________________________
56 AliEMCALQADataMaker::AliEMCALQADataMaker(const AliEMCALQADataMaker& qadm) :
60 SetName((const char*)qadm.GetName()) ;
61 SetTitle((const char*)qadm.GetTitle());
64 //__________________________________________________________________
65 AliEMCALQADataMaker& AliEMCALQADataMaker::operator = (const AliEMCALQADataMaker& qadm )
68 this->~AliEMCALQADataMaker();
69 new(this) AliEMCALQADataMaker(qadm);
73 //____________________________________________________________________________
74 void AliEMCALQADataMaker::EndOfDetectorCycle(AliQA::TASKINDEX task, TObjArray * list)
76 //Detector specific actions at end of cycle
78 AliQAChecker::Instance()->Run(AliQA::kEMCAL, task, list) ;
81 //____________________________________________________________________________
82 void AliEMCALQADataMaker::InitESDs()
84 //create ESDs histograms in ESDs subdir
86 TH1F * h0 = new TH1F("hEmcalESDs", "ESDs energy distribution in EMCAL", 100, 0., 100.) ;
89 TH1I * h1 = new TH1I("hEmcalESDsMul", "ESDs multiplicity distribution in EMCAL", 100, 0., 100) ;
94 //____________________________________________________________________________
95 void AliEMCALQADataMaker::InitHits()
97 // create Hits histograms in Hits subdir
98 TH1F * h0 = new TH1F("hEmcalHits", "Hits energy distribution in EMCAL", 100, 0., 100.) ;
100 Add2HitsList(h0, 0) ;
101 TH1I * h1 = new TH1I("hEmcalHitsMul", "Hits multiplicity distribution in EMCAL", 500, 0., 10000) ;
103 Add2HitsList(h1, 1) ;
106 //____________________________________________________________________________
107 void AliEMCALQADataMaker::InitDigits()
109 // create Digits histograms in Digits subdir
110 TH1I * h0 = new TH1I("hEmcalDigits", "Digits amplitude distribution in EMCAL", 500, 0, 5000) ;
112 Add2DigitsList(h0, 0) ;
113 TH1I * h1 = new TH1I("hEmcalDigitsMul", "Digits multiplicity distribution in EMCAL", 500, 0, 1000) ;
115 Add2DigitsList(h1, 1) ;
118 //____________________________________________________________________________
119 void AliEMCALQADataMaker::InitRecPoints()
121 // create Reconstructed Points histograms in RecPoints subdir
122 TH1F * h0 = new TH1F("hEmcalRecPoints", "RecPoints energy distribution in EMCAL", 100, 0., 100.) ;
124 Add2RecPointsList(h0, 0) ;
125 TH1I * h1 = new TH1I("hEmcalRecPointsMul", "RecPoints multiplicity distribution in EMCAL", 100, 0, 100) ;
127 Add2RecPointsList(h1, 1) ;
131 //____________________________________________________________________________
132 void AliEMCALQADataMaker::InitRaws()
134 AliInfo(Form("Raw QA infor for EMCAL not yet implemented"));
137 //____________________________________________________________________________
138 void AliEMCALQADataMaker::InitSDigits()
140 // create SDigits histograms in SDigits subdir
141 TH1F * h0 = new TH1F("hEmcalSDigits", "SDigits energy distribution in EMCAL", 100, 0., 100.) ;
143 Add2SDigitsList(h0, 0) ;
144 TH1I * h1 = new TH1I("hEmcalSDigitsMul", "SDigits multiplicity distribution in EMCAL", 500, 0, 10000) ;
146 Add2SDigitsList(h1, 1) ;
149 //____________________________________________________________________________
150 void AliEMCALQADataMaker::MakeESDs(AliESDEvent * esd)
152 // make QA data from ESDs
155 for ( Int_t index = 0; index < esd->GetNumberOfCaloClusters() ; index++ ) {
156 AliESDCaloCluster * clu = esd->GetCaloCluster(index) ;
157 if ( clu->IsEMCAL() ) {
158 GetESDsData(0)->Fill(clu->E()) ;
162 GetESDsData(1)->Fill(count) ;
165 //____________________________________________________________________________
166 void AliEMCALQADataMaker::MakeHits(TClonesArray * hits)
168 //make QA data from Hits
170 GetHitsData(1)->Fill(hits->GetEntriesFast()) ;
173 while ( (hit = dynamic_cast<AliEMCALHit *>(next())) ) {
174 GetHitsData(0)->Fill( hit->GetEnergy()) ;
178 //____________________________________________________________________________
179 void AliEMCALQADataMaker::MakeHits(TTree * hitTree)
181 // make QA data from Hit Tree
183 TClonesArray * hits = new TClonesArray("AliEMCALHit", 1000);
185 TBranch * branch = hitTree->GetBranch("EMCAL") ;
187 AliWarning("EMCAL branch in Hit Tree not found") ;
189 TClonesArray * tmp = new TClonesArray("AliEMCALHit", 1000) ;
190 branch->SetAddress(&tmp) ;
192 for (Int_t ientry = 0 ; ientry < branch->GetEntries() ; ientry++) {
193 branch->GetEntry(ientry) ;
194 for (Int_t ihit = 0 ; ihit < tmp->GetEntries() ; ihit++) {
195 AliEMCALHit * hit = dynamic_cast<AliEMCALHit *> (tmp->At(ihit)) ;
196 new((*hits)[index]) AliEMCALHit(*hit) ;
206 //____________________________________________________________________________
207 void AliEMCALQADataMaker::MakeDigits(TClonesArray * digits)
209 // makes data from Digits
211 GetDigitsData(1)->Fill(digits->GetEntriesFast()) ;
213 AliEMCALDigit * digit ;
214 while ( (digit = dynamic_cast<AliEMCALDigit *>(next())) ) {
215 GetDigitsData(0)->Fill( digit->GetAmp()) ;
219 //____________________________________________________________________________
220 void AliEMCALQADataMaker::MakeDigits(TTree * digitTree)
222 // makes data from Digit Tree
223 TClonesArray * digits = new TClonesArray("AliEMCALDigit", 1000) ;
225 TBranch * branch = digitTree->GetBranch("EMCAL") ;
227 AliWarning("EMCAL branch in Digit Tree not found") ;
229 branch->SetAddress(&digits) ;
230 branch->GetEntry(0) ;
235 //____________________________________________________________________________
236 void AliEMCALQADataMaker::MakeRaws(AliRawReader* rawReader)
238 //Raw QA info not yet implemented for EMCAL
242 //____________________________________________________________________________
243 void AliEMCALQADataMaker::MakeRecPoints(TTree * clustersTree)
245 // makes data from RecPoints
246 TBranch *emcbranch = clustersTree->GetBranch("EMCALECARP");
248 AliError("can't get the branch with the EMCAL clusters !");
251 TObjArray * emcrecpoints = new TObjArray(100) ;
252 emcbranch->SetAddress(&emcrecpoints);
253 emcbranch->GetEntry(0);
255 GetRecPointsData(1)->Fill(emcrecpoints->GetEntriesFast()) ;
256 TIter next(emcrecpoints) ;
257 AliEMCALRecPoint * rp ;
258 while ( (rp = dynamic_cast<AliEMCALRecPoint *>(next())) ) {
259 GetRecPointsData(0)->Fill( rp->GetEnergy()) ;
261 emcrecpoints->Delete();
266 //____________________________________________________________________________
267 void AliEMCALQADataMaker::MakeSDigits(TClonesArray * sdigits)
269 // makes data from SDigits
271 GetSDigitsData(1)->Fill(sdigits->GetEntriesFast()) ;
272 TIter next(sdigits) ;
273 AliEMCALDigit * sdigit ;
274 while ( (sdigit = dynamic_cast<AliEMCALDigit *>(next())) ) {
275 GetSDigitsData(0)->Fill( sdigit->GetAmp()) ;
279 //____________________________________________________________________________
280 void AliEMCALQADataMaker::MakeSDigits(TTree * sdigitTree)
282 // makes data from SDigit Tree
283 TClonesArray * sdigits = new TClonesArray("AliEMCALDigit", 1000) ;
285 TBranch * branch = sdigitTree->GetBranch("EMCAL") ;
287 AliWarning("EMCAL branch in SDigit Tree not found") ;
289 branch->SetAddress(&sdigits) ;
290 branch->GetEntry(0) ;
291 MakeSDigits(sdigits) ;
295 //____________________________________________________________________________
296 void AliEMCALQADataMaker::StartOfDetectorCycle()
298 //Detector specific actions at start of cycle