]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TPC/AliTPCCalibPulser.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Updated macro for ACORDE raw-data visualisation. We use now the ACORDE sym.names...
[u/mrichter/AliRoot.git] / TPC / AliTPCCalibPulser.cxx
CommitLineData
467a4b25 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/* $Id$ */
17
3cd27a08 18/////////////////////////////////////////////////////////////////////////////////////////
19// //
20// Implementation of the TPC pulser calibration //
21// //
22// Origin: Jens Wiechula, Marian Ivanov J.Wiechula@gsi.de, Marian.Ivanov@cern.ch //
23// //
24/////////////////////////////////////////////////////////////////////////////////////////
467a4b25 25/***************************************************************************
26 * Class Description *
27 ***************************************************************************
28
29 The AliTPCCalibPulser class is used to get calibration data concerning the FEE using
30 runs performed with the calibration pulser.
31
32 The information retrieved is
33 - Time0 differences
34 - Signal width differences
35 - Amplification variations
36
37 the seen differences arise from the manufacturing tolerances of the PASAs and are very small within
38 one chip and somewhat large between different chips.
39
40 Histograms:
41 For each ROC three TH2S histos 'Reference Histograms' (ROC channel vs. [Time0, signal width, Q sum]) is created when
42 it is filled for the first time (GetHisto[T0,RMS,Q](ROC,kTRUE)). The histos are stored in the
43 TObjArrays fHistoT0Array, fHistoRMSArray and fHistoQArray.
44
45
46 Working principle:
47 ------------------
48 Raw calibration pulser data is processed by calling one of the ProcessEvent(...) functions
49 (see below). These in the end call the Update(...) function.
50
51 - the Update(...) function:
52 In this function the array fPadSignal is filled with the adc signals between the specified range
53 fFirstTimeBin and fLastTimeBin for the current pad.
54 before going to the next pad the ProcessPad() function is called, which analyses the data for one pad
55 stored in fPadSignal.
56
57 - the ProcessPad() function:
58 Find Pedestal and Noise information
59 - use database information which has to be set by calling
60 SetPedestalDatabase(AliTPCCalPad *pedestalTPC, AliTPCCalPad *padNoiseTPC)
61 - if no information from the pedestal data base
62 is available the informaion is calculated on the fly ( see FindPedestal() function )
63
64 Find the Pulser signal information
65 - calculate mean = T0, RMS = signal width and Q sum in a range of -2+7 timebins around Q max
66 the Q sum is scaled by pad area
67 (see FindPulserSignal(...) function)
68
69 Fill a temprary array for the T0 information (GetPadTimesEvent(fCurrentSector,kTRUE)) (why see below)
70 Fill the Q sum and RMS values in the histograms (GetHisto[RMS,Q](ROC,kTRUE)),
71
72 At the end of each event the EndEvent() function is called
73
74 - the EndEvent() function:
75 calculate the mean T0 for each ROC and fill the Time0 histogram with Time0-<Time0 for ROC>
76 This is done to overcome syncronisation problems between the trigger and the fec clock.
77
78 After accumulating the desired statistics the Analyse() function has to be called.
79 - the Analyse() function
80 Whithin this function the mean values of T0, RMS, Q are calculated for each pad, using
81 the AliMathBase::GetCOG(...) function, and the calibration
82 storage classes (AliTPCCalROC) are filled for each ROC.
83 The calibration information is stored in the TObjArrays fCalRocArrayT0, fCalRocArrayRMS and
84 fCalRocArrayQ;
85
86
87
88 User interface for filling data:
89 --------------------------------
90
91 To Fill information one of the following functions can be used:
92
93 Bool_t ProcessEvent(eventHeaderStruct *event);
94 - process Date event
95 - use AliTPCRawReaderDate and call ProcessEvent(AliRawReader *rawReader)
96
97 Bool_t ProcessEvent(AliRawReader *rawReader);
98 - process AliRawReader event
99 - use AliTPCRawStream to loop over data and call ProcessEvent(AliTPCRawStream *rawStream)
100
101 Bool_t ProcessEvent(AliTPCRawStream *rawStream);
102 - process event from AliTPCRawStream
103 - call Update function for signal filling
104
105 Int_t Update(const Int_t isector, const Int_t iRow, const Int_t
106 iPad, const Int_t iTimeBin, const Float_t signal);
107 - directly fill signal information (sector, row, pad, time bin, pad)
108 to the reference histograms
109
110 It is also possible to merge two independently taken calibrations using the function
111
112 void Merge(AliTPCCalibPulser *sig)
113 - copy histograms in 'sig' if the do not exist in this instance
114 - Add histograms in 'sig' to the histograms in this instance if the allready exist
115 - After merging call Analyse again!
116
117
118
119 -- example: filling data using root raw data:
120 void fillSignal(Char_t *filename)
121 {
122 rawReader = new AliRawReaderRoot(fileName);
123 if ( !rawReader ) return;
124 AliTPCCalibPulser *calib = new AliTPCCalibPulser;
125 while (rawReader->NextEvent()){
126 calib->ProcessEvent(rawReader);
127 }
128 calib->Analyse();
129 calib->DumpToFile("SignalData.root");
130 delete rawReader;
131 delete calib;
132 }
133
134
135 What kind of information is stored and how to retrieve them:
136 ------------------------------------------------------------
137
138 - Accessing the 'Reference Histograms' (Time0, signal width and Q sum information pad by pad):
139
140 TH2F *GetHistoT0(Int_t sector);
141 TH2F *GetHistoRMS(Int_t sector);
142 TH2F *GetHistoQ(Int_t sector);
143
144 - Accessing the calibration storage objects:
145
146 AliTPCCalROC *GetCalRocT0(Int_t sector); // for the Time0 values
147 AliTPCCalROC *GetCalRocRMS(Int_t sector); // for the signal width values
148 AliTPCCalROC *GetCalRocQ(Int_t sector); // for the Q sum values
149
150 example for visualisation:
151 if the file "SignalData.root" was created using the above example one could do the following:
152
153 TFile fileSignal("SignalData.root")
154 AliTPCCalibPulser *sig = (AliTPCCalibPulser*)fileSignal->Get("AliTPCCalibPulser");
155 sig->GetCalRocT0(0)->Draw("colz");
156 sig->GetCalRocRMS(0)->Draw("colz");
157
158 or use the AliTPCCalPad functionality:
159 AliTPCCalPad padT0(ped->GetCalPadT0());
160 AliTPCCalPad padSigWidth(ped->GetCalPadRMS());
161 padT0->MakeHisto2D()->Draw("colz"); //Draw A-Side Time0 Information
162 padSigWidth->MakeHisto2D()->Draw("colz"); //Draw A-Side signal width Information
163*/
164
165
3cd27a08 166//Root includes
167#include <TObjArray.h>
168#include <TH1F.h>
169#include <TH2S.h>
170#include <TString.h>
171#include <TVectorF.h>
172#include <TMath.h>
173
174#include <TDirectory.h>
175#include <TSystem.h>
176#include <TROOT.h>
177#include <TFile.h>
178
179//AliRoot includes
180#include "AliRawReader.h"
181#include "AliRawReaderRoot.h"
182#include "AliRawReaderDate.h"
183#include "AliTPCRawStream.h"
08205ed7 184#include "AliTPCRawStreamFast.h"
3cd27a08 185#include "AliTPCCalROC.h"
186#include "AliTPCCalPad.h"
187#include "AliTPCROC.h"
188#include "AliTPCParam.h"
189#include "AliTPCCalibPulser.h"
190#include "AliTPCcalibDB.h"
191#include "AliMathBase.h"
d8faeea7 192#include "AliLog.h"
3cd27a08 193#include "TTreeStream.h"
467a4b25 194
3cd27a08 195//date
196#include "event.h"
197
198
199
200
201ClassImp(AliTPCCalibPulser)
202
203AliTPCCalibPulser::AliTPCCalibPulser() :
467a4b25 204 TObject(),
205 fFirstTimeBin(60),
206 fLastTimeBin(120),
207 fNbinsT0(200),
208 fXminT0(-2),
209 fXmaxT0(2),
210 fNbinsQ(200),
211 fXminQ(1),
212 fXmaxQ(40),
213 fNbinsRMS(100),
214 fXminRMS(0.1),
215 fXmaxRMS(5.1),
216 fLastSector(-1),
467a4b25 217 fROC(AliTPCROC::Instance()),
04049c81 218 fMapping(NULL),
467a4b25 219 fParam(new AliTPCParam),
220 fPedestalTPC(0x0),
221 fPadNoiseTPC(0x0),
d8faeea7 222 fOutliers(0x0),
467a4b25 223 fPedestalROC(0x0),
224 fPadNoiseROC(0x0),
225 fCalRocArrayT0(72),
226 fCalRocArrayQ(72),
227 fCalRocArrayRMS(72),
228 fCalRocArrayOutliers(72),
229 fHistoQArray(72),
230 fHistoT0Array(72),
231 fHistoRMSArray(72),
232 fPadTimesArrayEvent(72),
233 fPadQArrayEvent(72),
234 fPadRMSArrayEvent(72),
235 fPadPedestalArrayEvent(72),
236 fCurrentChannel(-1),
237 fCurrentSector(-1),
238 fCurrentRow(-1),
239 fMaxPadSignal(-1),
240 fMaxTimeBin(-1),
241 fPadSignal(1024),
242 fPadPedestal(0),
243 fPadNoise(0),
244 fVTime0Offset(72),
245 fVTime0OffsetCounter(72),
3cd27a08 246// fEvent(-1),
467a4b25 247 fDebugStreamer(0x0),
248 fDebugLevel(0)
249{
250 //
251 // AliTPCSignal default constructor
252 //
253
254}
255//_____________________________________________________________________
256AliTPCCalibPulser::AliTPCCalibPulser(const AliTPCCalibPulser &sig) :
257 TObject(sig),
258 fFirstTimeBin(sig.fFirstTimeBin),
259 fLastTimeBin(sig.fLastTimeBin),
260 fNbinsT0(sig.fNbinsT0),
261 fXminT0(sig.fXminT0),
262 fXmaxT0(sig.fXmaxT0),
263 fNbinsQ(sig.fNbinsQ),
264 fXminQ(sig.fXminQ),
265 fXmaxQ(sig.fXmaxQ),
266 fNbinsRMS(sig.fNbinsRMS),
267 fXminRMS(sig.fXminRMS),
268 fXmaxRMS(sig.fXmaxRMS),
269 fLastSector(-1),
467a4b25 270 fROC(AliTPCROC::Instance()),
04049c81 271 fMapping(NULL),
467a4b25 272 fParam(new AliTPCParam),
273 fPedestalTPC(0x0),
274 fPadNoiseTPC(0x0),
d8faeea7 275 fOutliers(0x0),
467a4b25 276 fPedestalROC(0x0),
277 fPadNoiseROC(0x0),
278 fCalRocArrayT0(72),
279 fCalRocArrayQ(72),
280 fCalRocArrayRMS(72),
281 fCalRocArrayOutliers(72),
282 fHistoQArray(72),
283 fHistoT0Array(72),
284 fHistoRMSArray(72),
285 fPadTimesArrayEvent(72),
286 fPadQArrayEvent(72),
287 fPadRMSArrayEvent(72),
288 fPadPedestalArrayEvent(72),
289 fCurrentChannel(-1),
290 fCurrentSector(-1),
291 fCurrentRow(-1),
292 fMaxPadSignal(-1),
293 fMaxTimeBin(-1),
294 fPadSignal(1024),
295 fPadPedestal(0),
296 fPadNoise(0),
297 fVTime0Offset(72),
298 fVTime0OffsetCounter(72),
3cd27a08 299// fEvent(-1),
467a4b25 300 fDebugStreamer(0x0),
301 fDebugLevel(sig.fDebugLevel)
302{
303 //
304 // AliTPCSignal default constructor
305 //
306
307 for (Int_t iSec = 0; iSec < 72; ++iSec){
308 const AliTPCCalROC *calQ = (AliTPCCalROC*)sig.fCalRocArrayQ.UncheckedAt(iSec);
309 const AliTPCCalROC *calT0 = (AliTPCCalROC*)sig.fCalRocArrayT0.UncheckedAt(iSec);
310 const AliTPCCalROC *calRMS = (AliTPCCalROC*)sig.fCalRocArrayRMS.UncheckedAt(iSec);
311 const AliTPCCalROC *calOut = (AliTPCCalROC*)sig.fCalRocArrayOutliers.UncheckedAt(iSec);
312
313 const TH2S *hQ = (TH2S*)sig.fHistoQArray.UncheckedAt(iSec);
314 const TH2S *hT0 = (TH2S*)sig.fHistoT0Array.UncheckedAt(iSec);
315 const TH2S *hRMS = (TH2S*)sig.fHistoRMSArray.UncheckedAt(iSec);
316
317 if ( calQ != 0x0 ) fCalRocArrayQ.AddAt(new AliTPCCalROC(*calQ), iSec);
318 if ( calT0 != 0x0 ) fCalRocArrayT0.AddAt(new AliTPCCalROC(*calT0), iSec);
319 if ( calRMS != 0x0 ) fCalRocArrayRMS.AddAt(new AliTPCCalROC(*calRMS), iSec);
320 if ( calOut != 0x0 ) fCalRocArrayOutliers.AddAt(new AliTPCCalROC(*calOut), iSec);
321
322 if ( hQ != 0x0 ){
323 TH2S *hNew = new TH2S(*hQ);
324 hNew->SetDirectory(0);
325 fHistoQArray.AddAt(hNew,iSec);
326 }
327 if ( hT0 != 0x0 ){
328 TH2S *hNew = new TH2S(*hT0);
329 hNew->SetDirectory(0);
330 fHistoQArray.AddAt(hNew,iSec);
331 }
332 if ( hRMS != 0x0 ){
333 TH2S *hNew = new TH2S(*hRMS);
334 hNew->SetDirectory(0);
335 fHistoQArray.AddAt(hNew,iSec);
336 }
337 }
338
339}
340//_____________________________________________________________________
341AliTPCCalibPulser& AliTPCCalibPulser::operator = (const AliTPCCalibPulser &source)
342{
343 //
344 // assignment operator
345 //
346 if (&source == this) return *this;
347 new (this) AliTPCCalibPulser(source);
348
349 return *this;
350}
351//_____________________________________________________________________
352AliTPCCalibPulser::~AliTPCCalibPulser()
353{
354 //
355 // destructor
356 //
357
d8faeea7 358 Reset();
359
360 if ( fDebugStreamer) delete fDebugStreamer;
361 delete fROC;
362 delete fParam;
363}
364void AliTPCCalibPulser::Reset()
365{
366 //
367 // Delete all information: Arrays, Histograms, CalRoc objects
368 //
467a4b25 369 fCalRocArrayT0.Delete();
370 fCalRocArrayQ.Delete();
371 fCalRocArrayRMS.Delete();
372 fCalRocArrayOutliers.Delete();
373
374 fHistoQArray.Delete();
375 fHistoT0Array.Delete();
376 fHistoRMSArray.Delete();
377
378 fPadTimesArrayEvent.Delete();
379 fPadQArrayEvent.Delete();
380 fPadRMSArrayEvent.Delete();
381 fPadPedestalArrayEvent.Delete();
467a4b25 382}
383//_____________________________________________________________________
3cd27a08 384Int_t AliTPCCalibPulser::Update(const Int_t icsector,
467a4b25 385 const Int_t icRow,
386 const Int_t icPad,
387 const Int_t icTimeBin,
388 const Float_t csignal)
389{
390 //
3cd27a08 391 // Signal filling methode on the fly pedestal and time offset correction if necessary.
467a4b25 392 // no extra analysis necessary. Assumes knowledge of the signal shape!
393 // assumes that it is looped over consecutive time bins of one pad
394 //
b401648b 395
396 if (icRow<0) return 0;
397 if (icPad<0) return 0;
398 if (icTimeBin<0) return 0;
467a4b25 399 if ( (icTimeBin>fLastTimeBin) || (icTimeBin<fFirstTimeBin) ) return 0;
400
d8faeea7 401 if ( icRow<0 || icPad<0 ){
402 AliWarning("Wrong Pad or Row number, skipping!");
403 return 0;
404 }
405
467a4b25 406 Int_t iChannel = fROC->GetRowIndexes(icsector)[icRow]+icPad; // global pad position in sector
407
408 //init first pad and sector in this event
409 if ( fCurrentChannel == -1 ) {
410 fCurrentChannel = iChannel;
411 fCurrentSector = icsector;
412 fCurrentRow = icRow;
413 }
414
415 //process last pad if we change to a new one
416 if ( iChannel != fCurrentChannel ){
417 ProcessPad();
418 fCurrentChannel = iChannel;
419 fCurrentSector = icsector;
420 fCurrentRow = icRow;
421 }
422
423 //fill signals for current pad
424 fPadSignal[icTimeBin]=csignal;
425 if ( csignal > fMaxPadSignal ){
426 fMaxPadSignal = csignal;
427 fMaxTimeBin = icTimeBin;
428 }
429 return 0;
430}
431//_____________________________________________________________________
432void AliTPCCalibPulser::FindPedestal(Float_t part)
433{
434 //
435 // find pedestal and noise for the current pad. Use either database or
436 // truncated mean with part*100%
437 //
438 Bool_t noPedestal = kTRUE;;
439 if (fPedestalTPC&&fPadNoiseTPC){
440 //use pedestal database
441 //only load new pedestals if the sector has changed
442 if ( fCurrentSector!=fLastSector ){
443 fPedestalROC = fPedestalTPC->GetCalROC(fCurrentSector);
444 fPadNoiseROC = fPadNoiseTPC->GetCalROC(fCurrentSector);
445 fLastSector=fCurrentSector;
446 }
447
448 if ( fPedestalROC&&fPadNoiseROC ){
449 fPadPedestal = fPedestalROC->GetValue(fCurrentChannel);
450 fPadNoise = fPadNoiseROC->GetValue(fCurrentChannel);
451 noPedestal = kFALSE;
452 }
453
454 }
455
456 //if we are not running with pedestal database, or for the current sector there is no information
457 //available, calculate the pedestal and noise on the fly
458 if ( noPedestal ) {
459 const Int_t kPedMax = 100; //maximum pedestal value
460 Float_t max = 0;
461 Float_t maxPos = 0;
462 Int_t median = -1;
463 Int_t count0 = 0;
464 Int_t count1 = 0;
465 //
466 Float_t padSignal=0;
467 //
468 UShort_t histo[kPedMax];
469 memset(histo,0,kPedMax*sizeof(UShort_t));
470
471 for (Int_t i=fFirstTimeBin; i<=fLastTimeBin; ++i){
472 padSignal = fPadSignal.GetMatrixArray()[i];
473 if (padSignal<=0) continue;
474 if (padSignal>max && i>10) {
475 max = padSignal;
476 maxPos = i;
477 }
478 if (padSignal>kPedMax-1) continue;
479 histo[Int_t(padSignal+0.5)]++;
480 count0++;
481 }
482 //
483 for (Int_t i=1; i<kPedMax; ++i){
484 if (count1<count0*0.5) median=i;
485 count1+=histo[i];
486 }
487 // truncated mean
488 //
d8faeea7 489 // what if by chance histo[median] == 0 ?!?
467a4b25 490 Float_t count=histo[median] ,mean=histo[median]*median, rms=histo[median]*median*median ;
491 //
492 for (Int_t idelta=1; idelta<10; ++idelta){
493 if (median-idelta<=0) continue;
494 if (median+idelta>kPedMax) continue;
495 if (count<part*count1){
496 count+=histo[median-idelta];
497 mean +=histo[median-idelta]*(median-idelta);
498 rms +=histo[median-idelta]*(median-idelta)*(median-idelta);
499 count+=histo[median+idelta];
500 mean +=histo[median+idelta]*(median+idelta);
501 rms +=histo[median+idelta]*(median+idelta)*(median+idelta);
502 }
503 }
504 fPadPedestal = 0;
505 fPadNoise = 0;
506 if ( count > 0 ) {
507 mean/=count;
508 rms = TMath::Sqrt(TMath::Abs(rms/count-mean*mean));
509 fPadPedestal = mean;
510 fPadNoise = rms;
511 }
512 }
513}
514//_____________________________________________________________________
515void AliTPCCalibPulser::FindPulserSignal(TVectorD &param, Float_t &qSum)
516{
517 //
518 // Find position, signal width and height of the CE signal (last signal)
519 // param[0] = Qmax, param[1] = mean time, param[2] = rms;
520 // maxima: array of local maxima of the pad signal use the one closest to the mean CE position
521 //
522
523 Float_t ceQmax =0, ceQsum=0, ceTime=0, ceRMS=0;
524 Int_t cemaxpos = fMaxTimeBin;
d8faeea7 525 Float_t ceSumThreshold = 10.*TMath::Max(fPadNoise,Float_t(1.)); // threshold for the signal sum
526 Float_t ceMaxThreshold = 5.*TMath::Max(fPadNoise,Float_t(1.)); // threshold for the signal max
467a4b25 527 const Int_t kCemin = 2; // range for the analysis of the ce signal +- channels from the peak
528 const Int_t kCemax = 7;
d8faeea7 529 param[0] = ceQmax;
530 param[1] = ceTime;
531 param[2] = ceRMS;
532 qSum = ceQsum;
467a4b25 533
d8faeea7 534 if (cemaxpos>0){
535 ceQmax = fPadSignal.GetMatrixArray()[cemaxpos]-fPadPedestal;
536 if ( ceQmax<ceMaxThreshold ) return;
467a4b25 537 for (Int_t i=cemaxpos-kCemin; i<cemaxpos+kCemax; ++i){
538 Float_t signal = fPadSignal.GetMatrixArray()[i]-fPadPedestal;
539 if ( (i>fFirstTimeBin) && (i<fLastTimeBin) && (signal>0) ){
540 ceTime+=signal*(i+0.5);
541 ceRMS +=signal*(i+0.5)*(i+0.5);
542 ceQsum+=signal;
543 }
544 }
545 }
d8faeea7 546 if (ceQsum>ceSumThreshold) {
467a4b25 547 ceTime/=ceQsum;
548 ceRMS = TMath::Sqrt(TMath::Abs(ceRMS/ceQsum-ceTime*ceTime));
d8faeea7 549 //only fill the Time0Offset if pad was not marked as an outlier!
550 if ( !fOutliers ){
551 fVTime0Offset.GetMatrixArray()[fCurrentSector]+=ceTime; // mean time for each sector
552 fVTime0OffsetCounter.GetMatrixArray()[fCurrentSector]++;
553 } else {
554 if ( !(fOutliers->GetCalROC(fCurrentSector)->GetValue(fCurrentChannel)) ){
555 fVTime0Offset.GetMatrixArray()[fCurrentSector]+=ceTime; // mean time for each sector
556 fVTime0OffsetCounter.GetMatrixArray()[fCurrentSector]++;
557 }
558 }
467a4b25 559
560 //Normalise Q to the pad area
561 Float_t norm = fParam->GetPadPitchWidth(fCurrentSector)*fParam->GetPadPitchLength(fCurrentSector,fCurrentRow);
562
563 ceQsum/=norm;
564 } else {
565 ceQmax=0;
566 ceTime=0;
567 ceRMS =0;
568 ceQsum=0;
569 }
570 param[0] = ceQmax;
571 param[1] = ceTime;
572 param[2] = ceRMS;
573 qSum = ceQsum;
574}
575//_____________________________________________________________________
3cd27a08 576void AliTPCCalibPulser::ProcessPad()
467a4b25 577{
578 //
579 // Process data of current pad
580 //
581
582 FindPedestal();
583 TVectorD param(3);
3cd27a08 584 Float_t qSum;
585 FindPulserSignal(param, qSum);
467a4b25 586
587 Double_t meanT = param[1];
588 Double_t sigmaT = param[2];
589
d8faeea7 590
467a4b25 591 //Fill Event T0 counter
592 (*GetPadTimesEvent(fCurrentSector,kTRUE)).GetMatrixArray()[fCurrentChannel] = meanT;
593
594 //Fill Q histogram
3cd27a08 595 GetHistoQ(fCurrentSector,kTRUE)->Fill( TMath::Sqrt(qSum), fCurrentChannel );
467a4b25 596
597 //Fill RMS histogram
598 GetHistoRMS(fCurrentSector,kTRUE)->Fill( sigmaT, fCurrentChannel );
599
600
601 //Fill debugging info
602 if ( fDebugLevel>0 ){
603 (*GetPadPedestalEvent(fCurrentSector,kTRUE))[fCurrentChannel]=fPadPedestal;
604 (*GetPadRMSEvent(fCurrentSector,kTRUE))[fCurrentChannel]=sigmaT;
3cd27a08 605 (*GetPadQEvent(fCurrentSector,kTRUE))[fCurrentChannel]=qSum;
467a4b25 606 }
467a4b25 607 ResetPad();
608}
609//_____________________________________________________________________
3cd27a08 610void AliTPCCalibPulser::EndEvent()
467a4b25 611{
612 //
613 // Process data of current event
614 //
615
616 //check if last pad has allready been processed, if not do so
617 if ( fMaxTimeBin>-1 ) ProcessPad();
618
619 //loop over all ROCs, fill Time0 histogram corrected for the mean Time0 of each ROC
620 for ( Int_t iSec = 0; iSec<72; ++iSec ){
621 TVectorF *vTimes = GetPadTimesEvent(iSec);
f0312c5d 622 if ( !vTimes || fVTime0OffsetCounter[iSec]==0 ) continue;
623 Float_t time0 = fVTime0Offset[iSec]/fVTime0OffsetCounter[iSec];
467a4b25 624 for ( UInt_t iChannel=0; iChannel<fROC->GetNChannels(iSec); ++iChannel ){
3cd27a08 625 Float_t time = (*vTimes).GetMatrixArray()[iChannel];
467a4b25 626
3cd27a08 627 GetHistoT0(iSec,kTRUE)->Fill( time-time0,iChannel );
d8faeea7 628 //GetHistoT0(iSec,kTRUE)->Fill( time,iChannel );
467a4b25 629
630
631 //Debug start
632 if ( fDebugLevel>0 ){
633 if ( !fDebugStreamer ) {
634 //debug stream
635 TDirectory *backup = gDirectory;
636 fDebugStreamer = new TTreeSRedirector("deb2.root");
637 if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer
638 }
639
640 Int_t row=0;
641 Int_t pad=0;
642 Int_t padc=0;
643
3cd27a08 644 Float_t q = (*GetPadQEvent(iSec)).GetMatrixArray()[iChannel];
645 Float_t rms = (*GetPadRMSEvent(iSec)).GetMatrixArray()[iChannel];
467a4b25 646
647 UInt_t channel=iChannel;
648 Int_t sector=iSec;
649
650 while ( channel > (fROC->GetRowIndexes(sector)[row]+fROC->GetNPads(sector,row)-1) ) row++;
651 pad = channel-fROC->GetRowIndexes(sector)[row];
652 padc = pad-(fROC->GetNPads(sector,row)/2);
653
654 TH1F *h1 = new TH1F(Form("hSignalD%d.%d.%d",sector,row,pad),
655 Form("hSignalD%d.%d.%d",sector,row,pad),
656 fLastTimeBin-fFirstTimeBin,
657 fFirstTimeBin,fLastTimeBin);
658 h1->SetDirectory(0);
659
660 for (Int_t i=fFirstTimeBin; i<fLastTimeBin+1; ++i)
661 h1->Fill(i,fPadSignal(i));
662
663 (*fDebugStreamer) << "DataPad" <<
3cd27a08 664// "Event=" << fEvent <<
467a4b25 665 "Sector="<< sector <<
666 "Row=" << row<<
667 "Pad=" << pad <<
668 "PadC=" << padc <<
669 "PadSec="<< channel <<
3cd27a08 670 "Time0=" << time0 <<
671 "Time=" << time <<
672 "RMS=" << rms <<
673 "Sum=" << q <<
467a4b25 674 "hist.=" << h1 <<
675 "\n";
676
677 delete h1;
678 }
679 //Debug end
467a4b25 680 }
681 }
08205ed7 682}
683//_____________________________________________________________________
684Bool_t AliTPCCalibPulser::ProcessEventFast(AliTPCRawStreamFast *rawStreamFast)
685{
686 //
687 // Event Processing loop - AliTPCRawStream
688 //
689 ResetEvent();
467a4b25 690
08205ed7 691 Bool_t withInput = kFALSE;
692
693 while ( rawStreamFast->NextDDL() ){
694 while ( rawStreamFast->NextChannel() ){
695 Int_t isector = rawStreamFast->GetSector(); // current sector
696 Int_t iRow = rawStreamFast->GetRow(); // current row
697 Int_t iPad = rawStreamFast->GetPad(); // current pad
08205ed7 698
699 while ( rawStreamFast->NextBunch() ){
5b549198 700 Int_t startTbin = (Int_t)rawStreamFast->GetStartTimeBin();
701 Int_t endTbin = (Int_t)rawStreamFast->GetEndTimeBin();
08205ed7 702 for (Int_t iTimeBin = startTbin; iTimeBin < endTbin; iTimeBin++){
703 Float_t signal=(Float_t)rawStreamFast->GetSignals()[iTimeBin-startTbin];
704 Update(isector,iRow,iPad,iTimeBin+1,signal);
705 withInput = kTRUE;
706 }
707 }
708 }
709 }
710 if (withInput){
711 EndEvent();
712 }
713 return withInput;
714}
715//_____________________________________________________________________
716Bool_t AliTPCCalibPulser::ProcessEventFast(AliRawReader *rawReader)
717{
718 //
719 // Event processing loop - AliRawReader
720 //
04049c81 721 AliTPCRawStreamFast *rawStreamFast = new AliTPCRawStreamFast(rawReader, (AliAltroMapping**)fMapping);
08205ed7 722 Bool_t res=ProcessEventFast(rawStreamFast);
723 delete rawStreamFast;
724 return res;
467a4b25 725}
726//_____________________________________________________________________
3cd27a08 727Bool_t AliTPCCalibPulser::ProcessEvent(AliTPCRawStream *rawStream)
467a4b25 728{
729 //
730 // Event Processing loop - AliTPCRawStream
731 //
732
467a4b25 733 ResetEvent();
734
735 Bool_t withInput = kFALSE;
736
737 while (rawStream->Next()) {
467a4b25 738 Int_t isector = rawStream->GetSector(); // current sector
739 Int_t iRow = rawStream->GetRow(); // current row
740 Int_t iPad = rawStream->GetPad(); // current pad
741 Int_t iTimeBin = rawStream->GetTime(); // current time bin
742 Float_t signal = rawStream->GetSignal(); // current ADC signal
743
744 Update(isector,iRow,iPad,iTimeBin,signal);
745 withInput = kTRUE;
746 }
467a4b25 747 if (withInput){
748 EndEvent();
749 }
467a4b25 750 return withInput;
751}
752//_____________________________________________________________________
753Bool_t AliTPCCalibPulser::ProcessEvent(AliRawReader *rawReader)
754{
755 //
756 // Event processing loop - AliRawReader
757 //
758
759
04049c81 760 AliTPCRawStream rawStream(rawReader, (AliAltroMapping**)fMapping);
467a4b25 761
762 rawReader->Select("TPC");
763
764 return ProcessEvent(&rawStream);
765}
766//_____________________________________________________________________
767Bool_t AliTPCCalibPulser::ProcessEvent(eventHeaderStruct *event)
768{
769 //
770 // Event processing loop - date event
771 //
772 AliRawReader *rawReader = new AliRawReaderDate((void*)event);
773 Bool_t result=ProcessEvent(rawReader);
774 delete rawReader;
775 return result;
776
777}
778//_____________________________________________________________________
3cd27a08 779TH2S* AliTPCCalibPulser::GetHisto(Int_t sector, TObjArray *arr,
467a4b25 780 Int_t nbinsY, Float_t ymin, Float_t ymax,
781 Char_t *type, Bool_t force)
782{
783 //
784 // return pointer to Q histogram
785 // if force is true create a new histogram if it doesn't exist allready
786 //
787 if ( !force || arr->UncheckedAt(sector) )
788 return (TH2S*)arr->UncheckedAt(sector);
789
790 // if we are forced and histogram doesn't yes exist create it
791 Char_t name[255], title[255];
792
793 sprintf(name,"hCalib%s%.2d",type,sector);
794 sprintf(title,"%s calibration histogram sector %.2d",type,sector);
795
796 // new histogram with Q calib information. One value for each pad!
797 TH2S* hist = new TH2S(name,title,
798 nbinsY, ymin, ymax,
799 fROC->GetNChannels(sector),0,fROC->GetNChannels(sector));
800 hist->SetDirectory(0);
801 arr->AddAt(hist,sector);
802 return hist;
803}
804//_____________________________________________________________________
3cd27a08 805TH2S* AliTPCCalibPulser::GetHistoT0(Int_t sector, Bool_t force)
467a4b25 806{
807 //
808 // return pointer to T0 histogram
809 // if force is true create a new histogram if it doesn't exist allready
810 //
811 TObjArray *arr = &fHistoT0Array;
812 return GetHisto(sector, arr, fNbinsT0, fXminT0, fXmaxT0, "T0", force);
813}
814//_____________________________________________________________________
3cd27a08 815TH2S* AliTPCCalibPulser::GetHistoQ(Int_t sector, Bool_t force)
467a4b25 816{
817 //
818 // return pointer to Q histogram
819 // if force is true create a new histogram if it doesn't exist allready
820 //
821 TObjArray *arr = &fHistoQArray;
822 return GetHisto(sector, arr, fNbinsQ, fXminQ, fXmaxQ, "Q", force);
823}
824//_____________________________________________________________________
3cd27a08 825TH2S* AliTPCCalibPulser::GetHistoRMS(Int_t sector, Bool_t force)
467a4b25 826{
827 //
828 // return pointer to Q histogram
829 // if force is true create a new histogram if it doesn't exist allready
830 //
831 TObjArray *arr = &fHistoRMSArray;
832 return GetHisto(sector, arr, fNbinsRMS, fXminRMS, fXmaxRMS, "RMS", force);
833}
834//_____________________________________________________________________
3cd27a08 835TVectorF* AliTPCCalibPulser::GetPadInfoEvent(Int_t sector, TObjArray *arr, Bool_t force)
467a4b25 836{
837 //
838 // return pointer to Pad Info from 'arr' for the current event and sector
839 // if force is true create it if it doesn't exist allready
840 //
841 if ( !force || arr->UncheckedAt(sector) )
842 return (TVectorF*)arr->UncheckedAt(sector);
843
844 TVectorF *vect = new TVectorF(fROC->GetNChannels(sector));
845 arr->AddAt(vect,sector);
846 return vect;
847}
848//_____________________________________________________________________
3cd27a08 849TVectorF* AliTPCCalibPulser::GetPadTimesEvent(Int_t sector, Bool_t force)
467a4b25 850{
851 //
852 // return pointer to Pad Times Array for the current event and sector
853 // if force is true create it if it doesn't exist allready
854 //
855 TObjArray *arr = &fPadTimesArrayEvent;
856 return GetPadInfoEvent(sector,arr,force);
857}
858//_____________________________________________________________________
3cd27a08 859TVectorF* AliTPCCalibPulser::GetPadQEvent(Int_t sector, Bool_t force)
467a4b25 860{
861 //
862 // return pointer to Pad Q Array for the current event and sector
863 // if force is true create it if it doesn't exist allready
864 // for debugging purposes only
865 //
866
867 TObjArray *arr = &fPadQArrayEvent;
868 return GetPadInfoEvent(sector,arr,force);
869}
870//_____________________________________________________________________
3cd27a08 871TVectorF* AliTPCCalibPulser::GetPadRMSEvent(Int_t sector, Bool_t force)
467a4b25 872{
873 //
874 // return pointer to Pad RMS Array for the current event and sector
875 // if force is true create it if it doesn't exist allready
876 // for debugging purposes only
877 //
878 TObjArray *arr = &fPadRMSArrayEvent;
879 return GetPadInfoEvent(sector,arr,force);
880}
881//_____________________________________________________________________
3cd27a08 882TVectorF* AliTPCCalibPulser::GetPadPedestalEvent(Int_t sector, Bool_t force)
467a4b25 883{
884 //
885 // return pointer to Pad RMS Array for the current event and sector
886 // if force is true create it if it doesn't exist allready
887 // for debugging purposes only
888 //
889 TObjArray *arr = &fPadPedestalArrayEvent;
890 return GetPadInfoEvent(sector,arr,force);
891}
892//_____________________________________________________________________
3cd27a08 893AliTPCCalROC* AliTPCCalibPulser::GetCalRoc(Int_t sector, TObjArray* arr, Bool_t force) const
467a4b25 894{
895 //
896 // return pointer to ROC Calibration
897 // if force is true create a new histogram if it doesn't exist allready
898 //
899 if ( !force || arr->UncheckedAt(sector) )
900 return (AliTPCCalROC*)arr->UncheckedAt(sector);
901
902 // if we are forced and histogram doesn't yes exist create it
903
904 // new AliTPCCalROC for T0 information. One value for each pad!
905 AliTPCCalROC *croc = new AliTPCCalROC(sector);
906 arr->AddAt(croc,sector);
907 return croc;
908}
909//_____________________________________________________________________
3cd27a08 910AliTPCCalROC* AliTPCCalibPulser::GetCalRocT0(Int_t sector, Bool_t force)
467a4b25 911{
912 //
913 // return pointer to Carge ROC Calibration
914 // if force is true create a new histogram if it doesn't exist allready
915 //
916 TObjArray *arr = &fCalRocArrayT0;
917 return GetCalRoc(sector, arr, force);
918}
919//_____________________________________________________________________
3cd27a08 920AliTPCCalROC* AliTPCCalibPulser::GetCalRocQ(Int_t sector, Bool_t force)
467a4b25 921{
922 //
923 // return pointer to T0 ROC Calibration
924 // if force is true create a new histogram if it doesn't exist allready
925 //
926 TObjArray *arr = &fCalRocArrayQ;
927 return GetCalRoc(sector, arr, force);
928}
929//_____________________________________________________________________
3cd27a08 930AliTPCCalROC* AliTPCCalibPulser::GetCalRocRMS(Int_t sector, Bool_t force)
467a4b25 931{
932 //
933 // return pointer to signal width ROC Calibration
934 // if force is true create a new histogram if it doesn't exist allready
935 //
936 TObjArray *arr = &fCalRocArrayRMS;
937 return GetCalRoc(sector, arr, force);
938}
939//_____________________________________________________________________
940AliTPCCalROC* AliTPCCalibPulser::GetCalRocOutliers(Int_t sector, Bool_t force)
941{
942 //
943 // return pointer to Outliers
944 // if force is true create a new histogram if it doesn't exist allready
945 //
946 TObjArray *arr = &fCalRocArrayOutliers;
947 return GetCalRoc(sector, arr, force);
948}
949//_____________________________________________________________________
3cd27a08 950void AliTPCCalibPulser::ResetEvent()
467a4b25 951{
952 //
953 // Reset global counters -- Should be called before each event is processed
954 //
955 fLastSector=-1;
956 fCurrentSector=-1;
957 fCurrentRow=-1;
958 fCurrentChannel=-1;
959
960 ResetPad();
961
962 fPadTimesArrayEvent.Delete();
963 fPadQArrayEvent.Delete();
964 fPadRMSArrayEvent.Delete();
965 fPadPedestalArrayEvent.Delete();
966
967 for ( Int_t i=0; i<72; ++i ){
968 fVTime0Offset[i]=0;
969 fVTime0OffsetCounter[i]=0;
970 }
971}
972//_____________________________________________________________________
3cd27a08 973void AliTPCCalibPulser::ResetPad()
467a4b25 974{
975 //
976 // Reset pad infos -- Should be called after a pad has been processed
977 //
978 for (Int_t i=fFirstTimeBin; i<fLastTimeBin+1; ++i)
979 fPadSignal[i] = 0;
980 fMaxTimeBin = -1;
981 fMaxPadSignal = -1;
982 fPadPedestal = -1;
983 fPadNoise = -1;
984}
985//_____________________________________________________________________
986void AliTPCCalibPulser::Merge(AliTPCCalibPulser *sig)
987{
988 //
989 // Merge reference histograms of sig to the current AliTPCCalibPulser
990 //
991
992 //merge histograms
993 for (Int_t iSec=0; iSec<72; ++iSec){
994 TH2S *hRefQmerge = sig->GetHistoQ(iSec);
995 TH2S *hRefT0merge = sig->GetHistoT0(iSec);
996 TH2S *hRefRMSmerge = sig->GetHistoRMS(iSec);
997
998
999 if ( hRefQmerge ){
1000 TDirectory *dir = hRefQmerge->GetDirectory(); hRefQmerge->SetDirectory(0);
1001 TH2S *hRefQ = GetHistoQ(iSec);
1002 if ( hRefQ ) hRefQ->Add(hRefQmerge);
1003 else {
1004 TH2S *hist = new TH2S(*hRefQmerge);
1005 hist->SetDirectory(0);
1006 fHistoQArray.AddAt(hist, iSec);
1007 }
1008 hRefQmerge->SetDirectory(dir);
1009 }
1010 if ( hRefT0merge ){
1011 TDirectory *dir = hRefT0merge->GetDirectory(); hRefT0merge->SetDirectory(0);
1012 TH2S *hRefT0 = GetHistoT0(iSec);
1013 if ( hRefT0 ) hRefT0->Add(hRefT0merge);
1014 else {
1015 TH2S *hist = new TH2S(*hRefT0merge);
1016 hist->SetDirectory(0);
1017 fHistoT0Array.AddAt(hist, iSec);
1018 }
1019 hRefT0merge->SetDirectory(dir);
1020 }
1021 if ( hRefRMSmerge ){
1022 TDirectory *dir = hRefRMSmerge->GetDirectory(); hRefRMSmerge->SetDirectory(0);
1023 TH2S *hRefRMS = GetHistoRMS(iSec);
1024 if ( hRefRMS ) hRefRMS->Add(hRefRMSmerge);
1025 else {
1026 TH2S *hist = new TH2S(*hRefRMSmerge);
1027 hist->SetDirectory(0);
1028 fHistoRMSArray.AddAt(hist, iSec);
1029 }
1030 hRefRMSmerge->SetDirectory(dir);
1031 }
1032
1033 }
1034}
1035//_____________________________________________________________________
3cd27a08 1036void AliTPCCalibPulser::Analyse()
467a4b25 1037{
1038 //
1039 // Calculate calibration constants
1040 //
1041
1042 TVectorD paramQ(3);
1043 TVectorD paramT0(3);
1044 TVectorD paramRMS(3);
1045 TMatrixD dummy(3,3);
1046
1047 for (Int_t iSec=0; iSec<72; ++iSec){
1048 TH2S *hT0 = GetHistoT0(iSec);
1049 if (!hT0 ) continue;
1050
1051 AliTPCCalROC *rocQ = GetCalRocQ (iSec,kTRUE);
1052 AliTPCCalROC *rocT0 = GetCalRocT0 (iSec,kTRUE);
1053 AliTPCCalROC *rocRMS = GetCalRocRMS(iSec,kTRUE);
1054 AliTPCCalROC *rocOut = GetCalRocOutliers(iSec,kTRUE);
1055
1056 TH2S *hQ = GetHistoQ(iSec);
1057 TH2S *hRMS = GetHistoRMS(iSec);
1058
3cd27a08 1059 Short_t *arrayhQ = hQ->GetArray();
1060 Short_t *arrayhT0 = hT0->GetArray();
1061 Short_t *arrayhRMS = hRMS->GetArray();
467a4b25 1062
1063 UInt_t nChannels = fROC->GetNChannels(iSec);
1064
1065 //debug
1066 Int_t row=0;
1067 Int_t pad=0;
1068 Int_t padc=0;
1069 //! debug
1070
1071 for (UInt_t iChannel=0; iChannel<nChannels; ++iChannel){
1072
467a4b25 1073 Float_t cogTime0 = -1000;
1074 Float_t cogQ = -1000;
1075 Float_t cogRMS = -1000;
1076 Float_t cogOut = 0;
1077
467a4b25 1078 Int_t offsetQ = (fNbinsQ+2)*(iChannel+1)+1;
1079 Int_t offsetT0 = (fNbinsT0+2)*(iChannel+1)+1;
1080 Int_t offsetRMS = (fNbinsRMS+2)*(iChannel+1)+1;
467a4b25 1081/*
3cd27a08 1082 AliMathBase::FitGaus(arrayhQ+offsetQ,fNbinsQ,fXminQ,fXmaxQ,&paramQ,&dummy);
1083 AliMathBase::FitGaus(arrayhT0+offsetT0,fNbinsT0,fXminT0,fXmaxT0,&paramT0,&dummy);
1084 AliMathBase::FitGaus(arrayhRMS+offsetRMS,fNbinsRMS,fXminRMS,fXmaxRMS,&paramRMS,&dummy);
467a4b25 1085 cogQ = paramQ[1];
1086 cogTime0 = paramT0[1];
1087 cogRMS = paramRMS[1];
1088*/
3cd27a08 1089 cogQ = AliMathBase::GetCOG(arrayhQ+offsetQ,fNbinsQ,fXminQ,fXmaxQ);
1090 cogTime0 = AliMathBase::GetCOG(arrayhT0+offsetT0,fNbinsT0,fXminT0,fXmaxT0);
1091 cogRMS = AliMathBase::GetCOG(arrayhRMS+offsetRMS,fNbinsRMS,fXminRMS,fXmaxRMS);
467a4b25 1092
467a4b25 1093 /*
1094 if ( (cogQ < ??) && (cogTime0 > ??) && (cogTime0<??) && ( cogRMS>??) ){
1095 cogOut = 1;
1096 cogTime0 = 0;
1097 cogQ = 0;
1098 cogRMS = 0;
1099 }
1100*/
1101 rocQ->SetValue(iChannel, cogQ*cogQ);
1102 rocT0->SetValue(iChannel, cogTime0);
1103 rocRMS->SetValue(iChannel, cogRMS);
1104 rocOut->SetValue(iChannel, cogOut);
1105
1106
1107 //debug
1108 if ( fDebugLevel > 0 ){
1109 if ( !fDebugStreamer ) {
1110 //debug stream
1111 TDirectory *backup = gDirectory;
1112 fDebugStreamer = new TTreeSRedirector("deb2.root");
1113 if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer
1114 }
1115
1116 while ( iChannel > (fROC->GetRowIndexes(iSec)[row]+fROC->GetNPads(iSec,row)-1) ) row++;
1117 pad = iChannel-fROC->GetRowIndexes(iSec)[row];
1118 padc = pad-(fROC->GetNPads(iSec,row)/2);
1119
1120 (*fDebugStreamer) << "DataEnd" <<
1121 "Sector=" << iSec <<
1122 "Pad=" << pad <<
1123 "PadC=" << padc <<
1124 "Row=" << row <<
1125 "PadSec=" << iChannel <<
1126 "Q=" << cogQ <<
1127 "T0=" << cogTime0 <<
1128 "RMS=" << cogRMS <<
1129 "\n";
1130 }
1131 //! debug
467a4b25 1132 }
467a4b25 1133 }
1134 delete fDebugStreamer;
1135 fDebugStreamer = 0x0;
1136}
1137//_____________________________________________________________________
1138void AliTPCCalibPulser::DumpToFile(const Char_t *filename, const Char_t *dir, Bool_t append)
1139{
1140 //
1141 // Write class to file
1142 //
1143
1144 TString sDir(dir);
1145 TString option;
1146
1147 if ( append )
1148 option = "update";
1149 else
1150 option = "recreate";
1151
1152 TDirectory *backup = gDirectory;
1153 TFile f(filename,option.Data());
1154 f.cd();
1155 if ( !sDir.IsNull() ){
1156 f.mkdir(sDir.Data());
1157 f.cd(sDir);
1158 }
1159 this->Write();
1160 f.Close();
1161
1162 if ( backup ) backup->cd();
1163}
1164//_____________________________________________________________________
1165//_________________________ Test Functions ___________________________
1166//_____________________________________________________________________
1167TObjArray* AliTPCCalibPulser::TestBinning()
1168{
1169 //
1170 // Function to test the binning of the reference histograms
1171 // type: T0, Q or RMS
1172 // mode: 0 - number of filled bins per channel
1173 // 1 - number of empty bins between filled bins in one ROC
1174 // returns TObjArray with the test histograms type*2+mode:
1175 // position 0 = T0,0 ; 1 = T0,1 ; 2 = Q,0 ...
1176
1177
1178 TObjArray *histArray = new TObjArray(6);
1179 const Char_t *type[] = {"T0","Q","RMS"};
1180 Int_t fNbins[3] = {fNbinsT0,fNbinsQ,fNbinsRMS};
1181
1182 for (Int_t itype = 0; itype<3; ++itype){
1183 for (Int_t imode=0; imode<2; ++imode){
1184 Int_t icount = itype*2+imode;
1185 histArray->AddAt(new TH1F(Form("hTestBinning%s%d",type[itype],imode),
1186 Form("Test Binning of '%s', mode - %d",type[itype],imode),
1187 72,0,72),
1188 icount);
1189 }
1190 }
1191
1192
1193 TH2S *hRef=0x0;
1194 Short_t *array=0x0;
1195 for (Int_t itype = 0; itype<3; ++itype){
1196 for (Int_t iSec=0; iSec<72; ++iSec){
1197 if ( itype == 0 ) hRef = GetHistoT0(iSec);
1198 if ( itype == 1 ) hRef = GetHistoQ(iSec);
1199 if ( itype == 2 ) hRef = GetHistoRMS(iSec);
1200 if ( hRef == 0x0 ) continue;
1201 array = (hRef->GetArray());
1202 UInt_t nChannels = fROC->GetNChannels(iSec);
1203
1204 Int_t nempty=0;
1205 for (UInt_t iChannel=0; iChannel<nChannels; ++iChannel){
1206 Int_t nfilled=0;
1207 Int_t offset = (fNbins[itype]+2)*(iChannel+1)+1;
1208 Int_t c1 = 0;
1209 Int_t c2 = 0;
1210 for (Int_t iBin=0; iBin<fNbins[itype]; ++iBin){
1211 if ( array[offset+iBin]>0 ) {
1212 nfilled++;
1213 if ( c1 && c2 ) nempty++;
1214 else c1 = 1;
1215 }
1216 else if ( c1 ) c2 = 1;
1217
1218 }
1219 ((TH1F*)histArray->At(itype*2))->Fill(nfilled);
1220 }
1221 ((TH1F*)histArray->At(itype*2+1))->Fill(iSec,nempty);
1222 }
1223 }
1224 return histArray;
1225}