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