HistoryBase.cc

Go to the documentation of this file.

#include <algorithm>

#include "FWCore/MessageLogger/interface/MessageLogger.h"

#include "SimTracker/TrackHistory/interface/HistoryBase.h"

void HistoryBase::traceRecoGenHistory(reco::GenParticle const *genParticle) {
  // fill the trace of reco::GenParticle for correct reading of the History
  // flags etc. This is called from the TrackingParticle->genParticle_begin()
  // which is a reco::GenParticle

  // Take onlt genParticles with a status() smaller than "depth_". Typically
  // depth is 2 and so you trace back untill you reach the hard partons see for
  // status():
  // https://twiki.cern.ch/twiki/bin/view/CMSPublic/WorkBookGenParticleCandidate#GenPCand
  if (genParticle->status() <= abs(depth_) && (genParticle->pdgId() < 88 || genParticle->pdgId() > 99)) {
    // If the particle is already in the history, it is looping and you should
    // stop
    if (recoGenParticleTrailHelper_.find(genParticle) != recoGenParticleTrailHelper_.end()) {
      return;
    }
    recoGenParticleTrail_.push_back(genParticle);
    recoGenParticleTrailHelper_.insert(genParticle);
    // Get the genParticle's mother and trace its history
    if (genParticle->mother() != nullptr) {
      traceRecoGenHistory((const reco::GenParticle *)genParticle->mother());
    }
  }
}

void HistoryBase::traceGenHistory(HepMC::GenParticle const *genParticle) {
  // Third stop criteria: status abs(depth_) particles after the hadronization.
  // The after hadronization is done by detecting the pdg_id pythia code from 88
  // to 99
  if (genParticle->status() <= abs(depth_) && (genParticle->pdg_id() < 88 || genParticle->pdg_id() > 99)) {
    genParticleTrail_.push_back(genParticle);
    // Get the producer vertex and trace it history
    traceGenHistory(genParticle->production_vertex());
  }
}

void HistoryBase::traceGenHistory(HepMC::GenVertex const *genVertex) {
  // Verify if has a vertex associated
  if (genVertex) {
    // Skip if already exist in the collection
    if (genVertexTrailHelper_.find(genVertex) != genVertexTrailHelper_.end())
      return;
    // Add vertex to the history
    genVertexTrail_.push_back(genVertex);
    genVertexTrailHelper_.insert(genVertex);
    // Verify if the vertex has incoming particles
    if (genVertex->particles_in_size())
      traceGenHistory(*(genVertex->particles_in_const_begin()));
  }
}

bool HistoryBase::traceSimHistory(TrackingParticleRef const &trackingParticle, int depth) {
  // first stop condition: if the required depth is reached
  if (depth == depth_ && depth_ >= 0)
    return true;

  // second stop condition: if a gen particle is associated to the TP
  if (!trackingParticle->genParticles().empty()) {
    LogDebug("TrackHistory") << "Particle " << trackingParticle->pdgId() << " has a GenParicle image." << std::endl;

    traceRecoGenHistory(&(**(trackingParticle->genParticle_begin())));
  }

  LogDebug("TrackHistory") << "No GenParticle image for " << trackingParticle->pdgId() << std::endl;

  // if no association between TP and genParticles is found: get a reference to
  // the TP's parent vertex and trace its history
  return traceSimHistory(trackingParticle->parentVertex(), depth);
}

bool HistoryBase::traceSimHistory(TrackingVertexRef const &trackingVertex, int depth) {
  // verify if the parent vertex exists
  if (trackingVertex.isNonnull()) {
    // save the vertex in the trail
    simVertexTrail_.push_back(trackingVertex);

    if (!trackingVertex->sourceTracks().empty()) {
      LogDebug("TrackHistory") << "Moving on to the parent particle." << std::endl;

      // select the original source in case of combined vertices
      bool flag = false;
      TrackingVertex::tp_iterator itd, its;

      for (its = trackingVertex->sourceTracks_begin(); its != trackingVertex->sourceTracks_end(); its++) {
        for (itd = trackingVertex->daughterTracks_begin(); itd != trackingVertex->daughterTracks_end(); itd++)
          if (itd != its) {
            flag = true;
            break;
          }
        if (flag)
          break;
      }

      if (!flag)
        return false;

      // verify if the new particle is not in the trail (looping partiles)
      if (std::find(simParticleTrail_.begin(), simParticleTrail_.end(), *its) != simParticleTrail_.end()) {
        LogDebug("TrackHistory") << "WARNING: Looping track found." << std::endl;
        return false;
      }

      // save particle in the trail
      simParticleTrail_.push_back(*its);
      return traceSimHistory(*its, --depth);
    } else if (!trackingVertex->genVertices().empty()) {
      // navigate over all the associated generated vertexes
      LogDebug("TrackHistory") << "Vertex has " << trackingVertex->genVertices().size() << "GenVertex image."
                               << std::endl;
      for (TrackingVertex::genv_iterator ivertex = trackingVertex->genVertices_begin();
           ivertex != trackingVertex->genVertices_end();
           ++ivertex)
        traceGenHistory(&(**(ivertex)));
      return true;
    } else {
      LogDebug("TrackHistory") << "WARNING: Source track for tracking vertex cannot be found." << std::endl;
    }
  } else {
    LogDebug("TrackHistory") << " WARNING: Vertex cannot be found.";
  }

  return false;
}