CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Types | Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes
pat::JetCorrFactorsProducer Class Reference

Produces a ValueMap between JetCorrFactors and the to the originating reco jets. More...

#include "PhysicsTools/PatAlgos/interface/JetCorrFactorsProducer.h"

Inheritance diagram for pat::JetCorrFactorsProducer:
edm::EDProducer edm::ProducerBase edm::ProductRegistryHelper

Public Types

typedef std::map
< JetCorrFactors::Flavor,
std::vector< std::string > > 
FlavorCorrLevelMap
 map of correction levels to different flavors More...
 
typedef edm::ValueMap
< pat::JetCorrFactors
JetCorrFactorsMap
 value map for JetCorrFactors (to be written into the event) More...
 
- Public Types inherited from edm::EDProducer
typedef EDProducer ModuleType
 
typedef WorkerT< EDProducerWorkerType
 
- Public Types inherited from edm::ProducerBase
typedef
ProductRegistryHelper::TypeLabelList 
TypeLabelList
 

Public Member Functions

 JetCorrFactorsProducer (const edm::ParameterSet &cfg)
 default constructor More...
 
virtual void produce (edm::Event &event, const edm::EventSetup &setup)
 everything that needs to be done per event More...
 
 ~JetCorrFactorsProducer ()
 default destructor More...
 
- Public Member Functions inherited from edm::EDProducer
 EDProducer ()
 
virtual ~EDProducer ()
 
- Public Member Functions inherited from edm::ProducerBase
 ProducerBase ()
 
void registerProducts (ProducerBase *, ProductRegistry *, ModuleDescription const &)
 
boost::function< void(const
BranchDescription &)> 
registrationCallback () const
 used by the fwk to register list of products More...
 
virtual ~ProducerBase ()
 

Static Public Member Functions

static void fillDescriptions (edm::ConfigurationDescriptions &descriptions)
 description of configuration file parameters More...
 
- Static Public Member Functions inherited from edm::EDProducer
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 

Private Member Functions

float evaluate (edm::View< reco::Jet >::const_iterator &jet, boost::shared_ptr< FactorizedJetCorrector > &corrector, int level)
 evaluate jet correction factor up to a given level More...
 
std::vector< std::string > expand (const std::vector< std::string > &levels, const JetCorrFactors::Flavor &flavor)
 
bool flavorDependent () const
 return true if the jec levels contain at least one flavor dependent correction level More...
 
int numberOf (const edm::Handle< std::vector< reco::Vertex > > &primaryVertices)
 determines the number of valid primary vertices for the standard L1Offset correction of JetMET More...
 
std::vector
< JetCorrectorParameters
params (const JetCorrectorParametersCollection &parameters, const std::vector< std::string > &levels) const
 return the jec parameters as input to the FactorizedJetCorrector for different flavors More...
 
std::string payload ()
 map jet algorithm to payload in DB More...
 

Private Attributes

bool emf_
 use electromagnetic fraction for jet energy corrections or not (will only have an effect for jets CaloJets) More...
 
std::string label_
 label of jec factors More...
 
FlavorCorrLevelMap levels_
 
std::string payload_
 label of payload More...
 
edm::InputTag primaryVertices_
 label for L1Offset primaryVertex collection More...
 
edm::InputTag rho_
 label for L1FastJet energy density parameter rho More...
 
edm::InputTag src_
 input jet collection More...
 
std::string type_
 type of flavor dependent JEC factors (only 'J' and 'T' are allowed) More...
 
bool useNPV_
 use the NPV and rho with the JEC? (used for L1Offset/L1FastJet and L1FastJet, resp.) More...
 
bool useRho_
 

Additional Inherited Members

- Protected Member Functions inherited from edm::EDProducer
CurrentProcessingContext const * currentContext () const
 
- Protected Member Functions inherited from edm::ProducerBase
template<class TProducer , class TMethod >
void callWhenNewProductsRegistered (TProducer *iProd, TMethod iMethod)
 

Detailed Description

Produces a ValueMap between JetCorrFactors and the to the originating reco jets.

The JetCorrFactorsProducer produces a set of correction factors, defined in the class pat::JetCorrFactors. This vector is linked to the originating reco jets through an edm::ValueMap. The initializing parameters of the module can be found in the recoLayer1/jetCorrFactors_cfi.py of the PatAlgos package. In the standard PAT workflow the module has to be run before the creation of the pat::Jet. The edm::ValueMap will then be embedded into the pat::Jet.

