fastjet::PseudoJet Class Reference

Class to contain pseudojets, including minimal information of use to to jet-clustering routines. More...

#include <PseudoJet.hh>

Inheritance diagram for fastjet::PseudoJet:

Inheritance graph
fastjet::ClusterSequenceActiveArea::GhostJet
[legend]
List of all members.

Public Types

enum  {
  X = 0, Y = 1, Z = 2, T = 3,
  NUM_COORDINATES = 4, SIZE = NUM_COORDINATES
}

Public Member Functions

 PseudoJet ()
 PseudoJet (const double px, const double py, const double pz, const double E)
 construct a pseudojet from explicit components
template<class L>
 PseudoJet (const L &some_four_vector)
 constructor from any object that has px,py,pz,E = some_four_vector[0--3],
double E () const
double e () const
double px () const
double py () const
double pz () const
const double phi () const
 returns phi (in the range 0..2pi)
const double phi_std () const
 returns phi in the range -pi..pi
const double phi_02pi () const
 returns phi in the range 0..2pi
double rap () const
 returns the rapidity or some large value when the rapidity is infinite
double rapidity () const
 the same as rap()
double pseudorapidity () const
 returns the pseudo-rapidity or some large value when the rapidity is infinite
double eta () const
double kt2 () const
 returns the squared transverse momentum
double perp2 () const
 returns the squared transverse momentum
double perp () const
 returns the scalar transverse momentum
double m2 () const
 returns the squared invariant mass // like CLHEP
double mperp2 () const
 returns the squared transverse mass = kt^2+m^2
double mperp () const
 returns the transverse mass = sqrt(kt^2+m^2)
double m () const
 specialization of the "reset" template for case where something is reset to a pseudojet -- it then takes the user and history indices from the psjet
double operator() (int i) const
 returns component i, where X==0, Y==1, Z==2, E==3
double operator[] (int i) const
 returns component i, where X==0, Y==1, Z==2, E==3
PseudoJetboost (const PseudoJet &prest)
 transform this jet (given in the rest frame of prest) into a jet in the lab frame;
PseudoJetunboost (const PseudoJet &prest)
 transform this jet (given in the rest frame of prest) into a jet in the lab frame;
const int & cluster_hist_index () const
 return the cluster_hist_index, intended to be used by clustering routines.
void set_cluster_hist_index (const int index)
 set the cluster_hist_index, intended to be used by clustering routines.
const int cluster_sequence_history_index () const
 alternative name for cluster_hist_index() [perhaps more meaningful]
