CMS 3D CMS Logo

CMSMidpointAlgorithm Class Reference

CMSMidpointAlgorithm is an algorithm for CMS jet reconstruction baded on the CDF midpoint algorithm. More...

#include <RecoJets/JetAlgorithms/interface/CMSMidpointAlgorithm.h>

List of all members.

Public Types

typedef std::vector< ProtoJet * > InternalCollection

Public Member Functions

 CMSMidpointAlgorithm (double fSeedThreshold, double fConeRadius, double fConeAreaFraction, int fMaxPairSize, int fMaxIterations, double fOverlapThreshold, int fDebugLevel)
 Constructor takes as input all the values of the algorithm that the user can change.
 CMSMidpointAlgorithm ()
 Default constructor.
void run (const JetReco::InputCollection &fInput, JetReco::OutputCollection *fOutput)
 Runs the algorithm and returns a list of caloJets.

Private Member Functions

void addClustersToPairs (const JetReco::InputCollection &fInput, std::vector< int > &testPair, std::vector< std::vector< int > > &pairs, std::vector< std::vector< bool > > &distanceOK, int maxClustersInPair)
 Add proto-jets to pairs, triplets, etc, prior to finding their midpoints.
void findStableConesFromMidPoints (const JetReco::InputCollection &fInput, InternalCollection *fOutput)
 Add to the list of proto-jets the list of midpoints.
void findStableConesFromSeeds (const JetReco::InputCollection &fInput, InternalCollection *fOutput)
 Find the list of proto-jets from the seeds.
void iterateCone (const JetReco::InputCollection &fInput, double startRapidity, double startPhi, double startPt, bool reduceConeSize, InternalCollection *fOutput)
 Iterate the proto-jet center until it is stable.
void splitAndMerge (const JetReco::InputCollection &fInput, InternalCollection *fProtoJets, JetReco::OutputCollection *fFinalJets)
 Split and merge the proto-jets, assigning each tower in the protojets to one and only one final jet.

Private Attributes

double theConeAreaFraction
double theConeRadius
int theDebugLevel
int theMaxIterations
int theMaxPairSize
double theOverlapThreshold
double theSeedThreshold


Detailed Description

CMSMidpointAlgorithm is an algorithm for CMS jet reconstruction baded on the CDF midpoint algorithm.

The algorithm is documented in the proceedings of the Physics at RUN II: QCD and Weak Boson Physics Workshop: hep-ex/0005012.

The algorithm runs off of generic Candidates

Author:
Robert M Harris, Fermilab
Version:
1st Version Feb. 4, 2005 Based on the CDF Midpoint Algorithm code by Matthias Toennesmann.

2nd Version Apr. 6, 2005 Modifications toward integration in new EDM.

3rd Version Oct. 19, 2005 Modified to work with real CaloTowers from Jeremy Mans

F.Ratnikov, Mar. 8, 2006. Work from Candidate

Id
CMSMidpointAlgorithm.h,v 1.8 2007/05/03 21:20:09 fedor Exp

Definition at line 27 of file CMSMidpointAlgorithm.h.


Member Typedef Documentation

typedef std::vector<ProtoJet*> CMSMidpointAlgorithm::InternalCollection

Definition at line 30 of file CMSMidpointAlgorithm.h.


Constructor & Destructor Documentation

CMSMidpointAlgorithm::CMSMidpointAlgorithm (  )  [inline]

Default constructor.

Definition at line 33 of file CMSMidpointAlgorithm.h.

00033                           :
00034     theSeedThreshold(3.0),
00035     theConeRadius(0.5),
00036     theConeAreaFraction(1.0),
00037     theMaxPairSize(2),
00038     theMaxIterations(100),
00039     theOverlapThreshold(0.75),
00040     theDebugLevel(0)
00041   { }

CMSMidpointAlgorithm::CMSMidpointAlgorithm ( double  fSeedThreshold,
double  fConeRadius,
double  fConeAreaFraction,
int  fMaxPairSize,
int  fMaxIterations,
double  fOverlapThreshold,
int  fDebugLevel 
) [inline]

Constructor takes as input all the values of the algorithm that the user can change.

