Skip to content

Object definitions

Baseline object definitions

We now move to the MyAnalysisName::ProcessEvent(AnalysisEvent *event) method, starting with the baseline object definitions. Baseline objects are selected with relatively loose identification criteria and are used for resolving possible reconstruction ambiguities.

For each object type, there is a corresponding getter method available from the TruthEvent class. For example, in order to get electrons, we need to call

1
virtual AnalysisObjects getElectrons(float ptCut,float etaCut,int isolation);

This method returns an std::vector of AnalysisObjects, which are subclasses of TLorentzVector with some convenience features.

Getters and setters

Have a look at TruthEvent.h for a list of all object getter methods. As you would expect, you can get e.g. jets using virtual AnalysisObjects getJets(float ptCut, float etaCut, int btag);.

Electrons

As can be seen from the getElectrons() method declaration, you can specify a p_\mathrm{T} as well as \eta requirement that should be applied to your electrons. In addition, you can specify the isolation and identification criteria. Looking at the AnalysisElectronID enum declared in AnalysisClass.h, you can see which working points are available:

 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
enum AnalysisElectronID {
  EVeryLooseLH = 1 << 0,
  ELooseLH = 1 << 1,
  EMediumLH = 1 << 2,
  ETightLH = 1 << 3,
  ELooseBLLH = 1 << 4,
  EIsoGradientLoose = 1 << 8,
  EIsoBoosted = 1 << 9,
  EIsoFixedCutTight = 1 << 10,
  EIsoLooseTrack = 1 << 11,
  EIsoLoose = 1 << 12,
  EIsoGradient = 1 << 13,
  EIsoFixedCutLoose = 1 << 14,
  EIsoFixedCutTightTrackOnly = 1 << 15,
  ED0Sigma5 = 1 << 16,
  EZ05mm = 1 << 17,
  EIsoFCTight = 1<< 18,
  EIsoFCTightTrackOnly = 1<< 19,
  EIsoTightTrackOnly = 1<<19,
  EIsoFCHighPtCaloOnly = 1 << 20,
  EIsoTightTrackOnly_FixedRad = 1 << 21,
  EIsoFCLoose = 1 << 22,
  EGood = EVeryLooseLH | ELooseLH | EMediumLH | ETightLH | ELooseBLLH | ED0Sigma5 | EZ05mm,
  EIsoGood = EGood | EIsoGradientLoose | EIsoBoosted | EIsoFixedCutTight | EIsoFCTight | EIsoLooseTrack | EIsoLoose | EIsoGradient | EIsoFixedCutLoose | EIsoFixedCutTightTrackOnly | EIsoFCTightTrackOnly | EIsoFCHighPtCaloOnly | EIsoTightTrackOnly_FixedRad | EIsoFCLoose
};

Effects of object ID working points

Note that the choice of object ID does not have an effect when running on unsmeared truth. Only when smearing is enabled, or when running over reco-level xAODs do the different ID working points have an effect.

Since this example follows the 1Lbb analysis, we will use baseline electrons that have a p_\mathrm{T} > 7 GeV , \eta < 2.47 and are identified using the LooseAndBLayerLLH, which corresponds to ELooseBLLH in SimpleAnalysis.

Thus, in order to get the right electrons, we have to call:

1
auto baselineElectrons = event->getElectrons(7.0, 2.47, ELooseBLLH);

Muons

The procedure is very similar for muons, where we use the virtual AnalysisObjects getMuons(float ptCut, float etaCut, int isolation) method. Baseline muons are required to have p_\mathrm{T} > 6 GeV , \eta < 2.70 and satisfy the Medium identification criterion. Additionally, cosmic (MuNotCosmic) and bad (MuQoPSignificance) muons are vetoed. The longitudinal impact parameter relative to the PV is required to satsify |z0sinθ| < 0.5 mm (MuZ05mm).

1
auto baselineMuons = event->getMuons(6.0, 2.70, MuMedium | MuNotCosmic | MuZ05mm | MuQoPSignificance);

Looking up available working points

For a list of all available working points, please consult the AnalysisMuonID enum in AnalysisClass.h.

Jets

We'll use baseline jets where p_\mathrm{T} > 20 GeV and \eta < 4.50. By now, you probably know how this works, so I'll leave this uncommented:

1
auto baselineJets = event->getJets(20.0, 4.5);

Missing transverse energy

Next, we will get the E_\mathrm{T}^\mathrm{miss} with:

1
2
auto met_Vect  = event->getMET(); //actually returns the 4-vector
float met = met_Vect.Et();

