CMS 3D CMS Logo

List of all members | Public Types | Public Member Functions | Private Member Functions | Private Attributes
JetDeltaRValueMapProducer< T, C > Class Template Reference
Inheritance diagram for JetDeltaRValueMapProducer< T, C >:
edm::stream::EDProducer<>

Public Types

typedef edm::ValueMap< float > JetValueMap
 
- Public Types inherited from edm::stream::EDProducer<>
typedef CacheContexts< T... > CacheTypes
 
typedef CacheTypes::GlobalCache GlobalCache
 
typedef AbilityChecker< T... > HasAbility
 
typedef CacheTypes::LuminosityBlockCache LuminosityBlockCache
 
typedef LuminosityBlockContextT< LuminosityBlockCache, RunCache, GlobalCacheLuminosityBlockContext
 
typedef CacheTypes::LuminosityBlockSummaryCache LuminosityBlockSummaryCache
 
typedef CacheTypes::RunCache RunCache
 
typedef RunContextT< RunCache, GlobalCacheRunContext
 
typedef CacheTypes::RunSummaryCache RunSummaryCache
 

Public Member Functions

 JetDeltaRValueMapProducer (edm::ParameterSet const &params)
 
 ~JetDeltaRValueMapProducer () override
 
- Public Member Functions inherited from edm::stream::EDProducer<>
 EDProducer ()=default
 
bool hasAbilityToProduceInLumis () const final
 
bool hasAbilityToProduceInRuns () const final
 

Private Member Functions

void produce (edm::Event &iEvent, const edm::EventSetup &iSetup) override
 

Private Attributes

const double distMax_
 
std::map< std::string, std::unique_ptr< const StringObjectFunction< C > > > evaluationMap_
 
const bool lazyParser_
 
const edm::EDGetTokenT< typename edm::View< C > > matchedToken_
 
bool multiValue_
 
const edm::EDGetTokenT< typename edm::View< T > > srcToken_
 
const std::string value_
 
const std::vector< std::string > valueLabels_
 
const std::vector< std::string > values_
 

Detailed Description

template<class T, class C = T>
class JetDeltaRValueMapProducer< T, C >

Definition at line 32 of file JetDeltaRValueMapProducer.cc.

Member Typedef Documentation

template<class T , class C = T>
typedef edm::ValueMap<float> JetDeltaRValueMapProducer< T, C >::JetValueMap

Definition at line 36 of file JetDeltaRValueMapProducer.cc.

Constructor & Destructor Documentation

template<class T , class C = T>
JetDeltaRValueMapProducer< T, C >::JetDeltaRValueMapProducer ( edm::ParameterSet const &  params)
inline

Definition at line 38 of file JetDeltaRValueMapProducer.cc.

References JetDeltaRValueMapProducer< T, C >::evaluationMap_, mps_fire::i, JetDeltaRValueMapProducer< T, C >::lazyParser_, JetDeltaRValueMapProducer< T, C >::multiValue_, JetDeltaRValueMapProducer< T, C >::value_, JetDeltaRValueMapProducer< T, C >::valueLabels_, and JetDeltaRValueMapProducer< T, C >::values_.

