File indexing completed on 2022-09-10 01:03:14
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012 #include "RecoMuon/MuonIdentification/plugins/MuonIdProducer.h"
0013
0014 #include "DataFormats/Math/interface/deltaPhi.h"
0015 #include "DataFormats/MuonDetId/interface/MuonSubdetId.h"
0016 #include "DataFormats/MuonDetId/interface/DTChamberId.h"
0017 #include "DataFormats/MuonDetId/interface/CSCDetId.h"
0018 #include "DataFormats/MuonDetId/interface/RPCDetId.h"
0019 #include "DataFormats/MuonDetId/interface/GEMDetId.h"
0020 #include "DataFormats/MuonDetId/interface/ME0DetId.h"
0021 #include "DataFormats/MuonReco/interface/MuonCocktails.h"
0022 #include "DataFormats/MuonReco/interface/MuonTime.h"
0023 #include "DataFormats/MuonReco/interface/MuonTimeExtra.h"
0024 #include "DataFormats/MuonReco/interface/MuonTimeExtraMap.h"
0025 #include "DataFormats/RecoCandidate/interface/IsoDeposit.h"
0026 #include "DataFormats/RecoCandidate/interface/IsoDepositFwd.h"
0027
0028 #include "PhysicsTools/IsolationAlgos/interface/IsoDepositExtractorFactory.h"
0029 #include "TrackingTools/Records/interface/TrackingComponentsRecord.h"
0030
0031 #include "RecoMuon/MuonIdentification/interface/MuonMesh.h"
0032 #include "RecoMuon/MuonIdentification/interface/MuonKinkFinder.h"
0033
0034 MuonIdProducer::MuonIdProducer(const edm::ParameterSet& iConfig)
0035 : geomTokenRun_(esConsumes<edm::Transition::BeginRun>()),
0036 propagatorToken_(esConsumes(edm::ESInputTag("", "SteppingHelixPropagatorAny"))) {
0037 LogTrace("MuonIdentification") << "RecoMuon/MuonIdProducer :: Constructor called";
0038
0039 produces<reco::MuonCollection>();
0040 produces<reco::CaloMuonCollection>();
0041 produces<reco::MuonTimeExtraMap>("combined");
0042 produces<reco::MuonTimeExtraMap>("dt");
0043 produces<reco::MuonTimeExtraMap>("csc");
0044
0045 minPt_ = iConfig.getParameter<double>("minPt");
0046 minP_ = iConfig.getParameter<double>("minP");
0047 minPCaloMuon_ = iConfig.getParameter<double>("minPCaloMuon");
0048 minNumberOfMatches_ = iConfig.getParameter<int>("minNumberOfMatches");
0049 addExtraSoftMuons_ = iConfig.getParameter<bool>("addExtraSoftMuons");
0050 maxAbsEta_ = iConfig.getParameter<double>("maxAbsEta");
0051 maxAbsDx_ = iConfig.getParameter<double>("maxAbsDx");
0052 maxAbsPullX2_ = iConfig.getParameter<double>("maxAbsPullX");
0053 maxAbsPullX2_ *= maxAbsPullX2_;
0054 maxAbsDy_ = iConfig.getParameter<double>("maxAbsDy");
0055 maxAbsPullY2_ = iConfig.getParameter<double>("maxAbsPullY");
0056 maxAbsPullY2_ *= maxAbsPullY2_;
0057 fillCaloCompatibility_ = iConfig.getParameter<bool>("fillCaloCompatibility");
0058 fillEnergy_ = iConfig.getParameter<bool>("fillEnergy");
0059 storeCrossedHcalRecHits_ = iConfig.getParameter<bool>("storeCrossedHcalRecHits");
0060 fillMatching_ = iConfig.getParameter<bool>("fillMatching");
0061 fillIsolation_ = iConfig.getParameter<bool>("fillIsolation");
0062 fillShowerDigis_ = iConfig.getParameter<bool>("fillShowerDigis");
0063 writeIsoDeposits_ = iConfig.getParameter<bool>("writeIsoDeposits");
0064 fillGlobalTrackQuality_ = iConfig.getParameter<bool>("fillGlobalTrackQuality");
0065 fillGlobalTrackRefits_ = iConfig.getParameter<bool>("fillGlobalTrackRefits");
0066 arbitrateTrackerMuons_ = iConfig.getParameter<bool>("arbitrateTrackerMuons");
0067 selectHighPurity_ = iConfig.getParameter<bool>("selectHighPurity");
0068
0069 fillTrackerKink_ = false;
0070 if (fillGlobalTrackQuality_)
0071 fillTrackerKink_ = iConfig.getParameter<bool>("fillTrackerKink");
0072
0073 ptThresholdToFillCandidateP4WithGlobalFit_ =
0074 iConfig.getParameter<double>("ptThresholdToFillCandidateP4WithGlobalFit");
0075 sigmaThresholdToFillCandidateP4WithGlobalFit_ =
0076 iConfig.getParameter<double>("sigmaThresholdToFillCandidateP4WithGlobalFit");
0077 caloCut_ = iConfig.getParameter<double>("minCaloCompatibility");
0078 arbClean_ = iConfig.getParameter<bool>("runArbitrationCleaner");
0079
0080
0081 const edm::ParameterSet parameters = iConfig.getParameter<edm::ParameterSet>("TrackAssociatorParameters");
0082 edm::ConsumesCollector iC = consumesCollector();
0083 parameters_.loadParameters(parameters, iC);
0084
0085
0086 edm::ParameterSet timingParameters = iConfig.getParameter<edm::ParameterSet>("TimingFillerParameters");
0087 theTimingFiller_ = std::make_unique<MuonTimingFiller>(timingParameters, consumesCollector());
0088
0089
0090 if (fillShowerDigis_ && fillMatching_) {
0091 edm::ParameterSet showerDigiParameters = iConfig.getParameter<edm::ParameterSet>("ShowerDigiFillerParameters");
0092 theShowerDigiFiller_ = std::make_unique<MuonShowerDigiFiller>(showerDigiParameters, consumesCollector());
0093 } else {
0094 theShowerDigiFiller_ = std::make_unique<MuonShowerDigiFiller>();
0095 }
0096
0097 if (fillCaloCompatibility_) {
0098
0099 const auto caloParams = iConfig.getParameter<edm::ParameterSet>("MuonCaloCompatibility");
0100 muonCaloCompatibility_.configure(caloParams);
0101 }
0102
0103 if (fillIsolation_) {
0104
0105 edm::ParameterSet caloExtractorPSet = iConfig.getParameter<edm::ParameterSet>("CaloExtractorPSet");
0106 std::string caloExtractorName = caloExtractorPSet.getParameter<std::string>("ComponentName");
0107 muIsoExtractorCalo_ =
0108 IsoDepositExtractorFactory::get()->create(caloExtractorName, caloExtractorPSet, consumesCollector());
0109
0110 edm::ParameterSet trackExtractorPSet = iConfig.getParameter<edm::ParameterSet>("TrackExtractorPSet");
0111 std::string trackExtractorName = trackExtractorPSet.getParameter<std::string>("ComponentName");
0112 muIsoExtractorTrack_ =
0113 IsoDepositExtractorFactory::get()->create(trackExtractorName, trackExtractorPSet, consumesCollector());
0114
0115 edm::ParameterSet jetExtractorPSet = iConfig.getParameter<edm::ParameterSet>("JetExtractorPSet");
0116 std::string jetExtractorName = jetExtractorPSet.getParameter<std::string>("ComponentName");
0117 muIsoExtractorJet_ =
0118 IsoDepositExtractorFactory::get()->create(jetExtractorName, jetExtractorPSet, consumesCollector());
0119 }
0120 if (fillIsolation_ && writeIsoDeposits_) {
0121 trackDepositName_ = iConfig.getParameter<std::string>("trackDepositName");
0122 produces<reco::IsoDepositMap>(trackDepositName_);
0123 ecalDepositName_ = iConfig.getParameter<std::string>("ecalDepositName");
0124 produces<reco::IsoDepositMap>(ecalDepositName_);
0125 hcalDepositName_ = iConfig.getParameter<std::string>("hcalDepositName");
0126 produces<reco::IsoDepositMap>(hcalDepositName_);
0127 hoDepositName_ = iConfig.getParameter<std::string>("hoDepositName");
0128 produces<reco::IsoDepositMap>(hoDepositName_);
0129 jetDepositName_ = iConfig.getParameter<std::string>("jetDepositName");
0130 produces<reco::IsoDepositMap>(jetDepositName_);
0131 }
0132
0133 inputCollectionLabels_ = iConfig.getParameter<std::vector<edm::InputTag> >("inputCollectionLabels");
0134 const auto inputCollectionTypes = iConfig.getParameter<std::vector<std::string> >("inputCollectionTypes");
0135 if (inputCollectionLabels_.size() != inputCollectionTypes.size())
0136 throw cms::Exception("ConfigurationError")
0137 << "Number of input collection labels is different from number of types. "
0138 << "For each collection label there should be exactly one collection type specified.";
0139 if (inputCollectionLabels_.size() > 7 || inputCollectionLabels_.empty())
0140 throw cms::Exception("ConfigurationError") << "Number of input collections should be from 1 to 7.";
0141
0142 debugWithTruthMatching_ = iConfig.getParameter<bool>("debugWithTruthMatching");
0143 if (debugWithTruthMatching_) {
0144 edm::LogWarning("MuonIdentification")
0145 << "========================================================================\n"
0146 << "Debugging mode with truth matching is turned on!!! Make sure you understand what you are doing!\n"
0147 << "========================================================================\n";
0148
0149 globalGeomToken_ = esConsumes();
0150 }
0151 if (fillGlobalTrackQuality_) {
0152 const auto& glbQualTag = iConfig.getParameter<edm::InputTag>("globalTrackQualityInputTag");
0153 glbQualToken_ = consumes<edm::ValueMap<reco::MuonQuality> >(glbQualTag);
0154 }
0155
0156 if (fillTrackerKink_) {
0157 trackerKinkFinder_ =
0158 std::make_unique<MuonKinkFinder>(iConfig.getParameter<edm::ParameterSet>("TrackerKinkFinderParameters"), iC);
0159 }
0160
0161 if (selectHighPurity_) {
0162 const auto& pvTag = iConfig.getParameter<edm::InputTag>("pvInputTag");
0163 pvToken_ = mayConsume<reco::VertexCollection>(pvTag);
0164 }
0165
0166
0167 meshAlgo_ = std::make_unique<MuonMesh>(iConfig.getParameter<edm::ParameterSet>("arbitrationCleanerOptions"));
0168
0169 edm::InputTag rpcHitTag("rpcRecHits");
0170 rpcHitToken_ = consumes<RPCRecHitCollection>(rpcHitTag);
0171
0172 edm::InputTag gemHitTag("gemRecHits");
0173 gemHitToken_ = consumes<GEMRecHitCollection>(gemHitTag);
0174
0175
0176 inputCollectionTypes_.resize(inputCollectionLabels_.size());
0177 for (unsigned int i = 0; i < inputCollectionLabels_.size(); ++i) {
0178 const auto inputLabel = inputCollectionLabels_[i];
0179 const auto inputType = ICTypes::toKey(inputCollectionTypes[i]);
0180
0181 if (inputType == ICTypes::INNER_TRACKS) {
0182 innerTrackCollectionToken_ = consumes<reco::TrackCollection>(inputLabel);
0183 } else if (inputType == ICTypes::OUTER_TRACKS) {
0184 outerTrackCollectionToken_ = consumes<reco::TrackCollection>(inputLabel);
0185 } else if (inputType == ICTypes::LINKS) {
0186 linkCollectionToken_ = consumes<reco::MuonTrackLinksCollection>(inputLabel);
0187 } else if (inputType == ICTypes::MUONS) {
0188 muonCollectionToken_ = consumes<reco::MuonCollection>(inputLabel);
0189 } else if (inputType == ICTypes::TEV_FIRSTHIT) {
0190 tpfmsCollectionToken_ = consumes<reco::TrackToTrackMap>(inputLabel);
0191 } else if (fillGlobalTrackRefits_ && inputType == ICTypes::TEV_PICKY) {
0192 pickyCollectionToken_ = consumes<reco::TrackToTrackMap>(inputLabel);
0193 } else if (fillGlobalTrackRefits_ && inputType == ICTypes::TEV_DYT) {
0194 dytCollectionToken_ = consumes<reco::TrackToTrackMap>(inputCollectionLabels_.at(i));
0195 }
0196
0197 inputCollectionTypes_[i] = inputType;
0198 }
0199 }
0200
0201 MuonIdProducer::~MuonIdProducer() {
0202
0203 }
0204
0205 void MuonIdProducer::init(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0206 innerTrackCollectionHandle_.clear();
0207 outerTrackCollectionHandle_.clear();
0208 linkCollectionHandle_.clear();
0209 muonCollectionHandle_.clear();
0210
0211 tpfmsCollectionHandle_.clear();
0212 pickyCollectionHandle_.clear();
0213 dytCollectionHandle_.clear();
0214
0215 trackAssociator_.setPropagator(&iSetup.getData(propagatorToken_));
0216
0217 if (fillTrackerKink_)
0218 trackerKinkFinder_->init(iSetup);
0219
0220 for (unsigned int i = 0; i < inputCollectionLabels_.size(); ++i) {
0221 const auto& inputLabel = inputCollectionLabels_[i];
0222 const auto inputType = inputCollectionTypes_[i];
0223 if (inputType == ICTypes::INNER_TRACKS) {
0224 iEvent.getByToken(innerTrackCollectionToken_, innerTrackCollectionHandle_);
0225 if (!innerTrackCollectionHandle_.isValid())
0226 throw cms::Exception("FatalError") << "Failed to get input track collection with label: " << inputLabel;
0227 LogTrace("MuonIdentification") << "Number of input inner tracks: " << innerTrackCollectionHandle_->size();
0228 } else if (inputType == ICTypes::OUTER_TRACKS) {
0229 iEvent.getByToken(outerTrackCollectionToken_, outerTrackCollectionHandle_);
0230 if (!outerTrackCollectionHandle_.isValid())
0231 throw cms::Exception("FatalError") << "Failed to get input track collection with label: " << inputLabel;
0232 LogTrace("MuonIdentification") << "Number of input outer tracks: " << outerTrackCollectionHandle_->size();
0233 } else if (inputType == ICTypes::LINKS) {
0234 iEvent.getByToken(linkCollectionToken_, linkCollectionHandle_);
0235 if (!linkCollectionHandle_.isValid())
0236 throw cms::Exception("FatalError") << "Failed to get input link collection with label: " << inputLabel;
0237 LogTrace("MuonIdentification") << "Number of input links: " << linkCollectionHandle_->size();
0238 } else if (inputType == ICTypes::MUONS) {
0239 iEvent.getByToken(muonCollectionToken_, muonCollectionHandle_);
0240 if (!muonCollectionHandle_.isValid())
0241 throw cms::Exception("FatalError") << "Failed to get input muon collection with label: " << inputLabel;
0242 LogTrace("MuonIdentification") << "Number of input muons: " << muonCollectionHandle_->size();
0243 } else if (fillGlobalTrackRefits_ && inputType == ICTypes::TEV_FIRSTHIT) {
0244 iEvent.getByToken(tpfmsCollectionToken_, tpfmsCollectionHandle_);
0245 if (!tpfmsCollectionHandle_.isValid())
0246 throw cms::Exception("FatalError") << "Failed to get input muon collection with label: " << inputLabel;
0247 LogTrace("MuonIdentification") << "Number of input muons: " << tpfmsCollectionHandle_->size();
0248 } else if (fillGlobalTrackRefits_ && inputType == ICTypes::TEV_PICKY) {
0249 iEvent.getByToken(pickyCollectionToken_, pickyCollectionHandle_);
0250 if (!pickyCollectionHandle_.isValid())
0251 throw cms::Exception("FatalError") << "Failed to get input muon collection with label: " << inputLabel;
0252 LogTrace("MuonIdentification") << "Number of input muons: " << pickyCollectionHandle_->size();
0253 } else if (fillGlobalTrackRefits_ && inputType == ICTypes::TEV_DYT) {
0254 iEvent.getByToken(dytCollectionToken_, dytCollectionHandle_);
0255 if (!dytCollectionHandle_.isValid())
0256 throw cms::Exception("FatalError") << "Failed to get input muon collection with label: " << inputLabel;
0257 LogTrace("MuonIdentification") << "Number of input muons: " << dytCollectionHandle_->size();
0258 } else
0259 throw cms::Exception("FatalError") << "Unknown input collection type: #" << ICTypes::toStr(inputType);
0260 }
0261
0262 iEvent.getByToken(rpcHitToken_, rpcHitHandle_);
0263 iEvent.getByToken(gemHitToken_, gemHitHandle_);
0264 if (fillGlobalTrackQuality_)
0265 iEvent.getByToken(glbQualToken_, glbQualHandle_);
0266 if (selectHighPurity_)
0267 iEvent.getByToken(pvToken_, pvHandle_);
0268 }
0269
0270 reco::Muon MuonIdProducer::makeMuon(edm::Event& iEvent,
0271 const edm::EventSetup& iSetup,
0272 const reco::TrackRef& track,
0273 MuonIdProducer::TrackType type) {
0274 LogTrace("MuonIdentification") << "Creating a muon from a track " << track.get()->pt()
0275 << " Pt (GeV), eta: " << track.get()->eta();
0276 reco::Muon aMuon(makeMuon(*(track.get())));
0277
0278 LogTrace("MuonIdentification") << "Muon created from a track ";
0279
0280 aMuon.setMuonTrack(type, track);
0281 aMuon.setBestTrack(type);
0282 aMuon.setTunePBestTrack(type);
0283
0284 LogTrace("MuonIdentification")
0285 << "Muon created from a track and setMuonBestTrack, setBestTrack and setTunePBestTrack called";
0286
0287 return aMuon;
0288 }
0289
0290 reco::CaloMuon MuonIdProducer::makeCaloMuon(const reco::Muon& muon) {
0291 LogTrace("MuonIdentification") << "Creating a CaloMuon from a Muon";
0292
0293 reco::CaloMuon aMuon;
0294 aMuon.setInnerTrack(muon.innerTrack());
0295
0296 if (muon.isEnergyValid())
0297 aMuon.setCalEnergy(muon.calEnergy());
0298
0299 if (fillCaloCompatibility_)
0300 aMuon.setCaloCompatibility(muonCaloCompatibility_.evaluate(muon));
0301 return aMuon;
0302 }
0303
0304 reco::Muon MuonIdProducer::makeMuon(const reco::MuonTrackLinks& links) {
0305 LogTrace("MuonIdentification") << "Creating a muon from a link to tracks object";
0306
0307 reco::Muon aMuon;
0308 reco::Muon::MuonTrackTypePair chosenTrack;
0309 reco::TrackRef tpfmsRef;
0310 reco::TrackRef pickyRef;
0311 reco::TrackRef dytRef;
0312 bool useSigmaSwitch = false;
0313
0314 if (tpfmsCollectionHandle_.isValid() && !tpfmsCollectionHandle_.failedToGet() && pickyCollectionHandle_.isValid() &&
0315 !pickyCollectionHandle_.failedToGet()) {
0316 tpfmsRef = muon::getTevRefitTrack(links.globalTrack(), *tpfmsCollectionHandle_);
0317 pickyRef = muon::getTevRefitTrack(links.globalTrack(), *pickyCollectionHandle_);
0318 dytRef = muon::getTevRefitTrack(links.globalTrack(), *dytCollectionHandle_);
0319
0320 if (tpfmsRef.isNull() && pickyRef.isNull() && dytRef.isNull()) {
0321 edm::LogWarning("MakeMuonWithTEV") << "Failed to get TEV refits, fall back to sigma switch.";
0322 useSigmaSwitch = true;
0323 }
0324 } else {
0325 useSigmaSwitch = true;
0326 }
0327
0328 if (useSigmaSwitch) {
0329 chosenTrack = muon::sigmaSwitch(links.globalTrack(),
0330 links.trackerTrack(),
0331 sigmaThresholdToFillCandidateP4WithGlobalFit_,
0332 ptThresholdToFillCandidateP4WithGlobalFit_);
0333 } else {
0334 chosenTrack = muon::tevOptimized(links.globalTrack(),
0335 links.trackerTrack(),
0336 tpfmsRef,
0337 pickyRef,
0338 dytRef,
0339 ptThresholdToFillCandidateP4WithGlobalFit_);
0340 }
0341 aMuon = makeMuon(*chosenTrack.first);
0342 aMuon.setInnerTrack(links.trackerTrack());
0343 aMuon.setOuterTrack(links.standAloneTrack());
0344 aMuon.setGlobalTrack(links.globalTrack());
0345 aMuon.setBestTrack(chosenTrack.second);
0346 aMuon.setTunePBestTrack(chosenTrack.second);
0347
0348 if (fillGlobalTrackRefits_) {
0349 if (tpfmsCollectionHandle_.isValid() && !tpfmsCollectionHandle_.failedToGet()) {
0350 reco::TrackToTrackMap::const_iterator it = tpfmsCollectionHandle_->find(links.globalTrack());
0351 if (it != tpfmsCollectionHandle_->end())
0352 aMuon.setMuonTrack(reco::Muon::TPFMS, (it->val));
0353 }
0354 if (pickyCollectionHandle_.isValid() && !pickyCollectionHandle_.failedToGet()) {
0355 reco::TrackToTrackMap::const_iterator it = pickyCollectionHandle_->find(links.globalTrack());
0356 if (it != pickyCollectionHandle_->end())
0357 aMuon.setMuonTrack(reco::Muon::Picky, (it->val));
0358 }
0359 if (dytCollectionHandle_.isValid() && !dytCollectionHandle_.failedToGet()) {
0360 reco::TrackToTrackMap::const_iterator it = dytCollectionHandle_->find(links.globalTrack());
0361 if (it != dytCollectionHandle_->end())
0362 aMuon.setMuonTrack(reco::Muon::DYT, (it->val));
0363 }
0364 }
0365 return aMuon;
0366 }
0367
0368 bool MuonIdProducer::isGoodTrack(const reco::Track& track) {
0369
0370 const double p = track.p();
0371 const double pt = track.pt();
0372 if (pt < minPt_ || (p < minP_ && p < minPCaloMuon_)) {
0373 LogTrace("MuonIdentification") << "Skipped low momentum track (Pt,P): " << pt << ", " << track.p() << " GeV";
0374 return false;
0375 }
0376
0377
0378 const double eta = track.eta();
0379 const double absEta = std::abs(eta);
0380 if (absEta > maxAbsEta_) {
0381 LogTrace("MuonIdentification") << "Skipped track with large pseudo rapidity (Eta: " << track.eta() << " )";
0382 return false;
0383 }
0384
0385 return true;
0386 }
0387
0388 unsigned int MuonIdProducer::chamberId(const DetId& id) {
0389 if (id.det() != DetId::Muon)
0390 return 0;
0391
0392 const auto subdetId = id.subdetId();
0393 if (subdetId == MuonSubdetId::DT) {
0394 return DTChamberId(id.rawId()).rawId();
0395 } else if (subdetId == MuonSubdetId::CSC) {
0396 return CSCDetId(id.rawId()).chamberId().rawId();
0397 }
0398
0399 return 0;
0400 }
0401
0402 int MuonIdProducer::overlap(const reco::Muon& muon, const reco::Track& track) {
0403 if (!muon.isMatchesValid() || track.extra().isNull() || track.extra()->recHitsSize() == 0)
0404 return 0;
0405
0406 int numberOfCommonDetIds = 0;
0407 const std::vector<reco::MuonChamberMatch>& matches(muon.matches());
0408 for (const auto& match : matches) {
0409 if (match.segmentMatches.empty())
0410 continue;
0411
0412 bool foundCommonDetId = false;
0413 for (auto hit = track.extra()->recHitsBegin(); hit != track.extra()->recHitsEnd(); ++hit) {
0414
0415
0416
0417
0418 if (chamberId((*hit)->geographicalId()) == match.id.rawId()) {
0419 foundCommonDetId = true;
0420 break;
0421 }
0422 }
0423 if (foundCommonDetId) {
0424 ++numberOfCommonDetIds;
0425 break;
0426 }
0427 }
0428 return numberOfCommonDetIds;
0429 }
0430
0431 void MuonIdProducer::beginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) {
0432 meshAlgo_->setCSCGeometry(&iSetup.getData(geomTokenRun_));
0433
0434 if (fillShowerDigis_ && fillMatching_)
0435 theShowerDigiFiller_->getES(iSetup);
0436 }
0437
0438 bool validateGlobalMuonPair(const reco::MuonTrackLinks& goodMuon, const reco::MuonTrackLinks& badMuon) {
0439 const int nHitsGood = goodMuon.globalTrack()->hitPattern().numberOfValidMuonHits();
0440 const int nHitsBad = badMuon.globalTrack()->hitPattern().numberOfValidMuonHits();
0441 if (std::min(nHitsGood, nHitsBad) > 10) {
0442 const double chi2Good = goodMuon.globalTrack()->normalizedChi2();
0443 const double chi2Bad = badMuon.globalTrack()->normalizedChi2();
0444 return (chi2Good <= chi2Bad);
0445 }
0446
0447 return (nHitsGood >= nHitsBad);
0448 }
0449
0450 void MuonIdProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0451 auto outputMuons = std::make_unique<reco::MuonCollection>();
0452 auto caloMuons = std::make_unique<reco::CaloMuonCollection>();
0453
0454 init(iEvent, iSetup);
0455
0456 if (fillShowerDigis_ && fillMatching_)
0457 theShowerDigiFiller_->getDigis(iEvent);
0458
0459
0460
0461
0462 if (muonCollectionHandle_.isValid()) {
0463 for (const auto& muon : *muonCollectionHandle_) {
0464 outputMuons->push_back(muon);
0465 }
0466 }
0467
0468
0469 if (linkCollectionHandle_.isValid()) {
0470 const auto nLink = linkCollectionHandle_->size();
0471 std::vector<bool> goodmuons(nLink, true);
0472 if (nLink > 1) {
0473
0474 for (unsigned int i = 0; i < nLink - 1; ++i) {
0475 const auto& iLink = linkCollectionHandle_->at(i);
0476 if (iLink.trackerTrack().isNull() || !checkLinks(&iLink))
0477 continue;
0478 for (unsigned int j = i + 1; j < nLink; ++j) {
0479 const auto& jLink = linkCollectionHandle_->at(j);
0480 if (!checkLinks(&jLink))
0481 continue;
0482 if (iLink.trackerTrack() == jLink.trackerTrack()) {
0483
0484
0485
0486
0487 if (validateGlobalMuonPair(iLink, jLink))
0488 goodmuons[j] = false;
0489 else
0490 goodmuons[i] = false;
0491 }
0492 }
0493 }
0494
0495 for (unsigned int i = 0; i < nLink - 1; ++i) {
0496 if (!goodmuons[i])
0497 continue;
0498 const auto& iLink = linkCollectionHandle_->at(i);
0499 if (iLink.standAloneTrack().isNull() || !checkLinks(&iLink))
0500 continue;
0501 for (unsigned int j = i + 1; j < nLink; ++j) {
0502 if (!goodmuons[j])
0503 continue;
0504 const auto& jLink = linkCollectionHandle_->at(j);
0505 if (!checkLinks(&jLink))
0506 continue;
0507 if (iLink.standAloneTrack() == jLink.standAloneTrack()) {
0508 if (validateGlobalMuonPair(iLink, jLink))
0509 goodmuons[j] = false;
0510 else
0511 goodmuons[i] = false;
0512 }
0513 }
0514 }
0515 }
0516 for (unsigned int i = 0; i < nLink; ++i) {
0517 if (!goodmuons[i])
0518 continue;
0519 const auto& iLink = linkCollectionHandle_->at(i);
0520 if (!checkLinks(&iLink))
0521 continue;
0522
0523 bool newMuon = true;
0524 for (const auto& muon : *outputMuons) {
0525 if (muon.track() == iLink.trackerTrack() && muon.standAloneMuon() == iLink.standAloneTrack() &&
0526 muon.combinedMuon() == iLink.globalTrack()) {
0527 newMuon = false;
0528 break;
0529 }
0530 }
0531 if (newMuon) {
0532 outputMuons->push_back(makeMuon(iLink));
0533 outputMuons->back().setType(reco::Muon::GlobalMuon | reco::Muon::StandAloneMuon);
0534 }
0535 }
0536 }
0537
0538
0539 if (innerTrackCollectionHandle_.isValid()) {
0540 LogTrace("MuonIdentification") << "Creating tracker muons";
0541 std::vector<TrackDetectorAssociator::Direction> directions1, directions2;
0542 directions1.push_back(TrackDetectorAssociator::InsideOut);
0543 directions1.push_back(TrackDetectorAssociator::OutsideIn);
0544 directions2.push_back(TrackDetectorAssociator::Any);
0545
0546 const GlobalTrackingGeometry* geometry = nullptr;
0547 if (debugWithTruthMatching_) {
0548 geometry = &iSetup.getData(globalGeomToken_);
0549 }
0550
0551 for (unsigned int i = 0; i < innerTrackCollectionHandle_->size(); ++i) {
0552 const reco::Track& track = innerTrackCollectionHandle_->at(i);
0553 if (!isGoodTrack(track))
0554 continue;
0555 if (selectHighPurity_ && !track.quality(reco::TrackBase::highPurity)) {
0556 const reco::VertexCollection* recoVertices = pvHandle_.product();
0557 if (!(*recoVertices)[0].isFake())
0558 continue;
0559 }
0560 const auto& trackRef = reco::TrackRef(innerTrackCollectionHandle_, i);
0561 bool splitTrack = false;
0562 if (track.extra().isAvailable() && TrackDetectorAssociator::crossedIP(track))
0563 splitTrack = true;
0564 const auto& directions = splitTrack ? directions1 : directions2;
0565 for (const auto direction : directions) {
0566
0567 reco::Muon trackerMuon(makeMuon(iEvent, iSetup, trackRef, reco::Muon::InnerTrack));
0568 fillMuonId(iEvent, iSetup, trackerMuon, direction);
0569
0570 if (debugWithTruthMatching_) {
0571
0572
0573 MuonIdTruthInfo::truthMatchMuon(iEvent, *geometry, trackerMuon);
0574 }
0575
0576
0577
0578
0579 bool newMuon = true;
0580 const bool goodTrackerMuon = isGoodTrackerMuon(trackerMuon);
0581 const bool goodRPCMuon = isGoodRPCMuon(trackerMuon);
0582 const bool goodGEMMuon = isGoodGEMMuon(trackerMuon);
0583 const bool goodME0Muon = isGoodME0Muon(trackerMuon);
0584 if (goodTrackerMuon)
0585 trackerMuon.setType(trackerMuon.type() | reco::Muon::TrackerMuon);
0586 if (goodRPCMuon)
0587 trackerMuon.setType(trackerMuon.type() | reco::Muon::RPCMuon);
0588 if (goodGEMMuon)
0589 trackerMuon.setType(trackerMuon.type() | reco::Muon::GEMMuon);
0590 if (goodME0Muon)
0591 trackerMuon.setType(trackerMuon.type() | reco::Muon::ME0Muon);
0592
0593 for (auto& muon : *outputMuons) {
0594 if (muon.innerTrack().get() == trackerMuon.innerTrack().get() &&
0595 std::abs(reco::deltaPhi(phiOfMuonInteractionRegion(muon), phiOfMuonInteractionRegion(trackerMuon))) <
0596 M_PI_2) {
0597 newMuon = false;
0598 muon.setMatches(trackerMuon.matches());
0599 if (trackerMuon.isTimeValid())
0600 muon.setTime(trackerMuon.time());
0601 if (trackerMuon.isEnergyValid())
0602 muon.setCalEnergy(trackerMuon.calEnergy());
0603 if (goodTrackerMuon)
0604 muon.setType(muon.type() | reco::Muon::TrackerMuon);
0605 if (goodRPCMuon)
0606 muon.setType(muon.type() | reco::Muon::RPCMuon);
0607 if (goodGEMMuon)
0608 muon.setType(muon.type() | reco::Muon::GEMMuon);
0609 if (goodME0Muon)
0610 muon.setType(muon.type() | reco::Muon::ME0Muon);
0611 LogTrace("MuonIdentification") << "Found a corresponding global muon. Set energy, matches and move on";
0612 break;
0613 }
0614 }
0615 if (newMuon) {
0616 if (goodTrackerMuon || goodRPCMuon || goodGEMMuon || goodME0Muon) {
0617 outputMuons->push_back(trackerMuon);
0618 } else {
0619 LogTrace("MuonIdentification") << "track failed minimal number of muon matches requirement";
0620 const reco::CaloMuon& caloMuon = makeCaloMuon(trackerMuon);
0621 if (isGoodCaloMuon(caloMuon))
0622 caloMuons->push_back(caloMuon);
0623 }
0624 }
0625 }
0626 }
0627 }
0628
0629
0630 if (outerTrackCollectionHandle_.isValid()) {
0631 LogTrace("MuonIdentification") << "Looking for new muons among stand alone muon tracks";
0632 for (unsigned int i = 0; i < outerTrackCollectionHandle_->size(); ++i) {
0633 const auto& outerTrack = outerTrackCollectionHandle_->at(i);
0634
0635
0636 bool newMuon = true;
0637 for (auto& muon : *outputMuons) {
0638 if (!muon.standAloneMuon().isNull()) {
0639
0640 if (muon.standAloneMuon().get() == &outerTrack ||
0641 (muon.standAloneMuon()->extra().isNonnull() &&
0642 muon.standAloneMuon()->extra().get() == outerTrack.extra().get())) {
0643 newMuon = false;
0644 break;
0645 }
0646 } else {
0647
0648
0649
0650
0651
0652
0653
0654 if (overlap(muon, outerTrack) > 0) {
0655 LogTrace("MuonIdentification") << "Found associated tracker muon. Set a reference and move on";
0656 newMuon = false;
0657 muon.setOuterTrack(reco::TrackRef(outerTrackCollectionHandle_, i));
0658 muon.setType(muon.type() | reco::Muon::StandAloneMuon);
0659 break;
0660 }
0661 }
0662 }
0663 if (newMuon) {
0664 LogTrace("MuonIdentification") << "No associated stand alone track is found. Making a muon";
0665 outputMuons->push_back(
0666 makeMuon(iEvent, iSetup, reco::TrackRef(outerTrackCollectionHandle_, i), reco::Muon::OuterTrack));
0667 outputMuons->back().setType(reco::Muon::StandAloneMuon);
0668 }
0669 }
0670 }
0671
0672 if (arbitrateTrackerMuons_) {
0673 fillArbitrationInfo(outputMuons.get());
0674 arbitrateMuons(outputMuons.get(), caloMuons.get());
0675 }
0676
0677 LogTrace("MuonIdentification") << "Dress up muons if it's necessary";
0678
0679 const int nMuons = outputMuons->size();
0680
0681 std::vector<reco::MuonTimeExtra> dtTimeColl(nMuons);
0682 std::vector<reco::MuonTimeExtra> cscTimeColl(nMuons);
0683 std::vector<reco::MuonTimeExtra> combinedTimeColl(nMuons);
0684 std::vector<reco::IsoDeposit> trackDepColl(nMuons);
0685 std::vector<reco::IsoDeposit> ecalDepColl(nMuons);
0686 std::vector<reco::IsoDeposit> hcalDepColl(nMuons);
0687 std::vector<reco::IsoDeposit> hoDepColl(nMuons);
0688 std::vector<reco::IsoDeposit> jetDepColl(nMuons);
0689
0690
0691 for (unsigned int i = 0; i < outputMuons->size(); ++i) {
0692 auto& muon = outputMuons->at(i);
0693
0694
0695 if ((fillMatching_ && !muon.isMatchesValid()) || (fillEnergy_ && !muon.isEnergyValid())) {
0696
0697
0698 if (muon.isStandAloneMuon()) {
0699 if (std::abs(reco::deltaPhi(phiOfMuonInteractionRegion(muon), muon.phi())) < M_PI_2) {
0700 fillMuonId(iEvent, iSetup, muon, TrackDetectorAssociator::InsideOut);
0701 } else {
0702 fillMuonId(iEvent, iSetup, muon, TrackDetectorAssociator::OutsideIn);
0703 }
0704 } else {
0705 LogTrace("MuonIdentification") << "THIS SHOULD NEVER HAPPEN";
0706 fillMuonId(iEvent, iSetup, muon);
0707 }
0708 }
0709
0710 if (fillGlobalTrackQuality_) {
0711
0712 fillGlbQuality(iEvent, iSetup, muon);
0713 }
0714 LogDebug("MuonIdentification") << "";
0715
0716 if (fillTrackerKink_) {
0717 fillTrackerKink(muon);
0718 }
0719
0720 if (fillCaloCompatibility_)
0721 muon.setCaloCompatibility(muonCaloCompatibility_.evaluate(muon));
0722
0723 if (fillIsolation_) {
0724 fillMuonIsolation(
0725 iEvent, iSetup, muon, trackDepColl[i], ecalDepColl[i], hcalDepColl[i], hoDepColl[i], jetDepColl[i]);
0726 }
0727
0728
0729 reco::MuonTime muonTime;
0730 reco::MuonTimeExtra dtTime;
0731 reco::MuonTimeExtra cscTime;
0732 reco::MuonTime rpcTime;
0733 reco::MuonTimeExtra combinedTime;
0734
0735 theTimingFiller_->fillTiming(muon, dtTime, cscTime, rpcTime, combinedTime, iEvent, iSetup);
0736
0737 muonTime.nDof = combinedTime.nDof();
0738 muonTime.timeAtIpInOut = combinedTime.timeAtIpInOut();
0739 muonTime.timeAtIpInOutErr = combinedTime.timeAtIpInOutErr();
0740 muonTime.timeAtIpOutIn = combinedTime.timeAtIpOutIn();
0741 muonTime.timeAtIpOutInErr = combinedTime.timeAtIpOutInErr();
0742
0743 muon.setTime(muonTime);
0744 muon.setRPCTime(rpcTime);
0745 dtTimeColl[i] = dtTime;
0746 cscTimeColl[i] = cscTime;
0747 combinedTimeColl[i] = combinedTime;
0748 }
0749
0750 LogTrace("MuonIdentification") << "number of muons produced: " << outputMuons->size();
0751 if (fillMatching_) {
0752 fillArbitrationInfo(outputMuons.get(), reco::Muon::TrackerMuon);
0753 fillArbitrationInfo(outputMuons.get(), reco::Muon::ME0Muon);
0754 fillArbitrationInfo(outputMuons.get(), reco::Muon::GEMMuon);
0755 }
0756 edm::OrphanHandle<reco::MuonCollection> muonHandle = iEvent.put(std::move(outputMuons));
0757
0758 auto fillMap = [](auto refH, auto& vec, edm::Event& ev, const std::string& cAl = "") {
0759 typedef edm::ValueMap<typename std::decay<decltype(vec)>::type::value_type> MapType;
0760 auto oMap = std::make_unique<MapType>();
0761 {
0762 typename MapType::Filler filler(*oMap);
0763 filler.insert(refH, vec.begin(), vec.end());
0764 vec.clear();
0765 filler.fill();
0766 }
0767 ev.put(std::move(oMap), cAl);
0768 };
0769 fillMap(muonHandle, combinedTimeColl, iEvent, "combined");
0770 fillMap(muonHandle, dtTimeColl, iEvent, "dt");
0771 fillMap(muonHandle, cscTimeColl, iEvent, "csc");
0772
0773 if (writeIsoDeposits_ && fillIsolation_) {
0774 fillMap(muonHandle, trackDepColl, iEvent, trackDepositName_);
0775 fillMap(muonHandle, ecalDepColl, iEvent, ecalDepositName_);
0776 fillMap(muonHandle, hcalDepColl, iEvent, hcalDepositName_);
0777 fillMap(muonHandle, hoDepColl, iEvent, hoDepositName_);
0778 fillMap(muonHandle, jetDepColl, iEvent, jetDepositName_);
0779 }
0780
0781 iEvent.put(std::move(caloMuons));
0782 }
0783
0784 bool MuonIdProducer::isGoodTrackerMuon(const reco::Muon& muon) {
0785 if (muon.track()->pt() < minPt_ || muon.track()->p() < minP_)
0786 return false;
0787 if (addExtraSoftMuons_ && muon.pt() < 5 && std::abs(muon.eta()) < 1.5 &&
0788 muon.numberOfMatches(reco::Muon::NoArbitration) >= 1)
0789 return true;
0790 return (muon.numberOfMatches(reco::Muon::NoArbitration) >= minNumberOfMatches_);
0791 }
0792
0793 bool MuonIdProducer::isGoodCaloMuon(const reco::CaloMuon& caloMuon) {
0794 if (!caloMuon.isCaloCompatibilityValid() || caloMuon.caloCompatibility() < caloCut_ || caloMuon.p() < minPCaloMuon_)
0795 return false;
0796 return true;
0797 }
0798
0799 bool MuonIdProducer::isGoodRPCMuon(const reco::Muon& muon) {
0800 if (muon.track()->pt() < minPt_ || muon.track()->p() < minP_)
0801 return false;
0802 if (addExtraSoftMuons_ && muon.pt() < 5 && std::abs(muon.eta()) < 1.5 &&
0803 muon.numberOfMatchedRPCLayers(reco::Muon::RPCHitAndTrackArbitration) > 1)
0804 return true;
0805 return (muon.numberOfMatchedRPCLayers(reco::Muon::RPCHitAndTrackArbitration) > minNumberOfMatches_);
0806 }
0807
0808 bool MuonIdProducer::isGoodGEMMuon(const reco::Muon& muon) {
0809
0810 if (!isGoodTrackerMuon(muon))
0811 return false;
0812 if (muon.track()->pt() < minPt_ || muon.track()->p() < minP_)
0813 return false;
0814 return (muon.numberOfMatches(reco::Muon::GEMSegmentAndTrackArbitration) +
0815 muon.numberOfMatches(reco::Muon::GEMHitAndTrackArbitration)) >= 1;
0816 }
0817
0818 bool MuonIdProducer::isGoodME0Muon(const reco::Muon& muon) {
0819
0820 if (muon.track()->p() < minP_)
0821 return false;
0822 return (muon.numberOfMatches(reco::Muon::ME0SegmentAndTrackArbitration) >= 1);
0823 }
0824
0825 void MuonIdProducer::fillMuonId(edm::Event& iEvent,
0826 const edm::EventSetup& iSetup,
0827 reco::Muon& aMuon,
0828 TrackDetectorAssociator::Direction direction) {
0829 LogTrace("MuonIdentification") << "RecoMuon/MuonIdProducer :: fillMuonId";
0830
0831
0832 const reco::Track* track = nullptr;
0833 if (aMuon.track().isNonnull())
0834 track = aMuon.track().get();
0835 else if (aMuon.standAloneMuon().isNonnull())
0836 track = aMuon.standAloneMuon().get();
0837 else
0838 throw cms::Exception("FatalError")
0839 << "Failed to fill muon id information for a muon with undefined references to tracks";
0840
0841 TrackDetMatchInfo info = trackAssociator_.associate(iEvent, iSetup, *track, parameters_, direction);
0842
0843 LogTrace("MuonIdentification") << "RecoMuon/MuonIdProducer :: fillMuonId :: fillEnergy = " << fillEnergy_;
0844
0845 if (fillEnergy_) {
0846 reco::MuonEnergy muonEnergy;
0847 muonEnergy.em = info.crossedEnergy(TrackDetMatchInfo::EcalRecHits);
0848 muonEnergy.had = info.crossedEnergy(TrackDetMatchInfo::HcalRecHits);
0849 muonEnergy.ho = info.crossedEnergy(TrackDetMatchInfo::HORecHits);
0850 muonEnergy.tower = info.crossedEnergy(TrackDetMatchInfo::TowerTotal);
0851 muonEnergy.emS9 = info.nXnEnergy(TrackDetMatchInfo::EcalRecHits, 1);
0852 muonEnergy.emS25 = info.nXnEnergy(TrackDetMatchInfo::EcalRecHits, 2);
0853 muonEnergy.hadS9 = info.nXnEnergy(TrackDetMatchInfo::HcalRecHits, 1);
0854 muonEnergy.hoS9 = info.nXnEnergy(TrackDetMatchInfo::HORecHits, 1);
0855 muonEnergy.towerS9 = info.nXnEnergy(TrackDetMatchInfo::TowerTotal, 1);
0856 if (storeCrossedHcalRecHits_) {
0857 muonEnergy.crossedHadRecHits.clear();
0858 for (auto hit : info.crossedHcalRecHits) {
0859 reco::HcalMuonRecHit mhit;
0860 mhit.energy = hit->energy();
0861 mhit.chi2 = hit->chi2();
0862 mhit.time = hit->time();
0863 mhit.detId = hit->id();
0864 muonEnergy.crossedHadRecHits.push_back(mhit);
0865 }
0866 }
0867 muonEnergy.ecal_position = info.trkGlobPosAtEcal;
0868 muonEnergy.hcal_position = info.trkGlobPosAtHcal;
0869 if (!info.crossedEcalIds.empty())
0870 muonEnergy.ecal_id = info.crossedEcalIds.front();
0871 if (!info.crossedHcalIds.empty())
0872 muonEnergy.hcal_id = info.crossedHcalIds.front();
0873
0874 DetId emMaxId = info.findMaxDeposition(TrackDetMatchInfo::EcalRecHits, 2);
0875 for (const auto& hit : info.ecalRecHits) {
0876 if (hit->id() != emMaxId)
0877 continue;
0878 muonEnergy.emMax = hit->energy();
0879 muonEnergy.ecal_time = hit->time();
0880 }
0881 DetId hadMaxId = info.findMaxDeposition(TrackDetMatchInfo::HcalRecHits, 1);
0882 for (const auto& hit : info.hcalRecHits) {
0883 if (hit->id() != hadMaxId)
0884 continue;
0885 muonEnergy.hadMax = hit->energy();
0886 muonEnergy.hcal_time = hit->time();
0887 }
0888 aMuon.setCalEnergy(muonEnergy);
0889 }
0890 if (!fillMatching_ && !aMuon.isTrackerMuon() && !aMuon.isRPCMuon() && !aMuon.isGEMMuon())
0891 return;
0892
0893
0894 LogTrace("MuonIdentification") << "RecoMuon/MuonIdProducer :: fillMuonId :: fill muon match info ";
0895 std::vector<reco::MuonChamberMatch> muonChamberMatches;
0896 unsigned int nubmerOfMatchesAccordingToTrackAssociator = 0;
0897 for (const auto& chamber : info.chambers) {
0898 if (chamber.id.subdetId() == MuonSubdetId::RPC && rpcHitHandle_.isValid())
0899 continue;
0900 if (chamber.id.subdetId() == MuonSubdetId::GEM && gemHitHandle_.isValid() &&
0901 GEMDetId(chamber.id.rawId()).station() != 0)
0902 continue;
0903 reco::MuonChamberMatch matchedChamber;
0904
0905 const auto& lErr = chamber.tState.localError();
0906 const auto& lPos = chamber.tState.localPosition();
0907 const auto& lDir = chamber.tState.localDirection();
0908
0909 const auto& localError = lErr.positionError();
0910 matchedChamber.x = lPos.x();
0911 matchedChamber.y = lPos.y();
0912 matchedChamber.xErr = sqrt(localError.xx());
0913 matchedChamber.yErr = sqrt(localError.yy());
0914
0915 matchedChamber.dXdZ = lDir.z() != 0 ? lDir.x() / lDir.z() : 9999;
0916 matchedChamber.dYdZ = lDir.z() != 0 ? lDir.y() / lDir.z() : 9999;
0917
0918 AlgebraicSymMatrix55 trajectoryCovMatrix = lErr.matrix();
0919 matchedChamber.dXdZErr = trajectoryCovMatrix(1, 1) > 0 ? sqrt(trajectoryCovMatrix(1, 1)) : 0;
0920 matchedChamber.dYdZErr = trajectoryCovMatrix(2, 2) > 0 ? sqrt(trajectoryCovMatrix(2, 2)) : 0;
0921
0922 matchedChamber.edgeX = chamber.localDistanceX;
0923 matchedChamber.edgeY = chamber.localDistanceY;
0924
0925 matchedChamber.id = chamber.id;
0926
0927 if (fillShowerDigis_ && fillMatching_) {
0928 theShowerDigiFiller_->fill(matchedChamber);
0929 } else {
0930 theShowerDigiFiller_->fillDefault(matchedChamber);
0931 }
0932
0933 if (!chamber.segments.empty())
0934 ++nubmerOfMatchesAccordingToTrackAssociator;
0935
0936
0937 for (const auto& segment : chamber.segments) {
0938 reco::MuonSegmentMatch matchedSegment;
0939 matchedSegment.x = segment.segmentLocalPosition.x();
0940 matchedSegment.y = segment.segmentLocalPosition.y();
0941 matchedSegment.dXdZ =
0942 segment.segmentLocalDirection.z() ? segment.segmentLocalDirection.x() / segment.segmentLocalDirection.z() : 0;
0943 matchedSegment.dYdZ =
0944 segment.segmentLocalDirection.z() ? segment.segmentLocalDirection.y() / segment.segmentLocalDirection.z() : 0;
0945 matchedSegment.xErr = segment.segmentLocalErrorXX > 0 ? sqrt(segment.segmentLocalErrorXX) : 0;
0946 matchedSegment.yErr = segment.segmentLocalErrorYY > 0 ? sqrt(segment.segmentLocalErrorYY) : 0;
0947 matchedSegment.dXdZErr = segment.segmentLocalErrorDxDz > 0 ? sqrt(segment.segmentLocalErrorDxDz) : 0;
0948 matchedSegment.dYdZErr = segment.segmentLocalErrorDyDz > 0 ? sqrt(segment.segmentLocalErrorDyDz) : 0;
0949 matchedSegment.t0 = segment.t0;
0950 matchedSegment.mask = 0;
0951 matchedSegment.dtSegmentRef = segment.dtSegmentRef;
0952 matchedSegment.cscSegmentRef = segment.cscSegmentRef;
0953 matchedSegment.gemSegmentRef = segment.gemSegmentRef;
0954 matchedSegment.me0SegmentRef = segment.me0SegmentRef;
0955 matchedSegment.hasZed_ = segment.hasZed;
0956 matchedSegment.hasPhi_ = segment.hasPhi;
0957
0958 bool matchedX = false;
0959 bool matchedY = false;
0960 LogTrace("MuonIdentification") << " matching local x, segment x: " << matchedSegment.x
0961 << ", chamber x: " << matchedChamber.x << ", max: " << maxAbsDx_;
0962 LogTrace("MuonIdentification") << " matching local y, segment y: " << matchedSegment.y
0963 << ", chamber y: " << matchedChamber.y << ", max: " << maxAbsDy_;
0964 const double matchedSegChDx = std::abs(matchedSegment.x - matchedChamber.x);
0965 const double matchedSegChDy = std::abs(matchedSegment.y - matchedChamber.y);
0966 if (matchedSegment.xErr > 0 && matchedChamber.xErr > 0)
0967 LogTrace("MuonIdentification") << " xpull: "
0968 << matchedSegChDx / std::sqrt(std::pow(matchedSegment.xErr, 2) +
0969 std::pow(matchedChamber.xErr, 2));
0970 if (matchedSegment.yErr > 0 && matchedChamber.yErr > 0)
0971 LogTrace("MuonIdentification") << " ypull: "
0972 << matchedSegChDy / std::sqrt(std::pow(matchedSegment.yErr, 2) +
0973 std::pow(matchedChamber.yErr, 2));
0974
0975 if (matchedSegChDx < maxAbsDx_)
0976 matchedX = true;
0977 else if (matchedSegment.xErr > 0 && matchedChamber.xErr > 0) {
0978 const double invMatchedSegChPullX2 = std::pow(matchedSegment.xErr, 2) + std::pow(matchedChamber.xErr, 2);
0979 if (matchedSegChDx * matchedSegChDx < maxAbsPullX2_ * invMatchedSegChPullX2)
0980 matchedX = true;
0981 }
0982 if (matchedSegChDy < maxAbsDy_)
0983 matchedY = true;
0984 else if (matchedSegment.yErr > 0 && matchedChamber.yErr > 0) {
0985 const double invMatchedSegChPullY2 = std::pow(matchedSegment.yErr, 2) + std::pow(matchedChamber.yErr, 2);
0986 if (matchedSegChDy * matchedSegChDy < maxAbsPullY2_ * invMatchedSegChPullY2)
0987 matchedY = true;
0988 }
0989 if (matchedX && matchedY) {
0990 if (matchedChamber.id.subdetId() == MuonSubdetId::ME0)
0991 matchedChamber.me0Matches.push_back(matchedSegment);
0992 else if (matchedChamber.id.subdetId() == MuonSubdetId::GEM)
0993 matchedChamber.gemMatches.push_back(matchedSegment);
0994 else
0995 matchedChamber.segmentMatches.push_back(matchedSegment);
0996 }
0997 }
0998 muonChamberMatches.push_back(matchedChamber);
0999 }
1000
1001
1002 LogTrace("MuonIdentification") << "RecoMuon/MuonIdProducer :: fillMuonId :: fill RPC info";
1003 if (rpcHitHandle_.isValid()) {
1004 for (const auto& chamber : info.chambers) {
1005 if (chamber.id.subdetId() != MuonSubdetId::RPC)
1006 continue;
1007 const auto& lErr = chamber.tState.localError();
1008 const auto& lPos = chamber.tState.localPosition();
1009 const auto& lDir = chamber.tState.localDirection();
1010
1011 reco::MuonChamberMatch matchedChamber;
1012
1013 LocalError localError = lErr.positionError();
1014 matchedChamber.x = lPos.x();
1015 matchedChamber.y = lPos.y();
1016 matchedChamber.xErr = sqrt(localError.xx());
1017 matchedChamber.yErr = sqrt(localError.yy());
1018
1019 matchedChamber.dXdZ = lDir.z() != 0 ? lDir.x() / lDir.z() : 9999;
1020 matchedChamber.dYdZ = lDir.z() != 0 ? lDir.y() / lDir.z() : 9999;
1021
1022 AlgebraicSymMatrix55 trajectoryCovMatrix = lErr.matrix();
1023 matchedChamber.dXdZErr = trajectoryCovMatrix(1, 1) > 0 ? sqrt(trajectoryCovMatrix(1, 1)) : 0;
1024 matchedChamber.dYdZErr = trajectoryCovMatrix(2, 2) > 0 ? sqrt(trajectoryCovMatrix(2, 2)) : 0;
1025
1026 matchedChamber.edgeX = chamber.localDistanceX;
1027 matchedChamber.edgeY = chamber.localDistanceY;
1028
1029 theShowerDigiFiller_->fillDefault(matchedChamber);
1030
1031 matchedChamber.id = chamber.id;
1032
1033 for (const auto& rpcRecHit : *rpcHitHandle_) {
1034 reco::MuonRPCHitMatch rpcHitMatch;
1035
1036 if (rpcRecHit.rawId() != chamber.id.rawId())
1037 continue;
1038
1039 rpcHitMatch.x = rpcRecHit.localPosition().x();
1040 rpcHitMatch.mask = 0;
1041 rpcHitMatch.bx = rpcRecHit.BunchX();
1042
1043 const double absDx = std::abs(rpcRecHit.localPosition().x() - chamber.tState.localPosition().x());
1044 if (absDx <= 20 or absDx * absDx <= 16 * localError.xx())
1045 matchedChamber.rpcMatches.push_back(rpcHitMatch);
1046 }
1047
1048 muonChamberMatches.push_back(matchedChamber);
1049 }
1050 }
1051
1052
1053 LogTrace("MuonIdentification") << "RecoMuon/MuonIdProducer :: fillMuonId :: fill GEM info";
1054 if (gemHitHandle_.isValid()) {
1055 for (const auto& chamber : info.chambers) {
1056
1057 if (chamber.id.subdetId() != MuonSubdetId::GEM || GEMDetId(chamber.id.rawId()).station() == 0 ||
1058 GEMDetId(chamber.id.rawId()).layer() == 0)
1059 continue;
1060 const auto& lErr = chamber.tState.localError();
1061 const auto& lPos = chamber.tState.localPosition();
1062 const auto& lDir = chamber.tState.localDirection();
1063
1064 reco::MuonChamberMatch matchedChamber;
1065
1066 LocalError localError = lErr.positionError();
1067 matchedChamber.x = lPos.x();
1068 matchedChamber.y = lPos.y();
1069 matchedChamber.xErr = sqrt(localError.xx());
1070 matchedChamber.yErr = sqrt(localError.yy());
1071
1072 matchedChamber.dXdZ = lDir.z() != 0 ? lDir.x() / lDir.z() : 9999;
1073 matchedChamber.dYdZ = lDir.z() != 0 ? lDir.y() / lDir.z() : 9999;
1074
1075 AlgebraicSymMatrix55 trajectoryCovMatrix = lErr.matrix();
1076 matchedChamber.dXdZErr = trajectoryCovMatrix(1, 1) > 0 ? sqrt(trajectoryCovMatrix(1, 1)) : 0;
1077 matchedChamber.dYdZErr = trajectoryCovMatrix(2, 2) > 0 ? sqrt(trajectoryCovMatrix(2, 2)) : 0;
1078
1079 matchedChamber.edgeX = chamber.localDistanceX;
1080 matchedChamber.edgeY = chamber.localDistanceY;
1081
1082 theShowerDigiFiller_->fillDefault(matchedChamber);
1083
1084 matchedChamber.id = chamber.id;
1085
1086 for (const auto& gemRecHit : *gemHitHandle_) {
1087 reco::MuonGEMHitMatch gemHitMatch;
1088
1089 if (GEMDetId(gemRecHit.gemId().region(),
1090 gemRecHit.gemId().ring(),
1091 gemRecHit.gemId().station(),
1092 gemRecHit.gemId().layer(),
1093 gemRecHit.gemId().chamber(),
1094 0)
1095 .rawId() != chamber.id.rawId())
1096 continue;
1097
1098 gemHitMatch.x = gemRecHit.localPosition().x();
1099 gemHitMatch.mask = 0;
1100 gemHitMatch.bx = gemRecHit.BunchX();
1101
1102 const double absDx = std::abs(gemRecHit.localPosition().x() - chamber.tState.localPosition().x());
1103 if (absDx <= 5 or absDx * absDx <= 16 * localError.xx())
1104 matchedChamber.gemHitMatches.push_back(gemHitMatch);
1105 }
1106
1107 muonChamberMatches.push_back(matchedChamber);
1108 }
1109 }
1110
1111 aMuon.setMatches(muonChamberMatches);
1112
1113 LogTrace("MuonIdentification") << "number of muon chambers: " << aMuon.matches().size() << "\n"
1114 << "number of chambers with segments according to the associator requirements: "
1115 << nubmerOfMatchesAccordingToTrackAssociator;
1116 LogTrace("MuonIdentification") << "number of segment matches with the producer requirements: "
1117 << aMuon.numberOfMatches(reco::Muon::NoArbitration);
1118
1119
1120 }
1121
1122 void MuonIdProducer::arbitrateMuons(reco::MuonCollection* muons, reco::CaloMuonCollection* caloMuons) {
1123 reco::Muon::ArbitrationType arbitration = reco::Muon::SegmentAndTrackArbitration;
1124
1125
1126 for (reco::MuonCollection::iterator muon = muons->begin(); muon != muons->end();) {
1127 if (muon->isTrackerMuon()) {
1128 if (muon->numberOfMatches(arbitration) < minNumberOfMatches_) {
1129
1130
1131
1132
1133 unsigned int mask = reco::Muon::TrackerMuon | reco::Muon::PFMuon | reco::Muon::GEMMuon;
1134 if ((muon->type() & (~mask)) == 0) {
1135 const reco::CaloMuon& caloMuon = makeCaloMuon(*muon);
1136 if (isGoodCaloMuon(caloMuon))
1137 caloMuons->push_back(caloMuon);
1138 muon = muons->erase(muon);
1139 continue;
1140 } else {
1141 muon->setType(muon->type() & (~(reco::Muon::TrackerMuon | reco::Muon::GEMMuon)));
1142 }
1143 }
1144 }
1145 muon++;
1146 }
1147 }
1148
1149 void MuonIdProducer::fillArbitrationInfo(reco::MuonCollection* pOutputMuons, unsigned int muonType) {
1150
1151
1152
1153 std::vector<std::pair<reco::MuonChamberMatch*, reco::MuonSegmentMatch*> > chamberPairs;
1154 std::vector<std::pair<reco::MuonChamberMatch*, reco::MuonSegmentMatch*> > stationPairs;
1155 std::vector<std::pair<reco::MuonChamberMatch*, reco::MuonSegmentMatch*> >
1156 arbitrationPairs;
1157
1158
1159 for (unsigned int muonIndex1 = 0; muonIndex1 < pOutputMuons->size(); ++muonIndex1) {
1160 auto& muon1 = pOutputMuons->at(muonIndex1);
1161
1162 for (auto& chamber1 : muon1.matches()) {
1163
1164 std::vector<reco::MuonSegmentMatch>* segmentMatches1 = getSegmentMatches(chamber1, muonType);
1165
1166 if (segmentMatches1->empty())
1167 continue;
1168 chamberPairs.clear();
1169
1170 for (auto& segment1 : *segmentMatches1) {
1171 chamberPairs.push_back(std::make_pair(&chamber1, &segment1));
1172 if (!segment1.isMask())
1173 {
1174 arbitrationPairs.clear();
1175 arbitrationPairs.push_back(std::make_pair(&chamber1, &segment1));
1176
1177
1178
1179 if (muon1.type() & muonType) {
1180
1181 for (unsigned int muonIndex2 = muonIndex1 + 1; muonIndex2 < pOutputMuons->size(); ++muonIndex2) {
1182 auto& muon2 = pOutputMuons->at(muonIndex2);
1183
1184 if (!(muon2.type() & muonType))
1185 continue;
1186
1187 for (auto& chamber2 : muon2.matches()) {
1188
1189 std::vector<reco::MuonSegmentMatch>* segmentMatches2 = getSegmentMatches(chamber2, muonType);
1190 for (auto& segment2 : *segmentMatches2) {
1191 if (segment2.isMask())
1192 continue;
1193 if (approxEqual(segment2.x, segment1.x) && approxEqual(segment2.y, segment1.y) &&
1194 approxEqual(segment2.dXdZ, segment1.dXdZ) && approxEqual(segment2.dYdZ, segment1.dYdZ) &&
1195 approxEqual(segment2.xErr, segment1.xErr) && approxEqual(segment2.yErr, segment1.yErr) &&
1196 approxEqual(segment2.dXdZErr, segment1.dXdZErr) &&
1197 approxEqual(segment2.dYdZErr, segment1.dYdZErr)) {
1198 arbitrationPairs.push_back(std::make_pair(&chamber2, &segment2));
1199 }
1200 }
1201 }
1202 }
1203 }
1204
1205
1206 if (arbitrationPairs.empty())
1207 continue;
1208 if (arbitrationPairs.size() == 1) {
1209 arbitrationPairs.front().second->setMask(reco::MuonSegmentMatch::BelongsToTrackByDRSlope);
1210 arbitrationPairs.front().second->setMask(reco::MuonSegmentMatch::BelongsToTrackByDXSlope);
1211 arbitrationPairs.front().second->setMask(reco::MuonSegmentMatch::BelongsToTrackByDR);
1212 arbitrationPairs.front().second->setMask(reco::MuonSegmentMatch::BelongsToTrackByDX);
1213 arbitrationPairs.front().second->setMask(reco::MuonSegmentMatch::Arbitrated);
1214 } else {
1215 sort(arbitrationPairs.begin(),
1216 arbitrationPairs.end(),
1217 SortMuonSegmentMatches(reco::MuonSegmentMatch::BelongsToTrackByDRSlope));
1218 arbitrationPairs.front().second->setMask(reco::MuonSegmentMatch::BelongsToTrackByDRSlope);
1219 sort(arbitrationPairs.begin(),
1220 arbitrationPairs.end(),
1221 SortMuonSegmentMatches(reco::MuonSegmentMatch::BelongsToTrackByDXSlope));
1222 arbitrationPairs.front().second->setMask(reco::MuonSegmentMatch::BelongsToTrackByDXSlope);
1223 sort(arbitrationPairs.begin(),
1224 arbitrationPairs.end(),
1225 SortMuonSegmentMatches(reco::MuonSegmentMatch::BelongsToTrackByDR));
1226 arbitrationPairs.front().second->setMask(reco::MuonSegmentMatch::BelongsToTrackByDR);
1227 sort(arbitrationPairs.begin(),
1228 arbitrationPairs.end(),
1229 SortMuonSegmentMatches(reco::MuonSegmentMatch::BelongsToTrackByDX));
1230 arbitrationPairs.front().second->setMask(reco::MuonSegmentMatch::BelongsToTrackByDX);
1231 for (auto& ap : arbitrationPairs) {
1232 ap.second->setMask(reco::MuonSegmentMatch::Arbitrated);
1233 }
1234 }
1235 }
1236
1237
1238 if (muonType == reco::Muon::TrackerMuon && chamber1.id.subdetId() == MuonSubdetId::CSC && arbClean_ &&
1239 CSCDetId(chamber1.id).ring() == 4) {
1240 for (auto& segment2 : chamber1.segmentMatches) {
1241 if (segment1.cscSegmentRef.isNull() || segment2.cscSegmentRef.isNull())
1242 continue;
1243 if (meshAlgo_->isDuplicateOf(segment1.cscSegmentRef, segment2.cscSegmentRef) &&
1244 (segment2.mask & 0x1e0000) && (segment1.mask & 0x1e0000)) {
1245 segment2.setMask(reco::MuonSegmentMatch::BelongsToTrackByME1aClean);
1246
1247 }
1248 }
1249 }
1250
1251 }
1252
1253
1254 if (chamberPairs.empty())
1255 continue;
1256 if (chamberPairs.size() == 1) {
1257 chamberPairs.front().second->setMask(reco::MuonSegmentMatch::BestInChamberByDRSlope);
1258 chamberPairs.front().second->setMask(reco::MuonSegmentMatch::BestInChamberByDXSlope);
1259 chamberPairs.front().second->setMask(reco::MuonSegmentMatch::BestInChamberByDR);
1260 chamberPairs.front().second->setMask(reco::MuonSegmentMatch::BestInChamberByDX);
1261 } else {
1262 sort(chamberPairs.begin(),
1263 chamberPairs.end(),
1264 SortMuonSegmentMatches(reco::MuonSegmentMatch::BestInChamberByDRSlope));
1265 chamberPairs.front().second->setMask(reco::MuonSegmentMatch::BestInChamberByDRSlope);
1266 sort(chamberPairs.begin(),
1267 chamberPairs.end(),
1268 SortMuonSegmentMatches(reco::MuonSegmentMatch::BestInChamberByDXSlope));
1269 chamberPairs.front().second->setMask(reco::MuonSegmentMatch::BestInChamberByDXSlope);
1270 sort(chamberPairs.begin(),
1271 chamberPairs.end(),
1272 SortMuonSegmentMatches(reco::MuonSegmentMatch::BestInChamberByDR));
1273 chamberPairs.front().second->setMask(reco::MuonSegmentMatch::BestInChamberByDR);
1274 sort(chamberPairs.begin(),
1275 chamberPairs.end(),
1276 SortMuonSegmentMatches(reco::MuonSegmentMatch::BestInChamberByDX));
1277 chamberPairs.front().second->setMask(reco::MuonSegmentMatch::BestInChamberByDX);
1278 }
1279 }
1280
1281
1282 for (int stationIndex = 1; stationIndex < 5; ++stationIndex) {
1283 for (int detectorIndex = 1; detectorIndex <= 5;
1284 ++detectorIndex)
1285 {
1286 stationPairs.clear();
1287
1288
1289 for (auto& chamber : muon1.matches()) {
1290 if (!(chamber.station() == stationIndex && chamber.detector() == detectorIndex))
1291 continue;
1292 std::vector<reco::MuonSegmentMatch>* segmentMatches = getSegmentMatches(chamber, muonType);
1293 if (segmentMatches->empty())
1294 continue;
1295
1296 for (auto& segment : *segmentMatches) {
1297 stationPairs.push_back(std::make_pair(&chamber, &segment));
1298 }
1299 }
1300
1301 if (stationPairs.empty())
1302 continue;
1303 if (stationPairs.size() == 1) {
1304 stationPairs.front().second->setMask(reco::MuonSegmentMatch::BestInStationByDRSlope);
1305 stationPairs.front().second->setMask(reco::MuonSegmentMatch::BestInStationByDXSlope);
1306 stationPairs.front().second->setMask(reco::MuonSegmentMatch::BestInStationByDR);
1307 stationPairs.front().second->setMask(reco::MuonSegmentMatch::BestInStationByDX);
1308 } else {
1309 sort(stationPairs.begin(),
1310 stationPairs.end(),
1311 SortMuonSegmentMatches(reco::MuonSegmentMatch::BestInStationByDRSlope));
1312 stationPairs.front().second->setMask(reco::MuonSegmentMatch::BestInStationByDRSlope);
1313 sort(stationPairs.begin(),
1314 stationPairs.end(),
1315 SortMuonSegmentMatches(reco::MuonSegmentMatch::BestInStationByDXSlope));
1316 stationPairs.front().second->setMask(reco::MuonSegmentMatch::BestInStationByDXSlope);
1317 sort(stationPairs.begin(),
1318 stationPairs.end(),
1319 SortMuonSegmentMatches(reco::MuonSegmentMatch::BestInStationByDR));
1320 stationPairs.front().second->setMask(reco::MuonSegmentMatch::BestInStationByDR);
1321 sort(stationPairs.begin(),
1322 stationPairs.end(),
1323 SortMuonSegmentMatches(reco::MuonSegmentMatch::BestInStationByDX));
1324 stationPairs.front().second->setMask(reco::MuonSegmentMatch::BestInStationByDX);
1325 }
1326 }
1327 }
1328
1329 }
1330
1331 if (arbClean_) {
1332
1333 meshAlgo_->clearMesh();
1334 meshAlgo_->runMesh(pOutputMuons);
1335 }
1336 }
1337
1338 void MuonIdProducer::fillMuonIsolation(edm::Event& iEvent,
1339 const edm::EventSetup& iSetup,
1340 reco::Muon& aMuon,
1341 reco::IsoDeposit& trackDep,
1342 reco::IsoDeposit& ecalDep,
1343 reco::IsoDeposit& hcalDep,
1344 reco::IsoDeposit& hoDep,
1345 reco::IsoDeposit& jetDep) {
1346 const reco::Track* track = nullptr;
1347 if (aMuon.track().isNonnull())
1348 track = aMuon.track().get();
1349 else if (aMuon.standAloneMuon().isNonnull())
1350 track = aMuon.standAloneMuon().get();
1351 else
1352 throw cms::Exception("FatalError")
1353 << "Failed to compute muon isolation information for a muon with undefined references to tracks";
1354
1355 reco::MuonIsolation isoR03, isoR05;
1356
1357
1358 reco::IsoDeposit depTrk = muIsoExtractorTrack_->deposit(iEvent, iSetup, *track);
1359 std::vector<reco::IsoDeposit> caloDeps = muIsoExtractorCalo_->deposits(iEvent, iSetup, *track);
1360 reco::IsoDeposit depJet = muIsoExtractorJet_->deposit(iEvent, iSetup, *track);
1361
1362 if (caloDeps.size() != 3) {
1363 LogTrace("MuonIdentification") << "Failed to fill vector of calorimeter isolation deposits!";
1364 return;
1365 }
1366
1367 reco::IsoDeposit depEcal = caloDeps.at(0);
1368 reco::IsoDeposit depHcal = caloDeps.at(1);
1369 reco::IsoDeposit depHo = caloDeps.at(2);
1370
1371
1372 if (writeIsoDeposits_) {
1373 trackDep = depTrk;
1374 ecalDep = depEcal;
1375 hcalDep = depHcal;
1376 hoDep = depHo;
1377 jetDep = depJet;
1378 }
1379
1380 isoR03.sumPt = depTrk.depositWithin(0.3);
1381 isoR03.emEt = depEcal.depositWithin(0.3);
1382 isoR03.hadEt = depHcal.depositWithin(0.3);
1383 isoR03.hoEt = depHo.depositWithin(0.3);
1384 isoR03.nTracks = depTrk.depositAndCountWithin(0.3).second;
1385 isoR03.nJets = depJet.depositAndCountWithin(0.3).second;
1386 isoR03.trackerVetoPt = depTrk.candEnergy();
1387 isoR03.emVetoEt = depEcal.candEnergy();
1388 isoR03.hadVetoEt = depHcal.candEnergy();
1389 isoR03.hoVetoEt = depHo.candEnergy();
1390
1391 isoR05.sumPt = depTrk.depositWithin(0.5);
1392 isoR05.emEt = depEcal.depositWithin(0.5);
1393 isoR05.hadEt = depHcal.depositWithin(0.5);
1394 isoR05.hoEt = depHo.depositWithin(0.5);
1395 isoR05.nTracks = depTrk.depositAndCountWithin(0.5).second;
1396 isoR05.nJets = depJet.depositAndCountWithin(0.5).second;
1397 isoR05.trackerVetoPt = depTrk.candEnergy();
1398 isoR05.emVetoEt = depEcal.candEnergy();
1399 isoR05.hadVetoEt = depHcal.candEnergy();
1400 isoR05.hoVetoEt = depHo.candEnergy();
1401
1402 aMuon.setIsolation(isoR03, isoR05);
1403 }
1404
1405 reco::Muon MuonIdProducer::makeMuon(const reco::Track& track) {
1406 const double energy = std::sqrt(track.p() * track.p() + 0.105658369 * 0.105658369);
1407 const math::XYZTLorentzVector p4(track.px(), track.py(), track.pz(), energy);
1408 return reco::Muon(track.charge(), p4, track.vertex());
1409 }
1410
1411 double MuonIdProducer::sectorPhi(const DetId& id) {
1412 double phi = 0;
1413 if (id.subdetId() == MuonSubdetId::DT) {
1414 DTChamberId muonId(id.rawId());
1415 if (muonId.sector() <= 12)
1416 phi = (muonId.sector() - 1) / 6. * M_PI;
1417 if (muonId.sector() == 13)
1418 phi = 3 / 6. * M_PI;
1419 if (muonId.sector() == 14)
1420 phi = 9 / 6. * M_PI;
1421 }
1422 if (id.subdetId() == MuonSubdetId::CSC) {
1423 CSCDetId muonId(id.rawId());
1424 phi = M_PI / 4 + (muonId.triggerSector() - 1) / 3. * M_PI;
1425 }
1426 if (phi > M_PI)
1427 phi -= 2 * M_PI;
1428 return phi;
1429 }
1430
1431 double MuonIdProducer::phiOfMuonInteractionRegion(const reco::Muon& muon) const {
1432 if (muon.isStandAloneMuon())
1433 return muon.standAloneMuon()->innerPosition().phi();
1434
1435 if (muon.matches().empty()) {
1436 if (muon.innerTrack().isAvailable() && muon.innerTrack()->extra().isAvailable())
1437 return muon.innerTrack()->outerPosition().phi();
1438 else
1439 return muon.phi();
1440 }
1441 return sectorPhi(muon.matches().at(0).id);
1442 }
1443
1444 void MuonIdProducer::fillGlbQuality(edm::Event& iEvent, const edm::EventSetup& iSetup, reco::Muon& aMuon) {
1445 if (aMuon.isGlobalMuon() && glbQualHandle_.isValid() && !glbQualHandle_.failedToGet()) {
1446 aMuon.setCombinedQuality((*glbQualHandle_)[aMuon.combinedMuon()]);
1447 }
1448
1449 LogDebug("MuonIdentification") << "tkChiVal " << aMuon.combinedQuality().trkRelChi2;
1450 }
1451
1452 void MuonIdProducer::fillTrackerKink(reco::Muon& aMuon) {
1453
1454 if (aMuon.innerTrack().isNull())
1455 return;
1456
1457 reco::MuonQuality quality = (aMuon.isQualityValid() ? aMuon.combinedQuality() : reco::MuonQuality());
1458
1459 const bool filled = trackerKinkFinder_->fillTrkKink(quality, *aMuon.innerTrack());
1460
1461 if (filled || aMuon.isQualityValid())
1462 aMuon.setCombinedQuality(quality);
1463 }
1464
1465 bool MuonIdProducer::checkLinks(const reco::MuonTrackLinks* links) const {
1466 const bool trackBAD = links->trackerTrack().isNull();
1467 const bool staBAD = links->standAloneTrack().isNull();
1468 const bool glbBAD = links->globalTrack().isNull();
1469 if (trackBAD || staBAD || glbBAD) {
1470 edm::LogWarning("muonIDbadLinks") << "Global muon links to constituent tracks are invalid: trkBad " << trackBAD
1471 << " standaloneBad " << staBAD << " globalBad " << glbBAD
1472 << ". There should be no such object. Muon is skipped.";
1473 return false;
1474 }
1475 return true;
1476 }
1477
1478 void MuonIdProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
1479 edm::ParameterSetDescription desc;
1480 desc.setAllowAnything();
1481
1482 desc.add<bool>("arbitrateTrackerMuons", false);
1483 desc.add<bool>("storeCrossedHcalRecHits", false);
1484 desc.add<bool>("fillShowerDigis", false);
1485 desc.ifValue(
1486 edm::ParameterDescription<bool>("selectHighPurity", false, true),
1487 true >> (edm::ParameterDescription<edm::InputTag>("pvInputTag", edm::InputTag("offlinePrimaryVertices"), true)) or
1488 false >> (edm::ParameterDescription<edm::InputTag>("pvInputTag", edm::InputTag(""), true)));
1489
1490 edm::ParameterSetDescription descTrkAsoPar;
1491 descTrkAsoPar.add<edm::InputTag>("GEMSegmentCollectionLabel", edm::InputTag("gemSegments"));
1492 descTrkAsoPar.add<edm::InputTag>("ME0SegmentCollectionLabel", edm::InputTag("me0Segments"));
1493 descTrkAsoPar.add<bool>("useGEM", false);
1494 descTrkAsoPar.add<bool>("useME0", false);
1495 descTrkAsoPar.setAllowAnything();
1496 desc.add<edm::ParameterSetDescription>("TrackAssociatorParameters", descTrkAsoPar);
1497
1498 edm::ParameterSetDescription descJet;
1499 descJet.setAllowAnything();
1500 descJet.add<edm::ParameterSetDescription>("TrackAssociatorParameters", descTrkAsoPar);
1501 desc.add<edm::ParameterSetDescription>("JetExtractorPSet", descJet);
1502
1503 edm::ParameterSetDescription descCalo;
1504 descCalo.setAllowAnything();
1505 descCalo.add<edm::ParameterSetDescription>("TrackAssociatorParameters", descTrkAsoPar);
1506 desc.add<edm::ParameterSetDescription>("CaloExtractorPSet", descCalo);
1507
1508 descriptions.addDefault(desc);
1509 }