From 8c2e4d128a23d3013f87c1bd9c9633c7ffef576d Mon Sep 17 00:00:00 2001
From: Adrien Matta <matta@lpccaen.in2p3.fr>
Date: Mon, 23 Aug 2021 11:07:22 +0200
Subject: [PATCH] * small changes in MUGAST analysis

---
 NPLib/Detectors/MDM/CMakeLists.txt   |  4 ++--
 NPLib/Detectors/Minos/CMakeLists.txt |  2 +-
 NPSimulation/Process/BeamReaction.cc |  2 ++
 Projects/MUGAST/Analysis.cxx         | 15 +++++++++++++--
 Projects/MUGAST/Analysis.h           |  2 ++
 5 files changed, 20 insertions(+), 5 deletions(-)

diff --git a/NPLib/Detectors/MDM/CMakeLists.txt b/NPLib/Detectors/MDM/CMakeLists.txt
index 99fc5615c..a3ad822b2 100644
--- a/NPLib/Detectors/MDM/CMakeLists.txt
+++ b/NPLib/Detectors/MDM/CMakeLists.txt
@@ -23,11 +23,11 @@ endif()
     add_definitions(-DUSE_RAYTRACE)
     enable_language(Fortran)
     set (CMAKE_Fortran_FLAGS "-O3 -finit-local-zero -falign-commons -fno-automatic -w")
-    add_library(NPMDM SHARED TMDMSpectra.cxx TMDMData.cxx TMDMPhysics.cxx TMDMPhysicsMinimizer TMDMDataDict.cxx TMDMPhysicsDict.cxx TMDMPhysicsMinimizerDict.cxx MDMTrace.cpp RAYTKIN1.F)
+    add_library(NPMDM SHARED TMDMSpectra.cxx TMDMData.cxx TMDMPhysics.cxx TMDMPhysicsMinimizer.cxx TMDMDataDict.cxx TMDMPhysicsDict.cxx TMDMPhysicsMinimizerDict.cxx MDMTrace.cpp RAYTKIN1.F)
  else()
 ## No fortran support, compile "fake" c-version
    message(STATUS "No Fortran support, disabling RAYTRACE in libMDM")
-   add_library(NPMDM SHARED TMDMSpectra.cxx TMDMData.cxx TMDMPhysics.cxx TMDMPhysicsMinimizer TMDMDataDict.cxx TMDMPhysicsDict.cxx TMDMPhysicsMinimizerDict.cxx MDMTrace.cpp)
+   add_library(NPMDM SHARED TMDMSpectra.cxx TMDMData.cxx TMDMPhysics.cxx TMDMPhysicsMinimizer.cxx TMDMDataDict.cxx TMDMPhysicsDict.cxx TMDMPhysicsMinimizerDict.cxx MDMTrace.cpp)
  endif()
 #else()
 #    message(STATUS "Fortran support only included with Unix Makefile generator, disabling RAYTRACE in libMDM.")
diff --git a/NPLib/Detectors/Minos/CMakeLists.txt b/NPLib/Detectors/Minos/CMakeLists.txt
index 60932dbce..d81fe157d 100644
--- a/NPLib/Detectors/Minos/CMakeLists.txt
+++ b/NPLib/Detectors/Minos/CMakeLists.txt
@@ -21,7 +21,7 @@ endif()
 
 if(libMinuit_FOUND)
 ## Link to Minuit library
-  add_library(NPMinos SHARED TMinosClust.cxx TMinosData.cxx TMinosPhysics.cxx TMinosResult oldTMinosPhysicsDict.cxx oldTMinosDataDict.cxx TMinosDataDict.cxx TMinosPhysicsDict.cxx TMinosClustDict.cxx TMinosResultDict.cxx MinosUtility.cxx) 
+  add_library(NPMinos SHARED TMinosClust.cxx TMinosData.cxx TMinosPhysics.cxx TMinosResult.cxx oldTMinosPhysicsDict.cxx oldTMinosDataDict.cxx TMinosDataDict.cxx TMinosPhysicsDict.cxx TMinosClustDict.cxx TMinosResultDict.cxx MinosUtility.cxx) 
   target_link_libraries(NPMinos ${ROOT_LIBRARIES} Minuit NPCore NPTrackReconstruction )
   install(FILES TMinosClust.h TMinosData.h TMinosPhysics.h TMinosResult.h MinosUtility.h DESTINATION ${CMAKE_INCLUDE_OUTPUT_DIRECTORY})
 else()
diff --git a/NPSimulation/Process/BeamReaction.cc b/NPSimulation/Process/BeamReaction.cc
index 2dbb46a31..0744ba1c0 100644
--- a/NPSimulation/Process/BeamReaction.cc
+++ b/NPSimulation/Process/BeamReaction.cc
@@ -37,6 +37,8 @@
 #include <Randomize.hh>
 #include <iostream>
 #include <string>
