PAT Examples: Jet Example

Contents

Introduction

This tutorial will be concerned with the performance of the Jet Energy Scale Calibration and the reconstruction of jets from different input collections. The following topics will be discussed:

  • Basic jet quantities for calorimeter based jets and particle flow jets.
  • Sensitivity of the invariant dijet mass on shifts of the jet energy scale.
  • The response of calorimeter based jets and particle flow jets at different jet energy correction levels.

We will make use of a standard ttbar input sample for our use cases. You can find the corresponding files on DBS here. We make use of a total of roughly 20k events. You will be guided through the exercises step by step. Questions will help you to probe your understanding of each section.

How to get the code

First of all connect to lxplus and go to some working directory. You can choose any directory, provided that you have enough space. You need ~5 MB of free disc space for this exercise. We recommend you to use your ~/scratch0 space. In case you don't have this (or do not even know what it is) check your quota typing fs lq and follow this link. If you don't have enough space, you may instead use the temporary space (/tmp/your_user_name), but be aware that this is lost once you log out of lxplus (or within something like a day). We will expect in the following that you have such a ~/scratch0 directory.

ssh lxplus
[ ... enter password ... ]
cd scratch0/

Create a directory for this exercise (to avoid interference with code from the other exercises).

mkdir jetExercise
cd jetExercise

Create a local release area and enter it.

cmsrel CMSSW_5_3_14
cd CMSSW_5_3_14/src 

The first command creates all directories needed in a local release area. Setting up the environment is done invoking the following script:

cmsenv

Sensitivity of basic quantities to shifts in the JES

In a first exercise we want explore a few basic quantities of jets and their sensitivity to shifts in the jet energy scale. For this we make use of a PAT tuple that has been pre-produced with the patTuple_standard_cfg.py configuration file as located in the test directory of the PatAlgos package. You can find this file in any CMSSW release. It contains the configuration to produce a standard PAT tuple. As a small extension we added a collection of AK5 jets based on particle flow candidates in addition to the same jet collection based on pure calorimeter information (which is the default collection in PAT). From this pre-processed PAT tuple you should expect the following event content:

[rwolf@tcx112]~/data/CMSSW_3_8_4/src% edmDumpEventContent patTuple.root 
edm::OwnVector >     "selectedPatJets"       "tagInfos"    "PAT."         
edm::OwnVector >     "selectedPatJetsAK5PF"    "tagInfos"    "PAT."         
vector                 "selectedPatJets"       "caloTowers"    "PAT."         
vector                 "selectedPatJetsAK5PF"    "caloTowers"    "PAT."         
vector             "cleanPatElectrons"     ""            "PAT."         
vector                  "cleanPatJets"          ""            "PAT."         
vector                  "cleanPatJetsAK5PF"     ""            "PAT."         
vector                  "patMETs"               ""            "PAT."         
vector                 "cleanPatMuons"         ""            "PAT."         
vector               "cleanPatPhotons"       ""            "PAT."         
vector                  "cleanPatTaus"          ""            "PAT."         
vector              "selectedPatJets"       "genJets"     "PAT."         
vector              "selectedPatJetsAK5PF"    "genJets"     "PAT."         
vector         "selectedPatJets"       "pfCandidates"    "PAT."         
vector         "selectedPatJetsAK5PF"    "pfCandidates"    "PAT." 

ALERT! Note: While the standard jet collection is retrieved via the cleanPatJets label, the additional jet collection is labels as cleanPatJetsAK5PF, where the postfix indicated the the jet algorithm and input collection. For the standard jet collection the AK5 algorithm and calorimeter towers as input are used.

To do the analysis you need to check out and compile the PatExamples package as given below:

addpkg PhysicsTools/PatExamples V00-04-16
scram b PhysicsTools/PatExamples

ALERT! Note: You might want to use more than one core to compile the package making use of the scram compiler flag -j followed by the number of cores you would like to make use of.

