GenerateTruthJets.py 5.91 KB
Newer Older
1
2
3
4
5
6
from AthenaCommon import Logging
jetlog = Logging.logging.getLogger("GenerateTruthJets")

# Attach jet algorithms
def PrepareTruthJetInputs(algseq):
    print "will prepare thruth jets"
7
8
9
10
#    if "JetTruthCopyAlg" in algseq:
    if hasattr( algseq, "JetTruthCopyAlg" ): 
        "JetTruthCopyAlg already present, no need to add again"
        return
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
    from AthenaCommon.AppMgr import ToolSvc
    from MCTruthClassifier import MCTruthClassifierConf
    from ParticleJetTools import ParticleJetToolsConf
    from JetRec import JetRecConf
    #from MCTruthClassifier.MCTruthClassifierConfig import firstSimCreatedBarcode

    ToolSvc += MCTruthClassifierConf.MCTruthClassifier("JetMCTruthClassifier",
                                                       barcodeG4Shift=1000000) # or 200k? Doesn't matter anyway?
                                                       #barcodeG4Shift=firstSimCreatedBarcode())

    ToolSvc += ParticleJetToolsConf.CopyTruthJetParticles("truthpartcopy",
                                                          OutputName="JetInputTruthParticles",
                                                          MCTruthClassifier=ToolSvc.JetMCTruthClassifier,
                                                          BarCodeFromMetadata=0)

    ToolSvc += ParticleJetToolsConf.CopyTruthJetParticles("truthpartcopywz",
                                                          OutputName="JetInputTruthParticlesNoWZ",
                                                          MCTruthClassifier=ToolSvc.JetMCTruthClassifier,
                                                          IncludePromptLeptons=False,
                                                          #IncludePromptPhotons=False,
                                                          BarCodeFromMetadata=0)

33
34
35
36
37
    ToolSvc += ParticleJetToolsConf.CopyTruthPartons("truthpartonscopy",
                                                     OutputName="TruthLabelPartons",
                                                     PtMin=5000)

    algseq += JetRecConf.JetAlgorithm("JetTruthCopyAlg", Tools=[ToolSvc.truthpartcopy,ToolSvc.truthpartcopywz,ToolSvc.truthpartonscopy])
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66

def ScheduleAntiKtTruthJets(jetradius,algseq,mods=""):
    jetcollname = 'AntiKt{0}{1}TruthJets'.format(int(jetradius*10),mods)
    if jetcollname in algseq: return
    from AthenaCommon.AppMgr import ToolSvc
    from JetRec import JetRecConf

    if hasattr(algseq, 'JetAlg'+jetcollname):
        jetlog.warning("Jet algorithm \"{0}\" already scheduled! Skipping...".format(jetcollname))
        return

    # Set up the PseudoJetGetter
    pjget = None
    if mods=="":
        if "truthget" in ToolSvc:
            pjget = ToolSvc.truthget
        else:
            pjget = JetRecConf.PseudoJetGetter("truthget",
                                               Label = "Truth",
                                               InputContainer = ToolSvc.truthpartcopy.OutputName,
                                               OutputContainer = "PseudoJetTruth",
                                               GhostScale = 0.0,
                                               SkipNegativeEnergy = True
                                               )
    elif mods=="WZ":
        if "truthwzget" in ToolSvc:
            pjget = ToolSvc.truthwzget
        else:
            pjget = JetRecConf.PseudoJetGetter("truthwzget",
67
                                               Label = "TruthWZ",
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
                                               InputContainer = ToolSvc.truthpartcopywz.OutputName,
                                               OutputContainer = "PseudoJetTruthWZ",
                                               GhostScale = 0.0,
                                               SkipNegativeEnergy = True,
                                               )        
    if pjget:
        ToolSvc += pjget
    else:
        jetlog.error("Unsupported jet collection \"{0}\" requested!".format(jetcollname))

    # Set up the jet builder (no area moments)
    if not "jbld" in ToolSvc:
        ToolSvc += JetRecConf.JetFromPseudojet("jbld")

    # Set up the jet finder
    finder = JetRecConf.JetFinder(jetcollname+"Finder",
                                  JetAlgorithm = "AntiKt",
                                  JetRadius = jetradius,
                                  JetBuilder = ToolSvc.jbld,
                                  GhostArea = 0.01,
                                  PtMin = 5000
                                  )
    ToolSvc += finder

92
93
94
95
    if "truthpartonget" in ToolSvc:
        truthpartonget = ToolSvc.truthpartonget
    else:
        truthpartonget = JetRecConf.PseudoJetGetter("truthpartonget",
96
                                                    Label = "GhostTruthPartons",
97
                                                    InputContainer = ToolSvc.truthpartonscopy.OutputName,
98
99
                                                    OutputContainer = "PseudoJetGhostTruthPartons",
                                                    GhostScale = 1e-40,
100
101
102
103
104
105
106
107
108
                                                    SkipNegativeEnergy = True,
                                                    )
        ToolSvc += truthpartonget


    from ParticleJetTools import ParticleJetToolsConf
    partontruthlabel = ParticleJetToolsConf.Analysis__JetPartonTruthLabel("partontruthlabel")
    ToolSvc += partontruthlabel

109
110
111
    #Now we setup a JetRecTool which will use the above JetFinder
    jetrectool = JetRecConf.JetRecTool(jetcollname+"Rec",
                                       JetFinder = finder,
112
113
114
                                       PseudoJetGetters = [pjget,truthpartonget],
                                       OutputContainer = jetcollname,
                                       JetModifiers = [partontruthlabel]
115
116
117
118
119
120
121
                                       )
    ToolSvc += jetrectool

    algseq += JetRecConf.JetAlgorithm(jetcollname+"Alg", Tools=[jetrectool])

def ScheduleAntiKtTruthWZJets(jetradius,algseq):
    ScheduleAntiKtTruthJets(jetradius,algseq,mods="WZ")