Jet Analysis in CMSSW 3_3_x

Complete: 3
Detailed Review status

Goals of this page

The purpose of this tutorial is to give an idea about what Jet objects are in CMSSW and how to work with them. CMSSW is in continuous development so the syntax for particular commands may change. To resolve this issue separate tutorial pages for different CMSSW releases are provided.
Other versions of this tutorial may be found on the primary Jet Reconstruction tutorial page WorkBookJetAnalysis. This is an updated version of the tutorial for previous CMSSW versions by Fedor Ratnikov and Robert Harris

Author: Hartmut Stadie

Contents

Set up your Environment

Set your runtime environment (shown for release 3_3_1):
cd CMSSW_3_3_1/src
cmsenv

No files need to be copied at this step - necessary files will be refered later.

Jet Reconstruction Basics

A few steps are necessary to walk from calorimeter digis (a.k.a. frames), which are essentially non-linear ADC counts, to Jet and MET objects:

Dataflow for Jet Reconstruction

  1. convert ADC counts to energy in one calorimeter cell
    • HcalSimpleReconstructor does this calculation for HCAL
    • ECAL does it in two steps
      • first produces uncalibrated hits: EcalWeightUncalibratedRecHitProducer
      • then calibrates hits: EcalAnalFitUncalibratedRecHitProducer
  2. combine ECAL and HCAL cells into projective towers corresponding to HCAL granularity
    • CaloTowersCreator does this operation
  3. convert CMS.CaloTowers into standard objects CaloTowerCandidateCreator (see more details in the next section)
  4. run clustering algorithm to produce Jets
    • Several basic algorithms are currently implemented for CMS
      • Midpoint Cone algorithm
      • Iterative Cone algorithm
      • KT algorithm
      • Seedless Infrared Safe cone algorithm
      • Anti-KT algorithm

Different Jet Types

Jet clusterization algorithms are very generic and may run on any set of 4-momenta. This feature is used to separate the jet algorithm itself from the physical nature of the jet constituents. Before running jet algorithm on the set of objects, these objects may need to be converted into standard format a.k.a. Candidates, if not inherited from Candidate already.

Four jet flavors are supported:

  • CMS.CaloJets - jets produced from CMS.CaloTowers as described above.
  • GenJets - jets produced from generator level MC particles.
  • PFJets - jets produced from Particle Flow candidates.
  • BasicJets - jets produced from arbitrary collection of Candidate inherited objects.

In this tutorial we use two kinds of input objects for Jet analysis: CMS.CaloJets and GenJets

Jets of any flavor always contain generic jet information. In addition some extra information specific to particular jet flavors is also available.

Generic Jet Properties

Every jet contains information about its 4-momentum. This is the quantity to use for generic kinematics analysis.

Specific Jet Properties

  • CaloJet also contains information about the fraction of energy deposited in a particular calorimeter compartment or group of compartments: ECAL, HCAL, HB(barrel), HO(outer), HE(endcap), HF(forward).

Note that for CMS.CaloJets, the functions p4().eta() and eta() both give you an estimate of the jet rapidity determined with the calorimeters alone. One can get better rapidity resolution (particularly for high energy jets) by instead taking the rapidity using functions physicsEtaQuick(float zPV) or physicsP4(float zPV). These correct the jet direction for the z-offset of the reconstructed primary vertex zPV from zero.

  • GenJet contains information about fraction from different types of generated particles, i.e. from charged hadrons.

  • PFJet conains information about energy contribution and number of contribuing different types of Particle Flow objects.

Jet Algorithms

Run Jet Reconstruction

The jet reconstruction procedure is run as part of the default reconstruction. In order to see final jets, one can simply run a full simulation and reconstruction step in order to see how CMSSW works. We simply use cmsDriver in order to:
  • Generating single pion events with fixed pT of about 50 GeV/c
  • GEANT simulation of detector response
  • Digitization - simulation of detector readout
  • Local calorimetry reconstruction
  • Combining hits into Calorimetry Towers
  • Running CMS.CaloJets reconstruction
  • Running GenJets reconstruction

In order to obtain such a python config file we run:

