CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_4/src/RecoJets/JetAlgorithms/src/CMSInsideOutAlgorithm.cc

Go to the documentation of this file.
00001 #include "RecoJets/JetAlgorithms/interface/CMSInsideOutAlgorithm.h"
00002 
00003 #include "RecoJets/JetAlgorithms/interface/CompoundPseudoJet.h"
00004 
00005 
00006 using namespace std;
00007 
00008 void
00009 CMSInsideOutAlgorithm::run(const std::vector<fastjet::PseudoJet>& fInput, std::vector<fastjet::PseudoJet> & fOutput)
00010 {
00011 
00012    //make a list of input objects
00013    list<fastjet::PseudoJet> input;
00014    for (std::vector<fastjet::PseudoJet>::const_iterator candIter  = fInput.begin();
00015         candIter != fInput.end(); 
00016         ++candIter) {
00017       input.push_back(*candIter);
00018    }
00019 
00020    while( !input.empty() && input.front().perp() > seedThresholdPt_ ) 
00021    {
00022       //get seed eta/phi
00023       double seedEta = input.front().eta();
00024       double seedPhi = input.front().phi();
00025 
00026       //find iterators to those objects that are in the max cone size
00027       list<inputListIter> maxCone;
00028       // add seed, then test elements after seed
00029       inputListIter iCand = input.begin();
00030       maxCone.push_back(iCand++);
00031       for(; iCand != input.end(); ++iCand)
00032       {
00033         const fastjet::PseudoJet& candidate = *iCand;
00034          if( reco::deltaR2(seedEta, candidate.eta(), seedPhi, candidate.phi()) < maxSizeSquared_ )
00035             maxCone.push_back(iCand);
00036       }
00037       //sort objects by increasing DR about the seed  directions
00038       maxCone.sort(ListIteratorLesserByDeltaR(seedEta, seedPhi));
00039       list<inputListIter>::const_iterator position = maxCone.begin(); 
00040       bool limitReached = false;
00041       double totalET    = (**position).perp();
00042       ++position;
00043       while(position != maxCone.end() && !limitReached)
00044       {
00045          const fastjet::PseudoJet& theCandidate = **position;
00046          double candidateET    = theCandidate.perp() + totalET; 
00047          double candDR2        = reco::deltaR2(seedEta, theCandidate.eta(), seedPhi, theCandidate.phi());
00048          if( candDR2 < minSizeSquared_ ||  candDR2*candidateET*candidateET < growthParameterSquared_ )
00049             totalET = candidateET;
00050          else
00051             limitReached = true;
00052          ++position;
00053       }
00054       //turn this into a final jet
00055       fastjet::PseudoJet final;
00056       for(list<inputListIter>::const_iterator iNewJet  = maxCone.begin();
00057                                               iNewJet != position;
00058                                             ++iNewJet)
00059       {
00060          final += **iNewJet;
00061          input.erase(*iNewJet);
00062       }
00063       fOutput.push_back(final);
00064    } // end loop over seeds
00065    GreaterByEtPseudoJet compJets;
00066    sort (fOutput.begin (), fOutput.end (), compJets);
00067 }
00068          
00069