Generate, digitize and reconstruct two signal events with decays Bd->J/Psi(mumu)Kshort

This tutorial shows you how to simulate, digitize and reconstruct signal events so that they can be analyzed in DaVinci. The versions of the code quoted are available in early 2009 and may not be the final version used for the main productions in 2009.

Using Gauss to simulate events

The slides that explain what the program Gauss does are here.

Make an executable version of Gauss available in your lxplus area.

SetupProject Gauss v36r0 --build-env 
getpack Sim/Gauss v36r0
cd Sim/Gauss/cmt
SetupProject Gauss v36r0 
gmake
cd ../options
emacs -nw Gauss-Job.py

Note: you only need to run SetupProject Package vNrP --build-env once to make the directory. However you need to run SetupProject Package vNrP every time you wish to start using the package.

Edit Gauss-Job.py to simulate the required events

Change the 5 in line

nEvts = 5
so that Gauss generates two events.

Change the name of the output file from GaussExample to something that will help you remember the content; for example

idFile = 'BsJPsiKs-2evt'

Save the file and quit emacs

Find out the event type for Bd->J/Psi(mumu)Kshort

Go into the EvtGen web page linked from that of Gauss and find out the event type for the decay you want to simulate. You will want to generate events with the decay products in the acceptance so pick one with DecProdCut in the name (the table is at the end). For Bd->J/Psi(mumu)Kshort this means you will pick either 11144101 or 11144103 depending if you want to simulate your event without or with CP violation.

Run Gauss and look at the monitoring output

To run Gauss you will use the command gaudirun.py and specify that you want

  • a standard Gauss jobs with the latest geometry
  • the event type you have chosen
  • what is specific to your job (i.e. number of events, name of output file,...)
You will do so providing the appropriate arguments to the command, as in this case
gaudirun.py $GAUSSOPTS/Gauss-2008.py $DECFILESROOT/options/11144101.opts $GAUSSOPTS/Gauss-Job.py

Wait while Gauss configures itself, Pythia, EvtGen and GEANT4, then generates two events.

You should see something like this:

Generation.SignalRepe...   INFO *************   Signal counters   ****************
Number of events for generator level cut, before : 2, after : 1
Efficiency of the generator level cut : 0.5 +/- 0.35355
Number of z-inverted events : 1

Number of events for generator particle level cut, before : 2, after : 1
Efficiency of the generator particle level cut : 0.5 +/- 0.35355

Number of signal B0 in sample : 2 [fraction : 1 +/- 0]
Number of signal B~0 in sample : 0 [fraction : 0 +/- 0]

Number of accepted B0 : 0 [fraction : nan +/- nan]
Number of accepted B+ : 0 [fraction : nan +/- nan]
Number of accepted Bs0 : 0 [fraction : nan +/- nan]
Number of accepted b-Baryon : 0 [fraction : nan +/- nan]
Number of accepted Bc+ : 0 [fraction : nan +/- nan]
Number of accepted anti-B0 : 2 [fraction : 1 +/- 0]
Number of accepted B- : 0 [fraction : 0 +/- 0]
Number of accepted anti-Bs0 : 0 [fraction : 0 +/- 0]
Number of accepted anti-b-Baryon : 0 [fraction : 0 +/- 0]
Number of accepted Bc- : 0 [fraction : 0 +/- 0]
Number of accepted (bb) : 0

Number of accepted D0 : 1 [fraction : 1 +/- 0]
Number of accepted D+ : 0 [fraction : 0 +/- 0]
Number of accepted Ds+ : 0 [fraction : 0 +/- 0]
Number of accepted c-Baryon : 0 [fraction : 0 +/- 0]
Number of accepted anti-D0 : 1 [fraction : 1 +/- 0]
Number of accepted D- : 0 [fraction : 0 +/- 0]
Number of accepted Ds- : 0 [fraction : 0 +/- 0]
Number of accepted anti-c-Baryon : 0 [fraction : 0 +/- 0]
Number of accepted (cc) : 0

Number of accepted B : 0 [fraction : 0 +/- 0]
Number of accepted B* : 2 [fraction : 1 +/- 0]
Number of accepted B** : 0 [fraction : 0 +/- 0]

Number of accepted D : 0 [fraction : 0 +/- 0]
Number of accepted D* : 2 [fraction : 1 +/- 0]
Number of accepted D** : 0 [fraction : 0 +/- 0]
showing the number of "signal" particle types created in each event. Note this is not the end of the log file and may have scrolled off the top of the screen. In general send the output of any program to a log file to be able to read it later.

Look at the histogram file produced called something like BsJPsiKs-2evt-histos.root and check there are hits in the VELO, RICH etc.

Digitize the two decays made in Gauss

Make an executable version of Boole in your lxplus area

SetupProject Boole v17r1 --build-env
getpack Digi/Boole v17r1
cd Digi/Boole/cmt
SetupProject Boole v17r1
gmake
cd ../options

Edit 2008-Files.py to setup Boole to digitize the events generated earlier

Change the input data file to the BsJPsiKs-2evt.sim file you generated earlier in Gauss. Do not bother with spillover for these events.

So make the following changes to 2008-Files.py

datasetname = 'BsJPsiKs-2evt'
EventSelector().Input = ["DATAFILE='PFN:$HOME/cmtuser/Gauss_v36r0/Sim/Gauss/v30r5/options/" + datasetName + ".sim' TYP='POOL_ROOTTREE' OPT='READ'"]

Run Boole and check the output

Run Boole with gaudirun.py Boole-2008.py 2008-Files.py

Look at the histograms produced by Boole in BsJPsiKs-2evt-histos.root, check there are entries in the histograms.

Use Brunel to reconstruct the events digitized with Boole

Use the same instructions as for Gauss and Boole to setup, getpack and compile Brunel v34r0, note the location of the Brunel package for getpack is Rec/Brunel.

In the options directory set the input file to the output produced by Boole, i.e. 11144103-2ev-20090109.digi, in 2008-Files.py.

Then run Brunel with gaudirun.py Brunel-2008-MC.py 2008-Files.py

Again check the output log file and the histograms produced.

-- DavidHutchcroft - 09 Jan 2009

Edit | Attach | Watch | Print version | History: r19 | r12 < r11 < r10 < r9 | Backlinks | Raw View | Raw edit | More topic actions...
Topic revision: r10 - 2009-02-05 - MarcoCattaneo
 
    • Cern Search Icon Cern Search
    • TWiki Search Icon TWiki Search
    • Google Search Icon Google Search

    LHCb All webs login

This site is powered by the TWiki collaboration platform Powered by PerlCopyright & 2008-2019 by the contributing authors. All material on this collaboration platform is the property of the contributing authors.
Ideas, requests, problems regarding TWiki? Send feedback