Reference Manual

Go to the documentation of this file.
00001 #include "RecoParticleFlow/PFProducer/interface/PFMuonAlgo.h"
00002 #include "DataFormats/ParticleFlowReco/interface/PFBlock.h"
00003 #include "DataFormats/ParticleFlowReco/interface/PFBlockElement.h"
00004 #include "DataFormats/ParticleFlowReco/interface/PFBlockElementTrack.h"
00005 #include "DataFormats/ParticleFlowReco/interface/PFBlockElementGsfTrack.h"
00006 #include "DataFormats/TrackReco/interface/Track.h"
00007 #include "DataFormats/MuonReco/interface/Muon.h"
00008 #include "DataFormats/MuonReco/interface/MuonSelectors.h"
00009 
00010 bool
00011 PFMuonAlgo::isMuon( const reco::PFBlockElement& elt ) {
00012 
00013   const reco::PFBlockElementTrack* eltTrack 
00014     = dynamic_cast<const reco::PFBlockElementTrack*>(&elt);
00015 
00016   assert ( eltTrack );
00017   reco::MuonRef muonRef = eltTrack->muonRef();
00018   
00019   return isMuon(muonRef);
00020 
00021 }
00022 
00023 bool
00024 PFMuonAlgo::isLooseMuon( const reco::PFBlockElement& elt ) {
00025 
00026   const reco::PFBlockElementTrack* eltTrack 
00027     = dynamic_cast<const reco::PFBlockElementTrack*>(&elt);
00028 
00029   assert ( eltTrack );
00030   reco::MuonRef muonRef = eltTrack->muonRef();
00031 
00032   return isLooseMuon(muonRef);
00033 
00034 }
00035 
00036 bool
00037 PFMuonAlgo::isGlobalTightMuon( const reco::PFBlockElement& elt ) {
00038 
00039   const reco::PFBlockElementTrack* eltTrack 
00040     = dynamic_cast<const reco::PFBlockElementTrack*>(&elt);
00041 
00042   assert ( eltTrack );
00043   reco::MuonRef muonRef = eltTrack->muonRef();
00044   
00045   return isGlobalTightMuon(muonRef);
00046 
00047 }
00048 
00049 bool
00050 PFMuonAlgo::isGlobalLooseMuon( const reco::PFBlockElement& elt ) {
00051 
00052   const reco::PFBlockElementTrack* eltTrack 
00053     = dynamic_cast<const reco::PFBlockElementTrack*>(&elt);
00054 
00055   assert ( eltTrack );
00056   reco::MuonRef muonRef = eltTrack->muonRef();
00057 
00058   return isGlobalLooseMuon(muonRef);
00059 
00060 }
00061 
00062 bool
00063 PFMuonAlgo::isTrackerTightMuon( const reco::PFBlockElement& elt ) {
00064 
00065   const reco::PFBlockElementTrack* eltTrack 
00066     = dynamic_cast<const reco::PFBlockElementTrack*>(&elt);
00067 
00068   assert ( eltTrack );
00069   reco::MuonRef muonRef = eltTrack->muonRef();
00070 
00071   return isTrackerTightMuon(muonRef);
00072 
00073 }
00074 
00075 bool
00076 PFMuonAlgo::isIsolatedMuon( const reco::PFBlockElement& elt ) {
00077 
00078   const reco::PFBlockElementTrack* eltTrack 
00079     = dynamic_cast<const reco::PFBlockElementTrack*>(&elt);
00080 
00081   assert ( eltTrack );
00082   reco::MuonRef muonRef = eltTrack->muonRef();
00083 
00084   return isIsolatedMuon(muonRef);
00085 
00086 }
00087 
00088 
00089 bool
00090 PFMuonAlgo::isMuon(const reco::MuonRef& muonRef ){
00091 
00092   return isGlobalTightMuon(muonRef) || isTrackerTightMuon(muonRef) || isIsolatedMuon(muonRef);
00093 }
00094 
00095 bool
00096 PFMuonAlgo::isLooseMuon(const reco::MuonRef& muonRef ){
00097 
00098   return isGlobalLooseMuon(muonRef) || isTrackerLooseMuon(muonRef);
00099 
00100 }
00101 
00102 bool
00103 PFMuonAlgo::isGlobalTightMuon( const reco::MuonRef& muonRef ) {
00104 
00105  if ( !muonRef.isNonnull() ) return false;
00106 
00107  if ( !muonRef->isGlobalMuon() ) return false;
00108  if ( !muonRef->isStandAloneMuon() ) return false;
00109  
00110  
00111  if ( muonRef->isTrackerMuon() ) { 
00112    
00113    bool result = muon::isGoodMuon(*muonRef,muon::GlobalMuonPromptTight);
00114    
00115    bool isTM2DCompatibilityTight =  muon::isGoodMuon(*muonRef,muon::TM2DCompatibilityTight);   
00116    int nMatches = muonRef->numberOfMatches();
00117    bool quality = nMatches > 2 || isTM2DCompatibilityTight;
00118    
00119    return result && quality;
00120    
00121  } else {
00122  
00123    reco::TrackRef standAloneMu = muonRef->standAloneMuon();
00124    
00125     // No tracker muon -> Request a perfect stand-alone muon, or an even better global muon
00126     bool result = false;
00127       
00128     // Check the quality of the stand-alone muon : 
00129     // good chi**2 and large number of hits and good pt error
00130     if ( ( standAloneMu->hitPattern().numberOfValidMuonDTHits() < 22 &&
00131            standAloneMu->hitPattern().numberOfValidMuonCSCHits() < 15 ) ||
00132          standAloneMu->normalizedChi2() > 10. || 
00133          standAloneMu->ptError()/standAloneMu->pt() > 0.20 ) {
00134       result = false;
00135     } else { 
00136       
00137       reco::TrackRef combinedMu = muonRef->combinedMuon();
00138       reco::TrackRef trackerMu = muonRef->track();
00139             
00140       // If the stand-alone muon is good, check the global muon
00141       if ( combinedMu->normalizedChi2() > standAloneMu->normalizedChi2() ) {
00142         // If the combined muon is worse than the stand-alone, it 
00143         // means that either the corresponding tracker track was not 
00144         // reconstructed, or that the sta muon comes from a late 
00145         // pion decay (hence with a momentum smaller than the track)
00146         // Take the stand-alone muon only if its momentum is larger
00147         // than that of the track
00148         result = standAloneMu->pt() > trackerMu->pt() ;
00149      } else { 
00150         // If the combined muon is better (and good enough), take the 
00151         // global muon
00152         result = 
00153           combinedMu->ptError()/combinedMu->pt() < 
00154           std::min(0.20,standAloneMu->ptError()/standAloneMu->pt());
00155       }
00156     }      
00157 
00158     return result;    
00159   }
00160 
00161   return false;
00162 
00163 }
00164 
00165 bool
00166 PFMuonAlgo::isTrackerTightMuon( const reco::MuonRef& muonRef ) {
00167 
00168   if ( !muonRef.isNonnull() ) return false;
00169     
00170   if(!muonRef->isTrackerMuon()) return false;
00171   
00172   reco::TrackRef trackerMu = muonRef->track();
00173   const reco::Track& track = *trackerMu;
00174   
00175   unsigned nTrackerHits =  track.hitPattern().numberOfValidTrackerHits();
00176   
00177   if(nTrackerHits<=12) return false;
00178   
00179   bool isAllArbitrated = muon::isGoodMuon(*muonRef,muon::AllArbitrated);
00180   
00181   bool isTM2DCompatibilityTight = muon::isGoodMuon(*muonRef,muon::TM2DCompatibilityTight);
00182   
00183   if(!isAllArbitrated || !isTM2DCompatibilityTight)  return false;
00184 
00185   if((trackerMu->ptError()/trackerMu->pt() > 0.10)){
00186     //std::cout<<" PT ERROR > 10 % "<< trackerMu->pt() <<std::endl;
00187     return false;
00188   }
00189   return true;
00190   
00191 }
00192 
00193 bool
00194 PFMuonAlgo::isGlobalLooseMuon( const reco::MuonRef& muonRef ) {
00195 
00196   if ( !muonRef.isNonnull() ) return false;
00197   if ( !muonRef->isGlobalMuon() ) return false;
00198   if ( !muonRef->isStandAloneMuon() ) return false;
00199   
00200   reco::TrackRef standAloneMu = muonRef->standAloneMuon();
00201   reco::TrackRef combinedMu = muonRef->combinedMuon();
00202   reco::TrackRef trackerMu = muonRef->track();
00203  
00204   unsigned nMuonHits =
00205     standAloneMu->hitPattern().numberOfValidMuonDTHits() +
00206     2*standAloneMu->hitPattern().numberOfValidMuonCSCHits();
00207     
00208   bool quality = false;
00209   
00210   if ( muonRef->isTrackerMuon() ){
00211 
00212     bool result = combinedMu->normalizedChi2() < 100.;
00213     
00214     bool laststation =
00215       muon::isGoodMuon(*muonRef,muon::TMLastStationAngTight);
00216         
00217     int nMatches = muonRef->numberOfMatches();    
00218     
00219     quality = laststation && nMuonHits > 12 && nMatches > 1;
00220 
00221     return result && quality;
00222     
00223   }
00224   else{
00225 
00226     // Check the quality of the stand-alone muon : 
00227     // good chi**2 and large number of hits and good pt error
00228     if (  nMuonHits <=15  ||
00229           standAloneMu->normalizedChi2() > 10. || 
00230           standAloneMu->ptError()/standAloneMu->pt() > 0.20 ) {
00231       quality = false;
00232     }
00233    else { 
00234       // If the stand-alone muon is good, check the global muon
00235       if ( combinedMu->normalizedChi2() > standAloneMu->normalizedChi2() ) {
00236         // If the combined muon is worse than the stand-alone, it 
00237         // means that either the corresponding tracker track was not 
00238         // reconstructed, or that the sta muon comes from a late 
00239         // pion decay (hence with a momentum smaller than the track)
00240         // Take the stand-alone muon only if its momentum is larger
00241         // than that of the track
00242 
00243         // Note that here we even take the standAlone if it has a smaller pT, in contrast to GlobalTight
00244         if(standAloneMu->pt() > trackerMu->pt() || combinedMu->normalizedChi2()<5.) quality =  true;
00245       }
00246       else { 
00247         // If the combined muon is better (and good enough), take the 
00248         // global muon
00249         if(combinedMu->ptError()/combinedMu->pt() < std::min(0.20,standAloneMu->ptError()/standAloneMu->pt())) 
00250           quality = true;
00251         
00252       }
00253    }         
00254   }
00255   
00256 
00257   return quality;
00258 
00259 }
00260 
00261 
00262 bool
00263 PFMuonAlgo::isTrackerLooseMuon( const reco::MuonRef& muonRef ) {
00264 
00265   if ( !muonRef.isNonnull() ) return false;
00266   if(!muonRef->isTrackerMuon()) return false;
00267 
00268   reco::TrackRef trackerMu = muonRef->track();
00269 
00270   if(trackerMu->ptError()/trackerMu->pt() > 0.20) return false;
00271 
00272   // this doesn't seem to be necessary on the small samples looked at, but keep it around as insurance
00273   if(trackerMu->pt()>20.) return false;
00274     
00275   bool isAllArbitrated = muon::isGoodMuon(*muonRef,muon::AllArbitrated);
00276   bool isTMLastStationAngTight = muon::isGoodMuon(*muonRef,muon::TMLastStationAngTight);
00277 
00278   bool quality = isAllArbitrated && isTMLastStationAngTight;
00279 
00280   return quality;
00281   
00282 }
00283 
00284 bool
00285 PFMuonAlgo::isIsolatedMuon( const reco::MuonRef& muonRef ){
00286 
00287 
00288   if ( !muonRef.isNonnull() ) return false;
00289   if ( !muonRef->isIsolationValid() ) return false;
00290   
00291   // Isolated Muons which are missed by standard cuts are nearly always global+tracker
00292   if ( !muonRef->isGlobalMuon() ) return false;
00293 
00294   // If it's not a tracker muon, only take it if there are valid muon hits
00295 
00296   reco::TrackRef standAloneMu = muonRef->standAloneMuon();
00297 
00298   if ( !muonRef->isTrackerMuon() ){
00299     if(standAloneMu->hitPattern().numberOfValidMuonDTHits() == 0 &&
00300        standAloneMu->hitPattern().numberOfValidMuonCSCHits() ==0) return false;
00301   }
00302   
00303   // for isolation, take the smallest pt available to reject fakes
00304 
00305   reco::TrackRef combinedMu = muonRef->combinedMuon();
00306   double smallestMuPt = combinedMu->pt();
00307   
00308   if(standAloneMu->pt()<smallestMuPt) smallestMuPt = standAloneMu->pt();
00309   
00310   if(muonRef->isTrackerMuon())
00311     {
00312       reco::TrackRef trackerMu = muonRef->track();
00313       if(trackerMu->pt() < smallestMuPt) smallestMuPt= trackerMu->pt();
00314     }
00315      
00316   double sumPtR03 = muonRef->isolationR03().sumPt;
00317   double emEtR03 = muonRef->isolationR03().emEt;
00318   double hadEtR03 = muonRef->isolationR03().hadEt;
00319   
00320   double relIso = (sumPtR03 + emEtR03 + hadEtR03)/smallestMuPt;
00321 
00322   if(relIso<0.1) return true;
00323   else return false;
00324 }
00325 
00326 bool 
00327 PFMuonAlgo::isTightMuonPOG(const reco::MuonRef& muonRef) {
00328 
00329   if(!muon::isGoodMuon(*muonRef,muon::GlobalMuonPromptTight)) return false;
00330 
00331   if(!muonRef->isTrackerMuon()) return false;
00332   
00333   if(muonRef->numberOfMatches()<2) return false;
00334   
00335   //const reco::TrackRef& combinedMuon = muonRef->combinedMuon();    
00336   const reco::TrackRef& combinedMuon = muonRef->globalTrack();    
00337   
00338   if(combinedMuon->hitPattern().numberOfValidTrackerHits()<11) return false;
00339   
00340   if(combinedMuon->hitPattern().numberOfValidPixelHits()==0) return false;
00341   
00342   if(combinedMuon->hitPattern().numberOfValidMuonHits()==0) return false;  
00343 
00344   return true;
00345 
00346 }
00347 
00348 void 
00349 PFMuonAlgo::printMuonProperties(const reco::MuonRef& muonRef){
00350   
00351   if ( !muonRef.isNonnull() ) return;
00352   
00353   bool isGL = muonRef->isGlobalMuon();
00354   bool isTR = muonRef->isTrackerMuon();
00355   bool isST = muonRef->isStandAloneMuon();
00356 
00357   std::cout<<" GL: "<<isGL<<" TR: "<<isTR<<" ST: "<<isST<<std::endl;
00358   std::cout<<" nMatches "<<muonRef->numberOfMatches()<<std::endl;
00359   
00360   if ( muonRef->isGlobalMuon() ){
00361     reco::TrackRef combinedMu = muonRef->combinedMuon();
00362     std::cout<<" GL,  pt: " << combinedMu->pt() 
00363         << " +/- " << combinedMu->ptError()/combinedMu->pt() 
00364              << " chi**2 GBL : " << combinedMu->normalizedChi2()<<std::endl;
00365     std::cout<< " Total Muon Hits : " << combinedMu->hitPattern().numberOfValidMuonHits()
00366              << "/" << combinedMu->hitPattern().numberOfLostMuonHits()
00367         << " DT Hits : " << combinedMu->hitPattern().numberOfValidMuonDTHits()
00368         << "/" << combinedMu->hitPattern().numberOfLostMuonDTHits()
00369         << " CSC Hits : " << combinedMu->hitPattern().numberOfValidMuonCSCHits()
00370         << "/" << combinedMu->hitPattern().numberOfLostMuonCSCHits()
00371         << " RPC Hits : " << combinedMu->hitPattern().numberOfValidMuonRPCHits()
00372              << "/" << combinedMu->hitPattern().numberOfLostMuonRPCHits()<<std::endl;
00373 
00374     std::cout<<"  # of Valid Tracker Hits "<<combinedMu->hitPattern().numberOfValidTrackerHits()<<std::endl;
00375     std::cout<<"  # of Valid Pixel Hits "<<combinedMu->hitPattern().numberOfValidPixelHits()<<std::endl;
00376   }
00377   if ( muonRef->isStandAloneMuon() ){
00378     reco::TrackRef standAloneMu = muonRef->standAloneMuon();
00379     std::cout<<" ST,  pt: " << standAloneMu->pt() 
00380         << " +/- " << standAloneMu->ptError()/standAloneMu->pt() 
00381         << " eta : " << standAloneMu->eta()  
00382         << " DT Hits : " << standAloneMu->hitPattern().numberOfValidMuonDTHits()
00383         << "/" << standAloneMu->hitPattern().numberOfLostMuonDTHits()
00384         << " CSC Hits : " << standAloneMu->hitPattern().numberOfValidMuonCSCHits()
00385         << "/" << standAloneMu->hitPattern().numberOfLostMuonCSCHits()
00386         << " RPC Hits : " << standAloneMu->hitPattern().numberOfValidMuonRPCHits()
00387         << "/" << standAloneMu->hitPattern().numberOfLostMuonRPCHits()
00388              << " chi**2 STA : " << standAloneMu->normalizedChi2()<<std::endl;
00389       }
00390 
00391 
00392   if ( muonRef->isTrackerMuon() ){
00393     reco::TrackRef trackerMu = muonRef->track();
00394     const reco::Track& track = *trackerMu;
00395     std::cout<<" TR,  pt: " << trackerMu->pt() 
00396         << " +/- " << trackerMu->ptError()/trackerMu->pt() 
00397              << " chi**2 TR : " << trackerMu->normalizedChi2()<<std::endl;    
00398     std::cout<<" nTrackerHits "<<track.hitPattern().numberOfValidTrackerHits()<<std::endl;
00399     std::cout<< "TMLastStationAngLoose               "
00400         << muon::isGoodMuon(*muonRef,muon::TMLastStationAngLoose) << std::endl       
00401         << "TMLastStationAngTight               "
00402         << muon::isGoodMuon(*muonRef,muon::TMLastStationAngTight) << std::endl          
00403         << "TMLastStationLoose               "
00404         << muon::isGoodMuon(*muonRef,muon::TMLastStationLoose) << std::endl       
00405         << "TMLastStationTight               "
00406         << muon::isGoodMuon(*muonRef,muon::TMLastStationTight) << std::endl          
00407         << "TMOneStationLoose                "
00408         << muon::isGoodMuon(*muonRef,muon::TMOneStationLoose) << std::endl       
00409         << "TMOneStationTight                "
00410         << muon::isGoodMuon(*muonRef,muon::TMOneStationTight) << std::endl       
00411         << "TMLastStationOptimizedLowPtLoose " 
00412         << muon::isGoodMuon(*muonRef,muon::TMLastStationOptimizedLowPtLoose) << std::endl
00413         << "TMLastStationOptimizedLowPtTight " 
00414         << muon::isGoodMuon(*muonRef,muon::TMLastStationOptimizedLowPtTight) << std::endl 
00415         << "TMLastStationOptimizedBarrelLowPtLoose " 
00416         << muon::isGoodMuon(*muonRef,muon::TMLastStationOptimizedBarrelLowPtLoose) << std::endl
00417         << "TMLastStationOptimizedBarrelLowPtTight " 
00418         << muon::isGoodMuon(*muonRef,muon::TMLastStationOptimizedBarrelLowPtTight) << std::endl 
00419         << std::endl;
00420 
00421   }
00422 
00423   std::cout<< "TM2DCompatibilityLoose           "
00424       << muon::isGoodMuon(*muonRef,muon::TM2DCompatibilityLoose) << std::endl 
00425       << "TM2DCompatibilityTight           "
00426            << muon::isGoodMuon(*muonRef,muon::TM2DCompatibilityTight) << std::endl;
00427 
00428 
00429 
00430   if (      muonRef->isGlobalMuon() &&  muonRef->isTrackerMuon() &&  muonRef->isStandAloneMuon() ){
00431     reco::TrackRef combinedMu = muonRef->combinedMuon();
00432     reco::TrackRef trackerMu = muonRef->track();
00433     reco::TrackRef standAloneMu = muonRef->standAloneMuon();
00434     
00435     double sigmaCombined = combinedMu->ptError()/(combinedMu->pt()*combinedMu->pt());    
00436     double sigmaTracker = trackerMu->ptError()/(trackerMu->pt()*trackerMu->pt());        
00437     double sigmaStandAlone = standAloneMu->ptError()/(standAloneMu->pt()*standAloneMu->pt());    
00438     
00439     bool combined = combinedMu->ptError()/combinedMu->pt() < 0.20;       
00440     bool tracker = trackerMu->ptError()/trackerMu->pt() < 0.20;          
00441     bool standAlone = standAloneMu->ptError()/standAloneMu->pt() < 0.20;         
00442   
00443     double delta1 =  combined && tracker ?       
00444       fabs(1./combinedMu->pt() -1./trackerMu->pt())      
00445       /sqrt(sigmaCombined*sigmaCombined + sigmaTracker*sigmaTracker) : 100.;     
00446     double delta2 = combined && standAlone ?     
00447       fabs(1./combinedMu->pt() -1./standAloneMu->pt())   
00448       /sqrt(sigmaCombined*sigmaCombined + sigmaStandAlone*sigmaStandAlone) : 100.;       
00449     double delta3 = standAlone && tracker ?      
00450       fabs(1./standAloneMu->pt() -1./trackerMu->pt())    
00451       /sqrt(sigmaStandAlone*sigmaStandAlone + sigmaTracker*sigmaTracker) : 100.;         
00452     
00453     double delta =       
00454       standAloneMu->hitPattern().numberOfValidMuonDTHits()+      
00455       standAloneMu->hitPattern().numberOfValidMuonCSCHits() > 0 ?        
00456       std::min(delta3,std::min(delta1,delta2)) : std::max(delta3,std::max(delta1,delta2));       
00457 
00458     std::cout << "delta = " << delta << " delta1 "<<delta1<<" delta2 "<<delta2<<" delta3 "<<delta3<<std::endl;   
00459     
00460     double ratio =       
00461       combinedMu->ptError()/combinedMu->pt()     
00462       / (trackerMu->ptError()/trackerMu->pt());          
00463     //if ( ratio > 2. && delta < 3. ) std::cout << "ALARM ! " << ratio << ", " << delta << std::endl;
00464     std::cout<<" ratio "<<ratio<<" combined mu pt "<<combinedMu->pt()<<std::endl;
00465     //bool quality3 =  ( combinedMu->pt() < 50. || ratio < 2. ) && delta <  3.;
00466 
00467 
00468   }
00469 
00470     double sumPtR03 = muonRef->isolationR03().sumPt;
00471     double emEtR03 = muonRef->isolationR03().emEt;
00472     double hadEtR03 = muonRef->isolationR03().hadEt;    
00473     double relIsoR03 = (sumPtR03 + emEtR03 + hadEtR03)/muonRef->pt();
00474     double sumPtR05 = muonRef->isolationR05().sumPt;
00475     double emEtR05 = muonRef->isolationR05().emEt;
00476     double hadEtR05 = muonRef->isolationR05().hadEt;    
00477     double relIsoR05 = (sumPtR05 + emEtR05 + hadEtR05)/muonRef->pt();
00478     std::cout<<" 0.3 Radion Rel Iso: "<<relIsoR03<<" sumPt "<<sumPtR03<<" emEt "<<emEtR03<<" hadEt "<<hadEtR03<<std::endl;
00479     std::cout<<" 0.5 Radion Rel Iso: "<<relIsoR05<<" sumPt "<<sumPtR05<<" emEt "<<emEtR05<<" hadEt "<<hadEtR05<<std::endl;
00480   return;
00481 
00482 }
CMSSW_4_4_3_patch1/src/RecoParticleFlow/PFProducer/src/PFMuonAlgo.cc