MC event weights

In case we want to use it, we are also getting the MC event weights from the xAOD::Event with:

1
auto weights = event->getMCWeights();
SimpleAnalysis automatically uses the first weight in the list of xAOD::Event MC weights, for filling into histograms, ntuple branches and computing acceptances. This behaviour can be changed through the command line.

Overlap removal

In order to resolve reconstruction ambiguities between electrons, muons and jets, an overlap removal procedure is applied on baseline objects. SimpleAnalysis allows easy overlap removal procedure through the AnalysisClass::overlapRemoval() method:

1
static AnalysisObjects overlapRemoval(const AnalysisObjects &cands, const AnalysisObjects &others, float deltaR, int passId=0);
This method returns all the AnalysisObjects &cands candidates that pass the passID criterion and do not overlap with any of the other AnalysisObjects &others AnalysisObjects within a distance of float deltaR.

Custom radius function

One can also pass a radius function to AnalysisClass::overlapRemoval() instead of a fixed float deltaR parameter. We will use this very soon.

Check your overlap removal, at least twice!

The following overlap removal procedure describes the default overlap removal that needs to be done in SimpleAnalysis when the analysis sticks to the default overlap removal procedure as defined in SUSYTools and the OverlapRemovalTool as of April 2020. This might be significantly different for your analysis, but it is still essential that you get this right when implementing your event selection into SimpleAnalysis.

Let's start by defining a variable cone size ΔR = min(0.4, 0.04 + 10 GeV/p_\mathrm{T}) that shrinks with increasing object p_\mathrm{T}:

1
2
3
auto radiusCalcLep = [] (const AnalysisObject& lep,const AnalysisObject&) {
    return (0.04 + 10/lep.Pt()) > 0.4 ? 0.4 : 0.04 + 10/lep.Pt();
  };
Usually, the first step in the overlap removal, is to remove electrons that share the same inner detector track with a muon. Since ID tracks are not included in the truth content, we will just remove electrons close to another muon with a fixed ΔR = 0.01.
1
baselineElectrons = overlapRemoval(baselineElectrons, baselineMuons, 0.01);
Next, jets within ΔR = 0.2 of an electron are removed, followed by removal of electrons close to jets using the shrinking ΔR cone:
1
2
baselineJets = overlapRemoval(baselineJets, baselineElectrons, 0.2);
baselineElectrons = overlapRemoval(baselineElectrons, baselineJets, radiusCalcLep);
Finally, the same procedure is repeated for muons and jets. Jets within ΔR = 0.2 of a muon are removed, except if the jet has more than 3 ID tracks. This is followed by the removal of muons close to a remaining jet using again the shrinking ΔR cone:
1
2
baselineJets = overlapRemoval(baselineJets, baselineMuons, 0.2, LessThan3Tracks);
baselineMuons = overlapRemoval(baselineMuons, baselineJets, radiusCalcLep);

Applying constraints on selections

Note how we use the float passId argument of the AnalysisClass::overlapRemoval() method to only apply the overlap removal to jets that satisfy the condition.

b-tag aware overlap removal

Keep in mind that any overlap removal involving jets might actually do different things for light jets than for b-jets. Check your analysis implementation!

Check your overlap removal!

Check what has been implemented specifically for your analysis! Overlap removal is a topic that is often not properly (or correctly) documented, so you might have to dig through the actual configuration and code of the central tools that you are using.

Signal object definitions

After having defined our baseline as well as loose (which is just baseline after overlap removal) objects, it is time to define the actual signal objects that are used as physics objects in the rest of the event selection. Signal objects are a subset of the baseline objects with tighter quality requirements and are used to define the search regions. SimpleAnalysis offers an easy way to select a subset of AnalysisObjects from another set of AnalysisObjects by using the AnalysisClass::filterObjects() method:

1
static AnalysisObjects filterObjects(const AnalysisObjects& cands, float ptCut, float etaCut=100., int id=0, unsigned int maxNum=10000);
This method returns a maxNum numbers of candidates AnalysisObjects& cands that satisfy a p_\mathrm{T}, \eta and identification requirement.

Electrons

Signal electrons have p_\mathrm{T} > 7 GeV, \eta < 2.47 and need to satisfy the TightLLH (ETightLH in SimpleAnalysis) identification working point as well as the FCLoose (EIsoFixedCutLoose in SimpleAnalysis) isolation working point. The transverse impact parameter d_0 needs to satisfy \vert d_0 / \sigma(d_0)\vert < 5 (ED0Sigma5). The longitudinal impact parameter relative to the PV is required to satsify z_0\sin\theta < 0.5 mm (EZ05mm).