Parameters:
fSeedThreshold,: minimum ET in GeV of an CaloTower that can seed a jet.
fTowerThreshold,: minimum ET in GeV of an CaloTower that is included in a jet
fConeRadius,: nominal radius of the jet in eta-phi space
fConeAreaFraction,: multiplier to reduce the search cone area during the iteration phase. introduced by CDF in 2002 to avoid energy loss due to proto-jet migration. Actively being discussed in 2005. A value of 1.0 gives the original run II algorithm in hep-ex/0005012. CDF value of 0.25 gives new search cone of 0.5 theConeRadius during iteration, but keeps the final cone at theConeRadius for the last iteration.
fMaxPairSize,: Maximum size of proto-jet pair, triplet, etc for defining midpoint. Both CDF and D0 use 2.
fMaxIterations,: Maximum number of iterations before finding a stable cone.
fOverlapThreshold,: When two proto-jets overlap, this is the merging threshold on the fraction of PT in the overlap region compared to the lower Pt Jet. If the overlapPt/lowerJetPt > theOverlapThreshold the 2 proto-jets will be merged into one final jet. if the overlapPt/lowerJetPt < theOverlapThreshold the towers in the two proto-jets will be seperated into two distinct sets of towers: two final jets. D0 has used 0.5, CDF has used both 0.5 and 0.75 in run 2, and 0.75 in run 1.
fDebugLevel,: integer level of diagnostic printout: 0 = no printout, 1 = minimal printout, 2 = pages per event.

Definition at line 65 of file CMSMidpointAlgorithm.h.

00066                                                                                                         : 
00067     theSeedThreshold(fSeedThreshold),
00068     theConeRadius(fConeRadius),
00069     theConeAreaFraction(fConeAreaFraction),
00070     theMaxPairSize(fMaxPairSize),
00071     theMaxIterations(fMaxIterations),
00072     theOverlapThreshold(fOverlapThreshold),
00073     theDebugLevel(fDebugLevel)
00074   { }


Member Function Documentation

void CMSMidpointAlgorithm::addClustersToPairs ( const JetReco::InputCollection fInput,
std::vector< int > &  testPair,
std::vector< std::vector< int > > &  pairs,
std::vector< std::vector< bool > > &  distanceOK,
int  maxClustersInPair 
) [private]

Add proto-jets to pairs, triplets, etc, prior to finding their midpoints.

Called by findStableConesFromMidPoints but public method to allow testing and studies.

Referenced by findStableConesFromMidPoints().

void CMSMidpointAlgorithm::findStableConesFromMidPoints ( const JetReco::InputCollection fInput,
InternalCollection fOutput 
) [private]

Add to the list of proto-jets the list of midpoints.

Called by run but public method to allow testing and studies.

Definition at line 236 of file CMSMidpointAlgorithm.cc.

References addClustersToPairs(), GenMuonPlsPt100GeV_cfg::cout, deltaR2(), lat::endl(), iterateCone(), mergeAndRegister_online::pairs, ProtoJet::phi(), python::multivaluedict::sort(), theConeRadius, theDebugLevel, theMaxPairSize, and ProtoJet::y().

Referenced by run().

