108 std::vector<std::pair<unsigned,double> > ordered_hits;
109 for(
unsigned i = 0;
i < input->size(); ++
i ) {
110 std::pair<unsigned,double> val = std::make_pair(
i,input->at(
i).energy());
111 auto pos = std::upper_bound(ordered_hits.begin(),ordered_hits.end(),
112 val, greaterByEnergy);
113 ordered_hits.insert(pos,val);
116 for(
const auto& idx_e : ordered_hits ) {
117 const unsigned i = idx_e.first;
118 if( !mask[i] )
continue;
120 int hitlayer = (int)rechit.
layer();
126 if( rechit.
energy() < clean._singleSpikeThresh )
continue;
127 const double rhenergy = rechit.
energy();
130 double surroundingEnergy = rechit.
energy();
131 double neighbourEnergy = 0.0;
132 double layerEnergy = 0.0;
134 if( !mask[neighbour.key()] )
continue;
135 const double sum = neighbour->energy();
136 surroundingEnergy += sum;
137 neighbourEnergy += sum;
138 layerEnergy += neighbour->energy();
143 const double fraction1 = surroundingEnergy/rhenergy;
146 const double f1Cut = ( clean._minS4S1_a*std::log10(rechit.
energy()) +
148 if( fraction1 < f1Cut ) {
152 std::pair<double,double> dcr = dCrack(phi,eta);
157 ( (aeta < 2.85 && dcrmin > 1.0) ||
158 (rhenergy > clean._eneThreshMod*clean._singleSpikeThresh &&
159 fraction1 < f1Cut/clean._fracThreshMod ) ) ) {
164 if( mask[i] && rhenergy > clean._doubleSpikeThresh ) {
166 double surroundingEnergyi = 0.0;
167 double enmax = -999.0;
171 if( !mask[neighbour.key()] )
continue;
172 const double nenergy = neighbour->energy();
173 surroundingEnergyi += nenergy;
174 if( nenergy > enmax ) {
176 mostEnergeticNeighbour = neighbour;
181 double surroundingEnergyj = 0.0;
183 mostEnergeticNeighbour->neighbours4();
186 surroundingEnergyj += neighbour->energy();
189 const double surroundingEnergyFraction =
190 (surroundingEnergyi+surroundingEnergyj) /
191 (rechit.
energy()+mostEnergeticNeighbour->energy()) - 1.;
192 if ( surroundingEnergyFraction < clean._doubleSpikeS6S2 ) {
196 std::pair<double,double> dcr = dCrack(phi,eta);
201 ( (aeta < 2.85 && dcrmin > 1.0) ||
202 (rhenergy > clean._eneThreshMod*clean._doubleSpikeThresh &&
203 surroundingEnergyFraction < clean._doubleSpikeS6S2/clean._fracThreshMod )
206 mask[mostEnergeticNeighbour.
key()] =
false;
void clean(const edm::Handle< reco::PFRecHitCollection > &input, std::vector< bool > &mask)
const PFRecHitRefVector & neighbours4() const
key_type key() const
Accessor for product key.
PFLayer::Layer layer() const
rechit layer
U second(std::pair< T, U > const &p)
Particle flow rechit (rechit + geometry and topology information). See clustering algorithm in PFClus...
Abs< T >::type abs(const T &t)
std::unordered_map< int, spike_cleaning > _thresholds
double energy() const
rechit energy
const REPPoint & positionREP() const