void set_cluster_sequence_history_index (const int index)
 alternative name for set_cluster_hist_index(.
const int & user_index () const
 return the user_index, intended to allow the user to "add" information
void set_user_index (const int index)
 set the user_index, intended to allow the user to "add" information
std::valarray< double > four_mom () const
 return a valarray containing the four-momentum (components 0-2 are 3-mom, component 3 is energy).
double kt_distance (const PseudoJet &other) const
 returns kt distance (R=1) between this jet and another
double plain_distance (const PseudoJet &other) const
 returns squared cylinder (rap-phi) distance between this jet and another
double squared_distance (const PseudoJet &other) const
 returns squared cylinder (rap-phi) distance between this jet and another
double beam_distance () const
 returns distance between this jet and the beam
void operator *= (double)
 multiply the jet's momentum by the coefficient
void operator/= (double)
 divide the jet's momentum by the coefficient
void operator+= (const PseudoJet &)
 add the other jet's momentum to this jet
void operator-= (const PseudoJet &)
 subtract the other jet's momentum from this jet
void reset (double px, double py, double pz, double E)
 reset the 4-momentum according to the supplied components and put the user and history indices back to their default values
void reset (const PseudoJet &psjet)
 reset the PseudoJet to be equal to psjet (including its indices); NB if the argument is derived from a PseudoJet then the "reset" used will be the templated version (which does not know about indices.
template<class L>
void reset (const L &some_four_vector)
 reset the 4-momentum according to the supplied generic 4-vector (accessible via indexing, [0]==px,.
template<>
 PseudoJet (const Cmomentum &four_vector)
 shortcut for converting siscone Cmomentum into PseudoJet

Private Member Functions

void _finish_init ()
 do standard end of initialisation
void _reset_indices ()
 set the indices to default values

Private Attributes

double _px
double _py
double _pz
double _E
double _phi
double _rap
double _kt2
int _cluster_hist_index
int _user_index

Detailed Description

Class to contain pseudojets, including minimal information of use to to jet-clustering routines.

Definition at line 52 of file PseudoJet.hh.


Member Enumeration Documentation

anonymous enum

Enumerator:
X 
Y 
Z 
T 
NUM_COORDINATES 
SIZE 

Definition at line 115 of file PseudoJet.hh.

00115 { X=0, Y=1, Z=2, T=3, NUM_COORDINATES=4, SIZE=NUM_COORDINATES };


Constructor & Destructor Documentation

fastjet::PseudoJet::PseudoJet (  )  [inline]

Definition at line 55 of file PseudoJet.hh.

Referenced by fastjet::operator+(), fastjet::operator-(), PseudoJet(), and fastjet::PtYPhiM().

00055 {};

fastjet::PseudoJet::PseudoJet ( const double  px,
const double  py,
const double  pz,
const double  E 
)

construct a pseudojet from explicit components

Definition at line 46 of file PseudoJet.cc.

References _E, _finish_init(), _px, _py, _pz, and _reset_indices().

00046                                                                                       {
00047   
00048   _E  = E ;
00049   _px = px;
00050   _py = py;
00051   _pz = pz;
00052 
00053   this->_finish_init();
00054 
00055   // some default values for the history and user indices
00056   _reset_indices();
00057 
00058 }

template<class L>
fastjet::PseudoJet::PseudoJet ( const L &  some_four_vector  )  [inline]

constructor from any object that has px,py,pz,E = some_four_vector[0--3],

Definition at line 267 of file PseudoJet.hh.

References _E, _finish_init(), _px, _py, _pz, and _reset_indices().

00267                                                                           {
00268 
00269   _px = some_four_vector[0];
00270   _py = some_four_vector[1];
00271   _pz = some_four_vector[2];
00272   _E  = some_four_vector[3];
00273   _finish_init();
00274   // some default values for these two indices
00275   _reset_indices();
00276 }

template<>
fastjet::PseudoJet::PseudoJet ( const Cmomentum &  four_vector  )  [inline]

shortcut for converting siscone Cmomentum into PseudoJet

Definition at line 55 of file SISConePlugin.cc.

References PseudoJet().

00055                                                              {
00056   (*this) = PseudoJet(four_vector.px,four_vector.py,four_vector.pz,
00057                       four_vector.E);
00058 }


Member Function Documentation

double fastjet::PseudoJet::E (  )  const [inline]

Definition at line 66 of file PseudoJet.hh.

Referenced by fastjet::ClusterSequenceActiveArea::_throw_unless_jets_have_same_perp_or_E(), boost(), fastjet::have_same_momentum(), fastjet::operator+(), fastjet::operator-(), fastjet::JetDefinition::DefaultRecombiner::preprocess(), fastjet::sorted_by_E(), and unboost().

00066 {return _E;}

double fastjet::PseudoJet::e (  )  const [inline]

Definition at line 67 of file PseudoJet.hh.

Referenced by operator()().

00067 {return _E;} // like CLHEP

double fastjet::PseudoJet::px (  )  const [inline]

Definition at line 68 of file PseudoJet.hh.

Referenced by _finish_init(), fastjet::ClusterSequenceActiveArea::_throw_unless_jets_have_same_perp_or_E(), boost(), fastjet::have_same_momentum(), operator()(), fastjet::operator+(), fastjet::operator-(), fastjet::JetDefinition::DefaultRecombiner::preprocess(), pseudorapidity(), and unboost().

00068 {return _px;}

double fastjet::PseudoJet::py (  )  const [inline]

Definition at line 69 of file PseudoJet.hh.

Referenced by _finish_init(), fastjet::ClusterSequenceActiveArea::_throw_unless_jets_have_same_perp_or_E(), boost(), fastjet::have_same_momentum(), operator()(), fastjet::operator+(), fastjet::operator-(), fastjet::JetDefinition::DefaultRecombiner::preprocess(), pseudorapidity(), and unboost().

00069 {return _py;}

double fastjet::PseudoJet::pz (  )  const [inline]

Definition at line 70 of file PseudoJet.hh.

Referenced by fastjet::ClusterSequenceActiveArea::_throw_unless_jets_have_same_perp_or_E(), boost(), fastjet::have_same_momentum(), operator()(), fastjet::operator+(), fastjet::operator-(), fastjet::JetDefinition::DefaultRecombiner::preprocess(), pseudorapidity(), and unboost().

00070 {return _pz;}

const double fastjet::PseudoJet::phi (  )  const [inline]

returns phi (in the range 0..2pi)

Definition at line 73 of file PseudoJet.hh.

Referenced by fastjet::RangeDefinition::is_in_range(), and fastjet::JetDefinition::DefaultRecombiner::recombine().

00073 {return phi_02pi();}

const double fastjet::PseudoJet::phi_std (  )  const [inline]

returns phi in the range -pi..pi

Definition at line 76 of file PseudoJet.hh.

References fastjet::pi, and twopi.

00076                                        {
00077     return _phi > pi ? _phi-twopi : _phi;}

const double fastjet::PseudoJet::phi_02pi (  )  const [inline]

returns phi in the range 0..2pi

Definition at line 80 of file PseudoJet.hh.

00080 {return _phi;}

double fastjet::PseudoJet::rap (  )  const [inline]

returns the rapidity or some large value when the rapidity is infinite

Definition at line 84 of file PseudoJet.hh.

Referenced by fastjet::RangeDefinition::is_in_range(), main(), fastjet::JetDefinition::DefaultRecombiner::recombine(), and fastjet::sorted_by_rapidity().

00084 {return _rap;}

double fastjet::PseudoJet::rapidity (  )  const [inline]

the same as rap()

Definition at line 87 of file PseudoJet.hh.

00087 {return _rap;} // like CLHEP

double fastjet::PseudoJet::pseudorapidity (  )  const

returns the pseudo-rapidity or some large value when the rapidity is infinite

Definition at line 136 of file PseudoJet.cc.

References fastjet::MaxRap, perp(), fastjet::pi, px(), py(), and pz().

00136                                        {
00137   if (px() == 0.0 && py() ==0.0) return MaxRap;
00138   if (pz() == 0.0) return 0.0;
00139 
00140   double theta = atan(perp()/pz());
00141   if (theta < 0) theta += pi;
00142   return -log(tan(theta/2));
00143 }

double fastjet::PseudoJet::eta (  )  const [inline]

Definition at line 92 of file PseudoJet.hh.

00092 {return pseudorapidity();}

double fastjet::PseudoJet::kt2 (  )  const [inline]

returns the squared transverse momentum

Definition at line 95 of file PseudoJet.hh.

Referenced by fastjet::ClusterSequence::jet_scale_for_algorithm(), and fastjet::sorted_by_pt().

00095 {return _kt2;}

double fastjet::PseudoJet::perp2 (  )  const [inline]

returns the squared transverse momentum

Definition at line 97 of file PseudoJet.hh.

Referenced by fastjet::ClusterSequenceActiveArea::_throw_unless_jets_have_same_perp_or_E(), fastjet::ClusterSequence::has_parents(), fastjet::ClusterSequence::inclusive_jets(), fastjet::JetDefinition::DefaultRecombiner::preprocess(), and fastjet::JetDefinition::DefaultRecombiner::recombine().

00097 {return _kt2;}  // like CLHEP

double fastjet::PseudoJet::perp (  )  const [inline]

returns the scalar transverse momentum

Definition at line 99 of file PseudoJet.hh.

Referenced by fastjet::ClusterSequenceAreaBase::get_median_rho_and_sigma(), fastjet::ClusterSequenceAreaBase::parabolic_pt_per_unit_area(), fastjet::ClusterSequenceActiveArea::parabolic_pt_per_unit_area(), pseudorapidity(), fastjet::ClusterSequenceActiveArea::pt_per_unit_area(), fastjet::JetDefinition::DefaultRecombiner::recombine(), fastjet::ClusterSequenceAreaBase::subtracted_jet(), and fastjet::ClusterSequenceAreaBase::subtracted_pt().

00099 {return sqrt(_kt2);}    // like CLHEP

double fastjet::PseudoJet::m2 (  )  const [inline]

returns the squared invariant mass // like CLHEP

Definition at line 101 of file PseudoJet.hh.

Referenced by _finish_init(), and m().

00101 {return (_E+_pz)*(_E-_pz)-_kt2;}    

double fastjet::PseudoJet::mperp2 (  )  const [inline]

returns the squared transverse mass = kt^2+m^2

Definition at line 103 of file PseudoJet.hh.

00103 {return (_E+_pz)*(_E-_pz);}

double fastjet::PseudoJet::mperp (  )  const [inline]

returns the transverse mass = sqrt(kt^2+m^2)

Definition at line 105 of file PseudoJet.hh.

00105 {return sqrt(std::abs(mperp2()));}

double fastjet::PseudoJet::m (  )  const [inline]

specialization of the "reset" template for case where something is reset to a pseudojet -- it then takes the user and history indices from the psjet

Definition at line 302 of file PseudoJet.hh.

References m2().

Referenced by boost(), and unboost().

00302                                  {
00303   double mm = m2();
00304   return mm < 0.0 ? -std::sqrt(-mm) : std::sqrt(mm);
00305 }

double fastjet::PseudoJet::operator() ( int  i  )  const

returns component i, where X==0, Y==1, Z==2, E==3

Definition at line 116 of file PseudoJet.cc.

References e(), px(), py(), pz(), T, X, Y, and Z.

00116                                           {
00117   switch(i) {
00118   case X:
00119     return px();
00120   case Y:
00121     return py();
00122   case Z:
00123     return pz();
00124   case T:
00125     return e();
00126   default:
00127     ostringstream err;
00128     err << "PseudoJet subscripting: bad index (" << i << ")";
00129     throw Error(err.str());
00130   }
00131   return 0.;
00132 }  

double fastjet::PseudoJet::operator[] ( int  i  )  const [inline]

returns component i, where X==0, Y==1, Z==2, E==3

Definition at line 112 of file PseudoJet.hh.

00112 { return (*this)(i); }; // this too

PseudoJet & fastjet::PseudoJet::boost ( const PseudoJet prest  ) 

transform this jet (given in the rest frame of prest) into a jet in the lab frame;

Definition at line 233 of file PseudoJet.cc.

References _E, _finish_init(), _px, _py, _pz, E(), m(), px(), py(), and pz().

00233                                                     {
00234   
00235   if (prest.px() == 0.0 && prest.py() == 0.0 && prest.pz() == 0.0) 
00236     return *this;
00237 
00238   double m = prest.m();
00239   assert(m != 0);
00240 
00241   double pf4  = (  px()*prest.px() + py()*prest.py()
00242                  + pz()*prest.pz() + E()*prest.E() )/m;
00243   double fn   = (pf4 + E()) / (prest.E() + m);
00244   _px +=  fn*prest.px();
00245   _py +=  fn*prest.py();
00246   _pz +=  fn*prest.pz();
00247   _E = pf4;
00248 
00249   _finish_init(); // we need to recalculate phi,rap,kt2
00250   return *this;
00251 }

PseudoJet & fastjet::PseudoJet::unboost ( const PseudoJet prest  ) 

transform this jet (given in the rest frame of prest) into a jet in the lab frame;

Definition at line 260 of file PseudoJet.cc.

References _E, _finish_init(), _px, _py, _pz, E(), m(), px(), py(), and pz().

00260                                                       {
00261   
00262   if (prest.px() == 0.0 && prest.py() == 0.0 && prest.pz() == 0.0) 
00263     return *this;
00264 
00265   double m = prest.m();
00266   assert(m != 0);
00267 
00268   double pf4  = ( -px()*prest.px() - py()*prest.py()
00269                  - pz()*prest.pz() + E()*prest.E() )/m;
00270   double fn   = (pf4 + E()) / (prest.E() + m);
00271   _px -=  fn*prest.px();
00272   _py -=  fn*prest.py();
00273   _pz -=  fn*prest.pz();
00274   _E = pf4;
00275 
00276   _finish_init(); // we need to recalculate phi,rap,kt2
00277   return *this;
00278 }

const int& fastjet::PseudoJet::cluster_hist_index (  )  const [inline]

return the cluster_hist_index, intended to be used by clustering routines.

Definition at line 127 of file PseudoJet.hh.

Referenced by fastjet::ClusterSequence::_potentially_valid(), fastjet::ClusterSequenceActiveAreaExplicitGhosts::area(), fastjet::ClusterSequenceVoronoiArea::area(), fastjet::ClusterSequenceActiveArea::area(), fastjet::ClusterSequenceActiveAreaExplicitGhosts::area_4vector(), fastjet::ClusterSequenceVoronoiArea::area_4vector(), fastjet::ClusterSequenceActiveArea::area_4vector(), fastjet::ClusterSequenceActiveArea::area_error(), fastjet::ClusterSequence::has_child(), fastjet::ClusterSequence::has_parents(), fastjet::ClusterSequence::has_partner(), fastjet::ClusterSequenceActiveAreaExplicitGhosts::is_pure_ghost(), and fastjet::ClusterSequence::object_in_jet().

00127 {return _cluster_hist_index;}

void fastjet::PseudoJet::set_cluster_hist_index ( const int  index  )  [inline]

set the cluster_hist_index, intended to be used by clustering routines.

Definition at line 129 of file PseudoJet.hh.

Referenced by _reset_indices().

00129 {_cluster_hist_index = index;}

const int fastjet::PseudoJet::cluster_sequence_history_index (  )  const [inline]

alternative name for cluster_hist_index() [perhaps more meaningful]

Definition at line 132 of file PseudoJet.hh.

00132                                                           {
00133     return cluster_hist_index();}

void fastjet::PseudoJet::set_cluster_sequence_history_index ( const int  index  )  [inline]

alternative name for set_cluster_hist_index(.

..) [perhaps more meaningful]

Definition at line 136 of file PseudoJet.hh.

00136                                                                   {
00137     set_cluster_hist_index(index);}

const int& fastjet::PseudoJet::user_index (  )  const [inline]

return the user_index, intended to allow the user to "add" information

Definition at line 141 of file PseudoJet.hh.

Referenced by fastjet::JetDefinition::DefaultRecombiner::preprocess().

00141 {return _user_index;}

void fastjet::PseudoJet::set_user_index ( const int  index  )  [inline]

set the user_index, intended to allow the user to "add" information

Definition at line 143 of file PseudoJet.hh.

Referenced by _reset_indices(), main(), fastjet::JetDefinition::DefaultRecombiner::preprocess(), and fastjet::JetDefinition::DefaultRecombiner::recombine().

00143 {_user_index = index;}

valarray< double > fastjet::PseudoJet::four_mom (  )  const

return a valarray containing the four-momentum (components 0-2 are 3-mom, component 3 is energy).

Definition at line 104 of file PseudoJet.cc.

References _E, _px, _py, and _pz.

00104                                            {
00105   valarray<double> mom(4);
00106   mom[0] = _px;
00107   mom[1] = _py;
00108   mom[2] = _pz;
00109   mom[3] = _E ;
00110   return mom;
00111 }

double fastjet::PseudoJet::kt_distance ( const PseudoJet other  )  const

returns kt distance (R=1) between this jet and another

Definition at line 301 of file PseudoJet.cc.

References _kt2, _phi, _rap, fastjet::pi, and fastjet::twopi.

00301                                                            {
00302   //double distance = min(this->kt2(), other.kt2());
00303   double distance = min(_kt2, other._kt2);
00304   double dphi = abs(_phi - other._phi);
00305   if (dphi > pi) {dphi = twopi - dphi;}
00306   double drap = _rap - other._rap;
00307   distance = distance * (dphi*dphi + drap*drap);
00308   return distance;
00309 }

double fastjet::PseudoJet::plain_distance ( const PseudoJet other  )  const

returns squared cylinder (rap-phi) distance between this jet and another

Definition at line 314 of file PseudoJet.cc.

References _phi, _rap, fastjet::pi, and fastjet::twopi.

00314                                                               {
00315   double dphi = abs(_phi - other._phi);
00316   if (dphi > pi) {dphi = twopi - dphi;}
00317   double drap = _rap - other._rap;
00318   return (dphi*dphi + drap*drap);
00319 }

double fastjet::PseudoJet::squared_distance ( const PseudoJet other  )  const [inline]

returns squared cylinder (rap-phi) distance between this jet and another

Definition at line 156 of file PseudoJet.hh.

00156                                                                 {
00157     return plain_distance(other);}

double fastjet::PseudoJet::beam_distance (  )  const [inline]

returns distance between this jet and the beam

Definition at line 165 of file PseudoJet.hh.

00165 {return _kt2;}

void fastjet::PseudoJet::operator *= ( double   ) 

multiply the jet's momentum by the coefficient

Definition at line 189 of file PseudoJet.cc.

References _E, _kt2, _px, _py, and _pz.

00189                                        {
00190   _px *= coeff;
00191   _py *= coeff;
00192   _pz *= coeff;
00193   _E  *= coeff;
00194   _kt2*= coeff*coeff;
00195   // phi and rap are unchanged
00196 }

void fastjet::PseudoJet::operator/= ( double   ) 

divide the jet's momentum by the coefficient

Definition at line 200 of file PseudoJet.cc.

00200                                        {
00201   (*this) *= 1.0/coeff;
00202 }

void fastjet::PseudoJet::operator+= ( const PseudoJet  ) 

add the other jet's momentum to this jet

Definition at line 207 of file PseudoJet.cc.

References _E, _finish_init(), _px, _py, and _pz.

00207                                                       {
00208   _px += other_jet._px;
00209   _py += other_jet._py;
00210   _pz += other_jet._pz;
00211   _E  += other_jet._E ;
00212   _finish_init(); // we need to recalculate phi,rap,kt2
00213 }

void fastjet::PseudoJet::operator-= ( const PseudoJet  ) 

subtract the other jet's momentum from this jet

Definition at line 218 of file PseudoJet.cc.

References _E, _finish_init(), _px, _py, and _pz.

00218                                                       {
00219   _px -= other_jet._px;
00220   _py -= other_jet._py;
00221   _pz -= other_jet._pz;
00222   _E  -= other_jet._E ;
00223   _finish_init(); // we need to recalculate phi,rap,kt2
00224 }

void fastjet::PseudoJet::reset ( double  px,
double  py,
double  pz,
double  E 
) [inline]

reset the 4-momentum according to the supplied components and put the user and history indices back to their default values

Definition at line 308 of file PseudoJet.hh.

References _E, _finish_init(), _px, _py, _pz, and _reset_indices().

00308                                                                       {
00309   _px = px;
00310   _py = py;
00311   _pz = pz;
00312   _E  = E;
00313   _finish_init();
00314   _reset_indices();
00315 }

void fastjet::PseudoJet::reset ( const PseudoJet psjet  )  [inline]

reset the PseudoJet to be equal to psjet (including its indices); NB if the argument is derived from a PseudoJet then the "reset" used will be the templated version (which does not know about indices.

..)

Definition at line 181 of file PseudoJet.hh.

00181                                              {
00182     (*this) = psjet;
00183   }

template<class L>
void fastjet::PseudoJet::reset ( const L &  some_four_vector  )  [inline]

reset the 4-momentum according to the supplied generic 4-vector (accessible via indexing, [0]==px,.

..[3]==E) and put the user and history indices back to their default values.

Definition at line 188 of file PseudoJet.hh.

00188                                                                    {
00189     reset(some_four_vector[0], some_four_vector[1],
00190           some_four_vector[2], some_four_vector[3]);
00191   }

void fastjet::PseudoJet::_finish_init (  )  [private]

do standard end of initialisation

Definition at line 63 of file PseudoJet.cc.

References _E, _kt2, _phi, _pz, _rap, m2(), fastjet::MaxRap, px(), py(), and fastjet::twopi.

Referenced by boost(), operator+=(), operator-=(), PseudoJet(), reset(), and unboost().

00063                               {
00064   _kt2 = this->px()*this->px() + this->py()*this->py();
00065   if (_kt2 == 0.0) {
00066     _phi = 0.0; } 
00067   else {
00068     _phi = atan2(this->py(),this->px());
00069   }
00070   if (_phi < 0.0) {_phi += twopi;}
00071   if (_phi >= twopi) {_phi -= twopi;} // can happen if phi=-|eps<1e-15|?
00072   if (this->E() == abs(this->pz()) && _kt2 == 0) {
00073     // Point has infinite rapidity -- convert that into a very large
00074     // number, but in such a way that different 0-pt momenta will have
00075     // different rapidities (so as to lift the degeneracy between
00076     // them) [this can be relevant at parton-level]
00077     double MaxRapHere = MaxRap + abs(this->pz());
00078     if (this->pz() >= 0.0) {_rap = MaxRapHere;} else {_rap = -MaxRapHere;}
00079   } else {
00080     // get the rapidity in a way that's modestly insensitive to roundoff
00081     // error when things pz,E are large (actually the best we can do without
00082     // explicit knowledge of mass)
00083     double effective_m2 = max(0.0,m2()); // force non tachyonic mass
00084     double E_plus_pz    = _E + abs(_pz); // the safer of p+, p-
00085     // p+/p- = (p+ p-) / (p-)^2 = (kt^2+m^2)/(p-)^2
00086     _rap = 0.5*log((_kt2 + effective_m2)/(E_plus_pz*E_plus_pz));
00087     if (_pz > 0) {_rap = - _rap;}
00088   }
00089 
00091   //if (this->E() != abs(this->pz())) {
00092   //  _rap = 0.5*log((this->E() + this->pz())/(this->E() - this->pz()));
00093   //    } else {
00094   //  // Overlapping points can give problems. Let's lift the degeneracy
00095   //  // in case of multiple 0-pT points (can be found at parton-level)
00096   //  double MaxRapHere = MaxRap + abs(this->pz());
00097   //  if (this->pz() >= 0.0) {_rap = MaxRapHere;} else {_rap = -MaxRapHere;}
00098   //}
00099 }

void fastjet::PseudoJet::_reset_indices (  )  [inline, private]

set the indices to default values

Definition at line 280 of file PseudoJet.hh.

References set_cluster_hist_index(), and set_user_index().

Referenced by PseudoJet(), and reset().

00280                                       { 
00281   set_cluster_hist_index(-1);
00282   set_user_index(-1);
00283 }


Member Data Documentation

double fastjet::PseudoJet::_px [private]

Definition at line 195 of file PseudoJet.hh.

Referenced by boost(), four_mom(), operator *=(), operator+=(), operator-=(), PseudoJet(), reset(), and unboost().

double fastjet::PseudoJet::_py [private]

Definition at line 195 of file PseudoJet.hh.

Referenced by boost(), four_mom(), operator *=(), operator+=(), operator-=(), PseudoJet(), reset(), and unboost().

double fastjet::PseudoJet::_pz [private]

Definition at line 195 of file PseudoJet.hh.

Referenced by _finish_init(), boost(), four_mom(), operator *=(), operator+=(), operator-=(), PseudoJet(), reset(), and unboost().

double fastjet::PseudoJet::_E [private]

Definition at line 195 of file PseudoJet.hh.

Referenced by _finish_init(), boost(), four_mom(), operator *=(), operator+=(), operator-=(), PseudoJet(), reset(), and unboost().

double fastjet::PseudoJet::_phi [private]

Definition at line 196 of file PseudoJet.hh.

Referenced by _finish_init(), kt_distance(), and plain_distance().

double fastjet::PseudoJet::_rap [private]

Definition at line 196 of file PseudoJet.hh.

Referenced by _finish_init(), kt_distance(), and plain_distance().

double fastjet::PseudoJet::_kt2 [private]

Definition at line 196 of file PseudoJet.hh.

Referenced by _finish_init(), kt_distance(), and operator *=().

int fastjet::PseudoJet::_cluster_hist_index [private]

Definition at line 197 of file PseudoJet.hh.

int fastjet::PseudoJet::_user_index [private]

Definition at line 197 of file PseudoJet.hh.


The documentation for this class was generated from the following files:
Generated on Tue Dec 18 17:05:53 2007 for fastjet by  doxygen 1.5.2