7.9 B-Tagging

Complete: 5
Detailed Review status

Goals of this page

This page introduces the principal algorithms used to tag b-jets and will help you decide which one will be most suitable for your physics analysis. It then explains how to access the most commonly used b tag information.

Please also note that from version 1.7.0 onwards the bTagging data formats have been frozen. Whatever you read on 1.7.X you can safely applied to any 2.X.Y and 3.X.Y versions (at least for what concerns data formats).

Expert users may also wish to consult the b tag chapter of the CMSSW Offline Guide for more detailed information, such as how to modify and rerun the b tagging algorithms.

Contents

b Tag Algorithms - which one to use ?

Several b tag algorithms have been implemented in CMSSW. Some exploit the long B hadron lifetime, others its semi-leptonic decay mode and others use kinematic variables related to the high B hadron mass and hard b fragmentation function.

All begin by associating good quality tracks with jets. By default, b tag algorithms are only run on iterativeCone5CMS.CaloJets (i.e. uncalibrated, iterative-cone jets produced with Delta(R) = 0.5) up to CMSSW_2_X_Y, and on AntiKt05 starting from CMSSW_3_X_Y. Instructions on how to tag alternative jet collections are given below.

In CMSSW, all the algorithms produce a b tag "discriminator" for each jet, on which one can cut more or less tightly, in order to distinguish b-jets from non-b jets. Although the definition of the discriminator varies between algorithms, it is always true that b-jets have more positive values of the discriminator. This page will help you choose a sensible cut value for the discriminator.

  • "Track Counting" algorithm: More... Close This is a very simple tag, exploiting the long lifetime of B hadrons. It calculates the signed impact parameter significance of all good tracks, and orders them by decreasing significance. Its b tag discriminator is defined as the significance of the N'th track. It comes in two variations for N = 2 (high efficiency) or N = 3 (high purity).
  • "Jet Probability" algorithm: More... Close This is a more sophisticated algorithm, also exploiting the long lifetime of B hadrons. Its b tag discriminator is equal to the negative logarithm of the confidence level that all the tracks in the jet are consistent with originating from the primary vertex. This confidence level is calculated from the signed impact parameter significances of all good tracks. It reads the resolution function on these from a database (DB). Indeed, we have two versions of this tagger: JetProbabilityBJetTags and JetBProbabilityBJetTags - the latter uses only the four most displaced tracks, matching the typical reconstructed multiplicity of a B decay vertex.
  • "Soft Muon" and "Soft Electron" algorithms: More... Close These two algorithms tag b jets by searching for the lepton from a semi-leptonic B decay, which typically has a large Pt_rel with respect to the jet axis. Their b tag discriminators are the output of neural nets based on the leptons Pt_rel, impact parameter significance and a few other variables. For each of these taggers, we want to have a simple one (cutting basically only on the presence and pT of the lepton), and a complex one using also jet quantities to compute a MVA analysis. This is already the case for Muons; Electron taggers need refinement and are not really usable up to 2.1.X.
  • "Simple Secondary Vertex" algorithms: More... Close These class of algorithms reconstructs the B decay vertex using an adaptive vertex finder, and then uses variables related to it, such as decay length significance to calculate its b tag discriminator. It has been found to be more robust to Tracker misalignment than the other lifetime-based tags. CMSSW releases <= 35X contain one version of the algorithm (simpleSecondaryVertexBJetTags). Starting from 36X two versions will be provided: simpleSecondaryVertexHighEffBJetTags (equivalent to the previous version) and simpleSecondaryVertexHighPurBJetTags (with increased purity due to a cut on the track multiplicity at the secondary vertex)
  • "Combined Secondary Vertex" algorithm: More... Close This sophisticated and complex tag exploits all known variables, which can distinguish b from non-b jets. Its goal is to provide optimal b tag performance, by combining information about impact parameter significance, the secondary vertex and jet kinematics. (Currently lepton information is not included). . The variables are combined using a likelihood ratio technique to compute the b tag discriminator. A variant of this tagger combines the variables using the Multivariant Analysis (MVA) tool.

