--- calo/flight/CaloEnergy/src/CaloEnergy.cpp	2009/07/13 16:10:23	1.3
+++ calo/flight/CaloEnergy/src/CaloEnergy.cpp	2009/08/11 13:28:27	1.5
@@ -35,7 +35,35 @@
   //
 }
 
+void CaloEnergy::UseCaloPreSampler(){
+  //
+  // use the presampler setting forcefitmode to 1000 means to force the DV routine to find the track inside the calorimeter using the "shower" approach developed for electrons
+  //
+  if ( !cp ) cp = new CaloPreSampler(L2);
+  cp->SplitInto(0,22);
+  cp->SetForceFitMode(1000);
+//   cp->UseTracker(false);
+//   cp->ForceCaloFit();
+//   cp->SetDebug(true);
+//   cp->Process();
+}
+
+
+void CaloEnergy::UseLongFit(){ 
+  fPl = 0; 
+  fLong = true;
+  if ( !clong ){
+    clong = new CaloLong(L2);
+    if ( cp ) clong->SetCaloLevel2Pointer(cp->GetLevel2Pointer());
+    clong->SplitInto(0,22);
+  };
+  //
+}
+
 void CaloEnergy::Set(){
+  //
+  // set default values, NB default conversion factor for energy is just very approximated!
+  //
   OBT = 0;
   PKT = 0;
   atime = 0;
@@ -48,9 +76,14 @@
   //
   debug = false;
   //
+  indep = false;
+  //
   fLong = false;
   fPl = 1;
   fRad = -1;
+  cp = NULL;
+  clong = NULL;
+  x0max = -1.;
   //
   this->DefineGeometry();
   fXosel =true;
@@ -75,6 +108,10 @@
   if ( section.Contains("YE") ) fYemin = plane;
 }
 
+void CaloEnergy::SetMinimumContainment(Int_t plane){
+  this->SetMinimumContainment("XEXOYEYO",plane);
+}
+
 Int_t CaloEnergy::GetMinimumContainment(TString section){
   section.ToUpper();
   if ( section.Contains("XO") ) return(fXomin);
@@ -93,6 +130,10 @@
   if ( section.Contains("YE") ) fConv_rye = conv;
 }
 
+void CaloEnergy::SetConversionFactor(Float_t conv){
+  this->SetConversionFactor("XEXOYEYO",conv);
+}
+
 Float_t CaloEnergy::GetConversionFactor(TString section){
   section.ToUpper();
   if ( section.Contains("XO") ) return(fConv_rxo);
@@ -112,92 +153,126 @@
   return(-1);
 }
 
