Analysis Starter Kit PAT Tutorial, 14May08

This is a page dedicated to the questions asked in advance for the Analysis Starter Kit and PAT tutorials on 14May08 at the CMS Run Plan Workshop at Fermilab.

The Indico webpage is here.

Usage of jet energy corrections (Ia Iashvili)

Jet energy corrections are applied to PAT objects by default. Currently we are using the (slightly out of date) prescription to use the MCJet corrections, and flavor-dependent corrections. We are looking for a volunteer to update to the latest Jet/MET group prescription.

In slightly more detail, we produce the jet energy corrections at layer 0 here

module jetCorrFactors = JetCorrFactorsProducer {
  InputTag jetSource = iterativeCone5CMS.CaloJets
  string defaultJetCorrector = "MCJetCorrectorIcone5"
  string udsJetCorrector = "L5FlavorJetCorrectorUds"
  string gluonJetCorrector = "L5FlavorJetCorrectorGluon"
  string cJetCorrector = "L5FlavorJetCorrectorC"
  string bJetCorrector = "L5FlavorJetCorrectorB"
}

These correction factors are then stored in a value map at layer 0.

At layer 1, all this information is collated and the access to the user is via the Jet class:

      /// return the correction factor to go to a non-calibrated jet
      JetCorrFactors jetCorrFactors() const;
      /// return the original non-calibrated jet
      JetType recJet() const;
      /// return the associated non-calibrated jet
      Jet noCorrJet() const;
      /// return the associated default-calibrated jet
      Jet defaultCorrJet() const;
      /// return the associated uds-calibrated jet
      Jet udsCorrJet() const;
      /// return the associated gluon-calibrated jet
      Jet gluCorrJet() const;
      /// return the associated c-calibrated jet
      Jet cCorrJet() const;
      /// return the associated b-calibrated jet
      Jet bCorrJet() const;
      /// return the jet calibrated according to the MC flavour truth
      Jet mcFlavCorrJet() const;
      /// return the jet calibrated with weights assuming W decay
      Jet wCorrJet() const;

As I've said, we're looking for someone to update to the latest recommendations from the Jet/MET group.

Dileptons from Z' (Dimitri Bourilkov)

I have written a simple CompositeCandidate plotting kit that can handle most cases where Candidate Combiners are used. This can be done entirely in config files, as shown here by this code snippet:


module CompositeKitDemo = CompositeKit {
    # debug file name
    string outputTextName     =  'CompositeKitKit_output.txt'
    # parameters to enable and disable
    string disable            =  ''
    string enable             =  ''
    # parameters to ntuplize
    string ntuplize           =  'all'
    # are we running over a CSA07 soup?
    bool   isCSA07Soup        =  false 
    # composite candidate collection name
    InputTag compositeCandTag = zToEE
    string description        = 'Z to e+e-'
    # Kinematic distributions for daughters
    double pt1                = 0
    double pt2                = 200
    double m1                 = 0
    double m2                 = 200
    # resonance mass range
    double resonanceM1        = 500
    double resonanceM2        = 1500
}


  ### H->ZZ ###
  module zToEE = NamedCandViewShallowCloneCombiner {
    string decay = "selectedLayer1Electrons@+ selectedLayer1Electrons@-"
    string cut = "0.0 < mass < 20000.0"
    string name = "zToEE"
    vstring roles = { 'electron1', 'electron2' }
  }

  # Filter on number of  desired composite candidates
  module compositeFilter = CandViewCountFilter{
    InputTag src = zToEE
    uint32 minNumber = 1
  }

  ### Output ###
  path p = {
    patLayer0,
    patLayer1,
    zToEE,
    compositeFilter,
    CompositeKitDemo
  }

There is nothing different from constructing Z->e+e- or Z'->e+e- except for the cuts on the mass, which I have made very loose here. Note, however, that the resonance mass must be set by hand. This is a bit of a kludge for the moment, when we develop the upgraded histogramming utilities this will be done in the config file in a more natural way.

zprime_to_ee.gif

The full prescription is to run over the config file here:

cmsRun CompositeKitDemo_Zprime.cfg

