BPHRecoBuilder.cc

Go to the documentation of this file.

/*
 *  See header file for a description of this class.
 *
 *  \author Paolo Ronchese INFN Padova
 *
 */

//-----------------------
// This Class' Header --
//-----------------------
#include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoBuilder.h"

//-------------------------------
// Collaborating Class Headers --
//-------------------------------
#include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoSelect.h"
#include "HeavyFlavorAnalysis/RecoDecay/interface/BPHMomentumSelect.h"
#include "HeavyFlavorAnalysis/RecoDecay/interface/BPHVertexSelect.h"
#include "HeavyFlavorAnalysis/RecoDecay/interface/BPHFitSelect.h"
#include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoCandidate.h"
#include "HeavyFlavorAnalysis/RecoDecay/interface/BPHTrackReference.h"
#include "DataFormats/RecoCandidate/interface/RecoCandidate.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"

//---------------
// C++ Headers --
//---------------
#include <iostream>
#include <map>

using namespace std;

//-------------------
// Initializations --
//-------------------

//----------------
// Constructors --
//----------------
BPHRecoBuilder::BPHRecoBuilder(const edm::EventSetup& es) : evSetup(&es), minPDiff(-1.0) {
  msList.reserve(5);
  vsList.reserve(5);
}

//--------------
// Destructor --
//--------------
BPHRecoBuilder::~BPHRecoBuilder() {
  int n = sourceList.size();
  while (n--) {
    delete sourceList[n]->collection;
    delete sourceList[n];
  }
  int m = srCompList.size();
  while (m--)
    delete srCompList[m];
  while (!compCollectList.empty()) {
    const vector<BPHRecoConstCandPtr>* cCollection = *compCollectList.begin();
    delete cCollection;
    compCollectList.erase(cCollection);
  }
}

//--------------
// Operations --
//--------------
BPHRecoBuilder::BPHGenericCollection* BPHRecoBuilder::createCollection(const vector<const reco::Candidate*>& candList,
                                                                       const string& list) {
  return new BPHSpecificCollection<vector<const reco::Candidate*> >(candList, list);
}

template <>
const reco::Candidate& BPHRecoBuilder::BPHSpecificCollection<vector<const reco::Candidate*> >::get(int i) const {
  return *(*cPtr)[i];
}

void BPHRecoBuilder::filter(const string& name, const BPHRecoSelect& sel) const {
  map<string, int>::const_iterator iter = sourceId.find(name);
  if (iter == sourceId.end())
    return;
  BPHRecoSource* rs = sourceList[iter->second];
  rs->selector.push_back(&sel);
  return;
}

void BPHRecoBuilder::filter(const string& name, const BPHMomentumSelect& sel) const {
  map<string, int>::const_iterator iter = srCompId.find(name);
  if (iter == sourceId.end())
    return;
  BPHCompSource* cs = srCompList[iter->second];
  cs->momSelector.push_back(&sel);
  return;
}

void BPHRecoBuilder::filter(const string& name, const BPHVertexSelect& sel) const {
  map<string, int>::const_iterator iter = srCompId.find(name);
  if (iter == sourceId.end())
    return;
  BPHCompSource* cs = srCompList[iter->second];
  cs->vtxSelector.push_back(&sel);
  return;
}

void BPHRecoBuilder::filter(const string& name, const BPHFitSelect& sel) const {
  map<string, int>::const_iterator iter = srCompId.find(name);
  if (iter == sourceId.end())
    return;
  BPHCompSource* cs = srCompList[iter->second];
  cs->fitSelector.push_back(&sel);
  return;
}

void BPHRecoBuilder::filter(const BPHMomentumSelect& sel) {
  msList.push_back(&sel);
  return;
}

void BPHRecoBuilder::filter(const BPHVertexSelect& sel) {
  vsList.push_back(&sel);
  return;
}

void BPHRecoBuilder::filter(const BPHFitSelect& sel) {
  fsList.push_back(&sel);
  return;
}

bool BPHRecoBuilder::accept(const BPHRecoCandidate& cand) const {
  int i;
  int n;
  n = msList.size();
  for (i = 0; i < n; ++i) {
    if (!msList[i]->accept(cand))
      return false;
  }
  n = vsList.size();
  for (i = 0; i < n; ++i) {
    if (!vsList[i]->accept(cand))
      return false;
  }
  n = fsList.size();
  for (i = 0; i < n; ++i) {
    if (!fsList[i]->accept(cand))
      return false;
  }
  return true;
}

