DeltaBetaWeights.cc

Go to the documentation of this file.

#include "CommonTools/ParticleFlow/plugins/DeltaBetaWeights.h"
 
DeltaBetaWeights::DeltaBetaWeights(const edm::ParameterSet& iConfig)
    : src_(iConfig.getParameter<edm::InputTag>("src")),
      pfCharged_(iConfig.getParameter<edm::InputTag>("chargedFromPV")),
      pfPU_(iConfig.getParameter<edm::InputTag>("chargedFromPU")) {
  produces<reco::PFCandidateCollection>();
 
  pfCharged_token = consumes<edm::View<reco::Candidate> >(pfCharged_);
  pfPU_token = consumes<edm::View<reco::Candidate> >(pfPU_);
  src_token = consumes<edm::View<reco::Candidate> >(src_);
 
  // pfCharged_token = consumes<reco::PFCandidateCollection>(pfCharged_);
  // pfPU_token = consumes<reco::PFCandidateCollection>(pfPU_);
  // src_token = consumes<reco::PFCandidateCollection>(src_);
}
 
DeltaBetaWeights::~DeltaBetaWeights() {
  // do anything here that needs to be done at desctruction time
  // (e.g. close files, deallocate resources etc.)
}
 
//
// member functions
//
 
// ------------ method called to produce the data  ------------
void DeltaBetaWeights::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
  using namespace edm;
 
  edm::Handle<edm::View<reco::Candidate> > pfCharged;
  edm::Handle<edm::View<reco::Candidate> > pfPU;
  edm::Handle<edm::View<reco::Candidate> > src;
 
  iEvent.getByToken(src_token, src);
  iEvent.getByToken(pfCharged_token, pfCharged);
  iEvent.getByToken(pfPU_token, pfPU);
 
  double sumNPU = .0;
  double sumPU = .0;
 
  std::unique_ptr<reco::PFCandidateCollection> out(new reco::PFCandidateCollection);
 
  for (const reco::Candidate& cand : *src) {
    if (cand.charge() != 0) {
      // this part of code should be executed only if input collection is not entirely composed of neutral candidates, i.e. never by default
      edm::LogWarning("DeltaBetaWeights")
          << "Trying to reweight charged particle... saving it to output collection without any change";
      out->emplace_back(cand.charge(), cand.p4(), reco::PFCandidate::ParticleType::X);
      (out->back()).setParticleType((out->back()).translatePdgIdToType(cand.pdgId()));
      continue;
    }
 
    sumNPU = 1.0;
    sumPU = 1.0;
    double eta = cand.eta();
    double phi = cand.phi();
    for (const reco::Candidate& chCand : *pfCharged) {
      double sum = (chCand.pt() * chCand.pt()) / (deltaR2(eta, phi, chCand.eta(), chCand.phi()));
      if (sum > 1.0)
        sumNPU *= sum;
    }
    sumNPU = 0.5 * log(sumNPU);
 
    for (const reco::Candidate& puCand : *pfPU) {
      double sum = (puCand.pt() * puCand.pt()) / (deltaR2(eta, phi, puCand.eta(), puCand.phi()));
      if (sum > 1.0)
        sumPU *= sum;
    }
    sumPU = 0.5 * log(sumPU);
 
    reco::PFCandidate neutral = reco::PFCandidate(cand.charge(), cand.p4(), reco::PFCandidate::ParticleType::X);
    neutral.setParticleType(neutral.translatePdgIdToType(cand.pdgId()));
    if (sumNPU + sumPU > 0)
      neutral.setP4(((sumNPU) / (sumNPU + sumPU)) * neutral.p4());
    out->push_back(neutral);
  }
 
  iEvent.put(std::move(out));
}