1 /**************************************************************************
2 * Copyright(c) 2007, 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 **************************************************************************/
18 // This class decodes the stream of ALTRO samples to extract
19 // the PHOS "digits" of current event. Uses fitting procedure
20 // to separate reasonable samples
23 // AliRawReader* rf = new AliRawReaderDate("2006run2211.raw");
24 // AliPHOSRawDecoder dc(rf);
25 // while (rf->NextEvent()) {
26 // dc.SubtractPedestals(kTRUE);
27 // while ( dc.NextDigit() ) {
28 // Int_t module = dc.GetModule();
29 // Int_t column = dc.GetColumn();
30 // Int_t row = dc.GetRow();
31 // Double_t amplitude = dc.GetEnergy();
32 // Double_t time = dc.GetTime();
33 // Bool_t IsLowGain = dc.IsLowGain();
38 // Author: Dmitri Peressounko using idea of Y.Kucheryaev
40 // --- ROOT system ---
51 // --- AliRoot header files ---
53 #include "AliPHOSRawDecoderv2.h"
54 #include "AliPHOSPulseGenerator.h"
57 ClassImp(AliPHOSRawDecoderv2)
59 //-----------------------------------------------------------------------------
60 AliPHOSRawDecoderv2::AliPHOSRawDecoderv2():AliPHOSRawDecoder(),
61 fNtimeSamples(25),fRMScut(11.)
63 //Default constructor.
72 //-----------------------------------------------------------------------------
73 AliPHOSRawDecoderv2::AliPHOSRawDecoderv2(AliRawReader* rawReader, AliAltroMapping **mapping):
74 AliPHOSRawDecoder(rawReader,mapping),
75 fNtimeSamples(25),fRMScut(11.)
77 //Construct a decoder object.
78 //Is is user responsibility to provide next raw event
79 //using AliRawReader::NextEvent().
88 //-----------------------------------------------------------------------------
89 AliPHOSRawDecoderv2::~AliPHOSRawDecoderv2()
95 //-----------------------------------------------------------------------------
96 AliPHOSRawDecoderv2::AliPHOSRawDecoderv2(const AliPHOSRawDecoderv2 &phosDecoder ):
97 AliPHOSRawDecoder(phosDecoder),
98 fNtimeSamples(25),fRMScut(11.)
101 fNtimeSamples=phosDecoder.fNtimeSamples ;
102 for(Int_t i=0; i<3;i++){
103 fLGpar[i]=phosDecoder.fLGpar[i] ;
104 fHGpar[i]=phosDecoder.fHGpar[i] ;
106 fRMScut=phosDecoder.fRMScut ;
109 //-----------------------------------------------------------------------------
110 AliPHOSRawDecoderv2& AliPHOSRawDecoderv2::operator = (const AliPHOSRawDecoderv2 &phosDecoder)
112 //Assignment operator.
114 fNtimeSamples=phosDecoder.fNtimeSamples ;
115 for(Int_t i=0; i<3;i++){
116 fLGpar[i]=phosDecoder.fLGpar[i] ;
117 fHGpar[i]=phosDecoder.fHGpar[i] ;
119 fRMScut=phosDecoder.fRMScut ;
123 //-----------------------------------------------------------------------------
124 Bool_t AliPHOSRawDecoderv2::NextDigit()
126 //Extract an energy deposited in the crystal,
127 //crystal' position (module,column,row),
128 //time and gain (high or low).
129 //First collects sample, then evaluates it and if it has
130 //reasonable shape, fits it with Gamma2 function and extracts
133 TCanvas * cs = (TCanvas*)gROOT->FindObjectAny("CSample") ;
135 cs = new TCanvas("CSample","CSample") ;
137 TH1D * h = (TH1D*)gROOT->FindObjectAny("hSample") ;
139 h=new TH1D("hSample","",200,0.,200.) ;
142 AliCaloRawStream* in = fCaloStream;
144 Int_t iBin = fSamples->GetSize();
146 Double_t pedMean = 0;
149 Double_t baseLine = 1.0;
150 const Int_t nPreSamples = 10;
153 while ( in->Next() ) {
155 // Fit the full sample
156 if(in->IsNewHWAddress() && iBin!=fSamples->GetSize()) {
158 Double_t pedestal =0. ;
161 pedestal = (Double_t)(pedMean/nPed);
165 for(Int_t i=0; i<fSamples->GetSize(); i++){
166 h->SetBinContent(i+1,fSamples->At(i)) ;
169 //calculate time and energy
172 for(Int_t i=iBin; i<fSamples->GetSize(); i++){
173 if(maxAmp<fSamples->At(i)){
175 maxAmp=fSamples->At(i) ;
178 if(maxBin==fSamples->GetSize()-1){//bad sample
183 fEnergy=Double_t(maxAmp)-pedestal ;
184 fOverflow =0 ; //look for plato on the top of sample
185 if(fEnergy>500 && //this is not fluctuation of soft sample
186 maxBin<fSamples->GetSize()-1 && fSamples->At(maxBin+1)==maxAmp){ //and there is a plato
191 // if(fRow==54 && fColumn==24){
192 // printf("fE=%f, ped=%f, row=%d, col=%d \n",fEnergy,pedestal,fRow,fColumn) ;
193 // if(fOverflow)printf(" Overflow \n") ;
194 // else printf("iBin=%d, maxBin=%d, maxAmp=%d,Amp(+1)=%d,Amp(-1)=%d \n",iBin,maxBin,maxAmp,fSamples->At(maxBin+1),fSamples->At(maxBin-1)) ;
202 return kTRUE ; //do not calculate energy and time for overflowed channels
204 if(fEnergy<baseLine){ //do not evaluate time, drop this sample
210 //else calculate time
227 for(Int_t i=iBin+1; i<fSamples->GetSize()&& cnts<fNtimeSamples; i++){
228 if(fSamples->At(i)<pedestal)
230 //Presently we do not use channels with overflow
231 // if(fOverflow && (fSamples->At(i)==maxAmp ||fSamples->At(i-1)==maxAmp)) //can not calculate time in overflow bin
233 if(fTimes->At(i)-fTimes->At(i-1)!=1) //do not use samples with non-consequtive points
235 Double_t de=fSamples->At(i)-pedestal ;
236 Double_t av = de+fSamples->At(i-1)-pedestal ;
237 if(av<=0.) //this is fluctuation around pedestal, scip
239 Double_t ds = fSamples->At(i)-fSamples->At(i-1) ;
240 Double_t ti = ds/av ; // calculate log. derivative
241 ti=a/(ti+b)-c*ti ; // and compare with parameterization
242 ti=Double_t(fTimes->At(i))-ti ;
243 Double_t wi = TMath::Abs(ds) ;
252 fQuality = tRMS/tW-fTime*fTime ;
254 //printf("t0=%f, RMS=%f, cut=%f \n",fTime,tRMS,fRMScut) ;
255 // if(tRMS>=fRMScut){ //bad sample
266 for(Int_t i=iBin+1; i<fSamples->GetSize()-1&&!isBad; i++){
267 if(fSamples->At(i)>fSamples->At(i-1)+5 && fSamples->At(i)>fSamples->At(i+1)+5) { //single jump
274 pedRMS=pedRMS/nPed-pedestal*pedestal ;
278 for(Int_t i=iBin+1; i<fSamples->GetSize()-1&&!isBad; i++){
279 if(fSamples->At(i)<pedestal-1)
286 if(fEnergy>10. && !isBad ){
287 printf("fE=%f, ped=%f, fQuality=%f, pedRMS=%f \n",fEnergy,pedestal,fQuality,pedRMS) ;
288 if(fOverflow)printf(" Overflow \n") ;
289 if(isBad)printf("bad") ;
290 // else printf("iBin=%d, maxBin=%d, maxAmp=%d,Amp(+1)=%d,Amp(-1)=%d \n",iBin,maxBin,maxAmp,fSamples->At(maxBin+1),fSamples->At(maxBin-1)) ;
298 return kTRUE ; //will scan further
301 fLowGainFlag = in->IsLowGain();
302 fModule = in->GetModule()+1;
303 fRow = in->GetRow() +1;
304 fColumn = in->GetColumn()+1;
307 // Fill array with samples
309 if(fPedSubtract && (in->GetTime() < nPreSamples)) {
310 pedMean += in->GetSignal();
311 pedRMS += in->GetSignal()*in->GetSignal();
314 fSamples->AddAt(in->GetSignal(),iBin);
315 fTimes->AddAt(in->GetTime(),iBin);