In CMSSW 1.5 - 1.6, the "Jet Probability" algorithm gives best performance and is consequently recommended for most physics analyses. In older versions of CMSSW, this is not available, and you should instead use the simpler, but fairly effective "Track Counting" algorithm. The "Soft Lepton" algorithms have several times lower efficiency due to the low B semi-leptonic branching ratio, but can sometimes be useful, particularly for early data taking or for measuring b tag performance with data. The "Combined Secondary Vertex" algorithm should ultimately give better performance than the Jet Probability algorithm, but only becomes competitive with it in CMSSW 1.7, and becomes the most powerful algorithm since CMSSW_2_X_Y (but beware it will be the most difficult to calibrate on data).

So the advice is

  • for initial analysis, stick to the TrackCounting;
  • when experience is gained on data, you can move to algorithms more complicated and performing.

b Tag AOD/RECO Data

The full data format is described here in the CMSSW Offline Guide. However, this section gives you an overview.

Each b tag algorithm produces a collection of JetTag and various kinds of TagInfo (Each collection has one entry for each jet, ordered as in the original jet collection). For most physics analyses, you will only need to access the JetTag object:

170_JetTag.png

This gives the b tag discriminator value and a references (link) to the jet. The JetTag collections are labelled by the b tag algorithm that produced them. In CMSSW 1.7, the most important of these labels are:

  • Track counting tag with N = 3: trackCountingHighPurBJetTags
  • Track counting tag with N = 2: trackCountingHighEffBJetTags
  • Jet probability tag: jetProbabilityBJetTags
  • Soft electron tag: softElectronBJetTags
  • Soft muon tag: softMuonBJetTags
  • Soft muon tag, not using muon impact parameter significance: softMuonNoIPBJetTags
  • Simple secondary vertex b tag: simpleSecondaryVertexBJetTags
  • Combined SV b tag using likelihood ratios: combinedSVBJetTags
  • Combined SV b tag using MVA: combinedSVMVABJetTags

For a full list of these labels in any CMSSW version, see the file RecoBTag/Configuration/python/RecoBTag_EventContent.cff .

The TagInfo contains the variables, such as track impact parameter significance, from which the b tag discriminator was calculated. It also gives you the jet and the tracks associated with it. There are several different types of TagInfo object, for the various b tag algorithms. The TrackIPTagInfo is shared by the Track Counting, Jet Probability and CombinedSV tags. Sophisticated b tags, such as the CombinedSV tag can use more than one TagInfo. (It also uses SecondaryVertexTagInfo).

170_TagInfo.png

The b tag algorithms also rely on the association of tracks with jets. This information is stored in JetTracksAssociation. However, it can be accessed by the JetTag or TagInfo classes, so you should not need to worry about this.

The above CMSSW 1.7 documentation is kept here and is valid for every later release (including CMSSW_2_X_Y and CMSSW_3_X_Y).

b Tag Analysis with bare ROOT, FWLite or EDAnalyzers

BEWARE - THIS DOES DOT WORK ANYMORE WITH CMSSW_3_X_Y; USE FWLITE INSTEAD

This section describes how to analyze b tag information with bare ROOT, FWLite and EDAnalyzers. The first two can be useful for a quick look at the data, whereas the latter is recommended for serious analysis. To follow the examples here, please find some t-tbar events with DBS. (For CMSSW 1.6, /store/mc/2007/10/20/RelVal-RelValTTbar-1192895175 /0000/00C41641-2A81-DC11-B6EA-0019DB29C620.root is available at CERN).

b Tag Analysis with bare ROOT

To quickly look at the b tag information inside a CMSSW dataset, it can be convenient to simply open the CMSSW dataset with bare ROOT. However, you will only be able to look at simple data members (e.g. float, int) inside classes, as illustrated below:

root rfio:///castor/cern.ch/cms//store/relval/CMSSW_2_1_4/RelValQCD_Pt_80_120/
GEN-SIM-DIGI-RAW-HLTDEBUG-RECO/IDEAL_V6_v1/0004/080739DB-816C-DD11-ADFC-001617E30F50.root

(note the need to add "rfio:///castor/cern.ch/cms/" to dataset name provided by DBS, in order to open it interatively).

