]>
Commit | Line | Data |
---|---|---|
045862db | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
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 | // Class for bad channels & bad runs identification | |
16 | // | |
17 | // This class is designed for the following tasks: | |
18 | // 1. find bad (dead/noisy/strange) cells on a per run basis; | |
19 | // 2. find the status/quality of analysed data (e.g. missing RCUs, run quality); | |
20 | // 3. find the extent of problems related to bad cells: required for | |
21 | // systematics estimation for a physical quantity related to a cluster. | |
22 | // | |
23 | // See also AliAnalysisTaskCellsQA. | |
24 | // Works for EMCAL and PHOS. Requires initialized geometry. | |
d30ac678 | 25 | // Class defaults are optimized for EMCAL. |
045862db | 26 | // |
27 | // Usage example for EMCAL: | |
28 | // | |
29 | // a) In LocalInit(): | |
30 | // | |
31 | // // to check the quality of analysed data and detect most of the problems: | |
32 | // AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1"); | |
33 | // cellsQA = new AliCaloCellsQA(10, AliCaloCellsQA::kEMCAL); // 10 supermodules | |
34 | // | |
35 | // // add this line for a full analysis (fills additional histograms): | |
36 | // cellsQA->ActivateFullAnalysis(); | |
37 | // | |
38 | // b) In UserCreateOutputObjects(): | |
39 | // | |
40 | // // not required, but suggested | |
41 | // cellsQA->InitSummaryHistograms(); | |
42 | // | |
43 | // | |
44 | // c) In UserExec(): | |
45 | // | |
46 | // AliVEvent *event = InputEvent(); | |
47 | // AliVCaloCells* cells = event->GetEMCALCells(); | |
48 | // | |
49 | // AliVVertex *vertex = (AliVVertex*) event->GetPrimaryVertex(); | |
50 | // Double_t vertexXYZ[3]; | |
51 | // vertex->GetXYZ(vertexXYZ); | |
52 | // | |
53 | // // clusters | |
f4dd0897 | 54 | // TObjArray clusArray; |
045862db | 55 | // for (Int_t i = 0; i < event->GetNumberOfCaloClusters(); i++) { |
56 | // AliVCluster *clus = event->GetCaloCluster(i); | |
57 | // | |
58 | // // filter clusters here, if necessary | |
0d81af58 | 59 | // // if (clus->E() < 0.3) continue; |
f4dd0897 | 60 | // // if (clus->GetNCells() < 2) continue; |
045862db | 61 | // |
f4dd0897 | 62 | // clusArray.Add(clus); |
045862db | 63 | // } |
64 | // | |
0d81af58 | 65 | // // apply afterburners, etc. |
045862db | 66 | // // ... |
67 | // | |
f4dd0897 | 68 | // cellsQA->Fill(event->GetRunNumber(), &clusArray, cells, vertexXYZ); |
045862db | 69 | // |
70 | // d) Do not forget to post data, where necessary: | |
71 | // PostData(1,cellsQA->GetListOfHistos()); | |
045862db | 72 | // |
73 | // TODO: add DCAL (up to 6 supermodules?) | |
74 | // | |
75 | //---- | |
76 | // Author: Olga Driga (SUBATECH) | |
77 | ||
78 | // --- ROOT system --- | |
79 | #include <TLorentzVector.h> | |
80 | ||
81 | // --- AliRoot header files --- | |
b028cd0a | 82 | #include <AliLog.h> |
045862db | 83 | #include <AliEMCALGeometry.h> |
84 | #include <AliPHOSGeometry.h> | |
85 | #include <AliCaloCellsQA.h> | |
86 | ||
d30ac678 | 87 | ClassImp(AliCaloCellsQA) |
045862db | 88 | |
89 | //_________________________________________________________________________ | |
90 | AliCaloCellsQA::AliCaloCellsQA() : | |
91 | fDetector(0), | |
92 | fNMods(0), | |
93 | fClusElowMin(0), | |
94 | fClusEhighMin(0), | |
95 | fPi0EClusMin(0), | |
96 | fkFullAnalysis(0), | |
0d81af58 | 97 | fNBinsECells(0), |
98 | fNBinsPi0Mass(0), | |
99 | fNBinsXNCellsInCluster(0), | |
100 | fNBinsYNCellsInCluster(0), | |
101 | fXMaxECells(0), | |
102 | fXMaxPi0Mass(0), | |
103 | fXMaxNCellsInCluster(0), | |
f4dd0897 | 104 | fRunNumbers(), |
105 | fNRuns(0), | |
d30ac678 | 106 | fRI(0), |
045862db | 107 | fAbsIdMin(0), |
108 | fAbsIdMax(0), | |
109 | fListOfHistos(0), | |
110 | fhNEventsProcessedPerRun(0), | |
111 | fhCellLocMaxNTimesInClusterElow(), | |
112 | fhCellLocMaxNTimesInClusterEhigh(), | |
113 | fhCellLocMaxETotalClusterElow(), | |
114 | fhCellLocMaxETotalClusterEhigh(), | |
115 | fhCellNonLocMaxNTimesInClusterElow(), | |
116 | fhCellNonLocMaxNTimesInClusterEhigh(), | |
117 | fhCellNonLocMaxETotalClusterElow(), | |
118 | fhCellNonLocMaxETotalClusterEhigh(), | |
119 | fhECells(), | |
120 | fhPi0Mass(), | |
121 | fhNCellsInCluster(), | |
122 | fhCellAmplitude(0), | |
123 | fhCellAmplitudeEhigh(0), | |
124 | fhCellAmplitudeNonLocMax(0), | |
125 | fhCellAmplitudeEhighNonLocMax(0), | |
126 | fhCellTime(0) | |
127 | { | |
d30ac678 | 128 | // Default constructor (for root I/O). Do not use it. |
045862db | 129 | |
d30ac678 | 130 | // tells the class to reinitialize its transient members |
131 | fRI = -1; | |
045862db | 132 | } |
133 | ||
134 | //_________________________________________________________________________ | |
135 | AliCaloCellsQA::AliCaloCellsQA(Int_t nmods, Int_t det, Int_t startRunNumber, Int_t endRunNumber) : | |
136 | fDetector(0), | |
137 | fNMods(0), | |
138 | fClusElowMin(0), | |
139 | fClusEhighMin(0), | |
140 | fPi0EClusMin(0), | |
141 | fkFullAnalysis(0), | |
0d81af58 | 142 | fNBinsECells(0), |
143 | fNBinsPi0Mass(0), | |
144 | fNBinsXNCellsInCluster(0), | |
145 | fNBinsYNCellsInCluster(0), | |
146 | fXMaxECells(0), | |
147 | fXMaxPi0Mass(0), | |
148 | fXMaxNCellsInCluster(0), | |
f4dd0897 | 149 | fRunNumbers(), |
150 | fNRuns(0), | |
d30ac678 | 151 | fRI(0), |
045862db | 152 | fAbsIdMin(0), |
153 | fAbsIdMax(0), | |
154 | fListOfHistos(0), | |
155 | fhNEventsProcessedPerRun(0), | |
156 | fhCellLocMaxNTimesInClusterElow(), | |
157 | fhCellLocMaxNTimesInClusterEhigh(), | |
158 | fhCellLocMaxETotalClusterElow(), | |
159 | fhCellLocMaxETotalClusterEhigh(), | |
160 | fhCellNonLocMaxNTimesInClusterElow(), | |
161 | fhCellNonLocMaxNTimesInClusterEhigh(), | |
162 | fhCellNonLocMaxETotalClusterElow(), | |
163 | fhCellNonLocMaxETotalClusterEhigh(), | |
164 | fhECells(), | |
165 | fhPi0Mass(), | |
166 | fhNCellsInCluster(), | |
167 | fhCellAmplitude(0), | |
168 | fhCellAmplitudeEhigh(0), | |
169 | fhCellAmplitudeNonLocMax(0), | |
170 | fhCellAmplitudeEhighNonLocMax(0), | |
171 | fhCellTime(0) | |
172 | { | |
173 | // Allows to set main class parameters. | |
174 | // | |
175 | // nmods -- maximum supermodule number + 1: | |
176 | // use 4 for EMCAL <= 2010; | |
177 | // use 4 for PHOS (PHOS numbers start from 1, not from zero); | |
178 | // use 10 for EMCAL >= 2011. | |
179 | // | |
180 | // det -- detector: AliCaloCellsQA::kEMCAL or AliCaloCellsQA::kPHOS. | |
181 | // Note: DCAL not yet implemented. | |
182 | // | |
183 | // startRunNumber, endRunNumber -- run range the analysis will run on | |
184 | // (to allow later merging); wide (but reasonable) range can be provided. | |
185 | ||
d30ac678 | 186 | fRI = -1; |
187 | ||
045862db | 188 | Init(nmods, det, startRunNumber, endRunNumber); |
189 | } | |
190 | ||
191 | //_________________________________________________________________________ | |
192 | AliCaloCellsQA::~AliCaloCellsQA() | |
193 | { | |
194 | delete fListOfHistos; | |
195 | } | |
196 | ||
197 | //_________________________________________________________________________ | |
198 | void AliCaloCellsQA::ActivateFullAnalysis() | |
199 | { | |
200 | // Activate initialization and filling of all the designed histograms. | |
201 | // Must be called before InitSummaryHistograms() and Fill() methods, if necessary. | |
202 | ||
203 | fkFullAnalysis = kTRUE; | |
204 | } | |
205 | ||
206 | //_________________________________________________________________________ | |
207 | void AliCaloCellsQA::InitSummaryHistograms(Int_t nbins, Double_t emax, | |
208 | Int_t nbinsh, Double_t emaxh, | |
209 | Int_t nbinst, Double_t tmin, Double_t tmax) | |
210 | { | |
211 | // Initialize summary (not per run) histograms. | |
0d81af58 | 212 | // Must be called after ActivateFullAnalysis(), before calling Fill() method, if necessary. |
045862db | 213 | // |
214 | // nbins, emax -- number of bins and maximum amplitude for fhCellAmplitude and fhCellAmplitudeNonLocMax; | |
215 | // nbinsh, emaxh -- number of bins and maximum amplitude for fhCellAmplitudeEhigh and fhCellAmplitudeEhighNonLocMax; | |
216 | // nbinst, tmin, tmax -- number of bins and minimum/maximum time for fhCellTime. | |
217 | ||
218 | // do not add histograms to the current directory | |
219 | Bool_t ads = TH1::AddDirectoryStatus(); | |
220 | TH1::AddDirectory(kFALSE); | |
221 | ||
222 | fhCellAmplitude = new TH2F("hCellAmplitude", | |
223 | "Amplitude distribution per cell", fAbsIdMax-fAbsIdMin,fAbsIdMin,fAbsIdMax, nbins,0,emax); | |
224 | fhCellAmplitude->SetXTitle("AbsId"); | |
225 | fhCellAmplitude->SetYTitle("Amplitude, GeV"); | |
226 | fhCellAmplitude->SetZTitle("Counts"); | |
227 | ||
228 | fhCellAmplitudeEhigh = new TH2F("hCellAmplitudeEhigh", | |
229 | "Amplitude distribution per cell, high energies", fAbsIdMax-fAbsIdMin,fAbsIdMin,fAbsIdMax, nbinsh,0,emaxh); | |
230 | fhCellAmplitudeEhigh->SetXTitle("AbsId"); | |
231 | fhCellAmplitudeEhigh->SetYTitle("Amplitude, GeV"); | |
232 | fhCellAmplitudeEhigh->SetZTitle("Counts"); | |
233 | ||
0d81af58 | 234 | fhCellTime = new TH2F("hCellTime", "Time distribution per cell", |
235 | fAbsIdMax-fAbsIdMin,fAbsIdMin,fAbsIdMax, nbinst,tmin,tmax); | |
236 | fhCellTime->SetXTitle("AbsId"); | |
237 | fhCellTime->SetYTitle("Time, microseconds"); | |
238 | fhCellTime->SetZTitle("Counts"); | |
239 | ||
045862db | 240 | fListOfHistos->Add(fhCellAmplitude); |
241 | fListOfHistos->Add(fhCellAmplitudeEhigh); | |
0d81af58 | 242 | fListOfHistos->Add(fhCellTime); |
045862db | 243 | |
244 | if (fkFullAnalysis) { | |
245 | fhCellAmplitudeNonLocMax = new TH2F("hCellAmplitudeNonLocMax", | |
246 | "Amplitude distribution per cell which is not a local maximum", | |
247 | fAbsIdMax-fAbsIdMin,fAbsIdMin,fAbsIdMax, nbins,0,emax); | |
248 | fhCellAmplitudeNonLocMax->SetXTitle("AbsId"); | |
249 | fhCellAmplitudeNonLocMax->SetYTitle("Amplitude, GeV"); | |
250 | fhCellAmplitudeNonLocMax->SetZTitle("Counts"); | |
251 | ||
252 | fhCellAmplitudeEhighNonLocMax = new TH2F("hCellAmplitudeEhighNonLocMax", | |
253 | "Amplitude distribution per cell which is not a local maximum, high energies", | |
254 | fAbsIdMax-fAbsIdMin,fAbsIdMin,fAbsIdMax, nbinsh,0,emaxh); | |
255 | fhCellAmplitudeEhighNonLocMax->SetXTitle("AbsId"); | |
256 | fhCellAmplitudeEhighNonLocMax->SetYTitle("Amplitude, GeV"); | |
257 | fhCellAmplitudeEhighNonLocMax->SetZTitle("Counts"); | |
258 | ||
045862db | 259 | fListOfHistos->Add(fhCellAmplitudeNonLocMax); |
260 | fListOfHistos->Add(fhCellAmplitudeEhighNonLocMax); | |
045862db | 261 | } |
262 | ||
263 | // return to the previous add directory status | |
264 | TH1::AddDirectory(ads); | |
265 | } | |
266 | ||
267 | //_________________________________________________________________________ | |
268 | void AliCaloCellsQA::Fill(Int_t runNumber, TObjArray *clusArray, AliVCaloCells *cells, Double_t vertexXYZ[3]) | |
269 | { | |
270 | // Does the job: fills histograms for current event. | |
271 | // | |
272 | // runNumber -- current run number; | |
273 | // clusArray -- array of clusters (TObjArray or TClonesArray); | |
274 | // cells -- EMCAL or PHOS cells; | |
275 | // vertexXYZ -- primary vertex position. | |
276 | ||
d30ac678 | 277 | InitTransientFindCurrentRun(runNumber); |
278 | ||
045862db | 279 | fhNEventsProcessedPerRun->Fill(runNumber); |
045862db | 280 | |
d30ac678 | 281 | FillCellsInCluster(clusArray, cells); |
282 | FillJustCells(cells); | |
283 | FillPi0Mass(clusArray, vertexXYZ); | |
045862db | 284 | } |
285 | ||
286 | //_________________________________________________________________________ | |
287 | void AliCaloCellsQA::SetClusterEnergyCuts(Double_t pi0EClusMin, Double_t elowMin, Double_t ehighMin) | |
288 | { | |
289 | // Set cuts for minimum cluster energies. | |
290 | // Must be called before calling Fill() method, if necessary. | |
291 | ||
292 | fClusElowMin = elowMin; | |
293 | fClusEhighMin = ehighMin; | |
294 | fPi0EClusMin = pi0EClusMin; | |
295 | } | |
296 | ||
297 | //_________________________________________________________________________ | |
298 | void AliCaloCellsQA::SetBinningParameters(Int_t nbins1, Int_t nbins2, Int_t nbins3x, Int_t nbins3y, | |
299 | Double_t xmax1, Double_t xmax2, Double_t xmax3) | |
300 | { | |
301 | // Set binning parameters for histograms hECells, hNCellsInCluster and hPi0Mass. | |
302 | // Must be called before Fill() method, if necessary. | |
303 | // | |
304 | // nbins1, xmax1 -- number of bins and maximum X axis value for hECells; | |
305 | // nbins2, xmax2 -- number of bins and maximum X axis value for hPi0Mass; | |
306 | // nbins3x, nbins3y, xmax3 -- number of bins in X and Y axes and maximum X axis value for hNCellsInCluster. | |
307 | ||
0d81af58 | 308 | fNBinsECells = nbins1; |
309 | fNBinsPi0Mass = nbins2; | |
310 | fNBinsXNCellsInCluster = nbins3x; | |
311 | fNBinsYNCellsInCluster = nbins3y; | |
312 | fXMaxECells = xmax1; | |
313 | fXMaxPi0Mass = xmax2; | |
314 | fXMaxNCellsInCluster = xmax3; | |
045862db | 315 | } |
316 | ||
317 | //_________________________________________________________________________ | |
318 | void AliCaloCellsQA::Init(Int_t nmods, Int_t det, Int_t startRunNumber, Int_t endRunNumber) | |
319 | { | |
320 | // Class initialization. | |
321 | // Defaults: fClusElowMin = 0.3GeV, fClusEhighMin = 1GeV, fPi0EClusMin = 0.5GeV, fkFullAnalysis = false. | |
322 | ||
323 | // check input (for design limitations only) | |
324 | if (det != kEMCAL && det != kPHOS) { | |
db05d873 | 325 | AliError("Wrong detector provided"); |
326 | AliInfo("I will use EMCAL"); | |
045862db | 327 | det = kEMCAL; |
328 | } | |
329 | if (nmods < 1 || nmods > 10) { | |
db05d873 | 330 | AliError("Wrong last supermodule number + 1 provided"); |
331 | AliInfo("I will use nmods = 10"); | |
045862db | 332 | nmods = 10; |
333 | } | |
334 | ||
335 | fDetector = det; | |
336 | fNMods = nmods; | |
337 | fkFullAnalysis = kFALSE; | |
338 | SetClusterEnergyCuts(); | |
339 | SetBinningParameters(); | |
340 | ||
341 | // minimum/maximum cell absId; | |
342 | // straightforward solution avoids complications | |
343 | fAbsIdMin = 0; | |
344 | fAbsIdMax = fNMods * 1152; | |
345 | if (fDetector == kPHOS) { | |
346 | fAbsIdMin = 1; | |
347 | fAbsIdMax = 1 + (fNMods-1)*3584; | |
348 | } | |
349 | ||
350 | fListOfHistos = new TObjArray; | |
351 | fListOfHistos->SetOwner(kTRUE); | |
352 | ||
353 | fhNEventsProcessedPerRun = new TH1D("hNEventsProcessedPerRun", | |
354 | "Number of processed events vs run number", endRunNumber - startRunNumber, startRunNumber, endRunNumber); | |
355 | fhNEventsProcessedPerRun->SetDirectory(0); | |
356 | fhNEventsProcessedPerRun->SetXTitle("Run number"); | |
357 | fhNEventsProcessedPerRun->SetYTitle("Events"); | |
358 | fListOfHistos->Add(fhNEventsProcessedPerRun); | |
045862db | 359 | } |
360 | ||
361 | //_________________________________________________________________________ | |
d30ac678 | 362 | void AliCaloCellsQA::InitTransientFindCurrentRun(Int_t runNumber) |
045862db | 363 | { |
d30ac678 | 364 | // Initialize transient members, add a new run if necessary. |
045862db | 365 | |
045862db | 366 | // try previous value ... |
d30ac678 | 367 | if (fRI >= 0 && fRunNumbers[fRI] == runNumber) return; |
045862db | 368 | |
369 | // ... or find current run index ... | |
f4dd0897 | 370 | for (fRI = 0; fRI < fNRuns; fRI++) |
d30ac678 | 371 | if (fRunNumbers[fRI] == runNumber) break; |
045862db | 372 | |
373 | // ... or add a new run | |
d30ac678 | 374 | if (fRI == fNRuns) { |
db05d873 | 375 | if (fNRuns >= 1000) AliFatal("Too many runs, how is this possible?"); |
045862db | 376 | |
d30ac678 | 377 | fRunNumbers[fNRuns] = runNumber; |
378 | InitHistosForRun(runNumber); | |
379 | fNRuns++; | |
380 | } | |
045862db | 381 | |
d30ac678 | 382 | // initialize transient class members; |
383 | // this happens once per run, i.e. not a big overhead | |
384 | InitTransientMembers(runNumber); | |
045862db | 385 | } |
386 | ||
387 | //_________________________________________________________________________ | |
d30ac678 | 388 | void AliCaloCellsQA::InitTransientMembers(Int_t run) |
389 | { | |
390 | // Initializes transient data members -- references to histograms | |
391 | // (e.g. in case the class was restored from a file); | |
392 | // run -- current run number. | |
393 | ||
394 | fhNEventsProcessedPerRun = (TH1D*) fListOfHistos->FindObject("hNEventsProcessedPerRun"); | |
395 | ||
396 | fhCellAmplitude = (TH2F*) fListOfHistos->FindObject("hCellAmplitude"); | |
397 | fhCellAmplitudeEhigh = (TH2F*) fListOfHistos->FindObject("hCellAmplitudeEhigh"); | |
398 | fhCellAmplitudeNonLocMax = (TH2F*) fListOfHistos->FindObject("hCellAmplitudeNonLocMax"); | |
399 | fhCellAmplitudeEhighNonLocMax = (TH2F*) fListOfHistos->FindObject("hCellAmplitudeEhighNonLocMax"); | |
400 | fhCellTime = (TH2F*) fListOfHistos->FindObject("hCellTime"); | |
401 | ||
402 | fhCellLocMaxNTimesInClusterElow = (TH1F*) fListOfHistos->FindObject(Form("run%i_hCellLocMaxNTimesInClusterElow",run)); | |
403 | fhCellLocMaxNTimesInClusterEhigh = (TH1F*) fListOfHistos->FindObject(Form("run%i_hCellLocMaxNTimesInClusterEhigh",run)); | |
404 | fhCellLocMaxETotalClusterElow = (TH1F*) fListOfHistos->FindObject(Form("run%i_hCellLocMaxETotalClusterElow",run)); | |
405 | fhCellLocMaxETotalClusterEhigh = (TH1F*) fListOfHistos->FindObject(Form("run%i_hCellLocMaxETotalClusterEhigh",run)); | |
406 | fhCellNonLocMaxNTimesInClusterElow = (TH1F*) fListOfHistos->FindObject(Form("run%i_hCellNonLocMaxNTimesInClusterElow",run)); | |
407 | fhCellNonLocMaxNTimesInClusterEhigh = (TH1F*) fListOfHistos->FindObject(Form("run%i_hCellNonLocMaxNTimesInClusterEhigh",run)); | |
408 | fhCellNonLocMaxETotalClusterElow = (TH1F*) fListOfHistos->FindObject(Form("run%i_hCellNonLocMaxETotalClusterElow",run)); | |
409 | fhCellNonLocMaxETotalClusterEhigh = (TH1F*) fListOfHistos->FindObject(Form("run%i_hCellNonLocMaxETotalClusterEhigh",run)); | |
410 | ||
411 | Int_t minsm = 0; | |
412 | if (fDetector == kPHOS) minsm = 1; | |
413 | ||
414 | // per supermodule histograms | |
415 | for (Int_t sm = minsm; sm < fNMods; sm++) { | |
416 | fhECells[sm] = (TH1F*) fListOfHistos->FindObject(Form("run%i_hECellsSM%i",run,sm)); | |
417 | fhNCellsInCluster[sm] = (TH2F*) fListOfHistos->FindObject(Form("run%i_hNCellsInClusterSM%i",run,sm)); | |
418 | ||
419 | for (Int_t sm2 = sm; sm2 < fNMods; sm2++) | |
420 | fhPi0Mass[sm][sm2] = (TH1F*) fListOfHistos->FindObject(Form("run%i_hPi0MassSM%iSM%i",run,sm,sm2)); | |
421 | } | |
422 | } | |
423 | ||
424 | //_________________________________________________________________________ | |
425 | void AliCaloCellsQA::InitHistosForRun(Int_t run) | |
045862db | 426 | { |
427 | // Initialize per run histograms for a new run number; | |
d30ac678 | 428 | // run -- run number. |
045862db | 429 | |
430 | // do not add histograms to the current directory | |
431 | Bool_t ads = TH1::AddDirectoryStatus(); | |
432 | TH1::AddDirectory(kFALSE); | |
433 | ||
d30ac678 | 434 | fhCellLocMaxNTimesInClusterElow = new TH1F(Form("run%i_hCellLocMaxNTimesInClusterElow",run), |
045862db | 435 | "Number of times cell was local maximum in a low energy cluster", fAbsIdMax-fAbsIdMin,fAbsIdMin,fAbsIdMax); |
d30ac678 | 436 | fhCellLocMaxNTimesInClusterElow->SetXTitle("AbsId"); |
437 | fhCellLocMaxNTimesInClusterElow->SetYTitle("Counts"); | |
045862db | 438 | |
d30ac678 | 439 | fhCellLocMaxNTimesInClusterEhigh = new TH1F(Form("run%i_hCellLocMaxNTimesInClusterEhigh",run), |
045862db | 440 | "Number of times cell was local maximum in a high energy cluster", fAbsIdMax-fAbsIdMin,fAbsIdMin,fAbsIdMax); |
d30ac678 | 441 | fhCellLocMaxNTimesInClusterEhigh->SetXTitle("AbsId"); |
442 | fhCellLocMaxNTimesInClusterEhigh->SetYTitle("Counts"); | |
045862db | 443 | |
d30ac678 | 444 | fhCellLocMaxETotalClusterElow = new TH1F(Form("run%i_hCellLocMaxETotalClusterElow",run), |
045862db | 445 | "Total cluster energy for local maximum cell, low energy", fAbsIdMax-fAbsIdMin,fAbsIdMin,fAbsIdMax); |
d30ac678 | 446 | fhCellLocMaxETotalClusterElow->SetXTitle("AbsId"); |
447 | fhCellLocMaxETotalClusterElow->SetYTitle("Energy"); | |
045862db | 448 | |
d30ac678 | 449 | fhCellLocMaxETotalClusterEhigh = new TH1F(Form("run%i_hCellLocMaxETotalClusterEhigh",run), |
045862db | 450 | "Total cluster energy for local maximum cell, high energy", fAbsIdMax-fAbsIdMin,fAbsIdMin,fAbsIdMax); |
d30ac678 | 451 | fhCellLocMaxETotalClusterEhigh->SetXTitle("AbsId"); |
452 | fhCellLocMaxETotalClusterEhigh->SetYTitle("Energy"); | |
045862db | 453 | |
d30ac678 | 454 | fListOfHistos->Add(fhCellLocMaxNTimesInClusterElow); |
455 | fListOfHistos->Add(fhCellLocMaxNTimesInClusterEhigh); | |
456 | fListOfHistos->Add(fhCellLocMaxETotalClusterElow); | |
457 | fListOfHistos->Add(fhCellLocMaxETotalClusterEhigh); | |
045862db | 458 | |
459 | ||
460 | if (fkFullAnalysis) { | |
d30ac678 | 461 | fhCellNonLocMaxNTimesInClusterElow = new TH1F(Form("run%i_hCellNonLocMaxNTimesInClusterElow",run), |
045862db | 462 | "Number of times cell wasn't local maximum in a low energy cluster", fAbsIdMax-fAbsIdMin,fAbsIdMin,fAbsIdMax); |
d30ac678 | 463 | fhCellNonLocMaxNTimesInClusterElow->SetXTitle("AbsId"); |
464 | fhCellNonLocMaxNTimesInClusterElow->SetYTitle("Counts"); | |
045862db | 465 | |
d30ac678 | 466 | fhCellNonLocMaxNTimesInClusterEhigh = new TH1F(Form("run%i_hCellNonLocMaxNTimesInClusterEhigh",run), |
045862db | 467 | "Number of times cell wasn't local maximum in a high energy cluster", fAbsIdMax-fAbsIdMin,fAbsIdMin,fAbsIdMax); |
d30ac678 | 468 | fhCellNonLocMaxNTimesInClusterEhigh->SetXTitle("AbsId"); |
469 | fhCellNonLocMaxNTimesInClusterEhigh->SetYTitle("Counts"); | |
045862db | 470 | |
d30ac678 | 471 | fhCellNonLocMaxETotalClusterElow = new TH1F(Form("run%i_hCellNonLocMaxETotalClusterElow",run), |
045862db | 472 | "Total cluster energy for not local maximum cell, low energy", fAbsIdMax-fAbsIdMin,fAbsIdMin,fAbsIdMax); |
d30ac678 | 473 | fhCellNonLocMaxETotalClusterElow->SetXTitle("AbsId"); |
474 | fhCellNonLocMaxETotalClusterElow->SetYTitle("Energy"); | |
045862db | 475 | |
d30ac678 | 476 | fhCellNonLocMaxETotalClusterEhigh = new TH1F(Form("run%i_hCellNonLocMaxETotalClusterEhigh",run), |
045862db | 477 | "Total cluster energy for not local maximum cell, high energy", fAbsIdMax-fAbsIdMin,fAbsIdMin,fAbsIdMax); |
d30ac678 | 478 | fhCellNonLocMaxETotalClusterEhigh->SetXTitle("AbsId"); |
479 | fhCellNonLocMaxETotalClusterEhigh->SetYTitle("Energy"); | |
045862db | 480 | |
d30ac678 | 481 | fListOfHistos->Add(fhCellNonLocMaxNTimesInClusterElow); |
482 | fListOfHistos->Add(fhCellNonLocMaxNTimesInClusterEhigh); | |
483 | fListOfHistos->Add(fhCellNonLocMaxETotalClusterElow); | |
484 | fListOfHistos->Add(fhCellNonLocMaxETotalClusterEhigh); | |
045862db | 485 | } |
486 | ||
487 | ||
488 | Int_t minsm = 0; | |
489 | if (fDetector == kPHOS) minsm = 1; | |
490 | ||
491 | // per supermodule histograms | |
492 | for (Int_t sm = minsm; sm < fNMods; sm++) { | |
d30ac678 | 493 | fhECells[sm] = new TH1F(Form("run%i_hECellsSM%i",run,sm), |
0d81af58 | 494 | "Cell amplitude distribution", fNBinsECells,0,fXMaxECells); |
d30ac678 | 495 | fhECells[sm]->SetXTitle("Amplitude, GeV"); |
496 | fhECells[sm]->SetYTitle("Number of cells"); | |
045862db | 497 | |
d30ac678 | 498 | fhNCellsInCluster[sm] = new TH2F(Form("run%i_hNCellsInClusterSM%i",run,sm), |
045862db | 499 | "Distrubution of number of cells in cluster vs cluster energy", |
0d81af58 | 500 | fNBinsXNCellsInCluster,0,fXMaxNCellsInCluster, fNBinsYNCellsInCluster,0,fNBinsYNCellsInCluster); |
d30ac678 | 501 | fhNCellsInCluster[sm]->SetXTitle("Energy, GeV"); |
502 | fhNCellsInCluster[sm]->SetYTitle("Number of cells"); | |
503 | fhNCellsInCluster[sm]->SetZTitle("Counts"); | |
045862db | 504 | |
d30ac678 | 505 | fListOfHistos->Add(fhECells[sm]); |
506 | fListOfHistos->Add(fhNCellsInCluster[sm]); | |
045862db | 507 | } |
508 | ||
509 | // pi0 mass spectrum | |
510 | for (Int_t sm = minsm; sm < fNMods; sm++) | |
511 | for (Int_t sm2 = sm; sm2 < fNMods; sm2++) { | |
d30ac678 | 512 | fhPi0Mass[sm][sm2] = new TH1F(Form("run%i_hPi0MassSM%iSM%i",run,sm,sm2), |
0d81af58 | 513 | "#pi^{0} mass spectrum", fNBinsPi0Mass,0,fXMaxPi0Mass); |
d30ac678 | 514 | fhPi0Mass[sm][sm2]->SetXTitle("M_{#gamma#gamma}, GeV"); |
515 | fhPi0Mass[sm][sm2]->SetYTitle("Counts"); | |
045862db | 516 | |
d30ac678 | 517 | fListOfHistos->Add(fhPi0Mass[sm][sm2]); |
045862db | 518 | } |
519 | ||
520 | // return to the previous add directory status | |
521 | TH1::AddDirectory(ads); | |
522 | } | |
523 | ||
524 | //_________________________________________________________________________ | |
d30ac678 | 525 | void AliCaloCellsQA::FillCellsInCluster(TObjArray *clusArray, AliVCaloCells *cells) |
045862db | 526 | { |
d30ac678 | 527 | // Fill histograms related to a cluster |
045862db | 528 | |
529 | Int_t sm; | |
530 | ||
531 | for (Int_t i = 0; i < clusArray->GetEntriesFast(); i++) | |
532 | { | |
533 | AliVCluster *clus = (AliVCluster*) clusArray->At(i); | |
534 | if ((sm = CheckClusterGetSM(clus)) < 0) continue; | |
535 | ||
d30ac678 | 536 | if (fhNCellsInCluster[sm]) |
537 | fhNCellsInCluster[sm]->Fill(clus->E(), clus->GetNCells()); | |
045862db | 538 | |
539 | if (clus->E() >= fClusElowMin) | |
540 | for (Int_t c = 0; c < clus->GetNCells(); c++) { | |
541 | Int_t absId = clus->GetCellAbsId(c); | |
542 | ||
543 | if (IsCellLocalMaximum(c, clus, cells)) {// local maximum | |
544 | if (clus->E() < fClusEhighMin) { | |
d30ac678 | 545 | fhCellLocMaxNTimesInClusterElow->Fill(absId); |
546 | fhCellLocMaxETotalClusterElow->Fill(absId, clus->E()); | |
045862db | 547 | } else { |
d30ac678 | 548 | fhCellLocMaxNTimesInClusterEhigh->Fill(absId); |
549 | fhCellLocMaxETotalClusterEhigh->Fill(absId, clus->E()); | |
045862db | 550 | } |
551 | } | |
552 | else if (fkFullAnalysis) {// not a local maximum | |
553 | if (clus->E() < fClusEhighMin) { | |
d30ac678 | 554 | fhCellNonLocMaxNTimesInClusterElow->Fill(absId); |
555 | fhCellNonLocMaxETotalClusterElow->Fill(absId, clus->E()); | |
045862db | 556 | } else { |
d30ac678 | 557 | fhCellNonLocMaxNTimesInClusterEhigh->Fill(absId); |
558 | fhCellNonLocMaxETotalClusterEhigh->Fill(absId, clus->E()); | |
045862db | 559 | } |
560 | } | |
561 | } // cells loop | |
562 | ||
563 | } // cluster loop | |
564 | ||
565 | } | |
566 | ||
567 | //_________________________________________________________________________ | |
d30ac678 | 568 | void AliCaloCellsQA::FillJustCells(AliVCaloCells *cells) |
045862db | 569 | { |
d30ac678 | 570 | // Fill cell histograms not related with a cluster |
045862db | 571 | |
572 | Short_t absId; | |
573 | Double_t amp, time; | |
574 | Int_t sm; | |
575 | ||
576 | for (Short_t c = 0; c < cells->GetNumberOfCells(); c++) { | |
577 | cells->GetCell(c, absId, amp, time); | |
578 | if ((sm = GetSM(absId)) < 0) continue; | |
579 | ||
d30ac678 | 580 | if (fhECells[sm]) fhECells[sm]->Fill(amp); |
045862db | 581 | |
0d81af58 | 582 | if (fhCellAmplitude) { // in case InitSummaryHistograms() was not called |
045862db | 583 | fhCellAmplitude->Fill(absId, amp); |
584 | fhCellAmplitudeEhigh->Fill(absId, amp); | |
0d81af58 | 585 | fhCellTime->Fill(absId, time); |
045862db | 586 | |
587 | // fill not a local maximum distributions | |
588 | if (fkFullAnalysis && !IsCellLocalMaximum(absId, cells)) { | |
589 | fhCellAmplitudeNonLocMax->Fill(absId, amp); | |
590 | fhCellAmplitudeEhighNonLocMax->Fill(absId, amp); | |
591 | } | |
592 | } | |
593 | } // cell loop | |
594 | ||
595 | } | |
596 | ||
597 | //_________________________________________________________________________ | |
d30ac678 | 598 | void AliCaloCellsQA::FillPi0Mass(TObjArray *clusArray, Double_t vertexXYZ[3]) |
045862db | 599 | { |
600 | // Fill gamma+gamma invariant mass histograms. | |
601 | // ri -- run index. | |
602 | ||
603 | Int_t sm1, sm2; | |
604 | TLorentzVector p1, p2, psum; | |
605 | ||
606 | // cluster loop | |
607 | for (Int_t i = 0; i < clusArray->GetEntriesFast(); i++) | |
608 | { | |
609 | AliVCluster *clus = (AliVCluster*) clusArray->At(i); | |
610 | if (clus->E() < fPi0EClusMin) continue; | |
611 | if ((sm1 = CheckClusterGetSM(clus)) < 0) continue; | |
612 | ||
613 | clus->GetMomentum(p1, vertexXYZ); | |
614 | ||
615 | // second cluster loop | |
616 | for (Int_t j = i+1; j < clusArray->GetEntriesFast(); j++) | |
617 | { | |
618 | AliVCluster *clus2 = (AliVCluster*) clusArray->At(j); | |
619 | if (clus2->E() < fPi0EClusMin) continue; | |
620 | if ((sm2 = CheckClusterGetSM(clus2)) < 0) continue; | |
621 | ||
622 | clus2->GetMomentum(p2, vertexXYZ); | |
623 | ||
624 | psum = p1 + p2; | |
625 | if (psum.M2() < 0) continue; | |
626 | ||
627 | // s1 <= s2 | |
628 | Int_t s1 = (sm1 <= sm2) ? sm1 : sm2; | |
629 | Int_t s2 = (sm1 <= sm2) ? sm2 : sm1; | |
630 | ||
d30ac678 | 631 | if (fhPi0Mass[s1][s2]) |
632 | fhPi0Mass[s1][s2]->Fill(psum.M()); | |
045862db | 633 | } // second cluster loop |
634 | } // cluster loop | |
635 | } | |
636 | ||
637 | //_________________________________________________________________________ | |
638 | Int_t AliCaloCellsQA::CheckClusterGetSM(AliVCluster* clus) | |
639 | { | |
640 | // Apply common cluster cuts and return supermodule number on success. | |
641 | // Return -1 if cuts not passed or an error occured. | |
642 | ||
643 | // check detector and number of cells in cluster | |
644 | if (fDetector == kEMCAL && !clus->IsEMCAL()) return -1; | |
645 | if (fDetector == kPHOS && !clus->IsPHOS()) return -1; | |
646 | if (clus->GetNCells() < 1) return -1; | |
647 | ||
648 | return GetSM(clus->GetCellAbsId(0)); | |
649 | } | |
650 | ||
651 | //_________________________________________________________________________ | |
652 | Int_t AliCaloCellsQA::GetSM(Int_t absId) | |
653 | { | |
654 | // Convert cell absId --> supermodule number. Return -1 in case of error. | |
655 | // Note: we use simple and straightforward way to find supermodule number. | |
656 | // This allows to avoid unnecessary external dependencies (on geometry). | |
657 | ||
658 | Int_t sm = -1; | |
659 | ||
660 | if (fDetector == kEMCAL) sm = absId/1152; | |
661 | if (fDetector == kPHOS) sm = 1 + (absId-1)/3584; | |
662 | ||
663 | // check for data corruption to avoid segfaults | |
664 | if (sm < 0 || sm > 9) { | |
db05d873 | 665 | AliError("Data corrupted"); |
045862db | 666 | return -1; |
667 | } | |
668 | ||
669 | return sm; | |
670 | } | |
671 | ||
672 | //____________________________________________________________ | |
673 | Bool_t AliCaloCellsQA::IsCellLocalMaximum(Int_t c, AliVCluster* clus, AliVCaloCells* cells) | |
674 | { | |
675 | // Returns true if cell is a local maximum in cluster (clusterizer dependent). | |
676 | // Cell fractions are taken into account. | |
677 | // | |
678 | // c -- cell number in the cluster. | |
679 | ||
680 | Int_t sm, eta, phi, sm2, eta2, phi2; | |
681 | ||
682 | Int_t absId = clus->GetCellAbsId(c); | |
683 | Double_t amp = cells->GetCellAmplitude(absId); | |
684 | if (clus->GetCellAmplitudeFraction(c) > 1e-5) amp *= clus->GetCellAmplitudeFraction(c); | |
685 | ||
686 | AbsIdToSMEtaPhi(absId, sm, eta, phi); | |
687 | ||
688 | // try to find a neighbour which has bigger amplitude | |
689 | for (Int_t c2 = 0; c2 < clus->GetNCells(); c2++) { | |
690 | if (c == c2) continue; | |
691 | ||
692 | Int_t absId2 = clus->GetCellAbsId(c2); | |
693 | Double_t amp2 = cells->GetCellAmplitude(absId2); | |
694 | if (clus->GetCellAmplitudeFraction(c2) > 1e-5) amp2 *= clus->GetCellAmplitudeFraction(c2); | |
695 | ||
696 | AbsIdToSMEtaPhi(absId2, sm2, eta2, phi2); | |
697 | if (sm != sm2) continue; | |
698 | ||
699 | Int_t deta = TMath::Abs(eta-eta2); | |
700 | Int_t dphi = TMath::Abs(phi-phi2); | |
701 | ||
702 | if ((deta == 0 && dphi == 1) || (deta == 1 && dphi == 0) || (deta == 1 && dphi == 1)) | |
703 | if (amp < amp2) return kFALSE; | |
704 | } // cell loop | |
705 | ||
706 | return kTRUE; | |
707 | } | |
708 | ||
709 | //____________________________________________________________ | |
710 | Bool_t AliCaloCellsQA::IsCellLocalMaximum(Int_t absId, AliVCaloCells* cells) | |
711 | { | |
712 | // Returns true if cell is a local maximum among cells (clusterizer independent). | |
713 | ||
714 | Int_t sm, eta, phi, sm2, eta2, phi2; | |
715 | ||
716 | Double_t amp = cells->GetCellAmplitude(absId); | |
717 | AbsIdToSMEtaPhi(absId, sm, eta, phi); | |
718 | ||
719 | // try to find a neighbour which has bigger amplitude | |
720 | for (Short_t c2 = 0; c2 < cells->GetNumberOfCells(); c2++) { | |
721 | Short_t absId2 = cells->GetCellNumber(c2); | |
722 | if (absId2 == absId) continue; | |
723 | ||
724 | AbsIdToSMEtaPhi(absId2, sm2, eta2, phi2); | |
725 | if (sm != sm2) continue; | |
726 | ||
727 | Int_t deta = TMath::Abs(eta-eta2); | |
728 | Int_t dphi = TMath::Abs(phi-phi2); | |
729 | ||
730 | if ((deta == 0 && dphi == 1) || (deta == 1 && dphi == 0) || (deta == 1 && dphi == 1)) | |
731 | if (amp < cells->GetAmplitude(c2)) | |
732 | return kFALSE; | |
733 | } // cell loop | |
734 | ||
735 | return kTRUE; | |
736 | } | |
737 | ||
738 | //____________________________________________________________ | |
739 | void AliCaloCellsQA::AbsIdToSMEtaPhi(Int_t absId, Int_t &sm, Int_t &eta, Int_t &phi) | |
740 | { | |
741 | // Converts absId --> (sm, eta, phi) for a cell. | |
742 | // Works both for EMCAL and for PHOS. | |
743 | // Geometry must be already initialized. | |
744 | ||
745 | // EMCAL | |
746 | if (fDetector == kEMCAL) { | |
f4dd0897 | 747 | AliEMCALGeometry *geomEMCAL = AliEMCALGeometry::GetInstance(); |
045862db | 748 | if (!geomEMCAL) |
db05d873 | 749 | AliFatal("EMCAL geometry is not initialized"); |
045862db | 750 | |
751 | Int_t nModule, nIphi, nIeta; | |
752 | geomEMCAL->GetCellIndex(absId, sm, nModule, nIphi, nIeta); | |
753 | geomEMCAL->GetCellPhiEtaIndexInSModule(sm, nModule, nIphi, nIeta, phi, eta); | |
754 | return; | |
755 | } | |
756 | ||
757 | // PHOS | |
758 | if (fDetector == kPHOS) { | |
f4dd0897 | 759 | AliPHOSGeometry *geomPHOS = AliPHOSGeometry::GetInstance(); |
045862db | 760 | if (!geomPHOS) |
db05d873 | 761 | AliFatal("PHOS geometry is not initialized"); |
045862db | 762 | |
763 | Int_t relid[4]; | |
764 | geomPHOS->AbsToRelNumbering(absId, relid); | |
765 | sm = relid[0]; | |
766 | eta = relid[2]; | |
767 | phi = relid[3]; | |
768 | } | |
769 | ||
770 | // DCAL | |
771 | // not implemented | |
772 | } |