]>
Commit | Line | Data |
---|---|---|
ee299369 | 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$ */ | |
ee299369 | 17 | |
e5bbbc4e | 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. | |
ee299369 | 27 | //*-- Author: Marco van Leeuwen (LBL) |
e5bbbc4e | 28 | |
ee299369 | 29 | #include "AliEMCALRawUtils.h" |
507751ce | 30 | #include <stdexcept> |
21cad85c | 31 | |
4fe71e02 | 32 | #include "TF1.h" |
33 | #include "TGraph.h" | |
6751f762 | 34 | #include <TRandom.h> |
e5bbbc4e | 35 | class TSystem; |
21cad85c | 36 | |
e5bbbc4e | 37 | class AliLog; |
72c58de0 | 38 | #include "AliRun.h" |
ee299369 | 39 | #include "AliRunLoader.h" |
e5bbbc4e | 40 | class AliCaloAltroMapping; |
ee299369 | 41 | #include "AliAltroBuffer.h" |
42 | #include "AliRawReader.h" | |
32cd4c24 | 43 | #include "AliCaloRawStreamV3.h" |
ee299369 | 44 | #include "AliDAQ.h" |
21cad85c | 45 | |
feedcab9 | 46 | #include "AliEMCALRecParam.h" |
ee299369 | 47 | #include "AliEMCALLoader.h" |
48 | #include "AliEMCALGeometry.h" | |
e5bbbc4e | 49 | class AliEMCALDigitizer; |
ee299369 | 50 | #include "AliEMCALDigit.h" |
916f1e76 | 51 | #include "AliEMCALRawDigit.h" |
20b636fc | 52 | #include "AliEMCAL.h" |
5e3106bc | 53 | #include "AliCaloCalibPedestal.h" |
9f467289 | 54 | #include "AliCaloFastAltroFitv0.h" |
c8603a2b | 55 | #include "AliCaloNeuralFit.h" |
16605c06 | 56 | #include "AliCaloBunchInfo.h" |
57 | #include "AliCaloFitResults.h" | |
7683df1d | 58 | #include "AliCaloRawAnalyzerFastFit.h" |
59 | #include "AliCaloRawAnalyzerNN.h" | |
16605c06 | 60 | #include "AliCaloRawAnalyzerLMS.h" |
61 | #include "AliCaloRawAnalyzerPeakFinder.h" | |
62 | #include "AliCaloRawAnalyzerCrude.h" | |
9f467289 | 63 | |
ee299369 | 64 | ClassImp(AliEMCALRawUtils) |
21cad85c | 65 | |
ee299369 | 66 | // Signal shape parameters |
89d338a6 | 67 | Int_t AliEMCALRawUtils::fgTimeBins = 256; // number of sampling bins of the raw RO signal (we typically use 15-50; theoretical max is 1k+) |
e5bbbc4e | 68 | Double_t AliEMCALRawUtils::fgTimeBinWidth = 100E-9 ; // each sample is 100 ns |
09974781 | 69 | Double_t AliEMCALRawUtils::fgTimeTrigger = 1.5E-6 ; // 15 time bins ~ 1.5 musec |
ee299369 | 70 | |
71 | // some digitization constants | |
72 | Int_t AliEMCALRawUtils::fgThreshold = 1; | |
73 | Int_t AliEMCALRawUtils::fgDDLPerSuperModule = 2; // 2 ddls per SuperModule | |
916f1e76 | 74 | Int_t AliEMCALRawUtils::fgPedestalValue = 0; // pedestal value for digits2raw, default generate ZS data |
e5bbbc4e | 75 | Double_t AliEMCALRawUtils::fgFEENoise = 3.; // 3 ADC channels of noise (sampled) |
ee299369 | 76 | |
16605c06 | 77 | AliEMCALRawUtils::AliEMCALRawUtils(fitAlgorithm fitAlgo) |
b4133f05 | 78 | : fHighLowGainFactor(0.), fOrder(0), fTau(0.), fNoiseThreshold(0), |
9f467289 | 79 | fNPedSamples(0), fGeom(0), fOption(""), |
46f1d25f | 80 | fRemoveBadChannels(kTRUE),fFittingAlgorithm(0),fUseFALTRO(kFALSE),fRawAnalyzer(0) |
8cb998bd | 81 | { |
b4133f05 | 82 | |
83 | //These are default parameters. | |
84 | //Can be re-set from without with setter functions | |
9f467289 | 85 | //Already set in the OCDB and passed via setter in the AliEMCALReconstructor |
46f1d25f | 86 | fHighLowGainFactor = 16. ; // Adjusted for a low gain range of 82 GeV (10 bits) |
87 | fOrder = 2; // Order of gamma fn | |
88 | fTau = 2.35; // in units of timebin, from CERN 2007 testbeam | |
89 | fNoiseThreshold = 3; // 3 ADC counts is approx. noise level | |
90 | fNPedSamples = 4; // Less than this value => likely pedestal samples | |
91 | fRemoveBadChannels = kFALSE; // Do not remove bad channels before fitting | |
92 | fUseFALTRO = kTRUE; // Get the trigger FALTRO information and pass it to digits. | |
4fe71e02 | 93 | SetFittingAlgorithm(fitAlgo); |
16605c06 | 94 | |
65bdc82f | 95 | //Get Mapping RCU files from the AliEMCALRecParam |
96 | const TObjArray* maps = AliEMCALRecParam::GetMappings(); | |
97 | if(!maps) AliFatal("Cannot retrieve ALTRO mappings!!"); | |
98 | ||
21cad85c | 99 | for(Int_t i = 0; i < 4; i++) { |
65bdc82f | 100 | fMapping[i] = (AliAltroMapping*)maps->At(i); |
101 | } | |
102 | ||
72c58de0 | 103 | //To make sure we match with the geometry in a simulation file, |
104 | //let's try to get it first. If not, take the default geometry | |
33c3c91a | 105 | AliRunLoader *rl = AliRunLoader::Instance(); |
916f1e76 | 106 | if (rl && rl->GetAliRun() && rl->GetAliRun()->GetDetector("EMCAL")) { |
72c58de0 | 107 | fGeom = dynamic_cast<AliEMCAL*>(rl->GetAliRun()->GetDetector("EMCAL"))->GetGeometry(); |
108 | } else { | |
109 | AliInfo(Form("Using default geometry in raw reco")); | |
937d0661 | 110 | fGeom = AliEMCALGeometry::GetInstance(AliEMCALGeometry::GetDefaultGeometryName()); |
65bdc82f | 111 | } |
112 | ||
72c58de0 | 113 | if(!fGeom) AliFatal(Form("Could not get geometry!")); |
114 | ||
65bdc82f | 115 | } |
116 | ||
117 | //____________________________________________________________________________ | |
16605c06 | 118 | AliEMCALRawUtils::AliEMCALRawUtils(AliEMCALGeometry *pGeometry, fitAlgorithm fitAlgo) |
5544799a | 119 | : fHighLowGainFactor(0.), fOrder(0), fTau(0.), fNoiseThreshold(0), |
9f467289 | 120 | fNPedSamples(0), fGeom(pGeometry), fOption(""), |
46f1d25f | 121 | fRemoveBadChannels(kTRUE),fFittingAlgorithm(0),fUseFALTRO(kFALSE),fRawAnalyzer() |
5544799a | 122 | { |
123 | // | |
124 | // Initialize with the given geometry - constructor required by HLT | |
125 | // HLT does not use/support AliRunLoader(s) instances | |
126 | // This is a minimum intervention solution | |
127 | // Comment by MPloskon@lbl.gov | |
128 | // | |
129 | ||
130 | //These are default parameters. | |
131 | //Can be re-set from without with setter functions | |
9f467289 | 132 | //Already set in the OCDB and passed via setter in the AliEMCALReconstructor |
46f1d25f | 133 | fHighLowGainFactor = 16. ; // adjusted for a low gain range of 82 GeV (10 bits) |
134 | fOrder = 2; // order of gamma fn | |
135 | fTau = 2.35; // in units of timebin, from CERN 2007 testbeam | |
136 | fNoiseThreshold = 3; // 3 ADC counts is approx. noise level | |
137 | fNPedSamples = 4; // Less than this value => likely pedestal samples | |
138 | fRemoveBadChannels = kFALSE; // Do not remove bad channels before fitting | |
139 | fUseFALTRO = kTRUE; // Get the trigger FALTRO information and pass it to digits. | |
4fe71e02 | 140 | SetFittingAlgorithm(fitAlgo); |
46f1d25f | 141 | |
5544799a | 142 | //Get Mapping RCU files from the AliEMCALRecParam |
143 | const TObjArray* maps = AliEMCALRecParam::GetMappings(); | |
144 | if(!maps) AliFatal("Cannot retrieve ALTRO mappings!!"); | |
145 | ||
21cad85c | 146 | for(Int_t i = 0; i < 4; i++) { |
5544799a | 147 | fMapping[i] = (AliAltroMapping*)maps->At(i); |
148 | } | |
149 | ||
150 | if(!fGeom) AliFatal(Form("Could not get geometry!")); | |
151 | ||
152 | } | |
153 | ||
154 | //____________________________________________________________________________ | |
65bdc82f | 155 | AliEMCALRawUtils::AliEMCALRawUtils(const AliEMCALRawUtils& rawU) |
156 | : TObject(), | |
157 | fHighLowGainFactor(rawU.fHighLowGainFactor), | |
b4133f05 | 158 | fOrder(rawU.fOrder), |
159 | fTau(rawU.fTau), | |
160 | fNoiseThreshold(rawU.fNoiseThreshold), | |
161 | fNPedSamples(rawU.fNPedSamples), | |
65bdc82f | 162 | fGeom(rawU.fGeom), |
9f467289 | 163 | fOption(rawU.fOption), |
164 | fRemoveBadChannels(rawU.fRemoveBadChannels), | |
16605c06 | 165 | fFittingAlgorithm(rawU.fFittingAlgorithm), |
46f1d25f | 166 | fUseFALTRO(rawU.fUseFALTRO), |
16605c06 | 167 | fRawAnalyzer(rawU.fRawAnalyzer) |
65bdc82f | 168 | { |
169 | //copy ctor | |
170 | fMapping[0] = rawU.fMapping[0]; | |
171 | fMapping[1] = rawU.fMapping[1]; | |
21cad85c | 172 | fMapping[2] = rawU.fMapping[2]; |
173 | fMapping[3] = rawU.fMapping[3]; | |
65bdc82f | 174 | } |
175 | ||
176 | //____________________________________________________________________________ | |
177 | AliEMCALRawUtils& AliEMCALRawUtils::operator =(const AliEMCALRawUtils &rawU) | |
178 | { | |
179 | //assignment operator | |
180 | ||
181 | if(this != &rawU) { | |
182 | fHighLowGainFactor = rawU.fHighLowGainFactor; | |
46f1d25f | 183 | fOrder = rawU.fOrder; |
184 | fTau = rawU.fTau; | |
185 | fNoiseThreshold = rawU.fNoiseThreshold; | |
186 | fNPedSamples = rawU.fNPedSamples; | |
187 | fGeom = rawU.fGeom; | |
188 | fOption = rawU.fOption; | |
9f467289 | 189 | fRemoveBadChannels = rawU.fRemoveBadChannels; |
190 | fFittingAlgorithm = rawU.fFittingAlgorithm; | |
46f1d25f | 191 | fUseFALTRO = rawU.fUseFALTRO; |
192 | fRawAnalyzer = rawU.fRawAnalyzer; | |
193 | fMapping[0] = rawU.fMapping[0]; | |
194 | fMapping[1] = rawU.fMapping[1]; | |
195 | fMapping[2] = rawU.fMapping[2]; | |
196 | fMapping[3] = rawU.fMapping[3]; | |
65bdc82f | 197 | } |
198 | ||
199 | return *this; | |
200 | ||
ee299369 | 201 | } |
65bdc82f | 202 | |
ee299369 | 203 | //____________________________________________________________________________ |
204 | AliEMCALRawUtils::~AliEMCALRawUtils() { | |
e5bbbc4e | 205 | //dtor |
65bdc82f | 206 | |
ee299369 | 207 | } |
65bdc82f | 208 | |
ee299369 | 209 | //____________________________________________________________________________ |
65bdc82f | 210 | void AliEMCALRawUtils::Digits2Raw() |
ee299369 | 211 | { |
212 | // convert digits of the current event to raw data | |
213 | ||
33c3c91a | 214 | AliRunLoader *rl = AliRunLoader::Instance(); |
ee299369 | 215 | AliEMCALLoader *loader = dynamic_cast<AliEMCALLoader*>(rl->GetDetectorLoader("EMCAL")); |
216 | ||
217 | // get the digits | |
218 | loader->LoadDigits("EMCAL"); | |
219 | loader->GetEvent(); | |
220 | TClonesArray* digits = loader->Digits() ; | |
221 | ||
222 | if (!digits) { | |
223 | Warning("Digits2Raw", "no digits found !"); | |
224 | return; | |
225 | } | |
65bdc82f | 226 | |
ee299369 | 227 | static const Int_t nDDL = 12*2; // 12 SM hardcoded for now. Buffers allocated dynamically, when needed, so just need an upper limit here |
228 | AliAltroBuffer* buffers[nDDL]; | |
229 | for (Int_t i=0; i < nDDL; i++) | |
230 | buffers[i] = 0; | |
231 | ||
e2c2134b | 232 | TArrayI adcValuesLow(fgTimeBins); |
233 | TArrayI adcValuesHigh(fgTimeBins); | |
ee299369 | 234 | |
ee299369 | 235 | // loop over digits (assume ordered digits) |
236 | for (Int_t iDigit = 0; iDigit < digits->GetEntries(); iDigit++) { | |
237 | AliEMCALDigit* digit = dynamic_cast<AliEMCALDigit *>(digits->At(iDigit)) ; | |
238 | if (digit->GetAmp() < fgThreshold) | |
239 | continue; | |
240 | ||
241 | //get cell indices | |
242 | Int_t nSM = 0; | |
243 | Int_t nIphi = 0; | |
244 | Int_t nIeta = 0; | |
245 | Int_t iphi = 0; | |
246 | Int_t ieta = 0; | |
247 | Int_t nModule = 0; | |
65bdc82f | 248 | fGeom->GetCellIndex(digit->GetId(), nSM, nModule, nIphi, nIeta); |
249 | fGeom->GetCellPhiEtaIndexInSModule(nSM, nModule, nIphi, nIeta,iphi, ieta) ; | |
ee299369 | 250 | |
21cad85c | 251 | //Check which is the RCU, 0 or 1, of the cell. |
ee299369 | 252 | Int_t iRCU = -111; |
253 | //RCU0 | |
254 | if (0<=iphi&&iphi<8) iRCU=0; // first cable row | |
255 | else if (8<=iphi&&iphi<16 && 0<=ieta&&ieta<24) iRCU=0; // first half; | |
256 | //second cable row | |
257 | //RCU1 | |
258 | else if(8<=iphi&&iphi<16 && 24<=ieta&&ieta<48) iRCU=1; // second half; | |
259 | //second cable row | |
260 | else if(16<=iphi&&iphi<24) iRCU=1; // third cable row | |
21cad85c | 261 | |
262 | if (nSM%2==1) iRCU = 1 - iRCU; // swap for odd=C side, to allow us to cable both sides the same | |
263 | ||
e36e3bcf | 264 | if (iRCU<0) |
265 | Fatal("Digits2Raw()","Non-existent RCU number: %d", iRCU); | |
ee299369 | 266 | |
267 | //Which DDL? | |
268 | Int_t iDDL = fgDDLPerSuperModule* nSM + iRCU; | |
269 | if (iDDL >= nDDL) | |
270 | Fatal("Digits2Raw()","Non-existent DDL board number: %d", iDDL); | |
271 | ||
272 | if (buffers[iDDL] == 0) { | |
273 | // open new file and write dummy header | |
274 | TString fileName = AliDAQ::DdlFileName("EMCAL",iDDL); | |
21cad85c | 275 | //Select mapping file RCU0A, RCU0C, RCU1A, RCU1C |
276 | Int_t iRCUside=iRCU+(nSM%2)*2; | |
277 | //iRCU=0 and even (0) SM -> RCU0A.data 0 | |
278 | //iRCU=1 and even (0) SM -> RCU1A.data 1 | |
279 | //iRCU=0 and odd (1) SM -> RCU0C.data 2 | |
280 | //iRCU=1 and odd (1) SM -> RCU1C.data 3 | |
281 | //cout<<" nSM "<<nSM<<"; iRCU "<<iRCU<<"; iRCUside "<<iRCUside<<endl; | |
282 | buffers[iDDL] = new AliAltroBuffer(fileName.Data(),fMapping[iRCUside]); | |
ee299369 | 283 | buffers[iDDL]->WriteDataHeader(kTRUE, kFALSE); //Dummy; |
284 | } | |
285 | ||
286 | // out of time range signal (?) | |
287 | if (digit->GetTimeR() > GetRawFormatTimeMax() ) { | |
288 | AliInfo("Signal is out of time range.\n"); | |
289 | buffers[iDDL]->FillBuffer((Int_t)digit->GetAmp()); | |
290 | buffers[iDDL]->FillBuffer(GetRawFormatTimeBins() ); // time bin | |
291 | buffers[iDDL]->FillBuffer(3); // bunch length | |
292 | buffers[iDDL]->WriteTrailer(3, ieta, iphi, nSM); // trailer | |
293 | // calculate the time response function | |
294 | } else { | |
e2c2134b | 295 | Bool_t lowgain = RawSampledResponse(digit->GetTimeR(), digit->GetAmp(), adcValuesHigh.GetArray(), adcValuesLow.GetArray()) ; |
ee299369 | 296 | if (lowgain) |
e2c2134b | 297 | buffers[iDDL]->WriteChannel(ieta, iphi, 0, GetRawFormatTimeBins(), adcValuesLow.GetArray(), fgThreshold); |
ee299369 | 298 | else |
e2c2134b | 299 | buffers[iDDL]->WriteChannel(ieta,iphi, 1, GetRawFormatTimeBins(), adcValuesHigh.GetArray(), fgThreshold); |
ee299369 | 300 | } |
301 | } | |
302 | ||
303 | // write headers and close files | |
304 | for (Int_t i=0; i < nDDL; i++) { | |
305 | if (buffers[i]) { | |
306 | buffers[i]->Flush(); | |
307 | buffers[i]->WriteDataHeader(kFALSE, kFALSE); | |
308 | delete buffers[i]; | |
309 | } | |
310 | } | |
65bdc82f | 311 | |
ee299369 | 312 | loader->UnloadDigits(); |
313 | } | |
314 | ||
315 | //____________________________________________________________________________ | |
916f1e76 | 316 | void AliEMCALRawUtils::Raw2Digits(AliRawReader* reader,TClonesArray *digitsArr, const AliCaloCalibPedestal* pedbadmap, TClonesArray *digitsTRG) |
ee299369 | 317 | { |
65bdc82f | 318 | // convert raw data of the current event to digits |
ee299369 | 319 | |
c47157cd | 320 | digitsArr->Clear(); |
ee299369 | 321 | |
c47157cd | 322 | if (!digitsArr) { |
ee299369 | 323 | Error("Raw2Digits", "no digits found !"); |
324 | return; | |
325 | } | |
326 | if (!reader) { | |
327 | Error("Raw2Digits", "no raw reader found !"); | |
328 | return; | |
329 | } | |
330 | ||
32cd4c24 | 331 | AliCaloRawStreamV3 in(reader,"EMCAL",fMapping); |
ee299369 | 332 | // Select EMCAL DDL's; |
7643e728 | 333 | reader->Select("EMCAL",0,43); // 43 = AliEMCALGeoParams::fgkLastAltroDDL |
feedcab9 | 334 | |
16605c06 | 335 | // fRawAnalyzer setup |
f57baa2d | 336 | fRawAnalyzer->SetNsampleCut(5); // requirement for fits to be done |
16605c06 | 337 | fRawAnalyzer->SetAmpCut(fNoiseThreshold); |
338 | fRawAnalyzer->SetFitArrayCut(fNoiseThreshold); | |
339 | fRawAnalyzer->SetIsZeroSuppressed(true); // TMP - should use stream->IsZeroSuppressed(), or altro cfg registers later | |
ee299369 | 340 | |
16605c06 | 341 | // channel info parameters |
ee299369 | 342 | Int_t lowGain = 0; |
e5bbbc4e | 343 | Int_t caloFlag = 0; // low, high gain, or TRU, or LED ref. |
ee299369 | 344 | |
32cd4c24 | 345 | // start loop over input stream |
346 | while (in.NextDDL()) { | |
916f1e76 | 347 | |
348 | // if ( in.GetDDLNumber() != 0 && in.GetDDLNumber() != 2 ) continue; | |
349 | ||
32cd4c24 | 350 | while (in.NextChannel()) { |
7643e728 | 351 | |
916f1e76 | 352 | /* |
353 | Int_t hhwAdd = in.GetHWAddress(); | |
354 | UShort_t iiBranch = ( hhwAdd >> 11 ) & 0x1; // 0/1 | |
355 | UShort_t iiFEC = ( hhwAdd >> 7 ) & 0xF; | |
356 | UShort_t iiChip = ( hhwAdd >> 4 ) & 0x7; | |
357 | UShort_t iiChannel = hhwAdd & 0xF; | |
358 | ||
359 | if ( !( iiBranch == 0 && iiFEC == 1 && iiChip == 3 && ( iiChannel >= 8 && iiChannel <= 15 ) ) && !( iiBranch == 1 && iiFEC == 0 && in.GetColumn() == 0 ) ) continue; | |
360 | */ | |
361 | ||
7643e728 | 362 | //Check if the signal is high or low gain and then do the fit, |
16605c06 | 363 | //if it is from TRU or LEDMon do not fit |
7643e728 | 364 | caloFlag = in.GetCaloFlag(); |
916f1e76 | 365 | // if (caloFlag != 0 && caloFlag != 1) continue; |
366 | if (caloFlag > 2) continue; // Work with ALTRO and FALTRO | |
367 | ||
368 | //Do not fit bad channels of ALTRO | |
369 | if(caloFlag < 2 && fRemoveBadChannels && pedbadmap->IsBadChannel(in.GetModule(),in.GetColumn(),in.GetRow())) { | |
5e3106bc | 370 | //printf("Tower from SM %d, column %d, row %d is BAD!!! Skip \n", in.GetModule(),in.GetColumn(),in.GetRow()); |
371 | continue; | |
372 | } | |
373 | ||
16605c06 | 374 | vector<AliCaloBunchInfo> bunchlist; |
32cd4c24 | 375 | while (in.NextBunch()) { |
16605c06 | 376 | bunchlist.push_back( AliCaloBunchInfo(in.GetStartTimeBin(), in.GetBunchLength(), in.GetSignals() ) ); |
377 | } // loop over bunches | |
7643e728 | 378 | |
916f1e76 | 379 | |
507751ce | 380 | if ( caloFlag < 2 ){ // ALTRO |
916f1e76 | 381 | |
507751ce | 382 | Float_t time = 0; |
383 | Float_t amp = 0; | |
384 | short timeEstimate = 0; | |
385 | Float_t ampEstimate = 0; | |
386 | Bool_t fitDone = kFALSE; | |
916f1e76 | 387 | |
7683df1d | 388 | if ( fFittingAlgorithm == kFastFit || fFittingAlgorithm == kNeuralNet || fFittingAlgorithm == kLMS || fFittingAlgorithm == kPeakFinder || fFittingAlgorithm == kCrude) { |
16605c06 | 389 | // all functionality to determine amp and time etc is encapsulated inside the Evaluate call for these methods |
390 | AliCaloFitResults fitResults = fRawAnalyzer->Evaluate( bunchlist, in.GetAltroCFG1(), in.GetAltroCFG2()); | |
391 | ||
392 | amp = fitResults.GetAmp(); | |
507751ce | 393 | time = fitResults.GetTime(); |
394 | timeEstimate = fitResults.GetMaxTimebin(); | |
395 | ampEstimate = fitResults.GetMaxSig(); | |
396 | if (fitResults.GetStatus() == AliCaloFitResults::kFitPar) { | |
397 | fitDone = kTRUE; | |
398 | } | |
16605c06 | 399 | } |
400 | else { // for the other methods we for now use the functionality of | |
401 | // AliCaloRawAnalyzer as well, to select samples and prepare for fits, | |
402 | // if it looks like there is something to fit | |
403 | ||
404 | // parameters init. | |
507751ce | 405 | Float_t pedEstimate = 0; |
16605c06 | 406 | short maxADC = 0; |
16605c06 | 407 | Int_t first = 0; |
408 | Int_t last = 0; | |
409 | Int_t bunchIndex = 0; | |
410 | // | |
411 | // The PreFitEvaluateSamples + later call to FitRaw will hopefully | |
412 | // be replaced by a single Evaluate call or so soon, like for the other | |
413 | // methods, but this should be good enough for evaluation of | |
414 | // the methods for now (Jan. 2010) | |
415 | // | |
416 | int nsamples = fRawAnalyzer->PreFitEvaluateSamples( bunchlist, in.GetAltroCFG1(), in.GetAltroCFG2(), bunchIndex, ampEstimate, maxADC, timeEstimate, pedEstimate, first, last); | |
7643e728 | 417 | |
16605c06 | 418 | if (ampEstimate > fNoiseThreshold) { // something worth looking at |
7643e728 | 419 | |
f57baa2d | 420 | time = timeEstimate; // maxrev in AliCaloRawAnalyzer speak; comes with an offset w.r.t. real timebin |
421 | Int_t timebinOffset = bunchlist.at(bunchIndex).GetStartBin() - (bunchlist.at(bunchIndex).GetLength()-1); | |
16605c06 | 422 | amp = ampEstimate; |
423 | ||
424 | if ( nsamples > 1 ) { // possibly something to fit | |
507751ce | 425 | FitRaw(first, last, amp, time, fitDone); |
f57baa2d | 426 | time += timebinOffset; |
507751ce | 427 | timeEstimate += timebinOffset; |
9f467289 | 428 | } |
16605c06 | 429 | |
16605c06 | 430 | } // ampEstimate check |
431 | } // method selection | |
507751ce | 432 | |
433 | if ( fitDone ) { // brief sanity check of fit results | |
434 | Float_t ampAsymm = (amp - ampEstimate)/(amp + ampEstimate); | |
435 | Float_t timeDiff = time - timeEstimate; | |
436 | if ( (TMath::Abs(ampAsymm) > 0.1) || (TMath::Abs(timeDiff) > 2) ) { | |
437 | // AliDebug(2,Form("Fit results amp %f time %f not consistent with expectations amp %f time %d", amp, time, ampEstimate, timeEstimate)); | |
438 | ||
439 | // for now just overwrite the fit results with the simple/initial estimate | |
440 | amp = ampEstimate; | |
441 | time = timeEstimate; | |
442 | fitDone = kFALSE; | |
443 | } | |
444 | } // fitDone | |
16605c06 | 445 | |
029fe7a2 | 446 | if (amp > fNoiseThreshold && amp<fgkRawSignalOverflow) { // something to be stored |
507751ce | 447 | if ( ! fitDone) { // smear ADC with +- 0.5 uniform (avoid discrete effects) |
448 | amp += (0.5 - gRandom->Rndm()); // Rndm generates a number in ]0,1] | |
449 | } | |
450 | ||
16605c06 | 451 | Int_t id = fGeom->GetAbsCellIdFromCellIndexes(in.GetModule(), in.GetRow(), in.GetColumn()) ; |
7643e728 | 452 | lowGain = in.IsLowGain(); |
453 | ||
16605c06 | 454 | // go from time-bin units to physical time fgtimetrigger |
455 | time = time * GetRawFormatTimeBinWidth(); // skip subtraction of fgTimeTrigger? | |
7643e728 | 456 | |
457 | AliDebug(2,Form("id %d lowGain %d amp %g", id, lowGain, amp)); | |
458 | // printf("Added tower: SM %d, row %d, column %d, amp %3.2f\n",in.GetModule(), in.GetRow(), in.GetColumn(),amp); | |
459 | // round off amplitude value to nearest integer | |
460 | AddDigit(digitsArr, id, lowGain, TMath::Nint(amp), time); | |
461 | } | |
462 | ||
916f1e76 | 463 | }//ALTRO |
46f1d25f | 464 | else if(fUseFALTRO) |
916f1e76 | 465 | {// Fake ALTRO |
466 | // if (maxTimeBin && gSig->GetN() > maxTimeBin + 10) gSig->Set(maxTimeBin + 10); // set actual max size of TGraph | |
916f1e76 | 467 | Int_t hwAdd = in.GetHWAddress(); |
468 | UShort_t iRCU = in.GetDDLNumber() % 2; // 0/1 | |
469 | UShort_t iBranch = ( hwAdd >> 11 ) & 0x1; // 0/1 | |
470 | ||
471 | // Now find TRU number | |
472 | Int_t itru = 3 * in.GetModule() + ( (iRCU << 1) | iBranch ) - 1; | |
473 | ||
474 | AliDebug(1,Form("Found TRG digit in TRU: %2d ADC: %2d",itru,in.GetColumn())); | |
475 | ||
476 | Int_t idtrg; | |
477 | ||
478 | Bool_t isOK = fGeom->GetAbsFastORIndexFromTRU(itru, in.GetColumn(), idtrg); | |
479 | ||
480 | Int_t timeSamples[256]; for (Int_t j=0;j<256;j++) timeSamples[j] = 0; | |
481 | Int_t nSamples = 0; | |
482 | ||
483 | for (std::vector<AliCaloBunchInfo>::iterator itVectorData = bunchlist.begin(); itVectorData != bunchlist.end(); itVectorData++) | |
484 | { | |
485 | AliCaloBunchInfo bunch = *(itVectorData); | |
486 | ||
487 | const UShort_t* sig = bunch.GetData(); | |
488 | Int_t startBin = bunch.GetStartBin(); | |
489 | ||
490 | for (Int_t iS = 0; iS < bunch.GetLength(); iS++) | |
491 | { | |
492 | Int_t time = startBin--; | |
493 | Int_t amp = sig[iS]; | |
494 | ||
495 | if ( amp ) timeSamples[nSamples++] = ( ( time << 12 ) & 0xFF000 ) | ( amp & 0xFFF ); | |
496 | } | |
497 | } | |
498 | ||
499 | if (nSamples && isOK) AddDigit(digitsTRG, idtrg, timeSamples, nSamples); | |
500 | }//Fake ALTRO | |
32cd4c24 | 501 | } // end while over channel |
502 | } //end while over DDL's, of input stream | |
16605c06 | 503 | |
ee299369 | 504 | return ; |
505 | } | |
506 | ||
916f1e76 | 507 | //____________________________________________________________________________ |
508 | void AliEMCALRawUtils::AddDigit(TClonesArray *digitsArr, Int_t id, Int_t timeSamples[], Int_t nSamples) | |
509 | { | |
510 | new((*digitsArr)[digitsArr->GetEntriesFast()]) AliEMCALRawDigit(id, timeSamples, nSamples); | |
511 | ||
512 | // Int_t idx = digitsArr->GetEntriesFast()-1; | |
513 | // AliEMCALRawDigit* d = (AliEMCALRawDigit*)digitsArr->At(idx); | |
514 | } | |
515 | ||
ee299369 | 516 | //____________________________________________________________________________ |
82cbdfca | 517 | void AliEMCALRawUtils::AddDigit(TClonesArray *digitsArr, Int_t id, Int_t lowGain, Int_t amp, Float_t time) { |
518 | // | |
519 | // Add a new digit. | |
520 | // This routine checks whether a digit exists already for this tower | |
521 | // and then decides whether to use the high or low gain info | |
522 | // | |
523 | // Called by Raw2Digits | |
524 | ||
525 | AliEMCALDigit *digit = 0, *tmpdigit = 0; | |
82cbdfca | 526 | TIter nextdigit(digitsArr); |
527 | while (digit == 0 && (tmpdigit = (AliEMCALDigit*) nextdigit())) { | |
528 | if (tmpdigit->GetId() == id) | |
529 | digit = tmpdigit; | |
530 | } | |
531 | ||
532 | if (!digit) { // no digit existed for this tower; create one | |
a7ec7165 | 533 | if (lowGain && amp > fgkOverflowCut) |
82cbdfca | 534 | amp = Int_t(fHighLowGainFactor * amp); |
535 | Int_t idigit = digitsArr->GetEntries(); | |
536 | new((*digitsArr)[idigit]) AliEMCALDigit( -1, -1, id, amp, time, idigit) ; | |
537 | } | |
538 | else { // a digit already exists, check range | |
b4133f05 | 539 | // (use high gain if signal < cut value, otherwise low gain) |
82cbdfca | 540 | if (lowGain) { // new digit is low gain |
b4133f05 | 541 | if (digit->GetAmp() > fgkOverflowCut) { // use if stored digit is out of range |
82cbdfca | 542 | digit->SetAmp(Int_t(fHighLowGainFactor * amp)); |
543 | digit->SetTime(time); | |
544 | } | |
545 | } | |
b4133f05 | 546 | else if (amp < fgkOverflowCut) { // new digit is high gain; use if not out of range |
82cbdfca | 547 | digit->SetAmp(amp); |
548 | digit->SetTime(time); | |
549 | } | |
550 | } | |
551 | } | |
552 | ||
553 | //____________________________________________________________________________ | |
507751ce | 554 | void AliEMCALRawUtils::FitRaw(const Int_t firstTimeBin, const Int_t lastTimeBin, Float_t & amp, Float_t & time, Bool_t & fitDone) const |
16605c06 | 555 | { // Fits the raw signal time distribution |
556 | ||
557 | //-------------------------------------------------- | |
558 | //Do the fit, different fitting algorithms available | |
559 | //-------------------------------------------------- | |
560 | int nsamples = lastTimeBin - firstTimeBin + 1; | |
507751ce | 561 | fitDone = kFALSE; |
ee299369 | 562 | |
16605c06 | 563 | switch(fFittingAlgorithm) { |
564 | case kStandard: | |
565 | { | |
7683df1d | 566 | if (nsamples < 3) { return; } // nothing much to fit |
16605c06 | 567 | //printf("Standard fitter \n"); |
7683df1d | 568 | |
16605c06 | 569 | // Create Graph to hold data we will fit |
7683df1d | 570 | TGraph *gSig = new TGraph( nsamples); |
571 | for (int i=0; i<nsamples; i++) { | |
572 | Int_t timebin = firstTimeBin + i; | |
f57baa2d | 573 | gSig->SetPoint(i, timebin, fRawAnalyzer->GetReversed(timebin)); |
7683df1d | 574 | } |
575 | ||
16605c06 | 576 | TF1 * signalF = new TF1("signal", RawResponseFunction, 0, GetRawFormatTimeBins(), 5); |
577 | signalF->SetParameters(10.,5.,fTau,fOrder,0.); //set all defaults once, just to be safe | |
578 | signalF->SetParNames("amp","t0","tau","N","ped"); | |
579 | signalF->FixParameter(2,fTau); // tau in units of time bin | |
580 | signalF->FixParameter(3,fOrder); // order | |
581 | signalF->FixParameter(4, 0); // pedestal should be subtracted when we get here | |
582 | signalF->SetParameter(1, time); | |
583 | signalF->SetParameter(0, amp); | |
507751ce | 584 | // set rather loose parameter limits |
585 | signalF->SetParLimits(0, 0.5*amp, 2*amp ); | |
586 | signalF->SetParLimits(1, time - 4, time + 4); | |
587 | ||
588 | try { | |
589 | gSig->Fit(signalF, "QROW"); // Note option 'W': equal errors on all points | |
590 | // assign fit results | |
591 | amp = signalF->GetParameter(0); | |
592 | time = signalF->GetParameter(1); | |
593 | ||
594 | // cross-check with ParabolaFit to see if the results make sense | |
595 | FitParabola(gSig, amp); // amp is possibly updated | |
596 | fitDone = kTRUE; | |
597 | } | |
598 | catch (const std::exception & e) { | |
599 | AliError( Form("TGraph Fit exception %s", e.what()) ); | |
600 | // stay with default amp and time in case of exception, i.e. no special action required | |
601 | fitDone = kFALSE; | |
602 | } | |
16605c06 | 603 | delete signalF; |
604 | ||
16605c06 | 605 | //printf("Std : Amp %f, time %g\n",amp, time); |
7683df1d | 606 | delete gSig; // delete TGraph |
16605c06 | 607 | |
608 | break; | |
609 | }//kStandard Fitter | |
610 | //---------------------------- | |
7683df1d | 611 | case kLogFit: |
16605c06 | 612 | { |
7683df1d | 613 | if (nsamples < 3) { return; } // nothing much to fit |
614 | //printf("LogFit \n"); | |
615 | ||
616 | // Create Graph to hold data we will fit | |
617 | TGraph *gSigLog = new TGraph( nsamples); | |
618 | for (int i=0; i<nsamples; i++) { | |
619 | Int_t timebin = firstTimeBin + i; | |
620 | gSigLog->SetPoint(timebin, timebin, TMath::Log(fRawAnalyzer->GetReversed(timebin) ) ); | |
7643e728 | 621 | } |
7683df1d | 622 | |
623 | TF1 * signalFLog = new TF1("signalLog", RawResponseFunctionLog, 0, GetRawFormatTimeBins(), 5); | |
624 | signalFLog->SetParameters(2.3, 5.,fTau,fOrder,0.); //set all defaults once, just to be safe | |
625 | signalFLog->SetParNames("amplog","t0","tau","N","ped"); | |
626 | signalFLog->FixParameter(2,fTau); // tau in units of time bin | |
627 | signalFLog->FixParameter(3,fOrder); // order | |
628 | signalFLog->FixParameter(4, 0); // pedestal should be subtracted when we get here | |
629 | signalFLog->SetParameter(1, time); | |
630 | if (amp>=1) { | |
631 | signalFLog->SetParameter(0, TMath::Log(amp)); | |
16605c06 | 632 | } |
7683df1d | 633 | |
634 | gSigLog->Fit(signalFLog, "QROW"); // Note option 'W': equal errors on all points | |
635 | ||
636 | // assign fit results | |
637 | Double_t amplog = signalFLog->GetParameter(0); //Not Amp, but Log of Amp | |
638 | amp = TMath::Exp(amplog); | |
639 | time = signalFLog->GetParameter(1); | |
507751ce | 640 | fitDone = kTRUE; |
7683df1d | 641 | |
642 | delete signalFLog; | |
643 | //printf("LogFit: Amp %f, time %g\n",amp, time); | |
644 | delete gSigLog; | |
16605c06 | 645 | break; |
7683df1d | 646 | } //kLogFit |
647 | //---------------------------- | |
648 | ||
16605c06 | 649 | //---------------------------- |
650 | }//switch fitting algorithms | |
fb070798 | 651 | |
16605c06 | 652 | return; |
653 | } | |
8cb998bd | 654 | |
16605c06 | 655 | //__________________________________________________________________ |
656 | void AliEMCALRawUtils::FitParabola(const TGraph *gSig, Float_t & amp) const | |
657 | { | |
658 | //BEG YS alternative methods to calculate the amplitude | |
659 | Double_t * ymx = gSig->GetX() ; | |
660 | Double_t * ymy = gSig->GetY() ; | |
661 | const Int_t kN = 3 ; | |
662 | Double_t ymMaxX[kN] = {0., 0., 0.} ; | |
663 | Double_t ymMaxY[kN] = {0., 0., 0.} ; | |
664 | Double_t ymax = 0. ; | |
665 | // find the maximum amplitude | |
666 | Int_t ymiMax = 0 ; | |
667 | for (Int_t ymi = 0; ymi < gSig->GetN(); ymi++) { | |
668 | if (ymy[ymi] > ymMaxY[0] ) { | |
669 | ymMaxY[0] = ymy[ymi] ; //<========== This is the maximum amplitude | |
670 | ymMaxX[0] = ymx[ymi] ; | |
671 | ymiMax = ymi ; | |
672 | } | |
673 | } | |
674 | // find the maximum by fitting a parabola through the max and the two adjacent samples | |
675 | if ( ymiMax < gSig->GetN()-1 && ymiMax > 0) { | |
676 | ymMaxY[1] = ymy[ymiMax+1] ; | |
677 | ymMaxY[2] = ymy[ymiMax-1] ; | |
678 | ymMaxX[1] = ymx[ymiMax+1] ; | |
679 | ymMaxX[2] = ymx[ymiMax-1] ; | |
680 | if (ymMaxY[0]*ymMaxY[1]*ymMaxY[2] > 0) { | |
681 | //fit a parabola through the 3 points y= a+bx+x*x*x | |
682 | Double_t sy = 0 ; | |
683 | Double_t sx = 0 ; | |
684 | Double_t sx2 = 0 ; | |
685 | Double_t sx3 = 0 ; | |
686 | Double_t sx4 = 0 ; | |
687 | Double_t sxy = 0 ; | |
688 | Double_t sx2y = 0 ; | |
689 | for (Int_t i = 0; i < kN ; i++) { | |
690 | sy += ymMaxY[i] ; | |
691 | sx += ymMaxX[i] ; | |
692 | sx2 += ymMaxX[i]*ymMaxX[i] ; | |
693 | sx3 += ymMaxX[i]*ymMaxX[i]*ymMaxX[i] ; | |
694 | sx4 += ymMaxX[i]*ymMaxX[i]*ymMaxX[i]*ymMaxX[i] ; | |
695 | sxy += ymMaxX[i]*ymMaxY[i] ; | |
696 | sx2y += ymMaxX[i]*ymMaxX[i]*ymMaxY[i] ; | |
697 | } | |
698 | Double_t cN = (sx2y*kN-sy*sx2)*(sx3*sx-sx2*sx2)-(sx2y*sx-sxy*sx2)*(sx3*kN-sx*sx2); | |
699 | Double_t cD = (sx4*kN-sx2*sx2)*(sx3*sx-sx2*sx2)-(sx4*sx-sx3*sx2)*(sx3*kN-sx*sx2) ; | |
700 | Double_t c = cN / cD ; | |
701 | Double_t b = ((sx2y*kN-sy*sx2)-c*(sx4*kN-sx2*sx2))/(sx3*kN-sx*sx2) ; | |
702 | Double_t a = (sy-b*sx-c*sx2)/kN ; | |
703 | Double_t xmax = -b/(2*c) ; | |
704 | ymax = a + b*xmax + c*xmax*xmax ;//<========== This is the maximum amplitude | |
029fe7a2 | 705 | amp = ymax; |
16605c06 | 706 | } |
707 | } | |
708 | ||
709 | Double_t diff = TMath::Abs(1-ymMaxY[0]/amp) ; | |
710 | if (diff > 0.1) | |
711 | amp = ymMaxY[0] ; | |
712 | //printf("Yves : Amp %f, time %g\n",amp, time); | |
713 | //END YS | |
ee299369 | 714 | return; |
715 | } | |
16605c06 | 716 | |
ee299369 | 717 | //__________________________________________________________________ |
718 | Double_t AliEMCALRawUtils::RawResponseFunction(Double_t *x, Double_t *par) | |
719 | { | |
8cb998bd | 720 | // Matches version used in 2007 beam test |
721 | // | |
ee299369 | 722 | // Shape of the electronics raw reponse: |
723 | // It is a semi-gaussian, 2nd order Gamma function of the general form | |
724 | // | |
7643e728 | 725 | // xx = (t - t0 + tau) / tau [xx is just a convenient help variable] |
726 | // F = A * (xx**N * exp( N * ( 1 - xx) ) for xx >= 0 | |
727 | // F = 0 for xx < 0 | |
ee299369 | 728 | // |
729 | // parameters: | |
8cb998bd | 730 | // A: par[0] // Amplitude = peak value |
731 | // t0: par[1] | |
732 | // tau: par[2] | |
733 | // N: par[3] | |
734 | // ped: par[4] | |
ee299369 | 735 | // |
736 | Double_t signal ; | |
8cb998bd | 737 | Double_t tau =par[2]; |
e5bbbc4e | 738 | Double_t n =par[3]; |
8cb998bd | 739 | Double_t ped = par[4]; |
740 | Double_t xx = ( x[0] - par[1] + tau ) / tau ; | |
ee299369 | 741 | |
5a056daa | 742 | if (xx <= 0) |
8cb998bd | 743 | signal = ped ; |
ee299369 | 744 | else { |
e5bbbc4e | 745 | signal = ped + par[0] * TMath::Power(xx , n) * TMath::Exp(n * (1 - xx )) ; |
ee299369 | 746 | } |
747 | return signal ; | |
748 | } | |
749 | ||
7683df1d | 750 | //__________________________________________________________________ |
751 | Double_t AliEMCALRawUtils::RawResponseFunctionLog(Double_t *x, Double_t *par) | |
752 | { | |
753 | // Matches version used in 2007 beam test | |
754 | // | |
755 | // Shape of the electronics raw reponse: | |
756 | // It is a semi-gaussian, 2nd order Gamma function of the general form | |
757 | // | |
758 | // xx = (t - t0 + tau) / tau [xx is just a convenient help variable] | |
759 | // F = A * (xx**N * exp( N * ( 1 - xx) ) for xx >= 0 | |
760 | // F = 0 for xx < 0 | |
761 | // | |
762 | // parameters: | |
763 | // Log[A]: par[0] // Amplitude = peak value | |
764 | // t0: par[1] | |
765 | // tau: par[2] | |
766 | // N: par[3] | |
767 | // ped: par[4] | |
768 | // | |
769 | Double_t signal ; | |
770 | Double_t tau =par[2]; | |
771 | Double_t n =par[3]; | |
772 | //Double_t ped = par[4]; // not used | |
773 | Double_t xx = ( x[0] - par[1] + tau ) / tau ; | |
774 | ||
775 | if (xx < 0) | |
776 | signal = par[0] - n*TMath::Log(TMath::Abs(xx)) + n * (1 - xx ) ; | |
777 | else { | |
778 | signal = par[0] + n*TMath::Log(xx) + n * (1 - xx ) ; | |
779 | } | |
780 | return signal ; | |
781 | } | |
782 | ||
ee299369 | 783 | //__________________________________________________________________ |
6751f762 | 784 | Bool_t AliEMCALRawUtils::RawSampledResponse(const Double_t dtime, const Double_t damp, Int_t * adcH, Int_t * adcL, const Int_t keyErr) const |
ee299369 | 785 | { |
786 | // for a start time dtime and an amplitude damp given by digit, | |
787 | // calculates the raw sampled response AliEMCAL::RawResponseFunction | |
788 | ||
ee299369 | 789 | Bool_t lowGain = kFALSE ; |
790 | ||
48a56166 | 791 | // A: par[0] // Amplitude = peak value |
792 | // t0: par[1] | |
793 | // tau: par[2] | |
794 | // N: par[3] | |
795 | // ped: par[4] | |
796 | ||
56e13066 | 797 | TF1 signalF("signal", RawResponseFunction, 0, GetRawFormatTimeBins(), 5); |
48a56166 | 798 | signalF.SetParameter(0, damp) ; |
56e13066 | 799 | signalF.SetParameter(1, (dtime + fgTimeTrigger)/fgTimeBinWidth) ; |
b4133f05 | 800 | signalF.SetParameter(2, fTau) ; |
801 | signalF.SetParameter(3, fOrder); | |
fe93d365 | 802 | signalF.SetParameter(4, fgPedestalValue); |
6751f762 | 803 | |
804 | Double_t signal=0.0, noise=0.0; | |
ee299369 | 805 | for (Int_t iTime = 0; iTime < GetRawFormatTimeBins(); iTime++) { |
916f1e76 | 806 | signal = signalF.Eval(iTime) ; |
807 | ||
7643e728 | 808 | // Next lines commeted for the moment but in principle it is not necessary to add |
ff10f540 | 809 | // extra noise since noise already added at the digits level. |
7643e728 | 810 | |
fe93d365 | 811 | //According to Terry Awes, 13-Apr-2008 |
812 | //add gaussian noise in quadrature to each sample | |
09974781 | 813 | //Double_t noise = gRandom->Gaus(0.,fgFEENoise); |
fe93d365 | 814 | //signal = sqrt(signal*signal + noise*noise); |
815 | ||
e2c2134b | 816 | // March 17,09 for fast fit simulations by Alexei Pavlinov. |
4fe71e02 | 817 | // Get from PHOS analysis. In some sense it is open questions. |
6751f762 | 818 | if(keyErr>0) { |
819 | noise = gRandom->Gaus(0.,fgFEENoise); | |
820 | signal += noise; | |
821 | } | |
822 | ||
ee299369 | 823 | adcH[iTime] = static_cast<Int_t>(signal + 0.5) ; |
b4133f05 | 824 | if ( adcH[iTime] > fgkRawSignalOverflow ){ // larger than 10 bits |
825 | adcH[iTime] = fgkRawSignalOverflow ; | |
ee299369 | 826 | lowGain = kTRUE ; |
827 | } | |
828 | ||
829 | signal /= fHighLowGainFactor; | |
830 | ||
831 | adcL[iTime] = static_cast<Int_t>(signal + 0.5) ; | |
b4133f05 | 832 | if ( adcL[iTime] > fgkRawSignalOverflow) // larger than 10 bits |
833 | adcL[iTime] = fgkRawSignalOverflow ; | |
ee299369 | 834 | } |
835 | return lowGain ; | |
836 | } | |
4fe71e02 | 837 | |
838 | //__________________________________________________________________ | |
839 | void AliEMCALRawUtils::SetFittingAlgorithm(Int_t fitAlgo) | |
840 | { | |
841 | //Set fitting algorithm and initialize it if this same algorithm was not set before. | |
916f1e76 | 842 | //printf("**** Set Algorithm , number %d ****\n",fitAlgo); |
843 | ||
4fe71e02 | 844 | if(fitAlgo == fFittingAlgorithm && fRawAnalyzer) { |
845 | //Do nothing, this same algorithm already set before. | |
846 | //printf("**** Algorithm already set before, number %d, %s ****\n",fitAlgo, fRawAnalyzer->GetName()); | |
847 | return; | |
848 | } | |
849 | //Initialize the requested algorithm | |
850 | if(fitAlgo != fFittingAlgorithm || !fRawAnalyzer) { | |
851 | //printf("**** Init Algorithm , number %d ****\n",fitAlgo); | |
852 | ||
853 | fFittingAlgorithm = fitAlgo; | |
854 | if (fRawAnalyzer) delete fRawAnalyzer; // delete prev. analyzer if existed. | |
855 | ||
856 | if (fitAlgo == kFastFit) { | |
857 | fRawAnalyzer = new AliCaloRawAnalyzerFastFit(); | |
858 | } | |
859 | else if (fitAlgo == kNeuralNet) { | |
860 | fRawAnalyzer = new AliCaloRawAnalyzerNN(); | |
861 | } | |
862 | else if (fitAlgo == kLMS) { | |
863 | fRawAnalyzer = new AliCaloRawAnalyzerLMS(); | |
864 | } | |
865 | else if (fitAlgo == kPeakFinder) { | |
866 | fRawAnalyzer = new AliCaloRawAnalyzerPeakFinder(); | |
867 | } | |
868 | else if (fitAlgo == kCrude) { | |
869 | fRawAnalyzer = new AliCaloRawAnalyzerCrude(); | |
870 | } | |
871 | else { | |
872 | fRawAnalyzer = new AliCaloRawAnalyzer(); | |
873 | } | |
874 | } | |
875 | ||
876 | } | |
877 | ||
878 |