# 21.4.6. Vertex fitting¶

## Introduction¶

In the broadest sense, we call vertex fitting a technique in which one uses prior knowledge on the nature of a decay to improve the measurement of its observables. The fits we are going to perform are of two main types:

• Geometric Fitting: We use the fit to determine the decay vertex of the particle. Usually this is done by fitting together the tracks of its charged decay products, which we know originate from a common point. Additional information could be available — for example, if the particle is short lived, we can improve this by adding an IP constraint, i.e. fit the beam spot together with the tracks. If there’s only one track, using the beam spot is the only way to obtain a vertex.

Warning

If no vertex fit is performed, the corresponding variables for the vertex position will not be filled.

• Kinematic Fitting: We use the fit to improve our knowledge of the particle kinematics. By default, composite particle kinematics are built off the decay products using 4-momentum conservation. If the particle we are reconstructing has a well defined mass (either stable, or a narrow resonance) it might make sense to apply a mass constraint to help reject combinatorial background.

Warning

If you apply a mass constraint, the invariant mass will be fixed to the nominal mass. This is problematic if you then want to use this variable, for example if you want to fit a peak. In that case, make sure you save the pre-fit mass separately.

Note

Several fitters exist. For this exercise we will focus on KFit which is the most basic one.

Exercise

Locate the documentation for vertex fitting functions and find KFit.

Hint

Use the search bar.

Solution

You can find it here: Vertex-fitting convenience functions.

## Basic Fitting¶

This lesson assumes you successfully reconstructed your $$B \to J/\Psi(\to e^+e^-)K_s(\to \pi^+\pi^+)$$ decay following the previous exercises. Now suppose you are interested in reconstructing the $$B$$ decay vertex position using a fit (for example, you’re trying to do a time-dependent CPV study).

Question

Which particles do you need to fit in order to reconstruct the $$B$$ vertex?

Hint

It can’t be the $$B$$ itself: out of its daughters, neither the $$J/\Psi$$ nor the $$K_s$$ are charged tracks.

You must fit the $$J/\Psi \to e^+e^-$$ vertex. The $$J/\Psi$$ is short-lived and therefore its vertex is a good approximation of the $$B$$. Meanwhile, the $$K_s$$ can decay several cm away from where it is produced.

Exercise

Call the fit function with the correct parameters and save the output. Include the true vertex position from MC for comparison.

Hint

Look up the variable collections for vertices. Don’t forget to import the vertex module!

Solution

import vertex
...
vertex.kFit("J/psi:ee", conf_level=0.0, path=main)
...
jpsi_ks_vars += vc.vertex + vc.mc_vertex


You can also set the confidence level to -1, which means failed fits will be included. The fit p-value is saved as part of the vc.vertex collection.

Exercise (optional)

Fit the $$K_s$$ as well. How does its flight length compare to the $$J/\Psi$$?

Exercise (optional)

Look up the documentation for TreeFitter and fit the whole $$B \to J/\Psi(\to e^+e^-)K_s(\to \pi^+\pi^+)$$ decay chain at once.

## Tag Vertex Fitting¶

Since $$B$$ mesons are produced in pairs, for every signal candidate we reconstruct, there will also be another (the “tag”) which is not explicitly reconstructed.

We might be interested in knowing the decay position of this meson without placing any requirements on its decay. This is done using the TagV module.

TagV performs a geometric fit over the tracks in the ROE to determine the tag decay vertex. However, not all tracks will necessarily come from the tag itself; consider for example our signal, where the pion tracks originate from a long lived kaon vertex. TagV is designed to iteratively downweight those tracks, ultimately excluding them from the fit.

Exercise

Locate the TagV documentation.

Solution

It’s in the same page as KFit.

Question

By default, TagV only uses tracks with PXD hits. Why?

Solution

Those tracks provide the best resolution close to the interaction point. As a bonus, this selection rejects tracks from displaced vertices.

Exercise

Call the TagV module and save the output.

Hint

In order to reinforce the fit, an IP constraint is applied to the TagV. If the signal is fully reconstructed, this condition can be relaxed along the signal $$B$$ flight direction.

Solution

vertex.TagV("B0", constraintType="tube", path=main)
...
b_vars += vc.tag_vertex + vc.mc_tag_vertex


## Conclusion and Plotting¶

Congratulations! Your steering file is ready! Time to run it and check the results.

Exercise

Run the steering file.

Solution

