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 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.
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.
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.
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.
46 Finally cout was changes to AliInfo.
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
59 In Update the RAW signals >= 3 ADC channels are stored in the 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
66 This is done sector by sector.
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.
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.
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
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.
92 #include <TDirectory.h>
98 #include "AliRawReader.h"
99 #include "AliRawReaderRoot.h"
100 #include "AliRawReaderDate.h"
101 #include "AliTPCRawStream.h"
102 #include "AliTPCRawStreamV3.h"
103 #include "AliTPCCalROC.h"
104 #include "AliTPCROC.h"
105 #include "AliMathBase.h"
106 #include "TTreeStream.h"
107 #include "AliTPCRawStreamFast.h"
111 #include "AliTPCCalPad.h"
112 #include "AliTPCPreprocessorOnline.h"
115 #include "AliTPCdataQA.h"
118 ClassImp(AliTPCdataQA)
120 AliTPCdataQA::AliTPCdataQA() : /*FOLD00*/
138 fHistQVsTimeSideA(0),
139 fHistQVsTimeSideC(0),
140 fHistQMaxVsTimeSideA(0),
141 fHistQMaxVsTimeSideC(0),
142 fHistOccupancyVsEvent(0),
143 fHistNclustersVsEvent(0),
146 fMaxEvents(500000), // Max events for event histograms
147 fEventsPerBin(1000), // Events per bin for event histograms
148 fSignalCounter(0), // Signal counter
149 fClusterCounter(0), // Cluster counter
158 // default constructor
162 //_____________________________________________________________________
163 AliTPCdataQA::AliTPCdataQA(AliRecoParam::EventSpecie_t es) :
181 fHistQVsTimeSideA(0),
182 fHistQVsTimeSideC(0),
183 fHistQMaxVsTimeSideA(0),
184 fHistQMaxVsTimeSideC(0),
185 fHistOccupancyVsEvent(0),
186 fHistNclustersVsEvent(0),
200 // ctor creating the histogram
201 char * name = Form("TPCRAW_%s", AliRecoParam::GetEventSpecieName(es)) ;
202 TH1F(name, name,100,0,100) ;
205 //_____________________________________________________________________
206 AliTPCdataQA::AliTPCdataQA(const AliTPCdataQA &ped) : /*FOLD00*/
208 fFirstTimeBin(ped.GetFirstTimeBin()),
209 fLastTimeBin(ped.GetLastTimeBin()),
210 fAdcMin(ped.GetAdcMin()),
211 fAdcMax(ped.GetAdcMax()),
225 fHistQVsTimeSideA(0),
226 fHistQVsTimeSideC(0),
227 fHistQMaxVsTimeSideA(0),
228 fHistQMaxVsTimeSideC(0),
229 fHistOccupancyVsEvent(0),
230 fHistNclustersVsEvent(0),
231 fEventCounter(ped.GetEventCounter()),
232 fIsAnalysed(ped.GetIsAnalysed()),
233 fMaxEvents(ped.GetMaxEvents()),
234 fEventsPerBin(ped.GetEventsPerBin()),
235 fSignalCounter(ped.GetSignalCounter()),
236 fClusterCounter(ped.GetClusterCounter()),
247 if(ped.GetNLocalMaxima())
248 fNLocalMaxima = new AliTPCCalPad(*ped.GetNLocalMaxima());
249 if(ped.GetMaxCharge())
250 fMaxCharge = new AliTPCCalPad(*ped.GetMaxCharge());
251 if(ped.GetMeanCharge())
252 fMeanCharge = new AliTPCCalPad(*ped.GetMeanCharge());
253 if(ped.GetNoThreshold())
254 fNoThreshold = new AliTPCCalPad(*ped.GetNoThreshold());
255 if(ped.GetNTimeBins())
256 fNTimeBins = new AliTPCCalPad(*ped.GetNTimeBins());
258 fNPads = new AliTPCCalPad(*ped.GetNPads());
259 if(ped.GetTimePosition())
260 fTimePosition = new AliTPCCalPad(*ped.GetTimePosition());
261 if(ped.GetOverThreshold10())
262 fOverThreshold10 = new AliTPCCalPad(*ped.GetOverThreshold10());
263 if(ped.GetOverThreshold20())
264 fOverThreshold20 = new AliTPCCalPad(*ped.GetOverThreshold20());
265 if(ped.GetOverThreshold30())
266 fOverThreshold30 = new AliTPCCalPad(*ped.GetOverThreshold30());
267 if(ped.GetHistQVsTimeSideA()) {
268 fHistQVsTimeSideA = new TProfile(*ped.GetHistQVsTimeSideA());
269 fHistQVsTimeSideA->SetDirectory(0);
271 if(ped.GetHistQVsTimeSideC()) {
272 fHistQVsTimeSideC = new TProfile(*ped.GetHistQVsTimeSideC());
273 fHistQVsTimeSideC->SetDirectory(0);
275 if(ped.GetHistQMaxVsTimeSideA()) {
276 fHistQMaxVsTimeSideA = new TProfile(*ped.GetHistQMaxVsTimeSideA());
277 fHistQMaxVsTimeSideA->SetDirectory(0);
279 if(ped.GetHistQMaxVsTimeSideC()) {
280 fHistQMaxVsTimeSideC = new TProfile(*ped.GetHistQMaxVsTimeSideC());
281 fHistQMaxVsTimeSideC->SetDirectory(0);
283 if(ped.GetHistOccupancyVsEventConst()) {
284 fHistOccupancyVsEvent = new TH1F(*ped.GetHistOccupancyVsEventConst());
285 fHistOccupancyVsEvent->SetDirectory(0);
287 if(ped.GetHistNclustersVsEventConst()) {
288 fHistNclustersVsEvent = new TH1F(*ped.GetHistNclustersVsEventConst());
289 fHistNclustersVsEvent->SetDirectory(0);
293 //_____________________________________________________________________
294 AliTPCdataQA::AliTPCdataQA(const TMap *config) : /*FOLD00*/
295 TH1F("TPCRAW","TPCRAW",100,0,100),
313 fHistQVsTimeSideA(0),
314 fHistQVsTimeSideC(0),
315 fHistQMaxVsTimeSideA(0),
316 fHistQMaxVsTimeSideC(0),
317 fHistOccupancyVsEvent(0),
318 fHistNclustersVsEvent(0),
333 // default constructor
335 if (config->GetValue("FirstTimeBin")) fFirstTimeBin = ((TObjString*)config->GetValue("FirstTimeBin"))->GetString().Atoi();
336 if (config->GetValue("LastTimeBin")) fLastTimeBin = ((TObjString*)config->GetValue("LastTimeBin"))->GetString().Atoi();
337 if (config->GetValue("AdcMin")) fAdcMin = ((TObjString*)config->GetValue("AdcMin"))->GetString().Atoi();
338 if (config->GetValue("AdcMax")) fAdcMax = ((TObjString*)config->GetValue("AdcMax"))->GetString().Atoi();
339 if (config->GetValue("MaxEvents")) fMaxEvents = ((TObjString*)config->GetValue("MaxEvents"))->GetString().Atoi();
340 if (config->GetValue("EventsPerBin")) fEventsPerBin = ((TObjString*)config->GetValue("EventsPerBin"))->GetString().Atoi();
343 //_____________________________________________________________________
344 AliTPCdataQA& AliTPCdataQA::operator = (const AliTPCdataQA &source)
347 // assignment operator
349 if (&source == this) return *this;
350 new (this) AliTPCdataQA(source);
356 //_____________________________________________________________________
357 AliTPCdataQA::~AliTPCdataQA() /*FOLD00*/
363 // do not delete fMapping, because we do not own it.
364 // do not delete fMapping, because we do not own it.
365 // do not delete fNoise and fPedestal, because we do not own them.
367 delete fNLocalMaxima;
373 delete fTimePosition;
374 delete fOverThreshold10;
375 delete fOverThreshold20;
376 delete fOverThreshold30;
377 delete fHistQVsTimeSideA;
378 delete fHistQVsTimeSideC;
379 delete fHistQMaxVsTimeSideA;
380 delete fHistQMaxVsTimeSideC;
381 delete fHistOccupancyVsEvent;
382 delete fHistNclustersVsEvent;
384 for (Int_t iRow = 0; iRow < fRowsMax; iRow++) {
385 delete [] fAllBins[iRow];
386 delete [] fAllSigBins[iRow];
389 delete [] fAllSigBins;
390 delete [] fAllNSigBins;
393 //_____________________________________________________________________
394 TH1F* AliTPCdataQA::GetHistOccupancyVsEvent()
397 // Create Occupancy vs event histogram
398 // (we create this histogram differently then the other histograms
399 // because this we want to be able to access and copy
400 // from AliTPCQAMakerRec before it normally would be created)
402 if(!fHistOccupancyVsEvent) {
404 Int_t nBins = fMaxEvents/fEventsPerBin;
405 fHistOccupancyVsEvent = new TH1F("hOccupancyVsEvent", "Occupancy vs event number (~time); Event number; Occupancy", nBins, 0, nBins*fEventsPerBin);
406 fHistOccupancyVsEvent->SetDirectory(0);
409 return fHistOccupancyVsEvent;
412 //_____________________________________________________________________
413 TH1F* AliTPCdataQA::GetHistNclustersVsEvent()
416 // Create Nclusters vs event histogram
417 // (we create this histogram differently then the other histograms
418 // because this we want to be able to access and copy
419 // from AliTPCQAMakerRec before it normally would be created)
421 if(!fHistNclustersVsEvent) {
423 Int_t nBins = fMaxEvents/fEventsPerBin;
424 fHistNclustersVsEvent = new TH1F("hNclustersVsEvent", "Nclusters vs event number (~time); Event number; Nclusters per event", nBins, 0, nBins*fEventsPerBin);
425 fHistNclustersVsEvent->SetDirectory(0);
428 return fHistNclustersVsEvent;
431 //_____________________________________________________________________
432 void AliTPCdataQA::UpdateEventHistograms()
434 // Update histograms that display occupancy and
435 // number of clusters as a function of number of
437 if (!fHistOccupancyVsEvent)
438 GetHistOccupancyVsEvent();
439 if (!fHistNclustersVsEvent)
440 GetHistNclustersVsEvent();
442 if(fEventCounter > fMaxEvents) {
444 // we have to expand the histogram to handle the larger number of
445 // events. The way it is done now is to double the range and the
446 // number of events per bin (so the number of histogram bins stays
451 // Change histogram limits
452 const Int_t nBins = fHistOccupancyVsEvent->GetXaxis()->GetNbins();
453 fHistOccupancyVsEvent->GetXaxis()->Set(nBins, fHistOccupancyVsEvent->GetXaxis()->GetNbins(), fMaxEvents);
454 fHistNclustersVsEvent->GetXaxis()->Set(nBins, fHistNclustersVsEvent->GetXaxis()->GetNbins(), fMaxEvents);
456 // Rebin the histogram
457 for(Int_t bin = 1; bin <= nBins; bin+=2) {
459 Int_t newBin = TMath::Nint(Float_t(bin+1)/2.0);
460 Float_t newContent = (fHistOccupancyVsEvent->GetBinContent(bin)+
461 fHistOccupancyVsEvent->GetBinContent(bin+1))/2.0;
462 fHistOccupancyVsEvent->SetBinContent(newBin, newContent);
464 newContent = (fHistNclustersVsEvent->GetBinContent(bin)+
465 fHistNclustersVsEvent->GetBinContent(bin+1))/2.0;
466 fHistNclustersVsEvent->SetBinContent(newBin, newContent);
469 // Set the remaining bins to 0
470 Int_t lastHalf = nBins/2 +1;
471 for(Int_t bin = lastHalf; bin <= nBins; bin++) {
473 fHistOccupancyVsEvent->SetBinContent(bin, 0);
474 fHistNclustersVsEvent->SetBinContent(bin, 0);
477 // In this case we should nut update but wait untill the new
478 // number of events per bin is reached!
482 const Int_t bin = TMath::Nint(Float_t(fEventCounter)/fEventsPerBin);
484 Float_t averageOccupancy =
485 Float_t(fSignalCounter)/fEventsPerBin/(fLastTimeBin - fFirstTimeBin +1.0)
486 / 570132.0; // 570,132 is number of pads
487 fHistOccupancyVsEvent->SetBinContent(bin, averageOccupancy);
490 Float_t averageNclusters =
491 Float_t(fClusterCounter)/fEventsPerBin;
492 fHistNclustersVsEvent->SetBinContent(bin, averageNclusters);
496 //_____________________________________________________________________
497 Bool_t AliTPCdataQA::ProcessEvent(AliTPCRawStreamV3 *const rawStreamV3)
500 // Event Processing loop - AliTPCRawStreamV3
502 Bool_t withInput = kFALSE;
504 Int_t lastSector = -1;
506 while ( rawStreamV3->NextDDL() ){
508 while ( rawStreamV3->NextChannel() ){
510 Int_t iSector = rawStreamV3->GetSector(); // current sector
511 Int_t iRow = rawStreamV3->GetRow(); // current row
512 Int_t iPad = rawStreamV3->GetPad(); // current pad
513 if (iRow<0 || iPad<0) continue;
514 // Call local maxima finder if the data is in a new sector
515 if(iSector != lastSector) {
518 FindLocalMaxima(lastSector);
521 lastSector = iSector;
525 while ( rawStreamV3->NextBunch() ){
527 Int_t startTbin = (Int_t)rawStreamV3->GetStartTimeBin();
528 Int_t bunchlength = (Int_t)rawStreamV3->GetBunchLength();
529 const UShort_t *sig = rawStreamV3->GetSignals();
531 for (Int_t iTimeBin = 0; iTimeBin<bunchlength; iTimeBin++){
532 Float_t signal=(Float_t)sig[iTimeBin];
533 nSignals += Update(iSector,iRow,iPad,startTbin--,signal);
540 if (lastSector>=0&&nSignals>0)
541 FindLocalMaxima(lastSector);
546 //_____________________________________________________________________
547 Bool_t AliTPCdataQA::ProcessEvent(AliRawReader *const rawReader)
550 // Event processing loop - AliRawReader
552 AliTPCRawStreamV3 rawStreamV3(rawReader, (AliAltroMapping**)fMapping);
553 Bool_t res=ProcessEvent(&rawStreamV3);
555 fEventCounter++; // only increment event counter if there is TPC data
557 if(fEventCounter%fEventsPerBin==0)
558 UpdateEventHistograms();
563 //_____________________________________________________________________
564 Bool_t AliTPCdataQA::ProcessEventFast(AliTPCRawStreamFast *const rawStreamFast)
567 // Event Processing loop - AliTPCRawStream
569 Bool_t withInput = kFALSE;
571 Int_t lastSector = -1;
573 while ( rawStreamFast->NextDDL() ){
574 while ( rawStreamFast->NextChannel() ){
576 Int_t iSector = rawStreamFast->GetSector(); // current sector
577 Int_t iRow = rawStreamFast->GetRow(); // current row
578 Int_t iPad = rawStreamFast->GetPad(); // current pad
579 // Call local maxima finder if the data is in a new sector
580 if(iSector != lastSector) {
583 FindLocalMaxima(lastSector);
586 lastSector = iSector;
590 while ( rawStreamFast->NextBunch() ){
591 Int_t startTbin = (Int_t)rawStreamFast->GetStartTimeBin();
592 Int_t endTbin = (Int_t)rawStreamFast->GetEndTimeBin();
594 for (Int_t iTimeBin = startTbin; iTimeBin < endTbin; iTimeBin++){
595 Float_t signal = rawStreamFast->GetSignals()[iTimeBin-startTbin];
596 nSignals += Update(iSector,iRow,iPad,iTimeBin+1,signal);
605 //_____________________________________________________________________
606 Bool_t AliTPCdataQA::ProcessEventFast(AliRawReader *const rawReader)
609 // Event processing loop - AliRawReader
611 AliTPCRawStreamFast rawStreamFast(rawReader, (AliAltroMapping**)fMapping);
612 Bool_t res=ProcessEventFast(&rawStreamFast);
614 fEventCounter++; // only increment event counter if there is TPC data
615 // otherwise Analyse (called in QA) fails
620 //_____________________________________________________________________
621 Bool_t AliTPCdataQA::ProcessEvent(AliTPCRawStream *const rawStream)
624 // Event Processing loop - AliTPCRawStream
627 Bool_t withInput = kFALSE;
629 Int_t lastSector = -1;
631 while (rawStream->Next()) {
633 Int_t iSector = rawStream->GetSector(); // current ROC
634 Int_t iRow = rawStream->GetRow(); // current row
635 Int_t iPad = rawStream->GetPad(); // current pad
636 Int_t iTimeBin = rawStream->GetTime(); // current time bin
637 Float_t signal = rawStream->GetSignal(); // current ADC signal
639 // Call local maxima finder if the data is in a new sector
640 if(iSector != lastSector) {
643 FindLocalMaxima(lastSector);
646 lastSector = iSector;
650 // Sometimes iRow==-1 if there is problems to read the data
652 nSignals += Update(iSector,iRow,iPad,iTimeBin,signal);
657 if (lastSector>=0&&nSignals>0)
658 FindLocalMaxima(lastSector);
664 //_____________________________________________________________________
665 Bool_t AliTPCdataQA::ProcessEventOld(AliRawReader *const rawReader)
668 // Event processing loop - AliRawReader
671 // if fMapping is NULL the rawstream will crate its own mapping
672 AliTPCRawStream rawStream(rawReader, (AliAltroMapping**)fMapping);
673 rawReader->Select("TPC");
674 Bool_t res = ProcessEvent(&rawStream);
677 fEventCounter++; // only increment event counter if there is TPC data
678 // otherwise Analyse (called in QA) fails
683 //_____________________________________________________________________
684 Bool_t AliTPCdataQA::ProcessEvent(eventHeaderStruct *const event)
687 // process date event
690 AliRawReaderDate rawReader((void*)event);
691 Bool_t result=ProcessEvent(&rawReader);
697 //_____________________________________________________________________
698 void AliTPCdataQA::DumpToFile(const Char_t *filename, const Char_t *dir, Bool_t append) /*FOLD00*/
701 // Write class to file
712 TDirectory *backup = gDirectory;
713 TFile f(filename,option.Data());
715 if ( !sDir.IsNull() ){
716 f.mkdir(sDir.Data());
722 if ( backup ) backup->cd();
726 //_____________________________________________________________________
727 Int_t AliTPCdataQA::Update(const Int_t iSector, /*FOLD00*/
730 const Int_t iTimeBin,
734 // Signal filling method
738 // Define the calibration objects the first time Update is called
739 // NB! This has to be done first even if the data is rejected by the time
740 // cut to make sure that the objects are available in Analyse
742 if (!fNLocalMaxima) fNLocalMaxima = new AliTPCCalPad("NLocalMaxima","NLocalMaxima");
743 if (!fMaxCharge) fMaxCharge = new AliTPCCalPad("MaxCharge","MaxCharge");
744 if (!fMeanCharge) fMeanCharge = new AliTPCCalPad("MeanCharge","MeanCharge");
745 if (!fNoThreshold) fNoThreshold = new AliTPCCalPad("NoThreshold","NoThreshold");
746 if (!fNTimeBins) fNTimeBins = new AliTPCCalPad("NTimeBins","NTimeBins");
747 if (!fNPads) fNPads = new AliTPCCalPad("NPads","NPads");
748 if (!fTimePosition) fTimePosition = new AliTPCCalPad("TimePosition","TimePosition");
749 if (!fOverThreshold10) fOverThreshold10 = new AliTPCCalPad("OverThreshold10","OverThreshold10");
750 if (!fOverThreshold20) fOverThreshold20 = new AliTPCCalPad("OverThreshold20","OverThreshold20");
751 if (!fOverThreshold30) fOverThreshold30 = new AliTPCCalPad("OverThreshold30","OverThreshold30");
752 if (!fHistQVsTimeSideA) {
753 fHistQVsTimeSideA = new TProfile("hQVsTimeSideA", "Q vs time (side A); Time [Timebin]; Q [ADC ch]", 100, 0, 1000);
754 fHistQVsTimeSideA->SetDirectory(0);
756 if (!fHistQVsTimeSideC) {
757 fHistQVsTimeSideC = new TProfile("hQVsTimeSideC", "Q vs time (side C); Time [Timebin]; Q [ADC ch]", 100, 0, 1000);
758 fHistQVsTimeSideC->SetDirectory(0);
760 if (!fHistQMaxVsTimeSideA) {
761 fHistQMaxVsTimeSideA = new TProfile("hQMaxVsTimeSideA", "Q_{MAX} vs time (side A); Time [Timebin]; Q_{MAX} [ADC ch]", 100, 0, 1000);
762 fHistQMaxVsTimeSideA->SetDirectory(0);
764 if (!fHistQMaxVsTimeSideC) {
765 fHistQMaxVsTimeSideC = new TProfile("hQMaxVsTimeSideC", "Q_{MAX} vs time (side C); Time [Timebin]; Q_{MAX} [ADC ch]", 100, 0, 1000);
766 fHistQMaxVsTimeSideC->SetDirectory(0);
769 // Make the arrays for expanding the data
775 // If Analyse has been previously called we need now to denormalize the data
776 // as more data is coming
778 if(fIsAnalysed == kTRUE) {
780 const Int_t nTimeBins = fLastTimeBin - fFirstTimeBin +1;
781 const Float_t denormalization = Float_t(fEventCounter * nTimeBins);
782 fNoThreshold->Multiply(denormalization);
784 fMeanCharge->Multiply(fNLocalMaxima);
785 fMaxCharge->Multiply(fNLocalMaxima);
786 fNTimeBins->Multiply(fNLocalMaxima);
787 fNPads->Multiply(fNLocalMaxima);
788 fTimePosition->Multiply(fNLocalMaxima);
789 fIsAnalysed = kFALSE;
795 if (iTimeBin<fFirstTimeBin) return 0;
796 if (iTimeBin>fLastTimeBin) return 0;
798 // if pedestal calibrations are loaded subtract pedestals
801 Float_t ped = fPedestal->GetCalROC(iSector)->GetValue(iRow, iPad);
802 // Don't use data from pads where pedestals are abnormally small or large
808 // In fNoThreshold we fill all data to estimate the ZS volume
809 Float_t count = fNoThreshold->GetCalROC(iSector)->GetValue(iRow, iPad);
810 fNoThreshold->GetCalROC(iSector)->SetValue(iRow, iPad,count+1);
812 // Require at least 3 ADC channels
816 // if noise calibrations are loaded require at least 3*sigmaNoise
819 Float_t noise = fNoise->GetCalROC(iSector)->GetValue(iRow, iPad);
821 if(signal < noise*3.0)
826 // This signal is ok and we store it in the cluster map
829 SetExpandDigit(iRow, iPad, iTimeBin, signal);
833 return 1; // signal was accepted
836 //_____________________________________________________________________
837 void AliTPCdataQA::FindLocalMaxima(const Int_t iSector)
840 // This method is called after the data from each sector has been
841 // exapanded into an array
842 // Loop over the signals and identify local maxima and fill the
843 // calibration objects with the information
846 Int_t nLocalMaxima = 0;
847 const Int_t maxTimeBin = fTimeBinsMax+4; // Used to step between neighboring pads
848 // Because we have tha pad-time data in a
851 for (Int_t iRow = 0; iRow < fRowsMax; iRow++) {
853 Float_t* allBins = fAllBins[iRow];
854 Int_t* sigBins = fAllSigBins[iRow];
855 const Int_t nSigBins = fAllNSigBins[iRow];
857 for (Int_t iSig = 0; iSig < nSigBins; iSig++) {
859 Int_t bin = sigBins[iSig];
860 Float_t *b = &allBins[bin];
863 // Now we check if this is a local maximum
868 // First check that the charge is bigger than the threshold
872 // Require at least one neighboring pad with signal
873 if (b[-maxTimeBin]+b[maxTimeBin]<=0) continue;
875 // Require at least one neighboring time bin with signal
876 if (b[-1]+b[1]<=0) continue;
879 // Check that this is a local maximum
880 // Note that the checking is done so that if 2 charges has the same
881 // qMax then only 1 cluster is generated
882 // (that is why there is BOTH > and >=)
884 if (b[-maxTimeBin] >= qMax) continue;
885 if (b[-1 ] >= qMax) continue;
886 if (b[+maxTimeBin] > qMax) continue;
887 if (b[+1 ] > qMax) continue;
888 if (b[-maxTimeBin-1] >= qMax) continue;
889 if (b[+maxTimeBin-1] >= qMax) continue;
890 if (b[+maxTimeBin+1] > qMax) continue;
891 if (b[-maxTimeBin+1] >= qMax) continue;
894 // Now we accept the local maximum and fill the calibration/data objects
898 Int_t iPad, iTimeBin;
899 GetPadAndTimeBin(bin, iPad, iTimeBin);
901 Float_t count = fNLocalMaxima->GetCalROC(iSector)->GetValue(iRow, iPad);
902 fNLocalMaxima->GetCalROC(iSector)->SetValue(iRow, iPad, count+1);
904 count = fTimePosition->GetCalROC(iSector)->GetValue(iRow, iPad);
905 fTimePosition->GetCalROC(iSector)->SetValue(iRow, iPad, count+iTimeBin);
907 Float_t charge = fMaxCharge->GetCalROC(iSector)->GetValue(iRow, iPad);
908 fMaxCharge->GetCalROC(iSector)->SetValue(iRow, iPad, charge + qMax);
911 count = fOverThreshold10->GetCalROC(iSector)->GetValue(iRow, iPad);
912 fOverThreshold10->GetCalROC(iSector)->SetValue(iRow, iPad, count+1);
915 count = fOverThreshold20->GetCalROC(iSector)->GetValue(iRow, iPad);
916 fOverThreshold20->GetCalROC(iSector)->SetValue(iRow, iPad, count+1);
919 count = fOverThreshold30->GetCalROC(iSector)->GetValue(iRow, iPad);
920 fOverThreshold30->GetCalROC(iSector)->SetValue(iRow, iPad, count+1);
924 // Calculate the total charge as the sum over the region:
932 // with qmax at the center C.
934 // The inner charge (i) we always add, but we only add the outer
935 // charge (o) if the neighboring inner bin (i) has a signal.
937 Int_t minP = 0, maxP = 0, minT = 0, maxT = 0;
939 for(Int_t i = -1; i<=1; i++) {
940 for(Int_t j = -1; j<=1; j++) {
945 Float_t charge1 = GetQ(b, i, j, maxTimeBin, minT, maxT, minP, maxP);
948 // see if the next neighbor is also above threshold
950 qTot += GetQ(b, 2*i, 2*j, maxTimeBin, minT, maxT, minP, maxP);
952 // we are in a diagonal corner
953 qTot += GetQ(b, i, 2*j, maxTimeBin, minT, maxT, minP, maxP);
954 qTot += GetQ(b, 2*i, j, maxTimeBin, minT, maxT, minP, maxP);
955 qTot += GetQ(b, 2*i, 2*j, maxTimeBin, minT, maxT, minP, maxP);
961 charge = fMeanCharge->GetCalROC(iSector)->GetValue(iRow, iPad);
962 fMeanCharge->GetCalROC(iSector)->SetValue(iRow, iPad, charge + qTot);
964 count = fNTimeBins->GetCalROC(iSector)->GetValue(iRow, iPad);
965 fNTimeBins->GetCalROC(iSector)->SetValue(iRow, iPad, count + maxT-minT+1);
967 count = fNPads->GetCalROC(iSector)->GetValue(iRow, iPad);
968 fNPads->GetCalROC(iSector)->SetValue(iRow, iPad, count + maxP-minP+1);
970 if((iSector%36)<18) { // A side
971 fHistQVsTimeSideA->Fill(iTimeBin, qTot);
972 fHistQMaxVsTimeSideA->Fill(iTimeBin, qMax);
974 fHistQVsTimeSideC->Fill(iTimeBin, qTot);
975 fHistQMaxVsTimeSideC->Fill(iTimeBin, qMax);
977 } // end loop over signals
978 } // end loop over rows
980 fClusterCounter += nLocalMaxima;
983 //_____________________________________________________________________
984 void AliTPCdataQA::Analyse()
987 // Calculate calibration constants
990 AliInfo("Analyse called");
992 if(fIsAnalysed == kTRUE) {
994 AliInfo("No new data since Analyse was called last time");
998 if(fEventCounter==0) {
1000 AliInfo("EventCounter == 0, Cannot analyse");
1004 Int_t nTimeBins = fLastTimeBin - fFirstTimeBin +1;
1005 AliInfo(Form("EventCounter: %d , TimeBins: %d", fEventCounter, nTimeBins));
1007 Float_t normalization = 1.0 / Float_t(fEventCounter * nTimeBins);
1008 fNoThreshold->Multiply(normalization);
1010 fMeanCharge->Divide(fNLocalMaxima);
1011 fMaxCharge->Divide(fNLocalMaxima);
1012 fNTimeBins->Divide(fNLocalMaxima);
1013 fNPads->Divide(fNLocalMaxima);
1014 fTimePosition->Divide(fNLocalMaxima);
1016 fIsAnalysed = kTRUE;
1020 //_____________________________________________________________________
1021 void AliTPCdataQA::MakeTree(const char *fname) const {
1023 // Export result to the tree -located in the file
1024 // This file can be analyzed using AliTPCCalibViewer
1026 AliTPCPreprocessorOnline preprocesor;
1028 if (fNLocalMaxima) preprocesor.AddComponent(fNLocalMaxima);
1029 if (fMaxCharge) preprocesor.AddComponent(fMaxCharge);
1030 if (fMeanCharge) preprocesor.AddComponent(fMeanCharge);
1031 if (fNoThreshold) preprocesor.AddComponent(fNoThreshold);
1032 if (fNTimeBins) preprocesor.AddComponent(fNTimeBins);
1033 if (fNPads) preprocesor.AddComponent(fNPads);
1034 if (fTimePosition) preprocesor.AddComponent(fTimePosition);
1035 if (fOverThreshold10) preprocesor.AddComponent(fOverThreshold10);
1036 if (fOverThreshold20) preprocesor.AddComponent(fOverThreshold20);
1037 if (fOverThreshold30) preprocesor.AddComponent(fOverThreshold30);
1039 preprocesor.DumpToFile(fname);
1043 //_____________________________________________________________________
1044 void AliTPCdataQA::MakeArrays(){
1046 // The arrays for expanding the raw data are defined and
1047 // som parameters are intialised
1049 AliTPCROC * roc = AliTPCROC::Instance();
1051 // To make the array big enough for all sectors we take
1052 // the dimensions from the outer row of an OROC (the last sector)
1054 fRowsMax = roc->GetNRows(roc->GetNSector()-1);
1055 fPadsMax = roc->GetNPads(roc->GetNSector()-1,fRowsMax-1);
1056 fTimeBinsMax = fLastTimeBin - fFirstTimeBin +1;
1059 // Since we have added 2 pads (TimeBins) before and after the real pads (TimeBins)
1060 // to make sure that we can always query the exanded table even when the
1061 // max is on the edge
1065 fAllBins = new Float_t*[fRowsMax];
1066 fAllSigBins = new Int_t*[fRowsMax];
1067 fAllNSigBins = new Int_t[fRowsMax];
1069 for (Int_t iRow = 0; iRow < fRowsMax; iRow++) {
1071 Int_t maxBin = (fTimeBinsMax+4)*(fPadsMax+4);
1072 fAllBins[iRow] = new Float_t[maxBin];
1073 memset(fAllBins[iRow],0,sizeof(Float_t)*maxBin); // set all values to 0
1074 fAllSigBins[iRow] = new Int_t[maxBin];
1075 fAllNSigBins[iRow] = 0;
1080 //_____________________________________________________________________
1081 void AliTPCdataQA::CleanArrays(){
1086 for (Int_t iRow = 0; iRow < fRowsMax; iRow++) {
1088 // To speed up the performance by a factor 2 on cosmic data (and
1089 // presumably pp data as well) where the ocupancy is low, the
1090 // memset is only called if there is more than 1000 signals for a
1091 // row (of the order 1% occupancy)
1092 if(fAllNSigBins[iRow]<1000) {
1094 Float_t* allBins = fAllBins[iRow];
1095 Int_t* sigBins = fAllSigBins[iRow];
1096 const Int_t nSignals = fAllNSigBins[iRow];
1097 for(Int_t i = 0; i < nSignals; i++)
1098 allBins[sigBins[i]]=0;
1101 Int_t maxBin = (fTimeBinsMax+4)*(fPadsMax+4);
1102 memset(fAllBins[iRow],0,sizeof(Float_t)*maxBin);
1105 fAllNSigBins[iRow]=0;
1109 //_____________________________________________________________________
1110 void AliTPCdataQA::GetPadAndTimeBin(Int_t bin, Int_t& iPad, Int_t& iTimeBin){
1112 // Return pad and timebin for a given bin
1115 // Int_t bin = iPad*(fTimeBinsMax+4)+iTimeBin;
1116 iTimeBin = bin%(fTimeBinsMax+4);
1117 iPad = (bin-iTimeBin)/(fTimeBinsMax+4);
1121 iTimeBin += fFirstTimeBin;
1123 R__ASSERT(iPad>=0 && iPad<=fPadsMax);
1124 R__ASSERT(iTimeBin>=fFirstTimeBin && iTimeBin<=fLastTimeBin);
1127 //_____________________________________________________________________
1128 void AliTPCdataQA::SetExpandDigit(const Int_t iRow, Int_t iPad,
1129 Int_t iTimeBin, const Float_t signal)
1134 R__ASSERT(iRow>=0 && iRow<fRowsMax);
1135 R__ASSERT(iPad>=0 && iPad<=fPadsMax);
1136 R__ASSERT(iTimeBin>=fFirstTimeBin && iTimeBin<=fLastTimeBin);
1138 iTimeBin -= fFirstTimeBin;
1142 Int_t bin = iPad*(fTimeBinsMax+4)+iTimeBin;
1143 fAllBins[iRow][bin] = signal;
1144 fAllSigBins[iRow][fAllNSigBins[iRow]] = bin;
1145 fAllNSigBins[iRow]++;
1148 Float_t AliTPCdataQA::GetQ(const Float_t* adcArray, const Int_t time,
1149 const Int_t pad, const Int_t maxTimeBins,
1150 Int_t& timeMin, Int_t& timeMax,
1151 Int_t& padMin, Int_t& padMax) const
1154 // This methods return the charge in the bin time+pad*maxTimeBins
1155 // If the charge is above 0 it also updates the padMin, padMax, timeMin
1156 // and timeMax if necessary
1158 Float_t charge = adcArray[time + pad*maxTimeBins];
1160 timeMin = TMath::Min(time, timeMin); timeMax = TMath::Max(time, timeMax);
1161 padMin = TMath::Min(pad, padMin); padMax = TMath::Max(pad, padMax);
1166 //______________________________________________________________________________
1167 void AliTPCdataQA::Streamer(TBuffer &xRuub)
1169 // Automatic schema evolution was first used from revision 4
1171 // http://root.cern.ch/root/roottalk/roottalk02/3207.html
1173 UInt_t xRuus, xRuuc;
1174 if (xRuub.IsReading()) {
1175 Version_t xRuuv = xRuub.ReadVersion(&xRuus, &xRuuc);
1176 //we use the automatic algorithm for class version > 3
1178 AliTPCdataQA::Class()->ReadBuffer(xRuub, this, xRuuv, xRuus,
1182 TH1F::Streamer(xRuub);
1183 xRuub >> fFirstTimeBin;
1184 xRuub >> fLastTimeBin;
1187 xRuub >> fNLocalMaxima;
1188 xRuub >> fMaxCharge;
1189 xRuub >> fMeanCharge;
1190 xRuub >> fNoThreshold;
1191 xRuub >> fNTimeBins;
1193 xRuub >> fTimePosition;
1194 xRuub >> fEventCounter;
1195 xRuub >> fIsAnalysed;
1196 xRuub.CheckByteCount(xRuus, xRuuc, AliTPCdataQA::IsA());
1198 AliTPCdataQA::Class()->WriteBuffer(xRuub,this);