1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
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
31 In Update the RAW signals >= 3 ADC channels are stored in the 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
38 This is done sector by sector.
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.
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.
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
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.
71 #include <TDirectory.h>
76 #include "AliRawReader.h"
77 #include "AliRawReaderRoot.h"
78 #include "AliRawReaderDate.h"
79 #include "AliTPCRawStream.h"
80 #include "AliTPCRawStreamV3.h"
81 #include "AliTPCCalROC.h"
82 #include "AliTPCROC.h"
83 #include "AliMathBase.h"
84 #include "TTreeStream.h"
85 #include "AliTPCRawStreamFast.h"
89 #include "AliTPCCalPad.h"
90 #include "AliTPCPreprocessorOnline.h"
93 #include "AliTPCdataQA.h"
96 ClassImp(AliTPCdataQA)
98 AliTPCdataQA::AliTPCdataQA() : /*FOLD00*/
116 fHistQVsTimeSideA(0),
117 fHistQVsTimeSideC(0),
118 fHistQMaxVsTimeSideA(0),
119 fHistQMaxVsTimeSideC(0),
120 fHistOccupancyVsEvent(0),
121 fHistNclustersVsEvent(0),
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
136 // default constructor
140 //_____________________________________________________________________
141 AliTPCdataQA::AliTPCdataQA(AliRecoParam::EventSpecie_t es) :
159 fHistQVsTimeSideA(0),
160 fHistQVsTimeSideC(0),
161 fHistQMaxVsTimeSideA(0),
162 fHistQMaxVsTimeSideC(0),
163 fHistOccupancyVsEvent(0),
164 fHistNclustersVsEvent(0),
178 // ctor creating the histogram
179 char * name = Form("TPCRAW_%s", AliRecoParam::GetEventSpecieName(es)) ;
180 TH1F(name, name,100,0,100) ;
183 //_____________________________________________________________________
184 AliTPCdataQA::AliTPCdataQA(const AliTPCdataQA &ped) : /*FOLD00*/
186 fFirstTimeBin(ped.GetFirstTimeBin()),
187 fLastTimeBin(ped.GetLastTimeBin()),
188 fAdcMin(ped.GetAdcMin()),
189 fAdcMax(ped.GetAdcMax()),
203 fHistQVsTimeSideA(0),
204 fHistQVsTimeSideC(0),
205 fHistQMaxVsTimeSideA(0),
206 fHistQMaxVsTimeSideC(0),
207 fHistOccupancyVsEvent(0),
208 fHistNclustersVsEvent(0),
209 fEventCounter(ped.GetEventCounter()),
210 fIsAnalysed(ped.GetIsAnalysed()),
211 fMaxEvents(ped.GetMaxEvents()),
212 fEventsPerBin(ped.GetEventsPerBin()),
213 fSignalCounter(ped.GetSignalCounter()),
214 fClusterCounter(ped.GetClusterCounter()),
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());
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());
245 if(ped.GetHistQVsTimeSideA()) {
246 fHistQVsTimeSideA = new TProfile(*ped.GetHistQVsTimeSideA());
247 fHistQVsTimeSideA->SetDirectory(0);
249 if(ped.GetHistQVsTimeSideC()) {
250 fHistQVsTimeSideC = new TProfile(*ped.GetHistQVsTimeSideC());
251 fHistQVsTimeSideC->SetDirectory(0);
253 if(ped.GetHistQMaxVsTimeSideA()) {
254 fHistQMaxVsTimeSideA = new TProfile(*ped.GetHistQMaxVsTimeSideA());
255 fHistQMaxVsTimeSideA->SetDirectory(0);
257 if(ped.GetHistQMaxVsTimeSideC()) {
258 fHistQMaxVsTimeSideC = new TProfile(*ped.GetHistQMaxVsTimeSideC());
259 fHistQMaxVsTimeSideC->SetDirectory(0);
261 if(ped.GetHistOccupancyVsEventConst()) {
262 fHistOccupancyVsEvent = new TH1F(*ped.GetHistOccupancyVsEventConst());
263 fHistOccupancyVsEvent->SetDirectory(0);
265 if(ped.GetHistNclustersVsEventConst()) {
266 fHistNclustersVsEvent = new TH1F(*ped.GetHistNclustersVsEventConst());
267 fHistNclustersVsEvent->SetDirectory(0);
271 //_____________________________________________________________________
272 AliTPCdataQA::AliTPCdataQA(const TMap *config) : /*FOLD00*/
273 TH1F("TPCRAW","TPCRAW",100,0,100),
291 fHistQVsTimeSideA(0),
292 fHistQVsTimeSideC(0),
293 fHistQMaxVsTimeSideA(0),
294 fHistQMaxVsTimeSideC(0),
295 fHistOccupancyVsEvent(0),
296 fHistNclustersVsEvent(0),
311 // default constructor
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();
317 if (config->GetValue("MaxEvents")) fMaxEvents = ((TObjString*)config->GetValue("MaxEvents"))->GetString().Atoi();
318 if (config->GetValue("EventsPerBin")) fAdcMax = ((TObjString*)config->GetValue("EventsPerBin"))->GetString().Atoi();
321 //_____________________________________________________________________
322 AliTPCdataQA& AliTPCdataQA::operator = (const AliTPCdataQA &source)
325 // assignment operator
327 if (&source == this) return *this;
328 new (this) AliTPCdataQA(source);
334 //_____________________________________________________________________
335 AliTPCdataQA::~AliTPCdataQA() /*FOLD00*/
341 // do not delete fMapping, because we do not own it.
342 // do not delete fMapping, because we do not own it.
343 // do not delete fNoise and fPedestal, because we do not own them.
345 delete fNLocalMaxima;
351 delete fTimePosition;
352 delete fOverThreshold10;
353 delete fOverThreshold20;
354 delete fOverThreshold30;
355 delete fHistQVsTimeSideA;
356 delete fHistQVsTimeSideC;
357 delete fHistQMaxVsTimeSideA;
358 delete fHistQMaxVsTimeSideC;
359 delete fHistOccupancyVsEvent;
360 delete fHistNclustersVsEvent;
362 for (Int_t iRow = 0; iRow < fRowsMax; iRow++) {
363 delete [] fAllBins[iRow];
364 delete [] fAllSigBins[iRow];
367 delete [] fAllSigBins;
368 delete [] fAllNSigBins;
371 //_____________________________________________________________________
372 TH1F* AliTPCdataQA::GetHistOccupancyVsEvent()
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)
380 if(!fHistOccupancyVsEvent) {
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);
388 return fHistOccupancyVsEvent;
391 //_____________________________________________________________________
392 TH1F* AliTPCdataQA::GetHistNclustersVsEvent()
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)
400 if(!fHistNclustersVsEvent) {
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);
408 return fHistNclustersVsEvent;
411 //_____________________________________________________________________
412 void AliTPCdataQA::UpdateEventHistograms()
414 // Update histograms that display occupancy and
415 // number of clusters as a function of number of
417 if (!fHistOccupancyVsEvent)
418 GetHistOccupancyVsEvent();
419 if (!fHistNclustersVsEvent)
420 GetHistNclustersVsEvent();
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);
429 Float_t averageNclusters =
430 Float_t(fClusterCounter)/fEventsPerBin;
431 if(fEventCounter<=fMaxEvents)
432 UpdateEventHisto(fHistNclustersVsEvent, averageNclusters);
436 //_____________________________________________________________________
437 void AliTPCdataQA::UpdateEventHisto(TH1F* hist, Float_t average)
439 // Do the actually updating of each histogram and
440 // change the visible range if needed
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);
447 if (hist->GetBinContent(bin)>0) {
448 AliError("Bin already filled. This should not happen.");
450 hist->SetBinContent(bin, average);
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));
459 //_____________________________________________________________________
460 Bool_t AliTPCdataQA::ProcessEvent(AliTPCRawStreamV3 *rawStreamV3)
463 // Event Processing loop - AliTPCRawStreamV3
465 Bool_t withInput = kFALSE;
467 Int_t lastSector = -1;
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) {
479 FindLocalMaxima(lastSector);
482 lastSector = iSector;
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();
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);
500 if (lastSector>=0&&nSignals>0)
501 FindLocalMaxima(lastSector);
506 //_____________________________________________________________________
507 Bool_t AliTPCdataQA::ProcessEvent(AliRawReader *rawReader)
510 // Event processing loop - AliRawReader
512 AliTPCRawStreamV3 rawStreamV3(rawReader, (AliAltroMapping**)fMapping);
513 Bool_t res=ProcessEvent(&rawStreamV3);
515 fEventCounter++; // only increment event counter if there is TPC data
517 if(fEventCounter%fEventsPerBin==0)
518 UpdateEventHistograms();
523 //_____________________________________________________________________
524 Bool_t AliTPCdataQA::ProcessEventFast(AliTPCRawStreamFast *rawStreamFast)
527 // Event Processing loop - AliTPCRawStream
529 Bool_t withInput = kFALSE;
531 Int_t lastSector = -1;
533 while ( rawStreamFast->NextDDL() ){
534 while ( rawStreamFast->NextChannel() ){
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) {
543 FindLocalMaxima(lastSector);
546 lastSector = iSector;
550 while ( rawStreamFast->NextBunch() ){
551 Int_t startTbin = (Int_t)rawStreamFast->GetStartTimeBin();
552 Int_t endTbin = (Int_t)rawStreamFast->GetEndTimeBin();
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);
565 //_____________________________________________________________________
566 Bool_t AliTPCdataQA::ProcessEventFast(AliRawReader *rawReader)
569 // Event processing loop - AliRawReader
571 AliTPCRawStreamFast rawStreamFast(rawReader, (AliAltroMapping**)fMapping);
572 Bool_t res=ProcessEventFast(&rawStreamFast);
574 fEventCounter++; // only increment event counter if there is TPC data
575 // otherwise Analyse (called in QA) fails
580 //_____________________________________________________________________
581 Bool_t AliTPCdataQA::ProcessEvent(AliTPCRawStream *rawStream)
584 // Event Processing loop - AliTPCRawStream
587 Bool_t withInput = kFALSE;
589 Int_t lastSector = -1;
591 while (rawStream->Next()) {
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
599 // Call local maxima finder if the data is in a new sector
600 if(iSector != lastSector) {
603 FindLocalMaxima(lastSector);
606 lastSector = iSector;
610 // Sometimes iRow==-1 if there is problems to read the data
612 nSignals += Update(iSector,iRow,iPad,iTimeBin,signal);
617 if (lastSector>=0&&nSignals>0)
618 FindLocalMaxima(lastSector);
624 //_____________________________________________________________________
625 Bool_t AliTPCdataQA::ProcessEventOld(AliRawReader *rawReader)
628 // Event processing loop - AliRawReader
631 // if fMapping is NULL the rawstream will crate its own mapping
632 AliTPCRawStream rawStream(rawReader, (AliAltroMapping**)fMapping);
633 rawReader->Select("TPC");
634 Bool_t res = ProcessEvent(&rawStream);
637 fEventCounter++; // only increment event counter if there is TPC data
638 // otherwise Analyse (called in QA) fails
643 //_____________________________________________________________________
644 Bool_t AliTPCdataQA::ProcessEvent(eventHeaderStruct *event)
647 // process date event
650 AliRawReaderDate rawReader((void*)event);
651 Bool_t result=ProcessEvent(&rawReader);
657 //_____________________________________________________________________
658 void AliTPCdataQA::DumpToFile(const Char_t *filename, const Char_t *dir, Bool_t append) /*FOLD00*/
661 // Write class to file
672 TDirectory *backup = gDirectory;
673 TFile f(filename,option.Data());
675 if ( !sDir.IsNull() ){
676 f.mkdir(sDir.Data());
682 if ( backup ) backup->cd();
686 //_____________________________________________________________________
687 Int_t AliTPCdataQA::Update(const Int_t iSector, /*FOLD00*/
690 const Int_t iTimeBin,
694 // Signal filling method
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
702 if (!fNLocalMaxima) fNLocalMaxima = new AliTPCCalPad("NLocalMaxima","NLocalMaxima");
703 if (!fMaxCharge) fMaxCharge = new AliTPCCalPad("MaxCharge","MaxCharge");
704 if (!fMeanCharge) fMeanCharge = new AliTPCCalPad("MeanCharge","MeanCharge");
705 if (!fNoThreshold) fNoThreshold = new AliTPCCalPad("NoThreshold","NoThreshold");
706 if (!fNTimeBins) fNTimeBins = new AliTPCCalPad("NTimeBins","NTimeBins");
707 if (!fNPads) fNPads = new AliTPCCalPad("NPads","NPads");
708 if (!fTimePosition) fTimePosition = new AliTPCCalPad("TimePosition","TimePosition");
709 if (!fOverThreshold10) fOverThreshold10 = new AliTPCCalPad("OverThreshold10","OverThreshold10");
710 if (!fOverThreshold20) fOverThreshold20 = new AliTPCCalPad("OverThreshold20","OverThreshold20");
711 if (!fOverThreshold30) fOverThreshold30 = new AliTPCCalPad("OverThreshold30","OverThreshold30");
712 if (!fHistQVsTimeSideA) {
713 fHistQVsTimeSideA = new TProfile("hQVsTimeSideA", "Q vs time (side A); Time [Timebin]; Q [ADC ch]", 100, 0, 1000);
714 fHistQVsTimeSideA->SetDirectory(0);
716 if (!fHistQVsTimeSideC) {
717 fHistQVsTimeSideC = new TProfile("hQVsTimeSideC", "Q vs time (side C); Time [Timebin]; Q [ADC ch]", 100, 0, 1000);
718 fHistQVsTimeSideC->SetDirectory(0);
720 if (!fHistQMaxVsTimeSideA) {
721 fHistQMaxVsTimeSideA = new TProfile("hQMaxVsTimeSideA", "Q_{MAX} vs time (side A); Time [Timebin]; Q_{MAX} [ADC ch]", 100, 0, 1000);
722 fHistQMaxVsTimeSideA->SetDirectory(0);
724 if (!fHistQMaxVsTimeSideC) {
725 fHistQMaxVsTimeSideC = new TProfile("hQMaxVsTimeSideC", "Q_{MAX} vs time (side C); Time [Timebin]; Q_{MAX} [ADC ch]", 100, 0, 1000);
726 fHistQMaxVsTimeSideC->SetDirectory(0);
729 // Make the arrays for expanding the data
735 // If Analyse has been previously called we need now to denormalize the data
736 // as more data is coming
738 if(fIsAnalysed == kTRUE) {
740 const Int_t nTimeBins = fLastTimeBin - fFirstTimeBin +1;
741 const Float_t denormalization = Float_t(fEventCounter * nTimeBins);
742 fNoThreshold->Multiply(denormalization);
744 fMeanCharge->Multiply(fNLocalMaxima);
745 fMaxCharge->Multiply(fNLocalMaxima);
746 fNTimeBins->Multiply(fNLocalMaxima);
747 fNPads->Multiply(fNLocalMaxima);
748 fTimePosition->Multiply(fNLocalMaxima);
749 fIsAnalysed = kFALSE;
755 if (iTimeBin<fFirstTimeBin) return 0;
756 if (iTimeBin>fLastTimeBin) return 0;
758 // if pedestal calibrations are loaded subtract pedestals
761 Float_t ped = fPedestal->GetCalROC(iSector)->GetValue(iRow, iPad);
762 // Don't use data from pads where pedestals are abnormally small or large
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);
772 // Require at least 3 ADC channels
776 // if noise calibrations are loaded require at least 3*sigmaNoise
779 Float_t noise = fNoise->GetCalROC(iSector)->GetValue(iRow, iPad);
781 if(signal < noise*3.0)
786 // This signal is ok and we store it in the cluster map
789 SetExpandDigit(iRow, iPad, iTimeBin, signal);
793 return 1; // signal was accepted
796 //_____________________________________________________________________
797 void AliTPCdataQA::FindLocalMaxima(const Int_t iSector)
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
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
811 for (Int_t iRow = 0; iRow < fRowsMax; iRow++) {
813 Float_t* allBins = fAllBins[iRow];
814 Int_t* sigBins = fAllSigBins[iRow];
815 const Int_t nSigBins = fAllNSigBins[iRow];
817 for (Int_t iSig = 0; iSig < nSigBins; iSig++) {
819 Int_t bin = sigBins[iSig];
820 Float_t *b = &allBins[bin];
823 // Now we check if this is a local maximum
828 // First check that the charge is bigger than the threshold
832 // Require at least one neighboring pad with signal
833 if (b[-maxTimeBin]+b[maxTimeBin]<=0) continue;
835 // Require at least one neighboring time bin with signal
836 if (b[-1]+b[1]<=0) continue;
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 >=)
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;
854 // Now we accept the local maximum and fill the calibration/data objects
858 Int_t iPad, iTimeBin;
859 GetPadAndTimeBin(bin, iPad, iTimeBin);
861 Float_t count = fNLocalMaxima->GetCalROC(iSector)->GetValue(iRow, iPad);
862 fNLocalMaxima->GetCalROC(iSector)->SetValue(iRow, iPad, count+1);
864 count = fTimePosition->GetCalROC(iSector)->GetValue(iRow, iPad);
865 fTimePosition->GetCalROC(iSector)->SetValue(iRow, iPad, count+iTimeBin);
867 Float_t charge = fMaxCharge->GetCalROC(iSector)->GetValue(iRow, iPad);
868 fMaxCharge->GetCalROC(iSector)->SetValue(iRow, iPad, charge + qMax);
871 count = fOverThreshold10->GetCalROC(iSector)->GetValue(iRow, iPad);
872 fOverThreshold10->GetCalROC(iSector)->SetValue(iRow, iPad, count+1);
875 count = fOverThreshold20->GetCalROC(iSector)->GetValue(iRow, iPad);
876 fOverThreshold20->GetCalROC(iSector)->SetValue(iRow, iPad, count+1);
879 count = fOverThreshold30->GetCalROC(iSector)->GetValue(iRow, iPad);
880 fOverThreshold30->GetCalROC(iSector)->SetValue(iRow, iPad, count+1);
884 // Calculate the total charge as the sum over the region:
892 // with qmax at the center C.
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.
897 Int_t minP = 0, maxP = 0, minT = 0, maxT = 0;
899 for(Int_t i = -1; i<=1; i++) {
900 for(Int_t j = -1; j<=1; j++) {
905 Float_t charge1 = GetQ(b, i, j, maxTimeBin, minT, maxT, minP, maxP);
908 // see if the next neighbor is also above threshold
910 qTot += GetQ(b, 2*i, 2*j, maxTimeBin, minT, maxT, minP, maxP);
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);
921 charge = fMeanCharge->GetCalROC(iSector)->GetValue(iRow, iPad);
922 fMeanCharge->GetCalROC(iSector)->SetValue(iRow, iPad, charge + qTot);
924 count = fNTimeBins->GetCalROC(iSector)->GetValue(iRow, iPad);
925 fNTimeBins->GetCalROC(iSector)->SetValue(iRow, iPad, count + maxT-minT+1);
927 count = fNPads->GetCalROC(iSector)->GetValue(iRow, iPad);
928 fNPads->GetCalROC(iSector)->SetValue(iRow, iPad, count + maxP-minP+1);
930 if((iSector%36)<18) { // A side
931 fHistQVsTimeSideA->Fill(iTimeBin, qTot);
932 fHistQMaxVsTimeSideA->Fill(iTimeBin, qMax);
934 fHistQVsTimeSideC->Fill(iTimeBin, qTot);
935 fHistQMaxVsTimeSideC->Fill(iTimeBin, qMax);
937 } // end loop over signals
938 } // end loop over rows
940 // cout << "Number of local maximas found: " << nLocalMaxima << endl;
941 fClusterCounter += nLocalMaxima;
944 //_____________________________________________________________________
945 void AliTPCdataQA::Analyse()
948 // Calculate calibration constants
951 cout << "Analyse called" << endl;
953 if(fIsAnalysed == kTRUE) {
955 cout << "No new data since Analyse was called last time" << endl;
959 if(fEventCounter==0) {
961 cout << "EventCounter == 0, Cannot analyse" << endl;
965 Int_t nTimeBins = fLastTimeBin - fFirstTimeBin +1;
966 cout << "EventCounter: " << fEventCounter << endl
967 << "TimeBins: " << nTimeBins << endl;
969 Float_t normalization = 1.0 / Float_t(fEventCounter * nTimeBins);
970 fNoThreshold->Multiply(normalization);
972 fMeanCharge->Divide(fNLocalMaxima);
973 fMaxCharge->Divide(fNLocalMaxima);
974 fNTimeBins->Divide(fNLocalMaxima);
975 fNPads->Divide(fNLocalMaxima);
976 fTimePosition->Divide(fNLocalMaxima);
982 //_____________________________________________________________________
983 void AliTPCdataQA::MakeTree(const char *fname){
985 // Export result to the tree -located in the file
986 // This file can be analyzed using AliTPCCalibViewer
988 AliTPCPreprocessorOnline preprocesor;
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);
1001 preprocesor.DumpToFile(fname);
1005 //_____________________________________________________________________
1006 void AliTPCdataQA::MakeArrays(){
1008 // The arrays for expanding the raw data are defined and
1009 // som parameters are intialised
1011 AliTPCROC * roc = AliTPCROC::Instance();
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)
1016 fRowsMax = roc->GetNRows(roc->GetNSector()-1);
1017 fPadsMax = roc->GetNPads(roc->GetNSector()-1,fRowsMax-1);
1018 fTimeBinsMax = fLastTimeBin - fFirstTimeBin +1;
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
1027 fAllBins = new Float_t*[fRowsMax];
1028 fAllSigBins = new Int_t*[fRowsMax];
1029 fAllNSigBins = new Int_t[fRowsMax];
1031 for (Int_t iRow = 0; iRow < fRowsMax; iRow++) {
1033 Int_t maxBin = (fTimeBinsMax+4)*(fPadsMax+4);
1034 fAllBins[iRow] = new Float_t[maxBin];
1035 memset(fAllBins[iRow],0,sizeof(Float_t)*maxBin); // set all values to 0
1036 fAllSigBins[iRow] = new Int_t[maxBin];
1037 fAllNSigBins[iRow] = 0;
1042 //_____________________________________________________________________
1043 void AliTPCdataQA::CleanArrays(){
1048 for (Int_t iRow = 0; iRow < fRowsMax; iRow++) {
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) {
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++)
1060 allBins[sigBins[i]]=0;
1063 Int_t maxBin = (fTimeBinsMax+4)*(fPadsMax+4);
1064 memset(fAllBins[iRow],0,sizeof(Float_t)*maxBin);
1067 fAllNSigBins[iRow]=0;
1071 //_____________________________________________________________________
1072 void AliTPCdataQA::GetPadAndTimeBin(Int_t bin, Int_t& iPad, Int_t& iTimeBin){
1074 // Return pad and timebin for a given bin
1077 // Int_t bin = iPad*(fTimeBinsMax+4)+iTimeBin;
1078 iTimeBin = bin%(fTimeBinsMax+4);
1079 iPad = (bin-iTimeBin)/(fTimeBinsMax+4);
1083 iTimeBin += fFirstTimeBin;
1085 R__ASSERT(iPad>=0 && iPad<=fPadsMax);
1086 R__ASSERT(iTimeBin>=fFirstTimeBin && iTimeBin<=fLastTimeBin);
1089 //_____________________________________________________________________
1090 void AliTPCdataQA::SetExpandDigit(const Int_t iRow, Int_t iPad,
1091 Int_t iTimeBin, const Float_t signal)
1096 R__ASSERT(iRow>=0 && iRow<fRowsMax);
1097 R__ASSERT(iPad>=0 && iPad<=fPadsMax);
1098 R__ASSERT(iTimeBin>=fFirstTimeBin && iTimeBin<=fLastTimeBin);
1100 iTimeBin -= fFirstTimeBin;
1104 Int_t bin = iPad*(fTimeBinsMax+4)+iTimeBin;
1105 fAllBins[iRow][bin] = signal;
1106 fAllSigBins[iRow][fAllNSigBins[iRow]] = bin;
1107 fAllNSigBins[iRow]++;
1110 Float_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)
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
1120 Float_t charge = adcArray[time + pad*maxTimeBins];
1122 timeMin = TMath::Min(time, timeMin); timeMax = TMath::Max(time, timeMax);
1123 padMin = TMath::Min(pad, padMin); padMax = TMath::Max(pad, padMax);
1128 //______________________________________________________________________________
1129 void AliTPCdataQA::Streamer(TBuffer &R__b)
1131 // Automatic schema evolution was first used from revision 4
1133 // http://root.cern.ch/root/roottalk/roottalk02/3207.html
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
1140 AliTPCdataQA::Class()->ReadBuffer(R__b, this, R__v, R__s,
1144 TH1F::Streamer(R__b);
1145 R__b >> fFirstTimeBin;
1146 R__b >> fLastTimeBin;
1149 R__b >> fNLocalMaxima;
1151 R__b >> fMeanCharge;
1152 R__b >> fNoThreshold;
1155 R__b >> fTimePosition;
1156 R__b >> fEventCounter;
1157 R__b >> fIsAnalysed;
1158 R__b.CheckByteCount(R__s, R__c, AliTPCdataQA::IsA());
1160 AliTPCdataQA::Class()->WriteBuffer(R__b,this);