cmsDriver.py SinglePiE50HCAL_cfi -s GEN,SIM,DIGI,L1,DIGI2RAW,RAW2DIGI,RECO --conditions=FrontierConditions_CMS.GlobalTag,MC_31X_V9::All -n 50 --no_exec
This should produce a config file called:
SinglePiE50HCAL_cfi_GEN_SIM_DIGI_L1_DIGI2RAW_RAW2DIGI_RECO_MC.py
To run the full chain we simply run
cmsRun SinglePiE50HCAL_cfi_GEN_SIM_DIGI_L1_DIGI2RAW_RAW2DIGI_RECO_MC.py
Note: you may want change settings like, the number of events, output filename, content, etc. in the config file. In order to learn about this, refer to the documentation of the python configs.

Jet Analysis

#CMS.DataSamples

Data sample

In following examples we will explore two data samples. One file is SinglePiE50HCAL_cfi_GEN_SIM_DIGI_L1_DIGI2RAW_RAW2DIGI_RECO.root produced on the previous step and containing 50 single pion events.

Bare ROOT

> root -l SinglePiE50HCAL_cfi_GEN_SIM_DIGI_L1_DIGI2RAW_RAW2DIGI_RECO.root
root [0] TBrowser b;

Travel to list of objects in the file: "ROOT Files" -> "SinglePiE50HCAL_cfi_GEN_SIM_DIGI_L1_DIGI2RAW_RAW2DIGI_RECO.root" -> "Events" Now you can browse EDM objects available in the file: Objects list,

for example CMS.CaloJets: CaloJet.

Note: clicking on the leaf is equivalent to executing "Draw" command from the ROOT prompt:

   root [2] Events->Draw("recoCaloJets_ak5CaloJets__RECO.obj.pt_");
 
produces transverse momentum distribution for reconstructed jets

FWLite Analysis

Using TTree::Draw

FWLite mode allows one to address not only members of objects, but also object's methods. See SWGuideEDMWithRoot for details. The precondition for using FWLite is to have regular CMSSW environment established, or at least CMSSW objects libraries available. Fortunatelly this is not a big restriction. The following two commands executed before openning data file enable autoloading of objects libraries:
   > root
   root [0] gSystem->Load("libFWCoreFWLite");
   root [1] AutoLibraryLoader::enable();
   root [2] TFile f ("SinglePiE50HCAL_cfi_GEN_SIM_DIGI_L1_DIGI2RAW_RAW2DIGI_RECO.root");
   root [3] TBrowser b;

Now Jet objects (like all other CMSSW objects) in the browser are accompanied by available functions, which are clickable and may be used in ROOT scripts.

  root [3] Events->Draw("recoCaloJets_ak5CaloJets__RECO.obj.energy()");
Corresponding methods are indicated by the exclamation leaf in the browser window.

Using FWLite::Event

Please see the FWLite part of the workbook for a nice example.

Access from Framework

The following example demonstrates use of an EDAnalyzer for CMS.CaloJets, GenJets and PFJets. It analyzes the simulated QCD data file: '/store/relval/CMSSW_3_3_0/RelValQCD_Pt_80_120/GEN-SIM-DIGI-RECO/STARTUP31X_V8_FastSim_Early10TeVCollision-v1/0001/0E3E1B64-56B7-DE11-A8E1-00304867924E.root'

Instructions to get and compile the code:

cmsrel CMSSW_3_3_4
cd CMSSW_3_3_4/src/
cmsenv
cvs co -r V00-04-09 RecoJets/JetAnalyzers
scram b
Then run the example:
cd RecoJets/JetAnalyzers/test
cmsRun JetPlotsExample.py
The analysis is done by the template class JetPlotsExample.cc which can handle CMS.CaloJets, GenJets or PFJets. It loops over a configurable number of leading jets in the event, makes some basic histograms, and writes them to a file. The jet reconstruction algorithm read from the event is specified by the label JetAlgorithm in the configuration file JetPlotsExample.py.

In this example, the iterativeCone5 algorithm is chosen, but other algorithms should work as well, provided jets from those algorithms are reconstructed in the event. Browsing the input data file in ROOT is a straightforward way to see what jet algorithms are currently available and what their label is (the second field Y in the branch name X_Y_Z). If the label for the algorithm changes in future releases then the user will need to change the label in the config file. As mentioned above, the user will need to change the input file names in the configuration file to point at whatever data file needs to be read.

