use the new AliCaloRawStreamV3 decoder (from Yuri), based on the new AliAltroRawStrea...
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALRawUtils.cxx
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
18 //_________________________________________________________________________
19 //  Utility Class for handling Raw data
20 //  Does all transitions from Digits to Raw and vice versa, 
21 //  for simu and reconstruction
22 //
23 //  Note: the current version is still simplified. Only 
24 //    one raw signal per digit is generated; either high-gain or low-gain
25 //    Need to add concurrent high and low-gain info in the future
26 //    No pedestal is added to the raw signal.
27 //*-- Author: Marco van Leeuwen (LBL)
28
29 #include "AliEMCALRawUtils.h"
30   
31 #include "TF1.h"
32 #include "TGraph.h"
33 class TSystem;
34   
35 class AliLog;
36 #include "AliRun.h"
37 #include "AliRunLoader.h"
38 class AliCaloAltroMapping;
39 #include "AliAltroBuffer.h"
40 #include "AliRawReader.h"
41 #include "AliCaloRawStreamV3.h"
42 #include "AliDAQ.h"
43   
44 #include "AliEMCALRecParam.h"
45 #include "AliEMCALLoader.h"
46 #include "AliEMCALGeometry.h"
47 class AliEMCALDigitizer;
48 #include "AliEMCALDigit.h"
49 #include "AliEMCAL.h"
50   
51 ClassImp(AliEMCALRawUtils)
52   
53 // Signal shape parameters
54 Double_t AliEMCALRawUtils::fgTimeBinWidth  = 100E-9 ; // each sample is 100 ns
55 Double_t AliEMCALRawUtils::fgTimeTrigger = 1.5E-6 ;   // 15 time bins ~ 1.5 musec
56
57 // some digitization constants
58 Int_t    AliEMCALRawUtils::fgThreshold = 1;
59 Int_t    AliEMCALRawUtils::fgDDLPerSuperModule = 2;  // 2 ddls per SuperModule
60 Int_t    AliEMCALRawUtils::fgPedestalValue = 32;     // pedestal value for digits2raw
61 Double_t AliEMCALRawUtils::fgFEENoise = 3.;          // 3 ADC channels of noise (sampled)
62
63 AliEMCALRawUtils::AliEMCALRawUtils()
64   : fHighLowGainFactor(0.), fOrder(0), fTau(0.), fNoiseThreshold(0),
65     fNPedSamples(0), fGeom(0), fOption("")
66 {
67
68   //These are default parameters.  
69   //Can be re-set from without with setter functions
70   fHighLowGainFactor = 16. ;          // adjusted for a low gain range of 82 GeV (10 bits) 
71   fOrder = 2;                         // order of gamma fn
72   fTau = 2.35;                        // in units of timebin, from CERN 2007 testbeam
73   fNoiseThreshold = 3;
74   fNPedSamples = 5;
75
76   //Get Mapping RCU files from the AliEMCALRecParam                                 
77   const TObjArray* maps = AliEMCALRecParam::GetMappings();
78   if(!maps) AliFatal("Cannot retrieve ALTRO mappings!!");
79
80   for(Int_t i = 0; i < 4; i++) {
81     fMapping[i] = (AliAltroMapping*)maps->At(i);
82   }
83
84   //To make sure we match with the geometry in a simulation file,
85   //let's try to get it first.  If not, take the default geometry
86   AliRunLoader *rl = AliRunLoader::Instance();
87   if(!rl) AliError("Cannot find RunLoader!");
88   if (rl->GetAliRun() && rl->GetAliRun()->GetDetector("EMCAL")) {
89     fGeom = dynamic_cast<AliEMCAL*>(rl->GetAliRun()->GetDetector("EMCAL"))->GetGeometry();
90   } else {
91     AliInfo(Form("Using default geometry in raw reco"));
92     fGeom =  AliEMCALGeometry::GetInstance(AliEMCALGeometry::GetDefaultGeometryName());
93   }
94
95   if(!fGeom) AliFatal(Form("Could not get geometry!"));
96
97 }
98
99 //____________________________________________________________________________
100 AliEMCALRawUtils::AliEMCALRawUtils(AliEMCALGeometry *pGeometry)
101   : fHighLowGainFactor(0.), fOrder(0), fTau(0.), fNoiseThreshold(0),
102     fNPedSamples(0), fGeom(pGeometry), fOption("")
103 {
104   //
105   // Initialize with the given geometry - constructor required by HLT
106   // HLT does not use/support AliRunLoader(s) instances
107   // This is a minimum intervention solution
108   // Comment by MPloskon@lbl.gov
109   //
110
111   //These are default parameters. 
112   //Can be re-set from without with setter functions 
113   fHighLowGainFactor = 16. ;          // adjusted for a low gain range of 82 GeV (10 bits)
114   fOrder = 2;                         // order of gamma fn
115   fTau = 2.35;                        // in units of timebin, from CERN 2007 testbeam
116   fNoiseThreshold = 3;
117   fNPedSamples = 5;
118
119   //Get Mapping RCU files from the AliEMCALRecParam
120   const TObjArray* maps = AliEMCALRecParam::GetMappings();
121   if(!maps) AliFatal("Cannot retrieve ALTRO mappings!!");
122
123   for(Int_t i = 0; i < 4; i++) {
124     fMapping[i] = (AliAltroMapping*)maps->At(i);
125   }
126
127   if(!fGeom) AliFatal(Form("Could not get geometry!"));
128
129 }
130
131 //____________________________________________________________________________
132 AliEMCALRawUtils::AliEMCALRawUtils(const AliEMCALRawUtils& rawU)
133   : TObject(),
134     fHighLowGainFactor(rawU.fHighLowGainFactor), 
135     fOrder(rawU.fOrder),
136     fTau(rawU.fTau),
137     fNoiseThreshold(rawU.fNoiseThreshold),
138     fNPedSamples(rawU.fNPedSamples),
139     fGeom(rawU.fGeom), 
140     fOption(rawU.fOption)
141 {
142   //copy ctor
143   fMapping[0] = rawU.fMapping[0];
144   fMapping[1] = rawU.fMapping[1];
145   fMapping[2] = rawU.fMapping[2];
146   fMapping[3] = rawU.fMapping[3];
147 }
148
149 //____________________________________________________________________________
150 AliEMCALRawUtils& AliEMCALRawUtils::operator =(const AliEMCALRawUtils &rawU)
151 {
152   //assignment operator
153
154   if(this != &rawU) {
155     fHighLowGainFactor = rawU.fHighLowGainFactor;
156     fOrder = rawU.fOrder;
157     fTau = rawU.fTau;
158     fNoiseThreshold = rawU.fNoiseThreshold;
159     fNPedSamples = rawU.fNPedSamples;
160     fGeom = rawU.fGeom;
161     fOption = rawU.fOption;
162     fMapping[0] = rawU.fMapping[0];
163     fMapping[1] = rawU.fMapping[1];
164     fMapping[2] = rawU.fMapping[2];
165     fMapping[3] = rawU.fMapping[3];
166   }
167
168   return *this;
169
170 }
171
172 //____________________________________________________________________________
173 AliEMCALRawUtils::~AliEMCALRawUtils() {
174   //dtor
175
176 }
177
178 //____________________________________________________________________________
179 void AliEMCALRawUtils::Digits2Raw()
180 {
181   // convert digits of the current event to raw data
182   
183   AliRunLoader *rl = AliRunLoader::Instance();
184   AliEMCALLoader *loader = dynamic_cast<AliEMCALLoader*>(rl->GetDetectorLoader("EMCAL"));
185
186   // get the digits
187   loader->LoadDigits("EMCAL");
188   loader->GetEvent();
189   TClonesArray* digits = loader->Digits() ;
190   
191   if (!digits) {
192     Warning("Digits2Raw", "no digits found !");
193     return;
194   }
195
196   static const Int_t nDDL = 12*2; // 12 SM hardcoded for now. Buffers allocated dynamically, when needed, so just need an upper limit here
197   AliAltroBuffer* buffers[nDDL];
198   for (Int_t i=0; i < nDDL; i++)
199     buffers[i] = 0;
200
201   Int_t adcValuesLow[fgkTimeBins];
202   Int_t adcValuesHigh[fgkTimeBins];
203
204   // loop over digits (assume ordered digits)
205   for (Int_t iDigit = 0; iDigit < digits->GetEntries(); iDigit++) {
206     AliEMCALDigit* digit = dynamic_cast<AliEMCALDigit *>(digits->At(iDigit)) ;
207     if (digit->GetAmp() < fgThreshold) 
208       continue;
209
210     //get cell indices
211     Int_t nSM = 0;
212     Int_t nIphi = 0;
213     Int_t nIeta = 0;
214     Int_t iphi = 0;
215     Int_t ieta = 0;
216     Int_t nModule = 0;
217     fGeom->GetCellIndex(digit->GetId(), nSM, nModule, nIphi, nIeta);
218     fGeom->GetCellPhiEtaIndexInSModule(nSM, nModule, nIphi, nIeta,iphi, ieta) ;
219     
220     //Check which is the RCU, 0 or 1, of the cell.
221     Int_t iRCU = -111;
222     //RCU0
223     if (0<=iphi&&iphi<8) iRCU=0; // first cable row
224     else if (8<=iphi&&iphi<16 && 0<=ieta&&ieta<24) iRCU=0; // first half; 
225     //second cable row
226     //RCU1
227     else if(8<=iphi&&iphi<16 && 24<=ieta&&ieta<48) iRCU=1; // second half; 
228     //second cable row
229     else if(16<=iphi&&iphi<24) iRCU=1; // third cable row
230
231     if (nSM%2==1) iRCU = 1 - iRCU; // swap for odd=C side, to allow us to cable both sides the same
232
233     if (iRCU<0) 
234       Fatal("Digits2Raw()","Non-existent RCU number: %d", iRCU);
235     
236     //Which DDL?
237     Int_t iDDL = fgDDLPerSuperModule* nSM + iRCU;
238     if (iDDL >= nDDL)
239       Fatal("Digits2Raw()","Non-existent DDL board number: %d", iDDL);
240
241     if (buffers[iDDL] == 0) {      
242       // open new file and write dummy header
243       TString fileName = AliDAQ::DdlFileName("EMCAL",iDDL);
244       //Select mapping file RCU0A, RCU0C, RCU1A, RCU1C
245       Int_t iRCUside=iRCU+(nSM%2)*2;
246       //iRCU=0 and even (0) SM -> RCU0A.data   0
247       //iRCU=1 and even (0) SM -> RCU1A.data   1
248       //iRCU=0 and odd  (1) SM -> RCU0C.data   2
249       //iRCU=1 and odd  (1) SM -> RCU1C.data   3
250       //cout<<" nSM "<<nSM<<"; iRCU "<<iRCU<<"; iRCUside "<<iRCUside<<endl;
251       buffers[iDDL] = new AliAltroBuffer(fileName.Data(),fMapping[iRCUside]);
252       buffers[iDDL]->WriteDataHeader(kTRUE, kFALSE);  //Dummy;
253     }
254     
255     // out of time range signal (?)
256     if (digit->GetTimeR() > GetRawFormatTimeMax() ) {
257       AliInfo("Signal is out of time range.\n");
258       buffers[iDDL]->FillBuffer((Int_t)digit->GetAmp());
259       buffers[iDDL]->FillBuffer(GetRawFormatTimeBins() );  // time bin
260       buffers[iDDL]->FillBuffer(3);          // bunch length      
261       buffers[iDDL]->WriteTrailer(3, ieta, iphi, nSM);  // trailer
262       // calculate the time response function
263     } else {
264       Bool_t lowgain = RawSampledResponse(digit->GetTimeR(), digit->GetAmp(), adcValuesHigh, adcValuesLow) ; 
265       if (lowgain) 
266         buffers[iDDL]->WriteChannel(ieta, iphi, 0, GetRawFormatTimeBins(), adcValuesLow, fgThreshold);
267       else 
268         buffers[iDDL]->WriteChannel(ieta,iphi, 1, GetRawFormatTimeBins(), adcValuesHigh, fgThreshold);
269     }
270   }
271   
272   // write headers and close files
273   for (Int_t i=0; i < nDDL; i++) {
274     if (buffers[i]) {
275       buffers[i]->Flush();
276       buffers[i]->WriteDataHeader(kFALSE, kFALSE);
277       delete buffers[i];
278     }
279   }
280
281   loader->UnloadDigits();
282 }
283
284 //____________________________________________________________________________
285 void AliEMCALRawUtils::Raw2Digits(AliRawReader* reader,TClonesArray *digitsArr)
286 {
287   // convert raw data of the current event to digits                                                                                     
288
289   digitsArr->Clear(); 
290
291   if (!digitsArr) {
292     Error("Raw2Digits", "no digits found !");
293     return;
294   }
295   if (!reader) {
296     Error("Raw2Digits", "no raw reader found !");
297     return;
298   }
299
300   AliCaloRawStreamV3 in(reader,"EMCAL",fMapping);
301   // Select EMCAL DDL's;
302   reader->Select("EMCAL");
303
304   //Updated fitting routine from 2007 beam test takes into account
305   //possibility of two peaks in data and selects first one for fitting
306   //Also sets some of the starting parameters based on the shape of the
307   //given raw signal being fit
308
309   TF1 * signalF = new TF1("signal", RawResponseFunction, 0, GetRawFormatTimeBins(), 5);
310   signalF->SetParameters(10.,0.,fTau,fOrder,5.); //set all defaults once, just to be safe
311   signalF->SetParNames("amp","t0","tau","N","ped");
312   signalF->SetParameter(2,fTau); // tau in units of time bin
313   signalF->SetParLimits(2,2,-1);
314   signalF->SetParameter(3,fOrder); // order
315   signalF->SetParLimits(3,2,-1);
316   
317   Int_t id =  -1;
318   Float_t time = 0. ; 
319   Float_t amp = 0. ; 
320   Int_t i = 0;
321   Int_t startBin = 0;
322
323   //Graph to hold data we will fit (should be converted to an array
324   //later to speed up processing
325   TGraph * gSig = new TGraph(GetRawFormatTimeBins()); 
326
327   Int_t lowGain = 0;
328   Int_t caloFlag = 0; // low, high gain, or TRU, or LED ref.
329
330   // start loop over input stream 
331   while (in.NextDDL()) {
332     while (in.NextChannel()) {
333       
334       id =  fGeom->GetAbsCellIdFromCellIndexes(in.GetModule(), in.GetRow(), in.GetColumn()) ;
335       caloFlag = in.GetCaloFlag();
336       lowGain = in.IsLowGain();
337
338       // There can be zero-suppression in the raw data, 
339       // so set up the TGraph in advance
340       for (i=0; i < GetRawFormatTimeBins(); i++) {
341         gSig->SetPoint(i, i , 0);
342       }
343       Int_t maxTime = 0;
344
345       while (in.NextBunch()) {
346         const UShort_t *sig = in.GetSignals();
347         startBin = in.GetStartTimeBin();
348
349         if (maxTime < in.GetStartTimeBin()) {
350           maxTime = in.GetStartTimeBin(); // timebins come in reverse order
351         }
352
353         if (maxTime < 0 || maxTime >= GetRawFormatTimeBins()) {
354           AliWarning(Form("Invalid time bin %d",maxTime));
355           maxTime = GetRawFormatTimeBins();
356         }
357
358         for (i = 0; i < in.GetBunchLength(); i++) {
359           time = startBin--;
360           gSig->SetPoint(time, time, sig[i]) ;
361         }
362       } // loop over bunches
363     
364       gSig->Set(maxTime+1);
365       FitRaw(gSig, signalF, amp, time) ; 
366     
367       if (caloFlag == 0 || caloFlag == 1) { // low gain or high gain 
368         if (amp > 0 && amp < 2000) {  //check both high and low end of
369         //result, 2000 is somewhat arbitrary - not nice with magic numbers in the code..
370           AliDebug(2,Form("id %d lowGain %d amp %g", id, lowGain, amp));
371         
372           AddDigit(digitsArr, id, lowGain, (Int_t)amp, time);
373         }
374         
375       }
376
377       // Reset graph
378       for (Int_t index = 0; index < gSig->GetN(); index++) {
379         gSig->SetPoint(index, index, 0) ;  
380       } 
381       // Reset starting parameters for fit function
382       signalF->SetParameters(10.,0.,fTau,fOrder,5.); //reset all defaults just to be safe
383
384    } // end while over channel   
385   } //end while over DDL's, of input stream 
386   
387   delete signalF ; 
388   delete gSig;
389   
390   return ; 
391 }
392
393 //____________________________________________________________________________ 
394 void AliEMCALRawUtils::AddDigit(TClonesArray *digitsArr, Int_t id, Int_t lowGain, Int_t amp, Float_t time) {
395   //
396   // Add a new digit. 
397   // This routine checks whether a digit exists already for this tower 
398   // and then decides whether to use the high or low gain info
399   //
400   // Called by Raw2Digits
401   
402   AliEMCALDigit *digit = 0, *tmpdigit = 0;
403   
404   TIter nextdigit(digitsArr);
405   while (digit == 0 && (tmpdigit = (AliEMCALDigit*) nextdigit())) {
406     if (tmpdigit->GetId() == id)
407       digit = tmpdigit;
408   }
409
410   if (!digit) { // no digit existed for this tower; create one
411     if (lowGain) 
412       amp = Int_t(fHighLowGainFactor * amp); 
413     Int_t idigit = digitsArr->GetEntries();
414     new((*digitsArr)[idigit]) AliEMCALDigit( -1, -1, id, amp, time, idigit) ;   
415   }
416   else { // a digit already exists, check range 
417          // (use high gain if signal < cut value, otherwise low gain)
418     if (lowGain) { // new digit is low gain
419       if (digit->GetAmp() > fgkOverflowCut) {  // use if stored digit is out of range
420         digit->SetAmp(Int_t(fHighLowGainFactor * amp));
421         digit->SetTime(time);
422       }
423     }
424     else if (amp < fgkOverflowCut) { // new digit is high gain; use if not out of range
425       digit->SetAmp(amp);
426       digit->SetTime(time);
427     }
428   }
429 }
430
431 //____________________________________________________________________________ 
432 void AliEMCALRawUtils::FitRaw(TGraph * gSig, TF1* signalF, Float_t & amp, Float_t & time) const 
433 {
434   // Fits the raw signal time distribution; from AliEMCALGetter 
435
436   amp = time = 0. ; 
437   Double_t ped = 0;
438   Int_t nPed = 0;
439
440   for (Int_t index = 0; index < fNPedSamples; index++) {
441     Double_t ttime, signal;
442     gSig->GetPoint(index, ttime, signal) ; 
443     if (signal > 0) {
444       ped += signal;
445       nPed++;
446     }
447   }
448
449   if (nPed > 0)
450     ped /= nPed;
451   else {
452     AliWarning("Could not determine pedestal");   
453     ped = 10; // put some small value as first guess
454   }
455
456   Int_t maxFound = 0;
457   Int_t iMax = 0;
458   Float_t max = -1;
459   Float_t maxFit = gSig->GetN();
460   Float_t minAfterSig = 9999;
461   Int_t tminAfterSig = gSig->GetN();
462   Int_t nPedAfterSig = 0;
463   Int_t plateauWidth = 0;
464   Int_t plateauStart = 9999;
465   Float_t cut = 0.3;
466
467   for (Int_t i=fNPedSamples; i < gSig->GetN(); i++) {
468     Double_t ttime, signal;
469     gSig->GetPoint(i, ttime, signal) ; 
470     if (!maxFound && signal > max) {
471       iMax = i;
472       max = signal;
473     }
474     else if ( max > ped + fNoiseThreshold ) {
475       maxFound = 1;
476       minAfterSig = signal;
477       tminAfterSig = i;
478     }
479     if (maxFound) {
480       if ( signal < minAfterSig) {
481         minAfterSig = signal;
482         tminAfterSig = i;
483       }
484       if (i > tminAfterSig + 5) {  // Two close peaks; end fit at minimum
485         maxFit = tminAfterSig;
486         break;
487       }
488       if ( signal < cut*max){   //stop fit at 30% amplitude(avoid the pulse shape falling edge)
489         maxFit = i;
490         break;
491       }
492       if ( signal < ped + fNoiseThreshold)
493         nPedAfterSig++;
494       if (nPedAfterSig >= 5) {  // include 5 pedestal bins after peak
495         maxFit = i;
496         break;
497       }
498     }
499     //Add check on plateau
500     if (signal >= fgkRawSignalOverflow - fNoiseThreshold) {
501       if(plateauWidth == 0) plateauStart = i;
502       plateauWidth++;
503     }
504   }
505
506   if(plateauWidth > 0) {
507     for(int j = 0; j < plateauWidth; j++) {
508       //Note, have to remove the same point N times because after each
509       //remove, the positions of all subsequent points have shifted down
510       gSig->RemovePoint(plateauStart);
511     }
512   }
513
514   if ( max - ped > fNoiseThreshold ) { // else its noise 
515     AliDebug(2,Form("Fitting max %d ped %d", max, ped));
516     signalF->SetRange(0,maxFit);
517
518     if(max-ped > 50) 
519       signalF->SetParLimits(2,1,3);
520
521     signalF->SetParameter(4, ped) ; 
522     signalF->SetParameter(1, iMax);
523     signalF->SetParameter(0, max);
524     
525     gSig->Fit(signalF, "QROW"); // Note option 'W': equal errors on all points
526     amp = signalF->GetParameter(0); 
527     time = signalF->GetParameter(1)*GetRawFormatTimeBinWidth() - fgTimeTrigger;
528   }
529   return;
530 }
531 //__________________________________________________________________
532 Double_t AliEMCALRawUtils::RawResponseFunction(Double_t *x, Double_t *par)
533 {
534   // Matches version used in 2007 beam test
535   //
536   // Shape of the electronics raw reponse:
537   // It is a semi-gaussian, 2nd order Gamma function of the general form
538   //
539   // t' = (t - t0 + tau) / tau
540   // F = A * t**N * exp( N * ( 1 - t) )   for t >= 0
541   // F = 0                                for t < 0 
542   //
543   // parameters:
544   // A:   par[0]   // Amplitude = peak value
545   // t0:  par[1]
546   // tau: par[2]
547   // N:   par[3]
548   // ped: par[4]
549   //
550   Double_t signal ;
551   Double_t tau =par[2];
552   Double_t n =par[3];
553   Double_t ped = par[4];
554   Double_t xx = ( x[0] - par[1] + tau ) / tau ;
555
556   if (xx <= 0) 
557     signal = ped ;  
558   else {  
559     signal = ped + par[0] * TMath::Power(xx , n) * TMath::Exp(n * (1 - xx )) ; 
560   }
561   return signal ;  
562 }
563
564 //__________________________________________________________________
565 Bool_t AliEMCALRawUtils::RawSampledResponse(
566 const Double_t dtime, const Double_t damp, Int_t * adcH, Int_t * adcL) const 
567 {
568   // for a start time dtime and an amplitude damp given by digit, 
569   // calculates the raw sampled response AliEMCAL::RawResponseFunction
570
571   Bool_t lowGain = kFALSE ; 
572
573   // A:   par[0]   // Amplitude = peak value
574   // t0:  par[1]                            
575   // tau: par[2]                            
576   // N:   par[3]                            
577   // ped: par[4]
578
579   TF1 signalF("signal", RawResponseFunction, 0, GetRawFormatTimeBins(), 5);
580   signalF.SetParameter(0, damp) ; 
581   signalF.SetParameter(1, (dtime + fgTimeTrigger)/fgTimeBinWidth) ; 
582   signalF.SetParameter(2, fTau) ; 
583   signalF.SetParameter(3, fOrder);
584   signalF.SetParameter(4, fgPedestalValue);
585
586   for (Int_t iTime = 0; iTime < GetRawFormatTimeBins(); iTime++) {
587     Double_t signal = signalF.Eval(iTime) ;     
588
589     //According to Terry Awes, 13-Apr-2008
590     //add gaussian noise in quadrature to each sample
591     //Double_t noise = gRandom->Gaus(0.,fgFEENoise);
592     //signal = sqrt(signal*signal + noise*noise);
593
594     adcH[iTime] =  static_cast<Int_t>(signal + 0.5) ;
595     if ( adcH[iTime] > fgkRawSignalOverflow ){  // larger than 10 bits 
596       adcH[iTime] = fgkRawSignalOverflow ;
597       lowGain = kTRUE ; 
598     }
599
600     signal /= fHighLowGainFactor;
601
602     adcL[iTime] =  static_cast<Int_t>(signal + 0.5) ;
603     if ( adcL[iTime] > fgkRawSignalOverflow)  // larger than 10 bits 
604       adcL[iTime] = fgkRawSignalOverflow ;
605   }
606   return lowGain ; 
607 }