You run the analysis starting cmsRun with the following configuration file:

rwolf@tcx112]~/data/CMSSW_3_8_4/src% cmsRun PhysicsTools/PatExamples/test/analyzePatShiftedJets_cfg.py
28-Sep-2010 22:10:00 CEST  Initiating request to open file file:patTuple_sec.root
28-Sep-2010 22:10:03 CEST  Successfully opened file file:patTuple_sec.root
Begin processing the 1st record. Run 1, Event 559002, LumiSection 1 at 28-Sep-2010 22:10:05 CEST
Begin processing the 2nd record. Run 1, Event 559008, LumiSection 1 at 28-Sep-2010 22:10:05 CEST
Begin processing the 3rd record. Run 1, Event 559010, LumiSection 1 at 28-Sep-2010 22:10:05 CEST
Begin processing the 4th record. Run 1, Event 559013, LumiSection 1 at 28-Sep-2010 22:10:05 CEST
Begin processing the 5th record. Run 1, Event 559014, LumiSection 1 at 28-Sep-2010 22:10:05 CEST
Begin processing the 6th record. Run 1, Event 559019, LumiSection 1 at 28-Sep-2010 22:10:05 CEST
Begin processing the 7th record. Run 1, Event 559020, LumiSection 1 at 28-Sep-2010 22:10:05 CEST
Begin processing the 8th record. Run 1, Event 559023, LumiSection 1 at 28-Sep-2010 22:10:05 CEST
Begin processing the 9th record. Run 1, Event 559024, LumiSection 1 at 28-Sep-2010 22:10:05 CEST
...

While your job is processing you might want to have a look into the configuration file for this exercise:

import FWCore.ParameterSet.Config as cms

process = cms.Process("EnergyShift")

## Declare input
from CMS.PhysicsTools.PatExamples.samplesDESY_cff import *

process.source = cms.Source("PoolSource",
  fileNames = ttbarJets
)

process.maxEvents = cms.untracked.PSet(
    input = cms.untracked.int32( -1 )
)

## Message logger configuration
process.load("FWCore.MessageLogger.MessageLogger_cfi")

## Configure jet energy scaling module
process.load("CMS.PhysicsTools.PatExamples.JetEnergyShift_cfi")
process.scaledJets.scaleFactor =  1.0

## Select good jets
from CMS.PhysicsTools.PatAlgos.selectionLayer1.jetSelector_cfi import *
process.goodJets = selectedPatJets.clone(
    src="scaledJets:cleanPatJets",
    cut = 'abs(eta) < 3 & pt > 30. &'
    'emEnergyFraction > 0.01       &'
    'jetID.fHPD < 0.98             &'
    'jetID.n90Hits > 1'    
    )

## Analyze jets
process.load("CMS.PhysicsTools.PatExamples.PatJetAnalyzer_cfi")
process.analyzePatJets.src = 'goodJets'

## Define output file
process.TFileService = cms.Service("TFileService",
  fileName = cms.string('analyzeEnergyShiftedJets.root')
)

process.p = cms.Path(
    process.scaledJets *
    process.goodJets   * 
    process.analyzePatJets
)

ALERT! Note: After the usual procedure of defining the input files and the number of events to be processed we configure a module scaledJets which is defined in the JetEnergyShift_cfi.py file. You can find this file in the python directory of the PatExamples package. Its are given below:

scaledJets = cms.EDProducer("JetEnergyShift",
    inputJets            = cms.InputTag("cleanPatJets"),
    inputMETs            = cms.InputTag("patMETs"),
    scaleFactor          = cms.double(1.0),
    jetPTThresholdForMET = cms.double(20.),
    jetEMLimitForMET     = cms.double(0.9)
)