00236                                                                                                                      {
00237   // We take the previous list of stable protojets from seeds as input and add to it
00238   // those from the midpoints between proto-jet pairs, triplets, etc.
00239   // distanceOK[i-1][j] = Is distance between stableCones i and j (i>j) less than 2*_theConeRadius?
00240   vector< vector<bool> > distanceOK;         // A vector of vectors
00241   distanceOK.resize(stableCones->size() - 1); // Set the outer vector size, num protojets - 1
00242   for(unsigned int nCluster1 = 1; nCluster1 < stableCones->size(); ++nCluster1){  // Loop over the protojets
00243     distanceOK[nCluster1 - 1].resize(nCluster1);               //Set inner vector size: monotonically increasing.
00244     const ProtoJet* cluster1 = (*stableCones)[nCluster1];
00245     for(unsigned int nCluster2 = 0; nCluster2 < nCluster1; ++nCluster2){         // Loop over the other proto-jets
00246       const ProtoJet* cluster2 = (*stableCones)[nCluster2];
00247       double dR2 = deltaR2 (cluster1->y(), cluster1->phi(), cluster2->y(), cluster2->phi());
00248       distanceOK[nCluster1 - 1][nCluster2] = dR2 < 4*theConeRadius*theConeRadius;
00249     }
00250   }
00251 
00252   // Find all pairs (triplets, ...) of stableCones which are less than 2*theConeRadius apart from each other.
00253   vector< vector<int> > pairs(0);  
00254   vector<int> testPair(0);
00255   int maxClustersInPair = theMaxPairSize;  // Set maximum proto-jets to a pair or a triplet, etc.
00256   if(!maxClustersInPair)maxClustersInPair = stableCones->size();  // If zero, then skys the limit!
00257   addClustersToPairs(fInput, testPair,pairs,distanceOK,maxClustersInPair);  // Make the pairs, triplets, etc.
00258   
00259   // Loop over all combinations. Calculate MidPoint. Make midPointClusters.
00260   bool reduceConeSize = false;  // Note that here we keep the iteration cone size fixed.
00261   for(vector<vector<int> >::const_iterator iPair = pairs.begin(); iPair != pairs.end(); ++iPair) {
00262     const vector<int> & Pair = *iPair;
00263     // Calculate rapidity, phi and energy of MidPoint.
00264     reco::Particle::LorentzVector midPoint (0,0,0,0);
00265     for(vector<int>::const_iterator iPairMember = Pair.begin(); iPairMember != Pair.end(); ++iPairMember) {
00266       midPoint += (*stableCones)[*iPairMember]->p4();
00267     }
00268     if(theDebugLevel>=2)cout << endl << "[CMSMidpointAlgorithm] midpoint " << iPair-pairs.begin() << ": y = " << midPoint.Rapidity() << ", phi=" << midPoint.Phi() <<
00269                           ", size=" << Pair.size() << endl;
00270     iterateCone(fInput, midPoint.Rapidity(),midPoint.Phi(),midPoint.e(),reduceConeSize,stableCones);
00271   }
00272   GreaterByPtPtr<ProtoJet> compJets;
00273   sort (stableCones->begin(), stableCones->end(), compJets);
00274 }

void CMSMidpointAlgorithm::findStableConesFromSeeds ( const JetReco::InputCollection fInput,
InternalCollection fOutput 
) [private]

Find the list of proto-jets from the seeds.

Called by run, but public method to allow testing and studies.

Definition at line 130 of file CMSMidpointAlgorithm.cc.

References GenMuonPlsPt100GeV_cfg::cout, lat::endl(), etOrderedCaloTowers(), i, iterateCone(), theDebugLevel, and theSeedThreshold.

Referenced by run().

00130                                                                                                               {
00131   // This dictates that the cone size will be reduced in the iterations procedure,
00132   // to prevent excessive cone movement, and then will be enlarged at the end of
00133   // the iteration procedure (ala CDF).
00134   bool reduceConeSize = true;  
00135   
00136   // Get the Seed Towers sorted by Et.  
00137   vector<InputItem> seedTowers = etOrderedCaloTowers(fInput, theSeedThreshold);  //This gets towers
00138 
00139   // Loop over all Seeds
00140   for(vector<InputItem>::const_iterator i = seedTowers.begin(); i != seedTowers.end(); ++i) {
00141     double seedEta = (*i)->eta ();
00142     double seedPhi = (*i)->phi (); 
00143 
00144     // Find stable cone from seed.  
00145     // This iterates the cone centroid, makes a proto-jet, and adds it to the list.
00146     if(theDebugLevel>=2) cout << endl << "[CMSMidpointAlgorithm] seed " << i-seedTowers.begin() << ": eta=" << seedEta << ", phi=" << seedPhi << endl;
00147     iterateCone(fInput, seedEta, seedPhi, 0, reduceConeSize, fOutput);
00148     
00149   }
00150 
00151 }

void CMSMidpointAlgorithm::iterateCone ( const JetReco::InputCollection fInput,
double  startRapidity,
double  startPhi,
double  startPt,
bool  reduceConeSize,
InternalCollection fOutput 
) [private]

Iterate the proto-jet center until it is stable.

Called by findStableConesFromSeeds and findStableConesFromMidPoints but public method to allow testing and studies.

Definition at line 155 of file CMSMidpointAlgorithm.cc.

References funct::abs(), GenMuonPlsPt100GeV_cfg::cout, ProtoJet::e(), lat::endl(), ProtoJet::getTowerList(), i, ProtoJet::p4(), ProtoJet::phi(), ProtoJet::pt(), ProtoJet::putTowers(), funct::sqrt(), t, theConeAreaFraction, theConeRadius, theDebugLevel, theMaxIterations, towersWithinCone(), and ProtoJet::y().

Referenced by findStableConesFromMidPoints(), and findStableConesFromSeeds().