38  :
39  srcToken_( consumes< typename edm::View<T> >( params.getParameter<edm::InputTag>("src") ) ),
40  matchedToken_( consumes< typename edm::View<C> >( params.getParameter<edm::InputTag>( "matched" ) ) ),
41  distMax_( params.getParameter<double>( "distMax" ) ),
42  value_( params.existsAs<std::string>("value") ? params.getParameter<std::string>("value") : "" ),
43  values_( params.existsAs<std::vector<std::string> >("values") ? params.getParameter<std::vector<std::string> >("values") : std::vector<std::string>() ),
44  valueLabels_( params.existsAs<std::vector<std::string> >("valueLabels") ? params.getParameter<std::vector<std::string> >("valueLabels") : std::vector<std::string>() ),
45  lazyParser_( params.existsAs<bool>("lazyParser") ? params.getParameter<bool>("lazyParser") : false ),
46  multiValue_(false)
47  {
48  if( value_!="" )
49  {
50  evaluationMap_.insert( std::make_pair( value_, std::unique_ptr<StringObjectFunction<C> >( new StringObjectFunction<C>( value_, lazyParser_ ) ) ) );
51  produces< JetValueMap >();
52  }
53 
54  if( !valueLabels_.empty() || !values_.empty() )
55  {
56  if( valueLabels_.size()==values_.size() )
57  {
58  multiValue_ = true;
59  for( size_t i=0; i<valueLabels_.size(); ++i)
60  {
61  evaluationMap_.insert( std::make_pair( valueLabels_[i], std::unique_ptr<StringObjectFunction<C> >( new StringObjectFunction<C>( values_[i], lazyParser_ ) ) ) );
62  produces< JetValueMap >(valueLabels_[i]);
63  }
64  }
65  else
66  edm::LogWarning("ValueLabelMismatch") << "The number of value labels does not match the number of values. Values will not be evaluated.";
67  }
68  }
const std::vector< std::string > valueLabels_
std::map< std::string, std::unique_ptr< const StringObjectFunction< C > > > evaluationMap_
const edm::EDGetTokenT< typename edm::View< C > > matchedToken_
const edm::EDGetTokenT< typename edm::View< T > > srcToken_
const std::vector< std::string > values_
template<class T , class C = T>
JetDeltaRValueMapProducer< T, C >::~JetDeltaRValueMapProducer ( )
inlineoverride

Definition at line 70 of file JetDeltaRValueMapProducer.cc.

70 {}

Member Function Documentation

template<class T , class C = T>
void JetDeltaRValueMapProducer< T, C >::produce ( edm::Event iEvent,
const edm::EventSetup iSetup 
)
inlineoverrideprivate

Definition at line 74 of file JetDeltaRValueMapProducer.cc.

References reco::deltaR2(), JetDeltaRValueMapProducer< T, C >::distMax_, JetDeltaRValueMapProducer< T, C >::evaluationMap_, funct::false, edm::helper::Filler< Map >::fill(), objects.autophobj::filler, edm::Event::getByToken(), mps_fire::i, edm::helper::Filler< Map >::insert(), JetDeltaRValueMapProducer< T, C >::matchedToken_, eostools::move(), JetDeltaRValueMapProducer< T, C >::multiValue_, edm::Event::put(), JetDeltaRValueMapProducer< T, C >::srcToken_, JetDeltaRValueMapProducer< T, C >::value_, JetDeltaRValueMapProducer< T, C >::valueLabels_, and MuonErrorMatrixValues_cff::values.

Referenced by JSONExport.JsonExport::export(), HTMLExport.HTMLExport::export(), and HTMLExport.HTMLExportStatic::export().

