Main Page | Namespace List | Class Hierarchy | Class List | Directories | File List | Namespace Members | Class Members | File Members

PseudoJet.cc File Reference

#include "fastjet/Error.hh"
#include "fastjet/PseudoJet.hh"
#include <valarray>
#include <iostream>
#include <sstream>
#include <cmath>
#include <algorithm>

Include dependency graph for PseudoJet.cc:

Go to the source code of this file.

Namespaces

namespace  fastjet

Functions

PseudoJet operator+ (const PseudoJet &jet1, const PseudoJet &jet2)
PseudoJet operator- (const PseudoJet &jet1, const PseudoJet &jet2)
PseudoJet operator * (double coeff, const PseudoJet &jet)
PseudoJet operator * (const PseudoJet &jet, double coeff)
PseudoJet operator/ (const PseudoJet &jet, double coeff)
bool have_same_momentum (const PseudoJet &jeta, const PseudoJet &jetb)
 returns true if the momenta of the two input jets are identical
PseudoJet PtYPhiM (double pt, double y, double phi, double m)
 return a pseudojet with the given pt, y, phi and mass
void sort_indices (vector< int > &indices, const vector< double > &values)
template<class T>
vector< T > objects_sorted_by_values (const vector< T > &objects, const vector< double > &values)
 given a vector of values with a one-to-one correspondence with the vector of objects, sort objects into an order such that the associated values would be in increasing order
vector< PseudoJet > sorted_by_pt (const vector< PseudoJet > &jets)
 return a vector of jets sorted into decreasing kt2
vector< PseudoJet > sorted_by_rapidity (const vector< PseudoJet > &jets)
 return a vector of jets sorted into increasing rapidity
vector< PseudoJet > sorted_by_E (const vector< PseudoJet > &jets)
 return a vector of jets sorted into decreasing energy


Function Documentation

bool fastjet::have_same_momentum const PseudoJet &  jeta,
const PseudoJet &  jetb
 

returns true if the momenta of the two input jets are identical

Definition at line 284 of file PseudoJet.cc.

References fastjet::PseudoJet::E(), fastjet::PseudoJet::px(), fastjet::PseudoJet::py(), and fastjet::PseudoJet::pz().

Referenced by fastjet::SISConePlugin::run_clustering().

00284                                                                         {
00285   return jeta.px() == jetb.px()
00286     &&   jeta.py() == jetb.py()
00287     &&   jeta.pz() == jetb.pz()
00288     &&   jeta.E()  == jetb.E();
00289 }

template<class T>
vector<T> fastjet::objects_sorted_by_values const vector< T > &  objects,
const vector< double > &  values
 

given a vector of values with a one-to-one correspondence with the vector of objects, sort objects into an order such that the associated values would be in increasing order

Definition at line 335 of file PseudoJet.cc.

References fastjet::sort_indices().

Referenced by fastjet::sorted_by_E(), fastjet::sorted_by_pt(), and fastjet::sorted_by_rapidity().

00337                                                       {
00338 
00339   assert(objects.size() == values.size());
00340 
00341   // get a vector of indices
00342   vector<int> indices(values.size());
00343   for (size_t i = 0; i < indices.size(); i++) {indices[i] = i;}
00344   
00345   // sort the indices
00346   sort_indices(indices, values);
00347   
00348   // copy the objects 
00349   vector<T> objects_sorted(objects.size());
00350   
00351   // place the objects in the correct order
00352   for (size_t i = 0; i < indices.size(); i++) {
00353     objects_sorted[i] = objects[indices[i]];
00354   }
00355 
00356   return objects_sorted;
00357 }

PseudoJet fastjet::operator * const PseudoJet &  jet,
double  coeff
 

Definition at line 178 of file PseudoJet.cc.

00178                                                           {
00179   return coeff*jet;
00180 } 

PseudoJet fastjet::operator * double  coeff,
const PseudoJet &  jet
 

Definition at line 168 of file PseudoJet.cc.

00168                                                           {
00169   //return PseudoJet(coeff*jet.four_mom());
00170   // the following code is hopefully more efficient
00171   PseudoJet coeff_times_jet(jet);
00172   coeff_times_jet *= coeff;
00173   return coeff_times_jet;
00174 } 

