CMSInsideOutAlgorithm.cc

Go to the documentation of this file.

#include "RecoJets/JetAlgorithms/interface/CMSInsideOutAlgorithm.h"
 
#include "RecoJets/JetAlgorithms/interface/CompoundPseudoJet.h"
 
using namespace std;
 
void CMSInsideOutAlgorithm::run(const std::vector<fastjet::PseudoJet>& fInput,
                                std::vector<fastjet::PseudoJet>& fOutput) {
  //make a list of input objects
  list<fastjet::PseudoJet> input;
  for (std::vector<fastjet::PseudoJet>::const_iterator candIter = fInput.begin(); candIter != fInput.end();
       ++candIter) {
    input.push_back(*candIter);
  }
 
  while (!input.empty() && input.front().perp() > seedThresholdPt_) {
    //get seed eta/phi
    double seedEta = input.front().eta();
    double seedPhi = input.front().phi();
 
    //find iterators to those objects that are in the max cone size
    list<inputListIter> maxCone;
    // add seed, then test elements after seed
    inputListIter iCand = input.begin();
    maxCone.push_back(iCand++);
    for (; iCand != input.end(); ++iCand) {
      const fastjet::PseudoJet& candidate = *iCand;
      if (reco::deltaR2(seedEta, candidate.eta(), seedPhi, candidate.phi()) < maxSizeSquared_)
        maxCone.push_back(iCand);
    }
    //sort objects by increasing DR about the seed  directions
    maxCone.sort(ListIteratorLesserByDeltaR(seedEta, seedPhi));
    list<inputListIter>::const_iterator position = maxCone.begin();
    bool limitReached = false;
    double totalET = (**position).perp();
    ++position;
    while (position != maxCone.end() && !limitReached) {
      const fastjet::PseudoJet& theCandidate = **position;
      double candidateET = theCandidate.perp() + totalET;
      double candDR2 = reco::deltaR2(seedEta, theCandidate.eta(), seedPhi, theCandidate.phi());
      if (candDR2 < minSizeSquared_ || candDR2 * candidateET * candidateET < growthParameterSquared_)
        totalET = candidateET;
      else
        limitReached = true;
      ++position;
    }
    //turn this into a final jet
    fastjet::PseudoJet final;
    for (list<inputListIter>::const_iterator iNewJet = maxCone.begin(); iNewJet != position; ++iNewJet) {
      final += **iNewJet;
      input.erase(*iNewJet);
    }
    fOutput.push_back(final);
  }  // end loop over seeds
  sort(fOutput.begin(), fOutput.end(), greaterByEtPseudoJet);
}