+Float_t CaloEnergy::GetEnergyAtMaxplane(TString section){
+  section.ToUpper();
+  if ( section.Contains("XO") ) return xomax_en;
+  if ( section.Contains("XE") ) return xemax_en;
+  if ( section.Contains("YO") ) return yomax_en;
+  if ( section.Contains("YE") ) return yemax_en;
+  return(-1);
+}
+
+Float_t CaloEnergy::GetMaxEnergy(TString section){
+  section.ToUpper();
+  if ( fLong ){ 
+    this->Process(section);    
+    return fXOen_maxplane;
+  } else {
+    if ( section.Contains("XO") ) return fXOen_maxplane;
+    if ( section.Contains("XE") ) return fXEen_maxplane;
+    if ( section.Contains("YO") ) return fYOen_maxplane;
+    if ( section.Contains("YE") ) return fYEen_maxplane;
+  };
+  return(-1);
+}
+
+Float_t CaloEnergy::GetMaxEnergy(){   
+  if ( fLong ){
+    if ( debug ) printf(" oh! call process! with asntr %s and sntr %s \n",asntr.Data(),sntr.Data());
+    this->Process(asntr);
+  };
+  return((fXEen_maxplane+fYOen_maxplane+fYEen_maxplane+fXOen_maxplane));
+}
+
 void CaloEnergy::DefineGeometry(){
   //
+  // Use CaloStrip to determine once the position of border strips for each section
+  //
+  //
   fM = 2. + 0.096; // real position from cbar
   fM1 = 2. - 0.122 - 0.096; // due to calculation of xe1 etc.
   if ( fM1 < 0. ) fM1 = 0.;
   //
   CaloStrip *cs = new CaloStrip(fSimu);
-
-  // vista y piano 0 striscia 0
+  //
+  // view y plane 0 strip 0
   cs->Set(1,0,0);
   xe1= cs->GetY();
-  z1= cs->GetZ();
-  // vista y piano 0 striscia 31
+  //  z1= cs->GetZ();
+  // view y plane 0 strip 31
   cs->Set(1,0,31);
   xe2= cs->GetY();
-  // vista y piano 0 striscia 32
+  // view y plane 0 strip 32
   cs->Set(1,0,32);
   xe3= cs->GetY();
-  // vista y piano 0 striscia 63
+  // view y plane 0 strip 63
   cs->Set(1,0,63);
   xe4= cs->GetY();
-  // vista y piano 0 striscia 64
+  // view y plane 0 strip 64
   cs->Set(1,0,64);
   xe5= cs->GetY();
-  // vista y piano 0 striscia 95
+  // view y plane 0 strip 95
   cs->Set(1,0,95);
   xe6= cs->GetY();
 
-  // vista x piano 0 striscia 0
+  // view x plane 0 strip 0
   cs->Set(0,0,0);
   yo1= cs->GetX();
-  z2= cs->GetZ();
-  // vista x piano 0 striscia 31
+  //  z2= cs->GetZ();
+  // view x plane 0 strip 31
   cs->Set(0,0,31);
   yo2= cs->GetX();
-  // vista x piano 0 striscia 32
+  // view x plane 0 strip 32
   cs->Set(0,0,32);
   yo3= cs->GetX();
-  // vista x piano 0 striscia 63
+  // view x plane 0 strip 63
   cs->Set(0,0,63);
   yo4= cs->GetX();
-  // vista x piano 0 striscia 64
+  // view x plane 0 strip 64
   cs->Set(0,0,64);
   yo5= cs->GetX();
-  // vista x piano 0 striscia 95
+  // view x plane 0 strip 95
   cs->Set(0,0,95);
   yo6= cs->GetX();
 
   
-  // vista y piano 1 striscia 0
+  // view y plane 1 strip 0
   cs->Set(1,1,0);
   xo1= cs->GetY();
-  z3= cs->GetZ();
-  // vista y piano 1 striscia 31
+  //  z3= cs->GetZ();
+  // view y plane 1 strip 31
   cs->Set(1,1,31);
   xo2= cs->GetY();
-  // vista y piano 1 striscia 32
+  // view y plane 1 strip 32
   cs->Set(1,1,32);
   xo3= cs->GetY();
-  // vista y piano 1 striscia 63
+  // view y plane 1 strip 63
   cs->Set(1,1,63);
   xo4= cs->GetY();
-  // vista y piano 1 striscia 64
+  // view y plane 1 strip 64
   cs->Set(1,1,64);
   xo5= cs->GetY();
-  // vista y piano 1 striscia 95
+  // view y plane 1 strip 95
   cs->Set(1,1,95);
   xo6= cs->GetY();
   
-  // vista x piano 1 striscia 0
+  // view x plane 1 strip 0
   cs->Set(0,1,0);
   ye1= cs->GetX();
-  z4= cs->GetZ();
-  // vista x piano 1 striscia 31
+  //  z4= cs->GetZ();
+  // view x plane 1 strip 31
   cs->Set(0,1,31);
   ye2= cs->GetX();
-  // vista x piano 1 striscia 32
+  // view x plane 1 strip 32
   cs->Set(0,1,32);
   ye3= cs->GetX();
-  // vista x piano 1 striscia 63
+  // view x plane 1 strip 63
   cs->Set(0,1,63);
   ye4= cs->GetX();
-  // vista x piano 1 striscia 64
+  // view x plane 1 strip 64
   cs->Set(0,1,64);
   ye5= cs->GetX();
-  // vista x piano 1 striscia 95
+  // view x plane 1 strip 95
   cs->Set(0,1,95);
   ye6= cs->GetX();
 
@@ -205,7 +280,7 @@
   for (Int_t p = 0; p<22; p ++){
     for (Int_t v = 0; v<2; v++ ){   
       cs->Set(v,p,0);
-      trk_z[p][v]= cs->GetZ();                  // << cooordinata Z piano
+      trk_z[p][v]= cs->GetZ();                  //  Z coord for each plane
     };
   };
   //
@@ -215,6 +290,8 @@
 
 void CaloEnergy::Clear(){
   //
+  // clear variables
+  //
   fPartsel = false;
   fSel = false;
   fXosel = false;
@@ -231,8 +308,13 @@
   fMax_planexo = 0;
   fMax_planexe = 0;
   fMax_planeyo = 0;
-  fMax_planeye = 0;
-  memset(enstrip,0,2*22*96*(sizeof(Float_t)));  
+  fMax_planeye = 0; 
+  xomax_en= 0.; 
+  xemax_en= 0.; 
+  yomax_en= 0.; 
+  yemax_en= 0.; 
+  //
+ memset(enstrip,0,2*22*96*(sizeof(Float_t)));  
   en = 0.;
   view = 0;
   plane = 0;
@@ -245,6 +327,10 @@
   fYeout = 0;
   fXoout = 0;
   fXeout = 0;
+  fXEen_maxplane = 0.;  
+  fXOen_maxplane = 0.;  
+  fYEen_maxplane = 0.;  
+  fYOen_maxplane = 0.;  
   memset(en_xep,0,11*sizeof(Float_t));
   memset(en_yep,0,11*sizeof(Float_t));
   memset(en_xop,0,11*sizeof(Float_t));
@@ -291,6 +377,7 @@
   printf(" energyxo:.............. %f \n",energyxo); 
   printf(" energyye:.............. %f \n",energyye); 
   printf(" energyyo:.............. %f \n",energyyo); 
+  printf(" x0max   :.............. %f \n",x0max); 
   printf(" debug   :.............. %i \n",debug); 
 
   printf("========================================================================\n");
@@ -299,11 +386,13 @@
 
 void CaloEnergy::Delete(){
   Clear();
+  delete this;
 }
 
 
 Bool_t CaloEnergy::IsInsideAcceptance(TString section){
   //
+  // check if the event is inside the acceptance of the given section(s)
   //  
   TString ntr = section;
   if ( !L2 ){
@@ -326,37 +415,64 @@
     newentry = true;
   };
   //
+  // if we have already called this method for this event and no input changed then return fSel and exit
+  //
   if ( !newentry ) return fSel;
   //
+  // process the event
+  //
   if ( debug ) printf(" ########## IsInsideAcceptance ######### \n");
   //
+  // clear variables
+  //
   this->Clear();
   //
   section.ToUpper();
   //
-  Int_t ss = Int_t(section.Contains("XO"))+Int_t(section.Contains("XE"))+Int_t(section.Contains("YO"))+Int_t(section.Contains("YE"));
-  if ( !ss ){
+  // Count the number of section(s) given as input
+  //
+  Int_t fNumSec = Int_t(section.Contains("XO"))+Int_t(section.Contains("XE"))+Int_t(section.Contains("YO"))+Int_t(section.Contains("YE"));
+  if ( !fNumSec ){
     printf(" ERROR: section must be XO or XE or YO or YE while it is %s \n",section.Data());
     return false;
   };
   //
+  // If the presampler object exists then use the presampler output instead of the level2 output
+  //
+  CaloLevel2 *cl2 = NULL;
+  if ( cp ){
+    cl2 = cp->GetCaloLevel2();
+  } else {
+    cl2 = L2->GetCaloLevel2();
+  };
+  //
+  // get the energy for every strip of the calorimeter
   //
   for (Int_t ch=0; ch< L2->GetCaloLevel1()->istrip; ch++){
     en = L2->GetCaloLevel1()->DecodeEstrip(ch,view,plane,strip);
     enstrip[view][plane][strip]=en;
   };
+  //
+  // sum energy plane by plane for each sections 
+  //
   for (Int_t i=0;i<11;i++){
     for(strip=0; strip<96; strip++) {
       if ( fRad < 0 ){
+	//
+	// run over all the strips of the plane
+	//
 	en_xep[i] += enstrip[1][2*i][strip];
 	en_yop[i] += enstrip[0][2*i][strip];
-	en_xop[i] += enstrip[1][2*i+1][strip];
-	en_yep[i] += enstrip[0][2*i+1][strip];
+	en_xop[i] += enstrip[1][(2*i)+1][strip];
+	en_yep[i] += enstrip[0][(2*i)+1][strip];
       } else {
-	if ( strip >= L2->GetCaloLevel2()->cibar[2*i][1]-1-fRad &&  strip <= L2->GetCaloLevel2()->cibar[2*i][1]-1+fRad ) en_xep[i] += enstrip[1][2*i][strip]; 
-	if ( strip >= L2->GetCaloLevel2()->cibar[2*i][0]-1-fRad &&  strip <= L2->GetCaloLevel2()->cibar[2*i][0]-1+fRad ) en_yop[i] += enstrip[0][2*i][strip]; 
-	if ( strip >= L2->GetCaloLevel2()->cibar[2*i+1][1]-1-fRad &&  strip <= L2->GetCaloLevel2()->cibar[2*i+1][1]-1+fRad ) en_xop[i] += enstrip[1][2*i+1][strip]; 
-	if ( strip >= L2->GetCaloLevel2()->cibar[2*i+1][0]-1-fRad &&  strip <= L2->GetCaloLevel2()->cibar[2*i+1][0]-1+fRad ) en_yep[i] += enstrip[0][2*i+1][strip]; 
+	//
+	// use only the strips inside a cylinder of given radius fRad
+	//
+	if ( strip >= cl2->cibar[2*i][1]-1-fRad &&  strip <= cl2->cibar[2*i][1]-1+fRad ) en_xep[i] += enstrip[1][2*i][strip]; 
+	if ( strip >= cl2->cibar[2*i][0]-1-fRad &&  strip <= cl2->cibar[2*i][0]-1+fRad ) en_yop[i] += enstrip[0][2*i][strip]; 
+	if ( strip >= cl2->cibar[(2*i)+1][1]-1-fRad &&  strip <= cl2->cibar[(2*i)+1][1]-1+fRad ) en_xop[i] += enstrip[1][(2*i)+1][strip]; 
+	if ( strip >= cl2->cibar[(2*i)+1][0]-1-fRad &&  strip <= cl2->cibar[(2*i)+1][0]-1+fRad ) en_yep[i] += enstrip[0][(2*i)+1][strip]; 
       };
     };
     energyxe += en_xep[i];
@@ -364,88 +480,218 @@
     energyxo += en_xop[i];
     energyye += en_yep[i];
   };
-  Float_t en_maxplane = 0.;  
-  Float_t max_en= 0.; 
   //
-  for (Int_t ipl =0; ipl < 11; ipl ++) {
-    if(en_xop[ipl] > max_en) {
-      max_en = en_xop[ipl];
-      fMax_planexo = ipl;
+  // Find the plane of maximum for each section
+  //
+  //
+  Int_t xen = 0;
+  Int_t yon = 0;
+  Int_t xon = 0;
+  Int_t yen = 0;
+  //
+  if ( section.Contains("XE") ){
+    yon++;
+    xon++;
+    yen++;
+    for (Int_t ipl =0; ipl < 11; ipl ++) {
+      if(en_xep[ipl] > xemax_en) {
+	xemax_en = en_xep[ipl];
+	fMax_planexe = ipl;
       };
+    };
+  };
+  //
+  if ( section.Contains("YO") ){
+    xon++;
+    yen++;
+    for (Int_t ipl =0; ipl < 11; ipl ++) {
+      if(en_yop[ipl] > yomax_en) {
+	yomax_en = en_yop[ipl];
+	fMax_planeyo = ipl;
+      };
+    };
   };
   //
-  max_en = 0.;
-  for (Int_t ipl =0; ipl < 11; ipl ++) {
-    if(en_xep[ipl] > max_en) {
-      max_en = en_xep[ipl];
-      fMax_planexe = ipl;
+  if ( section.Contains("XO") ){
+    yen++;
+    for (Int_t ipl =0; ipl < 11; ipl ++) {
+      if(en_xop[ipl] > xomax_en) {
+	xomax_en = en_xop[ipl];
+	fMax_planexo = ipl;
+      };
     };
   };
   //
-  max_en = 0.;
-  for (Int_t ipl =0; ipl < 11; ipl ++) {
-    if(en_yop[ipl] > max_en) {
-      max_en = en_yop[ipl];
-      fMax_planeyo = ipl;
+  if ( section.Contains("YE") ){
+    for (Int_t ipl =0; ipl < 11; ipl ++) {
+      if(en_yep[ipl] > yemax_en) {
+	yemax_en = en_yep[ipl];
+	fMax_planeye = ipl;
+      };
     };
   };
   //
-  max_en = 0.;
-  for (Int_t ipl =0; ipl < 11; ipl ++) {
-    if(en_yep[ipl] > max_en) {
-      max_en = en_yep[ipl];
-      fMax_planeye = ipl;
+  // the maximum is given externally as X0, convert it to plane and section
+  //
+  if ( x0max > 0. ){
+    if ( debug ) printf(" CALCULATE MAX PLANE GIVEN X0 ASSUMING PERPENDICULAR TRACK \n");
+    Int_t wpl = (Int_t)roundf(x0max/0.76);
+    Bool_t isY = false;
+    if ( ((x0max/0.76)-(Float_t)wpl) > 0. ) isY = true;
+    if ( !(wpl%2) ){
+      // 0, 2, 4, ...
+      if ( isY ){
+	yomax_en = 1000.;
+	xemax_en = 500.;
+	fMax_planeyo=wpl/2;
+	fMax_planexe=wpl/2;
+      } else {
+	yomax_en = 500.;
+	xemax_en = 1000.;
+	fMax_planeyo=wpl/2;
+	fMax_planexe=wpl/2;
+      };
+      xomax_en = 0.;
+      yemax_en = 0.;
+    } else {
+      // 1, 3, 5, ...
+      if ( isY ){
+	yemax_en = 1000.;
+	xomax_en = 500.;
+	fMax_planeye=(wpl-1)/2;
+	fMax_planexo=(wpl-1)/2;
+      } else {
+	yemax_en = 500.;
+	xomax_en = 1000.;
+	fMax_planeye=(wpl-1)/2;
+	fMax_planexo=(wpl-1)/2;
+      };
+      xemax_en = 0.;
+      yomax_en = 0.;
     };
+    if ( debug ) printf(" x0max %f wpl %i isY %i yomax_en %f xemax_en %f yemax_en %f xomax_en %f fMaxplane %i %i %i %i\n",x0max,wpl,isY,yomax_en,xemax_en,yemax_en,xomax_en,fMax_planeyo,fMax_planexe,fMax_planeye,fMax_planexo);
   };
   //
+  Int_t nPl = fPl;
+  //
+  // Set the maximum in case of coherent mode was selected
+  //
+  if ( !indep ){
+    nPl = 0;
+    if ( debug ) printf(" A: Check maximum, coherent mode: xoen %f yoen %f xeen %f yeen %f xomax %i yomax %i xemax %i yemax %i\n",xomax_en,yomax_en,xemax_en,yemax_en,fMax_planexo,fMax_planeyo,fMax_planexe,fMax_planeye);
+    Int_t nummod = 0;
+    Int_t numexpl = 0;
+    if ( xomax_en > xemax_en && xomax_en > yemax_en && xomax_en > yomax_en ){
+      //
+      // Section XO contains the maximum energy release per plane of the whole calorimeter
+      //
+      if ( debug ) printf(" XO is MAX %i %i %i %i\n",xen,yon,xon,yen);
+      //
+      // fMax_plane is the plane of maximum + number of additional dE/dx measurement counting planes from 0 to 43
+      //
+      fMax_plane = (fNumSec * fMax_planexo) +(Float_t)xon + fPl;
+      //
+      // nummod is the integer part of the number of modules in which the maximum is contained
+      //
+      nummod = (Int_t)(((Float_t)fMax_plane)/(Float_t)fNumSec);
+      //
+      // numexpl is the number of additional planes (0,1,2) inside the module
+      //
+      numexpl = (Int_t)((Float_t)fMax_plane-(Float_t)fNumSec*(Float_t)nummod);
+      //
+    };
+    if ( xemax_en > xomax_en && xemax_en > yemax_en && xemax_en > yomax_en ){ 
+      if ( debug ) printf(" XE is MAX %i %i %i %i\n",xen,yon,xon,yen);
+      fMax_plane = (fNumSec * fMax_planexe) +(Float_t)xen + fPl;
+      nummod = (Int_t)(((Float_t)fMax_plane)/(Float_t)fNumSec);
+      numexpl = (Int_t)((Float_t)fMax_plane-(Float_t)fNumSec*(Float_t)nummod);
+      //
+   };
+
+    if ( yemax_en > xomax_en && yemax_en > xemax_en && yemax_en > yomax_en ){
+      if ( debug ) printf(" YE is MAX %i %i %i %i\n",xen,yon,xon,yen);
+      fMax_plane = (fNumSec * fMax_planeye) +(Float_t)yen + fPl;
+      nummod = (Int_t)(((Float_t)fMax_plane)/(Float_t)fNumSec);
+      numexpl = (Int_t)((Float_t)fMax_plane-(Float_t)fNumSec*(Float_t)nummod);
+      //
+    };
+    if ( yomax_en > xemax_en && yomax_en > yemax_en && yomax_en > xomax_en ){
+      if ( debug ) printf(" YO is MAX %i %i %i %i\n",xen,yon,xon,yen);
+      fMax_plane = (fNumSec * fMax_planeyo) +(Float_t)yon + fPl;
+      nummod = (Int_t)(((Float_t)fMax_plane)/(Float_t)fNumSec);
+      numexpl = (Int_t)((Float_t)fMax_plane-(Float_t)fNumSec*(Float_t)nummod);
+      //
+    };
+    //
+    // find the plane up to which is necessary to integrate the energy for each section
+    //
+    Int_t a = 0;
+    Int_t b = 0;
+    Int_t c = 0;
+    if ( numexpl > xen ) a = 1;
+    if ( numexpl > yon ) b = 1;
+    if ( numexpl > xon ) c = 1;
+    fMax_planexe = nummod;
+    fMax_planeyo = nummod - 1 + a;
+    fMax_planexo = nummod - 1 + b;
+    fMax_planeye = nummod - 1 + c;       
+    if ( debug ) printf(" fMax_plane %f nummod %i numexpl %i a %i b %i c %i \n",fMax_plane,nummod,numexpl,a,b,c);
+    if ( debug ) printf(" DONE: Check maximum, coherent mode: xoen %f yoen %f xeen %f yeen %f xomax %i yomax %i xemax %i yemax %i\n",xomax_en,yomax_en,xemax_en,yemax_en,fMax_planexo,fMax_planeyo,fMax_planexe,fMax_planeye);
+  };
+  //
+  // for each plane of the calorimeter find the position of the track in the direction along the strip (where we do not have a measurement from the selected plane) by looking at the plane above/below of the other view and extrapolating the trajectory to the given plane
+  //
   Float_t  track_coordx[22][2];
   Float_t  track_coordy[22][2];
   //
   Float_t tgx_cl2;
   Float_t tgy_cl2;  
-  tgx_cl2 = L2->GetCaloLevel2()->tanx[0];
-  tgy_cl2 = L2->GetCaloLevel2()->tany[0];
+  tgx_cl2 = cl2->tanx[0];
+  tgy_cl2 = cl2->tany[0];
   //
   for (Int_t p=0; p<22; p++){    
-    track_coordy[p][1] = L2->GetCaloLevel2()->cbar[p][1]; 
-    track_coordx[p][1] = L2->GetCaloLevel2()->cbar[p][0] - fabs(trk_z[p][1]-trk_z[p][0])*tgx_cl2;
-    track_coordx[p][0] = L2->GetCaloLevel2()->cbar[p][0];
-    track_coordy[p][0] = L2->GetCaloLevel2()->cbar[p][1] - fabs(trk_z[p][1]-trk_z[p][0])*tgy_cl2;
+    track_coordy[p][1] = cl2->cbar[p][1]; 
+    track_coordx[p][1] = cl2->cbar[p][0] - fabs(trk_z[p][1]-trk_z[p][0])*tgx_cl2;
+    track_coordx[p][0] = cl2->cbar[p][0];
+    track_coordy[p][0] = cl2->cbar[p][1] - fabs(trk_z[p][1]-trk_z[p][0])*tgy_cl2;
     if ( debug ) printf(" p %i track_coordy[p][1] %f track_coordx[p][1] %f track_coordx[p][0] %f  track_coordy[p][0] %f \n",p,track_coordy[p][1],track_coordx[p][1],track_coordx[p][0],track_coordy[p][0]);
   };
   //
-  if ( debug ) printf(" acceptance ss %i tgx %f tgy %f\n",ss,tgx_cl2,tgy_cl2);
+  if ( debug ) printf(" acceptance fNumSec %i tgx %f tgy %f\n",fNumSec,tgx_cl2,tgy_cl2);
   //
   if ( section.Contains("XO") ){
+    //
+    // check event is inside XO acceptance
+    //
     for (Int_t i=0; i<11; i++) {
       if
 	(
-	 (((track_coordx[2*i+1][1]>=(-12.054+fM))&&
-	   (track_coordx[2*i+1][1]<=(-4.246-fM)))&&
-	  (((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo1 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo2 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo3 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo4 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo5 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo6 - fM1 ))))||
+	 (((track_coordx[(2*i)+1][1]>=(-12.054+fM))&&
+	   (track_coordx[(2*i)+1][1]<=(-4.246-fM)))&&
+	  (((cl2->cbar[(2*i)+1][1]>=xo1 + fM1)&&
+	    (cl2->cbar[(2*i)+1][1]<=xo2 - fM1 ))||
+	   ((cl2->cbar[(2*i)+1][1]>=xo3 + fM1)&&
+	    (cl2->cbar[(2*i)+1][1]<=xo4 - fM1 ))||
+	   ((cl2->cbar[(2*i)+1][1]>=xo5 + fM1)&&
+	    (cl2->cbar[(2*i)+1][1]<=xo6 - fM1 ))))||
 	 
-	 (((track_coordx[2*i+1][1]>=(-4.004+fM))&&
-	   (track_coordx[2*i+1][1]<=(3.804-fM)))&&
-	  (((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo1 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo2 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo3 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo4 - fM1))||
-	   ((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo5 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo6 - fM1 ))))||
+	 (((track_coordx[(2*i)+1][1]>=(-4.004+fM))&&
+	   (track_coordx[(2*i)+1][1]<=(3.804-fM)))&&
+	  (((cl2->cbar[(2*i)+1][1]>=xo1 + fM1)&&
+	    (cl2->cbar[(2*i)+1][1]<=xo2 - fM1 ))||
+	   ((cl2->cbar[(2*i)+1][1]>=xo3 + fM1)&&
+	    (cl2->cbar[(2*i)+1][1]<=xo4 - fM1))||
+	   ((cl2->cbar[(2*i)+1][1]>=xo5 + fM1)&&
+	    (cl2->cbar[(2*i)+1][1]<=xo6 - fM1 ))))||
 	 
-	 (((track_coordx[2*i+1][1]>=(4.046+fM))&&
-	   (track_coordx[2*i+1][1]<=(11.854-fM)))&&
-	  (((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo1 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo2 - fM1))||
-	   ((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo3 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo4 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i+1][1]>=xo5 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][1]<=xo6 - fM1 ))))
+	 (((track_coordx[(2*i)+1][1]>=(4.046+fM))&&
+	   (track_coordx[(2*i)+1][1]<=(11.854-fM)))&&
+	  (((cl2->cbar[(2*i)+1][1]>=xo1 + fM1)&&
+	    (cl2->cbar[(2*i)+1][1]<=xo2 - fM1))||
+	   ((cl2->cbar[(2*i)+1][1]>=xo3 + fM1)&&
+	    (cl2->cbar[(2*i)+1][1]<=xo4 - fM1 ))||
+	   ((cl2->cbar[(2*i)+1][1]>=xo5 + fM1)&&
+	    (cl2->cbar[(2*i)+1][1]<=xo6 - fM1 ))))
 	 ){
 	fXosel = true;
 	fXoout = i;
@@ -454,15 +700,20 @@
 	break;
       };
     };
-    if ( !fXosel && fXoout >= fXomin && fXoout >= (Int_t)(fMax_planexo+fPl) ){
-      if ( debug ) printf(" XO: this event is only partially contained: fXoout %i fXomin %i fMax_planexo + fPl %i \n",fXoout,fXomin,(Int_t)(fMax_planexo+fPl));
+    //
+    // if it goes out of the acceptance BUT the plane up to which we are integrating the energy is contained then the event is "partially" contained
+    //
+    if ( !fXosel && fXoout >= fXomin && fXoout >= (Int_t)(fMax_planexo+nPl) ){
+      if ( debug ) printf(" XO: this event is only partially contained: fXoout %i fXomin %i fMax_planexo + nPl %i \n",fXoout,fXomin,(Int_t)(fMax_planexo+nPl));
       fPartsel = true;
       fXosel = true;
     };
+    //
+    // event is contained (or partially contained) hence we can integrate energy up to the maximum and calculate the energy as measured by this section
+    //
     if ( fXosel ){
-      en_maxplane = 0.;  
-      for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planexo+fPl)) ;iplm++) en_maxplane += en_xop[iplm]; 
-      fEnergyxo = en_maxplane/fConv_rxo;
+      for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planexo+nPl)) ;iplm++) fXOen_maxplane += en_xop[iplm]; 
+      fEnergyxo = fXOen_maxplane/fConv_rxo;
     };
   };
   //
@@ -472,30 +723,30 @@
 	(
 	 (((track_coordx[2*i][1]>=(-11.854+fM))&& 
 	   (track_coordx[2*i][1]<=(-4.046-fM)))&&
-	  (((L2->GetCaloLevel2()->cbar[2*i][1]>=xe1 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][1]<=xe2 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i][1]>=xe3 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][1]<=xe4 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i][1]>=xe5 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][1]<=xe6 - fM1 ))))||
+	  (((cl2->cbar[2*i][1]>=xe1 + fM1)&&
+	    (cl2->cbar[2*i][1]<=xe2 - fM1 ))||
+	   ((cl2->cbar[2*i][1]>=xe3 + fM1)&&
+	    (cl2->cbar[2*i][1]<=xe4 - fM1 ))||
+	   ((cl2->cbar[2*i][1]>=xe5 + fM1)&&
+	    (cl2->cbar[2*i][1]<=xe6 - fM1 ))))||
 	 
 	 (((track_coordx[2*i][1]>=(-3.804+fM))&&
 	   (track_coordx[2*i][1]<=(4.004-fM)))&&
-	  (((L2->GetCaloLevel2()->cbar[2*i][1]>=xe1 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][1]<=xe2 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i][1]>=xe3 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][1]<=xe4 - fM1))||
-	   ((L2->GetCaloLevel2()->cbar[2*i][1]>=xe5 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][1]<=xe6 - fM1 ))))||
+	  (((cl2->cbar[2*i][1]>=xe1 + fM1)&&
+	    (cl2->cbar[2*i][1]<=xe2 - fM1 ))||
+	   ((cl2->cbar[2*i][1]>=xe3 + fM1)&&
+	    (cl2->cbar[2*i][1]<=xe4 - fM1))||
+	   ((cl2->cbar[2*i][1]>=xe5 + fM1)&&
+	    (cl2->cbar[2*i][1]<=xe6 - fM1 ))))||
 	 
 	 (((track_coordx[2*i][1]>=(4.246+fM))&&
 	   (track_coordx[2*i][1]<=(12.054-fM)))&&
-	  (((L2->GetCaloLevel2()->cbar[2*i][1]>=xe1 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][1]<=xe2 - fM1))||
-	   ((L2->GetCaloLevel2()->cbar[2*i][1]>=xe3 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][1]<=xe4 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i][1]>=xe5 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][1]<=xe6 - fM1 ))))
+	  (((cl2->cbar[2*i][1]>=xe1 + fM1)&&
+	    (cl2->cbar[2*i][1]<=xe2 - fM1))||
+	   ((cl2->cbar[2*i][1]>=xe3 + fM1)&&
+	    (cl2->cbar[2*i][1]<=xe4 - fM1 ))||
+	   ((cl2->cbar[2*i][1]>=xe5 + fM1)&&
+	    (cl2->cbar[2*i][1]<=xe6 - fM1 ))))
 	 ){
 	fXesel = true;
 	fXeout = i;
@@ -504,15 +755,14 @@
 	break;
       };
     };
