---

DurhamProtoJet.cpp

Go to the documentation of this file.

//
// Implementation of
// CLASS : DurhamProtoJet
// Object used to compute jets with Durham
//
// Written by J.P. Wellisch, 25. Jan 1999
// adapted to Alpha++ by C.Delaere
//

#include "DurhamProtoJet.h"
#include "QvecBase.h"
#include <float.h>

// Constructor (from QvecBase)

DurhamProtoJet::DurhamProtoJet(const QvecBase & alobj)
{
  if (alobj.A4V().e()==0)
    {
      cout << "Error: bad Eflw: E=0" << endl;
      exit(1);
    }
  theFourMomentum = alobj.A4V();
  theConstituants.push_back(&alobj);
}

// Operator += (depends on the scheme)
DurhamProtoJet & 
DurhamProtoJet::operator+=(const DurhamProtoJet & aJet)
{
  if (Scheme)
    {
      // recombination according to E-Scheme.
      theFourMomentum += aJet.get4Momentum();
    }
  else
    {
      // recombination according to E0-Scheme.
      double e = theFourMomentum.e() + aJet.get4Momentum().e();
      double p2i = (theFourMomentum + aJet.get4Momentum()).vect().mag();
      if (p2i <= pow(10,-14))
        {
          cout << " Error: Division by zero! " << endl;
          exit(1);
        }
      theFourMomentum.setVect(e/p2i*(theFourMomentum.vect()+aJet.get4Momentum().vect()));
      theFourMomentum.setE(e);
    }

  vector<const QvecBase *> aConstituantList;
  aConstituantList = aJet.getConstituants();
  vector<const QvecBase *>::iterator constitIter;
  for(constitIter=aConstituantList.begin(); 
      constitIter!=aConstituantList.end(); 
      constitIter++)
  {
    theConstituants.push_back(*constitIter);
  }
  
  return *this;
}  

// the Metric
double 
DurhamProtoJet::getYTimesE2(const DurhamProtoJet * aJet)
{
  double result;
  double e1 = theFourMomentum.e();
  double e2 = aJet->get4Momentum().e();
  double costh = theFourMomentum.vect() * aJet->get4Momentum().vect();
  costh /= theFourMomentum.vect().mag();
  costh /= aJet->get4Momentum().vect().mag();
  result = e1;
  if(e1>e2) result = e2;
  result *= result;
  result *= 2.*(1-costh);
  return result;
}

// find the Next Neighbour
vector<DurhamProtoJet>::iterator 
DurhamProtoJet::findNextNeighbour(const vector<DurhamProtoJet> & theProtoJets)
{
  vector<DurhamProtoJet>::const_iterator i = theProtoJets.begin();
  double minDist = DBL_MAX;
  vector<DurhamProtoJet>::const_iterator theNextNeighbour = 0;
  for (i = theProtoJets.begin(); i != theProtoJets.end(); i++)
  {
    if(i!=this)
    {
      double distance = getYTimesE2(i);
      if(distance<minDist)
      {
        minDist = distance;
        theNextNeighbour = i;
      }
    }
  }
  return (vector<DurhamProtoJet>::iterator)theNextNeighbour;
}

// the scheme (set by the algo, used by operator+=)
int DurhamProtoJet::Scheme=1;

---