Merge pull request #21821 from mariadalfonso/PRcleanupM2

HBHE: M2 small fixed

HLTrigger/Configuration/python/customizeHLTforCMSSW.py

@@ -17,6 +17,21 @@
 #                     pset.minGoodStripCharge = cms.PSet(refToPSet_ = cms.string('HLTSiStripClusterChargeCutNone'))
 #     return process
 
+def customiseFor21821(process):
+    for producer in producers_by_type(process, "HBHEPhase1Reconstructor"):
+        producer.algorithm.ts4Max = cms.vdouble(100., 20000., 30000)
+        del producer.algorithm.pedestalUpperLimit
+        del producer.algorithm.pedSigmaHPD
+        del producer.algorithm.pedSigmaSiPM
+        del producer.algorithm.noiseHPD
+        del producer.algorithm.noiseSiPM
+
+    for producer in producers_by_type(process, "HcalHitReconstructor"):
+        if hasattr(producer,"puCorrMethod"):
+            del producer.puCorrMethod
+
+    return process
+
 def customiseFor21664_forMahiOn(process):
     for producer in producers_by_type(process, "HBHEPhase1Reconstructor"):
         producer.algorithm.useMahi   = cms.bool(True)
@@ -66,5 +81,7 @@ def customizeHLTforCMSSW(process, menuType="GRun"):
 
     # add call to action function in proper order: newest last!
     # process = customiseFor12718(process)
-
+
+    process = customiseFor21821(process)
+
     return process

RecoLocalCalo/HcalRecAlgos/interface/HcalDeterministicFit.h

