11 double iMeanTime,
double iTimeSigmaHPD,
double iTimeSigmaSiPM,
12 const std::vector <int> &iActiveBXs,
int iNMaxItersMin,
int iNMaxItersNNLS,
13 double iDeltaChiSqThresh,
double iNnlsThresh) {
41 float& reconstructedEnergy,
42 float& reconstructedTime,
69 std::array<float,3> reconstructedVals {{ 0.0, -9999, -9999 }};
71 double tsTOT = 0, tstrig = 0;
82 double noisePhoto = 0;
84 noisePhoto =
sqrt((charge-ped)*channelData.
fcByPE());
91 nnlsWork_.
noiseTerms.coeffRef(iTS) = noiseADC*noiseADC + noisePhoto*noisePhoto + pedWidth*pedWidth;
93 tsTOT += (charge - ped)*channelData.
tsGain(0);
95 tstrig += (charge - ped)*channelData.
tsGain(0);
105 doFit(reconstructedVals,1);
107 doFit(reconstructedVals,0);
112 doFit(reconstructedVals,0);
117 reconstructedVals.at(0) = 0.;
118 reconstructedVals.at(1) = -9999.;
119 reconstructedVals.at(2) = -9999.;
122 reconstructedEnergy = reconstructedVals[0]*channelData.
tsGain(0);
123 reconstructedTime = reconstructedVals[1];
124 chi2 = reconstructedVals[2];
130 unsigned int bxSize=1;
146 for (
unsigned int iBX=0; iBX<bxSize; ++iBX) {
197 bool foundintime =
false;
198 unsigned int ipulseintime = 0;
210 correctedOutput.at(1) = arrivalTime;
211 correctedOutput.at(2) = chiSq;
219 double oldChiSq=9999;
220 double chiSq=oldChiSq;
234 double deltaChiSq = newChiSq - chiSq;
236 if (newChiSq==oldChiSq && newChiSq<chiSq) {
265 const double xx[4]={
t0, 1.0, 0.0, 3};
269 (*pfunctor_)(&xx[0]);
272 (*pfunctor_)(&xxm[0]);
275 (*pfunctor_)(&xxp[0]);
297 pulseCov(iTS,jTS) +=
tmp;
298 pulseCov(jTS,iTS) +=
tmp;
335 if (offset==0) itIndex=iBX;
354 for (
unsigned int iBX=0; iBX<npulse; ++iBX) {
379 const unsigned int nActive = npulse -
nnlsWork_.
nP;
382 Index idxwmaxprev = idxwmax;
383 double wmaxprev = wmax;
386 if (wmax<threshold || (idxwmax==idxwmaxprev && wmax==wmaxprev)) {
408 bool positive =
true;
421 for (
unsigned int ipulse=0; ipulse<
nnlsWork_.
nP; ++ipulse) {
425 if (ratio<minratio) {
427 minratioidx = ipulse;
515 using namespace Eigen;
519 Matrix<double,10,10>
temp = mat;
520 outvec.head<10>() = temp.ldlt().solve(invec.head<10>());
525 Matrix<double,9,9>
temp = mat.topLeftCorner<9,9>();
526 outvec.head<9>() = temp.ldlt().solve(invec.head<9>());
531 Matrix<double,8,8>
temp = mat.topLeftCorner<8,8>();
532 outvec.head<8>() = temp.ldlt().solve(invec.head<8>());
537 Matrix<double,7,7>
temp = mat.topLeftCorner<7,7>();
538 outvec.head<7>() = temp.ldlt().solve(invec.head<7>());
543 Matrix<double,6,6>
temp = mat.topLeftCorner<6,6>();
544 outvec.head<6>() = temp.ldlt().solve(invec.head<6>());
549 Matrix<double,5,5>
temp = mat.topLeftCorner<5,5>();
550 outvec.head<5>() = temp.ldlt().solve(invec.head<5>());
555 Matrix<double,4,4>
temp = mat.topLeftCorner<4,4>();
556 outvec.head<4>() = temp.ldlt().solve(invec.head<4>());
561 Matrix<double,3,3>
temp = mat.topLeftCorner<3,3>();
562 outvec.head<3>() = temp.ldlt().solve(invec.head<3>());
567 Matrix<double,2,2>
temp = mat.topLeftCorner<2,2>();
568 outvec.head<2>() = temp.ldlt().solve(invec.head<2>());
573 Matrix<double,1,1>
temp = mat.topLeftCorner<1,1>();
574 outvec.head<1>() = temp.ldlt().solve(invec.head<1>());
579 <<
"Weird number of pulses encountered in Mahi, module is configured incorrectly!";
void phase1Apply(const HBHEChannelInfo &channelData, float &reconstructedEnergy, float &reconstructedTime, bool &useTriple, float &chi2) const
SamplePulseMatrix invcovp
HcalTimeSlew::BiasSetting slewFlavor_
double tsGain(const unsigned ts) const
void nnlsConstrainParameter(Index minratioidx) const
MahiNnlsWorkspace nnlsWork_
double singlePulseShapeFunc(const double *x)
void solveSubmatrix(PulseMatrix &mat, PulseVector &invec, PulseVector &outvec, unsigned nP) const
double delay(double fC, BiasSetting bias=Medium) const
Returns the amount (ns) by which a pulse of the given number of fC will be delayed by the timeslew ef...
static constexpr int pedestalBX_
void resetPulseShapeTemplate(const HcalPulseShapes::Shape &ps)
SampleDecompLLT covDecomp
PulseVector ampvecpermtest
std::array< FullSampleVector, MaxPVSize > pulseShapeArray
void nnlsUnconstrainParameter(Index idxp) const
std::array< double, MaxSVSize > pulseM
double tsPedestal(const unsigned ts) const
SamplePulseMatrix pulseDerivMat
std::array< double, MaxSVSize > pulseN
double tsRawCharge(const unsigned ts) const
void onePulseMinimize() const
std::unique_ptr< FitterFuncs::PulseShapeFunctor > psfPtr_
Eigen::Matrix< double, FullSampleVectorSize, 1 > FullSampleVector
void swap(edm::DataFrameContainer &lhs, edm::DataFrameContainer &rhs)
void doFit(std::array< float, 3 > &correctedOutput, const int nbx) const
Eigen::Matrix< double, 1, 1 > SingleVector
Abs< T >::type abs(const T &t)
double calculateArrivalTime() const
SampleMatrix pedConstraint
std::vector< int > activeBXs_
std::array< FullSampleVector, MaxPVSize > pulseDerivArray
const unsigned int fullTSofInterest_
SampleMatrix covDecompLinv
void updatePulseShape(double itQ, FullSampleVector &pulseShape, FullSampleVector &pulseDeriv, FullSampleMatrix &pulseCov) const
void resetWorkspace() const
std::unique_ptr< ROOT::Math::Functor > pfunctor_
Eigen::Matrix< double, FullSampleVectorSize, FullSampleVectorSize > FullSampleMatrix
void resize(int bx, unsigned size)
double tsPedestalWidth(const unsigned ts) const
double calculateChiSq() const
unsigned int fullTSOffset
Eigen::Matrix< double, Eigen::Dynamic, 1, 0, PulseVectorSize, 1 > PulseVector
SamplePulseMatrix pulseMat
std::array< FullSampleMatrix, MaxPVSize > pulseCovArray
std::vector< std::vector< double > > tmp
void setParameters(bool iDynamicPed, double iTS4Thresh, double chiSqSwitch, bool iApplyTimeSlew, HcalTimeSlew::BiasSetting slewFlavor, double iMeanTime, double iTimeSigmaHPD, double iTimeSigmaSiPM, const std::vector< int > &iActiveBXs, int iNMaxItersMin, int iNMaxItersNNLS, double iDeltaChiSqThresh, double iNnlsThresh)
const HcalTimeSlew * hcalTimeSlewDelay_
const HcalPulseShapes::Shape * currentPulseShape_
unsigned nSamples() const
std::array< double, MaxSVSize > pulseP
Eigen::Matrix< double, 1, 1 > SingleMatrix
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic, 0, PulseVectorSize, PulseVectorSize > PulseMatrix
float tsDFcPerADC(const unsigned ts) const
void setPulseShapeTemplate(const HcalPulseShapes::Shape &ps, const HcalTimeSlew *hcalTimeSlewDelay)