diff --git a/HLTrigger/Configuration/python/customizeHLTforCMSSW.py b/HLTrigger/Configuration/python/customizeHLTforCMSSW.py index 7a71d4c06a3d0..ea543e716867c 100644 --- a/HLTrigger/Configuration/python/customizeHLTforCMSSW.py +++ b/HLTrigger/Configuration/python/customizeHLTforCMSSW.py @@ -188,6 +188,27 @@ def customiseForEcalTestPR22254thresholdC(process): +def customizeHLTForL3OIPR24267(process): + for seedproducer in producers_by_type(process, "TSGForOI"): + if "hltIterL3OISeedsFromL2Muons" == seedproducer.label(): + process.hltIterL3OISeedsFromL2Muons = cms.EDProducer("TSGForOIFromL2") + if "hltIterL3OISeedsFromL2MuonsOpenMu" == seedproducer.label(): + process.hltIterL3OISeedsFromL2MuonsOpenMu = cms.EDProducer("TSGForOIFromL2") + if "hltIterL3OISeedsFromL2MuonsNoVtx" == seedproducer.label(): + process.hltIterL3OISeedsFromL2MuonsNoVtx = cms.EDProducer("TSGForOIFromL2") + process.hltIterL3OISeedsFromL2MuonsNoVtx.src = cms.InputTag( 'hltL2Muons' ) + + + for trackproducer in producers_by_type(process, "CkfTrackCandidateMaker"): + if "hltIterL3OITrackCandidates" in trackproducer.label(): + trackproducer.reverseTrajectories =cms.bool(True) + + + return process + + + + # CMSSW version specific customizations def customizeHLTforCMSSW(process, menuType="GRun"): @@ -195,4 +216,6 @@ def customizeHLTforCMSSW(process, menuType="GRun"): # add call to action function in proper order: newest last! # process = customiseFor12718(process) + customizeHLTForL3OIPR24267(process) + return process diff --git a/RecoMuon/TrackerSeedGenerator/plugins/TSGForOIFromL2.cc b/RecoMuon/TrackerSeedGenerator/plugins/TSGForOIFromL2.cc new file mode 100644 index 0000000000000..0b288cd29e865 --- /dev/null +++ b/RecoMuon/TrackerSeedGenerator/plugins/TSGForOIFromL2.cc @@ -0,0 +1,461 @@ +/** + \class TSGForOIFromL2 + \brief Create L3MuonTrajectorySeeds from L2 Muons updated at vertex in an outside-in manner + \author Benjamin Radburn-Smith, Santiago Folgueras, Bibhuprasad Mahakud, Jan Frederik Schulte (Purdue University, West Lafayette, USA) + */ + +#include "RecoMuon/TrackerSeedGenerator/plugins/TSGForOIFromL2.h" +#include "DataFormats/TrackerCommon/interface/TrackerTopology.h" +#include "DataFormats/Math/interface/deltaR.h" + +#include + + +TSGForOIFromL2::TSGForOIFromL2(const edm::ParameterSet & iConfig) : + src_(consumes(iConfig.getParameter("src"))), + maxSeeds_(iConfig.getParameter("maxSeeds")), + maxHitlessSeeds_(iConfig.getParameter("maxHitlessSeeds")), + maxHitSeeds_(iConfig.getParameter("maxHitSeeds")), + numOfLayersToTry_(iConfig.getParameter("layersToTry")), + numOfHitsToTry_(iConfig.getParameter("hitsToTry")), + numL2ValidHitsCutAllEta_(iConfig.getParameter("numL2ValidHitsCutAllEta")), + numL2ValidHitsCutAllEndcap_(iConfig.getParameter("numL2ValidHitsCutAllEndcap")), + fixedErrorRescalingForHits_(iConfig.getParameter("fixedErrorRescaleFactorForHits")), + fixedErrorRescalingForHitless_(iConfig.getParameter("fixedErrorRescaleFactorForHitless")), + adjustErrorsDynamicallyForHits_(iConfig.getParameter("adjustErrorsDynamicallyForHits")), + adjustErrorsDynamicallyForHitless_(iConfig.getParameter("adjustErrorsDynamicallyForHitless")), + estimatorName_(iConfig.getParameter("estimator")), + minEtaForTEC_(iConfig.getParameter("minEtaForTEC")), + maxEtaForTOB_(iConfig.getParameter("maxEtaForTOB")), + useHitLessSeeds_(iConfig.getParameter("UseHitLessSeeds")), + updator_(new KFUpdator()), measurementTrackerTag_(consumes(iConfig.getParameter("MeasurementTrackerEvent"))), + pT1_(iConfig.getParameter("pT1")), + pT2_(iConfig.getParameter("pT2")), + pT3_(iConfig.getParameter("pT3")), + eta1_(iConfig.getParameter("eta1")), + eta2_(iConfig.getParameter("eta2")), + eta3_(iConfig.getParameter("eta3")), + eta4_(iConfig.getParameter("eta4")), + eta5_(iConfig.getParameter("eta5")), + eta6_(iConfig.getParameter("eta6")), + eta7_(iConfig.getParameter("eta7")), + SF1_(iConfig.getParameter("SF1")), + SF2_(iConfig.getParameter("SF2")), + SF3_(iConfig.getParameter("SF3")), + SF4_(iConfig.getParameter("SF4")), + SF5_(iConfig.getParameter("SF5")), + SF6_(iConfig.getParameter("SF6")), + tsosDiff1_(iConfig.getParameter("tsosDiff1")), + tsosDiff2_(iConfig.getParameter("tsosDiff2")), + propagatorName_(iConfig.getParameter("propagatorName")), + theCategory_(std::string("Muon|RecoMuon|TSGForOIFromL2")) +{ + produces >(); +} + +TSGForOIFromL2::~TSGForOIFromL2() { + +} + +// +// Produce seeds +// +void TSGForOIFromL2::produce(edm::StreamID sid, edm::Event& iEvent, const edm::EventSetup& iSetup) const { + + // Initialize variables + unsigned int numSeedsMade = 0; + unsigned int layerCount = 0; + unsigned int hitlessSeedsMadeIP = 0; + unsigned int hitlessSeedsMadeMuS = 0; + unsigned int hitSeedsMade = 0; + + // Surface used to make a TSOS at the PCA to the beamline + Plane::PlanePointer dummyPlane = Plane::build(Plane::PositionType(), Plane::RotationType()); + + // Read ESHandles + edm::Handle measurementTrackerH; + edm::ESHandle estimatorH; + edm::ESHandle magfieldH; + edm::ESHandle propagatorAlongH; + edm::ESHandle propagatorOppositeH; + edm::ESHandle tmpTkGeometryH; + edm::ESHandle geometryH; + + iSetup.get().get(magfieldH); + iSetup.get().get(propagatorName_, propagatorOppositeH); + iSetup.get().get(propagatorName_, propagatorAlongH); + iSetup.get().get(geometryH); + iSetup.get().get(tmpTkGeometryH); + iSetup.get().get(estimatorName_,estimatorH); + iEvent.getByToken(measurementTrackerTag_, measurementTrackerH); + + // Read L2 track collection + edm::Handle l2TrackCol; + iEvent.getByToken(src_, l2TrackCol); + + // The product + std::unique_ptr > result(new std::vector()); + + // Get vector of Detector layers + std::vector const& tob = measurementTrackerH->geometricSearchTracker()->tobLayers(); + std::vector const& tecPositive = tmpTkGeometryH->isThere(GeomDetEnumerators::P2OTEC) ? + measurementTrackerH->geometricSearchTracker()->posTidLayers() : + measurementTrackerH->geometricSearchTracker()->posTecLayers(); + std::vector const& tecNegative = tmpTkGeometryH->isThere(GeomDetEnumerators::P2OTEC) ? + measurementTrackerH->geometricSearchTracker()->negTidLayers() : + measurementTrackerH->geometricSearchTracker()->negTecLayers(); + + // Get suitable propagators + std::unique_ptr propagatorAlong = SetPropagationDirection(*propagatorAlongH,alongMomentum); + std::unique_ptr propagatorOpposite = SetPropagationDirection(*propagatorOppositeH,oppositeToMomentum); + + // Stepping Helix Propagator for propogation from muon system to tracker + edm::ESHandle SHPOpposite; + iSetup.get().get("hltESPSteppingHelixPropagatorOpposite", SHPOpposite); + + // Loop over the L2's and make seeds for all of them + LogTrace(theCategory_) << "TSGForOIFromL2::produce: Number of L2's: " << l2TrackCol->size(); + for (unsigned int l2TrackColIndex(0); l2TrackColIndex!=l2TrackCol->size(); ++l2TrackColIndex) { + + const reco::TrackRef l2(l2TrackCol, l2TrackColIndex); + + // Container of Seeds + std::vector out; + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2::produce: L2 muon pT, eta, phi --> " << l2->pt() << " , " << l2->eta() << " , " << l2->phi() << std::endl; + + FreeTrajectoryState fts = trajectoryStateTransform::initialFreeState(*l2, magfieldH.product()); + + dummyPlane->move(fts.position() - dummyPlane->position()); + TrajectoryStateOnSurface tsosAtIP = TrajectoryStateOnSurface(fts, *dummyPlane); + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2::produce: Created TSOSatIP: " << tsosAtIP << std::endl; + + // Get the TSOS on the innermost layer of the L2 + TrajectoryStateOnSurface tsosAtMuonSystem = trajectoryStateTransform::innerStateOnSurface(*l2, *geometryH, magfieldH.product()); + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2::produce: Created TSOSatMuonSystem: " << tsosAtMuonSystem <eta(); + double absL2muonEta=std::abs(L2muonEta); + bool useBoth = false; + if (outerTkStateInside.isValid() && outerTkStateOutside.isValid()) { + //following commented out variables dist1 (5 par compatibility of tsos at outertracker surface) + //dist2 (angle between two tsos) could further be explored in combination of L2 valid hits for seeding. So kept for + //future developers + //auto dist1 = match_Chi2(outerTkStateInside,outerTkStateOutside);//for future developers + //auto dist2 = deltaR(outerTkStateInside.globalMomentum(),outerTkStateOutside.globalMomentum());//for future developers + //if ((dist1 > tsosDiff1_ || dist2 > tsosDiff2_) && l2->numberOfValidHits() < 20) useBoth = true;//for future developers + if (l2->numberOfValidHits() < numL2ValidHitsCutAllEta_) useBoth = true; + if (l2->numberOfValidHits() < numL2ValidHitsCutAllEndcap_ && absL2muonEta > eta7_) useBoth = true; + if (absL2muonEta > eta1_ && absL2muonEta < eta1_) useBoth = true; + } + + numSeedsMade = 0; + hitlessSeedsMadeIP = 0; + hitlessSeedsMadeMuS = 0; + hitSeedsMade = 0; + + // calculate scale factors + double errorSFHits = (adjustErrorsDynamicallyForHits_? calculateSFFromL2(l2) : fixedErrorRescalingForHits_); + double errorSFHitless = (adjustErrorsDynamicallyForHitless_? calculateSFFromL2(l2) : fixedErrorRescalingForHitless_); + + // BARREL + if (absL2muonEta < maxEtaForTOB_) { + layerCount = 0; + for (auto it=tob.rbegin(); it!=tob.rend(); ++it) { + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2::produce: looping in TOB layer " << layerCount << std::endl; + if ( useHitLessSeeds_ && hitlessSeedsMadeIP < maxHitlessSeeds_ && numSeedsMade < maxSeeds_ ) + makeSeedsWithoutHits(**it, tsosAtIP, *(propagatorAlong.get()), estimatorH, errorSFHitless, hitlessSeedsMadeIP, numSeedsMade, out); + + // Do not create hitbased seeds in barrel region + if (absL2muonEta > 1.0 && hitSeedsMade < maxHitSeeds_ && numSeedsMade < maxSeeds_ ) + makeSeedsFromHits(**it, tsosAtIP, *(propagatorAlong.get()), estimatorH, measurementTrackerH, errorSFHits, hitSeedsMade, numSeedsMade, layerCount, out); + + if (useBoth) { + if ( useHitLessSeeds_ && hitlessSeedsMadeMuS < maxHitlessSeeds_ && numSeedsMade < maxSeeds_ ) + makeSeedsWithoutHits(**it, outerTkStateOutside, *(propagatorOpposite.get()), estimatorH, errorSFHitless, hitlessSeedsMadeMuS, numSeedsMade, out); + } + } + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2:::produce: NumSeedsMade = " << numSeedsMade << " , layerCount = " << layerCount << std::endl; + } + + // Reset number of seeds if in overlap region + if (absL2muonEta > minEtaForTEC_ && absL2muonEta < maxEtaForTOB_) { + numSeedsMade = 0; + hitlessSeedsMadeIP = 0; + hitlessSeedsMadeMuS = 0; + hitSeedsMade = 0; + } + + // ENDCAP+ + if (L2muonEta > minEtaForTEC_) { + layerCount = 0; + for (auto it=tecPositive.rbegin(); it!=tecPositive.rend(); ++it) { + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2::produce: looping in TEC+ layer " << layerCount << std::endl; + if ( useHitLessSeeds_ && hitlessSeedsMadeIP < maxHitlessSeeds_ && numSeedsMade < maxSeeds_ ) + makeSeedsWithoutHits(**it, tsosAtIP, *(propagatorAlong.get()), estimatorH, errorSFHitless, hitlessSeedsMadeIP, numSeedsMade, out); + + if (absL2muonEta > 1.0 && hitSeedsMade < maxHitSeeds_ && numSeedsMade < maxSeeds_ ) + makeSeedsFromHits(**it, tsosAtIP, *(propagatorAlong.get()), estimatorH, measurementTrackerH, errorSFHits, hitSeedsMade, numSeedsMade, layerCount, out); + + if (useBoth) { + if ( useHitLessSeeds_ && hitlessSeedsMadeMuS < maxHitlessSeeds_ && numSeedsMade < maxSeeds_ ) + makeSeedsWithoutHits(**it, outerTkStateOutside, *(propagatorOpposite.get()), estimatorH, errorSFHitless, hitlessSeedsMadeMuS, numSeedsMade, out); + } + } + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2:::produce: NumSeedsMade = " << numSeedsMade << " , layerCount = " << layerCount << std::endl; + } + + // ENDCAP- + if (L2muonEta < -minEtaForTEC_) { + layerCount = 0; + for (auto it=tecNegative.rbegin(); it!=tecNegative.rend(); ++it) { + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2::produce: looping in TEC- layer " << layerCount << std::endl; + if ( useHitLessSeeds_ && hitlessSeedsMadeIP < maxHitlessSeeds_ && numSeedsMade < maxSeeds_ ) + makeSeedsWithoutHits(**it, tsosAtIP, *(propagatorAlong.get()), estimatorH, errorSFHitless, hitlessSeedsMadeIP, numSeedsMade, out); + + if (absL2muonEta > 1.0 && hitSeedsMade < maxHitSeeds_ && numSeedsMade < maxSeeds_ ) + makeSeedsFromHits(**it, tsosAtIP, *(propagatorAlong.get()), estimatorH, measurementTrackerH, errorSFHits, hitSeedsMade, numSeedsMade, layerCount, out); + + if (useBoth) { + if ( useHitLessSeeds_ && hitlessSeedsMadeMuS < maxHitlessSeeds_ && numSeedsMade < maxSeeds_ ) + makeSeedsWithoutHits(**it, outerTkStateOutside, *(propagatorOpposite.get()), estimatorH, errorSFHitless, hitlessSeedsMadeMuS, numSeedsMade, out); + } + } + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2:::produce: NumSeedsMade = " << numSeedsMade << " , layerCount = " << layerCount << std::endl; + } + + for (std::vector::iterator it=out.begin(); it!=out.end(); ++it) { + result->push_back(*it); + } + + } // L2Collection + + edm::LogInfo(theCategory_) << "TSGForOIFromL2::produce: number of seeds made: " << result->size(); + + iEvent.put(std::move(result)); +} + + +// +// Create seeds without hits on a given layer (TOB or TEC) +// +void TSGForOIFromL2::makeSeedsWithoutHits( + const GeometricSearchDet& layer, + const TrajectoryStateOnSurface& tsos, + const Propagator& propagatorAlong, + edm::ESHandle& estimator, + double errorSF, + unsigned int& hitlessSeedsMade, + unsigned int& numSeedsMade, + std::vector& out) const { + + // create hitless seeds + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2::makeSeedsWithoutHits: Start hitless" << std::endl; + std::vector< GeometricSearchDet::DetWithState > dets; + layer.compatibleDetsV(tsos, propagatorAlong, *estimator, dets); + if (!dets.empty()) { + auto const& detOnLayer = dets.front().first; + auto const& tsosOnLayer = dets.front().second; + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2::makeSeedsWithoutHits: tsosOnLayer " << tsosOnLayer << std::endl; + if (!tsosOnLayer.isValid()) { + edm::LogInfo(theCategory_) << "ERROR!: Hitless TSOS is not valid!"; + } + else { + dets.front().second.rescaleError(errorSF); + PTrajectoryStateOnDet const& ptsod = trajectoryStateTransform::persistentState(tsosOnLayer,detOnLayer->geographicalId().rawId()); + TrajectorySeed::recHitContainer rHC; + out.push_back(TrajectorySeed(ptsod,rHC,oppositeToMomentum)); + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2::makeSeedsWithoutHits: TSOS (Hitless) done " << std::endl; + hitlessSeedsMade++; + numSeedsMade++; + } + } + +} + + +// +// Find hits on a given layer (TOB or TEC) and create seeds from updated TSOS with hit +// +void TSGForOIFromL2::makeSeedsFromHits( + const GeometricSearchDet& layer, + const TrajectoryStateOnSurface& tsos, + const Propagator& propagatorAlong, + edm::ESHandle& estimator, + edm::Handle& measurementTracker, + double errorSF, + unsigned int& hitSeedsMade, + unsigned int& numSeedsMade, + unsigned int& layerCount, + std::vector& out) const { + + if (layerCount > numOfLayersToTry_) return; + + // Error Rescaling + TrajectoryStateOnSurface onLayer(tsos); + onLayer.rescaleError(errorSF); + + std::vector< GeometricSearchDet::DetWithState > dets; + layer.compatibleDetsV(onLayer, propagatorAlong, *estimator, dets); + + // Find Measurements on each DetWithState + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2::makeSeedsFromHits: Find measurements on each detWithState " << dets.size() << std::endl; + std::vector meas; + for (std::vector::iterator it=dets.begin(); it!=dets.end(); ++it) { + MeasurementDetWithData det = measurementTracker->idToDet(it->first->geographicalId()); + if (det.isNull()) continue; + if (!it->second.isValid()) continue; // Skip if TSOS is not valid + + std::vector mymeas = det.fastMeasurements(it->second, onLayer, propagatorAlong, *estimator); // Second TSOS is not used + for (std::vector::const_iterator it2 = mymeas.begin(), ed2 = mymeas.end(); it2 != ed2; ++it2) { + if (it2->recHit()->isValid()) meas.push_back(*it2); // Only save those which are valid + } + } + + // Update TSOS using TMs after sorting, then create Trajectory Seed and put into vector + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2::makeSeedsFromHits: Update TSOS using TMs after sorting, then create Trajectory Seed, number of TM = " << meas.size() << std::endl; + std::sort(meas.begin(), meas.end(), TrajMeasLessEstim()); + + unsigned int found = 0; + for (std::vector::const_iterator it=meas.begin(); it!=meas.end(); ++it) { + TrajectoryStateOnSurface updatedTSOS = updator_->update(it->forwardPredictedState(), *it->recHit()); + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2::makeSeedsFromHits: TSOS for TM " << found << std::endl; + if (not updatedTSOS.isValid()) continue; + + edm::OwnVector seedHits; + seedHits.push_back(*it->recHit()->hit()); + PTrajectoryStateOnDet const& pstate = trajectoryStateTransform::persistentState(updatedTSOS, it->recHit()->geographicalId().rawId()); + TrajectorySeed seed(pstate, std::move(seedHits), oppositeToMomentum); + LogTrace("TSGForOIFromL2") << "TSGForOIFromL2::makeSeedsFromHits: Number of seedHits: " << seedHits.size() << std::endl; + out.push_back(seed); + found++; + numSeedsMade++; + hitSeedsMade++; + if (found == numOfHitsToTry_) break; + if (hitSeedsMade > maxHitSeeds_) return; + } + + if (found) layerCount++; + +} + + +// +// Calculate the dynamic error SF by analysing the L2 +// +double TSGForOIFromL2::calculateSFFromL2(const reco::TrackRef track) const { + + double theSF = 1.0; + // L2 direction vs pT blowup - as was previously done: + // Split into 4 pT ranges: eta()); + if (track->pt() <= pT1_) theSF = SF1_; + else if (track->pt() > pT1_ && track->pt() <= pT2_) { + if (abseta <= eta3_) theSF = SF3_; + else if (abseta > eta3_ && abseta <= eta6_) theSF = SF2_; + else if (abseta > eta6_) theSF = SF3_; + } + else if (track->pt() > pT2_ && track->pt() <= pT3_) { + if (abseta <= eta1_) theSF = SF6_; + else if (abseta > eta1_ && abseta <= eta2_) theSF = SF4_; + else if (abseta > eta2_ && abseta <= eta3_) theSF = SF6_; + else if (abseta > eta3_ && abseta <= eta4_) theSF = SF1_; + else if (abseta > eta4_ && abseta <= eta5_) theSF = SF1_; + else if (abseta > eta5_ ) theSF = SF5_; + } + else if (track->pt() > pT3_) { + if (abseta <= eta3_) theSF = SF5_; + else if (abseta > eta3_ && abseta <= eta4_) theSF = SF4_; + else if (abseta > eta4_ && abseta <= eta5_) theSF = SF4_; + else if (abseta > eta5_ ) theSF = SF5_; + } + + LogTrace(theCategory_) << "TSGForOIFromL2::calculateSFFromL2: SF has been calculated as: " << theSF; + + return theSF; + +} + +// +// calculate Chi^2 of two trajectory states +// +double TSGForOIFromL2::match_Chi2(const TrajectoryStateOnSurface& tsos1, + const TrajectoryStateOnSurface& tsos2) const { + + if ( !tsos1.isValid() || !tsos2.isValid() ) return -1.; + + AlgebraicVector5 v(tsos1.localParameters().vector() - tsos2.localParameters().vector()); + AlgebraicSymMatrix55 m(tsos1.localError().matrix() + tsos2.localError().matrix()); + + bool ierr = !m.Invert(); + + if ( ierr ) { + edm::LogInfo("TSGForOIFromL2") << "Error inverting covariance matrix"; + return -1; + } + + double est = ROOT::Math::Similarity(v,m); + + return est; + +} + + +// +// +// +void TSGForOIFromL2::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + + edm::ParameterSetDescription desc; + desc.add("src",edm::InputTag("hltL2Muons","UpdatedAtVtx")); + desc.add("layersToTry",2); + desc.add("fixedErrorRescaleFactorForHitless",2.0); + desc.add("hitsToTry",1); + desc.add("adjustErrorsDynamicallyForHits",false); + desc.add("adjustErrorsDynamicallyForHitless",true); + desc.add("MeasurementTrackerEvent",edm::InputTag("hltSiStripClusters")); + desc.add("UseHitLessSeeds",true); + desc.add("estimator","hltESPChi2MeasurementEstimator100"); + desc.add("maxEtaForTOB",1.8); + desc.add("minEtaForTEC",0.7); + desc.addUntracked("debug",false); + desc.add("fixedErrorRescaleFactorForHits",1.0); + desc.add("maxSeeds",20); + desc.add("maxHitlessSeeds",5); + desc.add("maxHitSeeds",1); + desc.add("numL2ValidHitsCutAllEta",20); + desc.add("numL2ValidHitsCutAllEndcap",30); + desc.add("pT1",13.0); + desc.add("pT2",30.0); + desc.add("pT3",70.0); + desc.add("eta1",0.2); + desc.add("eta2",0.3); + desc.add("eta3",1.0); + desc.add("eta4",1.2); + desc.add("eta5",1.6); + desc.add("eta6",1.4); + desc.add("eta7",2.1); + desc.add("SF1",3.0); + desc.add("SF2",4.0); + desc.add("SF3",5.0); + desc.add("SF4",7.0); + desc.add("SF5",10.0); + desc.add("SF6",2.0); + desc.add("tsosDiff1",0.2); + desc.add("tsosDiff2",0.02); + desc.add("propagatorName","PropagatorWithMaterialParabolicMf"); + descriptions.add("TSGForOIFromL2",desc); + +} + +DEFINE_FWK_MODULE(TSGForOIFromL2); diff --git a/RecoMuon/TrackerSeedGenerator/plugins/TSGForOIFromL2.h b/RecoMuon/TrackerSeedGenerator/plugins/TSGForOIFromL2.h new file mode 100644 index 0000000000000..a797fb2a1d8cd --- /dev/null +++ b/RecoMuon/TrackerSeedGenerator/plugins/TSGForOIFromL2.h @@ -0,0 +1,139 @@ +#ifndef RecoMuon_TrackerSeedGenerator_TSGForOIFromL2_H +#define RecoMuon_TrackerSeedGenerator_TSGForOIFromL2_H + +/** + \class TSGForOIFromL2 + \brief Create L3MuonTrajectorySeeds from L2 Muons updated at vertex in an outside-in manner + \author Benjamin Radburn-Smith, Santiago Folgueras, Bibhuprasad Mahakud, Jan Frederik Schulte (Purdue University, West Lafayette) + */ + +#include "DataFormats/TrackReco/interface/Track.h" +#include "FWCore/Framework/interface/global/EDProducer.h" +#include "FWCore/Framework/interface/Event.h" +#include "FWCore/Framework/interface/EventSetup.h" +#include "FWCore/Framework/interface/ESHandle.h" +#include "FWCore/Framework/interface/MakerMacros.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "Geometry/CommonDetUnit/interface/GlobalTrackingGeometry.h" +#include "Geometry/Records/interface/GlobalTrackingGeometryRecord.h" +#include "MagneticField/Engine/interface/MagneticField.h" +#include "MagneticField/Records/interface/IdealMagneticFieldRecord.h" +#include "RecoTracker/MeasurementDet/interface/MeasurementTrackerEvent.h" +#include "TrackingTools/GeomPropagators/interface/Propagator.h" +#include "TrackingTools/KalmanUpdators/interface/Chi2MeasurementEstimator.h" +#include "TrackingTools/KalmanUpdators/interface/KFUpdator.h" +#include "TrackingTools/MeasurementDet/interface/MeasurementDet.h" +#include "TrackingTools/PatternTools/interface/TrajMeasLessEstim.h" +#include "TrackingTools/Records/interface/TrackingComponentsRecord.h" +#include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h" +#include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h" +#include "TrackingTools/GeomPropagators/interface/StateOnTrackerBound.h" +#include "DataFormats/TrackerCommon/interface/TrackerTopology.h" + +class TSGForOIFromL2 : public edm::global::EDProducer<> { + + public: + + explicit TSGForOIFromL2(const edm::ParameterSet & iConfig); + ~TSGForOIFromL2() override; + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + void produce(edm::StreamID sid, edm::Event& iEvent, const edm::EventSetup& iSetup) const override; + + private: + + /// Labels for input collections + const edm::EDGetTokenT src_; + + /// Maximum number of seeds for each L2 + const unsigned int maxSeeds_; + + /// Maximum number of hitless seeds for each L2 + const unsigned int maxHitlessSeeds_; + + /// Maximum number of hitbased seeds for each L2 + const unsigned int maxHitSeeds_; + + /// How many layers to try + const unsigned int numOfLayersToTry_; + + /// How many hits to try per layer + const unsigned int numOfHitsToTry_; + + ///L2 valid hit cuts to decide seed creation by both states + const unsigned int numL2ValidHitsCutAllEta_; + const unsigned int numL2ValidHitsCutAllEndcap_; + + /// Rescale L2 parameter uncertainties (fixed error vs pT, eta) + const double fixedErrorRescalingForHits_; + const double fixedErrorRescalingForHitless_; + + /// Whether or not to use an automatically calculated scale-factor value + const bool adjustErrorsDynamicallyForHits_; + const bool adjustErrorsDynamicallyForHitless_; + + /// Estimator used to find dets and TrajectoryMeasurements + const std::string estimatorName_; + + /// Minimum eta value to activate searching in the TEC + const double minEtaForTEC_; + + /// Maximum eta value to activate searching in the TOB + const double maxEtaForTOB_; + + /// Switch ON (True) : use additional hits for seeds depending on the L2 properties (ignores numOfMaxSeeds_) + /// Switch OFF (False): the numOfMaxSeeds_ defines if we will use hitless (numOfMaxSeeds_==1) or hitless+hits (numOfMaxSeeds_>1) + const bool useHitLessSeeds_; + + /// KFUpdator defined in constructor + const std::unique_ptr updator_; + + const edm::EDGetTokenT measurementTrackerTag_; + + /// pT, eta ranges and scale factor values + const double pT1_,pT2_,pT3_; + const double eta1_,eta2_,eta3_,eta4_,eta5_,eta6_,eta7_; + const double SF1_,SF2_,SF3_,SF4_,SF5_,SF6_; + + /// Distance of L2 TSOSs before and after updated with vertex + const double tsosDiff1_; + const double tsosDiff2_; + + /// Counters and flags for the implementation + const std::string propagatorName_; + const std::string theCategory_; + + /// Create seeds without hits on a given layer (TOB or TEC) + void makeSeedsWithoutHits( + const GeometricSearchDet& layer, + const TrajectoryStateOnSurface& tsos, + const Propagator& propagatorAlong, + edm::ESHandle& estimator, + double errorSF, + unsigned int& hitlessSeedsMade, + unsigned int& numSeedsMade, + std::vector& out) const; + + /// Find hits on a given layer (TOB or TEC) and create seeds from updated TSOS with hit + void makeSeedsFromHits( + const GeometricSearchDet& layer, + const TrajectoryStateOnSurface& tsos, + const Propagator& propagatorAlong, + edm::ESHandle& estimator, + edm::Handle& measurementTracker, + double errorSF, + unsigned int& hitSeedsMade, + unsigned int& numSeedsMade, + unsigned int& layerCount, + std::vector& out) const; + + /// Calculate the dynamic error SF by analysing the L2 + double calculateSFFromL2(const reco::TrackRef track) const; + + /// Find compatability between two TSOSs + double match_Chi2(const TrajectoryStateOnSurface& tsos1, + const TrajectoryStateOnSurface& tsos2) const; + +}; + +#endif