From f28b567a8dc5d4c45ee13ae2c9f5bf8a2fcec44c Mon Sep 17 00:00:00 2001 From: Taeun_Kwon Date: Wed, 5 Aug 2020 11:23:42 -0400 Subject: [PATCH 1/9] adding noise correlation between adjacent time slices to noise covariance matrix --- .../HcalRecHit/interface/HBHEChannelInfo.h | 7 ++++++ DataFormats/HcalRecHit/src/classes_def.xml | 3 ++- DataFormats/HcalRecHit/test/test_hcal_reco.cu | 6 +++-- .../HcalRecAlgos/interface/MahiFit.h | 4 ++-- RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc | 22 +++++++++++++------ .../src/HBHEPhase1Reconstructor.cc | 2 ++ 6 files changed, 32 insertions(+), 12 deletions(-) diff --git a/DataFormats/HcalRecHit/interface/HBHEChannelInfo.h b/DataFormats/HcalRecHit/interface/HBHEChannelInfo.h index 7090944c580bd..9df1fc629f10a 100644 --- a/DataFormats/HcalRecHit/interface/HBHEChannelInfo.h +++ b/DataFormats/HcalRecHit/interface/HBHEChannelInfo.h @@ -28,6 +28,7 @@ class HBHEChannelInfo { darkCurrent_{0}, fcByPE_{0}, lambda_{0}, + noisecorr_{0}, riseTime_{0.f}, adc_{0}, dFcPerADC_{0.f}, @@ -51,6 +52,7 @@ class HBHEChannelInfo { darkCurrent_{0}, fcByPE_{0}, lambda_{0}, + noisecorr_{0}, riseTime_{0.f}, adc_{0}, dFcPerADC_{0.f}, @@ -73,6 +75,7 @@ class HBHEChannelInfo { darkCurrent_ = 0; fcByPE_ = 0; lambda_ = 0, dropped_ = true; + noisecorr_ = 0; hasLinkError_ = false; hasCapidError_ = false; } @@ -85,6 +88,7 @@ class HBHEChannelInfo { const double darkCurrent, const double fcByPE, const double lambda, + const double noisecorr, const bool linkError, const bool capidError, const bool dropThisChannel) { @@ -96,6 +100,7 @@ class HBHEChannelInfo { darkCurrent_ = darkCurrent; fcByPE_ = fcByPE; lambda_ = lambda, dropped_ = dropThisChannel; + noisecorr_ = noisecorr; hasLinkError_ = linkError; hasCapidError_ = capidError; } @@ -136,6 +141,7 @@ class HBHEChannelInfo { constexpr double darkCurrent() const { return darkCurrent_; } constexpr double fcByPE() const { return fcByPE_; } constexpr double lambda() const { return lambda_; } + constexpr double noisecorr() const { return noisecorr_; } constexpr bool isDropped() const { return dropped_; } constexpr bool hasLinkError() const { return hasLinkError_; } constexpr bool hasCapidError() const { return hasCapidError_; } @@ -254,6 +260,7 @@ class HBHEChannelInfo { double darkCurrent_; double fcByPE_; double lambda_; + double noisecorr_; // Signal rise time from TDC in ns (if provided) float riseTime_[MAXSAMPLES]; diff --git a/DataFormats/HcalRecHit/src/classes_def.xml b/DataFormats/HcalRecHit/src/classes_def.xml index 1288f7d8eec22..c6d35a85cd6ec 100644 --- a/DataFormats/HcalRecHit/src/classes_def.xml +++ b/DataFormats/HcalRecHit/src/classes_def.xml @@ -21,7 +21,8 @@ - + + diff --git a/DataFormats/HcalRecHit/test/test_hcal_reco.cu b/DataFormats/HcalRecHit/test/test_hcal_reco.cu index 5f5d39fc562a5..5557e0a203b02 100644 --- a/DataFormats/HcalRecHit/test/test_hcal_reco.cu +++ b/DataFormats/HcalRecHit/test/test_hcal_reco.cu @@ -23,7 +23,7 @@ __global__ void kernel_test_hcal_hfqie10info() { HFQIE10Info info; } __global__ void kernel_test_hcal_hbhechinfo(HBHEChannelInfo *other) { HBHEChannelInfo info{true, true}; - info.setChannelInfo(HcalDetId{0}, 10, 10, 10, 1, 2.0, 2.0, 2.0, false, false, false); + info.setChannelInfo(HcalDetId{0}, 10, 10, 10, 1, 2.0, 2.0, 2.0, 0.0, false, false, false); other->setChannelInfo(info.id(), info.recoShape(), info.nSamples(), @@ -32,6 +32,7 @@ __global__ void kernel_test_hcal_hbhechinfo(HBHEChannelInfo *other) { info.darkCurrent(), info.fcByPE(), info.lambda(), + info.noisecorr(), info.hasLinkError(), info.hasCapidError(), info.isDropped()); @@ -85,7 +86,7 @@ void test_hcal_hbhechinfo() { }; HBHEChannelInfo h_info, h_info_test{true, true}; - h_info_test.setChannelInfo(HcalDetId{0}, 10, 10, 10, 1, 2.0, 2.0, 2.0, false, false, false); + h_info_test.setChannelInfo(HcalDetId{0}, 10, 10, 10, 1, 2.0, 2.0, 2.0, 0.0, false, false, false); HBHEChannelInfo *d_info; cudaMalloc((void **)&d_info, sizeof(HBHEChannelInfo)); @@ -103,6 +104,7 @@ void test_hcal_hbhechinfo() { assert(h_info.darkCurrent() == h_info_test.darkCurrent()); assert(h_info.fcByPE() == h_info_test.fcByPE()); assert(h_info.lambda() == h_info_test.lambda()); + assert(h_info.noisecorr() == h_info_test.noisecorr()); assert(h_info.hasLinkError() == h_info_test.hasLinkError()); assert(h_info.hasCapidError() == h_info_test.hasCapidError()); diff --git a/RecoLocalCalo/HcalRecAlgos/interface/MahiFit.h b/RecoLocalCalo/HcalRecAlgos/interface/MahiFit.h index a176026b0fd03..6bd9aaa08cabc 100644 --- a/RecoLocalCalo/HcalRecAlgos/interface/MahiFit.h +++ b/RecoLocalCalo/HcalRecAlgos/interface/MahiFit.h @@ -116,7 +116,7 @@ class MahiFit { void phase1Debug(const HBHEChannelInfo& channelData, MahiDebugInfo& mdi) const; - void doFit(std::array& correctedOutput, const int nbx) const; + void doFit(std::array& correctedOutput, const int nbx, const double noisecorr) const; void setPulseShapeTemplate(int pulseShapeId, const HcalPulseShapes& ps, @@ -130,7 +130,7 @@ class MahiFit { private: typedef std::pair > ShapeWithId; - const float minimize() const; + const float minimize(const double noisecorr) const; void onePulseMinimize() const; void updateCov(const SampleMatrix& invCovMat) const; void resetPulseShapeTemplate(int pulseShapeId, const HcalPulseShapes& ps, unsigned int nSamples); diff --git a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc index 3cd7e5d733c20..c48445e34113c 100644 --- a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc +++ b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc @@ -78,10 +78,12 @@ void MahiFit::phase1Apply(const HBHEChannelInfo& channelData, if (iTS == nnlsWork_.tsOffset) tstrig += amplitude; } - + double noisecorr = channelData.noisecorr(); tsTOT *= channelData.tsGain(0); tstrig *= channelData.tsGain(0); + std::cout << noisecorr << std::endl; + useTriple = false; if (tstrig >= ts4Thresh_ && tsTOT > 0) { //Average pedestal width (for covariance matrix constraint) @@ -92,13 +94,13 @@ void MahiFit::phase1Apply(const HBHEChannelInfo& channelData, // only do pre-fit with 1 pulse if chiSq threshold is positive if (chiSqSwitch_ > 0) { - doFit(reconstructedVals, 1); + doFit(reconstructedVals, 1, noisecorr); if (reconstructedVals[2] > chiSqSwitch_) { - doFit(reconstructedVals, 0); //nbx=0 means use configured BXs + doFit(reconstructedVals, 0, noisecorr); //nbx=0 means use configured BXs useTriple = true; } } else { - doFit(reconstructedVals, 0); + doFit(reconstructedVals, 0, noisecorr); useTriple = true; } } else { @@ -112,7 +114,7 @@ void MahiFit::phase1Apply(const HBHEChannelInfo& channelData, chi2 = reconstructedVals[2]; } -void MahiFit::doFit(std::array& correctedOutput, int nbx) const { +void MahiFit::doFit(std::array& correctedOutput, int nbx, const double noisecorr) const { unsigned int bxSize = 1; if (nbx == 1) { @@ -175,7 +177,7 @@ void MahiFit::doFit(std::array& correctedOutput, int nbx) const { } } - const float chiSq = minimize(); + const float chiSq = minimize(noisecorr); bool foundintime = false; unsigned int ipulseintime = 0; @@ -203,7 +205,7 @@ void MahiFit::doFit(std::array& correctedOutput, int nbx) const { } } -const float MahiFit::minimize() const { +const float MahiFit::minimize(const double noisecorr) const { nnlsWork_.invcovp.setZero(nnlsWork_.tsSize, nnlsWork_.nPulseTot); nnlsWork_.ampVec.setZero(nnlsWork_.nPulseTot); @@ -211,6 +213,12 @@ const float MahiFit::minimize() const { invCovMat.setConstant(nnlsWork_.tsSize, nnlsWork_.tsSize, nnlsWork_.pedVal); invCovMat += nnlsWork_.noiseTerms.asDiagonal(); + //Add off-Diagonal components up to first order + for (unsigned int i = 1; i < nnlsWork_.tsSize; ++i) { + invCovMat(i - 1, i) += noisecorr * sqrt(nnlsWork_.noiseTerms.coeff(i - 1) * nnlsWork_.noiseTerms.coeff(i)); + invCovMat(i, i - 1) += noisecorr * sqrt(nnlsWork_.noiseTerms.coeff(i) * nnlsWork_.noiseTerms.coeff(i - 1)); + } + float oldChiSq = 9999; float chiSq = oldChiSq; diff --git a/RecoLocalCalo/HcalRecProducers/src/HBHEPhase1Reconstructor.cc b/RecoLocalCalo/HcalRecProducers/src/HBHEPhase1Reconstructor.cc index 5ebb60fce8c96..3074f7d401e49 100644 --- a/RecoLocalCalo/HcalRecProducers/src/HBHEPhase1Reconstructor.cc +++ b/RecoLocalCalo/HcalRecProducers/src/HBHEPhase1Reconstructor.cc @@ -492,6 +492,7 @@ void HBHEPhase1Reconstructor::processData(const Collection& coll, const double fcByPE = properties.siPMParameter->getFCByPE(); double darkCurrent = 0.; double lambda = 0.; + const double noisecorr = properties.siPMParameter->getauxi2(); if (!saveEffectivePeds || saveInfos_) { // needed for the dark current in the M2 in alternative of the effectivePed darkCurrent = properties.siPMParameter->getDarkCurrent(); @@ -571,6 +572,7 @@ void HBHEPhase1Reconstructor::processData(const Collection& coll, darkCurrent, fcByPE, lambda, + noisecorr, hwerr.first, hwerr.second, properties.taggedBadByDb || dropByZS); From 644125bb7d84f5154bb0e9614d6f546a075823cb Mon Sep 17 00:00:00 2001 From: TaeunKwon Date: Thu, 6 Aug 2020 00:32:36 +0900 Subject: [PATCH 2/9] Update HBHEChannelInfo.h correction on indentation --- DataFormats/HcalRecHit/interface/HBHEChannelInfo.h | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/DataFormats/HcalRecHit/interface/HBHEChannelInfo.h b/DataFormats/HcalRecHit/interface/HBHEChannelInfo.h index 9df1fc629f10a..3e84f08770532 100644 --- a/DataFormats/HcalRecHit/interface/HBHEChannelInfo.h +++ b/DataFormats/HcalRecHit/interface/HBHEChannelInfo.h @@ -28,7 +28,7 @@ class HBHEChannelInfo { darkCurrent_{0}, fcByPE_{0}, lambda_{0}, - noisecorr_{0}, + noisecorr_{0}, riseTime_{0.f}, adc_{0}, dFcPerADC_{0.f}, @@ -52,7 +52,7 @@ class HBHEChannelInfo { darkCurrent_{0}, fcByPE_{0}, lambda_{0}, - noisecorr_{0}, + noisecorr_{0}, riseTime_{0.f}, adc_{0}, dFcPerADC_{0.f}, From 38c4f8842d6b63fae5033ee671f7dc08a2ba4841 Mon Sep 17 00:00:00 2001 From: TaeunKwon Date: Thu, 6 Aug 2020 01:26:12 +0900 Subject: [PATCH 3/9] Update MahiFit.cc --- RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc index c48445e34113c..19ca09ca82184 100644 --- a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc +++ b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc @@ -78,7 +78,7 @@ void MahiFit::phase1Apply(const HBHEChannelInfo& channelData, if (iTS == nnlsWork_.tsOffset) tstrig += amplitude; } - double noisecorr = channelData.noisecorr(); + const double noisecorr = channelData.noisecorr(); tsTOT *= channelData.tsGain(0); tstrig *= channelData.tsGain(0); From 3bf261d4ebe0339cbc45a37400e2f47be33546c7 Mon Sep 17 00:00:00 2001 From: TaeunKwon Date: Thu, 6 Aug 2020 01:43:05 +0900 Subject: [PATCH 4/9] Update HBHEChannelInfo.h --- DataFormats/HcalRecHit/interface/HBHEChannelInfo.h | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/DataFormats/HcalRecHit/interface/HBHEChannelInfo.h b/DataFormats/HcalRecHit/interface/HBHEChannelInfo.h index 3e84f08770532..5b96c0abbfcd9 100644 --- a/DataFormats/HcalRecHit/interface/HBHEChannelInfo.h +++ b/DataFormats/HcalRecHit/interface/HBHEChannelInfo.h @@ -28,7 +28,7 @@ class HBHEChannelInfo { darkCurrent_{0}, fcByPE_{0}, lambda_{0}, - noisecorr_{0}, + noisecorr_{0}, riseTime_{0.f}, adc_{0}, dFcPerADC_{0.f}, @@ -52,7 +52,7 @@ class HBHEChannelInfo { darkCurrent_{0}, fcByPE_{0}, lambda_{0}, - noisecorr_{0}, + noisecorr_{0}, riseTime_{0.f}, adc_{0}, dFcPerADC_{0.f}, From eba3186fb60546d2bc84bd075f0f05fc31cbaad3 Mon Sep 17 00:00:00 2001 From: TaeunKwon Date: Thu, 6 Aug 2020 02:50:52 +0900 Subject: [PATCH 5/9] Update MahiFit.cc --- RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc | 2 -- 1 file changed, 2 deletions(-) diff --git a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc index 19ca09ca82184..b2a8a87c189a3 100644 --- a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc +++ b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc @@ -82,8 +82,6 @@ void MahiFit::phase1Apply(const HBHEChannelInfo& channelData, tsTOT *= channelData.tsGain(0); tstrig *= channelData.tsGain(0); - std::cout << noisecorr << std::endl; - useTriple = false; if (tstrig >= ts4Thresh_ && tsTOT > 0) { //Average pedestal width (for covariance matrix constraint) From 7185e1e1ddbd060ea3edc6ba886fb3522ea56952 Mon Sep 17 00:00:00 2001 From: TaeunKwon Date: Wed, 12 Aug 2020 01:33:45 +0900 Subject: [PATCH 6/9] Update RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc Co-authored-by: Slava Krutelyov --- RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc index b2a8a87c189a3..ef1fa7306460b 100644 --- a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc +++ b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc @@ -213,8 +213,9 @@ const float MahiFit::minimize(const double noisecorr) const { //Add off-Diagonal components up to first order for (unsigned int i = 1; i < nnlsWork_.tsSize; ++i) { - invCovMat(i - 1, i) += noisecorr * sqrt(nnlsWork_.noiseTerms.coeff(i - 1) * nnlsWork_.noiseTerms.coeff(i)); - invCovMat(i, i - 1) += noisecorr * sqrt(nnlsWork_.noiseTerms.coeff(i) * nnlsWork_.noiseTerms.coeff(i - 1)); + auto const noiseCorrTerm = noisecorr * sqrt(nnlsWork_.noiseTerms.coeff(i - 1) * nnlsWork_.noiseTerms.coeff(i)); + invCovMat(i - 1, i) += noiseCorrTerm; + invCovMat(i, i - 1) += noiseCorrTerm; } float oldChiSq = 9999; From 7347f447ccc269a22e0d3f4cd067f2aa5bf554d6 Mon Sep 17 00:00:00 2001 From: TaeunKwon Date: Wed, 12 Aug 2020 01:36:21 +0900 Subject: [PATCH 7/9] Update MahiFit.cc --- RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc index ef1fa7306460b..019f2fdcf8856 100644 --- a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc +++ b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc @@ -212,10 +212,12 @@ const float MahiFit::minimize(const double noisecorr) const { invCovMat += nnlsWork_.noiseTerms.asDiagonal(); //Add off-Diagonal components up to first order - for (unsigned int i = 1; i < nnlsWork_.tsSize; ++i) { - auto const noiseCorrTerm = noisecorr * sqrt(nnlsWork_.noiseTerms.coeff(i - 1) * nnlsWork_.noiseTerms.coeff(i)); - invCovMat(i - 1, i) += noiseCorrTerm; - invCovMat(i, i - 1) += noiseCorrTerm; + if (noisecorr > 0.){ + for (unsigned int i = 1; i < nnlsWork_.tsSize; ++i) { + auto const noiseCorrTerm = noisecorr * sqrt(nnlsWork_.noiseTerms.coeff(i - 1) * nnlsWork_.noiseTerms.coeff(i)); + invCovMat(i - 1, i) += noiseCorrTerm; + invCovMat(i, i - 1) += noiseCorrTerm; + } } float oldChiSq = 9999; From abefc0f05a530705ec177c6f63962c869aaf45ce Mon Sep 17 00:00:00 2001 From: TaeunKwon Date: Wed, 12 Aug 2020 01:49:16 +0900 Subject: [PATCH 8/9] Update MahiFit.cc --- RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc index 019f2fdcf8856..e5e5b731c6623 100644 --- a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc +++ b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc @@ -212,7 +212,7 @@ const float MahiFit::minimize(const double noisecorr) const { invCovMat += nnlsWork_.noiseTerms.asDiagonal(); //Add off-Diagonal components up to first order - if (noisecorr > 0.){ + if (noisecorr > 0.) { for (unsigned int i = 1; i < nnlsWork_.tsSize; ++i) { auto const noiseCorrTerm = noisecorr * sqrt(nnlsWork_.noiseTerms.coeff(i - 1) * nnlsWork_.noiseTerms.coeff(i)); invCovMat(i - 1, i) += noiseCorrTerm; From 600bf925ae358f6ec1f8954b15edea792ca12fa5 Mon Sep 17 00:00:00 2001 From: TaeunKwon Date: Wed, 12 Aug 2020 02:02:55 +0900 Subject: [PATCH 9/9] Update MahiFit.cc --- RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc index e5e5b731c6623..15775344d1696 100644 --- a/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc +++ b/RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc @@ -212,7 +212,7 @@ const float MahiFit::minimize(const double noisecorr) const { invCovMat += nnlsWork_.noiseTerms.asDiagonal(); //Add off-Diagonal components up to first order - if (noisecorr > 0.) { + if (noisecorr != 0) { for (unsigned int i = 1; i < nnlsWork_.tsSize; ++i) { auto const noiseCorrTerm = noisecorr * sqrt(nnlsWork_.noiseTerms.coeff(i - 1) * nnlsWork_.noiseTerms.coeff(i)); invCovMat(i - 1, i) += noiseCorrTerm;