Index: trunk/src/SmearUtil.cc
===================================================================
--- trunk/src/SmearUtil.cc	(revision 177)
+++ trunk/src/SmearUtil.cc	(revision 219)
@@ -13,15 +13,26 @@
 
 
-#include "interface/SmearUtil.h"
+#include "SmearUtil.h"
 #include "TRandom.h"
 
 #include <iostream>
+#include <fstream>
 #include <sstream>
-#include <fstream>
 #include <iomanip>
-
-
-
 using namespace std;
+
+
+ParticleUtil::ParticleUtil(const TLorentzVector &genMomentum, int pid) {
+  _pid=pid;
+  _e = genMomentum.E();
+  _px = genMomentum.Px();
+  _py = genMomentum.Py();
+  _pz = genMomentum.Pz();
+  _pt = genMomentum.Pt();
+
+  //_e, _px, _py, _pz, _pt;
+  //float _eta, _etaCalo, _phi, _phiCalo;
+  //int _pid;
+}
 
 //------------------------------------------------------------------------------
@@ -147,4 +158,218 @@
 }
 
+
+RESOLution::RESOLution(const RESOLution & DET) {
+  // Detector characteristics
+  CEN_max_tracker   = DET.CEN_max_tracker; 
+  CEN_max_calo_cen  = DET.CEN_max_calo_cen;
+  CEN_max_calo_fwd  = DET.CEN_max_calo_fwd;
+  CEN_max_mu        = DET.CEN_max_mu;
+  
+  // Energy resolution for electron/photon
+  ELG_Scen          = DET.ELG_Scen;
+  ELG_Ncen          = DET.ELG_Ncen;
+  ELG_Ccen          = DET.ELG_Ccen;
+  ELG_Cfwd          = DET.ELG_Cfwd;
+  ELG_Sfwd          = DET.ELG_Sfwd;
+  ELG_Nfwd          = DET.ELG_Nfwd;
+
+  // Energy resolution for hadrons in ecal/hcal/hf
+  HAD_Shcal         = DET.HAD_Shcal;
+  HAD_Nhcal         = DET.HAD_Nhcal;
+  HAD_Chcal         = DET.HAD_Chcal;
+  HAD_Shf           = DET.HAD_Shf;
+  HAD_Nhf           = DET.HAD_Nhf;
+  HAD_Chf           = DET.HAD_Chf;
+
+  // Muon smearing
+  MU_SmearPt        = DET.MU_SmearPt;
+
+  // Tracking efficiencies
+  TRACK_ptmin       = DET.TRACK_ptmin;
+  TRACK_eff         = DET.TRACK_eff;
+
+  // Calorimetric towers
+  TOWER_number      = DET.TOWER_number;
+  TOWER_eta_edges = new float[TOWER_number+1];
+  for(unsigned int i=0; i<TOWER_number +1; i++) TOWER_eta_edges[i] = DET.TOWER_eta_edges[i];
+
+  TOWER_dphi = new float[TOWER_number];
+  for(unsigned int i=0; i<TOWER_number; i++) TOWER_dphi[i] = DET.TOWER_dphi[i];
+
+  // Thresholds for reconstructed objetcs
+  PTCUT_elec      = DET.PTCUT_elec;
+  PTCUT_muon      = DET.PTCUT_muon;
+  PTCUT_jet       = DET.PTCUT_jet;
+  PTCUT_gamma     = DET.PTCUT_gamma;
+  PTCUT_taujet    = DET.PTCUT_taujet;
+
+  // General jet variable
+  JET_coneradius   = DET.JET_coneradius;
+  JET_jetalgo      = DET.JET_jetalgo;
+  JET_seed         = DET.JET_seed;
+
+  // Tagging definition
+  BTAG_b           = DET.BTAG_b;
+  BTAG_mistag_c    = DET.BTAG_mistag_c;
+  BTAG_mistag_l    = DET.BTAG_mistag_l;
+
+  // FLAGS
+  FLAG_bfield      = DET.FLAG_bfield;
+  FLAG_vfd         = DET.FLAG_vfd;
+  FLAG_trigger     = DET.FLAG_trigger;
+  FLAG_frog        = DET.FLAG_frog;
+
+  // In case BField propagation allowed
+  TRACK_radius      = DET.TRACK_radius;
+  TRACK_length      = DET.TRACK_length;
+  TRACK_bfield_x    = DET.TRACK_bfield_x;
+  TRACK_bfield_y    = DET.TRACK_bfield_y;
+  TRACK_bfield_z    = DET.TRACK_bfield_z;
+
+  // In case Very forward detectors allowed
+  VFD_min_calo_vfd  = DET.VFD_min_calo_vfd;
+  VFD_max_calo_vfd  = DET.VFD_max_calo_vfd;
+  VFD_min_zdc       = DET.VFD_min_zdc;
+  VFD_s_zdc         = DET.VFD_s_zdc;
+
+  RP_220_s          = DET.RP_220_s;
+  RP_220_x          = DET.RP_220_x;
+  RP_420_s          = DET.RP_420_s;
+  RP_420_x          = DET.RP_420_x;
+
+  // In case FROG event display allowed
+  NEvents_Frog      = DET.NEvents_Frog;
+
+  JET_overlap      = DET.JET_overlap;
+  // MidPoint algorithm definition
+  JET_M_coneareafraction = DET.JET_M_coneareafraction;
+  JET_M_maxpairsize      = DET.JET_M_maxpairsize;
+  JET_M_maxiterations    = DET.JET_M_maxiterations;
+  // Define Cone algorithm.
+  JET_C_adjacencycut     = DET.JET_C_adjacencycut;
+  JET_C_maxiterations    = DET.JET_C_maxiterations;
+  JET_C_iratch           = DET.JET_C_iratch;
+  //Define SISCone algorithm.
+  JET_S_npass            = DET.JET_S_npass;
+  JET_S_protojet_ptmin   = DET.JET_S_protojet_ptmin;
+
+  //For Tau-jet definition
+  TAU_energy_scone = DET.TAU_energy_scone;
+  TAU_track_scone  = DET.TAU_track_scone;
+  TAU_track_pt     = DET.TAU_track_pt;
+  TAU_energy_frac  = DET.TAU_energy_frac;
+
+  PT_QUARKS_MIN    = DET.PT_QUARKS_MIN;
+}
+
+RESOLution& RESOLution::operator=(const RESOLution& DET) {
+  if(this==&DET) return *this;
+  // Detector characteristics
+  CEN_max_tracker   = DET.CEN_max_tracker; 
+  CEN_max_calo_cen  = DET.CEN_max_calo_cen;
+  CEN_max_calo_fwd  = DET.CEN_max_calo_fwd;
+  CEN_max_mu        = DET.CEN_max_mu;
+  
+  // Energy resolution for electron/photon
+  ELG_Scen          = DET.ELG_Scen;
+  ELG_Ncen          = DET.ELG_Ncen;
+  ELG_Ccen          = DET.ELG_Ccen;
+  ELG_Cfwd          = DET.ELG_Cfwd;
+  ELG_Sfwd          = DET.ELG_Sfwd;
+  ELG_Nfwd          = DET.ELG_Nfwd;
+
+  // Energy resolution for hadrons in ecal/hcal/hf
+  HAD_Shcal         = DET.HAD_Shcal;
+  HAD_Nhcal         = DET.HAD_Nhcal;
+  HAD_Chcal         = DET.HAD_Chcal;
+  HAD_Shf           = DET.HAD_Shf;
+  HAD_Nhf           = DET.HAD_Nhf;
+  HAD_Chf           = DET.HAD_Chf;
+
+  // Muon smearing
+  MU_SmearPt        = DET.MU_SmearPt;
+
+  // Tracking efficiencies
+  TRACK_ptmin       = DET.TRACK_ptmin;
+  TRACK_eff         = DET.TRACK_eff;
+
+  // Calorimetric towers
+  TOWER_number      = DET.TOWER_number;
+  TOWER_eta_edges = new float[TOWER_number+1];
+  for(unsigned int i=0; i<TOWER_number +1; i++) TOWER_eta_edges[i] = DET.TOWER_eta_edges[i];
+
+  TOWER_dphi = new float[TOWER_number];
+  for(unsigned int i=0; i<TOWER_number; i++) TOWER_dphi[i] = DET.TOWER_dphi[i];
+
+  // Thresholds for reconstructed objetcs
+  PTCUT_elec      = DET.PTCUT_elec;
+  PTCUT_muon      = DET.PTCUT_muon;
+  PTCUT_jet       = DET.PTCUT_jet;
+  PTCUT_gamma     = DET.PTCUT_gamma;
+  PTCUT_taujet    = DET.PTCUT_taujet;
+
+  // General jet variable
+  JET_coneradius   = DET.JET_coneradius;
+  JET_jetalgo      = DET.JET_jetalgo;
+  JET_seed         = DET.JET_seed;
+
+  // Tagging definition
+  BTAG_b           = DET.BTAG_b;
+  BTAG_mistag_c    = DET.BTAG_mistag_c;
+  BTAG_mistag_l    = DET.BTAG_mistag_l;
+
+  // FLAGS
+  FLAG_bfield      = DET.FLAG_bfield;
+  FLAG_vfd         = DET.FLAG_vfd;
+  FLAG_trigger     = DET.FLAG_trigger;
+  FLAG_frog        = DET.FLAG_frog;
+
+  // In case BField propagation allowed
+  TRACK_radius      = DET.TRACK_radius;
+  TRACK_length      = DET.TRACK_length;
+  TRACK_bfield_x    = DET.TRACK_bfield_x;
+  TRACK_bfield_y    = DET.TRACK_bfield_y;
+  TRACK_bfield_z    = DET.TRACK_bfield_z;
+
+  // In case Very forward detectors allowed
+  VFD_min_calo_vfd  = DET.VFD_min_calo_vfd;
+  VFD_max_calo_vfd  = DET.VFD_max_calo_vfd;
+  VFD_min_zdc       = DET.VFD_min_zdc;
+  VFD_s_zdc         = DET.VFD_s_zdc;
+
+  RP_220_s          = DET.RP_220_s;
+  RP_220_x          = DET.RP_220_x;
+  RP_420_s          = DET.RP_420_s;
+  RP_420_x          = DET.RP_420_x;
+
+  // In case FROG event display allowed
+  NEvents_Frog      = DET.NEvents_Frog;
+
+  JET_overlap      = DET.JET_overlap;
+  // MidPoint algorithm definition
+  JET_M_coneareafraction = DET.JET_M_coneareafraction;
+  JET_M_maxpairsize      = DET.JET_M_maxpairsize;
+  JET_M_maxiterations    = DET.JET_M_maxiterations;
+  // Define Cone algorithm.
+  JET_C_adjacencycut     = DET.JET_C_adjacencycut;
+  JET_C_maxiterations    = DET.JET_C_maxiterations;
+  JET_C_iratch           = DET.JET_C_iratch;
+  //Define SISCone algorithm.
+  JET_S_npass            = DET.JET_S_npass;
+  JET_S_protojet_ptmin   = DET.JET_S_protojet_ptmin;
+
+  //For Tau-jet definition
+  TAU_energy_scone = DET.TAU_energy_scone;
+  TAU_track_scone  = DET.TAU_track_scone;
+  TAU_track_pt     = DET.TAU_track_pt;
+  TAU_energy_frac  = DET.TAU_energy_frac;
+
+  PT_QUARKS_MIN    = DET.PT_QUARKS_MIN;
+  return *this;
+}
+
+
+
+
 //------------------------------------------------------------------------------
 void RESOLution::ReadDataCard(const string datacard) {
@@ -155,5 +380,5 @@
   ifstream fichier_a_lire(datacard.c_str());
   if(!fichier_a_lire.good()) {
-    cout << datacard << "Datadard " << datacard << " not found, use default values" << endl;
+    cout <<"**        Datadard  not found, use default values                  **" << endl;
     return;
   }
@@ -252,5 +477,5 @@
 }
 