00160                                                                         {
00161   //  The workhorse of the algorithm.
00162   //  Throws a cone around a position and iterates the cone until the centroid and energy of the clsuter is stable.
00163   //  Uses a reduced cone size to prevent excessive cluster drift, ala CDF.
00164   //  Adds unique clusters to the list of proto-jets.
00165   //
00166   int nIterations = 0;
00167   bool keepJet = true;
00168   ProtoJet* trialCone = new ProtoJet;
00169   double iterationtheConeRadius = theConeRadius;
00170   if(reduceConeSize)iterationtheConeRadius *= sqrt(theConeAreaFraction);
00171   while(++nIterations <= theMaxIterations + 1 && keepJet){
00172     
00173     // Last iteration uses the full cone size.  Others use reduced cone size.
00174     if(nIterations == theMaxIterations + 1)iterationtheConeRadius = theConeRadius;
00175     
00176     //Add all towers in cone and over threshold to the cluster
00177     vector<InputItem> towersInSeedCluster 
00178       = towersWithinCone(fInput, startRapidity, startPhi, iterationtheConeRadius);
00179     if(theDebugLevel>=2)cout << "[CMSMidpointAlgorithm] iter=" << nIterations << ", towers=" <<towersInSeedCluster.size();    
00180     
00181     if(towersInSeedCluster.size()<1) {  // Just in case there is an empty cone.
00182       keepJet = false;
00183       if(theDebugLevel>=2)cout << endl;
00184     } else{
00185       //Put seed cluster into trial cone 
00186       trialCone->putTowers(towersInSeedCluster);
00187       double endRapidity = trialCone->y();
00188       double endPhi      = trialCone->phi();
00189       double endPt       = trialCone->pt();
00190       double endE        = trialCone->e();
00191       if(theDebugLevel>=2)cout << ", y=" << endRapidity << ", phi=" << endPhi << ", PT=" << endPt << ", constituents=" << trialCone->getTowerList().size () << endl;
00192       if(nIterations <= theMaxIterations){
00193         // Do we have a stable cone?
00194         if(abs(endRapidity-startRapidity)<.001 && 
00195            abs(endPhi-startPhi)<.001 && 
00196            abs(endE - startE) < .001) {
00197           nIterations = theMaxIterations;       // If cone size is reduced, then still one more iteration.
00198           if(!reduceConeSize) ++nIterations; // Otherwise, this is the last iteration.
00199         }
00200         else{
00201           // Another iteration.
00202           startRapidity = endRapidity;
00203           startPhi      = endPhi;
00204           startE        = endE;
00205         }
00206       } 
00207       if(nIterations==theMaxIterations+1) {
00208         for(vector<InputItem>::const_iterator i = towersInSeedCluster.begin(); i != towersInSeedCluster.end(); ++i) {
00209           InputItem t = *i;
00210           if(theDebugLevel>=2) cout << "[CMSMidpointAlgorithm] Tower " <<
00211                                  i-towersInSeedCluster.begin() << ": eta=" << t->eta() << 
00212                                  ", phi=" << t->phi() << ", ET="  << t->et() << endl;
00213         }
00214       }         
00215     }
00216   }
00217   
00218   if(keepJet){  // Our trial cone is now stable
00219     bool identical = false;
00220     // Loop over proto-jets and check that our trial cone is a unique proto-jet
00221     for (unsigned icone = 0; icone < stableCones->size(); icone++) {
00222       if (trialCone->p4() == (*stableCones) [icone]->p4()) identical = true;  // This proto-jet is not unique.
00223     }
00224     if(!identical){
00225       stableCones->push_back(trialCone);  // Save the unique proto-jets
00226       trialCone = 0;
00227       if(theDebugLevel>=2)cout << "[CMSMidpointAlgorithm] Unique Proto-Jet Saved" << endl;
00228     }    
00229   }
00230   delete trialCone;
00231 }

void CMSMidpointAlgorithm::run ( const JetReco::InputCollection fInput,
JetReco::OutputCollection fOutput 
)

Runs the algorithm and returns a list of caloJets.

The user declares the vector and calls this method.

Definition at line 92 of file CMSMidpointAlgorithm.cc.

References TestMuL1L2Filter_cff::cerr, GenMuonPlsPt100GeV_cfg::cout, lat::endl(), findStableConesFromMidPoints(), findStableConesFromSeeds(), index, splitAndMerge(), and theDebugLevel.

Referenced by FWLiteJetProducer::makeMidpointJets(), and cms::MidpointPilupSubtractionJetProducer::runAlgorithm().

