CMS 3D CMS Logo

Public Member Functions | Private Types | Private Attributes

tautools::RecoTauDistanceFromTruthPlugin Class Reference

Inheritance diagram for tautools::RecoTauDistanceFromTruthPlugin:
reco::tau::RecoTauCleanerPlugin reco::tau::RecoTauEventHolderPlugin reco::tau::RecoTauNamedPlugin

List of all members.

Public Member Functions

void beginEvent ()
double operator() (const reco::PFTauRef &) const
 RecoTauDistanceFromTruthPlugin (const edm::ParameterSet &pset)
virtual ~RecoTauDistanceFromTruthPlugin ()

Private Types

typedef edm::Association
< reco::GenJetCollection
GenJetAssociation

Private Attributes

edm::Handle< GenJetAssociationgenTauMatch_
edm::InputTag matchingSrc_

Detailed Description

Definition at line 16 of file RecoTauDistanceFromTruthPlugin.cc.


Member Typedef Documentation

Definition at line 24 of file RecoTauDistanceFromTruthPlugin.cc.


Constructor & Destructor Documentation

tautools::RecoTauDistanceFromTruthPlugin::RecoTauDistanceFromTruthPlugin ( const edm::ParameterSet pset)
virtual tautools::RecoTauDistanceFromTruthPlugin::~RecoTauDistanceFromTruthPlugin ( ) [inline, virtual]

Definition at line 19 of file RecoTauDistanceFromTruthPlugin.cc.

{}

Member Function Documentation

void tautools::RecoTauDistanceFromTruthPlugin::beginEvent ( ) [virtual]
double tautools::RecoTauDistanceFromTruthPlugin::operator() ( const reco::PFTauRef tauRef) const [virtual]

Implements reco::tau::RecoTauCleanerPlugin.

Definition at line 75 of file RecoTauDistanceFromTruthPlugin.cc.

References abs, edm::Ptr< T >::get(), edm::Ref< C, T, F >::get(), infinity, edm::Ref< C, T, F >::isNull(), reco::Candidate::pt(), and python::multivaluedict::sort().

                                                                                  {

  GenJetAssociation::reference_type truth = (*genTauMatch_)[tauRef];

  // Check if the matching exists, if not return +infinity
  if (truth.isNull())
    return std::numeric_limits<double>::infinity();

  // screw all this noise
  return std::abs(tauRef->pt() - truth->pt());

  typedef const reco::Candidate* CandPtr;
  typedef std::set<CandPtr> CandPtrSet;
  typedef std::vector<CandPtr> CandPtrVector;
  typedef std::list<CandPtr> CandPtrList;
  // Store all of our reco and truth candidates
  CandPtrList recoCands;
  CandPtrSet truthCandSet;

  BOOST_FOREACH(const reco::RecoTauPiZero& piZero,
      tauRef->signalPiZeroCandidates()) {
    recoCands.push_back(&piZero);
  }

  BOOST_FOREACH(const reco::PFCandidateRef& pfch,
      tauRef->signalPFChargedHadrCands()) {
    recoCands.push_back(pfch.get());
  }

  // Use a set to contain the truth pointers so we ensure that no pizeros
  // are entered twice.
  BOOST_FOREACH(const reco::CandidatePtr& ptr,
      truth->daughterPtrVector()) {
    // Get mother pi zeros in the case of gammas
    if (ptr->pdgId() == 22)
      truthCandSet.insert(ptr->mother());
    else
      truthCandSet.insert(ptr.get());
  }

  //Convert truth cands from set to vector so we can sort it.
  CandPtrVector truthCands(truthCandSet.begin(), truthCandSet.end());

  // Sort the truth candidates by descending pt
  std::sort(truthCands.begin(), truthCands.end(),
      boost::bind(&reco::Candidate::pt, _1) > boost::bind(&reco::Candidate::pt, _2));

  double quality = 0.0;
  BOOST_FOREACH(CandPtr truthCand, truthCands) {
    // Find the best reco match for this truth cand
    CandPtrList::iterator recoMatch = findBestMatch(truthCand,
        recoCands.begin(), recoCands.end());

    // Check if this truth cand is matched
    if (recoMatch != recoCands.end()) {
      // Add a penalty factor based on how different the reconstructed
      // is w.r.t. the true momentum
      quality += momentumDifference(truthCand, *recoMatch);
      // Remove this reco cand from future matches
      recoCands.erase(recoMatch);
    } else {
      // this truth cand was not matched!
      quality += truthCand->p4().P2();
    }
  }

  // Now add a penalty for the unmatched reco stuff
  BOOST_FOREACH(CandPtr recoCand, recoCands) {
    quality += recoCand->p4().P2();
  }

  return quality;
}

Member Data Documentation

Definition at line 25 of file RecoTauDistanceFromTruthPlugin.cc.

Referenced by beginEvent().

Definition at line 23 of file RecoTauDistanceFromTruthPlugin.cc.

Referenced by beginEvent(), and RecoTauDistanceFromTruthPlugin().