74  {
75 
77  iEvent.getByToken( srcToken_, h_jets1 );
79  iEvent.getByToken( matchedToken_, h_jets2 );
80 
81  std::vector<float> values( h_jets1->size(), -99999 );
82  std::map<std::string, std::vector<float> > valuesMap;
83  if( multiValue_ )
84  {
85  for( size_t i=0; i<valueLabels_.size(); ++i)
86  valuesMap.insert( std::make_pair( valueLabels_[i], std::vector<float>( h_jets1->size(), -99999 ) ) );
87  }
88  std::vector<bool> jets1_locks( h_jets1->size(), false );
89 
90  for ( typename edm::View<C>::const_iterator ibegin = h_jets2->begin(),
91  iend = h_jets2->end(), ijet = ibegin;
92  ijet != iend; ++ijet )
93  {
94  float matched_dR2 = 1e9;
95  int matched_index = -1;
96 
97  for ( typename edm::View<T>::const_iterator jbegin = h_jets1->begin(),
98  jend = h_jets1->end(), jjet = jbegin;
99  jjet != jend; ++jjet )
100  {
101  int index=jjet - jbegin;
102 
103  if( jets1_locks.at(index) ) continue; // skip jets that have already been matched
104 
105  float temp_dR2 = reco::deltaR2(ijet->eta(),ijet->phi(),jjet->eta(),jjet->phi());
106  if ( temp_dR2 < matched_dR2 )
107  {
108  matched_dR2 = temp_dR2;
109  matched_index = index;
110  }
111  }// end loop over src jets
112 
113  if( matched_index>=0 )
114  {
115  if ( matched_dR2 > distMax_*distMax_ )
116  edm::LogInfo("MatchedJetsFarApart") << "Matched jets separated by dR greater than distMax=" << distMax_;
117  else
118  {
119  jets1_locks.at(matched_index) = true;
120  if( value_!="" )
121  values.at(matched_index) = (*(evaluationMap_.at(value_)))(*ijet);
122  if( multiValue_ )
123  {
124  for( size_t i=0; i<valueLabels_.size(); ++i)
125  valuesMap.at(valueLabels_[i]).at(matched_index) = (*(evaluationMap_.at(valueLabels_[i])))(*ijet);
126  }
127  }
128  }
129  }// end loop over matched jets
130 
131  if( value_!="" )
132  {
133  std::unique_ptr< JetValueMap > jetValueMap ( new JetValueMap() );
134 
135  JetValueMap::Filler filler(*jetValueMap);
136  filler.insert(h_jets1, values.begin(), values.end());
137  filler.fill();
138 
139  // put in Event
140  iEvent.put(std::move(jetValueMap));
141  }
142  if( multiValue_ )
143  {
144  for( size_t i=0; i<valueLabels_.size(); ++i)
145  {
146  std::unique_ptr< JetValueMap > jetValueMap ( new JetValueMap() );
147 
148  JetValueMap::Filler filler(*jetValueMap);
149  filler.insert(h_jets1, valuesMap.at(valueLabels_[i]).begin(), valuesMap.at(valueLabels_[i]).end());
150  filler.fill();
151 
152  // put in Event
153  iEvent.put(std::move(jetValueMap), valueLabels_[i]);
154  }
155  }
156  }
edm::ValueMap< float > JetValueMap
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:137
const std::vector< std::string > valueLabels_
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:579
std::map< std::string, std::unique_ptr< const StringObjectFunction< C > > > evaluationMap_
const edm::EDGetTokenT< typename edm::View< C > > matchedToken_
const edm::EDGetTokenT< typename edm::View< T > > srcToken_
helper::Filler< ValueMap< float > > Filler
Definition: ValueMap.h:169
constexpr auto deltaR2(const T1 &t1, const T2 &t2) -> decltype(t1.eta())
Definition: deltaR.h:16
boost::indirect_iterator< typename seq_t::const_iterator > const_iterator
Definition: View.h:86
def move(src, dest)
Definition: eostools.py:511

Member Data Documentation

template<class T , class C = T>
const double JetDeltaRValueMapProducer< T, C >::distMax_
private
template<class T , class C = T>
std::map<std::string, std::unique_ptr<const StringObjectFunction<C> > > JetDeltaRValueMapProducer< T, C >::evaluationMap_
private
template<class T , class C = T>
const bool JetDeltaRValueMapProducer< T, C >::lazyParser_
private
template<class T , class C = T>
const edm::EDGetTokenT< typename edm::View<C> > JetDeltaRValueMapProducer< T, C >::matchedToken_
private
template<class T , class C = T>
bool JetDeltaRValueMapProducer< T, C >::multiValue_
private
template<class T , class C = T>
const edm::EDGetTokenT< typename edm::View<T> > JetDeltaRValueMapProducer< T, C >::srcToken_
private
template<class T , class C = T>
const std::string JetDeltaRValueMapProducer< T, C >::value_
private
template<class T , class C = T>
const std::vector<std::string> JetDeltaRValueMapProducer< T, C >::valueLabels_
private
template<class T , class C = T>
const std::vector<std::string> JetDeltaRValueMapProducer< T, C >::values_
private