Your steering file should look like this:

  1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 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 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 #!/usr/bin/env python3 import sys import basf2 as b2 import modularAnalysis as ma import stdV0s import flavorTagger as ft from variables import variables as vm # shorthand for VariableManager import variables.collections as vc import variables.utils as vu import vertex # get input file number from the command line filenumber = sys.argv[1] # set analysis global tag (needed for flavor tagging) b2.conditions.prepend_globaltag("analysis_tools_release-04-02") # create path main = b2.Path() # load input data from mdst/udst file ma.inputMdstList( environmentType="default", filelist=[b2.find_file(f"starterkit/2021/1111540100_eph3_BGx0_{filenumber}.root", "examples")], path=main, ) # fill final state particle lists ma.fillParticleList( "e+:uncorrected", "electronID > 0.1 and dr < 0.5 and abs(dz) < 2 and thetaInCDCAcceptance", path=main, ) stdV0s.stdKshorts(path=main) # apply Bremsstrahlung correction to electrons vm.addAlias( "goodFWDGamma", "passesCut(clusterReg == 1 and clusterE > 0.075)" ) vm.addAlias( "goodBRLGamma", "passesCut(clusterReg == 2 and clusterE > 0.05)" ) vm.addAlias( "goodBWDGamma", "passesCut(clusterReg == 3 and clusterE > 0.1)" ) vm.addAlias( "goodGamma", "passesCut(goodFWDGamma or goodBRLGamma or goodBWDGamma)" ) ma.fillParticleList("gamma:brems", "goodGamma", path=main) ma.correctBrems("e+:corrected", "e+:uncorrected", "gamma:brems", path=main) vm.addAlias("isBremsCorrected", "extraInfo(bremsCorrected)") # combine final state particles to form composite particles ma.reconstructDecay( "J/psi:ee -> e+:corrected e-:corrected ?addbrems", cut="dM < 0.11", path=main, ) # perform vertex fit of J/psi candidates vertex.kFit("J/psi:ee", conf_level=0.0, path=main) # combine J/psi and KS candidates to form B0 candidates ma.reconstructDecay( "B0 -> J/psi:ee K_S0:merged", cut="Mbc > 5.2 and abs(deltaE) < 0.15", path=main, ) # match reconstructed with MC particles ma.matchMCTruth("B0", path=main) # build the rest of the event ma.buildRestOfEvent("B0", fillWithMostLikely=True, path=main) track_based_cuts = "thetaInCDCAcceptance and pt > 0.075 and dr < 5 and abs(dz) < 10" ecl_based_cuts = "thetaInCDCAcceptance and E > 0.05" roe_mask = ("my_mask", track_based_cuts, ecl_based_cuts) ma.appendROEMasks("B0", [roe_mask], path=main) # call flavor tagging ft.flavorTagger(["B0"], path=main) # fit B vertex on the tag-side vertex.TagV("B0", fitAlgorithm="Rave", path=main) # perform best candidate selection b2.set_random_seed("Belle II StarterKit") ma.rankByHighest("B0", variable="random", numBest=1, path=main) # Create list of variables to save into the output file b_vars = [] standard_vars = vc.kinematics + vc.mc_kinematics + vc.mc_truth b_vars += vc.deltae_mbc b_vars += standard_vars # ROE variables roe_kinematics = ["roeE()", "roeM()", "roeP()", "roeMbc()", "roeDeltae()"] roe_multiplicities = [ "nROE_Charged()", "nROE_Photons()", "nROE_NeutralHadrons()", ] b_vars += roe_kinematics + roe_multiplicities # Let's also add a version of the ROE variables that includes the mask: for roe_variable in roe_kinematics + roe_multiplicities: # e.g. instead of 'roeE()' (no mask) we want 'roeE(my_mask)' roe_variable_with_mask = roe_variable.replace("()", "(my_mask)") b_vars.append(roe_variable_with_mask) b_vars += ft.flavor_tagging b_vars += vc.tag_vertex + vc.mc_tag_vertex # Variables for final states (electrons, positrons, pions) fs_vars = vc.pid + vc.track + vc.track_hits + standard_vars b_vars += vu.create_aliases_for_selected( fs_vars + ["isBremsCorrected"], "B0 -> [J/psi -> ^e+ ^e-] K_S0", prefix=["ep", "em"], ) b_vars += vu.create_aliases_for_selected( fs_vars, "B0 -> J/psi [K_S0 -> ^pi+ ^pi-]", prefix=["pip", "pim"] ) # Variables for J/Psi, KS jpsi_ks_vars = vc.inv_mass + standard_vars jpsi_ks_vars += vc.vertex + vc.mc_vertex b_vars += vu.create_aliases_for_selected(jpsi_ks_vars, "B0 -> ^J/psi ^K_S0") # Add the J/Psi mass calculated with uncorrected electrons: vm.addAlias( "Jpsi_M_uncorrected", "daughter(0, daughterCombination(M,0:0,1:0))" ) b_vars += ["Jpsi_M_uncorrected"] # Also add kinematic variables boosted to the center of mass frame (CMS) # for all particles cmskinematics = vu.create_aliases( vc.kinematics, "useCMSFrame({variable})", "CMS" ) b_vars += vu.create_aliases_for_selected( cmskinematics, "^B0 -> [^J/psi -> ^e+ ^e-] [^K_S0 -> ^pi+ ^pi-]" ) vm.addAlias( "withBremsCorrection", "passesCut(passesCut(ep_isBremsCorrected == 1) or passesCut(em_isBremsCorrected == 1))", ) b_vars += ["withBremsCorrection"] # Save variables to an output file (ntuple) ma.variablesToNtuple( "B0", variables=b_vars, filename="Bd2JpsiKS.root", treename="tree", path=main, ) # Start the event loop (actually start processing things) b2.process(main) # print out the summary print(b2.statistics) 