We have an input collection of jets (inputJets) and missing transverse momentum (inputMETs) a scale factor (scaleFactor) and a set of technical parameters, concerned with the type1 JES correction scheme for raw calorimeter based MET. We use this module to produce new collections of jets and MET from the corresponding input collections, where the jet energy scale has been shifted according to the value given by the parameter scaleFactor. This parameter should be interpreted as the relative shift of the jet energy scale (JES). Accordingly a value of scaleFactor=1.0 will not shift the the jet energy scale at all, while a value of scaleFactor=1.1 will shift the jet energy scale up by 10%. As the calorimeter based MET has been corrected for the JES it has to be modified as well to get things right. The output of the module consists of two collections:

  • scaledJets:cleanPatJets, which is the scaled jet collection.
  • scaleJets:patMETs, which is the recalculated MET collection.

ALERT! Note. the module label for both collections is the same (scaledJets). As the module produces two collections these are distinguished by an extra instance label, which corresponds to the module labels of the input collections. Most collections that you know about do not have an instance label, but just a module label. Here you have one of the rare cases, where an instance label is indeed necessary. When replacing input collections by collections, which carry an additional instance label you can appended the instance label to the module label string separated by a ':'. You can see such a parameter replacement in the following lines of the analyzePatShiftedJets_cfg.py file, where we pass on the shifted jet collection to a jet selector module to apply a typical selection:

from CMS.PhysicsTools.PatAlgos.selectionLayer1.jetSelector_cfi import *
process.goodJets = selectedPatJets.clone(
    src="scaledJets:cleanPatJets",
    cut = 'abs(eta) < 3 & pt > 30. &'
    'emEnergyFraction > 0.01       &'
    'jetID.fHPD < 0.98             &'
    'jetID.n90Hits > 1'    
    )

ALERT! Note: apart from the kinematics this selection corresponds to the tight jetID for calorimeter based jets. We also restrict ourselves to jets with pt>30 GeV in the central detector. From that point on we continue dealing with the goodJets collection that we feed into the EDAnalyzer module:

process.load("CMS.PhysicsTools.PatExamples.PatJetAnalyzer_cfi")
process.analyzePatJets.src = 'goodJets'

You can find the default configuration of this module in the PatJetAnalyzer_cfi.py file in the python directory of the PatExamples package:

analyzePatJets = cms.EDAnalyzer("PatJetAnalyzer",
  src  = cms.InputTag("cleanPatJets"),                              
  corrLevel = cms.string("abs")
)  

ALERT! Note: You see that the module expects two parameters:

  • src, which defines the jet collection to be analysed.
  • corrLevel, which is the level of jet energy corrections we want to apply to the jets.

You can find all available jet correction levels, which are supported by the CMS.JetMET PAG here. Per default all jet correction levels are stored with the PAT jet and the jet is corrected up to the L3Absolute jet energy correction level as recommended by the JetMET PAG.

You can find the implementation of the module in the PatJetAnalyzer.cc file in the plugins directory of the PatExamples package. The implementation is fairly simple. We will discuss the constructor first:

PatJetAnalyzer::PatJetAnalyzer(const edm::ParameterSet& cfg):
  corrLevel_(cfg.getParameter<std::string>("corrLevel")),
  jets_(cfg.getParameter<edm::InputTag>("src"))
{
  // register TFileService
  edm::Service<TFileService> fs;

  // jet multiplicity
  hists_["mult" ]=fs->make<TH1F>("mult" , "N_{Jet}"          ,   15,   0.,   15.);
  // jet pt (for all jets)
  hists_["pt"   ]=fs->make<TH1F>("pt"   , "p_{T}(Jet) [GeV]" ,   60,   0.,  300.);
  // jet eta (for all jets)
  hists_["eta"  ]=fs->make<TH1F>("eta"  , "#eta (Jet)"       ,   60,  -3.,    3.);
  // jet phi (for all jets)
  hists_["phi"  ]=fs->make<TH1F>("phi"  , "#phi (Jet)"       ,   60,  3.2,   3.2);
  // dijet mass (if available)
  hists_["mass" ]=fs->make<TH1F>("mass" , "M_{jj} [GeV]"     ,   50,   0.,  500.);
  // basic histograms for jet energy response
  for(unsigned int idx=0; idx<MAXBIN; ++idx){
    char buffer [10]; sprintf (buffer, "jes_%i", idx);
    char title  [50]; sprintf (title , "p_{T}^{rec}/p_{T}^{gen} [%i GeV - %i GeV]", (int)BINS[idx], (int)BINS[idx+1]);
    hists_[buffer]=fs->make<TH1F>(buffer, title,  100, 0., 2.);
  }  
}

