Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:26:45 

0001 
0002 
0003 // system include files
0004 #include <memory>
0005 #include <iostream>
0006 
0007 // user include files
0008 #include "FWCore/Framework/interface/Frameworkfwd.h"
0009 #include "FWCore/Framework/interface/global/EDFilter.h"
0010 
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/MakerMacros.h"
0013 
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 
0016 #include "DataFormats/MuonReco/interface/Muon.h"
0017 #include "DataFormats/MuonReco/interface/MuonFwd.h"
0018 #include "DataFormats/MuonReco/interface/MuonSelectors.h"
0019 #include "DataFormats/TrackReco/interface/Track.h"
0020 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0021 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
0022 
0023 //
0024 // class declaration
0025 //
0026 
0027 class MuonBadTrackFilter : public edm::global::EDFilter<> {
0028 public:
0029   explicit MuonBadTrackFilter(const edm::ParameterSet&);
0030   ~MuonBadTrackFilter() override;
0031 
0032 private:
0033   bool filter(edm::StreamID iID, edm::Event&, const edm::EventSetup&) const override;
0034   virtual std::string trackInfo(const reco::TrackRef& trackRef) const;
0035   virtual void printMuonProperties(const reco::MuonRef& muonRef) const;
0036 
0037   // ----------member data ---------------------------
0038 
0039   edm::EDGetTokenT<reco::PFCandidateCollection> tokenPFCandidates_;
0040   const bool taggingMode_;
0041   const double ptMin_;
0042   const double chi2Min_;
0043   const double p1_;
0044   const double p2_;
0045   const double p3_;
0046   const bool debug_;
0047 };
0048 
0049 //
0050 // constants, enums and typedefs
0051 //
0052 
0053 //
0054 // static data member definitions
0055 //
0056 
0057 //
0058 // constructors and destructor
0059 //
0060 MuonBadTrackFilter::MuonBadTrackFilter(const edm::ParameterSet& iConfig)
0061     : tokenPFCandidates_(consumes<reco::PFCandidateCollection>(iConfig.getParameter<edm::InputTag>("PFCandidates"))),
0062       taggingMode_(iConfig.getParameter<bool>("taggingMode")),
0063       ptMin_(iConfig.getParameter<double>("ptMin")),
0064       chi2Min_(iConfig.getParameter<double>("chi2Min")),
0065       p1_(iConfig.getParameter<double>("p1")),
0066       p2_(iConfig.getParameter<double>("p2")),
0067       p3_(iConfig.getParameter<double>("p3")),
0068       debug_(iConfig.getParameter<bool>("debug")) {
0069   produces<bool>();
0070 }
0071 
0072 MuonBadTrackFilter::~MuonBadTrackFilter() {}
0073 
0074 //
0075 // member functions
0076 //
0077 
0078 // ------------ method called on each new Event  ------------
0079 bool MuonBadTrackFilter::filter(edm::StreamID iID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0080   using namespace std;
0081   using namespace edm;
0082 
0083   Handle<reco::PFCandidateCollection> pfCandidates;
0084   iEvent.getByToken(tokenPFCandidates_, pfCandidates);
0085 
0086   bool foundBadTrack = false;
0087 
0088   for (unsigned i = 0; i < pfCandidates->size(); ++i) {
0089     const reco::PFCandidate& cand = (*pfCandidates)[i];
0090 
0091     if (std::abs(cand.pdgId()) != 13)
0092       continue;
0093 
0094     if (cand.pt() < ptMin_)
0095       continue;
0096 
0097     if (cand.muonRef().isNull())
0098       continue;
0099 
0100     const reco::MuonRef muon = cand.muonRef();
0101     if (debug_)
0102       printMuonProperties(muon);
0103 
0104     if (muon->muonBestTrack().isAvailable()) {
0105       if (muon->muonBestTrack()->hitPattern().numberOfValidMuonHits() == 0) {
0106         if (muon->globalTrack().isAvailable()) {
0107           if (muon->globalTrack()->normalizedChi2() > chi2Min_) {
0108             foundBadTrack = true;
0109             if (debug_)
0110               cout << "globalTrack numberOfValidMuonHits: " << muon->globalTrack()->hitPattern().numberOfValidMuonHits()
0111                    << " numberOfValidMuonCSCHits: " << muon->globalTrack()->hitPattern().numberOfValidMuonCSCHits()
0112                    << " numberOfValidMuonDTHits: " << muon->globalTrack()->hitPattern().numberOfValidMuonDTHits()
0113                    << " normalizedChi2: " << muon->globalTrack()->normalizedChi2() << endl;
0114             if (debug_)
0115               cout << "muonBestTrack numberOfValidMuonHits: "
0116                    << muon->muonBestTrack()->hitPattern().numberOfValidMuonHits()
0117                    << " numberOfValidMuonCSCHits: " << muon->muonBestTrack()->hitPattern().numberOfValidMuonCSCHits()
0118                    << " numberOfValidMuonDTHits: " << muon->muonBestTrack()->hitPattern().numberOfValidMuonDTHits()
0119                    << " normalizedChi2: " << muon->muonBestTrack()->normalizedChi2() << endl;
0120           }
0121         }
0122       }
0123     }
0124 
0125     // perform same check as for charged hadrons
0126     if (!cand.trackRef().isNull()) {
0127       const reco::TrackRef trackref = cand.trackRef();
0128       const double Pt = trackref->pt();
0129       const double DPt = trackref->ptError();
0130       const double P = trackref->p();
0131       const unsigned int LostHits = trackref->numberOfLostHits();
0132 
0133       if ((DPt / Pt) > (p1_ * sqrt(p2_ * p2_ / P + p3_ * p3_) / (1. + LostHits))) {
0134         foundBadTrack = true;
0135 
0136         if (debug_) {
0137           cout << cand << endl;
0138           cout << "muon \t"
0139                << "track pT = " << Pt << " +/- " << DPt;
0140           cout << endl;
0141         }
0142       }
0143     }
0144 
0145     // check if at least one track has good quality
0146     if (muon->innerTrack().isAvailable()) {
0147       const double P = muon->innerTrack()->p();
0148       const double DPt = muon->innerTrack()->ptError();
0149       if (P != 0) {
0150         if (debug_)
0151           cout << "innerTrack DPt/P: " << DPt / P << endl;
0152         if (DPt / P < 1) {
0153           if (debug_)
0154             cout << "innerTrack good" << endl;
0155           continue;
0156         }
0157       }
0158     }
0159     if (muon->pickyTrack().isAvailable()) {
0160       const double P = muon->pickyTrack()->p();
0161       const double DPt = muon->pickyTrack()->ptError();
0162       if (P != 0) {
0163         if (debug_)
0164           cout << "pickyTrack DPt/P: " << DPt / P << endl;
0165         if (DPt / P < 1) {
0166           if (debug_)
0167             cout << "pickyTrack good" << endl;
0168           continue;
0169         }
0170       }
0171     }
0172     if (muon->globalTrack().isAvailable()) {
0173       const double P = muon->globalTrack()->p();
0174       const double DPt = muon->globalTrack()->ptError();
0175       if (P != 0) {
0176         if (debug_)
0177           cout << "globalTrack DPt/P: " << DPt / P << endl;
0178         if (DPt / P < 1) {
0179           if (debug_)
0180             cout << "globalTrack good" << endl;
0181           continue;
0182         }
0183       }
0184     }
0185     if (muon->tpfmsTrack().isAvailable()) {
0186       const double P = muon->tpfmsTrack()->p();
0187       const double DPt = muon->tpfmsTrack()->ptError();
0188       if (P != 0) {
0189         if (debug_)
0190           cout << "tpfmsTrack DPt/P: " << DPt / P << endl;
0191         if (DPt / P < 1) {
0192           if (debug_)
0193             cout << "tpfmsTrack good" << endl;
0194           continue;
0195         }
0196       }
0197     }
0198     if (muon->dytTrack().isAvailable()) {
0199       const double P = muon->dytTrack()->p();
0200       const double DPt = muon->dytTrack()->ptError();
0201       if (P != 0) {
0202         if (debug_)
0203           cout << "dytTrack DPt/P: " << DPt / P << endl;
0204         if (DPt / P < 1) {
0205           if (debug_)
0206             cout << "dytTrack good" << endl;
0207           continue;
0208         }
0209       }
0210     }
0211     if (debug_)
0212       cout << "No tracks are good" << endl;
0213     foundBadTrack = true;
0214 
0215   }  // end loop over PF candidates
0216 
0217   bool pass = !foundBadTrack;
0218 
0219   iEvent.put(std::make_unique<bool>(pass));
0220 
0221   return taggingMode_ || pass;
0222 }
0223 
0224 std::string MuonBadTrackFilter::trackInfo(const reco::TrackRef& trackRef) const {
0225   std::ostringstream out;
0226 
0227   if (trackRef.isNull()) {
0228     out << "track ref not set";
0229   } else if (!trackRef.isAvailable()) {
0230     out << "track ref not available";
0231   } else {
0232     const reco::Track& track = *trackRef;
0233     out << "pt = " << track.pt() << " +- " << track.ptError() / track.pt() << " chi2 = " << track.normalizedChi2()
0234         << "; Muon Hits: " << track.hitPattern().numberOfValidMuonHits() << "/"
0235         << track.hitPattern().numberOfLostMuonHits() << " (DT: " << track.hitPattern().numberOfValidMuonDTHits() << "/"
0236         << track.hitPattern().numberOfLostMuonDTHits() << " CSC: " << track.hitPattern().numberOfValidMuonCSCHits()
0237         << "/" << track.hitPattern().numberOfLostMuonCSCHits()
0238         << " RPC: " << track.hitPattern().numberOfValidMuonRPCHits() << "/"
0239         << track.hitPattern().numberOfLostMuonRPCHits() << ")"
0240         << "; Valid inner hits:"
0241         << " TRK: " << track.hitPattern().numberOfValidTrackerHits()
0242         << " PIX: " << track.hitPattern().numberOfValidPixelHits();
0243   }
0244   return out.str();
0245 }
0246 
0247 void MuonBadTrackFilter::printMuonProperties(const reco::MuonRef& muonRef) const {
0248   if (!muonRef.isNonnull())
0249     return;
0250 
0251   bool isGL = muonRef->isGlobalMuon();
0252   bool isTR = muonRef->isTrackerMuon();
0253   bool isST = muonRef->isStandAloneMuon();
0254   bool isTPFMS = muonRef->tpfmsTrack().isNonnull() && muonRef->tpfmsTrack()->pt() > 0;
0255   bool isPicky = muonRef->pickyTrack().isNonnull() && muonRef->pickyTrack()->pt() > 0;
0256   bool isDyt = muonRef->dytTrack().isNonnull() && muonRef->dytTrack()->pt() > 0;
0257 
0258   reco::Muon::MuonTrackType tunePType = muonRef->tunePMuonBestTrackType();
0259   std::string tunePTypeStr;
0260   switch (tunePType) {
0261     case reco::Muon::InnerTrack:
0262       tunePTypeStr = "Inner";
0263       break;
0264     case reco::Muon::OuterTrack:
0265       tunePTypeStr = "Outer";
0266       break;
0267     case reco::Muon::CombinedTrack:
0268       tunePTypeStr = "Combined";
0269       break;
0270     case reco::Muon::TPFMS:
0271       tunePTypeStr = "TPFMS";
0272       break;
0273     case reco::Muon::Picky:
0274       tunePTypeStr = "Picky";
0275       break;
0276     case reco::Muon::DYT:
0277       tunePTypeStr = "DYT";
0278       break;
0279     default:
0280       tunePTypeStr = "unknow";
0281       break;
0282   }
0283 
0284   std::cout << "pt " << muonRef->pt() << " eta " << muonRef->eta() << " GL: " << isGL << " TR: " << isTR
0285             << " ST: " << isST << " TPFMS: " << isTPFMS << " Picky: " << isPicky << " DYT: " << isDyt
0286             << " TuneP: " << tunePTypeStr << " nMatches " << muonRef->numberOfMatches() << std::endl;
0287 
0288   if (isGL) {
0289     std::cout << "\tCombined " << trackInfo(muonRef->combinedMuon()) << std::endl;
0290     std::cout << "\tInner " << trackInfo(muonRef->innerTrack()) << std::endl;
0291   }
0292 
0293   if (isST) {
0294     std::cout << "\tOuter " << trackInfo(muonRef->standAloneMuon()) << std::endl;
0295   }
0296 
0297   if (isTR) {
0298     reco::TrackRef trackerMu = muonRef->innerTrack();
0299     // const reco::Track& track = *trackerMu;
0300     std::cout << "\tInner " << trackInfo(trackerMu) << std::endl;
0301     std::cout << "\t\tTMLastStationAngLoose               " << muon::isGoodMuon(*muonRef, muon::TMLastStationAngLoose)
0302               << std::endl
0303               << "\t\tTMLastStationAngTight               " << muon::isGoodMuon(*muonRef, muon::TMLastStationAngTight)
0304               << std::endl
0305               << "\t\tTMLastStationLoose               " << muon::isGoodMuon(*muonRef, muon::TMLastStationLoose)
0306               << std::endl
0307               << "\t\tTMLastStationTight               " << muon::isGoodMuon(*muonRef, muon::TMLastStationTight)
0308               << std::endl
0309               << "\t\tTMOneStationLoose                " << muon::isGoodMuon(*muonRef, muon::TMOneStationLoose)
0310               << std::endl
0311               << "\t\tTMOneStationTight                " << muon::isGoodMuon(*muonRef, muon::TMOneStationTight)
0312               << std::endl
0313               << "\t\tTMLastStationOptimizedLowPtLoose "
0314               << muon::isGoodMuon(*muonRef, muon::TMLastStationOptimizedLowPtLoose) << std::endl
0315               << "\t\tTMLastStationOptimizedLowPtTight "
0316               << muon::isGoodMuon(*muonRef, muon::TMLastStationOptimizedLowPtTight) << std::endl
0317               << "\t\tTMLastStationOptimizedBarrelLowPtLoose "
0318               << muon::isGoodMuon(*muonRef, muon::TMLastStationOptimizedBarrelLowPtLoose) << std::endl
0319               << "\t\tTMLastStationOptimizedBarrelLowPtTight "
0320               << muon::isGoodMuon(*muonRef, muon::TMLastStationOptimizedBarrelLowPtTight) << std::endl
0321               << std::endl;
0322   }
0323 
0324   if (isPicky) {
0325     std::cout << "\tPicky " << trackInfo(muonRef->pickyTrack()) << std::endl;
0326   }
0327 
0328   if (isDyt) {
0329     std::cout << "\tDyt " << trackInfo(muonRef->dytTrack()) << std::endl;
0330   }
0331 
0332   if (isTPFMS) {
0333     std::cout << "\tTPFMS " << trackInfo(muonRef->tpfmsTrack()) << std::endl;
0334   }
0335 
0336   std::cout << "TM2DCompatibilityLoose           " << muon::isGoodMuon(*muonRef, muon::TM2DCompatibilityLoose)
0337             << std::endl
0338             << "TM2DCompatibilityTight           " << muon::isGoodMuon(*muonRef, muon::TM2DCompatibilityTight)
0339             << std::endl;
0340 
0341   if (muonRef->isGlobalMuon() && muonRef->isTrackerMuon() && muonRef->isStandAloneMuon()) {
0342     reco::TrackRef combinedMu = muonRef->combinedMuon();
0343     reco::TrackRef trackerMu = muonRef->track();
0344     reco::TrackRef standAloneMu = muonRef->standAloneMuon();
0345 
0346     double sigmaCombined = combinedMu->ptError() / (combinedMu->pt() * combinedMu->pt());
0347     double sigmaTracker = trackerMu->ptError() / (trackerMu->pt() * trackerMu->pt());
0348     double sigmaStandAlone = standAloneMu->ptError() / (standAloneMu->pt() * standAloneMu->pt());
0349 
0350     bool combined = combinedMu->ptError() / combinedMu->pt() < 0.20;
0351     bool tracker = trackerMu->ptError() / trackerMu->pt() < 0.20;
0352     bool standAlone = standAloneMu->ptError() / standAloneMu->pt() < 0.20;
0353 
0354     double delta1 = combined && tracker ? fabs(1. / combinedMu->pt() - 1. / trackerMu->pt()) /
0355                                               sqrt(sigmaCombined * sigmaCombined + sigmaTracker * sigmaTracker)
0356                                         : 100.;
0357     double delta2 = combined && standAlone ? fabs(1. / combinedMu->pt() - 1. / standAloneMu->pt()) /
0358                                                  sqrt(sigmaCombined * sigmaCombined + sigmaStandAlone * sigmaStandAlone)
0359                                            : 100.;
0360 
0361     double delta3 = standAlone && tracker ? fabs(1. / standAloneMu->pt() - 1. / trackerMu->pt()) /
0362                                                 sqrt(sigmaStandAlone * sigmaStandAlone + sigmaTracker * sigmaTracker)
0363                                           : 100.;
0364 
0365     double delta =
0366         standAloneMu->hitPattern().numberOfValidMuonDTHits() + standAloneMu->hitPattern().numberOfValidMuonCSCHits() > 0
0367             ? std::min(delta3, std::min(delta1, delta2))
0368             : std::max(delta3, std::max(delta1, delta2));
0369 
0370     std::cout << "delta = " << delta << " delta1 " << delta1 << " delta2 " << delta2 << " delta3 " << delta3
0371               << std::endl;
0372 
0373     double ratio = combinedMu->ptError() / combinedMu->pt() / (trackerMu->ptError() / trackerMu->pt());
0374     //if ( ratio > 2. && delta < 3. ) std::cout << "ALARM ! " << ratio << ", " << delta << std::endl;
0375     std::cout << " ratio " << ratio << " combined mu pt " << combinedMu->pt() << std::endl;
0376     //bool quality3 =  ( combinedMu->pt() < 50. || ratio < 2. ) && delta <  3.;
0377   }
0378 
0379   double sumPtR03 = muonRef->isolationR03().sumPt;
0380   double emEtR03 = muonRef->isolationR03().emEt;
0381   double hadEtR03 = muonRef->isolationR03().hadEt;
0382   double relIsoR03 = (sumPtR03 + emEtR03 + hadEtR03) / muonRef->pt();
0383   double sumPtR05 = muonRef->isolationR05().sumPt;
0384   double emEtR05 = muonRef->isolationR05().emEt;
0385   double hadEtR05 = muonRef->isolationR05().hadEt;
0386   double relIsoR05 = (sumPtR05 + emEtR05 + hadEtR05) / muonRef->pt();
0387   std::cout << " 0.3 Rel Iso: " << relIsoR03 << " sumPt " << sumPtR03 << " emEt " << emEtR03 << " hadEt " << hadEtR03
0388             << std::endl;
0389   std::cout << " 0.5 Rel Iso: " << relIsoR05 << " sumPt " << sumPtR05 << " emEt " << emEtR05 << " hadEt " << hadEtR05
0390             << std::endl;
0391   return;
0392 }
0393 
0394 //define this as a plug-in
0395 DEFINE_FWK_MODULE(MuonBadTrackFilter);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!