393{
394 const EventContext& ctx = Gaudi::Hive::currentContext();
395 SG::ReadCondHandle<LArOnOffIdMapping> cablingHdl (
m_cablingKey, ctx);
396 const LArOnOffIdMapping*
cabling=(*cablingHdl);
397
399
400 const ILArPedestal* pedestal = nullptr;
401 const ILArNoise*
noise =
nullptr;
403 SG::ReadCondHandle<ILArNoise> noiseH (
m_noiseKey, ctx);
404 noise = noiseH.cptr();
405 }
406 else {
408 pedestal = pedH.cptr();
409 }
410
411 const CaloNoise* elecNoise = nullptr;
414 elecNoise = noiseH.cptr();
415 }
416
417 SG::ReadCondHandle<LArADC2MeV> adc2mev (
m_adc2mevKey, ctx);
418
419 FILE*
fp = fopen(
"calonoise.txt",
"w");
420
421 TBranch* b1 =
m_tree->GetBranch(
"luminosity");
422 TBranch* b2 =
m_tree->GetBranch(
"nevt");
423 TBranch* b3 =
m_tree->GetBranch(
"average");
424 TBranch* b4 =
m_tree->GetBranch(
"rms");
425 TBranch* b5 =
m_tree->GetBranch(
"nevt_good");
432 ATH_MSG_DEBUG (
" Number of entries in ntuple " << nentries );
433
434 std::vector<float> anoise;
435 std::vector<float> bnoise;
438
440 std::vector<float>
x;
441 std::vector<float>
y;
442 std::vector<float> ey;
443
445 b1->GetEntry(i);
446 b2->GetEntry(i);
447 b4->GetEntry(i);
448
453 }
454 }
455
457 if (
x.size()==1) anoise[
icell]=
y[0];
458 }
459 else {
460 HepMatrix
alpha(2,2);
462
463 for (
unsigned int i=0;
i<2;
i++) {
464 for (unsigned int j=0;j<2;j++) {
466 for (
unsigned int k=0;
k<
x.size();
k++) {
467 alpha[
i][j] += ((std::pow(
x[k],i))*(std::pow(
x[k],j))/(std::pow(ey[k],2)));
468 }
469 }
470 }
471 for (
unsigned int i=0;
i<2;
i++) {
473 for (
unsigned int k=0;
k<
x.size();
k++) {
474 beta[
i] += (
y[
k]*(std::pow(
x[k],i))/(std::pow(ey[k],2)));
475 }
476 }
477 HepVector
comp=solve(alpha,beta);
480
481 }
483 msg() << MSG::DEBUG <<
" cell " <<
icell <<
" lumi/noise ";
484 for (
unsigned int i=0;
i<
x.size();
i++) {
485 msg() << MSG::DEBUG <<
x[
i] <<
" " <<
y[
i] <<
" / ";
486 }
488 ATH_MSG_DEBUG (
" fitted a,b " << anoise[icell] <<
" " << bnoise[icell] );
489 }
490 }
491
492
493
495 if (anoise[icell]<3. ||
m_treeData->m_nevt_good[icell]==0) {
496 IdentifierHash idHash =
icell;
497 Identifier
id=
m_calo_id->cell_id(idHash);
500 Identifier regionId =
m_calo_id->region_id(
id);
502 int phimin =
m_calo_id->phi_min(regionId);
503 int phimax =
m_calo_id->phi_max(regionId);
504 int nring=0;
506 ATH_MSG_DEBUG (
" regionId,eta,phimin,phimax " << regionId <<
" " <<
eta <<
" " << phimin <<
" " << phimax );
511 msg() << MSG::DEBUG <<
" cell in ring " <<
m_calo_id->show_to_string(
id2) ;
512 IdentifierHash idHash2 =
m_calo_id->calo_cell_hash(
id2);
514 if (index>=0 && index<
m_ncell) {
516 if (anoise[index]>3. &&
m_treeData->m_nevt_good[index]>0) {
517 nring+=1;
519 }
520 }
521 }
522 }
523 if (nring>0) {
524 float patched_noise =
sum/((
float)(nring));
525 if (patched_noise>anoise[icell]) anoise[
icell] = patched_noise;
526 }
527 ATH_MSG_DEBUG(
" corrected noise nring, anoise[icell] " << nring <<
" " << anoise[icell] );
528 }
529 }
530 }
531
532
533
535 IdentifierHash idHash =
icell;
536 Identifier
id=
m_calo_id->cell_id(idHash);
537 HWIdentifier hwid=
cabling->createSignalChannelID(
id);
538 int subCalo;
539 IdentifierHash idSubHash =
m_calo_id->subcalo_cell_hash (idHash, subCalo);
540
541 int iCool=-1;
545 iCool=2;
546 else
547 iCool=1;
548 }
551 iCool=3;
552 else
553 iCool=0;
554 }
555
556 }
558 iCool=16;
559 }
561 iCool=32;
562 }
564 iCool=48;
565 }
566 int ii = (
int) (idSubHash);
567
568
569
576
577
578
581
582 float anoise_corr=anoise[
icell];
583
584 if (gain != gainref) {
585
586
587 float noise0=-1.;
589 else {
592 }
593 LArVectorProxy
594 polynom_adc2mev0 = adc2mev->ADC2MEV(id,gainref);
595 float adc2mev0=-1;
596 if (polynom_adc2mev0.size()>1) adc2mev0=polynom_adc2mev0[1];
597
598
599
600 float noise1=-1;
602 else {
605 }
606 LArVectorProxy
607 polynom_adc2mev1 = adc2mev->ADC2MEV(hwid,gain);
608 float adc2mev1=-1;
609 if (polynom_adc2mev1.size()>1) adc2mev1=polynom_adc2mev1[1];
610
611
612
613
614 if (noise0>0 && noise1>0 && adc2mev0>0 && adc2mev1>0.) {
615 anoise_corr = anoise[
icell]*noise1/noise0 * adc2mev1/adc2mev0;
616 }
617 }
618
619
620
621
622 if (elecNoise) {
623 float adb = elecNoise->
getNoise(
id,gain);
624
625
626 if (anoise_corr<1. && adb>1e-6) {
627 anoise_corr = adb;
628 ATH_MSG_WARNING (
" No noise found for cell: " <<
m_calo_id->show_to_string(
id) <<
" gain: " << gain <<
" Use reference value: " << adb );
629 }
630 if (adb>1e-6) {
631 float delta = std::fabs((anoise_corr-adb)/adb);
634 << " computed " << anoise_corr << " reference " << adb );
635 }
636 }
637
638
639 fprintf(fp,"%5d %5d %5d %8.3f %8.3f\n",iCool,ii,gain,anoise_corr,bnoise[icell]);
640 }
641 }
642
643 else {
650 float anoise_corr = anoise[
icell];
651 if (igain>0) {
652 }
653 fprintf(fp,"%5d %5d %5d %8.3f %8.3f\n",iCool,ii,gain,anoise_corr,bnoise[icell]);
654 }
655 }
656
657 }
658
659 fclose(fp);
660 return StatusCode::SUCCESS;
661}
Scalar eta() const
pseudorapidity method
Scalar phi() const
phi method
bool msgLvl(const MSG::Level lvl) const
SG::ReadCondHandleKey< ILArNoise > m_noiseKey
SG::ReadCondHandleKey< CaloNoise > m_elecNoiseKey
SG::ReadCondHandleKey< ILArPedestal > m_pedestalKey
SG::ReadCondHandleKey< LArOnOffIdMapping > m_cablingKey
SG::ReadCondHandleKey< LArADC2MeV > m_adc2mevKey
virtual float pedestalRMS(const HWIdentifier &id, int gain) const =0
access to RMS of Pedestal index by Identifier, and gain setting
bool is_valid() const
Check if id is in a valid state.