Jets corrected up to a given correction level can then be accessed via the pat::Jet member function correctedJet. For more details have a look into the class description of the pat::Jet.

ATTENTION: available options for flavor corrections are L5Flavor_gJ L7Parton_gJ gluon from dijets L5Flavor_qJ/_qT L7Parton_qJ/_qT quark from dijets/top L5Flavor_cJ/_cT L7Parton_cJ/_cT charm from dijets/top L5Flavor_bJ/_bT L7Parton_bJ/_bT beauty from dijets/top L7Parton_jJ/_tT mixture from dijets/top

where mixture refers to the flavor mixture as determined from the MC sample the flavor dependent corrections have been derived from. 'J' and 'T' stand for a typical dijet (ttbar) sample.

L1Offset corrections require the collection of offlinePrimaryVertices, which are supposed to be added as an additional optional parameter primaryVertices in the jetCorrFactors_cfi.py file.

L1FastJet corrections, which are an alternative to the standard L1Offset correction as recommended by the JetMET PAG the energy density parameter rho is supposed to be added as an additional optional parameter rho in the jetCorrFactors_cfi.py file.

NOTE: the mixed mode (mc input mixture from dijets/ttbar) only exists for parton level corrections. jJ and tT are not covered in this implementation of the JetCorrFactorsProducer there are no gluon corrections available from the top sample on the L7Parton level.

Definition at line 59 of file JetCorrFactorsProducer.h.

Member Typedef Documentation

typedef std::map<JetCorrFactors::Flavor, std::vector<std::string> > pat::JetCorrFactorsProducer::FlavorCorrLevelMap

map of correction levels to different flavors

Definition at line 64 of file JetCorrFactorsProducer.h.

value map for JetCorrFactors (to be written into the event)

Definition at line 62 of file JetCorrFactorsProducer.h.

Constructor & Destructor Documentation

JetCorrFactorsProducer::JetCorrFactorsProducer ( const edm::ParameterSet cfg)
explicit

default constructor

Definition at line 22 of file JetCorrFactorsProducer.cc.

References pat::JetCorrFactors::BOTTOM, pat::JetCorrFactors::CHARM, edm::hlt::Exception, edm::ParameterSet::existsAs(), expand(), spr::find(), edm::ParameterSet::getParameter(), pat::JetCorrFactors::GLUON, levels_, python.rootplot.argparse::message, pat::JetCorrFactors::NONE, primaryVertices_, rho_, pat::JetCorrFactors::UDS, useNPV_, and useRho_.

