Index: trunk/Delphes.cpp
===================================================================
--- trunk/Delphes.cpp	(revision 18)
+++ trunk/Delphes.cpp	(revision 19)
@@ -197,7 +197,56 @@
   beamline2->add(rp420_2);
 
-  
+  // we will have four jet definitions, and the first three will be
+  // plugins
+  fastjet::JetDefinition jet_def;
+  fastjet::JetDefinition::Plugin * plugins;
+
+  switch(DET->JETALGO) {
+  default:
+  case 1: {
+
+   // set up a CDF midpoint jet definition
+   #ifdef ENABLE_PLUGIN_CDFCONES
+   plugins = new fastjet::CDFJetCluPlugin(DET->C_SEEDTHRESHOLD,DET->CONERADIUS,DET->C_ADJACENCYCUT,DET->C_MAXITERATIONS,DET->C_IRATCH,DET->C_OVERLAPTHRESHOLD);
+   jet_def = fastjet::JetDefinition(plugins);
+   #else
+   plugins = NULL;
+   #endif
+   }
+   break;
+
+   case 2: {
+
+   // set up a CDF midpoint jet definition
+   #ifdef ENABLE_PLUGIN_CDFCONES
+   plugins = new fastjet::CDFMidPointPlugin (DET->M_SEEDTHRESHOLD,DET->CONERADIUS,DET->M_CONEAREAFRACTION,DET->M_MAXPAIRSIZE,DET->M_MAXPAIRSIZE,DET->C_OVERLAPTHRESHOLD);
+   jet_def = fastjet::JetDefinition(plugins);
+   #else
+   plugins = NULL;
+   #endif
+   }
+   break;
+   case 3: {
+   // set up a siscone jet definition
+   #ifdef ENABLE_PLUGIN_SISCONE
+   int npass = 0;               // do infinite number of passes
+   double protojet_ptmin = 0.0; // use all protojets
+   plugins = new fastjet::SISConePlugin (DET->CONERADIUS,DET->C_OVERLAPTHRESHOLD,npass, protojet_ptmin);
+   jet_def = fastjet::JetDefinition(plugins);
+   #else
+   plugins = NULL;
+   #endif
+   }
+   break;
+
+  case 4: {
+  jet_def = fastjet::JetDefinition(fastjet::kt_algorithm, DET->CONERADIUS);
+  //jet_defs[4] = fastjet::JetDefinition(fastjet::cambridge_algorithm,jet_radius);
+  //jet_defs[5] = fastjet::JetDefinition(fastjet::antikt_algorithm,jet_radius);
+  }
+  break;
+  }
  
-  // Loop over all events
+   // Loop over all events
   Long64_t entry, allEntries = treeReader->GetEntries();
   cout << "** Chain contains " << allEntries << " events" << endl;
@@ -409,56 +458,4 @@
       //*****************************
 
-      // we will have four jet definitions, and the first three will be
-      // plugins
-      fastjet::JetDefinition jet_def;
-      fastjet::JetDefinition::Plugin * plugins;
-      
-      switch(DET->JETALGO) {
-      default:
-      case 1: {
-	
-	// set up a CDF midpoint jet definition
-        #ifdef ENABLE_PLUGIN_CDFCONES
-	plugins = new fastjet::CDFJetCluPlugin(DET->C_SEEDTHRESHOLD,DET->CONERADIUS,DET->C_ADJACENCYCUT,DET->C_MAXITERATIONS,DET->C_IRATCH,DET->C_OVERLAPTHRESHOLD);
-	jet_def = fastjet::JetDefinition(plugins);
-        #else
-	plugins = NULL;
-        #endif
-	
-      } 
-	break;
-	
-      case 2: {
-	
-        // set up a CDF midpoint jet definition
-        #ifdef ENABLE_PLUGIN_CDFCONES
-        plugins = new fastjet::CDFMidPointPlugin (DET->M_SEEDTHRESHOLD,DET->CONERADIUS,DET->M_CONEAREAFRACTION,DET->M_MAXPAIRSIZE,DET->M_MAXPAIRSIZE,DET->C_OVERLAPTHRESHOLD);
-	jet_def = fastjet::JetDefinition(plugins);
-        #else
-	plugins = NULL;
-        #endif
-      }
-	break;
-	
-      case 3: {
-	// set up a siscone jet definition
-        #ifdef ENABLE_PLUGIN_SISCONE
-	int npass = 0;               // do infinite number of passes
-	double protojet_ptmin = 0.0; // use all protojets
-	plugins = new fastjet::SISConePlugin (DET->CONERADIUS,DET->C_OVERLAPTHRESHOLD,npass, protojet_ptmin);
-	jet_def = fastjet::JetDefinition(plugins);
-        #else
-        plugins = NULL;
-        #endif
-      }
-        break; 
-	
-      case 4: {
-        jet_def = fastjet::JetDefinition(fastjet::kt_algorithm, DET->CONERADIUS);
-        //jet_defs[4] = fastjet::JetDefinition(fastjet::cambridge_algorithm,jet_radius);
-        //jet_defs[5] = fastjet::JetDefinition(fastjet::antikt_algorithm,jet_radius);
-      }
-        break;
-      }
       // run the jet clustering with the above jet definition
       if(input_particles.size()!=0)
