File indexing completed on 2024-04-06 12:13:29
0001
0002 import FWCore.ParameterSet.Config as cms
0003
0004 process = cms.Process('GEN')
0005
0006
0007 process.load('Configuration.StandardSequences.Services_cff')
0008 process.load('SimGeneral.HepPDTESSource.pythiapdt_cfi')
0009 process.load('FWCore.MessageService.MessageLogger_cfi')
0010 process.load('Configuration.EventContent.EventContent_cff')
0011 process.load('SimGeneral.MixingModule.mixNoPU_cfi')
0012 process.load('Configuration.StandardSequences.GeometryRecoDB_cff')
0013 process.load('Configuration.StandardSequences.MagneticField_cff')
0014 process.load('Configuration.StandardSequences.Generator_cff')
0015 process.load('IOMC.EventVertexGenerators.VtxSmearedRealistic50ns13TeVCollision_cfi')
0016 process.load('GeneratorInterface.Core.genFilterSummary_cff')
0017 process.load('Configuration.StandardSequences.EndOfProcess_cff')
0018 process.load('Configuration.StandardSequences.FrontierConditions_GlobalTag_cff')
0019
0020 process.maxEvents = cms.untracked.PSet(
0021 input = cms.untracked.int32(200)
0022 )
0023
0024
0025 process.source = cms.Source("EmptySource")
0026
0027 process.options = cms.untracked.PSet()
0028
0029
0030
0031 process.RAWSIMoutput = cms.OutputModule("PoolOutputModule",
0032 SelectEvents = cms.untracked.PSet(
0033 SelectEvents = cms.vstring('generation_step')
0034 ),
0035 dataset = cms.untracked.PSet(
0036 dataTier = cms.untracked.string('GEN'),
0037 filterName = cms.untracked.string('')
0038 ),
0039 eventAutoFlushCompressedSize = cms.untracked.int32(5242880),
0040 fileName = cms.untracked.string('file:gen_lambdab.root'),
0041 outputCommands = process.RAWSIMEventContent.outputCommands,
0042 splitLevel = cms.untracked.int32(0)
0043 )
0044
0045
0046
0047
0048 process.genstepfilter.triggerConditions=cms.vstring("generation_step")
0049 from Configuration.AlCa.GlobalTag import GlobalTag
0050 process.GlobalTag = GlobalTag(process.GlobalTag, 'auto:run2_mc', '')
0051
0052 process.generator = cms.EDFilter("Pythia8GeneratorFilter",
0053 ExternalDecays = cms.PSet(
0054 EvtGen1 = cms.untracked.PSet(
0055 convertPythiaCodes = cms.untracked.bool(False),
0056 decay_table = cms.string('GeneratorInterface/EvtGenInterface/data/DECAY_NOLONGLIFE.DEC'),
0057 list_forced_decays = cms.vstring('MyLambda_b0','Myanti-Lambda_b0'),
0058 particle_property_file = cms.FileInPath('GeneratorInterface/EvtGenInterface/data/evt.pdl'),
0059 operates_on_particles = cms.vint32(5122),
0060 user_decay_embedded = cms.vstring(
0061 """
0062 Alias MyLambda_b0 Lambda_b0
0063 Alias Myanti-Lambda_b0 anti-Lambda_b0
0064 ChargeConj Myanti-Lambda_b0 MyLambda_b0
0065 #
0066 Decay MyLambda_b0
0067 1.000 p+ mu- anti-nu_mu PHOTOS Lb2plnuLCSR 1 1 1 1;
0068 #1.000 p+ mu- anti-nu_mu Lb2plnuLCSR 1 1 1 1;
0069 Enddecay
0070 CDecay Myanti-Lambda_b0
0071 #
0072 End
0073 """
0074 )
0075 ),
0076 parameterSets = cms.vstring('EvtGen1')
0077 ),
0078 PythiaParameters = cms.PSet(
0079 parameterSets = cms.vstring('pythia8CommonSettings',
0080 'pythia8CUEP8M1Settings',
0081 'processParameters'),
0082 processParameters = cms.vstring(
0083 'SoftQCD:nonDiffractive = on',
0084 'PTFilter:filter = on',
0085 'PTFilter:quarkToFilter = 5',
0086 'PTFilter:scaleToFilter = 1.0'),
0087 pythia8CUEP8M1Settings = cms.vstring('Tune:pp 14',
0088 'Tune:ee 7',
0089 'MultipartonInteractions:pT0Ref=2.4024',
0090 'MultipartonInteractions:ecmPow=0.25208',
0091 'MultipartonInteractions:expPow=1.6'),
0092 pythia8CommonSettings = cms.vstring('Tune:preferLHAPDF = 2',
0093 'Main:timesAllowErrors = 10000',
0094 'Check:epTolErr = 0.01',
0095 'Beams:setProductionScalesFromLHEF = off',
0096 'SLHA:minMassSM = 1000.',
0097 'ParticleDecays:limitTau0 = on',
0098 'ParticleDecays:tau0Max = 10',
0099 'ParticleDecays:allowPhotonRadiation = on')
0100 ),
0101 comEnergy = cms.double(13000.0),
0102 maxEventsToPrint = cms.untracked.int32(0),
0103 pythiaHepMCVerbosity = cms.untracked.bool(False),
0104 pythiaPylistVerbosity = cms.untracked.int32(0)
0105 )
0106
0107
0108 from GeneratorInterface.EvtGenInterface.EvtGenSetting_cff import *
0109 process.generator.PythiaParameters.processParameters.extend(EvtGenExtraParticles)
0110
0111 process.lbfilter = cms.EDFilter("PythiaDauVFilter",
0112 DaughterIDs = cms.untracked.vint32(2212, 13),
0113 MaxEta = cms.untracked.vdouble(2.5, 2.5),
0114 MinEta = cms.untracked.vdouble(-2.5, -2.5),
0115 MinPt = cms.untracked.vdouble(0.4, 1.5),
0116 MotherID = cms.untracked.int32(0),
0117 NumberDaughters = cms.untracked.int32(2),
0118 ParticleID = cms.untracked.int32(5122)
0119 )
0120
0121 process.ProductionFilterSequence = cms.Sequence(process.generator+process.lbfilter)
0122
0123
0124 process.generation_step = cms.Path(process.pgen)
0125 process.genfiltersummary_step = cms.EndPath(process.genFilterSummary)
0126 process.endjob_step = cms.EndPath(process.endOfProcess)
0127 process.RAWSIMoutput_step = cms.EndPath(process.RAWSIMoutput)
0128
0129
0130 process.schedule = cms.Schedule(process.generation_step,process.genfiltersummary_step,process.endjob_step,process.RAWSIMoutput_step)
0131
0132 for path in process.paths:
0133 getattr(process,path)._seq = process.ProductionFilterSequence * getattr(process,path)._seq
0134
0135 process.MessageLogger.cerr.FwkReport.reportEvery = 10
0136
0137