The particular example consists of three identical modules, each one of which is a specific instance of the template class, in order to handle CMS.CaloJets, GenJets and PFJets. Each module requires some configurable parameters: the jet algorithm (JetAlgorithm), the name of the output ROOT file (HistoFileName) and the number of leading jets in the event that will be analyzed (NJets). Apparently, the user can modify the path to ask for only one or two modules (e.g. if the user wants to analyze a real data file only CMS.CaloJets and PFJets will be available).

Jet Energy Corrections

Contribution of Konstantinos Kousouris

Introduction

You can find more information related to Jet Energy Corrections here.

The latest recommended corrections and global tags for data and MC simulation can be found at the link https://twiki.cern.ch/twiki/bin/view/CMS/JECDataMC

Instructions

How to connect to a global tag

process.load('Configuration.StandardSequences.Services_cff')
process.load('Configuration.StandardSequences.FrontierConditions_GlobalTag_cff')
process.GlobalTag.globaltag = 'GR_R_52_V9::All'

How to connect to a local SQLite file

The recommended way of accessing JEC is via a global tag, as described above. However, it is also possible to read them from an SQLite file, which is mainly exercised during the initial testing phase in preparation of a global tag.

In order to read conditions from an SQLite file, one should use module PoolDBESSource as follows:

from CondCore.DBCommon.CondDBSetup_cfi import CondDBSetup
process.jec = cms.ESSource('PoolDBESSource',
    CondDBSetup,
    connect = cms.string('sqlite:Fall15_V2_DATA.db'),
    toGet = cms.VPSet(
        cms.PSet(
            record = cms.string('JetCorrectionsRecord'),
            tag    = cms.string('JetCorrectorParametersCollection_Fall15_V2_DATA_AK4PFchs'),
            label  = cms.untracked.string('AK4PFchs')
        ),
        cms.PSet(
            record = cms.string('JetCorrectionsRecord'),
            tag    = cms.string('JetCorrectorParametersCollection_Fall15_V2_DATA_AK8PFPuppi'),
            label  = cms.untracked.string('AK8PFPuppi')
        ),
        # ...and so on for all jet types you need
    )
)

# Add an ESPrefer to override JEC that might be available from the global tag
process.es_prefer_jec = cms.ESPrefer('PoolDBESSource', 'jec')

You will need to adapt the above snippet to your specific case, providing the path to the SQLite file and listing all required jet types with appropriate tags. In order to resolve the location of the file with the FileInPath functionality, use prefix sqlite_fip:. The tags are specific to a particular SQLite file. Their list can be printed with the command conddb --db <dbfile.db> listTags. Although it is technically possible to provide an arbitrary string for the label attribute, it is recommended to stick to the jet type labels used in global tags:

  • AK4PF
  • AK4PFchs
  • AK4PFPuppi
  • AK8PF
  • AK8PFchs
  • AK8PFPuppi

How to get the JEC txt files from the database

For those who need to produce the corresponding JEC txt files (for standalone use for example), one can download these from the wiki JECDataMC. This is a preferred method. One can also reproduce the original txt files using the JetCorrectorDBReader module. Below is a complete cfg file. However, this method is error-prone: for GTs which are interleaved by IOVs, it produces text files corresponding to the first IOV encountered, which may not be what the user needs.

  import FWCore.ParameterSet.Config as cms
  process = cms.Process("jectxt")
  process.load('Configuration.StandardSequences.Services_cff')
  process.load('Configuration.StandardSequences.FrontierConditions_GlobalTag_condDBv2_cff')
  # define your favorite global tag
  process.GlobalTag.globaltag = cms.string('74X_mcRun2_asymptotic_v4')
  process.maxEvents = cms.untracked.PSet(input = cms.untracked.int32(1))
  process.source = cms.Source("EmptySource")
  process.readAK4PFchs    = cms.EDAnalyzer('JetCorrectorDBReader',  
        # below is the communication to the database 
        payloadName    = cms.untracked.string('AK4PFchs'),
        # this is used ONLY for the name of the printed txt files. You can use any name that you like, 
        # but it is recommended to use the GT name that you retrieved the files from.
        globalTag      = cms.untracked.string('74X_mcRun2_asymptotic_v4'),  
        printScreen    = cms.untracked.bool(False),
        createTextFile = cms.untracked.bool(True)
  )
  process.readAK8PFchs = process.readAK4PFchs.clone(payloadName = 'AK8PFchs')
  process.readAK4PF = process.readAK4PFchs.clone(payloadName = 'AK4PF')
  process.p = cms.Path(process.readAK4PFchs * process.readAK8PFchs * process.readAK4PF)