00093 {
00094   if (!fOutput) {
00095     std::cerr << "CMSMidpointAlgorithm::run-> ERROR: no output collection" << std::endl;
00096   }
00097   if(theDebugLevel>=1)cout << "[CMSMidpointAlgorithm] Num of input constituents = " << fInput.size() << endl;
00098   if(theDebugLevel>=2) {
00099     unsigned index = 0;
00100     for (; index < fInput.size (); index++) {
00101       cout << index << " consituent p/eta/phi: " 
00102            << fInput[index]->p() << '/'
00103            << fInput[index]->eta() << '/'
00104            << fInput[index]->phi() << std::endl;
00105     }
00106   }
00107   // Find proto-jets from the seeds.
00108   InternalCollection protoJets;         // Initialize working container
00109   // vector<ProtoJet> finalJets;         // Final proto-jets container
00110   findStableConesFromSeeds(fInput, &protoJets);   //Find the proto-jets from the seeds
00111   if(theDebugLevel>=1)cout << "[CMSMidpointAlgorithm] Num proto-jets from Seeds = " << protoJets.size() << endl;
00112 
00113   // Find proto-jets from the midpoints, and add them to the list
00114   if(protoJets.size()>0)findStableConesFromMidPoints(fInput, &protoJets);// Add midpoints
00115   if(theDebugLevel>=1)cout << "[CMSMidpointAlgorithm] Num proto-jets from Seeds and Midpoints = " << protoJets.size() << endl;
00116 
00117   // Split and merge the proto-jets, assigning each tower in the protojets to one and only one final jet.
00118   if(protoJets.size()>0)splitAndMerge(fInput, &protoJets, fOutput);    // Split and merge 
00119   if(theDebugLevel>=1)cout << "[CMSMidpointAlgorithm] Num final jets = " << fOutput->size() << endl;
00120 
00121   // Make the CaloJets from the final protojets
00122   //  MakeCaloJet(*theCtcp, finalJets, caloJets);
00123 
00124 }

void CMSMidpointAlgorithm::splitAndMerge ( const JetReco::InputCollection fInput,
InternalCollection fProtoJets,
JetReco::OutputCollection fFinalJets 
) [private]

Split and merge the proto-jets, assigning each tower in the protojets to one and only one final jet.

Called by run but public method to allow testing and studies.

Definition at line 311 of file CMSMidpointAlgorithm.cc.

References TestMuL1L2Filter_cff::cerr, GenMuonPlsPt100GeV_cfg::cout, deltaR2(), lat::endl(), find(), ProtoJet::getTowerList(), i, j, ProtoJet::p4(), ProtoJet::phi(), ProtoJet::pt(), ProtoJet::putTowers(), sameTower(), python::multivaluedict::sort(), theDebugLevel, theOverlapThreshold, and ProtoJet::y().

Referenced by run().