@@ -478,4 +475,5 @@
 	TLorentzVector JET;
 	JET.SetPxPyPzE(sorted_jets[i].px(),sorted_jets[i].py(),sorted_jets[i].pz(),sorted_jets[i].E());
+        cout<<"Jet.Pt() "<<JET.Pt()<<endl;
 	elementJet->Set(JET);
 	// b-jets
Index: trunk/Makefile
===================================================================
--- trunk/Makefile	(revision 18)
+++ trunk/Makefile	(revision 19)
@@ -56,5 +56,9 @@
 	interface/LHEFConverter.h \
 	interface/STDHEPConverter.h \
-	interface/SmearUtil.h
+	interface/SmearUtil.h \
+	Utilities/Fastjet/include/fastjet/PseudoJet.hh \
+	Utilities/Fastjet/include/fastjet/ClusterSequence.hh \
+	Utilities/Fastjet/include/fastjet/config.h \
+	interface/TreeClasses.h
 EXECUTABLE =  \
 	Delphes$(ExeSuf) \
@@ -68,9 +72,14 @@
 	Utilities/ExRootAnalysis/src/BlockClassesLinkDef.h \
 	Utilities/ExRootAnalysis/interface/BlockClasses.h
+tmp/src/TreeClassesDict.$(SrcSuf): \
+	src/TreeClassesLinkDef.h \
+	interface/TreeClasses.h
 DICT =  \
-	tmp/Utilities/ExRootAnalysis/src/BlockClassesDict.$(SrcSuf)
+	tmp/Utilities/ExRootAnalysis/src/BlockClassesDict.$(SrcSuf) \
+	tmp/src/TreeClassesDict.$(SrcSuf)
 
 DICT_OBJ =  \
-	tmp/Utilities/ExRootAnalysis/src/BlockClassesDict.$(ObjSuf)
+	tmp/Utilities/ExRootAnalysis/src/BlockClassesDict.$(ObjSuf) \
+	tmp/src/TreeClassesDict.$(ObjSuf)
 
 tmp/src/HEPEVTConverter.$(ObjSuf): \
@@ -101,4 +110,7 @@
 	src/SmearUtil.$(SrcSuf) \
 	interface/SmearUtil.h
+tmp/src/TreeClasses.$(ObjSuf): \
+	src/TreeClasses.$(SrcSuf) \
+	interface/TreeClasses.h
 tmp/Utilities/ExRootAnalysis/src/BlockClasses.$(ObjSuf): \
 	Utilities/ExRootAnalysis/src/BlockClasses.$(SrcSuf) \
@@ -281,4 +293,5 @@
 	tmp/src/STDHEPConverter.$(ObjSuf) \
 	tmp/src/SmearUtil.$(ObjSuf) \
+	tmp/src/TreeClasses.$(ObjSuf) \
 	tmp/Utilities/ExRootAnalysis/src/BlockClasses.$(ObjSuf) \
 	tmp/Utilities/ExRootAnalysis/src/ExRootProgressBar.$(ObjSuf) \
@@ -468,4 +481,9 @@
 	@touch $@
 
+interface/TreeClasses.h: \
+	Utilities/ExRootAnalysis/interface/BlockCompare.h \
+	Utilities/ExRootAnalysis/interface/BlockClasses.h
+	@touch $@
+
 interface/LHEFConverter.h: \
 	Utilities/ExRootAnalysis/interface/ExRootTreeBranch.h \
