Line data Source code
1 : // -*- mode: c++ -*-
2 : /**************************************************************************
3 : * This file is property of and copyright by *
4 : * the Relativistic Heavy Ion Group (RHIG), Yale University, US, 2009 *
5 : * *
6 : * Primary Author: Per Thomas Hille <perthomas.hille@yale.edu> *
7 : * *
8 : * Contributors are mentioned in the code where appropriate. *
9 : * Please report bugs to perthomas.hille@yale.edu *
10 : * *
11 : * Permission to use, copy, modify and distribute this software and its *
12 : * documentation strictly for non-commercial purposes is hereby granted *
13 : * without fee, provided that the above copyright notice appears in all *
14 : * copies and that both the copyright notice and this permission notice *
15 : * appear in the supporting documentation. The authors make no claims *
16 : * about the suitability of this software for any purpose. It is *
17 : * provided "as is" without express or implied warranty. *
18 : **************************************************************************/
19 :
20 : // The Peak-Finder algorithm
21 : // The amplitude is extracted as a
22 : // weighted sum of the samples using the
23 : // best possible weights.
24 : // The weights are calculated only once and the
25 : // actual extraction of amplitude and peak position
26 : // is done with a simple vector multiplication, allowing
27 : // extremely fast computations.
28 :
29 : #include "AliCaloRawAnalyzerPeakFinder.h"
30 : #include "AliCaloBunchInfo.h"
31 : #include "AliCaloFitResults.h"
32 : #include "TMath.h"
33 : #include "AliLog.h"
34 : #include "AliCDBEntry.h"
35 : #include "AliCDBManager.h"
36 : #include "TFile.h"
37 : #include "AliCaloPeakFinderVectors.h"
38 : #include <iostream>
39 :
40 : using namespace std;
41 :
42 42 : ClassImp( AliCaloRawAnalyzerPeakFinder )
43 :
44 :
45 0 : AliCaloRawAnalyzerPeakFinder::AliCaloRawAnalyzerPeakFinder() :AliCaloRawAnalyzer("Peak-Finder", "PF"),
46 0 : fPeakFinderVectors(0),
47 0 : fRunOnAlien(false),
48 0 : fIsInitialized(false)
49 0 : {
50 : // Ctor
51 0 : fAlgo= Algo::kPeakFinder;
52 0 : fPeakFinderVectors = new AliCaloPeakFinderVectors() ;
53 0 : ResetVectors();
54 0 : LoadVectorsOCDB();
55 0 : }
56 :
57 : void
58 : AliCaloRawAnalyzerPeakFinder::ResetVectors()
59 : {
60 : //As name implies
61 0 : for(int i=0; i < PF::MAXSTART; i++)
62 : {
63 0 : for(int j=0; j < PF::SAMPLERANGE; j++ )
64 : {
65 0 : for(int k=0; k < 100; k++ )
66 : {
67 0 : fPFAmpVectors[i][j][k] = 0;
68 0 : fPFTofVectors[i][j][k] = 0;
69 0 : fPFAmpVectorsCoarse[i][j][k] = 0;
70 0 : fPFTofVectorsCoarse[i][j][k] = 0;
71 : }
72 : }
73 : }
74 0 : }
75 :
76 :
77 : Double_t
78 : AliCaloRawAnalyzerPeakFinder::ScanCoarse(const Double_t *const array, const int length ) const
79 : {
80 : // Fisrt (coarce) estimate of Amplitude using the Peak-Finder.
81 : // The output of the first iteration is sued to select vectors
82 : // for the second iteration.
83 :
84 : Double_t tmpTof = 0;
85 : Double_t tmpAmp= 0;
86 :
87 0 : for(int i=0; i < length; i++)
88 : {
89 0 : tmpTof += fPFTofVectorsCoarse[0][length][i]*array[i];
90 0 : tmpAmp += fPFAmpVectorsCoarse[0][length][i]*array[i];
91 : }
92 :
93 0 : tmpTof = tmpTof / tmpAmp ;
94 0 : return tmpTof;
95 : }
96 :
97 :
98 : AliCaloFitResults
99 : AliCaloRawAnalyzerPeakFinder::Evaluate( const vector<AliCaloBunchInfo> &bunchvector, const UInt_t altrocfg1, const UInt_t altrocfg2 )
100 : {
101 : // Evaluation of amplitude and TOF
102 0 : if( fIsInitialized == false )
103 : {
104 0 : cout << __FILE__ << ":" << __LINE__ << "ERROR, peakfinder vectors not loaded" << endl;
105 0 : return AliCaloFitResults(kInvalid, kInvalid);
106 : }
107 :
108 : // Extracting the amplitude using the Peak-Finder algorithm
109 : // The amplitude is a weighted sum of the samples using
110 : // optimum weights.
111 :
112 0 : short maxampindex; //index of maximum amplitude
113 0 : short maxamp; //Maximum amplitude
114 0 : fAmp = 0;
115 0 : int index = SelectBunch( bunchvector, &maxampindex, &maxamp );
116 :
117 0 : if( index >= 0)
118 : {
119 0 : Float_t ped = ReverseAndSubtractPed( &(bunchvector.at(index)) , altrocfg1, altrocfg2, fReversed );
120 0 : Float_t maxf = TMath::MaxElement( bunchvector.at(index).GetLength(), fReversed );
121 0 : short timebinOffset = maxampindex - (bunchvector.at( index ).GetLength()-1);
122 0 : Float_t time = (timebinOffset*TIMEBINWITH)-fL1Phase;
123 0 : if( maxf < fAmpCut || maxamp > fOverflowCut ) // (maxamp - ped) > fOverflowCut = Close to saturation (use low gain then)
124 : {
125 0 : return AliCaloFitResults( maxamp, ped, Ret::kCrude, maxf, time, (int)time, 0, 0, Ret::kDummy);
126 : }
127 0 : else if ( maxf >= fAmpCut )
128 : {
129 0 : int first = 0;
130 0 : int last = 0;
131 0 : short maxrev = maxampindex - bunchvector.at(index).GetStartBin();
132 0 : SelectSubarray( fReversed, bunchvector.at(index).GetLength(), maxrev, &first, &last, fFitArrayCut);
133 0 : int nsamples = last - first;
134 :
135 0 : if( ( nsamples ) >= fNsampleCut ) // no if statement needed really; keep for readability
136 : {
137 0 : int startbin = bunchvector.at(index).GetStartBin();
138 0 : int n = last - first;
139 0 : int pfindex = n - fNsampleCut;
140 0 : pfindex = pfindex > PF::SAMPLERANGE ? PF::SAMPLERANGE : pfindex;
141 0 : int dt = maxampindex - startbin -2;
142 : int tmpindex = 0;
143 0 : Float_t tmptof = ScanCoarse( &fReversed[dt] , n );
144 :
145 0 : if( tmptof < -1 )
146 : {
147 : tmpindex = 0;
148 0 : }
149 : else
150 0 : if( tmptof > -1 && tmptof < 100 )
151 : {
152 : tmpindex = 1;
153 0 : }
154 : else
155 : {
156 : tmpindex = 2;
157 : }
158 :
159 : double tof = 0;
160 0 : for(int k=0; k < PF::SAMPLERANGE; k++ )
161 : {
162 0 : tof += fPFTofVectors[0][pfindex][k]*fReversed[ dt +k + tmpindex -1 ];
163 : }
164 0 : for( int i=0; i < PF::SAMPLERANGE; i++ )
165 : {
166 : {
167 :
168 0 : fAmp += fPFAmpVectors[0][pfindex][i]*fReversed[ dt +i +tmpindex -1 ];
169 : }
170 : }
171 :
172 0 : if( TMath::Abs( (maxf - fAmp )/maxf ) > 0.1 )
173 : {
174 0 : fAmp = maxf;
175 0 : }
176 :
177 0 : tof = timebinOffset - 0.01*tof/fAmp - fL1Phase/TIMEBINWITH; // clock
178 0 : Float_t chi2 = CalculateChi2(fAmp, tof-timebinOffset+maxrev, first, last);
179 0 : Int_t ndf = last - first - 1; // nsamples - 2
180 0 : Float_t toftime = (tof*TIMEBINWITH)-fL1Phase;
181 0 : return AliCaloFitResults( maxamp, ped , Ret::kFitPar, fAmp, toftime,
182 0 : (int)toftime, chi2, ndf,
183 0 : Ret::kDummy, AliCaloFitSubarray(index, maxrev, first, last) );
184 : }
185 : else
186 : {
187 0 : Float_t chi2 = CalculateChi2(maxf, maxrev, first, last);
188 0 : Int_t ndf = last - first - 1; // nsamples - 2
189 0 : return AliCaloFitResults( maxamp, ped , Ret::kCrude, maxf, time,
190 0 : (int)time, chi2, ndf, Ret::kDummy, AliCaloFitSubarray(index, maxrev, first, last) );
191 : }
192 0 : } // ampcut
193 0 : }
194 0 : return AliCaloFitResults(kInvalid, kInvalid);
195 0 : }
196 :
197 :
198 : void
199 : AliCaloRawAnalyzerPeakFinder::CopyVectors( const AliCaloPeakFinderVectors *const pfv )
200 : {
201 : // As name implies
202 0 : if ( pfv != 0)
203 : {
204 0 : for(int i = 0; i < PF::MAXSTART ; i++)
205 : {
206 0 : for( int j=0; j < PF::SAMPLERANGE; j++)
207 : {
208 0 : pfv->GetVector( i, j, fPFAmpVectors[i][j] , fPFTofVectors[i][j],
209 0 : fPFAmpVectorsCoarse[i][j] , fPFTofVectorsCoarse[i][j] );
210 :
211 0 : fPeakFinderVectors->SetVector( i, j, fPFAmpVectors[i][j], fPFTofVectors[i][j],
212 : fPFAmpVectorsCoarse[i][j], fPFTofVectorsCoarse[i][j] );
213 : }
214 : }
215 0 : }
216 : else
217 : {
218 0 : AliFatal( "pfv = ZERO !!!!!!!");
219 : }
220 0 : }
221 :
222 :
223 : void
224 : AliCaloRawAnalyzerPeakFinder::LoadVectorsOCDB()
225 : {
226 : //Loading of Peak-Finder vectors from the
227 : //Offline Condition Database (OCDB)
228 0 : AliCDBEntry* entry = AliCDBManager::Instance()->Get("EMCAL/Calib/PeakFinder/");
229 :
230 0 : if( entry != 0 )
231 : {
232 : //cout << __FILE__ << ":" << __LINE__ << ": Printing metadata !! " << endl;
233 : //entry->PrintMetaData();
234 0 : AliCaloPeakFinderVectors *pfv = (AliCaloPeakFinderVectors *)entry->GetObject();
235 0 : if( pfv == 0 )
236 : {
237 0 : cout << __FILE__ << ":" << __LINE__ << "_ ERRROR " << endl;
238 0 : }
239 0 : CopyVectors( pfv );
240 :
241 0 : if( pfv != 0 )
242 : {
243 0 : fIsInitialized = true;
244 0 : }
245 0 : }
246 0 : }
247 :
248 :
249 : void
250 : AliCaloRawAnalyzerPeakFinder::WriteRootFile() const
251 : { // Utility function to write Peak-Finder vectors to an root file
252 : // The output is used to create an OCDB entry.
253 0 : fPeakFinderVectors->PrintVectors();
254 0 : TFile *f = new TFile("peakfindervectors2.root", "recreate" );
255 0 : fPeakFinderVectors->Write();
256 0 : f->Close();
257 0 : delete f;
258 0 : }
259 :
260 :
261 : void
262 : AliCaloRawAnalyzerPeakFinder::PrintVectors()
263 : { // Utility function to write Peak-Finder vectors
264 0 : for(int i=0; i < 20; i++)
265 : {
266 0 : for( int j = 0; j < PF::MAXSTART; j ++ )
267 : {
268 0 : for( int k=0; k < PF::SAMPLERANGE; k++ )
269 : {
270 0 : cout << fPFAmpVectors[j][k][i] << "\t" ;
271 : }
272 : }
273 0 : cout << endl;
274 : }
275 0 : cout << __FILE__ << ":" << __LINE__ << ":.... DONE !!" << endl;
276 0 : }
|
