110 ATH_MSG_DEBUG(
"1) LArDigitContainer container size = " << digitContainer->size() );
112 ATH_MSG_DEBUG(
"2) LArDigitContainer container size = " << digitContainer->size() );
113 if( digitContainer->empty() ) {
115 return StatusCode::SUCCESS;
118 auto larRawChannelContainer = std::make_unique<LArRawChannelContainer>();
123 if (
detStore()->retrieve(larPedestal).isFailure()) {
125 ATH_MSG_DEBUG(
"No pedestal found in database. Use default values." );
133 return StatusCode::FAILURE;
148 int nMinEM(0),nMinHEC(0),nMinFCAL(0);
149 std::vector<float> fSumEM,fSumHEC,fSumFCAL;
150 std::vector<float> *pSum =
nullptr;
156 for(
int iloop=0;iloop<2;iloop++) {
161 for (
const LArDigit* digit : *digitContainer) {
164 const std::vector<short>& samples = digit->samples();
165 unsigned int nSamples = samples.size();
170 float thePedestal=-1;
172 float DBpedestal =larPedestal->
pedestal(chid,gain);
174 thePedestal=DBpedestal;
178 ATH_MSG_DEBUG(
"No pedestal found for this cell. Use default value " << thePedestal );
181 std::vector<float> mySamples;
182 mySamples.resize(samples.size());
187 for (
unsigned int i=0; i<samples.size(); i++ )
189 mySamples[i] = ((float)samples[i]) - thePedestal;
195 GainFactor = 9.8*9.8;
202 ATH_MSG_ERROR(
"Channel " << chid <<
"unknown gain: " << gain );
208 ATH_MSG_DEBUG( std::hex <<
" FebID / chid / channel = " << FebID <<
" / " << chid <<
" / " << channel );
216 unsigned int nAverage = 1;
219 const Identifier id = cabling->cnvToIdentifier(chid);
227 else if (
m_emId->is_lar_em(
id)) {
234 else if (
m_hecId->is_lar_hec(
id)) {
256 if ( nAverage > nSamples ) {
258 << nAverage <<
") is larger than total number of samples ("
259 << nSamples <<
")! adjusting nAverage ... " );
268 pSum->resize(nSamples,0);
270 for(
unsigned i=0;i<nSamples-nAverage+1;i++ ) {
271 for(
unsigned j=0;j<nAverage;j++ ) {
272 (*pSum)[i] += mySamples[i+j];
277 float maxADCPeak = 0.;
278 unsigned int iPeakSamp = 0;
279 float averagedADC = 0;
280 for(
unsigned i=0;i<nSamples;i++ ) {
281 if ( maxADCPeak < mySamples[i] )
283 maxADCPeak = mySamples[i];
286 if ( (
int)i >= nMin && i < nMin+nAverage )
287 averagedADC += mySamples[i];
289 averagedADC /= myScale;
291 bool CubicFailed =
false;
299 ADCPeak = maxADCPeak;
300 time = ((float)iPeakSamp) * 25.0 * nanosecond;
308 if(isFCAL &&
m_mode ==
"CUBIC" &&
315 const float invT[3][3]
321 if ( iPeakSamp <= 1 ) {
323 }
else if ( iPeakSamp >= nSamples - 1 ) {
331 float A[3] = {0, 0, 0};
334 for (
int ia = 0; ia < 3; ia++)
335 for (
int it = 0; it < 3; it++)
336 A[ia] += invT[ia][it] * mySamples[it0+it];
339 if ( not ( CubicFailed = (
A[2] == 0 ) ) ) {
340 dtmax = -1.0 *
A[1] / 2.0 /
A[2];
341 if ( ! ( CubicFailed = ( dtmax < 0 || dtmax > 2 ) ) ) {
344 for(
int ia = 0; ia < 3; ia++)
345 ADCPeak +=
A[ia] * pow(dtmax, ia);
351 float weightSum = 0.;
353 for(
int it=0; it<3; it++) {
354 timeSum += float(mySamples[it0+it]) * float(it);
355 weightSum += float(mySamples[it0+it]);
357 time = (float(it0) + timeSum / weightSum) * 25.0 * nanosecond;
364 }
else if(not isFCAL ||
m_mode !=
"CUBIC" ) {
367 if (
m_mode ==
"PARABOLA" &&
372 const Identifier id = cabling->cnvToIdentifier(chid);
374 if (
m_emId->is_em_barrel(
id)) {
375 layer=
m_emId->sampling(
id);
380 <<
" Cannot get offline identifier from online ID = " << chid );
385 ADCPeak = peak[0]-thePedestal;
386 if (peak.size()==2)time = peak[1];
389 ATH_MSG_DEBUG(
"No pic is computed from Parabola. Use scaled average of selected samples" );
390 ADCPeak = averagedADC;
394 ATH_MSG_FATAL(
"No parabola tool available ! Choose another mode" );
399 else if (
m_mode ==
"CUBIC" &&
403 const float invT[4][4]
405 { -1.83333, 3, -1.5, 0.333333},
407 {-0.166666, 0.5, -0.5, 0.166666} };
410 if ( iPeakSamp <= 1 ) {
412 }
else if ( iPeakSamp >= nSamples - 2 ) {
415 it0 = ( mySamples[iPeakSamp-2] > mySamples[iPeakSamp+2] )
421 float A[4] = {0, 0, 0, 0};
425 for (
int ia = 0; ia < 4; ia++)
426 for (
int it = 0; it < 4; it++)
427 A[ia] += invT[ia][it] * mySamples[it0+it];
430 disc =
A[2]*
A[2] - 3*
A[1]*
A[3];
431 if ( ! ( CubicFailed = ( disc < 0 ||
A[3] == 0 ) ) ) {
432 dtmax = (-
A[2]-std::sqrt(disc))/(
A[3]*3);
433 if ( ! ( CubicFailed = ( dtmax < 0 || dtmax > 3 ) ) ) {
434 time = (float(it0) + dtmax) * 25.0 * nanosecond;
435 for(
int ia = 0; ia < 4; ia++)
436 ADCPeak +=
A[ia] * pow(dtmax, ia);
444 if(
m_mode ==
"FIXED" || CubicFailed ) {
445 ADCPeak = averagedADC;
452 << (isEM?
"EM":(isHEC?
"HEC":(isFCAL?
"FCAL":
"none")))
465 float ADCPeakPower=ADCPeak;
469 if (ramp.size()>1 && ramp[1]<500 && ramp[1]>0) {
471 for (
unsigned i=1;i<ramp.size();i++)
472 {energy+=ramp[i]*ADCPeakPower;
474 ADCPeakPower*=ADCPeak;
479 ATH_MSG_DEBUG(
"No Ramp found for this cell. Use default values" );
482 float ADCtoMeV = 10.0;
490 const Identifier id = cabling->cnvToIdentifier(chid);
494 if (
m_emId->is_em_barrel(
id)) {
495 const int layer=
m_emId->sampling(
id);
498 if (layer==2 &&
eta<32)
499 ADCtoMeV=ADCtoMeV*(12./18.);
500 ATH_MSG_DEBUG(
" in EMB s="<< layer <<
", using ADCtoMeV = " << ADCtoMeV );
501 }
else if (
m_emId->is_em_endcap_inner(
id)) {
504 ATH_MSG_DEBUG(
" in EMEC inner s="<<
m_emId->sampling(
id)<<
", using ADCtoMeV = " << ADCtoMeV );
505 }
else if (
m_emId->is_em_endcap_outer(
id)) {
509 ATH_MSG_DEBUG(
" in EMEC outer s="<<
m_emId->sampling(
id)<<
", using ADCtoMeV = " << ADCtoMeV );
510 }
else if (
m_fcalId->is_lar_fcal(
id)) {
514 }
else if (
m_hecId->is_lar_hec(
id)) {
526 energy = ADCPeak * ADCtoMeV * GainFactor;
529 uint16_t iprovenance=0;
533 LArRawChannel larRawChannel(chid,(
int)energy,(
int)time,iquality,iprovenance, gain);
534 larRawChannelContainer->add(larRawChannel);
542 if ( !fSumEM.empty() ) {
543 for(i=0;i<fSumEM.size();i++) {
544 if (i == 0 || fSumEM[i] > tmpSum ) {
549 ATH_MSG_DEBUG(
"Found best EM window starting at sample <" << nMinEM <<
">" );
552 for(i=0;i<fSumHEC.size();i++) {
553 if (i == 0 || fSumHEC[i] > tmpSum ) {
558 if ( !fSumHEC.empty() ) {
559 ATH_MSG_DEBUG(
"Found best HEC window starting at sample <" << nMinHEC <<
">" );
562 for(i=0;i<fSumFCAL.size();i++) {
563 if (i == 0 || fSumFCAL[i] > tmpSum ) {
565 tmpSum = fSumFCAL[i];
568 if ( !fSumFCAL.empty() ) {
569 ATH_MSG_DEBUG(
"Found best FCAL window starting at sample <" << nMinFCAL <<
">" );
580 return StatusCode::SUCCESS;