@@ -497,4 +515,8 @@
 	@touch $@
 
+Utilities/ExRootAnalysis/interface/BlockClasses.h: \
+	Utilities/ExRootAnalysis/interface/BlockCompare.h
+	@touch $@
+
 interface/STDHEPConverter.h: \
 	Utilities/ExRootAnalysis/interface/BlockClasses.h \
@@ -502,8 +524,4 @@
 	Utilities/ExRootAnalysis/interface/LHEF.h \
 	interface/DataConverter.h
-	@touch $@
-
-Utilities/ExRootAnalysis/interface/BlockClasses.h: \
-	Utilities/ExRootAnalysis/interface/BlockCompare.h
 	@touch $@
 
Index: trunk/Resolutions.cpp
===================================================================
--- trunk/Resolutions.cpp	(revision 18)
+++ trunk/Resolutions.cpp	(revision 19)
@@ -19,8 +19,7 @@
 #include "Utilities/ExRootAnalysis/interface/ExRootTreeBranch.h"
 
-#include "Utilities/CDFCones/interface/JetCluAlgorithm.h"
-#include "Utilities/CDFCones/interface/MidPointAlgorithm.h"
-#include "Utilities/CDFCones/interface/PhysicsTower.h"
-#include "Utilities/CDFCones/interface/Cluster.h"
+#include "H_BeamParticle.h"
+#include "H_BeamLine.h"
+#include "H_RomanPot.h"
 
 #include "interface/DataConverter.h"
@@ -28,7 +27,25 @@
 #include "interface/LHEFConverter.h"
 #include "interface/STDHEPConverter.h"
+
 #include "interface/SmearUtil.h"
+#include "Utilities/Fastjet/include/fastjet/PseudoJet.hh"
+#include "Utilities/Fastjet/include/fastjet/ClusterSequence.hh"
+
+// get info on how fastjet was configured
+#include "Utilities/Fastjet/include/fastjet/config.h"
+
+// make sure we have what is needed
+#ifdef ENABLE_PLUGIN_SISCONE
+#  include "Utilities/Fastjet/plugins/SISCone/SISConePlugin.hh"
+#endif
+#ifdef ENABLE_PLUGIN_CDFCONES
+#  include "Utilities/Fastjet/plugins/CDFCones/CDFMidPointPlugin.hh"
+#  include "Utilities/Fastjet/plugins/CDFCones/CDFJetCluPlugin.hh"
+#endif
+
+#include<vector>
+#include<iostream>
+
 #include "interface/TreeClasses.h"
-
 using namespace std;
 
@@ -37,24 +54,20 @@
 
 
-void PairingJet(TLorentzVector &JETSm, TLorentzVector JET, vector<Cluster> jetsS)
+void PairingJet(TLorentzVector &JETSm, TLorentzVector JET,  vector<fastjet::PseudoJet> sorted_jetsS)
 {
   JETSm.SetPxPyPzE(0,0,0,0);
-  vector<Cluster>::iterator itJetS;
   float deltaRtest=5000;
-  LorentzVector jetMomentumS;
-  for(itJetS = jetsS.begin(); itJetS != jetsS.end(); ++itJetS) 
-     {
-       jetMomentumS = itJetS->fourVector;
-       TLorentzVector Att;
-       Att.SetPxPyPzE(jetMomentumS.px,jetMomentumS.py,jetMomentumS.pz,jetMomentumS.E);
-       if(DeltaR(JET.Phi(),JET.Eta(),Att.Phi(),Att.Eta()) < deltaRtest)
-            {
-              deltaRtest = DeltaR(JET.Phi(),JET.Eta(),Att.Phi(),Att.Eta());
-              if(deltaRtest < 0.25)
-                {
-                  JETSm = Att;
-                }
-            }
-     }
+  for (unsigned int i = 0; i < sorted_jetsS.size(); i++) {
+      TLorentzVector Att;
+      Att.SetPxPyPzE(sorted_jetsS[i].px(),sorted_jetsS[i].py(),sorted_jetsS[i].pz(),sorted_jetsS[i].E());
+      if(DeltaR(JET.Phi(),JET.Eta(),Att.Phi(),Att.Eta()) < deltaRtest)
+        {
+          deltaRtest = DeltaR(JET.Phi(),JET.Eta(),Att.Phi(),Att.Eta());
+          if(deltaRtest < 0.25)
+             {
+               JETSm = Att;
+             }
+        }
+    }
 }
 
