00001 #include "JadeProtoJet.h"
00002 #include <float.h>
00003
00004
00005
00006
00007
00008
00010
00023 template <class Type>
00024 AlephCollection<AlJet> AlephCollection<Type>::JadeJet(float theYCut, int Scheme, float Qelep) const
00025 {
00026
00027 JadeProtoJet::Scheme = Scheme;
00028
00029 AlephCollection<AlJet> result;
00030
00031 vector<JadeProtoJet> theProtoJets;
00032 const_iterator i = begin();
00033 double totalEnergy = 0.;
00034
00035 for (i = begin(); i != end(); i++)
00036 if (!((*i)->isLocked()))
00037 {
00038 totalEnergy += (*i)->A4V().e() ;
00039 JadeProtoJet aProtoJet(*(*i));
00040 theProtoJets.push_back(aProtoJet);
00041 }
00042 if(Qelep != 0.) totalEnergy=Qelep;
00043 int Np = theProtoJets.size();
00044 double totalE2 = totalEnergy*totalEnergy;
00045
00046
00047 if (Np < 2)
00048 {
00049 result.resize(0);
00050 fout <<
00051 "there is only one particle, therefore AJetJn cannot be performed"
00052 << endl;
00053 return result;
00054 }
00055
00056 double minY = 0;
00057 int iYcut = (int) theYCut;
00058
00059 while(((minY<theYCut) && (theYCut >= 0.)) ||
00060 (((int)theProtoJets.size() != -iYcut) && (theYCut < 0.)))
00061 {
00062 minY = DBL_MAX;
00063 vector<JadeProtoJet>::iterator aCandidate = 0;
00064 vector<JadeProtoJet>::iterator aNextNeighbour = 0;
00065 vector<JadeProtoJet>::iterator i1 = theProtoJets.begin();
00066 for(i1 = theProtoJets.begin(); i1!= theProtoJets.end(); i1++)
00067 {
00068 double y = DBL_MAX;
00069 vector<JadeProtoJet>::iterator aNeighbour = i1->findNextNeighbour(theProtoJets);
00070 if(aNeighbour != 0) y = i1->getYTimesE2(aNeighbour)/totalE2;
00071 if(y<minY)
00072 {
00073 minY = y;
00074 aCandidate = i1;
00075 aNextNeighbour = aNeighbour;
00076 }
00077 }
00078 if(((minY<theYCut) && (theYCut >= 0.)) ||
00079 (((int)theProtoJets.size() != -iYcut) && (theYCut < 0.)))
00080 {
00081 *aCandidate += *aNextNeighbour;
00082 theProtoJets.erase(aNextNeighbour);
00083 }
00084 }
00085
00086
00087
00088 vector<JadeProtoJet>::iterator aJetList;
00089 for(aJetList = theProtoJets.begin();
00090 aJetList != theProtoJets.end();
00091 aJetList++)
00092 result.push_back(*aJetList);
00093
00094 return result;
00095 }
00096