diff --git a/NPLib/Core/NPDetectorManager.cxx b/NPLib/Core/NPDetectorManager.cxx
index 0b000bc0b6c24a3aa33b0be6b732b639c41145a1..f8a1933c8f7c8712f0598ce13ce3dfc5042719de 100644
--- a/NPLib/Core/NPDetectorManager.cxx
+++ b/NPLib/Core/NPDetectorManager.cxx
@@ -413,7 +413,7 @@ void NPL::DetectorManager::InitThreadPool(){
 ////////////////////////////////////////////////////////////////////////////////
 void NPL::DetectorManager::StartThread(NPL::VDetector* det,unsigned int id){ 
   // Let the main thread start 
-  std::this_thread::sleep_for(std::chrono::milliseconds(100));
+  std::this_thread::sleep_for(std::chrono::milliseconds(250));
   while(true){
     // Do the job if possible
     if(m_Ready[id]){
diff --git a/NPLib/Detectors/Nebula/TNebulaData.h b/NPLib/Detectors/Nebula/TNebulaData.h
index 9fcaf75142c032ccc7768beff2b45d3d9a116811..54eaea6e3fd846bca0b4665c72960171d294e31d 100644
--- a/NPLib/Detectors/Nebula/TNebulaData.h
+++ b/NPLib/Detectors/Nebula/TNebulaData.h
@@ -100,24 +100,44 @@ class TNebulaData : public TObject {
     fNebula_Td_ID.push_back(ID);
     fNebula_Td_Time.push_back(Time);
     };//!
-/*
-    //////////////////////    GETTERS    ////////////////////////
-    // Energy
-    inline UShort_t GetMultEnergy() const
-      {return fNebula_E_DetectorNbr.size();}
-    inline UShort_t GetE_DetectorNbr(const unsigned int &i) const 
-      {return fNebula_E_DetectorNbr[i];}//!
-    inline Double_t Get_Energy(const unsigned int &i) const 
-      {return fNebula_Energy[i];}//!
 
+    //////////////////////    GETTERS    ////////////////////////
+    // MULT //
+    // Charge 
+    inline unsigned int GetChargeUpMult() const
+      {return fNebula_Qu_ID.size();};
+    // Time
+    inline unsigned int GetTimeUpMult() const
+      {return fNebula_Tu_ID.size();};
+    // Charge
+    inline unsigned int GetChargeDownMult() const
+      {return fNebula_Qd_ID.size();};
     // Time
-    inline UShort_t GetMultTime() const
-      {return fNebula_T_DetectorNbr.size();}
-    inline UShort_t GetT_DetectorNbr(const unsigned int &i) const 
-      {return fNebula_T_DetectorNbr[i];}//!
-    inline Double_t Get_Time(const unsigned int &i) const 
-      {return fNebula_Time[i];}//!
-*/
+    inline unsigned int GetTimeDownMult() const
+      {return fNebula_Td_ID.size();};
+
+    // Value // 
+    // Charge 
+    inline UShort_t GetChargeUpID(unsigned int& i) const
+      {return fNebula_Qu_ID[i];};
+    inline double GetChargeUp(unsigned int& i) const
+      {return fNebula_Qu_Charge[i];};
+    // Time 
+    inline UShort_t GetTimeUpID(unsigned int& i) const
+      {return fNebula_Tu_ID[i];};
+    inline double GetTimeUp(unsigned int& i) const
+      {return fNebula_Tu_Time[i];};
+    // Charge 
+    inline UShort_t GetChargeDownID(unsigned int& i) const
+      {return fNebula_Qd_ID[i];};
+    inline double GetChargeDown(unsigned int& i) const
+      {return fNebula_Qd_Charge[i];};
+    // Time 
+    inline UShort_t GetTimeDownID(unsigned int& i) const
+      {return fNebula_Td_ID[i];};
+    inline double GetTimeDown(unsigned int& i) const
+      {return fNebula_Td_Time[i];};
+
 
   //////////////////////////////////////////////////////////////
   // Required for ROOT dictionnary
diff --git a/NPLib/Detectors/Nebula/TNebulaPhysics.cxx b/NPLib/Detectors/Nebula/TNebulaPhysics.cxx
index 0664800e9952cbcba7c297f729a26294efc510e2..34490e56a8ec1ca22fabe50068b235002a83c631 100644
--- a/NPLib/Detectors/Nebula/TNebulaPhysics.cxx
+++ b/NPLib/Detectors/Nebula/TNebulaPhysics.cxx
@@ -35,28 +35,61 @@ using namespace std;
 #include "RootOutput.h"
 #include "NPDetectorFactory.h"
 #include "NPOptionManager.h"
-
+#include "NPSystemOfUnits.h"
 //   ROOT
 #include "TChain.h"
 
 ClassImp(TNebulaPhysics)
 
 
-///////////////////////////////////////////////////////////////////////////
+  ///////////////////////////////////////////////////////////////////////////
 TNebulaPhysics::TNebulaPhysics()
-   : m_EventData(new TNebulaData),
-     m_PreTreatedData(new TNebulaData),
-     m_EventPhysics(this),
-     m_Spectra(0),
-     m_E_RAW_Threshold(0), // adc channels
-     m_E_Threshold(0),     // MeV
-     m_NumberOfBars(0) {
-}
+  : m_EventData(new TNebulaData),
+  m_EventPhysics(this),
+  m_Spectra(0),
+  m_Q_RAW_Threshold(0), // adc channels
+  m_Q_Threshold(7),     // normal bars in MeV
+  m_V_Threshold(1),     // veto bars in MeV
+  m_NumberOfBars(0) {
+  }
 
 ///////////////////////////////////////////////////////////////////////////
 /// A usefull method to bundle all operation to add a detector
 void TNebulaPhysics::ReadXML(NPL::XmlParser xml){ 
-  m_NumberOfBars++;
+  std::vector<NPL::XML::block*> b = xml.GetAllBlocksWithName("NEBULA");  
+
+  for(unsigned int i = 0 ; i < b.size() ; i++){
+    m_NumberOfBars++;
+    unsigned int id = b[i]->AsInt("ID");
+
+    // position
+    PositionX[id] = b[i]->AsDouble("PosX"); 
+    PositionY[id] = b[i]->AsDouble("PosY"); 
+    PositionZ[id] = b[i]->AsDouble("PosZ"); 
+    // linear cal
+    aQu[id] = b[i]->AsDouble("QUCal");
+    bQu[id] = b[i]->AsDouble("QUPed");
+    aQd[id] = b[i]->AsDouble("QDCal");
+    bQd[id] = b[i]->AsDouble("QDPed");
+    aTu[id] = b[i]->AsDouble("TUCal");
+    bTu[id] = b[i]->AsDouble("TUOff");
+    aTd[id] = b[i]->AsDouble("TDCal");
+    bTd[id] = b[i]->AsDouble("TDOff");
+
+    // T average offset
+    avgT0[id] = b[i]->AsDouble("TAveOff");
+
+    // slew correction T= tcal +slwT/sqrt(Qcal)
+    slwTu[id] = b[i]->AsDouble("TUSlw");
+    slwTd[id] = b[i]->AsDouble("TDSlw");
+
+    // DT position cal
+    DTa[id] = b[i]->AsDouble("DTCal");//!
+    DTb[id] = b[i]->AsDouble("DTOff");//!
+
+
+  } 
+  cout << " -> " << m_NumberOfBars << " bars found" << endl;;
 } 
 
 ///////////////////////////////////////////////////////////////////////////
@@ -69,115 +102,109 @@ void TNebulaPhysics::BuildSimplePhysicalEvent() {
 ///////////////////////////////////////////////////////////////////////////
 void TNebulaPhysics::BuildPhysicalEvent() {
   // apply thresholds and calibration
-  PreTreat();
-/*
-  // match energy and time together
-  unsigned int mysizeE = m_PreTreatedData->GetMultEnergy();
-  unsigned int mysizeT = m_PreTreatedData->GetMultTime();
-  for (UShort_t e = 0; e < mysizeE ; e++) {
-    for (UShort_t t = 0; t < mysizeT ; t++) {
-      if (m_PreTreatedData->GetE_DetectorNbr(e) == m_PreTreatedData->GetT_DetectorNbr(t)) {
-        DetectorNumber.push_back(m_PreTreatedData->GetE_DetectorNbr(e));
-        Energy.push_back(m_PreTreatedData->Get_Energy(e));
-        Time.push_back(m_PreTreatedData->Get_Time(t));
-      }
-    }
-  }*/
-}
-
-///////////////////////////////////////////////////////////////////////////
-void TNebulaPhysics::PreTreat() {
-  // This method typically applies thresholds and calibrations
-  // Might test for disabled channels for more complex detector
-
-  // clear pre-treated object
-  ClearPreTreatedData();
 
   // instantiate CalibrationManager
   static CalibrationManager* Cal = CalibrationManager::getInstance();
-/*
-  // Energy
-  unsigned int mysize = m_EventData->GetMultEnergy();
-  for (UShort_t i = 0; i < mysize ; ++i) {
-    if (m_EventData->Get_Energy(i) > m_E_RAW_Threshold) {
-      Double_t Energy = Cal->ApplyCalibration("Nebula/ENERGY"+NPL::itoa(m_EventData->GetE_DetectorNbr(i)),m_EventData->Get_Energy(i));
-      if (Energy > m_E_Threshold) {
-        m_PreTreatedData->SetEnergy(m_EventData->GetE_DetectorNbr(i), Energy);
+  static double rawQup,calQup,rawQdown,calQdown,rawTup,calTup,rawTdown,calTdown,calQ,calT,Y;
+  static unsigned int ID;
+  // All vector size 
+  static unsigned int QUsize, QDsize, TUsize, TDsize ; 
+  QUsize = m_EventData->GetChargeUpMult();
+  QDsize = m_EventData->GetChargeDownMult();
+  TUsize = m_EventData->GetTimeUpMult();
+  TDsize = m_EventData->GetTimeDownMult();
+  static double threshold;
+  // loop on Qup
+  for (unsigned int qup = 0; qup < QUsize ; qup++) {
+
+    rawQup = m_EventData->GetChargeUp(qup);
+    rawTup=-1;
+    rawQdown=-1;
+    rawTdown=-1;
+    if (rawQup > m_Q_RAW_Threshold) {
+      ID = m_EventData->GetChargeUpID(qup);
+      if(ID<121)
+        threshold=m_Q_Threshold;
+      else
+        threshold=m_V_Threshold;
+
+      // look for associated Charge down
+      for(unsigned int qdown = 0 ; qdown < QDsize ; qdown++){
+        if(m_EventData->GetChargeDownID(qdown)==ID){
+          rawQdown=m_EventData->GetChargeDown(qdown); 
+          if(rawQdown > m_Q_RAW_Threshold){
+            // Look for the associate time 
+            for(unsigned int tdown = 0 ; tdown < TDsize; tdown++){
+              if(m_EventData->GetTimeDownID(qdown)==ID) {
+                rawTdown=m_EventData->GetTimeDown(qdown);
+                break;
+              }
+            }// TDown
+          }//if raw threshold down
+
+          break;
+        } //if match ID 
+
+      }// Qdwown 
+
+      if(rawTdown>0){ // Tdown is found, means Qdown as well
+        // look for Tup  
+        for(unsigned int tup = 0 ; tup < TUsize ; tup++){
+          if(m_EventData->GetTimeUpID(tup)==ID){
+            rawTup = m_EventData->GetTimeUp(tup);
+            break;
+          }
+        }
+      }
+      // Got everything, do the math
+      if(rawTup>0){
+        calQup=aQu[ID]*(rawQup-bQu[ID]);
+        calQdown=aQd[ID]*(rawQdown-bQd[ID]);
+        // cal T 
+        calTup=aTu[ID]*rawTup+bTu[ID];
+        // slew correction
+        calTup -= slwTu[ID]/sqrt(rawQup-bQu[ID]);
+
+        // cal T
+        calTdown=aTd[ID]*rawTdown+bTd[ID];
+
+        // slew correction
+        calTdown -= slwTd[ID]/sqrt(rawQdown-bQd[ID]);
+
+        // average value of Up and Down
+        calQ=sqrt(calQup*calQdown); 
+        if(calQ>threshold){
+          calT= (calTdown+calTup)*0.5+avgT0[ID]+Cal->GetPedestal("NEBULA_T_ID"+NPL::itoa(ID)); 
+          Y=(calTdown-calTup)*DTa[ID]+DTb[ID]+Cal->GetPedestal("NEBULA_Y_ID"+NPL::itoa(ID));
+
+          DetectorNumber.push_back(ID);
+          Charge.push_back(calQ);
+          TOF.push_back(calT);
+          PosY.push_back(Y+PositionY[ID]);
+          PosX.push_back(PositionX[ID]);
+          PosZ.push_back(PositionZ[ID]);
+
+          if(ID<120)
+            IsVeto.push_back(0);
+          else
+            IsVeto.push_back(1);
+
+        }
       }
-    }
-  }
 
-  // Time 
-  mysize = m_EventData->GetMultTime();
-  for (UShort_t i = 0; i < mysize; ++i) {
-    Double_t Time= Cal->ApplyCalibration("Nebula/TIME"+NPL::itoa(m_EventData->GetT_DetectorNbr(i)),m_EventData->Get_Time(i));
-    m_PreTreatedData->SetTime(m_EventData->GetT_DetectorNbr(i), Time);
-  }
-  */
+    }// if raw threshold up
+  } // Qup
 }
 
-
-
 ///////////////////////////////////////////////////////////////////////////
-void TNebulaPhysics::ReadAnalysisConfig() {
-  bool ReadingStatus = false;
-
-  // path to file
-  string FileName = "./configs/ConfigNebula.dat";
-
-  // open analysis config file
-  ifstream AnalysisConfigFile;
-  AnalysisConfigFile.open(FileName.c_str());
+void TNebulaPhysics::PreTreat() {
 
-  if (!AnalysisConfigFile.is_open()) {
-    cout << " No ConfigNebula.dat found: Default parameter loaded for Analayis " << FileName << endl;
-    return;
-  }
-  cout << " Loading user parameter for Analysis from ConfigNebula.dat " << endl;
-
-  // Save it in a TAsciiFile
-  TAsciiFile* asciiConfig = RootOutput::getInstance()->GetAsciiFileAnalysisConfig();
-  asciiConfig->AppendLine("%%% ConfigNebula.dat %%%");
-  asciiConfig->Append(FileName.c_str());
-  asciiConfig->AppendLine("");
-  // read analysis config file
-  string LineBuffer,DataBuffer,whatToDo;
-  while (!AnalysisConfigFile.eof()) {
-    // Pick-up next line
-    getline(AnalysisConfigFile, LineBuffer);
-
-    // search for "header"
-    string name = "ConfigNebula";
-    if (LineBuffer.compare(0, name.length(), name) == 0) 
-      ReadingStatus = true;
-
-    // loop on tokens and data
-    while (ReadingStatus ) {
-      whatToDo="";
-      AnalysisConfigFile >> whatToDo;
-
-      // Search for comment symbol (%)
-      if (whatToDo.compare(0, 1, "%") == 0) {
-        AnalysisConfigFile.ignore(numeric_limits<streamsize>::max(), '\n' );
-      }
+}
 
-      else if (whatToDo=="E_RAW_THRESHOLD") {
-        AnalysisConfigFile >> DataBuffer;
-        m_E_RAW_Threshold = atof(DataBuffer.c_str());
-        cout << whatToDo << " " << m_E_RAW_Threshold << endl;
-      }
 
-      else if (whatToDo=="E_THRESHOLD") {
-        AnalysisConfigFile >> DataBuffer;
-        m_E_Threshold = atof(DataBuffer.c_str());
-        cout << whatToDo << " " << m_E_Threshold << endl;
-      }
 
-      else {
-        ReadingStatus = false;
-      }
-    }
-  }
+///////////////////////////////////////////////////////////////////////////
+void TNebulaPhysics::ReadAnalysisConfig() {
 }
 
 
@@ -185,8 +212,12 @@ void TNebulaPhysics::ReadAnalysisConfig() {
 ///////////////////////////////////////////////////////////////////////////
 void TNebulaPhysics::Clear() {
   DetectorNumber.clear();
-  Energy.clear();
-  Time.clear();
+  Charge.clear();
+  TOF.clear();
+  PosY.clear();
+  PosX.clear();
+  PosZ.clear();
+  IsVeto.clear();
 }
 
 
@@ -227,7 +258,6 @@ void TNebulaPhysics::InitSpectra() {
 ///////////////////////////////////////////////////////////////////////////
 void TNebulaPhysics::FillSpectra() {
   m_Spectra -> FillRawSpectra(m_EventData);
-  m_Spectra -> FillPreTreatedSpectra(m_PreTreatedData);
   m_Spectra -> FillPhysicsSpectra(m_EventPhysics);
 }
 
@@ -268,8 +298,8 @@ void TNebulaPhysics::WriteSpectra() {
 void TNebulaPhysics::AddParameterToCalibrationManager() {
   CalibrationManager* Cal = CalibrationManager::getInstance();
   for (int i = 0; i < m_NumberOfBars; ++i) {
-    Cal->AddParameter("NEBULA_ID"+ NPL::itoa(i+1)+"_T");
-    Cal->AddParameter("NEBULA_ID"+ NPL::itoa(i+1)+"_Y");
+    Cal->AddParameter("NEBULA_T_ID"+ NPL::itoa(i+1));
+    Cal->AddParameter("NEBULA_Y_ID"+ NPL::itoa(i+1));
   }
 }
 
@@ -313,14 +343,14 @@ NPL::VDetector* TNebulaPhysics::Construct() {
 //            Registering the construct method to the factory                 //
 ////////////////////////////////////////////////////////////////////////////////
 extern "C"{
-class proxy_Nebula{
-  public:
-    proxy_Nebula(){
-      NPL::DetectorFactory::getInstance()->AddToken("NEBULA","Nebula");
-      NPL::DetectorFactory::getInstance()->AddDetector("NEBULA",TNebulaPhysics::Construct);
-    }
-};
+  class proxy_Nebula{
+    public:
+      proxy_Nebula(){
+        NPL::DetectorFactory::getInstance()->AddToken("NEBULA","Nebula");
+        NPL::DetectorFactory::getInstance()->AddDetector("NEBULA",TNebulaPhysics::Construct);
+      }
+  };
 
-proxy_Nebula p_Nebula;
+  proxy_Nebula p_Nebula;
 }
 
diff --git a/NPLib/Detectors/Nebula/TNebulaPhysics.h b/NPLib/Detectors/Nebula/TNebulaPhysics.h
index 224e49c224ec190b020b65782cc16fd082799cb4..be14226c6b1966f53693ec6f5b50e496267f6e9a 100644
--- a/NPLib/Detectors/Nebula/TNebulaPhysics.h
+++ b/NPLib/Detectors/Nebula/TNebulaPhysics.h
@@ -52,26 +52,64 @@ class TNebulaPhysics : public TObject, public NPL::VDetector {
     ~TNebulaPhysics() {};
 
 
-  //////////////////////////////////////////////////////////////
-  // Inherited from TObject and overriden to avoid warnings
+    //////////////////////////////////////////////////////////////
+    // Inherited from TObject and overriden to avoid warnings
   public: 
     void Clear();   
     void Clear(const Option_t*) {};
 
 
-  //////////////////////////////////////////////////////////////
-  // data obtained after BuildPhysicalEvent() and stored in
-  // output ROOT file
+    //////////////////////////////////////////////////////////////
+    // data obtained after BuildPhysicalEvent() and stored in
+    // output ROOT file
   public:
     vector<int>      DetectorNumber;
-    vector<double>   Energy;
-    vector<double>   Time;
+    vector<double>   Charge;
+    vector<double>   TOF;
+    vector<double>   PosY;
+    vector<double>   PosX;
+    vector<double>   PosZ;
+    vector<bool>     IsVeto;
 
-  /// A usefull method to bundle all operation to add a detector
-  void ReadXML(NPL::XmlParser); 
-  
-  //////////////////////////////////////////////////////////////
-  // methods inherited from the VDetector ABC class
+  public:
+    TVector3 GetPos(const unsigned int& i) const{
+      return TVector3(PosX[i],PosY[i],PosZ[i]);
+    }
+
+    // Return true if one veto fired
+    bool HasVeto(){
+      unsigned int size = IsVeto.size();
+      for(unsigned int i = 0 ; i < size ; i++){
+        if(IsVeto[i])
+          return true;
+      }
+      return false;
+    };
+
+    /////////// Get index of fastest neutron
+    int GetFirstHit(){
+      unsigned int size = TOF.size();
+      unsigned int index=0;
+
+      if(!size)
+        return -1;
+
+      double tof = TOF[0];
+      for(unsigned int i = 1 ; i < size ; i++){
+        if(tof<TOF[i]){
+          tof=TOF[i];
+          index=i;
+        }
+      }
+      return index;
+    };
+
+  public:
+    /// A usefull method to bundle all operation to add a detector
+    void ReadXML(NPL::XmlParser); 
+
+    //////////////////////////////////////////////////////////////
+    // methods inherited from the VDetector ABC class
   public:
     // read stream from ConfigFile to pick-up detector parameters
     void ReadConfiguration(NPL::InputParser);
@@ -126,40 +164,36 @@ class TNebulaPhysics : public TObject, public NPL::VDetector {
     void WriteSpectra();
 
 
-  //////////////////////////////////////////////////////////////
-  // specific methods to Nebula array
+    //////////////////////////////////////////////////////////////
+    // specific methods to Nebula array
   public:
     // remove bad channels, calibrate the data and apply thresholds
     void PreTreat();
 
-    // clear the pre-treated object
-    void ClearPreTreatedData()   {m_PreTreatedData->Clear();}
-
     // read the user configuration file. If no file is found, load standard one
     void ReadAnalysisConfig();
 
     // give and external TNebulaData object to TNebulaPhysics. 
     // needed for online analysis for example
     void SetRawDataPointer(TNebulaData* rawDataPointer) {m_EventData = rawDataPointer;}
-    
-  // objects are not written in the TTree
+
+    // objects are not written in the TTree
   private:
     TNebulaData*         m_EventData;        //!
-    TNebulaData*         m_PreTreatedData;   //!
     TNebulaPhysics*      m_EventPhysics;     //!
 
-  // getters for raw and pre-treated data object
+    // getters for raw and pre-treated data object
   public:
     TNebulaData* GetRawData()        const {return m_EventData;}
-    TNebulaData* GetPreTreatedData() const {return m_PreTreatedData;}
 
-  // parameters used in the analysis
+    // parameters used in the analysis
   private:
     // thresholds
-    double m_E_RAW_Threshold; //!
-    double m_E_Threshold;     //!
+    double m_Q_RAW_Threshold; //!
+    double m_Q_Threshold;     //!
+    double m_V_Threshold;     //!
 
-  // number of detectors
+    // number of detectors
   private:
     int m_NumberOfBars;  //!
 
@@ -168,16 +202,43 @@ class TNebulaPhysics : public TObject, public NPL::VDetector {
     std::map<unsigned int, bool> m_invertX;//!
     std::map<unsigned int, bool> m_invertY;//!
 
+  private: // xml calibration
+    // position
+    std::map<unsigned int , double > PositionX;//!
+    std::map<unsigned int , double > PositionY;//!
+    std::map<unsigned int , double > PositionZ;//!
+
+    // linear cal
+    std::map<unsigned int , double > aQu;//!
+    std::map<unsigned int , double > bQu;//!
+    std::map<unsigned int , double > aQd;//!
+    std::map<unsigned int , double > bQd;//!
+    std::map<unsigned int , double > aTu;//!
+    std::map<unsigned int , double > bTu;//!
+    std::map<unsigned int , double > aTd;//!
+    std::map<unsigned int , double > bTd;//!
+
+    // T average offset
+    std::map<unsigned int , double > avgT0;//!
+
+    // slew correction T= tcal +slwT/sqrt(Qcal)
+    std::map<unsigned int , double > slwTu;//!
+    std::map<unsigned int , double > slwTd;//!
+
+    // DT position cal
+    std::map<unsigned int , double > DTa;//!
+    std::map<unsigned int , double > DTb;//!
+
 
-  // spectra class
+    // spectra class
   private:
     TNebulaSpectra* m_Spectra; // !
 
-  // spectra getter
+    // spectra getter
   public:
     map<string, TH1*>   GetSpectra(); 
 
-  // Static constructor to be passed to the Detector Factory
+    // Static constructor to be passed to the Detector Factory
   public:
     static NPL::VDetector* Construct();
 
diff --git a/NPLib/Physics/NPParticle.cxx b/NPLib/Physics/NPParticle.cxx
index 66125ce2af4ed76602c43e7fc5e59900880fa170..72849bcaec006191de63d34a58c85266eea5fea0 100644
--- a/NPLib/Physics/NPParticle.cxx
+++ b/NPLib/Physics/NPParticle.cxx
@@ -510,7 +510,7 @@ void Particle::BetaToGamma(){
 
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 
 void Particle::BetaToVelocity(){
-  fVelocity = (c_light*1e6)*fBeta*1e-7;
+  fVelocity = c_light*fBeta;
 }
 
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 
diff --git a/NPLib/Physics/NPParticle.h b/NPLib/Physics/NPParticle.h
index 32594de241c9c0b31211136648c86fb71ed99855..e1aa9dd41a15db086eda13737fa02e8ff2759cd4 100644
--- a/NPLib/Physics/NPParticle.h
+++ b/NPLib/Physics/NPParticle.h
@@ -139,6 +139,7 @@ namespace NPL {
       EnergyToBeta();
       BetaToGamma();
       BetaToVelocity();}
+
     void SetExcitationEnergy(double Ex) {fExcitationEnergy=Ex;}
     void				SetBeta(double beta)					{fBeta = beta; BetaToGamma(); BetaToEnergy(); EnergyToTof(); EnergyToBrho();BetaToVelocity();}
       double GetEnergyCM(double EnergyLab, double ThetaLab, double PhiLab, double relativisticboost);
diff --git a/NPLib/TrackReconstruction/NPDCReconstructionMT.cxx b/NPLib/TrackReconstruction/NPDCReconstructionMT.cxx
index 73c443aeca852ef06f3b37fc1e67a563fe7a9eb4..390a6272281c105481431a3cdfa4a2620d90b9f0 100644
--- a/NPLib/TrackReconstruction/NPDCReconstructionMT.cxx
+++ b/NPLib/TrackReconstruction/NPDCReconstructionMT.cxx
@@ -229,7 +229,7 @@ void DCReconstructionMT::StartThread(unsigned int id){
   mini->SetPrintLevel(0);
 
   // Let the main thread start
-  std::this_thread::sleep_for(std::chrono::milliseconds(100));
+  std::this_thread::sleep_for(std::chrono::milliseconds(250));
   while(true){
     // Do the job if possible
     if(m_Ready[id]){
diff --git a/Projects/S034/Analysis.cxx b/Projects/S034/Analysis.cxx
index 544f76824ebf0db7a597f97ded69ad4c056fbacf..c0f4265c85afce31f86e495598866641eeba3483 100644
--- a/Projects/S034/Analysis.cxx
+++ b/Projects/S034/Analysis.cxx
@@ -28,6 +28,8 @@ using namespace std;
 #include"NPDetectorManager.h"
 #include"RootInput.h"
 #include"RootOutput.h"
+#include"NPParticle.h"
+
 ////////////////////////////////////////////////////////////////////////////////
 Analysis::Analysis(){
 }
@@ -38,11 +40,13 @@ Analysis::~Analysis(){
 ////////////////////////////////////////////////////////////////////////////////
 void Analysis::Init(){
   Minos= (TMinosPhysics*) m_DetectorManager->GetDetector("Minos");
-  //Nebula = (TNebulaPhysics*) m_DetectorManager->GetDetector("NEBULA");
+  Nebula = (TNebulaPhysics*) m_DetectorManager->GetDetector("NEBULA");
   BDC = (TSamuraiBDCPhysics*) m_DetectorManager->GetDetector("SAMURAIBDC");
   FDC0 = (TSamuraiFDC0Physics*) m_DetectorManager->GetDetector("SAMURAIFDC0");
   FDC2 = (TSamuraiFDC2Physics*) m_DetectorManager->GetDetector("SAMURAIFDC2");
   Hodo = (TSamuraiHodoscopePhysics*) m_DetectorManager->GetDetector("SAMURAIHOD");
+
+  FragmentTarget = NPL::EnergyLoss("He6_LH2.G4table","G4Table",1000 );
   m_field.LoadMap(30*deg,"field_map/180702-2,40T-3000.table.bin",10);
   m_field.SetFDC2Angle((59.930-90.0)*deg);
   m_field.SetFDC2R(FDC2->GetOffset().Z());
@@ -55,47 +59,84 @@ void Analysis::Init(){
 ////////////////////////////////////////////////////////////////////////////////
 void Analysis::TreatEvent(){
   Clear();
-  //cout << Trigger << " " ; 
+  
+  static NPL::Particle He6("6He");
+  static NPL::Particle He8("8He");
+  static NPL::Particle n("neutron");
+  static double mhe6 = He6.Mass();
+  static double mn   = n.Mass();
+  static double sumM= mhe6+mn;
+  static TLorentzVector LVHe6;
+  static TLorentzVector LVn;
+  He8.SetBeta(0.5152);
   Trigger=Trigger&0x00ff;
-  //cout << Trigger << endl;
-  // Compute Brho 
-  //
   if( FDC2->PosX>-1500 && FDC2->PosX<1000 
       && FDC2->PosY>-500 && FDC2->PosY<500 
       && FDC0->PosX>-80 && FDC0->PosX<80 
       && FDC0->PosY>-80 && FDC0->PosY<80 // both FDC ok
-      && (Minos->Tracks_P0.size()>1)) {  // p,pn or p,2p
+      && Minos->Tracks_P0.size()>0     ) {  // p,pn or p,2p
     // Compute ThetaX and PhiY using Minos vertex and FDC0 X
     // Check if both BDC are reconstructed
     TVector3 BDC1=BDC->GetPos(1);
     TVector3 BDC2=BDC->GetPos(2);
 
     if( BDC1.Z()!=-10000 && BDC2.Z()!=-10000){
-      TVector3 Vertex;
+      TVector3 Vertex,delta;
       TVector3 P1 = Minos->Tracks_P0[0]+Minos->Tracks_Dir[0];
-      TVector3 delta;
+
+      // Vertex in Samurai Reference frame 
       MinimumDistanceTwoLines(BDC1,BDC2, 
           Minos->Tracks_P0[0], P1,
           Vertex, delta) ;
+     
       TVector3 FDC0_Dir= FDC0->GetPos()-Vertex;
       FDC0_Dir=FDC0_Dir.Unit();
       TVector3 BDCDir=BDC2-BDC1;
       BDCDir=BDCDir.Unit();
       BDCDir*=(Vertex.Z()-BDC2.Z())/BDCDir.Z();
-      //cout << "BDC" << endl;
       BDCX=(BDC2+BDCDir).X();
       BDCY=(BDC2+BDCDir).Y();
-      // Z relative to Minos entrance
-      Z=Vertex.Z()+4657.39;
-      //double brho_param[6]={FDC0->PosX/*+1.77*/, FDC0->PosY, 0, 0, FDC2->PosX/*-252.55*/, FDC2->ThetaX};
 
-      if(FDC0_Dir.Z()>0.6){
+      // XYZ relative to Minos entrance
+      X=Vertex.X();
+      Y=Vertex.Y();
+      Z=Vertex.Z()+4650;
+
+      // Fragment analysis 
+      if(FDC0_Dir.Z()>0.6 && Z>0 && Z<150 && sqrt(X*X+Y*Y)<15){
         double FDC0_ThetaX = atan((FDC0->PosX-Vertex.X())/(1254.39-Z));
         double FDC0_PhiY   = atan((FDC0->PosY-Vertex.Y())/(1254.39-Z));
-        double brho_param[6]={FDC0->PosX, FDC0->PosY, tan(FDC0_ThetaX), tan(FDC0_PhiY), FDC2->PosX+252.416, FDC2->ThetaX};
-        BrhoP=r_fit(brho_param);
         Brho=m_field.FindBrho(FDC0->GetPos(),FDC0_Dir,FDC2->GetPos(),TVector3(0,0,1));
-        //    cout << Brho-BrhoP << endl;
+        He6.SetBrho(Brho);
+        double Energy = He6.GetEnergy();
+        if(Energy){
+          Energy=FragmentTarget.EvaluateInitialEnergy(Energy,150-Z,FDC0_Dir.Angle(TVector3(0,0,1)));
+          He6.SetKineticEnergy(Energy);
+          Beta_f=He6.GetBeta();
+          LVHe6.SetVectM(TVector3(0,0,0),mhe6);
+          LVHe6.Boost(Beta_f*FDC0_Dir.Unit());
+        }
+        
+      }
+      
+      // Neutron analysis
+      if(Nebula->Charge.size()>0 && !Nebula->HasVeto()){
+        // Index of first neutron hit
+        unsigned int first = Nebula->GetFirstHit();
+        TVector3 Pfirst = (Nebula->GetPos(first)-Vertex);
+        double L = Pfirst.Mag();
+        double TSBT= (Vertex.Z()+7377.56)/He8.GetVelocity();
+        double TOF = Nebula->TOF[first]-TSBT;
+        Beta_n = (L/TOF)/NPUNITS::c_light;
+        LVn.SetVectM(TVector3(0,0,0),mn);
+        LVn.Boost(Beta_n*Pfirst.Unit());
+      }
+      
+      // a full event
+      if(Beta_n&&Beta_f){
+        He6.SetEnergyImpulsion(LVHe6);
+        n.SetEnergyImpulsion(LVn);
+        Erel = (LVHe6+LVn).Mag()-sumM;
       }
       // Calib//////////////////////////////////////////////////////////////////
  /*     static int count=0;
@@ -132,1610 +173,45 @@ void Analysis::Clear(){
   Brho=-1000;
   BDCX=-1000;
   BDCY=-1000;
-  BrhoP=-1000;
   Beta_f=-1000;
-  Z=-1000;
+  Beta_n=-1000;
+  X=Y=Z=-1000;
 }
 
 ////////////////////////////////////////////////////////////////////////////////
 void Analysis::InitOutputBranch() {
   RootOutput::getInstance()->GetTree()->Branch("Brho",&Brho,"Brho/D");
-  RootOutput::getInstance()->GetTree()->Branch("BrhoP",&BrhoP,"BrhoP/D");
   RootOutput::getInstance()->GetTree()->Branch("BDCX",&BDCX,"BDCX/D");
   RootOutput::getInstance()->GetTree()->Branch("BDCY",&BDCY,"BDCY/D");
+  RootOutput::getInstance()->GetTree()->Branch("X",&X,"X/D");
+  RootOutput::getInstance()->GetTree()->Branch("Y",&Y,"Y/D");
   RootOutput::getInstance()->GetTree()->Branch("Z",&Z,"Z/D");
+  RootOutput::getInstance()->GetTree()->Branch("Erel",&Erel,"Erel/D");
   RootOutput::getInstance()->GetTree()->Branch("Beta_f",&Beta_f,"Beta_f/D");
+  RootOutput::getInstance()->GetTree()->Branch("Beta_n",&Beta_n,"Beta_n/D");
   RootOutput::getInstance()->GetTree()->Branch("Trigger",&Trigger,"Trigger/I");
 } 
 ////////////////////////////////////////////////////////////////////////////////
 void Analysis::InitInputBranch(){
   RootInput::getInstance()->GetChain()->SetBranchAddress("Trigger",&Trigger);
 }
-
-// -*- mode: c++ -*-
-// 
-// File simtree_energy.C generated by TMultiDimFit::MakeRealCode
-// on Fri May 28 00:25:23 2021
-// ROOT version 5.32/00
-//
-// This file contains the function 
-//
-//    double  e_fit(double *x); 
-//
-// For evaluating the parameterization obtained
-// from TMultiDimFit and the point x
-// 
-// See TMultiDimFit class documentation for more information 
-// 
-//
-// Static data variables
-//
-static int    e_gNVariables    = 6;
-static int    e_gNCoefficients = 76;
-static double e_gDMean         = 69.5201;
-// Assignment to mean vector.
-static double e_gXMean[] = {
-  0.00503955, -0.0144156, -0.101042, -0.00374735, 0.252452, -0.512868 };
-
-// Assignment to minimum vector.
-static double e_gXMin[] = {
-  -64.3965, -52.0748, -0.0335722, -0.0401226, -1065.23, -299.807 };
-
-// Assignment to maximum vector.
-static double e_gXMax[] = {
-  52.3459, 55.8258, 0.0432727, 0.04243, 1147.99, 931.885 };
-
-// Assignment to coefficients vector.
-static double e_gCoefficient[] = {
-  -871.253,
-  -1919.06,
-  -788.479,
-  2645.64,
-  443.115,
-  984.07,
-  243.005,
-  3552.25,
-  -157.556,
-  2188.2,
-  618.533,
-  1176.32,
-  1303.87,
-  -1056.05,
-  -772.684,
-  6.86646,
-  1343.46,
-  3390.39,
-  871.459,
-  7388.91,
-  -599.788,
-  -4023.35,
-  -299.575,
-  1670.13,
-  0.620256,
-  10.0888,
-  -0.557305,
-  -1.12306,
-  9.51285,
-  0.728805,
-  -0.620644,
-  -4.05522,
-  0.633792,
-  -11.1877,
-  -4.56183,
-  -5.71034,
-  -4.45376,
-  1.38591,
-  1.18557,
-  -4.92592,
-  -5.16675,
-  -4.66713,
-  -1506.82,
-  -306.44,
-  4269.09,
-  2292.62,
-  19.7398,
-  2612.26,
-  1682.7,
-  -1163.29,
-  -597.13,
-  -7838.7,
-  14404.4,
-  6607.08,
-  -532.085,
-  -317.801,
-  681.575,
-  -546.576,
-  -775.207,
-  -3803.59,
-  -1160.61,
-  0.0848953,
-  4.22188,
-  -1.96972,
-  0.794997,
-  -2260.46,
-  -1.21301,
-  -7409.46,
-  269.548,
-  -0.00925447,
-  2.77788,
-  1.95742,
-  0.453985,
-  -0.188488,
-  -0.675742,
-  -0.181284
-};
-
-// Assignment to error coefficients vector.
-static double e_gCoefficientRMS[] = {
-  66.6599,
-  130.954,
-  45.1508,
-  146.993,
-  49.7567,
-  42.2773,
-  8.59505,
-  193.537,
-  31.5244,
-  92.7045,
-  67.4399,
-  97.1825,
-  63.1346,
-  63.8809,
-  37.0527,
-  0.0875267,
-  99.2757,
-  205.578,
-  127.446,
-  361.809,
-  107.175,
-  235.193,
-  35.8634,
-  105.166,
-  2.78933,
-  2.0882,
-  3.14393,
-  3.06267,
-  8.25249,
-  0.616472,
-  1.47943,
-  2.11076,
-  0.930677,
-  7.39957,
-  2.26951,
-  2.50688,
-  0.95993,
-  2.42862,
-  0.393511,
-  0.458493,
-  2.0148,
-  2.57894,
-  72.1773,
-  61.4191,
-  180.529,
-  189.341,
-  7.78812,
-  193.422,
-  248.274,
-  208.754,
-  69.8817,
-  458.225,
-  704.917,
-  400.536,
-  30.129,
-  31.4199,
-  27.3927,
-  40.2961,
-  54.0649,
-  207.698,
-  79.9551,
-  9.4176,
-  1.89081,
-  0.422792,
-  1.61055,
-  155.773,
-  1.45066,
-  404.667,
-  36.5963,
-  0.729026,
-  2.55087,
-  0.355723,
-  1.83536,
-  0.547759,
-  2.82987,
-  2.05714
-};
-
-// Assignment to powers vector.
-// The powers are stored row-wise, that is
-//  p_ij = e_gPower[i * NVariables + j];
-static int    e_gPower[] = {
-  1,  1,  1,  1,  1,  1,
-  1,  1,  1,  1,  1,  2,
-  2,  1,  1,  1,  1,  1,
-  1,  1,  1,  1,  2,  1,
-  1,  1,  1,  2,  1,  1,
-  1,  1,  2,  1,  1,  1,
-  1,  1,  1,  1,  1,  3,
-  1,  1,  1,  1,  2,  2,
-  1,  3,  1,  1,  1,  1,
-  2,  1,  2,  1,  1,  1,
-  1,  1,  1,  2,  1,  2,
-  1,  2,  2,  1,  1,  1,
-  1,  1,  2,  1,  1,  2,
-  2,  1,  1,  1,  1,  2,
-  3,  1,  1,  1,  1,  1,
-  1,  1,  1,  1,  3,  1,
-  2,  1,  1,  2,  1,  1,
-  1,  1,  1,  2,  2,  1,
-  1,  2,  1,  2,  1,  1,
-  1,  1,  2,  1,  2,  1,
-  1,  1,  2,  2,  1,  1,
-  2,  1,  1,  1,  2,  1,
-  1,  1,  1,  3,  1,  1,
-  1,  1,  1,  1,  2,  3,
-  2,  1,  1,  3,  1,  1,
-  1,  2,  1,  3,  1,  1,
-  1,  1,  2,  3,  1,  1,
-  1,  1,  3,  2,  1,  1,
-  1,  2,  2,  2,  1,  1,
-  1,  3,  1,  1,  2,  1,
-  1,  3,  1,  2,  1,  1,
-  1,  2,  2,  1,  2,  1,
-  1,  1,  3,  1,  2,  1,
-  2,  2,  1,  2,  1,  1,
-  2,  1,  1,  2,  2,  1,
-  1,  2,  1,  2,  2,  1,
-  1,  1,  1,  3,  2,  1,
-  3,  1,  1,  2,  1,  1,
-  2,  1,  1,  1,  3,  1,
-  1,  1,  2,  1,  3,  1,
-  3,  2,  1,  1,  1,  1,
-  1,  2,  3,  1,  1,  1,
-  3,  1,  1,  1,  1,  2,
-  1,  3,  1,  1,  1,  2,
-  2,  1,  2,  1,  1,  2,
-  1,  2,  2,  1,  1,  2,
-  2,  2,  2,  1,  1,  1,
-  2,  1,  1,  2,  1,  2,
-  1,  2,  1,  2,  1,  2,
-  1,  1,  2,  2,  1,  2,
-  1,  1,  1,  3,  1,  2,
-  2,  1,  1,  1,  2,  2,
-  1,  1,  2,  1,  2,  2,
-  1,  1,  1,  2,  2,  2,
-  2,  1,  1,  1,  1,  3,
-  1,  2,  1,  1,  1,  3,
-  1,  1,  2,  1,  1,  3,
-  1,  2,  1,  1,  1,  1,
-  1,  2,  1,  1,  1,  2,
-  1,  2,  1,  1,  2,  1,
-  2,  2,  1,  1,  1,  1,
-  2,  1,  2,  2,  1,  1,
-  2,  2,  1,  1,  2,  1,
-  1,  1,  1,  2,  3,  1,
-  2,  1,  3,  1,  1,  1,
-  2,  2,  1,  1,  1,  2,
-  3,  1,  2,  1,  1,  1,
-  1,  2,  1,  1,  2,  2,
-  1,  1,  1,  2,  1,  3,
-  3,  1,  1,  1,  2,  1,
-  1,  1,  2,  2,  2,  1,
-  1,  2,  1,  1,  3,  1,
-  2,  3,  1,  1,  1,  1,
-  1,  1,  3,  1,  1,  1,
-  2,  1,  2,  1,  2,  1,
-  1,  3,  2,  1,  1,  1
-};
-
-// 
-// The function   double e_fit(double *x)
-// 
-double e_fit(double *x) {
-  double returnValue = e_gDMean;
-  int    i = 0, j = 0, k = 0;
-  for (i = 0; i < e_gNCoefficients ; i++) {
-    // Evaluate the ith term in the expansion
-    double term = e_gCoefficient[i];
-    for (j = 0; j < e_gNVariables; j++) {
-      // Evaluate the polynomial in the jth variable.
-      int power = e_gPower[e_gNVariables * i + j]; 
-      double p1 = 1, p2 = 0, p3 = 0, r = 0;
-      double v =  1 + 2. / (e_gXMax[j] - e_gXMin[j]) * (x[j] - e_gXMax[j]);
-      // what is the power to use!
-      switch(power) {
-        case 1: r = 1; break; 
-        case 2: r = v; break; 
-        default: 
-                p2 = v; 
-                for (k = 3; k <= power; k++) { 
-                  p3 = p2 * v;
-                  p3 = 2 * v * p2 - p1; 
-                  p1 = p2; p2 = p3; 
-                }
-                r = p3;
-      }
-      // multiply this term by the poly in the jth var
-      term *= r; 
-    }
-    // Add this term to the final result
-    returnValue += term;
-  }
-  return returnValue;
-}
-
-// EOF for simtree_energy.C
-// -*- mode: c++ -*-
-// 
-// File simtree_momentum.C generated by TMultiDimFit::MakeRealCode
-// on Fri May 28 00:25:23 2021
-// ROOT version 5.32/00
-//
-// This file contains the function 
-//
-//    double  p_fit(double *x); 
-//
-// For evaluating the parameterization obtained
-// from TMultiDimFit and the point x
-// 
-// See TMultiDimFit class documentation for more information 
-// 
-//
-// Static data variables
-//
-static int    p_gNVariables    = 6;
-static int    p_gNCoefficients = 76;
-static double p_gDMean         = 361.833;
-// Assignment to mean vector.
-static double p_gXMean[] = {
-  0.00503955, -0.0144156, -0.101042, -0.00374735, 0.252452, -0.512868 };
-
-// Assignment to minimum vector.
-static double p_gXMin[] = {
-  -64.3965, -52.0748, -0.0335722, -0.0401226, -1065.23, -299.807 };
-
-// Assignment to maximum vector.
-static double p_gXMax[] = {
-  52.3459, 55.8258, 0.0432727, 0.04243, 1147.99, 931.885 };
-
-// Assignment to coefficients vector.
-static double p_gCoefficient[] = {
-  -1679.5,
-  -3768.96,
-  -1773.35,
-  5927.56,
-  -1052.24,
-  786.483,
-  2247.73,
-  553.038,
-  7920.2,
-  -326.613,
-  1147.19,
-  3010.2,
-  -1520.59,
-  -2398.64,
-  12.0318,
-  2542.51,
-  7843.96,
-  1874.12,
-  17655.7,
-  -902.702,
-  -9667.32,
-  -2259.82,
-  -673.898,
-  3744.62,
-  23.4495,
-  -8.41591,
-  23.1301,
-  23.3445,
-  -0.850973,
-  -21.3182,
-  -5.51343,
-  -13.2501,
-  0.0420299,
-  4.8444,
-  4.76409,
-  -11.8105,
-  -4406.58,
-  -637.169,
-  8498.88,
-  4695.17,
-  4938.56,
-  3623.63,
-  -1738.8,
-  -1348.34,
-  34420.8,
-  15289.1,
-  -1177.71,
-  -587.346,
-  1529.61,
-  450.733,
-  4356.32,
-  2406.82,
-  0.828932,
-  -1494.72,
-  -7.21335,
-  -4.61121,
-  2.36389,
-  -9.47311,
-  9.7888,
-  -3.50166,
-  -12.2852,
-  -2912.85,
-  -0.277797,
-  42.1753,
-  -1.26653,
-  -18825.4,
-  -8755.82,
-  3.69156,
-  8.16369,
-  -7.45288,
-  -0.273504,
-  -9.0052,
-  -17056.8,
-  -0.145767,
-  -1.13422,
-  0.865558
-};
-
-// Assignment to error coefficients vector.
-static double p_gCoefficientRMS[] = {
-  1.69031e-11,
-  3.29223e-11,
-  6.40666e-11,
-  2.89085e-11,
-  6.58543e-11,
-  8.69785e-11,
-  7.0898e-11,
-  3.57517e-11,
-  5.62969e-11,
-  1.89724e-11,
-  1.69397e-10,
-  1.38078e-10,
-  1.28278e-10,
-  1.24785e-10,
-  2.44521e-11,
-  3.47315e-10,
-  1.44807e-10,
-  2.73924e-10,
-  1.18391e-10,
-  3.51888e-10,
-  1.21525e-10,
-  2.55786e-10,
-  1.82279e-11,
-  6.11644e-11,
-  7.41771e-11,
-  9.78353e-11,
-  1.16067e-09,
-  1.31561e-09,
-  3.16376e-11,
-  1.03428e-09,
-  6.01156e-10,
-  3.13162e-11,
-  9.7543e-11,
-  9.29464e-11,
-  8.80057e-11,
-  1.16199e-10,
-  4.98203e-10,
-  3.6952e-11,
-  4.19168e-10,
-  5.29444e-10,
-  6.7647e-10,
-  5.33585e-10,
-  6.85262e-10,
-  3.55029e-11,
-  2.30548e-10,
-  2.81963e-10,
-  1.35502e-10,
-  1.39257e-10,
-  1.49982e-10,
-  1.83996e-10,
-  2.15232e-10,
-  2.71861e-10,
-  1.88434e-11,
-  1.92304e-11,
-  7.69361e-11,
-  1.02682e-10,
-  3.52953e-10,
-  6.53796e-10,
-  5.36745e-10,
-  1.27112e-10,
-  7.41853e-11,
-  3.7455e-11,
-  7.64436e-11,
-  9.64104e-10,
-  8.7203e-11,
-  2.36674e-10,
-  1.34232e-10,
-  6.85812e-11,
-  5.67547e-10,
-  4.46781e-10,
-  3.2384e-11,
-  7.34435e-11,
-  2.61417e-10,
-  3.21076e-11,
-  8.11735e-11,
-  7.22006e-11
-};
-
-// Assignment to powers vector.
-// The powers are stored row-wise, that is
-//  p_ij = p_gPower[i * NVariables + j];
-static int    p_gPower[] = {
-  1,  1,  1,  1,  1,  1,
-  1,  1,  1,  1,  1,  2,
-  2,  1,  1,  1,  1,  1,
-  1,  1,  1,  1,  2,  1,
-  1,  2,  1,  1,  1,  1,
-  1,  1,  1,  2,  1,  1,
-  1,  1,  2,  1,  1,  1,
-  1,  1,  1,  1,  1,  3,
-  1,  1,  1,  1,  2,  2,
-  1,  3,  1,  1,  1,  1,
-  1,  1,  1,  2,  1,  2,
-  1,  1,  2,  1,  1,  2,
-  1,  2,  1,  1,  1,  2,
-  2,  1,  1,  1,  1,  2,
-  1,  1,  1,  1,  3,  1,
-  2,  1,  1,  2,  1,  1,
-  1,  1,  1,  2,  2,  1,
-  1,  2,  1,  2,  1,  1,
-  1,  1,  2,  1,  2,  1,
-  1,  1,  2,  2,  1,  1,
-  2,  1,  1,  1,  2,  1,
-  2,  2,  1,  1,  1,  1,
-  1,  1,  1,  3,  1,  1,
-  1,  1,  1,  1,  2,  3,
-  1,  2,  1,  3,  1,  1,
-  1,  1,  3,  2,  1,  1,
-  1,  2,  2,  2,  1,  1,
-  2,  1,  2,  2,  1,  1,
-  1,  3,  1,  1,  2,  1,
-  2,  2,  1,  2,  1,  1,
-  1,  2,  1,  2,  2,  1,
-  1,  1,  1,  3,  2,  1,
-  3,  1,  1,  2,  1,  1,
-  2,  1,  1,  1,  3,  1,
-  1,  2,  1,  1,  3,  1,
-  1,  1,  2,  1,  3,  1,
-  2,  2,  1,  1,  1,  2,
-  1,  3,  1,  1,  1,  2,
-  2,  1,  2,  1,  1,  2,
-  1,  2,  2,  1,  1,  2,
-  2,  1,  1,  2,  1,  2,
-  1,  2,  1,  2,  1,  2,
-  1,  1,  2,  2,  1,  2,
-  1,  1,  1,  3,  1,  2,
-  1,  1,  2,  1,  2,  2,
-  1,  1,  1,  2,  2,  2,
-  2,  1,  1,  1,  1,  3,
-  1,  2,  1,  1,  1,  3,
-  1,  1,  2,  1,  1,  3,
-  1,  1,  1,  2,  1,  3,
-  2,  1,  2,  1,  1,  1,
-  1,  2,  2,  1,  1,  1,
-  1,  1,  3,  1,  1,  1,
-  3,  1,  1,  1,  1,  1,
-  1,  1,  2,  3,  1,  1,
-  1,  3,  1,  2,  1,  1,
-  2,  1,  2,  1,  2,  1,
-  2,  1,  1,  2,  2,  1,
-  1,  1,  2,  2,  2,  1,
-  1,  1,  1,  2,  3,  1,
-  3,  2,  1,  1,  1,  1,
-  3,  1,  1,  1,  1,  2,
-  2,  1,  3,  1,  1,  1,
-  2,  2,  2,  1,  1,  1,
-  3,  1,  2,  1,  1,  1,
-  2,  1,  1,  1,  2,  2,
-  1,  2,  1,  1,  2,  1,
-  2,  1,  1,  3,  1,  1,
-  2,  2,  1,  1,  2,  1,
-  1,  2,  2,  1,  2,  1,
-  1,  1,  3,  1,  2,  1,
-  1,  2,  3,  1,  1,  1,
-  1,  2,  1,  1,  2,  2,
-  3,  1,  1,  1,  2,  1,
-  1,  3,  2,  1,  1,  1,
-  2,  3,  1,  1,  1,  1
-};
-
-// 
-// The function   double p_fit(double *x)
-// 
-double p_fit(double *x) {
-  double returnValue = p_gDMean;
-  int    i = 0, j = 0, k = 0;
-  for (i = 0; i < p_gNCoefficients ; i++) {
-    // Evaluate the ith term in the expansion
-    double term = p_gCoefficient[i];
-    for (j = 0; j < p_gNVariables; j++) {
-      // Evaluate the polynomial in the jth variable.
-      int power = p_gPower[p_gNVariables * i + j]; 
-      double p1 = 1, p2 = 0, p3 = 0, r = 0;
-      double v =  1 + 2. / (p_gXMax[j] - p_gXMin[j]) * (x[j] - p_gXMax[j]);
-      // what is the power to use!
-      switch(power) {
-        case 1: r = 1; break; 
-        case 2: r = v; break; 
-        default: 
-                p2 = v; 
-                for (k = 3; k <= power; k++) { 
-                  p3 = p2 * v;
-                  p3 = 2 * v * p2 - p1; 
-                  p1 = p2; p2 = p3; 
-                }
-                r = p3;
-      }
-      // multiply this term by the poly in the jth var
-      term *= r; 
-    }
-    // Add this term to the final result
-    returnValue += term;
-  }
-  return returnValue;
-}
-
-// EOF for simtree_momentum.C
-// -*- mode: c++ -*-
-// 
-// File simtree_rigidity.C generated by TMultiDimFit::MakeRealCode
-// on Fri May 28 00:25:23 2021
-// ROOT version 5.32/00
-//
-// This file contains the function 
-//
-//    double  r_fit(double *x); 
-//
-// For evaluating the parameterization obtained
-// from TMultiDimFit and the point x
-// 
-// See TMultiDimFit class documentation for more information 
-// 
-//
-// Static data variables
-//
-static int    r_gNVariables    = 6;
-static int    r_gNCoefficients = 76;
-static double r_gDMean         = 4.82778;
-// Assignment to mean vector.
-static double r_gXMean[] = {
-  0.00503955, -0.0144156, -0.101042, -0.00374735, 0.252452, -0.512868 };
-
-// Assignment to minimum vector.
-static double r_gXMin[] = {
-  -64.3965, -52.0748, -0.0335722, -0.0401226, -1065.23, -299.807 };
-
-// Assignment to maximum vector.
-static double r_gXMax[] = {
-  52.3459, 55.8258, 0.0432727, 0.04243, 1147.99, 931.885 };
-
-// Assignment to coefficients vector.
-static double r_gCoefficient[] = {
-  -22.4088,
-  -50.2876,
-  -23.6611,
-  79.0889,
-  -14.0396,
-  10.4937,
-  29.9905,
-  7.37895,
-  105.676,
-  -4.35786,
-  15.3065,
-  40.1638,
-  -20.2886,
-  -32.0041,
-  0.160535,
-  33.9237,
-  104.659,
-  25.0055,
-  235.572,
-  -12.0444,
-  -128.987,
-  -30.1518,
-  -8.99153,
-  49.9628,
-  0.312876,
-  -0.11229,
-  0.308614,
-  0.311475,
-  -0.0113542,
-  -0.28444,
-  -0.0735633,
-  -0.17679,
-  0.000560787,
-  0.0646367,
-  0.0635652,
-  -0.157583,
-  -58.795,
-  -8.50147,
-  113.397,
-  62.6455,
-  65.8931,
-  48.3485,
-  -23.2001,
-  -17.9904,
-  459.262,
-  203.996,
-  -15.7137,
-  -7.8367,
-  20.4089,
-  6.01393,
-  58.1245,
-  32.1131,
-  0.0110601,
-  -19.9434,
-  -0.0962446,
-  -0.0615253,
-  0.0315403,
-  -0.126396,
-  0.130608,
-  -0.0467211,
-  -0.163916,
-  -38.8649,
-  -0.00370653,
-  0.562727,
-  -0.0168988,
-  -251.179,
-  -116.825,
-  0.0492549,
-  0.108925,
-  -0.0994405,
-  -0.00364924,
-  -0.120152,
-  -227.582,
-  -0.0019449,
-  -0.0151334,
-  0.0115488
-};
-
-// Assignment to error coefficients vector.
-static double r_gCoefficientRMS[] = {
-  1.69031e-11,
-  3.29223e-11,
-  6.40666e-11,
-  2.89085e-11,
-  6.58543e-11,
-  8.69785e-11,
-  7.0898e-11,
-  3.57517e-11,
-  5.62969e-11,
-  1.89724e-11,
-  1.69397e-10,
-  1.38078e-10,
-  1.28278e-10,
-  1.24785e-10,
-  2.44521e-11,
-  3.47315e-10,
-  1.44807e-10,
-  2.73924e-10,
-  1.18391e-10,
-  3.51888e-10,
-  1.21525e-10,
-  2.55786e-10,
-  1.82279e-11,
-  6.11644e-11,
-  7.41771e-11,
-  9.78353e-11,
-  1.16067e-09,
-  1.31561e-09,
-  3.16376e-11,
-  1.03428e-09,
-  6.01156e-10,
-  3.13162e-11,
-  9.7543e-11,
-  9.29464e-11,
-  8.80057e-11,
-  1.16199e-10,
-  4.98203e-10,
-  3.6952e-11,
-  4.19168e-10,
-  5.29444e-10,
-  6.7647e-10,
-  5.33585e-10,
-  6.85262e-10,
-  3.55029e-11,
-  2.30548e-10,
-  2.81963e-10,
-  1.35502e-10,
-  1.39257e-10,
-  1.49982e-10,
-  1.83996e-10,
-  2.15232e-10,
-  2.71861e-10,
-  1.88434e-11,
-  1.92304e-11,
-  7.69361e-11,
-  1.02682e-10,
-  3.52953e-10,
-  6.53796e-10,
-  5.36745e-10,
-  1.27112e-10,
-  7.41853e-11,
-  3.7455e-11,
-  7.64436e-11,
-  9.64104e-10,
-  8.7203e-11,
-  2.36674e-10,
-  1.34232e-10,
-  6.85812e-11,
-  5.67547e-10,
-  4.46781e-10,
-  3.2384e-11,
-  7.34435e-11,
-  2.61417e-10,
-  3.21076e-11,
-  8.11735e-11,
-  7.22006e-11
-};
-
-// Assignment to powers vector.
-// The powers are stored row-wise, that is
-//  p_ij = r_gPower[i * NVariables + j];
-static int    r_gPower[] = {
-  1,  1,  1,  1,  1,  1,
-  1,  1,  1,  1,  1,  2,
-  2,  1,  1,  1,  1,  1,
-  1,  1,  1,  1,  2,  1,
-  1,  2,  1,  1,  1,  1,
-  1,  1,  1,  2,  1,  1,
-  1,  1,  2,  1,  1,  1,
-  1,  1,  1,  1,  1,  3,
-  1,  1,  1,  1,  2,  2,
-  1,  3,  1,  1,  1,  1,
-  1,  1,  1,  2,  1,  2,
-  1,  1,  2,  1,  1,  2,
-  1,  2,  1,  1,  1,  2,
-  2,  1,  1,  1,  1,  2,
-  1,  1,  1,  1,  3,  1,
-  2,  1,  1,  2,  1,  1,
-  1,  1,  1,  2,  2,  1,
-  1,  2,  1,  2,  1,  1,
-  1,  1,  2,  1,  2,  1,
-  1,  1,  2,  2,  1,  1,
-  2,  1,  1,  1,  2,  1,
-  2,  2,  1,  1,  1,  1,
-  1,  1,  1,  3,  1,  1,
-  1,  1,  1,  1,  2,  3,
-  1,  2,  1,  3,  1,  1,
-  1,  1,  3,  2,  1,  1,
-  1,  2,  2,  2,  1,  1,
-  2,  1,  2,  2,  1,  1,
-  1,  3,  1,  1,  2,  1,
-  2,  2,  1,  2,  1,  1,
-  1,  2,  1,  2,  2,  1,
-  1,  1,  1,  3,  2,  1,
-  3,  1,  1,  2,  1,  1,
-  2,  1,  1,  1,  3,  1,
-  1,  2,  1,  1,  3,  1,
-  1,  1,  2,  1,  3,  1,
-  2,  2,  1,  1,  1,  2,
-  1,  3,  1,  1,  1,  2,
-  2,  1,  2,  1,  1,  2,
-  1,  2,  2,  1,  1,  2,
-  2,  1,  1,  2,  1,  2,
-  1,  2,  1,  2,  1,  2,
-  1,  1,  2,  2,  1,  2,
-  1,  1,  1,  3,  1,  2,
-  1,  1,  2,  1,  2,  2,
-  1,  1,  1,  2,  2,  2,
-  2,  1,  1,  1,  1,  3,
-  1,  2,  1,  1,  1,  3,
-  1,  1,  2,  1,  1,  3,
-  1,  1,  1,  2,  1,  3,
-  2,  1,  2,  1,  1,  1,
-  1,  2,  2,  1,  1,  1,
-  1,  1,  3,  1,  1,  1,
-  3,  1,  1,  1,  1,  1,
-  1,  1,  2,  3,  1,  1,
-  1,  3,  1,  2,  1,  1,
-  2,  1,  2,  1,  2,  1,
-  2,  1,  1,  2,  2,  1,
-  1,  1,  2,  2,  2,  1,
-  1,  1,  1,  2,  3,  1,
-  3,  2,  1,  1,  1,  1,
-  3,  1,  1,  1,  1,  2,
-  2,  1,  3,  1,  1,  1,
-  2,  2,  2,  1,  1,  1,
-  3,  1,  2,  1,  1,  1,
-  2,  1,  1,  1,  2,  2,
-  1,  2,  1,  1,  2,  1,
-  2,  1,  1,  3,  1,  1,
-  2,  2,  1,  1,  2,  1,
-  1,  2,  2,  1,  2,  1,
-  1,  1,  3,  1,  2,  1,
-  1,  2,  3,  1,  1,  1,
-  1,  2,  1,  1,  2,  2,
-  3,  1,  1,  1,  2,  1,
-  1,  3,  2,  1,  1,  1,
-  2,  3,  1,  1,  1,  1
-};
-
-// 
-// The function   double r_fit(double *x)
-// 
-double r_fit(double *x) {
-  double returnValue = r_gDMean;
-  int    i = 0, j = 0, k = 0;
-  for (i = 0; i < r_gNCoefficients ; i++) {
-    // Evaluate the ith term in the expansion
-    double term = r_gCoefficient[i];
-    for (j = 0; j < r_gNVariables; j++) {
-      // Evaluate the polynomial in the jth variable.
-      int power = r_gPower[r_gNVariables * i + j]; 
-      double p1 = 1, p2 = 0, p3 = 0, r = 0;
-      double v =  1 + 2. / (r_gXMax[j] - r_gXMin[j]) * (x[j] - r_gXMax[j]);
-      // what is the power to use!
-      switch(power) {
-        case 1: r = 1; break; 
-        case 2: r = v; break; 
-        default: 
-                p2 = v; 
-                for (k = 3; k <= power; k++) { 
-                  p3 = p2 * v;
-                  p3 = 2 * v * p2 - p1; 
-                  p1 = p2; p2 = p3; 
-                }
-                r = p3;
-      }
-      // multiply this term by the poly in the jth var
-      term *= r; 
-    }
-    // Add this term to the final result
-    returnValue += term;
-  }
-  return returnValue;
-}
-
-// EOF for simtree_rigidity.C
-// -*- mode: c++ -*-
-// 
-// File simtree_length.C generated by TMultiDimFit::MakeRealCode
-// on Fri May 28 00:25:23 2021
-// ROOT version 5.32/00
-//
-// This file contains the function 
-//
-//    double  l_fit(double *x); 
-//
-// For evaluating the parameterization obtained
-// from TMultiDimFit and the point x
-// 
-// See TMultiDimFit class documentation for more information 
-// 
-//
-// Static data variables
-//
-static int    l_gNVariables    = 6;
-static int    l_gNCoefficients = 77;
-static double l_gDMean         = 9270.11;
-// Assignment to mean vector.
-static double l_gXMean[] = {
-  0.00503955, -0.0144156, -0.101042, -0.00374735, 0.252452, -0.512868 };
-
-// Assignment to minimum vector.
-static double l_gXMin[] = {
-  -64.3965, -52.0748, -0.0335722, -0.0401226, -1065.23, -299.807 };
-
-// Assignment to maximum vector.
-static double l_gXMax[] = {
-  52.3459, 55.8258, 0.0432727, 0.04243, 1147.99, 931.885 };
-
-// Assignment to coefficients vector.
-static double l_gCoefficient[] = {
-  -6590.16,
-  -13914.6,
-  -2835.19,
-  11286.2,
-  -2318.17,
-  1506.93,
-  2806.41,
-  1055.27,
-  15246.9,
-  -582.244,
-  14287.3,
-  2024.12,
-  5826.32,
-  2883.06,
-  -3211.64,
-  -2565.14,
-  -3508.49,
-  -3945.73,
-  81.3328,
-  6976.65,
-  15491.2,
-  3340.02,
-  32961.2,
-  -4257.58,
-  -18445.7,
-  -1204.58,
-  7213.73,
-  1.15164,
-  56.7961,
-  -6.55251,
-  -3.1322,
-  32.3364,
-  17.5881,
-  2.39586,
-  -17.4049,
-  4.99685,
-  -44.3394,
-  -18.7891,
-  -32.238,
-  -26.3708,
-  6.54775,
-  7.12538,
-  -7694.27,
-  -1133.08,
-  27850.4,
-  11364,
-  -5002.96,
-  13573.4,
-  6437.24,
-  -8275.75,
-  -2420.82,
-  -5.63532,
-  -35904.8,
-  64281.3,
-  30190.8,
-  -1440.27,
-  2620.77,
-  1005.06,
-  -17864.9,
-  -5198.67,
-  1.04599,
-  -1.61611,
-  -3.90405,
-  -15.7528,
-  -9.48741,
-  -25.0997,
-  -10129.3,
-  2.95363,
-  99.1311,
-  -34800.1,
-  -2335.5,
-  -8.51776,
-  11.7194,
-  9.32626,
-  -25.8219,
-  2.18449,
-  -0.917242
-};
-
-// Assignment to error coefficients vector.
-static double l_gCoefficientRMS[] = {
-  647.441,
-  1258.13,
-  378.249,
-  1228.53,
-  338.526,
-  439.055,
-  444.451,
-  68.299,
-  1618.83,
-  242.058,
-  929.339,
-  595.678,
-  754.718,
-  687.373,
-  454.652,
-  379.081,
-  542.258,
-  290.313,
-  1.10705,
-  744.15,
-  1703.8,
-  1010.19,
-  2870.46,
-  856.075,
-  1864.98,
-  296.741,
-  874.288,
-  33.2233,
-  23.4443,
-  37.7505,
-  37.4312,
-  100.661,
-  24.1279,
-  7.88562,
-  26.9518,
-  11.868,
-  89.8789,
-  28.9366,
-  32.0384,
-  12.1687,
-  29.4468,
-  4.91973,
-  565.232,
-  471.598,
-  1808.9,
-  1470.56,
-  738.43,
-  1449.61,
-  1967.42,
-  1666.8,
-  577.982,
-  17.9786,
-  3633.48,
-  5592.56,
-  3319.49,
-  257.515,
-  239.116,
-  312.322,
-  1701.52,
-  563.862,
-  9.27568,
-  114.437,
-  17.2737,
-  5.63337,
-  5.32805,
-  24.809,
-  1098.33,
-  25.1201,
-  96.0282,
-  3315.12,
-  242.887,
-  36.0287,
-  32.4886,
-  4.45102,
-  31.7885,
-  20.2812,
-  22.29
-};
-
-// Assignment to powers vector.
-// The powers are stored row-wise, that is
-//  p_ij = l_gPower[i * NVariables + j];
-static int    l_gPower[] = {
-  1,  1,  1,  1,  1,  1,
-  1,  1,  1,  1,  1,  2,
-  2,  1,  1,  1,  1,  1,
-  1,  1,  1,  1,  2,  1,
-  1,  2,  1,  1,  1,  1,
-  1,  1,  1,  2,  1,  1,
-  1,  1,  2,  1,  1,  1,
-  1,  1,  1,  1,  1,  3,
-  1,  1,  1,  1,  2,  2,
-  1,  3,  1,  1,  1,  1,
-  2,  1,  2,  1,  1,  1,
-  1,  1,  1,  2,  1,  2,
-  1,  2,  2,  1,  1,  1,
-  1,  1,  2,  1,  1,  2,
-  1,  2,  1,  1,  1,  2,
-  1,  1,  3,  1,  1,  1,
-  2,  1,  1,  1,  1,  2,
-  3,  1,  1,  1,  1,  1,
-  1,  1,  1,  1,  3,  1,
-  2,  1,  1,  2,  1,  1,
-  1,  1,  1,  2,  2,  1,
-  1,  2,  1,  2,  1,  1,
-  1,  1,  2,  1,  2,  1,
-  1,  1,  2,  2,  1,  1,
-  2,  1,  1,  1,  2,  1,
-  1,  1,  1,  3,  1,  1,
-  1,  1,  1,  1,  2,  3,
-  2,  1,  1,  3,  1,  1,
-  1,  2,  1,  3,  1,  1,
-  1,  1,  2,  3,  1,  1,
-  1,  1,  3,  2,  1,  1,
-  1,  2,  2,  2,  1,  1,
-  2,  2,  1,  1,  2,  1,
-  1,  3,  1,  1,  2,  1,
-  1,  2,  2,  1,  2,  1,
-  1,  1,  3,  1,  2,  1,
-  2,  2,  1,  2,  1,  1,
-  2,  1,  1,  2,  2,  1,
-  1,  2,  1,  2,  2,  1,
-  1,  1,  1,  3,  2,  1,
-  3,  1,  1,  2,  1,  1,
-  2,  1,  1,  1,  3,  1,
-  3,  1,  1,  1,  1,  2,
-  1,  3,  1,  1,  1,  2,
-  2,  1,  2,  1,  1,  2,
-  1,  2,  2,  1,  1,  2,
-  1,  1,  3,  1,  1,  2,
-  2,  1,  1,  2,  1,  2,
-  1,  2,  1,  2,  1,  2,
-  1,  1,  2,  2,  1,  2,
-  1,  1,  1,  3,  1,  2,
-  3,  1,  2,  1,  1,  1,
-  2,  1,  1,  1,  2,  2,
-  1,  1,  2,  1,  2,  2,
-  1,  1,  1,  2,  2,  2,
-  1,  2,  1,  1,  1,  3,
-  1,  1,  2,  1,  1,  3,
-  1,  1,  1,  2,  1,  3,
-  1,  2,  1,  1,  2,  1,
-  2,  2,  1,  1,  1,  1,
-  3,  1,  1,  1,  2,  1,
-  2,  1,  2,  2,  1,  1,
-  1,  3,  1,  2,  1,  1,
-  1,  1,  2,  1,  3,  1,
-  1,  1,  1,  2,  3,  1,
-  3,  2,  1,  1,  1,  1,
-  2,  2,  1,  1,  1,  2,
-  1,  3,  2,  1,  1,  1,
-  2,  2,  2,  1,  1,  1,
-  1,  2,  1,  1,  2,  2,
-  2,  1,  1,  1,  1,  3,
-  2,  1,  2,  1,  2,  1,
-  1,  1,  2,  2,  2,  1,
-  1,  2,  1,  1,  3,  1,
-  1,  2,  3,  1,  1,  1,
-  2,  1,  3,  1,  1,  1,
-  2,  3,  1,  1,  1,  1
-};
-
-// 
-// The function   double l_fit(double *x)
-// 
-double l_fit(double *x) {
-  double returnValue = l_gDMean;
-  int    i = 0, j = 0, k = 0;
-  for (i = 0; i < l_gNCoefficients ; i++) {
-    // Evaluate the ith term in the expansion
-    double term = l_gCoefficient[i];
-    for (j = 0; j < l_gNVariables; j++) {
-      // Evaluate the polynomial in the jth variable.
-      int power = l_gPower[l_gNVariables * i + j]; 
-      double p1 = 1, p2 = 0, p3 = 0, r = 0;
-      double v =  1 + 2. / (l_gXMax[j] - l_gXMin[j]) * (x[j] - l_gXMax[j]);
-      // what is the power to use!
-      switch(power) {
-        case 1: r = 1; break; 
-        case 2: r = v; break; 
-        default: 
-                p2 = v; 
-                for (k = 3; k <= power; k++) { 
-                  p3 = p2 * v;
-                  p3 = 2 * v * p2 - p1; 
-                  p1 = p2; p2 = p3; 
-                }
-                r = p3;
-      }
-      // multiply this term by the poly in the jth var
-      term *= r; 
-    }
-    // Add this term to the final result
-    returnValue += term;
-  }
-  return returnValue;
-}
-
-// EOF for simtree_length.C
-// -*- mode: c++ -*-
-// 
-// File simtree_tof.C generated by TMultiDimFit::MakeRealCode
-// on Fri May 28 00:25:23 2021
-// ROOT version 5.32/00
-//
-// This file contains the function 
-//
-//    double  t_fit(double *x); 
-//
-// For evaluating the parameterization obtained
-// from TMultiDimFit and the point x
-// 
-// See TMultiDimFit class documentation for more information 
-// 
-//
-// Static data variables
-//
-static int    t_gNVariables    = 6;
-static int    t_gNCoefficients = 77;
-static double t_gDMean         = 88.0547;
-// Assignment to mean vector.
-static double t_gXMean[] = {
-  0.00503955, -0.0144156, -0.101042, -0.00374735, 0.252452, -0.512868 };
-
-// Assignment to minimum vector.
-static double t_gXMin[] = {
-  -64.3965, -52.0748, -0.0335722, -0.0401226, -1065.23, -299.807 };
-
-// Assignment to maximum vector.
-static double t_gXMax[] = {
-  52.3459, 55.8258, 0.0432727, 0.04243, 1147.99, 931.885 };
-
-// Assignment to coefficients vector.
-static double t_gCoefficient[] = {
-  282.419,
-  600.69,
-  199.214,
-  -610.381,
-  16.8965,
-  -181.103,
-  -65.4683,
-  27.7219,
-  -586.193,
-  -221.165,
-  131.422,
-  257.64,
-  1.87553,
-  -199.501,
-  -2326.17,
-  1262.92,
-  88.6015,
-  -380.386,
-  -2.50251,
-  0.192506,
-  -0.0682083,
-  -0.0605175,
-  -0.00933975,
-  2.54256,
-  2.13194,
-  1.70986,
-  -0.734912,
-  0.725959,
-  0.566525,
-  1.41423,
-  262.324,
-  54.1268,
-  -1143.9,
-  255.544,
-  -5.4915,
-  -383.569,
-  158.012,
-  177.209,
-  2458.18,
-  2069.51,
-  -4536.21,
-  -1851.8,
-  144.418,
-  21.6852,
-  -147.382,
-  49.0783,
-  -791.882,
-  13.6185,
-  -269.147,
-  75.7361,
-  134.56,
-  -230.177,
-  -950.155,
-  1062.38,
-  81.3635,
-  -0.315919,
-  -1.87225,
-  -0.378093,
-  0.88221,
-  -0.437349,
-  -0.563023,
-  1.42746,
-  -525.101,
-  -447.633,
-  0.297019,
-  20.3115,
-  195.852,
-  0.424978,
-  -0.802113,
-  0.778162,
-  381.743,
-  -0.338946,
-  -0.109268,
-  -0.201413,
-  -0.0377476,
-  -0.150836,
-  0.0571845
-};
-
-// Assignment to error coefficients vector.
-static double t_gCoefficientRMS[] = {
-  51.6658,
-  100.306,
-  31.0942,
-  100.473,
-  36.5499,
-  36.1352,
-  5.48829,
-  19.3199,
-  74.8745,
-  55.5697,
-  31.1116,
-  44.4306,
-  0.110749,
-  81.8808,
-  231.073,
-  150.427,
-  24.4632,
-  71.4537,
-  2.1768,
-  3.60448,
-  0.792011,
-  1.63035,
-  3.6143,
-  8.62908,
-  3.21448,
-  1.21207,
-  0.490659,
-  0.55933,
-  0.532584,
-  2.39502,
-  45.341,
-  37.6402,
-  145.644,
-  60.5957,
-  9.27703,
-  159.431,
-  133.982,
-  47.623,
-  293.067,
-  271.101,
-  450.201,
-  272.399,
-  19.7108,
-  21.145,
-  19.4434,
-  28.0278,
-  132.355,
-  49.4246,
-  60.5858,
-  37.5262,
-  23.2998,
-  58.9152,
-  139.819,
-  139.147,
-  68.8468,
-  3.6068,
-  9.68357,
-  10.9918,
-  2.70767,
-  1.19191,
-  0.444182,
-  3.0704,
-  118.051,
-  114.732,
-  1.74076,
-  25.8244,
-  43.6475,
-  3.15977,
-  2.42206,
-  2.89927,
-  84.993,
-  3.25617,
-  1.97584,
-  2.41757,
-  0.930273,
-  2.8245,
-  2.14083
-};
-
-// Assignment to powers vector.
-// The powers are stored row-wise, that is
-//  p_ij = t_gPower[i * NVariables + j];
-static int    t_gPower[] = {
-  1,  1,  1,  1,  1,  1,
-  1,  1,  1,  1,  1,  2,
-  2,  1,  1,  1,  1,  1,
-  1,  1,  1,  1,  2,  1,
-  1,  1,  1,  2,  1,  1,
-  1,  1,  2,  1,  1,  1,
-  1,  1,  1,  1,  1,  3,
-  1,  3,  1,  1,  1,  1,
-  2,  1,  2,  1,  1,  1,
-  1,  1,  2,  1,  1,  2,
-  1,  1,  3,  1,  1,  1,
-  2,  1,  1,  1,  1,  2,
-  1,  1,  1,  1,  3,  1,
-  1,  2,  1,  2,  1,  1,
-  1,  1,  2,  1,  2,  1,
-  2,  1,  1,  1,  2,  1,
-  1,  1,  1,  3,  1,  1,
-  1,  1,  1,  1,  2,  3,
-  1,  2,  1,  3,  1,  1,
-  1,  1,  3,  2,  1,  1,
-  1,  3,  1,  1,  2,  1,
-  1,  3,  1,  2,  1,  1,
-  2,  1,  2,  1,  2,  1,
-  2,  2,  1,  2,  1,  1,
-  1,  2,  1,  2,  2,  1,
-  1,  1,  1,  3,  2,  1,
-  2,  1,  1,  1,  3,  1,
-  1,  1,  2,  1,  3,  1,
-  1,  1,  1,  2,  3,  1,
-  3,  2,  1,  1,  1,  1,
-  3,  1,  1,  1,  1,  2,
-  1,  3,  1,  1,  1,  2,
-  2,  1,  2,  1,  1,  2,
-  1,  1,  3,  1,  1,  2,
-  2,  2,  2,  1,  1,  1,
-  1,  2,  1,  2,  1,  2,
-  1,  1,  2,  2,  1,  2,
-  1,  1,  1,  3,  1,  2,
-  2,  1,  1,  1,  2,  2,
-  1,  2,  1,  1,  2,  2,
-  1,  1,  2,  1,  2,  2,
-  1,  1,  1,  2,  2,  2,
-  2,  1,  1,  1,  1,  3,
-  1,  2,  1,  1,  1,  3,
-  1,  1,  2,  1,  1,  3,
-  1,  2,  1,  1,  1,  1,
-  1,  1,  1,  1,  2,  2,
-  1,  1,  1,  2,  1,  2,
-  1,  2,  2,  1,  1,  1,
-  1,  2,  1,  1,  1,  2,
-  3,  1,  1,  1,  1,  1,
-  2,  1,  1,  2,  1,  1,
-  1,  1,  1,  2,  2,  1,
-  1,  2,  1,  1,  2,  1,
-  1,  1,  2,  2,  1,  1,
-  1,  1,  2,  3,  1,  1,
-  1,  2,  2,  2,  1,  1,
-  2,  1,  2,  2,  1,  1,
-  1,  2,  2,  1,  2,  1,
-  1,  1,  3,  1,  2,  1,
-  1,  2,  1,  1,  3,  1,
-  1,  2,  3,  1,  1,  1,
-  1,  2,  2,  1,  1,  2,
-  2,  1,  1,  2,  1,  2,
-  3,  1,  2,  1,  1,  1,
-  1,  1,  1,  2,  1,  3,
-  2,  2,  1,  1,  1,  1,
-  2,  1,  1,  3,  1,  1,
-  2,  2,  1,  1,  2,  1,
-  2,  1,  1,  2,  2,  1,
-  2,  2,  1,  1,  1,  2,
-  1,  1,  2,  2,  2,  1,
-  2,  1,  3,  1,  1,  1,
-  1,  3,  2,  1,  1,  1,
-  3,  1,  1,  1,  2,  1,
-  3,  1,  1,  2,  1,  1,
-  2,  3,  1,  1,  1,  1
-};
-
-// 
-// The function   double t_fit(double *x)
-// 
-double t_fit(double *x) {
-  double returnValue = t_gDMean;
-  int    i = 0, j = 0, k = 0;
-  for (i = 0; i < t_gNCoefficients ; i++) {
-    // Evaluate the ith term in the expansion
-    double term = t_gCoefficient[i];
-    for (j = 0; j < t_gNVariables; j++) {
-      // Evaluate the polynomial in the jth variable.
-      int power = t_gPower[t_gNVariables * i + j]; 
-      double p1 = 1, p2 = 0, p3 = 0, r = 0;
-      double v =  1 + 2. / (t_gXMax[j] - t_gXMin[j]) * (x[j] - t_gXMax[j]);
-      // what is the power to use!
-      switch(power) {
-        case 1: r = 1; break; 
-        case 2: r = v; break; 
-        default: 
-                p2 = v; 
-                for (k = 3; k <= power; k++) { 
-                  p3 = p2 * v;
-                  p3 = 2 * v * p2 - p1; 
-                  p1 = p2; p2 = p3; 
-                }
-                r = p3;
-      }
-      // multiply this term by the poly in the jth var
-      term *= r; 
-    }
-    // Add this term to the final result
-    returnValue += term;
+  ////////////////////////////////////////////////////////////////////////////////
+  //            Construct Method to be pass to the AnalysisFactory              //
+  ////////////////////////////////////////////////////////////////////////////////
+  NPL::VAnalysis* Analysis::Construct(){
+    return (NPL::VAnalysis*) new Analysis();
   }
-  return returnValue;
-}
 
-// EOF for simtree_tof.C
-////////////////////////////////////////////////////////////////////////////////
-//            Construct Method to be pass to the AnalysisFactory              //
-////////////////////////////////////////////////////////////////////////////////
-NPL::VAnalysis* Analysis::Construct(){
-  return (NPL::VAnalysis*) new Analysis();
-}
-
-////////////////////////////////////////////////////////////////////////////////
-//            Registering the construct method to the factory                 //
-////////////////////////////////////////////////////////////////////////////////
-extern "C"{
-  class proxy_analysis{
-    public:
-      proxy_analysis(){
-        NPL::AnalysisFactory::getInstance()->SetConstructor(Analysis::Construct);
-      }
-  };
+  ////////////////////////////////////////////////////////////////////////////////
+  //            Registering the construct method to the factory                 //
+  ////////////////////////////////////////////////////////////////////////////////
+  extern "C"{
+    class proxy_analysis{
+      public:
+        proxy_analysis(){
+          NPL::AnalysisFactory::getInstance()->SetConstructor(Analysis::Construct);
+        }
+    };
 
   proxy_analysis p_analysis;
 }
-
diff --git a/Projects/S034/Analysis.h b/Projects/S034/Analysis.h
index 7dfeb6e35196e4257c880fd6e5e278d3b75355af..89f158472a8143871f8b78ea63c0d5963d0d0786 100644
--- a/Projects/S034/Analysis.h
+++ b/Projects/S034/Analysis.h
@@ -22,6 +22,7 @@
  *****************************************************************************/
 
 #include"NPVAnalysis.h"
+#include"NPEnergyLoss.h"
 #include"TMinosPhysics.h"
 #include"TNebulaPhysics.h"
 #include"TSamuraiBDCPhysics.h"
@@ -50,18 +51,18 @@ class Analysis: public NPL::VAnalysis{
     TSamuraiFDC2Physics* FDC2;
     TSamuraiHodoscopePhysics* Hodo;
     SamuraiFieldMap m_field ;
-    ofstream file;
+//    ofstream file;
   private: // output variable
-    double Brho,BrhoP,BDCX,BDCY,Z;
+    double Brho,BDCX,BDCY,X,Y,Z,Erel;
     double Beta_f;
+    double Beta_n;
     int    Trigger;
+  private: // Energy loss table
+   NPL::EnergyLoss FragmentTarget ;
   public:
     void  Clear();
     void  InitOutputBranch();
     void  InitInputBranch();
-
 };
 
-// use for broh calculation
-double r_fit(double *x) ;
 #endif
diff --git a/Projects/S034/Calibration/Nebula/Nebula_T.txt b/Projects/S034/Calibration/Nebula/Nebula_T.txt
index d975f25e3cd53e9ba1a8ef65ad0439d733c3e6c7..e041de79d7d5107053887cc9d46add802b1bbaec 100644
--- a/Projects/S034/Calibration/Nebula/Nebula_T.txt
+++ b/Projects/S034/Calibration/Nebula/Nebula_T.txt
@@ -1,120 +1,120 @@
-NEBULA_ID1_T 19.2852
-NEBULA_ID2_T 19.2846
-NEBULA_ID3_T 19.4056
-NEBULA_ID4_T 19.1069
-NEBULA_ID5_T 19.4713
-NEBULA_ID6_T 19.2863
-NEBULA_ID7_T 19.2007
-NEBULA_ID8_T 19.3881
-NEBULA_ID9_T 19.1476
-NEBULA_ID10_T 19.0105
-NEBULA_ID11_T 19.1384
-NEBULA_ID12_T 19.1565
-NEBULA_ID13_T 18.9115
-NEBULA_ID14_T 19.0286
-NEBULA_ID15_T 19.4379
-NEBULA_ID16_T 18.9817
-NEBULA_ID17_T 18.8992
-NEBULA_ID18_T 18.6782
-NEBULA_ID19_T 19.0459
-NEBULA_ID20_T 18.8453
-NEBULA_ID21_T 18.7874
-NEBULA_ID22_T 18.9267
-NEBULA_ID23_T 18.6575
-NEBULA_ID24_T 18.8422
-NEBULA_ID25_T 19.072
-NEBULA_ID26_T 19.5044
-NEBULA_ID27_T 19.4049
-NEBULA_ID28_T 19.0738
-NEBULA_ID29_T 19.3946
-NEBULA_ID30_T 19.4726
-NEBULA_ID31_T 19.121
-NEBULA_ID32_T 19.1251
-NEBULA_ID33_T 19.2008
-NEBULA_ID34_T 18.9232
-NEBULA_ID35_T 18.8517
-NEBULA_ID36_T 18.9349
-NEBULA_ID37_T 19.2574
-NEBULA_ID38_T 18.7397
-NEBULA_ID39_T 19.0286
-NEBULA_ID40_T 18.5307
-NEBULA_ID41_T 19.006
-NEBULA_ID42_T 18.634
-NEBULA_ID43_T 19.2158
-NEBULA_ID44_T 18.8532
-NEBULA_ID45_T 18.8342
-NEBULA_ID46_T 18.9721
-NEBULA_ID47_T 18.968
-NEBULA_ID48_T 19.1271
-NEBULA_ID49_T 18.816
-NEBULA_ID50_T 18.952
-NEBULA_ID51_T 18.9897
-NEBULA_ID52_T 19.0074
-NEBULA_ID53_T 19.2299
-NEBULA_ID54_T 19.0185
-NEBULA_ID55_T 19.2308
-NEBULA_ID56_T 19.0939
-NEBULA_ID57_T 19.3677
-NEBULA_ID58_T 19.5316
-NEBULA_ID59_T 19.8612
-NEBULA_ID60_T 19.6356
-NEBULA_ID61_T 18.6837
-NEBULA_ID62_T 18.6049
-NEBULA_ID63_T 20.3056
-NEBULA_ID64_T 19.1702
-NEBULA_ID65_T 19.3322
-NEBULA_ID66_T 19.6835
-NEBULA_ID67_T 19.0867
-NEBULA_ID68_T 18.6978
-NEBULA_ID69_T 17.9037
-NEBULA_ID70_T 18.9826
-NEBULA_ID71_T 18.6267
-NEBULA_ID72_T 18.2024
-NEBULA_ID73_T 19.0313
-NEBULA_ID74_T 18.8189
-NEBULA_ID75_T 18.6591
-NEBULA_ID76_T 19.3614
-NEBULA_ID77_T 19.2587
-NEBULA_ID78_T 18.1176
-NEBULA_ID79_T 18.7292
-NEBULA_ID80_T 19.3139
-NEBULA_ID81_T 18.8809
-NEBULA_ID82_T 17.6727
-NEBULA_ID83_T 19.1507
-NEBULA_ID84_T 18.8959
-NEBULA_ID85_T 18.0162
-NEBULA_ID86_T 19.2333
-NEBULA_ID87_T 18.5977
-NEBULA_ID88_T 17.6151
-NEBULA_ID89_T 19.0762
-NEBULA_ID90_T 18.4975
-NEBULA_ID91_T 19.0267
-NEBULA_ID92_T 19.0031
-NEBULA_ID93_T 19.0938
-NEBULA_ID94_T 18.7937
-NEBULA_ID95_T 19.1556
-NEBULA_ID96_T 19.1339
-NEBULA_ID97_T 19.3789
-NEBULA_ID98_T 18.9279
-NEBULA_ID99_T 19.3174
-NEBULA_ID100_T 19.3137
-NEBULA_ID101_T 19.0579
-NEBULA_ID102_T 18.775
-NEBULA_ID103_T 19.3587
-NEBULA_ID104_T 18.874
-NEBULA_ID105_T 19.4789
-NEBULA_ID106_T 19.4282
-NEBULA_ID107_T 19.1013
-NEBULA_ID108_T 18.9991
-NEBULA_ID109_T 19.0751
-NEBULA_ID110_T 18.0316
-NEBULA_ID111_T 18.9544
-NEBULA_ID112_T 19.4328
-NEBULA_ID113_T 18.7139
-NEBULA_ID114_T 19.893
-NEBULA_ID115_T 18.3269
-NEBULA_ID116_T 18.545
-NEBULA_ID117_T 20.0339
-NEBULA_ID118_T 18.9662
-NEBULA_ID119_T 20.0821
-NEBULA_ID120_T 19.4934
+NEBULA_T_ID1 19.2852
+NEBULA_T_ID2 19.2846
+NEBULA_T_ID3 19.4056
+NEBULA_T_ID4 19.1069
+NEBULA_T_ID5 19.4713
+NEBULA_T_ID6 19.2863
+NEBULA_T_ID7 19.2007
+NEBULA_T_ID8 19.3881
+NEBULA_T_ID9 19.1476
+NEBULA_T_ID10 19.0105
+NEBULA_T_ID11 19.1384
+NEBULA_T_ID12 19.1565
+NEBULA_T_ID13 18.9115
+NEBULA_T_ID14 19.0286
+NEBULA_T_ID15 19.4379
+NEBULA_T_ID16 18.9817
+NEBULA_T_ID17 18.8992
+NEBULA_T_ID18 18.6782
+NEBULA_T_ID19 19.0459
+NEBULA_T_ID20 18.8453
+NEBULA_T_ID21 18.7874
+NEBULA_T_ID22 18.9267
+NEBULA_T_ID23 18.6575
+NEBULA_T_ID24 18.8422
+NEBULA_T_ID25 19.072
+NEBULA_T_ID26 19.5044
+NEBULA_T_ID27 19.4049
+NEBULA_T_ID28 19.0738
+NEBULA_T_ID29 19.3946
+NEBULA_T_ID30 19.4726
+NEBULA_T_ID31 19.121
+NEBULA_T_ID32 19.1251
+NEBULA_T_ID33 19.2008
+NEBULA_T_ID34 18.9232
+NEBULA_T_ID35 18.8517
+NEBULA_T_ID36 18.9349
+NEBULA_T_ID37 19.2574
+NEBULA_T_ID38 18.7397
+NEBULA_T_ID39 19.0286
+NEBULA_T_ID40 18.5307
+NEBULA_T_ID41 19.006
+NEBULA_T_ID42 18.634
+NEBULA_T_ID43 19.2158
+NEBULA_T_ID44 18.8532
+NEBULA_T_ID45 18.8342
+NEBULA_T_ID46 18.9721
+NEBULA_T_ID47 18.968
+NEBULA_T_ID48 19.1271
+NEBULA_T_ID49 18.816
+NEBULA_T_ID50 18.952
+NEBULA_T_ID51 18.9897
+NEBULA_T_ID52 19.0074
+NEBULA_T_ID53 19.2299
+NEBULA_T_ID54 19.0185
+NEBULA_T_ID55 19.2308
+NEBULA_T_ID56 19.0939
+NEBULA_T_ID57 19.3677
+NEBULA_T_ID58 19.5316
+NEBULA_T_ID59 19.8612
+NEBULA_T_ID60 19.6356
+NEBULA_T_ID61 18.6837
+NEBULA_T_ID62 18.6049
+NEBULA_T_ID63 20.3056
+NEBULA_T_ID64 19.1702
+NEBULA_T_ID65 19.3322
+NEBULA_T_ID66 19.6835
+NEBULA_T_ID67 19.0867
+NEBULA_T_ID68 18.6978
+NEBULA_T_ID69 17.9037
+NEBULA_T_ID70 18.9826
+NEBULA_T_ID71 18.6267
+NEBULA_T_ID72 18.2024
+NEBULA_T_ID73 19.0313
+NEBULA_T_ID74 18.8189
+NEBULA_T_ID75 18.6591
+NEBULA_T_ID76 19.3614
+NEBULA_T_ID77 19.2587
+NEBULA_T_ID78 18.1176
+NEBULA_T_ID79 18.7292
+NEBULA_T_ID80 19.3139
+NEBULA_T_ID81 18.8809
+NEBULA_T_ID82 17.6727
+NEBULA_T_ID83 19.1507
+NEBULA_T_ID84 18.8959
+NEBULA_T_ID85 18.0162
+NEBULA_T_ID86 19.2333
+NEBULA_T_ID87 18.5977
+NEBULA_T_ID88 17.6151
+NEBULA_T_ID89 19.0762
+NEBULA_T_ID90 18.4975
+NEBULA_T_ID91 19.0267
+NEBULA_T_ID92 19.0031
+NEBULA_T_ID93 19.0938
+NEBULA_T_ID94 18.7937
+NEBULA_T_ID95 19.1556
+NEBULA_T_ID96 19.1339
+NEBULA_T_ID97 19.3789
+NEBULA_T_ID98 18.9279
+NEBULA_T_ID99 19.3174
+NEBULA_T_ID100 19.3137
+NEBULA_T_ID101 19.0579
+NEBULA_T_ID102 18.775
+NEBULA_T_ID103 19.3587
+NEBULA_T_ID104 18.874
+NEBULA_T_ID105 19.4789
+NEBULA_T_ID106 19.4282
+NEBULA_T_ID107 19.1013
+NEBULA_T_ID108 18.9991
+NEBULA_T_ID109 19.0751
+NEBULA_T_ID110 18.0316
+NEBULA_T_ID111 18.9544
+NEBULA_T_ID112 19.4328
+NEBULA_T_ID113 18.7139
+NEBULA_T_ID114 19.893
+NEBULA_T_ID115 18.3269
+NEBULA_T_ID116 18.545
+NEBULA_T_ID117 20.0339
+NEBULA_T_ID118 18.9662
+NEBULA_T_ID119 20.0821
+NEBULA_T_ID120 19.4934
diff --git a/Projects/S034/Calibration/Nebula/Nebula_Y.txt b/Projects/S034/Calibration/Nebula/Nebula_Y.txt
index 95254e626efb70d3975f76dd45658cd0b07effc0..55ea73a36f96c0c69165b80a0b4ed79312435ef1 100644
--- a/Projects/S034/Calibration/Nebula/Nebula_Y.txt
+++ b/Projects/S034/Calibration/Nebula/Nebula_Y.txt
@@ -1,120 +1,120 @@
-NEBULA_ID1_Y 1.38638
-NEBULA_ID2_Y -3.8305
-NEBULA_ID3_Y 3.07288
-NEBULA_ID4_Y -4.23115
-NEBULA_ID5_Y -2.90823
-NEBULA_ID6_Y 2.06144
-NEBULA_ID7_Y -6.14942
-NEBULA_ID8_Y -4.20327
-NEBULA_ID9_Y 13.9879
-NEBULA_ID10_Y -3.06973
-NEBULA_ID11_Y -1.99039
-NEBULA_ID12_Y -4.86424
-NEBULA_ID13_Y 5.29633
-NEBULA_ID14_Y -0.662125
-NEBULA_ID15_Y 3.04338
-NEBULA_ID16_Y -4.29256
-NEBULA_ID17_Y -2.30441
-NEBULA_ID18_Y 0.160267
-NEBULA_ID19_Y -9.17392
-NEBULA_ID20_Y 7.05093
-NEBULA_ID21_Y 2.09011
-NEBULA_ID22_Y -6.41303
-NEBULA_ID23_Y 0.530808
-NEBULA_ID24_Y 4.12012
-NEBULA_ID25_Y 0.363813
-NEBULA_ID26_Y -4.18461
-NEBULA_ID27_Y -4.88908
-NEBULA_ID28_Y 1.10303
-NEBULA_ID29_Y 6.93376
-NEBULA_ID30_Y -0.585477
-NEBULA_ID31_Y 3.24728
-NEBULA_ID32_Y -0.673857
-NEBULA_ID33_Y -0.00853288
-NEBULA_ID34_Y 0.436468
-NEBULA_ID35_Y 2.1838
-NEBULA_ID36_Y 1.45213
-NEBULA_ID37_Y 3.2373
-NEBULA_ID38_Y -2.325
-NEBULA_ID39_Y -3.20755
-NEBULA_ID40_Y -1.26772
-NEBULA_ID41_Y 0.887659
-NEBULA_ID42_Y 2.37744
-NEBULA_ID43_Y -1.16329
-NEBULA_ID44_Y -3.6847
-NEBULA_ID45_Y 3.10103
-NEBULA_ID46_Y 2.02808
-NEBULA_ID47_Y -0.509373
-NEBULA_ID48_Y 6.99985
-NEBULA_ID49_Y -1.29247
-NEBULA_ID50_Y -0.304723
-NEBULA_ID51_Y 3.23443
-NEBULA_ID52_Y -5.04529
-NEBULA_ID53_Y -1.669
-NEBULA_ID54_Y 4.07673
-NEBULA_ID55_Y -2.16088
-NEBULA_ID56_Y 4.10968
-NEBULA_ID57_Y -5.66633
-NEBULA_ID58_Y 5.72757
-NEBULA_ID59_Y -0.284415
-NEBULA_ID60_Y 4.61204
-NEBULA_ID61_Y -11.215
-NEBULA_ID62_Y 2.20483
-NEBULA_ID63_Y 3.65855
-NEBULA_ID64_Y -6.00011
-NEBULA_ID65_Y -4.83289
-NEBULA_ID66_Y 4.86074
-NEBULA_ID67_Y -5.8572
-NEBULA_ID68_Y 4.79375
-NEBULA_ID69_Y -4.77623
-NEBULA_ID70_Y -2.30597
-NEBULA_ID71_Y -5.44358
-NEBULA_ID72_Y 5.76698
-NEBULA_ID73_Y -4.52093
-NEBULA_ID74_Y 0.669192
-NEBULA_ID75_Y 5.2536
-NEBULA_ID76_Y -3.06509
-NEBULA_ID77_Y -4.24188
-NEBULA_ID78_Y 2.02627
-NEBULA_ID79_Y -4.03842
-NEBULA_ID80_Y 3.1809
-NEBULA_ID81_Y -5.90454
-NEBULA_ID82_Y 1.81077
-NEBULA_ID83_Y -5.46129
-NEBULA_ID84_Y -1.16729
-NEBULA_ID85_Y 2.94024
-NEBULA_ID86_Y -3.34849
-NEBULA_ID87_Y -8.60917
-NEBULA_ID88_Y -1.84303
-NEBULA_ID89_Y 1.29457
-NEBULA_ID90_Y -1.06823
-NEBULA_ID91_Y 5.69821
-NEBULA_ID92_Y -5.7119
-NEBULA_ID93_Y 7.594
-NEBULA_ID94_Y 4.47444
-NEBULA_ID95_Y 3.8241
-NEBULA_ID96_Y -10.6418
-NEBULA_ID97_Y 14.1286
-NEBULA_ID98_Y -5.54334
-NEBULA_ID99_Y 5.50645
-NEBULA_ID100_Y -5.75486
-NEBULA_ID101_Y -3.05549
-NEBULA_ID102_Y 7.65411
-NEBULA_ID103_Y 0.353753
-NEBULA_ID104_Y -5.75768
-NEBULA_ID105_Y 4.14771
-NEBULA_ID106_Y -0.914708
-NEBULA_ID107_Y 3.08513
-NEBULA_ID108_Y 0.902261
-NEBULA_ID109_Y -0.553475
-NEBULA_ID110_Y 5.95105
-NEBULA_ID111_Y -2.47195
-NEBULA_ID112_Y -5.24212
-NEBULA_ID113_Y 13.0663
-NEBULA_ID114_Y -2.03506
-NEBULA_ID115_Y 5.80118
-NEBULA_ID116_Y -2.23628
-NEBULA_ID117_Y 6.53639
-NEBULA_ID118_Y 2.23348
-NEBULA_ID119_Y -6.16759
-NEBULA_ID120_Y 0.972052
+NEBULA_Y_ID1 1.38638
+NEBULA_Y_ID2 -3.8305
+NEBULA_Y_ID3 3.07288
+NEBULA_Y_ID4 -4.23115
+NEBULA_Y_ID5 -2.90823
+NEBULA_Y_ID6 2.06144
+NEBULA_Y_ID7 -6.14942
+NEBULA_Y_ID8 -4.20327
+NEBULA_Y_ID9 13.9879
+NEBULA_Y_ID10 -3.06973
+NEBULA_Y_ID11 -1.99039
+NEBULA_Y_ID12 -4.86424
+NEBULA_Y_ID13 5.29633
+NEBULA_Y_ID14 -0.662125
+NEBULA_Y_ID15 3.04338
+NEBULA_Y_ID16 -4.29256
+NEBULA_Y_ID17 -2.30441
+NEBULA_Y_ID18 0.160267
+NEBULA_Y_ID19 -9.17392
+NEBULA_Y_ID20 7.05093
+NEBULA_Y_ID21 2.09011
+NEBULA_Y_ID22 -6.41303
+NEBULA_Y_ID23 0.530808
+NEBULA_Y_ID24 4.12012
+NEBULA_Y_ID25 0.363813
+NEBULA_Y_ID26 -4.18461
+NEBULA_Y_ID27 -4.88908
+NEBULA_Y_ID28 1.10303
+NEBULA_Y_ID29 6.93376
+NEBULA_Y_ID30 -0.585477
+NEBULA_Y_ID31 3.24728
+NEBULA_Y_ID32 -0.673857
+NEBULA_Y_ID33 -0.00853288
+NEBULA_Y_ID34 0.436468
+NEBULA_Y_ID35 2.1838
+NEBULA_Y_ID36 1.45213
+NEBULA_Y_ID37 3.2373
+NEBULA_Y_ID38 -2.325
+NEBULA_Y_ID39 -3.20755
+NEBULA_Y_ID40 -1.26772
+NEBULA_Y_ID41 0.887659
+NEBULA_Y_ID42 2.37744
+NEBULA_Y_ID43 -1.16329
+NEBULA_Y_ID44 -3.6847
+NEBULA_Y_ID45 3.10103
+NEBULA_Y_ID46 2.02808
+NEBULA_Y_ID47 -0.509373
+NEBULA_Y_ID48 6.99985
+NEBULA_Y_ID49 -1.29247
+NEBULA_Y_ID50 -0.304723
+NEBULA_Y_ID51 3.23443
+NEBULA_Y_ID52 -5.04529
+NEBULA_Y_ID53 -1.669
+NEBULA_Y_ID54 4.07673
+NEBULA_Y_ID55 -2.16088
+NEBULA_Y_ID56 4.10968
+NEBULA_Y_ID57 -5.66633
+NEBULA_Y_ID58 5.72757
+NEBULA_Y_ID59 -0.284415
+NEBULA_Y_ID60 4.61204
+NEBULA_Y_ID61 -11.215
+NEBULA_Y_ID62 2.20483
+NEBULA_Y_ID63 3.65855
+NEBULA_Y_ID64 -6.00011
+NEBULA_Y_ID65 -4.83289
+NEBULA_Y_ID66 4.86074
+NEBULA_Y_ID67 -5.8572
+NEBULA_Y_ID68 4.79375
+NEBULA_Y_ID69 -4.77623
+NEBULA_Y_ID70 -2.30597
+NEBULA_Y_ID71 -5.44358
+NEBULA_Y_ID72 5.76698
+NEBULA_Y_ID73 -4.52093
+NEBULA_Y_ID74 0.669192
+NEBULA_Y_ID75 5.2536
+NEBULA_Y_ID76 -3.06509
+NEBULA_Y_ID77 -4.24188
+NEBULA_Y_ID78 2.02627
+NEBULA_Y_ID79 -4.03842
+NEBULA_Y_ID80 3.1809
+NEBULA_Y_ID81 -5.90454
+NEBULA_Y_ID82 1.81077
+NEBULA_Y_ID83 -5.46129
+NEBULA_Y_ID84 -1.16729
+NEBULA_Y_ID85 2.94024
+NEBULA_Y_ID86 -3.34849
+NEBULA_Y_ID87 -8.60917
+NEBULA_Y_ID88 -1.84303
+NEBULA_Y_ID89 1.29457
+NEBULA_Y_ID90 -1.06823
+NEBULA_Y_ID91 5.69821
+NEBULA_Y_ID92 -5.7119
+NEBULA_Y_ID93 7.594
+NEBULA_Y_ID94 4.47444
+NEBULA_Y_ID95 3.8241
+NEBULA_Y_ID96 -10.6418
+NEBULA_Y_ID97 14.1286
+NEBULA_Y_ID98 -5.54334
+NEBULA_Y_ID99 5.50645
+NEBULA_Y_ID100 -5.75486
+NEBULA_Y_ID101 -3.05549
+NEBULA_Y_ID102 7.65411
+NEBULA_Y_ID103 0.353753
+NEBULA_Y_ID104 -5.75768
+NEBULA_Y_ID105 4.14771
+NEBULA_Y_ID106 -0.914708
+NEBULA_Y_ID107 3.08513
+NEBULA_Y_ID108 0.902261
+NEBULA_Y_ID109 -0.553475
+NEBULA_Y_ID110 5.95105
+NEBULA_Y_ID111 -2.47195
+NEBULA_Y_ID112 -5.24212
+NEBULA_Y_ID113 13.0663
+NEBULA_Y_ID114 -2.03506
+NEBULA_Y_ID115 5.80118
+NEBULA_Y_ID116 -2.23628
+NEBULA_Y_ID117 6.53639
+NEBULA_Y_ID118 2.23348
+NEBULA_Y_ID119 -6.16759
+NEBULA_Y_ID120 0.972052
diff --git a/Projects/S034/Calibration/Nebula/convert.cxx b/Projects/S034/Calibration/Nebula/convert.cxx
index 79372a24fef45de76b1a4d3f9c7fe4117a3dfd6b..6be586f7d600a62409316bf27696c8276bcb4f38 100644
--- a/Projects/S034/Calibration/Nebula/convert.cxx
+++ b/Projects/S034/Calibration/Nebula/convert.cxx
@@ -7,7 +7,7 @@ void convert(){
  double val;
  getline(in,buffer); // ignore first line
  while(in >> id >> buffer){
-    out << "NEBULA_ID" << id << "_T " << buffer  << endl;
+    out << "NEBULA_T_ID" << id << " " << buffer  << endl;
  }
  out.close();
  in.close();
@@ -16,7 +16,7 @@ void convert(){
 
  getline(in,buffer); // ignore first line
  while(in >> id >> buffer){
-    out << "NEBULA_ID" << id << "_Y " << buffer  << endl;
+    out << "NEBULA_Y_ID" << id << " " << buffer  << endl;
  }
  out.close();
  in.close();
diff --git a/Projects/S034/db/s034_list.txt b/Projects/S034/db/s034_list.txt
index 2113fec8a1b74851c3b19af0503238d37be54e4f..1200b20f93f7f321c45bb7f1a5980996d2460d7c 100644
--- a/Projects/S034/db/s034_list.txt
+++ b/Projects/S034/db/s034_list.txt
@@ -7,3 +7,4 @@ db/SAMURAIFDC0_20200109.xml
 db/SAMURAIHOD_s034_all40mV_s037_20170702.xml
 db/SAMURAIBDC1.xml
 db/SAMURAIBDC2.xml
+db/SAMURAIPlastic.s034.2.xml
diff --git a/Projects/S034/s034.detector b/Projects/S034/s034.detector
index f0c1a6b7d82a2e99446b87da9344154b19482b9e..6ff1e04865f26125ca96b3b02f641e33eecf1001 100644
--- a/Projects/S034/s034.detector
+++ b/Projects/S034/s034.detector
@@ -9,8 +9,6 @@ SAMURAIBDC 1
 %%%%%%%%%%%%
 SAMURAIBDC 2
  XML= db/SAMURAIBDC2.xml
-% Offset= +0.25 +1.6 -5876.7 mm
-% Offset= 0 0 -5876.7 mm
  Offset= -0.0686464 0.294288 -5876.7 mm
  InvertX= 0 
  InvertY= 0