If the config file above is run, a number of txt files will be created in your local area. These can be used as inputs to your standalone code. The names of the txt files will have the format 74X_mcRun2_asymptotic_v4_XXXX_YYYY.txt, where XXXX is the correction level (e.g. L2Relative), and where YYYY is the payload name (jet algo).

Note: JEC uncertainty sources are not accessible from data base, but are provided as txt files at https://twiki.cern.ch/twiki/bin/view/CMS/JECUncertaintySources.

Applying the Jet Energy Corrections

UPDATED: Since CMSSW 7.6.X, the consumes interface is now mandatory. The jet corrector framework was not modular enough to enforce this, so a new one was developed. The old interface is deprecated (using JetCorrector::getJetCorrector and friends) and won't work anymore on CMSSW 7.6.X and above.

Sections below have been updated to reflect these changes, but you can find a more in-depth migration tutorial on this twiki page: https://twiki.cern.ch/twiki/bin/view/CMSPublic/FWMultithreadedConsumesJetCorrectorMigration

Jet Energy Corrections in FWLite

Simple jet corrections, which are dependent from jet properties only, may be directly accessed from ROOT scripts running in FWLite environment. The class that delivers the jet correction is called FactorizedJetCorrector and is located in the CondFormats/JetMETObjects package.

Step by step instructions for use in a ROOT macro.

  • Step1 (Load the FWLite libraries)
   gSystem->Load("libFWCoreFWLite.so");
   AutoLibraryLoader::enable(); 

  • Step2 (Construct a vector with the JetCorrectorParameters objects.)
   // Create the JetCorrectorParameter objects, the order does not matter.
   // YYYY is the first part of the txt files: usually the global tag from which they are retrieved
   JetCorrectorParameters *ResJetPar = new JetCorrectorParameters("YYYY_L2L3Residual_AK5PF.txt"); 
   JetCorrectorParameters *L3JetPar  = new JetCorrectorParameters("YYYY_L3Absolute_AK5PF.txt");
   JetCorrectorParameters *L2JetPar  = new JetCorrectorParameters("YYYY_L2Relative_AK5PF.txt");
   JetCorrectorParameters *L1JetPar  = new JetCorrectorParameters("YYYY_L1FastJet_AK5PF.txt");
   //  Load the JetCorrectorParameter objects into a vector, IMPORTANT: THE ORDER MATTERS HERE !!!! 
   vector<JetCorrectorParameters> vPar;
   vPar.push_back(*L1JetPar);
   vPar.push_back(*L2JetPar);
   vPar.push_back(*L3JetPar);
   vPar.push_back(*ResJetPar);
  • Step3 (Construct a FactorizedJetCorrector object)
   FactorizedJetCorrector *JetCorrector = new FactorizedJetCorrector(vPar)
  • Step4 (Set the necessary values, for every jet, inside the jet loop)
   JetCorrector->setJetEta(eta);
   JetCorrector->setJetPt(pt);
   JetCorrector->setJetA(area);
   JetCorrector->setRho(rho); 
  • Step5 (Get the correction factor or a vector of the individual correction factors)
   double correction = JetCorrector->getCorrection();
OR
   vector<double> factors = JetCorrector->getSubCorrections();

IMPORTANT: the getSubCorrections member function returns the vector of the subcorrections UP to the given level. For example in the example above, factors[0] is the L1 correction and factors[3] is the L1+L2+L3+Residual correction.