void BPHRecoBuilder::setMinPDiffererence(double pMin) {
  minPDiff = pMin;
  return;
}

vector<BPHRecoBuilder::ComponentSet> BPHRecoBuilder::build() const {
  daugMap.clear();
  compMap.clear();
  vector<ComponentSet> candList;
  ComponentSet compSet;
  build(candList, compSet, sourceList.begin(), sourceList.end(), srCompList.begin(), srCompList.end());
  return candList;
}

const edm::EventSetup* BPHRecoBuilder::eventSetup() const { return evSetup; }

bool BPHRecoBuilder::sameTrack(const reco::Candidate* lCand, const reco::Candidate* rCand, double minPDifference) {
  const reco::Track* lrcTrack = BPHTrackReference::getFromRC(*lCand);
  const reco::Track* rrcTrack = BPHTrackReference::getFromRC(*rCand);
  const reco::Track* lpfTrack = BPHTrackReference::getFromPF(*lCand);
  const reco::Track* rpfTrack = BPHTrackReference::getFromPF(*rCand);
  if ((lrcTrack != nullptr) && ((lrcTrack == rrcTrack) || (lrcTrack == rpfTrack)))
    return true;
  if ((lpfTrack != nullptr) && ((lpfTrack == rrcTrack) || (lpfTrack == rpfTrack)))
    return true;
  reco::Candidate::Vector pDiff = (lCand->momentum() - rCand->momentum());
  reco::Candidate::Vector pMean = (lCand->momentum() + rCand->momentum());
  double pDMod = pDiff.mag2();
  double pMMod = pMean.mag2();
  if (((pDMod / pMMod) < minPDifference) && (lCand->charge() == rCand->charge()))
    return true;
  return false;
}

void BPHRecoBuilder::add(const string& name, const BPHGenericCollection* collection, double mass, double msig) {
  BPHRecoSource* rs;
  if (sourceId.find(name) != sourceId.end()) {
    edm::LogPrint("TooManyParticles") << "BPHRecoBuilder::add: "
                                      << "Decay product already inserted with name " << name << " , skipped";
    return;
  }
  rs = new BPHRecoSource;
  rs->name = &sourceId.insert(make_pair(name, sourceList.size())).first->first;
  rs->collection = collection;
  rs->selector.reserve(5);
  rs->mass = mass;
  rs->msig = msig;
  sourceList.push_back(rs);
  return;
}

void BPHRecoBuilder::add(const string& name, const vector<BPHRecoConstCandPtr>& collection) {
  BPHCompSource* cs;
  if (srCompId.find(name) != srCompId.end()) {
    edm::LogPrint("TooManyParticles") << "BPHRecoBuilder::add: "
                                      << "Decay product already inserted with name " << name << " , skipped";
    return;
  }
  cs = new BPHCompSource;
  cs->name = &srCompId.insert(make_pair(name, srCompList.size())).first->first;
  cs->collection = &collection;
  srCompList.push_back(cs);
  return;
}