H->ZZ->ll + jj (Youn Roh)

As in the previous example, this can be handled entirely in the config file, as shown by this code snippet:



module CompositeKitDemo = CompositeKit {
    # debug file name
    string outputTextName     =  'CompositeKitKit_output.txt'
    # parameters to enable and disable
    string disable            =  ''
    string enable             =  ''
    # parameters to ntuplize
    string ntuplize           =  'all'
    # are we running over a CSA07 soup?
    bool   isCSA07Soup        =  false 
    # composite candidate collection name
    InputTag compositeCandTag = hToZZ
    string description        = 'Higgs to Z + Z'
    # Kinematic distributions for daughters
    double pt1                = 0
    double pt2                = 200
    double m1                 = 0
    double m2                 = 200
    # resonance mass range
    double resonanceM1        = 0
    double resonanceM2        = 200
}


  ### H->ZZ->lljj ###
  module zToMuMu = NamedCandViewShallowCloneCombiner {
    string decay = "selectedLayer1Muons@+ selectedLayer1Muons@-"
    string cut = "0.0 < mass < 20000.0"
    string name = "zToMuMu"
    vstring roles = { 'muon1', 'muon2' }
  }
  module zToQQ = NamedCandViewShallowCloneCombiner {
    string decay = "selectedLayer1Jets selectedLayer1Jets"
    string cut = "0.0 < mass < 20000.0"
    string name = "zToQQ"
    vstring roles = { 'q1', 'q2' }
  }

  module hToZZ = NamedCandViewCombiner {
    string decay = "zToMuMu zToQQ"
    string cut = "0.0 < mass < 20000.0"
    string name = "hToZZ"
    vstring roles = { 'Z1', 'Z2' }
  }

  # Filter on number of  desired composite candidates
  module compositeFilter = CandViewCountFilter{
    InputTag src = hToZZ
    uint32 minNumber = 1
  }

  ### Output ###
  path p = {
    patLayer0,
    patLayer1,
    zToMuMu,
    hToZZ,
    compositeFilter,
    CompositeKitDemo
  }

Unfortunately I was unable to locate data here at Fermilab that has H->ZZ->lljj, so I was unable to make plots for Higgs. I did find some ZZ inclusive, so I was able to make plots for those, but Z-> jet + jet requires a bit more care in the selection than a quick and dirty example, so it isn't so pretty.

z1_mass_lljj.gif z2_mass_lljj.gif

The full prescription is to run over the config file here:

cmsRun CompositeKitDemo_lljj.cfg

Z+jet and Gamma+jet (Pt-Jet1 / Pt(Z->mu+mu) and Pt-Jet1/Pt-photon) (Anwar Bhatti)

In order to make plots about jet-balancing, it will be necessary to make a new "kit" since the objects in question are not basic kinematic variables. We are currently working to make these plots without writing any code, but the new features are not yet ready.

In the meantime, I have demonstrated the simplicity of creating a "kit" that will handle this called ZGammaJetBalance. It will plot the requested information, and can be modified in different ways to produce more, using the first two plots as an example.

A code snippet is:


class ZGammaJetBalanceKit : public StarterKit 
{
public:
  explicit ZGammaJetBalanceKit(const edm::ParameterSet&);
  virtual ~ZGammaJetBalanceKit();
    
protected:
  virtual void beginJob(const edm::EventSetup&) ;
  virtual void produce( edm::Event&, const edm::EventSetup&);
  virtual void endJob() ;

  // ----------member data ---------------------------

  // Input info
  edm::Handle<std::vector<reco::NamedCompositeCandidate> > zHandle_;
  edm::InputTag                                            zHandleName_;

  // The plots to make
  pat::PhysVarHisto  *  ptJet1OverPtZ_;
  pat::PhysVarHisto  *  ptJet1OverPtGamma_;
  pat::HistoComposite *   zHistos_;

  // Adding plots to ntuple
  std::vector<pat::PhysVarHisto *>  zgammaNtVars_;  
};