PseudoJet fastjet::operator+ const PseudoJet &  jet1,
const PseudoJet &  jet2
 

Definition at line 148 of file PseudoJet.cc.

References fastjet::PseudoJet::E(), fastjet::PseudoJet::px(), fastjet::PseudoJet::py(), and fastjet::PseudoJet::pz().

00148                                                                      {
00149   //return PseudoJet(jet1.four_mom()+jet2.four_mom());
00150   return PseudoJet(jet1.px()+jet2.px(),
00151                    jet1.py()+jet2.py(),
00152                    jet1.pz()+jet2.pz(),
00153                    jet1.E() +jet2.E()  );
00154 } 

PseudoJet fastjet::operator- const PseudoJet &  jet1,
const PseudoJet &  jet2
 

Definition at line 158 of file PseudoJet.cc.

References fastjet::PseudoJet::E(), fastjet::PseudoJet::px(), fastjet::PseudoJet::py(), and fastjet::PseudoJet::pz().

00158                                                                      {
00159   //return PseudoJet(jet1.four_mom()-jet2.four_mom());
00160   return PseudoJet(jet1.px()-jet2.px(),
00161                    jet1.py()-jet2.py(),
00162                    jet1.pz()-jet2.pz(),
00163                    jet1.E() -jet2.E()  );
00164 } 

PseudoJet fastjet::operator/ const PseudoJet &  jet,
double  coeff
 

Definition at line 184 of file PseudoJet.cc.

00184                                                           {
00185   return (1.0/coeff)*jet;
00186 } 

PseudoJet fastjet::PtYPhiM double  pt,
double  y,
double  phi,
double  m = 0.0
 

return a pseudojet with the given pt, y, phi and mass

Definition at line 294 of file PseudoJet.cc.

00294                                                              {
00295   double ptm = sqrt(pt*pt+m*m);
00296   return PseudoJet(pt*cos(phi), pt*sin(phi), ptm*sinh(y), ptm*cosh(y));
00297 }

void fastjet::sort_indices vector< int > &  indices,
const vector< double > &  values
 

Definition at line 325 of file PseudoJet.cc.

Referenced by fastjet::objects_sorted_by_values().

00326                                                         {
00327   IndexedSortHelper index_sort_helper(&values);
00328   sort(indices.begin(), indices.end(), index_sort_helper);
00329 }

vector<PseudoJet> fastjet::sorted_by_E const vector< PseudoJet > &  jets  ) 
 

return a vector of jets sorted into decreasing energy

Definition at line 377 of file PseudoJet.cc.

References fastjet::objects_sorted_by_values().

Referenced by main().

00377                                                               {
00378   vector<double> energies(jets.size());
00379   for (size_t i = 0; i < jets.size(); i++) {energies[i] = -jets[i].E();}
00380   return objects_sorted_by_values(jets, energies);
00381 }

vector<PseudoJet> fastjet::sorted_by_pt const vector< PseudoJet > &  jets  ) 
 

return a vector of jets sorted into decreasing kt2

Definition at line 361 of file PseudoJet.cc.

References fastjet::objects_sorted_by_values().

Referenced by main(), print_jets(), and fastjet::ClusterSequenceAreaBase::subtracted_jets().

00361                                                                {
00362   vector<double> minus_kt2(jets.size());
00363   for (size_t i = 0; i < jets.size(); i++) {minus_kt2[i] = -jets[i].kt2();}
00364   return objects_sorted_by_values(jets, minus_kt2);
00365 }

vector<PseudoJet> fastjet::sorted_by_rapidity const vector< PseudoJet > &  jets  ) 
 

return a vector of jets sorted into increasing rapidity

Definition at line 369 of file PseudoJet.cc.

References fastjet::objects_sorted_by_values().

Referenced by main().

00369                                                                      {
00370   vector<double> rapidities(jets.size());
00371   for (size_t i = 0; i < jets.size(); i++) {rapidities[i] = jets[i].rap();}
00372   return objects_sorted_by_values(jets, rapidities);
00373 }


Generated on Fri Aug 15 13:45:44 2008 for fastjet by  doxygen 1.4.2