@@ -19,7 +19,8 @@ class HcalDeterministicFit {
   ~HcalDeterministicFit();
 
   enum FType {shapeLandau, shape205, shape206, shape207};
-  void init(HcalTimeSlew::ParaSource tsParam, HcalTimeSlew::BiasSetting bias, bool iApplyTimeSlew, PedestalSub pedSubFxn, double respCorr);
+
+  void init(HcalTimeSlew::ParaSource tsParam, HcalTimeSlew::BiasSetting bias, bool iApplyTimeSlew, double respCorr);
 
   void phase1Apply(const HBHEChannelInfo& channelData,
 		   float& reconstructedEnergy,

RecoLocalCalo/HcalRecAlgos/interface/HcalSimpleRecAlgo.h

@@ -18,11 +18,6 @@
 #include "CalibCalorimetry/HcalAlgos/interface/HcalPulseContainmentManager.h"
 #include "CondFormats/HcalObjects/interface/AbsOOTPileupCorrection.h"
 
-#include "RecoLocalCalo/HcalRecAlgos/interface/PulseShapeFitOOTPileupCorrection.h"
-#include "RecoLocalCalo/HcalRecAlgos/interface/HcalDeterministicFit.h"
-
-#include "RecoLocalCalo/HcalRecAlgos/interface/PedestalSub.h"
-
 #include "CalibCalorimetry/HcalAlgos/interface/HcalTimeSlew.h"
 
 /** \class HcalSimpleRecAlgo
@@ -88,13 +83,6 @@ class HcalSimpleRecAlgo {
   int puCorrMethod_;
 
   const HcalTimeSlew* hcalTimeSlew_delay_;
-
-  std::unique_ptr<PulseShapeFitOOTPileupCorrection> psFitOOTpuCorr_;
-
-  std::unique_ptr<PedestalSub> pedSubFxn_;
-
-  // S.Brandt Feb19 : Add a pointer to the HLT algo
-  std::unique_ptr<HcalDeterministicFit> hltOOTpuCorr_;
 };
 
 #endif

RecoLocalCalo/HcalRecAlgos/interface/PedestalSub.h

@@ -11,18 +11,11 @@ class PedestalSub
   PedestalSub();
   ~PedestalSub();
 
-  void init(int runCond, float threshold, float quantile);
-
   void calculate(const std::vector<double> & inputCharge, const std::vector<double> & inputPedestal, const std::vector<double> & inputNoise, std::vector<double> & corrCharge, int soi, int nSample) const;
 
   double getCorrection(const std::vector<double> & inputCharge, const std::vector<double> & inputPedestal, const std::vector<double> & inputNoise, int soi, int nSample) const;
 
 
- private:
-  float fThreshold;
-  float fQuantile;
-  float fCondition;
-
 };
 
 #endif 

RecoLocalCalo/HcalRecAlgos/interface/PulseShapeFitOOTPileupCorrection.h

@@ -42,7 +42,6 @@ class PulseShapeFitOOTPileupCorrection
 		     const std::vector<double> & its4Chi2, HcalTimeSlew::BiasSetting slewFlavor, int iFitTimes);
 
     const HcalPulseShapes::Shape* currentPulseShape_=nullptr;
-    void setChi2Term( bool isHPD );
 
     void setPulseShapeTemplate  (const HcalPulseShapes::Shape& ps, bool isHPD, unsigned nSamples);
     void resetPulseShapeTemplate(const HcalPulseShapes::Shape& ps, unsigned nSamples);
@@ -85,11 +84,6 @@ class PulseShapeFitOOTPileupCorrection
     double timeSigSiPM_;
     double pedMean_;
     double pedSig_;
-    double pedSigHPD_;
-    double pedSigSiPM_;
-    double noise_;    
-    double noiseHPD_;
-    double noiseSiPM_;
     HcalTimeSlew::BiasSetting slewFlavor_;    
 
     bool isCurrentChannelHPD_;

RecoLocalCalo/HcalRecAlgos/interface/PulseShapeFunctor.h

@@ -18,9 +18,9 @@ namespace FitterFuncs{
 
    class PulseShapeFunctor {
       public:
-     PulseShapeFunctor(const HcalPulseShapes::Shape& pulse,bool iPedestalConstraint, bool iTimeConstraint,bool iAddPulseJitter,bool iAddTimeSlew,
-		       double iPulseJitter,double iTimeMean,double iTimeSig,double iPedMean,double iPedSig,
-		       double iNoise, unsigned int nSamplesToFit);
+     PulseShapeFunctor(const HcalPulseShapes::Shape& pulse,bool iPedestalConstraint, bool iTimeConstraint,bool iAddPulseJitter,
+		       double iPulseJitter,double iTimeMean,double iPedMean,
+		       unsigned int nSamplesToFit);
      ~PulseShapeFunctor();
 
      double EvalPulse(const double *pars, unsigned int nPar);
@@ -35,6 +35,7 @@ namespace FitterFuncs{
      void setpsFitslew (double *slew  ){ for(int i=0; i<HcalConst::maxSamples; ++i) {psFit_slew [i] = slew [i]; } }
      double getSiPMDarkCurrent(double darkCurrent, double fcByPE, double lambda);
      void setinvertpedSig2(double x) { invertpedSig2_ = x; }
+     void setinverttimeSig2(double x) { inverttimeSig2_ = x; }
 
      double singlePulseShapeFunc( const double *x );
      double doublePulseShapeFunc( const double *x );
@@ -63,12 +64,10 @@ namespace FitterFuncs{
      double timeMean_;
      double timeSig_;
      double pedMean_;
-     double pedSig_;
-     double noise_;
      double timeShift_;
 
-     double inverttimeSig_, inverttimeSig2_;
-     double invertpedSig_, invertpedSig2_;
+     double inverttimeSig2_;
+     double invertpedSig2_;
      std::array<double,HcalConst::maxSamples> pulse_shape_;
      std::array<double,HcalConst::maxSamples> pulse_shape_sum_;
 

RecoLocalCalo/HcalRecAlgos/src/HcalDeterministicFit.cc

@@ -17,11 +17,10 @@ HcalDeterministicFit::HcalDeterministicFit() {
 HcalDeterministicFit::~HcalDeterministicFit() { 
 }
 
-void HcalDeterministicFit::init(HcalTimeSlew::ParaSource tsParam, HcalTimeSlew::BiasSetting bias, bool iApplyTimeSlew, PedestalSub pedSubFxn, double respCorr) {
+void HcalDeterministicFit::init(HcalTimeSlew::ParaSource tsParam, HcalTimeSlew::BiasSetting bias, bool iApplyTimeSlew, double respCorr) {
   fTimeSlew_=tsParam;
   fTimeSlewBias_=bias;
   applyTimeSlew_=iApplyTimeSlew;
-  fPedestalSubFxn_=pedSubFxn;
   frespCorr_=respCorr;
 }
 
@@ -130,7 +129,7 @@ void HcalDeterministicFit::phase1Apply(const HBHEChannelInfo& channelData,
   }
 
   float ch3,ch4,ch5, i3,n3,nn3, i4,n4,i5,n5;
-  ch3=0.f,ch4=0.f,ch5=0.f,i3=0.f,n3=0.f,nn3=0.f,i4=0.f,n4=0.f,i5=0.f,n5=0.f;
+  ch4=0.f,i3=0.f,n3=0.f,nn3=0.f,i4=0.f,n4=0.f,i5=0.f,n5=0.f;
 
   FType fType;
   if(channelData.hasTimeInfo() && channelData.recoShape()==205) fType = shape205;

RecoLocalCalo/HcalRecAlgos/src/HcalSimpleRecAlgo.cc

@@ -24,8 +24,6 @@ HcalSimpleRecAlgo::HcalSimpleRecAlgo(bool correctForTimeslew, bool correctForPul
 {  
   hcalTimeSlew_delay_ = nullptr;
   pulseCorr_ = std::make_unique<HcalPulseContainmentManager>(MaximumFractionalError);
-  pedSubFxn_ = std::make_unique<PedestalSub>();
-  hltOOTpuCorr_ = std::make_unique<HcalDeterministicFit>();
 }
 
 
@@ -278,8 +276,7 @@ namespace HcalSimpleRecAlgoImpl {
                      const int runnum, const bool useLeak,
                      const AbsOOTPileupCorrection* pileupCorrection,
                      const BunchXParameter* bxInfo, const unsigned lenInfo, 
-		     const int puCorrMethod, const PulseShapeFitOOTPileupCorrection * psFitOOTpuCorr,
-		     HcalDeterministicFit * hltOOTpuCorr, PedestalSub * hltPedSub, /* whatever don't know what to do with the pointer...*/
+		     const int puCorrMethod,
 		     const HcalTimeSlew* hcalTimeSlew_delay_)
   {
     double fc_ampl =0, ampl =0, uncorr_ampl =0, m3_ampl =0, maxA = -1.e300;
@@ -336,7 +333,7 @@ HORecHit HcalSimpleRecAlgo::reconstruct(const HODataFrame& digi, int first, int
 							   HcalTimeSlew::Slow,
                                                            runnum_, false, hoPileupCorr_.get(),
                                                            bunchCrossingInfo_, lenBunchCrossingInfo_, 
-							   puCorrMethod_, psFitOOTpuCorr_.get(),/*hlt*/hltOOTpuCorr_.get(),pedSubFxn_.get(),
+							   puCorrMethod_,
 							   hcalTimeSlew_delay_);
 }
 
@@ -348,7 +345,7 @@ HcalCalibRecHit HcalSimpleRecAlgo::reconstruct(const HcalCalibDataFrame& digi, i
 									 HcalTimeSlew::Fast,
                                                                          runnum_, false, nullptr,
                                                                          bunchCrossingInfo_, lenBunchCrossingInfo_, 
-									 puCorrMethod_, psFitOOTpuCorr_.get(),/*hlt*/hltOOTpuCorr_.get(),pedSubFxn_.get(),
+									 puCorrMethod_,
 									 hcalTimeSlew_delay_);
 }
 