-    if ( !fXesel && fXeout >= fXemin && fXeout >= (Int_t)(fMax_planexe+fPl) ){
-      if ( debug ) printf(" XE: this event is only partially contained: fXeout %i fXemin %i fMax_planexe + fPl %i \n",fXeout,fXemin,(Int_t)(fMax_planexe+fPl));
+    if ( !fXesel && fXeout >= fXemin && fXeout >= (Int_t)(fMax_planexe+nPl) ){
+      if ( debug ) printf(" XE: this event is only partially contained: fXeout %i fXemin %i fMax_planexe + nPl %i \n",fXeout,fXemin,(Int_t)(fMax_planexe+nPl));
       fPartsel = true;
       fXesel = true;
     };
     if ( fXesel ){
-      en_maxplane = 0.;  
-      for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planexe+fPl)) ;iplm++) en_maxplane += en_xep[iplm]; 
-      fEnergyxe = en_maxplane/fConv_rxe;
+      for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planexe+nPl)) ;iplm++) fXEen_maxplane += en_xep[iplm]; 
+      fEnergyxe = fXEen_maxplane/fConv_rxe;
     };
   };  
   //
@@ -520,32 +770,32 @@
     for (Int_t i=0; i<11; i++) {
       if
 	(
-	 (((track_coordy[2*i+1][0]>=(-12.154+fM))&& 
-	   (track_coordy[2*i+1][0]<=(-4.346-fM)))&&
-	  (((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye1 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye2 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye3 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye4 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye5 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye6 - fM1 ))))||
+	 (((track_coordy[(2*i)+1][0]>=(-12.154+fM))&& 
+	   (track_coordy[(2*i)+1][0]<=(-4.346-fM)))&&
+	  (((cl2->cbar[(2*i)+1][0]>=ye1 + fM1)&&
+	    (cl2->cbar[(2*i)+1][0]<=ye2 - fM1 ))||
+	   ((cl2->cbar[(2*i)+1][0]>=ye3 + fM1)&&
+	    (cl2->cbar[(2*i)+1][0]<=ye4 - fM1 ))||
+	   ((cl2->cbar[(2*i)+1][0]>=ye5 + fM1)&&
+	    (cl2->cbar[(2*i)+1][0]<=ye6 - fM1 ))))||
 	 
-	 (((track_coordy[2*i+1][0]>=(-4.104+fM))&&
-	   (track_coordy[2*i+1][0]<=(3.704-fM)))&&
-	  (((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye1 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye2 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye3 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye4 - fM1))||
-	   ((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye5 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye6 - fM1 ))))||
+	 (((track_coordy[(2*i)+1][0]>=(-4.104+fM))&&
+	   (track_coordy[(2*i)+1][0]<=(3.704-fM)))&&
+	  (((cl2->cbar[(2*i)+1][0]>=ye1 + fM1)&&
+	    (cl2->cbar[(2*i)+1][0]<=ye2 - fM1 ))||
+	   ((cl2->cbar[(2*i)+1][0]>=ye3 + fM1)&&
+	    (cl2->cbar[(2*i)+1][0]<=ye4 - fM1))||
+	   ((cl2->cbar[(2*i)+1][0]>=ye5 + fM1)&&
+	    (cl2->cbar[(2*i)+1][0]<=ye6 - fM1 ))))||
 	 
-	 (((track_coordy[2*i+1][0]>=(3.946+fM))&&
-	   (track_coordy[2*i+1][0]<=(11.754-fM)))&&
-	  (((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye1 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye2 - fM1))||
-	   ((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye3 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye4 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i+1][0]>=ye5 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i+1][0]<=ye6 - fM1 ))))
+	 (((track_coordy[(2*i)+1][0]>=(3.946+fM))&&
+	   (track_coordy[(2*i)+1][0]<=(11.754-fM)))&&
+	  (((cl2->cbar[(2*i)+1][0]>=ye1 + fM1)&&
+	    (cl2->cbar[(2*i)+1][0]<=ye2 - fM1))||
+	   ((cl2->cbar[(2*i)+1][0]>=ye3 + fM1)&&
+	    (cl2->cbar[(2*i)+1][0]<=ye4 - fM1 ))||
+	   ((cl2->cbar[(2*i)+1][0]>=ye5 + fM1)&&
+	    (cl2->cbar[(2*i)+1][0]<=ye6 - fM1 ))))
 	 ){
 	fYesel = true;
 	fYeout = i;
@@ -554,50 +804,49 @@
 	break;
       };
     };
