ATLAS Offline Software
Loading...
Searching...
No Matches
BPHY28 Namespace Reference

Functions

 BPHY28Kernel (flags)
 BPHY28Cfg (flags)

Variables

str BPHYDerivationName = "BPHY28"
str streamName = "StreamDAOD_BPHY28"
str MuMuContainerName = "BPHY28MuMuCandidates"
str BsPhiMuMuContainerName = "BPHY28BsKKMuMuCandidates"

Function Documentation

◆ BPHY28Cfg()

BPHY28.BPHY28Cfg ( flags)

Definition at line 204 of file BPHY28.py.

204def BPHY28Cfg(flags):
205 doLRT = flags.Tracking.doLargeD0
206 isSimulation = flags.Input.isMC
207 acc = BPHY28Kernel(flags)
208 from DerivationFrameworkCore.SlimmingHelper import SlimmingHelper
209 from OutputStreamAthenaPool.OutputStreamConfig import OutputStreamCfg
210 from xAODMetaDataCnv.InfileMetaDataConfig import SetupMetaDataForStreamCfg
211 BPHY28SlimmingHelper = SlimmingHelper("BPHY28SlimmingHelper", NamesAndTypes = flags.Input.TypedCollections, flags = flags)
212 from DerivationFrameworkBPhys.commonBPHYMethodsCfg import getDefaultAllVariables
213 AllVariables = getDefaultAllVariables()
214 StaticContent = []
215
216 # Needed for trigger objects
217 BPHY28SlimmingHelper.IncludeMuonTriggerContent = True
218 BPHY28SlimmingHelper.IncludeBPhysTriggerContent = True
219
220
221 AllVariables += ["PrimaryVertices"]
222 StaticContent += ["xAOD::VertexContainer#BPHY28RefittedPrimaryVertices"]
223 StaticContent += ["xAOD::VertexAuxContainer#BPHY28RefittedPrimaryVerticesAux."]
224
225
226
227 AllVariables += ["InDetTrackParticles", "InDetLargeD0TrackParticles"] if doLRT else ["InDetTrackParticles"]
228
229
232 AllVariables += ["CombinedMuonTrackParticles"]
233 AllVariables += ["ExtrapolatedMuonTrackParticles"]
234
235
236 AllVariables += ["Muons", "MuonsLRT"] if doLRT else ["Muons"]
237
238
239
240 StaticContent += ["xAOD::VertexContainer#%s" % MuMuContainerName]
241
242 StaticContent += ["xAOD::VertexAuxContainer#%sAux.-vxTrackAtVertex" % MuMuContainerName]
243
244 StaticContent += ["xAOD::VertexContainer#%s" % BsPhiMuMuContainerName]
245 StaticContent += ["xAOD::VertexAuxContainer#%sAux.-vxTrackAtVertex" % BsPhiMuMuContainerName]
246
247
248 # Tagging information (in addition to that already requested by usual algorithms)
249 AllVariables += [ "MuonSpectrometerTrackParticles" ]
250 SmartVar = ["Photons", "Electrons", "LRTElectrons"] if doLRT else ["Photons", "Electrons"] #[ tagJetCollections ]
251
252
253 # Truth information for MC only
254 if isSimulation:
255 AllVariables += ["TruthEvents","TruthParticles","TruthVertices","MuonTruthParticles", "egammaTruthParticles" ]
256
257
258 AllVariables = list(set(AllVariables)) # remove duplicates
259
260 BPHY28SlimmingHelper.AllVariables = AllVariables
261 BPHY28SlimmingHelper.StaticContent = StaticContent
262 BPHY28SlimmingHelper.SmartCollections = SmartVar
263 BPHY28ItemList = BPHY28SlimmingHelper.GetItemList()
264 acc.merge(OutputStreamCfg(flags, "DAOD_BPHY28", ItemList=BPHY28ItemList, AcceptAlgs=["BPHY28Kernel"]))
265 acc.merge(SetupMetaDataForStreamCfg(flags, "DAOD_BPHY28", AcceptAlgs=["BPHY28Kernel"], createMetadata=[MetadataCategory.CutFlowMetaData]))
266 acc.printConfig(withDetails=True, summariseProps=True, onlyComponents = [], printDefaults=True)
267 return acc
STL class.

◆ BPHY28Kernel()

BPHY28.BPHY28Kernel ( flags)

Definition at line 19 of file BPHY28.py.