Instructions for use in CRAB.

  • Step1 (Add the relative file path to the EDAnalyzer's configuration parameters in python)
   // YYYY is the first part of the txt files: usually the global tag from which they are retrieved
   L1File = cms.string("path/to/YYYY_L1FastJet_AK4PFchs.txt")
  • Step2 (Construct a vector with the JetCorrectorParameters objects.)
   // Create the JetCorrectorParameter objects, the order does not matter.
   JetCorrectorParameters *L1Parms;
   std::string l1file;
   l1file = edm::FileInPath(iConfig.getParameter<std::string> ("L1File")).fullPath();
   L1Parms = new JetCorrectorParameters(l1file);

   //  Load the JetCorrectorParameter objects into a vector, IMPORTANT: THE ORDER MATTERS HERE !!!! 
   vector<JetCorrectorParameters> vPar;
   vPar.push_back(*L1Parms);
  • Step3 (Continue using the ROOT macro procedure)

Jet Correction Service

The foundamental unit for the application of jet corrections in CMSSW is the Jet Correction Service. There is one service per level of correction plus a chain service that allows to apply several corrections sequentially. The typical structure of a service is the following (example with L2 Relative correction):

ak5CaloL2Relative = cms.ESSource(
    'LXXXCorrectionService',
    era = cms.string('Jec11V12'),
    section   = cms.string(''),
    level     = cms.string('L2Relative'),
    # the above 3 elements are needed only when the service is initialized from local txt files
    algorithm = cms.string('AK5PF'),
    # the 'algorithm' tag is also the name of the DB payload
    useCondDB = cms.untracked.bool(True)
    )
)

Modular application: corrected jet collection (UPDATED for CMSSW 7.6.X)

Once a jet correction service has been defined, as described above, the correction can be applied to create a corrected jet collection. More specifically, in every event all the jets are corrected and put into the event as a new collection. The corrected jet collection is sorted in decreasing pt. With this method, analysis must be done using the corrected jet collection. A typical module that produces the correction is the following:

process.load('JetMETCorrections.Configuration.JetCorrectors_cff')

process.ak4PFchsCorrectedJets   = cms.EDProducer('CorrectedPFJetProducer',
    src         = cms.InputTag('ak4PFchsJets'),
    correctors  = cms.VInputTag('ak4PFCHSL1FastL2L3Corrector')
    )

If you are not using unscheduled mode, you'll need to add process.ak4PFchsCorrectedJets to your path, as well the chain of correctors before your analyzer / producer. Relevant chain of correctors are process.ak4PFCHSL1FastL2L3CorrectorChain or process.ak4PFCHSL1FastL2L3ResidualCorrectorChain.

The corrected collection producer is templated to account for CaloJets, PFJets, JPTJets and TrackJets. The instances are named CorrectedCaloJetProducer, CorrectedPFJetProducer, CorrectedJPTJetProducer, CorrectedTrackJetProducer. In the above definition, the src is the name of the jet collection where the correction will be applied and the correctors is the name of the predefined corrections. You can find a list of all these corrections here.

Correcting jets "on-the-fly" (UPDATED for CMSSW 7.6.X)

You can also apply the correction jet by jet directly in your analysis code. Two parts are needed for this to work.

  • First, you need to load the necessary jet corrector you want to use. These are now independent modules which must be ran before your analyzers. This line is generally enough if you run in unscheduled mode.

process.load('JetMETCorrections.Configuration.JetCorrectors_cff')

You can find the complete list of all the predefined modules here. If you are not running in unscheduled mode, then you also need to add the sequence of correctors to your path. The complete list is available on the previous link, but relevant sequences would probably be process.ak4PFCHSL1FastL2L3CorrectorChain and process.ak4PFCHSL1FastL2L3ResidualCorrectorChain.

  • In your analyzer, you need to retrieve the jet corrector module. You first have to get a token to the module through the consume interface, and them retrieve the content of the token in the event loop.
    • The needed include file is:
#include "JetMETCorrections/JetCorrector/interface/JetCorrector.h"
    • The token you need:
edm::EDGetTokenT<reco::JetCorrector> mJetCorrector;
    • Retrieve the token (usually in your analyzer constructor). You'll need to replace CORRECTOR_LABEL by the name of the jet corrector you want to use (something like ak4PFCHSL1FastL2L3Corrector or ak4PFCHSL1FastL2L3ResidualCorrector):
mJetCorrector  = consumes<reco::JetCorrector>(edm::InputTag("CORRECTOR_LABEL"));
    • Inside your analyze or produce methode, retrieve the jet corrector
edm::Handle<reco::JetCorrector> corrector;
iEvent.getByToken(mJetCorrector, corrector);
    • Use it:
for (const auto& jet: jets) {
    double jec = corrector->correction(jet);
}

IMPORTANT: after applying corrections "on the fly", the resulting jets are NOT sorted by pt.

Re-correcting pat::Jets and jets from MiniAOD (UPDATED for CMSSW 7.6.4)

This recipe can be used to apply a new set of jet energy corrections, to the already corrected pat::Jets or to slimmedJets from MiniAOD, starting from CMSSW_7_4_0. The pat::Jets or slimmedJets already contain a set a jet energy correction. The following code reverts these corrections and applies the new ones, producing as output an new collection of pat::Jets which have the same content as the input jet, but new jet energy corrections applied. These modules can for example be placed in the configuration sequence before your ntuplizer.

There's two different recipes based on your CMSSW version. In both case, you first need to load a global tag or sqlite file file containing the jet energy correction if you want to see any difference. Instruction and recommendations for this are given on this TWiki page: WorkBookJetEnergyCorrections#JecGlobalTag, WorkBookJetEnergyCorrections#JecSqliteFile.

Before CMSSW 7.6.4

You need to add the following python code into your configuration. The new jets collection, with the updated JEC, is patJetsReapplyJEC

from PhysicsTools.PatAlgos.producersLayer1.jetUpdater_cff import patJetCorrFactorsUpdated
process.patJetCorrFactorsReapplyJEC = process.patJetCorrFactorsUpdated.clone(
  src = cms.InputTag("slimmedJets"),
  levels = ['L1FastJet', 
        'L2Relative', 
        'L3Absolute'],
  payload = 'AK4PFchs' ) # Make sure to choose the appropriate levels and payload here!

from PhysicsTools.PatAlgos.producersLayer1.jetUpdater_cff import patJetsUpdated
process.patJetsReapplyJEC = process.patJetsUpdated.clone(
  jetSource = cms.InputTag("slimmedJets"),
  jetCorrFactorsSource = cms.VInputTag(cms.InputTag("patJetCorrFactorsReapplyJEC"))
  )

process.p += cms.Sequence( process.patJetCorrFactorsReapplyJEC + process.patJetsReapplyJEC )

CMSSW 7.6.4 and above

Starting from CMSSW_7_6_4, the updateJetCollection function was added to PhysicsTools/PatAlgos/python/tools/jetTools.py which allows easy re-correction of pat::Jets as well as adding extra information to them. An example of its use can be found in PhysicsTools/PatAlgos/test/patTuple_updateJets_fromMiniAOD_cfg.py. The function is by default configured for pat::Jets stored in MiniAOD but can be used for any collection of pat::Jets with appropriately set input arguments.

If you just want to update the JEC in the MiniAOD, using the following python code is enough:

from PhysicsTools.PatAlgos.tools.jetTools import updateJetCollection

updateJetCollection(
   process,
   jetSource = cms.InputTag('slimmedJets'),
   labelName = 'UpdatedJEC',
   jetCorrections = ('AK4PFchs', cms.vstring(['L1FastJet', 'L2Relative', 'L3Absolute', 'L2L3Residual']), 'None')  # Update: Safe to always add 'L2L3Residual' as MC contains dummy L2L3Residual corrections (always set to 1)
)

The new jets collection, with the updated JEC, is updatedPatJetsUpdatedJEC, where by default no selection is applied (as in the standard PAT jet workflow). More generally, the final updated jet collection will be called updatedPatJets+labelName+postfix, with labelName and postfix set by default to an empty string.

If you are not using unscheduled mode (you really should using it), you need to add the following modules to your execution path (can be different if you change the labelName):

process.jecSequence = cms.Sequence(process.patJetCorrFactorsUpdatedJEC * process.updatedPatJetsUpdatedJEC)

Accessing the JEC uncertainties

There are two ways to access the JEC uncertainties:

A) from the data base, directly in an EDAnalyzer:

edm::ESHandle<JetCorrectorParametersCollection> JetCorParColl;
iSetup.get<JetCorrectionsRecord>().get("AK5PF",JetCorParColl); 
JetCorrectorParameters const & JetCorPar = (*JetCorParColl)["Uncertainty"];
JetCorrectionUncertainty *jecUnc = new JetCorrectionUncertainty(JetCorPar);

