Skip to content

Commit

Permalink
Merge pull request #4664 from ferencek/NjettinessAdderMadeMoreGeneric…
Browse files Browse the repository at this point in the history 
…_from-CMSSW_7_2_0_pre1

More generic implementation of NjettinessAdder
  • Loading branch information
@ktf
ktf committed Jul 18, 2014
2 parents 05547dc + a7c2001 commit 30c61a0
Show file tree
Hide file tree
Showing 3 changed files with 53 additions and 53 deletions.
31 changes: 19 additions & 12 deletions RecoJets/JetProducers/interface/NjettinessAdder.h
View file
Expand Up @@ -7,29 +7,36 @@
#include "FWCore/Framework/interface/EDProducer.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "FWCore/Utilities/interface/InputTag.h"
#include "DataFormats/PatCandidates/interface/Jet.h"
#include "DataFormats/JetReco/interface/Jet.h"
#include "DataFormats/Common/interface/ValueMap.h"

class NjettinessAdder : public edm::EDProducer {
public:
explicit NjettinessAdder(const edm::ParameterSet& iConfig) :
src_(iConfig.getParameter<edm::InputTag>("src")),
src_token_(consumes<edm::View<reco::PFJet>>(src_)),
cone_(iConfig.getParameter<double>("cone"))
explicit NjettinessAdder(const edm::ParameterSet& iConfig) :
src_(iConfig.getParameter<edm::InputTag>("src")),
src_token_(consumes<edm::View<reco::Jet>>(src_)),
cone_(iConfig.getParameter<double>("cone")),
Njets_(iConfig.getParameter<std::vector<unsigned> >("Njets"))
{
for ( std::vector<unsigned>::const_iterator n = Njets_.begin(); n != Njets_.end(); ++n )
{
produces<edm::ValueMap<float> >("tau1");
produces<edm::ValueMap<float> >("tau2");
produces<edm::ValueMap<float> >("tau3");
std::ostringstream tauN_str;
tauN_str << "tau" << *n;

produces<edm::ValueMap<float> >(tauN_str.str().c_str());
}
}

virtual ~NjettinessAdder() {}

void produce(edm::Event & iEvent, const edm::EventSetup & iSetup) ;
float getTau(int num,edm::Ptr<reco::PFJet> object) const;
float getTau(unsigned num, const edm::Ptr<reco::Jet> & object) const;

private:
edm::InputTag src_ ;
edm::EDGetTokenT<edm::View<reco::PFJet>> src_token_;
double cone_ ;
const edm::InputTag src_;
const edm::EDGetTokenT<edm::View<reco::Jet>> src_token_;
const double cone_ ;
const std::vector<unsigned> Njets_;
};

#endif
72 changes: 32 additions & 40 deletions RecoJets/JetProducers/plugins/NjettinessAdder.cc
View file
Expand Up @@ -5,57 +5,49 @@

