Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

More generic implementation of NjettinessAdder #4664

Merged
merged 3 commits into from Jul 18, 2014
Merged

More generic implementation of NjettinessAdder #4664

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
dc1d187
more generic implementation of NjettinessAdder so it works with any j…
ferencek Jul 16, 2014
994d68d
added LogWarning printout
ferencek Jul 16, 2014
a7c2001
further generalization for N's used in the tauN computation
ferencek Jul 16, 2014
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
31 changes: 19 additions & 12 deletions RecoJets/JetProducers/interface/NjettinessAdder.h
View file
Open in desktop
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
Open in desktop
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
Open in desktop
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)
)