Index: trunk/src/SmearUtil.cc
===================================================================
--- trunk/src/SmearUtil.cc	(revision 319)
+++ trunk/src/SmearUtil.cc	(revision 321)
@@ -102,7 +102,10 @@
   PTCUT_taujet    = 10.0;
 
+  // Isolation
   ISOL_PT         = 2.0;      //minimal pt of tracks for isolation criteria
   ISOL_Cone       = 0.5;      //Cone  for isolation criteria
-
+  ISOL_Calo_ET    = 1E99;     //minimal tower energy for isolation criteria. Default off = 1E99
+  ISOL_Calo_Cone  = 0.5;      //Cone for calorimetric isolation
+  ISOL_Calo_Grid  = 3;        //Grid size (N x N) for calorimetric isolation
 
   // General jet variable
@@ -228,6 +231,10 @@
   PTCUT_taujet    = DET.PTCUT_taujet;
 
-  ISOL_PT         = DET.ISOL_PT;      //minimal pt of tracks for isolation criteria
-  ISOL_Cone       = DET.ISOL_Cone;      //Cone  for isolation criteria
+  // Isolation
+  ISOL_PT         = DET.ISOL_PT;     // tracking isolation
+  ISOL_Cone       = DET.ISOL_Cone;  
+  ISOL_Calo_ET    = DET.ISOL_Calo_ET;  // calorimeter isolation, defaut off
+  ISOL_Calo_Cone  = DET.ISOL_Calo_Cone;
+  ISOL_Calo_Grid  = DET.ISOL_Calo_Grid; 
 
 
@@ -347,7 +354,10 @@
   PTCUT_taujet    = DET.PTCUT_taujet;
 
-  ISOL_PT         = DET.ISOL_PT;      //minimal pt of tracks for isolation criteria
-  ISOL_Cone       = DET.ISOL_Cone;      //Cone  for isolation criteria
-
+  // Isolation
+  ISOL_PT         = DET.ISOL_PT;       // tracking isolation
+  ISOL_Cone       = DET.ISOL_Cone;    
+  ISOL_Calo_ET    = DET.ISOL_Calo_ET;  // calorimeter isolation, defaut off
+  ISOL_Calo_Cone  = DET.ISOL_Calo_Cone;
+  ISOL_Calo_Grid  = DET.ISOL_Calo_Grid;
 
   // General jet variable
@@ -497,7 +507,9 @@
     else if(strstr(temp_string.c_str(),"PTCUT_taujet"))     {curstring >> varname >> value; PTCUT_taujet      = value;}
 
-    else if(strstr(temp_string.c_str(),"ISOL_PT"))          {curstring >> varname >> value; ISOL_PT        = value;}
-    else if(strstr(temp_string.c_str(),"ISOL_Cone"))        {curstring >> varname >> value; ISOL_Cone        = value;}
-
+    else if(strstr(temp_string.c_str(),"ISOL_PT"))          {curstring >> varname >> value; ISOL_PT           = value;}
+    else if(strstr(temp_string.c_str(),"ISOL_Cone"))        {curstring >> varname >> value; ISOL_Cone         = value;}
+    else if(strstr(temp_string.c_str(),"ISOL_Calo_ET"))     {curstring >> varname >> value; ISOL_Calo_ET      = value;}
+    else if(strstr(temp_string.c_str(),"ISOL_Calo_Cone"))   {curstring >> varname >> value; ISOL_Calo_Cone    = value;}
+    else if(strstr(temp_string.c_str(),"ISOL_Calo_Grid"))   {curstring >> varname >> ivalue; ISOL_Calo_Grid    = ivalue;}
 
     else if(strstr(temp_string.c_str(),"JET_coneradius"))   {curstring >> varname >> value; JET_coneradius    = value;}
@@ -787,4 +799,16 @@
   f_out << left << setw(40) <<"* Cone for isolation criteria: "<<""
         << left << setw(20) <<ISOL_Cone             <<""<< right << setw(10)<<"*"<<"\n";
