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"
41 #include "AliEMCALSDigitizer.h"
43 ClassImp(AliEMCALQADataMakerSim)
45 //____________________________________________________________________________
46 AliEMCALQADataMakerSim::AliEMCALQADataMakerSim() :
47 AliQADataMakerSim(AliQAv1::GetDetName(AliQAv1::kEMCAL), "EMCAL Quality Assurance Data Maker")
52 //____________________________________________________________________________
53 AliEMCALQADataMakerSim::AliEMCALQADataMakerSim(const AliEMCALQADataMakerSim& qadm) :
57 SetName((const char*)qadm.GetName()) ;
58 SetTitle((const char*)qadm.GetTitle());
61 //__________________________________________________________________
62 AliEMCALQADataMakerSim& AliEMCALQADataMakerSim::operator = (const AliEMCALQADataMakerSim& qadm )
65 this->~AliEMCALQADataMakerSim();
66 new(this) AliEMCALQADataMakerSim(qadm);
70 //____________________________________________________________________________
71 void AliEMCALQADataMakerSim::EndOfDetectorCycle(AliQAv1::TASKINDEX_t task, TObjArray ** list)
73 //Detector specific actions at end of cycle
75 ResetEventTrigClasses(); // reset triggers list to select all histos
76 AliQAChecker::Instance()->Run(AliQAv1::kEMCAL, task, list) ;
79 //____________________________________________________________________________
80 void AliEMCALQADataMakerSim::InitHits()
82 // create Hits histograms in Hits subdir
83 const Bool_t expert = kTRUE ;
84 const Bool_t image = kTRUE ;
86 TH1F * h0 = new TH1F("hEmcalHits", "Hits energy distribution in EMCAL;Energy [MeV];Counts", 200, 0., 2.) ; //GeV
88 Add2HitsList(h0, 0, !expert, image) ;
89 TH1I * h1 = new TH1I("hEmcalHitsMul", "Hits multiplicity distribution in EMCAL;# of Hits;Entries", 1000, 0, 10000) ;
91 Add2HitsList(h1, 1, !expert, image) ;
93 ClonePerTrigClass(AliQAv1::kHITS); // this should be the last line
96 //____________________________________________________________________________
97 void AliEMCALQADataMakerSim::InitDigits()
99 // create Digits histograms in Digits subdir
100 const Bool_t expert = kTRUE ;
101 const Bool_t image = kTRUE ;
103 TH1I * h0 = new TH1I("hEmcalDigits", "Digits amplitude distribution in EMCAL;Amplitude [ADC counts];Counts", 500, 0, 500) ;
105 Add2DigitsList(h0, 0, !expert, image) ;
106 TH1I * h1 = new TH1I("hEmcalDigitsMul", "Digits multiplicity distribution in EMCAL;# of Digits;Entries", 200, 0, 2000) ;
108 Add2DigitsList(h1, 1, !expert, image) ;
110 ClonePerTrigClass(AliQAv1::kDIGITS); // this should be the last line
113 //____________________________________________________________________________
114 void AliEMCALQADataMakerSim::InitSDigits()
116 // create SDigits histograms in SDigits subdir
117 const Bool_t expert = kTRUE ;
118 const Bool_t image = kTRUE ;
120 TH1F * h0 = new TH1F("hEmcalSDigits", "SDigits energy distribution in EMCAL;Energy [MeV];Counts", 200, 0., 20.) ;
122 Add2SDigitsList(h0, 0, !expert, image) ;
123 TH1I * h1 = new TH1I("hEmcalSDigitsMul", "SDigits multiplicity distribution in EMCAL;# of SDigits;Entries", 500, 0, 5000) ;
125 Add2SDigitsList(h1, 1, !expert, image) ;
127 ClonePerTrigClass(AliQAv1::kSDIGITS); // this should be the last line
130 //____________________________________________________________________________
131 void AliEMCALQADataMakerSim::MakeHits()
133 //make QA data from Hits
135 FillHitsData(1,fHitsArray->GetEntriesFast()) ;
136 TIter next(fHitsArray) ;
138 while ( (hit = dynamic_cast<AliEMCALHit *>(next())) ) {
139 FillHitsData(0, hit->GetEnergy()) ;
143 //____________________________________________________________________________
144 void AliEMCALQADataMakerSim::MakeHits(TTree * hitTree)
146 // make QA data from Hit Tree
148 fHitsArray->Clear() ;
150 fHitsArray = new TClonesArray("AliEMCALHit", 1000);
152 TBranch * branch = hitTree->GetBranch("EMCAL") ;
153 if ( ! branch ) { AliWarning("EMCAL branch in Hit Tree not found") ; return;}
155 branch->SetAddress(&fHitsArray) ;
156 for (Int_t ientry = 0 ; ientry < branch->GetEntries() ; ientry++) {
157 branch->GetEntry(ientry) ;
159 fHitsArray->Clear() ;
161 IncEvCountCycleHits();
162 IncEvCountTotalHits();
166 //____________________________________________________________________________
167 void AliEMCALQADataMakerSim::MakeDigits()
169 // makes data from Digits
171 FillDigitsData(1,fDigitsArray->GetEntriesFast()) ;
172 TIter next(fDigitsArray) ;
173 AliEMCALDigit * digit ;
174 while ( (digit = dynamic_cast<AliEMCALDigit *>(next())) ) {
175 FillDigitsData(0, digit->GetAmp()) ;
179 //____________________________________________________________________________
180 void AliEMCALQADataMakerSim::MakeDigits(TTree * digitTree)
182 // makes data from Digit Tree
185 fDigitsArray->Clear("C") ;
187 fDigitsArray = new TClonesArray("AliEMCALDigit", 1000) ;
189 TBranch * branch = digitTree->GetBranch("EMCAL") ;
190 if ( ! branch ) { AliWarning("EMCAL branch in Digit Tree not found") ; return; }
192 branch->SetAddress(&fDigitsArray) ;
193 branch->GetEntry(0) ;
196 IncEvCountCycleDigits();
197 IncEvCountTotalDigits();
200 //____________________________________________________________________________
201 void AliEMCALQADataMakerSim::MakeSDigits()
203 // makes data from SDigits
204 //Need a copy of the SDigitizer to calibrate the sdigit amplitude to
207 AliEMCALSDigitizer* sDigitizer = new AliEMCALSDigitizer();
209 FillSDigitsData(1,fSDigitsArray->GetEntriesFast()) ;
210 TIter next(fSDigitsArray) ;
211 AliEMCALDigit * sdigit ;
212 while ( (sdigit = dynamic_cast<AliEMCALDigit *>(next())) ) {
213 FillSDigitsData(0, sDigitizer->Calibrate(sdigit->GetAmp())) ;
218 //____________________________________________________________________________
219 void AliEMCALQADataMakerSim::MakeSDigits(TTree * sdigitTree)
221 // makes data from SDigit Tree
223 fSDigitsArray->Clear("C") ;
225 fSDigitsArray = new TClonesArray("AliEMCALDigit", 1000) ;
227 TBranch * branch = sdigitTree->GetBranch("EMCAL") ;
228 if ( ! branch ) { AliWarning("EMCAL branch in SDigit Tree not found"); return;}
230 branch->SetAddress(&fSDigitsArray);
234 IncEvCountCycleSDigits();
235 IncEvCountTotalSDigits();
239 //____________________________________________________________________________
240 void AliEMCALQADataMakerSim::StartOfDetectorCycle()
242 //Detector specific actions at start of cycle