AliRoot Core  edcc906 (edcc906)
AliCaloRawAnalyzerKStandard.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * This file is property of and copyright by *
3  * the Relativistic Heavy Ion Group (RHIG), Yale University, US, 2009 *
4  * *
5  * Primary Author: Per Thomas Hille <p.t.hille@fys.uio.no> *
6  * *
7  * Contributors are mentioned in the code where appropriate. *
8  * Please report bugs to p.t.hille@fys.uio.no *
9  * *
10  * Permission to use, copy, modify and distribute this software and its *
11  * documentation strictly for non-commercial purposes is hereby granted *
12  * without fee, provided that the above copyright notice appears in all *
13  * copies and that both the copyright notice and this permission notice *
14  * appear in the supporting documentation. The authors make no claims *
15  * about the suitability of this software for any purpose. It is *
16  * provided "as is" without express or implied warranty. *
17  **************************************************************************/
18 
19 // AliRoot/EMCal system
21 #include "AliCaloBunchInfo.h"
22 #include "AliCaloFitResults.h"
23 #include "AliLog.h"
24 
25 // Standard libraries
26 #include <stdexcept>
27 #include <iostream>
28 
29 // Root system
30 #include "TF1.h"
31 #include "TGraph.h"
32 #include "TRandom.h"
33 #include "TMath.h"
34 
35 #include "AliEMCALRawResponse.h"
36 
37 using namespace std;
38 
40 ClassImp( AliCaloRawAnalyzerKStandard ) ;
42 
45 //_______________________________________________________________________
47 {
49 }
50 
53 //_______________________________________________________________________
55 {
56  // delete fTf1;
57 }
58 
61 //_______________________________________________________________________
63 AliCaloRawAnalyzerKStandard::Evaluate( const vector<AliCaloBunchInfo> &bunchlist,
64  UInt_t altrocfg1, UInt_t altrocfg2 )
65 {
66  Float_t pedEstimate = 0;
67  short maxADC = 0;
68  Int_t first = 0;
69  Int_t last = 0;
70  Int_t bunchIndex = 0;
71  Float_t ampEstimate = 0;
72  short timeEstimate = 0;
73  Float_t time = 0;
74  Float_t amp=0;
75  Float_t chi2 = 0;
76  Int_t ndf = 0;
77  Bool_t fitDone = kFALSE;
78 
79  int nsamples = PreFitEvaluateSamples( bunchlist, altrocfg1, altrocfg2, bunchIndex, ampEstimate,
80  maxADC, timeEstimate, pedEstimate, first, last, (int)fAmpCut );
81 
82 
83  if (ampEstimate >= fAmpCut )
84  {
85  time = timeEstimate;
86  Int_t timebinOffset = bunchlist.at(bunchIndex).GetStartBin() - (bunchlist.at(bunchIndex).GetLength()-1);
87  amp = ampEstimate;
88 
89  if ( nsamples > 1 && maxADC< OVERFLOWCUT )
90  {
91  FitRaw(first, last, amp, time, chi2, fitDone);
92  time += timebinOffset;
93  timeEstimate += timebinOffset;
94  ndf = nsamples - 2;
95  }
96  }
97  if ( fitDone )
98  {
99  Float_t ampAsymm = (amp - ampEstimate)/(amp + ampEstimate);
100  Float_t timeDiff = time - timeEstimate;
101 
102  if ( (TMath::Abs(ampAsymm) > 0.1) || (TMath::Abs(timeDiff) > 2) )
103  {
104  amp = ampEstimate;
105  time = timeEstimate;
106  fitDone = kFALSE;
107  }
108  }
109  if (amp >= fAmpCut )
110  {
111  if ( ! fitDone)
112  {
113  amp += (0.5 - gRandom->Rndm());
114  }
115  time = time * TIMEBINWITH;
116  time -= fL1Phase;
117 
118  return AliCaloFitResults( -99, -99, fAlgo , amp, time,
119  (int)time, chi2, ndf, Ret::kDummy );
120  }
122 }
123 
124 
127 //_______________________________________________________________________
128 void
129  AliCaloRawAnalyzerKStandard::FitRaw( Int_t firstTimeBin, Int_t lastTimeBin,
130  Float_t & amp, Float_t & time, Float_t & chi2, Bool_t & fitDone) const
131 {
132  int nsamples = lastTimeBin - firstTimeBin + 1;
133  fitDone = kFALSE;
134  if ( nsamples < 3 ) return;
135 
136  TGraph *gSig = new TGraph(nsamples);
137 
138  for (int i=0; i<nsamples; i++)
139  {
140  Int_t timebin = firstTimeBin + i;
141  gSig->SetPoint(i, timebin, GetReversed(timebin));
142  }
143 
144  TF1 * signalF = new TF1("signal", AliEMCALRawResponse::RawResponseFunction, 0, TIMEBINS , 5);
145 
146  signalF->SetParameters(10.,5., TAU ,ORDER,0.); //set all defaults once, just to be safe
147  signalF->SetParNames("amp","t0","tau","N","ped");
148  signalF->FixParameter(2,TAU);
149  signalF->FixParameter(3,ORDER);
150  signalF->FixParameter(4, 0);
151  signalF->SetParameter(1, time);
152  signalF->SetParameter(0, amp);
153  signalF->SetParLimits(0, 0.5*amp, 2*amp );
154  signalF->SetParLimits(1, time - 4, time + 4);
155 
156  try
157  {
158  gSig->Fit(signalF, "QROW"); // Note option 'W': equal errors on all points
159  amp = signalF->GetParameter(0);
160  time = signalF->GetParameter(1);
161  chi2 = signalF->GetChisquare();
162  fitDone = kTRUE;
163  }
164  catch (const std::exception & e)
165  {
166  AliError( Form("TGraph Fit exception %s", e.what()) );
167  // stay with default amp and time in case of exception, i.e. no special action required
168  fitDone = kFALSE;
169  }
170 
171  delete signalF;
172  delete gSig; // delete TGraph
173  return;
174 }
175 
176 
Float_t fAmpCut
Max ADC - pedestal must be higher than this befor attemting to extract the amplitude.
virtual AliCaloFitResults Evaluate(const std::vector< AliCaloBunchInfo > &bunchvector, UInt_t altrocfg1, UInt_t altrocfg2)
Evaluation Amplitude and TOF.
const double TIMEBINWITH
each sample is 100 ns
Container class to hold results from fitting.
const double TAU
Approximate shaping time.
void FitRaw(Int_t firstTimeBin, Int_t lastTimeBin, Float_t &amp, Float_t &time, Float_t &chi2, Bool_t &fitDone) const
Fits the raw signal time distribution.
Double_t chi2
Definition: AnalyzeLaser.C:7
Raw data fitters base class.
Algo::fitAlgorithm fAlgo
Which algorithm to use.
Double_t GetReversed(const int i) const
int PreFitEvaluateSamples(const std::vector< AliCaloBunchInfo > &bunchvector, UInt_t altrocfg1, UInt_t altrocfg2, Int_t &index, Float_t &maxf, short &maxamp, short &maxampindex, Float_t &ped, int &first, int &last, int acut)
Method to do the selection of what should possibly be fitted.
const int ORDER
Order of shaping stages of the signal conditioning unit.
Raw data fitting: standard TMinuit fit.
Double_t fL1Phase
Phase of the ADC sampling clock relative to the LHC clock.
const int TIMEBINS
number of sampling bins of the raw RO signal (we typically use 15-50; max is 1k+) ...
#define AliError(message)
Definition: AliLog.h:591
static Double_t RawResponseFunction(Double_t *x, Double_t *par)
virtual ~AliCaloRawAnalyzerKStandard()
Destructor.
const int OVERFLOWCUT
sample overflow