-    if ( !fYesel && fYeout >= fYemin && fYeout >= (Int_t)(fMax_planeye+fPl) ){
-      if ( debug ) printf(" YE: this event is only partially contained: fYeout %i fYemin %i fMax_planeye + fPl %i \n",fYeout,fYemin,(Int_t)(fMax_planeye+fPl));
+    if ( !fYesel && fYeout >= fYemin && fYeout >= (Int_t)(fMax_planeye+nPl) ){
+      if ( debug ) printf(" YE: this event is only partially contained: fYeout %i fYemin %i fMax_planeye + nPl %i \n",fYeout,fYemin,(Int_t)(fMax_planeye+nPl));
       fPartsel = true;
       fYesel = true;
     };
     if ( fYesel ){
-      en_maxplane = 0.;  
-      for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planeye+fPl)) ;iplm++) en_maxplane += en_yep[iplm]; 
-      fEnergyye = en_maxplane/fConv_rye;
+      for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planeye+nPl)) ;iplm++) fYEen_maxplane += en_yep[iplm]; 
+      fEnergyye = fYEen_maxplane/fConv_rye;
     };
   };  
   //
   if ( section.Contains("YO") ){
     for (Int_t i=0; i<11; i++) {
-      if ( debug ) printf(" i %i track_coordy[2*i][0] %f L2->GetCaloLevel2()->cbar[2*i][0] %f \n",i,track_coordy[2*i][0],L2->GetCaloLevel2()->cbar[2*i][0]);
+      if ( debug ) printf(" i %i track_coordy[2*i][0] %f cl2->cbar[2*i][0] %f \n",i,track_coordy[2*i][0],cl2->cbar[2*i][0]);
       if ( debug ) printf(" i %i fm %f fm1 %f yo1 %g yo2 %f yo3 %f yo4 %f yo5 %f yo6 %f \n",i,fM,fM1,yo1,yo2,yo3,yo4,yo5,yo6);
       if
 	(
 	 (((track_coordy[2*i][0]>=(-11.954+fM))&& 
 	   (track_coordy[2*i][0]<=(-4.146-fM)))&&
-	  (((L2->GetCaloLevel2()->cbar[2*i][0]>=yo1 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][0]<=yo2 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i][0]>=yo3 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][0]<=yo4 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i][0]>=yo5 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][0]<=yo6 - fM1 ))))||
+	  (((cl2->cbar[2*i][0]>=yo1 + fM1)&&
+	    (cl2->cbar[2*i][0]<=yo2 - fM1 ))||
+	   ((cl2->cbar[2*i][0]>=yo3 + fM1)&&
+	    (cl2->cbar[2*i][0]<=yo4 - fM1 ))||
+	   ((cl2->cbar[2*i][0]>=yo5 + fM1)&&
+	    (cl2->cbar[2*i][0]<=yo6 - fM1 ))))||
 	 
 	 (((track_coordy[2*i][0]>=(-3.904+fM))&&
 	   (track_coordy[2*i][0]<=(+3.904-fM)))&&
-	  (((L2->GetCaloLevel2()->cbar[2*i][0]>=yo1 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][0]<=yo2 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i][0]>=yo3 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][0]<=yo4 - fM1))||
-	   ((L2->GetCaloLevel2()->cbar[2*i][0]>=yo5 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][0]<=yo6 - fM1 ))))||
+	  (((cl2->cbar[2*i][0]>=yo1 + fM1)&&
+	    (cl2->cbar[2*i][0]<=yo2 - fM1 ))||
+	   ((cl2->cbar[2*i][0]>=yo3 + fM1)&&
+	    (cl2->cbar[2*i][0]<=yo4 - fM1))||
+	   ((cl2->cbar[2*i][0]>=yo5 + fM1)&&
+	    (cl2->cbar[2*i][0]<=yo6 - fM1 ))))||
 	 
 	 (((track_coordy[2*i][0]>=(4.146+fM))&&
 	   (track_coordy[2*i][0]<=(11.954-fM)))&&
-	  (((L2->GetCaloLevel2()->cbar[2*i][0]>=yo1 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][0]<=yo2 - fM1))||
-	   ((L2->GetCaloLevel2()->cbar[2*i][0]>=yo3 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][0]<=yo4 - fM1 ))||
-	   ((L2->GetCaloLevel2()->cbar[2*i][0]>=yo5 + fM1)&&
-	    (L2->GetCaloLevel2()->cbar[2*i][0]<=yo6 - fM1 ))))
+	  (((cl2->cbar[2*i][0]>=yo1 + fM1)&&
+	    (cl2->cbar[2*i][0]<=yo2 - fM1))||
+	   ((cl2->cbar[2*i][0]>=yo3 + fM1)&&
+	    (cl2->cbar[2*i][0]<=yo4 - fM1 ))||
+	   ((cl2->cbar[2*i][0]>=yo5 + fM1)&&
+	    (cl2->cbar[2*i][0]<=yo6 - fM1 ))))
 	 ){
 	fYosel = true;
 	fYoout = i;
@@ -606,24 +855,45 @@
 	break;
       };
     };