void BPHRecoBuilder::build(vector<ComponentSet>& compList,
                           ComponentSet& compSet,
                           vector<BPHRecoSource*>::const_iterator r_iter,
                           vector<BPHRecoSource*>::const_iterator r_iend,
                           vector<BPHCompSource*>::const_iterator c_iter,
                           vector<BPHCompSource*>::const_iterator c_iend) const {
  if (r_iter == r_iend) {
    if (c_iter == c_iend) {
      compSet.compMap = compMap;
      compList.push_back(compSet);
      return;
    }
    BPHCompSource* source = *c_iter++;
    const vector<BPHRecoConstCandPtr>* collection = source->collection;
    vector<const BPHMomentumSelect*> momSelector = source->momSelector;
    vector<const BPHVertexSelect*> vtxSelector = source->vtxSelector;
    vector<const BPHFitSelect*> fitSelector = source->fitSelector;
    int i;
    int j;
    int n = collection->size();
    int m;
    bool skip;
    for (i = 0; i < n; ++i) {
      skip = false;
      BPHRecoConstCandPtr cand = collection->at(i);
      if (contained(compSet, cand))
        continue;
      m = momSelector.size();
      for (j = 0; j < m; ++j) {
        if (!momSelector[j]->accept(*cand))
          skip = true;
      }
      m = vtxSelector.size();
      for (j = 0; j < m; ++j) {
        if (!vtxSelector[j]->accept(*cand))
          skip = true;
      }
      m = fitSelector.size();
      for (j = 0; j < m; ++j) {
        if (!fitSelector[j]->accept(*cand))
          skip = true;
      }
      if (skip)
        continue;
      compMap[*source->name] = cand;
      build(compList, compSet, r_iter, r_iend, c_iter, c_iend);
      compMap.erase(*source->name);
    }
    return;
  }
  BPHRecoSource* source = *r_iter++;
  const BPHGenericCollection* collection = source->collection;
  vector<const BPHRecoSelect*>& selector = source->selector;
  int i;
  int j;
  int n = collection->size();
  int m = selector.size();
  bool skip;
  for (i = 0; i < n; ++i) {
    const reco::Candidate& cand = collection->get(i);
    if (contained(compSet, &cand))
      continue;
    skip = false;
    for (j = 0; j < m; ++j) {
      if (!selector[j]->accept(cand, this))
        skip = true;
    }
    if (skip)
      continue;
    BPHDecayMomentum::Component comp;
    comp.cand = &cand;
    comp.mass = source->mass;
    comp.msig = source->msig;
    comp.searchList = collection->searchList();
    compSet.daugMap[*source->name] = comp;
    daugMap[*source->name] = &cand;
    build(compList, compSet, r_iter, r_iend, c_iter, c_iend);
    daugMap.erase(*source->name);
    compSet.daugMap.erase(*source->name);
  }
  return;
}

bool BPHRecoBuilder::contained(ComponentSet& compSet, const reco::Candidate* cand) const {
  map<string, BPHDecayMomentum::Component>& dMap = compSet.daugMap;
  map<string, BPHDecayMomentum::Component>::const_iterator d_iter;
  map<string, BPHDecayMomentum::Component>::const_iterator d_iend = dMap.end();
  for (d_iter = dMap.begin(); d_iter != d_iend; ++d_iter) {
    const reco::Candidate* cChk = d_iter->second.cand;
    if (cand == cChk)
      return true;
    if (sameTrack(cand, cChk))
      return true;
  }
  return false;
}

bool BPHRecoBuilder::contained(ComponentSet& compSet, BPHRecoConstCandPtr cCand) const {
  map<string, BPHRecoConstCandPtr>::const_iterator c_iter;
  map<string, BPHRecoConstCandPtr>::const_iterator c_iend = compMap.end();
  const vector<const reco::Candidate*>& dCand = cCand->daughFull();
  int j;
  int m = dCand.size();
  int k;
  int l;
  for (j = 0; j < m; ++j) {
    const reco::Candidate* cand = cCand->originalReco(dCand[j]);
    map<string, BPHDecayMomentum::Component>& dMap = compSet.daugMap;
    map<string, BPHDecayMomentum::Component>::const_iterator d_iter;
    map<string, BPHDecayMomentum::Component>::const_iterator d_iend = dMap.end();
    for (d_iter = dMap.begin(); d_iter != d_iend; ++d_iter) {
      const reco::Candidate* cChk = d_iter->second.cand;
      if (cand == cChk)
        return true;
      if (sameTrack(cand, cChk))
        return true;
    }

    for (c_iter = compMap.begin(); c_iter != c_iend; ++c_iter) {
      const pair<string, BPHRecoConstCandPtr>& entry = *c_iter;
      BPHRecoConstCandPtr cCChk = entry.second;
      const vector<const reco::Candidate*>& dCChk = cCChk->daughFull();
      l = dCChk.size();
      for (k = 0; k < l; ++k) {
        const reco::Candidate* cChk = cCChk->originalReco(dCChk[k]);
        if (cand == cChk)
          return true;
        if (sameTrack(cand, cChk))
          return true;
      }
    }
  }

  return false;
}

bool BPHRecoBuilder::sameTrack(const reco::Candidate* lCand, const reco::Candidate* rCand) const {
  return sameTrack(lCand, rCand, minPDiff);
}