23 LogTrace(metname) <<
"Muon Trajectory Cleaner called" << endl;
25 TrajectoryContainer::iterator iter, jter;
26 Trajectory::DataContainer::const_iterator m1, m2;
28 if ( trajC.empty() )
return;
30 LogTrace(metname) <<
"Number of trajectories in the container: " <<trajC.size()<< endl;
36 map<int, vector<int> > seedmap;
41 vector<bool>
mask(trajC.size(),
true);
45 for ( iter = trajC.begin(); iter != trajC.end(); iter++ ) {
46 if ( !
mask[
i] ) { i++;
continue; }
47 if ( !(*iter)->isValid() || (*iter)->empty() ) {
mask[
i] =
false; i++;
continue; }
48 if(seedmap.count(i)==0) seedmap[i].push_back(i);
52 for ( jter = iter+1; jter != trajC.end(); jter++ ) {
53 if ( !
mask[j] ) { j++;
continue; }
54 if(seedmap.count(j)==0) seedmap[j].push_back(j);
57 for( m1 = meas1.begin(); m1 != meas1.end(); m1++ ) {
58 for( m2 = meas2.begin(); m2 != meas2.end(); m2++ ) {
59 if ( ( (*m1).recHit()->globalPosition() - (*m2).recHit()->globalPosition()).
mag()< 10
e-5 ) match++;
64 int nTotHits_i = (*iter)->measurements().size();
65 int nTotHits_j = (*jter)->measurements().size();
68 double chi2_dof_i = (*iter)->ndof() > 0 ? (*iter)->chiSquared()/(*iter)->ndof() : (*iter)->chiSquared()/1
e-10;
69 double chi2_dof_j = (*jter)->ndof() > 0 ? (*jter)->chiSquared()/(*jter)->ndof() : (*jter)->chiSquared()/1
e-10;
71 LogTrace(metname) <<
" MuonTrajectoryCleaner:";
72 LogTrace(metname) <<
" * trajC " << i
73 <<
" (pT="<<(*iter)->lastMeasurement().updatedState().globalMomentum().perp()
74 <<
" GeV) - chi2/nDOF = " << (*iter)->chiSquared() <<
"/" << (*iter)->ndof() <<
" = " << chi2_dof_i;
75 LogTrace(metname) <<
" - valid RH = " << (*iter)->foundHits() <<
" / total RH = " << nTotHits_i;
76 LogTrace(metname) <<
" * trajC " << j
77 <<
" (pT="<<(*jter)->lastMeasurement().updatedState().globalMomentum().perp()
78 <<
" GeV) - chi2/nDOF = " << (*jter)->chiSquared() <<
"/" << (*jter)->ndof() <<
" = " << chi2_dof_j;
79 LogTrace(metname) <<
" - valid RH = " << (*jter)->foundHits() <<
" / total RH = " << nTotHits_j;
82 int hit_diff = (*iter)->foundHits() - (*jter)->foundHits() ;
86 if (
abs(hit_diff) == 0 ) {
93 double fraction = (2.*
match)/((*iter)->measurements().size()+(*jter)->measurements().size());
97 if((*jter)->lastMeasurement().updatedState().globalMomentum().perp() <=
minPt) ++belowLimit;
98 if((*iter)->lastMeasurement().updatedState().globalMomentum().perp() <=
minPt) ++belowLimit;
100 if((*jter)->lastMeasurement().updatedState().globalMomentum().perp() >= dPt) ++above;
101 if((*iter)->lastMeasurement().updatedState().globalMomentum().perp() >= dPt) ++above;
104 if(fraction >= maxFraction && belowLimit == 1 && above == 1){
105 if((*iter)->lastMeasurement().updatedState().globalMomentum().perp() <
minPt){
108 seedmap[j].insert(seedmap[j].
end(), seedmap[i].
begin(), seedmap[i].
end());
110 LogTrace(metname) <<
"Trajectory # " << i <<
" (pT="<<(*iter)->lastMeasurement().updatedState().globalMomentum().perp()
111 <<
" GeV) rejected because it has too low pt";
115 seedmap[
i].insert(seedmap[i].
end(), seedmap[j].
begin(), seedmap[j].
end());
117 LogTrace(metname) <<
"Trajectory # " << j <<
" (pT="<<(*jter)->lastMeasurement().updatedState().globalMomentum().perp()
118 <<
" GeV) rejected because it has too low pt";
122 if (chi2_dof_i > chi2_dof_j) {
125 seedmap[j].insert(seedmap[j].
end(), seedmap[i].
begin(), seedmap[i].
end());
127 LogTrace(metname) <<
"Trajectory # " << i <<
" (pT="<<(*iter)->lastMeasurement().updatedState().globalMomentum().perp() <<
" GeV) rejected";
131 seedmap[
i].insert(seedmap[i].
end(), seedmap[j].
begin(), seedmap[j].
end());
133 LogTrace(metname) <<
"Trajectory # " << j <<
" (pT="<<(*jter)->lastMeasurement().updatedState().globalMomentum().perp() <<
" GeV) rejected";
138 if ( hit_diff < 0 ) {
141 seedmap[j].insert(seedmap[j].
end(), seedmap[i].
begin(), seedmap[i].
end());
143 LogTrace(metname) <<
"Trajectory # " << i <<
" (pT="<<(*iter)->lastMeasurement().updatedState().globalMomentum().perp() <<
" GeV) rejected";
147 seedmap[
i].insert(seedmap[i].
end(), seedmap[j].
begin(), seedmap[j].
end());
149 LogTrace(metname) <<
"Trajectory # " << j <<
" (pT="<<(*jter)->lastMeasurement().updatedState().globalMomentum().perp() <<
" GeV) rejected";
157 if(skipnext)
continue;
163 LogTrace(metname) <<
" Creating map between chosen seed and ghost seeds." << std::endl;
165 auto seedToSeedsMap = std::make_unique<L2SeedAssoc>();
166 if(!seeds->empty()) {
169 event.get(ptr.
id(), seedsHandle);
170 seedToSeedsMap = std::make_unique<L2SeedAssoc>(seedsHandle, seedsHandle);
175 for(
map<
int, vector<int> >::iterator itmap=seedmap.begin(); itmap!=seedmap.end(); ++itmap, ++seedcnt) {
179 int tmpsize=(*itmap).second.size();
180 for(
int cnt=0; cnt!=tmpsize; ++cnt) {
183 seedToSeedsMap->insert(tmpL2SeedRef1, tmpL2SeedRef2);
187 event.put(
std::move(seedToSeedsMap),
"");
192 for ( iter = trajC.begin(); iter != trajC.end(); iter++ ) {
194 result.push_back(*iter);
195 LogTrace(metname) <<
"Keep trajectory with pT = " << (*iter)->lastMeasurement().updatedState().globalMomentum().perp() <<
" GeV";
211 LogTrace(metname) <<
"Muon Trajectory Cleaner called" << endl;
213 if ( candC.size() < 2 )
return;
215 CandidateContainer::iterator iter, jter;
216 Trajectory::DataContainer::const_iterator m1, m2;
222 LogTrace(metname) <<
"Number of muon candidates in the container: " <<candC.size()<< endl;
231 vector<bool>
mask(candC.size(),
true);
235 for ( iter = candC.begin(); iter != candC.end(); iter++ ) {
236 if ( !
mask[
i] ) { i++;
continue; }
244 innerTSOS = (*iter)->trajectory()->firstMeasurement().updatedState();
247 innerTSOS = (*iter)->trajectory()->lastMeasurement().updatedState();
249 if ( !(innerTSOS.
isValid()) )
continue;
254 for ( jter = iter+1; jter != candC.end(); jter++ ) {
255 if ( !
mask[j] ) { j++;
continue; }
259 for ( m1 = meas1.begin(); m1 != meas1.end(); m1++ ) {
260 for ( m2 = meas2.begin(); m2 != meas2.end(); m2++ ) {
261 if ( (*m1).recHit()->isValid() && (*m2).recHit()->isValid() )
262 if ( ( (*m1).recHit()->globalPosition() - (*m2).recHit()->globalPosition()).
mag()< 10
e-5 ) match++;
267 <<
" MuonTrajectoryCleaner: candC " << i <<
" chi2/nRH = " 268 << (*iter)->trajectory()->chiSquared() <<
"/" << (*iter)->trajectory()->foundHits() <<
269 " vs trajC " << j <<
" chi2/nRH = " << (*jter)->trajectory()->chiSquared() <<
270 "/" << (*jter)->trajectory()->foundHits() <<
" Shared RecHits: " <<
match;
274 innerTSOS2 = (*jter)->trajectory()->firstMeasurement().updatedState();
277 innerTSOS2 = (*jter)->trajectory()->lastMeasurement().updatedState();
279 if ( !(innerTSOS2.
isValid()) )
continue;
285 float deta(fabs(eta1-eta2));
287 float dpt(fabs(pt1-pt2));
288 if ( dpt < deltaPt && deta < deltaEta && dphi < deltaPhi ) {
291 <<
" MuonTrajectoryCleaner: candC " << i<<
" and "<<j<<
" direction matched: " 296 bool hitsMatch = ((match > 0) && ((match/((*iter)->trajectory()->foundHits()) > 0.25) || (match/((*jter)->trajectory()->foundHits()) > 0.25))) ?
true :
false;
298 ( ( (*iter)->trackerTrack() == (*jter)->trackerTrack() ) &&
299 ( deltaR<double>((*iter)->muonTrack()->eta(),(*iter)->muonTrack()->phi(), (*jter)->muonTrack()->eta(),(*jter)->muonTrack()->phi()) < 0.2) );
302 if ( ( tracksMatch ) || (hitsMatch > 0) ) {
303 if ( (*iter)->trajectory()->foundHits() == (*jter)->trajectory()->foundHits() ) {
304 if ( (*iter)->trajectory()->chiSquared() > (*jter)->trajectory()->chiSquared() ) {
308 else mask[j] =
false;
311 if ( (*iter)->trajectory()->foundHits() < (*jter)->trajectory()->foundHits() ) {
315 else mask[j] =
false;
322 if(skipnext)
continue;
326 for ( iter = candC.begin(); iter != candC.end(); iter++ ) {
328 result.push_back(*iter);
330 delete (*iter)->trajectory();
331 delete (*iter)->trackerTrajectory();
void clean(TrajectoryContainer &muonTrajectories, edm::Event &evt, const edm::Handle< edm::View< TrajectorySeed > > &seeds)
Clean the trajectories container, erasing the (worst) clone trajectory.
const std::string metname
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
Geom::Phi< T > phi() const
static const double deltaEta
edm::AssociationMap< edm::OneToMany< L2MuonTrajectorySeedCollection, L2MuonTrajectorySeedCollection > > L2SeedAssoc
std::vector< TrajectoryMeasurement > DataContainer
MuonCandidate::CandidateContainer CandidateContainer
Abs< T >::type abs(const T &t)
MuonCandidate::TrajectoryContainer TrajectoryContainer
ProductID id() const
Accessor for product ID.
GlobalVector globalMomentum() const
std::pair< typename Association::data_type::first_type, double > match(Reference key, Association association, bool bestMatchByMaxValue)
Generic matching function.