CMS 3D CMS Logo

HBHENegativeEFilter.cc
Go to the documentation of this file.
1 #include <cmath>
2 #include <climits>
3 
5 
7  const std::vector<PiecewiseScalingPolynomial>& a1vec,
8  const std::vector<PiecewiseScalingPolynomial>& a2vec,
9  const std::vector<uint32_t>& iEtaLimits,
10  const std::vector<std::pair<double,double> >& cut,
11  const double minCharge,
12  const unsigned firstTimeSlice,
13  const unsigned lastTimeSlice)
14  : a1v_(a1vec),
15  a2v_(a2vec),
16  iEtaLimits_(iEtaLimits),
17  cut_(cut),
18  minCharge_(minCharge),
19  tFirst_(firstTimeSlice),
20  tLast_(lastTimeSlice)
21 {
22  if (!validate()) throw cms::Exception(
23  "Invalid HBHENegativeEFilter constructor arguments");
24 }
25 
27 {
28  if (cut_.empty())
29  return true;
30 
31  const std::size_t nLimits(iEtaLimits_.size());
32  if (nLimits >= static_cast<std::size_t>(UINT_MAX - 1U))
33  return false;
34  for (std::size_t i=1; i<nLimits; ++i)
35  if (!(iEtaLimits_[i-1] < iEtaLimits_[i]))
36  return false;
37 
38  if (a1v_.size() != nLimits + 1)
39  return false;
40  if (a2v_.size() != nLimits + 1)
41  return false;
42 
43  const std::size_t sz = cut_.size();
44  if (sz >= static_cast<std::size_t>(UINT_MAX - 1U))
45  return false;
46  for (std::size_t i=1; i<sz; ++i)
47  if (!(cut_[i-1U].first < cut_[i].first))
48  return false;
49 
50  if (tFirst_ < 2U)
51  return false;
52  if (!(tFirst_ <= tLast_))
53  return false;
54 
55  return true;
56 }
57 
59 {
60  if (cut_.empty() && r.cut_.empty())
61  return true;
62  else
63  return a1v_ == r.a1v_ &&
64  a2v_ == r.a2v_ &&
65  iEtaLimits_ == r.iEtaLimits_ &&
66  cut_ == r.cut_ &&
67  minCharge_ == r.minCharge_ &&
68  tFirst_ == r.tFirst_ &&
69  tLast_ == r.tLast_;
70 }
71 
73 {
74  const unsigned nLimits = iEtaLimits_.size();
75  unsigned which(0U);
76  if (nLimits)
77  {
78  const uint32_t uEta = std::abs(id.ieta());
79  const uint32_t* limits(&iEtaLimits_[0]);
80  for (; which<nLimits; ++which)
81  if (uEta < limits[which])
82  break;
83  }
84  return which;
85 }
86 
88  const double* ts, const unsigned lenTS) const
89 {
90  bool passes = true;
91  const unsigned sz = cut_.size();
92  if (sz)
93  {
94  double chargeInWindow = 0.0;
95  for (unsigned i=tFirst_; i<=tLast_ && i<lenTS; ++i)
96  chargeInWindow += ts[i];
97  if (chargeInWindow >= minCharge_)
98  {
99  // Figure out the cut value for this charge
100  const std::pair<double,double>* cut = &cut_[0];
101  double cutValue = cut[0].second;
102  if (sz > 1U)
103  {
104  // First point larger than charge
105  unsigned largerPoint = 0;
106  for (; cut[largerPoint].first <= chargeInWindow; ++largerPoint) {}
107 
108  // Constant extrapolation beyond min and max coords
109  if (largerPoint >= sz)
110  cutValue = cut[sz - 1U].second;
111  else if (largerPoint)
112  {
113  const double slope = (cut[largerPoint].second - cut[largerPoint-1U].second)/
114  (cut[largerPoint].first - cut[largerPoint-1U].first);
115  cutValue = cut[largerPoint-1U].second + slope*
116  (chargeInWindow - cut[largerPoint-1U].first);
117  }
118  }
119 
120  // Compare the modified time slices with the cut
121  const unsigned itaIdx = getEtaIndex(id);
122  const PiecewiseScalingPolynomial& a1(a1v_[itaIdx]);
123  const PiecewiseScalingPolynomial& a2(a2v_[itaIdx]);
124 
125  for (unsigned i=tFirst_; i<=tLast_ && i<lenTS && passes; ++i)
126  {
127  const double ecorr = ts[i] - a1(ts[i-1U]) - a2(ts[i-2U]);
128  passes = ecorr >= cutValue;
129  }
130  }
131  }
132  return passes;
133 }
std::vector< PiecewiseScalingPolynomial > a1v_
unsigned getEtaIndex(const HcalDetId &id) const
std::vector< PiecewiseScalingPolynomial > a2v_
static const double slope[3]
std::vector< std::pair< double, double > > cut_
std::vector< uint32_t > iEtaLimits_
bool checkPassFilter(const HcalDetId &id, const double *ts, unsigned lenTS) const
bool operator==(const HBHENegativeEFilter &r) const
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
def which(cmd)
Definition: eostools.py:336