B) from a local txt file

JetCorrectionUncertainty *jecUnc = new JetCorrectionUncertainty("Jec11_V12_Uncertainty_AK5PF.txt")

Note: JEC uncertainty sources (i.e. correlations, for advanced users) for final 2011 JEC are only available as text files at https://twiki.cern.ch/twiki/bin/view/CMS/JECUncertaintySources.

Please use for 2012 data the new uncertainty files that can be found at https://twiki.cern.ch/twiki/bin/view/CMS/JECUncertaintySources#2012_JEC.

In both cases above, the JetCorrectionUncertainty object must be initialized once but used inside the jet loop, by passing each time the jet pt (corrected) and eta:

for (PFJetCollection::const_iterator jet = jets->begin(); jet != jets->end(); jet++)  {
  ...............................
  jecUnc->setJetEta(eta);
  jecUnc->setJetPt(ptCor); // here you must use the CORRECTED jet pt
  double unc = jecUnc->getUncertainty(true);
  double ptCor_shifted = ptCor(1+shift*unc) ; // shift = +1(up), or -1(down)
}

Factorized Approach

IMPORTANT: NEW naming scheme for CMSSW_3_6_X and later. Calorimeter jets and Particle flow jets are supported, for all the default reconstructed jet algorithms in the official samples (IC5,KT4,KT6,AK5,AK7). In addition there are also available corrections on top of the JetPlusTrack algorithm and for TrackJets. The names of the corrected jet collections for the 7 TeV data and MC analysis with CMSSW_3_6_X are listed below:

  • kt4CaloJetsL2L3 (KT D=0.4 CaloJets corrected with L2+L3)
  • kt4PFJetsL2L3 (KT D=0.4 PFlow Jets corrected with L2+L3)
  • kt6CaloJetsL2L3 (KT D=0.6 CaloJets corrected with L2+L3)
  • kt6CaloJetsL2L3 (KT D=0.6 PFlow Jets corrected with L2+L3)
  • ak5CaloJetsL2L3 (anti-KT D=0.5 Calo Jets corrected with L2+L3)
  • ak5PFJetsL2L3 (anti-KT D=0.5 PFlow Jets corrected with L2+L3)
  • ak5JPTJetsL2L3 (anti-KT D=0.5 JPT Jets corrected with L2+L3)
  • ak5TrackJetsL2L3 (anti-KT D=0.5 Track Jets corrected with L2+L3)
  • ak7CaloJetsL2L3 (anti-KT D=0.7 Calo Jets corrected with L2+L3)
  • ak7PFJetsL2L3 (anti-KT D=0.7 PFlow Jets corrected with L2+L3)

Applying the Jet Energy Corrections

Instructions

Naming Scheme

Examples

Illustrated examples can be found in the RecoJets/JetAnalyzers package:

V00-04-09 RecoJets/JetAnalyzers

The default correction example which should be used by most users is runL2L3JetCorrectionExample_cfg.py in RecoJets/JetAnalyzers/test . Putting together all the above instructions, the full prescription for running this example is:

cmsrel CMSSW_3_3_4
cd CMSSW_3_3_4/src
cmsenv
cvs co -r V00-04-09 RecoJets/JetAnalyzers
scram b
cd RecoJets/JetAnalyzers/test
cmsRun runL2L3JetCorrectionExample_cfg.py

The above example runs on 100 events, corrects all the CMS.CaloJets in the events, histograms the two leading CMS.CaloJets, and writes out these histograms files:

  • CaloJetHisto.root contains histograms of the two leading CMS.CaloJets before L2+L3 corrrections.
  • CorJetHisto.root contains histograms of the two leading CMS.CaloJets after L2+L3 corrections.

Please see the Workbook pages on jet energy corrections for more details on the example and instructions for other CMSSW releases.

Jet Identification

The cuts to separate real jets from fakes are still evolving. Please look at the corresponding page for more details and the latest recommendations.

Advanced Topics

Jet Energy Corrections in FWLite

Simple jet corrections, that are corrections which are dependent from jet properties only, may be directly accessed from ROOT scripts running in FWLite environment. The class that delivers the jet correction is called CombinedJetCorrector and is located in the CondFormats/JetMETObjects package. In order to get the class, the following cvs checkout is necessary (works in the CMSSW_2_2_X releases and above):

