diff --git a/Projects/e793s/Analysis.cxx b/Projects/e793s/Analysis.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..b2d69a906e58b281938f0d3fd3741316d2344913
--- /dev/null
+++ b/Projects/e793s/Analysis.cxx
@@ -0,0 +1,659 @@
+/*****************************************************************************
+ * Copyright (C) 2009-2014    this file is part of the NPTool Project        *
+ *                                                                           *
+ * For the licensing terms see $NPTOOL/Licence/NPTool_Licence                *
+ * For the list of contributors see $NPTOOL/Licence/Contributors             *
+ *****************************************************************************/
+
+/*****************************************************************************
+ * Original Author: Adrien MATTA  contact address: a.matta@surrey.ac.uk      *
+ * Edited for e793s                                                          *
+ *                                                                           *
+ * Creation Date  :                                                          *
+ * Last update    : May 2021                                                 *
+ *---------------------------------------------------------------------------*
+ * Decription:                                                               *
+ * Class describing the property of an Analysis object                       *
+ *                                                                           *
+ *---------------------------------------------------------------------------*
+ * Comment:                                                                  *
+ *                                                                           *
+ *                                                                           *
+ *****************************************************************************/
+#include<iostream>
+#include<string>
+
+#include<fstream>
+
+using namespace std;
+
+#include "Analysis.h"
+#include "NPFunction.h"
+#include "NPAnalysisFactory.h"
+#include "NPDetectorManager.h"
+#include "NPOptionManager.h"
+
+////////////////////////////////////////////////////////////////////////////////
+Analysis::Analysis(){
+}
+////////////////////////////////////////////////////////////////////////////////
+Analysis::~Analysis(){
+}
+
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::Init() {
+  ///////////////////////////////////////////////////////////////////////////////  
+ 
+
+  agata_zShift=51*mm;
+  //BrhoRef=0.65;
+  
+  // initialize input and output branches
+  InitOutputBranch();
+  InitInputBranch();
+  
+  M2 = (TMust2Physics*)  m_DetectorManager -> GetDetector("M2Telescope");
+  MG = (TMugastPhysics*) m_DetectorManager -> GetDetector("Mugast");
+  ML = (TModularLeafPhysics*) m_DetectorManager -> GetDetector("ModularLeaf"); 
+//  CATS = (TCATSPhysics*) m_DetectorManager->GetDetector("CATSDetector");
+
+  // get reaction information
+  reaction.ReadConfigurationFile(NPOptionManager::getInstance()->GetReactionFile());
+  OriginalBeamEnergy = reaction.GetBeamEnergy();
+  // target thickness
+  TargetThickness = m_DetectorManager->GetTargetThickness();
+  string TargetMaterial = m_DetectorManager->GetTargetMaterial();
+  // energy losses
+  string light = NPL::ChangeNameToG4Standard(reaction.GetParticle3()->GetName());
+  string beam = NPL::ChangeNameToG4Standard(reaction.GetParticle1()->GetName());
+  LightTarget = NPL::EnergyLoss(light+"_"+TargetMaterial+".G4table","G4Table",100 );
+  LightAl = NPL::EnergyLoss(light+"_Al.G4table" ,"G4Table",100);
+  LightSi = NPL::EnergyLoss(light+"_Si.G4table" ,"G4Table",100);
+  BeamTargetELoss = NPL::EnergyLoss(beam+"_"+TargetMaterial+".G4table","G4Table",100);
+
+  FinalBeamEnergy = BeamTargetELoss.Slow(OriginalBeamEnergy, 0.5*TargetThickness, 0);
+  reaction.SetBeamEnergy(FinalBeamEnergy); 
+
+  cout << "Beam energy at mid-target: " << FinalBeamEnergy << endl;
+
+  // initialize various parameters
+  Rand = TRandom3();
+  ParticleMult = 0;
+  GammaMult = 0;
+  DetectorNumber = 0;
+  ThetaNormalTarget = 0;
+  ThetaM2Surface = 0;
+  ThetaMGSurface = 0;
+  Si_E_M2 = 0;
+  CsI_E_M2 = 0;
+  Energy = 0;
+  ThetaGDSurface = 0;
+  elab_tmp = 0;
+  thetalab_tmp = 0;
+  philab_tmp = 0;
+  ThetaGDSurface = 0;
+  X.clear();
+  Y.clear();
+  Z.clear();
+  dE = 0;
+  BeamDirection = TVector3(0,0,1);
+  //nbTrack=0;
+  //nbHits=0;
+  //count=0;
+  //AHeavy=reaction.GetNucleus4()->GetA();
+  AHeavy=reaction.GetParticle4()->GetA();
+  //ALight=reaction.GetNucleus3()->GetA(); 
+  ALight=reaction.GetParticle3()->GetA(); 
+  //MHeavy=reaction.GetNucleus4()->Mass();
+  MHeavy=reaction.GetParticle4()->Mass();
+  //MLight=reaction.GetNucleus3()->Mass();
+  MLight=reaction.GetParticle3()->Mass();
+
+  for(int i=0;i<GATCONF_SIZE;i++){ // loop over the bits
+    GATCONF_Counter[i] = 0 ; 
+  }
+
+}
+
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::TreatEvent(){
+  
+  // Reinitiate calculated variable
+  ReInitValue();
+  GATCONF_MASTER=ML->GetCalibratedValue("GATCONF_MASTER");
+  for(int i=0;i<GATCONF_SIZE;i++){ // loop over the bits
+    if(GATCONF_MASTER & (unsigned int)pow(2,i)){ // test if ith bit is on
+      GATCONF_Counter[i]++ ; // increment the array
+      //      cout <<  std::bitset<16>(GATCONF_MASTER) << " " << i+1 << endl;
+    }
+  }
+
+  double xbeam = 0.0; 
+  double ybeam = 0.0;
+
+  // double XTarget = CATS->GetPositionOnTarget().X();
+  // double YTarget = CATS->GetPositionOnTarget().Y();
+  TVector3 BeamDirection(xbeam,ybeam,1);
+  BeamImpact = TVector3(xbeam,ybeam,m_DetectorManager->GetTargetZ()); 
+
+  ParticleMult=M2->Si_E.size()+MG->DSSD_E.size();
+  //ParticleMult=M2->Si_E.size();
+  //  FinalBeamEnergy=BeamCD2.Slow(OriginalBeamEnergy,0.5*TargetThickness,BeamDirection.Angle(TVector(0,0,1)));
+  //reaction.SetBeamEnergy(FinalBeamEnergy);
+
+  ////////////////////////////////////////////////////////////////////////////
+  //////////////////////////////// LOOP on MUST2  ////////////////////////////
+  ////////////////////////////////////////////////////////////////////////////
+  unsigned int sizeM2 = M2->Si_E.size();
+  for(unsigned int countMust2 = 0; countMust2 < sizeM2; countMust2++){
+
+    /************************************************/
+    // Part 0 : Get the useful Data
+    // MUST2
+    int TelescopeNumber = M2->TelescopeNumber[countMust2];
+
+    /************************************************/
+    // Part 1 : Impact Angle
+    ThetaM2Surface = 0;
+    ThetaNormalTarget = 0;
+    thetalab_tmp = 0;
+    philab_tmp = 0;
+    TVector3 HitDirection = M2 -> GetPositionOfInteraction(countMust2) - BeamImpact ;
+    thetalab_tmp = HitDirection.Angle(BeamDirection);
+    philab_tmp = HitDirection.Phi();
+
+    X.push_back(M2->GetPositionOfInteraction(countMust2).X());
+    Y.push_back(M2->GetPositionOfInteraction(countMust2).Y());
+    Z.push_back(M2->GetPositionOfInteraction(countMust2).Z());
+
+    ThetaM2Surface = HitDirection.Angle(- M2->GetTelescopeNormal(countMust2) );
+    ThetaNormalTarget = HitDirection.Angle( TVector3(xbeam,ybeam,1) ) ;
+
+   // cout<<"Must2 Znormal:"<<M2 -> GetTelescopeNormal(countMust2).Z()<<endl;
+  //  cout<<"Must2 telescope:"<<M2->TelescopeNumber[countMust2]<<endl;
+
+    /************************************************/
+    // Part 2 : Impact Energy
+    Energy = 0;
+    elab_tmp = 0;
+    Si_E_M2 = M2->Si_E[countMust2];
+    CsI_E_M2= M2->CsI_E[countMust2];
+
+    // if CsI
+    if(CsI_E_M2>0 ){
+      // The energy in CsI is calculate form dE/dx Table because
+      Energy = CsI_E_M2;
+      Energy = LightAl.EvaluateInitialEnergy(Energy, 0.4*micrometer, ThetaM2Surface);
+      Energy+=Si_E_M2;
+    }
+
+    else
+      Energy = Si_E_M2;
+
+
+    RawEnergy.push_back(Energy); //CPx ADDITION
+
+
+    // Evaluate energy using the thickness
+    elab_tmp = LightAl.EvaluateInitialEnergy(Energy, 0.4*micrometer, ThetaM2Surface);
+    // Target Correction
+    elab_tmp = LightTarget.EvaluateInitialEnergy(elab_tmp, 0.5*TargetThickness, ThetaNormalTarget);
+    ELab.push_back(elab_tmp);
+
+    /************************************************/
+
+    /************************************************/
+    // Part 3 : Excitation Energy Calculation
+    Ex.push_back(reaction.ReconstructRelativistic(elab_tmp,thetalab_tmp));
+
+    Ecm.push_back(Energy*(AHeavy+ALight)/(4*AHeavy*cos(thetalab_tmp)*cos(thetalab_tmp)));
+    /************************************************/
+
+    /************************************************/
+    // Part 4 : Theta CM Calculation
+    
+    ThetaCM.push_back(reaction.EnergyLabToThetaCM(elab_tmp, thetalab_tmp)/deg);
+    /************************************************/
+
+    ThetaLab.push_back(thetalab_tmp/deg);
+    PhiLab.push_back(philab_tmp/deg);
+
+  }
+
+  ////////////////////////////////////////////////////////////////////////////
+  //////////////////////////////// LOOP on MUGAST ////////////////////////////
+  ////////////////////////////////////////////////////////////////////////////
+  unsigned int sizeMG = MG->DSSD_E.size();
+  for(unsigned int countMugast = 0; countMugast<sizeMG; countMugast++){
+    // Part 1 : Impact Angle
+    ThetaMGSurface = 0;
+    ThetaNormalTarget = 0;
+    thetalab_tmp = 0;
+    philab_tmp = 0;
+    TVector3 HitDirection = MG->GetPositionOfInteraction(countMugast) - BeamImpact;
+//    cout << HitDirection.X() << "\t" << HitDirection.Y() << "\t" << HitDirection.Z() << "\t\t" << MG->GetPositionOfInteraction(countMugast).X() << "\t" << MG->GetPositionOfInteraction(countMugast).Y() << "\t" << MG->GetPositionOfInteraction(countMugast).Z() << "\t\t" << BeamImpact.X() << "\t" << BeamImpact.Y() << "\t" << BeamImpact.Z() << endl;
+    //HitDirection = [Xp, Yp, Zp] - [Xb, Yb, Zb] = [Xd, Td, Zd] IN MY CODE!
+    thetalab_tmp = HitDirection.Angle(BeamDirection);
+    philab_tmp = HitDirection.Phi();
+
+    X.push_back( MG -> GetPositionOfInteraction(countMugast).X());
+    Y.push_back( MG -> GetPositionOfInteraction(countMugast).Y());
+    Z.push_back( MG -> GetPositionOfInteraction(countMugast).Z());
+
+    //ThetaMGSurface = HitDirection.Angle( TVector3(0,0,1) ) ;
+    ThetaMGSurface = HitDirection.Angle( MG -> GetTelescopeNormal(countMugast) );
+    ThetaNormalTarget = HitDirection.Angle( TVector3(xbeam,ybeam,-1) ) ;
+
+    // Part 2 : Impact Energy
+    Energy = elab_tmp = 0;
+    Energy = MG->GetEnergyDeposit(countMugast);
+
+    RawEnergy.push_back(Energy);
+
+    elab_tmp = LightAl.EvaluateInitialEnergy(Energy, 0.4*micrometer, ThetaMGSurface);
+   
+    elab_tmp = LightTarget.EvaluateInitialEnergy(elab_tmp, 0.5*TargetThickness, 0.);
+    //elab_tmp   = LightTarget.EvaluateInitialEnergy( elab_tmp ,TargetThickness*0.5, ThetaNormalTarget);
+    ELab.push_back(elab_tmp);
+
+    // Part 3 : Excitation Energy Calculation
+    Ex.push_back(reaction.ReconstructRelativistic(elab_tmp,thetalab_tmp));
+    Ecm.push_back(elab_tmp*(AHeavy+ALight)/(4*AHeavy*cos(thetalab_tmp)*cos(thetalab_tmp)));
+    // Part 4 : Theta CM Calculation
+    ThetaCM.push_back(reaction.EnergyLabToThetaCM(elab_tmp, thetalab_tmp)/deg);
+    ThetaLab.push_back(thetalab_tmp/deg);
+     
+    PhiLab.push_back(philab_tmp/deg);
+
+
+
+//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+/*
+cout << MG->TelescopeNumber[countMugast] << " "
+     << MG->GetPositionOfInteraction(countMugast).X() << " " 
+     << MG->GetPositionOfInteraction(countMugast).Y() << " " 
+     << MG->GetPositionOfInteraction(countMugast).Z() << " " 
+     << Energy <<  " " 
+     << endl; 
+*/
+
+///////////!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+
+
+
+
+
+
+    if(sizeMG==1){
+      MG_T = MG->DSSD_T[0];
+      MG_E = MG->DSSD_E[0];
+      MG_X = MG->DSSD_X[0];
+      MG_Y = MG->DSSD_Y[0];
+      MG_D = MG->TelescopeNumber[0];
+    }
+
+  }//end loop Mugast
+
+  ////////////////////////////////////////////////////////////////////////////
+  ///////////////////////////////// LOOP on AGATA ////////////////////////////
+  ////////////////////////////////////////////////////////////////////////////
+  // Agata by track
+  /*
+  for(int j=0; j<nbTrack; j++){ // all multiplicity
+    TLorentzVector GammaLV;
+    // Measured E
+    double Egamma=trackE[j]/1000.; // From keV to MeV 
+    // Gamma detection position
+    // TrackZ1 to be corrected there is a shift of +51mm
+    TVector3 GammaHit(trackX1[j],trackY1[j],trackZ1[j]+agata_zShift); 
+    // TVector3 GammaHit(trackX1[0],trackY1[0],trackZ1[0]); 
+    // Gamma Direction 
+    TVector3 GammaDirection = GammaHit-BeamImpact;
+    GammaDirection = GammaDirection.Unit();
+    // Beta from Two body kinematic
+    //TVector3 beta = reaction.GetEnergyImpulsionLab_4().BoostVector();
+    // Beta from the Beam mid target 
+    reaction.GetKinematicLine4();
+    TVector3 beta(0,0,-reaction.GetNucleus4()->GetBeta());
+
+    double ThetaGamma = GammaDirection.Angle(BeamDirection)/deg;
+    // Construct LV in lab frame
+    GammaLV.SetPx(Egamma*GammaDirection.X());
+    GammaLV.SetPy(Egamma*GammaDirection.Y());
+    GammaLV.SetPz(Egamma*GammaDirection.Z());
+    GammaLV.SetE(Egamma);
+    // Boost back in CM
+    GammaLV.Boost(beta);
+    // Get EDC
+    EDC.push_back(GammaLV.Energy());
+  }
+  */
+  // Agata add back is not always multiplicity 1 ?? NO, not necessarily!
+  for(int i=0; i<nbAdd; i++){
+  //if(nbAdd==1){
+    TLorentzVector GammaLV;
+    // Measured E
+    double Egamma = AddE[i]/1000.; // From keV to MeV 
+    // Gamma detection position
+    // TrackZ1 to be corrected there is a shift of +51mm
+    TVector3 GammaHit(AddX[i],AddY[i],AddZ[i]+agata_zShift); 
+    // TVector3 GammaHit(trackX1[0],trackY1[0],trackZ1[0]); 
+    // Gamma Direction 
+    TVector3 GammaDirection = GammaHit-BeamImpact;
+    GammaDirection = GammaDirection.Unit();
+    // Beta from Two body kinematic
+    //TVector3 beta = reaction.GetEnergyImpulsionLab_4().BoostVector();
+    // Beta from the Beam mid target 
+    reaction.GetKinematicLine4();
+    // TVector3 beta(0,0,-reaction.GetNucleus4()->GetBeta());
+    TVector3 beta(0,0,-0.1257);
+
+    double ThetaGamma = GammaDirection.Angle(BeamDirection)/deg;
+    // Construct LV in lab frame
+    GammaLV.SetPx(Egamma*GammaDirection.X());
+    GammaLV.SetPy(Egamma*GammaDirection.Y());
+    GammaLV.SetPz(Egamma*GammaDirection.Z());
+    GammaLV.SetE(Egamma);
+    // Boost back in CM
+    GammaLV.Boost(beta);
+    // Get EDC
+    AddBack_EDC.push_back(GammaLV.Energy());
+
+  }
+
+  float offset[21] = {0.,27.4,46.7,41.1,71.2,58.7,32.1,40.4,-46.3,15.1,9.9,-39.1,-39.8,-5.3,18.5,0.7,-28.5,-7,-22,-8.5,-8.5};
+  float slope[21] = {1.,0.96,0.93,0.93,0.89,0.91,0.95,0.93,1.06,0.97,0.99,1.05,1.04,1.01,0.97,1.,1.03,0.99,1.,1.02,1.02};
+
+  for(unsigned int i=0; i<MW_Nr; i++){
+    MWT[i] = -1500;
+    for(unsigned short j=0; j<21; j++)
+      if(MW_N[i]==j)
+        MWT[i] = (MW_T[i]-offset[j])/slope[j];
+  }
+
+}
+
+  ////////////////////////////////////////////////////////////////////////////////
+  void Analysis::End(){
+    cout << endl << "\e[1;32m GATCONF Statistics " << endl ;
+    for(int i=0;i<GATCONF_SIZE;i++){ // loop over the bits
+      cout << GATCONF_Label[i]  << "\t\t" << GATCONF_Counter[i] << endl ; 
+    }
+    cout << endl ;
+  }
+
+  ////////////////////////////////////////////////////////////////////////////////
+  void Analysis::InitOutputBranch(){
+    //RootOutput::getInstance()->GetTree()->Branch("Ex",&Ex,"Ex/D");
+    RootOutput::getInstance()->GetTree()->Branch("Ex",&Ex);
+    //RootOutput::getInstance()->GetTree()->Branch("EDC",&EDC,"EDC/D");
+    RootOutput::getInstance()->GetTree()->Branch("EDC",&EDC);
+    RootOutput::getInstance()->GetTree()->Branch("AddBack_EDC",&AddBack_EDC);
+    RootOutput::getInstance()->GetTree()->Branch("EAgata",&EAgata,"EAgata/D");
+    RootOutput::getInstance()->GetTree()->Branch("ELab",&ELab);
+    RootOutput::getInstance()->GetTree()->Branch("Ecm",&Ecm);
+    RootOutput::getInstance()->GetTree()->Branch("RawEnergy",&RawEnergy); // CPx ADDITION
+    RootOutput::getInstance()->GetTree()->Branch("ThetaLab",&ThetaLab);
+    RootOutput::getInstance()->GetTree()->Branch("PhiLab",&PhiLab);
+    RootOutput::getInstance()->GetTree()->Branch("ThetaCM",&ThetaCM);
+    RootOutput::getInstance()->GetTree()->Branch("Run",&Run,"Run/I");
+    RootOutput::getInstance()->GetTree()->Branch("X",&X);
+    RootOutput::getInstance()->GetTree()->Branch("Y",&Y);
+    RootOutput::getInstance()->GetTree()->Branch("Z",&Z);
+    RootOutput::getInstance()->GetTree()->Branch("dE",&dE,"dE/D");
+    RootOutput::getInstance()->GetTree()->Branch("MG_T",MG_T);
+    RootOutput::getInstance()->GetTree()->Branch("MG_E",MG_E);
+    RootOutput::getInstance()->GetTree()->Branch("MG_X",MG_X);
+    RootOutput::getInstance()->GetTree()->Branch("MG_Y",MG_Y);
+    RootOutput::getInstance()->GetTree()->Branch("MG_D",MG_D);
+
+    // Vamos 
+    RootOutput::getInstance()->GetTree()->Branch("LTS",&LTS,"LTS/l");
+/*    RootOutput::getInstance()->GetTree()->Branch("T_FPMW_CATS1",&T_FPMW_CATS1,"T_FPMW_CATS1/s");
+    RootOutput::getInstance()->GetTree()->Branch("T_FPMW_CATS2_C",&T_FPMW_CATS2_C,"T_FPMW_CATS2_C/F");
+    RootOutput::getInstance()->GetTree()->Branch("T_FPMW_HF",&T_FPMW_HF,"T_FPMW_CATS1/s");
+    RootOutput::getInstance()->GetTree()->Branch("T_FPMW_HF_C",&T_FPMW_HF_C,"T_FPMW_CATS1/F");
+    RootOutput::getInstance()->GetTree()->Branch ("IC",IC, "IC[12]/F");
+    RootOutput::getInstance()->GetTree()->Branch ("ICRawPUMult01",ICRawPUMult01, "ICRawPUMult01[7]/s");
+
+    RootOutput::getInstance()->GetTree()->Branch("T_CATS2_HF",&T_CATS2_HF,"T_CATS2_HF/s");
+    RootOutput::getInstance()->GetTree()->Branch("T_MUGAST_FPMW",&T_MUGAST_FPMW,"T_MUGAST_FPMW/s");
+    RootOutput::getInstance()->GetTree()->Branch("T_FPMW_CATS2",&T_FPMW_CATS2,"T_FPMW_CATS2/s");
+    RootOutput::getInstance()->GetTree()->Branch("T_FPMW_HF",&T_FPMW_HF,"T_FPMW_HF/s");
+
+    // DC0
+    RootOutput::getInstance()->GetTree()->Branch ("DC0_X", &DC0_X, "DC0_X/F");
+    RootOutput::getInstance()->GetTree()->Branch ("DC0_Y", &DC0_Y, "DC0_Y/F");
+    // DC1
+    RootOutput::getInstance()->GetTree()->Branch ("DC1_X", &DC1_X, "DC1_X/F");
+    RootOutput::getInstance()->GetTree()->Branch ("DC1_Y", &DC1_Y, "DC1_Y/F");
+    // DC2
+    RootOutput::getInstance()->GetTree()->Branch ("DC2_X", &DC2_X, "DC2_X/F");
+    RootOutput::getInstance()->GetTree()->Branch ("DC2_Y", &DC2_Y, "DC2_Y/F");
+    // DC3
+    RootOutput::getInstance()->GetTree()->Branch ("DC3_X", &DC3_X, "DC3_X/F");
+    RootOutput::getInstance()->GetTree()->Branch ("DC3_Y", &DC3_Y, "DC3_Y/F");
+
+    RootOutput::getInstance()->GetTree()->Branch ("Xf", &Xf, "Xf/F");
+    RootOutput::getInstance()->GetTree()->Branch ("Tf", &Tf, "Tf/F");
+
+    RootOutput::getInstance()->GetTree()->Branch ("Brho", &Brho, "Brho/F");
+    RootOutput::getInstance()->GetTree()->Branch ("Theta", &Theta, "Theta/F");
+    RootOutput::getInstance()->GetTree()->Branch ("Phi", &Phi, "Phi/F");
+    RootOutput::getInstance()->GetTree()->Branch ("Path", &Path, "Path/F");
+
+    RootOutput::getInstance()->GetTree()->Branch ("EWIRE_1_1", &EWIRE_1_1, "EWIRE_1_1/s");
+    RootOutput::getInstance()->GetTree()->Branch ("EWIRE_1_2", &EWIRE_1_2, "EWIRE_1_2/s");
+    RootOutput::getInstance()->GetTree()->Branch ("EWIRE_2_1", &EWIRE_2_1, "EWIRE_2_1/s");
+    RootOutput::getInstance()->GetTree()->Branch ("EWIRE_2_2", &EWIRE_2_2, "EWIRE_2_2/s");
+*/
+    RootOutput::getInstance()->GetTree()->Branch ("MW_Nr", &MW_Nr, "MW_Nr/I");
+    RootOutput::getInstance()->GetTree()->Branch ("MW_T", MW_T, "MW_T[MW_Nr]/F");
+    RootOutput::getInstance()->GetTree()->Branch ("MW_N", MW_N, "MW_N[MW_Nr]/s");
+    RootOutput::getInstance()->GetTree()->Branch ("MWT", MWT, "MWT[MW_Nr]/F");
+
+    RootOutput::getInstance()->GetTree()->Branch ("AGAVA_VAMOSTS", &AGAVA_VAMOSTS, "AGAVA_VAMOSTS/l");
+/*    RootOutput::getInstance()->GetTree()->Branch("mT",&mT,"mT/D");
+    RootOutput::getInstance()->GetTree()->Branch("mV",&mV,"mV/D");
+    RootOutput::getInstance()->GetTree()->Branch("mD",&mD,"mD/D");
+    RootOutput::getInstance()->GetTree()->Branch("mBeta",&mBeta,"mBeta/D");
+    RootOutput::getInstance()->GetTree()->Branch("mGamma",&mGamma,"mGamma/D");
+    RootOutput::getInstance()->GetTree()->Branch("mM_Q",&mM_Q,"mM_Q/D");
+    RootOutput::getInstance()->GetTree()->Branch("mM",&mM,"mM/D");
+    RootOutput::getInstance()->GetTree()->Branch("mE",&mE,"mE/D");
+    RootOutput::getInstance()->GetTree()->Branch("mdE",&mdE,"mdE/D");
+*/
+
+    // DC0->3_X/_Y Xf, Tf, Brho, Path, EWIRE_1_1 all, MW_Nr, MW_T, AGAVA_VAMOSTS, 
+    // Agata
+    // Time stamp of the agata trigger
+    /*
+    RootOutput::getInstance()->GetTree()->Branch("TStrack",&TStrack,"TStrack/l");
+
+    // Array of reconstructed tracks
+    RootOutput::getInstance()->GetTree()->Branch("nbTrack",&nbTrack,"nbTrack/I");
+    RootOutput::getInstance()->GetTree()->Branch("trackE",trackE,"trackE[nbTrack]/F");
+    RootOutput::getInstance()->GetTree()->Branch("trackX1",trackX1,"trackX1[nbTrack]/F");
+    RootOutput::getInstance()->GetTree()->Branch("trackY1",trackY1,"trackY1[nbTrack]/F");
+    RootOutput::getInstance()->GetTree()->Branch("trackZ1",trackZ1,"trackZ1[nbTrack]/F");
+    RootOutput::getInstance()->GetTree()->Branch("trackT",trackT,"trackT[nbTrack]/F");
+    RootOutput::getInstance()->GetTree()->Branch("trackCrystalID",trackCrystalID,"trackCrystalID[nbTrack]/I");
+    */
+    
+    //Agata Addback
+    RootOutput::getInstance()->GetTree()->Branch("nbAdd",&nbAdd,"nbAdd/I");
+    RootOutput::getInstance()->GetTree()->Branch("TSHit",&TSHit,"TSHit/l");
+    RootOutput::getInstance()->GetTree()->Branch("AddE",AddE,"AddE[nbAdd]/F");
+    RootOutput::getInstance()->GetTree()->Branch("AddX",AddX,"AddX[nbAdd]/F");
+    RootOutput::getInstance()->GetTree()->Branch("AddY",AddY,"AddY[nbAdd]/F");
+    RootOutput::getInstance()->GetTree()->Branch("AddZ",AddZ,"AddZ[nbAdd]/F");
+  }
+
+
+  ////////////////////////////////////////////////////////////////////////////////
+  void Analysis::InitInputBranch(){
+    SetBranchStatus();
+    // RootInput:: getInstance()->GetChain()->SetBranchAddress("GATCONF",&vGATCONF);
+    //
+    // Vamos
+    RootInput::getInstance()->GetChain()->SetBranchAddress("LTS",&LTS);
+/*    RootInput::getInstance()->GetChain()->SetBranchAddress("T_FPMW_CATS1",&T_FPMW_CATS1);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("T_FPMW_CATS2_C",&T_FPMW_CATS2_C);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("T_FPMW_HF_C",&T_FPMW_HF_C);
+
+    RootInput::getInstance()->GetChain()->SetBranchAddress("T_CATS2_HF",&T_CATS2_HF);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("T_MUGAST_FPMW",&T_MUGAST_FPMW);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("T_FPMW_CATS2",&T_FPMW_CATS2);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("T_FPMW_HF",&T_FPMW_HF);
+
+    RootInput::getInstance()->GetChain()->SetBranchAddress("Xf", &Xf);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("Tf", &Tf);
+
+    RootInput::getInstance()->GetChain()->SetBranchAddress("Brho", &Brho);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("Theta", &Theta);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("Phi", &Phi);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("Path", &Path);
+
+    RootInput::getInstance()->GetChain()->SetBranchAddress("EWIRE_1_1", &EWIRE_1_1);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("EWIRE_1_2", &EWIRE_1_2);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("EWIRE_2_1", &EWIRE_2_1);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("EWIRE_2_2", &EWIRE_2_2);
+*/
+    RootInput::getInstance()->GetChain()->SetBranchAddress("MWTVM", &MW_Nr);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("MWTV", MW_T);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("MWTVN", MW_N);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("AGAVA_VAMOSTS", &AGAVA_VAMOSTS);
+
+    // Agata
+     /*
+    RootInput::getInstance()->GetChain()->SetBranchAddress("TStrack",&TStrack);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("nbTrack",&nbTrack);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("trackE",trackE);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("trackX1",trackX1);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("trackY1",trackY1);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("trackZ1",trackZ1);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("trackT",trackT);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("trackCrystalID",trackCrystalID);
+    */
+    RootInput::getInstance()->GetChain()->SetBranchAddress("TSHit",&TSHit);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("nbAdd",&nbAdd);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("AddE",AddE);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("AddX",AddX);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("AddY",AddY);
+    RootInput::getInstance()->GetChain()->SetBranchAddress("AddZ",AddZ);
+  }
+
+  void Analysis::SetBranchStatus(){
+    // Set Branch status 
+    RootInput::getInstance()->GetChain()->SetBranchStatus("LTS",true);
+ /*   RootInput::getInstance()->GetChain()->SetBranchStatus("T_FPMW_CATS1",true);
+    RootInput::getInstance()->GetChain()->SetBranchStatus("T_FPMW_CATS2_C",true);
+    RootInput::getInstance()->GetChain()->SetBranchStatus("T_FPMW_HF",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("T_FPMW_HF_C",true);
+    
+    RootInput::getInstance()->GetChain()->SetBranchStatus("T_CATS2_HF",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("T_MUGAST_FPMW",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("T_FPMW_CATS2",true );
+
+    RootInput::getInstance()->GetChain()->SetBranchStatus("Xf",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("Tf",true );
+
+    RootInput::getInstance()->GetChain()->SetBranchStatus("Brho",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("Path",true );
+
+    RootInput::getInstance()->GetChain()->SetBranchStatus("EWIRE_1_1",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("EWIRE_1_2",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("EWIRE_2_1",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("EWIRE_2_2",true );
+*/
+    RootInput::getInstance()->GetChain()->SetBranchStatus("MW_Nr",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("MW_T",true  );
+
+    RootInput::getInstance()->GetChain()->SetBranchStatus("MWTV",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("MWTVN",true  );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("MWTVM",true  );
+
+    RootInput::getInstance()->GetChain()->SetBranchStatus("AGAVA_VAMOSTS",true  );
+    // Agata
+    /*
+    RootInput::getInstance()->GetChain()->SetBranchStatus("TStrack",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("nbTrack",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("trackE",true  );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("trackX1",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("trackY1",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("trackZ1",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("trackT",true  );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("trackCrystalID",true);
+    */
+    RootInput::getInstance()->GetChain()->SetBranchStatus("TSHit",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("nbAdd",true );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("AddE",true  );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("AddX",true  );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("AddY",true  );
+    RootInput::getInstance()->GetChain()->SetBranchStatus("AddZ",true  );
+  }
+
+  ////////////////////////////////////////////////////////////////////////////////
+  void Analysis::ReInitValue(){
+    Ex.clear();
+    Ecm.clear();
+    AddBack_EDC.clear();
+    EAgata = -1000;
+    ELab.clear();
+    RawEnergy.clear(); //CPx ADDITION
+    ThetaLab.clear();
+    PhiLab.clear();
+    ThetaCM.clear();
+    //relative_angle = -1000;
+    //sum_angle = -1000;
+    X.clear();
+    Y.clear();
+    Z.clear();
+    dE= -1000;
+    ParticleMult = 0;
+    GammaMult = 0;
+    DetectorNumber = 0;
+    ThetaNormalTarget = 0;
+    ThetaM2Surface = 0;
+    ThetaMGSurface = 0;
+    Si_E_M2 = 0;
+    CsI_E_M2 = 0;
+    Energy = 0;
+    ThetaGDSurface = 0;
+    BeamDirection = TVector3(0,0,1);
+    elab_tmp = 0;
+    thetalab_tmp = 0;
+    philab_tmp = 0;
+
+    MG_T=-1000;
+    MG_E=-1000;
+    MG_X=-1000;
+    MG_Y=-1000;
+    MG_D=-1000;
+
+    //EDC= -1000;
+    EDC.clear();
+  }
+
+  ////////////////////////////////////////////////////////////////////////////////
+  //            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);
+        }
+    };
+
+  proxy_analysis p_analysis;
+}
diff --git a/Projects/e793s/Analysis.h b/Projects/e793s/Analysis.h
new file mode 100644
index 0000000000000000000000000000000000000000..3dfb5fcf308dd31110300d5987a2b1d199255437
--- /dev/null
+++ b/Projects/e793s/Analysis.h
@@ -0,0 +1,219 @@
+#ifndef Analysis_h 
+#define Analysis_h
+/*****************************************************************************
+ * Copyright (C) 2009-2014    this file is part of the NPTool Project        *
+ *                                                                           *
+ * For the licensing terms see $NPTOOL/Licence/NPTool_Licence                *
+ * For the list of contributors see $NPTOOL/Licence/Contributors             *
+ *****************************************************************************/
+
+/*****************************************************************************
+ * Original Author: Adrien MATTA  contact address: a.matta@surrey.ac.uk      *
+ * Edited for e793s                                                          *
+ *                                                                           *
+ * Creation Date  :                                                          *
+ * Last update    : May 2021                                                 *
+ *---------------------------------------------------------------------------*
+ * Decription:                                                               *
+ * Class describing the property of an Analysis object                       *
+ *                                                                           *
+ *---------------------------------------------------------------------------*
+ * Comment:                                                                  *
+ *                                                                           *
+ *                                                                           *
+ *****************************************************************************/
+#include"NPVAnalysis.h"
+#include"NPEnergyLoss.h"
+#include"NPReaction.h"
+#include"RootOutput.h"
+#include"RootInput.h"
+#include "TMust2Physics.h"
+#include "TMugastPhysics.h"
+#include "TCATSPhysics.h"
+//#include "../../NPLib/Detectors/ModularLeaf/TModularLeafPhysics.h"
+#include "TModularLeafPhysics.h"
+#include <TRandom3.h>
+#include <TVector3.h>
+#include <TMath.h>
+#include <bitset>
+
+class Analysis: public NPL::VAnalysis{
+  public:
+    Analysis();
+    ~Analysis();
+
+  public: 
+    void Init();
+    void TreatEvent();
+    void End();
+
+    void InitOutputBranch();
+    void InitInputBranch();
+    void SetBranchStatus();
+    void CalculateVamos();
+    void ReInitValue();
+    static NPL::VAnalysis* Construct();
+
+  private:
+    unsigned int ParticleMult;
+    unsigned int GammaMult;
+    //double EDC;
+    std::vector<double> EDC;
+    vector<double> AddBack_EDC;
+    double EAgata;
+    std::vector<double> ELab;
+    std::vector<double> Ex;
+    std::vector<double> Ecm;
+    std::vector<double> RawEnergy; //CPx ADDITION
+    std::vector<double> ThetaLab;
+    std::vector<double> PhiLab;
+    std::vector<double> ThetaCM;
+    
+    NPL::Reaction reaction;
+    //	Energy loss table: the G4Table are generated by the simulation
+    NPL::EnergyLoss LightTarget;
+    NPL::EnergyLoss LightAl;
+    NPL::EnergyLoss LightSi;
+    NPL::EnergyLoss BeamTargetELoss;
+
+    double TargetThickness ;
+    double WindowsThickness;
+    // Beam Energy
+    double OriginalBeamEnergy ; // AMEV
+    double FinalBeamEnergy; 
+
+    // intermediate variable
+    TVector3 BeamDirection;
+    TVector3 BeamImpact;
+    TRandom3 Rand;
+    int Run;
+    int DetectorNumber;
+    double ThetaNormalTarget;
+    double ThetaM2Surface;
+    double ThetaMGSurface;
+    double elab_tmp;
+    double thetalab_tmp;
+    double philab_tmp;
+    double Si_E_M2;
+    double CsI_E_M2;
+    double Energy;
+#define GATCONF_SIZE 16
+    unsigned int GATCONF_Counter[GATCONF_SIZE];
+    string GATCONF_Label[GATCONF_SIZE] = { 
+      "MUVI1",
+      "MUVI2",
+      "MUVI3",
+      "",
+      "CATS2%",
+      "VAMOS%",
+      "MUST2%",
+      "",
+      "",
+      "VAMOS",
+      "ORAGATA", 
+    } ;
+    
+    double MG_T;
+    double MG_E;
+    int    MG_X;
+    int    MG_Y;
+    int    MG_D;
+
+    double ThetaGDSurface ;
+    std::vector<double> X ;
+    std::vector<double> Y ;
+    std::vector<double> Z ;
+    int AHeavy;
+    int ALight;
+    double MHeavy;
+    double MLight;
+
+    // Vamos Branches
+    unsigned long long int LTS;
+    unsigned short T_FPMW_CATS1; 
+    float T_FPMW_CATS2_C; 
+    unsigned short T_FPMW_HF; 
+    unsigned short T_CATS2_HF;
+    unsigned short T_MUGAST_FPMW;
+    unsigned short T_FPMW_CATS2;
+    float T_FPMW_HF_C; 
+    float IC[12]; 
+    unsigned short ICRawPUMult01[7]; 
+
+    float DC0_X;
+    float DC0_Y;
+
+    float DC1_X;
+    float DC1_Y;
+
+    float DC2_X;
+    float DC2_Y;
+
+    float DC3_X;
+    float DC3_Y;
+
+    float Xf;
+    float Tf;
+
+    float Brho;
+    float Theta;
+    float Phi;
+    float Path;
+
+    unsigned short EWIRE_1_1;
+    unsigned short EWIRE_1_2;
+    unsigned short EWIRE_2_1;
+    unsigned short EWIRE_2_2;
+
+    int MW_Nr;
+    float MW_T[1000];
+    float MWT[1000];
+    unsigned short MW_N[1000];
+    unsigned long long int AGAVA_VAMOSTS;
+
+    // Vamos Calculated
+    double BrhoRef;
+    double mT;    
+    double mV;   
+    double mD;    
+    double mBeta;
+    double mGamma;
+    double mM_Q;  
+    double mM;    
+    double mE;    
+    double mdE;   
+
+    // Agata branches
+    double agata_zShift;
+    /*
+    unsigned long long int TStrack;
+    int nbHits;
+    int nbTrack;
+    float *trackE  = new float(1000);
+    float *trackX1 = new float(1000);
+    float *trackY1 = new float(1000);
+    float *trackZ1 = new float(1000);
+    float *trackT  = new float(1000);
+    int *trackCrystalID = new int(1000);
+    */
+    int nbAdd;
+    unsigned long long int TSHit;
+    float AddE[1000] ;
+    float AddX[1000] ;
+    float AddY[1000] ;
+    float AddZ[1000] ;
+
+    //
+    double dE;
+    double dTheta;
+    // Branches and detectors
+    TMust2Physics* M2;
+    TMugastPhysics* MG;
+    //TCATSPhysics* CATS;
+    TModularLeafPhysics* ML;
+    unsigned int GATCONF_MASTER;
+
+    unsigned long long int count ;
+
+};
+#endif
diff --git a/Projects/e793s/CMakeLists.txt b/Projects/e793s/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..18d9ab375d9946eb246d297aacc7aee6813ff8a6
--- /dev/null
+++ b/Projects/e793s/CMakeLists.txt
@@ -0,0 +1,6 @@
+cmake_minimum_required (VERSION 2.8) 
+project(e793s)
+# Setting the policy to match Cmake version
+cmake_policy(VERSION ${CMAKE_MAJOR_VERSION}.${CMAKE_MINOR_VERSION})
+# include the default NPAnalysis cmake file
+include("../../NPLib/ressources/CMake/NPAnalysis.cmake")
diff --git a/Projects/e793s/Detector/mugast.detector b/Projects/e793s/Detector/mugast.detector
new file mode 100644
index 0000000000000000000000000000000000000000..e478739d677249516fff6af4184b237caa597ee9
--- /dev/null
+++ b/Projects/e793s/Detector/mugast.detector
@@ -0,0 +1,139 @@
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Target
+ THICKNESS= 3.2 micrometer
+ ANGLE= 0 deg
+ RADIUS=  24 mm
+ MATERIAL= CD2
+ X= 0 mm
+ Y= 0 mm
+ Z= 0 mm
+% CATS POSITION ARE FROM E748, so WRONG
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+CATSDetector
+ X1_Y1= -35.86  -35.26  -2658 mm
+ X28_Y1= 35.26  -35.26  -2658 mm
+ X1_Y28= -35.86  35.86  -2658 mm
+ X28_Y28= 35.26  35.86  -2658 mm
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+CATSDetector
+ X1_Y1= 34.86    -35.86   -2045 mm
+ X28_Y1= -36.26  -35.86   -2045 mm
+ X1_Y28= 34.86    35.26   -2045 mm
+ X28_Y28= -36.26  35.26   -2045 mm
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+
+%%%%%%% Telescope 1 %%%%%%%
+M2Telescope
+ X128_Y128=  115.88  9.61  154.54
+ X128_Y1=  104.8  101.89  125.09
+ X1_Y1=  14.55  102.4  160.63
+ X1_Y128=  25.63  10.12  190.08
+ SI=  1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 2 %%%%%%%
+M2Telescope
+ X128_Y128=  -11.23  102.42  160.87
+ X128_Y1=  -101.39  102.39  124.37
+ X1_Y1=  -113.17  10.36  153.56
+ X1_Y128=  -23.03  10.38  190.05
+ SI=  1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 3 %%%%%%%
+M2Telescope
+ X128_Y128=  -113.28  -12.52  153.32
+ X128_Y1=  -101.58  -104.77  124.82
+ X1_Y1=  -11.39  -104.58  161.48
+ X1_Y128=  -23.1  -12.34  189.98
+ SI= 1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 4 %%%%%%%
+M2Telescope
+ X128_Y128=  13.82  -104.92  160.72
+ X128_Y1=  104.3  -104.95  125.08
+ X1_Y1=  115.75  -12.73  153.76
+ X1_Y128=  25.23  -12.65  189.43
+ SI= 1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid 
+ DetectorNumber= 7
+ X128_Y128=  40.99  21.69  -98.21
+ X1_Y128=  23.34  39.63  -98.04
+ X1_Y1=  57.53  120.53  -32.58
+ X128_Y1=  122.2  55.73  -32.67
+
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid 
+ DetectorNumber= 3
+ X128_Y128=  -48.05  -23.23  -95.05
+ X1_Y128=  -24.83  -39.35  -99.33
+ X1_Y1=  -58.24  -119.46  -32.57
+ X128_Y1=  -122.86  -54.55  -32.41
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 1
+ X128_Y128=  -23.05  42.83  -99.66
+ X1_Y128=  -40.74  25.21  -99.65
+ X1_Y1=  -120.65  58.67  -32.39
+ X128_Y1=  -55.66  123.17  -32.55
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 5
+ X128_Y128=  21.68  -41.48  -100.19
+ X1_Y128=  39.55  -23.81  -100.11
+ X1_Y1=  119.55  -57.08  -33.23
+ X128_Y1=  54.81  -121.7  -33.12
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 4
+ X128_Y128=  -15.11  -44.73  -99.25
+ X1_Y128=  10.11  -44.71  -99.35
+ X1_Y1=  43.37  -125.07  -32.74
+ X128_Y1=  -48.18  -125.07  -32.45
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Square      
+ DetectorNumber= 10
+ X128_Y128=  111.36  97.74  89.82
+ X1_Y128=  147.47  9.83  90
+ X1_Y1=  147.7  9.88  -1.79
+ X128_Y1=  111.63  97.7  -2.02
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Square 
+ DetectorNumber= 9
+ X128_Y128= 107.49 -95.88 0.0
+ X1_Y128=  143.81 -8.21 0.0
+ X1_Y1= 143.81 -8.21 91.7 
+ X128_Y1= 107.49 -95.88 91.7
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Annular  
+ DetectorNumber= 11
+ Center=  0 0 -125 mm
+
+ModularLeaf 
+ DefaultValue= -5000
+ Leafs= GATCONF_MASTER GATCONF_SLAVE GATCONF_SLAVE2 CONFDEC_AGATA CONFDEC_VAMOS DATATRIG_CATS1 DATATRIG_CATS2 ADC1_9 CORREL_INVAMOS ADC_CATS_0
+ X= ADC1_9              ADC1_9               ADC_CATS_0
+ Y= CORREL_INVAMOS    ADC_CATS_0           CORREL_INVAMOS
+
diff --git a/Projects/e793s/Detector/mugast_Charlie05May_fixMG3.detector b/Projects/e793s/Detector/mugast_Charlie05May_fixMG3.detector
new file mode 100755
index 0000000000000000000000000000000000000000..570f57c679fb1594e63d37a748d278a930e7947c
--- /dev/null
+++ b/Projects/e793s/Detector/mugast_Charlie05May_fixMG3.detector
@@ -0,0 +1,161 @@
+%%%%%%%%%%Detector%%%%%%%%%%%%%%%%%%%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%1
+Target
+ THICKNESS= 4.76 micrometer
+ ANGLE= 0 deg
+ RADIUS= 10 mm
+ MATERIAL= CD2
+ X= 0
+ Y= 0
+ Z= 0
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%1
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%CATSDetector
+% X1_Y1=  35.56 -35.56  -2658 mm
+% X28_Y28= -35.56  35.56  -2658 mm
+% X1_Y28= 35.56  35.56  -2658 mm
+% X28_Y1=  -35.56 -35.56  -2658 mm
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%CATSDetector
+% X28_Y1=  -35.56  35.56   -2045 mm
+% X1_Y28= 35.56 -35.56   -2045 mm
+% X28_Y28= -35.56 -35.56   -2045 mm
+% X1_Y1=  35.56  35.56   -2045 mm
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%%%%%%% Telescope 1 %%%%%%%
+M2Telescope
+ X128_Y128=  115.88  9.61  154.54 mm
+ X128_Y1=  104.8  101.89  125.09 mm 
+ X1_Y1=  14.55  102.4  160.63 mm
+ X1_Y128=  25.63  10.12  190.08 mm
+ SI=  1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 2 %%%%%%%
+M2Telescope
+ X128_Y128=  -11.23  102.42  160.87 mm
+ X128_Y1=  -101.39  102.39  124.37 mm
+ X1_Y1=  -113.17  10.36  153.56 mm
+ X1_Y128=  -23.03  10.38  190.05 mm
+ SI=  1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 3 %%%%%%%
+M2Telescope
+ X128_Y128=  -113.28  -12.52  153.32 mm
+ X128_Y1=  -101.58  -104.77  124.82 mm
+ X1_Y1=  -11.39  -104.58  161.48 mm
+ X1_Y128=  -23.1  -12.34  189.98 mm
+ SI= 1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 4 %%%%%%%
+M2Telescope
+ X128_Y128=  13.82  -104.92  160.72 mm
+ X128_Y1=  104.3  -104.95  125.08 mm
+ X1_Y1=  115.75  -12.73  153.76 mm
+ X1_Y128=  25.23  -12.65  189.43 mm
+ SI= 1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 5 %%%%%%%
+M2Telescope
+ X128_Y128= 107.49 -95.88 0.0 mm
+ X1_Y128=  143.81 -8.21 0.0 mm
+ X1_Y1= 143.81 -8.21 91.7  mm
+ X128_Y1= 107.49 -95.88 91.7 mm
+ SI= 1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid 
+ DetectorNumber= 7
+ X128_Y128=  40.99  21.69  -98.21 mm
+ X1_Y128=  23.34  39.63  -98.04 mm
+ X1_Y1=  57.53  120.53  -32.58 mm
+ X128_Y1=  122.2  55.73  -32.67 mm
+
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid 
+ DetectorNumber= 3
+ X128_Y128=  -42.66  -21.68  -99.33 mm
+ X1_Y128=  -24.83  -39.35  -99.33 mm
+ X1_Y1=  -58.24  -119.46  -32.57 mm
+ X128_Y1=  -122.86  -54.55  -32.41 mm
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 1
+ X128_Y128=  -23.05  42.83  -99.66 mm
+ X1_Y128=  -40.74  25.21  -99.65 mm
+ X1_Y1=  -120.65  58.67  -32.39 mm
+ X128_Y1=  -55.66  123.17  -32.55 mm
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 5
+ X128_Y128=  21.68  -41.48  -100.19 mm
+ X1_Y128=  39.55  -23.81  -100.11 mm
+ X1_Y1=  119.55  -57.08  -33.23 mm
+ X128_Y1=  54.81  -121.7  -33.12 mm
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 4
+ X128_Y128=  -15.11  -44.73  -99.25 mm
+ X1_Y128=  10.11  -44.71  -99.35 mm
+ X1_Y1=  43.37  -125.07  -32.74 mm
+ X128_Y1=  -48.18  -125.07  -32.45 mm
+
+%%%%%%%%%%%%%%%%%%%%%POS TO UPDATE%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 2
+ X128_Y128=  -45.68  14.28  -98.74 mm
+ X1_Y128=  -45.68  -10.92  -98.74 mm
+ X1_Y1=  -125.65  -44.06  -31.63 mm
+ X128_Y1=  -125.65  47.42  -31.63 mm
+
+%%%%%%%%%%%%%%%%%%%%POS TO UPDATE%%%%%%%%%
+%Mugast Trapezoid      
+% DetectorNumber= 6
+% X128_Y128=  45.68  -14.28  -98.74 mm
+% X1_Y128=  45.68  10.92  -98.74 mm
+% X1_Y1=  125.65  44.06  -31.63 mm
+% X128_Y1=  125.65  -47.42  -31.63 mm
+
+%%%%%%%%%%%%%%%%%%%%POS TO UPDATE%%%%%%%%%
+%Mugast Trapezoid      
+% DetectorNumber= 8
+% X128_Y128=  -15.11  -44.73  -99.25 mm
+% X1_Y128=  10.11  -44.71  -99.35 mm
+% X1_Y1=  43.37  -125.07  -32.74 mm
+% X128_Y1=  -48.18  -125.07  -32.45 mm
+
+%Mugast Annular  
+% use Annular Micron only for the Annular S1 with micron std packaging
+%Mugast AnnularMicron 
+% DetectorNumber= 11
+% Center=  -0.3 1 -126.9 mm
+
+ModularLeaf 
+ DefaultValue= -5000
+ Leafs= GATCONF_MASTER CONFDEC_AGATA T_CATS2_HF T_MUGAST_CATS2 T_MUGAST_CATS2b E_AGATA T_FAG_CATS2 T_VAMOS_AGATA T_AGATA_CATS2 T_MUGAST_VAMOS T_MUGAST_AGATA T_VAMOS_CATS2 
+%%% Leafs= GATCONF_MASTER CONFDEC_AGATA CONFDEC_VAMOS T_CATS2_HF T_MUGAST_CATS2 T_MUGAST_CATS2b EVAMOS_AGATA E_AGATA T_FAG_CATS2 T_VAMOS_AGATA T_AGATA_CATS2 TVAMOS_MUGAST_HF T_MUGAST_VAMOS T_MUGAST_AGATA T_VAMOS_CATS2 ADC1_9 ADC1_10 ADC1_11 MW_00 MW_01 MW_02 MW_03 MW_04 MW_05 MW_06 MW_07 MW_08 MW_09 MW_10 MW_11 MW_12 MW_13 MW_14 MW_15 MW_16 MW_17 MW_18 MW_19 MW_21
+ X= EVAMOS_AGATA   
+ Y= E_AGATA       
+
diff --git a/Projects/e793s/Detector/mugast_Charlie24Mar.detector b/Projects/e793s/Detector/mugast_Charlie24Mar.detector
new file mode 100644
index 0000000000000000000000000000000000000000..6f6e67d604c8c39f9ce219bb5d86d304aeb979d0
--- /dev/null
+++ b/Projects/e793s/Detector/mugast_Charlie24Mar.detector
@@ -0,0 +1,163 @@
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% CPx: check what the actual values are for this!
+Target
+ THICKNESS= 3.2 micrometer
+ ANGLE= 0 deg
+ RADIUS=  24 mm
+ MATERIAL= CD2
+ X= 0 mm
+ Y= 0 mm
+ Z= 0 mm
+
+% CATS POSITION ARE FROM E748, so WRONG
+%% CPx: commented out the CATS detectors for now 
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%CATSDetector
+% X1_Y1= -35.86  -35.26  -2658 mm
+% X28_Y1= 35.26  -35.26  -2658 mm
+% X1_Y28= -35.86  35.86  -2658 mm
+% X28_Y28= 35.26  35.86  -2658 mm
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%CATSDetector
+% X1_Y1= 34.86    -35.86   -2045 mm
+% X28_Y1= -36.26  -35.86   -2045 mm
+% X1_Y28= 34.86    35.26   -2045 mm
+% X28_Y28= -36.26  35.26   -2045 mm
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%%%%%%% Telescope 1 %%%%%%%
+M2Telescope
+ X128_Y128=  115.88  9.61  154.54
+ X128_Y1=  104.8  101.89  125.09
+ X1_Y1=  14.55  102.4  160.63
+ X1_Y128=  25.63  10.12  190.08
+ SI=  1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 2 %%%%%%%
+M2Telescope
+ X128_Y128=  -11.23  102.42  160.87
+ X128_Y1=  -101.39  102.39  124.37
+ X1_Y1=  -113.17  10.36  153.56
+ X1_Y128=  -23.03  10.38  190.05
+ SI=  1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 3 %%%%%%%
+M2Telescope
+ X128_Y128=  -113.28  -12.52  153.32
+ X128_Y1=  -101.58  -104.77  124.82
+ X1_Y1=  -11.39  -104.58  161.48
+ X1_Y128=  -23.1  -12.34  189.98
+ SI= 1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 4 %%%%%%%
+M2Telescope
+ X128_Y128=  13.82  -104.92  160.72
+ X128_Y1=  104.3  -104.95  125.08
+ X1_Y1=  115.75  -12.73  153.76
+ X1_Y128=  25.23  -12.65  189.43
+ SI= 1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+
+%%%%%%% Telescope 5 %%%%%%%
+%% CPx: fifth MUGAST position???
+%% CPx: this is just pulled from old_mugast_e793s
+M2Telescope
+ X128_Y128= 107.49 -95.88 0.0 mm
+ X1_Y128=  143.81 -8.21 0.0 mm
+ X1_Y1= 143.81 -8.21 91.7  mm
+ X128_Y1= 107.49 -95.88 91.7 mm
+ SI= 1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid 
+ DetectorNumber= 7
+ X128_Y128=  40.99  21.69  -98.21
+ X1_Y128=  23.34  39.63  -98.04
+ X1_Y1=  57.53  120.53  -32.58
+ X128_Y1=  122.2  55.73  -32.67
+
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid 
+ DetectorNumber= 3
+ X128_Y128=  -48.05  -23.23  -95.05
+ X1_Y128=  -24.83  -39.35  -99.33
+ X1_Y1=  -58.24  -119.46  -32.57
+ X128_Y1=  -122.86  -54.55  -32.41
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 1
+ X128_Y128=  -23.05  42.83  -99.66
+ X1_Y128=  -40.74  25.21  -99.65
+ X1_Y1=  -120.65  58.67  -32.39
+ X128_Y1=  -55.66  123.17  -32.55
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 5
+ X128_Y128=  21.68  -41.48  -100.19
+ X1_Y128=  39.55  -23.81  -100.11
+ X1_Y1=  119.55  -57.08  -33.23
+ X128_Y1=  54.81  -121.7  -33.12
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 4
+ X128_Y128=  -15.11  -44.73  -99.25
+ X1_Y128=  10.11  -44.71  -99.35
+ X1_Y1=  43.37  -125.07  -32.74
+ X128_Y1=  -48.18  -125.07  -32.45
+
+%% CPx: SHOULD BE ANOTHER TRAPEZOIDAL!!
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% CPx: commented out the square detectors
+%Mugast Square      
+% DetectorNumber= 10
+% X128_Y128=  111.36  97.74  89.82
+% X1_Y128=  147.47  9.83  90
+% X1_Y1=  147.7  9.88  -1.79
+% X128_Y1=  111.63  97.7  -2.02
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%Mugast Square 
+% DetectorNumber= 9
+% X128_Y128= 107.49 -95.88 0.0
+% X1_Y128=  143.81 -8.21 0.0
+% X1_Y1= 143.81 -8.21 91.7 
+% X128_Y1= 107.49 -95.88 91.7
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% CPx: no annular in e793s
+%Mugast Annular  
+% DetectorNumber= 11
+% Center=  0 0 -125 mm
+
+ModularLeaf 
+ DefaultValue= -5000
+ Leafs= GATCONF_MASTER CONFDEC_AGATA CONFDEC_VAMOS T_CATS2_HF T_MUGAST_CATS2 T_MUGAST_CATS2b EVAMOS_AGATA E_AGATA T_FAG_CATS2 T_VAMOS_AGATA T_AGATA_CATS2 TVAMOS_MUGAST_HF T_MUGAST_VAMOS T_MUGAST_AGATA T_VAMOS_CATS2 ADC1_9 ADC1_10 ADC1_11 MW_00 MW_01 MW_02 MW_03 MW_04 MW_05 MW_06 MW_07 MW_08 MW_09 MW_10 MW_11 MW_12 MW_13 MW_14 MW_15 MW_16 MW_17 MW_18 MW_19 MW_21
+ X= EVAMOS_AGATA
+ Y= E_AGATA
+
+% Leafs= GATCONF_MASTER GATCONF_SLAVE GATCONF_SLAVE2 CONFDEC_AGATA CONFDEC_VAMOS DATATRIG_CATS1 DATATRIG_CATS2 ADC1_9 CORREL_INVAMOS ADC_CATS_0
+% X= ADC1_9              ADC1_9               ADC_CATS_0
+% Y= CORREL_INVAMOS    ADC_CATS_0           CORREL_INVAMOS
+
diff --git a/Projects/e793s/Detector/mugast_Charlie24Mar_2.detector b/Projects/e793s/Detector/mugast_Charlie24Mar_2.detector
new file mode 100755
index 0000000000000000000000000000000000000000..688b09d57ce711ce944df9943aab83324983b7ec
--- /dev/null
+++ b/Projects/e793s/Detector/mugast_Charlie24Mar_2.detector
@@ -0,0 +1,161 @@
+%%%%%%%%%%Detector%%%%%%%%%%%%%%%%%%%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%1
+Target
+ THICKNESS= 4.76 micrometer
+ ANGLE= 0 deg
+ RADIUS= 10 mm
+ MATERIAL= CD2
+ X= 0
+ Y= 0
+ Z= 0
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%1
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%CATSDetector
+% X1_Y1=  35.56 -35.56  -2658 mm
+% X28_Y28= -35.56  35.56  -2658 mm
+% X1_Y28= 35.56  35.56  -2658 mm
+% X28_Y1=  -35.56 -35.56  -2658 mm
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%CATSDetector
+% X28_Y1=  -35.56  35.56   -2045 mm
+% X1_Y28= 35.56 -35.56   -2045 mm
+% X28_Y28= -35.56 -35.56   -2045 mm
+% X1_Y1=  35.56  35.56   -2045 mm
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%%%%%%% Telescope 1 %%%%%%%
+M2Telescope
+ X128_Y128=  115.88  9.61  154.54 mm
+ X128_Y1=  104.8  101.89  125.09 mm 
+ X1_Y1=  14.55  102.4  160.63 mm
+ X1_Y128=  25.63  10.12  190.08 mm
+ SI=  1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 2 %%%%%%%
+M2Telescope
+ X128_Y128=  -11.23  102.42  160.87 mm
+ X128_Y1=  -101.39  102.39  124.37 mm
+ X1_Y1=  -113.17  10.36  153.56 mm
+ X1_Y128=  -23.03  10.38  190.05 mm
+ SI=  1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 3 %%%%%%%
+M2Telescope
+ X128_Y128=  -113.28  -12.52  153.32 mm
+ X128_Y1=  -101.58  -104.77  124.82 mm
+ X1_Y1=  -11.39  -104.58  161.48 mm
+ X1_Y128=  -23.1  -12.34  189.98 mm
+ SI= 1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 4 %%%%%%%
+M2Telescope
+ X128_Y128=  13.82  -104.92  160.72 mm
+ X128_Y1=  104.3  -104.95  125.08 mm
+ X1_Y1=  115.75  -12.73  153.76 mm
+ X1_Y128=  25.23  -12.65  189.43 mm
+ SI= 1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%% Telescope 5 %%%%%%%
+M2Telescope
+ X128_Y128= 107.49 -95.88 0.0 mm
+ X1_Y128=  143.81 -8.21 0.0 mm
+ X1_Y1= 143.81 -8.21 91.7  mm
+ X128_Y1= 107.49 -95.88 91.7 mm
+ SI= 1.00
+ SILI= 0.00
+ CSI= 1.00
+ VIS= all
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid 
+ DetectorNumber= 7
+ X128_Y128=  40.99  21.69  -98.21 mm
+ X1_Y128=  23.34  39.63  -98.04 mm
+ X1_Y1=  57.53  120.53  -32.58 mm
+ X128_Y1=  122.2  55.73  -32.67 mm
+
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid 
+ DetectorNumber= 3
+ X128_Y128=  -48.05  -23.23  -95.05 mm
+ X1_Y128=  -24.83  -39.35  -99.33 mm
+ X1_Y1=  -58.24  -119.46  -32.57 mm
+ X128_Y1=  -122.86  -54.55  -32.41 mm
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 1
+ X128_Y128=  -23.05  42.83  -99.66 mm
+ X1_Y128=  -40.74  25.21  -99.65 mm
+ X1_Y1=  -120.65  58.67  -32.39 mm
+ X128_Y1=  -55.66  123.17  -32.55 mm
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 5
+ X128_Y128=  21.68  -41.48  -100.19 mm
+ X1_Y128=  39.55  -23.81  -100.11 mm
+ X1_Y1=  119.55  -57.08  -33.23 mm
+ X128_Y1=  54.81  -121.7  -33.12 mm
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 4
+ X128_Y128=  -15.11  -44.73  -99.25 mm
+ X1_Y128=  10.11  -44.71  -99.35 mm
+ X1_Y1=  43.37  -125.07  -32.74 mm
+ X128_Y1=  -48.18  -125.07  -32.45 mm
+
+%%%%%%%%%%%%%%%%%%%%%POS TO UPDATE%%%%%%%%%
+Mugast Trapezoid      
+ DetectorNumber= 2
+ X128_Y128=  -45.68  14.28  -98.74 mm
+ X1_Y128=  -45.68  -10.92  -98.74 mm
+ X1_Y1=  -125.65  -44.06  -31.63 mm
+ X128_Y1=  -125.65  47.42  -31.63 mm
+
+%%%%%%%%%%%%%%%%%%%%POS TO UPDATE%%%%%%%%%
+%Mugast Trapezoid      
+% DetectorNumber= 6
+% X128_Y128=  45.68  -14.28  -98.74 mm
+% X1_Y128=  45.68  10.92  -98.74 mm
+% X1_Y1=  125.65  44.06  -31.63 mm
+% X128_Y1=  125.65  -47.42  -31.63 mm
+
+%%%%%%%%%%%%%%%%%%%%POS TO UPDATE%%%%%%%%%
+%Mugast Trapezoid      
+% DetectorNumber= 8
+% X128_Y128=  -15.11  -44.73  -99.25 mm
+% X1_Y128=  10.11  -44.71  -99.35 mm
+% X1_Y1=  43.37  -125.07  -32.74 mm
+% X128_Y1=  -48.18  -125.07  -32.45 mm
+
+%Mugast Annular  
+% use Annular Micron only for the Annular S1 with micron std packaging
+%Mugast AnnularMicron 
+% DetectorNumber= 11
+% Center=  -0.3 1 -126.9 mm
+
+ModularLeaf 
+ DefaultValue= -5000
+ Leafs= GATCONF_MASTER CONFDEC_AGATA T_CATS2_HF T_MUGAST_CATS2 T_MUGAST_CATS2b E_AGATA T_FAG_CATS2 T_VAMOS_AGATA T_AGATA_CATS2 T_MUGAST_VAMOS T_MUGAST_AGATA T_VAMOS_CATS2 
+%%% Leafs= GATCONF_MASTER CONFDEC_AGATA CONFDEC_VAMOS T_CATS2_HF T_MUGAST_CATS2 T_MUGAST_CATS2b EVAMOS_AGATA E_AGATA T_FAG_CATS2 T_VAMOS_AGATA T_AGATA_CATS2 TVAMOS_MUGAST_HF T_MUGAST_VAMOS T_MUGAST_AGATA T_VAMOS_CATS2 ADC1_9 ADC1_10 ADC1_11 MW_00 MW_01 MW_02 MW_03 MW_04 MW_05 MW_06 MW_07 MW_08 MW_09 MW_10 MW_11 MW_12 MW_13 MW_14 MW_15 MW_16 MW_17 MW_18 MW_19 MW_21
+ X= EVAMOS_AGATA   
+ Y= E_AGATA       
+
diff --git a/Projects/e793s/Reaction/47Kdp_0keV.reaction b/Projects/e793s/Reaction/47Kdp_0keV.reaction
new file mode 100755
index 0000000000000000000000000000000000000000..f3f1ac774900ccec92351d6ebea8bcbb9018befd
--- /dev/null
+++ b/Projects/e793s/Reaction/47Kdp_0keV.reaction
@@ -0,0 +1,30 @@
+%%%%%%%%%%%%%%%%%%%%%% E793S %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Beam
+  Particle= 47K 
+  ExcitationEnergy= 0 MeV
+  Energy= 362 MeV
+  SigmaEnergy= 0. MeV
+  SigmaThetaX= 0.0  deg
+  SigmaPhiY= 0.0    deg
+  SigmaX= 2 mm
+  SigmaY= 2 mm
+  MeanThetaX= 0 deg
+  MeanPhiY= 0 deg
+  MeanX= 0 mm
+  MeanY= 0 mm
+  %EnergyProfilePath=
+  %XThetaXProfilePath=
+  %YPhiYProfilePath=
+ 
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+TwoBodyReaction
+ Beam= 47K
+ Target= 2H
+ Light= 1H
+ Heavy= 48K
+ ExcitationEnergyLight= 0.0 MeV
+ ExcitationEnergyHeavy= 0.0 MeV
+ CrossSectionPath= flat.txt CSR
+ ShootLight= 1
+ ShootHeavy= 0
diff --git a/Projects/e793s/Reaction/47Kdp_Sim.reaction b/Projects/e793s/Reaction/47Kdp_Sim.reaction
new file mode 100755
index 0000000000000000000000000000000000000000..f3f1ac774900ccec92351d6ebea8bcbb9018befd
--- /dev/null
+++ b/Projects/e793s/Reaction/47Kdp_Sim.reaction
@@ -0,0 +1,30 @@
+%%%%%%%%%%%%%%%%%%%%%% E793S %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Beam
+  Particle= 47K 
+  ExcitationEnergy= 0 MeV
+  Energy= 362 MeV
+  SigmaEnergy= 0. MeV
+  SigmaThetaX= 0.0  deg
+  SigmaPhiY= 0.0    deg
+  SigmaX= 2 mm
+  SigmaY= 2 mm
+  MeanThetaX= 0 deg
+  MeanPhiY= 0 deg
+  MeanX= 0 mm
+  MeanY= 0 mm
+  %EnergyProfilePath=
+  %XThetaXProfilePath=
+  %YPhiYProfilePath=
+ 
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+TwoBodyReaction
+ Beam= 47K
+ Target= 2H
+ Light= 1H
+ Heavy= 48K
+ ExcitationEnergyLight= 0.0 MeV
+ ExcitationEnergyHeavy= 0.0 MeV
+ CrossSectionPath= flat.txt CSR
+ ShootLight= 1
+ ShootHeavy= 0
diff --git a/Projects/e793s/Reaction/47Kdt.reaction b/Projects/e793s/Reaction/47Kdt.reaction
new file mode 100755
index 0000000000000000000000000000000000000000..f22b82ace5945cf5c56a1777535e87a3940e0675
--- /dev/null
+++ b/Projects/e793s/Reaction/47Kdt.reaction
@@ -0,0 +1,30 @@
+%%%%%%%%%%%%%%%%%%%%%% E775S %%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Beam
+  Particle= 47K 
+  ExcitationEnergy= 0 MeV
+  Energy= 362 MeV
+  SigmaEnergy= 0. MeV
+  SigmaThetaX= 0.0  deg
+  SigmaPhiY= 0.0    deg
+  SigmaX= 2 mm
+  SigmaY= 2 mm
+  MeanThetaX= 0 deg
+  MeanPhiY= 0 deg
+  MeanX= 0 mm
+  MeanY= 0 mm
+  %EnergyProfilePath=
+  %XThetaXProfilePath=
+  %YPhiYProfilePath=
+ 
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+TwoBodyReaction
+ Beam= 47K
+ Target= 2H
+ Light= 3H
+ Heavy= 46K
+ ExcitationEnergyLight= 0.0 MeV
+ ExcitationEnergyHeavy= 0.0 MeV
+ CrossSectionPath= flat.txt CSR
+ ShootLight= 1
+ ShootHeavy= 0
diff --git a/Projects/e793s/Reaction/alpha.source b/Projects/e793s/Reaction/alpha.source
new file mode 100755
index 0000000000000000000000000000000000000000..b01a3a158a752d5b5f14d9fe656ae699f063b5ae
--- /dev/null
+++ b/Projects/e793s/Reaction/alpha.source
@@ -0,0 +1,17 @@
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%%%%%%%%%% An Isotropic Source to be used as EventGenerator %%%%%%%%%%%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%      Energy are given in MeV , Position in mm      % 
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Isotropic
+ EnergyLow= 3.5 MeV
+ EnergyHigh= 3.5 MeV
+ HalfOpenAngleMin= 0 deg
+ HalfOpenAngleMax= 180 deg
+ x0= 0 mm
+ y0= 0 mm
+ z0= 0 mm
+ Particle= alpha
+ ExcitationEnergy= 0 MeV
+
+% Supported particle type: proton, neutron, deuton, triton, He3 , alpha 
diff --git a/Projects/e793s/RunToTreat.txt b/Projects/e793s/RunToTreat.txt
new file mode 100644
index 0000000000000000000000000000000000000000..ec67499e059bacbdc96a5ddd41351eb0352de992
--- /dev/null
+++ b/Projects/e793s/RunToTreat.txt
@@ -0,0 +1,14 @@
+TTreeName 
+ TreeMaster
+RootFileName 
+ /home/charlottepaxman/StoreData/K48/run_0051/Tree_0000.root
+ /home/charlottepaxman/StoreData/K48/run_0052/Tree_0000.root
+ /home/charlottepaxman/StoreData/K48/run_0053/Tree_0000.root
+ /home/charlottepaxman/StoreData/K48/run_0054/Tree_0000.root
+ /home/charlottepaxman/StoreData/K48/run_0055/Tree_0000.root
+ /home/charlottepaxman/StoreData/K48/run_0056/Tree_0000.root
+ /home/charlottepaxman/StoreData/K48/run_0060/Tree_0000.root
+ /home/charlottepaxman/StoreData/K48/run_0061/Tree_0000.root
+ /home/charlottepaxman/StoreData/K48/run_0062/Tree_0000.root
+ /home/charlottepaxman/StoreData/K48/run_0063/Tree_0000.root
+ /home/charlottepaxman/StoreData/K48/run_0064/Tree_0000.root
diff --git a/Projects/e793s/configs/ConfigCATS.dat b/Projects/e793s/configs/ConfigCATS.dat
new file mode 100755
index 0000000000000000000000000000000000000000..65e5e0ff6e2bb19317b51210d5962c80ca52821c
--- /dev/null
+++ b/Projects/e793s/configs/ConfigCATS.dat
@@ -0,0 +1,44 @@
+ConfigCATS
+
+ZPROJ 0
+ZPROJ -2190
+
+% RECONSTRUCTION METHOD
+% 5 recosntructions are available :
+% FSECH : for fitted hyperbolic secante method 
+% ASECH : for analytic hyperbolic secante method
+% AGAUSS : for analytic Gaussian method
+% CENTROIDE : for a simple centroide method
+ RECONSTRUCTION_METHOD CATS1X ASECH
+ RECONSTRUCTION_METHOD CATS1Y ASECH
+ RECONSTRUCTION_METHOD CATS2X ASECH
+ RECONSTRUCTION_METHOD CATS2Y ASECH
+ SHIFT_CATS CATS1X 0
+ SHIFT_CATS CATS1Y 0
+ SHIFT_CATS CATS2X -6.03
+ SHIFT_CATS CATS2Y 0
+
+
+
+
+% INVERSION CATS2 11 and 12 inversed   
+%DISABLE_CHANNEL CATS2STRX27
+
+
+
+% CORRECTION (Correction3Points or Correction4Points : depends on the reconstruction method)
+% use Correction3Points for GAUSS, SECH or BAR3
+% use Correction4Points for BAR4
+%        CORRECTION_METHOD CATS1X Correction4Points
+%        CORRECTION_METHOD CATS1Y Correction4Points
+%        CORRECTION_METHOD CATS2X Correction4Points
+%        CORRECTION_METHOD CATS2Y Correction4Points
+
+% CORRECTION COEF
+% See TCATSPhyscis.cxx to see the correction formula where the coefficient of correction is applied. 
+ %       CORRECTION_COEF CATS1X 0.60
+ %       CORRECTION_COEF CATS1Y 0.60
+ %       CORRECTION_COEF CATS2X 0.60
+ %       CORRECTION_COEF CATS2Y 0.60
+
+
diff --git a/Projects/e793s/configs/ConfigMugast.dat b/Projects/e793s/configs/ConfigMugast.dat
new file mode 100755
index 0000000000000000000000000000000000000000..db0c1889238d943e901981ee9eaea2fe589a6372
--- /dev/null
+++ b/Projects/e793s/configs/ConfigMugast.dat
@@ -0,0 +1,112 @@
+ConfigMugast
+ TAKE_E_Y= 1 
+ TAKE_T_Y= 1 
+
+#Dead X strips
+ DISABLE_CHANNEL_X= 2 4
+
+ DISABLE_CHANNEL_X= 3 7
+ DISABLE_CHANNEL_X= 3 5
+ 
+ DISABLE_CHANNEL_X= 4 31
+
+ DISABLE_CHANNEL_X= 5 61
+ DISABLE_CHANNEL_X= 5 24
+ DISABLE_CHANNEL_X= 5 32
+ DISABLE_CHANNEL_X= 5 125
+ DISABLE_CHANNEL_X= 5 124
+ DISABLE_CHANNEL_X= 5 109
+ DISABLE_CHANNEL_X= 5 106
+ DISABLE_CHANNEL_X= 5 104
+ DISABLE_CHANNEL_X= 5 94
+ DISABLE_CHANNEL_X= 5 76
+ DISABLE_CHANNEL_X= 5 73
+ DISABLE_CHANNEL_X= 5 71
+ DISABLE_CHANNEL_X= 5 72
+
+ DISABLE_CHANNEL_X= 7 9
+ DISABLE_CHANNEL_X= 7 77
+ DISABLE_CHANNEL_X= 7 75
+ DISABLE_CHANNEL_X= 7 73
+ DISABLE_CHANNEL_X= 7 74
+ DISABLE_CHANNEL_X= 7 71
+ DISABLE_CHANNEL_X= 7 72
+ DISABLE_CHANNEL_X= 7 69
+ DISABLE_CHANNEL_X= 7 70
+ DISABLE_CHANNEL_X= 7 67
+ DISABLE_CHANNEL_X= 7 68
+ DISABLE_CHANNEL_X= 7 65
+
+#Dead Y strips
+ DISABLE_CHANNEL_Y= 1 85
+ DISABLE_CHANNEL_Y= 1 31
+ DISABLE_CHANNEL_Y= 1 36
+ DISABLE_CHANNEL_Y= 1 27
+ DISABLE_CHANNEL_Y= 1 35
+ DISABLE_CHANNEL_Y= 1 24
+ DISABLE_CHANNEL_Y= 1 33
+ DISABLE_CHANNEL_Y= 1 39
+ DISABLE_CHANNEL_Y= 1 37
+ DISABLE_CHANNEL_Y= 1 22
+ DISABLE_CHANNEL_Y= 1 46
+ DISABLE_CHANNEL_Y= 1 26
+ DISABLE_CHANNEL_Y= 1 13
+ DISABLE_CHANNEL_Y= 1 15
+
+ DISABLE_CHANNEL_Y= 2 39
+ DISABLE_CHANNEL_Y= 2 16
+ DISABLE_CHANNEL_Y= 2 62
+ 
+ DISABLE_CHANNEL_Y= 3 11
+
+ DISABLE_CHANNEL_Y= 4 100
+ DISABLE_CHANNEL_Y= 4 90
+ DISABLE_CHANNEL_Y= 4 102
+ DISABLE_CHANNEL_Y= 4 117
+ DISABLE_CHANNEL_Y= 4 105
+ DISABLE_CHANNEL_Y= 4 38
+ DISABLE_CHANNEL_Y= 4 36
+ DISABLE_CHANNEL_Y= 4 22
+ DISABLE_CHANNEL_Y= 4 50
+ DISABLE_CHANNEL_Y= 4 26
+ DISABLE_CHANNEL_Y= 4 23
+ DISABLE_CHANNEL_Y= 4 64
+ DISABLE_CHANNEL_Y= 4 4
+ DISABLE_CHANNEL_Y= 4 2
+ DISABLE_CHANNEL_Y= 4 55
+
+ DISABLE_CHANNEL_Y= 5 120
+ DISABLE_CHANNEL_Y= 5 122
+ DISABLE_CHANNEL_Y= 5 113
+ DISABLE_CHANNEL_Y= 5 111
+ DISABLE_CHANNEL_Y= 5 109
+ DISABLE_CHANNEL_Y= 5 93
+ DISABLE_CHANNEL_Y= 5 107
+ DISABLE_CHANNEL_Y= 5 105
+ DISABLE_CHANNEL_Y= 5 38
+ DISABLE_CHANNEL_Y= 5 36
+ DISABLE_CHANNEL_Y= 5 4
+
+ DISABLE_CHANNEL_Y= 7 126
+ DISABLE_CHANNEL_Y= 7 103
+ DISABLE_CHANNEL_Y= 7 128
+ DISABLE_CHANNEL_Y= 7 123
+ DISABLE_CHANNEL_Y= 7 121
+ DISABLE_CHANNEL_Y= 7 105
+ DISABLE_CHANNEL_Y= 7 38
+ DISABLE_CHANNEL_Y= 7 36
+ DISABLE_CHANNEL_Y= 7 35
+ DISABLE_CHANNEL_Y= 7 40
+ DISABLE_CHANNEL_Y= 7 43
+ DISABLE_CHANNEL_Y= 7 46
+ DISABLE_CHANNEL_Y= 7 54
+ DISABLE_CHANNEL_Y= 7 60
+ DISABLE_CHANNEL_Y= 7 63
+ DISABLE_CHANNEL_Y= 7 57
+
+ MAX_STRIP_MULTIPLICITY= 10
+ STRIP_ENERGY_MATCHING= 0.2 MeV
+ DSSD_X_E_RAW_THRESHOLD= 8250
+ DSSD_Y_E_RAW_THRESHOLD= 8100
+ DSSD_X_E_THRESHOLD= 1 MeV
+ DSSD_Y_E_THRESHOLD= 1 MeV
diff --git a/Projects/e793s/configs/ConfigMust2.dat b/Projects/e793s/configs/ConfigMust2.dat
new file mode 100755
index 0000000000000000000000000000000000000000..e117fa5f8a896ed8e1b0d6b0eec551dde21f2d33
--- /dev/null
+++ b/Projects/e793s/configs/ConfigMust2.dat
@@ -0,0 +1,45 @@
+ConfigMust2
+ MAX_STRIP_MULTIPLICITY 100
+ STRIP_ENERGY_MATCHING_NUMBER_OF_SIGMA 1
+ STRIP_ENERGY_MATCHING_SIGMA 0.2
+ % Disabling dead strips
+ DISABLE_CHANNEL MM1STRX3
+ DISABLE_CHANNEL MM1STRX48
+
+ DISABLE_CHANNEL MM1STRY105
+ DISABLE_CHANNEL MM1STRY119
+ DISABLE_CHANNEL MM1STRY123
+ DISABLE_CHANNEL MM1STRY126
+ 
+ DISABLE_CHANNEL MM2STRX7
+ DISABLE_CHANNEL MM2STRX12
+ DISABLE_CHANNEL MM2STRX44
+ DISABLE_CHANNEL MM2STRX68
+ DISABLE_CHANNEL MM2STRX111
+ DISABLE_CHANNEL MM2STRX121
+ DISABLE_CHANNEL MM2STRX126
+ 
+ DISABLE_CHANNEL MM2STRY1
+ DISABLE_CHANNEL MM2STRY2
+ DISABLE_CHANNEL MM2STRY3
+ DISABLE_CHANNEL MM2STRY4
+ DISABLE_CHANNEL MM2STRY19
+ DISABLE_CHANNEL MM2STRY20
+ DISABLE_CHANNEL MM2STRY24
+ DISABLE_CHANNEL MM2STRY25
+ DISABLE_CHANNEL MM2STRY26
+ DISABLE_CHANNEL MM2STRY108
+ DISABLE_CHANNEL MM2STRY108
+ 
+%keep filling this in........
+
+
+
+ SI_X_E_RAW_THRESHOLD 8200
+ SI_Y_E_RAW_THRESHOLD 8180
+ SI_X_E_THRESHOLD 1
+ SI_Y_E_THRESHOLD 1
+ CSI_E_RAW_THRESHOLD 8250
+ CSI_SIZE 34
+ TAKE_T_X
+ TAKE_E_Y
diff --git a/Projects/e793s/configs/OriginalConfigMugast.dat b/Projects/e793s/configs/OriginalConfigMugast.dat
new file mode 100755
index 0000000000000000000000000000000000000000..fea5c151c27a3893f4384ebf15ede62a8793b754
--- /dev/null
+++ b/Projects/e793s/configs/OriginalConfigMugast.dat
@@ -0,0 +1,14 @@
+ConfigMugast
+ TAKE_E_Y= 1 
+ TAKE_T_Y= 1 
+ %raw channel 86 corresponds to stripY 46
+ %DISABLE_CHANNEL_Y= 3 123
+ %DISABLE_CHANNEL_X= 7 63
+ %DISABLE_CHANNEL_X= 7 64
+ %DISABLE_CHANNEL_X= 7 65
+ MAX_STRIP_MULTIPLICITY= 10
+ STRIP_ENERGY_MATCHING= 0.2 MeV
+ DSSD_X_E_RAW_THRESHOLD= 8250
+ DSSD_Y_E_RAW_THRESHOLD= 8100
+ DSSD_X_E_THRESHOLD= 1 MeV
+ DSSD_Y_E_THRESHOLD= 1 MeV
diff --git a/Projects/e793s/configs/OriginalConfigMust2.dat b/Projects/e793s/configs/OriginalConfigMust2.dat
new file mode 100755
index 0000000000000000000000000000000000000000..30b74ed6d0f2c75f0ee7878138462babc81aadcb
--- /dev/null
+++ b/Projects/e793s/configs/OriginalConfigMust2.dat
@@ -0,0 +1,21 @@
+ConfigMust2
+ MAX_STRIP_MULTIPLICITY 100
+ STRIP_ENERGY_MATCHING_NUMBER_OF_SIGMA 1
+ STRIP_ENERGY_MATCHING_SIGMA 0.2
+ %DISABLE_CHANNEL MM1STRX12
+ %DISABLE_CHANNEL MM1STRX124
+ %DISABLE_CHANNEL MM2STRX12
+ %DISABLE_CHANNEL MM2STRX123
+ %DISABLE_CHANNEL MM2STRY123
+ %DISABLE_CHANNEL MM2STRX127
+ %DISABLE_CHANNEL MM3STRX12
+ %DISABLE_CHANNEL MM3STRX7
+ %DISABLE_CHANNEL MM4STRX12
+ SI_X_E_RAW_THRESHOLD 8200
+ SI_Y_E_RAW_THRESHOLD 8180
+ SI_X_E_THRESHOLD 1
+ SI_Y_E_THRESHOLD 1
+ CSI_E_RAW_THRESHOLD 8250
+ CSI_SIZE 34
+ TAKE_T_X
+ TAKE_E_Y
diff --git a/Projects/e793s/macro/BeamSpot/BeamSpot.C b/Projects/e793s/macro/BeamSpot/BeamSpot.C
new file mode 100755
index 0000000000000000000000000000000000000000..a53f6575e34b5b1405404f9a22f689af4e97c4db
--- /dev/null
+++ b/Projects/e793s/macro/BeamSpot/BeamSpot.C
@@ -0,0 +1,381 @@
+/***************************************************************************
+ *   Charlie Paxman (cp00474@surrey.ac.uk)                                 *
+ ***************************************************************************/
+
+//-------------------------------
+//C++
+#include <iostream>
+#include <fstream>
+#include <vector>
+#include <cmath>
+//Root
+#include <TVector3.h>
+//NPTool
+#include "NPEnergyLoss.h"
+#include "NPReaction.h"
+#include "NPSystemOfUnits.h"
+using namespace std;
+//-------------------------------
+
+// TO DO:
+// - Implement ROOT minimiser?
+// - Too much repeptitive code! clean up, chuck in functions
+// - Push project?
+
+void BeamSpot(){
+  vector <double> Xp, Yp, Zp, Ep;  //XYZE of detected particle. Constant.
+  vector <int>    DetNum;          //Telescope number of hit
+  vector <double> thetaLab, phiLab;
+  double ELab = 0.0, Ex = 0.0;
+  vector <double> Xd, Yd, Zd;      //Vector of particle direction. Calculated as Xp-Xb, Yp-Yb...
+  ifstream MugastDataFile;
+
+  /*** ITERATIVE GRID CONTROLS ***/
+  /***** pos varied as offset ****/
+  /**/ double xmin = +0.000;   /**/
+  /**/ double xmax = +0.100;   /**/
+  /**/ unsigned int xdiv = 10; /**/
+  /**/                         /**/
+  /**/ double ymin = -0.100;   /**/
+  /**/ double ymax = +0.100;   /**/
+  /**/ unsigned int ydiv = 10; /**/
+  /**/                         /**/
+  /**/ double zmin = -0.000;   /**/
+  /**/ double zmax = +0.000;   /**/
+  /**/ unsigned int zdiv =  1; /**/
+  /**/                         /**/
+  /***** thick varied as %ge *****/
+  /**/ unsigned int tmin = 7;  /**/
+  /**/ unsigned int tmax = 13; /**/
+  /**/ double tmult = 0.1;     /**/
+  /*******************************/
+
+  // Calculate size of iteratve steps
+  double xstp = (xmax-xmin)/ ((double) xdiv);
+  double ystp = (ymax-ymin)/ ((double) ydiv);
+  double zstp = (zmax-zmin)/ ((double) zdiv);
+  cout << "Xstp = " << xstp << "   Ystp = " << ystp << "   Zstp = " << zstp << endl;
+
+  // Vectors of the normal for each detector. UPDATE WITH NEW POSITIONS!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!1
+  TVector3 MugastNormal1{ -0.453915, +0.455463, -0.765842};
+  TVector3 MugastNormal2{ -0.642828, +0.000000, -0.766010};
+  TVector3 MugastNormal3{ -0.454594, -0.450670, -0.768271};
+  TVector3 MugastNormal4{ -0.002437, -0.638751, -0.769409};
+  TVector3 MugastNormal5{ +0.452429, -0.454575, -0.767248};
+  TVector3 MugastNormal7{ +0.443072, +0.443265, -0.779232};
+
+  // Read in NPTool values
+  NPL::EnergyLoss LightAl;
+  NPL::EnergyLoss LightTarget;
+  LightAl = NPL::EnergyLoss("proton_Al.G4table", "G4Table", 100);
+  LightTarget = NPL::EnergyLoss("proton_CD2.G4table", "G4Table", 100);
+  NPL::Reaction reaction;
+  reaction.ReadConfigurationFile("../../Reaction/47Kdp_0keV.reaction");
+
+  // Open and read event data file
+  //MugastDataFile.open("XYZE_gammaGated_Full.txt", ios::in);
+  MugastDataFile.open("XYZE_gammaGated_Run63.txt", ios::in);
+  //MugastDataFile.open("XYZE_writeRun62_May04.txt", ios::in);
+  //MugastDataFile.open("XYZE_writeRun63_May05.txt", ios::in);
+  if(!MugastDataFile){
+    cout << "ERROR: File not opened." << endl;
+  }
+  else{
+    int numTemp = 0;
+    double Xptemp = 0.0, Yptemp = 0.0, Zptemp = 0.0, Eptemp = 0.0; 
+    while(true){ // Exit through break command
+      MugastDataFile >> numTemp >> Xptemp >> Yptemp >> Zptemp >> Eptemp;
+      DetNum.push_back(numTemp);
+      Xp.push_back(Xptemp);
+      Yp.push_back(Yptemp);
+      Zp.push_back(Zptemp);
+      Ep.push_back(Eptemp);
+      
+      Xptemp = Yptemp = Zptemp = Eptemp = 0.0;
+      numTemp = 0;
+
+      if(MugastDataFile.eof()) {break;}
+    }
+  }
+
+  // Initilize various values and histograms
+  int size = Xp.size();
+  double tempTheta;
+  TH1F *tempHist = new TH1F("tempHist","All Detectors", 40,-1.0,1.0);
+  TH1F *tempMG1 = new TH1F("tempMG1","Individual MG#", 40,-1.0,1.0);
+  TH1F *tempMG2 = new TH1F("tempMG2","Individual MG#", 40,-1.0,1.0);
+  TH1F *tempMG3 = new TH1F("tempMG3","Individual MG#", 40,-1.0,1.0);
+  TH1F *tempMG4 = new TH1F("tempMG4","Individual MG#", 40,-1.0,1.0);
+  TH1F *tempMG5 = new TH1F("tempMG5","Individual MG#", 40,-1.0,1.0);
+  TH1F *tempMG7 = new TH1F("tempMG7","Individual MG#", 40,-1.0,1.0);
+
+  // ------------ ITERATE BEAM POSITION ------------ //
+  for(int thickiter=tmin; thickiter<tmax+1; thickiter++){
+    double TargetThickness = (thickiter * tmult) * 0.00476; //multiplier * original thickness
+    //cout << " Thickness @ " << TargetThickness << endl;
+    for(int xiter=0; xiter<xdiv+1; xiter++){
+      //cout << " X @ " << xmin+(xiter*xstp) << endl;
+      for(int yiter=0; yiter<ydiv+1; yiter++){
+        //cout << " Y @ " << ymin+(yiter*ystp) << endl;
+        for(int ziter=0; ziter<zdiv+1; ziter++){
+          //cout << " Z @ " << zmin+(ziter*zstp) << endl;
+
+  	  // Build beam spot vector for this iteration 
+          TVector3 beamSpot{ xmin+(xiter*xstp), ymin+(yiter*ystp), zmin+(ziter*zstp) }; 
+
+	  // Clean histograms
+	  tempHist->Reset();
+	  tempMG1->Reset();
+	  tempMG2->Reset();
+	  tempMG3->Reset();
+	  tempMG4->Reset();
+	  tempMG5->Reset();
+	  tempMG7->Reset();
+
+          // -------- ITERATE PARTICLE NUM -------- //
+          for(int i=0; i<size-1; i++){
+            // Build particle direction vector
+	    TVector3 particleDir{ Xp[i] - beamSpot.X(),   //Xd
+  	 	                  Yp[i] - beamSpot.Y(),   //Yd
+			          Zp[i] - beamSpot.Z() }; //Zd
+	  
+	    tempTheta = ELab = Ex = 0.0;
+
+	    switch(DetNum[i]){
+	      case 1:
+                tempTheta = particleDir.Angle(MugastNormal1);
+	        break;
+	      case 2:
+                tempTheta = particleDir.Angle(MugastNormal2);
+	        break;
+	      case 3:
+                tempTheta = particleDir.Angle(MugastNormal3);
+	        break;
+  	      case 4:
+                tempTheta = particleDir.Angle(MugastNormal4);
+	        break;
+	      case 5:
+                tempTheta = particleDir.Angle(MugastNormal5);
+	        break;
+	      case 7:
+                tempTheta = particleDir.Angle(MugastNormal7);
+	        break;
+	      default:
+	        cout << "ERROR! Invalid DetNum: " << DetNum[i] << " @" << i << endl;
+	        return; // Exit code
+	    }
+
+	    //micrometer defined in NPSystemOfUnits.h
+	    ELab = LightAl.EvaluateInitialEnergy(Ep[i], 0.4*micrometer, tempTheta); 
+	    ELab = LightTarget.EvaluateInitialEnergy(ELab, 0.5*TargetThickness, 0.);
+
+            // Change beam spot vector to beam direction vector
+	    beamSpot.SetZ(1.0);
+            Ex = reaction.ReconstructRelativistic( ELab, particleDir.Angle(beamSpot) );
+
+	    // Fill Ex histograms
+	    tempHist->Fill(Ex);
+	    switch(DetNum[i]){
+	      case 1:
+                tempMG1->Fill(Ex); 
+	        break;
+	      case 2:
+                tempMG2->Fill(Ex); 
+	        break;
+	      case 3:
+                tempMG3->Fill(Ex); 
+	        break;
+  	      case 4:
+                tempMG4->Fill(Ex); 
+	        break;
+	      case 5:
+                tempMG5->Fill(Ex); 
+	        break;
+	      case 7:
+                tempMG7->Fill(Ex); 
+	        break;
+	      default:
+	        cout << "ERROR! Invalid DetNum: " << DetNum[i] << " @" << i << endl;
+	        return; // Exit code
+	    }
+          }
+	  // ------ END OF PARTICLE ITERATION ----- //
+          // -------------------------------------- // 
+
+	  // Initilise gaussian variable arrays
+          double mean[7];
+          double meanErr[7];
+          double sigma[7];
+          double sigmaErr[7];
+
+	  // Draw Ex histograms
+          TCanvas *c1 = new TCanvas("c1","Ex Histograms",20,20,1600,800);
+          gStyle->SetOptStat(0);
+          c1->Divide(2,1,0.005,0.005,0);
+          c1->cd(1)->SetLeftMargin(0.15);
+          c1->cd(1)->SetBottomMargin(0.15);
+          gPad->SetTickx();
+          gPad->SetTicky();
+          c1->cd(2)->SetLeftMargin(0.15);
+          c1->cd(2)->SetBottomMargin(0.15);
+          gPad->SetTickx();
+          gPad->SetTicky();
+      
+          c1->cd(1);
+          tempMG1->SetMaximum(75.);
+          tempMG1->GetXaxis()->SetTitle("Ex [MeV]");
+          tempMG1->GetYaxis()->SetTitle("Counts");
+      
+          // ----- MG1 -----
+          tempMG1->SetLineColor(kRed);
+          tempMG1->SetFillStyle(3244);
+          tempMG1->SetFillColor(kRed);
+          tempMG1->Draw();
+      
+          tempMG1->Fit("gaus","WQ"); //add N to stop it drawing
+          TF1 *gaus = (TF1*)tempMG1->GetListOfFunctions()->FindObject("gaus");
+          mean[1] = gaus->GetParameter(1);
+          meanErr[1] = gaus->GetParError(1);
+          sigma[1] = gaus->GetParameter(2);
+          sigmaErr[1] = gaus->GetParError(2);
+          // ---------------
+          // ----- MG2 -----
+          tempMG2->SetLineColor(kOrange);
+          tempMG2->SetFillStyle(3244);
+          tempMG2->SetFillColor(kOrange);
+          tempMG2->Draw("same");
+
+          tempMG2->Fit("gaus","WQ"); //add N to stop it drawing
+          gaus = (TF1*)tempMG2->GetListOfFunctions()->FindObject("gaus");
+          mean[2] = gaus->GetParameter(1);
+          meanErr[2] = gaus->GetParError(1);
+          sigma[2] = gaus->GetParameter(2);
+          sigmaErr[2] = gaus->GetParError(2);
+          // ---------------
+          // ----- MG3 -----
+          tempMG3->SetLineColor(kGreen);
+          tempMG3->SetFillStyle(3344);
+          tempMG3->SetFillColor(kGreen);
+          tempMG3->Draw("same");
+
+          tempMG3->Fit("gaus","WQ"); //add N to stop it drawing
+          gaus = (TF1*)tempMG3->GetListOfFunctions()->FindObject("gaus");
+          mean[3] = gaus->GetParameter(1);
+          meanErr[3] = gaus->GetParError(1);
+          sigma[3] = gaus->GetParameter(2);
+          sigmaErr[3] = gaus->GetParError(2);
+          // ---------------
+          // ----- MG4 -----
+          tempMG4->SetLineColor(kTeal);
+          tempMG4->SetFillStyle(3444);
+          tempMG4->SetFillColor(kTeal);
+          tempMG4->Draw("same");
+
+	  tempMG4->Fit("gaus","WQ"); //add N to stop it drawing
+          gaus = (TF1*)tempMG4->GetListOfFunctions()->FindObject("gaus");
+          mean[4] = gaus->GetParameter(1);
+          meanErr[4] = gaus->GetParError(1);
+          sigma[4] = gaus->GetParameter(2);
+          sigmaErr[4] = gaus->GetParError(2);
+          // ---------------
+          // ----- MG5 -----
+	  tempMG5->SetLineColor(kBlue);
+          tempMG5->SetFillStyle(3544);
+          tempMG5->SetFillColor(kBlue);
+          tempMG5->Draw("same");
+
+          tempMG5->Fit("gaus","WQ"); //add N to stop it drawing
+          gaus = (TF1*)tempMG5->GetListOfFunctions()->FindObject("gaus");
+          mean[5] = gaus->GetParameter(1);
+          meanErr[5] = gaus->GetParError(1);
+          sigma[5] = gaus->GetParameter(2);
+          sigmaErr[5] = gaus->GetParError(2);
+          // ---------------
+	  // ----- MG7 -----
+	  tempMG7->SetLineColor(kViolet);
+          tempMG7->SetFillStyle(3644);
+          tempMG7->SetFillColor(kViolet);
+	  tempMG7->Draw("same");
+
+	  tempMG7->Fit("gaus","WQ"); //add N to stop it drawing
+	  gaus = (TF1*)tempMG7->GetListOfFunctions()->FindObject("gaus");
+          mean[6] = gaus->GetParameter(1);
+          meanErr[6] = gaus->GetParError(1);
+          sigma[6] = gaus->GetParameter(2);
+          sigmaErr[6] = gaus->GetParError(2);
+          // ---------------
+          // Format legend
+	  auto legend = new TLegend(0.15,0.7,0.35,0.9);
+          legend->AddEntry(tempMG1,"MUGAST 1","f");
+          legend->AddEntry(tempMG2,"MUGAST 2","f");
+          legend->AddEntry(tempMG3,"MUGAST 3","f");
+          legend->AddEntry(tempMG4,"MUGAST 4","f");
+          legend->AddEntry(tempMG5,"MUGAST 5","f");
+          legend->AddEntry(tempMG7,"MUGAST 7","f");
+          legend->Draw();
+          // ----- ALL -----
+	  c1->cd(2);
+          tempHist->GetXaxis()->SetTitle("Ex [MeV]");
+          tempHist->GetYaxis()->SetTitle("Counts");
+          tempHist->Draw();
+          tempHist->Fit("gaus", "WQ");
+
+	  gaus = (TF1*)tempHist->GetListOfFunctions()->FindObject("gaus");
+          mean[0] = gaus->GetParameter(1);
+          meanErr[0] = gaus->GetParError(1);
+          sigma[0] = gaus->GetParameter(2);
+          sigmaErr[0] = gaus->GetParError(2);
+          // ---------------
+
+	  // Caluclate metric
+          double metric[7];
+          for (int i=0; i<7; i++){
+            metric[i] = pow(mean[i]-0.143,2) + pow(0.5 * sigma[i],2);
+          }
+
+	  // Write values to screen
+          //cout << TargetThickness << "\t"
+	  //     << xmin+(xiter*xstp) << "\t"
+	  //     << ymin+(yiter*ystp) << "\t"
+	  //     << zmin+(ziter*zstp) << "\t\t"
+	  //     << metric[0] << "\t"
+          //     << metric[1] << "\t"
+          //     << metric[2] << "\t"
+          //     << metric[3] << "\t"
+          //     << metric[4] << "\t"
+          //     << metric[5] << "\t"
+          //     << metric[6] << endl;
+
+	  cout << TargetThickness << "\t"
+	       << xmin+(xiter*xstp) << "\t"
+	       << ymin+(yiter*ystp) << "\t"
+	       << zmin+(ziter*zstp) << "\t\t"
+	       << mean[0] << "\t" << sigma[0] << "\t"
+	       << mean[1] << "\t" << sigma[1] << "\t"
+	       << mean[2] << "\t" << sigma[2] << "\t"
+	       << mean[3] << "\t" << sigma[3] << "\t"
+	       << mean[4] << "\t" << sigma[4] << "\t"
+	       << mean[5] << "\t" << sigma[5] << "\t"
+	       << mean[6] << "\t" << sigma[6] << "\t"
+               << endl;
+
+
+	  // Write histograms to file
+          string fileOut = "./Histograms/targetIter";
+          fileOut += to_string(TargetThickness);
+          fileOut += "_x";
+          fileOut += to_string(xmin+(xiter*xstp));
+          fileOut += "_y";
+          fileOut += to_string(ymin+(yiter*ystp));
+          fileOut += "_z";
+          fileOut += to_string(zmin+(ziter*zstp));
+          fileOut += ".pdf";
+
+          c1->SaveAs(fileOut.c_str());
+        } //iterate Z
+      } //iterate Y
+    } //iterate X
+  } //iterate thickness
+  //cout << " --- COMPLETE --- " << endl;
+  // ----------------------------------------------- //
+}
diff --git a/Projects/e793s/macro/BeamSpot/EventReader.C b/Projects/e793s/macro/BeamSpot/EventReader.C
new file mode 100755
index 0000000000000000000000000000000000000000..ff1ec0d0fb2d2dab29e399d17321346a555f8b81
--- /dev/null
+++ b/Projects/e793s/macro/BeamSpot/EventReader.C
@@ -0,0 +1,95 @@
+/***************************************************************************
+ *   Charlie Paxman (cp00474@surrey.ac.uk)                                 *
+ ***************************************************************************/
+
+//-------------------------------
+//C++
+#include <iostream>
+#include <fstream>
+#include <vector>
+#include <cmath>
+//Root
+//#include <TVector3.h>
+//NPTool
+//#include "NPEnergyLoss.h"
+//#include "NPReaction.h"
+//#include "NPSystemOfUnits.h"
+using namespace std;
+//-------------------------------
+
+void EventReader(){
+
+  // Read ROOT file and pull the tree	
+  //auto DataFile = new TFile("../../../../Outputs/Analysis/47K_RawEnergyBranch_Run63_May11.root", "READ");
+  auto DataFile = new TFile("../../../../Outputs/Analysis/47K_RawEnergyBranch_Full_May11.root", "READ");
+  auto PhysTree = (TTree*) DataFile->FindObjectAny("PhysicsTree");
+
+  // Initilise the Mugast branch
+  auto Mugast = new TMugastPhysics();
+
+  // Initilise access variables for PhysTree
+  static double T_MUGAST_VAMOS;
+  static vector<double> *X, *Y, *Z, *RawEnergy, *AddBack_EDC;
+
+  // Pull PhysTree branches
+  auto Energy_Branch = PhysTree->GetBranch("RawEnergy");
+  auto Gamma_Branch = PhysTree->GetBranch("AddBack_EDC");
+  auto X_Branch = PhysTree->GetBranch("X");
+  auto Y_Branch = PhysTree->GetBranch("Y");
+  auto Z_Branch = PhysTree->GetBranch("Z");
+  auto MugVam_Branch = PhysTree->GetBranch("T_MUGAST_VAMOS");
+
+  // Set Mugast branch address
+  PhysTree->SetBranchAddress("Mugast",&Mugast);
+
+  // Set PhysTree variable addresses
+  Energy_Branch->SetAddress(&RawEnergy);
+  Gamma_Branch->SetAddress(&AddBack_EDC);
+  X_Branch->SetAddress(&X);
+  Y_Branch->SetAddress(&Y);
+  Z_Branch->SetAddress(&Z);
+  MugVam_Branch->SetAddress(&T_MUGAST_VAMOS);
+
+  // Build loop variables
+  unsigned int numEntries = PhysTree->GetEntries();
+  unsigned int multiplicity = 0;
+
+  // Open output file
+  ofstream outfile;
+  outfile.open("./XYZE_gammaGated_Full.txt");
+
+  // Loop on entries
+  for(unsigned int i=0; i<numEntries; i++){
+    PhysTree->GetEntry(i);
+    multiplicity = Mugast->TelescopeNumber.size(); 
+
+    // Loop on MUGAST multiplicity
+    for(int m=0; m<multiplicity; m++){
+
+      // Gate on Timing
+      if(abs(T_MUGAST_VAMOS-2777)<600){
+        int gammaMultip = AddBack_EDC->size();
+//        cout << "GAMMA MULT " << gammaMultip << endl;
+
+	// gamma stuff (fails otherwise??)
+        //if(gammaMultip!=0){
+        if(gammaMultip==1){
+        //for(int g; g<gammaMultip; g++){
+	
+	  // Gate on E_Gamma
+	  if(abs(AddBack_EDC->at(0)-0.143) < 0.005){
+	    outfile << Mugast->TelescopeNumber.at(m) << " " <<
+	               X->at(m) << " " <<
+                       Y->at(m) << " " <<
+	               Z->at(m) << " " <<
+	               RawEnergy->at(m) << endl;
+          }//if 0.143
+        //}//for g
+        }//if gamma
+      }//if timing
+    }//for m
+  }//for i
+  outfile.close();
+
+  cout << " --- COMPLETE --- " << endl;
+}
diff --git a/Projects/e793s/macro/DrawPlots.C b/Projects/e793s/macro/DrawPlots.C
new file mode 100644
index 0000000000000000000000000000000000000000..d8d416d71a4366d8803e1560eefcd971db9d6b88
--- /dev/null
+++ b/Projects/e793s/macro/DrawPlots.C
@@ -0,0 +1,457 @@
+//void AddTiStates(double E);
+#include "NPReaction.h"
+
+//TCutG* ETOF=NULL;
+//TCutG* EDE=NULL;
+TChain* chain=NULL ;
+
+char cond[1000];
+
+
+////////////////////////////////////////////////////////////////////////////////
+TChain* Chain(std::string TreeName, std::vector<std::string>& file, bool EventList){
+  TChain*  chain = new TChain(TreeName.c_str());
+  unsigned int size =file.size(); 
+    for(unsigned int i = 0 ; i < size ; i++){
+      cout << "Adding " << file[i] << endl;
+      chain->Add(file[i].c_str());
+     /* TFile* file = new TFile(file[i].c_str(),"READ");
+      if(EventList){
+        TEventList* current = (TEventList*) file->Get("GoodEntries");
+        if(!EL)
+          EL=current;
+        else{
+          
+          EL->Add(current);
+        }
+      }*/
+    } 
+  //chain->SetEntryListFile("$.root/"); 
+  return chain;
+}
+////////////////////////////////////////////////////////////////////////////////
+
+void LoadChainNP(){
+
+  vector<string> files;
+  //files.push_back("../../Outputs/Analysis/47K_Part_25Mar.root");
+  //files.push_back("../../Outputs/Analysis/47K_Full_24Mar.root");
+  //files.push_back("../../Outputs/Analysis/test_09Apr_newPhiAnalysis.root");
+//  files.push_back("../../Outputs/Analysis/47K_Full_23Apr.root");
+  files.push_back("../../../Outputs/Analysis/47K_RawEnergyBranch_Run63_May11.root");
+
+  chain = Chain("PhysicsTree",files,true);
+  //  chain = new TChain("PhysicsTree");
+  //  chain->Add("NPOutput/PhyAllCom2019_256_AL.root");
+  // chain->Add("NPOutput/PhyAllCom2019_257_AL.root");
+  //chain->Add("NPOutput/PhyAllCom2019_257.root");
+}
+//////////////////////////////////////////////////////////////////////////////////
+void LoadEventList(){
+}
+
+void plot_kine(NPL::Reaction r, double Ex,Color_t c,int w, int s){
+  r.SetExcitation4(Ex);
+
+  TGraph* g= r.GetKinematicLine3();
+  g->SetLineColor(c) ;
+  g->SetLineStyle(s) ;
+  g->SetLineWidth(w) ;
+
+  g->Draw("c");
+}
+////////////////////////////////////////////////////////////////////////////////
+void plot_state(double Ex,double max,Color_t c,int w, int s){
+  TLine* line = new TLine(Ex,0,Ex,max) ; 
+  line->SetLineColor(c) ;
+  line->SetLineStyle(s) ;
+  line->SetLineWidth(w) ;
+  line->Draw();
+
+}
+
+void AddTiStates(double E){
+  
+ NPL::Reaction Ti("47Ti(d,p)48Ti@362");
+   // Ti states
+ Ti.SetExcitationHeavy(E);
+ auto g = Ti.GetKinematicLine3();
+ g->SetLineWidth(1);
+ g->SetLineStyle(2);
+ g->Draw("c");
+
+}
+
+void DrawPlots(){
+
+  gStyle->SetOptStat("nei");
+  LoadChainNP();
+  
+  NPL::Reaction Kdp("47K(d,p)48K@362");
+  NPL::Reaction Kdt("47K(d,t)46K@362");
+  NPL::Reaction Kdd("47K(d,d)47K@362");
+  NPL::Reaction Kpp("47K(p,p)47K@362");
+  NPL::Reaction K12C12C("47K(12C,12C)47K@362");
+  NPL::Reaction Tidp("47Ti(d,p)48Ti@362");
+  NPL::Reaction Tidt("47Ti(d,t)46Ti@362");
+  NPL::Reaction Tidd("47Ti(d,d)47Ti@362");
+  NPL::Reaction Ti12C12C("47Ti(12C,12C)47Ti@362");
+
+
+  /*
+  TCanvas *cF = new TCanvas("cF","cF",1000,1000);
+  gPad->SetLogz();
+  chain->Draw("ELab:ThetaLab>>hcF(130,50,180,600,0,12)","abs(T_MUGAST_VAMOS-2777)<600","col");
+  TH2F* hcF = (TH2F*) gDirectory->Get("hcF");
+  hcF->GetXaxis()->SetTitle("#theta_{lab} [deg]");
+  hcF->GetYaxis()->SetTitle("E_{lab} [MeV]");
+  plot_kine(Kdp, 0, kBlack, 2, 9);
+  plot_kine(Kdp, 2, kRed, 2, 9);
+  plot_kine(Kdp, 4, kBlue, 2, 9);
+  plot_kine(Kdd, 0, kGreen+2, 2, 9);
+  plot_kine(Kpp, 0, kYellow, 2, 9);
+  */
+
+
+  /* MUGAST testing phi dependance */
+  //TCanvas *diagnose = new TCanvas("diagnose","diagnose",1000,1000);
+//  chain->Draw("PhiLab","abs(T_MUGAST_VAMOS-2777)<600 && (Mugast.TelescopeNumber==1 ||  Mugast.TelescopeNumber==2 || Mugast.TelescopeNumber==3 || Mugast.TelescopeNumber==4 || Mugast.TelescopeNumber==5 || Mugast.TelescopeNumber==7)","");
+//  chain->Draw("PhiLab>>h1","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==1","");
+//  chain->Draw("PhiLab>>h2","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==2","same");
+//  chain->Draw("PhiLab>>h3","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==3","same");
+//  chain->Draw("PhiLab>>h4","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==4","same");
+//  chain->Draw("PhiLab>>h5","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==5","same");
+//  chain->Draw("PhiLab>>h7","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==7","same");
+//  TH2F* histDiag = (TH2F*) gDirectory->Get("histDiag");
+  
+
+//  TH1F* h1 = (TH1F*) gDirectory->Get("h1");
+//  h1->SetTitle("Phi of MUGAST telescopes");
+//  h1->GetXaxis()->SetTitle("PhiLab [deg]");
+//  h1->GetYaxis()->SetTitle("Counts");
+//  h1->SetLineColor(kRed);
+  //h1->SetFillStyle(3001);
+  //h1->SetFillColor(kRed);
+//  TH1F* h2 = (TH1F*) gDirectory->Get("h2");
+//  h2->SetLineColor(kOrange);
+  //h2->SetFillStyle(3001);
+  //h2->SetFillColor(kOrange);
+//  TH1F* h3 = (TH1F*) gDirectory->Get("h3");
+//  h3->SetLineColor(kGreen);
+  //h3->SetFillStyle(3001);
+  //h3->SetFillColor(kGreen);
+//  TH1F* h4 = (TH1F*) gDirectory->Get("h4");
+//  h4->SetLineColor(kTeal);
+  //h4->SetFillStyle(3001);
+  //h4->SetFillColor(kTeal);
+//  TH1F* h5 = (TH1F*) gDirectory->Get("h5");
+//  h5->SetLineColor(kBlue);
+  //h5->SetFillStyle(3001);
+  //h5->SetFillColor(kBlue);
+//  TH1F* h7 = (TH1F*) gDirectory->Get("h7");
+//  h7->SetLineColor(kViolet);
+  //h7->SetFillStyle(3001);
+  //h7->SetFillColor(kViolet);
+//  h7->GetYaxis()->SetRangeUser(0.,4000.);
+//  h7->GetXaxis()->SetRangeUser(-180.,180.);
+
+  //auto legend = new TLegend(0.1,0.7,0.48,0.9);
+  //legend->AddEntry(h1,"MUGAST 1","f");
+  //legend->AddEntry(h2,"MUGAST 2","f");
+  //legend->AddEntry(h3,"MUGAST 3","f");
+  //legend->AddEntry(h4,"MUGAST 4","f");
+  //legend->AddEntry(h5,"MUGAST 5","f");
+  //legend->AddEntry(h7,"MUGAST 7","f");
+  //legend->Draw();
+
+  /* MUGAST good Phi-Ex histogram*/
+  //TCanvas *diagnose = new TCanvas("diagnose","diagnose",1000,1000);
+  //chain->Draw("Ex:PhiLab>>hist(180,-180,180,80,-2,6)","abs(T_MUGAST_VAMOS-2777)<600 && (Mugast.TelescopeNumber==1 || Mugast.TelescopeNumber==2 || Mugast.TelescopeNumber==3 || Mugast.TelescopeNumber==4 || Mugast.TelescopeNumber==5 || Mugast.TelescopeNumber==7 )","colz");
+
+  /* MUGAST testing theta dependance by detector */
+//  TCanvas *diagnose = new TCanvas("diagnose","diagnose",1000,1000);
+//  chain->Draw("Ex:ThetaLab>>histDiag(36,0,180,200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600","colz");
+//  chain->Draw("Ex:ThetaLab>>histDiag(36,0,180,200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==1","colz");
+//  chain->Draw("Ex:ThetaLab>>histDiag(36,0,180,200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==2","colz");
+//  chain->Draw("Ex:ThetaLab>>histDiag(36,0,180,200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==3","colz");
+//  chain->Draw("Ex:ThetaLab>>histDiag(36,0,180,200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==4","colz");
+//  chain->Draw("Ex:ThetaLab>>histDiag(36,0,180,200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==5","colz");
+//  chain->Draw("Ex:ThetaLab>>histDiag(36,0,180,200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==7","colz");
+//  TH2F* histDiag = (TH2F*) gDirectory->Get("histDiag");
+
+//  TCanvas *cG = new TCanvas("cG","cG",1000,1000);
+//  chain->Draw("Ex>>hcG(200,-3,7)","");
+//  chain->Draw("Ex>>hcG_T1(200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==1","");
+//  chain->Draw("Ex>>hcG_T2(200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==2","same");
+//  chain->Draw("Ex>>hcG_T3(200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==3","same");
+//  chain->Draw("Ex>>hcG_T4(200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==4","same");
+//  chain->Draw("Ex>>hcG_T5(200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==5","same");
+//  chain->Draw("Ex>>hcG_T7(200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==7","same");
+//  TH1F* hcG = (TH1F*) gDirectory->Get("hcG");
+//  hcG->GetXaxis()->SetTitle("E_{x} [MeV]");
+//  hcG->GetYaxis()->SetTitle("Counts / (50 keV)");
+//  hcG->SetLineColor(kBlack);
+//  TH1F* hcG_T1 = (TH1F*) gDirectory->Get("hcG_T1");
+//  hcG_T1->SetLineColor(kRed);
+//  hcG_T1->SetFillStyle(3001);
+//  hcG_T1->SetFillColor(kRed);
+//  TH1F* hcG_T2 = (TH1F*) gDirectory->Get("hcG_T2");
+//  hcG_T2->SetLineColor(kOrange);
+//  hcG_T2->SetFillStyle(3001);
+//  hcG_T2->SetFillColor(kOrange);
+//  TH1F* hcG_T3 = (TH1F*) gDirectory->Get("hcG_T3");
+//  hcG_T3->SetLineColor(kGreen);
+//  hcG_T3->SetFillStyle(3001);
+//  hcG_T3->SetFillColor(kGreen);
+//  TH1F* hcG_T4 = (TH1F*) gDirectory->Get("hcG_T4");
+//  hcG_T4->SetLineColor(kTeal);
+//  hcG_T4->SetFillStyle(3001);
+//  hcG_T4->SetFillColor(kTeal);
+//  TH1F* hcG_T5 = (TH1F*) gDirectory->Get("hcG_T5");
+//  hcG_T5->SetTitle("Misalignment of MUGAST telescopes");
+//  hcG_T5->GetXaxis()->SetTitle("E_{x} [MeV]");
+//  hcG_T5->GetYaxis()->SetTitle("Counts / (50 keV)");
+//  hcG_T5->SetLineColor(kBlue);
+//  hcG_T5->SetFillStyle(3001);
+//  hcG_T5->SetFillColor(kBlue);
+//  TH1F* hcG_T7 = (TH1F*) gDirectory->Get("hcG_T7");
+//  hcG_T7->SetLineColor(kViolet);
+//  hcG_T7->SetFillStyle(3001);
+//  hcG_T7->SetFillColor(kViolet);
+//  hcG_T5->GetYaxis()->SetRangeUser(0.,500.);
+//  hcG_T5->GetXaxis()->SetRangeUser(-3.,7.);
+
+
+
+
+  // MUGAST Displacement Ex histogram //
+  TCanvas *cG = new TCanvas("cG","cG",1000,1000);
+//  chain->Draw("Ex>>hcG(200,-3,7)","");
+  chain->Draw("Ex>>hcG_T1(200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==1","");
+  chain->Draw("Ex>>hcG_T2(200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==2","same");
+  chain->Draw("Ex>>hcG_T3(200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==3","same");
+  chain->Draw("Ex>>hcG_T4(200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==4","same");
+  chain->Draw("Ex>>hcG_T5(200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==5","same");
+  chain->Draw("Ex>>hcG_T7(200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==7","same");
+//  TH1F* hcG = (TH1F*) gDirectory->Get("hcG");
+//  hcG->GetXaxis()->SetTitle("E_{x} [MeV]");
+//  hcG->GetYaxis()->SetTitle("Counts / (50 keV)");
+//  hcG->SetLineColor(kBlack);
+  TH1F* hcG_T1 = (TH1F*) gDirectory->Get("hcG_T1");
+  hcG_T1->SetTitle("Misalignment of MUGAST telescopes");
+  hcG_T1->GetXaxis()->SetTitle("E_{x} [MeV]");
+  hcG_T1->GetYaxis()->SetTitle("Counts / (50 keV)");
+  hcG_T1->SetLineColor(kRed);
+  hcG_T1->SetFillStyle(3001);
+  hcG_T1->SetFillColor(kRed);
+  TH1F* hcG_T2 = (TH1F*) gDirectory->Get("hcG_T2");
+  hcG_T2->SetLineColor(kOrange);
+  hcG_T2->SetFillStyle(3001);
+  hcG_T2->SetFillColor(kOrange);
+  TH1F* hcG_T3 = (TH1F*) gDirectory->Get("hcG_T3");
+  hcG_T3->SetLineColor(kGreen);
+  hcG_T3->SetFillStyle(3001);
+  hcG_T3->SetFillColor(kGreen);
+  TH1F* hcG_T4 = (TH1F*) gDirectory->Get("hcG_T4");
+  hcG_T4->SetLineColor(kTeal);
+  hcG_T4->SetFillStyle(3001);
+  hcG_T4->SetFillColor(kTeal);
+  TH1F* hcG_T5 = (TH1F*) gDirectory->Get("hcG_T5");
+  hcG_T5->SetLineColor(kBlue);
+  hcG_T5->SetFillStyle(3001);
+  hcG_T5->SetFillColor(kBlue);
+  TH1F* hcG_T7 = (TH1F*) gDirectory->Get("hcG_T7");
+  hcG_T7->SetLineColor(kViolet);
+  hcG_T7->SetFillStyle(3001);
+  hcG_T7->SetFillColor(kViolet);
+  hcG_T7->GetYaxis()->SetRangeUser(0.,500.);
+
+  TLine *line = new TLine(0,0,0,500);
+  line->SetLineColor(kBlack);
+  line->Draw();
+
+  auto legend = new TLegend(0.1,0.7,0.48,0.9);
+  legend->AddEntry(hcG_T1,"MUGAST 1","f");
+  legend->AddEntry(hcG_T2,"MUGAST 2","f");
+  legend->AddEntry(hcG_T3,"MUGAST 3","f");
+  legend->AddEntry(hcG_T4,"MUGAST 4","f");
+  legend->AddEntry(hcG_T5,"MUGAST 5","f");
+  legend->AddEntry(hcG_T7,"MUGAST 7","f");
+  legend->Draw();
+
+
+/* ORIGINAL - DO NOT EDIT
+  TCanvas *cG = new TCanvas("cG","cG",1000,1000);
+  chain->Draw("Ex>>hcG(200,-3,7)","");
+  chain->Draw("Ex>>hcG2(200,-3,7)","abs(T_MUGAST_VAMOS-2777)<600","same");
+  TH1F* hcG = (TH1F*) gDirectory->Get("hcG");
+  hcG->GetXaxis()->SetTitle("E_{x} [MeV]");
+  hcG->GetYaxis()->SetTitle("Counts / (50 keV)");
+  hcG->SetLineColor(kBlack);
+  TH1F* hcG2 = (TH1F*) gDirectory->Get("hcG2");
+  hcG2->SetLineColor(kRed);
+*/
+
+  /* ORIGINAL 4-PLOT KINEMATIC SCREEN */
+/*  TCanvas *c0 = new TCanvas("c0", "Kinematics", 1000, 1000);
+  c0->Divide(2,2);
+  c0->cd(1);
+  gPad->SetLogz();
+  chain->Draw("ELab:ThetaLab>>hKine(200,0,180,600,0,12)","abs(T_MUGAST_VAMOS-2777)<600","col");
+  TH2F* hKine = (TH2F*) gDirectory->Get("hKine");
+  hKine->GetXaxis()->SetTitle("#theta_{lab} [deg]");
+  hKine->GetYaxis()->SetTitle("E_{lab} [MeV]");
+  plot_kine(Kdp, 0, kBlack, 2, 9);
+  plot_kine(Kdp, 2, kRed, 2, 9);
+  plot_kine(Kdp, 4, kBlue, 2, 9);
+  plot_kine(Kdd, 0, kGreen+2, 2, 9);
+  plot_kine(Kpp, 0, kYellow, 2, 9);
+  plot_kine(Kdt, 0, kMagenta, 2, 9);
+  plot_kine(K12C12C, 4, kCyan, 2, 9);
+  plot_kine(Tidp, 4, 42, 2, 9);
+  plot_kine(Tidt, 4, 42, 2, 9);
+  plot_kine(Tidd, 4, 42, 2, 9);
+  plot_kine(Ti12C12C, 4, 42, 2, 9);
+  
+  c0->cd(2);
+  chain->Draw("Ex>>hEx(300,-3,7)","abs(T_MUGAST_VAMOS-2777)<600");
+  TH1F* hEx = (TH1F*) gDirectory->Get("hEx");
+  hEx->GetXaxis()->SetTitle("E_{x} [MeV]");
+  hEx->GetYaxis()->SetTitle("Counts / (100 keV)");
+  hEx->SetLineColor(kBlack);
+  chain->Draw("Ex>>hEx_gateNoG(300,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && AddBack_EDC@.size()==0","same");
+  chain->Draw("Ex>>hEx_gate(300,-3,7)","abs(T_MUGAST_VAMOS-2777)<600 && abs(AddBack_EDC-0.143)< 0.01","same");
+  TH1F* hEx_gate = (TH1F*) gDirectory->Get("hEx_gate");
+  hEx_gate->SetLineColor(kRed);
+  auto AGATA_eff = new TF1("AGATA_eff","TMath::Exp([0]+[1]*TMath::Log(x)+[2]*pow(TMath::Log(x),2.0)+[3]*pow(TMath::Log(x),3.0)+[4]*pow(TMath::Log(x),4.0))",100,5000);
+  AGATA_eff->SetParameter(0,-7.84071e+00);
+  AGATA_eff->SetParameter(1, 6.44921e+00);
+  AGATA_eff->SetParameter(2, -1.42899e+00);
+  AGATA_eff->SetParameter(3, 1.37921e-01);
+  AGATA_eff->SetParameter(4, -5.23947e-03);
+  
+  //TH1F* hEx_gate_scaled = (TH1F*) hEx_gate->Clone("hEx_gate_scaled");
+  hEx_gate->Scale(1./(AGATA_eff->Eval(143)*0.01));
+  hEx_gate->Draw("histsame");
+  c0->Update();
+  double ymax = gPad->GetUymax();
+  plot_state(0, ymax, kBlack, 2, 9);
+  plot_state(2, ymax, kRed, 2, 9);
+  plot_state(3.8, ymax, kBlue, 2, 9);
+  plot_state(4.644, ymax, kGreen+2, 2, 1);
+  TLatex latex;
+  latex.SetTextAlign(13);
+  latex.SetTextSize(0.035);
+  latex.SetTextAngle(90);
+  latex.SetTextColor(kGreen+2);
+  latex.DrawLatex(4.8,0.6*ymax,"S_{n} = 4.64 MeV");
+
+  TLegend *legend = new TLegend(0.7,0.7,0.9,0.9);
+
+  c0->cd(3);
+  chain->Draw("AddBack_EDC:Ex>>hEgEx(300,-1,7,5000,0,5)","abs(T_MUGAST_VAMOS-2777)<600","col");
+*/
+/*  chain->Draw("AddBack_EDC:Ex>>hEgEx(150,0,6,1000,0,5)","abs(T_MUGAST_VAMOS-2777)<600 && Mugast.TelescopeNumber==5","col");
+  TH2F* hEgEx = (TH2F*) gDirectory->Get("hEgEx");
+  hEgEx->SetTitle("MUGAST#5 only");
+  hEgEx->GetXaxis()->SetTitle("E_{x} [MeV]");
+  hEgEx->GetYaxis()->SetTitle("E_{#gamma} [MeV]");
+  TLine *line = new TLine(0,0,6,6);
+  line->SetLineWidth(2);
+  line->SetLineColor(kRed);
+  line->Draw();
+*/
+/*
+  //  plot_state(0.143, ymax, kYellow, 2, 9);
+  
+  c0->cd(4);
+  chain->Draw("AddBack_EDC>>hEg(4000,0,4)","abs(T_MUGAST_VAMOS-2777)<600");
+  TH1F* hEg = (TH1F*) gDirectory->Get("hEg");
+  hEg->GetXaxis()->SetTitle("E_{#gamma} [MeV]");
+  hEg->GetYaxis()->SetTitle("Counts / (1 keV)");
+  hEg->SetLineColor(kBlack);
+  c0->Update();
+  ymax = gPad->GetUymax();
+  //plot_state(0.143, ymax, kMagenta, 2, 9);
+  //plot_state(0.279, ymax, kCyan, 2, 9);
+  //plot_state(1.863, ymax, kOrange, 2, 9);
+*/
+
+
+/*
+  TCanvas *gammaspec = new TCanvas("gammaSpec", "gammaSpec", 1000, 1000);
+  chain->Draw("AddBack_EDC>>hEg(4000,0,4)","abs(T_MUGAST_VAMOS-2777)<600");
+  TH1F* hEg = (TH1F*) gDirectory->Get("hEg");
+  hEg->GetXaxis()->SetTitle("E_{#gamma} [MeV]");
+  hEg->GetYaxis()->SetTitle("Counts / (1 keV)");
+  hEg->SetLineColor(kBlack);
+  //c0->Update();
+  //ymax = gPad->GetUymax();
+  //plot_state(0.143, ymax, kMagenta, 2, 9);
+  //plot_state(0.279, ymax, kCyan, 2, 9);
+  //plot_state(1.863, ymax, kOrange, 2, 9);
+*/
+
+
+
+
+
+  /*
+  TCanvas *c1 = new TCanvas("c1", "EDC with gates on Ex", 700, 1000);
+  c1->Divide(2,2);
+  c1->cd(1);
+  chain->Draw("AddBack_EDC>>hEg_0(4000,0,4)","abs(T_MUGAST_VAMOS-2777)<600 && abs(Ex)<0.5");
+  TH1F* hEg_0 = (TH1F*) gDirectory->Get("hEg_0");
+  hEg_0->GetXaxis()->SetTitle("E_{#gamma} [MeV]");
+  hEg_0->GetYaxis()->SetTitle("Counts / (1 keV)");
+  hEg_0->SetLineColor(kBlack);
+
+  c1->cd(1);
+  chain->Draw("AddBack_EDC>>hEg_1(4000,0,4)","abs(T_MUGAST_VAMOS-2777)<600 && abs(Ex-1)<0.5");
+  TH1F* hEg_1 = (TH1F*) gDirectory->Get("hEg_1");
+  hEg_1->GetXaxis()->SetTitle("E_{#gamma} [MeV]");
+  hEg_1->GetYaxis()->SetTitle("Counts / (1 keV)");
+  hEg_1->SetLineColor(kBlack);
+
+
+  c1->cd(2);
+  chain->Draw("AddBack_EDC>>hEg_2(4000,0,4)","abs(T_MUGAST_VAMOS-2777)<600 && abs(Ex-2)<0.5");
+  TH1F* hEg_2 = (TH1F*) gDirectory->Get("hEg_2");
+  hEg_2->GetXaxis()->SetTitle("E_{#gamma} [MeV]");
+  hEg_2->GetYaxis()->SetTitle("Counts / (1 keV)");
+  hEg_2->SetLineColor(kBlack);
+
+  c1->cd(3);
+  chain->Draw("AddBack_EDC>>hEg_3(4000,0,4)","abs(T_MUGAST_VAMOS-2777)<600 && abs(Ex-4)<0.5");
+  TH1F* hEg_3 = (TH1F*) gDirectory->Get("hEg_3");
+  hEg_3->GetXaxis()->SetTitle("E_{#gamma} [MeV]");
+  hEg_3->GetYaxis()->SetTitle("Counts / (1 keV)");
+  hEg_3->SetLineColor(kBlack);
+  */
+
+  // K states
+  /*auto gr = K.GetKinematicLine3();
+  gr->SetLineColor(kAzure+7);
+  gr->SetLineWidth(3);
+  gr->Draw("ac");
+  K.SetExcitationHeavy(4);
+  gr = K.GetKinematicLine3();
+  gr->SetLineColor(kAzure+7);
+  gr->SetLineWidth(2);
+  gr->SetLineStyle(1);
+  gr->Draw("c");
+
+  AddTiStates(0); 
+  AddTiStates(0.969); 
+  AddTiStates(2.2292); 
+  AddTiStates(2.419); 
+  AddTiStates(3.223); 
+  AddTiStates(3.332); 
+  AddTiStates(3.622); 
+  AddTiStates(4.388); 
+  AddTiStates(4.458); 
+  AddTiStates(4.719); 
+  AddTiStates(4.852); 
+  AddTiStates(5.151); 
+  */   
+}
diff --git a/Projects/e793s/macro/DrawSpectraRawCal.C b/Projects/e793s/macro/DrawSpectraRawCal.C
new file mode 100644
index 0000000000000000000000000000000000000000..04546d411f0bf14999dee847851bc8a60be4e2b9
--- /dev/null
+++ b/Projects/e793s/macro/DrawSpectraRawCal.C
@@ -0,0 +1,268 @@
+const int MUGAST_TELESCOPES = 6; /* Mugast: 1, 2, 3, 4, 5, 7 */
+const int MUST2_TELESCOPES = 5; 
+
+void DrawTimeRawSpectra(const char* input_file){
+
+  TFile f(input_file);
+  TDirectory *dir = (TDirectory*)(f.Get("ControlSpectra"));
+  dir->cd();
+  TDirectoryFile *dir_mugast = (TDirectoryFile*)(dir->Get("Mugast"));
+  dir_mugast->cd();
+  TH2F *h_mugast[MUGAST_TELESCOPES*2];
+  h_mugast[0] = (TH2F*) dir_mugast->Get("MG1_STRX_T_RAW");
+  h_mugast[1] = (TH2F*) dir_mugast->Get("MG2_STRX_T_RAW");
+  h_mugast[2] = (TH2F*) dir_mugast->Get("MG3_STRX_T_RAW");
+  h_mugast[3] = (TH2F*) dir_mugast->Get("MG4_STRX_T_RAW");
+  h_mugast[4] = (TH2F*) dir_mugast->Get("MG5_STRX_T_RAW");
+//  h_mugast[5] = (TH2F*) dir_mugast->Get("MG6_STRX_T_RAW");
+  h_mugast[5] = (TH2F*) dir_mugast->Get("MG7_STRX_T_RAW");
+//  h_mugast[7] = (TH2F*) dir_mugast->Get("MG9_STRX_T_RAW");
+//  h_mugast[8] = (TH2F*) dir_mugast->Get("MG11_STRX_T_RAW");
+  h_mugast[6] = (TH2F*) dir_mugast->Get("MG1_STRY_T_RAW");
+  h_mugast[7] = (TH2F*) dir_mugast->Get("MG2_STRY_T_RAW");
+  h_mugast[8] = (TH2F*) dir_mugast->Get("MG3_STRY_T_RAW");
+  h_mugast[9] = (TH2F*) dir_mugast->Get("MG4_STRY_T_RAW");
+  h_mugast[10] = (TH2F*) dir_mugast->Get("MG5_STRY_T_RAW");
+//  h_mugast[14] = (TH2F*) dir_mugast->Get("MG6_STRY_T_RAW");
+  h_mugast[11] = (TH2F*) dir_mugast->Get("MG7_STRY_T_RAW");
+//  h_mugast[16] = (TH2F*) dir_mugast->Get("MG9_STRY_T_RAW");
+//  h_mugast[17] = (TH2F*) dir_mugast->Get("MG11_STRY_T_RAW");
+
+  TCanvas *c0 = new TCanvas("c0", "Raw time spectra of Mugast (X e Y)", 1200, 900);
+  c0->Divide(4,3);
+  for(int i=0; i<MUGAST_TELESCOPES*2; i++){
+    h_mugast[i]->SetDirectory(0);
+    c0->cd(i+1); 
+    h_mugast[i]->Draw("colz");
+  }
+  dir_mugast->Close();
+  dir->cd();
+
+  TDirectoryFile *dir_must2 = (TDirectoryFile*)(dir->Get("Must2"));
+  dir_must2->cd();
+  TH2F *h_must[MUST2_TELESCOPES*2];
+  h_must[0] = (TH2F*) dir_must2->Get("MM1_STRX_T_RAW");
+  h_must[1] = (TH2F*) dir_must2->Get("MM2_STRX_T_RAW");
+  h_must[2] = (TH2F*) dir_must2->Get("MM3_STRX_T_RAW");
+  h_must[3] = (TH2F*) dir_must2->Get("MM4_STRX_T_RAW");
+  h_must[4] = (TH2F*) dir_must2->Get("MM5_STRX_T_RAW");
+  h_must[5] = (TH2F*) dir_must2->Get("MM1_STRY_T_RAW");
+  h_must[6] = (TH2F*) dir_must2->Get("MM2_STRY_T_RAW");
+  h_must[7] = (TH2F*) dir_must2->Get("MM3_STRY_T_RAW");
+  h_must[8] = (TH2F*) dir_must2->Get("MM4_STRY_T_RAW");
+  h_must[9] = (TH2F*) dir_must2->Get("MM5_STRY_T_RAW");
+
+  TCanvas *c1 = new TCanvas("c1", "Raw time spectra of MUST2 (X e Y)", 1200, 600);
+  c1->Divide(5,2);
+  for(int i=0; i<MUST2_TELESCOPES*2; i++){
+    h_must[i]->SetDirectory(0);
+    c1->cd(i+1); 
+    h_must[i]->Draw("colz");
+  }
+  dir_must2->Close();
+}
+
+void DrawTimeCalSpectra(const char* input_file){
+
+  TFile f(input_file);
+  TDirectory *dir = (TDirectory*)(f.Get("ControlSpectra"));
+  dir->cd();
+  TDirectoryFile *dir_mugast = (TDirectoryFile*)(dir->Get("Mugast"));
+  dir_mugast->cd();
+  TH2F *h_mugast[MUGAST_TELESCOPES*2];
+  h_mugast[0] = (TH2F*) dir_mugast->Get("MG1_STRX_T_CAL");
+  h_mugast[1] = (TH2F*) dir_mugast->Get("MG2_STRX_T_CAL");
+  h_mugast[2] = (TH2F*) dir_mugast->Get("MG3_STRX_T_CAL");
+  h_mugast[3] = (TH2F*) dir_mugast->Get("MG4_STRX_T_CAL");
+  h_mugast[4] = (TH2F*) dir_mugast->Get("MG5_STRX_T_CAL");
+//  h_mugast[5] = (TH2F*) dir_mugast->Get("MG6_STRX_T_CAL");
+  h_mugast[5] = (TH2F*) dir_mugast->Get("MG7_STRX_T_CAL");
+//  h_mugast[7] = (TH2F*) dir_mugast->Get("MG9_STRX_T_CAL");
+//  h_mugast[8] = (TH2F*) dir_mugast->Get("MG11_STRX_T_CAL");
+  h_mugast[6] = (TH2F*) dir_mugast->Get("MG1_STRY_T_CAL");
+  h_mugast[7] = (TH2F*) dir_mugast->Get("MG2_STRY_T_CAL");
+  h_mugast[8] = (TH2F*) dir_mugast->Get("MG3_STRY_T_CAL");
+  h_mugast[9] = (TH2F*) dir_mugast->Get("MG4_STRY_T_CAL");
+  h_mugast[10] = (TH2F*) dir_mugast->Get("MG5_STRY_T_CAL");
+//  h_mugast[14] = (TH2F*) dir_mugast->Get("MG6_STRY_T_CAL");
+  h_mugast[11] = (TH2F*) dir_mugast->Get("MG7_STRY_T_CAL");
+//  h_mugast[16] = (TH2F*) dir_mugast->Get("MG9_STRY_T_CAL");
+//  h_mugast[17] = (TH2F*) dir_mugast->Get("MG11_STRY_T_CAL");
+
+  TCanvas *c0 = new TCanvas("c0", "Calibrated time spectra of Mugast (X e Y)", 1200, 900);
+  c0->Divide(4,3);
+  for(int i=0; i<MUGAST_TELESCOPES*2; i++){
+    h_mugast[i]->SetDirectory(0);
+    c0->cd(i+1); 
+    h_mugast[i]->Draw("colz");
+  }
+  dir_mugast->Close();
+  dir->cd();
+
+  TDirectoryFile *dir_must2 = (TDirectoryFile*)(dir->Get("Must2"));
+  dir_must2->cd();
+  TH2F *h_must[MUST2_TELESCOPES*2];
+  h_must[0] = (TH2F*) dir_must2->Get("MM1_STRX_T_CAL");
+  h_must[1] = (TH2F*) dir_must2->Get("MM2_STRX_T_CAL");
+  h_must[2] = (TH2F*) dir_must2->Get("MM3_STRX_T_CAL");
+  h_must[3] = (TH2F*) dir_must2->Get("MM4_STRX_T_CAL");
+  h_must[4] = (TH2F*) dir_must2->Get("MM5_STRX_T_CAL");
+  h_must[5] = (TH2F*) dir_must2->Get("MM1_STRY_T_CAL");
+  h_must[6] = (TH2F*) dir_must2->Get("MM2_STRY_T_CAL");
+  h_must[7] = (TH2F*) dir_must2->Get("MM3_STRY_T_CAL");
+  h_must[8] = (TH2F*) dir_must2->Get("MM4_STRY_T_CAL");
+  h_must[9] = (TH2F*) dir_must2->Get("MM5_STRY_T_CAL");
+
+  TCanvas *c1 = new TCanvas("c1", "Calibrated time spectra of MUST2 (X e Y)", 1200, 600);
+  c1->Divide(5,2);
+  for(int i=0; i<MUST2_TELESCOPES*2; i++){
+    h_must[i]->SetDirectory(0);
+    c1->cd(i+1); 
+    h_must[i]->Draw("colz");
+  }
+  dir_must2->Close();
+}
+
+void DrawEnergyRawSpectra(const char* input_file){
+
+  TFile f(input_file);
+  TDirectory *dir = (TDirectory*)(f.Get("ControlSpectra"));
+  dir->cd();
+  TDirectoryFile *dir_mugast = (TDirectoryFile*)(dir->Get("Mugast"));
+  dir_mugast->cd();
+  int binning = 200;
+  TH2F *h_mugast[MUGAST_TELESCOPES*2];
+  h_mugast[0] = (TH2F*) dir_mugast->Get("MG1_STRX_E_RAW");
+  h_mugast[1] = (TH2F*) dir_mugast->Get("MG2_STRX_E_RAW");
+  h_mugast[2] = (TH2F*) dir_mugast->Get("MG3_STRX_E_RAW");
+  h_mugast[3] = (TH2F*) dir_mugast->Get("MG4_STRX_E_RAW");
+  h_mugast[4] = (TH2F*) dir_mugast->Get("MG5_STRX_E_RAW");
+//  h_mugast[5] = (TH2F*) dir_mugast->Get("MG6_STRX_E_RAW");
+  h_mugast[5] = (TH2F*) dir_mugast->Get("MG7_STRX_E_RAW");
+//  h_mugast[7] = (TH2F*) dir_mugast->Get("MG9_STRX_E_RAW");
+//  h_mugast[8] = (TH2F*) dir_mugast->Get("MG11_STRX_E_RAW");
+  h_mugast[6] = (TH2F*) dir_mugast->Get("MG1_STRY_E_RAW");
+  h_mugast[7] = (TH2F*) dir_mugast->Get("MG2_STRY_E_RAW");
+  h_mugast[8] = (TH2F*) dir_mugast->Get("MG3_STRY_E_RAW");
+  h_mugast[9] = (TH2F*) dir_mugast->Get("MG4_STRY_E_RAW");
+  h_mugast[10] = (TH2F*) dir_mugast->Get("MG5_STRY_E_RAW");
+//  h_mugast[14] = (TH2F*) dir_mugast->Get("MG6_STRY_E_RAW");
+  h_mugast[11] = (TH2F*) dir_mugast->Get("MG7_STRY_E_RAW");
+//  h_mugast[16] = (TH2F*) dir_mugast->Get("MG9_STRY_E_RAW");
+//  h_mugast[17] = (TH2F*) dir_mugast->Get("MG11_STRY_E_RAW");
+
+  TCanvas *c0 = new TCanvas("c0", "Raw energy spectra of Mugast (X e Y)", 1200, 600);
+  c0->Divide(4,3);
+  for(int i=0; i<MUGAST_TELESCOPES*2; i++){
+    h_mugast[i]->SetDirectory(0);
+    c0->cd(i+1);
+    if(i<MUGAST_TELESCOPES)
+      h_mugast[i]->GetYaxis()->SetRangeUser(8200,9000);
+    else
+      h_mugast[i]->GetYaxis()->SetRangeUser(7200,8200);
+    
+    h_mugast[i]->Draw("colz");
+    gPad->SetLogz();
+  }
+  dir_mugast->Close();
+  dir->cd();
+
+  TDirectoryFile *dir_must2 = (TDirectoryFile*)(dir->Get("Must2"));
+  dir_must2->cd();
+  binning = 10;
+  TH2F *h_must[MUST2_TELESCOPES*2];
+  h_must[0] = (TH2F*) dir_must2->Get("MM1_STRX_E_RAW");
+  h_must[1] = (TH2F*) dir_must2->Get("MM2_STRX_E_RAW");
+  h_must[2] = (TH2F*) dir_must2->Get("MM3_STRX_E_RAW");
+  h_must[3] = (TH2F*) dir_must2->Get("MM4_STRX_E_RAW");
+  h_must[4] = (TH2F*) dir_must2->Get("MM5_STRX_E_RAW");
+  h_must[5] = (TH2F*) dir_must2->Get("MM1_STRY_E_RAW");
+  h_must[6] = (TH2F*) dir_must2->Get("MM2_STRY_E_RAW");
+  h_must[7] = (TH2F*) dir_must2->Get("MM3_STRY_E_RAW");
+  h_must[8] = (TH2F*) dir_must2->Get("MM4_STRY_E_RAW");
+  h_must[9] = (TH2F*) dir_must2->Get("MM5_STRY_E_RAW");
+
+  TCanvas *c1 = new TCanvas("c1", "Raw energy spectra of MUST2 (X e Y)", 1200, 600);
+  c1->Divide(5,2);
+  for(int i=0; i<MUST2_TELESCOPES*2; i++){
+    h_must[i]->SetDirectory(0);
+    c1->cd(i+1); 
+    if(i<MUST2_TELESCOPES)
+      h_must[i]->GetYaxis()->SetRangeUser(8200,9000);
+    else
+      h_must[i]->GetYaxis()->SetRangeUser(7200,8200);
+ 
+    
+    h_must[i]->Draw("colz");
+    gPad->SetLogz();
+  }
+  dir_must2->Close();
+}
+
+void DrawEnergyCalSpectra(const char* input_file){
+
+  TFile f(input_file);
+  TDirectory *dir = (TDirectory*)(f.Get("ControlSpectra"));
+  dir->cd();
+  TDirectoryFile *dir_mugast = (TDirectoryFile*)(dir->Get("Mugast"));
+  dir_mugast->cd();
+  int binning = 200;
+  TH2F *h_mugast[MUGAST_TELESCOPES*2];
+
+  h_mugast[0] = (TH2F*) dir_mugast->Get("MG1_STRX_E_CAL");
+  h_mugast[1] = (TH2F*) dir_mugast->Get("MG2_STRX_E_CAL");
+  h_mugast[2] = (TH2F*) dir_mugast->Get("MG3_STRX_E_CAL");
+  h_mugast[3] = (TH2F*) dir_mugast->Get("MG4_STRX_E_CAL");
+  h_mugast[4] = (TH2F*) dir_mugast->Get("MG5_STRX_E_CAL");
+//  h_mugast[5] = (TH2F*) dir_mugast->Get("MG6_STRX_E_CAL");
+  h_mugast[5] = (TH2F*) dir_mugast->Get("MG7_STRX_E_CAL");
+//  h_mugast[7] = (TH2F*) dir_mugast->Get("MG9_STRX_E_CAL");
+//  h_mugast[8] = (TH2F*) dir_mugast->Get("MG11_STRX_E_CAL");
+  h_mugast[6] = (TH2F*) dir_mugast->Get("MG1_STRY_E_CAL");
+  h_mugast[7] = (TH2F*) dir_mugast->Get("MG2_STRY_E_CAL");
+  h_mugast[8] = (TH2F*) dir_mugast->Get("MG3_STRY_E_CAL");
+  h_mugast[9] = (TH2F*) dir_mugast->Get("MG4_STRY_E_CAL");
+  h_mugast[10] = (TH2F*) dir_mugast->Get("MG5_STRY_E_CAL");
+//  h_mugast[14] = (TH2F*) dir_mugast->Get("MG6_STRY_E_CAL");
+  h_mugast[11] = (TH2F*) dir_mugast->Get("MG7_STRY_E_CAL");
+//  h_mugast[16] = (TH2F*) dir_mugast->Get("MG9_STRY_E_CAL");
+//  h_mugast[17] = (TH2F*) dir_mugast->Get("MG11_STRY_E_CAL");
+
+  TCanvas *c0 = new TCanvas("c0", "Calibrated energy spectra of Mugast (X e Y)", 1200, 900);
+  c0->Divide(4,3);
+  for(int i=0; i<MUGAST_TELESCOPES*2; i++){
+    h_mugast[i]->SetDirectory(0);
+    h_mugast[i]->GetYaxis()->SetRange(4*binning, 7*binning);
+    c0->cd(i+1); 
+    h_mugast[i]->Draw("colz");
+    gPad->SetLogz();
+  }
+  dir_mugast->Close();
+  dir->cd();
+
+  TDirectoryFile *dir_must2 = (TDirectoryFile*)(dir->Get("Must2"));
+  dir_must2->cd();
+  binning = 10;
+  TH2F *h_must[MUST2_TELESCOPES*2];
+  h_must[0] = (TH2F*) dir_must2->Get("MM1_STRX_E_CAL");
+  h_must[1] = (TH2F*) dir_must2->Get("MM2_STRX_E_CAL");
+  h_must[2] = (TH2F*) dir_must2->Get("MM3_STRX_E_CAL");
+  h_must[3] = (TH2F*) dir_must2->Get("MM4_STRX_E_CAL");
+  h_must[4] = (TH2F*) dir_must2->Get("MM5_STRX_E_CAL");
+  h_must[5] = (TH2F*) dir_must2->Get("MM1_STRY_E_CAL");
+  h_must[6] = (TH2F*) dir_must2->Get("MM2_STRY_E_CAL");
+  h_must[7] = (TH2F*) dir_must2->Get("MM3_STRY_E_CAL");
+  h_must[8] = (TH2F*) dir_must2->Get("MM4_STRY_E_CAL");
+  h_must[9] = (TH2F*) dir_must2->Get("MM5_STRY_E_CAL");
+
+  TCanvas *c1 = new TCanvas("c1", "Calibrated energy spectra of MUST2 (X e Y)", 1200, 600);
+  c1->Divide(5,2);
+  for(int i=0; i<MUST2_TELESCOPES*2; i++){
+    h_must[i]->SetDirectory(0);
+    h_must[i]->GetYaxis()->SetRange(4*binning, 7*binning);
+    c1->cd(i+1); 
+    h_must[i]->Draw("colz");
+    gPad->SetLogz();
+  }
+  dir_must2->Close();
+}