RecoLocalCalo/HcalRecAlgos/src/HybridMinimizer.cc

@@ -190,7 +190,6 @@ bool HybridMinimizer::SetVariable(unsigned int ivar, const std::string & name, d
       std::string txtmsg("Wrong index used for the variable " + name);
       MN_INFO_MSG2("HybridMinimizer::SetVariable",txtmsg);  
       MN_INFO_VAL2("HybridMinimizer::SetVariable",minuit2Index);  
-      ivar = minuit2Index;
       return false;
    }
    fState.RemoveLimits(ivar);
@@ -320,17 +319,19 @@ bool HybridMinimizer::Minimize() {
    int strategyLevel = Strategy(); 
    fMinuitFCN->SetErrorDef(ErrorDef() );
 
+   /*
    if (PrintLevel() >=1) { 
       // print the real number of maxfcn used (defined in ModularFuncitonMinimizer)
       int maxfcn_used = maxfcn; 
       if (maxfcn_used == 0) { 
-         int nvar = fState.VariableParameters();
-         maxfcn_used = 200 + 100*nvar + 5*nvar*nvar;
+	int nvar = fState.VariableParameters();
+	maxfcn_used = 200 + 100*nvar + 5*nvar*nvar;
       }      
 //      std::cout << "HybridMinimizer: Minimize with max-calls " << maxfcn_used 
 //                << " convergence for edm < " << tol << " strategy " 
 //                << strategyLevel << std::endl; 
    }
+   */
 
    // internal minuit messages
    MnPrint::SetLevel(PrintLevel() );
@@ -701,17 +702,18 @@ bool HybridMinimizer::GetMinosError(unsigned int i, double & errLow, double & er
    // cut off too small tolerance (they are not needed)
    tol = std::max(tol, 0.01);
 
+   /*
    if (PrintLevel() >=1) { 
       // print the real number of maxfcn used (defined in MnMinos)
       int maxfcn_used = maxfcn; 
       if (maxfcn_used == 0) { 
-         int nvar = fState.VariableParameters();
-         maxfcn_used = 2*(nvar+1)*(200 + 100*nvar + 5*nvar*nvar);
+	int nvar = fState.VariableParameters();
+	maxfcn_used = 2*(nvar+1)*(200 + 100*nvar + 5*nvar*nvar);
       }
 //      std::cout << "HybridMinimizer::GetMinosError for parameter " << i << "  " << par_name
 //                << " using max-calls " << maxfcn_used << ", tolerance " << tol << std::endl; 
    }
-
+   */
 
    if (runLower) low = minos.Loval(i,maxfcn,tol);
    if (runUpper) up  = minos.Upval(i,maxfcn,tol);

RecoLocalCalo/HcalRecAlgos/src/MahiFit.cc

@@ -498,10 +498,11 @@ void MahiFit::resetPulseShapeTemplate(const HcalPulseShapes::Shape& ps) {
 
   // only the pulse shape itself from PulseShapeFunctor is used for Mahi
   // the uncertainty terms calculated inside PulseShapeFunctor are used for Method 2 only
-  psfPtr_.reset(new FitterFuncs::PulseShapeFunctor(ps,false,false,false,false,
-						   1,0,2.5,0,0.00065,1,10));
+  psfPtr_.reset(new FitterFuncs::PulseShapeFunctor(ps,false,false,false,
+						   1,0,0,10));
   pfunctor_ = std::unique_ptr<ROOT::Math::Functor>( new ROOT::Math::Functor(psfPtr_.get(),&FitterFuncs::PulseShapeFunctor::singlePulseShapeFunc, 3) );
 
+
 }
 
 void MahiFit::nnlsUnconstrainParameter(Index idxp) const {

RecoLocalCalo/HcalRecAlgos/src/PedestalSub.cc

@@ -5,17 +5,12 @@
 
 using namespace std;
 
-PedestalSub::PedestalSub() : fThreshold(2.7),fQuantile(0.0),fCondition(0){
+PedestalSub::PedestalSub() {
 }
 
 PedestalSub::~PedestalSub() { 
 }
 
-void PedestalSub::init(int runCond=0, float threshold=0.0, float quantile=0.0) {
-  fThreshold=threshold;
-  fQuantile=quantile;
-  fCondition=runCond;
-}
 
 void PedestalSub::calculate(const std::vector<double> & inputCharge, const std::vector<double> & inputPedestal, const std::vector<double> & inputNoise, std::vector<double> & corrCharge, int soi, int nTS) const {
 

RecoLocalCalo/HcalRecAlgos/src/PulseShapeFitOOTPileupCorrection.cc

@@ -9,8 +9,8 @@ PulseShapeFitOOTPileupCorrection::PulseShapeFitOOTPileupCorrection() : cntsetPul
 								       psfPtr_(nullptr), spfunctor_(nullptr), dpfunctor_(nullptr), tpfunctor_(nullptr),
 								       TSMin_(0), TSMax_(0), vts4Chi2_(0), pedestalConstraint_(false),
 								       timeConstraint_(false), addPulseJitter_(false), applyTimeSlew_(false),
-								       ts4Min_(0), vts4Max_(0), pulseJitter_(0), timeMean_(0), timeSig_(0), pedMean_(0), pedSig_(0),
-								       noise_(0) {
+								       ts4Min_(0), vts4Max_(0), pulseJitter_(0), timeMean_(0), timeSig_(0), pedMean_(0)
+								       {
    hybridfitter = new PSFitter::HybridMinimizer(PSFitter::HybridMinimizer::kMigrad);
    iniTimesArr = { {-100,-75,-50,-25,0,25,50,75,100,125} };
 }
@@ -25,13 +25,6 @@ double PulseShapeFitOOTPileupCorrection::getSiPMDarkCurrent(double darkCurrent,
   return sqrt(mu/pow(1-lambda,3)) * fcByPE;
 }
 
-void PulseShapeFitOOTPileupCorrection::setChi2Term( bool isHPD ) {
-
-  if(isHPD) timeSig_       = timeSigHPD_;
-  else timeSig_            = timeSigSiPM_;
-
-}
-
 
 void PulseShapeFitOOTPileupCorrection::setPUParams(bool   iPedestalConstraint, bool iTimeConstraint,bool iAddPulseJitter,
 						   bool iApplyTimeSlew,double iTS4Min, const std::vector<double> & iTS4Max,
@@ -68,17 +61,21 @@ void PulseShapeFitOOTPileupCorrection::setPulseShapeTemplate(const HcalPulseShap
 
   if (!(&ps == currentPulseShape_ && isHPD == isCurrentChannelHPD_))
     {
-      setChi2Term(isHPD);
       resetPulseShapeTemplate(ps,nSamples);
+
+      // redefine the inverttimeSig2
+      if(!isHPD) psfPtr_->setinverttimeSig2(1./(timeSigSiPM_*timeSigSiPM_));
+      else psfPtr_->setinverttimeSig2(1./(timeSigHPD_*timeSigHPD_));
+
       currentPulseShape_ = &ps;
       isCurrentChannelHPD_ = isHPD;
     }
 }
 
 void PulseShapeFitOOTPileupCorrection::resetPulseShapeTemplate(const HcalPulseShapes::Shape& ps, unsigned nSamples) {
    ++ cntsetPulseShape;
-   psfPtr_.reset(new FitterFuncs::PulseShapeFunctor(ps,pedestalConstraint_,timeConstraint_,addPulseJitter_,applyTimeSlew_,
-						    pulseJitter_,timeMean_,timeSig_,pedMean_,pedSig_,noise_,nSamples));
+   psfPtr_.reset(new FitterFuncs::PulseShapeFunctor(ps,pedestalConstraint_,timeConstraint_,addPulseJitter_,
+						    pulseJitter_,timeMean_,pedMean_,nSamples));
    spfunctor_    = std::unique_ptr<ROOT::Math::Functor>( new ROOT::Math::Functor(psfPtr_.get(),&FitterFuncs::PulseShapeFunctor::singlePulseShapeFunc, 3) );
    dpfunctor_    = std::unique_ptr<ROOT::Math::Functor>( new ROOT::Math::Functor(psfPtr_.get(),&FitterFuncs::PulseShapeFunctor::doublePulseShapeFunc, 5) );
    tpfunctor_    = std::unique_ptr<ROOT::Math::Functor>( new ROOT::Math::Functor(psfPtr_.get(),&FitterFuncs::PulseShapeFunctor::triplePulseShapeFunc, 7) );
@@ -306,7 +303,11 @@ void PulseShapeFitOOTPileupCorrection::phase1Apply(const HBHEChannelInfo& channe
   psfPtr_->setinvertpedSig2(1./(averagePedSig2GeV));
 
   if(channelData.hasTimeInfo()) { 
-    ts4Max_=vts4Max_[1]; ts4Chi2_=vts4Chi2_[1]; 
+
+    ts4Chi2_=vts4Chi2_[1];
+    if(channelData.id().depth()==1) ts4Max_=vts4Max_[1];
+    else ts4Max_=vts4Max_[2];
+
   } else {
     ts4Max_=vts4Max_[0]; ts4Chi2_=vts4Chi2_[0];
   }

RecoLocalCalo/HcalRecAlgos/src/PulseShapeFunctor.cc

@@ -8,9 +8,9 @@ namespace FitterFuncs{
 
   //Decalare the Pulse object take it in from Hcal and set some options
   PulseShapeFunctor::PulseShapeFunctor(const HcalPulseShapes::Shape& pulse,
-				       bool iPedestalConstraint, bool iTimeConstraint,bool iAddPulseJitter,bool iAddTimeSlew,
-				       double iPulseJitter,double iTimeMean,double iTimeSig,double iPedMean,double iPedSig,
-				       double iNoise,unsigned nSamplesToFit) :
+				       bool iPedestalConstraint, bool iTimeConstraint,bool iAddPulseJitter,
+				       double iPulseJitter,double iTimeMean,double iPedMean,
+				       unsigned nSamplesToFit) :
     cntNANinfit(0),
     acc25nsVec(HcalConst::maxPSshapeBin), diff25nsItvlVec(HcalConst::maxPSshapeBin),
     accVarLenIdxZEROVec(HcalConst::nsPerBX), diffVarItvlIdxZEROVec(HcalConst::nsPerBX), 
@@ -48,19 +48,13 @@ namespace FitterFuncs{
     timeConstraint_     = iTimeConstraint;
     addPulseJitter_     = iAddPulseJitter;
     pulseJitter_        = iPulseJitter*iPulseJitter;
+
+    // for M2
     timeMean_           = iTimeMean;
-    timeSig_            = iTimeSig;
     pedMean_            = iPedMean;
-    pedSig_             = iPedSig;
-    noise_              = iNoise;
     timeShift_          = 100.;
     timeShift_ += 12.5;//we are trying to get BX 
 
-    inverttimeSig_ = 1./timeSig_;
-    inverttimeSig2_ = inverttimeSig_*inverttimeSig_;
-    invertpedSig_ = 1./pedSig_;
-    invertpedSig2_ = invertpedSig_*invertpedSig_;
-
     nSamplesToFit_ = nSamplesToFit;
 
   }