cvs co -r V01-08-03 CondFormats/JetMETObjects

Once the package is checked out and compiled, the CombinedJetCorrector class will be loaded in a ROOT macro with the rest of the FWLite libraries.

  • Step by step instructions for use in a ROOT macro.
    • Step1: load the FWLite libraries:
   gSystem->Load("libFWCoreFWLite.so");
   AutoLibraryLoader::enable(); 
    • Step2: define an object of the class CombinedJetCorrector . There are 2 constructors:
   string Levels = "L2:L3";
   string Tags = "Summer08Redigi_L2Relative_IC5Calo:Summer08Redigi_L3Absolute_IC5Calo";
   CombinedJetCorrector *L2L3JetCorrector = new CombinedJetCorrector(Levels,Tags);
or
   string Levels = "L2:L3:L7";
   string Tags = "Summer08Redigi_L2Relative_IC5Calo:Summer08Redigi_L3Absolute_IC5Calo:L7Parton_IC5";
   string Options = "Parton:gJ";
   CombinedJetCorrector *L2L3L7JetCorrector = new CombinedJetCorrector(Levels,Tags,Options);
  • The string Levels defines which corrections will be applied. IMPORTANT: the order maters !!!!.
  • The string Tags defines the specific tags for each level of correction. The available tags can be found in CondFormats/JetMETObjects/data
  • The string Options defines the L5 Flavor option ("qJ","gJ","cJ","bJ","gT","qT","cT","bT") and/or the L7 Parton Option ("qJ","gJ","cJ","bJ","jJ","qT","cT","bT","tT"). If both L5 and L7 are requested:
   string Options = "Parton:gJ & Flavor:gJ"
    • Step3: get the correction factor or the vector with the individual corrections:
   double correction = L2L3JetCorrector->getCorrection(pt,eta,energy);
   vector<double> factors = L2L3JetCorrector->getSubCorrections(pt,eta,energy);
If the L4 EMF correction is requested, the above commands become:
   double correction = L2L3JetCorrector->getCorrection(pt,eta,energy,emf);
   vector<double> factors = L2L3JetCorrector->getSubCorrections(pt,eta,energy,emf);
IMPORTANT: the getSubCorrections member function returns the vector of the subcorrections UP to the given level. For example in the second constructor above (with L2, L3 & L7), factors[0] is the L2 correction, factors[1] is the L2+L3 correction and factors[2] is the L2+L3+L7 correction.

A full example ROOT macro is: CondFormats/JetMETObjects/test/TestCorrections.C

Jet Energy Corrections in PAT

Please see the corresponding page for information on PAT. Also see here for information on how to change the JEC in the PAT. Also see here for information on retrieving jet correction values from the pat::Jet.

Jet Energy Corrections using Jet+tracks algorithm

Contribution of Sasha Nikitenko and Olga Kodolova

For detailed information please look at the corresponding pages.

No permission to view CMS.JetPlusTracksCorrections

True Jet Flavour Definition

Contribution of Attilio Santocchia

Please see SWGuideBTagMCTools.

References

Review status

Reviewer/Editor and Date (copy from screen) Comments
HartmutStadie - 13-Oct-2009 start with copy from 210 version
HartmutStadie - 21-Dec-2009 make jec and jet id own topics, remove outdated examples

Responsible: HartmutStadie

-- HartmutStadie - 13-Oct-2009

Topic attachments
I Attachment History Action Size Date Who Comment
PNGpng piGun_CaloJet_fwlite.png r1 manage 49.6 K 2009-11-30 - 18:17 HartmutStadie  
PNGpng piGun_CaloJets.png r1 manage 44.4 K 2009-11-04 - 11:46 HartmutStadie  
PNGpng piGun_ptDist.png r1 manage 20.9 K 2009-11-04 - 11:48 HartmutStadie  
PNGpng pigunObjects.png r1 manage 48.4 K 2009-11-04 - 11:45 HartmutStadie  
Edit | Attach | Watch | Print version | History: r16 < r15 < r14 < r13 < r12 | Backlinks | Raw View | WYSIWYG | More topic actions
Topic revision: r16 - 2011-02-18 - BenjaminRadburnSmith
 
    • 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