]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TPC/AliTPCdataQA.cxx
D+ cuts updated (Renu)
[u/mrichter/AliRoot.git] / TPC / AliTPCdataQA.cxx
CommitLineData
0ffacf98 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
16
17/* $Id$ */
18
266f8637 19/*
1267cf3a 20 June 2010
21
22 This update should solve two problems mainly:
23 * The vs event histograms have been limited to a fixed size for the
24 DQM. The 500k seemed to be a big size but is no longer so, so we
25 need to dynamically expand the range. The non-trivial point is that
26 we also have to do it for the copy owned by AliTPCQADataMakerRec.
27 * The amoreGui now remembers the x-range of the first visualization so
28 the trick of setting the relevant event range as the histogram is
29 filled no longer works.
30
31 The fix is a bit crude but avoids creating a new histogram. Instead
32 the range is expanded (max events and events per bin is doubled) but
33 the number of bins is kept constant! In this way we can change just
34 double the max of the X-axis of the hist and rebin the data. The
35 same can easily be done for the copy owned by AliTPCQADataMakerRec.
36
37 CAUTION:
38 If we change the number of bins we could crash the whole system
39 because ROOT does not create space for extra bins! (but we do not do
40 this). In that way it is a crude solution.
41 The rebinning in the code only works for an even number of bins.
42
43 In addition to the above a bug in the reading of the config file was
44 also found and corrected. fAdcMax was set instead of fEventsPerBin.
45
46 Finally cout was changes to AliInfo.
47
266f8637 48 February 2008
49
50 The code has been heavily modified so that now the RAW data is
51 "expanded" for each sector and stored in a big signal array. Then a
52 simple version of the code in AliTPCclustererMI is used to identify
53 the local maxima and these are then used for QA. This gives a better
54 estimate of the charge (both max and total) and also limits the
55 effect of noise.
56
57 Implementation:
58
59 In Update the RAW signals >= 3 ADC channels are stored in the arrays.
60
61 There are 3 arrays:
62 Float_t** fAllBins 2d array [row][bin(pad, time)] ADC signal
63 Int_t** fAllSigBins 2d array [row][signal#] bin(with signal)
64 Int_t* fAllNSigBins; 1d array [row] Nsignals
65
66 This is done sector by sector.
67
68 When data from a new sector is encountered, the method
69 FindLocalMaxima is called on the data from the previous sector, and
70 the calibration/data objects are updated with the "cluster"
71 info. Finally the arrays are cleared.
72
73 The requirements for a local maxima is:
74 Charge in bin is >= 5 ADC channels.
75 Charge in bin is larger than all the 8 neighboring bins.
76 At least one of the two pad neighbors has a signal.
77 At least one of the two time neighbors has a signal.
78
79 Before accessing the data it is expected that the Analyse method is
80 called. This normalizes some of the data objects to per event or per
81 cluster.
82 If more data is passed to the class after Analyse has been called
83 the normalization is reversed and Analyse has to be called again.
84*/
85
0ffacf98 86
87//Root includes
88#include <TH1F.h>
89#include <TH2F.h>
90#include <TString.h>
91#include <TMath.h>
0ffacf98 92#include <TDirectory.h>
93#include <TFile.h>
266f8637 94#include <TError.h>
ac940b58 95#include <TMap.h>
6a50ff96 96#include <TProfile.h>
0ffacf98 97//AliRoot includes
98#include "AliRawReader.h"
99#include "AliRawReaderRoot.h"
100#include "AliRawReaderDate.h"
101#include "AliTPCRawStream.h"
0c25417d 102#include "AliTPCRawStreamV3.h"
0ffacf98 103#include "AliTPCCalROC.h"
104#include "AliTPCROC.h"
105#include "AliMathBase.h"
106#include "TTreeStream.h"
0ffacf98 107
108//date
109#include "event.h"
110#include "AliTPCCalPad.h"
258cd111 111#include "AliTPCPreprocessorOnline.h"
0ffacf98 112
113//header file
114#include "AliTPCdataQA.h"
ce0175fa 115#include "AliLog.h"
0ffacf98 116
0ffacf98 117ClassImp(AliTPCdataQA)
118
336156cc 119AliTPCdataQA::AliTPCdataQA() : /*FOLD00*/
0ffacf98 120 fFirstTimeBin(60),
121 fLastTimeBin(1000),
122 fAdcMin(1),
123 fAdcMax(100),
0ffacf98 124 fMapping(NULL),
f11b3071 125 fPedestal(0),
126 fNoise(0),
266f8637 127 fNLocalMaxima(0),
0ffacf98 128 fMaxCharge(0),
f11b3071 129 fMeanCharge(0),
130 fNoThreshold(0),
266f8637 131 fNTimeBins(0),
132 fNPads(0),
133 fTimePosition(0),
c94a79e1 134 fOverThreshold10(0),
135 fOverThreshold20(0),
136 fOverThreshold30(0),
23c9ab21 137 fHistQVsTimeSideA(0),
138 fHistQVsTimeSideC(0),
139 fHistQMaxVsTimeSideA(0),
140 fHistQMaxVsTimeSideC(0),
ce0175fa 141 fHistOccupancyVsEvent(0),
142 fHistNclustersVsEvent(0),
266f8637 143 fEventCounter(0),
144 fIsAnalysed(kFALSE),
ce0175fa 145 fMaxEvents(500000), // Max events for event histograms
146 fEventsPerBin(1000), // Events per bin for event histograms
147 fSignalCounter(0), // Signal counter
148 fClusterCounter(0), // Cluster counter
266f8637 149 fAllBins(0),
150 fAllSigBins(0),
151 fAllNSigBins(0),
152 fRowsMax(0),
153 fPadsMax(0),
154 fTimeBinsMax(0)
0ffacf98 155{
156 //
157 // default constructor
158 //
159}
160
57acd2d2 161//_____________________________________________________________________
162AliTPCdataQA::AliTPCdataQA(AliRecoParam::EventSpecie_t es) :
163fFirstTimeBin(60),
164fLastTimeBin(1000),
165fAdcMin(1),
166fAdcMax(100),
167fMapping(NULL),
168fPedestal(0),
169fNoise(0),
170fNLocalMaxima(0),
171fMaxCharge(0),
172fMeanCharge(0),
173fNoThreshold(0),
174fNTimeBins(0),
175fNPads(0),
176fTimePosition(0),
177fOverThreshold10(0),
178fOverThreshold20(0),
179fOverThreshold30(0),
23c9ab21 180fHistQVsTimeSideA(0),
181fHistQVsTimeSideC(0),
182fHistQMaxVsTimeSideA(0),
183fHistQMaxVsTimeSideC(0),
ce0175fa 184fHistOccupancyVsEvent(0),
185fHistNclustersVsEvent(0),
57acd2d2 186fEventCounter(0),
187fIsAnalysed(kFALSE),
ce0175fa 188fMaxEvents(500000),
189fEventsPerBin(1000),
190fSignalCounter(0),
191fClusterCounter(0),
57acd2d2 192fAllBins(0),
193fAllSigBins(0),
194fAllNSigBins(0),
195fRowsMax(0),
196fPadsMax(0),
197fTimeBinsMax(0)
198{
199// ctor creating the histogram
200 char * name = Form("TPCRAW_%s", AliRecoParam::GetEventSpecieName(es)) ;
201 TH1F(name, name,100,0,100) ;
202}
0ffacf98 203
204//_____________________________________________________________________
205AliTPCdataQA::AliTPCdataQA(const AliTPCdataQA &ped) : /*FOLD00*/
336156cc 206 TH1F(ped),
0ffacf98 207 fFirstTimeBin(ped.GetFirstTimeBin()),
208 fLastTimeBin(ped.GetLastTimeBin()),
209 fAdcMin(ped.GetAdcMin()),
210 fAdcMax(ped.GetAdcMax()),
266f8637 211 fMapping(NULL),
212 fPedestal(0),
213 fNoise(0),
214 fNLocalMaxima(0),
215 fMaxCharge(0),
216 fMeanCharge(0),
217 fNoThreshold(0),
266f8637 218 fNTimeBins(0),
219 fNPads(0),
220 fTimePosition(0),
c94a79e1 221 fOverThreshold10(0),
222 fOverThreshold20(0),
223 fOverThreshold30(0),
23c9ab21 224 fHistQVsTimeSideA(0),
225 fHistQVsTimeSideC(0),
226 fHistQMaxVsTimeSideA(0),
227 fHistQMaxVsTimeSideC(0),
ce0175fa 228 fHistOccupancyVsEvent(0),
229 fHistNclustersVsEvent(0),
266f8637 230 fEventCounter(ped.GetEventCounter()),
231 fIsAnalysed(ped.GetIsAnalysed()),
ce0175fa 232 fMaxEvents(ped.GetMaxEvents()),
233 fEventsPerBin(ped.GetEventsPerBin()),
234 fSignalCounter(ped.GetSignalCounter()),
235 fClusterCounter(ped.GetClusterCounter()),
266f8637 236 fAllBins(0),
237 fAllSigBins(0),
238 fAllNSigBins(0),
239 fRowsMax(0),
240 fPadsMax(0),
241 fTimeBinsMax(0)
0ffacf98 242{
243 //
244 // copy constructor
245 //
266f8637 246 if(ped.GetNLocalMaxima())
247 fNLocalMaxima = new AliTPCCalPad(*ped.GetNLocalMaxima());
248 if(ped.GetMaxCharge())
249 fMaxCharge = new AliTPCCalPad(*ped.GetMaxCharge());
250 if(ped.GetMeanCharge())
251 fMeanCharge = new AliTPCCalPad(*ped.GetMeanCharge());
252 if(ped.GetNoThreshold())
253 fNoThreshold = new AliTPCCalPad(*ped.GetNoThreshold());
254 if(ped.GetNTimeBins())
255 fNTimeBins = new AliTPCCalPad(*ped.GetNTimeBins());
256 if(ped.GetNPads())
257 fNPads = new AliTPCCalPad(*ped.GetNPads());
258 if(ped.GetTimePosition())
259 fTimePosition = new AliTPCCalPad(*ped.GetTimePosition());
260 if(ped.GetOverThreshold10())
261 fOverThreshold10 = new AliTPCCalPad(*ped.GetOverThreshold10());
262 if(ped.GetOverThreshold20())
263 fOverThreshold20 = new AliTPCCalPad(*ped.GetOverThreshold20());
264 if(ped.GetOverThreshold30())
265 fOverThreshold30 = new AliTPCCalPad(*ped.GetOverThreshold30());
1cb9ffdb 266 if(ped.GetHistQVsTimeSideA()) {
23c9ab21 267 fHistQVsTimeSideA = new TProfile(*ped.GetHistQVsTimeSideA());
1cb9ffdb 268 fHistQVsTimeSideA->SetDirectory(0);
269 }
270 if(ped.GetHistQVsTimeSideC()) {
23c9ab21 271 fHistQVsTimeSideC = new TProfile(*ped.GetHistQVsTimeSideC());
1cb9ffdb 272 fHistQVsTimeSideC->SetDirectory(0);
273 }
274 if(ped.GetHistQMaxVsTimeSideA()) {
23c9ab21 275 fHistQMaxVsTimeSideA = new TProfile(*ped.GetHistQMaxVsTimeSideA());
1cb9ffdb 276 fHistQMaxVsTimeSideA->SetDirectory(0);
277 }
278 if(ped.GetHistQMaxVsTimeSideC()) {
23c9ab21 279 fHistQMaxVsTimeSideC = new TProfile(*ped.GetHistQMaxVsTimeSideC());
1cb9ffdb 280 fHistQMaxVsTimeSideC->SetDirectory(0);
281 }
ce0175fa 282 if(ped.GetHistOccupancyVsEventConst()) {
283 fHistOccupancyVsEvent = new TH1F(*ped.GetHistOccupancyVsEventConst());
284 fHistOccupancyVsEvent->SetDirectory(0);
285 }
286 if(ped.GetHistNclustersVsEventConst()) {
287 fHistNclustersVsEvent = new TH1F(*ped.GetHistNclustersVsEventConst());
288 fHistNclustersVsEvent->SetDirectory(0);
289 }
0ffacf98 290}
291
ac940b58 292//_____________________________________________________________________
8ba97cc8 293AliTPCdataQA::AliTPCdataQA(const TMap *config) : /*FOLD00*/
ac940b58 294 TH1F("TPCRAW","TPCRAW",100,0,100),
295 fFirstTimeBin(60),
296 fLastTimeBin(1000),
297 fAdcMin(1),
298 fAdcMax(100),
299 fMapping(NULL),
300 fPedestal(0),
301 fNoise(0),
302 fNLocalMaxima(0),
303 fMaxCharge(0),
304 fMeanCharge(0),
305 fNoThreshold(0),
306 fNTimeBins(0),
307 fNPads(0),
308 fTimePosition(0),
309 fOverThreshold10(0),
310 fOverThreshold20(0),
311 fOverThreshold30(0),
23c9ab21 312 fHistQVsTimeSideA(0),
313 fHistQVsTimeSideC(0),
314 fHistQMaxVsTimeSideA(0),
315 fHistQMaxVsTimeSideC(0),
ce0175fa 316 fHistOccupancyVsEvent(0),
317 fHistNclustersVsEvent(0),
ac940b58 318 fEventCounter(0),
319 fIsAnalysed(kFALSE),
ce0175fa 320 fMaxEvents(500000),
321 fEventsPerBin(1000),
322 fSignalCounter(0),
323 fClusterCounter(0),
ac940b58 324 fAllBins(0),
325 fAllSigBins(0),
326 fAllNSigBins(0),
327 fRowsMax(0),
328 fPadsMax(0),
329 fTimeBinsMax(0)
330{
331 //
332 // default constructor
333 //
334 if (config->GetValue("FirstTimeBin")) fFirstTimeBin = ((TObjString*)config->GetValue("FirstTimeBin"))->GetString().Atoi();
335 if (config->GetValue("LastTimeBin")) fLastTimeBin = ((TObjString*)config->GetValue("LastTimeBin"))->GetString().Atoi();
336 if (config->GetValue("AdcMin")) fAdcMin = ((TObjString*)config->GetValue("AdcMin"))->GetString().Atoi();
337 if (config->GetValue("AdcMax")) fAdcMax = ((TObjString*)config->GetValue("AdcMax"))->GetString().Atoi();
ce0175fa 338 if (config->GetValue("MaxEvents")) fMaxEvents = ((TObjString*)config->GetValue("MaxEvents"))->GetString().Atoi();
1267cf3a 339 if (config->GetValue("EventsPerBin")) fEventsPerBin = ((TObjString*)config->GetValue("EventsPerBin"))->GetString().Atoi();
ac940b58 340}
0ffacf98 341
342//_____________________________________________________________________
343AliTPCdataQA& AliTPCdataQA::operator = (const AliTPCdataQA &source)
344{
345 //
346 // assignment operator
347 //
348 if (&source == this) return *this;
349 new (this) AliTPCdataQA(source);
350
351 return *this;
352}
353
354
355//_____________________________________________________________________
356AliTPCdataQA::~AliTPCdataQA() /*FOLD00*/
357{
358 //
359 // destructor
360 //
361
362 // do not delete fMapping, because we do not own it.
266f8637 363 // do not delete fMapping, because we do not own it.
364 // do not delete fNoise and fPedestal, because we do not own them.
365
366 delete fNLocalMaxima;
367 delete fMaxCharge;
368 delete fMeanCharge;
369 delete fNoThreshold;
370 delete fNTimeBins;
371 delete fNPads;
372 delete fTimePosition;
373 delete fOverThreshold10;
374 delete fOverThreshold20;
375 delete fOverThreshold30;
23c9ab21 376 delete fHistQVsTimeSideA;
377 delete fHistQVsTimeSideC;
378 delete fHistQMaxVsTimeSideA;
379 delete fHistQMaxVsTimeSideC;
ce0175fa 380 delete fHistOccupancyVsEvent;
381 delete fHistNclustersVsEvent;
266f8637 382
383 for (Int_t iRow = 0; iRow < fRowsMax; iRow++) {
384 delete [] fAllBins[iRow];
385 delete [] fAllSigBins[iRow];
386 }
387 delete [] fAllBins;
388 delete [] fAllSigBins;
389 delete [] fAllNSigBins;
0ffacf98 390}
ce0175fa 391
392//_____________________________________________________________________
393TH1F* AliTPCdataQA::GetHistOccupancyVsEvent()
394{
395 //
396 // Create Occupancy vs event histogram
397 // (we create this histogram differently then the other histograms
398 // because this we want to be able to access and copy
399 // from AliTPCQAMakerRec before it normally would be created)
400 //
401 if(!fHistOccupancyVsEvent) {
402
403 Int_t nBins = fMaxEvents/fEventsPerBin;
404 fHistOccupancyVsEvent = new TH1F("hOccupancyVsEvent", "Occupancy vs event number (~time); Event number; Occupancy", nBins, 0, nBins*fEventsPerBin);
405 fHistOccupancyVsEvent->SetDirectory(0);
ce0175fa 406 }
407
408 return fHistOccupancyVsEvent;
409}
410
411//_____________________________________________________________________
412TH1F* AliTPCdataQA::GetHistNclustersVsEvent()
413{
414 //
415 // Create Nclusters vs event histogram
416 // (we create this histogram differently then the other histograms
417 // because this we want to be able to access and copy
418 // from AliTPCQAMakerRec before it normally would be created)
419 //
420 if(!fHistNclustersVsEvent) {
421
422 Int_t nBins = fMaxEvents/fEventsPerBin;
423 fHistNclustersVsEvent = new TH1F("hNclustersVsEvent", "Nclusters vs event number (~time); Event number; Nclusters per event", nBins, 0, nBins*fEventsPerBin);
424 fHistNclustersVsEvent->SetDirectory(0);
ce0175fa 425 }
426
427 return fHistNclustersVsEvent;
428}
429
430//_____________________________________________________________________
431void AliTPCdataQA::UpdateEventHistograms()
432{
433 // Update histograms that display occupancy and
434 // number of clusters as a function of number of
435 // events
436 if (!fHistOccupancyVsEvent)
437 GetHistOccupancyVsEvent();
438 if (!fHistNclustersVsEvent)
439 GetHistNclustersVsEvent();
440
1267cf3a 441 if(fEventCounter > fMaxEvents) {
442
443 // we have to expand the histogram to handle the larger number of
444 // events. The way it is done now is to double the range and the
445 // number of events per bin (so the number of histogram bins stays
446 // constant)
447 fEventsPerBin *= 2;
448 fMaxEvents *= 2;
449
450 // Change histogram limits
451 const Int_t nBins = fHistOccupancyVsEvent->GetXaxis()->GetNbins();
452 fHistOccupancyVsEvent->GetXaxis()->Set(nBins, fHistOccupancyVsEvent->GetXaxis()->GetNbins(), fMaxEvents);
453 fHistNclustersVsEvent->GetXaxis()->Set(nBins, fHistNclustersVsEvent->GetXaxis()->GetNbins(), fMaxEvents);
454
455 // Rebin the histogram
456 for(Int_t bin = 1; bin <= nBins; bin+=2) {
457
458 Int_t newBin = TMath::Nint(Float_t(bin+1)/2.0);
459 Float_t newContent = (fHistOccupancyVsEvent->GetBinContent(bin)+
460 fHistOccupancyVsEvent->GetBinContent(bin+1))/2.0;
461 fHistOccupancyVsEvent->SetBinContent(newBin, newContent);
462
463 newContent = (fHistNclustersVsEvent->GetBinContent(bin)+
464 fHistNclustersVsEvent->GetBinContent(bin+1))/2.0;
465 fHistNclustersVsEvent->SetBinContent(newBin, newContent);
466 }
467
468 // Set the remaining bins to 0
469 Int_t lastHalf = nBins/2 +1;
470 for(Int_t bin = lastHalf; bin <= nBins; bin++) {
471
472 fHistOccupancyVsEvent->SetBinContent(bin, 0);
473 fHistNclustersVsEvent->SetBinContent(bin, 0);
474 }
475
476 // In this case we should nut update but wait untill the new
477 // number of events per bin is reached!
478 return;
479 }
480
481 const Int_t bin = TMath::Nint(Float_t(fEventCounter)/fEventsPerBin);
482
ce0175fa 483 Float_t averageOccupancy =
484 Float_t(fSignalCounter)/fEventsPerBin/(fLastTimeBin - fFirstTimeBin +1.0)
1267cf3a 485 / 570132.0; // 570,132 is number of pads
486 fHistOccupancyVsEvent->SetBinContent(bin, averageOccupancy);
ce0175fa 487 fSignalCounter = 0;
488
489 Float_t averageNclusters =
490 Float_t(fClusterCounter)/fEventsPerBin;
1267cf3a 491 fHistNclustersVsEvent->SetBinContent(bin, averageNclusters);
ce0175fa 492 fClusterCounter = 0;
493}
494
0c25417d 495//_____________________________________________________________________
6a50ff96 496Bool_t AliTPCdataQA::ProcessEvent(AliTPCRawStreamV3 *const rawStreamV3)
0c25417d 497{
498 //
499 // Event Processing loop - AliTPCRawStreamV3
500 //
501 Bool_t withInput = kFALSE;
502 Int_t nSignals = 0;
503 Int_t lastSector = -1;
504
505 while ( rawStreamV3->NextDDL() ){
1267cf3a 506
0c25417d 507 while ( rawStreamV3->NextChannel() ){
1267cf3a 508
0c25417d 509 Int_t iSector = rawStreamV3->GetSector(); // current sector
510 Int_t iRow = rawStreamV3->GetRow(); // current row
511 Int_t iPad = rawStreamV3->GetPad(); // current pad
512 if (iRow<0 || iPad<0) continue;
513 // Call local maxima finder if the data is in a new sector
514 if(iSector != lastSector) {
515
516 if(nSignals>0)
517 FindLocalMaxima(lastSector);
518
519 CleanArrays();
520 lastSector = iSector;
521 nSignals = 0;
522 }
523
524 while ( rawStreamV3->NextBunch() ){
1267cf3a 525
0c25417d 526 Int_t startTbin = (Int_t)rawStreamV3->GetStartTimeBin();
527 Int_t bunchlength = (Int_t)rawStreamV3->GetBunchLength();
528 const UShort_t *sig = rawStreamV3->GetSignals();
529
530 for (Int_t iTimeBin = 0; iTimeBin<bunchlength; iTimeBin++){
531 Float_t signal=(Float_t)sig[iTimeBin];
532 nSignals += Update(iSector,iRow,iPad,startTbin--,signal);
533 withInput = kTRUE;
534 }
535 }
536 }
537 }
1267cf3a 538
539 if (lastSector>=0&&nSignals>0)
0c25417d 540 FindLocalMaxima(lastSector);
541
542 return withInput;
543}
544
545//_____________________________________________________________________
6a50ff96 546Bool_t AliTPCdataQA::ProcessEvent(AliRawReader *const rawReader)
0c25417d 547{
548 //
549 // Event processing loop - AliRawReader
550 //
c75ba816 551 AliTPCRawStreamV3 rawStreamV3(rawReader, (AliAltroMapping**)fMapping);
552 Bool_t res=ProcessEvent(&rawStreamV3);
ce0175fa 553 if(res) {
0c25417d 554 fEventCounter++; // only increment event counter if there is TPC data
ce0175fa 555
1267cf3a 556 if(fEventCounter%fEventsPerBin==0)
ce0175fa 557 UpdateEventHistograms();
558 }
0c25417d 559 return res;
560}
0ffacf98 561
0ffacf98 562//_____________________________________________________________________
6a50ff96 563Bool_t AliTPCdataQA::ProcessEvent(AliTPCRawStream *const rawStream)
0ffacf98 564{
565 //
566 // Event Processing loop - AliTPCRawStream
567 //
568
0ffacf98 569 Bool_t withInput = kFALSE;
266f8637 570 Int_t nSignals = 0;
571 Int_t lastSector = -1;
0ffacf98 572
573 while (rawStream->Next()) {
574
575 Int_t iSector = rawStream->GetSector(); // current ROC
576 Int_t iRow = rawStream->GetRow(); // current row
577 Int_t iPad = rawStream->GetPad(); // current pad
578 Int_t iTimeBin = rawStream->GetTime(); // current time bin
579 Float_t signal = rawStream->GetSignal(); // current ADC signal
0c25417d 580
266f8637 581 // Call local maxima finder if the data is in a new sector
582 if(iSector != lastSector) {
583
584 if(nSignals>0)
0c25417d 585 FindLocalMaxima(lastSector);
266f8637 586
587 CleanArrays();
588 lastSector = iSector;
589 nSignals = 0;
590 }
0ffacf98 591
266f8637 592 // Sometimes iRow==-1 if there is problems to read the data
593 if(iRow>=0) {
594 nSignals += Update(iSector,iRow,iPad,iTimeBin,signal);
595 withInput = kTRUE;
596 }
0ffacf98 597 }
598
0c25417d 599 if (lastSector>=0&&nSignals>0)
600 FindLocalMaxima(lastSector);
601
0ffacf98 602 return withInput;
603}
604
605
606//_____________________________________________________________________
6a50ff96 607Bool_t AliTPCdataQA::ProcessEventOld(AliRawReader *const rawReader)
0ffacf98 608{
609 //
610 // Event processing loop - AliRawReader
611 //
612
613 // if fMapping is NULL the rawstream will crate its own mapping
614 AliTPCRawStream rawStream(rawReader, (AliAltroMapping**)fMapping);
615 rawReader->Select("TPC");
c75bf2f1 616 Bool_t res = ProcessEvent(&rawStream);
617
618 if(res)
619 fEventCounter++; // only increment event counter if there is TPC data
620 // otherwise Analyse (called in QA) fails
621 return res;
0ffacf98 622}
623
624
625//_____________________________________________________________________
6a50ff96 626Bool_t AliTPCdataQA::ProcessEvent(eventHeaderStruct *const event)
0ffacf98 627{
628 //
629 // process date event
630 //
631
c75ba816 632 AliRawReaderDate rawReader((void*)event);
633 Bool_t result=ProcessEvent(&rawReader);
0ffacf98 634 return result;
635}
636
637
638
639//_____________________________________________________________________
640void AliTPCdataQA::DumpToFile(const Char_t *filename, const Char_t *dir, Bool_t append) /*FOLD00*/
641{
642 //
643 // Write class to file
644 //
645
646 TString sDir(dir);
647 TString option;
648
649 if ( append )
650 option = "update";
651 else
652 option = "recreate";
653
654 TDirectory *backup = gDirectory;
655 TFile f(filename,option.Data());
656 f.cd();
657 if ( !sDir.IsNull() ){
658 f.mkdir(sDir.Data());
659 f.cd(sDir);
660 }
661 this->Write();
662 f.Close();
663
664 if ( backup ) backup->cd();
665}
666
667
668//_____________________________________________________________________
266f8637 669Int_t AliTPCdataQA::Update(const Int_t iSector, /*FOLD00*/
670 const Int_t iRow,
671 const Int_t iPad,
672 const Int_t iTimeBin,
673 Float_t signal)
0ffacf98 674{
675 //
676 // Signal filling method
677 //
266f8637 678
679 //
680 // Define the calibration objects the first time Update is called
681 // NB! This has to be done first even if the data is rejected by the time
682 // cut to make sure that the objects are available in Analyse
683 //
684 if (!fNLocalMaxima) fNLocalMaxima = new AliTPCCalPad("NLocalMaxima","NLocalMaxima");
0ffacf98 685 if (!fMaxCharge) fMaxCharge = new AliTPCCalPad("MaxCharge","MaxCharge");
f11b3071 686 if (!fMeanCharge) fMeanCharge = new AliTPCCalPad("MeanCharge","MeanCharge");
687 if (!fNoThreshold) fNoThreshold = new AliTPCCalPad("NoThreshold","NoThreshold");
266f8637 688 if (!fNTimeBins) fNTimeBins = new AliTPCCalPad("NTimeBins","NTimeBins");
689 if (!fNPads) fNPads = new AliTPCCalPad("NPads","NPads");
690 if (!fTimePosition) fTimePosition = new AliTPCCalPad("TimePosition","TimePosition");
0ffacf98 691 if (!fOverThreshold10) fOverThreshold10 = new AliTPCCalPad("OverThreshold10","OverThreshold10");
692 if (!fOverThreshold20) fOverThreshold20 = new AliTPCCalPad("OverThreshold20","OverThreshold20");
693 if (!fOverThreshold30) fOverThreshold30 = new AliTPCCalPad("OverThreshold30","OverThreshold30");
1cb9ffdb 694 if (!fHistQVsTimeSideA) {
23c9ab21 695 fHistQVsTimeSideA = new TProfile("hQVsTimeSideA", "Q vs time (side A); Time [Timebin]; Q [ADC ch]", 100, 0, 1000);
1cb9ffdb 696 fHistQVsTimeSideA->SetDirectory(0);
697 }
698 if (!fHistQVsTimeSideC) {
23c9ab21 699 fHistQVsTimeSideC = new TProfile("hQVsTimeSideC", "Q vs time (side C); Time [Timebin]; Q [ADC ch]", 100, 0, 1000);
1cb9ffdb 700 fHistQVsTimeSideC->SetDirectory(0);
701 }
702 if (!fHistQMaxVsTimeSideA) {
23c9ab21 703 fHistQMaxVsTimeSideA = new TProfile("hQMaxVsTimeSideA", "Q_{MAX} vs time (side A); Time [Timebin]; Q_{MAX} [ADC ch]", 100, 0, 1000);
1cb9ffdb 704 fHistQMaxVsTimeSideA->SetDirectory(0);
705 }
706 if (!fHistQMaxVsTimeSideC) {
23c9ab21 707 fHistQMaxVsTimeSideC = new TProfile("hQMaxVsTimeSideC", "Q_{MAX} vs time (side C); Time [Timebin]; Q_{MAX} [ADC ch]", 100, 0, 1000);
1cb9ffdb 708 fHistQMaxVsTimeSideC->SetDirectory(0);
709 }
710
266f8637 711 // Make the arrays for expanding the data
f11b3071 712
266f8637 713 if (!fAllBins)
714 MakeArrays();
715
716 //
717 // If Analyse has been previously called we need now to denormalize the data
718 // as more data is coming
719 //
720 if(fIsAnalysed == kTRUE) {
721
722 const Int_t nTimeBins = fLastTimeBin - fFirstTimeBin +1;
723 const Float_t denormalization = Float_t(fEventCounter * nTimeBins);
724 fNoThreshold->Multiply(denormalization);
725
726 fMeanCharge->Multiply(fNLocalMaxima);
727 fMaxCharge->Multiply(fNLocalMaxima);
728 fNTimeBins->Multiply(fNLocalMaxima);
729 fNPads->Multiply(fNLocalMaxima);
730 fTimePosition->Multiply(fNLocalMaxima);
731 fIsAnalysed = kFALSE;
732 }
f11b3071 733
266f8637 734 //
735 // TimeBin cut
736 //
737 if (iTimeBin<fFirstTimeBin) return 0;
738 if (iTimeBin>fLastTimeBin) return 0;
739
f11b3071 740 // if pedestal calibrations are loaded subtract pedestals
741 if(fPedestal) {
742
266f8637 743 Float_t ped = fPedestal->GetCalROC(iSector)->GetValue(iRow, iPad);
744 // Don't use data from pads where pedestals are abnormally small or large
745 if(ped<10 || ped>90)
f11b3071 746 return 0;
266f8637 747 signal -= ped;
f11b3071 748 }
266f8637 749
750 // In fNoThreshold we fill all data to estimate the ZS volume
751 Float_t count = fNoThreshold->GetCalROC(iSector)->GetValue(iRow, iPad);
752 fNoThreshold->GetCalROC(iSector)->SetValue(iRow, iPad,count+1);
753
f11b3071 754 // Require at least 3 ADC channels
266f8637 755 if (signal < 3.0)
f11b3071 756 return 0;
757
758 // if noise calibrations are loaded require at least 3*sigmaNoise
759 if(fNoise) {
266f8637 760
761 Float_t noise = fNoise->GetCalROC(iSector)->GetValue(iRow, iPad);
762
763 if(signal < noise*3.0)
f11b3071 764 return 0;
765 }
266f8637 766
f11b3071 767 //
266f8637 768 // This signal is ok and we store it in the cluster map
0ffacf98 769 //
f11b3071 770
266f8637 771 SetExpandDigit(iRow, iPad, iTimeBin, signal);
ce0175fa 772
773 fSignalCounter++;
266f8637 774
775 return 1; // signal was accepted
f11b3071 776}
266f8637 777
f11b3071 778//_____________________________________________________________________
266f8637 779void AliTPCdataQA::FindLocalMaxima(const Int_t iSector)
f11b3071 780{
781 //
266f8637 782 // This method is called after the data from each sector has been
783 // exapanded into an array
784 // Loop over the signals and identify local maxima and fill the
785 // calibration objects with the information
f11b3071 786 //
266f8637 787
788 Int_t nLocalMaxima = 0;
789 const Int_t maxTimeBin = fTimeBinsMax+4; // Used to step between neighboring pads
790 // Because we have tha pad-time data in a
791 // 1d array
792
793 for (Int_t iRow = 0; iRow < fRowsMax; iRow++) {
794
795 Float_t* allBins = fAllBins[iRow];
796 Int_t* sigBins = fAllSigBins[iRow];
797 const Int_t nSigBins = fAllNSigBins[iRow];
798
799 for (Int_t iSig = 0; iSig < nSigBins; iSig++) {
800
801 Int_t bin = sigBins[iSig];
802 Float_t *b = &allBins[bin];
803
804 //
805 // Now we check if this is a local maximum
806 //
807
808 Float_t qMax = b[0];
809
810 // First check that the charge is bigger than the threshold
811 if (qMax<5)
812 continue;
813
814 // Require at least one neighboring pad with signal
815 if (b[-maxTimeBin]+b[maxTimeBin]<=0) continue;
816
817 // Require at least one neighboring time bin with signal
818 if (b[-1]+b[1]<=0) continue;
819
820 //
821 // Check that this is a local maximum
822 // Note that the checking is done so that if 2 charges has the same
823 // qMax then only 1 cluster is generated
824 // (that is why there is BOTH > and >=)
825 //
826 if (b[-maxTimeBin] >= qMax) continue;
827 if (b[-1 ] >= qMax) continue;
828 if (b[+maxTimeBin] > qMax) continue;
829 if (b[+1 ] > qMax) continue;
830 if (b[-maxTimeBin-1] >= qMax) continue;
831 if (b[+maxTimeBin-1] >= qMax) continue;
832 if (b[+maxTimeBin+1] > qMax) continue;
833 if (b[-maxTimeBin+1] >= qMax) continue;
834
835 //
836 // Now we accept the local maximum and fill the calibration/data objects
837 //
838 nLocalMaxima++;
839
840 Int_t iPad, iTimeBin;
841 GetPadAndTimeBin(bin, iPad, iTimeBin);
842
843 Float_t count = fNLocalMaxima->GetCalROC(iSector)->GetValue(iRow, iPad);
844 fNLocalMaxima->GetCalROC(iSector)->SetValue(iRow, iPad, count+1);
845
846 count = fTimePosition->GetCalROC(iSector)->GetValue(iRow, iPad);
847 fTimePosition->GetCalROC(iSector)->SetValue(iRow, iPad, count+iTimeBin);
848
849 Float_t charge = fMaxCharge->GetCalROC(iSector)->GetValue(iRow, iPad);
850 fMaxCharge->GetCalROC(iSector)->SetValue(iRow, iPad, charge + qMax);
851
852 if(qMax>=10) {
853 count = fOverThreshold10->GetCalROC(iSector)->GetValue(iRow, iPad);
854 fOverThreshold10->GetCalROC(iSector)->SetValue(iRow, iPad, count+1);
855 }
856 if(qMax>=20) {
857 count = fOverThreshold20->GetCalROC(iSector)->GetValue(iRow, iPad);
858 fOverThreshold20->GetCalROC(iSector)->SetValue(iRow, iPad, count+1);
859 }
860 if(qMax>=30) {
861 count = fOverThreshold30->GetCalROC(iSector)->GetValue(iRow, iPad);
862 fOverThreshold30->GetCalROC(iSector)->SetValue(iRow, iPad, count+1);
863 }
864
865 //
866 // Calculate the total charge as the sum over the region:
867 //
868 // o o o o o
869 // o i i i o
870 // o i C i o
871 // o i i i o
872 // o o o o o
873 //
874 // with qmax at the center C.
875 //
876 // The inner charge (i) we always add, but we only add the outer
877 // charge (o) if the neighboring inner bin (i) has a signal.
878 //
879 Int_t minP = 0, maxP = 0, minT = 0, maxT = 0;
880 Float_t qTot = qMax;
881 for(Int_t i = -1; i<=1; i++) {
882 for(Int_t j = -1; j<=1; j++) {
883
884 if(i==0 && j==0)
885 continue;
886
77f88633 887 Float_t charge1 = GetQ(b, i, j, maxTimeBin, minT, maxT, minP, maxP);
888 qTot += charge1;
889 if(charge1>0) {
266f8637 890 // see if the next neighbor is also above threshold
891 if(i*j==0) {
892 qTot += GetQ(b, 2*i, 2*j, maxTimeBin, minT, maxT, minP, maxP);
893 } else {
894 // we are in a diagonal corner
895 qTot += GetQ(b, i, 2*j, maxTimeBin, minT, maxT, minP, maxP);
896 qTot += GetQ(b, 2*i, j, maxTimeBin, minT, maxT, minP, maxP);
897 qTot += GetQ(b, 2*i, 2*j, maxTimeBin, minT, maxT, minP, maxP);
898 }
899 }
900 }
901 }
902
903 charge = fMeanCharge->GetCalROC(iSector)->GetValue(iRow, iPad);
904 fMeanCharge->GetCalROC(iSector)->SetValue(iRow, iPad, charge + qTot);
905
906 count = fNTimeBins->GetCalROC(iSector)->GetValue(iRow, iPad);
907 fNTimeBins->GetCalROC(iSector)->SetValue(iRow, iPad, count + maxT-minT+1);
908
909 count = fNPads->GetCalROC(iSector)->GetValue(iRow, iPad);
910 fNPads->GetCalROC(iSector)->SetValue(iRow, iPad, count + maxP-minP+1);
911
23c9ab21 912 if((iSector%36)<18) { // A side
913 fHistQVsTimeSideA->Fill(iTimeBin, qTot);
914 fHistQMaxVsTimeSideA->Fill(iTimeBin, qMax);
915 } else {
916 fHistQVsTimeSideC->Fill(iTimeBin, qTot);
917 fHistQMaxVsTimeSideC->Fill(iTimeBin, qMax);
918 }
266f8637 919 } // end loop over signals
920 } // end loop over rows
f11b3071 921
ce0175fa 922 fClusterCounter += nLocalMaxima;
0ffacf98 923}
11ccf1c1 924
f11b3071 925//_____________________________________________________________________
926void AliTPCdataQA::Analyse()
927{
11ccf1c1 928 //
f11b3071 929 // Calculate calibration constants
11ccf1c1 930 //
f11b3071 931
1267cf3a 932 AliInfo("Analyse called");
f11b3071 933
266f8637 934 if(fIsAnalysed == kTRUE) {
935
1267cf3a 936 AliInfo("No new data since Analyse was called last time");
266f8637 937 return;
938 }
f11b3071 939
266f8637 940 if(fEventCounter==0) {
941
1267cf3a 942 AliInfo("EventCounter == 0, Cannot analyse");
f11b3071 943 return;
944 }
266f8637 945
f11b3071 946 Int_t nTimeBins = fLastTimeBin - fFirstTimeBin +1;
1267cf3a 947 AliInfo(Form("EventCounter: %d , TimeBins: %d", fEventCounter, nTimeBins));
f11b3071 948
f11b3071 949 Float_t normalization = 1.0 / Float_t(fEventCounter * nTimeBins);
266f8637 950 fNoThreshold->Multiply(normalization);
951
952 fMeanCharge->Divide(fNLocalMaxima);
953 fMaxCharge->Divide(fNLocalMaxima);
954 fNTimeBins->Divide(fNLocalMaxima);
955 fNPads->Divide(fNLocalMaxima);
956 fTimePosition->Divide(fNLocalMaxima);
957
958 fIsAnalysed = kTRUE;
11ccf1c1 959}
258cd111 960
961
266f8637 962//_____________________________________________________________________
6a50ff96 963void AliTPCdataQA::MakeTree(const char *fname) const {
258cd111 964 //
965 // Export result to the tree -located in the file
966 // This file can be analyzed using AliTPCCalibViewer
967 //
258cd111 968 AliTPCPreprocessorOnline preprocesor;
266f8637 969
970 if (fNLocalMaxima) preprocesor.AddComponent(fNLocalMaxima);
971 if (fMaxCharge) preprocesor.AddComponent(fMaxCharge);
972 if (fMeanCharge) preprocesor.AddComponent(fMeanCharge);
973 if (fNoThreshold) preprocesor.AddComponent(fNoThreshold);
974 if (fNTimeBins) preprocesor.AddComponent(fNTimeBins);
975 if (fNPads) preprocesor.AddComponent(fNPads);
976 if (fTimePosition) preprocesor.AddComponent(fTimePosition);
977 if (fOverThreshold10) preprocesor.AddComponent(fOverThreshold10);
978 if (fOverThreshold20) preprocesor.AddComponent(fOverThreshold20);
979 if (fOverThreshold30) preprocesor.AddComponent(fOverThreshold30);
980
258cd111 981 preprocesor.DumpToFile(fname);
982}
c322f08a 983
984
266f8637 985//_____________________________________________________________________
c322f08a 986void AliTPCdataQA::MakeArrays(){
987 //
266f8637 988 // The arrays for expanding the raw data are defined and
989 // som parameters are intialised
c322f08a 990 //
991 AliTPCROC * roc = AliTPCROC::Instance();
992 //
266f8637 993 // To make the array big enough for all sectors we take
994 // the dimensions from the outer row of an OROC (the last sector)
995 //
996 fRowsMax = roc->GetNRows(roc->GetNSector()-1);
997 fPadsMax = roc->GetNPads(roc->GetNSector()-1,fRowsMax-1);
998 fTimeBinsMax = fLastTimeBin - fFirstTimeBin +1;
999
1000 //
1001 // Since we have added 2 pads (TimeBins) before and after the real pads (TimeBins)
1002 // to make sure that we can always query the exanded table even when the
1003 // max is on the edge
1004 //
1005
c322f08a 1006
266f8637 1007 fAllBins = new Float_t*[fRowsMax];
1008 fAllSigBins = new Int_t*[fRowsMax];
1009 fAllNSigBins = new Int_t[fRowsMax];
1010
1011 for (Int_t iRow = 0; iRow < fRowsMax; iRow++) {
c322f08a 1012 //
266f8637 1013 Int_t maxBin = (fTimeBinsMax+4)*(fPadsMax+4);
c322f08a 1014 fAllBins[iRow] = new Float_t[maxBin];
266f8637 1015 memset(fAllBins[iRow],0,sizeof(Float_t)*maxBin); // set all values to 0
c322f08a 1016 fAllSigBins[iRow] = new Int_t[maxBin];
266f8637 1017 fAllNSigBins[iRow] = 0;
c322f08a 1018 }
1019}
1020
1021
266f8637 1022//_____________________________________________________________________
c322f08a 1023void AliTPCdataQA::CleanArrays(){
1024 //
1025 //
1026 //
266f8637 1027
1028 for (Int_t iRow = 0; iRow < fRowsMax; iRow++) {
42919b08 1029
1030 // To speed up the performance by a factor 2 on cosmic data (and
1031 // presumably pp data as well) where the ocupancy is low, the
1032 // memset is only called if there is more than 1000 signals for a
1033 // row (of the order 1% occupancy)
1034 if(fAllNSigBins[iRow]<1000) {
1035
1036 Float_t* allBins = fAllBins[iRow];
1037 Int_t* sigBins = fAllSigBins[iRow];
1038 const Int_t nSignals = fAllNSigBins[iRow];
1039 for(Int_t i = 0; i < nSignals; i++)
23c9ab21 1040 allBins[sigBins[i]]=0;
42919b08 1041 } else {
23c9ab21 1042
42919b08 1043 Int_t maxBin = (fTimeBinsMax+4)*(fPadsMax+4);
1044 memset(fAllBins[iRow],0,sizeof(Float_t)*maxBin);
1045 }
1046
c322f08a 1047 fAllNSigBins[iRow]=0;
1048 }
1049}
1050
266f8637 1051//_____________________________________________________________________
1052void AliTPCdataQA::GetPadAndTimeBin(Int_t bin, Int_t& iPad, Int_t& iTimeBin){
1053 //
1054 // Return pad and timebin for a given bin
c322f08a 1055 //
266f8637 1056
1057 // Int_t bin = iPad*(fTimeBinsMax+4)+iTimeBin;
1058 iTimeBin = bin%(fTimeBinsMax+4);
1059 iPad = (bin-iTimeBin)/(fTimeBinsMax+4);
1060
1061 iPad -= 2;
1062 iTimeBin -= 2;
1063 iTimeBin += fFirstTimeBin;
1064
1065 R__ASSERT(iPad>=0 && iPad<=fPadsMax);
1066 R__ASSERT(iTimeBin>=fFirstTimeBin && iTimeBin<=fLastTimeBin);
1067}
1068
1069//_____________________________________________________________________
1070void AliTPCdataQA::SetExpandDigit(const Int_t iRow, Int_t iPad,
1071 Int_t iTimeBin, const Float_t signal)
1072{
c322f08a 1073 //
266f8637 1074 //
c322f08a 1075 //
266f8637 1076 R__ASSERT(iRow>=0 && iRow<fRowsMax);
1077 R__ASSERT(iPad>=0 && iPad<=fPadsMax);
1078 R__ASSERT(iTimeBin>=fFirstTimeBin && iTimeBin<=fLastTimeBin);
1079
1080 iTimeBin -= fFirstTimeBin;
1081 iPad += 2;
1082 iTimeBin += 2;
c322f08a 1083
266f8637 1084 Int_t bin = iPad*(fTimeBinsMax+4)+iTimeBin;
1085 fAllBins[iRow][bin] = signal;
1086 fAllSigBins[iRow][fAllNSigBins[iRow]] = bin;
1087 fAllNSigBins[iRow]++;
1088}
1089
1090Float_t AliTPCdataQA::GetQ(const Float_t* adcArray, const Int_t time,
1091 const Int_t pad, const Int_t maxTimeBins,
1092 Int_t& timeMin, Int_t& timeMax,
6a50ff96 1093 Int_t& padMin, Int_t& padMax) const
266f8637 1094{
1095 //
1096 // This methods return the charge in the bin time+pad*maxTimeBins
1097 // If the charge is above 0 it also updates the padMin, padMax, timeMin
1098 // and timeMax if necessary
1099 //
1100 Float_t charge = adcArray[time + pad*maxTimeBins];
1101 if(charge > 0) {
1102 timeMin = TMath::Min(time, timeMin); timeMax = TMath::Max(time, timeMax);
1103 padMin = TMath::Min(pad, padMin); padMax = TMath::Max(pad, padMax);
1104 }
1105 return charge;
c322f08a 1106}
ce4b4255 1107
1108//______________________________________________________________________________
6a50ff96 1109void AliTPCdataQA::Streamer(TBuffer &xRuub)
ce4b4255 1110{
1111 // Automatic schema evolution was first used from revision 4
1112 // Code based on:
1113 // http://root.cern.ch/root/roottalk/roottalk02/3207.html
1114
6a50ff96 1115 UInt_t xRuus, xRuuc;
1116 if (xRuub.IsReading()) {
1117 Version_t xRuuv = xRuub.ReadVersion(&xRuus, &xRuuc);
ce4b4255 1118 //we use the automatic algorithm for class version > 3
6a50ff96 1119 if (xRuuv > 3) {
1120 AliTPCdataQA::Class()->ReadBuffer(xRuub, this, xRuuv, xRuus,
1121 xRuuc);
ce4b4255 1122 return;
1123 }
6a50ff96 1124 TH1F::Streamer(xRuub);
1125 xRuub >> fFirstTimeBin;
1126 xRuub >> fLastTimeBin;
1127 xRuub >> fAdcMin;
1128 xRuub >> fAdcMax;
1129 xRuub >> fNLocalMaxima;
1130 xRuub >> fMaxCharge;
1131 xRuub >> fMeanCharge;
1132 xRuub >> fNoThreshold;
1133 xRuub >> fNTimeBins;
1134 xRuub >> fNPads;
1135 xRuub >> fTimePosition;
1136 xRuub >> fEventCounter;
1137 xRuub >> fIsAnalysed;
1138 xRuub.CheckByteCount(xRuus, xRuuc, AliTPCdataQA::IsA());
ce4b4255 1139 } else {
6a50ff96 1140 AliTPCdataQA::Class()->WriteBuffer(xRuub,this);
ce4b4255 1141 }
1142}