+
+  if(ISOL_Calo_ET > 1E98) f_out<<"# No Calorimetric isolation applied                                  *"<<"\n";
+  else {
+   f_out << left << setw(40) <<"* Minimum ET for towers [GeV]: "<<""
+         << left << setw(20) <<ISOL_Calo_ET               <<""<< right << setw(10)<<"*"<<"\n";
+   f_out << left << setw(40) <<"* Cone for calorimetric isolation: "<<""
+         << left << setw(20) <<ISOL_Calo_Cone             <<""<< right << setw(10)<<"*"<<"\n";
+   f_out << left << setw(40) <<"* Grid size (NxN) for calorimetric isolation: "<<""
+         << left << setw(20) <<ISOL_Calo_Grid             <<""<< right << setw(10)<<"*"<<"\n";
+  }
+
+
   f_out<<"*                                                                    *"<<"\n";
   f_out<<"#***************                                                     *"<<"\n";
@@ -1072,21 +1096,48 @@
 //***********************Isolation criteria***********************
 //****************************************************************
-bool RESOLution::Isolation(const float phi, const float eta,const vector<TRootTracks> &tracks, const float& pt_second_track, const float& isolCone)
+bool RESOLution::Isolation(const D_Particle& part, const vector<TRootTracks> &tracks, const float& pt_second_track, const float& isolCone, float& ptiso )
 {
    bool isolated = false;
-   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_second_track)continue;
-       float genDeltaR = DeltaR(phi,eta,tracks[i].Phi,tracks[i].Eta); 
+   ptiso = 0; // sum of all track pt in isolation cone
+   float deltar=1E99; // Initial value; should be high; no further repercussion
+
+   // loop on all tracks, with p_t above threshold, close enough from the charged lepton
+   for(unsigned int i=0; i < tracks.size(); i++) {
+       if(tracks[i].PT < pt_second_track) continue; // ptcut on tracks
+       float genDeltaR = DeltaR(part.Phi(),part.Eta(),tracks[i].Phi,tracks[i].Eta); 
          if(
              (genDeltaR > deltar) ||
-             (genDeltaR==0)
+             (genDeltaR==0)  // rejets the track of the particle itself
            ) continue ;
-           deltar=genDeltaR;
+           deltar=genDeltaR;	// finds the closest track
+
+         // as long as (genDeltaR==0) is put above, the particle itself is not taken into account
+	 if( genDeltaR < ISOL_Cone) ptiso += tracks[i].PT; // dR cut on tracks
       }
    if(deltar > isolCone) isolated = true;
    return isolated;
+}
+
+// ******* Calorimetric isolation
+float RESOLution::CaloIsolation(const D_Particle& part,  const D_CaloTowerList & towers) {
+ // etrat, which is a percentage between 00 and 99. It is the ratio of the transverse energy 
+ // in a 3Ã3 grid surrounding the muon to the pT of the muon. For well-isolated muons, both ptiso and etrat will be small.
+   if(ISOL_Calo_ET>1E10) return UNDEFINED; // avoid doing anything unreasonable...
+   float etrat=0;
+   // available parameters: ISOL_Calo_ET , ISOL_Calo_Cone , 
+/*   for(unsigned int i=0; i < towers.size(); i++) {
+     if(towers[i].E > ISOL_Calo_ET) {
+        float genDeltaR = DeltaR(part.Phi(),part.Eta(),towers[i].getPhi(),towers[i].getEta());
+        if(genDeltaR < ISOL_Calo_Cone) {
+           ptiso += towers[i].getET();
+        }	    
+     }
+   } // loop on towers
+   ptiso -=
+*/
+  etrat = 100*etrat/part.Pt();
+  if(etrat<0) cout << "Error: negative etrat in CaloIsolation (" << etrat <<")\n";
+  else if(etrat>99) cout << "Error: etrat shoud be in [0;99] in CaloIsolation (" << etrat <<")\n";
+  return etrat; 
 }
 