ALERT! Note: at the very beginning you see the parameters src and corrLevel being passed on from the configuration file. Then histograms for the pt, eta and phi of the jets are booked, which will be filled for each jet in the collection, a histogram for the invariant dijet mass that will be calculated from the two leading jets and a set of histograms for the pt of the reconstructed jet over the pt of the associated generated jet (as reconstructed from stable hadrons at generator level).

In the analyze function you can see how the histograms are filled:

void
PatJetAnalyzer::analyze(const edm::Event& event, const edm::EventSetup& setup)
{
  // recieve jet collection label
  edm::Handle<edm::View<pat::Jet> > jets;
  event.getByLabel(jets_,jets);

  // loop jets
  for(edm::View<pat::Jet>::const_iterator jet=jets->begin(); jet!=jets->end(); ++jet){
    // fill basic kinematics
    fill( "pt" , jet->correctedJet(corrLevel_).pt());
    fill( "eta", jet->eta());
    fill( "phi", jet->phi());
    // basic plots for jet responds plot as a function of pt
    if( jet->genJet() ){
      double resp=jet->correctedJet(corrLevel_).pt()/jet->genJet()->pt();
      for(unsigned int idx=0; idx<MAXBIN; ++idx){
	if(BINS[idx]<=jet->genJet()->pt() && jet->genJet()->pt()<BINS[idx+1]){
	  char buffer [10]; sprintf (buffer, "jes_%i", idx);
	  fill( buffer, resp );
	}
      }
    }
  }
  // jet multiplicity
  fill( "mult" , jets->size());
  // invariane dijet mass
  if(jets->size()>2){ fill( "mass", ((*jets)[0].p4()+(*jets)[1].p4()).mass());}
}

ALERT! Note: this example shows the easy access both to the associated generator jet (if anassociation exists) as well as to the corresponding jet energy correction levels. The histograms for the jet energy response are filled in bins of the pt of the generated jet. These bins are defined as static const at the top of the file. Once your cmsRun job is finished open the analyzeEnergyShiftedJets.root file with a root Browser and answer the following questions:

Question Question a):
What is the mean value of pt and the invariant dijet mass w/o energy shift? How does the mean of the invariant dijet mass change if you redo the exercise and shift the jet energy scale by +5% or -5%?

Question Question b):
How does the number of all jets, that pass the goodJet selection change when applying the shifts (you can check this from the number of entries in the phi histogram)? How does the mean of the jet multiplicity plot change?

Question Question c):
How would you have to modify the analyer module to check the influence of the shift in the JES on the MET correction?

Jet Response for calorimeter based and particle flow jets

We will use the same EDAnalyzer module in a slightly different configuration to study the energy response for Uncorrected jets and for jets corrected up to the jet energy correction levels L2Relative and L3Absolute. To do this run the analyzePatJetEnergyScale_cfg.py configuration file. You can fin it in the test directory of the the PatExamples package:

