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 plots samples and qualify their quality according to their shape
21 // AliPHOSRawFitter *fitter=new AliPHOSRawFitter();
22 // fitter->SetChannelGeo(module,cellX,cellZ,caloFlag);
23 // fitter->SetCalibData(fgCalibData) ;
24 // fitter->Eval(sig,sigStart,sigLength);
25 // Double_t amplitude = fitter.GetEnergy();
26 // Double_t time = fitter.GetTime();
27 // Bool_t isLowGain = fitter.GetCaloFlag()==0;
29 // Author: Dmitri Peressounko (Oct.2008)
30 // Modified: Yuri Kharlov (Jul.2009)
32 // --- ROOT system ---
42 // --- AliRoot header files ---
44 #include "AliPHOSCalibData.h"
45 #include "AliPHOSRawFitterv2.h"
46 #include "AliPHOSPulseGenerator.h"
48 ClassImp(AliPHOSRawFitterv2)
50 //-----------------------------------------------------------------------------
51 AliPHOSRawFitterv2::AliPHOSRawFitterv2():
56 //Default constructor.
65 //-----------------------------------------------------------------------------
66 AliPHOSRawFitterv2::~AliPHOSRawFitterv2()
71 //-----------------------------------------------------------------------------
72 AliPHOSRawFitterv2::AliPHOSRawFitterv2(const AliPHOSRawFitterv2 &phosFitter ):
73 AliPHOSRawFitterv0(phosFitter),
78 fNtimeSamples=phosFitter.fNtimeSamples ;
79 for(Int_t i=0; i<3;i++){
80 fLGpar[i]=phosFitter.fLGpar[i] ;
81 fHGpar[i]=phosFitter.fHGpar[i] ;
83 fRMScut=phosFitter.fRMScut ;
86 //-----------------------------------------------------------------------------
87 AliPHOSRawFitterv2& AliPHOSRawFitterv2::operator = (const AliPHOSRawFitterv2 &phosFitter)
89 //Assignment operator.
91 fNtimeSamples=phosFitter.fNtimeSamples ;
92 for(Int_t i=0; i<3;i++){
93 fLGpar[i]=phosFitter.fLGpar[i] ;
94 fHGpar[i]=phosFitter.fHGpar[i] ;
96 fRMScut=phosFitter.fRMScut ;
100 //-----------------------------------------------------------------------------
101 Bool_t AliPHOSRawFitterv2::Eval(const UShort_t *signal, Int_t sigStart, Int_t sigLength)
103 //Extract an energy deposited in the crystal,
104 //crystal' position (module,column,row),
105 //time and gain (high or low).
106 //First collects sample, then evaluates it and if it has
107 //reasonable shape, fits it with Gamma2 function and extracts
115 const Float_t kBaseLine = 1.0;
116 const Int_t kPreSamples = 10;
121 UShort_t maxSample = 0;
124 TCanvas * cs = (TCanvas*)gROOT->FindObjectAny("CSample") ;
126 cs = new TCanvas("CSample","CSample") ;
128 TH1D * h = (TH1D*)gROOT->FindObjectAny("hSample") ;
129 if(!h) h = new TH1D("hSample","",200,0.,200.) ;
131 Double_t pedestal = 0;
132 for (Int_t i=0; i<sigLength; i++) {
135 pedMean += signal[i];
136 pedRMS += signal[i]*signal[i] ;
138 if(signal[i] > maxSample) maxSample = signal[i];
139 if(signal[i] == maxSample) nMax++;
140 h->SetBinContent(i+1,signal[i]) ;
142 fEnergy = (Double_t)maxSample;
143 if (maxSample > 900 && nMax > 2) fOverflow = kTRUE;
147 fPedestalRMS=(pedRMS - pedMean*pedMean/nPed)/nPed ;
148 if(fPedestalRMS > 0.)
149 fPedestalRMS = TMath::Sqrt(fPedestalRMS) ;
150 fEnergy -= (Double_t)(pedMean/nPed); // pedestal subtraction
156 //take pedestals from DB
157 pedestal = (Double_t) fAmpOffset ;
159 Float_t truePed = fCalibData->GetADCpedestalEmc(fModule, fCellZ, fCellX) ;
160 Int_t altroSettings = fCalibData->GetAltroOffsetEmc(fModule, fCellZ, fCellX) ;
161 pedestal += truePed - altroSettings ;
164 AliWarning(Form("Can not read data from OCDB")) ;
169 if (fEnergy < kBaseLine) {
179 Double_t a=0,b=0,c=0 ;
180 if(fCaloFlag == 0){ // Low gain
185 else if(fCaloFlag == 1){ // High gain
194 for(Int_t i=1; i<sigLength && cnts<fNtimeSamples; i++){
195 if(signal[i] < pedestal)
197 Double_t de = signal[i] - pedestal ;
198 Double_t av = signal[i-1] - pedestal + de;
199 if(av<=0.) //this is fluctuation around pedestal, skip it
201 Double_t ds = signal[i] - signal[i-1] ;
202 Double_t ti = ds/av ; // calculate log. derivative
203 ti = a/(ti+b)-c*ti ; // and compare with parameterization
205 Double_t wi = TMath::Abs(ds) ;
214 fQuality = tRMS/tW-fTime*fTime ;
223 for(Int_t i=1; i<sigLength-1&&!isBad; i++){
224 if(signal[i] > signal[i-1]+5 && signal[i] > signal[i+1]+5) { //single jump
231 if(fPedestalRMS > 0.1)
234 for(Int_t i=1; i<sigLength-1&&!isBad; i++){
235 if(signal[i] < pedestal-1)
239 if(fEnergy>10. && !isBad ){
240 printf("fE=%f, ped=%f, fQuality=%f, pedRMS=%f \n",fEnergy,pedestal,fQuality,pedRMS) ;
241 if(fOverflow)printf(" Overflow \n") ;
242 if(isBad)printf("bad") ;