-    if ( !fYosel && fYoout >= fYomin && fYoout >= (Int_t)(fMax_planeyo+fPl) ){
-      if ( debug ) printf(" YO: this event is only partially contained: fYoout %i fYomin %i fMax_planeyo + fPl %i \n",fYoout,fYomin,(Int_t)(fMax_planeyo+fPl));
+    if ( !fYosel && fYoout >= fYomin && fYoout >= (Int_t)(fMax_planeyo+nPl) ){
+      if ( debug ) printf(" YO: this event is only partially contained: fYoout %i fYomin %i fMax_planeyo + nPl %i \n",fYoout,fYomin,(Int_t)(fMax_planeyo+nPl));
       fPartsel = true;
       fYosel = true;
     };
     if ( fYosel ){
-      en_maxplane = 0.;  
-      for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planeyo+fPl)) ;iplm++) en_maxplane += en_yop[iplm]; 
-      fEnergyyo = en_maxplane/fConv_ryo;
+      for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planeyo+nPl)) ;iplm++) fYOen_maxplane += en_yop[iplm]; 
+      fEnergyyo = fYOen_maxplane/fConv_ryo;
     };
   };  
   //
-  fSel = ( fXesel || fYesel || fXosel || fYosel );
+  // Count the number of sections in which the event is contained
+  //
   fCount = (Float_t)((Int_t)fXesel+(Int_t)fXosel+(Int_t)fYesel+(Int_t)fYosel);