22  :
23  emf_(cfg.getParameter<bool>( "emf" )),
24  src_(cfg.getParameter<edm::InputTag>( "src" )),
25  type_ (cfg.getParameter<std::string>("flavorType")),
26  label_(cfg.getParameter<std::string>( "@module_label" )),
27  payload_( cfg.getParameter<std::string>("payload") ),
28  useNPV_(cfg.getParameter<bool>("useNPV")),
29  useRho_(cfg.getParameter<bool>("useRho"))
30 {
31 
32  std::vector<std::string> levels = cfg.getParameter<std::vector<std::string> >("levels");
33  // fill the std::map for levels_, which might be flavor dependent or not;
34  // flavor dependency is determined from the fact whether the std::string
35  // L5Flavor or L7Parton can be found in levels; if flavor dependent four
36  // vectors of strings will be filled into the map corresponding to GLUON,
37  // UDS, CHARM and BOTTOM (according to JetCorrFactors::Flavor), 'L5Flavor'
38  // and 'L7Parton' will be expanded accordingly; if not levels_ is filled
39  // with only one vector of strings according to NONE. This vector will be
40  // equivalent to the original vector of strings.
41  if(std::find(levels.begin(), levels.end(), "L5Flavor")!=levels.end() || std::find(levels.begin(), levels.end(), "L7Parton")!=levels.end()){
46  }
47  else{
48  levels_[JetCorrFactors::NONE ] = levels;
49  }
50  // if the std::string L1Offset can be found in levels an additional para-
51  // meter primaryVertices is needed, which should pass on the offline pri-
52  // mary vertex collection. The size of this collection is needed for the
53  // L1Offset correction.
54  if(useNPV_){
55  if(cfg.existsAs<edm::InputTag>("primaryVertices")){
56  primaryVertices_=cfg.getParameter<edm::InputTag>("primaryVertices");
57  }
58  else{
59  throw cms::Exception("No primaryVertices specified")
60  << "The configured correction levels contain an L1Offset or L1FastJet correction, \n"
61  << "which requires the number of offlinePrimaryVertices. Please specify this col- \n"
62  << "lection as additional optional parameter primaryVertices in the jetCorrFactors\n"
63  << "module. \n";
64  }
65  }
66  // if the std::string L1FastJet can be found in levels an additional
67  // parameter rho is needed, which should pass on the energy density
68  // parameter for the corresponding jet collection.
69  if(useRho_){
70  if(std::find(levels.begin(), levels.end(), "L1FastJet")!=levels.end()){
71  if(cfg.existsAs<edm::InputTag>("rho")){
72  rho_=cfg.getParameter<edm::InputTag>("rho");
73  }
74  else{
75  throw cms::Exception("No parameter rho specified")
76  << "The configured correction levels contain a L1FastJet correction, which re- \n"
77  << "quires the energy density parameter rho. Please specify this collection as \n"
78  << "additional optional parameter rho in the jetCorrFactors module. \n";
79  }
80  }
81  else{
82  edm::LogWarning message( "Parameter rho not used" );
83  message << "Module is configured to use the parameter rho, but but rho is only used \n"
84  << "for L1FastJet corrections at the moment. The configuration of levels \n"
85  << "does not contain L1FastJet corrections though, so rho will not be used \n"
86  << "throughout this module. \n";
87  }
88  }
89  produces<JetCorrFactorsMap>();
90 }
T getParameter(std::string const &) const
bool existsAs(std::string const &parameterName, bool trackiness=true) const
checks if a parameter exists as a given type
Definition: ParameterSet.h:187
std::string label_
label of jec factors
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:7
bool emf_
use electromagnetic fraction for jet energy corrections or not (will only have an effect for jets Cal...
std::string payload_
label of payload
std::vector< std::string > expand(const std::vector< std::string > &levels, const JetCorrFactors::Flavor &flavor)
edm::InputTag src_
input jet collection
edm::InputTag primaryVertices_
label for L1Offset primaryVertex collection
edm::InputTag rho_
label for L1FastJet energy density parameter rho
bool useNPV_
use the NPV and rho with the JEC? (used for L1Offset/L1FastJet and L1FastJet, resp.)
std::string type_
type of flavor dependent JEC factors (only &#39;J&#39; and &#39;T&#39; are allowed)
pat::JetCorrFactorsProducer::~JetCorrFactorsProducer ( )
inline

default destructor

Definition at line 70 of file JetCorrFactorsProducer.h.

70 {};

Member Function Documentation

float JetCorrFactorsProducer::evaluate ( edm::View< reco::Jet >::const_iterator &  jet,
boost::shared_ptr< FactorizedJetCorrector > &  corrector,
int  level 
)
private

evaluate jet correction factor up to a given level

Definition at line 130 of file JetCorrFactorsProducer.cc.

References emf_, reco::LeafCandidate::energy(), reco::LeafCandidate::eta(), reco::JPTJet::getCaloJetRef(), metsig::jet, testEve_cfg::level, p4, reco::LeafCandidate::phi(), and reco::LeafCandidate::pt().

Referenced by produce().

131 {
132  // add parameters for JPT corrections
133  const reco::JPTJet* jpt = dynamic_cast<reco::JPTJet const *>( &*jet );
134  if( jpt ){
135  TLorentzVector p4; p4.SetPtEtaPhiE(jpt->getCaloJetRef()->pt(), jpt->getCaloJetRef()->eta(), jpt->getCaloJetRef()->phi(), jpt->getCaloJetRef()->energy());
136  corrector->setJPTrawP4(p4);
137  }
138  corrector->setJetEta(jet->eta()); corrector->setJetPt(jet->pt()); corrector->setJetE(jet->energy());
139  if( emf_ && dynamic_cast<const reco::CaloJet*>(&*jet)){
140  corrector->setJetEMF(dynamic_cast<const reco::CaloJet*>(&*jet)->emEnergyFraction());
141  }
142  return corrector->getSubCorrections()[level];
143 }
const edm::RefToBase< reco::Jet > & getCaloJetRef() const
Definition: JPTJet.h:130
bool emf_
use electromagnetic fraction for jet energy corrections or not (will only have an effect for jets Cal...
virtual double eta() const
momentum pseudorapidity
virtual double energy() const
energy
double p4[4]
Definition: TauolaWrapper.h:92
Jets made from CaloJets corrected for ZSP and tracks.
Definition: JPTJet.h:29
virtual double pt() const
transverse momentum
tuple level
Definition: testEve_cfg.py:34
virtual double phi() const
momentum azimuthal angle
std::vector< std::string > JetCorrFactorsProducer::expand ( const std::vector< std::string > &  levels,
const JetCorrFactors::Flavor flavor 
)
private

return an expanded version of correction levels for different flavors; the result should be of type ['L2Relative', 'L3Absolute', 'L5FLavor_gJ', 'L7Parton_gJ']; L7Parton_gT will result in an empty string as this correction level is not available

Definition at line 93 of file JetCorrFactorsProducer.cc.

References pat::JetCorrFactors::BOTTOM, pat::JetCorrFactors::CHARM, pat::JetCorrFactors::GLUON, testEve_cfg::level, python.rootplot.argparse::message, type_, and pat::JetCorrFactors::UDS.

Referenced by JetCorrFactorsProducer().

94 {
95  std::vector<std::string> expand;
96  for(std::vector<std::string>::const_iterator level=levels.begin(); level!=levels.end(); ++level){
97  if((*level)=="L5Flavor" || (*level)=="L7Parton"){
98  if(flavor==JetCorrFactors::GLUON ){
99  if(*level=="L7Parton" && type_=="T"){
100  edm::LogWarning message( "L7Parton::GLUON not available" );
101  message << "Jet energy corrections requested for level: L7Parton and type: 'T'. \n"
102  << "For this combination there is no GLUON correction available. The \n"
103  << "correction for this flavor type will be taken from 'J'.";
104  }
105  expand.push_back(std::string(*level).append("_").append("g").append("J"));
106  }
107  if(flavor==JetCorrFactors::UDS ) expand.push_back(std::string(*level).append("_").append("q").append(type_));
108  if(flavor==JetCorrFactors::CHARM ) expand.push_back(std::string(*level).append("_").append("c").append(type_));
109  if(flavor==JetCorrFactors::BOTTOM) expand.push_back(std::string(*level).append("_").append("b").append(type_));
110  }
111  else{
112  expand.push_back(*level);
113  }
114  }
115  return expand;
116 }
std::vector< std::string > expand(const std::vector< std::string > &levels, const JetCorrFactors::Flavor &flavor)
tuple level
Definition: testEve_cfg.py:34
std::string type_
type of flavor dependent JEC factors (only &#39;J&#39; and &#39;T&#39; are allowed)
void JetCorrFactorsProducer::fillDescriptions ( edm::ConfigurationDescriptions descriptions)
static

description of configuration file parameters

Definition at line 262 of file JetCorrFactorsProducer.cc.

References edm::ConfigurationDescriptions::add(), and edm::ParameterSetDescription::add().

263 {
265  iDesc.add<bool>("emf", false);
266  iDesc.add<std::string>("flavorType", "J");
267  iDesc.add<edm::InputTag>("src", edm::InputTag("ak5CaloJets"));
268  iDesc.add<std::string>("payload", "AK5Calo");
269  iDesc.add<bool>("useNPV", true);
270  iDesc.add<edm::InputTag>("primaryVertices", edm::InputTag("offlinePrimaryVertices"));
271  iDesc.add<bool>("useRho", true);
272  iDesc.add<edm::InputTag>("rho", edm::InputTag("kt6PFJets", "rho"));
273 
274  std::vector<std::string> levels;
275  levels.push_back(std::string("L1Offset" ));
276  levels.push_back(std::string("L2Relative"));
277  levels.push_back(std::string("L3Absolute"));
278  levels.push_back(std::string("L5Flavor" ));
279  levels.push_back(std::string("L7Parton" ));
280  iDesc.add<std::vector<std::string> >("levels", levels);
281  descriptions.add("JetCorrFactorsProducer", iDesc);
282 }
ParameterDescriptionBase * add(U const &iLabel, T const &value)
void add(std::string const &label, ParameterSetDescription const &psetDescription)
bool pat::JetCorrFactorsProducer::flavorDependent ( ) const
inlineprivate

return true if the jec levels contain at least one flavor dependent correction level

Definition at line 78 of file JetCorrFactorsProducer.h.

References levels_.

Referenced by produce().

78 { return (levels_.size()>1); };
int pat::JetCorrFactorsProducer::numberOf ( const edm::Handle< std::vector< reco::Vertex > > &  primaryVertices)
inlineprivate

determines the number of valid primary vertices for the standard L1Offset correction of JetMET

Definition at line 124 of file JetCorrFactorsProducer.h.

Referenced by produce().

125  {
126  int npv=0;
127  for(std::vector<reco::Vertex>::const_iterator pv=primaryVertices->begin(); pv!=primaryVertices->end(); ++pv){
128  if(pv->ndof()>=4) ++npv;
129  }
130  return npv;
131  }
std::vector< JetCorrectorParameters > JetCorrFactorsProducer::params ( const JetCorrectorParametersCollection parameters,
const std::vector< std::string > &  levels 
) const
private

return the jec parameters as input to the FactorizedJetCorrector for different flavors

Definition at line 119 of file JetCorrFactorsProducer.cc.

References testEve_cfg::level, and JetCorrectorParametersCollection::push_back().

Referenced by produce().

120 {
121  std::vector<JetCorrectorParameters> params;
122  for(std::vector<std::string>::const_iterator level=levels.begin(); level!=levels.end(); ++level){
123  const JetCorrectorParameters& ip = parameters[*level]; //ip.printScreen();
124  params.push_back(ip);
125  }
126  return params;
127 }
void push_back(key_type i, value_type const &j, label_type const &flav="")
std::vector< JetCorrectorParameters > params(const JetCorrectorParametersCollection &parameters, const std::vector< std::string > &levels) const
return the jec parameters as input to the FactorizedJetCorrector for different flavors ...
tuple level
Definition: testEve_cfg.py:34
std::string JetCorrFactorsProducer::payload ( )
private

map jet algorithm to payload in DB

Definition at line 146 of file JetCorrFactorsProducer.cc.

References payload_.

Referenced by produce().

147 {
148  return payload_;
149 }
std::string payload_
label of payload
void JetCorrFactorsProducer::produce ( edm::Event event,
const edm::EventSetup setup 
)
virtual

everything that needs to be done per event

Implements edm::EDProducer.

Definition at line 152 of file JetCorrFactorsProducer.cc.

References evaluate(), edm::hlt::Exception, edm::helper::Filler< Map >::fill(), flavorDependent(), edm::EventSetup::get(), edm::Event::getByLabel(), edm::helper::Filler< Map >::insert(), patTestJEC_cfi::jec, metsig::jet, fwrapper::jets, edm::InputTag::label(), label_, levels_, numberOf(), Parameters::parameters, params(), payload(), primaryVertices_, rho, rho_, and src_.

153 {
154  // get jet collection from the event
156  event.getByLabel(src_, jets);
157 
158  // get primary vertices for L1Offset correction level if needed
159  edm::Handle<std::vector<reco::Vertex> > primaryVertices;
160  if(!primaryVertices_.label().empty()) event.getByLabel(primaryVertices_, primaryVertices);
161 
162  // get parameter rho for L1FastJet correction level if needed
163 
165  if(!rho_.label().empty()) event.getByLabel(rho_, rho);
166 
167  // retreive parameters from the DB
169  setup.get<JetCorrectionsRecord>().get(payload(), parameters);
170 
171  // initialize jet correctors
172  std::map<JetCorrFactors::Flavor, boost::shared_ptr<FactorizedJetCorrector> > corrector;
173  for(FlavorCorrLevelMap::const_iterator flavor=levels_.begin(); flavor!=levels_.end(); ++flavor){
174  corrector[flavor->first] = boost::shared_ptr<FactorizedJetCorrector>( new FactorizedJetCorrector(params(*parameters, flavor->second)) );
175  }
176 
177  // fill the jetCorrFactors
178  std::vector<JetCorrFactors> jcfs;
179  for(edm::View<reco::Jet>::const_iterator jet = jets->begin(); jet!=jets->end(); ++jet){
180  // the JetCorrFactors::CorrectionFactor is a std::pair<std::string, std::vector<float> >
181  // the string corresponds to the label of the correction level, the vector contains four
182  // floats if flavor dependent and one float else. Per construction jet energy corrections
183  // will be flavor independent up to the first flavor dependent correction and flavor de-
184  // pendent afterwards. The first correction level is predefined with label 'Uncorrected'.
185  // Per definition it is flavor independent. The correction factor is 1.
186  std::vector<JetCorrFactors::CorrectionFactor> jec;
187  jec.push_back(std::make_pair<std::string, std::vector<float> >(std::string("Uncorrected"), std::vector<float>(1, 1)));
188 
189  // pick the first element in the map (which could be the only one) and loop all jec
190  // levels listed for that element. If this is not the only element all jec levels, which
191  // are flavor independent will give the same correction factors until the first flavor
192  // dependent correction level is reached. So the first element is still a good choice.
193  FlavorCorrLevelMap::const_iterator corrLevel=levels_.begin();
194  if(corrLevel==levels_.end()){
195  throw cms::Exception("No JECFactors") << "You request to create a jetCorrFactors object with no JEC Levels indicated. \n"
196  << "This makes no sense, either you should correct this or drop the module from \n"
197  << "the sequence.";
198  }
199  for(unsigned int idx=0; idx<corrLevel->second.size(); ++idx){
200  bool flavorDependent=false;
201  std::vector<float> factors;
202  if(flavorDependent ||
203  corrLevel->second[idx].find("L5Flavor")!=std::string::npos ||
204  corrLevel->second[idx].find("L7Parton")!=std::string::npos){
205  flavorDependent=true;
206  // after the first encounter all subsequent correction levels are flavor dependent
207  for(FlavorCorrLevelMap::const_iterator flavor=corrLevel; flavor!=levels_.end(); ++flavor){
208  if(!primaryVertices_.label().empty()){
209  // if primaryVertices_ has a value the number of primary vertices needs to be
210  // specified
211  corrector.find(flavor->first)->second->setNPV(numberOf(primaryVertices));
212  }
213  if(!rho_.label().empty()){
214  // if rho_ has a value the energy density parameter rho and the jet area need
215  // to be specified
216  corrector.find(flavor->first)->second->setRho(*rho);
217  corrector.find(flavor->first)->second->setJetA(jet->jetArea());
218  }
219  factors.push_back(evaluate(jet, corrector.find(flavor->first)->second, idx));
220  }
221  }
222  else{
223  if(!primaryVertices_.label().empty()){
224  // if primaryVertices_ has a value the number of primary vertices needs to be
225  // specified
226  corrector.find(corrLevel->first)->second->setNPV(numberOf(primaryVertices));
227  }
228  if(!rho_.label().empty()){
229  // if rho_ has a value the energy density parameter rho and the jet area need
230  // to be specified
231  corrector.find(corrLevel->first)->second->setRho(*rho);
232  corrector.find(corrLevel->first)->second->setJetA(jet->jetArea());
233  }
234  factors.push_back(evaluate(jet, corrector.find(corrLevel->first)->second, idx));
235  }
236  // push back the set of JetCorrFactors: the first entry corresponds to the label
237  // of the correction level, which is taken from the first element in levels_. For
238  // L5Flavor and L7Parton the part including the first '_' indicating the flavor
239  // of the first element in levels_ is chopped of from the label to avoid confusion
240  // of the correction levels. The second parameter corresponds to the set of jec
241  // factors, which might be flavor dependent or not. In the default configuration
242  // the CorrectionFactor will look like this: 'Uncorrected': 1 ; 'L2Relative': x ;
243  // 'L3Absolute': x ; 'L5Flavor': v, x, y, z ; 'L7Parton': v, x, y, z
244  jec.push_back(std::make_pair((corrLevel->second[idx]).substr(0, (corrLevel->second[idx]).find("_")), factors));
245  }
246  // create the actual object with the scale factors we want the valuemap to refer to
247  // label_ corresponds to the label of the module instance
248  JetCorrFactors corrFactors(label_, jec);
249  jcfs.push_back(corrFactors);
250  }
251  // build the value map
252  std::auto_ptr<JetCorrFactorsMap> jetCorrsMap(new JetCorrFactorsMap());
253  JetCorrFactorsMap::Filler filler(*jetCorrsMap);
254  // jets and jetCorrs have their indices aligned by construction
255  filler.insert(jets, jcfs.begin(), jcfs.end());
256  filler.fill(); // do the actual filling
257  // put our produced stuff in the event
258  event.put(jetCorrsMap);
259 }
dictionary parameters
Definition: Parameters.py:2
boost::indirect_iterator< typename seq_t::const_iterator > const_iterator
Definition: View.h:81
std::string payload()
map jet algorithm to payload in DB
std::string label_
label of jec factors
Definition: DDAxes.h:10
bool flavorDependent() const
return true if the jec levels contain at least one flavor dependent correction level ...
int numberOf(const edm::Handle< std::vector< reco::Vertex > > &primaryVertices)
determines the number of valid primary vertices for the standard L1Offset correction of JetMET ...
Class for the storage of jet correction factors.
vector< PseudoJet > jets
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:356
helper::Filler< ValueMap< T > > Filler
Definition: ValueMap.h:157
edm::ValueMap< pat::JetCorrFactors > JetCorrFactorsMap
value map for JetCorrFactors (to be written into the event)
const T & get() const
Definition: EventSetup.h:55
std::vector< JetCorrectorParameters > params(const JetCorrectorParametersCollection &parameters, const std::vector< std::string > &levels) const
return the jec parameters as input to the FactorizedJetCorrector for different flavors ...
std::string const & label() const
Definition: InputTag.h:25
edm::InputTag src_
input jet collection
float evaluate(edm::View< reco::Jet >::const_iterator &jet, boost::shared_ptr< FactorizedJetCorrector > &corrector, int level)
evaluate jet correction factor up to a given level
edm::InputTag primaryVertices_
label for L1Offset primaryVertex collection
edm::InputTag rho_
label for L1FastJet energy density parameter rho

Member Data Documentation

bool pat::JetCorrFactorsProducer::emf_
private

use electromagnetic fraction for jet energy corrections or not (will only have an effect for jets CaloJets)

Definition at line 94 of file JetCorrFactorsProducer.h.

Referenced by evaluate().

std::string pat::JetCorrFactorsProducer::label_
private
FlavorCorrLevelMap pat::JetCorrFactorsProducer::levels_
private

jec levels for different flavors. In the default configuration this map would look like this: GLUON : 'L2Relative', 'L3Absolute', 'L5FLavor_jg', L7Parton_jg' UDS : 'L2Relative', 'L3Absolute', 'L5FLavor_jq', L7Parton_jq' CHARM : 'L2Relative', 'L3Absolute', 'L5FLavor_jc', L7Parton_jc' BOTTOM : 'L2Relative', 'L3Absolute', 'L5FLavor_jb', L7Parton_jb' or just like this: NONE : 'L2Relative', 'L3Absolute', 'L2L3Residual' per definition the vectors for all elements in this map should have the same size

Definition at line 120 of file JetCorrFactorsProducer.h.

Referenced by flavorDependent(), JetCorrFactorsProducer(), and produce().

std::string pat::JetCorrFactorsProducer::payload_
private

label of payload

Definition at line 102 of file JetCorrFactorsProducer.h.

Referenced by payload().

edm::InputTag pat::JetCorrFactorsProducer::primaryVertices_
private

label for L1Offset primaryVertex collection

Definition at line 104 of file JetCorrFactorsProducer.h.

Referenced by JetCorrFactorsProducer(), and produce().

edm::InputTag pat::JetCorrFactorsProducer::rho_
private

label for L1FastJet energy density parameter rho

Definition at line 106 of file JetCorrFactorsProducer.h.

Referenced by JetCorrFactorsProducer(), and produce().

edm::InputTag pat::JetCorrFactorsProducer::src_
private

input jet collection

Definition at line 96 of file JetCorrFactorsProducer.h.

Referenced by produce().

std::string pat::JetCorrFactorsProducer::type_
private
bool pat::JetCorrFactorsProducer::useNPV_
private

use the NPV and rho with the JEC? (used for L1Offset/L1FastJet and L1FastJet, resp.)

Definition at line 108 of file JetCorrFactorsProducer.h.

Referenced by JetCorrFactorsProducer().

bool pat::JetCorrFactorsProducer::useRho_
private

Definition at line 109 of file JetCorrFactorsProducer.h.

Referenced by JetCorrFactorsProducer().