14 : stripThreshold_(stripThreshold),
15 maxThreshold_(maxThreshold),
17 maxStripTime_(maxStripTime),
19 seedHitIetaMax_(seedHitIetaMax),
22 verboseLevel_(verboseLevel)
42 std::vector<HFRecHit> d1strip;
43 std::vector<HFRecHit> d2strip;
50 for (
auto& it : rec) {
53 if (it.id().depth() == 1) {
60 bool signStripIeta =
false;
63 int stripDepthMax = 0;
66 d1strip.push_back(d1max);
67 signStripIeta = std::signbit(d1max.
id().
ieta());
68 stripIphiMax = d1max.
id().
iphi();
69 stripIetaMax = d1max.
id().
ieta();
70 stripDepthMax = d1max.
id().
depth();
74 for (
auto& it : rec) {
77 if (it.id().depth() == 2) {
80 bool signIeta = std::signbit(it.id().ieta());
81 if (it.id().iphi() == stripIphiMax && signIeta == signStripIeta) {
91 if (d2max.
energy() > 0) d2strip.push_back(d2max);
97 signStripIeta = std::signbit(d2max.
id().
ieta());
98 stripIphiMax = d2max.
id().
iphi();
99 stripIetaMax = d2max.
id().
ieta();
100 stripDepthMax = d2max.
id().
depth();
103 std::stringstream ss;
106 ss <<
" MaxHit in Depth 1: ieta = " << d1max.
id().
ieta() <<
" iphi = " << stripIphiMax
107 <<
" energy = " << d1max.
energy() <<
" time = " << d1max.
time() << std::endl;
110 ss <<
" MaxHit in Depth 2: ieta = " << d2max.
id().
ieta() <<
" iphi = " << stripIphiMax
111 <<
" energy = " << d2max.
energy() <<
" time = " << d2max.
time() << std::endl;
118 for (
auto& it : rec) {
120 bool signIeta = std::signbit(it.id().ieta());
123 ss <<
" HF hit: ieta = " << it.id().ieta() <<
"\t iphi = " << it.id().iphi()
124 <<
"\t depth = " << it.id().depth() <<
"\t time = " << it.time() <<
"\t energy = " 125 << it.energy() <<
"\t flags = " << it.flags() << std::endl;
129 if (it.id().iphi() == stripIphiMax && signIeta == signStripIeta && it.time() <
maxStripTime_) {
130 if (it.id().depth() == 1) {
133 if (d1strip.empty()) {
135 d1strip.push_back(it);
139 for (
auto& it1 : d1strip) {
140 if (it.id().ieta() == it1.id().ieta() && it.energy() == it1.energy()) {
148 for (
auto& it1 : d1strip) {
150 if (
std::abs(it1.id().ieta() - it.id().ieta()) <=
gap_) {
151 d1strip.push_back(it);
156 else if (it.id().depth() == 2) {
159 if (d2strip.empty()) {
161 d2strip.push_back(it);
165 for (
auto& it1 : d2strip) {
166 if (it.id().ieta() == it1.id().ieta() && it.energy() == it1.energy()) {
175 for (
auto& it1 : d2strip) {
177 if (
std::abs(it1.id().ieta() - it.id().ieta()) <=
gap_) {
178 d2strip.push_back(it);
187 ss <<
" Lstrip1 = " << (
int)d1strip.size() <<
" (iphi = " << stripIphiMax
188 <<
") Lstrip2 = " << (
int)d2strip.size() << std::endl <<
" Strip1: ";
189 for (
auto& it1 : d1strip) {
190 ss << it1.energy() <<
" (" << it1.id().ieta() <<
") ";
192 ss << std::endl <<
" Strip2: ";
193 for (
auto& it1 : d2strip) {
194 ss << it1.energy() <<
" (" << it1.id().ieta() <<
") ";
209 int ietaMax1 = -1000;
210 for (
auto& it1 : d1strip) {
211 if (it1.id().ieta() < ietaMin1) ietaMin1 = it1.id().ieta();
212 if (it1.id().ieta() > ietaMax1) ietaMax1 = it1.id().ieta();
215 int ietaMax2 = -1000;
216 for (
auto& it1 : d2strip) {
217 if (it1.id().ieta() < ietaMin2) ietaMin2 = it1.id().ieta();
218 if (it1.id().ieta() > ietaMax2) ietaMax2 = it1.id().ieta();
222 int ietaMin = ietaMin1;
223 int ietaMax = ietaMax1;
225 if (ietaMin2 >= ietaMin1 && ietaMin2 <= ietaMax1) {
226 if (ietaMax2 > ietaMax1) ietaMax = ietaMax2;
227 }
else if (ietaMin2 <= ietaMin1 && ietaMax2 >= ietaMin1) {
228 if (ietaMin2 < ietaMin1) ietaMin = ietaMin2;
229 if (ietaMax2 > ietaMax1) ietaMax = ietaMax2;
234 for (
auto& it1 : d1strip) {
235 if (it1.id().ieta() < ietaMin || it1.id().ieta() > ietaMax)
continue;
236 eStrip += it1.energy();
238 for (
auto& it1 : d2strip) {
239 if (it1.id().ieta() < ietaMin || it1.id().ieta() > ietaMax)
continue;
240 eStrip += it1.energy();
244 ss <<
" ietaMin1 = " << ietaMin1 <<
" ietaMax1 = " << ietaMax1 << std::endl
245 <<
" ietaMin2 = " << ietaMin2 <<
" ietaMax2 = " << ietaMax2 << std::endl
246 <<
" Common strip: ietaMin = " << ietaMin <<
" ietaMax = " << ietaMax << std::endl;
250 double energyIphi1 = 0;
251 double energyIphi2 = 0;
252 int iphi1 = 0, iphi2 = 0;
257 for (
auto& it : rec) {
259 if (it.id().ieta() < ietaMin || it.id().ieta() > ietaMax)
continue;
261 bool neigh = myqual->
topo()->
decIPhi(
id, neighbour);
262 iphi1 = (neigh) ? neighbour.
iphi() : 0;
264 iphi2 = (neigh) ? neighbour.
iphi() : 0;
265 if (it.id().iphi() == iphi1) energyIphi1 += it.energy();
266 else if (it.id().iphi() == iphi2) energyIphi2 += it.energy();
269 double ratio1 = eStrip > 0 ? energyIphi1/eStrip : 0;
270 double ratio2 = eStrip > 0 ? energyIphi2/eStrip : 0;
273 ss <<
" iphi = " << stripIphiMax <<
" iphi1 = " << iphi1 <<
" iphi2 = " 274 << iphi2 << std::endl
275 <<
" Estrip = " << eStrip <<
" EnergyIphi1 = " << energyIphi1
276 <<
" Ratio = " << ratio1 << std::endl
277 <<
" " <<
" EnergyIphi2 = " << energyIphi2
278 <<
" Ratio = " << ratio2 << std::endl;
285 ss <<
" stripPass = false" << std::endl
286 <<
" mark hits in strips now " << std::endl;
290 for (
auto& it1 : d1strip) {
291 if (it1.id().ieta() < ietaMin || it1.id().ieta() > ietaMax)
continue;
293 if (hit != rec.end()) {
296 ss <<
" d1strip: marked hit = " << (*hit) << std::endl;
303 for (
auto& it1 : d2strip) {
304 if (it1.id().ieta() < ietaMin || it1.id().ieta() > ietaMax)
continue;
306 if (hit != rec.end()) {
309 ss <<
" d2strip: marked hit = " << (*hit) << std::endl;
320 ss <<
" stripPass = true" << std::endl;
330 return std::make_unique<HFStripFilter>(
347 desc.
add<
double>(
"stripThreshold", 40.0);
348 desc.
add<
double>(
"maxThreshold", 100.0);
349 desc.
add<
double>(
"timeMax", 6.0);
350 desc.
add<
double>(
"maxStripTime", 10.0);
351 desc.
add<
double>(
"wedgeCut", 0.05);
352 desc.
add<
int>(
"seedHitIetaMax", 35);
353 desc.
add<
int>(
"gap", 2);
354 desc.
add<
int>(
"lstrips", 2);
constexpr float energy() const
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
ParameterDescriptionBase * addUntracked(U const &iLabel, T const &value)
static edm::ParameterSetDescription fillDescription()
bool decIPhi(const HcalDetId &id, HcalDetId &neighbor) const
int depth() const
get the tower depth
bool incIPhi(const HcalDetId &id, HcalDetId &neighbor) const
void runFilter(HFRecHitCollection &rec, const HcalChannelQuality *myqual) const
constexpr float time() const
int ieta() const
get the cell ieta
Abs< T >::type abs(const T &t)
static std::unique_ptr< HFStripFilter > parseParameterSet(const edm::ParameterSet &ps)
ParameterDescriptionBase * add(U const &iLabel, T const &value)
std::vector< T >::iterator iterator
int iphi() const
get the cell iphi
HFStripFilter(double stripThreshold, double maxThreshold, double timeMax, double maxStripTime, double wedgeCut, int seedHitIetaMax, int gap, int lstrips, int verboseLevel)
constexpr void setEnergy(float energy)
const HcalTopology * topo() const