+#include "G4UserLimits.hh"
+
 ////////////////////////////////////////////////////////////////////////////////
 NPS::BeamReaction::BeamReaction(G4String modelName, G4Region* envelope)
   : G4VFastSimulationModel(modelName, envelope) {
diff --git a/Projects/MUGAST/Analysis.cxx b/Projects/MUGAST/Analysis.cxx
index 81e30e8c5..e5a58a006 100644
--- a/Projects/MUGAST/Analysis.cxx
+++ b/Projects/MUGAST/Analysis.cxx
@@ -135,7 +135,7 @@ void Analysis::TreatEvent() {
   BeamImpact = TVector3(XTarget,YTarget,0); 
   // determine beam energy for a randomized interaction point in target
   // 1% FWHM randominastion (E/100)/2.35
-  //reaction.SetBeamEnergy(Rand.Gaus(myInit->GetIncidentFinalKineticEnergy(),myInit->GetIncidentFinalKineticEnergy()/235));
+  //reaction.SetBeamEnergy(Rand.Gaus(ReactionConditions->GetIncidentFinalKineticEnergy(),ReactionConditions->GetIncidentFinalKineticEnergy()/235));
 
   ////////////////////////////////////////////////////////////////////////////
   //////////////////////////////// LOOP on MUST2  ////////////////////////////
@@ -191,6 +191,10 @@ void Analysis::TreatEvent() {
     /************************************************/
     // Part 3 : Excitation Energy Calculation
     Ex = reaction.ReconstructRelativistic( ELab , ThetaLab );
+    reaction.SetBeamEnergy(Initial->GetIncidentFinalKineticEnergy());
+    ExNoBeam=reaction.ReconstructRelativistic( ELab , ThetaLab );
+    reaction.SetBeamEnergy(FinalBeamEnergy);
+    ExNoProton = reaction.ReconstructRelativistic( ReactionConditions->GetKineticEnergy(0) , ReactionConditions->GetParticleDirection(0).Angle(TVector3(0,0,1)) );
     ThetaLab=ThetaLab/deg;
 
     /************************************************/
@@ -233,7 +237,11 @@ void Analysis::TreatEvent() {
 
     // Part 3 : Excitation Energy Calculation
     Ex = reaction.ReconstructRelativistic( ELab , ThetaLab );
-
+    reaction.SetBeamEnergy(Initial->GetIncidentFinalKineticEnergy()-0.005*Initial->GetIncidentFinalKineticEnergy());
+    ExNoBeam=reaction.ReconstructRelativistic( ELab , ThetaLab );
+    reaction.SetBeamEnergy(FinalBeamEnergy);
+    ExNoProton = reaction.ReconstructRelativistic( ReactionConditions->GetKineticEnergy(0) , ReactionConditions->GetParticleDirection(0).Angle(TVector3(0,0,1)) );
+    
 
     // Part 4 : Theta CM Calculation
     ThetaCM  = reaction.EnergyLabToThetaCM( ELab , ThetaLab)/deg;
@@ -273,6 +281,8 @@ void Analysis::End(){
 ////////////////////////////////////////////////////////////////////////////////
 void Analysis::InitOutputBranch() {
   RootOutput::getInstance()->GetTree()->Branch("Ex",&Ex,"Ex/D");
+  RootOutput::getInstance()->GetTree()->Branch("ExNoBeam",&ExNoBeam,"ExNoBeam/D");
+  RootOutput::getInstance()->GetTree()->Branch("ExNoProton",&ExNoProton,"ExNoProton/D");
   RootOutput::getInstance()->GetTree()->Branch("EDC",&Ex,"Ex/D");
   RootOutput::getInstance()->GetTree()->Branch("ELab",&ELab,"ELab/D");
   RootOutput::getInstance()->GetTree()->Branch("ThetaLab",&ThetaLab,"ThetaLab/D");
@@ -346,6 +356,7 @@ void Analysis::InitInputBranch(){
 ////////////////////////////////////////////////////////////////////////////////
 void Analysis::ReInitValue(){
   Ex = -1000 ;
+  ExNoBeam=ExNoProton=-1000;
   EDC= -1000;
   ELab = -1000;
   BeamEnergy = -1000;
diff --git a/Projects/MUGAST/Analysis.h b/Projects/MUGAST/Analysis.h
index f7ae5c274..a0881f4f3 100644
--- a/Projects/MUGAST/Analysis.h
+++ b/Projects/MUGAST/Analysis.h
@@ -52,6 +52,8 @@ class Analysis: public NPL::VAnalysis{
  
   private:
   double Ex;
+  double ExNoBeam;
+  double ExNoProton;
   double EDC;
   double ELab;
   double ThetaLab;
-- 
GitLab