You can now plot some relevant vertex variables. In general, the choice would depend on what you need for your analysis. A few examples would include:

• Vertex position in various coordinates, such as dz and dr.

• P-value of the fit.

• Resolution of the vertex fit ($$\sigma(x)/x$$) where x is each of the above variables.

• Pull ($$(x-x(MC)/\sigma(x)$$).

As an exercise we will focus on the first two.

Exercise

Plot the $$J/\Psi$$ vertex position and compare it with the true value. Plot the p-value distribution of the fit.

Hint: Variable names

You can either take another look at the variable collections that you included above, or you load your dataframe and then take a look at its columns print(list(df.columns)).

Hint: Plot ranges

Plotting was already discussed in The Rest of Event (ROE). For the sake of this exercise, remember we already set the minimum p-value of our fits to 0, so failed fits will not be included and you can plot it in the [0,1] interval.

Should you have changed that, failed fits will be included with a p-value of -1; in this case, make sure to change your plotting range accordingly to [-1,1].

Solution

  1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 import matplotlib.pyplot as plt from root_pandas import read_root plt.style.use("belle2") df = read_root("Bd2JpsiKS.root") m_bins = 50 # number of bins for the histograms of both plots # Z position fig, ax = plt.subplots(figsize=(8, 6)) m_range = [-0.1, 0.1] ax.set_xlim(left=-0.1, right=0.15) ax.hist(df["Jpsi_dz"], bins=m_bins, range=m_range, label=r"$J/\psi$ vertex") ax.hist( df["Jpsi_mcDecayVertexZ"], histtype="step", lw=2, color="black", linestyle="--", bins=m_bins, range=m_range, label=r"$J/\psi$ vertex(MC)", ) ax.set_xlabel("dz[cm]") ax.set_ylabel("Events") ax.legend() fig.savefig("vertex_jpsi_dz.svg") # P-value fig, ax = plt.subplots(figsize=(8, 6)) m_range = [0, 1] ax.set_xlim(left=-0.05, right=1.05) ax.hist( df["Jpsi_chiProb"], bins=m_bins, range=m_range, label=r"$J/\psi$ vertex", ) ax.set_yscale("log") # set a logarithmic scale in the y-axis ax.set_xlabel("p-value") ax.set_ylabel("Events") ax.legend() fig.savefig("vertex_pValue.svg") 

Fig. 21.25 Distribution of the fitted vertex position in Z

Fig. 21.26 Distribution of the fit p-values.

Exercises (optional)

• Compare the $$J/\Psi$$ and Tag vertex positions and show that they are both compatible with being $$B$$ vertices.

• If you’ve fit the $$K_s$$ vertex, compare its radial position with the $$J/\Psi$$. Is this what you expect?

Key points

• Use KFit to fit simple vertices.

• Think carefully which vertex you need to fit, and whether it will need additional constraints.

• Study the documentation if you need a different functionality, such as TreeFitter to fit complex trees.

• Use TagV to reconstruct a vertex from the ROE.

Stuck? We can help!

If you get stuck or have any questions to the online book material, the #starterkit-workshop channel in our chat is full of nice people who will provide fast help.

Refer to Collaborative Tools. for other places to get help if you have specific or detailed questions about your own analysis.

Improving things!

If you know how to do it, we recommend you to report bugs and other requests with JIRA. Make sure to use the documentation-training component of the Belle II Software project.

If you just want to give very quick feedback, use the last box “Quick feedback”.

Please make sure to be as precise as possible to make it easier for us to fix things! So for example:

• typos (where?)

• missing bits of information (what?)

• bugs (what did you do? what goes wrong?)

• too hard exercises (which one?)

• etc.

If you are familiar with git and want to create your first pull request for the software, take a look at How to contribute. We’d be happy to have you on the team!

Quick feedback!

Authors of this lesson

Francesco Tenchini