void NjettinessAdder::produce(edm::Event & iEvent, const edm::EventSetup & iSetup) {
// read input collection
edm::Handle<edm::View<reco::PFJet> > jets;
edm::Handle<edm::View<reco::Jet> > jets;
iEvent.getByToken(src_token_, jets);

// prepare room for output
std::vector<float> tau1; tau1.reserve(jets->size());
std::vector<float> tau2; tau2.reserve(jets->size());
std::vector<float> tau3; tau3.reserve(jets->size());
for ( std::vector<unsigned>::const_iterator n = Njets_.begin(); n != Njets_.end(); ++n )
{
std::ostringstream tauN_str;
tauN_str << "tau" << *n;

for ( typename edm::View<reco::PFJet>::const_iterator jetIt = jets->begin() ; jetIt != jets->end() ; ++jetIt ) {
reco::PFJet newCand(*jetIt);
edm::Ptr<reco::PFJet> jetPtr = jets->ptrAt(jetIt - jets->begin());
// prepare room for output
std::vector<float> tauN;
tauN.reserve(jets->size());

float t1=getTau(1, jetPtr );
float t2=getTau(2, jetPtr );
float t3=getTau(3, jetPtr );
for ( typename edm::View<reco::Jet>::const_iterator jetIt = jets->begin() ; jetIt != jets->end() ; ++jetIt ) {

tau1.push_back(t1);
tau2.push_back(t2);
tau3.push_back(t3);
}
edm::Ptr<reco::Jet> jetPtr = jets->ptrAt(jetIt - jets->begin());

std::auto_ptr<edm::ValueMap<float> > outT1(new edm::ValueMap<float>());
std::auto_ptr<edm::ValueMap<float> > outT2(new edm::ValueMap<float>());
std::auto_ptr<edm::ValueMap<float> > outT3(new edm::ValueMap<float>());
edm::ValueMap<float>::Filler fillerT1(*outT1);
edm::ValueMap<float>::Filler fillerT2(*outT2);
edm::ValueMap<float>::Filler fillerT3(*outT3);
fillerT1.insert(jets, tau1.begin(), tau1.end());
fillerT2.insert(jets, tau2.begin(), tau2.end());
fillerT3.insert(jets, tau3.begin(), tau3.end());
fillerT1.fill();
fillerT2.fill();
fillerT3.fill();

iEvent.put(outT1,"tau1");
iEvent.put(outT2,"tau2");
iEvent.put(outT3,"tau3");
float t=getTau( *n, jetPtr );

tauN.push_back(t);
}

std::auto_ptr<edm::ValueMap<float> > outT(new edm::ValueMap<float>());
edm::ValueMap<float>::Filler fillerT(*outT);
fillerT.insert(jets, tauN.begin(), tauN.end());
fillerT.fill();

iEvent.put(outT,tauN_str.str().c_str());
}
}

float NjettinessAdder::getTau(int num, edm::Ptr<reco::PFJet> object) const
float NjettinessAdder::getTau(unsigned num, const edm::Ptr<reco::Jet> & object) const
{
std::vector<const reco::PFCandidate*> all_particles;
for (unsigned k =0; k < object->getPFConstituents().size(); k++)
all_particles.push_back( object->getPFConstituent(k).get() );

std::vector<fastjet::PseudoJet> FJparticles;
for (unsigned particle = 0; particle < all_particles.size(); particle++){
const reco::PFCandidate *thisParticle = all_particles.at(particle);
FJparticles.push_back( fastjet::PseudoJet( thisParticle->px(), thisParticle->py(), thisParticle->pz(), thisParticle->energy() ) );
for (unsigned k = 0; k < object->numberOfDaughters(); ++k)
{
const reco::CandidatePtr & dp = object->daughterPtr(k);
if ( dp.isNonnull() && dp.isAvailable() )
FJparticles.push_back( fastjet::PseudoJet( dp->px(), dp->py(), dp->pz(), dp->energy() ) );
else
edm::LogWarning("MissingJetConstituent") << "Jet constituent required for N-subjettiness computation is missing!";
}
fastjet::contrib::NsubParameters paraNsub = fastjet::contrib::NsubParameters(1.0, cone_); //assume R=0.7 jet clusering used

fastjet::contrib::NsubParameters paraNsub = fastjet::contrib::NsubParameters(1.0, cone_);
fastjet::contrib::Njettiness routine(fastjet::contrib::Njettiness::onepass_kt_axes, paraNsub);
return routine.getTau(num, FJparticles);
}
Expand Down
3 changes: 2 additions & 1 deletion RecoJets/JetProducers/python/nJettinessAdder_cfi.py
View file
Expand Up @@ -2,5 +2,6 @@

Njettiness = cms.EDProducer("NjettinessAdder",
src=cms.InputTag("ca8PFJetsCHS"),
cone=cms.double(0.8)
cone=cms.double(0.8),
Njets=cms.vuint32(1,2,3)
)

0 comments on commit 30c61a0

Please sign in to comment.