1
auto signalElectrons = filterObjects(baselineElectrons, 7.0, 2.47, ETightLH | ED0Sigma5  | EZ05mm  | EIsoFixedCutLoose);

Muons

Signal muons have p_\mathrm{T} > 6 GeV, \eta < 2.70. Additionally, they need to fulfil \vert d_0 / \sigma(d_0)\vert < 3 (MuD0Sigma3) and z_0\sin\theta < 0.5 mm (MuZ05mm). Signal muons also need to satisfy the FCLoose isolation working point (MuIsoFixedCutLoose).

1
auto signalMuons = filterObjects(baselineMuons, 6.0, 2.7, MuD0Sigma3 | MuIsoFixedCutLoose);

Since we have now construced signal leptons, let's define two useful variables for later:

1
2
auto signalLeptons = signalElectrons + signalMuons;
auto baselineLeptons = baselineElectrons + baselineMuons;
Note how you can simply two lists of analysis objects by simple addition.

Jets

Signal jets need to fulfil p_\mathrm{T} > 30 GeV, \eta < 2.80. Jets with p_\mathrm{T} < 120 GeV need to be matched to the primary vertex by the Jet Vertex Tagger (JVT120Jet). Jets containing b-hadrons are identified using the MV2c10 tagger, using the 77% efficiency working point (BTag77MV2c10).

1
2
auto signalJets = filterObjects(baselineJets, 30.0, 2.80, JVT120Jet);
auto signalBJets = filterObjects(signalJets, 30.0, 2.8, BTag77MV2c10);

Instead of using small-R jets, you can also get large-R jets by using the reclusterJets() method:

1
auto largeRJets = reclusterJets(signalJets, 1.0, 30, 0.2, 0.05);
which uses the FastJet package for the clustering algorithm. Here, the parameters are the input jet collection, the radius of the reclustered large-R jets, the p_\mathrm{T} threshold as well as the two reclustering/trimming parameteres defining the radius-parameter of subjets as well as the minimum p_\mathrm{T} fraction.

Summary

We have now initialised our regions of interest, defined our baseline objects, performed the overlap removal procedure and finally defined the actual physics objects that are being used in the following. Summarising this, your MyAnalysisName::ProcessEvent() step should now 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
void MyAnalysisName::ProcessEvent(AnalysisEvent *event) {

  auto baselineElectrons = event->getElectrons(7.0, 2.47, ELooseBLLH);
  auto baselineMuons = event->getMuons(6.0, 2.70, MuMedium | MuNotCosmic | MuZ05mm | MuQoPSignificance);

  auto baselineJets = event->getJets(20.0, 4.5);

  auto met_Vect  = event->getMET(); //actually returns the 4-vector
  float met = met_Vect.Et();
  auto weights = event->getMCWeights();

  auto radiusCalcLep = [] (const AnalysisObject& lep,const AnalysisObject&) {
      return (0.04 + 10/lep.Pt()) > 0.4 ? 0.4 : 0.04 + 10/lep.Pt();
    };

  baselineElectrons = overlapRemoval(baselineElectrons, baselineMuons, 0.01);

  baselineJets = overlapRemoval(baselineJets, baselineElectrons, 0.2);
  baselineElectrons = overlapRemoval(baselineElectrons, baselineJets, radiusCalcLep);

  baselineJets = overlapRemoval(baselineJets, baselineMuons, 0.2, LessThan3Tracks);
  baselineMuons = overlapRemoval(baselineMuons, baselineJets, radiusCalcLep);

  auto signalElectrons = filterObjects(baselineElectrons, 7.0, 2.47, ETightLH | ED0Sigma5  | EZ05mm  | EIsoFixedCutLoose);
  auto signalMuons = filterObjects(baselineMuons, 6.0, 2.7, MuD0Sigma3 | MuIsoFixedCutLoose);

  auto signalLeptons = signalElectrons + signalMuons;
  auto baselineLeptons = baselineElectrons + baselineMuons;

  auto signalJets = filterObjects(baselineJets, 30.0, 2.80, JVT120Jet);
  auto signalBJets = filterObjects(signalJets, 30.0, 2.8, BTag77MV2c10);

  auto largeRJets = reclusterJets(signalJets, 1.0, 30, 0.2, 0.05);

  return;
}


Last update: February 16, 2021