366{
367 const EventContext& ctx = Gaudi::Hive::currentContext();
368 SG::ReadCondHandle<LArOnOffIdMapping> cablingHdl (
m_cablingKey, ctx);
369 const LArOnOffIdMapping*
cabling=(*cablingHdl);
370
372
373 const ILArPedestal* pedestal = nullptr;
374 const ILArNoise*
noise =
nullptr;
376 SG::ReadCondHandle<ILArNoise> noiseH (
m_noiseKey, ctx);
377 noise = noiseH.cptr();
378 }
379 else {
381 pedestal = pedH.cptr();
382 }
383
384 const CaloNoise* elecNoise = nullptr;
387 elecNoise = noiseH.cptr();
388 }
389
390 SG::ReadCondHandle<LArADC2MeV> adc2mev (
m_adc2mevKey, ctx);
391
392 FILE*
fp = fopen(
"calonoise.txt",
"w");
393
394 TBranch* b1 =
m_tree->GetBranch(
"luminosity");
395 TBranch* b2 =
m_tree->GetBranch(
"nevt");
396 TBranch* b3 =
m_tree->GetBranch(
"average");
397 TBranch* b4 =
m_tree->GetBranch(
"rms");
398 TBranch* b5 =
m_tree->GetBranch(
"nevt_good");
405 ATH_MSG_DEBUG (
" Number of entries in ntuple " << nentries );
406
407 std::vector<float> anoise;
408 std::vector<float> bnoise;
411
413 std::vector<float>
x;
414 std::vector<float>
y;
415 std::vector<float> ey;
416
418 b1->GetEntry(i);
419 b2->GetEntry(i);
420 b4->GetEntry(i);
421
426 }
427 }
428
430 if (
x.size()==1) anoise[
icell]=
y[0];
431 }
432 else {
433 HepMatrix
alpha(2,2);
435
436 for (
unsigned int i=0;
i<2;
i++) {
437 for (
unsigned int j=0;
j<2;
j++) {
439 for (
unsigned int k=0;
k<
x.size();
k++) {
440 alpha[
i][
j] += ((std::pow(
x[k],i))*(std::pow(
x[k],j))/(std::pow(ey[k],2)));
441 }
442 }
443 }
444 for (
unsigned int i=0;
i<2;
i++) {
446 for (
unsigned int k=0;
k<
x.size();
k++) {
447 beta[
i] += (
y[
k]*(std::pow(
x[k],i))/(std::pow(ey[k],2)));
448 }
449 }
450 HepVector
comp=solve(alpha,beta);
453
454 }
456 msg() << MSG::DEBUG <<
" cell " <<
icell <<
" lumi/noise ";
457 for (
unsigned int i=0;
i<
x.size();
i++) {
458 msg() << MSG::DEBUG <<
x[
i] <<
" " <<
y[
i] <<
" / ";
459 }
461 ATH_MSG_DEBUG (
" fitted a,b " << anoise[icell] <<
" " << bnoise[icell] );
462 }
463 }
464
465
466
468 if (anoise[icell]<3. ||
m_treeData->m_nevt_good[icell]==0) {
469 IdentifierHash idHash =
icell;
470 Identifier
id=
m_calo_id->cell_id(idHash);
473 Identifier regionId =
m_calo_id->region_id(
id);
475 int phimin =
m_calo_id->phi_min(regionId);
476 int phimax =
m_calo_id->phi_max(regionId);
477 int nring=0;
479 ATH_MSG_DEBUG (
" regionId,eta,phimin,phimax " << regionId <<
" " <<
eta <<
" " << phimin <<
" " << phimax );
484 msg() << MSG::DEBUG <<
" cell in ring " <<
m_calo_id->show_to_string(id2) ;
485 IdentifierHash idHash2 =
m_calo_id->calo_cell_hash(id2);
487 if (index>=0 && index<
m_ncell) {
489 if (anoise[index]>3. &&
m_treeData->m_nevt_good[index]>0) {
490 nring+=1;
492 }
493 }
494 }
495 }
496 if (nring>0) {
497 float patched_noise =
sum/((
float)(nring));
498 if (patched_noise>anoise[icell]) anoise[
icell] = patched_noise;
499 }
500 ATH_MSG_DEBUG(
" corrected noise nring, anoise[icell] " << nring <<
" " << anoise[icell] );
501 }
502 }
503 }
504
505
506
508 IdentifierHash idHash =
icell;
509 Identifier
id=
m_calo_id->cell_id(idHash);
510 HWIdentifier hwid=
cabling->createSignalChannelID(
id);
511 int subCalo;
512 IdentifierHash idSubHash =
m_calo_id->subcalo_cell_hash (idHash, subCalo);
513
514 int iCool=-1;
518 iCool=2;
519 else
520 iCool=1;
521 }
524 iCool=3;
525 else
526 iCool=0;
527 }
528
529 }
531 iCool=16;
532 }
534 iCool=32;
535 }
537 iCool=48;
538 }
539 int ii = (
int) (idSubHash);
540
541
542
549
550
551
554
555 float anoise_corr=anoise[
icell];
556
557 if (gain != gainref) {
558
559
560 float noise0=-1.;
562 else {
565 }
566 LArVectorProxy
567 polynom_adc2mev0 = adc2mev->ADC2MEV(id,gainref);
568 float adc2mev0=-1;
569 if (polynom_adc2mev0.size()>1) adc2mev0=polynom_adc2mev0[1];
570
571
572
573 float noise1=-1;
575 else {
578 }
579 LArVectorProxy
580 polynom_adc2mev1 = adc2mev->ADC2MEV(hwid,gain);
581 float adc2mev1=-1;
582 if (polynom_adc2mev1.size()>1) adc2mev1=polynom_adc2mev1[1];
583
584
585
586
587 if (noise0>0 && noise1>0 && adc2mev0>0 && adc2mev1>0.) {
588 anoise_corr = anoise[
icell]*noise1/noise0 * adc2mev1/adc2mev0;
589 }
590 }
591
592
593
594
595 if (elecNoise) {
596 float adb = elecNoise->
getNoise(
id,gain);
597
598
599 if (anoise_corr<1. && adb>1e-6) {
600 anoise_corr = adb;
601 ATH_MSG_WARNING (
" No noise found for cell: " <<
m_calo_id->show_to_string(
id) <<
" gain: " << gain <<
" Use reference value: " << adb );
602 }
603 if (adb>1e-6) {
604 float delta = std::fabs((anoise_corr-adb)/adb);
607 << " computed " << anoise_corr << " reference " << adb );
608 }
609 }
610
611
612 fprintf(fp,"%5d %5d %5d %8.3f %8.3f\n",iCool,ii,gain,anoise_corr,bnoise[icell]);
613 }
614 }
615
616 else {
623 float anoise_corr = anoise[
icell];
624 if (igain>0) {
625 }
626 fprintf(fp,"%5d %5d %5d %8.3f %8.3f\n",iCool,ii,gain,anoise_corr,bnoise[icell]);
627 }
628 }
629
630 }
631
632 fclose(fp);
633 return StatusCode::SUCCESS;
634}
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
Gaudi::Property< bool > m_doMC
SG::ReadCondHandleKey< LArOnOffIdMapping > m_cablingKey
SG::ReadCondHandleKey< LArADC2MeV > m_adc2mevKey
Gaudi::Property< int > m_nmin
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.
float j(const xAOD::IParticle &, const xAOD::TrackMeasurementValidation &hit, const Eigen::Matrix3d &jab_inv)