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>
32 // --- Standard library ---
34 // --- AliRoot header files ---
35 #include "AliESDCaloCluster.h"
37 #include "AliEMCALDigit.h"
38 #include "AliEMCALHit.h"
39 #include "AliEMCALQADataMakerSim.h"
40 #include "AliQAChecker.h"
42 ClassImp(AliEMCALQADataMakerSim)
44 //____________________________________________________________________________
45 AliEMCALQADataMakerSim::AliEMCALQADataMakerSim() :
46 AliQADataMakerSim(AliQA::GetDetName(AliQA::kEMCAL), "EMCAL Quality Assurance Data Maker")
51 //____________________________________________________________________________
52 AliEMCALQADataMakerSim::AliEMCALQADataMakerSim(const AliEMCALQADataMakerSim& qadm) :
56 SetName((const char*)qadm.GetName()) ;
57 SetTitle((const char*)qadm.GetTitle());
60 //__________________________________________________________________
61 AliEMCALQADataMakerSim& AliEMCALQADataMakerSim::operator = (const AliEMCALQADataMakerSim& qadm )
64 this->~AliEMCALQADataMakerSim();
65 new(this) AliEMCALQADataMakerSim(qadm);
69 //____________________________________________________________________________
70 void AliEMCALQADataMakerSim::EndOfDetectorCycle(AliQA::TASKINDEX_t task, TObjArray * list)
72 //Detector specific actions at end of cycle
74 AliQAChecker::Instance()->Run(AliQA::kEMCAL, task, list) ;
77 //____________________________________________________________________________
78 void AliEMCALQADataMakerSim::InitHits()
80 // create Hits histograms in Hits subdir
81 TH1F * h0 = new TH1F("hEmcalHits", "Hits energy distribution in EMCAL", 100, 0., 100.) ;
84 TH1I * h1 = new TH1I("hEmcalHitsMul", "Hits multiplicity distribution in EMCAL", 500, 0., 10000) ;
89 //____________________________________________________________________________
90 void AliEMCALQADataMakerSim::InitDigits()
92 // create Digits histograms in Digits subdir
93 TH1I * h0 = new TH1I("hEmcalDigits", "Digits amplitude distribution in EMCAL", 500, 0, 5000) ;
95 Add2DigitsList(h0, 0) ;
96 TH1I * h1 = new TH1I("hEmcalDigitsMul", "Digits multiplicity distribution in EMCAL", 500, 0, 1000) ;
98 Add2DigitsList(h1, 1) ;
101 //____________________________________________________________________________
102 void AliEMCALQADataMakerSim::InitSDigits()
104 // create SDigits histograms in SDigits subdir
105 TH1F * h0 = new TH1F("hEmcalSDigits", "SDigits energy distribution in EMCAL", 100, 0., 100.) ;
107 Add2SDigitsList(h0, 0) ;
108 TH1I * h1 = new TH1I("hEmcalSDigitsMul", "SDigits multiplicity distribution in EMCAL", 500, 0, 10000) ;
110 Add2SDigitsList(h1, 1) ;
113 //____________________________________________________________________________
114 void AliEMCALQADataMakerSim::MakeHits(TClonesArray * hits)
116 //make QA data from Hits
118 GetHitsData(1)->Fill(hits->GetEntriesFast()) ;
121 while ( (hit = dynamic_cast<AliEMCALHit *>(next())) ) {
122 GetHitsData(0)->Fill( hit->GetEnergy()) ;
127 //____________________________________________________________________________
128 void AliEMCALQADataMakerSim::MakeHits(TTree * hitTree)
130 // make QA data from Hit Tree
132 TClonesArray * hits = new TClonesArray("AliEMCALHit", 1000);
134 TBranch * branch = hitTree->GetBranch("EMCAL") ;
136 AliWarning("EMCAL branch in Hit Tree not found") ;
138 TClonesArray * tmp = new TClonesArray("AliEMCALHit", 1000) ;
139 branch->SetAddress(&tmp) ;
141 for (Int_t ientry = 0 ; ientry < branch->GetEntries() ; ientry++) {
142 branch->GetEntry(ientry) ;
143 for (Int_t ihit = 0 ; ihit < tmp->GetEntries() ; ihit++) {
144 AliEMCALHit * hit = dynamic_cast<AliEMCALHit *> (tmp->At(ihit)) ;
145 new((*hits)[index]) AliEMCALHit(*hit) ;
155 //____________________________________________________________________________
156 void AliEMCALQADataMakerSim::MakeDigits(TClonesArray * digits)
158 // makes data from Digits
160 GetDigitsData(1)->Fill(digits->GetEntriesFast()) ;
162 AliEMCALDigit * digit ;
163 while ( (digit = dynamic_cast<AliEMCALDigit *>(next())) ) {
164 GetDigitsData(0)->Fill( digit->GetAmp()) ;
169 //____________________________________________________________________________
170 void AliEMCALQADataMakerSim::MakeDigits(TTree * digitTree)
172 // makes data from Digit Tree
173 TClonesArray * digits = new TClonesArray("AliEMCALDigit", 1000) ;
175 TBranch * branch = digitTree->GetBranch("EMCAL") ;
177 AliWarning("EMCAL branch in Digit Tree not found") ;
179 branch->SetAddress(&digits) ;
180 branch->GetEntry(0) ;
186 //____________________________________________________________________________
187 void AliEMCALQADataMakerSim::MakeSDigits(TClonesArray * sdigits)
189 // makes data from SDigits
191 GetSDigitsData(1)->Fill(sdigits->GetEntriesFast()) ;
192 TIter next(sdigits) ;
193 AliEMCALDigit * sdigit ;
194 while ( (sdigit = dynamic_cast<AliEMCALDigit *>(next())) ) {
195 GetSDigitsData(0)->Fill( sdigit->GetAmp()) ;
200 //____________________________________________________________________________
201 void AliEMCALQADataMakerSim::MakeSDigits(TTree * sdigitTree)
203 // makes data from SDigit Tree
204 TClonesArray * sdigits = new TClonesArray("AliEMCALDigit", 1000) ;
206 TBranch * branch = sdigitTree->GetBranch("EMCAL") ;
208 AliWarning("EMCAL branch in SDigit Tree not found") ;
210 branch->SetAddress(&sdigits) ;
211 branch->GetEntry(0) ;
212 MakeSDigits(sdigits) ;
217 //____________________________________________________________________________
218 void AliEMCALQADataMakerSim::StartOfDetectorCycle()
220 //Detector specific actions at start of cycle