1 /**************************************************************************
2 * Copyright(c) 2004, 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 **************************************************************************/
15 // --- ROOT system ---
17 #include <TClonesArray.h>
22 // --- AliRoot header files ---
23 #include "AliESDEvent.h"
25 #include "AliFMDQADataMakerSim.h"
26 #include "AliFMDDigit.h"
27 #include "AliFMDHit.h"
28 #include "AliQAChecker.h"
29 #include "AliFMDParameters.h"
30 #include "AliFMDSDigit.h"
32 //_____________________________________________________________________
33 // This is the class that collects the QA data for the FMD during simulation.
34 // The following data types are picked up:
37 // The following data types are not supported (yet):
40 // Author : Hans Hjersing Dalsgaard, Niels Bohr Institute, hans.dalsgaard@cern.ch
41 //_____________________________________________________________________
43 ClassImp(AliFMDQADataMakerSim)
45 ; // This line is for Emacs - do not delete!
47 //_____________________________________________________________________
48 AliFMDQADataMakerSim::AliFMDQADataMakerSim()
49 : AliQADataMakerSim(AliQAv1::GetDetName(AliQAv1::kFMD),
50 "FMD Quality Assurance Data Maker"),
51 fSDigitsArray("AliFMDSDigit", 1000),
52 fDigitsArray("AliFMDDigit", 1000),
53 fHitsArray("AliFMDHit", 10)
59 //_____________________________________________________________________
60 AliFMDQADataMakerSim::AliFMDQADataMakerSim(const AliFMDQADataMakerSim& qadm)
61 : AliQADataMakerSim(),
62 fSDigitsArray(qadm.fSDigitsArray),
63 fDigitsArray(qadm.fDigitsArray),
64 fHitsArray(qadm.fHitsArray)
69 // qadm Object to copy from
72 //_____________________________________________________________________
73 AliFMDQADataMakerSim& AliFMDQADataMakerSim::operator = (const AliFMDQADataMakerSim& qadm )
75 fSDigitsArray = qadm.fSDigitsArray;
76 fDigitsArray = qadm.fDigitsArray;
77 fHitsArray = qadm.fHitsArray;
81 //_____________________________________________________________________
82 AliFMDQADataMakerSim::~AliFMDQADataMakerSim()
87 //_____________________________________________________________________
88 void AliFMDQADataMakerSim::EndOfDetectorCycle(AliQAv1::TASKINDEX_t task,
91 //Detector specific actions at end of cycle
93 AliLog::Message(5,"FMD: end of detector cycle",
94 "AliFMDQADataMakerSim","AliFMDQADataMakerSim",
95 "AliFMDQADataMakerSim::EndOfDetectorCycle",
96 "AliFMDQADataMakerSim.cxx",83);
97 AliQAChecker::Instance()->Run(AliQAv1::kFMD, task, list) ;
100 //_____________________________________________________________________
101 void AliFMDQADataMakerSim::InitSDigits()
103 // create SDigits histograms in SDigits subdir
104 const Bool_t expert = kTRUE ;
105 const Bool_t image = kTRUE ;
107 TH1I* hADCCounts = new TH1I("hADCCounts","Dist of ADC counts",1024,0,1024);
108 hADCCounts->SetXTitle("ADC counts");
109 Add2SDigitsList(hADCCounts, 0, !expert, image);
112 //____________________________________________________________________
113 void AliFMDQADataMakerSim::InitHits()
115 // create Digits histograms in Digits subdir
116 const Bool_t expert = kTRUE ;
117 const Bool_t image = kTRUE ;
119 TH1F* hEnergyOfHits = new TH1F("hEnergyOfHits","Energy distribution",100,0,3);
120 hEnergyOfHits->SetXTitle("Edep");
121 hEnergyOfHits->SetYTitle("Counts");
122 Add2HitsList(hEnergyOfHits, 0, !expert, image);
125 //_____________________________________________________________________
126 void AliFMDQADataMakerSim::InitDigits()
128 // create Digits histograms in Digits subdir
129 const Bool_t expert = kTRUE ;
130 const Bool_t image = kTRUE ;
132 TH1I* hADCCounts = new TH1I("hADCCounts","Dist of ADC counts",1024,0,1024);
133 hADCCounts->SetXTitle("ADC counts");
134 Add2DigitsList(hADCCounts, 0, !expert, image);
137 //_____________________________________________________________________
138 void AliFMDQADataMakerSim::MakeHits(TClonesArray * hits)
142 while ((hit = static_cast<AliFMDHit *>(next())))
143 GetHitsData(0)->Fill(hit->Edep()/hit->Length()*0.032);
146 //_____________________________________________________________________
147 void AliFMDQADataMakerSim::MakeHits(TTree * hitTree)
149 // make QA data from Hit Tree
153 TBranch * branch = hitTree->GetBranch("FMD") ;
155 AliWarning("FMD branch in Hit Tree not found") ;
159 TClonesArray* hitsAddress = &fHitsArray;
161 branch->SetAddress(&hitsAddress) ;
163 for (Int_t ientry = 0 ; ientry < branch->GetEntries() ; ientry++) {
164 branch->GetEntry(ientry);
165 MakeHits(hitsAddress); //tmp);
170 //_____________________________________________________________________
171 void AliFMDQADataMakerSim::MakeDigits(TClonesArray * digits)
173 // makes data from Digits
176 for(Int_t i = 0 ; i < digits->GetEntriesFast() ; i++) {
178 AliFMDDigit* digit = static_cast<AliFMDDigit*>(digits->At(i));
179 GetDigitsData(0)->Fill(digit->Counts());
183 //_____________________________________________________________________
184 void AliFMDQADataMakerSim::MakeDigits(TTree * digitTree)
187 fDigitsArray.Clear();
188 TBranch * branch = digitTree->GetBranch("FMD") ;
190 AliWarning("FMD branch in Digit Tree not found") ;
193 TClonesArray* digitAddress = &fDigitsArray;
194 branch->SetAddress(&digitAddress) ;
195 branch->GetEntry(0) ;
196 MakeDigits(digitAddress) ;
199 //_____________________________________________________________________
200 void AliFMDQADataMakerSim::MakeSDigits(TClonesArray * sdigits)
202 // makes data from Digits
205 for(Int_t i = 0 ; i < sdigits->GetEntriesFast() ; i++) {
207 AliFMDSDigit* sdigit = static_cast<AliFMDSDigit*>(sdigits->At(i));
208 GetSDigitsData(0)->Fill(sdigit->Counts());
212 //_____________________________________________________________________
213 void AliFMDQADataMakerSim::MakeSDigits(TTree * sdigitTree)
216 fSDigitsArray.Clear();
217 TBranch * branch = sdigitTree->GetBranch("FMD") ;
219 AliWarning("FMD branch in SDigit Tree not found") ;
222 TClonesArray* sdigitAddress = &fSDigitsArray;
223 branch->SetAddress(&sdigitAddress) ;
224 branch->GetEntry(0) ;
225 MakeSDigits(sdigitAddress) ;
228 //_____________________________________________________________________
229 void AliFMDQADataMakerSim::StartOfDetectorCycle()
233 //_____________________________________________________________________