-void RESOLution::Logfile(string LogName) { 
+void RESOLution::Logfile(const string& LogName) { 
   //void RESOLution::Logfile(string outputfilename) { 
   
@@ -617,5 +842,5 @@
   energyS = ((energyS1>0)?energyS1:0) + ((energyS2>0)?energyS2:0);
   } 
-  if(abs(hadron.Eta()) > CEN_max_calo_cen && fabs(hadron.Eta()) < CEN_max_calo_fwd){ 
+  if(fabs(hadron.Eta()) > CEN_max_calo_cen && fabs(hadron.Eta()) < CEN_max_calo_fwd){ 
     energyS = gRandom->Gaus(energy, sqrt(
 			    pow(HAD_Nhf,2) +
@@ -720,13 +945,13 @@
 //***********************Isolation criteria***********************
 //****************************************************************
-bool RESOLution::Isolation(Float_t phi,Float_t eta,const vector<TLorentzVector> &tracks,float PT_TRACK2)
+bool RESOLution::Isolation(const float phi, const float eta,const vector<TLorentzVector> &tracks, const float pt_second_track)
 {
    bool isolated = false;
-   Float_t deltar=5000.; // Initial value; should be high; no further repercussion
+   float deltar=5000.; // Initial value; should be high; no further repercussion
    // loop on all final charged particles, with p_t >2, close enough from the electron
    for(unsigned int i=0; i < tracks.size(); i++) 
       {
-       if(tracks[i].Pt() < PT_TRACK2)continue;
-       Float_t genDeltaR = DeltaR(phi,eta,tracks[i].Phi(),tracks[i].Eta());    // slower to evaluate
+       if(tracks[i].Pt() < pt_second_track)continue;
+       float genDeltaR = DeltaR(phi,eta,tracks[i].Phi(),tracks[i].Eta()); 
          if(
              (genDeltaR > deltar) ||
@@ -735,5 +960,5 @@
            deltar=genDeltaR;
       }
-   if(deltar > 0.5)isolated = true; // returns the closest distance
+   if(deltar > 0.5) isolated = true;
    return isolated;
 }
@@ -769,5 +994,7 @@
 float DeltaPhi(const float phi1, const float phi2) {
   float deltaphi=phi1-phi2; // in here, -PI < phi < PI
-  if(fabs(deltaphi) > PI) deltaphi=2.*PI-fabs(deltaphi);// put deltaphi between 0 and PI
+  if(fabs(deltaphi) > PI) {
+	deltaphi=2.*PI -fabs(deltaphi);// put deltaphi between 0 and PI
+	}
   else deltaphi=fabs(deltaphi);
   