rwolf@tcx112]~/data/CMSSW_3_8_4/src% cmsRun PhysicsTools/PatExamples/test/analyzePatJetEnergyScale_cfg.py 
28-Sep-2010 23:47:00 CEST  Initiating request to open file file:patTuple_sec.root
28-Sep-2010 23:47:03 CEST  Successfully opened file file:patTuple_sec.root
Begin processing the 1st record. Run 1, Event 559002, LumiSection 1 at 28-Sep-2010 23:47:05 CEST
Begin processing the 2nd record. Run 1, Event 559008, LumiSection 1 at 28-Sep-2010 23:47:05 CEST
Begin processing the 3rd record. Run 1, Event 559010, LumiSection 1 at 28-Sep-2010 23:47:05 CEST
Begin processing the 4th record. Run 1, Event 559013, LumiSection 1 at 28-Sep-2010 23:47:05 CEST
Begin processing the 5th record. Run 1, Event 559014, LumiSection 1 at 28-Sep-2010 23:47:05 CEST
Begin processing the 6th record. Run 1, Event 559019, LumiSection 1 at 28-Sep-2010 23:47:05 CEST
Begin processing the 7th record. Run 1, Event 559020, LumiSection 1 at 28-Sep-2010 23:47:06 CEST
Begin processing the 8th record. Run 1, Event 559023, LumiSection 1 at 28-Sep-2010 23:47:06 CEST
Begin processing the 9th record. Run 1, Event 559024, LumiSection 1 at 28-Sep-2010 23:47:06 CEST
...

We will have a closer look into the configuration file below:

import FWCore.ParameterSet.Config as cms

process = cms.Process("EnergyScale")

## Declare input
from CMS.PhysicsTools.PatExamples.samplesDESY_cff import *

process.source = cms.Source("PoolSource",
  fileNames = ttbarJets
)

process.maxEvents = cms.untracked.PSet(
    input = cms.untracked.int32( -1 )
)

## Message logger configuration
process.load("FWCore.MessageLogger.MessageLogger_cfi")

## Select good jets
from CMS.PhysicsTools.PatAlgos.selectionLayer1.jetSelector_cfi import *
## use this selection of calorimeter based jets
process.goodJets = selectedPatJets.clone(
    src="cleanPatJets",
    cut = 'abs(eta) < 3 & pt > 30. &'
    'emEnergyFraction > 0.01       &'
    'jetID.fHPD < 0.98             &'
    'jetID.n90Hits > 1'    
)

## use this selection of particle flow jets
#process. goodJets   = selectedPatJets.clone(
#    src = 'cleanPatJetsAK5PF',
#    cut = 'abs(eta) < 2.4 & pt > 30.          &'
#    'chargedHadronEnergyFraction > 0.0  &'
#    'neutralHadronEnergyFraction/corrFactor("raw") < 0.99 &'
#    'chargedEmEnergyFraction/corrFactor("raw")     < 0.99 &'
#    'neutralEmEnergyFraction/corrFactor("raw")     < 0.99 &'
#    'chargedMultiplicity > 0            &'
#    'nConstituents > 0'
#)

## Analyze jets
from CMS.PhysicsTools.PatExamples.PatJetAnalyzer_cfi import analyzePatJets
process.Uncorrected = analyzePatJets.clone(src = 'goodJets', corrLevel='raw')
process.L2Relative  = analyzePatJets.clone(src = 'goodJets', corrLevel='rel')
process.L3Absolute  = analyzePatJets.clone(src = 'goodJets', corrLevel='abs')

## Define output file
process.TFileService = cms.Service("TFileService",
  fileName = cms.string('analyzeJetEnergyScale.root')
)

process.p = cms.Path(
    process.goodJets    * 
    process.Uncorrected *
    process.L2Relative  *
    process.L3Absolute  
)

ALERT! Note: instead of a single instance of the analyzePatJets module we create three instances, labeled as:

  • Uncorrected
  • L2Relative
  • L3Absolute

We do this by cloning the module and replacing the parameter corrLevel. The input collection should be the goodJets collection in all cases.

Question Exercise a):
As soon as your cmsRun job is finished open root and run the macro monitorJetEnergyScale.C, which is located in the bin directory of the PatExamples package:

[rwolf@tcx112]~/data/CMSSW_3_8_4/src% root -l
root [0] .x PhysicsTools/PatExamples/bin/monitorJetEnergyScale.C+
Info in <TUnixSystem::ACLiC>: creating shared library /scratch/hh/lustre/cms/user/rwolf/CMSSW_3_8_4/src/./CMS.PhysicsTools/PatExamples/bin/monitorJetEnergyScale_C.so
/usr/bin/ld: skipping incompatible /usr/lib64/libm.so when searching for -lm
/usr/bin/ld: skipping incompatible /usr/lib64/libm.a when searching for -lm
/usr/bin/ld: skipping incompatible /usr/lib64/libc.so when searching for -lc
/usr/bin/ld: skipping incompatible /usr/lib64/libc.a when searching for -lc
<TCanvas::MakeDefCanvas>: created default TCanvas with name c1
 FCN=19.6573 FROM MIGRAD    STATUS=CONVERGED     145 CALLS         146 TOTAL
                     EDM=3.19386e-09    STRATEGY= 1      ERROR MATRIX ACCURATE 
  EXT PARAMETER                                   STEP         FIRST   
  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE 
   1  Constant     1.67527e+00   3.50393e-01   5.85318e-04   8.82597e-05
   2  Mean         5.85842e-01   1.67557e-01   3.27498e-04  -2.55278e-04
   3  Sigma        4.14015e-01   3.08051e-01   6.70345e-04  -1.38701e-04
 FCN=12.6336 FROM MIGRAD    STATUS=CONVERGED     136 CALLS         137 TOTAL
                     EDM=6.30454e-10    STRATEGY= 1      ERROR MATRIX ACCURATE 
...

This will result in a jet response plot for calorimeter based AK5 jets. You can have a closer look into the implementation of this macro to learn how these plots are determined from the basic histograms that have been booked and filled in the analyzePatJets module.

Question Exercise b):
Do the same exercise but replace the calorimeter based jets by particle flow based jets. Note that you should then apply a different jetID in you configration file. You can find the details below:

goodJets   = selectedPatJets.clone(
       src = 'cleanPatJetsAK5PF',
       cut = 'abs(eta) < 2.4 & pt > 30.          &'
                'chargedHadronEnergyFraction > 0.0  &'
                'neutralHadronEnergyFraction/corrFactor("raw") < 0.99 &'
                'chargedEmEnergyFraction/corrFactor("raw")     < 0.99 &'
                'neutralEmEnergyFraction/corrFactor("raw")     < 0.99 &'
                'chargedMultiplicity > 0            &'
                'nConstituents > 0'
)

This corresponds to the tight jetID for particle flow jets. Running the root macro should result in plots of this kind:

JetEnergyResponse_calo.png JetEnergyResponse_pflow.png

(left/upper) Jet energy response for calorimeter based jets. (right/lower) Jet energy response for particle flow based jets.

Question Question d):
Can you explain what you see in both plots? Can you explain where the difference of the two plots originates from?

Question Question e):
What would you have to change in the root macro to plot the resolution as a function of the pt of the generated jet?

Question Question f):
What would you have to change in the analyzePatJets module to plot the jet energy response as function of the eta of the generator jet and what would you expect for the L2Relative correction?

Review status

Reviewer/Editor and Date (copy from screen) Comments
RogerWolf - 29 Sept 2010 first version for German Physics School

Responsible: RogerWolf

Topic attachments
I Attachment History Action Size Date Who Comment
PNGpng JetEnergyResponse_calo.png r1 manage 17.3 K 2010-09-29 - 00:49 RogerWolf  
PNGpng JetEnergyResponse_pflow.png r1 manage 17.1 K 2010-09-29 - 00:51 RogerWolf  
Cascading Style Sheet filecss tutorial.css r1 manage 0.2 K 2010-09-28 - 21:59 RogerWolf  
Edit | Attach | Watch | Print version | History: r7 < r6 < r5 < r4 < r3 | Backlinks | Raw View | WYSIWYG | More topic actions
Topic revision: r7 - 2010-09-29 - RogerWolf
 
    • 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