19def BPHY28Kernel(flags):
20
21 # Lists for better code organization
22 augList = [] # List of active augmentation tools
23 skimList = [] # List of active skimming algorithms
24 thinList = [] # List of active thinning algorithms
25 thinTrkVtxList = [] # List of reconstructed candidates to use for the thinning of tracks from vertices
26
27
28 from DerivationFrameworkBPhys.commonBPHYMethodsCfg import (
29 BPHY_V0ToolCfg, BPHY_InDetDetailedTrackSelectorToolCfg,
30 BPHY_VertexPointEstimatorCfg, BPHY_TrkVKalVrtFitterCfg,
31 AugOriginalCountsCfg)
32 from JpsiUpsilonTools.JpsiUpsilonToolsConfig import PrimaryVertexRefittingToolCfg
33 acc = ComponentAccumulator()
34 isSimulation = flags.Input.isMC
35 V0Tools = acc.popToolsAndMerge(BPHY_V0ToolCfg(flags, BPHYDerivationName))
36 vkalvrt = acc.popToolsAndMerge(BPHY_TrkVKalVrtFitterCfg(flags, BPHYDerivationName)) # VKalVrt vertex fitter
37 acc.addPublicTool(vkalvrt)
38 acc.addPublicTool(V0Tools)
39 trackselect = acc.popToolsAndMerge(BPHY_InDetDetailedTrackSelectorToolCfg(flags, BPHYDerivationName))
40 acc.addPublicTool(trackselect)
41 vpest = acc.popToolsAndMerge(BPHY_VertexPointEstimatorCfg(flags, BPHYDerivationName))
42 acc.addPublicTool(vpest)
43 PVrefit = acc.popToolsAndMerge(PrimaryVertexRefittingToolCfg(flags))
44 acc.addPublicTool(PVrefit)
45
46 # LRT
47 doLRT = flags.Tracking.doLargeD0
48 if not doLRT : print("BPHY28: LRT tracks disabled")
49 mainMuonInput = "StdWithLRTMuons" if doLRT else "Muons"
50 mainIDInput = "InDetWithLRTTrackParticles" if doLRT else "InDetTrackParticles"
51 if doLRT:
52 from DerivationFrameworkLLP.LLPToolsConfig import LRTMuonMergerAlg
53 from AthenaConfiguration.Enums import LHCPeriod
54 acc.merge(LRTMuonMergerAlg( flags,
55 PromptMuonLocation = "Muons",
56 LRTMuonLocation = "MuonsLRT",
57 OutputMuonLocation = mainMuonInput,
58 CreateViewCollection = True,
59 UseRun3WP = flags.GeoModel.Run == LHCPeriod.Run3))
60 from DerivationFrameworkInDet.InDetToolsConfig import InDetLRTMergeCfg
61 acc.merge(InDetLRTMergeCfg(flags))
62 toRelink = ["InDetTrackParticles", "InDetLargeD0TrackParticles"] if doLRT else []
63 MuonReLink = [ "Muons", "MuonsLRT" ] if doLRT else []
64
65 BPHY28_AugOriginalCounts = acc.popToolsAndMerge(
66 AugOriginalCountsCfg(flags, name = "BPHY28_AugOriginalCounts"))
67 augList += [ BPHY28_AugOriginalCounts ]
68
69 BPHY28MuMuFinder = CompFactory.Analysis.JpsiFinder(
70 name = "BPHY28MuMuFinder",
71 muAndMu = True,
72 muAndTrack = False,
73 TrackAndTrack = False,
74 assumeDiMuons = True,
75 invMassUpper = 100000.0,
76 invMassLower = 0.0,
77 Chi2Cut = 100.,
78 oppChargesOnly = True,
79 combOnly = True,
80 atLeastOneComb = False,
81 useCombinedMeasurement = False, # Only takes effect if combOnly=True
82 muonCollectionKey = mainMuonInput,
83 TrackParticleCollection = mainIDInput,
84 useV0Fitter = False, # if False a TrkVertexFitterTool will be used
85 TrkVertexFitterTool = vkalvrt,
86 V0VertexFitterTool = None,
87 TrackSelectorTool = trackselect,
88 VertexPointEstimator = vpest,
89 useMCPCuts = False )
90
91 BPHY28MuMuSelectAndWrite = CompFactory.DerivationFramework.Reco_Vertex(name = "BPHY28MuMuSelectAndWrite",
92 VertexSearchTool = BPHY28MuMuFinder,
93 OutputVtxContainerName = MuMuContainerName,
94 PVContainerName = "PrimaryVertices",
95 V0Tools = V0Tools,
96 PVRefitter = PVrefit,
97 RefPVContainerName = "SHOULDNOTBEUSED",
98 RelinkTracks = toRelink,
99 RelinkMuons = MuonReLink,
100 DoVertexType = 7)
101 augList += [ BPHY28MuMuSelectAndWrite ]
102 thinTrkVtxList += [ MuMuContainerName ]
103
104
105 BPHY28_Select_Jpsi2mumu = CompFactory.DerivationFramework.Select_onia2mumu(
106 name = "BPHY28_Select_Jpsi2mumu",
107 HypothesisName = "Jpsi",
108 InputVtxContainerName = MuMuContainerName,
109 V0Tools = V0Tools,
110 VtxMassHypo = 3096.916,
111 MassMin = 2000.0,
112 MassMax = 3600.0,
113 Chi2Max = 100, Do3d = False,
114 DoVertexType = 7)
115 augList += [ BPHY28_Select_Jpsi2mumu ]
116
117
118 BPHY28BsKKMuMu = CompFactory.Analysis.JpsiPlus2Tracks(name = "BPHY28BsKKMuMu",
119 kaonkaonHypothesis = True,
120 pionpionHypothesis = False,
121 kaonpionHypothesis = False,
122 trkThresholdPt = 500.0,
123 trkMaxEta = 3.0,
124 BMassUpper = 5900.0,
125 BMassLower = 4900.0,
126 DiTrackMassUpper = 1220,
127 DiTrackMassLower = 820,
128 Chi2Cut = 100.0, # this is chi2/ndf cut
129 TrkQuadrupletMassUpper = 6000.0,
130 TrkQuadrupletMassLower = 4800.0,
131 JpsiContainerKey = MuMuContainerName,
132 TrackParticleCollection = mainIDInput,
133 MuonsUsedInJpsi = mainMuonInput,
134 TrkVertexFitterTool = vkalvrt,
135 TrackSelectorTool = trackselect,
136 UseMassConstraint = False)
137
138
139 BPHY28BsKKSelectAndWrite = CompFactory.DerivationFramework.Reco_Vertex(name = "BPHY28BsKKSelectAndWrite",
140 VertexSearchTool = BPHY28BsKKMuMu,
141 OutputVtxContainerName = BsPhiMuMuContainerName,
142 PVContainerName = "PrimaryVertices",
143 V0Tools = V0Tools,
144 PVRefitter = PVrefit,
145 RefPVContainerName = "BPHY28RefittedPrimaryVertices",
146 RefitPV = True, Do3d = False,
147 RelinkTracks = toRelink,
148 MaxPVrefit = 10000, DoVertexType = 7)
149 augList += [ BPHY28BsKKSelectAndWrite ]
150 thinTrkVtxList += [ BsPhiMuMuContainerName ]
151
152
153 BPHY28_Select_Bs2KKMuMu = CompFactory.DerivationFramework.Select_onia2mumu(
154 name = "BPHY28_Select_Bs2KKMuMu",
155 HypothesisName = "Bs",
156 InputVtxContainerName = BsPhiMuMuContainerName,
157 V0Tools = V0Tools,
158 TrkMasses = [105.658, 105.658, 493.677, 493.677],
159 VtxMassHypo = 5366.3,
160 MassMin = 4900.0,
161 MassMax = 5900.0, Do3d = False,
162 Chi2Max = 100)
163 augList += [ BPHY28_Select_Bs2KKMuMu ]
164
165
166 if not isSimulation: #Only Skim Data
167 from DerivationFrameworkTools.DerivationFrameworkToolsConfig import (
168 xAODStringSkimmingToolCfg)
169 BPHY28_SelectBsKKMuMuEvent = acc.getPrimaryAndMerge(xAODStringSkimmingToolCfg(
170 flags, name = "BPHY28_SelectBsKKMuMuEvent",
171 expression = "count(BPHY28BsKKMuMuCandidates.passed_Bs) > 0"))
172 skimList += [ BPHY28_SelectBsKKMuMuEvent ]
173
174
175
177 BPHY28_Thin_VtxTracks = CompFactory.DerivationFramework.Thin_vtxTrk( name = "BPHY28_Thin_VtxTracks",
178 StreamName = streamName,
179 TrackParticleContainerName = "InDetTrackParticles",
180 VertexContainerNames = thinTrkVtxList,
181 IgnoreFlags = True )
182 thinList += [ BPHY28_Thin_VtxTracks ]
183
184 # LRT ID tracks
185 if doLRT:
186 BPHY28_Thin_VtxTracks_LRT = CompFactory.DerivationFramework.Thin_vtxTrk( name = "BPHY28_Thin_VtxTracks_LRT",
187 StreamName = streamName,
188 TrackParticleContainerName = "InDetLargeD0TrackParticles",
189 VertexContainerNames = thinTrkVtxList,
190 IgnoreFlags = True )
191 thinList += [ BPHY28_Thin_VtxTracks_LRT ]
192
193
194 for t in augList + skimList + thinList : acc.addPublicTool(t)
195 acc.addEventAlgo(CompFactory.DerivationFramework.DerivationKernel("BPHY28Kernel",
196 AugmentationTools = augList,
197 #Only skim if not MC
198 SkimmingTools = skimList,
199 ThinningTools = thinList))
200
201 return acc
202
203
void print(char *figname, TCanvas *c1)

Variable Documentation

◆ BPHYDerivationName

str BPHY28.BPHYDerivationName = "BPHY28"

Definition at line 12 of file BPHY28.py.

◆ BsPhiMuMuContainerName

str BPHY28.BsPhiMuMuContainerName = "BPHY28BsKKMuMuCandidates"

Definition at line 16 of file BPHY28.py.

◆ MuMuContainerName

str BPHY28.MuMuContainerName = "BPHY28MuMuCandidates"

Definition at line 15 of file BPHY28.py.

◆ streamName

str BPHY28.streamName = "StreamDAOD_BPHY28"

Definition at line 13 of file BPHY28.py.