From f732d4ff385477fc08ec4de29c695e1034e3dd39 Mon Sep 17 00:00:00 2001 From: mmusich Date: Tue, 21 Feb 2023 12:00:20 +0100 Subject: [PATCH] add a fillDescriptions method for TrackerTrackHitFilter --- .../plugins/TrackerTrackHitFilter.cc | 39 ++++++++- .../python/TrackerTrackHitFilter_cff.py | 83 ++++++++++--------- 2 files changed, 81 insertions(+), 41 deletions(-) diff --git a/RecoTracker/FinalTrackSelectors/plugins/TrackerTrackHitFilter.cc b/RecoTracker/FinalTrackSelectors/plugins/TrackerTrackHitFilter.cc index 0e8eeffe7e672..1eda7ba7de73b 100644 --- a/RecoTracker/FinalTrackSelectors/plugins/TrackerTrackHitFilter.cc +++ b/RecoTracker/FinalTrackSelectors/plugins/TrackerTrackHitFilter.cc @@ -5,6 +5,8 @@ #include "FWCore/MessageLogger/interface/MessageLogger.h" #include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" #include "FWCore/Utilities/interface/InputTag.h" #include "FWCore/Utilities/interface/ESGetToken.h" @@ -85,6 +87,8 @@ namespace reco { std::vector &hits); void produceFromTrack(const edm::EventSetup &iSetup, const Track *itt, std::vector &hits); + static void fillDescriptions(edm::ConfigurationDescriptions &descriptions); + private: class Rule { public: @@ -300,7 +304,7 @@ namespace reco { stripFrontInvalidHits_(iConfig.getParameter("stripFrontInvalidHits")), stripBackInvalidHits_(iConfig.getParameter("stripBackInvalidHits")), stripAllInvalidHits_(iConfig.getParameter("stripAllInvalidHits")), - isPhase2_(iConfig.getUntrackedParameter("isPhase2", false)), + isPhase2_(iConfig.getParameter("isPhase2")), rejectBadStoNHits_(iConfig.getParameter("rejectBadStoNHits")), CMNSubtractionMode_(iConfig.getParameter("CMNSubtractionMode")), detsToIgnore_(iConfig.getParameter >("detsToIgnore")), @@ -990,6 +994,39 @@ namespace reco { return tTopo->side(id); } + // ------------ method fills 'descriptions' with the allowed parameters for the module ------------ + void TrackerTrackHitFilter::fillDescriptions(edm::ConfigurationDescriptions &descriptions) { + edm::ParameterSetDescription desc; + desc.setComment(""); + desc.add("src", edm::InputTag("generalTracks")); + desc.add("minimumHits", 3)->setComment("number of hits for refit"); + desc.add("replaceWithInactiveHits", false) + ->setComment( + " instead of removing hits replace them with inactive hits, so you still consider the multiple " + "scattering"); + desc.add("stripFrontInvalidHits", false) + ->setComment("strip invalid & inactive hits from any end of the track"); + desc.add("stripBackInvalidHits", false) + ->setComment("strip invalid & inactive hits from any end of the track"); + desc.add("stripAllInvalidHits", false)->setComment("dangerous to turn on, you might forget about MS"); + desc.add("isPhase2", false); + desc.add("rejectBadStoNHits", false); + desc.add("CMNSubtractionMode", std::string("Median"))->setComment("TT6"); + desc.add >("detsToIgnore", {}); + desc.add("useTrajectories", false); + desc.add("rejectLowAngleHits", false); + desc.add("TrackAngleCut", 0.25)->setComment("rad"); + desc.add("usePixelQualityFlag", false); + desc.add("PxlTemplateProbXYCut", 0.000125); + desc.add("PxlTemplateProbXYChargeCut", -99.); + desc.add >("PxlTemplateqBinCut", {0, 3}); + desc.add("PxlCorrClusterChargeCut", -999.0); + desc.add("tagOverlaps", false); + desc.add >("commands", {})->setComment("layers to remove"); + desc.add >("StoNcommands", {})->setComment("S/N cut per layer"); + descriptions.addWithDefaultLabel(desc); + } + } // namespace modules } // namespace reco diff --git a/RecoTracker/FinalTrackSelectors/python/TrackerTrackHitFilter_cff.py b/RecoTracker/FinalTrackSelectors/python/TrackerTrackHitFilter_cff.py index cce7b97dee595..1fd75c86247b9 100644 --- a/RecoTracker/FinalTrackSelectors/python/TrackerTrackHitFilter_cff.py +++ b/RecoTracker/FinalTrackSelectors/python/TrackerTrackHitFilter_cff.py @@ -1,43 +1,46 @@ import FWCore.ParameterSet.Config as cms -TrackerTrackHitFilter = cms.EDProducer("TrackerTrackHitFilter", - src = cms.InputTag("generalTracks"), - minimumHits =cms.uint32(3), ##min number of hits for refit - ## # layers to remove - commands = cms.vstring( - "drop PXB", "drop PXE" ### same works for TIB, TID, TOB, TEC, - #"drop TIB 3", ## you can also drop specific layers/wheel/disks - #"keep PXB 3", ## you can also 'keep' some layer after - ##having dropped the whole structure - ), - - ###list of individual detids to turn off, in addition to the structures above - detsToIgnore = cms.vuint32( ), - - ### what to do with invalid hits - replaceWithInactiveHits =cms.bool(False), ## instead of removing hits replace - ## them with inactive hits, so you still - ## consider the multiple scattering - stripFrontInvalidHits =cms.bool(False), ## strip invalid & inactive hits from - stripBackInvalidHits =cms.bool(False), ## any end of the track - - stripAllInvalidHits = cms.bool(False), ##not sure if it's better 'true' or 'false' - ## might be dangerous to turn on - ## as you will forget about MS +from RecoTracker.FinalTrackSelectors.trackerTrackHitFilter_cfi import trackerTrackHitFilter as _trackerTrackHitFilter +TrackerTrackHitFilter = _trackerTrackHitFilter.clone( + src = "generalTracks", + minimumHits = 3, ##min number of hits for refit + ## # layers to remove + commands = ["drop PXB", "drop PXE"], ### same works for TIB, TID, TOB, TEC, + # "drop TIB 3", ## you can also drop specific layers/wheel/disks + # "keep PXB 3", ## you can also 'keep' some layer after + # having dropped the whole structure - ### hit quality cuts - rejectBadStoNHits = cms.bool(False), - CMNSubtractionMode = cms.string("Median"), ## "TT6" - StoNcommands = cms.vstring( - "TIB 1.0 ", "TOB 1.0 999.0" - ), - useTrajectories=cms.bool(False), - rejectLowAngleHits=cms.bool(False), - TrackAngleCut=cms.double(0.25), ## in radians - tagOverlaps=cms.bool(False), - usePixelQualityFlag=cms.bool(False), - PxlTemplateProbXYCut=cms.double(0.000125), #recommended by experts - PxlTemplateProbXYChargeCut=cms.double(-99.), #recommended by experts - PxlTemplateqBinCut =cms.vint32(0, 3), #recommended by experts - PxlCorrClusterChargeCut = cms.double(-999.0) - )####end of module + ###list of individual detids to turn off, in addition to the structures above + detsToIgnore = [], + + ### what to do with invalid hits + replaceWithInactiveHits = False, ## instead of removing hits replace + ## them with inactive hits, so you still + ## consider the multiple scattering + + stripFrontInvalidHits = False, ## strip invalid & inactive hits from + stripBackInvalidHits = False, ## any end of the track + + stripAllInvalidHits = False, ## not sure if it's better 'true' or 'false' + ## might be dangerous to turn on + ## as you will forget about MS + + ### hit quality cuts + isPhase2 = False, + rejectBadStoNHits = False, + CMNSubtractionMode = "Median", ## "TT6" + StoNcommands = ["TIB 1.0 ", "TOB 1.0 999.0"], + useTrajectories = False, + rejectLowAngleHits = False, + TrackAngleCut = 0.25, ## in radians + tagOverlaps= False, + usePixelQualityFlag = False, + PxlTemplateProbXYCut = 0.000125, # recommended by experts + PxlTemplateProbXYChargeCut = -99., # recommended by experts + PxlTemplateqBinCut = [0, 3], # recommended by experts + PxlCorrClusterChargeCut = -999.0 +) #### end of module + +from Configuration.Eras.Modifier_phase2_tracker_cff import phase2_tracker +phase2_tracker.toModify(TrackerTrackHitFilter, + isPhase2 = True)