Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 12:46:08

0001 from __future__ import print_function
0002 import FWCore.ParameterSet.Config as cms
0003 import sys
0004 from enum import Enum
0005 from PhysicsTools.PatAlgos.patInputFiles_cff import filesRelValTTbarPileUpGENSIMRECO
0006 
0007 class RefitType(Enum):
0008      STANDARD = 1
0009      COMMON   = 2
0010  
0011 isDA = True
0012 isMC = True
0013 allFromGT = True
0014 applyBows = True
0015 applyExtraConditions = True
0016 theRefitter = RefitType.COMMON
0017 
0018 process = cms.Process("Demo") 
0019 
0020 ###################################################################
0021 # Event source and run selection
0022 ###################################################################
0023 process.source = cms.Source("PoolSource",
0024                             fileNames = filesRelValTTbarPileUpGENSIMRECO,
0025                             duplicateCheckMode = cms.untracked.string('checkAllFilesOpened')
0026                             )
0027 
0028 runboundary = 1
0029 process.source.firstRun = cms.untracked.uint32(int(runboundary))
0030 process.maxEvents = cms.untracked.PSet( input = cms.untracked.int32(10) )
0031 
0032 ###################################################################
0033 # JSON Filtering
0034 ###################################################################
0035 if isMC:
0036      print(">>>>>>>>>> testPVValidation_cfg.py: msg%-i: This is Simulation!")
0037      runboundary = 1
0038 else:
0039      print(">>>>>>>>>> testPVValidation_cfg.py: msg%-i: This is DATA!")
0040      import FWCore.PythonUtilities.LumiList as LumiList
0041      process.source.lumisToProcess = LumiList.LumiList(filename ='None').getVLuminosityBlockRange()
0042 
0043 ###################################################################
0044 # Messages
0045 ###################################################################
0046 process.load('FWCore.MessageService.MessageLogger_cfi')   
0047 process.MessageLogger.cerr.enable = False
0048 process.MessageLogger.PrimaryVertexValidation=dict()  
0049 process.MessageLogger.SplitVertexResolution=dict()
0050 process.MessageLogger.FilterOutLowPt=dict()  
0051 process.MessageLogger.cout = cms.untracked.PSet(
0052     enable = cms.untracked.bool(True),
0053     threshold = cms.untracked.string("INFO"),
0054     default   = cms.untracked.PSet(limit = cms.untracked.int32(0)),                       
0055     FwkReport = cms.untracked.PSet(limit = cms.untracked.int32(-1),
0056                                    reportEvery = cms.untracked.int32(1000)
0057                                    ),                                                      
0058     PrimaryVertexValidation = cms.untracked.PSet( limit = cms.untracked.int32(-1)),
0059     SplitVertexResolution   = cms.untracked.PSet( limit = cms.untracked.int32(-1)),
0060     FilterOutLowPt          = cms.untracked.PSet( limit = cms.untracked.int32(-1)),
0061     enableStatistics = cms.untracked.bool(True)
0062     )
0063 
0064 ####################################################################
0065 # Produce the Transient Track Record in the event
0066 ####################################################################
0067 process.load("TrackingTools.TransientTrack.TransientTrackBuilder_cfi")
0068 
0069 ####################################################################
0070 # Get the Magnetic Field
0071 ####################################################################
0072 process.load('Configuration.StandardSequences.MagneticField_cff')
0073 
0074 ###################################################################
0075 # Standard loads
0076 ###################################################################
0077 process.load("Configuration.Geometry.GeometryRecoDB_cff")
0078 
0079 ####################################################################
0080 # Get the BeamSpot
0081 ####################################################################
0082 process.load("RecoVertex.BeamSpotProducer.BeamSpot_cff")
0083 
0084 ####################################################################
0085 # Get the GlogalTag
0086 ####################################################################
0087 process.load("Configuration.StandardSequences.FrontierConditions_GlobalTag_cff")
0088 from Configuration.AlCa.GlobalTag import GlobalTag
0089 process.GlobalTag = GlobalTag(process.GlobalTag, 'auto:phase1_2017_realistic', '')
0090 
0091 if allFromGT:
0092      print(">>>>>>>>>> testPVValidation_cfg.py: msg%-i: All is taken from GT")
0093 else:
0094      ####################################################################
0095      # Get Alignment constants
0096      ####################################################################
0097      from CondCore.DBCommon.CondDBSetup_cfi import *
0098      process.trackerAlignment = cms.ESSource("PoolDBESSource",CondDBSetup,
0099                                              connect = cms.string('frontier://FrontierProd/CMS_CONDITIONS'),
0100                                              timetype = cms.string("runnumber"),
0101                                              toGet = cms.VPSet(cms.PSet(record = cms.string('TrackerAlignmentRcd'),
0102                                                                         tag = cms.string('TrackerAlignment_Upgrade2017_design_v4')
0103                                                                         )
0104                                                                )
0105                                              )
0106      process.es_prefer_trackerAlignment = cms.ESPrefer("PoolDBESSource", "trackerAlignment")
0107 
0108      ####################################################################
0109      # Get APE
0110      ####################################################################
0111      process.setAPE = cms.ESSource("PoolDBESSource",CondDBSetup,
0112                                    connect = cms.string('frontier://FrontierProd/CMS_CONDITIONS'),
0113                                    timetype = cms.string("runnumber"),
0114                                    toGet = cms.VPSet(cms.PSet(record = cms.string('TrackerAlignmentErrorExtendedRcd'),
0115                                                               tag = cms.string('TrackerAlignmentErrorsExtended_Upgrade2017_design_v0')
0116                                                               )
0117                                                      )
0118                                    )
0119      process.es_prefer_setAPE = cms.ESPrefer("PoolDBESSource", "setAPE")
0120 
0121      ####################################################################
0122      # Kinks and Bows (optional)
0123      ####################################################################
0124      if applyBows:
0125           print(">>>>>>>>>> testPVValidation_cfg.py: msg%-i: Applying TrackerSurfaceDeformations!")
0126           process.trackerBows = cms.ESSource("PoolDBESSource",CondDBSetup,
0127                                              connect = cms.string('frontier://FrontierProd/CMS_CONDITIONS'),
0128                                              toGet = cms.VPSet(cms.PSet(record = cms.string('TrackerSurfaceDeformationRcd'),
0129                                                                         tag = cms.string('TrackerSurfaceDeformations_zero')
0130                                                                         )
0131                                                                )
0132                                              )
0133           process.es_prefer_Bows = cms.ESPrefer("PoolDBESSource", "trackerBows")
0134      else:
0135           print(">>>>>>>>>> testPVValidation_cfg.py: msg%-i: MultiPVValidation: Not applying TrackerSurfaceDeformations!")
0136 
0137      ####################################################################
0138      # Extra corrections not included in the GT
0139      ####################################################################
0140      if applyExtraConditions:
0141 
0142           import CalibTracker.Configuration.Common.PoolDBESSource_cfi
0143           ##### END OF EXTRA CONDITIONS
0144  
0145      else:
0146           print(">>>>>>>>>> testPVValidation_cfg.py: msg%-i: Not applying extra calibration constants!")
0147      
0148 ####################################################################
0149 # Load and Configure event selection
0150 ####################################################################
0151 process.primaryVertexFilter = cms.EDFilter("VertexSelector",
0152                                              src = cms.InputTag("offlinePrimaryVertices"),
0153                                              cut = cms.string("!isFake && ndof > 4 && abs(z) <= 24 && position.Rho <= 2"),
0154                                              filter = cms.bool(True)
0155                                              )
0156 
0157 process.noscraping = cms.EDFilter("FilterOutScraping",
0158                                   applyfilter = cms.untracked.bool(True),
0159                                   src =  cms.untracked.InputTag("generalTracks"),
0160                                   debugOn = cms.untracked.bool(False),
0161                                   numtrack = cms.untracked.uint32(10),
0162                                   thresh = cms.untracked.double(0.25)
0163                                   )
0164 
0165 process.noslowpt = cms.EDFilter("FilterOutLowPt",
0166                                 applyfilter = cms.untracked.bool(True),
0167                                 src =  cms.untracked.InputTag("generalTracks"),
0168                                 debugOn = cms.untracked.bool(False),
0169                                 numtrack = cms.untracked.uint32(0),
0170                                 thresh = cms.untracked.int32(1),
0171                                 ptmin  = cms.untracked.double(3.),
0172                                 runControl = cms.untracked.bool(True),
0173                                 runControlNumber = cms.untracked.vuint32(int(runboundary))
0174                                 )
0175 
0176 if isMC:
0177      process.goodvertexSkim = cms.Sequence(process.noscraping)
0178 else:
0179      process.goodvertexSkim = cms.Sequence(process.primaryVertexFilter + process.noscraping + process.noslowpt)
0180 
0181 
0182 if(theRefitter == RefitType.COMMON):
0183 
0184      print(">>>>>>>>>> testPVValidation_cfg.py: msg%-i: using the common track selection and refit sequence!")          
0185      ####################################################################
0186      # Load and Configure Common Track Selection and refitting sequence
0187      ####################################################################
0188      import Alignment.CommonAlignment.tools.trackselectionRefitting as trackselRefit
0189      process.seqTrackselRefit = trackselRefit.getSequence(process, 'generalTracks',
0190                                                           isPVValidation=True, 
0191                                                           TTRHBuilder='WithAngleAndTemplate',
0192                                                           usePixelQualityFlag=True,
0193                                                           openMassWindow=False,
0194                                                           cosmicsDecoMode=True,
0195                                                           cosmicsZeroTesla=False,
0196                                                           momentumConstraint=None,
0197                                                           cosmicTrackSplitting=False,
0198                                                           use_d0cut=False,
0199                                                           )
0200      
0201 elif (theRefitter == RefitType.STANDARD):
0202 
0203      print(">>>>>>>>>> testPVValidation_cfg.py: msg%-i: using the standard single refit sequence!")          
0204      ####################################################################
0205      # Load and Configure Measurement Tracker Event
0206      # (needed in case NavigationSchool is set != '')
0207      ####################################################################
0208      # process.load("RecoTracker.MeasurementDet.MeasurementTrackerEventProducer_cfi") 
0209      # process.MeasurementTrackerEvent.pixelClusterProducer = 'generalTracks'
0210      # process.MeasurementTrackerEvent.stripClusterProducer = 'generalTracks'
0211      # process.MeasurementTrackerEvent.inactivePixelDetectorLabels = cms.VInputTag()
0212      # process.MeasurementTrackerEvent.inactiveStripDetectorLabels = cms.VInputTag()
0213 
0214      ####################################################################
0215      # Load and Configure TrackRefitter
0216      ####################################################################
0217      process.load("RecoTracker.TrackProducer.TrackRefitters_cff")
0218      import RecoTracker.TrackProducer.TrackRefitters_cff
0219      process.FinalTrackRefitter = RecoTracker.TrackProducer.TrackRefitter_cfi.TrackRefitter.clone()
0220      process.FinalTrackRefitter.src = "generalTracks"
0221      process.FinalTrackRefitter.TrajectoryInEvent = True
0222      process.FinalTrackRefitter.NavigationSchool = ''
0223      process.FinalTrackRefitter.TTRHBuilder = "WithAngleAndTemplate"
0224 
0225      ####################################################################
0226      # Sequence
0227      ####################################################################
0228      process.seqTrackselRefit = cms.Sequence(process.offlineBeamSpot*
0229                                              # in case NavigatioSchool is set !='' 
0230                                              #process.MeasurementTrackerEvent*
0231                                              process.FinalTrackRefitter)     
0232 
0233 ####################################################################
0234 # Output file
0235 ####################################################################
0236 process.TFileService = cms.Service("TFileService",
0237                                    fileName=cms.string("PVValidation_test_0.root")
0238                                   )
0239 
0240 ####################################################################
0241 # Imports of parameters
0242 ####################################################################
0243 from RecoVertex.PrimaryVertexProducer.OfflinePrimaryVertices_cfi import offlinePrimaryVertices
0244 ## modify the parameters which differ
0245 FilteringParams = offlinePrimaryVertices.TkFilterParameters.clone(
0246      maxNormalizedChi2 = 5.0,  # chi2ndof < 5
0247      maxD0Significance = 5.0,  # fake cut (requiring 1 PXB hit)
0248      maxEta = 5.0,             # as per recommendation in PR #18330
0249 )
0250 
0251 ## MM 04.05.2017 (use settings as in: https://github.com/cms-sw/cmssw/pull/18330)
0252 from RecoVertex.PrimaryVertexProducer.TkClusParameters_cff import DA_vectParameters
0253 DAClusterizationParams = DA_vectParameters.clone()
0254 
0255 GapClusterizationParams = cms.PSet(algorithm   = cms.string('gap'),
0256                                    TkGapClusParameters = cms.PSet(zSeparation = cms.double(0.2))  # 0.2 cm max separation betw. clusters
0257                                    )
0258 
0259 ####################################################################
0260 # Deterministic annealing clustering or Gap clustering
0261 ####################################################################
0262 def switchClusterizerParameters(da):
0263      if da:
0264           print(">>>>>>>>>> testPVValidation_cfg.py: msg%-i: Running DA Algorithm!")
0265           return DAClusterizationParams
0266      else:
0267           print(">>>>>>>>>> testPVValidation_cfg.py: msg%-i: Running GAP Algorithm!")
0268           return GapClusterizationParams
0269 
0270 ####################################################################
0271 # Configure the PVValidation Analyzer module
0272 ####################################################################
0273 process.PVValidation = cms.EDAnalyzer("PrimaryVertexValidation",
0274                                       TrackCollectionTag = cms.InputTag("FinalTrackRefitter"),
0275                                       VertexCollectionTag = cms.InputTag("offlinePrimaryVertices"),
0276                                       Debug = cms.bool(False),
0277                                       storeNtuple = cms.bool(False),
0278                                       useTracksFromRecoVtx = cms.bool(False),
0279                                       isLightNtuple = cms.bool(True),
0280                                       askFirstLayerHit = cms.bool(False),
0281                                       forceBeamSpot = cms.untracked.bool(False),
0282                                       probePt = cms.untracked.double(3.),
0283                                       minPt   = cms.untracked.double(1.),
0284                                       maxPt   = cms.untracked.double(30.),
0285                                       runControl = cms.untracked.bool(True),
0286                                       runControlNumber = cms.untracked.vuint32(int(runboundary)),
0287                                       TkFilterParameters = FilteringParams,
0288                                       TkClusParameters = switchClusterizerParameters(isDA)
0289                                       )
0290 
0291 ####################################################################
0292 # Path
0293 ####################################################################
0294 process.p = cms.Path(process.goodvertexSkim*
0295                      process.seqTrackselRefit*
0296                      process.PVValidation)
0297 
0298 ## PV refit part
0299 process.load("TrackingTools.TransientTrack.TransientTrackBuilder_cfi")
0300 process.offlinePrimaryVerticesFromRefittedTrks  = offlinePrimaryVertices.clone()
0301 process.offlinePrimaryVerticesFromRefittedTrks.TrackLabel                                       = cms.InputTag("FinalTrackRefitter")
0302 process.offlinePrimaryVerticesFromRefittedTrks.vertexCollections.maxDistanceToBeam              = 1
0303 process.offlinePrimaryVerticesFromRefittedTrks.TkFilterParameters.maxNormalizedChi2             = 20
0304 process.offlinePrimaryVerticesFromRefittedTrks.TkFilterParameters.minSiliconLayersWithHits      = 5
0305 process.offlinePrimaryVerticesFromRefittedTrks.TkFilterParameters.maxD0Significance             = 5.0
0306 # as it was prior to https://github.com/cms-sw/cmssw/commit/c8462ae4313b6be3bbce36e45373aa6e87253c59
0307 process.offlinePrimaryVerticesFromRefittedTrks.TkFilterParameters.maxD0Error                    = 1.0
0308 process.offlinePrimaryVerticesFromRefittedTrks.TkFilterParameters.maxDzError                    = 1.0
0309 process.offlinePrimaryVerticesFromRefittedTrks.TkFilterParameters.minPixelLayersWithHits        = 2
0310 
0311 ###################################################################
0312 # The trigger filter module
0313 ###################################################################
0314 from HLTrigger.HLTfilters.triggerResultsFilter_cfi import *
0315 process.HLTFilter = triggerResultsFilter.clone(
0316      triggerConditions = cms.vstring("HLT_ZeroBias_*"),
0317      #triggerConditions = cms.vstring("HLT_HT*"),
0318      hltResults = cms.InputTag( "TriggerResults", "", "HLT" ),
0319      l1tResults = cms.InputTag( "" ),
0320      throw = cms.bool(False)
0321 )
0322 
0323 ###################################################################
0324 # The analysis module
0325 ###################################################################
0326 process.myanalysis = cms.EDAnalyzer("GeneralPurposeTrackAnalyzer",
0327                                     TkTag  = cms.InputTag('FinalTrackRefitter'),
0328                                     isCosmics = cms.bool(False)
0329                                     )
0330 
0331 ###################################################################
0332 # The PV resolution module
0333 ###################################################################
0334 process.PrimaryVertexResolution = cms.EDAnalyzer('SplitVertexResolution',
0335                                                  storeNtuple         = cms.bool(True),
0336                                                  vtxCollection       = cms.InputTag("offlinePrimaryVerticesFromRefittedTrks"),
0337                                                  trackCollection     = cms.InputTag("FinalTrackRefitter"),
0338                                                  minVertexNdf        = cms.untracked.double(10.),
0339                                                  minVertexMeanWeight = cms.untracked.double(0.5),
0340                                                  runControl = cms.untracked.bool(True),
0341                                                  runControlNumber = cms.untracked.vuint32(int(runboundary))
0342                                                  )
0343 
0344 process.p2 = cms.Path(process.HLTFilter                               +
0345                       process.seqTrackselRefit                        +
0346                       process.offlinePrimaryVerticesFromRefittedTrks  +
0347                       process.PrimaryVertexResolution                 +
0348                       process.myanalysis
0349                       )