Using DaVinci to fill an nTuple

This page is designed to give a brief outline of how to use DaVinci to fill nTuples from data or Monte Carlo. For more details please refer to the DaVinci tutorials from LHCb: https://twiki.cern.ch/twiki/bin/view/LHCb/DaVinciTutorial My full options file is at the bottom of the page so you can use it for testing.

Setting up DaVinci

If you are working from eprexa, you first need to set up the LHCb environment:
 source /home/lhsoft/LbLogin.sh

Then set up the DaVinci project

 SetupProject DaVinci v36r0 

Note that you need to make sure the version of DaVinci you are using contains the version of Stripping you are using. For example, if you are using Stripping 20, you need to use a version of DaVinci higher than v32. If you don't specify a version, then the latest version available will be set up.

Using DecayTreeTuple

The tool used to fill nTuples is DecayTreeTuple. Create and options file (.py) and start and instance of DecayTreeTuple:
from Configurables import DecayTreeTuple
from DecayTreeTuple.Configuration import *
tuple=DecayTreeTuple() 

You then need to specify the input location for the tuple. This is where the output from the stripping line is stored in the file:

 tuple.Inputs=["/Event/Bhadron/Phys/B2XEtaB2etapKSLLLine/Particles"]
If you are producing and nTuple which is empty, it is likely that your input location is wrong.

You now need to add a list of tools to the tuple. This specifies exactly which information is filled in the nTuple. Below is an example. Please see the tutorials for a full list of TupleTools.

 tuple.ToolList +=  ["TupleToolGeometry" , "TupleToolKinematic", "TupleToolEventInfo"] 
You need to specify the decay you are looking at. ^ means that information from the TupleTools will be filled for that particle. By default the head particle will always be filled (B0).
tuple.Decay="[B0 -> ^(eta_prime -> ^pi+ ^pi- ^gamma) ^(KS0 -> ^pi+ ^pi-)]CC"
You can also specify some branches, so that different information can be filled for different particles. For example, if you want different information for the B0, use the branch below.
tuple.Branches={"B0":"B0:[B0 -> (eta_prime -> pi+ pi- gamma) (KS0 -> pi+ pi-)]CC" }
tuple.addTool( TupleToolDecay, name = "B0" )
In order to fill information for calculating trgger efficiencies, you will need the TISTOS tool. You need to specify VerboseL0,Hlt1,Hlt2=True to ensure that the information is filled properly. TriggerList should contain all the trigger lines you are interested in.
 from Configurables import TupleToolTISTOS
tistos = tuple.B0.addTupleTool(TupleToolTISTOS, name="TupleToolTISTOS")
tistos.VerboseL0=True
tistos.VerboseHlt1=True
tistos.VerboseHlt2=True
tistos.TriggerList=["L0PhotonDecision", "L0ElectronDecision", "Hlt1TrackAllL0Decision", "Hlt2Topo2BodyBBDTDecision"] 

Another useful tool is the EventTuple. This is similar to DecayTreeTuple, but will fill information for every event processed, not just the ones which pass your stripping line. This is useful for calculating efficiencies, and you only need to include one tool.

 etuple=EventTuple()
etuple.ToolList=["TupleToolEventInfo"] 

At the end of you options file, you need to specify the options for the DaVinci configurable. Your tuples need to be included in the UserAlgorithms, and you need to make sure the specified DataType? is correct for the data you are looking at. EvtMax is the number of events you want to process. -1 asks to process all events.

 DaVinci().UserAlgorithms+[tuple,etuple]
DaVinci().InputType='DST'
DaVinci().TupleFile="davinci.root"
DaVinci().HistogramFile="davinci-histos.root"
DaVinci().DataType="2012"
DaVinci().EvtMax=-1
DaVinci().PrintFreq=1000
DaVinci().MoniSequence=[tuple]
DaVinci().Simulation=False 

Running on Monte Carlo

DecayTreeTuple works slightly differently on Monte Carlo than data. Firstly you need to set the simulation flag in davinci to true, and then specify the database tags used to create the Monte Carlo. You can find these out from the bookkeeping.
DaVinci().Simulation=True
DaVinci().DDDBtag='MC11-20111102'
DaVinci().CondDBtag='sim-20111111-vc-mu100' 
The Input location for DecayTreeTuple seems to be slightly different for Monte Carlo. You should start at /Phys:
 tuple.Inputs=["/Phys/B2XEtaB2etapKSLLLine/Particles"]