+  //
+  if ( indep ){
+    //
+    // independent mode, average the energy measurement and max plane of the contained sections
+    //
+    fSel = ( fXesel || fYesel || fXosel || fYosel );
+    fMax_plane = (Float_t)(fMax_planeyo+fMax_planeye+fMax_planexo+fMax_planexe)/fCount;
+    fEnergy = (fEnergyxe+fEnergyyo+fEnergyye+fEnergyxo)/fCount;
+    //
+  } else {
+    //
+    // coherent mode, sum the energy [MIP] of the given sections and convert using fConv_rxo. **** NB: it is assumed that the conversion factor is unique and the method SetConvertionFactor(Float_t) has been used**** The event is selected only if it is contained in all the given sections
+    //
+    if ( fCount != fNumSec ){
+      fSel = false;
+    } else {
+      fSel = true;
+    };
+    fEnergy = (fXEen_maxplane+fYOen_maxplane+fYEen_maxplane+fXOen_maxplane)/fConv_rxo;
+  };
+  //
+  if ( debug ) printf("sel %i indep %i fMax_plane %f conv_r %f en_maxplane %f encalo %f  \n",fSel,indep,fMax_plane,fConv_rxo,fXOen_maxplane,fEnergy);
   if ( debug ) printf(" IsInside XE %i XO %i YE %i YO %i => SEL %i \n",fXesel,fXosel,fYesel,fYosel,fSel);