A code snippet is here:

  if ( jetHandle_->size() > 0 ) {

    // Get the first jet, sorted by pt
    pat::Jet const & ijet = jetHandle_->at(0);

    bool found = false;

    if ( muonHandle_->size() >= 2 || electronHandle_->size() >= 2 ) {
      // Get the composite candidates from upstream
      iEvent.getByLabel(zHandleName_,   zHandle_ );

      // Try Z + jet case first
      if ( zHandle_->size() > 0 ) {
   // Loop over the input Z collection
   vector<NamedCompositeCandidate>::const_iterator i = zHandle_->begin(),
     iend = zHandle_->end();
   for ( ; i != iend; ++i ) {
     // Fill Z histograms
     zHistos_->fill( *i );
     // Get the Pt of the Z
     double ptZ = i->pt();
     // Plot the Pt of the jet over Pt of the Z
     if ( ptZ > 0 ) {
       found = true;
       double f = ijet.pt() / ptZ;
       ptJet1OverPtZ_->fill( f );
     }
   }
    
      }
    }
...

The full analysis can be found here.

The plots are shown here:

ptj1_over_ptz.gif

ptj1_over_ptgamma.gif

The full prescription is to run over the config file here and here.

cmsRun ZGammaJetBalance.cfg
cmsRun ZGammaJetBalance_gamma.cfg

Accessing b-tagging information and flavor ID (Jochen Cammin)

Accessing the b-tagging information and the flavor ID is quite simple, the Layer 1 Jet object simply returns the TagInfo's in question:

      /// Get a tagInfo with the given name, or NULL if none is found. 
      /// You should omit the 'TagInfos' part from the label
      const reco::BaseTagInfo            * tagInfo(const std::string &label) const;
      /// Get a tagInfo with the given name and type or NULL if none is found. 
      /// If the label is empty or not specified, it returns the first tagInfo of that type (if any one exists)
      /// You should omit the 'TagInfos' part from the label
      const reco::TrackIPTagInfo         * tagInfoTrackIP(const std::string &label="") const;
      /// Get a tagInfo with the given name and type or NULL if none is found. 
      /// If the label is empty or not specified, it returns the first tagInfo of that type (if any one exists)
      /// You should omit the 'TagInfos' part from the label
      const reco::SoftLeptonTagInfo      * tagInfoSoftLepton(const std::string &label="") const;
      /// Get a tagInfo with the given name and type or NULL if none is found. 
      /// If the label is empty or not specified, it returns the first tagInfo of that type (if any one exists)
      /// You should omit the 'TagInfos' part from the label
      const reco::SecondaryVertexTagInfo * tagInfoSecondaryVertex(const std::string &label="") const;

Accessing the flavor is similarly trivial:

      /// return the matched generated parton
      const reco::Particle * genParton() const;
      /// return the matched generated jet
      const reco::GenJet * genJet() const;
      /// return the flavour of the parton underlying the jet
      int partonFlavour() const;

Accessing trigger information (Francisco Yumiceva)

Volker has a web page on the TWiki devoted to this:

https://twiki.cern.ch/twiki/bin/view/CMS/SWGuidePATLayer0#Trigger_matching

All PAT objects have a trigger match flag associated with them.

      /// trigger matches
      const std::vector<TriggerPrimitive> & triggerMatches() const;
      const std::vector<TriggerPrimitive> triggerMatchesByTrigger(const std::string & aTrig) const;
      const std::vector<TriggerPrimitive> triggerMatchesByHltFilter(const std::string & aFilt) const;
      const std::vector<TriggerPrimitive> triggerMatchesByL1Object(const std::string & aObj) const;

The TriggerPrimitive class is only available (currently) in the 1.6.x branch of the PAT.

Including default trigger matching is very easy, and by default they are turned on:

path p = {
    patLayer0_TriggerMatch , 
    patLayer1 ,
}

To turn them off, do

  replace allLayer1Electrons.addTrigMatch = false
  replace allLayer1Muons.addTrigMatch     = false
  replace allLayer1Jets.addTrigMatch      = false
  replace allLayer1METs.addTrigMatch      = false
  replace allLayer1Photons.addTrigMatch   = false
  replace allLayer1Taus.addTrigMatch      = false

path p = {
    patLayer0 , 
    patLayer1 ,
}

Then you can access the trigger primitives that are matched to the object in question.

If you want to match to a trigger that is not on by default, Volker has a web link here:

https://twiki.cern.ch/twiki/bin/view/CMS/SWGuidePATFAQs#How_can_I_implement_a_match_to_t

ttbar resonance production (Salvatore Rappoccio... hey, I'm allowed to request things too!)

I wanted to demonstrate how one would make a custom Layer 2 module and use the Starter kit to look at it.

As such, I have developed a BoostedTopProducer which will analyze boosted tops. The code can be obtained by

cvs co UserCode/srappocc/src/Analysis
mv UserCode/srappocc/src/Analysis .
scramv1 b
cd PhysicsTools/StarterKit/test
cmsRun CompositeKitDemo_CMS.BoostedTop.cfg

The output of this Layer 2 module is a collection of NameCompositeCandidates, which has the advantage of the full weight of the Physics Tools algorithms, selectors, and combiners. It also obviates any need for bookkeeping by the user. It is fully streamable, so can be written to and read from the event record trivially.

You can see the code to produce this here. This is a fairly involved example, so you can see how one would do various things with PAT objects in a complicated way.

This is run with the following snippet:


  ### PAT ###

  # PAT Layer 1 object production
  include "CMS.PhysicsTools/PatAlgos/data/patLayer0.cff"
  include "CMS.PhysicsTools/PatAlgos/data/patLayer1.cff"

  # turn off trigger matching
  replace allLayer1Electrons.addTrigMatch = false
  replace allLayer1Muons.addTrigMatch     = false
  replace allLayer1Jets.addTrigMatch      = false
  replace allLayer1METs.addTrigMatch      = false
  replace allLayer1Photons.addTrigMatch   = false
  replace allLayer1Taus.addTrigMatch      = false

  # PAT Layer 1 object selection
  replace minLayer1Jets.minNumber = 2
  replace countLayer1Leptons.minNumber = 1

  # TQAF Layer 2 for the ttbar semi-leptonic final state
  include "TopQuarkAnalysis/TopEventProducers/data/TtGenEvtProducer.cff"
  include "TopQuarkAnalysis/TopEventProducers/data/TtSemiEvtSolProducer.cfi"

  replace solutions.leptonFlavour = 'muon'

  # Boosted Top
  include "Analysis/BoostedTop/data/BoostedTopProducer.cfi"

  # Filter on number of boosted tops
  module boostedTopFilter = CandViewCountFilter{
    InputTag src = BoostedTopProducer
    uint32 minNumber = 1
  }

  # Starter kit
  include "CMS.PhysicsTools/StarterKit/test/CompositeKitDemo.cfi"

  replace CompositeKitDemo.compositeCandTag = BoostedTopProducer
  replace CompositeKitDemo.description = 'Boosted Top'
  replace CompositeKitDemo.isCSA07Soup = true


  service = TFileService {
    string fileName = "CompositeKit_CMS.BoostedTop.root"
  }


  ### Output ###
  path p = {
    #patPFTauTagging,                     # enable this on CSA07 AODSIMs (including the 1_5_2 ones)
    #pfRecoTauDiscriminationByIsolation,   # enable this instead on 1_6_10 RelVals
    csa07EventWeightProducer,            # enable this when you run on the CSA07 samples
#    patLayer0_TriggerMatch ,                           # use patLayer0_withoutPFTau if you can't get PFTaus working 
    patLayer0 ,
    patLayer1 ,
    makeGenEvt,
    solutions,
    BoostedTopProducer,
    boostedTopFilter,
    CompositeKitDemo
  }

A plot of the ttbar invariant mass for a few events of the Chowder sample is shown here. These are only W+jets events, because the ttbar events are further in on the file.

mttbar_boostedtop_chowder.gif

MC Truth matching (Dmitry Hits)

MC truth matching is totally trivial with the PAT:

   const Particle & mcmuon = muon->genLepton(); 
   const Particle & mcelectron = electron->genLepton();
   const Particle & genparton = jet->genParton();
   const GenJet & genjet = jet->genJet();

The details are that the MC matching is done at Layer 0, where an association map is produced for each layer 0 object and it's associated GenParticle.

2-d histograms (Dmitry Hits):

This feature is not yet available, but is coming soon.

Currently the best way to make multi-dimensional plots is to create an ntuple as demonstrated in the general tutorial, and then use that ntuple to make multiple dimensional plots.

lepton + jets with btagging (Cecilia Gerber):

This feature will be made simpler with the upcoming features of the Expression Histograms.

Currently there is no way to make "before and after" plots in a simple way, but this functionality will be added to the underlying histogram machinery in the coming weeks. There are two analyses (TtSemiEvtKit, and the Boosted Top example) which will demonstrate lepton+jets events, but the additional requirement of b-tagging to make all of these plots is not yet implemented, and will be in the coming weeks.

Addendum after Tutorial: I was able to find a solution within the current Starter Kit paradigm to make these plots. It involves making cuts using selectors and then adding multiple instances of the Starter Kit module after various stages of the passes.

So I implemented a filter to cut on the number of btags in the event:


struct BDiscriminatorSelector {
  BDiscriminatorSelector( std::string disc, double discCut ) : disc_(disc), discCut_(discCut) {}
  template<typename T>
  bool operator()( const T & t ) const { return t.bDiscriminator(disc_) > discCut_; }

private:
  std::string   disc_;
  double        discCut_;
};



namespace reco {
  namespace modules {
    
    template<>
      struct ParameterAdapter<BDiscriminatorSelector> {
      static BDiscriminatorSelector make( const edm::ParameterSet & cfg ) {
        return BDiscriminatorSelector( cfg.getParameter<std::string>( "disc" ),
                   cfg.getParameter<double> ( "discCut" ));
      }
    };

  }
}


namespace pat {

typedef ObjectCountFilter<
          std::vector<Jet>,
          BDiscriminatorSelector,
          MinNumberSelector
        > BDiscriminatorPatJetSelector;

}

This selector and accompanying filter will require minNumber of tags with a b discriminant disc above cut value discCut. The module config file looks like:


  ### ttbar semileptonic event kit ###

  # TtSemiEvtKit (layer 2) - calls LepJetMetKit by default (inherits from LepJetMetKit)
  module TtSemiEvtKitPreTag = TtSemiEvtKit {
   # These strings are currently not used
       string outputTextName    =  'TtSemiEvtKit_output.txt'
       string disable           =  'all'
       string enable            =  ''
       string ntuplize          =  'csa07Id,csa07Weight,lrJetCombProb,lrSignalEvtProb,kinFitProbChi2'
        bool isCSA07Soup         = true       
        InputTag EvtSolution     = solutions
  }
  module TtSemiEvtKitOneTag = TtSemiEvtKit {
   # These strings are currently not used
       string outputTextName    =  'TtSemiEvtKit_output.txt'
       string disable           =  'all'
       string enable            =  ''
       string ntuplize          =  'csa07Id,csa07Weight,lrJetCombProb,lrSignalEvtProb,kinFitProbChi2'
        bool isCSA07Soup         = true       
        InputTag EvtSolution     = solutions
  }
  module TtSemiEvtKitTwoTag = TtSemiEvtKit {
   # These strings are currently not used
       string outputTextName    =  'TtSemiEvtKit_output.txt'
       string disable           =  'all'
       string enable            =  ''
       string ntuplize          =  'csa07Id,csa07Weight,lrJetCombProb,lrSignalEvtProb,kinFitProbChi2'
        bool isCSA07Soup         = true       
        InputTag EvtSolution     = solutions
  }

  module TtSemiEvtKitThreeTag = TtSemiEvtKit {
   # These strings are currently not used
       string outputTextName    =  'TtSemiEvtKit_output.txt'
       string disable           =  'all'
       string enable            =  ''
       string ntuplize          =  'csa07Id,csa07Weight,lrJetCombProb,lrSignalEvtProb,kinFitProbChi2'
        bool isCSA07Soup         = true       
        InputTag EvtSolution     = solutions
  }

  module TtSemiEvtKitFourTag = TtSemiEvtKit {
   # These strings are currently not used
       string outputTextName    =  'TtSemiEvtKit_output.txt'
       string disable           =  'all'
       string enable            =  ''
       string ntuplize          =  'csa07Id,csa07Weight,lrJetCombProb,lrSignalEvtProb,kinFitProbChi2'
        bool isCSA07Soup         = true       
        InputTag EvtSolution     = solutions
  }

   ### Selection on number of btagged jets ###

  module btagFilter1 = BDiscriminatorPatJetSelector {
      InputTag src = selectedLayer1Jets
      string   disc   = 'trackCountingHighEffJetTags'
      double   discCut= 5.0
      uint32      minNumber= 1
  }
  module btagFilter2 = BDiscriminatorPatJetSelector {
      InputTag src = selectedLayer1Jets
      string   disc   = 'trackCountingHighEffJetTags'
      double   discCut= 5.0
      uint32      minNumber= 2
  }
  module btagFilter3 = BDiscriminatorPatJetSelector {
      InputTag src = selectedLayer1Jets
      string   disc   = 'trackCountingHighEffJetTags'
      double   discCut= 5.0
      uint32      minNumber= 3
  }
  module btagFilter4 = BDiscriminatorPatJetSelector {
      InputTag src = selectedLayer1Jets
      string   disc   = 'trackCountingHighEffJetTags'
      double   discCut= 5.0
      uint32      minNumber= 4
  }

  service = TFileService {
    string fileName = "TtSemiEvtKitHistos_1610.root"
  }

  



  ### Output ###
  path p = {
    #patPFTauTagging,                     # enable this on CSA07 AODSIMs (including the 1_5_2 ones)
    #pfRecoTauDiscriminationByIsolation,   # enable this instead on 1_6_10 RelVals
    csa07EventWeightProducer,            # enable this when you run on the CSA07 samples
    patLayer0_TriggerMatch ,                           # use patLayer0_withoutPFTau if you can't get PFTaus working 
    patLayer1 ,
    tqafLayer2_TtSemi,
    TtSemiEvtKitPreTag,
    btagFilter1,
    TtSemiEvtKitOneTag,
    btagFilter2,
    TtSemiEvtKitTwoTag,
    btagFilter3,
    TtSemiEvtKitThreeTag,
    btagFilter4,
    TtSemiEvtKitFourTag
  }

This will make five sets of plots, for

ntags >= 0, 1, 2, 3, 4
They are all stored in separate directories in the output file TtSemiEvtKitHistos_1610.root.

A full prescription is

cd $CMSSW_DIR/src
cvs co UserCode/srappocc/src/Analysis
mv UserCode/srappocc/src/Analysis .
scramv1 b
cd PhysicsTools/StarterKit/test
cmsRun TtSemiEvtKit_btag.cfg

Comments, questions and suggestions from the audience

-- SalvatoreRRappoccio - 14 May 2008

Topic attachments
I Attachment History Action Size Date Who Comment
GIFgif mttbar_boostedtop_chowder.gif r1 manage 6.9 K 2008-05-14 - 21:38 SalvatoreRRappoccio  
GIFgif ptj1_over_ptgamma.gif r1 manage 7.1 K 2008-05-14 - 21:36 SalvatoreRRappoccio  
GIFgif ptj1_over_ptz.gif r1 manage 6.4 K 2008-05-14 - 21:37 SalvatoreRRappoccio  
GIFgif z1_mass_lljj.gif r1 manage 7.2 K 2008-05-14 - 21:37 SalvatoreRRappoccio  
GIFgif z2_mass_lljj.gif r1 manage 7.0 K 2008-05-14 - 21:38 SalvatoreRRappoccio  
GIFgif zmass.gif r1 manage 7.5 K 2008-05-14 - 21:37 SalvatoreRRappoccio  
GIFgif zprime_to_ee.gif r1 manage 7.1 K 2008-05-14 - 21:38 SalvatoreRRappoccio  
Edit | Attach | Watch | Print version | History: r5 < r4 < r3 < r2 < r1 | Backlinks | Raw View | WYSIWYG | More topic actions
Topic revision: r5 - 2009-05-12 - 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