00312                                                                                                         {
00313   //
00314   // This can use quite a bit of optimization and simplification.
00315   //
00316   
00317   // Debugging
00318   if(theDebugLevel>=2){
00319     // Take a look at the proto-jets we were given
00320     int numProtojets = stableCones->size();
00321     for(int i = 0; i < numProtojets ; ++i){
00322       const ProtoJet* icone = (*stableCones)[i];
00323       int numTowers = icone->getTowerList().size();
00324       cout << endl << "[CMSMidpointAlgorithm] ProtoJet " << i << ": PT=" << (*stableCones)[i]->pt()
00325            << ", y="<< icone->y()
00326            << ", phi="<< icone->phi()  
00327            << ", ntow="<< numTowers << endl;
00328       ProtoJet::Constituents protojetTowers = icone->getTowerList(); 
00329       for(int j = 0; j < numTowers; ++j){
00330         cout << "[CMSMidpointAlgorithm] Tower " << j << ": ET=" << protojetTowers[j]->et() 
00331              << ", eta="<< protojetTowers[j]->eta()
00332              << ", phi="<<   protojetTowers[j]->phi() << endl;
00333       }
00334     }      
00335   }
00336   
00337   // Start of split and merge algorithm
00338   bool mergingNotFinished = true;
00339   while(mergingNotFinished){
00340     
00341     // Sort the stable cones (highest pt first).
00342     GreaterByPtPtr<ProtoJet> compJets;
00343     sort(stableCones->begin(),stableCones->end(),compJets);
00344     // clean removed clusters
00345     InternalCollection::const_iterator i = find (stableCones->begin(), stableCones->end(), (ProtoJet*)0);
00346     //    for (; i != stableCones->end(); ++i) std::cout << "CMSMidpointAlgorithm::splitAndMerge-> removing pointer " << *i << std::endl; 
00347     stableCones->erase (find (stableCones->begin(), stableCones->end(), (ProtoJet*)0), stableCones->end()); 
00348     
00349     // Start with the highest pt cone left in list. List changes in loop,
00350     // getting smaller with each iteration.
00351     InternalCollection::iterator stableConeIter1 = stableCones->begin();
00352     
00353     if(stableConeIter1 == stableCones->end()) {  // Stable cone list empty?
00354       mergingNotFinished = false;
00355     }
00356     else {
00357       ProtoJet* stableCone1 = *stableConeIter1;
00358       if (!stableCone1) {
00359         std::cerr << "CMSMidpointAlgorithm::splitAndMerge-> Error: stableCone1 should never be 0" << std::endl;
00360         continue;
00361       }
00362       bool coneNotModified = true;
00363       
00364       // Determine whether highest pt cone has an overlap with other stable cones.
00365       InternalCollection::iterator stableConeIter2 = stableConeIter1; 
00366       ++stableConeIter2;   // Iterator for 2nd highest pt cone, just past 1st cone.
00367       
00368       while(coneNotModified && stableConeIter2 != stableCones->end()){
00369         ProtoJet* stableCone2 = *stableConeIter2;
00370         if (!stableCone2) {
00371           std::cerr << "CMSMidpointAlgorithm::splitAndMerge-> Error: stableCone2 should never be 0" << std::endl;
00372           continue;
00373         }
00374         
00375         // Calculate overlap of the two cones.
00376         bool overlap = false;
00377         vector<InputItem> overlapTowers;  // Make a list to hold the overlap towers
00378         //cout << "1st cone num towers=" << stableCone1->getTowerList().size() << endl;
00379         //int numTowers1=0;
00380         
00381         //Loop over towers in higher Pt cone
00382         for(ProtoJet::Constituents::const_iterator towerIter1 = stableCone1->getTowerList().begin();
00383             towerIter1 != stableCone1->getTowerList().end();
00384             ++towerIter1){
00385           //cout << "1st cone tower " << numTowers1 << endl;
00386           //++numTowers1;
00387           bool isInCone2 = false;
00388           //cout << "2nd cone num towers=" << stableCone2->getTowerList().size() << endl;
00389           //int numTowers2=0;
00390           
00391           // Loop over towers in lower Pt cone
00392           for(ProtoJet::Constituents::const_iterator towerIter2 = stableCone2->getTowerList().begin();
00393               towerIter2 != stableCone2->getTowerList().end();
00394               ++towerIter2) {
00395             //cout << "2st cone tower " << numTowers2 << endl;
00396             //++numTowers2;
00397             
00398             // Check if towers are the same by checking for unique eta, phi and energy values.
00399             // Will want to replace this with checking the tower index when available.
00400             if(sameTower(*towerIter1, *towerIter2)) {
00401               isInCone2 = true;   //The tower is in both cones
00402               //cout << "Merging found overlap tower: eta=" << (*towerIter1)->eta << ", phi=" << (*towerIter1)->phi << " in 1st protojet "  << 
00403               //                             " and eta=" << (*towerIter2)->eta << ", phi=" << (*towerIter2)->phi << " in 2nd protojet "  << endl;
00404             }
00405           }
00406           if(isInCone2){
00407             overlapTowers.push_back(*towerIter1);  //Add tower in both cones to the overlap list
00408             overlap = true;
00409           }          
00410         }
00411         if(overlap){   
00412           // non-empty overlap.  Decide on splitting or merging.
00413           
00414           // Make a proto-jet with the overlap towers so we can calculate things for the overlap
00415           ProtoJet overlap;
00416           overlap.putTowers(overlapTowers);
00417           coneNotModified = false;
00418           
00419           // Compare the overlap pt with the overlap fractcion threshold times the lower jet pt.
00420           if(overlap.pt() >= theOverlapThreshold*stableCone2->pt()){
00421             
00422             // Merge the two cones.
00423             // Get a copy of the list of towers in higher Pt proto-jet 
00424             ProtoJet::Constituents stableCone1Towers = stableCone1->getTowerList(); 
00425             
00426             //Loop over the list of towers lower Pt jet
00427             for(ProtoJet::Constituents::const_iterator towerIter2 = stableCone2->getTowerList().begin();
00428                 towerIter2 != stableCone2->getTowerList().end();
00429                 ++towerIter2){
00430               bool isInOverlap = false;
00431               
00432               //Check if that tower is in the overlap region
00433               for(vector<InputItem>::iterator overlapTowerIter = overlapTowers.begin();
00434                   overlapTowerIter != overlapTowers.end();
00435                   ++overlapTowerIter){
00436                 // Check if towers are the same by checking for unique eta, phi and energy values.
00437                 // Will want to replace this with checking the tower index when available.
00438                 if(sameTower(*overlapTowerIter, *towerIter2))   isInOverlap = true;
00439               }
00440               //Add  non-overlap tower from proto-jet 2 into the proto-jet 1 tower list.
00441               if(!isInOverlap)stableCone1Towers.push_back(*towerIter2);  
00442             }
00443             
00444             if(theDebugLevel>=2)cout << endl << "[CMSMidpointAlgorithm] Merging: 1st Proto-jet grows: "
00445                                   " y=" << stableCone1->y() << 
00446                                   ", phi=" << stableCone1->phi() << 
00447                                   " increases from " << stableCone1->getTowerList().size() << 
00448                                   " to "  << stableCone1Towers.size() << " towers." << endl;
00449             
00450             // Put the new expanded list of towers into the first proto-jet
00451             stableCone1->putTowers(stableCone1Towers);
00452             
00453             if(theDebugLevel>=2)cout << "[CMSMidpointAlgorithm] Merging: 1st protojet now at y=" << stableCone1->y() <<
00454                                   ", phi=" << stableCone1->phi() << endl;
00455             
00456             if(theDebugLevel>=2)cout << "[CMSMidpointAlgorithm] Merging: 2nd Proto-jet removed:" 
00457                                   " y=" << stableCone2->y() << 
00458                                   ", phi=" << stableCone2->phi() << endl;
00459             
00460             // Remove the second proto-jet.
00461             delete *stableConeIter2;
00462             *stableConeIter2 = 0;
00463             
00464           }
00465           else{
00466             // Split the two proto-jets.
00467             
00468             // Create lists of towers to remove from each proto-jet
00469             vector<InputItem> removeFromCone1,removeFromCone2;
00470             
00471             // Which tower goes where?
00472             // Loop over the overlap towers
00473             for(vector<InputItem>::iterator towerIter = overlapTowers.begin();
00474                 towerIter != overlapTowers.end();
00475                 ++towerIter){
00476               double dR2Jet1 = deltaR2 ((*towerIter)->p4().Rapidity(), (*towerIter)->phi(), 
00477                                       stableCone1->y(), stableCone1->phi()); 
00478               // Calculate distance from proto-jet 2.
00479               double dR2Jet2 = deltaR2 ((*towerIter)->p4().Rapidity(), (*towerIter)->phi(),
00480                                       stableCone2->y(), stableCone2->phi()); 
00481               
00482               if(dR2Jet1 < dR2Jet2){
00483                 // Tower is closer to proto-jet 1. To be removed from proto-jet 2.
00484                 removeFromCone2.push_back(*towerIter);
00485               }
00486               else {
00487                 // Tower is closer to proto-jet 2. To be removed from proto-jet 1.
00488                 removeFromCone1.push_back(*towerIter);
00489               }
00490             }
00491             // Remove towers in the overlap region from the cones to which they have the larger distance.
00492             
00493             // Remove towers from proto-jet 1.
00494             vector<InputItem> towerList1 (stableCone1->getTowerList().begin(), stableCone1->getTowerList().end()); 
00495 
00496             // Loop over towers in remove list
00497             for(vector<InputItem>::iterator towerIter = removeFromCone1.begin();
00498                 towerIter != removeFromCone1.end();
00499                 ++towerIter) {
00500               // Loop over towers in protojet
00501               for(vector<InputItem>::iterator towerIter1 = towerList1.begin(); towerIter1 != towerList1.end(); ++towerIter1) {
00502                 // Check if they are equal
00503                 if(sameTower(*towerIter, *towerIter1)) {
00504                   // Remove the tower
00505                   towerList1.erase(towerIter1);
00506                   break;
00507                 }
00508               }
00509             }
00510             
00511             if(theDebugLevel>=2)cout << endl << "[CMSMidpointAlgorithm] Splitting: 1st Proto-jet  shrinks: y=" << 
00512                                   stableCone1->y() << 
00513                                   ", phi=" << stableCone1->phi() << 
00514                                   " decreases from" << stableCone1->getTowerList().size() << 
00515                                   " to "  << towerList1.size() << " towers." << endl;
00516             
00517             //Put the new reduced list of towers into proto-jet 1.
00518             stableCone1->putTowers(towerList1); 
00519             // Remove towers from cone 2.
00520             vector<InputItem> towerList2 (stableCone2->getTowerList().begin(), stableCone2->getTowerList().end()); 
00521             
00522             // Loop over towers in remove list
00523             for(vector<InputItem>::iterator towerIter = removeFromCone2.begin();
00524                 towerIter != removeFromCone2.end();
00525                 ++towerIter) {
00526               // Loop over towers in protojet
00527               for(vector<InputItem>::iterator towerIter2 = towerList2.begin(); towerIter2 != towerList2.end(); ++towerIter2){
00528                 // Check if they are equal
00529                 if(sameTower(*towerIter, *towerIter2)) {
00530                   // Remove the tower
00531                   towerList2.erase(towerIter2);
00532                   break;
00533                 }
00534               }
00535             }
00536             
00537             if(theDebugLevel>=2)cout << "[CMSMidpointAlgorithm] Splitting: 2nd Proto-jet shrinks: y=" <<
00538                                   stableCone2->y() << 
00539                                   ", phi=" << stableCone2->phi() << 
00540                                   " decreases from" << stableCone2->getTowerList().size() << 
00541                                   " to "  << towerList2.size() << " towers." << endl;
00542             
00543             //Put the new reduced list of towers into proto-jet 2.
00544             stableCone2->putTowers(towerList2); 
00545           }
00546         }
00547         else {
00548           if(theDebugLevel>=2)cout << endl << 
00549                                 "[CMSMidpointAlgorithm] no overlap between 1st protojet at  y=" << stableCone1->y() << 
00550                                 ", phi=" << stableCone1->phi() << 
00551                                 " and 2nd protojet at  y=" << stableCone2->y() << 
00552                                 ", phi=" << stableCone2->phi() <<endl;
00553         }
00554         
00555         ++stableConeIter2;  //Increment iterator to the next highest Pt protojet
00556       }
00557       if(coneNotModified){
00558         
00559         if(theDebugLevel>=2)cout << 
00560                               "[CMSMidpointAlgorithm] Saving: Proto-jet  at y=" << stableCone1->y() << 
00561                               ", phi=" << stableCone1->phi() <<  " has no overlap" << endl;
00562         
00563         fFinalJets->push_back(ProtoJet (stableCone1->p4(), stableCone1->getTowerList()));
00564         delete *stableConeIter1;
00565         *stableConeIter1 = 0;
00566       }
00567     }
00568   }
00569   
00570   GreaterByPt<ProtoJet> compJets;
00571   sort(fFinalJets->begin(),fFinalJets->end(),compJets);
00572 }