@@ -148,12 +161,36 @@
   DET->ReadDataCard(DetDatacard);
   
-  
   TLorentzVector genMomentum(0,0,0,0);
   LorentzVector jetMomentum;
-  vector<TLorentzVector> TrackCentral;  
-  vector<PhysicsTower> towersS;
-  vector<Cluster> jetsS; 
-  vector<PhysicsTower> towers;
-  vector<Cluster> jets;
+  vector<TLorentzVector> TrackCentral;
+
+  vector<fastjet::PseudoJet> input_particles;//for FastJet algorithm
+  vector<fastjet::PseudoJet> inclusive_jets;
+  vector<fastjet::PseudoJet> sorted_jets;
+
+  vector<fastjet::PseudoJet> input_particlesS;//for FastJet algorithm
+  vector<fastjet::PseudoJet> inclusive_jetsS;
+  vector<fastjet::PseudoJet> sorted_jetsS;
+
+  fastjet::JetDefinition jet_def;
+  fastjet::JetDefinition jet_defS;
+  fastjet::JetDefinition::Plugin * plugins;
+  fastjet::JetDefinition::Plugin * pluginsS;
+
+  // set up a CDF midpoint jet definition
+  #ifdef ENABLE_PLUGIN_CDFCONES
+  plugins = new fastjet::CDFJetCluPlugin(DET->C_SEEDTHRESHOLD,DET->CONERADIUS,DET->C_ADJACENCYCUT,DET->C_MAXITERATIONS,DET->C_IRATCH,DET->C_OVERLAPTHRESHOLD);
+  jet_def = fastjet::JetDefinition(plugins);
+  #else
+  plugins = NULL;
+  #endif
+
+  // set up a CDF midpoint jet definition
+  #ifdef ENABLE_PLUGIN_CDFCONES
+  pluginsS = new fastjet::CDFJetCluPlugin(2,DET->CONERADIUS,DET->C_ADJACENCYCUT,DET->C_MAXITERATIONS,DET->C_IRATCH,DET->C_OVERLAPTHRESHOLD);
+  jet_defS = fastjet::JetDefinition(pluginsS);
+  #else
+  pluginsS = NULL;
+  #endif
 
  
@@ -172,8 +209,11 @@
       itGen.Reset();
       TrackCentral.clear();
-      towers.clear();
-      towersS.clear();
       TSimpleArray<TRootGenParticle> NFCentralQ;
-
+      input_particles.clear();
+      inclusive_jets.clear();
+      sorted_jets.clear();
+      input_particlesS.clear();
+      inclusive_jetsS.clear();
+      sorted_jetsS.clear();
 
       // Loop over all particles in event
@@ -185,5 +225,8 @@
 	  float eta = fabs(particle->Eta);
 
-if(particle->Status == 1)towers.push_back(PhysicsTower(LorentzVector(genMomentum.Px(),genMomentum.Py(),genMomentum.Pz(), genMomentum.E())));
+          if(particle->Status == 1)
+             {
+               input_particles.push_back(fastjet::PseudoJet(genMomentum.Px(),genMomentum.Py(),genMomentum.Pz(), genMomentum.E()));
+             }
 
 	  // keeps only final particles, visible by the central detector, including the fiducial volume
@@ -224,12 +267,9 @@
 
 		if(pid != pMU && genMomentum.Pt() > DET->PT_TRACKS_MIN)