TBrowser b;

and then navigate in the TBrowser to the Event directory, where you can find the b tags: recoJetdmRefProdTofloatAssociationVector_jetProbabilityBJetTags__RECO, for example. The name is of the form `C++ namespace' `class name'_label__RECO', as shown in the following screen capture:

Picture_9.png

Double-clicking on recoJetdmRefProdTofloatAssociationVector_jetProbabilityBJetTags__RECO and then on recoJetdmRefProdTofloatAssociationVector_jetProbabilityBJetTags__RECO.obj lets one see the C++ data members of the JetTag object, as shown in this screen capture:

Picture_12.png

Double-clicking on recoJetdmRefProdTofloatAssociationVector_jetProbabilityBJetTags__RECO.obj.data shows the histogram of the b `probability' tag discriminator value:

Picture_11.png

The above CMSSW 1.7 documentation is kept here.

b Tag Analysis with CMSSW Framework Light (fwLite)

This is convenient for quickly looking at the data, with the bonus of being able to use the public member functions of any class stored in the event. You should first set up your CMSSW environment from inside your project area and start ROOT (we work in this example with CMSSW_3_3_5 and a RelVal TTbar file):

eval `scramv1 runtime -csh`
root

then make the CMSSW FWlite libraries available and open the dataset:

gSystem.Load("libFWCoreFWLite.so");
AutoLibraryLoader::enable();   

TFile * f= TFile::Open("/castor/cern.ch/cms/store/relval/CMSSW_3_3_5/RelValTTbar/GEN-SIM-RECO/
MC_31X_V9-v1/0008/ACFD3497-68DC-DE11-A974-00261894388A.root")

We can now set some aliases for B tag objects. (If you don't know the exact name of the B tag objects, check using TBrowser). While defining aliases is not necessary, it does make the subsequent analysis much simpler:

Events->SetAlias("bmu", 
"recoJetedmRefToBaseProdTofloatsAssociationVector_softMuonBJetTags__RECO.obj");

Functionality is somehow limited by the fact that an association map is involved, you can for example get the number of btagged jets as

Events->Draw("bmu.size()")

Screen_shot_2009-12-11_at_11.15.30.png

The above CMSSW 1.7 and beyond documentation is kept here.

b Tag Analysis with a CMSSW EDAnalyzer

Most physics analyses will use this method. You should start by creating your own generic EDAnalyzer and .cfg files, as described here earlier in the Workbook.

To look at b tag information in each event, you can add to your EDAnalyzer::analyze() code like this. :

#include "DataFormats/BTauReco/interface/JetTag.h"

# Get b tag information
edm::Handle<reco::JetTagCollection> bTagHandle;
iEvent.getByLabel("trackCountingHighEffBJetTags", bTagHandle);
const reco::JetTagCollection & bTags = *(bTagHandle.product());

# Loop over jets and study b tag info.
for (int i = 0; i != bTags.size(); ++i) {
   cout<<" Jet "<< i 
         <<" has b tag discriminator = "<<bTags[i].second
         << " and jet Pt = "<<bTags[i].first->pt()<<endl;
}

If you have a Jet, rather than a JetTag, and wish to know if it is b-tagged, there are several ways of doing so. One which always works is to perform angular matching between the Jet and the JetTag::jet(). (The match should be perfect if your JetCollection was used to produce the JetTagCollection).

Also add

 <use name=DataFormats/BTauReco> 
to the CMS.BuildFile.

The above CMSSW 1.7 and beyond documentation is kept here.

How to determine the True Flavour of a Jet.

Since CMSSW_2_X_Y the old way (see 1.6 chapter) of finding the flavour of a jet is deprecated, and we use GenParticles, which are available in any data tier up to (including) AODSIM.

To access the information, you have before to add in the cmsRun Path the sequence process.myPartons * process.AK5Flavour which can be defined via

process.load("CMS.PhysicsTools.JetMCAlgos.CaloJetsMCFlavour_cfi")  

Here is how to operate:

  • take a CMSSW version (let's say CMSSW_3_3_5)
  • prepare a cmsRun sequence which includes
process.load("CMS.PhysicsTools.JetMCAlgos.CaloJetsMCFlavour_cfi")  
and runs before your analyzer
process.myPartons * process.AK5Flavour

In the .cc code of the analyzer includes lines like

  edm::Handle<reco::JetFlavourMatchingCollection> jetMC;
  iEvent.getByLabel(jetMCSrc, jetMC);
to load the product. jetMCsrc is the collection you want to analyze, depending on the jet type. Default should be 'AK5byValAlgo'.

This product is a mapping between Reco::Jets and an int which is the flavour (please note it can be positive and negative, so to select b quark ask for std::abs(flavour)==5. Then, inside the loop in which you analyze jets, you can ask for the flavour. The easiest way is to produce event by event a searchable map, with lines like

typedef std::map<edm::RefToBase<reco::Jet>, unsigned int, JetRefCompare> FlavourMap;
FlavourMap flavours;
  for (reco::JetFlavourMatchingCollection::const_iterator iter = jetMC->begin();
       iter != jetMC->end(); iter++) {
    unsigned int fl = std::abs(iter->second.getFlavour());
    flavours.insert(FlavourMap::value_type(iter->first, fl));
  }

Then, when you have a RefToBase<reco::Jet> (for example from a view), you can ask for the flavour:

unsigned int myFlavour=0;
RefToBase<reco::Jet> aJet; // fill it from the collection you want to probe!
if (flavours.find (aJet) == flavours.end()) {
    std::cout <<" Cannot access flavour for this jet - not in the Map"<<std::endl;
} else {
   myFlavour = flavours[aJet];
}

Please note the value of myFlavour can be zero if there was not a positive identification of the jet flavour (for example, no partons inside the matching cone).

In summary the following things need to be changed in order to get the jet flavour.

  • Make changes in your MyEDAnalyzer.cc and it should look like this MyEDAnalyzer.cc

  • Make changes in the BuildFile and it should look like this CMS.BuildFile

Then do scram b and run it by doing cmsRun myedanalyzer_cfg.py

You should see something like this (click on Show result below):

Begin processing the 1st record. Run 1, Event 8753440, LumiSection 5674 at 18-May-2010 10:31:54 CEST
 Jet 0 has b tag discriminator = -100 and jet Pt = 1.51255
 Jet 1 has b tag discriminator = -100 and jet Pt = 1.12788
 Jet 2 has b tag discriminator = -100 and jet Pt = 1.00847
 Jet 0 has flavour = 21
 Jet 1 has flavour = 0
 Jet 2 has flavour = 0
Begin processing the 2nd record. Run 1, Event 8753441, LumiSection 5674 at 18-May-2010 10:31:54 CEST
 Jet 0 has b tag discriminator = -100 and jet Pt = 3.76868
 Jet 1 has b tag discriminator = -100 and jet Pt = 3.47249
 Jet 2 has b tag discriminator = -100 and jet Pt = 3.33955
 Jet 3 has b tag discriminator = -100 and jet Pt = 2.79301
 Jet 4 has b tag discriminator = -100 and jet Pt = 1.48905
 Jet 5 has b tag discriminator = -100 and jet Pt = 1.48498
 Jet 6 has b tag discriminator = -100 and jet Pt = 1.21182
 Jet 7 has b tag discriminator = -100 and jet Pt = 1.19971
 Jet 8 has b tag discriminator = -100 and jet Pt = 1.16671
 Jet 9 has b tag discriminator = -100 and jet Pt = 1.12592
 Jet 10 has b tag discriminator = -100 and jet Pt = 1.03013
 Jet 0 has flavour = 2
 Jet 1 has flavour = 21
 Jet 2 has flavour = 21
 Jet 3 has flavour = 0
 Jet 4 has flavour = 0
 Jet 5 has flavour = 0
 Jet 6 has flavour = 0
 Jet 7 has flavour = 0
 Jet 8 has flavour = 0
 Jet 9 has flavour = 0
 Jet 10 has flavour = 0
Begin processing the 3rd record. Run 1, Event 8753442, LumiSection 5674 at 18-May-2010 10:31:54 CEST
 Jet 0 has b tag discriminator = -2.35632 and jet Pt = 2.80738
 Jet 1 has b tag discriminator = -100 and jet Pt = 2.4376
 Jet 2 has b tag discriminator = -100 and jet Pt = 2.3545
 Jet 3 has b tag discriminator = -100 and jet Pt = 2.02223
 Jet 4 has b tag discriminator = -100 and jet Pt = 1.65683
 Jet 5 has b tag discriminator = -100 and jet Pt = 1.6328
 Jet 6 has b tag discriminator = -100 and jet Pt = 1.61783
 Jet 7 has b tag discriminator = -100 and jet Pt = 1.57096
 Jet 8 has b tag discriminator = -100 and jet Pt = 1.56623
 Jet 9 has b tag discriminator = -100 and jet Pt = 1.38087
 Jet 10 has b tag discriminator = -100 and jet Pt = 1.33699
 Jet 11 has b tag discriminator = -100 and jet Pt = 1.27284
 Jet 12 has b tag discriminator = -100 and jet Pt = 1.09305
 Jet 0 has flavour = 0
 Jet 1 has flavour = 0
 Jet 2 has flavour = 0
 Jet 3 has flavour = 21
 Jet 4 has flavour = 21
 Jet 5 has flavour = 0
 Jet 6 has flavour = 21
 Jet 7 has flavour = 0
 Jet 8 has flavour = 21
 Jet 9 has flavour = 21
 Jet 10 has flavour = 0
 Jet 11 has flavour = 0
 Jet 12 has flavour = 0

The above CMSSW 2_X_Y and beyond documentation is kept here

How to make Plots of b Tag Performance

Warning: Can't find topic CMSPublic.WorkBookBTagPerformance2XY

The above CMSSW 2.X and beyond documentation is kept here.

#Advancedtopics

Advanced topics

This section will cover additional topics related to b-tagging in CMSSW, in particular how to run a b-tag algorithm that was not included in the original data sample, or how to re-run some or all taggers with a different set of jets and/or tracks.

How to run a new tagger

Obtaining additional b-tag discriminators from the data stored in AOD is rather simple, as it only requires running the corresponding module. As an example, it will be shown how to run the softMuonByIP3d soft lepton tagger.

Add to your python configuration file the following snippet to take care of the general configuration (change the GlobalTag to match the dataset you are using):

# default configuration with frontier conditions
process.load("Configuration.StandardSequences.FrontierConditions_CMS.GlobalTag_cff")
process.load("Configuration.StandardSequences.MagneticField_cff")
process.load("Configuration.StandardSequences.Geometry_cff")
process.GlobalTag.globaltag = 'IDEAL_V9::All'

# b-tagging general configuration
process.load("RecoBTag.Configuration.RecoBTag_cff")

Then add the configuration for the extra tagger you want to run - n this example, softMuonByIP3d. Usually both an ESProducer and EDProducer are needed, the first taking care of the algorithm configuration, the second of the b-tag raw data to use:

# configure the softMuonByIP3d ESProducer and EDProducer
process.load("RecoBTag.SoftLepton.softLeptonByIP3dES_cfi")
process.load("RecoBTag.SoftLepton.softMuonByIP3dBJetTags_cfi")

Finally, add the producer to a path, and run it:

# run the b-tag sequence
process.btag  = cms.Path( process.softMuonByIP3dBJetTags )

Of course, remember to add the new module(s) before any analysys code that will use them.

How to run on a different jet collection

Running the b-tag algorithms on a different jet collection requires access to the full RECO data.

Note: an alternative is to use the b-tag information precomputed with the default jet collection, matching the jets from the new collection to those from the original one.

In order to rerun all the algorithms on a different jet (or track) collection, one needs to duplicate the original b-tag sequence. Only the EDProducers needs to be cloned.

import FWCore.ParameterSet.Config as cms

# default configuration with frontier conditions
process.load("Configuration.StandardSequences.MagneticField_cff")
process.load("Configuration.StandardSequences.Geometry_cff")
process.load("Configuration.StandardSequences.Reconstruction_cff")
process.load("Configuration.StandardSequences.FrontierConditions_CMS.GlobalTag_cff")
process.GlobalTag.globaltag = 'IDEAL_V9::All'

# b-tagging general configuration
process.load("RecoJets.JetAssociationProducers.ic5JetTracksAssociatorAtVertex_cfi")
process.load("RecoBTag.Configuration.RecoBTag_cff")

The first step is to clone the definition of the association between jets and tracks. This is where most of the customization will take place. As an example, here generalTracks and sisCone5CMS.CaloJets are used - to try different ones, change them to what you need throughout the configuration (check later the soft lepton part, too):

# create a new jets and tracks association
process.newJetTracksAssociatorAtVertex = process.ic5JetTracksAssociatorAtVertex.clone()
process.newJetTracksAssociatorAtVertex.jets = "sisCone5CMS.CaloJets"
process.newJetTracksAssociatorAtVertex.tracks = "generalTracks"

The one needs to clone the b-tag producers and instruct them to use this new collection. First the impact parameter-based b-tag:

# impact parameter b-tag
process.newImpactParameterTagInfos = process.impactParameterTagInfos.clone()
process.newImpactParameterTagInfos.jetTracks = "newJetTracksAssociatorAtVertex"
process.newTrackCountingHighEffBJetTags = process.trackCountingHighEffBJetTags.clone()
process.newTrackCountingHighEffBJetTags.tagInfos = cms.VInputTag( cms.InputTag("newImpactParameterTagInfos") )
process.newTrackCountingHighPurBJetTags = process.trackCountingHighPurBJetTags.clone()
process.newTrackCountingHighPurBJetTags.tagInfos = cms.VInputTag( cms.InputTag("newImpactParameterTagInfos") )
process.newJetProbabilityBJetTags = process.jetProbabilityBJetTags.clone()
process.newJetProbabilityBJetTags.tagInfos = cms.VInputTag( cms.InputTag("newImpactParameterTagInfos") )
process.newJetBProbabilityBJetTags = process.jetBProbabilityBJetTags.clone()
process.newJetBProbabilityBJetTags.tagInfos = cms.VInputTag( cms.InputTag("newImpactParameterTagInfos") )

Then the secondary vertex-based b-tag. Note that these producers inherit the jets and tracks they use from the impact parameter modules: # secondary vertex b-tag
process.newSecondaryVertexTagInfos = process.secondaryVertexTagInfos.clone()
process.newSecondaryVertexTagInfos.trackIPTagInfos = "newImpactParameterTagInfos"
process.newSimpleSecondaryVertexBJetTags = process.simpleSecondaryVertexBJetTags.clone()
process.newSimpleSecondaryVertexBJetTags.tagInfos = cms.VInputTag( cms.InputTag("newSecondaryVertexTagInfos") )
process.newCombinedSecondaryVertexBJetTags = process.combinedSecondaryVertexBJetTags.clone()
process.newCombinedSecondaryVertexBJetTags.tagInfos = cms.VInputTag( cms.InputTag("newImpactParameterTagInfos"), cms.InputTag("newSecondaryVertexTagInfos") )
process.newCombinedSecondaryVertexMVABJetTags = process.combinedSecondaryVertexMVABJetTags.clone()
process.newCombinedSecondaryVertexMVABJetTags.tagInfos = cms.VInputTag( cms.InputTag("newImpactParameterTagInfos"), cms.InputTag("newSecondaryVertexTagInfos") )

And the soft lepton b-tag. These producers will accept as input either the raw jets, or the association collection:

# soft electron b-tag
process.newSoftElectronTagInfos = process.softElectronTagInfos.clone()
process.newSoftElectronTagInfos.jets = "sisCone5CMS.CaloJets"
process.newSoftElectronBJetTags = process.softElectronBJetTags.clone()
process.newSoftElectronBJetTags.tagInfos = cms.VInputTag( cms.InputTag("newSoftElectronTagInfos") )

# soft muon b-tag
process.newSoftMuonTagInfos = process.softMuonTagInfos.clone()
process.newSoftMuonTagInfos.jets = "sisCone5CMS.CaloJets"
process.newSoftMuonBJetTags = process.softMuonBJetTags.clone()
process.newSoftMuonBJetTags.tagInfos = cms.VInputTag( cms.InputTag("newSoftMuonTagInfos") )
process.newSoftMuonByIP3dBJetTags = process.softMuonByIP3dBJetTags.clone()
process.newSoftMuonByIP3dBJetTags.tagInfos = cms.VInputTag( cms.InputTag("newSoftMuonTagInfos") )
process.newSoftMuonByPtBJetTags = process.softMuonByPtBJetTags.clone()
process.newSoftMuonByPtBJetTags.tagInfos = cms.VInputTag( cms.InputTag("newSoftMuonTagInfos") )

Finally, there needs to be a new path running all these modules

# prepare a path running the new modules
process.newJetTracksAssociator = cms.Sequence(
    process.newJetTracksAssociatorAtVertex
)

process.newJetBtaggingIP = cms.Sequence(
    process.newImpactParameterTagInfos * (
        process.newTrackCountingHighEffBJetTags +
        process.newTrackCountingHighPurBJetTags +
        process.newJetProbabilityBJetTags +
        process.newJetBProbabilityBJetTags
    )
)

process.newJetBtaggingSV = cms.Sequence(
    process.newImpactParameterTagInfos *
    process.newSecondaryVertexTagInfos * (
        process.newSimpleSecondaryVertexBJetTags +
        process.newCombinedSecondaryVertexBJetTags +
        process.newCombinedSecondaryVertexMVABJetTags
    )
)

process.newJetBtaggingEle = cms.Sequence(
    process.btagSoftElectrons *
    process.newSoftElectronTagInfos *
    process.newSoftElectronBJetTags
)

process.newJetBtaggingMu = cms.Sequence(
    process.newSoftMuonTagInfos * (
        process.newSoftMuonBJetTags +
        process.newSoftMuonByIP3dBJetTags +
        process.newSoftMuonByPtBJetTags
    )
)

process.newJetBtagging = cms.Sequence(
    process.newJetBtaggingIP +
    process.newJetBtaggingSV +
    process.newJetBtaggingEle +
    process.newJetBtaggingMu
)

process.newBtaggingPath = cms.Path(
    process.newJetTracksAssociator *
    process.newJetBtagging
)

Further Information

Review status

Reviewer/Editor and Date (copy from screen) Comments
JyothsnaK - 18-May-2010 Fixed the truth matching section to get the jet flavour
JyothsnaK - 17-Feb-2010 Fixed a typo in the EDAnalyzer section and also added in this section changes to be made to the BuildFile.
AndreaBocci - 29 Jan 2009 Added advanced topics
Main.Aresh - 27 Feb 2008 Changes in verbatim elements (in https://twiki.cern.ch/twiki/bin/edit/CMS/WorkBookBTagEdAnalyzer160, https://twiki.cern.ch/twiki/bin/edit/CMS/WorkBookBTagFWlite160) because of some lines too long for printable version
IanTomalin - 18 Oct 2007 Created from original WorkBookTauTagging page, with updates to reflect CMSSW 1.6, and simplifications to correspond more to what a new user would wish to know. More sophisticated documentation will be moved to Offline Guide.

Responsible: AlexanderSchmidt & LucaScodellaro

Last reviewed by: JyothsnaK - 18-May-2010

Topic attachments
I Attachment History Action Size Date Who Comment
PowerPointppt 160_btag_pictures.ppt r1 manage 303.5 K 2007-10-30 - 18:44 IanTomalin Powerpoint file used to produce all the pictures here
PowerPointppt 170_btag_pictures.ppt r2 r1 manage 23.5 K 2007-10-31 - 19:20 IanTomalin Powerpoint file used to produce all the pictures here
Edit | Attach | Watch | Print version | History: r33 < r32 < r31 < r30 < r29 | Backlinks | Raw View | WYSIWYG | More topic actions
Topic revision: r33 - 2014-01-24 - DanielBloch


ESSENTIALS

ADVANCED TOPICS


 
    • Cern Search Icon Cern Search
    • TWiki Search Icon TWiki Search
    • Google Search Icon Google Search

    CMSPublic 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