Member Data Documentation

double CMSMidpointAlgorithm::theConeAreaFraction [private]

Definition at line 109 of file CMSMidpointAlgorithm.h.

Referenced by iterateCone().

double CMSMidpointAlgorithm::theConeRadius [private]

Definition at line 108 of file CMSMidpointAlgorithm.h.

Referenced by findStableConesFromMidPoints(), and iterateCone().

int CMSMidpointAlgorithm::theDebugLevel [private]

Definition at line 113 of file CMSMidpointAlgorithm.h.

Referenced by findStableConesFromMidPoints(), findStableConesFromSeeds(), iterateCone(), run(), and splitAndMerge().

int CMSMidpointAlgorithm::theMaxIterations [private]

Definition at line 111 of file CMSMidpointAlgorithm.h.

Referenced by iterateCone().

int CMSMidpointAlgorithm::theMaxPairSize [private]

Definition at line 110 of file CMSMidpointAlgorithm.h.

Referenced by findStableConesFromMidPoints().

double CMSMidpointAlgorithm::theOverlapThreshold [private]

Definition at line 112 of file CMSMidpointAlgorithm.h.

Referenced by splitAndMerge().

double CMSMidpointAlgorithm::theSeedThreshold [private]

Definition at line 107 of file CMSMidpointAlgorithm.h.

Referenced by findStableConesFromSeeds().


The documentation for this class was generated from the following files:
Generated on Tue Jun 9 18:16:23 2009 for CMSSW by  doxygen 1.5.4