-                towersS.push_back
-                 (
-                   PhysicsTower
-                    (
-                      LorentzVector(genMomentum.Px(),genMomentum.Py(),genMomentum.Pz(), genMomentum.E())
-                    )
-                 );
-		// all final charged particles
+                  {
+                     input_particlesS.push_back(fastjet::PseudoJet(genMomentum.Px(),genMomentum.Py(),genMomentum.Pz(), genMomentum.E()));
+                  }
+		
+                // all final charged particles
 		if(
 		     ((rand()%100) < DET->TRACKING_EFF)  &&
@@ -251,29 +291,36 @@
 
       //*****************************
-      JetCluAlgorithm jetAlgoS(1,DET->CONERADIUS,DET->C_ADJACENCYCUT,DET->C_MAXITERATIONS,DET->C_IRATCH,DET->C_OVERLAPTHRESHOLD);
-      jetAlgoS.run(towersS, jetsS);
-
-      JetCluAlgorithm jetAlgo(0,DET->CONERADIUS,DET->C_ADJACENCYCUT,DET->C_MAXITERATIONS,DET->C_IRATCH,DET->C_OVERLAPTHRESHOLD);
-      jetAlgo.run(towers, jets);
-
-
-      vector<Cluster>::iterator itJet;
-      TLorentzVector JETSm(0,0,0,0);
-      for(itJet = jets.begin(); itJet != jets.end(); ++itJet) {
-          jetMomentum = itJet->fourVector;
-          TLorentzVector JET(0,0,0,0);
-          JET.SetPxPyPzE(jetMomentum.px,jetMomentum.py,jetMomentum.pz,jetMomentum.E);
-          PairingJet(JETSm,JET,jetsS);
-          if(JETSm.Pt()>3)
-            {
-              elementJet= (RESOLJET*) branchjet->NewEntry();
-              elementJet->NonSmearePT = JET.Pt();
-              elementJet->SmearePT = (JETSm.Pt()/JET.Pt());
-              /*cout<<"valeur obtenue "<<JETSm.Pt()/JET.Pt()<<endl;
-              cout<<" pt non smeare "<<JET.Pt()<<" phi "<<JET.Phi()<<" eta "<<JET.Eta()<<endl;
-              cout<<"pt smeare "<<JETSm.Pt()<<" phi "<<JETSm.Phi()<<" eta "<<JETSm.Eta()<<endl;
-              cout<<"************"<<endl;
+
+     double ptmin=1;
+     if(input_particles.size()!=0)
+        {
+          fastjet::ClusterSequence clust_seq(input_particles, jet_def);
+          inclusive_jets = clust_seq.inclusive_jets(ptmin);
+          sorted_jets = sorted_by_pt(inclusive_jets);
+        }
+
+     if(input_particlesS.size()!=0)
+        {
+          fastjet::ClusterSequence clust_seqS(input_particlesS, jet_defS);
+          inclusive_jetsS = clust_seqS.inclusive_jets(ptmin);
+          sorted_jetsS = sorted_by_pt(inclusive_jetsS);
+        }
+
+    TLorentzVector JETSm(0,0,0,0);
+    for (unsigned int i = 0; i < sorted_jets.size(); i++) {
+        TLorentzVector JET(0,0,0,0);
+        JET.SetPxPyPzE(sorted_jets[i].px(),sorted_jets[i].py(),sorted_jets[i].pz(),sorted_jets[i].E());
+        PairingJet(JETSm,JET,sorted_jetsS);
+        if(JETSm.Pt()>1)
+          {
+            elementJet= (RESOLJET*) branchjet->NewEntry();
+            elementJet->NonSmearePT = JET.Pt();
+            elementJet->SmearePT = (JETSm.Pt()/JET.Pt());
+            /*cout<<"valeur obtenue "<<JETSm.Pt()/JET.Pt()<<endl;
+            cout<<" pt non smeare "<<JET.Pt()<<" phi "<<JET.Phi()<<" eta "<<JET.Eta()<<endl;
+            cout<<"pt smeare "<<JETSm.Pt()<<" phi "<<JETSm.Phi()<<" eta "<<JETSm.Eta()<<endl;
+            cout<<"************"<<endl;
 */
-            }
+          }
 
 	} // for itJet : loop on all jets
Index: trunk/genMakefile.tcl
===================================================================
--- trunk/genMakefile.tcl	(revision 18)
+++ trunk/genMakefile.tcl	(revision 19)
@@ -232,5 +232,5 @@
 executableDeps {*.cpp} 
 
-dictDeps {DICT} {Utilities/ExRootAnalysis/src/*LinkDef.h}
+dictDeps {DICT} {Utilities/ExRootAnalysis/src/*LinkDef.h} {src/*LinkDef.h}
 
 sourceDeps {SOURCE} {src/*.cc} {Utilities/ExRootAnalysis/src/*.cc} {Utilities/Hector/src/*.cc} {Utilities/CDFCones/src/*cc} {Utilities/Fastjet/src/*.cc} {Utilities/Fastjet/plugins/CDFCones/*.cc} {Utilities/Fastjet/plugins/CDFCones/src/*.cc} {Utilities/Fastjet/plugins/SISCone/*.cc} {Utilities/Fastjet/plugins/SISCone/src/*.cc}