It is likely that the Monte Carlo you have has not been run with the correct stripping. You can run a different version of the Stripping with the following options:
    from StrippingConf.Configuration import StrippingConf, StrippingStream
    from StrippingSettings.Utils import strippingConfiguration
    from StrippingArchive.Utils import buildStreams
    from StrippingArchive import strippingArchive

    stripping="stripping20"
    config=strippingConfiguration(stripping)
    archive=strippingArchive(stripping)
    streams=buildStreams(stripping=config,archive=archive)

    MyStream = StrippingStream("MyStream")
    MyLines = ["StrippingB2XEtaB2etapKSLLLine]

    for stream in streams:
        for line in stream.lines:
            if line.name() in MyLines:
                MyStream.appendLines( [ line ])

    from Configurables import ProcStatusCheck
    filterBadEvents=ProcStatusCheck()

    sc=StrippingConf( Streams= [ MyStream ],
                      MaxCandidates = 2000,
                      AcceptBadEvents = False,
                      BadEventSelection = filterBadEvents)
   
   DaVinci().appendToMainSequence([sc.sequence()]) 

Running on a MicroDST?

MicroDST? (or uDST) is a file format designed to reduce the amount of data stored, and so reduce the strain of computing resources. There is certain information that is not available on a uDST, and some things are stored in a slightly different way, but the aim is to make it as close to processing a DST as possible. This means you can run DecayTreeTuple on uDSTs in the same way. There will be some TupleTools which aren't available, for example TupleToolTagging or TupleToolTrackIsolation. There are also some files which DecayTreeTuple will look for but don't exist, producing an error. You therefore need to kill these nodes:
from Configurables import EventNodeKiller
eventNodeKiller = EventNodeKiller('DAQkiller')
eventNodeKiller.Nodes = ['DAQ','pRec']
MySequencer.Members+=[eventNodeKiller] 

Full Options File

############################################ 
# DaVinci Options for filling ntuple with  #
# B0 -> KS0 (eta' -> pi pi gamma)          #
# for data or Monte Carlo                  #
############################################

from Gaudi.Configuration import *
from Configurables import DaVinci

#set your options for the job
simulation=True

#If its simulation, the monte carlo sample you are running on
sample='BdetapKs-pipig'
#If its data, the integrated luminosity from book-keeping
#2012
#lumi='1fb'
#2011: 480pb magup 600pb magdown
if(simulation):
    lumi="lumi"
else:
    lumi='1fb' #
mdst=False
dataType='2012'
magnet='magup'
stripping='stripping20r0p2'
sStream='Bhadron'
#Your stripping line (without the leading 'Stripping')
StrippingLine='B2XEtaB2etapKSDDLine'

head='B0'
daughter='etap'

#decay descriptor
#decay='[Lambda_b0 -> ^(eta_prime -> ^gamma ^pi+ ^pi-) ^(Lambda0 -> ^p+ ^pi-)]CC'
decay='B0 -> ^(eta_prime -> ^gamma ^pi+ ^pi-) ^(KS0 -> ^pi+ ^pi-)'

print 'decay: '+decay

from Configurables import GaudiSequencer
MySequencer = GaudiSequencer('Sequence')

#Setting the output file names, the location for the tuple and the database tags if MC
if(simulation):
    name="MC"+dataType+"-"+magnet+"-"+sample+"-"+StrippingLine+"-"+stripping+"-filtered_"+trigger
    location="Phys/"+StrippingLine+"/Particles"
        
    DaVinci().DDDBtag='Sim08-20130503-1'
    DaVinci().CondDBtag='Sim08-20130503-1-vc-mu100'
else:
    name="data"+dataType+"-"+magnet+"-"+sample+"-"+sStream+"-"+StrippingLine+"-"+stripping+"-"+lumi
    location="/Event/"+sStream+"/Phys/"+StrippingLine+"/Particles"

if(mdst):
    from Configurables import EventNodeKiller
    eventNodeKiller = EventNodeKiller('DAQkiller')
    eventNodeKiller.Nodes = ['DAQ','pRec']
    MySequencer.Members+=[eventNodeKiller]   

if(simulation):
    #######run stripping20 on MC11#################
    from StrippingConf.Configuration import StrippingConf, StrippingStream
    from StrippingSettings.Utils import strippingConfiguration
    from StrippingArchive.Utils import buildStreams
    from StrippingArchive import strippingArchive

    config=strippingConfiguration(stripping)
    archive=strippingArchive(stripping)
    streams=buildStreams(stripping=config,archive=archive)

    MyStream = StrippingStream("MyStream")
    MyLines = ["Stripping"+StrippingLine]

    for stream in streams:
        for line in stream.lines:
            if line.name() in MyLines:
                MyStream.appendLines( [ line ])

    from Configurables import ProcStatusCheck
    filterBadEvents=ProcStatusCheck()

    sc=StrippingConf( Streams= [ MyStream ],
                      MaxCandidates = 2000,
                      AcceptBadEvents = False,
                      BadEventSelection = filterBadEvents,
                      HDRLocation = "SomeNonExistingLocation") #only if running the same stripping again
                      
    DaVinci().appendToMainSequence([sc.sequence()])

#Creating the ntuple
from Configurables import DecayTreeTuple
from DecayTreeTuple.Configuration import *

tuple=DecayTreeTuple()
tuple.Inputs=[location]
tuple.ToolList += [
    "TupleToolGeometry"
    , "TupleToolDira"
    , "TupleToolAngles"
    , "TupleToolPid"
    , "TupleToolKinematic"
    , "TupleToolPropertime"
    , "TupleToolPrimaries"
    , "TupleToolEventInfo"
    , "TupleToolTrackInfo"
    , "TupleToolVtxIsoln"
    , "TupleToolPhotonInfo"
    , "TupleToolMCTruth"
    , "TupleToolMCBackgroundInfo"
    , "TupleToolTagging"
    ]
if(mdst!=True): tuple.ToolList+=["TupleToolTrackIsolation"]

tuple.Decay=decay

import re

brdecay=re.sub('\^','',decay)
print 'brdecay: '+brdecay


tuple.Branches={ "B0" : brdecay }

if(head=="B0"): 
elif(head=="Lb"):
    tuple.Branches={ "B0" : brdecay }
else:
    print "Unknown Head"
    
tuple.addTool(TupleToolDecay, name = "B0")

from Configurables import TupleToolDecayTreeFitter

tuple.B0.addTool(TupleToolDecayTreeFitter("PVFit"))
tuple.B0.PVFit.Verbose = True
tuple.B0.PVFit.constrainToOriginVertex = True

if(head=="B0"):
    if(daughter=="etap"):
        tuple.B0.PVFit.daughtersToConstrain = [ "eta_prime", "KS0" ]
    elif(daughter=="eta"):
        tuple.B0.PVFit.daughtersToConstrain = [ "eta", "KS0" ]
    else:
        Print("Unknown daughter")

elif(head=="Lb"):
    if(daughter=="etap"):
        tuple.B0.PVFit.daughtersToConstrain = [ "eta_prime", "Lambda0" ]
    elif(daughter=="eta"): 
        tuple.B0.PVFit.daughtersToConstrain = [ "eta", "Lambda0" ]
    else:
        print "Unknown Daughter"

else:
    print "Unkown Head"

tuple.B0.ToolList+=["TupleToolDecayTreeFitter/PVFit"]

from Configurables import TupleToolTISTOS
tistos = tuple.B0.addTupleTool(TupleToolTISTOS, name="TupleToolTISTOS")
tistos.VerboseL0=True
tistos.VerboseHlt1=True
tistos.VerboseHlt2=True
tistos.TriggerList=["L0PhotonDecision", "L0ElectronDecision", "Hlt1TrackPhotonDecision", "Hlt1TrackAllL0Decision","Hlt1TrackMuonDecision", "Hlt1TrackForwardPassThroughDecision","Hlt1TrackForwardPassThroughLooseDecision", "Hlt1SingleElectronNoIPDecision","L0HadronDecision","L0LocalPi0Decision","L0GlobalPi0Decision","L0MuonDecision","Hlt2Topo2BodyBBDTDecision","Hlt2Topo3BodyBBDTDecision","Hlt2Topo4BodyBBDTDecision", "Hlt2RadiativeTopoTrackTOSDecision", "Hlt2RadiativeTopoPhotonL0Decision","Hlt2TopoRad2BodyBBDTDecision","Hlt2TopoRad2plus1BodyBBDTDecision","Hlt2Topo2BodySimpleDecision","Hlt2Topo3BodySimpleDecision","Hlt2Topo4BodySimpleDecision"]

if(simulation):
    etuple=EventTuple()
    etuple.ToolList=["TupleToolEventInfo"]

    MySequencer.Members.append(etuple)
    
    from Configurables import MCDecayTreeTuple
    mctuple=MCDecayTreeTuple()
    mctuple.ToolList+=["MCTupleToolKinematic","MCTupleToolReconstructed","MCTupleToolHierarchy","MCTupleToolDecayType","MCTupleToolPID"]

    if(head=="B0"):
        decay2=re.sub('\^\(\KS0 \-\> \^pi\+ \^pi\-\)','^KS0',decay)
        mcdecay=re.sub('\)',' {,gamma}{,gamma}{,gamma}{,gamma})',decay2)      
        mcdecay = '['+mcdecay+']cc'
    elif(head=="Lb"):
        decay2=re.sub('\^\(Lambda0 \-\> \^p\+ \^pi\-\)','^Lambda0',decay)
        mcdecay=re.sub('\)',' {,gamma}{,gamma}{,gamma}{,gamma})',decay2) 

    print 'mcdecay= '+mcdecay

    mctuple.Decay=mcdecay    

MySequencer.Members.append(tuple)

if(mdst):
    DaVinci().InputType='MDST'
    DaVinci().UserAlgorithms+=[MySequencer]    
else:
    DaVinci().InputType='DST'
    DaVinci().UserAlgorithms+=[MySequencer]

DaVinci().TupleFile=name+".root"
DaVinci().HistogramFile=name+"-histos.root"
DaVinci().DataType=dataType
DaVinci().EvtMax=500
DaVinci().PrintFreq=1000
DaVinci().MoniSequence=[tuple]
DaVinci().Simulation=simulation

-- JamesMccarthy - 15 Nov 2012


This topic: General > Lhcb?  > DaVinci
Topic revision: r6 - 12 Sep 2014 - _47C_61UK_47O_61eScience_47OU_61Birmingham_47L_61ParticlePhysics_47CN_61james_32mccarthy
 
This site is powered by the TWiki collaboration platform Powered by Perl This site is powered by the TWiki collaboration platformCopyright © by the contributing authors. All material on this collaboration platform is the property of the contributing authors.
Ideas, requests, problems regarding TWiki? Send feedback