-  fMax_plane = (Float_t)(fMax_planeyo+fMax_planeye+fMax_planexo+fMax_planexe)/fCount;
-  fEnergy = (fEnergyxe+fEnergyyo+fEnergyye+fEnergyxo)/fCount;
-  if ( debug ) printf(" fMax_plane %f conv_r %f en_maxplane %f encalo %f  \n",fMax_plane,fConv_rxo,en_maxplane,fEnergy);
+  //
+  // finish exit
   //
   return fSel;
   // 
@@ -637,6 +907,8 @@
 
 void CaloEnergy::Process(TString section){
   //  
+  // process the event
+  //
   TString ntr = section;
   if ( !L2 ){
     printf(" ERROR: cannot find PamLevel2 object, use the correct constructor or check your program!\n");
@@ -658,25 +930,81 @@
     newentry = true;
   };
   //
+  // if we have already called this method for this event and no input changed then return fSel and exit
+  //
   if ( !newentry ) return;
   //
+  // process the event
+  //
   if ( debug ) printf(" Processing event at OBT %u PKT %u time %u section %s\n",OBT,PKT,atime,section.Data());
   //
-  //  this->Clear();
+  // check if the cylinder of integration can go out of the sensor given the frame which has been set (if we use all the calorimeter fRad is < 0 and the printout is suppressed)
   //
   if ( (fM1+0.122-0.244*(Float_t)fRad) < 0. ) printf("Error: (fM1+0.122-0.244*(Float_t)fRad) < 0. fM1 %f fRad %i %f \n",fM1,fRad,(fM1+0.122-0.244*(Float_t)fRad));
   //
   if ( fLong ){
-    // use long fit to measure energy
+    if ( debug ) printf(" ==================================================================> LONGITUDINAL FIT! \n"); 
+    //
+    // use long fit to measure energy 
+    //
+    if ( this->IsInsideAcceptance(section) ){
+      //
+      if ( debug ) printf(" ==================================================================> LONG INSIDE! \n"); 
+      //
+      Float_t myene[2][22];
+      memset(myene,0,(sizeof(Float_t))*2*22);
+      for (Int_t j=0; j<11; j++){
+	if ( section.Contains("XE") ) myene[1][2*j] = en_xep[j];
+	if ( section.Contains("YO") ) myene[0][2*j] = en_yop[j];
+	if ( section.Contains("XO") ) myene[1][(2*j)+1] = en_xop[j];
+	if ( section.Contains("YE") ) myene[0][(2*j)+1] = en_yep[j];	
+      };
+      clong->UnMaskSections();
+      if ( !(section.Contains("YE")) ) clong->MaskSection("YE");
+      if ( !(section.Contains("YO")) ) clong->MaskSection("YO");
+      if ( !(section.Contains("XO")) ) clong->MaskSection("XO");
+      if ( !(section.Contains("XE")) ) clong->MaskSection("XE");
+      clong->ForceNextFit();
+      clong->SetEnergies(myene);
+      if ( debug ){
+	clong->Fit(true);
+      } else {
+	clong->Fit();
+      };      
+      fXOen_maxplane = clong->Get_E0();
+      fYOen_maxplane  = 0.;
+      fYEen_maxplane = 0.; 
+      fXEen_maxplane = 0.;
+      fEnergy=fXOen_maxplane/fConv_rxo;
+      if ( fEnergy != fEnergy || clong->Get_fitresult() != 0 ) fEnergy = -1.;
+      //      if ( fEnergy != fEnergy ) fEnergy = 1.;
+      //
+    } else {
+      //
+      // if the event is not in the acceptance, return a negative energy.
+      //
+      if ( debug ) printf(" Outside acceptance \n");
+      fEnergy *= -1.;
+      //
+    };
+    //
   } else {
     //
+    // use the energy measurement
+    //
     if ( this->IsInsideAcceptance(section) ){
       //
+      // the event is good
+      //
       if ( debug ) printf(" XE %i XO %i YE %i YO %i \n",fXesel,fXosel,fYesel,fYosel);
       //
     } else {
+      //
+      // if the event is not in the acceptance, return a negative energy.
+      //
       if ( debug ) printf(" Outside acceptance \n");
-      fEnergy = -1.;
+      fEnergy *= -1.;
+      //
     };
   };
   //